© The Authors 2025

1 Introduction

According to the Λ cold dark matter cosmological model, galaxy assembly and evolution are mainly driven by hierarchical mergers, meaning that the halos of larger galaxies build up and grow via the accretion of smaller systems. Thus, we expect to see galaxies accompanied by a population of smaller satellites, as widely observed. The largest satellite of our Milky Way is the Large Magellanic Cloud (LMC), with a total mass of ~1-2 × 10¹¹ M_⊙ (e.g. Erkal et al. 2019, 2021; Vasiliev et al. 2021; Kacharov et al. 2024; Watkins et al. 2024). Given its high mass, which is only a factor of ten less than that of our Galaxy, the LMC would have brought its own system of dwarf galaxy satellites as it entered the halo of the Milky Way (e.g. Lynden-Bell 1976; D’Onghia & Lake 2008). Its most noticeable companion is the Small Magellanic Cloud (SMC). Both galaxies, the LMC and SMC, form a pair of interacting dwarf galaxies and are currently undergoing the early stages of a minor merger event. Beyond the SMC, a group of ultra-faint dwarf (UFD) galaxies is further expected to be associated with the LMC (Sales et al. 2011).

Using multi-epoch observations with the Hubble Space Telescope (HST), Kallivayalil et al. (2006b,a, 2013) measured precise proper motions (PMs) of the LMC and SMC, constraining, for the first time, their orbits around the Milky Way. The high tangential velocity of the galaxies indicates they are most likely on their first (e.g. Besla et al. 2007; Patel et al. 2017) or second (Vasiliev 2024) passage around our Galaxy. As a consequence, the system of LMC satellites should still be intact without being stripped by tidal forces of the Milky Way. Indeed, within the last several years, a numbr of UFD galaxies have been identified as satellites of the LMC (see e.g. Jethwa et al. 2016; Sales et al. 2017; Kallivayalil et al. 2018; Erkal & Belokurov 2020; Vasiliev 2024).

The ongoing interaction between the Clouds offers us the unique opportunity to witness the cosmological process of hierarchical structure formation close up. It is also natural to assume that the LMC has already experienced accretion events with smaller satellites in its past. Thus, reconstructing the interaction history of the LMC will provide us with valuable insights in the assembly of systems on the dwarf-galaxy scale. In retracing the evolution of our Galaxy, globular clusters (GCs) have proven to be a powerful tool (e.g. Searle & Zinn 1978; Forbes & Bridges 2010; Massari et al. 2019; Kruijssen et al. 2019). Compact stellar systems such as GCs are able to survive the merging event in which the host dwarf galaxy is dissolved by tidal forces. These clusters retain the kinematic signature of the merger event. Their characteristics, such as their chemical composition, reflect the properties of the environment they formed in, making them distinguishable from GCs formed in situ (e.g. Horta et al. 2020; Ceccarelli et al. 2024).

Our understanding of the assembly history of the Milky Way has benefitted from the precise astrometric measurements provided by the Gaia mission (Gaia Collaboration 2016, 2018a, 2021). With PMs making up the missing piece of kinematic information up until now, it has become possible to trace the orbits of GCs within the potential of the Milky Way (e.g. Gaia Collaboration 2018b; Baumgardt et al. 2019; Vasiliev & Baumgardt 2021). In combination with determinations of their chemical abundances and ages, most of the Galactic GC population can now be associated with different progenitors that form the present Milky Way (e.g. Massari et al. 2019; Myeong et al. 2019; Callingham et al. 2022; Malhan et al. 2022). Based on the likeliest properties of the progenitors, resulting from the characteristics of the relic GCs, Kruijssen et al. (2020) were able to construct the first merger tree of our Galaxy. However, chemo-dynamical information on its own is not sufficient to unambiguously disentangle the origin of each Milky Way GC (e.g. Minelli et al. 2021; Carretta & Bragaglia 2022). For this reason, Massari et al. (2023) initiated the Cluster Ages to Reconstruct the Milky Way Assembly (CARMA) project to determine in a homogeneous way the ages of the Galactic GCs with the aim to unambiguously assign each GC to a merging event or an insitu formation, based on the age-metallicity relations defined by the different progenitors.

With the availability of precise PM measurements even at the distance of the Clouds (e.g. Massari et al. 2021; Niederhofer et al. 2022; Niederhofer et al. 2024), we are now at the advent of studying in detail the chemo-dynamics of their star cluster systems, providing us with additional clues about the formation, assembly, and evolution of the galaxies, as well as the origins of their cluster population. Based on line-of-sight velocity measurements, the ancient (>10 Gyr) star clusters in the LMC seem to follow a flattened rotating disc-like structure (e.g. Freeman et al. 1983; Grocholski et al. 2006; Sharma et al. 2010), suggesting that the LMC does not harbour a population of halo GCs. Such a lack would be in strong contrast to the system of Milky Way GCs, which is composed of disc and halo clusters with different characteristics, raising the question whether the LMC star-cluster system assembled in a different way than the one of our Galaxy. Bekki (2007) discussed the implications of these observations and argued that given the halo of the LMC formed from accreted low-mass systems, akin the build-up of the Milky Way halo, the ancient star clusters and halo field stars should share the same kinematics. As a possible formation scenario, they proposed that there might exist a lower mass limit for a low-mass galaxy to form GCs. Thus, the LMC halo might have formed solely from systems below this threshold and the ancient star clusters arose in the early stages of the formation of the disc. Wagner-Kaiser et al. (2017) estimated the ages of a sample of six ancient LMC clusters, by determining the brightness of their main sequence turn-off points. Their results suggested that their sample of clusters being coeval to Galactic GCs, indicating they both formed at similar times.

Piatti et al. (2019) were the first to study the full 3D velocity structure of the old LMC star-cluster system by combining line-of-sight velocities with PM measurements from the second data release (DR2) from Gaia (Gaia Collaboration 2018b). In contrast to the findings based on 1D kinematics, the authors of the cited study claimed to have found evidence for two distinctive populations of clusters, one of which showing halo-like kinematics and spatial distributions. This result was later questioned by Bennet et al. (2022), who combined Gaia and HST data to measure the PMs of a sample of 32 young (<3 Gyr) and old (>10 Gyr) LMC star clusters. Including literature line-of-sight velocities and distances, the resulting dynamic structure they found does not show any indication of a bi-modality; furthermore, it is consistent with a single population with disc-like kinematics.

To investigate the origin of the LMC star-cluster system from a chemical point-of-view, Mucciarelli et al. (2021) carried out a detailed analysis of high-resolution spectra for a sample of 11 old LMC clusters. They discovered that NGC 2005 displays a peculiar chemical abundance pattern and argued this cluster likely formed in a low-mass galaxy that had been accreted by the LMC in a past merger event.

The inconclusive results regarding the nature of the old LMC star-cluster system illustrate the need for a homogeneous derivation and analysis of its chemo-dynamical properties. In this paper of the CARMA series, we employ our recently published HST astro-photometric catalogues (Niederhofer et al. 2024) to derive the PMs, ages, metallicities and distances of a sample of seven old LMC star clusters with available multi-epoch HST data in a robust and self-consistent way. This work represents the first successful attempt to combine dynamical information with accurate determinations of the cluster properties (age and metallicity) to investigate the formation of the LMC cluster system.

The paper is organised as follows. In Section 2, we describe the applied datasets and provide a brief outline of the photometric and astrometric reduction procedures. In Section 3, we describe the isochrone fitting routine and present the results for our sample of clusters. We derive the kinematics of the clusters within the frame of the LMC in Section 4 and discuss the results and its implications regarding the different clusters. In Section 5, we provide a summary of the paper and draw our final conclusions.

2 Data sets and reduction

For the seven clusters analysed in this study, we used multi-epoch HST observations taken with the Ultraviolet-Visible (UVIS) channel of the Wide Field Camera 3 (WFC3) and the Wide-Field Channel (WFC) of the Advanced Camera for Surveys (ACS). The astro-photometric catalogues of six of those clusters were presented in Niederhofer et al. (2024)¹. A detailed list of all observations used to create these catalogues can be found in the appendix of Niederhofer et al. (2024) and in the Mikulski Archive for Space Telescopes (MAST) under the following DOI: 10.17909/7d5e-s940. In addition to these six clusters, we also included the old LMC cluster Reticulum, which we found to have archival long time-baseline data, suitable for PM determinations. The list of observations for Reticulum is given in Table 1.

The photometric reduction is described in Niederhofer et al. (2024) and follows the prescriptions and methods presented in Bellini et al. (2017b, 2018). In brief, we performed our analysis on the un-resampled_flc images that have been corrected for the effects of imperfect charge transfer efficiency (CTE; see Anderson & Ryon 2018; Anderson et al. 2021). The reduction routine is a combination of a first- and second-pass photometric run, using hst1pass (Anderson 2022) and KS2 (see Sabbi et al. 2016; Bellini et al. 2017b, for details). To fit the PSFs of the stars, we used the focus-diverse PSF models², which we perturbed to adjust them to each exposure. We corrected the measured positions of the stars for geometric distortions, applying the solutions by Anderson & Bedin (2010) for ACS/WFC³ and Bellini & Bedin (2009) and Bellini et al. (2011) for WFC3/UVIS. Finally, we calibrated the instrumental magnitudes to the VEGA-mag system, as described in Bellini et al. (2017b).

We measured relative PMs in an iterative way following the techniques developed by Bellini et al. (2014, 2018); Libralato et al. (2018, 2022). We refer to these papers for a detailed description of the methods and outline the main steps below. We transformed the positions of the stars in the individual exposures to a common frame of reference by means of a general six-parameter transformation, using well-measured cluster members as reference objects. Thus, the PMs are relative to the bulk motions of the clusters. Within each iteration, this list of reference stars is refined, excluding sources that are not in agreement with the motion of the cluster (i.e. those that do not cluster around 0). The transformed positions of the stars as a function of observing epoch are then fitted with a straight line. The slopes of the fitted lines directly correspond to the PMs of the stars. We inspected the PMs and corrected them for spatially variable systematic effects, as described in Bellini et al. (2014). Finally, we transformed the relative PMs to an absolute scale using stars from the Gaia DR3 catalogue (Gaia Collaboration et al. 2023). Table 2 presents the absolute PMs of the star clusters, together with their RA and Dec positions and literature line-of-sight velocity measurements.

3 Isochrone fitting

In order to obtain homogeneous, self-consistent ages, metallic-ities, and distances of our sample of clusters, we employed the isochrone-fitting routine used in the CARMA project. This code was originally developed by Saracino et al. (2019) and updated and refined by Massari et al. (2023). In the following, we provide a brief overview of the employed method and refer to Massari et al. (2023) for a detailed description.

3.1 The isochrone-fitting methodology

The isochrone fitting routine developed within CARMA allows us to simultaneously infer the metallicity [M/H], reddening E(B – V), distance modulus (m – M), and age of the clusters. It is based on a Bayesian statistical framework and implements a Markov chain Monte Carlo (MCMC) sampling technique to explore the posterior probability distribution. The code performs a star-by-star comparison between an observed colourmagnitude diagram (CMD) in a given filter combination and a set of theoretical isochrone models, determining the associated log-likelihood. The log-likelihood function is the sum of two terms and is given by $$\mathrm{l}\mathrm{n}{\mathcal{L}}_{\mathrm{t}\mathrm{o}\mathrm{t}}=\mathrm{l}\mathrm{n}{\mathcal{L}}_{\mathrm{p}\mathrm{r}\mathrm{i}\mathrm{o}\mathrm{r}}+0.3\times \mathrm{l}\mathrm{n}{\mathcal{L}}_{\mathrm{f}\mathrm{i}\mathrm{t}}.$$(1)

The first term is the likelihood associated to the initial priors that we provide for the parameters [M/H], E(B – V), and (m – M). For each of the three parameters, we opted for a Gaussian form; thus, lnL_prior can be written as $$\mathrm{l}\mathrm{n}{\mathcal{L}}_{\mathrm{p}\mathrm{r}\mathrm{i}\mathrm{o}\mathrm{r}}=-0.5\times \sum _{i=1}^3\frac{(x_i-x_{i,\mathrm{p}\mathrm{r}\mathrm{i}\mathrm{o}\mathrm{r}}{)}^2}{x_{i,\mathrm{s}\mathrm{t}\mathrm{d}}^2},$$(2)

where x corresponds to [M/H], E(B – V), and (m – M), while x_prior and x_std are the provided priors and their associated uncertainties, respectively. The second term in Equation (1) (lnL_fit) gives the probability of the observed CMD for a given model. It can be written as $$\mathrm{l}\mathrm{n}{\mathcal{L}}_{\mathrm{f}\mathrm{i}\mathrm{t}}=-0.5\times \sum _{i=1}^N\frac{[\mathrm{m}\mathrm{i}\mathrm{n}(dis{t}_i){]}^2}{{\sigma }_i^2},$$(3)

where N is the total number of stars, min(dist_i) gives the minimum distance between star i and the model, and σ is the photometric uncertainty. In Equation (1), a factor of 0.3 is introduced in front of ln L_fit to prevent the code from getting stuck in local minima (see Massari et al. 2023).

To infer the best-fitting values for [M/H], E(B – V), (m -M), and age, along with its associated uncertainties, the code explores the 4D posterior probability space using the MCMC sampler emcee (Foreman-Mackey et al. 2013), which is a python implementation of the affine-invariant MCMC ensemble sampler (Goodman & Weare 2010).

As theoretical stellar models to which the observed CMDs are compared, the routine adopts the Bag of Stellar Tracks and Isochrones (BaSTI) database (Hidalgo et al. 2018; Pietrinferni et al. 2021). We used the newest release that also includes the effects of diffusive processes. As discussed in Massari et al. (2023) we made use of solar scaled models and fitted for the global metallicity [M/H], rather than the iron content [Fe/H]⁴. Since measurements of the α-element abundance do not exist for all clusters in our sample, this choice allows us to achieve homogeneous results that are not prone to any assumptions on the [α/Fe] content of the clusters.

To create the models that are fitted to the observed CMDs, we built a grid of theoretical BaSTI isochrones with ages ranging from 6 Gyr to 15 Gyr in steps of 100 Myr, along with metallic-ities going from [M/H] = 0.0 dex to [M/H] = −2.0 dex with a step of 0.01 dex. We further interpolated the models in age for a finer sampling. To apply interstellar reddening to the models, we determined the extinction coefficients, A_λ, in the different ACS/WFC and WFC3/UVIS filter bands, assuming a Cardelli et al. (1989) extinction law. To take into account the dependence of A_λ on the stellar spectral type, we created for each filter a grid of A_λ values over a large range of stellar effective temperatures and for E(B – V) values ranging between 0.0 mag and 1.0 mag, in steps of 0.1 mag. In the next step, we linearly interpolated this grid.

3.2 Selection criteria

To ensure the isochrone-fitting routine was run only on sources with the best measurements, we applied several photometric and astrometric selection criteria to our catalogues. To achieve this, we closely followed the prescriptions in Niederhofer et al. (2024) and applied the same quality cuts. Specifically, we used the following photometric diagnostic parameters, which are provided by the KS2 output: the quality-of-fit or QFIT parameter (which indicates how well a source was fitted by the PSF model), photometric rms error, the shape parameter RADXS (which indicates how extended a source is with respect to the PSF model, see Bedin et al. 2008), isolation parameter of a source (which gives the fraction of flux within the PSF fitting aperture that comes from neighbouring sources, before neighbour subtraction), and the fraction of good measurements of a source with respect to the total number of detections. The astrometric quality parameters are based on the PM measurements and include the reduced χ² of the PM fit, fraction of data points of a source actually used to determine its PM, and uncertainty in the PM measurement.

We further took advantage of the available kinematic information to minimise the contribution from field stars within our sample of well-measured stars. Although there is in most cases a considerable overlap between the PM distribution of cluster and field stars, we were able to remove a large fraction of the field stars based on their motions, given the larger velocity dispersion of the field stars. This selection is illustrated in Fig. 1 for the example of NGC 2210. Panel a shows the CMD of all well-measured stars in the field of NGC 2210. We selected (by hand) stars that follow the bulk of cluster stars in the m_F814W vs PM diagram (blue line in panel b). We note that here we exploit the relative PMs; thus, cluster stars have motions that are close to zero, whereas field stars show larger motions. Panel c shows the CMD of the PM-selected cluster stars. We further refined our selection by restricting it to stars that are within one effective radius of the clusters (as determined by Niederhofer et al. 2024). This final selection of cluster members is displayed in panel d of Fig. 1.

3.3 Differential reddening correction

Before running the isochrone fitting routine, we corrected the photometry of the seven clusters for the effects of differential reddening across the cluster field, applying the widely used techniques as described by Milone et al. (2012) and Bellini et al. (2017a). Briefly, we first selected for each cluster a sample of reference stars along the main sequence, excluding stars on the binary sequence. We then determined the fiducial line for these stars. For this, we divided the magnitude range covered by the reference stars into bins of 0.2 mag and calculated the sigmaclipped median colour within each bin. These median colours were then fit with a cubic spline model. For each reference star we subsequently determined its offset from the fiducial line in the direction of the reddening vector. As extinction coefficients in the different filter bands, we assumed the A_λ values, as explained in Section 3.1. Finally, we determined the reddening of each star as the 2.5σ-clipped median value of the closest 75 reference stars, excluding the target star itself. We note that this procedure only corrects for any variable reddening across the cluster; thus, the photometry is still affected by the mean global E(B – V). We found moderate reddening variations across the fields of Hodge 11, NGC 1841 (−0.02 mag ≲ δE(B -V)≲0.02 mag for both clusters), and NGC 1898 (−0.04 mag ≲ δE(B – V) ≲0.04 mag). We will use the corrected photometry of these clusters for the analysis. The remaining four clusters show negligible effects of differential reddening.

Fig. 1 Illustration of the selection of cluster member stars for NGC 2210. (a) The mF814W vs mF606W – mF814W CMD of well-measured stars in the field of the cluster for which PMs have been determined (b) 1D relative PMs as a function of the mF814W magnitude. The blue line (drawn by hand) follows our selection of likely cluster members (black dots) based on their motions. (c) CMD of PM-selected cluster members. (d) CMD of PM-selected cluster stars, including only stars within one effective radius from the cluster centre.

3.4 Results

For each cluster, we ran the code on two different CMDs, resulting from the combination of the available filters for each cluster. For Hodge 11, NGC 1841, NGC 2210, NGC 2257, and Reticulum, these are the (m_F814W vs m_F606W – m_F814W ) and (m_F606W vs m_F606W – m_F814W) CMDs. For NGC 1898, the (mm_F438W vs ^mm_F438W – m_F814W) and (^m_F814W vs m_F438W – m_F814W) CMDs and for NGC 2005 the (m_F435W vs m_F435W – m_F814W) and (m_F814W vs m_F435W – m_F814W) CMDs. We did not use observations in the F336W filter for the fitting, since the bolometric corrections show an increasing sensitivity on the exact chemical abundance pattern of the stars in the ultraviolet wavelengths (e.g. Cassisi et al. 2004; Pietrinferni et al. 2024); thus, they might have an impact on the final results when using only solar-scaled models. Although there are more Reticulum observations in more filter bands available, we opted to only use data in the ACS/WFC F606W and F814W filters to maintain consistency with the other clusters in our sample. Before fitting the isochrones, we had to tailor the input catalogues for the code to give reliable results. Specifically, we selected stars within each CMD that follow the cluster sequence, namely, excluding stars on the binary sequence, blue straggler stars, and left-over field interlopers. For this purpose, we first defined a magnitude bin of 0.1 mag, starting at the brightest magnitude in the catalogue. If this bin contained more than five stars, we determined the sigma-clipped median colour of the stars within this bin and kept all stars that are within 1.5-2.0 times (depending on the cluster) the median colour uncertainty in the bin. We then moved this bin by 0.001 mag to fainter magnitudes and repeated the process until we reached the faintest magnitude. Finally, we also defined a bright and faint limit to the list of stars. The main reason for this is to prevent giving stars on the main sequence too much weight in the fitting. Tests have shown that including stars that are within +2.5 mag and −3.5 mag from the main sequence turn-off yields the most robust results in terms of age and metallicity.

For each cluster, we set up Gaussian priors for the metallicity, distance modulus and reddening. For [M/H] we chose the spectroscopically determined mean values and uncertainties as given by Grocholski et al. (2006) and Mucciarelli et al. (2021), whereas for the distance modulus and E(B – V) values, we used the measurements provided by Wagner-Kaiser et al. (2017) and Milone et al. (2023). The spectroscopic metallicity measurements will provide the most robust prior information for the fitting and also break the age-metallicity degeneracy. The isochrone fitting routine does not require any priors on the ages to freely explore the parameter space. We created an ensemble of 100 walkers and ran the MCMC for a total of 1 000 steps, using the first 200 steps as burn-in phase. The inferred bestfit solutions and their associated uncertainties correspond to the 50th, 16th, and the 84th percentiles of the posterior distributions, respectively. In Appendix A, we provide, for each of the seven clusters, the specific CMD with the best-fitting isochrone model, along with the corresponding corner plots showing the posterior probability distribution.

The results of fitting process are reported in Table 3. The values provided here are the average values resulting from the individual fits to the two CMDs. The uncertainties associated with each value are derived accounting for the upper and lower limits of the inferred quantities from both fits. These are the values we adopted in our further analysis. We want to highlight that the values we present in this study (and the associated uncertainties) are to be intended in a relative sense, rather than in an absolute one.

As a test of consistency of our results, we compared our inferred values for age, metallicity, distance and reddening to values from the literature. The top-left panel in Fig. 2 shows the difference in cluster ages between our values and the ones from Milone et al. (2023), which are derived from isochrone fitting, as well. We note that there is a mean offset of ~0.6 Gyr, that most likely originates from the use of different theoretical isochrone models. While the BaSTI set of isochrones is used in this work, Milone et al. (2023) employed the database of Padova models (Marigo et al. 2017). However, since we are only interested in relative ages of the clusters, and not absolute ages, that offset is not an issue. More importantly, we found a good concordance in the relative ages, with a spread around the mean difference of only 0.2 Gyr.

In the top-right panel of Fig. 2, we present the comparison of our inferred metallicities with literature spectroscopic measurements from Grocholski et al. (2006), Song et al. (2021), and Mucciarelli et al. (2021), as well as with metallicity estimates based on periods and amplitudes of fundamental mode RR Lyrae stars (Sarajedini 2024). Since we determined the global metallicity [M/H] so that we could perform a proper comparison, we needed first to convert the measured literature [Fe/H] values to [M/H], taking into account the respective α-element abundances of the clusters. Salaris et al. (1993) showed that the global metallicity can be derived with the following relation, $$[\mathrm{M}/\mathrm{H}]=[\mathrm{F}\mathrm{e}/\mathrm{H}]+\mathrm{l}\mathrm{o}\mathrm{g}(0.694\times {10}^{[\alpha /\mathrm{F}\mathrm{e}]}+0.301).$$(4)

We note that here the coefficients have been slightly modified with respect to the original solution by Salaris et al. (1993), taking into account the reference solar mixture used in BaSTI models (see also Massari et al. 2023). We converted the [Fe/H] measurements from Mucciarelli et al. (2021) to [M/H], using their determined α-element abundances. Since other literature studies only measured the iron content for their samples of clusters, we determined an average α-abundance for the clusters in common by least-squares fitting Equation (4) for [α/Fe] using the [Fe/H] values from the different literature studies and the [M/H] inferred by us. We found mean α-abundances for the literature values of [α/Fe] = 0.32±0.05 dex. The shown error bars in the top-right panel of Fig. 2 are the combined uncertainties from our measurements, the literature values, as well as the uncertainties from the α-abundance estimates. The comparison with the literature [M/H] measurements shows a good agreement for most of the clusters, with an average difference between our estimates and the literature of −0.01 dex, with σ=0.09. For Hodge 11, NGC 1841, NGC 2005, NGC 2257, and Reticulum, the different measurements from the literature agree well within the uncertainties with our inferred values. For NGC 2210, our value is slightly higher than the measurements by Mucciarelli et al. (2021) and Sarajedini (2024). It is only the metallicity we inferred for NGC 1898 that shows a discrepancy of more than 2σ with respect to the value determined by Mucciarelli et al. (2021); however, it agrees within the 1σ uncertainties with the metallicity derived by Sarajedini (2024).

The bottom-left panel of Fig. 2 presents the comparison of our results for the distance modulus with the findings from Wagner-Kaiser et al. (2017) and Milone et al. (2023). In addition, Wagner-Kaiser et al. (2017) determined the distances to a sample of six old LMC clusters using two different methods. The first method is based on comparing the horizontal branch magnitudes to the ones of a sample of Galactic GCs with known distances. The second method involves fitting the colours and magnitudes of several reference subdwarfs to the main sequence fiducial lines of the LMC clusters. The comparison reveals an excellent agreement with the results from Milone et al. (2023) and the subdwarf-based distances from Wagner-Kaiser et al. (2017); whereas there are significant discrepancies from the values derived using the horizontal-branch brightness. These discrepancies are likely due to the number of assumptions involved when estimating the distances based on the horizontal branch brightness. As an additional consistency check on the distances and the sensitivity of the fitting on the assumed priors, we used literature RR Lyrae-based distances (from Walker 1989, 1990, 1992a,b; Kuehn et al. 2013; Cusano et al. 2021) for several of the clusters in our sample as new priors and performed the fitting again. The differences in distance modulus, compared to the original results, are all less than 0.01 mag (<0.2 kpc, at the distance of the LMC), suggesting the fitting is not very sensitive to the exact choice of the distance prior.

Finally, the bottom-right panel of Fig. 2 shows the comparison of our derived value for the interstellar reddening with measurements from Wagner-Kaiser et al. (2017) and Milone et al. (2023), as well as with values resulting from the extinction maps from Schlegel et al. (1998) and Skowron et al. (2021). Overall, we find a good agreement with the derived values for the clusters from the literature studies. Only for NGC 1841, Milone et al. (2023) derived a considerably smaller E(B – V) (about 0.05 mag lower) compared to our measurement. Furthermore, our reddenings match the values from the two extinction maps very well. The large error bars for NGC 1898 and NGC 2005 from the comparison with the Skowron et al. (2021) reddening map result from the high uncertainties within the central regions of the LMC that are present in their map.

The above consistency checks have demonstrated that our isochrone-fitting method is able to produce reliable results for old star clusters in the LMC. With these results in hand, we can first focus on the ages and metallicities of the clusters in our sample. Figure 3 shows the positions of the seven LMC clusters (blue symbols) in the age-metallicity space, compared to that of the nine Milky Way GCs selected as genuine members of the Gaia-Sausage-Enceladus (GSE, Helmi et al. 2018; Belokurov et al. 2018) merger event by Aguado-Agelet et al. (2025, shown as red symbols). We chose GSE for the comparison, since its estimated mass (~6 × 10¹⁰ M_⊙, Helmi et al. 2018; ~9 × 10¹⁰ M_⊙, Kruijssen et al. 2020 ) is comparable to the one of the LMC (1-2 × 10¹¹ M_⊙ Vasiliev 2023). The first feature that is immediately evident is that our sample of LMC clusters can be divided into two age groups. NGC 2005 and Hodge 11 with ages of 14.0 Gyr and 13.8 Gyr, respectively, are the two oldest clusters within our sample. The remaining five clusters (NGC 1841, NGC 1898, NGC 2210, NGC 2257, and Reticulum) are significantly younger (and have similar ages). Their mean age is 12.4 Gyr with a small dispersion of σ_t = 0.2 Gyr. GSE clusters describe two epochs of formation as well, but at clearly different times. The second remarkable feature is that the general trend of the age-metallicity relation (AMR) of the LMC clusters closely follows that shown by the GSE GCs. This is also evident when comparing the location of the LMC clusters with the parametric AMR model derived by Massari et al. (2019) for the GSE GCs (red solid line). The shape of this model ultimately depends on the mass of the progenitor and on its star-formation rate, assumed to be constant (Prantzos 2008). The similarity of the two observed AMRs points towards a similar early evolution of the two dwarf galaxies, as possibly indicated by the chemistry of their metal-poor populations as well (e.g. Matsuno et al. 2021). One small difference is evident based on the location of the two oldest LMC clusters, which appear to reside systematically above the theoretical model, even when considering the associated uncertainty (red dashed lines). Given the shape of the AMR (as described by Prantzos 2008), this might indicate that the mass (or the star-formation efficiency) of the LMC at the time of the formation of these clusters was slightly higher than that of GSE. We note here that, despite the similar [M/H] values, Hodge 11 and NGC 2005 are characterised by substantially different chemical compositions (Mucciarelli et al. 2021). NGC 2005 shows chemical abundance ratios in several elements that are in strong contrast to the ones of the other ancient LMC clusters. Thus, Hodge 11 and NGC 2005 formed from distinct gas clouds with different enrichment histories. We discuss the origin of NGC 2005 in more detail in Sect. 4.3. For the sake of comparison, in Fig. 3, we also show the AMR model used by Massari et al. (2019) to describe in situ MW GCs: the much larger early mass and star-formation rate of the MW manifest in a significantly different global trend (see the black solid line) compared to the other two dwarf galaxies (see also e.g. Marín-Franch et al. 2009; Leaman et al. 2013; Souza et al. 2024). The final striking feature is the unique position of NGC 1841 in the age-metallicity space. While the other clusters broadly follow the theoretical AMR model (suggesting they formed within the LMC), NGC 1841 has a significantly lower metallicity than the other clusters in the younger age group. Instead, its [M/H] value is consistent with that of NGC 2005, the oldest cluster in our sample. This anomalous position of NGC 1841 with respect to the other studied clusters prompts us to the speculation that NGC 1841 might not have formed within the LMC; rather, it originates from an environment with a lower star formation efficiency and, thus, a slower chemical enrichment rate.

As an independent test to verify the relative ages among the LMC star clusters in our sample, we employed the horizontal method (e.g. Vandenberg et al. 1990; Salaris & Weiss 1998), which is based on measuring the distance in a given colour between the main sequence turn-off point and the base of the red giant branch. This difference in colour for a given metallicity is dependent on the age of a cluster, with older clusters having smaller distances. Since there is a strong metallicity dependence in the horizontal method, we only performed direct comparisons between Hodge 11 and NGC 1841, as well as between NGC 2210, NGC 2257, and Reticulum. For NGC 2005 and NGC 1898, a direct comparison with the other clusters is not possible, since they have observations in different filters (F435W and F814W for NGC 2005, and F438W and F814W for NGC 1898). The tests verified that Hodge 11 is indeed older than NGC 1841 and, in addition, that NGC 2210, NGC 2257, and Reticulum are approximately coeval. In Appendix B, we illustrate the CMDs of Hodge 11, NGC 2005, and NGC 1841 together with isochrones of different age and metallicities.

Fig. 2 Comparison between the cluster parameters derived in this work and literature values. Top-left panel: differences between the measured cluster ages and ages derived by Milone et al. (2023). Top-right panel: differences between the inferred metallicity and the results from Grocholski et al. (2006), Mucciarelli et al. (2021), Song et al. (2021) and Sarajedini (2024). Bottom-left panel: differences between the distance modulus derived in this study and the measurements from Wagner-Kaiser et al. (2017) and Milone et al. (2023). Bottom-right panel: differences between the derived E(B – V) values and the reddenings from Wagner-Kaiser et al. (2017) and Milone et al. (2023), as well as the resulting values from the reddening maps from Schlegel et al. (1998) and Skowron et al. (2021).

Fig. 3 Resulting age-metallicity relation of the seven old LMC clusters studied in this work (blue symbols). For comparison, nine Milky Way GCs determined as genuine members of GSE (red symbols) are also shown. The red and black solid lines follow the parametric age-metallicity model derived by Massari et al. (2019) for GSE and the MW, respectively. The red dashed lines illustrate the associated uncertainties of the GSE model.

4 Kinematic structure of the old LMC clusters

To build on our hypothesis regarding the origin of NGC 1841 and to reach an unequivocal conclusion, we took into account the kinematics of the clusters within the LMC, as well. Combining the observed positions on the sky and 3D velocities from Table 2 with the determined distances from the isochrone fitting (Table 3), we have the full 6D phase-space information for the clusters, allowing us to examine the velocities and positions of the clusters within the LMC in detail, along with an estimation of their orbits within the potential of the galaxy.

4.1 Velocity and coordinate transformation

Before we can analyse the kinematics of the seven old LMC clusters in our sample and determine their orbits within the galaxy, we first need to transform their observed on-sky positions, distances, line-of-sight velocities, and PMs into the reference frame of the LMC. This frame is defined as a right-handed orthogonal coordinate system, centred on the dynamical centre of the galaxy. The X-Y plane is aligned with the disc of the LMC, where the X axis is along the line of nodes (the line where the LMC-disc plane intersects the sky plane). To project the observed motions and positions into the coordinate system of the LMC, we applied the formalism derived by van der Marel & Cioni (2001) and van der Marel et al. (2002). To solve the transformation equations, we adopted the following values for the position, orientation, and motion of the LMC disc: the coordinates of the LMC dynamical centre (α₀, δ₀) = (79.95°, −69.31°), the PM of the LMC centre-of-mass (μ_αcos(δ)₀, μ_δ,0) = (1.867, 0.314) mas yr⁻¹, the inclination angle i = 33.5°, and the angle of the line-of-nodes Θ = 129.8° (all taken from Niederhofer et al. 2022). Additionally, we adopted a distance to the LMC of 49.59 kpc (Pietrzynski et al. 2019) and a line-of-sight velocity of 262.2 km s⁻¹ (van der Marel et al. 2002). The transformed positions and velocities in the LMC frame, expressed in cylindrical coordinates, are presented in Table 4. The associated errors are calculated by propagating the measurement uncertainties throughout the transformation equations. We like to stress here that the exact values of the reported positions and velocities are specific to the adopted parameters of the LMC. In the literature there exist a broad variety of measurements of these parameters, employing different tracers and methods (see, e.g. Niederhofer et al. 2022, for an overview). However, any other choice of the adopted orientation, position and motion of the LMC would not have significantly affected our final results. We verified this by performing the dynamic analysis described below again assuming the parameters of the LMC as determined by van der Marel & Kallivayalil (2014).

4.2 Kinematics within the LMC

In Fig. 4, panels a and b present, for each cluster, the tangential (V_φ, panel a) and radial (V_R, panel b) velocity components, respectively, as a function of the cylindrical galactocentric radius, R. Also shown in panel a are two rotation curves of the LMC as determined by Kacharov et al. (2024), based on dynamical models. For the first model (orange line), the authors fitted axisymmetric Jeans dynamical models to stars in the Gaia DR3 catalogue with measured PM and line-of-sight velocities. In this model, it is assumed that the mass distribution of the LMC follows a pure dark-matter halo with a spherical Navarro, Frenk, & White (1997, NFW) mass profile. Their best-fit model gives a virial mass of the LMC of M₂₀₀ = 1.81 × 10¹¹ M_⊙ at a virial radius r₂₀₀ of 75 kpc, with an overdensity ∆_c = 200 relative to the critical density of the Universe. The second model (green dashed line) is based on a Schwarzschild orbit superposition method and is fitted to the line-of-sight velocity field of the LMC. In this model, Kacharov et al. (2024) also included the contribution from the stellar component of the galaxy, which is described as a triaxial bar and an axisymmetric disc. The shown model assumes a mass-to-light ratio M/L of 1.5 M_⊙/L_⊙ and a total luminosity of the LMC of 1.3 × 10⁹ L_⊙.

We can see from panels a and b of Fig. 4 that most clusters have tangential velocities that are smaller than the predicted model rotation curves and also have non-zero radial velocities, suggesting the clusters are on non-circular, elongated orbits. A particular case depicts NGC 1841, highlighted with a triangle symbol in all panels of Fig. 4. It has a small negative tangential velocity, which is almost consistent with being zero (V_φ = −7.5 ± 6.1 km s⁻¹) and a relatively large radial velocity (V_R = 41.9 ± 6.7 kms⁻¹), which leads us to the conclusion that NGC 1841 is on a highly eccentric, retrograde orbit.

Panel c of Fig. 4 shows the vertical, out-of-plane velocity component (V_z) of the clusters as a function of R. Studying the PMs of 15 old LMC star clusters using data from the Gaia DR2, Piatti et al. (2019) claimed the existence of two kinematically distinct populations of clusters, a disc, and a halo population. Their separation criterion was solely based on a sharp limit in |Vz|. According to their classification, NGC 1841, NGC 2005, and NGC 2210 would belong to the disc family of clusters, whereas Hodge 11, NGC 1898, NGC 2257, and Reticulum would be halo clusters (see Table 1 in Piatti et al. 2019). However, when we discuss the orbits of the clusters (see Section 4.3), it is clear that our data do not support this categorisation. Our measurements indicate that most clusters in our sample have absolute vertical velocities smaller than 40 kms⁻¹. The two remarkable exceptions with extreme out-of-plane velocities are NGC 1841 with V_z = 62.7 ± 5.4 kms⁻¹ and NGC 2210 with V_z = −85.7 ± 4.0 kms⁻¹. From their study of the kinematics of LMC star clusters, Bennet et al. (2022) determined an unusually high vertical velocity of NGC 2210 (−111 ± 8 kms⁻¹), as well. Their measured total velocity of NGC 2210 (~160 km s⁻¹ ) is close to the local escape velocity of the LMC (~166 kms⁻¹, as determined by Boubert et al. 2017, based on LMC runaway stars); thus they speculated that NGC 2210 might not physically be associated with the LMC and that, rather, it had originated from the Milky Way halo. Our measurements give a smaller total velocity of NGC 2210 of 143.6 ± 6.8 kms⁻¹ (mainly due to a closer distance to NGC 2210 and a smaller assumed line-of-sight velocity, compared to that employed by Bennet et al. 2022). Its position in the age-metallicity space further coincides with the other LMC clusters of similar ages. Thus, we conclude that NGC 2210 is likely to be a genuine member of the LMC.

Finally, panel d of Fig. 4 illustrates the vertical distance of the seven clusters from the disc plane of the LMC, as a function of the cylindrical radius R. Also in this plot, NGC 1841 stands out as an outlier with respect to the other clusters in the sample. While the current positions of most clusters are close to the LMC disc plane (|Z| <1.2 kpc), NGC 1841 is currently located more than 9 kpc above the plane of the galaxy, suggesting NGC 1841 is on a highly inclined orbit around the LMC.

Fig. 4 Outline of the clusters’ velocities and positions within the LMC. In each panel, NGC 1841 is highlighted with triangle symbol. Panel (a) shows the tangential velocity, Vφ, as a function of the radial distance from the LMC centre R. Model circular velocities resulting from a pure NFW profile (orange solid line) and a model composed of a dark matter halo and a stellar disc and bar (green dashed line) are also shown. The size of the stellar bar is indicated by the vertical dashed line. The radial velocity VR as a function of R is presented in panel (b). Panel (c) illustrates the out-of-plane velocity VZ as a function of R and panel (d) shows the vertical distance from the plane Z as a function of R.

4.3 Cluster orbits

The sum of the evidence collected above strongly suggests that NGC 1841 differs from the other studied LMC clusters also from a dynamical point of view. The measured velocities and positions of the clusters, however, only form a snapshot of the momentary state of the cluster system. To investigate the kinematic structure of the cluster system in a more descriptive way, we can model and examine the orbits of the clusters within the LMC. To trace the orbits, we used two of the LMC potentials derived by Kacharov et al. (2024) based on Jeans and Schwarzschild models, with the following considerations in mind. For clusters within the inner 2 kpc of the LMC (the region where the dynamics are dominated by the bar), we employed the multicomponent model that is composed of a spherical NFW dark matter halo, as well as an axisymmetric stellar disc and a triaxial stellar bar component. As pointed out by Kacharov et al. (2024), the pure NFW model fails to properly reproduce the stellar kinematics within the central parts of the LMC. We note here that their best-fit Schwarzschild model results in an un-physically low mass-to-light ratio of M/L = 0.3 M_⊙/L_⊙, and a virial mass of 4 × 10¹² M_⊙; this is more than a factor of 20 greater than what is now the widely accepted mass estimate of the LMC (1-2 × 10¹¹ M_⊙, see Vasiliev 2023, for a recent review). Thus, we opt to use their model with the mass-to-light ratio fixed at M/L = 1.5 M_⊙/L_⊙, which is a more realistic value for LMC-like galaxies (van der Marel et al. 2002). For radii larger than 2 kpc, we used the pure NFW mass profile, resulting from the Jeans modelling. Even though within about 2 kpc, both the pure NFW model and the composite models with different assumptions about the mass-to-light ratios produce similar results (see figure 17 in Kacharov et al. 2024), they start to deviate beyond 2 kpc. Thus, the pure NFW model depicts the best representation of the outer parts of the LMC. This is also backed up by independent measurements of the LMC potential. Watkins et al. (2024) used the sample of LMC clusters from Bennet et al. (2022) as tracers of the LMC potential, assuming a NFW mass density profile as well. Although the radial limits of the dynamical model fits from Kacharov et al. (2024) are only 6.2 kpc, and the dataset used by Watkins et al. (2024) extends out to about 13 kpc, the profile obtained by Watkins et al. (2024) is identical to the one determined by Kacharov et al. (2024). We want to emphasise here that the potential used here is a simplified (no spiral arms, central bar instead of off-centre bar) and static potential of an isolated (no tidal forces from the Milky Way and SMC) LMC galaxy. Thus, it reflects an idealised case and the interpretations below were made with these caveats in mind.

For the integration of the cluster orbits, we employed the python package galpy⁵ Bovy (2015). For each cluster, we integrated the orbit forward for 1 Gyr. To estimate the uncertainties in the orbital parameters of the clusters, we also created, for each cluster, 500 realisations of the initial conditions, by randomly drawing from Gaussian distributions centred around the measurements for their positions and velocities, and with standard deviations corresponding to the measurement uncertainties.

Figure 5 presents the resulting orbits of the clusters, together with 50 of these random realisations. Each panel shows a different projection of the 3D orbits: in the XY (top panel), XZ (middle panel), and the YZ plane (bottom panel). The current positions of the clusters are also indicated by coloured circles. It is immediately evident from the plots that NGC 1841 (orangecoloured lines) is on a peculiar orbit that is very different from those of the other clusters. With an inclination of 69 ± 6°, its orbital plane is highly inclined with respect to the LMC disc plane. Therefore, the cluster can reach a maximum vertical height of Z_max = 14.85 ± 1.29 kpc above the disc. Additionally, as already predicted from its measured velocity components, NGC 1841 is on a highly eccentric orbit, with an eccentricity of 0.90 ± 0.04. We also note that the Milky Way hosts several GCs with orbital properties similar to those of NGC 1841, where all of those clusters have an accreted origin. The combined properties of NGC 1841 (its unique position in the age-metallicity space and its peculiar kinematics) collected in this work form the consistent picture that NGC 1841 is most likely not a genuine LMC cluster. Since the LMC is massive enough to be accompanied by its own system of dwarf satellite galaxies, including the SMC and other, low-mass systems that have been identified within the last years (see, e.g. Erkal & Belokurov 2020; Patel et al. 2020; Vasiliev 2024), NGC 1841 might have been accreted by the LMC from a smaller companion during a minor merger event. Using E-MOSAICS simulations, Pfeffer et al. (2020) analysed the relation between the observed orbital parameters of accreted clusters and the properties of the progenitor host galaxies. They found that either recent merger events or low-mass progenitors preferentially place clusters at orbits with large apocentres. Thus, NGC 1841 might have been brought in during an early accretion event of a dwarf spheroidal galaxy. The fact that four of the in situ LMC clusters formed at similar times with comparable metallicities suggests that they emerged from a starburst event, which, potentially, could have been linked to that merger. NGC 1841 is about 0.5 Gyr older than these four clusters. Thus, if the merger that brought in NGC 1841 was indeed linked to this starburst event, then NGC 1841 has formed shortly before the merger. Alternatively, NGC 1841 might have been recently stripped from the SMC during a close encounter with the LMC.

The orbits of the other clusters display several interesting features and we discuss them below. Hodge 11, NGC 2257, and Reticulum (blue, brown, and pink lines, respectively, in Fig. 5) orbit the LMC within its disc plane. The orbits of NGC 2257 and Reticulum, with inclination angles of 3 ± 3° and 6 ± 12°, respectively, are consistent, within the errors, of being coplanar with the plane of the disc. Only the orbit of Hodge 11 shows a moderate inclination (18 ± 3°). Thus, it is evident that they follow clear disc-like kinematics. This is in contrast to the results from Piatti et al. (2019), who classified these three clusters as belonging to a halo population.

NGC 1898 and NGC 2005 (green and red lines, respectively, in Fig. 5) are situated close to the centre of the LMC, so their orbits are largely affected by the central bar structure of the galaxy. Both clusters seem to follow box-like orbits (Binney 1982), whereas the orbit of NGC 1898 is more extended (apocentre distance of 2.25 ± 0.38 kpc, with a maximum vertical distance of 1.30 ± 0.34 kpc). NGC 2005 follows a very elongated orbit (eccentricity of 0.89 ± 0.07) along the direction of the bar. Based on the peculiar chemical composition of NGC 2005, Mucciarelli et al. (2021) suggested that it once belonged to a small dwarf galaxy that has since been accreted by the LMC. From the inferred orbit, NGC 2005 appears to be contained within the very central parts of the LMC, with an apocentre distance of 1.18 ± 0.13 kpc and a maximum height above the plane of 0.57 ± 0.20 kpc. Thus, the position and dynamic properties of NGC 2005 measured in this study do not display any obvious indication that the cluster would have been accreted to the LMC (unless the cluster was accreted so long ago that it had the time to sink deep into the potential well of the LMC). Mucciarelli et al. (2021) also noted that the line-of-sight velocity of NGC 2005 is similar to other clusters in its vicinity. One possible explanation might be that the orbit of NGC 2005 was affected by dynamical friction, which caused the cluster to spiral to the central parts of the galaxy as it moved through the galactic disc. On the other hand, Piatti & Hirai (2023) argued that the chemical abundance pattern of NGC 2005 is not a univocal sign of ex situ formation and could instead be consistent with the cluster having formed within the LMC. In this case, NGC 2005 would constitute one of the oldest building blocks of the LMC, given its old age, which had formed near the centre of the potential well of the galaxy and its chemical composition might resemble that of the still inhomogeneous gas at the earliest time of the formation of the LMC. As already noted in Sect. 3.4, the specific abundances in NGC 2005 result in a [M/H] value that is compatible with the AMR of the in situ LMC clusters. Thus, the position of NGC 2005 in the age-metallicity plane is not sufficient for a confirmation or refutation of a potentially ex situ origin. Regardless of the nature of NGC 2005, the old age of the cluster demonstrates that peculiar chemical evolutions had already emerged at very early times.

Moving on from NGC 1841, we can see that NGC 2210 (purple-coloured lines in Fig. 5) follows an extended orbit that is significantly inclined with respect to the LMC disc, although the inclination is not as extreme as in the case of NGC 1841. We found an inclination angle of 37 ± 2° with a maximum height above the plane of Z_max = 6.01 ± 0.48 kpc and an apocentre distance of 11.02 ± 0.98 kpc. The current position of the cluster near its pericentre (4.29 ± 0.14 kpc) results in the rather high measured velocity. Given the kinematic characteristics of NGC 2210, it is tempting to interpret it as a halo cluster. However, as discussed above, it shares the same region in the age-metallicity space as the other clusters with disc-like kinematics and similar ages. Moreover, NGC 2210 is the only cluster in our sample that shows these types of orbital properties. Thus, at the moment, we cannot make any strong claim regarding its nature. An alternative interpretation might be that the orbit of NGC 2210 has been disturbed in the past by a close encounter with the SMC, tilting the orbital plane.

Fig. 5 Orbits of the seven clusters within the LMC, integrated forward for 1 Gyr. The three panels show projections in the XY- (top), XZ- (middle), and YZ-plane (bottom). The current positions of the clusters are indicated by a circle.

5 Summary and conclusions

In this work, we applied the isochrone-fitting algorithm developed within the CARMA project to a sample of seven old LMC star clusters to derive homogeneous and robust estimates for their distances, ages, and metallicities. In combination with their full 3D kinematic data, resulting from multi-epoch HST PM measurements and literature line-of-sight velocities, we determined the positions and velocity components of the clusters in the frame of the LMC galaxy. With this collection of information we explored any connection between the orbital properties of the clusters and their location in the age-metallicity space for the first time in a self-consistent way.

We found that the AMR of the LMC clusters is remarkably similar to that described by Milky Way GCs associated to the GSE merger event. This suggests that during their early evolution, the two galaxies share similar masses and levels of star formation efficiency. Two of the clusters in our sample, Hodge 11 and NGC 2005, are significantly older than the other five clusters, which are approximately coeval. Of this younger group of clusters, NGC 1841 stands out as a clear outlier in the age-metallicity space, owing to its low [M/H] value. The kinematic data further underline the peculiarity of NGC 1841. In contrast to the other clusters, NGC 1841 follows a highly eccentric orbit that is almost perpendicular to the plane of the LMC disc and reaches up to about 15 kpc above the galaxy. Based on the combined properties of NGC 1841, we concluded that this cluster must not have formed within the LMC, but had instead been accreted by the LMC from a smaller dwarf galaxy with a lower star formation efficiency, such as a dwarf spheroidal galaxy. There is also the possibility that NGC 1841 might have been stripped from the SMC during one of the last close encounters between the LMC and SMC.

We analysed the orbital properties of the other clusters as well and found that Hodge 11, NGC 2257, and Reticulum follow disc-like orbits. Our data suggest that NGC 2210 is likely bound to the LMC, but follows a more inclined orbit. Based on its orbit, this system might be a halo cluster, but since its position in the age-metallicity plane is consistent with those of clusters with disc-like kinematics, we speculated that its orbit might have been disturbed. Finally, we did not find any clear kinematic evidence that NGC 2005 has been accreted by the LMC. Instead, our data suggest the cluster, which is the oldest among those of the LMC, is constrained to the inner regions of the galaxy. We interpreted this result with two different scenarios: either NGC 2005 has an ex situ origin and it has been dragged to the centre by, for instance, dynamical friction, or NGC 2005 is a genuine member and formed within the early LMC, where its peculiar chemistry reflects the conditions at the time of formation of the galaxy.

Based on about half of the old LMC cluster population, our results have already shown that precise PM measurements used in combination with homogeneous determinations of the clusters’ distances, ages, and metallicities serve as powerful tools to assess the nature of the LMC cluster population. A natural follow-up of this work would be to expand our sample to the remaining ancient LMC clusters to study the full population. Since most of them have already one existing epoch of HST observations, such an expanded study could be done with minimal observational efforts in the future.

Acknowledgements

We thank the anonymous referee for constructive comments and suggestions that improved the quality of our paper. This research was funded by DLR grant 50 OR 2216. DM, and SC acknowledge financial support from PRIN-MIUR-22: CHRONOS: adjusting the clock(s) to unveil the CHRONOchemo-dynamical Structure of the Galaxy” (PI: S. Cassisi). A.M. acknowledges support from the project “LEGO – Reconstructing the building blocks of the Galaxy by chemical tagging” (P.I. A. Mucciarelli). granted by the Italian MUR through contract PRIN 2022LLP8TK_001. SS acknowledges funding from the European Union under the grant ERC-2022-AdG, “StarDance: the non-canonical evolution of stars in clusters”, Grant Agreement 101093572, PI: E. Pancino. Support for this work was provided by NASA through grants for program GO-16478 from the Space Telescope Science Institute (STScI), which is operated by the Association of Universities for Research in Astronomy (AURA), Inc., under NASA contract NAS5-26555. This work is based on observations made with the NASA/ESA Hubble Space Telescope, obtained from the Data Archive at the Space Telescope Science Institute. This work has made use of data from the European Space Agency (ESA) mission Gaia (https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing and Analysis Consortium (DPAC, https://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. This research made use of astropy⁶, a community-developed core python package for Astronomy (Astropy Collaboration 2013, 2018), iphython (Perez & Granger 2007), Jupyter Notebook (Kluyver et al. 2016), matplotlib (Hunter 2007), numpy (Harris et al. 2020) and scipy (Virtanen et al. 2020).

References

Appendix A Isochrone fitting results

Figs. A.1-A.7, illustrate the results of the isochrone fitting routine for the seven old LMC clusters analysed in this study. For each cluster, the upper two panels show their CMDs, resulting from the two different combinations of the available filters. Overplotted as a red solid line is the best-fitting isochrone model for that respective colour-magnitude space. Stars actually used for the fit are highlighted in green. The lower two panels show corner plots of the posterior distribution of the model parameters and their pairwise correlations.

Fig. A.1 Isochrone fitting results for Hodge 11. (a) Best-fit isochrone model in the mF814W vs mF606W — mF814W CMD. (b) Best-fit isochrone model in the mF606W vs mF606W – mF814W CMD. The red X-symbol shown in both CMDs corresponds to the location of the zero-age horizontal branch. (c): Corner plot of the posterior probability distributions of pairwise model parameters for the mF814W vs mF606W – mF814W CMD. The best-fit parameters are quoted in the labels. (d): Corner plot of the posterior probability distributions of pairwise model parameters for the mF606W vs mF606W – mF814W CMD. The best-fit parameters are quoted in the labels.

Fig. A.2 Same as Fig. A.1, but for NGC 1841.

Fig. A.3 Same as Fig. A.1, but now for NGC 1898 and for the mF438W vs mF438W – mF814W (left panels) and mF814W vs mF438W – mF814W (right panels) colour-magnitude combinations.

Fig. A.4 Same as Fig. A.1, but now for NGC 2005 and for the mF475W vs mF475W – mF475W (left panels) and mF814W vs mF475W – mF814W (right panels) colour-magnitude combinations.

Fig. A.5 Same as Fig. A.1, but now for NGC 2210.

Fig. A.6 Same as Fig. A.1, but now for NGC 2257.

Fig. A.7 Same as Fig. A.1, but now for Reticulum.

Appendix B Isochrone comparisons

In Fig. B.1, we illustrate the observed CMDs of Hodge 11, NGC 2005 and NGC 1841 together with isochrones of different ages and metallicities. In the left and middle panels, we compare the CMDs of Hodge 11 and NGC 2005, respectively, with isochrones of three different ages: one with the best-fit parameters (red curves) and two that have the same parameters as the best-fit one but are younger, 12.3 Gyr (green curve) and 12.8 Gyr (blue curve). For both clusters, the younger isochrones are able to fit the red-giant branch as well as the horizontal branch, but miss the main sequence turn-off, suggesting that Hodge 11 and NGC 2005 are indeed older. In the right panel, we compare for NGC 1841 the best-fitting isochrone (red curve) with an isochrone that resembles the old-metal-poor clusters, i.e. same metallicity as NGC 1841 but with an age similar to Hodge 11 (green line) and with an isochrone that resembles the younger metal-rich clusters in the disc, i.e. same age as NGC 1841 but higher metallicity (blue line). The old, metal-poor isochrone fits the red-giant branch and the horizontal branch of the cluster but misses the turn-off, suggesting a younger age. The metal-rich isochrone misses the turn-off as well as the slope of the red-giant branch, indicating that NGC 1841 is more metal-poor.

Fig. B.1 Observed CMDs of Hodge 11 (left panel), NGC 2005 (middle panel) and NGC 1841 (right panel), together with isochrones of different ages and metallicities (see text).

All Tables

All Figures

Fig. 1 Illustration of the selection of cluster member stars for NGC 2210. (a) The mF814W vs mF606W – mF814W CMD of well-measured stars in the field of the cluster for which PMs have been determined (b) 1D relative PMs as a function of the mF814W magnitude. The blue line (drawn by hand) follows our selection of likely cluster members (black dots) based on their motions. (c) CMD of PM-selected cluster members. (d) CMD of PM-selected cluster stars, including only stars within one effective radius from the cluster centre.

In the text

Fig. 2 Comparison between the cluster parameters derived in this work and literature values. Top-left panel: differences between the measured cluster ages and ages derived by Milone et al. (2023). Top-right panel: differences between the inferred metallicity and the results from Grocholski et al. (2006), Mucciarelli et al. (2021), Song et al. (2021) and Sarajedini (2024). Bottom-left panel: differences between the distance modulus derived in this study and the measurements from Wagner-Kaiser et al. (2017) and Milone et al. (2023). Bottom-right panel: differences between the derived E(B – V) values and the reddenings from Wagner-Kaiser et al. (2017) and Milone et al. (2023), as well as the resulting values from the reddening maps from Schlegel et al. (1998) and Skowron et al. (2021).

In the text

Fig. 3 Resulting age-metallicity relation of the seven old LMC clusters studied in this work (blue symbols). For comparison, nine Milky Way GCs determined as genuine members of GSE (red symbols) are also shown. The red and black solid lines follow the parametric age-metallicity model derived by Massari et al. (2019) for GSE and the MW, respectively. The red dashed lines illustrate the associated uncertainties of the GSE model.

In the text

Fig. 4 Outline of the clusters’ velocities and positions within the LMC. In each panel, NGC 1841 is highlighted with triangle symbol. Panel (a) shows the tangential velocity, Vφ, as a function of the radial distance from the LMC centre R. Model circular velocities resulting from a pure NFW profile (orange solid line) and a model composed of a dark matter halo and a stellar disc and bar (green dashed line) are also shown. The size of the stellar bar is indicated by the vertical dashed line. The radial velocity VR as a function of R is presented in panel (b). Panel (c) illustrates the out-of-plane velocity VZ as a function of R and panel (d) shows the vertical distance from the plane Z as a function of R.

In the text

	Fig. 5 Orbits of the seven clusters within the LMC, integrated forward for 1 Gyr. The three panels show projections in the XY- (top), XZ- (middle), and YZ-plane (bottom). The current positions of the clusters are indicated by a circle.
In the text

Fig. A.1 Isochrone fitting results for Hodge 11. (a) Best-fit isochrone model in the mF814W vs mF606W — mF814W CMD. (b) Best-fit isochrone model in the mF606W vs mF606W – mF814W CMD. The red X-symbol shown in both CMDs corresponds to the location of the zero-age horizontal branch. (c): Corner plot of the posterior probability distributions of pairwise model parameters for the mF814W vs mF606W – mF814W CMD. The best-fit parameters are quoted in the labels. (d): Corner plot of the posterior probability distributions of pairwise model parameters for the mF606W vs mF606W – mF814W CMD. The best-fit parameters are quoted in the labels.

In the text

	Fig. A.2 Same as Fig. A.1, but for NGC 1841.
In the text

	Fig. A.3 Same as Fig. A.1, but now for NGC 1898 and for the mF438W vs mF438W – mF814W (left panels) and mF814W vs mF438W – mF814W (right panels) colour-magnitude combinations.
In the text

	Fig. A.4 Same as Fig. A.1, but now for NGC 2005 and for the mF475W vs mF475W – mF475W (left panels) and mF814W vs mF475W – mF814W (right panels) colour-magnitude combinations.
In the text

	Fig. A.5 Same as Fig. A.1, but now for NGC 2210.
In the text

	Fig. A.6 Same as Fig. A.1, but now for NGC 2257.
In the text

	Fig. A.7 Same as Fig. A.1, but now for Reticulum.
In the text

	Fig. B.1 Observed CMDs of Hodge 11 (left panel), NGC 2005 (middle panel) and NGC 1841 (right panel), together with isochrones of different ages and metallicities (see text).
In the text