© The Authors 2026

1. Introduction

Observed correlations between host galaxy properties and central black hole activity (Peterson 2008; Kormendy & Ho 2013) suggest the existence of a mechanism linking nuclear-scale black hole behavior to galaxy-scale effects, the so-called active galactic nucleus (AGN) feedback, which is a key component in galaxy evolution models (Di Matteo et al. 2005; Hopkins & Elvis 2010), yet its exact workings and triggers remain unclear. Feedback from AGN is generally classified into two main modes: a radiative (or quasar) mode, in which radiation pressure from the accretion disk drives powerful multiphase winds, and a kinetic (or jet) mode, in which mechanical energy from relativistic jets interacts with the surrounding medium. These processes give rise to AGN-driven winds and jets that launch multiphase (i.e., highly ionized, warm ionized, neutral, and molecular) outflows, constituting a primary feedback mechanism supported by multiwavelength observations (e.g., Cicone et al. 2018; Fluetsch et al. 2019; Esposito et al. 2024). As in Zubovas & King (2012), we use “wind” to refer to the mildly relativistic (v ∼ 0.1c) ejection of accretion disk gas from the vicinity of the supermassive black hole resulting from Eddington accretion, and “outflow” for the large-scale nonrelativistic flows produced by the interaction between the wind and the galaxy’s ambient gas. These outflows may provide the connection between the black hole and its host galaxy required to reconcile theory with observations, carrying mass and energy out to larger (galactic) scales (Silk & Rees 1998; Hopkins et al. 2016 and Harrison et al. 2018, and references therein).

Some models (e.g., Faucher-Giguère & Quataert 2012; King & Pounds 2015) propose that feedback begins with a sub-relativistic wind launched near the accretion disk at velocities exceeding several 10³–10⁴ km s⁻¹ (v ∼ 0.01c − 0.1c), detected as highly ionized and highly blueshifted Fe K-shell absorption lines observed in X-rays as ultrafast outflows (UFOs; Cappi 2006; Tombesi et al. 2010, 2012). As these winds interact with the interstellar medium (ISM), they produce additional outflowing material, with lower ionization, observed in the optical (e.g., Harrison et al. 2014; Carniani et al. 2015; Marasco et al. 2020; Esposito et al. 2024) and molecular bands (Veilleux et al. 2013; Feruglio et al. 2015; Cicone et al. 2018; Lutz et al. 2020).

Zubovas & King (2012) presented a model that explains how the different outflow phases are connected and how efficiently the wind energy is transmitted to large-scale outflows. In their model, the transfer of wind energy strongly depends on how the wind interacts with the diffuse ISM of the host galaxy. The wind is hypersonic and, as it propagates through the ISM, it produces reverse and forward shocks separated by a contact discontinuity, which defines a shocked wind and a shocked ambient medium regions on the left and right sides of the contact discontinuity (see Fig. 1 in Zubovas & King 2012 and Fig. 7 in King & Pounds 2015). The shocked wind gas then acts as a piston, sweeping up the host ISM at a contact discontinuity that moves ahead of it. In brief, two possible outcomes arise depending on the cooling time of the shocked wind: (i) if the shocked wind cools on timescales short compared to the motion of the shock patterns, the shocked wind gas is compressed to high densities and radiates away almost all its kinetic energy. This corresponds to the “momentum-driven” regime, in which the outflow cannot reach large distances due to energy losses; (ii) if the shocked wind gas does not cool efficiently and instead expands as a hot bubble, the flow remains adiabatic and the wind retains its energy rate due to a momentum boost, which depends on the UFO velocity and the velocity of the gas phase considered (i.e., molecular or ionized). This corresponds to the “energy-conserving” regime, in which the sizes of the shocked wind and shocked ambient medium are much larger than those in the momentum-driven case and the energy is transferred up to larger distances. This latter model can explain the galaxy-wide outflows found, for instance, in quasar objects (e.g., Feruglio et al. 2010; Cicone et al. 2012; Feruglio et al. 2015).

Following King (2005) and Zubovas & King (2012, 2014), a multiphase outflow can be produced as the result of thermal instabilities in the outflow that would lead to a two-phase medium and the cool gas would become molecular and form stars. In the works by Richings & Faucher-Giguère (2018a,b), they showed, for the first time in their simulations, that the swept up gas from the ambient medium by the outflow is able to cool and form molecules within ∼1 Myr, and thus produce large molecular outflow rates (up to 140 M_⊙ yr⁻¹), supporting the in situ molecule formation after cooling and condensation of the hotter gas.

Thus, in principle, comparing the energetics of X-ray and molecular and/or ionized outflows can provide a test for the energy conservation and the strong coupling between the different phases of the wind (see broader discussion in Harrison & Ramos Almeida 2024). Currently, theoretical models and simulations do not include a standardized recipe for properly measuring AGN feedback (e.g., Nelson et al. 2019; Costa et al. 2020; Ward et al. 2024). Observational constraints generally come from a broad range of data covering various AGN types, wavelengths, and spatial and spectral resolutions (e.g., Harrison et al. 2018; Husemann et al. 2019; Wylezalek et al. 2020).

Spatially resolved observations, obtained via millimeter interferometry, optical or near-infrared integral field units (IFUs), are proving mostly adequate to determine the spatial scales and geometry of outflows (e.g., Cicone et al. 2014; Husemann et al. 2019; Mingozzi et al. 2019; Davies et al. 2020; Kakkad et al. 2020; Venturi et al. 2023; Harrison & Ramos Almeida 2024), at least in local AGNs. Among key findings, surveys such as the Galactic Activity, Torus, and Outflow Survey¹ (GATOS) that focus on hard X-ray selected Seyfert galaxies have revealed molecular outflows in Atacama Large Millimeter Array (ALMA) data for sources exhibiting significant nuclear cold gas deficiency (e.g., Alonso-Herrero et al. 2021; García-Burillo et al. 2021, 2024; Esposito et al. 2024). These results support the idea that outflows efficiently expel gas from the central regions. Regarding ionized gas, the MAGNUM survey (e.g., Mingozzi et al. 2019; Venturi et al. 2021), conducted with the Multi Unit Spectroscopic Explorer (MUSE) at the Very Large Telescope (VLT) on local Seyfert galaxies hosting low-power radio jets, has shown that besides a clear ionized outflow component aligned with ionization cones and radio jets along the galaxy disk, there exists an additional high velocity, kiloparsec-scale region of excited gas extending perpendicular to this direction. This feature likely results from the jet perturbing the gas in the disk. Recent James Webb Space Telescope (JWST) IFU observations have now uncovered similar phenomena in both high-redshift and local AGNs. For instance, in XID2028 (z = 1.59), JWST/NIRSpec revealed an expanding emission-line bubble likely driven by the interaction between the ISM and the radiation-driven outflow and low-power radio jet (Cresci et al. 2023). In the nearby Seyfert 1.5 galaxy NGC 7469, Armus et al. (2023), using the Medium Resolution Spectrometer (MRS) of the Mid-Infrared Instrument (MIRI) on board the JWST data, detected blueshifted high-ionization ([Mg V]) lines reaching velocities of ∼1700 km s⁻¹ co-spatial with the starburst ring. They found that the width of the broad emission and the broad-to-narrow line flux ratios correlate with ionization potential, suggesting the presence of a decelerating, stratified, AGN-driven outflow emerging from the nucleus (see also Robleto-Orús et al. 2021, who characterized the warm ionized outflows in NGC 7469 with MUSE).

Conversely, in recent years there have also been studies of AGN samples with a well-established outflow component in a single band. For instance, Salomé et al. (2023) reported on the molecular gas properties and higher star formation efficiency in a collection of ten narrow line Seyfert 1 (NLSy1) galaxies selected to have a solid detection of X-ray UFO. Among their many peculiarities, NLSy1 galaxies are known for their relatively small Hβ full width at half maximum (FWHM; i.e., < 2000 km s⁻¹; Goodrich 1989) compared to other Seyfert 1 (Sy1) galaxies, which suggests low black hole masses (M_BH ∼ 10⁶–10⁸ M_⊙; e.g., Peterson et al. 2000). They also typically exhibit strong FeII and weak [OIII] emission (Osterbrock & Pogge 1985), as well as pronounced X-ray variability and steep X-ray power-law continua (e.g., Panessa et al. 2011). Furthermore, these sources are characterized by significantly higher accretion rates compared to standard Seyfert galaxies (between L/L_Edd ∼ 0.1–1; Boroson & Green 1992), indicating an extreme central environment, which may suggest that NLSy1 galaxies are in the early phase of their evolution (i.e., Mathur et al. 2001; Véron-Cetty et al. 2001). X-ray UFOs with complex patterns have become fairly common in NLSy1 objects (Xu et al. 2024), a fact interpreted as a natural consequence of their high accretion rate, which enables the launching of disk winds via a radiatively driven mechanism (Matzeu et al. 2023).

There are also uncountable observational studies of individual AGNs showing multiband outflows that are greatly contributing to advance our understanding of the outflow phenomenology and the interplay of the different outflowing gas phases (e.g., Feruglio et al. 2015; Morganti et al. 2015; Fluetsch et al. 2021; Longinotti et al. 2023; Esposito et al. 2024). The study presented herein falls somewhat into this category. The radio-loud NLSy1 galaxy IRAS17020+4544² (IRAS17 hereafter) was the first AGN where a multicomponent UFO was identified using XMM-Newton observations (Longinotti et al. 2015), along with a complex pattern of slower X-ray velocity absorbers that were interpreted as the signature of the UFO shocking the ISM (Sanfrutos et al. 2018). Multiband follow-up studies have uncovered the presence of outflows in every single bandwidth: very long-baseline interferometry (VLBI) studies (Giroletti et al. 2017) revealed elongated synchrotron emission on ∼10 pc scales, possibly indicating that this low-power radio jet contributes to shock-driven interactions with the ISM. Longinotti et al. (2018) reported a powerful energy-conserving molecular outflow in single-dish observations obtained by the Large Millimeter Telescope (LMT). The UV study carried out with the Cosmic Origins Spectrograph (COS) on board the Hubble Space Telescope (HST) data by Mehdipour et al. (2023) identified Lyα absorption at 23 400 km s⁻¹ as the counterpart of the X-ray UFO, a rare finding interpreted as gas entrainment in the shocked wind. Interferometry observations with the NOrthern Extended Millimeter Array (NOEMA) provided the first spatially resolved view of the host galaxy of this interesting AGN: based on these data, Salomé et al. (2021) revealed an interaction with a satellite galaxy located up to ∼8 kpc north of the nucleus, and Longinotti et al. (2023) reported a biconical outflow structure in CO gas extending out to a radius of ≲3 kpc, with velocities ranging from 1000 to 1800 km s⁻¹, and fully consistent with energy-conserving regime, confirming the result obtained by the LMT.

Within the framework of the MATRIOSKA³ project, we aim to characterize the ionized gas phase of the outflowing component in IRAS17 using the Multi-Espectrógrafo en GTC de Alta Resolución para Astronomía (MEGARA; Gil de Paz et al. 2018; Carrasco et al. 2018) mounted on the Gran Telescopio Canarias (GTC). Our observations enabled a detailed spatially resolved (2D) kinematic study of the warm ionized gas, traced by the Hα-[NII]λλ6548,6583 and Hβ-[OIII]λλ4959,5007 complexes. This study has allowed us to assess the presence of the ionized outflowing gas, its connection with other phases (i.e., whether it still follows the previously reported energy-conserving regime for the molecular component), and its overall role in AGN feedback.

This paper is organized as follows. In Sect. 2 we present the observations and data reduction. In Sect. 3 we analyze the ionized gas kinematics. In Sect. 4 we present the kinematic results of the ionized gas phase traced by different emission lines (i.e., Hα and [OIII]), deriving key physical parameters of the outflow. In Sect. 6 we discuss our results, attempt to connect the different outflow phases, and assess the global impact of the outflow on the host galaxy. In Sect. 7 we summarize our main findings and present our conclusions. In Appendix A we provide details on the derivation of the visual extinction. In Appendix B we present additional kinematic maps derived for other emission lines observed in the different setups. Finally, in Appendix C we report new results on the ionized-gas emission from the companion galaxy, which is thought to be interacting with IRAS17.

Throughout this paper, we assume a Λ cold dark matter cosmology with H₀ = 70 km s⁻¹ Mpc⁻¹, Ω_M = 0.3, and Ω_Λ = 0.7. Using the E. L. Wright Cosmology calculator, which is based on the prescription given by Wright (2006), we used a a luminosity distance, D_L, of 274.3 Mpc, and an angular diameter distance, D_A = 243.5 Mpc, as derived in Salomé et al. (2021). This corresponds to a scale of 1.181 kpc/″ at a redshift of z = 0.0612.

2. Observations and data reduction

2.1. MEGARA/GTC observations

The data were gathered as part of two different programs and at different epochs (GTC8-20AMEX and GTC1-24AMEX; PI: Longinotti A.; see Table 1) using the MEGARA⁴ instrument. MEGARA is a high-resolution (HR) optical spectrograph mounted on the 10.4 m GTC telescope with two observing modes: (i) the IFU, or large compact bundle, mode, which covers a field of view (FoV) of 12.5 × 11.3 arcsec² (i.e., 13.6 × 14.5 kpc² at the distance of our object) with a spatial resolution of 0.62″, and (ii) the multi-object spectroscopy (MOS) mode, which allows to observed up to 100 objects (8 for the sky and 92 for the observations) in a region of 3.5 arcmin × 3.5 arcmin around the IFU bundle. Both the IFU and MOS capabilities of MEGARA provide intermediate-to-high spectral resolutions: R = λ/FWHM ∼ 6000, 12 000, and 20 000 for the low-resolution (LR), medium-resolution (MR), and high-resolution HR modes, respectively. In this work we used the MEGARA-IFU, which uses 567 fibers in a hexagonal tessellation. Additional 56 fibers are also used to simultaneously measure sky-background spectra with the observation of the scientific target to subtract them from the observations afterward. MEGARA has 18 different volume-phase holographic (VPH) transmission gratings that allow one to cover a range from 3650 to 9700 Å with the three modes.

For our study, we obtained both LR (LR-V and LR-R) and MR (MR-G) data, which allowed us to cover the spectral ranges 5140–7300 Å and 4960–5445 Å, respectively. In particular, the VPH gratings used have the following characteristics: (i) LR-V: wavelength range of 5140–6170 Å with an average resolution R = 6080 (i.e., 0.27 Å pix⁻¹); (ii) LR-R: wavelength range of 6100–7300 Å with an average resolution R = 6100 (i.e., 0.32 Å pix⁻¹); (iii) MR-G: wavelength range of 4970–5445 Å with an average resolution R = 12 035 (i.e., 0.13 Å pix⁻¹). These gratings cover relevant spectral features, such as the Hα-[NII] complex, the [OI]λλ6300,6363, and [SII]λλ6717,6731 doublets in the LR-R setup (see the nuclear spectrum in Fig. 1), while Hβ-[OIII]λλ4959,5007 complex in the LR-V and MR-G setups. Due to the lower signal-to-noise ratio (S/N) obtained in the previous MR-G observations, we recently requested additional MEGARA IFU data in the LR-V setup to spatially resolve and better characterize the ionized gas using the [OIII]λλ4959,5007 doublet emission. Indeed, although the LR-V data have lower spectral resolution than the MR-G setup and the Hβ emission line falls at the edge of the blue part of the spectrum (making it unusable), these data provide higher S/N per spaxel for the [OIII] doublet, enabling the detection of fainter components.

Fig. 1. Observed integrated spectrum of the inner region (R ∼ 1.9 kpc ≈ 1.6″) obtained with the LR-R setup. Several emission lines are detected, including the prominent Hα-[NII] complex, the fainter [SII] and [OI] doublets, as well as very faint [FeVII], FeII, and HeI lines.

2.2. Data reduction

MEGARA data were reduced following the MEGARA Data Reduction Pipeline (DRP; v0.12.0; Pascual et al. 2022), entirely written in a Python environment. The data reduction was performed according to the instructions provided in the MEGARA cookbook (i.e., Castillo-Morales et al. 2020⁵). We briefly describe the main steps of the data reduction process, which include: (1) the sky and bias subtraction; (2) the flat field correction; (3) the spectra tracing and extraction; (4) correction of fiber and pixel transmission; (5) wavelength calibration; and (6) flux calibration (using the spectrophotometric standard star observed each night). For an in-depth description of the data reduction, we refer the reader to the work by Chamorro-Cazorla et al. (2023).

We finally transformed the final raw-stacked spectra⁶ file to a standard integral-field spectroscopic (IFS) data cube by applying a regularization grid to obtain a resample spaxel size of 0.4″ (corresponding to a physical size of 0.472 kpc). The point spread function (PSF) of the MEGARA data cube can be estimated from the FWHM of the 2D profile brightness distribution of the standard star used for flux calibration (HR8634 and HR4963), which gives a FWHM of 1.2″ and 0.9″, for the grisms LR-R and MR-G, and LR-V, respectively.

3. Data analysis

3.1. Optical continuum maps

IRAS17 has been extensively studied at millimeter (Longinotti et al. 2018, 2023) and X-ray wavelengths (Longinotti et al. 2015; Sanfrutos et al. 2018), but it has received comparatively little attention in the optical regime (e.g., de Grijp et al. 1992). In Fig. 2 we show the continuum image obtained with the the Sloan Digital Sky Survey (SDSS) instrument in the r band, which was used as a reference for the astrometric alignment, since its spatial resolution is similar to that of the MEGARA data, allowing for an accurate registration of the features.

Fig. 2. Top: SDSS-r continuum image of IRAS17. The MEGARA-IFU FoV is indicated by a black rectangle. A foreground star (marked with a black star) is visible to the north along with the companion galaxy (indicated by a magenta square) detected by Salomé et al. (2021). The flux intensity map is shown in nanomaggy units (1 nanomaggy = 3.631 μJy). Bottom: Continuum image obtained with MEGARA-IFU within the wavelength range 6770–6840 Å. The flux intensity map is shown in janskys. North is up and east is to the left in both images. The continuum peak emission is marked with a black cross in both panels and corresponds to the AGN position, as supported by its coincidence with the optical AGN position measured by Gaia, the radio emission detected with e-MERLIN at 1.5 GHz (Longinotti et al. 2023), and the CO peak emission (Salomé et al. 2021).

The optical continuum image derived from the MEGARA-IFU data cube was extracted over the observed wavelength range 6770–6840 Å (Fig. 2, bottom panel). In the SDSS image, we can clearly identify the nuclear emission of IRAS17, its spiral arms, and the presence of a foreground star in the northern part of the FoV. Furthermore, a faint emission feature located just below the star, together forming a “comma-shaped” structure, is visible in the SDSS image, although it is only marginally detected in the MEGARA continuum map, where the contours appear elongated toward the southwest. We associate this feature with the newly identified companion galaxy reported by Salomé et al. (2021), which resembles a tidal-tail-like structure that appears to connect the companion with IRAS17, suggesting a past interaction between them.

Throughout this paper, we consider the peak of the continuum emission as the photometric center (AGN position) of the galaxy. This choice is supported by the spatial coincidence between the optical continuum peak and the optical AGN position measured by Gaia, the radio emission detected with e-MERLIN at 1.5 GHz (Longinotti et al. 2023), and the CO peak emission (Salomé et al. 2021). All these measurements are consistent within their respective uncertainties.

3.2. Multicomponent fit and map generation

In our analysis, we did not subtract the underlying stellar emission of the host galaxy in the MEGARA spectra, as the absence of typical stellar feature in the spectra (and the limited wavelength coverage in the MR-G setup) prevent effective modeling of the stellar continuum. As in Sy1 galaxies, we have minimal starlight contamination, since the continuum is largely dominated by the AGN. However, since the main goal of this work is to characterize the kinematics of the ionized outflow in the galaxy, this limitation does not impact our key results. As discussed in Bellocchi et al. (2019a), while continuum subtraction is important for accurately recovering line fluxes as well as the derivation of the equivalent widths, it is not crucial for kinematic measurements. Therefore, we directly fit the observed emission-line profiles to derive the kinematic properties.

The emission features in each spectral band, such as the Hα-[NII] complex in the LR-R grism and the Hβ-[OIII] complex in the MR-G and LR-V grisms, were analyzed using Gaussian profile fitting (see the example in Fig. 3). This was performed using the Interactive Data Language (IDL) routine MPFITEXPR (Markwardt 2009), which employs a Levenberg-Marquardt least-squares minimization to optimize the fit parameters in each spatial element (spaxel). To properly reproduce the observed line profiles in both complexes, we employed a multicomponent Gaussian fitting approach (e.g., Bellocchi et al. 2013, 2019b). This method allowed us to distinguish between a primary (narrow) component and a secondary (broader) one based on their velocity dispersions. In our case, the secondary emission component can even dominate in flux over the systemic one (i.e., Hα emission). In some spaxels, located in the innermost region, a tertiary component was required. This tertiary component is the broadest one, and it clearly traces a more extreme (i.e., fast) ionized outflow in both the Hα and [OIII] lines (see Sect. 4 for details). Thus, in each spaxel, when multiple Gaussian components are fit, the different components, i.e., primary (P), secondary (S) and tertiary (T), are labeled based on their velocity dispersion, from the narrowest to the broadest, with the narrowest component generally tracing the systemic emission. As a first step, we performed initial fits using a single Gaussian profile per line (hereafter referred to as the one-component, or 1c, model). However, in the central regions (R < 2 kpc), where line profiles are more complex, multicomponent fits (two or three Gaussians per line, 2c or 3c) were clearly required, first based on visual inspection, to adequately reproduce the observed features. We only considered spaxels with S/N > 4, a threshold that holds for each individual component. The S/N was computed by dividing the integrated flux in the velocity channels by the corresponding noise, estimated from the standard deviation of the data-model residuals of the fit around each line (within a range about 60 Å wide).

Fig. 3. Nuclear integrated spectra (within the PSF region) of the Hα-[NII]λλ6548,6583 (top) and Hβ-[OIII]λλ4959,5007 complexes (bottom), fit with three kinematic components. For the Balmer lines, the green, magenta, and blue curves represent the systemic, ILR, and broad (outflowing) components, respectively. For the [OIII] lines, the green, purple, and blue curves represent the systemic, intermediate (slow outflow), and broad (fast outflow) components.

To avoid overfitting, we adopted the ε criterion introduced by Cazzoli et al. (2022), which compares the standard deviation of a line-free continuum region (ε_cont) near the line of interest with the standard deviation of the residuals under the fit line profile (ε_line). For the Hα-[NII] complex, the residuals were computed across all three lines to properly account for blending in spaxels with broad profiles. In the case of the [OIII] doublet, we excluded the Hβ line and calculated the residuals using only the two [OIII] lines. We applied a threshold of ε_line/ε_cont ≥ 1.5 to justify the inclusion of additional components. This criterion, supported by visual inspection of asymmetries in the line profiles, motivated the use of a second component and, in some inner spaxels, a tertiary one, especially when residuals from the 2c fit remained significant (i.e., still applying the condition ε_line/ε_cont ≥ 1.5).

The line widths, σ, of each line is constrained to be broader than the instrumental resolution (σ_INS) and the observed velocity dispersion (σ_obs) was corrected for instrumental broadening to derive the intrinsic velocity dispersion, σ_line, by subtracting σ_INS in quadrature, following the relation ${\sigma }_{\mathrm{line}}=\sqrt{{\sigma }_{\mathrm{obs}}^2-{\sigma }_{\mathrm{INS}}^2}$ (with σ_INS ∼ 22.2 km s⁻¹ for the LR setups and ∼11.2 km s⁻¹ for the MR setup).

For each emission line, we derived flux intensity, velocity field⁷, and intrinsic velocity dispersion maps (see Figs. 7, 8, 9, 10, and 11 and the figures in Appendix B). A comprehensive analysis of the kinematic maps as well as their physical interpretation, is presented in the next sections.

3.2.1. LR-R setup for the Hα-[NII] complex

For the LR-R setup, we generally fit the Hα-[NII] complex with two kinematic components (i.e., primary and secondary), tying the velocity and velocity dispersion of Hα and both [NII] lines to minimize degeneracies and mitigate the effect of telluric absorption near [NII]6548⁸. For the [NII] doublets, we fixed the [NII]λλ6583/6548 flux ratios to ∼3 (e.g., Osterbrock & Ferland 2006). In the innermost region of the FoV, however, an additional intermediate-width Hα component was required to account for the emission from the intermediate-broad line region (ILR), as expected in a NLSy1 galaxy (see Fig. 3, top panel). In this region, we associated the ILR with the secondary component, whereas the broadest and most blueshifted component was identified as the ionized outflow (see Table 2).

3.2.2. MR and LR setups for the Hβ-[OIII] complex

For the MR-G and LR-V setups, we first subtracted the FeII emission template (see Fig. 4) by fitting a template for the multiple iron lines, adding an AGN power-law component (f_λ ∝ λ^β). We used the new template constructed by Park et al. (2022), based on a spectrum of the NLSy1 galaxy Mrk 493 obtained with the HST (e.g., Greene & Ho 2007). Compared to the canonical FeII template of I Zw 1 constructed by Véron-Cetty et al. (2004, see also Torres-Papaqui et al. 2020), Mrk 493 exhibits narrower broad-line widths, lower reddening, and a less extreme Eddington ratio, λ_Edd (i.e., λ_Edd ∼ 2.5 for I Zw 1 versus ∼0.5 for Mrk 493), more similar to that of IRAS17 (λ_Edd ∼ 0.7). The fitting method consists of scaling and velocity-broadening the FeII template to match the original MEGARA spectrum (see also Greene & Ho 2005 for details). Figure 4 shows an example of the FeII template fitting in the MR-G setup.

Fig. 4. FeII template subtraction from the observed Hβ-[OIII] complex at the position of the emission peak in the MR-G setup (within a single spaxel of 0.4″ × 0.4″). The observed spectrum is shown in light blue, the FeII template in orange, and the residuals (pure Hβ-[OIII] complex) in green. The vertical dashed red lines mark the iron emission lines in this wavelength range (i.e., rest-frame wavelengths: λFe II = 4923.92, 4993.35, and 5018.45 Å).

After subtracting the FeII emission template, we fit the Hβ-[OIII] complex independently in the MR-G setup⁹, motivated by the fact that these lines may trace gas with different kinematics (see also Cazzoli et al. 2022; Peralta de Arriba et al. 2023; Esposito et al. 2024; Hermosa Muñoz et al. 2024). In this approach, the kinematics of the Hβ components were allowed to vary independently from those of the [OIII] doublet, as tying their kinematics was previously tested and found to provide not very good fits to the three lines, because Hβ may trace partially different gas phases not present in [OIII], such as the ILR (see Fig. 3).

In particular, we modeled the [OIII]λλ4959,5007 doublet with three distinct components: a narrow primary component tracing the galactic rotation (systemic component), and intermediate and broad components. The kinematics of these lines were tied, and their flux ratio was fixed according to atomic physics (i.e., [OIII]λλ5007/4959∼3; Osterbrock & Ferland 2006).

Hβ, observed in the MR-G setup¹⁰, was also fit with three components following the Hα kinematic structure (i.e., a narrow primary component, an intermediate-width component, and a broad component associated with the outflow). Upper limits on the Hβ component widths were set to the average Hα widths of the corresponding components to ensure physically consistent broad-line properties (see Table 2).

4. Main results: Gas kinematics

As has recently been discussed by Salomé et al. (2021) and Longinotti et al. (2023), two sources can be identified within the MEGARA FoV alongside our target, IRAS17: (i) a foreground star to the north, clearly visible in all optical continuum images (Fig. 2); and (ii) a companion galaxy, located slightly to the northwest of the star, first detected via CO molecular gas emission by Salomé et al. (2021), which is not immediately apparent in the optical image presented here (Fig. 2). Our MEGARA/GTC data further reveal that the companion galaxy also exhibits Hα emission in the same region where the CO(1-0) emission was initially detected (Fig. 5; see the next section for details).

Fig. 5. Hα flux map of the primary component (red contours), with overlaid 12CO(1-0) contours of the systemic emission (white) from Salomé et al. (2021). The molecular and ionized gas in IRAS17 and its companion galaxy, traced on the northwest side by the CO contours, exhibit a remarkable spatial overlap.

Fig. 6. Hα velocity dispersion map of the primary component (white contours) with overlaid CO(1-0) contours of the molecular outflow (blue and red contours trace the approaching and receding molecular outflowing gas, respectively) from Longinotti et al. (2023).

4.1. Hα-[NII] complex

Within the LR-R setup, we detected the bright Hα and [NII] emission-line complex. The Hα emission was fit using three components in the inner region (see the spectra in Fig. 3, top), each corresponding to distinct kinematic features (see Fig. 7). At this wavelength, we identified the following kinematic components for IRAS17:

1. A primary component, corresponding to the systemic component of the galaxy, traces its rotation and shows a peculiar velocity dispersion map with a pronounced enhancement over a large area (i.e., a “butterfly” or “X-shaped” region), surrounded by lower-dispersion regions;

2. A secondary component, characterized by a complex velocity field. Within this component, we can identify an inner region with velocities close to the systemic value, possibly associated with the ILR (according to its velocity dispersion value; see Sect. 4.1.2), and a redshifted region toward the north, apparently connected to the companion galaxy. This secondary component also shows an enhanced velocity dispersion with a butterfly pattern (as for the systemic component) found around the ILR, while the companion galaxy shows a relatively low velocity dispersion;

3. A tertiary component, interpreted as the ionized outflow in Hα, based on its significant velocity shift relative to the systemic component (Δv = v_T − v_P ∼ −250 km s⁻¹) and very broad line widths (σ ∼ 900 km s⁻¹).

As was mentioned in the previous section, the kinematics of the [NII] doublet¹¹ were constrained to match those of the Hα line. The main difference between the [NII] and Hα panels (Figs. 7 and 8) is that the inner region of the secondary component, associated with the ILR, is absent in the forbidden-line tracer [NII].

Fig. 7. Kinematic maps for the ionized component as traced by the Hα line. From left to right, top to bottom: Flux intensity, velocity field and (intrinsic) velocity dispersion maps of the primary (systemic), secondary and tertiary components. The velocity maps of each component have been corrected for the systemic velocity, vsys, derived at the AGN position (i.e., at the continuum intensity peak, marked by a black cross in all the maps). Gray contours in the flux intensity systemic maps trace the R-band continuum emission, highlighting the position of the foreground star. The magenta circles in the primary and secondary maps mark two individual spaxels in the northeast (NE) and southwest (SW) directions within the butterfly region, from which the spectra shown in the bottom panel are extracted. The spectra are fit with two components, and the broad [NII] emission is detected at S/N > 4 (i.e., 5.5 in the NE and 4.1 in the SW). The light gray circle in the bottom-left panel is the MEGARA PSF (FWHM). For all the maps, north is up and east to the left.

Fig. 8. Similar to Fig. 7 but for the [NII]6583 line. The main difference between the [NII] and Hα (Fig. 7) panels is that the inner region of the [NII] secondary component, corresponding to the ILR in the Hα map, is absent in the forbidden [NII] line.

4.1.1. Hα primary (systemic) component

The primary (systemic) component covers almost the entire MEGARA FoV, extending primarily in the north-south direction. Its Hα emission highlights the base of the spiral arms extending toward the northwest and southeast, along with the bright central emission from our target. The Hα flux peaks at the same location as the R-band continuum emission, as shown in Fig. 7 (top-left panel). The continuum emission reveals the foreground star located to the north, very close to the Hα emission from the companion galaxy.

The velocity field exhibits the typical signature of a purely rotating disk (i.e., a “spider-like” pattern) with a velocity shear, v_shear, of ∼211 km s⁻¹ (similar to the rotational velocity derived for the molecular phase by Salomé et al. 2021; see our Table 3). The velocity shear is defined as half the difference between the median values of the lowest and highest 5% of the velocity distribution: $v_{\mathrm{shear}}-\frac{1}{2}(v_{max}^{5\%}-v_{min}^{5\%})$. The systemic (heliocentric) velocity, derived at the continuum peak position, is 18 331 km s⁻¹, corresponding to a redshift of z = 0.0611, in excellent agreement with the value derived from the molecular emission by Salomé et al. (2023, z = 0.0612; giving a redshift difference between these two tracers of 0.3 km s⁻¹) . At the position of the companion galaxy, the velocity field is slightly redshifted by ∼200 km s⁻¹, consistent with the velocity offset observed in the CO(1-0) emission (i.e., ∼210 km s⁻¹).

The velocity dispersion map reveals an X-shaped (or butterfly-shaped) structure, showing quite high velocity dispersion values (σ ∼ 200 km s⁻¹) in the central region (R ≲ 4 kpc) and perpendicular to the galaxy’s bar. The bar is known to be aligned in the north-south direction (see Ohta et al. 2007), and we suggest that the stream of velocity dispersion values with σ between 50 and 100 km s⁻¹ observed along the north-south direction may trace this structure (e.g., Hernandez et al. 2005; Erroz-Ferrer et al. 2015). The butterfly structure previously observed in the Low-Ionization Nuclear Emission-line Region (LINER) galaxy NGC 1052 by Cazzoli et al. (2022) (see also Hermosa Muñoz et al. 2024) may be associated with an outflow interacting with the surrounding ISM and/or the radio jet-ISM interaction, producing a region of enhanced velocity dispersion (i.e., turbulence). As shown in Fig. 6, the butterfly region is clearly spatially coincident with the molecular outflow traced in CO (i.e., Longinotti et al. 2023).

According to the observed (i.e., uncorrected for inclination) dynamical ratio, v/σ, between the velocity shear and the mean velocity dispersion (see Table 3), we found v/σ ∼ 2.3, which clearly indicates that the primary (systemic) component is dominated by rotation. Furthermore, we compare the Hα flux emission of the primary component with the molecular emission from the CO(1-0) transition observed with NOEMA (Salomé et al. 2023), as shown in Fig. 5. We found good agreement between the two at the positions of our source, IRAS17, and the companion galaxy. Although the interpretation of the Hα kinematic maps is not straightforward, our results suggest that the companion galaxy appears to extend over a larger area in Hα emission than in the molecular component. According to the Hα emission (see also Sect. 4.1.2), we identified an extended structure toward the southeast, resembling a tidal tail. This feature is characterized by redshifted velocities and low velocity dispersion values (σ ≲ 50 km s⁻¹; top-right panel in Fig. 7), may be connected to IRAS17.

4.1.2. Hα secondary component

The secondary component extends over a similar FoV as the primary component but is significantly brighter. Within this kinematic component, we can clearly distinguish the regions corresponding to the three sources: IRAS17, the foreground star, and the companion galaxy (see Fig. 7). The emission from the star has been removed from the map (see the gray flux contours in Fig. 7 toward the north). In this component, the emission from IRAS17 is the main contributor to the Hα flux, while the companion galaxy appears fainter than in the systemic flux map. Toward the north, an arm-like structure is visible, specially in the velocity field and velocity dispersion maps, likely associated with the companion galaxy. This structure is characterized by higher velocities (redshifted by ∼200 km s⁻¹) and lower velocity dispersion values (σ ∼ 100 km s⁻¹). Its contribution in the corresponding flux map is hard to discern.

The secondary component overall exhibits complex kinematics, with both the velocity field and velocity dispersion maps showing clear signs of turbulence and possible interaction effects. The width of this component is broader (> 250 km s⁻¹) than that of the systemic component (< 250 km s⁻¹), with a velocity shear of approximately 170 km s⁻¹ across the FoV (see Table 3).

However, this component is not characterized by rotation: no rotational velocity pattern is observed, and we find a v/σ of ∼0.5; it is thus dispersion-dominated. In the central region, we observed a median velocity very close to the systemic median velocity (i.e., with a velocity shift of ∼15 km s⁻¹). In contrast, the companion exhibits a (positive) velocity shift of ∼110 km s⁻¹ with respect to the systemic component. According to the features observed in the velocity dispersion map of the secondary component, we can identify three distinct regions:

• a.

  Lowσ (∼110 km s⁻¹): Located in the north, overlapping the companion’s position.

• b.

  Intermediateσ (∼280–380 km s⁻¹): Spanning a ∼28 arcsec² area (i.e., R ∼ 3.0″ ∼ 3.5 kpc; light green region in the map) around the central region, underlying the high-σ (i.e., butterfly) region.

• c.

  Highσ (∼600 km s⁻¹): Corresponding to the X-shaped (or butterfly) region, also seen in the primary component.

The intermediate-σ component (i.e., represented by the Gaussian shown in magenta in Fig. 3) is found in the region where the tertiary component (tracing the outflow) is added in the fit. In this area, the velocity field resembles that of the primary component, and the velocity dispersion values are indicative of an ILR. This interpretation is supported by the fact that the line widths are not as broad as those typically observed in Sy1 objects, where FWHM ≳ 2000 km s⁻¹ (corresponding to σ ≳ 850 km s⁻¹). In contrast, we observe a characteristic σ ∼ 280–380 km s⁻¹ (FWHM ∼ 700–900 km s⁻¹), consistent with typical values reported for ILRs (FWHM in the range 700–1200 kms⁻¹; Crenshaw & Kraemer 2007; Adhikari et al. 2016 and references therein). Furthermore, its velocity field closely matches the systemic velocity (Δv ∼ 0), reinforcing our interpretation. As presented in Adhikari et al. (2016), the ILR may represent an inner extension of the narrow line region.

In contrast, the high-σ region reveals significant turbulent motions surrounding the ILR, possibly induced by the outflow or representing the contribution from the ionization cones. Figure 7 (bottom panel) shows representative Hα-[NII] complex spectra extracted from two locations within the butterfly region, toward the northeast and southwest directions of the outflow. Both primary and very broad components are required to properly reproduce the total emission. In these regions, the broad Hα component dominates the flux emission.

Regarding the kinematics of the companion galaxy, the Hα emission in the secondary component does not align with that of the CO emission, appearing nearly perpendicular. Specifically, the CO emission forms a bridge aligned along the northwest-southeast direction, while the ionized gas appears to wrap around the molecular structure, moving toward the north-southwest and offset by approximately 40 degrees.

4.1.3. Hα tertiary (outflow) component

We identified the emission of the tertiary component with the ionized outflow due to its kinematic properties. It extends over a roughly circular area of approximately ∼6 arcsec² (observed radius of 1.4″ ∼ 1.65 kpc), corresponding to ∼8.3 kpc². This component is spatially resolved beyond the PSF region (see Sect. 5.1 for further details).

The velocity field shows a velocity shear of v_shear^T ∼ 140 km s⁻¹ and a median velocity dispersion of σ^T ∼ 956 km s⁻¹ (FWHM ∼ 2250 km s⁻¹). These velocity shifts are consistent with ionized outflows observed in local galaxies, whereas the velocity dispersion values are significantly higher than those typically found in local star-forming or AGN-dominated galaxies. Instead, they are more comparable to values observed in high-z sources, such as those in the SUPER survey and quasar samples (e.g., Kakkad et al. 2020, using [OIII] lines; FWHM = 600–3000 km s⁻¹). The velocity shift relative to the primary component, defined as Δv = v_T − v_P, ranges from −329.6 to −168.3 km s⁻¹, with a mean value of −252.5 km s⁻¹ (median: −243.9 km s⁻¹; see Table 3).

Although the geometry of the ionization cones and the outflow is difficult to define, under the assumption of a biconical outflow roughly perpendicular to the disk, the redshifted (receding) component may be partially obscured by the disk, while the blueshifted (approaching) component is less affected. In the case of IRAS17, the outflow velocity field shows blueshifted velocities with respect to the systemic component across the entire FoV, suggesting that the receding part is completely obscured by the disk. A few southeastern spaxels show the most extreme blueshifted velocities together with low velocity dispersions. This may indicate projection effects and obscuration, possibly associated with an interaction between the ionized outflow and denser ambient material in this region, which could lower the observed velocity dispersion.

In Sect. 5 we derive the physical and kinematic properties of the ionized outflow, such as its mass, kinetic energy, momentum, and velocity, v_out, based on this component. The main difference between the [NII] and Hα panels (Fig. 7) is that the inner region of the [NII] secondary component, corresponding to the ILR in the Hα map, is absent in the forbidden [NII] line.

4.2. Hβ-[OIII] complex in the MR-G and LR-V setups

Thanks to the MR-G and LR-V setups, we were able to characterize in detail the kinematics of the Hβ and [OIII] lines.In both setups, three kinematic components were required to fit the observed [OIII] emission (i.e., primary, secondary and tertiary) while three components were needed for the Hβ line in the MR-G setup, and only one component was considered in the LR-V (see Sect. 2.1 and Table 2). The secondary and tertiary components in [OIII] have been associated with a “modest” and a “fast” outflow, respectively, as the kinematics of the former are less extreme than that of the latter (see the following subsections).

As noted in Sect. 3.2.2, we linked the Hβ profile to the Hα kinematics¹², which enables us to characterize this fainter emission using three components: a narrow primary component, an intermediate secondary component associated with the ILR, and a broad tertiary feature associated with the outflow. Owing to the lower S/N of Hβ and [OIII] lines, their primary components are detected over a much smaller spatial extent than that of Hα.

In Figs. 9, 10, and 11 the kinematic maps of the Hβ-[OIII]5007 complex obtained using the MR-G and LR-V setups are shown. The [OIII] emission appears a bit more compact in the MR-G setup than in LR-V setup and this could be due to the longer exposure time (more than twice as long) used in the LR-V mode compared to MR-G, as well as better seeing conditions, which allowed us to better trace fainter emission found at larger radii. A proper characterization of the outflow extension in the different lines will be discussed in Sect. 5.1. Nevertheless, the main kinematic results are in agreement across both setups (see Tables 4 and 5) and three components were needed to properly fit the observed emission in both cases.

Fig. 9. Same as Fig. 7 but for the ionized component as traced by the Hβ line according to the MR-G setup.

Fig. 10. Same as Fig. 7 but for the ionized component as traced by the [OIII]5007 line according to the MR-G setup.

Fig. 11. Same as Fig. 7 but for the ionized component as traced by the Hβ (top) and [OIII]5007 (bottom) complex as derived in the LR-V setup.

As was mentioned in Sect. 3.2.2, the Hβ emission was modeled with a three-component fit in the MR-G setup (Fig. 9). In the LR-V setup, however, the line lies at the blue edge of the spectral range, close to the sensitivity limit, making the fit less reliable. Although some extended emission is detected in LR-V (but not in MR-G), we do not draw conclusions from it. Therefore, we adopt the kinematic results derived from the MR-G data, as they are more robust.

4.2.1. Hβ-[OIII] primary (systemic) components

The Hβ primary (systemic) component closely follows the global Hα kinematics and is consistent with the systemic component. The [OIII] primary (systemic) component extends in a smaller region than that covered by the Hα line (see Fig. 10 and Sect. 5.1 for further details), although it still overlaps with the IRAS17 emission, and the companion galaxy shows detectable emission. The main kinematic features are consistent across both the MR and LR setups, including a rotation-dominated velocity field and a velocity dispersion map that reveals a ring-like structure with enhanced σ values. However, we find a clear difference between the kinematic maps derived from the [OIII] and Hα emission lines, particularly in the velocity dispersion distribution. As shown in Fig. 12, the regions of enhanced σ traced by the [OIII] emission coincide with the edges of the Hα butterfly structure, surrounding the nuclear outflow.

Fig. 12. Hα velocity dispersion map of the primary component, with the [OIII] velocity dispersion contours (white) of the primary component overlaid, obtained from the LR-V setup.

4.2.2. Hβ-[OIII] secondary components

The secondary Hβ component follows the global Hα kinematics and is identified with the ILR, consistent with the Hα emission (se Table 2). The secondary component observed in [OIII] has been identified with an outflow. In particular, for this component in the MR-G mode, we derived a mean (median) velocity offset (with respect to the spatially resolved systemic velocity map) of Δv = −174.0 ± 9.0 km s⁻¹ (median: −157.6 km s⁻¹), and a mean (median) velocity dispersion of 223.4 ± 5.3 km s⁻¹ (median: 221.7 km s⁻¹). The less blueshifted¹³ region is located to the northwest, while the most blueshifted region lies to the southeast. The less blueshifted side displays a lower velocity dispersion, whereas the most blueshifted side exhibits higher σ, consistent with the interpretation that the latter suffers less extinction. This behavior is consistent with an orientation-dependent obscuration scenario: the most blueshifted component likely traces the near side of the outflow, which is less affected by extinction, allowing us to observe gas closer to the nucleus and with higher turbulence. In contrast, the less blueshifted components are plausibly associated with gas on the far side of the outflow (or partially embedded in the galactic disk), where line-of-sight obscuration preferentially suppresses emission from the most turbulent regions, leading to a lower observed velocity dispersion. These kinematic differences therefore reflect orientation and obscuration effects within a dusty circumnuclear environment. The same pattern is found for the [OIII] line observed in the LR-V setup, for which we also derived similar kinematic values (see Table 5). Both [OIII] outflows identified by the secondary components are resolved, according to the PSF values of each setup.

4.2.3. Hβ-[OIII] tertiary components

The Hβ tertiary component (detected in the MR-G setup) traces the global Hα kinematics and is likewise associated with the fast ionized outflow. In both setups, the tertiary component of the [OIII] line appears more compact than its secondary component, while exhibiting more extreme kinematics, consistent with a faster outflow. In particular, the [OIII] tertiary component is unresolved in the MR-G setup, whereas it is resolved in the LR-V setup (see Sect. 5.1). The kinematics of this outflow shows typical median values of velocity offset Δv = −540 km s⁻¹ and velocity dispersion σ ∼ 400 km s⁻¹ in the MR-G setup, as derived from the integrated emission within the PSF region. Therefore, for the resolved emission in the LR-V setup, we derived a median velocity offset of Δv ∼ −536 km s⁻¹, and a median velocity dispersion of ∼520 km s⁻¹ (see Table 5). The kinematic results obtained for both setups and emission lines are in mutual agreement and are reported in Tables 4 and 5.

5. Deriving the physical properties of the outflow(s)

5.1. Constraining the intrinsic size of the outflow(s)

We derived the size of the ionized outflow using the Hα and [OIII] lines by applying a fixed percentage of the peak flux to each outflow map, in order to avoid any dependence on the S/N of the individual lines. This method uses a flux threshold (e.g., 10% or 50% of the peak flux) to eliminate the S/N dependency of the size measurement. From these measurements, we computed an equivalent circular radius by converting the area enclosed by the isophote into the radius of a circle with the same area: R_eq=$\sqrt{\frac{\mathrm{Area}}{\pi }}.$ In our analysis, we adopted the 10% peak-flux isophote as the primary measure of the outflow size, as it captures the full spatial extent of the faint, extended emission and is less sensitive to flux concentration in the central regions compared to the 50% peak-flux isophote. We counted the number of pixels above this threshold and calculated the total area as Area = N_pix × pixscale². This approach considers only the number of pixels above the threshold, not the sum of their fluxes.

To derive the intrinsic size of the outflow, we accounted for the PSF size measured at 10% of the peak emission, given by PSF_10%  =  1.82 × PSF_FWHM. Using the 10% flux cutoff, we obtained half the PSF (PSF_10%/2) of 1.09″ for the LR-R and MR-G setups, and 0.82″ for the LR-V setup.

The intrinsic outflow sizes derived using this method for the Hα and [OIII] components are in good agreement, with values ≲1 kpc for both tracers. Specifically, for the Hα outflow we obtained an intrinsic size of 0.98 ± 0.32 kpc, while for [OIII] in the MR-G setup we measured 0.70 ± 0.41 kpc for the secondary component. The tertiary [OIII] component is unresolved in MR-G, for which we derived an upper limit to its size based on the PSF, < 1.29 kpc. In the LR-V setup, both components are resolved: for the secondary component, we measured a size of 0.70 ± 0.26 kpc, while for the tertiary component we obtained 0.47 ± 0.25 kpc. All these values have been used to derive the physical outflow parameters presented in Tables 6–10 (see Sect. 5.3).

5.2. Outflow electron density derivation

The electron density of ionized gas, particularly in both the systemic (disk) and outflowing components of galaxies, is a crucial parameter for understanding the physical conditions of the ionized ISM and for accurately determining the properties of outflows like its mass, mass outflow rate and kinetic energy. However, several works (i.e., Harrison et al. 2018; Kakkad et al. 2018) pointed out, the electron density, n_e, is the parameter with the largest uncertainties in these calculations. In literature, different methods can be used to derive this parameter (Davies et al. 2020). A typical method largely used is based on the [SII] doublet ratio, [SII]λ6717/λ6731 (hereafter the “[SII] method”; e.g., Osterbrock & Ferland 2006; Kakkad et al. 2018; Rodríguez del Pino et al. 2019; Venturi et al. 2021). This method suffers from saturation at high density and cannot probe high densities because collisional de-excitation dominates above 10⁴ cm⁻³: this ratio is thus only sensitive to densities in the range 50 < n_e < 5000 cm⁻³, giving lower electron densities by down to an order of magnitude than those obtained with other methods (see Davies et al. 2020). Thus, the values obtained from this method could lead to an overestimation of the mass of the outflow, as well as of the mass outflow rate. In our case, we could not isolate the outflow emission in the [SII]λλ6716,6731 doublet by fitting each line with a single-component (1c) model (see Fig. B.1, top panel). Therefore, as a first approximation, we used the flux ratio derived from the 1c fit (primary component) to calculate the electron density in each spaxel following the prescription of Sanders et al. (2015):

$$\begin{array}{c}\hfill n_{\mathrm{e}}=\frac{cR-ab}{a-R},\end{array}$$(1)

where R = [S II]λ6717/λ6731 and a = 0.4315, b = 2107, c = 627.1. According to this method, we derived a mean (median) electron density n_e = 266 ± 15 (245) cm⁻³ across the entire FoV (see the spatially resolved n_e map in Fig. 13, top), while in the inner region (i.e., the area covered by the outflow) we obtained a higher value of n_e ∼ 350 cm⁻³. The electron density profile is centrally peaked and decreases with increasing radius, although a few spaxels show relatively high values (> 500 cm⁻³), possibly due to lower S/N in the emission, which leads to more uncertain results. As we know, this method tends to underestimate the electron density because a significant fraction of the [SII] emission preferentially originates from the outer, less dense regions of the ionized outflow (see Kakkad et al. 2018). Furthermore, the lack of a detectable broad component in the [SII] emission lines prevents a reliable kinematic decomposition of this doublet, making the [SII]-based electron density estimate potentially unrepresentative of the outflowing gas. In particular, the derived value should be regarded as a lower limit, which would imply unrealistically high outflow masses if directly adopted. For this reason, we also derived the electron density of the outflow using the method described by Baron & Netzer (2019), who derived the n_e from the ionization parameter value, U (hereafter the “log U method”; see also Davies et al. 2020; Peralta de Arriba et al. 2023; Esposito et al. 2024). This method is more sensitive to the high-density, fully ionized regions, where Hα and [OIII] emission lines originate. It derives the n_e parameter using the ratio between the emission lines such as Hα/[NII],and [OIII]/Hβ, which is typically strong in AGNs. According to this method we can derive the ionization parameter, U, defined as the number of ionizing photons per atom: U = Q_H/(4πr²n_Hc), where Q_H is the rate of hydrogen-ionizing photons (s⁻¹ units), r is the distance from the ionizing source, n_H∼ n_e is the hydrogen density, and c is the speed of light. The electron density can be written as a function of the AGN luminosity (i.e., L_bol = 10^44.5 erg s⁻¹ from Longinotti et al. 2015), the distance from the AGN (r) and the ionization parameter (U):

$$\begin{array}{c}\hfill \begin{array}{c}\hfill n_{\mathrm{e}}\approx 3.2\left(\frac{L_{\mathrm{bol}}}{{10}^{45}\mathrm{e}rg{s}^{-1}}\right){\left(\frac{r}{1\mathrm{kpc}}\right)}^{-2}\left(\frac{1}{U}\right)\mathrm{c}{m}^{-3}.\end{array}\end{array}$$(2)

Fig. 13. Electron density (ne) maps in units of cm−3 derived using the [SII] (top) and log U bottom) methods. The [SII] method is applied to the total [SII]λλ6717,6731 line profile, while the log U map sums the flux contributions from the secondary and tertiary components. Typical mean and median values, together with the standard deviation and the median absolute deviation, respectively, are shown in the insets.

The U parameter is defined according to the following expression:

$$\begin{array}{cc}\hfill logU=& -3.766+0.191log\left(\frac{[\mathrm{O}III]}{\mathrm{H}\beta }\right)+0.778{log}^2\left(\frac{[\mathrm{O}III]}{\mathrm{H}\beta }\right)\hfill \\ \hfill & -0.251log\left(\frac{[\mathrm{N}II]}{\mathrm{H}\alpha }\right)+0.342{log}^2\left(\frac{[\mathrm{N}II]}{\mathrm{H}\alpha }\right).\hfill \end{array}$$(3)

In our case, we derived the ionization parameter, U, and the electron density, n_e, for the outflowing component by summing the flux contributions of the secondary and tertiary components. In particular, we used the secondary and tertiary components for the Hα, Hβ, and [OIII] lines, while for the [NII] line we only considered the secondary component associated with the outflow. Under these assumptions, we derived a median electron density, n_e ≈ 2500 cm⁻³ (Fig. 13, bottom), which we adopt for the derivation of the outflow parameters (see Tables 6–10). This value is roughly a factor of ∼10 larger than that derived using the [SII] method, and is in agreement with those found in the outflow component in nearby AGN systems (e.g., Rose et al. 2018; Singha et al. 2022). It is worth mentioning that, as noted in several works, the electron density can vary across the FoV; therefore, assuming a constant electron density can introduce an additional systematic uncertainty (e.g., Kakkad et al. 2018).

For comparison, we also derived the outflow parameters adopting an electron density of 500 cm⁻³, comparable to our [SII]-based estimate, and representative of local active galaxies (e.g., (ultra)luminous infrared galaxies, LINERs, and type 2 AGNs in the MAGNUM survey; Fluetsch et al. 2021; Venturi et al. 2021; Hermosa Muñoz et al. 2024). The classical [SII]λ6716/λ6731 ratio typically yields lower electron densities, characteristic of the diffuse ionized medium. In contrast, applying the Baron & Netzer (2019) method results in higher densities, as it combines multiple line diagnostics with photoionization modeling. The discrepancy reflects the different sensitivity regimes: the [SII] ratio becomes insensitive above a few ×10³ cm⁻³, tracing mainly extended, low-density gas, whereas the Baron & Netzer (2019) approach probes denser clumps that contribute significantly to the emission but are underrepresented in the [SII]-based estimate. Together, both measurements indicate a multiphase ionized medium, where compact, high-density knots coexist with more diffuse gas.

5.3. Deriving outflow parameters using the Hα line

To understand the impact of the outflow in this galaxy, we first need to derive its main properties, such as the mass of the outflow, M_out, the mass outflow rate, Ṁ_out, the kinetic energy, E_kin, the kinetic power, Ė_kin, and the momentum rate, Ṗ_out. These quantities were derived as in Fiore et al. (2017, see also Cresci et al. 2017 and Venturi et al. 2023, and references therein). We first derived the mass of the outflow from the extinction corrected luminosity as derived from the integrated total Hα flux of the tertiary component. This is a more robust estimation of the mass compared to that obtained using the [OIII] line because its luminosity does not depend on the gas metallicity or on the energy of the ionizing photons (e.g., Carniani et al. 2015; Venturi et al. 2021). To correct the Hα luminosity for the extinction, we used the Calzetti et al. (2000) attenuation law (assuming the typical R_V = 4.05 for starbursts) and an intrinsic ratio (Hα/Hβ)₀ = 3.1 (for an electron temperature of T_e = 10⁴ K; Osterbrock & Ferland 2006). We derived a (median) value of A_V = 2.73 mag for the outflow (see Appendix A for details). We then computed the mass of the ionized outflow using the relation

$$\begin{array}{c}\hfill M_{\mathrm{out}}=6.1\times {10}^8\left(\frac{L_{\mathrm{H}\alpha }}{{10}^{44}\text{erg}{\text{s}}^{-1}}\right)\left(\frac{500{\text{cm}}^{-3}}{n_{\mathrm{e}}}\right){\mathrm{M}}_{\odot },\end{array}$$(4)

where n_e was derived in the manner described in the previous section. We assumed two different values (i.e., n_e = 500 and 2500 cm⁻³; see Sect. 5.2 for details) and the outflow parameters obtained from using these values are shown in Table 6. In both cases, we obtained similar outflow masses, in the range ∼0.6–3.2 × 10⁷ M_⊙. M_out is the ionized gas outflow mass integrated within the outflow region.

A key parameter for the characterization of the outflow kinematics is the maximum outflow velocity, v_out, derived as in Rupke et al. (2005): $v_{\mathrm{out}}=|\mathrm{\Delta }v|+\frac{FWH{M}_{\mathrm{out}}}{2}=|\mathrm{\Delta }v|+1.18{\sigma }_{\mathrm{out}}$, where Δv = v_T − v_P is the maximum velocity difference between the outflowing (i.e., v_T) and primary (i.e., v_P) components, and σ_out is the median velocity dispersion of the outflow.

We then estimated the mass outflow rate, Ṁ_out, at a given radius, R_out, using the formula from Fiore et al. (2017), see also Lutz et al. 2020):

$$\begin{array}{c}\hfill {\dot{M}}_{\mathrm{out}}=C·\frac{M_{\mathrm{out}}{v}_{\mathrm{out}}}{R_{\mathrm{out}}},\end{array}$$(5)

where M_out was derived in Eq. (4) and R_out is the size of the outflow as derived in Sect. 5.1 (i.e., R(Hα) = 0.98 ± 0.32 kpc). The C factor depends on the adopted outflow history. A value of C = 1 assumes that the outflow is the result of a single explosive event in which clouds were ejected and continued to expand at a constant mass outflow rate until the present. In some cases, authors have assumed C = 3, which implies that the outflow is continuously replenished by clouds ejected from the galactic gaseous disk (see Lutz et al. 2020 for further details; see also Fiore et al. 2017; Hermosa Muñoz et al. 2024). This assumption corresponds to a constant volume density of the outflowing gas and a mass outflow rate decreasing to zero over time. However, this scenario is at odds with the presence of a UFO in our object, as the continuously replenished outflow implied by C = 3 assumes a slowly declining mass outflow rate fed by the galactic disk, whereas the UFO indicates a recent, highly energetic ejection directly from the AGN, inconsistent with a long-lived, disk-fed outflow. Therefore, in our case, we considered the same assumption used in Marasco et al. (2020) (C = 1) who assumed a constant mass outflow rate during the flow time, defined as R_out/v_out, which leads to a decreasing density of the outflowing gas. This approach is consistent with the “time-averaged thin-shell” model (e.g., Rupke et al. 2005), and is widely used to derive outflow rates in both ionized and neutral gas phases (e.g., Heckman et al. 2015; González-Alfonso et al. 2017; Catalán-Torrecilla et al. 2020; Marasco et al. 2020). This assumption (C = 1) also represents a more conservative approach for deriving the outflow parameters. Since the C factor is constant, it does not affect the main conclusions of this work, although we have taken the different approaches into account to compare our results with previous studies (see Longinotti et al. 2023). Under these assumptions, we derived a mass outflow rate between ∼10 and 50 M_⊙ yr⁻¹ using n_e = 2500 or 500 cm⁻³, respectively, as reported in Table 6.

We then computed the kinetic energy, E_kin, and the kinetic power, Ė_kin, of the outflow as in Rose et al. (2018, also Santoro et al. 2020):

$$\begin{array}{c}\hfill E_{\text{kin}}=\frac{1}{2}{\sigma }_{\mathrm{out}}^2{M}_{\mathrm{out}},\end{array}$$(6)

$$\begin{array}{c}\hfill {\dot{E}}_{\text{kin}}=\frac{{\dot{M}}_{\mathrm{out}}}{2}(v_{\mathrm{out}}^2+3{\sigma }_{\mathrm{out}}^2).\end{array}$$(7)

The corresponding momentum rate of the outflow was derived as Ṗ = Ṁ_out v_out, which yields a few ×10³⁵ cm g s⁻². When we compare the momentum rate of the outflow with the radiative momentum rate of the AGN, defined as Ṗ_rad = L_bol/c = 1.7 × 10³⁴ cm g s⁻² for our source, we obtained a ratio Ṗ_Hα/Ṗ_rad between ∼5 and 27 (depending on the assumed n_e value; see Table 6). Under our assumptions, we derive Ṗ_Hα/Ṗ_rad = 5.3 ± 3.2.

According to our results, we can now compare the momentum rate values of the outflow (i.e., Ṗ_out/Ṗ_rad) obtained for the different gas phases (molecular, ionized and highly ionized, derived, respectively, in Longinotti et al. 2023, this work and Longinotti et al. 2015) as a function of the outflow velocity (see Fig. 14). This plot is useful for comparing the energy transfer predicted by different models between UFOs and galaxy-scale outflows, under the assumption that the molecular phase carries most of the outflow mass (i.e., ${\dot{P}}_{\mathrm{out}}^{\mathrm{tot}}\sim {\dot{P}}_{\mathrm{out}}^{\mathrm{mol}}$; Tombesi et al. 2015; Feruglio et al. 2015).

Fig. 14. Force of the molecular, ionized, and X-ray phases of the outflow, Ṗout, in IRAS17, normalized to the source’s X-ray momentum rate, ṖUFO (i.e., ṖUFO/Ṗrad = 1.87 ± 0.62; see Longinotti et al. 2018), is shown as a function of the outflow velocities. The results for IRAS17 derived from NOEMA, MEGARA (i.e., Hα using ne = 2500 cm−3), and XMM-Newton are shown with blue circles. To compare results across the three phases, we used the same formula for vout as in Longinotti et al. (2023, i.e., Δv + 2 × σB) and normalized their Ṁout values assuming C = 1. Filled circles in orange and magenta represent the results derived for 2c component of the [OIII] emission (i.e., slow outflow) in the LR and MR setups, respectively. Two values are plotted for each setup (connected by a solid line), corresponding to the two assumed electron densities, ne = 500 and 2500 cm−3 (see text for details). Empty circles, following the same color code as for the 2c component, indicate the 3c component of the [OIII] emission (i.e., fast outflow) in both setups. The fast [OIII] outflow is spatially resolved in the LR setup (empty orange circle), while it remains unresolved in the MR setup (empty magenta circles), for which we show the two results (lower limits) derived assuming the two ne values. Results for three other well-known sources for which the energy-conserving regime has been observed, Mrk 231, IRAS 05189+2524, and IRAS11119+3257 (all ultraluminous infrared galaxies, shown using black, green, and red squares, respectively), are also included. Their data are extracted from Feruglio et al. (2015), Lutz et al. (2020), and Tombesi et al. (2015), respectively (see also Marasco et al. 2020 and references therein). The dashed blue and red lines represent the two linear fits marking the energy-conserving outflow regime: the former obtained using CO, Hα, and X-ray data, and the latter using these three points plus the [OIII] result from the LR-V setup (empty orange circle). The dotted gray line indicates the prediction for “momentum-conserving” outflows (i.e., Ṗout/ṖUFO ∼ 1).

To compare our results with those from previous works (i.e., Longinotti et al. 2023), we needed to homogenize the formulas used in this study with those adopted in earlier ones. In particular, Longinotti et al. (2023) used the formula proposed by Rupke & Veilleux (2013) for the derivation of the outflow velocity, which assumes v_out = Δv + 2 σ_out, while in this work we followed the definition by Rupke et al. (2005). The former gives larger outflow velocity than the latter used in this work. In the case of MEGARA data, we would have derived larger outflow velocity by a factor of ∼1.5. On the other hand, the C factor used by Longinotti et al. (2023) to derive the mass outflow rate is 3, whereas in our case we assumed C = 1. Both these assumptions increase the velocity and the energetics of the outflow, making the outflow to appear more extreme. Since we preferred to be more conservative in the derivation of the outflow physical parameters, we assumed the C = 1 as first assumption maintaining the (higher) v_out formula by Rupke & Veilleux (2013) for a direct comparison of the results.

We finally report the results by Longinotti et al. (2018) and Longinotti et al. (2023) adjusted according to the C = 1 assumption for a direct comparison with our results (i.e., Ṁ_out and Ṗ_out; see Table 6). The Hα outflow, as measured by MEGARA, confirms the presence of the energy-conserving regime for the ionized gas phase, previously identified for the molecular emission, in which a momentum boost is required as a function of the outflow velocity (i.e., $\dot{P}\propto v_{\mathrm{out}}^{-1}$). The ionized gas retains a significant fraction of the energy from the inner wind, suggesting that AGN feedback is dynamically important and can have a strong impact on galaxy evolution. The results are shown in Fig. 14.

5.4. Deriving outflow parameters from the [OIII] line

We also determined the mass of gas contained in the outflow using the [OIII] lines, which can be expressed with the relation (see Carniani et al. 2015; Fiore et al. 2017; Venturi et al. 2023)

$$\begin{array}{c}\hfill M_{\mathrm{out}}^{[\mathrm{OIII}]}=8.0\times {10}^7\left(\frac{L_{[\mathrm{OIII}]}}{{10}^{44}\text{erg}{\text{s}}^{-1}}\right)\left(\frac{500{\text{cm}}^{-3}}{<n_e>}\right)\frac{\mathcal{CF}}{{10}^{[\mathrm{O}/\mathrm{H}]-{[\mathrm{O}/\mathrm{H}]}_{\odot }}}{\mathrm{M}}_{\odot },\end{array}$$(8)

where < n_e> is the electron density averaged over the ionized outflow volume (i.e., < n_e> = ∫n_efdV/∫fdV) and 𝒞ℱ = < n_e> ²/< n_e²> is the so-called “condensation factor.” Under the simplifying assumption that all ionizing gas clouds have the same density, we set 𝒞ℱ = 1, thereby eliminating the dependence of the outflow mass on the filling factor of the emitting clouds. The [O/H] represents the oxygen abundance relative to the solar value, [O/H]_⊙. We assumed that the metallicity of the outflowing material is always assumed solar ([O/H] = [O/H]_⊙). The L_[OIII] is the outflow luminosity corrected for the extinction (see Appendix A). Since the [OIII] lines are highly sensitive to both temperature and ionization, they are poor tracers of the outflow mass. Therefore, this method can, at best, provide only an order-of-magnitude estimate of the outflow properties. Finally, in order to estimate the total outflow mass using [OIII] lines, we considered the approximation used by Fiore et al. (2017), for which $M_{\mathrm{tot}}^{\mathrm{out}}=3\times M_{\mathrm{out}}^{[\mathrm{OIII}]}$ (see also Peralta de Arriba et al. 2023). In order to compute the energetics of the ionized outflow, we considered the same equations used for the Hα line (see Eqs. 5, 6, 7). The outflow mass derived using the [OIII] line in the LR-V setup is higher than that obtained with the MR-G configuration by a factor of two. This difference can be attributed to the longer exposure time in LR-V (twice that used for the MR-G grism) and better seeing conditions, which result in higher sensitivity to lower surface brightness, mostly found in the outer regions of the galaxy (see Figs. 10 and 11). Nevertheless, the kinematic values derived for the outflow are consistent (within the errors) between the two setups.

In the [OIII] emission, we identified two distinct outflows. The secondary outflow component is less kinematically extreme yet more spatially extended (R_out ∼ 0.7 kpc) than its tertiary counterpart in both the MR-G and LR-V setups. The tertiary outflow component is compact (R_out ∼ 0.5 kpc) and kinematically extreme, characterized by large velocity shifts and broad line widths. It remains unresolved at the spatial resolution of the MR-G setup, consistent with the nuclear base of the outflow, but is resolved in the LR-V setup. In contrast, the lower-ionization lines show little evidence of outflowing gas: [NII] exhibits only a faint broad component, while [SII] shows none. This indicates that these lines are intrinsically weak in the highly ionized outflowing gas. Consequently, high-ionization lines such as [OIII] provide the clearest tracers of the outflow, and the limited broad emission in [NII] and [SII] reflects their low intrinsic flux rather than limitations in spectral decomposition on a spaxel-by-spaxel basis. Consequently, [OIII], along with Hα, provides the most robust tracer of the outflow kinematics in this source. In Tables 7–10 we present the results obtained for the [OIII] line with both the setups.

6. Discussion

6.1. Tracing stratification within the ionized gas phase

Our analysis of the MEGARA data show that the ionized outflow in IRAS17 is consistent with a multiphase, shock-driven wind, as predicted by theoretical models of AGN energy-conserving outflows. A compact, high-velocity phase is detected in both Hα (R_out ∼ 1.0 kpc and v_out ∼ 1460 km s⁻¹) and [OIII] emission (i.e., tertiary component in the LR-V setup with R_out ∼ 0.5 kpc, v_out ∼ 1240 km s⁻¹), tracing ionized gas compressed and accelerated by the nuclear wind. A secondary [OIII] component with more modest velocities (v_out ∼ 450 km s⁻¹ in the LR-V setup) at R_out ∼ 0.7 kpc traces a diffuse, radiatively cooled envelope that has become momentum-dominated. This configuration naturally reproduces a stratified structure, where the inner regions efficiently transfer kinetic energy into the dense ambient medium, while the outer layers represent a slower, cooled phase.

The Hα and [OIII] lines trace gas under different physical conditions within the outflow. Hα emission arises from denser ionized material at T ∼ 10⁴ K, where recombination is efficient, typically associated with photoionized clumps or shock-compressed regions. In contrast, [OIII] originates from more highly ionized and generally lower-density gas exposed to the AGN radiation field or fast shocks. Owing to its relatively high critical density (n_e ∼ 6.8 × 10⁵ cm⁻³), [OIII] is efficiently emitted in diffuse ionized regions, whereas in denser environments the line becomes suppressed by collisional de-excitation. These different excitation and density requirements naturally lead to distinct radial extents for the two tracers. The more compact [OIII] emission (R_out ∼ 0.5–0.7 kpc) originates where the gas remains sufficiently diffuse and strongly illuminated by the AGN. In contrast, Hα emission (R_out ∼ 1 kpc) is detected farther out because denser ionized clumps can survive and cool efficiently at larger distances, including regions where shocks compress the gas within the outflow.

In IRAS17, the detection of higher velocities in Hα suggests that this line is tracing denser or partially obscured clumps of ionized gas. These clumps remain visible in Hα, while the [OIII] emission is partially suppressed due to collisional de-excitation in dense regions and attenuation by dust. The velocities derived from both ionized lines are consistent with those measured for the molecular outflow with NOEMA (∼1000–1800 km s⁻¹ from Longinotti et al. 2023), supporting a coherent, multiphase outflow powered by the active nucleus (see the next section).

As found by Fiore et al. (2017) and other works (e.g., Harrison et al. 2018), the coupling efficiency Ė_out/L_bol can cover a large range of values, and depends on the gas phase considered. In particular, theoretical AGN feedback models predict a typical coupling efficiency of ∼5% between AGN power and the surrounding dense gas, corresponding to the molecular gas phase (Di Matteo et al. 2005; Hopkins & Elvis 2010; Zubovas & King 2012; Costa et al. 2014; see also Fig. 2 in Harrison et al. 2018, which compares theoretical predictions to observations). This coupling efficiency is required to reproduce the observed M_BH-σ relation in local galaxies. In contrast, a much lower coupling efficiency (< 1%, i.e., 0.01–0.1%) is typically observed in ionized outflows (e.g., Carniani et al. 2015; Bischetti et al. 2017; Musiimenta et al. 2023). Nevertheless, as shown in the simulations of Hopkins & Elvis (2010) and summarized in Harrison et al. (2018), even such low coupling efficiencies in ionized winds can still have a significant impact on the host galaxy.

In IRAS17, the comparison between the kinetic powers of the two [OIII] outflow (i.e., slow and fast) components indicates that both can be powered by the AGN. The fast component dominates the energetics, with Ė_fast = 0.8% × L_bol, while the slower component requires only a small fraction of this energy (Ė_slow = 0.04% × L_bol, corresponding to ∼5% × Ė_fast). The resulting energy ratio indicates that the AGN has sufficient power to drive both phases simultaneously according to simulations. Within this framework, the slow [OIII] outflow may also naturally arise from the interaction of the fast AGN-driven outflow with the surrounding ISM, through shocks and partial energy transfer.

The inferred kinetic coupling efficiencies for the slow and fast [OIII] outflows (i.e., Ė_out/L_bol = 0.04%–0.8%) together with the moderate momentum boosts (Ṗ_out/Ṗ_rad ∼ 0.3–3, respectively) indicate that the ionized outflows operate in a momentum-driven or weakly “energy-driven” regime, respectively. In this scenario, only a small fraction of the AGN radiative or mechanical power is converted into bulk kinetic energy of the ionized gas, while the momentum flux of the fast component slightly exceeds the single-scattering limit (Ṗ_out ∼ Ṗ_rad), likely due to shocks and multiple interactions within a clumpy ISM. The slow outflow shows a globally blueshifted kinematic pattern and might trace gas at lower velocity within the same ionized outflow due to projection effects, while the fast outflow traces outflowing gas directed more toward the line of sight, and therefore showing more extreme blueshifted velocities. Both the slow and fast [OIII] components can be considered as part of a single physical outflow, with kinematic differences arising from projection effects, orientation, and interactions with the surrounding ISM. The lower-velocity component, with a smaller momentum boost, can naturally arise from the interaction of the fast AGN-driven outflow with the surrounding ISM. Such efficiencies and momentum fluxes are commonly observed in ionized AGN-driven outflows and are consistent with theoretical expectations for inefficient coupling between the AGN and the ambient gas (e.g., Harrison et al. 2018).

6.2. Connecting the different outflow phases

In the model proposed by Faucher-Giguère & Quataert (2012), AGN-driven outflows generally operate in the energy-driven regime, in which the thermal energy of the shocked wind bubble is efficiently transferred to the swept-up ISM, powering the large-scale expansion (see also Zubovas & King 2012; King & Pounds 2015). In this scenario, the shocked wind itself remains largely adiabatic, with a cooling time much longer than the dynamical time (t_cool ≫ t_dyn), so it retains its thermal energy and drives the outflow expansion. The multiphase outflow forms primarily in the region affected by the forward shock, where the expanding hot bubble encounters and sweeps up the ambient ISM (i.e., Zubovas & King 2014; Richings & Faucher-Giguère 2018b). Here, the shocked ISM can cool efficiently if its radiative cooling time is shorter than the dynamical time (t_cool ≲ t_dyn), losing thermal energy and condensing into cooler, denser phases. This process produces a multiphase structure consisting of hot X-ray emitting plasma (from the shocked wind), warm ionized gas (traced by optical lines such as [OIII] and Hα), and cold molecular material (CO, OH, etc.) originating from the radiatively cooled shocked ISM. Hydrodynamic instabilities (Rayleigh-Taylor and Kelvin-Helmholtz) at the contact discontinuity between the shocked wind and shocked ISM (see Fig. 1 in Zubovas & King 2012) promote turbulent mixing, enhancing the coexistence of different temperature and density phases. In contrast, the reverse shock and the inner shocked wind region remain dominated by extremely hot (T ≳ 10⁹ K) gas, which is largely adiabatic and not directly involved in forming the multiphase material.

Fast molecular outflows may result from the entrainment of preexisting molecular clouds by a hot, high-velocity AGN wind, although simulations (e.g., Scannapieco & Brüggen 2015; Richings & Faucher-Giguère 2018a) show that this process is often inefficient due to cloud destruction by hydrodynamical instabilities and thermal evaporation. The destruction of these clouds would also naturally contribute to the observed deficit of molecular gas in the innermost regions of active nuclei (e.g., Alonso-Herrero et al. 2021; Esposito et al. 2024), while mass-loading the hot wind and promoting cooling instabilities.

Molecular gas can also form in situ through the cooling of previously shocked ionized or atomic material, as proposed by Richings & Faucher-Giguère (2018a,b). In these models, the swept-up gas behind the forward shock cools rapidly and may condense into a cold phase during the expansion of the AGN-driven wind. This picture is compatible with energy-conserving wind frameworks (e.g., King & Pounds 2015, Zubovas et al. 2022, 2024; Zubovas & Tartėnas 2025), in which the fast nuclear wind creates the conditions required for the gas to cool and form molecules within the outflow itself. Specifically: (1) the outflow begins as a hot, ionized, high-velocity wind launched close to the AGN; (2) shocks arise when this wind interacts with the surrounding ISM and the gas starts to cool; (3) in regions of higher density or enhanced thermal instabilities, the cooling becomes rapid enough to produce a cold phase; (4) thus, while preexisting molecular clouds may be destroyed early in the interaction, the molecular component observed at larger radii is likely formed directly within the outflow through rapid cooling and condensation. This in situ formation explains the presence of molecular gas observed at large distances from the nucleus and under conditions where molecules would otherwise be expected to be destroyed, requiring reformation processes. This is likely to happen in IRAS17: this scenario is favored by the spatial distribution of the ionized and molecular outflows (see Fig. 15), where the ionized outflow (R_ion ∼ 1 kpc and ∼0.5 kpc for the Hα and [OIII] lines) is found within the region of the molecular outflow (R_CO = 2.8 ± 0.3 kpc; Longinotti et al. 2023).

Fig. 15. Velocity field shift map, Δv (in km s−1), between the outflow and primary components. These maps are derived from the Hα tertiary component (left) and from the [OIII] line for the tertiary (resolved) component (right) using LR-V setup. A flux-emission threshold of 10% of the peak value has been applied to the outflow kinematic maps. The gray circle marks the PSF size at 10% of the flux peak (PSF10%; see text for details). The blue and red contours trace the approaching and receding molecular outflowing gas, respectively, as presented by Longinotti et al. (2023, see their Fig. 4). The ionized-gas emission nicely matches the region of the molecular outflow and shows blueshifted velocities (Δv < 0 km s−1; the minimum value is reported in each panel). The intrinsic outflow size for each component is also indicated in the corresponding panel.

Faucher-Giguère & Quataert (2012) show that the shocked ambient medium cools rapidly behind the forward shock, forming a thin, dense shell prone to instabilities and fragmentation. This shell is the site where cold gas can form within the outflow. This is also supported by the work of Nims et al. (2015) who claimed that most of the emission is produced by the forward shock driven into the ambient ISM.

As was mentioned in the previous section, the ionized outflow velocities are consistent with the molecular outflow velocities measured with NOEMA (average combined value 1280 ± 480 km s⁻¹ from Longinotti et al. 2023), indicating that both molecular and ionized gas belong to the same multiphase AGN-driven outflow. These cold, massive clouds can survive and propagate farther from the nucleus because they are less exposed to the AGN radiation field and possess greater inertia compared to the ionized gas. Thus, the observed radial stratification reflects the survival and entrainment conditions of each gas phase rather than a temporal evolutionary pathway (Veilleux et al. 2005, 2020). This structure is consistent with the multiphase nature of AGN-driven outflows: the ultrafast nuclear wind (T ∼ 10^7 − 8 K) shocks the surrounding ISM, generating hot gas that then cools and produces ionized phases (T ∼ 10⁴ K), traced by [OIII] and Hα. Cold molecular gas (i.e., H₂, as traced by CO, T ∼ 10–100 K) can arise either through radiative cooling of the shocked gas or via the entrainment of preexisting molecular material.

In an energy-conserving outflow scenario, theoretical models predict that, in general, the velocity of different gas phases should scale inversely with their density, with the fast nuclear wind transferring energy to progressively cooler and denser material as it shocks and expands into the ISM. Accordingly, molecular gas is often expected to move more slowly than the ionized phase, although faster velocities are possible if it is efficiently entrained by the outflow. However, in IRAS17, we found that the (average) velocity derived from the molecular CO emission is 1280 ± 480 km s⁻¹ (Longinotti et al. 2023), which is comparable to the velocities measured in the ionized gas traced by [OIII] of the fast outflow (i.e., 3c component) or Hα. This result suggests that the fast AGN wind may efficiently couple to the cold molecular gas, entrain it, or that projection effects are masking intrinsic differences between the phases. In an in situ formation scenario, the molecular gas would inherit the velocities of the shocked ionized gas from which it condensed. These observations support a scenario in which different gas phases are dynamically coupled, leading to similar outflow velocities despite their distinct physical and thermodynamic conditions.

Nowadays, most studies focus on the X-ray and molecular gas phases, while the ionized gas counterpart is often unavailable. This component is crucial for understanding the physical mechanisms driving the outflows. So far, several sources, including extreme quasars, appear to follow a momentum-conserving regime (Bonanomi et al. 2023; Baldini et al. 2024; Zanchettin et al. 2021, 2023, see also Marasco et al. 2020 and references therein), indicating that even powerful quasars can host momentum-driven outflows, in contrast to predictions from theoretical models (e.g., King & Pounds 2015; Zubovas et al. 2022). In particular, Marasco et al. (2020) and Tozzi et al. (2021) investigated the connection between nuclear X-ray UFOs and the ionized phase (Hα and [OIII]) of large-scale outflows in bright nearby and high-z quasars, respectively. They found that the ionized outflows of such quasars follow a momentum-conserving regime. Travascio et al. (2024) analyzed MUSE/VLT data of the radio-quiet quasar PDS 456 (z ∼ 0.2), a local analog of powerful Cosmic Noon quasars (see also Bischetti et al. 2019), which hosts an X-ray UFO and a clumpy molecular outflow extending to ∼5 kpc. The extended [OIII] ionized outflow shows a momentum rate consistent with an energy-conserving expansion. However, when the multiphase outflow (ionized + molecular) is considered, the total momentum boost indicates that the outflow is not fully in the energy-conserving regime.

On the other hand, Esposito et al. (2024), studied the nearby galaxy, NGC 5506, classified as NLSy1 and hosting an X-ray UFO, as in this work. They found that the molecular and ionized outflows (traced with ALMA and MEGARA/GTC data, respectively) follow a momentum-conserving regime. However, the ionized outflow mass was derived only using [OIII], rather than Hα, which would have allowed for a direct comparison with the [OIII] tracer.

Outflows found in nearby galaxies hosting an UFO in the energy-conserving regime are rare objects (i.e., IRAS11119+3257 in Tombesi et al. 2015, Mrk 231 in Feruglio et al. 2015, IRAS05189+2524 in Lutz et al. 2020, and the object studied in this work), and their nature warrants further investigation. Understanding their occurrence and physical conditions is essential to constrain feedback models and their impact on galaxy evolution.

6.3. Identifying the ionization mechanisms in the outflow: “Positive feedback” in action?

In this subsection, we present the application of different diagnostic diagrams, such as the Hα equivalent width versus [NII]/Hα diagram (WHAN; Cid Fernandes et al. 2011), the Hα equivalent width versus Hα velocity dispersion diagram (WHaD; Sánchez et al. 2024), and the Baldwin, Phillips & Terlevich (BPT) diagrams (Baldwin et al. 1981), to identify the dominant ionization mechanisms in different regions of the galaxy, such as the disk and the outflow.

6.3.1. WHAN diagram

To better understand the origin of the ionization source in the inner region of our target, we constructed the WHAN diagram (Cid Fernandes et al. 2011), which compares the equivalent width of the Hα line, EW(Hα), with the [NII]/Hα flux ratio. This diagnostic helps disentangle the contributions from star formation and AGN activity, differentiating Seyfert (strong AGN) from LINER (weak AGN) regions.

According to this diagram, galaxies can be classified into different categories¹⁴ based on their log([NII]/Hα) and EW(Hα) values: star-forming galaxies, weak and strong AGNs, and retired (passive) galaxies. Using our spatially resolved data, we derived WHAN diagram for each individual kinematic components (i.e., primary, secondary and outflow; see Peralta de Arriba et al. 2023), as shown in Figs. 16 (top panel) and 17. We clearly identified that the primary component is dominated by ionization from a strong AGN. The secondary component exhibits a significant contribution from star formation in the inner regions, while AGN activity becomes more prominent in the outer parts. This star-forming component (light green area in Fig. 16), overlaps with the butterfly region of enhanced velocity dispersion observed in the secondary Hα component (see Fig. 7). The tertiary component, which traces the outflow, is predominantly dominated by star formation throughout most of its spatial extent, with AGN-dominated ionization confined to a small region toward the south.

Fig. 16. Top: 2D maps corresponding to the WHAN regions according to Cid Fernandes et al. (2011) classification: star-forming galaxies (SF, in light green), weak AGN (wAGN, in cyan), strong AGN (sAGN, in dark blue), and passive and retired galaxies (in gray). Bottom: 2D maps corresponding to the WHaD regions according to Sánchez et al. (2024) classification: star-forming galaxies (SF, in blue), weak AGNs (wAGNs, in orange), strong AGNs (sAGNs, in red), and retired and unclassified galaxies (in black and gray).

Fig. 17. WHAN diagram for the primary, secondary and outflow components. The continuum emission has been corrected for the AGN fraction (i.e., fSF = 0.57) as described in Salomé et al. (2023) to take into account only the contribution from SF. The main results are not affected if this correction is not considered.

Deriving the continuum contribution associated with the different components to constrain the EW(Hα) parameter is challenging in IRAS17, as both star formation and AGN activity contribute to its evolution. To account for this, we corrected the total continuum emission to isolate only the star formation contribution. In the work by Salomé et al. (2023) using the CIGALE code, they found that the AGN contributes approximately 43% to the total continuum emission, leaving 57% attributable to star formation. This fraction was used to compute the EW(Hα) values shown in Fig. 17. If the AGN contribution had not been corrected for, the resulting EW(Hα) would represent a slightly lower value and thus a lower limit; however, this correction would not alter our conclusions. All these results suggest that the ionized outflow may be enhancing star formation by compressing the gas, thereby triggering positive feedback in IRAS17. However, based on our data, we cannot distinguish whether the observed star formation occurs within the outflow itself or is instead due to radiation from stars located in the disk plane (e.g., Maiolino et al. 2017; Gallagher et al. 2019).

Furthermore, we also note that the WHAN diagram cannot reliably distinguish between AGN excitation and shock-ionized gas, as both mechanisms can produce similarly high [NII]/Hα ratios. Additional diagnostics, such as BPT line ratios (i.e., Baldwin et al. 1981), could help in disentangling these contributions.

6.3.2. WHaD diagram

We also explored the WHaD diagram presented by Sánchez et al. (2024), which is used to distinguish between different ionization and kinematic regimes based on the equivalent width and velocity dispersion of the Hα line (i.e., EW(Hα)-σ(Hα)). Regions with high EW(Hα) and low σ(Hα) are typically associated with star formation due to young-massive OB-stars, tracing recent star forming activity, where the ionized gas is primarily influenced by thermal motions and turbulence driven by stellar processes. In contrast, areas with low EW(Hα) and high σ(Hα) often indicate AGN-driven outflows or (high-velocity) shock-ionized gas, where the kinematics are dominated by nonthermal processes and the contribution from young stars is minimal. Intermediate values may correspond to composite regions or partially resolved zones where multiple ionization sources coexist, i.e., “retired galaxies” (RGs), whose ionization is due to hot old low-mass evolved stars, such as post-AGB stars. Our results, shown in Figs. 16 (bottom panel) and 18, suggest that the Hα kinematics in IRAS17 are too extreme to be explained by star formation alone. Consequently, most of the FoV covered by the individual components falls into the AGN-dominated region of the diagram. This outcome is not in conflict with our previous findings, as both diagnostics are tracing different physical processes. On the one hand, the ionization indicators (e.g., the [NII]/Hα ratio) confirm the presence of star formation in the outflow, consistent with the enhanced star formation rate (SFR) derived from molecular gas (Salomé et al. 2023). On the other hand, the Hα kinematics reveal that pure star formation cannot account for the observed velocity dispersions, pointing to the additional influence of AGN activity or shocks. As stated in Sánchez et al. (2024), high-velocity shocks can also reproduce the typical values found in sources hosting an AGN. Overall, these findings indicate a complex interplay between AGN-driven shocks and star formation processes within the same regions, highlighting the role of feedback in shaping the evolution of the host galaxy.

Fig. 18. WHaD diagrams for the primary, secondary, and outflow components. The continuum emission has been corrected for the fSF = 0.57 from Salomé et al. (2023) to take into account only the contribution from SF. The main results are not affected if this correction is not considered.

6.3.3. BPT diagram(s) for the outflow

Using the [NII]-BPT diagrams (e.g., Baldwin et al. 1981; Kewley et al. 2001; Kauffmann et al. 2003; Schawinski et al. 2007) to investigate the excitation mechanism in the outflow component (i.e., the tertiary component for all the lines involved; Fig. 19), we find that the outflowing gas spans all three regions defined by the [NII]/Hα and [OIII]/Hβ ratios, corresponding to HII-like, composite, and Seyfert (Sy) emission. This indicates that the ionization close to the AGN is dominated by star formation, transitions through composite regions, and becomes AGN-dominated toward the southeast. When including fast radiative shocks modeled with the MAPPINGS III code (Allen et al. 2008), we observe significant overlap with all three regions, particularly in the Seyfert-dominated regime (Fig. 19, top). These results are consistent with those derived from the WHAN diagram and help better distinguish the ionization mechanisms near the AGN.

Fig. 19. Top: [NII]-BPT diagrams for the outflow component (i.e., tertiary component for all the lines involved). Shock+precursor grid models from Allen et al. (2008) are shown assuming Z = Z⊙ and different values of the shock velocities (v = 200 − 1000 km s−1) and magnetic field (B = 10−4–10 μG). The solid curve defines the theoretical upper bound for pure star formation from Kewley et al. (2001), while the dashed one in the [NII]-diagram is the Kauffmann et al. (2003) empirical classification. The dash-dotted line represents the demarcation line between Seyfert galaxies and shocks or LINERs from Schawinski et al. (2007). Bottom: 2D map derived applying the [NII]-BPT criteria for the outflowing component. The black cross identified the AGN position.

The observed ionization structure, with HII-like emission in the nucleus and a transition to AGN-dominated emission toward the southeast, is consistent with the extinction map. Low A_V in the center indicates a relatively dust-free line of sight, while higher A_V in the surrounding ring traces dense, dusty gas. These clouds can locally shield the compact AGN radiation, allowing stellar photoionization, which is more spatially distributed, to dominate in some regions, including the nucleus. However, some spaxels with elevated A_V still show AGN-dominated emission, reflecting anisotropic AGN illumination and the fact that A_V measures dust along our line of sight rather than the shielding between the AGN and the gas.

The presence of broad outflowing components and the overlap with fast shock models in the BPT diagrams indicate that the outflow likely drives shocks, which contribute alongside AGN and stellar photoionization to the observed line ratios and ionization structure. Overall, the ionization structure arises from the interplay of gas density, dust distribution, local shielding, shocks, and AGN geometry.

The BPT diagrams based on [OI] and [SII], often considered more reliable than [NII]-based BPTs (see discussion in Gallagher et al. 2019), could only be fit with a 1c model due to lower S/N in the emission. Nevertheless, they confirm the presence of SF-dominated ionization in the inner regions, supporting our interpretation.

Together, these findings indicate that the outflow simultaneously hosts star-forming regions and AGN-ionized gas, and may also drive shocks. This highlights the combined action of positive feedback (stimulating star formation) and energetic feedback (including gas ionization, shocks, and turbulence) in shaping the ionization and kinematics of the host galaxy.

6.4. Energy budget and feedback mechanisms of the ionized outflow

The kinetic coupling efficiency, ε, of the ionized outflow is a useful parameter to assess whether the AGN can transfer a significant fraction of its energy to the gas. Similarly, the momentum rate of the outflow, Ṗ_out, provides insight into the driving regime of the outflow and allows us to evaluate whether the different outflow phases (i.e., molecular and ionized) are dynamically linked.

6.4.1. Kinetic coupling efficiency and energy transfer

As is presented in Fiore et al. (2017), an important parameter essential for understanding the strength and impact of AGN feedback is the kinetic coupling efficiency defined as the ratio between the kinetic power of the outflow and the bolometric luminosity of the source: ε = Ė_out/L_bol. This parameter tells us how efficiently the energy output from the central source (an AGN or starburst) is transferred into the kinetic energy of the outflowing gas. It quantifies the fraction of AGN luminosity that goes into driving the outflow. A high efficiency (ε ≥ 1%) means the AGN is very effective at transferring its energy into feedback on the host galaxy, identifying a strong coupling. The outflow can have the potential to regulate star formation and impact galaxy evolution: this means that a significant fraction of the AGN’s energy is transferred to the gas, heating it and/or expelling it. This process can disperse the cold molecular gas, heat it, or prevent its collapse, thereby suppressing star formation (quenching). On the other hand, a low ε (≲1%) identifies a weak coupling, where the outflow does not strongly impact host ISM.

In our case, we found that ε ranges between ∼3%–15% for the Hα line (using the formula by Rupke et al. 2005 for the v_out and considering the two electron densities), while showing a lower value for the tertiary component of the [OIII] line (i.e., from the LR setup, ε ∼ 1%–4%)¹⁵, a value high enough to enable efficient energy transfer from the outflow to the ISM, potentially impacting the host galaxy (i.e., quenching star formation for ε > 1%; see Harrison et al. 2018). Moreover, strong turbulence, as seen in the Hα maps, with σ > 120 km s⁻¹, can inject additional energy into the gas, preventing it from cooling and collapsing under gravity. These findings suggest a scenario in which AGN activity can produce both negative (suppressing star formation) and positive feedback (potentially triggering star formation in some regions; see Sect. 6.3.1), acting on different regions and timescales across the host galaxy.

6.4.2. Momentum transfer and outflow driving mechanisms

Another important parameter is the ratio between the momentum rate of the outflow and that supplied by the AGN (i.e., Ṗ_out/Ṗ_rad), which is a powerful diagnostic of the driving mechanism and efficiency of AGN feedback. The AGN radiative momentum rate is defined as Ṗ_rad = L_bol/c = 1.7 × 10³⁴ g cm s⁻², with L_bol = −5.2 × 10⁴⁴ erg s⁻¹ (Giroletti et al. 2017; Longinotti et al. 2018). Values close to unity (≲1) suggest momentum-conserving outflows, likely driven by direct radiation pressure, while significantly higher values (a factor ≳10) indicate energy-conserving outflows, in which an AGN-driven wind efficiently transfers energy to the ISM via shocks. Ratios well below unity may imply inefficient coupling or that other processes, such as stellar winds or supernovae, are contributing to the observed outflow. Thus, this diagnostic helps distinguish between different feedback regimes and assess the AGNs impact on galaxy evolution. In our case, we found Ṗ_ion/Ṗ_rad ∼ 13 (see Table 6 and Fig. 14), which suggests an energy-conserving regime, where the AGN-driven wind shocks the ISM, forming a hot, adiabatic expanding bubble that accelerates the surrounding gas. Such momentum boosts are consistent with theoretical models of energy-driven feedback and imply that the AGN is efficiently transferring energy to the ISM, potentially regulating star formation through either negative or positive feedback (King & Pounds 2015 and reference therein). In energy-conserving outflows, the kinetic energy from the AGN is efficiently transferred through a hot wind that drives the cooler gas phases. This results in a coordinated evolution of all phases, forming a dynamically coupled, multiphase outflow. In contrast, momentum-conserving outflows may exhibit weaker coupling and more independent evolution of each phase due to rapid energy dissipation. In our source, the comparison of the momentum rates among phases reveals that the ionized and molecular components exceed the momentum rate of the UFO by factors of ∼5–27 and ∼20–30, respectively (see Table 6). These large momentum boosts reinforce the picture of an energy-conserving outflow where the AGN’s energy is efficiently deposited into the larger-scale ISM. This implies that the different phases of the outflow are dynamically linked, supporting the idea of a sustained and powerful AGN feedback capable of shaping the host galaxy evolution.

6.5. Gas depletion and feedback timescales in the molecular and ionized outflows: Witnessing “negative” feedback?

We compared the depletion times and the mass loading factor of the ionized and molecular gas phases in IRAS17. The former parameter is defined as the ratio between the outflow mass and the mass outflow rate, ${\tau }_{\mathrm{depl}}=\frac{M_{\mathrm{out}}}{{\dot{M}}_{\mathrm{out}}}$, and represents the time it would take for the outflow alone to remove all the available gas (if no new gas is accreted or formed). A short depletion time implies that the outflow can rapidly deplete the galaxy’s gas reservoir. The latter parameter is the ratio between the mass outflow rate and the SFR, η = Ṁ_out/SFR, and measures the strength or efficiency of the outflow relative to the galaxy’s star formation activity. A value of η > 1 indicates that the galaxy is ejecting gas at a rate higher than it is forming stars, suggesting strong feedback capable of regulating or even quenching star formation. Our results may be affected by several factors, such as the assumed wind geometry and the luminosity-to-mass conversion factors, which could introduce uncertainties of up to an order of magnitude (see also Herrera-Camus et al. 2020).

As presented in Salomé et al. (2023), NLSy1 galaxies, including IRAS17, clearly show very short depletion times, t_dep, SF = M(H₂)/SFR < 0.1 Gyr, compared to typical star-forming disks (i.e., ∼2 Gyr; see Leroy et al. 2013). On one hand, they proposed a nuanced view of AGN feedback: AGNs can initially promote star formation, particularly through mechanical compression from jets or outflows. However, this effect may be temporary, as AGN activity eventually depletes the cold gas reservoir, leading to long-term quenching. These findings fit within a “two-phase feedback” framework, where positive and negative AGN feedback coexist but act on different timescales. On the other hand, AGN feedback may be weak, ineffective, or simply not yet acting. This challenges the traditional view that AGNs always strongly suppress star formation in their host galaxies. Finally, as argued by Cresci & Maiolino (2018) and others, AGN feedback may not act instantaneously: instead, there could be a delay between AGN ignition and the onset of feedback effects, particularly those on larger scales, such as gas removal or heating.

For what concerns the depletion time of the molecular outflow, we derived a depletion time, ${\tau }_{\mathrm{depl}}({\mathrm{H}}_2)=\frac{M_{\mathrm{out}}({\mathrm{H}}_2)}{{\dot{M}}_{\mathrm{out}}({\mathrm{H}}_2)}$, in IRAS17 of ∼0.7 Myr from Longinotti et al. (2023), based on a molecular gas mass of ∼10⁸ M_⊙ and a mass outflow rate of 139 M_⊙ yr⁻¹. This value is shorter than the dynamical timescales of the outflow, calculated as t_dyn = R_out/v_out = 2.8 kpc/1280 km s⁻¹ ∼ 2.1 Myr, suggesting that the molecular gas can be largely expelled before the outflow propagates through the entire observed region. For the ionized gas, we computed a depletion time of ∼0.61 Myr, using the Hα ionized gas mass of 0.4 × 10⁷ M_⊙ and an outflow mass rate of 6.4 M_⊙ yr⁻¹. Its dynamical timescales is comparable, t_dyn ∼ 0.66 Myr (i.e., from the ratio 0.98 kpc/1461 km s⁻¹). The comparable values of τ_depl and t_dyn indicate that the outflow can efficiently remove the ionized gas from the host region within its dynamical timescale.

The fact that depletion times are comparable to or shorter than the dynamical times indicates high efficiency in gas removal from the galaxy’s central regions. In the molecular phase, depletion occurs on timescales shorter than the gas travel time, indicating that dense, star-forming gas is preferentially expelled, favoring negative feedback. In the ionized phase, the close match between depletion and dynamical times reinforces that the outflow is effective across gas phases. Combined with ionization diagnostics showing coexisting star formation, these results support a dual-feedback scenario: the outflow both removes a substantial gas reservoir (negative feedback) and compresses pockets of gas to trigger star formation (positive feedback).

Furthermore, the presence of line diagnostics indicating star formation in parts of the ionized gas suggests that the molecular and ionized phases are not independent. They are likely co-spatial or dynamically coupled: ionized winds may entrain molecular clouds or compress them, launching molecular gas outward or triggering localized star formation. In this scenario, while the molecular phase is more vulnerable to depletion, the ionized phase acts as the primary conduit for AGN energy, coupling it into the multiphase outflow and regulating the balance between suppression and triggering of star formation.

The mass loading factor, η, derived using a SFR of 26 M_⊙ yr⁻¹ (see Salomé et al. 2021), highlights the dominant role of the cold gas. We derived η_mol ∼ 5.3 for the molecular gas, indicating that the outflow is removing gas at more than five times the SFR, while the ionized gas shows a lower η_ion ∼ 0.3, insufficient by itself to regulate star formation. These results suggest a scenario in which AGN feedback initially affects the ionized and/or atomic gas, which can subsequently cool and form molecular clouds. The later expulsion of this molecular gas by the outflow represents the primary negative feedback mechanism, effectively suppressing star formation in the host galaxy.

In agreement with our results, Herrera-Camus et al. (2020) derived the mass loading factors for different gas phases in the nearby ultraluminous infrared galaxy IRAS F08572+3915. They found that the molecular phase dominates the mass loading compared to the atomic and ionized phases, and is the only phase where the gas ejection rate exceeds the rate of molecular gas consumption, yielding η_mol > 1. Overall, in their source, the molecular gas is the dominant component driving the outflow and regulating star formation, whereas the ionized phase contributes less to both the mass budget and feedback efficiency, as similarly observed for IRAS17 in this study.

6.6. Can the outflow(s) escape the galaxy?

To assess whether the ionized outflow can effectively remove gas from the galaxy, we modeled the kinematics of the primary (systemic) Hα emission using kinemetry (e.g., Krajnović et al. 2006; see also Bellocchi et al. 2012). This approach allowed us to determine the disk inclination and derive the intrinsic circular velocity, v_circ, which in turn was used to estimate the escape velocity, v_esc, as a function of radius. By comparing the measured outflow velocities from Hα with v_esc, we evaluated whether a significant fraction of the ionized gas can reach or exceed the escape velocity, and thus revealed the capability of the AGN-driven wind to expel gas from the galaxy.

Kinemetry models the observed line-of-sight velocity field by expanding it along concentric elliptical annuli using a Fourier series. The zeroth-order term, A₀, represents the systemic velocity, while the first-order terms, A₁ and B₁, describe the bulk rotational motion. The amplitude of the first harmonic, $k_1=\sqrt{A_1^2+B_1^2}$, provides a direct measure of the observed rotational velocity, v_obs, projected along the line of sight. Following Bellocchi et al. (2012), we employed a ten-term harmonic expansion in our analysis and adopted a covering factor (i.e., COVER) of 0.75¹⁶. The geometry of each annulus is defined by its position angle, Γ, and flattening, q = b/a, where a and b are the semimajor and semiminor axes, respectively. Under the assumption of a thin disk, the flattening is related to the disk inclination as cos i ∼ q. The intrinsic circular velocity is then obtained by correcting for inclination as $v_{\mathrm{circ}}(R)=k_1(R)/sini=k_1(R)/\sqrt{1-q^2}$. The resulting inclination-corrected rotation curve is used to estimate the escape velocity of the galaxy, v_esc, assuming a spherically symmetric potential modeled as a singular isothermal truncated sphere (e.g., Bischetti et al. 2019; Veilleux et al. 2020; Marconcini et al. 2025): $v_{\mathrm{esc}}(R)\approx \sqrt{2}{v}_{\mathrm{circ}}(R)$. This provides a first-order assessment of whether the ionized gas motions remain gravitationally bound.

In our case, we considered the last five data points of the first harmonic amplitude, k₁, which traces the observed rotational velocity and gives an observed circular velocity of 201.4 km s⁻¹, and we then corrected it for the mean inclination derived using such values, from which we derived a mean flattening of q ∼ 0.5; see Fig. 20). This yields an intrinsic circular velocity of v_circ ∼ 268.5 km s⁻¹. Assuming a simple spherical potential (Veilleux et al. 2020), this corresponds to an escape velocity of v_esc ∼ 380 km s⁻¹. The outflow velocity traced by Hα reaches ∼1450 km s⁻¹, with a comparable value derived from [OIII]. These velocities are well above v_esc (nearly four times higher), suggesting that a substantial fraction of the ionized gas is likely unbound (v_out ≫ v_esc). These results demonstrate that the outflow is sufficiently powerful to remove gas from the galaxy’s potential.

Fig. 20. Observed rotational curve parameter, k1 (top), and the flattening, q, of each ellipse (bottom), extracted at different radii using kinemetry. Error bars are smaller than the symbol size.

7. Summary and main conclusions

We have presented new GTC/MEGARA optical IFU observations of the nearby NLSy1 galaxy IRAS 17020+4544, one of the few known sources exhibiting both an X-ray UFO (Longinotti et al. 2015) and a molecular outflow (Longinotti et al. 2023) consistent with an energy-conserving regime. This galaxy also shows the presence of an elongated VLBI radio jet (Giroletti et al. 2017) that may be inducing shocks in the ISM.

Using MEGARA/GTC data in both LR (R ∼ 6000) and MR (R ∼ 12 000) modes, we characterized the ionized outflow in the galaxy’s central region by analyzing the Hα and [OIII] emission lines. The ionized outflow is confined within the more extended molecular outflow. By deriving its kinematic and energetic properties (velocity, mass, kinetic power, and momentum), we assessed the energy budget of the outflow and compared it with that of the molecular and X-ray phases.

The main results of this work are summarized below:

• The Hα and [OIII] emission line profiles are well fit with up to three Gaussian components, primarily detected in the central regions where the ionized outflow is located. We identify three components: (1) A primary (or narrow) component, present across the entire line emitting region, identified with the systemic component, characterized by relatively low velocity dispersions (σ_mean ∼ 90 − 100 km s⁻¹) and a clear rotational pattern. Peculiar features in the Hα velocity dispersion map, such as an X-shaped enhancement, may indicate turbulence driven by the outflow or by the jet-ISM interactions. (2) A secondary “broad” component of the Hα line covers roughly the same area as the systemic component and exhibits complex kinematic structures, with high velocity dispersions (σ_mean ∼ 380 km s⁻¹). In [OIII], this component traces an ionized outflow with modest velocities within an intrinsic radius of R_out ∼ 0.7 kpc. (3) A tertiary “very broad” component identifies the most extreme ionized outflow, present in both Hα and [OIII]. In Hα, its extension is R_out = 0.94 ± 0.32 kpc, while in the (LR-V) [OIII] emission it is more compact (R_out ∼ 0.5 kpc).

• The fast ionized outflow traced by the Hα and [OIII] lines (tertiary component) and the slower ionized outflow detected only in [OIII] (secondary component) are confined within the extent of the molecular outflow observed in CO (R_CO = 2.8 ± 0.3 kpc). The [OIII] outflow associated with the secondary component is resolved in both MR-G and LR-V setups, while the faster outflow is only resolved in the LR-V setup. The presence of multiple [OIII] velocity components indicates a stratified ionized gas structure: fast emission originates from dense, shock-compressed gas near the wind-ISM interface, whereas slower emission arises from more diffuse, partially cooled gas at larger radii. This stratified structure was previously observed in IRAS17 in the X-ray band, confirming a multiphase, multicomponent wind.

• We confirm that the ionized gas phase follows an energy-conserving regime, as has previously been found for the molecular phase (Longinotti et al. 2023). Both the Hα and fast [OIII] outflows exhibit large momentum boosts, consistent with efficient energy deposition from the AGN into the large-scale ISM. These results suggest that different gas phases are dynamically linked, indicating powerful AGN feedback capable of significantly influencing the host galaxy’s evolution.

• Diagnostic diagrams, depletion timescales, and mass-loading factors consistently point to a dual-feedback scenario. On the one hand, ionizing radiation detected in the outflow and in the butterfly region, together with the enhanced SFR inferred from the molecular gas, indicates localized, outflow-triggered star formation (positive feedback). This suggests that the outflow compresses the surrounding gas, triggering small-scale episodes of star formation in IRAS 17020+4544. On the other hand, the extreme kinematics of the ionized gas, which cannot be explained by stellar processes alone, together with the short depletion times and high molecular mass-loading factor (η> 1), reveal efficient removal or heating of the cold gas reservoir (negative feedback). The high velocities and energetics of both the ionized and molecular outflows indicate that the outflow is primarily AGN-driven, with only a minor contribution from star formation. While positive feedback may occur locally, our results show that negative feedback dominates globally and on longer timescales, effectively suppressing star formation across the host galaxy. This dual behavior underscores the complex and multifaceted role of AGN-driven winds in regulating galaxy evolution.

• Comparing the Hα outflow velocity with the escape velocity derived from its primary component, we found that the ionized outflow is capable of leaving the galaxy, since v_out ≫ v_esc. This result supports a scenario in which AGN-driven winds can deplete the host galaxy’s ISM and regulate star formation.

Our results shed light on the complex feedback processes occurring in this extreme and peculiar source. To better understand the interplay between different gas phases and spatial scales (i.e., UFOs, ionized outflows, and molecular outflows), studies of a larger, more complete sample of NLSy1 galaxies hosting X-ray UFOs are needed. Such a sample will allow us to test whether these systems are consistent with theoretical models of AGN feedback (e.g., Zubovas & King 2012; King & Pounds 2015; Zubovas et al. 2022).

Data availability

Reduced data cubes are available at the CDS via https://cdsarc.cds.unistra.fr/viz-bin/cat/J/A+A/709/A260

Acknowledgments

We thank the referee for a constructive report that helped improve the paper. EB, AGP, CCT, ACM acknowledge support from the Spanish grants PID2022-138621NB-I00 and PID2021-123417OB-I00, funded by MCIN/AEI/10.13039/501100011033/FEDER, EU. ALL acknowledges support from DGAPA-PAPIIT IA103625. YK acknowledges support form UNAM-PAPIIT grant 102023. This research has made use of the NASA/IPAC Extragalactic Database (NED), which is funded by the National Aeronautics and Space Administration and operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with the National Aeronautics and Space Administration. This paper is based on observations made with the Gran Telescopio Canarias (GTC), installed in the Spanish Observatorio del Roque de los Muchachos of the Instituto de Astrofśica de Canarias, in the island of La Palma. This work is based on data obtained with MEGARA instrument, funded by European Regional Development Funds (ERDF), trough Programa Operativo Canarias FEDER 2014-2020.

References

Appendix A: Deriving the extinction map, A_V

We calculated the visual extinction, A_V, from the Balmer decrement Hα/Hβ, assuming the Calzetti et al. (2000) attenuation law with R_V = 4.05, typical of starburst galaxies. An intrinsic ratio of Hα/Hβ = 3.1 and a fixed electron temperature of T_e = 10⁴ K were adopted. Based on the interpretation of the kinematic components presented in Tables 3, 4, and 5, the extinction map was derived using the Hβ flux from the MR-G data. This configuration allows us to resolve multiple kinematic components, and the Hβ line is well covered within the observed wavelength range. In contrast, in the LR-V setup only a single component could be fit, as the line lies near the blue edge of the spectrum, making its kinematics less reliable.

We derived the extinction of each component using the integrated flux maps of Hα and Hβ for the primary (systemic), secondary and tertiary (outflow) components. We also computed A_V using the total flux of the three components and by combining the Hα and Hβ fluxes of the secondary and tertiary components only (i.e., A_V(secondary+tertiary); see the bottom-right panel in Fig. A.1). The results show moderately different extinction values across components: the primary and secondary components display relatively high extinctions (A_V ∼2.8 and 3.5 mag, respectively), while the tertiary (outflowing) component exhibits a significantly lower value (A_V = 1.24 mag). When combining the secondary and tertiary components or using the total flux, we obtained values comparable to that of the primary component (A_V = 2.73 mag).

Fig. A.1. Extinction maps, AV, derived from the Hα and Hβ emission-line maps. From left to right, top to bottom: AV obtained from the total, primary, secondary, tertiary, and combined (secondary and tertiary) components. The black cross marks the continuum peak emission. The combined map (i.e., AV secondary+tertiary components) was finally adopted for estimating the outflow visual extinction, AV.

The primary and secondary components present similar, relatively high extinctions, suggesting that both are located behind comparable dust columns within the host galaxy or circumnuclear region. In contrast, the tertiary (outflowing) component shows much lower extinction, possibly indicating that the outflow lies on the near side of the dust distribution or has partially cleared its line of sight. However, the lower A_V of the outflow may also reflect the reduced S/N of its Hβ flux and uncertainties in line decomposition, continuum placement, or aperture coverage. Consequently, this value should be regarded as an indicative lower limit rather than a precise measurement.

Studies on dust extinction in outflows versus galactic disks yield mixed results. Some works (e.g., Rose et al. 2018; Mingozzi et al. 2019; Fluetsch et al. 2021) found that outflows are less affected by extinction (as in our case), while others (e.g., Holt et al. 2011; Villar Martín et al. 2014; Rodríguez del Pino et al. 2019) reported higher extinction in the outflows. Indeed, in some systems, outflows may be intrinsically dustier than the disk. Fluetsch et al. (2021) found a correlation between the difference in extinction between the outflow and the disk and the ionized outflow gas mass: the most massive ionized outflows tend to have more dust than the surrounding disk. In such cases, the preferential expulsion of dusty gas may enhance the extinction in the outflow while simultaneously lowering it in the disk. This provides further evidence that dust plays a key role in galactic outflows, particularly in the most massive ones. In particular, this finding supports models in which galactic winds are driven by radiation pressure on dusty clouds (e.g., Fabian 2012; Ishibashi et al. 2018).

All A_V maps reveal a clear ring-like structure with high extinction (A_V ≳ 4 mag), likely tracing regions of enhanced dust content. The inner region of the A_V total map, however, shows extinction values similar to those of the outflowing (tertiary) component. We finally adopt A_V = 2.73 mag as a representative extinction value for the ionized outflow.

Although the stellar continuum was not subtracted from our spectra, we acknowledge that this could affect the derivation of the extinction A_V, particularly the flux of the Hβ line, which is more strongly affected by underlying stellar absorption than Hα. This leads to an overestimated Hα/Hβ ratio and, consequently, to an artificially high extinction correction (A_V). In our case, we verified that using the extinction derived from the tertiary (outflowing) component (A_V = 1.24 mag), instead of the higher value obtained by combining the secondary and tertiary flux emissions (A_V = 2.73 mag), would reduce the estimated outflow momentum rate ratio (i.e., Ṗ_out/Ṗ_rad) from ∼13-61 to ∼4-20, depending on the adopted electron density (n_e = 2500 and 500 cm⁻³, respectively).

However, given the uncertainties on our measured quantities, the resulting values remain consistent within the errors and do not affect the main physical interpretation or conclusions of our analysis.

Appendix B: Additional kinematic maps

In this appendix we present the kinematic maps derived from the one component (1c) fits to the [SII] and [OI] doublets observed in the LR-R setup. The results are shown in Fig. B.1.

Fig. B.1. Same as Fig. 7 but for the [SII]6716,6731 (top) and [OI]6300,6363 (bottom) doublets.

For the [SII] doublet, both the velocity field and the velocity dispersion maps (Fig. B.1, top) follow the same kinematic pattern observed for the molecular gas traced by CO(1-0) (Salomé et al. 2021). An enhanced velocity dispersion is found in the inner region, particularly along the rotational axis, and a clear rotational pattern is also visible in the velocity field. Some spaxels in the companion galaxy also show [SII] emission, suggesting the presence of high-ionization regions, possibly associated with star-forming activity or shocks. AGN activity appears unlikely, given the absence of [OIII] emission in that region.

For the [OI] doublet (Fig. B.1, bottom), its emission is detected only in the nuclear region, where it follows the same rotational pattern seen in the [SII] maps, with a velocity dispersion distribution similar to that of the [SII] doublet.

Appendix C: Revealing ionized gas emission from the companion galaxy

Salomé et al. (2021), using NOEMA data, found that IRAS17 is interacting with a small companion galaxy located at a similar redshift. This interaction had not been previously identified, as IRAS17 appears to be an undisturbed barred spiral galaxy in optical images (Ohta et al. 2007). It remains unclear whether the companion flux emission corresponds to the comma-like (tidal tail) feature visible in the NOT optical image, located just near the foreground star, as shown in Fig. 2.

With our optical data, it is difficult to disentangle the contribution of the companion from that of IRAS17 using the Hα line alone. However, guided by the CO(1-0) systemic flux emission reported by Salomé et al. (2021), we identified the region likely associated with the companion and performed a two-component fit to the Hα spectra. The companion shows a slightly higher redshift than IRAS17, with z ∼ 0.0618 as derived by our MEGARA/GTC data in the Hα emission. A secondary Hα component is detected toward the southwest, redshifted relative to the systemic component and exhibiting a low velocity dispersion (σ ≲ 100 km s⁻¹), possibly indicating a connection with IRAS17 and a low level of turbulence (see the discussion in Sect. 4.1).

The companion also exhibits emission in other lines, such as [OIII], Hβ (detected only in the LR-V setup), and [SII]. A tentative BPT classification based on the [SII]/Hα ratio suggests that the companion is dominated by star formation. Based on our kinematic results, we propose that the companion may have interacted with IRAS17 in the past, and that we are now witnessing the residual effects of this merger in the ionized gas phase, which does not appear to be strongly perturbed. A more detailed investigation of this source is beyond the scope of this work.

All Tables

All Figures

	Fig. 1. Observed integrated spectrum of the inner region (R ∼ 1.9 kpc ≈ 1.6″) obtained with the LR-R setup. Several emission lines are detected, including the prominent Hα-[NII] complex, the fainter [SII] and [OI] doublets, as well as very faint [FeVII], FeII, and HeI lines.
In the text

Fig. 2. Top: SDSS-r continuum image of IRAS17. The MEGARA-IFU FoV is indicated by a black rectangle. A foreground star (marked with a black star) is visible to the north along with the companion galaxy (indicated by a magenta square) detected by Salomé et al. (2021). The flux intensity map is shown in nanomaggy units (1 nanomaggy = 3.631 μJy). Bottom: Continuum image obtained with MEGARA-IFU within the wavelength range 6770–6840 Å. The flux intensity map is shown in janskys. North is up and east is to the left in both images. The continuum peak emission is marked with a black cross in both panels and corresponds to the AGN position, as supported by its coincidence with the optical AGN position measured by Gaia, the radio emission detected with e-MERLIN at 1.5 GHz (Longinotti et al. 2023), and the CO peak emission (Salomé et al. 2021).

In the text

Fig. 3. Nuclear integrated spectra (within the PSF region) of the Hα-[NII]λλ6548,6583 (top) and Hβ-[OIII]λλ4959,5007 complexes (bottom), fit with three kinematic components. For the Balmer lines, the green, magenta, and blue curves represent the systemic, ILR, and broad (outflowing) components, respectively. For the [OIII] lines, the green, purple, and blue curves represent the systemic, intermediate (slow outflow), and broad (fast outflow) components.

In the text

Fig. 4. FeII template subtraction from the observed Hβ-[OIII] complex at the position of the emission peak in the MR-G setup (within a single spaxel of 0.4″ × 0.4″). The observed spectrum is shown in light blue, the FeII template in orange, and the residuals (pure Hβ-[OIII] complex) in green. The vertical dashed red lines mark the iron emission lines in this wavelength range (i.e., rest-frame wavelengths: λFe II = 4923.92, 4993.35, and 5018.45 Å).

In the text

	Fig. 5. Hα flux map of the primary component (red contours), with overlaid 12CO(1-0) contours of the systemic emission (white) from Salomé et al. (2021). The molecular and ionized gas in IRAS17 and its companion galaxy, traced on the northwest side by the CO contours, exhibit a remarkable spatial overlap.
In the text

	Fig. 6. Hα velocity dispersion map of the primary component (white contours) with overlaid CO(1-0) contours of the molecular outflow (blue and red contours trace the approaching and receding molecular outflowing gas, respectively) from Longinotti et al. (2023).
In the text

Fig. 7. Kinematic maps for the ionized component as traced by the Hα line. From left to right, top to bottom: Flux intensity, velocity field and (intrinsic) velocity dispersion maps of the primary (systemic), secondary and tertiary components. The velocity maps of each component have been corrected for the systemic velocity, vsys, derived at the AGN position (i.e., at the continuum intensity peak, marked by a black cross in all the maps). Gray contours in the flux intensity systemic maps trace the R-band continuum emission, highlighting the position of the foreground star. The magenta circles in the primary and secondary maps mark two individual spaxels in the northeast (NE) and southwest (SW) directions within the butterfly region, from which the spectra shown in the bottom panel are extracted. The spectra are fit with two components, and the broad [NII] emission is detected at S/N > 4 (i.e., 5.5 in the NE and 4.1 in the SW). The light gray circle in the bottom-left panel is the MEGARA PSF (FWHM). For all the maps, north is up and east to the left.

In the text

	Fig. 8. Similar to Fig. 7 but for the [NII]6583 line. The main difference between the [NII] and Hα (Fig. 7) panels is that the inner region of the [NII] secondary component, corresponding to the ILR in the Hα map, is absent in the forbidden [NII] line.
In the text

	Fig. 9. Same as Fig. 7 but for the ionized component as traced by the Hβ line according to the MR-G setup.
In the text

	Fig. 10. Same as Fig. 7 but for the ionized component as traced by the [OIII]5007 line according to the MR-G setup.
In the text

	Fig. 11. Same as Fig. 7 but for the ionized component as traced by the Hβ (top) and [OIII]5007 (bottom) complex as derived in the LR-V setup.
In the text

	Fig. 12. Hα velocity dispersion map of the primary component, with the [OIII] velocity dispersion contours (white) of the primary component overlaid, obtained from the LR-V setup.
In the text

Fig. 13. Electron density (ne) maps in units of cm−3 derived using the [SII] (top) and log U bottom) methods. The [SII] method is applied to the total [SII]λλ6717,6731 line profile, while the log U map sums the flux contributions from the secondary and tertiary components. Typical mean and median values, together with the standard deviation and the median absolute deviation, respectively, are shown in the insets.

In the text

Fig. 14. Force of the molecular, ionized, and X-ray phases of the outflow, Ṗout, in IRAS17, normalized to the source’s X-ray momentum rate, ṖUFO (i.e., ṖUFO/Ṗrad = 1.87 ± 0.62; see Longinotti et al. 2018), is shown as a function of the outflow velocities. The results for IRAS17 derived from NOEMA, MEGARA (i.e., Hα using ne = 2500 cm−3), and XMM-Newton are shown with blue circles. To compare results across the three phases, we used the same formula for vout as in Longinotti et al. (2023, i.e., Δv + 2 × σB) and normalized their Ṁout values assuming C = 1. Filled circles in orange and magenta represent the results derived for 2c component of the [OIII] emission (i.e., slow outflow) in the LR and MR setups, respectively. Two values are plotted for each setup (connected by a solid line), corresponding to the two assumed electron densities, ne = 500 and 2500 cm−3 (see text for details). Empty circles, following the same color code as for the 2c component, indicate the 3c component of the [OIII] emission (i.e., fast outflow) in both setups. The fast [OIII] outflow is spatially resolved in the LR setup (empty orange circle), while it remains unresolved in the MR setup (empty magenta circles), for which we show the two results (lower limits) derived assuming the two ne values. Results for three other well-known sources for which the energy-conserving regime has been observed, Mrk 231, IRAS 05189+2524, and IRAS11119+3257 (all ultraluminous infrared galaxies, shown using black, green, and red squares, respectively), are also included. Their data are extracted from Feruglio et al. (2015), Lutz et al. (2020), and Tombesi et al. (2015), respectively (see also Marasco et al. 2020 and references therein). The dashed blue and red lines represent the two linear fits marking the energy-conserving outflow regime: the former obtained using CO, Hα, and X-ray data, and the latter using these three points plus the [OIII] result from the LR-V setup (empty orange circle). The dotted gray line indicates the prediction for “momentum-conserving” outflows (i.e., Ṗout/ṖUFO ∼ 1).

In the text

Fig. 15. Velocity field shift map, Δv (in km s−1), between the outflow and primary components. These maps are derived from the Hα tertiary component (left) and from the [OIII] line for the tertiary (resolved) component (right) using LR-V setup. A flux-emission threshold of 10% of the peak value has been applied to the outflow kinematic maps. The gray circle marks the PSF size at 10% of the flux peak (PSF10%; see text for details). The blue and red contours trace the approaching and receding molecular outflowing gas, respectively, as presented by Longinotti et al. (2023, see their Fig. 4). The ionized-gas emission nicely matches the region of the molecular outflow and shows blueshifted velocities (Δv < 0 km s−1; the minimum value is reported in each panel). The intrinsic outflow size for each component is also indicated in the corresponding panel.

In the text

Fig. 16. Top: 2D maps corresponding to the WHAN regions according to Cid Fernandes et al. (2011) classification: star-forming galaxies (SF, in light green), weak AGN (wAGN, in cyan), strong AGN (sAGN, in dark blue), and passive and retired galaxies (in gray). Bottom: 2D maps corresponding to the WHaD regions according to Sánchez et al. (2024) classification: star-forming galaxies (SF, in blue), weak AGNs (wAGNs, in orange), strong AGNs (sAGNs, in red), and retired and unclassified galaxies (in black and gray).

In the text

	Fig. 17. WHAN diagram for the primary, secondary and outflow components. The continuum emission has been corrected for the AGN fraction (i.e., fSF = 0.57) as described in Salomé et al. (2023) to take into account only the contribution from SF. The main results are not affected if this correction is not considered.
In the text

	Fig. 18. WHaD diagrams for the primary, secondary, and outflow components. The continuum emission has been corrected for the fSF = 0.57 from Salomé et al. (2023) to take into account only the contribution from SF. The main results are not affected if this correction is not considered.
In the text

Fig. 19. Top: [NII]-BPT diagrams for the outflow component (i.e., tertiary component for all the lines involved). Shock+precursor grid models from Allen et al. (2008) are shown assuming Z = Z⊙ and different values of the shock velocities (v = 200 − 1000 km s−1) and magnetic field (B = 10−4–10 μG). The solid curve defines the theoretical upper bound for pure star formation from Kewley et al. (2001), while the dashed one in the [NII]-diagram is the Kauffmann et al. (2003) empirical classification. The dash-dotted line represents the demarcation line between Seyfert galaxies and shocks or LINERs from Schawinski et al. (2007). Bottom: 2D map derived applying the [NII]-BPT criteria for the outflowing component. The black cross identified the AGN position.

In the text

	Fig. 20. Observed rotational curve parameter, k1 (top), and the flattening, q, of each ellipse (bottom), extracted at different radii using kinemetry. Error bars are smaller than the symbol size.
In the text

	Fig. A.1. Extinction maps, AV, derived from the Hα and Hβ emission-line maps. From left to right, top to bottom: AV obtained from the total, primary, secondary, tertiary, and combined (secondary and tertiary) components. The black cross marks the continuum peak emission. The combined map (i.e., AV secondary+tertiary components) was finally adopted for estimating the outflow visual extinction, AV.
In the text

	Fig. B.1. Same as Fig. 7 but for the [SII]6716,6731 (top) and [OI]6300,6363 (bottom) doublets.
In the text