On the Connection between the Repeated X-Ray Quasiperiodic Oscillation and Warm Absorber in the Active Galaxy RE J1034+396

We conduct an in-depth spectral analysis of ∼1 Ms XMM-Newton data of the narrow line Seyfert 1 galaxy RE J1034+396. The long exposure ensures high spectral quality and provides us with a detailed look at the intrinsic absorption and emission features toward this target. Two warm-absorber (WA) components with different ionization states ( log(ξ/ergcms−1)∼4 and log(ξ/ergcms−1)∼2.5–3 ) are required to explain the intrinsic absorption features in the Reflection Grating Spectrometer spectra. The estimated outflow velocities are around −1400 km s−1 and −(100–300) km s−1 for the high- and low-ionization WA components, respectively. Both absorbers are located beyond the broad-line region and cannot significantly affect the host environment. We analyze the warm absorbers in different flux states. We also examine the 2007 May observation in the low and high phases of quasiperiodic oscillation (QPO). In contrast to previous analyses showing a negative correlation between the high-ionization WA and the QPO phase, we have found no such variation in this WA component. We discover a broad emission bump in the spectral range of ∼12–18 Å, covering the primary features of the high-ionization WA. This emission bump shows a dramatic change in different source states, and its intensity may positively correlate with the QPO phase. The absence of this emission bump in previous work may contribute to the suggested WA–QPO connection.


INTRODUCTION
Outflowing gas from Active Galactic Nuclei (AGNs) acts as a bridge that connects the central black hole (BH) and the host galaxy.It influences the surroundings of the BH and may regulate the environments of the host galaxy (see, e.g., King & Pounds 2015).Characterizing the physical properties of these outflowing winds is crucial for understanding their origins and getting knowl-edge on the growth of BH and its coevolution with the host (e.g., Kormendy & Ho 2013).
Warm absorbers (WAs) are a type of AGN outflow discovered through absorption lines and edges in the soft X-rays (e.g., Halpern 1984;Nicastro et al. 1999;Blustin et al. 2005;Steenbrugge et al. 2005;Krongold et al. 2007;McKernan et al. 2007;Laha et al. 2014).It has been reported that the WAs can be detected in around 50−65% of the nearby type 1 AGNs (Blustin et al. 2005;McKernan et al. 2007;Tombesi et al. 2013;Laha et al. 2014).Considering a typical lifetime of an AGN of 10 7 years, the high detection rate implies WAs to be long-lived and highly-covering outflowing clouds.The kinetic energy of the WAs is usually insufficient to generate significant feedback (e.g., Krongold et al. 2007;Reeves et al. 2009;Zhou et al. Ebrero et al. 2011;Gupta et al. 2013).However, WAs may still impact significantly the host surrounding environment during the entire lifetime (Blustin et al. 2005;Khanna et al. 2016).
The study of WAs remains incomplete ever since the first report of this outflow phenomenon (Halpern 1984).An important debate is on the launching mechanism, with several theoretical models proposed to explain wind generation.The widest accepted ones describe the WAs as thermally-evaporated winds (e.g., Krolik & Kriss 2001;Dorodnitsyn et al. 2008;Mizumoto et al. 2019), radiatively-driven winds (e.g., Proga & Kallman 2004;Dannen et al. 2019), or magneto-driven winds (e.g., Blandford & Payne 1982;Fukumura et al. 2010Fukumura et al. , 2018)).We also know very little about the wind location.The most frequently adopted method calculates the radial distance r through the definition of the ionization parameter.However, the degeneracy between n H and r makes the distance estimation ambiguous.
There are two approaches to determining the density of an outflowing absorber and thus breaking the n H − r degeneracy.The first is to analyze the response of the absorber to the source flux change and determine the recombination timescale of the wind plasma (e.g., Krongold et al. 2005Krongold et al. , 2007;;Kaastra et al. 2012;Khanna et al. 2016;Wang et al. 2022b), which inversely correlates with the gas density (Nicastro et al. 1999;Bottorff et al. 2000).The second approach is through analysis of the metastable absorption lines that are sensitive to the plasma density (e.g., Kaastra et al. 2004;Mao et al. 2017).However, accurately determining the plasma density using either method is challenging.Measuring the recombination timescale requires an apparent change in 1 Here we adopt the definition of ionization parameter following Tarter et al. (1969) as ξ = L ion /(n H r 2 ), where L ion is the 1 − 1000 Ryd source luminosity, n H the gas density, and r the radial distance between the wind and the central illuminating source.
the AGN luminosity, long exposures, and a high plasma density to ensure the recombination timescale is detectable (e.g., Krongold et al. 2007;Kaastra et al. 2012;Silva et al. 2016;Rogantini et al. 2022).As for the analysis of metastable lines, the insufficient effective area and spectral resolution of the current generation X-ray spectrometers make the lines hard to detect.RE J1034+396 is a narrow line Seyfert 1 galaxy at redshift z ∼ 0.0431 and is famous for a repeated one-hour signal of quasi-periodic oscillation (QPO) (Gierliński et al. 2008).The BH mass of this target is between ∼ 10 6 −10 7 M ⊙ (Czerny et al. 2016;Bian & Huang 2010;Chaudhury et al. 2018;Jin et al. 2021).In this work, we adopt a recently estimated mass of M BH = 3 × 10 6 M ⊙ (Chaudhury et al. 2018;Jin et al. 2021).
By analyzing the XMM-Newton EPIC-pn data observed in May 2007, Maitra & Miller (2010) identified an O viii absorption edge at ∼ 0.85 − 1.10 keV, presenting the first evidence for WA in this target.The hydrogen column density and the ionization parameter of the WA constrained by the oxygen edge were N H = (4 ± 1) × 10 21 cm −2 and log(ξ/erg cm s −1 ) = 3.2 +0.2 −0.1 , respectively.The authors also found that the oxygen edge disappeared when the QPO reached the high phase, suggesting the QPO is produced by the periodic obscuring of the WA clouds located at 15 R g , where R g is the gravitational radius.The WA was later confirmed by the narrow absorption lines of O viii, Fe xviii, Fe xix, etc, seen in the XMM-Newton RGS spectrum (Middleton et al. 2011), with N H = 2.23 × 10 21 cm −2 , log(ξ/erg cm s −1 ) = 2.7, and an outflow velocity on the order of 1000 km s −1 .This finding casts doubt on the obscuring scenario since the location of 15 R g is too small to explain the narrow broadening of the lines (∼ 400 − 1000 km s −2 ) seen in the high-resolution spectrum.Jin et al. (2021) studied the May-2007 and Oct-2018 spectra by a similar method applied in Maitra & Miller (2010) and found the equivalent width of the most significant WA feature (Fe xix at ∼ 0.9 keV) negatively correlates with the QPO phase.Therefore, they concluded that the WA is possibly correlated with the QPO phase.
RE J1034+396 is the only AGN for which a correlation between the WA and QPO is hinted.The potential existence of this correlation may open up a new avenue to study the nature of both the phenomena and the vicinity of the central BH.However, no detailed analysis of the WAs in this target has been presented.With more than 1 Ms XMM-Newton observations over the past dozen years, the high spectral quality makes RE J1034+396 one of the best targets to study the WA properties.
In this work, we aim to constrain the basic properties of WAs in RE J1034+396.We also inspect the spectral time variation at around 0.9 keV, which will benefit our understanding of the potential connection between WA and QPO.This paper is organized as follows.In section 2, we introduce the XMM-Newton observations and the data reduction procedures.In section 3, we describe the spectral fitting methods.The fitting result and discussions are presented in section 4. Throughout the paper, we adopt a flat ΛCDM model for all the calculations, with H 0 = 70 km s −1 Mpc −1 , Ω Λ = 0.70, and Ω M = 0.30.
2. OBSERVATION AND DATA RESUCTION RE J1034+396 was observed by XMM-Newton for more than 1 Ms before 2021 May 31.We use the Xray data of both the EPIC-pn and RGS observed in the small-window mode.This yielded 17 available observations, and Table 1 shows the detailed information of each observation.To accurately estimate the ionizing luminosity L ion , we also used the OM data in spectral fitting to constrain the broadband spectral energy distribution (SED).We only considered the UVW1 data of OM since the target was primarily observed by this filter, and the optical emission of RE J1034+396 is dominated by the starlight of the host galaxy (Bian & Huang 2010;Czerny et al. 2016;Jin et al. 2021).
We followed the standard procedures of the XMM-Newton Science Analysis System (SAS v21.0.0) to reduce all the data.The EPIC-pn data were processed using the task epproc.We filtered out bad pixels by setting 'FLAG = 0' and only adopted the single events ('PATTERN = 0') to reduce the pile-up effect.Source spectra were extracted by a 30 ′′ circular region centered on the target, and the background spectra were taken from a 30 ′′ source-free circle of the same chip.The 10.0 − 12.0 keV light curve was applied to filter the soft-proton flares with a threshold of RATE < 0.04 counts s −1 .As for the RGS, we used the task rgsproc to reduce the data.The light curve of CCD 9 of each RGS unit was adopted to exclude the bad intervals with rates above 0.1 counts s −1 .Only the first-order spectra were considered.Similar to Mehdipour et al. (2015), we processed the OM data by the standard task omichain.The derived UVW1 count rates in the source list were later converted into the standard OGIP format through the task om2pha.
Over the 12 years, the target shows two X-ray flux states with a maximum flux change of ∼ 50% at 0.3 − 10.0 keV spectral range.It has been reported The left side of the vertical line shows the RGS data, and the right side is the EPIC-pn data.We make an instrumental calibration between RGS and EPIC-pn, and the EPIC-pn spectra are re-scaled by a factor of 0.922 and 1.054 for the LFS and HFS, respectively.We re-bin the spectra for display purposes.(Bottom) Red and blue show the LFS and HFS light curves at 0.3−10 keV band overlapped with the UVW1 count rates in black dots.
that the QPO of RE J1034+396 can only be detected at the low flux state (Zoghbi & Fabian 2011;Alston et al. 2014).We, therefore, divided the spectra into two catalogs to track the WA properties in different flux states (see Figure 1).We also provide 1 − 4 keV light curves of the May-2007 and Oct-2018 observations for illustration purposes of the QPO in Appendix A. The highflux-state (HFS) spectrum consists of the two brightest observations (obs.ID 0561580201 and 0675440301) with a total net exposure of 77 ks and 48 ks for RGS and EPIC-pn, respectively.The other observations make up the low-flux-state (LFS) spectrum with a total net exposure of 1024 ks for RGS and 687 ks for EPIC-pn.To analyze the weak features with a high signal-to-noise ratio, we co-added all the spectra in each catalog using the tools epicspeccombine and rgscombine for EPIC-pn and RGS, respectively.Signal-to-noise ratio per resolution element at 20 Å is 26 and 11 for the stacked LFS and HFS RGS spectra, respectively.Despite the X-ray flux state, the UVW1 count rates show a small variation (≲ 11%) among the 10 observations when the filter is available.The amplitude of the UV variation is comparative to the 0.3−10 keV flux change of ∼ 10% among the LFS observations.Thus, all the UVW1 data were combined using 1/σ 2 as the weight for spectral fitting of both the flux states, where σ is the uncertainty of the UVW1 count rates.
The EPIC-pn, RGS, and OM spectra were analyzed together.We cross-calibrated the EPIC-pn and RGS to match their flux at a common wavelength band.The EPIC-pn flux was later rescaled by a factor of 0.922 and 1.054 for the co-added LFS and HFS spectra, respectively.In Figure 1, we show the stacked LFS and HFS spectra of RGS and EPIC-pn.The EPIC-pn spectra exhibit some wiggles at ≳ 7 keV spectral range, implying the presence of ultra-fast outflow (UFO, see Appendix B).The source has a very soft spectrum, and the variability between the two states is mainly due to the changes in the soft X-rays.

SPECTRAL ANALYSIS
In this paper, we use the SPEX package v3.07.03 (Kaastra et al. 1996(Kaastra et al. , 2020) ) to fit the time-averaged spectra and search for the optimal model by minimizing the C-statistic (Kaastra 2017).Solar abundance (Lodders et al. 2009) is adopted throughout the paper unless otherwise mentioned.Uncertainties are quoted at 1 σ significance range.To avoid oversampling the data, we re-binned the RGS spectra by a factor of three.The EPIC-pn data were re-binned optimally by the SPEX command obin (Kaastra & Bleeker 2016).

Continuum Modelling
The intrinsic SED consists of a disk blackbody component (dbb), a warm Comptonized disk component (comt), and a power-law tail (pow) (Middleton & Done 2010;Hu et al. 2014;Jin et al. 2021).The optical depth of the Comptonized disk was fixed at τ = 11 according to previous works (Done et al. 2012;Hu et al. 2014;Kaufman et al. 2017;Jin et al. 2021).The seed photons of the warm Comptionization were assumed to come from the inner disk emission (i.e., T dbb was linked to T 0 of comt).We adopted the plasma temperature of the Comptonized disk to constrain the low-energy cut-off for the power-law component, and the high-energy cut-off was set to be 300 keV (Gonzalez et al. 2018;Buhariwalla et al. 2020).We modeled the Galactic neutral absorption with the hot model (de Plaa et al. 2004;Steenbrugge et al. 2005).This model accounts for both absorption lines and edges under collisional ionization equilibrium.We let the gas temperature free to vary but fixed the hydrogen column density at 1.25 × 10 20 cm −2 according to the HI4PI survey (HI4PI Collaboration et al. 2016).
The best-fit continuum model gives C = 2293, with 1049 degrees of freedom (DoF).The dbb component was fixed at the best-fit value for further analysis.The neutral gas of the host galaxy was also considered at the beginning by adding another hot component at the source redshift.However, the improvement of the fitting is insignificant (∆C ∼ 0), and we did not find any neutral absorption lines like O i, O ii, and N i at the host redshift.The best-fit hydrogen column density is N H,host < 4 × 10 18 cm −2 .This value is much smaller than the previous results of N H,host ∼ (1−6)×10 20 cm −2 when only absorption edges were considered (e.g.Done et al. 2012;Jin et al. 2021).Therefore, we excluded this component due to its negligible contribution.

Absorption from Galactic Warm-hot Halo
We identify three absorption features centered at 21.618 Å, 28.773 Å, and 34.946 Å (see Figure 2), which are consistent with the rest wavelengths of O vii Heα, N vi Heα, and C v Heβ, respectively.These absorption lines come from the warm-hot halo (see, e.g., Tumlinson et al. 2017) of the Milky Way (MW) with a temperature of 10 5 −10 7 K.The ion fractions of all three species peak at ∼ 10 5 − 10 6 K under collisional ionization equilibrium (Gnat & Sternberg 2007), implying the warm phase (10 5 −10 6 K, e.g., Sembach et al. 2003 We applied the photoionization model pion (Miller et al. 2015;Mehdipour et al. 2016;Mao et al. 2018) to fit the features from the warm absorber.Pion is a selfconsistent model which calculates the thermal equilibrium, ionization balance, and absorption and emission spectrum of the photoionized gas.All the pion components can be simultaneously fitted with the spectral SED and other model components.The ionization balance is re-calculated in real-time during each iteration when the SED varies.Therefore, there is no requirement for a fixed prior SED as in the classical photoionization models.
Initially, we added one pion component to fit the WA features and set the WA properties to the previous results (Maitra & Miller 2010;Middleton et al. 2011;Jin et al. 2021): hydrogen column density of N H = 3 × 10 21 cm −2 , ionization parameter of log(ξ/erg cm s −1 ) = 3, outflow velocity of v out = −2000 km s −1 , and microscopic motion velocity (i.e., turbulence) of v mic = 300 km s −1 .Fitting this component improves the bestfit model by ∆C/∆DoF = 326/4.This WA component (WA 1 hereafter) takes responsibility for the highly ionized Fe lines at around ∼ 14 Å and part of the O viii line at 19.7 Å (observer frame), which is similar to previous findings (Middleton et al. 2011;Jin et al. 2021).However, absorption features with a lower ionization are not well-modeled by WA 1 , especially for the intrinsic O vii features shown in Figure 3.We thus adopted the second WA component with a lower ionization state (WA 2 hereafter).This process improved the fitting by ∆C/∆DoF = 120/4 (C/DoF = 1788/1039).In our modeling, the covering factor of both the WA components was assumed to be unity (full covering).The unobscured SED was adopted to ionize the WA 1 , while the leaked light from the WA 1 layer was assumed to be the ionizing SED for WA 2 .The two WA components help us recognize all the prominent absorption features in the LFS spectrum.However, some of the absorption features are not well described by the current model, implying additional absorption components.We regarded the third absorption component either as a WA modeled by a pion component or as the warm-hot halo of the host galaxy modeled by a hot component.However, the improvement of the fitting is limited in both cases.We found some emission features in the LFS spectrum overlap with absorption lines.Modeling these emission features may alleviate the imperfect fitting of the absorption lines.

Emission Features
The most significant emission feature in the residual spectrum is a broad emission bump (BEB) at around 12 − 16 Å as shown in the top panel of Figure 4.The BEB covers the primary features of WA 1 and influences the accurate estimation of the WA property.We find the BEB can be well-modeled by a broad Gaussian profile as described by the gaus model in SPEX.The best-fit statistic of applying a broad Gaussian is C/DoF = 1489/1036, with ∆C/∆DoF = 299/3.Resid-ual spectrum after adding the broad Gaussian is shown in the bottom panel of Figure 4. Absorption features of WA 1 , like Ne ix and Fe xvii-xix at around 14 Å, can then be well-modeled after introducing the broad Gaussian component.The BEB is hard to be explained by a photoionized emitting plasma (see Appendix C).Since the primary objective of this paper is to analyze the WA property, we applied this phenomenological model for further analysis.
Finally, we included a pion emission component to model the narrow emission lines from the warm emitter (WE).The unabsorbed SED is assumed to ionize the emission pion component (Mao et al. 2018).The velocity shift of the narrow emission lines is fairly slow, and we thus fixed the outflow velocity at 0 to reduce the fitting complexity.We introduced the macroscopic motion broadening to the emission pion component modeled by a Gaussian broadening model vgau.This broadening refers to the rotation around the black hole and is often degenerated with the microscopic motion velocity (Mao et al. 2018).We fixed the microscopic velocity at v mic = 100 km s −1 and let the emission covering factor (C em ) free to vary within a typical range of 0 − 0.1 (Mao et al. 2022).The narrow WE component improves the fitting by ∆C/∆DoF = 28/4, and the best-fit statistic of the final model is C/DoF = 1461/1032.

Stacked High-flux-state Spectrum
We applied the same continuum model used for the LFS spectrum (i.e., (dbb + comt + pow) * hot).The temperature of the Galactic cold gas was fixed at 3.3 eV according to the best-fit LFS result.The continuum model results in a fitting statistic of C/DoF = 1355/1015.We do not involve the second hot component for the Galactic warm-hot halo and the pion emission component for the warm emitter since the weak features cannot be identified in the low-quality spectrum.A broad Gaussian and two pion absorption components are still adopted to fit the BEB and WA features of the AGN, with the best-fit LFS parameters to be the initial values.This reduces the fitting statistic to C/DoF = 1242/1004 (∆C/∆DoF = 57/4 for WA 1 , ∆C/∆DoF = 16/4 for WA 2 , and ∆C/∆DoF = 40/3 for BEB).

RESULTS AND DISCUSSIONS
We show the optimal model parameters in Table 2. Figure 5 shows the best-fit SED, and Figures 6 and 7 are the best-fit LFS and HFS spectra, respectively.

Intrinsic SED
Between the two flux states, spectral variation focuses on the soft X-ray band and is due to the change of the   et al. 2021).This discrepancy is due to the difference in the disk blackbody model included in SPEX (Shakura & Sunyaev 1973) and XSPEC (see, e.g., Mitsuda et al. 1984).The bolometric luminosity of the AGN derived by the SED model is L bol ∼ 1.2 × 10 45 erg s −1 for both the flux states.This value is slightly higher than the previous estimation of 1.04 × 10 45 erg s −1 (Jin et al. 2012).Adopting a BH mass of M BH = 3 × 10 6 M ⊙ , the Eddington luminosity of RE J1034+396 is L Edd = 3.75× 10 44 erg s −1 using the empirical equation of L Edd = 1.25 × 10 38 × (M BH /M ⊙ ) (Rybicki & Lightman 1986).Thus, the Eddington ratio of RE J1034+396 is estimated to be λ Edd = L bol /L Edd = 3.2, with a range of λ Edd ∼ 1 − 10 due to the mass uncertainty (i.e., 10 6 − 10 7 M ⊙ ).

Absorption from the Milky Way
We have identified absorption features from the neutral and warm-hot gas components of our Galaxy in the RGS spectra.The best-fit temperature of the neutral component is 0.33 +0.08 −0.07 eV (3.8 +0.9 −0.8 × 10 3 K).This component absorbs the spectral continuum at > 10 Å areas and generates the absorption lines like O i, O ii, N i, etc.The warm-phase plasma dominates the Galactic warmhot halo toward the LOS of RE J1034+396.We find a best-fit temperature of 30.7 +3.5  −1.9 eV (3.6 +0.4 −0.2 × 10 5 K), and an optimal hydrogen column density of N H = 3.4 +0.8 −0.7 × 10 19 cm −2 for the warm-phase halo.As for Observed wavelength ( Å) the hot-phase halo (10 6 − 10 7 K, e.g., Fang et  We detect two WA components by analyzing the X-ray spectra of RE J1034+396.The highly ionized component WA 1 has an ionization parameter of log(ξ/erg cm s −1 ) ∼ 4 and an outflow velocity of v out ∼ −1400 km s −1 .The main absorption features of WA 1 are the Fe xvii-xxi lines at around 14 Å (∼ 0.9 keV, observer frame), indicating this is the WA component discovered in previous works (Maitra & Miller 2010;Middleton et al. 2011;Jin et al. 2021).WA 1 also contributes to the highly ionized lines like O vii-viii, and Ne ix-x, etc.The derived N H of WA 1 in this work is more than two times higher than that in Middleton et al. (2011).Though we get a consistent N H with Maitra & Miller (2010), the authors fitted the low-resolution EPIC spectrum and interpreted the Fe absorption lines as the O viii edge.From our high-resolution spectral analysis, our derived ionization parameter is higher than all the previous results with log(ξ/erg cm s −1 ) ∼ 3 (Maitra & Miller 2010;Middleton et al. 2011  WA 2 is a newly discovered component in this work with an outflow velocity of an order −100 km s −1 .The ionization parameter of this component is log(ξ/erg cm s −1 ) ∼ 2.5 − 3.0, and the hydrogen column density is N H ∼ 10 20 cm −2 .WA 2 contributes to the absorption lines with a moderate ionization state like O vii, O viii, C vi, N vii, etc. Comparing the WA 2 with WA 1 , we find the WA component with a higher ionization parameter also has higher hydrogen column density and outflow velocity.These correlations have been observed in many other WA targets (e.g., Tombesi et al. 2013;Laha et al. 2014).The two WA components in RE J1034+396 manifest a correlation of v out ∝ ξ 0.4 and v out ∝ ξ 1.2 in the LFS and HFS, respectively.The significant change on the index of v out − ξ correlation can not be explained by the pure case of radiatively-driven wind (v out ∝ ξ, King & Pounds 2003) or magneto-driven wind (v out ∝ ξ 0.5 , Fukumura et al. 2010).The derived v out −ξ relation is also far from the observed statistical relation of v ∝ ξ 0.12±0.03(Laha et al. 2014) and individual cases in Laha et al. (2016) and Wang et al. (2022a).This result implies a complex launching mechanism for the WAs (see, e.g., Laha et al. 2016;Wang et al. 2022a).Alternatively, the two WAs in RE J1034+396 may be separately generated at different times, locations, or with distinct launching mechanisms.

Radial Location and Mass Outflow Rate
We estimate the upper and lower distances of the WA winds to the central BH.The upper limit is a geometrical constraint that the thickness of the absorber can not exceed its distance to the central BH (i.e., ∆r/r ≤ 1).Combined with the definition of ionization parameter ξ = L ion /n H r 2 and N H = n H ∆r, we have r max = L ion /N H ξ. The lower boundary is obtained by assuming the outflow velocity is larger than the es- cape velocity at r. Thus, r min = 2GM BH /v 2 out , where M BH = 3 × 10 6 M ⊙ is the mass of the central BH, G is the gravitational constant.Substituting the best-fit quantities into the equations, the radial location of WA 1 is estimated to be 0.013 − 5.7 pc by the LFS spectrum and 0.014 − 9.1 pc by the HFS spectrum.As for the WA 2 , the corresponding values are 0.24 pc − 12.8 kpc and 5.3 pc − 1.6 kpc constrained by the LFS and HFS data, respectively.Considering the results of both the flux states, the radial location of WA 1 is estimated to be 0.014 pc < r WA1 < 5.7 pc (∼ 2 × 10 4 − 1 × 10 7 R g ).The location of WA 2 is constrained within 5.3 pc < r WA2 < 1.6 kpc (∼ 9 × 10 6 − 3 × 10 9 R g ).
To illustrate the WA location more intuitively, we roughly calculate the locations of the broad-line region (BLR), dusty torus, and narrow-line region (NLR).The size of BLR is estimated by r BLR = 39.08 × [λL λ (5100 Å)/(10 44 erg s −1 )] 0.518 light-days (Bentz et al. 2006).Adopting the luminosity derived in Bian & Huang (2010) of λL λ (5100 Å) = 5.6 × 10 42 erg s −1 , the size of BLR is 0.007 pc.The inner and outer radii of the torus are calculated according to Nenkova et al. (2008), where r in = 0.4 × (L bol /10 45 erg s −1 ) 0.5 × (1500 K/T d ) 2.6 pc, and r out < 30r in .Applying a dust sublimation temperature of T d = 1500 K, we have 0.44 pc < r dust < 13 pc.As for the NLR, its size is estimated by r NLR = 2.1(L [O III] /10 42 erg s −1 ) 0.52 kpc (Netzer et al. 2004).An O iii luminosity of 5.47×10 40 erg s −1 (Bian & Huang 2010) leads to a NLR size of 463 pc.We compare the distances of the WAs to those of BLR, torus, and NLR in Figure 8.Both WA components lie outside the BLR.WA 1 is located within the outer boundary of the torus, while WA 2 is likely to be located between the inner boundary of the torus and the NLR edge.
The mass outflow rate and the kinetic energy of the WA winds can then be calculated after the location estimation.For a uniform spherical outflow, the mass outflow rate can be estimated by Ṁout = 4πr 2 n H m p v out ≤ 4πrm p N H v out , where m p is the proton mass.Adopting an outflow velocity of −1400 km s −1 , a hydrogen column density of 10 21.75 cm −2 , and a maximum location of 5.7 pc for the WA 1 , the mass outflow rate of WA 1 is Ṁout < 4.6 M ⊙ yr −1 .As for the WA 2 , the mass outflow rate is Ṁout < 8.9 M ⊙ yr −1 estimated using v out = −300 km s −1 , N H = 10 20.38 cm −2 , and r max = 1.2 kpc.The kinetic energy of the WA wind is then ĖK = 1/2 Ṁout v 2 out < 2.9 × 10 42 erg s −1 for WA 1 and ĖK < 2.4 × 10 41 erg s −1 for WA 2 .Comparing the WA kinetic energy to the bolometric luminosity of RE J1034+396, we find ĖK /L bol < 0.24% for WA 1 and ĖK /L bol < 0.02% for WA 2 , respectively.The AGN feedback models predict a wind ĖK /L bol ratio of > 0.5% to be an efficient feedback (e.g., Hopkins & Elvis 2010).The low kinetic energy indicates that the WA winds in RE J1034+396 can not significantly affect the environment of the host galaxy.According to our best-fit models, WA 1 maintains its properties between the two flux states.In contrast, the ionization parameter of WA 2 seems to decrease when the source flux state becomes lower.It is against intuition that the outer plasma WA 2 ) responds to the flux change of the central illuminating source without any variation on the inner plasma (i.e., WA 1 ).However, the estimated WA properties of the HFS may be inaccurate due to the low-quality data.A longer exposure of the HFS is required to resolve the WA variation between the different flux states.
The previously discovered WA-QPO connection is based on the analysis of the absorption features at around 14 Å (Maitra & Miller 2010;Jin et al. 2021), indicating that it is the WA 1 potentially exhibits a connection with the QPO phase.It has been reported that the QPO in RE J1034+396 can only be found at the LFS (e.g., Zoghbi & Fabian 2011;Alston et al. 2014).However, our analysis suggests that there might be no difference in the properties of WA 1 between the two flux states, regardless of the existence of the QPO.This finding casts doubt on the previously discovered connection between the WA and QPO.
We discover a broad emission component centered at ∼ 14 − 15 Å in both states.The full width at half maxima of BEB is ∼ 0.14 − 0.24 keV (∼ 45, 000 − 90, 000 km s −1 ).The BEB shows significant variation in the strength, central wavelength, and broadening between the two states.We attempted to model the BEB in the stacked LFS spectrum using an emitting pion component (Appendix C).The result shows that the BEB could potentially be relativistically broadened features of highly ionized Ne or Fe.However, both explanations require a negligible abundance of oxygen and a very low abundance of metals like C and N.This kind of photoionized emitting plasma has never been discovered in the AGN vicinity.The nature of BEB is still unclear and needs a further check.
The BEB covers the spectral range of ∼ 12 − 18 Å as shown in Figures 6 and 7. Since the previous works were not aware of the BEB contribution, they made an inaccurate estimation on the continuum level at around 14 Å for the absorption line analysis.Considering the remarkable variation of BEB in the two flux states, it is likely that WA 1 does not change its properties in different QPO phases, but the BEB variation leads to the change in the measured equivalent width of the absorption feature.If this is indeed the case, the BEB flux at around 14 Å should have a positive correlation with the QPO phase.4.4.2.QPO-phase-resolved Spectra of the May-2007

Observation
We analyzed the QPO-phase-resolved spectra of the May-2007 observation, in which the WA-QPO connection was first reported.We used the 0.3 − 10 keV EPICpn light curve to determine the good time intervals of the high QPO phase (HQP) and low QPO phase (LQP) using the method applied in Maitra & Miller (2010).In short, we labeled the flux extremums in each QPO period, and the high QPO intervals were determined when the flux decreases 40% from the peaks to the nearby troughs, while the low QPO intervals were determined when flux increases 40% from the troughs to the nearby peaks.Similar to Maitra & Miller (2010), observational data behind 84 ks are excluded due to soft-proton flares.
We extracted the QPO-phase-resolved spectra of RGS and EPIC-pn using the derived phase intervals.The net exposures of the HQP and LQP are 32.6 ks and 29.3 ks, respectively.We re-binned the RGS spectra by a factor of six due to the poor data quality.The RGS and EPIC-pn spectra combined with the U V W 1 count rates were jointly analyzed.The spectral continuum was fitted by (pow+comt+dbb)*hot, with the Galactic neutral gas absorption fixed at the best-fit stacked LFS value.Absorption from WA 1 and emission from BEB were fitted by a pion component and a Gaussian profile, respectively.We did not model the features from Galactic warm-hot gas, WA 2 , and warm emitter due to their negligible contribution to the low-quality spectra (∆C ∼ 0).The best-fit model parameters and 10 − 18 Å RGS spectrum are shown in Table 3 and Figure 9, respectively.
The two-phase spectra exhibit a significant difference in the power-law component, with the HQP showing a steeper spectrum and a larger power-law normalization.This finding is consistent with previous timing analysis, which suggested the power-law component is likely to take responsibility for the periodic flux change (Middleton et al. 2009;Zoghbi & Fabian 2011).We detect the WA 1 in both QPO phases without any obvious difference, in contrast to the previous results showing a non-detection of WA 1 in the HQP spectrum (Maitra & Miller 2010;Jin et al. 2021).The BEB seems stronger (in normalization) and broader in the HQP than the LQP, implying a positive correlation between BEB intensity and the QPO phase.Therefore, the dramatic variation of BEB and the unawareness of BEB contribution in the early works may lead to the wrong impression of the WA-QPO connection.The correlation between BEB and QPO is not assertive due to poor data quality.A detailed analysis of the BEB nature and its timing properties will be presented in our future work on the target.

SUMMARY
In this work, we conduct a detailed analysis of WA winds in RE J1034+396 with more than 1 Ms XMM-Newton observations.We analyzed the properties of WA winds in low and high flux states as well as the May-2007 observation in low and high QPO phases.Our main findings are summarized as follows.
1. Two WA components are required to explain the intrinsic absorption features in the time-averaged RGS spectra.The highly ionized component, which has been discovered before (Maitra & Miller 2010;Middleton et al. 2011;Jin et al. 2021), has an ionization parameter of log(ξ/erg cm s −1 ) ∼ 4, an outflow velocity of around −1400 km s −1 , and contributes to the transitions like O vii-viii, Ne ix-x, and Fe xvii-xxi.The lower ionized component is newly discovered in this work.This component shows an outflow velocity of around −(100 − 300) km s −1 and a ionization parameter of log(ξ/erg cm s −1 ) ∼ 2.5 − 3. It primarily models the absorption features of O vii-viii, N vii, and C vi.
2. The highly ionized WA is likely to be located within the outer boundary of the torus but outside the BLR.While the low-ionization WA may be located between the inner boundary of the torus and the NLR.The kinetic energy of the wind is estimated to be < 0.24%L bol and < 0.02%L bol for the high-and low-ionization WAs, respectively, which is unlikely to produce a significant impact on the host galaxy.
3. We find no difference in the highly ionized WA between the two flux states and the two QPO phases, which is against the WA-QPO connection discovered in previous works.The ionization parameter of the low-ionization WA may positively correlate with the source flux.The suspected correlation may be biased by the poor data quality of the HFS and should be further explored with deeper observations.
4. We identify a broad emission bump at around 14 Å, which covers the primary features of the high-ionization WA.This component shows a significant variation between the two flux states, and its strength may positively correlate with the QPO phase.The dramatic variation of the broad emission bump may take responsibility for the misidentified WA-QPO connection in the early works.
5. The cold and warm-hot halos of the MW contribute to the absorption features seen in the RGS spectra.The best-fit temperatures of the two components are ∼ 3.8 × 10 3 K and ∼ 3.6 × 10 5 K, respectively.
We expect that our study will benefit future analyses on both the absorption and emission phenomena of the target.By conducting QPO-phase-resolved spectroscopy with more than 1 Ms XMM-Newton observations, future studies may provide a better understanding of the BEB nature and its potential correlation with the QPO, which will empower us to investigate the vicinity of the supermassive black hole with deeper insight.7)-( 8) ∆Cstat and ∆AICc with and without the line, and (9) signal-to-noise ratio of the line.∼ 9.0 keV in the stacked LFS spectrum.We adopted a power law to fit the spectral continuum and used three Gaussian profiles to fit the absorption lines.Table 4 lists the best-fit parameters of the absorption lines.We estimated the line significance using the corrected Akaike information criterion (∆AIC c ) between models with and without the interested line (Burnham et al. 2011).Two of the lines show moderate significance at 2 σ level, while the third is considered insignificant.We calculated the line velocity shift assuming an origin of Fe xxv Heα or Fe xxvi Lyα.The calculating result shows no overlapping of the line velocities, suggesting the three lines do not come from the same UFO component.All three lines are insignificant in the stacked HFS spectrum.The UFO is implied but cannot be confirmed by the current analysis.
C. MODELING THE BEB AS A PHOTOIONIZED PLASMA Centered at around 14 Å, the BEB may be relativistically broadened Ne or Fe lines, which possess strong features at this wavelength band.Based on the best-fit stacked LFS spectrum, we replaced the Gaussian profile with a pion emission component to model the BEB.When examining the Ne case, we fixed the Ne abundance at the solar value but let C, N, O, and Fe abundances free to vary.As for the Fe case, we fixed the Fe abundance to solar value, and the abundances of C, N, O, Ne, Mg, Al, Si, and S were treated as free parameters.The best-fit 10 − 18 Å spectra are displayed in Figure 12.
The best-fit Ne case results in a fitting statistic of C stat /DoF = 1470/1026.This explanation requires an ionization state of log(ξ/erg cm s −1 ) ∼ 1.5 with a velocity of ∼ 3000 km s −1 (inflow) and a broadening of ∼ 20000 km s −1 .
Zhou et al.The best-fit abundances of O and Fe are near zero, avoiding the generation of strong emission features unseen in the observed spectrum (see the orange curve in Figure 12 for an example).The best C and N abundances are around 0.1 solar.
The Fe case achieves a slightly better statistic of C stat /DoF = 1457/1022.In this case, the ionization state is similar to the high-ionization WA value of log(ξ/erg cm s −1 ) ∼ 4.0.This explanation requires a velocity of ∼ −9000 km s −1 (outflow) with a broadening of ∼ 30000 km s −1 .The best-fit abundances of N, O, Ne, Mg, Al, and S are near zero (0.012 solar abundance for oxygen), and C and Si abundances are 0.65 and 0.5 solar, respectively.The impact of oxygen abundance is less pronounced in this case, and setting the oxygen abundance to 0.3 solar only marginally increases the continuum level at around 19 Å.However, the best-fit oxygen abundance of 0.012 solar indicates that the Fe abundance is 83 times the oxygen abundance.
The weird abundance composition is unusual and has never been discovered before for the photoionized plasma.Moreover, the best-fit statistics of the two cases are worse than the original model using a simple Gaussian profile (C stat /DoF = 1461/1032).We did not find a photoionized explanation using a narrow component (i.e., multiple narrow Fe emission lines).Currently, the BEB nature remains uncertain and requires future investigations.

Figure 1 .
Figure1.(Top) Co-added 0.3 − 10.0 keV X-ray spectra of the low flux state (red) and the high flux state (blue).The left side of the vertical line shows the RGS data, and the right side is the EPIC-pn data.We make an instrumental calibration between RGS and EPIC-pn, and the EPIC-pn spectra are re-scaled by a factor of 0.922 and 1.054 for the LFS and HFS, respectively.We re-bin the spectra for display purposes.(Bottom) Red and blue show the LFS and HFS light curves at 0.3−10 keV band overlapped with the UVW1 count rates in black dots.

Figure 2 .
Figure 2. Absorption lines from Galactic warm-hot halo.The black curve is the observed data with 1 σ uncertainty shown in blue.The red curve is the best-fit hot model.The Galactic C v Heβ line is overlapped with the AGN absorption features.

Figure 4 .
Figure 4. Spectral residual without modeling the BEB (top) and modeling the BEB by a broad Gaussian component (bottom).The main absorption features of WA1 are labeled in red.Residuals of the absorption features disappear when BEB is introduced.

Figure 6 .
Figure 6.Time-averaged 10 − 36 Å LFS RGS spectrum with the best-fit model.The top panel of each segmental spectrum shows the transmission and emissivity of the absorption and emission components, respectively.Colors red and blue label the WA and WE components.Color green marks the Galactic absorption.The bottom panel shows the observed data in black, with 1 σ uncertainty in faint blue.All the prominent absorption features are labeled.

Figure 7 .
Figure 7. Same Figure 6, but for the time-averaged HFS RGS spectrum.

Figure 8 .
Figure 8. Distances of each structure away from the central BH.The blue and red lines show the nH −R solutions of WA1 and WA2 constrained by the LFS (solid) and HFS (dotted) spectra.The vertical lines mark the locations of BLR, torus, and NLR, respectively.

4. 4 .
Time Variation of WA and BEB 4.4.1.WA and BEB in the Two Flux States

Figure 11 .
Figure11.Stacked 4 − 10 keV EPIC-pn spectra of the LFS (left) and HFS (right), respectively.Blue points are the observed data.The orange curve shows the background level.We fit the spectral continuum by a power law.For the LFS spectrum, we also use three Gaussian profiles to fit the potential UFO features.The red curve is the best-fit model.The vertical green shadows in both panels highlight the potential UFO features seen in the stacked LFS spectrum.

Table 1 .
XMM-Newton observation log.The first two columns are the observational date and ID.The total and net exposures of EPIC-pn and RGS after filtering out the soft protons are listed in columns (3) to (5).Column (6) shows the 0.3 − 10.0 keV source count rates (correlated by the background) observed by the EPIC-pn.Columns (7) and (8) are the available OM filter and the U V W 1 count rates.The last column labels the X-ray flux state of each observation.At the bottom of the table, we show the exposure of each combined spectrum.
Intrinsic O vii and O viii absorption lines in the co-added LFS spectrum.The black curve shows the observed spectrum with 1 σ uncertainty in blue.The red-dashed and red-solid lines denote the best-fit modes with one and two WA components, respectively.
;Savage & Wakker 2009;Tumlinson et al. 2017;Qu et al. 2020) primarily contributes to the lines.Another hot model was used to fit the Galactic warm-hot lines, with the gas temperature a free parameter varied between 10 5 K and 10 7 K.The best-fit result, as shown in Figure2, improves the fitting by ∆C/∆DoF = 59/2 (C/DoF = 2234/1047).The local C v Heβ line is overlapped with the AGN intrinsic features.3.1.3.Warm Absorber

Table 3 .
Best-fit models of the QPO-phase-resolved spectra.
Figure 9. 10 − 18 Å RGS spectra of the LQP (left) and HQP (right) of the May-2007 observation.The top panel shows the model decompositions, with the BEB in blue and the high-ionization WA in red.The bottom panel shows the observed data in black and the 1σ uncertainty in blue.The red curve shows the best-fit model.

Table 4 .
Properties of the potential UFO lines in the stacked LFS spectrum.Columns are (1) observed line energy, (2) energy at the AGN rest frame, (3) full width at half maximum, (4) equivalent width, (5)-(6) outflow velocity if Fe xxv Heα line or Fe xxvi Lyα line, ( Figure 12.Stacked LFS BEB spectra fitted by the relativistically broadened Ne (left) or Fe (right) lines.The top panel only shows the decomposition of BEB and WA1 for clarity.The red curve in the spectrum panel is the best-fit model.The orange curve is the best model but with 0.3 times the solar oxygen abundance.