© C. Andonie et al. 2022

1. Introduction

X-ray emission is a universal characteristic of active galactic nuclei (AGN), thought to arise from inverse Compton scattering of optical-UV photons from the accretion disk by hot electrons in the corona (e.g., Haardt & Maraschi 1991). The intrinsic X-ray emission takes a power-law spectral form (f(E)∝E^−Γ, with typical photon indices of ⟨Γ⟩∼1.8–2.1; e.g., Nandra & Pounds 1994; Winter et al. 2009; Corral et al. 2011), but it can be modified due to an interaction with matter in the vicinity of the central supermassive black hole (SMBH). In particular, Compton scattering and photoelectric absorption of the primary X-ray continuum lead to two important features in the X-ray spectrum: the Fe Kα emission line and the so-called Compton-hump. By studying these reprocessed features, together known as the so-called AGN reflection component, we can infer the physical properties of the matter from which they originate, and hence probe the circumnuclear environments of central SMBHs.

The Fe Kα feature at 6.4 keV is produced by fluorescence processes related to the absorption of higher energy X-ray photons by neutral Fe atoms. Its spectral profile is generally comprised of broad and narrow components. The narrow component of the Fe Kα line (FWHM ≲ 10 000 km s⁻¹; e.g., Lubiński & Zdziarski 2001; Yaqoob & Padmanabhan 2004; Shu et al. 2010) is a ubiquitous spectral feature of AGN, and in a majority of cases the only component immediately visible, while the broad component is harder to pin down since it requires exceptional statistics and broad energy coverage to decouple the line from the underlying continuum and absorption components (e.g., Guainazzi et al. 2006; Marinucci et al. 2014). Nonetheless, when present, reverberation studies suggest that the broad component originates from a compact zone, only a few Schwarzschild radius (r_g) in extent, around the SMBH (e.g., Cackett et al. 2014), and hence it is strongly affected by Doppler and gravitational broadening (e.g., Mushotzky et al. 1995; Tanaka et al. 1995; Yaqoob et al. 1995). On the other hand, the narrow component has been thought to be produced somewhere amongst the outer accretion disk, the broad line region (BLR), and the torus clouds (e.g., Ghisellini et al. 1994; Krolik et al. 1994; Yaqoob et al. 1995), corresponding to months-to-year variability timescales; although, contributions from the more distant narrow line region (NLR) and ionization cone have been observed (e.g., Wang et al. 2011; Fabbiano et al. 2017). In practice, the line likely has contributions from all of the above structures.

One of the more straightforward narrow-line constraints that provided sufficient spectral resolution is the full width at half maximum (FWHM), which traces the spatially unresolved kinematics of the circumnuclear matter and hence can be used to estimate the average reprocessor location. Different studies have arrived at different conclusions with respect to the location of the Fe Kα emitting regions. Based on a sample of 14 bright AGN observed with the high energy grating (HEG) of the Chandra X-ray observatory (Weisskopf et al. 2000), Nandra (2006) found a lack of correlation between the Fe Kα FWHM and either the optically derived BLR line width or the SMBH mass, and they concluded that the Fe Kα core likely has a mix of contributions from the outer accretion disk, BLR, and torus, in differing proportions depending on the source, but it predominantly originates in regions outside the BLR, possibly near the inner edge of the torus. Shu et al. (2010) expanded upon this, using a sample of 36 nearby AGN with HEG spectra (27 with FWHM constraints), and arrived at similar conclusions. Later, Gandhi et al. (2015) analyzed 13 local type 1 AGN, also using HEG spectra, and found that the Fe Kα sizes, as estimated from the Fe Kα FWHM, appeared to be bounded by the dust sublimation radii (i.e., the inner wall of the torus), and they may predominantly originate in clouds associated with either the inner edge of the torus, the BLR, or even further inside. Among type 2 AGN, Shu et al. (2011) found no obvious differences compared to type 1 AGN, while Marinucci et al. (2012) presented a time, spectral, and imaging analysis of NGC 4945, showing that the reflecting structure is at a distance ≥30–50 pc, which is much larger than the typical torus scales.

Rapid X-ray continuum variability is commonly observed in unobscured, obscured (e.g., Uttley & McHardy 2005), and even some heavily obscured AGN (Puccetti et al. 2014) and suggests that the primary X-ray continuum emitting source (i.e., the corona) is produced in a compact zone very near to the SMBH (e.g., Mushotzky et al. 1993; De Marco et al. 2013). The X-ray light curve can be analyzed via the power-spectral density (PSD) function, which is typically characterized as a power-law of the form P_ν ∝ ν^α, where ν is the temporal frequency and α is the power-law slope (e.g., Green et al. 1993; Edelson & Nandra 1999; Vaughan et al. 2003b). Typical values for the power-law slope in AGN are α ∼ −1 at low frequencies, indicative of pink noise, and α ≲ −2 at high frequencies, indicative of red-noise (e.g., Vaughan et al. 2003b; Markowitz et al. 2003; McHardy et al. 2004, 2005; Uttley & McHardy 2005; Summons et al. 2007; Arévalo et al. 2008). The transition in the PSD between these two regimes is denoted as ν_B = 1/t_B, and is related to the characteristic X-ray variability timescales of the system (e.g., Markowitz et al. 2003; McHardy et al. 2006; González-Martín & Vaughan 2012).

Many studies have investigated the correlated variability from the X-ray continuum and the broad Fe Kα via reverberation mapping to constrain SMBH spin (see recent review by Uttley et al. 2014), but far fewer have investigated the relation between the variability of the X-ray continuum and the narrow Fe Kα line. The latter have focused on campaigns of individual nearby sources like MCG-6-30-15 (Iwasawa et al. 1996), MRK 509 (Ponti et al. 2013), NGC 2992 (Marinucci et al. 2018), NGC 4051 (Lamer et al. 2003), NGC 4151 (Zoghbi et al. 2019), and NGC 7314 (Yaqoob et al. 1996). In at least a few of these, the authors were able to place constraints on the location and size of the Fe Kα emitting region by studying the reaction of the Fe Kα line to continuum variations (e.g., Ponti et al. 2013; Zoghbi et al. 2019). Surprisingly, some studies found a tight correlation between the narrow Fe Kα line and the continuum on observational timescales of a few days, implying that the narrow Fe Kα component predominantly arises from regions interior to the BLR.

While timing and multiepoch spectral investigations can place important constraints on reflection close to a SMBH (e.g., light hours-to-years scales), a number of studies based on Chandra observations have also found Fe Kα emission extending out hundreds to thousands of light years in galaxies such as NGC 4151 (Wang et al. 2011), NGC 6240 (Wang et al. 2014), NGC 4945 (Marinucci et al. 2017), and ESO 428-G014 (e.g., Fabbiano et al. 2017). These scales are much larger than the putative size range of the dusty torus, which ALMA and IR reverberation studies have shown to be ≲10 pc (e.g., Gallimore et al. 2016; Imanishi et al. 2016; García-Burillo et al. 2016; Lyu & Rieke 2021). The fractional contribution from such highly extended reflection, however, is generally not dominant.

The goal of the present work is to enhance our understanding of the reflecting cloud structure in AGN. In particular, we aim to constrain the location(s) and size(s) where the narrow Fe Kα line is produced, and understand how such regions may vary among a large AGN sample, particularly as a function of various AGN and host galaxy properties. To this end, we carry out spectral, timing, and imaging analyses on a large ensemble of Chandra and XMM-Newton observations, investigating the FWHM of the Fe Kα line, its temporal properties both alone and as they relate to those of the X-ray continuum, and its potential spatial extent. The paper is organized as follows. In Sect. 2 we describe the observations and data reduction, in Sect. 3 we explain the X-ray the spectral fitting, while in Sect. 4 we present an analysis of the light curves. In Sect. 3.2 we analyze the Fe Kα line FWHM of the sample, in Sects. 4.1.2 and 4.2 we investigate the Fe Kα and continuum variability, and in Sect. 5 we outline our assessment of the Chandra images. We conclude with some discussion in Sect. 6 and a summary in Sect. 7.

2. Observations and data reduction

In the following subsections, we describe how we select the sample of bright, mostly local AGN, and the observations and data reduction for Chandra and XMM-Newton.

2.1. Sample selection

Our broad goal is to quantify the spectral and temporal properties of the Fe Kα line in AGN, and relate these to the variable X-ray continuum; we focus on the 2–10 keV continuum band in order to minimize spectral complexity associated with the soft excess and host contamination (e.g., Fabbiano 2006; Done et al. 2012). To start, we consider that the typical Fe Kα equivalent width of AGN ranges from ≈0.1–1 keV (e.g., Shu et al. 2010, 2011), which implies that ∼1–10% of the total photons in the 2–10 keV band will be Fe Kα photons assuming typical AGN spectra. Thus a clear requirement emerges such that the AGN be observed by a facility with a high 2–10 keV sensitivity. We therefore focus on observations from Chandra (Weisskopf et al. 2000) and XMM-Newton (Jansen et al. 2001), whereby even a relatively short 10 ks exposure of a typical AGN (Γ = 1.9 power-law) with a flux of f_2−10 keV = 4 × 10⁻¹² erg s⁻¹ cm⁻² yields ∼1700 (ACIS-I) and ∼5000 (pn) 2–10 keV photons, respectively. Such limits produce what we consider to be the bare minimum in terms of Fe Kα photon statistics (i.e., > 20–50 counts in the line) to enable spectroscopic constraints for typical observations (≥10–20 ks)¹. Clearly XMM-Newton observations provide better spectral and timing statistics, but suffer substantial background flaring and potential contamination from off-nuclear emission. On the other hand, Chandra observations can be more versatile since they offer higher spatial resolution to search for extended Fe Kα emission on ∼100-pc to kpc scales, and, when the HEG is deployed, sufficient spectral resolution to resolve the Fe Kα line.

To obtain the broadest possible sample, we adopt as our parent input sample the most recent 105-month Swift-Burst Alert Telescope (BAT) Survey (BASS, Oh et al. 2018), an all-sky survey in the ultra-hard X-ray band (14–195 keV), which provides a relatively unbiased AGN sample at least up to N_H ≲ 10²⁴ cm⁻² (Winter et al. 2009; Ricci et al. 2015; Koss et al. 2016). The 105-month Swift-BAT catalog is a uniform hard X-ray all-sky survey with a sensitivity of 8.4 × 10⁻¹² erg s⁻¹ cm⁻² over 90% of the sky in the 14–195 keV band. The survey catalogs 1632 hard X-ray sources, 947 of which are securely classified as AGN. We include one additional target in our sample, the well-known narrow-line Sy1 1H0707−495, which is relatively bright in the 2–10 keV band and is undetected in the BAT 105-month catalog, probably due to its very steep spectral index (e.g., Done et al. 2007; Boller et al. 2021). To limit our analysis to only those observations for which we have a high likelihood of constraining the Fe Kα line, we apply a flux cut of f_{14−195 keV} ≥ 10⁻¹¹ erg s⁻¹ cm⁻² or f_2−10 keV ≥ 4 × 10⁻¹² erg s⁻¹ cm⁻² (these are roughly equivalent for a Γ = 1.9 power-law); this resulted in the selection of 252 sources in the local universe (z < 0.1), and 28 more distant galaxies with redshifts between 0.1 and 0.56. In order to assess off-nuclear contamination and carry out both imaging and timing analysis (to assess long-term variability), we require a minimum of at least five Chandra observations; this reduces our final sample to 38 objects. However, for historical and technical reasons which we outline below, we reduce and extract Chandra spectra for all 280 sources that satisfy the flux cut previously described, in order to calibrate the flux of the annular spectra (see Appendix B for details). We complement our final sample of 38 sources with XMM-Newton pn observations when available. Table A.1 in the Appendix summarizes the observations analyzed in this work.

We stress that while this final sample is by no means complete, it spans a reasonable range of the parameter space to be considered representative of local hard X-ray selected AGN. The top panel of Fig. 1 shows the intrinsic X-ray luminosities and redshifts of our sample, while the bottom panel shows the X-ray photon index distribution. The X-ray luminosity and photon index values are taken from Ricci et al. (2017a), while redshifts and AGN properties used in this work are from Koss et al. (2017) and BASS DR2 (Koss et al. 2022), which extends the DR1 results of Koss et al. (2017) and Ricci et al. (2017a). A large majority of the sample reside in the local universe (z < 0.1) and have modest X-ray luminosities (L_2−10 keV < 10⁴⁴ erg s⁻¹). Additionally, the photon indices of our sample range from 1.3 to 3.13, with a median value of Γ = 1.82, which is broadly consistent with the range and median photon index of the nonblazar AGN in the BASS survey (Γ = 1.78 ± 0.01 Ricci et al. 2017a). The distributions of other AGN properties such as the line-of-sight column density [log(N_H/cm⁻²) = 20–25.9, with median of 21.5], AGN type (predominantly Seyfert types Sy1, Sy1.2, Sy1.5, Sy1.8 and Sy2), black hole mass [log(M_BH/M_⊙) = 6.13–9.83, with median 7.6] and Eddington ratio (ER = 3.4 × 10⁻⁷–1.13, with median 0.09) also span comparable ranges to the z < 0.85 nonblazar AGN in the BASS survey, and hence should be as representative (compare Figs. E.1 and E.2 with Fig. 13 of Koss et al. 2017).

Fig. 1. Some X-ray properties of our sample. Top figure: intrinsic 2–10 keV luminosity versus redshift of our sample. Magenta circles are radio-loud sources (see Table 1 for details), black circles are Compton-thick AGNs, green circles are dual AGNs, and gray circles are normal, radio-quiet AGNs. Bottom figure: distribution of X-ray photon indices of our sample.

2.2. Chandra X-ray Observatory

We downloaded 1001 observations with the Advanced CCD Imaging Spectrometer (ACIS), both the ACIS-S and ACIS-I CCD configurations, available from the Chandra archive associated with the 280 sources in our sample; we excluded the small number of High Resolution Camera observations, due to that instrument’s much lower spectral resolution and poorer high-energy sensitivity. Among the 1001 ACIS observations, 279 were acquired in the High-Energy Transmission Grating mode (Canizares et al. 2000, HETG). The HETG consists of two grating assemblies: the High Energy Grating (HEG, 0.8–10 keV) and the Medium Energy Grating (MEG, 0.4–5 keV); we only use the HEG spectra, since our focus is on the study of the 2–10 keV spectra, and in particular the region around the 6.4 keV Fe Kα line to make our analysis. The observations were acquired between 2000 and 2018.

The data reduction follows standard procedures with the CIAO software (v4.11) and calibration files (CALDB v4.8.3), using the chandra_repro script. The X-ray peaks in images of the same object are occasionally shifted up to ≈1″ from each other, or from the established optical centroid of the galaxy. To correct this, we manually choose the center of each observation and estimate the alignment offset with respect to the observations with higher exposure time, create a new aspect solution using the wcs_update task, and update the astrometry. We resort to manual alignment for two reasons. First, in some observations, there are not enough point sources in the field in common to match and align them automatically, and second, nuclear saturation (i.e., pileup) is very high in some observations, producing a hole in the center of the source which can confuse simple detection and alignment codes. Therefore, we apply new aspect solutions based on the different X-ray observations. Then, we remove background flares from the event files with the script deflare and subtract the readout streak, if present above the background, with the acisreadcorr task, to allow for more precise spectral and imaging analysis.

Finally, we reproject the events to a common tangent point using the reproject_events task and merge them for image analysis.

2.2.1. Spectral extraction

For both normal ACIS and zeroth-order HETG observations, we extract spectra with the specextract task, which creates spectra and responses files (ARF and RMF). We generate source spectra for each observation using both a 15 radius circular aperture and an annulus of 3″–5″, adopting an annular region of 20″–35″ aperture to estimate the background spectra. We mask all obvious off-axis point sources using the observation with the highest exposure time as our reference image. We use the wavdetect task to detect all the point sources in the field, setting the wavelet scale parameters to one and two pixels. By visual inspection, we confirmed that all point sources were detected and masked correctly. Then, we mask the off-axis point sources in the field using the default 3σ elliptical regions from the wavdetect output file before the spectral extraction². We use the circular extracted spectra for the observations unaffected by pileup and the annular spectra otherwise. We apply aperture corrections to both the 1.5″ circular spectrum and the annular spectra using the arfcorr task, which are applied through modified responses files.

For HETG observations, we extract 1st-order dispersed spectra using a 6-pixel width (≈3″) aperture. We first create a mask to delineate events in the HEG and MEG arms using the tg_create_mask task, and then we resolve the spectral orders using the tg_resolve_events task. Lastly, we generate spectra and response files using tgextract and mktgresp, respectively.

A key concern is that the spectra are not highly affected by off-axis sources or extended emission. As mentioned above, before extracting the spectra, we remove all contaminating sources in the field. For the nuclear 15 and grating spectra, we confirm that there is generally little contamination, due to the small aperture. However, for sources affected by extended emission (see Appendix B for further details), contamination can be larger and we do not use the 3″–5″ spectra in such cases.

2.2.2. Imaging analysis

To complement the spectral analysis, we also investigate the high angular resolution Chandra images to see whether some sources are spatially extended and compare this with the spectral results. To assess this, we extract radial profiles. Specifically, for each source, we merge all the event files with an off-axis angle < 2′ (to avoid observations affected by strong distortions to the point spread function, or PSF), using the task reproject_obs, which reprojects the observations to a common tangent point and then merges them. We mask off-nuclear point sources in the field, as well as dispersed photons related to the 1st-order spectrum for grating-mode observations. Then, we extract radial profiles for each source with dmextract, using sequential annuli, accounting for the masked area from each annulus.

2.3. XMM-Newton Observatory

When available, we complement our data with observations from the pn camera of XMM-Newton for 32 sources among our sample of AGN with more than five Chandra observations. We do not include observations from the EPIC-MOS cameras, as many of the observations are likely affected by pileup due to the high typical readout time (2.6 s). The pn camera has a full-frame time resolution of 73.3 ms per CCD (Jansen et al. 2001), such that the observations generally do not suffer much from pileup, making them well suited for variability analysis. When available, we use 118 spectra of 13 AGN, extracted from 30″ aperture provided by Tortosa et al. (in prep.). Additionally, we downloaded 172 60″-aperture pn spectra for 29 sources from the 4XMM–DR9 catalog (Webb et al. 2020). The 30″ aperture spectra are preferred, when available, to minimize host contamination. The procedure for the extraction and data reduction in Tortosa et al. (in prep) is described below. For details of the extraction and data reduction of sources in the 4XMM–DR9 catalog we refer to Webb et al. (2020).

The extraction of calibrated XMM-Newton spectra is performed by means of the Science Analysis System (SAS) software package (v.18.0.0). The XMM-Newton pn data are processed using the task epchain in order to obtain calibrated and concatenated event lists. XMM-Newton can also focus charged particles on the detection plane. Since these background flares have a large effect on the detected X-ray spectrum of all three EPIC cameras aboard XMM-Newton, it is imperative to remove these events during the data reduction process. To this end, the evselect command examines the count-rate of such events, selecting all single-pixel events (i.e., PATTERN==0) in the energy range sensitive to soft proton flares. Source and background spectra are extracted again using evselect, adopting a 30″ radius circular source region and a nearby source-free 50″ radius circular background region. The Redistribution Matrix Files (RMF) are generated using the command rmfgen, the Auxiliary Response Files (ARF) are generated with the command arfgen. When the input spectrum is multiplied by the ARF, the result will be the distribution of counts as would be seen by a detector with ideal resolution in energy. Then, the RMF is needed, in order to produce the final spectrum. All these files are compressed into a single file, easily readable by XSPEC (Arnaud 1996), using the SAS tool specgroup. EPIC spectra were binned in order to over-sample the instrumental resolution by at least a factor of 3 and to have no less than 30 counts in each background-subtracted spectral channel.

The 30″ and 60″ aperture XMM-Newton spectra are more likely to be affected by potential galaxy contamination than the nuclear Chandra spectra. We estimated the amount of contamination in the XMM-Newton spectra by computing the ratio between the encircled counts in a 15 aperture and the encircled counts in a 30″ and 60″ aperture of the merged Chandra images for each galaxy, after proper PSF correction. We found that the galaxy contamination for Centaurus A (Cen A) and Cygnus A dominates the spectra, contributing more than 60% of the total flux, and therefore the XMM-Newton spectral epochs for those AGN were not used. Four additional sources (IC 4329A, M51a, H1821+643 and 1H0707-495) have host galaxy contributions between 20–40% of the total flux, while for the rest of the sample the host galaxy contributions are less than 10–20%. We expect that the contribution of the host contamination will be roughly constant between the observations, and thus apply a correction to the measured flux by a constant factor. Thus, spatially resolved host contamination should be not have a significant impact on subsequent variability results.

3. X-ray spectral analysis

We fit the X-ray spectra of 652 Chandra and XMM-Newton observations in our sample. We describe our spectral fitting approach and the results for the Fe Kα line profile below.

3.1. Spectral fitting

We fit the spectra using the package PyXspec³, which is a python version of the popular standalone X-ray spectral fitting package XSPEC (version 12.10.1, Arnaud 1996). We apply a simple model to compute the fluxes of the Fe Kα line and the continuum between 2–10 keV, given by:

(1)

where the PHABS₁ and PHABS₂ components account for the Galactic and intrinsic AGN line-of-sight photoelectric absorption, respectively, the POWERLAW component corresponds to the continuum emission, and the ZGAUSS term reproduces the Fe Kα line emission at 6.4 keV. The Galactic absorption is fixed at the value obtained from the HI 4π survey (HI4PI Collaboration 2016). We leave free to vary the normalization and spectral index of POWERLAW, and line center (within 10% of its rest-frame energy, adopting redshifts fixed at their published values), the normalization and line width of ZGAUSS. The CPFLUX₁ and CPFLUX₂ components provide the unabsorbed photon fluxes of the continuum emission and the Fe Kα, respectively. To compute unabsorbed fluxes, we first freeze the free parameters after minimizing the fit and set the absorbing column densities to zero. Figure 2 shows an example of the fit for two observations of Cen A in instrumental (top panel) and unfolded (bottom panel) units. The model reproduces fairly well the 2–10 keV spectra, which are strongly attenuated by absorption, and yields a secure measure of the well-known emission line at 6.4 keV. We can see strong variability between the observations of about one order of magnitude for both the continuum and the line. Table A.1 lists the best-fitting parameters of the observations analyzed in this work. The fits are carried out with Cash statistics (Cash 1979).

Fig. 2. 2–10 keV spectra of the Type 2 AGN Centaurus A, extracted and aperture-corrected from 3″–5″ annuli in ACIS-I ObsIDs 7800 (black) and 8489 (red), fitted by the model described in Eq. (1). The top panel shows the spectra and model-to-data ratio in instrument units while the bottom panel shows the unfolded spectra. The spectra are attenuated by strong absorption. The sole emission line is the 6.4 keV Fe Kα line.

This model is simple but sufficient to measure the fluxes of unobscured sources. We confirmed this by replacing ZGAUSS by the PEXMON model, which reproduces the Compton neutral reflection with self-consistent Fe Kα, Fe Kβ and Ni Kα lines. We applied this model to three sources with different levels of obscuration: Circinus Galaxy (N_H = 4 × 10²⁴ cm⁻², Ricci et al. 2017a), Cen A (N_H = 1 × 10²³ cm⁻², Ricci et al. 2017a), and NGC 4051 (N_H ≲ 10²⁰ cm⁻², Ricci et al. 2017a). For the last two sources, the continuum and Fe Kα line fluxes are consistent at a 90% confidence level, which is expected since at N_H < 10²⁴ cm⁻² the spectrum is dominated by the transmitted component. For the Circinus Galaxy, however, the model with reflection predicts an unabsorbed continuum flux twice that of the simple model, although the Fe Kα fluxes are consistent. The PEXMON model predicts a higher column density and softer photon index, such that an increase in the continuum flux is expected. At N_H > 10²⁴ cm⁻² and energies lower than 10 keV, the spectrum is dominated by the reflected component, therefore it is not surprising that the PEXMON model predicts a different continuum flux. We do not consider this to be a problem, since the majority of our sources have N_H < 10²⁴ cm⁻² (see Fig. E.1c) and we take special care in the interpretation of the results for more heavily obscured sources.

To properly compare the fluxes measured from different spectra, we need to know whether they are affected by pileup, which occurs when two or more photons fall in the same pixel within a single readout “frame” and, in consequence, the information in that pixel is altered. If the source is very bright or the readout time is high, there is a higher probability to be affected. Normal CCD observations with the ACIS camera of Chandra are commonly affected by pileup given the sharp on-axis PSF and the nominal 3.2 s frame time. Chandra first-order grating spectra are much less affected by pileup since the photons are more dispersed. Therefore, if an observation is performed in HETG grating mode, we rely on the HEG 1st-order spectrum. If an observation is performed in normal CCD mode, we first estimate the amount of pileup to decide whether to use the 1.5″ circular aperture spectrum or the annular spectrum. The XMM-Newton pn observations are less affected by pileup due to the much shorter readout time of that instrument and the fact that the PSF is spread over many more pixels.

One issue we discovered during this process is a large inconsistency between the fluxes measured by the HEG first-order spectrum and the annular spectrum from the zeroth-order spectrum of the same data, indicating that the aperture correction performed by the task specextract is not accurate. Appendix B explains in detail how we calibrate an aperture correction factor between the HEG and 3″–5″ spectra, in order to incorporate the fluxes measured by the annular spectra in our analysis.

3.2. Fe Kα line FWHM

The Fe Kα line profile provides information about the kinematics of the reflecting clouds from which the line originates. We fit the Fe Kα FWHM for all sources with available Chandra HEG spectra, which provides the best spectral resolution among current observatories (FWHM ≈ 1860 km s⁻¹, or 39 eV at 6.4 keV); we note that the values reported in Table 2 are deconvolved FWHMs, as the line spread function information from the RMF is used during the fit. Among all the HEG observations, we only consider model fits with a lower limit different from zero and a finite upper limit at a 90% confidence. While this could bias our results against sources where the Fe Kα line is either faint or observed in a low state, we want to avoid fitting poorly detected lines where the profile could attempt to fit the underlying continuum, yielding inaccurate or overestimated line widths.

We estimate the radius of the narrow Fe Kα line emission (R_Fe Kα) assuming virial motion with the gravitational potential dominated by the SMBH of the emitting material as

(2)

where G is the gravitational constant and M_BH is the SMBH mass (e.g., Netzer 1990; Peterson et al. 2004). If the Fe Kα line originates in an outflow instead of the BLR or the dusty torus, the kinematics can no longer be modeled by virial motion. Nevertheless, for an outflow to reach velocities comparable to the typical υ_FWHM values of the narrow Fe Kα, large AGN luminosities are needed (log(L_bol/erg s⁻¹) > 45, Bischetti et al. 2017) and only a few sources of our sample are this powerful.

We compare the R_Fe Kα with the inner radii of the dusty torus, that is the dust sublimation radius, R_sub, given by Nenkova et al. (2008) as

(3)

where L_bol is the bolometric luminosity and T_sub is the dust sublimation temperature, generally assumed to be the sublimation temperature of graphite grains, T ≈ 1500 K (e.g., Kishimoto et al. 2007). We adopt a temperature uncertainty of ΔT = 500 K to compute R_sub, which broadly accounts for possible variations due to grain mineralogy, porosity and size, among others. We get the values of L_bol and M_BH from the DR2 of BASS (Koss et al. 2022). The values of L_bol and M_BH are tabulated in Table 1. We also compare R_Fe Kα to the optical BLR radius R_Hβ, inferred from Hβ reverberation studies of eleven sources (Bentz et al. 2009; Du et al. 2016; Zhang et al. 2019), as listed in Table 1.

We provide an estimate for the Fe Kα emission location for 24 out of 38 sources in our sample. Table 2 lists the R_Fe Kα, υ_FWHM and R_sub values. The left and right panels of Fig. 3 compare the location of the Fe Kα emitting regions to dust sublimation radius, R_sub, and R_Hβ, which we take to be the nominal optical BLR region radius. For 21 out of 24 sources, the bulk of the Fe Kα emission appears to arise from regions inside the dust sublimation radius, while for eight out of 11 sources, the Fe Kα emission appears to arise from regions near or a factor of several beyond R_Hβ. Thus, for most AGN, the narrow Fe Kα emission appears to originate primarily in the outer BLR. We cannot exclude that small portions of the narrow Fe Kα flux may arise from the outer accretion disk or the dusty torus, and additionally we note that a small minority (∼20%) of AGN diverge from this general behavior.

Fig. 3. Left: Fe Kα emitting radius (RFe Kα) versus dust sublimation radius (Rsub) in units of parsecs. Rsub values are estimated using Eq. (3) and RFe Kα values are estimated using Eq. (2). The gray shaded area represents RFe Kα ≥ Rsub. Right: Fe Kα emitting radius (RFe Kα) versus optical BLR radius (RHβ) in units of parsecs. RHβ is directly calculated from Hβ reverberation lags. The gray shaded area represents RFe Kα ≥ RHβ.

4. X-ray light curve analysis

We construct light curves for each source using the recalibrated and aperture-corrected Fe Kα line (F_Fe Kα) and the 2–10 keV continuum (F_2−10 keV) fluxes. Appendix Figs. D.1 and D.2 provide light curves for the continuum and Fe Kα line fluxes, respectively, for the entire sample.

In the following subsections, we characterize the variability of the light curves, estimate the size of the reprocessor from which the Fe Kα photons originate, and study possible correlations between the Fe Kα line and continuum variability.

4.1. Variability features

4.1.1. The variability probability P_var

To assess whether the light curves are variable or not, we compute P_var = P(χ²) (e.g., Paolillo et al. 2004; Lanzuisi et al. 2014; Sánchez et al. 2017), which is the probability that a χ² lower than that observed could occur by chance, for a nonvariable source. Here chi-squared (χ²) is calculated as:

(4)

where N_obs is the number of observations, x_i is the flux measured at each observation, is the mean flux amongst all observations in the light curve, and σ_err, i is the flux error. If a source is not variable, we expect that χ² ∼ N_obs − 1. The typical threshold to distinguish variable from nonvariable light curves is P_var = 0.95 (e.g., Paolillo et al. 2004; Papadakis et al. 2008), which indicates a 95% chance that the source is intrinsically variable, or, alternatively, a 5% chance that the variability observed is due to Poisson noise. We calculate P_var for both the continuum and Fe Kα light curves, as listed in Table 2.

4.1.2. The excess variance

The normalized excess variance (e.g., Edelson et al. 1990; Nandra et al. 1997; Vaughan et al. 2003a; Paolillo et al. 2004; Papadakis et al. 2008) is a quantitative measure of the variability amplitude of a light curve, defined as

(5)

Effectively, is the intrinsic variance of the light curve, normalized by the mean flux, producing a dimensionless quantifier that can be easily compared between objects of different brightness or light curves from different energy bands. The last term in Eq. (5) denotes the contribution of the observational noise to the total variance, which is subtracted in order to find the intrinsic contribution.

Low intrinsic variances compared to the Poisson noise can sometimes lead to negative values of the estimate, since the uncertainty in the Poisson noise can be larger than the difference in Eq. (5). We performed Monte Carlo simulations to estimate more accurately the contribution of the observational noise (second term in Eq. (5)) and its asymmetric uncertainties in all our light curves, in order to quantify their variability robustly, even in cases where measured excess variance is negative.

For this purpose, we perform flux randomization of the light curves by adding a Gaussian deviate to each light curve point, with σ equal to the error on the flux of each point. This procedure adds variance to the light curves, on average by an amount equal to the observational noise, as shown in Appendix H. In this way, each simulation has the intrinsic variance of the light curve and twice the observational noise. Subtracting the variance of the original light curve from the flux-randomized light curve produces one estimate of the observational noise. Repeating this process 1000 times for each light curve, we obtained the median variance produced by the observational noise and its 16% and 84% bounds. We compared the median and bounds of the resulting excess variances (i.e., total variance -noise estimate) to the excess variance and error formula expressed in equations 6, 8, 9, and 11 of Vaughan et al. (2003a), obtaining consistent results for most cases.

We consider light curves to be significantly variable if the lower 16% bound of the excess variance distribution is positive. We caution that this limit may result in a few nonvariable sources being misclassified as variable, but we adopt it nonetheless to improve the completeness of the variable sample at the cost of reducing its purity.

Ultimately, we want to compare the variability amplitudes of the continuum and Fe Kα line light curves, distinguishing between cases where the Fe Kα line variability is not significant because it is small compared to that of the continuum versus cases where the variability amplitudes are similar but the lower count rate in the Fe Kα line renders its variability insignificant. For this, we employ a similar set of simulations, with one realization of the noise for the continuum and another for the Fe Kα line light curve, to compile the distribution of the ratio between the excess variances of both. Table 2 shows the 50% percentile and 1σ uncertainties of the continuum (Col. 10) and Fe Kα (Col. 11) distributions of the excess variances of each light curve, as well as their ratio (Col. 12). These data are also plotted in Fig. 4 (panels a, b, and c). To illustrate how the ratio changes as a function of reprocessor size, the ratios corresponding to four simulated reprocessors with different diameters (d_r = 0,2,10,100 ld) and the same power-spectral bend timescale of 10 days are shown in Fig. 4c.

Fig. 4. Variability properties of the sample. Sources are ordered based on the ratio between the excess variances of the continuum and the Fe Kα line, . All four plots adopt the same color code as in Fig. 4a. Open black circles denote AGN without a special characteristic, magenta empty circles are radio-loud AGN, black filled circles are CTȦGN, green empty circles are dual AGN, and filled cyan circles are changing-look AGN. (a): normalized excess variance of the continuum light curves. (b): normalized excess variance of the Fe Kα line light curves. (c): the ratio between the normalized excess variance of the Fe Kα line and continuum light curves, where the four horizontal lines denote ratio for simulated reprocessors with diameters dr = 0,2,10,100 ld, assuming identical 10 days power-spectral bend timescales. (d): ratio between the normalized excess variance of the simulated and real continuum light curves, with unity denoted by the red line. The light red shaded region denotes up to a factor 2 error in our assumed normalization of the simulated power-spectral model. Magenta empty circles are radio-loud sources, black filled circles are Compton-thick sources, green empty circles represents dual AGN, and cyan empty circles are changing-look AGN. Mrk 3 is both Compton-thick and radio-loud, while NGC 1365 is both radio-loud and changing-look. Plotted values represent the median of each parameter while the error bars correspond to the 16% and 84% bounds of the normalized excess variance distributions.

Significant variability, denoted by a positive lower bound on the , is detected in the continuum of 37 out of 38 AGN; the only exception is 4C+29.30. Significant Fe Kα line variability is detected in 18 AGN: 4C+74.26, Cen A, Circinus Galaxy, IC4329A, MCG-6-30-15, MR2251-178, MRK 1040, MRK 1210, MRK 3, NGC 1068, NGC 1275, NGC 1365, NGC 2992, NGC 3783, NGC 4151, NGC 4388, NGC 6300 and NGC 7582. For the rest of the AGN, the lower bound on the excess variance is negative, implying that the Fe Kα line variability is consistent with observational noise, within uncertainties. The upper bound on the Fe Kα line variability in those AGN is still of interest, depending on how this bound compares to the variance of the continuum.

If the Fe Kα line flux tracks the fluctuations of the continuum flux, then we expect their variances to be related. The ratio of the variances should be similar to unity if the reflector is small compared to the timescale of the fluctuations, and smaller than 1 if the reflector is large. Therefore, upper bounds smaller than 1 in the ratio column of Table 2 can allow us to place a lower limit on the size of the reflector. This is true in the first two objects in Table 2. Even though their Fe Kα line variance is not significant, its upper limit is so far below the variance of the continuum that a lower limit can be placed on the size of the reflector, because if the reflector was any smaller, the Fe Kα line variability would be detectable. Conversely, lower bounds of the variance ratio larger than 0 can put upper limits on the reflector size. One or both of these limits are therefore measurable in many objects in the table. We describe this analysis in detail in Sect. 4.2.

The normalized excess variance of the continuum is generally well constrained in the sources listed in Table 2. The median values of the distributions range from 0.007 for 3C 445 to 3.59 for NGC 2992, as can be seen in Fig. 4a. This large difference can raise doubts about a common origin of the continuum variability in all sources. We know however, from dedicated monitoring campaigns of radio-quiet AGN, that the X-ray power spectra of different AGN is remarkably uniform in shape, but that the frequency ν_b of the only break in the power spectrum scales inversely with black hole mass and directly with accretion rate (McHardy et al. 2006). Since the excess variance equals the integral of the power spectrum over the timescales covered by the light curves, we can expect to measure significantly different values for objects of different SMBH mass, even if the light curves have similar lengths.

To test how unusual some of these excess variances are, we compute the expected variance for the respective SMBH masses and accretion rates. Since the variability is a stochastic process, and the underlying power-spectral shape is only realized on average, simulations provide a good way of estimating the possible range of variance measurements, for given SMBH parameters and light curve sampling pattern. Thus for each source we generate simulated light curves of a red-noise process following the method of Timmer & Koenig (1995). The underlying power-spectral shapes are bending power-laws, with a low-frequency slope α_L = −1 bending to a steeper slope at higher frequencies α_H at a characteristic timescale or break frequency ν_b:

(6)

where the normalization is A = 10⁻² s with frequencies measured in Hz. This means that at low frequencies, ν ≪ ν_b, the combination νP(ν) tends to this constant value of A. We adopt this as the ‘standard power-spectral density (PSD) model’. The variance scales linearly with the normalization A, such that by shifting A, the mean expected variance and its scatter shift by the same factor. Differences in A by a factor of 2 (higher or lower) are consistent with the majority of the monitored sample.

This shape has been shown to represent well the power spectrum of AGN X-ray light curves. When available, we use the high-frequency slopes and break frequencies measured for each object in Summons (2007) based on the long-term monitoring campaigns performed with the Rossi X-Ray Timing Explorer (RXTE) observatory. For the rest, we estimate the break frequency using the expression from González-Martín & Vaughan (2012):

(7)

where T_b is the power-spectral bend timescale (T_b = 1/ν_b), M_SMBH is the SMBH mass in 10⁶ M_⊙ units, and L_bol is the bolometric luminosity in 10⁴⁴ erg s⁻¹ units. The adopted bend frequencies ν_b are tabulated in Table 1. The high-frequency slopes α_H in Summons (2007; see also Vaughan et al. 2003a; Markowitz et al. 2003; McHardy et al. 2004, 2005; Uttley & McHardy 2005; Summons et al. 2007) range from 1.5 to 3.5. For sources in which α_H has not yet been measured, we adopt an intermediate value of α_H = 2. The power-spectral model used herein only considers one break timescale, which is the only one detectable with short-term light curves such as those used by González-Martín & Vaughan (2012). It is possible that the power spectra for our AGN sample feature a second break at lower frequencies, flattening further to a value of 0. Dedicated monitoring campaigns on a few sources such as MCG–6-30-15 (McHardy et al. 2004), NGC 4051 (McHardy et al. 2005) and NGC 3783 (Summons et al. 2007) do cover timescales similar to those covered in this work, and notably have not shown any second, lower-frequency break, so this consideration might not be relevant to the sources studied here. If there were a second break, we would observe lower variances than those predicted by the model, especially for sources with a higher break frequency, for which the second break could be at higher frequencies as well.

We generated 100 realizations of the continuum light curves following this underlying power-spectral model and the expected or measured bend timescales and high-frequency slopes (α_H) for each object. The simulations are run for at least 10 000 days or 3 times the length of the light curve, whichever is longest, and generated with a time step of 0.01 days or 1/100 times the bend timescale, whichever is smaller. Each simulation is then sampled in an identical manner to the corresponding real light curve, its variance is computed and the mean and root-mean-squared scatter of the resulting variances are recorded. This mean and scatter represent our expectation for the excess variance measured in the continuum light curves.

The ratio between the real and expected excess variances, σ_c/σ_c, sim, is plotted in Fig. 4d, where the error bars represent only the scatter of expected variances produced by the stochasticity of the intrinsic variability; no observational noise is added. Naively, we expect sources to cluster around σ_c/σ_c, sim ∼ 1 (solid red line); allowing for up to a factor 2 difference in the normalization of the power-spectral model (denoted by the light red shaded region), we find that a majority of the measured variances are consistent with this expectation. However, some sources remain significantly above and below this; we distinguish radio-loud sources (magenta empty circles), Compton-thick AGN (black filled circles), dual AGN (green empty circles), and changing look AGN (cyan filled circles), which tend to be outliers. The radio-loudness (RL) of the sample is computed using the 20 cm and 14–150 keV fluxes from Ricci et al. (2017a) as RL = log(f_20 cm/f_{14−150 keV}) (e.g., Terashima & Wilson 2003; Panessa & Giroletti 2013; Panessa et al. 2015), listed in column 6 of Table 1. We adopt a separation between radio-loud and quiet sources at RL = −4 following Panessa et al. (2015) and Ricci et al. (2017a). It is possible that the continuum (and Fe Kα) variability of radio-loud sources could be affected by beamed X-ray emission associated with the powerful jets, leading to stronger or more rapid variations (e.g., Ulrich et al. 1997; Chatterjee et al. 2008; Weaver et al. 2020); this could explain why 6/18 radio-loud sources in our sample have σ_c/σ_c, sim> 1. On the other hand, Compton-thick AGN should have reflection-dominated continua with little variation, and thus our simple X-ray spectral fitting approach would measure the flux of this relatively static component. Thus, it is not surprising that the three Compton-thick sources of our sample all have σ_c/σ_c, sim < 1. Appendix C contains notes on individual sources that do not seem to conform to the standard PSD model.

4.2. Estimating the size of the X-ray reflector

If we assume that the Fe Kα line emission is reprocessed from the same X-ray continuum we observe, we can expect the Fe Kα line flux to track the continuum fluctuations. The light curves of both, primary and reflected components, can still differ by light travel time effects, as the reflected light will travel on different and longer paths to the observer. The Fe Kα line light curve can therefore be delayed with respect to the continuum and can also be smoothed out, as variations on timescales shorter than the light crossing time of the reflector are damped.

The majority of the light curves obtained in this work are too sparsely sampled, however, to detect directly a delay or lag between the continuum and Fe Kα line fluctuations, as exemplified in Figs. D.1 and D.2. The potential reduction of the Fe Kα line variability amplitude with respect to the continuum, however, can still shed light on the size of the reflector, as larger reflectors will suppress a larger fraction of the intrinsic variance. We quantify the reduction in the observed variability amplitude by calculating the excess variance of the continuum () and Fe Kα line () light curves and taking their ratio as . The longer the light crossing time of the reflector, the smaller the expected value of var_r should be.

To demonstrate the effect of reflection on the light curves, we show in Fig. 5 a typical power spectrum of continuum X-ray variations, which are damped by a reprocessor at a distance of 50 light days; for simplicity we assume a simple thin spherical shell reprocessor, but acknowledge that more realistic distributions could be thicker, clumpier, and have a toroidal or ionization cone-like structure seen at a specific orientation with respect to the line-of-sight, all of which can impact delay times. The corresponding continuum and reflected Fe Kα light curves are shown in Fig. 6. The reflected light curve is differentially delayed (by design), which results in a clear drop in variability amplitude at high frequencies (or short timescales). If we observe these simulated light curves in a manner comparable to the real data, the sparse sampling would result in two effects: (1) the sampled light curves could appear uncorrelated, preventing a measurement of any lag, and (2) the reduction in variance will be less deterministic, as the reflected light curve can be sampled at times when it is more variable by chance. To illustrate this, Fig. 6 shows continuum (black dots) and Fe Kα line (cyan stars) flux measurements for a random selection of observing epochs. The sampled light curves are not obviously correlated, although the parent light curves are. We note that observational noise has not been added to these simulations; all the scatter is produced by the intrinsic fluctuations of the source.

Fig. 5. Power spectra of simulated X-ray light curves. The underlying power-spectral shape of the continuum was chosen as a bending powerlaw model with a bend timescale of 10 days (black) and a given realization of this model produces the power spectrum plotted in red. Reflecting the light curves over a spherical reflector of radius R = 50 light days suppresses the high frequency power and results in the model and realization plotted in cyan and blue, respectively.

Fig. 6. Example light curve realizations of the continuum (red) and reflected light (blue), which correspond to the power-spectral models plotted in Fig. 5. The reflected light curve reduces variability amplitude and is delayed by 50 days. Increasing the size of the reflector further reduces the amplitude of the reflected light curve and shifts it further forward in time. A random selection of observing epochs are plotted, representing measurements of the continuum (black dots) and Fe Kα line (cyan stars) fluxes. The sampled light curves are not obviously correlated, although the parent light curves are.

The variance of the continuum light curve in the above example, after sampling, is still larger than the variance of the reflected light curve, although the ratio of the variances depends on the sampling pattern. Intuitively, the degree of apparent randomness in the reflected light curve, and the scatter between the continuum and reflected variance ratios, increases as the time lag increases compared to the characteristic timescale of fluctuations, as quantified by the bend timescale in the power spectrum. The delay between continuum and reflected light curves, together with the sparse sampling, can explain part or all of the scatter seen in the flux-flux plots in Fig. D.3.

With this in mind, we examine the observed values of var_r in Fig. 4c. Somewhat surprisingly, we find that 11% (4/38) of sources lie > 1-σ above d_r = 1 ld (black continuous line), although all of these appear consistent with unity within 2-σ uncertainties. Similarly, 21% (8/38) of sources lie above the d_{r = 10 ld} line (∼0.009 pc, close to the lower grouping of sources in Fig. 3), but at 2-σ, only Cen A remains confidently above this value. In total, more than 50% have values consistent with unity, implying little damping of the continuum by the reflector (i.e., angle-averaged light-crossing timescales to the reflector are comparable to the continuum variability timescales), or alternatively that we are not observing the true continuum fluctuations (as is likely the case for the three Compton-thick AGN, which all have ratios consistent with unity). Notably, all the sources which show stronger than expected continuum excess variances in Fig. 4d consistently have var_r values below unity.

Although it is straightforward to estimate var_r for our fiducial reflector model as a function of the break frequency of the continuum light curve power spectrum and the light crossing time of the reflector, we expect a large scatter in measured values due to the small number of data points and the stochastic nature of the light curves. We therefore adopt a Monte Carlo approach, as described below, to place meaningful constraints on the size of the reflector, incorporating both the stochasticity of the intrinsic variations and the observational noise.

We employ the same setup as for simulating continuum light curves from a bending powerlaw underlying power spectrum, described in Sect. 4.1.2. To simulate reflected light curves, we multiply both the real and imaginary parts of the Fourier transform of the simulated continuum light curve by a sinc(ντ) function, which corresponds to a convolution with a top hat function of the light curve in time space. We then shifted the inverse Fourier-transformed curve forward in time by τ/2 days. This corresponds to reflection of the continuum by a spherical shell of radius R = cτ/2, leading to an immediate response from the front end of the reflector, followed by reflection from the rest of the shell until the light from the back end finally reaches the observer 2R/c = τ days after the start of the initial front-end response (e.g., Sect. 5.1 in Arévalo et al. 2009). This particular response function is chosen for simplicity, and reproduces the main characteristics of a reflected light curve, that is suppressed variability amplitude and average time delay. This is sufficient to estimate the size of the reflector, but not its geometry. Replacing the idealized thin shell reprocessor above (i.e., the sinc² multiplicative function on the power spectrum) with the average of three reprocessors of similar size (e.g., τ, τ × 1.5, τ/1.5) can approximate the response of a finite width reprocessor. We implemented this setup as well, noting that it made little difference on the resulting sizes. The sizes reported below correspond to the latter finite-width reprocessor case.

Transfer functions for different reprocessors, such as a thick shells with different radial matter distributions, and flat or flared disks are presented in Sect. 5.1 of Arévalo et al. (2009) for a given average delay τ between the direct and reprocessed light curves. All these functions have in common a flat response for frequencies below 1τ and consistently drop above it. Since we can only attempt to measure the decrease in variance of the Fe Kα line light curve compared to that of the continuum, our interest is to compare the integrals of the original and filtered power spectra, not their precise shape. Therefore, any response function that generally complies with the shape described above will suffice. We note however that reprocessors that are strongly asymmetric, such as clouds lying along only one axis or a compact region far from the nucleus are only partially captured by our approach. These reprocessors can have a smoothing effect proportional to their diameters, but a delay proportional to the average distance to the source. In these cases the reduction in Fe Kα line variance would give an estimate of the diameter of such reprocessors, but the delays, and therefore the allowed scatter between realizations of the variance ratio, would be underestimated.

For each source, we consider its spectral parameters and choose a sampling timescale dt that is related to its bending timescale as dt = T_b/100 or dt = 0.1, whichever was shortest, such that the sampling resolution is at least a 100th of the T_b scale. We resample the simulated continuum and Fe Kα light curves according to the observed epochs. Finally, we compute the ratio of the variances (). As expected, this ratio decreases as the light crossing time of the reflector increases. We run sets of 100 simulations and record the median and median-absolute-deviation (MAD) of var_r, sim for a range of values of the light crossing time. To determine the τ values which correspond to our estimated reprocessor size and its upper and lower limits, we adopt the following statistic:

(8)

where var_r, real is the median value of the var_r distribution and var_r, sim is the median value of the simulated ratios.

The τ values that return X = 0, 1, −1 correspond to the reprocessor size and its ±1σ uncertainties in light days, respectively. In order to constrain the size of the reprocessor, we first explore a range of τ values between 0.01–10 000 days, and then numerically look for the roots X = 0, 1, −1 using Newton’s method.

Table 2 summarizes our estimates for the size of the reprocessor (R_rep) derived from the values of τ at X = 0, 1, −1, as well as the ν_b value used in the simulations. The R_rep values are also shown in summary Fig. 12 (black squares), for comparison against the other reprocessor size estimates or limits derived from the spectral and imaging (see Sect. 5) analysis. The simulations provide estimates of the reprocessor size for 24 sources in our sample.

We note that this approach can only provide limits on the size of the reflector if a limit of the measured ratio is between 0 and 1. If the upper and lower limits of the measured ratio, considering errors from observational noise only, already cover the whole range of possible ratio values, then any reprocessor size is consistent with the data and no limits can be placed. This is the case for 14 of our sources, as can be seen in Fig. 4c. Conversely, if one or both limits on the measured ratio fall between 0 and 1, it is in principle possible to place limits on the size of the reflector. The additional scatter in possible ratio values produced by the red-noise nature and sparse sampling of the light curves can, however, extend the error bars beyond the 0–1 range and prevent a limit to be placed. This happened in the case of Pictor A. In addition, for sources where the lower limit on the measured ratio was above 1, an upper limit on the reprocessor size could still be placed in the case of NGC 1068 since the additional scatter in the ratio allowed a value of τ for which X = 1. For sources MCG-6-30-15, MR2251-178 and 4C+74.26, the smallest value of τ explored still produced a value of X > 1, so in these cases the upper limits on τ were set to the smallest value explored.

4.3. Correlations between the Fe Kα line and 2–10 keV continuum

A strong, positive correlation between Fe Kα line and 2–10 keV continuum fluxes, which is much simpler to assess compared to more complicated lag analyses, should indicate that the reflector lies in close proximity to the source of the X-ray continuum emission. Thus we explore potential correlations between the observed 2–10 keV continuum and Fe Kα line light curves for sources in our sample. As mentioned above, any lag in the Fe Kα line light curves, which is expected from reprocessing travel time delays, can reduce the apparent correlation between the observed light curves. The degree of loss of correlation is a function of the ratio between the light crossing time of the reprocessor and the characteristic timescales of fluctuations, as well as the geometry of the reprocessor.

For illustration, Fig. 7 shows a comparison between the fluxes from the sampled observing epochs of the simulated light curves in Fig. 6 (i.e., the dots and stars). In that simulation, the power-spectral bend timescale is 10 days, while the reprocessor diameter is 100 light days. The corresponding flux-flux points, normalized to their respective means, are plotted in blue open squares. The correlation is weak, with a large scatter, and the best-fitting linear regression yields a slope of 0.17 ± 0.19; the expected slope for a perfect correlation and equal variability amplitudes is 1. However, if the same continuum light curve is reprocessed by a smaller structure, the correlation improves and the slope approaches a value of 1. In the same figure we show the result of sampling the light curves on the same observing epochs but using a reprocessor diameter of 10 days (red stars) and 2 days (black dots). The best-fitting slopes for these smaller reprocessors are 0.46 ± 0.18 and 0.68 ± 0.12, respectively. We note that the loss of correlation in these sampled light curves is due entirely to the effect of the reflection, both smoothing and delaying the Fe Kα light curve, since no observational noise is added to these simulations. A key point here is that, at least for a symmetric reflector observed with a relatively sparse cadence, its light travel distance only needs to be a ∼10 times larger than the continuum source variability timescale for any correlation to be almost completely washed out.

Fig. 7. Comparison of simulated Fe Kα and 2–10 keV continuum mean-normalized fluxes. The continuum simulations use a power-spectral bend timescale of 10 days, while three different reprocessor diameters (dr) are chosen: 2 light days (black dots), 10 light days (red stars) and 100 light days (blue open squares). The simulated light curves are resampled using the observing epochs shown in Fig. 6. The different sets highlight how varying the reprocessor size can damp potential flux-flux correlations. Evidently, the larger the reprocessor, the larger the scatter in the flux-flux relation, and the flatter the best fitting linear regression (solid lines with corresponding colors).

To assess if the Fe Kα line and continuum light curves of our sources are correlated, we fit a line using the python package LINMIX⁴. This package is based on the code of Kelly (2007), which describes a Bayesian method to perform linear regressions with measurement errors in both variables. The code runs between 5000 and 100 000 steps of a Markov chain Monte Carlo to produce samples from the posterior distribution of the model parameters, given the data. We choose this code since it does not assume a particular distribution for the errors and incorporates nondetections, both of which are fundamental for this study since in several cases the fluxes of the Fe Kα line are upper limits. The flux errors obtained in the spectral analysis are mildly asymmetric, as the data often straddle the division between Poisson and Gaussian distributions. We symmetrize them adopting the average between upper and lower errors, in order to use them as input to LINMIX, but this does not strongly impact the results for most sources; for sources with well-constrained fluxes, the positive and negative errors only differ by a few percent (≲5%), while for poorly constrained fluxes, we adopt conservative (larger) errors or upper limits.

The LINMIX package requires at least five measurements, which all sources in our sample satisfy by definition. The best-fit linear regression results in 29 AGN with well-constrained slopes, and nine AGN with unconstrained slopes. Appendix Fig. D.3 compares the Fe Kα line and continuum fluxes, with the best-fit slope shown in the label. In each plot, the thin red lines correspond to the samples of the posterior distribution and the red thick line is the average posterior. Table 2 reports the values of the well-constrained slopes (m), the maximum light curve timespan (t_lc), and the P_var parameters obtained for the Fe Kα line and 2–10 keV continuum light curves.

We find three sources with a F_Fe Kα–F_2−10 keV slope consistent with one (i.e., m + 1σ > 1.0 and m − 2σ > 0.5). The most outstanding case is Cen A, with a slope m = 1.14 ± 0.20. The top plots of Fig. 8 show the light curves of Cen A for the continuum and Fe Kα line. The Fe Kα line light curve is strongly variable, with P_var, Fe = 1 and , and appears to perfectly track the continuum, with variations of up to 1 dex in five days. This suggests that Cen A has a compact reflector, very close to the source of the X-ray continuum. Similarly strong constraints are seen for H1821+643 and IC 4329A.

Fig. 8. Light curves for the 2–10 keV continuum (left) and Fe Kα line (right) fluxes for the galaxy Centaurus A (upper panels), 1H0707-495 (middle panels) and NGC 1275 (lower panels). The green and black stars denote fluxes measured from the Chandra HEG first-order spectra and Chandra ACIS 3″–5″ annular spectra, respectively, while yellow stars are fluxes measured from XMM-Newton pn observations. The flux units are erg cm−2 s−1 and the error plotted are 1σ confidence.

On the other hand, six objects have slopes consistent with zero (i.e., m − 1σ ≈ 0 and m + 2σ < 0.5); all have σ_Fe consistent with ≈0.0 as well, although only one shows no clear Fe Kα line variability (P_var, Fe < 0.5). This latter source is 1H0707-495, with a slope of m = −0.28 ± 0.27; its light curves are shown in the middle panels of Fig. 8. While the continuum light curve is clearly variable (P_var, c = 1, σ_c = 0.195 ± 0.005), the Fe Kα line light curve has a P_var = 0.37 and σ_Fe = −0.18–0.05, indicating no discernible line flux variability over a ∼19 year timescale. Such behavior suggests that the light crossing size between the X-ray continuum source and the reprocessor is substantially larger than the typical continuum variability timescale; based on the simulations carried out in the previous section, the reprocessor is ≳3 light days away, which is still consistent with BLR clouds. However, in a recent investigation, Boller et al. (2021) found that the X-ray spectrum of 1H0707-495 is dominated by relativistic reflection, and that the absorber is probably ionized, leading to weak Fe Kα variability. The other five sources are NGC 3783, NGC 4051, NGC 4151, NGC 5548, and NGC 7469.

A majority of the remaining AGN show some degree of variability of the Fe Kα emission (P_var, Fe > 0.9), but span a wide range in terms of slope uncertainties; some have well-determined intermediate slopes, while others have completely unconstrained slopes. The lower panels of Fig. 8 show the light curve for NGC 1275, which has a low but uncertain slope of m = 0.07 ± 0.31, although its Fe Kα line light curve appears variable, with P_var = 1 and σ_Fe = 0.029–0.83.

The analyses made in Sects. 4.1.2–4.3 indicate that different AGN have distinct morphological and structural distributions of reflecting clouds. Some sources have a compact reprocessor, whereby the Fe Kα flux tracks the continuum variations on timescales comparable to each observation exposure. For other sources the line varies but the correlation between the continuum and the line is weak due to the damping effects of a large or distant reprocessor. And finally, some sources show little variability in the line and do not have correlated fluxes, suggesting that the reflector is sufficiently large or distant, compared to the variability timescales of the continuum source, to wash out any reaction from the line.

4.4. Potential influence of relativistically broadened Fe Kα

As noted in Table 1, nearly half of the objects in our sample have been argued to have relativistically broadened Fe Kα emission in the literature. Thus, an important consideration for interpreting the above results is to understand the influence that any potential relativistically broadened Fe Kα emission will have on the variability measurements of either the continuum or the narrow Fe Kα line. We begin with a general caveat regarding the veracity of relativistically broadened Fe Kα detections in the literature, which remain generally controversial for the majority of objects for which a detection has been claimed. This is in part due to a combination of limited photon statistics (e.g., Brenneman & Reynolds 2009 argue that spectra with ∼10⁶ counts are generally required to confirm relativistically blurred components in < 10 keV spectra) and lack of good quality, simultaneous spectral constraints spanning both ∼0.5–8 keV and ∼8–100 keV, which are essential to lock down intrinisic spectral slopes and reflection fractions (these are often degenerate even in high-count ∼0.5–8 keV spectra due to potential combinations of neutral and ionized absorption). As a consequence, the constraints on the contributions of such relativistically blurred Fe Kα components remain poorly constrained (e.g., reflection fractions ranging from ≲0.01 to ∼20 Brenneman & Reynolds 2009).

Importantly, the continuum as we model it will automatically absorb the bulk of any relativistically broad/blurred line flux and variability. Since both components are strongly correlated, with reported lags on the order of minutes to hours (de Marco et al. 2009; Kara et al. 2016), our observation by observation analysis should not be strongly affected, given the typical length of the observations. Another consideration is whether the relativistically broadened Fe Kα profile peaks near 6.4 keV, and thus if a portion of that component could be assigned to the narrow Fe Kα lines that we fit, and lead to an enhancement in correlations. Hu et al. (2019) performed a stacking analysis of 193 RQ and 97 RL AGN, arguing that average broad line components are detected in both subsamples, although the average broad-line Fe Kα components are subdominant compared to the narrow component by factors of ∼4 for RQ and ∼2 for RL AGN in the 6.2–6.6 keV regime, respectively. Similarly, Falocco et al. (2014) fit the stacked spectra from 263 X-ray unabsorbed AGN, finding that narrow (“unresolved”) Fe Kα emission accounts for ≳70% of the total 6.2–6.6 keV line equivalent width. Given that our simple continuum model already incorporates for some of the broad Fe Kα in the 6.2–6.6 keV regime, we assume that the narrow line as we measure it suffers from < 10–20% contamination at most. Further considering the arguments presented in Appendix G, the Fe Kα line variances that we report should be related to the square of the aforementioned fractional contributions, and hence strongly dominated by the narrow component in all cases.

4.5. Variability properties compared with AGN and host galaxy properties

To gain further insight into the different behaviors of the Fe Kα line in our sample, we compare the σ_Fe/σ_c values and F_Fe Kα–F_2−10 keV slopes (for 29 sources in the sample with firm σ_Fe/σ_c and slopes measurements) with several AGN properties (see Figs. E.1 and E.2 of Appendix E). The line-of-sight column density and radio-loudness are taken from Ricci et al. (2017a), who derived the X-ray properties of the sources through X-ray spectral fitting. The SMBH mass, Eddington ratio, and Seyfert type are taken from Koss et al. (2017) and BASS DR2 (Koss et al. 2022). The SMBH masses were estimated using broad-line measurements assuming virial motion, or stellar velocity dispersion, assuming the M_BH − σ_* relation. To test possible correlations, we compute the Spearman correlation coefficient ρ with its p-value and also perform a linear regression between the variability features and the properties mentioned above. The Spearman coefficient varies between −1 and 1, with ρ = 0 and ρ = ±1 implying no correlation and a perfect correlation, respectively, and the p-value represents the probability of that the values are uncorrelated.

Figure 9 compares the slope and the central SMBH mass, which is the only set of properties where we find a weak correlation (ρ = 0.36, p-value = 0.05, and linear regression slope m = 0.22 ± 0.13). If the continuum variations are produced by the corona, we expect somewhat longer bend timescales for larger SMBHs following Eq. (7), presumably with some scatter due to the effect of spin on the innermost stable circular orbit. To further test whether the weak correlation is real, we perform a bootstrap analysis to estimate the confidence interval of the Spearman coefficient, finding that ρ ranges between 0.01 and 0.65 at a 95% confidence level. The wide confidence interval suggests that our results are compatible with both a positive and nonexistent correlation between the slope and SMBH mass, which does not allow us to draw a definitive conclusion. In contrast, Fig. E.1a compares σ_Fe/σ_c with central SMBH mass, where the linear regression and the Spearman test suggest no clear correlation (ρ = 0.21, p-value = 0.29, and slope m = −0.12 ± 0.2), implying that the observed damping does not depend on the SMBH mass. The disagreement between these two indicators, combined with the low strength of the relationship with the slope, suggests that the weak correlation in the former may not be real.

Fig. 9. Comparison between the FFe Kα–F2−10 keV slopes and the SMBH mass (MBH).

Similarly, we search for correlations between both σ_Fe/σ_c (Figs. E.1b, E.1c E.1d, and E.1e) and F_Fe Kα–F_2−10 keV slopes (Figs. E.2a, E.2b, E.2c and E.2d) versus the Eddington ratio, the line-of-sight column density, radio loudness and AGN type, respectively. Table 3 reports the linear regression slope, ρ, and p-value for each property. We might expect potential relations with any one of these parameters. For instance, the N_H and AGN type should trace the overall reflector geometry, which should imprint itself on the Fe Kα line properties. Likewise, the Eddington ratio has been linked to variations in N_H, potentially sculpting the inner few 10s of pc via radiative feedback (e.g., Ricci et al. 2017b). Finally, whether an AGN is radio-loud or not may be associated with certain physical conditions related to black hole spin and magnetic threading of the accretion disk, which could extend to the broader local environment. Among all of these possibilities, however, we find no clear trends.

Summarizing, we find a weak correlation between the F_Fe Kα–F_2−10 keV slopes with the SMBH mass, and no correlations with the Eddington ratio, column density, radio-loudness, or AGN type. Therefore, the only property that might affect the reaction of the Fe Kα line to the continuum variations appears to be the SMBH mass.

5. Imaging analysis with Chandra

To complement the spectral results above, we analyze the Chandra images to look for possible spatially extended Fe Kα emission in our sample. Extended Fe Kα emission has been observed in a handful of nearby, typically Compton-thick, AGN, and we want to examine whether sources where the Fe Kα flux reacts quickly to continuum variations are point-like, and sources with no statistical line variability may be more spatially extended.

To quantify this, we investigate the radial profiles of the AGN in our sample, using both the azimuthally averaged profile in comparison to the nominal Chandra PSF and by comparing averaged quadrant profiles to look for strong asymmetries. We analyze a rest-frame continuum-subtracted Fe Kα-only image, as well as the high-energy continuum in the rest-frame 5–6 keV band. The Fe Kα-only image is created by taking the rest-frame 6.2–6.5 keV band, which should capture the vast majority of 6.405 keV photons⁵, and subtracting the average between the rest-frame 5.9–6.2 and 6.5–6.8 keV bands, which should remove the underlying 6.2–6.5 keV continuum given that the ACIS effective area smoothly and linearly changes between these two bands. This scheme is not perfect, as the faint wings of the Fe Kα spectral profile will extend into the 5.9–6.2 and 6.5–6.8 keV bands, but should nonetheless yield the approximate spatial distribution of Fe Kα photons, is relatively straightforward to implement and most importantly does not suffer from the PSF calibration uncertainties at large annuli (see below). The rest-frame 5–6 keV band is adopted as a proxy for the broader 2–10 keV, because it will not be as strongly affected by absorption as lower energies, it is free of strong emission lines, and Chandra’s sensitivity is still relatively high here.

The Chandra PSF is predominantly a function of the energy and off-axis angle, being compact and roughly symmetric on-axis but broadening substantially both radially and azimuthally beyond off-axis angles ≳2′, with complex structure due to mirror (mis)alignment and aberrations, shadowing by mirror support struts, and increased scattering as a function of energy⁶. The observed PSF shape further depends on the number of counts, with pileup potentially impacting the innermost pixels and at least a few counts/pixel needed to fully sample the PSF wings and complex structure. This combination makes it extremely difficult to calibrate the full shape of Chandra’s PSF using in-flight observations of single bright targets⁷. As such, the PSF is typically simulated with numerical ray-trace calculations based upon a fiducial mirror model developed using preflight measurements. These models, however strongly underestimate the flux in the PSF wings at energies higher than 2 keV by factors of ∼2–3 beyond > 3″ (as we show in Sect. 5.2). Given such limitations, we describe our approach to analyze the Chandra images below.

5.1. SAOTRACE and MARX simulated PSFs

The Chandra PSF can be simulated using either SAOTrace⁸ or MARX (Davis et al. 2012), both of which model the on-orbit performance of the various operationals modes of Chandra for input source shapes and spectra. SAOTrace uses a more detailed physical model of the mirror geometry, which can be important for off-axis sources or to study the wings of the PSF out to several arcseconds, but is computationally very expensive for bright sources. MARX has as its default a slightly simplified description of the mirror, which is much faster to run and produces very similar overall results (differing by only a few percent in radial profiles). For these reasons, we use MARX (v5.4.0) for the simulations below. MARX takes as input the source spectrum, as well as parameters such as the sky position of the source AGN, exposure time of the observation, and grating type, to generate a simulated event list for the source.

To understand ability of MARX (and by extension SAOTrace) to reproduce faithfully the wings of the PSF, we compare it to a bright, point-like X-ray source which has been extensively observed by Chandra. For this purpose, we analyze the archive ACIS-S HETG data for the X-ray binary Hercules X-1 (HERX1, hereafter), following the procedures outlined in Sect. 2.2. Importantly, HERX1 lies well above the Galactic Plane (N_H, Gal = 1.5 × 10²⁰ cm⁻²; HI4PI Collaboration 2016), and thus should not show strong contamination at large angles due to dust scattering (e.g., Smith & Dwek 1998), which could complicate interpretation of extended features at soft energies. Specifically, we merge the images and grating spectra for ObsIDs 2749, 3821, 3822, 4585, 6149, and 6150, resulting in total exposure time of 166.6 ks. We model the first-order HEG and MEG spectra from each ObsID with a simple unabsorbed powerlaw model [Γ = −0.125, F_2−10 keV = (2.5)  ×  10⁻¹⁰ erg cm⁻² s⁻¹], which we provide as input to MARX. HERX1 is relatively bright, with a comparable number of counts to the brightest AGN in our sample, and therefore has a well-defined radial profile out to ∼200″ in the combined image. This high count rate, however, means the central few pixels are affected by pileup. For this reason, we restrict comparisons to outer radii, renormalizing the PSFs at 25; at this radius, the radial profile is always found to be declining, with counts per frame is < 0.01, well below the regime where pileup begins to occur, while the radial profile of grade 0 and 6 events remains flat and close to expected values⁹. We ran MARX using an aspect blur of 02¹⁰, and leaving pileup turned off, such that the MARX simulated images and radial profile yield intrinsic images and profiles.

After generating empirical and simulated images of HERX1, we extract radial profiles for both, which we compare in Fig. 10, in two energy ranges characteristic of the continuum and Fe Kα line. Both radial profiles show rough powerlaw declines. After normalizing the profiles at 25, we see that the real data fall below the simulation at small radii (as expected due to pileup), while at large radii the real data are systematically higher.

Fig. 10. Comparison between real and simulated radial profiles for the X-ray binary HERX1. Panels a,b: radial profiles for the 6.2–6.5 keV and 5–6 keV ranges, respectively. The profiles are renormalized at 25, to highlight the regime interior to which pileup effects may occur in the real data. The errors of profiles denote 99% confidence.

We expect that pileup should be a roughly symmetric effect, and hence we test the symmetry of the PSF for HERX1 on all scales by comparing the radial profiles in Fig. 11 for four distinct quadrants (q1, q2, q3, and q4) as shown in the NGC 1068 coordinate plane of Fig. F.4. Importantly, the quadrant profiles are consistent between 1″–80″ to within errors. We see significant differences in the 05 bin, which may be due to the known PSF artifact/spur at ∼06–08 caused by a misalignment in the HRMA mirrors¹¹. We also see differences beyond 80″, which are the result of the dispersed order spectra.

Fig. 11. Comparison of 6.2–6.5 keV radial profiles per quadrant (q1, q2, q3, and q4, as defined in Appendix Fig. F.4) for HERX1. The profiles are statistically consistent, indicating that the PSF appears to be approximately symmetric out to at least ∼100″ in the energy range which includes Fe Kα emission for our sample. The errors plotted are 99% confidence.

A fundamental point here is that if we model sources based on simulated PSFs, which is the only method provided by the Chandra X-ray Center, we would misinterpret point-like emission as extended emission. Thus, for the rest of our analysis, we adopt the radial profile of HERX1 as our empirical Chandra PSF, which we compare to our AGN sample.

5.2. Empirical PSF comparison

We first compare the 6.2–6.5 keV radial profiles of HERX1 to the 6.2–6.5 keV radial profiles for all the sources in our sample in Fig. F.1, to estimate whether any extended (Fe Kα+continuum) emission is observed. We restrict our radial profile comparisons with HERX1 to beyond 25; this only has a small impact in our analysis since the observations of most AGN in the sample are also affected by pileup in the central ∼2″, and thus are not reliable anyway. To make a proper comparison, we normalized each profile to the flux at 25, except for Cen A and 2MASXJ23444 where we normalized at 4″ to avoid their more severe pileup. Most of the sources have a profile consistent with that of HERX1, except for Cygnus A, H1821+643 and NGC 1275. For these latter sources, a visual inspection suggests contamination by jets, diffuse emission from the galaxy cluster which the galaxy are centrally embedded, or starburst emission.

Next, we compare the Fe Kα-only radial profiles (i.e., continuum-substracted) of our sample to that of HERX1 in Fig. F.2¹². The azimuthally averaged radial profiles of all the sources are formally consistent with that of HERX1 at 3-σ confidence, although some sources show ≈2-σ excesses beyond 25 in one or more bins compared with HERX1 (e.g., MRK 1040, NGC 4151, NGC 1068). The profiles of Cygnus A, H1821+643 and NGC 1275 have no counts after 25, since we could not estimate the underlying Fe Kα continuum. We exclude these sources from the remaining analysis since we cannot draw conclusions about the spatial extension of those sources.

To assess whether the Fe Kα emission is asymmetric in some cases, we compare the radial profile of each quadrant of the images. Figure F.3 shows the profiles divided by 25% of the renormalized HERX1 radial profile used in Fig. F.2, and the error bars are plotted at a 3-σ level. In several cases, when the Fe Kα emission is not very bright, we only obtain an upper limit for the flux at a given radius. The plots indicate that in most sources, the radial profiles are consistent between the quadrants, suggesting that the Fe Kα emission is symmetric. However, we clearly see extended emission in the Circinus Galaxy and NGC 1068. The Circinus Galaxy shows asymmetry around ∼1–5″, particularly in the E-W direction. As such, we rotate the quadrants by 45° to obtain a more robust measurement of the extent (see bottom panel of Fig. F.4). We find that quadrant 2 has more flux than the others, followed by quadrant 4. In the case of NGC 1068, quadrants 2 and 4 also have more flux than quadrants 1 and 3 extending out to ∼10″. This can be clearly appreciated in the Fe Kα-only images of Fig. F.4.

We estimate the physical extent of the Fe Kα emission for both sources by converting the angular size to physical distance. For the Circinus Galaxy (z = 0.0014), the Fe Kα emission extends up to 95 ± 15 pc, and for NGC 1068 (z = 0.0037) up to 795 ± 176 pc. The errors correspond to the 25% of the radii bin width. The rest of the sources do not show extended emission outside 25 (4″ in the cases of Cen A and 2MASXJ23444), which we can use as an upper limit for the Fe Kα emission. A summary of the constraints calculated from the imaging analysis are listed in Table 2.

We finally estimate the fractional contribution of the extended Fe Kα emission in Circinus and NGC 1068, by subtracting the Fe Kα counts in quadrants 1 and 3 (which show little or no extended emission) from quadrants 2 and 4 (which show extended emission), and compare this to the overall Fe Kα counts. We find that the extended emission contributes < 8.0 ± 0.9% and < 16.7 ± 7.7% of the total Fe Kα flux in Circinus Galaxy and NGC 1068, respectively. These values are strong upper limits to the real ones, since both images are affected by pileup. Although the value for NGC 1068 is higher, for both AGN, the extended emission represents a small portion of the total Fe Kα emission.

It is important to understand whether the extended emission seen in Circinus Galaxy and NGC 1068 is a result of some observational bias, rather than an intrinsic difference. In particular, the Circinus Galaxy and NGC 1068 have ∼1 Ms and ∼0.5 Ms of imaging data, respectively, which rank among the highest in the sample, compared to the majority of sources which only have ∼100 ks of data. Comparing the total number of Fe Kα counts among sources, these objects have factors of > 60 and > 5 more counts than the rest of the sample. They are also among the closest sources in the sample, such that only 15% and 30% of the sample would be able to resolve extended Fe Kα emission on similar angular scales. Finally, Circinus Galaxy and NGC 1068 are both Compton-thick AGN (two of the three in the sample), meaning we likely observe far less direct emission from the BLR, accretion disk and corona, and hence achieve better contrast with which to observe extended features. All of these factors likely play a part in producing radial profiles which exhibit extended emission for these two targets, and thus we do not rule out that higher quality data could find faint Fe Kα extension in other sample sources.

6. Constraints on the origin of the Fe Kα line

Comparing the spectral, timing and imaging results for a representative sample of 38 AGN with at least five Chandra observations, complemented by existing XMM-Newton observations, we find multiple Fe Kα location estimates for 37 of them. Figure 12 summarizes the different Fe Kα location constraints for our sample. Although there is a large dispersion among the constraints (with imaging being the least constraining due to fundamental instrument limitations), these three independent measurements show strong consistency for the majority of the sample, indicating that the bulk of the Fe Kα emission typically originates inside the dust sublimation radius, and thus is presumably associated with the BLR or outer accretion disk. As we discuss below, only a handful of AGN show clear discrepancies between the indicators (Circinus Galaxy, NGC 1068, NGC 4151, and NGC 5548) and only NGC 2992 exhibits indications that the bulk of Fe Kα arises from beyond the dust sublimation radius. Notably, several of these are among the best studied AGN in the local universe.

Fig. 12. Summary of the reprocessor size estimations computed through the spectral analysis (red stars/limits), the timing analysis (black squares/limits), and the imaging analysis (green lines/limits). The dust sublimation radius Rsub is plotted in blue as a reference. Black squares without errorbars are sources where the simulations could not provide upper or lower limits for the size of the reflector. The green shaded areas in Circinus and NGC 1068 represent the extent of the extended emission that we observe in the imaging analysis.

The spectral analysis of Sect. 3 implies that the bulk of the Fe Kα emission for 21 out of 24 sources arises from inside the dust sublimation radius, and that for eight out of 11 sources arises from outside the radius of the optical BLR, suggesting that most of the time the Fe Kα emitting clouds are associated with the BLR, and in a few cases with the outer accretion disk. This is consistent with several previous studies, which all rely on a subset of the HETG spectra used here. For instance, Shu et al. (2010, 2011) and Jiang et al. (2011) studied the Fe Kα line core for several type 1 and type 2 AGN, finding that a fraction of the Fe Kα line originates at factors of ∼0.7–2 and ∼0.7–11 times the radius of the optical BLR, for type 1 and type 2 AGN, respectively. Minezaki & Matsushita (2015) similarly conclude that the Fe Kα line lies between that of the broad Balmer emission lines and the dust reverberation radius for a similar sample of AGN. Alternatively, Gandhi et al. (2015) found that the Fe Kα radii are similar or smaller than the radius of the optical BLR, and concluded that the dust sublimation radius is an outer envelope to the Fe Kα bulk for type 1 AGN. Notably, all of the above studies, including ours, are based on the first-order HEG spectra, where a small fraction of the FWHM estimates are of the same order as the nominal spectral resolution. Liu (2016) demonstrated that the best-fit line widths from the second- and third-order spectra are systematically lower then the first-order spectral fits by ≈30%, although, in general, they remain consistent within the errors. In the case of NGC 1275, Hitomi Collaboration (2018) measured a line width of 500–1600 km s⁻¹ (90% confidence) with Hitomi, which a factor of 2–6 lower than the value we found, although consistent within errors. Even with such additional systematic errors, the values remain firmly inside of the dust sublimation radius, under the assumption of virial motion.

The timing analysis of Sect. 4 provides complementary information for the sample. The simple model considered in Sect. 4.1.2, that invokes a universal power-spectral shape with a single break time-scale scaling with mass and accretion rate, accounts for much but not all of the differences in continuum variance between the sources. This can be observed in Fig. 4d, which shows that for some sources the continuum is more quiescent than expected (σ_c/σ_c, sim ≪ 1), while in others it is significantly more variable (σ_c/σ_c, sim ≫ 1). Notable deviations are seen in the Compton-thick sources, which show lower variability than the simple model predicts, consistent with the fact that the direct coronal emission is fully obscured in the observed band and variations might arise from slight changes in the very heavy obscuration (e.g., Marinucci et al. 2016). On the other hand, the one changing-look AGN in the sample (NGC 1365) varies more than expected probably due to well-documented large changes in obscuration which appear to act in addition to intrinsic variations of the continuum, which remains directly visible (e.g., Rivers et al. 2015). Additionally, radio-loud sources show more scatter in σ_c/σ_c, sim than radio-quiet ones, possibly due to different levels of jet contribution, none of which is included in the model.

The simulations described in Sect. 4.2 provide estimates of the X-ray reprocessor size for 24 sources. Ten have upper limits on the reprocessor size of < 1 pc, suggesting that the whole Fe Kα flux arises in regions smaller than R_sub, the fiducial inner wall of the torus, and in a few cases perhaps even from regions smaller than the BLR. For five sources, the upper limit on the reprocessor size is further out, on parsec scales, leaving open the possibility that a portion of the Fe Kα flux could arise in the dusty torus. For the other eight sources, the simulations only provide a lower limit or a single value (without lower and upper limits) of the reprocessor size. In some of the latter cases, this is due to the fact that the light curves have very few observations (≲7–9) to detect any correlation between the continuum and Fe Kα fluxes, or the sampling of the observations in not constraining. In other cases, the Fe Kα line flux does not show variability during the observations, possibly because the reflector is too large compared to the continuum variability timescales, and the simulations are not able to find an upper limit for the reflector.

For Circinus Galaxy, the simulations predict a small reprocessor, with an upper limit of 0.15 pc. This is contrary to several investigations (e.g, Marinucci et al. 2013; Uematsu et al. 2021), which have predicted a larger reflector. This inconsistency could arise from the fact that we are not be adequately modeling the X-ray continuum from this source, as it is heavily obscured. In the cases of NGC 4151 and NGC 5548, the estimated reprocessor sizes appear to be inconsistent with the FWHM-based sizes, with the latter being smaller. This may suggest that the Fe Kα emission has a more extended configuration than the bulk of the sample or that the FWHM sizes (and hence velocity profiles) are not purely virial in nature.

Based on our variability analyses, we single out Cen A and IC 4329A, where significant Fe Kα variability is detected on timescales of days and estimate upper limits for the reprocessor size of 50 and 98 light days, respectively. Fukazawa et al. (2016) studied the narrow Fe Kα emission of both sources using a sample of multiple observations of Suzaku (six observations per source) and XMM-Newton (two observations for Cen A) over periods of 4 and 9 years, respectively, and they only detected significant variability on timescales of 1000–2000 days. The main difference between these results and our could be related to the number of observations used in each study. We analyze nine observations of IC 4329A and 39 observations of Cen A over a timespan of ∼18 years, which are sensitive to significant variability on shorter (and longer) timescales.

Fukazawa et al. (2016) also analyzed four observations of NGC 4151, finding significant variability on scales of ∼1000 days. By contrast, Zoghbi et al. (2019) studied two decades of XMM-Newton observations of NGC 4151 and found a time delay between the Fe Kα and continuum variability of days. In Sect. 4, we find a lower limit for the reprocessor of ∼550 light days for NGC 4151, suggesting that the bulk of the Fe Kα emission mainly arises from or inside the torus, although the velocity of the line suggests that a portion of the line originates on scales of 40 light days from the center. Our study, by design (i.e., only analyzing spectroscopy on observational exposure timescales), will not be very sensitive to rapid variability of the continuum or Fe Kα and avoids highly complex spectral fits, which may be why we do not more robustly confirm the results of Zoghbi et al. (2019).

Many of our results appear consistent with past findings in smaller samples. One notable exception is the high incidence of radio-loud AGN showing stronger excess variance compared to expectations for a standard disk+corona model. Relativistic beaming coupled with variations in jet or corona structure as a function of time could lead to a high level of associated continuum variability (e.g., Ulrich et al. 1997; Chatterjee et al. 2008; Weaver et al. 2020), even if the jet component is not completely dominant in the X-rays (e.g., Cowperthwaite & Reynolds 2012; Lohfink et al. 2013). In large samples of radio-loud AGN (e.g., MacLeod et al. 2010; Zhu et al. 2021), the X-ray emission still appears to be dominated by the corona, although the variability is lower compared to radio-quiet AGNs due to contaminated from a jetted component. This could explain why four radio-loud sources in our sample are less variable than expected by the standard PSD model, but does not explain the stronger excess variance seen in others.

In Sect. 4.5, we examine possible correlations between the observed X-ray continuum and Fe Kα fluxes. Overall, in most of the cases the F_Fe Kα–F_2−10 keV slope is consistent with zero or has a large scatter. The two sources for which we were able to detect strong correlations between the Fe Kα to continuum variations are IC4329A and Cen A. These two outstanding sources show Fe Kα flux variations on timescales of a few days, again suggesting a small reflector, which is confirmed by the simulations, as mentioned above. The case of Cen A, shown in Fig. 8, is particularly striking, with ≈1 dex continuum variations on timescales of weeks to months, and ≈1.5 dex Fe Kα variations on similar timescales. This level of variation is significantly higher than predicted by the simulations, implying a nonstandard source of variability (perhaps associated with a radio-jet launching region or a highly variable and luminous X-ray emitting knot in the jet at larger distances). IC4329A potentially demonstrates similar levels of variation, although this source only exhibits a strong downward variation during a single epoch, and otherwise shows fairly typical variability behavior.

An intriguing result is the case of NGC 2992, which is the only source in the sample where the spectral and timing analyses estimate a Fe Kα emission size larger than the dust sublimation radius. The spectral analysis finds that for this source υ_FWHM = 1990 ± 672 km s⁻¹, which is consistent with the value found by Murphy et al. (2017), who also analyzed the Chandra-HETG spectrum of NGC 2992. Although we find correlated variability between the intrinsic X-ray continuum and the Fe Kα line (0.24 ± 0.06), our timing analysis locates the reflector in parsec scales (2.2 pc < R_rep < 7.1 pc), which is expected since the Fe Kα line responds to the continuum changes on timescales of years. Weaver et al. (1996) reached a similar conclusion by performing a timing analysis of 16 years of ASCA observations and estimated that the reprocessor is ∼ 3.2 pc away from the X-ray source.

We perform imaging analyses on the Chandra ACIS images for 35 sources in our sample with sufficiently clean images to allow the removal of off-nuclear point-like sources and analyze their radial profiles. We find that two of the sources with the highest signal-to-noise data of the sample show extended emission outside of 2.5″, while we can only provide upper limits on the reflecting region for the rest of the sample. We cannot discard the possibility that faint extended emission might be detected in the rest of the sources if comparably high signal-to-noise data could be obtained. For the Circinus Galaxy and NGC 1068, we find that the Fe Kα emission extends up to 95 and 795 pc, respectively. Extended emission in the Circinus Galaxy was already reported by Marinucci et al. (2013) on scales of tens of parsec, based on radial profile analyses; by Arévalo et al. (2014) who found an overall broadening of the iron line in the Chandra HEG spectrum due to spatial extension along the dispersion direction; and by Andonie et al. (2022) who found that the Fe Kα line width systematically increases as the spectral aperture increases. Similarly, Bauer et al. (2015) found that a portion of the Fe Kα emission is extended up to hundreds of parsec scales in NGC 1068.

Recently, Yi et al. (2021) found that NGC 4388 shows kpc-scale Fe Kα extended emission by comparing the radial profiles of the Chandra ACIS-S PSF with the one of NGC 4388 in the 6.2–6.7 keV band. On the contrary, our analysis finds that NGC 4388 does not present Fe Kα emission outside of ∼450 pc. Our analysis takes a very different approach, as we construct an empirical PSF that we compare to the radial profiles, as well as consider excesses for both continuum+Fe Kα and Fe Kα-only emission (i.e., after subtracting the continuum). On the other hand, Guainazzi et al. (2012) concluded that the Fe Kα emission in Mrk 3 extends up to ≃300 pc by comparing a Chandra image with the PSF. Due to the nature of our imaging study, and in particular the inherent pileup of our empirical PSF, we only probe extended emission beyond 2″5, which for Mrk 3 (z = 0.0135) is equivalent to 717 pc; thus cannot confirm nor refute the possibility of extended emission for this object.

The spectral, timing, and imaging analysis suggest that the Fe Kα emission predominantly arises from a region inside the dust sublimation radius, probably in the BLR, or perhaps in some cases in the outer accretion disk. Typical spectral energy distributions of AGN indicate that a large fraction of the primary X-ray and UV emission is reprocessed by gas in the molecular torus, with the dust sublimation radius only setting a lower limit for the location of the dust. Thus, we might naively expect that the X-ray reflection would occur in the molecular torus. However, significant cold X-ray absorption is found to occur inside this radius, implying that neutral gas torus coexists with the BLR (e.g., Davies et al. 2015). Moreover, Ichikawa et al. (2019) studied a sample of 606 AGN in the BASS survey and demonstrated that an extra neutral gas component, in addition to the gas in the dusty torus, is needed to reproduce the X-ray reflection in AGN. Our results, in addition to the studies mentioned above plus many others (e.g., Risaliti et al. 2007, 2011; Maiolino et al. 2010), provide robust evidence for the presence of abundant dust-free gas inside the dust sublimation radius, which appears to be giving rise to the bulk of reflection in most local AGN. On the other hand, our imaging analysis finds modest Fe Kα extended emission out to ∼0.1-kpc scales for two Compton-thick sources in our sample, suggesting that heavily-obscured reflection; probably in the Compton-thin (N_H ∼ 10²³ cm⁻²) or Compton-thick regime, has a complex structure, with secondary/sub dominant reflection likely arising from the molecular torus, the ionization cone, and potentially broader scale gas in the host galaxy.

7. Conclusions

In this work, we analyzed observations for 38 objects from Chandra and XMM-Newton observatories to constrain the neutral reflecting material in AGN, invoking spectral, timing and imaging analyses of Fe Kα line emission. We fit the X-ray spectra of all the observations available of each source, measuring the Fe Kα line FWHM and flux, and the 2–10 keV continuum flux. Then, we created light curves for the Fe Kα line and the continuum and used them to provide an estimate of the reflector size through simulations that compute the delay between the Fe Kα line and continuum variability. Finally, we perform an imaging analysis using Chandra-ACIS images, to assess whether the sources of our sample are spatially extended.

In total, we are able to provide estimates of the reflector size for 37 out of 38 AGN. We find:

1. For 21 out of 24 sources with a Fe Kα FWHM measurement, under the assumption of virialized motion, we estimate that the bulk of the Fe Kα originates inside the dust sublimation radius, probably in the BLR or outer accretion disk.

2. 37 sources show significant variability in the continuum, but only 18 AGN show variability in the Fe Kα flux. The light curve simulations provide reflector size estimates for 24 sources. Of these, six have continuum excess variances that are ≳2-σ above expectations based on simulated power-spectral models, hinting at a different origin for the variability for those sources; meanwhile four lie ≲2-σ below expectations, implying significant continuum damping. The simple model invoking a universal power-spectral shape accounts for much but not all of the differences in continuum variance between the sources.

3. Nearly half of the sample show observed Fe Kα excess variances equal to or higher (possibly by factors of up to ∼3–10 in a few cases) than the continuum ones, implying very little damping by any reprocessor geometry. Further specialized simulations, analyses, and possibly dedicated observations are required to determine plausible geometries.

4. Despite the wide range of variability properties, our constraints on the Fe Kα photon reprocessor sizes confirm the picture from the FWHMs, whereby for ≈83% of the systems the Fe Kα emitting regions are consistent with being inside R_sub, albeit with looser size limits bounded between ∼10⁻³–10 pc. Also, for eight out of 11 sources, the Fe Kα emission appears to originate from a radius consistent with or larger than R_Hβ. One interesting outlier is NGC 2992, whose FWHM and timing constraints place the Fe Kα-emitting region a factor of ≳5 − ∼10 times that of R_sub.

5. For most of the sources, we do not detect a clear correlation between the observed Fe Kα line and continuum fluxes, although two outstanding objects, Cen A and IC 4329A, show strong ≈0.5–1.5 dex Fe Kα variability on timescales of weeks to months. We find no significant trends between the variability features of the Fe Kα line and AGN/host galaxy properties (e.g., SMBH mass, Eddington ratio, AGN type, radio-loudness and column density).

6. Finally, we find that the Fe Kα radial profiles for 33 of 35 AGN are consistent with being point-like beyond 2.5″, and derive loose upper limits for the reflecting clouds using that radius (∼6 to 7096 pc). In the cases of the Circinus Galaxy and NGC 1068, we observe Fe Kα emission that extends out to 95 ± 15 and 795 ± 176 pc, respectively, suggesting the existence of reflecting clouds on scales of hundreds of parsecs. The extended emission in these sources is by no means dominant, but demonstrates that overall heavily-obscured reflection in AGN has a broad (physical) range of secondary contributions to the total Fe Kα emission.

As a main conclusion, we confirm that the Fe Kα emission in the vast majority of AGN appears to arise from regions smaller than and presumably inside R_sub, yet equal to or beyond R_Hβ, and thus is associated either with either the outer BLR or accretion disk. That being said, the wide variety of continuum and Fe Kα variability properties are not easily accommodated by a universal scenario for the production of the bulk of the Fe Kα emission among local AGNs. More detailed analyses and future observational campaigns are likely required to make further sense of the variability properties hinted at in this study.

In particular, future observations dedicated to modeling the resolved Fe Kα line profiles of bright AGN with the calorimeters on XRISM (XRISM Science Team 2020) and the variability properties thereof, should thus improve our understanding of how reflection arises from the BLR, outer accretion disk or jet base regions. Likewise, observations being acquired over the next several years with eROSITA (Predehl et al. 2021) should provide more systematic constraints on intermediate timescales albeit at relatively low signal-to-noise for bright AGN, while dedicated observing campaigns and archival studies of outliers such as Cen A and IC 4329A should be pursued. Ultimately, dedicated campaigns with Athena (Barret et al. 2018), given its high sensitivity and spectral resolution, will place the best constraints on the location and origin of the ubiquitous Fe Kα emission which emanates from AGN.

As part of this work, we additionally detected an inconsistency between the observed and simulated Chandra PSFs beyond ∼3″–20″, where the simulations systematically underpredict the real data. This is notable both in direct PSF comparisons, normalized at 2″5, between observed and simulated sources, as well in indirect comparisons of fluxes derived from HEG and ACIS annular spectra. We derived a calibration constant between the HEG and 3″–5″ annular spectra of 0.424 ± 0.00507, where we expect a value of 1 if the nominal PSF correction is accurate. Our findings imply that any aperture corrections and flux limits based on simulated Chandra PSFs are likely to be systematically underestimated by at least a few to several percent.

Acknowledgments

This project has received funding from The European Union’s Horizon 2020 research and innovation programme under the Marie Sklodowska-Curie grant agreement No 860744 (CA), ANID – Millennium Science Initiative Program – ICN12_009 (CA, FEB), CATA-Basal – AFB-170002 (CA, FEB, ET), CATA-Puente – ACE210002 (FEB, ET), CATA2 – FB210003 (FEB, ET), FONDECYT Regular – 1190818 (ET, FEB) and 1200495 (CA, FEB, ET), PIA ACT172033 (ET, PA), NCN19_058 TITANS (PA, ET), the Max-Planck Society through a Max-Planck Partner Group with the university of Valparaíso (PA, RC), Fondecyt Iniciacion grant 11190831 (CR), ANID BASAL project FB210003 (CR), the Science and Technology Facilities Council for support through grant code ST/T000244/1 (DMA), and FONDECYT Postdoctorado for the project n. 3190213 (AT).

References

Appendix A: Observational details and best-fitting parameters

Table A.1 provides a summary of each observation used in this paper, together with the best-fitting parameters. Here we only show a random number of entries as an example. The full version of the table is available in electronic format at the CDS via http://cdsweb.u-strasbg.fr/cgi-bin/qcat?J/A+A/ and also reports the fluxes uncertainties.

Appendix B: Pileup and annular spectra recalibration

We cannot use saturation-affected data to study flux variability, and thus we need to know whether each ACIS observation is affected by pileup. A good way to estimate the amount of pileup of an observation is to compute the count rate per frame time of the observation, since it has a relation with the pileup fraction¹³. When the count rate per frame is less than 0.1, the pileup fraction is < 5%. We used this threshold to separate saturated from unsaturated observations.

To compute the total expected count rate per frame time of each observation, we used the 3″–5″ annular spectra for all the ACIS observations, and the following expression:

(B.1)

where the factor 33.3 is the effective aperture correction, which is the inverse of the fractional power encircled in a 3″–5″ annulus, which is approximated as 0.03 based on the digitization and integration of the 1 keV encircled energy profile of Fig. 6.1 from the Chandra Proposer’s Observatory Guide.¹⁴ We assume that the 1 keV value provides a rough average over the entire relevant Chandra energy range.

If the observation is affected by pileup, we adopted the aperture-corrected 3″–5″ annular spectrum, if not, we used the aperture-corrected 1.5″ circular spectrum.

We discovered an inconsistency when comparing fluxes given by the 1.5″ aperture and 3″–5″ annular spectra, suggesting an inconsistency in the aperture correction for the annular spectrum. Given that the central aperture spectra could be affected by saturation, we performed a consistency check between the 3″–5″ annular spectrum with the HEG spectrum in the 2-10 keV energy range, for all the observations performed in grating mode. To do this, we used PyXspec to fit the following simple model:

(B.2)

The description of each component is the same as in Section 3, aside from the CONSTANT component, which is added to compute the calibration constant between the spectra. The constant value was fixed to 1 for the grating spectrum and let free to vary for the annular spectrum.

The upper panel of Figure B.1 shows the calibration constants calculated for all the BASS sources with HEG spectra. We obtained an evident consistency for 43/56 galaxies. For 13 sources, higher values of the calibration constant are found, and therefore we did not use the 3″–5″ annular spectra for those cases. A possible explanation for these outliers is the presence of extended emission, either from the AGN or circumnuclear contamination. Some of the sources with a constant higher than one are already claimed to present extended emission, such as the Circinus Galaxy (Arévalo et al. 2014), NGC 1068 (Bauer et al. 2015) and NGC 4945 (Marinucci et al. 2017). A deeper analysis of the rest of the sources, beyond the scope of the present work, is needed to confirm the possibility of extended, but for now this is the most probable explanation.

Fig. B.1. Calibration constants between the Chandra HEG first-order and ACIS/HETG 0th order 3″–5″ annular spectra. The upper panel shows all the calibration constants, and the lower panel shows the calibration constants consistent with the "floor" found at 0.424 ± 0.005, marked by a red dashed line.

The lower panel of Fig. B.1 shows the values of the calibration constants for the 43 consistent sources. If the aperture corrections of both the HETG and annular spectra are accurate, we would expect a constant consistent with unity. However, the average value for the constant is 0.424 ± 0.005, which is highly constrained. Presumably the few percent of remaining variation in the lower panel of Fig. B.1 is due to differences in underlying spectral shape amongst the sources.

The CIAO aperture corrections are based on a theoretical model of the Chandra point spread function (PSF), and hence the deviation from unity for the above constant implies a likely problem with the wings of the model. This inconsistency between the theoretical/model and real PSFs is discussed on the Chandra X-ray Center website, in $ 4.2.3 of the Proposer’s Observatory Guide¹⁵ and the ChaRT webpages,¹⁶. In particular, uncertainties of 30–50% are thought to exist at high energies, consistent with our findings.

We use the average value of the constant to rescale the fluxes obtained for all the annular spectra. We further compare these against the 15 aperture fluxes obtained from the non-piled-up spectra, and confirm consistency within statistical errors. We revisit the issue of the PSF wings in the analysis of the Chandra images, presented in §5.

Appendix C: Notes on individual sources

• 3C 120: this radio-loud sources has σ_c/σ_c, sim < 1, suggesting that the X-ray emission of this object is not jet-dominated but coronal-dominated, as predicted by Zhu et al. (2021), which also found that radio-loud objects coronal-dominated are 20% less variable than radio-quiet objects. Additionally, Lohfink et al. (2013) argue that 3C 120 may have a truncated disk, depleted by the ejection of a new superluminal knot.

• MRK 273: ultra-luminous infrared galaxy with dual obscured nuclei separated by ≈1″. Our extraction and spectral model does not attempt to distinguish between the NE and SW components, which have distinct spectral shapes and variability properties (Iwasawa et al. 2018; Liu et al. 2019); while the more variable obscured SW component dominates the < 10 keV emission, the Compton-thick NE component dominates above ∼12 keV. This added complexity may lead to its moderately high σ_c/σ_c, sim value.

• NGC 1365: this source has a σ_c/σ_c, sim > 1, not following the typical radio-quiet power spectrum. It is a well-known changing-look AGN, with N_H and covering factor variations between ∼10²²–10²⁴ cm⁻² and ∼0.2–1.0, respectively (e.g., Rivers et al. 2015); while our simple spectral model accounts for the variable N_H, it does not incorporate covering factor changes, which appear to increase the continuum variability over expectations from the X-ray corona.

• NGC 2992: this is the only source in the sample with a reflector size estimate from the simulations and Fe Kα FWHM consistently larger than the dust sublimation radius. Also, σ_c/σ_c, sim ∼ 50, being much more variable than expected for a radio-quiet AGN. This intriguing AGN has optically changed from Seyfert 1.5 to Seyfert 2 (Trippe et al. 2008), and also has transitioned from a passive state to an active state in the X-ray band (Gilli et al. 2000). Long- and short-term variability with large amplitude (≳1 order of magnitude) has already been reported for this source by many studies (e.g., Murphy et al. 2007; Guolo et al. 2021; Ghosh & Pal 2021).

• NGC 7469: this is the only source that does not reproduce the well-known radio-quiet power spectrum which seems to defy explanation, which by most accounts appears to be a fairly typical unobscured AGN, yet appears to have the lowest measure of σ_c/σ_c, sim on record. Based on Fig. D.1, NGC 7469’s 2–10 keV continuum flux appears to vary by a factor of ≈1.4 on 1-year timescales, but remains within that window over the full extent of the light curve. The measurements from 2015-2016 are fully consistent with the Swift XRT light curve published by Behar et al. (2020), who find a factor of ≈2 variation in 2–10 keV count rate on 80-day timescales. The only "peculiarity" is that the AGN in NGC 7469 appears to show two UV to X-ray delay timescales of 0.37 and −3.5 days (e.g., Pahari et al. 2020). It could be that there is unresolved emission contaminating the aperture, damping the expected variability, although its hard X-ray luminosity is sufficiently high that any contaminant would have to be very atypically bright (and constant) itself.

Appendix D: Light curves

Figures D.1 and D.2 show the continuum and Fe Kα light curves, respectively, for all of the sources in our sample. Machine-readable files for each source will be provided upon publication.

Fig. D.1. Light curves for the 2–10 keV continuum flux (F2−10 keV) for all sources in the sample. In each subplot, green stars denote Chandra HEG spectra, black stars denote Chandra ACIS 3″-5″ annular spectra, blue stars denote 1.5″ circular Chandra ACIS spectra, and yellow stars denote XMM-Newton pn spectra. The units of the flux are erg cm−2 s−1.

Fig. D.2. Light curves for the Fe Kα line flux (FFe Kα) for all sources in the sample. In each subplot, green stars denote Chandra HEG spectra, black stars denote Chandra ACIS 3″-5″ annular spectra, blue stars denote 1.5″ circular Chandra ACIS spectra, and yellow stars denote XMM-Newton pn spectra. The units of the flux are erg cm−2 s−1.

Fig. D.3. Comparison between Fe Kα line (F6.4 keV) and 2–10 keV continuum (F2−10 keV) fluxes, for all sources in our sample. In each subplot, green stars denote Chandra HEG spectra, black stars denote Chandra ACIS 3″-5″ annular spectra, blue stars denote 1.5″ circular Chandra ACIS spectra, and yellow stars denote XMM-Newton pn spectra. The thick red line shows the best-fit linear regression, while thin red lines represent fits from Monte Carlo resampling to demonstrate uncertainties; the best-fit slope and error are indicated in the legend of each plot. Fits to sources with fewer than five measurements in both fluxes are often poorly constrained. The values for the slopes are reported in the label of each plot. The flux units are erg cm−2 s−1.

Appendix E: Variability versus AGN and host galaxy properties

Figures E.1 and E.2 compare the properties of the continuum and Fe Kα variability to various AGN and host galaxy properties.

Fig. E.1. Comparison of the σFe Kα/σc slopes with respect to (a) SMBH MBH, (b) Eddington ratio, (c) line-of-sight column density, (d) radio-loudness, and (e) AGN type. For the first five plots, the thick red line denotes the best-fit linear regression, while thin red lines denote fits from Monte Carlo resampling to estimate uncertainties; the best-fit slope and error are indicated in the legend of each plot.

Fig. E.2. Comparison of the FFe Kα–F2−10 keV slopes with respect to (a) Eddington ratio, (b) line-of-sight column density, (c) radio-loudness, and (d) AGN type. For the first three plots, the thick red line denotes the best-fit linear regression, while thin red lines denote fits from Monte Carlo resampling to estimate uncertainties; the best-fit slope and error are indicated in the legend of each plot.

Appendix F: Radial profiles

Figure F.2 shows a comparison of the continuum-subtracted Fe Kα radial profiles for each source in our sample with the continuum-subtracted Fe Kα radial profile of HERX1, which serves as our nominal PSF. The HERX1 radial profile is renormalized to match each source light curve at 25, as pixels inside this radius are potentially affected by pileup both in the sample sources and HERX1; we renormalized at 4″ in the case of Cen A, as the pileup is worse for that source. In a few cases, the sample AGN are embedded in strong, diffuse cluster emission, and thus we crop the radial profiles at the point where the cluster background starts to dominate.

Similarly, Fig. F.1 shows a comparison of the 6.2–6.5 keV band (Fe Kα line and continuum) for each source in our sample and the 6.2–6.5 keV profile of HERX1, renormalized at 25.

Fig. F.1. Comparison between the 6.2–6.5 keV radial profile of HERX1 (black) and all the sources in the sample (red). The green dotted line on each plot represents the renormalization radius (i.e., pileup-free radius), which is 25 for all the sources except Cen A and 2MASXJ23444, where we adopt 4″. The errors plotted are 99.987% confidence.

Fig. F.2. Comparison between the Fe Kα radial profile of HERX1 (black) and all the sources in the sample (red). The green dotted line on each plot represents the renormalization radius (i.e., pileup-free radius), which is 25 for all the sources except Cen A and 2MASXJ23444, where we adopt 4″. The errors plotted are 99.987% confidence.

Finally, Fig. F.3 compares the continuum-subtracted Fe Kα radial profiles split into four quadrants for all of the sources in our sample. The profiles have been divided by 25% of the renormlized HERX1 radial profile used in Fig. F.2, such that azimuthally symmetric profiles consistent with the PSF should have values consistent with 1 at all radii.

Fig. F.3. Comparison of continuum-subtracted Fe Kα radial profiles per quadrant (q1, q2, q3, and q4, as defined in Appendix Fig. F.4) for all the sources in the sample. The profiles have been divided by 25% of the renormlized HERX1 radial profile used in Fig. F.2, such that azimuthally symmetric profiles consistent with the PSF should have values consistent with 1 at all radii. The errors plotted are 99.987% confidence. The radial bins of the different quadrants have been shifted slightly so that the error bars are easier to compare.

Fig. F.4. Fe Kα images of the central 18″ × 18″ of NGC 1068 (top) and the central 8″ × 8″ of Circinus Galaxy (bottom). The quadrants are denoted in the image.

Appendix G: Contribution from multiple reflectors

The size of the reflector, which is estimated through the relative variance of the Fe Kα light curve compared to the continuum light curve, can be viewed as a weighted average of the distances of different reflecting structures. This method uses the ratio between the flux-normalized variances of the intrinsic (continuum) and reflected (Fe Kα) light curves to estimate the size of the reflector. Therefore when comparing the contribution of different reflectors, the key is to understand how they all affect the measured variance of the reflected light curve.

Assuming only two different reflecting regions for simplicity of the argument, we can construct the Fe Kα light curve as the sum of two reflected light curves R₁ and R₂. The flux-normalized variance of this light curve will be

(G.1)

(G.2)

(G.3)

(G.4)

(G.5)

i.e., the weighted average of the flux-normalized variances of each reflected light curve, where the weighting factor is the square of the fraction of Fe Kα flux contributed by each reflector, plus a mixed term. Separating the reflected light curves into their mean fluxes and their deviations around the mean δ_1, 2, the part within parentheses of the mixed term (G.5) can be expressed as

(G.6)

(G.7)

(G.8)

(G.9)

(G.10)

(G.11)

If the reflectors are sufficiently different, the deviations around the mean fluxes of their reflected light curves will be shifted in time by different delay times and therefore will tend to be uncorrelated. In this case, the mixed term will tend to vanish and the variance of the Fe Kα light curve tends toward the weighted average discussed above.

Appendix H: Poisson noise estimate using simulations

As described in the main text we used flux-randomization to estimate the contribution of observational noise to the variance. The variance of the flux-randomized light curve can be viewed as

(H.1)

where lc is the flux-randomized light curve, s is the intrinsic variability signal, e are the deviations of each flux measurement due to observation error and n are the flux deviates from the flux randomization. We have avoided the explicit time dependence of all these quantities for simplicity. Assuming that the intrinsic signal, the sign of the observational error and of the additional noise are uncorrelated and that the last two are distributed around zero, the cross terms and the expectation value of n and e can be canceled, which results in

(H.2)

(H.3)

i.e., the intrinsic variance, plus the variance of the noise, , and the variance of the additional noise, . Therefore, taking the difference between the variance of the flux-randomized light curve and the original light curve results in an estimate of the variance of the noise, i.e., . Repeating the flux-randomization many times then produces the distribution of possible contributions of observational noise to the variance, for the particular flux and error values of a given light curve.

All Tables

All Figures

	Fig. 1. Some X-ray properties of our sample. Top figure: intrinsic 2–10 keV luminosity versus redshift of our sample. Magenta circles are radio-loud sources (see Table 1 for details), black circles are Compton-thick AGNs, green circles are dual AGNs, and gray circles are normal, radio-quiet AGNs. Bottom figure: distribution of X-ray photon indices of our sample.
In the text

Fig. 2. 2–10 keV spectra of the Type 2 AGN Centaurus A, extracted and aperture-corrected from 3″–5″ annuli in ACIS-I ObsIDs 7800 (black) and 8489 (red), fitted by the model described in Eq. (1). The top panel shows the spectra and model-to-data ratio in instrument units while the bottom panel shows the unfolded spectra. The spectra are attenuated by strong absorption. The sole emission line is the 6.4 keV Fe Kα line.

In the text

Fig. 3. Left: Fe Kα emitting radius (RFe Kα) versus dust sublimation radius (Rsub) in units of parsecs. Rsub values are estimated using Eq. (3) and RFe Kα values are estimated using Eq. (2). The gray shaded area represents RFe Kα ≥ Rsub. Right: Fe Kα emitting radius (RFe Kα) versus optical BLR radius (RHβ) in units of parsecs. RHβ is directly calculated from Hβ reverberation lags. The gray shaded area represents RFe Kα ≥ RHβ.

In the text

Fig. 4. Variability properties of the sample. Sources are ordered based on the ratio between the excess variances of the continuum and the Fe Kα line, . All four plots adopt the same color code as in Fig. 4a. Open black circles denote AGN without a special characteristic, magenta empty circles are radio-loud AGN, black filled circles are CTȦGN, green empty circles are dual AGN, and filled cyan circles are changing-look AGN. (a): normalized excess variance of the continuum light curves. (b): normalized excess variance of the Fe Kα line light curves. (c): the ratio between the normalized excess variance of the Fe Kα line and continuum light curves, where the four horizontal lines denote ratio for simulated reprocessors with diameters dr = 0,2,10,100 ld, assuming identical 10 days power-spectral bend timescales. (d): ratio between the normalized excess variance of the simulated and real continuum light curves, with unity denoted by the red line. The light red shaded region denotes up to a factor 2 error in our assumed normalization of the simulated power-spectral model. Magenta empty circles are radio-loud sources, black filled circles are Compton-thick sources, green empty circles represents dual AGN, and cyan empty circles are changing-look AGN. Mrk 3 is both Compton-thick and radio-loud, while NGC 1365 is both radio-loud and changing-look. Plotted values represent the median of each parameter while the error bars correspond to the 16% and 84% bounds of the normalized excess variance distributions.

In the text

Fig. 5. Power spectra of simulated X-ray light curves. The underlying power-spectral shape of the continuum was chosen as a bending powerlaw model with a bend timescale of 10 days (black) and a given realization of this model produces the power spectrum plotted in red. Reflecting the light curves over a spherical reflector of radius R = 50 light days suppresses the high frequency power and results in the model and realization plotted in cyan and blue, respectively.

In the text

Fig. 6. Example light curve realizations of the continuum (red) and reflected light (blue), which correspond to the power-spectral models plotted in Fig. 5. The reflected light curve reduces variability amplitude and is delayed by 50 days. Increasing the size of the reflector further reduces the amplitude of the reflected light curve and shifts it further forward in time. A random selection of observing epochs are plotted, representing measurements of the continuum (black dots) and Fe Kα line (cyan stars) fluxes. The sampled light curves are not obviously correlated, although the parent light curves are.

In the text

Fig. 7. Comparison of simulated Fe Kα and 2–10 keV continuum mean-normalized fluxes. The continuum simulations use a power-spectral bend timescale of 10 days, while three different reprocessor diameters (dr) are chosen: 2 light days (black dots), 10 light days (red stars) and 100 light days (blue open squares). The simulated light curves are resampled using the observing epochs shown in Fig. 6. The different sets highlight how varying the reprocessor size can damp potential flux-flux correlations. Evidently, the larger the reprocessor, the larger the scatter in the flux-flux relation, and the flatter the best fitting linear regression (solid lines with corresponding colors).

In the text

Fig. 8. Light curves for the 2–10 keV continuum (left) and Fe Kα line (right) fluxes for the galaxy Centaurus A (upper panels), 1H0707-495 (middle panels) and NGC 1275 (lower panels). The green and black stars denote fluxes measured from the Chandra HEG first-order spectra and Chandra ACIS 3″–5″ annular spectra, respectively, while yellow stars are fluxes measured from XMM-Newton pn observations. The flux units are erg cm−2 s−1 and the error plotted are 1σ confidence.

In the text

	Fig. 9. Comparison between the FFe Kα–F2−10 keV slopes and the SMBH mass (MBH).
In the text

	Fig. 10. Comparison between real and simulated radial profiles for the X-ray binary HERX1. Panels a,b: radial profiles for the 6.2–6.5 keV and 5–6 keV ranges, respectively. The profiles are renormalized at 25, to highlight the regime interior to which pileup effects may occur in the real data. The errors of profiles denote 99% confidence.
In the text

	Fig. 11. Comparison of 6.2–6.5 keV radial profiles per quadrant (q1, q2, q3, and q4, as defined in Appendix Fig. F.4) for HERX1. The profiles are statistically consistent, indicating that the PSF appears to be approximately symmetric out to at least ∼100″ in the energy range which includes Fe Kα emission for our sample. The errors plotted are 99% confidence.
In the text

Fig. 12. Summary of the reprocessor size estimations computed through the spectral analysis (red stars/limits), the timing analysis (black squares/limits), and the imaging analysis (green lines/limits). The dust sublimation radius Rsub is plotted in blue as a reference. Black squares without errorbars are sources where the simulations could not provide upper or lower limits for the size of the reflector. The green shaded areas in Circinus and NGC 1068 represent the extent of the extended emission that we observe in the imaging analysis.

In the text

	Fig. B.1. Calibration constants between the Chandra HEG first-order and ACIS/HETG 0th order 3″–5″ annular spectra. The upper panel shows all the calibration constants, and the lower panel shows the calibration constants consistent with the "floor" found at 0.424 ± 0.005, marked by a red dashed line.
In the text

	Fig. D.1. Light curves for the 2–10 keV continuum flux (F2−10 keV) for all sources in the sample. In each subplot, green stars denote Chandra HEG spectra, black stars denote Chandra ACIS 3″-5″ annular spectra, blue stars denote 1.5″ circular Chandra ACIS spectra, and yellow stars denote XMM-Newton pn spectra. The units of the flux are erg cm−2 s−1.
In the text

	Fig. D.2. Light curves for the Fe Kα line flux (FFe Kα) for all sources in the sample. In each subplot, green stars denote Chandra HEG spectra, black stars denote Chandra ACIS 3″-5″ annular spectra, blue stars denote 1.5″ circular Chandra ACIS spectra, and yellow stars denote XMM-Newton pn spectra. The units of the flux are erg cm−2 s−1.
In the text

Fig. D.3. Comparison between Fe Kα line (F6.4 keV) and 2–10 keV continuum (F2−10 keV) fluxes, for all sources in our sample. In each subplot, green stars denote Chandra HEG spectra, black stars denote Chandra ACIS 3″-5″ annular spectra, blue stars denote 1.5″ circular Chandra ACIS spectra, and yellow stars denote XMM-Newton pn spectra. The thick red line shows the best-fit linear regression, while thin red lines represent fits from Monte Carlo resampling to demonstrate uncertainties; the best-fit slope and error are indicated in the legend of each plot. Fits to sources with fewer than five measurements in both fluxes are often poorly constrained. The values for the slopes are reported in the label of each plot. The flux units are erg cm−2 s−1.

In the text

Fig. E.1. Comparison of the σFe Kα/σc slopes with respect to (a) SMBH MBH, (b) Eddington ratio, (c) line-of-sight column density, (d) radio-loudness, and (e) AGN type. For the first five plots, the thick red line denotes the best-fit linear regression, while thin red lines denote fits from Monte Carlo resampling to estimate uncertainties; the best-fit slope and error are indicated in the legend of each plot.

In the text

Fig. E.2. Comparison of the FFe Kα–F2−10 keV slopes with respect to (a) Eddington ratio, (b) line-of-sight column density, (c) radio-loudness, and (d) AGN type. For the first three plots, the thick red line denotes the best-fit linear regression, while thin red lines denote fits from Monte Carlo resampling to estimate uncertainties; the best-fit slope and error are indicated in the legend of each plot.

In the text

	Fig. F.1. Comparison between the 6.2–6.5 keV radial profile of HERX1 (black) and all the sources in the sample (red). The green dotted line on each plot represents the renormalization radius (i.e., pileup-free radius), which is 25 for all the sources except Cen A and 2MASXJ23444, where we adopt 4″. The errors plotted are 99.987% confidence.
In the text

	Fig. F.2. Comparison between the Fe Kα radial profile of HERX1 (black) and all the sources in the sample (red). The green dotted line on each plot represents the renormalization radius (i.e., pileup-free radius), which is 25 for all the sources except Cen A and 2MASXJ23444, where we adopt 4″. The errors plotted are 99.987% confidence.
In the text

Fig. F.3. Comparison of continuum-subtracted Fe Kα radial profiles per quadrant (q1, q2, q3, and q4, as defined in Appendix Fig. F.4) for all the sources in the sample. The profiles have been divided by 25% of the renormlized HERX1 radial profile used in Fig. F.2, such that azimuthally symmetric profiles consistent with the PSF should have values consistent with 1 at all radii. The errors plotted are 99.987% confidence. The radial bins of the different quadrants have been shifted slightly so that the error bars are easier to compare.

In the text

	Fig. F.4. Fe Kα images of the central 18″ × 18″ of NGC 1068 (top) and the central 8″ × 8″ of Circinus Galaxy (bottom). The quadrants are denoted in the image.
In the text