© The Authors 2025

1 Introduction

It is essential to investigate the earliest and most deeply embedded stages of star formation to understand (proto)stellar evolution, the properties of surrounding protoplanetary disks, and the formation of embedded planets. A protostellar system in the Class 0/I stage during low-mass star formation is still surrounded by a large mass reservoir that allows further mass growth (Evans et al. 2009; Dunham et al. 2014). Understanding how the envelope material feeds protostellar systems requires multiscale observations on scales from the large envelope to the small disk (Pineda et al. 2023).

A spherically symmetric collapse of an isolated core was first considered in early theoretical models of the formation of a solartype star (e.g., Larson 1969; Penston 1969; Shu 1977). Sensitive observations with a high spatial resolution from the infrared to radio wavelengths over the past decades have revealed that star formation is a complex process in which the environment and asymmetries have to be taken into account, such as filamentary structures in molecular clouds (André et al. 2014; Hacar et al. 2023), multiplicity (Offner et al. 2023), protoplanetary disks (Miotello et al. 2023), chemical diversity (Ceccarelli et al. 2023), and streamers (Pineda et al. 2020). We use the definition that a streamer is an asymmetric structure with kinematically confirmed infall motions. More recent models and simulations that for example include turbulence and magnetic fields also show asymmetries at all scales relevant for star formation (e.g., Walch et al. 2009; Offner et al. 2010; Ballesteros-Paredes et al. 2015; Seifried et al. 2015; Chira et al. 2018; Kuffmeier et al. 2019; Kuznetsova et al. 2022; Kuffmeier 2024).

Asymmetric infalling streamers have been identified in interferometric data in various molecular lines, such as HC₃N (Pineda et al. 2020; Valdivia-Mena et al. 2024), CO (Ginski et al. 2021), ¹³CO (Gupta et al. 2024), C¹⁸O (Yen et al. 2014; Flores et al. 2023; Kido et al. 2023; Mercimek et al. 2023; Cacciapuoti et al. 2024; Kido et al. 2025), HCO⁺ (Yen et al. 2019; Garufi et al. 2022), C₂H (Tanious et al. 2024), H₂CO (Valdivia-Mena et al. 2022; Podio et al. 2024), and CH₃ CN (Fernández-López et al. 2023), and deuterated species such as DCN (Hsieh et al. 2023; Gieser et al. 2024, Cortes et al. subm.). The origin and effects of these infalling structures, including mass infall rates and potential shocks at the landing sites in the disk, have yet to be thoroughly investigated. Most importantly, the connection between streamers and the surrounding envelope material remains unclear because the investigation is challenged by the limited field of view (FOV) and the insufficient line sensitivity of current interferometers.

Recently, a double infalling streamer structure, referred to as the bridge, has been found toward the protostar IRS3A, which is located in L1448N (Gieser et al. 2024). The L1448N starforming complex is located in the Perseus molecular cloud at a distance of 288 ± 6 ± 13 pc (with statistical and systemic uncertainties, respectively, Zucker et al. 2018). This region contains three Class 0/I protostellar systems (IRS3A, IRS3B, and IRS3C) that host a single, triple, and binary system, respectively (Tobin et al. 2018). While IRS3B and IRS3C are Class 0 systems, it is unclear whether IRS3A is in the Class 0 or I stage (e.g., Ciardi et al. 2003). In the following, we therefore refer to all systems as Class 0/I. Continuum observations with a high angular resolution with the Atacama Large Millimeter/submillimeter Array (ALMA) have revealed a ring and spiral-arm structures in the disks of IRS3A and IRS3B, respectively, while the two protostars in the IRS3C binary system each have a compact disk (Tobin et al. 2016; Reynolds et al. 2021, 2024). An overview of the dust continuum distribution of L1448N is presented in Fig. 1. It also includes the orientation of the outflows in the systems (Lee et al. 2016; Reynolds et al. 2021; Gieser et al. 2024; Dunham et al. 2024). The outflow in the IRS3B multiple system can be associated with the protostar IRS3B-C (Reynolds et al. 2021), while currently available data for IRS3C prevent us from distinguishing which protostar in the binary drives the outflow (Tobin et al. 2018).

The three Class 0/I systems of L1448N are target regions of the survey called PROtostars & DIsks: Global Evolution (PRODIGE). This large program with the Northern Extended Millimeter Array (NOEMA) is part of the Max-Planck IRAM Observatory program (MIOP). The NOEMA continuum and molecular line observations at 1.4 mm probe spatial scales down to 300 au at the distance of the Perseus molecular cloud (1″ angular resolution) and trace the environment of the protostel-lar systems. PRODIGE observations of DCN (3 - 2) revealed an extended gas bridge that surrounds systems IRS3A and IRS3B. The bridge is infalling onto IRS3A (Gieser et al. 2024). The gas bridge was even larger in other molecular line tracers such as C¹⁸O (2 - 1) and H₂CO (3_0,3 - 2_0,2). The spatial filtering in the NOEMA observations due to missing short-baseline data prevented us from studying the connection of the gas bridge to the environment and from estimating masses and infall rates of the streamers, however.

To complement the NOEMA data of PRODIGE, we are carrying out the IRAM 30 m telescope large program called “quantifying the connection from the envelope to the disk with MIOP-PRODIGE” (ANTIHEROES). The sensitive single-dish maps allow us to study the large-scale properties of the molecular gas that surrounds the protostellar systems from a few 10 000–3600 au (12″). The combination of the interferometric and single-dish observations allow us to recover the flux on all scales, and thus, we are able to reliably estimate molecular column densities and kinematic properties to 300 au scales for the entire PRODIGE sample.

In this work, we analyze the connection of infalling streamers to the larger-scale material from 36 000 (120″) to 300 au (1″) scales based on a combination of IRAM 30 m telescope and NOEMA molecular line observations to demonstrate the potential of the combined data. The paper is organized as follows: Section 2 contains a description of the PRODIGE and ANTIHEROES data reduction, and the combination of both datasets. The analysis of the molecular line data is presented in Sect. 3, and our results are further discussed in Sect. 4. Our conclusions are summarized in Sect. 5.

Fig. 1 Continuum images toward L1448N. The left panel shows large-scale 1.1 mm continuum emission of L1448 (Enoch et al. 2006). The beam size is shown in the bottom left corner. The map size of the IRAM 30 m observations is highlighted by the dashed gray rectangle. The right panel shows the NOEMA 1.4 mm continuum of L1448N with contour levels at 5, 25, 60, 120, 200 times σcont with σcont = 0.78 mJy beam-1. The triangles and white labels mark the positions of individual protostars (Tobin et al. 2016), and the pink cross marks the recently detected source IRS3-ALMA (Looney et al. 2025). The bipolar outflow directions launched by the three protostellar systems are marked by red and blue arrows, including the blueshifted outflow lobe that is launched by L1448-mm toward the southeast. The synthesized beam of the NOEMA continuum data is shown in the bottom left corner. The scale bars in the two panels are shown in the top left corner.

2 Observations

In this section, we explain the data reduction of the NOEMA and IRAM 30 m datasets and the merging of the data. The molecular line properties (quantum numbers, rest-frequency ν, upper energy level E_u/k_B, and Einstein A-coefficient A_ul) of all transitions we analyzed, and the angular resolution and noise properties of the observations, are summarized in Table A.1.

2.1 NOEMA continuum and molecular line data

The PRODIGE project (program id: L19MB, PIs: P. Caselli and Th. Henning) covers L1448N with NOEMA in the 1.4 mm band. Three individual pointings toward the protostellar systems IRS3A, IRS3B, and IRS3C were obtained from December 2019 to October 2022 in the C- and D-array configurations. The projected baselines range from 15.9 up to 402 m. The frequency coverage of the wideband data is 214.7-222.8 GHz in the lower sideband (LSB) and 230.2-238.3 GHz in the upper sideband (USB) at a channel width of 2 MHz (≈2.7 km s^-1 at 1.4 mm). A total of 39 high spectral resolution units with smaller bandwidths (narrowband data) were placed within the LSB and USB frequency ranges with a channel width of 62.5 kHz (≈0.09 km s^-1 at 1.4 mm). A summary of all molecular transitions covered by the narrowband data that we analyzed is presented in Table A.1.

A detailed description of the NOEMA data-calibration procedure using CLIC¹, including phase self-calibration using MAPPING, of the IRS3A and IRS3B PRODIGE data is presented in Sect. 2.1 in Gieser et al. (2024). For the purpose of this work, we applied the same method to the additional IRS3C pointing. The three NOEMA pointings overlap sufficiently within their primary beam for us to combine all pointings into a mosaic (as described in Sect. 2.2 in Gieser et al. 2024).

Using MAPPING, we imaged the continuum and SO₂ NOEMA data with robust weighting (robust parameter of one) and natural weighting, respectively, and cleaned them with the Hogbom algorithm (Högbom 1974) because the emission is compact (Figs. 1 and 2). Because SO₂ is not detected in the IRAM 30 m data (Fig. 2), we only used the NOEMA data for this molecular transition, as only noise would be added in the data merging. For all other transitions we analyzed (Table A.1), extended emission is detected in the IRAM 30 m data. It is therefore crucial for these transitions to merge the two datasets to recover missing flux from the spatial filtering. The remaining molecular line data were first merged with the IRAM 30 m data, imaged with natural weighting, and then cleaned with the SDI algorithm (Steer et al. 1984), which takes the extended emission into account (as described further in Sect. 2.3). We applied primary beam correction with a primary beam size of ≈22" at 1.4 mm to all data. To increase the signal-to-noise ratio (S/N) of the spectral line data, we rebinned the narrowband data to a common channel width of 0.1 km s^-1.

The continuum mosaic data were cleaned down to three times the theoretical noise with a support mask around the emission of the three protostellar systems. The synthesized beam (θ_maj × θ_min) of the continuum data is 0.″91 × 0.″73 with a position angle (PA) of 18°. The continuum noise level, σ_cont, is 0.78 mJy beam^-1 as estimated from an emission-free area. The SO₂ NOEMA data were cleaned down to the theoretical noise and no support mask. The imaging and cleaning of the remaining molecular line data that were first merged with the IRAM 30 m data is further presented in Sect. 2.3.

2.2 Line data from the IRAM 30 m telescope

To avoid missing extended line emission that may be resolved out by the interferometer, we combined the NOEMA data with complementary zero-spacing spectral line observations from the IRAM 30 m telescope using the Eight MIxer Receiver (EMIR), which covers the same frequency range as the NOEMA data (program code 096-23, PI: C. Gieser). These observations of L1448N are the pilot project for the ongoing ANTIHEROES large program (program code 162-24, PI: C. Gieser), which includes the remaining protostellar systems of the entire PRODIGE sample. The fast Fourier transform spectrometer (FTS) backend provides a channel width of 50 kHz (~0.07kms^-1).

The IRAM 30 m observations of L1448N were carried out in the on-the-fly (OTF) mode in February and April 2024. The central position of the map is α_J2000 = 3^h25^m36^s and δ_J2000 = 30°45′30″. The OTF maps have a size of 120″ × 120″, and the region was observed in zigzag scanning mode (scanning speed 5″ s^-1) in right ascension and declination. The reference offset position (α_J2000 = 3^h26^m40^s and δ_J2000 = 30°50′0″) was selected based on the nondetection of ¹³CO emission in the COordinated Molecular Probe Line Extinction Thermal Emission (COMPLETE) survey data (Ridge et al. 2006). Pointing and focus settings were checked and corrected every ≈1 and ≈2h, respectively, on nearby and bright sources.

The single-dish data were reduced with the CLASS package. In individual spectra, baseline subtraction was carried out by fitting a first-order polynomial around the emission lines. The data were converted from antenna temperature ($T_{\mathrm{A}}^{\ast }$) into main-beam temperature (T_MB) using $T_{\mathrm{M}\mathrm{B}}=\frac{F_{\mathrm{e}\mathrm{f}\mathrm{f}}}{B_{\mathrm{e}\mathrm{f}\mathrm{f}}}\times T_{\mathrm{A}}^{\ast }$ with a forward efficiency F_eff = 0.94 and a main-beam efficiency B_eff = 0.63². The single-dish spectral line data were also rebinned to 0.1 km s^-1 to match the spectral grid of the NOEMA data. Standalone IRAM 30 m data cubes were created using xy_map with a pixel scale of 3″. The half-power beam width (HPBW) is ≈12″ at 1.4 mm, which corresponds to 3600 au at the distance of L1448N.

2.3 Data merging of the line data

We merged the interferometric and single-dish spectral line data using the uv_short task in MAPPING. The merged line data were imaged with natural weighting and cleaned without a support mask down to the theoretical noise using the SDI algorithm. As a second threshold, we set a maximum number of 5000 clean components per channel to avoid overcleaning channels with very extended emission. The pixel scale was set to 0.15″, and the angular resolution was ≈1″ (300 au at the distance of L1448N).

The channel width δv, beam size, and noise level of the IRAM 30 m (θ, σ_line,30m) and merged or NOEMA (θ_maj × θ_min (PA), σ_line) spectral line data are summarized for all transitions in Table A.1. The line sensitivity was ≈100 mK and ≈0.4 K for the IRAM 30 m and merged data, respectively. This dataset provides sensitive molecular line observations from scales of 36 000 au (120″) to 300 au (1″), which is crucial for understanding the origin of molecular gas, such as streamers, that feeds the protostellar systems. The importance of merging the IRAM 30 m data with the NOEMA observations is further evaluated in Appendix B, where we compare the spatial distribution and individual spectra of DCN (3 - 2) and C¹⁸O (2 - 1). The comparison shows that more than 90% is filtered out in the PRODIGE C¹⁸O (2 - 1) NOEMA data compared to the merged ANTIHEROES-PRODIGE data. Most transitions (CO, C¹⁸O, SO, DCN, and c-C₃HC₂) we studied show extended emission. It is therefore crucial to include short-spacing data to properly analyze the kinematic properties and measure reliable fluxes.

Fig. 2 Line-integrated intensity map of C18O (2 - 1), SO (65 - 54), and SO2 (42,2 - 31,3) of the IRAM 30 m data (left) and merged or NOEMA data (right). The line-integrated intensity is shown in color in all panels. The contours and arrows are the same as in Fig. 1. The synthesized beam of the line and continuum data is shown in the bottom left and right corners, respectively. The scale bars are shown in the top left corner in all panels. In the right panels, distinct structures (bridge and stick) seen in molecular emission are indicated by dashed ellipses.

3 Results

We present our analysis of the standalone single-dish data, and the combined NOEMA + IRAM 30 m molecular line data. In Sect. 3.1 we separate the complex line profiles of C¹⁸O and SO₂ in order to extract infalling streamers and connect their kinematic properties to the larger-scale environment. In Sect. 3.2 we use c-C₃H₂ to estimate the temperature in the L1448N region. In Sect. 3.3 we infer the source velocity and central mass of the IRS3C system. A streamline model is applied in Sect. 3.4 to the C¹⁸O and SO₂ streamers found toward IRS3C. The mass and infall rates of the C¹⁸O and SO₂ streamers of IRS3C and DCN streamers of IRS3A (Gieser et al. 2024) are inferred in Sect. 3.5.

3.1 Kinematic properties of the molecular gas

3.1.1 Spatial morphology of C¹⁸O

The line-integrated (v=1-11 km s^-1) intensity maps of the C¹⁸O (2 - 1) transition of the IRAM 30 m and merged data are presented in Fig. 2. In the IRAM 30 m map, C¹⁸O is brightest toward the locations of the three protostellar systems with a NW-SE elongation. In the merged data, this elongated feature shows even more substructure. The bridge structure discovered previously in the PRODIGE data in which IRS3A and IRS3B are embedded is clearly detected (Gieser et al. 2024).

North of IRS3A there is an elongated feature that we refer to as the stick in the following. This structure is nearby the IRS3-ALMA source that was recently detected in millimeter continuum by Looney et al. (2025). Although these authors argued that this source is most likely a background galaxy, the stick structure might be associated with the source given the spatial vicinity. IRS3-ALMA might therefore also be a very faint and young protostar or a prestellar core associated with L1448N. It is beyond the scope of this work to further investigate the nature of IRS3-ALMA.

Elongated structures surround the IRS3C system. The feature connecting IRS3A/B to IRS3C was previously reported by Volgenau et al. (2006) in C¹⁸O (1 - 0) emission (their Fig. 1). Based on the orientation along the outflow direction of IRS3C, the morphology indicates that C¹⁸O traces the remnant envelope that is carved out by the outflow.

3.1.2 Spatial morphology and channel maps of SO and SO₂

The line-integrated (v = 2-6 km s^-1) intensity maps of the SO (6₅ - 5₄) and SO₂ (4_2,2 - 3₁,₃) transitions of the IRAM 30 m and merged and NOEMA data, respectively, are presented in Fig. 2. In the single-dish data, SO₂ is not detected, while the large-scale distribution of SO peaks toward IRS3C, but also east and west of the L1448N region, where no known protostars are located. In the NOEMA data, however, SO₂ is bright toward IRS3C, with a bright arc toward the south of the protostellar system. This elongated arc structure is also detected in bright SO emission.

The channel maps from 2.4 to 5.5 km s^-1 toward IRS3C for the SO (6₅ - 5₄) and SO₂ (4_2,2 - 3_1,3) transitions are presented in Fig. 3. For comparison, the red- and blueshifted morphology of the outflow is highlighted using CO (2 - 1) merged data of ANTIHEROES-PRODIGE (Table A.1). There is compact SO and SO₂ emission at high velocities toward the IRS3C proto-stellar system that rotates southwest (blueshifted) to northeast (redshifted), which is consistent with high angular resolution C¹⁸O data that trace the circumbinary material (Tobin et al. 2018). The redshifted side is strongly blended by extended emission, however: In both tracers, a blob from 3.6 km s^-1 to 4.6 km s^-1 lies north of the IRS3C continuum and the elongated arc structure toward the southwest (marked in the 4.2 km s^-1 channel in Fig. 3). In the SO channel map, extended SO emission also lies close to the region, with velocities between 4.3 and 4.8 km s^-1. The emission of the SO blob disappears at 4.7 and 4.8 km s^-1, which might be due to absorption against the extended component because it is still detected in the SO₂ channel maps. At higher velocities, the blob again becomes bright. The extended emission of SO is known to be abundant in the core envelope because it was also detected in starless and prestel-lar cores (Swade 1989; Spezzano et al. 2017; Vastel et al. 2018; Harju et al. 2020; Lattanzi et al. 2020; Rodríguez-Baras et al. 2021). Although SO is brighter and morphologically more complex, we performed the following kinematic analysis on SO₂ because it lacks additional extended emission.

3.1.3 Gaussian decomposition and clustering

In the following, we analyze the kinematic properties of the molecular gas, which has a complex spatial morphology (Fig. 2) in C¹⁸O and in sulfur-bearing species (SO and SO₂). We study the nature of asymmetric infalling streamers and the connection to the surrounding envelope material in more detail. The dynamics of the infalling gas bridge toward IRS3A have been analyzed in detail using DCN (3 - 2) emission by Gieser et al. (2024) (see also Fig. B.2). While the previous NOEMA-only data could not be used to analyze brighter and more extended molecular emission, such as C¹⁸O, the complementary ANTIHEROES data allow us a detailed kinematic analysis over a large area. We applied the same method as presented in Gieser et al. (2024) to the C¹⁸O (2 - 1) single-dish and merged data, and the SO₂ (4_2,2 - 3_1,3) NOEMA data to infer the large-scale kinematic properties of the L1448N envelope and the smaller-scale properties toward IRS3C.

To obtain the kinematic properties, we fit up to three Gaussian velocity components with the pyspeckit package (Ginsburg et al. 2022) to each position in which the peak S/N in that spectrum was higher than 5. For further details of the method, we refer to Gieser et al. (2024), and a summary is presented in Appendix C.

Following the method by Gieser et al. (2024), we used the tool called density-based spatial clustering of applications with noise (DBSCAN) of the Python package scikit-learn (Pedregosa et al. 2011) to extract individual cluster from the dataset of the extracted velocity components. The input was the position (x,y) and peak velocity 3 of all velocity components after quality assessment (Appendix C). We note that the line width and amplitude were not used as an input, but were used to confirm that the clustering results provided reliable results not only in the velocity, but also in line width and amplitude maps. The main clusters relevant for this work are presented in Fig. 4, and the remaining clusters as well as data points discarded by DBSCAN are shown in Fig. C.2.

On large scales such as traced by the IRAM 30 m data, C¹⁸O has a velocity gradient from east (redshifted, C¹⁸O (IRAM 30 m) cluster 1) to west (blueshifted, C¹⁸O (IRAM 30 m) cluster 2), indicating rotation of the large-scale envelope material (left panel in Fig. 4). The protostellar systems are oriented along the rotation axis of this large-scale motion, that is, they are located at the change in the velocities from red- to blueshifted. This large-scale motion is consistent with the rotation of the disks (Tobin et al. 2018; Reynolds et al. 2021). The infalling DCN bridge also shows the same velocity pattern, with the northern redshifted streamer in the east and the southern blueshifted streamer in the west (Gieser et al. 2024). The agreement with the large-scale C¹⁸O velocity map thus implies that the infalling material from the gas bridge originates from the large-scale environment. We further discuss this result in Sect. 4.2. Toward the protostellar systems, the line width also increases, with ∆v up to 1.5 km s^-1, while in the surrounding environment ∆v <0.5km s^-1. An increase like this in the line width along the streamer toward the protostellar system is found toward other streamers such as IRS3A (Gieser et al. 2024) and Per-emb-50 (Valdivia-Mena et al. 2022).

The large-scale envelope rotation is also seen in the higher angular resolution C¹⁸O dataset (C¹⁸O (merged) cluster 1). The redshifted side of the infalling bridge toward IRS3A can also be extracted in the merged C¹⁸O dataset (C¹⁸O (merged) cluster 6). We note that due to the large number of data points and the complex line profiles close to the infalling bridge, it is not possible to extract the blue side of the IRS3A bridge in C¹⁸O, even though it is evident in the integrated intensity map (Fig. 2). In addition to the rotation of the envelope, an extended emission feature surrounds IRS3C (C¹⁸O (merged) cluster 5), where velocities become more blueshifted as the system is approached from the south. The DCN bridge around IRS3A showed a very similar behavior, that is, the velocities become more red- and blueshifted closer to the protostar. This is evidence for infall motions.

The NOEMA SO₂ clustering results extract the arc structure toward IRS3C (SO₂ (NOEMA) cluster 1) where the emission becomes more redshifted as the system is approached. In higher-resolution ALMA data, Artur de la Villarmois et al. (2023) detected this extended arc in SO and SO₂ emission (with transitions at higher upper energy levels of 48-95 K). This suggests that along the arc itself, the gas is denser and/or hotter than in the rest of the gas. The C¹⁸O and SO₂ asymmetric structures (C¹⁸O (merged) cluster 5 and SO₂ (NOEMA) cluster 1) we extracted from this clustering method and their infall motions onto IRS3C are further characterized in Sect. 3.4.

Fig. 3 Channel maps of SO 65 - 54 (top, merged data) and SO2 42,2 - 31,3 (bottom, NOEMA-only data) toward L1448 IRS3C. The green color map shows the line intensity with the velocity of the channel in km s-1 labeled in the bottom right. The black contours are the same as in Fig. 1. The dashed red and solid blue contours indicate the CO (2 - 1) integrated intensity of the red- and blueshifted outflow launched by IRS3C. The integration ranges are 5.9-27 km s-1 and 2-3.4 km s-1 (narrow range to avoid foreground contamination) with contour levels at 0.1, 0.3, 0.5, and 0.7 and 0.5, 0.7, and 0.9 times the peak-integrated intensity, respectively. In the top left panel, the black scale bar in the top right corner indicates a spatial scale of 1 000 au. The synthesized beam of the line and continuum data are shown in the bottom left and right corners, respectively. In each panel, the scale is set from a minimum brightness temperature of -0.3 K to a maximum of 7 K. The IRS3C source velocity is 3.8 km -1 (Sect. 3.3).

Fig. 4 Main clusters extracted with DBSCAN. The velocity, peak intensity, and line width are shown in color in the top, central, and bottom panels, respectively. The contours and arrows are the NOEMA 1.4 mm continuum with the same levels as in Fig. 1. In cluster 1 of each analyzed dataset, the synthesized beam of the line and continuum data is shown in the bottom left and right corner, respectively, and scale bars are shown in the top left corner. The remaining clusters are shown in Fig. C.2.

3.2 Temperature distribution in L1448N derived from c-C₃H₂

While hot corino sources in the PRODIGE sample show a plethora of emission lines by complex organic molecules (Hsieh et al. 2024; Busch et al. 2025), the PRODIGE spectral setup typically only covers one transition of smaller molecules. One of the exceptions is c-C₃H₂, for which three transitions are detected in the PRODIGE+30m data in the L1448N region (Table A.1). The velocity-integrated intensity maps (integrated between 3 and 7km s^-1) of the IRAM 30 m and merged data are presented in Fig. D.1.

The emission of the three c-C₃ H₂ transitions (Fig. D.1) is spatially extended in the 30 m data, where it embeds the three protostellar systems, similar to the distribution of C¹⁸O. In the high-resolution NOEMA+30m data, the gas bridge, the stick, and the elongated envelope surrounding IRS3C are traced. The extended emission of c-C₃H₂ suggests that L1448N is a source of warm carbon-chain chemistry (WCCC; even though c-C₃H₂ is the cyclic isomer), as suggested previously by Sakai et al. (2009). The authors found a high abundance (>10^-9) of C₄H derived from single-dish observations and proposed that L1448N is a candidate WCCC source.

The underlying excitation conditions, such as rotation temperature and column density, can be inferred by fitting the radiative transfer equation of the rotational transitions. We used the eXtended CASA Line Analysis Software Suite (XCLASS, Möller et al. 2017, version 1.4.3) to fit the observed IRAM 30 m and merged the spectral line data of c-C₃ H₂. We applied the same method to the two datasets as explained in detail in Appendix D. Although the spectral line profiles of c-C₃H₂ show multiple components in some locations (e.g., close to the protostars; Fig. C.1), we chose to fit the observed line profiles with only one velocity component in order to obtain the average excitation conditions in the region. The resulting best-fit parameter maps (T_rot, N, ∆v, v) of c-C₃H₂ are presented in Fig. 5.

The mean excitation temperature in the region is 10 and 14 K in the 30 m and merged data, respectively, indicating that the molecular gas is still cold in the L1448N region over a large area. The merged data shows areas with enhanced temperatures, however, for example, toward the IRS3B disk with temperatures higher than 30 K. The triple system that is embedded in the disk efficiently heats the surrounding gas. The northern envelope material of IRS3C has higher temperatures as well (25 K), which is gas that is heated along the cavity walls of the outflow.

The mean column density in the two datasets is similar, at ≈2×10¹³ cm^-2, but the merged data spatially resolve regions with a higher column density toward the protostellar systems and the surrounding envelope material. Notably, the stick is located very close to IRS3-ALMA (a projected distance of ≈1″, which corresponds to ≈300 au) and also shows higher temperatures and column densities in the merged data. Although we modeled only one c-C₃H₂ emission component, the line width and velocity maps are similar to those extracted from the more sophisticated kinematic analysis of C¹⁸O in Sect. 3.1. Thus, C¹⁸O and c-C₃H₂ trace the same molecular gas components (a spectrum comparing the line profile of molecular analyzed in this work toward all three protostellar systems is shown in Fig. C.1). In Sect. 3.5 we use the temperature determined by the XCLASS fitting of c-C₃H₂ in order to estimate the mass of infalling streamers and the L1448N envelope.

3.3 Systemic velocity and mass of IRS3C

In order to test whether the structures detected toward IRS3C in C¹⁸ O and SO₂ (Fig. 4) are due to gravitational infall using streamline models (Pineda et al. 2020, Sect. 3.4), the precise source velocity and the central mass of IRS3C have first to be determined. We estimated the source velocity using DCN (3 - 2) emission following the approach by Gieser et al. (2024). DCN is only weakly detected toward the continuum peak position of IRS3C (Fig. 6), but the spectral line profile is not as complex as other species, such as C¹⁸O and SO₂ (Fig. C.1). We fit two Gaussian velocity components to the DCN (3 - 2) spectral profile and assigned the brighter component to the IRS3C source velocity, at v_LSR = 3.81 ± 0.05 km s^-1 ■ Despite its low S/N, the second component was included in the fit to avoid biasing the fitted velocity of the brighter component. The peak of the brighter component coincides with the peak emission of other molecular transitions (Fig. C.1) and therefore is a good estimate of the IRS3C source velocity.

The central mass can be estimated from fitting the position velocity (PV) diagram of molecular lines that trace the disk rotation. The C¹⁸O (2 - 1) transition was used for IRS3A (Gieser et al. 2024), but the emission around the IRS3C system is too complex and receives a significant distribution from the surrounding envelope material (Fig. 2). We used the SO₂ (4_2,2 - 3_1,3) transition to estimate the central mass, but because it is faint, we were not able to fit the PV data (Fig. 7). A match by eye of the observations can still provide an estimate of the mass (see also, e.g., Reynolds et al. 2021; Valdivia-Mena et al. 2022). Given the low S/N of the SO₂ PV diagram, more sophisticated PV diagram modeling including rotation and infall is needed, but this is beyond the scope of this work. When we extracted the PV diagram using the pvextractor package (Ginsburg et al. 2015), the blueshifted SO₂ velocities appeared to be roughly compatible with an inclination-corrected mass of ≈1 M_⊙ (Fig. 7). The position angle and inclination of the disk were taken from Reynolds et al. (2024). For comparison, we also show Keplerian models for the inclination-corrected masses of 0.5 M_⊙ and 2 M_⊙.

Fig. 5 Parameter maps of c-C3 H2 (rotation temperature, column density, FWHM line width, and peak velocity) of the IRAM 30 m (left) and merged (right) data. The contours and arrows are the same as in Fig. 1. The beams of the line and continuum data are shown in the bottom left and right corners, respectively. The scale bars are indicated in the top left corner.

Fig. 6 Spectrum of DCN (3 - 2) extracted from the IRS3C continuum peak position The observed spectrum is shown in black, and a two-Gaussian velocity componen fi s hown in bue and red. The dashed purple line marks the 3σline level. The source velocity estimated from the red componen is indicated by the vertical dashed line.

Fig. 7 Position velocity diagram of SO2 (42,2 - 31,3) extracted along the disk orientation of IRS3C. The black contours show the 4, 8, and 16σline levels. The red line shows the profile for a Keplerian disk of 1 M⊙ inclined by 50°, and for comparison, the dotted cyan and green lines show profiles for 0.5 M⊙ and 2 M⊙, respectively. A spectral and spatial resolution element (0.1 km s-1 ×283 au) is indicated in the bottom left corner.

Fig. 8 Streamline models toward IRS3C. In the left panels, the line-integrated intensity and the velocities of the structure are shown in gray and color, respectively. The trajectory of the best-match streamline model is shown in green. The contours are the NOEMA 1.4 mm continuum with the same levels as in Fig. 1. The beam of the line and continuum data are shown in the bottom left and right corners, respectively. A scale bar is shown in the top right corner. In the right panels, we indicate the velocity distribution extracted from the areas by dotted ellipses (left) in red or blue, and the best-match streamline model is shown in green. The source velocity of IRS3C (3.8 km s-1) is indicated by the horizontal dashed line.

3.4 Streamline model

We combined the observed velocity profiles of the C¹⁸O and SO₂ asymmetric structures (C¹⁸O (merged) cluster 5 and SO₂ (NOEMA) cluster 1) in Fig. 4) that showed signs of infall motions toward IRS3C, and the determined source velocity and central mass (Sect. 3.3). We tested with this whether the velocities were consistent with gravitational infall toward IRS3C using a streamline model. The analytical solutions (Mendoza et al. 2009) and the applications to infalling streamers are described in detail in Pineda et al. (2020) and are implemented in the Python package velocity_tools³. In short, the streamline model follows the trajectory of a particle in the rotating envelope and is infalling toward the central gravitational sink, that is, the protostellar system.

We ran a grid of models in which we varied the radial separation r₀, initial position (θ₀, φ₀), radial velocity v₀ = v(r₀), and angular velocity (Ω₀) of the particle. We set the rotation of the system to the observed rotation of the circumbinary envelope (Tobin et al. 2018). We note that the angular resolution is not sufficient to spatially resolve the disks of each binary, and we therefore treated the system as one. This also implies that we were unable to determine whether the streamers feed only one of the protostars in the binary or both.

The observed velocity profiles as a function of projected distance from the protostar are shown in the right panel in Fig. 8. The velocities on the blue side of the C¹⁸O structure and the SO₂ structure become more blue- and redshifted, respectively, as the protostars are approached. For the redshifted side of the C¹⁸O structure, however, the velocity becomes more redshifted with distance from the protostar, which is opposite to what is expected from an infalling structure. Because only the line-of-sight velocity component can be observed, however, it might still be infalling. Figure 8 shows the best-match trajectory projected along the plane of the sky on the left and the observed and streamline velocity profile for the C¹⁸O and SO₂ streamers on the right.

For the blue C¹⁸O and SO₂ infalling structures, the trajectory and velocity profile in the streamline model match very well. For C¹⁸O, no cluster matches the SO₂ streamer velocity profile. This might be due to the high spectral complexity of C¹⁸O line profiles toward the protostellar systems (Fig. C.1). IRS3C is the first example that shows asymmetric infall from the remnant envelope (C¹⁸O) with an additional streamer (bright in SO and SO₂) that cannot be directly connected to the envelope material.

3.5 Streamer mass and infall rates

As reported by Gieser et al. (2024), both the red and blue side of the DCN bridge showed infall motions toward IRS3A. Because no short-spacing data were available at the time, no mass infall rates were determined. A comparison of the NOEMA-only and merged DCN data is presented in Appendix B. The DCN bridge surrounding IRS3A is not affected by spatial filtering, and we therefore use the kinematic results by Gieser et al. (2024) below to estimate the mass and infall rates of the streamers toward IRS3A. For IRS3C, we estimate the mass and infall rates of the C¹⁸O and SO₂ streamers found in this work (Fig. 8).

Molecular column density (N) maps of the C¹⁸O and DCN streamers in IRS3C and IRS3A, respectively, were calculated using an optically thin approximation (Eq. (80) in Mangum & Shirley 2015) within the Python package molecular_columns⁴. For the excitation temperature, we used T_rot from the c-C₃H₂ map, assuming T — 10 K in all pixels without an available c-C₃H₂ temperature estimate. The integrated intensity was calculated using the Gaussian-fit results (intensity and line width) of the corresponding cluster. The abundance ratio X of C¹⁸O/H₂ and DCN/H₂ were assumed to be 1.8 × 10^-7 (Frerking et al. 1982; Wilson & Rood 1994) and 2.4 × 10^-10 (Hsieh et al. 2023), respectively. To compare the streamer masses with the total envelope mass traced by the IRAM 30 m C¹⁸O data, we applied the same method as for the column density calculation using the C¹⁸O-integrated intensity map (Fig. 2) and the temperature from the c-C₃H₂ single-dish map assuming T — 10 K in all pixels without an available c-C₃H₂ temperature estimate.

A column density map of the SO₂ streamer was created using XCLASS (a description of the tool is presented in Sect. 3.2). We fixed the source size and rotation temperature parameters to θ_source = 2″ and T_rot = T_rot(c-C₃H₂) and T_rot=10K otherwise). The line width and velocity were also fixed according to the results from the clustering (SO₂ (NOEMA) cluster 1, Fig. 4). The column density was then the remaining free parameter. The abundance ratio X of SO₂/H₂ was assumed to be 10^-10 (Artur de la Villarmois et al. 2023; Miranzo-Pastor et al. 2025). The column density maps of the streamers and the large-scale envelope are presented in Fig. E.1.

The total mass per pixel was calculated as N × (∆s)² × m(H₂)/X, with the pixel scale ∆s. The total mass of the streamers (M_streamer) and envelope (M_env) is the sum of all pixels. The freefall timescale was calculated as $t_{\mathrm{f}\mathrm{f}}=\sqrt{R_{\mathrm{s}\mathrm{t}\mathrm{r}\mathrm{e}\mathrm{a}\mathrm{m}\mathrm{e}\mathrm{r}}^3/(G\times M_{\ast })}$, with the streamer length R_streamer = r₀ determined from the streamline model (Sect. 3.4). We neglected the contribution of the initial velocity v₀ of the best-match streamline model. The mass M_* of the central protostellar system was estimated from the PV diagram of C¹⁸O and SO₂ for IRS3A (Gieser et al. 2024) and IRS3C (this work, Sect. 3.3), respectively. The mass-infall rate Ṁ_streamer is Ṁ_streamer = ^M_streamer/^t_ff ∙ The total mass of the envelope is 16.4 M_⊙, which is consistent with recent estimates from dust continuum emission (≈15.9M_⊙ Murillo et al. 2024). The streamer properties are summarized in Table 1 and are further discussed in the next section.

The uncertainties of the streamer and envelope mass and the infall rate calculation are dominated by the uncertainties of the abundance ratios. We used average values at the Galactic distance of the Perseus star-forming region for the C¹⁸O/H₂ ratio. For DCN/H₂, we used estimates from the nearby protostar SVS13A, which is also located in the same molecular cloud. Local chemical differences might cause a significant deviation from these values, however. In the case of SO₂/H₂, Artur de la Villarmois et al. (2023) measured a ratio of 1.4 ± 0.5 × 10^-10 toward L1448 IRS3C within an aperture with a diameter of 0.64″, which is consistent with recent estimates from PRODIGE toward the source (1.3 × 10^-11-1.2 × 10^-10, Miranzo-Pastor et al. 2025). Along the streamer, however, this ratio might differ as well. We therefore assumed that the streamer and envelope masses and infall rates are uncertain by a factor of a few.

4 Discussion

The L1448N region consisting of three protostellar systems shows multiple infalling streamers toward IRS3A and IRS3C as revealed by several molecules. A sketch of the streamers in comparison to the bipolar outflow directions and rotation of the disks or circumbinary material in the case of IRS3C is presented in Fig. 9. In the following, we discuss the properties of the streamers and their connection to the larger-scale environment.

4.1 Streamer properties

An overview of the IRS3A and IRS3C properties and the results of M_streamer, t_ff, and M_streamer are summarized in Table 1. The streamers might provide an additional mass of ≈1 and ≈8-17% toward IRS3A and IRS3C, respectively, compared to the central mass M_* of the systems, assuming that all of the streamer mass was accreted onto the protostar. The mass-infall rates are about 10^-6 M_⊙ yr^-1 for IRS3A and considerably higher (5-18×10^-6 M_⊙ yr^-1) for the younger system IRS3C.

While the protostellar mass-accretion rate of IRS3A is 5 × 10^-7 M_⊙ yr^-1 (derived from the luminosity and stellar mass, Reynolds et al. 2021), the infall rates are higher by a factor of two. From mid-infrared imaging, the envelope mass infall rate was estimated to be higher by one order of magnitude at 10^-5 M_⊙ yr^-1 (Ciardi et al. 2003), whereas IRS3B was not detected in the MIR.

The IRS3C infalling streamer detected in SO₂ and SO might indicate a moderate shock and higher temperatures that are caused by the higher infall rate as the material is infalling from the larger-scale environment to the compact envelope and disk system. ALMA observations of SO₂ and SO transitions with higher upper-energy levels also trace the streamer structure (Artur de la Villarmois et al. 2023). This change in the chemical composition where sulfur-bearing species become more abundant was found previously toward the transition zone from the envelope to disk regime toward the protostar IRAS 043681+2557 (Sakai et al. 2014). In most cases in the literature, streamers were found serendipitously in various molecular tracers, and their chemical composition in comparison to the surrounding material was therefore not studied in depth. Because not all streamers originate from the natal envelope, it is not necessary for them to share physical and chemical conditions with the envelope. This analysis is now possible with a large sample of Class 0/I protostellar systems in the Perseus star-forming region using data of ANTIHEROES-PRODIGE.

The infall rates estimated toward other Class 0/I systems are similar; for instance, 10^-6 M_⊙ yr^-1 (Per-emb-2 and Per-emb-50, Pineda et al. 2020; Valdivia-Mena et al. 2022), 6 × 10^-6 M_⊙ yr^-1 (HH212, Lee et al. 2006), 3-6 × 10^-6 M_⊙ yr^-1 (L1527, Flores-Rivera et al. 2021). We note that the mass of the infalling streamers is not high enough to maintain these high infall rates for a long time, and the infall rate is therefore expected to vary with time. This agrees with simulations by Kuffmeier et al. (2019), who showed that structures like this that surround and feed protostellar systems exist up to a few 10 kyr. The IRAM 30 m molecular line maps revealed, however, that a large mass reservoir of ≈ 16 M_⊙ still surrounds the L1448N protostellar systems. It might supply even more material to the protostellar scales in the future.

A similar study that compared the large-scale envelope and streamer material was recently carried out by Kido et al. (2025). The Class 0 protostellar system IRAS 16544-1604 is being fed by three streamers that were detected in C¹⁸O emission with masses of 1-4×10^-3 M_⊙ and mass-infall rates of 1-5×10^-8 M_⊙ yr^-1. Using a Keplerian rotation model, the authors estimated the central mass to be 0.14 M_⊙, which is much lower than for IRS3A and IRS3C (Table 1). The authors estimated a mass of 0.5 M_⊙ in the envelope, which is lower by a factor of 30 than in L1448N. The larger mass reservoir can explain the higher mass-infall rates in the L1448N systems compared to the lower-mass IRAS 165441604 protostar. In agreement with Kido et al. (2025), we found that the streamer masses in L1448N are much lower (<1%) than the overall envelope mass (16.4 M_⊙), which implies that a single streamer event does not deliver a significant amount of mass from the envelope to the smaller scales. A comparison of the streamer mass with the central mass of the protostellar systems (Table 1) shows, however, that a single streamer might significantly impact the masses of the protostellar systems. This might have implications on the disk angular momentum and stability and even on the accretion rate onto the protostar.

Fig. 9 Sketch of L1448N (expanded from Fig. 10 in Gieser et al. 2024).

4.2 Kinematic properties from the envelope to disk scales

We have analyzed the kinematic properties of molecular gas from ≈36 000 au scales to 300 au. We found that the velocity profile of the two DCN streamers in IRS3A (Gieser et al. 2024) and the C¹⁸O streamer in IRS3C (Fig. 8, this work) agree with the large-scale rotation of the L1448N region (C¹⁸O 30m clusters 1 and 2 and C¹⁸O merged cluster 1 in Fig. 4) with blue- and redshifted velocities toward the west and east, respectively. This implies that these streamers are infalling structures that originate from the natal large-scale envelope material. The streamline model explicitly assumes that the particle trajectory originates from the larger-scale rotating envelope (Sect. 3.4). By using the combination of large-scale single-dish maps and sensitive NOEMA observations, we confirmed for the first time that this can be the case for observed streamers. When we compare our results with high-resolution ALMA observations resolving the disk rotation (Tobin et al. 2018; Reynolds et al. 2021), the rotational motion is coherent from the large-scale envelope to the infalling streamers and down to the disk scales. This implies that the protostellar systems have not formed in isolation, but originated from the same envelope material.

The SO₂ streamer that falls onto IRS3C (which is also bright in SO emission; Fig. 3) does not follow the large-scale rotational pattern (Fig. 8), however. The emission of SO and SO₂ is known to trace both accretion and ejection processes (Artur de la Villarmois et al. 2023), and SO is also prominent in the hot inner envelope (Tychoniec et al. 2021). Since the SO₂ line width is not broadened (Fig. 4), the emission is likely not caused by a strong shock. It is more likely that the SO₂ streamer is affected by the hotter inner region of IRS3C, and additional interactions with the outflow material might cause a change in the chemical composition with more abundant sulfur-bearing species. This streamer might originate from a cloudlet-capture event, as was proposed for the streamer in DG Tau (Hanawa et al. 2024). A case of three infalling streamers has been found toward the source IRAS 04239+2436 (Lee et al. 2023), where the spiral structures of the streamer were detected in SO emission, and one in addition in SO₂. By comparing the observational results with magnetohydrodynamic simulations, the authors found that the streamers were consistent with gravitational interactions between the turbulent infalling envelope material and the central triple system and that the emission of S-bearing species originated from low-velocity shocks from these interactions. IRS3C is a binary system, so that this might also explain the SO₂ streamer. The streamers traced by C¹⁸O and DCN in the region, on the other hand, are not bright in SO and SO₂. This might imply a lower density and lower temperatures along these streamers, and shocks might therefore not play an important role. In contrast, the infall rate of the SO/SO₂ streamer is higher. This might cause shocks through interaction with the envelope material, and potentially, even outflows.

The morphology of the large-scale IRAM 30 m map of SO is different from that of C¹⁸O (Fig. 2). These differences between molecular emission following the dust-continuum peak and molecular emission of SO and SO₂ showing asymmetries with respect to the dust continuum can be evident already in the prestellar stage, as was observed toward the prestellar core L1544 (Spezzano et al. 2017). The SO₂ streamer might originate from this gas reservoir. To investigate this potential connection, deep observations of sulfur-bearing species with a larger FOV are required.

5 Conclusions

We analyzed the kinematic properties of the L1448N region, which is located in the Perseus molecular cloud. To do this, we combined single-dish observations with the IRAM 30 m telescope (ANTIHEROES) and interferometric observations with NOEMA (PRODIGE) that trace scales from 36 000 to 300 au.

1. The C¹⁸O emission traces the large-scale envelope material in the single-dish maps, where we spectrally resolved large-scale rotation in which the three protostellar systems (IRS3A, IRS3B, and IRS3C) are embedded. The higher resolution of the combined data resolves highly substructured material.

2. Comparing NOEMA-only and merged (NOEMA+30m) C¹⁸O data, we found that more than 90% of the flux is filtered out in the interferometric data. With the exception of one SO₂ transition, all other lines we analyzed were heavily affected by spatial filtering, with more than 50% missing flux. Not only the total flux, but also the line profiles are severely affected. The ANTIHEROES data therefore provide invaluable information of all molecular lines with extended emission of the entire PRODIGE sample.

3. Temperature maps using the c-C₃H₂ lines indicate a low temperatures (<20 K) on average in the region in the singledish and merged data. Higher temperatures are found in the merged data toward the location of the embedded protostars and toward the cavity walls of IRS3C.

4. Toward IRS3C, we found two infalling streamers that were detected in the C¹⁸O and sulfur-bearing species (SO, SO₂). The velocity profiles and trajectories of the infalling streamers found toward IRS3C are consistent with gravitational infall.

5. The large-scale rotation of the L1448N envelope material traced by C¹⁸O (2 - 1) is coherent down to the streamer kinematics, as is the disk rotation of the three protostel-lar systems. The infalling streamers in C¹⁸O in IRS3C (this work) and DCN in IRS3A (Gieser et al. 2024) accelerate toward the protostellar systems. These infalling streamers are consistent with remnant envelope material that is carved out by the bipolar outflows. The SO₂ streamer of IRS3C shows a different velocity pattern and might originate from an additional nearby gas reservoir that is detected in SO in the single-dish maps.

6. The infall rates of the streamers that fall onto IRS3A and IRS3C are 1 × 10^-6 M_⊙ yr^-1 and 5-18 × 10^-6 M_⊙ yr^-1, respectively. Thus, IRS3C, which is expected to be younger based on its bolometric temperature, accretes at higher accretion rates than IRS3A.

We found that streamers in young Class 0/I systems are connected to the remnant infalling envelope. In the same system, asymmetric gas inflows with an origin that is not consistent with the overall envelope motion can occur at the same time. The ongoing ANTIHEROES program for the first time enables a concise characterization of the origin and chemical composition of infalling material in order to better constrain the nature and impact of streamers in a sample of Class 0/I systems in the Perseus star-forming region. Further sensitive observations with a high angular resolution of molecular line tracers are then necessary to spatially resolve the connection to and effect of the infalling streamers on the disk properties.

Acknowledgements

The authors thank the anonymous referee for their constructive report that helped improve the clarity of this paper. The authors thank the IRAM staff at the NOEMA observatory and IRAM 30 m telescope for their support in the observations and data calibration. This work is based on observations carried out under project numbers L19MB and 096-23 with the IRAM NOEMA Interferometer and IRAM 30 m telescope, respectively. IRAM is supported by INSU/CNRS (France), MPG (Germany) and IGN (Spain). This project is co-funded by the European Union (ERC, SUL4LIFE, grant agreement No101096293). A.F. also thanks project PID2022-137980NB-I00 funded by the Spanish Ministry of Science and Innovation/State Agency of Research MCIN/AEI/10.13039/501100011033 and by “ERDF A way of making Europe”. D.S. has received funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (PROTOPLANETS, grant agreement No. 101002188).

References

Appendix A Observation properties

Table A.1 summarizes the transition properties and observation parameters of all molecular lines analyzed in this work. The data calibration and imaging procedures are described in detail in Sect. 2.

Appendix B Comparison of DCN (3 - 2) and C¹⁸O (2 - 1) NOEMA-only and merged data products

Fig. B.1 Spectra of DCN 3 - 2 (left) and C18O 2 - 1 (right) extracted from the continuum peak positions of IRS3A, IRS3B, and IRS3C. In black the NOEMA only spectrum is shown, while in blue the merged (NOEMA+30m) data are presented. The source vLSR is indicated by the vertical dashed line.

Figure B.1 shows a comparison between DCN (3 - 2) and C¹⁸O (2 - 1) spectra extracted from the continuum peak positions of IRS3A, IRS3B and IRS3C using the NOEMA-only data and the merged data. Line-integrated intensity maps of both transitions comparing the single-dish, merged, and NOEMA-only data are shown in Fig. B.2.

The spatial distribution of DCN is mostly compact tracing the gas bridge toward IRS3A and IRS3B (Gieser et al. 2024), but the additional IRAM 30m short spacing information reveals weak extended emission around IRS3C. This extended emission surrounding IRS3C accounts for in total 57% missing flux in DCN (3 - 2) considering all areas with emission > 5σ in the line-integrated intensity maps (Fig. B.2). The kinematic properties of the gas bridge analyzed in Gieser et al. (2024) using the DCN (3 - 2) transition is not affected by spatial filtering with similar intensities measured in the merged and NOEMA data toward IRS3A and IRS3B.

The emission of C¹⁸O (2 - 1) however is very extended being detected in the entire single-dish map. Comparing the merged and NOEMA-only integrated intensity maps stresses that while compact structures are recovered by NOEMA, most of the extended emission is not. In the C¹⁸O (2 - 1) NOEMA data the missing flux is in total 93%. While the integrated intensity maps of the NOEMA data have large areas with negative values, adding additional 30m data removes these artifacts completely.

In Table A.1 the percentage of missing flux is summarized for all transitions with extended emission used in this work. The amount of missing flux depends heavily on the transition and varies between 57% and 98%. Both the total flux Fig. B.2) and velocity profiles (Fig. B.1) of the C¹⁸O (2 - 1) are different in the merged and NOEMA-only data. Hence, the NOEMA-only data should not be used in the case of molecular emission with extended emission when the analysis relies total flux measurements and/or kinematic properties.

Fig. B.2 Line-integrated intensity map of DCN (3 - 2) and C18O (2 - 1) of the IRAM 30m data (left) and merged data (center), and NOEMA data (right). Red contours highlight emission of -5σ (dotted) and 5σ (solid). If no contours are shown, the emission is > 5σ in the entire FOV. Pink contours and arrows are the same as in Fig. 1. The synthesized beam of the line and continuum data is shown in the bottom left and right corner, respectively. Scale bars are shown in the top left corner.

Appendix C Details of the kinematic analysis

Figure C.1 shows spectra extracted from IRS3A, IRS3B, and IRS3C continuum peak positions comparing the line profiles of DCN (3 - 2), c-C₃H₂ (6_0,6 - 5_1,5), C¹⁸O (2 - 1), SO (6₅ — 5₄), and SO₂ (4_2,2 - 3_1,3). The strong line absorption of DCN and C¹⁸O toward IRS3B is due to absorption against the bright continuum background emission. The high spectral resolution of the molecular line data allows a detailed analysis of the kinematic properties of the gas conducted in Sect. 3.1. Here, a detailed description of the methodology based on the analysis of IRS3A (Gieser et al. 2024) is presented.

Fig. C.1 Molecular line spectra extracted from the continuum peak positions of IRS3A (left), IRS3B (center), and IRS3C (right). The source vLSR is indicated by the vertical dashed line. For SO2 the NOEMA-only data, while for the remaining transitions the merged (NOEMA+IRAM 30m) data is shown (Table A.1).

In Sect. 3.1 we perform a Gaussian decomposition of C¹⁸O (2-1) and SO₂ (4_2,2 - 3_1,3) where up to three velocity components are fitted in each spectrum. Given the increase of noise toward the edges of the FOV, the noise (per 0.1 km s^-1 channel) is estimated in each spectrum individually. The minimum amplitude of each Gaussian velocity component is set to at least 5× the noise and the minimum full-width at half maximum (FWHM) line width is set to 0.1 km s^-1 . We discard all pixels with a defined maximum noise value σ_line,max to mask noisy areas at the edges of the FOV.

We evaluate the goodness of the one-, two-, and three-component fits in each spectrum based on the Akaike Information Criterion (AIC). The best-fit is selected based on the lowest AIC value, but only if the difference to the second lowest value is higher than ∆AIC, defined for each transition. We perform a second quality assessment where we discard all velocity components which have relative uncertainties in any of the fit parameters (peak velocity v, FWHM line width ∆v, amplitude I_peak) higher than a defined percentage. The quality assessment of the fitted Gaussian velocity components is summarized in Table C.1 for all analyzed datasets.

Given the extended emission and complex spectral line profiles of C¹⁸O in both the single-dish and merged data, our criteria are more conservative compared to SO₂, in order to only extract the most reliable velocity components. We note that in the C¹⁸O merged data due to our conservative approach many velocity components are discarded after the quality assessment, even though line emission may be bright, but in many cases line profiles are too complex to fit them properly. This is especially the case toward the complex bridge region surrounding IRS3A and IRS3B where in addition to the infalling gas bridge, there is a significant contribution from outflows and rotation of the compact envelopes and disks (Fig. C.1).

The dataset with up to three velocity components per spectrum is sorted using the clustering algorithm DBSCAN. The clustering in DBSCAN is set by two parameters, the neighborhood parameter eps (ε) and the minimum number of data points within a cluster minsamples. For all datasets we check the parameter space of these two parameters and find a combination of the two that 1) creates smooth velocity maps in each cluster that also produce coherent amplitude and line width maps, 2) has a low number of discarded points, 3) is able to disentangle the velocity components at a certain position, so averaging of the data points within a single pixel is minimized. The input parameters and clustering results are summarized in Table C.1. Additional clusters as well as the discarded data points extracted with DBSCAN in Sect. 3.1 are summarized in Fig. C.2.

For example, large-scale redshifted emission in C¹⁸O single-dish data (C¹⁸O (IRAM 30m) cluster 3) is most likely connected to a background cloud component because there is no known outflow in this direction. The disk rotation of IRS3C is extracted in SO₂ (NOEMA) cluster 2. The redshifted SO₂ blob (SO₂ (NOEMA) cluster 3) toward the northwest might be connected to the arc structure (SO₂ (NOEMA) cluster 1), that is, the arc is wrapping around behind the source where it is infalling. In the merged C¹⁸O dataset, many redshifted blobs are extracted in four clusters, oriented along the redshifted outflow lobe launched by IRS3C (C¹⁸O (merged) clusters 3, 4, 7, and 8). Blueshifted outflow components are seen as well for IRS3B (C¹⁸O (merged) cluster 2) and IRS3C (C¹⁸O (merged) cluster 9).

Fig. C.2 The same as Fig. 4 presenting the main clusters, while here remaining clusters extracted with DBSCAN are shown.

Appendix D Excitation conditions of c-C₃H₂

Integrated intensity maps of all three detected c-C₃ H₂ transitions are shown in Fig. D.1. In Sect. 3.2 we use the XCLASS tool to derive the rotation temperature and column density maps. The detailed setup of the line fitting with XCLASS is described in the following.

Since the noise is not homogeneous and in general increases toward the edges of the FOV, we evaluate the noise in each pixel individually. The noise is computed from the average noise calculated in line-free velocity ranges on the blue- and redshifted sides of each emission line. All pixels with noise values higher than 0.15 K and 0.4 K in the 30m and merged data, respectively, are flagged which correspond to areas at the edge of the FOV. We then compute S/N maps of all transitions and only spectra where the S/N is larger than 6 for all three transitions are fitted with XCLASS.

The myXCLASSFit function solves the radiative transfer equation of an isothermal source for all selected transitions of a molecule assuming local thermal equilibrium (LTE) conditions. The critical densities for the three c-C₃H₂ transitions range between 4 - 5 × 10⁷ cm^-2 (using collisional rates at 30 K from the Leiden Atomic and Molecular Database, van der Tak et al. 2020). In a spherical symmetric envelope approximation, at distances of several hundreds of au from the protostar, densities can reach such high densities (Gerin et al. 2017). At larger distances, the underlying kinetic temperature might be slightly underestimated from this approach. The fit parameter set includes the source size θ_source, rotation temperature T_rot, column density N, FWHM line width ∆v, and peak velocity v.

To reduce the fit parameter set, we fix the source size θ_source to 12″ and 2″ for all 30 m and merged data, respectively. This implies a beam filling factor of 1, which is reasonable since the emission is extended for all c-C₃H₂ transitions (Fig. D.1). The best-fit model (based on the lowest reduced χ²) is determined by using an algorithm chain starting with the Genetic algorithm (100 iterations) followed by the Levenberg-Marquart algorithm (50 iterations) optimizing the parameter set toward global and local minima, respectively.

Fig. D.1 Same as Fig. 2, but for the c-C3 H2 transitions.

Appendix E Column density maps of the envelope and the streamers

Column density maps of the L1448N envelope and the streamers infalling onto IRS3A and IRS3C are presented in Fig. E.1 and are used to estimate the masses in Sect. 3.5.

Fig. E.1 Column density maps of the L1448N and the streamers. Contours and arrows are the same as in Fig. 1. The synthesized beam of the line and continuum data is shown in the bottom left and right, respectively. A scale bar is indicated in the top left corner. The IRS3A DCN streamers were analyzed in Gieser et al. (2024) and IRS3C C18O and SO2 streamers were analyzed in this work.

All Tables

All Figures

Fig. 1 Continuum images toward L1448N. The left panel shows large-scale 1.1 mm continuum emission of L1448 (Enoch et al. 2006). The beam size is shown in the bottom left corner. The map size of the IRAM 30 m observations is highlighted by the dashed gray rectangle. The right panel shows the NOEMA 1.4 mm continuum of L1448N with contour levels at 5, 25, 60, 120, 200 times σcont with σcont = 0.78 mJy beam-1. The triangles and white labels mark the positions of individual protostars (Tobin et al. 2016), and the pink cross marks the recently detected source IRS3-ALMA (Looney et al. 2025). The bipolar outflow directions launched by the three protostellar systems are marked by red and blue arrows, including the blueshifted outflow lobe that is launched by L1448-mm toward the southeast. The synthesized beam of the NOEMA continuum data is shown in the bottom left corner. The scale bars in the two panels are shown in the top left corner.

In the text

Fig. 2 Line-integrated intensity map of C18O (2 - 1), SO (65 - 54), and SO2 (42,2 - 31,3) of the IRAM 30 m data (left) and merged or NOEMA data (right). The line-integrated intensity is shown in color in all panels. The contours and arrows are the same as in Fig. 1. The synthesized beam of the line and continuum data is shown in the bottom left and right corners, respectively. The scale bars are shown in the top left corner in all panels. In the right panels, distinct structures (bridge and stick) seen in molecular emission are indicated by dashed ellipses.

In the text

Fig. 3 Channel maps of SO 65 - 54 (top, merged data) and SO2 42,2 - 31,3 (bottom, NOEMA-only data) toward L1448 IRS3C. The green color map shows the line intensity with the velocity of the channel in km s-1 labeled in the bottom right. The black contours are the same as in Fig. 1. The dashed red and solid blue contours indicate the CO (2 - 1) integrated intensity of the red- and blueshifted outflow launched by IRS3C. The integration ranges are 5.9-27 km s-1 and 2-3.4 km s-1 (narrow range to avoid foreground contamination) with contour levels at 0.1, 0.3, 0.5, and 0.7 and 0.5, 0.7, and 0.9 times the peak-integrated intensity, respectively. In the top left panel, the black scale bar in the top right corner indicates a spatial scale of 1 000 au. The synthesized beam of the line and continuum data are shown in the bottom left and right corners, respectively. In each panel, the scale is set from a minimum brightness temperature of -0.3 K to a maximum of 7 K. The IRS3C source velocity is 3.8 km -1 (Sect. 3.3).

In the text

Fig. 4 Main clusters extracted with DBSCAN. The velocity, peak intensity, and line width are shown in color in the top, central, and bottom panels, respectively. The contours and arrows are the NOEMA 1.4 mm continuum with the same levels as in Fig. 1. In cluster 1 of each analyzed dataset, the synthesized beam of the line and continuum data is shown in the bottom left and right corner, respectively, and scale bars are shown in the top left corner. The remaining clusters are shown in Fig. C.2.

In the text

	Fig. 5 Parameter maps of c-C3 H2 (rotation temperature, column density, FWHM line width, and peak velocity) of the IRAM 30 m (left) and merged (right) data. The contours and arrows are the same as in Fig. 1. The beams of the line and continuum data are shown in the bottom left and right corners, respectively. The scale bars are indicated in the top left corner.
In the text

	Fig. 6 Spectrum of DCN (3 - 2) extracted from the IRS3C continuum peak position The observed spectrum is shown in black, and a two-Gaussian velocity componen fi s hown in bue and red. The dashed purple line marks the 3σline level. The source velocity estimated from the red componen is indicated by the vertical dashed line.
In the text

Fig. 7 Position velocity diagram of SO2 (42,2 - 31,3) extracted along the disk orientation of IRS3C. The black contours show the 4, 8, and 16σline levels. The red line shows the profile for a Keplerian disk of 1 M⊙ inclined by 50°, and for comparison, the dotted cyan and green lines show profiles for 0.5 M⊙ and 2 M⊙, respectively. A spectral and spatial resolution element (0.1 km s-1 ×283 au) is indicated in the bottom left corner.

In the text

Fig. 8 Streamline models toward IRS3C. In the left panels, the line-integrated intensity and the velocities of the structure are shown in gray and color, respectively. The trajectory of the best-match streamline model is shown in green. The contours are the NOEMA 1.4 mm continuum with the same levels as in Fig. 1. The beam of the line and continuum data are shown in the bottom left and right corners, respectively. A scale bar is shown in the top right corner. In the right panels, we indicate the velocity distribution extracted from the areas by dotted ellipses (left) in red or blue, and the best-match streamline model is shown in green. The source velocity of IRS3C (3.8 km s-1) is indicated by the horizontal dashed line.

In the text

	Fig. 9 Sketch of L1448N (expanded from Fig. 10 in Gieser et al. 2024).
In the text

	Fig. B.1 Spectra of DCN 3 - 2 (left) and C18O 2 - 1 (right) extracted from the continuum peak positions of IRS3A, IRS3B, and IRS3C. In black the NOEMA only spectrum is shown, while in blue the merged (NOEMA+30m) data are presented. The source vLSR is indicated by the vertical dashed line.
In the text

Fig. B.2 Line-integrated intensity map of DCN (3 - 2) and C18O (2 - 1) of the IRAM 30m data (left) and merged data (center), and NOEMA data (right). Red contours highlight emission of -5σ (dotted) and 5σ (solid). If no contours are shown, the emission is > 5σ in the entire FOV. Pink contours and arrows are the same as in Fig. 1. The synthesized beam of the line and continuum data is shown in the bottom left and right corner, respectively. Scale bars are shown in the top left corner.

In the text

	Fig. C.1 Molecular line spectra extracted from the continuum peak positions of IRS3A (left), IRS3B (center), and IRS3C (right). The source vLSR is indicated by the vertical dashed line. For SO2 the NOEMA-only data, while for the remaining transitions the merged (NOEMA+IRAM 30m) data is shown (Table A.1).
In the text

	Fig. C.2 The same as Fig. 4 presenting the main clusters, while here remaining clusters extracted with DBSCAN are shown.
In the text

	Fig. D.1 Same as Fig. 2, but for the c-C3 H2 transitions.
In the text

	Fig. E.1 Column density maps of the L1448N and the streamers. Contours and arrows are the same as in Fig. 1. The synthesized beam of the line and continuum data is shown in the bottom left and right, respectively. A scale bar is indicated in the top left corner. The IRS3A DCN streamers were analyzed in Gieser et al. (2024) and IRS3C C18O and SO2 streamers were analyzed in this work.
In the text