Spectral Evolution of the Ultraluminous X-Ray Sources M82 X-1 and X-2

Ultraluminous X-ray sources (ULXs), first discovered by the Einstein X-ray observatory (Giacconi et al. 1979; Fabbiano 1989), are observed as bright sources of X-rays that appear to exceed the Eddington limit of the typical 10 M_⊙ black holes found in our own Galaxy. ULXs are not coincidental with the nuclei of their host galaxies, so they are not be powered by the supermassive black holes that power active galactic nuclei (AGNs), although ∼25% of candidate ULXs are likely to be background AGNs (e.g., Walton et al. 2011; Sutton et al. 2015; Earnshaw et al. 2019). At first, ULXs were promising intermediate-mass black hole (IMBH) candidates, as larger black holes can radiate at higher luminosities due to the Eddington limit scaling with mass (e.g., Colbert & Mushotzky 1999; Miller et al. 2003). But as more observations were made, the data did not appear consistent with this scenario, and instead, lower-mass black holes accreting at super-Eddington luminosities were favored for most sources (e.g., Mizuno et al. 1999), with a few IMBH candidates still remaining (e.g., ESO 243-49 HLX-1, Farrell et al. 2009).

M82 (the "cigar galaxy," or NGC 3034; Watson et al. 1984) hosts two well-known ULXs that, as they are separated by only 5'' on the sky, were only first resolved by Chandra (Matsumoto et al. 2001). The history of studies of M82 X-1 (CXOU J095550.2+694047), the brightest ULX in M82, embodies the above narrative. It has long been considered one of the best IMBH candidates because of its high luminosity, which can reach ∼10⁴¹ erg s⁻¹ (e.g., Ptak & Griffiths 1999; Rephaeli & Gruber 2002; Kaaret et al. 2006); detection of low-frequency quasi-periodic oscillations (QPOs) in the power spectrum (54 mHz, Strohmayer & Mushotzky 2003; Dewangan et al. 2006; Mucciarelli et al. 2006), indicative of a compact, unbeamed source; as well as twin-peaked QPOs at 3.3 and 5.1 Hz, which lead to a mass estimate of 400 M_⊙ using scaling laws between the QPOs frequencies and mass used for stellar-mass black holes (Pasham et al. 2014). Additionally, Feng & Kaaret (2010) (hereafter F10) observed the source with XMM-Newton and Chandra over the course of a flaring episode, and fitted the spectra with the standard thin accretion disk model. They found that the luminosity of the disk, L, scaled with inner temperature as L ∝ T⁴, which is expected from a thin accretion disk with a constant inner radius. From this, they inferred a black hole mass in the range 300–810 M_⊙, assuming that the black hole is rapidly spinning to avoid extreme violations of the Eddington limit, therefore adding support to the IMBH scenario.

The luminosity argument no longer stands, however, as the ULX NGC 5907 ULX1, which also reaches ∼10⁴¹ erg s⁻¹ (Sutton et al. 2013; Fürst et al. 2017), and was once considered an IMBH candidate, was discovered to be powered by a neutron star with only 1–2 M_⊙ (Israel et al. 2017) from the detection of coherent pulsations. The mass measurement from the twin QPOs is ambiguous too, as it is not known where these originate, and it is not clear if the scaling relationship used extends to the IMBH range. In Brightman et al. (2016b, hereafter B16), we presented a combined spectral analysis of X-1, also during a flaring episode, using simultaneous observations with Chandra, NuSTAR and Swift. With broader-band data than Feng & Kaaret (2010), we found that the temperature profile as a function of disk radius (T(r) ∝ r^−p) is significantly flatter than expected for a standard thin accretion disk as implied by F10, and instead characteristic of a slim disk that is expected at high accretion rates. As only one observation was analyzed, the L ∝ T⁴ relationship could not be tested. Nevertheless, the mass estimates inferred from the inner disk radius for this model imply a stellar-remnant black hole (${M}_{\mathrm{BH}}={26}_{-6}^{+9}\,{M}_{\odot }$), when assuming zero spin, or an IMBH (${M}_{\mathrm{BH}}={125}_{-30}^{+45}\,{M}_{\odot }$) when assuming maximal spin.

M82 X-2 (CXOM82 J095551.1+694045) is typically the second brightest X-ray source in the galaxy, and was the first ultraluminous X-ray pulsar (ULXP) discovered (Bachetti et al. 2014) using NuSTAR (Harrison et al. 2013). With only a 5'' angular separation from X-1, studying the spectral properties of this source in detail has been limited to Chandra observations. We conducted a study of the spectral and temporal properties of M82 X-2 in Brightman et al. (2016a), finding that the source's luminosity varies over two orders of magnitude over the range ${10}^{38}\mbox{--}{10}^{40}$ erg s⁻¹. Using timing analyses, we were able to isolate the pulsed emission from this source with NuSTAR, finding that it was well described by a cutoff power law. In a follow-up study, we found evidence that the variations in luminosity are modulated on a ∼60 days period, and that as the orbit of the neutron star and its companion is known to be 2.5 days, the ∼60 days period must be super-orbital in origin (Brightman et al. 2019).

In this paper, we present analysis of a systematic monitoring campaign on M82 by Chandra that took place in 2016. The primary goals of this campaign were to perform a temporal analysis of X-1 and X-2 to search for orbital and super-orbital modulations; to perform spectroscopic studies of the ULXs; and to study the nature of the other binary systems in M82. We report the temporal analysis in Brightman et al. (2019), and here we focus on the second objective: to present the most comprehensive X-ray spectral analysis of the two ULXs, M82 X-1 and X-2, to date. We combine simultaneous observations with Chandra to spatially resolve the two sources from each other and the other sources in M82, and simultaneous observations with NuSTAR to extend the spectral coverage up to 79 keV. We also present results from Swift monitoring observations of M82 that have been ongoing since 2012. Throughout this paper, we assume a distance to M82 of 3.3 Mpc (Foley et al. 2014), which is inferred from the light curve of SN2014J.

All Chandra and NuSTAR observations studied here were taken simultaneously, or quasi-simultaneously (within 24 hr of each other). Table 1 provides a description of the observational data. The following sections describe the individual observations and data reduction.

2.1. Chandra

As the angular separation of X-1 and X-2 is only 5'', only Chandra (Weisskopf 1999) can spatially resolve the emission from these two sources. The majority of the Chandra data analyzed here were taken during 2016 (Cycle 17) as part of a large program aimed at systematic monitoring of binaries in M82. The program consisted of 12 individual observations taken at ∼monthly intervals, but only eight having simultaneous NuSTAR observations were used in this work. We additionally use two observations taken in 2015 that also have a simultaneous NuSTAR observation. Full details are listed in Table 1. All 2016 observations were taken with ACIS-I at the optical axis using a 1/8 subarray of pixels on chip I1 or I3, depending on the roll angle. The ULXs at the center of M82 were placed 35 off-axis to smear out the point-spread function (PSF) enough to reduce the effects of pileup, but not so much as to cause significant blending of the PSFs. The subarray of pixels was used to decrease the readout time of the detector to 0.4 s, further reducing the effects of pileup.

We proceeded to extract the Chandra spectra in the same way as described in B16, using the ciao (v4.7, CALDB v4.6.5) tool specextract. For the point sources, spectra were extracted from elliptical regions drawn by eye to encompass the shape of the off-axis Chandra PSF. For X-1, we used an ellipse with a semimajor axis of 2''–3'' and a semiminor axis of 1''–2''. For X-2, the major and minor axes were 2'' and 1'', respectively. A small rectangular region close by was used for background subtraction. Figure 1 shows examples of the extraction regions used.

Zoom In Zoom Out Reset image size

We used the ciao tool pileup_map to give an indication of the level of pileup in each observation. The output, which is in counts per frame, ranges from 0.06 to 0.35 at maximum, which occurs at the position of X-1. These numbers correspond to pileup fractions of  <5% for  <0.1 counts per frame to  >10% for  >0.2 counts per frame. For three observations, obsIDs 17678, 18065, and 18067, X-1 has a pileup fraction of  >10%. We will discuss the potential effects of this on our results later in the paper.

For X-ray emission from M82 that does not come from X-1 or X-2, but that contributes to the emission seen by NuSTAR, such as from the fainter point sources and the diffuse emission, we extract spectra from a 49'' radius circular region centered on X-1, but with the X-1 and X-2 regions as described above masked out. A larger background region external to the galaxy was extracted to assess the background coming from the cosmic X-ray and particle backgrounds.

2.2. NuSTAR

2.3. Swift

Swift conducted monitoring of M82 with a typical cadence of a few days between 2012 and 2018. As this monitoring ran contemporaneously with our Chandra and NuSTAR observations, the well-sampled light curve provides context for our study. A total of 113 observations, mostly consisting of obsIDs 00032503099–154 and 00092202001–051, have been made of the galaxy over the 2015–2016 period which we use here to calculate a long-term light curve.

We calculate the fluxes via spectral fitting. We use the heasoft (v 6.16) tool xselect to filter events from a 49'' radius region centered on the ULXs and to extract the spectrum. This extraction region encloses all sources of X-ray emission in the galaxy. Background events were extracted from a nearby circular region of the same size. We group the spectra with a minimum of one count per bin using the heasoft tool grppha. We conduct spectral fitting in the range 0.2–10 keV. We fit the spectra with a simple power law subjected to absorption intrinsic to M82 at z = 0.00067 (zwabs∗powerlaw in xspec) with the Cash statistic (Cash 1979), which uses a Poisson likelihood function and is hence most suitable for low numbers of counts per bin. From this model we calculate the observed flux in the 0.5–8 keV range, equivalent to the Chandra band. The light curve is presented in Figure 2.

Zoom In Zoom Out Reset image size

The Swift light curve shows that our observations took place during a period of flaring activity from M82, which we found to be due to increased activity from X-1 (B16). The first Chandra/NuSTAR observations took place before the increase in activity, and the remaining observations tracked the activity over the next two years.

All Chandra and NuSTAR spectra were grouped with a minimum of 20 counts per bin using the heasoft tool grppha. Spectral fitting was carried out using xspec v12.8.2 (Arnaud 1996) and the ${\chi }^{2}$ statistic was used for spectral fitting to background subtracted spectra. All uncertainties quoted are 90%. We present the Chandra spectra of X-1, X-2, and the additional X-ray emission from M82, and the NuSTAR spectra of the entire galaxy for each of the 10 observational epochs in Figure 3.

Zoom In Zoom Out Reset image size

While the presence of X-2 in a bright state introduced some ambiguity to the results on the spectral analysis of X-1 in B16, in four observations here, Chandra obsIDs 18062, 18065, 18069, and 18070, X-2 was observed in a low state, allowing an unambiguous view of X-1 for the first time using NuSTAR. For these four observations, we neglect the emission from X-2 in our spectral fits, whereas for the rest of the observations we included it.

We proceed to conduct a spectral analysis for all 10 observations, where we fitted the spectra in the same way as described in B16. A cross-calibration constant was applied to each spectral data set to allow for absolute differences in normalization, and allowed to vary by ±10% (Madsen et al. 2015). For the diffuse emission from M82, we use a combination of three absorbed zwabs∗apec models with the temperatures and abundances fixed to the values found in B16. We allow the normalizations, both with respect to B16 and to each other, to vary here to account for small differences in the detector responses. For the spectrum of X-1, we use the zwabs∗diskpbb model, where zwabs is a redshifted neutral absorption component and diskpbb is a model representing emission from a multicolour accretion disk, with a variable radial temperature profile. This model combination was found to best represent the emission from X-1 with regards to other disk models in B16. Additionally, as X-1 is a bright source, and despite measures taken to reduce pileup, the source still suffers from pileup. We account for this in spectral fitting using the pileup convolution model, with frame time set to 0.4 s (Davis 2001).

In the top panel of Figure 4, we show the residuals to this set of models which reveal strong residuals above 10 keV, especially in obsID 18062, which is reminiscent of the hard tail seen in the NuSTAR spectra of several other ULXs such as Holmberg II X-1 (e.g., Walton et al. 2015), including those already identified as neutron star accretors such as NGC 7793 P13 (Walton et al. 2018a) and NGC 5907 ULX (Fürst et al. 2017; Walton et al. 2018b). This was also seen in our analysis in B16, but as X-2 was bright during that observation, there was ambiguity regarding its origin. Here it is clear that it originates from X-1.

Zoom In Zoom Out Reset image size

We test two models to fit this component, simpl, which describes the power-law emission from the Comptonization of a disk spectrum (Steiner et al. 2009), and a cutoffpl model used to model the pulsed emission from ULXPs (Brightman et al. 2016a; Walton et al. 2018a, 2018b). In Walton et al. (2018b), where no pulsations from a source had been detected, as is the case for M82 X-1, Γ was fixed at 0.5 for the cutoffpl model which is the average for the pulsed emission from ULXPs. Therefore we do the same when using the cutoffpl model and fix Γ = 0.5. The energy of the cutoff was allowed to vary. We find that the simpl model provides a better description of the data, with ${\rm{\Delta }}{\chi }^{2}$ = 70–200 for the addition of two free parameters from simpl and ${\rm{\Delta }}{\chi }^{2}$ = 3–90 for the addition of two free parameters from cutoffpl as shown in the middle and bottom panels of Figure 4. We do not find a significant improvement for the cutoffpl model if we allow Γ to vary. Therefore we use the simpl model to fit the high-energy emission from X-1 for the full data set.

We proceed to fit the Chandra and NuSTAR spectra of all 10 observational epochs with the pileup∗zwabs∗diskpbb∗simpl model combination to describe the emission from X-1. For the six observations where X-2 is bright, we model the emission from this source with an absorbed cutoff power-law model, zwabs∗cutoffpl, which was used in B16 to model the pulsed emission. Figure 5 presents the unfolded spectra with the different model components shown. Table 2 presents the best-fit spectral parameters for X-1, and Table 3 presents the best-fit spectral parameters for X-2.

Zoom In Zoom Out Reset image size

Table 2.  Spectral Fitting Results for X-1

Chandra	${N}_{{\rm{H}}}$	Tin	p	log10 norm	Γ	${f}_{\mathrm{scatt}}$	${\chi }^{2}$/DoF	${F}_{{\rm{X}}}$	${L}_{{\rm{X}}}$
ObsID	(1022 cm−2)	(keV)						(10−11 erg cm−2 s−1)	(1040 erg s−1)
17578	${0.98}_{-0.33}^{+0.15}$	${2.8}_{-1.9}^{+2.7}$	${0.59}_{-0.03}^{+0.13}$	$-{2.6}_{-1.2}^{+2.2}$	${3.0}_{-2.0}^{+2.0}$	${0.98}_{-0.93}^{+0.02}$	941.5/ 840	${0.97}_{-0.04}^{+0.04}$	${1.3}_{-0.1}^{+0.1}$
17678	${1.35}_{-0.08}^{+0.17}$	${1.7}_{-0.4}^{+0.4}$	${0.64}_{-0.03}^{+0.08}$	$-{0.8}_{-0.4}^{+0.4}$	${3.7}_{-0.5}^{+0.7}$	${0.73}_{-0.40}^{+0.27}$	1984.5/1259	${4.26}_{-0.13}^{+0.19}$	${5.6}_{-0.2}^{+0.3}$
18062	${1.81}_{-0.13}^{+0.13}$	${2.3}_{-0.7}^{+0.8}$	${0.52}_{-0.02}^{+0.03}$	$-{1.8}_{-0.5}^{+0.6}$	${3.3}_{-0.8}^{+0.4}$	${0.60}_{-0.47}^{+0.40}$	2070.3/1246	${3.73}_{-0.15}^{+0.08}$	${4.9}_{-0.2}^{+0.1}$
18063	${1.25}_{-0.11}^{+0.07}$	${6.7}_{-3.1}^{+0.6}$	${0.57}_{-0.01}^{+0.01}$	$-{3.7}_{-0.2}^{+0.7}$	${4.5}_{-4.5}^{+0.5}$	${0.00}_{-0.00}^{+0.26}$	1445.8/1282	${1.82}_{-0.05}^{+0.05}$	${2.4}_{-0.1}^{+0.1}$
18064	${1.09}_{-0.15}^{+0.16}$	${4.5}_{-2.3}^{+1.7}$	${0.58}_{-0.03}^{+0.04}$	$-{3.3}_{-0.6}^{+0.9}$	${2.7}_{-2.7}^{+2.4}$	${0.23}_{-0.23}^{+0.77}$	1486.0/1066	${1.04}_{-0.03}^{+0.03}$	${1.4}_{-0.0}^{+0.0}$
18065	${1.25}_{-0.08}^{+0.08}$	${2.3}_{-0.3}^{+0.7}$	${0.57}_{-0.02}^{+0.02}$	$-{1.6}_{-0.4}^{+0.2}$	${3.7}_{-0.6}^{+1.3}$	${0.92}_{-0.60}^{+0.08}$	2069.2/1404	${4.21}_{-0.10}^{+0.16}$	${5.5}_{-0.1}^{+0.2}$
18067	${1.38}_{-0.10}^{+0.11}$	${2.0}_{-0.6}^{+0.5}$	${0.56}_{-0.02}^{+0.03}$	$-{1.4}_{-0.2}^{+0.5}$	${3.5}_{-0.7}^{+1.1}$	${0.58}_{-0.36}^{+0.42}$	2410.8/1567	${3.10}_{-0.11}^{+0.12}$	${4.0}_{-0.1}^{+0.2}$
18068	${0.94}_{-0.14}^{+0.12}$	${3.8}_{-1.7}^{+0.9}$	${0.64}_{-0.02}^{+0.06}$	$-{2.7}_{-0.4}^{+1.0}$	${2.4}_{-1.4}^{+1.1}$	${0.31}_{-0.24}^{+0.58}$	1502.7/1366	${0.94}_{-0.05}^{+0.07}$	${1.2}_{-0.1}^{+0.1}$
18069	${1.25}_{-0.09}^{+0.10}$	${3.6}_{-1.0}^{+1.3}$	${0.56}_{-0.02}^{+0.02}$	$-{2.7}_{-0.3}^{+0.5}$	${2.7}_{-1.4}^{+1.1}$	${0.28}_{-0.26}^{+0.33}$	1529.0/1275	${2.31}_{-0.06}^{+0.09}$	${3.0}_{-0.1}^{+0.1}$
18070	${1.30}_{-0.11}^{+0.10}$	${3.1}_{-0.5}^{+0.6}$	${0.56}_{-0.02}^{+0.02}$	$-{2.5}_{-0.1}^{+0.3}$	${4.0}_{-1.4}^{+1.0}$	${0.91}_{-0.86}^{+0.09}$	1598.2/1252	${2.02}_{-0.05}^{+0.05}$	${2.6}_{-0.1}^{+0.1}$

Note. Best-fit parameters for the diskpbb model fitted to X-1. Fluxes and luminosities are given in the 0.5−30 keV range, and are corrected for absorption and pileup. Luminosities are calculated assuming a distance of 3.3 Mpc to M82.

Download table as:  ASCIITypeset image

Table 3.  Spectral Fitting Results for X-2

Chandra	${N}_{{\rm{H}}}$	Γ	${E}_{{\rm{C}}}$	log10 norm	${\chi }^{2}$/DoF	${F}_{{\rm{X}}}$	${L}_{{\rm{X}}}$
ObsID	(1022 cm−2)		(keV)			(10−11 erg cm−2 s−1)	(1040 erg s−1)
17578	${1.7}_{-1.4}^{+2.4}$	$-{1.3}_{-3.6}^{+3.4}$	${1.5}_{-0.4}^{+22}$	$-{3.6}_{-0.5}^{+0.9}$	899.9/ 840	${0.4}_{-0.1}^{+0.3}$	${0.5}_{-0.1}^{+0.4}$
17678	${3.8}_{-1.1}^{+1.2}$	${1.7}_{-0.7}^{+0.7}$	${15}_{-16}^{+77}$	$-{2.7}_{-0.4}^{+0.6}$	1378.4/1259	${1.4}_{-0.3}^{+0.5}$	${1.8}_{-0.4}^{+0.6}$
18063	${1.9}_{-0.4}^{+0.4}$	$-{1.2}_{-0.4}^{+0.4}$	${2.0}_{-0.3}^{+0.4}$	$-{3.6}_{-0.1}^{+0.2}$	1370.8/1282	${1.0}_{-0.1}^{+0.1}$	${1.2}_{-0.1}^{+0.1}$
18064	${2.3}_{-1.0}^{+1.3}$	$-{1.3}_{-1.6}^{+1.3}$	${1.4}_{-0.4}^{+0.9}$	$-{3.4}_{-0.4}^{+0.4}$	1233.9/1066	${0.5}_{-0.1}^{+0.2}$	${0.6}_{-0.1}^{+0.3}$
18067	${3.0}_{-0.6}^{+0.7}$	${1.3}_{-0.6}^{+0.5}$	${16}_{-10}^{+25}$	$-{3.0}_{-0.2}^{+0.3}$	1725.4/1567	${1.0}_{-0.2}^{+0.1}$	${1.3}_{-0.2}^{+0.2}$
18068	${0.8}_{-0.8}^{+0.2}$	$-{2.9}_{-2.0}^{+6.9}$	${1.2}_{-0.1}^{+0.0}$	$-{4.0}_{-0.0}^{+0.0}$	1478.3/1366	${0.7}_{-0.0}^{+0.0}$	${1.0}_{-0.0}^{+0.1}$

Note. Best-fit parameters for the cutoffpl model fitted to X-2. Fluxes and luminosities are given in the 0.5–30 keV range, and are corrected for absorption. Luminosities are calculated assuming a distance of 3.3 Mpc to M82.

Download table as:  ASCIITypeset image

To calculate the intrinsic luminosity of each source, we use the cflux model component in xspec placed after the absorption component with the normalization of the main model fixed. We calculate the intrinsic flux over the range 0.5–30 keV, and present these values with their uncertainties in Tables 2 and 3. We calculate the intrinsic luminosities over these ranges assuming a distance to M82 of 3.3 Mpc (Foley et al. 2014).

Our main goal is to explore the spectral evolution of M82 X-1. For each pair of spectral parameters, we compute Spearman's rank correlation, using the idl tool r_correlate.pro, to assess the presence of any correlation and its significance. This assesses how well the relationship between two variables can be described using a monotonic function. The results from this test, being the rank correlation coefficient, ρ, and the two-sided significance of its deviation from zero, p, are presented in Table 4.

Table 4.  Spearman's Rank Correlation Results for X-1

	${L}_{{\rm{X}}}$	${N}_{{\rm{H}}}$	Tin	p	log10norm	Γ
${N}_{{\rm{H}}}$	0.76
	0.0111
Tin	−0.76	−0.56
	0.011	0.090
p	−0.33	−0.62	0.03
	0.347	0.054	0.934
log10norm	0.79	0.62	−0.99	−0.10
	0.006	0.054	0.000	0.777
Γ	0.48	0.58	−0.22	−0.36	0.28
	0.162	0.082	0.533	0.310	0.425
${f}_{\mathrm{scatt}}$	0.28	0.09	−0.67	0.05	0.65	0.21
	0.425	0.803	0.033	0.881	0.043	0.043

Note. For each pair of parameters we list the rank correlation coefficient (top) and the two-sided significance of its deviation from zero (bottom).

Download table as:  ASCIITypeset image

The significance is a value in the interval [0.0, 1.0] and a small value indicates a significant correlation. The commonly used threshold to judge that a correlation is significant is p < 0.05. We find that this criterion is met for the following pairs of parameters: ${L}_{{\rm{X}}}$ and ${N}_{{\rm{H}}}$; ${L}_{{\rm{X}}}$ and T_in; ${L}_{{\rm{X}}}$ and log₁₀ norm; T_in and log₁₀ norm; ${L}_{{\rm{X}}}$ and ${f}_{\mathrm{scatt}}$; T_in and ${f}_{\mathrm{scatt}}$; and Γ and ${f}_{\mathrm{scatt}}$.

However, for many parameters, especially those of the simpl model, the uncertainties are large, which the correlation test does not account for. First, to account for this, we conduct 1000 Monte Carlo simulations where we randomly draw a value from the 90% confidence interval for each pair of parameters. We calculate from how many of the 1000 simulations we find a correlation with a p-value less than 0.05. For all correlations involving parameters of the simpl model, we find that in less than 25% of the simulations, a p-value less than 0.05 is recovered. Therefore, we do not have confidence in these correlations being real, due to the large uncertainties in the parameters.

Furthermore, for the correlations where the uncertainties in the parameters are smaller, some of these correlations can be expected due to degeneracies between parameters. To determine if these correlations are driven by degeneracies, we explore the two-dimensional ${\chi }^{2}$ space around the best fit using the xspec command steppar, and overplot the 3σ contours on Figure 6. We do this only for one observation, Chandra obsID 18069, as it is computationally expensive, but this should be adequate to reveal any spectral degeneracies. This data set is a typical observation where X-2 is at low fluxes, pileup for X-1 is low, and the measured parameters are in the middle of the distributions.

Zoom In Zoom Out Reset image size

The contours show that there is slight degeneracy between ${L}_{{\rm{X}}}$ and ${N}_{{\rm{H}}}$, but that does not appear to be strong enough to induce the correlation. There is no apparent degeneracy between ${L}_{{\rm{X}}}$ and T_in. A correlation between ${L}_{{\rm{X}}}$ and log₁₀ norm is expected, and a clear degeneracy between T_in and log₁₀ norm is seen.

As the background is higher in the NuSTAR data than in the Chandra data, and the signal to noise lower, we investigated whether binning the NuSTAR spectra with more counts would affect our results. We grouped the NuSTAR spectra with a minimum of 60 counts per bin rather than 20, and refit the spectra. We looked at the temperature of the diskpbb component, which is most sensitive to the NuSTAR data and the one of the key parameters involved in our results. We found the average difference in temperature between the stronger binning and the weaker one to be −0.02 keV, which is much smaller than the typical uncertainty on the parameter. We therefore conclude that the spectral binning does not affect our results.

Taking into account the uncertainties in parameters and degeneracies between them, we can only say with confidence that there is a correlation and therefore physical link between the neutral column density and the X-ray luminosity, and the the inner disk temperature and X-ray luminosity.

We test the apparent dependence of ${N}_{{\rm{H}}}$ on ${L}_{{\rm{X}}}$ further by fixing ${N}_{{\rm{H}}}$ at the approximate mean value of 1.3 × 10²² cm⁻² in all fits and note the difference in ${\chi }^{2}$. In all cases, ${\chi }^{2}$ is worse or the same as when ${N}_{{\rm{H}}}$ is a free parameter. For Chandra obsID 18062, ${L}_{{\rm{X}}}$ is high and the ${N}_{{\rm{H}}}$ measured is particularly high at 1.8 × 10²² cm⁻². When ${N}_{{\rm{H}}}$ is fixed at the mean value, the difference in ${\chi }^{2}$ is 80. For obsID 18068, ${L}_{{\rm{X}}}$ is low and ${N}_{{\rm{H}}}$ = 9 × 10²¹ cm⁻². When fixed at 1.3 × 10²² cm⁻², the ${\chi }^{2}$ increases by 150, supporting the finding that ${N}_{{\rm{H}}}$ does indeed depend on ${L}_{{\rm{X}}}$.

The anticorrelation between ${L}_{{\rm{X}}}$ and T_in that we find here is in contrast with the correlation found between these two parameters for M82 X-1 by Feng & Kaaret (2010), where the apparent ${L}_{{\rm{X}}}$ ∝ T⁴ relationship lead the authors to conclude that the source was observed in the thermal dominant state. To determine the exponent on the relationship, we take the logarithm of both ${L}_{{\rm{X}}}$ and T_in, such that log₁₀${L}_{{\rm{X}}}$ = αlog₁₀T_in+β and conduct a linear regression analysis. We use the idl tool linmix_err.pro, which takes into account uncertainties in both parameters. We find that log₁₀${L}_{{\rm{X}}}$ = (−1.47 ± 1.00)log₁₀T_in + (41.16 ± 0.47), where the uncertainties are 1σ. Therefore, the exponent is −1.47 ± 1.00, which is a >5σ difference from the results from F10, where the exponent was implied to be 4, although the uncertainty in this parameter was not listed. If we run the same linear regression analysis on the results presented in F10, we find that the relationship is essentially unconstrained.

The key differences between our analysis and that of F10 are that they did not have accompanying NuSTAR data, so their analysis was restricted to the narrow bandpass of 0.7–7 keV. The spectral models used are also different, whereby F10 used the diskbb model, which is related to diskpbb that we use when p is fixed at 0.75. However, as we showed in B16, and confirm here, the data are inconsistent with the p = 0.75 that describes a geometrically thin accretion disk. Furthermore, without NuSTAR data, the high-energy excess that we detect and fit with simpl was not observable by F10. Finally, the uncertainties related to degeneracies in the pileup model are reduced when data from an instrument without pileup such as NuSTAR is used.

We investigate what effect pileup has on our results by exploring the dependence of ${L}_{{\rm{X}}}$, ${N}_{{\rm{H}}}$ and T_in on the pileup model parameter, α. We find that none of these parameters show any dependence on α, leading us to conclude that this model component does not drive our results. Furthermore, as mentioned earlier, for three observations, odsID 17678, 18065, and 18067, the pileup fraction for X-1 is greater than 10%, and therefore the pileup model used in the spectral fitting may not be able to reliably account for this effect. We test the dependence of our results on these observations by removing them from our analysis. In doing so, we still find a significant correlation between ${L}_{{\rm{X}}}$ and ${N}_{{\rm{H}}}$, and indeed the significance increases, but for ${L}_{{\rm{X}}}$ and T_in the correlation is no longer significant. A fit with a linear relationship reveals log₁₀${L}_{{\rm{X}}}$ = (−1.23 ± 2.29)log₁₀T_in+(41.03 ± 1.28) and therefore we can only rule out a ${L}_{{\rm{X}}}$ ∝ T⁴ dependency at ∼2σ.

Finally, we note that some of our spectral fits have large ${\chi }^{2}$ relative to the number of degrees of freedom (DoF), indicating an unacceptable fit. This is likely due to the fits being quite complex, with many DoFs and data sets being fitted simultaneously. However, we note the presence of a potential excess of counts in the NuSTAR data at 3–4 keV that has been found in bright X-ray binaries and is a known calibration issue. We find that if we ignore data from NuSTAR below 5 keV, the fits improve and the result of this is to systematically reduce the temperature of the diskpbb component. This reduction is consistent with the uncertainties on this component, however, and does not alter our result that there is an anticorrelation between ${L}_{{\rm{X}}}$ and T_in because the effect is systematic across all observations.

These findings therefore rule out a thermal state for sub-Eddington accretion and therefore do not support M82 X-1 as an IMBH candidate.

We conduct the same analysis for X-1 on X-2, where the spectral parameters are plotted against each other in Figure 7 and the correlation analysis is shown in Table 5. Here our analysis suggests correlations between ${N}_{{\rm{H}}}$ and Γ, E_C and Γ, ${N}_{{\rm{H}}}$ and log₁₀norm, and Γ and log₁₀norm. However, the ${\chi }^{2}$ contours show that it is likely that strong degeneracies between these parameters drive these apparent correlations (Figure 7). We therefore do not find evidence for any significant spectral evolution in X-2.

Zoom In Zoom Out Reset image size

Table 5.  Spearman's Rank Correlation Results for X-2

	${L}_{{\rm{X}}}$	${N}_{{\rm{H}}}$	EC	Γ
${N}_{{\rm{H}}}$	0.71
	0.111
EC	0.71	0.77
	0.111	0.072
Γ	0.77	0.83	0.94
	0.072	0.042	0.005
log10norm	0.71	1.00	0.77	0.83
	0.111	0.000	0.072	0.042

Note. For each pair of parameters, we list the rank correlation coefficient (top) and the two-sided significance of its deviation from zero (bottom).

Download table as:  ASCIITypeset image

While the emission from X-2 is difficult to disentangle from the other sources of emission in M82, we were able to isolate the pulsed emission in the NuSTAR band from this source in Brightman et al. (2016a). We found that the pulsed emission is best fit by a power law with a high-energy cutoff, where Γ = 0.6 ± 0.3 and ${E}_{{\rm{C}}}={14}_{-3}^{+5}\,\mathrm{keV}$. In Figure 7, we show the parameters of the pulsed emission as a separate data point for comparison. We see that the values for Γ and E_C from our broadband fits are consistent with the pulsed emission when X-2 is at its highest luminosities, ${L}_{{\rm{X}}}$ > 10⁴⁰ erg s⁻¹, indicating that at these times the pulsations are most likely to be detected. Bachetti et al. (2019) have recently detected pulsations again from a NuSTAR observation taken on 2016 September 10, which unfortunately did not have any simultaneous Chandra observations, so we did not include it in our analysis here.

6.1. M82 X-1 as an Intermediate-mass Black Hole Candidate

6.2. M82 X-1 as a Super-Eddington Accreting Stellar Remnant

We conducted a comprehensive investigation into the spectral evolution of the ultraluminous X-ray sources M82 X-1 and X-2 using 10 simultaneous Chandra and NuSTAR observations. The Chandra data allowed us to spatially resolve the sources, separated by only 5'', while the NuSTAR data allowed a broadband X-ray spectral analysis. We found that for X-1, the luminosity of the disk scales with the inner temperature as L ∝ T^−3/2, which is contrary to previous findings of a L ∝ T⁴ scaling that supported a standard accretion disk powering the system. We furthermore find evidence that the neutral column density of the material in the line of sight increases with ${L}_{{\rm{X}}}$, perhaps due to an increased mass outflow with accretion rate. For X-2, we do not find any significant spectral evolution, but we can constrain the spectral parameters, and find the broadband emission to be consistent with the pulsed emission at the highest X-ray luminosities.

This research has made use of data obtained with NuSTAR, a project led by the California Institute of Technology, managed by the Jet Propulsion Laboratory, and funded by NASA. We thank the NuSTAR Operations, Software and Calibration teams for support with the execution and analysis of these observations. This research has made use of the NuSTAR Data Analysis Software (NuSTARDAS) jointly developed by the ASI Science Data Center (ASDC, Italy) and the California Institute of Technology (USA). Support for this work was provided by the National Aeronautics and Space Administration through Chandra Award Number GO6-17080X issued by the Chandra X-ray Center, which is operated by the Smithsonian Astrophysical Observatory for and on behalf of the National Aeronautics Space Administration under contract NAS8-03060. This research has made use of software provided by the Chandra X-ray Center (CXC) in the application package ciao. We also acknowledge the use of public data from the Swift data archive. D.J.W. acknowledges support from an STFC Ernest Rutherford Fellowship.

Facilities: Chandra (ACIS) - , NuSTAR - , Swift (XRT). -