# The Reliability of Accretion Disk Inclination Derived from X-ray
Spectroscopy of Active Galaxies

Rong Du School of Physics, Peking University, Beijing 100871, China
Department of Astronomy, School of Physics, Peking University, Beijing 100871,
China  Yuanze Ding Division of Physics, Mathematics and Astronomy, California
Institute of Technology, Pasadena, CA 91125, USA  Luis C. Ho Department of
Astronomy, School of Physics, Peking University, Beijing 100871, China  Kavli
Institute for Astronomy and Astrophysics, Peking University, Beijing 100871,
China  Ruancun Li Department of Astronomy, School of Physics, Peking
University, Beijing 100871, China  Kavli Institute for Astronomy and
Astrophysics, Peking University, Beijing 100871, China

###### Abstract

The inclination angle of substructures in active galaxies gives insights into
physical components from scales of the vicinity of the central black hole to
the entire host galaxy. We use the self-consistent reflection spectral model
RELXILL to measure the inclination of the inner region of accretion disks with
broadband ($0.3-78\,\rm keV$) X-ray observations, systematically studying the
reliability of this methodology. To test the capability of the model to return
statistically consistent results, we analyze multi-epoch, joint XMM-Newton and
NuSTAR data of the narrow-line Seyfert 1 galaxy I Zwicky 1 and the broad-line
radio galaxy 3C 382, which exhibit different degrees of spectral complexity
and reflection features. As expected, we find that adding more data for
analysis narrows the confidence interval and that multi-epoch, joint
observations return optimal measurements; however, even single-epoch data can
be well-fitted if the reflection component is sufficiently dominant. Mock
spectra are used to test the capability of RELXILL to recover input parameters
from typical single-epoch, joint observations. We find that inclination is
well-recovered at 90% confidence, with improved constraints at higher
reflection fraction and higher inclination. Higher iron abundance and corona
temperature tighten the constraints as well, but the effect is not as
significant as a higher reflection fraction. The spin, however, have little
effect in reflection-based inclination measurements. We conclude that
broadband reflection spectroscopy can reliably measure inner accretion disk
inclination.

Accretion (14), Active galaxies (17), X-ray active galactic nuclei (2035),
Supermassive black holes (1663)

##  1 Introduction

An active galaxy comprises axisymmetric structures on many scales. On the
largest dimensions, the principal axes of the host galaxy of the active
galactic nucleus (AGN) are anchored by the global angular momentum
distribution of the stars, which in general also defines the orientation of
the large-scale interstellar medium. At the opposite extreme, the accretion
disk around the central black hole (BH) establishes its orientation based on
the angular momentum of its latest fueling episode (e.g., Scheuer & Feiler,
1996; King et al., 2005; Li et al., 2015). Launched by AGNs, the relativistic
radio jet that extends from the vicinity of the central BH to possibly beyond
the size of the host galaxy is believed to align with the spin axis of the BH
(Blandford & Znajek, 1977; however, see Natarajan & Pringle, 1998). In the
unified model of AGNs, an axisymmetric dusty torus along the equatorial plane
of the nucleus obscures the optical or ultraviolet emission from the inner
regions (Antonucci, 1993). The structure and inclination of the torus can be
constrained from detailed studies of the infrared spectral energy distribution
(e.g., Zhuang et al., 2018).

The fuelling process of AGNs is expected to affect the alignment of the inner
accretion disk with respect to the stellar disk of the host galaxy (e.g.,
Hopkins et al., 2012). However, in practice, the orientation on the smallest
scales is observed to be poorly aligned with the host galaxy on larger scales
(e.g., Kinney et al., 2000; see review in Section 3.3.3 of Kormendy & Ho,
2013). No preferential alignment exists between AGN jets and the minor axis of
elliptical (Birkinshaw & Davies, 1985) or disk (Wu et al., 2022) galaxies. The
circumnuclear disks traced by water megamasers also appear misaligned with
respect to the host galaxy on global scales (e.g., Greenhill et al., 2009;
Pjanka et al., 2017). Nonetheless, after discarding a minority of highly
misaligned sources, Middleton et al. (2016) report a loose one-to-one
correlation, significant at the $3\sigma$ level, between the inclination of
the galactic stellar disk and the inclination of the AGN accretion disk
derived from reflection spectroscopy, hinting on the possibility that AGNs are
fed by different mechanisms.

Radiation generated from the accretion disk photoionizes the broad-line region
(e.g., Krolik, 1998), whose size and structure can be probed by reverberation
mapping experiments (Blandford & McKee, 1982). Dynamical modeling of the
broad-line region reveals that its structure is usually axisymmetric (disk-
like; e.g., Pancoast et al., 2014; Li et al., 2018; Williams et al., 2018),
although its formation mechanism is still not firmly established. Hypotheses
linking accretion disk winds to the origin of the broad-line region have been
discussed for years (e.g., Emmering et al., 1992; Murray et al., 1995; Czerny
et al., 2017), with recent attention focusing on the scenario of a failed
dusty outflow (e.g., Czerny & Hryniewicz, 2011; Czerny et al., 2016; Baskin &
Laor, 2018). Cross-matching the inclination of the accretion disk and the
inclination of the broad-line region inferred from dynamical modeling can
offer valuable insights into the physical relationship between these two
fundamental components of the AGN central engine.

Measuring the inner accretion disk inclination itself plays a crucial role in
understanding the environment around supermassive BHs and the processes of
accretion and outflow. Theoretically, if we accept general relativity
(Einstein, 1916), the no-hair conjecture (e.g., Misner et al., 1973, pp.
875–877) limits parameters describing the configuration of every Kerr-Newman
BH (Kerr, 1963; Newman et al., 1965) to its mass $M_{\mathrm{BH}}$, charge
$Q$111Charged BHs are usually not considered in an astrophysical context
(however, see e.g., Komissarov, 2022)., and angular momentum $\mathbf{J}$,
which can be decomposed further to a dimensionless spin parameter
$a_{\ast}\coloneqq|\mathbf{J}|c/GM_{\mathrm{BH}}^{2}$ and two angular
components, including the orientation of the BH. The inclination
$\theta_{\mathrm{disk}}$ of the inner accretion disk in AGNs, defined as the
viewing angle of the system with respect to its normal, is believed to align
with the orientation of the central BH because of Lense-Thirring precession
(the Bardeen-Petterson effect: Bardeen & Petterson, 1975; Pringle, 1992;
Papaloizou & Lin, 1995). On a more pragmatic level, knowing the inner disk
inclination helps to improve estimates of the intrinsic bolometric luminosity
of the AGN, assuming that the inner disk aligns with the outer disk emitting
in other wavelengths222This is not always true (e.g., see simulations of Liska
et al., 2021).. Disk inclination measurements also help determine the
intrinsic values of projected quantities, such as radial velocities. Precise
measurements of inclinations would open up the possibility of employing
inclination-dependent theoretical models in data analysis. As an instance in
which X-ray spectroscopy utilizes the forward-folding technique to minimize
the difference between the observed spectra and an instrument-convolved model,
fitting the thermal blackbody continuum of a relativistic, optically thick,
geometrically thin accretion disk (Shakura & Sunyaev, 1973; extended in the
relativistic regime by Novikov & Thorne, 1973) with a weak X-ray corona gives
loose constraints on BH spin, whose accuracy would be improved with knowledge
of the inclination (e.g., Czerny et al., 2011; Done et al., 2013; Reynolds,
2021).

Various efforts have been made to measure the inclination angle of AGN
accretion disks. For sources that have a jet, an approximate orientation of
the beam, which is equal to the inner disk inclination, can be readily
estimated from the radio core dominance parameter (e.g., Ghisellini et al.,
1993; Wills & Brotherton, 1995) via very-long-baseline interferometry
measurements of the fastest part of the jet. Other possible approaches include
fitting the profile of double-peaked broad H$\alpha$ emission lines with a
relativistic disk model (Eracleous & Halpern, 1994; Storchi-Bergmann et al.,
1995; e.g., see Storchi-Bergmann et al., 1997 for NGC 1097 and Ho et al., 2000
for NGC 4450), studying $\mathrm{H_{2}O}$ megamaser disks (e.g., in NGC 4258;
Herrnstein et al., 1999), reproducing near-infrared and mid-infrared
photometric and interferometric observations (e.g., for NGC 1068; Hönig et
al., 2007), and measuring the configuration of circumnuclear gas with high-
spatial resolution, near-infrared spectroscopy (e.g., for NGC 3227, NGC 4151,
and NGC 7469; Hicks & Malkan, 2008). In works that assume the unified model of
AGNs, the inner disk inclination is deemed to roughly equal the inclination of
other substructures, for example, the narrow-line region (Fischer et al.,
2013, 2014) and broad-line region (Wu & Han, 2001; Zhang & Wu, 2002) clouds.

Another approach to derive key disk properties of AGNs, including inclination,
is by fitting the reflection components of their broadband X-ray spectrum
(e.g., Fabian et al., 1989). Empirically, the X-ray emission in AGNs shows
similar properties, leading to the conventional view that electrons from a hot
corona in the vicinity of the BH inverse Compton (1923) scatter thermal
optical/ultraviolet photons from the accretion disk into the X-rays, creating
a power-law continuum with an exponential cutoff (Thorne & Price, 1975; Haardt
& Maraschi, 1991, 1993; Dove et al., 1997; Belmont et al., 2008). Apart from
this power-law continuum, observations of AGNs also identify two distinct
features commonly depicted as the consequence of coronal X-rays reflected by a
standard thin disk (Novikov & Thorne, 1973; Shakura & Sunyaev, 1973): one is
the Fe K$\alpha$ line at $6-7\,\rm keV$ (Fabian et al., 1989; Tanaka et al.,
1995), the other the Compton hump peaking at $20-40\,\rm keV$ (Guilbert &
Rees, 1988; Lightman & White, 1988; Matt et al., 1991). The former originates
from X-ray fluorescence, while the latter traces continuum emission Compton-
scattered into the shape of a hump, with the range determined by photoelectric
absorption on the red edge and Compton recoil, Klein & Nishina (1929) cross-
section, and cutoff temperature on the blue side. As the irradiation process
mainly occurs in the inner disk region, gravitational redshift, relativistic
beaming, and the Doppler effect modify the emission lines, which stretch the
low-energy wing and create a sharper blueshifted peak. In addition, the
relativistically smeared reflection fluorescent lines in the soft band may be
a possible origin of the soft excess (e.g., Fabian et al., 2002). On measuring
the inclination, the key physical process is the broadening of the Fe
K$\alpha$ profile, which offers information on the velocity of the emitting
fluid element as a function of disk inclination and the radius of the
emitters. The red wing of the Fe line, used to measure the BH spin (e.g.,
Reynolds, 2014, 2021), is sensitive to the inner radius of the emitters,
whereas the blue edge is sensitive to the inclination owing to Doppler
boosting (Rees, 1966). Consequently, the inclination in disk reflection models
is degenerate with the spin and the emissivity profile. The iron abundance of
the disk, to the extent that it modifies the strength of the Fe line, also
comes into play. Moreover, according to Fabian et al. (2015), broadband X-ray
spectroscopy can place stringent constraints on the geometry and radiative
mechanisms of the corona. Thus, the measurement of inclination can be
influenced by the coronal temperature and the overall proportion of the
reflection component (reflection fraction, $R_{f}$), which affect the
spectrum’s flux, the shape of the Compton hump, and the relative strength of
the Fe line.

Numerical models have been developed for calculating reflection spectra.
PEXRAV (Magdziarz & Zdziarski, 1995), REFLIONX (Ross & Fabian, 2005), and
XILLVER (García & Kallman, 2010; García et al., 2011, 2013) specialize in
calculating the intrinsic spectrum. For relativistic smearing, kernels like
DISKLINE (Fabian et al., 1989), LAOR (Laor, 1991), KYRLINE (Dovčiak et al.,
2004), KERRDISK (Brenneman & Reynolds, 2006), and RELLINE (Dauser et al.,
2010, 2013) have been readily applied. Advanced reflection models that
calculate broadband relativistic reflection continuum include the KY package
(Dovčiak et al., 2004, 2022), REFLKERR (Niedźwiecki et al., 2019), RELTRANS
(Ingram et al., 2019; Mastroserio et al., 2021), and RELXILL (Dauser et al.,
2014; García et al., 2014). As a self-consistent, angle-dependent scheme,
RELXILL is developed by point-wise convolving XILLVER with RELLINE. The former
model derives the reflected radiation field from a thin slab using atomic
calculations from XSTAR (Bautista & Kallman, 2001; Kallman & Bautista, 2001),
and the latter relativistically blurs the output of XILLVER using a general
relativistic ray-tracing technique.

Our larger aim, the subject of forthcoming works in this series, is to explore
the connection between the inner accretion disk inclination and other facets
of the AGN and its host galaxy. As an initial step, this paper focuses on
ascertaining the degree to which the RELXILL reflection model can deliver
reliable measurements of the inner disk inclination under conditions that
mimic observations of interest. Systematic uncertainties need to be understood
given the complexity of the model and the degeneracy of model parameters. Our
method is explained in Section 2. Section 3 describes the tests and results on
observed spectra, and Section 4 presents the analysis of mock spectra. Section
5 discusses the sensitivity of the inclination measurements to the reflection
fraction, iron abundance, coronal temperature, model selection, and some
technical concerns about binning. Conclusions appear in Section 6. We adopt a
cosmology with $H_{0}=70\,\rm km\,s^{-1}\,Mpc^{-1}$, $\Omega_{\Lambda}=0.73$
and $\Omega_{m}=0.27$.

##  2 Methodology

###  2.1 Systematic Inclination Measurement

We systematically measure inner accretion disk inclination by analyzing X-ray
spectra with RELXILL as a baseline model. Instead of rigidly fitting the
spectra by brute-force computation, we exploit the flexibility of RELXILL to
evaluate judiciously and systematically models with an increasing number of
emission or absorption components (e.g., soft excess and ionized absorbers).
The fit jointly solves for $\theta_{\mathrm{disk}}$ with other parameters,
including the spin $a_{\ast}$, the iron abundance $A_{\mathrm{Fe}}$, the
power-law photon index $\Gamma$, the ionization state of the inner disk $\xi$,
the electron temperature of the corona $kT_{e}$, the reflection fraction
$R_{f}$, and the radial emissivity profile of the primary source (i.e., the
corona). We implement RelxillCp in the RELXILL family to account
simultaneously for part of the soft excess, the broad Fe line, and the Compton
hump. This choice is motivated by the thermal Comptonization model (NTHCOMP;
Zdziarski et al., 1996; Życki et al., 1999) of the corona in RelxillCp, which
is more physical than the cutoff power-law model employed in, for example,
RELXILL and RelxillD.

To achieve a wide energy range from 0.3 to 78 keV, we combine the high-quality
spectra from the pn charge-coupled device (CCD; Strüder et al., 2001) and the
Metal Oxide Semiconductor (MOS; Turner et al., 2001) of the European Photon
Imaging Camera (EPIC) detectors onboard the X-ray Multi-Mirror Mission (XMM-
Newton; Jansen et al., 2001) to capture complex emission and obscuration
features in the softer bands, together with spectra from the Nuclear
Spectroscopy Telescope Array (NuSTAR; Harrison et al., 2013) mission, which
extends the bandpass coverage to include the Compton hump and offers more
information on the Fe line. For basic discussions, we always take advantage of
XMM-Newton (pn, MOS1, and MOS2) to provide the high-resolution soft X-ray
spectrum of the soft excess and complicated absorption components, which are
crucial in broadband fitting.

We fit the spectra with RelxillCp in the RELXILL family (version 2.2), as
implemented in Xspec (12.12.1; Arnaud, 1996). During the fit, we make the
underlying assumption that the “disk-related” parameters
($\theta_{\mathrm{disk}}$, $a_{\ast}$, and $A_{\mathrm{Fe}}$), except for
ionization $\xi$, do not vary on our timescales of interest, but that the
“corona-related” parameters ($\Gamma$, $kT_{e}$, $R_{f}$, and the emissivity
parameters, including the two indices Index1 and Index2, as well as the broken
radius $R_{\mathrm{br}}$), can vary freely. The fixed parameters are tied
among the spectra from different epochs; the rest are untied. For simplicity,
we lump $\xi$ together with the set of corona-related parameters hereinafter,
although it is actually a disk parameter. Our baseline model is a combination
of Galactic absorption, a blackbody soft excess, and RelxillCp333In the
terminology of Xspec:
constant*TBabs*(zBBody+RelxillCp). For different combinations of data, we
adopt different technical strategies.

We divide the data and corresponding strategy into five cases:

  1. 1.

Single-epoch XMM-Newton EPIC. The lack of higher energy observation is the
main motivation to use this method. All relevant parameters are set free and
fit jointly.

  2. 2.

Single-epoch XMM-Newton observation and a non-simultaneous NuSTAR observation.
We fix $\theta_{\mathrm{disk}}$, $a_{\ast}$, and $A_{\mathrm{Fe}}$, and fit
the data jointly while setting $\Gamma$, $\xi$, $kT_{\mathrm{e}}$, $R_{f}$,
and the emissivity parameters (Index1, Index2, $R_{\mathrm{br}}$) free.

  3. 3.

Single-epoch XMM-Newton observation and a simultaneous NuSTAR observation. In
this fairly ideal case, one that is commonly adopted in recent reflection
modeling, we tie both the disk-related and corona-related parameters and set
free the cross-calibration constant and normalization parameter.

  4. 4.

Multi-epoch, simultaneous XMM-Newton and NuSTAR observations, together with
additional non-simultaneous NuSTAR data. Where there are multiple epochs of
simultaneous XMM-Newton and NuSTAR observations, adding additional NuSTAR
data, even if not contemporaneous, may be helpful to improve the fit of the
Compton hump because the count rate above 20 keV is relatively low for typical
Seyfert 1 galaxies. Similar to previous cases, we tie relevant parameters
within each simultaneous epoch and leave free the cross-calibration terms. The
additional data are fitted with extra corona-related parameters, except that
the disk-related parameters are tied for all spectra.

  5. 5.

Multi-epoch, simultaneous XMM-Newton and NuSTAR data. We assume that all
parameters in the same epoch are the same but may vary between different
epochs.

In all five cases, the blackbody temperature $kT$ is a free parameter but
assumed constant within all the spectral epochs, while the Galactic absorption
column density is fixed to the value retrieved from the Leiden/Argentine/Bonn
Galactic H I Survey (Kalberla et al., 2005). Note that within each
simultaneous epoch, the parameters for the spectra from different detectors
are tied. Cross-calibration constants ($\sim 1$) account for the differences
between instruments, and normalization parameters are used to compensate for
differences in integrated flux and exposure time. Upon achieving the best-fit
parameters, we generate a posterior parameter distribution to determine their
confidence intervals.

###  2.2 Reliability Tests

In this section, we address the challenges inherent in the RELXILL methodology
and outline our design for reliability tests. Despite its self-consistent
integration of physical parameters, the model’s high dimensionality and
nonlinearity pose challenges, particularly in retrieving global minima during
spectral fitting. This complexity introduces intricacies in interpreting
results, compounded by the multitude of parameter degeneracies involved.

To establish the reliability of our methodology, we systematically discuss
three aspects:

  1. 1.

Statistical consistency: We scrutinize the consistency of model fits across
datasets with varying quality and quantity, with a primary focus on
statistical errors during observed spectra fitting. By analyzing model
performance under diverse data conditions, we identify cases that minimize
errors. Insights and technical strategies to enhance inclination measurements
are provided based on broadband fits to the observed spectra of I Zwicky 1 and
3C 382 under various combinations.

  2. 2.

Computational robustness: Assessing the robustness and accuracy of parameter
retrieval, we acknowledge potential biases introduced by the effectiveness of
the model coding and the theoretical explainability of the model. Parameter
retrieval is evaluated through tests on mock spectra simulated with known
input parameters, offering insights into model performance under idealized
conditions. Our emphasis in this work diverges from previous studies adopting
similar approaches (Bonson & Gallo, 2016; Choudhury et al., 2017; Kammoun et
al., 2018).

  3. 3.

Theoretical generalizability: We explore the model’s capacity to explain
observed data and underlying physics across different scenarios. Although
generalizability primarily relates to the model’s development, observers play
a crucial role in selecting an appropriate model. We provide a detailed
rationale for choosing RelxillCp in Section 5.6.

For a complex model such as this, parameter space exploration and error
calculation are non-trivial. We employ a Markov chain Monte Carlo (MCMC)
algorithm (Metropolis et al., 1953) to generate the posterior distribution of
baseline parameters. We check for chain convergence and use the 90% confidence
interval to document the results of spectral fits, facilitating the assessment
of the consistency of $\theta_{\mathrm{disk}}$ measurements across various
circumstances.

##  3 Fitting Observed Spectra

Table 1: Epochs of I Zwicky 1 Observations Epoch | Observatory | Obs. ID | Effective Exposure
---|---|---|---
|  |  | (ks)
$a$ | NuSTAR | 60501030002 | 53.3
$b$ |  | XMM-Newton
---
NuSTAR
| 0851990101
---
60501030002
| 50.3
---
38.2
$c$ | NuSTAR | 60501030002 | 40.2
$d$ |  | XMM-Newton
---
NuSTAR
| 0851990201
---
60501030002
| 47.2
---
35.9
$e$ | NuSTAR | 60501030002 | 63.9
Table 2: Epochs of 3C 382 Observations Epoch | Observatory | Obs. ID | Effective Exposure
---|---|---|---
|  |  | (ks)
$a$ |  | XMM-Newton
---
NuSTAR
| 0790600101
---
60202015002
| 18.8
---
23.0
$b$ |  | XMM-Newton
---
NuSTAR
| 0790600201
---
60202015004
| 19.0
---
24.6
$c$ |  | XMM-Newton
---
NuSTAR
| 0790600301
---
60202015006
| 21.7
---
20.8
$d$ |  | XMM-Newton
---
NuSTAR
| 0790600401
---
60202015008
| 19.0
---
21.7
$e$ |  | XMM-Newton
---
NuSTAR
| 0790600501
---
60202015010
| 16.4
---
21.0
Figure 1: (a) The XMM-Newton and NuSTAR broadband spectral fit of I Zwicky 1,
showing the best fit that contains the Galactic absorbed reflection model,
warm absorbers (WAs), and ultra-fast outflows (UFOs). For clarity, only the
spectra of epoch $d$ are shown. Note that the spectra of all five epochs are
fitted simultaneously. The residuals in terms of sigmas with error bars of
size one are shown for various models: (b) the phenomenological redshifted
power-law model, (c) RelxillCp, (d) RelxillCp plus a redshifted blackbody
component to represent a mild soft excess, and (e) the best-fit model, with
ionized absorbers added to the model of panel (c). Figure 2: (a) The XMM-
Newton and NuSTAR broadband spectral fit of 3C 382, showing the best fit that
contains the baseline model and a warm absorber (WA). For clarity, only the
spectra of epoch $c$ are shown. For plotting purposes, the XMM-Newton data
have been rebinned to $\rm S/N>30$, while the NuSTAR data have been rebinned
to $\rm S/N>10$. The spectra of five epochs are fitted simultaneously. The
residuals in terms of sigmas with error bars of size one are shown for various
models: (b) the phenomenological redshifted power-law model, (c) RelxillCp,
(d) RelxillCp plus a redshifted blackbody component to represent a mild soft
excess, and (e) the best-fit model, with a WA added to the previous model.
Table 3: Best-fit Parameters for the Continuum of I Zwicky 1

Epoch | $\theta_{\mathrm{disk}}$ | $a_{\ast}$ | $A_{\mathrm{Fe}}$ | $R_{f}$ | $kT_{e}$ | $\Gamma$ | Index1 | Index2 | $R_{\mathrm{br}}$ | $\log\xi$
---|---|---|---|---|---|---|---|---|---|---
| (∘) |  |  |  | (keV) |  |  |  | ($r_{g}$) | ($\mathrm{erg~{}cm~{}s^{-1}}$)
(1) | (2) | (3) | (4) | (5) | (6) | (7) | (8) | (9) | (10) | (11)
$a$ |  |  |  | $>8.95$ | $49.2_{-8.8}^{+9.8}$ | $2.10_{-0.07}^{+0.06}$ | $7.64_{-0.63}^{+0.77}$ | $3.00_{-0.37}^{+0.28}$ | $2.8\pm 0.4$ | $3.27_{-0.10}^{+0.14}$
$b$ |  |  |  | $1.93_{-0.31}^{+0.42}$ | $85.9_{-23.1}^{+22.8}$ | $2.13\pm 0.01$ | $>8.63$ | $2.96_{-0.48}^{+0.36}$ | $2.4_{-0.2}^{+0.5}$ | $3.47_{-0.11}^{+0.09}$
$c$ | $44.3_{-2.4}^{+2.8}$ | $>0.978$ | $3.01\pm 0.34$ | $3.46\pm 0.33$ | $293_{-40}^{+65}$ | $2.09_{-0.07}^{+0.06}$ | $8.81_{-0.77}^{+0.67}$ | $3.00_{-0.25}^{+0.31}$ | $2.0_{-0.6}^{+0.3}$ | $3.47_{-0.16}^{+0.21}$
$d$ |  |  |  | $3.66_{-0.37}^{+0.69}$ | $17.5_{-1.4}^{+1.2}$ | $1.88\pm 0.01$ | $>8.83$ | $2.98_{-0.34}^{+0.23}$ | $2.5_{-0.2}^{+0.3}$ | $3.37\pm 0.08$
$e$ |  |  |  | $2.10_{-0.49}^{+0.72}$ | $134_{-10}^{+6}$ | $2.15_{-0.06}^{+0.05}$ | $6.16_{-1.13}^{+1.03}$ | $2.95_{-0.63}^{+0.33}$ | $2.1\pm 0.3$ | $3.31_{-0.12}^{+0.17}$

Note. — Col. (1): Epoch. Col. (2): Inclination angle of the inner accretion
disk. Col. (3): Dimensionless BH spin. Col. (4): Iron abundance. Col. (5):
Reflection fraction. Col. (6): Temperature of the corona. Col. (7): Photon
index. Col. (8): Power-law index 1 for the broken power-law disk emissivity.
Col. (9): Power-law index 2 for the broken power-law disk emissivity. Col.
(10): Break radius for the broken power-law disk emissivity, in units of the
gravitational radius, $r_{g}=GM/c^{2}$. Col. (11): Ionization parameter of the
accretion disk. The spin, inclination, and iron abundance are fitted
simultaneously for the five epochs, instead of only for epoch $c$.

Table 4: Best-fit Parameters for the Continuum of 3C 382

Epoch | $\theta_{\mathrm{disk}}$ | $a_{\ast}$ | $A_{\mathrm{Fe}}$ | $R_{f}$ | $kT_{e}$ | $\Gamma$ | Index1 | Index2 | $R_{\mathrm{br}}$ | $\log\xi$
---|---|---|---|---|---|---|---|---|---|---
| (∘) |  |  |  | (keV) |  |  |  | ($r_{g}$) | ($\mathrm{erg~{}cm~{}s^{-1}}$)
(1) | (2) | (3) | (4) | (5) | (6) | (7) | (8) | (9) | (10) | (11)
$a$ |  |  |  | $0.24_{-0.06}^{+0.05}$ | $22.1_{-3.2}^{+4.9}$ | $1.75\pm 0.01$ | $5.15_{-0.53}^{+0.43}$ | $0.11_{-0.02}^{+0.05}$ | $23.7_{-3.9}^{+5.2}$ | $3.10\pm 0.22$
$b$ |  |  |  | $0.34_{-0.07}^{+0.08}$ | $15.6_{-1.5}^{+2.3}$ | $1.72\pm 0.01$ | $5.65_{-0.42}^{+0.46}$ | $1.32_{-0.28}^{+0.33}$ | $10.4_{-2.2}^{+2.1}$ | $3.44\pm 0.13$
$c$ | $41.8_{-4.1}^{+2.7}$ | $>0.992$ | $>8.19$ | $0.49_{-0.13}^{+0.10}$ | $23.0_{-3.0}^{+5.8}$ | $1.79\pm 0.01$ | $7.79_{-1.65}^{+1.26}$ | $3.00_{-0.29}^{+0.40}$ | $2.7_{-0.5}^{+0.4}$ | $2.90_{-0.10}^{+0.12}$
$d$ |  |  |  | $0.50_{-0.12}^{+0.07}$ | $16.8_{-1.8}^{+3.0}$ | $1.78\pm 0.01$ | $5.86_{-0.41}^{+0.35}$ | $1.06_{-0.14}^{+0.15}$ | $12.0_{-1.8}^{+2.5}$ | $2.79_{-0.10}^{+0.12}$
$e$ |  |  |  | $0.55_{-0.10}^{+0.05}$ | $29.7_{-5.2}^{+7.5}$ | $1.80_{-0.02}^{+0.01}$ | $8.38_{-0.80}^{+0.51}$ | $1.89_{-0.25}^{+0.34}$ | $4.6_{-0.5}^{+0.6}$ | $3.02_{-0.16}^{+0.14}$

Note. — Col. (1): Epoch. Col. (2): Inclination angle of the inner accretion
disk. Col. (3): Dimensionless BH spin. Col. (4): Iron abundance. Col. (5):
Reflection fraction. Col. (6): Temperature of the corona. Col. (7): Photon
index. Col. (8): Power-law index 1 for the broken power-law disk emissivity.
Col. (9): Power-law index 2 for the broken power-law disk emissivity. Col.
(10): Break radius for the broken power-law disk emissivity, in units of the
gravitational radius, $r_{g}=GM/c^{2}$. Col. (11): Ionization parameter of the
accretion disk. The spin, inclination, and iron abundance are fitted
simultaneously for the five epochs, instead of only for epoch $c$.

###  3.1 The Sources

The prototypical narrow-line Seyfert 1 galaxy I Zwicky 1 has a complicated
soft X-ray absorption structure and a strong Fe K$\alpha$ line (Ding et al.,
2022), which is thought to be evidence of strong relativistic disk reflection
(however, see Reeves & Braito, 2019). The warm absorbers (WAs) in the soft
X-rays, a combination of fast outflows and corresponding shocked circumnuclear
gas (Ding et al., 2022), evolve over long timescales. Changes in ionization
may be driven by the different densities of the clumps crossing the observer’s
line-of-sight, in which the “skin” layer facing the source has a higher
ionization state (Costantini et al., 2007; Silva et al., 2018). I Zwicky 1 was
observed once in 2020 by NuSTAR in a $\sim 5$ day-long exposure and twice by
XMM-Newton, each with a $\sim 70\,\rm ks$ exposure (PI: D. R. Wilkins).
Because of the strong variability, we divide the NuSTAR data into five epochs
according to their overlap with the XMM-Newton observations (as in Ding et
al., 2022, see our Table 1).

The broad-line radio galaxy 3C 382 exhibits a WA, a soft excess, and an Fe
K$\alpha$ line, but no Compton hump (Torresi et al., 2010; Ballantyne et al.,
2014). According to Ursini et al. (2018), the soft excess comes from a warm
Comptonization component with a photon index $\Gamma\approx 2.4-2.5$, a
temperature of 0.6 keV, and a corresponding optical depth of $\sim 20$. These
parameters are all in general agreement with those in radio-quiet Seyfert
galaxies (e.g., Jin et al., 2012; Petrucci et al., 2018). The Fe K$\alpha$
line is broad in two epochs and possibly broad in the other three, which,
together with the lack of a Compton hump, seemingly undermines the reflection
model previously applied to explain the properties of the source (e.g.,
Sambruna et al., 2011). However, given that the lower limit of the coronal
temperature may be as low as $20-40\,\rm keV$ and the reflection fraction is
merely $\sim 0.1$ (Ursini et al., 2018), the reflection model is intrinsically
poorly constrained, and additional information is needed to test it. For our
purposes, the broadening of the iron line is a sufficient indicator of
reflection-based inclination measurements. Compared to I Zwicky 1, 3C 382 has
weaker reflection features and less complex outflow structures. It was
observed simultaneously five times in 2016 (PI: F. Ursini; Table 2).

###  3.2 Data Reduction

For XMM-Newton data, we extract the light curve and spectrum from the pn and
MOS of the EPIC detectors with the System Analysis Software (version 20.0.0)
and the calibration files released on 31 March 2022. For each bundle of
observational data, we use epproc and emproc to retrieve its event lists,
evselect to extract the spectrum of the source and background, and arfgen and
rmfgen to generate, respectively, the ancillary files and redistribution
matrices. Checking for pile-up and background flares reveals that no
corrections are needed for either source. The spectra are extracted from a
circular region with a 40′′ diameter around the source. In the case of the pn
observations, a source-free polygon region excluding serendipitous sources and
chip edges is chosen for background extraction, whereas for MOS observations
the background spectrum is extracted from a 300′′-diameter circular region
within the range of the outer CCD. For I Zwicky 1, we combine the spectra from
the MOS1 and MOS2 detectors with the task epiccombine to enhance the signal-
to-noise ratio (S/N). We rebin the spectra of 3C 382 with specgroup so as to
not oversample the full width at half-maximum resolution by more than a factor
of 3, and to ensure that each spectral bin contains at least 25 counts.

The NuSTAR data are reduced using NuSTARDAS version 2.1.2, with event lists
retrieved from the two Focal Plane Module (FPM) detectors processed by the
standard NUPIPELINE, which is integrated into the HEASoft bundle (6.30.1). We
use calibration files from CALDB version 20220802. A source-free
300′′-diameter circle near the source region is used to measure the
background. For I Zwicky 1, we extract the spectra and light curve with a
small circular region of 30′′ diameter centered on the point source, as
recommended for faint sources, and then combine the data from FPMA and FPMB
following Wilkins et al. (2021). For 3C 382, to stay consistent with XMM-
Newton procedures, we also use a circular region of 40′′ diameter centered on
the point source to extract the spectrum. The source spectra of 3C 382 are
further rebinned with ftgrouppha in HEASoft such that each spectral bin
contains at least 25 counts.

Table 5: Best-fit Parameters for Ionized Absorbers of I Zwicky 1

Epoch | IA | $N_{\mathrm{H}}$ | $\log\xi$ | $v/c$
---|---|---|---|---
|  | ($10^{20}\,\rm cm^{-2}$) | ($\mathrm{erg~{}cm~{}s^{-1}}$) |
(1) | (2) | (3) | (4) | (5)
$b$ | UFO1 | $2.09^{+0.25}_{-0.31}$ | $-1.72^{+0.14}_{-0.12}$ | $0.34\pm 0.01$
UFO2 | $7.07^{+0.87}_{-0.82}$ | $0.85^{+0.12}_{-0.08}$ | $0.112^{+0.0008}_{-0.0013}$
WA1 | $2.46\pm 0.29$ | $1.82^{+0.21}_{-0.27}$ | $0.035\pm 0.004$
WA2 | $6.16^{+0.37}_{-0.33}$ | $0.695\pm 0.063$ | $0.044\pm 0.001$
$c$ | WA3 | $4790^{+590}_{-410}$ | $3.44^{+0.16}_{-0.36}$ | $0.25^{+0.02}_{-0.03}$
$d$ | UFO1 | $2.45\pm 0.21$ | $-0.711_{-0.037}^{+0.035}$ | $0.276_{-0.010}^{+0.008}$
UFO2 | $3.44^{+0.32}_{-0.37}$ | $0.283^{+0.026}_{-0.025}$ | $0.124^{+0.009}_{-0.006}$
WA1 | $1.38^{+0.11}_{-0.07}$ | $1.034^{+0.047}_{-0.051}$ | $0.036^{+0.003}_{-0.002}$
WA2 | $2.42_{-0.21}^{+0.32}$ | $1.21_{-0.08}^{+0.11}$ | $0.046\pm 0.001$
WA3 | $495^{+47}_{-39}$ | $3.64^{+0.14}_{-0.32}$ | $0.22^{+0.02}_{-0.03}$
$e$ | WA3 | $2270\pm 180$ | $3.94^{+0.36}_{-0.46}$ | $0.27^{+0.02}_{-0.04}$

Note. — Col. (1): Epoch. Col. (2): Ionized absorber (IA) component, including
two ultra-fast outflows (UFO1, UFO2) and three warm absorbers (WA1, WA2, WA3).
Col. (3): Hydrogen column density. Col. (4): Ionization parameter. Col. (5):
Outflow velocity divided by the speed of light.

Table 6: Best-fit Parameters for Ionized Absorbers of 3C 382

Epoch | $N_{\mathrm{H}}$ | $\log\xi$ | $v/c$
---|---|---|---
| ($10^{20}\,\rm cm^{-2}$) | ($\mathrm{erg~{}cm~{}s^{-1}}$) |
(1) | (2) | (3) | (4)
$a$ | $16.1_{-3.2}^{+2.8}$ | $1.92\pm 0.15$ | $-0.132_{-0.019}^{+0.009}$
$b$ | $357_{-79}^{+80}$ | $2.88_{-0.04}^{+0.03}$ | $0.007_{-0.006}^{+0.002}$
$c$ | $23.5_{-4.9}^{+4.7}$ | $1.68_{-0.08}^{+0.10}$ | $-0.125_{-0.010}^{+0.019}$
$d$ | $55.2_{-10.8}^{+25.3}$ | $2.51_{-0.11}^{+0.13}$ | $-0.013_{-0.008}^{+0.010}$
$e$ | $25.9_{-4.5}^{+3.2}$ | $2.02_{-0.07}^{+0.0.05}$ | $0.006\pm 0.008$

Note. — Col. (1): Epoch. Col. (2): Hydrogen column density. Col. (3):
Ionization parameter. Col. (4): Outflow velocity divided by the speed of
light.

Table 7: Luminosity, Flux, and Best-fit Parameters for Cross-calibration
Constants of I Zwicky 1

Epoch | $L_{\mathrm{0.1-200\,keV}}$ | $F_{\mathrm{0.3-2\,keV}}$ | $F_{\mathrm{2-10\,keV}}$ | $F_{\textsc{RELXILL}}$ | constant
---|---|---|---|---|---
| ($10^{44}\,\rm erg~{}s^{-1}$) | ($10^{-11}\,\rm erg~{}cm^{-2}~{}s^{-1}$) | ($10^{-11}\,\rm erg~{}cm^{-2}~{}s^{-1}$) | ($10^{-11}\,\rm erg~{}cm^{-2}~{}s^{-1}$) | MOS/PN
(1) | (2) | (3) | (4) | (5) | (6)
$a$ | $0.71\pm 0.02$ | $0.06\pm 0.01$ | $0.39\pm 0.01$ | $0.78\pm 0.02$ |
$b$ | $1.01\pm 0.01$ | $0.53\pm 0.01$ | $0.48\pm 0.01$ | $1.03\pm 0.01$ | $1.02\pm 0.01$
$c$ | $1.19\pm 0.03$ | $0.09\pm 0.01$ | $0.62\pm 0.01$ | $1.26\pm 0.03$ |
$d$ | $0.94\pm 0.01$ | $0.40\pm 0.01$ | $0.51\pm 0.01$ | $0.96\pm 0.01$ | $1.03\pm 0.01$
$e$ | $1.06\pm 0.02$ | $0.09\pm 0.01$ | $0.53\pm 0.01$ | $1.12\pm 0.02$ |

Note. — Col. (1): Epoch. Col. (2): $0.1-200\,\rm keV$ luminosity of the model.
Col. (3): $0.3-2\,\rm keV$ flux of the data. Col. (4): $2-10\,\rm keV$ flux of
the data. Col. (5): $0.3-78\,\rm keV$ flux of the RelxillCp component. Col.
(6): Cross-calibration constants of epochs $b$ and $d$; we set the constant of
PN and FPM to 1 and show the results of the MOS spectra.

Table 8: Luminosity, Flux, and Best-fit Parameters for Cross-calibration
Constants of 3C 382

Epoch | $L_{\mathrm{0.1-200\,keV}}$ | $F_{\mathrm{0.3-2\,keV}}$ | $F_{\mathrm{2-10\,keV}}$ | $F_{\textsc{RELXILL}}$ | constant
---|---|---|---|---|---
| ($10^{44}\,\rm erg~{}s^{-1}$) | ($10^{-11}\,\rm erg~{}cm^{-2}~{}s^{-1}$) | ($10^{-11}\,\rm erg~{}cm^{-2}~{}s^{-1}$) | ($10^{-11}\,\rm erg~{}cm^{-2}~{}s^{-1}$) | MOS1,2/PN; FPMB/A
(1) | (2) | (3) | (4) | (5) | (6)
$a$ | $9.35\pm 0.12$ | $0.41\pm 0.01$ | $3.96\pm 0.05$ | $11.69\pm 0.15$ |  $1.00\pm 0.01$, $1.02\pm 0.01$; $1.00\pm 0.03$
$b$ | $9.73\pm 0.11$ | $0.47\pm 0.01$ | $4.54\pm 0.05$ | $12.60\pm 0.15$ |  $1.03\pm 0.01$, $0.95\pm 0.01$; $1.01\pm 0.03$
$c$ | $10.35\pm 0.13$ | $0.49\pm 0.01$ | $4.61\pm 0.06$ | $13.03\pm 0.16$ |  $1.02\pm 0.01$, $1.03\pm 0.01$; $1.01\pm 0.03$
$d$ | $9.11\pm 0.12$ | $0.45\pm 0.01$ | $4.24\pm 0.05$ | $11.80\pm 0.15$ |  $1.03\pm 0.01$, $1.02\pm 0.01$; $1.01\pm 0.02$
$e$ | $10.98\pm 0.13$ | $0.52\pm 0.01$ | $4.82\pm 0.06$ | $13.52\pm 0.16$ |  $1.02\pm 0.01$, $1.02\pm 0.01$; $1.04\pm 0.03$

Note. — Col. (1): Epoch. Col. (2): $0.1-200\,\rm keV$ luminosity of the model.
Col. (3): $0.3-2\,\rm keV$ flux of the data. Col. (4): $2-10\,\rm keV$ flux of
the data. Col. (5): $0.3-78\,\rm keV$ flux of the RelxillCp component. Col.
(6): Cross-calibration constants; we set the constant of PN and FPMA to 1, and
arrange the results in the sequence of MOS1, MOS2, and FPMB.

###  3.3 Broadband Model Fitting

Based on the data rebinning scheme, we apply the Cash (1979) statistics
($C$-stat) for I Zwicky 1 during statistics minimization, and the $\chi^{2}$
statistics for 3C 382. I Zwicky 1 falls in case 4 and 3C 382 in case 5
discussed in Section 2.1. All the epochs of 3C 382 and epochs $b$ and $d$ of I
Zwicky 1 belong to case 3. We produce a case 4 fit for I Zwicky 1 and a case 5
fit for 3C 382. We adopt the conventional scheme by starting from a
phenomenological power-law model and gradually modify the model components.
Because X-ray spectra are folded, fitting with a simple model first allows us
to get an initial overall sketch of the spectra.

We start by fitting the data with a redshifted power-law model (zPowerlw)
modified by Galactic absorption (TBabs; see Wilms et al., 2000). The broadband
spectra of both sources show deviations from the phenomenological model, which
require the addition of ionized absorbers at $\sim$0.5 keV, a broad Fe
K$\alpha$ line, and a high-energy Compton hump, motivating us to use
RelxillCp. A mild excess in the residuals of the soft band prompts us to add a
redshifted blackbody component (zBBody), but we do not consider intrinsic
absorbers as in previous works (Ursini et al., 2018; Ding et al., 2022). As
discussed for cases 4 and 5, corona-related parameters are untied for non-
simultaneous observations, while disk-related parameters are tied.
Furthermore, we fix parameters of the accretion disk that cannot be derived
directly from our analysis, including its inner radius ($R_{\mathrm{in}}$,
radius of the innermost stable circular orbit), outer radius
($R_{\mathrm{out}}=400\,r_{g}$, with $r_{g}\coloneqq GM/c^{2}$), and density
$n=10^{15}\,\mathrm{cm^{-3}}$. The broadband model fitting procedure is
illustrated in Figure 1 for I Zwicky 1 and in Figure 2 for 3C 382444Figure 2b
shows a discrepancy between XMM–Newton and NuSTAR spectra ($\Delta\Gamma\sim
0.1$), as also reported in Ursini et al. (2018). According to Fürst et al.
(2016) and Ponti et al. (2018), NuSTAR sometimes measures a higher flux than
XMM-Newton, especially between 3 and 5 keV, which may lead to steeper FPM
spectra. In our case the difference in fluxes is present in the entire
overlapping range of 3–10 keV..

We remark on a technical issue related to the blackbody component in best-fit
models. Figure 1e presents the best-fit model for I Zwicky 1 as
constant*TBabs*XSTAR*RelxillCp. In contrast to a direct combination of XSTAR
tables and the baseline model (with zBBody), we exclude the blackbody
component because the successful modeling of complex emission and absorption
features in the soft band with XSTAR tables renders the addition of a zBBody
inconsequential. In this case, if a blackbody is added to the best-fit model,
the reduced $C$-stat changes negligibly ($\Delta
C\text{-stat}/\Delta\mathrm{DOF}=-0.1/-2$). Our best-fit model for I Zwick 1
aligns with that of Ding et al. (2022). Yet for 3C 382 the blackbody component
is crucial. As a test, removing the epoch-independent redshifted blackbody
component leads to a significant reduction in fit statistics
($\Delta\chi^{2}/\Delta\mathrm{DOF}=-186.8/-2$). Ursini et al. (2018) also
emphasized the inadequacy of pure reflection models in capturing the soft
excess. However, this epoch-independent treatment might include biases, which
will be addressed in Section 5.7.

Ionized absorbers exist in both sources. Since we focus on interpreting the
reflection model, the absorption components are fitted following Ding et al.
(2022) for I Zwicky 1 and Torresi et al. (2010) for 3C 382. Specifically, we
fit the ultra-fast outflows and WAs in I Zwicky 1 with five XSTAR absorbers
and five XSTAR additive tables, adopting the parameters of Ding et al. (2022).
We fit the WA in 3C 382 with multiplicative XSTAR tables. The reason for not
adopting the two-state model in Ursini et al. (2018) is that our analysis does
not include Reflection Grating Spectrometer (RGS; den Herder et al., 2001)
data and may not be able to resolve the two states. In terms of fit statistics
over degrees of freedom (DOF), the single-state absorber model returns
$\chi^{2}/\mathrm{DOF}=5682.6/5004=1.136$, whereas the double-state model
yields
$\chi^{2}/\mathrm{DOF}=5641.4/4989=1.131\>(\Delta\chi^{2}/\Delta\mathrm{DOF}=-41.2/-15)$.
The improvement is not significant.

We generate the posterior probability distribution for each parameter using
MCMC implemented in XSPEC. Adopting the Goodman & Weare (2010) algorithm, we
use 200 walkers with a chain of $6\times 10^{6}$ elements ($3\times 10^{4}$
iterations). The proposal distribution is generated from a Gaussian
distribution around the best-fit values. The chain convergence is tested using
graphs to monitor the evolution of parameters and $\chi^{2}$ with regard to
MCMC steps, and with the integrated autocorrelation time $\tau_{f}$ (Goodman &
Weare, 2010) estimated by the function in EMCEE (Foreman-Mackey et al., 2013).
We ensure that the MCMC chains are longer than $\sim 1000\,\tau_{f}$, as
suggested by Sokal (1996). The first $6\times 10^{5}$ elements are rejected as
the burn-in phase prior to storing the chain, in order that the chain forgets
the initialization point and samples in the equilibrium state. The results of
reflection continuum fitting of I Zwicky 1 and 3C 382 are given in Tables 3
and 4, respectively. The best-fit parameters of the ionized absorbers are
shown in Table 5 for I Zwicky 1 (also see Ding et al., 2022, their Table 2),
and in Table 6 for 3C 382.

Last but not least, Tables 7 and 8 provide the calculated luminosity and
fluxes. Cross-calibration constants $\sim 1$ yield good agreement between the
MOS detectors and PN on XMM-Newton, as well as between the FPMA and FPMB on
NuSTAR.

Figure 3: A scatter plot showing inclination measurements from (a) different
combinations of I Zwicky 1 data groups (cases 1–5) and (b) two simultaneous
epochs ($b$ and $d$). The areas shaded in gray show the confidence range for
$\theta_{\mathrm{disk}}$ derived from broadband model fitting (Section 3.3).
The two dashed lines give the average of the upper and lower edges of the 90%
error bars of the two epochs. Figure 4: A scatter plot showing inclination
measurements from (a) different combinations of 3C 382 data groups (cases 1–5)
and (b) five simultaneous epochs ($a-e$). The areas shaded in gray show the
confidence range for $\theta_{\mathrm{disk}}$ derived from broadband model
fitting (Section 3.3). The two dashed lines give the average of the upper and
lower edges of the 90% error bars of the five epochs.

###  3.4 Inclination with Different Combinations

We rearrange the observed data to mimic the different circumstances (cases)
discussed in Section 2.1. For I Zwicky 1, the result from epoch $b$ is
selected to present a case 3 fit. We produce a case 1 fit with the XMM-Newton
observation from epoch $b$, a case 2 fit with the XMM-Newton spectra from
epoch $b$ and NuSTAR spectra from epoch $d$, and a case 5 fit with epochs $b$
and $d$. The broadband model fitting in the previous section is a case 4 fit.
For 3C 382, the result from epoch $c$ is selected to represent a case 3 fit.
We conduct a case 1 fit with the XMM-Newton observation from epoch $c$, a case
2 fit with the XMM-Newton spectra from epoch $c$ and NuSTAR spectra from epoch
$b$, and a case 4 fit with epochs $b$ and $c$ and NuSTAR data from epochs $a$,
$d$, and $e$. The broadband model fitting in the previous section is a case 5
fit. Additionally, we produce case 3 fits for all available simultaneous
epochs of both sources. Errors are calculated with MCMC, and chain convergence
is assured. We show the inclination measurements of the different data
combinations for I Zwicky 1 in Figure 3, while the results for 3C 382 are
given in Figure 4.

The inclination measurements of each source are internally consistent, in view
of the overlap of the confidence intervals. The consistency of measurements
for all the different cases and different epochs confirms the overall
reliability of our method. Comparing the results from different cases, it is
clear that the more abundant the data, the smaller the errors, which means
that more information provided by the broadband spectra does lead to a better
fit, hence a more precise measurement. Moreover, if the several groups of
spectra are simultaneous, the dimension of the model is reduced, which further
diminishes the flexibility and complexity of the model and improves the fits.
Examining the results from different epochs reveals that the measured value
from the broadband model fitting in Section 3.3 lies within the average value
of the five simultaneous epochs. Inclination measurements for both sources can
be obtained consistently despite the differences in their outflow structure.
Yet, the accuracy of the measurement can be improved if the reflection
features in the source spectra are more dominant. For I Zwicky 1 ($R_{f}>1$),
the inclination measurements from different cases and different epochs only
vary a little with confidence intervals mostly overlapping, especially for the
two simultaneous epochs, whereas the inclination measurements of 3C 382
($R_{f}\lesssim 0.5$) fluctuate between different combinations of data. We
note that combining more data improves the fits more significantly for 3C 382
than for I Zwicky 1. That said, although case 5 is intuitively the most ideal
one among the five cases, the naïve expectation that more data give better
performance does not necessarily hold. We have shown that applying RELXILL to
single-epoch data that exhibit strong reflection signal (I Zwicky 1) suffices
to offer satisfactory results with errors of only several degrees. Yet, for
data with weak reflection features like those in the single epochs of 3C 382,
the fits can yield dubious results with measured inclination angles varying by
tens of degrees.

In summary: inner accretion disk inclinations can be measured more
consistently if the reflection component is more dominant in the spectra or if
more observations are combined and fitted simultaneously.

##  4 Fitting Simulated Spectra

###  4.1 Simulator Settings

We produce mock spectra with the model constant*(zBBody+RelxillCp) using the
fakeit method in XSPEC, as implemented in the Python interface PyXspec v.2.1.0
(Arnaud, 2016). To mimic simultaneous, single-epoch observations (case 3) of
nearby AGNs, each group of spectra is generated in the $0.3-10\,\rm keV$ band
of XMM-Newton/EPIC and the $3-78\,\rm keV$ band of NuSTAR, using response
files from the epoch $c$ observations of 3C 382 and their respective ancillary
files. The source exposure and background exposure time are all 25 ks. As the
background flux of the simulated NuSTAR spectra overwhelms the source flux at
$\sim 50\,\rm keV$, we ignore energy bins above $50\,\rm keV$ for the mock
FPMA/B spectra. The mock spectra are generated by setting $R_{\mathrm{in}}$ to
the radius of the innermost stable circular orbit,
$R_{\mathrm{out}}=400r_{g}$, and $n=10^{15}\,\mathrm{cm^{-3}}$. The blackbody
temperature is set to 0.1 keV, the redshift to $z=0.05557$ (van den Bosch et
al., 2015), and other parameters are initialized as shown in Table 9.

We study all combinations in Table 9, creating 13,860 distinct simulations.
The range for $\theta_{\mathrm{disk}}$ covers nearly face-on to nearly edge-on
systems, which have a high degree of relativistic smearing. The input values
of the iron abundance $A_{\mathrm{Fe}}$ are solar to super-solar to ensure
prominent Fe features and to represent the physical conditions in most
observations (García et al., 2013). The photon index $\Gamma$ and ionization
$\xi$ lie within the physical range in AGNs. The thermal temperature $kT_{e}$
of the coronal electrons shapes the high-energy cutoff. Although the physical
nature of the corona currently cannot be derived from first principles (for
progress in magnetohydrodynamics simulations, see Kinch et al. 2016, 2020 and
Jiang et al. 2019), a firm upper limit on $kT_{e}$ is given by the lower limit
of $\sim$511 keV imposed by the temperature of electron-positron pair
production. Observations typically find coronal temperatures of $\sim 150\,\rm
keV$ (e.g., Ricci et al., 2017), and the cutoff energy set in several models
is as high as $300\,\rm keV$ to even $\sim 500\,\rm keV$ (e.g., Baloković et
al., 2020; Kamraj et al., 2022).

The radial emission profile, which is degenerate with
$\theta_{\mathrm{disk}}$, is intricately linked to the selection of coronal
models—a topic we delve into in Section 5.6. While our primary focus for the
mock tests does not involve a comprehensive exploration of the impact of
radial emissivity parameters (Index1, Index2, $R_{\mathrm{br}}$) on
inclination measurements, we assign them reasonable values and conduct tests
to evaluate their successful recoveries.

The spin parameter influences the entire simulation grid. We first select a
specific value of $a_{\ast}$ for preliminary discussions on other parameters,
and subsequently extend the coverage from $a_{\ast}=-0.998$ to $+0.998$.
Theoretically, higher values of $a_{\ast}$ usually allow for a larger grid for
$R_{f}$ and $\theta_{\mathrm{disk}}$ (see Section 5.4). In the reflection
modeling of actual observations, reported spin values often lean toward higher
values ($a_{\ast}>0.9$; Reynolds, 2021). Therefore, exploring the parameter
space with a higher spin is more informative. Moreover, the degeneracies
intertwined within $a_{\ast}$ that affect the retrieval of other parameters
are fortunately disentangled with higher spin, as accurate retrieval of the
spin is typically achieved for values above $0.8$ in mock analyses (Bonson &
Gallo, 2016; Choudhury et al., 2017; Kammoun et al., 2018). Thus, we initially
set $a_{\ast}$ to its theoretical maximum of 0.998 to mitigate the influence
of spin degeneracy and study this more flexible slice of the parameter space.
In Section 5.4, we extend the discussion to other values of $a_{\ast}$ (Table
9) to assess the effect of spin degeneracy on the precision of inclination
measurement.

Table 9: Parameter Values of the Simulation Grid Parameter | Input Range
---|---
$\theta_{\mathrm{disk}}$ | [3, 15, 30, 45, 60, 75, 87]∘
$A_{\mathrm{Fe}}$ | [1, 3, 5, 7, 9] $A_{\odot}$
$\Gamma$ | 2
$\log\xi$ | 3
$kT_{e}$ | [30, 150, 300] keV
Index1 | 7
Index2 | 3
$R_{\mathrm{br}}$ | 10 $r_{g}$
$a_{\ast}$ | [$\pm 0.998$, $\pm 0.9$, $\pm 0.7$, $\pm 0.5$, $\pm 0.3$, $0$]
$R_{f}$ | [0, 0.1, 0.3, 0.5, 0.7, 0.9, 1, 3, 5, 7, 9, 10]
Figure 5: (a) An example of the simulated case 3 spectra, with the model in
black lines. (b) The $C$-stat residuals between the generated spectra and the
model. The input values are $\theta_{\mathrm{disk}}=45^{\circ}$,
$a_{\ast}=0.998$, $A_{\mathrm{Fe}}=3$, $R_{f}=1$, and $kT_{e}=\rm 150\,keV$.
For plotting purposes, the XMM-Newton data have been rebinned to $\rm S/N>30$,
and the NuSTAR data have been rebinned to $\rm S/N>10$. Figure 6: Results of
the MCMC analysis for the relevant parameters of the spectra shown in Figure
5. The red lines correspond to the measured values, and the blue lines give
the input values of the simulation. In general, all parameters are well
recovered. Figure 7: Inclination measurements for simulated spectra with solar
iron abundance ($A_{\mathrm{Fe}}=1$) and mild coronae
($kT_{e}=30\,\mathrm{keV}$) in the maximum Kerr case ($a_{\ast}=0.998$). The
colored regions depict the 90% confidence intervals for different input
reflection fractions $R_{f}$. The dashed line depicts the expected values of
$\theta_{\mathrm{disk}}$.

###  4.2 Mock Analysis

Retrieving results from mock spectra can be challenging. Since we know the
input parameters in advance, it is technically equivalent to the circumstance
that we have already achieved best-fit values when fitting the observed
spectra. Nevertheless, knowing the inputs is an idealized situation, for in
practice we have little confidence that best-fit parameters are actual model
parameters (global minima), although prior knowledge of the physical
properties of the source may serve as a guide. Hence, we need to calculate
properly the posterior distribution in mock analysis, such that the sampler
explores the parameter space as much as possible. However, for our model with
much flexibility and complexity, the sampler can get stuck easily in local
minima if we start at a random point very far away from the actual values
(e.g., Bonson & Gallo, 2016). Therefore, it is problematic to initialize the
sampler at random values or to start with the known model values.

We analyze the spectra under technical settings of case 3, where basically all
parameters in Table 9 and the blackbody temperature are free. As in Section
3.3, errors on the parameters are calculated using the Goodman-Weare MCMC
algorithm in XSPEC. Following Choudhury et al. (2017), we set a Gaussian
proposal distribution initialized at $\mathbf{X}+3\mathbf{\sigma_{X}}$, where
$\mathbf{X}$ is the input model parameter and $\mathbf{\sigma_{X}}$ is the
standard deviation estimated from the covariance matrix of a preliminary fit.
For boundary inputs $\theta_{\mathrm{disk}}=[3^{\circ},87^{\circ}]$ and
$R_{f}=[0,10]$, the covariance information is changed to a diagonal matrix
constructed from 100 times the step length of the parameter. For each group of
the mock spectra, we use 200 walkers for the chain and sample to a total of
$1.2\times 10^{6}$ elements (6000 iterations) with the first $1.2\times
10^{5}$ rejected. As an example of the mock and error calculation, we show a
group of spectra in Figure 5 and its corresponding MCMC corner in Figure 6.

We check the parameter recovery of the mock spectra by seeing whether the
“measured” (fitted) values are close to the input values, and whether the
error range includes the input. Figure 7 compares the measured versus input
values of inclination, with colored regions depicting 90% confidence
intervals: the dispersion suggests that the input values are generally
recoverable. There is little dependence on the iron abundance or coronal
temperature. Since our simulation grid is five-dimensional, the figures in the
following sections are provided to show the trend of parameter recovery
regardless of the specific parameter choice in one or two dimensions.
Additionally, larger reflection fraction and inclination both lead to tighter
parameter constraints, as to be expected because either large $R_{f}$ or large
$\theta_{\mathrm{disk}}$ increases the relativistic smearing of the intrinsic
spectra. Dauser et al. (2014) find a similar trend for the reflection model
based on the lamp-post geometry (Matt et al., 1991; Martocchia & Matt, 1996;
Reynolds & Begelman, 1997; Miniutti & Fabian, 2004). The confidence interval
shrinks most markedly when $R_{f}$ varies from 0.3 to 1, with less dramatic
but still notable improvements for $R_{f}=5-10$. Not surprisingly, the
completely reflection-free case ($R_{f}=0$) can barely measure
$\theta_{\mathrm{disk}}$; the posterior distribution (blue region) is very
close to a uniform distribution.

Figure 8: The dependence of inner disk inclination $\theta_{\mathrm{disk}}$ on
input iron abundance $A_{e}$ and reflection fraction $R_{f}$, for fixed
$kT_{e}=30\rm\,keV$, $\theta_{\mathrm{disk}}=45^{\circ}$, and
$a_{\ast}=0.998$. Panel (a) shows the measurements of $\theta_{\mathrm{disk}}$
and the 90% confidence intervals, and panel (b) the absolute offset in
$\theta_{\mathrm{disk}}$. Figure 9: The dependence of the relative offset in
inner disk inclination ($\sigma$ refers to the average statistical error for
each point in each panel) on coronal temperature $kT_{e}$ and input (a)
reflection fraction $R_{f}$ (with solar $A_{\mathrm{Fe}}$) and (b) iron
abundance $A_{e}$ (with $R_{f}=5$), for fixed
$\theta_{\mathrm{disk}}=45^{\circ}$ and $a_{\ast}=0.998$.

##  5 Discussion

###  5.1 Sensitivity to Reflection Fraction

Envisioning the corona as a razor-thin layer that hugs the entire disk, the
RelxillCp model defines the reflection fraction as

$R_{f}\coloneqq\frac{f_{\mathrm{disk}}}{f_{\mathrm{\infty}}},$ (1)

where $f_{\mathrm{disk}}$ denotes the fraction of photons hitting the disk and
$f_{\mathrm{\infty}}$ the fraction of photons escaping to infinity (Dauser et
al., 2016) in the rest-frame of the corona. Given that the corona precisely
overlays the disk with a large covering fraction, we anticipate the reflection
fraction to be $\sim 1$. Therefore, in our spectral fits, $R_{f}$ is not
intended to serve as an indicator of the geometry or relativistic light-
bending; instead, it serves as an empirical representation of the relative
strength of the reflection component, encompassing the soft excess, the Fe
K$\alpha$ feature, and the Compton hump. Although this parameter lacks
physical significance, it is nevertheless important because it provides
insight into the proportion of the reflected emission
($\frac{R_{f}}{1+R_{f}}$). This, in turn, helps us to estimate the systematic
error and lends confidence in applying the reflection model.

Referring back to the mock results, in a typical case of $R_{f}=1$ (e.g.,
Figure 7), the relative error can be as large as 50% to 90% when
$\theta_{\mathrm{disk}}$ decreases to $\lesssim 60^{\circ}$. However, for
larger $\theta_{\mathrm{disk}}$, the error is $\lesssim 10^{\circ}$. When
$R_{f}\gtrsim 5$, the proportion of reflection is already larger than 83%, and
further increasing $R_{f}$ does not significantly raise this proportion,
resulting in negligible changes to the confidence interval. On the other hand,
for $R_{f}\lesssim 0.3$, where the reflection signal is weak ($\lesssim
23\%$), robust inclination measurements are challenging to achieve. The impact
of reflection fraction is also seen in Section 3.3, where the results of
$\theta_{\mathrm{disk}}$ from different simultaneous epochs of I Zwicky 1 are
more consistent compared to 3C 382. Indeed, even a single-epoch XMM-Newton
observation of I Zwicky 1 offers a fair constraint on inclination. Without
loss of generality, we define a spectrum as reflection-dominant when
$R_{f}\gtrsim 5$ and reflection-subordinate when $R_{f}\lesssim 0.3$. In
general, it is usually easier to achieve tighter constraints on inner disk
inclination for sources with higher $R_{f}$, and, of course, larger
inclination, because the observed spectrum of more edge-on systems naturally
has a stronger reflection signal. However, it is also worth noting that
combining more groups of data can counterbalance the poor constraints of
reflection-subordinate sources, as examplified by the case of 3C 382.

###  5.2 Sensitivity to Iron Abundance

Figure 8 examines the effect of iron abundance on the measurement of inner
disk inclination. Large $A_{\mathrm{Fe}}$ increases the strength of the whole
Fe profile (García et al., 2013). In our simulations, both the absolute value
of inclination offset and the confidence interval of $\theta_{\mathrm{disk}}$
generally decrease with a larger input value of $A_{\mathrm{Fe}}$. Consistent
with the results of previous sections, this trend holds as well for smaller
$R_{f}$, although the constraints on $\theta_{\mathrm{disk}}$ are
significantly looser because reflection spectroscopy is intrinsically less
reliable for reflection-subordinate cases. In contrast with $R_{f}$, the
influence of $A_{\mathrm{Fe}}$ on the accuracy of spectrum fitting is not as
significant, as evidenced by comparing Figure 7 and Figure 8b, or by comparing
the effect of $R_{f}$ and $A_{\mathrm{Fe}}$ in Figure 8 itself. On the one
hand, the inclination results are rather reliable already for sources with
solar or mildly super-solar abundance; on the other hand, it implies that the
model is not that sensitive to $A_{\mathrm{Fe}}$ compared to $R_{f}$. This,
aside from the support for super-solar abundance in AGN reflection modeling
(Fabian, 2006), also explains why the broadband fitting for 3C 382 returns
such a huge value of $A_{\mathrm{Fe}}$ with a large error even though we
simultaneously fit five epochs. We cannot use $A_{\mathrm{Fe}}$ as a basis for
judging the robustness of reflection models.

###  5.3 Sensitivity to Coronal Temperature

Figure 9 documents the effect of coronal temperature on the measurement of
inner disk inclination. The inclination offset in terms of error generally
decreases with larger input values of $kT_{e}$. Consistent with previous
sections, this trend holds regardless of the values of $R_{f}$ and
$A_{\mathrm{Fe}}$, although the effect from $R_{f}$ is greater than that from
$A_{\mathrm{Fe}}$. Similar to the case of $A_{\mathrm{Fe}}$, the inclinations
are fairly reliable for sources with $kT_{e}=30\,\rm keV$. In all cases,
however, $|\theta_{\mathrm{disk}}|/\sigma<1$, meaning that deviations between
the measured and input inclinations all lie within the $1\sigma$ confidence.

Figure 10: Spin measurements for simulated spectra with different values of
(a) reflection fraction $R_{f}$ and fixed $\theta_{\mathrm{disk}}=45^{\circ}$,
and (b) different values of $\theta_{\mathrm{disk}}$ and fixed $R_{f}=1$. In
all panels, we fix $A_{\mathrm{Fe}}=1$ and $kT_{e}=30\,\rm keV$. Figure 11:
Panel (a) shows the dependence of the measurements of $\theta_{\mathrm{disk}}$
on spin $a_{\ast}$ with $R_{f}=1$; the colored regions depict the 90%
confidence intervals and the dashed line the expected values. Panel (b)
illustrates the relative offset in $\theta_{\mathrm{disk}}$ with respect to
different $a_{\ast}$ and $R_{f}$. In all panels, we fix $A_{\mathrm{Fe}}=1$
and $kT_{e}=30\,\rm keV$.

###  5.4 Spin Degeneracy

Both the spin and the inner disk inclination contribute to the broadening of
the Fe K$\alpha$ emission line, with the red edge of the profile being more
sensitive to $a_{\ast}$ and the blue edge more sensitive to
$\theta_{\mathrm{disk}}$. Different combinations of $a_{\ast}$ and
$\theta_{\mathrm{disk}}$ can yield similar Fe profiles, making them
challenging to distinguish during spectral fitting. This inherent degeneracy
is further complicated by the coupling of both parameters with the entire
model. Given the intricate interplay between $a_{\ast}$,
$\theta_{\mathrm{disk}}$, and other model parameters, our mock tests are
conducted with specific considerations for two corona-related parameters
($\xi$ and $\Gamma$), disk geometry, and, once again, the physical background
of $R_{f}$:

  1. 1.

While $\xi$ is better measured with higher values, incorrect measurements of
$\xi$, followed by $A_{\mathrm{Fe}}$, $\Gamma$, and $\theta_{\mathrm{disk}}$,
pose challenges for spin retrieval (Choudhury et al., 2017). We use
$\log\xi=3$ as the input and find that the 90% error of $\xi$ is $\lesssim
300\,\rm erg~{}cm~{}s^{-1}$ for all mock spectra.

  2. 2.

According to Choudhury et al. (2017), underestimating $\xi$ corresponds to
overestimating $\Gamma$, and vice versa. Lower $\xi$ shifts the spectra toward
harder X-rays, as larger $\Gamma$ produces less soft X-ray emission. Here,
$\xi$ influences recombination and photoionization strength, directly shaping
the recombination continuum. We generate spectra with $\xi=1000\,\rm
erg~{}cm~{}s^{-1}$ as an input value that can be acceptably retrieved,
reducing its influence on $a_{\ast}$ and the overall model fitting. As the
power-law exponent in the spectral model, $\Gamma$ should be easier to
independently constrain in our broadband fits.

  3. 3.

The value of $a_{\ast}$ sets a physical lower limit for the inner radius of
the accretion disk with a negative correlation, which is the main effect of
the spin in observations (Dauser et al., 2010). Besides, changing $a_{\ast}$
modifies the geodesics in ray tracing, which has a more complicated influence
on observed spectra. For a disk geometry, $a_{\ast}=0.998$ supports a disk
reaching $\sim 1.24\,r_{g}$, while $a_{\ast}=-0.998$ pushes the innermost
stable circular orbit to $9\,r_{g}$. We set the broken radius to $10\,r_{g}$
to avoid physically impossible circumstances of
$R_{\mathrm{br}}<R_{\mathrm{in}}$.

  4. 4.

The $R_{\mathrm{in}}-a_{\ast}$ relation also affects the reflection fraction.
The maximal possible $R_{f}$ declines for smaller $a_{\ast}$ or larger
$R_{\mathrm{in}}$ (Dauser et al., 2014). Although our phenomenological $R_{f}$
is not directly constrained by $a_{\ast}$, some extreme combinations in our
simulation grid are rare in real observations.

  5. 5.

Aside from the spin, measured $R_{f}$ is further affected by
$\theta_{\mathrm{disk}}$ because the observer only receives flux traveling at
a certain orientation with respect to the system. According to Suebsuwong et
al. (2006), extremely large $R_{f}$ ($\gtrsim 10$) is possible only for
systems with high inclinations ($\theta_{\mathrm{disk}}\gtrsim 70^{\circ}$).
Still, since $R_{f}$ in RelxillCp is not physical, a few extreme combinations
are possible in mocks but do not necessarily exist. The only motivation for
including those points in the parameter space is to make the discussion more
complete computationally.

In conclusion, we provide rationale for selecting the specific values for
$\xi$, $\Gamma$, and $R_{\mathrm{br}}$ summarized in Table 9. Additionally,
our previous analyses indicate that $A_{\mathrm{Fe}}$ and $kT_{e}$ have a
lesser impact on $\theta_{\mathrm{disk}}$ derivation compared to $R_{f}$.
Consequently, the ensuing discussions will focus on the interplay between spin
and inclination measurements and the influence of the reflection fraction.

Like the case for inclination, the measurement of $a_{\ast}$ is influenced by
its intrinsic value. All three preceding studies on mock spectra of RELXILL
focused on the reliability of spin measurements and observed a more accurate
recovery of $a_{\ast}$ with higher $R_{f}$ and larger $a_{\ast}$. This general
trend is confirmed in our analysis, as depicted Figure 10a. A larger
$a_{\ast}$ allows the inner edge of the disk to approach the BH more closely,
resulting in more pronounced redshifting of photons and the formation of a
broader Fe line that is easier to detect. Figure 10b examines the effect of
inclination on spin measurement, with colored regions depicting 90% confidence
intervals. Generally, larger inclinations lead to less constrained spins
because the differences in line broadening caused by $a_{\ast}$ are more
evident at smaller $\theta_{\mathrm{disk}}$ (e.g., see Dauser et al., 2010).
The inclination measurements are affected minimally by the spin (Figures 11a
and 11b).

Considering the entire simulation grid of 13,860 spectra, 85% exhibit an
offset between the input and measured inclination smaller than $10^{\circ}$,
77% have an inclination offset smaller than $5^{\circ}$, and 62% have an
absolute inclination offset smaller than $5^{\circ}$, with their 90%
confidence interval covering the input inclination. All reduced fit statistics
are below 1.1.

###  5.5 Comparison with Previous Work

In this section, we compare both observed and simulated spectra analyses with
previous work. Overall, eight X-ray studies have been reported for I Zwicky 1
since the advent of XMM-Newton. For 3C 382, four works have been reported with
XMM-Newton observations. Three mock analyses have been published for the
RELXILL family, all focusing on the recovery of $a_{\ast}$.

For I Zwicky 1, Gallo et al. (2004) discovered a broad iron line, a hard X-ray
flare, and spectral hardening during the flare with the $\sim 20\,\rm ks$ XMM-
Newton observation taken in 2002. Costantini et al. (2007) discovered two WAs
in soft X-rays from the RGS data, which was further confirmed by Silva et al.
(2018) from their two XMM-Newton observations in 2015. Based on the same
observations in our work, Wilkins et al. (2021) and Wilkins et al. (2022)
suggested that both the broad Fe K$\alpha$ line extending to 3 keV and the
Compton hump around 25 keV were well explained with the disk reflection model.
Using a lamp-post coronal geometry, an ionization gradient, and a twice-broken
power-law emissivity, they obtained a corona height of $h>14\,r_{g}$ during
the flares and $h\approx 3.7\,r_{g}$ before and after flares. Ding et al.
(2022) continued with the idea of reflection and showed that a single broken
power-law emissivity implemented by RelxillCp was sufficient, and that an
ionization gradient was unnecessary once the absorbers were added correctly
(see their Figure 2 and Section 3.2). While most of the reflection parameters
($\Gamma$, $kT_{e}$, $R_{f}$) indicated the same evolutionary trend as in
Wilkins et al. (2022), their result for the spin ($a_{\ast}>0.973$) is
consistent with that ($a_{\ast}>0.94$) of Wilkins et al. (2022); their
inclination angle $\theta_{\mathrm{disk}}\approx 42.3^{\circ}$ is $\sim
10^{\circ}$ smaller; and their $R_{f}\gtrsim 1$ is an order of magnitude
bigger, which was argued to be the consequence of applying a different ray-
tracing kernel. Modeling the reflection continuum and ionized absorbers as in
Ding et al. (2022) yields generally similar results. However, we provide
measurements for Index2 ($\sim 3$) while they fixed it to 3. Including
ultraviolet data from XMM-Newton observations, Rogantini et al. (2022) fit the
spectral energy distribution via both a frequentist and a Bayesian method and
further acknowledged the components found in previous works. With refl in SPEX
(Kaastra et al., 1996) assuming solar metal abundance and zero ionization,
they reported a reflection scale factor (equivalent to $R_{f}$) of $0.43$.

Reeves & Braito (2019) analyzed the 2015 XMM-Newton observations (3–10 keV)
with four different models consecutively: a P Cygni profile (Done et al.,
2007), a photoionization continuum, a reflection model (XILLVER$\ast$KDBLUR),
and a fast ($\sim c$), wide-angle disk wind model (Sim et al., 2008). They
claimed that the broad Fe K profile and the deep 9 keV absorption can be
described by a disk wind, although the reflection component could not be
excluded. Their reflection modeling returned an ionization $\xi\approx
10^{3}$, a reflection fraction $\sim 0.6$, a flat single power-law emissivity
with an index around 2.2, and an inclination angle even 2–20 degrees larger
than that of Wilkins et al. (2022). However, their spectra excluded the hard
X-rays, and their model did not have spin as a parameter due to the use of
another ray-tracing kernel. On modeling the outflowing wind, they set the wind
location to the escape radius of the BH and obtained an accretion disk
inclination angle $49.7^{\circ}$ as an indicator of blueshifted absorption
features.

For 3C 382, Torresi et al. (2010) resolved a WA in the RGS spectra from the
2008 observation, which was the first observed WA in a broad-line radio
galaxy. Sambruna et al. (2011) analyzed the EPIC spectrum from the same
observation using parameters derived from spectral fits of joint Suzaku and
Swift/BAT data and discovered both a broad and narrow Fe K$\alpha$ line, a
hump above 10 keV, a soft excess, and a WA. They argued for a two-ionization
reflection with a highly ionized component reproducing the broad emission line
and a mildly ionized component accounting for the narrow emission line, the
hump, and the soft excess. Both reflection components contributed to $\sim
10\%$ of the total continuum. Using two XILLVER models, they reported
$\log\xi=1.54$ and $2.93$. Due to low S/N, they could not constrain the spin
with both RELCONV and RELLINE. However, for an assumed emissivity index of
2.5, the inner and outer radii of the reflection region were around
$R_{\mathrm{in}}\approx 10\,r_{g}$ and $r_{\mathrm{out}}\approx 25\,r_{g}$,
respectively, and the inner accretion disk inclination is around 30 degrees.
With the simultaneous NuSTAR and Swift observations in 2012 and the NuSTAR
observation in 2013, Ballantyne et al. (2014) found a weak Fe line in the
NuSTAR spectra and argued that it originated from the outer regions ($\gtrsim
50\,r_{g}$) of an outflowing corona and ionized accretion disk. In their
attempts to utilize RELXILL, they reported an upper limit for the reflection
fraction ($R_{f}\lesssim 0.1$) and a lower limit for the inclination angle
($\theta_{\mathrm{disk}}>25^{\circ}$), which was insensitive to the data.
Based on the same observations in our work, Ursini et al. (2018) jointly fit
the PN and NuSTAR spectra and argued that the line broadening can be explained
by instrumental effects or reflection on the torus or the broad-line region.
With RELXILL, they reported a very small reflection fraction consistent with
that estimated by Ballantyne et al. (2014), and yet an upper limit of
$15^{\circ}$ for the inclination angle. They also calculated the intrinsic
line width for all five epochs: $\sigma<0.2$, $<0.35$, $0.5^{+0.3}_{-0.2}$,
$0.13^{+0.09}_{-0.08}$, and $<0.2$, respectively.

Bonson & Gallo (2016) generated over 4000 high-quality mock Seyfert 1 spectra
to study spin recovery with RELXILL, assuming reflection fractions of 1 and 5.
They discovered that most parameters were overestimated with spectral coverage
limited to 2.5–10.0 keV and that extending the data to 70 keV improved the
measurements. However, the model was found insensitive to Index1, and
$a_{\ast}$ was well constrained only for large values of $a_{\ast}$ (e.g.,
with error $\pm 0.10$ for input $a_{\ast}>0.8$). Choudhury et al. (2017)
tested the recovery of $a_{\ast}$ and $\xi$ in the lamp-post flavor of
RELXILL–RelxillLp–with 5800 NuSTAR simulations of a bright source. They found
that both parameters are well-recovered at 90% confidence with improving
constraints at higher $R_{f}$ and $a_{\ast}$, and lower source height $h$.
Besides, they stressed the importance of choosing a reasonable initial point
for spectral fitting and spanning the parameter space carefully for nonlinear-
behaving parameters (e.g., $\xi$) to avoid the fit being trapped in a local
minimum. Kammoun et al. (2018) generated 60 mock spectra with more complex
components including warm and neutral absorbers, relativistic (RelxillLp) and
distant reflection, and thermal emission. Their blind fit was produced
starting with selecting models, which emulated the typical procedure of
observed spectral fitting. They reported that neither the absorbers nor the
input $a_{\ast}$ affected the fits significantly, although $h$ was important.
Compared with our mock analysis, however, both Bonson & Gallo (2016) and
Choudhury et al. (2017) excluded the soft band. Tests in Kammoun et al.
(2018), although more complete in spectral coverage, included a plethora of
spectral models beyond our main focus. More importantly, all three works
generate mock spectra with a long exposure time to show the upper limitations
of model performance with high counts and high-S/N data. By contrast, our
analysis tests the model performance with data settings more closely related
to real observations. In addition, Choudhury et al. (2017) and Kammoun et al.
(2018) adopted the lamp-post geometry, whereas we do not assume any particular
coronal configuration.

###  5.6 Model Selection

In practice, selecting a model to interpret the X-ray observations of AGNs is
challenging because their typically complicated spectra can be fit by several
candidate models. The comparisons in the previous section clearly illustrate
the statistical and potential systematic uncertainties introduced during model
selection (aspect 3 in Section 2.2). For I Zwicky 1, the presence of obvious
line broadening and a high-energy hump may seem to provide exclusive evidence
for the reflection model, but the wind model cannot be fully excluded. In the
case of 3C 382, a weak reflection signal leaves ample room for alternative
models, and the results of reflection models are significantly shifted based
on data quality and quantity (Section 3.4), as well as the strategy of
parameter selection in fitting. Even within the framework of reflection, the
choice of the emissivity profile associated with the coronal model also
affects the measurement of inclination and other parameters. Therefore,
current reflection modeling cannot circumvent case-by-case studies to fit
observed spectra. For reflection-subordinate sources, combining multi-epoch
broadband spectra and carefully exploring the parameter space becomes crucial,
as it can at least minimize the statistical error.

Within the framework of reflection, we choose RelxillCp over RELXILL or
RelxillD due to its thermal Comptonization mechanism of the corona. The main
difference among the models in the RELXILL family lies in the configuration of
the corona, which might influence inclination measurements because it
determines the radial emissivity profile of the system. In practice, the
radial emissivity profile can be parameterized for various geometries (e.g.,
Dauser et al., 2013). For RelxillCp, the emissivity is given as
$r^{-\mathrm{Index}_{1}}$ between $R_{\mathrm{in}}$ and $R_{\mathrm{br}}$ and
$r^{-\mathrm{Index}_{2}}$ between $R_{\mathrm{br}}$ and $R_{\mathrm{out}}$.
The fact that these empirical parameterizations often find a high Index1
($\gtrsim 5$; see also Section 3.3) has motivated the lamp-post geometry–an
isotropically irradiating source located in the vicinity and on the rotational
axis of the BH–which has been widely employed in recent years to interpret the
X-ray spectra of AGNs (e.g., Wilkins & Fabian, 2011; Dauser et al., 2012;
Miller et al., 2015; Beuchert et al., 2017). This leads to the advent of the
RelxillLp series. We stress, however, that the configuration of the corona
remains unclear and difficult to distinguish from spectral analysis (e.g.,
Tortosa et al., 2018). Further determination of coronal geometry usually
requires additional information, for instance from X-ray reverberation mapping
(e.g., Reynolds et al., 1999; Fabian et al., 2009; Kara et al., 2016;
Caballero-García et al., 2020) or X-ray spectropolarimetry (e.g., Tagliacozzo
et al., 2023; but see Ingram et al., 2023). Comparing well-understood,
theoretical emissivity profiles to the profile determined from the Fe line,
and applying constraints from reverberation lags, can help establish the
geometry (Wilkins & Fabian, 2012), but this is beyond the scope of this work.

Table 10: Comparison of Best-fit Results of I Zwicky 1 Using RelxillCp and
RelxillLpCp

Flavor |  $\theta_{\mathrm{disk}}$ (∘) | $a_{\ast}$ |  $C$-stat | DOF | $C_{\mathrm{red}}$ | $k$ | AIC
---|---|---|---|---|---|---|---
(1) | (2) | (3) | (4) | (5) | (6) | (7) | (8)
-Cp | $44.3_{-2.4}^{+2.8}$ | $>0.978$ | 926 | 829 | 1.117 | 64 | 1054
-LpCp | $44.1_{-6.4}^{+8.2}$ | $>0.964$ | 1010 | 824 | 1.226 | 69 | 1148

Note. — Col. (1): Flavor of the model used. Col. (2): Inclination angle of the
inner accretion disk. Col. (3): Dimensionless BH spin. Col. (4): Cash
statistics. Col. (5): Degree of freedom. Col. (6): Reduced $C$-stat. Col. (7):
Number of parameters. Col. (8): Akaike Information Criterion.

Aside from the arguments in coronal geometry, RelxillLpCp includes another
feature compared to RelxillCp: the ionization gradient. Photoionized by the
irradiation from the corona, the accretion disk receives a flux as a
decreasing function of the radius, leading to a radial distribution of the
ionization parameter,

$\xi\coloneqq\frac{4\pi F(r)}{n},$ (2)

where $n$ is the number density of the plasma receiving an ionizing flux $F$.
In RELXILL, the parameter $n$ (or $n_{\mathrm{H}}$) usually denotes the
hydrogen number density, following XSTAR. However, the $n$ in Equation (2)
should be the electron number density, equal to 1.2 times the hydrogen density
since most H atoms are ionized, and $F$ becomes the net integrated flux in the
1–1000 Ry energy range ($F_{\rm X}$; Tarter et al., 1969; García et al.,
2013). Thus, the ionization parameter in RelxillCp is $\xi=4\pi F_{\rm
X}/1.2n_{\mathrm{H}}$. Svoboda et al. (2012) argued that the radial ionization
profile should be considered to derive the very steep emissivity profiles
(Index${}_{1}\geq 7$), which is degenerate with both the spin and the
inclination.

For 3C 382, it is pointless to debate whether we need to include additional
features because reflection modeling suffers from a low reflection fraction.
But it is worthwhile to see if the ionization gradient affects the inclination
angle for typical sources like I Zwicky 1. We conduct model comparison via the
likelihood ratio ($\Lambda$) and the Akaike (1974) information criterion
(AIC). Because we use $C$-stat for spectral fitting, the likelihood-ratio test
is equivalent to the $\chi^{2}$ hypothesis test (Wilks, 1938, 1963),

$-2\log\Lambda=\Delta C\to\chi^{2}_{\Delta k},$ (3)

where $k$ is the number of model parameters. We apply AIC to compare the
ability of models to properly describe the observed spectrum without
significant over-fitting and under-fitting. The AIC is defined by

$\mathrm{AIC}\coloneqq-2\log\mathcal{L}+2k=C+2k,$ (4)

which not only includes the likelihood function but also punishes the model
with more parameters. We calculate AIC for our fit of I Zwicky 1 using
RelxillCp as the null model, and we produce alternative model fitting using
RelxillLpCp with coronal parameters initialized according to the best-fit
values in Wilkins et al. (2022). During the fit, the density of the disk and
the ionization gradient are considered as corona-related parameters, meaning
that they are free to vary across spectral epochs. The radial density profile
is not considered because its effect is incomparable with that of irradiation.
We report the inclination, spin, and other values related to statistical tests
in Table 10. As the inclination and spin measurements are not significantly
affected, the reduced $C$-stat for RelxillLpCp is slightly larger than for
RelxillCp. A likelihood-ratio test cannot reject the null model, and the
alternative model is disfavored by $\Delta\mathrm{AIC}=+94$. Therefore, the
ionization gradient is unnecessary for I Zwicky 1, which is consistent with
Ding et al. (2022). In addition, returning radiation is also included in
RelxillLpCp (Dauser et al., 2022), which is tested jointly. Without a notable
effect on $\theta_{\mathrm{disk}}$, we have observed stronger reflection
($R_{f}>0.5$ to $>4.3$) and a lower coronal height ($h\lesssim 5\,r_{g}$ to
$\lesssim 2\,r_{g}$) comparing to Wilkins et al. (2022) ($R_{f}<0.5$;
$h>14\,r_{g}$ to $\gtrsim 3\,r_{g}$), even though all parameters evolve in the
same trend.

For generality and our purpose of inclination measurements, it is readily
acceptable to adopt RelxillCp to use its two emissivity indices (Index1 and
Index2) and the broken radius ($R_{\mathrm{br}}$), which implicitly has taken
into consideration various coronal geometries without further assumptions.
However, it is far more difficult to interpret than simply applying a lamp-
post geometry with which many of the model parameters are directly calculated
(e.g., the relation between $R_{f}$ and $\theta_{\mathrm{disk}}$ is well-
studied in Dauser et al., 2014). Moreover, the model is insensitive to the
inner emissivity Index1 (Bonson & Gallo, 2016), rendering it challenging to be
constrained in a blind fit. We suggest assigning proper initial values for the
emissivity parameters and always providing the fitted emissivity along with
the inclination measurements for reproducibility (as in Section 3.3).

Table 11: Comparison of Best-fit Results of 3C 382 with Different Treatments
of the Soft Band

Fit |  $\theta_{\mathrm{disk}}$ (∘) | $a_{\ast}$ |  $C$-stat | DOF | $C_{\mathrm{red}}$ | $k$ | AIC
---|---|---|---|---|---|---|---
(1) | (2) | (3) | (4) | (5) | (6) | (7) | (8)
I | $41.8_{-4.1}^{+2.7}$ | $>0.992$ | 5682.7 | 5004 | 1.136 | 68 | 5819
II | $41.8_{-4.0}^{+3.2}$ | $>0.992$ | 5682.7 | 4996 | 1.137 | 72 | 5827
III | — | — | — | — | — | — | —
IV | $33.2_{-6.5}^{+6.7}$ | $>0.996$ | 4677.5 | 4238 | 1.104 | — | —

Note. — Col. (1): Spectral fit. Col. (2): Inclination angle of the inner
accretion disk. Col. (3): Dimensionless BH spin. Col. (4): Cash statistics.
Col. (5): Degree of freedom. Col. (6): Reduced $C$-stat. Col. (7): Number of
parameters. Col. (8): Akaike Information Criterion. In view of Fit III
yielding no results, a comparative analysis of the AIC between Fit III and IV
is omitted.

###  5.7 Effect of Soft Excess Interpretation

The interpretation of the soft excess has been a subject of ongoing debate
(e.g., Crummy et al., 2006; Done et al., 2012; Jin et al., 2017; García et
al., 2019; Gliozzi & Williams, 2020). The soft excess is, at least in part, a
reflection feature that influences the determination of parameters such as
$a_{\ast}$ (e.g., Walton et al., 2013), ionization, and the emissivity
profile. However, the soft band often introduces complications with complex
absorption and scattering phenomena, mimicking the relativistic effects
associated with reflection. Consequently, distinguishing whether the soft
excess is due to reflection-based relativistic smearing or purely reprocessing
poses a challenge. This ambiguity is often addressed using high-S/N broadband
spectra (e.g., Guainazzi et al., 2006; Risaliti et al., 2013; Kammoun et al.,
2018). In Ding et al. (2022) and our sections on observed spectra, the utility
of multi-epoch, broadband, reflection-dominant spectra with well-modeled
absorption reduces the complexities.

Using 3C 382 as an example, we then delve into a detailed discussion of the
phenomenological modeling of the soft excess with little influence of the WA
since WAs are well-modeled in Section 3.3. As the RelxillCp model inadequately
describes the soft excess in 3C 382, we choose to model it using an additional
phenomenological redshifted blackbody with a temperature tied among spectral
epochs ($kT=8.7_{-0.3}^{+0.3}\times 10^{-2}\,\rm keV$). Despite statistical
robustness, this introduces systematic uncertainty into the spin, ionization,
and emissivity profile. The blackbody temperature may vary between epochs, and
alternatively, adopting the zBBody model may introduce a selection bias. We
conduct four spectral fits to elucidate the effect of our zBBody
interpretation of the soft excess:

  1. I.

The spectral fitting in Section 3.

  2. II.

Untying blackbody temperatures between epochs.

  3. III.

Ignoring spectral bins with energies lower than 2 keV and directly producing
the fit to check for the existence of a selection bias.

  4. IV.

Ignoring spectral bins with energies lower than 2 keV and deleting the
absorber table components in the soft band. We then directly start an MCMC
error calculation with the best-fit parameters from our original fit (Fit I).

Results are presented in Table 11. The inclination $\theta_{\mathrm{disk}}$
and spin $a_{\ast}$ show no significant impact in Fit II. Despite introducing
four additional model parameters to account for the temperature, the $C$-stat
remains unchanged ($\Delta C\text{-stat}=0.02$). The likelihood-ratio test
cannot reject the null model (Fit I), and this scenario is slightly disfavored
with $\Delta\mathrm{AIC}=+8$. The blackbody temperatures across the five
epochs align with the tied value in Fit I
($|kT_{\mathrm{a,b,c,d,e}}-0.087|<0.002$). Fit III yields unattainable results
when information in the soft band is completely lost. The values for
inclination ($\theta_{\mathrm{disk}}\leq 3^{\circ}$), spin ($a_{\ast}=0.998$),
Fe abundance ($A_{\mathrm{Fe}}=10$), and most emissivity indices (7 out of 10
Indexes(1,2) $\to 10,0$) are all pegged to their boundaries. The reflection
fraction ($R_{f}\lesssim 0.29$) is notably smaller than the values in Table 4.
Therefore, addressing the bias of using a specific soft excess model is
challenging through a comparison of fits with broadband data and pure hard
X-ray data. In Fit IV, we abandon pursuing a fit and instead start error
calculation on the known ”best-fit” parameters directly. As expected,
$a_{\ast}$, $A_{\mathrm{Fe}}$, and $\Gamma$ generally align with Fit I,
exhibiting overlapping 90% error. However, $R_{f}$ is smaller ($<0.25$), and
ionization is larger ($\log\xi\gtrsim 3.2$). The larger $\xi$ possibly shifts
$\theta_{\mathrm{disk}}$ to a smaller value. Notably, RELXILL still fits the
harder energies well in this scenario.

In summary, the epoch-independent blackbody effectively models the soft
excess. We emphasize the importance of incorporating information from the soft
band to enhance the precision of the reflection model fit.

###  5.8 Data Rebinning

Finally, we address a crucial technical aspect that can introduce biases in
spectrum analysis—data rebinning. This process involves redividing and
recombining data into a new set of bins. In spectroscopy, the decision to
rebin over energy channels occurs after data reduction. While this reduces raw
spectrum information, introduces potential errors, and increases the risk of
creating artifacts, it is scientifically and statistically justified.

Handling broadband data from diverse instruments necessitates rebinning when
assigning equal weights to all instruments. Analyzing spectra with
significantly different photon counts, such as combining optical and X-ray
spectra or XMM-Newton and NuSTAR spectra with varying X-ray emissions, can be
challenging. High-count spectra would dominate the fit if the raw data were
used without proper rebinning555Combining spectra from different instruments
is generally discouraged due to differences in response matrices, which may
introduce artifacts. Nevertheless, specific science applications may still
benefit from this approach (e.g., Wilkins et al. 2021; Ding et al. 2022; see
our data deduction for I Zwicky 1 in Section 3.2).. Moreover, even a single
spectrum from one instrument may lose important information, like emission
lines, if data bins are too narrow (e.g., Kaastra, 1999; Kaastra & Bleeker,
2016).

For statistical convenience, a common practice in X-ray observations is to
rebin spectra so that each energy bin contains a minimum of 25 counts. This
ensures approximately Gaussian data, facilitating interpretation with standard
$\chi^{2}$ statistics. However, blindly adopting this approximation introduces
bias in fitting (e.g., Humphrey et al., 2009; Choudhury et al., 2017). Poisson
likelihood and Cash statistics are available for unbinned spectra in event-
based X-ray observations. Our analysis of mock spectra (Section 4.2) uses Cash
statistics. The debate on statistical choices originates from frequentist
statistics, and introducing Bayesian approaches to X-ray spectrum analysis can
alleviate it. The nested sampling method (Skilling, 2004), widely embraced in
cosmology (e.g., Planck Collaboration et al., 2020) and gravitational wave
studies (e.g., Abbott et al., 2019), could be beneficial (see, e.g., Buchner
et al., 2014). This method offers advantages like model selection with the
Bayes factor and robust parameter estimation. The nested sampling method is
particularly useful for handling multiple modes/maxima in the posterior
parameter distribution of reflection models, challenging for typical MCMC
algorithms.

Rebinning introduces computational challenges, especially for instruments with
high resolution or when the forward-folding process is repeated frequently in
spectral fitting with many free parameters (see detailed discussion in Section
2 of Kaastra & Bleeker, 2016). This problem is magnified in upcoming high-
resolution X-ray missions like XRISM (Tashiro et al., 2018) and Athena (Nandra
et al., 2013). While not rebinning the simulated spectra might extend MCMC
analysis time, it helps avoid potential biases introduced during data
reduction. Given that the simulations only involve a broad Fe line and lack
complex emission lines, we believe oversampling does not compromise important
information.

In summary, the choice of data binning involves a compromise between
scientific and computational considerations. Balancing this tension often
means rebining spectra to avoid oversampling the instrumental resolution by a
factor of 3 (1.5 times the Nyquist rate; see Nyquist 1928; Shannon 1949). The
application of systematic approaches like that proposed by Kaastra & Bleeker
(2016), considering both S/N and instrumental resolution, needs validation
through tests on mock spectra to confirm its ability to reduce computational
time while preserving critical spectral properties, especially complex
reprocessing features. This remains a topic for future investigation.

##  6 Conclusions

We propose a systematic method to measure the inner accretion disk inclination
of AGNs with broadband ($0.3-78\,\rm keV$) X-ray reflection spectroscopy.
Using RelxillCp from the self-consistent reflection model family RELXILL,
without assuming any particular coronal geometry, we test the reliability of
measuring inclination by analyzing joint XMM-Newton and NuSTAR spectra of the
narrow-line Seyfert 1 galaxy I Zwicky 1 and the broad-line radio galaxy 3C
382, alongside 13,860 simulated spectra that mimic joint XMM-Newton and NuSTAR
observations of 3C 382.

Statistically consistent measurement of inclination can be achieved with
RelxillCp. More accurate measurements can be obtained if multiple observations
are combined and fitted simultaneously, and if the spectrum is more
reflection-dominant. Multi-epoch, simultaneous data offer the best constraints
on inner accretion disk inclination, but as long as prominent reflection
features exist within the relevant instrumental coverage, inclination
measurements can still be obtained, even for a single-epoch XMM-Newton
observation. The presence of ionizing absorbers, albeit complicating the
spectrum analysis process, in principle does not affect the inclination
measurements. In the future, a systematic investigation of the impact of warm
absorbers on X-ray spectral analysis would be worthwhile. At the current
stage, case-by-case treatments of the soft band should be employed to help the
fit converge and reduce systematic errors by successful modeling of the
absorption.

Our mock tests indicate that systems with higher inclinations and larger
reflection fractions recover the input inclinations more reliably. Systems
with weaker reflection ($R_{f}\lesssim 0.3$) yield greater scatter in the
inclination recovery, whereas for reflection-dominant systems ($R_{f}\gtrsim
5$) the precision improves mildly with a larger reflection fraction. The
reflection fraction can be treated as a qualitative indicator of the
systematic bias of applying RelxillCp. Higher iron abundance and coronal
temperature help tighten the constraints on inner accretion disk inclination,
but the improvements are not as significant as the reflection fraction. The
spin of the BH, although degenerate with inclination, does not affect
significantly the measurement of inclination.

## Acknowledgements

LCH was supported by the National Science Foundation of China (11721303,
11991052, 12011540375, 12233001), the National Key R&D Program of China
(2022YFF0503401), and the China Manned Space Project (CMS-CSST-2021-A04, CMS-
CSST-2021-A06). We thank the anonymous referee for helpful suggestions. RD is
grateful to Cosimo Bambi, James F. Steiner, and Javier A. García for useful
suggestions in the early stages of the paper, Kishalay Choudhury for
explaining the details of Choudhury et al. (2017), Elias S. Kammoun for
explaining Kammoun et al. (2018) and sharing their script on resolution
binning, and Shuo Feng for advice on writing.

All figures in this paper are produced with SciencePlots (Garrett, 2021).
Other computational platforms and packages used include: AstroPy (Astropy
Collaboration et al., 2013, 2018, 2022), Corner (Foreman-Mackey, 2016), Numpy
(van der Walt et al., 2011; Harris et al., 2020), and Matplotlib (Hunter,
2007).

## References

  * Abbott et al. (2019) Abbott, B. P., Abbott, R., Abbott, T. D., et al. 2019, ApJ, 882, L24
  * Akaike (1974) Akaike, H. 1974, IEEE Transactions on Automatic Control, 19, 716
  * Antonucci (1993) Antonucci, R. 1993, ARA&A, 31, 473
  * Arnaud (1996) Arnaud, K. A. 1996, in ASP Conf. Ser., Vol. 101, Astronomical Data Analysis Software and Systems V, ed. G. H. Jacoby & J. Barnes (San Francisco, CA: ASP), 17
  * Arnaud (2016) Arnaud, K. A. 2016, HEAD, 115.02
  * Astropy Collaboration et al. (2022) Astropy Collaboration, Price-Whelan, A. M., Lim, P. L., et al. 2022, ApJ, 935, 167
  * Astropy Collaboration et al. (2018) Astropy Collaboration, Price-Whelan, A. M., Sipőcz, B. M., et al. 2018, AJ, 156, 123
  * Astropy Collaboration et al. (2013) Astropy Collaboration, Robitaille, T. P., Tollerud, E. J., et al. 2013, A&A, 558, A33
  * Ballantyne et al. (2014) Ballantyne, D. R., Bollenbacher, J. M., Brenneman, L. W., et al. 2014, ApJ, 794, 62
  * Baloković et al. (2020) Baloković, M., Harrison, F. A., Madejski, G., et al. 2020, ApJ, 905, 41.
  * Bardeen & Petterson (1975) Bardeen, J. M., & Petterson, J. A. 1975, ApJ, 195, L65
  * Baskin & Laor (2018) Baskin, A., & Laor, A. 2018, MNRAS, 474, 1970
  * Bautista & Kallman (2001) Bautista, M. A., & Kallman, T. R. 2001, ApJS, 134, 139
  * Belmont et al. (2008) Belmont, R., Malzac, J., & Marcowith, A. 2008, A&A, 491, 617
  * Beuchert et al. (2017) Beuchert, T., Markowitz, A. G., Dauser, T., et al. 2017, A&A, 603, A50
  * Birkinshaw & Davies (1985) Birkinshaw, M. & Davies, R. L. 1985, ApJ, 291, 32
  * Blandford & McKee (1982) Blandford, R. D., & McKee, C. F. 1982, ApJ, 255, 419
  * Blandford & Znajek (1977) Blandford, R. D., & Znajek, R. L. 1977, MNRAS, 179, 433
  * Bonson & Gallo (2016) Bonson, K., & Gallo, L. C. 2016, MNRAS, 458, 1927
  * Brenneman & Reynolds (2006) Brenneman, L. W., & Reynolds, C. S. 2006, ApJ, 652, 1028
  * Buchner et al. (2014) Buchner, J., Georgakakis, A., Nandra, K., et al. 2014, A&A, 564, A125
  * Caballero-García et al. (2020) Caballero-García, M. D., Papadakis, I. E., Dovčiak, M., et al. 2020, MNRAS, 498, 3184
  * Cash (1979) Cash, W. 1979, ApJ, 228, 939
  * Choudhury et al. (2017) Choudhury, K., García, J. A., Steiner, J. F., & Bambi, C. 2017, ApJ, 851, 57
  * Compton (1923) Compton, A. H. 1923, Physical Review, 21, 483
  * Costantini et al. (2007) Costantini, E., Gallo, L. C., Brandt, W. N., Fabian, A. C., & Boller, T. 2007, MNRAS, 378, 873
  * Crummy et al. (2006) Crummy, J., Fabian, A. C., Gallo, L., et al. 2006, MNRAS, 365, 1067
  * Czerny et al. (2016) Czerny, B., Du, P., Wang, J.-M., & Karas, V. 2016, ApJ, 832, 15
  * Czerny & Hryniewicz (2011) Czerny, B., & Hryniewicz, K. 2011, A&A, 525, L8
  * Czerny et al. (2011) Czerny, B., Hryniewicz, K., Nikołajuk, M., & Sadowski, A. 2011, MNRAS, 415, 2942
  * Czerny et al. (2017) Czerny, B., Li, Y.-R., Hryniewicz, K., et al. 2017, ApJ, 846, 154
  * Dauser et al. (2022) Dauser, T., García, J. A., Joyce, A., et al. 2022, MNRAS, 514, 3965
  * Dauser et al. (2014) Dauser, T., Garcia, J., Parker, M. L., Fabian, A. C., & Wilms, J. 2014, MNRAS, 444, L100
  * Dauser et al. (2016) Dauser, T., García, J., Walton, D. J., et al. 2016, A&A, 590, A76
  * Dauser et al. (2013) Dauser, T., Garcia, J., Wilms, J., et al. 2013, MNRAS, 430, 1694
  * Dauser et al. (2012) Dauser, T., Svoboda, J., Schartel, N., et al. 2012, MNRAS, 422, 1914
  * Dauser et al. (2010) Dauser, T., Wilms, J., Reynolds, C. S., & Brenneman, L. W. 2010, MNRAS, 409, 1534
  * den Herder et al. (2001) den Herder, J. W., Brinkman, A. C., Kahn, S. M., et al. 2001, A&A, 365, L7
  * Ding et al. (2022) Ding, Y., Li, R., Ho, L. C., & Ricci, C. 2022, ApJ, 931, 77
  * Done et al. (2012) Done, C., Davis, S. W., Jin, C., et al. 2012, MNRAS, 420, 1848
  * Done et al. (2013) Done, C., Jin, C., Middleton, M., & Ward, M. 2013, MNRAS, 434, 1955
  * Done et al. (2007) Done, C., Sobolewska, M. A., Gierliński, M., et al. 2007, MNRAS, 374, L15
  * Dovčiak et al. (2004) Dovčiak, M., Karas, V., & Yaqoob, T. 2004, ApJS, 153, 205
  * Dovčiak et al. (2022) Dovčiak, M., Papadakis, I. E., Kammoun, E. S., et al. 2022, A&A, 661, A135
  * Dove et al. (1997) Dove, J. B., Wilms, J., Maisack, M., & Begelman, M. C. 1997, ApJ, 487, 759
  * Einstein (1916) Einstein, A. 1916, Annalen der Physik, 354, 769
  * Emmering et al. (1992) Emmering, R. T., Blandford, R. D., & Shlosman, I. 1992, ApJ, 385, 460
  * Eracleous & Halpern (1994) Eracleous, M., & Halpern, J. P. 1994, ApJS, 90, 1
  * Fabian (2006) Fabian, A. C. 2006, Astronomische Nachrichten, 327, 943
  * Fabian et al. (2015) Fabian, A. C., Lohfink, A., Kara, E., et al. 2015, MNRAS, 451, 4375
  * Fabian et al. (1989) Fabian, A. C., Rees, M. J., Stella, L., & White, N. E. 1989, MNRAS, 238, 729
  * Fabian et al. (2002) Fabian, A. C., Vaughan, S., Nandra, K., et al. 2002, MNRAS, 335, L1
  * Fabian et al. (2009) Fabian, A. C., Zoghbi, A., Ross, R. R., et al. 2009, Nature, 459, 540
  * Fischer et al. (2013) Fischer, T. C., Crenshaw, D. M., Kraemer, S. B., & Schmitt, H. R. 2013, ApJS, 209, 1
  * Fischer et al. (2014) Fischer, T. C., Crenshaw, D. M., Kraemer, S. B., Schmitt, H. R., & Turner, T. J. 2014, ApJ, 785, 25
  * Foreman-Mackey (2016) Foreman-Mackey, D. 2016, The Journal of Open Source Software, 1, 24
  * Foreman-Mackey et al. (2013) Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306
  * Fürst et al. (2016) Fürst, F., Müller, C., Madsen, K. K., et al. 2016, ApJ, 819, 150
  * Gallo et al. (2004) Gallo, L. C., Boller, T., Brandt, W. N., et al. 2004, A&A, 417, 29
  * García et al. (2014) García, J., Dauser, T., Lohfink, A., et al. 2014, ApJ, 782, 76
  * García et al. (2013) García, J., Dauser, T., Reynolds, C. S., et al. 2013, ApJ, 768, 146
  * García & Kallman (2010) García, J., & Kallman, T. R. 2010, ApJ, 718, 695
  * García et al. (2011) García, J., Kallman, T. R., & Mushotzky, R. F. 2011, ApJ, 731, 131
  * García et al. (2019) García, J. A., Kara, E., Walton, D., et al. 2019, ApJ, 871, 88
  * Garrett (2021) Garrett, J. 2021, SciencePlots (v1.0.9), v.1.0.9, Zenodo, doi:10.5281/zenodo.5512926
  * Ghisellini et al. (1993) Ghisellini, G., Padovani, P., Celotti, A., & Maraschi, L. 1993, ApJ, 407, 65
  * Gliozzi & Williams (2020) Gliozzi, M. & Williams, J. K. 2020, MNRAS, 491, 532
  * Goodman & Weare (2010) Goodman, J., & Weare, J. 2010, Communications in Applied Mathematics and Computational Science, 5, 65
  * Greenhill et al. (2009) Greenhill, L. J., Kondratko, P. T., Moran, J. M., et al. 2009, ApJ, 707, 787
  * Guainazzi et al. (2006) Guainazzi, M., Bianchi, S., & Dovčiak, M. 2006, Astronomische Nachrichten, 327, 1032
  * Guilbert & Rees (1988) Guilbert, P. W., & Rees, M. J. 1988, MNRAS, 233, 475
  * Haardt & Maraschi (1991) Haardt, F., & Maraschi, L. 1991, ApJ, 380, L51
  * Haardt & Maraschi (1993) Haardt, F., & Maraschi, L. 1993, ApJ, 413, 507
  * Harris et al. (2020) Harris, C. R., Millman, K. J., van der Walt, S. J., et al. 2020, Nature, 585, 357
  * Harrison et al. (2013) Harrison, F. A., Craig, W. W., Christensen, F. E., et al. 2013, ApJ, 770, 103
  * Herrnstein et al. (1999) Herrnstein, J. R., Moran, J. M., Greenhill, L. J., et al. 1999, Nature, 400, 539
  * Hicks & Malkan (2008) Hicks, E. K. S., & Malkan, M. A. 2008, ApJS, 174, 31
  * Ho et al. (2000) Ho, L. C., Rudnick, G., Rix, H.-W., et al. 2000, ApJ, 541, 120
  * Hönig et al. (2007) Hönig, S. F., Beckert, T., Ohnaka, K., & Weigelt, G. 2007, in The Central Engine of Active Galactic Nuclei, ed. L. C. Ho & J.-M. Wang (San Francisco, CA: ASP), 487
  * Hopkins et al. (2012) Hopkins, P. F., Hernquist, L., Hayward, C. C., & Narayanan, D. 2012, MNRAS, 425, 1121
  * Humphrey et al. (2009) Humphrey, P. J., Liu, W., & Buote, D. A. 2009, ApJ, 693, 822
  * Hunter (2007) Hunter, J. D. 2007, Computing in Science and Engineering, 9, 90
  * Ingram et al. (2023) Ingram, A., Ewing, M., Marinucci, A., et al. 2023, MNRAS, 525, 5437
  * Ingram et al. (2019) Ingram, A., Mastroserio, G., Dauser, T., et al. 2019, MNRAS, 488, 324
  * Jansen et al. (2001) Jansen, F., Lumb, D., Altieri, B., et al. 2001, A&A, 365, L1
  * Jiang et al. (2019) Jiang, Y.-F., Blaes, O., Stone, J. M., & Davis, S. W. 2019, ApJ, 885, 144
  * Jin et al. (2017) Jin, C., Done, C., & Ward, M. 2017, MNRAS, 468, 3663
  * Jin et al. (2012) Jin, C., Ward, M., Done, C., & Gelbord, J. 2012, MNRAS, 420, 1825
  * Kaastra et al. (1996) Kaastra, J. S., Mewe, R., & Nieuwenhuijzen, H. 1996, UV and X-ray Spectroscopy of Astrophysical and Laboratory Plasmas, 411
  * Kaastra (1999) Kaastra, J. S. 1999, in X-ray Spectroscopy in Astrophysics, Vol. 520, ed. J. van Paradijs & J. A. M. Bleeker (Berlin and Heidelberg, Germany: Springer), 269
  * Kaastra & Bleeker (2016) Kaastra, J. S., & Bleeker, J. A. M. 2016, A&A, 587, A151
  * Kalberla et al. (2005) Kalberla, P. M. W., Burton, W. B., Hartmann, D., et al. 2005, A&A, 440, 775
  * Kallman & Bautista (2001) Kallman, T., & Bautista, M. 2001, ApJS, 133, 221
  * Kammoun et al. (2018) Kammoun, E. S., Nardini, E., & Risaliti, G. 2018, A&A, 614, A44
  * Kamraj et al. (2022) Kamraj, N., Brightman, M., Harrison, F. A., et al. 2022, ApJ, 927, 42.
  * Kara et al. (2016) Kara, E., Alston, W. N., Fabian, A. C., et al. 2016, MNRAS, 462, 511
  * Kerr (1963) Kerr, R. P. 1963, Phys. Rev. Lett., 11, 237
  * Kinch et al. (2020) Kinch, B. E., Noble, S. C., Schnittman, J. D., & Krolik, J. H. 2020, ApJ, 904, 117
  * Kinch et al. (2016) Kinch, B. E., Schnittman, J. D., Kallman, T. R., & Krolik, J. H. 2016, ApJ, 826, 52
  * King et al. (2005) King, A. R., Lubow, S. H., Ogilvie, G. I., et al. 2005, MNRAS, 363, 49
  * Kinney et al. (2000) Kinney, A. L., Schmitt, H. R., Clarke, C. J., et al. 2000, ApJ, 537, 152
  * Klein & Nishina (1929) Klein, O., & Nishina, T. 1929, Zeitschrift fur Physik, 52, 853
  * Komissarov (2022) Komissarov, S. S. 2022, MNRAS, 512, 2798
  * Kormendy & Ho (2013) Kormendy, J., & Ho, L. C. 2013, ARA&A, 51, 511
  * Krolik (1998) Krolik, J. H. 1998, Active Galactic Nuclei (Princeton: Princeton Univ. Press)
  * Laor (1991) Laor, A. 1991, ApJ, 376, 90
  * Li et al. (2018) Li, Y.-R., Songsheng, Y.-Y., Qiu, J., et al. 2018, ApJ, 869, 137
  * Li et al. (2015) Li, Y.-R., Wang, J.-M., Cheng, C., et al. 2015, ApJ, 804, 45
  * Lightman & White (1988) Lightman, A. P., & White, T. R. 1988, ApJ, 335, 57
  * Liska et al. (2021) Liska, M., Hesp, C., Tchekhovskoy, A., et al. 2021, MNRAS, 507, 983
  * Magdziarz & Zdziarski (1995) Magdziarz, P., & Zdziarski, A. A. 1995, MNRAS, 273, 837
  * Martocchia & Matt (1996) Martocchia, A., & Matt, G. 1996, MNRAS, 282, L53
  * Mastroserio et al. (2021) Mastroserio, G., Ingram, A., Wang, J., et al. 2021, MNRAS, 507, 55
  * Matt et al. (1991) Matt, G., Perola, G. C., & Piro, L. 1991, A&A, 247, 25
  * Metropolis et al. (1953) Metropolis, N., Rosenbluth, A. W., Rosenbluth, M. N., Teller, A. H., & Teller, E. 1953, J. Chem. Phys., 21, 1087
  * Middleton et al. (2016) Middleton, M. J., Parker, M. L., Reynolds, C. S., Fabian, A. C., & Lohfink, A. M. 2016, MNRAS, 457, 1568
  * Miller et al. (2015) Miller, J. M., Tomsick, J. A., Bachetti, M., et al. 2015, ApJ, 799, L6
  * Miniutti & Fabian (2004) Miniutti, G., & Fabian, A. C. 2004, MNRAS, 349, 1435
  * Misner et al. (1973) Misner, C. W., Thorne, K. S., & Wheeler, J. A. 1973, Gravitation (San Francisco, CA: W. H. Freeman and Company)
  * Murray et al. (1995) Murray, N., Chiang, J., Grossman, S. A., et al. 1995, ApJ, 451, 498
  * Nandra et al. (2013) Nandra, K., Barret, D., Barcons, X., et al. 2013, arXiv:1306.2307
  * Natarajan & Pringle (1998) Natarajan, P., & Pringle, J. E. 1998, ApJ, 506, L97
  * Newman et al. (1965) Newman, E. T., Couch, E., Chinnapared, K., et al. 1965, Journal of Mathematical Physics, 6, 918
  * Niedźwiecki et al. (2019) Niedźwiecki, A., Szanecki, M., & Zdziarski, A. A. 2019, MNRAS, 485, 2942
  * Novikov & Thorne (1973) Novikov, I. D., & Thorne, K. S. 1973, in Black Holes (Les Astres Occlus), ed. C. DeWitt & B. DeWitt (New York: Gordon and Breach), 343
  * Nyquist (1928) Nyquist, H. 1928, Physical Review, 32, 110
  * Pancoast et al. (2014) Pancoast, A., Brewer, B. J., Treu, T., et al. 2014, MNRAS, 445, 3073
  * Papaloizou & Lin (1995) Papaloizou, J. C. B., & Lin, D. N. C. 1995, ApJ, 438, 841
  * Petrucci et al. (2018) Petrucci, P. O., Ursini, F., De Rosa, A., et al. 2018, A&A, 611, A59
  * Pjanka et al. (2017) Pjanka, P., Greene, J. E., Seth, A. C., et al. 2017, ApJ, 844, 165
  * Planck Collaboration et al. (2020) Planck Collaboration, Akrami, Y., Arroja, F., et al. 2020, A&A, 641, A10
  * Ponti et al. (2018) Ponti, G., Bianchi, S., Muñoz-Darias, T., et al. 2018, MNRAS, 473, 2304
  * Pringle (1992) Pringle, J. E. 1992, MNRAS, 258, 811
  * Rees (1966) Rees, M. J. 1966, Nature, 211, 468
  * Reeves & Braito (2019) Reeves, J. N., & Braito, V. 2019, ApJ, 884, 80
  * Reynolds (2014) Reynolds, C. S. 2014, Space Sci. Rev., 183, 277
  * Reynolds (2021) Reynolds, C. S. 2021, ARA&A, 59, 117
  * Reynolds & Begelman (1997) Reynolds, C. S., & Begelman, M. C. 1997, ApJ, 488, 109
  * Reynolds et al. (1999) Reynolds, C. S., Young, A. J., Begelman, M. C., & Fabian, A. C. 1999, ApJ, 514, 164
  * Ricci et al. (2017) Ricci, C., Trakhtenbrot, B., Koss, M. J., et al. 2017, ApJS, 233, 17.
  * Risaliti et al. (2013) Risaliti, G., Harrison, F. A., Madsen, K. K., et al. 2013, Nature, 494, 449
  * Rogantini et al. (2022) Rogantini, D., Costantini, E., Gallo, L. C., et al. 2022, MNRAS, 516, 5171
  * Ross & Fabian (2005) Ross, R. R., & Fabian, A. C. 2005, MNRAS, 358, 211
  * Sambruna et al. (2011) Sambruna, R. M., Tombesi, F., Reeves, J. N., et al. 2011, ApJ, 734, 105
  * Scheuer & Feiler (1996) Scheuer, P. A. G. & Feiler, R. 1996, MNRAS, 282, 291
  * Shakura & Sunyaev (1973) Shakura, N. I., & Sunyaev, R. A. 1973, A&A, 24, 337
  * Shannon (1949) Shannon, C. E. 1949, IEEE Proceedings, 37, 10
  * Sharma (2017) Sharma, S. 2017, ARA&A, 55, 213
  * Silva et al. (2018) Silva, C. V., Costantini, E., Giustini, M., et al. 2018, MNRAS, 480, 2334
  * Sim et al. (2008) Sim, S. A., Long, K. S., Miller, L., et al. 2008, MNRAS, 388, 611
  * Skilling (2004) Skilling, J. 2004, in AIP Conf. Proc., Vol. 735, Bayesian Inference and Maximum Entropy Methods Science and Engineering: 24th International Workshop on Bayesian Inference and Maximum Entropy Methods in Science and Engineering, ed. R. Fischer, R. Preuss, & U. V. Toussaint (Melville, NY: AIP), 395
  * Sokal (1996) Sokal, A. D. 1996, Nuclear Physics B Proceedings Supplements, 47, 172
  * Storchi-Bergmann et al. (1995) Storchi-Bergmann, T., Eracleous, M., Livio, M., et al. 1995, ApJ, 443, 617
  * Storchi-Bergmann et al. (1997) Storchi-Bergmann, T., Eracleous, M., Teresa Ruiz, M., et al. 1997, ApJ, 489, 87
  * Strüder et al. (2001) Strüder, L., Briel, U., Dennerl, K., et al. 2001, A&A, 365, L18
  * Suebsuwong et al. (2006) Suebsuwong, T., Malzac, J., Jourdain, E., et al. 2006, A&A, 453, 773
  * Svoboda et al. (2012) Svoboda, J., Dovčiak, M., Goosmann, R. W., et al. 2012, A&A, 545, A106
  * Tagliacozzo et al. (2023) Tagliacozzo, D., Marinucci, A., Ursini, F., et al. 2023, MNRAS, 525, 4735
  * Tanaka et al. (1995) Tanaka, Y., Nandra, K., Fabian, A. C., et al. 1995, Nature, 375, 659
  * Tarter et al. (1969) Tarter, C. B., Tucker, W. H., & Salpeter, E. E. 1969, ApJ, 156, 943
  * Tashiro et al. (2018) Tashiro, M., Maejima, H., Toda, K., et al. 2018, in Proc. SPIE, Vol. 10699, Space Telescopes and Instrumentation 2018: Ultraviolet to Gamma Ray, ed. J.-W. A. den Herder, S. Nikzad, & K. Nakazawa, 1069922
  * Thorne & Price (1975) Thorne, K. S., & Price, R. H. 1975, ApJ, 195, L101
  * Torresi et al. (2010) Torresi, E., Grandi, P., Longinotti, A. L., et al. 2010, MNRAS, 401, L10
  * Tortosa et al. (2018) Tortosa, A., Bianchi, S., Marinucci, A., Matt, G., & Petrucci, P. O. 2018, A&A, 614, A37
  * Turner et al. (2001) Turner, M. J. L., Abbey, A., Arnaud, M., et al. 2001, A&A, 365, L27
  * Ursini et al. (2018) Ursini, F., Petrucci, P. O., Matt, G., et al. 2018, MNRAS, 478, 2663
  * van den Bosch et al. (2015) van den Bosch, R. C. E., Gebhardt, K., Gültekin, K., Yıldırım, A., & Walsh, J. L. 2015, ApJS, 218, 10
  * van der Walt et al. (2011) van der Walt, S., Colbert, S. C., & Varoquaux, G. 2011, Computing in Science and Engineering, 13, 22
  * Walton et al. (2013) Walton, D. J., Nardini, E., Fabian, A. C., et al. 2013, MNRAS, 428, 2901
  * Wilkins & Fabian (2011) Wilkins, D. R., & Fabian, A. C. 2011, MNRAS, 414, 1269
  * Wilkins & Fabian (2012) Wilkins, D. R., & Fabian, A. C. 2012, MNRAS, 424, 1284
  * Wilkins et al. (2021) Wilkins, D. R., Gallo, L. C., Costantini, E., Brandt, W. N., & Blandford, R. D. 2021, Nature, 595, 657
  * Wilkins et al. (2022) Wilkins, D. R., Gallo, L. C., Costantini, E., et al. 2022, MNRAS, 512, 761
  * Wilks (1938) Wilks, S. S. 1938, Ann. Math. Stat., 9, 60
  * Wilks (1963) Wilks, S. S. 1963, Mathematical Statistics (Princeton: Princeton University Press), Chapter 13
  * Williams et al. (2018) Williams, P. R., Pancoast, A., Treu, T., et al. 2018, ApJ, 866, 75
  * Wills & Brotherton (1995) Wills, B. J., & Brotherton, M. S. 1995, ApJ, 448, L81
  * Wilms et al. (2000) Wilms, J., Allen, A., & McCray, R. 2000, ApJ, 542, 914
  * Wu & Han (2001) Wu, X.-B., & Han, J. L. 2001, ApJ, 561, L59
  * Wu et al. (2022) Wu, Z., Ho, L. C., & Zhuang, M.-Y. 2022, ApJ, 941, 95
  * Zdziarski et al. (1996) Zdziarski, A. A., Johnson, W. N., & Magdziarz, P. 1996, MNRAS, 283, 193
  * Zhang & Wu (2002) Zhang, T.-Z., & Wu, X.-B. 2002, Chinese J. Astron. Astrophys., 2, 487
  * Zhuang et al. (2018) Zhuang, M.-Y., Ho, L. C., & Shangguan, J. 2018, ApJ, 862, 118
  * Życki et al. (1999) Życki, P. T., Done, C., & Smith, D. A. 1999, MNRAS, 309, 561