© ESO 2021

1. Introduction

Pulsars in the Galactic centre of the Milky Way, in particular if found orbiting a stellar mass or intermediate-mass black hole, or the central supermassive black hole (SMBH) Sagittarius A* (hereafter Sgr A*), could allow direct accurate measurements of black hole properties such as mass, spin, or even quadrupole moment, and potentially enable the most stringent tests of general relativity and alternative theories of gravity to date (Wex & Kopeikin 1999; Kramer et al. 2004; Pfahl & Loeb 2004; Liu et al. 2012; Psaltis et al. 2016). In addition, Galactic centre pulsars would enable precision measurements of the interstellar medium (ISM) properties along the line of sight and at the central galactic region, such as electron density, scattering, or magnetic field (e.g., Shannon & Johnston 2013; Eatough et al. 2013a; Spitler et al. 2014; Bower et al. 2014; Dexter et al. 2017). Moreover, discovering a population of pulsars in the nuclear cluster of the Milky Way could help us understand the enigmatic and complex star formation history and evolution in the region (for a review, see e.g., Genzel et al. 2010). Such pulsar-based measurements, in particular those related to Sgr A*, would furthermore strongly complement the ongoing efforts to understand black holes, their environments, and the theory of gravity through the study of orbital dynamics (Gravity Collaboration 2018, 2020; Do et al. 2019) or through event-horizon-scale imaging (Goddi et al. 2017; Event Horizon Telescope Collaboration 2019a,b).

Because of the scientific potential of pulsars located in the Galactic centre, significant efforts have been undertaken to try to find pulsars in the surroundings of Sgr A* (e.g., Johnston et al. 1995, 2006; Kramer et al. 2000; Morris et al. 2002; Klein 2005; Deneva et al. 2009; Macquart et al. 2010; Bates et al. 2011; Wharton et al. 2012a; Eatough et al. 2013b; Siemion et al. 2013). Interestingly, all attempts to find a population of pulsars within the inner ∼20 arcmin of the Milky Way have been unsuccessful so far: only one radio magnetar, PSR J1745−2900, was found located ≈2.4 arcsec (≈0.1 pc in projection) from Sgr A* (Kennea et al. 2013; Mori et al. 2013; Rea et al. 2013; Eatough et al. 2013a; Shannon & Johnston 2013). Unfortunately, PSR J1745−2900 is still too far from the SMBH to enable the desired gravity tests, although it allowed us to measure the gas properties close to Sgr A*, showing that the gas is dense, turbulent, highly magnetised, and varies on relatively short timescales (Eatough et al. 2013a; Desvignes et al. 2018).

It is puzzling that no more pulsars are found in the central part of the Galaxy. Massive stars populate the region with pulsars progenitors (for a review, see e.g., Figer 2009), and supernova remnants and pulsar wind nebulae candidates have been identified (Wang et al. 2006; Muno et al. 2008; Zhang et al. 2020), as well as stellar mass black holes (Hailey et al. 2018; Generozov et al. 2018). The scenario of a large population of Galactic centre pulsars is also supported by the Galactic centre excess of gamma-ray emission, which can be explained by an undetected large population of millisecond pulsars (MSPs; e.g., Abazajian 2011; Gordon & Macias 2014; Brandt & Kocsis 2015). Theoretical estimates of the number of observable pulsars in the inner part of the Galaxy vary considerably between a few to several hundreds (Pfahl & Loeb 2004; Faucher-Giguère & Loeb 2011; Wharton et al. 2012b; Dexter & O’Leary 2014; Chennamangalam & Lorimer 2014; Rajwade et al. 2018; Schödel et al. 2020), showing the complexity of the Galactic centre stellar population, evolution, and its environment.

Different theories trying to explain the lack of detections have been proposed. For example, a deficit of normal pulsars in the Galactic centre with an unusual population composed mainly of magnetars (which rarely emit in radio, and whose emission may be intermittent and therefore difficult to detect; Dexter & O’Leary 2014), or millisecond pulsars (which might tend to be less luminous than normal pulsars and thus beyond the reach of our currently achievable sensitivities; Bailes et al. 1997; Kramer et al. 1998; Burgay et al. 2013; Macquart & Kanekar 2015). Geodetic precession of pulsars due to massive companions could move their emission beams in and out our line of sight, making them undetectable during long periods of time (e.g., Kramer et al. 1998; Perera et al. 2010), and eclipses may temporarily obscure the radio emission in certain binary systems (e.g., Freire 2005). Similarly, for tight binary pulsars, certain epochs might face orbital ranges with characteristics that are difficult to detect by searching algorithms (e.g., Eatough et al. 2021). Even dark matter has been considered, which could produce a rapid collapse of the neutron stars into black holes through accretion, making the life of pulsars unusually short in the inner Galactic centre region (Bramante & Linden 2014).

Although some of the previously outlined reasons are plausible, the scattering of the pulsed signals by the dense turbulent gas in the direction of the Galactic centre is currently widely accepted as the most probable cause for the non-detections (Lazio & Cordes 1998; Cordes & Lazio 2002). Nevertheless, Spitler et al. (2014) and Bower et al. (2014) showed that the scattering in the direction of PSR J1745−2900, only ≈2.4 arcsec from Sgr A*, is much weaker than predicted by theoretical models (see also Wucknitz 2014). A possible explanation for the low scattering measured for PSR J1745−2900 and the paucity of detections of other Galactic centre pulsars is that scattering in the Galactic centre region may be complex and strongly dependent on the line of sight because of multiple screens. This situation has also been found in other directions (e.g., Roy 2013; Schnitzeler et al. 2016; Dexter et al. 2017).

The ISM effects, and in particular, scattering, are strongly inverse dependent on the observing frequency (see e.g., Lorimer & Kramer 2004). The detrimental effect of the scattering that may prevent the detection of the Galactic centre pulsars can therefore be significantly alleviated by increasing the frequency of the radio observations. The drawback of this approach is that pulsars, which usually are steep spectral sources (⟨α⟩= − 1.8 ± 0.2, for S ∝ ν^α, Maron et al. 2000), become extremely weak at high frequencies. However, some pulsars have shown flat spectral indices (α > −1.0, see e.g., the PSRCAT database, Manchester et al. 2005) that keep the radio emission strong even at very high radio frequencies. Interestingly, certain pulsar emission models furthermore predict that incoherent emission may take over in pulsar spectra between radio and infrared wavelengths, potentially making the emission brighter at millimetre wavelengths than the power-law extrapolation from centimetre wavelengths (see Michel 1978; Kramer et al. 1996, 1997). The highest radio frequency at which pulsations from a normal pulsar were detected is 138 GHz (λ = 2.17 mm) (Torne 2017; Torne et al., in prep.), and radio magnetars have been detected up to ∼300 GHz (Torne et al. 2017, 2020b). This supports the hypothesis that we can try to discover pulsars using very high radio frequencies.

Because of the lack of pulsar discoveries and in order to overcome the ISM effects, in particular, the strong scattering, targeted surveys at the Galactic centre have progressively increased the radio frequency of the observations. The highest radio frequencies used to date for Galactic centre surveys were around 20 GHz (Siemion et al. 2013; Eatough et al. 2013b). In this work, we present a new targeted pulsar survey at the Galactic centre that for the first time uses short millimetre wavelengths between 3 and 2 mm (ν ∼ 84–156 GHz). The paper is structured as follows. The observations are introduced in Sect. 2, and the data analysis is explained in detail in Sect. 3. Section 4 presents and discusses our results, and Sect. 5 includes a summary and our conclusions.

2. Observations

The observations were made with the Institut de Radioastronomie Millimétrique (IRAM) 30m telescope on Pico Veleta, Spain, during different campaigns aimed primarily at monitoring the Galactic centre magnetar PSR J1745−2900 between December 2016 and May 2018 (under project numbers 159-16, 039-17, and 145-17). The receiver we used was the Superconductor-Insulator-Superconductor (SIS) Eight MIxer Receiver (EMIR, Carter et al. 2012), tuned to different frequencies in the 3.4, 2.1, 1.3, and 1.0 mm bands (from ∼80 to ∼300 GHz) at different epochs. For this pulsar survey, only the observations from the 3.4 and 2.1 mm bands (central frequencies 86 and 102, and 138 and 154 GHz, hereafter referred to as 3 and 2 mm bands) are analysed. This selection was made because the telescope is more efficient at these bands¹, and pulsars should be brighter here than at even higher frequencies. Additionally, we minimise the effect on signal quality and sensitivity by atmospheric effects, in particular, the low-frequency noise and signal absorption, which are magnified by the low telescope elevations that are necessary to observe the Galactic centre from the IRAM 30m telescope (elevation ≲24 deg).

The backend used to record the signal from EMIR was the Broad-Band Continuum backend (BBC), which provided one time series with a resolution of 100 μs for each frequency band and polarisation. This combination of EMIR and the BBC resulted in an effective non-contiguous bandwidth of ∼32 GHz (four bands of ∼8 GHz centred at the above-mentioned central frequencies of 86, 102, 138, and 154 GHz). With this set-up, each observation produced eight time series: from four frequency bands times two linear polarisations.

The observations at each epoch typically consisted of sessions of three to four hours, in which continuous tracks of PSR J1745−2900 of 10–20 min were alternated with calibration measurements. Tracking a source continuously is an observing strategy that differs from typical observing methods in millimetre astronomy, where the beam is usually switched on and off the target source to remove undesired instrumental and atmospheric effects. The reason for using tracking when observing pulsars is mainly twofold: to keep the coherence of the time series, which is not ensured if the telescope uses beam-switching techniques; and to gain sensitivity because with beam-switching techniques half of the time is spent off source (losing of sensitivity), with an additional of loss when the off is subtracted.

The beam size of the IRAM 30m telescope is ≈29, 24, 18, and 16 arcsec at 86, 102, 138, and 154 GHz, respectively. Sgr A* is 2.4 ± 0.3 arcsec from PSR J1745−2900 (Rea et al. 2013), therefore it is included close to the centre of the beams for all frequencies. In the worst case, that is, for the smallest beam sizes at the 2 mm band, the 2.4 arcsec separation from Sgr A* would reduce the sensitivity by ≲20%, reaching negligible reductions of sensitivity for the 3 mm beams. At a distance to the Galactic centre of d_GC =  8.18 kpc (Gravity Collaboration 2019), the projected diameter around PSR J1745−2900 is 1.16, 0.96, 0.72, and 0.64 pc for each of the four central observing frequencies.

The total number of analysed epochs between December 2016 and May 2018 is 28. The observations amount to 62.2 h covering the Galactic centre. The individual epochs are presented in Table 1.

3. Data analysis

The data analysis was split into four parts: calibration, data cleaning and preparation, pulsar searching, and pre-processing tests with synthetic pulsar signal injection and sensitivity analysis. We describe each of them in the following subsections.

3.1. Calibration

The first step of the analysis consisted of calibrating the data by calculating conversion factors from counts to Jansky and multiplying the time series correspondingly. This was made individually for each frequency band and each polarisation channel. The values required to calculate the receiver temperature and telescope gain were computed from a hot-cold calibration on loads of known physical temperature. The opacity of the atmosphere was calculated through a short measurement on sky 600 arcsec away in right ascension from the science target, using the information from the telescope weather station and the model called atmopsheric transmission at microwaves (ATM; Pardo et al. 2001) to derive the equivalent sky temperature. With this, we calculated the calibration temperature, T_cal, which was used to convert from counts into , that is, antenna temperature outside the atmosphere (see Kramer 1997). The scaling factor from to Jansky scale, , is obtained from the observatory efficiency tables². The factor is not given for all frequencies, therefore a second-order polynomial was fit to those provided, assuming an error of 10 and 15% for the 3 and 2 mm bands, respectively. We obtain , 6.0, 6.3, and 6.5 Jy K⁻¹ for 86, 102, 138, and 154 GHz. Finally, a correction to the gain due to gain dependence on the telescope elevation was applied (Peñalver 2012). The calibration measurements were made every 10–20 min to follow possible gain variations and opacity changes through observations. The final calibration factor applied to each time series was calculated for the middle time of each observation by interpolating linearly between the calibration results just before and after each scan.

3.2. Data cleaning and preparation

We used a receiver designed for spectroscopic observations, for which long-term stability was not a requirement. In our fast-sampled long integrations³, the data therefore occasionally showed short period of instabilities with an increased noise amplitude towards negative values and narrow drop-offs of the signal intensity. In addition, the data showed an excess of power weighted towards the low-frequency part of the Fourier spectrum (hereafter referred to as red noise). Red noise mainly arises from the atmospheric opacity variations during the observations. After the data were calibrated, we therefore applied a filtering to subtract or reduce the instrumental effects and the red noise.

In a first step, the intensity drop-offs were clipped by a moving a window of 3.5 s over the time series, substituting the samples that were below −5σ with the median value of the samples within the window. The drop-offs occurred only a few times per scan, and their filtering does not affect our sensitivity. This running-clip filtering also partly reduced the short-term instabilities that show increased negative amplitudes in the noise when they exceed −5σ. Next, a Fourier transform was carried out, and a number of frequency bins were zapped (i.e. the real and imaginary parts were set to the average value and zero, respectively) to filter a strong 1 Hz signal caused by the cryo generator cycle and the 50 Hz power mains. Five and nine harmonics of the two signals were also zapped. We then returned to the time domain by applying an inverse Fourier transform and applied a filtering to reduce the red noise, which consists of a moving window of 3 s that subtracts a fitted polynomial of first order. This resulted in time series that were flat and did almost not affect the potential pulsations, whether inside the noise or even potential single pulses. However, we note that this filtering scheme may cause very slowly rotating pulsars to appear below the noise mean level after folding. Our tests showed that only periods above P ≳ 4 s may be significantly affected, which reduces our sensitivity to long-period pulsars both in the periodicity search and when folding to produce candidate plots. Finally, we also note that the noise from the EMIR receiver shows locally generated periodic signals of strong intensity that proved difficult to model and filter because their amplitude and frequency structures vary slightly in time and can in addition vary with each new receiver tuning. We used an adaptive filtering in the Fourier domain to partly subtract these signals and discuss its potential effect in our pulsar detection capabilities in Sect. 4.

After the calibration and the cleaning, the two polarisation channels of each frequency band were summed, producing the total intensity time series per frequency band. At this stage, we had four time series per scan (one for each frequency band). From these total intensity time series, three different datasets were created: One set combines the two frequency bands at 3 mm (86+102 GHz), one set combines the two frequency bands at 2 mm (138+154 GHz), and a third set combines the four frequency bands (86+102+138+154 GHz), which we name 3+2 mm. The combination was made by averaging together the total intensity time series of the different frequency bands. We remark that the time series of the different frequencies are aligned in time, that is, the start time for all of them is exactly the same, thus enabling this simple averaging. We discuss the reasons for not focusing on the theoretically most sensitive 3+2 mm dataset alone, which includes the effect of uncorrected interstellar dispersion and technical motivations, in Sect. 4.

Lastly, as the observations were split into scans of 10–20 min, we coherently concatenated together all the scans for each dataset per epoch. The concatenation was made by filling the gaps with the median value (which was zero after the red noise filtering applied before). The total length of the concatenated time series varied for different epochs, with typical values around 3–4 h.

Figure 1 shows an example of data before and after the calibration, data cleaning, and preparation prior to searching. Table 1 summarises the observations and presents the data length before and after concatenation of scans, and the average system temperature per epoch and frequency band.

Fig. 1. Visualisation of an example time series of the data from EMIR plus BBC before (upper panels) and after (bottom panels) the processing presented in Sects. 3.1 and 3.2. The processing includes calibration, data cleaning, and preparation prior to the pulsar searching. Upper left panel: raw data as recorded by the BBC during an observation consisting of five scans of 20 min each during a total time span of about 2 h. The atmospheric opacity variations and the sporadic intensity drop-offs are apparent. The gaps in between scans show intervals when the telescope was not observing the Galactic centre, e.g., in the course of calibration scans, or pointing and focus adjustments. Upper right panel: histogram with the distribution of the count levels through the observation, showing that it is not ideal to use the data directly in this state for statistics-based algorithms such as the pulsar searching. Bottom left and bottom right panels: same time series after the processing. The slow level variations and the intensity drop-offs are corrected for. Furthermore, the intensity scale at this stage is calibrated in Jansky, which is necessary in order to correctly combine the different polarisations and frequency bands. Bottom left panel: slight increase in the standard deviation of the noise toward the end of the observation. This is because the Galactic centre decreases in elevation on sky, which forces the telescope to observe through a larger amount of atmosphere (or airmass), thus increasing the system temperature. The data shown correspond to the 102 GHz frequency band from 2017 September 28. This epoch was chosen to illustrate a case with a high opacity and thus a high system temperature, which typically is more challenging to clean.

3.3. Pulsar searching

We performed the search with a pipeline based on the PRESTO software version 2.1⁴ (Ransom 2001). We searched the three datasets (3, 2 and 3+2 mm) of each epoch individually. The pipeline consists of the following searching algorithm.

In a first step, we analyse the data for interfering signals using rfifind, which tries to identify chunks of bad data by comparing global and block-size statistics. rfifind produces a file with a mask that can later be used to flag specific blocks of data identified as having poor statistics. This automatic analysis did not always work; it often masked 100% of the data. We attribute this large masking to the somewhat peculiar statistics of the noise from EMIR, which can deviate from well-behaved Gaussian noise, especially when integrated without off-source-position subtraction (see Sect. 2). Whenever the masked fraction was above 30%, we passed a filtered version of the data to rfifind that we created using the rednoise filter of PRESTO on each scan before concatenating them. The rednoise filter acts in the Fourier domain by subtracting the slope of the power series using adaptive window sizes. A normalisation is also applied by rednoise so that the standard deviation of the powers in the filtered series equals one. This filtering proved very effective in reducing the undesired periodic instrumental signals in the data and improving the data statistics (see Fig. 2). However, rednoise can also partially filter potential pulsar signals when their spin frequency is very low or very high. For this reason, we did not use it to remove the red noise during data preparation, and opted for the running-fit filter (see Sect. 3.2). The program rfifind worked better in this rednoise-filtered version of the data, although the analysis of certain epochs still resulted in unacceptably high masking. If the masking fraction was still above 30% after the rednoise filter was applied, we did not allow rfifind to automatically compute a bad-data mask. The typical percentage of the automatically masked fraction of data after this routine is about 1−3%. After these steps, we manually passed the intervals of the gaps between scans that where filled with the median value to rfifind so that they were excluded during data analysis.

Fig. 2. Example of the ability and importance of the PRESTO rednoise filter to reduce the EMIR locally generated signals that interfere with the search algorithm. The plots show a blind detection of the Galactic centre magnetar PSR J1745−2900 in one epoch in which the pulsar was relatively weak (2017 May 18, mJy). Left panel: candidate plot from prepfold when the raw data are folded and rednoise is not applied. Right panel: same candidate when the data are cleaned by rednoise prior to folding. The filter does subtract part of the pulsar signal, as apparent by the narrow profile and the dip in the profile below the mean noise level, but the improvement compared to a non-filtered version is so substantial that the pulsar may easily have been missed if the filtering had not been applied. The main panels in each plot show (a) the folded profile twice, (b) the candidate and observation details including name, telescope, epoch, time sampling, statistics, source coordinates, period and period derivatives (in topocentric and barycentric references), and binary parameters if applicable, (c) waterfall plot of folded signal intensity vs. time, (d) accumulated reduced χ2 of the integrated profile vs. a model of noise alone, and (e) two-dimensional waterfall plot of the reduced χ2 of the candidate profile vs. noise as a function of folding period and period derivative. The white gaps in the phase vs. time panel are intervals in which the telescope did not observe the source, e.g., during calibration, pointing, or focus, or for some epochs, when observations were taking place at different frequency bands that were not used for the pulsar search (see Sect. 2).

A barycentred time series was then produced with prepsubband using the previously created rfifind mask to substitute bad data blocks with median values, and without dispersion correction. The time series was then Fourier transformed, the rednoise filter applied to the Fourier series, and a number of persistent periodic instrumental signals from a database created in pre-processing tests were removed by zapping their corresponding frequency bins. The resulting Fourier series was searched for pulsar candidates with accelsearch, a routine that detects significant periodicities in the Fourier domain, including harmonic summing and a template-matched algorithm to recover signals with a Doppler shift (as occurs in pulsars that are accelerated when they orbit companions).

Two passes of accelsearch were made. The first pass set the PRESTO zmax parameter to zero. This zmax parameter is the maximum width in Fourier bins of the templates that are used to recover accelerated signals in the Fourier domain (Ransom et al. 2002). Thus, the first pass is sensitive to isolated or very lowly accelerated pulsars. In order to be sensitive to potential binary systems, including pulsars orbiting stellar and intermediate-mass black holes and the central SMBH, the zmax parameter of accelsearch was set to 1200 in a second pass. This is the maximum currently allowed by the software. This high zmax parameter is motivated by the requirement that as many harmonics as possible need to be recovered to increase the sensitivity, in combination with the long total integration times, in particular after concatenating the scans per epoch (see Table 1). We remark that the drift in the Fourier domain (ṙ) due to the frequency shift (ḟ) by the Doppler effect of an accelerated pulsar has a quadratic dependence on integration time (t_int), ṙ = ḟ(t_int)² (Ransom 2001). In both passes we allowed for summing up to 32 harmonics, and set the sigma threshold to 2.0 and 3.0 for the zmax = 0 and 1200, respectively. The candidates from the two passes were then sifted to harmonically remove related candidates, duplicates, and candidates with sifting significance < 2σ. In addition, we only required candidates to be detected with one harmonic to pass the sifting filter. Although the sigma thresholds may look low, pre-processing tests conducted to find optimum parameters for the pipeline found that setting the significance limits this low enabled detections of weak pulsar signals in the data (see Sect. 3.5) while still producing a manageable number of total candidates. The number of candidates produced in each epoch is shown in Table 1. Typically, there were well below 200 candidates per observation.

The last step in the pipeline folds the data and produces plots to visually analyse the final candidates. This was done with prepfold. We carried out the folding twice for each candidate. At first, we folded the concatenated filterbank that is passed to the pipeline. Then we folded a concatenated filterbank that is filtered by the rednoise filter with the same parameters. The reason for folding twice is that the rednoise-filtered data tend to show a significantly cleaner profile and candidate signal (i.e. it is effective in removing parts of the local periodic interfering signals, as mentioned earlier), but the filtering may also subtract part of the candidate signal (see Fig. 2). Thus, by producing and reviewing the two folded candidate plots (original and filtered), we avoid missing a potential real pulsar by a filtering that is too strong or insufficient.

Finally, we selected a number of observations made under good weather conditions, at different epochs, and searched them with a modified version of the pipeline that included searching in jerk space (i.e. derivative in acceleration). Everything was analogous in these cases except for the parameters for the acceleration search of accelsearch, which were zmax = 300 and wmax = 900. The parameter wmax controls the maximum size of the matched templates in the dimension (Andersen & Ransom 2018). A visual summary of the data analysis is shown in Fig. 3.

Fig. 3. Flowchart of the data analysis as described in Sect. 3. The algorithm is applied to each observing epoch individually. The steps inside the dashed block are applied in parallel to each single-polarisation time series of each frequency band. When this block finalises, the algorithm continues to reduce the data, creating total intensity time series and averaging frequency bands to create the three datasets 3, 2, and 3+2 mm (see Sect. 3.2). These datasets then enter the searching algorithm as summarised in the right column of the diagram.

3.4. Sensitivity analysis

To estimate the survey sensitivity, we used the radiometer equation adapted for pulsar observations to calculate a minimum detectable flux density (see e.g., Lorimer & Kramer 2004),

(1)

where (S/N)_min is the required signal-to-noise ratio for a detection, T_sys is the telescope system temperature, G is the telescope gain⁵, n_p is the number of polarisations, t_int is the integration time, Δν is the instantaneous observing bandwidth, and W_eff and P are the effective width of the pulse and spin period of the pulsar, respectively. Taking into account the effects of scattering, dispersion, and instrumentation on the pulse widths, W_eff is calculated as the sum in quadrature of the intrinsic width of the pulse profile W_int, the characteristic scattering time τ_s, the dispersion smearing across one frequency channel Δt_DM, and the time resolution of the digitised data δt (e.g., Manchester et al. 1996),

(2)

Given the relatively varied dataset presented in this work, with observations taken at different frequency bands and weather conditions (and therefore with different T_sys) and of different total integration time, we selected representative values of observations of good quality for our sensitivity analysis, that is, those obtained when the integration time was long and the weather was good. We focus on the scenario with the averaged 3+2 mm dataset (i.e. a central observing frequency of 120 GHz with a bandwidth of 32 GHz) as it is considered the most sensitive due to the large bandwidth.

Inside Eq. (1), two parameters depend on the particular pulsar observed: the spin period P, and the pulse width, W_eff. These two parameters can be expressed combined in terms of duty cycle, δ = W_eff/P. Typical duty cycles are about 5 and 20% for normal and millisecond pulsars, respectively (see Kramer et al. 1998, and references therein). We assumed a main scenario with δ = 0.1, that is, a duty cycle of 10%. For completeness and to illustrate how the duty cycle and uncorrected smearing by dispersion affect the results, we discuss scenarios with δ = 0.05 and δ = 0.4 and include the effect of the uncorrected smearing in Sect. 4. For comparison and to place the sensitivity of the IRAM 30m telescope into context, we considered a potential and realistic survey with the Atacama Large Millimeter Array (ALMA)⁶. As ALMA has a better visibility of the Galactic centre, we assumed a total on-source integration time of 5 h, with the maximum bandwidth set to 8 GHz, limited by the current capability of the recording backends (Remijan et al. 2020). Table 2 summarises the parameters we used to calculate the minimum detectable flux densities, resulting in mJy and mJy for the IRAM 30m and for ALMA, respectively.

The S_min value can be used to directly compare the sensitivity with those of other surveys, but to better estimate the potential of this survey of detecting pulsars in the Galactic centre, we simulated a scenario in which we populated the Galactic centre region that is covered by our telescope beam with pulsars resembling those that are known to exist in the Milky Way. Then we calculated the percentage of these presumed Galactic centre pulsars that emit above the minimum detectable threshold. Those pulsars emitting above the threshold could have been discovered by the survey.

The putative Galactic centre population was derived using the Online Pulsar Catalog PSRCAT⁷ version 1.62 (Manchester et al. 2005). We extracted the pseudo-luminosity of those pulsars for which this information is available⁸. In total, 2125 pulsars out of 2800 (∼76%) have pseudo-luminosity information in our PSRCAT dataset. However, the reported pseudo-luminosity is always for frequencies below those used in this survey; mainly at 400 (∼8.5%) and 1400 MHz (∼67%). To asses any effect from possible spectral gigahertz-peaked turn-overs (see e.g., Kijak et al. 2017), we verified that the results did not vary when only the 1400 MHz pseudo-luminosities were used or the 400 and 1400 MHz information was included. We extrapolated the pseudo-luminosity to 120 GHz (the central frequency of the 3+2 mm dataset) with the following method: If a pulsar had a measured spectral index available in PSRCAT, we used this spectral index for the extrapolation. If no spectral index was available for a pulsar, we drew a value from a normal distribution with mean value and standard deviation calculated from the 788 pulsars (i.e. ∼28% of the total) in PSRCAT with known spectral indices, separating the MSP population⁹ from the rest. The mean spectral index obtained for the MSP population is ⟨α_MSP⟩= − 1.95 ± 0.54, and for the rest, it is ⟨α⟩= − 1.76 ± 0.79 (see Fig. 4, cf. Maron et al. 2000; Jankowski et al. 2018). In both cases we assumed a single power law for the extrapolation, that is, S_ν ∝ ν^α. The population simulation was repeated 5000 times to derive an error figure from the standard deviation. We note that the mean spectral indices derived for the two populations are statistically equivalent because they deviate by less than 1σ. However, our simulations of pulsar luminosities at millimetre wavelengths show that for similar average values, the standard deviation of the mean spectral index plays a fundamental role in the number of pulsars that may posses a sufficiently flat spectrum to be detectable¹⁰. As a consequence, the MSP population, showing a smaller standard deviation in their spectral indices than the rest, tend to appear at a lower luminosity at higher radio frequencies. This property may affect the discovery potential of MSPs in our survey, adding to another effect that complicates the detection of fast MSPs due to uncorrected interstellar dispersion, as explained in Sect. 4.

Fig. 4. Histograms of spectral index distribution for the populations of MSPs (P  ≤  30 ms) and normal population (P > 30 ms) from the database PSRCAT v1.62 (Manchester et al. 2005). A normal distribution is fit in each case, resulting in ⟨αMSP⟩= − 1.95 ± 0.54 for the MSP-only population, and ⟨α⟩= − 1.76 ± 0.79 for the rest.

Finally, we calculated the detection limit of the survey by converting our minimum detectable flux density into a pseudo-luminosity limit at the Galactic centre multiplying S_min by the square of the distance to the Galactic centre (see e.g., Lorimer & Kramer 2004), d_GC = 8.18 kpc (Gravity Collaboration 2019),

(3)

With the luminosity limit at the Galactic centre, we computed the number of detectable pulsars by the survey from the simulated population as those where the pseudo-luminosity at the observing frequency was higher than the limit (see Sect. 4).

3.5. Mock pulsar injection

The theoretical limit presented in Sect. 3.4 assumes Gaussian noise from the instruments and does not take the potential effect of undesired interfering signals, red noise, of the filtering applied to reduce them, or the limitation of the searching algorithms to recover highly accelerated pulsars into account. To try to account for these factors and evaluate the sensitivity thresholds, we carried out an independent preliminary analysis for which we injected six synthetic pulsar signals into the real data for a large number of epochs and performed the searching process. We manually adjusted the flux density, S, of the injected mock signals until they reached the detection limit, defined as the flux density of the synthetic signals for which the pulsars are detected by the pipeline in good weather and long integration epochs, but missed in those epochs when the conditions are not optimal. To produce the mock pulsar signals, we used two codes: for the isolated pulsars, the PRESTO injectpsr.py routine, and for the pulsars orbiting companions, a custom version of fake from the software package SIGPROC¹¹. The pulsed-averaged flux density of the synthetic pulsars was calculated by folding a noiseless time series with the pulsar signal and measuring S as the integrated profile divided by the spin period. We took the masked fractions in the searched data into account, which mainly account for the gaps in between scans (see Sect. 3.2 and Table 1). The masks have different sizes for different observations, and we used the median value of the fraction of filled gaps over the total time after concatenation (32%) to calculate representative values of S. The properties of the injected pulsars and the empirically derived flux density values are summarised in Table 3.

The mock data injection tests also allowed us to validate our data cleaning steps and to evaluate the limits of the algorithms in recovering highly accelerated signals in our data and optimise the significance thresholds both for accelsearch and the candidate-sifting step. In addition, these pre-processing tests, in particular when simulating tight orbits around black holes, served as a confirmation of the ability of the pipeline to recover the mock pulsars, and therefore similar real pulsar would they exist in the data¹². Finally, we established the need of a very high value for the parameter zmax of the PRESTO accelsearch for the most extreme systems (see Sect. 4, also Eatough et al. 2021).

4. Results and discussion

The search of the presented dataset produced a total of 5431 pulsar candidates, yielding 10862 plots that were individually reviewed¹³. The most prominent candidate repeatedly is the Galactic centre magnetar PSR J1745−2900, which is detected at a high significance in most epochs. In addition to PSR J1745−2900, a number of potential pulsar candidates were identified, but most showed a relatively low significance, and all appeared only once. Without a high-significance and consistently repeated detection of a given candidate periodicity, we considered that no new pulsars were distinctly discovered in this survey. Because a detection of any of these candidates at lower frequencies, where they could be brighter and where information of the dispersion measure can be obtained, could work as a confirmation of a potential real new pulsar, we add the properties of our 20 best candidates in Table 4 and show the 4 best candidate plots in Fig. 5. These candidates were chosen from a combination of detection significance and their pulsar-like characteristics in the candidate plots after folding.

Fig. 5. Four best candidates from the IRAM 30m pulsar survey as evaluated from a combination of detection significance and the characteristics of the candidate after folding. The main properties on which our selection is based are persistence of the signal in time, significance in the period-period derivative panel, σ value, and the profile shape. A list with the properties of these four and the other selected best candidates is given in Table 4. Upper left panel: candidate with a close period to one of the synthetic injected pulsars, but we remark that the candidate arises from data with no mock pulsar injected. Bottom right panel: the candidate has a relatively close period to PSR J1745−2900, but it seems to be harmonically unrelated. A description of the different panels of the plots is presented in the caption of Fig. 2.

When we assume that pulsars exist in the Galactic centre, two main reasons explain the lack of new discoveries. On one hand, the decreased flux density of pulsars at short millimetre wavelengths due to their typically steep spectral index. On the other hand, the current sensitivity of the IRAM 30m telescope, which even though it is very high for a facility operating in this frequency regime, is relatively low when compared to 100 m class telescopes and may be currently insufficient to detect the weak pulsar signals. Our analysis of the survey sensitivity following Sect. 3 results in a representative minimum detectable flux density of Jy. The same limit expressed in terms of pseudo-luminosity at the distance of the Galactic centre (following Eq. (3)) is . For a putative Galactic centre pulsar population resembling that of the Milky Way, only about 3.5 ± 0.3% of the pulsars could be detected with this limit. The error figure of 0.3% shows the 1σ statistical standard deviation related to the population simulation (see Sect. 3.4). However, when mock pulsar signals are used to evaluate a sensitivity limit, the result is even lower. The analysis shows that our detection threshold can sometimes be higher by up to a factor ∼10 than the theoretical one (see the S column of Table 3 and the position of the synthetic pulsars with respect to the theoretical limit in Fig. 6). With these higher thresholds, our survey would be sensitive to only about 2% of isolated pulsars and about only 1% of those in highly accelerated binary systems. We furthermore remark that our coverage to MSPs is close to zero. The latter fact arises from the simulation showing that under the assumptions on spectral index we made, the number of MSPs possessing a flat spectrum that therefore are bright enough to be detected at our observing frequencies is much lower than for the remaining pulsar population.

Fig. 6. Sensitivity to pulsars in the Galactic centre by this survey (continuous red line) and a hypothetical survey with ALMA at similar frequencies (dashed light blue line). The black dots represent the pseudo-luminosities for the known population of pulsars extrapolated to 120 GHz (see Sect. 3.4 for details), simulating a population of pulsars in the Galactic centre. The dots above each corresponding threshold line are pulsars that may theoretically be detected by the surveys. The sensitivity lines are calculated assuming a duty cycle for the minimum detectable flux density of δ = 0.1. We assumed a hyper-strong scattering towards the Galactic centre following the NE2001 model, s (Cordes & Lazio 2002) to consider a worst-case scenario. The effect of scattering is negligible (see text). The dominant interstellar effect is the uncorrected dispersion smearing for the IRAM 30m data, which is apparent at the very end of the short spin periods as a raising limit. The dispersion smearing effect is not affecting the ALMA limit because we consider that the ALMA data will enable the correction for dispersion. The red stars indicate the synthetic pulsar signals injected for the sensitivity tests. The three circles around the red stars with P ≃ 1.92, 55.22, and 141.87 ms mark the three pulsars orbiting black holes that were injected to test the sensitivity to highly accelerated pulsar systems. The orange diamond shows the average pseudo-luminosity of the Galactic centre magnetar PSR J1745−2900, and the blue triangle and green square mark the 120 GHz pseudo-luminosities of the Vela pulsar (PSR B0833−45) and PSR B0355+54, respectively.

We relate the decreased sensitivity shown by the mock pulsar analysis to the characteristics of the noise from the EMIR receiver, and for the binary systems, to an additional component from the limitation of the algorithms in recovering pulsar signals of highly accelerated systems. For example, for the synthetic binary in a 30 h orbit of a 10 M_⊙ black hole, even with a zmax = 1200, accelsearch can only recover eight harmonics in our datasets. For the fast-spinning MSP (P = 1.92 ms) in a six-month orbit of Sgr A* only one harmonic is recovered. We note, nevertheless, that the limitations of accelsearch occur only for extremely accelerated systems and are largely augmented by our long total integration times after concatenating the scans filling the gaps in between them (see Sect. 3). For a more detailed discussion of the limits for accelerated pulsar recovery with PRESTO, including effects such as jerk or eccentricity, see Eatough et al. (2021), Liu et al. (2021).

Furthermore, the effect of the duty cycle of the potential pulsars to be detected is non-negligible (see Eq. (1)). For instance, a shorter duty cycle of δ = 0.05 would increase the theoretical percentage of average detectable population to ∼4.4%, and a duty cycle δ = 0.4 would lower it to ∼1.9%. For comparison, our simulated survey with ALMA (for δ = 0.1) would theoretically cover about 9.6% of the population (11.1% if the derived sensitivity from QA2 is used), showing one path towards higher-sensitivity millimetre-wavelength surveys in the future. Figure 6 illustrates the sensitivity thresholds over the simulated Galactic centre population, including our synthetic injected pulsars, and three real pulsars for reference.

We remark that we may have an additional component of sensitivity loss to the MSP population. Firstly, MSPs in binary systems are more challenging to detect with the applied searching algorithms because the drifts in frequency space due to Doppler effects are larger (and therefore are more difficult to correct, see Eatough et al. 2021) than for slow pulsars. Secondly, our inability to correct for interstellar dispersion with a continuum backend could translate into a lower sensitivity to fast-spinning MSPs (with periods P ≲ 3 ms) due to pulse smearing. This last issue is one reason why three separated datasets were searched (3, 2, and 3+2 mm, see Sect. 3.2). The dispersion delay for the 3 mm-only data (between 104 to 84 GHz) for a dispersion measure (DM), DM = 2000 pc cm⁻³, is Δt_DM ≈  0.37 ms. For 2 mm-only (156–136 GHz), the dispersion delay would be Δt_DM ≈  0.11 ms. In the case of 3+2 mm (156–84 GHz), a DM = 2000 pc cm⁻³ would produce a frequency-dependent delay of Δt_DM ≈  0.83 ms. In this latest case, the dispersion is appreciably larger than the sampling interval and can therefore smear a sufficiently narrow pulse. We graphically show in Fig. 6 the corresponding potential decreased sensitivity to MSPs by including the effect of uncorrected dispersive smearing in the IRAM 30m data through Eqs. (2) and (1), assuming a DM = 2000 pc cm⁻³. The DM = 2000 pc cm⁻³ was chosen based on the DM measured for the closest pulsar to Sgr A* known today: PSR J1745−2900 with DM = 1778 pc cm⁻³ (Eatough et al. 2013a). Even so, the DM closer to Sgr A* could be even higher for example by the contribution of a denser gas or even an accretion disk (Cordes & Lazio 2002), which would reduce our ability to detect short-period MSPs in this survey even further.

To the previously outlined challenges of detecting fast-spinning MSPs, we add that the probability is high that due to the very high stellar densities of the Galactic centre (and thus a much higher probability of stellar interactions), a good portion of the pulsars populating the region could be MSPs (e.g., Wharton et al. 2012b; Macquart & Kanekar 2015). A scenario in which the Galactic centre pulsar population is dominated by MSPs could therefore explain the lack of discoveries in our survey. Interestingly, the majority of the best candidates yielded by the IRAM 30m survey are MSPs, with P < 10 ms. The limitations discussed above make us cautious to consider these MSP candidates as real, but the results might at the same time be indicative of an MSP-dominant population in the Galactic centre.

In addition to optimising the potential effect of uncorrected interstellar dispersion, there is a second reason for searching the three frequency datasets separately. It relates to the EMIR locally generated periodic signals that are more prominent in the 3 mm band than in the 2 mm band¹⁴. As a result, when we combined the 3 and 2 mm bands, we increased the theoretical sensitivity by increasing the bandwidth, but we added these undesired signals as well. On the other hand, the 2 mm band of EMIR is more stable and produces cleaner noise, so that the Fourier spectrum is also cleaner. Despite the effects of the larger and somewhat less-clean noise at 3 mm, we verified during the mock injected signal analysis and with the blind detections of PSR J1745−2900 that the combination of the 3 and 2 mm bands can increase the detection capability for weak pulsars with relatively flat spectra. In summary, the 3 mm dataset is potentially better suited to detect pulsars with not-so-flat spectral indices. The 2 mm band shows a better-behaved noise, but it is not optimal because the frequency is higher and pulsars tend to be dimmer. The 3+2 mm data are finally the most sensitive with their larger combined bandwidth, but shows the undesired stronger interference from the combined frequency bands. To minimise the potential adverse effects and maximise the chances of detecting pulsars in our full dataset, we searched the three datasets individually.

Some of the locally generated signals are periodic and polluted certain regions of the Fourier series. The rednoise filter of PRESTO  reduced them in part, but their intensity was sometimes so strong that we considered that parts of the Fourier series cannot be used to detect pulsars. Consequently, we might lose pulsars even with luminosities above our detection thresholds. This would nonetheless only occur if a pulsar rotational frequency lies exactly within a bin in the power spectrum that is affected by the local periodic signals. We estimate that ∼3.5% on average of the Fourier bins are affected. Although this is not too significant, it is an additional cause for potential sensitivity loss.

A final challenge to mention is the effect of the red noise. It has been demonstrated that pulsar surveys lose sensitivity to long-period pulsars when red noise is present in the data (Lazarus et al. 2015). Our IRAM 30m data, and generally all pulsar observations at very high radio frequencies where variable opacity from the atmosphere contributes significantly to the system noise, show a considerable amount of red noise. The running-fit filter prior to searching has been developed to be effective at reducing this red noise while at the same time not affecting significantly potential pulsar signals. This has been verified by optimising detections of PSR J1745−2900, with a spin period P ≃ 3.76 s. We therefore successfully reduced the effect of the red noise in our survey, but we cannot eliminate it completely. We thus consider our sensitivity to very long-period pulsars (P ≳ 4 s) reduced further than the limits presented above.

Despite the limitations and challenges discussed earlier, the method we followed here is capable of detecting new pulsars. This is demonstrated by previous successful detections of pulsars with the EMIR receiver at the IRAM 30m telescope that covered the same frequency ranges as we used in this survey (Torne et al. 2015, 2017, 2020a), the injection of simulated pulsar signals correctly detected in our pipeline, and the consistent blind detection of the magnetar PSR J1745−2900 by our searching algorithm. To place our detection thresholds further into context, PSR J1745−2900 belongs to the 0.5% most luminous pulsars at 120 GHz (Torne et al. 2015, 2017), and PSR B0355+54 (one of the brightest pulsars with one of the flattest spectra known) belongs to the 10% most luminous population at 120 GHz (Morris et al. 1997, Torne et al. in prep.). The Vela pulsar (PSR B0833−45) is very bright at 1.4 GHz, but its spectral index of α ≃  −1.5 (taken for the component with the flatter spectrum, Keith et al. 2011) places its luminosity at ∼ 100 GHz (Liu et al. 2019) at a level below even the ALMA threshold. These examples show that we only reach the highest-luminosity population in this survey, and at the same time, that we are highly biased to flat-spectrum pulsars when we search in the short-millimetre regime. It would therefore be reasonable to believe that there may still be pulsars in the Galactic centre that have not yet been detected simply because they are dim at short-millimetre wavelengths and are beyond the reach of our current sensitivity.

The survey presented here has two remarkable advantages over previous surveys carried out at typical pulsar-observing wavelengths in the centimetre regime. In the first place, the scattering is negligible. Equations (2) and (1) show how the temporal scattering would normally affect the sensitivity to detect pulsars, although the effect can be more complex, for instance affecting the periodicity-search algorithms (Kramer et al. 2000; Macquart et al. 2010). Even if we consider a worst-case scenario and assume a hyper-strong scattering to the Galactic centre with a scattering time of 2000 s at 1 GHz (NE2001 model, Cordes & Lazio 2002), the scattering in our lowest frequency range would be τ_s ≈ 2000 ⋅ 86⁻⁴ = 36 μs (see e.g., Lorimer & Kramer 2004). This scattering is negligible even for the fastest-spinning pulsars known. For comparison, in the same situation, a survey at an observing frequency of 2 GHz (λ = 15 cm) would be exposed to a scattering time of τ_s = 125 s, and even at 5 GHz (λ = 6 cm), this would be τ_s = 3.2 s. These long scattering times would disrupt the detection of pulsations at centimetre wavelengths. Even in the case of much lower scattering towards the Galactic centre, the population of fast-spinning MSPs may still be beyond the reach of surveys below ∼5 GHz (Spitler et al. 2014). The second advantage of the survey is the very low dispersion at millimetre wavelengths. This means that this parameter space does not require an intensive search. This situation allows us to save computing time that can be spent surveying the acceleration parameter space. Thus, this first short-millimetre-wavelength survey with the IRAM 30m telescope covers the parameter space of highly scattered and highly accelerated pulsars in the Galactic centre region for the first time (cf. e.g., Johnston et al. 2006; Deneva et al. 2009; Macquart et al. 2010; Bates et al. 2011; Wharton et al. 2012a; Eatough et al. 2013b; Siemion et al. 2013), and it is virtually unbiased to population coverage regardless of the scattering characteristics in this direction.

We note that the small effect of interstellar dispersion that allows us to search using a continuum backend and to expand our search space in acceleration becomes a great disadvantage in terms of discriminating interference or local signals from celestial ones. On one hand, the dispersion correction in a pulsar-searching algorithm at low radio frequencies smears local broadband interference, reducing its effect in the periodicity search, especially at high DM values such as are expected in the Galactic centre. At short millimetre wavelengths such as we used in the presented survey, we greatly limit this ability to suppress local broadband signals, and the local interference enters the search algorithm at almost full strength. This has two main consequences. First, as discussed earlier, the periodic interference may contaminate specific parts of the Fourier series, with the potential to hinder the detection of a pulsar with an equal or close spin period to the period of the interfering signal. On the other hand, the interference will increase the number of candidates produced by the search, sometimes very significantly, with the consequent additional difficulty in the candidate sifting and review steps. Related to this issue, the lack of spectral information of the pulsar candidates in our search (from using a continuum backend) affects our capability to discern broad- and narrow-band signals. It was therefore not possible to identify pulsar-resembling candidates produced by narrow-band locally generated signals.

Finally, it is remarkable that this survey covers a long time span, with repeated observations. This is an advantage compared to single-epoch or few-epoch surveys, and makes the search robust against time-varying factors that can prevent the detection of pulsars. One of the factors is weather (atmospheric opacity and turbulence), which can affect the sensitivity of short-millimetre observations greatly if, when only a few epochs are observed, the weather is not optimal. Other important factors whose effects are diminished by the repeated observations are the potential geodetic precession of pulsars in orbit with massive objects due to curved space-time, which can move the radio beams in and out of our line of sight (e.g., Desvignes et al. 2019), or eclipses that can temporarily hide the pulsar radio emission (e.g., Lyne et al. 1990). Lastly, the repeated observations at different epochs increase the opportunities of observing an extreme binary system in parts of the orbit in which the acceleration effects are less detrimental for the signal recovery in the searching algorithm, for example where the jerk is low (see e.g., Eatough et al. 2021).

When we look into the future and focus on short-millimetre wavelengths, one logical way to improve the surveys is by increasing the sensitivity by using larger telescopes such as the Large Millimeter Telescope (LMT), the NOrthern Extended Millimetre Array (NOEMA)¹⁵, or ALMA. From the instrumental side, other methods for increasing the sensitivity include the extension of the instantaneous bandwidth, the use of different receiver technology such as Kinetic Inductance Detectors (Torne et al. 2020b), and the use of fast spectrometers to enable dedispersion and obtaining the spectral information of the candidates. Other improvements may arise from the searching algorithms, for example by applying the fast-folding algorithm (FFA, see e.g., Cameron et al. 2017) or extending the acceleration and jerk parameter spaces surveyed. In particular at very high radio frequencies, a more advanced data cleaning applying smart filters to reduce the undesired potential interfering signals and the red noise while not affecting potential pulsar signals can also help to increase our sensitivity. Surveying using the frequency range of ν ≃ 30−50 GHz (λ ≃ 10−6 mm) would also be compelling because the scattering will still be largely reduced but pulsars are more luminous on average. To take advantage of long observing campaigns, searching in an incoherent or even a coherent stack of a number of epochs promises to significantly improve the sensitivity (Eatough et al. 2013b; Pan et al. 2016; Lentati et al. 2018). Finally, another procedure to try to find new pulsars and other signals such as those from rotating radio transients (RRATs, see e.g., Keane & McLaughlin 2011) or fast radio bursts (FRBs, see e.g., Cordes & Chatterjee 2019) is through the search for single pulsations or transient emission. We did not apply this technique in our survey because in our experience, an insufficient dispersion effect and frequency information to filter local interfering transient signals prevented an efficient transient-like candidate classification.

In spite of challenges such as the large distance, potential strong scattering, unknown population, and sensitivity limitations, the undeniable importance of finding pulsars in the Galactic centre and in particular, those orbiting Sgr A*, justify additional and continued efforts to survey this region. Surveys in the short-millimetre regime are novel and need to be further improved, but they offer a new tool in our quest to uncover the hidden population of pulsars in the centre of the Milky Way and behind extremely scattered regions.

5. Summary and conclusions

We presented a targeted pulsar survey of the Galactic centre that for the first time used short-millimetre wavelengths. After the analysis of a total of 62.2 hours of observations in 28 different epochs in a period spanning almost 1.5 years, no clear detections of new pulsars were found. The non-detections can be explained because the sensitivity of the observations is limited. It was estimated that only about 2 and 1% of a hypothetical Galactic centre population located close to Sgr A* would be detected for isolated and highly accelerated binary pulsars, respectively. For the particular case of Galactic centre MSPs, the current sensitivity appears insufficient to expect detections at short-millimetre wavelengths under our assumptions for population simulations.

The main reasons for the limited sensitivity to pulsars in the Galactic centre are a combination of the use of very high observing frequencies to overcome the scattering with the typical steep spectral index of pulsars in the radio band, making pulsars extremely dim at millimetre wavelengths, the large distance to the centre of the galaxy, and the relatively low sensitivity of the IRAM 30m telescope mainly due to its dish size (when compared to 100 m class telescopes that are generally used for pulsar science). We also found an additional loss of sensitivity from other factors. Firstly, from a considerable amount of red noise in the data, typical from slow atmospheric opacity variations during the observations. The excess of red noise mainly affects the detection of long-period pulsars. This red noise effect was reduced by our filtering steps prior to the search by using a running-fit filter. Nonetheless, the detrimental effect of red noise cannot be fully eliminated. Second, periodic locally generated signals generated inside the EMIR receiver pollute parts of the Fourier series, which diminishes our detection capability for certain particular periodicities. Although this only affects a small percentage of our Fourier series, we successfully used the PRESTO rednoise filter to significantly reduce these undesired signals, but for certain regions of the Fourier spectrum, we may still lose real pulsars if the locally generated signals have periods equal to the potential pulsar spin period. Finally, we were unable to correct for dispersive smearing by using a continuum backend. This had two main drawbacks: an additional loss in sensitivity to fast-spinning MSPs, and the impossibility of distinguishing locally generated periodic or interfering signals mimicking pulsars from potential real pulsar signals due to the lack of DM measurement. In addition, we were unable to differentiate pulsar candidates produced from narrow- and broad-band signals, making the candidate classification less effective. For this last obstacle, future observations with a very high frequency and temporal resolution, or at a lower observing frequency, may be needed to fully confirm a candidate by the detection of a dispersion measure different from zero and the confirmation of the broad-band nature of the emission as expected for radio pulsars.

In contrast, the principal advantages of searching for pulsars using short-millimetre wavelengths are a reduction of the scattering to negligible levels, being effectively unbiased in population coverage, and the possibility of searching a much larger acceleration and jerk parameter space for a given computing time in exchange for the very small (or negligible) space of dispersion measures to be blindly searched. Despite the sensitivity and the challenging data reduction, this pulsar survey with the IRAM 30m telescope is the first to cover a parameter space that is basically unaffected by temporal scattering and thus has the potential of detecting fast-spinning MSPs in the Galactic centre if the scattering times are long in this direction. The survey results reveal that no fundamental problem prevents us from carrying out pulsar surveys using a standard pulsar-searching software at these very high radio frequencies, and therefore motivates other similar surveys, conceivably with more sensitive instruments.

The interest of searching the IRAM 30m data is proven by the detections in our blind-search pipeline of the magnetar PSR J1745−2900 and a number of mock pulsar signals, including a pulsar-stellar black hole system and two pulsars in close orbit with Sgr A*. Our tests show that the survey was capable of detecting new pulsars if they were sufficiently bright, and the non-detections rule out that other bright radio magnetars similar to PSR J1745−2900 or highly luminous pulsars in the frequency ranges covered and with beams pointing towards us exist in the inner region of the Galactic centre, at least during the time span covered by the survey.

Similar pulsar surveys with improved sensitivity and involving larger telescopes, such as the LMT, NOEMA, and ALMA are encouraged. In the future, the use of new instrumentation with better performance in pulsar observations at millimetre observatories will help us to further constrain the characteristics of the potential pulsar population living inside the inner parsec of our galaxy.

Acknowledgments

The authors thank the anonymous referee for the quick review of the manuscript. PT wishes to thank Scott Ransom for his help to clarify details on the PRESTO routines, Ciriaco Goddi and Sergio Dzib for useful comments on the draft, and the staff of the IRAM 30 m telescope for their great support during the observations. This work is based on observations carried out with the IRAM 30 m telescope. IRAM is supported by INSU/CNRS (France), MPG (Germany) and IGN (Spain). Financial support by the European Research Council for the ERC SynergyGrant BlackHoleCam (ERC-2013-SyG, GrantAgreement no. 610058) is gratefully acknowledged. RE is supported by a “FAST Fellowship” under the “Cultivation Project for FAST Scientific Payoff and Research Achievement of the Center for Astronomical Mega-Science, Chinese Academy of Sciences (CAMS-CAS)”. This research has made use of NASA’s Astrophysics Data System.

References

All Tables

All Figures

Fig. 1. Visualisation of an example time series of the data from EMIR plus BBC before (upper panels) and after (bottom panels) the processing presented in Sects. 3.1 and 3.2. The processing includes calibration, data cleaning, and preparation prior to the pulsar searching. Upper left panel: raw data as recorded by the BBC during an observation consisting of five scans of 20 min each during a total time span of about 2 h. The atmospheric opacity variations and the sporadic intensity drop-offs are apparent. The gaps in between scans show intervals when the telescope was not observing the Galactic centre, e.g., in the course of calibration scans, or pointing and focus adjustments. Upper right panel: histogram with the distribution of the count levels through the observation, showing that it is not ideal to use the data directly in this state for statistics-based algorithms such as the pulsar searching. Bottom left and bottom right panels: same time series after the processing. The slow level variations and the intensity drop-offs are corrected for. Furthermore, the intensity scale at this stage is calibrated in Jansky, which is necessary in order to correctly combine the different polarisations and frequency bands. Bottom left panel: slight increase in the standard deviation of the noise toward the end of the observation. This is because the Galactic centre decreases in elevation on sky, which forces the telescope to observe through a larger amount of atmosphere (or airmass), thus increasing the system temperature. The data shown correspond to the 102 GHz frequency band from 2017 September 28. This epoch was chosen to illustrate a case with a high opacity and thus a high system temperature, which typically is more challenging to clean.

In the text

Fig. 2. Example of the ability and importance of the PRESTO rednoise filter to reduce the EMIR locally generated signals that interfere with the search algorithm. The plots show a blind detection of the Galactic centre magnetar PSR J1745−2900 in one epoch in which the pulsar was relatively weak (2017 May 18, mJy). Left panel: candidate plot from prepfold when the raw data are folded and rednoise is not applied. Right panel: same candidate when the data are cleaned by rednoise prior to folding. The filter does subtract part of the pulsar signal, as apparent by the narrow profile and the dip in the profile below the mean noise level, but the improvement compared to a non-filtered version is so substantial that the pulsar may easily have been missed if the filtering had not been applied. The main panels in each plot show (a) the folded profile twice, (b) the candidate and observation details including name, telescope, epoch, time sampling, statistics, source coordinates, period and period derivatives (in topocentric and barycentric references), and binary parameters if applicable, (c) waterfall plot of folded signal intensity vs. time, (d) accumulated reduced χ2 of the integrated profile vs. a model of noise alone, and (e) two-dimensional waterfall plot of the reduced χ2 of the candidate profile vs. noise as a function of folding period and period derivative. The white gaps in the phase vs. time panel are intervals in which the telescope did not observe the source, e.g., during calibration, pointing, or focus, or for some epochs, when observations were taking place at different frequency bands that were not used for the pulsar search (see Sect. 2).

In the text

Fig. 3. Flowchart of the data analysis as described in Sect. 3. The algorithm is applied to each observing epoch individually. The steps inside the dashed block are applied in parallel to each single-polarisation time series of each frequency band. When this block finalises, the algorithm continues to reduce the data, creating total intensity time series and averaging frequency bands to create the three datasets 3, 2, and 3+2 mm (see Sect. 3.2). These datasets then enter the searching algorithm as summarised in the right column of the diagram.

In the text

	Fig. 4. Histograms of spectral index distribution for the populations of MSPs (P  ≤  30 ms) and normal population (P > 30 ms) from the database PSRCAT v1.62 (Manchester et al. 2005). A normal distribution is fit in each case, resulting in ⟨αMSP⟩= − 1.95 ± 0.54 for the MSP-only population, and ⟨α⟩= − 1.76 ± 0.79 for the rest.
In the text

Fig. 5. Four best candidates from the IRAM 30m pulsar survey as evaluated from a combination of detection significance and the characteristics of the candidate after folding. The main properties on which our selection is based are persistence of the signal in time, significance in the period-period derivative panel, σ value, and the profile shape. A list with the properties of these four and the other selected best candidates is given in Table 4. Upper left panel: candidate with a close period to one of the synthetic injected pulsars, but we remark that the candidate arises from data with no mock pulsar injected. Bottom right panel: the candidate has a relatively close period to PSR J1745−2900, but it seems to be harmonically unrelated. A description of the different panels of the plots is presented in the caption of Fig. 2.

In the text

Fig. 6. Sensitivity to pulsars in the Galactic centre by this survey (continuous red line) and a hypothetical survey with ALMA at similar frequencies (dashed light blue line). The black dots represent the pseudo-luminosities for the known population of pulsars extrapolated to 120 GHz (see Sect. 3.4 for details), simulating a population of pulsars in the Galactic centre. The dots above each corresponding threshold line are pulsars that may theoretically be detected by the surveys. The sensitivity lines are calculated assuming a duty cycle for the minimum detectable flux density of δ = 0.1. We assumed a hyper-strong scattering towards the Galactic centre following the NE2001 model, s (Cordes & Lazio 2002) to consider a worst-case scenario. The effect of scattering is negligible (see text). The dominant interstellar effect is the uncorrected dispersion smearing for the IRAM 30m data, which is apparent at the very end of the short spin periods as a raising limit. The dispersion smearing effect is not affecting the ALMA limit because we consider that the ALMA data will enable the correction for dispersion. The red stars indicate the synthetic pulsar signals injected for the sensitivity tests. The three circles around the red stars with P ≃ 1.92, 55.22, and 141.87 ms mark the three pulsars orbiting black holes that were injected to test the sensitivity to highly accelerated pulsar systems. The orange diamond shows the average pseudo-luminosity of the Galactic centre magnetar PSR J1745−2900, and the blue triangle and green square mark the 120 GHz pseudo-luminosities of the Vela pulsar (PSR B0833−45) and PSR B0355+54, respectively.

In the text