© The Authors 2025

1 Introduction

Studying and understanding the evolution of protoplanetary disks from very embedded Class 0 objects to Class II disks is fundamental to unraveling the processes of planet formation. Protoplanetary disks emerge at very early evolutionary stages (Jørgensen et al. 2009; Harsono et al. 2014; Maury et al. 2019; Tobin et al. 2020), with evidence that suggests that planet formation is already underway during these initial phases (Harsono et al. 2018; Tychoniec et al. 2020). Although many of these early-stage disks appear smooth in current observations (Hsieh et al. 2024), substructures such as rings, gaps, and spirals could already be forming but remain undetected or hidden due to resolution limitations, optical depth effects (Ohashi et al. 2023), and temperature effects (Nazari et al. 2025). Even though studying planet formation in Class 0/I sources remains a challenge, primarily because of the optically thick dust continuum emission, it is still possible to study the chemical composition of these young disks, which can provide crucial insights into their evolution.

In this context, numerous Class 0/I protostars have been identified with hot, dense environments known as hot corinos, which are analogs of the hot cores found in high-mass star-forming regions (Ceccarelli 2004). Many of these sources have been found to contain high abundances of complex organic molecules (COMs), species with six or more atoms (Herbst & van Dishoeck 2009; Caselli & Ceccarelli 2012; Ceccarelli 2023). Most of these studies, however, have been conducted at low spatial resolution, limiting our ability to resolve the detailed structure and stratification of molecular species within these young disks. With the unprecedented sensitivity and resolution of the Atacama Large Millimeter/submillimeter Array (ALMA), we are now able to investigate in more detail the molecular distribution of these hot corinos.

As of today, only a handful of Class 0 hot corino sources have been resolved at high angular resolution (∼0.1″): IRAS 162932422 A (Maureira et al. 2022, 2020), SVS13-A (Bianchi et al. 2022), IRAS2AA (Jørgensen et al. 2022), B335 (Okoda et al. 2022; Diaz-Rodriguez et al. 2022; Hsieh et al. 2025), IRAS4A2 (Frediani et al. 2025), HH212 (Lee et al. 2022), and L1448-mm (Nazari et al. 2024). Interestingly, in HH 212, Lee et al. (2022) find molecular stratification, which was attributed to the binding energies of the molecules (sublimation temperature), as the temperature is expected to increase closer to the host star. Studying these young protoplanetary disks, where molecules sublimate at specific temperatures in different regions, can provide insights into where planet formation and dust growth may already be occurring, which we aim to identify in IRAS4A2.

NGC1333 IRAS4A2 was the second hot corino source ever identified (Bottinelli et al. 2004). It is located at a distance of 293 pc (Zucker et al. 2018) in the NGC1333 region of Perseus, and it is part of a binary Class 0 system, separated by 1″.8 from its primary component, IRAS4A1 (Tobin et al. 2018). IRAS4A2 has been extensively studied in the literature (Choi et al. 2010; Taquet et al. 2015; López-Sepulcre et al. 2017; De Simone et al. 2017; Chuang et al. 2021; De Simone et al. 2024; Chahine et al. 2024b); it has a rich inventory of COMs and strong absorption features that partially affect its emission (De Simone et al. 2020b). IRAS4A1 and IRAS4A2 share a common envelope with a total mass of approximately 8 M_⊙ (Maury et al. 2019), which contributes significantly to the millimeter continuum emission. The inner regions of IRAS4A2 are believed to be optically thick at millimeter wavelengths (Galametz et al. 2019; Ko et al. 2020). The system hosts well-characterized bipolar outflows observed in various tracers, including CO, SiO, and SO (Choi et al. 2011; Santangelo et al. 2015; Ching et al. 2016; Chuang et al. 2021). An S-shaped morphology has been observed in the IRAS4A2 outflow, along with an inverted velocity gradient in the inner regions compared to the large-scale outflow, as traced by SO (Santangelo et al. 2015; Chuang et al. 2021). This inner velocity gradient is consistent with that traced by ammonia emission (Choi et al. 2010). Altogether, these features underscore the complexity and rich structure of this particular deeply embedded Class 0 protostellar system.

In this work, we present a detailed study of the IRAS4A2 hot corino, from its large-scale outflows to its innermost compact emission. In Sect. 2, we describe the observations and the process to obtain the spectral line emission. Section 3 discusses our line identification method, along with the techniques used to create the images of the outflows and the extended and compact emission. Finally, in Sect. 4 we explore the possible spatial distribution of these molecules and their implications for IRAS4A2 in the context of dust, disk evolution, and planet formation.

2 Observations

The spatially resolved line observations for IRAS4A2 were obtained using ALMA. Specifically, Band 4 and Band 6 data were gathered under the project code 2018.1.00510.S (PI: James Di Francesco). For a more comprehensive understanding of the observations, including details about data characteristics, calibration, self-calibration, and imaging of the continuum emission, refer to Guerra-Alvarado et al. (2024).

Line spectral windows (SPWs) were separated from the continuum ones, and in this instance, no channel averaging was applied. The self-calibration tables were derived from the continuum analysis and subsequently used in the line SPWs. In Band 6, observations were carried out using two spectral settings, which we refer to as O6B1 and O6B2 for simplicity. O6B1 has a rest frequency of 267.557 GHz and encompasses two SPWs, each featuring 960 channels. These SPWs have central frequencies of 267.568 GHz and 267.569 GHz, respectively. O6B2, on the other hand, has a rest frequency of 265.886 GHz with four SPWs, each with 960 channels and central frequencies of 265.874 GHz, 265.894 GHz, 265.896 GHz, and 265.898 GHz, respectively. In total, the Band 6 data of O6B1 and O6B2 cover frequency ranges from 265.8214 GHz to 265.9385 GHz and from 267.4925 GHz to 267.6103 GHz, respectively.

The Band 4 line data have a rest frequency of 140.839 GHz with two SPWs, each containing 1920 channels with a central frequency of 140.795 GHz and a frequency range between 140.7319 and 140.8512 GHz. Across all observations, the system velocity was assumed to be 6.95 km/s (Di Francesco et al. 2001), which was later used for the analysis.

The final data cubes were produced using the tclean task in CASA, with specmode set to “cube” and the de-convolver to “Multiscale.” We used scale values of 0, 2.5, 10, and 20 times the pixel size, which was fixed at 0.02″ for all line cubes, approximately ten times smaller than the beam size. The spectral resolution is 0.5 km/s for Band 6 observations and 1 km/s for Band 4. These values were selected to balance sensitivity with maintaining the highest possible spectral resolution. Briggs weighting with a robust parameter of 0.5 was adopted for both datasets, as it was found to be the optimal choice based on the cubes images and also in the continuum images presented in Guerra-Alvarado et al. (2024). It also provided the best resolution against other robust values for determining the spatial distribution of the molecular emission, without blurring or missing any large-scale structures or weaker lines. This approach was confirmed through comparison with images made via natural weighting. The final characteristics are shown in Table 1 with typical beam sizes of 0.1″ and rms values of ∼1 mJy/beam per 1 km/s channel.

Finally, for the short baseline data, continuum subtraction was performed prior to imaging using the uvcontsub task in CASA, resulting in three short baseline cubes. For the concatenated data, which has higher resolution and greater sensitivity to small features, determining the continuum level was more challenging. To enhance the accuracy of continuum subtraction and better determine the continuum emission in the image cubes, the STATCONT tool (Sánchez-Monge et al. 2018) was employed. STATCONT utilizes a “corrected sigma-clipping algorithm” for determining continuum emission in line-rich sources like IRAS4A2. This method was applied to all concatenated cubes in the study, resulting in the creation of three other cubes. We used only the short-baseline data to account for the most extended emission, to study the large-scale outflows, and the concatenated (short plus large baseline) data to study the emission surrounding IRAS4A2.

3 Results

Our data revealed emission lines from SO₂ 13_{3 11}-13_{2 12}, OCS 22-21, HCN 3-2, H₂CO 2_{1 2}-1_{1 1} and HCO⁺ 3-2, along with several COMs with upper energy levels of 20 to 487 K (see Table B.1). These emission lines exhibit a variety of spatial components and kinematic structures. In the following section, we first describe the identification of the emission lines around IRAS4A2, we then characterize the spatial distribution of the integrated emission and analyze the kinematics, beginning with the outflow, followed by the extended emission, and concluding with the compact COM emission.

3.1 Analysis and line detection with CASSIS

We conducted an exploration of the spectral emission windows within the concatenated cubes and identified several emission lines originating from a small region near the center of IRAS4A2. Across our cubes, we observed some strong absorption features toward a small region inside the continuum emission. To mitigate this, spectra from all of our data were extracted within CARTA (Cube Analysis and Rendering Tool for Astronomy; Comrie et al. 2018) by selecting a pixel as close as possible to the continuum peak (∼0.08″, 23 au separation), but still outside the region affected by absorption, located at Dec = 31:13:31.9791986 and RA = 3:29:10.4336922. Subsequently, these spectra were converted into brightness temperature and processed for analysis using the CASSIS spectral analysis tool version 5.1.1 (Vastel et al. 2015). Following the methodology outlined in Nazari et al. (2022), we assumed local thermodynamic equilibrium conditions and employed a fitting-by-eye approach. Potential transitions were investigated within CASSIS, utilizing the Cologne Database for Molecular Spectroscopy (CDMS; Müller et al. 2001) and the Jet Propulsion Laboratory (JPL) database (Pickett et al. 1998).

Commonly detected molecules and COMs were prioritized in the search, with upper energy state levels and Einstein coefficient limits set to E_up ≤ 800 K and A_ul ≥ 9 × 10⁻⁶ s⁻¹. Since we extracted the spectra from single pixels, the source size parameter in CASSIS was uniformly set to 1000 for all lines, i.e., unity beam filling factors were assumed. The full width at half maximum (FWHM) of each line was adjusted to fit the spectra better, and the V_lsrk was uniformly set to 6.95 km/s for all spectral cubes.

Following the criteria from Pickett et al. (1998), a molecule was considered detected if at least three lines were found at 3σ otherwise, it was categorized as a “tentative detection.” Figure A.1 displays the observed spectra, where lines were identified, from our cubes in IRAS4A2 alongside our CASSIS model of the lines. Additionally, Table B.1 provides detailed information on all identified molecular lines and their transitions, and Table B.2 summarizes their extent and physical properties, including excitation temperatures and column densities.

Many COMs show little to no absorption toward IRAS4A2, whereas in IRAS4A1 their emission appears completely absorbed by the dust in our data. This difference may reflect their evolutionary stages, with IRAS4A1 being younger (<10³ yr; Santangelo et al. 2015) and IRAS4A2 older (>10⁴ yr; Choi et al. 2010). Recent studies, such as De Simone et al. (2020a), López-Sepulcre et al. (2017), and Guerra-Alvarado et al. (2024), have highlighted the differences in optical depth between the two sources, and some COM emission toward IRAS4A1 has even been observed in absorption against a hot, optically thick dust background Sahu et al. (2019), which may explain the lack of emission in our data of A1. Despite the challenge in observing emission from IRAS4A1, analyzing the diverse species emanating from IRAS4A2 offers valuable insights into the early stages of Class 0 objects and the beginnings of protoplanetary disk formation.

Our IRAS4A2 CASSIS model confirmed the detection of four COMs: C₂H₃CN (vinyl cyanide), CH₂(OH)CHO (glycolaldehyde ), CH₃OCHO (methyl formate), and CH₃C¹⁵N (methyl cyanide). While we were unable to fully fit some of the emission around IRAS4A2, nearly all lines were covered within our frequency range. Tentative detections were also made: the aGa′ glycol line, while providing temperature and column density estimates, was consistently blended with other molecules. For C₂H₅¹³CN, only one detection was made, but two tentative features are consistent only with warmer and denser gas than included in our fit. Finally, the higher frequency line of CH₃NCO is partially cut off in the frequency range, leading to underestimations of the associated excitation temperature and column densities. Despite column densities comparable with those found in López-Sepulcre et al. (2017) for different COMs, the high dust opacity in hot corinos results in underestimations of the abundances (De Simone et al. 2020a). Regardless of these limitations, this work presents the first detection of C₂H₃CN and aGg′-glycol around IRAS4A2.

3.2 Outflows, extended emission, and compact emission in IRAS4A2

3.2.1 Outflows

We employed the code bettermoments for collapsing the concatenated spectral cube, generating moment maps of the outflows of the IRAS4A system (both components) and the emission lines around only IRAS4A2 (Teague & Foreman-Mackey 2018). For a comprehensive view of the large-scale outflows, moment 0 and moment 1 images were generated for HCN, HCO⁺, and H₂CO using the short baseline data. We present these in Fig. 1, providing information about the spatial distribution and velocity structures of these outflows.

These outflows extend to approximately 18″ (∼5300 AU). We resolve both outflows in all lines, which allows for a clear distinction between the emission that originates from IRAS4A1 and IRAS4A2. We observed the distinctive S-shaped pattern (observed in the moment 0 maps) of the IRAS4A2 outflow, consistent with previous findings by Santangelo et al. (2015) and Chuang et al. (2021). This shape is attributed to a misalignment between the initial core angular momentum vector and the magnetic field. It is suggested that this misalignment can happen if the Hall-induced magnetic field flips the angular momentum vector between the flattened envelope and the circumstellar envelope (Zhao et al. 2020; Wurster et al. 2016).

Fig. 1 Integrated emission and velocity maps of the three main outflow tracers in our spectra - HCN, H2CO, and HCO+ - generated using the short baseline line data for IRAS4A. The emission from the outflows of both sources, IRAS4A1 and IRAS4A2, is resolved and can be distinguished separately. The dust continuum emission of IRAS4A1 and IRAS4A2 is outlined with small white contours showing 5 σ and 10 σ levels. The beam is drawn in the lower-left part of each panel, 0.28″ × 0.17″ for HCN and HCO+ and 0.44″ × 0.28″ for H2CO. The red and blue arrows in the middle-left panel indicate the outflow directions of IRAS4A1 (left) and IRAS4A2 (right), respectively (Chahine et al. 2024a).

3.2.2 Extended emission on ~100 au scales

Four lines were identified as tracers of more extended emission in the concatenated data, SO₂, HCN, H₂CO, and OCS, each exhibiting its own distinct morphology. In Fig. 2 we present the moment 0 and 1 maps of these molecules, with white contours outlining the continuum dust emission at 1.3 mm. All molecular emission images of IRAS4A2 are centered at RA = 3:29:10.43, Dec = 1:13:32.01. For these specific spectral lines, we obtained position-velocity (P-V) diagrams to determine the nature of the observed emission, whether it was driven by infall motion, Keplerian rotation from a disk, or a rotating infalling flattened envelope. To analyze this, we used CARTA. We extracted the emission from a cut perpendicular to the outflow axis of IRAS4A2 (29°, Chahine et al. 2024a), extending from the northeast to the southwest. We analyzed the H₂CO, and OCS image by selecting 24 points along the edge of the emission by eye and fitted them using a nonlinear least squares fitting procedure and the following models: For infall motion, we adopted

$$v_{\text{infall}}(r)=\sqrt{2}\cdot \sqrt{\frac{GM}{|r|}},$$(1)

For Keplerian rotation, we used

$$v_{\text{Keplerian}}(r)=\sqrt{\frac{GM}{|r|}},$$(2)

and for a rotating infalling flattened envelope, following the methodology in Sakai et al. (2014) and Lee et al. (2017), we adopted

$$v_r(r)=-\sqrt{\frac{2GM}{r}-\frac{l^2}{r^2}}$$(3)

and

$$v_{\varphi }(r)=\frac{l}{r},$$(4)

where G is the gravitational constant, M is the mass of the central protostar, l is the specific angular momentum, and r is the distance from the central protostar. For the infall and Keplerian models, the central protostellar mass was treated as a free parameter. In the rotating infall model, both the protostellar mass and the specific angular momentum l were fitted.

For SO₂ and HCN, the emission is more complex than what our simple models can fully explain. Therefore, we used the average of the fits from H₂CO and OCS to examine whether any regions of these molecules’ emission remain consistent with these models.

In Fig. 3, we present the P-V diagrams for the four extended molecules, overlaid with Keplerian, infall, and rotating infalling envelope models. While H₂CO and OCS tend to yield relatively low protostellar mass estimates, their values are generally consistent with each other (see the mass labels in Fig. 3). They also show similar angular momentum (l) and centrifugal barrier (CB) radii, which are summarized in Table 2. These similarities provide encouraging evidence that H₂CO and OCS originate from the same physical component, given their comparable dynamical structures.

There are, however, notable differences in the best-fit models for OCS and H₂CO, as well as in the overlay comparisons for HCN and SO₂, which may indicate distinct kinematic structures or physical origins. For OCS and H₂CO, the rotating infalling flattened envelope model reproduces the emission reasonably well. However, in the case of H₂CO, the model fails to capture the low-velocity emission, which could suggest that a purely Keplerian model is more appropriate. In contrast, for HCN and SO₂, both the Keplerian and pure infall models reproduce certain features of the observed emission patterns. Nevertheless, significant emission that is not well described by any of the models considered remains. Some possible contributing factors here are: absorption in certain velocity channels, that the molecules could be tracing multiple kinematic components, or that the molecular lines fitted are blended with other molecules, particularly for SO₂, which is blended with CH₃OCHO. Contamination from CH₃OCHO in the SO₂ line is difficult to quantify, as it appears to be completely blended in velocity. The expected peak of CH₃OCHO emission lies near 3 km/s in the SO₂ P-V diagram; however, it is located so close to the SO₂ and OCS features that it becomes challenging to distinguish where one line ends and the other begins, let alone to accurately assess the contribution from each component.

Notably, the emission that fits between the chosen models aligns well with previous observations of NH₃ (Choi et al. 2010). In contrast, the emission in other quadrants in HCN and SO₂ that cannot be reproduced by any of the models shows a distinct change in the velocity gradient. This shift coincides with the outflow position, as seen in the moment 1 maps in Fig. 2. Similar structures have been reported in other tracers associated with outflows, including SiO (De Simone et al. 2024), H₂CO (Su et al. 2019), and CO and SO (Chuang et al. 2021), suggesting that these molecular lines are likely tracing multiple dynamical components.

Finally, by examining the moment 1 maps, we find that the velocity fields in the inner regions of HCN, H₂CO, OCS, and SO₂ are misaligned (flipped) with those of the outflow (see Fig. 2). This inversion of the velocity maps has been previously observed in various disks before, such as VLA 1623B (Ohashi et al. 2022), BHR71 (Tobin et al. 2019), and IRAS 16293-2422 (Zapata et al. 2013). Multiple origins for this inversion have been proposed, including interaction with the outflow, hydrodynamical interaction with the circumbinary disk, accretion shocks, and, more recently, material accelerated by the disk wind as it crosses the plane of the sky (De Simone et al. 2024). The exact cause of the inverted velocity gradient, however, remains incompletely understood and requires further investigation.

Fig. 2 Integrated emission and velocity maps of SO2, OCS, H2CO, and HCN created using the long+short baseline data. These molecular lines trace a larger extended emission in IRAS4A2. We identified these features visually from the spectra as the most intense and spatially extended emission lines. The white and black contours represent the 6 to 60σ levels of the continuum emission from the ALMA Band 6 and Band 4 data, highlighting the underlying dust distribution. The synthesized beam is shown in the lower-left corner of all panels, and the dashed green lines indicate the cuts used to extract the P-V diagrams perpendicular to the outflow axis.

Fig. 3 P-V diagrams of SO2, OCS, H2CO, and HCN extracted along a cut perpendicular to the jet axis. For OCS and H2CO, we fitted the emission using three kinematic models: pink lines represent the infall motion fit, blue lines the Keplerian rotation fit, and green lines the rotating infalling envelope model. The labels give the mass of the central protostar derived from each fit, and the negative offset marks the northwest direction. For SO2 and HCN, we overlay the model with the average of the fit parameters obtained from the models in OCS and H2CO for comparison, as their emission is more complicated.

3.2.3 Compact emission on sub-100 au scales

Several COMs, like CH₃OCHO, CH₃CDO, C₂H₃CN,CH₃NCO, CH₃¹⁸OH, and more, show much more compact emission, on scales of 20-50 au with a range of morphologies (see Fig. C.1). As an illustration, Fig. 4 displays the moment 0 maps generated by bettermoments for three representative molecules that were not blended with any other line emission. These maps showcase three different morphologies traced within the IRAS4A2 system, indicating that the compact COM emission may arise from different regions like the disk, a warm inner envelope, or both (Tychoniec et al. 2021). However, with the current data, a definitive distinction cannot be made. In addition, Fig. 4 includes the P-V diagram of CH₃CDO and CH₃OCHO, taken along a cut perpendicular to the outflow axis, and the moment 1 map of CH₃OCHO showing the velocity structure of this particular molecule.

To investigate the spatial distribution of molecular emission in IRAS4A2, we analyzed the emission radius using aperture photometry on the moment 0 maps. We applied a signal-tonoise ratio threshold of >5, masking all values below this limit, and then determined the 68% flux radius of the emission, r₆₈. We applied the methodology described in Appendix D of Podio et al. (2021), where we extract a cut perpendicular to the outflow axis going through the continuum peak to measure the spatial extent of the emission. We then assumed the measured r₆₈ corresponds to the FWHM/2 of the emission profile and deconvolved this value from the beam size to correct for beam smearing and avoid biases due to intensity variations across different molecular lines. Figure 5 shows the resulting radius plotted against the excitation temperature derived from CASSIS for molecules with available data. Using the same approach, Fig. 5 also shows the measured radius of each molecule plotted against its respective upper energy levels.

It is important to note that many observed emission lines are blended with others, complicating the estimation of their emission radii and reducing accuracy. Furthermore, some molecular emission exhibits low signal-to-noise ratios or remains unresolved at the current spatial resolution. In these cases, we fixed the FWHM to the minimum beam size to avoid overestimating their spatial extent. Despite these concerns, we chose to proceed with the analysis, as even blended lines, when observed at high resolution, can still provide valuable insights into the chemical structure of young Class 0 objects like IRAS4A2. For a more precise determination of temperature and radius, multiple transitions of unblended molecular emission lines are required, as demonstrated in Frediani et al. (2025).

We note from Fig. 5 that the estimated temperatures are broadly consistent with those reported by Frediani et al. (2025). However, one peculiar result arises for C₂H₅¹³CN (ethyl cyanide, propionitrile). Given its chemical similarity, we would expect this molecule to exhibit comparable behavior to C₂H₃CN (Vinyl cyanide). While the two molecules display similar emitting radii, their derived excitation temperatures differ significantly. In the case of C₂H₅¹³CN, two of the transitions used in the analysis are blended. It is therefore possible that we underestimated the intensity contribution of these blended transitions, which may have resulted in an underestimation of the excitation temperature. Furthermore, we acknowledge that C₂H₅¹³CN is an isotopologue, unlike C₂H₃CN, which could also contribute to the observed differences and should be considered in the interpretation.

Fig. 4 Top panels: integrated intensities of the three representative molecular lines that were not blended with any other, showcasing the three visually identified regions they trace within the IRAS4A2 disk: C2H3CN, CH3CDO, and CH3OCHO. The white contours show the 30 σ to 60 σ levels of the continuum emission, and the synthesized beams are indicated in the lower-left corner of each image. Middle panels: P-V diagrams of CH3CDO and CH3OCHO, extracted along a cut perpendicular to the outflow axis. The dashed green curve represents the Keplerian rotation model from H2CO, while white contours indicate emission levels from 1.5 σ to 3 σ, based on the Band 6 line data. Bottom panel: moment 1 (M1) map of CH3OCHO, indicating the velocity structure across the emitting region.

Fig. 5 Left panel: deconvolved radii of the emission as a function of excitation temperature for the confirmed detected molecules in IRAS4A2. The error bars represent the minimum and maximum radii observed for all molecular lines, while the markers indicate the median radius values. All molecules have at least one blended transition, except for methyl isocyanate, which is highlighted with a black hat on top of its marker to indicate the absence of blending. Right panel: emission radius plotted against the upper state energy levels. Different colors and symbols represent emission from different molecules, as highlighted in the rightmost panel.

4 Discussion

4.1 Spatial scales of the dust and gas of IRAS4A2

In IRAS4A2, the distributions of emission from COMs and other molecules appear to trace distinct physical processes. Identifying these at high resolution can reveal key regions where planet-forming and disk-forming processes occur, as well as the sublimation temperatures that dictate where certain molecules transition from solid to gas. In particular, around IRAS4A2, molecular emission extends from the most inner regions (∼10 au) to large-scale outflows (∼5300 au).

4.1.1 Extended emission and its structure in IRAS4A2

The extended emission in Fig. 2 exhibits different morphologies, with some complexity. First, H₂CO is a well-established tracer of the gas disk in young stellar objects (van’t Hoff et al. 2020; Tychoniec et al. 2021). For CH₃OCHO + SO, the emission would appear to originate from the same region as H₂CO, as they trace the same spatial extent, but it could also be tracing a potential disk wind (De Simone et al. 2024). Meanwhile, OCS appears to be an effective tracer of small grains, especially following the continuum emission bridge previously observed between IRAS4A1 and IRAS4A2 (as shown by the contours in Fig. 2 and in Guerra-Alvarado et al. 2024). For both OCS and SO_2,v=0, it has been suggested that their extended emission in IRAS4A could originate from interactions between the outflow and the surrounding material (Taquet et al. 2020). The integrated emission maps further support the idea that both the grains and OCS are influenced by the S-shaped morphology of the IRAS4A2 outflow. A similar trend is observed for HCN, although in this case, its emission more clearly traces the shape of the disk continuum in the inner regions, suggesting it is simultaneously tracing both the outflow and the gas component of the disk.

On the other hand, H₂CO and SO₂ do not fully align with shape of the continuum emission. They are more compact and do not reach the outflow. This difference may arise from the fact that OCS transitions to the gas phase at lower temperatures than H₂CO and SO₂. Moreover, smaller dust grains, less affected by radial drift, tend to remain distributed over more extended colder regions, making OCS a suitable tracer of regions where only smaller dust grains are present. In contrast, the other two lines require higher temperatures, needing to be closer to the star to reach their sublimation temperatures, and therefore trace more evolved structures. Given this difference, one might expect that OCS and HCN trace a region with small dust grains, while H₂CO and SO₂ likely trace more inner structures, possibly those related to the gas disk. As shown in Fig. 3, the kinematics of H₂CO and OCS appear more similar compared to the other molecules, suggesting they are tracing the same or closely related physical components. However, the moment 1 map of OCS differs in shape from those of the other molecules, particularly at larger offsets. In the inner regions, all molecules display a clearer velocity gradient, but the deviation observed in OCS at greater distances may indicate that its velocity map is influenced by additional components traced by this molecule. Additionally, it is important to note that the fit for H₂CO is not ideal, and much of the negative-velocity emission is not well reproduced. The interpretation for this line therefore remains subject to uncertainties. The main issues probably come from the poor S/N at the current resolution and the absorption features present in some spectral channels. While there is not much we can do about the absorption, future observations with higher S/N and improved spectral resolution would help to obtain a more reliable fit for this particular molecule.

The OCS, HCN, and H₂CO emission exhibits significant absorption by the dust continuum opacity in the inner regions of the IRAS4A2 disk, with HCO⁺ being completely absorbed. It is also worth noting that some of this molecular emission may be subject to self-absorption due to infalling material along the line of sight, which could either contribute to or even account for the observed absorption features, rather than dust absorption alone. In contrast, SO₂ shows little or no absorption compared to the other extended molecules, suggesting it traces higher layers of the source, possibly the envelope. This behavior could indicate that HCN and H₂CO emission is likely tracing a more middle region, while HCO⁺ emission, due to its strong absorption, may trace regions closer to the mid-plane of the source, where dust absorption from the shrinking background might still be significant. This stratification of molecular emission was also observed in HH212 by Lee et al. (2022), in IRAS4A2 by De Simone et al. (2024), and is consistent with findings in HD 163296, where emission from the same molecules were found to trace a similar vertical order (Paneque-Carreño et al. 2023).

In summary, our analysis highlights clear differences in the behavior of these molecules. For example, H₂CO and SO₂ show more compact emission that does not extend into the outflow. H₂CO and OCS share a similar P-V diagram model fit and exhibit aligned masses. Moreover, all extended molecules (OCS, HCN, H₂CO, HCO⁺, and SO₂) display distinct absorption patterns, and both HCN and OCS moment 1 maps suggest the presence of two separate components. Together, these differences, and in some cases, shared behaviors, point to the possibility that the molecules are tracing multiple physical components simultaneously, even on very small scales (<20 au), with variations in vertical distribution across radii and emission extent. However, further observations are required to confirm which physical components are being traced and to disentangle possible chemical and structural stratification within IRAS4A2.

4.1.2 Compact emission and the origin of the dust absorption

As previously discussed, the compact COMs in Fig. 4 also reveal different morphologies. Some molecules, such as C₂H₃CN and CH₃¹⁸OH, exhibit extremely compact emission within radii smaller than 20 au, well below the beam size (0.065″). This indicates that their emission is still very unresolved, and, combined with the low signal-to-noise, makes accurate size estimation particularly challenging.

Despite the uncertainties associated with signal-to-noise and the effects of absorption, the Keplerian model derived from H₂CO provides a reasonable fit to the CH₃CDO and CH₃OCHO emission shown in Fig. 4. Moreover, the CH₃OCHO moment 1 map reveals a clearer velocity gradient perpendicular to the outflow axis, consistent with the rotation seen in the flattened envelope. Combined with the fact that neither emission is aligned with the outflow direction, this suggests that these molecules likely serve as reliable tracers of the disk itself. Further investigation of these two molecules may provide valuable information about the embedded disk and the early conditions enabling planet formation.

Some of these compact regions appear to trace the innermost part of the IRAS4A2 system, as shown in Fig. 5, especially given that emission from a couple of molecules remains unresolved. Moreover, low- and high-energy COM transitions probe a wide range of regions across the system. The most compact emission among these low upper state energy level lines could likely originate from colder regions in higher layers of the disk if the temperatures are lower there than in the mid-plane (Zamponi et al. 2021), and if they are effectively shielded from the protostar’s radiation. Further supporting this interpretation is the observation that some compact and extended molecular emission still exhibits absorption features around IRAS4A2, whereas many others do not. This contrast implies that the emitting and sublimation regions of certain molecules lie in higher and warmer layers of the disk structure, where the dust optically thick absorption appears to be now limited only to the innermost regions with lower scale height.

In Class 0/I objects, the molecular emission often seems to be absorbed over large spatial scales when having an optically thick dust background at submillimeter wavelengths (Galván-Madrid et al. 2018; Galametz et al. 2019). These absorption features have been observed and studied in other star-forming regions as well. For example, Lee et al. (2019) suggest that the absorption and reduction of line emission is due to additional heating sources, which could be altering the temperature structure in the disk. This view comes from observations of systems experiencing accretion outbursts, similar to what has been observed in V8 8 3 Ori (Alarcón et al. 2024; Yamato et al. 2024). In the context of young stellar objects, simulations by Baillié & Charnoz (2014) suggest that viscous heating could significantly impact the disk temperature beyond 20 au, affecting both the mid-plane and higher layers of the disk. Moreover, various hot spots have been identified in other protostellar sources, such as IRAS 162932422 A (Maureira et al. 2022). While the exact origin of these hot spots remains unclear, accretion shocks, particularly those associated with infalling material along filamentary structures (streamers) from the extended envelope, are a primary candidate for generating such localized heating. As highlighted in Maureira et al. (2022), analyzing and identifying these diverse heating mechanisms is crucial for accurately modeling Class 0, and even Class I, protoplanetary disks as these systems have shown signs of late-stage accretion shocks (e.g., Garufi et al. 2022).

While additional heating processes may help explain some of the absorption features observed in IRAS4A2, another important factor to consider is the effect of continuum over-subtraction. It is crucial to distinguish between true absorption, which occurs when lower-excitation material lies in front of or above a warmer, dust continuum (Su et al. 2019), and over-subtraction, an effect that arises when the optically thick continuum and the overlying optically thick gas have similar excitation temperatures, leading to absorption in the line profiles after subtraction of the continuum (see Nazari et al. 2024). In the case of IRAS4A2, although several portions of the spectrum exhibit negative features, some due to true absorption and some to noise, many are most likely the result of over-subtraction. Regardless of the underlying mechanism, this behavior is particularly clearly observed in the inner regions of IRAS4A2, very close to the central star. This behavior could suggest an evolutionary track, as previously mentioned, where strong absorption is now confined to the innermost regions for some molecules, in contrast to IRAS4A1 and some other regions, where much of its emission is still absorbed. This difference indicates that the optically thick dust background has dissipated in the outer regions, has become optically thin in most areas, or has shrunk and moved closer to the star and the mid-plane (as large dust grains are affected by settling and rapid radial drift). While these findings are intriguing, the absorbed emission may obscure critical details, preventing us from seeing the full picture with current observational capabilities. Furthermore, drawing conclusions about high-temperature and more compact molecular tracers remains challenging due to the limitations of both the CASSIS model and our spectral data (see Sect. 3), which explains the lack of a clear correlation between the excitation temperature and the radius of the emission.

4.2 Accretion shocks

OCS and SO₂ emission exhibits strong emission peaks that are offset from the peak of the continuum emission (see Fig. 2). This observed emission pattern could be linked to shocks, as discussed in van Gelder et al. (2021), where sulfur-bearing molecules serve as tracers for material flowing from the envelope into the disk, likely originating from filaments formed by what are commonly referred to as streamers. The formation of SO and SO₂ typically occurs toward the end of a shock, suggesting that their emission is tracing such processes.

The origin of shocks in protostellar systems could be attributed to various mechanisms, including outflows (Persson et al. 2012), accretion shocks (Sakai et al. 2014; Garufi et al. 2022), or disk winds (Tabone et al. 2017). Notably, Artur de la Villarmois et al. (2022) identified SO₂ as the most reliable tracer of accretion shocks, while Maureira et al. (2022) suggested that such accretion shocks could be responsible for generating hot dust spots. Furthermore, OCS and SO₂ emission appears too extended to be solely associated with the envelope (Harsono et al. 2021).

In IRAS4A2, shocks are a plausible mechanism behind the observed emission of OCS and SO₂ . Identifying a single, dominant origin for these shocks remains challenging, however, as they can arise in various contexts within the system. A key distinction lies between accretion shocks, which occur when infalling material decelerates upon impacting the disk surface, generating localized heating, and outflow-induced shocks, which result from the interaction between high-velocity outflows and the surrounding envelope, often producing more spatially extended emission with broader line profiles. Given the observations shown in Fig. 2, we propose that the compact, localized intensity peaks seen in OCS and SO₂ integrated emission maps are most consistent with accretion shocks. Meanwhile, the more extended and broader emission features may instead be associated with outflow-driven shocks (Taquet et al. 2020), reinforcing the idea that emission from these molecules could be tracing multiple components of the system.

It is difficult to determine from the P-V diagrams whether or not SO₂ exhibits the symmetric Keplerian signature expected from a disk wind (Tabone et al. 2017). If the absorption feature in the OCS P-V diagram is ignored, however, a pattern emerges, resembling the expected kinematics (broad and slower than the jet) of a rotating disk wind. This could explain why OCS appears to trace the smallest grains and suggests that it is simultaneously probing multiple components, such as accretion shocks and the disk wind. Nonetheless, determining the precise origin of these emissions remains challenging, and it is still possible that these lines trace a different component than assumed here.

Fig. 6 Normalized molecular abundance ratios relative to methanol for the species detected in this work. The values are compared to those observed in both high-mass star-forming regions - AFGL4176 (B0gelund et al. 2019) and Sgr B2 (Belloche et al. 2017) - and a low-mass star-forming region, IRAS 16293B (Jørgensen et al. 2016).

4.3 High-mass star-forming regions versus IRAS4A2

To understand the origin of hot corinos, and in particular the one associated with IRAS4A2, it is valuable to compare their chemical composition with those of hot cores of other star-forming regions and identify both similarities and differences. Hot corinos have been proposed as low-mass analogs of the hot cores observed in high-mass star-forming regions, sharing comparable chemical signatures despite differences in mass and environmental conditions. In Fig. 6, we compare the molecular abundances relative to methanol for the species detected in IRAS4A2 with those reported in both high-mass (B0gelund et al. 2019; Belloche et al. 2017) and low-mass star-forming regions (Jørgensen et al. 2016). The column density of methanol was adopted from Frediani et al. (2025), using the value that does not account for dust absorption effects for consistency with this work. When necessary, isotopologue abundances were adjusted to reflect the main species using solar system isotopic ratios from Wilson & Rood (1994).

From Fig. 6, IRAS4A2 appears to exhibit molecular abundances similar to those of high-mass star-forming regions such as AFGL4176 and Sgr B2, at least for the molecules where direct comparisons are possible. In the case of C₂H₅CN, where we predict a high molecular abundance, we suspect the CASSIS fitting is unreliable because only one unblended line is clearly detected in our spectra. The remaining lines are either heavily blended or faint, introducing significant uncertainty in the total molecular abundance of this molecule.

When comparing the molecular inventory of IRAS4A2 with that of the high-mass star-forming region AFGL 4176 in B0gelund et al. (2019), we find that the two regions share many molecular species, though AFGL 4176 exhibits a greater number of transitions and lines, probably due to its higher mass. Interestingly, the only molecules not identified in AFGL 4176 are NCC(O)NH₂ and CH₃NCO. CH₃NCO, however, has been detected in other high-mass star-forming regions, including Orion KL and Sgr B2(N) (Cernicharo et al. 2016; Belloche et al. 2017), and a tentative detection of NCC(O)NH₂ has recently been reported in Sgr B2(N) (Li et al. 2024). This absence is not entirely surprising, as previous studies have noted that AFGL 4176 is more comparable to the low-mass protobinary system IRAS 16293 (Jørgensen et al. 2016) than to other high-mass sources like Sgr B2(N) and Orion KL. In contrast, IRAS4A2 appears more aligned with Sgr B2(N), supporting the idea that hot corinos are low-mass analogs of the hot cores in high-mass star-forming regions, sharing similar complex chemistry but on a smaller scale.

Regarding the spatial distribution of these molecules, we observe several similarities with AFGL 4176 B0gelund et al. (2019). First, emission from transitions with lower upper-state energies appears to be more spatially extended than that from higher upper-state energy transitions. Additionally, emission from sulfur-bearing molecules is among the most extended compared to emission from other molecules, which aligns with detections in AFGL 4176 and our spatially resolved observations in IRAS4A2. Finally, we also observe a reverse trend with C₂H₃CN, like in AFGL 4176, with some higher upper-state energy transitions corresponding to larger radii. The reason for this reversal is still unclear and needs more investigation.

5 Summary and conclusions

We have studied the hot corino of IRAS4A2 at high resolution using ALMA, analyzing different spatial scales, from large-scale outflow-like emission down to a small warm inner envelope or disk (∼20 au) emission. Our findings are as follows:

• We confirm the detection of four COMs around IRAS4A2: C₂H₃CN, CH₂(OH)CHO, CH₃OCHO, and CH₃C¹⁵N, along with tentative detections of aGa′-glycol, C₂H₅¹³CN, and CH₃NCO. This is the first time C₂H₃CN and aGg′-glycol are detected toward IRAS4A2;

• Molecular emission appears to trace distinct regions within the IRAS4A2 system; however, determining their exact origins remains challenging with the current data. Several species with small spatial extents and high excitation temperatures are likely associated with the innermost regions of IRAS4A2. In contrast, molecules such as CH₃CDO and CH₃OCHO may serve as effective tracers of the disk structure. Further studies are necessary to determine the potential of these molecules for tracing the disk. If confirmed, however, these molecules could serve as key tracers for studying the disk during its earliest stages, when planet formation may have already begun;

• Given the inconsistencies between absorption features, kinematic signatures (e.g., P-V diagrams), and other observed properties, it is likely that extended emission is tracing multiple components simultaneously. These components include the outflow, envelope, disk, and warm inner envelope, and possibly accretion shocks, outflow shocks, and the disk wind, reflecting the complex and dynamic environment of this young protostellar system;

• The emission peaks observed in the OCS and SO₂ moment 0 maps are more consistent with signatures typically associated with shocks, particularly accretion shocks, which may be responsible for heating the dust in IRAS4A2. This hot dust could explain the prominent absorption features detected across the innermost regions of the spectra. In contrast, the more extended emission coming from the same molecules is better aligned with characteristics of outflow-driven shocks, suggesting that multiple shock mechanisms are contributing to the observed molecular emission;

• A comparison with high-mass star-forming regions such as AFGL 4176, Sgr B2(N), and Orion KL shows that these sources share most of the same molecular species. This commonality suggests that hot corinos do not exhibit significant differences in molecular composition regardless of mass or environment. Furthermore, the spatial distribution of these molecules appears to be similar across different sources, but more studies at high resolution are needed to confirm this;

• We resolved the large-scale outflows (∼5300 au) of HCN, HCO⁺, and H₂CO, distinguishing the outflows from IRAS4A1 and IRAS4A2. We confirm the previously reported S-shaped morphology and the inverted velocity gradient observed in this region between the outflow emission and other components of IRAS4A2.

Acknowledgements

We thank J. Alejandro López-Vázquez, Annaëlle Maury, and Marta de Simone for their helpful discussions and insights. We acknowledge assistance from Allegro, the European ALMA Regional Centre node in the Netherlands. This paper makes use of the following ALMA data: ADS/JAO.ALMA#2018.1.00510.S. ALMA is a partnership of ESO (representing its member states), NSF (USA), and NINS (Japan), together with NRC (Canada), MOST and ASIAA (Taiwan), and KASI (Republic of Korea), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO and NAOJ. P.N. acknowledges support from the ESO and IAU Gruber Foundation Fellowships.

References

Appendix A Spectral analysis and CASSIS modeling of the IRAS4A2 region

Fig. A.1 CASSIS models of the spectral lines in IRAS4A2 (pink) compared with the observed spectra (blue) at frequencies where line emission is detected. The molecules are labeled in black.

Appendix B IRAS4A2 molecules in spectra

Appendix C IRAS4A2 moment 0 maps

Fig. C.1 Integrated emission (moment 0) maps of all the remaining molecular lines identified in the ALMA Band 6 and Band 4 spectra of IRAS4A2.

All Tables

All Figures

Fig. 1 Integrated emission and velocity maps of the three main outflow tracers in our spectra - HCN, H2CO, and HCO+ - generated using the short baseline line data for IRAS4A. The emission from the outflows of both sources, IRAS4A1 and IRAS4A2, is resolved and can be distinguished separately. The dust continuum emission of IRAS4A1 and IRAS4A2 is outlined with small white contours showing 5 σ and 10 σ levels. The beam is drawn in the lower-left part of each panel, 0.28″ × 0.17″ for HCN and HCO+ and 0.44″ × 0.28″ for H2CO. The red and blue arrows in the middle-left panel indicate the outflow directions of IRAS4A1 (left) and IRAS4A2 (right), respectively (Chahine et al. 2024a).

In the text

Fig. 2 Integrated emission and velocity maps of SO2, OCS, H2CO, and HCN created using the long+short baseline data. These molecular lines trace a larger extended emission in IRAS4A2. We identified these features visually from the spectra as the most intense and spatially extended emission lines. The white and black contours represent the 6 to 60σ levels of the continuum emission from the ALMA Band 6 and Band 4 data, highlighting the underlying dust distribution. The synthesized beam is shown in the lower-left corner of all panels, and the dashed green lines indicate the cuts used to extract the P-V diagrams perpendicular to the outflow axis.

In the text

Fig. 3 P-V diagrams of SO2, OCS, H2CO, and HCN extracted along a cut perpendicular to the jet axis. For OCS and H2CO, we fitted the emission using three kinematic models: pink lines represent the infall motion fit, blue lines the Keplerian rotation fit, and green lines the rotating infalling envelope model. The labels give the mass of the central protostar derived from each fit, and the negative offset marks the northwest direction. For SO2 and HCN, we overlay the model with the average of the fit parameters obtained from the models in OCS and H2CO for comparison, as their emission is more complicated.

In the text

Fig. 4 Top panels: integrated intensities of the three representative molecular lines that were not blended with any other, showcasing the three visually identified regions they trace within the IRAS4A2 disk: C2H3CN, CH3CDO, and CH3OCHO. The white contours show the 30 σ to 60 σ levels of the continuum emission, and the synthesized beams are indicated in the lower-left corner of each image. Middle panels: P-V diagrams of CH3CDO and CH3OCHO, extracted along a cut perpendicular to the outflow axis. The dashed green curve represents the Keplerian rotation model from H2CO, while white contours indicate emission levels from 1.5 σ to 3 σ, based on the Band 6 line data. Bottom panel: moment 1 (M1) map of CH3OCHO, indicating the velocity structure across the emitting region.

In the text

Fig. 5 Left panel: deconvolved radii of the emission as a function of excitation temperature for the confirmed detected molecules in IRAS4A2. The error bars represent the minimum and maximum radii observed for all molecular lines, while the markers indicate the median radius values. All molecules have at least one blended transition, except for methyl isocyanate, which is highlighted with a black hat on top of its marker to indicate the absence of blending. Right panel: emission radius plotted against the upper state energy levels. Different colors and symbols represent emission from different molecules, as highlighted in the rightmost panel.

In the text

	Fig. 6 Normalized molecular abundance ratios relative to methanol for the species detected in this work. The values are compared to those observed in both high-mass star-forming regions - AFGL4176 (B0gelund et al. 2019) and Sgr B2 (Belloche et al. 2017) - and a low-mass star-forming region, IRAS 16293B (Jørgensen et al. 2016).
In the text

	Fig. A.1 CASSIS models of the spectral lines in IRAS4A2 (pink) compared with the observed spectra (blue) at frequencies where line emission is detected. The molecules are labeled in black.
In the text

	Fig. C.1 Integrated emission (moment 0) maps of all the remaining molecular lines identified in the ALMA Band 6 and Band 4 spectra of IRAS4A2.
In the text