https://orcid.org/0000-0002-0273-218X
ABSTRACT
In this paper, we present an updated version of our model (KYNXiltr) which considers thermal reverberation of a standard Novikov–Thorne accretion disc illuminated by an X-ray point-like source. Previously, the model considered only two cases of black hole spins, and assumed a colour correction factor fcol = 2.4. Now, we extend the model to any spin value and fcol. In addition, we consider two scenarios of powering the X-ray corona, either via accretion, or external to the accretion disc. We use KYNXiltr to fit the observed time lags obtained from intense monitoring of four local Seyfert galaxies (NGC 5548, NGC 4593, Mrk 817, and Fairall 9). We consider various combinations of black hole spin, colour correction, corona height, and fraction of accretion power transferred to the corona. The model fits well the overall time-lag spectrum in these sources (for a large parameter space). For NGC 4593 only, we detect a significant excess of delays in the U band. The contribution of the diffuse BLR emission in the time-lag spectrum of this source is significant. It is possible to reduce the large best-fitting parameter space by combining the results with additional information, such as the observed Eddington ratio and average X-ray luminosity. We also provide an update to the analytic expression provided by Kammoun et al., for an X-ray source that is not powered by the accretion process, which can be used for any value of fcol, and for two values of the black hole spin (0 and 0.998).
1 INTRODUCTION
The current paradigm assumes that active galactic nuclei (AGNs) are powered by the accretion of matter onto a supermassive black hole (SMBH) from a geometrically thin and optically thick disc. Thermal UV photons emerge from the disc and are Compton up scattered in a medium of hot electrons, known as the corona, and emitted in X-rays (e.g. Lightman & White 1988; Haardt 1993). Assuming an isotropic emission, part of theses photons are directly emitted in the direction of the observer. The other part will illuminate back the accretion disc, where they will be partially reflected in X-rays (e.g. George & Fabian 1991; Matt, Perola & Piro 1991) and partially absorbed and then re-emitted in the UV/optical range in the form of thermalized emission. As the X-ray source varies in time, correlated variability should be detected in the UV/optical band with a time-lag depending on the wavelength (e.g. Cackett, Horne & Winkler 2007).
Intense, multiwavelength, monitoring campaign using space-based and ground-based telescopes confirmed the presence of correlation between AGN light curves in various UV/optical bands where the time delay increases as function of wavelength (e.g. Edelson et al. 2015; Fausnaugh et al. 2016; Cackett et al. 2018; McHardy et al. 2018; Edelson et al. 2019; Cackett et al. 2020; Hernández Santisteban et al. 2020; Pahari et al. 2020; Kara et al. 2021; Vincentelli et al. 2021, 2022; Donnan et al. 2023; Kara et al. 2023; McHardy et al. 2023). In many of these cases modelling the data assuming thermal reverberation was able to explain the shape of the time lags as function of wavelength (time-lag spectra), while underpredicting its amplitude (see e.g. Fausnaugh et al. 2016). However, the model thermal reverberation time lags were calculated using basic approximations like the viscous heating and the assumption that X-rays contribute equal amounts of energy to the disc at all radii. In addition, these models estimate the time lags by simply considering the light-travel delay radius emitting light at a wavelength λ, without considering the source height and general relativity effects due to the presence of the central black hole (BH), while the radius is computed by assuming Wien’s displacement law. This is a law which can predict the peak wavelength (λpeak) assuming a certain blackbody temperature, but it cannot give the correct temperature given λ (obviously the temperature would be different if instead of λ one would consider the frequency ν). In addition, the previous models would not consider BH spins other than zero. Despite these short comings, results from model fits to the observed time lags were accepted as an indication that the discs are larger than expected, requiring accretion rates that are larger than the ones inferred from multiwavelength analysis and broad-band spectral-energy distribution (SED) fitting.
In recent studies, Kammoun, Papadakis & Dovčiak (2019); Kammoun et al. (2021b, hereafter K21b) presented KYNXiltr, a model able to compute the response functions (Ψ) of Novikov–Thorne (NT) accretion discs (Novikov & Thorne 1973) illuminated by a lamp-post X-ray corona, taking into account all general relativity and disc ionization effects. These responses were used then to estimate the average time delay at a given wavelength (λ). K21b studied the effect of various AGN parameters on the time lags, and presented an analytic prescription that can be used to model the observed time-lag spectra. This was then used by Kammoun, Papadakis & Dovčiak (2021a, hereafter K21a) to successfully model the time-lag spectra in seven AGNs. The estimated accretion rate required to fit the time-lag spectra were in good agreement with literature. Panagiotou et al. (2020, 2022) used the same model to compute UV/optical power spectral densities (PSD) and were able to model the PSD from the long monitoring campaign of NGC 5548. Recently, Dovčiak et al. (2022) presented a new model (KYNSED) of the SED around accreting black holes considering the X-ray irradiation of an NT accretion disc by a lamp-post corona. KYNSED considers the case where the X-ray luminosity is independent of the accretion power of the disc, and the case where the X-ray luminosity is equal to the accretion power generated in the disc within a radius, rtransf. The X-ray luminosity is then parametrized by the ratio of the accretion power within rtransf to the total accretion power, Ltransf/Ldisc.
In this work, we present a new code KYNXiltr.1 which can be used to fit the observed time-lag spectra. There are a few differences between the new code and the approach of K21b. While the analytical time-lag equations presented by K21b assume that the X-ray luminosity is not part of the power that is liberated by the accretion process in the disc. The new code can compute time lags under this assumption but also under the assumption of the X-ray luminosity being equal to the accretion power generated within a certain radius, rtransf. Moreover, the model time lags in K21b were computed assuming a colour correction factor of fcol = 2.4 for the disc emission and that the disc extends from the innermost stable circular orbit (ISCO) to a (fixed) outer radius of |$R_{\rm out} = 10^4\, r_{\rm g}$|. In the new code, time lags can be computed for any colour correction factor and for any disc outer radius. Finally, K21b considered only two spins, namely a spin of 1 and zero, while the new code we present in this work can compute time lags for any spin value, from 0 up to 1.
In Section 2, we describe how the new code works, and we present a discussion of the time-lag dependence on the new parameters we introduce in our model. In Section 3, we introduce the AGN sample that we use to fit their observed time lags. In Section 4 we present the fitting technique we use and the results. Finally, we discuss the results and summarize our work in Section 5.
2 MODELLING THE X-RAY THERMAL REVERBERATION OF THE DISC
Similar to K21b, we assume a lamp-post X-ray corona at a given height (h) on the rotational axis of the black hole, which emits isotropically in its rest frame a power-law spectrum of the form |$f_X(t)=N(t)E^{-\Gamma } \rm {exp}\left(-E/E_{cut}\right)$|, where Γ is assumed to be constant. We need to determine the disc response in order to model the time lags when the disc is illuminated by variable X-rays. We provide below a short description of how this is done (a detailed description is given by Kammoun et al. 2021b).
Part of the X-ray flux received by the disc will be reflected and re-emitted in X-rays (this is the ‘disc reflection component’), and part of it will be absorbed. The absorbed X-rays will thermalize in the disc. They will act as an extra source of heating hence the local disc temperature, and the disc emission, will increase. As a result, the total disc flux, as observed by a distant observer, will vary with time because the X-ray flash first illuminates the inner disc and then propagates to the outer parts.
The disc response function, Ψ(λ, tobs), is equal to the flux that the disc emits due to X-ray heating at wavelength λ, and at time tobs, as measured by a distant observer. We identify all disc elements that brighten up at time tobs, for the observer, we compute the sum of their flux (say Ftot(λ, tobs)), we subtract the sum of their intrinsic disc flux, say |$F_{\rm NT, t_{obs}}(\lambda)$| (NT stands for the disc flux calculated following Novikov & Thorne 1973), and then the disc response is defined as
Once the disc response is known, then the disc flux that a distant observer will detect in the UV/optical bands, when the disc is constantly being illuminated by variable X-rays, will be given by
where Fobs(λ, t) is the flux at wavelength λ emitted by the whole disc at time t, and FNT(λt), is the flux emitted at λ by an NT disc.
The response function determines, to a large extent, the cross-correlation between the X-ray and the UV/optical variations. In fact, the centroid of the response function, should be representative of the maximum peak in the cross-correlation function (CCF) that is measured between the X-ray and the UV/optical light curves. The centroid is defined as follows:
Like K21b, KYNXiltr uses the equation above to compute model time lags. However, there are some differences between K21b and the new code, which we describe below.
K21b used the observed 2–10 keV luminosity of the corona (in Eddington units), LXobs, Edd(t), to normalize the response function in equation 1. In general, LX in this equation can be any quantity representative of the X-ray corona luminosity. In this work, we follow Dovčiak et al. (2022) and we use the total luminosity of the corona to normalize the disc response function. We parametrize LX as the ratio of the total X-ray luminosity over the accretion power, Ltransf/Ldisc. If Ltransf/Ldisc is positive, then LX is equal to the power that is released by the accretion process below a radius, which is transferred to the corona by an unknown physical mechanism. The total LX is this case is equal to (Ltransf/Ldisc)Ldisc, where Ldisc is set by the accretion rate, |$\dot{M}$|, as: |$L_{\rm disc}=\eta \dot{M}c^2$| (η being the accretion efficiency). A negative value of Ltransf/Ldisc would mean that the power given to the corona is external to the power released by the accretion power in the disc. In this case, Ltransf/Ldisc can be larger than unity, contrary to the case when Ltransf/Ldisc > 0, when it cannot be larger than 1 (in fact, the code does not allow the use of a value larger than 0.9 for Ltransf/Ldisc.)
The disc fluxes in equation 1 are modelled by a colour-corrected blackbody of the form
where Iν is the specific intensity, k is Boltzmann’s constant, h is Planck’s constant, and T is the disc temperature. In the case of |$F_{\rm NT,t_{obs}}(\lambda)$| the temperature of the disc elements that an observer detects at tobs is set according to NT. In the case of Ftot(λ, tobs) we use the new disc temperature, which is computed after adding to the local disc flux the X-ray flux absorbed by the disc. For both terms though, the emitted spectrum is not equal to that of a blackbody. It is well known that electron scattering plays an important role and can lead to deviations from blackbody emission. In general, transfer of energy between electrons and photons that are Compton scattered enforces a Wien tail at the high-energy end of the spectrum. The resulting spectra can be approximately modelled by the equation above, where fcol is the multiplicative factor by which spectral features are shifted to higher energies, while fcol−4 keeps the frequency integrated flux fixed. K21b assumed that fcol = 2.4, following Ross, Fabian & Mineshige (1992). The new code can provide time lags for any fcol value, including the prescription of Done et al. (2012), when fcol = −1.
An additional improvement of KYNXiltr when compared to the analytic expressions of K21b is related with the outer disc radius, Rout. Although K21b had studied in detail the effects of the disc outer radius to the time lags (see their Section 3.5), the analytic expressions presented by K21b were computed assuming an outer radius fixed at |$10^4\, r_{\rm g}$|. The new code can compute time lags from any outer disc radius, as this value may affect significantly the observed time lags (see for example McHardy et al. 2023). Finally, K21b presented analytic expressions for the time lags only in the case of a non-rotating and a maximally rotating black hole, while KYNXiltr can compute time lags for any BH spin value.
We investigated the effects of the BH mass, spin, accretion rate, corona height, disc inclination, and photon index on the response functions and the time lags, using KYNXiltr. The effects of these parameters are identical to the ones described in detail by K21b. In the following, we present only the effects of Ltransf/Ldisc and fcol on thermal reverberation time lags, as these are the only two parameters that were not studied by K21b.
2.1 Dependence of the time lags on Ltransf/Ldisc
Fig.1 shows the effect of changing Ltransf/Ldisc on the response functions at 2000 Å and 10000 Å, and the time-lag spectra for negative and positive values of Ltransf/Ldisc (upper and middle panels, respectively), and for different values of fcol (bottom panels). These simulations are performed for a* = 0.7, |$M_{\rm BH} = 5\times 10^7\, \rm M_\odot$|, |$h = 10\, \rm r_g$|, |$\dot{m}_{\rm Edd} = 0.05$|, Γ = 2, |$E_{\rm cut} = 300\, \rm keV$|, |$R_{\rm out} = 5000 \, \rm r_g$|. In the case of different values of Ltransf/Ldisc, we fixed fcol at 1. In the case of different values of fcol we fixed Ltransf/Ldisc at 0.5.
Response functions at 2000 Å and 10000 Å (left and middle panels, respectively) and time-lag spectra (right panels) for different negative values Ltransf/Ldisc (top), positive values of Ltransf/Ldisc (centre), and fcol (bottom). The fiducial parameters are assumed to be a* = 0.7, |$M_{\rm BH} = 5\times 10^7\, \rm M_\odot$|, |$h = 10\, \rm r_g$|, |$\dot{m}_{\rm Edd} = 0.05$|, Γ = 2, |$E_{\rm cut} = 300\, \rm keV$|, |$R_{\rm out} = 5000 \, \rm r_g$|. In the case of different values of Ltransf/Ldisc, we fixed fcol at 1. In the case of different values of fcol, we fixed Ltransf/Ldisc at 0.5.
In the case of negative Ltransf/Ldisc (which is the case in K21b) the effect is similar to changing the X-ray luminosity in K21b. Our results confirm the non-linearity of the disc response. This effect is clearer in fig. 15 of K21b, as for this case, we consider a wider range of X-ray luminosity to better highlight this effect. The observed |$2\!-\!10\, \rm keV$| luminosity considered in K21b ranges between 0.001 and 0.5 times the Eddington luminosity. However, in this work, given the considered values Ltransf/Ldisc, this translates into observed |$2\!-\!10\, \rm keV$| luminosity between 0.0004 and 0.003 times the Eddington luminosity. We recall that the disc response is normalized to the observed X-ray luminosity (see equation 3 in K21b). In this case, if the response of the disc scales linearly with the X-ray luminosity, we would expect the response functions for the various X-ray luminosity to overlap (in all bands), which is not the case. The response functions (at all wavelengths) decrease in amplitude and broaden in time as Ltransf/Ldisc increases. This is due to the fact that the thermalized flux does not contribute to the response in each waveband in the same way at all times. In addition, the increase in luminosity leads to an increase in the ionization state of the disc, especially in the inner regions, which will affect the reverberated flux. As a result, the increase in Ltransf/Ldisc(for the negative values) leads to an increase in the time lags at all wavelengths. These effects are detailed in Section 3.6 of K21b.
In the case of positive values of Ltransf/Ldisc, i.e. when the X-rays are powered by extracting energy from the inner parts of the accretion disc, early times are dominated by the reprocessed emission from illuminating the innermost regions of the disc, within which the power is extracted. Contrary to the case of the externally powered X-ray source, the emission from this part is broader and brighter for larger Ltransf/Ldisc. As Ltransf/Ldisc increases the radius rtransf within which the accretion power is transferred to the X-ray source increases which leads to a broader disc response from that region. In addition, the X-ray luminosity increases (assuming the same height) which leads to an increase in the amplitude of the reprocessed emission. As the time passes, we detect the emission from the outer parts of the disc that gets fainter with time. The peak due to the illumination of the innermost parts of the disc shifts the centroid of the response function towards smaller times, which leads to a decrease in the time-lag as Ltransf/Ldisc increases which, again, is contrary to what happens in the case of the externally illuminated X-ray source. This effect is amplified at longer wavelengths.
We note that the response function for the case of negative and positive Ltransf/Ldisc are identical at times when we observe the thermal reverberation from large radii (r > rtransf). At smaller times, the amplitude of Ψ is larger in the case when the X-ray source is powered by the accretion process. This difference is due to the fact that the response function is estimated as the difference between the total observed flux and the intrinsic NT disc flux. In the case of accretion powered X-ray source, the intrinsic disc emission is zero, which leads to larger value of Ψ compared to the case where the X-ray corona is not powered by accretion. This difference in Ψ at short time scales shifts the time lags to smaller values in the case where the corona is powered by the accretion process for a given value of Ltransf/Ldisc. The difference in time lags between the two cases increases as the X-ray luminosity increases, to reach |$\sim 60~{{\ \rm per\ cent}}$|. It is also worth noting that the range of time lags for different values of Ltransf/Ldisc is larger in the case of a corona powered by the accretion process.
2.2 Dependence of the time lags on fcol
The bottom panel of Fig. 1 shows the effect of changing fcol on the response functions and the time-lag spectra. We fixed all the parameters to the same values as before, considered Ltransf/Ldisc = 0.5, and we varied fcol between 1 and 2.5. The response functions decrease in amplitude and get broader for larger values of fcol. This effect is similar to an increase in |$\dot{m}/\dot{m}_{\rm Edd}$|. In fact, as fcol increases the observed temperature at each disc element increases. In addition, the total observed flux is re-normalized by |$f_{\rm col}^4$| in order to conserve the total emitted flux. For these reasons, as fcol increases the response function (considered as the difference between the total observed flux and the NT flux at a given wavelength) decreases in amplitude. Thus, the thermalized flux at a given wavelength will be emitted from outer regions of the disc and thus will last longer. This explains the fact that the response functions last longer (i.e. are broader) for larger fcol. As a result, the time lags increase as fcol increases.
3 SAMPLE
In order to demonstrate how well the new code works in practice, we chose to fit the time lags of NGC 4593, NGC 5548, Fairall 9, and Mrk 817. The time lags for these sources have been determined using long, densely sampled light curves at many wavelengths. Table 1 list the characteristics of these sources, together with references for the time-lag spectra we used. The data of NGC 4593 are taken from Cackett et al. (2018), which includes data from the HST and Swift monitoring. We note that these results are in agreement with a previous analysis of the Swift data by McHardy et al. (2018). The NGC 5548, Mrk 817, and Fairall 9 data are taken from Fausnaugh et al. (2016), Kara et al. (2021), and Hernández Santisteban et al. (2020), respectively. The Eddington ratios, λEdd = Lbol/LEdd, are listed in the fourth column of Table 1, and are based on measurements reported in the aforementioned papers (Lbol and LEdd are the bolometric and Eddington luminosities, respectively). We will accept those as indicators of the actual accretion rate of the source, normalized to the Eddington limit, although this assumption is not straight forward (see below). As for the X-ray properties (photon index and 2–10 keV luminosity) of the sources, we considered the values reported in K21a for NGC 5548 and NGC 4593. For Mrk 817, we considered the values reported by Kara et al. (2021) for the XMM-Newton observation which coincides with a state where the source is close to its mean flux during the monitoring campaign. As for Fairall 9, we followed the approach of K21a and we fitted the spectrum of the source by considering time intervals in which the source was close to its average flux during the monitoring campaign. The spectrum of the source was extracted using the automatic Swift/XRT generator2 (Evans et al. 2009). We fit the spectrum assuming a power-law model plus a reflection (xillver; García et al. 2013; García et al. 2016), considering only Galactic absorption. The model can be written in the XSPEC parlance as follows:
We fixed the Galactic column density to |$N_{\rm H} = 2.86\times 10^{20}\, \rm cm^{-2}$| (HI4PI Collaboration et al. 2016), the inclination of XILLVER to 20○, and the iron abundance to the solar value. The model resulted in a statistically good fit (|$\rm \chi ^2/dof = 218/197$|) with a best-fitting photon index Γ = 1.90 ± 0.03, an ionization parameter of the reflecting medium of |$\log \left(\xi /\rm erg\, cm\, s^{-1}\right) = 1.3_{-0.3}^{+0.2}$|, and an intrinsic X-ray luminosity of |$\rm \log L_{2\!-\!10} = 43.99 \pm 0.01$|. The values of the photon index and X-ray luminosity for all the four sources are reported in the last two columns of Table 1, and will be used in the rest of this analysis.
The sources in our sample. Luminosity distance, DL, BH mass, MBH, the Eddington ratio (λEdd), the photon index (Γ), and the 2–10 keV luminosity obtained from the X-ray spectral analysis (see text for details).
| Source | |$D_{\rm L}^{^a}$| | |$M_{\rm BH}^{^{b}}$| | |$\lambda _{\rm Edd}^{^{c}}$| | Γ | |$L_{\rm 2\!-\!10}^{^d}$| |
|---|---|---|---|---|---|
| (Mpc) | |$(10^7\, M_\odot$|) | (|$10^{43}\, \rm cgs$|) | |||
| NGC 4593 | 38.8 | |$0.76^{+ 0.16}_{-0.16}$| | 0.08 | 1.74 | 0.6 ± 0.2 |
| NGC 5548 | 80.1 | |$7.0 ^{+2.0}_{-2.8}$| | 0.05 | 1.70 | 2.5 ± 0.7 |
| Fairall 9 | 209 | |$19.9_{-4.6}^{+3.9}$| | 0.03 | 1.90 | 9.7 ± 1.9 |
| Mrk 817 | 138.7 | |$3.8^{+0.6}_{-0.6}$| | 0.20 | 1.90 | 1.9 ± 1.3 |
| Source | |$D_{\rm L}^{^a}$| | |$M_{\rm BH}^{^{b}}$| | |$\lambda _{\rm Edd}^{^{c}}$| | Γ | |$L_{\rm 2\!-\!10}^{^d}$| |
|---|---|---|---|---|---|
| (Mpc) | |$(10^7\, M_\odot$|) | (|$10^{43}\, \rm cgs$|) | |||
| NGC 4593 | 38.8 | |$0.76^{+ 0.16}_{-0.16}$| | 0.08 | 1.74 | 0.6 ± 0.2 |
| NGC 5548 | 80.1 | |$7.0 ^{+2.0}_{-2.8}$| | 0.05 | 1.70 | 2.5 ± 0.7 |
| Fairall 9 | 209 | |$19.9_{-4.6}^{+3.9}$| | 0.03 | 1.90 | 9.7 ± 1.9 |
| Mrk 817 | 138.7 | |$3.8^{+0.6}_{-0.6}$| | 0.20 | 1.90 | 1.9 ± 1.3 |
aWe computed DL using the source redshift, assuming a flat Universe with: |$H_0 = 67.8\, \rm km\, s^{-1}\, Mpc^{-1}$|, ΩΛ = 0.7, and ΩM = 0.3.
b MBH are from ‘The AGN Black Hole Mass Database’ (Bentz & Katz 2015) except for NGC 5548 we use the estimate by Horne et al. (2021); we use the estimates when considering all emission lines.
cReferences to the Eddington rations: Cackett et al. (2018) for NGC 4593, Fausnaugh et al. (2016) for NGC 5548, Vasudevan & Fabian (2007) for Fairall 9, and Kara et al. (2021) for Mrk 817.
dThe uncertainty on the 2–10 keV luminosity represent the scatter around the mean in the X-ray light curves of each of the sources.
The sources in our sample. Luminosity distance, DL, BH mass, MBH, the Eddington ratio (λEdd), the photon index (Γ), and the 2–10 keV luminosity obtained from the X-ray spectral analysis (see text for details).
| Source | |$D_{\rm L}^{^a}$| | |$M_{\rm BH}^{^{b}}$| | |$\lambda _{\rm Edd}^{^{c}}$| | Γ | |$L_{\rm 2\!-\!10}^{^d}$| |
|---|---|---|---|---|---|
| (Mpc) | |$(10^7\, M_\odot$|) | (|$10^{43}\, \rm cgs$|) | |||
| NGC 4593 | 38.8 | |$0.76^{+ 0.16}_{-0.16}$| | 0.08 | 1.74 | 0.6 ± 0.2 |
| NGC 5548 | 80.1 | |$7.0 ^{+2.0}_{-2.8}$| | 0.05 | 1.70 | 2.5 ± 0.7 |
| Fairall 9 | 209 | |$19.9_{-4.6}^{+3.9}$| | 0.03 | 1.90 | 9.7 ± 1.9 |
| Mrk 817 | 138.7 | |$3.8^{+0.6}_{-0.6}$| | 0.20 | 1.90 | 1.9 ± 1.3 |
| Source | |$D_{\rm L}^{^a}$| | |$M_{\rm BH}^{^{b}}$| | |$\lambda _{\rm Edd}^{^{c}}$| | Γ | |$L_{\rm 2\!-\!10}^{^d}$| |
|---|---|---|---|---|---|
| (Mpc) | |$(10^7\, M_\odot$|) | (|$10^{43}\, \rm cgs$|) | |||
| NGC 4593 | 38.8 | |$0.76^{+ 0.16}_{-0.16}$| | 0.08 | 1.74 | 0.6 ± 0.2 |
| NGC 5548 | 80.1 | |$7.0 ^{+2.0}_{-2.8}$| | 0.05 | 1.70 | 2.5 ± 0.7 |
| Fairall 9 | 209 | |$19.9_{-4.6}^{+3.9}$| | 0.03 | 1.90 | 9.7 ± 1.9 |
| Mrk 817 | 138.7 | |$3.8^{+0.6}_{-0.6}$| | 0.20 | 1.90 | 1.9 ± 1.3 |
aWe computed DL using the source redshift, assuming a flat Universe with: |$H_0 = 67.8\, \rm km\, s^{-1}\, Mpc^{-1}$|, ΩΛ = 0.7, and ΩM = 0.3.
b MBH are from ‘The AGN Black Hole Mass Database’ (Bentz & Katz 2015) except for NGC 5548 we use the estimate by Horne et al. (2021); we use the estimates when considering all emission lines.
cReferences to the Eddington rations: Cackett et al. (2018) for NGC 4593, Fausnaugh et al. (2016) for NGC 5548, Vasudevan & Fabian (2007) for Fairall 9, and Kara et al. (2021) for Mrk 817.
dThe uncertainty on the 2–10 keV luminosity represent the scatter around the mean in the X-ray light curves of each of the sources.
We note that for all sources we used the centroid values of the ICCF as reported in the corresponding papers. For NGC 4593 and Mrk 817, we changed the reference wavelength to be the smallest one (1150 Å and 1180 Å, respectively) in order to avoid any contribution from the Balmer jump.
4 TIME-LAG FITTING
X-ray to UV/optical time lags depend on many physical parameters of the system, such as the BH mass, spin, accretion rate, height of the corona, inclination, photon index, Ltransf/Ldisc, and Rout. Furthermore, one may assume an X-ray source which is powered by the accretion process, or by an external source of power. Hence, it is not straight forward to fit the observed time lags in practice. We believe it is not possible to fit the time-lag spectra by letting all parameters free to vary. In fact, some of the model parameters are degenerate (i.e. Ltransf/Ldisc, |$\dot{m}/\dot{m}_{\rm Edd}$|, and h can all affect the time lags in similar ways), and in some cases, the number of parameters may even be larger than the number of points in the observed time-lag spectra. We present below a possible approach in fitting the observed time lags in the four AGN in our sample. Throughout this analysis, we will assume that a part of the accretion power that is released below a certain radius is transferred to the corona (by an unknown mechanism), hence, we will assume a positive Ltransf/Ldisc.
4.1 Fitting method
In order to avoid the possible degeneracy between the various parameters and to speed up the fitting process, we proceeded as follows. For each of the sources, we fixed the BH mass and the photon index to the values given in Table 1. We also fixed the inclination to 45○. Then we considered three values of the BH spin (0., 0.7, and 0.998), three values of Ltransf/Ldisc|$(0.1, 0,5, \, \rm and\, 0.9)$|, and three values of fcol (1, 1.7, and 2.5). For each combination of these values we fixed the height at 10 values between 3 rg and 96 rg fitting only for |$\dot{m}/\dot{m}_{\rm Edd}$|. This results in 270 different combinations of the four model parameters (a*, h, Ltransf/Ldisc, fcol), with |$\dot{m}/\dot{m}_{\rm Edd}$| being the only free parameter of the fit. While the data quality does not allow for a direct fitting to be applied, leaving all parameters free, our approach can be used to exclude specific parts of the parameter range, and to also constrain the parameters of interest when including independent information, as will be shown below.
We fitted the time lags using a combination of two functions from the scipy.optimize library: curve_fit3 and basinhopping4 The first function solves a nonlinear least-squares problem with bounds on the variables. This method was chosen in order to reduce the fitting time comparing to other fitting scheme in the same library. However, using this method does not avoid the risk of falling into a local minimum. To overcome this problem, we use the basinhopping algorithm (Wales & Doye 1997), which first minimizes χ2, then randomly chooses a new starting point depending on the minimization result and launches a new minimization starting from this new point. After that, the results of both minimizations are compared based on their χ2 values. The best result is kept and the algorithm repeats the above starting from the point it has selected. The number of such cycles can be specified by the user as a parameter of this function. For our fitting tool, we modified the basinhopping method from SciPy to use curve_fit as the optimization method. We found that, since we are fitting for only one parameter, the fit converges to a point independent of the starting point after two basinhopping cycles, which is assumed to be the global minimum.
4.2 Results
We fit the time-lag spectra of each of the sources as described above. We did not consider time lags measurements between 2000–4000 Å in NGC 4593. According to (Cackett et al. 2018), the time lags spectrum of NGC 4593 shows a clear excess around the 3646Å Balmer jump, which could imply that diffuse emission from gas in broad-line region (BLR) may contribute significantly to the observed time lags. Edelson et al. (2015) also noticed that the U-band time-lag measurement of NGC 5548 was larger than the measurements in the surrounding bands, which could also be attributed to the Balmer continuum emission. Similar comments were made by Hernández Santisteban et al. (2020) in the case of the Fairall 9 time-lag observations. For that reason, we did not consider the measurements in the Swift/U and ground-based u bands in NGC 5548 and Fairall 9. We considered all time-lag measurements in the case of Mrk817.
For Fairall 9 and NGC 5548 we fixed the outer radius of the disc to |$R_{\rm out}=5000\, r_{\rm g}$|. However, for NGC 4593 and Mrk 817, this value of Rout underestimated the lag at longer wavelength. This is due to the fact that these two sources have lower BH masses and higher accretion rates compared to the former two sources. Thus, their discs are hotter, and a larger value of Rout is needed to fit the time-lag spectra. So we fixed it at 10000 rg for NGC 4593 and Mrk 817.
Figs A1–A4 of Appendix A show the time-lag spectra fitted using the various combinations of spin, height, fcol, and Ltransf/Ldisc considered in this work by letting |$\dot{m}/\dot{m}_{\rm Edd}$| as a free parameter. As it can be seen from these figures, a wide range of parameters can provide a good fit to the data (for different mdot values), for all of the sources. However, not all parameters are accepted. For example, the case a* = 0, Ltransf/Ldisc = 0.9, and fcol = 1 in Fig. A2 (bottom panel) clearly shows that the low height values do not provide a good fit to the time-lag spectra.
The range of the best-fitting |$\dot{m}/\dot{m}_{\rm Edd}$| and height depends strongly on the assumed combination of (a*, Ltransf/Ldisc, fcol). This is illustrated in Fig. 2, where we show the best-fitting |$\dot{m}/\dot{m}_{\rm Edd}$| as a function of the corona height for the various fcol values we consider, assuming a* = 0 and Ltransf/Ldisc = 0.5. The colour code in this figure shows the χ2 value obtained in each best-fitting (for all sources). The plots in this figure clearly show that various combinations of |$\dot{m}/\dot{m}_{\rm Edd}$|and X-ray source height can fit the observed time lags, in all sources, equally well.
Accretion rate versus coronal height for all of the sources assuming a* = 0 and Ltransf/Ldisc = 0.5 for fcol = 1, 1.7, and 2.5 (circles, squares, and triangles, respectively). The colour bars correspond to the value of χ2 obtained for each fit (the degrees of freedom for NGC 5548, NGC 4593, Fairall 9, and Mrk 817 are 12, 14, 6, and 12, respectively). The shaded areas correspond to the uncertainty by a factor of two around the Eddington ratio obtained from the literature for each of the sources.
Additional constraints can be put on the resulting best-fitting parameters by checking which values of the best-fitting |$\dot{m}/\dot{m}_{\rm Edd}$| are in agreement with the λEdd values listed in Table 1. As we mentioned above, we use λEdd as a measure of the accretion rate in each source, however this is highly uncertain. For example, the bolometric luminosity is usually based on measuring the observed luminosity in a given spectral band and then applying a bolometric correction factor. However, the observed luminosity may not be representative of the mean luminosity in this band, and the correction factor is quite uncertain. Furthermore, converting from Lbol/LEdd to |$\dot{m}/\dot{m}_{\rm Edd}$| requires a knowledge of the inclination of the system (which we do not know). Thus, we consider a conservative uncertainty on the observed Lbol/LEdd by a factor of ∼2, to account for all the aforementioned factors.
The horizontal grey shaded areas in Fig. 2 indicate the λEdd values listed in Table 1 together with the assumed uncertainty, as explained above. For the specific combination of a* = 0 and Ltransf/Ldisc = 0.5, Fig. 2 shows that fcol = 1 is not an accepted value, as all best-fitting values of |$\dot{m}/\dot{m}_{\rm Edd}$| are above the grey area, for all sources. On the other hand, there are quite a few heights for which the best-fitting |$\dot{m}/\dot{m}_{\rm Edd}$| lies within the grey area when fcol = 1.7 and 2.5. In some cases, the range of the heights which results in |$\dot{m}/\dot{m}_{\rm Edd}$| within the grey area is relatively narrow. See for example the bottom left panel in Fig. 2, for Fairall 9. The range of the accepted heights is quite narrow for fcol = 1.7, while fcol = 2.5 puts a strong low limit on the corona height at 10 rg. We conclude that, although many model parameters can result in model time lags which can fit the observations well, when we consider the (assumed) accretion rate for each source, we can still rule out some model parameter combinations.
An additional selection can be performed to reduce further the best-fitting parameter space. For a certain set of parameters, the code computes a posteriori the expected observed 2–10 keV luminosity. This quantity can be directly compared to the observed, average 2–10 keV luminosity, as listed in Table 1. Fig. 3 shows the best-fitting |$\dot{m}/\dot{m}_{\rm Edd}$| and height values (upper and bottom panels, for each source), as a function of the model L2–10 in units of the Eddington luminosity LEdd (hereafter L2-10, fit/LEdd). In this figure, we show only the mdot and heights that lead to a statistically acceptable fit of the time-lag spectra (i.e. pnull > 0.01). Each of the columns in this figure corresponds to a fixed spin value. We show in different colours the results for various Ltransf/Ldisc, and in different symbols the results for various fcol. The vertical lines indicate the observed 2–10 keV luminosity divided by LEdd (L2-10, obs/LEdd). The vertical grey shaded area indicates the 1σ scatter around the mean as estimated from the Swift /XRT light curves. The horizontal solid lines indicate the λEdd values listed in Table 1, while the green shaded area shows the uncertainty associated with it (as explained above).
Best-fitting |$\dot{m}/\dot{m}_{\rm Edd}$| and h versus L2-10, fit/LEdd (top and bottom rows, respectively) for NGC 4593, NGC 5548, and Fairall 9 (top to bottom). The results are shown for a* = 0, 0.7, and 0.998 (left to right) for Ltransf/Ldisc = 0.1, 0.5, and 0.9 (red, blue, and black, respectively) and fcol = 1, 1.7, and 2.5 (circles, triangles, and squares, respectively). The vertical lines and the shaded grey areas represent the average value of L2-10, obs/LEdd and the corresponding 1σ uncertainty. The horizontal lines and the shaded green areas represent the values of the Eddington ratio obtained from the literature and the corresponding factor of 2 uncertainty, respectively best-fitting results for Mrk 817.
This figure shows clearly that, when we consider the additional constraints of the observed X-ray luminosity, then the accepted parameter space of the best-fitting |$\dot{m}/\dot{m}_{\rm Edd}$| and height values is further reduced. For example, the left-hand side, top panels in Fig. 3 show the best-fitting |$\dot{m}/\dot{m}_{\rm Edd}$| and height values in the case of NGC 4593, when we assume a BH with zero spin. The uppermost panel indicates that only the model parameters indicated by the black triangles can fit the time lags and, at the same time, can also result in accretion rates and 2–10 keV band luminosity which are in agreement with the observations. Triangles imply an fcol of 1.7, while the black colour indicates an Ltransf/Ldisc of 0.9. The panel just below the uppermost left panel, show that only heights of |$\sim 20-70\, r_{\rm g}$| are consistent with both the |$\dot{m}/\dot{m}_{\rm Edd}$| and X-ray luminosity constraints. The red circles in the same panel indicate that the combination of fcol = 1, Ltransf/Ldisc = 0.1 and a source height above |$\sim 5\, r_{\rm g}$| could also fit the observed time lags well and predict the observed X-ray luminosity, but the necessary |$\dot{m}/\dot{m}_{\rm Edd}$| value would be ∼8 times larger than the λEdd value listed in Table 1, so we do not consider this as a probable combination of model parameters for this source.
Following the same procedure, we selected the best-fitting model parameters that can provide a good fit to the observed time lags and, at the same time, can also predict an accretion rate and 2–10 keV luminosity which are consistent with λEdd and the observed luminosity values listed in Table 1, within their uncertainties (defined as explained above). Our final, best-fitting model parameters for the four sources are listed in Table 2. The last column in this Table lists the best-fitting χ2 values together with the number of degrees of freedom (dof). All models listed in this table provide an acceptable fit to the data. For example, even in the case of the Ltransf/Ldisc = 0.9, a* = 0, fcol = 2.5 model fit to the NGC 5548 data, which resulted in the largest |$\rm \chi ^2/dof$| ratio, the null hypothesis probability is 0.025, i.e. larger than 0.01).
Range of the best-fitting parameters obtained by fitting the time-lag spectra of each of the four sources.
| Ltransf/Ldisc | a* | fcol | |$h\, (r_{\rm g})$| | |$\rm \chi ^2_{min}/ dof$| |
|---|---|---|---|---|
| NGC 4593 | ||||
| 0.5 | 0.7 | 1.7 | [5.6; 42.7] | 18.6/14 |
| 0.998 | 1.7 | ≥64 | 19.9/14 | |
| 0.9 | 0 | 1.7 | [28.5; 64] | |$\rm 16.8/14$| |
| 0.7 | 1.7 | ≥42.7 | |$\rm 18.4/14$| | |
| 0.998 | 1.7 | ≥64 | 20.3/14 | |
| 0.998 | 2.5 | [5.6; 8.4] | 18.5/14 | |
| NGC 5548 | ||||
| 0.5 | 0 | 1.7 | ≥8.4 | 20.4/12 |
| 0.7 | 1.7 | ≥64 | 19.7/12 | |
| 0.998 | 2.5 | [5.6; 28.5] | 20.3/12 | |
| 0.9 | 0 | 2.5 | [18.9; 28.5] | 23.4/12 |
| 0.998 | 2.5 | [18.9; 42.7] | 20.5/12 | |
| Fairall 9 | ||||
| 0.9 | 0 | 1.7 | ≥64 | 5.4/6 |
| 0 | 2.5 | 8.4 | 5.5/6 | |
| 0.7 | 1.7 | ≥64 | 5.1/6 | |
| 0.998 | 1.7 | 96 | 4.7/6 | |
| Mrk 817 | ||||
| 0.1 | 0 | 2.5 | ≥4 | 11.2/12 |
| 0.7 | 2.5 | ≥64 | 11.2/12 | |
| Ltransf/Ldisc | a* | fcol | |$h\, (r_{\rm g})$| | |$\rm \chi ^2_{min}/ dof$| |
|---|---|---|---|---|
| NGC 4593 | ||||
| 0.5 | 0.7 | 1.7 | [5.6; 42.7] | 18.6/14 |
| 0.998 | 1.7 | ≥64 | 19.9/14 | |
| 0.9 | 0 | 1.7 | [28.5; 64] | |$\rm 16.8/14$| |
| 0.7 | 1.7 | ≥42.7 | |$\rm 18.4/14$| | |
| 0.998 | 1.7 | ≥64 | 20.3/14 | |
| 0.998 | 2.5 | [5.6; 8.4] | 18.5/14 | |
| NGC 5548 | ||||
| 0.5 | 0 | 1.7 | ≥8.4 | 20.4/12 |
| 0.7 | 1.7 | ≥64 | 19.7/12 | |
| 0.998 | 2.5 | [5.6; 28.5] | 20.3/12 | |
| 0.9 | 0 | 2.5 | [18.9; 28.5] | 23.4/12 |
| 0.998 | 2.5 | [18.9; 42.7] | 20.5/12 | |
| Fairall 9 | ||||
| 0.9 | 0 | 1.7 | ≥64 | 5.4/6 |
| 0 | 2.5 | 8.4 | 5.5/6 | |
| 0.7 | 1.7 | ≥64 | 5.1/6 | |
| 0.998 | 1.7 | 96 | 4.7/6 | |
| Mrk 817 | ||||
| 0.1 | 0 | 2.5 | ≥4 | 11.2/12 |
| 0.7 | 2.5 | ≥64 | 11.2/12 | |
Range of the best-fitting parameters obtained by fitting the time-lag spectra of each of the four sources.
| Ltransf/Ldisc | a* | fcol | |$h\, (r_{\rm g})$| | |$\rm \chi ^2_{min}/ dof$| |
|---|---|---|---|---|
| NGC 4593 | ||||
| 0.5 | 0.7 | 1.7 | [5.6; 42.7] | 18.6/14 |
| 0.998 | 1.7 | ≥64 | 19.9/14 | |
| 0.9 | 0 | 1.7 | [28.5; 64] | |$\rm 16.8/14$| |
| 0.7 | 1.7 | ≥42.7 | |$\rm 18.4/14$| | |
| 0.998 | 1.7 | ≥64 | 20.3/14 | |
| 0.998 | 2.5 | [5.6; 8.4] | 18.5/14 | |
| NGC 5548 | ||||
| 0.5 | 0 | 1.7 | ≥8.4 | 20.4/12 |
| 0.7 | 1.7 | ≥64 | 19.7/12 | |
| 0.998 | 2.5 | [5.6; 28.5] | 20.3/12 | |
| 0.9 | 0 | 2.5 | [18.9; 28.5] | 23.4/12 |
| 0.998 | 2.5 | [18.9; 42.7] | 20.5/12 | |
| Fairall 9 | ||||
| 0.9 | 0 | 1.7 | ≥64 | 5.4/6 |
| 0 | 2.5 | 8.4 | 5.5/6 | |
| 0.7 | 1.7 | ≥64 | 5.1/6 | |
| 0.998 | 1.7 | 96 | 4.7/6 | |
| Mrk 817 | ||||
| 0.1 | 0 | 2.5 | ≥4 | 11.2/12 |
| 0.7 | 2.5 | ≥64 | 11.2/12 | |
| Ltransf/Ldisc | a* | fcol | |$h\, (r_{\rm g})$| | |$\rm \chi ^2_{min}/ dof$| |
|---|---|---|---|---|
| NGC 4593 | ||||
| 0.5 | 0.7 | 1.7 | [5.6; 42.7] | 18.6/14 |
| 0.998 | 1.7 | ≥64 | 19.9/14 | |
| 0.9 | 0 | 1.7 | [28.5; 64] | |$\rm 16.8/14$| |
| 0.7 | 1.7 | ≥42.7 | |$\rm 18.4/14$| | |
| 0.998 | 1.7 | ≥64 | 20.3/14 | |
| 0.998 | 2.5 | [5.6; 8.4] | 18.5/14 | |
| NGC 5548 | ||||
| 0.5 | 0 | 1.7 | ≥8.4 | 20.4/12 |
| 0.7 | 1.7 | ≥64 | 19.7/12 | |
| 0.998 | 2.5 | [5.6; 28.5] | 20.3/12 | |
| 0.9 | 0 | 2.5 | [18.9; 28.5] | 23.4/12 |
| 0.998 | 2.5 | [18.9; 42.7] | 20.5/12 | |
| Fairall 9 | ||||
| 0.9 | 0 | 1.7 | ≥64 | 5.4/6 |
| 0 | 2.5 | 8.4 | 5.5/6 | |
| 0.7 | 1.7 | ≥64 | 5.1/6 | |
| 0.998 | 1.7 | 96 | 4.7/6 | |
| Mrk 817 | ||||
| 0.1 | 0 | 2.5 | ≥4 | 11.2/12 |
| 0.7 | 2.5 | ≥64 | 11.2/12 | |
If we accept the λEdd values listed in Table 1 as reliable estimates of the accretion rate in these objects, our results indicate that fcol cannot be equal to one, i.e. the emission from the accretion disc in these objects is not simply equal to the flux emitted by multitemperature blackbodies, where the temperature depends on radius according to the Novikov & Thorne (1973) and Shakura & Sunyaev (1973) prescriptions. In addition, our results show that, if we assume that the X-ray source is powered by the accretion process, then our results indicate that a fraction of the accretion power released in the accretion disc that is larger than 50 per cent should be transferred to the X-ray corona. The exception is Mrk 817, which is the object with the highest accretion rate among the four sources. Our results indicate that the X-ray corona luminosity should be smaller than 50 per cent of the accretion power in this source. This source is also the only one for which we cannot find a good fit to the data for a* = 0.998. All spins we consider can provide good fits to the data in the other three sources. Furthermore, the accepted range of the X-ray source height is rather broad. In many cases we find solutions where the source height is large (this is mainly the case with Fairall 9), but we also get solutions where the X-ray source height is as small as |$\sim 5-6\, r_{\rm g}$|.
NGC 5548 and NGC 4593 are the only sources in common with K21a. For NGC 5548, K21a reported |$h \in [23\, r_{\rm g}, 58\, r_{\rm g}]\, ([33\, r_{\rm g}, 78\, r_{\rm g}])$| and |$\dot{m}/\dot{m}_{\rm Edd}\le 0.008\, (0.05)$| for |$a^\ast = 0\, (a^\ast = 1)$|, at the 1σ confidence level. Considering fcol = 2.5 (comparable to the assumptions of K21b), we find a good agreement with the results of K21a for high spin values. We find a larger value of |$\dot{m}/\dot{m}_{\rm Edd}$|, while the ranges of height are in agreement, for the non-spinning BH case. The larger |$\dot{m}/\dot{m}_{\rm Edd}$| value is due to the fact that we consider the case of an X-ray source powered via accretion. This results into smaller time lags compared to an X-ray source powered externally, which was the case in K21a. To compensate for this, the model fits the data but with a larger |$\dot{m}/\dot{m}_{\rm Edd}$|. As for NGC 4593, K21a report |$h \le 25\, r_{\rm g}\, (23\, r_{\rm g})$| and |$\dot{m}/\dot{m}_{\rm Edd}\in [0.006, 0.016]\, ([0.06, 0.22])$| for |$a^\ast = 0\, (a^\ast = 1)$|, at the 1σ confidence level. For the non-spinning case, we do not find any solution that fits the time-lag spectra. This is comparable to the results of K21a, as the values of |$\dot{m}/\dot{m}_{\rm Edd}$| are smaller than the observed Eddington ratio. The results for a maximally spinning black hole are in agreement with K21a. We note that for both sources, the 1σ parameter range from the current fit is significantly smaller when compared to K21a.
5 DISCUSSION AND CONCLUSIONS
We present a new code which can be used to fit the time-lag spectra that has been observed the last few years in many AGN, under the assumption of disc thermal reverberation, due to X-ray illumination. Our work extends the work presented by K21b. These authors presented analytic functions for the time lags in the case when the source that powers the X-ray corona is not associated with the accretion power and the BH spin is either zero or one. They also assumed a colour correction factor of fcol = 2.4, and an infinite outer disc radius. The new code is based on the work of K21b, as it assumes a point-like X-ray source illuminating an NT accretion disc in the lamp-post geometry. However, we also take into account the recent work of D22. Consequently, the new code offers significant improvements when compared with the analytic functions of K21b in many ways: a) it can be used both in the case when the X-ray source is powered by a source which is not associated with the accretion process (as in K21b) but also in the case when the X-ray luminosity is a fraction of the accretion process (which is, somehow, transferred to the X-ray corona) b) it can be used to fit the data for any BH spin, from zero to 0.998, c) any fcol value, and d) any outer disc radius. We note that, the analytic functions of K21b, as well as the new code, can be used to fit time lags by keeping the BH mass as a free parameter. As we have already argued, one should try to minimize the number of free parameters when fitting the time lags, however, recent work by Pozo Nuñez et al. (2019) suggests that black hole mass measurements in AGN may be underestimated due to the unknown BLR geometry. In this case, one can leave the BH mass as a free parameter when fitting the time lags in order to investigate this possibility.
We studied in detail the dependence of the time lags on Ltransf/Ldisc in the case when the X-ray source is powered by an external source and in the case when it is powered by the accretion process. We find that, for the same Ltransf/Ldisc, the time lags are smaller in the latter case. This is because, in this case, the (relative) contribution of the inner disc emission to the observed flux in each optical/UV band increases with respect to the case of an externally powered X-ray corona. As for fcol, we find strong effects on the time lags that are similar to the effects of the accretion rate. The time lags increase with larger fcol, just like the time lags increase with increasing accretion rate. Therefore, fcol and |$\dot{m}/\dot{m}_{\rm Edd}$| (at least) should be degenerate, when it comes to the determination of either parameter from fits to the observed time lags. It is also worth noting that in this work we use the values of Eddington ratio reported in the literature, which do not take into consideration the effect of intrinsic reddening that may be important (see e.g. Gaskell et al. 2023). If intrinsic reddening is significant, this will affect the value of the observed Eddington ratio. This should not alter the shape of the time lags, however it may affect the values of the best-fitting parameters.
We used the code to fit the observed time lags in four AGN. We chose these objects mainly because their time lags have been determined in many wavebands, from the far-UV up to ∼8000–9000 Å, thus they are the best time-lag spectra to test theoretical models. To fit the data, we fix the BH mass and outer disc radius, and then we consider a large number of fcol, Ltransf/Ldisc, and spin values, and we fit the data by letting just |$\dot{m}/\dot{m}_{\rm Edd}$| and the X-ray source height to be free parameters. This approach is necessary in order to fit the data because, as we have already mentioned, many model parameters affect the time lags in similar ways, hence introducing various degrees of degeneracy between the parameters. We find that the time lags of all sources are well fitted by the model for a large range of model parameters (see for example the best-fitting models plotted in Figs A1–A4 in appendix). By introducing further observational constraints such as the observed λEdd and the observed |$2\!-\!10\, \rm keV$| luminosity, we are able to reduce the range of parameters that can explain the time-lag spectra and be consistent with the other constrains as well is reduced.
For all sources, we find values of fcol that are greater than 1. We emphasize that this result depends on whether the observed λEddvalues listed in Table 1 are indeed representative of the accretion rate in these objects or not. If the intrinsic accretion rate is larger in these objects, then fcol could be closer to unity. Our results are in agreement with previous works such as Ross et al. (1992); Shimura & Takahara (1995); Done et al. (2012); Davis & El-Abd (2019), who indicated the need for modifications to the standard blackbody accretion disc model, in the form of a colour temperature corrected blackbody. In particular, Davis & El-Abd (2019) found a moderate variation of fcol between 1.4–2, for accretion rates between 0.01 and 1 of the Eddington rate. Their equation (10) presents an analytic approximation of fcol as a function of the BH mass and accretion rate for a* = 0 and 0.9. We used this equation and we found that fcol factors of ∼1.6(1.7), 1.55(1.65), 1.4(1.5) and 1.73(1.83) for NGC 4593, NGC 5548, Fairall 9 and Mrk 817, respectively, in the case of a* = 0(0.9), when assuming that α = 0.1. Our results, based on time-lag modelling, predict slightly larger values of fcol, but we did not investigate a dense range of fcol values in our work. We suspect that models with fcol values equal to the ones reported above will almost certainly give good fits to the observed time lags. However, Davis & El-Abd (2019) considered a non X-ray illuminated accretion disc. In fact, for X-ray illuminated discs, fcol may be slightly different than the ones reported by Davis & El-Abd (2019). Given this uncertainty, we believe our measurements are in very good agreement with the predictions of Davis & El-Abd (2019), and this result indicates that colour correction factors may be necessary when fitting broadband AGN SEDs, as well as timing results.
We suggest to always use KYNXiltr when fitting the observed time-lag spectra, as it can cover both scenarios of the externally and the internally heated X-ray corona. In fact, in this way, we may get indirect evidence regarding the source of power that heats the X-ray corona. For example, if the internally powered X-ray corona models do not fit well the time-lag spectra, even for the maximum allowed value of Ltransf/Ldisc = 0.9, this could be an indication that the X-ray corona is powered by other mechanisms (most probably associated with transfer of power from the vicinity of the BH), as long as the model can fit the data in the case of the externally heated corona with a larger Ltransf/Ldisc. The analytic expressions of K21b are valid in the case of an externally heated corona, as long as a* = 0 or 1, and the outer disc radius is larger than 10000 rg, and fcol = 2.4. The restriction of the two spin values may not be very strong, as it seems rather unlikely that it would be possible to derive a spin parameter estimate with a small error with the currently observed time-lag spectra, given the large number of physical parameters and their inter-dependencies. We provide updated versions of the analytic expressions in Appendix B, which take into account the model dependence on fcol. We also discuss ways one could get a rough estimate of time lags in the case of internally powered X-ray coronæ, using the analytic time-lag equations of K21b.
When fitting the time lags, we omitted data points in the U-band (NGC 5548 and Fairall 9) while in the case of NGC 4593, between 2000–4000 Å. This is because the time-lag measurements in these wavelengths may be affected by line and continuum emission from gas in the BLR. This issue was already noticed by many authors in the past. Recently, Netzer (2022) suggested that the total lag-spectrum, and its normalization, could be due to diffuse emission from radiation pressure supported clouds in a BLR with a covering factor of about 0.2. Their results were based on the assumption of X-ray disc thermal reverberation time lags which are significantly different than the ones we present here. For example, the continuum time lags shown by the dashed line in Fig. 1 of Netzer (2022) are very different than our best-fitting models presented in Figs A1–A4.
Our model fits show that the majority of the observed time lags can be fully explained when self-consistently modelling the X-ray irradiation of an NT accretion disc, for reasonable physical parameters. We fit the Mrk 817 time-lag spectrum very well without omitting any points in the U-band. All the observed time lags in this object can be fully accounted by X-ray, thermal reverberation. In the case of Fairall 9, the best-list models listed in Table 2 can provide an acceptable fit to the full time-lag spectrum even if we consider the U-band measurements (χ2 increases from 5–5.5 to 15–15.5/8 dof, for all the models listed in Table 2, which implies pnull ∼ 0.05–0.06). The U-band bump is more pronounced in the case of NGC 5548, where pnull becomes less than 0.01 when we consider the U-band time lags as well, for some of the best-fitting models listed in Table 2. However, we do not detect the 5000–10, 000 Å feature (close to the wavelength of the Paschen jump at 8200 Å), as suggested by Netzer (2022), in the best-fitting residuals of NGC 5548 and Fairall 9.
The X-ray reverberation time lags cannot account for the U-band bump in the NGC 4593 time lags. A bump at around 6000–8000 Å is also observed in this object, although it is not as pronounced as the one in the U-band. We therefore conclude that, at least in some AGN, the U-band time lags cannot be fully explained if we assume only X-rays illumination of the disc. It is unclear why the time lags due to diffuse emission in the BLR may be present in one object and not in others. It may depend on the geometry of the BLR, the way the central source illuminates the clouds etc. Investigating this issue is not straightforward. In any case, we plan to investigate in the future whether we can update our code by including the contribution of the time lags from the diffuse BLR emission.
Our work shows that considering only time lags is not enough to directly constrain all the different parameters of the corona/disc configuration. In order to entirely lift any degeneracy between the various model parameters, and find a unique set of parameters that can explain the observed data, it is also required to fit the observed broad-band X-ray to UV/Optical SED (Dovčiak et al. 2022) and the power spectra (Panagiotou et al. 2020, 2022). This will be invoked in follow-up publications.
ACKNOWLEDGEMENTS
We thank the anonymous referee for their helpful comments. ESK and LR acknowledge financial support from the Centre National d’Etudes Spatiales (CNES) from which part of this work was completed. IEP acknowledges support from the European Union’s Horizon 2020 Framework Programme under the AHEAD2020 project (grant agreement n. 871158). This work made use of data supplied by the UK Swift Science Data Centre at the University of Leicester. This work makes use of Matplotlib (Hunter 2007), NumPy (Harris et al. 2020), and SciPy (Virtanen et al. 2020).
DATA AVAILABILITY
The data for the time-lag spectra are all published in previous works cited in Table 1. The code is publicly available at https://projects.asu.cas.cz/dovciak/kynxiltr.
Footnotes
The code is publicly available at https://projects.asu.cas.cz/dovciak/kynxiltr. We present also a detailed description of how the code can be used to perform simulations of time lags or fit observed time-lag spectra.
References
APPENDIX A: TIME-LAG SPECTRA FITS
Fitting the time-lag spectra of NGC 4593 for different values of Ltransf/Ldisc. We assume a* = 0, 0.7, and 0.998 (top, middle, and bottom rows, respectively) and fcol = 1, 1.7, and 2.5 (left to right columns). The colour maps correspond to the different values of the coronal height. Open circles show time-lag measurements that we did not consider when fitting the data, as they may be affected by emission from the much larger (than the disc) BLR.
Same as Fig. A1 but for NGC 5548.
Same as Fig. A1 but for Fairall 9.
Same as Fig. A1 but for Mrk 817.
APPENDIX B: UPDATING THE ANALYTIC PRESCRIPTION
K21b derived an analytic expression of the time-lag spectra as a function of the BH mass, corona height, observed |$2\!-\!10\, \rm keV$| luminosity, and the accretion rate, for spins of 0 and 0.998. The resulting functions were computed by assuming that the X-ray source is powered via an external process that is not linked to accretion (equivalent to a negative value of Ltransf/Ldisc in the current code) and for fcol = 2.4. In this section, we provide a new term that can be included in the K21b equations to take into account the use of fcol values different than 2.4. In addition, we compare the time lags when we assume an in internally and an externally powered X-ray corona, and we comment on how the analytical functions could be scaled to provide time-lag estimates in the former case.
B1 The dependence of time lags on fcol
In this section, we investigate the dependency of the time lags on fcol, in the case of externally powered corona. We performed a set of simulations assuming the same fiducial parameters as K21b, but for a range of fcol between 1 and 2.5, for a* = 0 and 0.998. We assumed, in all of them, Ltransf/Ldisc = −0.5. We remind the reader that the analytic expression presented by K21b was estimated using fcol = 2.4. Fig. B1 shows the ratio of the time-lag for a given value of fcol divided by the time-lag for fcol = 2.4 as a function of wavelength for a* = 0 and 0.998 (left and right panels, respectively). This figure shows that the time lags increase as fcol increases, by the same factor at all wavelengths (i.e. the ratios are constant as function of wavelength). In addition, these ratios are constant with spin. Fig. B2 shows the ratio of the time lags for fcol≠ 2.4) over the time lags for fcol=2.4. The solid line shows that a straight line gives a very good fit to the data (both for spin 0 and 0.998). This result allows us to update the analytic expression derived by the K21b as follows:
where τ(λ) and τ(λR) are the time lags between X-rays and the bands at λ and λR, respectively (λ1950 = λ/1950 Å, and λR, 1950 = λR/1950Å). In the equation above, h10 is the height of the lamp-post in units of 10 rg, M7 is the BH mass in units of 107M⊙, |$\dot{m}_{0.05}$| is the accretion rate in units of 5 per cent of the Eddington luminosity, LX, 0.01 is the observed, 2–10 keV luminosity in units of 0.01 of the Eddington luminosity. Functions |$A(h_{10}), B(h_{10}), f_1(\dot{m}_{0.05}),$| and f2(LX, 0.01), are given by equations 9–12 and 13–16, for spin 0 and 1, in K21b. As for the fcol, the term,
is the same for both a* = 0 and 0.998. This is the equation for the best-fitting line plotted in Fig. B2.
Ratio of the time lags over the time lag assuming fcol = 2.4 as a function of wavelength for a* = 0 and 0.998 (left and right panels, respectively). The colour code corresponds to the values of fcol.
Average ratio of the time lags over the time lag assuming fcol = 2.4 as a function of fcol, a* = 0 and 0.998 (filled circles and empty squares, respectively). The solid black line corresponds to the best-fitting linear relation to the data points.
B2 Effects of Ltransf/Ldisc
We repeated the simulations as K21b by choosing the same fiducial parameters of |$M_{\rm BH}, \dot{m}/\dot{m}_{\rm Edd}, h, \theta , \Gamma ,$| and Ecut. We also considered two values of the spin, 0 and 0.998, and two values of fcol being 1 and 2.4.
Fig.B3 shows the ratios of the time lags for external power corona divided by the time lags for internal power, for the different combinations of spin and fcol. In all cases, the time lags of externally powered corona are larger than the ones of internally powered corona, as already mentioned in Section 2.1. The difference in shape and amplitude depends on the chosen values of spin and fcol. In the case of a* = 0.998, the time lags difference is smaller than 20 per cent, for any value of fcol. Most probably, such differences cannot be detected with the currently available time-lag spectra. We conclude that equation (B1) above works well, both in the case of externally and the case of the internally powered corona, in the case of high spins.
Ratio of the time lag assuming external power of the X-ray source over time lag assuming an X-ray source powered by the accretion flow as a function of wavelength for a* = 0 and 0.998 (top and bottom panels, respectively) and fcol = 1 and 2.4 (left and right panels, respectively). The colour codes correspond to the values of |Ltransf/Ldisc|.
The upper panels in Fig. B3 show that for Ltransf/Ldisc values smaller than ∼0.4, the time lags difference is smaller than ∼20 per cent. We therefore conclude that, even in the case of zero spin, equation (B1) should provide the correct time lags, irrespective of how the X-ray corona is powered, as long as Ltransf/Ldisc is smaller than ∼0.4. However, this is not the case when Ltransf/Ldisc increases. In this case, the shape of the time-lag spectra changes, and the differences depend on fcol as well. In general, we can see that the time lags can change by a factor of up to ∼1.6 at long wavelengths, when fcol = 1, and by the same factor at small wavelengths, when fcol = 2.4 (and Ltransf/Ldisc = 0.9).
We cannot find a straightforward prescription to add into equation (B1) that could take into account the assumption of an externally or internally powered corona that will be valid for Ltransf/Ldisc values larger than ∼0.4. For that reason, we recommend the use of the code we present in this paper to model the observed time lags, under the desired assumption regarding the source of power for the X-ray corona.
