# X-ray view of the 2021 outburst of SXP 15.6: Constraints on the binary orbit and magnetic field of the neutron star

Georgios Vasilopoulos, Gaurav K. Jaisawal, Chandreyee Maitra, Frank Haberl, Pierre Maggi, Apostolos Sotirios Karaferias

► To cite this version:

Georgios Vasilopoulos, Gaurav K. Jaisawal, Chandreyee Maitra, Frank Haberl, Pierre Maggi, et al. X-ray view of the 2021 outburst of SXP 15.6: Constraints on the binary orbit and magnetic field of the neutron star. Astronomy and Astrophysics - A&A, 2022, 664, 10.1051/0004-6361/202243909. insu-03782489

HAL Id: insu-03782489

https://insu.hal.science/insu-03782489

Submitted on 21 Sep 2022

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
# X-ray view of the 2021 outburst of SXP 15.6: Constraints on the binary orbit and magnetic field of the neutron star

G. Vasilopoulos<fnref n="1" />, G. K. Jaisawal<fnref n="2" />, C. Maitra<fnref n="3" />, F. Haberl<fnref n="3" />, P. Maggi<fnref n="1" />, and A. S. Karaferias<fnref n="4" />

Received 28 April 2022 / Accepted 13 June 2022

## ABSTRACT

**Aims.** We conducted a spectral and temporal analysis of X-ray data from the Be X-ray binary pulsar SXP 15.6 located in the Small Magellanic Cloud based on *NuSTAR*, *NICER*, and *Swift* observations during the 2021 outburst.

**Methods.** We present the broadband X-ray spectra of the system based on simultaneous *NuSTAR* and *NICER* observations for the first time. Moreover, we used monitoring data to study the spectral and temporal properties of the system during the outburst.

**Results.** Comparison of the evolution of the 2021 outburst with archival data reveals a consistent pattern of variability, with multiple peaks occurring at time intervals similar to the orbital period of the system (~36 d). Our spectral analysis indicates that most of the energy is released at high energies above 10 keV, while we found no cyclotron absorption line in the spectrum. Analysing of the spectral evolution during the outburst, we find that the spectrum is softer when brighter, which in turn reveals that the system is probably in the super-critical regime in which the accretion column is formed. This places an upper limit on the magnetic field of the system of about $7 \times 10^{11}$ G. The spin-evolution of the neutron star (NS) during the outburst is consistent with an NS with a low magnetic field (~5 $\times 10^{11}$ G), while there is evident orbital modulation that we modelled, and we derived the orbital parameters. We found the orbit to have a moderate eccentricity of ~0.3. Our estimates of the magnetic field are consistent with the lack of an electron cyclotron resonance scattering feature in the broadband X-ray spectrum.

**Key words.** stars: neutron – X-rays: individuals: SXP 15.6 – Magellanic Clouds – X-rays: binaries

## 1. Introduction

Accreting X-ray pulsars (XRs) in high-mass X-ray binaries (HMXBs) are of key importance in the study of accretion and binary evolution. The majority of XRs are found in the so-called Be X-ray binaries (BeXRBs; see *Reig 2011*, for a review), where mass transfer from the donor to the neutron star (NS) occurs through a slow-moving equatorial disk (i.e., decrease disk). A plethora of information about the physical properties of the systems may be acquired during outburst. Type I outbursts occur during the NS periastron passage, while major outbursts with an X-ray luminosity $L_X > 10^{37}$ erg s$^{-1}$ are quite rare and occur on timescales of years to decades (Okazaki et al. 2013; Martin et al. 2014; Martin & Franchini 2021). These rare outbursts enable the study of the broadband X-ray spectral shape and the search for cyclotron resonant scattering features (CRSFs), which offer the only tool for directly measuring the magnetic field at the NS surface. During outbursts, such systems have been found to show both evolution in the continuum spectra and – in some cases – CRSFs (Jaisawal & Naik 2016; Maitra et al. 2018; Staubert et al. 2019), yielding diagnostic data for the state of accretion (Becker & Wolff 2007; Becker et al. 2012; Postnov et al. 2015). Although the Magellanic Clouds (MCs), especially the Small Magellanic Cloud (SMC), are populous in BeXRBs (Haberl & Sturm 2016), the presence of a CRSF is known only for a handful of systems, mainly owing to the lack of coverage in the hard X-ray band during its bright state.

Moreover, X-ray pulsars in the MCs offer quite favourable conditions as they are placed at known distances (in contrast to most Galactic systems) and have only low Galactic absorption. This allows a precise determination of their luminosities. In addition, the knowledge of the spin period and the magnetic field is crucial for constraining the accretion torque models and examine whether most of the BeXRBs are in spin equilibrium. Even if the magnetic field cannot be directly measured via CRSF, the study of the spin evolution can deliver constraints or indirect estimates of the magnetic field of the NS.

SXP 15.6 (also known as XMMU J004855.5–734946) is a BeXRB pulsar in the SMC (Vasilopoulos et al. 2017b). The orbital period of the system has been proposed to be 36.432 d based on optical monitoring data from OGLE (McBride et al. 2017), while a more recent analysis of OGLE data derived an optical period of 36.411 d (Coe et al. 2022). The spin period of the NS was detected in 2016 from *Chandra* observations, and since then, no further strong outburst was witnessed until 2021. On 2021 November 19, a strong outburst was detected by *Swift/XRT* (Coe et al. 2021) at a luminosity level of $2 \times 10^{37}$ erg s$^{-1}$ (0.3–10 keV) for a distance of 62 kpc (Graczyk et al. 2014). Follow-up *NICER* and *NuSTAR* target of opportunity (ToO) observations enabled us to study the spin evolution of the NS as well as its broadband spectral properties. Based on these data, we report estimates of the fundamental properties of the system such as the orbital parameters and the magnetic field of the NS.

<footnote n="1">¹ Université de Strasbourg, CNRS, Observatoire astronomique de Strasbourg, UMR 7550, 67000 Strasbourg, France
e-mail: georgios.vasilopoulos@astro.unistra.fr</footnote>

<footnote n="2">² National Space Institute, Technical University of Denmark, Elektrovej 327-328, 2800 Lyngby, Denmark</footnote>

<footnote n="3">³ Max-Planck-Institut für extraterrestrische Physik, Gießenbachstraße 1, 85748 Garching, Germany</footnote>

<footnote n="4">⁴ Department of Physics, National and Kapodistrian University of Athens, University Campus Zografos, 15783, Athens, Greece</footnote>
## 2. X-ray view of SXP 15.6

### 2.1. 2021 outburst

The 2021 outburst of SXP 15.6 was reported by Coe et al. (2021) as a result of monitoring from the *Swift* SMC Survey (S-CUBED; Kennea et al. 2018). S-CUBED is a high-cadence (one to two weeks) shallow X-ray survey of the SMC that consists of ~140 tiled pointings covering the optical extent of the SMC. The survey has enabled the early detection of some of the brightest outbursts in the MCs (e.g., SMC X-3; Kennea et al. 2016). Following the announcement of the outburst, NICER monitoring observations were performed, while NuSTAR observed the system near its peak flux with a ToO observation. In the following paragraphs we provide information for the NICER and NuSTAR data that were collected during the outburst, as well as Chandra archival data that were used for comparative studies.

### 2.2. Data analysis

#### 2.2.1. NuSTAR

The Nuclear Spectroscopic Telescope Array (NuSTAR) mission carries the first focusing high-energy X-ray telescope in orbit operating in the band from 3 to 79 keV (Harrison et al. 2013). NuSTAR observed the system with a 42 ks DDT observation (obsid: 90701339002) on 2019 November 26 (MJD 59544.40–59545.25). NuSTAR data were analysed with version 1.8.0 of the NuSTAR data analysis software (DAS) and instrumental calibration files from CalDB v20220301. The data were calibrated using the standard settings on the NUPIPELINE script, reducing internal high-energy background, and screening for passages through the South Atlantic Anomaly. We used the NUPRODUCTS script to extract phase-averaged spectra for source and background regions (60'' radius) for each of the two focal plane modules (FPMA/B). Finally, we performed barycentric corrections to the event times of arrival using the satellite's orbital ephemeris files.

#### 2.2.2. NICER

The NICER X-ray Timing Instrument (XTI; Gendreau et al. 2012, 2016) is a non-imaging, soft X-ray telescope on board the International Space Station. The XTI consists of an array of 56 co-aligned concentrator optics (52 currently active) with a field of view of ~30 arcmin² in the sky. Each unit is associated with a silicon drift detector (Prigozhin et al. 2012), operating in the 0.2–12 keV band, yielding a ~100 ns time resolution and spectral resolution of ~85 eV at 1 keV.

For the current study, we analysed NICER data obtained between MJD 59535 and MJD 59605. Data were reduced using HEASOFT version 6.29, NICER DAS version 2020-04-23_V007a, and the calibration database (CALDB) version v20210707. For the analysis, we selected good-time intervals with the nimaketime script using standard options. After inspecting the resulting light curves in different bands, we identified increased flaring activity due to background contamination. To mitigate the effects of the background, we altered some of the standard filtering parameters. We opted for ISS not in the South Atlantic Anomaly region, a source elevation >15° above the Earth limb (>30° above the bright Earth), and a magnetic cut-off rigidity (COR_SAX) > 2.0 GeV c⁻<fnref n="1" />. Unfortunately, a large fraction of the obtained data – especially during the lower flux states – was not useful due to enhanced background activity during the monitoring period and was filtered out. Finally, for the timing analy-

sis, we performed barycentric corrections to the event times of arrival using the barycorr tool and the JPL DE405 planetary ephemeris.

Because NICER is not an imaging instrument, the X-ray background is calculated indirectly. For systems in the direction of the MCs, the 3C50 tool (Remillard et al. 2022) method is optimal (see Treiber et al. 2021). This approach uses a number of background proxies in the NICER data to define the basis states of the background database.

#### 2.2.3. Swift

X-ray monitoring observations of SXP 15.6 have been obtained by the Neil Gehrels Swift Observatory (Swift; Gehrels et al. 2004) X-ray Telescope (XRT; Burrows et al. 2005). All archival XRT data were retrieved through the UK Swift science data centre<fnref n="1" /> and were analysed using standard procedures (Evans et al. 2007, 2009).

#### 2.2.4. Chandra

SXP 15.6 was observed by Chandra (obsid:18885) in July 2016 (MJD 57575.34) with an exposure time of 25 ks. Analysis of these data has revealed the presence of a single-peaked pulse profile at $L_X = 4 \times 10^{36}$ erg s⁻<fnref n="1" /> in the 0.3–10.0 keV band (Vasilopoulos et al. 2017b). In the current study, we used the same data to create pulse profiles for comparison with the pulse profiles during the 2021 outburst. Data reduction was performed with the CIAO v4.13 software (Fruscione et al. 2006) using standard options through chandra_repro script. Source events were extracted from a 5'' region.

## 3. Results

### 3.1. Broadband spectral properties

Spectral analysis was performed using XSPEC v12.10.1f (Arnaud 1996). Photoelectric absorption by the interstellar gas was modelled by the tbabs component in XSPEC, with solar abundances set according to Wilms et al. (2000) and atomic cross sections from Verner et al. (1996). To fit the NICER spectra, we used PG-statistics, which implements Cash-statistics (Cash 1979), with a non-Poisson background model. For the modelling of NuSTAR spectra (3.0–78 keV), we used Cash-statistics. We note that NuSTAR spectra are dominated by background above 40 keV and thus lack the sensitivity to detect faint features at high energies.

It is typical for BeXRBs in the MCs to model the X-ray absorption with a combination of two components to account for Galactic absorption and intrinsic absorption within the SMC and around the binary (e.g., Vasilopoulos et al. 2013, 2016, 2017a). The first component was fixed to the Galactic foreground value of $4 \times 10^{20}$ cm⁻<fnref n="1" /> (Dickey & Lockman 1990). The second component was left free to account for the absorption near the source or within the SMC, thus elemental abundances were fixed at 0.2 solar (Russell & Dopita 1992). In the spectral modelling, we determine whether there is evidence of absorption in addition to the Galactic absorption.

In BeXRBs, hard X-ray spectra originate from the accretion column and can be fitted by a phenomenological power-law-like shape with an exponential high-energy cut-off above (or of about) 10 keV (e.g., Müller et al. 2013; Sturm et al. 2014;

<footnote n="1">¹ http://www.swift.ac.uk/user_objects/</footnote>
<figure><img src="image_2.png" /><figcaption>Fig. 1. Broadband spectrum of SXP 15.6. NuSTAR spectra are shown in black, and two NICER snapshots are shown in red and cyan. The upper panel presents the best-fit model composed of an absorbed cut-off power law and a soft thermal component. The ratio plots show the comparison of the tested models and the data: (b) Absorbed cut-off power law, (c) absorbed cut-off power law plus disk blackbody, and (d) cut-off power law with partial coverage. The model parameters are listed in Table 1.</figcaption></figure>

Jaisawal et al. 2018; Vasilopoulos et al. 2020). In many cases, BeXRB spectra show residuals at soft energies. These residuals are often referred to as a ‘soft-excess’ with a physical origin that is attributed to one mechanism or a combination of mechanisms such as emission from the accretion disk, emission from the NS surface, hot plasma around the magnetosphere, or partial absorption from material around the NS (Hickox et al. 2004). In this study we investigate the broadband spectrum for signatures of these mechanisms.

NuSTAR data were obtained quasi-simultaneously (i.e., less than a day apart) with two NICER visits (obsid: 4202430107 and 4202430108). An absorbed cut-off power-law model sufficiently describes the NuSTAR spectrum. However, fitting the combined NICER and NuSTAR data, we see significant structure in the residuals caused by the soft excess (see Fig. 1). To eliminate the residual structure, we need to add either a partially covering absorber or a soft thermal component to the model (see the lower panels in Fig. 1). For the thermal component, we used a disk blackbody (diskbb in XSPEC). We also included a model with a combined partial coverage on a continuum composed of a cut-off power law and a thermal component. In all our tests, the column density of the SMC intrinsic absorption was unconstrained and tends to zero. For the partial covering model, the Tuebingen-Boulder ISM model (tbpcf in xspec) only provided spectral multiplicative components for standard solar abundances because a first-order approximation for the SMC abundance $N_H$ values obtained by this model should be increased by a factor of ~5 compared to Galactic ones. The best-fit parameters of all tested models are shown in Table 1. Based on the derived parameters, we also computed the mass-

accretion rate corresponding to the bolometric luminosity by assuming that all gravitational energy is converted into radiation<fnref n="2" /> (e.g., Campana et al. 2018; Frank et al. 2002). We note that for all the tested models, the best fit yields an acceptable fit statistics with a reduced $\chi^2$ lower than one. To estimate the uncertainties, we implemented a Markov chain Monte Carlo (MCMC) model through XSPEC. We used the Goodman–Weare algorithm with 20 walkers and a total length of 50 000. For the initial burn-in phase, we needed 30 000 steps before the chain reached equilibrium. We then generated parameter errors (90% confidence) based on the chain values. To further test the goodness of the fit, we also simulated spectra based on the MCMC chain parameters. We found that for all except for the simplest models (i.e., absorbed cut-off power law), only a small number of the simulated spectra had a better fit statistics than the real spectra. Thus we should be at the limit of our capabilities in testing more complicated spectral models. We finally note that we found no evidence of an Fe K$_\alpha$ line in the spectrum nor any broad absorption feature that is consistent with a CRSF.

## 3.2. Long-term light curves and outburst evolution

Monitoring data in the soft X-ray band with NICER enables us to study the evolution of the 2021 outburst and to compare it with Chandra data from 2016. NICER monitoring data are of sufficient statistical quality to allow us to extract enough counts and perform spectral modelling. The typical Swift/XRT exposure within one day is about 1000–2000 s, but the effective area of the detector is significantly smaller than that of NICER. Thus we used XRT data to only estimate average count rates for each XRT data set.

The 20 individual NICER spectra were fitted with an empirical absorbed power-law model in the 0.5–8 keV band. We also attempted to fit the spectra with a power law with a cut-off, but the cut-off always converged to very high values (i.e., above 100 keV). This was expected because the cut-off seen in the broad spectra is well above the upper bound of NICER spectra. For the spectral modelling, we used two absorption components as described above. We found the absorption to be consistent with the fixed foreground value of $4 \times 10^{20} \text{ cm}^{-1}$, while the power-law photon index had a mean value of $\Gamma = 0.97$, with evidence of the spectrum becoming softer when brighter at the brightest phase of the outburst. In Fig. 2 we plot the evolution of the spectral parameters and the flux for a 40-day interval (only one observation exists after MJD 59580).

After analysing the NICER and Swift/XRT monitoring data as well as the broadband spectra, we computed the evolution of the bolometric luminosity during the 2021 outburst and compared it with archival data. We converted the XRT count rates into 0.3–10.0 keV fluxes using the average spectral parameters inferred from the NICER spectral fits. A conversion factor of $3.16 \times 10^{37} \text{ erg cm}^{-2} \text{ s}^{-1}/(\text{c/s})$ was used for all XRT data, while errors were estimated based on the count rate uncertainties. The broadband unabsorbed $L_X$ was estimated from the broadband spectra. Most of the energy is emitted above 10 keV because the ratio of the broadband (0.3–80 keV) to narrow band (0.3–10.0 keV) luminosity was ~3.5. The 2021 X-ray light curve is shown in Fig. 3. In the same figure we overplot the 2016 XRT monitoring data, time-shifted so that the main peaks of each outburst match.

Monitoring data can also be used to investigate whether the system has entered the super-critical regime, in which

<footnote n="2">² I.e. $L_X = G\dot{M}/R \approx 0.2 \dot{M}c^2(M/1.4M_\odot)(R/10 \text{ km})^{-1}$.</footnote>
<table><caption>**Table 1.** Best-fit parameters of the spectral empirical models.</caption><thead><tr><th>Model</th><th>Parameter</th><th>cPL</th><th>cPL+dBB</th><th>PC*cPL</th><th>PC*(cPL+dBB)</th><th>PC*(cPL+dBB)</th><th>Units</th></tr><tr><th></th><th></th><th>Value</th><th></th><th></th><th></th><th></th><th></th></tr></thead><tbody><tr><td></td><td>N<sub>H</sub> Gal <sup>(a)</sup></td><td>4 (fixed)</td><td>4 (fixed)</td><td>4 (fixed)</td><td>4 (fixed)</td><td>4 (fixed)</td><td>10<sup>20</sup> cm<sup>-2</sup></td></tr><tr><td></td><td>N<sub>H</sub> SMC <sup>(a)</sup></td><td>→ 0</td><td>→ 0</td><td>→ 0</td><td>→ 0</td><td>&lt;3</td><td>10<sup>20</sup> cm<sup>-2</sup></td></tr><tr><td>TBpcf</td><td>N<sub>H</sub></td><td>—</td><td>—</td><td>53<sup>+7</sup><sub>-6</sub></td><td>47<sup>+12</sup><sub>-10</sub></td><td>53<sup>+7</sup><sub>-5</sub></td><td>10<sup>22</sup> cm<sup>-2</sup></td></tr><tr><td></td><td>cov. frac.</td><td>—</td><td>—</td><td>45 ± 3</td><td>32<sup>+8</sup><sub>-7</sub></td><td>45 ± 3</td><td>%</td></tr><tr><td>dBB <sup>(b)</sup></td><td>kT<sub>BB</sub></td><td>—</td><td>0.44<sup>+0.014</sup><sub>-0.06</sub></td><td>—</td><td>0.39<sup>+0.07</sup><sub>-0.08</sub></td><td>0.07 (fix)</td><td>keV</td></tr><tr><td></td><td>Norm<sub>BB</sub></td><td>—</td><td>6.3–6.0<sup>+5.3</sup><sub>-5.3</sub></td><td>—</td><td>6.7<sup>+5</sup><sub>-2.5</sub></td><td>2100<sup>+23000</sup><sub>-1400</sub></td><td>sin θ (R<sub>BB</sub>/D<sub>10</sub>)<sup>2</sup></td></tr><tr><td>cPL <sup>(c)</sup></td><td>Γ</td><td>0.85±0.02</td><td>0.49<sup>+0.05</sup><sub>-0.03</sub></td><td>1.05 ± 0.03</td><td>0.82<sup>+0.11</sup><sub>-0.1</sub></td><td>1.04<sup>+0.04</sup><sub>-0.03</sub></td><td>—</td></tr><tr><td></td><td>E<sub>c</sub></td><td>20.6 ± 0.9</td><td>12.8<sup>+0.7</sup><sub>-0.7</sub></td><td>18.6<sup>+1.1</sup><sub>-0.9</sub></td><td>15.6<sup>+1.7</sup><sub>-1.2</sub></td><td>18.6<sup>+1.1</sup><sub>-0.9</sub></td><td>keV</td></tr><tr><td></td><td>Norm <sup>(c)</sup></td><td>101 ± 2</td><td>89.8<sup>+1.7</sup><sub>-1.6</sub></td><td>109.4<sup>+1.8</sup><sub>-1.8</sub></td><td>100<sup>+5</sup><sub>-4</sub></td><td>109.4<sup>+2.1</sup><sub>-1.8</sub></td><td>10<sup>-12</sup> erg cm<sup>-2</sup> s<sup>-1</sup></td></tr><tr><td colspan="8">Other information</td></tr><tr><td></td><td>Total fit stat./d.o.f.</td><td>1399.56/1035</td><td>954.11/1033</td><td>940.51/1033</td><td>923.36/1031</td><td>940.47/1032</td><td></td></tr><tr><td></td><td>Goodness <sup>(d)</sup></td><td>100%</td><td>4%</td><td>1%</td><td>0%</td><td>4%</td><td></td></tr><tr><td></td><td>L<sub>X</sub> <sup>(e)</sup></td><td>4.67 ± 0.08</td><td>4.13 ± 0.08</td><td>5.03 ± 0.08</td><td>4.6 ± 0.2</td><td>5.03 ± 0.10</td><td>10<sup>37</sup> erg s<sup>-1</sup></td></tr><tr><td></td><td>M̄ <sup>(e)</sup></td><td>2.60 ± 0.03</td><td>2.30 ± 0.04</td><td>2.80 ± 0.04</td><td>2.55 ± 0.13</td><td>2.8 ± 0.05</td><td>10<sup>17</sup> g s<sup>-1</sup></td></tr></tbody></table>

**Notes.** <sup>(a)</sup>Galactic absorption was fixed to this value (see text for details). <sup>(b)</sup>The disk blackbody (diskbb in xspec) radius may be estimated from the normalisation of the model and a distance of 62 kpc (i.e., D<sub>10</sub> = 6.2) assuming a disk inclination angle θ. <sup>(c)</sup>Cut-off power law (cut-offpl in xspec), where the normalisation is the unabsorbed X-ray flux in the 0.3–80 keV band. <sup>(d)</sup>Based on XSPEC simulations, denotes the percentage of simulated spectra that when fitted with the same model yield lower test statistics than the data. <sup>(e)</sup>Unabsorbed X-ray luminosity (0.3–80 keV) for a distance of 62 kpc for the cut-off power-law component. Mass-accretion rate onto the NS, assuming L<sub>X</sub> = 0.2M<sub>c</sub>s<sup>-2</sup>.

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

**Fig. 2.** Spectral results of NICER monitoring data obtained during the 2021 outburst of SXP 15.6. The bottom panel shows the unabsorbed flux computed in the 0.3–10.0 keV band. The red shaded region indicates the epoch of the NuSTAR ToO. No strong evidence for a change in the column density of the absorbing material (top panel) is seen, but we found evidence of spectral evolution at higher fluxes with the spectral shape (characterised by Γ, middle panel) becoming softer when brighter above 2×10<sup>-11</sup> erg cm<sup>-2</sup> s<sup>-1</sup> (0.3–10.0 keV).

the accretion column has been formed above the NS surface (Becker et al. 2012). The simplest proxy for this transition

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

**Fig. 3.** X-ray light curves of the 2016 and 2021 outbursts. Luminosities are absorption corrected (0.3–80.0 keV) and computed for the SMC distance. Archival data from 2016 (first point at MJD 57547) are shifted in time to match the 2021 luminosity peaks. Vertical green lines mark the 36.411 d optical period phased on the first X-ray maximum. For comparison we also show the shifted epoch of the Chandra ToO that falls around the same orbital phase as the NuSTAR ToO.

is the change in hardness of the spectrum with intensity (Reig & Nespoli 2013). Following the nomenclature of Reig & Nespoli (2013), for low luminosities, the spectra of many BeXRB pulsars appear to be harder when brighter (i.e., so-called horizontal branch), while above a critical limit, the systems enter the diagonal branch, where they appear to be softer when brighter. For the intensity, we used the bolometrically corrected L<sub>X</sub>, while for the colour proxy, we used the power-law photon index from NICER and the Swift/XRT hardness ratios. We define the hardness ratios (HR) as HR = (H − S)/(H + S), where [H, S] is the count rate in a specific hard and soft energy band. In Fig. 4 we plot the intensity–colour diagram of SXP 15.6 from 2021 monitoring data. There is evidence that the system has
<figure><img src="image_5.png" /><figcaption>**Fig. 4.** Hardness intensity diagram of the 2021 outburst as monitored by NICER and Swift/XRT. The upper panel shows the results of the NICER spectral analysis, and the lower panels show the results based on HR. There is evidence that the source has entered the DB, thus the system appears to be close to or above the critical regime.</figcaption></figure>

entered the diagonal branch and appears to be close to or above the critical limit for accretion column formation.

### 3.3. Temporal properties – pulse profiles

To search for a periodic signal, we used the epoch-folding Z-search method (Buccheri et al. 1983) implemented through HENdrics command-line scripts and Stingray (Huppenkothen et al. 2019). For NuSTAR data, we searched for a periodic signal in the 3–40.0 keV range (35-960 PI channel). Our final estimate of the spin period and its uncertainties was based on the time of arrival (ToA) method (e.g., Tsygankov et al. 2020). We first used HENdrics to derive a most probable period, then we estimated ToAs of individual pulses for 16 intervals, and we finally used PINT<fnref n="3" /> (Luo et al. 2021). From the above we derived a period of 15.6395 ± 0.0004 s for the 2021 NuSTAR data. The reported period for the 2016 data was 15.6398 ± 0.0009 s and is consistent within the uncertainties with the newly derived period. For consistency, we implemented the ToA procedure to estimate a period for the 2016 Chandra data. Results of the ToA method are presented in Table 2. This yielded a period of 15.6396 ± 0.0014 s. All tests indicated that the period of the NS has remained unchanged within uncertainties for more than five years.

The strength of the periodic modulation is typically quantified through the root-mean-squared (rms) pulsed fraction (PF). This is given by

$$PF_{rms} = \frac{\left(\sum_{j=1}^{N} (R_j - \bar{R})^2 / N\right)^{1/2}}{\bar{R}}, \quad (1)$$

<table><caption>**Table 2.** Spin periods SXP 15.6.</caption><thead><tr><th>Obsid</th><th>MJD<sup>(†)</sup> (d)</th><th>Δt (d)</th><th>P (s)</th><th>#ToAs<sup>(‡)</sup></th></tr></thead><tbody><tr><td colspan="5"><em>NICER</em></td></tr><tr><td>4202430101</td><td>59538.50</td><td>0.8</td><td>15.6398 ± 0.0004</td><td>9</td></tr><tr><td>4202430102</td><td>59539.12</td><td>0.3</td><td>15.6394 ± 0.0019</td><td>5</td></tr><tr><td>4202430103</td><td>59540.51</td><td>0.7</td><td>15.6395 ± 0.0005</td><td>6</td></tr><tr><td>4202430104</td><td>59541.44</td><td>0.9</td><td>15.6394 ± 0.0005</td><td>9</td></tr><tr><td>4202430106</td><td>59543.48</td><td>0.8</td><td>15.6393 ± 0.0009</td><td>5</td></tr><tr><td>4202430107</td><td>59544.34</td><td>0.7</td><td>15.6394 ± 0.0007</td><td>5</td></tr><tr><td>4202430109</td><td>59550.70</td><td>0.5</td><td>15.6414 ± 0.0009</td><td>8</td></tr><tr><td>4202430110</td><td>59551.61</td><td>0.5</td><td>15.6424 ± 0.0013</td><td>4</td></tr><tr><td>4202430111</td><td>59552.45</td><td>0.7</td><td>15.6430 ± 0.0007</td><td>5</td></tr><tr><td>4202430112</td><td>59553.04</td><td>0.1</td><td>15.6438 ± 0.0010</td><td>EF</td></tr><tr><td>4202430114</td><td>59556.72</td><td>0.03</td><td>15.6444 ± 0.003</td><td>EF</td></tr><tr><td>4202430117</td><td>59560.77</td><td>0.1</td><td>15.650 ± 0.007</td><td>3</td></tr><tr><td>4202430118</td><td>59561.58</td><td>0.1</td><td>15.65 ± 0.004</td><td>3</td></tr><tr><td>4202430119</td><td>59563.27</td><td>0.1</td><td>15.645 ± 0.007</td><td>3</td></tr><tr><td>4202430121</td><td>59571.45</td><td>0.9</td><td>15.6378 ± 0.0006</td><td>9</td></tr><tr><td colspan="5"><em>NuSTAR</em></td></tr><tr><td>90701339002</td><td>59544.8</td><td>0.8</td><td>15.6395 ± 0.0004</td><td>16</td></tr></tbody></table>

**Notes.** (†) Middle epoch of observation. (‡) Number of ToAs used to refine the spin period. For observations for which this was not possible, we used periods based on epoch folding.

where *N* is the number of phase bins, $R_j$ is the background-subtracted count rate in the *j*th phase bin, and $\bar{R}$ is the average count rate in all bins (e.g., Wilson-Hodge et al. 2018). We used this definition to estimate the PF<sub>rms</sub> in different NuSTAR energy bands. We found no significant change in the PF within different energy bands up to 20 keV. However, above 20 keV, the PF increases and the pulsations become almost twice as strong at the highest energies. We also note that at above 40 keV, the contribution of background photons is ~50% of the net counts. The increasing PF with energy is typical for accreting pulsars and is attributed to hard photons that are emitted from the sides of the accretion column. They are more beamed than soft photons (e.g., Lutovinov & Tsygankov 2009).

In Fig. 5 we present the folded pulse profiles for different energy bands covering the full NuSTAR energy range. We opted to also show the soft energy band (1.6–5.0 keV) in order to compare with NICER and Chandra pulse profiles. In Fig. A.1 we present the folded pulse profiles from all NICER observations for which a period could be estimated.

### 3.4. Spin evolution

To investigate the spin evolution during the outburst, we implemented the same method on the NICER data between 0.5 and 8.0 keV. We first computed the most probable period by epoch folding, and then we refined the period and its uncertainty based on the ToA of individual pulses. Given the shape of the pulse profile, the template used for ToAs can drastically vary from one observation to the next. Thus, we used one universal template for all observations to characterise the off-phase of the pulse. The method was successful when two conditions occurred. Firstly, we need more than three NICER snapshots to be performed within one day, and secondly, the total number of counts must be high enough to obtain meaningful pulse profiles. For all other snapshot phases, it was challenging to connect the ToAs and even impossible at times due to multiple peaks

<footnote n="3">³ https://github.com/nanograv/pint/</footnote>
<figure><img src="image_6.png" /><figcaption>Fig. 5. Pulse profiles of SXP 15.6 from the NuSTAR observation. Different energy ranges are used to track the changes with energy. Shaded horizontal regions denote the limit of the statistical significance variability above the constant level hypothesis.</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>Fig. 6. Period evolution of SXP 15.6 within the 2021 outburst. The predicted spin-up rate based on the observed luminosity and three values of the magnetic field are plotted with dashed lines.</figcaption></figure>

in the periodogram with similar intensity. This problem is often encountered in slow pulsars that are observed with gaps (e.g., Zolotukhin et al. 2017; Vasilopoulos et al. 2018a).

The period evolution is shown in Fig. 6. The overall trend seems linear, although a large gap occurs in the data. Although due to the sampling variability due to orbital Doppler shifts is visible, the span of the NICER points with good timing solutions

is ~33 d and is comparable to the ~36.4 d optical period. Thus the secular evolution between the first and last point alone should not be affected much by orbital effects. With this assumption, we found an average spin-up of $\dot{P} = (7.0 \pm 2.5) \times 10^{-10} \text{ s s}^{-1}$ (or $\dot{\nu} = (2.8 \pm 1.0) \times 10^{-12} \text{ Hz s}^{-1}$). This value should be the approximate intrinsic spin-up due to accretion. Alternatively, the intrinsic spin up of the NS can be calculated (see Vasilopoulos et al. 2019, 2020, for method) due to the mass accretion rate as derived by the observed bolometric $L_X$ (see Fig. 3).

Here are the basic steps for our calculation. We assumed mass transfer from a Keplerian disk, thus the induced torque due to accretion alone is $N_{acc} \approx \dot{M} \sqrt{GM_{NS}R_M}$. The total torque can be expressed in the form of $N_{tot} = n(\omega_{fast})N_{acc}$, where $n(\omega_{fast})$ is a dimensionless function that accounts for the coupling of the magnetic field lines to the accretion disk (for details, see Wang 1995; Parfrey et al. 2016). The spin-up rate of the NS is then given by

$$ \dot{\nu} = \frac{n(\omega_{fast})}{2\pi I_{NS}} \dot{M} \sqrt{GM_{NS}R_M}, \qquad (2) $$

where $I_{NS} \approx (1 - 1.7) \times 10^{45} \text{ g cm}^2$ is the moment of inertia of the NS (e.g., Steiner et al. 2015).

To model the intrinsic spin evolution due to accretion, we just need to numerically solve Eq. (2) in time assuming a constant magnetic field strength. We used 1000 time steps in this span over the NICER monitoring. For each time step, the mass accretion was estimated by interpolating the observed flux in the 0.3–10 keV (see Fig. 2). Then we converted the flux into bolometric luminosity assuming the spectral parameters of the broadband spectra. Bolometric luminosity was then used as a proxy for mass accretion (i.e., $L_X = GM\dot{M}/R_M$). For $n(\omega_{fast})$, we followed the Wang (1995) model (see their Eq. (19)). For all calculations we adopted standard NS parameters (i.e. 12 km radius, $1.4 M_{\odot}$ mass and $I_{NS} = 1.3 \times 10^{45} \text{ g cm}^2$). We repeated this process for various magnetic field strengths, while the results are shown in Fig. 6. It is evident that the observed secular spin up is consistent with a magnetic field strength close to $3 \times 10^{11} \text{ G}$, while very low (<$10^{11} \text{ G}$) or high (>$10^{12} \text{ G}$) B values seem to be inconsistent with observations. Because the first and last pointing are separated by about one orbital period, we can neglect any orbital effects in this first-order approximation. However, we investigate any effects in the next section.

### 3.5. Orbital evolution

The observed spin evolution shows significant variation from the expected evolution caused by the intrinsic spin-up that is due to accretion. The remaining residuals may be due to orbital modulation. Keplerian orbits are described with five orbital elements: the orbital period ($P_{orb}$), the orbital eccentricity (e), the argument of periastron of the stellar orbit (ω), the velocity semi-projected axis ($a \sin i$ in light seconds), and for the orbital phase, we finally used the time of a mean longitude of 90 deg (i.e., $T_{\pi/2}$). The orbital modulation is often modelled with the intrinsic spin-up due to accretion, as described in the previous section (see Sect. 3.4), for Fermi/GBM pulsars (Sugizaki et al. 2017). However, most of the GBM monitored pulsars are Galactic sources that are monitored for extensive periods. To model our data set, the parameter space must be properly modelled in order to identify degeneracies between the model parameters. Thus, to fit the model to the data, we implemented a nested sampling algorithm for a Bayesian parameter estimation and estimated posterior distributions for the parameters of the standard accretion
<figure><img src="image_8.png" /><figcaption>**Fig. 7.** Modelling of the spin period evolution of SXP 15.6 based on orbital modulation and intrinsic spin-up. The most probable solution is shown as the dashed line, and the dotted black line shows the Doppler shifts due to orbital modulation alone, without the spin-up due to accretion. A family of orbits drawn from the sample of the posterior distribution of parameters from the Bayesian modelling is shown as red lines. The shaded region is the 99% quantile of the models. Although Swift/XRT data present hits of periodic modulation, the uncertainties are quite high and cannot constrain our model. For visualisation purposes, we plot one of these XRT measurements.</figcaption></figure>

torque models and binary orbital parameters. In terms of statistical treatment, similar methods have been used to model radial velocity curves from binary systems (Fulton et al. 2018).

To derive the posterior probability distributions and the Bayesian evidence, we used the nested sampling MC algorithm MLFriends (Skilling et al. 2004; Buchner 2019) that employs the ULTRANEST<fnref n="4" /> package (Buchner 2021). An outline of this method with applications to accreting pulsars will be presented by Karaferias et al. (2022).

Because of the gaps in the NICER monitoring and the high background, good-timing data exist only for the first 20 days of the monitoring. In an effort to improve our data set, we also searched for pulsations in Swift/XRT data. However, because of the low statistics and small number of XRT snapshots within a day, typical period uncertainties were about 0.01–0.05 s (see Fig. 7). They therefore offer little information for our study. To constrain our model parameters, we fixed the orbital period to 36.411 days and limited the magnetic field strength (i.e., log (B[G]) ∈ [10, 13]). With these assumptions and for typical parameters of the NS (12 km and 1.4 $M_⊙$), we estimated the posterior distribution for our other model parameters. These parameters are listed in Table 3, while in Fig. A.2 we show the corner plot of posterior distributions for our solution. In Fig. 7 we show a sample of 100 random orbital models from the posterior distribution together with the most probable model.

## 4. Discussion

The latest outburst of SXP 15.6 started in late November 2021, while the system remained in a bright state until early 2022. The bolometric luminosity during the event reached $7 \times 10^{37}$ erg s$^{-1}$, which presents the brightest stage ever observed for the system. This luminosity translates into ~30% of the Eddington limit for a typical NS (i.e., 12 km and 1.4 $M_⊙$).

The combined Swift/XRT and NICER light curve of the 2021 outburst reveals a complex structure with two peaks separated by the orbital period. The increased flux at later epochs indicates

<table><caption>**Table 3.** Orbital parameters of the binary SXP 15.6.</caption><thead><tr><th>Parameter</th><th>Prior</th><th>Result</th><th>Units</th></tr></thead><tbody><tr><td>P<sub>orb</sub></td><td>36.411 (fixed)</td><td></td><td>d</td></tr><tr><td>e</td><td>0–0.99</td><td>0.29<sup>+0.15</sup><sub>-0.16</sub></td><td>–</td></tr><tr><td>ω</td><td>0–360</td><td>20<sup>+26</sup><sub>-27</sub></td><td>°</td></tr><tr><td>a sin i</td><td>50–400</td><td>162<sup>+50</sup><sub>-40</sub></td><td>1-s</td></tr><tr><td>T<sub>π/2</sub></td><td>59 520–59 560</td><td>59537.9<sup>+1.2</sup><sub>-1.2</sub></td><td>MJD (d)</td></tr></tbody></table>

a third peak, but not enough monitoring data were collected to further explore this behaviour. The general structure of the event matches the behaviour seen in 2016 very well (see Fig. 3), with minor differences in the relative intensity of the three peaks. The mismatch of the X-ray and optical period could be attributed to a precessing decretion disk around the Be star (for an observational example and a theoretical application, see Treiber et al. 2021; Martin & Franchini 2021). This precession can cause an evolving period between outbursts if the disk is moving retrograde to the NS orbit. Nevertheless, the self-similarity of the 2016 and 2021 outbursts with three peaks is quite intriguing and reveals that the Be disk and the NS geometrical configuration behave in a repetitive manner.

The NS spin evolution shows evidence of orbital and intrinsic spin-up due to accretion. Because the first and last day of the NICER observations (MJD 59535-75 interval) consisted of a large number of snapshots, we were able to constrain the intrinsic spin-up and found a value of $\dot{v} = (2.8 \pm 1.0) \times 10^{12}$ Hz s$^{-1}$. To improve this estimate, we simultaneously modelled the orbital modulation and intrinsic spin-up with a Bayesian approach. We found that the spin-up is consistent with an NS pulsar with a magnetic field strength of $5 \times 10^{11}$ G with an uncertainty of a factor of ~3. Moreover, the orbit has a moderate eccentricity with a value of ~0.29, while values up to 0.6 cannot be statistically excluded. Interestingly, the corner plot in Fig. A.2 shows that lower magnetic fields are apparently favoured by more eccentric orbits, and the uncertainties in some parameters are still high. From the orbital parameters (i.e., $\omega$ and $T_{\pi/2}$) we can also estimate the epoch of periastron $T_{per}$, which is found to be MJD 59531 ± 4. The value of $T_{per}$ seems to match or lead by a few days the X-ray maxima observed in the X-ray light curve (see Fig. 3). Figure 7 interestingly shows that quite a few orbits with higher eccentricity are still statistically acceptable and cannot be excluded without better data. Future independent measurements that could constrain the orbital parameters would help us revisit the system and tighten the constraints on the magnetic field strength.

The value of the magnetic field that we estimate can be compared with the estimates of Coe et al. (2022), who found a value of $3.7 \times 10^{12}$ G, assuming the source is in spin equilibrium. For an NS pulsar to rotate near spin-equilibrium, constant mass accretion is assumed and the disk inner radius is required to rotate with about the same angular velocity as the NS. However, there are a few caveats in these assumptions. Because BeXRBs are extremely variable systems, the assumption of steady accretion does not hold, and for systems rotating near equilibrium, torque transfer is a non-linear problem. Moreover, rotation near equilibrium does not mean that the spin-up rate is zero. Under the assumption of steady accretion, even some of the systems with the highest spin-up rate can be argued to be near equilibrium (Pan et al. 2022). In addition, finding the secular spin evolution to be very small does not necessarily mean that the NS is in equilibrium. The spin evolution

<footnote n="4">⁴ https://johannesbuchner.github.io/UltraNest/</footnote>
is a dynamical problem in which the accretion duty cycle plays an equally important role in determining any secular spin change between two epochs. Estimates of the magnetic field assuming spin-equilibrium are therefore subject to large systematic uncertainties and should be considered much less accurate than our estimates from torque modelling.

The broadband spectra from the combined NuSTAR and NICER observations lack any significant line features that could be associated with emission from hot plasma (e.g., an Fe Kα line) or a CRSF. The high-energy part of the spectrum is well explained by a power law with a cut-of,f in agreement with other BeXRBs. However, the soft part of the spectrum requires additional components in order to be explained. We tested whether the addition of a thermal component or a partial absorber could improve the quality of the fit. Either or both of these components would result in an acceptable fit.

The thermal component is particularly interesting because it can be interpreted as the inner region of an accretion disk. An analytical form for the magnetospheric radius can be expressed in the following form (Frank et al. 2002; Campana et al. 2018):

$$R_M = 0.5 \left( \frac{\mu^4}{2GM\dot{M}^2} \right)^{1/7}, \qquad (3)$$

where M is the NS mass, $\mu = BR^3/2$ is the magnetic dipole moment, with R the NS radius and B the NS magnetic field strength at the magnetic poles. For typical parameters of the NS (12 km and 1.4 M⊙) and for the observed mass-accretion rate and assuming a polar magnetic field strength of $5 \times 10^{11}$ G, the disk radius should be about 600 km. However, our best-fit spectral model yields a much smaller inner disk radius (i.e., ~15–20 km, see Table 1), which is comparable to the NS radius. This radius would translate into a B value of $10^9$ G, which is unrealistically small for a BeXRB pulsar. This value is also at odds with the results of the torque modelling because Fig. 7 shows that a magnetic field lower than $10^{11}$ G would not be able to explain the spin evolution over the 30-day period, unless a quite eccentric orbit is assumed. Because the disk could be seen edge on, another way to estimate its size is from its temperature. By adopting the mass-accretion rate from the pulsed continuum and for a temperature of 0.4 keV (from the fit), the inner radius of a standard disk would be about 60–70 km, depending on the spectral hardening parameter (see Kubota et al. 1998; Zimmerman et al. 2005). This is still significantly smaller than the size estimated from torque modelling, which predicts a disk size of about 600 km for a magnetic field of $5 \times 10^{11}$ G and $\dot{M} \sim 2.5 \times 10^{17}$ g s$^{-1}$. For these parameters, the disk temperature should be ~70 eV. As seen in the last column of Table 1, a model with a disk with a temperature fixed at 70 eV can sufficiently explain the observed spectra. The disk size for this temperature is about 200–900 km.

With respect to the partial coverage model, the rationale behind this model is that the line of sight between the observer and the source is not a line in a mathematical sense. Because of the extent of the emitting region, the geometrical problem is better described as a superposition of multiple line of sights. The absorber present in the vicinity of the source can be imagined as a collection of dense clouds that lie inside the magnetosphere or engulf the binary itself. Partial covering manifests itself as broad humps or bumps below 10 keV depending on the column density of the partial absorber. For example for the Galactic pulsars GX 304-1 and Her X-1, the partial absorber has typical values of $N_H \sim 10^{-70} \times 10^{22}$ cm$^{-1}$ and a covering fraction of 30–40% depending on the selection of the continuum model (e.g., Endo et al. 2000; Asami et al. 2014; Jaisawal et al. 2016). These

values are quite similar to the values we find in SXP 15.6 and other BeXRB systems in the MCs (e.g., Vasilopoulos et al. 2018b).

This discussion demonstrates the complexity of the soft excess, which in our case cannot be modelled without overfitting the data. Thus we caution about the spectral parameters for the soft excess. For these empirical models, a physical interpretation is difficult because an acceptable goodness of fit and almost flat spectral residuals can be obtained for a wide range of parameters.

We can use the results of the broadband spectroscopy to convert observed fluxes in the 0.3–10.0 keV band into bolometric luminosities. The bolometric correction is ~3.5–4, where its uncertainty is due to the model selection because partial covering models generally give higher unabsorbed flux. This is slightly higher than what is observed in accreting pulsars, since the study of nine frequently observed Galactic systems has shown that 37–62% of the energy is released in the 0.5–10.0 keV band, when the total Luminosity is below the Eddington limit (Anastasopoulou et al. 2022). The value we find is consistent (within $2\sigma$) with the expected value in the observed luminosity range as determined by Anastasopoulou et al. (2022). The bolometric correction is important if we would like to track transitions in the spectral properties associated with changes in the accretion column (Becker et al. 2012). Following Reig & Nespoli (2013), we can argue that the spectral behaviour during the 2021 monitoring is consistent with the super-critical regime. This means that we are above the critical luminosity where the accretion column has started to form (see also Postnov et al. 2015). Following Becker et al. (2012) this luminosity is given by

$$L_{\text{crit}} = \left( \frac{B}{0.688 \times 10^{12} \text{ G}} \right)^{16/15} \times 10^{37} \text{ erg s}^{-1}, \qquad (4)$$

which holds for typical parameters (see Eq. (32) of Becker et al. 2012, for more details). From the colour intensity diagram (see Fig. 4) we can set an upper limit to the critical transition and thus the magnetic field of the NS. For $L_{\text{crit}}$ of ~$10^{37}$ erg s$^{-1}$ we found an upper limit for the magnetic field of about $0.7 \times 10^{12}$ G. This quantitative estimate seems to agree with the lack of any evident CRSF in the broadband spectrum. Such a feature would appear at energies ($E_{\text{CRSF}}$) and be related to the magnetic field through the following relation:

$$B_{\text{CRSF}} = (1+z) \left( \frac{E_{\text{CRSF}}}{11.57 \text{ keV}} \right) \times 10^{12} \text{ G}, \qquad (5)$$

where $z \sim 0.3$ is the gravitational redshift from the NS, which is related to the NS compactness (e.g., Riley et al. 2019). The upper limit for the critical luminosity would thus translate into an upper limit for the $E_{\text{CRSF}}$ of ~6 keV. Because the spectrum at lower energies is complex and there is soft excess, it is thus not surprising that we did not detect any CRSF in the spectrum. This is also consistent with the fact that no accreting pulsar has a detected electron CRSF below 8 keV, and claims of such detections in some systems still wait to be confirmed (Staubert et al. 2019).

The pulse shape of the pulse profile is single peaked (Figs. 5 and A.1). At higher energies (>40 keV), the pulse profile is sharper, with a triangular shape, while it appears to be broader in the soft band, with weak evidence of a secondary peak. We note the similarity of the pulse profile in the soft NuSTAR band (see Fig. 5) to the profiles obtained from the NICER monitoring (see Fig. A.1). A single-peaked pulse profile at what otherwise appears to be close to or above the critical limit for the super-critical regime is quite unusual because complex profiles are the norm for these bright BeXRBs (e.g., Epili et al. 2017;
Koliopanos & Vasilopoulos 2018). Further investigation of the pulse profiles with physical models (e.g., Cappallo et al. 2017; Mushtukov et al. 2018) could provide further information of this somewhat puzzling feature.

## 5. Conclusions

We have analysed broadband spectra from the 2021 outburst of SXP 15.6. We did not identify any CRSF that could provide a direct measurement of the NS magnetic field. Nevertheless, the lack of such a feature does not exclude its presence because a CRSF could also be quite weak or depend on the pulse phase (e.g., Tiengo et al. 2013). In the broadband spectra, we found evidence of a soft spectral component that could be associated with an accretion disk, but its parameters are poorly constrained. An alternative explanation is that the soft excess is a result of partial absorption, and we do not favour one model over the other. The evolution of the spectral properties during the 2021 outburst is consistent with an accretion column above the NS and with the system accreting close to or above the critical limit. Finally, we did not measure any secular evolution of the spin period of the pulsar between 2016 *Chandra* and 2021 *NuSTAR* observations, which is consistent with the findings of Coe et al. (2022). Nevertheless, there is evident modulation in the *NICER* monitoring data that is consistent with orbital motion of the binary. Modelling of the orbital motion and intrinsic spin-up during the outburst enabled us to constrain the magnetic field strength and the orbital parameters of the system. All the derived quantitative and qualitative results consistently provide indirect constraints on the NS magnetic field strength. Because a cyclotron line is lacking, the most reliable measurement comes from the intrinsic spin up, based on which we find a value of $5 \times 10^{11}$ G with an uncertainty of about a factor of ~3.

*Acknowledgements.* The authors would like to thank the anonymous referee for the constructive report that helped to greatly improve the manuscript. GV would like to thank P.S. Ray for his comments, suggestions and advice on using PINT. This work was supported by NASA through the *NICER* mission and the Astrophysics Explorers Program. Facilities: *NICER*, *NuSTAR*. We acknowledge the use of public data from the *Swift* data archive.

## References

Anastasopoulou, K., Zezas, A., Steiner, J. F., & Reig, P. 2022, MNRAS, 513, 1400

Arnaud, K. A. 1996, in Astronomical Data Analysis Software and Systems V, eds. G. H. Jacoby, & J. Barnes, ASP Conf. Ser., 101, 17

Asami, F., Enoto, T., Iwakiri, W., et al. 2014, PASJ, 66, 44

Becker, P. A., & Wolff, M. T. 2007, ApJ, 654, 435

Becker, P. A., Klochkov, D., Schönherr, G., et al. 2012, A&A, 544, A123

Buccheri, R., Bennett, K., Bignami, G. F., et al. 1983, A&A, 128, 245

Buchner, J. 2019, PASP, 131, 108005

Buchner, J. 2021, J. Open Source Softw., 6, 3001

Burrows, D. N., Hill, J. E., Nousek, J. A., et al. 2005, Space Sci. Rev., 120, 165

Campana, S., Stella, L., Mereghetti, S., & de Martino, D. 2018, A&A, 610, A46

Cappallo, R., Laycock, S. G. T., & Christodoulou, D. M. 2017, PASP, 129, 124201

Cash, W. 1979, ApJ, 228, 939

Coe, M. J., Evans, P. A., Kennea, J. A., et al. 2021, The ATel, 15054, 1

Coe, M. J., Monageng, I. M., Kennea, J. A., et al. 2022, MNRAS, 513, 5567

Dickey, J. M., & Lockman, F. J. 1990, ARA&A, 28, 215

Endo, T., Nagase, F., & Mihara, T. 2000, PASJ, 52, 223

Epili, P., Naik, S., Jaisawal, G. K., & Gupta, S. 2017, MNRAS, 472, 3455

Evans, P. A., Beardmore, A. P., Page, K. L., et al. 2007, A&A, 469, 379

Evans, P. A., Beardmore, A. P., Page, K. L., et al. 2009, MNRAS, 397, 1177

Frank, J., King, A., & Raine, D. J. 2002, Accretion Power in Astrophysics, 3rd edn. (Cambridge: Cambridge University Press)

Fruscione, A., McDowell, J. C., Allen, G. E., et al. 2006, Proc. SPIE, 6270, 62701V

Fulton, B. J., Petigura, E. A., Blunt, S., & Sinukoff, E. 2018, PASP, 130, 044504

Gehrels, N., Chincarini, G., Giommi, P., et al. 2004, ApJ, 611, 1005

Gendreau, K. C., Arzoumanian, Z., & Okajima, T. 2012, in The Neutron star Interior Composition ExploreR (NICER): an Explorer Mission of Opportunity for Soft x-ray Timing Spectroscopy, SPIE Conf. Ser., 8443, 844313

Gendreau, K. C., Arzoumanian, Z., Adkins, P. W., et al. 2016, in The Neutron star Interior Composition Explorer (NICER): Design and Development, SPIE Conf. Ser., 9905, 99051H

Graczyk, D., Pietrzyński, G., Thompson, I. B., et al. 2014, ApJ, 780, 59

Haberl, F., & Sturm, R. 2016, A&A, 586, A81

Harrison, F. A., Craig, W. W., Christensen, F. E., et al. 2013, ApJ, 770, 103

Hickox, R. C., Narayan, R., & Kallman, T. R. 2004, ApJ, 614, 881

Huppenkothen, D., Bachetti, M., Stevens, A. L., et al. 2019, ApJ, 881, 39

Jaisawal, G. K., & Naik, S. 2016, MNRAS, 461, L97

Jaisawal, G. K., Naik, S., & Chenevez, J. 2018, MNRAS, 474, 4432

Jaisawal, G. K., Naik, S., & Epili, P. 2016, MNRAS, 457, 2749

Karaferias, A.S., Vasilopoulos, G., Petropoulou, M., et al. 2022, MNRAS, submitted

Kennea, J. A., Coe, M. J., Evans, P. A., et al. 2016, The ATel, 9362, 1

Kennea, J. A., Coe, M. J., Evans, P. A., Waters, J., & Jasko, R. E. 2018, ApJ, 868, 47

Koliopanos, F., & Vasilopoulos, G. 2018, A&A, 614, A23

Kubota, A., Tanaka, Y., Makishima, K., et al. 1998, PASJ, 50, 667

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

Lutovinov, A. A., & Tsygankov, S. S. 2009, Astron. Lett., 35, 433

Maitra, C., Paul, B., Haberl, F., & Vasilopoulos, G. 2018, MNRAS, 480, L136

Martin, R. G., & Franchini, A. 2021, ApJ, 922, L37

Martin, R. G., Nixon, C., Armitage, P. J., Lubow, S. H., & Price, D. J. 2014, ApJ, 790, L34

McBride, V. A., González-Galán, A., Bird, A. J., et al. 2017, MNRAS, 467, 1526

Müller, S., Ferrigno, C., Kühnel, M., et al. 2013, A&A, 551, A6

Mushtukov, A. A., Verhagen, P. A., Tsygankov, S. S., et al. 2018, MNRAS, 474, 5425

Okazaki, A. T., Hayasaki, K., & Moritani, Y. 2013, PASJ, 65, 41

Pan, Y. Y., Li, Z. S., Zhang, C. M., & Zhong, J. X. 2022, MNRAS, 513, 6219

Parfrey, K., Spitkovsky, A., & Beloborodov, A. M. 2016, ApJ, 822, 33

Postnov, K. A., Gornostaev, M. I., Klochkov, D., et al. 2015, MNRAS, 452, 1601

Prigozhin, G., Gendreau, K., Foster, R., et al. 2012, in Characterization of the Silicon Drift Detector for NICER instrument, SPIE Conf. Ser., 8453, 845318

Reig, P. 2011, Ap&SS, 332, 1

Reig, P., & Nespoli, E. 2013, A&A, 551, A1

Remillard, R. A., Loewenstein, M., Steiner, J. F., et al. 2022, AJ, 163, 130

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

Russell, S. C., & Dopita, M. A. 1992, ApJ, 384, 508

Skilling, J. 2004, in American Institute of Physics Conference Series, eds. R. Fischer, R. Preuss, & U. V. Toussaint, 735, 395

Staubert, R., Trümper, J., Kendziorra, E., et al. 2019, A&A, 622, A61

Steiner, A. W., Gandolfi, S., Fattoyev, F. J., & Newton, W. G. 2015, Phys. Rev. C, 91, 015804

Sturm, R., Haberl, F., Vasilopoulos, G., et al. 2014, MNRAS, 444, 3571

Sugizaki, M., Mihara, T., Nakajima, M., & Makishima, K. 2017, PASJ, 69, 100

Tiengo, A., Esposito, P., Mereghetti, S., et al. 2013, Nature, 500, 312

Treiber, H., Vasilopoulos, G., Bailyn, C. D., et al. 2021, MNRAS, 503, 6187

Tsygankov, S. S., Doroshenko, V., Mushtukov, A. A., et al. 2020, A&A, 637, A33

Vasilopoulos, G., Maggi, P., Haberl, F., et al. 2013, A&A, 558, A74

Vasilopoulos, G., Haberl, F., Delvaux, C., Sturm, R., & Udalski, A. 2016, MNRAS, 461, 1875

Vasilopoulos, G., Haberl, F., & Maggi, P. 2017a, MNRAS, 470, 1971

Vasilopoulos, G., Zezas, A., Antoniou, V., & Haberl, F. 2017b, MNRAS, 470, 4354

Vasilopoulos, G., Haberl, F., Carpano, S., & Maitra, C. 2018a, A&A, 620, L12

Vasilopoulos, G., Maitra, C., Haberl, F., Hatzidimitriou, D., & Petropoulou, M. 2018b, MNRAS, 475, 220

Vasilopoulos, G., Petropoulou, M., Koliopanos, F., et al. 2019, MNRAS, 488, 5225

Vasilopoulos, G., Ray, P. S., Gendreau, K. C., et al. 2020, MNRAS, 494, 5350

Verner, D. A., Ferland, G. J., Korista, K. T., & Yakovlev, D. G. 1996, ApJ, 465, 487

Wang, Y.-M. 1995, ApJ, 449, L153

Wilms, J., Allen, A., & McCray, R. 2000, ApJ, 542, 914

Wilson-Hodge, C. A., Malacaria, C., Jenke, P. A., et al. 2018, ApJ, 863, 9

Zimmerman, E. R., Narayan, R., McClintock, J. E., & Miller, J. M. 2005, ApJ, 618, 832

Zolotukhin, I. Y., Bachetti, M., Sartore, N., Chilingarian, I. V., & Webb, N. A. 2017, ApJ, 839, 125
## Appendix A: Additional figures and tables

<table><caption>Table A.1. Power-law fit to NICER data.</caption><thead><tr><th>OBSID<br>-</th><th>MJD<br>d</th><th>N<sub>H</sub><sup>a</sup><br>×10<sup>22</sup> cm<sup>-2</sup></th><th>Γ<br>-</th><th>F<sub>X</sub> (0.3-10.0 keV)<br>10<sup>-11</sup> erg cm<sup>-2</sup> s<sup>-1</sup></th><th>χ<sup>2</sup><sub>red</sub><br>-</th><th>DOF<br>-</th></tr></thead><tbody><tr><td>4202430101</td><td>59538.5</td><td>0.0263<sup>+0.016</sup><sub>-0.016</sub></td><td>1.17<sup>+0.03</sup><sub>-0.03</sub></td><td>3.24<sup>+0.07</sup><sub>-0.07</sub></td><td>1.15</td><td>424</td></tr><tr><td>4202430102</td><td>59539.1</td><td>&lt;0.05</td><td>1.03<sup>+0.02</sup><sub>-0.02</sub></td><td>3.67<sup>+0.09</sup><sub>-0.09</sub></td><td>1.22</td><td>351</td></tr><tr><td>4202430103</td><td>59540.5</td><td>&lt;0.018</td><td>1.12<sup>+0.03</sup><sub>-0.02</sub></td><td>3.19<sup>+0.08</sup><sub>-0.08</sub></td><td>1.11</td><td>379</td></tr><tr><td>4202430104</td><td>59541.5</td><td>&lt;0.007</td><td>1.06<sup>+0.02</sup><sub>-0.02</sub></td><td>3.29<sup>+0.06</sup><sub>-0.05</sub></td><td>1.24</td><td>397</td></tr><tr><td>4202430105</td><td>59542.5</td><td>&lt;0.4</td><td>1.14<sup>+0.07</sup><sub>-0.05</sub></td><td>2.40<sup>+0.12</sup><sub>-0.13</sub></td><td>1.16</td><td>154</td></tr><tr><td>4202430106</td><td>59543.6</td><td>&lt;0.1</td><td>1.06<sup>+0.02</sup><sub>-0.02</sub></td><td>2.66<sup>+0.05</sup><sub>-0.07</sub></td><td>1.08</td><td>328</td></tr><tr><td>4202430107</td><td>59544.4</td><td>&lt;0.015</td><td>1.12<sup>+0.03</sup><sub>-0.02</sub></td><td>2.61<sup>+0.05</sup><sub>-0.07</sub></td><td>0.95</td><td>364</td></tr><tr><td>4202430108</td><td>59545.5</td><td>&lt;0.03</td><td>1.02<sup>+0.08</sup><sub>-0.08</sub></td><td>2.77<sup>+0.18</sup><sub>-0.14</sub></td><td>0.86</td><td>42</td></tr><tr><td>4202430109</td><td>59550.7</td><td>&lt;0.006</td><td>0.87<sup>+0.02</sup><sub>-0.02</sub></td><td>2.19<sup>+0.05</sup><sub>-0.05</sub></td><td>1.18</td><td>404</td></tr><tr><td>4202430110</td><td>59551.7</td><td>&lt;0.007</td><td>0.69<sup>+0.04</sup><sub>-0.04</sub></td><td>2.12<sup>+0.10</sup><sub>-0.10</sub></td><td>1.14</td><td>269</td></tr><tr><td>4202430111</td><td>59552.5</td><td>&lt;0.003</td><td>0.88<sup>+0.02</sup><sub>-0.02</sub></td><td>2.03<sup>+0.06</sup><sub>-0.06</sub></td><td>1.19</td><td>355</td></tr><tr><td>4202430112</td><td>59553.1</td><td>&lt;0.11</td><td>0.91<sup>+0.12</sup><sub>-0.10</sub></td><td>1.28<sup>+0.13</sup><sub>-0.13</sub></td><td>0.90</td><td>46</td></tr><tr><td>4202430113</td><td>59555.5</td><td>&lt;0.11</td><td>0.92<sup>+0.10</sup><sub>-0.09</sub></td><td>1.12<sup>+0.09</sup><sub>-0.09</sub></td><td>1.03</td><td>103</td></tr><tr><td>4202430114</td><td>59556.7</td><td>&lt;0.04</td><td>0.69<sup>+0.07</sup><sub>-0.07</sub></td><td>1.42<sup>+0.09</sup><sub>-0.13</sub></td><td>0.97</td><td>96</td></tr><tr><td>4202430116</td><td>59559.5</td><td>&lt;0.07</td><td>0.99<sup>+0.08</sup><sub>-0.07</sub></td><td>2.24<sup>+0.14</sup><sub>-0.14</sub></td><td>1.03</td><td>157</td></tr><tr><td>4202430117</td><td>59560.8</td><td>&lt;0.06</td><td>0.92<sup>+0.10</sup><sub>-0.07</sub></td><td>1.87<sup>+0.14</sup><sub>-0.16</sub></td><td>0.87</td><td>115</td></tr><tr><td>4202430118</td><td>59561.6</td><td>&lt;0.10</td><td>0.90<sup>+0.09</sup><sub>-0.09</sub></td><td>1.45<sup>+0.12</sup><sub>-0.11</sub></td><td>1.09</td><td>116</td></tr><tr><td>4202430119</td><td>59563.3</td><td>&lt;0.008</td><td>0.97<sup>+0.04</sup><sub>-0.04</sub></td><td>1.98<sup>+0.10</sup><sub>-0.09</sub></td><td>1.11</td><td>171</td></tr><tr><td>4202430121</td><td>59571.5</td><td>&lt;0.13</td><td>0.74<sup>+0.17</sup><sub>-0.16</sub></td><td>2.14<sup>+0.4</sup><sub>-0.3</sub></td><td>1.41</td><td>15</td></tr><tr><td>4202430141</td><td>59602.6</td><td>&lt;0.07</td><td>1.09<sup>+0.11</sup><sub>-0.08</sub></td><td>2.15<sup>+0.18</sup><sub>-0.18</sub></td><td>1.05</td><td>50</td></tr></tbody></table>

<footnote>(a) Column density intrinsic to the SMC and the source. Most values are not well constrained and are consistent with an upper limit. For the fit, the Galactic column density was fixed to a value of 4×10<sup>20</sup> cm<sup>-2</sup>.</footnote>
<figure><img src="image_9.png" /></figure>

**Fig. A.1.** Pulse profiles (0.5–8.0 keV) of all NICER observations where significant pulsations were detectable. The labels within each panel denote the last three digits of the NICER obsid number, i.e. 4202430XXX. All profiles are normalised to the average intensity, and the minimum is shifted to zero phase. Shaded horizontal regions denote the limit of the statistically significant variability above the constant level hypothesis. For comparison, the 2016 pulse profile as measured by Chandra is plotted with a dotted grey line in each panel.
<figure><img src="image_10.png" /><figcaption>**Fig. A.2.** Corner plot for SXP 15.6 using the model described in the text. We plot the logarithm of the polar magnetic field strength *B* in G, the eccentricity (*e*), the longitude of periastron in degrees (*ω*), the projected semi-major axis in light seconds (*a sin i*), the time of a mean longitude of 90 degrees *T*<sub>π/2</sub> and the pulsar frequency *F*<sub>0</sub> in mHz. The reference epoch for *F*<sub>0</sub> is the start of the NICER monitoring, i.e. MJD 59538.5.</figcaption></figure>