# The hot corino-like chemistry of four FUor-like protostars

F. Cruz-Sáenz de Miera, A. Coutens, Á. Kóspál, P. Ábrahám, A. Dutrey, S. Guilloteau

► To cite this version:

F. Cruz-Sáenz de Miera, A. Coutens, Á. Kóspál, P. Ábrahám, A. Dutrey, et al.. The hot corino-like chemistry of four FUor-like protostars. Astronomy & Astrophysics - A&A, 2025, 696, pp.A18. 10.1051/0004-6361/202451487. [hal-05016887](https://hal.science/hal-05016887v1)

HAL Id: hal-05016887

https://hal.science/hal-05016887v1

Submitted on 2 Apr 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished 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
# The hot corino-like chemistry of four FUor-like protostars

F. Cruz-Sáenz de Miera<fnref n="1" />,<fnref n="2" />,*, A. Coutens<fnref n="1" />, Á. Kóspál<fnref n="2" />,<fnref n="3" />,<fnref n="4" />, P. Ábrahám<fnref n="2" />,<fnref n="3" />,<fnref n="5" />, A. Dutrey<fnref n="6" />, and S. Guilloteau<fnref n="6" />

Received 12 July 2024 / Accepted 26 February 2025

## ABSTRACT

**Context.** Compared to Class 0 protostars, the higher densities and lower temperatures of the disk midplanes of Class I young stellar objects (YSOs) limit the detectability of complex organic molecules (COMs). The elevated luminosities of eruptive YSOs increase disk temperatures causing frozen molecules to sublimate, thus making them easier to detect.

**Aims.** Our aim is to investigate the chemical composition of four FUor-like Class I YSOs (L1551 IRS 5, Haro 5a IRS, V346 Nor, and OO Ser) and to compare their abundances of COMs with other YSOs in the literature.

**Methods.** We searched for COMs line emission in ALMA Band 6 observations. We used the CASSIS software to determine their column densities (N) and excitation temperatures (T<sub>ex</sub>) assuming local thermodynamical equilibrium.

**Results.** We detected 249 transitions from 12 COMs. In L1551 IRS 5 we identified CH<sub>3</sub>OH, <sup>13</sup>CH<sub>3</sub>OH, CH<sub>3</sub><sup>18</sup>OH, CH<sub>2</sub>DOH, CH<sub>3</sub>CHO, CH<sub>3</sub>OCH<sub>3</sub>, CH<sub>3</sub>OCHO, CH<sub>3</sub>COCH<sub>3</sub>, C<sub>2</sub>H<sub>5</sub>OH, C<sub>2</sub>H<sub>4</sub>CN, <sup>13</sup>CH<sub>3</sub>CN, and CH<sub>3</sub>C<sup>15</sup>N. Haro 5a IRS and OO Ser have emission from CH<sub>3</sub>OH, CH<sub>3</sub>CHO, CH<sub>3</sub>OCH<sub>3</sub>, and CH<sub>3</sub>OCHO. CH<sub>3</sub>COCH<sub>3</sub> is also detected in OO Ser. In V346 Nor we found CH<sub>3</sub>OH, CH<sub>2</sub>DOH, CH<sub>3</sub>CHO, CH<sub>3</sub>OCH<sub>3</sub>, CH<sub>3</sub>OCHO, and C<sub>2</sub>H<sub>5</sub>CN. The emission of COMs is compact in all targets. The analysis indicates that their temperatures are above 100 K. The abundance ratios of COMs derived for these eruptive YSOs, and for other protostars in the literature, span several orders of magnitude without any clear differentiation between the eruptive and quiescent YSOs. The column density of the main isotopolog of CH<sub>3</sub>OH should not be used as a reference, as most of the lines are optically thick.

**Conclusions.** The hot and compact emission of COMs indicates that the four FUor-like targets are hot corino-like. Spectral studies of these objects can be useful for investigating the complex organic chemistry at later evolutionary stages than the usual Class 0 stage.

**Key words.** astrochemistry – protoplanetary disks – stars: formation – stars: protostars – ISM: molecules

## 1. Introduction

Recent results have proposed that planets start forming in the early stages of protostellar evolution when the young stellar objects (YSOs) are still deeply embedded in their envelopes (e.g., Drążkowska et al. 2023). The chemical composition of these planets is expected to be inherited from the protoplanetary disks in which they are forming. Part of this composition is dictated by the presence of complex organic molecules (COMs). These are molecules composed of at least six atoms, of which at least one is carbon; they are the precursors of more complex pre-biotic species (Herbst & van Dishoeck 2009). During the early stages of the star formation process, when the temperatures are more elevated due to the higher mass accretion rates. COMs are more easily detected as they are still in the gas phase. In the case of low-mass protostars, the compact regions in which COMs are abundant in the gas phase are called hot corinos (Ceccarelli 2004). Hot corinos have been detected in both Class 0 and Class I sources (e.g., Cazaux et al. 2003; Bottinelli et al. 2004; Bergner et al. 2019; Bianchi et al. 2019). However, the more evolved Class I protostars have lower COM detection rates (Artur de la Villarmois et al. 2019; Yang et al. 2022). This

lack of detections can be attributed to a combination of a lower temperature, a smaller hot region, and a lower envelope density (van't Hoff et al. 2020; Hsu et al. 2023).

FU Orionis-type stars (FUors) are young stellar objects that experience a sudden increase in their mass accretion rate of 2 to 5 orders of magnitude, reaching values as high as $10^{-4} M_{\odot} year^{-1}$. Fischer et al. (2022) presented the most recent review of eruptive young stars, including FUors. This surge of the mass accretion rate is reflected as an increase in the accretion luminosity such that the bolometric luminosity of the protostar is then dominated by it. During this state of elevated accretion luminosity, the temperatures within ~10 au from the protostar rise to >100 K (e.g., FU Ori and V883 Ori; Lykou et al. 2022; Bourdarot et al. 2023; Alarcón et al. 2024), so that the protostellar disks and the envelope undergo detectable changes in their physical, mineralogical and chemical properties (e.g., Molyarova et al. 2018; Tobin et al. 2023; Zwicky et al. 2024). In particular, the increase in temperature causes the desorption of molecules that were frozen on the dust grains in the system.

Theoretical works have studied the effects of these outbursts, including the estimation of the changes in the abundances and distribution of molecules throughout the FUor envelopes and disks. Taquet et al. (2016) studied the changes in the abundances of dimethyl ether (CH<sub>3</sub>OCH<sub>3</sub>) and methyl formate (CH<sub>3</sub>OCHO).

<footnote n="1">¹ Institut de Recherche en Astrophysique et Planétologie, Université de Toulouse, OMP, CNRS, CNES, 9 av. du Colonel Roche, 31028 Toulouse Cedex 4, France</footnote>

<footnote n="2">² Konkoly Observatory, HUN-REN Research Centre for Astronomy and Earth Sciences, MTA Centre of Excellence, Konkoly-Thege Miklós út 15–17, 1121 Budapest, Hungary</footnote>

<footnote n="3">³ Institute of Physics and Astronomy, ELTE Eötvös Loránd University, Pázmány Péter sétány 1/A, 1117 Budapest, Hungary</footnote>

<footnote n="4">⁴ Max Planck Institute for Astronomy, Königstuhl 17, 69117 Heidelberg, Germany</footnote>

<footnote n="5">⁵ Institute for Astronomy (IfA), University of Vienna, Türkenschanzstrasse 17, 1180 Vienna, Austria</footnote>

<footnote n="6">⁶ Laboratoire d'astrophysique de Bordeaux, Univ. Bordeaux, CNRS, B18N, allée Geoffroy Saint-Hilaire, 33615 Pessac, France</footnote>

<footnote>* Corresponding author; fernando.cruz-saenz@irap.omp.eu</footnote>
out to envelope scales, and found that during a FUor-type outburst the abundances of these molecules with respect to methanol ($CH_3OH$) should be a factor of a few higher than pre-or post-outburst. Molyarova et al. (2018) modeled the changes in the protostellar disk during and after a FUor-type outburst, and found that the abundances of several molecules (including some COMs) can increase by ~3–8 orders of magnitude during the outburst and remain elevated between 20 and 120 years after the protostar has returned to its quiescent state. Zwicky et al. (2024) carried out a similar study to estimate the impact of these powerful accretion events and found a subset of molecules (e.g., $N_2H^+$, $C^{18}O$, $H_2CO$ and HCN) that can be used as tracers of previous outbursts.

Interferometric observations have resulted in detections of multiple COMs in the disks of five eruptive young stars: SVS13-A (Bianchi et al. 2017a, 2019, 2022b), V883 Ori (van 't Hoff et al. 2018; Lee et al. 2019; Yamato et al. 2024), L1551 IRS 5 (Bianchi et al. 2020), HOPS 373W (Lee et al. 2023), and B335 (Lee et al. 2025). SVS13-A (also known as V512 Per) is a Class I YSO that underwent a $\Delta K \sim 1.5$ mag brightening in 1990 and has been fading since then (Eisloffel et al. 1991; Hodapp & Chini 2014). Based on the amplitude of its brightening and on its optical/near-infrared spectrum, it was first classified as an EX Lupi-type outbursting young star (EXor). However, the length of the accretion outburst (>30 years) is vastly longer than what is expected in EXors (≤1 year), and instead its duration is more comparable to that of FUor-type outbursts. The COMs detected around SVS13-A are $CH_3CHO$, $CH_3OCHO$, $CH_3OCH_3$, $C_2H_5OH$, $NH_2CHO$, $CH_3CN$, $CH_2DCN$, and $HCOCH_2OH$ (De Simone et al. 2017; Bianchi et al. 2019, 2022a,b). V883 Ori is a Class I young eruptive star located in Orion, and is classified as FUor-like. This classification means that the brightening due to the accretion outburst was not detected, but the protostar presents the near-infrared spectroscopic signatures typical of FUors (Connelley & Reipurth 2018). The molecules detected in V883 Ori are $CH_3OH$, $CH_2DOH$, $^{13}CH_3OH$, $CH_3CHO$, $^{13}CH_3CHO$, $CH_3OCHO$, $CH_3O^{13}CHO$, $^{13}CH_3OCHO$, $CH_2DOCHO$, $CH_3OCDO$, $CH_3OCH_3$, $CH_3COCH_3$, $CH_3CN$, $c-C_2H_4O$, $t-HCOOH$, $t-C_2H_3CHO$, $s-C_2H_5CHO$, and $CH_3SH$ (van 't Hoff et al. 2018; Lee et al. 2019; Yamato et al. 2024). HOPS 373W is a young protostar in Orion that experienced a brief outburst in 2019 and began another accretion event in 2020 (Yoon et al. 2022). Lee et al. (2023) observed HOPS 373W with Band 7 of the Atacama Large Millimeter Array (ALMA). They identified multiple COMs including $CH_3OH$ and several deuterated forms of it, $CH_3CHO$, $CH_3OCHO$, $C_2H_5OH$, $CH_3CN$, and $HCOCH_2OH$. B335 is a Class I protostar that experienced a burst between 2010/2012 and 2023, detected only at mid- and far-infrared wavelengths. Lee et al. (2025) examined four epochs of ALMA data to study the changes in the column densities of $^{13}CH_3OH$, $CH_2DOH$, $CH_3OCHO$, $C_2H_5OH$, $HC(O)NH_2$, and $DC(O)NH_2$. L1551 IRS 5 is one of the FUor-type eruptive young stars that we focus on during this paper and is described in detail below.

In this paper, we present a chemical study of four FUor-like eruptive young stars, L1551 IRS 5, Haro 5a IRS, OO Ser, and V346 Nor, using ALMA Band 6 observations. L1551 IRS 5 is a Class I protobinary system located in the Taurus star-forming region at a distance of 147 pc from the Sun with a bolometric luminosity ($L_{bol}$) of 40 $L_⊙$ (Liseau et al. 2005). The northern star is more massive and luminous than the southern one (Liseau et al. 2005). The separation between the two protostars is 50 au. Their circumstellar disks have radii <10 au;

both disks are more massive than those of typical young stellar objects (Lim et al. 2016; Cruz-Sáenz de Miera et al. 2019; Kóspál et al. 2021). L1551 IRS 5 also contains a ~140 au circumbinary ring which is less massive than the two protostellar disks and has a slightly different inclination angle (Lim et al. 2016; Cruz-Sáenz de Miera et al. 2019; Takakuwa et al. 2020). Using ALMA Band 6 observations, Bianchi et al. (2020) detected three COMs toward this target: methanol ($CH_3OH$, $^{13}CH_3OH$, and $CH_2DOH$), methyl formate ($HCOOCH_3$), and ethanol ($CH_3CH_2OH$). Their beam size was not sufficient to resolve the two components of the system. However, using the methanol emission, they found differences in the systemic velocities (~3 km s⁻¹) along the length of the system in the north-to-south direction. Andreu et al. (2023) observed HDO and $H_2^{18}O$ in L1551 IRS 5, and while their Northern Extended Millimetre Array (NOEMA) observations also did not resolve the protobinary, they showed a comparable difference in velocities. Takakuwa et al. (2020) resolved the protobinary using ALMA Band 7 observations of $C^{18}O$. They confirmed that the difference in velocities is due to a combination of the difference in radial velocities between the two protostars and the velocity gradient due to the disk rotations. Mercimek et al. (2022) observed L1551 IRS 5 using the Institut de radioastronomie millimétrique (IRAM) 30 m antenna, and covering two spectral windows: one between 214.5 and 222.2 GHz, and another between 230.2 and 238.0 GHz (i.e., ~1 mm in wavelength). They identified 33 different species in the spectrum of L1551 IRS 5, including six COMs: $CH_3OH$, $CH_2DOH$, $CH_3CN$, $CH_3CCH$, $CH_3CHO$, and $CH_3OCHO$. Marchand et al. (2024) also used the 30 m telescope to observe L1551 IRS 5 at 2 mm and 3 mm. The authors identified emission lines of 75 different species, including ten COMs: $CH_3OH$, $CH_2DOH$, $CH_3CHO$, $CH_3OCH_3$, $CH_3OCHO$, $CH_3CN$, $l-C_4H_2$, $CH_3CCH$, $CH_2DCCH$, and $HC_5N$. Although their angular resolution was not enough to resolve the protostellar disks from the envelope, they found that some species have two components at two different temperatures (10 K and >30 K) that trace the cold envelope and the warmer emission close to the protostars.

For the other three FUor-like sources, this is the first study delving into the chemistry of their disks. Haro 5a IRS is a Class I protostar at a distance of 391 pc in Orion with $L_{bol} = 50 L_⊙$ (Reipurth & Aspin 1997), and it has a weaker companion with a separation of 0.70" to the southwest (Tobin et al. 2020; Kóspál et al. 2021). Haro 5a IRS and its companion have dust disk radii of 39 au and 34 au (Kóspál et al. 2021). V346 Nor is a Class 0/I in the Norma star-forming region with a distance of 700 pc with $L_{bol} = 160 L_⊙$ (Kóspál et al. 2017), whose outburst began before 1978. In 2010 V346 Nor dimmed suddenly; however, after a few months in quiescence, it began to brighten again, indicating the onset of a new outburst (Kraus et al. 2016; Kóspál et al. 2020). OO Ser is a Class I protostar located in the Serpens NW star-forming region, with a distance of 438 pc. Its outburst began sometime between 1991 and 1994, and it started dimming after reaching its peak brightness in late 1995 (Kóspál et al. 2007; Hodapp et al. 2012). Kóspál et al. (2007) calculated the pre-outburst and outburst bolometric luminosities of OO Ser, which, after we adjusted them to the assumed distance of 438 pc, are 9 and 70 $L_⊙$, respectively. ALMA observations revealed that its dust disk has a radius of 26 au (Kóspál et al. 2021).

We present the ALMA observations and the data reduction in Sect. 2. The resulting spectra and moment maps are shown in Sect. 3. In Sects. 4 and 5, we present the analysis and the discussion. Finally, we present our conclusions in Sect. 6.
<figure><img src="image_2.png" /><figcaption>**Fig. 1.** Continuum maps of the four FUor sources at 1.33 mm. In all panels the cross markers indicate the continuum peaks, i.e., the positions of the protostars, and the square markers indicate the pixels used to extract the spectra. In the case of L1551 IRS 5, the diamond markers indicate the positions of the alternative pixels used to examine potential spatial changes of the abundances (see Sect. 5.5). Panel a shows the map of L1551 IRS 5 centered on the position between the two protostars. The contours are at 8, 18, 44, 105, and 250σ with σ = 0.37 mJy beam⁻¹. Panel b is the map of Haro 5a IRS with contours at 5, 12, 30, 75, and 186σ with σ = 0.22 mJy beam⁻¹. Panel c is V346 Nor with contours at 5, 10, 23, 49, and 106σ with σ = 0.05 mJy beam⁻¹. Panel d shows OO Ser with contours at 5, 12, 30, 75, and 186σ with σ = 0.22 mJy beam⁻¹.</figcaption></figure>

## 2. Observations

The four FUor-like sources were observed as part of the ALMA program 2016.1.00209.S (PI: Takami) in Cycle 4. The four targets were observed using ALMA's 7m array plus two extended configurations (C40-3 and C40-6) resulting in baseline ranges between 8 m and ~3650 m. The four targets had similar integration times. For L1551 IRS 5 these were 4.6, 3.6, and 11.3 minutes for the 7m, C40-3, and C40-6 arrays, respectively. In the case of Haro 5a IRS these were 4.1, 3.1, and 9.3 minutes, respectively, with the 7m and C40-6 observations being repeated due to weather conditions. For V346 Nor and OO Ser the integration times were 4.1, 3.1, and 9.8 minutes for the three arrays. In the case of OO Ser, the 7m array observations were also repeated. The dates of the observations and other details can be found in Table 1 of Kóspál et al. (2021). Data were taken in Band 6, with three spectral windows centered at the ¹²CO, ¹³CO, and C¹⁸O J = 2 − 1 lines, and two spectral windows of 1875 MHz bandwidth each, centered at 216.9 GHz and 232.2 GHz. These latter two spectral windows are the focus of this work. They were taken at a spectral resolution of ~0.490 MHz, i.e., 0.677 km s⁻¹ and 0.633 km s⁻¹.

The data sets were calibrated using Common Astronomy Software Applications (CASA) v.5.1.1 (McMullin et al. 2007). We carried out self-calibration for both the visibility phases and amplitudes to improve the signal-to-noise ratio (S/N), and we applied these solutions to the continuum-subtracted measurement sets. The common baselines between the three different observing configurations were consistent so we concatenated them into a single measurement set. This concatenated measurement set includes observations using the full baseline range and thus recovers the emission at most spatial scales.

The continuum-subtracted data cubes were created using tclean in CASA with a Briggs scheme and a robust value of 0.5, thus generating maps with a compromise between sensitivity and angular resolution. Both data cubes were made using the same image size, pixel size, and beam shape that was obtained using the spectral window with the shortest frequencies. The beam shapes for each target and the resulting rms for each spectral window can be found in Table 1.

<table><caption>**Table 1.** Beam shapes and position angles for each target, and the rms calculated for each spectral window and a spectral resolution of 0.490 MHz.</caption><thead><tr><th>Target</th><th>Beam size (′′)</th><th>Beam PA (°)</th><th>rms<sub>216</sub> (K)</th><th>rms<sub>232</sub> (K)</th></tr></thead><tbody><tr><td>L1551 IRS 5</td><td>0.204 × 0.189</td><td>-60.8</td><td>1.35</td><td>1.27</td></tr><tr><td>Haro 5a IRS</td><td>0.278 × 0.141</td><td>-78.2</td><td>1.35</td><td>1.27</td></tr><tr><td>V346 Nor</td><td>0.171 × 0.116</td><td>+87.3</td><td>1.34</td><td>1.30</td></tr><tr><td>OO Ser</td><td>0.221 × 0.152</td><td>-68.3</td><td>1.35</td><td>1.27</td></tr></tbody></table>

## 3. Results

In this section we present the extraction of the spectra used in our analysis, the identification of the emission lines, and discuss the morphology of the emission of the COMs and the smaller molecules.

### 3.1. Spectra

In their analysis of the continuum, Kóspál et al. (2021) found that the inner disks of our targets are optically thick. Thus, line emission is expected to be weaker close to the peak continuum position. We verified this by extracting spectra from the pixel that matches the position of the peak in the continuum maps, i.e., the position of the protostar, and at some random pixels separated from the center. We found weaker line emission in the spectra extracted on the pixel closest to the protostar. Thus, we decided to extract the spectra from off-center positions where we could find the strongest line emission. We chose these positions by extracting spectra using a pixel-size aperture placed at multiple positions throughout the disks, and finding the pixels where the line emission was stronger in both spectral windows. Our decision to use an individual pixel to extract the spectra, instead of a larger aperture to boost the S/N, was to minimize the blending caused by the distribution of velocities found in a larger aperture. The continuum maps showing the positions of the protostars and our apertures are shown in Fig. 1. The coordinates of all the protostars and the ones of the selected positions for the extraction of spectra can be found in Table 2. For all targets,
<table><caption>**Table 2.** Positions of the protostars and the pixels used to extract the spectra, and the separation between the two positions.</caption><thead><tr><th>Target</th><th>Protostar coordinates</th><th>Spectra coordinates</th><th>Separation (′′)</th><th>v<sub>sys</sub> (km s<sup>-1</sup>)</th></tr></thead><tbody><tr><td>L1551 IRS 5 N</td><td>04:31:34.1623 +18:08:04.717</td><td>04:31:34.1575 +18:08:04.793</td><td>0.102</td><td>7.56</td></tr><tr><td>L1551 IRS 5 S</td><td>04:31:34.1664 +18:08:04.356</td><td>04:31:34.1669 +18:08:04.253</td><td>0.103</td><td>5.16</td></tr><tr><td>Haro 5a IRS</td><td>05:35:26.5594 –05:03:55.130</td><td>05:35:26.5616 –05:03:55.064</td><td>0.074</td><td>14.06</td></tr><tr><td>V346 Nor</td><td>16:32:32.1971 –44:55:30.769</td><td>16:32:32.1951 –44:55:30.769</td><td>0.021</td><td>-4.58</td></tr><tr><td>OO Ser</td><td>18:29:49.1515 +01:16:19.676</td><td>18:29:49.1535 +01:16:19.646</td><td>0.042</td><td>8.22</td></tr></tbody></table>

**Notes.** The $v_{LSR}$ at the coordinates where we extracted the spectra is measured on the bright H$_2$S 2(2,0)–2(1,1) transition at 216.71 GHz detected in these observations.

the separation between the peak of the continuum and where we extracted the spectra is less than the beam size.

For the line identification we used the CASSIS<fnref n="1" /> software, which includes the Cologne Database for Molecular Spectroscopy<fnref n="2" /> (CDMS; Müller et al. 2001, 2005), the Jet Propulsion Laboratory (JPL; Pickett et al. 1998) molecular database<fnref n="3" />, and the Lille Spectroscopic Database<fnref n="4" /> (LSD). We used CDMS for most of the species. For CHD$_2$OH, CH$_3$COCH$_3$, and CH$_3$OCHO, we utilized the JPL database, and for CH$_3$CHO we used LSD. The spectroscopic references for the different molecules are: CH$_3$OH (Pickett et al. 1981; Sastry et al. 1984; Herbst et al. 1984; Anderson et al. 1990; Odashima et al. 1995; Müller et al. 2004; Xu et al. 2008), $^{13}$CH$_3$OH (Xu et al. 1996; Xu & Lovas 1997), CH$_3$$^{18}$OH (Fisher et al. 2007), CH$_2$DOH (Pearson et al. 2012), CH$_3$CHO (Smirnov et al. 2014), CH$_3$OCH$_3$ (Endres et al. 2009, and references therein), CH$_3$OCHO (Ilyushin et al. 2009, and references therein), CH$_3$COCH$_3$ (Oldag & Sutter 1992; Groner et al. 2002, and references therein), C$_2$H$_5$OH (Pearson et al. 1995, 1996, 2008; Müller et al. 2016), C$_2$H$_5$CN (Johnson et al. 1977; Boucher et al. 1980; Pearson et al. 1994; Fukuyama et al. 1996; Brauer et al. 2009, and references therein), and $^{13}$CH$_3$CN and CH$_3$$^{C15}$N (Müller et al. 2009, and references therein).

Among our targets, L1551 IRS 5 is the one with the most detected species, which we suspect is due to a combination of its proximity and its luminosity. In its northern disk we identified simple organics (H$_2$CO, D$_2$CO, DCN), sulfur-based species (H$_2$S, OC$^{33}$S, SO$_2$), and 12 COMs (including different isotopologs): CH$_3$OH, $^{13}$CH$_3$OH, CH$_3$$^{18}$OH, CH$_2$DOH, CH$_3$CHO, CH$_3$OCH$_3$, CH$_3$OCHO, CH$_3$COCH$_3$, C$_2$H$_5$OH, C$_2$H$_5$CN, $^{13}$CH$_3$CN, and CH$_3$$^{C15}$N. The L1551 IRS 5 N spectrum includes 249 transitions from these 12 COMs. The line emission from the southern disk of L1551 IRS 5 is weaker, and we only found emission from a subset of the molecules found in the northern disk (CH$_3$OH, CH$_3$OCH$_3$, CH$_3$OCHO, $^{13}$CH$_3$CN, and H$_2$S) for a total of 37 COM transitions. In the case of Haro 5a IRS the continuum absorption is stronger than in the other targets decreasing the amount of species detected: CH$_3$OH, CH$_3$CHO, CH$_3$OCH$_3$, CH$_3$OCHO, H$_2$S, and DCN. The spectra of Haro 5a IRS also include 46 transitions from the COMS, and show tentative emission from H$_2$CO and OC$^{33}$S but their low S/N prevented us from gaining any information about them. The continuum in V346 Nor and OO Ser did not show significant absorption. However, comparing these two targets

to L1551 IRS 5 N, we detected fewer COMs and fewer transitions in V346 Nor (6 COMs and 70 transitions) and OO Ser (6 COMs and 51 transitions). In V346 Nor we detected CH$_3$OH, CH$_2$DOH, CH$_3$CHO, CH$_3$OCH$_3$, CH$_3$OCHO, C$_2$H$_5$CN, H$_2$S, DCN, H$_2$CO, D$_2$CO, and SO$_2$. Finally, in OO Ser we detected CH$_3$OH, CH$_3$CHO, CH$_3$OCH$_3$, CH$_3$OCHO, CH$_3$COCH$_3$, H$_2$S, DCN, H$_2$CO, D$_2$CO, and SO$_2$.

## 3.2. Spatial distribution of COMs

In Fig. 2, we present the Moment 1 maps of CH$_3$OH for each of our targets, which were created with a $3\sigma$ cutoff, and in this subsection we discuss their morphologies.

Six transitions of CH$_3$OH are strongly detected in L1551 IRS 5 N, V346 Nor, and OO Ser, and they show a velocity distribution that suggests disk rotation. The southern disk of L1551 IRS 5 is quite faint and its systemic velocity is bluer than the limits of the color scale in the plot, thus this disk is not seen in the Moment maps. Finally, the Haro 5a IRS maps show few pixels with notable emission, which is due to the continuum absorption of its disk. In the case of L1551 IRS 5, the transitions with lower $E_{up}$ values are more extended toward the east, which indicates a difference in the physical conditions of the material of the circumbinary disk or the envelope. We suspect this extended emission traces material of the circumbinary disk or envelope that has been heated by the jets emanating from L1551 IRS 5 (Fridlund & Liseau 1998). However, the other three sources are also embedded and two of them (Haro 5a IRS and V346 Nor) drive jets or outflows (Kóspál et al. 2017; Cruz-Sáenz de Miera et al. 2023) yet we did not find a similar extended emission next to them. The three disks have comparable inclinations (Kóspál et al. 2021), which rules out a different perspective as an explanation for this lack of detection. Instead, it can be a matter of lack of sensitivity and angular resolution because L1551 IRS 5 is 3 and 5 times closer to us than V346 Nor and Haro 5a IRS, respectively.

We estimated the size of the CH$_3$OH emission from the disk by using the *imfit* tool in CASA to fit a 2D Gaussian to the Moment 0 maps of the 232.78 GHz transition of each target. In the case of L1551 IRS 5 N, the transitions with $E_{up} < 400$ K extend beyond the dust disk toward the west of the protostar, so we chose to only use the transition at 232.783 GHz ($E_{up} = 446.5$ K) to measure the extent of CH$_3$OH in all protostars. Part of the output of *imfit* is the semi-major full width at half maximum (FWHM), which we converted to a Gaussian $\sigma$ by dividing the FWHM by $2\sqrt{2\ln 2}$, and we used this $\sigma$ as a proxy for the CH$_3$OH disk radius.

For L1551 IRS 5 North, the disk has a deconvolved major-axis FWHM of $0.521 \pm 0.108''$, which gives a CH$_3$OH radius of

<footnote n="1">¹ CASSIS has been developed by IRAP-UPS/CNRS (Vastel et al. 2015). https://cassis.irap.omp.eu/</footnote>

<footnote n="2">² https://cdms.astro.uni-koeln.de/</footnote>

<footnote n="3">³ https://spec.jpl.nasa.gov/</footnote>

<footnote n="4">⁴ https://lsd.univ-lille.fr/</footnote>
<figure><img src="image_3.png" /><figcaption>**Fig. 2.** Moment 1 maps of the different CH₃OH transitions found in our targets (from top to bottom in increasing *E*<sub>up</sub>). From left to right are L1551 IRS 5, Haro 5a IRS, V346 Nor, and OO Ser. The color scale of each target has been centered on the systemic velocity of the protostar (see Table 2), and is limited to ±3 km s⁻¹ with respect to this velocity. The contours trace the continuum emission and are at the same levels as in Fig. 1. The star markers indicate the positions of the protostars. The beams are shown in the bottom left of each panel. The southern disk of L1551 IRS 5 is not seen in the maps as it is quite faint and its systemic velocity is bluer than the limits of the color scale of the plot.</figcaption></figure>

30 au, an inclination of 46±17°, and a position angle (PA) of 0.9±22.5°. This gas disk radius is a factor of 3 larger than the emission from the continuum at 1.3 mm (Cruz-Sáenz de Miera et al. 2019). The CH₃OH emission from southern disk was not resolved in most of the transitions, indicating that it has an extent of <12 au, i.e., comparable to the extent of the disk in the continuum. The continuum absorption of the Haro 5a IRS disk caused that its CH₃OH maps only show a few pixels with significant

emission, thus we were unable to compare its extent with respect to the dust. For V346 Nor the best-fit values are a deconvolved major-axis FWHM of 0.177±0.017'' (resulting in a disk radius of 53±5 au), an inclination of 24±17°, and a PA of 96±82°. Therefore, the geometrical properties of the CH₃OH disk in V346 Nor are comparable, within the uncertainties, to the results obtained for the continuum emission via Gaussian fitting and via radiative transfer by Kóspál et al. (2021). The disk of OO Ser was
barely resolved, and the fit resulted in a deconvolved major-axis FWHM of 0.155±0.068", an inclination of 20±145°, and a PA of 144±87°. The CH₃OH disk radius of OO Ser is 29±12 au which is the same as the dust radius found by Kóspál et al. (2021), and its inclination and PA are also in agreement with those from the dust disk.

In Fig. 3, we show the Moment 0 maps for different representative transitions of other COMs found in our targets. In L1551 IRS 5, these COMs trace the emission of the northern disk, while the southern disk is only barely detected. We measured the separations between the peaks in the Moment 0 maps of the COMs and the peaks of the continuum, and compared them to the semi-major axis of the beam of each map. These separation are ≤25% the size of the beam for both V346 Nor and OO Ser, ~30–40% for Haro 5a IRS, and ~45–80% for L1551 IRS 5 N. We consider these separations to be an indication of the COM emission being off-center only for L1551 IRS 5 N and Haro 5a IRS, and we corroborated this via a visual inspection. In the case of Haro 5a IRS, the emission peaks off-center, similar to CH₃OH, and its extent is much smaller than the dust disk. For V346 Nor the emission is found close to the peak of the continuum, tracing an area similar to CH₃OH. Finally, in the maps of OO Ser, most of the COMs have emission close to the peak of the continuum. Indeed, CH₃CHO, CH₃OCH₃ and CH₃OCHO have morphologies similar to CH₃OH but CH₃COCH₃ and C₂H₅CN peak off-center.

Figure 4 contains the velocity maps (moment 1) of the same representative transitions as Fig. 3, which we use to examine which molecules in which targets show evidence of rotation. The maps of L1551 IRS 5 show disk rotation for all the transitions in a southeast-northwest direction, comparable to the position angle of the disk found in the continuum (Cruz-Sáenz de Miera et al. 2019; Kóspál et al. 2021). The Haro 5a IRS disk only shows rotation in the maps of CH₃OCH₃ and CH₃OCHO, which has an orientation also in agreement with the position angle estimated from the continuum (Kóspál et al. 2021). For the other detected COM, CH₃CHO, only the redshifted side is detected. The remaining COMs were not detected and thus are not seen in the moment maps. The redshifted emission seen in the C₂H₅CN map is from a nearby CH₃CHO transition. In the case of V346 Nor, only CH₃OCH₃ shows the velocity gradient expected from rotation. For OO Ser, CH₃CHO, CH₃OCH₃ and CH₃OCHO indicate rotation.

### 3.3. Spatial distribution of smaller molecules

In Fig. 5, we show the Moment 1 maps for the small molecules that we examined. These maps were constructed with a $4\sigma$ cut-off. In this figure the color scales for Haro 5a IRS, V346 Nor, and OO Ser were centered around their systemic velocities (see Table 2). For L1551 IRS 5, the color scale was centered with respect to the average velocity of its two disks, which matches the systemic velocity obtained from APEX observations of ¹³CO J=3–2 (Cruz-Sáenz de Miera et al. 2023).

As was the case for the COMs, the emission is stronger and most extended in L1551 IRS 5. H₂S is the most extended and traces the circumbinary material. This molecule shows high velocity emission in the outer parts of the protostellar disks (i.e. the northern side of the northern disk and the southern side of the southern disk) that aligns to where Takakuwa et al. (2020) suggested the circumbinary spiral arms make contact with both protostars. The distributions of DCN, H₂CO, and D₂CO are similar to that of the lower $E_{up}$ transitions of CH₃OH, including the extended emission toward the western side of the northern

disk. However, DCN also shows emission in the outer sides of both disks, at the same positions as in H₂S. Meanwhile, the two remaining sulfur species, OC³S and SO₂, extend beyond the continuum emission.

Similar to the case of the COMs, the dust continuum in Haro 5a IRS absorbs the emission from most of these small molecules, and only H₂S and H₂CO show emission with hints of disk rotation aligned with the dust disk (Kóspál et al. 2021).

V346 Nor shows emission in all these small molecules except for OC³S. The small molecules extend beyond the COMs, and H₂S and DCN extend beyond the continuum emission, and H₂CO, D₂CO and SO₂ have emission with an extent comparable to the continuum. The distribution of velocities from these molecules is different from that of the COMs, which suggests contamination from the surrounded envelope.

Finally, the emission maps of OO Ser show detections of H₂S, DCN, D₂CO, and SO₂. These four species have barely extended emission and their velocity profiles are consistent with rotation similar to the COMs. H₂CO has emission in a few pixels close to the peak of the continuum. H₂S in OO Ser shows a blueshifted arm extending toward the south of the protostar; however, we lack the information to learn whether this is an embedded companion or signs of the surrounding envelope feeding material toward the FUor-like target similar to the recently reported streamer in FU Ori (Hales et al. 2024).

## 4. Analysis

After identifying the detected molecules, we determined the column densities ($N$) and excitation temperatures ($T_{ex}$) of the different species assuming local thermodynamic equilibrium (LTE), similar to what has been done for other hot corinos (e.g., Jørgensen et al. 2016; Coutens et al. 2016, 2018). Andreu et al. (2023) checked the LTE assumption for water isotopologs in L1551 IRS 5 and found consistent results using both LTE and non-LTE methods. Thus, we assume that the other three targets have densities comparable to L1551 IRS 5 so that the LTE approximation is also valid for them. We used the Markov chain Monte Carlo (MCMC) approach provided in CASSIS for our analysis. This consists in generating synthetic spectra then comparing them to the observed data, and finding the set of parameters which results in the minimum $\chi^2$. This MCMC approach enabled us to explore a wide parameter space. Each synthetic spectrum is created using a combination of column density ($N$), excitation temperature ($T_{ex}$), systemic velocity ($v_{LSR}$), FWHM of the line, and source size ($s$). A $\chi^2$ is computed for all the transitions selected for the comparison (meaning lines that are not blended and with sufficient S/N). The uncertainties we consider for this $\chi^2$ are the rms of each spectral window (see Table 1) and the nominal flux calibration uncertainty of ALMA (10%; ALMA Partnership 2016). Using this methodology instead of a rotational diagram approach allows us to verify that all the transitions in the full spectral window (not just those used in the fit) match the observational data. In addition, this approach of comparing synthetic spectra to observations does not assume that all the lines are optically thin, and instead the use of moderately optically thick lines can help constrain the source size of the analyzed species.

Table A.1 lists all the identified transitions We indicated the lines that are detected in each of the sources with diamond symbols. A filled diamond means that the transition is used in the fitting procedure, and the transitions with open markers are ignored due to either blending or very high optical depths.
<figure><img src="image_4.png" /><figcaption>**Fig. 3.** Integrated emission maps of a representative line for seven species on each of the FUor-type young stars. Each map has been centered around the peak of the continuum of the protostar. From top to bottom, the *E*<sub>up</sub> values of the transitions are 108.3, 80.9, 111.5, 115.5, 175.3, 157.7, and 142.4 K. The contours trace the continuum emission and are at the same levels as in Fig. 1.</figcaption></figure>
<figure><img src="image_5.png" /><figcaption>**Fig. 4.** Velocity field map (Moment 1) of the same representative transitions as in Fig. 3. The maps were generated with a $4\sigma$ clipping to minimize the noise in the resulting image. The velocities are clipped at $\pm 3 \text{ km s}^{-1}$. The velocities were centered on the $v_{\text{LSR}}$ of each target (see Table 2). The contours trace the continuum emission and are at the same levels as in Fig. 1.</figcaption></figure>

### 4.1. L1551 IRS 5 North

Our first MCMC runs for this disk consisted of letting the column density, excitation temperature, the systemic velocity, the FWHM of the line, and source size be free parameters for each species. The parameter space that we explored with MCMC for all targets was $10^{12} \text{ cm}^{-2} \le N \le 10^{19} \text{ cm}^{-2}$, $10 \text{ K} \le T_{\text{ex}} \le 300 \text{ K}$, $7.0 \text{ km s}^{-1} \le v_{\text{LSR}} \le 12.5 \text{ km s}^{-1}$, $1 \text{ km s}^{-1} \le \text{FWHM} \le 10 \text{ km s}^{-1}$,

and $0.01'' \le s \le 1.0''$. $N$, $T_{\text{ex}}$ and $s$ were explored in logarithmic scale, while linear scales were used for $v_{\text{LSR}}$ and FWHM. We ran the MCMC for 50 000 iterations and dropped the first 45 000 as a burn-in phase. Afterward we determined the best-fit parameters and their uncertainties using the 0.16, 0.50, and 0.84 percentiles. We found that shortening the burn-in phase did not cause changes in the values of the best-fit parameters.
<figure><img src="image_6.png" /><figcaption>**Fig. 5.** Moment 1 maps of the small molecules found in our targets. From left to right are L1551 IRS 5, Haro 5a IRS, OO Ser, and V346 Nor. The color scale of L1551 IRS 5 is centered around the average velocity of its two disks, while for the remainder of the targets the color scales have been centered on the systemic velocity of each protostar. All the color scales are limited to ±3 km s⁻¹. The contours trace the continuum emission and are at the same levels as in Fig. 1. The beams are shown in the bottom left corner of each panel.</figcaption></figure>

In our fitting procedure for CH₂DOH, we followed the recommendations provided by the JPL database and only included the lines with $K_a < 10$ in the upper transition. Furthermore, based on the analysis by Oyama et al. (2023), we only kept the CH₂DOH transitions that were present in the JPL database and in their list, and those with $J < 26$.

After obtaining a first version of the best-fit spectrum for all species, we incorporated it into a reference spectrum, which was used to double-check that the lines used in the fitting process were not blended. The full list of best-fit values and their statistical uncertainties for this first set of MCMC runs can be

found in Table 3. We emphasize that these uncertainties are purely statistical and should be considered as underestimated. The subtraction of the baseline is known to be challenging for sources with very rich line spectra. Among the sources studied here, L1551 IRS 5 North is the one with the highest number of lines. We consequently used that source to estimate how an error on the baseline subtraction could affect the column densities. We assumed that our continuum level was overestimated by 1 rms so we added 1 rms to the spectrum and performed a new set of fits. In this additional fitting process, we fixed $s$, $v_{LSR}$, and $T_{ex}$ to the values obtained when all parameters were initially free (i.e., top
<table><caption>**Table 3.** Best-fit parameters obtained for the L1551 IRS 5 northern disk using two different approaches.</caption><thead><tr><th>Species</th><th>n<sub>lines</sub></th><th>E<sub>up</sub> (K)</th><th>N (cm<sup>-2</sup>)</th><th>T<sub>ex</sub> (K)</th><th>FWHM (km s<sup>-1</sup>)</th><th>v<sub>LSR</sub> (km s<sup>-1</sup>)</th><th>s (")</th></tr></thead><tbody><tr><td colspan="8">COMs with free parameters</td></tr><tr><td>CH<sub>3</sub>OH</td><td>18</td><td>56–946</td><td>(1.72<sup>+0.08</sup><sub>-0.08</sub>) × 10<sup>19</sup></td><td>294.90<sup>+3.41</sup><sub>-4.40</sub></td><td>1.97<sup>+0.03</sup><sub>-0.03</sub></td><td>9.59<sup>+0.02</sup><sub>-0.01</sub></td><td>0.05<sup>+0.01</sup><sub>-0.01</sub></td></tr><tr><td>CH<sub>3</sub>OH<sup>a</sup></td><td>7</td><td>657–946</td><td>(1.11<sup>+0.17</sup><sub>-0.12</sub>) × 10<sup>20</sup></td><td>140.30<sup>+3.90</sup><sub>-4.60</sub></td><td>1.92<sup>+0.06</sup><sub>-0.04</sub></td><td>9.77<sup>+0.03</sup><sub>-0.03</sub></td><td>0.05<sup>+0.01</sup><sub>-0.01</sub></td></tr><tr><td><sup>13</sup>CH<sub>3</sub>OH</td><td>3</td><td>162–254</td><td>(1.22<sup>+0.25</sup><sub>-0.10</sub>) × 10<sup>17</sup></td><td>128.40<sup>+9.20</sup><sub>-13.00</sub></td><td>1.98<sup>+0.12</sup><sub>-0.09</sub></td><td>9.72<sup>+0.04</sup><sub>-0.04</sub></td><td>0.25<sup>+0.03</sup><sub>-0.06</sub></td></tr><tr><td>CH<sub>3</sub><sup>18</sup>OH</td><td>9</td><td>33–96</td><td>(1.12<sup>+0.16</sup><sub>-0.10</sub>) × 10<sup>16</sup></td><td>91.31<sup>+7.98</sup><sub>-8.76</sub></td><td>2.56<sup>+0.21</sup><sub>-0.16</sub></td><td>9.75<sup>+0.05</sup><sub>-0.05</sub></td><td>0.42<sup>+0.04</sup><sub>-0.04</sub></td></tr><tr><td>CH<sub>2</sub>DOH</td><td>9</td><td>49–352</td><td>(6.57<sup>+0.57</sup><sub>-0.54</sub>) × 10<sup>18</sup></td><td>145.60<sup>+4.90</sup><sub>-5.30</sub></td><td>1.94<sup>+0.13</sup><sub>-0.13</sub></td><td>10.06<sup>+0.06</sup><sub>-0.06</sub></td><td>0.05<sup>+0.01</sup><sub>-0.01</sub></td></tr><tr><td>CH<sub>2</sub>DOH<sup>b</sup></td><td>6</td><td>49–352</td><td>(5.25<sup>+0.90</sup><sub>-0.72</sub>) × 10<sup>18</sup></td><td>210.30<sup>+24.70</sup><sub>-12.50</sub></td><td>1.78<sup>+0.05</sup><sub>-0.04</sub></td><td>9.98<sup>+0.03</sup><sub>-0.02</sub></td><td>0.06<sup>+0.01</sup><sub>-0.01</sub></td></tr><tr><td>CH<sub>3</sub>CHO</td><td>27</td><td>6–387</td><td>(1.12<sup>+0.23</sup><sub>-0.15</sub>) × 10<sup>17</sup></td><td>212.30<sup>+8.40</sup><sub>-6.30</sub></td><td>2.28<sup>+0.06</sup><sub>-0.05</sub></td><td>9.80<sup>+0.02</sup><sub>-0.02</sub></td><td>0.07<sup>+0.01</sup><sub>-0.01</sub></td></tr><tr><td>CH<sub>3</sub>OCH<sub>3</sub></td><td>8</td><td>81–253</td><td>(9.57<sup>+0.65</sup><sub>-0.51</sub>) × 10<sup>16</sup></td><td>184.60<sup>+6.20</sup><sub>-5.90</sub></td><td>2.59<sup>+0.06</sup><sub>-0.06</sub></td><td>9.31<sup>+0.03</sup><sub>-0.03</sub></td><td>0.26<sup>+0.01</sup><sub>-0.01</sub></td></tr><tr><td>CH<sub>3</sub>OCHO</td><td>34</td><td>106–412</td><td>(3.58<sup>+0.30</sup><sub>-0.22</sub>) × 10<sup>17</sup></td><td>162.50<sup>+5.40</sup><sub>-8.48</sub></td><td>2.11<sup>+0.05</sup><sub>-0.04</sub></td><td>9.55<sup>+0.02</sup><sub>-0.02</sub></td><td>0.09<sup>+0.01</sup><sub>-0.01</sub></td></tr><tr><td>CH<sub>3</sub>COCH<sub>3</sub></td><td>13</td><td>111–116</td><td>(6.55<sup>+5.46</sup><sub>-3.35</sub>) × 10<sup>16</sup></td><td>104.60<sup>+10.10</sup><sub>-10.40</sub></td><td>2.38<sup>+0.26</sup><sub>-0.22</sub></td><td>9.60<sup>+0.07</sup><sub>-0.08</sub></td><td>0.05<sup>+0.02</sup><sub>-0.01</sub></td></tr><tr><td>C<sub>2</sub>H<sub>5</sub>OH</td><td>24</td><td>24–464</td><td>(5.24<sup>+1.52</sup><sub>-2.04</sub>) × 10<sup>17</sup></td><td>190.70<sup>+19.60</sup><sub>-23.30</sub></td><td>2.23<sup>+0.09</sup><sub>-0.10</sub></td><td>9.95<sup>+0.03</sup><sub>-0.03</sub></td><td>0.07<sup>+0.02</sup><sub>-0.01</sub></td></tr><tr><td>C<sub>2</sub>H<sub>5</sub>CN</td><td>20</td><td>147–368</td><td>(2.43<sup>+0.07</sup><sub>-0.07</sub>) × 10<sup>15</sup></td><td>190.60<sup>+5.30</sup><sub>-5.20</sub></td><td>2.48<sup>+0.07</sup><sub>-0.08</sub></td><td>9.74<sup>+0.03</sup><sub>-0.03</sub></td><td>0.42<sup>+0.02</sup><sub>-0.02</sub></td></tr><tr><td><sup>13</sup>CH<sub>3</sub>CN</td><td>5</td><td>78–257</td><td>(4.71<sup>+1.73</sup><sub>-1.30</sub>) × 10<sup>15</sup></td><td>170.60<sup>+15.50</sup><sub>-18.39</sub></td><td>2.11<sup>+0.11</sup><sub>-0.09</sub></td><td>9.62<sup>+0.03</sup><sub>-0.04</sub></td><td>0.06<sup>+0.01</sup><sub>-0.01</sub></td></tr><tr><td>CH<sub>3</sub>C<sup>15</sup>N</td><td>3</td><td>78–335</td><td>(4.64<sup>+0.95</sup><sub>-0.78</sub>) × 10<sup>14</sup></td><td>242.55<sup>+37.25</sup><sub>-35.53</sub></td><td>2.77<sup>+0.83</sup><sub>-0.38</sub></td><td>9.27<sup>+0.17</sup><sub>-0.18</sub></td><td>0.43<sup>+0.04</sup><sub>-0.12</sub></td></tr><tr><td colspan="8">COMs with T<sub>ex</sub> = 177.6 K, FWHM = 2.17 km s<sup>-1</sup>, v<sub>LSR</sub> = 9.67 km s<sup>-1</sup> and s = 0.08"</td></tr><tr><td>CH<sub>3</sub>OH</td><td>18</td><td>56–946</td><td>(4.66<sup>+0.16</sup><sub>-0.18</sub>) × 10<sup>18</sup></td><td></td><td></td><td></td><td></td></tr><tr><td>CH<sub>3</sub>OH<sup>a</sup></td><td>7</td><td>657–946</td><td>(1.60<sup>+0.06</sup><sub>-0.07</sub>) × 10<sup>19</sup></td><td></td><td></td><td></td><td></td></tr><tr><td><sup>13</sup>CH<sub>3</sub>OH</td><td>3</td><td>162–254</td><td>(4.97<sup>+0.30</sup><sub>-0.27</sub>) × 10<sup>17</sup></td><td></td><td></td><td></td><td></td></tr><tr><td>CH<sub>3</sub><sup>18</sup>OH</td><td>9</td><td>33–96</td><td>(9.94<sup>+0.49</sup><sub>-0.49</sub>) × 10<sup>16</sup></td><td></td><td></td><td></td><td></td></tr><tr><td>CH<sub>2</sub>DOH</td><td>9</td><td>49–352</td><td>(1.08<sup>+0.07</sup><sub>-0.06</sub>) × 10<sup>18</sup></td><td></td><td></td><td></td><td></td></tr><tr><td>CH<sub>2</sub>DOH<sup>b</sup></td><td>6</td><td>49–352</td><td>(2.81<sup>+0.19</sup><sub>-0.17</sub>) × 10<sup>18</sup></td><td></td><td></td><td></td><td></td></tr><tr><td>CH<sub>3</sub>CHO</td><td>27</td><td>65–387</td><td>(8.04<sup>+0.18</sup><sub>-0.17</sub>) × 10<sup>16</sup></td><td></td><td></td><td></td><td></td></tr><tr><td>CH<sub>3</sub>OCH<sub>3</sub></td><td>8</td><td>81–253</td><td>(3.41<sup>+0.11</sup><sub>-0.10</sub>) × 10<sup>17</sup></td><td></td><td></td><td></td><td></td></tr><tr><td>CH<sub>3</sub>OCHO</td><td>34</td><td>106–412</td><td>(5.05<sup>+0.11</sup><sub>-0.11</sub>) × 10<sup>17</sup></td><td></td><td></td><td></td><td></td></tr><tr><td>CH<sub>3</sub>COCH<sub>3</sub></td><td>13</td><td>111–116</td><td>(1.10<sup>+0.07</sup><sub>-0.07</sub>) × 10<sup>17</sup></td><td></td><td></td><td></td><td></td></tr><tr><td>C<sub>2</sub>H<sub>5</sub>OH</td><td>24</td><td>24–464</td><td>(3.39<sup>+0.12</sup><sub>-0.12</sub>) × 10<sup>17</sup></td><td></td><td></td><td></td><td></td></tr><tr><td>C<sub>2</sub>H<sub>5</sub>CN</td><td>20</td><td>147–368</td><td>(8.84<sup>+0.22</sup><sub>-0.23</sub>) × 10<sup>15</sup></td><td></td><td></td><td></td><td></td></tr><tr><td><sup>13</sup>CH<sub>3</sub>CN</td><td>5</td><td>78–257</td><td>(2.73<sup>+0.11</sup><sub>-0.12</sub>) × 10<sup>15</sup></td><td></td><td></td><td></td><td></td></tr><tr><td>CH<sub>3</sub>C<sup>15</sup>N</td><td>3</td><td>78–335</td><td>(1.15<sup>+0.16</sup><sub>-0.15</sub>) × 10<sup>15</sup></td><td></td><td></td><td></td><td></td></tr></tbody></table>

**Notes.** The top part of the table shows the results from the runs with all parameters left free, and the bottom part when fixing the T<sub>ex</sub>, FWHM, v<sub>LSR</sub>, and s to the median values obtained from the all-free runs. The uncertainties listed in each column only represent the statistical uncertainties from the fitting procedure. For CH<sub>3</sub>OH we include the results obtained using only the high excitation lines (flagged with <sup>a</sup>), and for CH<sub>2</sub>DOH we include the results obtained using only the lines with the lowest opacities (flagged with <sup>b</sup>). The second and third columns are the number of transitions used in the fit and their E<sub>up</sub> range.

panel of Table 3). We then determined the best-fit column densities and FWHM, and found that the values of N were higher by a factor between 1.2 to 1.8, depending on the COM. As the spectra of the other sources are less line crowded, this factor is expected to be an upper limit for the other sources.

The observed spectrum of CH<sub>3</sub>OH in L1551 IRS 5 North includes transitions with 56 K ≲ E<sub>up</sub> ≲ 946 K. The best-fit model is dominated by the lower-excitation lines with E<sub>up</sub> ≲ 335 K (see Fig. B.1) and the fluxes of the high excitation lines are underestimated. As the low-excitation lines are commonly highly optically thick in hot corinos (e.g., Jørgensen et al. 2016; Bianchi et al. 2020), we ran our fitting again using only lines

with E<sub>up</sub> ≳ 657 K. After comparing the results from both fits, we found that s, FWHM and the v<sub>LSR</sub> did not change significantly. However, the best-fit values of N and T<sub>ex</sub> are different. This high-excitation fitting results in a higher N and lower T<sub>ex</sub> than the one including all the CH<sub>3</sub>OH lines (see Table 3). Using highly optically thick lines could lead to an underestimation of the column density. So by using less optically thick high-excitation lines only, our best-fit N should be more accurate. To test this, we computed the column densities of CH<sub>3</sub>OH using its less abundant isotopologs, <sup>13</sup>CH<sub>3</sub>OH and CH<sub>3</sub><sup>18</sup>OH, and assuming standard isotopic ratios from the literature: [<sup>12</sup>C]/[<sup>13</sup>C] = 68 (Milam et al. 2005) and [<sup>16</sup>O]/[<sup>18</sup>O] = 560
<table><caption>**Table 4.** Similar to Table 3, but for the southern disk of L1551 IRS 5 assuming a source size of 0.08".</caption><thead><tr><th>Species</th><th>n<sub>lines</sub></th><th>E<sub>up</sub> (K)</th><th>N (cm<sup>-2</sup>)</th><th>T<sub>ex</sub> (K)</th><th>FWHM (km s<sup>-1</sup>)</th><th>v<sub>LSR</sub> (km s<sup>-1</sup>)</th></tr></thead><tbody><tr><td colspan="7">COMs with free parameters</td></tr><tr><td>CH<sub>3</sub>OH</td><td>8</td><td>56–746</td><td>(3.99<sup>+0.24</sup><sub>-0.24</sub>) × 10<sup>17</sup></td><td>193.90<sup>+10.90</sup><sub>-11.20</sub></td><td>2.23<sup>+0.12</sup><sub>-0.09</sub></td><td>2.81<sup>+0.05</sup><sub>-0.04</sub></td></tr><tr><td>CH<sub>3</sub>OCH<sub>3</sub></td><td>4</td><td>81</td><td>(2.92<sup>+0.35</sup><sub>-0.26</sub>) × 10<sup>16</sup></td><td>88.80<sup>+10.30</sup><sub>-7.14</sub></td><td>2.65<sup>+0.21</sup><sub>-0.21</sub></td><td>3.45<sup>+0.06</sup><sub>-0.07</sub></td></tr><tr><td>CH<sub>3</sub>OCHO</td><td>15</td><td>106–312</td><td>(5.00<sup>+0.57</sup><sub>-0.62</sub>) × 10<sup>16</sup></td><td>200.25<sup>+12.27</sup><sub>-12.15</sub></td><td>1.79<sup>+0.18</sup><sub>-0.18</sub></td><td>3.52<sup>+0.10</sup><sub>-0.11</sub></td></tr><tr><td><sup>13</sup>CH<sub>3</sub>CN</td><td>4</td><td>78–142</td><td>(4.30<sup>+0.69</sup><sub>-0.52</sub>) × 10<sup>14</sup></td><td>83.25<sup>+29.88</sup><sub>-12.65</sub></td><td>2.96<sup>+0.27</sup><sub>-0.31</sub></td><td>3.49<sup>+0.16</sup><sub>-0.17</sub></td></tr><tr><td colspan="7">COMs with T<sub>ex</sub> = 193.9 K, FWHM = 2.23 km s<sup>-1</sup> and v<sub>LSR</sub> = 3.47 km s<sup>-1</sup></td></tr><tr><td>CH<sub>3</sub>OCH<sub>3</sub></td><td>4</td><td>81</td><td>(6.86<sup>+0.39</sup><sub>-0.36</sub>) × 10<sup>16</sup></td><td></td><td></td><td></td></tr><tr><td>CH<sub>3</sub>OCHO</td><td>15</td><td>106–312</td><td>(5.58<sup>+0.21</sup><sub>-0.21</sub>) × 10<sup>16</sup></td><td></td><td></td><td></td></tr><tr><td><sup>13</sup>CH<sub>3</sub>CN</td><td>4</td><td>78–142</td><td>(5.87<sup>+0.61</sup><sub>-0.61</sub>) × 10<sup>14</sup></td><td></td><td></td><td></td></tr></tbody></table>

**Notes.** In the case of a single value of $E_{up}$, when $n_{lines} > 1$ it means that those lines are all within 1 K of each other.

(Wilson & Rood 1994). For this test, we fixed the $T_{ex}$, FWHM, $v_{LSR}$ and $s$ to the median values of the best fits found for all COMs (assuming all the CH$_3$OH isotopologs and the COMs come from the same region), which left the column density as the only free parameter. When determining the medians, we excluded the results of CH$_3$OH (and CH$_2$DOH) that were obtained using optically thick lines. Here we ran the MCMC for 10 000 iterations and dropped the first 5000 as a burn-in phase. We multiplied the best-fit $N$ by the aforementioned isotopic ratios and computed methanol column density values of $3.38 \times 10^{19}$ cm$^{-2}$ and $5.57 \times 10^{19}$ cm$^{-2}$ for the $^{13}$C and $^{18}$O isotopologs, respectively. These are in reasonable agreement with the fit of the high-excitation lines ($1.60 \times 10^{19}$ cm$^{-2}$) by factors of 2–3, and are an order of magnitude higher than the value obtained when using the optically thick lines ($4.66 \times 10^{18}$ cm$^{-2}$). Therefore, for the remainder of the paper, we only utilize the column density obtained using the high-excitation lines. The result of this high-excitation fit is marked with an $^a$ in Table 3.

A similar problem is seen for CH$_2$DOH. The high opacity lines are well fitted while the fluxes of the optically thin lines are underestimated. So we ran the MCMC again only including the transitions with the lowest opacities. The fit obtained when including the optically thick lines resulted in a higher $T_{ex}$ (210 K vs 146 K) and a lower $N$ ($5.25 \times 10^{18}$ cm$^{-2}$ vs $6.57 \times 10^{18}$ cm$^{-2}$). The results are included in Table 3, labeled with a $^b$. Just as for methanol, for the remainder of the paper we work with the fit found using only the low opacity lines.

For the other COMs, we also ran the second batch of MCMC runs where we fixed the $T_{ex}$, $v_{LSR}$, FWHM and $s$ parameters to the medians of the best-fit values previously determined assuming free parameters. The models are of comparable quality with respect to the observational data as the ones obtained using free parameters (see Appendix B). For the majority of the lines, the differences between the data (black line) and the total model (green line) are below the $3\sigma$ level (shaded area). The full list of best-fit values and their statistical uncertainties for both batches of MCMC runs can be found in Table 3.

## 4.2. L1551 IRS 5 South

The southern disk of L1551 IRS 5 has fewer and weaker lines than its northern counterpart. Due to a lack of constraints, we assumed for these MCMC runs the same source size as the northern disk (0.08"). The first set of MCMC runs were with

$N$, $T_{ex}$, $v_{LSR}$ and FWHM left as free parameters. The parameter space explored for these variables was the same as for the northern disk, except for $v_{LSR}$ which we limited to within 2 and 5 km s$^{-1}$.

We computed the medians of this first set of MCMC runs and found $T_{ex} = 193.90$ K, FWHM = $2.23$ km s$^{-1}$ and $v_{LSR} = 3.47$ km s$^{-1}$. For the computation of the medians, we excluded the results for CH$_3$OCH$_3$ because its best-fit $T_{ex}$ value cannot be constrained due to the four lines detected in this species having the same $E_{up}$ value. The median values of $T_{ex}$ and FWHM are the same as those found independently for CH$_3$OH. However, the best-fit value of $v_{LSR}$ found for CH$_3$OH is $0.66$ km s$^{-1}$ bluer than the median $v_{LSR}$, which is large enough to cause a bad fit of CH$_3$OH. It is unclear whether the discrepancy between $v_{LSR}$ values is real or an effect caused by absorption at redshifted velocities. We carried out a second set of MCMC runs for CH$_3$OCH$_3$, CH$_3$OCHO and $^{13}$CH$_3$CN by fixing $T_{ex}$, $v_{LSR}$ and FWHM to the medians found from the first run. This second set of column densities will be the ones used to compute the abundances with respect to methanol (see below) even though CH$_3$OH has a different $v_{LSR}$ than the other three COMs. These abundance ratios could then be quite uncertain.

The best-fit values and their statistical uncertainties for both MCMC runs, and the median values used for the second set of runs, are presented in Table 4. Similar to its northern counterpart, the models in L1551 IRS 5 S (Appendix C) match the data, with the uncertainties being less than the $3\sigma$ level.

## 4.3. Haro 5a IRS

In the case of Haro 5a IRS, we were not able to constrain the source size. There is a degeneracy between the column densities and the source size. Consequently, we decided to set the source size to 0.1", which is comparable to the size of the emission seen in the moment maps (see Figs. 2 and 3). We also had to fix the FWHM at $2.5$ km s$^{-1}$ in order to decrease the number of free parameters and better fit the low S/N lines. This FWHM value was measured with CASSIS by fitting a Gaussian profile on the lines of the different COMs separately. The parameter space of $N$ and $T_{ex}$ were the same as for the two previous disks, and for the $v_{LSR}$, it was between 13 and 15 km s$^{-1}$.

As for the two previous disks, a second set of MCMC runs were done using the median values of $T_{ex}$ and $v_{LSR}$. However, after a visual inspection of the quality of these follow-up fits, we
<table><caption>**Table 5.** Best-fit parameters obtained for Haro 5a IRS assuming a source size of 0.1".</caption><thead><tr><th>Species</th><th>n<sub>lines</sub></th><th>E<sub>up</sub> (K)</th><th>N (cm<sup>-2</sup>)</th><th>T<sub>ex</sub> (K)</th><th>FWHM (km s<sup>-1</sup>)</th><th>v<sub>LSR</sub> (km s<sup>-1</sup>)</th></tr></thead><tbody><tr><td colspan="7">COMs with free parameters and FWHM = 2.50 km s<sup>-1</sup></td></tr><tr><td>CH<sub>3</sub>OH</td><td>5</td><td>56–374</td><td>(9.88<sup>+1.19</sup><sub>-1.13</sub>) × 10<sup>16</sup></td><td>225.10<sup>+21.90</sup><sub>-23.77</sub></td><td></td><td>13.76<sup>+0.25</sup><sub>-0.23</sub></td></tr><tr><td>CH<sub>3</sub>CHO</td><td>13</td><td>65–129</td><td>(6.36<sup>+0.53</sup><sub>-0.41</sub>) × 10<sup>15</sup></td><td>110.60<sup>+6.64</sup><sub>-7.40</sub></td><td></td><td>13.42<sup>+0.07</sup><sub>-0.06</sub></td></tr><tr><td>CH<sub>3</sub>OCH<sub>3</sub></td><td>4</td><td>81</td><td>(2.30<sup>+0.48</sup><sub>-0.21</sub>) × 10<sup>16</sup></td><td>127.50<sup>+24.60</sup><sub>-7.70</sub></td><td></td><td>13.71<sup>+0.10</sup><sub>-0.04</sub></td></tr><tr><td>CH<sub>3</sub>OCHO</td><td>3</td><td>109–109</td><td>(1.67<sup>+0.60</sup><sub>-0.20</sub>) × 10<sup>16</sup></td><td>87.50<sup>+28.75</sup><sub>-18.70</sub></td><td></td><td>13.72<sup>+0.11</sup><sub>-0.12</sub></td></tr><tr><td colspan="7">COMs with T<sub>ex</sub> = 225.10 K, FWHM = 2.50 km s<sup>-1</sup> and v<sub>LSR</sub> = 13.76 km s<sup>-1</sup></td></tr><tr><td>CH<sub>3</sub>OH</td><td>5</td><td>56–374</td><td>(9.89<sup>+0.64</sup><sub>-0.70</sub>) × 10<sup>16</sup></td><td></td><td></td><td></td></tr><tr><td>CH<sub>3</sub>CHO</td><td>13</td><td>65–129</td><td>(1.34<sup>+0.07</sup><sub>-0.08</sub>) × 10<sup>16</sup></td><td></td><td></td><td></td></tr><tr><td>CH<sub>3</sub>OCH<sub>3</sub></td><td>4</td><td>81</td><td>(5.54<sup>+0.29</sup><sub>-0.30</sub>) × 10<sup>16</sup></td><td></td><td></td><td></td></tr><tr><td>CH<sub>3</sub>OCHO</td><td>3</td><td>109–109</td><td>(3.71<sup>+0.28</sup><sub>-0.27</sub>) × 10<sup>16</sup></td><td></td><td></td><td></td></tr></tbody></table>

**Notes.** In the case of a single value of $E_{up}$, when $n_{lines} > 1$ it means that those lines are all within 1 K of each other.

found that the models were not a good fit for all the molecules. We ran a third set of MCMC fits using the $T_{ex}$ and $v_{LSR}$ from the CH$_3$OH fit, and found that the resulting models were a slightly better fit than the ones using the medians. Thus, we kept the models obtained using the CH$_3$OH constraints. We examined the residuals between the data and the best-fit model and found that the residual of all lines used for the fit fall within $\pm 3\sigma$. The best-fit values and their statistical uncertainties are in Table 5, and the best-fit spectra are shown in Appendix D.

## 4.4. V346 Nor

For V346 Nor, the source size was also not constrained. So we assumed the size determined from the integrated emission map of CH$_3$OH, 0.177" (see Sect. 3.2). We first ran the MCMC with $N$, $T_{ex}$, $v_{LSR}$ and FWHM as free parameters. The parameter space for $v_{LSR}$ was from $-5.0$ km s$^{-1}$ to $-2.0$ km s$^{-1}$ while for the other parameters it was the same as for the other disks. The CH$_3$CHO lines appear to experience absorption at redshifted velocities, causing the best-fit $v_{LSR}$ value to be 0.3 km s$^{-1}$ lower compared to the other COMs. The best-fit values and their statistical uncertainties are in Table 6, and the best-fit spectra are shown in Appendix E.

After fitting the COMs with all free parameters, we ran a secondary fitting using the median values obtained from it for the second set of runs. Just as the case of Haro 5a IRS, the models obtained using the median values was not satisfactory, and we ran a third set of fits using the results of the CH$_3$OH fit. We carried out a visual comparison between the models obtained using the medians and the models obtained with the values found for CH$_3$OH. In the case of CH$_3$OH, we find that the model computed using the medians underestimates the fluxes of three lines by factors between 15% and 50%. For CH$_3$CHO, the medians model results in six transitions with fluxes $\le 25\%$ higher than the observed spectrum, whereas the lines from the CH$_3$OH model fit the data in shape and scale. Thus, the CH$_3$OH model is better at describing CH$_3$CHO. In CH$_3$OCH$_3$ and CH$_3$OCHO, the two models represent these molecules equally well. For $C_2H_5CN$, the CH$_3$OH model slightly underestimates the flux of one line by $\sim 25\%$, while the medians model overestimates one line by $\sim 10\%$. Due to the higher S/N of the CH$_3$OH and CH$_3$CHO lines, we decided to keep the column densities obtained using the CH$_3$OH parameters for the remaining of this paper.

We noted that 18 lines in this source have line profiles different from the synthetic ones (see Appendix E). Indeed, the majority of the CH$_3$OH lines in V346 Nor (Fig. E.1) have flatter line profiles. For the lines with strong emission (panels a, b, c, d, f, g, and j) it could indicate that they are slightly optically thick (which suggests that our column density is a lower limit, as discussed below). For the other weaker lines, they are too close to the noise level to accurately discuss their line profiles. Some CH$_3$CHO lines also differ from the synthetic spectrum, featuring a narrow blueshifted peak (see panels a, b, f, j, and m of Fig. E.4). However, the angular and spectral resolutions of the observations are not sufficient to identify what causes these different line profiles.

## 4.5. OO Ser

Similarly to V346 Nor and Haro 5a IRS, the source size was not constrained and we assumed the one obtained from the Moment 0 maps, 0.155". Most of the lines in OO Ser appear to be made of two components, with the secondary component being blueshifted and $\le 60\%$ weaker (depending on the line) than the main component. We attempted to add a second component to the fitting procedure but, because some of the intensities of this possible secondary component are close to the noise level, we could not find consistent results during the runs. To try to circumvent this issue, we focused on the strongest component. We also fixed the values of FWHM and $v_{LSR}$ to 2.0 km s$^{-1}$ and 9.25 km s$^{-1}$, respectively. Therefore, $N$ and $T_{ex}$ were the only free parameters, covering the same parameter space as the runs for the previous targets. In a second step, we fitted the remaining COMs with $N$ being the sole free parameter and $T_{ex}$ fixed to the median of the values from in the previous run. The best-fit values and their statistical uncertainties are in Table 7, and the best-fit spectra are shown in Appendix F. The dominant peaks in the line profiles of OO Ser appear to be well described with one component. With the exception of panel c in Fig. F.4, the residuals between the data and the best-fit models are all within the $3\sigma$ level.

## 4.6. Small molecules

For the small molecules, we only have one transition of each. To derive column densities, we consequently need to assume a
<table><caption>**Table 6.** Best-fit parameters obtained for V346 Nor assuming a source size of 0.177".</caption><thead><tr><th>Species</th><th>n<sub>lines</sub></th><th>E<sub>up</sub> (K)</th><th>N (cm<sup>-2</sup>)</th><th>T<sub>ex</sub> (K)</th><th>FWHM (km s<sup>-1</sup>)</th><th>v<sub>LSR</sub> (km s<sup>-1</sup>)</th></tr></thead><tbody><tr><td colspan="7">COMs with free parameters</td></tr><tr><td>CH<sub>3</sub>OH</td><td>13</td><td>56–835</td><td>(2.06<sup>+0.04</sup><sub>-0.04</sub>) × 10<sup>18</sup></td><td>210.70<sup>+3.90</sup><sub>-3.00</sub></td><td>5.43<sup>+0.10</sup><sub>-0.06</sub></td><td>-4.66<sup>+0.05</sup><sub>-0.01</sub></td></tr><tr><td>CH<sub>2</sub>DOH</td><td>3</td><td>113–167</td><td>(5.98<sup>+0.92</sup><sub>-1.34</sub>) × 10<sup>17</sup></td><td>245.40<sup>+36.37</sup><sub>-53.10</sub></td><td>3.91<sup>+0.16</sup><sub>-0.16</sub></td><td>-4.71<sup>+0.08</sup><sub>-0.09</sub></td></tr><tr><td>CH<sub>3</sub>CHO</td><td>17</td><td>65–357</td><td>(1.71<sup>+0.15</sup><sub>-0.14</sub>) × 10<sup>16</sup></td><td>123.30<sup>+13.30</sup><sub>-12.70</sub></td><td>5.04<sup>+0.24</sup><sub>-0.18</sub></td><td>-5.46<sup>+0.11</sup><sub>-0.10</sub></td></tr><tr><td>CH<sub>3</sub>OCH<sub>3</sub></td><td>8</td><td>81–253</td><td>(7.28<sup>+0.73</sup><sub>-0.67</sub>) × 10<sup>16</sup></td><td>152.00<sup>+9.30</sup><sub>-11.30</sub></td><td>5.37<sup>+0.11</sup><sub>-0.15</sub></td><td>-4.73<sup>+0.03</sup><sub>-0.04</sub></td></tr><tr><td>CH<sub>3</sub>OCHO</td><td>7</td><td>106–312</td><td>(6.12<sup>+0.49</sup><sub>-0.35</sub>) × 10<sup>16</sup></td><td>156.40<sup>+11.80</sup><sub>-8.10</sub></td><td>4.99<sup>+0.18</sup><sub>-0.33</sub></td><td>-4.69<sup>+0.08</sup><sub>-0.13</sub></td></tr><tr><td>C<sub>2</sub>H<sub>5</sub>CN</td><td>9</td><td>147–311</td><td>(2.98<sup>+0.18</sup><sub>-0.17</sub>) × 10<sup>15</sup></td><td>161.35<sup>+25.10</sup><sub>-14.65</sub></td><td>4.87<sup>+0.38</sup><sub>-0.20</sub></td><td>-3.73<sup>+0.13</sup><sub>-0.12</sub></td></tr><tr><td colspan="7">COMs with T<sub>ex</sub> = 210.7 K, FWHM = 5.43 km s<sup>-1</sup> and v<sub>LSR</sub> = -4.66 km s<sup>-1</sup></td></tr><tr><td>CH<sub>3</sub>OH</td><td>13</td><td>56–835</td><td>(2.04<sup>+0.03</sup><sub>-0.03</sub>) × 10<sup>18</sup></td><td></td><td></td><td></td></tr><tr><td>CH<sub>2</sub>DOH</td><td>3</td><td>113–167</td><td>(3.40<sup>+0.03</sup><sub>-0.03</sub>) × 10<sup>17</sup></td><td></td><td></td><td></td></tr><tr><td>CH<sub>3</sub>CHO</td><td>17</td><td>65–347</td><td>(2.17<sup>+0.05</sup><sub>-0.05</sub>) × 10<sup>16</sup></td><td></td><td></td><td></td></tr><tr><td>CH<sub>3</sub>OCH<sub>3</sub></td><td>8</td><td>81–253</td><td>(1.10<sup>+0.02</sup><sub>-0.02</sub>) × 10<sup>17</sup></td><td></td><td></td><td></td></tr><tr><td>CH<sub>3</sub>OCHO</td><td>7</td><td>106–312</td><td>(7.68<sup>+0.20</sup><sub>-0.16</sub>) × 10<sup>16</sup></td><td></td><td></td><td></td></tr><tr><td>C<sub>2</sub>H<sub>5</sub>CN</td><td>9</td><td>147–311</td><td>(2.85<sup>+0.07</sup><sub>-0.07</sub>) × 10<sup>15</sup></td><td></td><td></td><td></td></tr></tbody></table>

**Notes.** In the case of a single value of $E_{up}$, when $n_{lines} > 1$ it means that those lines are all within 1 K of each other.

<table><caption>**Table 7.** Best-fit parameters obtained for OO Ser assuming a source size of 0.155".</caption><thead><tr><th>Species</th><th>n<sub>lines</sub></th><th>E<sub>up</sub> (K)</th><th>N (cm<sup>-2</sup>)</th><th>T<sub>ex</sub> (K)</th><th>FWHM (km s<sup>-1</sup>)</th><th>v<sub>LSR</sub> (km s<sup>-1</sup>)</th></tr></thead><tbody><tr><td colspan="7">COMs with FWHM = 2.0 km s<sup>-1</sup> and v<sub>LSR</sub> = 9.25 km s<sup>-1</sup></td></tr><tr><td>CH<sub>3</sub>OH</td><td>8</td><td>56–746</td><td>(2.52<sup>+0.11</sup><sub>-0.11</sub>) × 10<sup>17</sup></td><td>158.70<sup>+9.40</sup><sub>-5.44</sub></td><td></td><td></td></tr><tr><td>CH<sub>3</sub>CHO</td><td>14</td><td>65–351</td><td>(6.39<sup>+0.64</sup><sub>-0.46</sub>) × 10<sup>15</sup></td><td>142.40<sup>+15.00</sup><sub>-11.80</sub></td><td></td><td></td></tr><tr><td>CH<sub>3</sub>OCH<sub>3</sub></td><td>4</td><td>81</td><td>(3.00<sup>+0.70</sup><sub>-0.68</sub>) × 10<sup>16</sup></td><td>174.35<sup>+22.31</sup><sub>-25.25</sub></td><td></td><td></td></tr><tr><td>CH<sub>3</sub>OCHO</td><td>7</td><td>106–312</td><td>(3.89<sup>+0.29</sup><sub>-0.18</sub>) × 10<sup>16</sup></td><td>150.90<sup>+11.80</sup><sub>-6.30</sub></td><td></td><td></td></tr><tr><td>CH<sub>3</sub>COCH<sub>3</sub></td><td>8</td><td>115</td><td>(1.95<sup>+0.73</sup><sub>-0.59</sub>) × 10<sup>16</sup></td><td>151.60<sup>+20.70</sup><sub>-19.57</sub></td><td></td><td></td></tr><tr><td colspan="7">COMs with T<sub>ex</sub> = 151.6 K, FWHM = 2.0 km s<sup>-1</sup> and v<sub>LSR</sub> = 9.25 km s<sup>-1</sup></td></tr><tr><td>CH<sub>3</sub>OH</td><td>8</td><td>56–746</td><td>(2.54<sup>+0.08</sup><sub>-0.08</sub>) × 10<sup>17</sup></td><td></td><td></td><td></td></tr><tr><td>CH<sub>3</sub>CHO</td><td>14</td><td>65–351</td><td>(6.58<sup>+0.28</sup><sub>-0.32</sub>) × 10<sup>15</sup></td><td></td><td></td><td></td></tr><tr><td>CH<sub>3</sub>OCH<sub>3</sub></td><td>4</td><td>81</td><td>(2.35<sup>+0.11</sup><sub>-0.10</sub>) × 10<sup>16</sup></td><td></td><td></td><td></td></tr><tr><td>CH<sub>3</sub>OCHO</td><td>7</td><td>106–312</td><td>(3.17<sup>+0.15</sup><sub>-0.15</sub>) × 10<sup>16</sup></td><td></td><td></td><td></td></tr><tr><td>CH<sub>3</sub>COCH<sub>3</sub></td><td>8</td><td>115</td><td>(1.97<sup>+0.12</sup><sub>-0.12</sub>) × 10<sup>16</sup></td><td></td><td></td><td></td></tr></tbody></table>

**Notes.** In the case of a single value of $E_{up}$, when $n_{lines} > 1$ it means that those lines are all within 1 K of each other.

source size and $T_{ex}$. Thus, we estimated their column densities following two approaches.

First, for each target, we fixed the source size and $T_{ex}$ to the values obtained for the COMs, and we let $N$, $v_{LSR}$ and FWHM be free parameters. We then ran an MCMC as for the previous fits. The best-fit values and uncertainties of all targets are reported in the second, third, and fourth columns of Table 8.

Second, we assumed a source size of 2", so that there is no beam dilution, and considered two $T_{ex}$ values (50 K and 300 K) to provide a range of possible column densities for these molecules. Here we fixed the $v_{LSR}$ and FWHM to the values

obtained in the first approach, and computed another set of MCMC runs. The best-fit values and uncertainties obtained for these two temperatures are shown in columns 5 and 6 of Table 8.

## 5. Discussion

In this section we discuss the results for each individual FUor-like target, and explore whether they have the qualities necessary to be defined as a hot corino: compact and warm regions that exhibit COMs. We complement our discussion by comparing our FUor-like targets with other YSOs.
<table><caption>**Table 8.** Best-fit parameters obtained for the small molecules.</caption><thead><tr><th rowspan="2">Species</th><th colspan="3">s and T<sub>ex</sub> from COMs</th><th colspan="2">s=2'' and T<sub>ex</sub>=50 K and 300 K</th></tr><tr><th>N (cm<sup>-2</sup>)</th><th>FWHM (km s<sup>-1</sup>)</th><th>v<sub>LSR</sub> (km s<sup>-1</sup>)</th><th>N<sub>50K</sub> (cm<sup>-2</sup>)</th><th>N<sub>300K</sub> (cm<sup>-2</sup>)</th></tr></thead><tbody><tr><td colspan="6" style="text-align: center; font-weight: bold;">L1551 IRS 5 North</td></tr><tr><td>H<sub>2</sub>S</td><td>(3.01<sup>+0.20</sup><sub>-0.18</sub>) × 10<sup>16</sup></td><td>2.61<sup>+0.09</sup><sub>-0.09</sub></td><td>9.38<sup>+0.05</sup><sub>-0.06</sub></td><td>(1.96<sup>+0.31</sup><sub>-0.22</sub>) × 10<sup>16</sup></td><td>(3.81<sup>+0.22</sup><sub>-0.24</sub>) × 10<sup>16</sup></td></tr><tr><td>DCN</td><td>(5.00<sup>+0.34</sup><sub>-0.35</sub>) × 10<sup>14</sup></td><td>2.61<sup>+0.15</sup><sub>-0.12</sub></td><td>9.68<sup>+0.06</sup><sub>-0.06</sub></td><td>(1.96<sup>+0.18</sup><sub>-0.15</sub>) × 10<sup>14</sup></td><td>(6.53<sup>+0.42</sup><sub>-0.35</sub>) × 10<sup>14</sup></td></tr><tr><td>H<sub>2</sub>CO</td><td>(8.83<sup>+0.62</sup><sub>-0.62</sub>) × 10<sup>16</sup></td><td>2.30<sup>+0.11</sup><sub>-0.10</sub></td><td>9.43<sup>+0.06</sup><sub>-0.06</sub></td><td>(1.38<sup>+0.12</sup><sub>-0.11</sub>) × 10<sup>17</sup></td><td>(8.54<sup>+0.61</sup><sub>-0.53</sub>) × 10<sup>16</sup></td></tr><tr><td>D<sub>2</sub>CO</td><td>(2.77<sup>+0.18</sup><sub>-0.20</sub>) × 10<sup>15</sup></td><td>2.90<sup>+0.14</sup><sub>-0.19</sub></td><td>9.35<sup>+0.07</sup><sub>-0.09</sub></td><td>(6.57<sup>+0.53</sup><sub>-0.52</sub>) × 10<sup>14</sup></td><td>(4.91<sup>+0.32</sup><sub>-0.33</sub>) × 10<sup>15</sup></td></tr><tr><td>OC<sup>33</sup>S</td><td>(2.65<sup>+0.40</sup><sub>-0.35</sub>) × 10<sup>15</sup></td><td>2.04<sup>+0.23</sup><sub>-0.20</sub></td><td>9.35<sup>+0.12</sup><sub>-0.13</sub></td><td>(2.31<sup>+0.35</sup><sub>-0.33</sub>) × 10<sup>15</sup></td><td>(2.52<sup>+0.40</sup><sub>-0.33</sub>) × 10<sup>15</sup></td></tr><tr><td>SO<sub>2</sub></td><td>(1.24<sup>+0.16</sup><sub>-0.14</sub>) × 10<sup>16</sup></td><td>2.42<sup>+0.30</sup><sub>-0.26</sub></td><td>9.10<sup>+0.13</sup><sub>-0.12</sub></td><td>(4.47<sup>+0.56</sup><sub>-0.51</sub>) × 10<sup>16</sup></td><td>(1.04<sup>+0.12</sup><sub>-0.11</sub>) × 10<sup>16</sup></td></tr><tr><td colspan="6" style="text-align: center; font-weight: bold;">L1551 IRS 5 South</td></tr><tr><td>H<sub>2</sub>S</td><td>(10.00<sup>+0.76</sup><sub>-0.78</sub>) × 10<sup>16</sup></td><td>3.08<sup>+0.17</sup><sub>-0.14</sub></td><td>4.12<sup>+0.07</sup><sub>-0.06</sub></td><td>(1.10<sup>+0.09</sup><sub>-0.09</sub>) × 10<sup>16</sup></td><td>(2.85<sup>+0.17</sup><sub>-0.17</sub>) × 10<sup>16</sup></td></tr><tr><td colspan="6" style="text-align: center; font-weight: bold;">Haro 5a IRS</td></tr><tr><td>H<sub>2</sub>S</td><td>(1.21<sup>+0.20</sup><sub>-0.14</sub>) × 10<sup>16</sup></td><td>4.09<sup>+0.73</sup><sub>-0.58</sub></td><td>13.63<sup>+0.21</sup><sub>-0.24</sub></td><td>(1.51<sup>+0.20</sup><sub>-0.18</sub>) × 10<sup>15</sup></td><td>(5.26<sup>+1.35</sup><sub>-0.81</sub>) × 10<sup>15</sup></td></tr><tr><td>DCN</td><td>(1.59<sup>+0.32</sup><sub>-0.28</sub>) × 10<sup>14</sup></td><td>2.96<sup>+1.23</sup><sub>-0.41</sub></td><td>14.00<sup>+0.27</sup><sub>-0.36</sub></td><td>(1.38<sup>+0.35</sup><sub>-0.27</sub>) × 10<sup>13</sup></td><td>(6.29<sup>+1.19</sup><sub>-1.00</sub>) × 10<sup>13</sup></td></tr><tr><td colspan="6" style="text-align: center; font-weight: bold;">V346 Nor</td></tr><tr><td>H<sub>2</sub>S</td><td>(7.89<sup>+0.38</sup><sub>-0.36</sub>) × 10<sup>16</sup></td><td>5.24<sup>+0.15</sup><sub>-0.14</sub></td><td>-4.33<sup>+0.09</sup><sub>-0.09</sub></td><td>(3.35<sup>+0.35</sup><sub>-0.28</sub>) × 10<sup>16</sup></td><td>(6.56<sup>+0.29</sup><sub>-0.30</sub>) × 10<sup>16</sup></td></tr><tr><td>DCN</td><td>(2.94<sup>+0.17</sup><sub>-0.17</sub>) × 10<sup>15</sup></td><td>6.47<sup>+0.59</sup><sub>-0.41</sub></td><td>-4.14<sup>+0.19</sup><sub>-0.19</sub></td><td>(1.41<sup>+0.48</sup><sub>-0.25</sub>) × 10<sup>15</sup></td><td>(2.57<sup>+0.17</sup><sub>-0.17</sub>) × 10<sup>15</sup></td></tr><tr><td>H<sub>2</sub>CO</td><td>(1.98<sup>+0.13</sup><sub>-0.12</sub>) × 10<sup>17</sup></td><td>6.57<sup>+0.41</sup><sub>-0.31</sub></td><td>-4.52<sup>+0.19</sup><sub>-0.16</sub></td><td>(2.23<sup>+0.17</sup><sub>-0.16</sub>) × 10<sup>17</sup></td><td>(1.53<sup>+0.10</sup><sub>-0.10</sub>) × 10<sup>17</sup></td></tr><tr><td>D<sub>2</sub>CO</td><td>(1.12<sup>+0.09</sup><sub>-0.09</sub>) × 10<sup>16</sup></td><td>5.48<sup>+0.68</sup><sub>-0.58</sub></td><td>-4.54<sup>+0.18</sup><sub>-0.21</sub></td><td>(1.57<sup>+0.12</sup><sub>-0.11</sub>) × 10<sup>15</sup></td><td>(1.25<sup>+0.09</sup><sub>-0.09</sub>) × 10<sup>16</sup></td></tr><tr><td>SO<sub>2</sub></td><td>(9.66<sup>+0.54</sup><sub>-0.56</sub>) × 10<sup>16</sup></td><td>6.95<sup>+0.36</sup><sub>-0.54</sub></td><td>-3.72<sup>+0.26</sup><sub>-0.15</sub></td><td>(3.53<sup>+0.24</sup><sub>-0.23</sub>) × 10<sup>17</sup></td><td>(6.70<sup>+0.35</sup><sub>-0.36</sub>) × 10<sup>16</sup></td></tr><tr><td colspan="6" style="text-align: center; font-weight: bold;">OO Ser</td></tr><tr><td>H<sub>2</sub>S</td><td>(1.14<sup>+0.12</sup><sub>-0.11</sub>) × 10<sup>16</sup></td><td>2.79<sup>+0.32</sup><sub>-0.23</sub></td><td>8.32<sup>+0.08</sup><sub>-0.12</sub></td><td>(4.03<sup>+0.49</sup><sub>-0.43</sub>) × 10<sup>15</sup></td><td>(1.25<sup>+0.14</sup><sub>-0.13</sub>) × 10<sup>16</sup></td></tr><tr><td>DCN</td><td>(3.51<sup>+0.27</sup><sub>-0.26</sub>) × 10<sup>14</sup></td><td>2.52<sup>+0.17</sup><sub>-0.14</sub></td><td>9.51<sup>+0.06</sup><sub>-0.07</sub></td><td>(9.41<sup>+0.81</sup><sub>-0.82</sub>) × 10<sup>13</sup></td><td>(3.82<sup>+0.32</sup><sub>-0.31</sub>) × 10<sup>14</sup></td></tr><tr><td>H<sub>2</sub>CO</td><td>(1.34<sup>+0.09</sup><sub>-0.09</sub>) × 10<sup>16</sup></td><td>2.76<sup>+0.15</sup><sub>-0.21</sub></td><td>8.88<sup>+0.09</sup><sub>-0.09</sub></td><td>(1.40<sup>+0.10</sup><sub>-0.09</sub>) × 10<sup>16</sup></td><td>(1.11<sup>+0.08</sup><sub>-0.08</sub>) × 10<sup>16</sup></td></tr><tr><td>D<sub>2</sub>CO</td><td>(1.30<sup>+0.16</sup><sub>-0.17</sub>) × 10<sup>15</sup></td><td>2.72<sup>+0.73</sup><sub>-0.37</sub></td><td>9.72<sup>+0.10</sup><sub>-0.11</sub></td><td>(2.16<sup>+0.29</sup><sub>-0.27</sub>) × 10<sup>14</sup></td><td>(1.82<sup>+0.22</sup><sub>-0.22</sub>) × 10<sup>15</sup></td></tr><tr><td>SO<sub>2</sub></td><td>(1.59<sup>+0.16</sup><sub>-0.17</sub>) × 10<sup>16</sup></td><td>2.57<sup>+0.39</sup><sub>-0.31</sub></td><td>9.43<sup>+0.15</sup><sub>-0.11</sub></td><td>(4.62<sup>+0.59</sup><sub>-0.54</sub>) × 10<sup>16</sup></td><td>(1.06<sup>+0.13</sup><sub>-0.12</sub>) × 10<sup>16</sup></td></tr></tbody></table>

**Notes.** The second, third, and fourth columns were obtained using the same source size and $T_{ex}$ as for the COMs. The fifth and sixth columns show the best-fit column densities assuming a source size of 2'', the FWHM and $v_{LSR}$ of the third and fourth columns, and $T_{ex}$ values of 50 K (fifth column) and 300 K (sixth column).

<table><caption>**Table 9.** Abundances with respect to CH<sub>3</sub>OH (N/N<sub>CH<sub>3</sub>OH</sub>). The lines of CH<sub>3</sub>OH are probably optically thick, and thus the listed abundances with respect to it are considered as upper limits.</caption><thead><tr><th>Species</th><th colspan="2">L1551 IRS 5 N</th><th>L1551 IRS 5 S</th><th>Haro 5a IRS</th><th>V346 Nor</th><th>OO Ser</th></tr><tr><th>Isotopolog</th><th>CH<sub>3</sub>OH<sup>(a)</sup></th><th><sup>13</sup>CH<sub>3</sub>OH<sup>(b)</sup></th><th>CH<sub>3</sub><sup>18</sup>OH<sup>(c)</sup></th><th>CH<sub>3</sub>OH</th><th>CH<sub>3</sub>OH</th><th>CH<sub>3</sub>OH</th></tr></thead><tbody><tr><td>CH<sub>3</sub>CHO</td><td>&lt;5.02 × 10<sup>-3</sup></td><td>2.38 × 10<sup>-3</sup></td><td>1.44 × 10<sup>-3</sup></td><td>—</td><td>&lt;1.35 × 10<sup>-1</sup></td><td>&lt;2.59 × 10<sup>-2</sup></td></tr><tr><td>CH<sub>3</sub>OCH<sub>3</sub></td><td>&lt;2.13 × 10<sup>-2</sup></td><td>1.01 × 10<sup>-2</sup></td><td>6.13 × 10<sup>-3</sup></td><td>&lt;1.71 × 10<sup>-1</sup></td><td>&lt;5.60 × 10<sup>-1</sup></td><td>&lt;9.25 × 10<sup>-2</sup></td></tr><tr><td>CH<sub>3</sub>OCHO</td><td>&lt;3.16 × 10<sup>-2</sup></td><td>1.49 × 10<sup>-2</sup></td><td>9.07 × 10<sup>-3</sup></td><td>&lt;1.40 × 10<sup>-1</sup></td><td>&lt;3.75 × 10<sup>-1</sup></td><td>&lt;1.25 × 10<sup>-1</sup></td></tr><tr><td>CH<sub>3</sub>COCH<sub>3</sub></td><td>&lt;6.88 × 10<sup>-3</sup></td><td>3.25 × 10<sup>-3</sup></td><td>1.98 × 10<sup>-3</sup></td><td>—</td><td>—</td><td>&lt;7.76 × 10<sup>-2</sup></td></tr><tr><td>C<sub>2</sub>H<sub>5</sub>OH</td><td>&lt;2.12 × 10<sup>-2</sup></td><td>1.00 × 10<sup>-2</sup></td><td>6.09 × 10<sup>-3</sup></td><td>—</td><td>—</td><td>—</td></tr><tr><td>C<sub>2</sub>H<sub>5</sub>CN</td><td>&lt;5.53 × 10<sup>-4</sup></td><td>2.62 × 10<sup>-4</sup></td><td>1.59 × 10<sup>-4</sup></td><td>—</td><td>&lt;1.40 × 10<sup>-3</sup></td><td>—</td></tr><tr><td>CH<sub>3</sub>CN<sup>(d)</sup></td><td>&lt;1.16 × 10<sup>-2</sup></td><td>5.49 × 10<sup>-3</sup></td><td>3.34 × 10<sup>-3</sup></td><td>&lt;1.00 × 10<sup>-1</sup></td><td>—</td><td>—</td></tr></tbody></table>

**Notes.** <sup>(a)</sup>Using the column density obtained from the high-excitation fit. <sup>(b)</sup>Assuming an isotopic ratio of [¹²C]/[¹³C] = 68. <sup>(c)</sup>Assuming an isotopic ratio of [¹⁶O]/[¹⁸O] = 560. <sup>(d)</sup>The column density of CH<sub>3</sub>CN was calculated from our best-fit N of ¹³CH<sub>3</sub>CN and assuming an isotopic ratio of [¹²C]/[¹³C] = 68.
## 5.1. L1551 IRS 5

We found that both disks of L1551 IRS 5 have emission lines from multiple COMs in the inner areas of their disks, and their excitation temperatures are above 100 K. Thus, both disks have the requirements to classify as hot corinos. Bianchi et al. (2020) had already classified the northern disk of L1551 IRS 5 as a hot corino and suggested that the southern component might also be one. Here we confirm that both disks have hot corino-like properties.

Bianchi et al. (2020) ran a non-LTE fitting procedure for L1551 IRS 5 and found a lower limit of $1 \times 10^{19} \text{ cm}^{-2}$ for the *N* of CH$_3$OH with a $T_{\text{ex}}$ of 100 K, a FWHM of 3.5 km s$^{-1}$, and a source size of 0.15". While our column densities for the northern disk are in the same order of magnitude, a more quantitative comparison between our results and theirs is not trivial due to the different excitation temperatures, source sizes and FWHMs assumed. In the case of CH$_2$DOH, Bianchi et al. (2020) fixed the source size to be the same size as for the main isotopolog of methanol, and their best-fit values of $T_{\text{ex}}$ and *N* are 88 K and $\ge 5 \times 10^{17} \text{ cm}^{-2}$, respectively. Their $T_{\text{ex}}$ is lower than the value we obtained (210.3 K) and to the median of all the COMs (177.6 K). Both of our best-fit values of *N* ($5.25 \times 10^{18}$ and $2.81 \times 10^{18} \text{ cm}^{-2}$) agree with their lower limit. For C$_2$H$_5$OH, Bianchi et al. (2020) fixed the temperature and source size to their values obtained from methanol and obtained a best-fit *N* that is approximately 2–3 times lower than our values, which is in agreement with our results using a smaller source size (0.08"). Lastly, their estimate of the column density of CH$_3$OCHO is comparable to our best-fit values.

We used the best-fit column densities obtained using the median $T_{\text{ex}}$ of 177.6 K to calculate the deuteration ratio, D/H, of methanol for the northern disk of L1551 IRS 5. The *N* ratio between CH$_2$DOH ($2.81 \times 10^{18} \text{ cm}^{-2}$) and CH$_3$OH ($1.60 \times 10^{19} \text{ cm}^{-2}$) is $1.8 \times 10^{-1}$, which gives a statistical D/H ratio of CH$_3$OH of 5.9%. The column densities of CH$_3$OH from the $^{13}$C and $^{18}$O isotopologs ($3.38 \times 10^{19} \text{ cm}^{-2}$ and $5.57 \times 10^{19} \text{ cm}^{-2}$, respectively) lead to isotopic ratios of $8.2 \times 10^{-2}$ and $5.1 \times 10^{-2}$, and statistical D/H ratios of 2.8% and 1.7%. These isotopic ratios are one or two orders of magnitude higher than what was found for HDO/H$_2$O ($1.6 - 2.1 \times 10^{-3}$; Andreu et al. 2023). However, finding a higher D/H ratio for COMs than for water is common in low-mass protostars (e.g., Persson et al. 2013; Jørgensen et al. 2016, for IRAS 16293-2422), which is explained by the earlier formation of water in the star formation process (e.g., Ceccarelli et al. 2014; Furuya et al. 2016; Jensen et al. 2021). Our calculated D/H fractionation is in the same order to what was found for CH$_2$DOH in IRAS 16293 A (2.8%; Manigand et al. 2020), IRAS 16293 B (2.4%; Jørgensen et al. 2018), IRAS 2A and IRAS 4A (5.8% and 3.7%, respectively; Taquet et al. 2019), and Ser-emb 11 (4.0%; Martín-Doménech et al. 2021). At the same time, it is higher than what was found for SVS13-A (0.05–0.24%; Bianchi et al. 2017a), HH 212 (0.8%; Bianchi et al. 2017b), L483 (0.8%; Jacobsen et al. 2019) and for V883 Ori (0.2–0.7%; Lee et al. 2019; Yamato et al. 2024).

## 5.2. Haro 5a IRS

The line emission from Haro 5a IRS appears to be compact (see Figs. 2 and 3). Taking this with the presence of four COMs and $T_{\text{ex}}$ above 100 K, it suggests that Haro 5a IRS is also similar to a hot corino. Haro 5a IRS has a companion 0.7" to the southwest, and we did not find line emission from it.

Compared to the CH$_3$OH maps of the other FUor-like targets (Fig. 2), Haro 5a IRS is the only one that does not show a disk-like velocity distribution of the gas. Indeed, only the transitions with $E_{\text{up}}$ values from 165.3 K to 190.4 K show emission in more than a few pixels compared to the other three shown transitions. We think that most of the line emission is hidden by the optically thick dust disk.

## 5.3. V346 Nor

Similar to the other targets, V346 Nor can be considered as a hot corino-like given its compact emission of multiple COMs with $T_{\text{ex}}$ higher than 100 K. We used the best-fit column densities obtained using the median $T_{\text{ex}}$ of 210.7 K to calculate the deuteration ratio, D/H, of methanol. The *N* ratio between CH$_2$DOH ($3.40 \times 10^{17} \text{ cm}^{-2}$) and CH$_3$OH ($2.04 \times 10^{18} \text{ cm}^{-2}$) is $1.7 \times 10^{-1}$, which gives a statistical D/H ratio of CH$_3$OH of 5.6%. These numbers are the same as those found for L1551 IRS 5 North when using the main methanol isotopolog. However, the possible underestimation of the column density of CH$_3$OH and the low quality of the CH$_2$DOH model make this value highly uncertain.

## 5.4. OO Ser

The COM emission is compact and comparable with the size of the continuum disk. While the $T_{\text{ex}}$ values of this target are high enough to be considered a hot corino, its median $T_{\text{ex}}$ (151.6 K) is the lowest among the other targets, which have $T_{\text{ex}}$ values between 177 K and 225 K. The lower temperature might be due to the outburst in this source finishing ~20 years ago. Indeed, Molyarova et al. (2018) showed that after the outburst finishes, different COMs stay in the gas phase for different lengths of time (see their Fig. 5). An extreme example of this is the case of CH$_3$OCHO, which, in some models, had enhanced gas-phase column densities for more than 100 years after the outburst finished.

## 5.5. Impact of the chosen pixels

As discussed in Sect. 3, the spectra analyzed in this study are taken at a pixel position which maximizes the number of detected lines and species. However, the moment maps in Fig. 3 indicate that different COMs could exhibit different spatial distributions, which suggests that abundance ratios relative to CH$_3$OH may vary across the disk. We consequently propose to check it for L1551 IRS 5 North.

For L1551 IRS 5 South, Haro 5a IRS, V346 Nor, and OO Ser we cannot explore this as the emission is not extended enough. At alternative positions located one beam away from the primary pixels, we do not detect most of the COMs.

In contrast, the spectral lines at other positions of L1551 IRS 5 North were strong enough to conduct further analysis. We examined three additional positions in this FUor, at the northeast (NE), southeast (SE), and southwest (SW) sides of the disk, with our original pixel located on the northwest (NW). See Fig. 1. The spectra at these three positions present weaker lines, leading to non-detections of CH$_3$OCH$_3$, CH$_3$COCH$_3$, and CH$_3$$^{C15}$N. For the detected COMs at these positions (CH$_2$DOH, CH$_3$CHO, CH$_3$OCHO, C$_2$H$_5$OH, C$_2$H$_5$CN, and CH$_3$$^{C15}$N), we conducted line fitting by fixing the source size to the same size as for the NW position (0.08", see Table 3).

We followed the same procedure as in Sect. 4 to obtain the column densities at these alternative positions. We used the lines of CH$_3$$^{18}$OH to estimate the column density of methanol. We
<table><caption>**Table 10.** Abundance ratios with respect to CH<sub>3</sub>OH between different positions in L1551 IRS 5 North.</caption><thead><tr><th rowspan="2">Species</th><th>NW</th><th>NW</th><th>NW</th><th>NE</th><th>NE</th><th>SE</th></tr><tr><th>NE</th><th>SE</th><th>SW</th><th>SE</th><th>SW</th><th>SW</th></tr></thead><tbody><tr><td>CH<sub>2</sub>DOH</td><td>4.22</td><td>5.00</td><td>3.43</td><td>1.18</td><td>0.81</td><td>0.69</td></tr><tr><td>CH<sub>3</sub>CHO</td><td>0.99</td><td>1.60</td><td>0.93</td><td>1.62</td><td>0.94</td><td>0.58</td></tr><tr><td>CH<sub>3</sub>OCHO</td><td>1.04</td><td>0.67</td><td>0.79</td><td>0.64</td><td>0.76</td><td>1.18</td></tr><tr><td>C<sub>2</sub>H<sub>5</sub>OH</td><td>0.81</td><td>0.62</td><td>0.79</td><td>0.76</td><td>0.98</td><td>1.29</td></tr><tr><td>C<sub>2</sub>H<sub>5</sub>CN</td><td>0.53</td><td>0.38</td><td>0.46</td><td>0.71</td><td>0.87</td><td>1.23</td></tr><tr><td>CH<sub>3</sub>C<sup>15</sup>N</td><td>0.68</td><td>0.67</td><td>0.68</td><td>0.97</td><td>1.00</td><td>1.02</td></tr></tbody></table>

estimated the abundance ratios for the COMs detected at the three additional positions and compared them with those from the original fit. The results are shown in Table 10.

We found that, aside from CH<sub>2</sub>DOH, relative abundances to methanol were consistent across the four positions, within a factor of 3. For CH<sub>2</sub>DOH, the original position showed higher abundances compared to the three other positions by factors of 3–5, possibly indicating an actual enhancement in the column density of this molecule. Finally, the D/H statistical ratio at these three alternative positions (0.3–0.5%) would be comparable to what was found for SVS13-A, HH212, L483, and V883 Ori in Sect. 5.1.

## 5.6. Comparison with other YSOs

Molyarova et al. (2018) modeled the changes in the chemical composition of a protoplanetary disk due to the effects of a FUor-like luminosity outburst. They found that some molecules are more sensitive to the outburst than others, meaning that for some of them their abundances with respect to H<sub>2</sub> increase by at least two orders of magnitude. Four of the COMs we detected are among the ones analyzed by Molyarova et al. (2018): CH<sub>3</sub>OH, CH<sub>3</sub>CHO, CH<sub>3</sub>OCH<sub>3</sub>, and CH<sub>3</sub>OCHO. They found that the abundances with respect to H<sub>2</sub> of these species increase by at least three orders of magnitude during a FUor-type outburst. Their estimated abundance increases (their Table 1) suggest that for CH<sub>3</sub>CHO, CH<sub>3</sub>OCH<sub>3</sub> and CH<sub>3</sub>OCHO, their abundances with respect to CH<sub>3</sub>OH should increase one and two orders of magnitude during an outburst. Furthermore, Molyarova et al. (2018) found that CH<sub>3</sub>OH freezes out as soon as the outburst finishes. However, in their models with medium sized grains the abundance of CH<sub>3</sub>OH with respect to H<sub>2</sub> experiences a second enhancement after the return to quiescence. The models with large dust grains did not exhibit this behavior. Their models for CH<sub>3</sub>CHO and CH<sub>3</sub>OCHO show that they have enhanced abundances even 20–300 years after the outburst (depending on the species and the model; see their Fig. 5). Therefore, for large dust grains, the abundances of these two COMs with respect to CH<sub>3</sub>OH would be higher in the post-outburst stage than during the outburst (see also Taquet et al. 2016). In contrast, for medium dust grains, the abundances would be similar during and after the outburst.

We compare the chemical composition of these FUor-like eruptive young stars with other Class 0/I objects including quiescent sources from the literature by determining the abundance ratios of COMs with respect to CH<sub>3</sub>OH. The goal is to determine whether the presence of bursts could affect the molecular composition and search for differences in the abundances with respect to methanol that were predicted by Molyarova et al. (2018). The comparison plot is shown in Fig. 6. We limited

this comparison only to sources that have been observed with interferometers, as the single-dish observations can cause significant beam-dilution of the emission line fluxes and are more susceptible to envelope or large scale shock contamination. Our comparison was done against IRAS 16293-2422 (Jørgensen et al. 2016; Manigand et al. 2020), V883 Ori (Lee et al. 2019), L483 (Jacobsen et al. 2019), the CALYPSO survey (Belloche et al. 2020), SVS13-A (Bianchi et al. 2022b), HOPS 108 (Chahine et al. 2022), the PEACHES survey (Yang et al. 2022), and HOPS 373SW (Lee et al. 2023). IRAS 16293-2422, also a protobinary, L1483, HOPS 108, and B335 are Class 0 protostars. As mentioned in Sect. 1, V883 Ori, SVS13-A, HOPS 373SW, and B335 are outbursting protostars. The near-infrared spectrum of V883 Ori is similar to the prototypical FUors but its photometric brightening was not detected, thus it is classified as FUor-like (Connelley & Reipurth 2018). SVS13 is a protobinary that went into outburst in 1990 and has been fading since then. HOPS 373SW began experiencing a powerful outburst in 2020 (Yoon et al. 2022). B335 began its burst between 2010 and 2012, and it lasted until 2023. Lee et al. (2025) found that the line intensities of the COMs were higher during the two epochs during the burst and the one at the end of the event compared to the epoch before the burst. Lee et al. (2025) interpreted this as an indication of an enhancement of the COM column densities due to the sublimation brought by the outburst, and that the COMs had not returned to their pre-burst abundances months after the object returned to quiescence. The PEACHES survey is composed of 50 Class 0/I sources in Perseus. From the CALYPSO survey, we included into our analysis only the five sources with estimated abundances (NGC1333-IRAS2A1, Serps-MM18a, L1448-C, SerpM-S68N, and L1157), which all happen to be Class 0 objects.

For the four eruptive young stars of this work, we used the abundances calculated with respect to CH<sub>3</sub>OH. We also compared the abundances of L1551 IRS 5 North calculated using <sup>13</sup>CH<sub>3</sub>OH and CH<sub>3</sub><sup>18</sup>OH. In the case of HOPS 108 and the two protostars in IRAS 16293, the abundances were computed using CH<sub>3</sub><sup>18</sup>OH and an isotopic ratio [<sup>16</sup>O]/[<sup>18</sup>O] of 560. For V883 Ori (Lee et al. 2019), SVS13-A (Bianchi et al. 2022b), and HOPS 373SW (Lee et al. 2023), their respective authors calculated their abundances using <sup>13</sup>CH<sub>3</sub>OH and an isotopic ratio [<sup>12</sup>C]/[<sup>13</sup>C] of 60. We computed the column density ratios of B335 using the <sup>13</sup>CH<sub>3</sub>OH column densities reported by Lee et al. (2025), and an isotopic ratio [<sup>12</sup>C]/[<sup>13</sup>C] of 68. For the PEACHES survey (Yang et al. 2022), the abundances were calculated with respect to CH<sub>3</sub>OH.

In Fig. 6, we see that the abundance ratios of COMs with respect to CH<sub>3</sub>OH vary depending on the used isotopolog. On average, the ratios increase with the abundance of the methanol isotopolog. While it is clearly shown in the COMs abundances of L1551 IRS 5 North for which we also have the <sup>13</sup>C and <sup>18</sup>O isotopologs of methanol, this effect can also be seen along the sources from the literature. This trend strongly suggests that ratios using the main isotopolog of methanol should be considered with extreme caution. Indeed, most of the lines of the CH<sub>3</sub>OH main isotopolog are optically thick. Studies based on the main isotopolog only consequently leads to overestimates of the column density ratios with respect to CH<sub>3</sub>OH. Comparison using lower abundant isotopologs of CH<sub>3</sub>OH should be favored when possible.

In the northern disk of L1551 IRS 5, we also find that the abundance ratios calculated using <sup>13</sup>CH<sub>3</sub>OH are higher than those computed using CH<sub>3</sub><sup>18</sup>OH. This suggests that opacity effects could also affect the column densities of <sup>13</sup>CH<sub>3</sub>OH if the
<figure><img src="image_7.png" /><figcaption>**Fig. 6.** Comparison of abundances between our four targets and other young stellar objects. The panels from left to right show the abundance ratios with respect to methanol based on CH<sub>3</sub><sup>18</sup>OH, <sup>13</sup>CH<sub>3</sub>OH, and CH<sub>3</sub>OH. The sources considered to be or that have been in outburst (L1551 IRS 5 North, Haro 5a IRS, V346 Nor, OO Ser, VLA4A/B, HOPS 373SW, V883 Ori) are slightly shifted to the right, while the quiescent sources are shifted to the left. In the case of B335 we plotted two epochs, one from 2014 representing the quiescent state and another from 2016 when the source was bursting. We consider the abundance ratios calculated with respect to the main isotopolog of CH<sub>3</sub>OH as upper limits.</figcaption></figure>

assumed isotopic ratios of C and O are correct. Bianchi et al. (2020) suggested that the <sup>13</sup>CH<sub>3</sub>OH lines could be optically thick in this source too. So it is possible that it is also the case in other protostars and some abundances of COMs determined relative to <sup>13</sup>CH<sub>3</sub>OH in the literature may also represent upper limits.

The range of abundances found using <sup>13</sup>CH<sub>3</sub>OH and CH<sub>3</sub><sup>18</sup>OH spans one or two orders of magnitude depending on the species. However, the number of sources available for comparison is very small for a valuable comparison. Regarding the values of IRAS16293 A taken from Manigand et al. (2020), it should be noted that the position chosen to model the species could affect the abundances ratios, as some species (e.g., CH<sub>3</sub>CHO) seem to be more compact than others. The abundance ratios given for this source could consequently be different toward a more central position, but blending issues prevented us from constraining them precisely. The low CH<sub>3</sub>CHO/CH<sub>3</sub>OH ratio seen in Fig. 6 for IRAS16293 A may then be underestimated.

To circumvent the optical thickness issues in CH<sub>3</sub>OH, we also computed abundances with respect to the other COMs. The results are shown in Fig. 7, which do not include L1557 as only CH<sub>3</sub>OH and CH<sub>3</sub>OCHO were detected in it. These abundances also range between one or two orders of magnitude. IRAS16293 A is a common outlier when using CH<sub>3</sub>CHO, due to the aforementioned dependence on the position where the spectrum was extracted. HOPS 108 is also separated from the rest of the sample in a few ratios (CH<sub>3</sub>OCHO/CH<sub>3</sub>CHO, CH<sub>3</sub>OCHO/CH<sub>3</sub>OCH<sub>3</sub>).

Regardless of the COM used to calculate the abundance ratios, no trend is seen between the outbursting and quiescent sources but the sample is too small for any firm conclusion. In order to properly characterize the diversity of chemical composition in the Class 0/I stages of the star formation process, it is necessary to increase the number of protostars (both quiescent and outbursting) with high-quality measurements of multiple COMs, including rarer isotopologs of CH<sub>3</sub>OH.

## 6. Conclusions

We analyzed ALMA Band 6 observations of four FUor-type eruptive young stars, and found emission from multiple COMs and other simple species in each target. The emission of the COMs is generally compact and comparable to the size of the dust disk. We used CASSIS to run a line-fitting procedure assuming LTE, obtaining a best-fit column density *N* and excitation temperature *T*<sub>ex</sub> for each molecule in each protostar. Our best fit resulted in most COMs having *T*<sub>ex</sub> values above 100 K, indicating warm emission. Therefore, we propose that these four FUor-type protostars are hot corinos.

One of these sources, L1551 IRS 5, is a protobinary with COMs emission in both the northern and southern disks. Earlier works have proposed the northern disk to be the one currently experiencing the outburst, which aligns with our results of it having stronger emission from more COMs than the southern disk, and with higher *T*<sub>ex</sub> values. Haro 5a IRS has the weakest line emission, although its luminosity and distance are comparable
<figure><img src="image_8.png" /><figcaption>**Fig. 7.** Abundances of the four FUors and other young stellar objects with respect to different COMs. Here the symbols are also slightly shifted to the left or right if the source is considered to be quiescent or outbursting, respectively.</figcaption></figure>

to those of OO Ser. We expect this to be due to high dust optical depths. V346 Nor is the source where we detected the second-largest amount of species. OO Ser has the lowest $T_{ex}$ values among the outbursting sources, but still high enough to be considered as a hot corino. This could be due to the fact that this protostar experienced a 11-year outburst that finished in 2006 (Kóspál et al. 2007).

We calculated the abundances of COMs with respect to CH₃OH for each target (see Table 9). In the case of the northern disk of L1551 IRS 5, several high excitation lines of CH₃OH were detected. Their opacities are lower than those of the low excitation lines, which are clearly optically thick. The best-fit *N* calculated using these high excitation lines is in agreement, within a factor of 2–3, with the value obtained for ¹³CH₃OH and CH₃¹⁸OH assuming standard isotopic ratios. The ¹³CH₃OH lines of L1551 IRS 5 North could be slightly optically thick, if the assumed isotopic ratios are correct. For the remaining FUors, the CH₃OH emission lines should be considered optically thick and the abundance ratios of COMs as upper limits. Comparisons of COMs abundance ratios with respect to CH₃OH using the rarer ¹⁸O or ¹³C isotopologs should be favored.

The abundance ratios of COMs in L1551 IRS 5 are compared to other values in the literature. The sample of sources with ratios using the ¹⁸O or ¹³C isotopologs of methanol is small, which limits the comparison. The range of values spans 1–2 orders of magnitude for each species (CH₃CHO, CH₃OCH₃, CH₃OCHO, CH₃COCH₃, C₂H₅OH, C₂H₅CN). No trend is seen in the outbursting sources versus the more quiescent hot corinos.

The detection of COMs does not only occur in outbursting stars. However, it is clear that FUor-type eruptive young stars

could be used as experiments to study the chemistry of protostars that otherwise would not be luminous enough to have significant COM detection. Comparing different sources observed using different instrumentation or even spectral configurations is not trivial due to the difference in optical depths or range of $E_{up}$ values analyzed. Therefore, surveys with deep integrations and broad spectral bandwidth that covers multiple transitions from different isotopologs is the best option to study the chemistry of protostars and to carry out comparisons between them.

## Data availability

The appendices are available at https://zenodo.org/records/14967612.

*Acknowledgements.* We are grateful to the anonymous referee for their constructive comments, which significantly enhanced the quality of this paper. This study is part of a project that has received funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation programme (Grant agreements no. 949278, Chemtrip, and no. 716155, SACRED). This work was also supported by the NKFIH NKKP grant ADVANCED 149943 and the NKFIH excellence grant TKP2021-NKTA-64.Project no.149943 has been implemented with the support provided by the Ministry of Culture and Innovation of Hungary from the National Research, Development and Innovation Fund, financed under the NKKP ADVANCED funding scheme. This paper makes use of the following ALMA data: ADS/JAO.ALMA#2016.1.00209.S. ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan), together with NRC (Canada) and NSC and ASIAA (Taiwan) and KASI (Republic of Korea), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO and NAOJ. On behalf of the SACRED project we thank for the usage of MTA Cloud (https://cloud.csfk.org/) that significantly helped us achieving the results published in this paper.
References

Alarcón, F., Casassus, S., Lyra, W., Pérez, S., & Cieza, L. 2024, MNRAS, 527, 9655

ALMA Partnership, Asayama, S., Biggs, A., de Gregorio, I., et al. 2016, ALMA Cycle 4 Technical Handbook, https://almascience.eso.org/documents-and-tools/cycle4/alfa-technical-handbook

Anderson, T., De Lucia, F., & Herbst, E. 1990, ApJS, 72, 797

Andreu, A., Coutens, A., Cruz-Sáenz de Miera, F., et al. 2023, A&A, 677, L17

Artur de la Villarmois, E., Jørgensen, J. K., Kristensen, L. E., et al. 2019, A&A, 626, A71

Belloche, A., Maury, A. J., Maret, S., et al. 2020, A&A, 635, A198

Bergner, J. B., Martín-Doménech, R., Öberg, K. I., et al. 2019, ACS Earth Space Chem., 3, 1564

Bianchi, E., Codella, C., Ceccarelli, C., et al. 2017a, MNRAS, 467, 3011

Bianchi, E., Codella, C., Ceccarelli, C., et al. 2017b, A&A, 606, L7

Bianchi, E., Codella, C., Ceccarelli, C., et al. 2019, MNRAS, 483, 1850

Bianchi, E., Chandler, C. J., Ceccarelli, C., et al. 2020, MNRAS, 498, L87

Bianchi, E., Ceccarelli, C., Codella, C., et al. 2022a, A&A, 662, A103

Bianchi, E., López-Sepulcre, A., Ceccarelli, C., et al. 2022b, ApJ, 928, L3

Bottinelli, S., Ceccarelli, C., Neri, R., et al. 2004, ApJ, 617, L69

Boucher, D., Dubrulle, A., Demaison, J., & Dreizler, H. 1980, Z. Naturfor. A, 35, 1136

Bourdarot, G., Berger, J. P., Lesur, G., et al. 2023, A&A, 676, A124

Brauer, C. S., Pearson, J. C., Drouin, B. J., & Yu, S. 2009, ApJS, 184, 133

Cazaux, S., Tielens, A. G. G. M., Ceccarelli, C., et al. 2003, ApJ, 593, L51

Ceccarelli, C. 2004, in Astronomical Society of the Pacific Conference Series, 323, Star Formation in the Interstellar Medium: In Honor of David Hollenbach, eds. D. Johnstone, F. C. Adams, D. N. C. Lin, D. A. Neufeeld, & E. C. Ostriker, 195

Ceccarelli, C., Caselli, P., Bockelée-Morvan, D., et al. 2014, in Protostars and Planets VI, eds. H. Beuther, R. S. Klessen, C. P. Dullemond, & T. Henning, 859

Chahine, L., López-Sepulcre, A., Neri, R., et al. 2022, A&A, 657, A78

Connelley, M. S., & Reipurth, B. 2018, ApJ, 861, 145

Coutens, A., Jørgensen, J. K., van der Wiel, M. H. D., et al. 2016, A&A, 590, L6

Coutens, A., Willis, E. R., Garrod, R. T., et al. 2018, A&A, 612, A107

Cruz-Sáenz de Miera, F., Kóspál, Á., Ábrahám, P., Liu, H. B., & Takami, M. 2019, ApJ, 882, L4

Cruz-Sáenz de Miera, F., Kóspál, Á., Ábrahám, P., et al. 2023, ApJ, 945, 80

De Simone, M., Codella, C., Testi, L., et al. 2017, A&A, 599, A121

Drążkowska, J., Bitsch, B., Lambrechts, M., et al. 2023, in Astronomical Society of the Pacific Conference Series, 534, Protostars and Planets VII, eds. S. Inutsuka, Y. Aikawa, T. Muto, K. Tomida, & M. Tamura, 717

Eisloeffel, J., Guenther, E., Hessman, F. V., et al. 1991, ApJ, 383, L19

Endres, C. P., Drouin, B. J., Pearson, J. C., et al. 2009, A&A, 504, 635

Fisher, J., Paciga, G., Xu, L.-H., et al. 2007, J. Mol. Spectrosc., 245, 7

Fischer, W. J., Hillenbrand, L. A., Herczeg, G. J., et al. 2022, arXiv e-prints [arXiv:2203.11257]

Fridlund, C. V. M., & Liseau, R. 1998, ApJ, 499, L75

Fukuyama, Y., Odashima, H., Takagi, K., & Tsunekawa, S. 1996, ApJS, 104, 329

Furuya, K., van Dishoeck, E. F., & Aikawa, Y. 2016, A&A, 586, A127

Groner, P., Albert, S., Herbst, E., et al. 2002, ApJS, 142, 145

Hales, A. S., Gupta, A., Ruíz-Rodríguez, D., et al. 2024, ApJ, 966, 96

Herbst, E., & van Dishoeck, E. F. 2009, ARA&A, 47, 427

Herbst, E., Messer, J. K., De Lucia, F. C., & Helminger, P. 1984, J. Mol. Spectrosc., 108, 42

Hodapp, K. W., & Chini, R. 2014, ApJ, 794, 169

Hodapp, K. W., Chini, R., Watermann, R., & Lemke, R. 2012, ApJ, 744, 56

Hsu, S.-Y., Liu, S.-Y., Johnstone, D., et al. 2023, ApJ, 956, 120

Ilyushin, V., Kryvda, A., & Alekseev, E. 2009, J. Mol. Spectrosc., 255, 32

Jacobsen, S. K., Jørgensen, J. K., Di Francesco, J., et al. 2019, A&A, 629, A29

Jensen, S. S., Jørgensen, J. K., Kristensen, L. E., et al. 2021, A&A, 650, A172

Johnson, D. R., Lovas, F. J., Gottlieb, C. A., et al. 1977, ApJ, 218, 370

Jørgensen, J. K., van der Wiel, M. H. D., Coutens, A., et al. 2016, A&A, 595, A117

Jørgensen, J. K., Müller, H. S. P., Calcutt, H., et al. 2018, A&A, 620, A170

Kóspál, Á., Ábrahám, P., Prusti, T., et al. 2007, A&A, 470, 211

Kóspál, Á., Ábrahám, P., Csengeri, T., et al. 2017, ApJ, 843, 45

Kóspál, Á., Szabó, Z. M., Ábrahám, P., et al. 2020, ApJ, 889, 148

Kóspál, Á., Cruz-Sáenz de Miera, F., White, J. A., et al. 2021, ApJS, 256, 30

Kraus, S., Caratti o Garatti, A., Garcia-Lopez, R., et al. 2016, MNRAS, 462, L61

Lee, J.-E., Lee, S., Baek, G., et al. 2019, Nat. Astron., 3, 314

Lee, J.-E., Baek, G., Lee, S., et al. 2023, ApJ, 956, 43

Lee, J.-E., Evans, N. J., Baek, G., et al. 2025, ApJ, 978, L3

Lim, J., Yeung, P. K. H., Hanawa, T., et al. 2016, ApJ, 826, 153

Liseau, R., Fridlund, C. V. M., & Larsson, B. 2005, ApJ, 619, 959

Lykou, F., Ábrahám, P., Chen, L., et al. 2022, A&A, 663, A86

Manigand, S., Jørgensen, J. K., Calcutt, H., et al. 2020, A&A, 635, A48

Marchand, P., Coutens, A., Scigliuto, J., et al. 2024, A&A, 687, A195

Martín-Doménech, R., Bergner, J. B., Öberg, K. I., et al. 2021, ApJ, 923, 155

McMullin, J. P., Waters, B., Schiebel, D., Young, W., & Golap, K. 2007, in Astronomical Society of the Pacific Conference Series, 376, Astronomical Data Analysis Software and Systems XVI, eds. R. A. Shaw, F. Hill, & D. J. Bell, 127

Mercimek, S., Codella, C., Podio, L., et al. 2022, A&A, 659, A67

Milam, S. N., Savage, C., Brewster, M. A., Ziurys, L. M., & Wyckoff, S. 2005, ApJ, 634, 1126

Molyarova, T., Akimkin, V., Semenov, D., et al. 2018, ApJ, 866, 46

Müller, H. S. P., Thorwirth, S., Roth, D. A., & Winnewisser, G. 2001, A&A, 370, L49

Müller, H. S. P., Menten, K. M., & Mäder, H. 2004, A&A, 428, 1019

Müller, H. S. P., Schlöder, F., Stutzki, J., & Winnewisser, G. 2005, J. Mol. Struct., 742, 215

Müller, H. S. P., Drouin, B. J., & Pearson, J. C. 2009, A&A, 506, 1487

Müller, H. S. P., Belloche, A., Xu, L.-H., et al. 2016, A&A, 587, A92

Odashima, H., Matsushima, F., Nagai, K., Tsunekawa, S., & Takagi, K. 1995, J. Mol. Spectrosc., 173, 404

Oldag, F., & Sutter, D. H. 1992, Z. Naturfor. A, 47, 527

Oyama, T., Ohno, Y., Tamanai, A., et al. 2023, ApJ, 957, 4

Pearson, J. C., Sastry, K. V. L. N., Herbst, E., & De Lucia, F. C. 1994, ApJS, 93, 589

Pearson, J. C., Sastry, K. V. L. N., Winnewisser, M., Herbst, E., & De Lucia, F. C. 1995, J. Phys. Chem. Ref. Data, 24, 1

Pearson, J. C., Sastry, K. V. L. N., Herbst, E., & De Lucia, F. C. 1996, J. Mol. Spectrosc., 175, 246

Pearson, J. C., Brauer, C. S., & Drouin, B. J. 2008, J. Mol. Spectrosc., 251, 394

Pearson, J. C., Yu, S., & Drouin, B. J. 2012, J. Mol. Spectrosc., 280, 119

Persson, M. V., Jørgensen, J. K., & van Dishoeck, E. F. 2013, A&A, 549, L3

Pickett, H. M., Cohen, E. A., Brinza, D. E., & Schaefer, M. M. 1981, J. Mol. Spectrosc., 89, 542

Pickett, H. M., Poynter, R. L., Cohen, E. A., et al. 1998, J. Quant. Spec. Radiat. Transf., 60, 883

Reipurth, B., & Aspin, C. 1997, AJ, 114, 2700

Sastry, K. V. L. N., Lees, R. M., & De Lucia, F. C. 1984, J. Mol. Spectrosc., 103, 486

Smirnov, I. A., Alekseev, E. A., Ilyushin, V. V., et al. 2014, J. Mol. Spectrosc., 295, 44

Takakuwa, S., Saigo, K., Matsumoto, T., et al. 2020, ApJ, 898, 10

Taquet, V., Wirström, E. S., & Charnley, S. B. 2016, ApJ, 821, 46

Taquet, V., Bianchi, E., Codella, C., et al. 2019, A&A, 632, A19

Tobin, J. J., Sheehan, P. D., Megeath, S. T., et al. 2020, ApJ, 890, 130

Tobin, J. J., van't Hoff, M. L. R., Leemker, M., et al. 2023, Nature, 615, 227

van't Hoff, M. L. R., Tobin, J. J., Trapman, L., et al. 2018, ApJ, 864, L23

van't Hoff, M. L. R., Harsono, D., Tobin, J. J., et al. 2020, ApJ, 901, 166

Vastel, C., Bottinelli, S., Caux, E., Glorian, J. M., & Boiziot, M. 2015, in SF2A-2015: Proceedings of the Annual meeting of the French Society of Astronomy and Astrophysics, 313

Wilson, T. L., & Rood, R. 1994, ARA&A, 32, 191

Xu, L.-H., & Lovas, F. J. 1997, J. Phys. Chem. Ref. Data, 26, 17

Xu, L.-H., Walsh, M. S., & Lees, R. M. 1996, J. Mol. Spectrosc., 179, 269

Xu, L.-H., Fisher, J., Lees, R. M., et al. 2008, J. Mol. Spectrosc., 251, 305

Yamato, Y., Notsu, S., Aikawa, Y., et al. 2024, AJ, 167, 66

Yang, Y.-L., Evans, N. J., Karska, A., et al. 2022, ApJ, 925, 93

Yoon, S.-Y., Herczeg, G. J., Lee, J.-E., et al. 2022, ApJ, 929, 60

Zwicky, L., Molyarova, T., Akimkin, V., et al. 2024, MNRAS, 527, 7652