© ESO 2022

1. Introduction

Many massive galaxies tend to host super-massive black holes (SMBHs) at their centers (e.g., Genzel & Townes 1987; Genzel et al. 1996; Eckart et al. 2001; Kormendy & Ho 2013; Gaspari et al. 2019). These SMBH usually grow through accreting gaseous material and are visible as “Active Galactic Nuclei” (AGN). Whether or not these AGN have an important contribution in driving the evolution of their host galaxies had been an unanswered question until AGN feedback was shown to be essential for matching the galaxy luminosity function in models to observations.

Benson et al. (2003) examined the physics determining the shape of the galaxy luminosity function in a cold dark matter (CDM) universe with a high baryonic content (Ω_b = 0.02). They stated that although the cooling and condensation of gas in a merging hierarchy of dark matter halos could be the dominant process behind galaxy formation, that is not sufficient to explain the sharp cutoff of the galaxy luminosity function and scarcity of highly luminous galaxies. They stated that the sharp cutoff of the galaxy luminosity function can be achieved if the cooling of gas in these massive halos could be efficiently suppressed by (1) the reheating of cold disk gas; (2) energy injection from supernovae and stellar wind heating up and expanding the diffuse gas halo; or (3) complete expulsion of cold gas from the disk by highly energetic superwinds. They also added that superwinds can only expel gas if the formation of central SMBHs and the generation of the associated energy result in limiting the star formation in their host galaxies. Bower et al. (2006) and Croton et al. (2006) agree that the observed scarcity of very bright and massive galaxies can be explained by taking into account the energy released from AGN, which could effectively prevent the gas from cooling to form the unobserved, highly luminous galaxies.

These results motivated other theoretical studies to analyze the effect of AGN feedback on galaxy growth (Churazov et al. 2005; Somerville et al. 2008; McCarthy et al. 2010; Hopkins et al. 2016; Schaye et al. 2015; Gaspari & Sądowski 2017). Most of the galaxy evolution models involved a form of energy injection by the AGN in which the AGN drives galactic scale outflows and hence expels the gas out of their host-galaxies, prohibiting BH growth. This is known as the radiative feedback mode (e.g., Benson et al. 2003; King 2003; Di Matteo et al. 2005; Hopkins et al. 2006; Hopkins & Elvis 2010; Booth & Schaye 2010). The relativistic jets from the AGN could also heat up the surrounding environment and obstruct star formation. This is known as the mechanical feedback mode (e.g., McNamara & Nulsen 2007; Dubois et al. 2010; Bîrzan et al. 2012; Bower et al. 2012). In this feedback mode, the accretion of material onto the SMBH leads to the production of highly relativistic, plasma jets (Best & Heckman 2012). The mechanical feedback is suspected to be responsible for heating of the intra-cluster medium (ICM) and quenching of cooling flows in the centers of clusters (e.g., McNamara & Nulsen 2007; Gaspari et al. 2020). Even though star formation processes such as supernovae and stellar-winds drive galactic scale outflows (e.g., Heckman et al. 1984; Lehnert & Heckman 1996; Swinbank et al. 2009; Roberts-Borsani & Saintonge 2019), AGN are necessary to drive the outflows of observed velocities > 1000 km s⁻¹ (Benson et al. 2003; McCarthy et al. 2011).

In recent years, X-ray and UV spectroscopy demonstrated that a considerable fraction of the highly accreting AGN could drive ultrafast (velocity ∼0.1c) outflows (UFOs) close to their accretion disks (Blustin et al. 2003; Reeves et al. 2003; Ganguly & Brotherton 2008; Tombesi et al. 2013; Serafinelli et al. 2019). But, in order to drive the host galaxy’s evolution, one key requirement is that the AGN needs to drive galactic scale outflows as mentioned in several galaxy evolution models (Granato et al. 2004; Erb et al. 2006; Oppenheimer & Davé 2008; Hopkins & Elvis 2010; Lilly et al. 2013). The optical, highly ionized [O III] line is commonly used to trace the signs of AGN ionization over galaxy scales. As [O III] is a forbidden line, it cannot trace the dense broad emission line region (BLR) kinematics, but the kinematics in the narrow line region (NLR) which extends from pc to few tens of kpc scale (e.g., Husemann et al. 2016, 2019a). Therefore, any kinematic disturbance in the NLR would be well reflected in the [O III] line shape over the galactic scales.

Crenshaw & Kraemer (2000) suggested that the outflow away from the nucleus of NGC 1068 could follow a biconical kinematic model. They found that the similar amplitudes of the blueshifted and redshifted curves on the northeast side in the position-velocity diagram. This infers that the axis of the bicone is in vicinity to the plane of the sky. The lack of low radial velocities where the curves peak could indicate that that the bicone is evacuated along its axis.

Bae & Woo (2016) also found that a simple biconical model of radial outflow is sufficient to explain the general trend of observed outflow velocities. In addition, turbulent motions can also contribute to the broadening of the spectral lines (e.g., Tremblay et al. 2018; Simionescu et al. 2019; Gaspari et al. 2020).

Mullaney et al. (2013) used the SDSS spectroscopic data to study what drives the kinematics of the kilo-parsec scale [O III]λ5007 emitting gas. They analyzed a sample of 24 264 z < 0.4 optically selected type 1 and type 2 AGN and performed a multicomponent fitting to the optical emission-line profiles to estimate the [O III]λ5007 kinematics. For type 2 AGN, the line of sight is toward the central engine, intercepting the equatorial dust. On the other hand, the central engine will be directly visible for a type 1 AGN (Urry & Padovani 1995). The transition between the type 1 to type 2 usually occurs between inclination angle 45 − 60° (Marin 2014). Mullaney et al. (2013) reported that both type 1 and type 2 AGN exhibit asymmetry in [O III] line profile, which they interpreted as evidence of outflowing ionized gas. They suggested that the ionized outflows are prevalent in both type 1 and type 2 AGN and are driven by compact radio jets. Previous studies by Heckman et al. (1981), Whittle (1992), Blundell & Beasley (1998), Thean et al. (2001), Jarvis et al. (2019), Molyneux et al. (2019) also concluded that the compact radio sources could be responsible as the launching mechanism of these ionized outflows. Woo et al. (2016) analyzed a sample of ∼39 000 type 2 AGN at z < 0.3 and showed that ∼45% of type 2 AGN at z < 0.3 show signs of outflows in [O III]. They reported an increase in [O III] velocity dispersion as AGN luminosity and Eddington ratio increase. Woo et al. (2016) concluded that the radiation pressure from the AGN accretion disk could launch the ionized outflows.

Although the one-dimensional spectroscopic data provide insight into the ionized gas kinematics, they cannot provide any information about the spatial extent or structure of these outflows. This is where spatially resolved spectroscopy has significant advantage. Over the last two decades, both long-slit and integral field unit (IFU) observations of AGN have found evidence of high-velocity and kinematically disturbed ionized gas over several kpc scales (e.g., McCarthy et al. 1996; Colina et al. 1999; Swinbank et al. 2005; Westmoquette et al. 2012; Riffel et al. 2013; Diniz et al. 2015; Humire et al. 2018; Powell et al. 2018; Slater et al. 2019; Kakkad et al. 2020; Scholtz et al. 2020). Previous work (Liu et al. 2013a; Harrison et al. 2014; McElroy et al. 2015) found that the AGN activity resulted in high-velocity and disturbed ionized gas (velocity dispersion up to 500 km s⁻¹) extended to ∼10 kpc in their targets; and the majority (70%) of their AGN (z < 0.2) population exhibited highly ionized [O III] regions which extended to kpc scales with high velocity (∼510–1100 km s⁻¹) ionized gas.

Villar-Martin et al. (2016) investigated the similarity of the [O III] light distributions and point-spread functions (PSF) with very large telescope (VLT) FORS2 observations of type 2 AGN to show that the outflowing ionized [O III] intensity distribution follows a point-like profile. Husemann et al. (2016) reported that the bright point-like emission from the center could scatter out and confuse the extended flux. This is known as beam smearing. Recent studies have reported sub-kpc scale outflows (e.g., Villar-Martin et al. 2016; Tadhunter et al. 2018; Kawaguchi et al. 2018; Bellocchi et al. 2019; Baron & Netzer 2019). Therefore, the spatial extensions of these outflows still remains as an open question. Fortunately, if the [O III] wing surface-brightness profile mimics the emission profile of a theoretical point source, it is possible to perform a 2D modelling (Gadotti 2008) to estimate the center of the [O III] wing emission as well as the location of the central-engine (SMBH). Furthermore, one can achieve the precision of a hundredth of a pixel based on high spatial resolution data using spectro-astrometry (for a discussion, see Bailey 1998; Whelan & Garcia 2008; Gnerucci et al. 2010, 2011). Therefore, using type 1 QSOs, we aim to compare the [O III] light distribution to that of the PSF (the light distribution of the Hβ broad line region, BLR), and if they are similar, to determine how far off is the region of the maximum [O III] wing emission from the SMBH. Similar methods have been applied by Kakkad et al. (2020), Santoro et al. (2018).

Moreover, the energetics associated with the ionized outflow is strongly dependent upon the value of the electron densities in the NLR (which can take values from 10–10⁴ cm⁻³ (Nesvadba et al. 2006; Liu et al. 2013b; Harrison et al. 2014). The electron density is responsible for the largest fraction of systematic uncertainty when estimating the outflow-energetics for warm, ionized gas outflows (Müller-Sánchez et al. 2011; Kakkad et al. 2016; Perna et al. 2017). It therefore remains unclear how the kinetic power of the outflowing gas couples with AGN properties or the star-formation in the host galaxy. Hence, it is important to understand how the electron density varies on the different scales in AGN host galaxies and how it is affected by the warm, ionized outflows. Liu et al. (2013b) assumed an electron density of n_e ∼ 100 cm⁻³ at 1 kpc away and n_e < 10 cm⁻³ at 6 kpc away from the nucleus for the outflowing ionized gas. Studies by Bennert et al. (2006), Sharp & Bland-Hawthorn (2010), Cresci et al. (2015), Freitas et al. (2018) have reported that the electron density in the central region can go up to a few 10³ cm⁻³, whereas it drops to ∼10² cm⁻³ at ∼1.2 kpc away from the nucleus. Recent spatially resolved studies by Kakkad et al. (2018), Baron et al. (2020) and X-shooter observations from the study by Davies et al. (2020) reported a similar decrease for the electron densities of the outflowing gas with increasing distance from the central region. A systematic estimation of electron densities is hence required to properly quantify the outflow energetics.

In this paper we explore optical IFU observations of a large number of nearby, luminous type 1 QSOs from the close AGN reference survey (CARS, Husemann et al. 2017, 2019b, 2022). We focus on characterizing the light distribution for the asymmetric [O III] wing component with respect to the PSF. Thereby we can systematically verify whether such putative outflows are unresolved (at the limit of our observations) or possibly extended over kpc scales. In addition, we also measure electron densities associated with the outflowing wing component based on the [S II] λλ6717, 6731 doublet to test previous assumptions incorporated into outflow energetics calculations.

Throughout this paper we adopt the standard ΛCDM cosmology with H₀ = 70 km s⁻¹ Mpc⁻¹, Ω_m = 0.3, and Ω_Λ = 0.7.

2. The IFU data set

For this study we use the IFU data obtained as part of CARS. CARS represents an ambitious spatially resolved multi-wavelength survey of a representative sample of the most luminous type 1 AGN at redshifts 0.01 < z < 0.06. The CARS targets are drawn from the Hamburg-ESO Survey (HES, Wisotzki et al. 2000), which is a purely flux-limited AGN catalog based on B-band photometry and slitless spectroscopy. A random subsample of HES AGN in this redshift range have follow-up in the cold gas phase (Bertram et al. 2007; König et al. 2009) which defines the CARS targets. It is important to note that none of the sample selection criteria are related to ionized gas outflow properties in any sense, which makes this sample unbiased and fully representative for the outflow properties of luminous local AGN. All details of the sample can be found in Bertram et al. (2007), Husemann et al. (2022). IFU observations have been obtained with multi unit spectroscopic explorer (MUSE) (Bacon et al. 2010, 2014) for 36 objects (programs 094.B-0345(A) and 095.B-0015(A), PI: B. Husemann), visible multi object spectrograph (VIMOS) (Le Fèvre et al. 2003) for HE 1310−1051 and HE 1338−1423 (program 083.B-0801(A), PI: K. Jahnke), and the potsdam multi-aperture spectrophotometer (PMAS, Roth et al. 2005) at the Calar Alto Observatory for HE 0853−0126 and HE 0949−0122 (program H18-3.5-010, PI: B. Husemann).

The MUSE data were reduced by the standard ESO pipeline for the instrument (Weilbacher et al. 2020) while the VIMOS and PMAS data were reduced with the Py3D fibre-fed IFU reduction package developed initially to reduce the data of the CALIFA survey (Husemann et al. 2013). For more technical details on the IFUs and the data reduction process we refer the reader to Husemann et al. (2019b, 2022). Apart from the basic CARS IFU data processing as used by Husemann et al. (2019b) and explained in Husemann et al. (2022), a stellar continuum and emission-line modelling of the 3D IFU data is performed after QSO-subtraction using PYPARADISE (Husemann et al. 2016, 2019b). For the analysis presented in this paper we make use of the stellar continuum model datacubes. The original MUSE and QSO-subtracted data will be available at CARS data release 1 (DR1). In our analysis we neglect the CARS target HE 0021−1810 due to its poor data quality comprising just a single short MUSE exposure, which is not sufficient for robust measurements. Furthermore, PMAS IFU has a sampling size of 1″ which is clearly not as precise as MUSE/VIMOS. Therefore, we excluded the 2 AGN from the PMAS IFU in this paper.

3. Analysis and results

3.1. Modelling the [O III] line shape

A primary signature of ionized gas outflows from AGN are the commonly observed broad and blue-shifted wings in the [O III] emission lines. In order to extract the integrated [O III] emission line profile parameters, we co-added all the spectra within a 3″ diameter around the AGN location. As our sample consists of only type 1 AGN, the broad Hβ and Fe IIλλ4923,5018 lines from the broad line region (BLR) are blended with the narrow emission lines of Hβ and [O III] λλ4959, 5007, respectively. We modeled the spectral region from 4750 Å to 5090 Å (rest-frame wavelength) simultaneously to deblend all broad and narrow emission lines in this region. This is a standard process for the analysis of type 1 AGN spectra (e.g., Jin et al. 2012; Mullaney et al. 2013; Husemann et al. 2016, 2019b).

Here, we assumed a simple superposition of Gaussian profiles to model all lines in the spectrum. We adopt 1–2 Gaussian components to model each of the narrow [O III] λλ4959, 5007 doublet and the narrow Hβ component. The first Gaussian component corresponds to the core and the second Gaussian component corresponds to the broad, usually blue-shifted, wing of the lines which might be potentially attributed to ionized gas outflows. We further deployed 1 to 2 Gaussian components for each of the BLR Hβ and Fe IIλλ4923, 5018 emission lines depending upon their complexities. Finally, we approximated the local AGN continuum with a first-order polynomial given that the considered wavelength range is small.

Additionally, we always assume that the absolute velocity and the intrinsic velocity dispersion of the narrow [O III] λλ4959, 5007 and Hβ components are coupled, because they likely originate from the same physical region with similar kinematics. Similar coupling is assumed for the BLR components of the Hβ and Fe IIλλ4923, 5018 emission lines. We further put two constraints on the line fluxes in addition to our previous assumptions of emission line kinematics. One is that the [O III] λλ4959, 5007 has a line flux ratio of 1/3 (Storey & Zeippen 2000). Secondly, we found that the line flux ratio of Fe IIλ4923 and Fe IIλ5018 was nearly constant across the sample with a mean ratio to 0.81; hence we used 0.81 as Fe IIλ4923/Fe IIλ5018 to constrain the line fluxes of the Fe II doublet.

We independently fit the brightest AGN spaxel and the integrated 3″ spectra using a nonlinear Levenberg-Marquardt algorithm. For the error estimation, we create 100 different spectra by fluctuating the original spectrum with the noise spectrum and re-model the spectra to obtain the associated uncertainties for each parameter as the standard deviation of the repeated modelling. For the two Gaussian components of [O III], we only consider the second, broad Gaussian component when it has a signal-to-noise ratio of S/N > 5 and the fractional error in the velocity dispersion in each of the Gaussian components is < 1. All sources in our sample satisfy these criteria. An example of the best-fit AGN spectral model is shown for HE 1126−0407 in Fig. 1.

Fig. 1. Example of the multicomponent modelling for the spectrum of HE 0345+0056. Upper panel: brightest spaxel spectrum of HE 0345+0056 and the best-fit spectral model. The blue green lines refer to the narrow core and the blue lines refer to the broad wing components of Hβ + [O III]. The gray lines refer to the individual emission lines of the BLR Hβ and Fe II lines. The narrow width of the BLR component representing Hβ is due to the fact that HE 0345+0056 is a narrow line Seyfert 1 (NLS1), where the broad line region tends to show narrower Balmer lines. These narrow BLR components are common in these objects (e.g., Mathur et al. 2000; Mullaney & Ward 2008). Lower panel: residual spectrum normalized by the error spectrum. The dark gray shaded area highlights the ±3σ limit and the red dashed line provides the reference to 0.

From the best-fit emission-line model of the central 3″ × 3″ spectra, we report the [O III] line fluxes, their velocity difference Δv = v_wing − v_core and the [O III] line width (W₈₀) in Table 1. Any negative Δv will represent a blue-wing and a positive Δv will denote a red-wing [O III] asymmetry. W₈₀ is defined as the velocity width of the emission-line that contains 80% of the line flux. If v₁₀ and v₉₀ are the velocities at the 10th and 90th percentiles, then W₈₀ = v₉₀ − v₁₀ (Harrison et al. 2014). In addition, we also report the total BLR Hβ flux (f_HβBLR).

In Fig. 2, we compare the [O III] line width described by W₈₀ with two often used bolometric luminosity indicators, the total [O III] luminosity and the broad Hβ luminosity. The CARS sample is consistent with the bulk population of AGN as probed by the SDSS (Mullaney et al. 2013) and is therefore a representative subsample of the low-redshift AGN population in this parameter space. Several other IFU studies of SDSS-selected AGN instead focused on more extreme systems.

Fig. 2. [O III] line width described by W80 against AGN bolometric [O III] luminosity (left panel) and BLR Hβ luminosity (right panel). The yellow circles denote the type 2 AGN sample from Liu et al. (2013a), the blue diamonds denote the type 2 AGN sample by Harrison et al. (2014), the lime triangles denote the type 2 AGN sample by Sun et al. (2017) and the red stars denote our CARS sample of type 1 AGN. The dense scatter plot corresponds to the overall SDSS sample of both type 1 and type 2 (left panel) and type 1 AGN (right panel) from the ALPAKA library published by Mullaney et al. (2013).

The asymmetry in the [O III] line-shape does not necessarily correspond to an outflow: turbulent motions (e.g., Gaspari et al. 2018) and simple gravity may also drive the broadening and related red/blue wings in the [O III] line-profile toward the center below certain velocity limits. Nevertheless, the velocity offsets are seen at the galaxy center co-located with the AGN, which indicates that this asymmetric feature in the [O III] line profile is potentially related to an AGN-driven outflow as the dominant process. Previous detailed IFU studies of ionized gas outflows (e.g., Liu et al. 2013a,b; Harrison et al. 2014; McElroy et al. 2015; Sun et al. 2017) targeted rare subsamples at the edges of the distribution, with significantly higher W₈₀ and L_[O III] values. Those are therefore not representative for the overall AGN population and their results cannot be necessarily extrapolated to be common properties of the overall AGN population.

3.2. Spatial distribution of the [O III] wing component

The kinematics associated with the [O III] broad wings in the ionized gas can be easily obtained from the AGN spectra for large samples as described in the previous section. However, the spatial distribution and the size of the wings can only be determined from spatially resolved spectroscopy. In recent years, numerous IFU studies mapped the surface brightness and the kinematics of the ionized gas mainly around type 2 AGN. While several studies reported ionized gas outflows on kpc scales (e.g., Liu et al. 2013a; Harrison et al. 2014; McElroy et al. 2015; Wylezalek et al. 2020), recent work highlighted that ionized gas outflows can also be quite compact on 10–100 pc scales (e.g., Villar-Martin et al. 2016; Tadhunter et al. 2018; Storchi-Bergmann & Schnorr-Müller 2019). Those diskrepancies may be caused either by different AGN sample characteristics or in the treatment of the data.

Due to the high gas density and large ionizing flux close to an AGN, the [O III] emission is usually centrally concentrated and can outshine the extended [O III] core emission from the NLR on larger scale depending on seeing and redshift (e.g., Husemann et al. 2016). Hence, it is important to characterize the point-spread function (PSF) for an observation to unambiguously separate compact from extended [O III] emission. The main advantage of using type 1 AGN is that the broad Balmer emission lines (i.e., Hα, Hβ) allows us to reconstuct the PSF, with the assumption that the BLR emission is intrinsically unresolved (e.g., Jahnke et al. 2004).

We take advantage of the type 1 AGN in CARS to accurately characterize the PSF. In the following, we describe the analysis used to test whether the broad wing [O III] kinematic components are compact or extended at the limits of our seeing-limited spatial resolution and constrain the actual location using spectroastrometry. Therefore, we start with the hypothesis that the light distribution of the [O III] wing component is point like and centrally located. Our hypothesis will be diskarded in case the light distribution is spatially extended implying large-scale ionized gas outflows.

3.2.1. Constructing the 2D light profiles of emission lines

Following our hypothesis that the light distribution associated with all the emission lines is compact and that a possible extended component is faint with respect to its central counterpart, we fixed the emission line kinematics and radial velocities to the ones derived from the central brightest pixel modeling as described in Sect. 1. Only the line fluxes and the continuum shape parameters were allowed to vary within the entire 3″ around the AGN position. For VIMOS, we allowed the window to be 6″ as the pixel size and spatial resolution is worse than MUSE. Thereafter, the line fluxes for each spaxel were obtained using a nonlinear Levenberg-Marquardt algorithm and we reconstructed the 2D intensity maps for each line separately. We further created 100 different data cubes by fluctuating the main cube with the associated error cubes and repeated the measurements 100 times to obtain proper error maps for the emission line fluxes implicitly taking cross-talk between the lines into account. The resulting emission line maps are exemplarily shown for HE 0040−1105 in the top panels of Fig. 3.

Fig. 3. Central 3″ intensity maps associated with different spectral lines. The leftmost map describes the light distribution from Hβ broad emission-line region, which is the PSF. The bottom-right map describes the light distribution from blue-shifted ionized [O III] wing; and the top-right map describes light-distribution from the bright, narrow core of [O III].

3.2.2. Modeling the Point-Spread Function

We assume a 2D Moffat function (Moffat 1969) to model the PSF following the work by Gadotti (2008). We use a nonlinear Levenberg-Marquardt algorithm to model the observed light profiles. The reference PSF is derived by modelling the BLR Hβ light distribution, where all Moffat model parameters (flux amplitude, centers, PSF shape parameters, and ellipticity) are free parameters.

3.3. Classification of compact and extended outflows

If the spatial distribution of the [O III] wing is extended, it will not follow a PSF. We assume an extended spatial distribution for the [O III] wing component if the deviation from the PSF model is 5σ. While the [O III] emission line profile and gas kinematics can change over kpc scales, such changes would be rather smooth as each kinematic component itself necessarily follows a PSF in the spatial domain. Those kinematic changes would therefore translate into amplitude changes of the wing and core components even at fixed kinematics. Those amplitude changes will manifest themselves as deviations from a pure PSF. Here, we do not consider any [O III] kinematic information beyond the central region and distinct line shape components on several kpc-scale would be missed. However, the main purpose of our approach is to characterize the 2D spatial distribution of the [O III] wing emission close to the nucleus within the characteristic radius implicitly set by multi-object fibre spectroscopy, like SDSS.

We modeled the 2D light distributions of the broad wing (and also the narrow core) of [O III] by varying the position and amplitude normalization of the initially determined Moffat PSF. All other parameters were kept fixed because the PSF shape should not significantly change within the considered narrow wavelength range. We further fluctuated the 2D intensity maps with the error maps and repeated the entire process 100 times to estimate the uncertainty associated with the measurement of the Moffat parameters. Afterwards, we constructed the residual maps normalized by the errors (R_norm) as shown in Fig. 4.

Fig. 4. Results of our spectro-astrometric analysis for HE 934+0119 from the MUSE IFU observations. We have mapped the surface brightness within the central 3″ for all spectral components in the brightest spectrum assuming fixed kinematics. Here we only show the measure light distribution from the BLR Hβ (left), broad [O III] wing (middle), and narrow [O III] core (right) component. Those 2D light distributions are fitted with a 2D Moffat profile. The red cross indicates the Moffat center for the best-fit elliptical Moffat model to the surface brightness map of the BLR Hβ. The green and red crosses indicate the Moffat centers for the best-fit elliptical Moffat models to the surface brightness maps of the [O III] wing and core, in which the amplitude and position of the Moffat model were varied only. The histograms show the distributions of the normalized residuals (residual/error) for both the core and wing components in comparison to the BLR Hβ one.

We distinguish between unresolved and resolved [O III] wing emission based on the normalized residuals of the PSF fit. As the unresolved emission will be dominant close to the nucleus, we focus on the central 2″ region around the AGN for the MUSE sources, and central 4.5″ region around the AGN for the VIMOS sources given the larger PSF and poor sampling. In order to quantify the similarity between the normalized residual distributions of BLR Hβ and [O III] wing, we perform a Levene’s test (Levene 1960) within the aforementioned region. The Levene’s test checks only the difference between two distributions, but not their relative widths. The p-value from Levene’s test alone is not necessarily sufficient to conclude whether the [O III] wing emission is extended or not. Consequently, we further calculate the standard deviations of the normalized residuals distributions and compare σ_[O III] and σ_Hβ within the central region around the nucleus. When σ_{[O III:Hβ]} = σ_[O III]/σ_Hβ < 1, we consider the [O III] wing to be unresolved, because the residual distribution for the [O III] wing is consistent with the PSF at the given noise level. This assumption only holds if the 2D Moffat PSF model is a good model for the empirical PSF as inferred from the BLR Hβ 2D light distribution. In combination with the Levene’s test we define two necessary criteria for unresolved [O III] wing emission as follows:

1. p > 0.05

2. σ_{[O III:Hβ]} < 1 if p ≤ 0.05.

In Fig. 4, we show the surface brightness profiles for BLR Hβ, [O III] wing and cores for HE 09340+0119, along with their Moffat models and the normalized residuals. We find that the distribution of the normalized residuals for the [O III] wing are fully consistent with that of the BLR Hβ given a p-value of 0.27 for the Levene’s test. It confirms that the outflow appears to be point-like at the limit of the MUSE data. Therefore, the BLR Hβ and [O III] wing components are consistent with a point source with a slight offset in position by ∼27 pc. Comparing the normalized residuals of the [O III] core component to that of the BLR Hβ we report a p-value ∼10⁻⁷, indicating a completely different distribution than a point source.

We report 23 sources with compact [O III] wing (unresolved or point-like) and 13 sources with extended [O III] wing emission possibly extended to kpc scale (see Fig. 5). Furthermore, we calculate the total [O III] wing flux from its light distributions or surface brightness map (F_aperture) and the Moffat model (F_Moffat) within 2″ region around the Moffat center. F_aperture gives us an idea how much emission from the broad wing component is contained within the central 2″ region; whereas, F_Moffat demonstrates how much emission would be contained if the emission was entirely due to unresolved or point-like emission. We define the flux ratio of the [O III] wing emission within the central 2″ to the [O III] wing point-like emission as r_[O III] = F_aperture/F_Moffat. The advantage of using this parameter is that it gives an indication of how significant the emission is in a region compared to a pure point source assumption, where r_[O III] = 1 for an ideal point-like source. For these 23 AGN, we calculate how far is the region of the region of maximum [O III] wing emission from the AGN central engine (S_[O III]). The distance is a lower limit of the size of the outflowing region, while we cannot exclude that the outflow extends beyond this distance. We present the characterization of the [O III] wing emission for all CARS sources in Table 2.

Fig. 5. Results from our analysis: Out of 36 CARS targets, 23 agree with a point-like broad [O III] wing emission profile, whereas the other 13 targets have broad [O III] wing emission profile that is not point-like. We also show one source from each of compact (HE 1011−0403) and extended (HE 1248−1356) [O III] wing emitting sources. We further plot the green contours showing the [O III] wing emission-profiles, along with the black contours of the BLR Hβ (PSF) emission profile.

In Fig. 6, we plot the [O III] luminosity against the flux ratio, r_[O III]. For the 23 sources with unresolved [O III] wing emission, the ratio is between 0.9–1.05. 3 sources (HE 1107−0813, HE 2128−0221 and HE 2211−3903) appear to have less emission within their surface brightness map than the Moffat model (see Appendix D). This is caused by the fact that the surface brightness in every region does not mimic a 2D Moffat profile, which introduces the difference in total flux values. But the surface brightness profile is so close to the Moffat profile that the residual emission fluctuates within the noise range.

Fig. 6. r[O III] against [O III] luminosity for the sources with compact and extended [O III] wing emission. The dashed line represent the region where Faperture = FMoffat.

The 13 sources with extended [O III] wing emission show diverse emission profiles. HE 0212−0059 and HE 1353−1917 have around 15–40% more flux, clearly showing the extended [O III] emission that cannot be described by the Moffat model. 3 sources (HE 0021−1819, HE 0119−0118 and HE 0429−0247) with extended wing emission seem to have < 10% [O III] flux within the 2″ diameter around their peak emission. We note that these sources have the [O III] emission toward a certain direction in contrast to a purely spherically symmetric emission profile. While estimating the total flux for these sources, the excess [O III] wing emission compared to the Moffat PSF in some regions is often compensated by the other region which has less flux than the Moffat PSF. Therefore, the total flux from the [O III] wing emission is within 10% of a pure point-like emission.

On the other hand, the 4 out of 13 sources with extended [O III] wing emission have up to 10% less flux in the central 2″ compared to a purely unresolved or point-like emission profile. This feature is most prominently seen in HE 1338−1423, which exhibits a large bi-polar outflow (Husemann et al. 2014). The outflow in HE 1338−1423 is extended in the southeast and northwest direction (see Appendix D). Although, the central 3″ diameter region has less flux than a pure point source (which explains the flux deficiency), there is an excess [O III] wing emission beyond the central 3″ visible till the edge of the central 4″ window suggesting that the ionized gas might be extended to several kpc scales.

3.4. Systematic error estimation for offset measurements

Previously, we reported that 23 AGN in our sample are consistent with spatially unresolved [O III] wing emission. The 2D Moffat centroid for the [O III] wing emission was offset from the central engine by < 100 pc, with a median offset of ∼27 pc. Although we estimated the measurement error for the 2D offset by repeating the process several times with fluctuated data, we also need to consider the systematic uncertainty of positional measurements. The systematics include other instrument effects like systematic flat-fielding residuals, geometrical distortion of the CCD, or detector noise, which would artificially change the flux distribution at different wavelength and thereby affect the measured flux-weighted centroids.

In order to estimate the systematic positional uncertainty, we construct an artificial source based on the AGN spectrum model and the best-fit 2D PSF model of the BLR Hβ and add it in random locations within the MUSE FoV. We excluded regions with (a) a star, and (b) the presence of bright AGN continuum and verify if S/N < 5 in all spaxels within a 3″ window.

We put the artificial source on 300 arbitrary locations within IFU FoV following the aforementioned criteria and carried out a spectro-astrometric analysis to determine the angular offset between the known and measured offset for the outflowing ionized [O III]. The standard deviation of the angular offset distribution gives us an indication of systematic error in estimating the offsets, which we present in Table 2.

3.5. Estimation of the electron densities

With the sensitivity of MUSE we can robustly detect the optical emission lines [S II]λλ6716,6731 ([S II] doublet). The doublet lines occur due to two different energy levels with very similar excitation energies. Their relative fluxes are dependent upon the electron density that occupies these energy levels. In this scenario, any degeneracy that could arise because of different ions is removed (Osterbrock & Ferland 2006). One key issue is that, the measurement of the electron density is temperature dependent. But, the uncertainty that will be introduced with the assumption of an electron temperature, T_e = 10⁴ K, will be lower than the measurement uncertainty of the electron densities (Sanders et al. 2016). Additionally, the [S II] doublet is sensitive only to a limited range of electron densities, which is well within the typical range of electron densities observed in the Extended Narrow Line Region (ENLR). The ENLR is extended over a few kpc scales where densities in the diffuse ISM are at the level of 100 cm⁻³ (e.g., Harrison et al. 2014; Husemann et al. 2016).

To infer the electron densities of the core and wing components we deblend the [S II] doublet accordingly. We model the doublet within the rest-frame wavelength range from 6680 Å to 6750 Å which covers any possible blue-shifted component. Because the [S II] has a lower S/N than other emission lines (i.e., [O III], [N II]), it is difficult to obtain enough signal from a single pixel. The majority of the outflows in our sample agree with an unresolved source, so that the majority of the emission is confined within the central 1″ diameter around the AGN. Therefore, we co-add all the spectra of 3 × 3 spaxels around the AGN to improve the S/N. Similar to our Hβ + [O III] + Fe II spectral modelling, we model the [S II] doublet with kinematically coupled multi-Gaussian components. For the underlying continuum we adopt a 2nd order polynomial to capture the tail of the Hα BLR emission. Uncertainties are again estimated by fitting 100 mock spectra fluctuated within their noise and by calculating the standard deviation of the resulting model parameter distributions. An example of the best-fit AGN spectral model for HE 0253−1641 along with the data and the residual is shown in Fig. 7.

Fig. 7. Example of the multicomponent modelling for the spectrum of HE 0253−1641. Upper panel: integrated 3 × 3 spaxel-spectrum of HE 0253−1641 and the best-fit spectral model. The green lines refer to the narrow core and the blue lines refer to the broad wing components of [S II]. Lower panel: residual spectrum normalized by the error spectrum. The dark gray shaded area highlights the ±3σ limit and the red dashed line provides the reference to 0.

We use the PyNeb package (Luridiana et al. 2015) and assume a gas temperature, T = 10⁴ K for all calculations to convert the [S II] doublet flux ratio to electron densities. Since the wing component of the [S II] doublet has much lower S/N, we can only reliably determine the electron density for 9 sources among which 6 sources have unresolved [O III] wing emission and 3 other sources have extended [O III] wing emission. We present our finding in Table 3.

Following the work of Osterbrock & Ferland (2006) it is evident that as the electron density decreases, the flux ratio of the [S II] doublet approaches the upper limit of 1.45, whereas if the electron density increases, the flux ratio of [S II] doublet approaches the lower limit of 0.45. For the core component we find a flux ratio of ∼1.0 indicating an electron density ∼500 cm⁻³.

We derive a flux ratio of ∼0.7 for the wing component. Such a flux ratio corresponds to a much higher electron density > 10³ cm⁻³. Hence our results suggest that the electron densities for the wing [S II] component are systematically higher than the bright, narrow core component.

4. Discussion

4.1. Origin of majority of the broad [O III] wing emission

4.1.1. The range in the extent of turbulent circum-nuclear ionized gas

Numerous studies used HST/ACS imaging and spatially resolved spectroscopy to quantify the extent of turbulent ionized gas potentially related to AGN-driven outflows. Whereas some studies report outflows extending over several kpcs (e.g., Liu et al. 2013a,b; Harrison et al. 2014; McElroy et al. 2015; Karouzos et al. 2016a; Kang & Woo 2018), others have reported sub-kpc scale outflows (e.g., Villar-Martin et al. 2016; Tadhunter et al. 2018; Kawaguchi et al. 2018; Bellocchi et al. 2019; Baron & Netzer 2019). This paper focuses on the location of the [O III] wing emission as an indication of turbulent gas motion driven either by outflows (e.g., Crenshaw et al. 2010; Zhang 2021) or other processes (e.g., Tremblay et al. 2018; Gaspari et al. 2018). We use the term outflow to describe this turbulent gas region. Although we cannot physically confirm for individual objects that AGN outflows are the main driver of the turbulence, the outflow interpretation is motivated by other studies. Crenshaw et al. (2010) found that the NLR gas has a velocity offset > 50 km s⁻¹ in the [O III] line-shapes; where 35% of the quasars exhibit blue wings and 6% of the quasars show red wings. They concluded that these red and blue wings could be due to the mass outflow in ionized gas and dust extinction in the inner disk. Venturi et al. (2018), Singha et al. (2021) stated that the broad [O III] wings could indicate outflow and be related to AGN activities.

In order to determine the location of the [O III] wing emission, we compare the 2D light distribution of the [O III] wing with the empricially constrained PSF from the BLR emission approximated by a 2D Moffat analytic function. When the distribution is consistent with the Moffat PSF, we can determine the spatial locations of [O III] wing emission with high precision using spectro-astrometry. The location of this turbulent gas region might be considered as the current base of an outflow. When the [O III] wing emission is spatially extended, estimating a centroid is meaningless and requires a careful PSF debelending to recover the instrinsic distribution and kinematics. This is beyond the scope of the paper and will follow-up on this in future work.

For the compact outflows we can compare the results for our unbiased type 1 AGN sample with previous studies. Different methods have been used to quantify outflow size. We study the difference in inferred sizes by comparing our spectro-astrometric results with other methods. Here, we use the approach taken by Karouzos et al. (2016a) as an example for various other works based on seeing-limited observations and showcase the order of magnitude difference in outflow size/location that can result for the [O III] wing emission when the impact of beam smearing is not properly considered.

In Fig. 8, we plot the offset between the region of the maximum [O III] wing emission (S_[O III]) against L_[O III] for the unresolved outflows cases in our sample. In comparison, we also show measurements of the outflow sizes R_[O III] as reported from past works. It is important to note that S_[O III] is not exactly the same parameter as R_[O III] as defined, e,g., in Liu et al. (2013a), Harrison et al. (2014), McElroy et al. (2015), Karouzos et al. (2016a), Kang & Woo (2018) because S_[O III] estimates the distance between the maximum [O III] wing emitting region and the nucleus; whereas R_[O III] estimates the extent of the region dominated by outflowing ionized gas. Still, it is useful to compare both quantities as it highlights the impact of beam smearing on ground-based observations. Here, we also compute R_[O III] following the approach of Karouzos et al. (2016) for the compact outflow cases. The resulting R_[O III] is order of magnitude larger than S_[O III]. Intriguingly, the kpc-scale values appear fully consistent with the mild L_[O III] dependence despite the fact that we know the emission is unresolved. An explanation for this could be that past studies did not check for consistency with the PSF and therefore could not distinguish unresolved from resolved emission components. Below we diskuss a number of different approaches that have been used to estimate outflow sizes in the literature in the context of the spectro-astrometric method used in this paper.

Fig. 8. Spatial extension of the [O III] wing emission vs. [O III] luminosity for our samples along with outflow sizes obtained by other studies in the literature at similar or higher [O III] luminosity range that find outflows of several kpc scales. For the sources with compact [O III] wing emission, we assume the offset sizes obtained by spectro-astrometry to indicate the how far is the peak of the emission from the nucleus. As previous studies considered the [O III] wing to indicate outflows (e.g., Crenshaw et al. 2010; Zhang 2021), we further plot the outflow sizes measured using the approach taken by Karouzos et al. (2016a), Kang & Woo (2018). We connect the outflow-sizes obtained using the method of Karouzos et al. (2016a) and the offset-sizes using spectro-astrometry for the CARS sources with compact [O III] wing emission to quantify the difference in determining spatial extension using these two methods.

Liu et al. (2013a) analyzed a sample of 11 type 2 radio-quiet QSOs at z ∼ 0.5 from the sample of Reyes et al. (2008). They reported that [O III]/Hβ ∼ 10 in the central region until it reaches a ‘break radius’ (R_br) where it starts to decrease, reaching the ‘matter bounded regime’. Liu et al. (2013b) used R_br as the outflow radius (D) in order to calculate the outflow energetics. They estimated R_br to be ∼4–11 kpc for all the QSOs. Liu et al. (2013a) obtained a flux-calibrated [O III] and Hβ surface brightness map by collapsing the continuum subtracted IFU data cube over the [O III] and Hβ wavelength range. The drawback of their approach is that it could contain not only the contribution from the outflowing ionized gas but also from non-outflowing photo-ionized gas on several kpc scales. To visualize the effect of their prescription, we perform a similar analysis on HE 0040−1105 – one of the AGN with clearly compact [O III] wing emission. We notice that the [O III]/Hβ∼ 10 up to 2 kpc where it sharply falls of to [O III]/Hβ∼ 6 and stays almost constant up to 8 kpc where it falls to even lower value. Contrary to this kpc scale extension, we show that the emission profile of the [O III] wing follows the PSF very closely; and the peak of maximum [O III] wing emission is ∼100 pc away from the nucleus which is about a magnitude lower than the ‘break radius’.

Harrison et al. (2014) analyzed a sample of 16 z < 0.2 type 2 AGN. They reported the high velocity (W₈₀ > 600 km s⁻¹) ionized gas are prevalent over a distance, D₆₀₀, which they defined as the maximum projected spatial extent of the broad [O III] wing emission profile, where the outflows dominate the gas kinematics. They found that all of their AGN showed signatures of spherical outflows and bi-polar super-bubbles and are extended to ∼6–16 kpc. Harrison et al. (2014) concluded that the ionized outflows are not only ubiquitous in luminous (L_[O III] > 5 × 10⁴¹ erg s⁻¹) type 2 AGN, but also are extended over several kpc scales in ≥70% (3σ confidence) of cases. However, they selected QSOs with significantly broader [O III] line profiles than the bulk of the local AGN population as shown in Fig. 2. Hence, their sample is not comparable to our sample given the selection criteria and their respective scopes.

McElroy et al. (2015) studied a sample of 17 type 2 nearby (mean z ∼ 0.08), luminous QSOs (average L_[O III] ∼ 6 × 10⁴² erg s⁻¹). They defined the spatial extent of the outflow to be the extent of the region in which two-Gaussian components are needed to fit the spaxels. They reported that the high velocity (W_80, max ∼ 400–1600 km s⁻¹ with mean ∼790 ± 90 km s⁻¹) ionized gas is indicative of the AGN-driven outflows. They also reported that these outflows are extended to 3–16 kpc scale around the nucleus. The downside of their approach is that if the emission profile of [O III] closely follows the PSF, then the kinematics will roughly stay constant close (a few kpc) to the nucleus. For HE 0040−1105, we find the W₈₀ stays almost consistently ∼400 km s⁻¹ in the central 2.5 kpc (in diameter) region, although the [O III] wing emission clearly originates from a region only 100 pc away from the nucleus. This yields again an order of magnitude overestimation compared to the results obtained from our spectro-astrometric analysis. McElroy et al. (2015) approximated the outflow size as the extent of the region within which the [O III] spectra in all spaxels needed two-Gaussian components to be fitted. Again, if the emission profile closely follows the PSF, observing the spectra at several kpc away from the nucleus does not indicate that the kinematically disturbed gas has already traveled that distance. For example, in HE 0040−1105, the surface brightness maps in Fig. 3 indicates that the emission of both the core and wing components extends ∼2.5 kpc away from nucleus; although the [O III] wing emission originates from a region only ∼100 pc away from the nucleus as inferred from the PSF comparison analysis. We find that the outflow sizes obtained for the sources with point-like [O III] wing emission is significantly larger than the offsets we obtained using the spectro-astrometric analysis.

Karouzos et al. (2016a) analyzed a sample of 6 local (z < 0.1), luminous (average L_[O III] > 10⁴² erg s⁻¹) type 2 AGN. They reported the wing component to have a large velocity (< 600 km s⁻¹), reaching velocity dispersion < 800 km s⁻¹, suggestive of AGN driven outflows. In order to estimate the extents of the outflows, they used the effective radius (r_eff), defined as the radius within which half of the total [O III] wing emission flux is contained. They assumed that the outflow sizes are to be the distance at which- (1) the [O III] wing component’s velocity and velocity dispersion become similar to that of the [O III] core component and (2) the extreme kinematics of the [O III] wing component is prevalent. Their findings suggest that the outflow sizes are roughly ∼2–3 r_eff (1.3–2.1 kpc). Following the method by Karouzos et al. (2016a) we estimate an outflow size of ∼2.5 kpc for HE 0040−1105.

Kang & Woo (2018) analyzed a sample of 23, z < 0.2, luminous (L_[O III] > 3 × 10⁴⁰ erg s⁻¹) type 2 AGN. They traced the radial decrease of the [O III] velocity dispersion and estimated the size of the warm, ionized gas outflows. Kang & Woo (2018) defined the outflow sizes to be the distance at which [O III] velocity dispersion becomes same as the host galaxy’s stellar velocity dispersion. In order to take into account beam smearing, they subtracted the seeing sizes from the kinematically derived sizes of the outflowing regions and estimated the size of the outflows to be between 0.60–7.45 kpc. Furthermore, they obtained a relation between [O III] luminosity and outflow radius as-

$$\begin{array}{cc}& log\left(\frac{R_{\mathrm{out}}}{\mathrm{kpc}}\right)=(0.279\pm 0.035)\times log\left(\frac{L_{[\text{O}\text{III}]}}{{10}^{42}{\mathrm{erg}\mathrm{s}}^{-1}}\right)\hfill \\ & +(0.427\pm 0.025).\hfill \end{array}$$(1)

Following the size-luminosity relation from the work of Kang & Woo (2018) we estimate an outflow size of ∼3 kpc for HE 0040−1105. This overestimation occurs as the emission from the unresolved [O III] wing covers ∼3 kpc although the region of emission is only offset by ∼100 pc from the nucleus. Although they subtracted the seeing from the kinematically derived outflowing sizes, they could not compare the [O III] wing emission to the PSF because their samples consist of type 2 AGN. Subtracting the seeing sizes from the total sizes of the region, does not ensure that the beam smearing effect has been completely taken out from the measurements particularly for spatially unresolved emission.

Other work uses analytic biconical models based on the kinematics of double-peaked [O III] to directly constrain the geometry and energetics of outflows. Nevin et al. (2018) used follow-up longslit observations to investigate the outflows in sample of 18 AGN from SDSS at z < 0.1 with L_bol ∼ 10⁴³ − 10⁴⁵ erg s⁻¹. The used a parameter called the turnover radius, at which the outflowing gas starts to decelerate in the best fit biconical model, as an indicator of the size of the outflows. They found that 2 out of 18 AGN have outflow sizes < 1 kpc and the other 16 AGN have outflow sizes ranging from 1.1–6.5 kpc, where the outflows retain enough kinetic energy that they could drive a two-staged feedback process in their host galaxies.

The major drawback of these studies is that the outflow sizes were calculated from an averaged property of the two-Gaussian models, from the extent of the emission of the wing component, or from an analytic model based on longslit data. They did not test the [O III] emission profile against the PSF as their sample consisted of type 2 AGN. Their results are valid for the AGN in which the [O III] emission profile does not follow the PSF. For the AGN, where the [O III] emission profile follow the PSF, beam smearing (Husemann et al. 2016) could make the emission profiles appear to be extended to several kpc scales, although in reality the region of maximum [O III] emission could be just a ∼ a few hundred pc away from the nucleus. Therefore, if the emission profile is not tested against that of the PSF, a significant overestimation in outflow sizes could occur. Furthermore, they studied highly luminous quasars as seen from Fig. 2. The findings from those studies will not account for the general nearby quasar population as outlier samples of AGN could have extreme properties.

The situation changes for various studies which either use high-angular resolution HST observations or apply PSF deblending techniques to infer outflow sizes. In Fig. 9, we plot the offset sizes for the sources in our sample with compact [O III] wing components against L_[O III], along with previous studies finding compact outflows, usually < 1 kpc (Tadhunter et al. 2018; Fischer et al. 2018). These studies used high resolution HST long-slit spectroscopic (STIS) and imaging (ACS and WFCP2) data in order to investigate the AGN driven outflow sizes and gas kinematics. Furthermore, these studies took account of the dominant point like emission from the PSF. Their outflow sizes are in much better agreement with the sizes reported in this paper and support even the compact outflows found by our approach. Below we briefly diskuss some studies individually to diskuss their approaches.

Fig. 9. Spatial extension of the [O III] wing emission vs [O III] luminosity for our samples along with outflow sizes obtained by Villar-Martin et al. (2016), Tadhunter et al. (2018), Fischer et al. (2018) at similar or higher [O III] luminosity range that find compact outflows. The lower and upper limits of the error bars indicate the lower and upper limits of the outflow sizes obtained by Villar-Martin et al. (2016).

Villar-Martin et al. (2016) analyzed 18 luminous type 2 AGNs at 0.3 < z < 0.6 using VLT-FORS2 spectroscopic data to investigate the presence of the extended ionized gas outflows. They reported that 15 out of 18 AGN in their sample show a signature of ionized outflows. They took account of the beam smearing effect and stated that among the 15 AGN with ionized outflows, the [O III] wing emission profile follows the emission profile of a standard star very closely, and hence indicate unresolved [O III] wing emission. They report a lower limit of outflow sizes ∼several × 100 pc in these AGN. In 3 other AGN, where the [O III] wing emission did not follow that of the standard star (or the seeing), they report sizes < 1 − 2 kpc. They also estimated an upper limit of the outflow sizes to be ∼2 kpc for the majority of the AGN. The only exceptions were SDSS J0955+03 and SDSS J0903+02 where they reported upper limits of up to 3.5 and 10 kpc respectively. They did not explore the location of the maximum [O III] wing emission, and therefore could not provide an estimation of the offset from the nucleus.

Tadhunter et al. (2018) investigated the AGN-driven outflows in 8 ULIRGs at z < 0.15 with optical AGN nuclei with narrow-band [O III]λ5007 imaging data from the deep Hubble Space Telescope (HST ACS) observations. They further complimented the HST observations with the ground-based spectroscopic observations from William Herschel Telescope (WHT). Tadhunter et al. (2018) found that the outflowing regions dominated by the warm, ionized gas are usually compact for the majority of the AGN in their sample. They found that the ionized gas outflows are barely resolved at the resolution of HST (65–120 pc) in 3 objects. HST spatially resolves the outflows in 4 other objects but their flux weighted mean radii (R_[O III]) fall in the range 0.65 < R_[O III] < 1.2 kpc. In one other object (Mrk 273), they found conclusive evidence for a large scale outflow, which has a maximum extent of R_[O III] ∼ 5 kpc. To account for the beam smearing effect, they estimated the FWHM for the star (FWHM_star) and the [O III] emission profile (FWHM_[O III]) by fitting 2D Gaussians to the HST/ACS images and then subtracted FWHM_star from FWHM_[O III] in quadrature to estimate R_[O III]. We see that, in almost half of the objects the outflows are unresolved with size of ∼100 pc which is very similar to what we report.

Mingozzi et al. (2019) studied the interstellar medium (ISM) properties of the disk and the outflowing ionized gas in the central regions of 9 nearby (D_L < 50 Mpc) Seyfert galaxies using VLT/MUSE observations. They stated that the redshifted and the blueshifted wings of the emission lines are outflows whereas the narrow core represents the disk of the host galaxies. Although they did not take into account the beam smearing effect, which could affect the analysis of 2 objects (NGC 1068 and NGC 1386), the spatial extension of other kpc scale outflows are unlikely to be affected by it. Based on their [O III]λ5007 flux contours of the red and blue wings, the outflow sizes exhibit a large range from 200 pc to 3.5 kpc. However, the velocity dispersion of the wing component is much smaller than usually seen in the AGN spectra so that their results are not fully comparable with our study. It is possible that the turbulence associated with the small velocity dispersion has a different origin than an AGN-driven outflow.

Fischer et al. (2018) investigated the extent and kinematics of their AGN-ionized gas outflows of 12 nearby (z < 0.12), luminous (L_[O III] > 3 × 10⁴² erg s⁻¹) type 2 QSOs with long-slit spectroscopic and imaging data from HST (STIS, ACS and WFCP2). They reported that, the region where the emission-lines exhibit high velocity offsets from the systematic velocities (> 300 km s⁻¹), are likely affected by the AGN driven outflows. They estimated the maximum projected distance of outflowing kinematics/outflow radius to be 0.15–1.5 kpc (average outflow radius ∼600 pc). They reported that only 2 of QSOs have outflow radius beyond 1 kpc, while 10 other QSOs have radii < 600 pc.

Bellocchi et al. (2019) investigated the ionized gas outflows in 6 type 2 QSOs with observations from the optical imager and long-slit spectrograph (OSIRIS) mounted on the 10.4m Gran Telescopio Canarias Spanish telescope (GTC). They identified ionized outflows in 4 out of 6 type 2 QSOs and reported that the ionized outflows are spatially unresolved in 2 QSOs and compact in 1 QSO. The extent of the compact ionized outflow (R) is, R = 0.8 ± 0.3 kpc. In only 1 QSO (SDSS 0741+3020), Bellocchi et al. (2019) detected a large outflow extending to at least ∼4 kpc from the central engine.

Baron & Netzer (2019) investigated 2377 type 2 AGN showing visible [O III] wing emission from the parent sample of 24 264 AGN following Mullaney et al. (2013). They estimated the mass-weighted outflow or wind location according to the SED fitting method (Netzer 2013) and estimated that the winds are located at 0.1–0.5 kpc away from the AGN.

Recently, Zhang (2021) analyzed a sample of SDSS spectra of 2621 Type-1 AGN and 1987 Type-2 AGN, with reliable [O III] wing component detections. They found different properties of [O III] wing components, which could be confirmed between the Type-1 and Type-2 AGN: statistically lower broad [O III] luminosities and statistically lower flux ratio of the wing to the core [O III] component in the Type-2 AGN, considering necessary effects. Their results indicate stronger obscuration on the [O III] wing components in the Type-2 AGN because [O III] wing emission regions are nearer to the nucleus, under the framework of the Unified model for AGN. They estimated that the distance between broad [O III] emission regions to central BHs [O III] wing emitting region is, R_B3 ∼ R_sub × tan(θ/2) where R_sub is the dust sublimation radius and θ is the opening angle for type 2 AGN. They assumed a global mean value of the opening angle 25° and estimated R_B3 ∼ 20 × R_BLR. If we assume R_BLR = 0.12 pc as per Sturm et al. (2018), the method prescribed by Zhang (2021) will estimate the [O III] wing location to be ∼2 pc away from the nucleus. Additionally, Zhang (2021) stated that the red and blue [O III] wings could be robust signs of central outflows. Consequently, Zhang (2021) concluded that their results provide evidence for obscured central outflows in Type-2 AGN. But the drawback of their methods was that they estimated the location of [O III] wings using SDSS 1D spectra They did not have the IFU data, that we have in the CARS survey. We can accurately estimate, using spectro-astrometric analysis, how far from the nucleus those maximum [O III] wing emitting regions are located, which is impossible with fibre spectroscopy. Therefore, R_B3 might be associated with significant uncertainties.

As per our results, for 23 out of 36 sources we measure a median projected offset of ∼27 pc between the peak of [O III] wing emitting region and the nucleus. Now, the [O III] wing emission we have studied in this work could be from highly excited ionized gas clouds. There might be other [O III] wing emitting gas clouds, whose emissions could be fainter than the ones we have seen so far. Those fainter [O III] wing emitting gas clouds could be extended to kpc scales. The fact that Zhang (2021) found very compact [O III] wing emitting regions, itself indicates that the past work showing AGN driven outflows extended to kpc scale, likely was concerned with the fainter [O III] wing component in Type 2 AGN, which is not obscured by dust.

4.1.2. Effect of AGN inclination angle

As we are only able to measure projected sizes, the classification of compact and extended outflows might be affected by the inclination of the AGN with respect to our line of sight. We consider the effect of inclination on the outflow-size vs. luminosity relation provided by Kang & Woo (2018) for inclination angles of 5°, 20° and 45° as shown in Fig. 8. In these cases, 0° inclination would indicate that an observer will not be facing the obscuring torus structure and will be straight looking into the pole of the AGN (face-on), whereas 90° inclination would indicate that an observer will be facing the obscuring torus structure and will be straight looking into torus (edge on). As the inclination angle decreases, the projected outflow size also decreases. At an inclination of 5°, the size-luminosity relation would be matching with our compact [O III] wing component offsets. The question is whether a random sample of type 1 AGN would contain so many objects at such low inclinations.

As the CARS sample consists of 39 type 1 AGN, we investigate the inclination distribution of a similar number of AGN to see whether inclination could alter the spatial offsets significantly or not. In Fig. 10, we plot the probability distribution function and one possible random sample of 39 type 1 and 39 type 2 AGN assuming an underlying bi-conical geometry. The sample number is equivalent to the number of AGN in the CARS sample. As the bolometric luminosity in our CARS sample is between 6 × 10⁴¹–3 × 10⁴⁵ erg s⁻¹ (Husemann et al. 2017, 2022), we estimate that the half-opening angle is around 60° using the work of Merloni et al. (2014). Hence, anything below 60° would be seen as a type 1 AGN and anything beyond will be viewed as type 2.

Fig. 10. An example of possible inclination distributions of 39 random samples of AGN with inclination bin-sizes of 5°. We see 39 type 1 AGN randomly distributed between inclination angle between 0–60°; and 39 type 2 AGN randomly distributed between inclination angle 60–90°. The mean inclination angle for type 1 AGN ∼(40 ± 2)°, and for type 2 AGN it is ∼(75 ± 2)°. We over-plot the probability distribution functions of the inclination angle for type 1 and type 2 AGN as the dashed lines.

We compute the probabilities (see Appendix. B) for the detailed derivation) of seeing a type 1 and type 2 AGN at any inclination angle (i) and find that, the mean inclination angle for a sample is (40 ± 2)° for a sample of 39 type 1 AGN and (75 ± 2)° for a sample of 39 type 2 AGN. For our type 1 AGN sample, this mean inclination angle will increase the projected size of the [O III] wing component only by a factor of ∼1.5. In Fig. 10, we notice that the number of AGN with low inclinations rapidly decreases. More specifically, we find that the number of type 1 AGN (N_type 1) with inclination angles between 0–10° is at most 3, with 99.5% confidence. This is in contrast to the 23 compact [O III] wing sources in our sample which is much larger than can be predicted by the random type 1 AGN sample of the same sample size. We already get a mean inclination and an error bar from simulating the distributions of 39 type 1 AGN. The mean inclination is model dependent and could be very different in case the simple ionization cone model is not valid. Therefore, having a larger sample will not alter the mean inclination angle drastically as it is model dependent.

Therefore, we conclude that inclination alone cannot explain the large number of targets with compact [O III] wing locations < 100 pc away from the nucleus in our sample. This notion is further supported by other studies that report compact or spatially-unresolved broad [O III] on sub-kpc scales even in type 2 AGN (e.g., Villar-Martin et al. 2016; Tadhunter et al. 2018), which would be even harder to explain by inclination effects.

4.1.3. Relation to the electron density

Several previous studies investigating AGN-driven outflows assumed uniform electron densities in the outflowing medium (e.g., Veilleux et al. 2005; Cano-Díaz et al. 2012; Genzel et al. 2014; McElroy et al. 2015; Husemann et al. 2016; Kakkad et al. 2016; Leung et al. 2017); but the downside is that the electron density could vary spatially and hence assuming a constant electron density would introduce further systematic uncertainty (Kakkad et al. 2018). Bennert et al. (2006) estimated the electron density as a function of radius from the nucleus in a sample of low redshift type 1 and type 2 AGN, and reported that the electron density decreases from ∼10³ cm⁻³ in the central region to ∼10² cm⁻³ at ∼1.2 kpc away from the nucleus. A similar decrease in electron density from the nucleus starting at > 10³ cm⁻³ toward the outskirts was also found in AGN with outflows by Sharp & Bland-Hawthorn (2010), Cresci et al. (2015), and Freitas et al. (2018). Recently, Kakkad et al. (2018), Baron et al. (2020), and Davies et al. (2020) measured the electron density for the broad wing component of [S II]λλ6716, 6731 and reported radially decreasing electron densities with similar amplitude as previous works. Furthermore, other studies involving AGN have clearly demonstrated that the outflowing gas is denser than the ambient, non outflowing gas in luminous AGN (e.g., Holt et al. 2011; Villar Martín et al. 2014; Villar-Martin et al. 2016). Theoretical/numerical studies (e.g., Gaspari et al. 2017) find analogous radial decrease in the electron gas density in simulations of Chaotic Cold Accretion (CCA) feeding and feedback. Thus, this appears to be a natural outcome of galaxy evolution. Overall, we would expect the electron densities of the compact and extended [O III] wing emitting sources to be similar if the broad [O III] wing emission intrinsically originates from kpc scales and if our classification of compact and extended is just an observational bias due to projection effects or sensitivity limitations.

In Fig. 11, we compare the mean electron densities for 7 of our sources with 4 point-like (compact) and 3 extended spatial [O III] wing component distributions. We computed the electron densities within the central 0.6″ (< 500 pc) for those 7 sources which have sufficient S/N to separate the wing and core component even in the faint [S II] line. We find that the median electron density for the sources with a compact [O III] wing regions is n_e = 1915 ± 920 cm⁻³. For sources with extended regions [O III] wing, the median electron density is n_e = 503 ± 307 cm⁻³, which is a factor of 4 lower. A two sample Kolmogorov-Smirnov test returns a p-value = 0.007, which confirms that the two samples of electron densities are significantly different. Such a bimodality in electron density would be unexpected if the emission of the wing component originated from similar spatial scales within the galaxy. This dichotomy in electron densities therefore indicates that the [O III] wing emission for the extended sources is physically much more extended than that of the compact sources. In agreement with the observations and simulations quoted above, the significantly higher electron densities support our measurements that the compact wing regions originate much closer to the nucleus than the extended ones. This observation further strengthens our conclusion that projection effects due to different inclination angles cannot explain the measured spatial distribution of [O III] and that other parameters need to play an important role.

Fig. 11. Bar graphs of electron-densities for the compact and extended outflows. On the left, the bar represents average electron density for the sources with compact outflows; and on the right, the bar represents average electron density for the sources with extended outflows.

4.2. Investigating the AGN luminosity dependence of [O III] wing extent

Assuming that the [O III] wing is clearly related to an AGN-driven wind, previous studies found a clear trend of increasing outflow size with AGN luminosity. Liu et al. (2013a), Kang & Woo (2018) reported a correlation between the [O III] luminosity and the outflow sizes. As per their finding, the ionized outflow size would increase if the AGN had higher [O III] luminosity. Zakamska & Greene (2014), Wylezalek et al. (2020) reported that the AGN luminosity is primarily responsible for the launching and detection of the outflows.

Based on a comparison of the compact versus extended outflows identified for the CARS sources, we find that AGN luminosity is not the main driver of the dichotomity. In Fig. 12, we show the histograms of the BLR Hβ (L_{BLR Hβ}) luminosity for the sources with compact and extended outflows. For the sources with compact [O III] wing emission, L_{BLR Hβ} tends to cover the range of 2 × 10⁴⁰–6 × 10⁴² erg s⁻¹ whereas for the sources with extended outflows the range is wider (6 × 10³⁹–6 × 10⁴² erg s⁻¹). The mean luminosity for the sources with both compact and extended [O III] wing emission is 2.7 × 10⁴¹ ergs⁻¹. This is in agreement with the data from Mullaney et al. (2013) where the highest number of type 1 AGN is located at L_{BLR Hβ} = 3 × 10⁴¹ erg s⁻¹.

Fig. 12. Distribution of the BLR Hβ luminosity for QSOs with compact and extended ionized outflows. The black vertical dashed lines represent the mean BLR Hβ luminosities for the sources with compact [O III] wing emission and the brown vertical dashed lines represent the mean BLR Hβ luminosities for the sources with extended [O III] wing emission.

The sample size in Fig. 12 imposes limits on the minimum difference between the samples that can be diagnosed as statistically significant. For the given sample size, we use a two-sample Anderson-Darling (AD) test (Anderson & Darling 1952), which results in a p-value of 0.99 indicating that the compact and extended outflow samples in Fig. 12 are consistent with being drawn from the same parent distribution.

However, an offset in AGN bolometric luminosity between the two samples may be hidden due to the systematic uncertainties of the bolometric correction from L_{BLR Hβ} to L_bol. Here, we assume a systematic uncertainty of 0.3 dex. This is motivated by the fact that L_{BLR Hβ} has been reported to be tightly correlated to L_5100 Å (Greene & Ho 2005), while Richards et al. (2006) estimated the mean and standard deviation in their B-band bolometric correction factor to be 10.4 and 2.5, resulting in L_bol ≈ (8 − 12) × λL_L5100 Å. We simulate the samples from Gaussian distributions with a standard deviation of 0.3 dex with a range of offsets between the mean values of the two distributions. The AD test shows that at a separation of > 0.4 dex in luminosity we can be sure that the samples are different at a 5% confidence level. Hence, there could be a small difference in luminosity that we cannot diagnose. On the other hand, if we use the outflow size vs [O III] luminosity relation from the work of Kang & Woo (2018), we find that the bolometric luminosity (combining with the work from Heckman & Best 2014) difference between AGN with outflow sizes = 100 pc and 1 kpc will be 2.34 dex. This difference in bolometric luminosity is much larger than what we report here.

4.3. Dependence of the mass outflow rates on measurements methods

In the last decade, several optical IFU studies investigated the warm ionized phase of AGN driven outflows, by tracing the rest-frame optical emission line [O III] (e.g., Liu et al. 2013b; Mullaney et al. 2013; Harrison et al. 2014; Karouzos et al. 2016a,b; Woo et al. 2016). But the outflow sizes (e.g., Husemann et al. 2016; Villar-Martin et al. 2016; Harrison et al. 2018; Tadhunter et al. 2018) and the electron densities (e.g., Sharp & Bland-Hawthorn 2010; Cresci et al. 2015; Freitas et al. 2018; Kakkad et al. 2018; Baron et al. 2020) are still uncertain. This imposes a challenge when deriving the outflow energetics, which strongly depends upon of the outflow kinematics, spatial extent and electron density.

A major advantage of CARS is that, the sample consists of nearby (0.01 < z < 0.06) and luminous (L_[O III] ∼ 9 × 10³⁹ − 3 × 10⁴² erg s⁻¹) AGN, which well represents the general type 1 AGN population at low redshifts. Because of their proximity, the CARS host galaxies can be studied at sub-kpc scales, while still probing an important part of the local AGN luminosity function (Husemann et al. 2017, 2019b, 2022).

IFU data allows us to directly constrain critical parameters for outflow energetics, such as outflow radius, electron density, outflow velocity and extinction. The issue is that the IFU studies only observe a small sample, whereas fibre spectroscopy covers a large sample. On the other hand, fibre spectroscopy does not probe the outflows radius and the electron density at the outflow radius. Therefore, it is crucial to put our results in perspective with the previous spectroscopic studies of outflow energetics to investigate how much the values of the energetics parameters could change. Assuming a case B recombination (Nesvadba et al. 2006, 2011), with an electron temperature, T_e = 10⁴ K (Osterbrock & Ferland 2006), and adopting an intrinsic line ratio of [O III]λ5007 / Hβ = 3.8,

$$\frac{M_{\mathrm{g}\mathrm{a}\mathrm{s}}}{4\times {10}^7{M}_{\mathrm{s}}}=\left(\frac{L_{[\text{O} \text{III}], \mathrm{w}\mathrm{i}\mathrm{n}\mathrm{g}}}{{10}^{43}{\mathrm{e}\mathrm{r}\mathrm{g}\mathrm{s}}^{-1}}\right){\left(\frac{n_{\mathrm{e}, \mathrm{w}\mathrm{i}\mathrm{n}\mathrm{g}}}{100{\mathrm{c}\mathrm{m}}^{-3}}\right)}^{-1}$$(2)

where M_s is the solar mass, L_{[O III],wing} is [O III] wing luminosity and n_{e,  wing} is the electron density of the gas emitting [O III] wing.

Rakshit & Woo (2018) used the SDSS sample to generate a large sample of type 1 AGN. However, estimating the outflow energetics requires the outflow-size within which the energetics is to be calculated. Additionally, outflow energetics values are highly dependent upon the assumption of outflow geometry. Rakshit & Woo (2018) assumed a bi-conical outflow geometry and ‘Case B’ recombination (Crenshaw et al. 2010; Bae & Woo 2016). As their work was confined to a SDSS 1D spectroscopic sample of ∼5000 type 1 AGNs at z < 0.3, they estimated the outflow radius (D) using the outflow-size luminosity relation from Kang & Woo (2018). They approximated the outflow velocity as

$$\begin{array}{c}\hfill v_{\mathrm{out}}={({{\sigma }_{[\text{O}\text{III}]}}^2+{V_{[\text{O}\text{III}]}}^2)}^{1/2}\end{array}$$(3)

where V_[O III] is the velocity shift, defined as the velocity offset of the first moment of the single/double Gaussian [O III]λ5007 line profile from the systematic velocity. σ_[O III] is the second moment of the single/double Gaussian [O III]λ5007 line profile. Rakshit & Woo (2018) reported that 1526 objects in their sample overlap with the sample of Rakshit et al. (2017) with a median n_e = 272 cm⁻³. Following the work of Maiolino et al. (2012), Bae et al. (2017), Rakshit & Woo (2018) estimated the mass outflow rates to be:

$$\begin{array}{c}\hfill {\dot{M}}_{\mathrm{out}}=3M_{\mathrm{out}}\frac{v_{\mathrm{out}}}{D}.\end{array}$$(4)

We investigate how the energetics will change if the outflows are predominantly compact and unrelated to AGN luminosity.

Although, the measurement of the [O III] line shape are very accurate, their interpretation to determine the correct outflow velocity is difficult as different studies use different parameters as outflow velocities. Therefore, we aim to focus on the dependence of mass outflow rates on outflow radius (D) and electron densities n_e. The energy and momentum outflow rates could be obtained by multiplying the outflow velocity term to the mass outflow rates. Therefore, we entirely focus on the mass outflow rates for an initial comparison of assumptions.

In Fig. 13, we plot the mass outflow rates against the BLR Hβ luminosity (L_{Hβ BLR}) for the type 1 AGN of the samples of Mullaney et al. (2013). We only choose AGN with L_Hβ BLR > 4 × 10³⁷ erg s⁻¹. In the left panel, we calculate the outflow radius, D from their respective [O III] luminosities (Kang & Woo 2018) and use n_e = 272 cm⁻³ (Rakshit et al. 2017; Rakshit & Woo 2018). We also use the scikitlearn package (Pedregosa et al. 2011) to perform a linear regression between Ṁ_out and L_{Hβ BLR}. We obtain –

$$\begin{array}{c}\hfill \log ({\dot{M}}_{\mathrm{out}})=0.35\log (L_{\mathrm{H}}\beta BLR)-15.07.\end{array}$$(5)

Fig. 13. BLR Hβ luminosity against the mass outflow rates for the type 1 AGN from Mullaney et al. (2013). Left panel: we use the outflow size vs. [O III] luminosity relation from the work of Kang & Woo (2018) to calculate the outflow radius (D). We also use the median electron density ne = 272 cm−3 as reported by Rakshit & Woo (2018). The green dashed line indicates the linear regression line of the mass outflow rates and the BLR Hβ luminosity. Mid panel: we assume that all of the outflows in this sample of type 1 AGN of Mullaney et al. (2013) to have point like [O III] wing emission; in which case the median offset between the AGN nucleus and the [O III] wing Moffat center is ∼27 pc. Therefore, we use 27 pc as the outflow radius in this case. We keep the electron density fixed to ne = 272 cm−3 as reported by Rakshit & Woo (2018). Right panel: we use ne = 1915 cm−3 which is the median electron density for the AGN with compact [O III] wing emission. We use an outflow radius of 27 pc which is about the median offset between the Moffat center of the [O III] wing and the nucleus. We over-plot the linear regression line of the mass outflow rates and the BLR Hβ luminosity in all panels to showcase the difference in mass outflow rates in the approach used in all the panels.

In the middle panel, we assume that all the AGN in this particular subsample have spatially-unresolved [O III] wing emission. In our CARS sample, the median offset of the Moffat center of the [O III] wing from the nucleus is ∼27 pc. Therefore, we assume the outflow radius, D = 27 pc. We keep the electron density fixed at 272 cm⁻³ to highlight the sole impact of a small outflow size. We report 30 times increase in mass outflow rate on average compared to the left panel.

In the right panel, we assume the outflow radius, D = 27 pc similar to the middle panel. This time we also change the electron density to see the cumulative effect of changing the D and n_e, compared to the left panel. The AGN with compact [O III] wing emission in our CARS sample have an median electron density n_e = 1915 ± 920 cm⁻³. Therefore, we use n_e ∼ 1900 cm⁻³ to see how the mass outflow rate changes. We report an average ∼4 times increase in the mass-outflow rate, which is about 9 times less than the middle panel. Hence, we see a net increase in mass outflow rate compared to the method by Rakshit & Woo (2018) when assuming compact outflows, which should be frequent in this luminosity range according to our studies. While we cannot reconstruct the actual mass outflow rate distribution for the SDSS sample, we predict that many local AGN have more energetic outflows than previously measured if they are also more centrally concentrated.

Outflows from AGN have been observed in various gas phases – from highly ionized X-ray and UV absorption lines (e.g., Blustin et al. 2003; Reeves et al. 2003), ionized optical emission lines (e.g., Mullaney et al. 2013; Harrison et al. 2014; Husemann et al. 2016; Perna et al. 2017), to atomic and molecular emission and absorption lines (e.g., Feruglio et al. 2010; Rupke & Veilleux 2013; Veilleux et al. 2013). Some theoretical studies have also confirmed the presence of the multi-phase nature of these outflows (Zubovas & King 2012; Gaspari & Sądowski 2017). The nature of these multi-phase outflows still remains largely unconstrained, as different gas phases have been found to show different outflow velocities, mass, and covering factors (e.g., Fiore et al. 2017; Cicone et al. 2018; Veilleux et al. 2020). Additionally, the warm ionized gas phase traced by the optical emission lines only represents a small fraction of the mass of the multi-phase outflows. Rupke & Veilleux (2013) observed a higher mass in the neutral gas outflows than in the ionized gas outflows from ULIRGs. Therefore, other gas phases of these multi-phase outflows need to be investigated in order to get a detailed view of their energetics.

5. Summary and conclusions

We used optical IFU observations from CARS, covering the [O III] λλ4959,5007 and Hβ emission lines to investigate the spatial location of the [O III] blue wing component. Our CARS targets are selected from the Hamburg-ESO survey, a sample of luminous QSOs with 12.5 < B < 17.5 (Wisotzki et al. 2000) and representative of the bulk of the local AGN population (10⁴⁰ erg s⁻¹ < L_[O III] < 3 × 10⁴² erg s⁻¹) at z < 0.06. All CARS sources show asymmetries in [O III]]λλ4959, 5007 as a potential sign of ionized outflows. In this paper we use a spectro-astrometric analysis to characterize the [O III] wing emission. Our main results are:

• The [O III] wing emission profiles in 23 out of 36 AGN are consistent with the PSF profiles with peak positions offset from the nucleus by < 100 pc. The other 13 sources display a clearly extended [O III] wing spatial emission profile. Therefore, the spatial morphology of the [O III] wing emission is diverse and can be compact or extended in an unbiased luminous AGN sample, where both are likely to occur.

• The electron density in the compact [O III] wing regions (n_e ∼ 1900 cm⁻³) is nearly a magnitude higher than the one in the extended [O III] wing regions (n_e ∼ 500 cm⁻³).

• The sources with compact and extended [O III] wings follow the same distribution in terms of AGN (BLR Hβ) luminosity. The inclination correction that would be required to interpret the compact outflows as extended kpc size outflows is inconsistent with the distribution of inclinations expected for our sample.

• Inclination angle and AGN bolometric luminosity alone cannot explain the large number of compact [O III] wing emission in CARS sources.

• The mass outflow rates for type 1 AGN calculated over kpc scales for L_{BLR Hβ} = 10³⁸ erg s⁻¹ to 3 × 10⁴² erg s⁻¹ are Ṁ_out = 0.01 − 1.6 M_⊙ yr⁻¹.

• For a given electron density, the mass outflow rates increase up to 30 times if estimated near the nucleus (D = 27 pc) compared to the mass outflow rates calculated over kpc scales.

• The mass outflow rates increase only marginally if estimated near the nucleus (D = 27 pc) assuming with high electron density (n_e ∼ 1900 cm⁻³), compared to the mass outflow rates calculated over kpc scales, assuming low electron density (n_e ∼ 272 cm⁻³).

In this paper, we were able to recover compact sub-kpc scale out to extended kpc scale ionized outflows among a homogeneous sample of 36 luminous AGN observed with seeing-limited MUSE and VIMOS IFU observations. Strikingly, AGN luminosity is not the key parameter for the size of circum-nuclear outflow regions as previously proposed from various seeing-limited IFU studies (Liu et al. 2013b; Kang & Woo 2018; Wylezalek et al. 2020). Therefore, future studies should investigate not only correlations with AGN luminosity, but also with other properties/mechanisms, such as radio jets, that could potentially drive these AGN outflows. Indeed, mechanical (and not radiative) AGN feedback has been shown to be the leading process in shaping the evolution of galaxies over the long term (e.g., Santoro et al. 2018; Husemann et al. 2019b; Gaspari et al. 2020 for reviews). Generally, other factors such as the inclination and duty cycle could play a potential role in setting the apparent size and statistics of these outflows.

Therefore, a lot of things need to be explored in the future with CARS. Clearly, we will focus on measuring the sizes for the sources with extended [O III] wing emission, and to properly obtain their outflow energetics. Additionally, we will accurately estimate the electron density associated with these outflows by using the additional spectra from newly acquired, deep VLT/X-SHOOTER observations where MUSE data are currently not deep enough. Moreover, we will investigate the impact of radio jets on outflows by comparing the radio jet sizes and properties from dedicated VLA observations; and systematically check the presence of any obvious relation between the compact outflows and the compact radio jets in radio-quiet QSOs. Lastly, we will also explore the multi-phase nature of such outflows by combining cold, neutral and ionized gas outflows observations (e.g., Husemann et al. 2019b) for a large fraction of CARS targets.

Acknowledgments

We thank the referee, Dr. Montserrat Villar Martin, for her valuable and thoughtful comments which significantly improved the quality of the manuscript. The work of MS was supported in part by the University of Manitoba Faculty of Science Graduate Fellowship (Cangene Award), and by the University of Manitoba Graduate Enhancement of Tri-Council Stipends (GETS) program. BH is greatful to financial support by the DFG grant GE625/17-1 and DLR grant 50OR1911. BH and ISP are additionally supported by the DAAD travel grant 57509925. TU acknowledges funding by the Competitive Fund of the Leibniz Association through grants SAW-2013-AIP-4, SAW-2015-AIP-2 and SAW-2016-IPHT-2. JS is supported by the international Gemini Observatory, a program of NSF’s NOIRLab, which is managed by the Association of Universities for Research in Astronomy (AURA) under a cooperative agreement with the National Science Foundation, on behalf of the Gemini partnership of Argentina, Brazil, Canada, Chile, the Republic of Korea, and the United States of America. MS, CO and SB acknowledge partial support from the Natural Sciences and Engineering Research Council (NSERC) of Canada. MG acknowledges partial support by NASA Chandra GO8-19104X/GO9-20114X and HST GO-15890.020-A. TR is supported by the Science and Technology Facilities Council (STFC) through grant ST/R504725/1. TAD acknowledges support from the UK Science and Technology Facilities Council through grant ST/S00033X/1. BAT was supported by the Harvard Future Faculty Leaders Postdoctoral Fellowship. Based on observations collected at the European Southern Observatory under ESO programme(s) 083.B-0801, 094.B-0345(A) and 095.B-0015(A). GRT acknowledges support from NASA through grant numbers HST-GO-15411 and HST-GO-15440 from the Space Telescope Science Institute, which is operated by AURA, Inc., under NASA contract NAS 5-26555. MPT acknowledges financial support from the State Agency for Research of the Spanish MCIU through the “Center of Excellence Severo Ochoa” award to the Instituto de Astrofísica de Andalucía (SEV-2017-0709) and through grant PGC2018-098915-B-C21 (MCI/AEI/FEDER, UE). The Science, Technology and Facilities Council is acknowledged by JN for support through the Consolidated Grant Cosmology and Astrophysics at Portsmouth, ST/S000550/1. Based on observations collected at the Centro Astronómico Hispano-Alemán (CAHA) at Calar Alto, operated jointly by Junta de Andalucía and Consejo Superior de Investigaciones Científicas (IAA-CSIC). This research made use of Astropy, a community-developed core Python package for Astronomy (Astropy Collaboration 2013, 2018), SciPy (Virtanen et al. 2020), NumPy (van der Walt et al. 2011); and the plotting packages Matplotlib (Hunter 2007), CMasher (van der Velden 2020).

References

Appendix A: Moffat model formulation

We assumed a 2D Moffat function to model the PSF following the work by Gadotti (2008)

$$\begin{array}{c}\hfill f(r)=a[1+(\frac{r}{\alpha }{)}^2{]}^{-\beta }\end{array}$$(A.1)

where a is the amplitude of the PSF, α denotes the radius of the PSF, and β controls the radial shape of the intensity profile. Because the PSF might not be exactly circular, we used an elliptical Moffat model which is mathematically described as

$$\begin{array}{c}\hfill f(x,y)=a[1+A{(x-x_0)}^2+B{(y-y_0)}^2+C(x-x_0)(y-y_0){]}^{-\beta }\end{array}$$(A.2)

with A, B and C as follows

$$\begin{array}{c}\hfill A={\left(\frac{cos\varphi }{{\alpha }_1}\right)}^2+{\left(\frac{sin\varphi }{{\alpha }_2}\right)}^2\end{array}$$(A.3)

$$\begin{array}{c}\hfill B={\left(\frac{sin\varphi }{{\alpha }_1}\right)}^2+{\left(\frac{cos\varphi }{{\alpha }_2}\right)}^2\end{array}$$(A.4)

$$\begin{array}{c}\hfill C=2sin\varphi cos\varphi [{\left(\frac{1}{{\alpha }_1}\right)}^2-{\left(\frac{1}{{\alpha }_2}\right)}^2]\end{array}$$(A.5)

Where α₁ and α₂ are the minor and major axes and the ϕ is the position angle. The eccentricity (e) of the ellipse could be defined as,

$$\begin{array}{c}\hfill e=\sqrt{1-\frac{{\alpha }_1^2}{{\alpha }_2^2}}\end{array}$$(A.6)

Like an ellipse, there are two axes eventuating two FWHMs which depend on ϕ. In this case, the FWHM changes due to its ellipticity.

$$\begin{array}{c}\hfill {\mathrm{FWHM}}_1=2{\alpha }_1\sqrt{2^{1/\beta }-1}\end{array}$$(A.7)

$$\begin{array}{c}\hfill {\mathrm{FWHM}}_2=2{\alpha }_2\sqrt{2^{1/\beta }-1}\end{array}$$(A.8)

Appendix B: The inclination angle - AGN probability relation

We assume a differential area of a hemisphere which is

$$\begin{array}{c}\hfill dS=r^{2}sin(\theta )d\theta d\varphi \end{array}$$(B.1)

Now, the azimuthal angle ϕ goes from 0 to 2π and θ goes from θ_max to θ_min.

Therefore, integrating the area of subtended by the part of hemisphere is,

$$\begin{array}{c}\hfill S=r^2{\displaystyle {\int }_{{\theta }_{min}}^{{\theta }_{max}}\sin (\theta )d\theta {\int }_0^{2\pi }d\varphi =2\pi r^2[\cos ({\theta }_{min})-\cos ({\theta }_{max})]}\end{array}$$(B.2)

As the area of a hemisphere is 2πr², therefore the probability of an AGN to be within angle θ_min to θ_max is,

$$\begin{array}{c}\hfill P=\frac{2\pi r^2[\cos ({\theta }_{\min })-\cos ({\theta }_{\max })]}{2\pi r^2}=\cos ({\theta }_{\min })-\cos ({\theta }_{\max })\end{array}$$(B.3)

For type 1 AGN, the inclination angle spans the range 0 to 60 degrees with a inclination bin of 5 degrees. Therefore, the total probability of an AGN to be within these inclination angles should be = 1.

If the normalizing constant is A

$$\begin{array}{c}\hfill A\times \cos (0°)-\cos (60°)=1\end{array}$$(B.4)

Which yields, A = 2

Therefore, for type 1 AGN, therefore the probability of an AGN to be within angle θ_min to θ_max is,

$$\begin{array}{c}\hfill P=2[\cos ({\theta }_{\min })-\cos ({\theta }_{\max })]\end{array}$$(B.5)

Similarly, for type 2 AGN, therefore the probability of an AGN to be within angle θ_min to θ_max is,

$$\begin{array}{c}\hfill P=2[\mathit{cos}({\theta }_{\min })-\mathit{cos}({\theta }_{\max })]\end{array}$$(B.6)

Appendix C: Spectral Modelling

In this section, we list all the spectral modeling of the AGN in our analysis. We first present the Hβ + [O III]+ Fe II modeling for the integrated spectra from all the spaxels in the central 3″ region (comparable to SDSS). Next we present the [S II] modeling for 9 AGN where detected the visible broad wing (6 from MUSE and 2 from VIMOS).

Fig. C.1. AGN spectral modeling for all the objects. Figures are identical to the example case in Fig. 1.

Fig. C.1. continued.

Fig. C.1. continued.

Fig. C.1. continued.

Fig. C.1. continued.

Fig. C.1. continued.

Fig. C.2. AGN spectral modeling of the [S II] doublet for all the objects where a deblending into wing and core component is possible. Figures are identical to the example case in Fig. 7.

Fig. C.2. continued.

Appendix D: Flux maps

In this section, we list all the 2D surface brightness maps Hβ + [O III], their 2D Moffat modelling and the normalized residuals. For the MUSE AGN, the 2D maps were performed in the central 3″ region; whereas, for the VIMOS sources, it covers central 6″ region.

Fig. D.1. Results of the spectro-astrometric analysis for all remaining objects. Figures are identical to the example case in Fig. 4.

Fig. D.1. continued.

Fig. D.1. continued.

Fig. D.1. continued.

Fig. D.1. continued.

Fig. D.1. continued.

Fig. D.1. continued.

Fig. D.1. continued.

Fig. D.1. continued.

Fig. D.1. continued.

Fig. D.1. continued.

Fig. D.1. continued.

Fig. D.1. continued.

Fig. D.1. continued.

Fig. D.1. continued.

Fig. D.1. continued.

Fig. D.1. continued.

Fig. D.1. continued.

All Tables

All Figures

Fig. 1. Example of the multicomponent modelling for the spectrum of HE 0345+0056. Upper panel: brightest spaxel spectrum of HE 0345+0056 and the best-fit spectral model. The blue green lines refer to the narrow core and the blue lines refer to the broad wing components of Hβ + [O III]. The gray lines refer to the individual emission lines of the BLR Hβ and Fe II lines. The narrow width of the BLR component representing Hβ is due to the fact that HE 0345+0056 is a narrow line Seyfert 1 (NLS1), where the broad line region tends to show narrower Balmer lines. These narrow BLR components are common in these objects (e.g., Mathur et al. 2000; Mullaney & Ward 2008). Lower panel: residual spectrum normalized by the error spectrum. The dark gray shaded area highlights the ±3σ limit and the red dashed line provides the reference to 0.

In the text

Fig. 2. [O III] line width described by W80 against AGN bolometric [O III] luminosity (left panel) and BLR Hβ luminosity (right panel). The yellow circles denote the type 2 AGN sample from Liu et al. (2013a), the blue diamonds denote the type 2 AGN sample by Harrison et al. (2014), the lime triangles denote the type 2 AGN sample by Sun et al. (2017) and the red stars denote our CARS sample of type 1 AGN. The dense scatter plot corresponds to the overall SDSS sample of both type 1 and type 2 (left panel) and type 1 AGN (right panel) from the ALPAKA library published by Mullaney et al. (2013).

In the text

	Fig. 3. Central 3″ intensity maps associated with different spectral lines. The leftmost map describes the light distribution from Hβ broad emission-line region, which is the PSF. The bottom-right map describes the light distribution from blue-shifted ionized [O III] wing; and the top-right map describes light-distribution from the bright, narrow core of [O III].
In the text

Fig. 4. Results of our spectro-astrometric analysis for HE 934+0119 from the MUSE IFU observations. We have mapped the surface brightness within the central 3″ for all spectral components in the brightest spectrum assuming fixed kinematics. Here we only show the measure light distribution from the BLR Hβ (left), broad [O III] wing (middle), and narrow [O III] core (right) component. Those 2D light distributions are fitted with a 2D Moffat profile. The red cross indicates the Moffat center for the best-fit elliptical Moffat model to the surface brightness map of the BLR Hβ. The green and red crosses indicate the Moffat centers for the best-fit elliptical Moffat models to the surface brightness maps of the [O III] wing and core, in which the amplitude and position of the Moffat model were varied only. The histograms show the distributions of the normalized residuals (residual/error) for both the core and wing components in comparison to the BLR Hβ one.

In the text

Fig. 5. Results from our analysis: Out of 36 CARS targets, 23 agree with a point-like broad [O III] wing emission profile, whereas the other 13 targets have broad [O III] wing emission profile that is not point-like. We also show one source from each of compact (HE 1011−0403) and extended (HE 1248−1356) [O III] wing emitting sources. We further plot the green contours showing the [O III] wing emission-profiles, along with the black contours of the BLR Hβ (PSF) emission profile.

In the text

	Fig. 6. r[O III] against [O III] luminosity for the sources with compact and extended [O III] wing emission. The dashed line represent the region where Faperture = FMoffat.
In the text

Fig. 7. Example of the multicomponent modelling for the spectrum of HE 0253−1641. Upper panel: integrated 3 × 3 spaxel-spectrum of HE 0253−1641 and the best-fit spectral model. The green lines refer to the narrow core and the blue lines refer to the broad wing components of [S II]. Lower panel: residual spectrum normalized by the error spectrum. The dark gray shaded area highlights the ±3σ limit and the red dashed line provides the reference to 0.

In the text

Fig. 8. Spatial extension of the [O III] wing emission vs. [O III] luminosity for our samples along with outflow sizes obtained by other studies in the literature at similar or higher [O III] luminosity range that find outflows of several kpc scales. For the sources with compact [O III] wing emission, we assume the offset sizes obtained by spectro-astrometry to indicate the how far is the peak of the emission from the nucleus. As previous studies considered the [O III] wing to indicate outflows (e.g., Crenshaw et al. 2010; Zhang 2021), we further plot the outflow sizes measured using the approach taken by Karouzos et al. (2016a), Kang & Woo (2018). We connect the outflow-sizes obtained using the method of Karouzos et al. (2016a) and the offset-sizes using spectro-astrometry for the CARS sources with compact [O III] wing emission to quantify the difference in determining spatial extension using these two methods.

In the text

Fig. 9. Spatial extension of the [O III] wing emission vs [O III] luminosity for our samples along with outflow sizes obtained by Villar-Martin et al. (2016), Tadhunter et al. (2018), Fischer et al. (2018) at similar or higher [O III] luminosity range that find compact outflows. The lower and upper limits of the error bars indicate the lower and upper limits of the outflow sizes obtained by Villar-Martin et al. (2016).

In the text

Fig. 10. An example of possible inclination distributions of 39 random samples of AGN with inclination bin-sizes of 5°. We see 39 type 1 AGN randomly distributed between inclination angle between 0–60°; and 39 type 2 AGN randomly distributed between inclination angle 60–90°. The mean inclination angle for type 1 AGN ∼(40 ± 2)°, and for type 2 AGN it is ∼(75 ± 2)°. We over-plot the probability distribution functions of the inclination angle for type 1 and type 2 AGN as the dashed lines.

In the text

	Fig. 11. Bar graphs of electron-densities for the compact and extended outflows. On the left, the bar represents average electron density for the sources with compact outflows; and on the right, the bar represents average electron density for the sources with extended outflows.
In the text

	Fig. 12. Distribution of the BLR Hβ luminosity for QSOs with compact and extended ionized outflows. The black vertical dashed lines represent the mean BLR Hβ luminosities for the sources with compact [O III] wing emission and the brown vertical dashed lines represent the mean BLR Hβ luminosities for the sources with extended [O III] wing emission.
In the text

Fig. 13. BLR Hβ luminosity against the mass outflow rates for the type 1 AGN from Mullaney et al. (2013). Left panel: we use the outflow size vs. [O III] luminosity relation from the work of Kang & Woo (2018) to calculate the outflow radius (D). We also use the median electron density ne = 272 cm−3 as reported by Rakshit & Woo (2018). The green dashed line indicates the linear regression line of the mass outflow rates and the BLR Hβ luminosity. Mid panel: we assume that all of the outflows in this sample of type 1 AGN of Mullaney et al. (2013) to have point like [O III] wing emission; in which case the median offset between the AGN nucleus and the [O III] wing Moffat center is ∼27 pc. Therefore, we use 27 pc as the outflow radius in this case. We keep the electron density fixed to ne = 272 cm−3 as reported by Rakshit & Woo (2018). Right panel: we use ne = 1915 cm−3 which is the median electron density for the AGN with compact [O III] wing emission. We use an outflow radius of 27 pc which is about the median offset between the Moffat center of the [O III] wing and the nucleus. We over-plot the linear regression line of the mass outflow rates and the BLR Hβ luminosity in all panels to showcase the difference in mass outflow rates in the approach used in all the panels.

In the text

	Fig. C.1. AGN spectral modeling for all the objects. Figures are identical to the example case in Fig. 1.
In the text

	Fig. C.1. continued.
In the text

	Fig. C.1. continued.
In the text

	Fig. C.1. continued.
In the text

	Fig. C.1. continued.
In the text

	Fig. C.1. continued.
In the text

	Fig. C.2. AGN spectral modeling of the [S II] doublet for all the objects where a deblending into wing and core component is possible. Figures are identical to the example case in Fig. 7.
In the text

	Fig. C.2. continued.
In the text

	Fig. D.1. Results of the spectro-astrometric analysis for all remaining objects. Figures are identical to the example case in Fig. 4.
In the text

	Fig. D.1. continued.
In the text

	Fig. D.1. continued.
In the text

	Fig. D.1. continued.
In the text

	Fig. D.1. continued.
In the text

	Fig. D.1. continued.
In the text

	Fig. D.1. continued.
In the text

	Fig. D.1. continued.
In the text

	Fig. D.1. continued.
In the text

	Fig. D.1. continued.
In the text

	Fig. D.1. continued.
In the text

	Fig. D.1. continued.
In the text

	Fig. D.1. continued.
In the text

	Fig. D.1. continued.
In the text

	Fig. D.1. continued.
In the text

	Fig. D.1. continued.
In the text

	Fig. D.1. continued.
In the text

	Fig. D.1. continued.
In the text