# Self-consistent Analysis of Stellar Clusters: An Application to HST Data of the Halo Globular Cluster NGC 6752

S. O. Souza, L. O. Kerber, B. Barbuy, A. Pérez-Villegas, R. A. P. Oliveira, D.

Nardiello

► To cite this version:

S. O. Souza, L. O. Kerber, B. Barbuy, A. Pérez-Villegas, R. A. P. Oliveira, et al. Self-consistent Analysis of Stellar Clusters: An Application to HST Data of the Halo Globular Cluster NGC 6752. The Astrophysical Journal, 2020, 890, 10.3847/1538-4357/ab6a0f. insu-03667135

HAL Id: insu-03667135

https://insu.hal.science/insu-03667135v1

Submitted on 7 Aug 2025

**HAL** is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published 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
<figure><img src="image_2.png" /></figure>

# Self-consistent Analysis of Stellar Clusters: An Application to HST Data of the Halo Globular Cluster NGC 6752

S. O. Souza<fnref n="1" />, L. O. Kerber<fnref n="1" />,<fnref n="2" />, B. Barbuy<fnref n="1" />, A. Pérez-Villegas<fnref n="1" />, R. A. P. Oliveira<fnref n="1" />, and D. Nardiello<fnref n="3" />,<fnref n="4" />,<fnref n="5" />

Received 2019 August 31; revised 2020 January 6; accepted 2020 January 8; published 2020 February 11

## Abstract

The Bayesian isochrone fitting using the Markov chain Monte Carlo algorithm is applied, to derive the probability distribution of the parameters age, metallicity, reddening, and absolute distance modulus. We introduce the SIRIUS code by means of simulated color–magnitude diagrams (CMDs), including the analysis of multiple stellar populations (MPs). The population tagging is applied from the red giant branch to the bottom of the main sequence. Through sanity checks using synthetic Hubble Space Telescope CMDs of globular clusters we verify the code reliability in the context of simple and MPs. In such tests, the formal uncertainties in age or age difference, metallicity, reddening, and absolute distance modulus can reach 400 Myr, 0.03 dex, 0.01 mag, and 0.03 mag, respectively. We apply the method to analyze NGC 6752, using Dartmouth stellar evolutionary models. Assuming a single stellar population, we derive an age of 13.7 ± 0.5 Gyr and a distance of $d⊙ = 4.11 ± 0.08$ kpc, with the latter in agreement within $3σ$ with the inverse Gaia parallax. In the analysis of the MPs, three populations are clearly identified. From the Chromosome Map and UV/Optical two-color diagrams inspection, we found a fraction of stars of 25 ± 5, 46 ± 7, and 29 ± 5 %, for the first, second, and third generations, respectively. These fractions are in good agreement with the literature. An age difference of 500 ± 410 Myr between the first and the third generation is found, with the uncertainty decreasing to 400 Myr when the helium enhancement is taken into account.

Unified Astronomy Thesaurus concepts: Bayesian statistics (1900); Globular star clusters (656); Open star clusters (1160); Hertzsprung Russell diagram (725)

## 1. Introduction

The study of stellar clusters has implications in a wide variety of astrophysical topics, which includes star formation, stellar evolution and nucleosynthesis, stellar dynamics, Galactic structure, and galaxy formation and evolution (e.g., VandenBerg et al. 2013; Barbuy et al. 2018).

With the advent of space-based telescopes, in particular the Hubble Space Telescope (HST) and more recently the Gaia Data Release 2 (DR2, Gaia Collaboration et al. 2018a), as well as multi-object and high-resolution spectrographs, a wealth of high-quality and spatially resolved data have been collected for Milky Way globular and open clusters (GCs and OCs), and for stellar clusters in neighboring galaxies. Combined with sophisticated analysis, these data have opened an unprecedented opportunity for very accurate physical parameter derivation.

Milky Way globular clusters (GCs) formed during the early stages of the Galaxy formation (e.g., VandenBerg et al. 2013; Barbuy et al. 2018) are studied in the present work.

The phenomenon of multiple stellar populations (MPs) was observed for the first time by Osborn (1971) from CN-band strengths, but at the time this was not identified as due to the presence of two stellar populations. Later, MPs were clearly revealed by e.g., Lee et al. (1999); Bedin et al. (2004); Piotto et al. (2005); Milone et al. (2017), and hints on self-enrichment to explain abundance variations within a GC were discussed by Gratton et al. (2004). Evidence of MPs from spectroscopic work was reviewed by Carretta (2019, and references therein). The photometric counterpart of the CN anomaly is detectable in

the ultraviolet (UV) filters (Piotto et al. 2015; Lee 2019). These filters are sensitive to C, N, and O abundances, allowing to disentangle the different stellar populations (Piotto et al. 2015).

With the purpose of correlating the cluster age with the presence of MPs, Martocchia et al. (2018, 2019) analyzed a sample of Magellanic Clouds (MCs) and Milky Way clusters. They estimated the N abundance spread in the color–magnitude diagram (CMD), which is an indicator of the presence of MPs, and found that clusters older than ~2 Gyr host MPs, while those younger than this age show no evidence of the spread in N abundance. On the other hand, it is known that the presence of MPs is related to the mass of the cluster (Milone et al. 2017). For this reason, age cannot be the only parameter to constrain the presence of MPs. This fact is evident for the case of Berkeley 39 (Martocchia et al. 2018) and Lindsay 38 (Martocchia et al. 2019), both having an age of ~6.5 Gyr, without showing N abundance spread. Another counterexample was given by Lagioia et al. (2019), having found that the GC Terzan 7 is consistent with a single stellar population (SSP), despite a relatively old age and high mass. Therefore, the study of MPs helps understanding the formation and evolution of stellar systems in general.

Isochrone fitting to CMDs has been extensively used to obtain the star cluster properties age, distance modulus, and reddening. Previously, a visual method known as “chi-by-eye” was usually employed to fit theoretical isochrones to CMDs. Later on, to benefit from improved data quality and to extract physical parameters with meaningful uncertainties, several statistical isochrone fitting techniques were developed, most of

<footnote n="1">¹Universidade de São Paulo, IAG, Rua do Matão 1226, Cidade Universitária, São Paulo 05508-900, Brazil; stefano.souza@usp.br</footnote>

<footnote n="2">²Universidade Estadual de Santa Cruz, Rodovia Jorge Amado km 16, Ilhéus 45662-000, Brazil</footnote>

<footnote n="3">³Dipartimento di Fisica e Astronomia “Galileo Galilei,” Università di Padova, Vicolo dell’Osservatorio 3, Padova, I-35122, Italy</footnote>

<footnote n="4">⁴Istituto Nazionale di Astrofisica—Osservatorio Astronomico di Padova, Vicolo dell'Osservatorio 5, Padova, I-35122, Italy</footnote>

<footnote n="5">⁵Aix Marseille Université, CNRS, CNES, LAM, Marseille, France</footnote>
<figure><img src="image_3.png" /><figcaption>**Figure 1.** SIRIUS flow-chart shows the steps to perform the isochrone fitting.</figcaption></figure>

them based on $\chi^2$, maximum likelihood statistics, or Bayesian approach (Kerber & Santiago 2005; Naylor & Jeffries 2006; von Hippel et al. 2006; Hernandez & Valls-Gabaud 2008; Monteiro et al. 2010). In almost all these developments, synthetic CMDs are employed for validation of the methods.

The Bayesian approach has the advantage of being able to get distributions and to explore the information a priori about the data or models. Recent examples of isochrone fitting codes using Bayesian inference are ASteCA (Perren et al. 2015) and BASE-9 (Stenning et al. 2016), where the latter allows analysis of MPs to derive their difference on the helium content (Y). Ramírez-Siordia et al. (2019) also applied Bayes’ theorem to a Monte Carlo method to get the posterior distributions of the same parameters as BASE-9, neglecting helium enhancements. They applied their software to the scarce stellar populations of ultra-faint dwarf galaxies and Large Magellanic Cloud star clusters.

In the present work, we carry out a detailed analysis of CMDs assuming both cases of clusters as SSPs and MPs. With this purpose, we developed the code named SIRIUS,<fnref n="6" /> standing for Statistical Inference of physical parameters of single and multiple populations in Stellar clusters, to extract information on a stellar cluster from its CMDs. The SIRIUS code was applied to analyze NGC 6752, with data from the HST UV-Legacy Survey of Galactic GCs (Piotto et al. 2015). Gratton et al. (2003) obtained for this halo GC an age of $13.4 \pm 1.1$ and Carretta et al. (2012) found three distinct stellar populations (Milone et al. 2013). Whereas the precision in parameter derivation from CMDs has been improving, it is also important to stress that a new era is now open: the age difference between stellar populations in a GC can give us a better understanding on its formation.

This work is organized as follows. In Section 2 the SIRIUS code is described in detail. Experiments to check the validity of the method and analysis of sources of uncertainties are presented in Section 3. An application to HST data of the halo GC NGC 6752 is presented in Section 4. Conclusions are drawn in Section 5.

## 2. The SIRIUS Code

This section gives a detailed description of the SIRIUS code, built to carry out isochrone fitting to CMDs, following the flow-chart presented in Figure 1.

### 2.1. CMD Data

SIRIUS was designed to analyze stellar clusters, applied here both to synthetic data and to observed data. SIRIUS has already been successfully applied to derive the parameters of two bulge GCs. For HP 1, a multi-band ($K_S$ and $J$ from Gemini-GSAOI+GeMS, and F606W from HST-Advanced Camera for Survey; ACS) isochrone fitting was applied (Kerber et al. 2019). For ESO 456-SC38, HST photometry in the filters F606W from ACS and F110W from WFC3, and FORS2@VLT photometry in V and I were used (Ortolani et al. 2019). These studies confirmed that HP 1 and ESO 456-SC38 are among the oldest GCs in the Milky Way, with an age of ~12.8 Gyr.

SIRIUS can create synthetic CMDs using the following method. The Monte Carlo algorithm is used to generate random data from a given probability distribution, and can be applied to describe many physical systems. In the case of CMDs of stellar clusters the main probability distribution of the system is the initial mass function (IMF), here adopted to be the Kroupa IMF (Kroupa 2001). The method to generate a sample of data similar to a stellar cluster is called as Synthetic CMD (Kerber et al. 2007). Points are randomly generated and interpolated in mass within theoretical points of isochrones. From an error function, these random points are dispersed by Gaussian distributions to simulate the spread seen in observed CMDs.

### 2.2. Stellar Evolutionary Models and Parameter Space

The library of isochrones adopted include two sets of stellar evolutionary models: Dartmouth Stellar Evolutionary Database (DSED; Dotter et al. 2008) and A Bag of Stellar Tracks and Isochrones (BaSTI; Pietrinferni et al. 2006). We perform linear

<footnote n="6">⁶ The code is available upon request to the authors.</footnote>
<figure><img src="image_4.png" /><figcaption>Figure 2. Graphical explanation of how the main five parameters change the morphology and position of the isochrone. The first panel shows the variation due to changes in $E(B-V)$, the second in $(m-M)_0$, the third in age, the fourth in [Fe/H], and the last one in Y.</figcaption></figure>

regressions to interpolate the isochrones in steps of 0.1 Gyr in age in the range of 10.0–15.0 Gyr, and 0.01 dex in [Fe/H] in the range of −2.00 < [Fe/H] < 0.00.<fnref n="7" /> It is relevant to mention that the range and step size of age we adopted here are consistent with the context of Galactic GCs. For the case of younger stellar clusters, e.g., MCs clusters, the age range should allow ages below 10 Gyr, and the step size should be narrower than the value used here.

The simple $\chi^2$ isochrone fitting procedures do not necessarily represent a physical interpretation of a GC CMD. Since the best fit is the isochrone that appears most similar to the CMD, many combinations of the parameters can be found as the best fit (minimum $\chi^2$; D’Antona et al. 2018).

The morphology of the isochrone depends on the age, reddening, absolute distance modulus, metallicity, and helium abundance. Figure 2 illustrates the effects on the shape of isochrones, due to the change in each of these parameters. The reddening $E(B-V)$ changes the location of the isochrone in the diagonal direction because it contributes to the apparent distance modulus $(m-M)_\lambda$ and reddening $E(\lambda_1 - \lambda_2)$, without varying the morphology of the isochrone (first panel). For high values of reddening, a second-order correction, from the effective temperature (e.g., Ortolani et al. 2017; Kerber et al. 2019), has to be taken into account in the isochrone fitting. A vertical displacement is the result of a change in distance modulus $(m-M)_0$ (second panel). Age $\tau$ affects essentially the position of the turn-off point (TO; third panel). The metallicity [Fe/H] has a complex effect on the isochrone, but more strikingly by changing the slope of the red giant branch (RGB), with a sub-giant branch (SGB) and RGB steeper toward lower metallicities (fourth panel of Figure 2). A variation in Y changes the slope of the SGB and the location of the TO, shifting the isochrone to the bluer region of the CMD (last panel). A review on the interpretation of CMDs in terms of stellar evolution models can be found in Gallart et al. (2005).

## 2.3. Bayesian Statistics: Isochrone Fitting

The Bayesian statistics is based on Bayes’ Theorem. The probability that two events (M and D) are true, at the same time, according to a null hypothesis H is given by the product probability law:

$$P(M, D|H) = P(M|D, H) \times P(D|H),$$

where $P(M|D, H)$ represents the probability of M to be true if D is true as well according to H, and $P(D|H)$ is the probability of D following H. The opposite is also valid:

$$P(D, M|H) = P(D|M, H) \times P(M|H).$$

From the hypothesis of the conditional probability of M and D to be the same as D and M, results in Bayes’ theorem:

$$P(M|D) = \frac{P(D|M) \times P(M)}{P(D)},$$

where, in our case, the evolutionary model is represented by M and the data by D.

The posterior distributions $P(M|D)$ are the distributions a posteriori of the model (M) and will give the distributions for each parameter. On the right-hand $P(M)$ are the prior distributions that give the information a priori about the model. The priors are distributions that constrain the parameters with the physical information.

Assuming that stars are distributed in color and magnitude following a Gaussian distribution and disconsidering the dependence of color with magnitude, the likelihood is given by:

$$P(D|M) = \prod_{i}^{N} \prod_{j}^{M} e^{-\varphi_{\text{color}}^{2}} \cdot e^{-\varphi_{\text{Mag}}^{2}},$$

where N is the total number of the analyzed stars and M is the number of points in the isochrone. The $\varphi^2$ is defined as, for

<footnote n="7">⁷ The usual notation [Fe/H] = log(Fe/H)star−log(Fe/H)⊙ is adopted.</footnote>
example:

$$ \varphi_{\text{color}_{i,j}}^2 = \frac{1}{2} \left( \frac{\text{color}_i^{\text{obs}} - \text{color}_j^{\text{iso}}}{S_i + \sigma_i^{\text{Cor}}} \right)^2, $$

where S represents the entropy term of likelihood. This term is responsible for smoothing the region of highest spread and number of stars. The $S_i$, $|color_i^{\text{obs}} - \xi_f|$, is calculated for each star by comparison with the fiducial color $\xi_f$, which is defined as the median color for a bin of magnitude centered on the magnitude of the i-th star.

The maximum likelihood $\mathcal{L}$ corresponds to a maximization of the likelihood function in the parameter space. It is given by (in logarithm form):

$$ \mathcal{L} = \max \left\{ - \sum_{i=1}^{N} \sum_{j=1}^{M} \left[ \varphi_{\text{color}_{i,j}}^2 + \varphi_{\text{Mag}_{i,j}}^2 \right] \right\}. $$

Since the exponential function can reach high values quickly, it is convenient to work with Bayes' theorem in the logarithmic form:

$$ \ln P(M|D) = \ln P(M) + \mathcal{L}. $$

**Priors**—The prior distributions ($P(M)$) are the main difference between the Bayesian and the frequentist statistics. These distributions impose constraints on the free parameters, restricting the set of parameters to be explored. In an isochrone fitting, these priors reflect the physical constraints, such as: (a) the upper age limit as the age of the universe (Planck Collaboration et al. 2016); (b) the metallicity values taken from high-resolution spectroscopy; (c) distances constrained and primordial He content from RR Lyrae mean magnitudes, for example; and (d) non-negative reddening values.

**Marginalization**—In order to explore the parameter space as a whole and to get the posterior distributions of each parameter, we applied Bayes' theorem with the Metropolis-Hastings (MH) algorithm (Metropolis et al. 1953; Hastings 1970). The method is basically an exclusion iterative algorithm, built first to solve problems of statistical physics. The MH method compares the random probabilities trying to reach the minimum energy state, which justifies that we can neglect the normalization term of Bayes' law. The final result of MH is a chain with n energies for m states that is known as the Markov chain. For the applications with random distributions, which means Monte Carlo methods, the result from the MH algorithm is called the Markov chain Monte Carlo (MCMC; Hogg & Foreman-Mackey 2018). To get the probability distributions of the parameters, the marginalization is executed by the integral:

$$ \mathcal{P}(\phi) = \int \mathcal{L}(\phi) \times p(\phi) d\phi, $$

where ($\phi$) represents the parameter space. To perform the marginalization from MH algorithm and MCMC method, we employed the Python library `emcee` (Foreman-Mackey et al. 2013).

## 2.4. Multiple Stellar Populations in GCs

Before carrying on the analysis of MPs, in this section we describe the separation of stellar populations in the CMDs. The stellar population tagging allows us to distinguish the first (1G) and second (2G) generation stars (and subsequent ones) from a

given CMD. Figure 3 shows the procedure we follow to separate the stellar populations in each region of the created synthetic CMD with $\Delta\tau_{1G,2G} = 0.50$ Gyr. We adopted a Dartmouth (DSED) isochrone with $[Fe/H] = -1.26$, $E(B-V) = 0.18$, $(m-M)_0 = 14.38$, and $\tau = 13.0$ Gyr.

In Milone et al. (2013) the pseudo-color C was defined, with the purpose to maximize the separation among MPs on the CMD. Piotto et al. (2005) have shown the power of HST UV filters F275W, F336W, and F438W to separate the MPs. F275W is sensitive to OH and F438W to CN and CH. For these filters, the 1G stars are fainter than the 2G because the latter are oxygen- and carbon-poorer than the 2G ones. For the filter F336W, which is sensitive to NH, the 1G stars are brighter than the 2G stars, given the fact that the 2G stars are nitrogen-richer. Note that stronger lines lead to larger opacity, and lower brightness. For these reasons, the color (F275W-F438W) inverts the stellar populations on the CMD with respect to the color (F336W-F438W). In that color, the 2G stars seem to be redder than the 1G stars (Piotto et al. 2005, their Figure 2).

**Chromosome maps (RGB and MS)**—Milone et al. (2017) describe the method of MP separation using chromosome maps based on combinations of UV HST filters. Lee (2019) used UVB data to distinguish MPs, and reviewed methods discussed earlier. To construct the chromosome map diagrams, we adopt the method presented in Milone et al. (2017) that is briefly described below. For the CMDs $m_{F814W}$ versus $C_{F275W,F336W,F438W}$ and $m_{F814W}$ versus $(m_{F275W} - m_{F814W})$, the red and blue fiducial lines are defined by 96th and 4th percentiles, respectively. The top- and bottom-middle panels of Figure 3 show the red and blue fiducial lines enclosing the RGB and MS stars, respectively. The axis of chromosome map are the relative distance between each stars and the fiducial lines, defined by:

$$ \Delta_C F275W,F336W,F438W = \frac{C_r - C}{C_r - C_b}, $$

$$ \Delta_{F275W,F814W} = \frac{G - G_r}{G_r - G_b}, $$

where the indices r and b refer to the red and blue fiducial lines, respectively. The color G represents $m_{F275W} - m_{F814W}$.

The diagram $\Delta_C F275W,F336W,F438W$ versus $\Delta_{F275W,F814W}$ quantifies the color distance of each star to the blue and red envelopes, so that the $\Delta$-value is closer to zero as the star is closer to the red envelope. The right panels of Figure 3 show the final chromosome maps for the RGB (top) and MS (bottom), respectively, for the synthetic CMD.

Some modifications on the identification of the MPs were implemented in the original method from Milone et al. (2017), in order to preserve a uniformity in the MPs separation for the three evolutionary stages (MS, SGB, RGB). The identification of the MPs is done using the Gaussian Mixture Models (GMM), that is a non-supervised machine learning algorithm, which searches to fit K Gaussian distributions to a sample of N data. The fit comes from the basic equation of Bayes' theorem:

$$ G(x) = \sum_{i=1}^{K} \phi_i \times \mathcal{N}(x | \mu_i, \sigma_i), $$

where $\mathcal{N}(x | \mu_i, \sigma_i)$ represents the i-th Gaussian distribution with mean of $\mu_i$ and standard deviation of $\sigma_i$. This algorithm was adopted from the python library `Scikit-learn` (Pedregosa et al. 2011).
<figure><img src="image_5.png" /><figcaption>**Figure 3.** MP separation and population tagging applied to synthetic data with $\Delta\tau = 0.5$ Gyr. Left panel shows the pseudo-color C, which gives a pronounced MP separation. Middle panels show the procedure we apply to separate the stellar populations, from top to bottom are the RGB, SGB, and MS stars, respectively. Right panels show the stars identified to belong to the 1G and 2G.</figcaption></figure>

We here assume two subclasses for GMM in a two-dimensional plane. Then, each star is classified as 1G or 2G according to the strength of the two Gaussian distributions on that point of the chromosome map. The separation between the two populations includes clear members of both, but as well stars in the limiting intersection, that can contaminate each other samples. This analysis can be improved by increasing the number of subdivisions in GMM to select the bona-fide stars of each stellar populations, as in Milone et al. (2018).

*Two-color diagrams (SGB)*—Since the SGB sequence, depending on the adopted filter and the metallicity of the cluster, could be nearly horizontal and their MPs could appear mixed, the chromosome maps are not effective with these stars. Therefore, we applied a conventional two-color diagram $m_{F336W} - m_{F438W}$ versus $m_{F275W} - m_{F336W}$, as described in Nardiello et al. (2015b). In order to apply the GMM procedure (same as described in the previous section), $\Delta_1$ and $\Delta_2$ are the axes that were normalized and then rotated counterclockwise by an angle of 45°. The method is graphically represented in Figure 3 (middle panels).

## 2.5. Age Difference $\Delta\tau$

The origin of the 2G (and subsequent populations) stars is a major challenge in the MP analyses. Most scenarios trying to

explain MP formation predict an age difference ($\Delta\tau$) between the first and the later populations (Bastian & Lardo 2018). For example, the scenario of asymptotic giant branch (AGB) stars polluting the second and subsequent populations, predicts a difference around 100 Myr (D'Antona et al. 2016), up to 200–700 Myr from the delay of X-ray binaries (Renzini 2013; Renzini et al. 2015). Another scenario is that of the supermassive stars (SMS). MPs can be formed from multiple bursts of SMSs with intervals of a few Myr (Gieles et al. 2018). Another possibility are the fast rotating massive stars (FRMSs) that would enrich the interstellar medium in about 40 Myr (Decressin et al. 2007; Krause et al. 2013). Therefore, the age difference between the first and next populations is an important parameter to give hints to their plausible origin.

From our population tagging method, we can analyze separately each stellar population from their CMDs. To perform the isochrone fitting in the context of MPs we developed a hierarchical algorithm to estimate the $\Delta\tau$ between the first and subsequent populations. The hierarchical algorithm considers the stars as a SSP first, and subsequently each stellar population. For a SSP we leave all parameters free. In the context of MPs, it is expected that the age of a SSP is a weighted average age of each stellar population. Consequently, for the example of two stellar populations, the ages could be
<table><caption>**Table 1**
Input Parameters for the Construction of the Synthetic Catalogs</caption><thead><tr><th>Parameter</th><th>No-spread</th><th>Spread</th></tr></thead><tbody><tr><td>Evolutionary Model</td><td>DSED</td><td>DSED</td></tr><tr><td>N<sub>total</sub></td><td>260</td><td>10,000</td></tr><tr><td>τ<sub>SSP</sub> (Gyr)</td><td>13.0</td><td>13.0</td></tr><tr><td>Δτ (Gyr)</td><td>...</td><td>0.1, 0.5, 1.5</td></tr><tr><td>[Fe/H] (dex)</td><td>-1.26</td><td>-1.26</td></tr><tr><td>E(B - V)</td><td>0.18</td><td>0.18</td></tr><tr><td>(m-M)<sub>0</sub></td><td>14.38</td><td>14.38</td></tr><tr><td>f<sub>bin</sub></td><td>...</td><td>0.30</td></tr><tr><td>q<sub>min</sub></td><td>...</td><td>0.60</td></tr><tr><td>N<sub>1G</sub>/N<sub>total</sub></td><td>1.000</td><td>0.360</td></tr></tbody></table>

derived from:

$$ \tau_{1G} = \tau_{SSP} + \Delta\tau \times (N_{1G}/N_{\text{total}}), \\ \tau_{2G} = \tau_{SSP} - \Delta\tau \times (1 - N_{1G}/N_{\text{total}}). $$

The hierarchical method fits the first population and applies the constraints of distance, reddening, and metallicity to the second (or subsequent) one(s). Hence, the procedure to compute the $\Delta\tau$ turns out simply to be $\Delta\tau = \tau_{1G} - \tau_{2G}$. This procedure considers that 1G stars were formed earlier than others, which is logical when our objective is to estimate a $\Delta\tau$. The likelihood of hierarchical procedure $\ln P(M|D)$ takes into account the constraints of a stellar cluster as a whole. For example, all stars must have the same values of distance and must be influenced by interstellar dust in the same way. Therefore, the likelihood of 1G ($\mathcal{L}(1G)$) and NG ($\mathcal{L}(NG)$) are dependent on the likelihood of SSP ($\mathcal{L}_{SSP}$). The total likelihood $\ln P(M|D)$ is a linear combination of the priors and the likelihood of each stellar population with influence of SSP parameters:

$$ \ln P(M|D) = \ln P(M) + \sum_{i=1}^{N} [\mathcal{L}([i]G)_{SSP} + \ln(f_{[i]G})], $$

where $f_{[i]G}$ represents the fraction of stars that belong to the i-th population. A similar likelihood based on MPs and weighted by the fraction of stars is applied in Ramírez-Siordia et al. (2019).

Here, we are adopting that the 1G stars have primordial helium content ($Y$), which is consistent with the literature (Bastian & Lardo 2018). Wagner-Kaiser et al. (2016) performed a bayesian isochrone fitting, in the context of MPs, for a sample of 30 GCs. Differently from the present work, they fitted the value of $Y$ for the 1G stars, resulting in some cases in a high content of $Y_{1G} \sim 0.30$. They also assumed the same age for both analyzed stellar populations. In contrast, we are interested in finding if there is an age difference between the stellar populations. Even though our approach is similar to the one applied in Wagner-Kaiser et al. (2016), the methods are based on different assumptions.

### 3. Controlled Experiment

In this section, we test the reliability of our analysis by using synthetic CMDs. First, we constructed a synthetic CMD using an error function obtained from the atlas extracted by Nardiello et al. (2018) from the data of the HST UV-Legacy Survey of Galactic GCs (Piotto et al. 2015), allowing us to simulate MPs with the synthetic data. The stellar evolutionary model adopted was the DSED isochrone with $Z \sim 0.002$ with $[\alpha/Fe] = +0.4$, and age of 13.0 Gyr, as reported in Table 1, corresponding to

typical values of moderately metal-poor bulge GCs (e.g., Kerber et al. 2018, 2019). We simulated the CMD of a cluster with a total number of 10,000 stars ($N_{\text{total}}$) that host 36% of 1G stars with an age of 13.0 Gyr and 64% of 2G stars 0.5 Gyr younger than 1G stars. We considered a fraction of binaries ($f_{\text{bin}}$) of 30% and a minimum mass ratio ($q_{\min}$) of 0.60. Resulting CMDs combining the different available filters are shown in Figure 4.

#### 3.1. Sources of Uncertainty

In our method, during the isochrone fitting, we compute the likelihood star-by-star. To keep the high performance of MCMC, we imposed a range in magnitudes based on stellar evolutionary models. The third panel of Figure 2 shows that there is no significant difference regarding the age for the $\sim 3$ magnitudes brighter than the TO. For this reason, we do not take into account stars above this limit in the likelihood calculation.

The faintest stars are limited to the completeness limit, meaning that the number of faint stars depends on the photometric depth. There are no differences between the isochrones in the databases employed in SIRIUS for the faintest stars ($\sim 2$ magnitudes below the TO), therefore the fit does not depend on the faintest stars. Ramírez-Siordia et al. (2019) presented an analysis considering the faintest stars. They concluded that the effect of faintest stars only increases the uncertainties without changing the mode of distribution, since the isochrones do not seem to be different for the faintest stars, as shown in Figure 2 (third panel).

As regards to binary stars, their magnitudes represent the combination of fluxes from the two companion stars. Since the magnitude is the logarithm of the stellar flux, for a binary system with two stars of the same mass, the magnitude of this system corresponds to the magnitude of one star subtracted by $2.5 \times \log(2) \sim 0.75$ (Kerber et al. 2002, 2007). The decrement in magnitude tends to have the binary stars to be brighter and redder on the CMD. To reduce the effect of binary systems during the isochrone fitting, SIRIUS takes into account only the stars within $3\sigma$ from the fiducial line of the CMD.

The standard BaSTI isochrones overestimate ages by $\sim 0.80$ Gyr, with respect to DSED isochrones. The main reason for this discrepancy is that BaSTI isochrones do not include atomic diffusion in the calculations, among other differences in basic physics. Whereas the solar alpha-to-iron more complete models, including atomic diffusion, are already available in Hidalgo et al. (2018), the available alpha-enhanced models taking this effect into account are not yet available.

#### 3.2. Sanity Check

In the optical wavelengths some filters are more sensitive to some properties than others. For the NIR filters the effect of interstellar medium extinction is considerably lower than for the UV filters. Also, a color combining filters with a small bandwidth is more suitable to observe the structures on the CMD. Therefore, the combination of magnitudes and colors on the CMD is very important regarding the information that is expected to be obtained from isochrone fitting. In order to estimate the effect of the choice of color we performed the isochrone fitting using ten different colors, without spreading the stars, combining the five HST filters available in the UV-Legacy survey of GCs (Piotto et al. 2015).
<figure><img src="image_6.png" /><figcaption>**Figure 4.** CMDs for the synthetic data using a DSED isochrone with age =13.0 Gyr, [Fe/H] = -1.26, E(B-V) = 0.18, (m-M)$_0$ = 14.38, Δτ = 0.50 Gyr and fraction of 1G stars ($N_{1G}/N_{total}$) = 0.360, generated from HST filters. All available combinations of filters are shown.</figcaption></figure>

First, we perform the fit considering the SSP without taking into account the photometric spread of stars. The DSED isochrones are here fitted to the synthetic no-spread catalog data (Table 1) with the purpose of checking if the input parameters of the synthetic CMD are recovered. For this test, we adopted uniform distribution priors for all parameters. The range of values we used are: for age, between 10 and 15 Gyr; for the metallicity, between 0.00 and -2.00 dex; for reddening, between 0.0 and 1.0 mag; and for the distance modulus, between 12.0 and 16.0 mag. Figure 5 shows the behavior of the parameter space as a function of color. It can be observed that the age is the most sensitive parameter to the filters, whereas the other parameters vary only slightly with the choice of filters. For color 8 (third lower panel in Figure 4), which is equivalent to B-V, there is a strong effect on the age, whereas for color 6 (first lower panel in Figure 4) the parameters are closer to the original ones. Color 10 ($m_{F606W}$ - $m_{F814W}$, last lower panel in Figure 4), is also close to the input values and has small uncertainties due to its lowest reddening-dependency. Therefore, for our analysis, we chose color 10.

Second, to verify the sensitivity of the method, we simulate real data through synthetic CMDs to perform the isochrone

fitting, taking into account a spread of stars, and assuming Gaussian priors centered on the parameters given in Table 1 (spread). In Figure 6, we show the isochrone fitting for the synthetic CMD with Δτ = 0.50 Gyr, assuming that it is SSP (left panel) and MPs (right panel). We employ the corner plots to present the posterior distributions (Figure 7). They show the N parameter space in a 2D representation, where it is possible to see the correlations between the parameters. As the best value for each parameter we adopted the mode of the distributions. For the confidence interval, we selected the 16th and 84th percentile of the distributions that give us the values inside 1σ from the mode. The top-left panel in Figure 7 shows the corner-plot for the DSED SSP isochrone fitting. Figure 7, in the top-right, bottom-left, and bottom-right panels show the results for the age derivation in the context of MPs using DSED.

Even though the spread of stars changes the visual aspect of the CMD, the parameters obtained from the isochrone fitting given in Table 2 for SSP and MPs, are both in good agreement with the input values from Table 1. In conclusion, in this section we were able to describe the approach and check the validity of SIRIUS in the context of MPs.
<figure><img src="image_7.png" /><figcaption>**Figure 5.** Sanity check with no-spread data, the parameter space as a function of color. The posterior distributions of each parameter for the ten combinations of *HST* filters of the UV-Legacy survey of GCs (Piotto et al. 2015). DSED isochrones are adopted. The numbers represent each color.</figcaption></figure>

<table><thead><tr><th></th><th>COLOURS</th></tr></thead><tbody><tr><td>1)</td><td>m<sub>F275W</sub> - m<sub>F336W</sub></td></tr><tr><td>2)</td><td>m<sub>F275W</sub> - m<sub>F438W</sub></td></tr><tr><td>3)</td><td>m<sub>F275W</sub> - m<sub>F606W</sub></td></tr><tr><td>4)</td><td>m<sub>F275W</sub> - m<sub>F814W</sub></td></tr><tr><td>5)</td><td>m<sub>F336W</sub> - m<sub>F438W</sub></td></tr><tr><td>6)</td><td>m<sub>F336W</sub> - m<sub>F606W</sub></td></tr><tr><td>7)</td><td>m<sub>F336W</sub> - m<sub>F814W</sub></td></tr><tr><td>8)</td><td>m<sub>F438W</sub> - m<sub>F606W</sub></td></tr><tr><td>9)</td><td>m<sub>F438W</sub> - m<sub>F814W</sub></td></tr><tr><td>10)</td><td>m<sub>F606W</sub> - m<sub>F814W</sub></td></tr></tbody></table>

<figure><img src="image_8.png" /><figcaption>**Figure 6.** Sanity check with spread data, isochrone fitting for the synthetic CMD considering SSP (left) and MPs (right) for DSED isochrones. The gray dots are discarded for the fit.</figcaption></figure>

## 4. Application to the Halo Globular Cluster NGC 6752

*HST* photometric data for NGC 6752 in the UV filters within the UV-Legacy Survey GO-13297 (PI. G. Piotto), and in the optical within GO-10775 (PI. A. Sarajedini) are used. These programs made available data in the UV filters F275W, F336W, and F438W from the Wide Field Camera 3 (WFC3), and the optical filters F606W and F814W from the WFC of the ACS. The newly reduced catalogs presented in Nardiello et al. (2018) are used.

NGC 6752 is a halo cluster, located at $l = 336.49^\circ$, $b = -25.63^\circ$, with a distance from the Sun $d_\odot = 4.0$ kpc (Harris 1996, edition 2010).<fnref n="8" /> A metallicity of [Fe/H] = −1.48 ± 0.07 dex was derived by Gratton et al. (2005) from high-resolution spectroscopy (R = 40,000) of seven stars near the red giant branch bump. Gratton et al. (2003) and VandenBerg et al. (2013) obtained an age of 12.50 ± 0.25 Gyr and 13.4 ± 1.1 Gyr, respectively. Carretta et al. (2012) identified three stellar populations based on three values of abundances of O, Na, Mg, Al, and Si elements that are sensitive to stellar populations in GCs, denominated as first (P), intermediate (I), and extreme (E) populations. Milone et al. (2013) gave the first photometric evidence of three stellar populations by using *HST* data. Nardiello et al. (2015a), using FORS2/VLT data, have observed the split of the MS of NGC 6752 using UBI filters, and calculated the radial distribution of the populations and the difference in helium between the 1G and 2G stars. Milone et al. (2019) confirmed the existence of three stellar populations from NIR photometric data on MS stars. Cordoni et al. (2020) analyzed the kinematics of the P and E populations of NGC 6752, and they found that there is no difference in rotation between the two stellar populations.

In order to separate the populations P, I, and E (hereafter 1G, 2G, and 3G), the number of components on GMM were increased to three for the RGB and SGB, and to four for the MS. The classification of 1G, 2G, and 3G stars is in agreement with Milone et al. (2013), since a clear distinction of three stellar populations can be verified in Figure 8. Milone et al. (2013) derived the mass fraction of each population to be of ~25, ~45, and ~30%, respectively. We found a fraction of stars of 25 ± 5%, 46 ± 7%, and 29 ± 5% for the 1G, 2G, and 3G, respectively, in excellent agreement with Milone et al. (2013).

<footnote n="8">⁸ www.physics.mcmaster.ca/~harris/mwgc.dat</footnote>
<figure><img src="image_9.png" /><figcaption>**Figure 7.** Sanity check two, corner plots using DSED isochrones, relating physical parameters. Top-left panel: results of the sanity check applied to a synthetic SSP CMD where Monte Carlo spread of data is implemented, with a $\Delta\tau = 0.50$ Gyr. Other panels: 1G and 2G combined for $\Delta\tau = 0.10$ Gyr (top right), $\Delta\tau = 0.50$ Gyr (bottom left), and $\Delta\tau = 1.50$ Gyr (bottom right).</figcaption></figure>

In the following the analysis of NGC 6752 is restricted to DSED isochrones. The procedure starts with the isochrone fitting assuming the CMD to consist of a SSP, and the method is subsequently applied to the MPs. In order to carry out the isochrone fitting, we employed the same CMD $m_{F606W}$ versus ($m_{F606W} - m_{F814W}$) used for the synthetic data. In the left panel of Figure 8 is shown the CMD of NGC 6752 including all stars as a SSP. The value of $[Fe/H] = -1.48$ dex was used as prior through Gaussian distribution with standard deviation of 0.07. A prior in distance was applied with the value of apparent distance modulus ($m - M)_V = 13.26 \pm 0.08$ taken from Gratton et al. (2003). The results of SSP isochrone fitting are shown in Table 3

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

and Figures 9 and 10. The SSP age derivation of $13.7 \pm 0.5$ Gyr is in good agreement with Gratton et al. (2003), that obtained $13.4 \pm 1.1$ Gyr, and with the Bayesian technique from Wagner-Kaiser et al. (2017) that resulted in an age of $13.202^{+0.174}_{-0.152}$ Gyr. The parallax from Gaia DR2 (Gaia Collaboration et al. 2018b) for the NGC 6752, $\bar{\omega} = 0.2610 \pm 0.0011$ mas, corrected by the zero-point of $-0.03$ mas given by Lindegren et al. (2018), gives a heliocentric distance of $3.85 \pm 0.02$ kpc. Considering NGC 6752 as a SSP, the derived distance is $4.11 \pm 0.08$ kpc, in agreement within $3\sigma$ with Gaia DR2.

The metallicity estimated from SSP isochrone fitting, $[Fe/H] = -1.49^{+0.05}_{-0.05}$, was fixed for the MPs approach. The
<figure><img src="image_11.png" /><figcaption>**Figure 8.** MPs in NGC 6752. Left panel: SSP; middle panel: same as left panel, but color-identified stars; right panel: pseudo-color showing the clear separation of three stellar populations.</figcaption></figure>

<table><caption>**Table 2**

Sanity Check with Spread Data, Results Summarized for Synthetic Data in SSP Context and MPs</caption><thead><tr><th>Sanity Check</th><th>N<sub>1G</sub>/N<sub>Tot</sub></th><th>Model</th><th>τ<sub>SSP</sub><br>(Gyr)</th><th>Δτ<sub>I<sub>G,2G</sub></sub><br>(Gyr)</th><th>[Fe/H]<br>(dex)</th><th>E(B - V)<br>(mag)</th><th>(m - M)<sub>0</sub><br>(mag)</th></tr></thead><tbody><tr><td rowspan="2">SSP</td><td rowspan="2">...</td><td>DSED</td><td>12.70<sup>+0.36</sup><sub>-0.37</sub></td><td>...</td><td>-1.26<sup>+0.03</sup><sub>-0.03</sub></td><td>0.18<sup>+0.01</sup><sub>-0.01</sub></td><td>14.35<sup>+0.03</sup><sub>-0.03</sub></td></tr><tr><td>BaSTI</td><td>13.80<sup>+0.61</sup><sub>-0.61</sub></td><td>...</td><td>-1.26<sup>+0.03</sup><sub>-0.03</sub></td><td>0.18<sup>+0.01</sup><sub>-0.01</sub></td><td>14.30<sup>+0.04</sup><sub>-0.03</sub></td></tr><tr><td rowspan="2">MPs Δτ = 0.10 Gyr</td><td rowspan="2">0.377 ± 0.011</td><td>DSED</td><td>...</td><td>0.11<sup>+0.36</sup><sub>-0.38</sub></td><td>-1.26<sup>+0.02</sup><sub>-0.03</sub></td><td>0.18<sup>+0.01</sup><sub>-0.01</sub></td><td>14.38<sup>+0.03</sup><sub>-0.03</sub></td></tr><tr><td>BaSTI</td><td>...</td><td>0.19<sup>+0.49</sup><sub>-0.49</sub></td><td>-1.26<sup>+0.03</sup><sub>-0.03</sub></td><td>0.18<sup>+0.01</sup><sub>-0.01</sub></td><td>14.33<sup>+0.03</sup><sub>-0.03</sub></td></tr><tr><td rowspan="2">MPs Δτ = 0.50 Gyr</td><td rowspan="2">0.370 ± 0.012</td><td>DSED</td><td>...</td><td>0.41<sup>+0.43</sup><sub>-0.37</sub></td><td>-1.26<sup>+0.03</sup><sub>-0.02</sub></td><td>0.18<sup>+0.01</sup><sub>-0.01</sub></td><td>14.38<sup>+0.03</sup><sub>-0.03</sub></td></tr><tr><td>BaSTI</td><td>...</td><td>0.51<sup>+0.54</sup><sub>-0.54</sub></td><td>-1.26<sup>+0.02</sup><sub>-0.02</sub></td><td>0.18<sup>+0.01</sup><sub>-0.01</sub></td><td>14.33<sup>+0.03</sup><sub>-0.03</sub></td></tr><tr><td rowspan="2">MPs Δτ = 1.50 Gyr</td><td rowspan="2">0.339 ± 0.008</td><td>DSED</td><td>...</td><td>1.20<sup>+0.44</sup><sub>-0.38</sub></td><td>-1.26<sup>+0.02</sup><sub>-0.03</sub></td><td>0.18<sup>+0.01</sup><sub>-0.01</sub></td><td>14.37<sup>+0.03</sup><sub>-0.03</sub></td></tr><tr><td>BaSTI</td><td>...</td><td>1.47<sup>+0.53</sup><sub>-0.46</sub></td><td>-1.26<sup>+0.03</sup><sub>-0.02</sub></td><td>0.18<sup>+0.01</sup><sub>-0.01</sub></td><td>14.35<sup>+0.03</sup><sub>-0.03</sub></td></tr></tbody></table>

metallicity can be fixed because no [Fe/H] variation is detected in this cluster.

To derive the age difference between the stellar populations, the hierarchical likelihood described in Section 2.5 with $N = 3$ is applied. The fit is carried out simultaneously to 1G, 2G, and 3G. First, we consider the primordial helium content value for all populations. In a second run, we assume a helium enhancement by a type of polluter star, changing the amount of helium for each generation, according to values computed by Milone et al. (2019): $δY_{I<sub>G,2G</sub>} = 0.010$ and $δY_{I<sub>G,3G</sub>} = 0.042$ for the 2G, and 3G, respectively (Figures 8, 10, and Table 3). We assumed the helium enhancement values from Milone et al. (2019) since they were derived using the same DSED stellar evolutionary models employed here, therefore there is compatibility. For the metallicity of NGC 6752, the corresponding canonical helium content in the DSED isochrones is 0.247, which was associated to 1G. The 2G and 3G helium

contents were assumed to be of 0.257 and 0.289, adopting the $δY$ values from Milone et al. (2019).

Table 3 and Figure 11 provide the results of isochrone fitting to the MPs. The derived distances using canonical helium and helium enhanced are $4.13 ± 0.06$ and $4.11 ± 0.08$ kpc, respectively. The latter distance determination is in agreement with the distance from the inverse Gaia DR2 parallax (Gaia Collaboration et al. 2018b; see above). We derive age differences of $Δτ_{I<sub>G,2G</sub>} = 300 ± 400$ Myr, and $Δτ_{I<sub>G,3G</sub>} = 500 ± 400$ Myr, relative to the age of 1G stars, considering that there is no helium enhancement within the GC. However, taking into account the GC helium enhancement see Milone et al. (2019), and noting that the method fits the three stellar populations simultaneously, the 1G is less old (even if its He is still canonical), and the age differences are of $Δτ_{I<sub>G,2G</sub>} = 200 ± 400$ Myr, and $Δτ_{I<sub>G,3G</sub>} = 500 ± 400$ Myr. These results
<figure><img src="image_12.png" /><figcaption>**Figure 9.** Results for the SSP analysis of NGC 6752. Left panel: CMD with the result from isochrone fitting, green line is the most probable solution, and the blue strip is the solutions within $1\sigma$. Right panel: the posterior distributions.</figcaption></figure>

<table><caption>**Table 3**
Results of Isochrone Fitting for NGC 6752 in SSP Context and MPs</caption><thead><tr><th></th><th>Y</th><th>τ<br>(Gyr)</th><th>Δτ<sub>1G,2G</sub><br>(Gyr)</th><th>Δτ<sub>2G,3G</sub><br>(Gyr)</th><th>Δτ<sub>1G,3G</sub><br>(Gyr)</th><th>[Fe/H]<br>(dex)</th><th>E(B-V)</th><th>(m - M)<sub>0</sub></th><th>(m - M)<sub>V</sub></th><th>d<sub>☉</sub><br>(kpc)</th></tr></thead><tbody><tr><td>SSP</td><td>Y(Z)<sup>a</sup></td><td>13.70<sup>+0.61</sup><sub>-0.38</sub></td><td>...</td><td>...</td><td>...</td><td>-1.49<sup>+0.05</sup><sub>-0.05</sub></td><td>0.03<sup>+0.02</sup><sub>-0.02</sub></td><td>13.07<sup>+0.03</sup><sub>-0.03</sub></td><td>13.16<sup>+0.07</sup><sub>-0.07</sub></td><td>4.11 ± 0.08</td></tr><tr><td colspan="12">MPs with Y canonical</td></tr><tr><td>1G</td><td>0.247</td><td>13.80<sup>+0.45</sup><sub>-0.40</sub></td><td>0.30<sup>+0.42</sup><sub>-0.39</sub></td><td>0.20<sup>+0.39</sup><sub>-0.38</sub></td><td>0.50<sup>+0.43</sup><sub>-0.39</sub></td><td>-1.49<sup>b</sup></td><td>0.04<sup>+0.01</sup><sub>-0.01</sub></td><td>13.08<sup>+0.02</sup><sub>-0.02</sub></td><td>13.20<sup>+0.03</sup><sub>-0.03</sub></td><td>4.13 ± 0.06</td></tr><tr><td>2G</td><td>0.247</td><td>13.50<sup>+0.39</sup><sub>-0.38</sub></td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td></tr><tr><td>3G</td><td>0.247</td><td>13.30<sup>+0.39</sup><sub>-0.38</sub></td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td></tr><tr><td colspan="12">MPs with Y enhancement</td></tr><tr><td>1G</td><td>0.247</td><td>13.50<sup>+0.39</sup><sub>-0.42</sub></td><td>0.20<sup>+0.38</sup><sub>-0.41</sub></td><td>0.30<sup>+0.37</sup><sub>-0.41</sub></td><td>0.50<sup>+0.38</sup><sub>-0.42</sub></td><td>-1.49<sup>b</sup></td><td>0.04<sup>+0.01</sup><sub>-0.01</sub></td><td>13.07<sup>+0.03</sup><sub>-0.03</sub></td><td>13.19<sup>+0.03</sup><sub>-0.03</sub></td><td>4.11 ± 0.08</td></tr><tr><td>2G</td><td>0.257</td><td>13.20<sup>+0.39</sup><sub>-0.41</sub></td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td></tr><tr><td>3G</td><td>0.289</td><td>13.00<sup>+0.41</sup><sub>-0.41</sub></td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td><td>...</td></tr></tbody></table>

**Notes.**

<sup>a</sup> Y as function of Z, defined by: $0.245 + 1.5 \times Z$.

<sup>b</sup> Fixed value from the SSP isochrone fitting.

could give hints on the possible mechanism of GC internal pollution.

It is interesting to note that, for the He enhanced populations, the result is similar to those with no He enhancement. Assuming the primordial helium for the 1G, 2G, and 3G stars, the $\chi^2$ values are 0.10, 0.13, and 0.12, respectively, resulting in a total value of 0.35. For He enhanced isochrones, the values of $\chi^2$ are 0.09, 0.14, and 0.11, for the 1G, 2G, and 3G stars, respectively and with a total of 0.34. Therefore, the fitting using He enhanced isochrones are similarly well-fit.

Even though the uncertainties on the age derivation do not take into account the differences between the stellar evolutionary models, our uncertainty determinations are of the same order of magnitude as those by Monty et al. (2018). Given that we did not propagate the uncertainties from the grid size of the parameter space, the uncertainties given here are the formal

errors from MCMC algorithm and they are larger than the ones reported by Wagner-Kaiser et al. (2017).

**5. Conclusions**

We have developed the SIRIUS code to extract the maximum information from CMDs of stellar clusters, through a detailed analysis. SIRIUS was tested in terms of synthetic data. High precision parameter derivations were obtained with sanity checks that demonstrate the good performance of the code. Small fluctuations of the solutions were found in terms of the choice of CMD colors, relative to the input parameters of the synthetic data (Figure 5). Applying a Monte Carlo spread of stars, these fluctuations increase somewhat, as can be seen in Table 2. In any case, the solution obtained is within the uncertainties and limited because of the grid resolution in the parameter space.
<figure><img src="image_13.png" /><figcaption>**Figure 10.** Corner plots for NGC 6752. Left panel: simultaneous fitting of the three stellar populations, adopting canonical helium abundance; right panel: same as in left panel, but taking into account helium abundance differences.</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>**Figure 11.** Isochrone fitting for NGC 6752. Left panel: MPs all together. Second to fourth panels: isochrone fitting to 1G, 2G, and 3G. Upper panels: Canonical helium. Lower panels: enhanced helium. The strips are the solutions within $1\sigma$.</figcaption></figure>
The SIRIUS code is applied to analyze the halo globular cluster NGC 6752 of metallicity [Fe/H] ≈ −1.49. Three stellar populations are identified, confirming previous findings by Carretta et al. (2012) from spectroscopy, and Milone et al. (2019) from photometry. The age derivation of the three stellar populations, taking into account He abundance differences from Milone et al. (2019), results to be of 200/300 ± 400 Myr between 1G and 2G and between 2G and 3G. This points to a possible interpretation of having the same mechanism producing 2G, and later the 3G.

Many authors have extensively discussed the probable candidates to produce the chemical abundance patterns of second (and subsequent) stellar populations from self-enrichment of the cluster. The main candidates are the AGB stars, and SMS, in both cases through their winds, as well as FRMSs (Decressin et al. 2007; Krause et al. 2013). All of them predict an age difference between the stellar populations.

In conclusion, given the uncertainties in the models of pollution, and the uncertainties in the age difference derived from the CMDs, it is not possible to firmly indicate a scenario for the formation of a second stellar population. The age differences derived for NGC 6752 could be compatible with the AGB scenario if only the best value determinations are taken into account. However, considering the uncertainties, the results could be compatible with all scenarios regarding the origin of MPs (SMS and FRMS), even those with no age difference. Further analyses of age differences of MPs are of great interest. In particular, within the HST Legacy survey collaboration, Nardiello et al. (2015b) derived the relative age of NGC 6352 MPs from χ² minimization isochrone fitting, assuming each of them as SSPs, and Oliveira et al. (2020) apply the methods described here to derive the ages for seven bulge globular clusters and their MPs.

We acknowledge the anonymous referee for a detailed review and helpful suggestions, which allowed us to improve the manuscript. S.O.S. acknowledges the FAPESP PhD fellowship 2018/22044-3. L.O.K. and B.B. acknowledge partial financial support from FAPESP, CNPq, and CAPES—Finance Code 001. A.P.V. acknowledges the FAPESP postdoctoral fellowship no. 2017/15893-1. R.A.P.O. acknowledges the FAPESP postdoctoral fellowship no. 2018/22181-0. D.N. acknowledges partial support by the Università degli Studi di Padova Progetto di Ateneo BIRD178590. A.P.V. and S.O.S. acknowledge the DGAPA-PAPIIT grant IG100319. Based on observations with the NASA/ESA Hubble Space Telescope, obtained at the Space Telescope Science Institute, which is operated by AURA, Inc., under NASA contract NAS 5-26555.

## ORCID iDs

S. O. Souza https://orcid.org/0000-0001-8052-969X

L. O. Kerber https://orcid.org/0000-0002-7435-8748

B. Barbuy https://orcid.org/0000-0001-9264-4417

A. Pérez-Villegas https://orcid.org/0000-0002-5974-3998

R. A. P. Oliveira https://orcid.org/0000-0002-4778-9243

D. Nardiello https://orcid.org/0000-0003-1149-3659

## References

Barbuy, B., Chiappini, C., & Gerhard, O. 2018, ARA&A, 56, 223

Bastian, N., & Lardo, C. 2018, ARA&A, 56, 83

Bedin, L. R., Piotto, G., Anderson, J., et al. 2004, ApJL, 605, L125

Carretta, E. 2019, A&A, 624, A24

Carretta, E., Bragaglia, A., Gratton, R. G., Lucatello, S., & D'Orazi, V. 2012, ApJL, 750, L14

Cordoni, G., Milone, A. P., Mastrobuono-Battisti, A., et al. 2020, ApJ, 889, 18

D'Antona, F., Caloi, V., & Tailo, M. 2018, NatAs, 2, 270

D'Antona, F., Vesperini, E., D'Ercole, A., et al. 2016, MNRAS, 458, 2122

Decressin, T., Meynet, G., Charbonnel, C., Prantzos, N., & Ekström, S. 2007, A&A, 464, 1029

Dotter, A., Chaboyer, B., Jevremović, D., et al. 2008, ApJS, 178, 89

Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306

Gaia Collaboration, Brown, A. G. A., Vallenari, A., et al. 2018a, A&A, 616, A1

Gaia Collaboration, Helmi, A., van Leeuwen, F., et al. 2018b, A&A, 616, A12

Gallart, C., Zoccali, M., & Aparicio, A. 2005, ARA&A, 43, 387

Gieles, M., Charbonnel, C., Krause, M. G. H., et al. 2018, MNRAS, 478, 2461

Gratton, R., Sneden, C., & Carretta, E. 2004, ARA&A, 42, 385

Gratton, R. G., Bragaglia, A., Carretta, E., et al. 2003, A&A, 408, 529

Gratton, R. G., Bragaglia, A., Carretta, E., et al. 2005, A&A, 440, 901

Harris, W. E. 1996, AJ, 112, 1487

Hastings, W. K. 1970, Biometrika, 57, 97

Hernandez, X., & Valls-Gabaud, D. 2008, MNRAS, 383, 1603

Hidalgo, S. L., Pietrinferni, A., Cassisi, S., et al. 2018, ApJ, 856, 125

Hogg, D. W., & Foreman-Mackey, D. 2018, ApJS, 236, 11

Kerber, L. O., Libralato, M., Souza, S. O., et al. 2019, MNRAS, 484, 5530

Kerber, L. O., Nardiello, D., Ortolani, S., et al. 2018, ApJ, 853, 15

Kerber, L. O., & Santiago, B. X. 2005, A&A, 435, 77

Kerber, L. O., Santiago, B. X., & Brocato, E. 2007, A&A, 462, 139

Kerber, L. O., Santiago, B. X., Castro, R., & Valls-Gabaud, D. 2002, A&A, 390, 121

Krause, M., Charbonnel, C., Decressin, T., Meynet, G., & Prantzos, N. 2013, A&A, 552, A121

Kroupa, P. 2001, MNRAS, 322, 231

Lagioia, E. P., Milone, A. P., Marino, A. F., Cordoni, G., & Tailo, M. 2019, AJ, 158, 202

Lee, J.-W. 2019, ApJ, 883, 166

Lee, Y. W., Joo, J. M., Sohn, Y. J., et al. 1999, Natur, 402, 55

Lindegren, L., Hernández, J., Bombrun, A., et al. 2018, A&A, 616, A2

Martocchia, S., Cabrera-Ziri, I., Lardo, C., et al. 2018, MNRAS, 473, 2688

Martocchia, S., Dalessandro, E., Lardo, C., et al. 2019, MNRAS, 487, 5324

Metropolis, N., Rosenbluth, A. W., Rosenbluth, M. N., Teller, A. H., & Teller, E. 1953, JChPh, 21, 1087

Milone, A. P., Marino, A. F., Bedin, L. R., et al. 2019, MNRAS, 484, 4046

Milone, A. P., Marino, A. F., Piotto, G., et al. 2013, ApJ, 767, 120

Milone, A. P., Marino, A. F., Renzini, A., et al. 2018, MNRAS, 481, 5098

Milone, A. P., Piotto, G., Renzini, A., et al. 2017, MNRAS, 464, 3636

Monteiro, H., Dias, W. S., & Caetano, T. C. 2010, A&A, 516, A2

Monty, S., Puzia, T. H., Miller, B. W., et al. 2018, ApJ, 865, 160

Nardiello, D., Libralato, M., Piotto, G., et al. 2018, MNRAS, 481, 3382

Nardiello, D., Milone, A. P., Piotto, G., et al. 2015a, A&A, 573, A70

Nardiello, D., Piotto, G., Milone, A. P., et al. 2015b, MNRAS, 451, 312

Naylor, T., & Jeffries, R. D. 2006, MNRAS, 373, 1251

Oliveira, R. A. P., Souza, S. O., Kerber, L. O., et al. 2020, arXiv:2001.08611

Ortolani, S., Cassisi, S., & Salaris, M. 2017, Galax, 5, 28

Ortolani, S., Held, E. V., Nardiello, D., et al. 2019, A&A, 627, A145

Osborn, W. 1971, Obs, 91, 223

Pedregosa, F., Varoquaux, G., Gramfort, A., et al. 2011, Journal of Machine Learning Research, 12, 2825

Perren, G. I., Vázquez, R. A., & Piatti, A. E. 2015, A&A, 576, A6

Pietrinferni, A., Cassisi, S., Salaris, M., & Castelli, F. 2006, ApJ, 642, 797

Piotto, G., Milone, A. P., Bedin, L. R., et al. 2015, AJ, 149, 91

Piotto, G., Villanova, S., Bedin, L. R., et al. 2005, ApJ, 621, 777

Planck Collaboration, Ade, P. A. R., Aghanim, N., et al. 2016, A&A, 594, A13

Ramírez-Siordia, V. H., Bruzual, G., Cervantes Sodi, B., & Bitsakis, T. 2019, MNRAS, 486, 5567

Renzini, A. 2013, MmSAI, 84, 162

Renzini, A., D'Antona, F., Cassisi, S., et al. 2015, MNRAS, 454, 4197

Stenning, D. C., Wagner-Kaiser, R., Robinson, E., et al. 2016, ApJ, 826, 41

VandenBerg, D. A., Brogaard, K., Leaman, R., & Casagrand, E. L. 2013, ApJ, 775, 134

von Hippel, T., Jefferys, W. H., Scott, J., et al. 2006, ApJ, 645, 1436

Wagner-Kaiser, R., Sarajedini, A., von Hippel, T., et al. 2017, MNRAS, 468, 1038

Wagner-Kaiser, R., Stenning, D. C., Sarajedini, A., et al. 2016, MNRAS, 463, 3768