# NICER Detection of Thermal X-Ray Pulsations from the Massive Millisecond Pulsars PSR J0740+6620 and PSR J1614-2230

Michael T. Wolff, Sebastien Guillot, Slavko Bogdanov, Paul S. Ray, Matthew Kerr, Zaven Arzoumanian, Keith C. Gendreau, M. Coleman Miller, Alexander J. Dittmann, Wynn C.G. Ho, et al.

► To cite this version:

Michael T. Wolff, Sebastien Guillot, Slavko Bogdanov, Paul S. Ray, Matthew Kerr, et al.. NICER Detection of Thermal X-Ray Pulsations from the Massive Millisecond Pulsars PSR J0740+6620 and PSR J1614-2230. The Astrophysical journal letters, 2021, 918 (2), pp.L26. 10.3847/2041-8213/ac158e hal-03243446

HAL Id: hal-03243446

https://hal.science/hal-03243446

Submitted on 7 Dec 2021

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 Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# NICER Detection of Thermal X-Ray Pulsations from the Massive Millisecond Pulsars PSR J0740+6620 and PSR J1614-2230

Michael T. Wolff¹, Sebastien Guillot²,³, Slavko Bogdanov⁴, Paul S. Ray¹, Matthew Kerr¹, Zaven Arzoumanian⁵,

Keith C. Gendreau⁵, M. Coleman Miller⁶, Alexander J. Dittmann⁶, Wynn C. G. Ho⁷, Lucas Guillemot⁸,⁹, Ismael Cognard⁸,⁹, Gilles Theureau⁸,⁹,¹⁰, and Kent S. Wood¹¹

¹Space Science Division, U.S. Naval Research Laboratory, Washington, DC 20375-5352, USA; Michael.Wolff@nrl.navy.mil

²IRAP, CNRS, 9 avenue du Colonel Roche, BP 44346, F-31028 Toulouse Cedex 4, France

³Université de Toulouse, CNES, UPS-OMP, F-31028 Toulouse, France

⁴Columbia Astrophysics Laboratory, Columbia University, 550 West 120th Street, New York, NY 10027, USA

⁵Astrophysics Science Division, NASA's Goddard Space Flight Center, Greenbelt, MD 20771, USA

⁶Department of Astronomy and Joint Space-Science Institute, University of Maryland, College Park, MD 20742-2421, USA

⁷Department of Physics and Astronomy, Haverford College, 370 Lancaster Avenue, Haverford, PA 19041, USA

⁸Laboratoire de Physique et Chimie de l'Environnement et de l'Espace, LPC2E, CNRS-Université d'Orléans, F-45071 Orléans, France

⁹Station de Radioastronomie de Nançay, Observatoire de Paris, PLS Research University, CNRS/INSU, F-18330 Nançay, France

¹⁰LUTH, Observatoire de Paris, PSL Research University, CNRS, Université Paris Diderot, Sorbonne Paris Cité, F-92195 Meudon, France

¹¹Technology Service Corporation, 6304 Potomac Avenue, Alexandria, VA 22307, USA

Received 2021 April 16; revised 2021 July 16; accepted 2021 July 17; published 2021 September 8

## Abstract

We report the detection of X-ray pulsations from the rotation-powered millisecond-period pulsars PSR J0740+6620 and PSR J1614-2230, two of the most massive neutron stars known, using observations with the Neutron Star Interior Composition Explorer (NICER). We also analyze X-ray Multi-Mirror Mission (XMM-Newton) data for both pulsars to obtain their time-averaged fluxes and study their respective X-ray fields. PSR J0740+6620 exhibits a broad double-peaked profile with a separation of ~0.4 in phase. PSR J1614-2230, on the other hand, has a broad single-peak profile. We show the NICER detections of X-ray pulsations for both pulsars and also discuss the phase relationship to their radio pulsations. The XMM-Newton X-ray spectra of both pulsars shows they are thermally dominated but in the case of PSR J1614-2230 a weak nonthermal high energy tail appears to be present in the spectrum. The thermally dominated spectra along with broad modulations for both pulsars are indicative of thermal radiation from one or more small regions of the stellar surface. For PSR J0740+6620, this paper documents the data reduction performed to obtain the pulsation detection and prepare for pulse light curve modeling analysis.

*Unified Astronomy Thesaurus concepts: Millisecond pulsars (1062); Neutron stars (1108); Rotation powered pulsars (1408); X-ray binary stars (1811)*

## 1. Introduction

Rotation-powered millisecond pulsars (MSPs) are usually old neutron stars (NSs) with spin periods in the ~1–30 ms range. Their fast rotation is thought to result from the transfer of angular momentum of accreted matter from a close binary companion during a so-called “recycling” phase (Alpar et al. 1982; Radhakrishnan & Srinivasan 1982). MSPs can produce emission from radio to gamma-rays through particle acceleration and radiation processes in the pulsar magnetosphere. However, in the X-ray band substantial emission can originate from hot spots on the NS surface, and such emission can be an important probe of NS properties. This paper will focus on the X-ray emission from two pulsars found to be massive NSs.

Radio pulse timing observations of MSPs have resulted in precisely measured NS masses near or above $2 M_\odot$ for a number of objects to date. Shapiro-delay measurements for PSR J0740+6620, initially reported by Cromartie et al. (2020), yielded $M_{NS} = 2.14_{-0.09}^{+0.10} M_\odot$, and Fonseca et al. (2021) used additional data to infer $M_{NS} = 2.08 \pm 0.07 M_\odot$. For PSR J1614-2230, Demorest et al. (2010) found $M_{NS} = 1.97 \pm 0.04 M_\odot$, later refined

by Arzoumanian et al. (2018) to $M_{NS} = 1.908 \pm 0.016 M_\odot$. A potential third Shapiro-delay pulsar is PSR J1811-2405 for which Ng et al. (2020) found $M_{NS} = 2.0_{-0.5}^{+0.8} M_\odot$. Radial velocity measurements in the optical band have also identified potential high NS mass candidates. For PSR J0348+0432, Antoniadis et al. (2013) found $M_{NS} = 2.01 \pm 0.04 M_\odot$. Romani et al. (2021) measured $M_{NS} = 2.13 \pm 0.04 M_\odot$ for the PSR J1810+1744. Kandel & Romani (2020) reported $M_{NS} = 2.18 \pm 0.09 M_\odot$ for the NS in the black widow system PSR J1959+2048 and $M_{NS} = 2.28_{-0.09}^{+0.10} M_\odot$ for the NS in the redback binary system PSR J2215+5135 (see also Linares et al. 2018), although these last two mass measurements are somewhat model dependent. Considering these various NS mass measurements, Romani et al. (2021) explore the possibility that a substantial fraction of the NS population may exist with high masses near $2.0 M_\odot$ and above. How NSs can acquire such high masses is still an unresolved question, however. Depending on the evolutionary pathway of a given MSP, an appreciable amount of mass can be accreted from the companion star: for example, Tauris et al. (2011) found that PSR J1614-2230 could have accreted between ~0.05 $M_\odot$ and ~0.45 $M_\odot$. On the other hand, Cognard et al. (2017) argued that the moderately massive NS in PSR J2222-0137 ($M_{NS} = 1.75 \pm 0.06 M_\odot$) likely acquired its mass at birth and not via accretion.

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

Original content from this work may be used under the terms of the [Creative Commons Attribution 4.0 licence](https://creativecommons.org/licenses/by/4.0/). Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.
Such massive NSs are of great interest for understanding the NS mass–radius relation and the equation of state (EOS) of cold, dense matter. High-mass NSs enable us to constrain the very uncertain EOS at the highest densities and connect with the EOS at lower densities obtained from less massive NSs. Furthermore, each theoretical model of NS structure predicts a particular maximum mass (Oppenheimer & Volkoff 1939; see also, e.g., Lattimer 2012), above which black hole formation occurs. For example, GW190814 indicates the merger of a black hole and a second compact object with mass $2.6 M_\odot$ (Abbott et al. 2020). If the second object is also a black hole, then the maximum mass of a (slowly rotating) NS is between $2 M_\odot$ and $2.6 M_\odot$, and there is either a small gap or no gap between the highest-mass NS and lowest-mass black hole. On the other hand, if knowledge gained from $2 M_\odot$ NSs indicates an EOS that can support NSs of much higher mass, then objects such as the secondary of GW190814 could be an NS (see Abbott et al. 2018).

Even more stringent constraints on the dense matter EOS can be obtained if the mass measurement is combined with a radius estimate from modeling of the pulsed surface thermal X-ray radiation (see, e.g., Pavlov & Zavlin 1997; Miller & Lamb 1998; Bogdanov et al. 2007; Psaltis and Özel 2014). Obtaining NS radius measurements using this technique is one of the primary scientific goals of the Neutron star Interior Composition Explorer (NICER) mission. The first such constraints were presented in Miller et al. (2019) and Riley et al. (2019) based on a deep NICER exposure of the nearby isolated MSP PSR J0030+0451 (Bogdanov et al. 2019).

NICER observations of massive MSPs have the potential to provide unique constraints on the properties of cold, catalyzed matter beyond nuclear density even when very X-ray faint, as discussed in depth by Miller (2016) for the case of PSR J1614−2230. Because more compact stars tend to have lower fractional modulation amplitudes, all else being equal, measurement of the thermal X-ray pulsation amplitude will provide an upper limit on the compactness $GM/Rc^2$ and thus a lower limit on the radius given the precisely determined mass of this star. A lower limit to the radius of a high-mass NS, if that lower limit exceeds 10–11 km, would provide important input to nuclear theories. An upper limit on the stellar radius cannot be obtained in the same way. For example, moving the spot center toward a rotational pole can reduce the amplitude by an arbitrary factor regardless of the size of the star. A robust minimum radius can be inferred, using assumptions that otherwise maximize the overall amplitude (including the assumptions that the spot is point like and centered on the rotational equator). This is contingent, however, on being able to place constraints on any weak nonthermal pulsar emission that could potentially bias the derived radius constraints (Guillot et al. 2016).

PSR J0740+6620 is an MSP with $P_{\text{spin}} = 2.89$ ms discovered in the course of the Green Bank Northern Celestial Cap pulsar survey (Lynch et al. 2018). It is in a nearly circular 4.8 day binary orbit, possibly with an ultra-cool white dwarf companion (Beronya et al. 2019). Recently, Fonseca et al. (2021) obtained an updated mass measurement of $M_{\text{NS}} = 2.08 \pm 0.07 M_\odot$ (at 68.3% confidence) for PSR J0740+6620, from the observed radio Shapiro delay, making it the NS with the largest precisely established mass (first reported by Cromartie et al. 2020). Fonseca et al. (2021) also reported an updated distance measurement from

timing parallax and the Shklovskii effect of $d = 1.14_{-0.15}^{+0.17}$ kpc (68.3% confidence).

PSR J1614−2230 was discovered in a Parkes radio search targeting unidentified EGRET sources (Hessels et al. 2005; Crawford et al. 2006). It is an MSP with $P_{\text{spin}} = 3.15$ ms bound to a massive white dwarf companion in a $P_b = 8.7$ d orbit. The pulsar is important for the NICER mission because it is one of the most massive NSs known, with the most recent measurement of $M = 1.908 \pm 0.016 M_\odot$ (Arzoumanian et al. 2018). This precise measurement resulted from the detection of a strong signature of Shapiro delay in the binary (Demorest et al. 2010). PSR J1614−2230 has a parallax of $1.5 \pm 0.1$ mas, so it lies at a parallax distance of $d = 670_{-40}^{+50}$ pc (Arzoumanian et al. 2018). X-ray pulsations from this pulsar were first claimed at the ~4σ confidence level based on X-ray Multi-Mirror Mission (XMM-Newton) observations (Pancrazi et al. 2012).

Initial exploratory X-ray observations of both pulsars with NICER, XMM-Newton, and Swift showed that they are very faint sources. We have undertaken a systematic survey with NICER of the sample of nearby rotation-powered MSPs to detect and characterize their pulsed X-ray radiation, with a focus on finding thermal X-ray pulsations that are desirable for constraints on the NS mass–radius relation and the dense matter EOS. The NICER observations of eight other MSPs were presented in Ray et al. (2019) and Guillot et al. (2019).

In the present paper, we describe NICER observations of the two nearby binary MSP systems hosting massive NSs, PSR J0740+6620 and PSR J1614−2230, and report the detection of thermal X-ray pulsations from both. The paper is organized as follows. In Section 2 we describe the X-ray observations and the data reduction procedures. Because the photon event data for PSR J0740+6620 have been made available for the NS radius inference analyses of Miller et al. (2021) and Riley et al. (2021), the data extraction procedure is more stringent for that pulsar than for PSR J1614−2230. In Section 3 we present our findings from the X-ray observations together with the known radio properties of these two pulsars. In Section 4 we show how faint the pulsars are by discussing their XMM-Newton spectra. In Section 5, we discuss issues that are raised either by the observations themselves or the follow-on inference analyses. Finally, in Section 6 we summarize our principal conclusions.

## 2. Observations

The NICER (see Gendreau & Arzoumanian 2017) is a NASA Astrophysics Explorer mission that has been operating on the International Space Station (ISS) since 2017 June. The NICER X-ray Timing Instrument (XTI) is an array of 52 active silicon drift detectors each paired with concentrator optics that are nominally sensitive to photons in the 0.2–12 keV range and have absolute timing precision of ~100 nanoseconds rms for tagging of photon detection times (Prigozhin et al. 2016). The event data for both sources were processed using HEASOFT v6.28 and the NICER-specific package NICERDAS v6, together with the NICER calibration files v20200202. However, somewhat different criteria for extracting NICER event data were used for the two pulsars. Furthermore, the XMM-Newton data were obtained under a number of different circumstances, as described below.
<table><caption>**Table 1**
NICER Observation Summary for PSR J0740+6620</caption><thead><tr><th>ObsID<br>Range</th><th>Dates<br>Range</th><th>Exposure<br>Raw (ks)</th><th>Exposure<br>Filtered (ks)</th></tr></thead><tbody><tr><td>1031020101–1031020128</td><td>2018-09-21–2018-10-29</td><td>322.1</td><td>302.8</td></tr><tr><td>2031020101–2031020257</td><td>2019-04-05–2020-02-29</td><td>1142.2</td><td>983.5</td></tr><tr><td>3031020201–3031020244</td><td>2020-03-01–2020-04-17</td><td>403.2</td><td>360.1</td></tr></tbody></table>

## 2.1. NICER

PSR J0740+6620 was observed by NICER for 1867.5 ks between 2018 September 21 and 2020 April 17, while PSR J1614−2230 was observed by NICER for 852.5 ks over the period spanning 2017 July 4 to 2021 March 1. Both pulsars are challenging to observe because of their faintness, although PSR J1614−2230 is somewhat brighter than PSR J0740+6620. Because both are faint, everything possible has been done to minimize backgrounds, including making observations primarily during ISS orbital night. Special care is warranted for PSR J0740+6620 because data sets are being used to obtain an NS radius constraint from the light curve modeling of the pulsed component (which is derived from NICER) relative to the total phase-averaged flux (which is derived from XMM-Newton imaging). Furthermore, there are nearby sources of comparable brightness that the NICER concentrating detector optics do not exclude and their contributions must be estimated. A summary of NICER observations of our target pulsars can be found in Tables 1 and 2.

For PSR J0740+6620, all events included in our and subsequent analyses were obtained when the pulsar was at least 80° from the Sun. This minimized the effects of scattered solar X-rays and detector noise resulting from photons whose energies are just below the NICER energy range. We extract the event data in an ObsID by ObsID manner. We first filter out event data outside the 0.25–3.0 keV energy range. We also filter the observations keeping only those times that have a cosmic-ray cutoff rigidity greater than 2.0 GeV $c^{-1}$ and a space weather Kp index less than 5.0. We further require that the observations by NICER include the full complement of 52 detectors active during the observations. We exclude the detector with DET_ID 34 since this detector tends to be noisy. This ensures that the event data passed to the light curve analysts have an internally consistent set of 51 out of 52 available detectors and thus utilize a constant effective area configuration.

We also demand that NICER is always pointed to better than 54" of the source in order to include the data. The NICER vignetting function is known to be relatively flat-topped within 2' of the instrument boresight so the pointing accuracy filter removes rare, large pointing errors as a possible significant source of brightness variation. In 2020 August maintenance on the individual detectors was begun and this involved cycling individual detectors off for a period of annealing. Thus, for our initial investigation of PSR J0740+6620, utilizing event data obtained after the detector maintenance was started would have created a NICER XTI effective area inconsistency across this date boundary, and we choose to avoid this complication. The pulsar phase for each event is calculated with the photo-nphase tool of the PINT software package (Luo et al. 2021)

using the timing solutions for each pulsar described in Section 2.3. The end result here is a set of good time intervals (GTIs) with event data and low background noise.

For PSR J1614−2230, because we are not yet at the stage of performing a full light curve analysis to obtain an estimated NS radius, we slightly relaxed our event acceptance criteria. Namely, we employ a Kp index constraint of 4.0 and a Sun angle constraint of 60°. Furthermore, some of the NICER detectors suffered from a brief single event upset (SEU) on 2019 July 8, that caused them to develop anomalous time stamps until those detectors could be reset. We exclude the three days of NICER event data between the SEU occurrence and the reset. Because PSR J1614−2230 is somewhat brighter than PSR J0740+6620 this change in event acceptance criteria was found to be the best trade-off between integration time accumulation and background noise.

The NICER non-imaging field of view (FOV) with a ~3' radius does not resolve X-ray sources near the targets and thus such sources can contribute to the effective background counts for the target source. Therefore, the true X-ray phase-averaged count rate for faint pulsars may be difficult to accurately determine. In order to partially mitigate this issue, we make use of XMM-Newton observations to better determine the X-ray fluxes of the pulsars and investigate the presence of other X-ray sources in the NICER FOV.

## 2.2. XMM-Newton

For PSR J0740+6620, we acquired three separate XMM-Newton data sets on 2019 October 26 (ObsID 0851181601), 2019 October 28 (ObsID 0851181401), and 2019 November 1 (ObsID 0851181501) through the Director's Discretionary Time program. The three EPIC instruments were operated in "Full Frame" mode with the "Thin" optical blocking filters in place. Due to the long instrument event readout times in "Full Frame" mode (73.4 ms for EPIC-pn and 2.6 s for EPIC-MOS), the data cannot be used for a pulse timing analysis. As these observations occurred close to the telescope perigee passage, a large fraction of the exposures were affected by intense particle flaring. These are filtered out using the standard procedure, with special care to exclude times of high background flaring. For the EPIC-MOS cameras, we extracted events with $PATTERN <= 12$, and for the EPIC-pn camera we extracted events with $PATTERN <= 4$. For both instruments we excluded time ranges with significantly higher than average count rates in both the full 0.2–12 keV band and the high energy bands (10–12 keV for EPIC-MOS and 10–15 keV for EPIC-pn). This procedure resulted in the exposure times listed in Table 3.

For PSR J1614−2230, there are four separate observations of the pulsar in the XMM-Newton HEASARC archive. We utilize the archival XMM-Newton observations (ObsID: 0404790101) originally reported by Pancrazi et al. (2012) and shown in Table 4.

## 2.3. Radio Pulsar Timing

To search for pulsations, we require pulsar timing models that provide sufficient precision over the full NICER observation span to assign pulse phases to the NICER data without needing any additional trials over folding periods. These models are provided by ground-based pulsar timing programs at radio frequencies.
<table><caption>**Table 2**
NICER Observation Summary for PSR J1614−2230</caption><thead><tr><th>ObsID<br>Range</th><th>Dates<br>Range</th><th>Exposure<br>Raw (ks)</th><th>Exposure<br>Filtered (ks)</th></tr></thead><tbody><tr><td>0060310101–0060310109</td><td>2017-07-04–2017-07-16</td><td>17.2</td><td>15.2</td></tr><tr><td>1060310101–1060310257</td><td>2017-07-18–2019-03-01</td><td>410.3</td><td>332.3</td></tr><tr><td>2060310201–2060310264</td><td>2019-03-06–2019-09-12</td><td>143.7</td><td>126.8</td></tr><tr><td>3060310201–3060310298</td><td>2020-04-02–2021-03-01</td><td>281.3</td><td>260.0</td></tr></tbody></table>

<table><caption>**Table 3**
XMM-Newton Observation Summary for PSR J0740+6620</caption><thead><tr><th>Instrument</th><th>ObsID</th><th>Date</th><th>Exposure Raw (ks)</th><th>Exposure Filtered (ks)</th></tr></thead><tbody><tr><td rowspan="3">XMM-MOS1</td><td>0851181601</td><td>2019-10-26</td><td>6.55</td><td>6.00</td></tr><tr><td>0851181401</td><td>2019-10-28</td><td>8.69</td><td>7.00</td></tr><tr><td>0851181501</td><td>2019-11-01</td><td>11.45</td><td>8.70</td></tr><tr><td rowspan="3">XMM-MOS2</td><td>0851181601</td><td>2019-10-26</td><td>9.33</td><td>7.90</td></tr><tr><td>0851181401</td><td>2019-10-28</td><td>8.68</td><td>6.20</td></tr><tr><td>0851181501</td><td>2019-11-01</td><td>11.44</td><td>9.00</td></tr><tr><td rowspan="3">XMM-pn</td><td>0851181601</td><td>2019-10-26</td><td>6.42</td><td>4.95</td></tr><tr><td>0851181401</td><td>2019-10-28</td><td>6.00</td><td>3.50</td></tr><tr><td>0851181501</td><td>2019-11-01</td><td>8.42</td><td>1.90</td></tr></tbody></table>

**Note.** All XMM-Newton exposures were obtained in “Full-frame” mode. The “Raw” exposure times report the value of the “LIVETIME” FITS keyword. The “Filtered” exposure times report the sum of the GTIs after the filtering described in Section 2.2.

<table><caption>**Table 4**
XMM-Newton Observation Summary for PSR J1614−2230</caption><thead><tr><th>Instrument</th><th>ObsID</th><th>Date</th><th>Exposure Raw (ks)</th><th>Exposure Filtered (ks)</th></tr></thead><tbody><tr><td>XMM-MOS1</td><td>0404790101</td><td>2007-02-08</td><td>45.6</td><td>35.4</td></tr><tr><td>XMM-MOS2</td><td>0404790101</td><td>2007-02-08</td><td>45.7</td><td>34.7</td></tr><tr><td>XMM-pn</td><td>0404790101</td><td>2007-02-08</td><td>38.7</td><td>24.2</td></tr></tbody></table>

**Note.** All XMM-Newton exposures were obtained in “Full-frame” mode. The “Raw” exposure times report the value of the “LIVETIME” FITS keyword. The “Filtered” exposure times report the sum of the GTIs after the filtering described in Section 2.2.

PSR J0740+6620 is regularly observed as part of the NANOGrav pulsar timing program (Alam et al. 2021) with the Green Bank Telescope (GBT), and also by the Canadian Hydrogen Intensity Mapping Experiment (CHIME). Fonseca et al. (2021) have performed an analysis of these data, extending the work done by Cromartie et al. (2020). Because this covers the entire span of the NICER data, we use that timing model to assign pulse phases for the present work.

For PSR J1614−2230, we constructed a timing model using a merged data set from the NANOGrav 12.5 yr release (Alam et al. 2021) and the Nançay radio telescope (NRT). The NRT data were obtained from MJD 57,424 (2016 February 6) through 58,921 (2020 March 13). The NRT backends and data processing procedure are described in Guillemot et al. (2016). We analyzed this combined data set using TEMPO2 and fit a pulsar timing model containing astrometric, rotational, and binary parameters for the pulsar. The model included a

<table><caption>**Table 5**
Ephemeris of PSR J1614−2230</caption><thead><tr><th>Parameter</th><th>Value</th></tr></thead><tbody><tr><td>R.A. (J2000)</td><td>16:14:36.507995(7)</td></tr><tr><td>decl. (J2000)</td><td>-22:30:31.3329(5)</td></tr><tr><td>Proper motion in R.A. (mas yr<sup>-1</sup>)</td><td>3.87(3)</td></tr><tr><td>Proper motion in Decl. (mas yr<sup>-1</sup>)</td><td>-32.3(2)</td></tr><tr><td>Parallax, π (mas)</td><td>1.58(4)</td></tr><tr><td>Epoch of position (MJD)</td><td>56,823</td></tr><tr><td>Spin frequency, ν (Hz)</td><td>317.3789418919128(3)</td></tr><tr><td>Spin frequency derivative, ˙ν (s<sup>-2</sup>)</td><td>-9.69490(4) × 10<sup>-16</sup></td></tr><tr><td>Epoch of period (MJD)</td><td>56,823</td></tr><tr><td>Dispersion measure (DM) (pc cm<sup>-3</sup>)</td><td>34.48759(3)</td></tr><tr><td>DM derivative, DM1</td><td>-0.000302(6)</td></tr><tr><td>DM second derivative<sup>a</sup>, DM2</td><td>-1.75(4) × 10<sup>-4</sup></td></tr><tr><td>Epoch of DM (MJD)</td><td>56,823</td></tr><tr><td>Solar wind density at 1 AU (NE_SW; cm<sup>-3</sup>)</td><td>7.33(16)</td></tr><tr><td>Binary model</td><td>ELL1</td></tr><tr><td>Binary orbital period, P<sub>b</sub> (days)</td><td>8.6866194173(2)</td></tr><tr><td>First derivative of orbital period, ˙P<sub>b</sub></td><td>1.28(5) × 10<sup>-12</sup></td></tr><tr><td>Projected semimajor axis of orbit, x (lt-s)</td><td>11.29119744(4)</td></tr><tr><td>Epoch of ascending node passage, T<sub>asc</sub> (MJD)</td><td>52,331.17011020(5)</td></tr><tr><td>First Laplace parameter, sin ω</td><td>9.7(6) × 10<sup>-8</sup></td></tr><tr><td>Second Laplace parameter cos ω</td><td>-1.345(4) × 10<sup>-6</sup></td></tr><tr><td>Companion mass, M<sub>c</sub>(M<sub>☉</sub>)</td><td>0.4943(19)</td></tr><tr><td>Sine of orbital inclination, sin i</td><td>0.999896(4)</td></tr><tr><td>Terrestrial time standard (CLK)</td><td>TT(BIPM2017)</td></tr><tr><td>Units of barycentric time</td><td>TDB</td></tr><tr><td>Solar System ephemeris</td><td>DE421</td></tr></tbody></table>

**Notes.** Digits in parentheses represent the 1σ uncertainty on the last quoted digit of a parameter value.

<sup>a</sup> The default definition of DM2 (and higher terms) in TEMPO2 changed in 2020 June. We use the corrected Taylor series representation here.

time-variable dispersion measure (DM) using a second-order polynomial plus a solar wind model. The resulting pulse timing solution is presented in Table 5.

### 3. Pulsation Detection

After standard processing and filtering of the NICER data as described above, we isolate the events from the 0.25–3.0 keV range and search for pulsations. Because for both pulsars we have established radio ephemerides we ignore trials factors and simply assign pulse phases as described above. For PSR J0740+6620 a simple $Z_2^2$ test yields a test statistic of $Z_2^2 = 251.19$, corresponding to a detection significance of $15.35\sigma$ for 574405 events in the 0.31–1.22 keV range. We take advantage of the fact that the filtered NICER count rate (about $0.6–1.0 \text{ s}^{-1}$) is dominated by the background, with the pulsar contributing only a few percent to the total count rate. We can use the total count rate in each GTI as
<figure><img src="image_4.png" /><figcaption>(a)</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>(b)</figcaption></figure>

**Figure 1.** NICER pulse profiles of (a) PSR J0740+6620 in the energy range 0.30–1.18 keV, and of (b) PSR J1614−2230 in the energy range 0.59–2.20 keV, folded with their respective radio spin ephemerides. Two rotational cycles are shown for clarity. The blue dashed lines correspond to the radio profiles obtained with the NRT in the 1.4 GHz band.

a proxy for the effective background count rate and employ the GTI sorting method described in Guillot et al. (2019). The GTI sorting method breaks each original interval into shorter segments and adds each segment into the evaluated time series, making sure that the overall pulsation significance increases with each addition. The lowest count rate segments are included first, and we add in segments that have increasing count rates (thus higher background) so long as the detection significance continues to increase. Once the background noise becomes sufficiently large that adding further GTIs from the original distribution decreases rather than increases the detection significance, we stop the procedure. When we apply this method to PSR J0740+6620 we increase the overall pulsation detection significance to $\sigma = 15.46$ for 521,004 events in the 0.31–1.18 keV range. The GTI selection is particularly effective for faint pulsars with count rates $\sim 0.01 \text{ s}^{-1}$ such as PSR J0740+6620. The event lists and GTIs are then saved in both text and FITS file format so that they can be utilized in subsequent processing, in particular the light curve modeling of Miller et al. (2021) and Riley et al. (2021).

For PSR J1614−2230 we confirm the pulsations originally found by Pancrazi et al. (2012) using XMM-Newton timing-mode data. Searching the originally extracted 734.3 ks of event data for this pulsar in the 0.25–3.00 keV range detects the pulsations with significance $Z_2^2 = 72.89$ giving $\sigma = 7.81$ for 335,614 events in the 0.59–2.20 keV range. When we optimize the GTIs as above, our detection significance rises to $\sigma = 8.31$ but the total time is cut somewhat to 644.9 ks, the number of included events is likewise reduced to 267,122, but the optimal energy range is the same at 0.59–2.20 keV.

Figure 1 shows the folded NICER pulse profiles of PSRs J0740+6620 and PSR J1614−2230 in the energy bands where the maximum pulse detection significance is found as determined by the $Z_2^2$-test (de Jager et al. 1989; de Jager & Büsching 2010), together with their radio pulse profiles aligned in rotational phase according to the absolute pulse arrival times provided by the radio timing model of each pulsar.

## 4. Spectral Analysis

We have extracted and analyzed the XMM-Newton spectra for PSR J0740+6620 and PSR J1614−2230. PSR J1614−2230

is the brighter of the two, and we have a much longer total duration of observation for this source than for PSR J0740+6620, but we have sufficient data to generate analyzable spectra for both pulsars.

First, for PSR J0740+6620, the XMM-Newton data are in three relatively short duration observations, as shown in Table 3. We utilize the *eregionanalyse* utility of the SAS analysis suite to analyze the pulsar point source in the X-ray fields of each ObsID, and we note that any other sources in the field are more than 1' away from the pulsar (see Figure 2). Based on this analysis, we choose a 20'' radius circle (79% energy encirclement fraction for EPIC-MOS, 81% energy encirclement fraction for EPIC-pn) around the pulsar radio position to extract spectral events. An indication of the faintness of the pulsar can be seen by considering the number of extracted events from the longest duration ObsID (0851181501), which yields only 57 and 53 pulsar events in the 0.2–12.0 keV range from the MOS1 and MOS2, respectively, and 25 events in the 0.2–15.0 keV range from EPIC-pn. The total number of events extracted across all three ObsIDs from MOS1, MOS2, and EPIC-pn are 100, 105, and 155, respectively. With such few photon events, it is not practical to produce spectra from each individual ObsID for each individual instrument. Thus we combine (using the FTOOLS utility *addascaspec*) the extracted spectra for each instrument and ObsID so that we generate the spectrum for the pulsar from each instrument averaged across the three ObsIDs. Finally, we group the counts into energy bins with at least 10 counts in each instrument spectrum. Ten counts per spectral bin does not fully justify a Gaussian $\chi^2$ statistical treatment of the likelihood but it does preserve sufficient numbers of channels to give us the approximate spectral shape and flux so we elect to use $\chi^2$ as a goodness of fit statistic. We generated response matrices using the SAS commands *arfgen* and *rmfgen*.

The count rate for the pulsar is very low, such that we only fit a simple spectral model to the data. We employ an absorbed blackbody model with $TBabs \times bbodyrad$ in XSPEC to model the pulsar spectrum in the 0.3–5.0 keV energy range. Allowing the hydrogen column density to vary in the fit results in a poorly constrained value of $N_H$, so we fix $N_H$ at a value based on the observed DM ($DM = 14.96 \text{ pc cm}^{-3}$, Cromartie et al. 2020;
<figure><img src="image_6.png" /><figcaption>**Figure 2.** Left panel: XMM-Newton combined MOS1 and MOS2 0.5–3.0 keV image of PSR J0740+6620 and nearby sources. The green cross shows the radio position of the pulsar and the blue crosses show the positions of the field X-ray sources found by the XMM-Newton Pipeline Processing System (PPS) for ObsID 0851181501 data. The inner and outer white circles are 2' and 4' away from the pulsar where the NICER vignetting function drops quickly from 0.90 to 0.11, respectively. The AGN source SDSS J074115.14+662234.9 is 3.57 to the northeast of PSR J0740+6620. At this position SDSS J074115.14+662234.9 falls at the 21% position on the NICER vignetting function. Right panel: XMM-Newton combined MOS1 and MOS2 0.5–3.0 keV image of PSR J1614–2230 and nearby sources. The green cross shows the radio position of the pulsar and the inner and outer white circles are 2' and 4' away from the pulsar, respectively. The blue crosses are nearby sources from the Chandra source catalog (Evans et al. 2010).</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>**Figure 3.** XMM-Newton spectra of (a) PSR J0740+6620 and (b) PSR J1614–2230. All three instruments are shown in each panel with MOS1 in black, MOS2 in red, and the EPIC-pn in green. Even though we have relatively few events to work with, the thermal character of the emission is evident. For PSRs J0740+6620, we do not find compelling evidence for a nonthermal spectral component above 2 keV (see the text). For PSR J1614–2230, we confirm the finding of Pancrazi et al. (2012) of a nonthermal spectral component above 2 keV. For the PSR J1614–2230 spectrum the dotted lines indicate the best-fit blackbody component at low energies and the best-fit power-law component at higher energies, respectively.</figcaption></figure>

Fonseca et al. 2021) for the pulsar (see He et al. 2013) of 4.5 × 10²⁰ cm⁻². The fitted value of the blackbody temperature is kT = 0.140{0.017}^+0.020_-^0.017keV with χ² = 20.25 for 22 DOF. The PSR J0740+6620 spectra from all three instruments can be seen in Figure 3. The total flux from the model is Fₓ = 1.12^{+0.18}_{-0.53} × 10⁻¹⁴ erg cm⁻² s⁻¹ (90% confidence) in the 0.5–8.0 keV range. Employing this model in WebPIMMS<fnref n="12" /> to convert to NICER count rate we find that NICER should observe 0.0128 s⁻¹ from PSR J0740+6620 in the 0.5–8.0 keV range.

We also investigated the possible presence of a nonthermal X-ray power-law component at energies above the thermal

emission in the PSR J0740+6620 spectrum. We employed a simple absorbed blackbody plus power-law model with TBabs×(bbodyrad+powerlaw). However, we could not constrain either the power-law normalization or the spectral index in an acceptable fit. When the power-law normalization and spectral index are allowed to float free in the fit the power-law component completely overwhelms the blackbody component. When we try imposing a fixed spectral index we find that the model converges to one where the power-law normalization is so small that the nonthermal component contributes negligible flux to the model. If a nonthermal tail to the X-ray spectrum is present these observations are not sensitive enough to detect it.

<footnote n="12">¹² https://heasarc.gsfc.nasa.gov/cgi-bin/Tools/w3pimms/w3pimms.pl</footnote>
For PSR J1614−2230, the pulsar is somewhat brighter than for PSR J0740+6620. We utilize the one ObsID that we have focused on for this pulsar and extract events from all three XMM-Newton instruments from a 17" radius region around the pulsar for both MOS instruments and an 18" region for the EPIC-pn. This yields 175 events (MOS1), 170 events (MOS2) in the range 0.2–12.0 keV, and 462 events (EPIC-pn) in the range 0.2–15.0 keV. The reason we use a somewhat smaller extraction region for this pulsar than for PSR J0740+6620 is that a faint source (CXO J161435.6-223104) lies only ~35" away from the pulsar position and we want to prevent contamination of the pulsar X-ray spectrum. Finally, again we group the counts into energy bins with at least 10 counts in each instrument spectrum. The spectrum of PSR J1614−2230 can be seen in Figure 3. We use the model TBabs×(bbodyrad + powerlaw). We note that the $N_H$ value is again poorly constrained by the fit so we fix it at $N_H = 1.03 \times 10^{21} \text{cm}^{-2}$ as indicated by the pulsar DM ($DM = 34.49 \text{ pc cm}^{-3}$) with the He et al. (2013) conversion. However, the fitted blackbody temperature is somewhat higher than in the previous case at $kT = 0.176_{-0.015}^{+0.016}$. The power-law component is required to generate a good fit to the data, but the power-law index is poorly constrained so we fix $\Gamma$ at 1.6 using as guidance here the detailed X-ray spectral fits of PSR J0437−4715 reported by Guillot et al. (2016). This $\Gamma$ value was found by Guillot et al. for the NuSTAR energy range above the thermal emission and it is consistent with the spectral index ranges for the Guillot et al. best thermal plus power-law models. This $\Gamma$ value is also well within the uncertainties of the spectral index measurement obtained by Pancrazi et al. (2012) for PSR J1614−2230. With the inclusion of the power-law component we obtain a $\chi^2$ of 62.09 for 60 DOF. The total flux from the model is $F_x = 3.02_{-0.61}^{+0.45} \times 10^{-14} \text{ erg cm}^{-2} \text{ s}^{-1}$ in the 0.5–8.0 keV range. Employing the same procedure as before to convert to NICER count rate we find that NICER should observe $0.0218 \text{ s}^{-1}$ from PSR J1614−2230 in the 0.5–8.0 keV range (but see the discussion below).

## 5. Discussion

The NICER timing observations presented here reveal PSR J0740+6620 and PSR J1614−2230 as pulsed millisecond-variable X-ray sources with phase-broadened pulsations. The fact that the pulsations are optimized in the 0.31–1.18 keV energy range for PSR J0740+6620 and 0.59–2.20 keV for PSR J1614−2230 suggest a thermal character for the pulsations. Both NSs have masses near $2 M_\odot$ and thus are important X-ray sources to elucidate the cold matter EOS near the maximum NS mass. PSR J0740+6620 shows a double pulse light curve in the NICER energy range, while PSR J1614−2230 shows one broad pulse. The relative phasing of the radio and X-rays shows that for PSR J0740+6620 the radio pulses are aligned with the thermal X-ray pulses and that the X-ray pulses are separated by roughly 0.4 in phase. In the case of PSR J1614−2230 the dominant radio pulse occurs near the end of the broad single thermal pulse, although there is a very small second pulse in the radio light curve. Both pulsars are established as emitting from radio to γ-ray energies, indicating that nonthermal emission from the magnetosphere is present. We find evidence of weak nonthermal emission in the spectrum of PSR J1614−2230, and the fact that the detection of pulsations extends up to 2.29 keV in the NICER data opens the possibility that, with the accumulation of more data,

pulsations may extend to higher energies. For PSR J0740+6620, the thermal character of the X-ray spectrum, our finding pulsations only up to 1.18 keV, and the general faintness of the source suggests that if any nonthermal pulsations are to be found, more sensitive observations in the X-ray range above 3 keV will be necessary.

The goals of these NICER and XMM-Newton observations are to establish the existence of X-ray pulsations from both pulsars and to estimate the time-averaged flux of each pulsar. The XMM-Newton measurements help constrain the total flux from the NS, enabling a better estimate of the pulsed emission fraction in the NICER data, and thus the estimate of the NS compactness. A correct pulsed NICER flux from the pulsar, if paired with an erroneously high pulse-averaged flux from XMM-Newton, would yield a decreased pulsed fraction and thus an incorrect lower NS radius. Conversely, a spuriously low XMM-Newton flux would yield a larger NS radius. Because of the high NS masses in these systems, even a solid lower limit in the radius of each pulsar would be valuable for constraining the NS EOS because the cores of such massive NSs contain higher-density matter. An upper limit on the compactness, along with the known NS masses, can yield a lower limit on the pulsar radius.

These NICER data sets for PSR J0740+6620 and PSR J1614−2230 are very quiet in that the background count rates are strongly minimized. Also, for PSR J0740+6620 the effective area as a function of time should be very constant because we always utilized the events from the same 51 of 52 active detectors.

The XMM-Newton observations of the field of PSR J0740+6620 enable us to estimate the expected count rates of the likely AGN in the field 3'57 north and east from the pulsar. In fact, investigation of this X-ray source shows that there are really two extragalactic sources within 3''16 of one another (SDSS J074114.62+662235.2 and SDSS J074115.14+662234.9, from the Sloan Digital Sky Survey Data Release 9: Ahn et al. 2012). Because this angular distance is similar to the angular resolution of the EPIC-MOS and EPIC-pn instruments the XMM-Newton observations we obtained cannot easily distinguish if the X-rays are coming from one or both of these extragalactic sources. For SDSS J074115.14+662234.9, we extracted the spectrum and obtained a simple absorbed power-law fit with absorbing column $N_H = 2.3 \times 10^{21} \text{ cm}^2$ and spectral index $\Gamma \sim 2.71$. This gives a NICER count rate of ~0.031 c s⁻¹ for SDSS J074115.14+662234.9 in the 0.5–8.0 keV energy range. However, because this AGN (at 3'57 away) is at the 21% level of the NICER vignetting curve (Okajima et al. 2016) as a function of angular radius from the XTI look direction, the AGN only adds ~0.50 of the pulsar count rate when NICER is pointed directly at PSR J0740+6620 in the 0.5–8.0 keV energy range. Other weak sources in the field will also contribute some counts because the NICER concentrator optics will cause some photons from these sources to fall sufficiently close to the center of the detectors to be included in the background.

Chandra observations of the field of PSR J1614−2230 allowed Evans et al. (2010) to estimate the 0.5–7.0 keV fluxes of field sources near this pulsar. We have already incorporated the field source positions into Figure 2 to show their proximity to the pulsar. By using the fluxes in the Evans et al. catalog with an assumed Crab-like spectral model ($N_H = 3.5 \times 10^{21} \text{ cm}^2$, spectral index $\Gamma = 2.05$) we can estimate the NICER count rate in the 0.5–8.0 keV energy band for all the sources,
including the pulsar. For example, we noted above that there is a faint source roughly ~35" from PSR J1614-2230 that we sought to avoid in our spectral work by using a small extraction radius. Based on our simple model for the field source spectra, this close-by field source will account for ~11% of the total NICER count rate when NICER is pointed right at the pulsar. The NICER vignetting fraction for sources 4' away from the pointing direction falls to ~0.11 (see Okajima et al. 2016). If we use the simple Crab model for the spectra of all field sources within 4' of PSR J1614-2230 we expect to observe a total of ~0.066 counts per second above background, including the pulsar, in the 0.5-8.0 keV energy band. Furthermore, when NICER is pointed at the pulsar, PSR J1614-2230 itself will give about ~31% of the total count rate and the field sources will give the other ~69% of the observed counts. We emphasize, however, that this estimate is very strongly spectral model dependent. Using various values of the interstellar absorption and/or a thermal model for the field source spectra can change the calculated rates from the field sources. Nevertheless, the general result here is that when NICER is pointed at PSR J1614-2230 somewhat less than half of the observed count rate in the 0.5-8.0 keV range should be from the pulsar itself with the rest coming from field X-ray sources.

## 6. Conclusions

The event data from NICER and XMM-Newton we described here for PSR J0740+6620 were used by Miller et al. (2021), Raaijmakers et al. (2021), and Riley et al. (2021) to investigate the radius and EOS of the pulsar in this system. For PSR J0740+6620 the NICER observations have shown X-ray pulsations with two broad thermal X-ray peaks separated in phase by about ~0.4. In the case of PSR J1614-2230, we have confirmed the X-ray pulsations originally reported by Pancrazi et al. (2012) and find one broad thermal X-ray peak in the phase diagram. The spectra of both pulsars are dominated by thermal emission but we also confirm the presence of a nonthermal spectral component found by Pancrazi et al. (2012) for PSR J1614-2230. We show the fields of these pulsars in XMM-Newton images and find that both pulsars are faint, and nearby X-ray sources will contribute significantly to their respective NICER concentrator count rates.

We gratefully acknowledge discussions with James Lattimer, Sharon Morsink, Thomas Riley, and Anna Watts. We also thank an anonymous referee whose comments helped improve the manuscript. This work was supported in part by NASA through the NICER mission and the Astrophysics Explorers Program. NICER science team members at NRL are supported by NASA under Interagency Agreement NNG20OB08A. S.B. was funded in part by NASA grants NNX17AC28G and 80NSSC20K0275. M.C.M. thanks the Radboud Excellence Initiative for supporting his stay at Radboud University. M.C. M. and A.J.D. are also supported in part by NASA ADAP grant 80NSSC21K0649. W.C.G.H. acknowledges support through grant 80NSSC20K0278 from NASA. This research has made use of data and/or software provided by the High Energy Astrophysics Science Archive Research Center (HEASARC), which is a service of the Astrophysics Science Division at NASA/GSFC and the High Energy Astrophysics Division of the Smithsonian Astrophysical Observatory. We acknowledge use of NASA's Astrophysics Data System (ADS) bibliographic services and the ArXiv. The authors thank the

XMM-Newton observatory, specifically the help of Felix Fuerst, for the Director's Discretionary Time observations of PSR J0740+6620.

**Facilities:** NICER, XMM-Newton.

**Software:** astropy (Astropy Collaboration et al. 2013, ascl:1304.002), PINT (ascl:1902.007, https://github.com/nanograv/pint), HEAsoft (ascl:1408.004, https://heasarc.nasa.gov/lheasoft/), DS9 (Joye & Mandel 2003), ascl:0003.002), emcee (ascl:1303.002, https://github.com/dfm/emcee), and the XMM-Newton Scientific Analysis System Scientific Analysis System (SAS, ascl:1404.004, https://www.cosmos.esa.int/web/xmm-newton/sas).

## ORCID iDs

Michael T. Wolff https://orcid.org/0000-0002-4013-5650

Sebastien Guillot https://orcid.org/0000-0002-6449-106X

Slavko Bogdanov https://orcid.org/0000-0002-9870-2742

Paul S. Ray https://orcid.org/0000-0002-5297-5278

Matthew Kerr https://orcid.org/0000-0002-0893-4073

Keith C. Gendreau https://orcid.org/0000-0001-7115-2819

M. Coleman Miller https://orcid.org/0000-0002-2666-728X

Alexander J. Dittmann https://orcid.org/0000-0001-6157-6722

Wynn C. G. Ho https://orcid.org/0000-0002-6089-6836

Lucas Guillemot https://orcid.org/0000-0002-9049-8716

Ismael Cognard https://orcid.org/0000-0002-1775-9692

## References

Abbott, B. P., Abbott, R., Abbott, T. D., et al. 2018, PhRvL, 121, 161101

Abbott, R., Abbott, T. D., Abraham, S., et al. 2020, ApJL, 896, L44

Ahn, C. P., Alexandroff, R., Allende Prieto, C., et al. 2012, ApJS, 203, 21

Alam, M. F., Arzoumanian, Z., Baker, P. T., et al. 2021, ApJS, 252, 4

Alpar, M. A., Cheng, A. F., Ruderman, M. A., & Shaham, J. 1982, Natur, 300, 728

Antoniadis, J., Freire, P. C. C., Wex, N., et al. 2013, Sci, 340, 448

Arzoumanian, Z., Brazier, A., Burke-Spolaor, S., et al. 2018, ApJS, 235, 37

Astropy Collaboration, Robitaille, T. P., & Tollerud, E. J. 2013, A&A, 558, A33

Beronya, D. M., Karpova, A. V., Kirichenko, A. Y., et al. 2019, MNRAS, 485, 3715

Bogdanov, S., Guillot, S., Ray, P., et al. 2019, ApJL, 887, L25

Bogdanov, S., Rybicki, G. B., & Grindlay, J. E. 2007, ApJ, 670, 668

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

Crawford, F., Roberts, M. S. E., Hessels, J. W. T., et al. 2006, ApJ, 652, 1499

Cromartie, H. T., Fonseca, E., Ransom, S. M., et al. 2020, NatAs, 4, 72

de Jager, O. C., & Büsching, I. 2010, A&A, 517, L9

de Jager, O. C., Raubenheimer, B. C., & Swanepoel, J. W. H. 1989, A&A, 221, 180

Demorest, P. B., Pennucci, T., Ransom, S. M., Roberts, M. S. E., & Hessels, J. W. T. 2010, Natur, 467, 1081

Evans, I. N., Primini, F. A., Glotfelty, K. J., et al. 2010, ApJS, 189, 37

Fonseca, E., Cromartie, H. T., Pennucci, T. T., et al. 2021, ApJL, 915, L12

Gendreau, K., & Arzoumanian, Z. 2017, NatAs, 1, 895

Guillemot, L., Smith, D. A., Laffon, H., et al. 2016, A&A, 587, A109

Guillot, S., Kaspi, V. M., Archibald, R. F., et al. 2016, MNRAS, 463, 2612

Guillot, S., Kerr, M., Ray, P. S., et al. 2019, ApJL, 887, L27

He, C., Ng, C.-Y., & Kaspi, V. M. 2013, ApJ, 768, 64

Hessels, J., Ransom, S., Roberts, M., et al. 2005, in ASP Conf. Ser. 328, Binary Radio Pulsars, ed. F. A. Rasio & I. H. Stairs (San Francisco, CA: ASP), 395

Joye, W. A., & Mandel, E. 2003, in ASP Ser. 295, Astronomical Data Analysis Software and Systems XII, ed. H. E. Payne, R. I. Jedrzejewski, & R. N. Hook (San Francisco, CA: ASP), 489

Kandel, D., & Romani, R. W. 2020, ApJ, 892, 101

Lattimer, J. M. 2012, ARNPS, 62, 485

Linares, M., Shahbaz, T., & Casares, J. 2018, ApJ, 859, 54

Luo, J., Ransom, S., Demorest, P., et al. 2021, ApJ, 911, 45

Lynch, R. S., Swiggum, J. K., Kondratiev, V. I., et al. 2018, ApJ, 859, 93

Miller, M. C. 2016, ApJ, 822, 27

Miller, M. C., & Lamb, F. K. 1998, ApJL, 499, L37
Miller, M. C., Lamb, F. K., Dittmann, A. J., et al. 2019, ApJL, 887, L24

Miller, M. C., Lamb, F. K., Dittmann, A. J., et al. 2021, ApJL, in press

Ng, C., Guillemot, L., Freire, P. C. C., et al. 2020, MNRAS, 493, 1261

Okajima, T., Soong, Y., Balsamo, E. R., et al. 2016, Proc. SPIE, 9905, 99054X

Oppenheimer, J. R., & Volkoff, G. M. 1939, PhRv, 55, 374

Pancrazi, B., Webb, N. A., Becker, W., et al. 2012, A&A, 544, A108

Pavlov, G. G., & Zavlin, V. E. 1997, ApJL, 490, L91

Prigozhin, G., Gendreau, K., Doty, J. P., et al. 2016, Proc. SPIE, 9905, 99051I

Psaltis, D., & Özel, F. 2014, ApJ, 792, 87

Raaijmakers, G., Greif, S. K., Riley, T. E., et al. 2021, ApJL, in press

Radhakrishnan, V., & Srinivasan, G. 1982, CSci, 51, 1096

Ray, P. S., Guillot, S., Ransom, S. M., et al. 2019, ApJL, 878, L22

Riley, T. E., Watts, A. L., Bogdanov, S., et al. 2019, ApJL, 887, L21

Riley, T. E., Watts, A. L., Ray, P. S., et al. 2021, ApJL, in press

Romani, R. W., Kandel, D., Filippenko, A. V., Brink, T. G., & Zheng, W. 2021, ApJL, 908, L46

Tauris, T. M., Langer, N., & Kramer, M. 2011, MNRAS, 416, 2130