ALMA-IMF

Rodrigo H. Álvarez-Gutiérrez, Amelia M. Stutz, Nicolás A. Sandoval-Garrido, Fabien Louvet, Frédérique Motte, Roberto Galván-Madrid, Nichol Cunningham, Patricio Sanhueza, Mélisse Bonfand, Sylvain Bontemps, et al.

► To cite this version:

Rodrigo H. Álvarez-Gutiérrez, Amelia M. Stutz, Nicolás A. Sandoval-Garrido, Fabien Louvet, Frédérique Motte, et al. ALMA-IMF. Astronomy & Astrophysics - A&A, 2024, 689, pp.A74. 10.1051/0004-6361/202450321. hal-04783179

HAL Id: hal-04783179

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

Submitted on 15 Nov 2024

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

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

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

Distributed under a Creative Commons Attribution 4.0 International License
ALMA-IMF

XIII. N$_{2}$H$^{+}$ kinematic analysis of the intermediate protocluster G353.41

R. H. Álvarez-Gutiérrez<fnref n="1" />, A. M. Stutz<fnref n="1" />, N. Sandoval-Garrido<fnref n="1" />, F. Louvet², F. Motte², R. Galván-Madrid³,
N. Cunningham²,⁴, P. Sanhueza⁵,⁶, M. Bonfand⁷, S. Bontemps⁸, A. Gusdorf⁹,¹⁰, A. Ginsburg¹¹,
T. Csengeri⁸, S. D. Reyes<fnref n="1" />,¹², J. Salinas<fnref n="1" />, T. Baug¹³, L. Bronfman¹⁴, G. Busquet¹⁵,¹⁶,¹⁷
D. J. Díaz-González³, M. Fernandez-Lopez¹⁸, A. Guzmán¹⁹, A. Koley<fnref n="1" />, H.-L. Liu²⁰, F. A. Olguin²¹,
M. Valeille-Manet⁸, and F. Wyrowski²²

(Affiliations can be found after the references)

Received 10 April 2024 / Accepted 12 June 2024

ABSTRACT

The ALMA-IMF Large Program provides multi-tracer observations of 15 Galactic massive protoclusters at a matched sensitivity and spatial resolution. We focus on the dense gas kinematics of the G353.41 protocluster traced by N<sub>2</sub>H<sup>+</sup> (1&ndash;0), with a spatial resolution of ~0.02 pc. G353.41, at a distance of ~2 kpc, is embedded in a larger-scale (~8 pc) filament and has a mass of ~2.5 × 10<sup>3</sup> M<sub>☉</sub> within 1.3 × 1.3 pc<sup>2</sup>. We extracted the N<sub>2</sub>H<sup>+</sup> (1&ndash;0) isolated line component and decomposed it by fitting up to three Gaussian velocity components. This allows us to identify velocity structures that are either muddled or impossible to identify in the traditional position-velocity diagram. We identify multiple velocity gradients on large (~1 pc) and small scales (~0.2 pc). We find good agreement between the N<sub>2</sub>H<sup>+</sup> velocities and the previously reported DCN core velocities, suggesting that cores are kinematically coupled with the dense gas in which they form. We have measured nine converging “V-shaped” velocity gradients (VGs) (~20 km s<sup>-1</sup> pc<sup>-1</sup>) that are well resolved (sizes ~0.1 pc), mostly located in filaments, which are sometimes associated with cores near their point of convergence. We interpret these V-shapes as inflowing gas feeding the regions near cores (the immediate sites of star formation). We estimated the timescales associated with V-shapes as VG<sup>-1</sup>, and we interpret them as inflow timescales. The average inflow timescale is ~67 kyr, or about twice the free-fall time of cores in the same area (~33 kyr) but substantially shorter than protostar lifetime estimates (~0.5 Myr). We derived mass accretion rates in the range of (0.35&ndash;8.77) × 10<sup>-4</sup> M<sub>☉</sub> yr<sup>-1</sup>. This feeding might lead to further filament collapse and the formation of new cores. We suggest that the protocluster is collapsing on large scales, but the velocity signature of collapse is slow compared to pure free-fall. Thus, these data are consistent with a comparatively slow global protocluster contraction under gravity, and faster core formation within, suggesting the formation of multiple generations of stars over the protocluster’s lifetime.

**Key words.** ISM: clouds – ISM: kinematics and dynamics – ISM: molecules – ISM: structure

**1. Introduction**

While star clusters have been studied extensively over many
decades at comparatively short wavelengths, their precursors,
protoclusters, have not been studied in depth until recently. Pro-
toclusters (or embedded clusters) are the gas-dominated mater-
nal environments in which star clusters are born and whose
stellar constituents will ultimately populate the field of our
Galaxy. Protoclusters are distinct entities from star clusters. Both
are defined as relatively compact configurations in which the
gravity is strong enough to influence the dynamics of their con-
stituents. But in the latter, there is little to no gas, and the gravity
of the cluster is dominated by the stars themselves. In protoclus-
ters, in contrast, gravity is dominated by the cold gas in which
the stars themselves are forming (Stutz & Gould 2016; Csengeri
et al. 2017; Stutz 2018; Motte et al. 2018). Protoclusters are
more accessible now than ever before thanks to ALMA and its
exquisitely high-resolution interferometric millimeter-wave data
tracing the cold gas in which the stars form (Sanhueza et al. 2019;
Liu et al. 2020; Motte et al. 2022). Inside protoclusters, we wit-
ness the ongoing conversion of gas into compact and extremely
dense stars, a process mediated by gas filaments (Stutz 2018;
González Lobos & Stutz 2019; Álvarez-Gutiérrez et al. 2021)

feeding gas structures called “cores” (André et al. 2010; Stutz
& Kainulainen 2015; Kuznetsova et al. 2015, 2018). Cores are
compact gas mass concentrations, often defined to be of a size
matching the resolution limit of the observations. In this case,
we define cores to be in the order ~2 kau, for reasons described
below.

In this paper, we focus on the G353.41 protocluster (see Fig. 1), and in particular, on the dense gas kinematics observationally accessible from the protocluster scale (2.9 pc²) to the core scale. We trace this dense and cold gas using the N₂H⁺ (1–0) line observed with ALMA. Given that N₂H⁺ is detected at high column densities (N(H₂) ≥ 10²² cm⁻²; Tafalla et al. 2021), we gain access to the inner dense gas “skeleton” of the protocluster structure, free from confusion induced by lower-density gas. Meanwhile, ALMA permits us to obtain the resolution needed to trace structures down to the core scales at which individual or small numbers of stars may be forming.

The ALMA-IMF Large Program<fnref n="1" /> (LP) maps 15 dense,
neighbby (2–5.5 kpc), and massive (2–32 × 10³ M☉) Milky
Way protoclusters down to ~2 kau scales (Motte et al. 2022),
at a matched spatial resolution. ALMA-IMF provides a large

<footnote n="1">¹ Proposal ID 2017.1.01355.L, PIs: Motte, Ginsburg, Louvet, Sanhueza.</footnote>
<figure><img src="image_2.png" /><figcaption>**Fig. 1.** Composite image of G353: IRAC 3.6 µm (in blue), 4.5 µm (green), and 5.8 µm (red). We indicate the ALMA-IMF N<sub>2</sub>H<sup>+</sup> (1–0) coverage with a light green contour. We highlight ATLASGAL emission (870 µm) at 40 mJy beam<sup>-1</sup> with the gray contour, corresponding roughly to a Herschel-derived N(H) of ~5.5 10<sup>22</sup> cm<sup>-2</sup>.</figcaption></figure>

protocluster sample in order to test the universality of the stellar initial mass function (IMF) (Bastian et al. 2010; Offner et al. 2014). The ALMA-IMF LP also provides a vast catalog of molecular lines, in bands 3 (2.6–3.6 mm) and 6 (1.1–1.4 mm). This rich molecular treasure trove allows for a detailed kinematic characterization of the gas, protostellar cores, and young stellar objects (YSOs) present in these protoclusters. The current publicly available ALMA-IMF data include, but are not limited to, continuum maps (Ginsburg et al. 2022; Díaz-González et al. 2023), 12 m data cubes of all spectral windows (Cunningham et al. 2023), core catalogs (Pouteau et al. 2023; Louvet et al. 2024), and hot core and outflow catalogs (Cunningham et al. 2023; Nony et al. 2023; Towner et al. 2024; Armante et al. 2024; Bonfand et al. 2024, Valeille-Manet et al. in prep). The data products derived from the ALMA-IMF LP allow us to constrain the different star-forming environments, in which we can analyze column densities, temperatures, outflow masses, core properties, and multi-tracer gas kinematics. This approach offers a thorough characterization of the processes taking place in these regions.

Motte et al. (2022) present a method of classifying these 15 protoclusters based on their evolutionary stage, assuming that they exhibit more H II regions as they evolve. They take into account the flux ratio between the 1 mm to 3 mm continuum maps ($S_{1.3\,\text{mm}}/S_{3\,\text{mm}}^{\text{cloud}}$), and the free-free emission at the frequency of H$\alpha$ ($\Sigma_{H41\alpha}^{\text{free-free}}$). They find that as protoclusters evolve, $S_{1.3\,\text{mm}}/S_{3\,\text{mm}}^{\text{cloud}}$ decreases, while $\Sigma_{H41\alpha}^{\text{free-free}}$ increases (Motte et al. 2022, see their Fig. 3). Using these constraints, they group their 15 protocluster as being in a young, intermediate, or evolved evolutionary state. Out of these 15 regions, we analyze the G353.41 protocluster (hereafter G353). In Fig. 1, we indicate the ALMA-IMF N<sub>2</sub>H<sup>+</sup> (1–0) coverage of G353 (centered at α,δ (J2000) = 17:30:26.28, −34:41:49.7) and its parent filament (dark lane traced by ATLASGAL 870 µm emission; Schuller et al. 2009) with light green and gray contours, respectively.

Motte et al. (2022) classify this protocluster as being in an intermediate evolutionary state, located at ~2 kpc, and hosting a total mass of 2.5 × 10<sup>3</sup> M<sub>☉</sub>. They describe G353 as isolated, without obvious interaction with massive nearby stellar clusters. Using moment maps derived from the N<sub>2</sub>H<sup>+</sup> (1–0) 12 m dataset, they suggest the presence of multiple velocity components, indicating a complex velocity field. They propose that G353 is composed of filaments interacting at the central hub. As is presented in Bonfand et al. (2024), this region is an outlier in the ALMA-IMF hot core sample. Only one weak, low-mass (<2 M<sub>☉</sub>) compact methyl formate source is detected and it lacks strong emission from complex organic molecules. They state that this protocluster is in a chemically poor stage, in which further characterization of this region is required.

The N<sub>2</sub>H<sup>+</sup> (1–0) transition (ν = 93.173809 GHz), given its high critical density, n<sub>crit</sub> = 2 × 10<sup>5</sup> cm<sup>-3</sup> (Ungerechts et al. 1997), allows us to access the dense gas kinematics present in the innermost parts of star-forming regions (Caselli et al. 2002a; Bergin et al. 2002; Tafalla et al. 2004; Lippok et al. 2013; Storm et al. 2014; Hacar et al. 2018; Chen et al. 2019; González Lobos & Stutz 2019; Álvarez-Gutiérrez et al. 2021). The J = 1→0 transition presents seven hyperfine components (Cazzoli et al. 1985; Caselli et al. 1995, 2002a). The kinematic analysis of this complex emission can be simplified by considering only the well-separated isolated component (93.17631 GHz; F<sub>1</sub>, F = 0, 1 → 1, 2; Cazzoli et al. 1985). Such simplification is convenient to study the complex velocity fields found at the center of filaments. These regions present the densest environments for star formation, usually presenting multiple, blended velocity components, in which the velocity distributions exhibit twists, turns, spirals, and wave-like patterns (Csengeri et al. 2011; Fernández-López et al. 2014; Stutz & Gould 2016; Liu et al. 2019; González Lobos & Stutz 2019; Henshaw et al. 2020; Álvarez-Gutiérrez et al. 2021; Sanhueza et al. 2021; Redaelli et al. 2022; Olguin et al. 2023). Recent techniques, such as the intensity-weighted position-velocity (PV) diagrams (González Lobos & Stutz 2019; Álvarez-Gutiérrez et al. 2021), allow us to characterize processes such as infall, outflow, or rotation present in these environments, in which high spatial and spectral resolution studies open a window onto the small-scale gas kinematics of star-forming regions. In addition to the PV diagrams, we can create position-position-velocity (PPV) diagrams, in order to identify coherent structures that might be both spatially and kinematically associated (Chen et al. 2019; Henshaw et al. 2019, 2020; Sanhueza et al. 2021; Redaelli et al. 2022).

In this paper we investigate the N<sub>2</sub>H<sup>+</sup> dense gas kinematics of G353 from large (protocluster) to small (cores) scales. In Sect. 2, we present the data. In Sect. 3 we introduce our N<sub>2</sub>H<sup>+</sup> isolated extraction procedure. In Sect. 4, we model and decompose the multiple velocity components found in the N<sub>2</sub>H<sup>+</sup> isolated component spectra. In Sect. 5, we show our gas kinematic analysis, from protocluster to core scales. In Sect. 6, we show that G353 might be under gravitational collapse at small and large scales. In Sect. 7, we estimate mass accretion rates for multiple velocity gradients (VGs) characterized in our N<sub>2</sub>H<sup>+</sup> data. We discuss our results in Sect. 8, and we present our summary and conclusions in Sect. 9.

## 2. Data

### 2.1. ALMA-IMF data

We make use of the N<sub>2</sub>H<sup>+</sup> (1–0) 12 m, 7 m, and Total Power observations described in Motte et al. (2022) for our analysis,
<table><caption>**Table 1.** Relevant parameters of the N₂H⁺ 7m+12m imaging.</caption><thead><tr><th>Field size</th><th>Pixel scale</th><th>Beam size</th><th>BPA (°)</th><th><sup>(a)</sup>RMS (K)</th><th>Channel width (km s<sup>-1</sup>)</th><th>RMS velocity range (km s<sup>-1</sup>)</th><th><sup>(b)</sup>V<sub>LSR</sub> (km s<sup>-1</sup>)</th></tr></thead><tbody><tr><td>176" × 172"</td><td>0.72"</td><td>1.96" × 2.29"</td><td>80.19</td><td>0.37</td><td>0.23</td><td>[-43 ; -32], [0 ; +7]</td><td>-17</td></tr><tr><td>1.72 pc × 1.67 pc</td><td>1.44 kau</td><td>~4 kau × 4.6 kau</td><td></td><td></td><td></td><td></td><td></td></tr></tbody></table>

Notes. (a) RMS value at the peak of the RMS distribution. (b) Obtained from Motte et al. (2022).

providing robust uv plane coverage. We image the combination of the N₂H⁺ 7 m and 12 m (from now on called “7m+12m”) measurement set of G353, using the publicly available imaging scripts from the ALMA-IMF github repository<fnref n="2" />. These data are corrected by the primary beam response pattern. Due to the missing large-scale emission, we find that near the V<sub>LSR</sub> of the protocluster (−17 km s<sup>-1</sup>; Wienen et al. 2015; Motte et al. 2022) some subregions in the 7m+12m cube present deep negative artifacts (“negative bowls”). To cover all possible uv scales, we combine the N₂H⁺ 7m+12m continuum-subtracted cube with the Total Power observations from the ALMA-IMF LP. We use the feather<fnref n="3" /> task from CASA 5.6.0. With this combination, we were able to recover the missing flux, seen as negative bowls, present in the interferometric-only data. We produce a fully combined, multi-scale, feathered dataset which we use for our dense gas kinematic analysis.

To estimate and subtract the continuum emission present in the 7 m+12 m cube, we use the imcontsub<fnref n="4" /> CASA task. We select the emission-free channels between −43 km s<sup>-1</sup> and −33 km s<sup>-1</sup>, and set the polynomial degree of the continuum fit (fitorder) to zero. We list relevant final image parameters in Table 1, such as the field size, pixel scale, beam size, root-mean-square noise (RMS), and channel width.

We use 12 m datacubes from Cunningham et al. (2023) to compare the shock tracers SiO (5–4) and CO (2–1) to our N₂H⁺ kinematic analysis. We use DCN and N₂H⁺ data to determine core velocities (Sect. 4.2). To determine total masses in specific regions we use the N(H₂) map from Díaz-González et al. (2023).

## 2.2. Core properties from published catalogs

We use the cores catalog<fnref n="5" /> from Louvet et al. (2024). These cores were identified using the getsf algorithm, specialized in source extraction on regions with complex filamentary structures (Men'shchikov 2021). This procedure was done using the 1.3 mm continuum maps, smoothed at a common resolution of ~2700 au, obtaining a total of 45 sources for G353. We also use the DCN core velocities (15 sources, Cunningham et al. 2023) and the SiO outflow catalog (16 sources, Towner et al. 2024) in order to look for correlation between the N₂H⁺ gas kinematics and cores/outflows position and properties. It is worth mentioning that, within a radius of 0.3 pc from the center of G353 (Motte et al. 2022), we find 60% of the 1.3 mm cores (27 sources), and ~70% of the cores with DCN velocities and SiO outflows (11 sources from each catalog). Of these 11 outflows, 7 are “red” 3 are “blue” (monopolar), and 1 is “bipolar” (Towner et al. 2024). The presence of these sources might imply a complex velocity

field in this region, given that cores and outflows disturb the kinematics of the surrounding gas.

## 3. N₂H⁺ isolated component analysis and filament extraction

The N₂H⁺ (1–0) transition is characterized by its hyperfine emission composed by seven components (Caselli et al. 1995, see their Fig. 1). We present an ideal example of N₂H⁺ emission in Fig. 2, panel d. In this work, we refer to the triplet of hyperfine components that present the highest intensities as the main N₂H⁺ components, located at the center of the line emission at v<sub>rest</sub> = 93.173806 GHz. We refer to the most blueshifted hyperfine component as the isolated component, at v<sub>rest</sub> = 93.17631 GHz, shifted by ~8 km s<sup>-1</sup> relative to the main N₂H⁺ component (see Table 1 from Cazzoli et al. 1985).

### 3.1. Isolated component identification and extraction

We developed an algorithm to extract only the isolated hyperfine component from every pixel in the feathered datacube. This is in order to reduce the complexity of our data, given that it may contain multiple velocity components in addition to the hyperfine line emission. Considering that the N₂H⁺ emission moves in velocity across the protocluster, our approach is to find the velocity where the emission of the isolated component ends and remove the rest of the line emission. We also preserve the emission-free channels, at low (−43 km s<sup>-1</sup> to −31.5 km s<sup>-1</sup>) and high (0.7 km s<sup>-1</sup> to 6.7 km s<sup>-1</sup>) velocities, to improve future RMS estimations if needed. It should be noted that in the procedures described below, we use find Peaks<fnref n="6" /> to detect peaks and valleys in the different spectra.

Our extraction approach is separated into two procedures, for low and for high signal-to-noise ratio (S/N) data (see Fig. 3, and text below). In order to determine which data have low or high S/N, we obtained the mean spectrum over all the spatial pixels of the cube, which serves as a guide to determine the velocity at the “mean dip” (V<sub>mean dip</sub> = −22 km s<sup>-1</sup>, dashed gray line in Fig. 2 panel “a”). This velocity represents the mean location of the intensity valley between the isolated and the main components of the N₂H⁺ emission. We defined ΔV<sub>mean</sub> = 3.2 km s<sup>-1</sup> as the difference between V<sub>mean dip</sub> and the velocity at the peak of the mean isolated component V<sub>mean peak</sub> (dashed black line in Fig. 2 panel “a”), used in our velocity guesses for the high S/N extraction procedure (see below).

To create a S/N map of the isolated component, we first measured the RMS noise in emission-free channels (−43 km s<sup>-1</sup> to −31.5 km s<sup>-1</sup>), and the peak intensity in the channels range in which the mean isolated component is located (−43 km s<sup>-1</sup> to V<sub>mean dip</sub>). This approach allows us to exclude the

<footnote n="6">⁶ https://docs.scipy.org/doc/scipy/reference/generated/scipy.signal.find_peak.html</footnote>

<footnote n="2">² https://github.com/ALMA-IMF/reduction</footnote>

<footnote n="3">³ https://casadocs.readthedocs.io/en/stable/api/tt/casatasks.imaging.feather.html</footnote>

<footnote n="4">⁴ https://casadocs.readthedocs.io/en/stable/api/tt/casatasks_analysis.imcontsub.html</footnote>

<footnote n="5">⁵ Available at www.almaimf.com</footnote>
<figure><img src="image_3.png" /></figure>

**Fig. 2.** N<sub>2</sub>H<sup>+</sup> spectra involved in our isolated component extraction procedure. Panel a: normalized average N<sub>2</sub>H<sup>+</sup> spectrum (solid black line) over the entire region. We show the location of the mean peak of the isolated component (dashed black line) and the “mean dip” (dashed gray line), see text. Middle and right panels: normalized example spectra (within a pixel) of the N<sub>2</sub>H<sup>+</sup> isolated velocity component extraction procedure (see Sect. 3). We show G353 N<sub>2</sub>H<sup>+</sup> spectra with solid black lines and the extracted isolated component, along with emission-free channels, with dashed red lines. Panel d: expected N<sub>2</sub>H<sup>+</sup> emission for an excitation temperature of 15 K, an opacity of 1, velocity centroid of −17 km s<sup>−1</sup>, and a line width of 0.2 km s<sup>−1</sup>. We see the seven hyperfine components characteristic of this tracer, in which the most blueshifted corresponds to the isolated component. To derive this emission we use “n2hp_vtau” model from PySpecKit. In panels b, c, e, and f) we indicate V<sub>mean dip</sub> with a dashed gray line. We present data with S/N < 5 in panel b, in which we make a rough extraction based on the V<sub>mean dip</sub>. We show data with S/N ≥ 5 in panels c, e, and f, presenting clear single, double, and triple N<sub>2</sub>H<sup>+</sup> isolated velocity components, respectively. In these examples we represent the selected velocity guess that separates the isolated component emission from the main line emission with dashed blue lines. The offset positions (Δl, Δb) of the spectra in panels b, c, e, and f are (0.68 pc, 0.27 pc), (0.12 pc, −0.28 pc), (−0.08 pc, 0.47 pc), and (0.16 pc, −0.07 pc), respectively. These offsets are estimated relative to the center of the region (see Sect. 1).

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

**Fig. 3.** G353 N<sub>2</sub>H<sup>+</sup> isolated component S/N map. The white contour indicates the location of data with an isolated component *S/N* ≥ 5. We show the location of the 1.3 mm cores presented in Louvet et al. (2024) with black ellipses and are located in regions with *S/N* ≥ 15. We indicate the beam size of these data with a gray ellipse at the bottom left corner. Outside the S/N contour we make a rough extraction of the isolated component (see text). For data inside the S/N contours, we implement a procedure based on detection of peaks and valleys, to individually extract high (≥5) S/N isolated components (see Sect. 3).

emission of the main line components. We encountered spurious emission at the edges of the S/N map. We adopted the procedure from Towner et al. (2024) by using the image processing techniques implemented by binary_erosion<sup>7</sup> (1 iteration) and binaryPropagation<sup>8</sup> to clean the data for further analysis. binary_erosion allows us to remove the spurious emission in the outskirt of the map, although this approach also removes high S/N edges of our protocluster. Then, we used binaryPropagation on the cleaned S/N map, using the original S/N map mask, to restore only the protocluster edges. To test our cleaning approach we computed the total integrated intensity using the Python package SpectralCube<sup>9</sup> in the range of -31.5 km s<sup>-1</sup> to V<sub>mean dip</sub> using the original and cleaned S/N mask. We estimated that the removed spurious emission accounts for ~2% of the total integrated intensity for data with S/N > 5.

In Fig. 3, we show the N<sub>2</sub>H<sup>+</sup> isolated component S/N map, where at S/N values ≥5 we capture the cloud emission while excluding noise (white contour). We set our isolated component S/N threshold to 5, in order to use one of the two extraction procedures (see below). In this section we refer to high (low) S/N spectrum if its isolated component S/N ≥ 5 (< 5). For low S/N

<footnote n="7">⁷ https://docs.scipy.org/doc/scipy/reference/generated/scipy.ndimage-binary_erosion.html</footnote>

<footnote n="8">⁸ https://docs.scipy.org/doc/scipy/reference/generated/scipy.ndimage binaryPropagation.html</footnote>

<footnote n="9">⁹ https://github.com/radio-astro-tools/spectral-cube</footnote>
spectra, we extracted all the channels in the velocity range from −43 km s⁻¹ up until V<sub>mean dip</sub> (panel “b” in Fig. 2). We take this approach given that for low S/N data we cannot identify peaks in the N₂H⁺ emission in a reliable manner. For high S/N spectra the extraction procedure consists of creating different velocity guesses that represent the location of the intensity valley, similar to the definition of V<sub>mean dip</sub> (Fig. 2). Then, we selected a velocity guess based on its associated weight (see description below). This approach is described in detail here:

- First, we implemented a rolling average along each spectrum. This is in order to smooth over intensity bumps that might result in false positives for the detection of peaks and valleys. For this procedure, we averaged considering two channels before and after each velocity.

- After smoothing, for each spectrum we identified the isolated component peak using find_peaks. We called the velocity associated with this peak V<sub>isolated component</sub>. In the case of multiple velocity components it represents the most blueshifted one. We found the intensity valley between the isolated component and the N₂H⁺ main line emission by inverting the spectrum and finding the first peak which is the inverted intensity valley. We defined the associated velocity to this intensity valley as V<sub>first minima</sub>.

- We created three velocity guesses based on the properties of each spectrum in our cube (see points below). These are the 1st guess: V<sub>isolated component</sub> + ΔV<sub>mean</sub>. 2nd guess: V<sub>first minima</sub>, in the case of multiple isolated components this guess might incorrectly capture the intensity valley after the first isolated component. In that case, the other guesses are needed for a reliable isolated component extraction. 3rd guess: V<sub>mean dip</sub> + ΔV<sub>mean</sub> to provide a velocity cut further away from the V<sub>mean dip</sub>. This guess is mainly useful in the case where multiple isolated components cover a velocity range larger than the one probed by the other two guesses.

- From each velocity guess we estimated two parameters to later decide which one to use. One is the absolute value of its associated intensity “I<sub>i</sub>” (i.e., intensity at the guess velocity), and the other is distance in velocity “dV<sub>i</sub>” to the mean dip. The i subscript represents the guess associated with these parameters. We save the parameters of each guess in the lists “I” and “dV”.

- We normalized these lists by their minimum value ensuring that the guess with the smallest “I<sub>i</sub>” and “dV<sub>i</sub>” will have a weight (w) of 1, defined in Eq. (1). We do not encounter divergences in this normalization given these parameters are not exactly zero.

$$w = (I_{\text{norm}} \times 0.2 + dV_{\text{norm}} \times 0.8)^{-1}, \quad (1)$$

in which the “norm” subscript indicates that the parameter list is normalized by dividing it by its minimum value.

- By visual inspection we considered that we obtain good extraction results when the weight is mostly dependent on dV and in a minor part on I. This is reflected by the 0.2 and 0.8 factors multiplying I<sub>norm</sub> and dV<sub>norm</sub>, respectively, in the definition of “w” in Eq. (1).

- We chose the guess with the weight closest to unity.

- Similarly as for S/N < 5, we extracted the spectra from −43 km s⁻¹ up until the velocity of the chosen guess, and preserved the emission-free channels from 0.7 km s⁻¹ to 6.7 km s⁻¹.

Various examples of N₂H⁺ spectra and isolated hyperfine component extraction are shown in Fig. 2, in which we can see

<figure><img src="image_5.png" /><figcaption>Fig. 4. Moment 0 map of the extracted N₂H⁺ isolated component emission. We use the FilFinder Python package to identify the main filamentary structure present in G353 (see Appendix A). We identify three filaments (F1, F2, and F3; green lines) converging toward the central hub. The location of most of the 1.3 mm cores (red ellipses), projected in the POS, lie on top of the spine of these filaments, especially in the hub.</figcaption></figure>

spectra containing one (panel “c”), two (panel “e”), and three (panel “f”) velocity components, all well extracted by our procedure. In Stutz et al. (in prep) this approach is generalized to all ALMA-IMF regions for N₂H⁺, providing reliable results.

## 3.2. Filamentary identification

We use the FilFinder Python package (Koch & Rosolowsky 2015) in order to detect the most prominent filaments in this region (green lines in Fig. 4). The procedure, including the parameters we used for the filamentary identification, is presented in Appendix A. In Fig. 4, we indicate the detected filaments with green lines, on top of the moment zero map of the multiple N₂H⁺ isolated components. We see that G353 is a hub-filament system (HFS), composed by three main filaments converging at its center. The HFSs are a characteristic feature of early stages of star formation, in which gas flows through the filaments toward the central hub triggering star formation (Myers 2009; Galván-Madrid et al. 2010; Busquet et al. 2013; Galván-Madrid et al. 2013; Peretto et al. 2014; Kumar et al. 2022; Zhou et al. 2022). We see that in the plane of the sky (POS) most of the 1.3 mm cores (red ellipses) are located on top of the filaments. This spatial agreement between filaments and protostellar cores is consistent with filamentary fragmentation (André et al. 2010; Busquet et al. 2013; Stutz & Kainulainen 2015; Kuznetsova et al. 2015, 2018).

# 4. N₂H⁺ isolated component velocity decomposition

In Fig. 2, we see that clear multiple isolated velocity components are present in our dataset. In this section we present our approach
<figure><img src="image_6.png" /><figcaption>**Fig. 5.** Spatial distribution of the modeled N<sub>2</sub>H<sup>+</sup> isolated velocity components. We indicate the main filament structure with green lines (see Fig. 4). In blue, green, and red we indicate the first, second, and third velocity components, respectively. We indicate the beam size of these data with a gray ellipse at the bottom left corner. The emission of the first and second components is more extended and intense than the third, most red-shifted velocity component. The 1.3 mm cores match regions with high integrated intensity, mostly traced by the first and second velocity components.</figcaption></figure>

regarding the line modeling of these multiple velocity compo-
nents. We also show how we use these models to determine core
velocities.

## 4.1. Modeling of the isolated component emission

To characterize the complex dense-gas kinematics traced by N<sub>2</sub>H<sup>+</sup> we followed the method in Álvarez-Gutiérrez et al. (2021), and we used the spectroscopic toolkit PySpecKit (Ginsburg & Mirocha 2011; Ginsburg et al. 2022) to model and decompose the isolated component emission. PySpecKit adjusts a fixed number of components set by the user, based on visual inspection of the data we imposed three velocity components to every spectrum and then remove false positives (see below). Given the kinematic complexity of the data and cursory inspection of the spectra, a simpler analysis with only two components contradicts the data. In essence, three components is the simplest possible choice, given the data. While this might fail for a small number of spectra that could require ≥ 4 velocity components, the residuals indicate that this could occur in a severe minority of cases, and hence more components is not warranted given the S/N and resolution of this particular dataset. To improve the convergence of PySpecKit, we created a set of ranges for the parameters that define each of the three Gaussian velocity components, namely the peak intensity, central velocity, and velocity dispersion. After testing different parameter ranges, we set the intensity range between 1.76 K (four times the mean RMS) and 30 K, the velocity centroid from -30 km s<sup>-1</sup> to -20 km s<sup>-1</sup>, and the velocity dispersion from 0.22 km s<sup>-1</sup> to 1 km s<sup>-1</sup>.

From the results using the ranges defined above, we notice
that some modeled components do not fit any emission. These
fits are the result of imposing to the fitter a fixed number of com-
ponents, given these spectra can be better represented by one
or two velocity components. In these fits there is no uncertainty
estimation for both the peak intensity and velocity dispersion.
Based on these two criteria we removed those velocity fits from
the modeled cube. With this cleaning approach we are left with
spectra characterized by one (~34%), two (~53%), and three
(~13%) Gaussian velocity components. We present the Gaussian
fits of the high S/N spectra from Fig. 2 in Appendix B.

In Fig. 5, we show the spatial distribution of the multiple Gaussian velocity components. In gray we indicate the main filamentary structure in the region (see Sect. 3.2). The first and second velocity components, in blue and green, respectively, present most of the high intensity emission and they also spatially dominate over the third, most red-shifted component. Both the first and second components trace mostly the filaments F1 and F3 from Fig. 4, where most of the 1.3 mm cores are located. The position of these cores coincide with high integrated intensity regions in these isolated velocity components. The most redshifted component is compact and less intense compared to the first and second velocity components. This velocity distribution is located mostly along the filament F2 and the central hub (see Fig. 4). In Fig. 6, we present the number of Gaussian velocity components for each spectrum, where we highlight that:

- Most of the N<sub>2</sub>H<sup>+</sup> data presents emission characterized by two velocity components.

- Most of the spectra described by three velocity components are located in the innermost parts of the region.

- Most of the cores (black ellipses; Louvet et al. 2024) are located in regions with spectra presenting two to three Gaussian velocity components, indicating kinematic complexity even at ~4 kau (N<sub>2</sub>H<sup>+</sup> spatial resolution).

- Single velocity component spectra are located preferentially in the outskirt of the protocluster.

In Fig. 7, we show the histogram of the fitted velocity centroid of each Gaussian velocity component. The peaks of these distributions are located at -27, -24.7, and -23.3 km s⁻¹, respectively, well-separated in velocity. From hereafter we refer to these distribution as blue, green, and red, respectively. Most of the velocity components appear to be associated with the blue and green distributions. For consistency with the different tracers used in further analysis, we shift the isolated component velocities by +8 km s⁻¹, to the reference frame of the main line components of N₂H⁺ (Cazzoli et al. 1985).

## 4.2. DCN and N<sub>2</sub>H<sup>+</sup> derived core velocities

In this section, our goal is to increase the sample of core
velocities from the already published DCN catalog, aiming to
<figure><img src="image_7.png" /><figcaption>**Fig. 6.** Spatial distribution of the modeled N<sub>2</sub>H<sup>+</sup> isolated component spectra. These models are composed by up to three Gaussian velocity components. The 1.3 mm cores and beam size are the same as in Fig. 3. Most of the 1.3 mm cores are located in regions with spectra presenting two to three velocity components (see Sect. 4.1).</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>**Fig. 7.** Normalized velocity centroid distributions of each N<sub>2</sub>H<sup>+</sup> Gaussian velocity component. The velocities at the peak of the distributions are –26.9 km s<sup>-1</sup>, –24.7 km s<sup>-1</sup>, and –23.3 km s<sup>-1</sup> for component 1 (blue), 2 (green), and 3 (red), respectively (see Sect. 4.1).</figcaption></figure>

explore all the potential in these types of dataset. Given the relatively high $n_{crit}$ of DCN (3–2) (~10<sup>7</sup> cm<sup>-3</sup>) compared to N<sub>2</sub>H<sup>+</sup> (1–0) (2 × 10<sup>5</sup> cm<sup>-3</sup>), DCN (3–2) is known to coincide well with continuum peaks associated with cores (Liu et al. 2015; Cunningham et al. 2016; Minh et al. 2018), while N<sub>2</sub>H<sup>+</sup> is characterized by tracing the dense gas at the innermost parts of star-forming regions (Fernández-López et al. 2014; Hacar et al. 2018; González Lobos & Stutz 2019).

In Cunningham et al. (2023), they use ALMA-IMF 12 m observations of DCN (3–2) to study cores kinematics. They apply line emission fits for the DCN spectra inside the 1.3 mm cores from Louvet et al. (2024). For this procedure they determine core velocities in all ALMA-IMF targets. They classify as DCN single and complex core velocities, spectra that can be fitted with one or multiple Gaussian velocity components, respectively. Due to a global conservative S/N threshold the DCN fitting process missed the velocity estimation of some cores. For G353 only 15 out of the 45 cores present DCN velocity fits.

We used the ALMA-IMF DCN 12 m data from Cunningham et al. (2023), which presents a velocity resolution of ~0.34 km s<sup>-1</sup>. For each DCN core velocity described by a single component (Cunningham et al. 2023) we compared the emission of the DCN and modeled N<sub>2</sub>H<sup>+</sup> isolated spectra. We found an average velocity offset between the DCN peak and the closest N<sub>2</sub>H<sup>+</sup> isolated component peak of ~0.65 km s<sup>-1</sup>, less than two DCN channel widths. We used this approach for the remaining 30 cores, in order to determine their DCN velocities.

Here we estimated the RMS of the DCN data in emission-free channels in the range of –42 km s<sup>-1</sup> to –25 km s<sup>-1</sup> and we obtain the S/N map by dividing the peak intensity by the RMS. For the procedure below we only used DCN spectra with S/N > 3. Next, we extracted the average DCN and modeled N<sub>2</sub>H<sup>+</sup> isolated component spectrum of these 30 cores. We identified the N<sub>2</sub>H<sup>+</sup> isolated velocity component closest to the DCN peak within three DCN channel widths. We found that 11 out of these 30 cores present DCN with S/N > 3 close to one N<sub>2</sub>H<sup>+</sup> velocity component. Here, we defined the velocity of these cores as the velocity where the DCN emission peaks. On average, these cores have a velocity offset between these two tracers less than 0.8 km s<sup>-1</sup> (<2.5 DCN channels), similar to the results obtained for the 15 cores with DCN velocities from Cunningham et al. (2023), and they present an average velocity offset of 0.38 km s<sup>-1</sup>. Throughout this paper we refer to these cores as “DCN & N<sub>2</sub>H<sup>+</sup> cores” given they are derived from the comparison of these two tracers. In Appendix C we present two examples of the comparison of the DCN and N<sub>2</sub>H<sup>+</sup> spectra in cores in which we see clear agreement between these tracers (Fig. C.1). In Table C.1 we show these obtained core velocities from our comparison between DCN & N<sub>2</sub>H<sup>+</sup> spectra, complementing the DCN catalog from Cunningham et al. (2023).

# 5. Analysis of position-velocity diagrams

## 5.1. Traditional position-velocity diagram

We start by analyzing the “traditional” PV diagram shown in Fig. 8. We created this diagram by taking the total intensity along the Galactic longitude, in which $\Delta b$ indicates the distance in parsec relative to the center of G353, assuming a distance to the protocluster of 2 kpc (Motte et al. 2022). We see general agreement between the DCN core velocities and the N<sub>2</sub>H<sup>+</sup> velocity distribution. This suggests that most of the cores are still kinematically coupled to the dense gas in which they formed. As is presented in Sect. 4.2, the DCN and N<sub>2</sub>H<sup>+</sup> velocities match within 0.8 km s<sup>-1</sup> (<2.5 DCN channels).

Regarding the dense gas velocity distribution, in Fig. 8 we see a velocity spread of ~8 km s<sup>-1</sup> in the subregion between $\Delta b$ ~ -0.3 pc to 0.1 pc. Most of the intensity on this diagram is located at the upper part of this subregion, at $\Delta b$ ± 0.1 pc. This spread is also present in the PV diagram along $\Delta b$ and $\Delta l$ shown
<figure><img src="image_9.png" /></figure>

**Fig. 8.** “Traditional” PV diagram of the N₂H⁺ modeled isolated components, created by collapsing the *l* coordinate. Δb indicates the distance along *b* in pc, relative to the center of G353, assuming a distance of 2 kpc (Motte et al. 2022). The colormap indicates the total intensity along *l*. With fuchsia and black crosses we show the 1.3 mm cores with single and complex DCN velocities detections, respectively (Cunningham et al. 2023). With dark cyan “×” markers we show the 1.3 mm cores with velocities derived from DCN and N₂H⁺ data (see Sect. 4.2). The size of the markers indicate relative mass (Louvet et al. 2024). Black contours indicate total intensities at 40, 160, 280, and 400 K. We see a large-scale velocity spread (ΔV ~8 km s⁻¹) around Δb ~ -0.3 pc – 0.13 pc (see also Sect. 5.2). We show the major axis of the beam and the channel width with a white rectangle at the bottom left corner.

in the top right and bottom left panel of Fig. 9. We explore the
possible origin of this structure in Sect. 6.

5.2. *Intensity-weighted PV diagrams*

In the top left panel of Fig. 9 we show the spatial distribution of the fitted Gaussian velocity components (see Sect. 4.1). The blue, green, and red color maps indicate the integrated intensity of the first, second, and third velocity components of the N₂H⁺ spectra, respectively. We note that the spatial overlap between any of these components is presented in Fig. 6.

As is seen in Fig. 8, the traditional PV diagram provides
information on the dynamics on the large, protocluster-scale,
environment. Meanwhile, the intensity-weighted PV diagram
(Fig. 9), in which the color of each point indicates its integrated
intensity, highlights the small core-scale kinematics. Similarly
as in González Lobos & Stutz (2019) and Álvarez-Gutiérrez
et al. (2021), from the isolated component line decomposition
(Sect. 4.1), we derived the integrated intensity and velocity
centroid for each Gaussian velocity component. Using these
parameters we created intensity-weighted PV diagrams along the
*b* and *l* coordinates. We present these N₂H⁺ PV diagrams in the

bottom left and top right panels of Fig. 9. The key features on the
position-position (PP) and on the top right PV diagram, are:

- The agreement between the DCN core velocities and the overall N₂H⁺ PV structures suggests that cores are still kinematically coupled to the dense gas in which they formed.

- We see at least nine clear and prominent V-shaped velocity gradients (see Appendix D), across all velocity components. The orientation of these V-shape, pointing to the left/right (top right panel) or up/down (bottom left panel), follow no clear preference.

- In some cases, the vertex of these V-shapes is close spatially and in velocity to the location of cores.

- In the plane of the sky (POS), all three velocity components overlap in most of the region.

- This technique recovers the large-scale velocity spread present in Fig. 8 and highlights small-scale structures.

- The most prominent V-shape is located at (Δb, V) = (-0.14 pc, -20.5 km s⁻¹), between two 1.3 mm cores with DCN detections (see Sect. 6).

For a better visualization of the 3D structure of these V-shapes
we provide an interactive 3D PPV diagram [online](#).

5.3. *Velocity gradients*

In this section, we focus on the most prominent blue V-
shape (Fig. 9, top right panel). In Fig. 10, we show this
velocity distribution, named V-shape “C”, in detail. Given
the visual linearity of the VGs composing the V-shape, we
applied a linear fit to these distributions in order to charac-
terize them. For these fits we considered data only above an
integrated intensity threshold of 8 K km s⁻¹ and 3 K km s⁻¹,
for the upper and lower gradient, respectively. We removed
data not related with the VG, clustered in the ranges of
(Δb, V) ~ (-0.025 – 0.04 pc, –19.5 – –18.5 km s⁻¹), which
lie just outside the filament hosting this V-shape on the POS.
Additionally, we weighted each point based on their integrated
intensity to make our fits more robust. The slopes of the linear
fits represent the VGs in km s⁻¹ pc⁻¹. These linear fits follow the
VGs distribution and these are somewhat asymmetric, the upper
gradient is slightly shallower than the bottom gradient. Given
the unknown inclination angle (θ) of these structures relative
to the POS, the observed VG is just a fraction of the original
VG. These are related as VG = VG<sub>original</sub> · sin(θ). These VGs
present values between ~13 to ~18 km s⁻¹ pc⁻¹ (see Fig. 10).
Additionally, we estimated the center of this V-shape as the
velocity-weighted mean position of the points composing this
structure. With this approach the position of the points closest
to the V-shape apex present more weight, obtaining the center
of this V-shape at (l, b) = (353.4135°, –0.3657°). This position
is located between “core 2” and “core 3”, both of them having
DCN velocity fits (Cunningham et al. 2023). These cores present
masses of 20.7 and 9.4 M<sub>☉</sub>, respectively. We inspected the core
catalog from Louvet et al. (2024) and concluded that the location
of this V-shape do not coincide with any core.

In the left panel of Fig. 11, we show the integrated intensity
of the multiple modeled N₂H⁺ isolated components. With col-
ored boxes we show the areas where we create the different PV
diagrams presented on the right panel. These boxes are centered
at the main blue V-shape, matching the area of this V-shaped
structure (see Fig. 12). We show that the overall structure in PV
space is conserved at different angles, excluding the possibility
of this velocity feature being the result of projection effects.
<figure><img src="image_10.png" /><figcaption>**Fig. 9.** N₂H⁺ intensity-weighted position-position and PV diagrams of G353. Top left: spatial distribution of the fitted N₂H⁺ Gaussian isolated velocity components (blue, green, and red, see Sect. 4.1). Ellipses indicate the location of the 1.3 mm continuum cores (Louvet et al. 2024). Orange indicates cores with no DCN detections. Fuchsia and black represent cores with single and complex DCN velocities (Cunningham et al. 2023). DCN & N₂H⁺ cores are indicated with cyan. We show the beam size with a gray ellipse in the bottom left corner. Top right and bottom left: intensity-weighted PV diagrams along the *b* and *l* coordinates, respectively. For the 1.3 mm core velocities we use the same colors and markers convention from Fig. 8. For reference we indicate with a red arrow, in both the top right and bottom left panels, a VG of 10 km s⁻¹ pc⁻¹ corresponding to a timescale of ~0.1 Myr. We see multiple V-shapes near the location of cores across all velocities in the PV diagrams, more prominently in the top right panel (see Fig. D.1). The most prominent V-shape is located in the blue component, at (V, Δb) ~ (−20.5km s⁻¹, −0.14 pc) (see Sect. 5.2). We provide an interactive 3D PPV diagram online.</figcaption></figure>

For this V-shape we used its composing VGs to derive timescales as $t_{VG} = 1/\text{VG}$, similar to the procedure for a rotating filament presented in Álvarez-Gutiérrez et al. (2021), and we suggest these can be interpreted as inflow timescales. The $t_{VG}$ values for this V-shape range between ~50–70 kyr. These timescales are short compared to the ~0.21 Myr free fall time ($t_{ff}$) of the protocluster (Motte et al. 2022), and a few times larger than the $t_{ff}$ of nearby cores (~20 kyr, within 0.1 pc of this V-shape). To determine the cores $t_{ff}$, we use the 1.3 mm core masses from Louvet et al. (2024).

We characterized eight more N₂H⁺ V-shaped structures. In Fig. 13, we show the position of these V-shapes in the POS. Previous studies have detected VGs along filaments, toward a converging point (Peretto et al. 2014; Pan et al. 2024; Rawat et al. 2024). In the case of G353 we see that instead of detecting a single V-shape at the intersection of its filaments, the hub

appears to be fragmented into multiple, small-scale V-shaped VGs. Only V-shapes “A”, “D”, and “E” are outside of the hub, with V-shape “D” located on top of filament “F3”. This indicates that the V-shaped structures are not exclusive to the central parts of HFS, but are also present in comparatively isolated regions. We note that within a ~beam size from the apex of V-shape “B” (see Fig. D.2), the continuum core “7” (~6 M☉) is located.

In Henshaw et al. (2014), they propose two scenarios that might produce these V-shaped VGs (see their Fig. 12). One scenario suggests that gas in a filament is flowing toward a denser region (infall), while the other scenario suggests that a protostellar outflow moves the dense gas located in its vicinity. To analyze the different dynamical processes present in this region, we used the ALMA-IMF 12 m data of the shock, outflow tracer SiO (5–4), from Cunningham et al. (2023). From those data we created its intensity-weighted PV diagram, presented in Fig. E.1
<figure><img src="image_11.png" /><figcaption>Fig. 10. Zoomed-in version of the top right panel of Fig. 9, centered at the prominent blue V-shape (“C”) located at ($\Delta b, V$) = (−0.14 pc, −20.5 km s⁻¹). We indicate the major axis of the cores with vertical lines in fuchsia (DCN single) and cyan (N₂H⁺), similarly we represent the major axis of the beam with a vertical orange line. We apply linear fits to the upper and lower distributions, represented by darker points. These points are selected based on an integrated intensity threshold (see Sect. 5.3). We weight each point by their integrated intensity and derive VGs from the slope of these linear fits. The range of the obtained VGs is ~13–18 km s⁻¹ pc⁻¹. We defined the timescale associated with the VG as $t_{VG} = VG^{-1}$, being in the range of ~50–70 kyr. We show eight more well characterized V-shapes in Appendix D.</figcaption></figure>

(see Appendix E for more details). We note that there is almost no SiO emission nor outflow sources at the location of the N₂H⁺ V-shape (see Fig. 12). The ΔV within 10'' from this main blue V-shape is ~40 km s⁻¹, while for the whole SiO dataset is ~80km s⁻¹, showing a clear difference in the SiO ΔV and the core velocities. Furthermore, the SiO ΔV is ~10 times larger than the one of N₂H⁺. This difference in traced velocities implies that these two molecules trace vastly different physical phenomena. We suggest that the small velocity range probed by N₂H⁺ indicates that the VGs can be considered as infall signatures (see Sect. 6). We discuss the possible morphology of the filaments hosting V-shapes in Sect. 8.

## 6. G353 as a collapsing region

We used the SiO (5–4) and CO (2–1) data from Towner et al. (2024) and Cunningham et al. (2023) to identify possible outflows near the V-shape presented in Fig. 10. For CO we measured the RMS in the emission-free velocity range of 145–286 km s⁻¹. We used only CO data with $S/N > 3$ for our analysis. The cleaning of the SiO data is described in Appendix E. In Fig. 12, with a fuchsia “X” we indicate the center of the V-shape from Sect. 5.3, located between two 1.3 mm cores 2 and 3 (black ellipses), which present DCN velocity detections. From these diagrams we see there is neither SiO nor CO outflow detection at the location of the main blue V-shape ($\Delta l, \Delta b$) = (0 pc, 0 pc). In the right panel we show the position of the data composing this V-shape, where the velocity peaks toward the center of this velocity feature.

We derived the mass-weighted mean position of the two cores (2 and 3) closest to the V-shape to determine their barycenter (yellow “X” in Fig. 12). We found that the difference between the center of the V-shape and the barycenter of these cores is ~0.3” (~600 au), well below the beam size of the N₂H⁺ data. A

similar offset is also present in the intensity and velocity profiles along filaments from ATOMS data (Zhou et al. 2022, see their Fig. 6). This small spatial offset might suggest that the gas flowing in the V-shape is produced by the gravitational pull toward the barycenter of cores 2 & 3, where the N₂H⁺ radial velocities peak. This interpretation is similar to the one proposed in Zhou et al. (2023) in their kinematic analysis of the G333 complex. They describe the V-shaped VGs as the result of gas funneling from the molecular cloud to clumps which is then funneled onto cores (see their Fig. 9) consistent with gravitational acceleration.

In Fig. 14, we show the mean spectrum of different tracers at the central position of the main blue V-shape. These spectra were measured over a circular region with diameter equal to the major axis of the N₂H⁺ beam (2.28”; ~0.02 pc). This circular area results in a coverage of 1.14 times the N₂H⁺ beam, and ~5.6 times the beams of the H₂CO, DCN, and H₂¹³CO data. We see N₂H⁺ and H₂CO present double component spectrum with asymmetric peaks. Between these peaks we detect DCN and H₂¹³CO emission. The asymmetric spectrum present in N₂H⁺ and H₂CO is consistent with the “blue asymmetry” spectral feature, usually interpreted as infall signature, suggesting that this region is under gravitational collapse (e.g., Anglada et al. 1987; Mardones et al. 1997; Lee et al. 1999, 2001; Smith et al. 2012). Based on the idea that the V-shapes are the result of flowing gas along filaments toward denser regions, the blue-asymmetry detected at the center of the main blue V-shape suggests that gravitational collapse is taking place at the apex of the V-shaped structure.

Regarding large scales, in the traditional PV diagram presented in Sect. 5.1 (see Fig. 8), we see a clear velocity spread around $\Delta b \sim -0.1$ pc, also present in the top right panel of Fig. 9. Below we compare this velocity spread with the velocity distribution produced by infall, in which the gas velocities increase as the distance to the center of infall (“r”) decreases:

$$V_{\text{infall}} = -\sqrt{\frac{2G M}{r}}. \qquad (2)$$

For this comparison, we modeled a sphere with a total mass of 150 M☉, a radius of 0.5 pc, and a power law density profile described by:

$$\rho(r) = \rho_0 \left(\frac{r}{\text{pc}}\right)^{-\gamma}, \quad \gamma = 5.65, \quad \rho_0 = 6.1 \times 10^{-5} \frac{M_{\oplus}}{\text{pc}^3}. \qquad (3)$$

We provide the derivation of $\rho(r)$ in Appendix F. $\gamma = 5.65$ was determined by visual inspection by comparing the obtained radial velocities of the model (see below), at different $\gamma$ values, with the overall shape of the PV distribution.

We then estimated the infall velocity of each point, based on the cumulative mass distribution (“M”) at any given distance to the center (Eq. (2)). We obtain the radial component of the infalling velocities as:

$$V_r = V_{\text{infall}} \times \cos(\arctan(X/Z)), \qquad (4)$$

in which X represents the horizontal coordinate in the POS, while Z represents the (unobserved) depth of the sphere. The spatial coordinates for this model range from −0.5 pc to 0.5 pc.

In Fig. 15, we show the coverage of the PV distribution from our model, described in Eqs. (2)–(4), with a solid white line. We find good agreement between the PV distributions of our approach and the data. The PV distribution of our infall model is consistent with previous work that provide the expected
<figure><img src="image_12.png" /><figcaption>**Fig. 11.** Left panel: integrated intensity map of the modeled N<sub>2</sub>H<sup>+</sup> isolated component within a region of radius of 0.1 pc, centered at V-shape “C”. With boxes in shades of green and yellow, we indicate the different paths taken to create the PV diagrams shown on the right panel. The area covered by these four PVs matches the extent of this V-shape (see Fig. 12). The 1.3 mm cores are indicated using the same convention from Fig. 9. For the cores within the radius of 0.1 pc we indicate their IDs in black. Right panel: PV diagrams associated with the colored rectangles in the left panel. At the bottom right corner of each sub panel, with colored values, we indicate the angle (counter-clockwise) of each PV path. The PV paths match the V-shape coverage (right panel of Fig. 12). We see the overall structure of this V-shape persists at different angles, indicating that these structures are not a result of projection in the POS.</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>**Fig. 12.** Multi-tracer diagram at the location of the main blue V-shape of G353 (Sect. 5.3). The different panels show the step by step construction of the final plot (right panel). Left panel: the black to white background shows the 1.3 mm continuum emission from Díaz-González et al. (2023). With black ellipses we show the 1.3 mm continuum cores from Louvet et al. (2024). With black text we indicate the IDs of the cores closest to the center of the V-shape, marked with a fuchsia “x”. We indicate the barycenter of cores 2 & 3 (see Fig. 11) with a yellow “x”. Filled contours represent the CO (2–1) emission in the velocity ranges of –50 km s<sup>-1</sup> to –15 km s<sup>-1</sup> (cyan) and 15 km s<sup>-1</sup> to 50 km s<sup>-1</sup> (beige), relative to the V<sub>LSR</sub> = –17 km s<sup>-1</sup> of G353. With a white contour we indicate the log(N(H<sub>2</sub>) cm<sup>2</sup>) = 23.3. Middle panel: in addition to the left panel, we include the SiO (5–4) emission with blue and red contours in the same (negative and positive) velocity range as with CO. The pink contour indicates the integrated intensity of the most blueshifted modeled N<sub>2</sub>H<sup>+</sup> Gaussian velocity component (see Fig. 9; blue distribution), at a value of 7 K km s<sup>-1</sup>. Right panel: with open circles we show the location of the data composing the main N<sub>2</sub>H<sup>+</sup> blue V-shape. The colors indicate their velocity centroid, and their (increasing) size indicates how close the gas velocities are to the velocity apex of the V-shape. This VG seems to converge to the barycenter of cores 2 and 3, and it is oriented along filament “F3”.</figcaption></figure>

PV distributions for spherical protostellar envelopes under infall (Tobin et al. 2012). At large scales, we interpret the agreement between the PV diagrams of our model and the data as proto-cluster scale collapse due to gravitational contraction. It is worth noting that the inferred mass from our model is ~5.5 times lower than the one derived from the N(H<sub>2</sub>) map (Díaz-González et al. 2023). We speculate that a model considering complex processes such as turbulence, radiative transfer, and magnetic fields might

solve this mass discrepancy while still matching the observed PV distribution.

## 7. Mass accretion rates in the V-shaped structure

Based on the idea that the V-shapes are a result of gas flowing toward cores, in this section we provide estimates of their mass accretion rates ($\dot{M}_{in}$) for N<sub>2</sub>H<sup>+</sup> and H<sub>2</sub>.
<figure><img src="image_14.png" /><figcaption>Fig. 13. V-shapes location on the plane of the sky. We indicate the position of each V-shape with colored points and their ID with black text. With red lines we show the filamentary structure identified in Sect. 3.2. In the background we show the integrated intensity map of all three velocity components, similarly as in Fig. 9. We see most of the V-shapes are located at the hub.</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>Fig. 14. Mean spectra within a 1.14" (~0.01 pc) radius around the location of the main blue V-shape (pink "x" in Fig. 12). Both N₂H⁺ and H₂CO show blue asymmetry, known to characterize infall motions. The difference in velocity between the two N₂H⁺ peaks is ~2.5 km s⁻¹.</figcaption></figure>

### 7.1. H₂ mass accretion rate

Here we describe our approach for V-shape “C”, and to ensure that we estimate the V-shape M(H₂) on the area of this structure, we used the CASA task imregrid to obtain the continuum-derived N(H₂) map from Díaz-González et al. (2023) at the resolution of the N₂H⁺ data. We determined that in this V-shape the total N(H₂) ~ 1.17 × 10²⁶ cm⁻² in an area of 0.013 pc². Here, we derived a M(H₂) map using Eq. (5):

$$M(\text{H}_2) = 2 \times N(\text{H}_2) \times \text{area}_{\text{pixel}} \times \text{m}_{\text{proton}}, \quad (5)$$

from this M(H₂) map we considered only the points that are part of the V-shape. To determine the mass associated with flowing A74, page 12 of 23

<figure><img src="image_16.png" /><figcaption>Fig. 15. PV coverage of a gravitationally collapsing sphere. The white contour represents the coverage of the synthetic radial velocities derived from this model (Sect. 6). The background and cores are the same as in Fig. 8. For the modeled sphere we set its total mass to 150 M☉, within a radius of 0.5 pc.</figcaption></figure>

motions we subtracted the core masses (from Louvet et al. 2024) that are located inside this V-shape. We note here that this mass map is an upper limit given that we do not apply a background correction. We obtained a total of M(H₂) ~ 53 M☉. Considering t<sub>VG mean</sub> used in Eq. (8), we derived the M<sub>in</sub> (H₂) as:

$$\dot{M}_{\text{in}}(\text{H}_2) = \frac{M(\text{H}_2)}{t_{\text{VG mean}}} = 8.22 \times 10^{-4} \text{ M}_{\odot} \text{ yr}^{-1}. \quad (6)$$

We used the procedure described above to estimate the M<sub>in</sub> (H₂) of the other eight V-shapes shown in Fig. D.2. We include these values in Table D.1. The average M<sub>in</sub> (H₂) of these V-shapes is 3.4 × 10⁻⁴ M☉ yr⁻¹. We note that V-shapes “H”, “C”, “F”, and “B” present the largest M<sub>in</sub> (H₂) and they are located near or at the convergence point of the filaments (see Fig. 13).

For comparison, using the core masses from Louvet et al. (2024) we estimate the free-fall time of all 45 1.3 mm cores and their mass accretion rates. These values present large scattering, ranging from (0.07–25) × 10⁻⁴ M☉ yr⁻¹, with 28 of these cores presenting M<sub>in</sub> (H₂) < 10⁻⁴ M☉ yr⁻¹. For cores 2 & 3, the average M<sub>in</sub> (H₂) is 15.5 × 10⁻⁴ M☉ yr⁻¹, about twice the M<sub>in</sub> (H₂) of the main blue V-shape, located between these two cores.

### 7.2. N₂H⁺ mass accretion rate and relative abundance

To derive the N₂H⁺ M<sub>in</sub> associated with the main blue V-shape (Fig. 10), we need to estimate its N₂H⁺ mass. For this procedure we used PySpecKit with the n2hp_vtau fitter, to fit the full N₂H⁺ line. The fitted parameters (see below) allow us to derive the N(N₂H⁺). The V-shape structure contains 77, 143, and 25 N₂H⁺ spectra with one, two, and three velocity components, respectively. The bluest velocities in the three velocity component spectra accounts for ~2% of the total number of velocities
<table><caption>**Table 2.** N₂H⁺ full line fitting parameter ranges.</caption><tr><td>Excitation temperature (K)</td><td>(T<sub>ex</sub>)</td><td>2.73–80</td></tr><tr><td>Opacity</td><td>(τ)</td><td>0.01–40</td></tr><tr><td>Centroid velocity (km s<sup>-1</sup>)</td><td>(v)</td><td>-25–15</td></tr><tr><td>Velocity dispersion (km s<sup>-1</sup>)</td><td>(Δv)</td><td>0.20–3</td></tr></table>

**Notes.** Ranges used for the full N₂H⁺ line fitting.

in this V-shape. For this reason, we modeled the full N₂H⁺ hyperfine line structure with one and two velocity components. In Table 2 we list the parameters and ranges used for this procedure.

After obtaining the modeled N₂H⁺ cube, we removed modeled components that presented $\tau_{RMS}/\tau > 0.3$, in which $\tau$ and $\tau_{RMS}$ represent the estimated opacity and its associated error, respectively. This criterion is to ensure that we use reliable fitted parameters to determine our N(N₂H⁺) values. For the fitting of two N₂H⁺ components we only analyzed the most blueshifted component. The resulting opacities follow a log-normal distribution with a peak at $\tau = 0.13$.

We derived the N(N₂H⁺) of the V-shape by using Eq. (7) (Caselli et al. 2002b):

$$N(\text{N}_2\text{H}^+) = \frac{8\pi^{3/2}\Delta v}{2\sqrt{\ln 2}\lambda^3 A} \frac{g_l}{g_u} \frac{\tau}{1 - \exp(-hv/kT_{\text{ex}})} \frac{Q_{\text{rot}}}{g_l \exp(-E_l/kT_{\text{ex}})}, \quad (7)$$

in which $\tau$, $T_{\text{ex}}$, and $\Delta v$ are the opacity, excitation temperature, and velocity dispersion, respectively, obtained from the full line fitting. The Planck and Boltzmann constants are represented by $h$ and $k$, respectively, $\nu$ and $\lambda$ are the frequency and wavelength of N₂H⁺, $A$ is the Einstein coefficient of the N₂H⁺ (1–0) transition, $g_l$ and $g_u$ are the statistical weights of the lower and upper energy levels, $E_l$ is the energy of the lower level, $Q_{\text{rot}}$ is the partition function estimated using the excitation temperature of the full N₂H⁺ fits (see Eq. (A2) from Caselli et al. 2002b).

From the above procedure, inside the V-shaped structure, we get a total N(N₂H⁺) = 5.24 × 10¹⁵ cm⁻² and a total M(N₂H⁺) = 5.9 × 10⁻⁸ M☉, within its extent of 0.013 pc². We used the average timescale of the VGs from Sect. 5.3 (Fig. 10), $t_{\text{VG mean}} = 64.5$ kyr, to determine the $\dot{M}_{\text{in}}$ (N₂H⁺) as:

$$\dot{M}_{\text{in}}(\text{N}_2\text{H}^+) = \frac{M(\text{N}_2\text{H}^+)}{t_{\text{VG mean}}} = 9.1 \times 10^{-13} \text{ M}_\odot \text{ yr}^{-1}. \quad (8)$$

The $\dot{M}_{\text{in}}(\text{N}_2\text{H}^+)$ estimate (and $\dot{M}_{\text{in}}(\text{H}_2)$ below) should be multiplied by $\sin(\theta)$, in order to account for the unknown inclination angle ($\theta$) of the protocluster/filaments relative to the POS.

For the main blue V-shape (Fig. 10), we derived the N₂H⁺ relative abundance X(N₂H⁺), using the N₂H⁺ and H₂ column densities obtained above, as:

$$X(\text{N}_2\text{H}^+) = \frac{N(\text{N}_2\text{H}^+)}{N(\text{H}_2)} = 4.8 \times 10^{-11}. \quad (9)$$

The X(N₂H⁺) value obtained above appears lower than typical estimates in massive Galactic star-forming regions reported in different works ([1.6–3.8] × 10⁻¹⁰; Caselli et al. 2002a; Henshaw et al. 2014, Sandoval-Garrido et al. in prep.). We consider it a good agreement taking into account the uncertainties of the involved measurements (i.e., column density estimates). For a comparison of the N₂H⁺ relative abundance between regions composing the V-shapes and the rest of protocluster we would require to model the full N₂H⁺ line emission in the whole region.

# 8. Discussion

## 8.1. V-shaped velocity gradients in the literature

The V-shaped VGs described in this work have been detected across multiple Galactic star-forming system. Stutz & Gould (2016) introduced the Slingshot mechanism in the Integral Shaped Filament (ISF) located in Orion A. They show undulations of the region in both position and velocity, suggesting that these features appear to be ejecting protostars (see their Fig. 2). Furthermore, Stutz (2018) characterize a standing wave in the neighborhood of the ISF, consistent with the Slingshot mechanism. It is possible that the undulations in the works above might result in the observed V-shaped structures seen in different studies (see below). González Lobos & Stutz (2019) identify six evenly spaced (every ~0.44 pc) velocity peaks along the spine of the ISF in Orion A. They suggest that this periodicity is consistent with the wave-like perturbation in the gas caused by the Slingshot mechanism. In Álvarez-Gutiérrez et al. (2021) they analyze the L1482 filament located in the California Molecular Cloud. In all of the analyzed tracers there is a clear velocity peak in their north region (length ~1.8 pc, mass ~10³ M☉). This subregion contains a higher gas density and higher number of YSOs compared to the more quiescent south part.

While the two regions described above are considered nearby (both at D ≲ 500 pc), more distant regions also present these velocity features which we list below. Zhou et al. (2022) study the velocity profiles along filaments from the ATOMS survey (Liu et al. 2020). The median mass of their sources is ~1.4 × 10³ M☉ with a median length of the filaments at ~1.35 pc. By analyzing the H¹³CO⁺ (1–0) emission they find converging VGs along filaments (see their Figs. 6 and 10), which they also detected using simulations from Gómez & Vázquez-Semadeni (2014). These VGs at scales comparables to the V-shapes presented here are consistent with our VGs estimates (see their Figs. 7 and 8). In Zhou et al. (2023) they analyzed H¹³CO (2–1) APEX/LASMA data of the G333 complex. They identify multiple V-shaped VG (see their Fig. 7) which they describe as the PV projection of a funneling structure in PPV space (see their Fig. 9). The origin of this structure is due to material inflowing toward the central hub and also due to gravitational contraction of star-forming clouds or clumps.

Redaelli et al. (2022) use ALMA N₂H⁺ (1–0) isolated component data of the high-mass (5200 M☉) clump AGAL014.492-00.139 identifying multiple coherent structures in PPV space (trees “B” and “G”; right panel of their Figs. 7 and 9). These are characterized by multiple undulations, and possible V-shaped VGs. For their “G” PV distribution, they suggest that one scenario is where the dense gas is flowing along the filament (of length ~0.2 pc) from protostar “p3” toward the protostar “p2”. This motion has an $\dot{M}_{\text{in}} = 2.2 \times 10^{-4} \text{ M}_\odot \text{ yr}^{-1}$, being in the range of the $\dot{M}_{\text{in}}$ we derive for our VGs (see Table D.1).

In Rawat et al. (2024) they analyze H¹³CO(1–0) data obtained with the Purple Mountain Observatory, as part of the Milky Way Imaging Scroll Painting survey. They detect a V-shaped VG (see their Fig. 14) along the ridge of the G148.24+00.41 (G148) cloud. This V-shape peaks toward the dense clump at the center of this region, possibly indicating gas inflow along their filaments F2 and F6 toward the hub. The length of the V-shape in G148 is ~15 pc, while our most prominent V-shape (Fig. 10) is ~0.2 pc. Also the masses and lengths of their identified filaments are (1.3 to 6.9) × 10³ M☉ and 14 to 38 pc, respectively, large compared to the total mass (2.5 × 10³ M☉) and extent of G353 (~1.2 pc). This difference in probed lengths and masses might be reflected by their mean VG ~0.05 km s⁻¹ pc⁻¹, ~2.5 orders of magnitude
smaller than our VGs. This is consistent with the analysis pre-
sented in Zhou et al. (2022, 2023), in which they observe an
inverse relation between the spatial scale of a region and their
VGs (see their Fig. 8).

In Pan et al. (2024) using APEX C¹⁸O (2–1) data of the fil-
amentary cloud G034.43+00.24 (G34) they identify converging
VGs of lengths ~1 pc toward the “middle ridge” (see their Fig. 3,
top panel). They interpret these VG as gas flowing from the fil-
aments onto dense clumps, located at the center of G34. These
VGs of their southern and northern filaments are in the range
of ~0.3–0.4 km s⁻¹ pc⁻¹, and they estimate the total mass inflow
rate toward the middle ridge as ~5.5 × 10⁻⁴ M☉ yr⁻¹, similar to
our Ṁin estimates.

Current work by Sandoval-Garrido et al. (in prep.) in G351.77 (intermediate protocluster, located at 2 kpc; Motte et al. 2022; Reyes-Reyes et al. 2024) use a similar analysis as we present in this work, in which they identify multiple V-shaped velocity structures. In Salinas et al. (in prep) they analyze the kinematics of the evolved protocluster G012.80 (located at 2.4 kpc; Motte et al. 2022), in which they implement similar techniques and find velocity signatures of filamentary rotation.

8.2. Filamentary 3D morphology

V-shaped VGs appear to be a generic feature across a wide range
of star-forming environments, probing VGs with differences of
up to ~2 orders of magnitude in spatial scales ranging from 0.1 to
~10 pc. Despite being commonly detected in recent studies, it is
still not clear how they are produced. Henshaw et al. (2014) high-
lights the degeneracy regarding the opposite interpretations of
these V-shaped velocity structures. They suggest that these VGs
can be a signature of gas flows along kinked filaments toward
a core located at their convergence point. From our analysis
regarding the most prominent V-shape (see Sects. 5.3 and 6) we
see that no core is located at its apex, although cores 2 and 3 are
within ~0.05 pc. Also the spatial offset between the center of the
V-shape with the barycenter between these two cores is ~600 au.
This is consistent with the idea of small-scale gravitational col-
lapse within the protocluster, similar to clump decoupling from
their parent molecular cloud (Peretto et al. 2023). Based on this,
we conclude that cores may be located in the vicinity of the
velocity apex, and not necessarily on top of it. This gas flows
toward denser regions may result in the formation of high-mass
cores in later stages during the evolution of the protocluster.

Regarding the kinked morphology of the regions hosting V-
shapes, one scenario regarding magnetized shocks is presented
in Inoue & Fukui (2013) and Inoue et al. (2018). They use mag-
netohydrodynamics simulations to characterize the interaction of
molecular clouds and a magnetized shock produced by a cloud-
cloud collision. They find that the shock layer decelerates as
it collides with denser regions. This deceleration reshapes the
shock layer to be oblique, leading to the formation of kinked fila-
ments and converging flows, which are oriented toward the apex
of these filaments. They predict that magnetic fields present in
the region should be perpendicular to these filaments and bend
with the shock around the filament (Inoue et al. 2018, see their
Fig. 3). In Bonne et al. (2020) and Bonne et al. (2023) they pro-
pose that this scenario takes place in the Musca and the DR21
filaments. In both of these regions they detect V-shaped VGs
which they suggest are the result of cloud-cloud collisions bend-
ing the magnetic field (Bonne et al. 2020, see their Figs. 22 and
23). Further observations of magnetic field polarization in the
POS, along with information along the line of sight is required
to evaluate these models.

Another possibility is that these kinked structures could be
caused by mechanisms such as the Slingshot. This mechanism
proposes a standing wave, longitudinal gravitational instabilities,
or large-scale oscillations possibly caused by a possibly helical
magnetic field morphology, causing ejections of protostars and
protoclusters from their maternal filament (Stutz & Gould 2016;
Stutz 2018; Stutz et al. 2018; Liu et al. 2019).

A different interpretation is that they are the product of out-
flowing material coming from a forming protostar interacting
with the surrounding dense gas (see their Fig. 12). To shed
light into this degeneracy in G353 we compare the N₂H⁺ (dense
gas tracer) radial velocities and SiO (shock/outflow tracer) as
a proxy for energies. The velocity range ΔV covered by the
N₂H⁺ emission is ~8 km s⁻¹, while for SiO is ~80 km s⁻¹. Given
the difference in probed velocities between SiO and N₂H⁺ (see
Sect. 5.3) and the analysis presented in Sect. 6 we suggest that
the V-shapes in G353 are a signature of infall.

The multiple VGs that conform the V-shapes present in G353
have values of ~8 to ~31 km s⁻¹ pc⁻¹, with timescales ranging
from ~35–173 kyr, and Ṁin values (~0.4–9) × 10⁻⁴ M☉ yr⁻¹.
These values are similar to VGs in other regions with compara-
ble sizes and masses. In G353 it is likely that these VGs are the
result of dense gas moving through filaments, possibly increas-
ing the density of the central regions, shaping the overall velocity
field at large and small scales, and leading to a further increase
in the core population and star formation activity.

8.3. Timescales and mass accretion rates

One important aspect of the V-shapes that is still not well under-
stood is the timescale associated with the VGs (t<sub>VG</sub> = VG<sup>-1</sup>). It
is not clear if nor how these timescales determine core formation
lifetimes or impact the star formation environment in general. In
our sample of V-shapes the timescales are in tens of kyrs with
the average value of ~67 kyr, ~ two times the cores t<sub>ff</sub>, while the
t<sub>ff</sub> of the whole protocluster is ~0.21 Myr. In Rawat et al. (2024)
they estimate the longitudinal collapse timescales for their
filaments, being in the range of 5–15 Myr. Using their derived
VGs we estimate their associated timescales to be between
~16 and ~50 Myr, ~ one to two orders of magnitude larger than
our small-scale V-shapes timescales. We suggest that the VG
timescales might serve as an upper limit for filamentary collapse
timescales. In Zhou et al. (2022) they determine gas accretion
times as a function of the lengths of their filaments, assuming
that the VGs produced by gas inflow (see their Fig. 11). At
filament lengths comparable to our V-shapes (~0.1 pc) their
gas accretion timescales are on the order of our estimates (see
Table D.1).

8.4. Depletion timescales

It is also interesting to consider the mass accretion rates mea-
sured here compared to the available protocluster mass reservoir
to explore implications for the duration of the gas dominated
phase. The total mass accretion rate of our V-shaped structures
is Ṁin, Tot = 3 × 10⁻³ M☉ yr⁻¹ (see Table D.1). Considering the
total mass (ṀTot) of G353 as a mass reservoir, we estimated
the depletion timescale (t<sub>dep</sub>) as the time needed to fully con-
sume the gas. Here we assumed that Ṁin, Tot is representative of
flows feeding gas onto cores. We estimated t<sub>dep</sub> = ṀTot / Ṁin, Tot,
in which the total mass of the region is 2500 Ṁ (Motte et al.
2022). We obtained a t<sub>dep</sub> = 0.8 Myr, of similar magnitude but
about four times larger than the t<sub>ff</sub> of the protocluster. Consid-
ering that our estimate of Ṁin, Tot is certainly a lower limit (see
discussion above), the actual value of t<sub>dep</sub> is likely to be
shorter, so closer to the free-fall time estimate. Given that the
protocluster does not appear to be in a state of free-fall (see Sect. 6) but instead undergoing comparatively slow gravitational contraction, the similarity in these relatively crude estimates seems remarkable. While we do not yet have an explanation for why relatively good match in timescales, it would seem to indicate that protocluster evolution may be a self-regulating process. Larger samples and similar analysis will test this hypothesis.

Moreover, the approximate concordance of $t_{\text{dep}}$ and $t_{\text{ff}}$ may indicate that the “phase transition” of protocluster gas mass being converted into stellar mass could contribute a relevant “negative pressure” counteracting effects of, for example, feedback over the lifetime of the protocluster.

## 9. Summary and conclusions

We characterize the complex dense gas kinematics of G353 using ALMA-IMF LP observations. The data used in this paper mainly consist of the fully combined $N_2H^+$ data cube, but we also include 1.3 mm continuum cores and DCN cores velocity catalogs, SiO 12 m observations, and a $N(H_2)$ 1.3 mm continuum-derived map. We summarize our main results below:

- With our $N_2H^+$ isolated component modeling, we find that most of the 1.3 mm cores are located in regions with two to three velocity components. This indicates kinematic complexity down to ~4 kau scales;

- We increase the number of cores with a $V_{\text{LSR}}$ estimate in this region by further examining the DCN emission and comparing it with the $N_2H^+$ data extracted toward the core positions. We find that 11 cores that were previously undetected in the DCN background-subtracted fitting from Cunningham et al. (2023) are identified with our method. With this approach we increase our core velocities sample from 15 to 26, accounting for ~58% of the total 45 1.3 mm continuum cores. These are presented in Table C.1;

- We show that the traditional PV diagram highlights large, protocluster scale kinematics. In contrast, the intensity-weighted PV diagram allows us access to the small, core scale dynamics (see Figs. 8 and 9);

- From the PV diagrams, we see the DCN core velocities are in agreement with the $N_2H^+$ velocity distribution (within a few DCN channel widths). This suggests coupling between the cores and the dense gas in which they formed;

- In the intensity-weighted PVs we see clear V-shaped velocity structures, composed by two linear VGs converging into a common point. These VGs are present across all $N_2H^+$ velocity components. Some of them are near the location of cores in both position and velocity (see Sect. 5.2);

- We successfully characterize nine V-shaped VGs well detected in our $N_2H^+$ data (see Figs. 13 and D.2). These structures are located mostly at the center of the protocluster, where three filaments converge;

- V-shape “C” (see Figs. 10–12) is the most prominent across our sample. It is centered between cores 2 and 3, two of the most massive cores in this region;

- We estimate the barycenter of cores 2 & 3, presenting an offset relative to the center of the V-shape of ~0.3" (~600 au) well below the beam size of our $N_2H^+$ data;

- For V-shape “B” we find that core “7”, with a mass of $6 M_\odot$, is located within a ~beam size from its apex;

- We suggest that the dense gas is flowing along the filament, producing the V-shaped structure toward the derived barycentre;

- We characterize the VGs composing our sample of V-shapes by applying linear fits to these distributions. We estimate timescales associated with the VGs as $t_{\text{VG}} = \text{VG}^{-1}$. These

timescales are between ~35 to ~170 kyr, with an average of ~67 kyr. These values are short compared to the $t_{\text{ff}}$ of the protocluster (~0.21 Myr), and ~2 times larger than the cores average $t_{\text{ff}}$ (~32 kyr);

- We suggest that at small scales the $N_2H^+$ V-shaped structures indicate gas motions along filaments, toward denser regions. Thus we interpret $t_{\text{VG}}$ as inflow timescales;

- Using an $H_2$ mass map and the V-shapes mean timescales, we derive $H_2$ mass accretion rates of (0.35 to 8.77) × $10^{-4} M_\odot \text{ yr}^{-1}$, consistent with previous studies on regions that present gas flows along filaments toward denser object or regions, such as protostars and clumps. Moreover, V-shapes “H”, “C”, “F”, and “B” present the largest $\dot{M}_{\text{in}}$ ($H_2$) and they are located near or at the convergence point of the filaments (see Fig. 13);

- In SiO, the PV structure covers a velocity range ($\Delta V$) of ~80 km s$^{-1}$, while for $N_2H^+$ $\Delta V$ is ~8 km s$^{-1}$. This difference suggests that $N_2H^+$ is tracing infall, a less energetic processes compared to SiO, a shock and outflow tracer;

- We model the protocluster as a gravitationally collapsing sphere. The derived radial velocities are consistent with the large-scale morphology of the traditional PV diagram. This agreement suggests that at large scales the G353 protocluster is undergoing gravitational contraction.

Overall, it is imperative to replicate the kinematic analysis presented in this work in the remaining ALMA-IMF fields and other Galactic star-forming regions. By increasing the sample of analyzed fields we might find correlations between evolutionary state (young, intermediate, or evolved; see Motte et al. 2022), star formation activity, cores and outflow population properties, and their velocity field. This approach will allow us to better describe the kinematic processes taking place in this early stage of star formation.

## Data availability

Interactive figure is available as supplementary material at: https://www.aanda.org

**Acknowledgements.** The authors thank the referee for helpful comments that improved the text. This paper makes use of the following ALMA data: ADS/JAO.ALMA#2017.1.01355.L. 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. We thank Elena Redaelli, Diego R. Matus Carrillo, and Vineet Rawat for very helpful discussions. R.A. gratefully acknowledges support from ANID Beca Doctorado Nacional 21200897. A.S. gratefully acknowledges support by the Fondecyt Regular (project code 1220610) and ANID BASAL project FB210003. F.L. acknowledges the support of the Marie Curie Action of the European Union (project MagiKStar, Grant agreement number 841276) F.M. acknowledges support from the French Agence Nationale de la Recherche (ANR) under reference ANR-20-CE31-009, of the Programme National de Physique Stellaire and Physique et Chimie du Milieu Interstellaire (PNPS and PCMI) of CNRS/INSU (with INC/INP/IN2P3). R.G.M. and D.D.G. acknowledge support from UNAM-PAPIIT project IN108822 and from CONACyT Ciencia de Frontera project ID 86372. F.M., F.L., and N.C. acknowledge support from the European Research Council (ERC) via the ERC Synergy Grant ECOGAL (grant 855130). N.C. acknowledges funding from the ERC under the European Union's Horizon 2020 research. P.S. was partially supported by a Grant-in-Aid for Scientific Research (KAKENHI Number JP22H01271 and JP23H01221) of JSPS. M.B. is a postdoctoral fellow in the University of Virginia's VICO collaboration and is funded by grants from the NASA Astrophysics Theory Program (grant number 80NSSC18K0558) and the NSF Astronomy & Astrophysics program (grant number 2206516). A.G. acknowledges support from the NSF under grants AAG 2008101 and CAREER 2142300. T.Cs. has received financial support from the French State in the framework of the IdEx Université de Bordeaux Investments for the future Program. S.D.R. acknowledges the funding and support of ANID-Subdirección de Capital Humano Magíster/Nacional/2021-22211000. T.B. acknowledges the support from S. N. Bose National Centre for Basic Sciences under the
Department of Science and Technology, Govt. of India. G.B. acknowledges
financial support from the grants PID2020-117710GB-100 and CEX2019-000918
funded by MCIN/AEI/10.13039/501100011033. A.K. and L.B. gratefully
acknowledge support from ANID BASAL project FB210003. F.O. acknowledges
the support of the Ministry of Science and Technology of Taiwan, projects No.
109-2112-M-007-008-, 110-2112-M-007-023-, and 110-2112-M-007-034-.

References

Álvarez-Gutiérrez, R. H., Stutz, A. M., Law, C. Y., et al. 2021, ApJ, 908, 86

André, P., Men'shchikov, A., Bontemps, S., et al. 2010, A&A, 518, L102

Anglada, G., Rodriguez, L. F., Canto, J., Estalella, R., & Lopez, R. 1987, A&A,
186, 280

Armante, M., Gusdorf, A., Louvet, F., et al. 2024, A&A, 686, A122

Bastian, N., Covey, K. R., & Meyer, M. R. 2010, ARA&A, 48, 339

Bergin, E. A., Alves, J., Huard, T., & Lada, C. J. 2002, ApJ, 570, L101

Bonfand, M., Csengeri, T., Bontemps, S., et al. 2024, A&A, 687, A163

Bonne, L., Bontemps, S., Schneider, N., et al. 2020, A&A, 644, A27

Bonne, L., Bontemps, S., Schneider, N., et al. 2023, ApJ, 951, 39

Busquet, G., Zhang, Q., Palau, A., et al. 2013, ApJ, 764, L26

Caselli, P., Myers, P. C., & Thaddeus, P. 1995, ApJ, 455, L77

Caselli, P., Myers, P. C., & Thaddeus, P. 1995, ApJ, 455, L77

Caselli, P., Benson, P. J., Myers, P. C., & Tafalla, M. 2002a, ApJ, 572, 238

Caselli, P., Walmsley, C. M., Zucconi, A., et al. 2002b, ApJ, 565, 344

Cazzoli, G., Corbelli, G., Degli Esposti, C., & Favero, P. 1985, Chem. Phys. Lett.,
118, 164

Chen, H.-R. V., Zhang, Q., Wright, M. C. H., et al. 2019, ApJ, 875, 24

Csengeri, T., Bontemps, S., Schneider, N., Motte, F., & Dib, S. 2011, A&A, 527, A135

Csengeri, T., Bontemps, S., Wyrowski, F., et al. 2017, A&A, 601, A60

Cunningham, N., Lumsden, S. L., Cyganowski, C. J., Maud, L. T., & Purcell, C. 2016, MNRAS, 458, 1742

Cunningham, N., Ginsburg, A., Galván-Madrid, R., et al. 2023, A&A, 678, A194

Díaz-González, D. J., Galván-Madrid, R., Ginsburg, A., et al. 2023, ApJS, 269, 55

Fernández-López, M., Arce, H. G., Looney, L., et al. 2014, ApJ, 790, L19

Galván-Madrid, R., Zhang, Q., Keto, E., et al. 2010, ApJ, 725, 17

Galván-Madrid, R., Liu, H. B., Zhang, Z. Y., et al. 2013, ApJ, 779, 121

Ginsburg, A., & Mirocha, J. 2011, PySpecKit: Python Spectroscopic Toolkit, Ver.
0.1.23, Astrophysics Source Code Library [record ascl:1109.0001]

Ginsburg, A., Sokolov, V., de Val-Borro, M., et al. 2022, AJ, 163, 291

Gómez, G. C., & Vázquez-Semadeni, E. 2014, ApJ, 791, 124

González Lobos, V., & Stutz, A. M. 2019, MNRAS, 489, 4771

Hacar, A., Tafalla, M., Forbrich, J., et al. 2018, A&A, 610, A77

Henshaw, J. D., Caselli, P., Fontani, F., Jiménez-Serra, I., & Tan, J. C. 2014, MNRAS, 440, 2860

Henshaw, J. D., Ginsburg, A., Haworth, T. J., et al. 2019, MNRAS, 485, 2457

Henshaw, J. D., Kruijssen, J. M. D., Longmore, S. N., et al. 2020, Nat. Astron.,
4, 1064

Inoue, T., & Fukui, Y. 2013, ApJ, 774, L31

Inoue, T., Hennebelle, P., Fukui, Y., et al. 2018, PASJ, 70, S53

Koch, E. W., & Rosolowsky, E. W. 2015, MNRAS, 452, 3435

Kumar, M. S. N., Arzoumanian, D., Men'shchikov, A., et al. 2022, A&A, 658, A114

Kuznetsova, A., Hartmann, L., & Ballesteros-Paredes, J. 2015, ApJ, 815, 27

Kuznetsova, A., Hartmann, L., & Ballesteros-Paredes, J. 2018, MNRAS, 473, 2372

Lee, C. W., Myers, P. C., & Tafalla, M. 1999, ApJ, 526, 788

Lee, C. W., Myers, P. C., & Tafalla, M. 2001, ApJS, 136, 703

Lippok, N., Launhardt, R., Semenov, D., et al. 2013, A&A, 560, A41

Liu, H. B., Galván-Madrid, R., Jiménez-Serra, I., et al. 2015, ApJ, 804, 37

Liu, H.-L., Stutz, A., & Yuan, J.-H. 2019, MNRAS, 487, 1259

Liu, T., Evans, N. J., Kim, K.-T., et al. 2020, MNRAS, 496, 2790

Louvet, F., Sanhueza, P., Stutz, A., et al. 2024, A&A, in press, https://doi.org/10.1051/0004-6361/202345986

Mardones, D., Myers, P. C., Tafalla, M., et al. 1997, ApJ, 489, 719

Men'shchikov, A. 2021, A&A, 649, A89

Minh, Y. C., Liu, H. B., Galvań-Madrid, R., et al. 2018, ApJ, 864, 102

Motte, F., Bontemps, S., & Louvet, F. 2018, ARA&A, 56, 41

Motte, F., Bontemps, S., Csengeri, T., et al. 2022, A&A, 662, A8

Myers, P. C. 2009, ApJ, 700, 1609

Nony, T., Galván-Madrid, R., Motte, F., et al. 2023, A&A, 674, A75

Offner, S. S. R., Clark, P. C., Hennebelle, P., et al. 2014, in Protostars and Planets VI, eds. H. Beuther, R. S. Klessen, C. P. Dullemond, & T. Henning, 53

Olguin, F. A., Sanhueza, P., Chen, H.-R. V., et al. 2023, ApJ, 959, L31

Pan, S., Liu, H.-L., & Qin, S.-L. 2024, ApJ, 960, 76

Peretto, N., Fuller, G. A., André, P., et al. 2014, A&A, 561, A83

Peretto, N., Rigby, A. J., Louvet, F., et al. 2023, MNRAS, 525, 2935

Pouteau, Y., Motte, F., Nony, T., et al. 2023, A&A, 674, A76

Rawat, V., Samal, M. R., Walker, D. L., et al. 2024, MNRAS, 528, 2199

Redaelli, E., Bovino, S., Sanhueza, P., et al. 2022, ApJ, 936, 169

Reyes-Reyes, S. D., Stutz, A. M., Megeath, S. T., et al. 2024, MNRAS, 529, 2220

Sanhueza, P., Contreras, Y., Wu, B., et al. 2019, ApJ, 886, 102

Sanhueza, P., Girart, J. M., Padovani, M., et al. 2021, ApJ, 915, L10

Schuller, F., Menten, K. M., Contreras, Y., et al. 2009, A&A, 504, 415

Smith, R. J., Shetty, R., Stutz, A. M., & Klessen, R. S. 2012, ApJ, 750, 64

Storm, S., Mundy, L. G., Fernández-López, M., et al. 2014, ApJ, 794, 165

Stutz, A. M. 2018, MNRAS, 473, 4890

Stutz, A. M., & Gould, A. 2016, A&A, 590, A2

Stutz, A. M., & Kainulainen, J. 2015, A&A, 577, L6

Stutz, A. M., Gonzalez-Lobos, V. I., & Gould, A. 2018, arXiv e-prints
[arXiv:1807.11496]

Tafalla, M., Myers, P. C., Caselli, P., & Walmsley, C. M. 2004, A&A, 416, 191

Tafalla, M., Usero, A., & Hacar, A. 2021, A&A, 646, A97

Tobin, J. J., Hartmann, L., Bergin, E., et al. 2012, ApJ, 748, 16

Towner, A. P. M., Ginsburg, A., Dell'Ova, P., et al. 2024, ApJ, 960, 48

Ungerechts, H., Bergin, E. A., Goldsmith, P. F., et al. 1997, ApJ, 482, 245

Wienen, M., Wyrowski, F., Menten, K. M., et al. 2015, A&A, 579, A91

Zhou, J.-W., Liu, T., Evans, N. J., et al. 2022, MNRAS, 514, 6038

Zhou, J. W., Wyrowski, F., Neupane, S., et al. 2023, A&A, 676, A69

<footnote n="1">¹ Departamento de Astronomía, Universidad de Concepción, Casilla 160-C, Concepción, Chile
e-mail: rodralvarezz@gmail.com</footnote>

<footnote n="2">² Univ. Grenoble Alpes, CNRS, IPAG, 38000 Grenoble, France</footnote>

<footnote n="3">³ Instituto de Radioastronomía y Astrofísica, Universidad Nacional Autónoma de México, Morelia, Michoacán 58089, Mexico</footnote>

<footnote n="4">⁴ SKA Observatory, Jodrell Bank, Lower Withington, Macclesfield SK11 9FT, UK</footnote>

<footnote n="5">⁵ National Astronomical Observatory of Japan, 2-21-1 Osawa, Mitaka, Tokyo 181-8588, Japan</footnote>

<footnote n="6">⁶ Astronomical Science Program, The Graduate University for Advanced Studies, SOKENDAI, 2-21-1 Osawa, Mitaka, Tokyo 181-8588, Japan</footnote>

<footnote n="7">⁷ Departments of Astronomy, University of Virginia, Charlottesville,
VA 22904, USA</footnote>

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

<footnote n="9">⁹ Laboratoire de Physique de l'École Normale Supérieure, ENS, Université PSL, CNRS, Sorbonne Université, Université Paris Cité, 75005 Paris, France</footnote>

<footnote n="10">¹⁰ Observatoire de Paris, PSL University, Sorbonne Université, LERMA, 75014, Paris, France</footnote>

<footnote n="11">¹¹ Department of Astronomy, University of Florida, P.O. Box 112055, Gainesville, FL 32611, USA</footnote>

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

<footnote n="13">¹³ S. N. Bose National Centre for Basic Sciences, Sector-III, Salt Lake, Kolkata 700106, India</footnote>

<footnote n="14">¹⁴ Astronomy Department, Universidad de Chile, Camino El Observatorio 1515, Las Condes, Santiago, Chile</footnote>

<footnote n="15">¹⁵ Departament de Física Quàntica i Astrofísica (FQA), Universitat de Barcelona (UB), Martí i Franquès 1, 08028 Barcelona, Catalonia, Spain</footnote>

<footnote n="16">¹⁶ Institut de Ciències del Cosmos (ICCUB), Universitat de Barcelona, Martí i Franquès, 1, 08028 Barcelona, Catalonia, Spain</footnote>

<footnote n="17">¹⁷ Institut d'Estudis Espacials de Catalunya (IEEC), Gran Capità, 2-4, 08034 Barcelona, Catalonia, Spain</footnote>

<footnote n="18">¹⁸ Instituto Argentino de Radioastronomía (CCT-La Plata, CONICET; CICPBA), C.C. No. 5, 1894, Villa Elisa, Buenos Aires, Argentina</footnote>

<footnote n="19">¹⁹ Joint Alma Observatory (JAO), Alonso de Córdova 3107, Vitacura, Santiago, Chile</footnote>

<footnote n="20">²⁰ School of Physics and Astronomy, Yunnan University, Kunming, 650091, PR China</footnote>

<footnote n="21">²¹ Institute of Astronomy and Department of Physics, National Tsing Hua University, Hsinchu 30013, Taiwan</footnote>

<footnote n="22">²² Max-Planck-Institut für Radioastronomie, Auf dem Hügel 69, 53121 Bonn, Germany</footnote>
# Appendix A: Filamentary identification with FilFinder

Here we describe the procedure to identify the main filamentary structures presented in Sect. 3.2 using FilFinder (Koch & Rosolowsky 2015).

For this approach we used the moment 0 map of the extracted N₂H⁺ isolated components that present a S/N ≥ 5. To estimate the moment 0 map we used the moment task from the SpectralCube Python package, within the velocity range of -31.5 km s⁻¹ to 0 km s⁻¹. As part of the pre-processing of the moment 0 map before the filamentary detection, we decreased the contrast in the image by using the preprocess_image task and its argument flatten_percent set to 90. Now, in order to indicate to FilFinder the area in which to identify filaments, we used the subtask create_mask with the following parameters: glob_thresh: 4.5 K km s⁻¹, size_thresh: 0.25 pc², smooth_size: 0.12 pc, border_masking: False, fill_hole_size: 0.013 pc². The resulting mask is presented in Fig. A.1 with a white contour.

Then, we obtained the skeletons of the mask by using medskel. The derived structures are presented with red and green lines in Fig. A.1. Given we are only interested in the large-scale filaments, we used analyze_skeletons in order to “prune” the small-scale structures. For this pruning we used branch_thresh: 0.3 pc, prune_criteria: 'length', max_prune_iter: 0. This approach results in removing the small filaments (red lines in Fig. A.1) from the original skeleton and to obtain the main filamentary structure in G353 (green lines in Fig. A.1).

<figure><img src="image_17.png" /><figcaption>Fig. A.1: FilFinder filamentary identification. The background indicates the moment 0 map of the extracted N₂H⁺ isolated components. The white contour shows the area where FilFinder identifies multiple filamentary structures (red and green lines). We remove the small-scale structures (in red) by “pruning” the skeleton structure from medskel, obtaining the main filaments of G353. We represent these filaments with green lines.</figcaption></figure>

# Appendix B: Examples of the isolated components fitting

In Sect. 3.1 we decomposed the multiple isolated component emission using PySpecKit. In Fig. B.1, we present the results of the Gaussian fitting for the high S/N spectra shown in Fig. 2 (panels c, e, and f).
<figure><img src="image_18.png" /></figure>

Fig. B.1: Gaussian velocity fits of the extracted N<sub>2</sub>H<sup>+</sup> isolated components. In black we show the high S/N isolated components from panels c), e), and f) in Fig. 2. The individual Gaussian components and the obtained model are represented with dashed blue and solid red lines, respectively. On the right side of each panel we indicate the peak intensity (I), the velocity centroid (v), and velocity dispersion (σ) of each Gaussian component. The notations 1<sup>st</sup>, 2<sup>nd</sup>, and 3<sup>rd</sup> indicate the Gaussian velocity components from left to right.

**Appendix C: DCN and N<sub>2</sub>H<sup>+</sup> derived core velocity**

In Table C.1 we provide the 1.3 mm core velocities obtained from
DCN & N<sub>2</sub>H<sup>+</sup> data (see Sect. 4.2), complementing the published
DCN core velocities catalog from Cunningham et al. (2023).
In the last column we indicate the number of N<sub>2</sub>H<sup>+</sup> velocity
components detected in these cores.

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

Fig. C.1: DCN and N<sub>2</sub>H<sup>+</sup> normalized mean spectrum of cores 35 (top) and 46 (bottom). We show the multiple N<sub>2</sub>H<sup>+</sup> isolated velocity components with blue, red, and green colors. We present the DCN emission in black. We see a match between the DCN emission and one of the N<sub>2</sub>H<sup>+</sup> velocity components. We determine the N<sub>2</sub>H<sup>+</sup> velocity for 11 cores with no DCN velocity fits. These are listed in Table C.1.
<table><caption>Table C.1: 1.3 mm core catalog of DCN & N₂H⁺ velocities.</caption><thead><tr><th>Core Number</th><th>RA [°]</th><th>DEC [°]</th><th>F<sub>A</sub> ["]</th><th>F<sub>B</sub> ["]</th><th>PA [°]</th><th>Mass [M<sub>☉</sub>]</th><th>V<sub>LSR</sub> [km s<sup>-1</sup>]</th><th>Type</th><th>Number of N<sub>2</sub>H<sup>+</sup> components</th></tr></thead><tbody><tr><td>2</td><td>262.6165032</td><td>-34.6955865</td><td>1.98</td><td>1.59</td><td>64.00</td><td>20.7</td><td>-20.45±0.065</td><td>DCN, Single</td><td>2</td></tr><tr><td>3</td><td>262.6184156</td><td>-34.6965240</td><td>2.59</td><td>1.79</td><td>146.00</td><td>9.4</td><td>-16.48±0.040</td><td>DCN, Single</td><td>2</td></tr><tr><td>4</td><td>262.6103159</td><td>-34.6932659</td><td>1.56</td><td>1.46</td><td>104.00</td><td>5.2</td><td>-16.53±0.155</td><td>DCN, Single</td><td>3</td></tr><tr><td>5</td><td>262.6101515</td><td>-34.6960014</td><td>2.03</td><td>1.75</td><td>79.00</td><td>16.0</td><td>-16.94±0.087</td><td>DCN, Single</td><td>3</td></tr><tr><td>6</td><td>262.6049155</td><td>-34.6934384</td><td>1.60</td><td>1.48</td><td>129.00</td><td>4.9</td><td>-18.68±0.126</td><td>DCN, Single</td><td>2</td></tr><tr><td>7</td><td>262.6137738</td><td>-34.6947298</td><td>1.63</td><td>1.27</td><td>80.99</td><td>6.0</td><td>-19.79</td><td>DCN & N<sub>2</sub>H<sup>+</sup></td><td>3</td></tr><tr><td>8</td><td>262.6039531</td><td>-34.6936374</td><td>2.26</td><td>1.56</td><td>97.32</td><td>6.6</td><td>—</td><td>—</td><td>—</td></tr><tr><td>9</td><td>262.6192359</td><td>-34.6903650</td><td>2.69</td><td>1.96</td><td>62.96</td><td>3.7</td><td>-16.09</td><td>DCN & N<sub>2</sub>H<sup>+</sup></td><td>3</td></tr><tr><td>11</td><td>262.6243189</td><td>-34.6880780</td><td>2.12</td><td>1.95</td><td>172.70</td><td>2.8</td><td>-16.09</td><td>DCN & N<sub>2</sub>H<sup>+</sup></td><td>3</td></tr><tr><td>12</td><td>262.6072148</td><td>-34.6969795</td><td>2.98</td><td>2.02</td><td>86.00</td><td>10.3</td><td>-17.83±0.019</td><td>DCN, Single</td><td>3</td></tr><tr><td>13</td><td>262.6078228</td><td>-34.6996836</td><td>1.91</td><td>1.48</td><td>153.90</td><td>2.7</td><td>—</td><td>—</td><td>—</td></tr><tr><td>14</td><td>262.6147937</td><td>-34.6946762</td><td>2.00</td><td>1.58</td><td>124.00</td><td>6.2</td><td>-14.36±0.090</td><td>DCN, Single</td><td>3</td></tr><tr><td>15</td><td>262.6107433</td><td>-34.6964412</td><td>1.96</td><td>1.59</td><td>93.00</td><td>6.6</td><td>-14.78±0.087</td><td>DCN, Complex</td><td>—</td></tr><tr><td>16</td><td>262.6215941</td><td>-34.6989408</td><td>2.68</td><td>2.49</td><td>19.37</td><td>1.5</td><td>—</td><td>—</td><td>—</td></tr><tr><td>17</td><td>262.5954514</td><td>-34.6916168</td><td>1.57</td><td>1.42</td><td>87.50</td><td>0.9</td><td>—</td><td>—</td><td>—</td></tr><tr><td>18</td><td>262.5927434</td><td>-34.7052494</td><td>1.97</td><td>1.56</td><td>89.26</td><td>0.8</td><td>—</td><td>—</td><td>—</td></tr><tr><td>19</td><td>262.6064012</td><td>-34.7019756</td><td>1.55</td><td>1.20</td><td>57.74</td><td>0.5</td><td>—</td><td>—</td><td>—</td></tr><tr><td>20</td><td>262.6111096</td><td>-34.6932787</td><td>1.67</td><td>1.63</td><td>178.00</td><td>1.3</td><td>—</td><td>—</td><td>—</td></tr><tr><td>21</td><td>262.6131910</td><td>-34.6939495</td><td>2.10</td><td>1.48</td><td>108.00</td><td>2.4</td><td>-18.56±0.063</td><td>DCN, Single</td><td>3</td></tr><tr><td>22</td><td>262.6118441</td><td>-34.6946150</td><td>1.89</td><td>1.59</td><td>96.65</td><td>1.7</td><td>-16.76</td><td>DCN & N<sub>2</sub>H<sup>+</sup></td><td>3</td></tr><tr><td>23</td><td>262.6028175</td><td>-34.6925438</td><td>1.84</td><td>1.31</td><td>112.70</td><td>0.8</td><td>—</td><td>—</td><td>—</td></tr><tr><td>24</td><td>262.6198349</td><td>-34.6960383</td><td>1.88</td><td>1.72</td><td>76.00</td><td>0.8</td><td>-19.50±0.089</td><td>DCN, Single</td><td>2</td></tr><tr><td>25</td><td>262.6155222</td><td>-34.6952591</td><td>1.87</td><td>1.20</td><td>137.30</td><td>2.5</td><td>-18.78</td><td>DCN & N<sub>2</sub>H<sup>+</sup></td><td>2</td></tr><tr><td>26</td><td>262.6143434</td><td>-34.6917027</td><td>1.49</td><td>1.24</td><td>137.30</td><td>0.8</td><td>—</td><td>—</td><td>—</td></tr><tr><td>27</td><td>262.6000802</td><td>-34.6910324</td><td>3.39</td><td>2.51</td><td>48.16</td><td>1.8</td><td>—</td><td>—</td><td>—</td></tr><tr><td>28</td><td>262.6253977</td><td>-34.6999713</td><td>2.56</td><td>1.87</td><td>39.31</td><td>0.8</td><td>—</td><td>—</td><td>—</td></tr><tr><td>29</td><td>262.6133074</td><td>-34.6919187</td><td>1.67</td><td>1.26</td><td>76.67</td><td>0.7</td><td>—</td><td>—</td><td>—</td></tr><tr><td>30</td><td>262.6114686</td><td>-34.6962602</td><td>1.52</td><td>1.41</td><td>30.00</td><td>2.3</td><td>-12.81±0.074</td><td>DCN, Single</td><td>2</td></tr><tr><td>31</td><td>262.6096651</td><td>-34.6925680</td><td>1.72</td><td>1.34</td><td>119.10</td><td>0.7</td><td>—</td><td>—</td><td>—</td></tr><tr><td>32</td><td>262.6094126</td><td>-34.6910985</td><td>2.44</td><td>2.02</td><td>59.22</td><td>2.0</td><td>—</td><td>—</td><td>—</td></tr><tr><td>33</td><td>262.6287106</td><td>-34.6862068</td><td>2.32</td><td>2.00</td><td>16.64</td><td>1.7</td><td>—</td><td>—</td><td>—</td></tr><tr><td>34</td><td>262.5982914</td><td>-34.6919006</td><td>1.81</td><td>1.36</td><td>139.00</td><td>0.8</td><td>-18.77±0.080</td><td>DCN, Single</td><td>2</td></tr><tr><td>35</td><td>262.6142011</td><td>-34.6940134</td><td>2.31</td><td>1.81</td><td>111.00</td><td>3.4</td><td>-16.09</td><td>DCN & N<sub>2</sub>H<sup>+</sup></td><td>3</td></tr><tr><td>36</td><td>262.6202758</td><td>-34.7001995</td><td>2.19</td><td>1.88</td><td>141.00</td><td>0.8</td><td>—</td><td>—</td><td>—</td></tr><tr><td>37</td><td>262.6010398</td><td>-34.6950114</td><td>1.74</td><td>1.33</td><td>75.00</td><td>0.8</td><td>-17.83±0.053</td><td>DCN, Single</td><td>3</td></tr><tr><td>38</td><td>262.5971034</td><td>-34.6920396</td><td>2.01</td><td>1.61</td><td>104.80</td><td>0.7</td><td>—</td><td>—</td><td>—</td></tr><tr><td>39</td><td>262.6054437</td><td>-34.6963773</td><td>2.05</td><td>1.88</td><td>158.50</td><td>1.4</td><td>-16.43</td><td>DCN & N<sub>2</sub>H<sup>+</sup></td><td>3</td></tr><tr><td>40</td><td>262.5917454</td><td>-34.6897316</td><td>1.85</td><td>1.50</td><td>92.93</td><td>0.6</td><td>—</td><td>—</td><td>—</td></tr><tr><td>41</td><td>262.6095777</td><td>-34.6983259</td><td>2.23</td><td>1.83</td><td>72.51</td><td>1.3</td><td>—</td><td>—</td><td>—</td></tr><tr><td>42</td><td>262.5975992</td><td>-34.6876666</td><td>1.51</td><td>1.22</td><td>24.46</td><td>0.3</td><td>—</td><td>—</td><td>—</td></tr><tr><td>43</td><td>262.6035166</td><td>-34.6966807</td><td>2.48</td><td>1.78</td><td>95.00</td><td>1.2</td><td>-17.33±0.062</td><td>DCN, Single</td><td>3</td></tr><tr><td>44</td><td>262.6030115</td><td>-34.6956424</td><td>1.76</td><td>1.55</td><td>96.00</td><td>0.5</td><td>-16.99±0.091</td><td>DCN, Single</td><td>3</td></tr><tr><td>45</td><td>262.6143008</td><td>-34.6909376</td><td>3.40</td><td>2.75</td><td>94.93</td><td>2.6</td><td>—</td><td>—</td><td>—</td></tr><tr><td>46</td><td>262.6187648</td><td>-34.6912377</td><td>3.10</td><td>2.15</td><td>40.67</td><td>0.9</td><td>-15.08</td><td>DCN & N<sub>2</sub>H<sup>+</sup></td><td>3</td></tr><tr><td>47</td><td>262.6178453</td><td>-34.6919943</td><td>3.09</td><td>2.48</td><td>45.43</td><td>1.6</td><td>—</td><td>—</td><td>—</td></tr></tbody></table>

**Notes.** Velocities of the 1.3 mm continuum-derived cores. The column “Type” indicates if the core velocity is determined by a single or complex DCN spectra (Cunningham et al. 2023), or using both DCN & N₂H⁺ data (this work). In the last column indicate the number of N₂H⁺ velocity components inside these cores. For completeness we include the properties of the 1.3 mm cores with no velocity determinations. These cores present a “-” mark in the last three columns. The core catalog from Louvet et al. (2024) does not consider core “1” nor core “10” given their selection criteria reject cores with free-free contamination or cores undetected at 3 mm.
## Appendix D: V-shaped structures

In Sect. 5.3, we characterized the most prominent V-shaped structure we detect in Fig. 9. We repeated this process for the other eight different V-shaped structures, including the linear fits to the velocity gradients. In Fig. D.1, we indicate the V-shapes location in PV space with dark points, arrows, and their IDs. In Fig. D.2, we present individual close-ups for each V-shape. In Table D.1 we list their VGs, timescales, $H_2$ masses, and mass accretion rates.

Here we list a few clarifications due to projection effects seen in these V-shaped structures:

- In position-position space, only V-shape “B” presents a core within a ~beam size from its apex. This continuum core, core “7”, has a mass of ~6 $M_\odot$.

- For V-shape “A”, the 1.3 mm core with DCN single velocity component, located at the apex of this V-shape, is not spatially related with it.

- In V-shapes “G” and “H” we see the same 1.3 mm cores with $N_2H^+$ velocities. These V-shapes are not the same distribution. They are overlapped in PV space and spatially separated by ~10".

- We improve the clarity of V-shape “B” by rotating the data in PP space by 80° counter-clockwise. We apply this process for V-shapes “E”, “G”, and “H” with an angle of 33° clockwise.

- V-shapes G and H overlap in PV space but these are structures spatially separated.

<figure><img src="image_20.png" /><figcaption>Fig. D.1: V-shapes location in PV space. We highlight the V-shapes listed in Table D.1 with black points and indicate them with red arrows and their ID. The core velocities and the $N_2H^+$ velocity distributions follow the same definitions from the top right panel in Fig. 9. V-shapes G and H overlap in PV space but these structures are spatially separated (left panel).</figcaption></figure>

<table><caption>Table D.1: Characterized V-shaped structures.</caption><thead><tr><th>V-shape ID</th><th><i>l</i><br>[°]</th><th><i>b</i><br>[°]</th><th><i>M</i>(H<sub>2</sub>)<br>[<i>M</i><sub>☉</sub>]</th><th>Upper / lower VG<br>[km s<sup>-1</sup> pc<sup>-1</sup>]</th><th>Mean VG<br>[km s<sup>-1</sup> pc<sup>-1</sup>]</th><th>Upper / lower <i>t</i><sub>VG</sub><br>[kyr]</th><th><i>t</i><sub>VG mean</sub><br>[kyr]</th><th><i>M</i><sub>in</sub>(H<sub>2</sub>)<br>[10<sup>-4</sup> <i>M</i><sub>☉</sub> yr<sup>-1</sup>]</th></tr></thead><tbody><tr><td>A</td><td>353.3981</td><td>-0.3506</td><td>8.13</td><td>20.85 / 17.04</td><td>18.95</td><td>46.89 / 57.39</td><td>52.14</td><td>1.56</td></tr><tr><td>B</td><td>353.4127</td><td>-0.3632</td><td>13.29</td><td>25.34 / 17.22</td><td>21.28</td><td>38.59 / 56.77</td><td>47.68</td><td>2.79</td></tr><tr><td>C</td><td>353.4135</td><td>-0.3657</td><td>53.02</td><td>17.69 / 13.26</td><td>15.48</td><td>55.28 / 73.75</td><td>64.52</td><td>8.22</td></tr><tr><td>D</td><td>353.4133</td><td>-0.3727</td><td>7.96</td><td>12.85 / 21.22</td><td>17.15</td><td>76.07 / 46.07</td><td>61.07</td><td>1.31</td></tr><tr><td>E</td><td>353.4096</td><td>-0.3521</td><td>6.27</td><td>12.47 / 3.63</td><td>8.05</td><td>78.38 / 269.57</td><td>173.98</td><td>0.36</td></tr><tr><td>F</td><td>353.4128</td><td>-0.3604</td><td>27.24</td><td>15.75 / 15.61</td><td>15.68</td><td>62.09 / 62.64</td><td>62.37</td><td>4.37</td></tr><tr><td>G</td><td>353.4110</td><td>-0.3630</td><td>3.24</td><td>22.59 / 24.46</td><td>23.53</td><td>43.29 / 39.97</td><td>41.63</td><td>0.79</td></tr><tr><td>H</td><td>353.4140</td><td>-0.3627</td><td>31.46</td><td>21.28 / 39.68</td><td>30.48</td><td>45.94 / 24.64</td><td>35.29</td><td>8.91</td></tr><tr><td>I</td><td>353.4091</td><td>-0.3595</td><td>16.13</td><td>22.34 / 11.22</td><td>16.78</td><td>43.77 / 87.17</td><td>65.47</td><td>2.46</td></tr></tbody></table>

**Notes.** Properties of the nine, well characterized, V-shaped structures identified in our $N_2H^+$ data. The coordinates indicate the position of the velocity apex of each V-shape. V-shape “C” represents the structure analyzed in Sect. 5.3. We show the PV distribution of these structures in Fig. D.2.
<figure><img src="image_21.png" /><figcaption>Fig. D.2: V-shaped structures listed in Table D.1, with the exception of “C” shown in Fig. 10. We indicate the “V-shape ID” from Table D.1 at the top/bottom left corner of each plot. The colors of the distributions, DCN and DCN & N₂H⁺ derived core velocities, and beam size follow the same color and marker convention from Fig. 9. See Appendix D for clarifications regarding projection effects on these diagrams.</figcaption></figure>
## Appendix E: SiO Intensity-weighted position-velocity diagram

To create the SiO intensity-weighted PV diagram, first, we removed most of the noisy spectra by considering data with S/N ≥ 2.5. Then, we estimated the integrated intensity and velocity centroid at each pixel. We found improvements in our cleaning by using only spectra with integrated intensity ≥ 4 K km s⁻¹. Using the coordinate, integrated-intensity, and velocity centroid of each spectrum, we create the SiO intensity-weighted PV diagrams we show in E.1.

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

Fig. E.1: ALMA-IMF 12 m SiO equivalent of Fig. 9 using data from Cunningham et al. (2023). For the cores, we use the same marker and color convention from Fig. 9. With filled blue, red, and ‘red+blue’ circles we represent the SiO outflow candidates (Towner et al. 2024). With red arrows we indicate a VG = 400 km s⁻¹ pc⁻¹ corresponding to a timescale $t_{VG} = 2.5$ kyr. The velocity range (ΔV) covered by the SiO emission is ~80km s⁻¹, about 10 times the velocity range traced by N₂H⁺. This velocity difference suggests that SiO is tracing processes (outflows) ~100 times more energetic ($e_k = \Delta V/2$) than N₂H⁺ (possibly infall).
## Appendix F: G353 power law density profile

Here we provide the derivation of the density profile used for the gravitationally collapsing sphere. We assumed a power law density profile defined as:

$$ \rho(r) = \rho_0 \left(\frac{r}{\text{pc}}\right)^{-\gamma}, \qquad (\text{F.1}) $$

in which $\gamma = 5.65$ provides a good fit to the edges of the PV distribution seen in Fig. 8.

To determine the value of $\rho_0$, we integrated this expression in a sphere (Eq. (F.2)), with $r_{min} < r < 0.5$ pc. Based on different tests, probing total masses from $50 - 10^3 M_\odot$ and $\gamma$ values from 2 - 6, we set the total mass of the sphere to 150 $M_\odot$. We defined $r_{min} \sim 0.007$ pc which corresponds to the pixel size of the $N_2H^+$ data at a distance of 2 kpc.

$$ M_{enc}(r = 0.5 \text{ pc}) = 4\pi\rho_0 \int_{r_{min}}^{0.5 \text{ pc}} \left(\frac{r}{\text{pc}}\right)^{-\gamma} r^2 dr \qquad (\text{F.2}) $$

$$ = 4\pi\rho_0 \frac{r^{3-\gamma}}{3-\gamma} \Bigg|_{r=r_{min}}^{r=0.5 \text{ pc}} \text{pc}^{\gamma} \qquad (\text{F.3}) $$

$$ = \frac{4\pi\rho_0}{3-\gamma} (0.5^{3-\gamma} - r_{min}^{3-\gamma}) \text{pc}^{\gamma}, \qquad (\text{F.4}) $$

in which $M_{enc}(r = 0.5 \text{ pc}) = 150 M_\odot$, $r_{min} = 7 \times 10^{-3}$ pc, and $\gamma = 5.65$, from Eq. (F.4), we obtained:

$$ \rightarrow \rho_0 = 6.1 \times 10^{-5} \frac{M_\odot}{\text{pc}^3}. \qquad (\text{F.5}) $$