© The Authors 2026

1 Introduction

The Milky Way (MW) galaxy is surrounded by numerous self-gravitating stellar systems orbiting within its halo. These satellites are often categorized into two distinct groups: globular clusters (GCs) and dwarf galaxies. The main distinction between these two families is their dark matter (DM) content: dwarf galaxies are considered to be DM-dominated (e.g., Willman & Strader 2012), and GCs are considered to be devoid of DM (e.g., Bradford et al. 2011; Ibata et al. 2013).

Milky Way GCs are old (with ages around 12 Gyr), luminous (visual magnitudes M_V ~ −7), typically populated by tens of thousands to millions of stars, have very high stellar densities (Harris 1996; Bahramian et al. 2013), and are characterized as compact objects with half-light radii of a few tens of parsecs (see Gratton et al. 2019 for a review of GCs). Given the small spread in age and metallicity of their stars, they are often considered to be close approximations of simple stellar populations, although there is a great deal of evidence showing otherwise (see the review by Milone & Marino 2022).

Due to their ages, GCs serve as key witnesses to early galaxy formation and evolution, acting as valuable tracers of the structure of galaxies, including the MW (e.g., Searle & Zinn 1978; Forbes & Bridges 2010; Massari et al. 2019; Kruijssen et al. 2020; Massari et al. 2023). However, there is still no clear consensus regarding their formation mechanism. For example, Peebles (1984) first proposed that GCs could be formed within DM mini-halos, while alternative formation scenarios that do not involve DM have also been suggested (e.g., Conroy et al. 2011; Naoz & Narayan 2014; Lake et al. 2023, 2025). In fact, recent cosmological simulations of realistic mechanisms have shown that both scenarios are possible (Trenti et al. 2015; Kimm et al. 2016; Ricotti et al. 2016; Keller et al. 2020; Ma et al. 2020), initiating a debate on the possibility that GCs originate via multiple formation channels. Observationally, one of the main reasons why the understanding of GC formation remains elusive is that most of the consequences of these formation channels are exhibited in the outskirts of the GCs, where stellar densities are extremely low and it therefore becomes difficult to distinguish between bona fide GC stars and MW foreground and background stars.

For some time now, it has been reported that some of the outer halo GCs of the MW possess stars located beyond their King tidal radius (r_t), i.e., the observable radius beyond which stars are gravitationally unbound to the cluster (e.g., Carballo-Bello et al. 2011; Sanna et al. 2012; Marino et al. 2014; Bellazzini et al. 2014; Kuzma et al. 2018; Zhang et al. 2022; Yang et al. 2022; Kuzma et al. 2025). These stellar populations, when present, are generally attributed to one of two phenomena: (1) the formation of tidal tails, which are symmetric streams originating from the cluster due to the tidal pull of the MW (see, e.g., Odenkirchen et al. 2001 for Palomar 5), or (2) they correspond to a diffuse stellar envelope surrounding the cluster symmetrically at distances greater than expected based on their King tidal radii (see, e.g., Kuzma et al. 2016 for NGC 7089). In the latter case, it is currently unclear where these stars originate and how they manage to remain at the periphery of the cluster. A possible explanation involves GCs forming within DM halos, which allows DM to prevent stars from escaping the system’s self-gravitational potential. This scenario has been explored over the years by other authors, and models of GCs with DM halos have been proposed (e.g., Moore 1996; Mashchenko & Sills 2005; Ibata et al. 2013; Vitral & Boldrini 2022; Carlberg & Grillmair 2022). In particular, Peñarrubia et al. (2017) predict that the presence of DM imprints discernible effects on the kinematics and photometric properties of the GCs. Specifically, clusters with higher DM content tend to be more extended, with an inflated radial velocity dispersion profile at large radii and an outer density profile characterized by a power law with index γ > –3, which asymptotically approaches –3 as the DM mass increases.

In this vein, one interesting case study of outer halo GCs is NGC 5824. This cluster, first discovered in 1925 (Innes 1925), is located 32.1 kpc from the Sun; it is the second brightest of the outer halo clusters, with M_V = –9.3 ± 0.04 (Muñoz et al. 2018a), and has a mass of ≈10⁶ M_⊙ (Yuan et al. 2022). One particular point of interest is the fact that its density profile is not well described by a King profile (Carballo-Bello et al. 2011) but instead seems well fitted by a power law (Sanna et al. 2014). In addition, several studies report a large extension (Kuzma et al. 2018; Yang et al. 2022). Muñoz et al. (2018b) determined that the spatial extent of NGC 5824 exceeds 30 times its Sérsic effective radius (r_e,s) and that NGC 5824 has the highest concentration parameter (r_t/r_c = 172, with r_c being its core radius) of the whole sample of the MW outer halo GCs, despite having a relatively small effective radius.

In this study, we focused on testing the predictions of the Peñarrubia model for the outer slopes of the density profile of NGC 5824 by comparing it with the observed profile based on deep photometry and proper motion data. In addition, we analyzed the luminosity function (LF) of NGC 5824 to determine whether there is a diffuse stellar envelope and thus whether the primary assumption of the Peñarrubia et al. (2017) model is satisfied. We used the public photometry obtained by Muñoz et al. (2018a) and the Gaia Data Release 3 (DR3) catalog, and supplementary Dark Energy Camera (DECam) photometry from Kuzma et al. (2018). Furthermore, we contrasted the density profile behavior of NGC 5824 with that of NGC 2419, the largest and brightest of the outer halo clusters with M_V = –9.35 ± 0.03 and r_e,s ~ 25.7 ± 0.2 pc, r_t = 227 pc (Muñoz et al. 2018a), to provide a reference point for a cluster with a lower concentration parameter (r_t/r_c = 50). The contrast between these two clusters is evident not only due to their significantly different concentration parameters, but also because they represent two extremes in terms of the Sérsic index and exhibit clearly different central surface brightnesses. Specifically, NGC 5824 has a Sérsic index of n = 3.82 ± 0.05 and μ_V,0 = 11.15 ± 0.08 mag arcsec⁻², whereas NGC2419 has n = 1.71 ± 0.02 and μV₀ = 18.83 ± 0.05 mag arcsec⁻² (Marchi-Lasch et al. 2019).

The paper is organized as follows. In Sect. 2 we summarize the data acquisition and data reduction process presented by Muñoz et al. (2018a). Section 3 describes our methodology, which involves using multiple independent datasets to construct density profiles and fit King and power-law models for comparison with the prediction of Peñarrubia et al. (2017). Additionally, we analyze LFs to assess the cluster’s spatial symmetry and extent. Our results are presented in Sect. 4. In Sect. 5 we discuss the implications of our findings by comparing them with the literature, and finally in Sect. 6 we present our conclusions.

2 Data and reduction

In this study, we analyzed the photometric and proper motion data of the GC NGC5824 (α₀ = 225.9943deg, δ₀ = –33.0685 deg) and compared them with the predictions made by the Peñarrubia et al. (2017) model. In addition, as a control case, we analyzed the photometric properties of NGC 2419 (α₀ = 114.5354deg, δ₀ =+38.8819 deg), the brightest and one of the most distant outer halo clusters. The data used for both clusters come from the survey of the MW outer halo satellites carried out by Muñoz et al. (2018b). This survey used the MegaCam imager on both the Magellan Clay telescope and the Canada-France-Hawaii Telescope (CFHT), providing deep imaging in the g and r bands over a wide field. The typical limiting magnitudes for point sources at a confidence level of 5σ are $g_{lim}\simeq 25.6$ and r_lim ∝ 25.3 in AB magnitudes. On the other hand, NGC2419 was observed with CFHT. The MegaCam on Clay covers an area of 0.4 deg × 0.4 deg on the sky, while on CFHT it covers 1 deg × 1 deg.

The MegaCam data were preprocessed following the procedure specified in Muñoz et al. (2018b), which can be briefly described as follows: The initial astrometric and photometric solutions were refined. For the astrometric calibration of the GCs in this paper, the solutions were improved using the Gaia DR1 catalog, achieving root-mean-square uncertainties of ∼0″.04-0″.06. Point-spread function photometry was then carried out using DAOPHOT/ALLSTAR and ALLFRAME packages. Finally, photometric calibration was performed by determining zero points and color terms using Sloan Digital Sky Survey objects that overlapped with the survey, either as direct calibrators or as secondary standards.

Additionally, to clean the catalog of unresolved background galaxies and artifacts, we carried out a cut in the DAOPHOT morphological parameter sharp. This parameter is a measure of how much broader or narrower an object’s light profile is compared to that of a point source. In essence, sharp reflects how widespread or narrow the point spread function is. We applied cuts following upward and downward exponential functions, as shown in the top and bottom panels of Fig. 1 for NGC 5824 and NGC 2419, respectively.

In the case of NGC 5824, to improve our ability to detect member stars extending to the outskirts of the cluster, we complemented these data with G - and G_rp -band photometry and proper motions from the Gaia DR3 archive¹ (Gaia Collaboration 2023). This additional data cover a circular region of 60′ from the cluster center, with a limiting magnitude of G_lim = 20.7 mag (Hodgkin et al. 2021). We cleaned this catalog by selecting objects with a probability of being a single star, but not a white dwarf, greater than 0.8, according to the Gaia DR3 Data Source Classification System (Delchambre et al. 2023). In addition, we incorporated DECam photometry in the g and i band from Kuzma et al. (2018), who traced NGC 5824 stars out to 230 pc by analyzing its 2D density map. For this dataset, we restricted the covering region up to 35′ from the cluster center.

3 Methods

Given that in the outer regions of the cluster the stellar densities are much lower, making it difficult to assess the real extent of the cluster stars, we computed the radial density profiles through four different tests. Each test corresponds to a different approach for identifying cluster stars, from which radial surface density profiles were constructed to perform King and power-law fits.

The first two tests focus on main sequence (MS), red giant branch (RGB), and blue horizontal branch (BHB) stars, using MegaCam g- and r-band photometry. The third test targets RGB and BHB stars, and is based on Gaia DR3 G - and G_RP-band photometry. Additionally, proper motion constraints in RA and Dec from Gaia DR3 are incorporated in the second and third tests to further refine cluster membership. Finally, the fourth test uses DECam g- and i-band data, also targeting MS, RGB, and BHB stars. A summary of the characteristics of each test is presented in Table 1.

Independently of these tests, we constructed LFs using MegaCam and DECam g-band photometry. Incorporating DECam data allows us to probe a larger spatial area than MegaCam, enabling the use of greater radial distances to establish the cluster’s maximum detectable extent.

Fig. 1 sharp as a function of radius for NGC 5824 (top panel) and NGC 2419 (bottom panel). Gray dots represent all the detections within the MegaCam field. Blue dots indicate upward and downward exponentials applied as the sharp filter. Black dots denote all the detections that survived the sharp filtering. Purple dots highlight the stars selected for Test 1, which include MS, red giant, and BHB stars.

3.1 Star membership identification

For the first test, we used the photometric catalog of Clay-MegaCam from Muñoz et al. (2018b) and constructed g versus g - r color-magnitude diagrams (CMDs). We used these CMDs to select MS and RGB stars with the highest probability of belonging to NGC 5824. For the MS selection, we considered stars with 21.5 ≤ g ≤ 24. This upper magnitude limit was chosen because, first, at this magnitude, the completeness level is ≈90% (see Fig. 2 from Muñoz et al. 2010 for a reference), and second, it minimizes the contamination from unresolved galaxies, which begin to dominate beyond g = 22, while maintaining low photometric errors (see Fig. 2). To refine the MS selection, we incorporated a dynamic maximum Euclidean distance in the CMD space, d_dyn, from the best-fit isochrone to the MS, that is, a distance that increases to fainter magnitudes, which follows the equation $$d_{dyn}=d_0+k|g-g_{min}|,$$(1)

where d₀ is the base minimum distance, measured in mag, k is a dimensionless scale factor, and g is the magnitude value in the g band. After some experimentation, we adopted the following values: d₀ = 0.11, k = 0.0095, and g_min = 21.5. The values chosen here, and for subsequent star selections, best trace the observed features of the sequences.

The isochrones used in this test and hereafter were obtained from the Dartmouth Stellar Evolution Database² (Dotter et al. 2008). Additionally, note that the adopted isochrone parameters are used only to facilitate the GC’s star selection, rather than to determine any physical properties of the cluster. Particularly, for this test we used the CFHT-MegaCam ugriz photometric system, with age = 13Gyr, [Fe/H] = –1.91, ${\mathrm{d}}_{\odot }=42\mathrm{k}\mathrm{p}\mathrm{c}$ and E(g – r) = 0.25.

We considered RGB stars as those with 13.6 ≤ g ≤ 21.5 and selected them with the same dynamical distance criterion from Eq. (1), but with different parameters: d₀ = 0.04, k = 0.0027, and g_min = 13.6. We used the same isochrone, but with a different extinction value of E(g – r) = 0.24. We also included stars in the BHB region of the CMD, a region typically with very low MW background contamination. BHB stars were selected by defining a specific CMD region that encloses the BHB locus, with boundaries in color g - r and in magnitude g sufficiently narrow to minimize contamination. Specifically, stars with 18.38 ≤ g ≤ 18.9 and –0.09 ≤ g - r ≤ 0.18 were classified as BHB stars.

Thus, for this first test, we considered stars from the MS, RGB and BHB. The left panel of Fig. 3 shows the cluster’s CMD and the stars selected for this test. To evaluate the sensitivity of the results to the selection criteria, we applied different constraints on the sharp parameter and MS star selection.

First, we applied a more restrictive sharp filter while maintaining the original MS selection shown in Fig. 3. Then, we compared two extreme cases, one minimizing and the other increasing photometric errors, through a variation in the depth of the MS selection between g = 23 and g = 25 while keeping the original sharp filtering illustrated in Fig. 1. Finally, the original MS depth (g = 24) and sharp filter (Fig. 1) were fixed but with adjusted MS width, testing both narrower (using d₀ = 0.07) and wider (using d₀ = 0.15) selections than that shown in Fig. 3.

For the second test, we combined data from the Clay catalog with data from Gaia DR3 (Gaia Collaboration 2023) to aid the identification of distant cluster stars using proper motion and CMD properties. We first selected stars from the Clay-MegaCam catalog with g < 21.5 that have Gaia DR3 proper motion measurements in both right ascension (μ^*_α)³ and declination (μ_δ).

This selection was made using the TOPCAT⁴ software, where we matched stars from Clay-Megacam with the ones in the Gaia DR3 catalog by position (right ascension and declination), allowing a maximum error of 1" between the positions. In contrast to the first test, we did not include the MS stars since the Gaia DR3 limiting magnitude is G_lim = 20.7 mag (Hodgkin et al. 2021), which is brighter than for Clay-MegaCam data.

From this subset, the first filter used to select NGC 5824 star candidates was based on the CMD positions of the stars. We considered the stars from the BHB and RGB of the CMD as shown in the left panel of Fig. 3. For this selection, we adopted the same methodology as for Test 1.

We subsequently used the proper motion information to further refine the filter used to select cluster stars. Since the BHB region of the CMD is relatively clean of MW stars, we used this group of stars to provide a reliable representation of NGC 5824’s proper motion. Figure 3 (upper right panel) shows a well-defined area in proper motion space, further confirming their status as representative stars.

The lower right panel of Fig. 3 shows the proper motion of RGB stars (marked in red in the CMD). The region in proper motion space is similar but with a wider distribution, likely due to the presence of MW foreground stars and the fact that the RGB extends to fainter magnitudes than the BHB and therefore the proper motion errors for these stars are larger. To refine the selection of NGC 5824 star candidates, we defined a proper motion box (using the BHB star candidates region as reference) to be: μ^*_α ∈ [−1.65, 0.65] mas yr⁻¹ and μ_δ ∈ [−2.75, −1.8] mas yr⁻¹, which is in agreement with the systemic proper motion found by Yang et al. (2022).

In the third test, we used both photometry and proper motions from the Gaia DR3 catalog only, up to 60′ from the cluster’s center. The procedure to select GC stars is as follows: First we filtered the stars by proper motion within the same specified range detailed above. Then, we selected BHB and RGB stars according to their CMD location and their proximity to the best-fit isochrone, respectively. For this isochrone, we used the Gaia DR2 Revised photometric system, with age = 13Gyr, [Fe/H] = –1.5, d₀ = 35kpc and E(g – r) = 0.11. To facilitate the RGB selection, we adopted the same dynamical distance criterion of Eq. (1) used in Test 1 and 2, using d₀ = 0.07, k = 0.005 and g_min = 14.8.

Finally, for the fourth test, we used DECam photometry from Kuzma et al. (2018) covering up to 35′ from the cluster’s center, which corresponds to 2.5 times its King tidal radius, according to the measurement by Muñoz et al. (2018b). The star selection procedure is similar to Test 1. The isochrone used was in the DECam photometric system, with age = 13 Gyr, [Fe/H] = –1.91, ${\mathrm{d}}_{\odot }=32.1\mathrm{k}\mathrm{p}\mathrm{c}$ and distinct extinction values for MS and RGB selections, E(g – r) = 0.12 and E(g – r) = 0.06, respectively. MS stars, with 21.05 ≤ g ≤ 22.3, and RGB stars, with 17 ≤ g ≤ 21.04, were selected incorporating the dynamical distance criterion, using d₀ = 0.09, k = 0.03 and g_min = 21.05 for MS stars, and d₀ = 0.035, k = 0.002 and g_min = 17 for RGB stars.

For the case of NGC 2419, used as a control GC, we followed a procedure similar to that used in the first test (Clay-MegaCam catalog). We constructed the NGC 2419 CMD using the sharp -filtered stars from Fig. 1. Then, we selected the stars located in the BHB and RGB, the latter with 18 ≤ g ≤ 23.38, selected within a fixed distance of 0.05 mag from the best-fit isochrone, and stars in the MS, with 23.4 ≤ g ≤ 25, within a dynamic distance d_dyn from the isochrone as specified in Eq. (1) using d₀ = 0.013, k = 0.005 and g_min = 23.4. The isochrone is also in the CFHT-MegaCam ugriz photometric system, as for Tests 1 and 2, but with age = 12Gyr and [Fe/H] = –1.5, ${\mathrm{d}}_{\odot }=98\mathrm{k}\mathrm{p}\mathrm{c}$ and E(g – r) = 0.04. We extended the MS selection to fainter magnitudes than for NGC 5824 because adopting a brighter limit would trace a limited portion of the sequence. Therefore, this choice is motivated by a trade-off between magnitude depth and minimizing photometric errors, which remain low at magnitude g = 25 as illustrated in lower panels of Fig. 2.

Fig. 2 Photometric errors as a function of MegaCam magnitudes for NGC 5824 (top panels) and NGC 2419 (bottom panels). Left panels show the g band and right panels the r band. The dashed vertical red line indicates magnitude 24 for NGC 5 824 and 25 for NGC 2419 in both bands.

Fig. 3 Left: CMD of NGC 5824 with Clay photometry. Black dots correspond to the full sharp-cleaned data. Light blue dots are BHB candidates. Red dots are RGB candidates. Blue and yellow dots correspond to the isochrone. Right: NGC 5824 proper motion (PM) diagrams. Gray dots correspond to the full data with PM measurements. Light blue dots (upper panel) are BHB candidates. Red dots (lower panel) are RGB candidates.

3.2 Density profiles

To construct the radial density profiles, we spatially divided the previously selected GC stars into several rings, centered on the GC coordinates, using a geometric progression, and computed the surface number density in each ring. We subtracted the average background surface density contribution Σ_b from each surface density value of the GC. This background is determined in different ways depending on the area covered by the selected stars.

In the case of NGC 5824, for the first and second tests, we determined the background by averaging the surface density measured in rectangular areas at the corners of the cluster’s field. The idea is to use an area as far as possible from the cluster’s center, in order to include as few cluster stars as possible (assuming they extend up to the edge of our covered area). For the third and fourth tests, considering the area covered extends out to a radius of 60′ and 35′, respectively, we calculated the background using the density of the maximum concentric circular annuli. For NGC 2419, we used the same method as for the first and second tests for NGC 5824.

An important consideration is that the photometry used in this study suffers from overcrowding effects at the innermost regions of the clusters. For this reason, we used the surface brightness profiles of NGC 5824 and NGC 2419 from Trager et al. (1995, hereafter Trager) to complement our profiles at small radii. To ensure consistency in the units, we converted this brightness density profile to a number density profile, i.e., converting measurements from flux units to star-count units. A cubic spline interpolation between Trager data and Clay/CFHT-MegaCam, Gaia DR3, or DECam data, depending on the case, was conducted to facilitate the conversion. With the interpolated surface number density values, the following conversion from flux to star-counts was used (Carballo-Bello et al. 2011): $$\rho =c-0.4\mu ,$$(2)

where ρ represents the interpolated Trager number surface density values, μ is the associated Trager flux surface density value, and c is a constant.

Instead of providing statistical uncertainties, Trager subjectively assigned to each data point i a weight w_i ∈ [0,1], reflecting the quality of the measurement. These weights, however, are not clearly related to conventional statistical error bars. To estimate realistic uncertainties for the surface brightness values δμ_i, where i represents each Trager data point, the methodology proposed by McLaughlin & van der Marel (2005) was implemented, which assumes that the errors are inversely proportional to the assigned weights, i.e., δμ_i ∝ 1/w_i. This approach provides a statistically and physically consistent calibration of the uncertainties for each cluster. The proportionality constant for the case of NGC 5824 is 0.106. Therefore, the individual uncertainties for each Trager density data point are then computed as δμ_i = 0.106/w_i.

Finally, we fitted a King profile to the number density profiles for each test of NGC 5824 and NGC 2419. This King profile follows the following equation: $${\mathrm{\Sigma }}_{King}(R)={\mathrm{\Sigma }}_0[{(1+\frac{R^2}{R_c^2})}^{-\frac{1}{2}}-{(1+\frac{R_t^2}{R_c^2})}^{-\frac{1}{2}}],$$(3)

where Σ₀ is the central surface density of the cluster, R_c is its core radius and R_t is its King tidal radius.

For the parameter estimation procedure, we used the emcee library minimizing a χ² and exploring the parameter space following the rules of Markov chain Monte Carlo (MCMC). In the end, we obtained a corner plot that shows the distribution of each fitted parameter and the relations between them. These plots are shown in Figs. 4, A.1, and B.1 for Tests 1, 2, and 3, respectively.

To characterize the outskirts of the clusters region, we fitted a power-law profile only to the points corresponding to the outer region of the profiles for the four different tests. This profile follows the equation $${\mathrm{\Sigma }}_{pl}(R)=\sigma R^{\gamma },$$(4)

where γ is the power-law index that represents the slope of the profile, and σ is a normalization factor. Additionally, the χ² values for both the King (${\chi }_{King}^2$) and power-law (${\chi }_{pw}^2$) models were computed for each test to evaluate which model better fits the data.

3.3 g-Band luminosity function

To estimate the approximate extent of NGC 5824 and to assess the presence of a diffuse stellar envelope, we analyzed cumulative g-band LFs, tracing the turn-off point and MS of the cluster across different annuli. These LFs were compared with the normalized background cumulative LF to assess the radial distance from the cluster center at which we can reliably detect cluster members above the background level. For this analysis, we separately tested MegaCam and DECam photometry.

The selection of stars used to construct the cumulative LFs follows a color range in (g - r) for MegaCam and (g - i) for DECam, and a g magnitude filter as shown in Fig. 5. Specifically, MegaCam stars follow 18 ≤ g ≤ 24 and 0.3 ≤ g - r ≤ 0.5, while DECam stars follow 18 ≤ g ≤ 22.34 and 0.3 ≤ g - i ≤ 0.5. The chosen magnitude range is defined to ensure that we are accurately tracing the background, the turn-off and the MS of the cluster, while the chosen color range maximizes the contrast between cluster and background stars.

The MegaCam and DECam backgrounds were defined consistently with the approach used in the density profile section: for MegaCam, four squares located in the corners of the Clay-MegaCam field of view were used, while for DECam, the maximum concentric circular annulus was employed. Background stars were further divided into four quadrants. We then constructed and characterized the LFs of both the total background and its individual quadrants.

Using the selected cluster stars, we constructed the cumulative LFs for each annulus, each with a size of 3′ for the MegaCam case (see Fig. 7), and 5′ for the DECam case (see Fig. 9). To compare the cluster LF at different annuli with the background LF, the background LF was normalized at g = 20, as there is no MS at this magnitude in either MegaCam or DECam photometry. We defined a contrast parameter (C) to compare the cluster and background cumulative LFs, defined as the ratio between the maximum LF value of the cluster ($N_c^{max}$) and the background ($N_b^{max}$): $$C=\frac{N_c^{max}}{N_b^{max}}.$$(5)

To assess the symmetry of the cluster, we analyzed the LF of the outermost ring that remained distinguishable from the background (see Fig. 8 for the MegaCam case and Fig. 10 for DECam case). This ring was divided into four quadrants: the upper-right, upper-left, lower-right, and lower-left quadrants. As a result of this division, the number of cluster stars decreased; therefore, we adopted a larger number of bins for constructing these LFs. To maintain consistency, the background LF was also constructed with the same number of bins. It is important to note that the background stars were not divided into quadrants, instead they represent the same overall background.

Fig. 4 NGC 5824 results for Test 1. Top left: CMD and a star-count map of the cluster. Green, red, and blue dots correspond to the MS, RGB, and BHB star candidates, respectively, and the gray circle represents the best.fit Rt parameter for this test. Note that MS candidates are omitted from the star-count map because they occupy the entire MegaCam field. Top right: Observed surface number density profile, constructed using MS, RGB, and BHB stars, along with data from Trager et al. (1995) for the cluster’s innermost region. Bottom left and bottom right: corner plots obtained for the King profile and power-law fits, respectively. The dashed black lines correspond to percentiles 0.16, 0.5, and 0.83, in that order, from left to right. Rc and Rt are expressed in arcminutes. Σ0 is expressed in arcmin−2. γ is dimensionless, and σ is a normalization factor of the power-law model.

Fig. 5 CMDs of NGC 5824 based on MegaCam photometry (left) and DECam photometry (right). In the MegaCam CMD, black dots correspond to the full sharp-cleaned data, while in the DECam CMD, black dots correspond to stars within 1′ ≤ R ≤ 35′, where R is the angular distance from the cluster center. The yellow shaded regions indicate the color and magnitude selection used to construct the cumulative LFs.

Fig. 6 Background LFs of NGC 5824 based on MegaCam (upper panel) and Kuzma et al. (2018) DECam data (lower panel). Each panel contains a star-count map (left subpanel), the cumulative LF of the entire background (top right), and the cumulative LFs for each corner of the background stars, as indicated in the star-count map (bottom right). In the left subpanels, black dots correspond to the full sharp-cleaned data (only for the MegaCam data). Colored dots represent the color filtered stars used as background. In the right subpanels, the solid black line corresponds to the total background divided by four. The shaded areas in both panels and the black error bars in the bottom panels represent the corresponding Poisson errors, √Nstars.

Fig. 7 MegaCam cumulative LFs of NGC 5824 for different rings. The green line and the black line represent the cluster and the background, respectively. The shaded area corresponds to the Poisson error bars ($\sqrt{N_{stars}}$). C is the contrast parameter.

Fig. 8 MegaCam cumulative LFs for quadrants of the ring 13′ < R < 16′. The colors have the same meanings as in Fig. 7.

4 Results

The star selection and density profile results for NGC 5824 are presented in Figs. 4, A.1, B.1, and C.1 for Tests 1 through 4, respectively, and in Fig. D.1 for NGC 2419. For each test, the best-fitting parameters reported are the King limiting radius, R_t, and the slope γ. For NGC5824, R_t ~ 21′ for Test 1, R_t ~ 17′ for Test 2, R_t ~ 29′ for Test 3, and R_t ~ 21′ for Test 4. The power-law index remains consistent across all tests, with an average γ = –2.5 ± 0.1. In the case of NGC 2419, R_t ~ 10′, and the slope is steeper, with γ ~ −4.5. Table 2 summarizes the bestfit parameters for all tests. For both clusters, the ${\chi }_{pw}^2$ values are consistently lower than ${\chi }_{King}^2$, meaning that a power-law profile better describes the outer region of the cluster.

The evaluation of the sensitivity for Test 1 showed that, regardless of the star selection criteria, the King and power-law profile fits remain stable, with only minor variations that do not affect the main conclusions. Specifically, the power-law index (γ) varies between –2.7 and –2.4, consistent with the statement by Peñarrubia et al. (2017). The consistency of power-law fits across all tests further supports the robustness of these results, as they are derived from independent datasets.

Notably, in Test 2, NGC 5824 stars are detected out to a radius of 147 pc (∼15′.7 in the MegaCam field). This extended stellar distribution is clearly visible in the star-count map shown in the upper-left panel of Fig. A.1.

Finally, with respect to the NGC 5824 g-band LFs, the MegaCam LFs for rings 1′ < R < 4′ to 16′ < R < 18 are shown in Fig. 7, with LFs of the quadrants for the ring 13′ < R < 16′ displayed in Fig. 8. Similarly, the DECam LFs for rings 5′ < R < 10′ to 30′ < R < 35 are shown in Fig. 9, with LFs of the quadrants of the ring 15′ < R < 20′ presented in Fig. 10.

Fig. 9 DECam cumulative LFs of NGC 5824 for different rings. The colors have the same meanings as in Fig. 7. The shaded area corresponds to the Poisson error bars ($\sqrt{N_{stars}}$). C is the contrast parameter.

Fig. 10 DECam cumulative LFs for quadrants of the ring 15′ < R < 20′. The colors have the same meanings as in Fig. 7.

5 Analysis and discussion

As seen in multiple studies (e.g., McConnachie 2012; Muñoz et al. 2018b), GCs typically do not extend beyond ∼100 pc, although exceptions have been observed. For example, NGC 1851 has an extended stellar envelope of ∼240 pc (Olszewski et al. 2009; Kuzma et al. 2018; Zhang et al. 2022). Using the MegaCam data, we find that the LF analysis of NGC 5824 detects MS stars out to the edge of the field, at ∼16′ from its center (equivalent to ∼150 pc), as illustrated in Fig. 7. In Fig. 8, we use the MegaCam data to show a deeper analysis of the cluster’s outermost region where the MS is still detectable, from which we can conclude that this region is tracing a symmetric structure. We note that the increase in the contrast parameter C in the east panels of this outermost ring could be attributed to the spatial variation in the background shown in Fig. 6.

On the other hand, Kuzma et al. (2018) detected cluster stars out to 230 pc at a 3σ level of significance from 2D distribution imaging. Since MegaCam data have a limited radial coverage compared to Kuzma et al. (2018), we incorporated DECam data from their study into our analysis to test whether stars beyond MegaCam coverage can be detected using the LF methodology. According to Fig. 9, cluster stars are clearly detected up to at least 20′, while the ring between 20′ and 25′ suggests a marginal presence of stars beyond 20′, as the cluster signal is slightly above the background in this region, consistent with the detection limit of ∼24′ reported by Kuzma et al. (2018).

The King profile fits displayed in the top-right panels of Figs. 4, A.1, and B.1 (corresponding to Test 1, 2, and 3, respectively) reveal an absence of a clear tidal cutoff. This absence may suggest that the MW has not yet truncated its stellar distribution. For comparison, we computed the cluster’s Jacobi radius following Eq. (1) from (Burkert 1997), with the orbital parameters of the clusters taken from Bajkova & Bobylev (2021). The Jacobi radius corresponds to the theoretical radius beyond which stars are unbound to the cluster. We assumed a Galactocentric distance of 25.6 kpc for NGC 5824 and used the enclosed Galactic mass at this radius as described by Eq. (3) in Burkert (1997), whose result is consistent with more recent MW mass profiles (see Fig. 11 from Shen et al. 2022). We estimated the Jacobi radius to be ∼19′, in line with the fitted King limiting radii. The fact that cluster stars are detected beyond these radii reinforces the notion that NGC 5824 is unusually extended, presenting a “halo” of stars beyond its expected tidal cutoff radius.

One possibility to explain these observations is to recall the Peñarrubia et al. (2017) model, where the authors simulate GCs with varying amounts of DM, ranging from none, intermediate, and higher DM content. The evolution of the particles in the simulations is analyzed as they interact with the gravitational field of the system (stars+DM). The study reveals a clear trend in which the response of the density profile depends significantly on the amount of DM. For the zero DM case, the cluster develops steeper outer density profiles, while those with increasing DM content become more spatially extended (see Fig. 1 from Peñarrubia et al. 2017) and exhibit progressively shallower density profile slopes asymptotically approaching a value of —3, without showing a clear tidal cutoff (see the upper panels of Fig. 2 from Peñarrubia et al. 2017).

Figures 4, A.1, B.1, and C.1 show that, at large distances, a power-law profile with a slope of γ~  -2.6 is more appropriate than a King profile, which is consistent with the predictions of Peñarrubia et al. (2017). In contrast, NGC 2419 displays a much steeper slope of γ~ −4.5 (see Fig. D.1). This steep profile is in line with the no DM scenario of Peñarrubia et al. (2017). Moreover, the fact that power-law profiles systematically outperform King models in the outer regions, as in this case, is not new. Carballo-Bello et al. (2011) already reported this behavior for several clusters, attributing it to the improved ability of modern surveys to detect stars in their low-density outskirts. Altogether, we argue that the findings of this work suggest that NGC 5824 is a strong candidate to be embedded in a DM halo, which could account for its unusual extension and photometric properties. We note that there are other potential explanations for this observed extension. For example, if NGC 5824 were approaching its pericenter, the Jacobi radius would shrink as the GC approaches the Galactic center, making stars that were initially bound now unbound. However, according to Bajkova & Bobylev (2021), NGC 5824’s pericenter is ∼ 14.26 kpc, and currently the cluster is located at 25.6 kpc from the galactic center, not quite near the perigalacticon.

In a contrasting view, Yang et al. (2022) proposed the presence of tidal debris in the vicinity of NGC 5824. Notably, they claim the existence of a trailing tail observed to extend ∼50 de g from the cluster’s center. Additionally, the authors propose that the Triangulum stream serves as a component of the leading tail associated with NGC 5824. They compare the leading tail of a simulated stream of NGC 5824 with Triangulum stream stars (see their Fig. 3). Although this comparison is tantalizing, the phase space distribution of the mock stars spans approximately 15 de g in RA, which implies that the model is not precisely constrained. The authors argue that Triangulum stars are consistent with being associated with NGC 5824 according to their metallicity and CMD positions. However, we consider that their match-filter map (see their Fig. 8), which is presented as the detection of the trailing tail of the cluster, shows partial detections of a structure that we cannot directly associate with the cluster as a distinctive structure indicative of a NGC 5824 stream. More importantly, in their discussion, the authors acknowledge a significant inconsistency between Triangulum stream candidate stars and NGC 5824 stars, i.e., according to their Fig. 10, most of the proposed members of the supposed leading tail stars lie at nearly twice the distance of the model prediction, corresponding to a difference of approximately 1.5 magnitudes. For all of the above, in our view the analysis presented in Yang et al. (2022), the existence and detection of tidal tails is plausible but not definitive to explain the extended stellar envelope of the cluster.

Although several mechanisms have been proposed to explain how stars in GCs reach large radii (e.g., Weatherford et al. 2023), including evaporation, binary interactions, and supernova explosions, we did not attempt to identify the origin of such configurations. Instead, our analysis focused on the empirical fact that these stars exist and remain gravitationally bound to the cluster center at large distances, which we interpret as indicative of the presence of DM.

Based on the results presented above, we argue that NGC 5824 is a good candidate to further investigate the presence of DM. Strong support for this hypothesis requires high-quality radial velocity profile data (see Figs. 2 and 3 from Peñarrubia et al. 2017), which are not yet available. Future studies addressing this problem could be key to resolving this question.

If the results of future observations were to confirm that DM is indeed present, one might ask whether NGC 5824 should be classified instead as a dwarf galaxy, but the answer is not straightforward. Other factors, most notably NGC 5824’s chemically homogeneous stellar population, support the classification as a GC, leaving the door open to a reconsideration of what the fundamental boundaries between GCs and dwarf galaxies are.

6 Conclusions

In this study, we analyzed the surface density profile of the GC NGC 5824 using observations from the MegaCam imager at the Magellan Clay telescope in the g and r bands, with a limiting magnitude of g_lim ≃ 25.6. These data were complemented by the Trager catalog for the inner regions of the cluster, and by Gaia DR3 G and G_rp bands for the outskirts. We also used Gaia DR3 proper motions to constrain cluster stars. This combination of deep and homogeneous imaging with high-precision astrometry enhanced our ability to probe the extended stellar profile of NGC 5824. Additionally, we incorporated DECam data from Kuzma et al. (2018), who report a cluster detection up to 230 pc, which is highly unusual for a GC.

Our primary goal was to assess the model predictions from Peñarrubia et al. (2017) regarding the outer slopes of GCs with or without DM content. We also examined the cumulative LF to assess the cluster’s spatial symmetry and extent. As a control case, we also analyzed the outer halo GC NGC 2419, for which we also had observations from the MegaCam imager at the CFHT.

Previous studies of the GC NGC 5824 indicate that the cluster presents multiple unusual characteristics, including an extended stellar distribution and the absence of an observed King tidal cutoff. Our results align with these findings, further affirming the unique nature of NGC 5824 in this regard. From our results shown in Figs. 4, A.1, B.1, and C.1, we find that in the outermost regions of NGC 5824, a power-law profile is more appropriate than a King profile. Moreover, we estimate an average slope of –2.6 ± 0.1 based on several different tests using three independent datasets, which is consistent with the predictions. Additionally, the cumulative LFs displayed in Figs. 7 and 8 for MegaCam data, as well as in Figs. 9 and 10 for DECam data from Kuzma et al. (2018), further corroborate the cluster’s spatial symmetry and extended structure.

Based on these findings, we conclude that the density profile of NGC 5824 is unusual and corroborate the predictions made by Peñarrubia et al. (2017). The absence of truncation in its density profile supports the idea that it still retains a DM component. However, kinematic data are necessary to confirm this. Future spectroscopic observations of NGC 5824 stars in the outskirts will provide radial velocities, which will be critical for distinguishing whether the extended structure is due to DM or other dynamical effects.

Acknowledgements

P.B.D., B.M, R.R.M. and V.S. gratefully acknowledges support by the Agencia Nacional de Investigación y Desarrollo (ANID) Programa de Financiamiento Basal para Centros Científicos y Tecnológicos de Excelencia (BASAL) project FB210003. JAC-B acknowledges support from FONDECYT Regular N 1220083. P.B.K. is supported by a UKRI Future Leaders Fellowship (MR/S018859/1). This paper is based on observations obtained with the MegaPrime/MegaCam, a joint project of the Canada-France-Hawaii Telescope (CFHT) and CEA/DAPNIA, at CFHT, which is operated by the National Research Council (NRC) of Canada, the Institut National des Science de l’Univers of the Centre National de la Recherche Scientifique (CNRS) of France, and the University of Hawaii. This project used data obtained with the Dark Energy Camera (DECam), which was constructed by the Dark Energy Survey (DES) collaborating institutions: Argonne National Lab, University of California Santa Cruz, University of Cambridge, Centro de Investigaciones Energeticas, Medioambientales y Tecnologicas-Madrid, University of Chicago, University College London, DES-Brazil consortium, University of Edinburgh, ETH-Zurich, University of Illinois at Urbana-Champaign, Institut de Ciencies de l’Espai, Institut de Fisica d’Altes Energies, Lawrence Berkeley National Lab, Ludwig-Maximilians Universitat, University of Michigan, National Optical Astronomy Observatory, University of Nottingham, Ohio State University, University of Pennsylvania, University of Portsmouth, SLAC National Lab, Stanford University, University of Sussex, and Texas A&M University. Funding for DES, including DECam, has been provided by the U.S. Department of Energy, National Science Foundation, Ministry of Education and Science (Spain), Science and Technology Facilities Council (UK), Higher Education Funding Council (England), National Center for Supercomputing Applications, Kavli Institute for Cosmological Physics, Financiadora de Estudos e Projetos, Fundaçâo Carlos Chagas Filho de Amparo a Pesquisa, Conselho Nacional de Desenvolvimento Científico e Tecnológico and the Ministério da Ciência e Tecnologia (Brazil), the German Research Foundation-sponsored cluster of excellence “Origin and Structure of the Universe” and the DES collaborating institutions.

References

Appendix A Test 2 results

Fig. A.1 Same as Fig. 4 but for Test 2, which uses MegaCam photometry and Gaia proper motions. The top-left panel includes a proper motion diagram where gray dots represent the full dataset with PM measurements, and red and blue dots correspond to the RGB and BHB stars candidates, respectively, which are also shown in the CMD and star-count map.

Appendix B Test 3 results

Fig. B.1 Same as Fig. A.1 but for Test 3, which uses Gaia DR3 photometry and proper motion.

Appendix C Test 4 results

Fig. C.1 Same as Fig. 4 but for Test 4, which uses DECam photometry from Kuzma et al. (2018).

Appendix D NGC2419 results

Fig. D.1 Same as Fig. 4 but for NGC 2419. For the King profile fitting, we used density values up to the background level, corresponding to ≈9′ from the cluster center. This is motivated by the fact that density values below the background level exhibit high signal-to-noise ratios, causing the MCMC algorithm to misinterpret and include them in the fit, leading to an inadequate result.

All Tables

All Figures

Fig. 1 sharp as a function of radius for NGC 5824 (top panel) and NGC 2419 (bottom panel). Gray dots represent all the detections within the MegaCam field. Blue dots indicate upward and downward exponentials applied as the sharp filter. Black dots denote all the detections that survived the sharp filtering. Purple dots highlight the stars selected for Test 1, which include MS, red giant, and BHB stars.

In the text

	Fig. 2 Photometric errors as a function of MegaCam magnitudes for NGC 5824 (top panels) and NGC 2419 (bottom panels). Left panels show the g band and right panels the r band. The dashed vertical red line indicates magnitude 24 for NGC 5 824 and 25 for NGC 2419 in both bands.
In the text

Fig. 3 Left: CMD of NGC 5824 with Clay photometry. Black dots correspond to the full sharp-cleaned data. Light blue dots are BHB candidates. Red dots are RGB candidates. Blue and yellow dots correspond to the isochrone. Right: NGC 5824 proper motion (PM) diagrams. Gray dots correspond to the full data with PM measurements. Light blue dots (upper panel) are BHB candidates. Red dots (lower panel) are RGB candidates.

In the text

Fig. 4 NGC 5824 results for Test 1. Top left: CMD and a star-count map of the cluster. Green, red, and blue dots correspond to the MS, RGB, and BHB star candidates, respectively, and the gray circle represents the best.fit Rt parameter for this test. Note that MS candidates are omitted from the star-count map because they occupy the entire MegaCam field. Top right: Observed surface number density profile, constructed using MS, RGB, and BHB stars, along with data from Trager et al. (1995) for the cluster’s innermost region. Bottom left and bottom right: corner plots obtained for the King profile and power-law fits, respectively. The dashed black lines correspond to percentiles 0.16, 0.5, and 0.83, in that order, from left to right. Rc and Rt are expressed in arcminutes. Σ0 is expressed in arcmin−2. γ is dimensionless, and σ is a normalization factor of the power-law model.

In the text

Fig. 5 CMDs of NGC 5824 based on MegaCam photometry (left) and DECam photometry (right). In the MegaCam CMD, black dots correspond to the full sharp-cleaned data, while in the DECam CMD, black dots correspond to stars within 1′ ≤ R ≤ 35′, where R is the angular distance from the cluster center. The yellow shaded regions indicate the color and magnitude selection used to construct the cumulative LFs.

In the text

Fig. 6 Background LFs of NGC 5824 based on MegaCam (upper panel) and Kuzma et al. (2018) DECam data (lower panel). Each panel contains a star-count map (left subpanel), the cumulative LF of the entire background (top right), and the cumulative LFs for each corner of the background stars, as indicated in the star-count map (bottom right). In the left subpanels, black dots correspond to the full sharp-cleaned data (only for the MegaCam data). Colored dots represent the color filtered stars used as background. In the right subpanels, the solid black line corresponds to the total background divided by four. The shaded areas in both panels and the black error bars in the bottom panels represent the corresponding Poisson errors, √Nstars.

In the text

	Fig. 7 MegaCam cumulative LFs of NGC 5824 for different rings. The green line and the black line represent the cluster and the background, respectively. The shaded area corresponds to the Poisson error bars ($\sqrt{N_{stars}}$). C is the contrast parameter.
In the text

	Fig. 8 MegaCam cumulative LFs for quadrants of the ring 13′ < R < 16′. The colors have the same meanings as in Fig. 7.
In the text

	Fig. 9 DECam cumulative LFs of NGC 5824 for different rings. The colors have the same meanings as in Fig. 7. The shaded area corresponds to the Poisson error bars ($\sqrt{N_{stars}}$). C is the contrast parameter.
In the text

	Fig. 10 DECam cumulative LFs for quadrants of the ring 15′ < R < 20′. The colors have the same meanings as in Fig. 7.
In the text

	Fig. A.1 Same as Fig. 4 but for Test 2, which uses MegaCam photometry and Gaia proper motions. The top-left panel includes a proper motion diagram where gray dots represent the full dataset with PM measurements, and red and blue dots correspond to the RGB and BHB stars candidates, respectively, which are also shown in the CMD and star-count map.
In the text

	Fig. B.1 Same as Fig. A.1 but for Test 3, which uses Gaia DR3 photometry and proper motion.
In the text

	Fig. C.1 Same as Fig. 4 but for Test 4, which uses DECam photometry from Kuzma et al. (2018).
In the text

	Fig. D.1 Same as Fig. 4 but for NGC 2419. For the King profile fitting, we used density values up to the background level, corresponding to ≈9′ from the cluster center. This is motivated by the fact that density values below the background level exhibit high signal-to-noise ratios, causing the MCMC algorithm to misinterpret and include them in the fit, leading to an inadequate result.
In the text