# Physical Constraints on the Extended Interstellar Medium of the z = 6.42 Quasar J1148+5251: [C II]158 µm, [N II]205 µm, and [O I]146 µm Observations

Romain A. Meyer, Fabian Walter, Claudia Cicone, Pierre Cox, Roberto Decarli, Roberto Neri, Mladen Novak, Antonio Pensabene, Dominik Riechers,

Axel Weiss

► To cite this version:

Romain A. Meyer, Fabian Walter, Claudia Cicone, Pierre Cox, Roberto Decarli, et al.. Physical Constraints on the Extended Interstellar Medium of the z = 6.42 Quasar J1148+5251: [C II]158 µm, [N II]205 µm, and [O I]146 µm Observations. The Astrophysical Journal, 2022, 927, 10.3847/1538-4357/ac4e94. insu-03748292

HAL Id: insu-03748292

https://insu.hal.science/insu-03748292

Submitted on 9 Aug 2022

**HAL** is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire **HAL**, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
OPEN ACCESS

<figure><img src="image_2.png" /></figure>

CrossMark

# Physical Constraints on the Extended Interstellar Medium of the $z=6.42$ Quasar J1148+5251: [C II]$_{158}$ µm, [N II]$_{205}$ µm, and [O I]$_{146}$ µm Observations

Romain A. Meyer¹, Fabian Walter¹, Claudia Cicone², Pierre Cox³, Roberto Decarli⁴, Roberto Neri⁵, Mladen Novak¹, Antonio Pensabene⁴, Dominik Riechers⁶, and Axel Weiss⁷

¹Max Planck Institute for Astronomy, Königstuhl 17, D-D-69117, Heidelberg, Germany; meyer@mpia.de

²Institute of Theoretical Astrophysics, University of Oslo, P.O. Box 1029, Blindern, NO-0315 Oslo, Norway

³Sorbonne Université, UPMC Université Paris 6 and CNRS, UMR 7095, Institut d'Astrophysique de Paris, 98bis Boulevard Arago, F-75014 Paris, France

⁴INAF—Osservatorio di Astrofisica e Scienza dello Spazio di Bologna, via Gobetti 93/3, I-40129, Bologna, Italy

⁵Institut de Radioastronomie Millimétrique (IRAM), 300 Rue de la Piscine, F-38400 Saint-Martin-d'Hères, France

⁶I. Physikalisches Institut, Universität zu Köln, Zülpicher Strasse 77, D-50937 Köln, Germany

⁷Max-Planck-Institut für Radioastronomie, Auf dem Hügel 71, D-53121 Bonn, Germany

Received 2021 July 6; revised 2021 October 6; accepted 2022 January 20; published 2022 March 11

## Abstract

We report new Northern Extended Millimeter Array observations of the [C II]$_{158}$ µm, [N II]$_{205}$ µm, and [O I]$_{146}$ µm atomic fine structure lines (FSLs) and dust continuum emission of J1148+5251, a $z = 6.42$ quasar, which probe the physical properties of its interstellar medium (ISM). The radially averaged [C II]$_{158}$ µm and dust continuum emission have similar extensions (up to $\theta = 2.51_{-0.25}^{+0.46}$ arcsec, corresponding to $r = 9.8_{-2.1}^{+3.3}$ kpc, accounting for beam convolution), confirming that J1148+5251 is the quasar with the largest [C II]$_{158}$ µm-emitting reservoir known at these epochs. Moreover, if the [C II]$_{158}$ µm emission is examined only along its NE–SW axis, a significant excess ($>5.8\sigma$) of [C II]$_{158}$ µm emission (with respect to the dust) is detected. The new wide-bandwidth observations enable us to accurately constrain the continuum emission, and do not statistically require the presence of broad [C II]$_{158}$ µm line wings that were reported in previous studies. We also report the first detection of the [O I]$_{146}$ µm and (tentatively) [N II]$_{205}$ µm emission lines in J1148+5251. Using FSL ratios of the [C II]$_{158}$ µm, [N II]$_{205}$ µm, [O I]$_{146}$ µm, and previously measured [C I]$_{369}$ µm emission lines, we show that J1148+5251 has similar ISM conditions compared to lower-redshift (ultra)luminous infrared galaxies. CLOUDY modeling of the FSL ratios excludes X-ray-dominated regions and favors photodissociation regions as the origin of the FSL emission. We find that a high radiation field ($10^{3.5-4.5} G_0$), a high gas density ($n \approx 10^{3.5-4.5} \text{ cm}^{-3}$), and an H I column density of $10^{23} \text{ cm}^{-2}$ reproduce the observed FSL ratios well.

Unified Astronomy Thesaurus concepts: Quasars (1319); High-redshift galaxies (734); Interstellar medium (847); Interferometry (808); Millimeter astronomy (1061)

## 1. Introduction

Luminous quasar activity is a key process of galaxy evolution. Indeed, massive outflows driven by the radiation pressure generated by the accretion of gas onto the central supermassive black hole (SMBH) or so-called “Active Galactic Nuclei (AGNs) feedback” are invoked in most models of galaxy formation to clear massive galaxies of their gas and quench star formation (e.g., Silk & Rees 1998; Di Matteo et al. 2005; Springel et al. 2005; King 2010; Costa et al. 2014; Ishibashi & Fabian 2015; Richardson et al. 2016; Negri & Volonteri 2017; Oppenheimer et al. 2020; Koudmani et al. 2021). Luminous quasars are most interesting at high redshift in particular, when they probe the early phase of coevolution between the first galaxies and their central black hole. The advent of large optical and infrared surveys has enabled the discovery of quasars up to $z \sim 7.5$ (Bañados et al. 2018; Yang et al. 2020; Wang et al. 2021), with several hundreds at $z > 6$ (see Bosman 2020 for an up-to-date list). These early quasars harbor SMBHs with $M_* \gtrsim 10^8 M_\odot$ and accrete gas at or near the Eddington limit for most of their life, challenging models of SMBH formation and growth (e.g., De Rosa et al. 2014; Mazzucchelli et al. 2017;

Bañados et al. 2018; Wang et al. 2021). Because the bright quasar light outshines that of the host in the optical and near-infrared, the galaxies hosting early luminous quasars remained relatively mysterious until the advent of modern (sub)millimeter observatories.

Numerous observations of $z > 6$ quasars targeting the bright far-infrared (FIR) [C II]$_{158}$ µm emission line have revealed their host galaxies to be infrared luminous, dusty, and actively forming stars with estimated rates of $10^2-10^3 M_\odot \text{yr}^{-1}$ (e.g., Walter et al. 2003, 2009a; Maiolino et al. 2005, 2012; Wang et al. 2013, 2019; Bañados et al. 2015; Decarli et al. 2018; Venemans et al. 2018, 2020; Novak et al. 2019, 2020; Yang et al. 2020). [C II]$_{158}$ µm kinematics also show that most quasar hosts are massive galaxies ($M_* \sim 10^{10} M_\odot$) displaying a variety of morphologies, such as stable disks, bulge-dominated galaxies, and mergers with nearby companions (Wang et al. 2013, 2019; Shao et al. 2017, 2019; Decarli et al. 2019a, 2019b; Neeleman et al. 2019, 2021). The reports of broad [C II]$_{158}$ µm line wings in the $z = 6.42$ quasar J1148+5251 (Maiolino et al. 2012; Cicone et al. 2015) spurred the use of the [C II]$_{158}$ µm emission line to identify quasar outflow signatures in the early universe. Such features have remained rare, however, and stacking analyses have led to contradictory results (Bischetti et al. 2019; Novak et al. 2020). Recently, Izumi et al. (2021a, 2021b) reported broad [C II]$_{158}$ µm line wings in two low-luminosity quasars at $z = 6.72$ and $z = 7.07$.

<figure><img src="image_3.png" /></figure>

Original content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.
Two decades of [C II]158 µm and CO studies have shown that early luminous quasars provide an unparalleled observational window into the physics of the earliest (and most massive) galaxies in the universe. However, multiline studies using lines other than [C II]158 µm and CO transitions have been much rarer until now. Since different fine structure lines (FSLs) trace different gas densities and excitation levels, only in combination can they probe the ionized and neutral atomic gas phases, and the excitation source(s) of the gas (see, e.g., Carilli & Walter 2013 for a review). Besides [C II]158 µm, potential atomic FSLs of interest include [N II]122 µm, [N II]205 µm, [O I]63 µm, [O I]145 µm, [O III]88 µm, and [C I]369 µm, all accessible at $z \sim 6$ with (sub)millimeter arrays such as the Atacama Large Millimeter Array or the Northern Extended Millimeter Array (NOEMA). Moreover, these lines have been observed with Herschel in large samples of local (ultra)luminous infrared galaxies—or (U)LIRGs—which can be readily compared to $z \sim 6$ quasar hosts (e.g., Díaz-Santos et al. 2017; Herrera-Camus et al. 2018a).

Detections of FSLs other than [C II]158 µm in $z > 6$ quasars are still relatively recent (Walter et al. 2018; Hashimoto et al. 2019; Novak et al. 2019; Li et al. 2020), and a complex picture is emerging from these first results. Emission lines probing the neutral phase ([O I]145 µm and [C I]369 µm) show good agreement between the line ratios and line-to-FIR ratios of distant quasars and local (U)LIRGs (Novak et al. 2019; Li et al. 2020). While Novak et al. (2019) and Li et al. (2020) report potentially high [O III]88 µm to [C II]158 µm ratios (~2) in J1342+0928 and J2310+1855, respectively, this is not the case for the quasar J2100–1715 and its companion galaxy (Walter et al. 2018), which have ratios similar to the average of the local population of LIRGS and AGNs (~0.1–1; Herrera-Camus et al. 2018a). These first results are difficult to interpret, however, since the origins of FSLs are not always clearly determined, and can be linked to different phases (as is the case for the [C II]158 µm emission line). Clearly, more FIR multiline studies of $z > 6$ quasars are needed to understand the interstellar medium (ISM) of their host galaxies.

SDSS J1148+5251 is one of the earliest high-redshift quasars discovered in the Sloan Digital Sky Survey (SDSS; Fan et al. 2003; $z = 6.4189$), and harbors a $3 \times 10^9 M_{\odot}$ SMBH (Willott et al. 2003). Being the redshift record holder for many years after its discovery, it was extensively observed with the Very Large Array and the IRAM Plateau de Bure Interferometer (PdBI), and was the first object detected in CO and [C II]158 µm at $z > 5$ (Bertoldi et al. 2003b, 2003a; Walter et al. 2003, 2004, 2009a, 2009b; Maiolino et al. 2005; Riechers et al. 2009). These pioneering studies probed the host galaxy star formation rate (SFR), dust, and ISM properties. Additionally, Maiolino et al. (2012) and Cicone et al. (2015) reported the presence of a broad [C II]158 µm emission ($\sigma_v = 900 \text{ km s}^{-1}$) component in the PdBI data, suggesting the presence of an outflow as well as spatially extended [C II]158 µm emission (up to $r \sim 30$ kpc). In this paper, we return to J1148+5251 with a new set of NOEMA observations, capitalizing on larger bandwidths and more antennas, thus improving on the image fidelity as compared to earlier PdBI observations. The new observations targeted atomic FSLs ([O I]146 µm, [N II]205 µm) and other molecular (CO, H$_2$O) rotational transitions. The aim of the observations was to dissect the ISM phases without relying on assumptions about the origin of [C II]158 µm that can come from both the ionized and neutral phases. Indeed, [O I]146 µm traces exclusively the neutral phase/photodissociated regions (PDRs), whereas [N II]205 µm

traces the ionized/H II regions. This set of observations is complemented with earlier [C I]369 µm data (Riechers et al. 2009) that trace the neutral/molecular gas. Thanks to the wide spectral coverage of the new NOEMA correlator *PolyFix* and the upgraded NOEMA array, these observations achieved a high fidelity that resulted in deep [C II]158 µm observations and tight constraints on the underlying dust continuum.

The structure of this paper is as follows. We present in Sections 2 and 3 the new observations of the [C II]158 µm, [N II]205 µm, and [O I]146 µm emission lines in J1148+5251, as well as the FIR continuum observed between 200 and 280 GHz. We focus on the [C II]158 µm emission line to reassess the evidence for a broad velocity component and investigate its spatial extension in Section 4. In Section 5, we derive ISM properties from the strength of the atomic FSLs observed, before concluding our study in Section 6. Throughout this paper, we assume a concordance cosmology with $H_0 = 70 \text{ km s}^{-1}\text{Mpc}^{-1}$, $\Omega_M = 0.3$, and $\Omega_\Lambda = 0.7$. At the redshift of the target ($z = 6.42$), 1" corresponds to 5.62 proper kpc.

## 2. Observations and Data Reduction

We have observed the $z = 6.4189$ quasar J1148+5251 using NOEMA (Project ID: w17ex001/w17ex001; PI: F. Walter). The pointing and phase center of our observations were chosen to correspond to the quasar position in the optical SDSS imaging (R.A. = 11:48:16.64, decl. = +52:51:50.32). The observations included two spectral setups taking advantage of the new *PolyFix* correlator, covering simultaneously two 7.744 GHz wide sidebands.<fnref n="8" /> The first spectral setup was centered at 267 GHz, such that the lower sideband (255–263 GHz) covers the redshifted [C II]158 µm emission with one sideband, while the upper sideband (271–279 GHz) covers the [O I]146 µm emission. The second setup was centered at 208 GHz to cover the [N II]205 µm in the lower sideband (196–204 GHz). The setups also covered two high-J CO (14–13 and 13–12) and H$_2$O rotational transitions (5$_{23}$–4$_{32}$ and 3$_{22}$–3$_{13}$), which will be discussed in future works. The observations were executed between 2017 December and 2018 May. The [C II]158 µm and [O I]146 µm setup was mostly observed in configuration 9D, except for two tracks using eight antennas (with baselines ranging from 24 to 176 m), for a total observing time of 18.7 hr. Data was (remotely) reduced at IRAM Grenoble using the CLIC package within the GILDAS framework (jan2021a version).<fnref n="9" /> We reach an rms noise of 0.64(0.88) mJy beam$^{-1}$ in 50 km s$^{-1}$ channels, and the synthesized beam FWHM size is 1"97 × 1"59 (1"83 × 1"51) for the [C II]158 µm ([O I]146 µm) line observations. The [N II]205 µm and CO lines were observed with eight antennas for a total of 15.7 hr, with baselines ranging from 24 to 176 m. For the [N II]205 µm line, the noise rms is 0.78 mJy beam$^{-1}$ in 50 km s$^{-1}$ channels, and the synthesized beam size 1"79 × 1"51. For the continuum, the synthesized beam size achieved is 1"78 × 1"51, 1"65 × 1"46, 1"91 × 1"71, and 1"88 × 1"52 at 200, 212, 259, and 272 GHz, respectively.

Imaging and cleaning was performed using the latest version of MAPPING/GILDAS (jan2021a). The dirty maps were obtained from the visibilities without tapering and using natural weighting. The data were not primary beam-corrected.<fnref n="10" />

<footnote n="8">⁸ https://www.iram-institute.org/EN/content-page-96-7-56-96-0-0.html</footnote>

<footnote n="9">⁹ https://www.iram.fr/~gildas/dist/</footnote>

<footnote n="10">¹⁰ As the source is located at the phase center, even at the edge of the aperture considered in the analysis ($r = 3''$) the correction is minimal (1.06) and does not impact the results.</footnote>
<figure><img src="image_4.png" /><figcaption>**Figure 1.** Dust continuum emission maps integrated over the ~7.68 GHz effective bandwidth of the Polyfix correlator for each sideband, with central frequencies 200, 216, 259, and 274 GHz. The dotted circle centered on J1148+5251 indicates the aperture radius of 3" adopted throughout this paper for flux density estimates. An additional source, reported previously by Leipski et al. (2010) and Cicone et al. (2015), is visible to the NW in the three higher-frequency bands at ~10" (the black dotted circle). The contours are logarithmic (−4, −2, 2, 4, 8, 16, 32)σ (rms). The color scaling is log-linear, the threshold being at 3σ.</figcaption></figure>

Cleaning was performed down to $2\sigma$ (where $\sigma$ is the rms noise in the dirty map) using a circular clean region of radius $r = 5''$. The reason for choosing such a wide radius was the earlier reports of extended [C II]$_{158 \mu m}$ emission (Maiolino et al. 2012; Cicone et al. 2015). An additional clean region with radius $r = 2''$ was added on the NW source reported by Leipski et al. (2010) and Cicone et al. (2015), which is also detected in the NOEMA data.<fnref n="11" /> The final products were created using the following procedure. First, data cubes with 50 km s⁻¹ channels were produced to search for significant emission lines. The [C II]$_{158 \mu m}$, [N II]$_{205 \mu m}$, [O I]$_{146 \mu m}$, H₂O (523–432 and 322–313, $\nu_{rest} = 258.7$ GHz), CO(14–13), and CO(13–12) emission lines were fitted with a single-Gaussian profile in order to estimate their FWHM. Continuum maps were created using all channels at least 1.25 × FWHM away from each emission line. The continuum, determined from the line-free channels using an order 1 interpolation (GILDAS UV_BASE-LINE routine) over the ~7.6 GHz sidebands, was subtracted in the uv plane to create continuum-subtracted cubes.

In order to determine the significance of emission lines, velocity-integrated emission line maps (“line maps”) were created by integrating channels over 1.2 times the FWHM of the [C II]$_{158 \mu m}$ line (i.e., 482 km s⁻¹) at the redshifted frequency of the line. We use the [C II]$_{158 \mu m}$ redshift ($z = 6.4189$) to determine the redshifted frequency of all lines. Such maps, assuming the line is Gaussian, contain by definition 84% of the total flux<fnref n="12" /> (for a short derivation, see Appendix A of Novak et al. 2020). All of the total line fluxes measured from the line maps and reported in this paper account for this effect. Additionally, all continuum and line fluxes in this paper were computed using the residual scaling method (e.g., Jorsater & van Moorsel 1995; Walter & Brinks 1999; Walter et al. 2008; Novak et al. 2019).

In order to determine the aperture needed to recover most of the flux of the emission lines and the dust continuum, a curve of growth approach was adopted. We show in Appendix A that all

<caption>**Table 1**</caption>

Continuum Measurements at 200–274 GHz Extracted from the Continuum Images (Figure 1, Effective Bandwidth ~7.68 GHz (Masking the Edges and Central Baseband Gap) and Integrated in a Circular Aperture with Radius $r = 3''$)

<table><thead><tr><th>ν<sub>obs</sub> (GHz)</th><th>S<sub>v</sub> (mJy)</th><th>Rms (mJy)</th><th>Beam Size</th></tr></thead><tbody><tr><td>200</td><td>1.78 ± 0.23</td><td>0.06</td><td>1"78 × 1"51</td></tr><tr><td>216</td><td>2.33 ± 0.21</td><td>0.05</td><td>1"65 × 1"46</td></tr><tr><td>259</td><td>4.64 ± 0.26</td><td>0.08</td><td>1"91 × 1"71</td></tr><tr><td>274</td><td>5.85 ± 0.29</td><td>0.08</td><td>1"88 × 1"52</td></tr></tbody></table>

line and continuum fluxes reach a maximum or plateau at an aperture radius $r = 3''$, which corresponds to 16.9 kpc at $z = 6.42$. A nominal aperture radius of 3" is thus adopted throughout the paper. The line subtraction procedure described above was repeated using the $r = 3''$ aperture to obtain the final data cubes. Additionally, we investigate the use of multiscale cleaning in Appendix B. We conclude that the [C II]$_{158 \mu m}$ emission and the 259–272 GHz continuum are better recovered using multiscale cleaning, and we therefore use multiscale cleaning with Gildas/MAPPING for these data throughout the paper.

## 3. Results

### 3.1. Dust Continuum Emission

We first present the FIR continuum maps in Figure 1. The FIR continuum is clearly detected in all four sidebands. The measured continuum flux densities are tabulated in Table 1. Due to the upgraded bandwidth of the NOEMA PolyFix correlator, these continuum measurements have higher sensitivity than previous observations at ~260 GHz (Walter et al. 2009a; Maiolino et al. 2012; Cicone et al. 2015). The new continuum flux density at 259 GHz (4.64 ± 0.26 mJy) is in good agreement<fnref n="13" /> with the earlier PdBI measurement from Walter et al. (2009a) and the 1.2 mm continuum measurement (5.0 ± 0.6 mJy; Bertoldi et al. 2003a). We combine our new continuum measurements with previous literature results at different frequencies (Bertoldi et al. 2003b; Walter et al. 2003; Robson et al. 2004; Riechers et al. 2009; Leipski et al. 2010;

<footnote n="11">¹¹ We found no emission lines for this source in any of the spectral setups used in this work.</footnote>

<footnote n="12">¹² We have checked that this correction holds for the [C II]$_{158 \mu m}$ emission, where the flux recovered in a (−1500, +1500) km s⁻¹ velocity-integrated map is $S_{line}\Delta v = 10.5 \pm 0.92$ Jy km s⁻¹, and the 1/0.84-corrected flux in the nominal 1.2 × FWHM$_{[C II]} = 482$ km s⁻¹ channel is $S_{line}\Delta v = 10.4 \pm 0.5$ Jy km s⁻¹, where $S_{line}$ is integrated in a circular aperture with $r = 3''$. Indeed, the choice of a 1.2 × FWHM-wide channel maximizes the signal-to-noise ratio without flux losses.</footnote>

<footnote n="13">¹³ We have checked that there is no continuum offset between the two sidebands by imaging the calibrators (1150+497 and 1216+487) for every track and the stacked data.</footnote>
<figure><img src="image_5.png" /><figcaption>**Figure 2.** Continuum measurements from the literature and this work with the best-fit modified blackbody dust emission model. The effect of the CMB is accounted for as prescribed by Da Cunha et al. (2013), and the data points at $ν > 10^3$ GHz are not used for the fit. The FIR luminosity, integrated over the shaded orange area, is not significantly affected by this choice.</figcaption></figure>

Gallerani et al. 2014) to fit the FIR spectral energy distribution (SED).

Assuming optically thin dust emission at $λ > 40$ µm (e.g., Beelen et al. 2006), we use a modified blackbody model for the dust emission and correct both for contrast and cosmic microwave background (CMB) heating, as prescribed by Da Cunha et al. (2013). The dust mass is derived assuming an opacity $κ_{ν_{rest}} = κ_{ν_0}(ν_{rest}/ν_0)^β$, with $ν_0 = c/(125~µm)$ and $κ_0 = 2.64~m^2 kg^{-1}$, following Dunne et al. (2003), with $β$ being the dust spectral emissivity index. Our purpose is primarily to measure the FIR luminosity to constrain the SFR in J1148+5251. Therefore, we omit data points at $ν_{rest} > 1000$ GHz ($λ_{rest} \lesssim 125~µm$), where contamination by the quasar nonthermal and torus emission becomes significant (e.g., Leipski et al. 2010, 2014). The dust SED model uses three free parameters (total dust mass $M_d$, dust emissivity index $β$, and dust temperature $T_d$) and is fitted using a Markov Chain Monte Carlo with the `emcee` package (Foreman-Mackey et al. 2013). The resulting best-fit and observational constraints are shown in Figure 2, and the posterior probability distribution of the dust SED parameters is displayed in Appendix C. The median dust mass is $3.2 \times 10^8 M_\odot$, the median dust SED index is $β = 1.77_{-0.26}^{+0.30}$, and the median dust temperature is moderately high ($T_d = 51.3$ K), in agreement with earlier studies (e.g., Beelen et al. 2006; Leipski et al. 2010; Cicone et al. 2015), which is not surprising, considering that most of the constraining power comes from observations at lower and higher frequencies than those reported in this work.

We integrate the modified blackbody to derive the total infrared (IR; 8–1000 µm) and FIR (42.5–122.5 µm; e.g., Helou et al. 1985) luminosities. The total IR luminosity is $L_{IR} = (20.9 \pm 6.8) \times 10^{12} L_\odot$ and the FIR luminosity $L_{FIR} = (13.4 \pm 2.4) \times 10^{12} L_\odot$, in agreement with earlier studies of J1148+5251. Assuming that dust heating is dominated by young stars, the IR luminosity can be converted to an SFR using the Kennicutt (1998) and Kennicutt & Evans (2012) conversions, giving SFR = (1830 ± 595) – (2090 ± 680)$M_\odot$ yr$^{-1}$, respectively. This is in agreement with earlier studies that found that J1148+5251 is in an intense starburst phase (Maiolino et al. 2005, 2012; Walter et al. 2009a; Cicone et al. 2015).

## 3.2. FSL Detections

We present the line maps for the [C II]$_{158~µm}$, [O I]$_{146~µm}$, and [N II]$_{205~µm}$ emission in Figure 3. For comparison, we also show maps of the same channel width (482 km s$^{-1}$) as the main line map, but offset by ±482 km s$^{-1}$ from the line emission, to visually assess the robustness of our detections. The spectra of the [C II]$_{158~µm}$ line are presented in Figure 5 (see also Appendix D), while the [O I]$_{146~µm}$ and [N II]$_{205~µm}$ spectra are shown in Figure 4.

We detect [C II]$_{158~µm}$ and [O I]$_{146~µm}$ at $42σ$ and $5.3σ$ (where the signal-to-noise ratio (S/N) is calculated using the peak surface brightness in the line maps and the pixel rms level). [N II]$_{205~µm}$ is only marginally detected ($3.7σ$), and its peak emission is potentially offset from the rest-frame UV, dust, [C II]$_{158~µm}$, and [O I]$_{146~µm}$ emission. In the central pixel of the [C II]$_{158~µm}$ emission, [N II]$_{205~µm}$ is formally undetected. We defer further discussion of the [N II]$_{205~µm}$ emission line to Section 5.2. In our subsequent analysis of the ISM of J1148+5251, we focus on the extended aperture-integrated fluxes. All of the fluxes and derived line luminosities are tabulated in Table 2. For comparison, we have also measured the flux density of the [C II]$_{158~µm}$ emission by fitting the aperture-integrated spectrum with a Gaussian, and we show a comparison between our work and earlier studies for various aperture sizes in Appendix D. In summary, we find a good agreement between the line map flux estimates and the spectrum best-fit values, as well as good agreement between previous studies and the data presented here, as long as the “core” ($<3''$ emission, with no velocity offsets) [C II]$_{158~µm}$ emission is considered.

# 4. The Spatial and Velocity Structure of the [C II]$_{158~µm}$ Emission in J1148+5251

## 4.1. Spectral Analysis

Figure 5 shows the total and continuum-subtracted [C II]$_{158~µm}$ spectra (extracted in an $r = 3''$ aperture), as well as single- and double-Gaussian fits to the data. We found an emission feature at 262.2 GHz that was subsequently masked to determine and subtract the continuum. An emission line map of the feature does not reveal any convincing emission, and we are not aware of any strong molecular or ionized line that could correspond to this redshifted frequency. We checked that masking this feature or not marking it does not impact any of the results that follow. We find that a single Gaussian describes the [C II]$_{158~µm}$ emission well, both in the total and continuum-subtracted spectra, and that the single-Gaussian residuals are consistent with the rms noise (Figure 5). For the total spectrum, the best-fit single-Gaussian and continuum model has a monochromatic continuum flux density of $S_ν = 4.6 \pm 0.2$ mJy, [C II]$_{158~µm}$ observed frequency $ν_{[C II]obs} = 256.158 \pm 0.008$ GHz, [C II]$_{158~µm}$ FWHM = $408_{-23}^{+2}$ km s$^{-1}$, and integrated flux $I_ν = 10.0 \pm 0.8$ Jy km s$^{-1}$. The best-fit double-Gaussian and continuum model has a continuum flux density of $S_ν = 4.5 \pm 0.2$ mJy, $ν_{[C II]obs}^{narrow} = 256.159 \pm 0.009$ GHz, [C II]$_{158~µm}$ FWHM$narrow = 373_{-36}^{+36}$ km s$^{-1}$, and narrow-component integrated flux $I_ν^{narrow} = 8.7 \pm 1.0$ Jy km s$^{-1}$. The broad component has a central frequency $ν_{[C II]obs}^{broad} = 256.19 \pm 0.20$ GHz, [C II]$_{158~µm}$ FWHM$broad = 1330 \pm 820$ km s$^{-1}$, and integrated flux $I_ν^{broad} = 2.9 \pm 3.0$ Jy km s$^{-1}$. The $χ^2$ difference between the single- and double-Gaussian models is only $Δχ^2 = 3.66$ for three additional degrees of freedom, which corresponds to a $p$-value of $p = 0.7$ or evidence at the $0.52σ$ level for a broad component. Equivalently, the Bayesian Information Criterion (BIC) difference is $BIC_1$ –
<figure><img src="image_6.png" /><figcaption>**Figure 3.** Line maps, velocity-integrated over $1.2 \times \text{FWHM}_{[C\,II]} = 482 \text{ km s}^{-1}$ at the expected frequencies of [C II]158 µm, [N II]205 µm, and [O I]146 µm, and assuming the [C II]158 µm redshift. For each line, we also plot two additional collapsed maps, centered at ±482 km s⁻¹ away from the line emission. For [O I]146 µm, only one of the adjacent maps is empty, as the emission is at the edge of the band. On the central [C II]158 µm map, we plot the quadrants used for the spatial extension analysis in Section 4.2. The contours are logarithmic (−4, −2, 2, 4, 8, 16, 32)σ (rms).</figcaption></figure>

$BIC_2 = -12.0$, where a lower $BIC$ implies a better model, and we adopt a threshold $\Delta BIC_{12} > 10$ for strong significance to prefer the more complex model (e.g., Kass & Raftery 1995). We also perform an Anderson-Darling test on the residuals to find any deviations from normally distributed residuals (with a variance equal to the measured rms squared). We find no significant deviation, for either fit, at the $p > 0.15$ level. We thus do not find evidence for a broad spectral [C II]158 µm component in J1148+5251 based on the new NOEMA data alone.

For the continuum-subtracted spectrum, the best-fit single- and double-Gaussian models yield consistent results within the uncertainties (for more details, including the impact of different aperture sizes, see Appendix D), and the residuals are nearly identical (Figure 5, second row, lower right). The $\chi^2$ improvement obtained with a double Gaussian is $\Delta\chi^2 = 2.24$, which for an increase of three free parameters gives a $p$-value of 0.46 for rejecting the simple model in favor of a double Gaussian, in agreement with the analysis of the total spectrum. The BIC
<figure><img src="image_7.png" /></figure>

<figure><img src="image_8.png" /><figcaption>**Figure 4.** Continuum-subtracted spectrum of the [N II]$_{205}$ µm (upper panel) and [O I]$_{146}$ µm (lower panel). The spectra (black lines, shaded blue) are extracted in *r* = 1" apertures centered on the peak of the emission in the velocity-integrated maps (Figure 3) using residual scaling (see Section 2). The 1σ noise level per channel for an *r* = 1" aperture is shown in red, and a single-Gaussian fit is shown in orange (the best-fit parameters and uncertainties are in the upper right corners).</figcaption></figure>

difference is BIC₁ − BIC₂ = − 14.0, e.g., the single-Gaussian model is strongly favored, and an Anderson–Darling test on the residuals also concludes that they are consistent with the measured noise at the *p* > 0.15 level.

We note that, regardless of the aperture chosen or the use of residual scaling, the addition of a broad [C II]$_{158}$ µm component to the spectral fit is not preferred by the NOEMA data (see Appendix D). We have checked that the continuum fluxes in the lower- and higher-frequency sidebands are consistent within ~2%. This ensures that the continuum of J1148+5251, determined from the full ~7.6 GHz sideband, is consistent with that determined solely from the single ~3.8 GHz baseband containing the [C II]$_{158}$ µm line. We also show in Appendix E that the choice of the line masking width and the flux calibration

<table><caption>**Table 2**
Atomic FSL Flux Measurements from the Line Maps
($\Delta\nu = 1.2 \times \text{FWHM}_{\text{C II}} = 482 \text{ km s}^{-1}$), Flux-corrected by 1/0.84</caption><thead><tr><th></th><th>[N II]205 µm</th><th>[C II]158 µm</th><th>[O I]146 µm</th></tr></thead><tbody><tr><td>3" Aperture</td><td></td><td></td><td></td></tr><tr><td>S/N</td><td>3.7</td><td>42.0</td><td>5.4</td></tr><tr><td>S<sub>line</sub>Δν [Jy km s<sup>-1</sup>]</td><td>0.5 ± 0.3</td><td>10.2 ± 0.5</td><td>1.0 ± 0.6</td></tr><tr><td>L<sub>line</sub> [10<sup>9</sup>L<sub>☉</sub>]</td><td>0.4 ± 0.2</td><td>10.6 ± 0.5</td><td>1.1 ± 0.7</td></tr><tr><td>L'<sub>line</sub>[10<sup>9</sup> K km s<sup>-1</sup> pc<sup>-2</sup>]</td><td>4.4 ± 2.1</td><td>48.1 ± 2.5</td><td>4.1 ± 2.4</td></tr><tr><td>Peak [C II] pixel</td><td></td><td></td><td></td></tr><tr><td>S/N</td><td>1.1</td><td>42.0</td><td>5.4</td></tr><tr><td>S<sub>line</sub>Δν [Jy km s<sup>-1</sup>]</td><td>&lt;0.12</td><td>6.4 ± 0.2</td><td>1.0 ± 0.2</td></tr><tr><td>L<sub>line</sub> [10<sup>9</sup>L<sub>☉</sub>]</td><td>&lt;0.10</td><td>6.7 ± 0.2</td><td>1.2 ± 0.2</td></tr><tr><td>L'<sub>line</sub>[10<sup>9</sup> K km s<sup>-1</sup> pc<sup>-2</sup>]</td><td>&lt;1.0</td><td>30.4 ± 0.7</td><td>4.1 ± 0.8</td></tr></tbody></table>

**Note.** Upper limits are given at the 3σ level. The first three lines give the luminosity for the total integrated flux and the last three for the peak [C II]$_{158}$ µm surface brightness. All of the integrated fluxes are derived by applying residual scaling correction, and the luminosities are computed following the definitions of Solomon et al. (1997).

differences between the two ~3.8 GHz basebands do not impact the recovered [C II]$_{158}$ µm line.

The presence of a broad [C II]$_{158}$ µm wing reported in previous works (Maiolino et al. 2012; Cicone et al. 2015) and its absence in the NOEMA data seems puzzling at first. Indeed, Cicone et al. (2015) reach a sensitivity in the [C II]$_{158}$ µm line comparable to ours (0.46 mJy beam⁻¹ versus 0.45 mJy beam⁻¹ per 100 km s⁻¹ channel) and have a similar angular resolution (1".3 × 1".2 versus 1".97 × 1".59). However, their observations were performed with a smaller number of antennas (six, corresponding to 15 baselines, whereas the new observations were done with eight to nine, corresponding to 28–36 baselines), using a correlator that provided a spectral coverage (~3.6 GHz) twice as narrow as that of the new observations (~7.7 GHz) around the [C II]$_{158}$ µm line. Hence, the discrepancy cannot stem from sensitivity or resolution issues, and we expect that, due to our slightly larger beam, larger number of antennas, and higher imaging fidelity, we would be more sensitive, if anything, to extended [C II]$_{158}$ µm structures in J1148+5251 compared to earlier data. The wider spectral coverage of the new Polyfix also enables us to better constrain the dust continuum emission and the exact shape of a potential broad [C II]$_{158}$ µm emission line.

In Appendix F, we discuss in detail the origins of these discrepancies, which can be ascribed to a combination of (i) continuum subtraction and (ii) residual scaling effects. We show that when taking these into account, the previous PdBI and new NOEMA data are consistent. We present the combined data in Appendix F, and find that, using the same approach as used here, broad [C II]$_{158}$ µm line wings are only tentatively recovered in the merged PdBI and NOEMA data set. We finally note that the statistical methodology used here to assess the presence of an outflow is arguably quite conservative. Indeed, with the sensitivity and spectral coverage of the data at hand, only extremely powerful outflows would produce a broad line component whose significance is such that a ξ² or BIC criterion would prefer a double-component fit over a single-component one. In low-redshift studies, it is common to use multiple diagnostics (PV diagrams and other outflow tracers) to assume the presence of an outflow first, and then fit the spectral profile with two components (see, e.g., the discussion in Cicone et al. 2014). In any case, these
<figure><img src="image_9.png" /></figure>

<figure><img src="image_10.png" /></figure>

<figure><img src="image_11.png" /></figure>

<figure><img src="image_12.png" /></figure>

**Figure 5.** The total (top panels) and continuum-subtracted (bottom panels) spectra of the [C II]158 µm line (black), extracted in a $r = 3''$ aperture using residual scaling and one- and two-component Gaussian fits (the orange line and the dashed-dotted blue line, respectively). The thin dashed dark blue lines show the two components of the double-Gaussian model. Both models are fitted to the entire frequency range, with the exception of the overlap between the two sidebands that is coincident with an atmospheric absorption feature (hatched in gray) and the 262 GHz feature (see the text). The right panels show a zoomed-in version of the left plots and the residuals of the two models. The shaded red areas in the bottom of the right panels show the 1σ noise level, computed using the rms noise in each channel in mJy beam⁻¹ with the sigma-clipping scaled accordingly for the number of beams in the $r = 3''$ aperture. The median noise level at $ν_{obs} > 263$ GHz is 2.36 ± 0.14 mJy, whereas a direct rms estimate from the continuum-subtracted spectrum at the same frequencies gives 2.00 ± 0.15 mJy. We find no evidence for a broad [C II]158 µm emission line in the new NOEMA data (see the text for details).

new NOEMA data, thanks to the improved continuum subtraction and analysis methodology, indicate a much less prominent broad component than previously reported, and so rule out the presence of an extremely strong high-velocity outflow. We have, however, provided the two-Gaussian best-fit results above, in case future observations provide ancillary evidence in favor of an outflow.

## 4.2. Spatially Extended [C II]158 µm Emission

The [C II]158 µm emission is spatially extended (Figure 3). In Figure 6, we show the dust and [C II]158 µm radial profiles to investigate their spatial extension. The dust continuum emission (at ~259 GHz) is as extended (within 1σ errors) as the [C II]158 µm emission (up to $θ = 2.51_{-0.25}^{+0.46}$ arcsec, at the 3σ level; see Figure 6, leftmost panel), suggesting that both trace star-

forming material with physical conditions similar to those that we normally ascribe to the ISM of galaxies, rather than an extended [C II]158 µm halo<fnref n="14" /> or outflow. This corresponds to a radius of $r = 9.8_{-2.1}^{+3.3}$ kpc, accounting for beam convolution. This is at least twice as large as any of the 27 quasars observed in [C II]158 µm by Venemans et al. (2020), confirming earlier reports (Maiolino et al. 2012; Cicone et al. 2015) that J1148+5251 is an outlier in terms of [C II]158 µm emission extension. This result holds as well if the dust extension is measured from other spectral setups (e.g., Figure 1).

<footnote n="14">¹⁴ We consider here and hereafter a [C II]158 µm halo to be significant emission beyond that of the dust continuum emission (e.g., Fujimoto et al. 2019, 2020; Herrera-Camus et al. 2021).</footnote>
<figure><img src="image_13.png" /><figcaption>**Figure 6.** [C II]158 µm (red) and dust continuum emission (blue) profiles for J1148+5251, which both extend up to $\theta = 2.51^{+0.46}_{-0.25}$ arcsec at the 3$\sigma$ level, corresponding to $r = 9.8^{+3.3}_{-1.1}$ (accounting for beam convolution). The shaded areas show the corresponding 1$\sigma$ noise level. A point-source profile is shown in black and the dirty (before cleaning and residual scaling) profile is shown in dashed-dotted red. The three panels show, from left to right: the radial profile integrated over the full $2\pi$, only in the NW and SE quadrants, and finally only in the NE and SE quadrants (e.g., along the extended emission seen in Figure 3). In the latter case, the [C II]158 µm emission is more extended than the FIR continuum (5.8$\sigma$ significance).</figcaption></figure>

The extension of [C II]158 µm is asymmetric, however, with a prominent NE-SW axis (Figure 3, although this is less pronounced with a larger channel width; see Figure 8, first panel). In the second and third panels of Figure 6, we compare the dust continuum and the [C II]158 µm radial surface brightness profile for the radially averaged case, the NE-SW axis, and the NW-SE axis (see Figure 3). While we find no evidence for an extended [C II]158 µm halo when averaging radially (or along the NW-SE axis), along the NE-SW axis the [C II]158 µm emission is significantly more extended than the dust continuum (5.8$\sigma$).

We finally explore the extension of the [C II]158 µm emission and the dust continuum directly in the uv plane. To that end, we use the UV_FIT and UV_CIRCLE routines in GILDAS/MAPPING to fit the visibilities with a point source and 2D Gaussian emission model, and then bin the modeled and observed visibilities radially. We plot the real part of the visibilities against the uv radius in Figure 7. In such plots, a point source gives a constant flux density at all uv radii, while a Gaussian emission model yields a Gaussian profile centered at $r=0$. We find good agreement between the observed visibilities and a composite emission model comprising a point source and a 2D Gaussian. We fit both a circular Gaussian and an elliptical Gaussian to the dust and [C II]158 µm continuum visibilities. For the [C II]158 µm emission, the best-fit model gives a corrected velocity-integrated flux of $4.2 \pm 0.5$ Jy km s$^{-1}$ in the point-source component, and $5.0 \pm 0.5$ Jy km s$^{-1}$ in the extended Gaussian, in good agreement with the flux derived from the line map (Table 2) and the fitted spectrum (see Appendix D). The elliptical Gaussian models yield velocity-integrated fluxes of $5.0 \pm 0.3$ Jy km s$^{-1}$ and $5.5 \pm 0.4$ Jy km s$^{-1}$ for the point-source and resolved components, respectively. In either case, we are consistent with the large fraction of emission coming from the unresolved component reported in Cicone et al. (2015), but not with the overall flux, which is possibly due to continuum subtraction differences, as discussed in Appendix F. Additionally, the unresolved flux is consistent with that measured at higher resolution (Walter et al. 2009a), suggesting the difference in total flux is due to the extended component being resolved out in high-resolution configurations.

<figure><img src="image_14.png" /><figcaption>**Figure 7.** Observed (black) and modeled (blue and orange) real parts of the visibilities for the [C II]158 µm emission (integrated over $\Delta v = 482$ km s$^{-1}$) and dust continuum in J1148+5251. The model points are slightly offset for presentation purposes. The difference between the two models (blue: point source + circular Gaussian; orange: point source + elliptical Gaussian) is not statistically significant ($p$-value = 0.83; e.g., 0.96$\sigma$ significance).</figcaption></figure>
<figure><img src="image_15.png" /><figcaption>**Figure 8.** From left to right: a map of the [C II]158 µm emission, averaged over the velocity ranges (–1400, 1200) km s⁻¹, (–1400, 300) km s⁻¹, and (+400, +1200) km s⁻¹. These velocity ranges are similar to those presented in Cicone et al. (2015), who reported 6σ [C II]158 µm clumps at 2''–4'' from the source center in the blue-/redshifted wings of the [C II]158 µm emission. Such spatially and spectrally offset emission is not recovered in the NOEMA data. The 1σ rms values in the above three maps are 0.27, 0.14, and 0.15 mJy beam⁻¹. This figure shows contours in logarithmic increments (–4, –2, 2, 4, 8, 16, 32)σ (rms).</figcaption></figure>

The FWHM of the [C II]158 µm Gaussian component is FWHM = 1''6 ± 0''2(9 ± 1 kpc), or minor/major axes $a = 3''6 ± 0''3(20 ± 2 kpc)$ and $b = 1''8 ± 0''2(10 ± 1 kpc)$ for the elliptical model. This is in agreement with the scale of the 3σ extension directly measured on the cleaned image. The elliptical Gaussian has a major/minor axis difference and angle (PA = 53 ± 4 deg) in agreement with the asymmetric [C II]158 µm emission discussed above. However, a likelihood ratio test does not prefer the elliptical Gaussian model over the circular one (p-value = 0.83; e.g., 0.96σ significance). We do not find any difference between the elliptical Gaussian and circular Gaussian models for the dust emission (see Figure 7). The best-fit model for the dust emission gives a flux density of 1.0 ± 0.7 mJy and 3.3 ± 0.7 mJy for the unresolved/ resolved components, respectively, with the sum in agreement with the measurement of the aperture-integrated flux density (Table 1). The dust uv profile (Figure 7, bottom panel) does not show clear evidence of a flattening at large uv distances that is characteristic of a point source. The Gaussian component has an FWHM of $r = 0''7 ± 0''1(3.9 ± 0.6 kpc)$.

Cicone et al. (2015) reported a complex velocity structure of the [C II]158 µm outflow, with emission clumps detected up to $r \sim 30$ kpc and $\sim 1000$ km s⁻¹ offsets from the central emission, particularly visible in an emission line map integrated over a (–1400, 1200) km s⁻¹ interval (we show in Appendix F how we recover these structures in the previous PdBI data, although at a lower significance level). Such features are not recovered in the new NOEMA data when averaging over a similar velocity interval (Figure 8, left panel; see also Appendix F). The [C II]158 µm emission is mostly concentrated within 3'' (~11 kpc), and no significant emission is recovered beyond that radius.<fnref n="15" /> We further show in Appendix G that no significant velocity structure is detected in J1148+5251 in the NOEMA data, neither by inspecting the channel maps (with channel width $\Delta v = 120$ km s⁻¹) nor by means of a kinematical analysis.

# 5. ISM Properties

## 5.1. [C II]158 µm Emission

In the following analysis, we adopt luminosities measured from the line maps. We measure a [C II]158 µm line luminosity $L_{[C II]} = (10.6 \pm 0.5) \times 10^8 L_{\odot} (L'_{[C II]} = (48.1 \pm 2.5) \times 10^9 \text{ K km s}^{-1}\text{pc}^{-2})$. Following the [C II]–SFR relation of De Looze et al. (2014) for high-redshift ($z > 0.5$) galaxies, we find SFR = 2041 ± 111 $M_{\odot}$ yr⁻¹, in good agreement with the SFR inferred from the continuum luminosities. Using the De Looze et al. (2014) calibration for AGNs does lower the result (SFR = 857 ± 35 $M_{\odot}$ yr⁻¹), and that for ULIRGS increases the inferred SFR (SFR = 5562 ± 257 $M_{\odot}$ yr⁻¹). This agreement between the [C II]158 µm and the infrared luminosity-based SFR is also observed in other high-redshift quasars (e.g., Novak et al. 2019; Venemans et al. 2020) where the FIR luminosity is smaller. We caution that the SFR–[C II]158 µm relations used above are dependent on the FIR luminosity and the FIR surface brightness of the galaxy populations used to calibrate them (e.g., Díaz-Santos et al. 2017; Herrera-Camus et al. 2018b). These SFR values should only be used to compare with other high-redshift sources where the FIR is not well constrained. For J1148+5251, the [C II]/FIR ratio is $\approx (7.8 \pm 1.3) \times 10^{-4}$. We thus use the appropriate “high $\Sigma_{\text{FIR}}$” [C II]158 µm–SFR calibration from Herrera-Camus et al. (2018b) and find SFR = 3100 ± 1600 $M_{\odot}$ yr⁻¹, in agreement with the other estimates presented above.

Following Weiß et al. (2003, 2005), in the optically thin limit, the total mass of ions X emitting photons with rest-frame frequency $\nu_{ij}$ stemming from a transition between the i and j levels can be derived from the observed luminosity as:

$$M_X = m_X \frac{8\pi k \nu_{ij}^2}{hc^3 A_{ij}} Q(T_{\text{ex}}) \frac{1}{g_i} e^{-T_{\text{ex}}/T_i} L'_{[X_{i \to j}]}, \quad (1)$$

where k is the Boltzmann constant, c is the speed of light, h is the Plank constant, $m_X$ is the mass of a single ion X, $Q(T_{\text{ex}}) = \sum_i g_i e^{-T_i/T_{\text{ex}}}$ is the partition function of the species, with $g_i$ being the statistical weight of level i, $T_i$ being the energy of (above-ground) level i, $T_{\text{ex}}$ being the excitation temperature of the i→j transition, and $L'_{[X_{ij}]}$ being the observed integrated source brightness temperature of the line in K km s⁻¹ pc⁻². Note

<footnote n="15">¹⁵ Following the analysis of Cicone et al. (2015), we also plot the blue- and redshifted components ((–1400, –300) km s⁻¹ and (+400, +1200) km s⁻¹) of the [C II]158 µm emission, and find only marginal detections (3.0 and 3.4σ) that are not spatially offset (Figure 8).</footnote>
that Equation (1) neglects heating from the CMB, which is negligible for the [O I]146 µm, [C II]158 µm, and [N II]205 µm transitions, but not [C I]369 µm at z ~ 6.4. For C+, Equation (1) reduces to:

$$M_{C^+}/M_{\odot} = 2.92 \times 10^{-4} Q(T_{\text{ex}}) \frac{1}{4} e^{-91.2/T_{\text{ex}}} L'_{[\text{C II}]}, \quad (2)$$

where $Q(T_{\text{ex}}) = 2 + 4e^{91.2/T_{\text{ex}}}$ is the [C II]158 µm partition function. The optically thin limit assumption, while widespread in the literature, is uncertain at higher redshift, since optical depth measurements, although pointing to moderate values, are scarce (e.g., Neri et al. 2014; Gullberg et al. 2015). Lagache et al. (2018) and Vallini et al. (2015) note that in the optically thick limit, only emission from PDRs would reach the observer. As we will show in Section 5.4, most of the [C II]158 µm emission in J1148+5251 indeed comes from PDRs, and hence we do not expect the optically thin limit assumption to affect our results. Assuming an excitation temperature $T_{\text{ex}} = 50$ K (from the CO modeling; Riechers et al. 2009), we find $M_{C^+} = (5.8 \pm 0.3) \times 10^7 M_{\odot}$. This is five times the neutral carbon mass ($M_C = 1.1 \times 10^7 M_{\odot}$) measured by Riechers et al. (2009).

## 5.2. [N II]205 µm Emission

We report an [N II]205 µm marginal detection in J1148+5251 at $S/N = 3.7$ and an integrated ($r = 3''$) luminosity $L_{[N II]} = (0.4 \pm 0.2) \times 10^9 L_{\odot}$. This is in slight tension with the $3\sigma$ limit ($< 0.4 \times 10^9 L_{\odot}$) reported by Walter et al. (2009b) using the IRAM 30 m telescope. Nonetheless, we would expect better image fidelity with the new NOEMA facility. We consider the [N II]205 µm detection marginal ($3.7\sigma$), and urge caution in interpretations that rely upon it.

The marginal ($1''/59/8.9$ kpc ~ 1 beam) offset between [N II]205 µm and [C II]158 µm in J1148+5251 could be of interest. On the one hand, spatial offsets between low- and high-ionization lines have been reported in other high-redshift quasars, most often between [C II]158 µm and [O III]88 µm (e.g., Novak et al. 2019). Spatial offsets have also been predicted in theoretical simulations (Katz et al. 2017, 2019), where they arise from different gas phases with different temperatures and densities. Although nitrogen and carbon have a similar ionization level, Katz et al. (2019) show that the [N II]205 µm and [C II]158 µm can arise from different gas phases, with [C II]158 µm originating in lower-temperature and higher-density regions,<fnref n="16" /> which could explain the offset of [N II]205 µm toward the outskirts of J1148+5251. On the other hand, offsets often seen at low resolution in high-redshift galaxies can disappear in higher-resolution data once fainter emission components are detected (see, e.g., HZ10 in Pavesi et al. 2016, 2019). Following Ferkinhoff et al. (2010, 2011), we derive the minimum $H^+$ for the [N II]205 µm luminosity observed. To do so, we assume high densities and a high temperature (as found around O and B stars), such that all nitrogen in the H II regions is ionized. We use Equation (1) to derive the mass of $N^+$ in J1148+5251 from the observed luminosity, with $A_{10} = 2.1 \times 10^{-6}$ being the Einstein coefficient of the $^3P_1 \rightarrow ^3P_0$ transition, $g_1 = 3$ being the statistical weight of the $^3P_1$ emitting level, $\nu_{10} = 1461.1$ GHz being the rest-frame frequency, and $g_t \simeq 9$ being the partition function. We can then derive the minimum $H^+$ mass by

assuming the upper limit on the ionized nitrogen to ionized hydrogen ratio (i.e., $\chi(N^+) = N^+/H^+$) to be the total nitrogen abundance ratio $\chi(N)$, such that

$$M(H^+) \ge M_{N^+} - \frac{m_N}{m_H \chi(N^+)} = \frac{L_{[N II]205 \mu m}}{\frac{g_1 A_{10} h \nu_{10}}{g_t}} \frac{m_H}{\chi(N^+)}. \quad (3)$$

We adopt the abundance value for H II regions $\xi(N) = 9.3 \times 10^{-5}$ from Savage & Sembach (1996). Hence, we estimate an ionized hydrogen mass $M(H^+) \ge (2.1 \pm 1.0) \times 10^9 M_{\odot}$ ($2\sigma$ level). Using the H$_2$ gas mass from the CO luminosity (Riechers et al. 2009), the ionized to molecular gas ratio is $M(H^+)/M(H_2) > 0.1$. This is significantly higher than what is found by Ferkinhoff et al. 2011 (based on the compilation by Brauher et al. 2008) for local galaxies using the [N II] 122 µm line ($<0.01$). We caution here again that the [N II]205 µm line is marginal.

The [C II]/FIR ratio ($\simeq (7.8 \pm 1.3) \times 10^{-4}$), [N II]/FIR ratio ($\sim (2.9 \pm 1.5) \times 10^{-5}$), and dust temperature ($53 \pm 8$ K) are all in agreement with the ratio trends (extrapolated to the high temperature observed in J1148+5251) observed in local ULIRGS (Díaz-Santos et al. 2017) and high-redshift galaxies or quasars (De Breuck et al. 2019; Novak et al. 2019; Pavesi et al. 2019; Li et al. 2020; Pensabene et al. 2021; see Figure 9). We note that the [C II]/[N II] luminosity ratio ($10.6 \pm 3.7$) is slightly low for the dense PDRs/X-ray-dominated regions (XDRs) expected around high-redshift quasars (e.g., Decarli et al. 2014; Pavesi et al. 2019), though the uncertainties are large given the marginal detection in the [N II]205 µm line. Finally, we also test the Zhao et al. (2016), [N II]205 µm-SFR scaling relation derived from local ULIRGS, which gives SFR = 794 ± 274 $M_{\odot}$ yr$^{-1}$, lower than the [C II]158 µm and FIR-derived values. We note, however, that [N II]205 µm is a poor tracer of the SFR in intense star-forming regions, such as J1148+5251, due to its low critical density (e.g., Herrera-Camus et al. 2016; Zhao et al. 2016).

## 5.3. [O I]146 µm Emission

The peak of the [O I]146 µm emission is located at the same position as the [C II]158 µm emission. The [O I]146 µm/FIR ratio is $1.3 \times 10^{-5}$, in good agreement with that of local ULIRGS with comparable FIR surface flux density (Herrera-Camus et al. 2018a). The [O I]146 µm/[C II]158 µm line luminosity ratio $L_{[O I]}/L_{[C II]} = 0.10 \pm 0.07$ is typical of that observed in high-redshift submillimeter galaxies and quasars (e.g., De Breuck et al. 2019; Yang et al. 2019; Li et al. 2020; Lee et al. 2021), and is at the higher end of the range spanned by local AGNs (e.g., Fernández-Ontiveros et al. 2016).

Again using Equation (1) for the [O I]146 µm $^3P_0 \rightarrow ^3P_1$ transition ($T_0 = 329$ K, $g_0 = 1$, $Q(T_{\text{ex}}) = 5 + 3e^{-329/T_{\text{ex}}} + e^{-228/T_{\text{ex}}}$), we derive a neutral oxygen mass

$$M_{\text{O}}/M_{\odot} = 6.19 \times 10^{-5} Q(T_{\text{ex}}) e^{329/T_{\text{ex}}} L'_{[\text{O I}]}. \quad (4)$$

Assuming an excitation temperature $T_{\text{ex}} = 50$ K (from the CO modeling; Riechers et al. 2009), we find $M_{\text{O}} = (9 \pm 5) \times 10^8 M_{\odot}$. This estimate is extremely sensitive to the excitation temperature, e.g., ranging from $0.7 \times 10^7 M_{\odot}$ at 200 K to $8 \times 10^8 M_{\odot}$ at $T_{\text{ex}} = 50$ K. Note that the above expression only applies in the optically thin limit, which is probably not the case for [O I]146 µm, therefore it underestimates the neutral oxygen mass.

<footnote n="16">¹⁶ Although a minor fraction of [C II]158 µm traces the same H II regions as [N II]205 µm; see Section 5.4.</footnote>
<figure><img src="image_16.png" /><figcaption>**Figure 9.** [C II]$_{158}$ µm and [N II]$_{205}$ µm line deficits of local (U)LIRGS (shaded gray; Díaz-Santos et al. 2017), the high-redshift galaxy SPT0418-47 (the purple triangle; De Breuck et al. 2019), and high-redshift quasars (the orange pentagon, dark orange star, and dark red circles; Novak et al. 2019; Li et al. 2020; Pensabene et al. 2021, respectively) and companion galaxies (the dark red circles; Pensabene et al. 2021). Our measurements for J1148+5251 are shown via the red squares.</figcaption></figure>

## 5.4. CLOUDY Modeling of the FSL Ratios

We now use the FSL ratios to determine the physical properties of the ISM of J1148+5251. In order to do so, we make use of CLOUDY (Ferland et al. 2017), a spectral synthesis code designed to simulate the spectra of astrophysical plasmas used to study both low- and high-redshift galaxies’ submillimeter lines. The grid of models used in this work was generated for studying both high-redshift quasar hosts and their companions, by Pensabene et al. (2021), to which we refer for further details. To summarize briefly, the model grid includes both PDR and XDR predictions for total hydrogen column densities $N_H/[\text{cm}^{-2}] = 10^{23}-10^{24}$, total hydrogen number density $1 < \log n/[\text{cm}^{-3}] < 6$, and local far-UV radiation field $2 < \log G/[G_0] < 6$ in units of Habing flux (for the PDRs) or X-ray flux $-2 < \log F_X/[\text{erg s}^{-1} \text{cm}^{-2}] < 2$ (for the XDRs). For each combination of parameters, the fluxes of the various lines of interest are predicted, and the ISM conditions can be determined by comparing the observed line ratios.

Whereas [C II]$_{158}$ µm can originate from both the neutral and ionized gas phase, [O I]$_{146}$ µm and [C I]$_{369}$ µm originate solely in the neutral gas phase/PDRs. It is therefore crucial to determine which fraction of the [C II]$_{158}$ µm emission in J1148+5251 comes from the neutral gas. [N II]$_{205}$ µm and [C II]$_{158}$ µm have very close critical densities in ionized media, therefore the [C II]$_{158}$ µm to [N II]$_{205}$ µm ratio is primarily a function of the $N^+/C^+$ abundance ratio (e.g., Oberst et al. 2006).

Photoionization models (e.g., Oberst et al. 2006; Pavesi et al. 2016; Croxall et al. 2017) predict a relatively constant ratio (2.5–3) for [N II]$_{205}$ µm/[C II]$_{158}$ µm for a large range of electron densities. For consistency with the existing literature, we adopt [N II]$_{205}$ µm/[C II]$_{158}$ µm = 3. Any deviation from that ratio can be interpreted as [C II]$_{158}$ µm emission from the neutral phase. The fraction of [C II]$_{158}$ µm emission originating in the neutral phase is therefore

$$f([\text{C II}], \text{neutral}) \approx 1 - 3 \frac{L_{[\text{N II}]_{205 \text{ µm}}}}{L_{[\text{C II}]_{158 \text{ µm}}}} \quad (5)$$

The ratio of the luminosities measured in this work for J1148+5251 is $27 \pm 13$ for the aperture-integrated ($r=3''$) flux and $>100$ ($2\sigma$ level) at the peak of the [C II]$_{158}$ µm emission. In either case, most of the [C II]$_{158}$ µm emission (88%–97%) comes from the neutral phase.

Given the resolution of the NOEMA data, we use the $r=3''$ aperture-integrated line luminosities to compute line luminosity ratios and compare them to the grid of CLOUDY models. In order to properly include the [C I]$_{369}$ µm measurements from Riechers et al. (2009), we have repeated the line map procedure detailed in Section 2 on the original [C I]$_{369}$ µm data. We produce a line map with width $\Delta v = 482 \text{ km s}^{-1}$, and integrate the flux up to $r=3''$, finding $S_{[\text{C I}]_{369} \text{ µm}}\Delta v = 0.45 \pm 0.16 \text{ Jy km s}^{-1}$ and a line luminosity $L_{[\text{C I}]_{369} \text{ µm}} = (0.20 \pm 0.07) \times 10^9 L_{\odot}$. The new luminosity value is in agreement, within $1.4\sigma$, with the previously published value of $(0.10 \pm 0.02) \times 10^9 L_{\odot}$ (Riechers et al. 2009). The difference between the two measurements can be explained by the larger aperture, and is comparable to that seen for [C II]$_{158}$ µm (see Appendix A).

We show the radiation field and density predictions of CLOUDY for the observed FSL ratios in Figures 10 and 11. We find that the different line ratios are not perfectly reproduced for a single density and radiation field by the chosen grid of models. However, the [C II]$_{158}$ µm/[C I]$_{369}$ µm luminosity ratio ($47 \pm 3$) excludes XDR models, which cannot reproduce such high ratios (the maximum possible being ~15; Figure 10; e.g., Venemans et al. 2017a, 2017b; Novak et al. 2019; Pensabene et al. 2021). Besides, our XDR model grid cannot reproduce all of the FSL to FIR luminosity ratios with a single hydrogen number density and X-ray flux. We note that this result, based on the analysis of the FSLs, is in tension with Gallerani et al. (2014), who concluded that an XDR component is needed in J1148+5251 to explain their reported CO(17–16) line. However, as noted by the authors, this line is potentially contaminated by a nearby $OH^+$ emission line, and thus additional observations of high-$J$ CO lines are needed to clarify the situation.

The PDR models in agreement with the observed ratios (except [O I]$_{146}$ µm/[C I]$_{369}$ µm) have an H I column density of $N_{H1} = 10^{23} \text{ cm}^{-2}$, high radiation fields ($10^{3.5-4.5} G_0$), and moderate hydrogen number densities ($n \sim 10^{3.5-4.5} \text{ cm}^{-3}$), commensurable with other studies of high-redshift quasars (e.g., Novak et al. 2019; Pensabene et al. 2021). The model grids do not reproduce the observed [N II]$_{205}$ µm/[C II]$_{158}$ µm ratio, since it was not created to model the ionized phase (H II regions) traced by [N II]$_{205}$ µm emission. The [O I]$_{146}$ µm/[C I]$_{369}$ µm discrepancy could stem from the different quality and resolution of the data, as well as the low S/N of both lines. Another possibility is that a fraction of the [C II]$_{158}$ µm emission could come from collisional excitation in the cold neutral medium, rather than solely from the PDRs, as is the case for [O I]$_{146}$ µm and [C I]$_{369}$ µm, leading to a small discrepancy between the ratios
<figure><img src="image_17.png" /><figcaption>**Figure 10.** ISM constraints for J1148+5251 from the FSL ratios predicted by CLOUDY for XDRs with hydrogen column densities $N_H/[\text{cm}^{-2}] = 10^{23}$ (top panel) and $N_H/[\text{cm}^{-2}] = 10^{24}$ (lower panel). The observed ratio constraints from the aperture-integrated fluxes are plotted in solid lines, and the $\pm 1\sigma$ values in the shaded areas of the same color. The ratios absent from the plots (such as the $[\text{C II}]_{158\,\mu\text{m}}$ ratios) are not reproduced by the XDR grid of models.</figcaption></figure>

involving the $[\text{C II}]_{158\,\mu\text{m}}$ line and others in our PDR-only models. This hypothesis is supported by the spatial extent of the $[\text{C II}]_{158\,\mu\text{m}}$ line emission (Section 4.2), which is beyond what is normally observed for the disk component of $z \sim 6$ massive star-forming galaxies and other quasar hosts (e.g., Venemans et al. 2020). Deeper and higher-resolution observations of the inner and outer $[\text{C II}]_{158\,\mu\text{m}}$-emitting regions, combining the presence or absence of $[\text{O I}]_{146\,\mu\text{m}}$, $[\text{N II}]_{205\,\mu\text{m}}$, and $[\text{C I}]_{369\,\mu\text{m}}$, would be required to further investigate the ISM of J1148+5251.

## 6. Conclusions

We report new NOEMA observations of the $z = 6.42$ quasar J1148+5251 in the atomic FSLs of $[\text{C II}]_{158\,\mu\text{m}}$, $[\text{O I}]_{146\,\mu\text{m}}$, and $[\text{N II}]_{205\,\mu\text{m}}$, and the underlying dust continuum emission. The high-fidelity data, together with the large instantaneous bandwidths that NOEMA provides, enabled us to revisit the properties of the quasar’s host galaxy and derive the physical

<figure><img src="image_18.png" /><figcaption>**Figure 11.** ISM constraints from the FSL ratios predicted by CLOUDY for PDRs with hydrogen column densities $N_H/[\text{cm}^{-2}] = 10^{23}$ (top panel) and $N_H/[\text{cm}^{-2}] = 10^{24}$ (lower panel). The observed ratio constraints from the aperture-integrated fluxes are plotted in solid lines, and the $\pm 1\sigma$ values in the shaded areas of the same color.</figcaption></figure>

conditions of its ISM. The main conclusions of this paper are as follows:

1. A uv plane analysis confirms the presence of an extended $[\text{C II}]_{158\,\mu\text{m}}$ and dust emission component (FWHM = 1.6 $\pm$ 0.2 (9 $\pm$ 1 kpc)), accounting for ~50%–60% of the total $[\text{C II}]_{158\,\mu\text{m}}$ and most of the dust emission, in agreement with earlier studies (Maiolino et al. 2012; Cicone et al. 2015). J1148+5251 thus remains an outlier with very extended $[\text{C II}]_{158\,\mu\text{m}}$ and dust emission compared to other $z > 6$ quasars (e.g., Venemans et al. 2020).

2. The $[\text{C II}]_{158\,\mu\text{m}}$ emission shows a similar spatial extent as the FIR dust when averaged radially (up to $\theta = 2.51_{-0.25}^{+0.46}$ arcsec, corresponding to $r = 9.8_{-2.1}^{+3.3}$ physical kpc, accounting for the beam convolution). However, if the $[\text{C II}]_{158\,\mu\text{m}}$ emission is examined only along its NE–SW axis, a significant $(5.8\sigma)$ excess $[\text{C II}]_{158\,\mu\text{m}}$ emission (w.r.t. to the dust) is detected.
3. The [C II]158 µm line profile can be fitted with a single Gaussian with a FWHM of 408 ± 23 km s⁻¹. The new NOEMA data has enabled us to reassess the significance of the broad [C II]158 µm line wings reported in previous studies. We find that an additional broad [C II]158 µm component is not statistically required to describe the NOEMA data (or the merged NOEMA and previous PdBI data).

4. We report the detection of [O I]146 µm and [N II]205 µm (tentatively) in J1148+5251. Using various empirical relations, we report a C⁺ mass of M<sub>C⁺</sub> = (2.7 ± 0.1)×10<sup>7</sup>M<sub>☉</sub>, an oxygen mass of M<sub>O</sub> = (12 ± 4) × 10<sup>8</sup>M<sub>☉</sub>, and a lower limit on the ionized hydrogen mass M<sub>H⁺</sub>> (4.6 ± 0.16) × 10<sup>8</sup>M<sub>☉</sub> (2σ level).

5. The FSL line ratios are consistent with the trends observed in local (U)LIRGs. We find that a large fraction (~90%) of the [C II]158 µm emission originates in the neutral phase (PDR) of the gas.

6. We have compared CLOUDY models to the observed FSL ratios in J1148+5251. The [C II]158 µm/[C I]369 µm ratio and the multiple FSL to FIR ratios exclude XDR models in favor of PDRs. We find good agreement for models that have a high radiation field (10<sup>3.5-4.5</sup> G₀), moderate hydrogen number densities (n ≈ 10<sup>3.5-4.5</sup> cm<sup>-3</sup>), and H I column density N<sub>H I</sub> = 10<sup>23</sup> cm<sup>-2</sup>.

The results described here highlight the importance of large instantaneous bandwidths when observing high-redshift quasars (or galaxies) to search for weak extended emission of atomic or molecular lines. Our findings enabled a renewed view on the host galaxy of the J1148+5251 quasar, shedding light on the feedback activity and providing new constraints on the excitation conditions of its ISM that appear similar to what is found in local ULIRGs. Higher-angular-resolution and higher-sensitivity data, in particular for the [C I]369 µm emission line, would be required to put additional constraints on the gas density and the ionization source of J1148+5251, and further explore its properties.

The authors thank the anonymous referee for useful suggestions that improved the manuscript. R.A.M., F.W., and M.N. acknowledge support from the ERC Advanced Grant 740246 (Cosmic_Gas). D.R. acknowledges support from the National Science Foundation under grant Nos. AST-1614213 and AST-1910107. D.R. also acknowledges support from the Alexander von Humboldt Foundation through a Humboldt Research Fellowship for Experienced Researchers.

This work is based on observations carried out with the IRAM
Northern Extended Millimeter Array (NOEMA). IRAM is
supported by INSU/CNRS (France), MPG (Germany), and IGN
(Spain). The authors would like to wholeheartedly thank the
staff at IRAM for their help with the reduction and analysis
software GILDAS; and in particular M. Krips, V. de Souza, C.
Herrera, and J.-M. Winters.

Facility: NOEMA (IRAM).

Software: astropy (The Astropy Collaboration et al. 2018), CLOUDY (Ferland et al. 2017), Numpy (Harris et al. 2020), Scipy (Virtanen et al. 2020), Matplotlib (Hunter 2007), Interferopy (Boogaard, Meyer & Novak 2021).

Appendix A
Curve of Growth Analysis

In this appendix, we investigate what is the size of aperture
that is necessary to recover most of the continuum or [C II]158 µm
line fluxes. Figure 12 shows the line flux density of [C II]158 µm,
[O I]146 µm, and [N II]205 µm, as well as the continuum flux
density, as a function of increasing aperture radius. We find that
all of the line/continuum fluxes reach a maximum or plateau at
an aperture radius *r* = 3'', which corresponds to 16.9 kpc at
*z* = 6.42. The [C II]158 µm flux presents tentative evidence for
additional flux up to 4'' (~1 Jy km s⁻¹), but this is within the 1σ
errors. Note that at large radii, where there is no more cleaned
flux, residual scaling can become numerically unstable.
Throughout this paper, an aperture of *r* = 3'' is therefore adopted
for measurements, unless specified otherwise.

<figure><img src="image_19.png" /><figcaption>**Figure 12.** FSL flux densities (left panel) and integrated continuum flux (right panel) as a function of aperture radius. The solid and dashed lines show the fluxes with and without residual scaling correction (see Section 2), respectively. The final aperture radius of 3'' was chosen to encompass all of the [C II]158 µm emission and the continuum at the higher frequencies.</figcaption></figure>
## Appendix B
### Multiscale and Högbom Cleaning Methods

In this appendix, we briefly detail our experiments with different cleaning methods for our interferometric data. Högbom cleaning (Högbom 1974) is one of the standard methods for cleaning interferometric data. It relies on iteratively finding peaks in the data, and subtracting the dirty beam at that location, until the residuals reach a desired level. It is particularly efficient for point sources, but struggles with large-scale emission, which it tries to reconstruct using a multitude of point sources. In that case, the so-called multiscale algorithms, which convolve the beam with various Gaussians to subtract larger scales, are preferable (e.g., Wakker & Schwarz 1988).

Figure 13 shows the clean map, the dirty map, and the residuals for the Högbom and multiscale clean on the [C II]158 µm

map (integrated over 482 km s⁻¹). Clearly, the Högbom clean residuals show a flat excess of 2σ flux filling the 3" aperture, which would not be expected if the source was a single (or a limited number of) point source(s). On the contrary, the residuals of the multiscale algorithm are closer to zero on average. Therefore, we chose the multiscale algorithm for all of the [C II]158 µm-derived quantities and images in this paper.

Figure 14 shows the clean, dirty, and residual maps for the continuum maps and the other FSL maps for a Högbom clean. For the lower-frequency continuum maps, the [O I]146 µm maps, and the [N II]205 µm maps, the residuals are well-behaved and do not require the use of multiscale cleaning. For the 259 and 274 GHz continuum maps, some residuals are seen, and the multiscale algorithm is adopted for those in the paper, albeit changing only the final continuum flux by ≲2%.

<figure><img src="image_20.png" /><figcaption>**Figure 13.** Clean, dirty, and residual maps of the [C II]158 µm emission for two different cleaning algorithms: Högbom (the first row) and multiscale (the second row). Both are cleaned down to 2σ, where σ is the rms noise of the dirty map.</figcaption></figure>
<figure><img src="image_21.png" /><figcaption>**Figure 14.** Clean, dirty, and residual maps for a Högbom clean to $2\sigma$ of the 200 GHz continuum, [N II]$_{205\,\mu m}$ emission, the 216, 259, and 274 GHz continua, and finally the [O I]$_{146\,\mu m}$ emission.</figcaption></figure>

## Appendix C
### Dust SED Parameter Posterior Distribution

We present in Figure 15 the posterior distribution of the dust SED parameters fitted in Section 2. The median dust properties

derived are consistent with the existing literature on high-redshift quasars (e.g., Venemans et al. 2020).

<figure><img src="image_22.png" /><figcaption>**Figure 15.** The dust SED-fitting posterior distribution of the dust temperature *T*, dust mass *M*, and the dust spectral emissivity index *β*.</figcaption></figure>
## Appendix D
[C II]158 µm Flux Density

In Table 3, we compare the [C II]158 µm fluxes measured from the extracted aperture-integrated [C II]158 µm spectra for various aperture sizes between previous studies and this work. This work's spectra and best-fit Gaussians are displayed in Figure 16 for circular apertures with radii $r = 1''$, $2''$, and $4''$, and a central-pixel-only spectrum. The $r = 3''$ case has already been presented in the main text (Figure 5). For Cicone et al. (2015), we report

their best-fit parameters and fluxes for the “narrow” (–200, 200) km s⁻¹ component defined in their study. We find good agreement for all apertures between this work and previous studies. This implies that most of the excess flux in larger apertures ($r \gtrsim 3''$) described in Cicone et al. (2015) is due to the reported presence of blue-/redshifted components with large spatial offsets that are not recovered in this work, as discussed in Section 4.

<figure><img src="image_23.png" /><figcaption>**Figure 16.** Continuum-subtracted aperture-integrated [C II]158 µm spectra (black) for circular apertures with radii $r = 1''$, $2''$, and $4''$. The best-fit Gaussian is shown in orange, and the best-fit parameters are displayed in the upper right corners of each plot.</figcaption></figure>

<table><caption>**Table 3**
[C II]158 µm Line Fluxes Measured from a Gaussian Fit to the Aperture-integrated Spectra in This Work and Previous Studies</caption><thead><tr><th rowspan="2"></th><th colspan="2">r = 1"</th><th colspan="2">r = 2"</th><th colspan="2">r = 3"</th><th colspan="2">r = 4"</th></tr><tr><th>Cicone et al.<br>(2015)</th><th>This work</th><th>Cicone et al.<br>(2015)</th><th>This work</th><th>Maiolino et al.<br>(2012)</th><th>This work</th><th>Cicone et al.<br>(2015)</th><th>This work</th></tr></thead><tbody><tr><td>σ<sub>v</sub> [km s<sup>-1</sup>]</td><td>146 ± 11</td><td>158 ± 6</td><td>148 ± 16</td><td>162 ± 4</td><td>150 ± 20</td><td>171 ± 10</td><td>150 ± 20</td><td>175 ± 16</td></tr><tr><td>S<sub>peak</sub> [mJy]</td><td>14.5 ± 0.9</td><td>9.41 ± 0.32</td><td>30 ± 3</td><td>20.1 ± 0.7</td><td>23 ± 2</td><td>23.5 ± 1.2</td><td>34 ± 4</td><td>24.7 ± 2.0</td></tr><tr><td>I<sub>ν</sub> [Jy km s<sup>-1</sup>]</td><td>5.3 ± 0.5</td><td>3.7 ± 0.2</td><td>11.0 ± 1.5</td><td>8.1 ± 0.4</td><td>14 ± 3</td><td>10.0 ± 0.8</td><td>13 ± 3</td><td>11.0 ± 1.0</td></tr></tbody></table>

**Note.** For each aperture radius and study, we give the best-fit velocity width $\sigma_v$, the peak line flux density $S_{\text{peak}}$, and the integrated flux $I_{\nu}$. For the Cicone et al. (2015) values, only the narrow-component results ([C II]158 µm emission integrated between (–200, 200) km s⁻¹) are reported. The Walter et al. (2009a) best-fit values are: $\sigma_v = 122 \pm 12$ km s⁻¹, $S_{\text{peak}} = 12.7 \pm 1.1$ mJy, and $I_{\nu} = 3.9 \pm 0.3$ Jy km s⁻¹. Similarly, the Maiolino et al. (2005) IRAM 30 m measurement gives: $\sigma_v = 149 \pm 21$ km s⁻¹, $S_{\text{peak}} = 11.8$ mJy, and $I_{\nu} = 4.1 \pm 0.5$ Jy km s⁻¹.
## Appendix E

### Continuum Subtraction and Masking of the [C II]$_{158}$ µm Line

In this appendix, we provide details about the continuum subtraction procedures, focusing on the [C II]$_{158}$ µm spectral setup. Figure 17 shows the aperture-integrated [C II]$_{158}$ µm

spectrum for various half-width masking regions, ranging from 4.5 × FWHM([C II]) = 1733 km s⁻¹ to 0.5 × FWHM ([C II]) = 193 km s⁻¹. In all cases, we fit single- and double-Gaussian models to the [C II]$_{158}$ µm emission (as in Section 4.1)

<figure><img src="image_24.png" /></figure>

<figure><img src="image_25.png" /></figure>

<figure><img src="image_26.png" /></figure>

<figure><img src="image_27.png" /></figure>

<figure><img src="image_28.png" /></figure>

<figure><img src="image_29.png" /><figcaption>**Figure 17.** Spectra of the [C II]$_{158}$ µm line for different continuum subtractions. The line half-width masking region is indicated in the title of each panel. The colors and symbols are the same as in Figure 5.</figcaption></figure>
<figure><img src="image_30.png" /><figcaption>**Figure 18.** Aperture-integrated continuum flux density at 259 GHz as a function of the line masking half-width. The measurement has converged at the adopted masking width of 1.25 × FWHM([C II]).</figcaption></figure>

and find no evidence to reject the single-Gaussian model in favor of a double-Gaussian emission profile. Finally, Figure 18 shows the dust continuum measurement as a function of the masking half-width. The continuum is predictably higher for small masking regions, but has converged at the masking half-width adopted in this paper (1.25 × FWHM([C II]).

# Appendix F

## Previous PdBI Data and Merge with the New NOEMA Data

In this appendix, we analyze the previous PdBI observations of J1148+5251 published in Cicone et al. (2015). We do not recalibrate their data, and we follow our analysis method as presented in Sections 2 and 3. Based on the previous PdBI data, we show the total spectrum extracted in a $r = 3''$ aperture, applying residual scaling corrections (e.g., Jorsater & van Moorsel 1995; Walter & Brinks 1999; Walter et al. 2008; Novak et al. 2019), in Figure 19 (where we also show the new NOEMA spectrum for comparison). Clearly, the two data sets are compatible, and a single-Gaussian+continuum fit gives similar continuum levels and [C II]$_{158~\mu m}$ peak fluxes within the ~10%-20% amplitude calibration errors. We note that the noise in the previous PdBI spectrum is not uniform, as different data sets with different frequency coverages and antenna configurations were stitched together (Cicone et al. 2015). This leads to increased noise in the ranges $\nu \lesssim 254.7$ GHz and $\nu \gtrsim 258.0$ GHz. The best-fit continuum to the total spectrum in $r = 3''$ (4.6 $\pm$ 0.3 mJy) is perfectly consistent with the new measurement (4.5 $\pm$ 0.2 mJy) and the best-fit dust SED (Section 3.1),

and in tension with the 3.3 mJy continuum value (at 256 GHz) published in Cicone et al. (2015), derived from a best-fit model in the UV plane to line-free channels. We therefore infer that this tension is probably due to an unfortunate choice of continuum channels or an unstable UV plane-fitting routine in the old pipeline.

We then show in Figure 20 different spectra extracted from the previous PdBI data set with and without residual scaling. It can be clearly seen that the continuum flux at the edges of the band drops significantly in larger apertures when residual scaling is not used. This subtle effect is due to the fact that the continuum is only detected at the $\lesssim 2\sigma$ level in most line-free channels (e.g., at the edges of the band). When cleaning down to a typical $2\sigma$ threshold, there is in fact no or little flux above that threshold, and therefore no or few clean components are added to the final “clean” map (which is the sum of the clean Gaussian components and the residual map that does not have “clean beam” but “dirty beam” units) for channels that only contain continuum emission. As a consequence, these channels are dominated by the residual map that imprints the dirty beam pattern. With large apertures, e.g., $r > 3''$ in this case, one starts to integrate over negative sidelobes of the synthesized beam, decreasing the integrated flux significantly. This effect does not play a dominant role where the [C II]$_{158~\mu m}$ line is present, since in those channels the flux is dominated by the actual “clean” components. This can then lead, when combined with an undersubtracted continuum, to the illusion of a strong broad component in the final spectrum. However, with the residual scaling approach taken into account, this broad component is less significant (although still tentatively detected; see Figure 21), and the results become compatible with the new NOEMA data.

We now comment on the [C II]$_{158~\mu m}$ structure reported in Cicone et al. (2015). We show in Figure 22 the collapsed [C II]$_{158~\mu m}$ channels in the range of ($-1400, +1200$) km s$^{-1}$ based on the previous PdBI data only. We find a similar structure as previously reported, but at a different significance level (in our case, not exceeding the $3\sigma$ level in the extended regions). We therefore speculate that the rms might have been underestimated by the earlier GILDAS/MAPPING software release, which would explain the numerous $-6$ and $-3\sigma$ regions in the previously published [C II]$_{158~\mu m}$ map.

We conclude this appendix by presenting the results of our analysis of the merged PdBI and NOEMA [C II]$_{158~\mu m}$ spectra. We show in Figure 23 the total and continuum-subtracted spectra of the merged [C II]$_{158~\mu m}$ data. Although the statistical analysis of the spectrum extracted from this combined data set continues to prefer a single-component fit, there is evidence for the detection of a small flux excess in the redshifted [C II]$_{158~\mu m}$ line wing that might indicate the presence of weak outflow or an unresolved companion as in J0305–3150 (see Venemans et al. 2017a, 2019).
<figure><img src="image_31.png" /><figcaption>**Figure 19.** The total spectrum (black) of the [C II]158 µm line, extracted in *r* = 3" apertures, using residual scaling, in the previous PdBI data (left) and the new NOEMA data (right). The single-Gaussian (+continuum) fits to both data sets are shown in orange, and the parameters are given in the upper left corners. The best-fit continuum (constant) and single-Gaussian model parameters are consistent between the two data sets. We note that the increases in noise (shown as the red lines) in the PdBI data at low (<254.7 GHz) and high (>258 GHz) frequencies are due to the way in which the data were frequency-stitched at the time.</figcaption></figure>

<figure><img src="image_32.png" /></figure>
<figure><img src="image_33.png" /></figure>

<figure><img src="image_34.png" /></figure>

<figure><img src="image_35.png" /></figure>

<figure><img src="image_36.png" /><figcaption>**Figure 20.** The total spectra of the [C II]$_{158}$ µm line, with (first row) and without (second row) residual scaling, extracted in *r* = 2'' and *r* = 4'' apertures (the first and second columns, respectively), in the PdBI data. The single-Gaussian (+continuum) fits to both data sets are shown in orange. In the absence of residual scaling, the continuum flux drops away from the [C II]$_{158}$ µm line in the larger apertures (see the text for details). Even though no continuum was subtracted in the UV plane, the continuum flux approaches ~0 mJy at the band edges in the *r* = 4'' aperture (without residual scaling), creating the illusion of a broad emission line.</figcaption></figure>
<figure><img src="image_37.png" /><figcaption>**Figure 21.** The total spectrum of the [C II]158 µm line (black), presented as in Figure 5, but using the previous PdBI data, following the methodology used in this paper. Both the single-Gaussian and double-Gaussian fits are consistent with those performed on the NOEMA data only (see Figure 5). With ΔBIC12 = 3.0 < 10, a broad [C II]158 µm component is only tentatively detected.</figcaption></figure>

<figure><img src="image_38.png" /><figcaption>**Figure 22.** [C II]158 µm line map integrated over (−1400, +1200) km s⁻¹ using the previous PdBI data. The (dashed) black contours are at the (−2, 2, 4, 8)σ level, where σ is the rms computed directly from the map with σ-clipping. The extended [C II]158 µm structure reported in Cicone et al. (2015) is recovered, albeit at a lower significance level (see the main text for details).</figcaption></figure>
<figure><img src="image_39.png" /></figure>

<figure><img src="image_40.png" /><figcaption>**Figure 23.** The total spectrum of the [C II]158 µm line (black), presented as in Figure 5, but using the merged PdBI and NOEMA data set. Both the single-Gaussian and double-Gaussian fits are consistent with those performed on the NOEMA data only (see Figure 5). With $\Delta\text{BIC}_{12} = 4.7 < 10$, a broad [C II]158 µm component is only tentatively detected in the merged PdBI and NOEMA data set.</figcaption></figure>

## Appendix G
### Channel and Moment Maps of the [C II]158 µm Emission

In this appendix, we present additional visualizations of the [C II]158 µm emission velocity structure in J1148+5251. First, we present channel maps in Figure 24 with a channel width of 117 km s⁻¹. No significant emission ($>3\sigma$) is detected at large radii or at velocity offsets $>400$ km s⁻¹ from the peak of the

[C II]158 µm line (see Cicone et al. 2015). Second, we present the integrated flux, mean velocity, and velocity dispersion (so-called moment maps) in Figure 25. These maps are generated using *Qubefit* (Neeleman et al. 2020) using all [C II]158 µm voxels detected at S/N $>3$ and standard parameters. We find no kinematic evidence for a bulge-dominated dispersion or a rotating disk model.

<figure><img src="image_41.png" /><figcaption>**Figure 24.** Channel map of the [C II]158 µm line of J1148+5251, in channels of 117 km s⁻¹. The contours are logarithmic (−8, −4, −2, 2, 4, 8, 16, 32)σ (rms). The color scaling is log-linear, the threshold being at 3σ (rms).</figcaption></figure>
<figure><img src="image_42.png" /><figcaption>**Figure 25.** Moment maps of the [C II]158 µm emission. The integrated velocity and velocity dispersion are only shown in pixels at the 3σ level. The absence of any velocity structure is expected, given that the 3σ detected area covers only ~4 beams.</figcaption></figure>

## ORCID iDs

Romain A. Meyer https://orcid.org/0000-0001-5492-4522

Fabian Walter https://orcid.org/0000-0003-4793-7880

Claudia Cicone https://orcid.org/0000-0003-0522-6941

Pierre Cox https://orcid.org/0000-0003-2027-8221

Roberto Decarli https://orcid.org/0000-0002-2662-8803

Roberto Neri https://orcid.org/0000-0002-7176-4046

Mladen Novak https://orcid.org/0000-0001-8695-825X

Antonio Pensabene https://orcid.org/0000-0001-9815-4953

Dominik Riechers https://orcid.org/0000-0001-9585-1462

Axel Weiss https://orcid.org/0000-0003-4678-3939

## References

* Bañados, E., Decarli, R., Walter, F., et al. 2015, ApJL, 805, L8

* Bañados, E., Venemans, B. P., Mazzucchelli, C., et al. 2018, Natur, 553, 473

* Beelen, A., Cox, P., Benford, D. J., et al. 2006, ApJ, 642, 694

* Bertoldi, F., Carilli, C. L., Cox, P., et al. 2003a, A&A, 406, L55

* Bertoldi, F., Cox, P., Neri, R., et al. 2003b, A&A, 409, L47

* Bischetti, M., Maiolino, R., Carniani, S., et al. 2019, A&A, 630, A59

* Boogaard, L., Meyer, R. A., & Novak, M. 2021, Interferopy: Analyzing Data Cubes from Radio-to-submm Observations, Zenodo, doi:10.5281/ZENODO.5775603

* Bosman, S. E. I. 2020, All $z>5.7$ quasars currently known v1.15 Zenodo, doi:10.5281/zenodo.6039724

* Brauher, J. R., Dale, D. A., & Helou, G. 2008, ApJS, 178, 280

* Carilli, C., & Walter, F. 2013, ARA&A, 51, 105

* Cicone, C., Maiolino, R., Gallerani, S., et al. 2015, A&A, 574, 14

* Cicone, C., Maiolino, R., Sturm, E., et al. 2014, A&A, 562, A21

* Costa, T., Sijacki, D., & Haehnelt, M. G. 2014, MNRAS, 444, 2355

* Croxall, K. V., Smith, J. D., Pellegrini, E., et al. 2017, ApJ, 845, 96

* Da Cunha, E., Groves, B., Walter, F., et al. 2013, ApJ, 766, 13

* De Breuck, C., Weiß, A., Béthermin, M., et al. 2019, A&A, 631, A167

* De Looze, I., Cormier, D., Lebouteiller, V., et al. 2014, A&A, 568, A62

* De Rosa, G., Venemans, B. P., Decarli, R., et al. 2014, ApJ, 790, 145

* Decarli, R., Dotti, M., Bañados, E., et al. 2019a, ApJ, 880, 157

* Decarli, R., Mignoli, M., Gilli, R., et al. 2019b, A&A, 631, L10

* Decarli, R., Walter, F., Carilli, C., et al. 2014, ApJ, 782, 78

* Decarli, R., Walter, F., Venemans, B. P., et al. 2018, ApJ, 854, 97

* Di Matteo, T., Springel, V., & Iernquist, L. 2005, Natur, 433, 604

* Díaz-Santos, T., Armus, L., Charmandaris, V., et al. 2017, ApJ, 846, 32

* Dunne, L., Eales, S. A., & Edmunds, M. G. 2003, MNRAS, 341, 589

* Fan, X., Strauss, M. A., Schneider, D. P., et al. 2003, AJ, 125, 1649

* Ferkinhoff, C., Brisbin, D., Nikola, T., et al. 2011, ApJL, 740, 29

* Ferkinhoff, C., Hailey-Dunsheath, S., Nikola, T., et al. 2010, ApJL, 714, 147

* Ferland, G. J., Chatzikos, M., Guzmán, F., et al. 2017, The 2017 Release of Cloudy, Zenodo, doi:10.5281/zenodo.4110791

* Fernández-Ontiveros, J. A., Spinoglio, L., Pereira-Santaella, M., et al. 2016, ApJS, 226, 19

* Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306

* Fujimoto, S., Ouchi, M., Ferrara, A., et al. 2019, ApJ, 887, 107

* Fujimoto, S., Silverman, J. D., Bethermin, M., et al. 2020, ApJ, 900, 1

* Gallerani, S., Ferrara, A., Neri, R., & Maiolino, R. 2014, MNRAS, 445, 2848

* Gullberg, B., De Breuck, C., Vieira, J. D., et al. 2015, MNRAS, 449, 2883

* Harris, C. R., Jarrod Millman, K., van der Walt, S. J., et al. 2020, Natur, 585, 357

* Hashimoto, T., Inoue, A. K., Tamura, Y., et al. 2019, PASJ, 71, 71

* Helou, G., Soifer, B. T., & Rowan-Robinson, M. 1985, ApJL, 298, L7

* Herrera-Camus, R., Bolatto, A., Smith, J. D., et al. 2016, ApJ, 826, 175

* Herrera-Camus, R., Schreiber, N. F., Genzel, R., et al. 2021, A&A, 649, A31

* Herrera-Camus, R., Sturm, E., Graciá-Carpio, J., et al. 2018a, ApJ, 861, 94

* Herrera-Camus, R., Sturm, E., Graciá-Carpio, J., et al. 2018b, ApJ, 861, 95

* Högbom, J. A. 1974, A&AS, 15, 417

* Hunter, J. D. 2007, CSE, 9, 90

* Ishibashi, W., & Fabian, A. C. 2015, MNRAS, 451, 93

* Izumi, T., Matsuoka, Y., Fujimoto, S., et al. 2021a, ApJ, 914, 36

* Izumi, T., Onoue, M., Matsuoka, Y., et al. 2021b, ApJ, 908, 235

* Jorsater, S., & van Moorsel, G. A. 1995, AJ, 110, 2037

* Kass, R. E., & Raftery, A. E. 1995, J. Am. Stat. Assoc., 90, 773

* Katz, H., Kimm, T., Haehnelt, M. G., et al. 2019, MNRAS, 483, 1029

* Katz, H., Kimm, T., Sijacki, D., & Haehnelt, M. G. 2017, MNRAS, 468, 4831

* Kennicutt, R. C. 1998, ARA&A, 36, 189

* Kennicutt, R. C., & Evans, N. J. 2012, ARA&A, 50, 531

* King, A. R. 2010, MNRAS Lett., 408, L95

* Koudmani, S., Henden, N. A., & Sijacki, D. 2021, MNRAS, 503, 3568

* Lagache, G., Cousin, M., & Chatzikos, M. 2018, A&A, 609, 130

* Lee, M. M., Nagao, T., De Breuck, C., et al. 2021, ApJ, 913, 41

* Leipski, C., Meisenheimer, K., Klaas, U., et al. 2010, A&A, 518, L34

* Leipski, C., Meisenheimer, K., Walter, F., et al. 2014, ApJ, 785, 154

* Li, J., Wang, R., Cox, P., et al. 2020, ApJ, 900, 131

* Maiolino, R., Cox, P., Caselli, P., et al. 2005, A&A, 440, 51

* Maiolino, R., Gallerani, S., Neri, R., et al. 2012, MNRAS, 425, L66

* Mazzucchelli, C., Bañados, E., Venemans, B. P., et al. 2017, ApJ, 849, 91

* Neeleman, M., Bañados, E., Walter, F., et al. 2019, ApJ, 882, 10

* Neeleman, M., Novak, M., Venemans, B. P., et al. 2021, ApJ, 911, 141

* Neeleman, M., Prochaska, J. X., Kanekar, N., et al. 2020, qubefit: MCMC kinematic modeling, Astrophysics Source Code Library, ascl:2005.013

* Negri, A., & Volonteri, M. 2017, MNRAS, 467, 3475

* Neri, R., Downes, D., Cox, P., & Walter, F. 2014, A&A, 562, 35

* Novak, M., Bañados, E., Decarli, R., et al. 2019, ApJ, 881, 63

* Novak, M., Venemans, B. P., Walter, F., et al. 2020, ApJ, 904, 131

* Oberst, T. E., Parshley, S. C., Stacey, G. J., et al. 2006, ApJL, 652, L125

* Oppenheimer, B. D., Davies, J. J., Crain, R. A., et al. 2020, MNRAS, 491, 2939

* Pavesi, R., Riechers, D. A., Capak, P. L., et al. 2016, ApJ, 832, 151

* Pavesi, R., Riechers, D. A., Faisst, A. L., Stacey, G. J., & Capak, P. L. 2019, ApJ, 882, 168

* Pensabene, A., Decarli, R., Bañados, E., et al. 2021, A&A, 652, A66

* Richardson, M. L. A., Scannapieco, E., Devriendt, J., et al. 2016, ApJ, 825, 83

* Riechers, D. A., Walter, F., Bertoldi, F., et al. 2009, ApJ, 703, 1338

* Robson, I., Priddey, R. S., Isaak, K. G., & McMahon, R. G. 2004, MNRAS, 351, L29

* Savage, B. D., & Sembach, K. R. 1996, ARA&A, 34, 279

* Shao, Y., Wang, R., Carilli, C. L., et al. 2019, ApJ, 876, 99

* Shao, Y., Wang, R., Jones, G. C., et al. 2017, ApJ, 845, 138

* Silk, J., & Rees, M. J. 1998, A&A, 331, L1

* Solomon, P. M., Downes, D., Radford, S. J. E., & Barrett, J. W. 1997, ApJ, 478, 144

* Springel, V., Di Matteo, T., & Hernquist, L. 2005, MNRAS, 361, 776

* The Astropy Collaboration, Price-Whelan, A. A. M., Sipocz, B. M., et al. 2018, AJ, 156, 19
Vallini, L., Gallerani, S., Ferrara, A., Pallottini, A., & Yue, B. 2015, ApJ, 813, 36

Venemans, B., Neeleman, M., Walter, F., et al. 2019, ApJL, 874, L30

Venemans, B. P., Decarli, R., Walter, F., et al. 2018, ApJ, 866, 159

Venemans, B. P., Walter, F., Decarli, R., et al. 2017a, ApJL, 851, L8

Venemans, B. P., Walter, F., Decarli, R., et al. 2017b, ApJ, 837, 146

Venemans, B. P., Walter, F., Neeleman, M., et al. 2020, ApJ, 904, 130

Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, NatMe, 17, 261

Wakker, B. P., & Schwarz, U. J. 1988, A&A, 200, 312

Walter, F., Bertoidl, F., Carilli, C., et al. 2003, Natur, 424, 406

Walter, F., & Brinks, E. 1999, AJ, 118, 273

Walter, F., Brinks, E., de Blok, W. J. G., et al. 2008, AJ, 136, 2563

Walter, F., Carilli, C., Bertoldi, F., et al. 2004, ApJL, 615, L17

Walter, F., Riechers, D., Cox, P., et al. 2009a, Natur, 457, 699

Walter, F., Riechers, D., Novak, M., et al. 2018, ApJL, 869, L22

Walter, F., Wei, A., Riechers, D. A., et al. 2009b, ApJ, 691, L1

Wang, F., Yang, J., Fan, X., et al. 2021, ApJL, 907, L1

Wang, R., Shao, Y., Carilli, C. L., et al. 2019, ApJ, 887, 40

Wang, R., Wagg, J., Carilli, C. L., et al. 2013, ApJ, 773, 44

Weiß, A., Henkel, C., Downes, D., & Walter, F. 2003, A&A, 409, L41

Weiß, A., Downes, D., Henkel, C., & Walter, F. 2005, A&A, 429, L25

Willott, C. J., McLure, R. J., & Jarvis, M. J. 2003, ApJL, 587, L15

Yang, J., Venemans, B., Wang, F., et al. 2019, ApJ, 880, 153

Yang, J., Wang, F., Fan, X., et al. 2020, ApJL, 897, L14

Zhao, Y., Lu, N., Xu, C. K., et al. 2016, ApJ, 819, 69