# Modelling annual scintillation arc variations in PSR J1643-1224 using the Large European Array for Pulsars

G. Mall, R.A. Main, J. Antoniadis, C.G. Bassa, M. Burgay, S. Chen, Ismaël Cognard, R. Concu, A. Corongiu, M. Gaikwad, et al.

► To cite this version:

G. Mall, R.A. Main, J. Antoniadis, C.G. Bassa, M. Burgay, et al., Modelling annual scintillation arc variations in PSR J1643-1224 using the Large European Array for Pulsars. Monthly Notices of the Royal Astronomical Society, 2022, 511 (1), pp.1104-1114. 10.1093/mnras/stac096. hal-03545020

HAL Id: hal-03545020

https://hal.science/hal-03545020

Submitted on 9 Mar 2023

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

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

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

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# Modelling annual scintillation arc variations in PSR J1643–1224 using the Large European Array for Pulsars

G. Mall,<sup>1,2,3*</sup> R. A. Main,<sup>1</sup> J. Antoniadis,<sup>1,4,5</sup> C. G. Bassa<sup>1</sup>,<sup>6</sup> M. Burgay<sup>1</sup>,<sup>7</sup> S. Chen<sup>1</sup>,<sup>8,9,10</sup> I. Cognard,<sup>8,10</sup> R. Concu,<sup>7</sup> A. Corongiu,<sup>7</sup> M. Gaikwad,<sup>1</sup> H. Hu<sup>1</sup>,<sup>1</sup> G. H. Janssen,<sup>6,11</sup> R. Karuppusamy,<sup>1</sup> M. Kramer<sup>1</sup>,<sup>1,12</sup> K. J. Lee,<sup>13</sup> K. Liu,<sup>1</sup> J. W. McKee<sup>1</sup>,<sup>2</sup> A. Melis,<sup>7</sup> M. B. Mickaliger,<sup>12</sup> D. Perrodin,<sup>7</sup> M. Pilia,<sup>7</sup> A. Possenti,<sup>7,14</sup> D. J. Reardon<sup>1</sup>,<sup>15,16</sup> S. A. Sanidas,<sup>12,17</sup> T. Sprenger<sup>1</sup>,<sup>1</sup> B. W. Stappers<sup>1</sup>,<sup>12</sup> L. Wang<sup>1</sup>,<sup>12,18</sup> O. Wucknitz<sup>1</sup> and W. W. Zhu<sup>18</sup>

Accepted 2022 January 7. Received 2022 January 6; in original form 2021 July 23

## ABSTRACT

In this work, we study variations in the parabolic scintillation arcs of the binary millisecond pulsar PSR J1643–1224 over five years using the Large European Array for Pulsars (LEAP). The two-dimensional (2D) power spectrum of scintillation, called the secondary spectrum, often shows a parabolic distribution of power, where the arc curvature encodes the relative velocities and distances of the pulsar, ionized interstellar medium, and Earth. We observe a clear parabolic scintillation arc, which varies in curvature throughout the year. The distribution of power in the secondary spectra is inconsistent with a single scattering screen, which is fully 1D or entirely isotropic. We fit the observed arc curvature variations with two models: an isotropic scattering screen and a model with two independent 1D screens. We measure the distance to the scattering screen to be in the range 114–223 pc, depending on the model, consistent with the known distance of the foreground large-diameter H II region Sh 2-27 (112 ± 17 pc), suggesting that it is the dominant source of scattering. We obtain only weak constraints on the pulsar’s orbital inclination and longitude of ascending node, since the scintillation pattern is not very sensitive to the pulsar’s motion and the screen is much closer to the Earth than the pulsar. More measurements of this kind – where scattering screens can be associated with foreground objects – will help to inform the origins and distribution of scattering screens within our galaxy.

**Key words:** pulsars: general – pulsars: individual: PSR J1643–1224 – H II regions.

## 1 INTRODUCTION

Pulsars are remarkably stable clocks. This property has proved them to be incredibly successful laboratories for testing the predictions of general relativity using pulsar timing (Taylor, Fowler & McCulloch

1979; Kramer et al. 2006). Pulsar timing arrays (PTAs) aim to detect gravitational waves (GWs) by monitoring many millisecond pulsars (MSPs) over time to measure a spatially correlated signal in their timing residuals (Hellings & Downs 1983). The primary PTAs to date are the European Pulsar Timing Array (EPTA; van Haasteren et al. 2011), which combines data from different European telescopes; the North American Nanohertz Observatory for Gravitational Waves (NANOGrav; Demorest et al. 2013); the Parkes Pulsar Timing Array

© 2022 The Author(s).

Published by Oxford University Press on behalf of Royal Astronomical Society. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted reuse, distribution, and reproduction in any medium, provided the original work is properly cited.

<footnote n="1"><sup>1</sup>Max-Planck-Institut für Radioastronomie, Auf dem Hügel 69, D-53121 Bonn, Germany</footnote>

<footnote n="2"><sup>2</sup>Canadian Institute for Theoretical Astrophysics, University of Toronto, 60 St. George Street, Toronto, ON M5S 3H8, Canada</footnote>

<footnote n="3"><sup>3</sup>Department of Physics, University of Toronto, 60 St. George Street, Toronto, ON M5S 1A7, Canada</footnote>

<footnote n="4"><sup>4</sup>Institute of Astrophysics, FORTH, Department of Physics, University of Crete, Voutes, University Campus, GR-71003 Heraklion, Greece</footnote>

<footnote n="5"><sup>5</sup>Argelander Institut für Astronomie, Auf dem Hügel 71, D-53121 Bonn, Germany</footnote>

<footnote n="6"><sup>6</sup>ASTRON, the Netherlands Institute for Radio Astronomy, Oude Hoogeveensledijk 4, NL-7991 PD Dwingeloo, the Netherlands</footnote>

<footnote n="7"><sup>7</sup>INAF – Osservatorio Astronomico di Cagliari, via della Scienza 5, I-09047 Selargius (CA), Italy</footnote>

<footnote n="8"><sup>8</sup>Station de Radioastronomie de Nançay, Observatoire de Paris, PSL Research University, CNRS/INSU, F-18330 Nançay, France</footnote>

<footnote n="9"><sup>9</sup>FEMTO-ST, Department of Time and Frequency, UBFC and CNRS, F-25030 Besançon, France</footnote>

<footnote n="10"><sup>10</sup>Laboratoire de Physique et Chimie de l'Environnement et de l'Espace LPC2E CNRS-Université d'Orléans, F-45071 Orléans, France</footnote>

<footnote n="11"><sup>11</sup>Department of Astrophysics/IMAPP, Radboud University, P.O. Box 9010, NL-6500 GL Nijmegen, the Netherlands</footnote>

<footnote n="12"><sup>12</sup>Jodrell Bank Centre for Astrophysics, School of Physics and Astronomy, The University of Manchester, Manchester M13 9PL, UK</footnote>

<footnote n="13"><sup>13</sup>Kavli Institute for Astronomy and Astrophysics, Peking University, Beijing 100871, P. R. China</footnote>

<footnote n="14"><sup>14</sup>Dipartimento di Fisica, Università di Cagliari, S.P. Monserrato-Sestu Km 0,700, I-09042 Monserrato (CA), Italy</footnote>

<footnote n="15"><sup>15</sup>Centre for Astrophysics and Supercomputing, Swinburne University of Technology, P.O. Box 218, Hawthorn, VIC 3122, Australia</footnote>

<footnote n="16"><sup>16</sup>Australian Research Council Centre of Excellence for Gravitational Wave Discovery (OzGrav)</footnote>

<footnote n="17"><sup>17</sup>Anton Pannekoek Institute for Astronomy, University of Amsterdam, Science Park 904, NL-1098 XH Amsterdam, the Netherlands</footnote>

<footnote n="18"><sup>18</sup>National Astronomical Observatories, Chinese Academy of Sciences, A20 Datun Rd, Chaoyang District, Beijing 100012, P. R. China</footnote>

<footnote>*E-mail: gmall@mpifr-bonn.mpg.de</footnote>
(PPTA; Manchester et al. 2013); and the International Pulsar Timing Array (IPTA; Hobbs et al. 2010), which is a collaboration between the aforementioned individual PTAs.

Recently, several PTAs reported a detection of a common red-noise process in their 12.5 yr data set (Arzoumanian et al. 2020; Chen et al. 2021; Goncharov et al. 2021), but did not observe the significant spatial correlation needed to claim a GW detection, further explored by Goncharov et al. (2021), who argue that a signal of this type can arise from pulsars with independent red-noise properties. As PTAs may be nearing the detection of a GW signal, we need to better understand all sources of correlated noise in timing residuals; a large contributing factor is the ionized interstellar medium (IISM), which introduces time-varying, chromatic variations in the electron column density and multipath propagation.

The effects of multipath propagation can often be seen through scintillation, a pattern in time and frequency caused by interference between different deflected images of the pulsar. Scintillation is now commonly studied using the secondary spectrum – the two-dimensional (2D) power spectrum of scintillation – where a single dominant scattering screen results in a parabolic distribution of power (Stinebring et al. 2001; Walker et al. 2004; Cordes et al. 2006). While pulsar timing is primarily sensitive to changes in radial motion, the rate of scintillation (or equivalently, the scintillation time-scale) depends on the velocity on the plane of the sky. Measurements of variable scintillation rate can then be used to obtain an additional constraint on the pulsar’s orbit (Lyne 1984; Rickett et al. 2014; Reardon et al. 2019).

Measurements of the scintillation time-scale are dependent on the distribution of power along the scattering screen, while the curvatures of scintillation arcs are far less model-dependent. Annual and orbital variations in arc curvature can be used to measure the screen distance and geometry, and precisely measure orbital inclinations and longitude of ascending node, as shown in 16 years of scintillation arc measurements of PSR J0437-4715 by Reardon et al. (2020).

In this paper, we study PSR J1643-1224, a 4.622 ms period pulsar in a 147 d binary orbit with a white dwarf companion, which is observed as part of all aforementioned PTAs. We summarize the relevant theory of scintillation arcs needed for our paper in Section 2, discuss our observations and data reduction in Section 3, and interpret our arc curvature measurements in Section 4. In Section 5, we describe the models used for the arc curvature variations, and we present our results in Section 6. Finally, we discuss the ramifications of our results in Section 7.

## 2 BACKGROUND ON SCINTILLATION

Pulsar scintillation is caused by deflection of pulsar signals by inhomogeneities in the electron densities in the IISM between the pulsar and observer. These deflections create multiple images, which interfere with each other and produce an interference pattern, which changes with time due to the relative velocities between the pulsar, the IISM, and the Earth. The dynamic spectrum $I(f, t)$ shows the observed intensity as a function of frequency $f$ and time $t$. The squared modulus of the 2D Fourier transform of the dynamic spectrum, $I_2 = |\tilde{I}(f_D, \tau)|^2$, is called the secondary spectrum, where the $\tilde{I}$ denotes a Fourier transform. The secondary spectrum expresses the power as a function of the Doppler rate $f_D$ and geometric time delay $\tau$ between each pair of interfering images (Stinebring et al. 2001; Walker et al. 2004; Cordes et al. 2006).

A notable feature of secondary spectra is parabolic arcs (and sometimes inverted arclets), which imply the presence of a dominant and often anisotropic scattering screen between the pulsar and

observer (Stinebring et al. 2001; Walker et al. 2004; Cordes et al. 2006). The arc curvature $\eta$ of a parabolic arc at a central observing wavelength $\lambda$ is given by

$$ \eta = \frac{d_{\text{eff}} \lambda^2}{2 c \mathbf{v}_{\text{eff}}^2 \cos^2 \alpha}, \quad (1) $$

where $c$ is the speed of light, $\alpha$ is the angle between the anisotropy axis of the screen on the plane of the sky and the effective velocity $\mathbf{v}_{\text{eff}}$, which depends on the velocities of the pulsar $\mathbf{v}_{\text{psr}}$, the IISM $\mathbf{v}_{\text{scr}}$, and the Earth $\mathbf{v}_{\oplus}$ perpendicular to the line of sight

$$ \mathbf{v}_{\text{eff}} = \left(\frac{1}{s} - 1\right) \mathbf{v}_{\text{psr}} + \mathbf{v}_{\oplus} - \frac{1}{s} \mathbf{v}_{\text{scr}}, \quad (2) $$

and the effective distance $d_{\text{eff}}$ is given by

$$ d_{\text{eff}} = \left(\frac{1}{s} - 1\right) d_{\text{psr}}, \quad (3) $$

where $d_{\text{psr}}$ and $d_{\text{scr}}$ are the distance to the pulsar and screen, respectively, and $s \equiv 1 - d_{\text{scr}}/d_{\text{psr}}$.

The variation of the observed arc curvature with time then depends on the distance, geometry, and velocity of the scattering screen, as well as the distance and velocity of the pulsar. These properties will be used in Section 5.

## 3 OBSERVATIONS AND DATA

The data description and reduction in this work are largely the same as in Main et al. (2020); in this section, we reiterate the important points, and specific reduction parameters for PSR J1643-1224.

### 3.1 LEAP data

The Large European Array for Pulsars (LEAP) is an experiment designed to increase the sensitivity of pulsar timing observations, by coherently combining signals of the five largest European telescopes. These telescopes are the Effelsberg Telescope, the Nançay Radio Telescope, the Sardinia Radio Telescope, the Westerbork Synthesis Radio Telescope, and the Lovell Telescope at Jodrell Bank. The data from the five telescopes are coherently added, and the resulting signal-to-noise ratio (S/N) is the linear sum of the S/N from the individual telescopes (Bassa et al. 2016). Combining these dishes results in an effective aperture equivalent a 195 m diameter circular dish.

LEAP has been observing more than 20 MSPs monthly since 2012, at a frequency band centred on 1396 MHz with a bandwidth of 128 MHz, divided into contiguous 16 MHz sub-bands. LEAP observes with whichever telescopes are available, and the baseband data are correlated and coherently added in software at the Jodrell Bank Observatory (Smits et al. 2017). The coherently combined baseband data are stored on magnetic tapes and can be retroactively processed to generate pulse-profile data at arbitrary time and/or frequency resolution. The high sensitivity and flexible data product have allowed LEAP to carry out more than just timing analyses, such as studies of MSP single pulses (Liu et al. 2016; McKee et al. 2019) and scintillation properties (Main et al. 2020). Note that observations in 2012 had short observing times of ~10 min to allow for periodic scans on a phase calibrator, while from 2013 onwards, the typical observing times were extended to 30–60 min. As such, in this work, we select all observations of PSR J1643-1224 from 2013 until 2018, beyond which this pulsar was no longer observed as part of the regular monthly LEAP programme.
## 3.2 Creating dynamic and secondary spectra

For each 16 MHz sub-band, we use DSPSR (van Straten & Bailes 2011) to fold the coherently added baseband data into 10 s time bins, 16 phase bins, and 2048 frequency channels of width 7.8125 kHz. These are combined in frequency using PSRADD from PSRCHIVE (Hotan, van Straten & Manchester 2004; van Straten, Demorest & Oslowski 2012) to create the final folded spectrum. The small channel width of 7.8125 kHz allows us to resolve time delays caused by scintillation up to 64 µs. We sum polarizations to form total intensity, resulting in a 3D data cube with dimensions of time, frequency, and phase.

Before creating the dynamic spectrum, we flag and mask sub-integrations influenced by radio frequency interference, and remove the influence of the bandpass. We sum over time and frequency to form the pulse profile; the bottom half of the profile is selected as the off-pulse region. For every sub-integration, we compute the standard deviation in the off-pulse region, and values more than 5σ greater than the mean of the off-pulse region are masked. The off-pulse region, rather than the full phase window, is used so as to not inadvertently mask bright sub-integrations caused by scintillation maxima. The data cube is then divided by the time average of the off-gates to remove the effects of the bandpass. The cleaned time- and frequency-averaged pulse profile is used as a template. We use the template to weight each phase bin, and then sum over phase to create the dynamic spectrum $I(t, f)$.

The arc curvature changes as a function of frequency (see equation 1); approaches to deal with this include Fourier transforming over $\lambda$ instead of frequency (Reardon et al. 2020) or over a time axis scaled by frequency (Sprenger et al. 2021). In our case, the fractional bandwidth is small, so we compute the secondary spectrum $|\dot{I}(f_D, \tau)|^2$ directly as the squared amplitude of the 2D Fourier transform of the dynamic spectrum $I(t, f)$.

# 4 SECONDARY SPECTRA

## 4.1 Interpretation of the observed arcs

Fig. 1 shows the secondary spectra created for each of our observations, obtained using the methods described in Section 3.3. We see clear scintillation arcs, varying annually. Orbital variation would be clearly observable by comparing arcs on the same month across several years; the 147 d pulsar orbital period is coincidentally close to a 2/5 of a year, effectively causing a 2:5 orbital resonance. After one year, the pulsar will have moved over 2.5 orbits, while after two years, the pulsar will have moved over 5 orbits, and the orbital motion and Earth's motion will then be aligned and anti-aligned on alternating years. However, we do not see a strong effect of the orbital motion on our data, with the arc curvature of a given month appearing similar at every year.

An immediately apparent feature of the secondary spectra is how their distribution of power varies throughout the year. Months September–February show a clear arc indicative of a dominant anisotropic scattering screen. However, months March–July do not show discernible arcs, but rather a more diffuse distribution of power across $f_D$. In these months, it often appears that the power at low and high time delays follows parabola of different curvatures.

The observed secondary spectra cannot be explained either through a single 1D screen or a fully isotropic screen. A 1D screen may explain the clean arcs, but would collapse to a line on $f_D = 0$ when $v_{\text{eff}}^2 \cos^2 \alpha = 0$, while an isotropic screen would likely not show such clear arcs and would have a qualitatively similar distribution

of power when seen at different angles of $v_{\text{eff}}$. To fully explain our observations, we need either a second misaligned screen or for the dominant screen to be elongated with an axial ratio $\gtrsim 2$. We consider both of the possible models in our analysis in the later sections. The possible screen geometries and their resultant secondary spectra are explored further using simple illustrative simulations in Appendix A.

## 4.2 Measuring arc curvatures

To measure arc curvatures, we first average the secondary spectrum to 256 sub-samples of width 0.5 µs along the $\tau$-axis. For each $\tau$ sub-sample, we then fit a double-peaked Gaussian curve as

$$ \log_{10} (|I(f_D, \tau_i)|^2) = a \exp \left( \frac{-(f_D - \mu_{f_D,i})^2}{2\sigma^2} \right) + b \exp \left( \frac{-(f_D + \mu_{f_D,i})^2}{2\sigma^2} \right). \quad (4) $$

While not a physically motivated choice, this approximation is useful for finding the power centroid for each $\tau_i$, used to fit for arc curvatures. We take only values with $a$ or $b$ greater than four times the rms of the background noise for each $I(f_D, \tau_i)$ and remove points that converged to an anomalously large ($\sigma > 5$ mHz) or small ($\sigma < 0.1$ mHz) Gaussian width. We then have a series of independent data points $\tau_i$, and dependent data points $\mu_{f_D,i}$ with uncertainties, which we fit with $|f_D| = \sqrt{\tau/\eta}$. The proportionality constant $1/\sqrt{\eta}$ from the square root fitting is directly proportional to $|v_{\text{eff}}|$.

As described in Section 4.1, at certain months of the year, we see wide and diffuse arcs at low time delays, which cannot be explained by a single 1D or isotropic screen. In such cases, the secondary spectrum may not follow a single parabola. To account for the presence of a second screen, or a secondary screen axis, we restrict our fit to low time delays (taken as $\tau \lesssim 6$ µs). The wide parabola at low time delays represents the screen with the highest projected velocity – either a second scattering screen or the points along the axis of motion in a 2D screen. Attempting to measure a curvature of the points at high time delays is more difficult and may lead to a biased measurement depending on the screen model. An example fit is shown in Fig. 2, showing a case where there is a clear and dominant arc (top panels), and other case where the curvature at low time delays does not match the curvature at high time delays (bottom panels). The transition between the behaviour at low and high values of $\tau$ happens at $\sim 6-9$ µs, motivating our choice of the $\tau$ cutoff in our fitting.

# 5 MODELLING THE VARYING ARC CURVATURE

From equations (1) and (2), the changing velocity from the Earth's and pulsar's orbit results in arc curvature variations. A model of the arc curvature will include properties of the pulsar, specifically the distance and proper motion, which are already constrained through pulsar timing, and unknown values of orbital inclination ($i$) and longitude of ascending node ($\Omega$). We take measurements of the pulsar's distance ($d_{\text{psr}} = 0.85 \pm 0.35$ kpc) from the EPTA timing results of Desvignes et al. (2016). During modelling, we fix proper motions to their timing values, but allow distance to vary, using the timing value as a Gaussian prior.

The distance measurement $d_{\text{psr}} = 0.85 \pm 0.35$ kpc comes directly from the parallax measurement ($\pi = 1.17 \pm 0.26$ mas from Desvignes et al. 2016) and is consistent with the recent PPTA distance (Reardon et al. 2021) of $\pi = 0.82 \pm 0.17$ mas, $d_{\text{psr}} = 1.2_{-0.3}^{+0.4}$ kpc.
<figure><img src="image_3.png" /><figcaption>**Figure 1.** Secondary spectra showing parabolic arcs for all of observations of PSR J1643−1224 used in this paper. Rows correspond to months (Jan–Dec) and columns correspond to years (2013–2018). The blank entries in our data correspond to missing observations. We can see variations in arc curvature by comparing the parabolas. This variation is due to annual and orbital variations from the velocity of Earth and the pulsar, respectively, as described in Section 4.1.</figcaption></figure>
<figure><img src="image_4.png" /><figcaption>**Figure 2.** Dynamic and secondary spectra, with curvature fits for 2014 February (left-hand panel) and 2013 July (right-hand panel) observations. In each panel: Left-hand images: dynamic spectra zoomed into a bandwidth of 4 MHz to discern the fine-scale of scintillation. Dynamic spectra were binned to 30 s for plotting purposes. Top right-hand images: corresponding secondary spectra, with the same scaling as in Fig. 1. Bottom right-hand sub-plots: Fitted points $\mu_{fD,i}$ in slices $\tau_i$, as described in Section 4.2, from which $\eta$ is measured. 2014 February shows a distinct arc well described by one screen and associated diagonal striping in the dynamic spectra. In comparison, 2013 July shows power distributed in a way inconsistent with a single 1D or isotropic screen, as can be seen in the different behaviour at low and high $\tau$.</figcaption></figure>

The quoted distance from Desvignes et al. (2016) is $d_{psf} = 0.76 \pm 0.16$ kpc, which applies the Lutz-Kelker bias correction including a constraint from the pulsar’s luminosity estimate. This reduces the error by more than a factor of 2 compared to the direct parallax distance – we adopt the above value and errors to be conservative. The proper motion values are $\mu_{\alpha} = 6.04 \pm 0.04$ mas yr$^{-1}$, $\mu_{\delta} = 4.07 \pm 0.15$ mas yr$^{-1}$.

For a given screen model, we must also include the screen distance ($d_{scr}$), and parameters related to the screen velocities and geometry, which we explicitly describe in the following sections. The function to compute arc curvatures from the the pulsar’s orbital motion and Earth’s known velocity was taken from SCINTOOLS<fnref n="1" /> described in Reardon et al. (2020), slightly modified to fit for $1/\sqrt{\eta} \propto v_{eff}$ instead of fitting $\eta$ directly, and adding a two-screen model as described below. Our Markov chain Monte Carlo (MCMC) fits are performed using the EMCEE implementation in LMFIT (Newville et al. 2014).

## 5.1 Isotropic model

In an isotropic screen, $\cos^2\alpha = 1$ for any orientation of $v_{eff}$, and the arc curvature depends only on the magnitude of $v_{eff}$ on the plane of the sky. A full model of the screen needs only three parameters, the screen distance, the two components of the 2D screen velocity on the plane of the sky.

As discussed in Section 4.1 and Appendix A, the observed secondary spectra cannot be fully described by a single isotropic screen, needing either multiple screens or an elongated screen with axial ratio $\gtrsim 2$. However, in measuring the arc curvatures at low enough time delays, we are always measuring the magnitude of $v_{eff}$ on the plane of the sky, and the model for the arc curvature variations of any 2D screen is equivalent to the isotropic case. This was one of the primary motivations for restricting time delays to $\tau < 6$ µs, below the visible transition. A fit to the full power distribution of

the secondary spectrum for an elongated 2D screen would require at least two more parameters: the axial ratio and orientation of the screen.

## 5.2 Anisotropic two-screen model

For the case of a 1D screen model, it is only possible to measure $v_{eff}$ parallel to the screen. The observed secondary spectra can be qualitatively produced through the existence of two screens (Section 4.1 and Appendix A). In this model, we have two screens, each with a separate distance $d_{s1}$ and $d_{s2}$, angle $\Psi_1$ and $\Psi_2$, and velocity along the screen’s axis $V_{s1,\Psi_1}$ and $V_{s2,\Psi_2}$.

The measured curvature at low time delays is a measure of the screen with maximum projected effective velocity at any given time (or more precisely, the screen with minimum arc curvature at any given time). For any set of model parameters, there are model predictions of $\eta$ arising from both screens, and our final model that is fit takes the minimum curvature between these models at any given time.

## 5.3 Treatment of uncertainties

In Section 4.2, we described our measurements of the arc curvature $\eta$ and the formal statistical uncertainties $\delta\eta$. However, the formal uncertainties $\delta\eta$ may be underestimated due to unmodelled systematic errors, which could arise from, for example, asymmetric power distribution of the arcs, unresolved arclets, or contribution to the arc curvature from a dimension perpendicular to the primary screen axis. Underestimated errors will lead to biases in the final posterior distributions. To address this issue, we take the approach of using ‘EQUAD’ and ‘EFAC’ values typically used in pulsar timing, which describe the corrected errors as

$$ \delta\eta_{corr} = \sqrt{(EFAC \times \delta\eta)^2 + EQUAD^2}. \qquad (5) $$

In a grid of EQUAD, EFAC values (EQUAD ranging from 0 to 2, EFAC ranging from 0.5 to 5 in 100 steps), we perform a Kolmogorov–

<footnote n="1">¹https://github.com/danielreardon/scintools/</footnote>
<table><caption>**Table 1.** The main table shows the resulting values from the fit, $\chi^2$, degrees of freedom (dof), and Bayesian information criterion (BIC) values for both models. The parameters are defined in Section 5.</caption><thead><tr><th>Parameters</th><th>Isotropic</th><th>Anisotropic</th></tr></thead><tbody><tr><td>d<sub>psr</sub> (kpc)</td><td>1.0 ± 0.3</td><td>1.2 ± 0.3<sup>a</sup></td></tr><tr><td>d<sub>scr</sub> (kpc)</td><td>0.2</td><td>—</td></tr><tr><td>V<sub>scr</sub>, α (km s<sup>-1</sup>)</td><td>7.6 ± 0.8</td><td>—</td></tr><tr><td>V<sub>scr</sub>, δ (km s<sup>-1</sup>)</td><td>-11 ± 1</td><td>—</td></tr><tr><td>d<sub>s1</sub> (kpc)</td><td>—</td><td>0.13</td></tr><tr><td>V<sub>s1</sub>, ψ<sub>1</sub> (km s<sup>-1</sup>)</td><td>—</td><td>10 ± 1</td></tr><tr><td>Ψ<sub>1</sub> (°)</td><td>—</td><td>150 ± 4</td></tr><tr><td>d<sub>s2</sub> (kpc)</td><td>—</td><td>0.3 ± 0.1</td></tr><tr><td>V<sub>s2</sub>, ψ<sub>2</sub> (km s<sup>-1</sup>)</td><td>—</td><td>-6 ± 3</td></tr><tr><td>Ψ<sub>2</sub> (°)</td><td>—</td><td>31 ± 9</td></tr><tr><td>i (°)</td><td>135 ± 19</td><td>95 ± 32<sup>b</sup></td></tr><tr><td>Ω (°)</td><td>278 ± 17</td><td>231 ± 109<sup>b</sup></td></tr><tr><td>χ<sup>2</sup></td><td>45.9</td><td>37.3</td></tr><tr><td>dof</td><td>35</td><td>32</td></tr><tr><td>BIC</td><td>23.2</td><td>25.9</td></tr></tbody></table>

Note. <sup>a</sup>As mentioned in Section 5, we take a Gaussian prior of $d_{psr}=0.85 \pm 0.35$ kpc. <sup>b</sup>For the anisotropic model, the solution for i and $\Omega$ converges to two local minima at $i=68^{\circ} \pm 16^{\circ}$, $\Omega=101^{\circ} \pm 44^{\circ}$ and $i=117^{\circ} \pm 17^{\circ}$, $\Omega=284^{\circ} \pm 43^{\circ}$ (see Fig. B2).

Smirnov (KS) test on the scaled residuals to test how well they are described by a standard normal distribution. We find a maximum value at EFAC = 3.0, EQUAD = 0.36; we adopt these values to correct the errors on $1/\sqrt{\eta}$ before performing MCMC fits.

# 6 DISCUSSION

## 6.1 Results of model fitting

The best-fitting parameters to the variations in arc curvature are listed in Table 1, and the data alongside the model fits are shown in Fig. 3. The bottom panels show the isolated effects of the annual and orbital variation, after subtracting the model velocity of either the pulsar or the Earth from the data. Figs B1 and B2 in Appendix B show the full posterior probability distributions for the two models.

The isotropic and anisotropic models yield $\chi^2$ values of 45.9 and 37.3, respectively, suggesting that the anisotropic model is a better fit to the data. However, the anisotropic model has more free parameters, so we compute the Bayesian information criterion (BIC) for both, resulting in BIC values of 23.2 and 25.9 for the isotropic and anisotropic models, respectively – by this criterion, the isotropic model is slightly preferred to describe the variations in arc curvature.

The scattering screens in both models lie closer to the Earth than the pulsar, so the effective velocity $v_{eff}$ of the system is more sensitive to the motion of the Earth than to that of the pulsar, resulting in relatively poor constraints on orbital parameters. Despite this, the isotropic model clearly gives $i > 90^{\circ}$, resolving the sense of the orbit, while the anisotropic two-screen model finds two local solutions for i and $\Omega$. In addition, the values of i and $\Omega$ are consistent with recent PPTA constraints from the annual-orbital parallax of PSR J1643-1224 (Reardon et al. 2021, fig. 3). As can be seen in Figs B1 and B2, we observe strong degeneracies between pulsar distance $d_{psr}$ and screen parameters, so obtaining accurate measurements is highly dependent on the accuracy of the priors on distance. As there has been slight tension between PTA measurements, with the previous IPTA value appearing lower (Verbiest et al. 2016) and the recent PPTA value appearing higher

(Reardon et al. 2021) than the EPTA value (Desvignes et al. 2016), we adopted a conservative distance prior, as described in Section 5.

## 6.2 Screen association with Sh 2-27

The H II region Sh 2-27, behind which PSR J1643-1224 lies, is ionized by the O-star $\zeta$ Oph at a distance of 112 ± 3 pc (Ocker, Cordes & Chatterjee 2020) and has an inferred diameter of 0.034 kpc assuming spherical symmetry (Harvey-Smith, Madsen & Gaensler 2011). The distance to the screen is then 0.112 ± 0.003 ± 0.017 kpc, with uncertainties of the distance and the H II region's radius.

Sh 2-27 is commonly assumed to be the major contributing region for the scattering of PSR J1643-1224. In addition, several features have been observed in PSR J1643-1224, which may be attributed to propagation. A variation in the flux density of PSR J1643-1224 was observed from 1997 to 2000, interpreted as an extreme scattering event and modelled as arising from a single ionizing cloud (Maitia, Lestrade & Cognard 2003). Between 2010 and 2016, timing residuals for PSR J1643-1224 were measured by Shannon et al. (2016) using the Parkes radio telescope, and in 2015 the pulsar displayed a timing event when a new component of emission suddenly appeared in its pulse profile. The L-band (1–2 GHz) pulse profile for 2015 March showed significant time of arrival (TOA) shifts ($\approx 10$ µs), with larger deviations at 3 GHz, and no observed shift at 600 MHz. Due to the opposite expected scaling of the time delays with frequency, Shannon et al. (2016) concluded this event to be of magnetospheric in origin, intrinsic to the pulsar. However, a following study by NANOGrav suggests that the pulse variation may be caused by propagation, given the qualitative resemblance of the additional profile component with an echo (Brook et al. 2018). We, however, do not notice any significant, qualitative change in the secondary spectra spanning this period (see Fig. 1). A quantitative analysis of the scattering time delays measured for this pulsar from LEAP data will be addressed in a future paper.

In our isotropic model, we measured the distance to the scattering screen to be $d_{scr} = 0.208 \pm 0.015$ kpc, while in our two-screen anisotropic model, we measured the distance to the dominant scattering screen to be $d_{s1} = 0.129 \pm 0.015$ kpc, with the secondary screen at $d_{s2} = 0.34 \pm 0.09$ kpc. Depending on our choice of modelling, the scattering screen is consistent with being entirely within (in the anisotropic case) or near the boundary of (isotropic case) Sh 2-27, strongly suggesting that the scattering resulting in the observed scintillation arcs is associated with Sh 2-27.

Harvey-Smith et al. (2011) estimated the maximum possible mean electron density inside Sh 2-27 to be $n_{e,H II} = 10.1 \pm 1.1 \text{ cm}^{-3}$. Using this value, Ocker et al. (2020) estimated the Dispersion Measure (DM) contribution of Sh 2-27 to the pulsar’s total observed DM to be between DM<sub>H II</sub> = 34.4 ± 4.5 and 56.2 ± 0.4 pc cm<sup>-3</sup>, comprising at least half of the total DM $\approx 62.4$ pc cm<sup>-3</sup> [the range is due to slight tension between the cited IPTA parallax (Verbiest et al. 2016) compared to the EPTA (Desvignes et al. 2016) and PPTA (Reardon et al. 2016) parallax measurements]. Given the large electron density, scattering within Sh 2-27 is almost inevitable, and the observed secondary spectra may even be the result of several scattering screens within the H II region.

We note that other studies have found associations with foreground sources. The most clear associations are for pulsars within supernova remnants – the screen in the Crab nebula has been used to probe pulsar separation of giant pulse emission regions (Main et al. 2021), and Yao et al. (2021) measured the distance between PSR J0538+2817 and the supernova remnant shell using scintillation arcs, suggesting
<figure><img src="image_5.png" /><figcaption>**Figure 3.** Measurements and modelling of the variable arc curvature of PSR J1643−1224. The left-hand plots show the isotropic 2D screen model, and the right-hand plots show the model with two independent 1D screens, as described in Section 5. Top panels show $1/\sqrt{\eta}$ versus time, with data points in black and model fit in red. The bottom sub-plots show the variation of $1/\sqrt{\eta}$ with time of year (left-hand plot) and variation with orbital phase (right-hand plot). The data points and model corresponding to the second screen are shown in grey and blue, respectively.</figcaption></figure>

a 3D spin–velocity alignment of the pulsar. Dexter et al. (2017) used the Very Long Baseline Array (VLBA) + Very Large Array (VLA) to measure angular broadening of several pulsars and were able to associate scattering screens of three sources with foreground H II regions, and three sources within the Carina–Sagittarius spiral arm of the Milky Way. Gupta, Rickett & Lyne (1994) found a persistent slope in the autocorrelation function (ACF) of scintillation of PSR B1642−03, which they postulate arises from the edge of Sh 2-27 (the same H II region studied in this paper). Bhat & Gupta (2002) used scintillation from a sample of pulsars, finding enhanced scattering in many pulsars, likely caused by the Loop I bubble. Reardon et al. (2020) measured the distance to two scattering screens towards PSR J0437−4715, which could be near the edges of the Local Bubble.

## 7 CONCLUSIONS

We presented and modelled five years of variable scintillation arcs of PSR J1643−1224. Depending on the choice of screen model, the distance to the dominant scattering screen $d_{scr}$ is found to be within ~120–200 pc, likely associated with the foreground H II region Sh 2-27. We additionally measure an orbital inclination and longitude of ascending node, which are consistent between our two models, but poorly constrained because the scattering screen is much closer to the Earth than to the pulsar. Generically, if we can associate scattering screens with known astrophysical objects – particularly sources less extended than Sh 2-27 – then scintillation arc modelling may allow for an independent determination of precise pulsar distances.

We find that both an isotropic screen model and the two-screen model can reproduce our observed variable arc curvature, with a slight preference to the isotropic model, given its fewer free parameters. However, the appearance of the power distributions in our secondary spectra suggests that the scattering screen cannot be fully isotropic (see Appendix A for more details). Therefore, to fully describe this scattering screen, we would need at least two more parameters: its degree of anisotropy and its orientation. A full modelling of the observed secondary spectra may distinguish between models, as the data are not uniquely described by a single measurement of an arc curvature. Additionally, sensitive future

observations of PSR J1643−1224 at higher frequencies could reveal multiple, sharper arcs, which could be used to inform the choice of model. Secondary spectra show a well-defined arc at higher frequencies because the thickness of the arc is strongly dependent on frequency (Stinebring, Rickett & Ocker 2019).

If the scattering of a pulsar is dominated by a single, highly anisotropic scattering screen, we may be able to employ the $\theta-\theta$ method, which transforms the secondary spectra variables $\tau$ and $f_D$ to angular coordinates and shows the secondary spectra with parallel linear features, which are more convenient to interpret than parabolic arcs. Sprenger et al. (2021) introduces and describes the transformation, and Baker et al. (2022) show how it can be used to measure precise arc curvatures. However, when multiple or 2D scattering screens are present, as is the case for our observations, there may be biases, as $\theta-\theta$ assumes highly anisotropic scattering screens. The potential biases will first need to be explored in detail using simulations.

Many of our observations of secondary spectra display clear asymmetry in their power distributions. While not investigated here, such behaviour may be related to a potential local gradient in DM in the scattering screen along the direction of $v_{eff}$ (e.g. Rickett et al. 2014). In addition, the secondary spectrum shows the geometric time delay between interfering images, and under some assumptions, it can be used to estimate the total time delay due to multipath propagation, which may improve timing (Hemberger & Stinebring 2008; Main et al. 2020). An investigation of the link between DM, scattering time-scale, and scintillation arc variability for PSR J1643−1224, along with other LEAP sources, will be the focus of future work.

## ACKNOWLEDGEMENTS

We thank the anonymous referee for their comments, which improved this work. GM thanks C. R. H. Walker and Rik van Lieshout for detailed and valuable suggestions. This work was supported by the ERC Advanced Grant ‘LEAP’, Grant Agreement Number 227947 (PI: M. Kramer). JA acknowledges support by the Stavros Niarchos Foundation (SNF) and the Hellenic Foundation for Research and Innovation (H.F.R.I.) under the 2nd Call of ‘Science and Society’
Action Always strive for excellence – ‘Theodoros Papazoglou’ (Project Number: 01431). KL is supported by the European Research Council for the ERC Synergy Grant BlackHoleCam under contract no. 610058. JWM is a CITA Postdoctoral Fellow. This work was supported by the Natural Sciences and Engineering Research Council of Canada (NSERC) (funding reference #CITA 490888-16). TS is a member of the International Max Planck Research School for Astronomy and Astrophysics at the Universities of Bonn and Cologne.

The European Pulsar Timing Array (EPTA) is a collaboration between European Institutes, namely ASTRON (NL), INAF/Osservatorio Astronomico di Cagliari (IT), the Max-Planck-Institut für Radioastronomie (GER), Nançay/Paris Observatory (FRA), the University of Manchester (UK), the University of Birmingham (UK), the University of Cambridge (UK), and the University of Bielefeld (GER), with an aim to provide high-precision pulsar timing to work towards the direct detection of low-frequency gravitational waves. The Effelsberg 100 m telescope is operated by the Max-Planck-Institut für Radioastronomie. Pulsar research at the Jodrell Bank Centre for Astrophysics and the observations using the Lovell Telescope are supported by a consolidated grant from the STFC in the UK. The Westerbork Synthesis Radio Telescope is operated by the Netherlands Foundation for Radio Astronomy, ASTRON, with support from NWO. The Nançay Radio Observatory is operated by the Paris Observatory, associated with the French Centre National de la Recherche Scientifique. The Sardinia Radio Telescope (SRT) is funded by the Department of Universities and Research (MIUR), the Italian Space Agency (ASI), and the Autonomous Region of Sardinia (RAS), and is operated as a National Facility by the National Institute for Astrophysics (INAF).

## DATA AVAILABILITY

The timing data used in this paper shall be shared on reasonable request to the corresponding author. Institute for Radio Astronomy, Bonn.

## REFERENCES

Arzoumanian Z. et al., 2020, ApJ, 905, L34

Baker D., Brisken W., van Kerkwijk M. H., Main R., Pen U.-L., Sprenger T., Wucknitz O., 2022, MNRAS, 510, 4573

Bassa C. G. et al., 2016, MNRAS, 456, 2196

Bhat N. R., Gupta Y., 2002, ApJ, 567, 342

Brook P. R. et al., 2018, ApJ, 868, 122

Chen S. et al., 2021, MNRAS, 508, 4970

Cordes J. M., Rickett B. J., Stinebring D. R., Coles W. A., 2006, ApJ, 637, 346

Demorest P. B. et al., 2013, ApJ, 762, 94

Desvignes G. et al., 2016, MNRAS, 458, 3341

Dexter J. et al., 2017, MNRAS, 471, 3563

Foreman-Mackey D., 2016, Journal of Open Source Software, 1, 24

Goncharov B. et al., 2021, ApJ, 917, L19

Gupta Y., Rickett B. J., Lyne A. G., 1994, MNRAS, 269, 1035

Harvey-Smith L., Madsen G. J., Gaensler B. M., 2011, ApJ, 736, 83

Hellings R. W., Downs G. S., 1983, ApJ, 265, L39

Hemberger D. A., Stinebring D. R., 2008, ApJ, 674, L37

Hobbs G. et al., 2010, Class. Quantum Gravity, 27, 084013

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

Kramer M. et al., 2006, Science, 314, 97

Liu K. et al., 2016, MNRAS, 463, 3239

Lyne A. G., 1984, Nature, 310, 300

McKee J. W. et al., 2019, MNRAS, 483, 4784

Main R. A. et al., 2020, MNRAS, 499, 1468

Main R., Lin R., van Kerkwijk M. H., Pen U.-L., Rudnitskii A. G., Popov M. V., Soglasnov V. A., Lyutikov M., 2021, ApJ, 915, 65

Maitia V., Lestrade J. F., Cognard I., 2003, ApJ, 582, 972

Manchester R. et al., 2013, PASA, 30, e017

Newville M., Stensitzki T., Allen D. B., Ingargiola A., 2014, LMFIT: Non-Linear Least-Square Minimization and Curve-Fitting for Python, https://zenodo.org/record/11813

Ocker S. K., Cordes J. M., Chatterjee S., 2020, ApJ, 897, 124

Reardon D. J. et al., 2016, MNRAS, 455, 1751

Reardon D. J., Coles W. A., Hobbs G., Ord S., Kerr M., Bailes M., Bhat N. D. R., Venkatraman Krishnan V., 2019, MNRAS, 485, 4389

Reardon D. J. et al., 2020, ApJ, 904, 104

Reardon D. J. et al., 2021, MNRAS, 507, 2137

Rickett B. J. et al., 2014, ApJ, 787, 161

Shannon R. M. et al., 2016, ApJ, 828, L1

Smits R. et al., 2017, Astron. Comput., 19, 66

Sprenger T., Wucknitz O., Main R., Baker D., Brisken W., 2021, MNRAS, 500, 1114

Stinebring D., McLaughlin M., Cordes J., Becker K., Goodman J. E., Kramer M., Scheckard J., Smith C., 2001, ApJ, 549, L97

Stinebring D. R., Rickett B. J., Ocker S. K., 2019, ApJ, 870, 82

Taylor J. H., Fowler L., McCulloch P., 1979, Nature, 277, 437

van Haasteren R. et al., 2011, MNRAS, 414, 3117

van Straten W., Bailes M., 2011, PASA, 28, 1

van Straten W., Demorest P., Oslowski S., 2012, Astron. Res. Technol., 9, 237

Verbiest J. et al., 2016, MNRAS, 458, 1267

Walker M. A., Melrose D. B., Stinebring D., Zhang C., 2004, MNRAS, 354, 43

Yao J. et al., 2021, Nat. Astron., 5, 788

## APPENDIX A: SIMULATIONS OF DIFFERENT SCREEN GEOMETRIES

In Section 4.1 and Fig. 1, we note that the power distribution we see cannot be reproduced with a single 1D screen or a perfectly isotropic screen, and we adopt the possibility of multiple screens, or an elongated 2D screen. Here, we performed simulations to illustrate the different models, with the same simulation code used in Baker et al. (2022).

The simulation uses a set of image positions along a thin screen, treating each image as a stationary phase point with a random amplitude and phase at each point. The combination of dispersive and geometric delays remains constant along the screen at a reference frequency. For a grid of time and frequency values, the electric field at the observer is computed as the coherent summation of each point,

$$E(t_i, f_j) = \sum_k \mu_k e^{i\phi_k}. \qquad (A1)$$

The relative geometric delays of the images change over time, producing the time-variable electric field. The dynamic spectrum is then calculated from the amplitude squared of the electric field, and the secondary spectrum as the squared modulus of the 2D Fast fourier transform (FFT) of the dynamic spectrum.

We perform these simulations for the four geometries: a 1D screen, a statistically isotropic screen, a 2D screen with a 2:1 axial ratio, and two misaligned 1D screens at different distances. In each case, we simulate 200 random stationary phase points. In Fig. A1, we show the image distribution of each simulation, and secondary spectra corresponding to two different values of $v_{\text{eff}}$. For the case of two screens, all pairs of images are allowed to interfere, while multiply deflected paths are not considered.

Qualitatively, a 1D screen collapses to the $f_D = 0$ when $v_{\text{eff}} \to 0$, inconsistent with what we see for PSR J1643-1224. The other two
<figure><img src="image_6.png" /><figcaption>**Figure A1.** Simulations of different screen geometries, seen at two different effective velocity vectors. Grey-scale images show the wavefield $|E(f_t, \tau)|^2$, while the coloured images show the secondary spectra $|I(f_t, \tau)|^2$. Top: a perfectly 1D screen, second row: a statistically isotropic 2D screen, third row: a 2D screen with a 2:1 axial ratio, and bottom: two misaligned 1D screens. Details of the simulations are given in Appendix A. Qualitatively, neither the perfectly 1D screen nor the fully isotropic screen fits the distribution of power observed for PSR J1643−1224 throughout the year, with the 1D screen collapsing too thin when $v_{\text{eff}} \to 0$, and the isotropic screen never showing clear arcs. An anisotropic 2D screen or two 1D screens can produce clear primary arcs and a dominant lower curvature at low time delays depending on the orientation of $v_{\text{eff}}$.</figcaption></figure>

models can qualitatively reproduce our results; when the velocity is aligned with the dominant axis (either the elongated axis in 2D or the dominant screen with two 1D screens), we can observe a dominant scintillation arc, which is only slightly smeared. When the velocity is misaligned with the dominant axis, the curvature at low time delays is dominated by the velocity projected along the second axis (proportional to $v_{\text{eff}}$ for a 2D screen, $|v_{\text{eff}}|$ parallel to the second screen for two screens), while the power at high time delays is caused by images along the dominant axis interfering with all other images.

## APPENDIX B: POSTERIOR DISTRIBUTIONS

Figs B1 and B2 show the posterior distributions between different parameters. In these figures, blue lines indicate mean values and contours show $1\sigma$, $2\sigma$, and $3\sigma$ confidence levels. Fig. B1 shows the distributions resulting from the isotropic model discussed in Section 5.1, and Fig. B2 shows the distributions resulting from the anisotropic two-screen model discussed in Section 5.2. In both of these figures, we observe high covariances between pulsar's and screen's distance and velocity.
<figure><img src="image_7.png" /></figure>

**Figure B1.** Posterior probability distributions for the orbital and scintillation parameters of J1643−1224, obtained using an isotropic model fit and an MCMC sampler. The blue lines indicate the mean values and the contours show $1\sigma$, $2\sigma$, and $3\sigma$ confidence levels. The plot shows significant covariances between the parameter pairs ($d_{scr}$, $V_{scr,\alpha}$), ($d_{scr}$, $V_{scr,\delta}$), ($d$, $V_{scr,\delta}$), and ($d$, $d_{scr}$) with correlation coefficients of 0.83, −0.86, −0.78, and 0.76, respectively. This figure was created with corner.py (Foreman-Mackey 2016).
<figure><img src="image_8.png" /><figcaption>**Figure B2.** Posterior probability distributions for the orbital and scintillation parameters of J1643−1224, obtained using an anisotropic two-screen model fit and an MCMC sampler. The blue lines indicate the mean values and the contours show $1\sigma$, $2\sigma$, and $3\sigma$ confidence levels. There are two local maxima in $i$ and $\Omega$, so this model cannot unambiguously determine the sense of orbit. The plot shows significant covariances between the parameter pairs ($V_{s2,\psi_2}$, $\psi_2$) and ($d_{s1}$, $V_{s1,\psi_1}$) with correlation coefficients of 0.91 and 0.79, respectively. This figure was created with corner.py (Foreman-Mackey 2016).</figcaption></figure>

This paper has been typeset from a TeX/LATEX file prepared by the author.