© The Authors 2026

1. Introduction

The gas surrounding galaxies plays a crucial role in regulating their formation and evolution and tracks the cycle and exchange of baryons between galaxies and their environments (e.g., reviews by Tumlinson et al. 2017; Péroux & Howk 2020). Observing the cosmic gas from the filaments of the intergalactic medium (IGM) to the circumgalactic medium (CGM) provides a comprehensive view of the buildup of structures, especially around the epoch of the peak of the cosmic star-formation rate density (Madau & Dickinson 2014).

Direct imaging of diffuse gas surrounding galaxies not only reveals its properties but also serves as a powerful diagnostic of the gas-galaxy connection. However, the diffuse nature of this low-density gas makes its detection challenging. Thanks to the advent of integral-field spectrographs such as the Multi Unit Spectroscopic Explorer (MUSE, Bacon et al. 2010) at the Very Large Telescope (VLT), the H I gas has become observable via its Lyαλ1216 emission (see the reviews by Cantalupo et al. 2017; Ouchi et al. 2020) down to surface brightness (SB) levels of ≈ 10⁻²⁰ − 10⁻¹⁸ erg s⁻¹ cm⁻² arcsec⁻². Systematic campaigns around quasar fields at z ≳ 3 have unveiled the ubiquitous presence of Lyα gaseous nebulae extended over ≳ 100 kpc (e.g., Steidel et al. 2000; Cantalupo et al. 2014; Hennawi et al. 2015; Borisova et al. 2016; Arrigoni Battaia et al. 2019; Fossati et al. 2021) and filamentary structures connecting overdensities of active galactic nuclei (AGNs) and star-forming galaxies (Umehata et al. 2019; Bacon et al. 2021; Banerjee et al. 2025; Tornotti et al. 2025a). Although being a powerful tracer of cool diffuse gas, the Lyα emission suffers from strong radiative transfer effects, hampering the precise determination of the systemic redshift of the quasar host galaxies, and a reliable interpretation of the gas kinematics, which is essential to resolve the gas flows in the CGM. Detecting emission lines other than the Lyα is therefore crucial.

At z ≳ 3, the C IVλ1548 and the nonresonant He IIλ1640 transitions can be observed at optical wavelengths (e.g., Villar-Martin et al. 1997; Arrigoni Battaia et al. 2015b,a; Cai et al. 2016, 2017; Marino et al. 2019; Marques-Chaves et al. 2019). The combination of such emission lines allows the intrinsic properties of the gas to be constrained, such as its density, clumpiness, and total mass, as well as the powering mechanism (e.g., Cantalupo et al. 2019; Pezzulli & Cantalupo 2019). Although sparse, the study of such extended nebulae at z ≈ 2 − 4 points to a picture in which the gas surrounding galaxies is the result of coexisting several processes and components including multiple intersecting filaments (Herenz et al. 2020), inflows of gas (Zhang et al. 2023b,a), interacting clouds (Jiménez-Andrade et al. 2023), and galaxies (Sabhlok et al. 2024).

A view of the cosmic gas complementary to that provided by rest-frame UV lines of the ionized phases comes from (sub)millimeter surveys, especially those conducted with the Atacama Large Millimeter Array (ALMA). These provide information on the dust and the cold molecular gas content, and allow a robust determination of the quasar host galaxy systemic redshift (e.g., Ginolfi et al. 2017; Decarli et al. 2019a; Li et al. 2021; Emonts et al. 2023; Li et al. 2023; Sulzenauer et al. 2025). Furthermore, especially in overdense regions, such observations at longer wavelengths have revealed the presence of a population of massive, star-forming, and heavily dust-obscured galaxies that were missed in optical data (e.g., Hodge et al. 2013; Miller et al. 2018; Oteo et al. 2018; Jin et al. 2021; Umehata et al. 2015, 2018, 2019; Pensabene et al. 2024, 2025), completing the view of the gas-galaxy connection and the effects of the large-scale environment on different galaxy populations.

In this work, we focus on the MUSE Quasar Nebula 04 (MQN04) surrounding the quasar Q0055-269 at z ≈ 3.66. This source was originally included in the blind survey for extended Lyα nebulae conducted by Borisova et al. (2016), built on a MUSE observational campaign around 17 radio-quiet quasars at 3.0 < z < 3.9 (see Sect. 2.1). They found that this system hosts one of the brightest Lyα-emitting nebulae (L_Lyα = 4 × 10⁴⁴ erg s⁻¹) discovered so far at z ≳ 3, extending over ≈180 kpc, which also exhibits exceptionally bright C IV and He II emission spread across ≳ 50 kpc (see also Guo et al. 2020). In particular, all the tracers revealed the presence of two compact regions located ≈2.5″ away from the quasar. These properties make MQN04 a remarkable laboratory for investigating the structure, ionization, and kinematics of the CGM around quasars. Despite these advances, key aspects of this system remain uncertain. In particular, the lack of a precise systemic redshift for the quasar has so far limited our ability to interpret the kinematics of the surrounding gas. Moreover, the spatial relationship between the extended emission, possible companion galaxies or AGNs, and the cold gas phase traced by millimeter observations remains unexplored. To address these questions, we combined new ALMA Band 3 observations with existing MUSE data, providing the first measurement of the quasar systemic redshift and enabling a direct comparison between the ionized and molecular gas phases. Together with new high-resolution MUSE Narrow-Field Mode (NFM) observations, this dataset allowed us to (i) reliably trace the gas kinematics of the nebula through nonresonant He II emission, (ii) search for companion galaxies or AGNs within the inner few kiloparsecs, and (iii) investigate the morphology and clumpiness of the CGM on small scales. These complementary observations offer an unprecedented opportunity to understand the interplay between quasar activity, gas flows, and the assembly of massive galaxies at high redshift.

The structure of the paper is as follows: we describe the details and reduction processes of the data in Sect. 2. In Sect. 3, we describe the identification of galaxies in the field. In Sect. 4, we present the detection and properties of the extended gaseous nebulae, and a characterization of the large-scale galaxy environment. In Sect. 5, we discuss possible scenarios that can be invoked to explain the observed emission. Finally, we summarize the key results of this work and draw our conclusions in Sect. 6. Throughout this work, we adopt a standard Λ cold dark matter cosmology with H₀ = 67.7 km s⁻¹ Mpc⁻¹, Ω_m = 0.310, and Ω_Λ = 1 − Ω_m from Planck Collaboration VI (2020), make use of AB magnitude, report distances in physical units, and quote 3σ limits unless stated otherwise.

2. Observations and data reduction

2.1. MUSE observations

We analyze MUSE observations¹ centered on the quasar Q0055-269 at z ≊ 3.66. These have been collected as part of the MUSE Guaranteed Time Observations (GTO; program IDs 094.A-0131(B) and 096.A-0222(A), PI J. Schaye) and are included in the MUSE Quasar-field Blind Emitters Survey (MUSEQuBES; Muzahid et al. 2020). The data consists of a combination of 14 exposures taken in WFM-NOAO-N mode and resulting in an integrated exposure time of 10h on source. To mitigate the impact of bad pixels and cosmic rays, each exposure was shifted by a small offset and rotated by 90° relative to the previous one. The MUSE field-of-view (FoV) spans 1′×1′ with a sampling of 0.2″ per pixel and covers the optical wavelength range 4800–9300 Å with a sampling of 1.25 Å and a spectral resolution ranging from R ≈ 1800 to ≈ 3600. We refer to Muzahid et al. (2021) for details about the data reduction, while only a brief description of the key steps is given here. First, the raw data were processed using the ESO MUSE pipeline (Weilbacher et al. 2020, version 1.6) to perform bias subtraction, flat fielding, twilight corrections, wavelength calibration, and corrections for telluric absorption. To improve the quality of the final datacube, the products of this stage of data reduction were post-processed using the CUBEXTRACTOR package (CUBEX hereafter, version 1.6, (see de Beer et al. 2023 and Cantalupo et al., in prep.). The subroutines CUBEF IX and CUBESHARP were employed to perform self-calibration and flux-conserving sky subtraction, respectively (see e.g., Borisova et al. 2016; Bacon et al. 2017; Marino et al. 2018; Cantalupo et al. 2019). Finally, the individual exposures were combined into a fully reduced datacube using the CUBECOMBINE tool. We calibrated the astrometry of the datacube using stars from Gaia Collaboration (2023). In the final datacubes, we obtain a Lyα 1σ surface brightness limit of ≈ 1.5 × 10⁻¹⁹ erg s⁻¹ cm⁻² arcsec⁻² (corresponding to ≈ 4.2 × 10⁻²² erg s⁻¹ cm⁻² pixels⁻²) in a single layer of 1.25 Å within circular apertures of radius r = 1″ (corresponding to 5 pixels).

We also performed MUSE NFM-AO-N observations of the central quasar (PID 109.232M, PI: A. Travascio). The key steps of the data reduction are the same as those detailed for the WFM observations above, except for the treatment of the sky. While in the case of WFM observations, the sky frame was generated directly from the OBJECT exposures, for the NFM data we created the sky frame from the pixel tables of additional exposures taken in empty regions of the sky. The observing sequence consists of 250 s of sky exposures every 2300 s of OBJECT ones. The reduced and combined data result in an integrated on-source exposure time of 2.3 h, spanning a 7.4″ × 7.4″ FoV with a sampling of 0.025″ per pixel. The spectral coverage is the same as for the WFM observations, except for the range 5800−6000 Å which is masked due to the contamination of the NaD lines of the notch filter required for the laser guide star in AO mode. Due to the absence of bright stars in the small FoV, we calibrated the astrometry of the combined datacube by aligning the centroid of the QSO emission to the ALMA detection of its CO(4–3) emission line (described in Sect. 3.1). For the NFM data, we obtain in circular apertures of radius r = 0.125″ (corresponding to 5 pixels) a Lyα 1σ SB limit in a single layer of ≈ 1.4 × 10⁻¹⁷ erg s⁻¹ cm⁻² arcsec⁻² (corresponding to ≈8.7 × 10⁻²¹ erg s⁻¹ cm⁻² pixels⁻²). Finally, we characterized the point spread function (PSF) of the MUSE observations as described in Sect. 4.1. The effective PSF full width at half maximum (FWHM) estimated from the central quasar is $\approx 0\stackrel{″}{.}9$ and $\approx 0\stackrel{″}{.}12$ in the WFM and NFM data, respectively.

2.2. ALMA observations

We employ data from ALMA Cycle 11 program 2024.1.00499.S (PI: A. Pensabene). The observations used ALMA band 3 and were designed to target the rotational transition of carbon monoxide CO(4–3) around the quasar redshift of z ≈ 3.66 (rest-frame frequency 461.041 GHz) as well as the underlying 3 mm dust continuum (rest-frame wavelengths 650 μm) with a medium angular resolution of $\approx 0\stackrel{″}{.}35$ in the quasar field. The half power beam width of the 12 m ALMA antenna is $58\stackrel{″}{.}8$ at the reference frequency of 99.00 GHz. The frequency setup consists of two adjacent 4-bit 0.938 GHz-wide spectral windows (SPWs) in the lower-side band (LSB) covering the CO(4–3) line within ≈ ± 2780 km s⁻¹ around the quasar systemic redshift, and two 2-bit 1.875 GHz-wide SPWs placed in the upper-side band (USB) for continuum detection. The native spectral resolution of the acquired data is ≊3.9 MHz channel⁻¹ (≊ 11.8 km s⁻¹ channel⁻¹ at 99.00 GHz). The observations were taken in five execution blocks during the period 2-5 April 2025 employing a total on-source exposure time of ≈4.07 h, and utilizing a number of ALMA main array antennas ranging from 42 to 47 with baselines within 15.1 − 2516.9 m yielding a naturally weighted synthesized beam FWHM of $0\stackrel{″}{.}54\times 0\stackrel{″}{.}40$. The weather conditions during the executions provided a mean precipitable water vapor ranging between 0.9 and 1.9 mm achieving a nominal sensitivity of 48 μJy beam⁻¹ over a bandwidth of 33.190 MHz (≊100 km s⁻¹) for line observations and 6 μJy beam⁻¹ over the aggregate continuum bandwidth of 4.363 GHz.

We calibrated the data using the Common Astronomy Software Application (CASA, McMullin et al. 2007; Hunter et al. 2023) by running the calibration pipeline scriptForPI provided along with the measurement using CASA version 6.6.1. No self-calibration was performed on the dataset. We imaged the ALMA band 3 visibilities by using the CASA tclean task. We first obtained a “dirty” datacube (specmode=“cube”) with 25 km s⁻¹ channel width, and a continuum image (specmode=“mfs”) by Fourier transforming raw visibilities (niter = 0) adopting a natural weighting scheme to maximize the sensitivity per beam. We set the pixel scale to $0\stackrel{″}{.}08$ to achieve Nyquist sampling of the beam minor axis with at least 5 pixels. Additionally, in order to improve the surface brightness sensitivity for faint extended emission, we obtained a “dirty” continuum image and cube applying a Gaussian taper in the uv plane with FWHM_uv ≊ 182 kλ and setting the pixel scale to $0\stackrel{″}{.}2$, yielding a synthesized beam FWHM of $1\stackrel{″}{.}36\times 1\stackrel{″}{.}27$. Finally, we obtained “cleaned” data by running the tclean task down to a signal-to-noise threshold of nsigma = 1.0 placing masks around the observable emission of the sources (see Sect. 3).

2.3. Ancillary archival datasets

We complement our observations with archival imaging and spectroscopic data of the quasar and its surrounding field, as described below.

2.3.1. HST imaging

The Hubble Space Telescope (HST) observed the field around the Q0055-269 on 2024 June 25 (program ID: 17483; PI: R. Dutta) using the F160W filter on the Wide Field Camera 3 (WFC3) infrared (IR) detector with a total integration time of 1809 s split into three frames with dithering offsets as in the WFC3-IR-DITHER-LINE-3PT template. The data reduction was performed using the standard HST pipeline, in particular the DrizzlePac version 3.9.1 (Hoffmann et al. 2021; Fruchter et al. 2010). The pixel scale of the final image is $0\stackrel{″}{.}128{\mathrm{pix}}^{-1}$. At the redshift of the structure, z ≈ 3.66, the pivot wavelength of the filter (λ_ref = 15370 Å) samples the rest-frame UV (329.8 nm) emission. The data described here can be obtained from the MAST archive at². Finally, we obtained a background-subtracted HST/WFC3 F160W image by using photutils Python package (Bradley et al. 2023) with sigma-clipping statistics to remove the sources, and SExtractionBackground as an estimator of the background. The astrometry of the final product is registered to Gaia Collaboration (2023).

2.3.2. UVES spectrum of the quasar

High-resolution spectroscopy of the quasar Q0055-269 taken with the Ultraviolet and Visual Echelle Spectrograph (UVES, Dekker et al. 2000) at the VLT was already publicly available³ from the first data release of the Spectral Quasar Absorption Database (SQUAD) project (Murphy et al. 2019). We retrieved the coadded and continuum-normalized UVES spectrum of the quasar. This consists of a combination of 24 exposures (programme IDs 092.A-0011(A) and 65.O-0296(A), PI: J. Schaye, S. D’Odorico) resulting in a total exposure time t_exp = 89 245 s and covering the wavelength range 3139.73−9766.69 Å with a resolution of R ≈ 45 000. This data provides a complementary view of the quasar environment in absorption, which we later compare with the results obtained in emission from MUSE.

3. Identification of galaxies in the field

Leveraging the broad wavelength coverage of our data, we searched within the field to identify different galaxy populations and explore how they relate to the gas observed in emission.

3.1. ALMA-detected galaxies

We inspected ALMA datasets and identified the CO(4–3) line and 3 mm continuum emission of the quasar host galaxy. Additionally, we detected two other sources (hereafter named MQN04-QC and AG1) in the field showing emission lines at ≈ − 1170 km s⁻¹ and ≈ − 2200 km s⁻¹ relative to the quasar systemic redshift, and marginally detected in the 3 mm continuum. In particular, we note that MQN04-QC is detected at 3 mm with F_3mm = 19 ± 7 μJy, corresponding to a signal-to-noise ratio of S/N ≈ 2.5. The spatial coincidence of the continuum emission peak with that of the CO(4–3) line emission strengthens the reliability of this detection (see Fig 1).

Fig. 1. Spatial distribution of galaxies in the MQN04 field detected within |Δv ≤ 100 km s−1|, the ALMA-identified MQN04-QC galaxy at Δv ≈ −1172 km s−1 of the quasar systemic redshift, and the Lyα, C IV, and He II-emitting gas. Left panel: HST/F160W image. The star-forming galaxies detected in MUSE data are indicated by red and blue squares, depending on their line-of-sight velocity relative to the quasar’s systemic redshift as determined from the CO(4–3) line emission. The quasar host and the MQN04-QC companion galaxy are shown as circles. The contours indicate the detected extent of the Lyα (gold), C IV (cyan), and He II (magenta) nebulae (see Sect. 4). Central panel: QSO PSF-subtracted cutout of the HST/F160W image. The dashed lines mark the position of the quasar host and MQN04-QC galaxy as identified from ALMA observations. The region masked for the PSF modeling is also shown as a hatched circle (see Sect. 3.1). Right panel: ALMA CO(4–3) line-velocity integrated map and 3 mm continuum. The gray (CO(4–3) line emission) and black (3 mm continuum) contours correspond to [ − 2, 2.5n]σ (dotted and solid for the negative and positive levels, respectively), where n ≥ 1 is an integer, and σ is the RMS noise. The ALMA synthesized beam is shown as a yellow ellipse in the bottom left corner.

We searched for counterparts of these two ALMA-detected galaxies in the HST/WFC3 F160W image and the MUSE datacube. We detect stellar continuum and an emission line in the spectrum of AG1, which we identified as [O II] at z ≈ 1.314 using Marz tool (Hinton et al. 2016). Therefore, the emission line detected in ALMA corresponds to CO(2–1). MQN04-QC remains instead undetected in MUSE observations. We therefore searched for rest-frame UV emission in the HST/WFC3 F160W image at the location of MQN04-QC. To do so, we first subtracted the PSF of the QSO from the image. To model the PSF, we used one star in the field with similar observed flux to the quasar as a template, then performed a fit using MCMC ensemble sampler emcee using 30 walkers, and 500 steps discarding the first 100 as a burning phase. To optimize the subtraction of the PSF spikes and wings, we masked the quasar emission within the central $0\stackrel{″}{.}8$. In addition, we excluded outer regions beyond $2\stackrel{″}{.}8$ where the presence of additional sources may introduce biases in the fit. We used as free parameters the flux rescaling factor and the position of the PSF centroid, assuming uniform priors. The resulting PSF-subtracted image is shown in Fig. 1. We do not detect any counterparts of MQN04-QC or any additional sources⁴. The absence of any detected rest-frame UV emission associated with this galaxy can be explained by the presence of a large amount of dust obscuring its stellar light, as commonly observed in high-redshift submillimeter galaxies (SMGs, see, e.g., Casey et al. 2014). This may suggest that this source is a dust-obscured quasar companion at z ≈ 3.65. Given the lack of additional redshift information, the observed line may be interpreted either as a lower-J (J_up ≤ 3) or J_up ≥ 5 CO lines, the frequency of which enter in the LSB of the ALMA observations from redshift z ≲ 2.5 and z ≳ 4.8, respectively. However, based on the CO luminosity functions from Decarli et al. (2019b) and Boogaard et al. (2023), the expected number of galaxy interlopers in the volumes probed by the current ALMA observations is < 0.2. For the subsequent analysis, we therefore assume MQN04-QC to be located at z = 3.65.

We analyzed the CO(4–3) emission lines of the quasar host and MQN04-QC by extracting their spectra from the naturally weighted ALMA band 3 datacube. To this purpose, we applied the 2-σ photometry method (e.g., Béthermin et al. 2020; Pensabene et al. 2024) which enabled us to recover the line flux with an optimal S/N while preventing significant flux losses: (i) we extracted a spectrum from the source peak pixel as determined by the 3 mm dust continuum image, we corrected for the primary beam (PB) response, and fitted the spectrum with a Gaussian component for the line profile and a constant for the continuum. (ii) We obtained the line-velocity integrated map (0th moment) by summing together all the channels within ±2σ from the fitted line centroid of the initial spectrum. (iii) We extracted a new spectrum by summing the signal from all the corresponding pixels in the 0th moment map enclosed within the S/N ≥ 2 isophote. We converted the spectrum in flux density by dividing it by the synthesized beam area (in pixel units) and computed the associated uncertainties by multiplying the peak flux uncertainty (in mJy beam⁻¹) by the square root of the independent beam elements over the extracting region. (iv) We repeated the line fit on this new spectrum as described above, and we iterated the procedure until convergence. Hence, we performed an accurate fit of the resulting line profile by minimizing the residuals between the data and the model using the MCMC ensemble sampler emcee Python package with 100 walkers, 5000 steps, and discarding the first 1000 as a burning phase. We adopted Gaussian uncertainties in the definition of the likelihood and flat priors on the free parameters based on the last fitting iteration as described above. Finally, we derived the line luminosity as (e.g., Solomon et al. 1997): $L_{\mathrm{CO}}[L_{\odot }]=1.04\times {10}^{-3}{F}_{\mathrm{CO}}{\nu }_{\mathrm{obs}}{D}_{\mathrm{L}}^2$, $L{\prime }_{\mathrm{CO}}[{\mathrm{K}\mathrm{km}\mathrm{s}}^{-1}{\mathrm{pc}}^2]=3.25\times {10}^7{F}_{\mathrm{CO}}{D}_{\mathrm{L}}^2{(1+z)}^{-3}{\nu }_{\mathrm{obs}}^{-2}$, where F_CO is the velocity-integrated line flux in units of Jy km s⁻¹, ν_obs = ν_rest/(1 + z) is the observed central frequency of the line in GHz at the redshift z, and D_L(z) is the luminosity distance in Mpc. The two luminosity measurements are related via $L_{\mathrm{CO}}=3\times {10}^{-11}{\nu }_{\mathrm{rest}}^{3}L{\prime }_{\mathrm{CO}}$.

In addition to line measurements, we recovered the faint (observed-frame) 3 mm continuum emission of the sources by applying a single-beam photometry measuring the continuum source flux peak in the $\approx 1\stackrel{″}{.}3$ tapered image. In Fig. 1, we show the velocity-integrated CO(4–3) line maps of the quasar host and the foreground companion galaxy MQN04-QC, as well as the HST/WFC3 image in the F160W filter. In Fig. 2 we show the extracted CO line profile of the sources along with the best fit model. We report our measurements in Table 1.

Fig. 2. Signal-to-noise spectra extracted from the regions identified in the line-emitting nebulae. Top panels: Lyα (orange), C IV (blue), and He II (red) emission lines extracted from the MUSE-WFM data within the regions E (left), B (middle), and W (right). The Lyα S/N is rescaled by ×0.3 for visualization purposes. Bottom-left panels: He II emission line S/N spectra extracted from the MUSE-NFM data within the clumps C1 (top) and C2 (bottom). The shaded band marks spectral regions contaminated by skylines. Bottom-center panels: CO(4–3) emission line profiles of the quasar host (top) and the MQN04-QC galaxy (bottom). The data are shown in yellow, the single Gaussian best-fit curve is in black with 1σ uncertainty shaded in gray. The flux is in units of mJy. Bottom-right panel: H I absorption line detected in the (continuum-normalized) UVES spectrum of the quasar. The green line indicates the best-fit Voigt profile. In all the panels, the line-of-sight velocity is relative to the systemic redshift of the quasar host galaxy as determined from its CO(4–3) emission line. The green dashed line and shaded area indicate the centroid and the FWHM of the H I absorption line, respectively.

We estimated an order-of-magnitude indication of the dust and molecular gas mass of MQN04-QC and the quasar host galaxy from the observed 3-mm flux density and the CO(4–3) line luminosity. For this purpose, we modeled the dust emission with a single-temperature modified blackbody. In the case of MQN04-QC, we adopted a dust emissivity index of 1.8, and assumed a dust temperature of T_dust = 32 K, as observed in submillimeter galaxies (SMGs) at z ∼ 3 (e.g., Dudzevičiūtė et al. 2020). Under these hypotheses, we obtain a dust mass ranging within M_dust ∼ (0.8 − 1.6)×10⁸ M_⊙, taking into account the 3 mm flux density uncertainties. We then estimated the gas mass by assuming a CO(4–3)-to-CO(1–0) luminosity conversion factor of r₄₁ = 0.46 from Carilli & Walter (2013) and a CO(1-0) luminosity-to-H₂ molecular gas mass conversion factor of α_CO = 1.7 M_⊙ (k km s⁻¹ pc²)⁻¹ measured in a sample of local ultra-luminous infrared galaxies (Montoya Arroyave et al. 2023, see also Bolatto et al. 2013 for a review). With these assumptions, we obtain $M_{\mathrm{gas}}={\alpha }_{\mathrm{CO}}{r}_{41}^{-1}L{\prime }_{\mathrm{CO}}\sim 2\times {10}^{10}{M}_{\odot }$, and therefore a gas-to-dust ratio of δ_dg ∼ 95 − 190. These estimates are in broad agreement with typical properties of SMGs at z = 2 − 4 see, e.g., Swinbank et al. 2014; Dudzevičiūtė et al. 2020; Pantoni et al. 2021). With the same approach, we obtained estimates of the gas and dust mass for the quasar host galaxy. In this case, according to expectations for high-redshift quasar hosts (see, e.g., Downes & Solomon 1998; Beelen et al. 2006; Carilli & Walter 2013; Leipski et al. 2013, 2014), we assumed T_dust = 45 K, r₄₁ = 0.87, α_CO = 0.8 M_⊙ (K km s⁻¹ pc²)⁻¹. We obtain M_dust ∼ (1.3 − 2)×10⁸ M_⊙, M_gas ∼ 5 × 10⁹ M_⊙, and a gas-to-dust ratio of δ_dg ∼ 20 − 40. Finally, under the assumption of dispersion dominated systems, we derived a rough estimate of the virial dynamical mass enclosed within a radius R of the two galaxies using the equation M_dyn(< R)≈(3/2) σ² R/G (e.g., Decarli et al. 2018; Pensabene et al. 2025). With σ as the CO linewidth and R as half of the beam size $\approx 0\stackrel{″}{.}25\approx 1.8\mathrm{kpc}$, considered an upper limit for the unresolved galaxy size, we obtain M_dyn(< 1.8 kpc) = (3.1 ± 1.2)×10¹⁰ M_⊙ and (1.6 ± 0.3)×10¹⁰ M_⊙ for the quasar host galaxy and MQN04-QC, respectively. These values are in agreement with the gas masses inferred from the CO line luminosity. In particular, the gas mass of MQN04-QC appears to be comparable to its dynamical mass, therefore suggesting that this galaxy is very rich in molecular gas. However, we stress that the evaluation of the various mass budgets provided above relies on assumptions that cannot be verified with the limited available information, and therefore should only be taken as order-of-magnitude estimates.

3.2. MUSE continuum-detected galaxies

Additionally, we performed a blind search for continuum-detected star-forming galaxies in the MUSE-WFM datacube following the approach described by Galbiati et al. (2025). Briefly, we collapsed the datacube along the wavelength axis, producing a white-light image used for the detection. We then extracted the continuum-emitting sources by running SExtractor (Bertin & Arnouts 1996, version 2.24.2) with THRESHOLD = 1.5 weighted on the local RMS background level to make it uniform across the FoV, using MINAREA = 3 pixels, and default deblending parameters. The spectra of each source were then extracted within an aperture of radius r = 2.5 × R_Kron, where R_Kron is the Kron radius, chosen to collect the total flux of the sources. Finally, we determined the spectroscopic redshifts using the Marz tool⁵ based on Lyαλ1216Å emission line and, if visible, looking for Si IIλ1260, λ1526Å, O I+Si IIλ1303Å, C IIλ1334Å, Si IVλ1393, 1402Å, and C IVλ1548, 1550Å. With this procedure, we identify five galaxies within |Δv_QSO ≤ 1000 km s⁻¹| relative to the central quasar⁶ and with r-band magnitudes m_r ≲ 27 mag (see Table A.1). In Fig. 1, we show the location of all the ALMA- and MUSE-selected galaxies, together with the maximum extent of the detected nebular emission (see Sect. 4). We note that additional analysis is required to obtain a complete census of galaxies in the field, since emission line dominated galaxies, such as fainter Lyα emitters (LAEs), could be missed by our current search for continuum-bright sources. However, the presence of bright and extended Lyα emission in the FoV prevents an unambiguous identification of LAEs around the quasar in the absence of continuum-detected counterparts and spectral features other than their Lyα line.

4. Detection of Lyα, C IV, and He II extended emission

4.1. The gas in emission: Extraction of the nebulae

We searched the MUSE datacube for extended Lyαλ1216Å, C IVλ1548Å, and He IIλ1640Å emission around the central quasar in the WFM observations as well as, with higher resolution, for the central part of the extended He II emission in the NFM ones⁷. To this end, we applied additional processing steps using CUBEX, following the approach adopted in previous studies of quasar nebulae at similar redshifts (e.g., Borisova et al. 2016; Arrigoni Battaia et al. 2019; Cantalupo et al. 2019; Fossati et al. 2021). The first step consists of subtracting the quasar PSF with the CUBEPSFSUB routine in CUBEX. This tool constructs a narrowband image centered on the quasar. From this, a PSF image is generated and rescaled to match the quasar specific flux at each channel within the central 1 × 1 arcsec² and 0.25 × 0.25 arcsec² region for the WFM and NFM observations, respectively. To mitigate the effect of artifacts and cosmic rays, we computed the rescaling factor using an iterative sigma-clipped mean with a symmetric ±3σ rejection threshold around the median, iterated until convergence. Then, a circular cutout of radius 25 and 50 pixels for the WFM and NFM data, respectively (i.e., 5″ and 1.25″, corresponding to ≈5× and ≈10× the PSF FWHM) was extracted from the rescaled PSF image and subtracted from the datacube. This procedure was performed for each wavelength layer of the cube. We first ran CUBEPSFSUB, setting the spectral filter width (zPSFsize) to 600 layers to produce a PSF-subtracted datacube that we visually inspected to identify the spectral regions containing the line emission of the nebulae. We then ran the routine on the original datacube applying masks to such regions and using zPSFsize = 150, corresponding to 187.5 Å. We estimated the FWHM of the PSF of the central quasar by computing a mean radial profile, obtaining $\approx 0\stackrel{″}{.}9$ and $\approx 0\stackrel{″}{.}12$ in the WFM and NFM, respectively. After the PSF subtraction, we subtracted the background and any continuum emission using the CUBEBKGSUB tool (see Borisova et al. 2016). Briefly, the spectrum of each spaxel was rebinned into spectral bins and median-filtered to estimate the continuum, which was then subtracted from the data.

We then produced optimally extracted images of the Lyα, C IV, and He II emission around the quasar from the WFM data and of the He II emission from the NFM ones. The steps we performed are the following: (i) we smoothed each layer of the data and of the variance cubes using a two-dimensional Gaussian kernel with standard deviations σ_x = σ_y = 1.5 pixel (corresponding to $0\stackrel{″}{.}3$ and $0\stackrel{″}{.}04$ for the WFM and NFM data, respectively). (ii) We searched for all connected voxels above a signal-to-noise threshold of S/N = 3.5 using a friends-of-friends algorithm with a connectivity of 26. We set the minimum size of a structure to be extracted to 2000 voxels. (iii) We collapsed the datacube by adding together all the selected voxels along the spectral axes. We note that the results are robust against small variations in the adopted parameters. In Fig. 3 we report the results.

Fig. 3. Optimally extracted images of the extended Lyα λ1215 Å (top-left panel), C IV λ1548 Å (top-right panel), and He II λ1460 Å (MUSE-WFM and NFM data, bottom-left and bottom-right panels, respectively) nebulae detected around the quasar Q0055-269 in the MQN04 field. The contours are matched in flux density and indicate 2n × s0 SB levels with n ≥ 1 integer, and s0 = 0.5 × 10−18 and 16 × 10−18 erg s−1 cm−2 arcsec−2 for the WFM and NFM, respectively. The layer of the datacube at the quasar redshift has been added as a reference noise background for each transition. The green, magenta, and blue ellipses delineate the selected E (east), B (bridge), and W (west) regions identified based on the He II emission detected in the MUSE-WFM data, while the white and orange ellipses enclose the two clumps (C1 and C2, respectively) revealed by the MUSE-NFM observations. The gray circle indicates the location of the quasar host as determined by ALMA observations, and the size is rescaled to match the PSF FWHM. The blue circle marks the position of MQN04-QC. The contours of the He II nebula detected in the MUSE-WFM data (white), the quasar CO(4–3) emission (green, 2n × RMS, with n ≥ 1), as well as the ALMA synthesized beam, are overlaid on the MUSE-NFM data in the bottom-right panel.

In the WFM datacube we detect bright Lyα, C IV, and He II emission with SB ≥ 10^−18.3 erg s⁻¹ cm⁻² arcsec⁻². The Lyα emission is the most extended and appears elongated on scales of ≈130 kpc (or ≈600 comoving kpc; ckpc) while the C IV and He II nebulae have sizes of ≈55 and ≈45 kpc (or ≈250, and 210 ckpc), respectively, as measured at surface brightness levels of SB = 10^−18.3 erg s⁻¹ cm⁻² arcsec⁻². In the NFM data, we find bright and clumpy He II emission extending eastward up to projected distances of ≈ 3.7 kpc (or ≈ 20 ckpc) from the quasar, which was entirely hidden under the quasar PSF in the WFM observations.

4.2. The gas in absorption

To obtain a more complete view of the gaseous environment around Q0055-269, we complemented the emission-line analysis with a search for absorption features in the quasar spectrum obtained with UVES at the redshift of MQN04. As a first step, we converted the wavelength of the UVES spectrum of the quasar from air to vacuum and then visually inspected the continuum-normalized spectrum searching for H I absorption as well as proximate C IV and Si IV absorption-line systems. As previously reported by Murphy et al. (2019), no damped Lyα (DLA) absorbers are present along the quasar sightline. However, we identified an H I absorption system on the blue side of the quasar Lyα emission showing Lyβλ1025Å and Lyαλ1216Å. We do not detect any metal lines associated with this absorber. We fitted Voigt profiles to the Lyβ and Lyα features simultaneously to mitigate the uncertainties arising due to the saturation of the latter. To do so, we used the Bayesian code MC-ALF (Longobardi et al. 2023), following the approach of Lofthouse et al. (2023) and Galbiati et al. (2023). The free parameters are (i) Doppler parameter (b); (ii) column density (N_H I); (iii) redshift; (iv) a multiplicative constant ranging from 0.98 to 1.02 to account for small uncertainties in the continuum normalization; (v) spectral resolution, allowed to vary within ± 1 km s⁻¹ around the nominal value (≈ 7 km s⁻¹ on average). From this analysis, we determine that the H I absorber lies at z ≈ 3.654 (corresponding to Δv_QSO ≈ −590 km s⁻¹ relative to the quasar along the line of sight), has a column density of log(N_H I/cm⁻²) = 14.64 ± 0.10 and a Doppler parameter of b = 25.05 ± 0.58 km s⁻¹. The H I column density is much lower than that expected in the self-shielding limit (log(N_H I/cm⁻²) ≈ 17.2, see, e.g., Rahmati et al. 2013), which indicates that the detected absorber is associated with a gaseous medium that is highly ionized. In Fig. 2, we report the H I absorber identified in the quasar spectrum along with the best-fit model. We note that additional absorption features are present next to the H I absorber. However, the lack of additional lines of the Lyman series associated with these features prevents us from unambiguously identifying them as H I.

4.3. 3D multiphase gaseous structure around the quasar

While Lyα and C IV are the brightest UV/optical emission lines often detected in nebulae around quasars, the nonresonant nature of the He II line emission provides us with a key tracer of the intrinsic properties of the extended, diffuse, and ionized gas around the quasar. We produced two-dimensional surface brightness maps by collapsing the datacube around the Lyα, C IV, and He II emission lines as described in Sect. 4.1. In the maps obtained from the MUSE-WFM observations, the He II emission appears to be distributed in two structures, respectively located on the east (E) and west (W) sides with respect to the quasar location, as well as one diffuse emission “bridge” connecting these two regions (B), as shown in Fig. 3. Similar morphological features can also be recognized in the C IV nebula, and tentatively in the Lyα map. However, the resonant scattering of the Lyα photons severely affects the observed spatial distribution of the extended emission, making these structures less apparent. Similarly, we extracted an optimal image of the He II emission from the MUSE-NFM data, which provides a view of the inner ≲5 kpc around the quasar. We identified two clumps (C1 and C2) that are located near, but with a measurable displacement, the peak of the He II emission in the E region. This spatial separation between E and C1/C2 suggests that the emission on large scales is dominated by a diffuse gas component rather than the clumpy structures resolved on kpc scales in the NFM data.

We quantified the velocity offsets among the various regions and characterized the profile of the emission lines by extracting the spectra from the elliptical apertures that encompass the main morphological features identified in the He II nebulae and represented in Fig. 3. In Fig. 2 we show the S/N spectra of the Lyα, C IV and He II line emission from the E, W, and B regions as well as those of the He II emission extracted from the C1 and C2 clumps. The emission lines are robustly detected at S/N ≳ 5 − 8 in all these regions, whereas region E suffers from increased noise as a result of its proximity to the quasar PSF. We derived the characteristics of the He II emission in each region identified in the WFM observations by performing a fit of the line profile with a single Gaussian using emcee with 50 walkers, and 2500 steps, discarding the first 250 as a burning phase. We report the results in Table 2. The He II profiles exhibit a similar line width of FWHM ≈ 500 km s⁻¹ throughout all three emitting regions. In addition, the view provided by the observation of the gas in emission is complemented by the H I detected in absorption in the quasar sightline (see Sect. 4.2). We found that the H I absorber lies ≲100 km s⁻¹ from the W emitting region. However, we note that the low-column density of this absorber suggests that the observed H I gas is likely arising from the IGM rather than tracing local neutral gas coexisting with the ionized emission. Therefore, the redshift alignment with the W region could plausibly arise by chance. In Sect. 5, we discuss how different scenarios can account for the spatial and kinematical structures revealed by the observed He II emission.

The two bright clumps identified in the NFM observations provide us with additional information on the E region and its connection with the immediate surroundings of the quasar. The He II emission lines extracted from C1 and C2 are broad (FWHM ≳ 500 km s⁻¹) and show a double peaked profile with a minimum within ≲50 km s⁻¹ of the redshift of the quasar. When comparing these He II line profiles with that of the E region extracted from the WFM data, we find that the line centroid is shifted by ≈ − 200 km s⁻¹ with respect to those measured in C1 and C2. Interestingly, the He II line luminosity in E exceeds by ≈60% the total luminosity in the combined C1+C2 regions (see Table 2). Part of the observed difference may be driven by the different PSFs and surface brightness sensitivities of the WFM and NFM observations (see Sect. 2.1), which affect the apparent flux distribution. However, if the excess flux observed in region E was produced by additional compact clumps, these would be expected to appear in the NFM data with luminosities comparable to those of C1 and C2. The absence of such compact structures in the NFM observations therefore suggests that a significant fraction of the He II emission in E arises from a diffuse gas component on scales larger than $\approx 0\stackrel{″}{.}3$ (or ≈2.2 kpc at z = 3.66).

We then used the nonresonant He II line emission detected in the WFM and NFM data to investigate the 3D structures and kinematics of the gas surrounding the quasar. A kinematical analysis of the MQN04 nebulae was previously performed by Guo et al. (2020), although limited by the lack of information on the quasar systemic redshift. In the present work, thanks to our ALMA CO(4–3) line observations (see Sect. 3) and the high-resolution MUSE-NFM data, we can now accurately study, for the first time, the observed CGM gas kinematics down to ≈1 kpc relative to the location of the central quasar host galaxy. In Fig. 4 we show the channel maps of the He II line emission, the velocity shift relative to the quasar systemic redshift (1st moment), as well as a position-velocity diagram (PVD) extracted from the S/N cube along a polygonal path of 1″ width covering the brightest emitting regions. The two-dimensional maps of the velocity shift and velocity dispersion (2nd moment) are instead shown in Fig. B.1 and described in Appendix B. The E and W regions are clearly identifiable not only in the SB maps but also in velocity space, as well as the diffuse gas bridge, labeled B, connecting them. Remarkably, in all these regions, the He II emission appears to be blueshifted with respect to the quasar systemic redshift. Specifically, the He II arising from the W region spans velocities ranging in the interval ≈[−800, −400] km s⁻¹, with the peak of the observed emission located at sky-projected distances of ≈18 kpc. The range of observed line-of-sight velocities in this region is also consistent with that of the H I absorber detected along the quasar sightline. This could be indicative that the emitting gas observed in W is part of a more extended structure residing in the foreground of the quasar. On the other hand, in the E region, the He II appears to be closer to the quasar both in projection (5 kpc) and in velocity (≈[−400, 0] km s⁻¹). This may suggest that at least part of the observed emission arises from the interstellar medium (ISM) of the quasar host galaxy. However, the 1st moment map reveals that the W and E regions appear to be kinematically linked by a smooth velocity gradient. These conclusions are further supported by the inspection of the PVD and the channel maps. Altogether, the combined spatial and kinematical properties of the system suggest that the observed emission traces a physically connected structure rather than multiple unrelated structures. Finally, we note the presence of a low-S/N emission extending further eastward and aligned in redshift with the E region, which could instead represent an additional distinct component.

Fig. 4. Kinematic analysis of the He II extended emission extracted from the MUSE-WFM data. Left panels: channel maps of the He II emission obtained from the PSF and continuum-subtracted datacube by collapsing the voxels within windows of 100 km s−1. The contours refer to He II SB levels as in Fig. 3. Central panel: line-of-sight velocity map (1st moment) of the He II emission. The position of the quasar is indicated by the gray circle with a size corresponding to the QSO PSF FWHM. Right panel: position-velocity diagram extracted along the polygonal slit shown in the central panel. The boundaries of the E, B, and W region crossed by the slit are reported with vertical lines. The dotted and solid contours indicate S/N = 4 and 7 levels, respectively. The spectral and spatial resolution is represented as a rectangle in the bottom-left corner.

A similar analysis conducted on lines other than the He II reveals that the gas kinematics traced by the Lyα and C IV emission is more complex and does not show a clear velocity gradient across the nebulae. However, except for the E region, the overall velocity shifts are consistent with those measured from the He II line (see Fig. 2). This indicates that in B and W regions radiative transfer effects seem to only moderately affect the ability of resonant lines to trace the gas kinematics, as expected in a highly ionized medium.

4.4. Spatially resolved line ratios

The simultaneous detections of Lyα, C IV, and He II emission lines provide powerful diagnostics to get insights into the origin and physical properties of the extended gas nebulae surrounding the quasar such as the ionization parameter, gas density and metallicity (e.g., Villar-Martin et al. 1997; Arrigoni Battaia et al. 2015b,a; Prescott et al. 2015; Cantalupo et al. 2019; Marino et al. 2019; Jiménez-Andrade et al. 2023). Thanks to the exquisite S/N of the observations presented in this work, we can obtain spatially resolved line ratios in the MQN04 nebula. To this purpose, we maximized the recovery of faint He II and C IV emission by applying a resampling scheme along the spatial axes of the smoothed datacube using the adaptive two-dimensional Voronoi binning technique (Cappellari & Copin 2003). In this approach, we required each cell to achieve S/N = 8, 9 thresholds, for the He II and C IV, respectively, chosen to preserve the native sampling in the regions of interest while maximizing the detection of faint emission in low S/N pixels. To obtain consistent and optimal line ratio maps, we resampled the datacubes using the Voronoi binning obtained for the least extended nebula of the ratios. We report the results in Fig. 5 and in Table 2.

Fig. 5. Line ratios in the MQN04 nebula: C IV/Lyα (left panel), C IV/He II (central panel), and He II/Lyα (right panel). The line ratios are computed in Voronoi-resampled datacubes to maximize the signal-to-noise of the faintest regions. The gray contour in each panel represents the SB mask of the optimally extracted image, as in Fig. 3, of the least extended nebula included in the line ratio. The colored ellipses mark the same regions identified in Fig. 3. The gray circle indicates the location of the quasar host as determined by ALMA observations and is sized as the PSF FWHM.

Interestingly, all of the line ratios are significantly different between the E and W regions, indicating that different physical properties and/or ionization conditions are in place. In particular, following Cantalupo et al. (2019), assuming that the observed emission is powered by photoionization and the medium is characterized by a single density distribution, the observed difference between the He II/Lyα ratios indicates that the E and W regions are located at different physical distances. Under the same hypotheses, the He II emission in the E region is enhanced due to its proximity to the quasar ionization source. Furthermore, the high C IV/Lyα and C IV/He II line ratio observed in W region may indicate significant metal enrichment of the gas (e.g., Arrigoni Battaia et al. 2015a; Marques-Chaves et al. 2019; Kolwa et al. 2019). However, the determination of the metallicity of this region is challenging given the lack of any metal absorption features at the redshift of the H I aligned with W in the UVES spectrum of the quasar. In Sect. 5, we discuss various scenarios to interpret these results.

4.5. Overdensity of SF galaxies and large-scale environment

Direct observations of the cosmic gas in emission, particularly in regions surrounding quasars (e.g., Cai et al. 2018; Tornotti et al. 2025b) and hosting dense concentrations of star-forming galaxies (Cai et al. 2016; Umehata et al. 2019; Pensabene et al. 2024; Banerjee et al. 2025; Galbiati et al. 2025; Tornotti et al. 2025a), offer powerful means to probe its large-scale distribution and explore its connection to galaxy evolution. We inspected the MUSE data to conduct a census of the population of star-forming galaxies embedded in the gaseous structures observed in emission, as described in Sect. 3. We identified 5 galaxies with spectroscopic redshifts within |Δv_QSO| ≲ 500 km s⁻¹ (median Δv_QSO ≈ −68 km s⁻¹) from the systemic redshift of the central quasar. As shown in Fig. 1, three of these galaxies are concentrated within ≈ 500 ckpc and embedded in the Lyα nebula. The remaining two, although lying at larger projected distances from the quasar, appear to be aligned along the direction in which the nebula extends. This suggests that filaments of gas may continue further along this axis and may not be detected in emission in the current data due to their lower density or because they are not efficiently illuminated by the radiation of the quasar or additional undetected AGNs.

We are interested in comparing this system with other known overdensities found at similar redshifts and proven to be connected by filaments of gas observed in Lyα emission. We compared MQN04 field with two other structures for which MUSE data with similar exposure times are available, namely the protoclusters SSA22 at z ≈ 3.08 (Umehata et al. 2019) and MQN01 at z ≈ 3.25 (Pensabene et al. 2024; Galbiati et al. 2025). For these fields, we repeated the same procedure described in Sect. 3.2 and selected all the star-forming galaxies with R-band magnitudes m_R ≤ 27 mag, which is the completeness limit of the samples based on the luminosity functions in Galbiati et al. (2025). We then measured the number density of these galaxies within velocities |Δv| ≤ 1000 km s⁻¹ and projected separation R ≤ 20″ (defined to match the FoV of MQN04) from the central quasar in MQN01 field. In the case of SSA22, since there are no bright quasars in the structure, we took as the fiducial value of the number density the median of the distribution obtained by bootstrapping over the position and redshift of the individual galaxies, taken one by one as the center of the structure. We chose the same velocity interval as MQN01, i.e., |Δv| ≤ 1000 km s⁻¹, and projected separation R ≤ 2′, calibrated to encompass at least five neighbors since the galaxies are more sparse in the FoV. As a control field sample, we took the one described by Galbiati et al. (2025) assembled from fields where quasars are not present at z ≈ 3 that were observed with MUSE and have depth similar to MQN04 data.

As a result, we measure the following comoving number densities: ρ_MQN04 = (0.041 ± 0.002) cMpc⁻³, ρ_MQN01 = (0.056 ± 0.007) cMpc⁻³, ρ_SSA22 = (0.015 ± 0.002) cMpc⁻³, and ρ_field = (0.0011 ± 0.0001) cMpc⁻³. This implies overdensities relative to the field of δ_MQN04 = 41 ± 4, δ_SSA22 = 15 ± 2, and δ_MQN01 = 56 ± 9. We note, however, that the larger projected radius adopted for SSA22 results in a significantly larger surveyed volume compared to the other structures. This difference may dilute the measured galaxy number density if the structure is centrally concentrated, and therefore the overdensities should be compared with this caveat in mind. We found that MQN04 is substantially more overdense than SSA22, though not as extreme as MQN01. In MQN04 field, the galaxies appear spatially concentrated around the quasar and lie within small velocity separations. However, none show AGN-like emission lines (e.g., N Vλ1240Å, C IVλ1548, 1550Å, He IIλ1640Å). Assuming the same AGN fraction observed in MQN01 for galaxies with M_★ ≳ 10^9.5 M_⊙ (f_AGN ≈ 0.4; Travascio et al. 2025), we would expect roughly two AGN in MQN04, yet none are detected. Furthermore, the only SMG detected with ALMA, MQN04-QC, lies relatively far from the quasar in both redshift and projected separation. This contrasts with MQN01, where a significant AGN overdensity is clustered within an area of 2′×2′ around the central quasar (Travascio et al. 2025), and a massive submillimeter companion lies only Δv ≈ − 300 km s⁻¹ and 10 kpc from the quasar (Pensabene et al. 2025). Interestingly, none of these features is present in SSA22.

Overall, the three structures appear to trace different environments and may occupy distinct locations within the cosmic web. SSA22 is relatively diffuse, MQN04 is denser but still AGN-poor, and MQN01 is highly overdense and AGN-rich, with the latter showing the strongest indication of tracing a massive node (Cantalupo et al., in prep.). Alternatively, the differences among them may reflect observations of structures at different evolutionary stages (e.g., Shimakawa et al. 2018).

5. Discussion

The extended gas directly imaged via multiple emission lines in the MQN04 field provides powerful insights into the complex and multiphase structure of the CGM gas. In particular, the morphology of the nonresonant He II emission suggests the presence of multiple components that can be spatially and kinematically connected. In what follows, we address the plausibility of three different scenarios as possible explanations for the observed characteristics of the diffuse gas.

5.1. Scenario A: Inflow of a gas stream

We first investigate the possibility that the observed He II line emission arises from inflows of cool gas streams that accrete onto the quasar host galaxy from intergalactic scales, as predicted by hydrodynamical simulations (e.g., Fardal et al. 2001; Kereš et al. 2005; Dekel et al. 2009; van de Voort et al. 2012; Waterval et al. 2025). Since a large portion of the gas within these streams can exhibit high angular momentum, the accreting gas may appear as an inspiraling structure (e.g., Danovich et al. 2015; Stewart et al. 2017; Arrigoni Battaia et al. 2018; Wang et al. 2022; Zhang et al. 2023a,b) and produce velocity gradients as observed in MQN04 (see Fig. 4).

In this framework, the velocity structure of the gas can be interpreted as follows. The blueshifted region W would be a stream of gas infalling toward the quasar host from the far side, while the lower line-of-sight velocity in the bridge B would be due to a change in direction as a result of the inspiraling motion. Finally, region E and clumps C1/C2, close to the QSO systemic redshift, can represent the last leg of the accretion path, where the inflowing gas merges with the ISM of the central galaxy. We note that diffuse blueshifted emission further eastward of region E cannot be explained by the same stream originating in region W and should rather be attributed to a separate accretion flow. Within this scenario, the absence of any extended redshifted emission would suggest that neither stream makes a full turn around the galaxy, indicative either of relatively low angular momentum, or that some angular momentum has been lost during the inflow, possibly due to interactions with the hot phase of the CGM (see, e.g., the semi-analytical models from Afruni et al. 2021; Lan & Mo 2019). In the scenario described above, the H I absorption could be seen as a Lyα forest absorber that is physically unrelated to the emitting gas and appears aligned with W in velocity just for coincidence.

Instead, if the emitting gas in W is physically associated with the H I absorber, as similar velocities suggest, then it must be located on the near side. The observed blueshift can be reconciled with an inflow scenario only if the gas extends over intergalactic scales where the Hubble flow dominates its relative motion. In this case, assuming that the inflow velocity of the gas is small compared to the Hubble flow, the inflow should be a contiguous filament of gas extending at least ≈ 1.6 Mpc mainly aligned along the line-of-sight. However, the observed high surface brightness of the He II emission cannot be explained at such large distances from the quasar without the presence of additional multiple AGN with a hard ionizing spectrum (e.g., Humphrey et al. 2019; Pezzulli & Cantalupo 2019). While a possible source of such ionizing radiation may be identified in MQN04-QC, assuming it hosts an obscured AGN, it still lies too far away both in redshift space and in projection from the He II emitting gas. Also, no evidence for additional AGN appears in the data. This scenario would also require the presence of some redshifted emission arising from distances smaller than the turnaround radius (≈4 R_vir, e.g., de Beer et al. 2023), where the inflow velocity dominates over the Hubble flow, which is not observed. Although several aspects of these scenarios cannot be confirmed with the existing observations, and despite the limitations discussed above, the inflow interpretation remains a plausible explanation for the observed properties.

5.2. Scenario B: Multiphase ionized outflow

Another interpretation of blueshifted emission and H I absorption could be an outflow. A standard axisymmetric biconical outflow scenario can be ruled out, unless the redshifted components are completely obscured by dust (e.g., Villar-Martín et al. 2011). However, there is no indication from our deep ALMA observations of the presence of dust on CGM scales. An outflow with more complex geometry (non-biconical) cannot be excluded. We note that the spatial variations in the He II/Lyα line ratios (see Sect. 4.4) could indicate that E and W (as well as B) regions are located at different physical distances from the quasar. This pattern could arise if outflows propagate through media of unequal resistance, caused either by the dynamical pressure of denser gas, such as cold inflows, or by turbulent viscosity (e.g., Gabor & Bournaud 2014). Alternatively, the observed regions might trace a single outflow whose cones have varying orientations or were launched at different epochs. In addition, a hint of a blueshifted excess in the CO(4–3) emission line profile of the quasar host at ≈ − 300 km s⁻¹ (see Fig. 2) might be associated with the molecular phase of an outflow. However, we do not observe CO(4–3) emission co-spatial with the extended He II emission. This might suggest that molecules are dissociated within the outflow on large scales. Deep ALMA data combined with IFS observations with JWST/NIRSpec targeting Hα emission to map the ionized gas at high resolution may help further explore this scenario.

5.3. Scenario C: Stripping from a galaxy encounter

Alternatively to the previous scenarios, we consider the possibility of tidal stripping due to a close encounter between a gas-rich galaxy and the quasar host. To test this scenario, we compute a lower limit on the quasar halo mass assuming that the close encounter occurs with the following simplifying conditions: (i) the observed velocity of the gas in W, the brightest and extended region, corresponds to the galaxy velocity; (ii) the observed line-of-sight velocity of the gas in W (≈ − 600 km s⁻¹) is entirely due to its peculiar motion relative to the quasar systemic redshift; (iii) the gas is moving with velocity fully aligned with the line of sight, and (iv) it is indicative of the velocity of the galaxy at the pericenter of its orbit. By equating the observed velocity of the galaxy with the escape velocity⁸, we obtain that a minimum halo mass of $M_{\mathrm{halo},\min }^{200}>8\times {10}^{11}{\mathrm{M}}_{\odot }$ is required to accelerate a test mass up to 600 km s⁻¹ at the pericenter distance of ≈18.7 kpc (see Table 2). The estimated lower limit is consistent with the typical halo mass of quasars at cosmic noon (e.g., Trainor & Steidel 2012; Eftekharzadeh et al. 2015; de Beer et al. 2023; Pizzati et al. 2024). To complete the picture, the region B and the H I absorber can trace different portions of the tidal tail arising from the quasar-galaxy interaction, while E and the region extending further eastward can represent gas stripped during a previous close passage.

Alternatively, the large amount of diffuse, extended gas detected around the quasar can also be explained by the presence of a hot medium, such as a proto-intracluster medium (ICM) (e.g., Di Mascolo et al. 2023) that is efficient in stripping the gas from a galaxy passing through. The large overdensity of star-forming galaxies clustered around the quasar (see Sect. 4.5) supports the hypothesis of a massive quasar host halo, as detected, for instance, in MQN01 (Travascio et al. 2025).

While the galaxy encounter scenario currently offers a plausible explanation for the observed extended He II emission, confirming the existence of galaxies physically interacting with the quasar, and eventually unveiling the hypothesized hot CGM, requires deeper, higher-resolution, multiwavelength imaging and spectroscopic observations. We emphasize that the scenarios discussed above are not exhaustive: more complex configurations remain possible, including, for example, an asymmetric outflow triggered by a major merger event, or an asymmetric inflow of gas previously stripped from a companion galaxy. We also note that all the scenarios described above may need to be revised to account for the possibility that the quasar host does not trace the center of mass of the system and its redshift measurement is affected by its proper motion (e.g., Pensabene et al. 2025). Finally, it remains to be understood what is the origin of the more extended portions of the Lyα nebula currently undetected in He II.

6. Summary and conclusions

In this work, we presented ALMA and MUSE-NFM observations toward the field of the quasar Q0055-269 at z ≈ 3.66 and the surrounding MQN04, one of the brightest Lyα nebulae known at high redshift, which is also detected in C IV and He II emission extending over ≳ 50 kpc (Borisova et al. 2016; Guo et al. 2020). Combined with the view of the ionized gas phase offered by previous MUSE-WFM observations, these new data allow us to resolve the kinematics, morphology, and clumpiness of the CGM gas flows on scales ranging from tens of kiloparsecs down to #x2248; 1 kpc away from the quasar.

• We detected both the (observed-frame) 3 mm continuum and the CO(4–3) line emission from the quasar host galaxy, providing the first measurement of its systemic redshift and allowing us to accurately study the kinematics of the surrounding ionized gas relative to the quasar. We also detected the emission of the CO(4–3) line emission from a second source, MQN04-QC, offset by Δv_QSO ≈ −1100 km s⁻¹ from the quasar. The non-detection of any stellar continuum in UV observations suggests that MQN04-QC is a heavily dust-obscured galaxy.

• The combined analysis of the Lyα, C IV, and He II lines reveals emission extending over ≈ 50 kpc in projection. We distinguish two main bright structures on the east and west sides of the quasar host, which appear kinematically connected by a bridge of diffuse gas. The emission line profiles throughout these regions are broad (FWHM ≳ 500 km s⁻¹) and blueshifted by ≈ 0 − 800 km s⁻¹ relative to the quasar systemic redshift. Interestingly, the most blueshifted component aligns in velocity with an intervening H I absorber with a column density N_H I = 10^14.64 cm⁻² detected in the quasar spectrum, thus suggesting that the gas is highly ionized. As expected for resonant lines, Lyα and C IV show smoother spatial and velocity distributions than nonresonant He II yet, with the exception of the central region, they do not exhibit a systematic velocity shift relative to He II. This indicates that these lines reliably trace the average ionized-gas kinematics across most of the nebula, likely because quasar illumination keeps the CGM highly ionized and radiative-transfer effects relatively modest. On smaller scales, the MUSE-NFM observations uncover two additional compact He II-emitting clumps located ≈ 5 kpc away from the quasar and at a similar redshift. Similarly to the emission on larger scales, these structures are characterized by broad (FWHM ≳ 500 km s⁻¹) and complex He II line profiles. We explored three possible scenarios to explain the observed He II emission: an inflowing gas stream, an outflow, and gas stripping induced by a close galaxy encounter or by a passage of a galaxy through a hot CGM.

• We conducted a census of the galaxy population and identified five star-forming galaxies clustered within |Δv_QSO| ≲ 1000 km s⁻¹ around the central quasar. Three of them are embedded in the Lyα nebula, while the remaining two are located at larger distances, but aligned with the orientation of the nebula, and we speculate that they may also be connected with fainter yet undetected Lyα emission. We measure an overdensity of δ ≈ 41. Compared to other overdensities at similar redshifts, this structure resembles known protoclusters (e.g., MQN01) but lacks enhanced AGN activity or an excess of submillimeter galaxies. This may indicate that such systems trace different regions of the cosmic web or are captured at distinct stages of their evolution.

Follow-up observations with ALMA achieving higher sensitivity and with JWST/NIRSpec are crucial to map the cold molecular gas distribution and detect additional nonresonant lines, such as Hα, to uncover both the mechanisms powering the extended emission and any additional AGN or massive quiescent galaxies that remain undetected in the current data. Our joint analysis of spatially resolved kinematics traced by the nonresonant He II emission and secure systemic redshift determination via the CO(4–3) line shows that MQN04 fits within the broader picture in which gas flows from large to small scales are shaped by galaxy interactions and by the interplay between feeding and feedback. Most importantly, for the first time, we unambiguously measure the presence of a strong systematic velocity shift between the CGM kinematics and the quasar host systemic redshift, demonstrating that the CGM in this system is highly asymmetric. Overall, systems such as MQN04 are ideal laboratories for studying the baryon cycle in massive halos at high redshift and establishing how the interplay between the kinematics of different gas phases and the large-scale environment drives the coevolution of galaxies and the cosmic web across cosmic time.

Acknowledgments

We thank the anonymous referees for their valuable suggestions, which helped improve this work. We thank Mr. Cosimo Marconcini and Dr. Giacomo Venturi for their contribution and insightful discussion on the outflow scenario. This project was supported by Progetto FARE 2020 Svelare i nodi massicci della CosmicWeb ID 2021-NAZ-0326/PER, the European Research Council (ERC) Consolidator Grant 864361 (CosmicWeb), and by Fondazione Cariplo grant no. 2020-0902. AP acknowledges the support from the Independent Research Fund Denmark (DFF) under grant 3120-00043B. This paper makes use of the following ALMA data: ADS/JAO.ALMA#2024.1.00499.S ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan), together with NRC (Canada), NSTC 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 work is also based on observations collected at the European Southern Observatory under ESO programmes IDs 094.A-0131(B), 096.A-0222(A), 109.232M and data obtained from the ESO Science Archive Facility with DOI under https://doi.eso.org/10.18727/archive/42. This research made use of Astropy (http://www.astropy.org), a community-developed core Python package for Astronomy (Astropy Collaboration 2013, 2018), NumPy (Harris et al. 2020), SciPy (Virtanen et al. 2020), Matplotlib (Hunter 2007).

References

Appendix A: MUSE-selected galaxies in the field

In Table A.1, we list the galaxies detected by MUSE in the MQN04 field within |Δv_QSO| ≤ 1000 km s⁻¹.

Appendix B: Velocity field and velocity distribution of the He II nebula

In Fig. B.1, we present He II velocity shift (1st moment) and dispersion (2nd moment) extracted from MUSE WFM and NFM datacubes. The He II emitting gas appears blueshifted up to ≈ − 700 km s⁻¹ with respect to the quasar systemic redshift determined by the ALMA CO(4–3) line detection, differently from what is reported by Guo et al. (2020) who employed as a reference frame for the line-of-sight velocities the redshift derived from the Lyα emission line profile of the diffuse nebula. The observed velocity dispersion reaches values of ≳300 km s⁻¹ around regions E and B, consistent with the findings of Guo et al. (2020). In addition, NFM data reveal similarly high velocity dispersion values in the regions C1 and C2, within the inner few kiloparsecs from the quasar host galaxy. Such an increase of the He II gas velocity dispersion in the proximity of the quasar might be associated with higher turbulence caused by merger activity, outflows, or gas inflows into the gravitational potential well of the quasar halo (see Sect. 5).

Fig. B.1. The observed He II gas velocity shift relative to the quasar host galaxy systemic redshift (1st moment, top panels), and velocity dispersion (2nd moment, bottom panels) as derived from the MUSE WFM and NFM data (left and right panels, respectively). The contours refer to the He II SB levels as in Fig. 3. The regions identified in Sect. 4.3 are also reported.

All Tables

All Figures

Fig. 1. Spatial distribution of galaxies in the MQN04 field detected within |Δv ≤ 100 km s−1|, the ALMA-identified MQN04-QC galaxy at Δv ≈ −1172 km s−1 of the quasar systemic redshift, and the Lyα, C IV, and He II-emitting gas. Left panel: HST/F160W image. The star-forming galaxies detected in MUSE data are indicated by red and blue squares, depending on their line-of-sight velocity relative to the quasar’s systemic redshift as determined from the CO(4–3) line emission. The quasar host and the MQN04-QC companion galaxy are shown as circles. The contours indicate the detected extent of the Lyα (gold), C IV (cyan), and He II (magenta) nebulae (see Sect. 4). Central panel: QSO PSF-subtracted cutout of the HST/F160W image. The dashed lines mark the position of the quasar host and MQN04-QC galaxy as identified from ALMA observations. The region masked for the PSF modeling is also shown as a hatched circle (see Sect. 3.1). Right panel: ALMA CO(4–3) line-velocity integrated map and 3 mm continuum. The gray (CO(4–3) line emission) and black (3 mm continuum) contours correspond to [ − 2, 2.5n]σ (dotted and solid for the negative and positive levels, respectively), where n ≥ 1 is an integer, and σ is the RMS noise. The ALMA synthesized beam is shown as a yellow ellipse in the bottom left corner.

In the text

Fig. 2. Signal-to-noise spectra extracted from the regions identified in the line-emitting nebulae. Top panels: Lyα (orange), C IV (blue), and He II (red) emission lines extracted from the MUSE-WFM data within the regions E (left), B (middle), and W (right). The Lyα S/N is rescaled by ×0.3 for visualization purposes. Bottom-left panels: He II emission line S/N spectra extracted from the MUSE-NFM data within the clumps C1 (top) and C2 (bottom). The shaded band marks spectral regions contaminated by skylines. Bottom-center panels: CO(4–3) emission line profiles of the quasar host (top) and the MQN04-QC galaxy (bottom). The data are shown in yellow, the single Gaussian best-fit curve is in black with 1σ uncertainty shaded in gray. The flux is in units of mJy. Bottom-right panel: H I absorption line detected in the (continuum-normalized) UVES spectrum of the quasar. The green line indicates the best-fit Voigt profile. In all the panels, the line-of-sight velocity is relative to the systemic redshift of the quasar host galaxy as determined from its CO(4–3) emission line. The green dashed line and shaded area indicate the centroid and the FWHM of the H I absorption line, respectively.

In the text

Fig. 3. Optimally extracted images of the extended Lyα λ1215 Å (top-left panel), C IV λ1548 Å (top-right panel), and He II λ1460 Å (MUSE-WFM and NFM data, bottom-left and bottom-right panels, respectively) nebulae detected around the quasar Q0055-269 in the MQN04 field. The contours are matched in flux density and indicate 2n × s0 SB levels with n ≥ 1 integer, and s0 = 0.5 × 10−18 and 16 × 10−18 erg s−1 cm−2 arcsec−2 for the WFM and NFM, respectively. The layer of the datacube at the quasar redshift has been added as a reference noise background for each transition. The green, magenta, and blue ellipses delineate the selected E (east), B (bridge), and W (west) regions identified based on the He II emission detected in the MUSE-WFM data, while the white and orange ellipses enclose the two clumps (C1 and C2, respectively) revealed by the MUSE-NFM observations. The gray circle indicates the location of the quasar host as determined by ALMA observations, and the size is rescaled to match the PSF FWHM. The blue circle marks the position of MQN04-QC. The contours of the He II nebula detected in the MUSE-WFM data (white), the quasar CO(4–3) emission (green, 2n × RMS, with n ≥ 1), as well as the ALMA synthesized beam, are overlaid on the MUSE-NFM data in the bottom-right panel.

In the text

Fig. 4. Kinematic analysis of the He II extended emission extracted from the MUSE-WFM data. Left panels: channel maps of the He II emission obtained from the PSF and continuum-subtracted datacube by collapsing the voxels within windows of 100 km s−1. The contours refer to He II SB levels as in Fig. 3. Central panel: line-of-sight velocity map (1st moment) of the He II emission. The position of the quasar is indicated by the gray circle with a size corresponding to the QSO PSF FWHM. Right panel: position-velocity diagram extracted along the polygonal slit shown in the central panel. The boundaries of the E, B, and W region crossed by the slit are reported with vertical lines. The dotted and solid contours indicate S/N = 4 and 7 levels, respectively. The spectral and spatial resolution is represented as a rectangle in the bottom-left corner.

In the text

Fig. 5. Line ratios in the MQN04 nebula: C IV/Lyα (left panel), C IV/He II (central panel), and He II/Lyα (right panel). The line ratios are computed in Voronoi-resampled datacubes to maximize the signal-to-noise of the faintest regions. The gray contour in each panel represents the SB mask of the optimally extracted image, as in Fig. 3, of the least extended nebula included in the line ratio. The colored ellipses mark the same regions identified in Fig. 3. The gray circle indicates the location of the quasar host as determined by ALMA observations and is sized as the PSF FWHM.

In the text

	Fig. B.1. The observed He II gas velocity shift relative to the quasar host galaxy systemic redshift (1st moment, top panels), and velocity dispersion (2nd moment, bottom panels) as derived from the MUSE WFM and NFM data (left and right panels, respectively). The contours refer to the He II SB levels as in Fig. 3. The regions identified in Sect. 4.3 are also reported.
In the text