Articles

NuSTAR REVEALS AN INTRINSICALLY X-RAY WEAK BROAD ABSORPTION LINE QUASAR IN THE ULTRALUMINOUS INFRARED GALAXY MARKARIAN 231

, , , , , , , , , , , , , , , , , , , , , , , , , , and

Published 2014 March 21 © 2014. The American Astronomical Society. All rights reserved.
, , Citation Stacy H. Teng et al 2014 ApJ 785 19 DOI 10.1088/0004-637X/785/1/19

0004-637X/785/1/19

ABSTRACT

We present high-energy (3–30 keV) NuSTAR observations of the nearest quasar, the ultraluminous infrared galaxy (ULIRG) Markarian 231 (Mrk 231), supplemented with new and simultaneous low-energy (0.5–8 keV) data from Chandra. The source was detected, though at much fainter levels than previously reported, likely due to contamination in the large apertures of previous non-focusing hard X-ray telescopes. The full band (0.5–30 keV) X-ray spectrum suggests the active galactic nucleus (AGN) in Mrk 231 is absorbed by a patchy and Compton-thin ($N_{\rm H} \sim 1.2^{+0.3}_{-0.3}\times 10^{23}$ cm−2) column. The intrinsic X-ray luminosity (L0.5 − 30 keV ∼ 1.0 × 1043 erg s−1) is extremely weak relative to the bolometric luminosity where the 2–10 keV to bolometric luminosity ratio is ∼0.03% compared to the typical values of 2%–15%. Additionally, Mrk 231 has a low X-ray-to-optical power law slope (αOX ∼ −1.7). It is a local example of a low-ionization broad absorption line quasar that is intrinsically X-ray weak. The weak ionizing continuum may explain the lack of mid-infrared [O iv], [Ne v], and [Ne vi] fine-structure emission lines which are present in sources with otherwise similar AGN properties. We argue that the intrinsic X-ray weakness may be a result of the super-Eddington accretion occurring in the nucleus of this ULIRG, and may also be naturally related to the powerful wind event seen in Mrk 231, a merger remnant escaping from its dusty cocoon.

Export citation and abstract BibTeX RIS

1. INTRODUCTION

The presence of outflows in luminous quasars and active galactic nuclei (AGNs) is most evident in broad absorption line (BAL) systems (Lynds 1967). Approximately 20% of quasars show BAL features at some level (e.g., Gibson et al. 2009), understood to be the result of an observational sightline through AGN-driven wind material (e.g., Weymann et al. 1991; Ogle et al. 1999; DiPompeo et al. 2013). In a common model for BAL quasars (e.g., Murray et al. 1995; Proga et al. 2000), the wind is launched from the accretion disk approximately 0.01–0.1 pc from the black hole and radiatively driven by ultraviolet (UV) line pressure. However, X-ray emission from typical quasars would overionize the gas, quenching the line-driving mechanism and thus preventing the appearance of BAL features. Therefore, this so-called disk-wind model requires that BAL quasars have lower levels of X-ray emission, which could either be intrinsic or be due to shielding very close to the black hole, near the base of the wind. Indeed, BAL systems are very often observed to be X-ray-faint (e.g., Gallagher et al. 2006; Gibson et al. 2009; Wu et al. 2010), though observations have yet to provide a clear picture of the origin of this X-ray faintness. Absorption of a nominal X-ray continuum has been established in several cases (e.g., Gallagher et al. 2002; Grupe et al. 2003; Shemmer et al. 2005; Giustini et al. 2008), but there may be diversity among the population. Following the terminology of other studies of X-ray emission in BAL quasars (e.g., Luo et al. 2013 and references therein), the term "X-ray weakness or faintness" refers to the observed X-ray emission being significantly lower than expected given the optical–UV emission. This may be caused by a normal X-ray continuum modified by absorption from gas and dust. The term "intrinsic X-ray weakness" refers to an innate property of the AGN, which is one possible cause of the observed X-ray faintness.

In order to investigate whether intrinsic X-ray weakness or shielding gas is the most likely explanation for the observed X-ray faintness of BAL quasars, sensitive observations above 10 keV are highly desirable as they can distinguish between the presence of Compton-thick (NH ⩾ 1.5 × 1024 cm−2) shielding material and a less-obscured, intrinsically X-ray weak nucleus. The Nuclear Spectroscopic Telescope Array (NuSTAR; Harrison et al. 2013) is the first focusing high-energy X-ray telescope in orbit, and provides unprecedented sensitivity and angular resolution for high-energy, or hard (3–79 keV) X-ray photons. Recently, NuSTAR observed two of the optically brightest BAL quasars known, PG 1004+130 and PG 1700+518, with the goal of investigating the origin of their X-ray faintness (Luo et al. 2013). Unfortunately, low count rates prevented a definitive answer, but they conclusively showed that any absorption present must be Compton-thick. Luo et al. (2013) were also able to demonstrate that about 17%–40% of BAL quasars are intrinsically X-ray weak by stacking Chandra observations of z ≈ 1.5–3 BAL quasars in the Large Bright Quasar Survey.

Confirmation of intrinsically X-ray weak AGN in some BAL quasars would have important consequences regarding the origin of BAL features. Instead of invoking shielding gas with an arbitrarily thick column density, intrinsic X-ray weakness would favor a scenario in which the BAL features are coupled to abnormally faint coronal X-ray emission, perhaps due to some mechanism that (periodically) quenches the X-ray corona (see Section 4.2 of Luo et al. 2013 for discussion). To this end, we have obtained an extremely sensitive, broadband (0.5–30 keV) X-ray spectrum of Markarian 231 (Mrk 231), the nearest BAL quasar, using Chandra and NuSTAR.

This paper is organized as follows. Section 2 details the multi-wavelength properties of Mrk 231. Section 3 presents the new X-ray observations obtained for this study. Section 4 discusses the previously reported hard X-ray spectra of Mrk 231, which are inconsistent with our new data. Section 5 presents our modeling of the new broadband X-ray spectrum. Section 6 discusses our favored scenario, that Mrk 231 is indeed intrinsically X-ray weak; Section 7 summarizes our results. Throughout this paper, we adopt H0 = 71 km s−1 Mpc−1, ΩM = 0.27, and Λ = 0.73 (Hinshaw et al. 2009). Luminosities taken from the literature have been recalculated for our assumed cosmology.

2. THE MULTIWAVELENGTH PROPERTIES OF MRK 231

One of the first BAL quasars observed by NuSTAR was Mrk 231, a galaxy with Type 1 optical spectral classification (e.g., Osterbrock 1978). It is the nearest quasar with a redshift of 0.042 (Boksenberg et al. 1977). At its distance of only 183 Mpc, 1'' subtends 0.82 kpc. It is also one of the most infrared-luminous galaxies in the local universe (LIRL8 − 1000 μm = 3.6 × 1012L; Sanders et al. 2003) and is therefore also considered an ultraluminous infrared galaxy (ULIRG). Mrk 231 is a merger remnant, containing both an intense starburst and a luminous AGN. The starburst has a star formation rate (SFR) of ∼140 M yr−1 (Rupke & Veilleux 2013), estimated from its total infrared luminosity. However, there is considerable uncertainty in this value because of large scatters in the luminosity-to-SFR relations (see Section 5.2 for details). Modeling of the 1–1000 μm spectral energy distribution (SED) by Farrah et al. (2003) found an SFR of up to 450 M yr−1.

Since the SED of ULIRGs is dominated by emission in the infrared band, the bolometric luminosity of ULIRGs is assumed to be 1.15 times the total infrared luminosity (Lbol = 1.15 × LIR; Kim & Sanders 1998). The fractional contribution of the AGN and the starburst to the bolometric luminosity in composite sources like ULIRGs is very difficult to determine. There is ambiguity in the literature concerning the exact AGN contribution (between 30% and 70%) to the bolometric luminosity (∼1.5 × 1046 erg s−1) for Mrk 231. Modeling of the far-infrared SED by Farrah et al. (2003) implies ∼30% contribution from the AGN. More recently, Veilleux et al. (2009b) used six independent methods of measuring the AGN contribution to the bolometric luminosity in their sample of ULIRGs and quasars. These methods included fine-structure line ratios, the equivalent width (EW) of the 7.7 μm polycyclic aromatic hydrocarbon feature, and mid-infrared to far-infrared continuum flux ratios. For Mrk 231, these methods found the AGN contributes 64%–74%, for an average of ∼70%, to the bolometric luminosity. For simplicity, we will refer to the AGN portion of the bolometric luminosity (1.1 × 1046 erg s−1) as Lbol, AGN assuming a 70% AGN contribution to the bolometric luminosity.

Mrk 231 is probably the best local example of quasar-mode feedback, where the active nucleus is thought to be driving a powerful kiloparsec-scale outflow with velocities up to 1400 km s−1 (e.g., Rupke & Veilleux 2013). This wide-angle wind is multi-phased, containing ionized (Rupke & Veilleux 2013), neutral (Rupke & Veilleux 2013), and molecular outflow phases (Feruglio et al. 2010; Fischer et al. 2010; González-Alfonso et al. 2014). In the optical and UV, a BAL system is detected significantly in Mrk 231 in several transitions, including low-ionization species (e.g., Na i D, He i, Ca ii, Mg ii, Mg i, and Fe ii) and weakly in C iv (e.g., Adams & Weedman 1972; Gallagher et al. 2002, 2005; Veilleux et al. 2013). Thus, Mrk 231 is categorized as a rare (∼2% of the BAL quasar population) iron low-ionization BAL (FeLoBAL) quasar (Dai et al. 2012).

The structure of Mrk 231 is complex; it contains tidal tails, multiple molecular disks, and jets in addition to the host galaxy and the AGN. In Figure 1, we show a schematic of the structure of Mrk 231 and our assumed viewing geometry toward the nucleus described in this section. Although considered to be a radio-quiet object, Mrk 231 has a weak radio jet that extends from the parsec scale to the kiloparsec scale (Carilli et al. 1998; Ulvestad et al. 1999). On the kiloparsec scale, the jet becomes misaligned with the inner parsec-scale jet (Carilli et al. 1998; Ulvestad et al. 1999). The kiloparsec jet is thought to be providing additional acceleration to parts of the neutral outflow detected in spatially resolved maps in Na i D (Rupke & Veilleux 2013); single-dish observations in H i have also detected what is thought to be a weak jet-induced outflow (Morganti 2011; Teng et al. 2013). CO observations of Mrk 231 by Downes & Solomon (1998) revealed two distinct molecular disks. The denser inner molecular disk, extending from 75 to 1100 pc (≲ 1farcs3), is nearly face-on. This molecular disk is likely the site of intense star formation (see Section 5.2.3). Assuming an average dust temperature in the inner molecular disk of ∼200 K, heated by the AGN, the equivalent intrinsic column density derived from the Downes & Solomon (1998) measurements is NH ∼ 2.4 × 1023 cm−2. Given the orientation of the parsec-scale radio jet and assuming it is orthogonal to the accretion disk, the accretion disk has a different inclination than the inner molecular disk. The outer molecular disk extends out to a radius of ∼3'', overlapping with the starburst arc observed to have multiple UV and optical star-forming knots (Surace et al. 1998).

Figure 1.

Figure 1. Schematic showing the complexity of Mrk 231's structure and our assumed viewing geometry based on multiwavelength data. The host galaxy is face-on on the plane of the sky with optical emission extending out to a radius of 8 kpc (9farcs8; e.g., Rupke & Veilleux 2011). The thin (23 pc) molecular disk is tilted by 10° with the denser inner region extending from ∼75 pc (0farcs1) to 1.1 kpc (1farcs3) and the outer molecular disk extending out to 2.5 kpc (3farcs0; Downes & Solomon 1998). The accretion disk is assumed to be orthogonal to the parsec-scale jet, which is assumed to be inclined by 65° (see Section 5.2.2). In this figure, north is up.

Standard image High-resolution image

In the Chandra band, Mrk 231 is dominated by a point source, and has extended soft X-ray emission from the host galaxy (Gallagher et al. 2005). The far-UV is largely unobscured and unpolarized, but the near-UV and optical spectra (∼3000–4000 Å) appear to be dust reddened by a screen of AV ∼ 7 mag (e.g., Veilleux et al. 2013). Historic X-ray studies of Mrk 231 have found its spectrum to be relatively flat (Γ ∼ 1.2) suggestive of a reflection-dominated spectrum (Maloney & Reynolds 2000; Gallagher et al. 2002, 2005; Ptak et al. 2003). However, the low Fe Kα EW of ∼0.2 keV is difficult to reconcile with a significant reflection component. Maloney & Reynolds (2000) suggested that a large scattering fraction can dilute the Fe line EW. Braito et al. (2004, hereafter B04) and Piconcelli et al. (2013, hereafter P13) have both claimed detections above 10 keV that imply Compton-thick absorption. B04 measured an observed 0.5–10 keV luminosity of ∼3 × 1042 erg s−1; after correcting for absorption and removing the starburst component, the AGN luminosity in the 2–10 keV band is ∼5 × 1043 erg s−1, about 0.5% of the bolometric luminosity. Saez et al. (2012) found no significant variations in the intrinsic AGN luminosity in the 0.5–8 keV band. This is in agreement with P13 whose analysis suggests Compton-thick partial covering absorbers are responsible for the apparent flux variability below 10 keV seen on the time scale of years (e.g., Gallagher et al. 2002; Ptak et al. 2003). Changes in the absorber are likely responsible for the apparent variability observed in Mrk 231. These multiwavelength data imply complex geometries for the structure of Mrk 231.

3. OBSERVATIONS AND DATA REDUCTION

Mrk 231 was observed by NuSTAR for a total of 44 ks starting on UT 2012 August 27 (ObsID: 60002025002). In order to check for hard X-ray variability, it was observed again for a total of 29 ks starting on UT 2013 May 9 (ObsID: 60002025004). The first of these observations was scheduled to coincide with part of a 400 ks Chandra program with ACIS-S (PI: Veilleux; ObsID:13947).

The NuSTAR data were reduced with the pipeline software NuSTARDAS version 0.11.1 and CALDB version 20130509 with the standard corrections applied. After calibration and removal of bad time intervals, the final exposure is 41.1 and 28.6 ks on source for the first and second set of observations, respectively. For the first epoch of NuSTAR data, the absolute astrometric offset relative to the Very Large Array position of Mrk 231 (Becker et al. 1995) is ∼11'' and 10'' for focal plane module (FPM) A and B, respectively. Similarly, for the second epoch, the absolute astrometric offset is ∼8'' and 13'' for FPMA and B, respectively. The field of view (FOV) of each FPM is 12' by 12' (Harrison et al. 2013). On average, the first epoch of observations is off-axis by ∼4'. The second epoch is off-axis by an average of ∼2'.

The Chandra data were reduced using the standard ACIS-S scripts in CIAO version 4.5 and CALDB version 4.5.6. The Chandra data reduction is fully described in Veilleux et al. (2014). The Chandra data contemporaneous with the first epoch of NuSTAR observations were extracted to complement our analysis of the hard X-ray data. The strictly simultaneous good time intervals were extracted using the FTOOL "mgtime." There is no obvious temporal variability detected during either the NuSTAR or Chandra exposure. The two epochs of NuSTAR exposures also show no significant variability (see Section 5.1).

To be consistent between the two telescopes, source spectra were extracted using circles with 1' radii for both the NuSTAR and Chandra data centered at the X-ray peak. The Chandra data confirm that there are no contaminating point sources within the 1' source extraction regions. The Chandra background spectrum was extracted in the standard way, using a nearby circular extraction region where no obvious X-ray source is present. The NuSTAR background extraction is more complicated due to the inhomogeneity of the background throughout each of the two detectors, caused by stray light incompletely blocked by the aperture stop (hereafter referred to as the aperture background). A simulated total background that combines both the aperture and instrumental backgrounds was produced, following the method of background simulation detailed in Wik et al. (2014). Briefly, the background models for the NuSTAR exposures were normalized to "blank sky" observations of the COSMOS and the Extended Chandra Deep Field-South observations accumulated by NuSTAR through 2013 May. For each epoch, a background spectrum was simulated from the model and then scaled to match the size of the source extraction region. Conservatively, we estimate that the broadband systematic uncertainties in the derived background spectrum are ≲5%. In Figure 2, we show our first epoch of NuSTAR data before and after subtraction of the aperture background which dominates at ≲20 keV.

Figure 2.

Figure 2. First epoch of NuSTAR FPMA observations of Mrk 231 smoothed with a 16 pixel radius Gaussian kernel. Mrk 231 is clearly detected up to 20 keV. We show the raw data in the 3–20 keV band on the left and data after the subtraction of the aperture background in the 10–20 keV band to demonstrate the quality of our background model. In each panel, the white circle represents the 1' radius source extraction region and the red square marks the 12' × 12' FOV of NuSTAR. The color bar represents the total counts per pixel in log scale.

Standard image High-resolution image

Spectral analysis was performed using HEASoft version 6.13. The Chandra and NuSTAR source spectra were binned using the FTOOL "grppha" such that, in each bin, the total number of counts are four times the noise. In the first epoch of NuSTAR observations, FPMA and FPMB detected 1012 and 1015 counts, respectively, from Mrk 231 in the 3–30 keV band. Similarly, 819 and 771 counts were detected by FPMA and FPMB, respectively, in the second epoch. When modeling the spectra, an additional constant factor, typically on the order of a few percent, is applied to account for the cross-normalization between FPMA and B and between FPMA and Chandra ACIS-S. We assumed the Wilms et al. (2000) abundances and the Verner et al. (1996) photoelectric cross sections in our spectral modeling with XSPEC. The column density due to Galactic absorption is assumed to be NH, Gal = 1.26 × 1020 cm−2 (Dickey & Lockman 1990). All errors quoted in this paper are at the 90% confidence level (Δχ2 = 2.706 for a single parameter).

4. PREVIOUS HARD X-RAY OBSERVATIONS OF MRK 231

Previous detections of Mrk 231 in the hard X-ray (>10 keV) have been claimed using the non-imaging BeppoSAX Phoswich Detector System (PDS) (half-power diameter (HPD) ∼ 1fdg4; B04) and Suzaku Positive Intrinsic Negative detector (PIN detector) (HPD of 34' by 34' with a total square FOV ∼ 65farcm5 on each side; P13). The PDS and PIN 15–30 keV fluxes reported for Mrk 231 are 3.0 and 3.7 × 10−12 erg s−1 cm−2, respectively, in those publications. B04 concluded through their modeling of the broadband 0.5–50 keV spectrum that Mrk 231 is a Compton-thick quasar with a column density of ∼1.8–2.6 × 1024 cm−2. P13 supported this conclusion using their more recent Suzaku observations, but a comparison of the B04 data with their own suggested Mrk 231 is variable below 10 keV. They suggested that the variability is consistent with changes in the partial covering fraction of the Compton-thick absorber.

More recently, Koss et al. (2013) reported a faint detection (∼4σ) of Mrk 231 from the 70 month Swift Burst Alert Telescope (BAT) survey (Baumgartner et al. 2013). The BAT has a nominal point spread function (PSF) of ∼17', smaller than both Suzaku and BeppoSAX. By assuming a power law model (Γ = 1.9), they estimate the 14–195 keV BAT flux of Mrk 231 to be ∼5 × 10−12 erg s−1 cm−2, equivalent to ∼1.2 × 10−12 erg s−1 cm−2 in the 15–30 keV band, approximately half of the values measured by B04 and P13. For a random source to be considered a significant detection, the BAT team requires a flux limit of 4.8σ (Baumgartner et al. 2013). This blind detection limit is above the peak BAT signal of Mrk 231. Given the prior knowledge of the source position, the detection limit can be reduced to 3σ because knowing the position of an astronomical source reduces the likelihood that the peak signal at that location is random noise (Koss et al. 2013). However, when the full BAT band image is divided into eight sub-band images, the position of the BAT peak changes from one band to another by several arcminutes, suggesting a likely spurious detection (W. Baumgartner 2013, private communication). The BAT flux of Mrk 231, which is at least a factor of two lower than the BeppoSAX and Suzaku detections, suggests that part of the hard X-ray flux measured by BeppoSAX and Suzaku has been resolved out by the smaller PSF of BAT or is variable. However, the BAT detection has significant uncertainties which can be addressed by the NuSTAR observations.

5. MODELING THE X-RAY SPECTRUM OF MRK 231

5.1. The NuSTAR Data

The spectra from each epoch of NuSTAR observations are well-fit by a single power law model modified by only Galactic absorption and a narrow Gaussian to represent the Fe K emission. For the first epoch, the best-fit $\Gamma = 1.09^{+0.12}_{-0.11}$2 = 75.9 for 68 degrees of freedom). The observed 15–30 keV flux is 1.5 ± 0.2 × 10−12 erg s−1 cm−2, approximately the value measured by Koss et al. (2013) despite their steeper assumed spectral slope. The EW of the Fe K emission line at ∼6.5 keV is 0.29$^{+0.12}_{-0.15}$ keV. The second epoch spectrum has a similar best-fit Γ of 1.33 ± 0.14 (χ2 = 64.6 for 63 degrees of freedom) with an observed 15–30 keV flux of 1.2 ± 0.2 × 10−12 erg s−1 cm−2. Fe K emission is insignificant in the second epoch which has poorer photon statistics due to the shorter exposure. Although the two epochs of NuSTAR data are separated by about seven months, there is no obvious sign of variability above 10 keV. Therefore, to improve statistics, we created a combined spectrum of the two epochs of data using the "addascaspec" tool for each FPM. The summed spectrum is used for the remainder of this paper, and is best fit by a power law model with Γ = 1.22 ± 0.08 (χ2 = 105.5 for 108 degrees of freedom). This model implies a 15–30 keV flux of 1.4 ± 0.1 × 10−12 erg s−1 cm−2.

The derived 15–30 keV observed flux from the NuSTAR data are approximately a third to half of the values measured by B04 and P13; the NuSTAR spectrum is also significantly harder than those measured by B04 and P13. Comparing the NuSTAR best-fit model to the BeppoSAX PDS and Suzaku PIN data, we find that the PDS and PIN detections are brighter and have steeper slopes (Figure 3). The PIN data are only detected at ∼5% above the total background. We were unable to check the accuracy of the "tuned" non-X-ray background provided by the Suzaku team since Earth occulted observations are unavailable for these data. Taking into account the ∼3% systematic error of the "tuned" background,23 the lower limit of the PIN data shows an even steeper slope than the nominal value. The change in spectral shape could be due to either contamination or variability. The latter is unlikely because of the lack of variability in our NuSTAR epochs separated by seven months and the lack of intrinsic AGN variability measured in historic X-ray data (Gallagher et al. 2005; Saez et al. 2012). Unless the 10–30 keV X-ray spectrum of Mrk 231 hardened in the 14 months between the PIN and the first set of NuSTAR data, yet remained unchanged in the decade between the PDS and PIN observations, we conclude that the previous detections are contaminated by field sources outside of the NuSTAR FOV.

Figure 3.

Figure 3. Left: an image of the Mrk 231 field taken by XMM-Newton's EPIC-pn in 0.5–10 keV (ObsID: 0081340201). The data were minimally processed and smoothed with a four-pixel radius Gaussian and clearly show the richness of the field surrounding Mrk 231. The NuSTAR FOV of the first (white) and second (green) epochs are shown with the Mrk 231 spectral extraction region marked by the dashed cyan circle. The Suzaku PIN HPD (34' by 34') is also shown in magenta. The color bar represents log intensity in counts per pixel. Multiple sources within the PIN FOV may have contaminated the PIN and BeppoSAX PDS (HPD ∼ 1fdg4) data. Right: unfolded combined NuSTAR FPMA/B (red/blue) spectra of Mrk 231 with PIN data from P13 (black), and PDS data from B04 (green). The PIN data are nominally detected at only ∼5% above the total background and the black dashed lines approximate the upper and lower limits of the PIN detection when the ∼3% systematic error of the background is folded in. The model shown is a simple power law best fit to the NuSTAR data. The PIN spectrum (binned to 3σ) has a steeper slope than that found by NuSTAR. However, the PDS data (also binned to 3σ) is not inconsistent with either the NuSTAR or PIN spectrum.

Standard image High-resolution image

Although not discussed in B04 or P13, contamination from nearby sources can be problematic for detectors with limited angular resolution such as the PDS or a large single-pixel FOV such as the PIN. The PDS and PIN data match in flux and spectral index (Figure 3), implying that if a contaminant exists for both sets of observations, then it must lie within the PIN FOV. Taking advantage of the larger FOV of XMM-Newton EPIC-pn over that of Chandra ACIS-S, we show the 17 ks 0.5–10 keV image of the Mrk 231 field taken by XMM-Newton in 2001 (ObsID: 0081340201) in Figure 3. The field surrounding Mrk 231 is richly populated with point sources. Mrk 231 dominates the field in the 0.5–10 keV band, so it is reasonable to assume that it also dominates the field above 10 keV as other authors have done. It is unlikely that the source of contamination is a single object brighter than Mrk 231. However, with so many field sources, some may even be outside of the 27' by 26' XMM-Newton FOV; it is likely that the integrated spectrum of multiple sources fainter than Mrk 231 add to the PIN and PDS detections. A number of point sources detected by XMM-Newton within the NuSTAR FOVs are weakly detected by NuSTAR at 10–20 keV (see Figure 2) and there may be brighter sources outside the NuSTAR FOV contributing to the PIN detection.

While we cannot definitively rule out the possibility that the change in Mrk 231 flux above 10 keV is due to intrinsic variability, we conclude that the BeppoSAX and PIN detections were probably due to contamination for the following reasons.

  • 1.  
    The PIN detection was at only 5.1% above the background between 15 and 40 keV. At the 90% confidence level, the reproducibility of the total PIN background in this energy range is ∼3%.23 Therefore, about 5% of the time, the background is higher than that modeled, so the PIN detection could be consistent with a non-detection. Similarly, Mrk 231 was detected by BeppoSAX PDS at a rate of 0.066 ± 0.022 counts s−1 (Braito et al. 2004). The full band systematic error for PDS is 0.020 ± 0.015 counts s−1.24 Thus, the PDS detection is marginal and could also be consistent with a non-detection.
  • 2.  
    Using the hard band (2.5–8 keV) log N − log S (e.g., Kim et al. 2007) derived from random fields from all over the sky and assuming a flux S equal to the hard band flux of Mrk 231, there is a 33% chance a background source is of that flux or brighter within 1 deg2 (FOV of the PIN). Since the PIN response falls off linearly with off-axis angle, to produce the assumed detected flux, it is most likely that the source has to be brighter than S. Conservatively, within the central 0.25 deg2, the chances of contamination by a source of flux ⩾2S is better than 3%. Because the PIN's lack of imaging capability cannot rule out serendipitous sources, we expect the chance of a source falling within the PIN FOV of 1 deg2 is better than 10%. This probability increases when we consider sources fainter than 2S contaminating the field and detectors with larger FOV such as that of BeppoSAX. The 2.5–8 keV log N − log S may not match perfectly with the 10–30 keV log N − log S, but the only log N − log S currently available above 10 keV is from the BAT survey which only covers the brightest sources.
  • 3.  
    As stated above, the NuSTAR 10–30 keV flux is consistent with the BAT results integrated over nearly five years (Koss et al. 2013). Since the temporal coverage is so long, we can assume that the BAT value is the typical average flux of Mrk 231. The derived BeppoSAX and PIN fluxes are more than twice the BAT limit. Therefore, it would be a remarkable coincidence if the AGN flared twice only at times observed by detectors with much larger PSFs than NuSTAR.
  • 4.  
    If the AGN flared when observed by BeppoSAX and PIN but not when observed by NuSTAR, the lack of variability below 10 keV or change in the Fe line EW is puzzling (Gallagher et al. 2005; Saez et al. 2012). It would require a contrived scenario to have a Compton-thick cloud present to obscure the hard X-ray continuum only when the 10–30 keV flux is in the high state so that continuum variability is not detected below 10 keV or without change to the Fe line EW.

Therefore, we conclude that contamination is the simplest explanation given all the circumstances.

5.2. The Broadband Chandra and NuSTAR Spectrum

Our Chandra and NuSTAR observations of Mrk 231 comprise the most sensitive broadband X-ray spectrum of this object to date. In this section, we explore various models in order to constrain the three-dimensional geometry of the AGN and its associated absorber system as well as the intrinsic luminosity of the AGN. As stated previously, the NuSTAR-only spectra are well-fit by a single power law model. With the Chandra data anchoring the lower energies, we modeled the full 0.5–30 keV spectrum (χ2 = 314.7 for 261 degrees of freedom) with a MEKAL gas temperature of 0.82 ± 0.04 keV and a very flat power law ($\Gamma = 1.20^{+0.09}_{-0.08}$) modified by a small column ($N_{\rm H} = 2.9^{+0.7}_{-0.5} \times 10^{22}$ cm−2). There is no telltale sign of a Compton-reflection hump above 10 keV that would be indicative of significant reflection; this is corroborated by the weak Fe K line (EW ∼ 0.30 keV). If Mrk 231 were to contain a normal AGN with a typical photon index (e.g., 1.5 < Γ < 2.2; Nandra & Pounds 1994; Reeves & Turner 2000), then its apparent flat X-ray spectrum would require more complex models and geometries.

Since Mrk 231 is unresolved by the NuSTAR beam, star formation is responsible for a non-negligible fraction of the total X-ray emission detected by NuSTAR. For this reason, our modeling will include a starburst and an AGN component. The starburst component consists of a MEKAL model representing the hot gas associated with the starburst and a power law model with Γ = 1.1 and cutoff energy of 10 keV representing the non-thermal component of high mass X-ray binaries (HMXBs; e.g., Nagase 1989; Persic & Raphaeli 2002; Lutovinov et al. 2005; B04). The AGN component is represented by different versions of an absorbed power law model.

In the following subsections, we present a detailed discussion of our spectral modeling of the broadband X-ray spectra. We first applied simple geometry models such as those presented by other authors (e.g., B04; P13). We then used the more sophisticated MYTorus model (Murphy & Yaqoob 2009) to model the AGN emission. In the process, we find that it is difficult to constrain the various starburst and AGN components independently; these models result in poor constraints on the HMXB component. The derived HMXB luminosities from SFR–LX scaling relations are too high given the infrared-derived SFR (Lehmer et al. 2010; Mineo et al. 2012b) due to degeneracies in the spectral models. Exploiting the spatial resolution of Chandra, we separated out the X-ray emission into nuclear and host galaxy components. By fixing the host galaxy component, we are able to robustly measure the AGN emission. We conclude that the AGN is extremely under luminous in the X-rays.

5.2.1. Simple Geometries

First, we modeled the broadband spectrum from the 1' aperture with the AGN component being a power law modified by a single partial covering absorber. The AGN photon index (ΓAGN) tended toward a flatter slope than the canonical value of 1.8 (ΓAGN ∼ 1.4). We tried fixing ΓAGN to 1.8, but the resulting fit is a poor model of the spectrum which underestimates the flux above ∼10 keV. Thus, we allowed ΓAGN to be free. The best-fit values are listed in Column 2 of Table 1 and the spectrum is shown in Figure 4. The derived gas temperature of ∼0.81 keV is consistent with those previously measured in ULIRGs (e.g., Ptak et al. 2003; Teng & Veilleux 2010). Assuming a total SFR of ∼140 M yr−1 derived from the total infrared luminosity of the non-AGN component, the expected X-ray luminosity of the starburst is derived using the SFR–LX relations found by Mineo et al. (2012a, 2012b). The luminosity we find for the hot gas component, ∼1.5 × 1041 erg s−1, is consistent with the expected value of ∼7.3 × 1040 erg s−1 to within the 0.34 dex scatter of the Mineo et al. (2012b) relation. However, the 0.5–8 keV luminosity of the HMXB power law component is nearly a factor of nine lower than the expected value of ∼3.6 × 1041 erg s−1 based on the Mineo et al. (2012a) relation. We assumed no obscuration for this component, but an average column of ∼3 × 1024 cm−2 is required in order for the derived HMXB luminosity to be consistent with prediction. It is improbable that such a significant column is present over the whole galaxy. The absorption-corrected 2–10 keV AGN luminosity only accounts for ∼0.05% of Lbol, AGN, far lower than the typical value of ∼2% in luminous quasars (Elvis et al. 1994).

Figure 4.

Figure 4. Chandra (black) and NuSTAR FPMA/B (red/blue) spectra of Mrk 231 shown with various models. (A) Single power law modified by a single partial covering neutral absorber. (B) MYTorus without independent constraints on the starburst component. (C) MYTorus model with the X-ray emission from the star formation of the host galaxy constrained. (D) The best-fit MYTorus model simultaneously fit to 2003 (green) and 2012 (black) Chandra data. All models result in a similar quality of fit, though as described in Section 5, several of the models are deemed improbable because of required galaxy-scale Compton-thick columns or atypical SFR-to-X-ray luminosity values.

Standard image High-resolution image

Table 1. Best-fit Parameters

Model Partial Covering MYTorusa MYTorusa Comment on
Parameter Absorber (2012) (2003+2012) Parameter
(1) (2) (3) (4) (5)
B/A 0.97$^{+0.07}_{-0.07}$ 0.97$^{+0.07}_{-0.07}$ 0.97$^{+0.07}_{-0.07}$ FPMB-A cross-normalization
CXO/A 1.00$^{+0.10}_{-0.10}$ 0.97$^{+0.11}_{-0.10}$ 0.90$^{+0.09}_{-0.08}$ Chandra-FPMA cross-normalization
kT (keV) 0.82$^{+0.07}_{-0.09}$ 0.25$^{+0.05}_{-0.06}$, 0.93$^{+0.11}_{-0.10}$ 0.26$^{+0.06}_{-0.05}$, 0.87$^{+0.09}_{-0.08}$ MEKAL gas temperature from the starburst
Abs. 1 (× 1022 cm−2)b 6.57$^{+1.82}_{-1.40}$ 11.2$^{+2.9}_{-2.4}$ 19.4$^{+5.7}_{-4.4}$c, 9.5$^{+2.3}_{-1.9}$d Neutral absorber 1
cf 1 0.82$^{+0.12}_{-0.03}$ 1 (f) 1 (f) Covering factor 1
Abs. 2 (× 1022 cm−2)e  ⋅⋅⋅  24 (f) 24 (f) Neutral absorber 2
ΓHMXB 1.1 (f) 1.1 (f) 1.1 (f) HMXB cutoff power law index with cutoff energy at 10 keV
ΓAGN 1.39$^{+0.10}_{-0.11}$ 1.40$^{+0.04}_{...}$ 1.40$^{+0.03}_{...}$ AGN power law index (MYTorus lower limit fixed at 1.4)
Inc (°)  ⋅⋅⋅  65 (f) 65 (f) Inclination angle
Eline(keV) 6.63$^{+0.07}_{-0.06}$ 6.67$^{+0.03}_{-0.10}$ 6.66$^{+0.04}_{-0.04}$ Narrow ionized Fe line energy (Fe Kα self consistent within MYTorus)
EWline(keV) 0.168$^{+0.093}_{-0.095}$ 0.146$^{+0.121}_{-0.075}$ 0.166$^{+0.119}_{-0.070}$ Ionized Fe line equivalent width
Const. (C-thin)  ⋅⋅⋅  0.19$^{+0.04}_{-0.03}$ 0.23$^{+0.04}_{-0.03}$c, 0.20$^{+0.04}_{-0.04}$d Compton-thin fraction
f0.5 − 2 (× 10−13 erg s−1 cm−2) $1.07^{+1.76}_{-0.03}$ 1.16$^{+0.05}_{-0.10}$ $1.19^{+0.04}_{-0.10}$c, $1.16^{+0.04}_{-0.06}$d Observed 0.5–2 keV flux
f2 − 10 (× 10−13 erg s−1 cm−2) $9.42^{+3.59}_{-0.89}$ 9.20$^{+0.31}_{-2.30}$ $7.69^{+0.20}_{-1.75}$c, $8.72^{+0.27}_{-1.94}$d Observed 2–10 keV flux
f10 − 30 (× 10−12 erg s−1 cm−2) $1.82^{+0.22}_{-0.37}$ 1.77$^{+0.01}_{-0.62}$ $1.75^{+0.01}_{-0.56}$ Observed 10–30 keV flux
LMEKAL (erg s−1) 1.51 × 1041 2.30 × 1041 2.35 × 1041 Intrinsic MEKAL 0.5–30 keV luminosity
LHMXB (erg s−1) 8.53 × 1039 8.57 × 1041 8.57 × 1041 Intrinsic HMXB 0.5–30 keV luminosity
L0 (erg s−1) 1.40 × 1043 1.13 × 1043 1.01 × 1043 Intrinsic AGN 0.5–30 keV luminosity
L0 (2–10 keV) (erg s−1) 4.90 × 1042 3.94 × 1042 3.84 × 1042 Intrinsic AGN 2–10 keV luminosity
$L_{\rm C{\rm -}thin}$ (erg s−1)  ⋅⋅⋅  2.18 × 1042 2.38 × 1042c, 2.23 × 1042d Compton-thin AGN 0.5–30 keV luminosity
LAGN/Lbol, AGN (%) 0.05 0.03 0.03 2–10 keV X-ray-to-bolometric luminosity ratio for the AGN
χ2/d.o.f. 258.3/258 249.7/258 390.2/380 Goodness of fit

Notes. Column 1: model parameter for each fit. Column 2: a single partial covering absorber model. Column 3: MYTorus fits to the contemporaneous Chandra and NuSTAR data. Column 4: same as Column 3, but fitting 2003 and 2012 Chandra data simultaneously. Column 5: comments on model parameter. aBest-fit model: Const.×NH, Galactic(MEKAL1+MEKAL2+AbsnuclearHMXB ×cutoffPLnuclearHMXB+Line(6.7 keV)+MYTorus×PLAGN+Const.$_{\rm C{\rm -}thin}\times$PLAGN). bAbsorbing column of the torus in the MYTorus fits. c2003 values. d2012 values. eAbsorbing column of the HMXB component in the MYTorus fits.

Download table as:  ASCIITypeset image

Following Gallagher et al. (2005) and P13, who modified the B04 model of Mrk 231 to explain the spectral variability at <10 keV, we applied a double partial covering absorber model to the new broadband data. Such a model is a poor fit to the data. Since there is no sign of Compton-reflection in the NuSTAR spectrum, the reflection component is statistically unnecessary (Δχ2 = 3.7 for 2 degrees of freedom). The need for the reflection component by P13 to fit their Suzaku data is possibly due to the likely contamination discussed in Section 4. Therefore, we do not discuss this model further.

5.2.2. Toroidal Geometries

We also implemented a toroidal model that self-consistently accounts for the transmitted and scattered AGN emission in a torus geometry in the Compton-thick regime (MYTorus; Murphy & Yaqoob 2009). While the model assumes only neutral material, it is also relevant for partially ionized material since Compton scattering dominates in the Compton-thick regime. The MYTorus model is a table model that contains three main components: the direct (or zeroth order) AGN component, the Compton-scattered component, and some fluorescence emission line components for self-consistent modeling. An additional Compton-thin power law component can be added to represent "leaked" emission due to a partial covering absorber. The MYTorus table models are relevant for Γ in the range of 1.4–2.5 and NH up to 1025 cm−2.

The MYTorus model defaults to a half-opening angle of 60° which corresponds to a covering factor of 0.5. Therefore, for an inclination angle of 60°, the sightline to the center skirts the edge of the absorber. In the modeling, the strength of the Fe K line provides the main constraints on the inclination angle and the column density. Since the Fe K complex in Mrk 231 is so weak, we held the inclination angle fixed at 65° in order to derive a robust value for the column density. The assumption of 65° is based on several observational factors including the position angle of the inner parsec-scale jet (Ulvestad et al. 1999), which does not contribute significantly to the X-ray emission (Gallagher et al. 2002). Additionally, the profile of the Lyα line suggests the accretion disk must have an inclined geometry (inclination ≳45°; Veilleux et al. 2013). The partial covering result of P13 suggests that the line of sight through the torus must be near the "fuzzy" edge of the torus (Murphy & Yaqoob 2009). The robustness of our results based on the assumed inclination of 65° is discussed below in Section 5.2.3.

By allowing all the relevant model parameters to be free, the best-fit model requires a Compton-thin component, which is akin to the "leaked" direct emission in a partial covering absorber model favored by P13. The best-fit model suggests that the AGN is border-line Compton-thick ($N_{\rm H} \sim 2.1^{+1.7}_{-1.7}\times 10^{24}$ cm−2). However, the starburst element is poorly constrained. The derived luminosity for the HMXB component is unphysical (more than 10 times too luminous) with respect to the SFR-to-X-ray luminosity scaling relations assuming an infrared-derived SFR of 140 M yr−1 (Lehmer et al. 2010; Mineo et al. 2012b).

To better constrain the HMXB contribution to the broadband X-ray spectrum, we took advantage of Chandra's excellent spatial resolution and excised the inner 1'' (0.82 kpc) emission from the spectrum. Using this host-only Chandra spectrum, we measured the expected contribution of the starburst to the X-ray emission separate from that of the AGN. Accounting for the wings of the ACIS PSF (∼20% of the total counts lie beyond 1''), we find no significant contribution from HMXBs in the host spectrum (e.g., Veilleux et al. 2014). The host-only Chandra spectrum is best fit by two MEKAL components with temperatures at 0.25$^{+0.05}_{-0.06}$ and 0.93$^{+0.11}_{-0.10}$ keV. The requirement of a two-temperature MEKAL component is consistent with the two temperatures derived for the extended X-ray emission by Gallagher et al. (2002). The best-fit spectrum is shown in Figure 4. Contrary to the results from when the starburst component is allowed to vary, Mrk 231 is, in fact, Compton-thin ($N_{\rm H} \sim 1.3^{+0.03}_{-0.02}\times 10^{23}$ cm−2), but with ΓAGN of 1.40–1.47. Additionally, compared to Lbol, AGN, the intrinsic 2–10 keV X-ray luminosity of Mrk 231 is only 0.04%. The modeling of the X-ray emission from circumnuclear star formation is discussed in the following subsection.

5.2.3. The Preferred Model and the Variability of Mrk 231

From the previous sections, it is clear that accurate accounting of the starburst contribution to the X-ray luminosity is key to finding a robust measurement of the intrinsic AGN luminosity. Due to the lack of spatial resolution and the effects of obscuration within 1'' of the nucleus, it is difficult to separate the starburst and the AGN emission in the X-rays. Spatially resolved integral field unit maps of Mrk 231 in Hα from Rupke & Veilleux (2013) place poor constraints on the SFR in the nuclear region because of heavy extinction. Similar maps in Paα taken recently with adaptive optics show that the majority of star formation is taking place within the inner 0.5 kpc (0farcs6) diameter of the nucleus (D. Rupke 2013, private communication).

Consistent with the Paα constraints, CO measurements from Downes & Solomon (1998) imply that most of the molecular material for star formation is actually concentrated within the inner 1farcs3 (1.1 kpc). Assuming the IRAM 30 m single-dish CO(1–0) flux of 97 Jy km s−1 is emitted by a molecular disk with minor and major axes of 300 and 340 pc, respectively, and a height of 23 pc derived from CO rotation curves (Downes & Solomon 1998), the CO(1–0) flux implies a nuclear SFR in the range of ∼100–350 M yr−1 using the molecular gas density to SFR density relation of Kennicutt (1998). For these calculations, the CO-to-H2 mass conversion (αCO) is assumed to be 3–7 M (K km s−1 pc2)−1 based on the two-phase local velocity gradient models of the circumnuclear molecular gas by Papadopoulos et al. (2012). Given that the host-only Chandra spectrum shows no significant contribution from HMXBs, we will assume all of the X-ray luminosity from HMXBs arises from circumnuclear star formation within 1''. The CO flux also implies a column density of NH ∼ 2.4 × 1023 cm−2. Therefore, the X-ray emission from HMXBs associated with the nuclear starburst would be heavily absorbed and thus not contribute to the X-ray spectrum of Mrk 231 below 2 keV, as discussed by Gallagher et al. (2002).

To account for the nuclear starburst, we added a new obscured HMXB component to our modeling. In addition, deeper Chandra data show an ionized Fe line at ∼6.7 keV, likely arising from the circumnuclear starburst (e.g., Strickland & Heckman 2007). This line is not self-consistently accounted for by MYTorus, thus we include an additional Gaussian feature to model this line. To recap, the current model now includes a two-component MEKAL model, an Fe line at 6.7 keV, and a heavily obscured HMXB component for the nuclear star formation (SFR of 140 M yr−1), and a leaky MYTorus component for the AGN emission (see Table 1 notes for an equation form of this model). For the nuclear starburst component, we assumed values that are consistent with the CO measurements where the absorbing column is assumed to be 2.4 × 1023 cm−2 and the unabsorbed 0.5–8 keV flux is consistent with the Mineo et al. (2012a) relation for SFR of 140 M yr−1. All other components in the model are allowed to vary.

The best-fit values are listed in Column 3 of Table 1 and the best-fit spectrum (shown with the model components) is in Figure 5. The EW of the neutral Fe line compared to the zeroth order power law model is ∼0.04 keV and that of the ionized Fe line is ∼0.17 keV. If the 6.7 keV line comes from thermal emission of a starburst, then it should be accompanied by a bremsstrahlung continuum. However, an additional bremsstrahlung continuum component to the nuclear starburst model contributes very little to the overall fit (Δχ2 ∼ 0.6). This bremsstrahlung component is assumed to have a peak temperature of 4.5 keV in order to produce the line at 6.7 keV with an unabsorbed 0.5–2 keV luminosity equal to the expected value for a 140 M yr−1 SFR and modified by the same column as the HMXB component. The 6.7 keV line is likely associated with the HMXBs. The ionized line seems uncommon in high resolution low-mass X-ray binary spectra, but is also rare in HMXBs where only 20% were detected in the high spectral resolution survey of 10 HMXBs conducted by Torrejon et al. (2010). For the HMXB Cen X-3, its 6.7 keV line has EW of 0.02–0.4 keV (Naik et al. 2011). Assuming this is typical of HMXBs, because few high resolution measurements are available, and diluting the EW of the line by 80% to account for non-detections, the expected EW from HMXBs would be between 0.004–0.1 keV, approximately consistent with the EW range we measured in Mrk 231 (0.06–0.25 keV). The slightly higher EW measured in Mrk 231 can be partially accounted for by the AGN.

Figure 5.

Figure 5. Unfolded Chandra and NuSTAR spectra shown with the best-fit model components (dashed lines) with the total model in cyan. The starburst MEKAL and HMXB components are in green. The AGN component includes Compton-thin AGN emission unaffected by the absorber (gray), the "zeroth" order direct emission from the AGN transmitted through the absorber (purple), the scattered AGN emission through the absorber (orange), and the Fe K emission lines modeled self-consistently within MYTorus as well as the 6.7 keV Fe line (magenta).

Standard image High-resolution image

The derived ΓAGN of 1.4 is still flatter than the usual range, but not inconsistent with these values. Fixing the value to 1.8 results in a poorer fit (Δχ2 = 42.9 for a single parameter) and underestimates the flux above 10 keV. The intrinsic AGN luminosity remains weak, at only ∼0.03% of the bolometric luminosity. The line-of-sight column density is Compton-thin with $N_{\rm H}=1.1^{+0.3}_{-0.2}\times 10^{23}$ cm−2. A Compton-thin AGN is consistent with the weak Fe line emission (EW ∼ 0.2 keV).

The inclination angle of the torus has so far been fixed at 65°; this inclination angle implies a sightline through the edge of the obscuring torus. We also tested geometries where the sightline does not go through the torus at less than 60° (the opening angle assumed by the MYTorus model). These resulted in poor fits where the model underestimates the emission above 10 keV and overestimates the emission below 5 keV. If the inclination angle is fixed at 90°, the model parameters and fit statistics are approximately the same as those of the best-fit model listed in Column 3 of Table 1 except for the column density. For a geometry looking through the thickest part of the torus, the column density is patchy and Compton-thin ($N_{\rm H} = 1.1^{+0.1}_{-0.1}\times 10^{22}$ cm−2).

The inclination angle constraint in MYTorus is mainly driven by the strength of the Fe line and the shape of the hard X-ray spectrum. Since the Fe line observed in Mrk 231 is so weak and there is no evidence of a Compton-reflection hump, the MYTorus model cannot fully constrain both the inclination angle and the covering factor. Following Yaqoob (2012), we used the MYTorus model to constrain the line-of-sight column density and place an upper limit on the average column density by decoupling the scattered AGN component from that of the direct zeroth order AGN component. By assuming the best-fit model, we find that the line-of-sight column density is ∼5.8 × 1022 cm−2 and the upper limit to the average column density is ∼4.0 × 1024 cm−2. The decoupled column densities are in agreement with the limits we derived from our spectral modeling assuming a toroidal geometry and a sightline cutting through the torus.

Both Gallagher et al. (2005) and P13 suggested that the observed flux of Mrk 231 below 2 keV is variable, but Saez et al. (2012) noted that the intrinsic AGN luminosity of Mrk 231 does not vary. We re-reduced a 40 ks Chandra observation from 2003 (ObsID: 04028) presented in Gallagher et al. (2005) with the same calibration as our 2012 Chandra data so that the time-dependent ACIS efficiency is taken into account. Compared to the 2003 Chandra data, the observed 0.5–2 keV and 2–8 keV fluxes in 2012 have changed by $+3^{+4}_{-45}$% and $-18^{+3}_{-19}$%, respectively. We tested the robustness of our MYTorus results by modeling the NuSTAR data simultaneously with both the 2003 and 2012 sets of Chandra data. Assuming that the intrinsic luminosity of the AGN and the starburst component have not changed as discussed in Sections 4 and 5.1, we allowed the AGN column density and the AGN Compton-thin luminosity fraction to vary between the 2003 and 2012 data sets in order to explain the observed spectral variability. The best-fit results are listed in Column 4 of Table 1 and shown in Figure 4. The spectral variability is well-explained by patchy absorption varying between 1.2–2.0 × 1023 cm−2. This model implies an intrinsic 2–10 keV to bolometric luminosity ratio of 0.03%; thus, Mrk 231 is intrinsically X-ray weak.

6. THE X-RAY WEAK AGN IN MRK 231

6.1. The Flat ΓAGN of Mrk 231

The physically motivated MYTorus model appears to provide an excellent description of our broadband X-ray spectrum of Mrk 231. The best-fit model implies a Compton-thin absorbing column (NH ∼ 1.2 × 1023 cm−2). Although ΓAGN of 1.4 is not very much lower than the typical range seen in AGN and the upper limit is within the expected range, the derived value is an artificial constraint as it is limited by the lower bound (1.4) of the MYTorus table model. It is possible that the true value of ΓAGN is even lower than 1.4. As discussed in Section 2, the flatness of the Mrk 231 spectrum below 10 keV has been assumed to be due to a strong reflection component, and therefore the absorbing column is Compton-thick. It is not until our NuSTAR observations that the lack of a distinct reflection hump above 10 keV is definitively shown.

One explanation for the flat ΓAGN is that Mrk 231 is so obscured that even MYTorus, with an NH upper bound of 1025 cm−2, is ineffective in fully describing the broadband spectrum because the actual column is above the upper bound. This is unlikely. If the intrinsic AGN luminosity were higher, then in order to produce the observed spectrum, a thicker absorbing column is required. This in turn implies that more of the direct AGN emission is Compton scattered by the torus. Since the observed soft X-ray spectrum is almost completely accounted for by the starburst component, there is little leeway in the observed spectrum to accommodate more scattered emission below 2 keV. If the absorbing column is above 1025 cm−2 and has a high covering factor (≈4π), then it is unlikely that a bright optical nucleus would be detected. It is also possible that the assumed half-opening angle of MYTorus (60°) is incorrect. We examined this possibility by applying a version of the MYTorus model25 available to us that assumes a half-opening angle of 37° to our data. The model parameters from this test are consistent with those derived from the standard MYTorus model and a sightline going through the torus.

Another possibility is that the flat ΓAGN is an intrinsic property of the AGN in Mrk 231. In ULIRG-type sources, black hole growth occurs most rapidly after coalescence and during the obscured phase (e.g., Martínez-Sansigre et al. 2005; Teng & Veilleux 2010). The high accretion rate of the black hole also drives feedback which disperses the obscuring material. In most cases, a steep ΓAGN is associated with high Eddington ratios (e.g., Shemmer et al. 2008; Brightman et al. 2013). In this case, perhaps a disturbed corona, a product of the ongoing high accretion rate of the AGN and powerful outflow, produces a flatter than typical photon index. Galactic X-ray binaries in the "ultra soft state" are characterized by a high Eddington ratio but low X-ray luminosity (e.g., Done et al. 2007). Mrk 231 being an AGN in a rare ultra soft state is consistent with the lack of intrinsic variability and the low X-ray-to-optical power law slope (αOX, see below). However, ΓAGN for the ultra soft state would be ∼2 (e.g., Done et al. 2007), steeper than our measurement. A flat Γ (<1.4) implies the lack of low-energy seed photons that can be Compton up-scattered to higher energies (e.g., Zdziarski 1988; Fabian et al. 1988); this does not seem likely for Mrk 231 given its UV and optical properties. Furthermore, given the super-Eddington accretion rate measured for Mrk 231 (see Section 6.2), a standard scenario of an advection-dominated accretion flow (e.g., Narayan & Yi 1994, 1995), where the flat ΓAGN mimics a thermal bremsstrahlung spectrum (e.g., Pozdnyakov et al. 1983), is unlikely.

It is also possible that Mrk 231 is in fact a radio-loud object. The radio-loud BAL quasar PG 1004+130 is thought to be a jetted source that is likely to be intrinsically X-ray weak (Luo et al. 2013). X-ray spectra of radio-loud sources can have relatively flat slopes (e.g., NGC 4261; Sambruna et al. 2003). If the AGN corona is radiatively inefficient and is feeding a jet, then the X-ray emission may be dominated by inverse Compton, resulting in a flat ΓAGN (e.g., Markoff et al. 2005). Recent Very Long Baseline Array monitoring of Mrk 231 by Reynolds et al. (2013) discovered a blazar-like radio flare in Mrk 231. The first epoch of our NuSTAR data was taken before the most recent radio flare while the second epoch was taken soon after the peak of the flare. If the flat ΓAGN is related to the launching of a jet or a flare, then one would likely expect changes in the X-ray slope as flares come and go. We find no significant changes in ΓAGN between the two epochs of NuSTAR data nor significant variability in ΓAGN and X-ray flux in historic X-ray data of Mrk 231. In addition, the expected X-ray flux from synchrotron self-Compton derived from the radio flux density of the flare is only ∼10−21 erg s−1 cm−2, far below the measured X-ray flux. Therefore, it is unlikely that the flares contribute significantly to the observed X-ray properties of Mrk 231.

6.2. Mrk 231 is an Outlier among Quasars

Regardless of the reason behind the flat photon index of Mrk 231, the absorption-corrected luminosities from all of our modeling imply that the total intrinsic AGN luminosity for Mrk 231 in 0.5–30 keV is ∼1.0 × 1043 erg s−1. Compared to Lbol, AGN of 1.1 × 1046 erg s−1, the 2–10 keV X-ray luminosity is 0.03–0.05% of the AGN bolometric luminosity. Even if we consider that the AGN only contributes to ∼30% of the bolometric luminosity as derived by the modeling of the far-infrared SED (Lbol, AGN = 4.7 × 1045 erg s−1; Farrah et al. 2003), the ratio between 2–10 keV luminosity and Lbol, AGN is still ∼0.08%. The derived AGN fraction is much lower than the value expected for this type of source. For comparison, Seyfert 1 galaxies and radio-quiet quasars sampled by Elvis et al. (1994) for which X-ray data were available have X-ray to bolometric ratios of ∼2 to 15%, with the most luminous objects typically having the lowest ratios. However, this does not seem to be unusual for ULIRG-type sources. Teng & Veilleux (2010) found that Seyfert 1-like ULIRGs have mean 2–10 keV to bolometric ratios of ∼0.1% compared to the mean of ∼3% for PG quasars.

Recent studies by Vasudevan & Fabian (2007, 2009) and Lusso et al. (2010, 2012) have independently confirmed that AGN with the highest Eddington ratios require the highest bolometric corrections. Assuming the dynamically derived black hole mass of $1.7^{+4.0}_{-1.2}\times 10^7$M using stellar velocity dispersions measured from near-infrared spectra (Dasyra et al. 2006)26 and Lbol, AGN of 1.1 × 1046 erg s−1, the Eddington ratio of the AGN in Mrk 231 is 5.0$^{+11.3}_{-3.5}$, assuming isotropic emission. The 2–10 keV to bolometric luminosity ratio for an AGN with such an Eddington ratio is expected to be ∼0.3% (Lusso et al. 2010, 2012) to ∼0.7% (Vasudevan & Fabian 2009). However, the correlation has a large dispersion and lacks sampling at the very X-ray weak end. Therefore, it is difficult to assess whether the significant X-ray weakness can be explained simply by super-Eddington accretion.

Although rare, Mrk 231 is not the only intrinsically X-ray weak AGN known. The X-ray-to-optical power law slope (αOX; Figure 6) of Mrk 231 is not as weak as that of the known X-ray weak AGN PHL 1811, an X-ray variable AGN with a steep photon index (ΓAGN = 2.3; Leighly et al. 2007a, 2007b). X-ray stacking analysis of a sample of radio-quiet high-redshift (z ∼ 2.2) analogs of PHL 1811 shows that, on average, they have Γeff ∼ 1.1 (Wu et al. 2011). These sources appear to have large blueshifted C iv absorption line EWs, likely due to winds dominating their broad emission line regions (Wu et al. 2011). Since PHL 1811 itself has a soft spectrum of ΓAGN = 2.3, its analogs are plausibly expected to have similarly soft X-ray spectra. Therefore, the flat effective X-ray photon indices measured from the analogs are interpreted as indications of heavy obscuration and probably even reflection-dominated spectra. To place Mrk 231 in the context of other LoBAL quasars and X-ray weak AGN, we redshifted our best-fit model of Mrk 231 to the redshift of the other radio-quiet LoBAL quasar, PG 1700+518, observed by NuSTAR (z = 0.292; Luo et al. 2013). The effective photon index of the redshifted spectrum (Γeff ∼ 0.9) is similar to that measured in PG 1700+518 (Γeff ∼ 0.5). Mrk 231 seems to be a local example of an X-ray faint LoBAL (17–40% of which are intrinsically X-ray weak; Luo et al. 2013) that innately possesses a weak X-ray continuum.

Figure 6.

Figure 6. αOX vs. the 2500 Å monochromatic luminosity for Mrk 231. The 2500 Å flux of Mrk 231 has been corrected for dust reddening based on modeling of the UV continuum by Veilleux et al. (2013). The range of αOX for Mrk 231 plotted here represent the observed total 2 keV flux (−1.92) and the absorption-corrected direct AGN 2 keV flux (−1.71). Also included are the X-ray weak AGN PHL 1811 (Leighly et al. 2007a, 2007b) and BAL quasar PG 1700+518 (Luo et al. 2013). The gray points are normal quasars from Steffen et al. (2006) and the gray line is their αOX–optical luminosity relation for radio-quiet quasars. The dotted lines mark deviations of −0.50 and −0.80 from the Steffen et al. (2006) relation. The X-ray weakness of Mrk 231 is distinct compared to the "normal" AGN studied by Steffen et al. (2006).

Standard image High-resolution image

The X-ray weakness of the ionizing continuum is also manifested in the lack of [O iv] (25.89 μm), [Ne v] (14.32, 24.32 μm), and [Ne vi] (7.65 μm) line emission from Mrk 231 (Armus et al. 2007; Farrah et al. 2007; Veilleux et al. 2009b). These lines arise in the narrow line region, located above the torus, and are thus less likely to be affected by obscuration. The [O iv], [Ne v] and [Ne vi] lines require sufficient high-energy photons (>77 eV) to be produced. The weak continuum of Mrk 231 becomes more evident if we compare it to another local (z = 0.043) AGN ULIRG, IRAS F05189–2524. Both ULIRGs are in the same merger stage (Veilleux et al. 2009a), have approximately the same dynamically derived black hole masses (Dasyra et al. 2006), have nearly identical infrared SEDs (Armus et al. 2007), and have almost the same Lbol, AGN (Veilleux et al. 2009b). However, they have very different derived Eddington ratios (∼1 versus ∼5). Given the uncertainties in the measured black hole masses and Eddington ratios, this is only a qualitative comparison between Mrk 231 and IRAS F05189–2524. The comparison is valid because the black hole masses were measured using the same method and thus have similar systemic measurement errors. IRAS F05189–2524 has strong [O iv], [Ne v], and [Ne vi] emission lines and a strong intrinsic X-ray continuum (Ptak et al. 2003; Teng et al. 2009) while Mrk 231 is the opposite. Although we only have a limited sample, the distinction between these two sources may explain the large scatter of [O iv]-to-bolometric luminosity relations observed in nearby AGN (e.g., Melendez et al. 2008; Rigby et al. 2009; Weaver et al. 2010). The intrinsic X-ray weakness of Mrk 231 is likely associated with the super-Eddington accretion rate of the black hole that also drives its powerful outflows (Rupke & Veilleux 2013).

7. CONCLUSIONS

We have obtained the most sensitive broadband (0.5–30 keV) spectrum of Mrk 231, the closest BAL quasar and an ULIRG, to date. NuSTAR detects Mrk 231 at a fainter level above 10 keV than previous observations with poorer angular resolution by BeppoSAX PDS and Suzaku PIN. Two epochs of NuSTAR data separated by 7 months show no evidence of hard X-ray variability so contamination is the likely cause of the discrepancy between NuSTAR and previous, non-focusing hard X-ray observations.

Analysis of contemporaneous Chandra and NuSTAR broadband spectra find significant X-ray emission from a powerful circumnuclear starburst. The direct AGN emission is absorbed and scattered by a patchy torus that is Compton-thin with $N_{\rm H}\sim 1.2^{+0.3}_{-0.3}\times 10^{23}$ cm−2. The Compton-thin absorption is compatible with the apparent weak Fe K emission (EW ∼ 0.2 keV).

Mrk 231 appears intrinsically X-ray weak. The absorption-corrected 2–10 keV luminosity is ≲ 0.1% of the AGN bolometric luminosity with αOX ∼ −1.7. A reason for the low X-ray luminosity could be that the AGN is in a rare ultra soft state similar to those seen in X-ray binaries, though this is inconsistent with the measured ΓAGN. It is unclear why the AGN is intrinsically X-ray weak, but there are not very many known examples of such sources as points of comparison. The weak ionizing continuum explains the lack of mid-infrared [O iv], [N v], and [Ne vi] fine-structure emission lines which are present in sources with otherwise similar AGN properties. The intrinsic X-ray weakness may be a result of the super-Eddington accretion occurring in the nucleus of the ULIRG. These results on Mrk 231 are consistent with its being a merger remnant emerging from its dusty cocoon through a powerful wind event.

We are grateful to the anonymous referee for providing useful comments which improved our manuscript. We thank Wayne Baumgartner, Bret Lehmer, Richard Mushotzky, Jeremy Schnittman, Tahir Yaqoob, and Andreas Zezas for useful discussions. We would also like to thank Lee Armus who provided useful comments in the early planning phase of the NuSTAR ULIRG program. We also thank Roberto Maiolino, David Rupke, and Eckhard Sturm who are co-investigators of the Chandra program. This work was supported under NASA Contract No. NNG08FD60C, and made use of data from the NuSTAR mission, a project led by the California Institute of Technology, managed by the Jet Propulsion Laboratory, and funded by the National Aeronautics and Space Administration. 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). The scientific results reported in this article are based in part on observations made by the Chandra X-Ray Observatory and data obtained from the Chandra Data Archive published previously in cited articles. This work, in part, made use of observations obtained with XMM-Newton, an ESA science mission with instruments and contributions directly funded by ESA Member States and the USA (NASA). We made use of the NASA/IPAC Extragalactic Database (NED), which is operated by the Jet Propulsion Laboratory, Caltech, under contract with NASA. S.H.T. is supported by a NASA Postdoctoral Program (NPP) Fellowship. W.N.B. and B.L. acknowledge support by California Institute of Technology (Caltech) NuSTAR subcontract 44A-1092750 and NASA ADP grant NNX10AC99G. F.E.B. acknowledges support from Basal-CATA (PFB-06/2007) and CONICYT-Chile (under grants FONDECYT 1101024 and Anillo ACT1101). A.C. acknowledges support from ASI/INAF grant I/037/12/0-011/13. P.G. acknowledges support from STFC grant reference ST/J003697/1.

Facilities: NuSTAR - The NuSTAR (Nuclear Spectroscopic Telescope Array) mission, CXO - Chandra X-ray Observatory satellite

Footnotes

  • 23 
  • 24 

    http://heasarc.gsfc.nasa.gov/docs/sax/abc/saxabc/saxabc.html#SECTION00052300000000000000. The typical systematic error in the PDS count rate is given by a 1998 technical report by M. Guainazzi and A. Matteuzzi that is no longer available publicly. The authors compared PDS spectra from a number of blank fields to derive the residual count rate. Since the original document is no longer available online, the error is assumed to be at 90% confidence level similar to the confidence quoted on other errors in the current document and is typical of X-ray studies.

  • 25 

    These models were provided by T. Yaqoob with new calculations only relevant to the continuum (i.e., the transmitted and Compton-scattered components), so the strength of the Fe line is not taken into account. This version of the model may not be reliable below ∼3 keV.

  • 26 

    There is some uncertainty to the derived black hole mass and therefore the Eddington ratio. Using Hβ line widths, Kawakatu et al. (2007) found the black hole in Mrk 231 is 8.7 × 107M with an uncertainty that is within a factor of ∼3. Therefore, the Hβ value is broadly consistent with the Dasyra et al. (2006) value that we have assumed. The Hβ spectra were taken using a 3'' aperture, so there may be contamination from the host. Thus, we consider the Hβ-derived black hole mass an upper limit, leading to an Eddington ratio lower limit of ∼1.1.

Please wait… references are loading.
10.1088/0004-637X/785/1/19