Issue 
A&A
Volume 674, June 2023



Article Number  A71  
Number of page(s)  13  
Section  Stellar structure and evolution  
DOI  https://doi.org/10.1051/00046361/202346228  
Published online  02 June 2023 
Radio timing constraints on the mass of the binary pulsar PSR J1528−3146
^{1}
Laboratoire de Physique et Chimie de l’Environnement et de l’Espace, Université d’Orléans/CNRS, 45071 Orléans Cedex 02, France
email: anais.berthereau@cnrsorleans.fr
^{2}
Observatoire Radioastronomique de Nançay, Observatoire de Paris, Université PSL, Université d’Orléans, CNRS, 18330 Nançay, France
^{3}
MaxPlanckInstitut für Radioastronomie, Auf dem Hügel 69, 53121 Bonn, Germany
^{4}
LUTH, Observatoire de Paris, Université PSL, Université Paris Cité, CNRS, 92195 Meudon, France
^{5}
Centre for Astrophysics and Supercomputing, Swinburne University of Technology, PO Box 218 Hawthorn, Vic 3122, Australia
^{6}
ARC Centre of Excellence for Gravitational Wave Discovery, OzGrav, Mail H29, Swinburne University of Technology, PO Box 218 Hawthorn, VIC 3122, Australia
^{7}
Australia Telescope National Facility, CSIRO, Space and Astronomy, PO Box 76 Epping, NSW 1710, Australia
Received:
23
February
2023
Accepted:
12
April
2023
Context. PSR J1528−3146 is a 60.8 ms pulsar orbiting a heavy white dwarf (WD) companion, with an orbital period of 3.18 d. The pulsar was discovered in the early 2000 s in a survey at 1.4 GHz of intermediate Galactic latitudes conducted with the Parkes radio telescope. The initial timing analysis of PSR J1528−3146, using data recorded from 2001 and 2004, did not reveal any relativistic perturbations to the orbit of the pulsar or to the propagation of its pulses. However, with an orbital eccentricity of ∼0.0002 and a large companion mass on the order of 1 M_{⊙}, this system has been deemed likely to exhibit measurable perturbations.
Aims. This work is aimed at characterizing the pulsar’s astrometric, spin, and orbital parameters by analyzing timing measurements conducted at the Parkes, MeerKAT, and Nançay radio telescopes over nearly two decades. The measurement of postKeplerian perturbations to the pulsar’s orbit can be used to constrain the masses of the two component stars of the binary and, in turn, to offer insights into the history of the system.
Methods. We analyzed timing data from the Parkes, MeerKAT, and Nançay radio telescopes collected over about 16 yr, obtaining a precise rotation ephemeris for PSR J1528−3146. A Bayesian analysis of the timing data was carried out to constrain the masses of the two components and the orientation of the orbit. We further analyzed the polarization properties of the pulsar to constrain the orientation of the magnetic axis and of the line of sight with respect to the spin axis.
Results. We measured a significant rate of advance of periastron, for the first time, and we set constraints on the Shapiro delay in the system and on the rate of change of the projected semimajor axis of the pulsar’s orbit. The Bayesian analysis yielded measurements for the pulsar and companion masses of M_{p} = 1.61_{−0.13}^{+0.14} M_{⊙} and M_{c} = 1.33_{−0.07}^{+0.08} M_{⊙} (68% C.L.), respectively, confirming that the companion is indeed massive. This companion mass as well as other characteristics of PSR J1528−3146 indicate that this pulsar is very similar to PSR J2222−0137, a 32.8 ms pulsar orbiting a WD whose heavy mass (∼1.32 M_{⊙}) has been considered unique among pulsarWD systems until now. Our measurements suggest common evolutionary scenarios for PSRs J1528−3146 and J2222−0137.
Key words: ephemerides / pulsars: individual: J1528–3146
© The Authors 2023
Open Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.
1. Introduction
Pulsars are highly magnetized, rapidly rotating neutron stars that emit periodic trains of radio pulses. The pulsar population can be divided into two main groups: “normal” pulsars and “millisecond” pulsars (MSPs). The latter are defined as pulsars with small rotational periods (P ≲ 30 ms) and very small period increases (Ṗ ≲ 10^{−17} s s^{−1}). In most cases, MSPs are found to be members of binary systems, and are thought to have been “recycled” by the accretion of matter and transfer of angular momentum from their binary companion, spinning up their rotation to millisecond periods (BisnovatyiKogan & Komberg 1974; Alpar et al. 1982). In some cases, this accretion process can stop before the pulsar gets fully recycled, leading to socalled “mildly recycled” pulsars with rotational periods between 20 and 100 ms. This process happens mainly if the companion stars are more massive: such stars evolve more rapidly and therefore any accretion episodes will generally be much shorter.
PSR J1528−3146 is a member of the latter subclass of mildlyrecycled pulsars. This P ∼ 60.8 ms pulsar was discovered in a survey conducted at 1.4 GHz with the Parkes radio telescope (also known as ‘Murriyang’) in Australia, as part of an extension of the Swinburne Intermediate Latitude survey, which discovered a total of 26 pulsars (Jacoby et al. 2007). Timing measurements conducted from 2001 to 2004 revealed that the pulsar is in a binary system with an orbital period of P_{b} ∼ 3.18 days, a projected semimajor axis x ≡ a sin i/c ∼ 11.4 s (where a is the semimajor axis of the pulsar’s orbit, i is the orbital inclination, and c is the speed of light), and an orbital eccentricity e ∼ 0.0002 which is large for systems with comparable orbital periods (e.g., Berezina et al. 2017). The orbital period and semimajor axis imply a relatively massive companion: assuming a pulsar mass of 1.4 M_{⊙}, the mass function of 0.159 M_{⊙} indicates that the minimum companion mass (calculated for an edgeon orbit, i.e., with an inclination of 90°) is ∼0.94 M_{⊙}. With such properties the system can be classified as an intermediatemass binary pulsar (IMBP, see e.g., Camilo 1996; Tauris et al. 2012); the large mass of the companion suggests that it must be a CO or ONeMg white dwarf (WD). Optical observations of the field around PSR J1528−3146 with the 6.5 m Baade telescope at the Magellan Observatory led to the detection of a potential counterpart, with an Rband magnitude of about 24.2 and a Bband magnitude of about 24.5 (Jacoby et al. 2006). If the counterpart is indeed associated with the J1528−3146 system, then the optical emission very likely originates from the WD.
The initial timing study by Jacoby et al. (2007) over approximately three years did not reveal any relativistic perturbations to the pulsar’s orbit, such as the rate of periastron advance. Nevertheless, the nonnegligible orbital eccentricity and large companion mass make it a promising system for measuring such relativistic perturbations, that can be described in pulsar timing analyses by socalled “postKeplerian” (PK) parameters (e.g., Damour & Taylor 1992).
Two decades after the discovery of PSR J1528−3146, the accumulation of highquality pulsar timing data by the Nançay radio telescope (NRT) in France and the advent of highly sensitive newgeneration instruments such as the MeerKAT 64dish array in South Africa (Jonas 2009) and the ultrawideband low receiver at the Parkes telescope (UWL, see Hobbs et al. 2020) motivated a new study of this pulsar and of its timing properties. In particular, a new timing study of PSR J1528−3146 should enable more precise measurements of its parameters and could yield first detections of PK parameters in this system and possibly measurements of the masses of the pulsar and of its companion.
In this article we report on the results from the analysis of a ∼16 yr timing data set on PSR J1528−3146 with the Parkes, NRT, and MeerKAT radio telescopes. The article is organized as follows: in Sect. 2, we describe the details of the observations of PSR J1528−3146 and the extraction of pulsar timing data from these observations. In Sect. 3, we present an analysis of the linear polarization properties of the pulsar and the implications of the results in terms of emission geometry. In Sect. 4, we present the main timing results from our study and the constraints on the mass of the pulsar and of its companion. We finally summarize our results and discuss on the main implications from this study in Sect. 5.
2. Data
2.1. Observations
The NRT dataset consists of observations carried out with the BerkeleyOrléansNançay (BON) backend until August 2011 and the NUPPI (a version of the Green Bank Ultimate Pulsar Processing Instrument, GUPPI, made for the NRT) backend afterwards.
The BON observations of PSR J1528−3146 started in November 2006, and were carried out at a central frequency of 1398 MHz. These observations were coherently dedispersed and carried out over a frequency bandwidth of 64 MHz until July 2008, and over a bandwidth of 128 MHz after that date.
The NUPPI instrument, which became the primary pulsar timing instrument at Nançay after August 2011, increased the bandwidth of NRT pulsar observations to 512 MHz that are also coherently dedispersed. NUPPI observations were conducted at a central frequency of 1484 MHz. In this study we considered NUPPI data recorded until May 2022. Although pulsar observations with the BON backend were maintained for another few years after 2011 for testing purposes, we discarded BON observations made simultaneously with NUPPI observations, as they represented duplicated astrophysical information.
The Parkes radio telescope dataset included the timing data presented in Jacoby et al. (2007), recorded between September 2001 and September 2004, and more recent data taken with the UWL recorded from October 2019 to January 2021. The former observations were conducted with the 21cm multibeam receiver (StaveleySmith et al. 1996) at a central frequency of 1390 MHz, using the 2 × 0.5 × 512 MHz analogue filterbank backend (AFB) and were thus incoherently dedispersed. Observations with the UWL receiver were carried out using the MEDUSA backend, that performed coherent dedispersion over a much larger bandwidth of 3328 MHz, centered at 2368 MHz.
The MeerKAT dataset we used for the study is comprised of observations made between March 2019 and April 2021, as part of the RelBin program of the MeerTime large science project. The RelBin program aims at improving the measurement and detecting new relativistic effects in known pulsars, increasing the number of pulsars with mass measurements, and testing gravity theories, with the excellent timing capabilities of MeerKAT (Kramer et al. 2021a). Pulsar observing setups, verifications and early science results with MeerKAT are presented in Bailes et al. (2020). Observations were calibrated in polarization as described in Serylak et al. (2021). MeerKAT observations were conducted at a central frequency of 1284 MHz, and covered a frequency bandwidth of about 800 MHz. While the bulk of the observations had a typical duration of about 30 min, an intensive orbital campaign was conducted around 2019 December 12 (MJD 58829). On that day, the system was observed for about four hours to cover a wide range of orbital phases centered on superior conjunction in order to significantly increase the sensitivity of our subsequent timing analyses to a faint Shapiro delay (Shapiro 1964) signal that had already been detected in Nançay data. Additionally, three other 30 min observations of the pulsar were conducted within three days centered on MJD 58829, in order to sample orbital phase ranges away from superior conjunction.
In Fig. 1 we show the MeerKAT pulse profile for PSR J1528−3146 from the latter observation, which represented the highest S/N observation in our total dataset. The main radio component in the pulse profile is at rotational phase ∼0.54 in this graph. As can be seen from Fig. 1, the unprecedented sensitivity of MeerKAT observations has enabled us to detect a previously unseen weak component separated from the main pulse by approximately half a rotation. The presence of this interpulse component provides useful information about the pulsar’s emission geometry, as discussed further in Sect. 3.
Fig. 1. MeerKAT pulse profiles at 1284 MHz for PSR J1528−3146, from a ∼4h observation at MJD 58829, corresponding to the highest S/N observation in our dataset. Left: total intensity as a function of observing frequency and rotational phase. Right: integrated pulse profile, with the main radio pulse centered at rotational phase ∼0.54 and a previouslyunknown weak radio component at phase ∼0.03, indicated with an arrow. The inset at the top left shows a zoomin of the pulse profile around the weak radio component. 
2.2. Data analysis
All the data reduction steps for the datasets described above were done using the PSRCHIVE analysis library (Hotan et al. 2004). Observations were cleaned of radio frequency interference (RFI) and also polarizationcalibrated.
The typically ∼1 h long Nançay observations with the BON and NUPPI backends were split into subintegrations of 10 min each. For the BON dataset, these subintegrations were split into two frequency subbands of 32 MHz each for the data recorded before July 2008, and two subbands of 64 MHz for the data taken after the bandwidth increase. Subintegrations for the NUPPI data were divided into eight frequency subbands of 64 MHz each. For the BON and NUPPI data, we combined 10 observations with the highest signaltonoise ratios (S/N) that we then smoothed to produce noisefree standard profiles at 1.4 GHz with the same frequency resolution. Times of arrival (TOAs) to be used in the timing analysis (described in Sect. 4) were derived from the Nançay data using the Matrix Template Matching (MTM) method of the pat routine of PSRCHIVE (van Straten 2006).
The TOAs were extracted from the “historic” Parkes dataset by comparing observations with a high S/N template profile, as described in Jacoby et al. (2007). Each observation was split into two frequency subbands of 128 MHz, so that two TOAs were extracted from each observation. On the other hand, the very large frequency bandwidth covered by Parkes observations with the UWL prompted us to consider a different approach for extracting TOAs. Instead of splitting the very large bandwidth into several subbands (which often results in subbands containing no clear detections at the highest frequencies), we used the wideband template matching technique implemented in the PulsePortraiture library^{1} (Pennucci et al. 2014) and extracted one TOA for the entire bandwidth for every 10 min of observation.
Finally, the MeerKAT observations were divided into eight frequency subbands 97 MHz (resp. 107 MHz) for observations with a 759 MHz (resp. 856 MHz) bandwidth and the subintegrations were grouped into segments of 10 min. For each resulting element of time and frequency, we extracted a TOA using the standard Fourier Domain Markov chain Monte Carlo (FDM) method of pat, comparing the profiles with a smoothed version of the highest S/N MeerKAT observation with the same number of frequency subbands. A summary of the TOA dataset produced for the timing analysis is presented in Table 1.
Summary of observation parameters and TOA properties for each of the datasets considered in the timing study.
3. Pulse polarization analysis
As mentioned in Sect. 2, the unprecedented sensitivity of MeerKAT radio observations of PSR J1528−3146 enabled us to detect a weak profile component that is well separated from the main emission pulse, for the first time. This newly detected weak component can be seen in the polar representation of the pulse profile from the MeerKAT observation at MJD 58829 shown in Fig. 2. The dashed lines show the position of this profile component and that of the same position shifted by half a rotation. The phase location was determined from a fit of the pulse profile displayed in Fig. 1 with a Lorentzian line for the highest peak and two Gaussian lines for the other two components, on top of a constant baseline. Denoting ϕ_{M, 1} and ϕ_{M, 2} as the phases of the two components of the main pulse, and ϕ_{IP} that of the weak component, we found ϕ_{M, 1} = 0.469 (2), ϕ_{M, 2} = 0.54410 (2), and ϕ_{IP} = 0.033 (8) from the fit of the pulse profile. The new component is therefore separated from the main emission cone by almost exactly half a rotation and it can thus be considered as interpulse emission, likely originating from the magnetic pole opposite to the one responsible for the main emission pulse. The presence of this interpulse component indicates that the angle between the pulsar’s spin axis and its magnetic axis, α, must be close to 90°.
Fig. 2. Polar representation of the MeerKAT pulse profile at 1284 MHz displayed in Fig. 1, with radio flux densities shown in logarithmic scale. The black dashed line shows the position of the interpulse component and the red dashed line shows the same position shifted by 180°. 
Tighter constraints on the value of α and that of ζ, namely, the angle between the pulsar’s spin axis and the observer’s line of sight, may be placed by fitting the variation of the polarization position angle (PPA) Ψ of the linearly polarized component of the pulsar as a function of the rotational phase, with a rotating vector model (RVM, Radhakrishnan & Cooke 1969). The position angles Ψ are related to the Q and U Stokes parameters of the pulsar’s emission as . In the following, we use the modified form introduced in Johnston & Kramer (2019):
In the above expression, ϕ_{0} denotes the pulse phase for which Ψ = Ψ_{0} (assuming Δ = 0), and the term Δ can take account of a potential difference in emission heights between the two magnetic poles of the pulsar. We used PPA values as a function of rotational phase extracted from the 1284 MHz MeerKAT observation at MJD 58829. We first formed a timeaveraged and dedispersed version of the observation, keeping the initial frequency resolution of 928 channels over the total bandwidth of 776 MHz. We used the “Bayesian LombScargle Periodogram” (BGLSP) method of RMcalc^{2} to fit for the pulsar’s rotation measure, RM, which quantifies the amount of Faraday rotation occurring along the line of sight caused by the presence of magnetic field in the interstellar medium. Denoting λ as the radio wavelength, we have:
Our BGLSP analysis yielded an RM value of −19 (1) rad m^{−1}, consistent with the value published in Kramer et al. (2021a). We finally corrected the values of Ψ across the frequency bandwidth of the MeerKAT observation for the abovedescribed Faraday rotation term and built a frequencyaveraged pulse profile for PSR J1528−3146 with the four Stokes parameters. Using this profile, we computed the linearly and circularly polarized intensity as shown as red and blue lines in Fig. 3a (top panel), together with the total intensity as a black line. We then performed several RVM fits to the derived PPA data, shown in the bottom panel of Fig. 3a. The fits were obtained using the UltraNest sampler by Buchner (2021), using different priors for the angles in different runs. In each case, we used a standard expression, ln(ℒ) = − χ^{2}/2, as a likelihood function. In the first fit, we assumed our random α and ζ angles to form spin and magnetic field vectors that are distributed uniformly over a unitsphere, while Δ = 0. The prior for Ψ_{0} was a uniform distribution over [ − 90° ,90° ] and, similarly, ϕ_{0} had a uniform prior over [0, 1]. A corner plot of the obtained posterior distributions of the joint RVM model parameters, with the offdiagonal elements representing the correlations between parameters and the diagonal elements denoting the marginalized histograms is shown in Fig. 3b. The PPA values of the discovered interpulse component significantly help to constrain the derived parameter range. In addition to a solution at (α, ζ) = (∼103°, ∼92°), a second solution (less favoured by the loglikelihood statistic) is seen for (α, ζ) = (∼80°, ∼72°). Either case suggests an orthogonal rotator, as expected when an interpulse is present. The solid black line in the lower panel of Fig. 3a is the derived RVM corresponding to the preferred solution (with the gray band indicating 68% uncertainties). The dashed black lines show the RVM shifted by ±90°. If this fit is accurate, it indeed suggests that the main pulse and interpulse are emitting in different orthogonal emission states, as is found quite frequently in interpulsepulsars (see e.g., Johnston & Kramer 2019).
Fig. 3. Rotating vector model analysis of PSR J1528−3146. Left: top panel shows the polarization profile of PSR J1528−3146 as observed at 1284 MHz. The black line represents total power, the red shows linearly and the blue line circularly polarized intensity, respectively. Bottom panel displays the derived PPA values and a number of RVM fits. The black lines represents the resulting RVM for an unconstrained fit covering the whole parameter space. The dashed lines indicate vertical shifts of the RVM by ±90° to represent possible orthogonal polarized modes. The red and blue curves are RVM fits for a restricted range of ζ, forcing it to be consistent with values of the orbital inclination angle as derived from timing in this paper. The magenta line displays a possible shift in the emission height between the main pulse and the interpulse. See text for details. Right: corner plot of the posterior distributions of the joint RVM model parameters, with the offdiagonal elements representing the correlations between parameters and the diagonal elements denoting the marginalized histograms. The posteriors shown here are derived for the unconstrained fit leading to the black line on the left. 
As we discuss further in Sect. 4.3, the Bayesian analysis of the timing data suggests an orbital inclination angle of either i ∼ 57° or i ∼ 123°. For a fully recycled pulsar, we can expect the spin axis of the pulsar to be aligned with the orbital angular momentum vector, so that i ≈ ζ, and also Ψ_{0} ≈ Ω, where Ω is the longitude of the ascending node (see Sect. 4.3). The usage of polarization data and RVM fits to determine the orbital configuration was discussed recently, for instance, by Kramer et al. (2021a) or Guo et al. (2021). Clearly, the above ζvalues are not consistent with the results from timing. In a second and third fit, we therefore applied uniform priors for ζ restricted to a 6° interval centred on the respective inclination angle values. The red curves display the results for ζ ∼ 57° (α = 69 (1)°, ζ = 59 (1)°, Ψ_{0} = −23 (1)°), while the blue lines are for ζ ∼ 123° (α = 128 (1)°, ζ = 120 (1)°, and Ψ_{0} = −10 (1)°), respectively. The latter case, results in an “outer” line of sight, whereas the former is an “inner” one (see e.g., Lorimer & Kramer 2012). It is clear that the red curve, ζ ∼ 57°, is a better fit to the data, but it is still far from ideal.
We can improve the fit by allowing for Δ ≠ 0. This was done in our last fit, the result of which is shown as the magenta line in the PPA panel, resulting in Δ = 0.07 (2) and a significant difference in emission height: (α = 67 (2)°, ζ = 59 (2)°, Ψ_{0} = −25 (1)°). Regardless of a nonzero Δ or otherwise, overall, the solution for ζ ∼ 57° describes the data much better than the solution with i = 123°. It is interesting to also study the values for Ω derived from timing (see Sect. 4 and Table 2) and listed in Table 3. Our value derived for an RVM with ζ close to 57°, Ψ_{0} = −23° is consistent with a value for Ω in the range derived for Solution 4 in Table 3. Nevertheless, the overall quality of the fit when forcing ζ to lie close to the derived orbital inclination angle is clearly limited. An unconstrained fit prefers a value much closer to 90°. As shown in Sect. 4.3, an orbital inclination value close to 90° is excluded by the timing analysis of PSR J1528−3146. Whether the applicability of the RVM is at fault here or whether the assumption of aligned spinvectors is wrong, these are interesting questions that merit further investigation and discussion.
Parameters for PSR J1528−3146 from the timing analysis, and derived quantities.
Details of the parameters used in the χ^{2} sampling.
4. Timing analysis
The timing analysis was carried out by analyzing the topocentric TOA data presented in Sect. 2, using the TEMPO2 pulsar timing package (Hobbs et al. 2006). The measured TOAs were first converted to Barycentric Coordinate Time (TCB), taking the known clock corrections for the different telescopes and using the JPL DE438 solar system ephemeris^{3}. The orbital motion of the pulsar was described using the DDH orbital model, which uses the orthometric parameterization of the Shapiro delay (Freire & Wex 2010). This parameterization is wellsuited for systems with low orbital inclinations (which is the case of PSR J1528−3146, as shown later in this section). The reason is that it describes the Shapiro delay using two parameters: the orthometric amplitude (h_{3}) and the orthometric ratio (ς), which are, in such systems, much less correlated than the range (r) and shape (s) parameters of the DD orbital model (Damour & Deruelle 1985, 1986).
We used TEMPO2 to fit for the pulsar’s sky position and proper motion, its spin parameters (spin frequency and first and second time derivatives), orbital parameters, and PK parameters describing perturbations to the classical Keplerian orbit. The TOAs from the different datasets were aligned with each other using time offsets, which were fitted for as part of the timing analysis. The dispersion measure (DM) and its variations were modeled as a piecewiseconstant function, using the “DMX” model. The TOA dataset was divided into segments of 100 days, and constant DM values were fitted in each of the intervals containing timing data. In a first iteration of the analysis, we obtained a preliminary timing model for PSR J1528−3146 which we then used to determine “EFAC” factors for each of the individual datasets, namely, scaling factors applied to the TOA uncertainties such that the individual reduced χ^{2} values become equal to unity. The EFAC factors we determined through this procedure are listed in Table 1. As can be seen from the table, the EFACs we obtained were all close to 1. A final iteration of the timing analysis with scaled TOA uncertainties yielded the parameters for PSR J1528−3146 listed in Table 2. The differences between measured TOAs and those predicted by the bestfit timing model (the socalled timing residuals) are presented in Figs. 4 and 5. The residuals are normallydistributed around 0, demonstrating that the parameters in Table 2 describe the timing data appropriately. Figure 6 shows the bestfit piecewiseconsant DM values obtained from the timing analysis. The DM displays longterm variations around the value of 18.163 (6) measured by Jacoby et al. (2007).
Fig. 4. Timing residuals as a function of time (top) and orbital phase (bottom), with the datasets considered in our analysis shifted by different amounts for ease of comparison. In the bottom plot, the orbital phase corresponding to superior conjunction (i.e., the moment of the orbit when the pulsar is on the opposite side of the companion) was shifted to occur at 0.25. 
Fig. 5. Timing residuals as a function of time (top) and orbital phase (bottom), with no time offsets applied to the timing residuals for the different TOA datasets. No trends are visible in the timing residuals plotted as a function of orbital phase or as a function of time, indicating that the timing model describes the data appropriately. 
Fig. 6. Bestfit piecewiseconstant DM values from the timing analysis (see Sect. 4). We used segments of 100 days and excluded time intervals that contained no timing data. 
The analysis of ∼16 yr of timing data, which includes data of unprecedented quality from the MeerKAT telescope, enabled us to measure the pulsar’s proper motion and several PK parameters for the first time. In the following subsections, we present these measurements in detail and discuss their implications with respect to PSR J1528−3146.
4.1. Proper motion, distance and spin parameters
The timing analysis yielded a precise proper motion measurement. We find the total transverse proper motion to be 3.9 (2) mas yr^{−1}, with a position angle of Θ_{μ} = 200.8 (1.8)°. Given that 0° indicates north in the “observer’s” convention we used, the position angle of the proper motion of PSR J1528−3146 indicates a southward motion. Our analysis did not enable us to measure the timing parallax, which otherwise provides a direct constraint on the distance. However, the NE2001 (Cordes & Lazio 2002) and YMW16 (Yao et al. 2017) Galactic freeelectron density models both place the pulsar at a distance of d = 0.8 (2) kpc, assuming conservative uncertainty levels of 20%; we use this distance in the estimates given below.
With these proper motion and distance values, we can estimate the intrinsic spindown rate of the pulsar, Ṗ_{int}, as:
where P is the rotational period of the pulsar and Ṗ is its apparent spindown rate listed in Table 2, μ_{T} is the total transverse proper motion of the pulsar, d its distance, c is the speed of light, and a represents the difference of the accelerations of the pulsar and the Solar System in the gravitational field of the Galaxy, projected along the line of sight from the Solar System to the pulsar.
The term represents the contribution from the Shlovskii effect (Shklovskii 1970) and it is relatively small for PSR J1528−3146 due to its relatively small distance and transverse proper motion. The last term, a/c, can be calculated using the GalPot model of gravitational potential of the Galaxy (McMillan 2017; Dehnen & Binney 1998); it is also small given the small distance of the pulsar. The latter two terms combined represent a total kinematic contribution to the spin derivative ΔṖ = 1.2(3) × 10^{−21} s s^{−1} that is small compared to the apparent spindown rate Ṗ = 2.4850 (3) × 10^{−19} s s^{−1}, so that our estimate for the pulsar’s intrinsic spin period derivative Ṗ_{int} = Ṗ − ΔṖ is very close to Ṗ (see Table 2).
The derived intrinsic spindown rate was then used to determine the pulsar’s spindown power, Ė, magnetic field at the stellar surface, B_{surf}, and characteristic age, τ_{c}. Although PSR J1528−3146 is located at a small distance to the Earth, its spindown power, Ė ∼ 4.3 × 10^{31} erg s^{−1}, is very low compared to that of known γray pulsars (see e.g., Abdo et al. 2013), and no detection of γray pulsations has been reported to date for this pulsar. The nearest object in the Fermi Large Area Telescope (LAT) 12Year Point Source Catalog (4FGLDR3, Abdollahi et al. 2022) is located 45′ away from the pulsar, meaning that PSR J1528−3146 does not, in fact, have any counterpart in the Fermi LAT catalog of γray sources.
4.2. Keplerian and postkeplerian parameters
4.2.1. Rate of advance of periastron
From our timing analysis, we measured yr^{−1}. In the absence of a third object, the observed advance of periastron is given by Lorimer & Kramer (2012):
where the first term is caused by the leading order relativistic effect, which, assuming GR, only depends on the total mass of the system, M_{tot} (Robertson 1938; Taylor & Weisberg 1982):
where s is an exact quantity, the solar mass parameter (, Prša et al. 2016) in time units. Assuming , we find a total mass of M_{tot} = 2.8 (2) M_{⊙}. This total mass constraint is displayed graphically in Fig. 7 (red lines) and forms the basis of the comparison with the estimated mass from the Bayesian analysis in Sect. 4.3.
Fig. 7. Constraints on the mass of the companion of PSR J1528−3146 as a function of the cosine of the orbital inclination (left), and as a function of the mass of the pulsar (right). In the lefthand plot, the gray region is excluded by the requirement that the pulsar mass must be greater than zero, and the narrow orange region near cos i = 0 is excluded by the marginally significant constraints on ẋ. On the righthand side of the figure, mass values in the gray region are excluded by the mass function measurement. The solid red lines indicate the contours of the 1σ uncertainty regions for , the 1σ uncertainty regions for h_{3} (resp., ς) are displayed as solid blue (resp., dashed blue) lines; these parameters are measured in the DDH orbital model. The solid contours include 1σ and 2σ equivalent percentiles for the joint posterior pdfs for cos i and M_{c} (left) and M_{c} − M_{c} (right). The marginalized posterior pdfs for cos i, M_{c} and M_{p} are displayed for each plot axis. In the righthand panel, the probability density increases from right to left. 
We now go on to evaluate the above assumption that . The second term in Eq. (4) is caused by kinematic effects and is given by Kopeikin (1996):
which is at most 1.2 × 10^{−6} ° yr^{−1}, given the constraints on i and Ω from the Bayesian analysis presented in Sect. 4.3. This is four orders of magnitude smaller than the observed advance of periastron. The last term, , comes from the spinorbit coupling. This term is currently only measurable in the case of the double pulsar, PSR J0737−3039A/B (Kramer et al. 2021b). Given PSR J1528−3146’s much wider orbit compared to those of PSRs J0737−3039A and B, this term about two orders of magnitude smaller. Therefore, the initial assumption is warranted.
4.2.2. Shapiro delay
The DDH solution uses the orthometric parameterization of Freire & Wex (2010) to describe the Shapiro delay. The orthometric amplitude (h_{3}) and the orthometric ratio (ς) parameters can be written as a function of the Shapiro range (r) and shape (s) of the DD model, as follows:
where . From our timing, we obtain significant measurements of both parameters: h_{3} = 1.2 (1) μs and ς = 0.44 (8). These values can thus be used to constrain the individual masses of the components of the system. In Fig. 7, these constraints are displayed as blue solid (respectively, dashed) lines for h_{3} (resp., ς).
4.2.3. Variation of the orbital period
Our timing analysis of PSR J1528−3146 using the DDH model did not yield a significant measurement of the variation of the system’s orbital period, Ṗ_{b}. The timing model presented in Table 2 therefore does not include a measurement of Ṗ_{b}. Fitting for this parameter leads to Ṗ_{b,obs} = −0.7 (2.2) × 10^{−14} s s^{−1}, which is consistent with zero given the present uncertainties.
The observable orbital period change can be written as follows (Lorimer & Kramer 2012):
where the first term corresponds to the intrinsic orbital period variation caused by gravitational wave emission. Assuming the masses for the pulsar and its companion derived from the Bayesian analysis presented in Sect. 4.3, we estimate this term to be about −4.8 × 10^{−15} s s^{−1} (Peters 1964), namely, an order of magnitude smaller than the uncertainties.
The second term is the orbital period change caused by mass loss of the system:
where G denotes Newton’s gravitational constant, and I ∼ 10^{45} g cm^{2} is the moment of inertia of the pulsar. Assuming M_{tot} = 2.97 M_{⊙} (see Sect. 4.3), we find that the expected contribution from mass loss to the orbital period change is four orders of magnitude smaller than Ṗ_{b,obs}, and is thus negligible. Similarly, the third term in Eq. (8), which is due to tidal torques, can be neglected since the companion star is not deformed and is not interacting with the pulsar.
Finally, the last term corresponds to the kinematic contribution, that includes the Shklosvkii effect and the Galactic differential acceleration. We estimate it to be Ṗ_{b,kin} = 5.4 (1.4) × 10^{−15} s s^{−1}. This is of the same order as Ṗ_{b,GW}, so that the two terms dominate in the calculation of Ṗ_{b,obs}. However, they mostly cancel each other: we find a prediction for Ṗ_{b} of ∼0.7 (1.4)×10^{−15} s s^{−1}. Our measurement of the apparent orbital period variation is compatible with the predicted Ṗ_{b}, but its uncertainty is still several times larger than that prediction.
4.2.4. Variation of the projected semimajor axis
The projected semimajor axis x and its secular change ẋ have been measured using the DDH model (see Table 2) to be ẋ = 2.6 (1.9) × 10^{−15} lts s^{−1}. Assuming the absence of a third body in the system, the change in x can be written as the sum of various contributions (Lorimer & Kramer 2012):
A detailed evaluation of all these terms was done by Guo et al. (2021) for PSR J2222−0137, whose orbital properties are similar to those of PSR J1528−3146. Guo et al. (2021) concluded that all terms except for the first one are negligible. Similar calculations for PSR J1528−3146 enable us to conclude that the same applies to this system. The leading term in Eq. (10) describes the secular change in ẋ caused by the proper motion of the system, as follows (modified from Kopeikin 1996):
We find that ẋ_{PM} is at most 3.5 × 10^{−15} lts s^{−1}, comparable to ẋ_{obs}. This means that the latter measurement already constrains the orbital orientation of the system, in particular, the value of sin(Θ_{μ} − Ω). These constraints are depicted graphically in Fig. 8.
Fig. 8. Constraints on the longitude of the ascending node, Ω, as a function of the cosine of the orbital inclination. The gray shaded regions are excluded by the mass function and the estimate for the total mass of the binary system. The horizontal dashed orange line indicates the position angle of the proper motion, Θ_{μ}. The constraints given by ẋ_{obs} and its 1σ uncertainty are displayed as solid orange lines, the constraints given by the measurement of ς are given by the blue dashed lines. These parameters are as measured with the DDH model. The grayscale maps represent the joint posterior pdf for cos i and Ω, with darker shades corresponding to higher probabilities. The marginal plots display the marginalized posterior pdfs for cos i and Ω. In the righthand panel, the probability density increases from right to left. 
4.3. Bayesian analysis
In order to estimate the physical parameters of the system (pulsar mass, M_{p}, companion mass, M_{c}, orbital inclination angle, i, and longitude of the ascending node, Ω) while selfconsistently taking all available constraints into account, we performed a Bayesian analysis of the timing data. Following the method described in Stovall et al. (2019), we constructed a χ^{2} map of the threedimensional (3D) parameter space of cos i − Ω − M_{tot}. For randomly aligned orbital angular momenta and systems with unknown total mass, these parameters have a priori uniform probability densities. These uniform distributions are our Bayesian priors.
The orbital model we used to take the kinematic effects into account was the T2 model implemented in TEMPO2 (Edwards et al. 2006), which not only takes into account the secular effects on ẋ and (Kopeikin 1996), but also the effects caused by the varying line of sight to the system caused by the orbital motion of the Earth, such as the annual orbital parallax (Kopeikin 1995). The latter terms are generally very small and are not included in any other other orbital models. In some cases the annual orbital parallax allows for an elimination of the degeneracy between i and 180° −i (see e.g., Stovall et al. 2019). The T2 model uses Ω and i as parameters, calculating all kinematic effects internally from these two numbers, which are specified in the “observer’s convention”, where Ω is measured counterclockwise from north through east and i is the angle between the orbital angular momentum and the line of sight from the pulsar to the Earth.
For each value of the cos i − Ω − M_{tot} grid, we fixed Ω and i and then derived the PK parameters M_{c}, γ and from the total mass, the inclination and the mass function, using the relevant GR equations (from e.g., Lorimer & Kramer 2012). We then fit a solution with these fixed parameters to the timing data, while allowed all the other timing parameters to vary. When TEMPO2 is done minimizing the weighted sum of the squares of all timing residuals, we record the χ^{2} value obtained for that point of the cos i − Ω − M_{tot} grid.
We performed the analysis in two distinct cos i regions instead of a continuous range to reduce the computational cost: from −0.8 to −0.4 for the first region and from 0.4 to 0.8 for the second. We explored total mass M_{tot} values between 2 and 4 M_{⊙}. The ranges and step sizes for the three parameters are given in Table 3. The resulting threedimensional χ^{2} map was then used to derive the 3D probability density map using the Bayesian likelihood function (Splaver et al. 2002):
where is the minimum χ^{2} value over the entire map. From this threedimensional joint posterior probability density function (PDF), we derived 2D joint posterior PDFs for the M_{c} − cos i and M_{p} − M_{c} planes (see main panels in Fig. 7), for the cos i − Ω plane (see the main panel of Fig. 8), and marginalized posterior pdfs for cos i, Ω, M_{tot}, M_{p} and M_{c}; most of these appear in the lateral panels of Figs. 7 and 8.
The marginalized posterior pdf of cos i has two maxima, at −0.55 and 0.55. As can be seen from Table 3, these solutions have virtually the same χ^{2} and therefore the same associated probabilities. The preferred value for cos i of −0.55 corresponds to an inclination of , the other solution corresponds to . The fact that the two maxima have such similar probabilities indicates that the annual orbital parallax is not detected in the considered dataset. As noted in Sect. 3, the polarimetry suggests that the latter may be the preferred solution.
Figure 8 shows four distinct solutions in the cos i − Ω space, which are listed in Table 3. Again, none of the solutions is strongly preferred over the other, as can be seen from Fig. 8. Furthermore, there is a plateau of high probability between these peaks. The highprobability regions are well described by the ẋ and ς from the timing analysis, the broad distribution of probability can thus be understood to result from the relatively low precision of ẋ_{obs}.
For the masses we obtained, from the marginalized posterior pdfs, , , and ; here, the values are the medians of the pdfs and the uncertainties include 68.3% of the total probability around these medians. If we include 95.44% of the probability around the medians, then we obtain , , and . The value of M_{tot} matches that estimated based on the measurement, of 2.8 (2) M_{⊙}.
5. Discussion and prospects
In this work, we describe the polarimetric and timing analysis of our radio observations of PSR J1528−3146. This effort has resulted in substantially improved physical parameters for this system as well as vastly improved knowledge of the radio emission properties for the pulsar.
The mass measurements for the PSR J1528−3146 system are remarkable for several reasons: these values depend mostly on the fact that even with an orbital eccentricity as low as 2.1 × 10^{−4}, we can measure the slow periastron advance of 0.057° yr^{−1}; furthermore, with a low orbital inclination of about 57°, we were able to measure the Shapiro delay with enough precision to be useful, at least in conjunction with . Another feature of the system is that the TOA rms of 9.2 μs represents only 0.00015 of the spin period. The system is thus reminiscent of PSR J2222−0137 (Boyles et al. 2013), but with less pronounced characteristics: the latter pulsar spins faster (P = 32.8 ms), its orbital period is slightly smaller (2.44 days), and its orbit is slightly more eccentric (3.8 × 10^{−4}). On the other hand, PSR J1528−3146 has a magnetic field that is five times larger. The similarities between the two systems do suggest a close evolutionary relation.
For PSR J2222−0137, the high orbital inclination (i = 85.27 (4)°), the better timing precision and the slightly more eccentric, faster orbit allow more precise measurements of the Shapiro delay and . This results in precise mass measurements: M_{p} = 1.831 (10) M_{⊙} and M_{c} = 1.319 (4) M_{⊙}; the high total mass M_{tot} = 3.150 (14) M_{⊙} implies that it is the most massive double degenerate binary known in the Galaxy (Guo et al. 2021). Furthermore, as already highlighted by Cognard et al. (2017), the amount of material accreted by PSR J2222−0137 is very small, which implies that the mass of that pulsar is basically its birth mass. This shows that NSs are born with a wide range of masses.
For the same reasons, the mass of PSR J1528−3146 is also very likely its birth mass – its slower spin and higher magnetic field indicate it has been even less recycled than PSR J2222−0137, and therefore it is likely to have accreted less mass. Thus, the precise measurement of the mass of this pulsar will be important – at the very least, it will add to the distribution of known NS birth masses.
The mass of the companion WD is also interesting. First, it could be the most massive WD around any pulsar measured to date. Second, because, as found by McKee et al. (2020), the masses of WD stars orbiting pulsars do seem to clump around a few preferred values. For He WD companions around MSPs, the measured masses increase with orbital period (Tauris & Savonije 1999), reaching a maximum of about 0.4 M_{⊙}. Then, for several systems with slower spin periods (tens of ms), there is a second clump in companion masses, between 0.7 and 0.9 M_{⊙}, and possibly narrower. Finally, PSR J2222−0137 makes a third clump with companion masses around 1.3 M_{⊙}. In this respect, the current mass measurement of the WD companion to PSR J1528−3146 is already important because it confirms the existence of this third clump, and highlights the general pattern. There are very few WD masses outside these ranges: in two cases, they orbit pulsars that are not recycled, and in very few other cases (e.g., PSR J1614−2230, Crawford et al. 2006; Arzoumanian et al. 2018; Shamohammadi et al. 2023), the system is thought to have evolved through case A Roche lobe overflow (Tauris et al. 2011), which is a likely explanation for their very short spin periods.
In the near future, the precision of and ẋ values for this system will keep improving relatively fast, with the uncertainties decreasing faster than the nominal T^{−3/2} rate (where T is the timing baseline; see Damour & Taylor 1992), because of the recent improvement in timing precision achieved by MeerKAT. A better will improve the mass measurements, and along with better ẋ will greatly improve the constraints on the orbital orientation of the system. The measurement of Ṗ_{b} will improve faster, with uncertainties decreasing faster than the nominal rate of T^{−5/2} – again, because of the timing precision of MeerKAT.
The improved Ṗ_{b} that can be achieved in the near future has the potential for a sensitive search for dipolar gravitational wave (DGW) emission, similar to that done with PSR J2222−0137 (Guo et al. 2021). The nondetection of DGW emission in the latter system was important for ruling out the phenomenon of spontaneous scalarization (Zhao et al. 2022), which was predicted by a class of scalartensor gravity theories of Damour & EspositoFarèse (1993). This happened because of the specific mass of PSR J2222−0137, which is in the socalled “mass gap” of spontaneous scalarization (Shibata et al. 2014), which extends from 1.6 to 1.9 M_{⊙} (Shao et al. 2017). Whether a similar test for PSR J1528−3146 is useful or not will depend on the mass of the pulsar, which is not yet known with sufficient precision: the lower half of the uncertainty region is below the mass gap, any value between the median and the 2σ upper limit (1.89 M_{⊙}) would place the pulsar in that gap. Thus, continued timing of this pulsar will be essential in order to improve the mass measurements and determine its potential as a gravitational laboratory.
Acknowledgments
The MeerKAT telescope is operated by the South African Radio Astronomy Observatory, which is a facility of the National Research Foundation, an agency of the Department of Science and Innovation. SARAO acknowledges the ongoing advice and calibration of GPS systems by the National Metrology Institute of South Africa (NMISA) and the time space reference systems department department of the Paris Observatory. MeerTime data is housed on the OzSTAR supercomputer at Swinburne University of Technology. The Parkes radio telescope is funded by the Commonwealth of Australia for operation as a National Facility managed by CSIRO. We acknowledge the Wiradjuri people as the traditional owners of the Observatory site. This research has made extensive use of NASAs Astrophysics Data System (https://ui.adsabs.harvard.edu/) and includes archived data obtained through the CSIRO Data Access Portal (http://data.csiro.au). Parts of this research were conducted by the Australian Research Council Centre of Excellence for Gravitational Wave Discovery (OzGrav), through project number CE170100004 and the Laureate fellowship number FL150100148. 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 Energies” (PNHE) of CNRS/INSU, France. The authors acknowledge the constant support from the MaxPlanckSociety. M.B. acknowledges support through ARC grant CE170100004 (OzGrav).
References
 Abdo, A. A., Ajello, M., Allafort, A., et al. 2013, ApJS, 208, 17 [Google Scholar]
 Abdollahi, S., Acero, F., Baldini, L., et al. 2022, ApJS, 260, 53 [NASA ADS] [CrossRef] [Google Scholar]
 Alpar, M. A., Cheng, A. F., Ruderman, M. A., & Shaham, J. 1982, Nature, 300, 728 [NASA ADS] [CrossRef] [Google Scholar]
 Arzoumanian, Z., Brazier, A., BurkeSpolaor, S., et al. 2018, ApJS, 235, 37 [CrossRef] [Google Scholar]
 Bailes, M., Jameson, A., Abbate, F., et al. 2020, PASA, 37, e028 [NASA ADS] [CrossRef] [Google Scholar]
 Berezina, M., Champion, D. J., Freire, P. C. C., et al. 2017, MNRAS, 470, 4421 [NASA ADS] [CrossRef] [Google Scholar]
 BisnovatyiKogan, G. S., & Komberg, B. V. 1974, Soviet Astron., 18, 217 [NASA ADS] [Google Scholar]
 Boyles, J., Lynch, R. S., Ransom, S. M., et al. 2013, ApJ, 763, 80 [NASA ADS] [CrossRef] [Google Scholar]
 Buchner, J. 2021, J. Open Source Softw., 6, 3001 [CrossRef] [Google Scholar]
 Camilo, F. 1996, in IAU Colloq. 160: Pulsars: Problems and Progress, eds. S. Johnston, M. A. Walker, & M. Bailes, ASP Conf. Ser., 105, 539 [NASA ADS] [Google Scholar]
 Cognard, I., Freire, P. C. C., Guillemot, L., et al. 2017, ApJ, 844, 128 [NASA ADS] [CrossRef] [Google Scholar]
 Cordes, J. M., & Lazio, T. J. W. 2002, ArXiv eprints [arXiv:astroph/0207156] [Google Scholar]
 Crawford, F., Roberts, M. S. E., Hessels, J. W. T., et al. 2006, ApJ, 652, 1499 [NASA ADS] [CrossRef] [Google Scholar]
 Damour, T., & Deruelle, N. 1985, Ann. Inst. Henri Poincaré Phys. Théor, 43, 107 [Google Scholar]
 Damour, T., & Deruelle, N. 1986, Ann. Inst. Henri Poincaré Phys. Théor, 44, 263 [Google Scholar]
 Damour, T., & EspositoFarèse, G. 1993, Phys. Rev. Lett., 70, 2220 [NASA ADS] [CrossRef] [Google Scholar]
 Damour, T., & Taylor, J. H. 1992, Phys. Rev. D, 45, 1840 [NASA ADS] [CrossRef] [Google Scholar]
 Dehnen, W., & Binney, J. 1998, MNRAS, 294, 429 [Google Scholar]
 Edwards, R. T., Hobbs, G. B., & Manchester, R. N. 2006, MNRAS, 372, 1549 [Google Scholar]
 Freire, P. C. C., & Wex, N. 2010, MNRAS, 409, 199 [NASA ADS] [CrossRef] [Google Scholar]
 Guo, Y. J., Freire, P. C. C., Guillemot, L., et al. 2021, A&A, 654, A16 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Hobbs, G. B., Edwards, R. T., & Manchester, R. N. 2006, MNRAS, 369, 655 [Google Scholar]
 Hobbs, G., Manchester, R. N., Dunning, A., et al. 2020, PASA, 37, e012 [NASA ADS] [CrossRef] [Google Scholar]
 Hotan, A. W., van Straten, W., & Manchester, R. N. 2004, PASA, 21, 302 [Google Scholar]
 Jacoby, B. A., Chakrabarty, D., van Kerkwijk, M. H., Kulkarni, S. R., & Kaplan, D. L. 2006, ApJ, 640, L183 [NASA ADS] [CrossRef] [Google Scholar]
 Jacoby, B. A., Bailes, M., Ord, S. M., Knight, H. S., & Hotan, A. W. 2007, ApJ, 656, 408 [NASA ADS] [CrossRef] [Google Scholar]
 Johnston, S., & Kramer, M. 2019, MNRAS, 490, 4565 [NASA ADS] [CrossRef] [Google Scholar]
 Jonas, J. L. 2009, IEEE Proc., 97, 1522 [Google Scholar]
 Kopeikin, S. M. 1995, ApJ, 439, L5 [NASA ADS] [CrossRef] [Google Scholar]
 Kopeikin, S. M. 1996, ApJ, 467, L93 [NASA ADS] [CrossRef] [Google Scholar]
 Kramer, M., Stairs, I. H., Venkatraman Krishnan, V., et al. 2021a, MNRAS, 504, 2094 [CrossRef] [Google Scholar]
 Kramer, M., Stairs, I. H., Manchester, R. N., et al. 2021b, Phys. Rev. X, 11, 041050 [NASA ADS] [Google Scholar]
 Lorimer, D. R., & Kramer, M. 2012, Handbook of Pulsar Astronomy (Cambridge: Cambridge University Press) [Google Scholar]
 McKee, J. W., Freire, P. C. C., Berezina, M., et al. 2020, MNRAS, 499, 4082 [NASA ADS] [CrossRef] [Google Scholar]
 McMillan, P. J. 2017, MNRAS, 465, 76 [Google Scholar]
 Pennucci, T. T., Demorest, P. B., & Ransom, S. M. 2014, ApJ, 790, 93 [NASA ADS] [CrossRef] [Google Scholar]
 Peters, P. C. 1964, Phys. Rev., 136, 1224 [Google Scholar]
 Prša, A., Harmanec, P., Torres, G., et al. 2016, AJ, 152, 41 [Google Scholar]
 Radhakrishnan, V., & Cooke, D. J. 1969, Astrophys. Lett., 3, 225 [NASA ADS] [Google Scholar]
 Robertson, H. P. 1938, Ann. Math., 39, 101 [Google Scholar]
 Serylak, M., Johnston, S., Kramer, M., et al. 2021, MNRAS, 505, 4483 [NASA ADS] [CrossRef] [Google Scholar]
 Shamohammadi, M., Bailes, M., Freire, P. C. C., et al. 2023, MNRAS, 520, 1789 [NASA ADS] [CrossRef] [Google Scholar]
 Shao, L., Sennett, N., Buonanno, A., Kramer, M., & Wex, N. 2017, Phys. Rev. X, 7, 041025 [Google Scholar]
 Shapiro, I. I. 1964, Phys. Rev. Lett., 13, 789 [Google Scholar]
 Shibata, M., Taniguchi, K., Okawa, H., & Buonanno, A. 2014, Phys. Rev. D, 89, 084005 [NASA ADS] [CrossRef] [Google Scholar]
 Shklovskii, I. S. 1970, Soviet Astron., 13, 562 [NASA ADS] [Google Scholar]
 Splaver, E. M., Nice, D. J., Arzoumanian, Z., et al. 2002, ApJ, 581, 509 [NASA ADS] [CrossRef] [Google Scholar]
 StaveleySmith, L., Wilson, W. E., Bird, T. S., et al. 1996, PASA, 13, 243 [NASA ADS] [CrossRef] [Google Scholar]
 Stovall, K., Freire, P. C. C., Antoniadis, J., et al. 2019, ApJ, 870, 74 [NASA ADS] [CrossRef] [Google Scholar]
 Tauris, T. M., & Savonije, G. J. 1999, A&A, 350, 928 [NASA ADS] [Google Scholar]
 Tauris, T. M., Langer, N., & Kramer, M. 2011, MNRAS, 416, 2130 [NASA ADS] [CrossRef] [Google Scholar]
 Tauris, T. M., Langer, N., & Kramer, M. 2012, MNRAS, 425, 1601 [NASA ADS] [CrossRef] [Google Scholar]
 Taylor, J. H., & Weisberg, J. M. 1982, ApJ, 253, 908 [NASA ADS] [CrossRef] [Google Scholar]
 van Straten, W. 2006, ApJ, 642, 1004 [NASA ADS] [CrossRef] [Google Scholar]
 Yao, J. M., Manchester, R. N., & Wang, N. 2017, ApJ, 835, 29 [NASA ADS] [CrossRef] [Google Scholar]
 Zhao, J., Freire, P. C. C., Kramer, M., Shao, L., & Wex, N. 2022, Class. Quant. Grav., 39, 11LT01 [Google Scholar]
All Tables
Summary of observation parameters and TOA properties for each of the datasets considered in the timing study.
All Figures
Fig. 1. MeerKAT pulse profiles at 1284 MHz for PSR J1528−3146, from a ∼4h observation at MJD 58829, corresponding to the highest S/N observation in our dataset. Left: total intensity as a function of observing frequency and rotational phase. Right: integrated pulse profile, with the main radio pulse centered at rotational phase ∼0.54 and a previouslyunknown weak radio component at phase ∼0.03, indicated with an arrow. The inset at the top left shows a zoomin of the pulse profile around the weak radio component. 

In the text 
Fig. 2. Polar representation of the MeerKAT pulse profile at 1284 MHz displayed in Fig. 1, with radio flux densities shown in logarithmic scale. The black dashed line shows the position of the interpulse component and the red dashed line shows the same position shifted by 180°. 

In the text 
Fig. 3. Rotating vector model analysis of PSR J1528−3146. Left: top panel shows the polarization profile of PSR J1528−3146 as observed at 1284 MHz. The black line represents total power, the red shows linearly and the blue line circularly polarized intensity, respectively. Bottom panel displays the derived PPA values and a number of RVM fits. The black lines represents the resulting RVM for an unconstrained fit covering the whole parameter space. The dashed lines indicate vertical shifts of the RVM by ±90° to represent possible orthogonal polarized modes. The red and blue curves are RVM fits for a restricted range of ζ, forcing it to be consistent with values of the orbital inclination angle as derived from timing in this paper. The magenta line displays a possible shift in the emission height between the main pulse and the interpulse. See text for details. Right: corner plot of the posterior distributions of the joint RVM model parameters, with the offdiagonal elements representing the correlations between parameters and the diagonal elements denoting the marginalized histograms. The posteriors shown here are derived for the unconstrained fit leading to the black line on the left. 

In the text 
Fig. 4. Timing residuals as a function of time (top) and orbital phase (bottom), with the datasets considered in our analysis shifted by different amounts for ease of comparison. In the bottom plot, the orbital phase corresponding to superior conjunction (i.e., the moment of the orbit when the pulsar is on the opposite side of the companion) was shifted to occur at 0.25. 

In the text 
Fig. 5. Timing residuals as a function of time (top) and orbital phase (bottom), with no time offsets applied to the timing residuals for the different TOA datasets. No trends are visible in the timing residuals plotted as a function of orbital phase or as a function of time, indicating that the timing model describes the data appropriately. 

In the text 
Fig. 6. Bestfit piecewiseconstant DM values from the timing analysis (see Sect. 4). We used segments of 100 days and excluded time intervals that contained no timing data. 

In the text 
Fig. 7. Constraints on the mass of the companion of PSR J1528−3146 as a function of the cosine of the orbital inclination (left), and as a function of the mass of the pulsar (right). In the lefthand plot, the gray region is excluded by the requirement that the pulsar mass must be greater than zero, and the narrow orange region near cos i = 0 is excluded by the marginally significant constraints on ẋ. On the righthand side of the figure, mass values in the gray region are excluded by the mass function measurement. The solid red lines indicate the contours of the 1σ uncertainty regions for , the 1σ uncertainty regions for h_{3} (resp., ς) are displayed as solid blue (resp., dashed blue) lines; these parameters are measured in the DDH orbital model. The solid contours include 1σ and 2σ equivalent percentiles for the joint posterior pdfs for cos i and M_{c} (left) and M_{c} − M_{c} (right). The marginalized posterior pdfs for cos i, M_{c} and M_{p} are displayed for each plot axis. In the righthand panel, the probability density increases from right to left. 

In the text 
Fig. 8. Constraints on the longitude of the ascending node, Ω, as a function of the cosine of the orbital inclination. The gray shaded regions are excluded by the mass function and the estimate for the total mass of the binary system. The horizontal dashed orange line indicates the position angle of the proper motion, Θ_{μ}. The constraints given by ẋ_{obs} and its 1σ uncertainty are displayed as solid orange lines, the constraints given by the measurement of ς are given by the blue dashed lines. These parameters are as measured with the DDH model. The grayscale maps represent the joint posterior pdf for cos i and Ω, with darker shades corresponding to higher probabilities. The marginal plots display the marginalized posterior pdfs for cos i and Ω. In the righthand panel, the probability density increases from right to left. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.