© The Authors 2023

1. Introduction

Galaxy formation and evolution is a complex process involving several different physical phenomena acting simultaneously on different physical and temporal scales. Gas is a key player in this picture, feeding star formation and the accretion onto the central supermassive black hole (SMBH), and is in turn affected by feedback mechanisms. The feedback can manifest through powerful winds that are able to blow away the gas from the centre of the galaxy, quenching star formation and starving the BH of fuel (Veilleux et al. 2020). Active galactic nucleus (AGN) feedback processes play a fundamental role in galaxy growth and evolution; they are thought to be at the origin of the M_BH − σ_⋆ relation (King 2010; Silk & Rees 1998) and to prevent the overgrowth of massive galaxies (Bower et al. 2012; Hopkins et al. 2014).

In the hard X-ray spectrum, hot (T ∼ 10⁶ − 10⁷ K) ultra-fast outflows (UFOs) have been observed in ∼40% of the bright nearby local AGN population (Tombesi et al. 2010, 2012; Gofford et al. 2013, 2015). These winds, developed from the AGN accretion disk (≤1 pc), are observed through the detection of blueshifted (velocities up to v ∼ 0.3c) absorption lines associated with highly ionised iron transitions in the hard X-ray spectrum (Reeves et al. 2003; Tombesi et al. 2010).

Massive galaxy-scale cold (T_kin ∼ 10 − 100 K) molecular outflows with velocities between hundreds of and a few thousand km s⁻¹ have been observed in the last decade (Feruglio et al. 2010; Fischer et al. 2010). These winds can be detected by P Cygni profiles of the OH molecule in the far-IR regime (Fischer et al. 2010; Feruglio et al. 2010; Sturm et al. 2011; Veilleux et al. 2013) as well as blue- and redshifted high-velocity wings in the CO, HCN, or HCO⁺ profiles using interferometric observations in the millimetre band (Aalto et al. 2012; Cicone et al. 2012, 2014).

The theoretical model that is usually invoked to explain large-scale outflows launched by AGNs is the blast-wave scenario (Silk & Rees 1998; King 2010; Faucher-Giguère & Quataert 2012). According to this model, a nuclear wind arises from the accretion disk of an AGN and impacts on the interstellar medium (ISM) producing a forward shock and a reverse shock. The forward shock propagates through the unperturbed ISM producing a large-scale outflow. This outflow could be either energy- or momentum-driven, depending on whether cooling of the reverse shock is efficient. If it is, the energy is conserved and outflow propagates adiabatically (energy-driven), showing a momentum boost with respect to the X-ray wind (${\dot{p}}_{\mathrm{out}}^{\mathrm{ISM}}/{\dot{p}}_{\mathrm{out}}^{\mathrm{X}-\mathrm{ray}}\gg 1$; Faucher-Giguère & Quataert 2012). Otherwise, the energy is radiated away and only momentum is transferred to the ISM (${\dot{p}}_{\mathrm{out}}^{\mathrm{ISM}}/{\dot{p}}_{\mathrm{out}}^{\mathrm{X}-\mathrm{ray}}\sim 1$; King 2010).

Which model is most favoured by observations is a highly debated question. Simultaneous observations of X-ray winds and large-scale outflows are needed to test their predictions. The momentum rate versus wind velocity diagram is a widely used tool to visualise and compare the properties of different outflows (Feruglio et al. 2015; Tombesi et al. 2015). Smith et al. (2019) recently summarised the momentum rate versus the wind velocity for a sample of ten objects with observed X-ray UFOs and large-scale galactic outflows (see also Fig. 16, this work). Some of these sources, such as the luminous quasar PDS 456 (Bischetti et al. 2019) and the ultra-luminous infrared galaxy (ULIRG¹) IRAS F11119+3257 (Tombesi et al. 2015; Veilleux et al. 2017), seem to favour the momentum-driven scenario, while other objects, such as the ULIRG Mrk 231 (Feruglio et al. 2015) and the Seyfert 1 galaxy IRAS 17020+4544 (Longinotti et al. 2015, 2018), show large-scale outflows whose momentum rate is boosted compared to the X-ray wind. Finally, other sources, such as the multiple-lensed quasar SDSS J1353+1138 (Tozzi et al. 2021), do not appear to favour either of the AGN feedback models. Overall, as clearly drawn by Smith et al. (2019), the picture is much more complex than expected from AGN blast-wave feedback models.

Furthermore, most sources in the sample explored by Smith et al. (2019) are ULIRGs. These objects have an intense star formation, whose contribution to feedback processes is hard to distinguish from the AGN contribution. Testing the prediction of the blast-wave scenario in galaxies with a more moderate star formation activity is necessary to overcome this issue. The work by Sirressi et al. (2019) on the local Seyfert 2 galaxy MCG-03-58-007, with a star formation rate (SFR) ∼20 M_⊙ yr⁻¹, was a first step towards this direction. These authors detected a compact H₂ component that, if interpreted as an outflow, would present a momentum rate equal to ∼40% of that of the X-ray UFO. Our study on the LIRG and galaxy pair IRAS 05054+1718 (also known as CGCG 468-002) also fits into this context. The main target of this work is the western component of the pair (hereafter NED01), a local LIRG hosting a Seyfert 1.9 nucleus. The source shows a moderate SFR of 5–10 M_⊙ yr⁻¹ (De Looze et al. 2014; Pereira-Santaella et al. 2015) and hosts a powerful X-ray wind (Ballo et al. 2015), being a suitable candidate to test the AGN feedback scenario reducing the possible contamination from star formation-driven outflows.

Our aim is to investigate the presence of a large-scale molecular outflow in NED01 and in the companion NED02, by studying the distribution and the kinematics of the molecular gas. The latter is the phase of the ISM that is most tightly connected to star formation (Wong & Blitz 2002; Bigiel et al. 2008; Leroy et al. 2008), as stars form primarily in molecular clouds (e.g. Lada et al. 2010; André et al. 2014). We use carbon monoxide (CO) as a tracer of cold molecular hydrogen (H₂), because it is the second most abundant molecule after H₂, and its rotational transitions are easily observable at submillimetre–millimetre wavelengths. The lowest J levels of CO can be easily excited at molecular cloud temperatures (T ∼ 10 K, Omont 2007) and so the CO(1–0) and CO(2–1) lines can be used to estimate the total molecular gas mass in galaxies. Since these lines are optically thick at typical molecular cloud conditions, their luminosity is not proportional to the H₂ gas column density, but a CO(1–0)-to-H₂ conversion factor (α_CO) needs to be assumed. The estimate of this value is not straightforward as it depends on the physical state of the gas and needs to be calibrated using multiple molecular line tracers, which are often difficult to detect in extragalactic sources. For the Milky Way ISM and for normal star-forming galaxies, a value of α_CO = 4.3 M_⊙ (K km s⁻¹ pc²) is widely accepted (Bolatto et al. 2013). For different ISM environments, such as the massive molecular outflows discovered in local starbursts and AGN host galaxies (see review by Veilleux et al. 2020), the α_CO parameter is very poorly constrained. In this work we assume that molecular outflows have an α_CO = 2.1 ± 1.2 M_⊙ (K km s⁻¹ pc²), which is the value measured by Cicone et al. (2018) on the molecular outflow of the well-studied local ULIRG NGC 6240. We use α_CO = 4.3 ± 1.3 M_⊙ (K km s⁻¹ pc²), a value typically assumed when treating the molecular ISM of isolated galaxies like the Milky Way (Bolatto et al. 2013), to evaluate the molecular mass of the galaxy disk.

This paper is organised as follows. The selected targets are presented in Sect. 2. In Sect. 3 we describe the new high-sensitivity Atacama Large Millimetre/submillimetre Array (ALMA) CO(1–0) observations used in this work. The analysis performed on the data is reported in Sect. 4. In Sect. 5 we model the kinematics of the CO(1–0) emission using the 3D-Based Analysis of Rotating Object via Line Observations (3D-BAROLO) software for the two targets. In Sects. 6 and 7 we present the observations and the analysis of the new X-ray datasets, and in Sect. 8 we derive the energetics of the X-ray wind. The interpretation of the results is discussed in Sect. 9, where we test different hypothesis in the AGN-driven feedback scenario. In Sect. 10 we summarise our results and conclusions.

Throughout the paper we adopt a standard ΛCDM cosmological model with H₀ = 67.8 km s⁻¹ Mpc⁻¹, Ω_Λ = 0.692, and Ω_M = 0.308 (Planck Collaboration XIII 2016). At the distance of IRAS 05054+1718 NED01 (z = 0.0178, revised in this work), the physical scale is 0.373 kpc arcsec⁻¹.

2. Target description

The western and eastern pair members, in this work indicated respectively as NED01 and NED02 (see also Pereira-Santaella et al. 2015), have a projected distance of ∼29.6″ ∼ 11 kpc. Figure 1 shows the ALMA CO(1–0) contours overlayed onto a g-band image from the Pan-STARR Survey 1 (Chambers et al. 2016). According to Stierwalt et al. (2013), the system is in an early merger stage after a first encounter between the two galaxies since their disks are still symmetric, but show signs of tidal tales.

Fig. 1. ALMA CO(1–0) map (black contours) overlayed onto the g-band optical image from the Pan-STARR Survey 1. The ALMA CO(1–0) emission was averaged over a spectral range corresponding to v = [ − 490, +230] km s−1 with respect to the systemic redshift of NED01, which includes the CO(1–0) emission from both members of the galaxy pair. The contours correspond to the (3, 6, 9, 12, 24, 50) × σRMS levels, with σRMS = 0.2 mJy beam−1 being the average rms of the ALMA CO(1–0) map (not corrected for the primary beam).

The western galaxy NED01, at z = 0.0178 ± 0.0004², hosts a Seyfert 1.9 nucleus, and it is classified as a LIRG (log(L_{IR(8 − 1000 μm)}/L_⊙) = 10.6, Pereira-Santaella et al. 2015). Because of the presence of the AGN, its SFR is not well constrained in the literature, with values ranging between 5 M_⊙ yr⁻¹ (De Looze et al. 2014; Pereira-Santaella et al. 2015) and ∼10 M_⊙ yr⁻¹ (Howell et al. 2010). Based on the ratio of the SFR to the BH accretion rate (log(SFR/ṁ_BH) ∼ 2), obtained from the [Ne II]15.56 μm and [O III]λ5007 gas velocity dispersion, the stellar velocity dispersion, and the 8–1000 μm IR-luminosity, Alonso-Herrero et al. (2013) suggested that NED01 is transitioning from a H II-dominated to a Seyfert-dominated LIRG.

NED01 represents an interesting case study for the effects of AGN feedback on galaxies. Ballo et al. (2015) detected a deep absorption trough at E ∼ 7.8 keV (2.1σ significance) in its Swift-XRT (X-ray telescope) spectrum, which has been interpreted as a highly ionised (log ξ ∼ 3 erg cm⁻² s⁻¹), high column density (N_H ∼ 10²³ cm⁻²), and ultra-fast (v_out = (0.11 ± 0.03)c) disk wind.

The companion NED02 is also a LIRG (L_IR = 10¹¹ L_⊙, Pereira-Santaella et al. 2015) and has a measured redshift of z = 0.016812 ± 0.000003³. NED02 was classified as a composite galaxy according to the BPT classification by Pereira-Santaella et al. (2015), but no evidence for AGN emission has been detected to date (see e.g. Alonso-Herrero et al. 2012). The SFR estimates for this galaxy range between SFR_{(1 − 10) Myr} ∼ 15 M_⊙ yr⁻¹ (Pereira-Santaella et al. 2015) and SFR ∼ 20 M_⊙ yr⁻¹ (Howell et al. 2010; De Looze et al. 2014).

3. ALMA CO(1–0) observations

The ALMA Band 3 (84.0–116.0 GHz) observations of IRAS 05054+1718 were carried out in Cycle 5 (Project code: 2016.1.00694.S, PI: P. Severgnini). The primary target was NED01 (corresponding to the phase centre of the interferometric dataset), but the field of view and spectral bandwidth of the data also cover the CO(1–0) emission from the companion NED02, and so we include the latter in our analysis. We only use data from the two scheduling blocks that have passed the ALMA data quality assurance (QA0), which are also the only datasets that were delivered to the PI, with observing dates 5 and 6 March 2017. According to the QA0 report, the total observing time including overheads for the two combined valid execution blocks was 100 min, and the total time on target was 62 min. The 40 ALMA 12 m antennas were arranged in the most compact configuration (C40-1), with baselines ranging from 14 m to 310 m. The precipitable water vapour (PWV) varied from 3 mm to 8mm, wind speed was 3.3 m s⁻¹, and humidity ∼50%. The quasar J0423-0120 was used for flux calibration, J0510+1800 was instead used for band-pass response, phase calibration, and pointing, and both sources in addition to the main target were used for atmospheric calibration and radiometric phase correction.

We employed four spectral windows, two for each side band of the ALMA correlator. Two adjacent high-resolution, 1.875 GHz wide spectral windows (960 channels each, channel width of 1953.13 kHz, corresponding to 5.2 km s⁻¹) were centred at sky frequencies of 113.179 GHz and 111.438 GHz in order to sample both the CO(1–0) line and the N = 1 spin-doublet transition of CN, which have rest-frame frequencies of ${\nu }_{\mathrm{CO}(1-0)}^{\mathrm{rest}}=115.2712$ GHz and ${\nu }_{\mathrm{CN}(1-0)}^{\mathrm{rest}}\simeq 113.4910$ GHz⁴. Two additional low-resolution 2 GHz wide spectral windows (128 channels, 15.625 MHz channel width, corresponding to ∼50 km s⁻¹) were centred at sky frequencies of 101.190 GHz and 99.387 GHz to probe the 3 mm continuum. In this work, we focus on the CO(1–0) line data, and we postpone the analysis of the CN(1–0) line to a future publication. Through the spectral line modelling described in Sect. 4.1, we found the CO(1–0) lines of NED01 and NED02 to be respectively centred at (sky) frequencies of ν = 113.2578 GHz and ν = 113.3651 GHz, which we used to compute new estimates of the systemic redshift of the two galaxies. Except for Sect. 4.1, where we worked with the initial datacubes not corrected for the right redshift, the rest of the analysis presented in this paper was performed on two separate datacubes (one for NED01 and one for NED02) corrected for their new CO-based redshift estimates.

The data were calibrated by running the version 5.4.0 of the Common Astronomy Software Applications (CASA) package calibration pipeline (McMullin et al. 2007). For the cleaning and other analysis steps we used CASA software version 5.6.1-8. We combined the measurement sets of the two execution blocks by using the task concat, after having pre-selected with split the CO(1–0) line spectral windows relevant to our target. An analysis of the continuum at 100 GHz, conducted using the two line-free spectral windows (see further details in Appendix A), shows a clear detection of both NED01 and NED02, with respective continuum peak flux densities equal to $S_{\mathrm{cont}}^{\mathrm{NED}01}=0.71\pm 0.03$ mJy beam⁻¹ and $S_{\mathrm{cont}}^{\mathrm{NED}02}=2.26\pm 0.05$ mJy beam⁻¹. For this reason, before proceeding with the analysis of the CO kinematics, we subtracted the continuum from the CO(1–0) spectral windows in the uv visibility plane using the task uvcontsub. We selected a zeroth-order polynomial fit to the continuum channels adjacent to the CO(1–0) line in the 112.26–113.0 GHz and 113.57–114.12 GHz sky frequency ranges (corresponding to v ∈ [ − 2650, −690] km s⁻¹ and v ∈ [820, 2280] km s⁻¹ with respect to the CO(1–0) line centre). We then worked exclusively on the continuum-subtracted CO(1–0) line data.

The cleaning procedure for modelling the true sky brightness distribution of the source out of the uv visibility data was performed using the task tclean, by selecting the automasking algorithm (auto-multithresh parameter), which creates a different mask for every channel, minimising negative sidelobes. We used Briggs weighting with the robust parameter set equal to zero, and a cell size of 0.2″. The synthesised beam size of the resulting cleaned datacube changes slightly with spectral channels, with a median value of 2.59″ × 2.25″, corresponding to an average spatial resolution of 0.97 kpc × 0.84 kpc. We adopted the native spectral resolution of 5.17 km s⁻¹ for the channel size. We selected a cleaning threshold equal to our first estimate of the average line rms of 2.5 mJy beam⁻¹ per channel. In order to account for the bias of the primary beam (PB) pattern on the image, we divided the cleaned datacube by the PB response using the task impbcor.

The mean rms CO(1–0) line sensitivity for source NED01 as measured in the cleaned and PB-corrected cube is 2.2 mJy beam⁻¹ per 5.2 km s⁻¹ channel. This value was calculated with the task imstat, by selecting the central 40″ portion of the field of view, and so it is adequate to characterise the noise fluctuations of the CO(1–0) data for NED01. We also verified that the noise follows a Gaussian distribution, so we adopted the mean rms value across the whole CO(1–0) spectral range on which our analysis is focused. For the companion NED02 instead, given its proximity to the edge of the PB’s FWHM, the CO(1–0) line sensitivity is lower, with an average 1σ rms value of 4.0 mJy beam⁻¹ per 5.2 km s⁻¹ channel. For this source, the 3D-BAROLO kinematic analysis (presented in Sect. 5.2) will be conducted on a portion of the datacube centred on NED02 and not corrected for the PB. Table 1 summarises the main observational parameters.

4. Analysis of the CO(1–0) line emission

4.1. CO-based redshift estimates

4.1.1. NED01

The continuum-subtracted CO(1–0) spectrum of NED01, reported in the left panels of Fig. 2, shows a double-peaked emission line centred at ν = 113.26 GHz. This is slightly different from the (redshifted) CO(1–0) central frequency of 113.189 GHz that was expected from a previous redshift estimate of this source (z = 0.0184), derived from the heliocentric velocity (error-weighted average of the optical and radio velocities) reported in the HyperLeda catalogue (see Makarov et al. 2014). To refine the redshift estimate of NED01 we performed two spectral fits, one with a single-Gaussian function and one with two Gaussians, both displayed in the top left panel of Fig. 2, with results listed in Table 2. We computed a new systemic redshift using the central frequency value derived from the single-Gaussian fit (ν = 113.2578 ± 0.0009 GHz). However, the source shows a double-peaked profile with clear asymmetries. To take this into account, we conservatively assigned to such CO-based redshift an uncertainty equal to the average frequency difference between the two CO line peaks (as measured from the double-Gaussian fit) and the single-Gaussian fit peak frequency value. We obtained z = 0.0178 ± 0.0004, which is offset by v ∼ −210 km s⁻¹ with respect to the previously known redshift. We further checked that the new redshift estimate matches with the kinematic centre of NED01’s host galaxy disk. The CO(1–0) position-velocity (PV) diagrams displayed in Fig. 3 were computed in CASA from slit-like apertures with sizes of ∼30″ and ∼9″, aligned with the major and minor axes of the CO(1–0) disk of NED01 respectively (see kinematic modelling reported in Sect 5). The new redshift estimate (obtained through the spectral analysis described above) is shown as a black cross at the centre of the two PV diagrams, hence confirming the correspondence with the rotational centre of the CO(1–0) disk.

Fig. 2. ALMA CO(1–0) continuum-subtracted spectra of the interacting galaxy pair IRAS 05054+1718: NED01 (left panels) and NED02 (right panels). The top and bottom panels display the same spectral data, but with different units on the x-axes. In particular, the top panels report the CO(1–0) flux density as a function of observed (sky) frequency, and were used to refine the systemic redshift estimates for the two galaxies (see Sect. 4.1). The bottom panels show the spectra corrected for redshift, where the CO(1–0) flux density is reported as a function of optical velocity along the line of sight. The CO(1–0) line spectrum of NED01 was extracted from an elliptical aperture maximising the CO(1–0) flux, with size 23″ × 14″, centred at RA = 05h08m19.858s, Dec = +17°21′45.898″. The CO(1–0) spectrum of NED02 was extracted from an 10″ × 12″ elliptical aperture centred at RA = 05h08m21.212s, Dec = +17°22′08.660″. The best-fit spectral parameters are listed in Table 2.

Fig. 3. CO(1–0) PV diagrams of NED01, extracted from a datacube where velocities are calculated with respect to the previously known redshift of the source (z = 0.0184). The PV diagrams confirm that our new redshift estimate of z = 0.0178 ± 0.0004, indicated with a black cross, closely matches the kinematic centre of the CO(1–0) source. The upper panel shows the PV diagram obtained from a slit-like aperture along the axis connecting the blue and red peaks of the CO(1–0) emission, with a size of ∼30″ and a position angle of 95° (measured anti-clockwise from the north direction). The PV diagram shown in the bottom panel was computed from a slit-like aperture with a size of ∼9″ orthogonal to the previous one.

4.1.2. NED02

A similar analysis aimed at refining the systemic redshift of the host galaxy was performed on the companion source, NED02. The total continuum-subtracted CO(1–0) spectrum of NED02, displayed in the right panel of Fig. 2, presents a single peak. We modelled it using a single-Gaussian function, whose best-fit parameters are reported in Table 2. The previous redshift estimate for this source, z = 0.016842, derived from the heliocentric velocity reported in the 2MASS Redshift Survey (Huchra et al. 2012), would have produced a CO(1–0) emission line peaked at 113.1385 GHz. Instead, our spectral analysis of the new ALMA CO(1–0) observations of NED02 shows that the CO(1–0) line is centred at a higher frequency of 113.3651 ± 0.0003 GHz. We used this value and its associated uncertainty to compute a new systemic redshift of NED02 equal to z = 0.016812 ± 0.000003, which is offset by ∼50 km s⁻¹ with respect to the previously known redshift.

4.2. CO(1–0) morphology

The redshift-corrected CO(1–0) emission line spectra of NED01 and NED02, plotted as a function of line-of-sight velocity, are shown in the bottom panels of Fig. 2. The corresponding best-fit spectral parameters are listed in Table 2. By using the CO(1–0) line fluxes obtained from the spectral fits, we computed the total CO(1–0) line luminosities of the two galaxies, and from these estimated their molecular line masses, by adopting a standard CO-to-H₂ conversion factor of α_CO = 4.3 ± 1.3 M_⊙ K km s⁻¹ pc². These values are reported in Table 3.

Figure 4 displays the zeroth, first, and second moment maps of the source NED01, computed by applying an intensity threshold of 0.01 Jy beam⁻¹ and within a box 30″ in size centred on the AGN. The velocity map shows clearly that the CO(1–0) gas kinematics in NED01 is dominated by ordered rotation. The blueshifted emission arises from the western side of the galaxy, and the redshifted emission from the eastern side, following the typical pattern of a rotating disk. Therefore, from Fig. 4 we can already infer that if a molecular outflow is present, it does not seem to impact the bulk of the CO(1–0) kinematics in this source. This was also evident from the CO(1–0) spectrum of NED01, which presents a clear double peak and no evidence for very high-velocity wings, which are typical of extreme molecular outflows detected in some local (U)LIRGs (see e.g. Cicone et al. 2014, and several references in Veilleux et al. 2020).

Fig. 4. Moment maps of NED01 computed within a box of size 30″ and including pixels above a CO(1–0) intensity threshold of 0.01 Jy beam−1 and within v ∈ [ − 500, 500] km s−1. Left: Intensity (moment 0) map, with contours plotted at [0.1, 0.3, 1, 3, 5, 8] Jy beam−1 km s−1. Centre: Velocity (moment 1) map, with contours plotted every 30 km s−1; the letters indicate the three CO(1–0) components discussed in the main text, where A indicates the main rotating CO(1–0) disk. Right: Velocity dispersion (moment 2) map, with contours plotted at intervals of 30 km s−1.

In addition to the central molecular disk (component A), the CO moment maps of NED01 reveal two apparently disconnected CO(1–0) emitting structures in this galaxy, labelled B and C in Fig. 4. Component B is ∼3.3 kpc and component C is ∼2.3 kpc from the centre of NED01’s main disk. Despite their apparent offset from the disk in the moment maps shown in Fig. 4, which results from the adoption of a sensitivity threshold, these two components are physically linked to the central disk, and also follow the same velocity pattern. This will be confirmed by the kinematic modelling with the BBarolo software presented in the next section. The CO(1–0) spectral line profiles of components B and C are single-peaked and narrow. Component B is redshifted, centred at a velocity of v ≃ 150 km s⁻¹, with a velocity dispersion of σ_v ≃ 13 km s⁻¹, and entrains a CO(1–0) flux of 0.65 Jy km s⁻¹, which, if adopting the same α_CO as above, corresponds to ∼4 × 10⁷ M_⊙ of molecular hydrogen gas. Component C is only slightly redshifted (v ∼ 9 km s⁻¹), and can be modelled with a single-Gaussian with σ_v ≃ 30 km s⁻¹ and an integrated CO(1–0) flux of 2 Jy km s⁻¹, corresponding to 1.3 × 10⁸ M_⊙. Figure 1 shows that components B and C of the CO(1–0) emission from NED01 (detected respectively at 3 and 10σ in the CO(1–0) channel map shown as contours in Fig. 1) do not correspond to any significant sub-structure in the optical continuum; instead, they overlap with diffuse lower surface brightness stellar light. We can rule out the hypothesis that the CO components B and C are not distinguishable in the g-band Pan-STARR image because of high dust extinction, since in this case we would expect to detect them in the ALMA 3 mm continuum map, which is not the case (see Appendix A). We therefore suggest that the CO(1–0) components B and C are ISM substructures of the main galaxy disk, possibly tracing a spiral arm.

The companion galaxy NED02 is located close to the edge of the PB, and so the moment maps, reported in Fig. 5, were computed from the datacube not corrected for the PB. The CO(1–0) emission in this edge-on galaxy (see Fig. 1) appears more compact than in NED01. The spider diagram shows a velocity gradient skewed towards blueshifted velocities, tracing a disturbed rotating disk. The modelling presented in the next section supports this interpretation.

Fig. 5. Moment maps of NED02 computed from the non-PB-corrected datacube within a box of size 14″ centred at RA(ICRS) = 05:08:21.108; Dec(ICRS) = 17.22.08.489, including pixels above a CO(1–0) intensity threshold of 0.009 Jy beam−1 and within v ∈ [ − 500, 500] km s−1. Left: Intensity (moment 0) map, with contours plotted at [0.1, 1, 5, 8, 11] Jy beam−1 km s−1. Centre: Velocity (moment 1) map, with contours plotted every 30 km s−1. Right: Velocity dispersion (moment 2) map, with contours plotted at intervals of 20 km s−1.

5. Modelling of the CO(1–0) kinematics

In this section we analyse the gas kinematics in IRAS 05054+1718. As shown by Fig. 4, the bulk of the CO(1–0) emission from NED01 traces a rotating molecular disk, while the disk of NED02 appears more disturbed (Fig. 5). In order to study the presence of CO(1–0) components that are not participating in the rotation and may trace the effect of AGN feedback on the large-scale ISM, we first model the disk kinematics, and then study any residual emission, similar to Sirressi et al. (2019). For completeness, we perform the same analysis on the companion galaxy NED02, even though we do not have any evidence for the presence of an AGN in this source. A uniform and common analysis of the CO(1–0) gas kinematics of both members of the galaxy pair can give us insights into the role of galaxy interactions in shaping the cold ISM kinematics.

5.1. Modelling of the CO(1–0) disk in NED01

We model the disk rotation in NED01 using the 3D-Based Analysis of Rotating Object via Line Observations (3D-BAROLO, Di Teodoro & Fraternali 2015, also known as BBarolo). BBarolo identifies the set of geometrical and kinematic parameters that best fit the rotating gaseous disk observations, and uses these parameters to produce a mock datacube of the best-fit model. An important assumption of the BBarolo model, to consider when interpreting the results, is the hypothesis that the line emission is distributed in a geometrically thin disk whose kinematics is dominated by pure rotational motion. In other words, the model does not consider the presence of additional outflow components, tidal tails, or other features that do not belong to the main rotating disk. The disk model is built up by combining concentric rings with a user-defined width, and the comparison between models and data is performed ring by ring.

As shown by Fig. 4, the CO(1–0) line emission from NED01 has a complex morphology. In order to study the effect of adding the extended CO(1–0) line components in the BBarolo modelling, we run BBarolo on three regions with different sizes, and compare the results. The regions are as follows (see Fig. B.1): the Small region, which includes only component A, with size = 10″ × 8″ (3.7 × 3.0 kpc); the Medium region, which includes components A and C, with size = 16″ × 15″ (6.0 × 5.6 kpc); and the big region, which includes all three components, with size = 35″ × 35″ (13 × 13 kpc).

We produced the CO(1–0) datacubes that need to be given as an input to the software by cropping the continuum-subtracted ALMA CO(1–0) (clean) datacube with the CASA task imsubimage and then exporting it into a FITS file using the exportfits task. We selected a spectral range corresponding to v = [ − 600, +600] km s⁻¹. For each region, we also produced a residual map by subtracting the BBarolo disk model from the input datacube.

As parameters of the BBarolo model we kept the inclination (ι, i.e. the angle of the disk with respect to the line of sight) and the position angle (PA) of the molecular gas disk fixed at all radii. Using the ALMA moment 1 map shown in Fig. B.1, we estimated PA = 95° ±10° and ι = 45° ±15°. The inclination was inferred from the ratio of the minor to the major axis of the molecular disk. The fitting was performed considering only pixels with a S/N > 2.3 to enable inclusion of the low-S/N extended features.

In the following we present the results obtained by using the Big region, which includes all three CO(1–0) emission line components (A, B, and C) observed in the moment 1 map of NED01, while the results from the other two regions are reported in Appendix B. We focus on the Big region because, once it has been subtracted from the ALMA datacube, this is the fit that produces the lowest residual CO(1–0) emission. Additionally, the analysis of the Big region fit allows us to study the global CO(1–0) emission from NED01 as this region includes also components B and C. We therefore believe that the fit on the Big region is the most conservative one for our goal, which is to study residual emission whose kinematics is inconsistent with a rotating disk.

We constructed the BBarolo model disk by using seven rings of width 2″ (0.7 kpc), covering in total a region of 28″ in diameter (10.4 kpc) around NED01. The results of this fit are shown in Figures 6–8. Although the rotating disk clearly provides a very good fit for the bulk of CO(1–0) emission in and around NED01, there are statistically significant residuals. We note that the moment maps and the PV diagrams produced by BBarolo (Figs. 6 and 7) are not optimised to visualise residual emission, because they use a colour scale that is fine-tuned to enhance the match between data and model. The bottom rows of the channel maps in Fig. 8, displaying the residual CO(1–0) emission, show that positive (green contours) and negative (grey contours) residuals with similar significance (∼3σ) are present in each channel map with velocity |v|> 200 km s⁻¹, consistent with low-level residual rotation that is not accounted for by the fit. However, at blueshifted velocities, within a range −290 < v < −200 km s⁻¹, we detected only positive residuals at much higher significance (S/N ≥ 15), showing a compact morphology.

Fig. 6. Moment maps showing the comparison between the data and the BBarolo best fit performed on the Big region around NED01. Shown (from left to right) are data, BBarolo model, and residual emission. The rows show (rom top to bottom) the moment 0 map (intensity), moment 1 map (line of sight velocity), and moment 2 map (velocity dispersion). The black cross indicates the rotating disk centre found by the fit, and the dashed black line shows the major axis of the disk model.

Fig. 7. Position velocity diagrams showing the comparison between the data and the BBarolo best-fit model performed on the Big region around NED01. The red contours show the disk model, and they are overplotted onto the observed ALMA CO(1–0) line data (displayed in grey, with blue contours at 3, 6, 9, 15, 20, 30, 40, and 50σ). The upper panel shows the PV diagram extracted along the major axis. The yellow data points indicate the best-fit rotational velocity of each ring. The lower panel shows the PV diagram computed along the minor axis.

Fig. 8. CO channel maps computed by integrating the CO(1–0) emission in channels of Δv = 42 km s−1, showing the comparison between data and BBarolo best-fit disk model within the Big region around NED01. In each of the three panels, the top row (blue contours) displays the data, the central row the model (red contours), and the bottom row the residual emission (green contours). The contours are plotted at 3, 6, 9, 15, 20, 30, 40, 50σ. Symmetric negative contours are plotted in grey. The yellow cross indicates the centre of the disk model.

Figure 9 displays the CO(1–0) residual emission in NED01 after subtracting the best-fit disk model computed by BBarolo. The left panel of Fig. 9 shows a contour map of such compact residual structure, obtained integrating the CO(1–0) residual flux within the velocity ranges where it is most prominent, that is −290 ≲ v[km s⁻¹] ≲ − 200. This structure can be enclosed within a box of 8.4″ × 8.4″ (∼3 kpc, see black box overlayed on the map). We find that a similar residual, compact, and blueshifted structure is common to all three fits performed with BBarolo, and it is actually minimised in this one performed on the Big region compared to those using the smaller regions. The spectrum reported in the central panel of Fig. 9, extracted from the black squared aperture in the left panel, shows that the blueshifted CO(1–0) residual emission is spectrally resolved into two emission peaks, which we modelled with double-Gaussian functions (see best-fit parameters in Table 4). There is also a negative feature at the systemic velocity, common to all three BBarolo fits, which is probably an indication of overfitting of the rotational structure. By mapping separately the two blueshifted spectral components, we find that they trace two apparently independent, spatially offset structures, as shown in the right panel of Fig. 9. Here we plotted with different colours the first component of this residual CO(1–0) emission, which is the most blueshifted one integrated between v ∈ [ − 400, −200] km s⁻¹ (blue contours), and the second component, which was integrated between v ∈ [ − 200, 0] km s⁻¹ and displayed using the black contours. The first component peaks closer to the dynamic centre of the molecular disk, at a distance of R = 0.7 ± 0.2 kpc, while the second component peaks further away, at a distance of R = 1.5 ± 0.2 kpc. Table 4 lists for each of the two components their spectral best-fit Gaussian parameters and corresponding CO(1–0) flux, CO(1–0) luminosity, molecular gas mass, and their physical distance from the best-fit BBarolo disk model.

Fig. 9. Residual CO(1–0) emission in NED01 after subtracting the best-fit disk model computed by BBarolo in the Big region. Left: Map of the blueshifted CO(1–0) residual emission integrated in the range −290 ≲ v[km s−1] ≲ − 200. The blue contours are plotted at 3, 6, 9, 15, and 20σ. The black cross indicates the centre of the best-fit molecular disk. The blue cross indicates the peak of the CO(1–0) residual emission. The residual structure can be enclosed within a box of 8.4″ × 8.4″ (∼3 kpc, see black box overlayed on the map). Centre: Spectrum of the CO(1–0) residual emission extracted from an aperture corresponding to the black box shown in the left panel. The two CO(1–0) peaks are fit with double-Gaussian functions whose parameters are reported in Table 4. Right: Map of the two blueshifted spectral components visible in the central panel, plotted with different colours. The first component (blue contours) is the most blueshifted spectral peaks, integrated within v ∈ [ − 400, −200] km s−1; the second component (in black) corresponds to the secondary peak, integrated between v ∈ [ − 200, 0] km s−1. Contours are plotted at 3, 6, 9, 15, and 20σ.

5.2. Modelling of the CO(1–0) disk in NED02

Compared to NED01, the CO(1–0) moment maps of NED02 in Fig. 5 show a more disturbed, higher inclination, molecular disk. The CO(1–0) emission detected by ALMA in this source is also more compact in size, and hence its morphology and kinematics are more affected by beam smearing effects. We used BBarolo to model this CO(1–0) disk, following a procedure similar to NED01. However, given the absence of significant extended features, we did not deem it necessary for NED02 to run BBarolo on regions of different sizes. The input 3D dataset for the BBarolo modelling was obtained by cropping the cleaned (non-PB-corrected) continuum-subtracted ALMA datacube around NED02 using a box size of 27.4″ × 25.4″ (i.e. 9.4 × 8.7 kpc), and selecting a velocity range of v ∈ [ − 1000, +1000] km s⁻¹, computed with respect to the systemic redshift of NED02. From the ALMA observations, we estimated a position angle of PA = 120° ±10° and an inclination of ι = 70° ±15°. We applied a S/N threshold of S/N > 3 for each pixel included in the modelling. We used five rings, each of 2″ width (0.7 kpc).

Figures 10, 11, and B.2 show the moment maps, PV diagrams, and channel maps of NED02 produced by the fitting code, comparing the data with the best-fit model of a rotating disk. Although BBarolo can reproduce the global rotation pattern, we detect significant residuals, especially at blueshifted velocities, similar to what was found for NED01. The map of such residual CO(1–0) emission around NED02, integrated over the spectral range v ∈ [ − 170, 20] km s⁻¹, is shown in the left panel of Fig. 12. The box shown on the map enclosing the 3σ level CO(1–0) contours has a size of 8″ × 8″ (∼3 × 3 kpc). The spectrum of the CO(1–0) residual emission extracted from this region (central panel of Fig. 12) shows two peaks. The two spectral features peak at different positions (right panel of Fig. 12), where the CO(1–0) residual emission components integrated between v ∈ (− 350, −120) km s⁻¹ (blue contours) and v ∈ (− 120, 0) km s⁻¹ (black contours) are plotted separately. The two components do not overlap on the map, suggesting they trace physically distinct structures. Table 5 lists the best-fit spectral parameters of these two CO(1–0) residual features, their corresponding fluxes, luminosities, molecular gas mass estimates, and their distance from the CO(1–0) rotation centre of NED02.

Fig. 10. CO(1–0) moment maps of the companion galaxy NED02 showing the comparison between the data and the best-fit molecular disk model found by BBarolo.

Fig. 11. Position velocity diagrams of NED02 showing the comparison between the data (blue contours, plotted at 3, 6, 9, 15, 20, 30, 40, and 50σ) and the BBarolo best-fit rotating disk model (red contours). The upper panel shows the PV diagram extracted along the major axis, and the lower panel shows the one extracted along the minor axis of rotation.

Fig. 12. Residual CO(1–0) emission in NED02 after subtracting the best-fit rotating disk model computed by BBarolo. Left: Map of the CO(1–0) residual emission integrated within v ∈ [ − 170, 20] km s−1. Contours are plotted at 3, 6, 9, 15, and 20σ. The black cross indicates the dynamic centre of the rotating disk and the blue cross shows the centroid of the residual emission. The residual structure can be enclosed within a box of 8″ × 8″ (∼3 kpc, see black box overlayed on the map). Centre: Spectrum of the CO(1–0) residual emission in NED02, extracted from the squared aperture shown in the left panel. The results of the spectral fitting with Gaussian functions are given in Table 5. Right: Map of the two blueshifted spectral components visible in the residual spectrum shown in the central panel. The first component, integrated within v ∈ [ − 350, −120] km s−1, is displayed using blue contours. The second component, integrated within v ∈ [ − 120, 0] km s−1 is shown in black. Contours are plotted at 3, 6, 9, 15, and 20σ.

6. X-ray observations

6.1. Observations and data reduction

IRAS 05054+1718 was observed three times in the X-ray band. In 2012 it was observed by Nuclear Spectroscopic Telescope Array (NuSTAR, Harrison et al. 2013) for a total of ∼16 ksec; in 2014 it was the target of a monitoring programme (totalling ∼72 ksec) performed with the Neil Gehrels Swift Observatory (hereafter Swift, Gehrels et al. 2004); and more recently we obtained a simultaneous deep observation with XMM-Newton and NuSTAR (PI: V. Braito, 2021). The results of the Swift observations were published by Ballo et al. (2015), while the 2012 NuSTAR data are reported in Ricci et al. 2017. According to these works, the X-ray spectrum of NED01 is that of a moderately absorbed Seyfert 2 galaxy (with photon index Γ ∼ 1.7, column density N_H ∼ (1 − 2) × 10²² cm⁻², and flux in the 2–10 keV regime F_{2 − 10 keV} = (0.6 − 1)×10⁻¹¹ erg cm⁻² s⁻¹) with the NuSTAR 3–10 keV flux in 2012 being a factor 1.6 brighter than the Swift value. Ballo et al. (2015) first reported the presence of an absorption trough at E ∼ 7.8 keV (at 2.1σ), which was interpreted as the presence of highly ionised wind, outflowing with a velocity of v_out ∼ 0.1 c.

In this work, we re-analyse all X-ray spectra collected for NED01; we re-reduced the XMM-Newton and NuSTAR data, while for the Swift monitoring we consider the spectrum extracted by Ballo et al. (2015). In Table 6 we report the summary of the XMM-Newton and NuSTAR observations, while the list of the Swift observations can be found in Table 1 of Ballo et al. (2015). Although the dominant X-ray source is NED01, the inspection of the soft (0.3–1.5 keV) XMM images revealed a faint X-ray source at the expected location of the NED02 as well as possible extended emission. In particular, for NED02 we detected ∼200 X-ray counts in the 0.3–1.5 keV band with XMM. Assuming that they are due to thermal emission (kT = 0.6 ± 0.2 keV), as generally seen in star-forming galaxies, we derived F_{(0.5 − 2) keV} ∼ 2 × 10⁻¹⁴ erg cm⁻² s⁻¹ and L_{(0.5 − 2) keV} ∼ 2 × 10⁴⁰ erg s⁻¹. At higher energies, we cannot assess if there is any X-ray emission from NED02 that could be indicative of a weak or highly obscured AGN because it falls within the wings of the PSF of the bright companion.

6.2. XMM-Newton

We processed and cleaned the XMM-Newton data using the Science Analysis Software (SAS ver. 18.0.0) and the resulting spectra were analysed using standard software packages (FTOOLS ver. 6.30.1, XSPEC ver. 12.11; Arnaud 1996). The XMM-Newton-EPIC instruments operated in full-frame mode and with the thin filter. We first filtered the EPIC data for high background, which only moderately affected the observations. The EPIC-pn, MOS1, and MOS2 source spectra were extracted using a circular region with a radius of 20″, while for the background we adopted two circular regions with a radius of 20″. The response matrices and the ancillary response files at the source position were generated using the SAS tasks arfgen and rmfgen and the latest calibration available. After checking for consistency we combined the spectra from each of the individual MOS detectors into a single spectrum. Both the pn and MOS spectra were then binned to at least 100 counts per bin.

6.3. NuSTAR

NuSTAR observed IRAS 05054+1718 twice; the second observation was coordinated with XMM-Newton, starting at the same time but extending for 100 ksec beyond the XMM-Newton observation (see Table 6). We reduced the data following the standard procedure using the HEASOFT task NUPIPELINE (version 0.4.9) of the NuSTAR Data Analysis Software (NUSTARDAS, ver. 1.8.0). We used the calibration files released with the CALDB version 20220426 and applied the standard screening criteria, where we filtered for the passages through the SAS setting the mode to ‘optimised’ in NUCALSAA. For each of the Focal Plane Modules (FPMA and FPMB) the source spectra were extracted from a circular region with a radius of 50″, while the background spectra were extracted from two circular regions with a 50″ radius located on the same detector. Since the second set of observations partially overlap with the XMM-Newton exposure, we also extracted light curves in the 7–10 keV and 10–30 keV band from the same regions using the NUPRODUCTS task to check whether we could use the averaged spectra when performing a joint fit with the XMM spectra. The FPMA and FPMB background-subtracted light curves were then combined into a single curve. The inspection of these light curves revealed that in the last 100 ksec, not covered by XMM-Newton, NED01 was slightly brighter; we thus created a good time intervals (GTI) file corresponding to the strictly simultaneous part of the observation. This GTI file was then used to extract source and background spectra and the corresponding response files.

7. X-ray spectral analysis

The spectral fits were performed with XSPEC (ver. 12.11) and in all the models we included the galactic absorption in the direction of NED01 (N_H = 1.93 × 10²¹ cm⁻², HI4PI Collaboration 2016), which was modelled with the Tuebingen–Boulder absorption model (tbabs component in XSPEC, Wilms et al. 2000). The source spectra were binned to have at least 20 counts in each energy bin for the Swift data, and 100 counts for the EPIC-pn, the EPIC-MOS, and the NuSTAR spectra. For the portion of the NuSTAR spectra that is strictly simultaneous with XMM we adopted a lower binning of 50 counts. We employed χ² statistics and errors are quoted at the 90% confidence level for one interesting parameter. All the outflow velocities are relativistically corrected.

In Fig. 13 we show all the X-ray spectra collected for NED01, obtained by unfolding the data against a power-law model with photon index Γ = 2. It is clear that the Swift and the 2012 NuSTAR observations (black and red data points, respectively) caught NED01 in a similar bright state, while the most recent XMM-Newton and NuSTAR observations found the source in a much fainter state (blue and light blue data points, respectively). We therefore fit separately the 2012–2014 and 2021 epochs. It is also noticeable that the hard (E > 10 keV) X-ray emission measured in the first NuSTAR observation is at the same level as the averaged Swift-BAT spectrum from the 70-month survey⁵ suggesting the bright one as the most frequent state (grey data points in Fig. 13).

Fig. 13. Broad-band (rest frame 0.3–60 keV) X-ray spectra of all the observations of IRAS 05054+1718. The EFE spectra are obtained unfolding the data against a simple power-law model with Γ = 2. The 2012 NuSTAR observation is shown in red, the 2014 Swift is in black, while the spectra obtained with the coordinated XMM-Newton and NuSTAR observations are shown in dark and light blue, respectively. We also included the averaged Swift-BAT spectrum from the 70-month survey (Baumgartner et al. 2013). The residuals obtained from the joint fit for the Swift and NuSTAR spectra using a typical model for a moderately obscured Seyfert 2 as NED01 are shown in the inset. The model, composed of an absorbed primary power-law component, a scattered component, with the same photon index (Γ), and a reflected component, is able to reproduce the broad-band emission, but it leaves some residuals in absorption (see Sect. 7.1).

7.1. 2012 NuSTAR and 2014 Swift observations

We performed a joint fit for the Swift and NuSTAR spectra and considered the 0.3–10 keV and the 3.5–60 keV data for the Swift and FPM data, respectively. We first tested the baseline continuum model as reported in Ballo et al. (2015), which is a typical model for a moderately obscured Seyfert 2 as NED01. The model is composed of an absorbed primary power-law component, a scattered component, with the same Γ, and a reflected component (modelled with the XILLVER component in XSPEC). The reflected component represents the emission produced by a distant material, such as the putative torus. We tied only the photon index. We allowed the N_H of the neutral absorber, and the normalisation of the primary power-law and of the reflection components to differ.

This model provides a reasonable fit (χ² = 343.9/293 degrees of freedom) and requires a standard photon index (Γ = 1.81 ± 0.06). The column density of the neutral absorber is N_H,SW = (2.4 ± 0.2) × 10²² cm⁻² and N_H,NU = (4.4 ± 1.9) × 10²² cm⁻² for the Swift and the NuSTAR spectra, respectively. The NuSTAR 3–10 keV flux is a factor of 1.6 higher than the Swift value. Overall, this model is able to reproduce the broad-band emission, but it leaves some residuals in absorption that are shown in the inset of Fig. 13. We thus included in the model two Gaussian absorption lines, and the fit improved by Δχ²/Δν = 19/6. One feature is detected in the Swift spectrum at E = 7.8 ± 0.3 keV (Δχ² = 6) and the second at E = 8.6 ± 0.2 keV (Δχ² = 13) in the NuSTAR data. The equivalent widths (EWs) of the absorption troughs are 220 ± 150 eV and 190 ± 100 eV in the Swift and NuSTAR data, respectively.

To self-consistently account for these features, we replaced the two Gaussian absorption lines with a multiplicative grid of photoionised absorbers generated with the XSTAR photoionisation code (Kallman et al. 2004). Since the lines appear to be broad (σ_SW ∼ 0.2 keV, σ_NU = 0.3 ± 0.2 keV), we chose a grid that was generated with a high turbulence velocity (v_turb = 5000 km s⁻¹) and with column density in the range N_H = 10²¹ – 3 × 10²⁴ cm⁻² and ionisation⁶ log ξ = 1 − 6. We assumed that the ionised absorber has the same ionisation, but allowed the N_H and velocity to be different in the two epochs. The inclusion of the fast-outflowing and highly ionised (log ξ = 3.4 ± 0.5) absorber improves the fit by Δχ²/ν = 20/5. For the Swift observation we found $N_{\mathrm{H},\mathrm{SW}}=1.3_{-0.9}^{+1.7}\times {10}^{23}$ cm⁻² and $v_{\mathrm{out}}/c=-0.{12}_{-0.05}^{+0.02}$, while for the NuSTAR spectrum we found N_H, NU = 1.1 ± 0.6 × 10²³ cm⁻² and v_out/c = −0.23 ± 0.02. Our independent analysis confirms the results reported by Ballo et al. (2015), and also suggests that the wind in IRAS 05054+1718 could be variable or have multiple velocity components, as seen for other X-ray disk winds (e.g. PDS 456, Matzeu et al. 2017; Reeves et al. 2018; IRAS 13224-3809, Parker et al. 2018; MCG-03-58-007, Braito et al. 2022; HS0810+2554, Chartas et al. 2016; and IRAS 17020+4544, Longinotti et al. 2015). The best-fit model is given in Table 7, while the spectra with their respective best-fit models are shown in Fig. 14.

Fig. 14. Spectra and best-fit model (blue line) for the Swift (black data points) and 2012 NuSTAR (red) observations. The continuum model includes an absorbed power-law component and a reflected component (grey dotted lines). The reflected component (modelled with XILLVER) also includes the expected Fe Kα emission line at 6.4 keV. The NH of the neutral absorber, and the normalisation of the primary power law and of the reflected components are all allowed to vary. The model also includes an ionised outflowing absorber, for which the same ionisation is assumed, but the column densities and velocity are allowed to differ in the two spectra. The Swift data were rebinned for plotting purposes.

7.2. 2021 simultaneous XMM-Newton and NuSTAR observations

The 2021 XMM and NuSTAR observation caught NED01 in a faint state (see Fig. 13), where the observed 2–10 keV X-ray flux dropped by a factor of ∼4 with respect to the 2012 and 2014 observations. At this lower flux level the NuSTAR spectrum becomes background-dominated above 40 keV. We thus limited our spectral analysis to the 0.3–10 keV energy band for the XMM data, and the 3.5–40 keV band for the NuSTAR spectrum. Figure 13 also reveals that the average NuSTAR spectrum (light blue) lies well above the EPIC-pn spectrum (blue data points). If we apply the best-fit continuum model derived in the previous section, we find a cross-normalisation between XMM and NuSTAR of C = 1.36 ± 0.05, which is well above the current uncertainties in the cross normalisation between the two observatories. We therefore considered only the portion of the NuSTAR observation that is strictly simultaneous with XMM, where we found C = 1.17 ± 0.06. We note that the variation is mainly a flux change between the simultaneous spectrum and the last 100 ksec exposure, with no clear spectral variations. We applied the same baseline continuum model to the new observation, which describes well the overall spectral shape and already provides a reasonable fit (χ² = 255.7/205 d.o.f.) to the spectra with no evident residuals. The broad-band spectra and the best-fit models are presented in Fig. 15. The best-fit parameters are reported in Table 7, where we note that the faint state is caused by a drop in the primary emission and not by an increase in absorption, as is often found for obscured AGN. Here the primary emission drops by a factor of ∼4 with respect to the 2014 Swift spectrum and by an order of magnitude with respect to the 2012 NuSTAR observation.

Fig. 15. Spectra and best-fit model (red line) for the 2021 XMM (blue data points) and NuSTAR (light blue) observation. The continuum model includes an absorbed power-law component and a reflected component (grey dotted lines).

Although there are no residuals that can be associated with the disk wind observed in the past observations, we estimated the upper limits on the column density for a putative ionised absorber outflowing at v_out = 0.12 c or v_out = 0.23 c for the same ionisation parameter as found previously, of log ξ = 3.4. We found that they are of the order of ∼3 × 10²² cm⁻² and ∼4 × 10²² cm⁻² for the slower and the faster phase, respectively. The upper limits are only marginally below the 90% lower value of N_H estimated in the past observations.

8. The X-ray disk wind energetics of NED01

Here we compute a first-order estimate of the energetics of the disk wind detected in NED01. The values were estimated using only the 2012 NuSTAR and the 2014 Swift data where the wind was clearly detected.

Following the same arguments presented in many works on disk winds (e.g. Nardini et al. 2015; Reeves & Braito 2019; Braito et al. 2021), we estimate the mass outflow rate of the X-ray wind using the equation

$$\begin{array}{c}\hfill {\dot{M}}_{\mathrm{out}}^{\mathrm{X}-\mathrm{ray}}=\mathrm{\Omega }\mathrm{\mu }{\mathrm{m}}_{\mathrm{p}}{v}_{\mathrm{out}}{R}_{\mathrm{w}}{N}_{\mathrm{H}},\end{array}$$(1)

which assumes a biconical geometry for the flow (Krongold et al. 2007). In Eq. (1), μ is a constant factor set to μ = n_H/n_e = 1.2 for solar abundances, Ω is the wind solid angle, R_w is the disk wind launching radius, N_H and v_out are the column density and the velocity of the disk wind. The main uncertainties in Eq. (1) are the wind opening angle (Ω) and the launch radius R_w. We assume that the wind subtends Ω/4π = 0.5. This is justified by the systematic searches of ultra-fast disk winds in bright nearby AGN (Gofford et al. 2013; Tombesi et al. 2010), which resulted in a detection rate of about 40%, thus suggesting that the disk winds have a wide opening angle. A lower limit on the launching radius can be derived assuming that the wind is launched at its escape radius R_min = 2 GM_BH/v². We note that by adopting this value for R_w, we obtain the most conservative estimate of the mass outflow rate and energetics (see Gofford et al. 2015; Tombesi et al. 2012). Since the main uncertainty in the R_min derivation is the black hole mass of NED01, we decided to normalise the mass outflow rate by the Eddington rate,

$$\begin{array}{c}\hfill {\dot{M}}_{\mathrm{Edd}}=\frac{4\pi G{m}_{\mathrm{p}}{M}_{\mathrm{BH}}}{{\sigma }_{\mathrm{T}}\eta c},\end{array}$$(2)

where σ_T is the Thomson cross section and η = 0.1 is the accretion efficiency. Thus, combining Eqs. (1) and (2), and substituting for R_min, we obtain

$$\begin{array}{c}\hfill \frac{{\dot{M}}_{\mathrm{out}}^{\mathrm{X}-\mathrm{ray}}}{{\dot{M}}_{\mathrm{Edd}}}=2\frac{\mathrm{\Omega }}{4\pi }\mu N_{\mathrm{H}}{\sigma }_{\mathrm{T}}\eta {\left(\frac{v_{\mathrm{out}}}{c}\right)}^{-1}.\end{array}$$(3)

The X-ray wind kinetic power (${\dot{E}}_{\mathrm{kin}}^{\mathrm{X}-\mathrm{ray}}=1/2{\dot{M}}_{\mathrm{out}}{v}_{\mathrm{out}}^2$) normalised by the Eddington luminosity (L_Edd = η Ṁ_Edd c²) is

$$\begin{array}{c}\hfill \frac{{\dot{E}}_{\mathrm{kin}}^{\mathrm{X}-\mathrm{ray}}}{L_{\mathrm{Edd}}}=\frac{\mathrm{\Omega }}{4\pi }\mu N_{\mathrm{H}}{\sigma }_{\mathrm{T}}\frac{v_{\mathrm{out}}}{c}.\end{array}$$(4)

For the nuclear wind detected by Swift, we measure $N_H=1.3_{-0.9}^{+1.7}\times {10}^{23}$ cm⁻² and $v_{\mathrm{out}}/c=-0.{12}_{-0.05}^{+0.02}$. The (Swift) mass outflow rate is thus ${\dot{M}}_{\mathrm{out}}^{\mathrm{X}-\mathrm{ray}}/{\dot{M}}_{\mathrm{Edd}}\sim 0.09$ and the wind kinetic power is then ${\dot{E}}_{\mathrm{out}}^{\mathrm{X}-\mathrm{ray}}\sim 0.6$% of Eddington. Conversely, for the 2012 NuSTAR spectrum we derive ${\dot{M}}_{\mathrm{out}}^{\mathrm{X}-\mathrm{ray}}/{\dot{M}}_{\mathrm{Edd}}\sim 0.04$ and ${\dot{E}}_{\mathrm{out}}^{\mathrm{X}-\mathrm{ray}}\sim 1$% of Eddington. Considering a BH mass of ∼10⁸ M_⊙, as estimated by Alonso-Herrero et al. (2013) in NED01, these values correspond to ${\dot{E}}_{\mathrm{out}}^{\mathrm{Swift}}\sim 8\times {10}^{43}$ erg s⁻¹ and ${\dot{E}}_{\mathrm{out}}^{\mathrm{Nustar}}\sim 1.3\times {10}^{44}$ erg s⁻¹, where the higher value measured in the 2012 observation is driven by the higher velocity of the wind. We note that these estimates are clearly affected by large errors; only considering the 90% errors on the column density of the wind, we derive that the kinetic power of the X-ray wind detected in the Swift observations could range between 2.5 × 10⁴³ erg s⁻¹ and 1.7 × 10⁴⁴ erg s⁻¹, while for the NuSTAR wind the range is 0.6 − 2 × 10⁴⁴ erg s⁻¹. The momentum rate of the X-ray wind is ${\dot{p}}_{\mathrm{out}}^{\mathrm{X}-\mathrm{ray}}\sim (4\pm 2)\times {10}^{34}$ g cm s⁻² from NuSTAR, which is consistent with the value obtained from the Swift data, although the latter has a larger uncertainty.

9. Discussion

Since the large amount of energy released by AGN and the associated feedback processes affect the ISM on different physical scales, studying outflows in different ISM phases using multi-wavelength data is the only way to probe the AGN–ISM interplay. A multi-phase analysis is the only method that leads to a complete description of the outflows, allowing us to determine their physical parameters (e.g. their full extent and mass) and to obtain an overview of their driving mechanisms to test their theoretical models. The goal of our ALMA CO(1–0) observations was to investigate the presence of a high-velocity molecular outflow in NED01, which is expected to be experiencing strong AGN feedback since it hosts a powerful X-ray wind. At first look, the ALMA CO(1–0) data did not show any evidence for a massive powerful molecular outflow in this source, and instead presented clear evidence for CO(1–0) rotating disks in both NED01 and its companion galaxy NED02. Our BBarolo modelling confirmed that the bulk of the CO(1–0) emissions from the two galaxies can be fitted using a rotating disk model. We then focused on studying the presence of any CO(1–0) residual emission, not accounted for by the disk modelling, following a methodology commonly adopted in the literature (e.g. Sirressi et al. 2019; Bewketu Belete et al. 2021; Ramos Almeida et al. 2022).

By subtracting the best-fit BBarolo disk model from the ALMA datacube, we detected CO(1–0) residual emission at the ∼10% level around NED01, for all three models constructed by applying BBarolo on regions of different sizes (see also Appendix B). The fit that leaves the least residuals in NED01, with flux corresponding to ∼8% of its total CO(1–0) emission, is the one described in Sect. 5.1 and performed on the Big region. The subtraction of this fit leaves mainly blueshifted residual emission within v ∈ (− 400, −100) km s⁻¹, whose contours are enclosed in a squared box of 3 kpc in size, which is rather compact compared to the main rotating disk that has a diameter of ∼10 kpc. Quite strikingly, blueshifted residual CO(1–0) emission was detected consistently in the data after subtracting all three BBarolo best-fit disk models, with similar morphology and spectral properties, and so we believe it traces structures that may be truly disconnected from the main disk and whose origin is worth exploring further. On the other hand, redshifted residual emission was not detected in the fit performed on the Big region around NED01, although it was detected in the other two fits performed on smaller regions, which did not properly account for the extended redshifted CO(1–0) components (labelled B and C in the central panel of Fig. 4).

For NED02, after subtracting the best-fit BBarolo disk model from the ALMA data, we also detected residual CO(1–0) emission, corresponding to ∼7% of the total line flux from this galaxy. Similar to NED01, the residuals are blueshifted with respect to the CO(1–0) systemic velocity and have a compact morphology, being enclosed within the central 3 kpc of the source. In the following, we discuss possible interpretations for the blueshifted CO(1–0) residual emission detected in the IRAS 05054+1718 galaxy pair.

9.1. Hypothesis I: Molecular outflows

The first hypothesis we consider is that the CO(1–0) residuals trace a molecular outflow, which was the scenario we aimed to test with the acquisition of these new ALMA data.

Under the hypothesis that the CO(1–0) residual emission in NED01 traces a molecular outflow, we used the best-fit Gaussian parameters reported in Table 4 to estimate its energetics, by summing the contribution from the two blueshifted CO(1–0) spectral components. For each spectral feature, we set v_out equal to the Gaussian central velocity, the outflow radius R as the projected distance from the disk centre (see Tables 4 and 5 for NED01 and NED02, respectively), and compute the dynamical timescale of the outflow as ${\tau }_{\mathrm{dyn}}^{\mathrm{mol}}=R/v_{\mathrm{out}}$. Using the molecular gas mass inferred from the CO(1–0) luminosity reported in Table 4, we computed the mass-loss rate (${\dot{M}}_{\mathrm{out}}^{\mathrm{mol}}=M_{\mathrm{out}}^{\mathrm{mol}}/{\tau }_{\mathrm{dyn}}$), momentum rate ($v{\dot{M}}_{\mathrm{out}}^{\mathrm{mol}}$) and kinetic power ($0.5{\dot{M}}_{\mathrm{out}}^{\mathrm{mol}}{v}^2$) for each of the two CO residual spectral features, and report these numbers in Table 8. A similar exercise was performed for NED02, with results listed in Table 8.

The total properties of the putative outflow in NED01 were computed by summing the contribution of all the residual components obtaining a total mass-loss rate of ${\dot{M}}_{\mathrm{out}}^{\mathrm{mol}}=\sum _i{M}_{\mathrm{out},i}^{\mathrm{mol}}=19\pm 14$ M_⊙ yr⁻¹ (including a systematic uncertainty on the α_CO factor), a total momentum rate of ${\dot{p}}_{\mathrm{out}}^{\mathrm{mol}}=\sum _i{v}_i{\dot{M}}_{\mathrm{out},i}^{\mathrm{mol}}=(2.7\pm 2.4)\times {10}^{34}$ g cm s⁻¹, and a total kinetic power of ${\dot{E}}_{\mathrm{out}}^{\mathrm{mol}}=\sum _{i}0.5{\dot{M}}_{\mathrm{out},i}^{\mathrm{mol}}{v}_i^2=(3.2\pm 2.2)\times {10}^{41}$ erg s⁻¹. For NED02, we get a putative outflow mass-loss rate of 20 ± 16 M_⊙ yr⁻¹. However, a molecular outflow, if present, cannot be brighter than the residuals, implying all computed values as upper limits on the outflow energetics parameters.

In a blast-wave AGN feedback scenario (see models by Faucher-Giguère & Quataert 2012; Fabian 1999), the energetics of the molecular outflow detected on galactic scales should be compared with that of the X-ray disk wind, since the former is theoretically expected to be triggered by the latter. According to Ballo et al. (2015), and to the new analysis based on the 2012 NuSTAR data (see Sect. 6), the X-ray wind in NED01 has a velocity of v ∼ 0.2c with corresponding X-ray wind outflow and momentum rates values, although highly uncertain⁷, of ∼0.13 M_⊙ yr⁻¹ and ∼4 × 10³⁴ g cm s⁻¹, respectively. The kinetic power is estimated to range between 2.5 × 10⁴³ erg s⁻¹ and 1.7 × 10⁴⁴ erg s⁻¹ and 0.6 − 2 × 10⁴³ erg s⁻¹ for the Swift and the NuSTAR winds, respectively. This is two orders of magnitude higher than that inferred for the putative molecular outflow in this source. The corresponding ratio of the momentum rates of the molecular and X-ray winds is ${\dot{p}}_{\mathrm{out}}^{\mathrm{mol}}/{\dot{p}}_{\mathrm{out}}^{\mathrm{X}-\mathrm{ray}}\sim 0.67$, consistent with unity within the large uncertainties. Even in the (unlikely) scenario that all of the CO(1–0) residual emission detected in NED01 can be ascribed to a molecular outflow, there is no evidence for a boost in the momentum rate of the molecular outflow with respect to the X-ray wind, which would be an expectation of a ‘classic’ energy-driven outflow model (Faucher-Giguère & Quataert 2012). We can robustly conclude that the current data do not require any momentum rate boost for the large-scale outflow.

The position of NED01 in the outflow momentum rate ratio versus velocity diagram proposed by Feruglio et al. (2015), Tombesi et al. (2015) is shown in brown in Fig. 16, compared to other galaxies studied in the literature. One of these sources is the Seyfert 2 galaxy MCG-03-58-007, whose ALMA CO(1–0) observations were analysed by Sirressi et al. (2019). With a similar methodology to this work, based on the subtraction from ALMA CO(1–0) data of a rotating disk model obtained with the BBarolo code, Sirressi et al. (2019) detected compact CO(1–0) residual emission within ∼2 kpc of the AGN, presenting two symmetric components at blueshifted velocities and one at redshifted velocities (v = ±170 km s⁻¹). These authors cautiously interpreted such CO(1–0) residuals as possible evidence for either a compact molecular outflow or a kiloparsec-scale rotating structure disconnected from the main disk. In the outflow scenario, the momentum rate of the putative molecular outflow in MCG-03-58-007 would be ${\dot{p}}_{\mathrm{out}}^{\mathrm{mol}}\sim 6\times {10}^{34}$ g cm s⁻², corresponding to 40% of the momentum rate of the X-ray UFO, consistent with our new findings on IRAS 05054+1718 NED01. Overall, Fig. 16 shows that the interplay between nuclear disk winds and large-scale molecular outflows is more complex than predicted by theoretical blast-wave feedback models, and the data suggest a range of values for the outflow momentum rates on small and large scales.

Fig. 16. Outflow momentum rate plotted against its velocity for twelve objects with both ultra-fast winds and large-scale galactic outflows, adapted from Smith et al. (2019). Different colors and shapes identify different objects and outflows. Solid error bars are used for values with upper and lower errors calculated, dotted bars if only a range of value was available, and arrows for limits. Values reported in this plot taken from this work, Bischetti et al. (2019), Chartas et al. (2020), Cicone et al. (2014), Feruglio et al. (2015, 2017), Fluetsch et al. (2019), García-Burillo et al. (2014), González-Alfonso et al. (2017), Luminari et al. (2018), Longinotti et al. (2018), Marasco et al. (2020), Mizumoto et al. (2019), Longinotti et al. (2015), Rupke et al. (2017), Sirressi et al. (2019), Smith et al. (2019), Tombesi et al. (2015, 2017), Tozzi et al. (2021), Veilleux et al. (2017).

In order to directly compare observational values with predictions of the blast-wave model, Fig. 17 shows the ratio of the outflow rate to the X-ray UFO momentum rate, and the ratio of the estimated values of the related energy- and momentum-driven regimes. Following Marasco et al. (2020), Tozzi et al. (2021), we extended this analysis including our targets and some sources from Mizumoto et al. (2019). Since Marasco et al. (2020), Tozzi et al. (2021) adopted α_CO = 0.8 ± 1.2 M_⊙ (K km s⁻¹ pc²)⁻¹ (Downes & Solomon 1998; Carilli & Walter 2013; Bolatto et al. 2013) in their analysis, we re-derived the energetics applying this conversion factor, obtaining ${\dot{p}}_{\mathrm{out}}^{\mathrm{mol}}/{\dot{p}}_{\mathrm{out}}^{\mathrm{X}-\mathrm{ray}}\sim 0.3$. Overall, the energetics of the outflows appears to be consistent in most of the targets (15 out of 17 objects) with either a momentum-driven (12) or an energy-driven (3) scenario. There are only two cases where the outflow energetics seems completely unrelated to the nuclear source. The momentum ratio of the quasar SDSS J1353+1138 is ∼100 times smaller than the momentum-driven prediction. According to Tozzi et al. (2021), this low value may be associated with a massive molecular outflow not considered in that work or to a high variability among the quasar activity. The other exception is the Seyfert 1 galaxy IRAS 17020+4544 studied by Longinotti et al. (2015, 2018), with its extremely high molecular outflow momentum rate, possibly related to an uncertain estimate of ${\dot{p}}_{\mathrm{out}}^{\mathrm{mol}}$ or again to the high AGN variability. Overall, according to Marasco et al. (2020), Fig. 17 shows that the blast-wave scenario either in a momentum- or an energy-driven regime appears to describe the interplay between nuclear winds and large-scale outflows, although the majority of measurements suggests the absence of a significant momentum boost.

Fig. 17. Ratio of the large-scale outflow rate to the UFO momentum rate (adapted from Marasco et al. 2020; Tozzi et al. 2021). Different shapes identify different gas phases of the winds. The theoretical prediction for momentum-driven (Ṗout/ṖUFO = 1) and energy-driven winds are shown by a dashed line and orange rectangles, respectively. Empty and filled squares indicate Hα and [OIII]-based measurements. Values reported in this plot are taken from this work, Bischetti et al. (2019), Braito et al. (2018), Chartas et al. (2009), Cicone et al. (2014), Feruglio et al. (2015, 2017), García-Burillo et al. (2014), González-Alfonso et al. (2017), Longinotti et al. (2015, 2018), Lutz et al. (2020), Marasco et al. (2020), Mizumoto et al. (2019), Reeves & Braito (2019), Rupke et al. (2017), Sirressi et al. (2019), Smith et al. (2019), Tombesi et al. (2015), Tozzi et al. (2021), Veilleux et al. (2017).

Many of the local galaxies with a molecular outflow detection shown in Fig. 16 (but also see the review by Veilleux et al. 2020) are intensely star forming, and so their large-scale outflows may be driven by star formation. In comparison to the extreme ULIRGs that have been the focus of previous investigations, the SFRs of both members of the IRAS 05054+1718 galaxy pair are rather modest. From its total CO(1–0) line luminosity (see Table 3), by assuming an $L_{\mathrm{CO}(1-0)}^{\prime }$–SFR relation typical of local star-forming galaxies (see e.g. Cicone et al. 2017), we estimate a SFR_NED01($L_{\mathrm{CO}}^{\prime }$) = 2.5 ± 0.3 M_⊙ yr⁻¹ for NED01, which is lower than the range reported in the literature obtained from the IR luminosity and SED modelling (SFR = 5 − 10 M_⊙ yr⁻¹; see Sect. 2), possibly because of AGN contamination affecting photometric measurements of this source. For NED02, using the same procedure, we estimate from the total CO(1–0) luminosity SFR_NED02($L_{\mathrm{CO}}^{\prime }$) = 3.6 ± 0.5 M_⊙ yr⁻¹, which is also lower than the values reported in the literature based on IR luminosity. Since the presence of an AGN has not been confirmed in this source, the discrepancy could be due to a particularly high star formation efficiency in this source (SFE=SFR/M_mol) or to the presence of additional molecular gas undetected by our ALMA CO(1–0) line observations. As shown in Fig. 1, the CO(1–0) emission from NED02 is rather compact compared to its optical extent. According to Veilleux et al. (2005), assuming a maximum coupling efficiency of ∼10%, star formation feedback in galaxies with the SFRs inferred for NED01 and NED02 would drive weak outflows, with mass loss rates of Ṁ_out < 1 M_⊙ yr⁻¹, which would be undetectable in our data and as not consistent with the blueshifted CO residuals investigated here.

9.2. Hypothesis II: Compact rotating structures

An alternative hypothesis that can account for the residual emission is the presence of a compact structure, similar to the circumnuclear disks (CNDs) detected in some local AGNs (see e.g. García-Burillo et al. 2014; Combes et al. 2019), whose rotation departs from the rotating molecular disk modelled by BBarolo. Although this was a viable hypothesis for the CO(1–0) residual components detected in MCG-03-58-007 by Sirressi et al. (2019), it seems an unlikely scenario in the case of the galaxy pair IRAS 05054+1718 studied in this work. The CO(1–0) residuals in both NED01 and NED02 have a complex, asymmetrical, spectral profile, with multiple peaks at blueshifted velocities, and lack a reliable redshifted counterpart.

9.3. Hypothesis III: ISM disk asymmetries

Finally, the blueshifted CO residuals detected in both NED01 and NED02 could be due to other types of asymmetries of the molecular disks that are not properly captured by the BBarolo model. The negative residuals seen around the systemic velocity in both sources would also be consistent with this hypothesis. As discussed by Di Teodoro & Fraternali (2015), although BBarolo manages quite well minor asymmetries such as the disk warp in NGC 5055, the code always favours the model that generates the lowest residuals, which necessarily results in an average rotation curve that may follow one side of the disk more closely than the other. In the case of our target, the ongoing gravitational interaction between the two members of the pair could have already slightly perturbed their ISM, hence producing perturbations of the molecular disk that are not properly captured by BBarolo, causing the residual emission. We note that in this case, the residuals obtained after subtracting the best-fit BBarolo model may not correspond directly to the perturbed components of the ISM, as it could be that the model is forcing the fit to adapt better to the perturbations than to the rest of the (unperturbed) disk, hence producing artificial residuals with spectral and spatial characteristics that are unrelated to the underlying physical process. Therefore, this third scenario can only be tested using higher resolution data.

10. Summary and conclusions

In this study we presented new high-sensitivity ALMA CO(1–0) observations of the interacting galaxy pair IRAS 05054+1718 (projected distance between the components ∼29.6″ ∼ 11 kpc).

We investigated the presence of a galactic-scale molecular outflow in the main target NED01 (z = 0.0178 ± 0.0004, this work) and its possible relation with the powerful X-ray wind detected in its nucleus (see Sect. 7 and Ballo et al. 2015). Given its moderate SFR (∼5 M_⊙/yr, De Looze et al. 2014; Pereira-Santaella et al. 2015), NED01 is a suitable object to study AGN feedback mechanisms and the blast-wave scenario, usually tested in environments that are extremely bright in the IR and with high SFRs where the AGN contribution is hard to isolate from the stellar contribution. The companion NED02 (z = 0.016812 ± 0.000003, this work) was included in the field of view of our compact-configuration ALMA data. This source is a LIRG, classified as composite by Pereira-Santaella et al. (2015), but there is no clear evidence of the presence of an AGN (Alonso-Herrero et al. 2012).

We investigated the distribution and kinematics of the CO(1–0) emission at a spatial resolution of ∼1 kpc, which appears dominated by the rotating molecular gas disk for both targets. We modelled the disk using 3D-BAROLO and studied the residual emission to identify structures deviating from the rotation. We detected residual emission (S/N ≥ 15) at v ∈ (− 400, −100) km s⁻¹ corresponding to the 8% of its total CO(1–0) emission in the main target NED01, characterised by a compact morphology and extending within 3 kpc of the galaxy centre. We also detected a similar blueshifted residual CO(1–0) emission (∼7% of the total CO(1–0) flux from this galaxy) in the companion NED02. A molecular outflow in NED01, if present, cannot be brighter than this residual, implying an upper limit on its mass-loss rate of 19 ± 14 M_⊙ yr⁻¹, on its momentum rate of (2.7 ± 2.4)×10³⁴ g cm s⁻¹, and on its kinetic power of (3.2 ± 2.2)×10⁴¹ erg s⁻¹. For NED02, we obtain a putative outflow mass-loss rate of 20 ± 16 M_⊙ yr⁻¹.

Regarding the X-ray UFO, the analysis performed on the 2012 NuSTAR data confirmed the results already found by Ballo et al. (2015) regarding the presence of a highly ionised and high-velocity outflowing disk wind in the NED01 nucleus. Although no disk wind features have been detected in the more recent XMM-Newton and NuSTAR data, the upper limits estimated for the column density of a putative ionised absorber outflowing are only marginally below the 90% lower value of N_H estimated in the past observations. This suggests that a small change in the opacity of the wind can explain a non-detection of the wind at the low-flux state during 2021. This, combined with the factor of 4 lower intrinsic flux of the 2021 X-ray observations, can account for the non-detection in the most recent data.

As for the companion NED02, we derive a faint soft (0.5–2 keV) flux of the order of ∼2 × 10⁻¹⁴ erg cm⁻² s⁻¹, corresponding to L_{(0.5 − 2) keV} ∼ 2 × 10⁴⁰ erg s⁻¹. According to the L_X − SFR relations derived from several surveys of LIRGs and ULIRGs (see e.g. Cicone et al. 2017) this would correspond to a SFR ∼ 7 M_⊙ yr⁻¹. At higher energies, only X-ray observations with the high spatial resolution, such as that offered by the Chandra observatory, will allow us to properly separate the X-ray emission of NED02 from the brighter companion NED01 and reach a definitive answer on the possible presence of a weak and/or obscured AGN in NED02.

We tested the AGN feedback blast-wave scenario in NED01 by comparing the energetics of the UFO and the putative molecular outflow. The results we obtained are not consistent with an energy-driven scenario since there is no evidence for a boost in the momentum rate of the molecular outflow with respect to the X-ray wind. We also considered the hypothesis of an outflow triggered by star formation activity. The SFR of NED01, equal to ${\mathrm{SFR}}_{\mathrm{NED}01}({L\prime }_{\mathrm{CO}})=2.5\pm$ 0.3 M_⊙ s⁻¹, refined using its CO(1–0) luminosity and a local tight L_CO − SFR relation, is too low to account for a putative star formation-driven molecular outflow.

Since the presence of an AGN has not been confirmed in the companion NED02, the hypothesis of AGN feedback cannot apply to this case. Following the hypothesis of a star formation-driven outflow, we estimated CO(1–0) luminosity in NED02 as ${\mathrm{SFR}}_{\mathrm{NED}01}({L\prime }_{\mathrm{CO}})=3.6\pm$ 0.5 M_⊙ s⁻¹, lower than the values reported in the literature based on IR luminosity. As for the main target, the estimated SFR is too low to account for the observed residual emission.

Because of the ambiguous interpretation of the CO residual emission obtained after subtracting the best-fit disks modelled by BBarolo, the residuals observed in NED01 and NED02 may have a different explanation from that of a molecular outflow. We considered the alternative hypothesis that the detected residual CO emission is due to an additional compact structure similar to a circumnuclear disk whose rotation cannot be associated with the central molecular disk. This hypothesis seems unlikely since the residuals detected in NED01 and NED02 show a complex asymmetrical profile without a redshifted counterpart. We also considered a third hypothesis, according to which the residuals may be due to disk asymmetries that are not properly captured by the BBarolo model. In this case the position and spectral properties of such residuals may not be highly informative of the underlying physical process responsible for the disk perturbations, and only higher resolution data could reveal their origin.

Our study, together with other recent works in the literature (e.g. Smith et al. 2019; Sirressi et al. 2019; Veilleux et al. 2020), has increased the sample of AGN hosting an ultra-fast X-ray wind with observational constraints on the presence and properties of a molecular outflow, on which it is possible to test the predictions of the AGN-driven feedback scenario. The increased statistics of observations appears to complicate the scenario of AGN-driven feedback mechanisms more and more, drawing a much more complex picture than expected. Furthermore, the high sensitivity to extended structures enabled by the use of a compact array configuration has revealed important information on the close environment of the target galaxy. These data have indeed allowed us to study NED01 simultaneously with its nearby companion NED02, and to investigate the low surface brightness components of their respective molecular disks, hence to compare their molecular ISM morphologies. Our understanding of diffuse low surface brightness and extended molecular gas components due to feedback processes and gravitational interactions will be transformed by a facility such as the Atacama Large Aperture Submillimetre Telescope (AtLAST)⁸.

Acknowledgments

This paper makes use of the following ALMA data: ADS/JAO.ALMA#2016.1.00694.S. ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan), together with NRC (Canada), MOST and ASIAA (Taiwan), and KASI (Republic of Korea), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO and NAOJ. This project has received funding from the European Union’s Horizon 2020 research and innovation programme under grant agreement No 951815 (AtLAST). FB acknowledges funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (Grant agreement No. 851435). CV acknowledges funding from PRIN MIUR 2017PH3WAT (‘Black hole winds and the baryon life cycle of galaxies’). CC, PS, VB, CV, RDC, LB and MD acknowledge financial contributions from Bando Ricerca Fondamentale INAF 2022 Large Grant “Dual and binary supermassive black holes in the multi-messenger era: from galaxy mergers to gravitational waves” and from the agreement ASI-INAF n.2017-14-H.O. VB acknowledges financial support through the NASA grant 80NSSC22K0220. VB and JR acknowledge support from NASA grant 80NSSC22K0474.

References

Appendix A: ALMA map of the 3mm continuum

The ALMA interferometric map of the 3 mm continuum emission from the target is shown in Figure A.1. The map was obtained by combining the data (total bandwidth Δν = 4 GHz) from the two low-resolution spectral windows centred at 101.190 GHz and 99.387 GHz, which are not contaminated by line emission. A 2D Gaussian fit performed on the PB-corrected map shown in Fig. A.1 delivers best-fit peak flux densities of $S_{\mathit{cont}}^{NED01}=0.71\pm 0.03$ mJy beam⁻¹ and $S_{\mathit{cont}}^{NED02}=2.26\pm 0.05$ mJy beam⁻¹, indicating that the 3 mm continuum of NED02 is approximately three times brighter than NED01. The fit also provides information about the size and position of the two sources. For NED01, we get a peak centre at RA(ICRS) = 05:08:19.713 and Dec(ICRS) = 17.21.48.11, with a deconvolved size of 3.0″ × 2.4″ (PA = 84.5 deg). For NED02, the ALMA continuum peaks at RA(ICRS) = 05:08:21.211 and Dec(ICRS) = 17.22.08.40 (PA = 98.4 deg), with a deconvolved size of 3.0″ × 2.5″. The source-integrated continuum flux obtained through the fit is S_3mm = 0.76 mJy for NED01 and S_3mm = 2.4 mJy for NED02. These results indicate that the 3 mm continuum emission from the targets is only marginally resolved at the resolution of our ALMA data, and much more compact than the CO(1–0) line emission. This is further shown by Figure A.2, which displays the overlay of the ALMA 3 mm continuum contours on the CO(1–0) line emission map.

Fig. A.1. ALMA 3 mm continuum map computed at an average frequency of 100 GHz. The map was obtained by averaging the data over a bandwidth of 4 GHz. Contours are plotted at (− 3, 3, 6, 9, 12, 18)×σcont, with σcont = 0.03 mJy/beam being the average 1σ rms noise level within the central 10″ portion of the FoV.

Fig. A.2. Overlay of ALMA 3 mm continuum contours (in black) on the CO(1–0) line map. This highlights the compactness of the continuum emission compared to the line emission in both members of the IRAS 05054+1718 galaxy pair. Neither the continuum contours nor the CO(1–0) line map is corrected for the PB, to ease the visualisation of the companion galaxy NED02.

Appendix B: Alternative BBarolo CO(1–0) disk models for NED01

The Small region, shown in Fig. B.1, includes only the central molecular disk of NED01 (i.e. component A in Fig. 4), whereas the Medium region also embeds the southern extension (component C) of the CO(1–0) emission from NED01. We ran the BBarolo disk modelling for NED01 on these two additional regions, and the results are reported in Fig. B.3. For the Small region we used three rings, each 2″ in width, while for the Medium region we increased the number of rings to seven. The rotation patterns resulting from these two fits does not differ significantly from the fit shown in the main text, which was performed by selecting the Big region around NED01. However, these two fits, performed on smaller areas, produce more significant residuals, at both blue- and redshifted velocities, as shown in Fig. B.4. The best-fit parameters of the CO(1–0) residuals obtained in these two fits performed on the Small and Medium regions are listed in Table B.1.

Fig. B.1. Velocity map of NED01 (see also Fig. 4), with the three regions used for the BBarolo analysis overplotted.

Fig. B.2. CO channel maps of NED02, computed by integrating CO(1–0) emission in channels of dv = 47 km/s. This shows the comparison between the ALMA data and the BBarolo disk model. The contours are plotted at 3, 6, 9, 15, 20, 30, 40, 50σ.

Fig. B.3. CO(1–0) moment maps showing comparison between data and BBarolo best-fit molecular disk model performed on the Small (upper panel) and Medium (lower panel) regions around NED01.

Fig. B.4. CO(1–0) residual emission around NED01 obtained after subtracting the best-fit disk model computed by BBarolo on the Small region (top panels) and on the Medium region (bottom panels). Top left: Map of the CO(1–0) residual emission integrated within v ∈ [170, 220] km s−1 (red contours) and v ∈ [ − 290, −200] km s−1 (blue contours). Contours are plotted at 3, 6, 9, 15, 20σ. The black cross indicates the centre of the best-fit disk model, while the red and blue crosses indicate the centroids of the red- and blueshifted residuals, respectively. Top centre: Spectrum of the CO(1–0) residual emission extracted from the black box (size 9.2″×6.4″) in the left panel. Top right: Imaging of the blueshifted CO(1–0) residual components identified in the spectrum. Bottom left: Map of the CO(1–0) residual emission integrated within v ∈ [170, 220] km s−1 (red contours) and v ∈ [ − 290, −120] km s−1 (blue contours). Bottom centre: Spectrum of the CO(1–0) residual emission extracted from the black box (size of 11.2″×7.2″) in the left panel. Bottom right: Imaging of the blueshifted CO(1–0) spectral components identified in the central panel.

Quite strikingly, the spectral profile and morphology of the blueshifted residuals are consistent with those produced by the fit performed on the Big region, shown in the main body of the paper shown in Sect. 5.1. This is best visualised in the left panel of Fig. B.5, which compares on the same map the blueshifted CO(1–0) residual emission obtained from all three fits performed on NED01. Our analysis focuses on these blueshifted residuals, which are persistent in all three fits. In addition, the two fits reported here produce residuals at receding velocities, east of the NED01 galaxy centre, which were not detected in the Big region fit. The right panel of Fig. B.5 compares the morphology of this redshifted CO(1–0) residual emission obtained in the Small and Medium region fits. The contours from the two fits overlap only partially. These redshifted residuals must be due to the additional extended CO(1–0) features that are not captured in these smaller regions, but were instead well modelled by the fit performed on the Big region, where they displayed redshifted velocities on the eastern side of the galaxy consistent with components that follow the same pattern as the main disk rotation.

Fig. B.5. Comparison between CO(1–0) residual emission revealed in NED01 after subtracting the different BBarolo disk models. The green, red, and blue contours show respectively the CO(1–0) residuals obtained in the Big, Medium, and Small regions fits, and they are plotted at 5, 10, 15, 20σ. Left panel: Blueshifted residuals. Right panel: Redshifted residuals.

All Tables

All Figures

Fig. 1. ALMA CO(1–0) map (black contours) overlayed onto the g-band optical image from the Pan-STARR Survey 1. The ALMA CO(1–0) emission was averaged over a spectral range corresponding to v = [ − 490, +230] km s−1 with respect to the systemic redshift of NED01, which includes the CO(1–0) emission from both members of the galaxy pair. The contours correspond to the (3, 6, 9, 12, 24, 50) × σRMS levels, with σRMS = 0.2 mJy beam−1 being the average rms of the ALMA CO(1–0) map (not corrected for the primary beam).

In the text

Fig. 2. ALMA CO(1–0) continuum-subtracted spectra of the interacting galaxy pair IRAS 05054+1718: NED01 (left panels) and NED02 (right panels). The top and bottom panels display the same spectral data, but with different units on the x-axes. In particular, the top panels report the CO(1–0) flux density as a function of observed (sky) frequency, and were used to refine the systemic redshift estimates for the two galaxies (see Sect. 4.1). The bottom panels show the spectra corrected for redshift, where the CO(1–0) flux density is reported as a function of optical velocity along the line of sight. The CO(1–0) line spectrum of NED01 was extracted from an elliptical aperture maximising the CO(1–0) flux, with size 23″ × 14″, centred at RA = 05h08m19.858s, Dec = +17°21′45.898″. The CO(1–0) spectrum of NED02 was extracted from an 10″ × 12″ elliptical aperture centred at RA = 05h08m21.212s, Dec = +17°22′08.660″. The best-fit spectral parameters are listed in Table 2.

In the text

Fig. 3. CO(1–0) PV diagrams of NED01, extracted from a datacube where velocities are calculated with respect to the previously known redshift of the source (z = 0.0184). The PV diagrams confirm that our new redshift estimate of z = 0.0178 ± 0.0004, indicated with a black cross, closely matches the kinematic centre of the CO(1–0) source. The upper panel shows the PV diagram obtained from a slit-like aperture along the axis connecting the blue and red peaks of the CO(1–0) emission, with a size of ∼30″ and a position angle of 95° (measured anti-clockwise from the north direction). The PV diagram shown in the bottom panel was computed from a slit-like aperture with a size of ∼9″ orthogonal to the previous one.

In the text

Fig. 4. Moment maps of NED01 computed within a box of size 30″ and including pixels above a CO(1–0) intensity threshold of 0.01 Jy beam−1 and within v ∈ [ − 500, 500] km s−1. Left: Intensity (moment 0) map, with contours plotted at [0.1, 0.3, 1, 3, 5, 8] Jy beam−1 km s−1. Centre: Velocity (moment 1) map, with contours plotted every 30 km s−1; the letters indicate the three CO(1–0) components discussed in the main text, where A indicates the main rotating CO(1–0) disk. Right: Velocity dispersion (moment 2) map, with contours plotted at intervals of 30 km s−1.

In the text

Fig. 5. Moment maps of NED02 computed from the non-PB-corrected datacube within a box of size 14″ centred at RA(ICRS) = 05:08:21.108; Dec(ICRS) = 17.22.08.489, including pixels above a CO(1–0) intensity threshold of 0.009 Jy beam−1 and within v ∈ [ − 500, 500] km s−1. Left: Intensity (moment 0) map, with contours plotted at [0.1, 1, 5, 8, 11] Jy beam−1 km s−1. Centre: Velocity (moment 1) map, with contours plotted every 30 km s−1. Right: Velocity dispersion (moment 2) map, with contours plotted at intervals of 20 km s−1.

In the text

Fig. 6. Moment maps showing the comparison between the data and the BBarolo best fit performed on the Big region around NED01. Shown (from left to right) are data, BBarolo model, and residual emission. The rows show (rom top to bottom) the moment 0 map (intensity), moment 1 map (line of sight velocity), and moment 2 map (velocity dispersion). The black cross indicates the rotating disk centre found by the fit, and the dashed black line shows the major axis of the disk model.

In the text

Fig. 7. Position velocity diagrams showing the comparison between the data and the BBarolo best-fit model performed on the Big region around NED01. The red contours show the disk model, and they are overplotted onto the observed ALMA CO(1–0) line data (displayed in grey, with blue contours at 3, 6, 9, 15, 20, 30, 40, and 50σ). The upper panel shows the PV diagram extracted along the major axis. The yellow data points indicate the best-fit rotational velocity of each ring. The lower panel shows the PV diagram computed along the minor axis.

In the text

Fig. 8. CO channel maps computed by integrating the CO(1–0) emission in channels of Δv = 42 km s−1, showing the comparison between data and BBarolo best-fit disk model within the Big region around NED01. In each of the three panels, the top row (blue contours) displays the data, the central row the model (red contours), and the bottom row the residual emission (green contours). The contours are plotted at 3, 6, 9, 15, 20, 30, 40, 50σ. Symmetric negative contours are plotted in grey. The yellow cross indicates the centre of the disk model.

In the text

Fig. 9. Residual CO(1–0) emission in NED01 after subtracting the best-fit disk model computed by BBarolo in the Big region. Left: Map of the blueshifted CO(1–0) residual emission integrated in the range −290 ≲ v[km s−1] ≲ − 200. The blue contours are plotted at 3, 6, 9, 15, and 20σ. The black cross indicates the centre of the best-fit molecular disk. The blue cross indicates the peak of the CO(1–0) residual emission. The residual structure can be enclosed within a box of 8.4″ × 8.4″ (∼3 kpc, see black box overlayed on the map). Centre: Spectrum of the CO(1–0) residual emission extracted from an aperture corresponding to the black box shown in the left panel. The two CO(1–0) peaks are fit with double-Gaussian functions whose parameters are reported in Table 4. Right: Map of the two blueshifted spectral components visible in the central panel, plotted with different colours. The first component (blue contours) is the most blueshifted spectral peaks, integrated within v ∈ [ − 400, −200] km s−1; the second component (in black) corresponds to the secondary peak, integrated between v ∈ [ − 200, 0] km s−1. Contours are plotted at 3, 6, 9, 15, and 20σ.

In the text

	Fig. 10. CO(1–0) moment maps of the companion galaxy NED02 showing the comparison between the data and the best-fit molecular disk model found by BBarolo.
In the text

	Fig. 11. Position velocity diagrams of NED02 showing the comparison between the data (blue contours, plotted at 3, 6, 9, 15, 20, 30, 40, and 50σ) and the BBarolo best-fit rotating disk model (red contours). The upper panel shows the PV diagram extracted along the major axis, and the lower panel shows the one extracted along the minor axis of rotation.
In the text

Fig. 12. Residual CO(1–0) emission in NED02 after subtracting the best-fit rotating disk model computed by BBarolo. Left: Map of the CO(1–0) residual emission integrated within v ∈ [ − 170, 20] km s−1. Contours are plotted at 3, 6, 9, 15, and 20σ. The black cross indicates the dynamic centre of the rotating disk and the blue cross shows the centroid of the residual emission. The residual structure can be enclosed within a box of 8″ × 8″ (∼3 kpc, see black box overlayed on the map). Centre: Spectrum of the CO(1–0) residual emission in NED02, extracted from the squared aperture shown in the left panel. The results of the spectral fitting with Gaussian functions are given in Table 5. Right: Map of the two blueshifted spectral components visible in the residual spectrum shown in the central panel. The first component, integrated within v ∈ [ − 350, −120] km s−1, is displayed using blue contours. The second component, integrated within v ∈ [ − 120, 0] km s−1 is shown in black. Contours are plotted at 3, 6, 9, 15, and 20σ.

In the text

Fig. 13. Broad-band (rest frame 0.3–60 keV) X-ray spectra of all the observations of IRAS 05054+1718. The EFE spectra are obtained unfolding the data against a simple power-law model with Γ = 2. The 2012 NuSTAR observation is shown in red, the 2014 Swift is in black, while the spectra obtained with the coordinated XMM-Newton and NuSTAR observations are shown in dark and light blue, respectively. We also included the averaged Swift-BAT spectrum from the 70-month survey (Baumgartner et al. 2013). The residuals obtained from the joint fit for the Swift and NuSTAR spectra using a typical model for a moderately obscured Seyfert 2 as NED01 are shown in the inset. The model, composed of an absorbed primary power-law component, a scattered component, with the same photon index (Γ), and a reflected component, is able to reproduce the broad-band emission, but it leaves some residuals in absorption (see Sect. 7.1).

In the text

Fig. 14. Spectra and best-fit model (blue line) for the Swift (black data points) and 2012 NuSTAR (red) observations. The continuum model includes an absorbed power-law component and a reflected component (grey dotted lines). The reflected component (modelled with XILLVER) also includes the expected Fe Kα emission line at 6.4 keV. The NH of the neutral absorber, and the normalisation of the primary power law and of the reflected components are all allowed to vary. The model also includes an ionised outflowing absorber, for which the same ionisation is assumed, but the column densities and velocity are allowed to differ in the two spectra. The Swift data were rebinned for plotting purposes.

In the text

	Fig. 15. Spectra and best-fit model (red line) for the 2021 XMM (blue data points) and NuSTAR (light blue) observation. The continuum model includes an absorbed power-law component and a reflected component (grey dotted lines).
In the text

Fig. 16. Outflow momentum rate plotted against its velocity for twelve objects with both ultra-fast winds and large-scale galactic outflows, adapted from Smith et al. (2019). Different colors and shapes identify different objects and outflows. Solid error bars are used for values with upper and lower errors calculated, dotted bars if only a range of value was available, and arrows for limits. Values reported in this plot taken from this work, Bischetti et al. (2019), Chartas et al. (2020), Cicone et al. (2014), Feruglio et al. (2015, 2017), Fluetsch et al. (2019), García-Burillo et al. (2014), González-Alfonso et al. (2017), Luminari et al. (2018), Longinotti et al. (2018), Marasco et al. (2020), Mizumoto et al. (2019), Longinotti et al. (2015), Rupke et al. (2017), Sirressi et al. (2019), Smith et al. (2019), Tombesi et al. (2015, 2017), Tozzi et al. (2021), Veilleux et al. (2017).

In the text

Fig. 17. Ratio of the large-scale outflow rate to the UFO momentum rate (adapted from Marasco et al. 2020; Tozzi et al. 2021). Different shapes identify different gas phases of the winds. The theoretical prediction for momentum-driven (Ṗout/ṖUFO = 1) and energy-driven winds are shown by a dashed line and orange rectangles, respectively. Empty and filled squares indicate Hα and [OIII]-based measurements. Values reported in this plot are taken from this work, Bischetti et al. (2019), Braito et al. (2018), Chartas et al. (2009), Cicone et al. (2014), Feruglio et al. (2015, 2017), García-Burillo et al. (2014), González-Alfonso et al. (2017), Longinotti et al. (2015, 2018), Lutz et al. (2020), Marasco et al. (2020), Mizumoto et al. (2019), Reeves & Braito (2019), Rupke et al. (2017), Sirressi et al. (2019), Smith et al. (2019), Tombesi et al. (2015), Tozzi et al. (2021), Veilleux et al. (2017).

In the text

	Fig. A.1. ALMA 3 mm continuum map computed at an average frequency of 100 GHz. The map was obtained by averaging the data over a bandwidth of 4 GHz. Contours are plotted at (− 3, 3, 6, 9, 12, 18)×σcont, with σcont = 0.03 mJy/beam being the average 1σ rms noise level within the central 10″ portion of the FoV.
In the text

	Fig. A.2. Overlay of ALMA 3 mm continuum contours (in black) on the CO(1–0) line map. This highlights the compactness of the continuum emission compared to the line emission in both members of the IRAS 05054+1718 galaxy pair. Neither the continuum contours nor the CO(1–0) line map is corrected for the PB, to ease the visualisation of the companion galaxy NED02.
In the text

	Fig. B.1. Velocity map of NED01 (see also Fig. 4), with the three regions used for the BBarolo analysis overplotted.
In the text

	Fig. B.2. CO channel maps of NED02, computed by integrating CO(1–0) emission in channels of dv = 47 km/s. This shows the comparison between the ALMA data and the BBarolo disk model. The contours are plotted at 3, 6, 9, 15, 20, 30, 40, 50σ.
In the text

	Fig. B.3. CO(1–0) moment maps showing comparison between data and BBarolo best-fit molecular disk model performed on the Small (upper panel) and Medium (lower panel) regions around NED01.
In the text

Fig. B.4. CO(1–0) residual emission around NED01 obtained after subtracting the best-fit disk model computed by BBarolo on the Small region (top panels) and on the Medium region (bottom panels). Top left: Map of the CO(1–0) residual emission integrated within v ∈ [170, 220] km s−1 (red contours) and v ∈ [ − 290, −200] km s−1 (blue contours). Contours are plotted at 3, 6, 9, 15, 20σ. The black cross indicates the centre of the best-fit disk model, while the red and blue crosses indicate the centroids of the red- and blueshifted residuals, respectively. Top centre: Spectrum of the CO(1–0) residual emission extracted from the black box (size 9.2″×6.4″) in the left panel. Top right: Imaging of the blueshifted CO(1–0) residual components identified in the spectrum. Bottom left: Map of the CO(1–0) residual emission integrated within v ∈ [170, 220] km s−1 (red contours) and v ∈ [ − 290, −120] km s−1 (blue contours). Bottom centre: Spectrum of the CO(1–0) residual emission extracted from the black box (size of 11.2″×7.2″) in the left panel. Bottom right: Imaging of the blueshifted CO(1–0) spectral components identified in the central panel.

In the text

	Fig. B.5. Comparison between CO(1–0) residual emission revealed in NED01 after subtracting the different BBarolo disk models. The green, red, and blue contours show respectively the CO(1–0) residuals obtained in the Big, Medium, and Small regions fits, and they are plotted at 5, 10, 15, 20σ. Left panel: Blueshifted residuals. Right panel: Redshifted residuals.
In the text