© The Authors 2024

1 Introduction

In 1955, Edwin Salpeter reported that the high-mass tail of the initial mass function of stars, the IMF, could be represented by a power law of the form $\frac{\text{d}N}{ \text{d}M}\propto M^{-\alpha }$ with α = 2.35. Since then, numerous studies reported its universality (see the review by Bastian et al. 2010) and studied its origin (e.g. Hennebelle & Chabrier 2008). To investigate the origin of the IMF, astronomers studied the fragmentation of molecular clouds in the form of small substructures (≲0.03 pc; the cores) and the link between the core mass function (CMF) and the IMF (e.g. Motte et al. 1998; Fiorellino et al. 2021). Because of observational limitations, these studies focused mostly on nearby clouds (≈140–400 pc) until 2018. These clouds mainly form low- and intermediate-mass stars. All studies reported CMF slopes at the high-mass end that were compatible with the slope of the IMF, suggesting that the IMF directly inherits its shape from the CMF (e.g. Alves et al. 2007; Könyves et al. 2010).

In 2018, taking advantage of the angular resolution and sensitivity of the Atacama Large Millimetre/submillimetre Array (ALMA), Motte et al. (2018) studied the high-mass protocluster W43-MM1 at a distance of ≃5.5 kpc. They reported a CMF with a high-mass tail that was flatter than the canonical IMF, with α=1.96±0.13. Since then, other teams have reported similar results in other high-mass protoclusters, but often used singlepointing observations, thus narrow fields that do not image parsec-scale clouds (e.g. Liu et al. 2018; Cheng et al. 2018; Kong 2019; O’Neill et al. 2021). These results cast doubt on the universality of the IMF, on a direct link between the CMF and the IMF, or on both. However, the comparison of core samples between these studies, which are captured at different evolutionary stages, becomes complex or impossible when the observations are carried out at disparate physical resolutions and/or sensitivities. Thus, it became imperative to obtain a sample of high-mass protoclusters that are representative of the Milky Way and were observed homogeneously.

To directly address this, we present in this article the results of the ALMA-IMF Large Program (PIs: F. Motte, A. Ginsburg, F. Louvet, P. Sanhueza). ALMA-IMF observed 15 high-mass protoclusters at comparable sensitivity and physical resolution. The main driver of the protocluster selection is a mass criterion that first chose protoclusters with inner parts of at least 500 M_⊙, as described in detail in Sect. 2.1 of Motte et al. (2022).

A second criterion was the distance of the protoclusters. Not too close by (>2 kpc) to permit imaging a large area with mosaic observations, and not too far (<6 kpc) to prevent excessive integration time. Our final sample consists of 15 protoclusters that were observed with a continuum sensitivity below 60 mK at 1.3 mm, and at similar physical resolution (≈2 kau).

In this article, we determine and analyse the combined CMF of all 15 protoclusters from ALMA-IMF. We show the 1.3 and 3 mm continuum images in Sect. 2; they were first presented in Ginsburg et al. (2022). We describe the source extraction in Sect. 3. We describe the source selection and core mass derivation in Sect. 4, paying special attention to free-free emission and to the core temperatures. We present the global CMF for all 15 fields. We discuss our results regarding the current knowledge about the CMF in Sect. 6, and we summarise our conclusions in Sect. 7.

2 Observations

The data were taken from the ALMA-IMF large programme (Project ID: 2017.1.01355.L), entitled ALMA transforms our view of the origin of stellar masses. Motte et al. (2022), hereafter referred to as Paper I, described the project, its choice of targets, datasets, objectives, and first results. Ginsburg et al. (2022), hereafter referred to as Paper II, described the reduction pipeline for the continuum emission maps at 1.3 mm and 3 mm in the ALMA-IMF sample. Paper II presented two versions of the continuum maps: the BSENS and the CLEANEST. The BSENS maps were constructed using the full range of the continuum spectral windows, which may include significant molecular emission, especially toward hot cores. The CLEANEST maps were constructed by flagging all channels contaminated by molecular emission (see Fig. 3 in Paper II). Because we wish to assemble the most reliable sample of cores, we used the CLEANEST continuum maps, which are publicly available on the ALMA-IMF website¹. Table 1 lists the 15 regions investigated by ALMA-IMF together with the observational characteristics taken from Papers I and II that are relevant to the present paper. We refer to Paper I for an in-depth description of the 15 regions.

3 Source extraction

3.1 Extraction tools

To extract sources in the continuum emission maps, we used two extraction tools: getsf (v210414, Men’shchikov 2021b) and GExt2D (v210603, Bontemps et al., in prep.). This approach permitted us to assess the biases inherent to our source extraction algorithms.

The multi-scale multi-wavelength source and filament extraction method getsf supersedes the previous algorithms getsources, getfilaments, and getimages (Men’shchikov et al. 2012; Men’shchikov 2013, 2017, respectively). The getsf method separates the structural components of sources, filaments, and background in the observations into independent images of the different components. Source and filament extraction are performed with their respective images, where the contributions of the other components are largely removed. For more details on getsf, illustrations of its applications, quantitative extraction performance, and comparisons with getsources, we refer to Men’shchikov (2021a,b).

GExt2D is based on the same idea as the cutex tool (Molinari et al. 2011). As a first step, it uses the second derivative of the maps to identify the locations of maximum curvature as indicative of compact sources. It starts from the strongest curvatures down to fainter and shallower fluctuations close to those expected from measured local noise. The second step involves source boundary identification and removal of background emission as the interpolated values in the map inside the boundaries. The third step corresponds to an automated Gaussian fitting at each location of compact sources after subtracting the derived background.

Fig. 1 Scatter plots of the source size (black dots) extracted with getsf before (left) and after (right) the maps were smoothed to a homogeneous spatial resolution of 2700 au. The blue triangles display the beam sizes. The orange diamonds display the mean sizes of sources. In the left panel, the y-axis is the linear scale (beam × distance) of the maps in Band 6, and the orange line shows the fit to the mean source size distribution for each field. In the right panel, the orange line shows the mean core size of all 15 protoclusters, and the red crosses illustrate the cores exceeding 25 M⊙.

3.2 Source extraction at native angular resolutions

With both getsf and GExt2D, the source extraction strategy has the same two steps. Firstly, we detected source candidates in the continuum images at 1.3 mm (Band 6) without primary-beam correction. These images have a homogeneous noise level in the entire field. This choice permits the detection of sources with a constant signal-to-noise ratio threshold. Secondly, we measured the source fluxes and sizes in the continuum maps at 1.3 mm and 3 mm (Band 3) after primary-beam correction. In this way, the fluxes we obtained account for any attenuation resulting from being farther from the phase centre of the fields. This strategy replicates that applied in Pouteau et al. (2022), hereafter referred to as Paper III, dedicated to the W43-MM2 and W43-MM3 protoclusters. Considering all 15 protoclusters together, we extracted 820 sources with getsf, and 930 with GExt2D. GExt2D detects the sources in the observed image, while getsf decomposes the image first and requires detecting a source in successive scales to validate its detection. Hence, getsf is more conservative and extracts fewer sources than GExt2D.

The distances of the protoclusters range from 2.5 to 5.5 kpc (see Table 1 in Paper I). The ALMA-IMF large programme aimed to image the 15 protoclusters with the same spatial resolution of 2100 au, requesting the corresponding angular resolution for each field. Despite this effort, the available antenna configurations at the observatory led to variations in the spacial resolutions by up to a factor of two (see Table 1 and Paper I). These variations in physical scale impact the size and mass of the sources that can be extracted. Louvet et al. (2021) performed convolution tests on projected density cubes obtained from hydro-dynamical simulations and on column density maps created from Herschel observations and showed that the sizes and integrated fluxes of the extracted sources closely follow the size of the linear scale. The left panel of Fig. 1 shows a linear dependence between the mean source size and the linear resolution. To obtain a homogeneous sample of sources, the data must be smoothed to a matched spatial resolution.

3.3 Source extraction at matched physical resolution

To extract a spatially homogeneous core sample, we smoothed the 1.3 and 3 mm maps to the lowest native resolution of ≃2700 au. We convolved each map with a Gaussian kernel whose width, θ_conv, was calculated as $${\theta }_{\text{conv}}=\sqrt{{\left(\frac{2700\text{au}}{D}\right)}^2-{\theta }_{\text{native}}^2},$$(1)

where D is the distance to the target expressed in parsecs, and θnative is the native angular resolution of the map. Table 2 lists the angular resolution of each field after smoothing. We present the smoothed continuum emission maps at 1.3 mm in Fig. 2.

We adopted the same strategy as explained in Sect. 3.2 to extract sources in the smoothed maps. In total, we retrieved about 680 sources with getsf and 1020 with GExt2D ². Table 2 lists the number of sources extracted per field with both algorithms, and Fig. 2 shows the sources extracted with getsf³. The right panel of Fig. 1 shows the sizes of the sources we extracted with getsf for each field. The mean size is similar for each field with L = 3 840 ± 270 au⁴.

4 Results

4.1 Getsf versus GExt2D extraction methods

In total, we retrieved 677 sources with getsf, and 1020 with GExt2D. The agreement between these extractions is good because ≈80% of the sources found by getsf are also found by GExt2D (Table 2). To compare the source flux measurements, we constructed the source flux functions (SFFs) and complementary cumulative SFFs for both algorithms (Fig. 3). The first function shows that the source fluxes range between 10 mJy and 10Jy for both algorithms, with a similar median value at ≃12mJy. It also shows that the additional sources found by GExt2D are low- and intermediate-flux sources from ≈30mJy to 200 mJy. The complementary SFF permitted us to estimate the slope of the high-flux end of the SFF. We fitted the cumulative SFF with a power law from the median value of the samples, that is, 12 mJy, up to the flux at which the number of cores is below ten⁵. The fitted high-flux slopes are similar for both catalogues, with ζ′ ≃ 0.85. These exponents, adjusted onto the complementary cumulative SFFs, correspond to power laws with exponents ζ = ζ′ + 1 ≃ 1.85. This analysis proves that the overall statistics of the source flux are independent of the two extraction algorithms we employed. In the following, we choose to focus on the extractions performed by getsf since it is more conservative than GExt2D, which is confirmed by the large fraction of getsf sources confirmed with GExt2D. When we compare the sources extracted at the native angular resolutions (see Sect. 3.2) to the sources extracted in the smoothed maps (see Sect. 3.3), the agreement between the getsf catalogues is better than that between GExt2D catalogues. About 95% of the sources extracted in the smoothed maps were also extracted at the native angular resolutions with the getsf algorithm. This fraction drops to ≃75% in the case of the GExt2D extractions.

Fig. 2 Fifteen ALMA-IMF protoclusters as traced by their continuum emission at 1.3 mm at a matched spatial resolution of 2700 au. The ellipses locate all the cores found by getsf. The name of the field is indicated in the top left corner of each panel, and the beam size is shown in the bottom left corner. For each protocluster, the centre coordinates are those specified in Table 1. We show the remaining fields in Fig. A.1.

Fig. 3 Source flux distributions for the getsf and GExt2D source extractions. The left panel displays the SFFs of the sources extracted by getsf in grey and by GExt2D in blue. The vertical blue and grey lines indicate the median values for each source catalogues, which differ by ≃5%. The right panel displays the complementary cumulative SFFs from the getsf source extraction (in grey) and from the GExt2D source extraction (in blue). The grey and blue lines display the best fits from linear regressions on the complementary cumulative SFFs from the median value at 12 mJy up to the flux at which the number of cores is smaller than 10. The corresponding ζ power-law indexes are indicated in the top right corner.

4.2 From source flux to core mass

4.2.1 Removing contaminated sources

The integrated fluxes comprise thermal dust emission and potentially free-free emission for the more evolved sources. In the following, we filter the sources that are arguably contaminated by free-free emission. To do this, we used the integrated flux measurements at 1.3 and 3 mm, when available, to compute their spectral index. To correct for the source size differences between 1.3 and 3 mm extractions, we linearly adjusted the integrated fluxes at 3mm: $S_{3 \text{mm}}^{\text{int}*}=S_{3 \text{mm}}^{\text{int}}\times \left(\frac{{\Theta }_A^{1.3 \text{mm}}\times {\Theta }_B^{1.3 \text{mm}}}{{\Theta }_A^{3 \text{mm}}\times {\Theta }_B^{3 \text{mm}}}\right)$, where Θ_A and Θ_B are the major and minor axes of the sources at 1.3 or 3 mm, respectively. This operation corrected the 3 mm fluxes by 17% on average. This flux rescaling replicates the method we used in Paper III. This rescaling works for an optically thick emission, and an optically thin emission for an isothermal protostellar envelope with a density profile ρ(r) ∝ r⁻², where r is the radius of the envelope. We then computed the spectral index as $$\gamma =\frac{\log \left(S_{1.3\text{mm}}^{\text{int}}/S_{3\text{mm}}^{\text{int}\ast }\right)}{\log \left(v_{1.3\text{mm}}/v_{3\text{mm}}\right)},$$(2)

where v_1.3mm and v_3mm are the central frequencies of the ALMA Band 6 and Band 3, respectively (see Table 1). We removed 68 sources with γ < 2, which are presumably contaminated by free-free emission. These sources are represented by the pink ellipses in Fig. 4.

To investigate the ≃500 sources that lack measurements at 3 mm, we built the spectral index map for each field using Eq. (2) on a pixel-by-pixel basis by replacing $S_{1.3 \text{mm}}^{\text{int}}\text{ and }{S}_{3 \text{mm}}^{\text{int}}$ with the pixel intensities. We point out that the angular resolutions and pixel sizes of the images at 1.3 mm equal that of the images at 3 mm. We show the spectral index maps in Fig. 4. The red areas (γ > 2) represent pixels with thermal dust emission, those in green (0 < γ < 2) pixels contaminated by free-free emission, and those in blue (γ < 0) pixels dominated by free-free emission. We manually selected all the sources lying over free-free emission areas (γ < 2; see yellow ellipses in Fig. 4). For the sources that are undetected at 3 mm, we estimated their flux, $S_{3 \text{mm}}^{\text{upper}}$, using the integrated flux in the aperture corresponding to the footprint of the source detected at 1.3 mm. We emphasise that although getsf did not detect these 1.3 mm sources as compact sources at 3 mm, they do have a background emission above 3σ at 3 mm by definition of the spectral index maps. We used this upper limit to compute an upper limit on their spectral index as $$\gamma =\frac{\log \left(S_{1.3\text{mm}}^{\text{int}}/S_{3\text{mm}}^{\text{upper}}\right)}{\log \left(v_{1.3\text{mm}}/v_{3\text{mm}}\right)},$$(3)

and we rejected the 16 sources with a spectral index below one, those for which free-free emissions could substantially contribute to the emission.

This filtering rejected 84 of the 677 sources extracted by getsf. The degree of rejection depended on the evolutionary stages of the regions: ≃0.5% in the young regions, ≃10% in the intermediate regions, and ≃25% in the evolved regions. Table 3 lists the number of sources that were rejected when these two levels of filtering were applied for each field.

Fig. 4 Maps of the spectral index in the 15 ALMA-IMF fields. We chose the colours such that red traces dust thermal emission, blue traces free-free emission, and green likely traces a mixture of dust and free-free emissions. The maps display the spectral index value only where the 3 mm intensity exceeds 3 σ. The black ellipses display the sources extracted by getsf. We show the remaining fields in Fig. B.1.

4.2.2 Core temperatures

Following the pilot study by Motte et al. (2018), we used the Bayesian procedure PPMAP (Marsh et al. 2015; Dell’Ova et al. 2024) to build the temperature maps of the ALMA-IMF protoclusters. The procedure uses several continuum emission maps to compute a cube of column densities as a function of dust temperature. Along with the ALMA data in Band 6 decontaminated from free-free emission (Galván-Madrid et al. 2024), our PPMAP execution takes the following maps as input:

• SOFIA/HAWC+ data at 214 µm with an angular resolution of 19.0″ (only for G012.80, G351.77, W51-E, and W51-IRS2, Vaillancourt 2016; Pillai & Simplifi Team 2023).

• Herschel/PACS and Herschel/SPIRE data at 70, 160, 250, 350, and 500 µm with angular resolutions of 5.6, 10.7, 17.6, 23.9, and 35.2″, respectively (Molinari et al. 2010; Motte et al. 2010).

• APEX/SABOCA data at 350 µm with an angular resolution of 7.8″ (Lin et al. 2019).

• APEX/LABOCA data at 870 µm with an angular resolution of 19.2″ (Schuller et al. 2009).

We refer to Dell’Ova et al. (2024) for an in-depth description of the method. Finally, we obtained a unique temperature for each core at an angular resolution of 2.5″. The derived temperatures vary from 19.4 to 62.8 K (see tables published on CDS). The error associated with the PPMAP derivation is ≈5 K; to account for potential systematic contributions, we adopted a 25% error on the PPMAP temperature estimates⁶.

We adopted the PPMAP temperatures for all cores except for the hot-core candidates, for which we adopted the method proposed in Bonfand et al. (2024): We cross-correlated the position of continuum sources (see Sect. 3.3) with methyl formate (CH₃ OCHO) emission maps that were observed by ALMA-IMF. Methyl formate forms at the surface of dust grains in lukewarm (30–40 K) regions and is then released in the gas phase when the temperature reaches ≈100K (e.g. Garrod et al. 2009). Methyl formate can be used as a proxy to trace regions where heating is present and dust surface products have started to sublimate. Following the method discussed in Bonfand et al. (2024), we set 100±50 K to the 49 sources whose position corresponded to the peaks of extended methyl formate emission. These sources are classified as hot-core candidates.

Additionally, we set 300±100K to sources n°1, n°2, and n°1 in W51-IRS2, W51-E, and G327.29 respectively, following their detailed modelling by Ginsburg et al. (2017), Goddi et al. (2020), and the adopted temperature in Bonfand et al. (2024), respectively. These three sources are associated with strong and extended emission structures in the methyl formate emission maps⁷.

4.2.3 Core mass and boundedness

We converted the measured integrated flux at 1.3 mm, $S_{1.3\text{mm}}^{\text{int}}$, into a mass using the formula presented in Paper III, which includes a first-order correction for the optical thickness, $$M_{\text{core}}=-{\text{Ω}}_{\text{beam}}^{1.3\text{mm}}\frac{D^2}{{\kappa }_{1.3\text{mm}}}\frac{S_{1.3\text{mm}}^{\text{int}}}{S_{1.3\text{mm}}^{\text{peak}}}\ln \left(1-\frac{S_{1.3\text{mm}}^{\text{peak}}}{{\text{Ω}}_{\text{beam}}^{1.3\text{mm}}B(T,v)}\right)\text{,}$$(4)

where D is the distance to the target, ^K1.3_mm = 0.1 (v/1000GHz)^β cm²g⁻¹ is the dust opacity per unit mass (dust + gas), with an opacity index β = 1.5 typical of cold and dense environments (Ossenkopf & Henning 1994), and B(T, v) is the Planck function computed at the representative frequency of the observations in Band 6 with the dust temperature T (see Table 1 for central frequencies and distances, and the tables published at the CDS for the source temperature).

Errors on the masses mostly arise from the uncertainties on the opacity index β and on the dust temperature, T. The masses would be a factor of two lower for β = 2 instead of β = 1.5. As for the temperature, a 25% difference leads to a mass shift by ≲40%.

To address the boundedness of cores, we computed the ratio M_BE/M_core, where M_core is the core mass and M_BE is the mass of the critical Bonnor-Ebert sphere (Bonnor 1956), whose size matches that of the core: $M_{\text{BE}}=2.4\times R\times {\sigma }_{\text{th}}^2/G$. Here, R is the equivalent core radius, estimated as $R={\left(a\times b^2\right)}^{\frac{1}{3}}\times D$⁸, where a and b are the major and minor axis of the source ellipse, respectively. Moreover, σ_th is the thermal broadening of lines, as ${\sigma }_{\text{th}}=\sqrt{\gamma k_{\text{B}}T/\left(\mu \times m_{\text{p}}\right)}$, where γ = 1 is the adiabatic index (isothermal), k_B is the Boltzmann constant, µ = 2.4 is the mean molecular weight per free particle, and m_p is the mass of a proton. We considered the sources to be gravitationally bound when M_BE/M_source < 2 (see e.g. Louvet et al. 2021). As listed in Table 3, this filter excludes no sources. We note that this is a first-order check of the boundedness of cores. A more accurate determination would require computing the equilibrium of each source taking into account its turbulence, external pressure, and magnetic field support in addition to the thermal support. Unfortunately, we lack all these pieces of information at the moment.

In Fig. 2, we show the thermal-dust cores for each field. The tables published at the CDS list all the sources detected by the getsf extraction algorithm in each field. The first group of sources in each table corresponds to thermal dust cores that are gravitationally bound, and the second group corresponds to sources whose fluxes are arguably contaminated by free-free emission. For each source, the last column indicates whether it was also detected by GExt2D, which is true for 80% of the extracted sources.

5 Core mass function

5.1 Completeness tests

In order to draw a coherent sample gathering cores from the 15 ALMA-IMF protoclusters, we need to ascertain that the core extractions are complete at the selected lower-mass limit. To do this, we performed completeness tests in each protocluster: We injected synthetic cores on top of their background emission (i.e. the emission subtracted from each source). The synthetic cores were injected in the form of 2D circular Gaussian with an FWHM corresponding to 3840 au, which is the mean core size (see Sect. 3.3). We injected these synthetic sources randomly provided the centres of the cores were separated by more than 2.5^′′ . We forbid the injection of sources within 5^′′ of the border of the maps, where the noise increases due to the primary-beam correction (see Fig. 2). The flux of the synthetic sources ranges from equivalent masses of 0.5 to 5 M_⊙ in all fields for gas at 20 K, except in W51-E, where we adopted a mass range of 1 to 10 M_⊙. The fluxes were equally split into ten flux bins. From these fluxes, we computed the mass of the synthetic cores when adopting the mean temperature of the cores extracted in the corresponding region. We then extracted the synthetic cores with the getsf method (see Sect. 3.1). To obtain good statistics in each mass bin, we repeated this procedure four times per field. Cumulating the four draws, we obtained ≃3800 cores in total in each field (or ≃380 per bin), and we probed ≃80% of the background area. We plot the percentage of synthetic cores we extracted as a function of their mass in the supplementary plots published on Zenodo. Following Paper III, we considered that the completeness limit was reached when 90% of the synthetic cores were extracted. The mass completeness varies from 0.80 to 1.64 M_⊙ with a mean value of 1.34 M_⊙ and a standard deviation of 0.2 M_⊙. These numbers exclude W51-E, for which we find a completeness as high as 3.9 M_⊙. We excluded W51-E from the analysis of the global CMF (see Sect. 6.1), and we restricted the core samples to cores exceeding 1.64 M_⊙ to fit the high-mass tail of the CMF.

5.2 Fitting core mass functions

5.2.1 Method

The representation of a distribution in the form of histograms, as in Fig. 3, further fitted through linear regression, may give inaccurate results (Clauset et al. 2009). In addition, the representation in the form of a CMF in log-log scales prohibits us from estimating the fit uncertainty since the noise is not Gaussian. Moreover, the choice of the bin width induces a free parameter that thwarts any attempts to estimate the uncertainty. The representation in the form of a cumulative CMF, however, prevents us from estimating the uncertainty because the data are not independent. To circumvent these issues, we chose to fit the high-mass end of the CMF using maximum likelihood estimates (MLE). We followed the procedure presented by Clauset et al. (2009) and used its implementation in the Python package powerlaw presented by Alstott et al. (2014). When we assume that the CMF can be represented by a power law p(x) = Cx^−α, where C is a constant, there must be some lower-mass value, x_min , from which the power-law fit is accurate and that prevents the divergence of p when x → 0. After normalisation (imposing that ${{{\displaystyle \int }}^{\text{}}}_{x_{\min }^{}}^{+\infty }p(x)dx=1$), this reads $$p(x)=\frac{\alpha -1}{x_{\min }}{\left(\frac{x}{x_{\min }}\right)}^{-\alpha }.$$(5)

When x_min is known, the MLE gives an estimate of the exponent of the power law as $$\alpha =1+n{\left[{\displaystyle \sum }_{i=1}^n\ln \frac{x_i}{x_{\min }}\right]}^{-1},$$(6)

where x_i is the mass of the core i (provided x_i ≥ x_min), and the uncertainty on α reads $$\sigma =\frac{\alpha -1}{\sqrt{n}}+O\left(\frac{1}{n}\right),$$(7)

where O is the mathematical notation for “not negligible in front of”.

The originality of the method by Clauset et al. (2009) is that it proposes a method for determining the best parameters x_min and α using the Kolmogorov-Smirnov (KS) distance. First, the method computes α for x_min taking successively each mass value of the core sample. Then it computes the maximal distance (the KS distance, KS_D) between the observational CMF for all elements whose mass exceeds x_min and the synthetic CMF following a probability distribution function as defined by Eq. (5). The best value for the parameters x_min and α is the parameter that minimises the KS_D. Figure 5 shows the evolution of α and KS_D for all x_min values. In our case, we selected the highest value between the x_min that minimised the KS_D and our completeness limit (see Sect. 5.1). The first minimum in the KS_D was met for x_min =1.40 M_⊙, below our completeness level at 1.64 M_⊙. We therefore selected x_min = 1.64 M_⊙, which corresponds to a core sample of 330 elements (excluding W51-E). From this subsample, we obtained a power-law index fit α = 2.11 ± 0.06⁹. A sample of 330 cores above the completeness limit conveys a statistically conclusive sample. Clauset et al. (2009) reported stable index values when the sample exceeded ≈50 elements above the completeness limits. Consistently, Louvet et al. (2021) reported stable index values when the sample exceeded ≈40 elements above the completeness limits.

Fig. 5 Convergence towards the best fit of a power law through the MLE and KSD method for the global getsf bound core catalogue without the free-free sources. The left vertical dashed line highlights the first minimal KSD (in blue), occurring at xmin = 1.40 M⊙. The right vertical dashed line highlights the completeness limits at xmin = 1.64 M⊙. The orange curve shows the exponent of the power-law fit of the high-mass tail of the CMF as a function of xmin.

5.2.2 Uncertainty of the fit

As pointed out in Sect. 4.2.3, errors in the mass estimates arise primarily from the uncertainties on the opacity index and core temperature. To study how these uncertainties affect the slope of the CMF, we performed a Monte Carlo simulation on our 330 cores by allowing their opacity index, their temperature, and their source flux to vary simultaneously. For the source flux, we used the Gaussian error associated with the flux measurement. Koen & Kondlo (2009) showed that including these errors on the fluxes flattens the distribution. However, since our relative uncertainty diminishes for brighter sources, we expect this effect to be marginal at the high-mass end of our CMF. In our Monte Carlo simulation, the opacity index β can take random normal values with a mean value of 1.5 (Ossenkopf & Henning 1994) and a standard deviation of 0.2 for each core. In parallel, we allowed the core temperature to take random normal values with the mean and standard deviation as described in Sect. 4.2.2. We show the resulting cumulative CMFs generated from 10³ trials in Fig. 6. Each draw was fitted via MLE (see Sect. 5.2.1), fixing x_min to 1.64 M_⊙. We obtain a mean exponent value for the slope of the CMF of 1.97±0.06, where the uncertainty is the quadratic sum of the statistical uncertainty from Eq. (7) (σ ≃ 0.06) and from the uncertainties on the temperature, opacity index, and flux measurement uncertainties (σ ≃ 0.02). We show in the supplementary plots published on Zenodo the results of the 10³ trials for a power-law fit.

Fig. 6 The grey curves show the complementary cumulative CMFs of the 103 core samples obtained by including flux measurement uncertainties, varying core temperatures, and applying opacity index variations on a per-regions basis (see Sect. 5.2.2), overlaid with the average fit by a power law (in blue) and the best fit with a power-law index 2.35 (in red).

5.2.3 Comparison with the Salpeter slope

The ALMA-IMF large programme aims to test whether the high-mass slope of the CMF differs from the high-mass slope of the canonical IMF (α = 2.35 when fitted by a power law; Salpeter 1955). To determine whether our sample of 330 cores (with M > x_min = 1.64 M_⊙) permits such a claim, we selected 330 sources from a perfect mass distribution with an exponent of 2.35 and fitted the slope of its high-mass tail through the MLE. We repeated the operation 10⁵ times. Figure 7 shows the probability of retrieving an exponent α when the parental distribution has a slope α = 2.35. With 330 cores, the probability of retrieving a slope compatible with ALMA-IMF, α = 1.97 ± 0.06, is lower than 1% (σ ≈ 2.4). Therefore, we report that the CMF in the protoclusters studied by ALMA-IMF is flatter than and cannot reasonably be reconciled with the Salpeter-IMF slope to a 2.4σ level¹⁰.

Fig. 7 Power-law exponents α retrieved when sorting 330 cores from a synthetic core sample whose CMF has a power-law exponent of 1.97 (in blue) and 2.35 (in orange).

6 Discussion

6.1 Global top-heavy core mass function

In Sects. 5.2.1 and 5.2.2, we presented the best fit to the high-mass tail of the CMF by a power law. The best fit yields α = 1.97 ± 0.06. This result, based on 330 cores, confirms the many investigations conducted in recent years that reported top-heavy CMFs in high-mass protoclusters (e.g. Motte et al. 2018; Cheng et al. 2018; Sanhueza et al. 2019; Kong 2019; Moser et al. 2020). We show that this slope cannot be reconciled, at the 2.4σ level, with the canonical IMF slope (α = 2.35, see Sect. 5.2.3) even after the inclusion of all known uncertainties.

Observing top-heavy CMFs in high-mass protoclusters contrasts with nearby low-mass star-forming regions where a high-mass tail compatible with Salpeter is consistently found (see e.g. Table 6 of Fiorellino et al. 2021 for all CMF exponents from the Herschel Gould Belt Survey).

6.2 From the core mass function to the initial mass function

To link the CMF to the IMF, one must assume that cores constitute the mass reservoir for the accretion of protostars. One also assumes a fragmentation cascade and a mass transfer efficiency ϵ from the cores to the protostars. The common assumption for protostars is that one core will give birth to one star (and therefore assumes no sub-fragmentation) and that ϵ remains constant regardless of the core mass. If we adopt these assumptions, our massive protoclusters should engender a top-heavy IMF with a high-mass slope equal to that of the CMF: α ≃ 1.97. However, Bontemps et al. (2010), Louvet et al. (2014), and Csengeri et al. (2017) showed that ϵ depends almost linearly on the local density. The massive cores (>25 M_⊙) are spread equally over the observed size range and are therefore denser on average than low- and intermediate-mass cores (see the red crosses in Fig. 1). This suggests that massive cores will convert gas into stars more efficiently than low- and intermediate-mass cores, which causes a flatter IMF than the parental CMF. In the past ten years, some examples of top-heavy IMFs have been reported in young massive clusters in and beyond our Galaxy (see for example Lu et al. 2013; Schneider et al. 2018). As a consequence, we speculate that the ALMA-IMF protoclusters might be the precursors of star clusters whose IMF is top-heavy. From another perspective, we now more regularly observe star-forming regions with a top-heavy CMF than star clusters with a top-heavy IMF. To reconcile these two observables, we hypothesise that some of the star-forming regions hosting a top-heavy CMF will evolve to produce a Salpeter-like IMF.

6.3 Limitations and perspectives

We recognise two limitations in our study that we describe below.

• Free-free emission: In Sect. 4.2, we rejected all sources that are arguably contaminated by free-free emission from our core sample. We finally rejected 84 sources (see Table 3). Galván-Madrid et al. (2024) reported estimates of the free-free contribution based on the H41α recombination lines collected by ALMA-IMF. This characterisation will allow us to correct and re-introduce these sources in the core sample. These sources, with free-free emission, arise from massive young stellar objects and would steepen our CMF even further. Therefore, our conclusions should remain unchanged by the inclusion of the free-free contaminated cores.

• Sub-fragmentation: We probed the cores at a spatial resolution of 2700 au. Recent observations in the high-mass protostellar core G335 showed that a binary system, and perhaps even a triple system, could take place below 1000 au (Olguin et al. 2021, 2022). Similar results were obtained by Izquierdo et al. (2018) in W33A. However, other works reported no fragmentation in G336.01-0.82 (Olguin et al. 2023) or HH80-81 (Girart et al. 2018). The multiplicity in the ALMA-IMF fields is therefore currently uncertain. A straightforward path would consist of conducting observations at an even higher angular resolution to reach a spatial resolution of ≃200 au, the typical scale of solar-type proto-planetary disks (e.g. Louvet et al. 2018). The main difficulty consists of obtaining high sensitivities with an optically thin tracer at the scale of protoclusters.

7 Conclusions

We presented the core catalogues of the 15 high-mass protoclusters observed by the ALMA-IMF large programme. At a homogeneous sensitivity and spatial resolution (2700 au), we collected about 680 sources. Rejecting sources arguably contaminated by free-free emission, we analysed a core sample of ≃600 elements. We performed completeness tests and found the matched physical resolution core catalogues to be complete down to 1.64 M_⊙, with the exception of the W51-E protocluster, whose cores were discarded from the CMF analysis. In total, we analysed 330 cores with masses above the completeness limit of 1.64 M_⊙ and with little to no free-free contamination.

We fitted the high-mass slope of the core mass function (CMF) with the maximum likelihood estimate technique and found a best-fit power-law probability distribution function (PDF) $\frac{\text{d}N}{\text{d}M}\propto M^{-\alpha }$ with α = 1.97 ± 0.06. This exponent is flatter than and cannot be reconciled with the Salpeter-IMF slope of α ≃ 2.35. We confirm that the CMF in a representative sample of high-mass Galactic protoclusters is shallower than the Salpeter slope. We suggest that these massive protoclusters will give birth Ginsburg, A., Goddi, C., Kruijssen, J. M. D., et al. 2017, ApJ, 842, 92 to top-heavy stellar clusters. Alternatively, in order to reconcile our results with the Universal IMF, the CMF needs to evolve and become Salpeter-like in later stages of the cluster formation.

Together with this paper, we provide the core catalogues for the 15 protoclusters at both the native (1300 to 2700 au) and smoothed (2700 au) linear resolutions. These catalogues include the sources position, size (in arcsec and in au), peak, and integrated fluxes at both 1.3 and 3 mm, their temperature estimates from PPMAP, and their estimated mass.

Data availability

Tables D.1–D.24 are available at the CDS via anonymous ftp to cdsarc.cds.unistra.fr (130.79.128.5) or via https://cdsarc.cds.unistra.fr/viz-bin/cat/J/A+A/690/A33

Acknowledgements

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. This paper also uses the following ALMA data: ADS/JAO. ALMA#2013.1.01365.S and ADS/JAO.ALMA#2015.1.01273.S. FL acknowledges support by the Marie Curie Action of the European Union (project MagiKStar, Grant agreement number 841276). PS was supported by a Grant-in-Aid for Scientific Research (KAKENHI Number JP22H01271 and JP23H01221) of the Japan Society for the Promotion of Science (JSPS). P.S. and H.-L.L. gratefully acknowledge the support from the NAOJ Visiting Fellow Program to visit the National Astronomical Observatory of Japan in 2019, February. FM, YP, NC, and BT acknowledge the support of the European Research Council (ERC) via the ERC Synergy Grant ECOGAL (grant 855130) and from the French Agence Nationale de la Recherche (ANR) through the project COSMHIC (ANR-20-CE31-0009). AG acknowledges support from the NSF under grants AST 2008101 and CAREER 2142300. SB acknowledges support from the French Agence Nationale de la Recherche (ANR) through the project GENESIS (ANR-16-CE92-0035-01). T. Cs. has received financial support from the French State in the framework of the IdEx Université de Bordeaux Investments for the Future Program. RG-M and TN acknowledge support from UNAM-PAPIIT project IN108822 and from CONACyT Ciencia de Frontera project ID: 86372. AS gratefully acknowledges support by the Fondecyt Regular (project code 1220610), and ANID BASAL project FB210003. LB gratefully acknowledges support by the ANID BASAL projects ACE210002 and FB210003. GB acknowledges support from the PID2020-117710GB-I00 grant funded by MCIN/ AEI /10.13039/501100011033. GB and ALS acknowledge funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme, for the Project “The Dawn of Organic Chemistry” (DOC), grant agreement No 741002. TB acknowledges the support from S. N. Bose National Centre for Basic Sciences under the Department of Science and Technology, Govt. of India. MB 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). The project leading to this publication has received support from ORP, which is funded by the European Union’s Horizon 2020 research and innovation programme under grant agreement No 101004719.

Appendix A Continuum emission maps

Fig. A.1 Continuum emission maps (continued from Fig. 2).

Fig. A.1 Continuum emission maps (continued).

Fig. A.1 Continuum emission maps (continued).

Appendix B Spectral index maps

Fig. B.1 Spectral index maps (continued from Fig. 4). The ellipses in pink highlight sources automatically flagged due to their spectral index; those in yellow display the sources whose spectral index were further checked manually (see Eq. 2).

Fig. B.1 Spectral index maps (continued).

Fig. B.1 Spectral index maps (continued).

References

All Tables

All Figures

Fig. 1 Scatter plots of the source size (black dots) extracted with getsf before (left) and after (right) the maps were smoothed to a homogeneous spatial resolution of 2700 au. The blue triangles display the beam sizes. The orange diamonds display the mean sizes of sources. In the left panel, the y-axis is the linear scale (beam × distance) of the maps in Band 6, and the orange line shows the fit to the mean source size distribution for each field. In the right panel, the orange line shows the mean core size of all 15 protoclusters, and the red crosses illustrate the cores exceeding 25 M⊙.

In the text

Fig. 2 Fifteen ALMA-IMF protoclusters as traced by their continuum emission at 1.3 mm at a matched spatial resolution of 2700 au. The ellipses locate all the cores found by getsf. The name of the field is indicated in the top left corner of each panel, and the beam size is shown in the bottom left corner. For each protocluster, the centre coordinates are those specified in Table 1. We show the remaining fields in Fig. A.1.

In the text

Fig. 3 Source flux distributions for the getsf and GExt2D source extractions. The left panel displays the SFFs of the sources extracted by getsf in grey and by GExt2D in blue. The vertical blue and grey lines indicate the median values for each source catalogues, which differ by ≃5%. The right panel displays the complementary cumulative SFFs from the getsf source extraction (in grey) and from the GExt2D source extraction (in blue). The grey and blue lines display the best fits from linear regressions on the complementary cumulative SFFs from the median value at 12 mJy up to the flux at which the number of cores is smaller than 10. The corresponding ζ power-law indexes are indicated in the top right corner.

In the text

Fig. 4 Maps of the spectral index in the 15 ALMA-IMF fields. We chose the colours such that red traces dust thermal emission, blue traces free-free emission, and green likely traces a mixture of dust and free-free emissions. The maps display the spectral index value only where the 3 mm intensity exceeds 3 σ. The black ellipses display the sources extracted by getsf. We show the remaining fields in Fig. B.1.

In the text

Fig. 5 Convergence towards the best fit of a power law through the MLE and KSD method for the global getsf bound core catalogue without the free-free sources. The left vertical dashed line highlights the first minimal KSD (in blue), occurring at xmin = 1.40 M⊙. The right vertical dashed line highlights the completeness limits at xmin = 1.64 M⊙. The orange curve shows the exponent of the power-law fit of the high-mass tail of the CMF as a function of xmin.

In the text

	Fig. 6 The grey curves show the complementary cumulative CMFs of the 103 core samples obtained by including flux measurement uncertainties, varying core temperatures, and applying opacity index variations on a per-regions basis (see Sect. 5.2.2), overlaid with the average fit by a power law (in blue) and the best fit with a power-law index 2.35 (in red).
In the text

	Fig. 7 Power-law exponents α retrieved when sorting 330 cores from a synthetic core sample whose CMF has a power-law exponent of 1.97 (in blue) and 2.35 (in orange).
In the text

	Fig. A.1 Continuum emission maps (continued from Fig. 2).
In the text

	Fig. A.1 Continuum emission maps (continued).
In the text

	Fig. A.1 Continuum emission maps (continued).
In the text

	Fig. B.1 Spectral index maps (continued from Fig. 4). The ellipses in pink highlight sources automatically flagged due to their spectral index; those in yellow display the sources whose spectral index were further checked manually (see Eq. 2).
In the text

	Fig. B.1 Spectral index maps (continued).
In the text

	Fig. B.1 Spectral index maps (continued).
In the text