Origin of the system of globular clusters in the Milky Way

TL;DR

This study uses Gaia data to analyze the origins of Milky Way globular clusters, linking their properties to known merger events and distinguishing between in situ formation and accretion origins.

Contribution

It combines Gaia kinematic data with cluster ages to identify the likely formation scenarios and merger associations of globular clusters in the Milky Way.

Findings

Approximately 40% of clusters formed in situ.

35% of clusters are associated with known merger events.

Age-metallicity relations differ by progenitor, supporting association accuracy.

Abstract

$Context$. The assembly history experienced by the Milky Way is currently being unveiled thanks to the data provided by the $Gaia$ mission. It is likely that the globular cluster system of our Galaxy has followed a similarly intricate formation path. $Aims$. To constrain this formation path, we explore the link between the globular clusters and the known merging events that the Milky Way has experienced. $Methods$. To this end, we combined the kinematic information provided by $Gaia$ for almost all Galactic clusters, with the largest sample of cluster ages available after carefully correcting for systematic errors. To identify clusters with a common origin we analysed their dynamical properties, particularly in the space of integrals of motion. $Results$. We find that about 40% of the clusters likely formed in situ. A similarly large fraction, 35%, appear to be possibly associated to known merger events, in particular to $Gaia$-Enceladus (19%), the Sagittarius dwarf galaxy (5%), the progenitor of the Helmi streams (6%), and to the Sequoia galaxy (5%), although some uncertainty remains due to the degree of overlap in their dynamical characteristics. Of the remaining clusters, 16% are tentatively associated to a group with high binding energy, while the rest are all on loosely bound orbits and likely have a more heterogeneous origin. The resulting age-metallicity relations are remarkably tight and differ in their detailed properties depending on the progenitor, providing further confidence on the associations made. $Conclusions$. We provide a table listing the likely associations. Improved kinematic data by future Gaia data releases and especially a larger, systematic error-free sample of cluster ages would help to further solidify our conclusions.

Tables

Table 1. Table 1: List of Galactic GCs and associated progenitors. M-D stands for Main-Disc; M-B stands for Main-Bulge; G-E stands for Gaia -Enceladus; Sag stands for Sagittarius dwarf; H99 stands for Helmi Streams; Seq stands for Sequoia galaxy; L-E stands for unassociated Low-Energy; H-E stands for unassociated High-Energy. Finally XXX mark clusters with no available kinematics. Most of these are bulge GCs too extincted to be observed by Gaia.

GC	Progenitor	GC	Progenitor	GC	Progenitor	GC	Progenitor
NGC 104	M-D	NGC 5927	M-D	HP 1	M-B	GLIMPSE 2	XXX
NGC 288	G-E	NGC 5946	L-E	NGC 6362	M-D	NGC 6584	H-E
NGC 362	G-E	BH 176	M-D	Liller 1	XXX	NGC 6624	M-B
Whiting 1	Sag	NGC 5986	L-E	NGC 6380	M-B	NGC 6626	M-B
NGC 1261	G-E	Lynga 7	M-D	Terzan 1	M-B	NGC 6638	M-B
Pal 1	M-D	Pal 14	H-E	Ton 2	L-E	NGC 6637	M-B
AM 1	H-E	NGC 6093	L-E	NGC 6388	M-B	NGC 6642	M-B
Eridanus	H-E	NGC 6121	L-E	NGC 6402	L-E	NGC 6652	M-B
Pal 2	G-E?	NGC 6101	Seq/G-E	NGC 6401	L-E	NGC 6656	M-D
NGC 1851	G-E	NGC 6144	L-E	NGC 6397	M-D	Pal 8	M-D
NGC 1904	G-E	NGC 6139	L-E	Pal 6	L-E	NGC 6681	L-E
NGC 2298	G-E	Terzan 3	M-D	NGC 6426	H-E	GLIMPSE 1	XXX
NGC 2419	Sag	NGC 6171	M-B	Djorg 1	G-E	NGC 6712	L-E
Ko 2	XXX	1636-283	M-B	Terzan 5	M-B	NGC 6715	Sag
Pyxis	H-E	NGC 6205	G-E	NGC 6440	M-B	NGC 6717	M-B
NGC 2808	G-E	NGC 6229	G-E	NGC 6441	L-E	NGC 6723	M-B
E 3	H99?	NGC 6218	M-D	Terzan 6	M-B	NGC 6749	M-D
Pal 3	H-E	FSR 1735	L-E	NGC 6453	L-E	NGC 6752	M-D
NGC 3201	Seq/G-E	NGC 6235	G-E	UKS 1	XXX	NGC 6760	M-D
Pal 4	H-E	NGC 6254	L-E	NGC 6496	M-D	NGC 6779	G-E
Ko 1	XXX	NGC 6256	L-E	Terzan 9	M-B	Terzan 7	Sag
NGC 4147	G-E	Pal 15	G-E?	Djorg 2	M-B	Pal 10	M-D
NGC 4372	M-D	NGC 6266	M-B	NGC 6517	L-E	Arp 2	Sag
Rup 106	H99?	NGC 6273	L-E	Terzan10	G-E	NGC 6809	L-E
NGC 4590	H99	NGC 6284	G-E	NGC 6522	M-B	Terzan 8	Sag
NGC 4833	G-E	NGC 6287	L-E	NGC 6535	L-E/Seq	Pal 11	M-D
NGC 5024	H99	NGC 6293	M-B	NGC 6528	M-B	NGC 6838	M-D
NGC 5053	H99	NGC 6304	M-B	NGC 6539	M-B	NGC 6864	G-E
NGC 5139	G-E/Seq	NGC 6316	M-B	NGC 6540	M-B	NGC 6934	H-E
NGC 5272	H99	NGC 6341	G-E	NGC 6544	L-E	NGC 6981	H99
NGC 5286	G-E	NGC 6325	M-B	NGC 6541	L-E	NGC 7006	Seq
AM 4	XXX	NGC 6333	L-E	2MS-GC01	XXX	NGC 7078	M-D
NGC 5466	Seq	NGC 6342	M-B	ESO-SC06	G-E	NGC 7089	G-E
NGC 5634	H99/G-E	NGC 6356	M-D	NGC 6553	M-B	NGC 7099	G-E
NGC 5694	H-E	NGC 6355	M-B	2MS-GC02	XXX	Pal 12	Sag
IC 4499	Seq	NGC 6352	M-D	NGC 6558	M-B	Pal 13	Seq
NGC 5824	Sag	IC 1257	G-E	IC 1276	M-D	NGC 7492	G-E
Pal 5	H99?	Terzan 2	M-B	Terzan12	M-D	Crater	H-E
NGC 5897	G-E	NGC 6366	M-D	NGC 6569	M-B	FSR 1716	M-D
NGC 5904	H99/G-E	Terzan 4	M-B	BH 261	M-B	FSR 1758	Seq

Equations

$$Z(t)=-p\ln \mu (t)=-p\ln \frac{t_f-t}{t_f-t_i}\mathrm{for}t\ge t_i,$$

Full text

*Aims. *To constrain this formation path, we explore the link between the globular clusters and the known merging events that the Milky Way has experienced.

*Methods. *To this end, we combined the kinematic information provided by Gaia for almost all Galactic clusters, with the largest sample of cluster ages available after carefully correcting for systematic errors. To identify clusters with a common origin we analysed their dynamical properties, particularly in the space of integrals of motion.

*Results. *We find that about 40% of the clusters likely formed in situ. A similarly large fraction, 35%, appear to be possibly associated to known merger events, in particular to Gaia-Enceladus (19%), the Sagittarius dwarf galaxy (5%), the progenitor of the Helmi streams (6%), and to the Sequoia galaxy (5%), although some uncertainty remains due to the degree of overlap in their dynamical characteristics. Of the remaining clusters, 16% are tentatively associated to a group with high binding energy, while the rest are all on loosely bound orbits and likely have a more heterogeneous origin. The resulting age–metallicity relations are remarkably tight and differ in their detailed properties depending on the progenitor, providing further confidence on the associations made.

*Conclusions. *We provide a table listing the likely associations. Improved kinematic data by future Gaia data releases and especially a larger, systematic error-free sample of cluster ages would help to further solidify our conclusions.

1 Introduction

According to the $\Lambda$CDM cosmological paradigm, structure formation proceeds bottom-up, as small structures merge together to build up the larger galaxies we observe today. The Milky Way is a prime example of this formation mechanism, as first demonstrated by the discovery of the Sagittarius dwarf spheroidal galaxy in the process of disruption (Ibata et al., 1994), then in halo stellar streams crossing the solar neighbourhood (Helmi et al., 1999), and more recently by the discovery of stellar debris from Gaia-Enceladus, revealing the last significant merger experienced by our Galaxy (Helmi et al. 2018, see also Belokurov et al. 2018).

As a natural result of such events, not only field stars but also globular clusters (GCs) may have been accreted (Peñarrubia et al. 2009). Starting with Searle & Zinn (1978) there has been a quest to understand which of the approximately 150 GCs hosted by the Galaxy actually formed in situ and which formed in different progenitors that were only later accreted. Recently, the availability of precise relative ages (with formal errors of $\lesssim 1$ Gyr; e.g. Marín-Franch et al., 2009; VandenBerg et al., 2013) and homogeneous metallicity measurements (Carretta et al. 2009) led to the discovery that the age–metallicity relation (AMR) of Galactic GCs is bifurcated (Marín-Franch et al. 2009; Forbes & Bridges 2010; Leaman et al. 2013). Although limited, kinematic information (e.g. Dinescu et al. 1997, 1999; Massari et al. 2013) nonetheless helped to reveal that the metal-poor branch of young GCs have halo-like kinematics (and are therefore more likely to be accreted), whereas GCs in the young and metal-rich branch have disc-like kinematics, and are consistent with having formed in situ (see also Recio-Blanco 2018).

With the advent of the second data release (DR2) of the Gaia mission (Gaia Collaboration et al. 2018), we have, for the first time, full six-dimensional phase space information for almost all of the Galactic GCs (Gaia Collaboration et al., 2018; Vasiliev, 2019a). Therefore, it is timely to revisit the origin of the Galactic GC system. The goal of this Letter is to use this information to provide a more complete picture of which GCs formed outside our Galaxy and in which progenitor (among those currently known or yet to be discovered). The main result of this analysis is given in Table 1, which lists all the Galactic GCs and the progenitors they have likely been associated with. These associations may be seen to reflect the best of our current understanding, although some are only tentative. We may expect to be able to draw firmer conclusions with improved GC age datasets and larger samples of field stars with full-phase space kinematics.

2 The dataset: dynamical properties, ages, and metallicities

We put together a dataset of $151$ clusters with full 6D phase-space information of Galactic GCs known based on the compilations by Gaia Collaboration et al. (2018) and Vasiliev (2019a) (for more details, see Appendix A.1). We used the AGAMA package (Vasiliev 2019b) with the McMillan (2017) potential to compute GCs orbital parameters like the apocenter (apo), maximum height from the disc ($Z_{\rm max}$), and eccentricity (ecc). We also computed the orbital circularity as circ = $L_{Z}/L_{Z,{\rm circ}}$, where $L_{Z,{\rm circ}}$ is the angular momentum of a circular orbit with the cluster energy, which thus takes extreme values +1 or -1 for co-planar circular prograde or retrograde orbits, respectively.

In this work we study the AMR for a subsample of GCs that has a homogeneous set of ages and metallicities. This sample includes the catalogue by VandenBerg et al. (2013), who provide absolute ages and uncertainty for 55 GCs, and objects from the compilation by Forbes & Bridges (2010) who gathered relative age estimates from Salaris & Weiss (1998), De Angeli et al. (2005), Marín-Franch et al. (2009), which we add after checking for systematic offsets with respect to VandenBerg et al. (2013). Our final sample consists of 69 GCs with ages and metallicities (see Appendix A.2 for details and a discussion of the care needed when handling GC ages).

3 Assignment of clusters

Figure 1 shows the AMR for the clusters in our sample, colour-coded according to various dynamical properties, namely apo, circ, $Z_{\rm max}$ and ecc. It is immediately clear that the clusters located on the young and metal-rich branch of the AMR are dynamically different from those on the young and metal-poor branch. Young and metal-rich GCs typically do not reach high altitudes above the Galactic plane ($Z_{\rm max}$), have smaller apocentres, and tend to have lower eccentricities. As already recognised in the literature (e.g. Leaman et al. 2013), these are the clusters formed in situ, either in the disc or in the bulge, in what we hereafter refer to as the Main Progenitor. For the first time we can recognise from Fig. 1 some old metal-poor GCs with orbital properties characteristic of the young metal-rich branch, which thus would also belong to the Main Progenitor.

3.1 In situ clusters

We therefore use the dynamical properties of the GCs that populate the young and metal-rich branch of the AMR of Fig. 1, and which were likely born in the Milky Way, to define simple criteria to identify Main Progenitor clusters:

• •

Bulge clusters: those placed on highly bound orbits (with apo ¡ 3.5 kpc). There are 36 GCs selected in this way111Nuclear clusters of dwarf galaxies accreted long ago could also end up on bulge-like orbits due to dynamical friction. The current uncertainty on the proper motions of bulge GCs (they are typically highly extincted and only few stars are detected by Gaia), and the lack of an age estimate for most of them, prevents us from investigating this possibility..

• •

Disc clusters satisfy: $i$) $Z_{\rm max}<5$ kpc and $ii$) ${\sf circ}>0.5$. While this does not guarantee a “pure” disc sample, and there may be a small amount of contamination, the effectiveness of these criteria is supported by the fact that the vast majority of the selected clusters tend to describe an AMR qualitatively similar to that found by Leaman et al. (2013), except for two clusters located on the young and metal-poor branch (NGC6235 and NGC6254). The dynamical parameters of these clusters are in the extremes of those characteristic of the Main Progenitor, and therefore we exclude them. This leaves a total of 26 Disc clusters.

We note that for [Fe/H] $<-1.5$, Main Progenitor clusters are older than the average. The $62$ in situ clusters are listed in Table1 as Main-Disc (M-D) or Main-Bulge (M-B).

3.2 Accreted clusters

We now analyse the remaining clusters looking for a common association with the progenitors of known merger events experienced by the Milky Way. To do so, we investigate the integrals of motion (IOM) space defined by $E$, $L_{Z}$, and $L_{perp}$, the latter being the angular momentum component perpendicular to $L_{Z}$, which, despite not being fully conserved in an axisymmetric potential like that of the Milky Way, still helps in discriminating groups of stars (or clusters) with similar origin (Helmi & de Zeeuw 2000), also if the potential has varied with time (see Peñarrubia et al., 2006; Gómez et al., 2013). In particular, we use the (known) extent of each progenitor stellar debris in the IOM space to provisionally identify associated GCs.

3.2.1 The Sagittarius dwarf spheroidal galaxy

The Sagittarius dwarf spheroidal galaxy constituted the first discovery of a merger with the Galaxy (Ibata et al., 1994). By exploiting numerical models that accurately reproduce the position and radial velocity of stars belonging to the Sagittarius streams, Law & Majewski (2010) provided a list of candidate GCs that could have been associated to the dwarf. This list has been recently refined by adding the information on their proper motions as measured from HST and Gaia observations (Massari et al. 2017; Sohn et al. 2018), and currently includes 6 GCs, namely M54, Arp 2, Pal 12, Terzan 7, Terzan 8, and NGC2419.

These six clusters describe a well-defined subgroup in IOM space: $i$) $3700\!<\!L_{perp}\!<\!6200$ km/s kpc and $ii$) $0<L_{Z}<3000$ km/s kpc, as seen in Fig. 2. Two more clusters are found in this region of IOM space: NGC 5824 and Whiting 1, which had previously been tentatively associated with the dwarf by Bellazzini et al. (2003) and Law & Majewski (2010), respectively.

3.2.2 The progenitor of the Helmi streams

Recently, Koppelman et al. (2019) used Gaia data to characterise the progenitor of the Helmi streams (hereafter H99, Helmi et al., 1999). According to these authors and based on their dynamical properties (and a comparison with a numerical simulation that best reproduces the properties of the streams), the seven GCs shown in orange in Fig. 2 could be associated to this object. Interestingly, these seven clusters were shown to follow a tight and low-normalisation AMR.

To explore whether or not additional members could exist, we started from the dynamical criteria suggested in their work, and revisited the location of the selection boundaries while requiring consistency with the AMR. The following criteria: $i$) $350<L_{Z}<3000$ km/s kpc, $ii$) $1000<L_{perp}<3200$ km/s kpc, and $iii$) $E<-1.0\times 10^{5}$ km2/s2 seem to be the most appropriate. Allowing for lower values of $L_{Z}$ leads to the inclusion of very old clusters (not consistent with the AMR of the core members) whereas increasing the upper limit leads to including two disc clusters lacking an age estimate (Pal 1 and BH 176). Increasing the $E$ limit would mean including a cluster with ${\sf apo}$ ¿100 kpc, whereas the typical value for the core members is $\sim 30$ kpc. Finally, the limits on $L_{perp}$ are given by Sagittarius clusters on one side and by the consistency with the AMR on the other.

Out of the ten GCs that would be associated to H99 according to the above criteria, there is no age information for three of them: Rup 106, E3, and Pal 5, and therefore we consider them to be tentative members (orange open symbols in Fig. 2). We note that Pal 5 and E3 have the lowest ecc ($\sim 0.2$) in the set, with E3 having the lowest $Z_{\rm max}$ ($\sim 7$ kpc), and Rup 106 the largest apo ($\sim 34$ kpc). Also in comparison to the H99 stars (see Koppelman et al., 2019, for details), E3 and Pal 5 have more extreme orbital properties, with Rup 106 being on a looser orbit. However, this does not necessarily preclude membership since GCs are expected to be less bound than the stars.

3.2.3 Gaia-Enceladus

To look for GCs associated to Gaia-Enceladus (G-E hereafter, Helmi et al., 2018), we directly compare the distribution in IOM space of GCs to that of field stars with 6D kinematics from Gaia, as shown in Fig. 3.

Based on this comparison, we associate clusters to G-E according to the following criteria: $i$) $-800<L_{Z}<620$ km/s kpc, $ii$) $-1.86\times 10^{5}<E<-0.9\times 10^{5}$ km2/s2, and $iii$) $L_{perp}<3500$ km/s kpc. This selection associates 28 GCs to G-E. With the exception of NGC 7492 (its apo being $\sim 28$ kpc), all of them have apocenter values $<25$ kpc, as reported for G-E stars (see Deason et al., 2018). The resulting AMR is remarkably tight (see below), and this tightness can be used to explore the energy boundaries. By decreasing the lower limit of $E$, a very old globular cluster enters the selection which is significantly off the AMR described by the other members. This energy limit is slightly higher than that used by Myeong et al. (2018a, who adopt the same Galactic potential as ours) to define the clusters belonging to the progenitor of the same accretion event, that they call Gaia-Sausage (see also Belokurov et al. 2018). Moving the upper limit to $E=-1.1\times 10^{5}$ km2/s2 excludes Pal 2 and Pal 15, which we therefore consider to be tentative members (open symbols in Fig. 2).

Some of the associated clusters are located in regions of the IOM space that are shared by stellar debris of different progenitors. Therefore, the associations are more uncertain and deserve further discussion (see e.g. Borsato et al., 2019). For example, two clusters (NGC5904 and NGC5634) lie near the region occupied by H99 debris. They have ecc $\sim 0.8$, and their $Z_{\rm max}$ and apo are somewhat larger but not inconsistent with those typical of G-E clusters. While NGC5634 has no age estimate, the age of NGC5904 ($11.5$ Gyr old) is consistent with both the AMRs. We therefore consider them as tentative members to both progenitors.

Also of particular interest is the overdensity of stars seen in Fig. 3 at $L_{Z}\sim-3000$ km/s kpc and $E\sim-10^{5}$ km2/s2, which has been associated to G-E by Helmi et al. (2018) because of its resemblance to a feature seen in numerical simulations of a merger event with similar characteristics to G-E (Villalobos & Helmi 2008). Two GCs (NGC3201 and NGC6101) are located in this region of IOM space, and therefore we mark them as tentative members but discuss them further in the following section.

With our selection criteria, $\omega$-Cen222for which we used the metallicity of the most metal-poor and oldest population, see Bellini et al. 2017. is the cluster with the highest binding energy among those associated to G-E. This is in good agreement with the idea that this cluster is in reality the remnant of the nuclear star cluster of an accreted dwarf (e.g. Bekki & Freeman 2003) as suggested by its peculiar chemistry.

After these considerations, we are left with $26$ possible members of G-E (and $6$ tentative ones). Although this number is large, it is consistent within the scatter with the relation between the number of GCs and host halo mass (van Dokkum et al., 2017), given the mass estimate of $6\times 10^{10}$ M*⊙* from Helmi et al. (2018).

3.2.4 Sequoia

Recently, Myeong et al. (2019) proposed the existence of merger debris from a galaxy named Sequoia that would have been accreted about 9 Gyr ago. Based on clustering algorithms performed over a scaled action space, these authors found five GCs likely associated to this system. We find seven GCs to be possibly associated when considering a selection box in $E$ and $L_{Z}$ corresponding to the stars from Sequoia according to Myeong et al. (2019), namely $-3700<L_{Z}<-850$ km/s kpc and $-1.5\times 10^{5}<E<-0.7\times 10^{5}$ km2/s2. Three of these GCs are in common, namely FSR 1758, NGC 3201, and NGC 6101. The other four (IC 4499, NGC 5466, NGC 7006 and Pal 13) were excluded by Myeong et al. because of their slightly larger eccentricity ($\langle$ecc$\rangle\sim 0.75$, compared to their initial estimate of $\sim 0.6$). Three clusters have known ages and these follow a low-normalisation AMR similar to the H99 GCs, which is consistent with the low stellar mass estimated for Sequoia (Myeong et al. 2019).

As mentioned earlier, the Sequioa IOM selection has some overlap with the arch-like overdensity ascribed to G-E debris by Helmi et al. (2018). Discerning which is the actual progenitor of NGC3201 and NGC6101 is therefore not possible, and so in Table 1 we link them to both systems. Nonetheless, there may be a slight preference for Sequoia given their ages and metallicities. On the other hand, Myeong et al. (2019) associated $\omega$-Cen and NGC6535 to Sequoia because of their location in IOM space. However, these two GCs follow a much higher AMR, typical of clusters from more massive progenitors. For this reason, we prioritise the association of $\omega$-Cen to G-E and of NGC6535 to one of the groups described below, though we acknowledge both the interpretations in Table 1.

3.2.5 The remaining clusters

We were not able to associate 36 of the 151 GCs with full phase-space information to known merger events. From their distribution in the IOM space (Fig. 2), it is clear that at least a significant fraction of them (25) could tentatively be part of a structure at low energy, with $E<-1.86\times 10^{5}$ km2/s2, low $L_{perp}$, and with $L_{Z}\sim 0$ km/s kpc (pink symbols in Fig. 2); we label these L-$E$.

The remaining 11 GCs all have high energy ($E>-1.5\times 10^{5}$ km2/s2, in cyan in Fig. 2), but span a very large range in $L_{Z}$ and $L_{perp}$. Therefore, they cannot have a common origin. Most likely instead, they have been accreted from different low-mass progenitors which have not contributed debris (field stars) to the Solar vicinity (as otherwise we would have identified corresponding overdensities in Fig. 3). For convenience only, we use a single label for all these objects (H-$E$, for high energy) in Table 1. Upcoming datasets, especially of field stars with full phase-space information across the Galaxy, could be key to understanding their origin given their heterogeneous properties.

3.3 Age–metallicity relation

The various panels in Fig. 4 show the AMR of the clusters associated to the different structures discussed so far, colour-coded as in Fig. 2. Although consistency with the AMR of each group is checked after the IOM selection, it is quite remarkable that the dynamical identification of associations of GC results in AMRs that are all well-defined and depict different shapes or amplitudes.

The clusters of the Main progenitor constitute the largest group and have the highest normalisation, that is the most metal-poor oldest clusters were born in the Galaxy itself. The G-E AMR is remarkably tight and has a high normalisation, though as expected this is not as high as that of the Main progenitor. Similarly, the L-$E$ group depicts a reasonably coherent AMR, with a high normalisation, which seems even higher than that of G-E clusters, thus possibly suggesting the presence of yet-to-be-discovered merger debris located preferentially in the Galactic bulge and originating in a more massive object. On the other hand, the AMR of the H99 members is remarkably low in terms of normalisation, and this is consistent with the fact that this progenitor is less massive ($M_{\star}\sim 10^{8}$ M*⊙*, Koppelman et al. 2019).

We can describe the various AMRs with a leaky-box chemical evolution model (Prantzos, 2008; Leaman et al., 2013; Boecker et al., 2019), where the metallicity of the system evolves as

[TABLE]

where $\mu(t)=M_{g}(t)/M_{g}(0)$ is the gas fraction, and $p$ the (effective) yield. We obtained this expression by assuming a constant star formation rate starting at time $t_{i}$ after the Big Bang and ending at time $t_{f}$, which we take to be the time of accretion (which is constrained by estimates in the literature and which we took to range from 3.2 Gyr for Sequioa to 5.7 Gyr for Sagittarius). For the yield $p$ we assumed the dependence on $M_{\star}$ derived by Dekel & Woo (2003) for star-forming dwarf galaxies (see also Prantzos, 2008). The time $t_{i}$ is a free parameter which we varied for each progenitor to obtain a reasonable description of the observed points. The only constraint we apply is that more massive progenitors should start forming stars earlier, which led to $t_{i}$ values in the range of 0.5 Gyr (for the Main progenitor) to 1.1 Gyr (for Sequioa). The resulting curves for each progenitor are shown in the panels of Fig. 4 with green solid lines. We stress that these curves do not represent fits, but merely show that a simple leaky box chemical evolution model is reasonably adequate to describe the AMRs found for each set of clusters associated to the individual progenitors. The scatter around each curve (typically $<0.5$ Gyr) is qualitatively consistent with the individual age uncertainties. This might indicate that for the clusters lacking an age uncertainty, our assumed value of $0.75$ Gyr (see Appendix A.2) was too conservative, especially when considering relative ages rather than absolute ones.

4 Summary and Conclusions

In this Letter we exploit the complete kinematic information for 151 Galactic GCs, in combination with metallicity and homogeneous age estimates for a subset of 69 GCs. Our goal is to elucidate which GCs formed in situ and which could have been accreted, associating the latter to a particular progenitor based on their dynamical properties and, where needed, on the shape of the AMR.

We found that 62 GCs likely formed in the Milky Way, and we separated them into disc and bulge clusters based on their orbital parameters. Among the accreted clusters, we assessed their possible associations with the progenitor of four known merger events: Gaia-Enceladus, the Sagittarius dwarf, the Helmi streams, and the Sequoia galaxy. We identified $26$ (and an additional 6 tentative) GCs associated to Gaia-Enceladus. This large number as well as the high normalisation of their AMRs are consistent with Gaia-Enceladus being the most massive among these four objects. According to our findings, $\omega$-Centauri would be its nuclear star cluster. We further identified eight clusters associated with the Sagittarius dwarf, possibly ten clusters to the progenitor of the Helmi streams, and a plausible seven to Sequoia. Despite not being very populated, the AMRs of these two groups are consistent with the lower literature estimates for their mass. There is an inherent uncertainty to these assignments, because debris from different progenitors partly overlaps in IOM space, as in the case of Sequoia and G-E, and to a lesser degree for G-E and the Helmi streams.

Interestingly, the specific frequency per unit of galaxy mass (TN) as defined by Zepf, & Whitmore (1993) for each progenitor closely follows the relation reported in Zaritsky et al. (2016) for a sample of galaxies in the same stellar mass range. For example, we find $\log T_{N}$ = [1.8, 2.0, 1.2, 2.1] for G-E, the Helmi streams, Sgr, and Sequioa, while the expected values using Zaritsky et al. (2016) are respectively, [1.8, 2.2, 1.8, 2.4], which is well within the observed scatter of the relation ($\sigma(\log T_{N})\sim 0.7$). The values of $\log T_{N}$ were computed using stellar mass estimates of [5, 1, 5]$\times 10^{8}M_{\odot}$ for G-E, the Helmi streams, and Sagittarius. These were derived with different methods, using chemistry in the first case (Helmi et al. 2018) and N-body models for the latter two (Koppelman et al. 2019; Dierickx, & Loeb 2017). In the case of Sequioa, we used the mass estimated by Myeong et al. (2019) on the basis of the specific frequency, and it is therefore less surprising to find good agreement (although our estimate uses a different number of associated clusters).

The 36 clusters that we have not associated to known debris can be split in two groups based on orbital energy. While the class of GCs with low binding energy is very heterogeneous and likely has several sites of origin, the low-energy (highly bound) group with 25 tentative members is relatively highly clustered in its dynamical properties and shows a reasonably tight and high-normalisation AMR, possibly suggesting the presence of debris towards the Galactic bulge from a large hitherto unknown galaxy. This finding is consistent with the conclusion that another significant accretion event is required to explain the overall AMR of Galactic GCs; for example, a merger with a “Kraken”-like galaxy (Kruijssen et al., 2019). We note however, that most of the clusters reported in Kruijssen et al. (2019) as possible members of Kraken are not dynamically coherent, since three are in common with the H99 GCs, seven with G-E, two with Sequoia, and one with the Main Progenitor. Therefore, although the L-$E$ group is Kraken-like, the associated clusters are different. It transpires that taking into account the dynamical properties is fundamental to establishing the origin of the different GCs of our Galaxy.

The next data release of the Gaia mission will provide improved astrometry and photometry for all the GCs, as well as for a much larger sample halo stars. This will be crucial to achieving a complete and accurate sample of GCs with absolute ages. Moreover it will lead to a better understanding of the debris of the known progenitors, and possibly to the discovery of new ones. The combination of these factors will result in significant progress in the field and allow to better pin down tentative associations and possibly also to assess under which conditions the different clusters formed, such as for example whether formation was prior to or during the different merger events.

Acknowledgements.

We thank the anonymous referee for comments and suggestions which improved the quality of our paper. DM, HHK and AH acknowledge financial support from a Vici grant from NWO. This work has made use of data from the European Space Agency (ESA) mission *Gaia *(http://www.cosmos.esa.int/gaia), processed by the *Gaia *Data Processing and Analysis Consortium (DPAC, http://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has been provided by national institutions, in particular the institutions participating in the *Gaia *Multilateral Agreement.

Appendix A Details of the sample of globular clusters

A.1 Kinematical properties

To put together the sample of GCs, we started from the 75 GCs analysed by Gaia Collaboration et al. (2018), who combined the Gaia measured proper motions with distances and line-of-sight velocities available from the compilation by Harris (1996, 2010 edition). We then added the data for the remaining GCs from Vasiliev (2019a), who determined the 6D coordinates combining Gaia measurements with line-of-sight velocities also from Baumgardt et al. (2019).

We then transformed the observed measurements of the kinematics of the resulting catalogue of $151$ clusters to the Galactocentric reference frame. To this end we assumed a Local Standard of Rest velocity $V_{\rm LSR}=232.8$ km/s (McMillan, 2017), a solar motion $(U,V,W)=(11.1,12.24,7.25)$ km/s (Schönrich et al., 2010), and a distance from the Sun to the Galactic centre of $R_{0}=8.2$ (McMillan, 2017).

A.2 Homogeneous cluster ages

Many methods have been applied to determine the absolute age of a GC. Photometric errors, poor calibration, and uncertainties on the cluster distance and reddening however, can affect such age estimates. To overcome these, it has often been preferred to determine relative ages (e.g. Buonanno et al. 1989; Bono et al. 2010; Massari et al. 2016), although these then need to be calibrated to some absolute scale (e.g. Marín-Franch et al. 2009). This explains why available age compilations of GCs in the literature are so heterogeneous, and therefore dangerous to blindly combine together as different methods can result in systematic differences that amount to several gigayears.

Figure 5 highlights the need for care when combining different clusters age datatsets. The top panel shows the difference in metallity between the VandenBerg et al. (2013) estimates based on the spectroscopic scale of Carretta et al. (2009), in comparison to those of Forbes & Bridges (2010). This difference is constant for [Fe/H]$\lesssim-1.1$, above which it rises steeply with metallicity to about +0.5 dex. The bottom panel of Fig. 5 shows the effect of this trend on age, namely that the clusters with [Fe/H] $\gtrsim-1.1$ (red symbols) appear to be systematically older by $\sim 2$ Gyr. When excluding these clusters, the mean difference between the age estimates is $\Delta t$ = 0.08 Gyr, with a spread of $\sigma_{t}$=0.75 Gyr (r.m.s = 0.14 Gyr), and thus fairly consistent with zero. Therefore, we only consider clusters from the Forbes & Bridges (2010) sample with [Fe/H] $\lesssim-1.1$, and assign an uncertainty to their age estimates (which lack errors) that equals the observed spread around the mean difference ($\sigma_{t}=0.75$ Gyr). As for metallicities, we adopted the scale from Carretta et al. (2009).

We also studied the estimates reported in Rosenberg et al. (1999), Dotter et al. (2011), Roediger et al. (2014), but found either poorly constrained values (with ages as old as 15 Gyr and very large uncertainties), or no new entries. Moreover, we decided not to include age estimates performed on single objects because of the impossibility of having systematic effects under control. The only exception is NGC5634, estimated to be as old as NGC4372 by Bellazzini et al. (2002), and for which we used NGC4372 age estimate by VandenBerg et al. (2013).