A publishing partnership

The following article is Open access

Search for Continuous Gravitational Waves from Known Pulsars in the First Part of the Fourth LIGO-Virgo-KAGRA Observing Run

, , , , , , , , ,

Published 2025 April 10 • © 2025. The Author(s). Published by the American Astronomical Society.
, , Citation A. G. Abac et al 2025 ApJ 983 99DOI 10.3847/1538-4357/adb3a0

Download Article PDF
DownloadArticle ePub

You need an eReader or compatible software to experience the benefits of the ePub3 file format.

Article metrics

403 Total downloads
0 Video abstract views

Share this article

0004-637X/983/2/99

Abstract

Continuous gravitational waves (CWs) emission from neutron stars carries information about their internal structure and equation of state, and it can provide tests of general relativity. We present a search for CWs from a set of 45 known pulsars in the first part of the fourth LIGO–Virgo–KAGRA observing run, known as O4a. We conducted a targeted search for each pulsar using three independent analysis methods considering single-harmonic and dual-harmonic emission models. We find no evidence of a CW signal in O4a data for both models and set upper limits on the signal amplitude and on the ellipticity, which quantifies the asymmetry in the neutron star mass distribution. For the single-harmonic emission model, 29 targets have the upper limit on the amplitude below the theoretical spin-down limit. The lowest upper limit on the amplitude is 6.4 × 10−27 for the young energetic pulsar J0537−6910, while the lowest constraint on the ellipticity is 8.8 × 10−9 for the bright nearby millisecond pulsar J0437−4715. Additionally, for a subset of 16 targets, we performed a narrowband search that is more robust regarding the emission model, with no evidence of a signal. We also found no evidence of nonstandard polarizations as predicted by the Brans–Dicke theory.

Export citation and abstractBibTeXRIS

Previous article in issue
Next article in issue

Original content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.

1. Introduction

Since their discovery in 1967, pulsars have been crucial in advancing our understanding of fundamental physics. These extremely dense and compact objects possess strong magnetic fields and rotate rapidly, emitting beams of electromagnetic (EM) radiation from hot spots near the poles or from the higher magnetosphere (A. Philippov & M. Kramer 2022). EM observations across various wavelengths (radio, X-rays, and gamma rays) have provided detailed insights into pulsar properties, allowing precise measurements of pulsar parameters. Given their stability and predictability, pulsars present an excellent opportunity for the search of continuous gravitational waves (CWs) in the LIGO–Virgo–KAGRA (LVK) data.

In contrast to transient gravitational waves (GWs) emitted by binary black hole (and neutron star) mergers (R. Abbott et al. 2021a, 2021b, 2023), CWs have yet to be observed. These signals should be nearly monochromatic, with amplitude and frequency exhibiting small variations over year-long timescales.

CWs are expected from a time-varying nonaxisymmetric mass distribution in rotating neutron stars (M. Zimmermann & E. Szedenits 1979). This could be the result of strain in the elastic crust (G. Ushomirsky et al. 2000), accretion from a companion star (L. Bildsten 1998; A. Melatos & D. J. B. Payne 2005; F. Gittins & N. Andersson 2021), or a strong inner magnetic field (S. Bonazzola & E. Gourgoulhon 1996; C. Cutler 2002). Alternatively, the deformation could be caused by fluid oscillations, such as those due to r-modes (N. Andersson 1998; J. L. Friedman & S. M. Morsink 1998). However, the sources of CWs are likely to possess smaller mass quadrupoles compared to the sources of transient GWs, resulting in a weaker signal. Therefore, for signals to be detectable with current detectors, we need to consider nearby sources integrating long stretches (months' or years' worth) of detector data.

Observation of CWs from a neutron star would yield crucial insights into the star's structure and its equation of state (B. Haskell & M. Bejger 2023; F. Gittins 2024). Moreover, the form of the signal can be employed to test general relativity by measuring (or constraining) the presence of nonstandard polarizations (M. Isi et al. 2017; B. P. Abbott et al. 2019a). Possible mechanisms of CW emission from neutron stars are discussed in greater detail by K. Glampedakis & L. Gualtieri (2018) and K. Riles (2023).

The primary challenge for CW searches lies in accurately accounting for the various modulations that affect the signal received on the Earth. These include the Doppler effect due to the Earth's motion, the pulsar's rotational evolution (i.e., the slow spin-down), and relativistic effects. EM observations provide accurate measurements of the sky position and rotation parameters that allow us to predict and correct these modulations, thereby enhancing the search sensitivity.

Based on the knowledge of the source parameters, different strategies can be used to search for CW signals in LVK data (see K. Wette 2023 for a complete review of search methods). Targeted searches, the primary subject of this paper, aim to detect CWs from known pulsars, whose timing solutions can be calculated from known rotation phases and spin-down rates. Targeted searches use full-coherent methods that integrate data over long observation times, maintaining phase coherence over time. By assuming that the GW phase evolution follows the EM solution, the parameter space can be reduced to the unknown signal amplitude and polarization parameters. This assumption is relaxed in narrowband searches, which are performed in a narrow band around the frequency and spin-down rate (e.g., B. P. Abbott et al. 2017a, 2019b; R. Abbott et al. 2022b). However, because these searches decrease the sensitivity and increase the computational cost, they are often performed on fewer targets.

Although several CW searches have been conducted in recent years, targeting both isolated pulsars and those in binary systems, so far none of these searches have produced evidence of CWs (e.g., B. P. Abbott et al. 2017b, 2019c; L. Nieder et al. 2019, 2020; R. Abbott et al. 2020, 2021c, 2022a; A. Ashok et al. 2021). In the absence of a signal, these searches set upper limits on the GW amplitude and on the ellipticity, the physical parameter that quantifies the asymmetry in the mass distribution. In many cases, these limits are more stringent than (i.e., said to have "surpassed") the so-called spin-down limit, the theoretical limit calculated by assuming 100% of each pulsar's spin-down luminosity to be radiated through GWs. In the most recent targeted search (R. Abbott et al. 2021c, 2022a) considering data from the third observing run O3, 24 pulsars surpassed their spin-down limits, including the Crab and Vela pulsars, J0537−6910, and two millisecond pulsars, J0437−4715 and J0711−6830. Additionally, searches for an r-mode emission are described in B. Rajbhandari et al. (2021) for the Crab and in L. Fesik & M. A. Papa (2020), B. Rajbhandari et al. (2021), and R. Abbott et al. (2021d) for J0537−6910. Continuous improvements in detector sensitivity and data analysis techniques are progressively enhancing our ability to detect these faint signals.

In this paper, we present a targeted search for CWs from a set of 45 known pulsars, considering LIGO data from the initial part of the most recent (fourth) LVK observing run, known as O4a run. Pulsar selection is based on the available EM observations (see Section 4.2) and on the anticipated sensitivity for targeted searches near or below the spin-down limit. Considering two different emission models, we find no evidence of CW signals in the data, and we set upper limits on the amplitude and on the ellipticity for each target. Additionally, we perform a narrowband search for a subset of 16 targets and a search for nonstandard polarizations as predicted by the Brans–Dicke theory (C. Brans & R. H. Dicke 1961).

The paper is structured as follows. In Section 2, we briefly describe the expected signal for the emission models that we considered. The data analysis methods used in the paper are discussed in Section 3, while in Section 4 we describe the EM and GW data used for the analysis. The results are summarized and discussed in detail in Sections 5 and 6. Conclusions are reported in Section 7.

2. Signal Model

2.1. Standard Signal

For the targeted search, we assume that the GW signal is locked to the rotational phase of the pulsar obtained through EM observations. For an isolated triaxial star, rotating steadily about one of its principal axes of inertia, the GW emission frequency is twice the pulsar's spin frequency, frot. For the single-harmonic emission model, we search for signals at 2frot. However, there are additional mechanisms that could emit GWs at other frequencies. A superfluid component beneath the crust, rotating with a spin axis misaligned to the star's rotation axis, would produce an additional emission at the rotation frequency, so that overall we have a dual-harmonic emission at both once and twice the rotation frequency (D. I. Jones 2010). This would not impact the EM signature. Therefore, a dual-harmonic search is performed at both frot and 2frot. Additionally, a single-harmonic narrowband search is performed around 2frot, allowing for the possibility of a difference in rotation rate between the pulsation-producing magnetosphere and the part of the star responsible for the CW emission.

For the general dual-harmonic emission model, the signals h21 and h22 at frot and 2frot can be defined as (M. Pitkin et al. 2015)

where C21 and C22 are the dimensionless constants that give the component amplitudes, the angles (α, δ) are the R.A. and decl. of the source, the angles (ι, ψ) describe the orientation of the source's spin axis with respect to the observer in terms of inclination and polarization, and are phase angles at a defined epoch, and Φ(t) is the rotational phase of the source. The antenna functions and describe how the two polarization components (plus and cross) are projected onto the detector. These waveforms are detailed in D. I. Jones (2010) and used in R. Abbott et al. (2022a).

For the triaxial star described earlier, which only emits GWs at 2frot, C21 in Equation (1) is 0, which leaves only Equation (2), which contains C22. The amplitude h0 is defined as the amplitude of the circularly polarized signal observable for a source directly above or below the plane of the detector with its spin axis pointed directly toward or away from the detector. It can be calculated as

where d is the distance of the source. The equatorial ellipticity ε for the triaxial star emitting GWs at only 2frot is defined as

where Ixx, Iyy, and Izz are the source's principal moments of inertia, with the star rotating about the z-axis. From the ellipticity, the pulsar's mass quadrupole, Q22, can be calculated using (B. J. Owen 2005)

The spin-down limit of a source is given by

where is the rotation frequency derivative, or spin-down rate. This limit is the maximum GW amplitude allowed, assuming all the lost rotational energy of the star is due to conversion into GW energy. It should be noted that there are two types of spin-down rates: observed, which can be affected by the transverse velocity of the source (e.g., the Shklovskii effect, described in I. S. Shklovskii 1970), and intrinsic. Therefore, where possible, the intrinsic spin-down rate is used in calculating the spin-down limit.

2.2. Nonstandard Polarization Signal

In this paper, similar to the analysis in R. Abbott et al. (2022a), we search for GWs with polarizations predicted by the Brans–Dicke modification of general relativity (GR). Brans–Dicke theory includes two tensor polarizations, like in GR, and an additional scalar polarization. The dominant scalar radiation stems from the time-dependent dipole moment D. The dipole radiation occurs at the pulsar's rotational frequency frot. Assuming the dipole moment is along the x-axis, the amplitude of the signal is given by

where ζ is the parameter of the Brans–Dicke theory (see P. Verma 2021 for details).

3. Methods

In this section, we describe the data analysis methods used in this work: three independent pipelines for the targeted searches, one pipeline for the narrowband searches, and one pipeline for the nonstandard polarization searches. We use three targeted pipelines to compare independent results as our methods rely on different statistical approaches (Bayesian or frequentist) and on different preprocessings and handling of nonstationary or non-Gaussian noise disturbances in the data.

3.1. Time-domain Bayesian Method

The Continuous (gravitational) Wave Inference in Python (CWInPy) package is used to perform the Bayesian analysis (M. Pitkin 2022) following the method described in B. P. Abbott et al. (2019c) and summarized here.

First, a complex, slowly evolving heterodyne is used to remove the phase evolution of the source. This includes corrections for the relative motion of the source with respect to the detector and relativistic effects (R. J. Dupuis & G. Woan 2005). Then, a low-pass antialiasing filter is applied to the data to remove the upper sideband produced from the heterodyne and limit the possibility of disturbances from spectral lines. Next, the data is down-sampled, centered about the expected signal frequency, which has been shifted to 0 Hz by the heterodyne. For the dual-harmonic search, this method is repeated so that time series centered at both frot and 2frot are obtained.

There are several unknown signal parameters, with the amplitude being of primary interest. Bayesian inference is used to estimate these parameters as well as the evidence for the signal model. The priors used are the same as those detailed in Appendix 2 of B. P. Abbott et al. (2017b) except for the amplitude priors, for which we use flat priors with an upper cutoff much higher than the detector sensitivity. This value is 1.0 × 10−21 for all pulsars. The Bayesian stochastic sampling algorithm used is dynesty (J. Skilling 2004, 2006), as wrapped with bilby (G. Ashton et al. 2019), with 1024 live points (the number of points drawn from the prior). In the absence of a signal, we calculate 95% credible upper bounds derived from the posterior probability distributions.

3.1.1. Restricted Priors

For some pulsars, there is sufficient information to restrict our uninformative prior assumptions on their inclination and polarization angles, for example if we have EM observations of their pulsar wind nebulae or, in the case of J1952+3252, from proper motion measurements and observations of Hα lobes bracketing the bow shock (C.-Y. Ng & R. W. Romani 2004). In these cases, the parameter estimation is repeated and used along with the results from the original uninformed priors. Table 1 shows the pulsars for which we could restrict the priors and the values used for the restrictions. These values were obtained using the data from Table 2 of C.-Y. Ng & R. W. Romani (2008) and the methods described in Appendix B of B. P. Abbott et al. (2017b). Each prior range is assumed to be Gaussian about the given mean and standard deviation. The two values for ι are to incorporate the unknown rotation direction in the search by using a bimodal distribution. The additional ι2 is simply π − ι1 radians.

Table 1. Pulsars for Which We Can Restrict Their Orientation Priors Using EM Observations

PulsarΨ (rad)ι1 (rad)ι2 (rad)
J0205+64491.5760  ±  0.00781.5896  ±  0.02191.5519  ±  0.0219
J0534+2200 (Crab)2.1844  ±  0.00161.0850  ±  0.01492.0566  ±  0.0149
J0537−69102.2864  ±  0.03831.6197  ±  0.01651.5219  ±  0.0165
J0540−69192.5150  ±  0.01441.6214  ±  0.01061.5202  ±  0.0106
J0835−4510 (Vela)2.2799  ±  0.00151.1048  ±  0.01052.0368  ±  0.0105
J1952+32520.2007  ±  0.1501......
J2021+36510.7854  ±  0.02501.3788  ±  0.03901.7628  ±  0.0390
J2229+61141.7977  ±  0.04540.8029  ±  0.11002.3387  ±  0.1100

Note. Ψ is the polarization angle, and ι is the inclination angle, with ι2 = π − ι1 in radiant (rad).

Download table as:  ASCIITypeset image

3.1.2. Glitches

Glitches are transient events, where the normally stable pulsar spin suddenly increases in both rotation frequency and spin-down rate (C. M. Espinoza et al. 2011; M. Yu et al. 2013; A. Basu et al. 2022). These events are most common in younger, nonrecycled pulsars, with rarer glitches seen in some millisecond pulsars (I. Cognard & D. C. Backer 2004; J. W. McKee et al. 2016). Some searches (e.g., D. Keitel et al. 2019; R. Abbott et al. 2022b; A. G. Abac et al. 2024a) look for transient GWs in the aftermath of glitches.

When a glitch occurs in a pulsar during the course of a GW observing run, we assume that the GW phase is affected in the same way as the EM phase. However, since the phase offset is unknown at the time of the glitch, it is incorporated into the parameter inference for the Bayesian method. This adjustment is not necessary for pulsars that glitch before or after the observation period. Only two of the pulsars in our sample experienced glitches during the length of O4a: J0537−6910, with an epoch of ∼MJD 60223, and J0540−6919, with an epoch of ∼MJD 60150 (see, e.g., C. M. Espinoza et al. 2024; Y. Tuo et al. 2024).

3.2. 5-vector Targeted Pipeline

The 5-vector method is a frequentist data analysis approach (P. Astone et al. 2010), which has been used in the last decade in several LVK searches (B. P. Abbott et al. 2017b, 2019c; R. Abbott et al. 2020, 2022a). The 5-vector targeted pipeline relies on the Band-Sampled Data (BSD; O. J. Piccinni et al. 2018) framework, i.e., a database of subsampled complex time-domain files (so-called BSD files) that covers 10 Hz and spans 1 month of the original data set that allows us the computational cost of the analysis to be reduced.

To correct for signal phase modulations due to, for example, pulsar spin-down and the Doppler effect, a heterodyne method is used. Recently (L. D'Onofrio et al. 2023), the 5-vector method has also been applied to pulsars in binary systems with the implementation of Doppler correction due to orbital motion (P. Leaci & R. Prix 2015; P. Leaci et al. 2017). After the spin-down and Doppler/relativistic corrections, there is a residual modulation at the detector due to the Earth sidereal motion that splits the frequency of the expected CW signal.

Assuming the single-harmonic emission scenario, the central peak at twice the pulsar rotation frequency has four sidebands whose distance to the central peak is ±1, ± 2f⊕, where f⊕ is the Earth's sidereal rotation frequency. A 5-vector consists of a complex array with five components corresponding to the Fourier transform of the detector's data (data 5-vector) or of the antenna patterns (template 5-vectors for the plus and cross component) at the five signal frequencies. The 5-vector method defines two matched filters in the frequency domain, defined by the normalized scalar product between the data 5-vector and the two template 5-vectors. To extend the analysis to n detectors, the data 5-vector from each detector are combined together with coefficients that take into account the detector's sensitivity and observation time to construct the 5n-vectors (L. D'Onofrio et al. 2024). The two matched filters are linearly combined to define a detection statistic (L. D'Onofrio et al. 2023). To assess the significance of a candidate, a p-value is computed from the noise distribution of the detection statistic using off-source frequencies as the noise background in the analyzed frequency bands. In case of no detection, a mixed Bayesian-frequentist upper limit procedure (described in B. P. Abbott et al. 2019c) is used to set the upper limit on the amplitude, assuming a uniform prior and considering informative priors on the polarization parameters, if present (see Table 1).

It is not straightforward to generalize the 5-vector method to the dual-harmonic emission model due to the complexity of the expected CW signal. In this case, an additional analysis considering the emission at only the rotation frequency is performed; this would be a good approximation if the star's moment of inertia tensor is biaxial, with a small misalignment angle between its symmetry axis and the rotation axis (see D. I. Jones 2010).

3.3.  -/-/-statistic Method

The time-domain -/-/-statistic method utilizes the statistic derived in P. Jaranowski et al. (1998) and the statistic derived in P. Jaranowski & A. Królak (2010). The input data for this analysis are the heterodyned data used in time-domain Bayesian analysis. The statistic is employed when the amplitude, phase, and polarization of the signal are unknown, whereas the statistic is applied when only the amplitude and phase are unknown, and the polarization is known (as described in Section 3.1.1). These methods have been used in several analyses of LIGO and Virgo data (J. Abadie et al. 2011; J. Aasi et al. 2014; B. P. Abbott et al. 2017b; R. Abbott et al. 2020, 2022a). Additionally, the method incorporates the statistic developed in P. Verma (2021) to search for dipole radiation in the Brans–Dicke theory. The -statistic search was first performed in the LIGO and Virgo analysis presented in R. Abbott et al. (2022a).

The three statistics are derived by calculating the maximum likelihood estimators of the signal's constant amplitude parameters. This is done by maximizing the likelihood function and then substituting the amplitude values with their respective maximum likelihood estimators. As a result, we obtain statistics that are independent of the amplitudes. In this method, a signal is detected if the value of the , , or statistic exceeds a certain threshold corresponding to an acceptable false-alarm probability. We consider a false-alarm probability of 1% for a signal to be deemed significant. The , , and statistics are computed for each detector and each interglitch period separately. The results from different detectors or interglitch periods are then combined incoherently by summing the respective statistics. When the values of the statistics are not statistically significant, we set an upper limit with a frequentist approach on the amplitude of the GW signal.

3.4. 5n-vector Narrowband Pipeline

The 5n-vector narrowband pipeline makes use of the 5n-vector as in P. Astone et al. (2014) and follows the same principle of the method described in Section 3.2.

While the former searches for a CW tightly locked to the EM emission, this assumption is here relaxed, searching for CWs in a narrow frequency and spin-down range around twice the values inferred from EM observations, namely

with δ = 10−3 (R. Abbott et al. 2022b) and an analogous expression for .

The method makes use of the Short Fourier Data Base (SFDB; P. Astone et al. 2005), which is a collection of short-duration (2048 s) fast Fourier transforms overlapped by half.

For every target, data are Doppler corrected in the time domain with a nonuniform resampling that is independent of the CW frequency, then they are subsampled at 1 Hz. At this point, the time series are match-filtered (using 5-vectors as for the targeted search) to estimate the two CW polarizations using a template bank in the f − ḟ space. We build the template bank considering bin widths equal to the inverse of the time series duration for the frequency, while we use its inverse squared for the spin-down grid. Higher-order spin-down terms, if provided in the ephemerides, are fixed at twice the rotation terms to track the GW frequency evolution over time, without exploring any additional template (P. Astone et al. 2014).

Then, the matched filter results from different detectors are coherently combined (S. Mastrogiovanni et al. 2017) to evaluate the detection statistic. From the collection of statistic values, we order them along the frequency axis to select the maximum every 10−4 Hz, marginalizing over the spin-down. We use a p-value threshold (i.e., the tail probability of the noise-only distribution) of 1% to determine whether a selected point has to be considered an outlier. The noise-only distribution is inferred from the tail of the histogrammed nonmaxima values, which is fitted with an exponential (A. Singhal et al. 2019). We use here a threshold of 1%, similar to previous searches (e.g., R. Abbott et al. 2022b), after taking into account the trial factor.

If CW detection is not claimed, we set the 95% confidence level upper limit through software-injection campaigns.

4. Data Sets

4.1. GW Data Set

The considered data set is the first part of the fourth observing run, known as O4a, of the LIGO Livingston (L1) and LIGO Hanford (H1) detectors.342 The O4a run took place between 2023 May 24 15:00:00 UTC and ended January 16 2024 16:00:00 UTC. The duty factors for L1 and H1 were 69.0% and 67.5%, respectively. The Virgo detector has not been considered since it joined the O4 run on 2024 April 10 while the KAGRA detector is planned to join O4 by the end of the run. For a description of the upgrades to the Advanced LIGO (E. Capote et al. 2025), Advanced Virgo, and KAGRA detectors in preparation for the O4 run, we refer to Appendix A in A. G. Abac et al. (2024b).

The LIGO detectors are calibrated using photon radiation pressure actuation, where an amplitude-modulated laser beam is directed onto the end test masses, causing a known change in the arm length from the equilibrium position (S. Karki et al. 2016; A. D. Viets et al. 2018). For the O4a data used in this analysis, the worst 1σ calibration uncertainty is within 10% in amplitude, and 10° in phase, over the range 10–2000 Hz. The uncertainty at specific frequencies or times can be significantly smaller.

The data set used underwent cleaning processes depending on the data framework used by each pipeline (we refer to S. Soni et al. (2024) for more details on data quality for CW searches). For the Bayesian pipeline, outliers are removed before the heterodyne using a median-absolute-deviation method as described in Chapter 3 of B. Iglewicz & D. Hoaglin (1993) with a threshold of 3.5. Concerning the 5-vector targeted pipeline, two cleaning steps are applied. First, short-duration time-domain disturbances are identified and subtracted from the data when the SFDB—from which BSD files are built—are created (this cleaning step is shared with the 5-vector narrowband search, which also starts from the SFDB; P. Astone et al. 2005). A second cleaning step is applied to the BSD files (O. J. Piccinni et al. 2018), removing large time-domain outliers, which were not visible in the full band time series, with quadratic values larger than 10 times the quadratic sum of the medians of the data's real and imaginary parts (computed over nonzero samples). The -, -, or -statistic method performs further cleaning of the fine heterodyne data through the Grubbs test (see Appendix D of J. Abadie et al. 2011).

4.2. EM Data Set

The timing solutions used as inputs for the GW search were produced with EM data in the gamma-ray, X-ray, and radio wavelengths. The gamma-ray timing solutions were obtained from Fermi-Large Area Telescope (LAT; W. B. Atwood et al. 2009); X-ray timing solutions were obtained from Chandra (M. C. Weisskopf et al. 2002) and the Neutron Star Interior Composition Explorer (NICER; K. C. Gendreau et al. 2016), while the radio timing solutions were obtained from the Nançay Radio Telescope (L. Guillemot et al. 2023), the Jodrell Bank Observatory (JBO), Instituto Argentino de Radioastronomía (IAR; G. Gancio et al. 2020), the Mount Pleasant Radio Observatory (D. R. Lewis et al. 2003), the Five-hundred-meter Aperture Spherical Telescope (FAST; R. Smits et al. 2009), and the Canadian Hydrogen Intensity Mapping Experiment (CHIME; M. Amiri et al. 2021).

For the radio-emitting pulsars, we used pazi or rfifind tasks in PSRCHIVE (W. van Straten et al. 2012) and PRESTO (S. Ransom 2011) packages, respectively, to mitigate radio-frequency interferences. Next, we folded observations with prepfold or fold_psrfits tasks in PRESTO and PSRFITS UTILS (A. W. Hotan et al. 2004) packages. We then cross-correlated the folded profiles with a noise-free template profile with a high signal-to-noise ratio to obtain the times of arrivals (ToAs). Next, we select ToAs during the course of the O4a run, so the solutions are valid for this time range. We use TEMPO (D. Nice et al. 2015), TEMPO2 (R. T. Edwards et al. 2006; G. B. Hobbs et al. 2006; G. Hobbs et al. 2009), or PINT (J. Luo et al. 2019; J. Luo et al. 2021) to characterize the rotation of each pulsar by fitting the ToAs to a Taylor series expansion:

where t0 is the reference epoch, ϕ0 is the phase at t0, and frot, and are the rotation frequency of the pulsar, and its first and second derivatives, respectively. If higher-order derivatives are measured, we also include the corresponding terms in the Taylor expansion.

During a glitch, the rotation frequency abruptly increases. This glitch-induced alteration in the rotational phase can be taken into account in the timing model as (M. Yu et al. 2013)

The uncertainty on the glitch epoch tg is counteracted here by Δϕ, and the step changes in frot, , and at tg are represented by , , and . Finally, denotes temporal frequency increases that decay in days.

The ToA-fitting process also provides the astrometric parameters of each pulsar and the orbital parameters for binary pulsars (D. R. Lorimer & M. Kramer 2004). Uncertainties in the values of the pulsar and orbital parameters derived from fitting the ToAs are not taken into account in targeted/narrowband searches.

For many pulsars, their distances (G. Hobbs et al. 2004) are based on the observed dispersion measure using the Galactic electron density distribution model YMW16 (J. M. Yao et al. 2017). The uncertainties in these measurements can be as large as a factor of 2. For other pulsars, the distance can be determined by measuring the parallax with the timing solution (R. Smits et al. 2011) or, if the pulsar is in a binary system, the orbital period derivative (J. P. W. Verbiest et al. 2008). The first method usually results in an uncertainty ranging from 5% to 50% (M. Shamohammadi et al. 2024). The second one offers significantly higher accuracy, achieving uncertainties as low as 0.1% (J. P. W. Verbiest et al. 2008; D. J. Reardon et al. 2016). Other methods such as very long baseline interferometry (R. Lin et al. 2023) were also used to determine pulsar distances.

For this analysis, we selected pulsars with a rotation frequency close to or greater than 10Hz, to lie in the bandwidth of the LIGO detectors, and with an expected targeted search sensitivity for the strain amplitude within a factor 3 of the spin-down limit (see Figure 1 for the pulsars' sky locations and Table 2 for the list of analyzed pulsars). Out of the 45 pulsars analyzed in this work, 11 pulsars belong to binary systems, and there are 10 millisecond pulsars with frequencies higher than 100 Hz.

Figure 1. Refer to the following caption and surrounding text.

Figure 1. Sky location in equatorial coordinates of the analyzed targets. All targets are selected for the targeted search, while the triangles are the targets also analyzed with the narrowband search.

Download figure:

Standard image High-resolution image

Table 2. Results for the Targeted Search on the Set of 45 Known Pulsars for the Three Considered Pipelines Described in Section 3

Pulsar Namefrot Distance Analysis ε95% Statistic@Statistic#
(J2000)(Hz)(s s−1)(kpc) Method  (kg m2)   l=2,m=1,2l=2,m=2
J0030+0451α205.5−4.2 × 10−160.3 (a)3.5 × 10−27Bayesian1.0 × 10−261.9 × 10−81.5 × 10302.91.1 × 10−265.2 × 10−27−5.1−10
      statistic1.4 × 10−262.7 × 10−8J0537 2.1 × 10304.11.6 × 10−266.5 × 10−270.470.04
     5n-vector9.3 × 10−271.7 × 10−81.3 × 10302.65.8 × 10−27⋯0.610.038
J0058−7218β45.94−6.1 × 10−1159.70 (b)1.56 × 10−26Bayesian8.8 × 10−275.9 × 10−54.5 × 10330.562.4 × 10−264.4 × 10−27−5.3−10
      statistic4.9 × 10−273.3 × 10−52.5 × 10330.313.5 × 10−266.2 × 10−270.320.86
     5n-vector6.5 × 10−274.4 × 10−53.4 × 10330.419.4 × 10−27⋯10.80
J0117+5914γ9.86−5.69 × 10−131.77 (c)1.1 × 10−25Bayesian1.7 × 10−257.4 × 10−45.8 × 10341.67.7 × 10−238.8 × 10−26−3.8−4.9
      statistic2.3 × 10−251.0 × 10−37.8 × 10342.28.4 × 10−231.2 × 10−250.390.06
     5n-vector3.2 × 10−251.4 × 10−31.1 × 10352.9⋯⋯⋯0.87
J0205+6449γ15.22−4.49 × 10−113.20 (d)4.33 × 10−25Bayesian3.2(4.2) × 10−261.0(1.4) × 10−48.0(10) × 10330.073(0.096)5.6(4.8) × 10−251.5(2.0) × 10−26−4.7(-4.5)-8.3(-8.3)
      statistic2.4(4.4) × 10−260.8(1.5) × 10−46.0(10) × 10330.055(0.10)1.3(0.56) × 10−241.1(2.2) × 10−260.24(0.44)0.75(0.47)
     5n-vector2.2(3.4) × 10−260.72(1.1) × 10−45.6(8.6) × 10330.051(0.078)(5.7) × 10−25⋯0.860.84
J0437−4715δ173.69−4.15 × 10−160.16 (e)7.95 × 10−27Bayesian7.2 × 10−278.8 × 10−96.8 × 10290.901.1 × 10−263.4 × 10−27-5.4-10
      statistic9.0 × 10−271.1 × 10−88.5 × 10291.11.5 × 10−264.4 × 10−270.310.26
     5n-vector5.8 × 10−277.1 × 10−95.5 × 10290.325.1 × 10−27⋯0.750.38
J0534+2200γ29.95−3.78 × 10−102.00 (f)1.43 × 10−24Bayesian1.1(0.9) × 10−265.9(5.0) × 10−64.6(3.9) × 10320.0078(0.0067)6.6(6.0) × 10−265.4(4.5) × 10−27−5.1(-5.2)-9.5(-9.2)
      statistic1.5(1.2) × 10−268.0(6.3) × 10−66.3(4.9) × 10320.011(0.0085)9.4(7.0) × 10−268.8(6.1) × 10−270.19(0.18)0.30(0.36)
     5n-vector1.1(1.1) × 10−265.7(5.8) × 10−64.3(4.4) × 10320.0074(0.0075)2.6(2.6) × 10−26⋯0.970.44
J0537−6910β62.03−1.99 × 10−1049.70 (g)2.91 × 10−26Bayesian6.4(8.9) × 10−272.0(2.7) × 10−51.5(2.1) × 10330.22(0.31)2.4(1.4) × 10−263.0(4.7) × 10−27−5.5(-5.3)-10.2(-10.4)
      statistic8.8(4.5) × 10−272.7(1.4) × 10−52.0(1.0) × 10330.29(0.15)3.2(9.0) × 10−262.0(2.9) × 10−270.24(0.25)1.00(0.92)
     5n-vector0.79(1.2) × 10−262.4(3.7) × 10−51.9(2.8) × 10330.27(0.41)0.98(1.5) × 10−26⋯0.450.34
J0540−6919β19.77−1.87 × 10−1049.70 (g)4.99 × 10−26Bayesian2.7(3.4) × 10−268.1(10) × 10−46.3(7.9) × 10340.54(0.69)1.9(1.6) × 10−251.4(1.7) × 10−26−4.6(-4.2)-8.5(-8.4)
      statistic2.9(2.9) × 10−268.6(8.6) × 10−46.7(6.7) × 10340.58(0.58)2.9(8.8) × 10−251.4(1.5) × 10−260.53(0.34)0.32(0.45)
     5n-vector1.6(2.4) × 10−264.8(7.3) × 10−43.7(5.6) × 10340.27(0.41)3.7(5.7) × 10−25⋯0.660.76
J0614−3329α317.59−1.76 × 10−150.63 (h)3.01 × 10−27Bayesian1.1 × 10−261.6 × 10−81.2 × 10303.68.9 × 10−275.4 × 10−27−5.1-10
      statistic1.1 × 10−263.1 × 10−82.4 × 10303.61.0 × 10−265.3 × 10−270.660.31
     5n-vector8.0 × 10−271.2 × 10−89.1 × 10292.77.3 × 10−27⋯0.150.21
J0737−3039Aα44.05−3.41 × 10−151.10 (i)6.45 × 10−27Bayesian8.9 × 10−271.2 × 10−69.2 × 10311.41.8 × 10−264.3 × 10−27−5.3-10
      statistic7.5 × 10−271.0 × 10−67.6 × 10311.21.7 × 10−256.4 × 10−270.990.87
     5n-vector6.2 × 10−278.3 × 10−76.4 × 10310.961.4 × 10−26⋯0.490.92
J0835−4510epsilon11.19−1.57 × 10−110.28 (i)3.41 × 10−24Bayesian9.0(7.7) × 10−264.7(4.1) × 10−53.7(3.1) × 10330.026(0.023)2.9(2.4) × 10−244.1(3.7) × 10−26−4.1(-4.1)-7.1(-7.1)
      statistic8.2(8.2) × 10−264.4(4.3) × 10−53.4(3.3) × 10330.024(0.024)3.5(2.3) × 10−244.0(3.9) × 10−260.65(0.81)0.45(0.19)
     5n-vector8.5(8.7) × 10−264.5(4.6) × 10−53.5(3.6) × 10330.025(0.025)2.7(2.8) × 10−24⋯0.660.17
J1231−1411α271.45−5.92 × 10−160.42 (c)2.83 × 10−27Bayesian1.2 × 10−261.6 × 10−81.2 × 10304.11.2 × 10−265.9 × 10−27−4.6-9.9
      statistic1.2 × 10−261.6 × 10−81.2 × 10304.11.7 × 10−265.9 × 10−270.140.14
     5n-vector6.6 × 10−278.9 × 10−96.9 × 10291.45.4 × 10−27⋯0.680.51
J1412+7922β17.18−9.72 × 10−133.30 (j)5.81 × 10−26Bayesian3.1 × 10−268.2 × 10−56.3 × 10330.539.9 × 10−251.4 × 10−26−4.7-8.0
      statistic2.5 × 10−266.6 × 10−55.1 × 10330.431.2 × 10−241.3 × 10−260.170.57
     5n-vector2.3 × 10−266.4 × 10−54.9 × 10330.44.0 × 10−25⋯10.44
J1537+1155α26.38−1.65 × 10−150.93 (k)6.82 × 10−27Bayesian1.3 × 10−264.0 × 10−63.1 × 10321.88.5 × 10−265.6 × 10−27−5.2-9.6
      statistic1.3 × 10−264.0 × 10−63.1 × 10321.87.1 × 10−266.4 × 10−270.620.91
     5n-vector1.6 × 10−265.0 × 10−63.9 × 10322.33.3 × 10−26⋯0.900.22
J1623−2631γ90.29−5.26 × 10−151.85 (l)3.33 × 10−27Bayesian8.5 × 10−274.6 × 10−73.5 × 10312.61.4 × 10−264.1 × 10−27−5.2-10
      statistic1.0 × 10−265.4 × 10−74.1 × 10313.11.2 × 10−265.2 × 10−270.790.19
     5n-vector5.6 × 10−273.0 × 10−72.3 × 10311.26.8 × 10−27⋯0.780.68
J1719−1438α172.71−2.22 × 10−160.34 (c)2.72 × 10−27Bayesian8.7 × 10−272.3 × 10−81.8 × 10303.21.3 × 10−263.7 × 10−27−5.3-10
      statistic7.4 × 10−272.0 × 10−81.6 × 10302.82.0 × 10−263.5 × 10−270.200.51
     5n-vector7.4 × 10−272.0 × 10−81.5 × 10302.67.0 × 10−27⋯0.330.14
J1744−1134α245.43−4.34 × 10−160.40 (e)2.72 × 10−27Bayesian9.2 × 10−271.4 × 10−81.1 × 10303.41.1 × 10−264.4 × 10−27−5.0-10
      statistic1.4 × 10−262.1 × 10−81.7 × 10305.21.2 × 10−267.1 × 10−270.340.02
     5n-vector5.7 × 10−278.8 × 10−96.8 × 10291.98.8 × 10−27⋯0.110.69
J1745−0952α51.61−2.3 × 10−160.23 (c)7.5 × 10−27Bayesian1.1 × 10−262.1 × 10−71.6 × 10311.42.6 × 10−265.0 × 10−27−5.3-9.7
      statistic8.6 × 10−271.6 × 10−71.3 × 10311.41.6 × 10−264.4 × 10−270.660.62
     5n-vector1.2 × 10−262.5 × 10−71.9 × 10311.61.1 × 10−26⋯0.910.023
J1756−2251γ35.14−1.26 × 10−150.73 (m)6.6 × 10−27Bayesian1.1 × 10−261.5 × 10−61.2 × 10321.65.6 × 10−265.4 × 10−27−5.1-9.4
      statistic1.6 × 10−262.2 × 10−61.7 × 10322.36.3 × 10−267.1 × 10−270.170.22
     5n-vector1.0 × 10−261.5 × 10−61.1 × 10321.63.1 × 10−26⋯0.0310.31
J1809−1917γ12.08−3.73 × 10−123.27 (c)1.37 × 10−25Bayesian8.1 × 10−264.3 × 10−43.3 × 10340.593.0 × 10−243.9 × 10−26−4.2-7.0
      statistic7.9 × 10−264.7 × 10−43.6 × 10340.654.3 × 10−243.9 × 10−260.520.05
     5n-vector4.0 × 10−262.1 × 10−41.6 × 10340.282.7 × 10−24⋯0.0790.64
J1811−1925β15.46−1.05 × 10−115.00 (n)1.33 × 10−25Bayesian2.6 × 10−261.3 × 10−41.0 × 10340.207.3 × 10−251.3 × 10−26−4.8-8.1
      statistic3.3 × 10−261.7 × 10−41.3 × 10340.254.3 × 10−251.7 × 10−260.670.85
     5n-vector3.3 × 10−261.6 × 10−41.2 × 10340.246.1 × 10−25⋯0.230.23
J1813−1246κ20.80−7.6 × 10−122.63 (o)1.85 × 10−25Bayesian1.5 × 10−262.2 × 10−51.7 × 10330.0821.4 × 10−257.8 × 10−27−5.0-9.1
      statistic3.1 × 10−264.7 × 10−53.4 × 10330.0551.3 × 10−251.5 × 10−260.720.07
     5n-vector1.6 × 10−262.3 × 10−51.8 × 10330.0877.7 × 10−26⋯0.730.50
J1813−1749θ22.35−6.34 × 10−116.15 (c)2.21 × 10−25Bayesian2.4 × 10−266.9 × 10−55.3 × 10330.118.8 × 10−261.2 × 10−26−4.4-8.7
      statistic3.1 × 10−268.9 × 10−56.8 × 10330.145.3 × 10−261.5 × 10−260.910.02
     5n-vector1.8 × 10−265.3 × 10−54.1 × 10330.0836.6 × 10−26⋯0.530.12
J1823−3021Aγ183.82−1.14 × 10−138.02 (l)2.5 × 10−27Bayesian7.5 × 10−274.2 × 10−73.3 × 10313.01.6 × 10−264.0 × 10−27−5.2-9.7
      statistic4.9 × 10−272.7 × 10−72.2 × 10312.02.3 × 10−263.5 × 10−270.120.78
     5n-vector4.5 × 10−272.5 × 10−71.9 × 10311.86.5 × 10−27⋯0.400.91
J1824−2452Aα327.41−1.73 × 10−135.37 (l)3.45 × 10−27Bayesian8.1 × 10−279.6 × 10−87.4 × 10302.41.5 × 10−263.7 × 10−27−5.2-10.0
      statistic7.0 × 10−278.3 × 10−86.4 × 10302.42.2 × 10−263.3 × 10−270.100.70
     5n-vector6.3 × 10−277.4 × 10−85.8 × 10301.87.2 × 10−27⋯0.170.64
J1826−1334γ9.85−7.31 × 10−123.61 (c)1.93 × 10−25Bayesian1.5 × 10−251.3 × 10−31.0 × 10350.794.9 × 10−237.0 × 10−26−3.8-5.7
      statistic1.5 × 10−251.3 × 10−31.0 × 10350.793.7 × 10−237.6 × 10−260.420.39
     5n-vector4.3 × 10−253.7 × 10−32.9 × 10352.2⋯⋯⋯0.42
J1828−1101γ13.88−2.85 × 10−124.77 (c)7.67 × 10−26Bayesian4.2 × 10−262.5 × 10−41.9 × 10340.555.3 × 10−242.1 × 10−26−4.6-7.4
      statistic3.5 × 10−262.1 × 10−41.6 × 10340.463.1 × 10−242.1 × 10−260.980.75
     5n-vector5.9 × 10−263.5 × 10−42.7 × 10340.783.8 × 10−24⋯0.660.19
J1831−0952γ14.87−1.84 × 10−123.68 (c)7.7 × 10−26Bayesian4.2 × 10−261.6 × 10−41.3 × 10340.546.4 × 10−251.9 × 10−26−4.6-8.0
      statistic3.2 × 10−261.2 × 10−40.99 × 10340.419.1 × 10−252.1 × 10−260.440.64
     5n-vector2.9 × 10−261.2 × 10−49.0 × 10330.385.2 × 10−25⋯0.850.77
J1833−0827γ11.72−1.26 × 10−124.50 (i)5.88 × 10−26Bayesian6.2 × 10−264.8 × 10−43.7 × 10341.02.7 × 10−243.0 × 10−26−4.4-7.2
      statistic5.4 × 10−264.2 × 10−43.2 × 10340.873.5 × 10−242.8 × 10−260.400.66
     5n-vector6.4 × 10−264.9 × 10−43.8 × 10341.12.4 × 10−24⋯0.430.35
J1837−0604γ10.38−4.84 × 10−124.78 (c)1.15 × 10−25Bayesian7.8 × 10−268.2 × 10−46.3 × 10340.685.6 × 10−243.6 × 10−26−4.4-7.0
      statistic8.8 × 10−269.3 × 10−47.1 × 10340.771.2 × 10−243.6 × 10−270.970.95
     5n-vector9.3 × 10−269.8 × 10−47.6 × 10340.804.7 × 10−24⋯0.980.43
J1838−0655β14.18−9.9 × 10−126.60 (p)1.02 × 10−25Bayesian6.6 × 10−265.1 × 10−44.0 × 10340.651.3 × 10−243.4 × 10−26−3.8-6.7
      statistic7.1 × 10−265.5 × 10−44.3 × 10340.702.6 × 10−243.6 × 10−260.150.04
     5n-vector4.2 × 10−263.3 × 10−42.5 × 10340.418.3 × 10−25⋯0.690.29
J1849−0001β25.96−9.54 × 10−127.00 (q)6.98 × 10−26Bayesian1.4 × 10−263.3 × 10−52.6 × 10330.191.3 × 10−256.0 × 10−27−5.2-9.1
      statistic1.7 × 10−264.0 × 10−53.2 × 10330.231.4 × 10−258.1 × 10−270.180.47
     5n-vector1.4 × 10−263.4 × 10−52.6 × 10330.26.3 × 10−26⋯0.0780.77
J1856+0245γ12.36−9.49 × 10−126.32 (c)1.12 × 10−25Bayesian5.8 × 10−265.7 × 10−44.4 × 10340.521.3 × 10−242.9 × 10−26−4.4-7.8
      statistic5.3 × 10−265.2 × 10−44.0 × 10340.481.8 × 10−242.5 × 10−260.690.66
     5n-vector5.2 × 10−265.1 × 10−44.0 × 10340.471.4 × 10−24⋯0.850.40
J1913+1011γ27.85−2.62 × 10−124.61 (c)5.36 × 10−26Bayesian1.1 × 10−261.6 × 10−51.2 × 10330.216.5 × 10−265.2 × 10−27−5.2-9.7
      statistic1.2 × 10−261.7 × 10−51.3 × 10330.234.0 × 10−265.8 × 10−270.960.60
     5n-vector1.2 × 10−261.8 × 10−51.4 × 10330.233.6 × 10−26⋯0.890.35
J1925+1720γ13.22−1.83 × 10−125.05 (c)5.94 × 10−26Bayesian5.9 × 10−264.1 × 10−43.1 × 10341.02.2 × 10−242.8 × 10−26−4.2-6.2
      statistic5.4 × 10−263.8 × 10−42.8 × 10340.924.2 × 10−242.9 × 10−260.030.33
     5n-vector6.8 × 10−264.6 × 10−43.6 × 10341.11.5 × 10−24⋯0.290.014
J1935+2025γ12.48−9.47 × 10−124.60 (c)1.53 × 10−25Bayesian5.5 × 10−263.8 × 10−43.0 × 10340.361.5 × 10−242.6 × 10−26−4.6-7.8
      statistic4.8 × 10−263.3 × 10−42.6 × 10340.317.7 × 10−252.3 × 10−260.980.69
     5n-vector3.7 × 10−262.6 × 10−42.0 × 10340.241.6 × 10−24⋯0.570.81
J1952+3252γ25.30−3.74 × 10−123.00 (i)1.03 × 10−25Bayesian1.1(1.0) × 10−261.2(1.1) × 10−59.2(8.2) × 10320.10(0.09)1.9(68) × 10−254.2(2.8) × 10−27−5.3(-5.4)-8.6(-8.2)
      statistic1.0 × 10−261.1 × 10−58.4 × 10320.0912.2 × 10−254.4 × 10−270.030.77
     5n-vector1.1(1.1) × 10−261.2(1.2) × 10−59.1(9.5) × 10320.1(0.11)5.1(5.1) × 10−26⋯0.540.72
J2016+3711ζ19.68−2.81 × 10−116.10 (r)1.58 × 10−25Bayesian2.0 × 10−267.4 × 10−55.7 × 10330.131.6 × 10−258.8 × 10−27−5.0-9.1
      statistic3.2 × 10−261.2 × 10−49.1 × 10330.211.2 × 10−251.7 × 10−260.840.07
     5n-vector3.0 × 10−262.4 × 10−41.9 × 10340.413.8 × 10−25⋯0.620.020
J2021+3651γ9.64−8.89 × 10−121.80 (s)4.3 × 10−25Bayesian2.3(2.2) × 10−251.0(1.0) × 10−38.0(7.7) × 10340.53(0.51)6.0(5.0) × 10−231.1(1.1) × 10−25−3.1(-3.2)-4.9(-5.2)
      statistic2.2(2.1) × 10−250.96(0.95) × 10−37.7(7.4) × 10340.51(0.49)2.0(1.3) × 10−221.1(1.1) × 10−250.09(0.65)0.08(0.13)
     5n-vector2.8(4.0) × 10−251.3(1.8) × 10−31(1.4) × 10350.66(0.93)⋯⋯⋯1
J2022+3842β20.59−3.65 × 10−1110.00 (t)1.07 × 10−25Bayesian1.9 × 10−261.0 × 10−48.1 × 10330.181.4 × 10−258.6 × 10−27−5.0-9.1
      statistic1.7 × 10−268.9 × 10−57.2 × 10330.168.0 × 10−269.0 × 10−270.900.60
     5n-vector1.4 × 10−267.5 × 10−55.8 × 10330.136.3 × 10−26⋯10.77
J2043+2740γ10.40−1.37 × 10−131.48 (c)6.26 × 10−26Bayesian1.1 × 10−253.7 × 10−42.8 × 10341.81.2 × 10−235.2 × 10−26−4.1-6.0
      statistic1.0 × 10−253.4 × 10−42.5 × 10341.61.3 × 10−234.4 × 10−260.420.78
     5n-vector7.9 × 10−262.5 × 10−42.0 × 10341.26.9 × 10−24⋯0.260.65
J2124−3358α202.79−2.94 × 10−160.41 (e)2.37 × 10−27Bayesian9.2 × 10−272.2 × 10−81.7 × 10303.99.9 × 10−274.8 × 10−27−5.0-10
      statistic9.7 × 10−272.3 × 10−81.8 × 10304.19.5 × 10−274.8 × 10−270.780.23
     5n-vector5.4 × 10−271.3 × 10−89.8 × 10291.35.1 × 10−27⋯0.730.62
J2214+3000η320.59−1.31 × 10−150.60 (u)2.71 × 10−27Bayesian9.0 × 10−271.2 × 10−89.6 × 10293.38.6 × 10−274.6 × 10−27−5.3-10
      statistic4.5 × 10−270.60 × 10−84.8 × 10291.79.0 × 10−275.8 × 10−270.810.87
     5n-vector5.4 × 10−277.4 × 10−95.7 × 10291.87.5 × 10−27⋯0.110.84
J2222−0137α30.47−4.99 × 10−160.27 (v)1.22 × 10−26Bayesian1.3 × 10−269.0 × 10−76.9 × 10311.14.1 × 10−266.4 × 10−27−5.1-9.7
      statistic2.2 × 10−261.5 × 10−61.2 × 10321.95.0 × 10−261.1 × 10−260.930.06
     5n-vector1.3 × 10−269.1 × 10−77.0 × 10313.34.9 × 10−26⋯0.0160.15
J2229+6114γ19.36−2.9 × 10−113.00 (w)3.29 × 10−25Bayesian1.5(0.9) × 10−262.8(1.8) × 10−52.2(1.4) × 10330.045(0.028)1.8(1.6) × 10−256.8(4.6) × 10−27−5.1(-5.3)-9.2(-9.4)
      statistic1.5(0.63) × 10−262.9(1.2) × 10−52.2(0.98) × 10330.046(0.032)2.9(2.5) × 10−254.7(2.5) × 10−270.24(0.41)0.99(0.99)
     5n-vector1.3(0.99) × 10−262.4(1.9) × 10−51.8(1.4) × 10330.037(0.029)3.9(3.2) × 10−27⋯0.750.88

Note. The last two columns refer to the significance of the data against the noise hypothesis (@ for the dual-harmonic emission model, # for the single-harmonic emission model). For the Bayesian method, the columns show the base-10 logarithm of the Bayesian odds, comparing coherent signal model modes to incoherent signal models. For the F-statistic and the 5n-vector methods, the columns show the p-value (for the 5n-vector method, considering a signal at just the l = 2 , m = 1 mode for the dual-harmonic emission model). The following is a list of references for pulsar ephemeris data used in this analysis: Nançay: α; NICER: β; JBO: γ; IAR: δ; Hobart: epsilon; FAST: ζ; CHIME: η; Chandra: θ; Fermi-LAT: κ.References. The following is a list of references for pulsar distances and intrinsic period derivatives, and they should be consulted for information on the associated uncertainties on these quantities: (a) H. Ding et al. (2023), (b) J. Storm et al. (2004), (c) J. M. Yao et al. (2017), (d) D. A. Roberts et al. (1993), (e) D. J. Reardon et al. (2016), (f) V. Trimble (1968), (g) A. R. Walker (2012), (h) C. G. Bassa et al. (2016), (i) J. P. W. Verbiest et al. (2012), (j) S. Mereghetti et al. (2021), (k) H. Ding et al. (2021), (l) H. Baumgardt & E. Vasiliev (2021), (m) R. D. Ferdman et al. (2014), (n) D. A. Green et al. (1988), (o) D. F. Torres et al. (2019), (p) L. C.-C. Lin et al. (2009), (q) H. E. S. S. Collaboration et al. (2018), (r) Q.-C. Liu et al. (2024), (s) A. Kirichenko et al. (2015), (t) Z. Arzoumanian et al. (2011), (u) L. Guillemot et al. (2016), (v) Y. J. Guo et al. (2021), (w) J. P. Halpern et al. (2001).

Download table as:  ASCIITypeset images: 1 2 3 4

5. Results

5.1. Targeted Searches

We found no statistical evidence of a CW signal in the O4a data for any of the analyzed targets. In this section, we present the results of the targeted search conducted using three different analysis methods across the full set of 45 pulsars.

The results are shown in Table 2. The 95% upper limit is given for the single-harmonic search along with the mass quadrupole and ellipticity ε95% upper limits calculated using the distance listed in the table and a fiducial moment of inertia Izz = 1038 kg m2. Uncertainties on these parameters are not taken into account, and for reference, we report the used best-fit distance values in Table 2 provided by the EM observation. However, for pulsars that did not surpass their spin-down limits, these Q22 and epsilon upper limits are unphysical since they would lead to spin-down rates that are greater than their measured values. From the upper limit on the amplitude, we also compute the spin-down ratio as with defined in Equation (6). The upper limits for the dual-harmonic search are included as and . Finally, for the Bayesian method, the odds of a coherent signal versus incoherent noise are given for both the single and dual-harmonic searches. For the -statistic and 5-vector method, to assess the statistical significance of a candidate and quantify the consistency with the assumption of just noise, we report the p-value.

For the two glitching pulsars, J0537−6910 and J0540−6919, the Bayesian results are produced when incorporating an additional phase offset in the parameter inference, while for the -statistic and 5-vector method, an incoherent approach is used summing the statistics from the interglitch periods. In cases with sufficient observations of the pulsar wind nebulae, results using restricted priors of inclination and polarization angles are listed in parentheses.

Figure 2 shows the upper limits from the Bayesian analysis for the single-harmonic search against an estimate of the sensitivity of the search using both detectors during O4a. The results for each pulsar are compared with the corresponding spin-down limit. The results from this analysis for each pulsar are represented by the blue stars, with their corresponding spin-down limit shown by the gray triangles at the same frequency. The sensitivity curve is shown as a pink line. Some highlighted results for individual pulsars include the Crab pulsar (J0534+2200), which had the lowest spin-down ratio of 0.00783, the Vela pulsar (J0835−4510), J2021+3651, which had the highest odds of coherent signal versus incoherent noise with −3.1, J0537−6910, which had the most constraining amplitude upper limit of 6.38 × 10−27, and J0437−4715, which had the most constraining ellipticity upper limit of 8.8 × 10−9. The distribution of spin-down ratios for these results is shown in Figure 3 with 29 targets that surpass the spin-down limit and the remaining targets, which all have a spin-down ratio below 5.

Figure 2. Refer to the following caption and surrounding text.

Figure 2. Upper limits on h0 for the 45 pulsars in this analysis using the time-domain Bayesian method and considering the single-harmonic emission model. The blue stars show 95% credible upper limits on the amplitudes of h0. Gray triangles represent the spin-down limits for each pulsar (based on the distance measurement stated in Table 2 and assuming the canonical moment of inertia). The pink curve gives an estimate of the expected strain sensitivity of both detectors combined during the course of O4a. The upper limits from the other two pipelines are broadly consistent, as shown in Table 2.

Download figure:

Standard image High-resolution image
Figure 3. Refer to the following caption and surrounding text.

Figure 3. A histogram of the spin-down ratio considering the single-harmonic emission model for 45 pulsars from the Bayesian analysis.

Download figure:

Standard image High-resolution image

In Figure 4, the ellipticity ε95% and mass quadrupole upper limits are plotted against the GW frequency and compared with the corresponding spin-down limits for the ellipticity. The contours of equal characteristic age have been calculated using , which can be derived with the assumption that GW emission alone is driving the spin-down.

Figure 4. Refer to the following caption and surrounding text.

Figure 4. 95% credible upper limits on ellipticity ε95% and mass quadrupole for all 45 pulsars using the Bayesian analysis method and considering the single-harmonic emission model. The upper limits for each pulsar are represented by blue circles while their spin-down limits are shown as gray triangles. Also included are purple contour lines of equal characteristic age , assuming that GW emission alone is causing spin-down. Only the results for pulsars that surpassed their spin-down are physically meaningful. The histogram on the right shows the distribution of the ellipticities obtained from the results (blue filled bars) and in the spin-down limit (gray bars).

Download figure:

Standard image High-resolution image

As shown in Table 2, there is broad agreement among the different pipelines, despite these pipelines being largely independent, and the statistical procedures used to derive the upper limits are different. The data frameworks and preprocessing procedures used by each pipeline account for the differences found in the upper limit results.

5.2. Narrowband Results

In this section, we detail the results obtained with the narrowband pipeline, presented in Section 3.4. The search did not highlight, for any of the 16 considered targets, outliers with a false-alarm probability (FAP) < 10−2 after taking into account the trial factor.

For the narrowband search, we consider only those targets with above the expected sensitivity that is typically worse by a factor of 2 (P. Astone et al. 2014) than that of targeted pipelines due to the trial factor. In this way, we selected 16 pulsars out of which eight have not been analyzed in O3 (R. Abbott et al. 2022b). Our data set includes the two pulsars that glitched, J0537–6910 and J0540–6919. Similar to O3, for these pulsars we split O4a data into two segments that exclude from the analysis the period around [tg-1 d, tg+2 d], with tg the glitch epoch. The segments are then analyzed independently.

The search did not highlight any statistically significant outlier since the measured p-values are well above the threshold set by a FAP of 10−2 corrected for the trial factor. In Table 3, we report for each target the lowest p-value found during the analysis and the threshold.

Table 3. Upper Limits on the Strain Amplitude Set with the Narrowband Search for Each of the Considered Targets

Pulsar NameNtrials Lowest p-valueThreshold p-value
(J2000)×106×10−26 ×10−5×10−11
J0205+6449909.600.221.541.10
J0534+220013834.160.030.770.07
J0537-69102983.521.210.130.33
J0537-6910 pg1193.751.29 0.84
J0540-69191114.232.450.448.95
J0540-6919 pg2487.381.270.020.40
J0835-45102321.530.068.504.18
J1811-1925239.710.720.454.29
J1813-1246245.810.314.074.03
J1813-17491825.170.230.040.55
J1838-06552112.981.272.934.68
J1913+1011154.140.771.406.45
J1935+20251617.841.170.635.94
J1952+3252184.240.411.065.55
J2016+3711736.310.401.371.37
J2021+36511237.110.860.257.73
J2022+3842995.800.540.251.01
J2229+6114776.190.190.621.28

Note. As a reference, we also list the number of templates, Ntrials, in the f − ḟ plane. We also report the lowest p-value found in the analysis with the corresponding threshold set after correcting a FAP of 10−2 for the trial factor. The nomenclature "pg" identifies the postglitch analysis.

Download table as:  ASCIITypeset image

In the absence of any detections, we have calculated upper limits at the 95% confidence level for each of the analyzed targets. Our results are listed in Table 3 and shown in Figure 5, comparing them with the expected sensitivity.

Figure 5. Refer to the following caption and surrounding text.

Figure 5. Expected sensitivity of the narrowband search using the O4a (shaded pink region) data set from the two LIGO detectors considering the single-harmonic emission model. The curve is compared with the spin-down limits (triangles) and the upper limits (stars) averaged over all the 10−4 Hz bands for each source. Upper limits below the spin-down limit are highlighted with orange circles.

Download figure:

Standard image High-resolution image

5.3. Brans–Dicke Theory

Table 4 shows the results for the analyses of 45 pulsars using the statistic to search for dipole radiation predicted by the Brans–Dicke theory. No outliers have been found in the analysis, and we set upper limits on the expected amplitude defined in Equation (7). The upper limits in parentheses show the results using informative priors on the polarization parameters for the pulsars in Table 1.

Table 4. Limits on GW Amplitude from Dipole Radiation in the Brans–Dicke Theory

Pulsar Namefrot Distance FAP
(J2000)(Hz)(Hz s−1)(kpc)  
J0030+0451α205.53−4.23 × 10−160.33a8.77 × 10−270.98
J0058−7218β45.94−6.1 × 10−1159.70b3.02 × 10−260.73
J0117+5914γ9.86−5.69 × 10−131.77c2.03 × 10−231
J0205+6449γ15.22−4.49 × 10−113.20d3.78(2.57) × 10−250.99(0.79)
J0437−4715δ173.69−4.15 × 10−160.16e8.13 × 10−270.98
J0534+2200γ29.95−3.78 × 10−102.00f1.37(0.86) × 10−250.46(0.99)
J0537−6910β62.03−1.99 × 10−1049.70g1.84(1.11) × 10−260.98(0.60)
J0540−6919β19.77−1.87 × 10−1049.70g2.41(1.44) × 10−250.92(0.36)
J0614−3329α317.59−1.76 × 10−150.63h1.78 × 10−260.56
J0737−3039Aα44.05−3.41 × 10−151.10i4.08 × 10−260.83
J0835−4510epsilon11.19−1.57 × 10−110.28i4.12(2.73) × 10−240.80(0.33)
J1231−1411α271.45−5.92 × 10−160.42c3.21 × 10−261
J1412+7922β17.18−9.72 × 10−133.30j5.56 × 10−260.93
J1537+1155α26.38−1.65 × 10−150.93k5.93 × 10−260.97
J1623−2631γ90.29−5.26 × 10−151.85l3.10 × 10−250.77
J1719−1438α172.71−2.22 × 10−160.34c5.44 × 10−271
J1744−1134α245.43−4.34 × 10−160.40e1.60 × 10−260.92
J1745−0952α51.61−2.3 × 10−160.23c4.30 × 10−260.81
J1756−2251γ35.14−1.26 × 10−150.73m5.29 × 10−260.64
J1809−1917γ12.08−3.73 × 10−123.27c3.25 × 10−240.90
J1811−1925β15.46−1.05 × 10−115.00n6.20 × 10−250.98
J1813−1246β20.80−7.6 × 10−122.63o3.60 × 10−250.44
J1813−1749β22.35−6.34 × 10−116.15c1.09 × 10−250.97
J1823−3021Aγ183.82−1.14 × 10−138.02l1.84 × 10−260.58
J1824−2452Aα327.41−1.73 × 10−135.37l2.08 × 10−260.49
J1826−1334γ9.85−7.31 × 10−123.61c3.69 × 10−250.98
J1828−1101γ13.88−2.85 × 10−124.77c3.68 × 10−241
J1831−0952γ14.87−1.84 × 10−123.68c1.13 × 10−240.87
J1833−0827γ11.72−1.26 × 10−124.50i3.79 × 10−240.82
J1837−0604γ10.38−4.84 × 10−124.78c9.77 × 10−240.91
J1838−0655β14.18−9.9 × 10−126.60p2.91 × 10−240.51
J1849−0001β25.96−9.54 × 10−127.00q1.70 × 10−250.34
J1856+0245γ12.36−9.49 × 10−126.32c7.32 × 10−251
J1913+1011γ27.85−2.62 × 10−124.61c8.43 × 10−260.88
J1925+1720γ13.22−1.83 × 10−125.05c4.00 × 10−240.70
J1935+2025γ12.48−9.47 × 10−124.60c1.57 × 10−240.95
J1952+3252γ25.30−3.74 × 10−123.00i1.01 × 10−250.78
J2016+3711ζ19.68−2.81 × 10−116.10r1.91 × 10−250.90
J2021+3651γ9.64−8.89 × 10−121.80s1.32(0.61) × 10−220.39(0.08)
J2022+3842β20.59−3.65 × 10−1110.00t7.44 × 10−240.67
J2043+2740γ10.40−1.37 × 10−131.48c3.93 × 10−240.97
J2214+3000η320.59−1.31 × 10−150.60u5.68 × 10−271
J2222−0137α30.47−4.99 × 10−160.27v5.21 × 10−260.96
J2229+6114γ19.36−2.9 × 10−113.00w2.31(1.61) × 10−250.86(0.25)

Note. For references and other notes, see Table 2. Values in parentheses are those produced using the restricted orientation priors described in Section 3.1.1. The last column shows the FAP for a signal, assuming that the value has a χ2 distribution with two degrees of freedom.

Download table as:  ASCIITypeset image

6. Discussion

In this section, we discuss the results in Table 2. Motivated by the comparable results among the three pipelines, we consider the Bayesian pipeline as a reference. As described in the previous section, we have no evidence of a CW signal in any of the searches we conducted.

We compare the O4a results with previous targeted searches by the LVK Collaboration (B. P. Abbott et al. 2017b, 2019c; R. Abbott et al. 2020, 2022a), considering the first three observing runs. The ratio between the O4a upper limits on h0 and C21 and the corresponding upper limits set in previous searches is shown in Figure 6.

Figure 6. Refer to the following caption and surrounding text.

Figure 6. Blue unfilled stars show the ratio between the O4a h0 upper limits for the analyzed targets (excluding the glitching pulsars), assuming the single-harmonic model divided by the corresponding h0 upper limits in R. Abbott et al. (2022a) for the Bayesian method as a function of the corresponding frequency at twice the rotation frequency (red circles refer instead to the C21 parameter at the rotation frequency assuming the dual-harmonic model). Blue filled stars show the h0 upper limit ratios considering the targets (J0205+6449, J0737−3039A, J1813−1246, J1831−0952, J1837−0604) analyzed using data from the second observing run O2 (B. P. Abbott et al. 2019c) and from the first observing run O1 (blue asterisk for J1826−1334; B. P. Abbott et al. 2017b).

Download figure:

Standard image High-resolution image

Thirty-four pulsars out of the 45 considered targets in Table 2 have been already analyzed in the joint O2+O3 analysis (R. Abbott et al. 2022a). Overall, the corresponding upper limits on the GW amplitude are comparable for the single-harmonic search, with some targets showing better results in O4a and some targets having worse results than those in R. Abbott et al. (2022a); see Figure 6. This is expected since the targeted search sensitivity of the O2+O3 data set is comparable to the O4a one, except at very low frequencies. The targeted search sensitivity can be expressed in terms of minimum detectable amplitude (L. D'Onofrio et al. 2025), hmin. For a multidetector analysis considering n detectors and averaging over the sky position and polarization parameters,

where the factor C ≃ 11 (the exact value depending on the considered pipeline), while Ti and Si are, respectively, the effective observation time and the average power spectral density (PSD) for the ith detector. For the O4a targeted search sensitivity (with an observation time of approximately 1.3 × 107 s for both detectors), see the pink curve in Figure 2. The O4a targeted searches have a sensitivity depth of around 500; here, Sh is the power spectrum taking the harmonic mean of the data over time and over detectors, and the upper limit on the pulsar amplitude (K. Wette 2023; B. Behnke et al. 2015; C. Dreissigacker et al. 2018).

The O4a PSDs for the two LIGO detectors are generally better, by almost a factor of 1.5–2, compared to the corresponding O3 PSD343 depending on the considered frequency band. However, the effective observation time for O4a is reduced by a factor of approximately 1.6, which diminishes the benefit of the improved detector sensitivity in O4a. As a result, we expect the sensitivity of the O4a search to be comparable to that of the combined O2+O3 search. Upper limits on C21, on the other hand, are lower on average than the corresponding O2+O3 results due to a better search sensitivity at frequencies below 20 Hz.

For the remaining targets, five pulsars (J0205+6449, J0737−3039A, J1813−1246, J1831−0952, and J1837−0604) have been analyzed in the O2 search (B. P. Abbott et al. 2019c), and J1826−1334 has been analyzed in the O1 search (B. P. Abbott et al. 2017b). For these pulsars, we have a clear improvement in the upper limits, as shown in Figure 6. The remaining targets (J0058−7218, J1811−1925, J2016+3711, J2021+3651, J2022+3842) have not been analyzed in recent targeted searches, and we surpass the spin-down limit for all these targets.

Many studies have been dedicated to illustrate how a future successful detection of CW will provide a wealth of information about neutron stars (see, e.g., M. Sieniawska & D. I. Jones 2022; N. Lu et al. 2023), and even help to constrain the nuclear equation of state (see, e.g., A. Idrisy et al. 2015; S. Ghosh 2023; S. Ghosh et al. 2023). It is interesting to note that, with improving sensitivity of the searches with each observing run, even nondetection of a CW signal sets more and more stringent upper limits on the ellipticity and possible sources of nonaxisymmetric deformations in rotating neutron stars. This may lead to a better understanding of properties of the crust, internal magnetic fields, and accretion physics (L. Bildsten 1998; A. Melatos & D. J. B. Payne 2005; R. Ciolfi & L. Rezzolla 2013), and even rule out certain scenarios related to r-modes or limit their maximum saturation amplitudes (R. Abbott et al. 2021d).

Theoretical estimates of the maximum mountain sizes that an elastically deformed neutron star can sustain are subject to significant uncertainties, with estimates for the ellipticity epsilon ranging from ∼10−6 for conventional neutron stars, to as large as ∼10−3 for stars with exotic solid phases (see, e.g., G. Ushomirsky et al. 2000; B. J. Owen 2005; B. Haskell et al. 2007; N. K. Johnson-McDaniel & B. J. Owen 2013; F. Gittins & N. Andersson 2021; J. A. Morales & C. J. Horowitz 2022). Comparison with the results given in Figure 4 and Table 2 show that our observationally obtained upper limits overlap with these ranges, confirming we are continuing to push into the regime of astrophysical interest. Estimates for magnetically induced ellipticities are similarly uncertain (see, e.g., B. Haskell et al. 2008; K. Glampedakis et al. 2012; K. Fujisawa et al. 2022). Nevertheless, to give a concrete example, S. Dall'Osso & R. Perna (2017) have suggested that the apparent gradual increase in the angle between the spin axis and magnetic axis of the Crab pulsar provides evidence for a magnetically induced ellipticity epsilon ∼ (3–10) × 10−6 (see, however, S. K. Lander & D. I. Jones 2018). This to be compared with our upper limit of epsilon ≈ 6 × 10−6 for the Crab, again confirming we are probing a regime of astrophysical interest.

Theoretical estimates of the minimum mountain sizes are presented in G. Woan et al. (2018), which provides population-based evidence for millisecond pulsars having a minimum ellipticity of epsilon ≈ 10−9.

We stress that our upper limits are subject to the uncertainties from the detector calibration, as described in Section 4.1, as well as statistical uncertainties that are dependent on the particular analysis method.

The narrowband results in Table 3 and Figure 5 show no evidence of a CW signal for the considered subset of pulsars. No outlier was found for any of the targets. Out of the 16 analyzed targets, 12 searches report an upper limit below the corresponding spin-down limit (see Table 3), ranging from a factor of 1.16 for J2021+3651 to 33 for J0534+2200. As for other methods, the targets analyzed with O3 data (R. Abbott et al. 2022b) report upper limits comparable to those in Table 3. We stress that the narrowband search sensitivity is worse by at least a factor of 2 compared to the targeted search sensitivity as it depends on the number of templates explored for each target (S. Mastrogiovanni et al. 2017).

The search for non-GR polarizations as predicted by the Brans–Dicke theory shows no evidence of dipole radiation. The most constraining upper limit for dipole radiation is obtained for the millisecond pulsar J1719−1438. Together with results from the O3 targeted search analysis (R. Abbott et al. 2022a), the obtained upper limits constitute the first constraints on dipole radiation from pulsar GW observations.

7. Conclusion

In this work, we present a search for CW signals from a set of 45 known pulsars using O4a data from the two LIGO detectors. Pulsars are chosen considering an expected sensitivity for the amplitude below or slightly above the theoretical spin-down limit with a rotation frequency close to or greater than 10 Hz. EM observations were employed to constrain the pulsars' sky positions and rotational parameters covering the O4a data period.

We performed a targeted search utilizing three independent data analysis methods and two different emission models. No evidence of a CW signal was found for any of the targets. The upper limit results show that 29 targets surpass the theoretical spin-down limit. For 11 of the 45 pulsars not analyzed in the last LVK targeted search, we have a notable improvement in detection sensitivity compared to previous searches. For these targets, we surpass or equal the theoretical spin-down limit for the single-harmonic emission model. We also have, on average, an improvement in the upper limits for the low-frequency component of the dual-harmonic search for all analyzed pulsars. For the remaining targets, the O4a upper limits are comparable to the results of the joint O2+O3 analysis described in R. Abbott et al. (2022b), which considered data with lower sensitivity but a longer observation time.

We also conducted a narrowband search for 16 pulsars and a search for non-GR polarization as predicted by the Brans–Dicke theory. No evidence of a CW signal was found in any of these searches.

The analysis of the full O4 data set will improve the sensitivity of targeted/narrowband searches for some of the pulsars analyzed in R. Abbott et al. (2022b) and here, including the Crab and Vela pulsars.

Acknowledgments

This material is based upon work supported by NSF's LIGO Laboratory, which is a major facility fully funded by the National Science Foundation. The authors also gratefully acknowledge the support of the Science and Technology Facilities Council (STFC) of the United Kingdom, the Max-Planck-Society (MPS), and the State of Niedersachsen/Germany for support of the construction of Advanced LIGO and construction and operation of the GEO 600 detector. Additional support for Advanced LIGO was provided by the Australian Research Council. The authors gratefully acknowledge the Italian Istituto Nazionale di Fisica Nucleare (INFN), the French Centre National de la Recherche Scientifique (CNRS), and the Netherlands Organization for Scientific Research (NWO) for the construction and operation of the Virgo detector and the creation and support of the EGO consortium. The authors also gratefully acknowledge research support from these agencies as well as by the Council of Scientific and Industrial Research of India, the Department of Science and Technology, India, the Science and Engineering Research Board (SERB), India, the Ministry of Human Resource Development, India, the Spanish Agencia Estatal de Investigación (AEI), the Spanish Ministerio de Ciencia, Innovación, y Universidades, the European Union NextGenerationEU/PRTR (PRTR-C17.I1), the ICSC—CentroNazionale di Ricerca in High Performance Computing, Big Data, and Quantum Computing, funded by the European Union NextGenerationEU, the Comunitat Autonòma de les Illes Balears through the Conselleria d'Educació i Universitats, the Conselleria d'Innovació, Universitats, Ciència i Societat Digital de la Generalitat Valenciana, and the CERCA Programme Generalitat de Catalunya, Spain, the Polish National Agency for Academic Exchange, the National Science Centre of Poland and the European Union—European Regional Development Fund; the Foundation for Polish Science (FNP), the Polish Ministry of Science and Higher Education, the Swiss National Science Foundation (SNSF), the Russian Science Foundation, the European Commission, the European Social Funds (ESF), the European Regional Development Funds (ERDF), the Royal Society, the Scottish Funding Council, the Scottish Universities Physics Alliance, the Hungarian Scientific Research Fund (OTKA), the French Lyon Institute of Origins (LIO), the Belgian Fonds de la Recherche Scientifique (FRS-FNRS), Actions de Recherche Concertées (ARC) and Fonds Wetenschappelijk Onderzoek—Vlaanderen (FWO), Belgium, the Paris Île-de-France Region, the National Research, Development and Innovation Office of Hungary (NKFIH), the National Research Foundation of Korea, the Natural Sciences and Engineering Research Council of Canada (NSERC), the Canadian Foundation for Innovation (CFI), the Brazilian Ministry of Science, Technology, and Innovations, the International Center for Theoretical Physics South American Institute for Fundamental Research (ICTP-SAIFR), the Research Grants Council of Hong Kong, the National Natural Science Foundation of China (NSFC), the Israel Science Foundation (ISF), the US-Israel Binational Science Fund (BSF), the Leverhulme Trust, the Research Corporation, the National Science and Technology Council (NSTC), Taiwan, the United States Department of Energy, and the Kavli Foundation. The authors gratefully acknowledge the support of the NSF, STFC, INFN, and CNRS for provision of computational resources.

This work was supported by MEXT, the JSPS Leading-edge Research Infrastructure Program, JSPS Grant-in-Aid for Specially Promoted Research 26000005, JSPS Grant-in-Aid for Scientific Research on Innovative Areas 2905: JP17H06358, JP17H06361, and JP17H06364, JSPS Core-to-Core Program A, Advanced Research Networks, JSPS Grants-in-Aid for Scientific Research (S) 17H06133 and 20H05639, JSPS Grant-in-Aid for Transformative Research Areas (A) 20A203: JP20H05854, the joint research program of the Institute for Cosmic Ray Research, the University of Tokyo, the National Research Foundation (NRF), the Computing Infrastructure Project of Global Science experimental Data hub Center (GSDC) at KISTI, the Korea Astronomy and Space Science Institute (KASI), the Ministry of Science and ICT (MSIT) in Korea, Academia Sinica (AS), the AS Grid Center (ASGC) and the National Science and Technology Council (NSTC) in Taiwan under grants including the Rising Star Program and Science Vanguard Research Program, the Advanced Technology Center (ATC) of NAOJ, and the Mechanical Engineering Center of KEK.

C.M.E. acknowledges support from ANID/FONDECYT, grant 1211964. S. Guillot acknowledges the support of the CNES. W.C.G.H. acknowledges support through grant 80NSSC22K1305 from NASA. This work is supported by NASA through the NICER mission and the Astrophysics Explorers Program and uses data and software provided by the High Energy Astrophysics Science Archive Research Center (HEASARC), which is a service of the Astrophysics Science Division at NASA/GSFC and High Energy Astrophysics Division of the Smithsonian Astrophysical Observatory.

The Nançay radio Observatory is operated by the Paris Observatory, associated with the French Centre National de la Recherche Scientifique (CNRS), and partially supported by the Region Centre in France. We acknowledge financial support from "Programme National de Cosmologie and Galaxies" (PNCG), and "Programme National Hautes Energies" (PNHE) funded by CNRS/INSU-IN2P3-INP, CEA, and CNES, France.

A.B.P. is a Banting Fellow, a McGill Space Institute (MSI) Fellow, and a Fonds de Recherche du Quebec—Nature et Technologies (FRQNT) postdoctoral fellow.

E. Fonseca is supported by the NSF grant AST-2407399.

The activities at the Instituto Argentino de Radioastronomía (IAR) are supported by the national agency CONICET, the Province of Buenos Aires agency CIC, and the National University of La Plata (UNLP).

Software: The 5-vector method is based on the BSD framework (O. J. Piccinni et al. 2018) while the narrowband method makes use of the SFDB framework (P. Astone et al. 2005), both of which are based on the Virgo Rome Snag software. The Bayesian analysis uses CWInPy (M. Pitkin 2022), which uses dynesty (J. Skilling 2004, 2006) within bilby (G. Ashton et al. 2019). Plots are produced using matplotlib (J. D. Hunter 2007). Many pulsar ephemerides are produced with Tempo (D. Nice et al. 2015) and Tempo2 (G. B. Hobbs et al. 2006).

Footnotes

  • 342 

    We considered L1:GDS–CALIB_STRAIN_CLEAN_AR and H1:GDS–CALIB_STRAIN_CLEAN_AR frame channels with CAT1 vetoes for L1 and H1, respectively.

  • 343 

    For this estimation, we only consider the O3 run for simplicity since it dominates the combined data sets.

10.3847/1538-4357/adb3a0
undefined