# Improving pulsar polarization and timing measurements with the Nançay Radio Telescope

Lucas Guillemot, Ismaël Cognard, Willem van Straten, Gilles Theureau, Eric Gérard

► To cite this version:

Lucas Guillemot, Ismaël Cognard, Willem van Straten, Gilles Theureau, Eric Gérard. Improving pulsar polarization and timing measurements with the Nançay Radio Telescope. Astronomy & Astrophysics - A&A, 2023, 678, pp.A79. 10.1051/0004-6361/202347018. hal-04189901

HAL Id: hal-04189901

https://hal.science/hal-04189901v1

Submitted on 9 Nov 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Improving pulsar polarization and timing measurements with the Nançay Radio Telescope

L. Guillemot<fnref n="1" />,<fnref n="2" />,✉, I. Cognard<fnref n="1" />,<fnref n="2" />✉, W. van Straten<fnref n="3" />,<fnref n="4" />✉, G. Theureau<fnref n="1" />,<fnref n="2" />,<fnref n="5" />✉, and E. Gérard<fnref n="6" />

Received 26 May 2023 / Accepted 2 August 2023

## ABSTRACT

**Context.** Accurate polarimetric calibration of the radio pulse profiles from pulsars is crucial for studying their radiation properties at these wavelengths. Additionally, inaccurate calibration can distort recorded pulse profiles, introducing noise in time of arrival (TOA) data and thus degrading pulsar timing analyses. One method for determining the full polarimetric response of a given telescope is to conduct observations of bright polarized pulsars over wide ranges of parallactic angles, to sample different orientations of their polarization angle and in turn determine the cross-couplings between polarization feeds.

**Aims.** The Nançay decimetric Radio Telescope (NRT) is a 94 m equivalent meridian telescope, capable of tracking a given pulsar for approximately 1 h around transit. The NRT therefore cannot sample wide ranges of parallactic angles when observing a given pulsar, so until late 2019 the polarimetric calibration of 1.4 GHz pulsar observations with the NRT was rudimentary. We therefore aimed to develop a method for improving the calibration of NRT observations, overcoming the above-mentioned limitation. Ultimately, our goal was to improve the quality of NRT pulsar timing, with better calibrated pulsar pulse profiles.

**Methods.** In November 2019, we began conducting regular observations of the bright and highly linearly polarized pulsar PSR J0742−2822, in a special observing mode in which the feed horn rotates by ~180° over the course of the 1 h observation, mimicking wide parallactic angle variations and in principle enabling us to determine the polarimetric response of the NRT at 1.4 GHz. In addition, we assessed the quality of the NRT timing of a selection of millisecond pulsars (MSPs), namely, J1730−2304, J1744−1134, and J1857+0953, with conventional TOAs extracted from total intensity pulse profiles, and TOAs extracted with the Matrix Template Matching (MTM) technique, designed to compensate for putative polarimetric calibration errors.

**Results.** From the analysis of the rotating horn observations of PSR J0742−2822 we could determine the cross-couplings between the polarization feeds and also constrain the Stokes parameters of the noise diode signal, which prior to this work was erroneously assumed to be ideal and was used as the only reference source for the calibration of pulsar observations. The improved polarimetric response of the NRT as determined from these observations was applied to observations of a selection of MSPs with published polarimetric properties. We find that the new polarimetric profiles and polarization position angles are consistent with previous findings, unlike NRT polarimetric results obtained with the previously used method of calibration. The analysis of the timing data shows that the new calibration method improves the quality of the timing, and the MTM method proves very effective at reducing noise from imperfect calibration. For pulsars with sufficient degrees of polarization, the MTM method appears to be the preferred method of extracting TOAs from NRT observations.

**Key words.** polarization – pulsars: general – pulsars: individual: J0742−2822 – pulsars: individual: J1730−2304 – pulsars: individual: J1744−1134 – pulsars: individual: J1857+0943

## 1. Introduction

Pulsars are rapidly rotating, highly magnetized neutron stars that emit beams of electromagnetic radiation. These beams are swept across the sky as they rotate, so that the emission appears to be “pulsed” to distant observers whose lines of sight are crossed by the beams. Pulsars can be used as tools to probe a variety of astrophysical topics, via the “pulsar timing” technique. Pulsar timing consists of measuring pulse times of arrival (TOAs), which correspond to epochs at which a fiducial phase of the

pulsar’s periodic signal is detected at a telescope. The TOAs of a given pulsar encode much information about its spin properties, its orbital properties if the pulsar belongs to a multiple system, information about the intervening interstellar medium, and other characteristics of the pulsar and its environment (see, for instance, Taylor 1992; Stairs 2003; Lorimer & Kramer 2004, for reviews of pulsar timing and its applications). These properties can be determined by minimizing the differences between the TOAs and the predictions from a model for the pulsar encoding the above-listed information, the so-called “timing residuals”.

<footnote n="1">¹ Laboratoire de Physique et Chimie de l'Environnement et de l'Espace, Université d'Orléans/CNRS, 45071 Orléans Cedex 02, France
e-mail: lucas.guillemot@cnrs-orleans.fr</footnote>

<footnote n="2">² Observatoire Radioastronomique de Nançay, Observatoire de Paris, Université PSL, Université d'Orléans, CNRS, 18330 Nançay, France</footnote>

<footnote n="3">³ Institute for Radio Astronomy & Space Research, Auckland University of Technology, Private Bag 92006, Auckland 1142, New Zealand</footnote>

<footnote n="4">⁴ Manly Astrophysics, 15/41-42 East Esplanade, Manly, NSW 2095, Australia</footnote>

<footnote n="5">⁵ LUTH, Observatoire de Paris, Université PSL, Université Paris-Cité, CNRS, 92195 Meudon, France</footnote>

<footnote n="6">⁶ GEPI, Observatoire de Paris, Université PSL, Université Paris-Cité, CNRS, 92195 Meudon, France</footnote>
For some pulsars with high rotational stability, timing residuals with standard deviations on the order of 100 ns can be achieved (see e.g., Liu et al. 2020).

The Nançay Radio Telescope (NRT) is a meridian transit-type telescope of the Kraus/Ohio State design (see e.g., Kraus 1960, 1966). It consists of a tiltable flat primary mirror (with dimensions 200 × 40 m), a fixed spherical secondary mirror (300 × 35 m) separated from the primary mirror by 460 m, and a receiver carriage that follows the focal point during observations by moving on a 80 m train track. The NRT is equivalent to a parabolic dish with a diameter of about 94 m, and can track objects with declinations above ~ -39° for approximately 1 h around transit. It is equipped with two cryogenically cooled receivers that respectively cover the [1.1; 1.8] GHz and [1.7; 3.5] GHz frequency ranges, allowing for measurement of the four Stokes parameters over those ranges. Since August 2011, pulsar observations with the NRT have been conducted using the NUPPI backend, a version of the Green Bank Ultimate Pulsar Processing Instrument used at the Green Bank Telescope (DuPlain et al. 2008), and designed for the NRT. The NUPPI instrument (for further details, see e.g., Desvignes et al. 2011; Cognard et al. 2013) uses Graphics Processing Units (GPUs) to coherently dedisperse and fold the dual linear polarization signals from the receiver in real time, over a total bandwidth of 512 MHz. The 512 MHz bandwidth is split into 128 channels of 4 MHz each, and in practice the bulk of pulsar observations with the NUPPI backend are conducted at a central frequency of 1484 MHz, to exploit as much of the frequency interval covered by the 1.1-1.8 GHz receiver as possible. A smaller fraction of the observations is conducted with the high-frequency receiver at a central frequency of 2539 MHz, or at other central frequencies, depending on the pulsar or science objective. Brief descriptions of the instruments used for conducting pulsar timing observations with the NRT prior to NUPPI can be found in Desvignes et al. (2016).

Pulsar timing data recorded with the NUPPI backend feature prominently in several recent high-precision pulsar timing studies. Examples include the timing of PSR J0737-3039A in the so-called "double pulsar" system, enabling highly precise tests of General Relativity (GR) in the strong-field regime (Kramer et al. 2021), tests of the strong equivalence principle of GR via the timing of PSR J0337+1715 in a stellar triple system (Voisin et al. 2020), the measurement of the mass of PSR J2222-0137 in a massive neutron star - white dwarf binary system (Cognard et al. 2017; Guo et al. 2021), the timing of the eclipsing black widow pulsar PSR J2055+3829 discovered at Nançay (Guillemot et al. 2019), or searches for low-frequency gravitational waves using Pulsar Timing Arrays (PTAs; see e.g., Chen et al. 2021). All of these studies analyzed TOAs that had been extracted from NUPPI data by cross-correlating a standard pulse profile for the pulsar of interest with individual observations, partially averaged in frequency and/or time. A comprehensive review of the TOA extraction process and of sources of measurement uncertainties can be found in e.g., Verbiest & Shaifullah (2018), and a recent comparison of TOA creation practices using data from the NRT and other telescopes can be found in Wang et al. (2022). One crucial assumption made when extracting TOAs from timing data by cross-correlating a standard (or template) pulse profile with observations of a given pulsar is that the observations resemble the standard profile; that is, they are scaled versions of the template profile on top of varying white noise levels. Distortions of the observed pulse profile caused by time-varying instrumental artifacts can cause individual observations to differ from the template profile, affecting the accuracy of the TOAs

and introducing systematic timing errors in the TOA dataset. As a consequence, the individual pulsar timing observations need to be carefully calibrated in polarization. Moreover, accurate polarization properties of a given pulsar can provide a useful insight into the physical processes at work in its magnetosphere (see e.g., Philippov & Kramer 2022).

As will be presented in further detail in Sect. 2, until recently the procedure followed for calibrating pulsar observations with the NRT was limited. With seemingly correct pulsar timing results with the NRT (as e.g., determined in the context of multi-telescope pulsar timing studies), and a relative paucity of pulsar polarization studies with this telescope (see for instance Theureau et al. 2011; Dyks et al. 2016; Desvignes et al. 2022, for examples of NRT pulsar polarization measurements), the validity of the former polarimetric calibration procedure had not been studied in detail prior to this work. However, detailed comparisons of NRT polarization profiles with reference results obtained at other telescopes, and significant modifications in the observed polarization properties of pulsars concurrent with instrumental changes in 2019, prompted us to revisit the procedure for calibrating pulsar observations with the NRT, and the procedure for extracting TOAs from NRT data. In Sect. 2 we present efforts to improve the polarimetric calibration procedure. In Sect. 3, we present the 1.4 GHz NUPPI data archive, and assess the impact of the improved calibration procedure on the quality of the timing, investigating a novel use of the Matrix Template Matching technique (MTM; see van Straten 2006) with NRT data. This technique uses the timing information in all four Stokes components of the pulsar signal while compensating for residual polarimetric calibration errors when extracting TOAs from observations. Finally, we present a summary of our findings and some prospects from this study in Sect. 4.

## 2. Observations and polarimetric calibration

### 2.1. First order polarimetric calibration of NUPPI observations

In August 2011, the NUPPI backend became the primary pulsar instrument in operation at the NRT, replacing the Berkeley-Orléans-Nançay (BON) pulsar backend. As mentioned in the introduction, the majority of NRT pulsar timing observations are conducted at a central frequency of 1.484 GHz with the low-frequency receiver. In this work we focus on these observations, leaving observations with previous backends or with the high-frequency receiver for future work.

At the beginning of each pulsar observation with the NRT, a 3.33 Hz noise diode injects a polarized reference signal into the receiver feed horn for ten seconds. An example of a noise diode observation is displayed in Fig. 1. As can be seen from the figure, the noise diode used at the NRT has three different states: an “OFF” state, in which the diode is inactive, and two “ON” states. In the first “ON” state (i.e., the “ON1” state), the noise diode theoretically injects a 100% linearly polarized signal that illuminates both feeds equally, which corresponds to 100% Stokes *U* polarization (note that in Fig. 1 the noise diode observation is not calibrated, so the “ON1” state appears to contain little emission along Stokes *U*, and strong circularly polarized emission). Under the latter assumption, the signal from the noise diode can thus be used to perform a first-order calibration of the complex gains of the two polarization feeds, as implemented in *e.g.*, the SINGLEAXIS calibration method of the PSRCHIVE software library (Hotan et al. 2004) used for manipulating the pulsar data presented throughout this article. The “ON2” state,
on the other hand, theoretically consists of a 100% linearly polarized signal along Stokes −U and is not used by PSRCHIVE in the calibration procedure.

In ten seconds, the noise diode is typically detected in each 4 MHz channel with a signal-to-noise ratio (S/N) of 2000. Throughout this article, S/N values are calculated using Eq. (7.1) of Lorimer & Kramer (2004),

$$S/N = \frac{1}{\sigma_{\text{Off}} \sqrt{W_{\text{eq}}}} \sum_{i=1}^{N_{\text{bins}}} (p_i - p_{\text{Off}}), \qquad (1)$$

where $p_i$ denotes the amplitude of the i-th bin, $\sigma_{\text{Off}}$ and $p_{\text{Off}}$ represent the off-pulse standard deviation and mean amplitude, and $W_{\text{eq}}$ is the equivalent width of a top-hat pulse with the same area and peak height as the profile, expressed in number of bins. With an S/N of 2000, and given the characteristics of the NRT (see Table 6 of Wu et al. 2018), we estimate, using Eq. (7.10) of Lorimer & Kramer (2004), that the noise diode increases the receiver temperature by ~20% of the system temperature. This significant temperature increase enables accurate determination of the calibration parameters.

A comprehensive description of the formalism employed by PSRCHIVE for the calibration of pulsar observations can be found in van Straten (2004). In a nutshell, the SINGLEAXIS method of PSRCHIVE assumes that the polarization feeds are perfectly orthogonal and that the reference signal from the noise diode is 100% linearly polarized and illuminates the two feeds equally and in phase, as should be the case under the assumption made above. Under these hypotheses, the 2 × 2 complex-valued Jones matrix $J_{\text{SA}}$ describing the response of the instrument can be written as:

$$J_{\text{SA}} = G \mathbf{B}_{\hat{q}}(\gamma) \mathbf{R}_{\hat{q}}(\varphi), \qquad (2)$$

where G is the absolute gain, $\gamma$ denotes the differential gain between the feeds, $\varphi$ is the differential phase, and $\mathbf{B}_{\hat{q}}$ and $\mathbf{R}_{\hat{q}}$ denote Lorentz boost and rotation transformations (Britton 2000). Under the assumptions made here, the Jones matrix has the following simple form:

$$J_{\text{SA}} = G \begin{pmatrix} e^{\gamma+i\varphi} & 0 \\ 0 & e^{-(\gamma+i\varphi)} \end{pmatrix}. \qquad (3)$$

The measured coherency matrix, $\rho'_{\text{ref}}$, of the reference signal from the noise diode is then given by the polarization measurement equation:

$$\rho'_{\text{ref}} = J_{\text{SA}} \rho_{\text{ref}} J_{\text{SA}}^{\dagger}, \qquad (4)$$

where $\rho_{\text{ref}}$ is the coherency matrix of the input ideal reference signal, in our case:

$$\rho_{\text{ref}} = \frac{C_0}{2} \begin{pmatrix} 1 & 1 \\ 1 & 1 \end{pmatrix}, \qquad (5)$$

with $C_0$ the flux density of the noise diode at the considered frequency. The SINGLEAXIS method uses Eq. (4) and measurements of $\rho'_{\text{ref}}$ at different frequencies to determine G, $\gamma$, and $\varphi$ across the bandwidth recorded by NUPPI. Figure 2 shows the measured G, $\gamma$, and $\varphi$ parameters as a function of frequency, as obtained from the analysis of the noise diode observation plotted in Fig. 1.

The SINGLEAXIS Jones matrix characterizing the instrumental response determined from the analysis of the noise diode

<figure><img src="image_2.png" /><figcaption>Fig. 1. Example of a noise diode observation conducted prior to a pulsar observation with the NRT. The black line represents the total intensity, I, of the signal as a function of pulse phase (the noise diode signal is periodic, with a frequency of 3.33 Hz), and the red, magenta, and blue lines correspond respectively to the Q, U, and V Stokes parameters. Note that the I, Q, U, and V Stokes parameters are uncalibrated in this figure. The data were extracted from a 10-s observation of the noise diode with the NUPPI backend, conducted on MJD 59368 prior to an observation of PSR J0742−2822 at the central frequency of 1.484 GHz. The data shown in this figure were taken from the 46th 4-MHz channel recorded by NUPPI, corresponding to the frequency of 1.410 GHz.</figcaption></figure>

<figure><img src="image_3.png" /><figcaption>Fig. 2. First order calibration parameters as determined from the SINGLEAXIS analysis of the noise diode observation shown in Fig. 1, as a function of frequency. Top: Differential phase, $\varphi$. Middle: Differential gain, $\gamma$. Bottom: Absolute gain, G, expressed in units of the square root of the flux density of the noise diode in the considered frequency channel, $c_0 = \sqrt{C_0}$.</figcaption></figure>
observation can then be used to perform the polarimetric calibration of the pulsar observation conducted subsequently. The Jones matrix to be used for this transformation is $J_{SA} R_{\theta}(\Phi)$, where $J_{SA}$ is the Jones matrix describing the instrumental response and $\Phi$ is the time-varying parallactic angle, given in the case of the NRT by:

$$ \sin(\Phi) = \sin(\delta) \times \sin(\text{HA}), \qquad (6) $$

where $\delta$ and HA are the declination and the (time-varying) hour angle of the source.

In the left-hand panels of Figs. 3–6, we show polarimetric profiles at 1.4 GHz for a selection of MSPs observed with NUPPI and calibrated using the SINGLEAXIS method of PSRCHIVE. For each pulsar we selected a high S/N observation, which was then calibrated using the first order calibration parameters (gains, differential gains, and differential phases) determined from the analysis of the noise diode observation conducted prior to it. The data were cleaned of radio-frequency interference (RFI) using the SURGICAL method of the COASTGUARD pulsar analysis library (Lazarus et al. 2016). The original number of 2048 profile bins was reduced by a factor of up to eight to improve readability, and the data were corrected for Faraday rotation across the band caused by magnetic fields along the line-of-sight to the pulsars, using the Rotation Measure (RM) values from Dai et al. (2015). The average baseline levels for the I, Q, U, and V Stokes parameters were set to zero. Stokes parameters in the figures are in accordance with the conventions described in van Straten et al. (2010). The average parallactic angles of all the considered observations were within ±1°, except for that of J1730–2304, for which the average parallactic angle was −1.8°. Additional details of the selected observations can be found in Table 1. The right-hand panels of Figs. 3–6 show “reference” polarimetric profiles from Dai et al. (2015), obtained by analyzing 20 cm data from the Parkes telescope in Australia. These profiles were obtained by summing up multiple observations conducted at a central frequency of 1.369 GHz, and over a bandwidth of 256 MHz, i.e., from 1.241 to 1.497 GHz. We discarded NUPPI frequency channels outside of this frequency range to prevent potential pulse profile evolution with frequency from altering the integrated NUPPI profiles, and therefore facilitate a comparison with the Dai et al. (2015) results. Unlike NUPPI data, the 20 cm Parkes data analyzed in Dai et al. (2015) were not coherently dedispersed. As a result, the profiles are affected by dispersion smearing, indicated by markers at the top right. Details on dispersion smearing and the amount of smearing for each of the considered pulsars are given in Dai et al. (2015). In a few cases, such as those of PSRs J0613–0200, J1600–3053, and J1909–3744, narrow profile structures visible in the NUPPI profiles can be smeared out in the Parkes profiles, due to this effect. Finally, in Figs. 7 and 8 we show plots of the position angles (PAs), $\Psi$, of the linear polarization as a function of rotational phase for the same pulsars, calculated as:

$$ \Psi = \frac{1}{2} \arctan \left( \frac{U}{Q} \right). \qquad (7) $$

We note that the PAs shown in Figs. 7 and 8 were obtained by correcting Faraday rotation to infinite frequency<fnref n="1" />, hence the fact that the “reference” PAs are different from those shown in Dai et al. (2015).

It is clear, from Figs. 3–8, that the polarimetric profiles and linear polarization angles derived from NUPPI data calibrated

using the SINGLEAXIS method are highly inconsistent with the reference results from Dai et al. (2015). The Parkes Multibeam receiver, used to collect the data presented in the latter article, is also equipped with a noise diode injecting a 100% linearly polarized signal, which is used to determine differential gains and phases prior to the pulsar observations. However, in addition to the noise diode observations, observations of the nearby, highly polarized MSP J0437–4715 (Johnston et al. 1993) are regularly conducted over wide ranges of parallactic angles in order to determine cross-couplings and polarization leakages between the feeds, as described in van Straten (2004). This model, developed further in Sect. 2.2, provides a more accurate estimate of the instrumental response. Polarimetric profiles and PA measurements from NUPPI data calibrated using the SINGLEAXIS method are therefore incorrect. This polarimetric calibration procedure was used until an improved calibration scheme (see next section) was developed in late 2019.

## 2.2. Improving the polarimetric calibration of NUPPI observations

As demonstrated in Sect. 2.1, a simple polarimetric calibration of NUPPI data that assumes perfectly orthogonally polarized feeds and an ideal, 100% linearly polarized noise diode signal illuminating the two feeds equally and in phase does not produce satisfactory results. Cross-coupling between the two feeds, and a potentially non-ideal noise diode signal, can alter polarimetric results. In order to determine these second-order calibration corrections, a number of studies have used observations of bright, broadband-emitting polarized pulsars over wide ranges of parallactic angles (see e.g., Stinebring et al. 1984; Xilouris 1991; Johnston 2002). Doing so, different orientations of the angle of polarization of the pulsar are sampled during the observation, and cross-coupling of the polarization feeds can in turn be determined. The NRT, however, can only observe sources for approximately 1 h around transit. From Eq. (6), it can be seen that in the most favorable situation, corresponding to $\delta = 90^{\circ}$, the parallactic angle varies by only 15° in an hour of observation, insufficient for appropriately sampling pulsar polarization angles. However, the discrepancies between NUPPI polarimetric profiles obtained with a simple SINGLEAXIS calibration with published results (such as the comparisons made in Sect. 2.1) prompted us to explore strategies for circumventing the above-mentioned limitation of the telescope. The solution found in late 2019 is to conduct observations of bright, polarized pulsars in a special observational mode, in which the feed horn, which is held fixed during normal observations, is made to rotate by 180° during the ~1-h observation, i.e., at a rate of ~3° per minute, close to the maximum horn rotation rate allowed by the instrument. Doing so effectively mimics a variation of the parallactic angle, thereby enabling the NRT to sample much wider ranges of pulsar polarization angles during the observation.

The framework for simultaneously determining the full polarimetric response of the receiver and the actual Stokes parameters of the reference signal from the analysis of such observations over wide parallactic angle ranges, or, as in our case, over wide horn orientation ranges, is presented in van Straten (2004) and is implemented in the RECEPTION calibration model of PSRCHIVE. In this model, the Jones matrix representing the instrumental response is given by

$$ J_R = G B_{\hat{q}}(\gamma) R_{\hat{q}}(\varphi) C. \qquad (8) $$

<footnote n="1">¹ With the “--aux_rm” option of the PSRCHIVE tool “pam”.</footnote>
<figure><img src="image_4.png" /><figcaption>**Fig. 3.** Polarimetric profiles for PSRs J0613−0200, J1022+1001, and J1024−0719. In each panel, the black line represents the total intensity (Stokes parameter *I*), the red and magenta lines correspond to the Stokes parameters *Q* and *U* describing the linear polarization, and the blue line is the circular polarization (Stokes parameter *V*). For each pulsar, the left-hand panel shows NUPPI pulse profiles calibrated with the SINGLEAXIS method as described in Sect. 2.1. The middle panel shows NUPPI polarimetric profiles calibrated with the improved scheme presented in Sect. 2.2. Reference results from Dai et al. (2015) are displayed in the right-hand panel. In the latter panel, the marker at the top right indicates the dispersion smearing resulting from incoherent dedispersion of the Parkes data. The displayed pulse phase ranges were restricted to intervals with significant emission features, to facilitate comparisons. The profiles were normalized to the maximum value of the total intensity. See Table 1 for details on the NUPPI observations presented in this figure.</figcaption></figure>
<figure><img src="image_5.png" /><figcaption>Fig. 4. Same as Fig. 3, for PSRs J1600–3053, J1643–1224, and J1713+0747.</figcaption></figure>

In this equation, the term C, which represents the response of a receiver with non-ideal feeds, is given by

$$C = \delta_0 S(\theta_0, \epsilon_0) + \delta_1 S(\theta_1, \epsilon_1), \quad (9)$$

where $\delta_i$ is the $2 \times 2$ selection matrix:

$$\delta_i = \begin{pmatrix} \delta_{0i} & 0 \\ 0 & \delta_{1i} \end{pmatrix}, \quad (10)$$

with $\delta_{ij}$ the Kronecker delta. In Eq. (9), $\epsilon$ and $\theta$ represent the ellipticities and orientations of the two feeds, and $S(\theta, \epsilon) =$

$R_{\tilde{u}}(\epsilon) R_{\tilde{v}}(\theta)$. The model describing the noise diode observation is then $\rho'_{\text{ref}} = J_R \rho_{\text{ref}} J_R^{\dagger}$, where $\rho_{\text{ref}}$ no longer corresponds to the ideal reference signal given in Eq. (5), and the model representing the pulsar observation is

$$\rho'_{\text{psr}}(\Phi') = J_R R_{\tilde{v}}(\Phi') \rho_{\text{psr}} R_{\tilde{v}}^{\dagger}(\Phi') J_R^{\dagger}, \quad (11)$$

with $\Phi' = \Phi + \alpha$ where $\alpha$ is the orientation of the horn at a given time, and $\rho_{\text{psr}}$ is the coherency matrix that is intrinsic to the observed pulsar signal. Given $\rho'_{\text{ref}}$ and $\rho'_{\text{psr}}$ for different
<figure><img src="image_6.png" /><figcaption>Fig. 5. Same as Fig. 3, for PSRs J1730−2304, J1744−1134, and J1857+0943.</figcaption></figure>

values of $\Phi'$, and assuming $\rho'_{psr}$ is constant in all time samples, the Reception model can solve for the actual Stokes parameters of the reference noise diode signal, and for the values of $G$, $\gamma$, and $\phi$ and of the ellipticities, $\epsilon$, and orientations, $\theta$, of the two feeds. In practice, multiple on-pulse phase bins from a given pulsar can be included as input source polarizations, to increase the number of constraints for the fit. Finally, as described in Ord et al. (2004), the polarimetric response determined from the above-described analysis at a reference epoch can be used to calibrate

regular pulsar observations using the following Jones matrix:

$$J = G' B_{\hat{q}}(\gamma') R_{\hat{q}}(\phi') J_0. \qquad (12)$$

In the above expression, $J_0$ is the Jones matrix representing the instrumental response at the reference epoch, and $G'$, $\gamma'$, and $\phi'$ denote the changes in absolute gain, differential gain, and differential phase since that reference epoch, as determined from a SINGLEAXIS analysis of a new observation of the noise
<figure><img src="image_7.png" /><figcaption>**Fig. 6.** Same as Fig. 3, for PSRs J1909−3744, J2124−3358, and J2145−0750.</figcaption></figure>

diode, using the actual Stokes parameters of the diode determined from the Reception analysis. This assumes that the polarization properties of the noise diode have not changed since the reference epoch, and similarly, the model described in Eq. (12) assumes that the ellipticities and orientations also have not varied.

We conducted the first observation of a bright polarized pulsar with feed horn rotation on MJD 58816 (November 29,

2019). The observed pulsar was PSR J0742−2822 (B0740−28; Bonsignori-Facondi et al. 1973), a relatively bright (flux density at 1.4 GHz of 26 ± 2 mJy; see Jankowski et al. 2018), slow-period (P ~ 0.167 s) pulsar with a very strong linear polarization component: Johnston & Kerr (2018) measured a percentage of linear polarization (L/I) of ~90% for this pulsar, making it an ideal target for our calibration observations. Since MJD 58816, PSR J0742−2822 has been observed with NUPPI in calibration
<table><caption>**Table 1.** Properties of the 1.4 GHz NUPPI observations presented in Figs. 3–6.</caption><thead><tr><th>Pulsar</th><th>Discovery article</th><th>P (ms)</th><th>DM (pc cm<sup>-3</sup>)</th><th>RM (rad m<sup>-2</sup>)</th><th>Epoch (MJD)</th><th>Duration (s)</th><th>N<sub>bins</sub></th><th>Orig. S/N</th><th>New S/N</th></tr></thead><tbody><tr><td>J0613-0200</td><td>Lorimer et al. (1995)</td><td>3.06</td><td>38.77</td><td>9.7</td><td>59346</td><td>3192.3</td><td>256</td><td>279.1</td><td>316.6 (+13.4%)</td></tr><tr><td>J1022+1001</td><td>Camilo et al. (1996)</td><td>16.45</td><td>10.27</td><td>-0.6</td><td>59339</td><td>3931.7</td><td>2048</td><td>6898.1</td><td>6384.0 (-7.5%)</td></tr><tr><td>J1024-0719</td><td>Bailes et al. (1997)</td><td>5.16</td><td>6.49</td><td>-8.2</td><td>59323</td><td>3901.1</td><td>512</td><td>349.1</td><td>454.6 (+30.2%)</td></tr><tr><td>J1600-3053</td><td>Jacoby et al. (2007)</td><td>3.60</td><td>52.33</td><td>-15.5</td><td>59338</td><td>3467.7</td><td>1024</td><td>794.8</td><td>841.5 (+5.9%)</td></tr><tr><td>J1643-1224</td><td>Lorimer et al. (1995)</td><td>4.62</td><td>62.40</td><td>-308.1</td><td>59400</td><td>1897.0</td><td>512</td><td>561.1</td><td>579.3 (+3.2%)</td></tr><tr><td>J1713+0747</td><td>Foster et al. (1993)</td><td>4.57</td><td>15.99</td><td>8.4</td><td>59384</td><td>2830.2</td><td>2048</td><td>4736.7</td><td>4966.0 (+4.8%)</td></tr><tr><td>J1730-2304</td><td>Lorimer et al. (1995)</td><td>8.12</td><td>9.63</td><td>-7.2</td><td>59368</td><td>2080.6</td><td>512</td><td>484.4</td><td>509.3 (+5.1%)</td></tr><tr><td>J1744-1134</td><td>Bailes et al. (1997)</td><td>4.07</td><td>3.14</td><td>-1.6</td><td>59347</td><td>1902.1</td><td>1024</td><td>1236.4</td><td>1510.6 (+22.2%)</td></tr><tr><td>J1857+0943</td><td>Segelstein et al. (1986)</td><td>5.36</td><td>13.30</td><td>16.4</td><td>59401</td><td>3039.3</td><td>512</td><td>1025.0</td><td>1018.3 (-0.7%)</td></tr><tr><td>J1909-3744</td><td>Jacoby et al. (2003)</td><td>2.95</td><td>10.39</td><td>-6.6</td><td>59359</td><td>3707.4</td><td>2048</td><td>1242.1</td><td>1418.1 (+14.2%)</td></tr><tr><td>J2124-3358</td><td>Bailes et al. (1997)</td><td>4.93</td><td>4.60</td><td>-5.0</td><td>59365</td><td>2738.4</td><td>512</td><td>439.7</td><td>425.8 (-3.2%)</td></tr><tr><td>J2145-0750</td><td>Bailes et al. (1994)</td><td>16.05</td><td>9.00</td><td>-1.3</td><td>59338</td><td>3013.8</td><td>2048</td><td>3004.3</td><td>3133.6 (+4.3%)</td></tr></tbody></table>

**Notes.** The first five columns give the pulsar names and discovery articles, their rotational periods (P), and Dispersion Measure (DM) values, and the Rotation Measure (RM) values from Dai et al. (2015) used for correcting the NUPPI data for Faraday rotation. The following columns give the observation epochs and durations, the numbers of pulse phase bins of the NUPPI profiles displayed in Figs. 3–6, and the signal-to-noise (S/N) ratios of the profiles, when calibrating the data using the SINGLEAXIS method of PSRCHIVE ("Orig. S/N") and with the improved calibration method presented in Sect. 2.2 ("New S/N"). See Sect. 2.2 and Eq. (1) for details on the calculation of the S/N values.

mode on a regular basis, with one to two observations per month on average. We here present the results of the analysis of a calibration mode observation of PSR J0742-2822 conducted on MJD 59368 (June 3, 2021), close in time to the regular mode NUPPI observations presented in Figs. 3–6. We integrated the ~1-h observation in time to form one sub-integration per minute of observation. The original frequency resolution of 128 channels of 4 MHz each was kept. The original number of phase bins per rotation of 2048 was reduced by a factor of four in order to increase the S/N of the pulsar in the individual on-pulse phase bins and time sub-integrations. We then used the PSRCHIVE tool "pcm" to carry out the Reception model analysis of the observation of PSR J0742-2822 and of the noise diode conducted prior to the pulsar observation, selecting 32 on-pulse phase bins as sources of input pulsar polarization data for the analysis. As explained in detail in Appendix B of van Straten (2004), solutions to Eq. (11) are degenerate under commutation, and two assumptions need to be made to constrain the mixing between Stokes I and V on the one hand and between Q and U on the other hand. In our analysis we first assumed that the feeds have equal ellipticities, which implies that the mixing between I and V is zero. In that case, the circular polarization of the reference signal from the noise diode is allowed to vary in the fit. Second, we assumed that the misorientation of the first polarization feed is zero. Under this hypothesis, the Q parameter of the noise diode signal is also allowed to vary in the fit. Finally, the data were normalized by the total invariant interval (Britton 2000) in order to limit apparent variations of the gain caused by scintillation of the pulsar across the observation.

In Fig. 9, we show a plot of the Stokes parameters as a function of time for the MJD 58816 observation of J0742-2822, along with the modeled Stokes parameters represented as dashed lines. The data shown in this figure were taken from the phase bin in which the pulsar is brightest. Stokes parameters Q and U, which describe the linear component of the polarized signal, display clear sinusoidal variations, as expected given that the feed horn rotated during the observation. The modeled Stokes parameters generally represent the actual measurements well,

with slight discrepancies (particularly in the Q and U parameters) observed during the first half of the observation. These discrepancies, which suggest possible time-dependence of the instrumental response during observations, will be investigated further in future work.

The gain, differential gain, and differential phase values, as well as the ellipticities and orientation parameters as determined from the Reception analysis, are plotted in Fig. 10. The Stokes parameters of the noise diode plotted as a function of channel index are shown in Fig. 11. The analysis finds significant non-orthogonality in most frequency channels, confirming the presence of cross-coupling between the polarization feeds. The gain and differential phase parameters are generally similar to those shown in Fig. 1, determined from a simple SINGLEAXIS analysis of the noise diode observation. Differential gains, however, strongly differ from those plotted in the latter figure in most channels. As can be seen from Fig. 11, the analysis confirms that the noise diode signal mainly consists of linear emission along the U Stokes parameter, albeit with significant Stokes Q and V emission in most frequency channels, indicating that the noise diode signal is also not ideal. The ellipticities and orientations of the feeds, as well as the Stokes parameters of the noise diode plotted in Figs. 10 and 11, need to be taken with a grain of salt since, as mentioned in the previous paragraph, strong assumptions about the non-orthogonality parameters needed to be made to constrain the mixing between I and V, and between Q and U. The analysis nevertheless demonstrates the existence of significant cross-coupling between the feeds and that the noise diode is most likely non-ideal.

The improved calibration should in principle lead to more accurate polarimetric results. We repeated the analysis of a selection of 1.4 GHz observations of MSPs with NUPPI presented in Sect. 2.1, this time using improved calibration parameters obtained from the analysis of observations of PSR J0742-2822 with feed horn rotation. For each of the MSP observations listed in Table 1, we used the PSRCHIVE command "pac" to calibrate the data using the second order calibration parameters from the closest-in-time observation of PSR J0742-2822 and the noise
<figure><img src="image_8.png" /><figcaption>**Fig. 7.** Position angle of the linear polarization as a function of pulse phase, for PSRs J0613−0200, J1022+1001, J1024−0719, J1600−3053, J1643−1224, and J1713+0747. The position angles, $\Psi$, were calculated using the *Q* and *U* Stokes parameters plotted in Figs. 3–6, as $\Psi = \frac{1}{2} \arctan\left(\frac{U}{Q}\right)$. Results from NUPPI observations calibrated with the SINGLEAXIS method are represented as orange squares. Those obtained with the improved calibration scheme are shown as green circles. Reference results from Dai et al. (2015) are shown as blue crosses. For each pulsar, the pulse profile is shown in light gray for reference.</figcaption></figure>

diode observation conducted before the MSP observation. Equation (12) was then used by *pac* to perform the calibration. The new polarimetric profiles are displayed in Figs. 3–6, and the new PA measurements are shown in Figs. 7 and 8. Additionally, in Table 1 we list the profile S/N values as obtained when calibrating the data with the simple SINGLEAXIS calibration scheme, and with the improved method. The polarimetric profiles formed with the improved calibration method are much closer to the Dai et al. (2015) profiles than those obtained with the generic calibration method, as can be seen from Figs. 3–6. Similarly, the new PA values are very different from those obtained with the SINGLEAXIS calibration method (see, for instance, the results for PSR J1713+0747 in Fig. 7 or PSR J1857+0943 in Fig. 8 for particularly discrepant results), and are also very close to the Dai et al. (2015) measurements. The latter results and those

obtained with the improved calibration method are, in all cases highlighted in Figs. 7 and 8, qualitatively consistent. The results are, however, not identical: some differences between the PA values are for example seen in the leading edge of the first peak of PSR J1022+1001 or in the main peak of J1713+0747. It should be noted that the NUPPI observation of PSR J1713+0747 was carried out about 60 days after the pulsar underwent an abrupt profile change on MJD 59321 (Xu et al. 2021; Singha et al. 2021; Jennings et al. 2022), while the reference data were taken several years before the event. As a consequence, the NUPPI and Parkes profiles are not expected to be identical. Residual calibration imperfections, or inaccurate RM values at the epochs of the NUPPI observations, could also produce discrepancies. However, the new NUPPI measurements generally agree well with the reference results, indicating that NUPPI polarimetric
<figure><img src="image_9.png" /><figcaption>**Fig. 8.** Same as Fig. 7, for PSRs J1730−2304, J1744−1134, J1857+0943, J1909−3744, J2124−3358, and J2145−0750.</figcaption></figure>

profiles calibrated with the improved scheme can be confidently used for studying the polarization emission properties of pulsars or determining RM values along their lines-of-sight.

As can be seen from Table 1, the new S/N values are in almost all cases larger than those of the profiles calibrated with the SINGLEAXIS method, with increases of about 20% in the case of PSR J1744−1134, or 30% for PSR J1024−0719. One notable exception appears to be that of PSR J1022+1001, for which the new S/N value is significantly smaller with the improved calibration method, and for which the total intensity profile as obtained with the SINGLEAXIS calibration method is closer to the reference profile than that obtained with the improved calibration. As mentioned earlier, we assumed that the mixing between Stokes *I* and *V* is zero when determining the polarimetric response. In PSR J1022+1001, the maximum of the total intensity profile coincides with the maximum of Stokes *V*, in absolute value. We could thus expect significant

differences in Stokes *I* in parts of the profile where |*V*| is large. It should also be noted that the pulse profile of PSR J1022+1001 is known to exhibit significant time variations (see e.g., Camilo et al. 1996; van Straten 2013; Padmanabh et al. 2021) over timescales of years, and could therefore differ in the NUPPI and Parkes observations considered here. Besides, the Stokes *Q*, *U*, and *V* parameters as determined with the improved scheme are more consistent with the reference results. We therefore conclude that the new polarimetric profiles for PSR J1022+1001 shown in Fig. 3 are a better representation of the true polarimetric properties of this pulsar, despite the lower S/N value.

The new differential gain values shown in Fig. 10 are significantly different from those derived from the analysis of the noise diode observation alone. S/N values are calculated using the total intensity profiles, which are formed by adding data from the two polarization feeds. Under the assumption of ideal feeds recording independent and identically distributed (i.i.d.) data, the noise
<figure><img src="image_10.png" /><figcaption>**Fig. 9.** Measured and modeled Stokes parameters for PSR J0742−2822 as a function of time and position angle $\Phi' = \Phi + \alpha$, where $\Phi$ is the parallactic angle and $\alpha$ is the orientation of the horn at a given time. The pulsar was observed at 1484 MHz with the NUPPI backend, in an observing mode where the feed horn is made to rotate by 3° per minute, during the ~1-h observation. The modeled Stokes parameters, represented as dashed lines, were determined by analyzing the data with the Reception model of PSRCHIVE, enabling the determination of the polarimetric response of the telescope. The data points shown in this plot correspond to the Stokes parameters in the phase bin in which the pulsar is brightest, with 512 bins per rotation in this case. See Sect. 2.2 for the description of the analysis.</figcaption></figure>

levels in each feed have equal variance, and the S/N is proportional to the square root of the number of polarizations summed (two in this case). However, if data are calibrated using inaccurate differential gains and unaccounted cross-couplings then the i.i.d. assumption no longer applies, thereby reducing the effective number of polarizations added to less than two. We argue that the larger S/N values in Table 1 are a consequence of the fact that the improved calibration, and particularly the substantial corrections to the differential gain parameter over large portions of the observed bandwidth, restored the i.i.d. assumption in most cases, resulting in larger effective sample sizes. Foster et al. (2015) studied the influence of inaccurate data calibration on pulsar observations and pulsar timing, based on simulations, finding that inaccurate calibration results in lower observed S/N than the ideal value. In their simulations, however, noise was added after the instrumental polarization transformation, whereas in our case the majority of the system noise was added before the differential gain and was thus subjected to the dominant boost transformation. Our conclusions are therefore consistent with and complementary to those of Foster et al. (2015).

In general, the improved calibration method led to larger S/N values for the observations we analyzed. As a consequence, NUPPI TOAs extracted from these observations should have smaller uncertainties. Furthermore, regular NRT observations of PSR J0742−2822 with horn rotation enable us to determine the polarimetric response of the instrument on a regular basis. Therefore, we also expect NUPPI TOAs to be more accurate, in addition to being more precise, as pulse profile distortions from inaccurate calibration should in principle be mitigated. Effects of the improved polarimetric calibration on NUPPI timing results are discussed in detail in the following section.

## 3. Pulsar timing with NUPPI

### 3.1. The 1.4 GHz NUPPI timing archive

To assess the influence of the improved calibration on the quality of NUPPI timing data, we consider three of the MSPs analyzed in Sect. 2: PSRs J1730−2304, J1744−1134, and J1857+0943, all three of them being well-timed PTA pulsars observed with high cadence since the beginning of NUPPI operations. The rationale for selecting these three MSPs is the following: of the MSPs considered in Sect. 2, PSR J1744−1134 is the most polarized pulsar and should thus be particularly susceptible to calibration inaccuracies. It is also among the best targets used by PTAs to search for gravitational waves (see e.g., Chen et al. 2021). Conversely, PSR J1857+0943 has the lowest degree of polarization in the sample according to Dai et al. (2015), and thus represents an interesting test case. Finally, we chose to consider PSR J1730−2304, for which van Straten (2013) predicted that the MTM technique (see Sect. 3.2) should determine TOAs with much (~30%) lower uncertainties than techniques that extract TOAs from total intensity profiles. In this subsection we present the 1.4 GHz NUPPI timing archive using PSR J1744−1134 as an illustrative example, since this pulsar has been observed with high cadence with the NRT since NUPPI was commissioned, and with previous backends also.

In Fig. 12 we show timing residuals for PSR J1744−1134 as a function of time between August 2011 and December 2022, where the residuals correspond to the differences between the measured TOAs and the predictions from an accurate timing model for the pulsar. From the figure it can be seen that pulsar timing observations have been conducted at a relatively stable cadence between August 2011 and December 2022, with
<figure><img src="image_11.png" /><figcaption>**Fig. 10.** Best-fit calibration parameters as determined from a joint analysis of an observation of PSR J0742−2822 with the feed horn rotating by 180° across the ~1-h observation, and of an observation of the noise diode conducted prior to the pulsar observation. The three bottom panels are the same as in Fig. 2. The blue points represent the best-fit results from this analysis, and the yellow squares show the results from the SINGLEAXIS analysis of the noise diode observation only, as in Fig. 2, for comparison. The two top panels show the feed ellipticities (assumed to be equal in our analysis, see Sect. 2.2 for details), and the orientations ($\theta_0$ is assumed to be 0 in the analysis).</figcaption></figure>

intervals of weeks to months when no timing observations were carried out. We here provide a list of the main interruptions, and what prevented observations from being carried out during those time intervals:

a) Painting activity in July 2013 and an intense observational campaign on the magnetar at a similar right ascension J1745−2900 (Eatough et al. 2013) prevented J1744−1134 from being monitored normally between April and July 2013.

b) A fire in the electrical cabinet of the focal carriage and subsequent repairs prevented the NRT from observing between late May and early July 2015.

c) Due to the failure of an amplifier of the low-frequency receiver in June 2016, no 1.4 GHz observations could be carried out for several weeks. Incidentally, the failure of a pump in the hydraulic system of the focal carriage in September 2016 made it unable to move and track sources.

<figure><img src="image_12.png" /><figcaption>**Fig. 11.** Q, U, and V Stokes parameters of the noise diode signal as a function of frequency, expressed in units of the flux density of the diode in the considered channel, $C_0$. The results plotted in this figure were determined by analyzing the noise diode observation shown in Fig. 1, as part of the procedure for improving the calibration of NUPPI data described in Sect. 2.2.</figcaption></figure>

The latter issue required heavy repairs, and normal observations resumed in January 2017.

d) Issues with the motion of the focal carriage along the track prevented NRT observations between August and November 2017.

e) The NRT was inactive for a few weeks in August 2018 due to electrical problems.

f) The failure of a pre-amplifier, which was replaced after a few weeks, prevented NRT observations with the low-frequency receiver from happening between late March and late May 2019.

g) A full repair of the focal track in early 2021 prevented the NRT from operating normally for a few weeks.

h) Electrical problems on the site caused the NRT to be inactive from late December 2022 to mid-February 2023.

The above-mentioned interruptions are highlighted in the top panel of Fig. 12.

## 3.2. Timing data preparation

We selected NUPPI observations of PSRs J1730−2304, J1744−1134, and J1857+0943 at 1.4 GHz, and cleaned the data of RFI using the SURGICAL method of the CoastGuard package. We then used the cleaned observations to prepare three distinct datasets for each pulsar. In the first dataset, the observations were left uncalibrated. In the second, observations were calibrated with the SINGLEAXIS method described in Sect. 2.1, i.e., using calibration parameters obtained from the analysis of the noise diode observations conducted prior to each of the pulsar
<figure><img src="image_13.png" /><figcaption>**Fig. 12.** NUPPI timing data on PSR J1744−1134. Top: Timing residuals (i.e., differences between measured pulse arrival times and arrival times predicted by a model for the pulsar) as a function of time for the millisecond pulsar J1744−1134 observed with the Nançay Radio Telescope at the central frequency of 1.484 GHz using the NUPPI backend. The 512 MHz of frequency bandwidth recorded by NUPPI were split into eight sub-bands of 64 MHz each, and the observations were integrated in time. Pulse arrival times were extracted from the timing data using the MTM method, as further described in Sect. 3. Time intervals when the NRT was inactive or when the low-frequency receiver was unavailable are highlighted with arrows (for interruptions shorter than 40 days) and shaded areas (when intervals were longer than 40 days). See Sect. 3.1 for details on the individual interruptions. Bottom: Average number of 1.4 GHz NUPPI observations of PSR J1744−1134 per month as a function of time. The average cadence was estimated by counting the number of observations within a time window of 100-d centered on the considered date, and normalizing the obtained number to 30 days.</figcaption></figure>

observations, and assuming that the noise diode signal and the polarization feeds are ideal. In the third dataset, NUPPI observations were calibrated using the improved method described in Sect. 2.2. In practice, the calibration parameters from the closest-in-time rotating horn observation of PSR J0742−2822 were used, setting a maximum threshold for the time separation of 100 days. The first observation of J0742−2822 with feed horn rotation was conducted on MJD 58816; therefore, for NUPPI observations of PSRs J1730−2304, J1744−1134, and J1857+0943 conducted before MJD 58716, a generic SINGLEAXIS calibration was applied, as in the second dataset. Finally, all observations in the three datasets were integrated in time and frequency, forming eight frequency sub-bands of 64 MHz each.

In Fig. 13, we compare the S/N values of the total intensity pulse profiles from the second and third datasets. As expected, between MJDs 55800 and ~58 700 the S/N values are identical, the calibration procedure being the same. However, it is clear from the figure that the improved calibration procedure generally leads to significantly larger S/N values as a consequence of accurate modeling of cross-couplings and polarization leakages, and particularly so for the highly polarized pulsar J1744−1134, with median increases of 4%, 20%, and 2%, respectively, for PSRs J1730−2304, J1744−1134, and J1857+0943.

In standard pulsar timing analyses, TOAs are extracted from the data by comparing total intensity profiles and a template profile for the pulsar at the considered frequency, using a cross-correlation algorithm. Wang et al. (2022) compared different

methods for generating the template profiles, and found that for bright pulsars such as the ones considered here, templates formed by smoothing the addition of a selection of high S/N observations of the pulsar were preferred over other standard template preparation methods. Wang et al. (2022) further recommended the usage of the “FDM” method, as implemented in the PSRCHIVE tool “pat” for extracting TOAs from the observations. In the FDM method, the phase-gradient of the cross-power spectrum of the template profile and the observation are fitted in the Fourier domain (Taylor 1992), as in the PGS method of pat, to determine the TOAs, but with TOA uncertainties estimated from a Monte Carlo Analysis. The FDM method was found by Wang et al. (2022) to produce TOAs and TOA uncertainty estimates that are generally more reliable than in the other tested methods. We adopted the recommendations of this systematic study of TOA creation procedures, and generated template profiles for each of the three datasets on the three pulsars by adding the total intensity profiles for the ten highest S/N observations made after MJD 58800 and then smoothing the summed profiles. The time selection was used to ensure that the template profile for the third dataset be formed by adding observations calibrated with the same procedure. TOAs for the three datasets were then extracted by comparing the observations with the respective template profiles, using the FDM method of pat.

As mentioned above, calibration observations of PSR J0742−2822 have been conducted at Nançay since
<figure><img src="image_14.png" /><figcaption>Fig. 13. Ratios of S/N and TOA uncertainty values for PSRs J1730−2304, J1744−1134, and J1857+0943 as a function of time. The values are derived from NUPPI pulse profiles observed at 1.4 GHz and calibrated using the SINGLEAXIS method (“no corr.”) and the improved method described in Sect. 2.2 (“corr.”). For each pulsar, the top panel (blue points) shows the ratios of the S/N values and the bottom panel (orange points) shows the ratios of the TOA uncertainties. Individual observations were fully integrated in time, and integrated in frequency to form eight sub-bands. TOAs were extracted from the observations using the MTM technique. See Sect. 3.2 for additional details about the data preparation.</figcaption></figure>

MJD 58816. We therefore do not have measurements of the full polarimetric response of the NRT before that date. In addition, as mentioned in Sect. 3.1, various instrumental changes, such as the replacement of a pre-amplifier in 2019, occurred over the time span of the NUPPI dataset. Some of these instrumental changes likely modified the polarimetric response of the NRT, potentially leading to varying total intensity pulse profiles at different epochs. In addition to extracting TOAs from total intensity profiles using the FDM method, we thus also investigated the use of the MTM technique presented in van Straten (2006) and implemented in pat, to extract TOAs from the NUPPI data. Unlike the FDM method, the MTM technique uses all four Stokes parameters and models a transformation between the template profile and the considered observation, to reduce timing errors caused by instrumental distortions of the pulse profile. In this approach, observed polarization profiles represented by the coherency matrix, $\rho'$, are related to the template profile for the pulsar, $\rho_0$, by the equation:

$$ \rho'(\Phi) = J \rho_0(\Phi - \Delta\Phi) J^\dagger + \rho_{\text{DC}} + \rho_{\text{noise}}(\Phi), \quad (13) $$

where $\Phi$ is the pulse phase, $\Delta\Phi$ is the phase shift, $\rho_{\text{DC}}$ and $\rho_{\text{noise}}$ are the coherency matrices for the DC offset between the profiles and for the noise, and $J$ is the polarimetric transformation between $\rho_0$ and $\rho'$. By comparing the four Stokes parameters of the observations and of the template profile in the Fourier domain, the MTM method determines the $\Delta\Phi$ parameter and the Jones matrix, $J$. We generated template profiles for each of the three datasets on the three MSPs by summing the ten highest S/N observations carried out after MJD 58 800 and smoothing the added profiles, this time keeping all four Stokes parameters. The PSRCHIVE command pat was finally used to extract TOAs using these template profiles and the MTM algorithm. In the bottom panel of Fig. 13 we compare the uncertainties of TOAs extracted with the MTM method from the observations from the second and third datasets. Before MJD ~58 700, the ratios are centered on one, but are not necessarily equal to one due to the fact that, although the procedure followed for calibrating the observations is the same for these epochs, the template profiles used for extracting the TOAs are different between the two datasets. After MJD ~58 700, TOA uncertainties for the
third dataset are generally smaller, as expected, given that the S/N values are generally higher with the improved calibration scheme.

After the above-described TOA preparation steps, we obtained six different TOA datasets for each pulsar: TOAs extracted from uncalibrated NUPPI total intensity profiles with the FDM algorithm (we will refer to these TOA datasets as “Uncal. FDM”), TOAs derived from uncalibrated NUPPI polarimetric profiles with the MTM algorithm (“Uncal. MTM”), TOAs obtained after a simple SINGLEAXIS calibration with the FDM algorithm (“Cal. FDM”) and with the MTM method (“Cal. MTM”), and finally TOAs obtained after application of the improved calibration method with the FDM algorithm (“Cal. + corr. FDM”) and with the MTM method (“Cal. + corr. MTM”). Figure 14 presents histograms of the TOA uncertainties for the six different TOA datasets on each MSP. One can see that, as expected, the improved polarimetric calibration scheme leads to the smallest TOA uncertainties among the MTM and FDM datasets, while the largest TOA uncertainties are obtained for the uncalibrated data. Additionally, for the “Cal. + corr.” datasets, we find that the ratios of the median TOA uncertainties estimated by the MTM and FDM algorithms are respectively 0.68, 1.49, and 0.91 for PSRs J1730−2304, J1744−1134, and J1857+0943. These values are very close to the theoretical expectations listed in Table 1 of van Straten (2013), of 0.71, 1.56, and 0.89, respectively. For PSR J1744−1134, which is almost 100% linearly polarized (see Dai et al. 2015), MTM finds larger TOA uncertainties due to large correlations between the phase shift and the parameters describing the polarization transformation. As can be seen from Table 1 of van Straten (2013), MTM is expected to find smaller TOA uncertainties in most pulsars, because of their lower fractional polarizations.

### 3.3. Timing results

The “Cal. + corr. MTM” TOA datasets were built from better calibrated observations than the other datasets, and the TOAs were extracted using the MTM algorithm, which models putative profile distortions caused by calibration artifacts at the same time as it determines the TOAs. Of the different TOA datasets described in Sect. 3.2, they thus likely represent our best datasets for timing analyses of PSRs J1730−2304, J1744−1134, and J1857+0943. We analyzed the three “Cal. + corr. MTM” datasets using the TEMPO2 pulsar timing package (Hobbs et al. 2006). Topocentric TOAs were converted to Barycentric Coordinate Time (TCB) using the DE440 solar system ephemeris described in Park et al. (2021). The timing models were purposefully kept simple: for each pulsar we fitted for the sky coordinates, the proper motion and timing parallax, the pulsar’s spin frequency and its first two time derivatives, the orbital parameters in the case of PSR J1857+0943, and for the DM and its first two time derivatives. That is, we did not model any spin frequency or DM red noise components that could absorb timing instabilities caused by calibration issues. Similarly, we also did not model putative systematic errors in TOAs and/or in their uncertainties using EFAC and EQUAD parameters (see e.g., Lentati et al. 2014). The best-fit ephemerides are consistent with previously published results on these pulsars (e.g., Desvignes et al. 2016). In all cases the second spin frequency derivatives were significant, indicating the presence of spin noise in these pulsars.

The timing residuals for PSR J1744−1134 as a function of time for the best-fit model and the “Cal. + corr. MTM” TOA dataset are plotted in Fig. 12. The weighted rms residual is

<figure><img src="image_15.png" /><figcaption>Fig. 14. Histograms of the TOA uncertainties, with TOAs extracted from uncalibrated NUPPI data (“Uncal.”), data calibrated using the SINGLEAXIS method (“Cal.”) and data calibrated with the improved calibration method described in Sect. 2.2 (“Cal. + corr.”). Three of the TOA datasets were derived from the NUPPI data using the FDM algorithm; the other three were extracted using the MTM technique.</figcaption></figure>

~0.628 µs, confirming PSR J1744−1134 as one of the most precisely timed pulsars known currently. The reduced $\chi^2$ value is ~2, indicating that the timing data could be whitened further than what is done by our simple ephemeris for this pulsar. We used the best-fit model for each pulsar to determine timing residuals for the five other TOA datasets, and compared the timing residuals with those for the “Cal. + corr. MTM” dataset. In Fig. 15 we show the differences in the timing residuals (i.e., the differences between the residuals for the various TOA datasets, and those for the “Cal. + corr. MTM” dataset) for PSR J1744−1134. As expected, the largest differences are
<figure><img src="image_16.png" /><figcaption>**Fig. 15.** Differences between timing residuals of PSR J1744−1134 for the TOA datasets described in Sect. 3.2, and those for the “Cal. + corr. MTM” TOA dataset, as a function of time. In all cases, the residuals were determined using an ephemeris for the pulsar obtained by analyzing the “Cal. + corr. MTM” dataset with TEMPO2. Colored symbols represent the residual differences, and black dots show averages over 30-d sliding windows.</figcaption></figure>

seen for the “Uncal. FDM” dataset (see panel a), with significant structures in the residual differences and many data points well away from zero. These differences reflect changes in the properties of the instrument that strongly influence the total intensity pulse profile of J1744−1134 at a given epoch, and can in particular make it inconsistent with the template profile used for extracting the TOAs. Interestingly, much less structure is seen for the “Uncal. MTM” TOAs (panel b), demonstrating that the MTM algorithm manages to compensate for many of the profile changes induced by instrumental artifacts, albeit not completely, likely due to bandwidth depolarization (van Straten 2002) in the 8 × 64 MHz frequency sub-bands. Much smaller residual differences are observed in the case of the “Cal. FDM” dataset (panel c) than for the “Uncal. FDM” TOAs. However, a systematic offset is detected between the data taken before and after the pre-amplifier change at the NRT in 2019, showing that the total intensity pulse profiles are inconsistent in these two parts of the NUPPI archive. This offset is compensated for by the MTM algorithm, as can be seen from panel d of Fig. 15. In the latter case, the residual differences are in all cases close to zero, indicating that for this pulsar the MTM algorithm finds similar results, independently of the calibration scheme. Finally, it can be seen from *e* that the above-mentioned offset is not seen in data prepared with the improved polarimetric calibration scheme. However, as in panels a and c, significant structure is seen in the residual differences, corresponding to cases where the MTM algorithm found differences between the observed pulse profiles and the template profile, which it compensated for. For the other two MSPs, J1730−1134 and J1857+0943, differences

between timing residuals display the same behavior, albeit with more noise than in J1744−1134 due to the fact that they are less precisely timed objects.

Strictly speaking, although significant structure is seen in panel e, we cannot tell from Fig. 15 if one of the two TOA datasets (“Cal. + corr. FDM” or “Cal. + corr. MTM”) is a better dataset than the other. We therefore conducted independent TEMPO2 analyses of all TOA datasets for the three pulsars to obtain best-fit timing models for each dataset and in turn determine weighted rms and reduced χ² values. These parameters are listed in Table 2. As expected, of the three TOA datasets extracted from total intensity pulse profiles with the FDM algorithm, the lowest weighted rms residual and reduced χ² values are obtained for the dataset that is based on the best calibration method. The same conclusion also applies to TOA datasets derived using the MTM algorithm, although in this case, the “Cal. MTM” and “Cal. + corr. MTM” datasets have similar weighted rms residuals and reduced χ² values. In Fig. 16 we show the distributions of the reduced χ² values of the residual profiles from template matching for each pulsar and TOA dataset. The median goodness-of-fit values are in all cases close to unity, indicating that template matching fits are generally good. However, for all three pulsars, the ranges between the 5th and 95th percentiles are significantly larger for the uncalibrated data, and in the case of the calibrated data the ranges for MTM TOAs are systematically smaller than for FDM TOAs. Profile residuals therefore seem to provide a metric for assessing the correctness of the polarization calibration in this case. We note finally that in all three MSPs the weighted rms residuals and the
<table><caption>**Table 2.** Weighted rms residual values and reduced $\chi^2$ values from individual timing analyses of the six TOA datasets described in Sect. 3.2 with TEMPO2.</caption><thead><tr><th>TOA dataset</th><th>Uncal. FDM</th><th>Cal. FDM</th><th>Cal. + corr. FDM</th><th>Uncal. MTM</th><th>Cal. MTM</th><th>Cal. + corr. MTM</th></tr></thead><tbody><tr><td colspan="7" style="text-align: center;">J1730-2304</td></tr><tr><td>Weighted rms residual (µs)</td><td>2.213</td><td>1.710</td><td>1.604</td><td>1.090</td><td>0.974</td><td>0.960</td></tr><tr><td>Reduced χ²</td><td>4.091</td><td>2.426</td><td>2.259</td><td>2.149</td><td>1.839</td><td>1.843</td></tr><tr><td colspan="7" style="text-align: center;">J1744-1134</td></tr><tr><td>Weighted rms residual (µs)</td><td>1.612</td><td>1.042</td><td>0.705</td><td>0.905</td><td>0.635</td><td>0.628</td></tr><tr><td>Reduced χ²</td><td>27.951</td><td>10.034</td><td>5.648</td><td>4.480</td><td>1.895</td><td>2.041</td></tr><tr><td colspan="7" style="text-align: center;">J1857+0943</td></tr><tr><td>Weighted rms residual (µs)</td><td>1.364</td><td>1.359</td><td>1.280</td><td>1.272</td><td>1.166</td><td>1.120</td></tr><tr><td>Reduced χ²</td><td>2.001</td><td>1.958</td><td>1.903</td><td>2.203</td><td>1.898</td><td>1.817</td></tr></tbody></table>

**Notes.** The TOA datasets for PSRs J1730–2304, J1744–1134, and J1857+0943, respectively, comprised 2379, 2085, and 1612 TOAs.

reduced $\chi^2$ values from the timing analyses are higher for the “Cal. + corr. FDM” TOA datasets than for their MTM counterparts. This analysis confirms our initial assumption that the “Cal. + corr. MTM” datasets, which are based on better calibrated observations and include TOAs extracted with the MTM algorithm to account for potential changes of the instrumental response with time, are our best datasets for conducting the timing analysis of PSRs J1730–2304, J1744–1134, and J1857+0943 with NUPPI. Additionally, the analysis independently confirms the conclusion from Foster et al. (2015), based on simulations, that the MTM method produces better timing than methods that use total intensity profiles. van Straten (2006) noted that the MTM method may not perform well when the analyzed pulsar has a low degree of polarization or the estimated phase shift is highly covariant with the MTM model parameters that describe the unknown polarimetric transformation. In these situations where the FDM method or another method using the total intensity pulse profiles would be preferred over the MTM algorithm, our analysis has shown that better calibrated data lead to better timing datasets, as expected.

## 4. Summary and prospects

Until recently, pulsar data recorded with the NUPPI backend of the NRT at 1.4 GHz were polarization-calibrated using the observations of a noise diode that injects a polarized reference signal into the feeds, assuming that both the instrumental response and the noise diode signal are ideal. In late 2019, we began observing the bright, highly linearly polarized pulsar J0742–2822 at a regular cadence, in a special observing mode in which the feed horn, which is fixed during normal pulsar observations, rotates by ~180° across the typically 1-h observation. Doing so enables us to more fully determine the polarimetric response of the telescope at the epochs of these observations as well as the actual Stokes parameters of the noise diode signal. We found significant cross-couplings between the polarization feeds and found that the noise diode does not emit 100% linearly polarized emission along the Stokes U parameter, as was assumed prior to this work. We applied these improved calibration parameters to the observations of a selection of MSPs with published polarimetric results, and found our new results to be consistent with those found by Dai et al. (2015) using the Parkes radio telescope. Finally, we analyzed NUPPI timing data on the

MSPs J1730–2304, J1744–1134, and J1857+0943, and found that TOAs extracted from total intensity pulse profiles obtained with the improved calibration procedure have lower residual rms and reduced $\chi^2$ values than TOAs extracted from uncalibrated data or from data calibrated with the generic method used prior to this work. We also found that the MTM method for extracting TOAs from polarimetric profiles while compensating for potential calibration errors gives excellent results for these three MSPs, with TOA residuals that have lower rms and $\chi^2$ values than those obtained with TOAs extracted from total intensity profiles. It should be noted that, even though calibration errors are expected to affect the pulse profiles of pulsars differently depending on their polarization properties, they introduce time-correlated signals in the timing data of different pulsars, and thus contribute noise to searches for gravitational waves by PTAs. Correcting polarization calibration errors is thus critical to these searches.

The polarimetric response of the NRT at the epochs of the rotating horn observations and the properties of the reference source were determined by solving Eq. (11). As explained in detail in van Straten (2004), solutions to this equation are degenerate under commutation, and assumptions need to be made to constrain the mixing between *I* and *V* and between *Q* and *U*. We assumed that the two feeds have equal ellipticities and that the misorientation of the first receptor is zero. A possibility for determining the individual receptor ellipticities would be to combine the rotating horn observations with observations of standard flux candles. In future work we will carry out such observations of standard candles, close in time to rotating horn observations, to determine the individual ellipticities across the bandwidth covered by NUPPI. Additionally, some of the properties of the NRT, such as the gain of the telescope, are known to vary with the declination of the observed source (see e.g., Theureau et al. 2005). In this work we have calibrated observations of MSPs at different sky locations using instrumental responses determined from observations of J0742–2822, which may or may not be optimal, given the potential dependency of some of the calibration parameters on declination or other pulsar parameters. This, and the possible variations of the instrumental response during observations (as noted in Sect. 2.2), will also be investigated in future work. Pulsar observations conducted in a mode in which appropriate feed horn rotation would maintain a fixed parallactic angle across the observations may be informative in this regard.
<figure><img src="image_17.png" /><figcaption>**Fig. 16.** Distributions of the reduced $\chi^2$ values of the residual profiles from template matching for each TOA, for the three pulsars and six dataset types considered in Sect. 3.3. Orange lines represent the median goodness-of-fit values, boxes encompass the first to third quartiles, whiskers encompass the 5th and 95th percentiles, and gray dots represent outliers.</figcaption></figure>

As mentioned in Sect. 2, the first pulsar observations with the NUPPI backend were conducted in August 2011, and regular observations of J0742–2822 with feed horn rotation began in late 2019. The improved calibration solutions derived from these later data cannot be directly applied to earlier observations because the response of the NRT likely evolved with the various instrumental changes that occurred since 2011. However, using the well-calibrated observations of pulsars observed since late 2019 as polarized reference sources, it is possible to model the instrumental response over the first eight years of NUPPI observations using Measurement Equation Template Matching (METM; van Straten 2013). The improved calibration solutions for the NUPPI observation archive can be expected to lead to improved timing results, which is of great interest for pulsar timing studies and PTA analyses.

Observations of PSR J0742–2822 with feed horn rotation are conducted at Nançay to overcome the fact that the NRT cannot sample wide ranges of parallactic angles when observing a given pulsar. The work presented in this article is thus relevant for any telescope affected by a similar limitation. The results obtained with the improved calibration procedure demonstrate that accurate polarization calibration is essential for any pulsar polarimetric or timing study, and the need for accurate calibration will get more and more acute for the increasingly sensitive radio telescopes of the future.

*Acknowledgements.* We would like to thank the anonymous referee for very constructive and valuable comments that helped to improve the manuscript. The Nançay Radio Observatory is operated by the Paris Observatory, associated with the French Centre National de la Recherche Scientifique (CNRS). We acknowledge financial support from the “Programme National de Cosmologie et Galaxies” (PNCG) and “Programme National Hautes Énergies” (PNHE) of CNRS/INSU, France. Throughout this work, we made extensive use of the

ATNF Pulsar Catalogue database (Manchester et al. 2005), available at: https://www.atnf.csiro.au/research/pulsar/psrcat/. Part of this research has also made use of the EPN Database of Pulsar Profiles maintained by the University of Manchester, available at: http://www.jodrellbank.manchester.ac.uk/research/pulsar/Resources/epn/. The published data displayed in Figs. 3–6 were obtained from the EPN database under the Creative Commons Attribution 4.0 International licence.

## References

Bailes, M., Harrison, P. A., Lorimer, D. R., et al. 1994, ApJ, 425, L41

Bailes, M., Johnston, S., Bell, J. F., et al. 1997, ApJ, 481, 386

Bonsignori-Facondi, S. R., Salter, C. J., & Sutton, J. M. 1973, A&A, 27, 67

Britton, M. C. 2000, ApJ, 532, 1240

Camilo, F., Nice, D. J., Shrauner, J. A., & Taylor, J. H. 1996, ApJ, 469, 819

Chen, S., Caballero, R. N., Guo, Y. J., et al. 2021, MNRAS, 508, 4970

Cognard, I., Theureau, G., Guillemot, L., et al. 2013, in SF2A-2013: Proceedings of the Annual meeting of the French Society of Astronomy and Astrophysics, eds. L. Cambresy, F. Martins, E. Nuss, & A. Palacios, 327

Cognard, I., Freire, P. C. C., Guillemot, L., et al. 2017, ApJ, 844, 128

Dai, S., Hobbs, G., Manchester, R. N., et al. 2015, MNRAS, 449, 3223

Desvignes, G., Barott, W. C., Cognard, I., Lespagnol, P., & Theureau, G. 2011, AIP Conf. Ser., 1357, 349

Desvignes, G., Caballero, R. N., Lentati, L., et al. 2016, MNRAS, 458, 3341

Desvignes, G., Cognard, I., Smith, D. A., et al. 2022, A&A, 667, A79

DuPlain, R., Ransom, S., Demorest, P., et al. 2008, SPIE Conf. Ser., 7019, 70191D

Dyks, J., Serylak, M., Oslowski, S., et al. 2016, A&A, 593, A83

Eatough, R. P., Falcke, H., Karuppusamy, R., et al. 2013, Nature, 501, 391

Foster, R. S., Wolszczan, A., & Camilo, F. 1993, ApJ, 410, L91

Foster, G., Karastergiou, A., Paulin, R., et al. 2015, MNRAS, 453, 1489

Guillemot, L., Octau, F., Cognard, I., et al. 2019, A&A, 629, A92

Guo, Y. J., Freire, P. C. C., Guillemot, L., et al. 2021, A&A, 654, A16

Hobbs, G. B., Edwards, R. T., & Manchester, R. N. 2006, MNRAS, 369, 655

Hotan, A. W., van Straten, W., & Manchester, R. N. 2004, PASA, 21, 302

Jacoby, B. A., Bailes, M., van Kerkwijk, M. H., et al. 2003, ApJ, 599, L99
Jacoby, B. A., Bailes, M., Ord, S. M., Knight, H. S., & Hotan, A. W. 2007, ApJ, 656, 408

Jankowski, F., van Straten, W., Keane, E. F., et al. 2018, MNRAS, 473, 4436

Jennings, R. J., Cordes, J. M., Chatterjee, S., et al. 2022, ApJ, submitted [arXiv:2210.12266]

Johnston, S. 2002, PASA, 19, 277

Johnston, S., & Kerr, M. 2018, MNRAS, 474, 4629

Johnston, S., Lorimer, D. R., Harrison, P. A., et al. 1993, Nature, 361, 613

Kramer, M., Stairs, I. H., Manchester, R. N., et al. 2021, Phys. Rev. X, 11, 041050

Kraus, J. D. 1960, AJ, 65, 54

Kraus, J. D. 1966, Radio astronomy (New York: McGraw-Hill)

Lazarus, P., Karuppusamy, R., Graikou, E., et al. 2016, MNRAS, 458, 868

Lentati, L., Alexander, P., Hobson, M. P., et al. 2014, MNRAS, 437, 3004

Liu, K., Guillemot, L., Istrate, A. G., et al. 2020, MNRAS, 499, 2276

Lorimer, D. R., & Kramer, M. 2004, Handbook of Pulsar Astronomy (Cambridge, UK: Cambridge University Press), 4

Lorimer, D. R., Nicastro, L., Lyne, A. G., et al. 1995, ApJ, 439, 933

Manchester, R. N., Hobbs, G. B., Teoh, A., & Hobbs, M. 2005, AJ, 129, 1993

Ord, S. M., van Straten, W., Hotan, A. W., & Bailes, M. 2004, MNRAS, 352, 804

Padmanabh, P. V., Barr, E. D., Champion, D. J., et al. 2021, MNRAS, 500, 1178

Park, R. S., Folkner, W. M., Williams, J. G., & Boggs, D. H. 2021, AJ, 161, 105

Philippov, A., & Kramer, M. 2022, ARA&A, 60, 495

Segelstein, D. J., Rawley, L. A., Stinebring, D. R., Fruchter, A. S., & Taylor, J. H. 1986, Nature, 322, 714

Singha, J., Surnis, M. P., Joshi, B. C., et al. 2021, MNRAS, 507, L57

Stairs, I. H. 2003, Living Rev. Relativ., 6, 5

Stinebring, D. R., Cordes, J. M., Rankin, J. M., Weisberg, J. M., & Boriakoff, V. 1984, ApJS, 55, 247

Taylor, J. H. 1992, Philos. Trans. Royal Soc. A, 341, 117

Theureau, G., Coudreau, N., Hallet, N., et al. 2005, A&A, 430, 373

Theureau, G., Parent, D., Cognard, I., et al. 2011, A&A, 525, A94

van Straten, W. 2002, ApJ, 568, 436

van Straten, W. 2004, ApJS, 152, 129

van Straten, W. 2006, ApJ, 642, 1004

van Straten, W. 2013, ApJS, 204, 13

van Straten, W., Manchester, R. N., Johnston, S., & Reynolds, J. E. 2010, PASA, 27, 104

Verbiest, J. P. W. & Shaifullah, G. M. 2018, Class. Quant. Gravity, 35, 133001

Voisin, G., Cognard, I., Freire, P. C. C., et al. 2020, A&A, 638, A24

Wang, J., Shaifullah, G. M., Verbiest, J. P. W., et al. 2022, A&A, 658, A181

Wu, J., Clark, C. J., Pletsch, H. J., et al. 2018, ApJ, 854, 99

Xilouris, K. M. 1991, A&A, 248, 323

Xu, H., Huang, Y. X., Burgay, M., et al. 2021, ATel, 14642, 1