Estimating the spin of the black hole candidate MAXI J1659-152 with the X-ray continuum-fitting method

The PCA spectra are extracted according to the standard procedures described in RXTE Cook Book, based on the New PCA Tools⁴⁴4https://heasarc.gsfc.nasa.gov/docs/xte/recipes2/Overview.html with HEAsoft⁵⁵5https://heasarc.gsfc.nasa.gov/docs/software/heasoft/download.html v6.28. The whole procedure uses the PCA calibration files⁶⁶6https://heasarc.gsfc.nasa.gov/docs/heasarc/caldb/caldb_supported_missions.html v20200515. We select the data from the top layer of the best-calibrated detector PCU2 in Standard 2 mode. A dead time correction is taken into account. The latest bright-source background model is used to generate background spectra. After generating standard data products, we also apply calibration correction pcacorr ⁷⁷7http://www.srl.caltech.edu/personnel/javier/crabcorr/index.html which can reduce the systematic error that accounts for the uncertainties in the instrumental responses (García et al. 2014). We adopt a $0.1\%$ systematic error. In this work, we choose the energy band of 3.0-45.0 keV. Following our stringent screening requirements on $l$ and $f_{\mathrm{sc}}\lesssim 25\%$, we select 9 spectra (SP1-SP9). To account for the significant flux-normalization differences between X-ray missions, we follow the previous work (Steiner et al. 2010) in attempt to standardize the calibration using the Crab as a reference. We opt for observations of Crab that are close to the target dates for SP1-SP9, and adopt $\Gamma=2.1$, $N=9.7\text{photons}s^{-1}\text{keV}^{-1}$ (Toor & Seward 1974) as the standard values. For spectra SP1-SP8, we determine that the normalization-correction coefficient $C_{\mathrm{TS}}$ is $1.128$ and the slope difference $\Delta\Gamma_{\mathrm{TS}}$ is $0.021$. For SP9, $C_{\mathrm{TS}}$ and $\Delta\Gamma_{\mathrm{TS}}$ are 1.123, 0.018, respectively. This standardization is implemented in XSPEC via crabcor. We list the detailed information about SP1-SP9 in Table 1. It is worth noting that MAXI J1659 is a bright source, whose peak flux reaches about 300 mCrab (Kalamkar et al. 2011). Jahoda et al. (2006) put forward that the threshold of pile-up significance of RXTE/PCA is 10,000 $\mathrm{cts}~{}\mathrm{s}^{-1}~{}\mathrm{PCU}^{-1}$. In other words, below this threshold, the photon pile-up is not significant. SP1-SP9 are all far below this threshold.

With RXTE/PCA coverage beginning at 3 keV, the hydrogen column density ($N_{\mathrm{H}}$) can not be well constrained. So we turn to XMM-Newton data which covers the low-energy band in which $N_{\mathrm{H}}$ is most prominent. We employ EPIC-pn timing-mode data during observation (ObsID 0656780601) which spanned $\sim$23 ks during the X-ray outburst. Data reduction is carried out using the Science Analysis System (SAS)⁸⁸8https://www.cosmos.esa.int/web/xmm-newton/download-and-install-sas v18.0 with the latest calibration files. We follow the standard procedures and remove the influence caused by pile-up. In addition, we include $0.5\%$ systematic error as suggested by the XMM-Newton team⁹⁹9https://heasarc.gsfc.nasa.gov/docs/xmm/sl/epic/image/sas_cl.html. Crab calibration is also made for XMM-Newton data. Implementing TBabs*(diskbb+powerlaw), eventually, we get a $N_{\mathrm{H}}=3.22\times 10^{21}\mathrm{cm}^{-2}$ via fitting the spectrum over the energy range of 0.7-12.0 keV with a reduced chi-square $\chi_{v}^{2}=1.179~{}(2451.8/2080.0)$. In subsequent work, all $N_{\mathrm{H}}$ are fixed to this value.

As shown in Figure 1, MAXI/GSC has been monitoring MAXI J1659 for up to eleven years. The horizontal axis represents MJD and the vertical axis represents the flux in 2.0-20.0 keV detected by MAXI/GSC. As can be seen from Figure 1, MAXI J1659 has exhibited just one major outburst and shown no other signs of activity. As displayed in Figure 2, the hardness-intensity diagram (HID) shows a classical “q”-like shape, which is a typical for an outbursting black hole LMXB.

Firstly, we present the non-relativistic case, which consists of multiple blackbody disk components diskbb (Mitsuda et al. 1984; Makishima et al. 1986). The composite non-relativistic model crabcor*TBabs*(simpl*diskbb) is used to fit the data for the 9 screened RXTE/PCA spectral. The model crabcor and model simpl have been mentioned earlier (see Section 2 and Section 1 respectively). The photoionization cross-sections of the interstellar medium (ISM) in model TBabs are based on Verner et al. (1996) and the abundances are based on Wilms et al. (2000), and we set $N_{\mathrm{H}}$ to $3.22\times 10^{21}\mathrm{cm}^{-2}$ (see Section 2). The best-fitting results of SP1-SP9 are listed in Table 2. As can be seen, these 9 spectra are well fitted and the reduced chi-square $\chi_{\nu}^{2}$ of SP1-SP9 are basically maintained in the vicinity of 1. Figure 3 shows the fitting of SP1 and SP4 as representatives. We can discover that the model fitted well without any significant residual.

We next move to the fully-relativistic accretion-disk model kerrbb2 (McClintock et al. 2006). kerrbb2 merges two different disk models, bhspec and kerrbb. Specifically, bhspec is used to determine the value of the spectral hardening factor $f\equiv T_{\mathrm{col}}/T_{\mathrm{eff}}$ (also known as the color correction factor, Davis et al. 2005), while kerrbb employs ray-tracing computations to model the disk (Li et al. 2005). We first use bhspec¹⁰¹⁰10http://people.virginia.edu/~swd8g/xspec.html to compute spectral-hardening look-up tables according to two representative values of the viscosity parameter ($\alpha=0.1$ and $\alpha=0.01$). It is worth noting that the calculations with the model BHSPEC are made for spins ranging from -1 to 1. During kerrbb2 fitting, the $f$-table is read in and used to automatically set $f$’s value; otherwise the performance is identical to the original kerrbb. The relativistic model can be expressed as crabcor*TBabs*(simpl*kerrbb2). By default we adopt $\alpha=0.1$ and note that $\alpha=0.01$ returns slightly larger values for spin.

We adopt the revised system parameters from improved by Torres et al. (2021), and consider inclination angles of $70^{\circ}$ and $80^{\circ}$. The mass is restricted to $5.7\pm 1.8~{}\mathrm{M}_{\odot}$ (1$\sigma$) and $4.9\pm 1.6~{}\mathrm{M}_{\odot}$ (1$\sigma$) respectively. Following Torres et al. (2021), when the distance is $D=6\pm 2$\text  kpc, the center value of the measured mass falls within the range of 4.3-6.0$~{}\mathrm{M}_{\odot}$ calculated based on formula 5 by Yamaoka et al. (2012). In addition, this distance is also consistent with Kuulkers et al. (2013). Therefore, we adopt the distance of $D=6\pm 2$$\text{~{}kpc}$ in this paper.

The self-irradiation of the disk (rflag=1) and the effect of limb-darkening (lflag=1) are taken into account. And we set the torque at the inner boundary of the disk to be zero (eta=0). For both bracketing sets of system parameters: $M=5.7\pm 1.8~{}\mathrm{M}_{\odot}$ (1$\sigma$), $i=70.0^{\circ}$, $D=6\pm 2$ kpc (hereafter group 1) and $M=4.9\pm 1.6~{}\mathrm{M}_{\odot}$ (1$\sigma$), $i=80.0^{\circ}$, $D=6\pm 2$ kpc (hereafter group 2), we show our results in Table 3. Above the double horizontal lines is group 1, below the double horizontal lines is group 2. Group 1 (70^∘) would have MAXI J1659 with a moderate positive spin, whereas group 2 (80^∘) finds an extreme retrograde spin as most likely.

Although in this paper, we focus upon results for $\alpha=0.1$ as our default, for the alternative case of $\alpha=0.01$ we have performed the same calculations. Those results are in Table 4. The errors that appear in Table 3 and 4 are only due to the statistical uncertainties estimated by XSPEC for $90\%$ confidence. Next, we will discuss the error from uncertainties of input parameters: $M$, $i$, $D$ in detail.

The dominant source of uncertainty in CF is the uncertainty of the input parameters: $M$, $i$, $D$. In order to determine those uncertainties, Monte Carlo (MC) methods are often used (see, e.g., Liu et al. 2008; Gou et al. 2009). We follow these works, but differ in that, unlike those, here there are two distinct sets of input system parameters (group 1 and group 2). We first consider these groups separately.

For each individual spectrum in SP1-SP9, assuming that randomly generated data points are independent of each other and obey the Gaussian distribution, we have generated 3000 data points of M and D respectively. $i=70.0^{\circ}$ ($i=80.0^{\circ}$), and these 3000 data points constituted 3000 data sets as input parameters. We next compute the look-up tables of $f$ for these 3000 data sets. Lastly, we use the composite model crabcor*TBabs*(simpl*kerrbb2) (see Section 2) to fit 3000 data sets to obtain the histogram of spin, so as to determine the errors. The MC-determined error analysis of each spectrum is shown in Figure 4. The histograms of $a_{*}$ for SP1-SP9 is shown in Figure 5. As shown in Figure 5, group 1 with higher-mass and lower inclination yields spin $a_{*}=0.21_{-0.20}^{+0.14}$ (1$\sigma$); group 2 gives $a_{*}=-0.79_{-0.20}^{+0.31}$ (1$\sigma$).

We now combine groups 1 and 2 in order to establish a net constraint on spin. We are interested in assessing the lower and upper bounds on spin which are mostly constrained by groups 2 and 1, respectively. We mark the appropriate 1-$\sigma$ limit for each spectrum in Figure 4, and for the composite result in Figure 5. In total, we find that the spin of MAXI J1659 is poorly constrained, with an allowable range $-1<a_{*}\lesssim 0.35$ in 1 sigma interval. The 90% upper limit on spin is 0.44 from group 1.

In order to assess whether the hydrogen column density significantly affects the spin results, we explore varying the value of $N_{\mathrm{H}}$ from 0.322 to 0.2, 0.4, 0.5 in units of $10^{22}\mathrm{cm}^{-2}$. These alternates are chosen as round numbers covering the range of values given by Kennea et al. (2011). The fitting is systematically affected for all the spectra in the same way. We illustrate this difference by showing such results for SP1 in Table 5. When $N_{\mathrm{H}}$ increases from $2\times 10^{21}\mathrm{cm}^{-2}$ to $5\times 10^{21}\mathrm{cm}^{-2}$, for group 1, $a_{*}$ varies from 0.210 to 0.134, with $\Delta a_{*}$ equaling to 0.076. And for group 2, $a_{*}$ changes from -0.834 to -0.957, and $\Delta a_{*}$ is 0.123. As expected, a slight increase of $N_{\mathrm{H}}$, causes the inferred spin decrease. Compared to the dynamical sources of uncertainty in hand, changing the value of $N_{\mathrm{H}}$ has very minor effect on the final spin results.

In addition to fitting each spectrum individually with the relativistic model, we also consider the case that all the spectra are fitted simultaneously and looked into its effect. We linked the spin parameters among all 9 spectra, and letting the other fit parameters free. For the convenience of comparison, we also listed the results of simultaneous fits in Table 3 and 4 to group 1 and 2, and use a single horizontal line as a distinction. As a showcase, the simultaneous fit to all the spectra is also shown in Figure 6 (a). As we can see, the spin ranging from $-1<a_{*}\lesssim 0.33$ (1$\sigma$) obtained from MC result (Figure 6 (b)), and it is clear that the results are fully consistent with ones obtained from individual spectral fits (Figure 5), which is expected. Therefore, the effect of joint fits is negligible and our work uses the fit results from the individual spectra as our primary results.

Rout et al. (2020) construct a wider range of system parameter ($M$, $i$, $D$) grids based on an older set of measurement results. They found a simultaneous data set with XMM/EPIC-pn and RXTE/PCA overlapping on September 28, 2010. Fitting the two data sets together, they simultaneously adopt the X-ray reflection fitting and continuum-fitting method. On this basis, they use meaningful values of mass accretion rate $\dot{M}$ to constrain the spin of MAXI J1659. In their results of the X-ray continuum-fitting show that the lower limit of the spin is pegged at -0.998, while the upper limit is 0.4. However, it should be noted that the spectra in Rout et al. (2020) do not meet the $f_{\mathrm{sc}}\lesssim 25\%$ criterion. Different from Rout et al. (2020), we firstly screen the spectra with $f_{\mathrm{sc}}\lesssim 25\%$ for the X-ray continuum-fitting and made crabcor correction. In addition, we use the relativistic model kerrbb2 that allows the spectral hardening factor $f$ to change, which will be closer to the actual situation. And we do not ignore the influence of self-irradiation in kerrbb2. Based on the use of different spectra, different models and updated dynamical parameters, we find a spin ranging from $-1<a_{*}\lesssim 0.35$ (1$\sigma$). As reported by Rout et al. (2020), we also rule out extreme prograde spin.

To verify that a stronger hard tail does not introduce bias, we also check for possible correlation between the scattering fraction in the spectrum and the spin. We take SP1 as an example, and find for any individual spectrum, the contour map between black hole spin and $f_{\mathrm{sc}}$ shows a strong degenerate relation (see subgraph (a) of Figure 7). However, we find that as an ensemble which encompasses a range of $f_{\mathrm{sc}}$ values, the spin does not exhibit correlation (see subgraph (b) of Figure 7). Accordingly, we conclude that $f_{\mathrm{sc}}$ does not affect our spin results.

King & Kolb (1999) considered the effect of prograde accretion in changing the black hole mass and its spin (See their Figure 3). We also consider the influence of retrograde accretion in changing the black hole mass and spin in Appendix A. As stated above, Kuulkers et al. (2013) found that the companion of MAXI J1659 had an initial mass of about $1.5~{}\mathrm{M}_{\odot}$, which evolved to its current mass in about 4.6-5.7 billion years. Let us make a simple estimate. We assume that MAXI J1659 is a retrograde black hole with an extreme spin, approaching $a_{*}=-1$. And we suppose the initial mass of the black hole is $4~{}\mathrm{M}_{\odot}$. The outburst-averaged dimensionless mass accretion rate ($\dot{m}=\dot{M}/\dot{M}_{\mathrm{Edd}}$) of the standard thin disk model is usually in the range of 0.01-0.3 (Narayan & Yi 1995; Ohsuga et al. 2002), nonetheless, if we take into account that only one 7-month outburst has been observed from the source in the last 25 years, the time-averaged accretion rate $\dot{m}=2\times 10^{-4}$ would be more realistic. When $\dot{m}=2\times 10^{-4}$, the Eddington accretion rate ($\dot{M}_{\mathrm{Edd}}$) of a black hole with $4~{}\mathrm{M}_{\odot}$ is about $10^{-8}~{}\mathrm{M}_{\odot}~{}\mathrm{yr}^{-1}$. According to Equation A1, in 4.6-5.7 billion years, the accreted mass onto MAXI J1659 should fall within the range 0.0092-0.0114$~{}\mathrm{M}_{\odot}$ as $\dot{m}=$$2\times 10^{-4}$. In other words, $\Delta M/M_{\rm i}$ (defined in Appendix A) should fall within the range of 0.0023-0.00285. This obviously is insufficient to rule out the possibility that MAXI J1659 has an extreme negative spin (see Figure 8, two of the red dash lines represent 0.0023 and 0.00285, respectively).