© ESO 2018

1 Introduction

Recent observations suggest that some ultraluminous X-ray (ULX) sources contain a neutron star in their center; a typical coherent pulsation from magnetized pulsar was detected in the cases of M82X-2 (Bachetti et al. 2014), NGC 5907 ULX1 (Israel et al. 2017a), and NGC 7793 P13 (Fürst et al. 2016; Israel et al. 2017b). The mass of the central object may equal 1.4 solar masses, which implies that the accretion disk has to be highly super-Eddington. Such a consequence directly results from the high values of the observed X-ray luminosities from those systems, which always exceed the theoretical maximum for spherical infall onto a stellar-mass black hole (Roberts 2007; Liu et al. 2013). The broadband spectra of those sources obtained simultaneously by XMM-Newton and NuSTAR X-ray satellites are, in most cases, inadequately described by only one multi-temperature disk component. A second spectral component, such as thermal Comptonization or a blackbody from the neutron star surface, is needed to fully explain the spectral shape (Walton et al. 2013, 2018).

A theoretical explanation for the broad spectral shape was given by Różańska et al. (2017) in the case of low mass X-ray binaries (LMXB), where we presented spectral shape from a non-spherical system that contains a neutron star with an accretion disk. Our single model component contains emission from both regions with the effects of mutual occultation, which gives the proper model normalization. The final spectrum depends on the viewing angle over the whole energy range, and for the assumed emitting surface it is inversely proportional to the square of the distance to the source.

In this paper, we reduce archival broadband data of three sources: NGC 7793 P13, NGC 5907 ULX1, and Circinus ULX5 (hereafter P13, ULX1, and ULX5). The first two systems are two of the three neutron star ULXs known at present. The observations were taken with the XMM-Newton (Jansen et al. 2001) and NuSTAR (Harrison et al. 2013) telescopes simultaneously. Spectral analyses were recently done by Walton et al. (2018) for P13, by Israel et al. (2017b) for ULX1, and by Walton et al. (2013) for ULX5, where their final models contained several components. In this paper we show that our single model of non-spherical emission from a neutron star with the multi-color blackbody emission from an accretion disk fits the broadband spectra of P13, ULX1, and ULX5 perfectly. The normalization of our model constrains the distance to the source. In the case of P13, the resulting distance perfectly agrees with the distance determination based on the Cepheid method to the hosting galaxy NGC 7793 (Pietrzyński et al. 2010). Furthermore, the distance to ULX5 from our method also agrees with the distance found with the radial velocity method to Circinus galaxy (Koribalski et al. 2004).

2 Single model of emission and parameters

We assume that a ULX source is a non-spherical system containing a neutron star with an accretion disk around it. The global disk model is not specified here, and the total disk emission is parametrized by multi-color blackbody shape. This approach is correct since we do not know the exact accretion efficiency and we are dealing with objects where standard accretion does not account for the high luminosity. For this object geometry, the infinitesimal energy $d{\mathcal{F}}_{\nu }$ measured by a distant observer is defined as $d{\mathcal{F}}_{\nu }=I_{\nu }d\omega$, where dω is a solid angle in steradian (sr; Mihalas 1978). This formula is applicable for the flat space, as we assume in this paper. When the emitter is located close to a black hole, both general and special relativistic corrections should be taken into account, as in Fabian et al. (1989, Eq. A4). Integrating the above formula over the solid angle subtended by the source, we obtained the energy dependent intensity ${\mathcal{F}}_{\nu }$ as seen by the observer. This quantity diminishes with increasing distance, and it is not an intrinsic property of the source. In the case of non-spherical systems the observed intensity should be integrated over the emitting region individually depending on the source geometry. In this paper, we consider an emitting region with cylindrical symmetry.

In the case of emission from the whole system, i.e., a neutron star with an accretion disk around it, we have the contribution from different emitting parts (for details, see Różańska et al. 2017) and the final observed energy dependent intensity directed to the observer is computed analytically as $$\begin{array}{lll}{\mathcal{F}}_{\nu ,\text{All}}\hfill & =\hfill & {\left(\frac{1}{D}\right)}^2[\pi R_{\text{NS}}^2\left({\displaystyle {\int }_0^1{I}_{\nu }}\mu \text{d}\mu +{\displaystyle {\int }_{\cos {\theta }^{}\prime }^1{I}_{\nu }}\mu \text{d}\mu \right)\hfill \\ \hfill & \hfill & +2({\displaystyle {\int }_0^{R_{\text{NS}}}{I}_{\nu }}\sin {\theta }^{}\prime \sqrt{R_{\text{NS}}^2-x^2}\text{d}x\hfill \\ \hfill & \hfill & -{\displaystyle {\int }_0^{R_{\text{NS}}\sin {\theta }^{}\prime }{I}_{\nu }}\sqrt{R_{\text{NS}}^2{\sin }^2{\theta }^{}\prime -x^2}\text{d}x)\hfill \\ \hfill & \hfill & +\pi \sin {\theta }^{}\prime \left({\displaystyle {\int }_{R_{\text{in}}}^{R_{\text{out}}}{I}_{\nu }}R\text{d}R+{\displaystyle {\int }_{R_{\text{boost}}}^{R_{\text{out}}}{I}_{\nu }}R\text{d}R\right)],\hfill \end{array}$$(1)

where R_NS is the neutron star radius, D is the distance to the system, and I_ν is the specific intensity emitted from the source surface towards the observer; μ = cos θ, where θ is the angle between the direction of the light beam and the normal to the neutron star surface; θ′ is the viewing angle related to the disk inclination as i = 90° − θ′ (see Fig. 3 in Różańska et al. 2017); and x is the variable of integration.

The first two parts of Eq. (1) correspond to the emission from the neutron star surface, taking into account the occultation by the disk, while the second two parts describe the multi-color blackbody disk emission with eventual occultation by the neutron star. Due to mutual occultation, one part of the disk is integrated over the radius from the innermost stable orbit R_in all the way up to the outer radius R_out, while the second part of the disk is integrated from R_boost = R_NS∕ sin θ′ to the outer radius R_out, where R_boost is the radius up to which the neutron star covers the inner disk. The above formula was derived analytically by Różańska et al. (2017), where we have demonstrated that the broadband emitted spectrum from a non-spherical system depends on the viewing angle, since multi-color blackbody disk emission and neutron star emission both change with angle.

We assume that the neutron star radiation equals the blackbody intensity at the given effective temperature I_ν = B_ν(T_eff,NS). Furthermore, we assume that the emission at different disk radii equals the local Planck function I_ν = B_ν(T_eff(R)), with the effective temperature given by a standard multi-color blackbody formula, $T_{\text{eff}}(R)=T_{\text{in}}{(R/R_{\text{in}})}^{-p}$, where T_in is the inner disk temperatureand exponent p equals 3/4 for the standard disk. Following an earlier approach to the spectral fitting of ULXs (Walton et al. 2013; Bachetti et al. 2014; Fürst et al. 2016; Israel et al. 2017b; Walton et al. 2018), we do not connect here the inner disk temperature with an accretion rate since this relation strongly depends on the assumed global disk model and the accretion efficiency. Since the latter quantity is not well known and since we expect super-Eddington accretion in ULXs, we continue disk parametrization by the innerdisk temperature. Nevertheless, we note here that our model is useful for systems where the angle dependent specific intensity is given as the result of the radiative transfer calculations (Madej 1989, 1991; Hubeny et al. 2001; Davis et al. 2005; Różańska et al. 2011). We plan to implement atmospheric models in a future work.

For the purpose of comparing our model to the observed X-ray spectrum of ULX source P13, we constructed a grid of models for arbitrarily assumed parameters. A non-rotating neutron star has a canonical mass 1.4 M_⊙, radius 12 km, and 11 various effective temperatures, ranging from 2 × 10⁶ K to 4 × 10⁷ K. The disk local multi-color emission was computed assuming 11 inner disk temperatures ranging from 2.24 × 10⁶ to 1.26 × 10⁷ K. The grid of θ′ angles spans from 10° to 90°. The lowest value of this angle corresponds to an almost edge-on disk, and the highest value to a face-on disk.

For each disk model, we calculated a multi-color blackbody spectrum from R_in to R_out in the range 3–1000 R_Schw, where the R_Schw = 2GM_NS∕c². The inner disk radius, and therefore the inner disk temperature, can change due to the high value of the magnetic field, the strong boundary layer, and when the relativistic corrections are taken into account. We plan to include them in a future paper together with the full ray tracing procedure.

Our model normalized by the distance of D = 10 kpc was prepared as a table model in the Flexible Image Transport System (FITS) format (Wells et al. 1981), and for the purposes of this paper we call it neutron star multi-color blackbody (nsmcbb). The final table of total emission from the system is expressed in terms of four parameters, whose values will be determined during the fitting procedure: neutron star effective temperature T_eff,NS, inner disk temperature T_in, viewing angle θ′, and normalization N. Since the normalization of our model is proportional to 1∕D², the value of this normalization, which is determined from the fitting procedure, is the direct indicator of the distance to the system, which equals $10/\text{}\text{}\sqrt{N}$ kpc. Below we show that for P13 and ULX5 the distance derived from our fitting is in excellent agreement with independent distance measurement.

3 P13, ULX1, ULX5, and their X-ray observations

P13 is a variable ULX source that is thought to harbor a stellar black hole of mass less than 15 M_⊙ (Motch et al. 2014). However, recent studies by Fürst et al. (2016) and Israel et al. (2017b) have discovered that P13 hosts an accreting neutron starwith a spin period of 0.42 s. The source is a part of the binary system, where it circulates around B91a star of 18–23 M_⊙, with a period of 64 days (Motch et al. 2014). Previously, P13 was reported to reach luminosities from ~ 2 × 10³⁹, detected in 1979 by the Einstein satellite, to ~10⁴⁰ erg s⁻¹, reported recently with joint XMM-Newton and NuSTAR data by Walton et al. (2018). The source is unobscured with a galactic warm absorption value of N_H = 9.60 ± 0.01 × 10²⁰ cm⁻² (Israel et al. 2017b). On the other hand, Walton et al. (2018) divided this absorption between galactic absorption, with a value of neutral hydrogen estimated towards the object of N_H = 1.2 × 10²⁰ cm⁻² (Kalberla et al. 2005), and intrinsic source absorption resulting from the fit, N_H = 8 ± 1 × 10²⁰ cm⁻². P13 is located in the NGC 7793 galaxy which is a part of the Sculptor Group. The first distance measurement to the source was done by Karachentsev et al. (2003) as a distance to the galaxy itself and was estimated to be 3.91 ± 0.41 Mpc. Further distance derivation was done within the ARAUCARIA project (Gieren et al. 2005) where for the first time Cepheid variables were detected in the Sculptor Group. The distance to the NGC 7793 galaxy with the Cepheid method was found to be 3.4 ± 0.17 Mpc (Pietrzyński et al. 2010).

ULX1 is located in the NGC 5907 galaxy and is a part of binary system; its period is 78 days. Recent studies by Israel et al. (2017b) have discovered that ULX1 harbors a neutron star at its center, with spin period ranging from 1.13 s to 1.43 s. The first estimations of the distance from Earth amounted to 13 Mpc (Tully et al. 2009); however, later studies discovered that ULX1 may be located even farther form Earth, with its distance totaling 17 Mpc (Tully et al. 2016). Such a great distance results in peak luminosities reaching values up to 6 × 10⁴⁰ erg s⁻¹, for the distance 13 Mpc (Fürst et al. 2017), which is many times higher then the luminosity for a standard neutron star. The typical values of galactic warm absorption toward the ULX1 given by those authors range from 5.7 to 8.5 × 10²¹ cm⁻².

ULX5 is located in the outskirts of the Circinus galaxy, with its distance estimated to be 2.78 Mpc (Koribalski et al. 2004; distance to the galaxy via the radial velocity method). The first observations of ULX5 at hard X-rays were performed and analyzed by Walton et al. (2013). Depending on the fitted model, the values of galactic warm absorption toward the source were from 1 to 9.5 × 10²¹ cm⁻². Its peak luminosity was calculated to reach a value of 1.6 × 10⁴⁰ erg s⁻¹. This luminosity was calculated for higher distance of 4 Mpc measured with the Tully–Fisher method when Circinus galaxy was discovered (Freeman et al. 1977).

In this work we made use of coordinated NuSTAR (OBSID: 80201010002, 80001042002, 60002039002) and XMM-Newton (OBSID: 0781800101, 0729561301, 0701981001) observations of three sources: P13, ULX1, and ULX5. The observations were taken on 2016 May 20, 2014 July 09, and 2013 February 03, respectively. Their exposure times were as follows: for P13, 118 ks for NuSTAR and 22/46 ks for the EPIC-pn/EPIC-MOS cameras on board the XMM-Newton X-ray telescope; for ULX1, 65 ks and 38 ks/43 ks, and for ULX5, 54 ks and 24 ks/47 ks for the same instruments.

We reduced the NuSTAR data from both focal plane modules, FPMA and FPMB, by using NuSTAR Data Analysis Software in accordance with guidelines provided in the NuSTAR Data Analysis Software Guide (v1.9.2)¹. Calibration files were taken from the actual database CALDB v20170817 throughout the whole process. The NUPIPELINE tool was used in order to produce filtered event files, with standard filtering applied. We utilized the NUPRODUCTS tool to extract source products and instrumental response files from circular regions of radii 60″ for P13, 40″ for ULX1, and 55″ for ULX5. Background products were obtained from regions that were four times larger, located on the same detector chip as the sources. In addition to the basic science data we also included spacecraft science data following the guidelines outlined in the NuSTAR Data Analysis Software Guide (see also Walton et al. 2016). Therefore, we were able to increase our exposure time by approximately 10% and maximize the signal-to-noise (S/N) ratio. Additional data obtained by this action were merged with standard scientific data by running the FTOOLS task, ADDASCASPEC. Since we intend to fit broadband data, small errors in responses generated by this tool do not affect our results significantly. The final spectra, averaged over exposure time, were binned to at least 20 counts per bin for the data fitting process. The extracted normalized counts are presented in Fig. 1.

The XMM-Newton data were reduced with the XMM-Newton Science Analysis System (SAS) v16.0.0², following standard guidelines outlined in the science analysis threads. We generated the calibrated and concatenated event lists by running the EPCHAIN task for EPIC-pn and the EMCHAIN task for both EPIC-MOS modules. Next, we created an event file and then subjected it to filtering for background flaring with the help of SAS tasks evselect and tabgtigen. To produce the final event files with spectra, the source products were extracted from circular regions of radii: 36″ for P13, 36″/45″ for ULX1 (PN/MOS), and 22.5″/27″ for ULX5 (PN/MOS), and background products were extracted from regions that were nine times larger on the same CCD chip. All regions were extracted in so as to avoid the CCD borders. In evselect, we used the filters FLAG == 0 PATTERN < 4 for EPIC-pn and FLAG == 0 PATTERN < 12 for the EPIC-MOS cameras. The appropriate response files and ancillary files were generated using the SAS commands rmfgen and arfgen, respectively. Finally, we combined the spectra from both EPIC-MOS detectors by running the FTOOLS task ADDASCASPEC. The final EPIC-pn spectrum was binned to at least 10 counts per energy bin, while spectra from the two MOS cameras were binned to at least 20 counts per energy bin for the data fitting process. The extracted normalized counts are presented in Fig. 1.

Fig. 1 Normalized counts from all detectors used in our spectral fitting analysis for P13, ULX1, and ULX5. Black and red crosses correspond to the XMM-Newton detectors EPIC-pn and EPIC-MOS. Green and blue crosses are data from NuSTAR FPMA and FPMB, respectively. Black solid lines are the best fit models.

Fig. 2 Upperpanels: unfolded photon spectra from all detectors used in our spectral fitting analysis. Lower panels: ratio i.e. data divided by model. All colors have the same meaning as in Fig. 1.

4 Spectral analysis

We performed spectral fitting of P13, ULX1, and ULX5 data with a single model of emission nsmcbb, presented in Sect. 2 with the use of XSPEC fitting package, version 12.9.0³. In case of each source, all sets of data were fitted simultaneously with nsmcbb model multiplied by galactic absorption model tbnew. All values of metal abundances in the tbnew model were frozen. Only hydrogen column density was fitted in this model. Table 1 presents all fitted parameters and uncertainties given at the 90% confidence level.

The quality of fits is excellent since the data are very good, giving 1245, 859, and 762 degrees of freedom (d.o.f.) for each source, respectively. The reducedstatistics are χ²∕d.o.f. = 1.08 for P13, 1.01 for ULX1, and 1.14 for ULX5, indicating that our single model of emission from the non-spherical system perfectly agrees with the data of each source, which we illustrate in Fig. 2. Unfolded spectra from all detectors are broad and correctly agree with the models up to 20 keV. There is a small deviation for higher energies resulting in the value of the ratio (data/model) reaching 4, but this feature is always present in broadband spectral analysis of ULX sources. Even in the case of multi-component spectral fitting of P13 by Walton et al. (2018), this ratio is between 2 and 3.5, which may suggest that the hard energy tail is still not well detected by us, and is also not well modeled.

The fitted values of galactic absorption are consistent with previous estimations given in Sect. 3, and since those values are quite low they do not affect other parameters. The neutron star effective temperatures are high, giving rise to the hard energy bump in the observed spectrum. It is clearly seen in the unfolded energy spectrum plotted in E * F_E versus E in Fig. 3. For P13 and ULX5, the resulting inclination angles suggest that the whole systems are observed edge-on, and the neutron star is strongly covered by the accretion disk. The case is different for ULX1, whose orientation is almost face-on, and therefore occultation is not very strong. For all sources considered in this paper, the integration over such emitting area fully explains the shape of the X-ray observed spectrum.

Since the emission from the non-spherical region fits observations, we can calculate distances to all sources from the model normalization, listed in Table 2. In the case of ULX1, the previous distance estimations are uncertain, and give two values 13 and 17 Mpc. The distance to ULX1 resulting from our modeling is a factor of two lower: $D={6.55}_{-0.81}^{+0.69}$ Mpc.

In the case of ULX5, the situation is much better since the distance to the source obtained by our method $D={2.60}_{-0.03}^{+0.05}$ Mpc agrees with the value of 2.78 Mpc obtained from the radial velocity method by Koribalski et al. (2004) (assuming Hubble constant H₀ = 75 km s⁻¹ Mpc). However, the best agreement we found for P13, the obtained value of $D={3.41}_{-0.10}^{+0.11}$ Mpc from our modeling fully agrees with two earlier distance estimations to the Sculptor Group, 3.91 ± 0.41 Mpc (Karachentsev et al. 2003), and to the host galaxy NGC 7793 from Cepheid method, –3.4 ± 0.17 Mpc (Pietrzyński et al. 2010). It shows that our model is correct and the proper emission from several regions with mutual occultation should be taken into account during broadband data analysis of accreting systems.

Fig. 3 Unfolded energy spectrum from all detectors used in our spectral fitting analysis. The E * FE quantity is plotted to show the maximum emission from the hard energy tail, which is associated with the emission. Black and red crosses correspond to the XMM-Newton detectors EPIC-pn and EPIC-MOS. Green and blue crosses are data from NuSTAR FPMA and FPMB, respectively. Black solid lines are the best fit models.

5 Discussion

In this paper we showed that the broadband spectra of ULX sources, P13 in NGC 7793, ULX1 in NGC 5907, and ULX5 in Circinus galaxy, are well fitted by a single model component. The fit statistics are excellent: χ² ∕d.o.f. = 1.08 for P13, 1.01 for ULX1, and 1.14 for ULX5. The model is the result of integration over a non-spherical and non-uniform emitting region. The best fit model indicates that all ULXs have a hot neutron star in their centers, with the disk with high inner temperatures. For each source the inclination of the whole system was determined. Our model does not exclude the presence of hot regions, such as hot corona, in the center since we can observe a weak deviation of the model from the data at 20 keV. Nevertheless, we aimed to show the effect of non-spherical emission in the simplest way, and we claim that it fully works for the ULXs analyzed in this paper.

Our results clearly show that the integration over true intensity emitted by a given surface should be applied to explain emission from non-spherical systems. When this is done, we obtain a new possibility to explain broadband spectra via a single model component. Next, we can derive the distance to the source from the model normalization. We derived the distance to P13, $D={3.41}_{-0.10}^{+0.11}$ Mpc, which is in very good agreement with the Cepheids distance 3.4 ± 17 Mpc (Pietrzyński et al. 2010). In the case of ULX1, the distance resulting from our model is $D={6.55}_{-0.81}^{+0.69}$ Mpc, a factor of two lower than previous estimates (Tully et al. 2009, 2016), but both previous estimations reported two different numbers (13 and 17 Mpc) for the distance to this source. The distance to ULX5 obtained by our method $D={2.60}_{-0.03}^{+0.05}$ Mpc is in agreement with the value of 2.78 Mpc found by Koribalski et al. (2004) using radial velocity measurements.

There is a more general result of our analysis. Any additional soft X-ray bump, which is very often observed in X-ray spectra of accreting objects, is usually fitted by separate model components, either disk or blackbody. Here, we have proven that this problem can be solved when we integrate emission over the disk with a hot inner source. Such a central hot source may be for instance a neutron star or a hot corona, and may be partially attenuated by the disk.

Acknowledgments

This research was supported by Polish National Science Center grants No. 2015/17/B/ST9/03422, 2015/18/M/ST9/00541, and 2016/21/N/ST9/03311. This research has made use of data obtained with NuSTAR, a project led by Caltech, funded by NASA, and managed by NASA/JPL, and has utilized the NUSTARDAS software package, jointly developed by ASDC (Italy) and Caltech (USA). This work has also made use of data obtained with XMM-Newton directly funded by ESA Member States.

References

All Tables

All Figures

	Fig. 1 Normalized counts from all detectors used in our spectral fitting analysis for P13, ULX1, and ULX5. Black and red crosses correspond to the XMM-Newton detectors EPIC-pn and EPIC-MOS. Green and blue crosses are data from NuSTAR FPMA and FPMB, respectively. Black solid lines are the best fit models.
In the text

	Fig. 2 Upperpanels: unfolded photon spectra from all detectors used in our spectral fitting analysis. Lower panels: ratio i.e. data divided by model. All colors have the same meaning as in Fig. 1.
In the text

Fig. 3 Unfolded energy spectrum from all detectors used in our spectral fitting analysis. The E * FE quantity is plotted to show the maximum emission from the hard energy tail, which is associated with the emission. Black and red crosses correspond to the XMM-Newton detectors EPIC-pn and EPIC-MOS. Green and blue crosses are data from NuSTAR FPMA and FPMB, respectively. Black solid lines are the best fit models.

In the text