Discovery of Extended Tidal Tails around the Globular Cluster Palomar 13

Stellar streams—remnants of tidally disrupted star clusters and satellite galaxies—provide information about the formation history and dark-matter distribution of the Milky Way (e.g., Johnston 1998; Helmi & White 1999; Bonaca & Hogg 2018). Dynamically cold streams from disrupting star clusters are also extremely sensitive to gravitational perturbations from massive substructures (e.g., Erkal & Belokurov 2015; Bonaca et al. 2019). The study of stellar streams is thus a promising avenue for studying the distribution of dark matter at subgalactic scales.

From our perspective in the Galaxy, the discovery and study of stellar streams requires wide-field sky surveys (e.g., Rockosi et al. 2002; Majewski et al. 2003). The population of known streams has increased substantially due to deep, well-calibrated photometric surveys such as the Dark Energy Survey (DES; DES Collaboration 2016; see Shipp et al. 2018), and our ability to characterize stream motions has been recently revolutionized by the Gaia mission (Gaia Collaboration et al. 2018; see, e.g., Malhan et al. 2018; Price-Whelan & Bonaca 2018; Shipp et al. 2019) and coordinated spectroscopic follow-up (e.g., the ${S}^{5}$ survey; Li et al. 2019).

The known population of stellar streams now comprises ∼70 candidates distributed throughout the Milky Way stellar halo.⁶ Most stellar streams likely originated from globular clusters, but the vast majority have no known surviving progenitors. In contrast, many clusters have extended extra-tidal features (Belokurov et al. 2006; Niederste-Ostholt et al. 2010; Sollima et al. 2011; Myeong et al. 2017; Sollima 2020; Dalessandro et al. 2015; Kuzma et al. 2018; Röser & Schilbach 2019), but very few appear to host extended thin, dense tidal tails (Rockosi et al. 2002; Grillmair 2019; Ibata et al. 2019). It is currently unclear how this discrepancy connects to the surface brightness evolution of tidal tails (e.g., Balbinot & Gieles 2018).

A well-known example of a thin stream is the globular cluster Palomar 5 (Pal 5) and its tidal tails (e.g., Rockosi et al. 2002), which span ∼25° on the sky (e.g., Bonaca et al. 2020) and contain a total stellar mass comparable to the surviving cluster (e.g., Ibata et al. 2017). The Pal 5 system is an archetype of stellar stream formation and has provided dynamical constraints on the dark-matter distribution around the Milky Way (e.g., Küpper et al. 2015). However, recent studies have found that the Pal 5 orbit and stream may be perturbed by the time-dependent influence of the Galactic bar (Pearson et al. 2017) and/or massive substructures (Erkal et al. 2017). These realizations motivate the study of other globular cluster streams as a way to disentangle the complex dynamical phenomena affecting streams in the inner Galactic halo.

In order to increase the sample of stellar streams with identified progenitors, we have initiated a systematic search for tidal debris structures using photometric data from the DESI Legacy Imaging Surveys (Dey et al. 2019). A first result from this search is the detection of extended tidal structures associated with the globular cluster Palomar 13 (Pal 13). Pal 13 is a low-luminosity globular cluster (M_V ∼ −2.8; Bradford et al. 2011) at a Galactocentric distance of ∼25 kpc (Côté et al. 2002). Several studies have previously suggested that Pal 13 is undergoing tidal disruption because of the spatial distribution of its blue straggler population (Siegel et al. 2001), its large radial velocity dispersion (Côté et al. 2002), and its extended radial profile (${r}_{1/2}=1\buildrel{\,\prime}\over{.} 27\pm 0\buildrel{\,\prime}\over{.} 16;$ Bradford et al. 2011). More recently, it was shown that Pal 13 member stars display a significant proper-motion scatter (using Gaia DR2 astrometry; Yepez et al. 2019), and its stellar population extends to almost twice the estimated Jacobi radius of ∼5'–10' (using photometry from the Dark Energy Camera Legacy Survey, DECaLS; Piatti & Fernández-Trincado 2020). In this paper, we analyze the DECaLS data and present evidence for a linear debris structure aligned with the proper motion and orbit of the Pal 13 cluster, extending ∼5° in either direction from the cluster center.

We perform our search using DECaLS data distributed in data release 8 (DR8) of the DESI Legacy Imaging Surveys (Dey et al. 2019). Source detection is performed with a point-source function (PSF) and spectral energy distribution (SED) matched-filter detection on the stacked images, with a 6σ detection limit. Morphological fitting and photometry are performed with the Tractor code initialized at the positions of detected sources (Lang et al. 2016), and we consider sources with TYPE = ''PSF'' to be stars. We perform our search on the g-, r-, and z-band data, and we require that sources are detected in all three bands. In addition, we require that sources pass quality cuts: (ANYMASK_G==0) & (ANYMASK_R==0) & (ANYMASK_Z==0) and (FRACFLUX_G < 0.05) & (FRACFLUX_R < 0.05) & (FRACFLUX_Z < 0.05). We correct the measured fluxes for interstellar extinction using the provided extinction values (which make use of Schlafly & Finkbeiner 2011).

We perform an unweighted matched-filter search to this data in color–magnitude (g versus g − r) space (e.g., Rockosi et al. 2002), following the procedure outlined in Section 3.1 of Shipp et al. (2018). We select stars around a synthetic isochrone derived from the stellar evolution models provided by the Dartmouth Stellar Evolution Database (Dotter et al. 2008). These isochrones are accessed programmatically through ugali, a codebase developed for the discovery and analysis of resolved stellar systems in DES (Bechtol et al. 2015; Drlica-Wagner et al. 2015).⁷ For the Pal 13 search, we select stars consistent with an isochrone of age 13.5 Gyr and metallicity Z = 0.0001 ($[\mathrm{Fe}/{\rm{H}}]=-2.2$). The selection is defined as in Equation (4) of Shipp et al. (2018), with the selection width parameters determined empirically by comparison to old, metal-poor stellar populations in the data, such as globular clusters and dwarf galaxies (our initial search is not specifically tailored to Pal 13). We limit our isochrone selection to an absolute magnitude of M_g > 3.5 to select stars lying along the higher signal-to-noise main sequence. We limit the observed magnitude range of our matched filter to g < 23 due to the decreased uniformity of the data at fainter magnitudes.

We scan our isochrone filter in distance modulus from 15 < m − M < 19 (10 < D < 63 kpc) in steps of 0.1 mag. At each distance modulus step, we fit a fifth-order polynomial to the filtered stellar density as a function of ${\alpha }_{2000},{\delta }_{2000}$ and smooth the data by applying a Gaussian filter with a smoothing kernel of 025. We subtract the background polynomial model from this smoothed data and visually inspect the residual stellar density maps as shown in Figure 1. We find several known features in the stellar density, including the Triangulum stream (Bonaca et al. 2012) and the Sagittarius stream (Majewski et al. 2003). The most prominent unknown linear feature is coincident with the globular cluster Pal 13, centered at (${\alpha }_{2000}$, ${\delta }_{2000}$ = 3467, 128).

In order to further examine the linear feature around Pal 13, we first determine optimized isochrone filter parameters (tuned to the stellar population of Pal 13) by eye via comparison to the Pal 13 color–magnitude diagram in the left-hand panel of Figure 2. We do not attempt to fit physically significant isochrone parameters to the cluster, and instead empirically determine a filter to most effectively select the Pal 13 signal (outlined in cyan in all panels of Figure 2). We choose to broaden our filter in the blue direction in order to select stars at the blue edge of the main sequence, which typically provide the highest signal-to-noise in matched-filter searches. We adopt isochrone parameters: $m-M=16.8$, $\tau$ = 13.5 Gyr and Z = 0.0001, where $\tau$ refers to the age of the stellar population. These parameters were selected based on the Pal 13 cluster color–magnitude diagram. However, we note that these parameters also appear to maximize the signal-to-noise of the tidal tails. Due to the width of the isochrone filter, the tidal tails are visible at a wide range of isochrone parameters, but we do find that significant adjustments to the age or metallicity noticeably reduce the visibility of the tidal tails in Figure 1.

Zoom In Zoom Out Reset image size

The selection region is also defined by a magnitude broadening (Δμ = 0.5), an asymmetric color broadening, C_1,2 = (−0.02, 1.0), and a multiplicative factor on the magnitude uncertainty (E = 2), which describes the error-dependent spread in color. All parameters are described in greater detail in Section 3.1 of Shipp et al. (2018).

The linear features detected by our isochrone filter extend ∼5° in each direction from the center of Pal 13. We estimate the end points to be (${\alpha }_{2000}$, ${\delta }_{2000}$) = (−9.8, 18.2)° and (−15.7, 8.9)° and find that the track of the tails is well matched by the great circle connecting these end points. We define a new coordinate system such that the end points lie along the x-axis and the cluster is positioned at (ϕ₁, ϕ₂) = (0, 0)°. The coordinate system is defined by the rotation matrix

$$\begin{eqnarray}R=\left[\begin{array}{ccc}0.94906836 & -0.22453560 & 0.22102719\\ -0.06325861 & 0.55143610 & 0.83181523\\ -0.30865450 & -0.80343138 & 0.50914675\end{array}\right],\end{eqnarray} \tag{ 1 }$$

which transforms from International Celestial Reference System (ICRS) to stream-aligned coordinates. Throughout the text, we consider the on-stream region to be between the end points (−45 < ϕ₁ < 64) and −025 < ϕ₂ < 025, excluding the region within 04 of the globular cluster. We define the radius of the cluster selection to be 04. This area was defined by determining the radius in the direction perpendicular to the tidal tails beyond which there is no visible contribution by the cluster signal in the color–magnitude diagram.

The region around Pal 13 is complicated by a linear feature in the interstellar dust adjacent to the globular cluster (right panel of Figure 1). Imperfect reddening corrections can lead to artificial color shifts that could increase or decrease the number of stars passing our isochrone filter. Indeed, we see a deficit of stars to the northwest of Pal 13 coinciding with this dust feature. A more subtle concern is that the underdensity of stars associated with the dust feature could conspire with the rapidly falling foreground stellar density to manifest as an apparent linear overdensity of stars directly adjacent to the dust lane.⁸

To distinguish a spurious density variation from the presence of a distant, metal-poor stellar population associated with Pal 13, we examine the significance of the tidal tails as a function of distance and compare to neighboring, off-stream regions. If the tails are an artifact of the foreground population, they would have the same distance dependence as the neighboring regions, while tidal features associated with Pal 13 would be more prominent at larger distances. In the top two panels of Figure 3, we show the normalized number of stars along the tidal tails, and in three adjacent equal-area regions above and below the stream. We count the number of stars passing the isochrone selection at a given distance modulus in each of the spatial regions. The counts in the off-stream regions are then normalized to the same maximum value as the on-stream region, in order to highlight the difference in shape between the curves. The number of stars falls off more slowly with distance modulus in the on-stream region than in any of the off-stream regions, including those overlapping with the dust feature visible in Figure 1. In the lower panel of Figure 3, we show the ratio of the normalized number of stars in the on-stream region and a region offset by +2°. We produce 20 realizations of these on-stream and off-stream regions, by selecting regions at the same Galactic latitude, but shifted in Galactic longitude within the DECaLS footprint. We exclude regions that overlap with the Sagittarius stream. In the lower panel of Figure 3, the solid black line shows the ratio between the on-stream region and the off-stream region at the location of Pal 13. The colored lines represent the ratios of counts in each of the other regions. We consider these ratios to follow a distribution corresponding to the ratio of two Poisson distributions. Given this assumption, the dashed black line represents the mean of the distribution, and the shaded region represents the ±1σ range. At the distance modulus of Pal 13, the ratio of counts in the Pal 13 on- and off-stream regions is a 4.2σ outlier. We determine that the candidate tidal tails have a different distance dependence than the foreground stars, and that this difference is significantly larger than that of other regions at the same Galactic latitude. We therefore conclude that the tidal tails are more likely associated with the distant stellar population of Pal 13 than with the foreground stellar population.

Zoom In Zoom Out Reset image size

In addition, we note the presence of a second feature extending from (${\alpha }_{2000}$, ${\delta }_{2000}$) = (−62, 199) to (−40, 260), which we refer to as the north–east arc (dotted red line in Figure 1). The feature is disconnected from the observed Pal 13 tidal tails; however, the orientation of the two features bears some resemblance to known broken streams, such as the ATLAS and Aliqa Uma stream (Li et al. 2019). Therefore, although we find no concrete evidence of an association between this north–east arc and Pal 13, this feature may be worth investigating further.

3.1. Color–Magnitude Diagram

3.2. RR Lyrae Stars

We search for RR Lyrae stars (RRLs) associated with the cluster and tidal tails using the Pan-STARRS1 (PS1) RRL catalog (Sesar et al. 2017). Following the same procedure as in Price-Whelan et al. (2019), we perform the search using the 61,795 RRLs selected by Sesar et al. (2017) as bona fide based on their classification scores, plus the 6459 non-bona fide PS1 stars also identified as RRLs in the Gaia Specific Objects Study catalog. We adopt photometric distances for the RRLs reported by Sesar et al. (2017), which are based on the PS1 (dust-corrected) i-band magnitudes. We also correct for the systematic offset in the distance of RRc stars present in the data published in Sesar et al. (2017) by recomputing their distances using their fundamentalized periods (see Price-Whelan et al. 2019 for full details). The resulting distances have a precision of ∼3% (Sesar et al. 2017). We retrieve astrometric information for these RRLs by cross-matching to the Gaia DR2 catalog with a 1'' tolerance.

Pal 13 is known to have four type-ab RRLs (V1–V4; Rosino 1957), which were recently confirmed as members based on Gaia proper motions (Yepez et al. 2019). Stars V1, V3, and V4 are present in the PS1 RRL catalog, but V2 is missing from both the PS1 and Gaia DR2 RRL catalogs (likely due to crowding, as it is within 3'' of the cluster center). Although this star was not identified as an RRL in Gaia DR2, it is present in the main catalog table gaia source where its proper motion is available, which allowed Yepez et al. (2019) to confirm its membership to the cluster. Since it is missing from the PS1 catalog, we will not include it in our distance analysis to avoid its G- or V-band inferred distance introducing potential systematics due to differences in the period–luminosity–metallicity (PLZ) calibrations.

We select RRLs with distances (17–35 kpc), positions ($| {\phi }_{2}| \lt 1^\circ$), and a prograde proper motions (corrected for the solar reflex motion; ${\mu }_{{\phi }_{1}}\gt 0$) consistent with Pal 13. This selection yields a total of eight type-ab RRLs: the three stars in the cluster plus another five in the tails. Of these five, two are likely contaminants as they have proper motions $({\mu }_{{\phi }_{1}},{\mu }_{{\phi }_{2}})$ close to zero, consistent with the background distribution. The other three RRLs have ${\mu }_{{\phi }_{1}}$ similar to the cluster RRLs. Our search does not reveal any significant overdensity of RRLs near the north–east arc.

The six likely members (three in the cluster, three in the tails) are shown as large markers in the upper panel of Figure 4. The arrows represent the proper motion of each star, and the color shows the distance offset in kiloparseconds from the cluster at $({\phi }_{1},{\phi }_{2})=(0,0)^\circ$. The Gaia DR2 source_ids for the RRLs are 2814894112367752192, 2814893910504601600, 2811888666052959232 (cluster) and 2712246494232725120, 2815432636842520960, 2705495870795198464 (tails).

Zoom In Zoom Out Reset image size

Using the three RRLs in the cluster, we compute a mean distance of 23.6 ± 0.2 kpc. We have not corrected for metallicity since the expected systematic error will be negligible (<0.3%), as the cluster almost matches the metallicity assumed by Sesar et al. (2017) in the PLZ relation to compute RR Lyrae distances. Our distance estimate is in excellent agreement with the previous estimate of 23.67 ± 0.57 kpc based on RRLs (Yepez et al. 2019) and is slightly shorter than previous isochrone distances ${24.3}_{-1.1}^{+1.2}$ kpc (Côté et al. 2002). Although well within the random uncertainties, this offset is accounted for by the differences in the assumed metallicity (−0.4 dex lower) and extinction law coefficients (they used the coefficients from Schlegel et al. 1998 available at the time; our distances are based on the new ones from Schlafly & Finkbeiner 2011).

We use the number of RRLs N_RR found to provide a rough estimation of the cluster's initial total luminosity, following the procedure described in Mateu et al. (2018), which uses the $\mathrm{log}{N}_{{RR}}\mbox{--}{M}_{V}$ relation observed for globular clusters and dwarf galaxies. The observed number of 7(3) type-ab RRLs in total (tails) yields an estimate of ${L}_{V}={5.1}_{-3.4}^{+9.7}\times {10}^{3}\,{L}_{\odot }$ (${1.4}_{-0.9}^{+2.6}\times {10}^{3}\,{L}_{\odot }$), slightly higher than the cluster luminosity estimated by Bradford et al. (2011; ${L}_{V}={1.1}_{-0.3}^{+0.5}\times {10}^{3}\,{L}_{\odot }$), but consistent within the uncertainties. Assuming a stellar mass-to-light ratio of M/L = 2.4 estimated for the cluster by Bradford et al. (2011), we estimate the stellar mass to be ${1.2}_{-0.8}^{+2.3}\times {10}^{4}\,{M}_{\odot }$ in total and ${3.3}_{-2.2}^{+6.3}\times {10}^{3}\,{M}_{\odot }$ for the tails. Our findings for the total luminosity and stellar mass therefore support the claim made by Siegel et al. (2001) and Bradford et al. (2011) that the abnormally large blue straggler population observed in the cluster suggests a higher initial mass and, therefore, luminosity.

3.3. Stream Model

The cluster Pal 13 is thought to be on a highly eccentric orbit (e.g., Siegel et al. 2001; Vasiliev 2019), but satellites that disrupt on very radial orbits are expected to form more diffuse tidal debris structures rather than linear streams (e.g., Helmi & White 1999). Is it then surprising that the tidal tails of Pal 13 appear stream-like? We use the sky track of the tidal tails determined here and the kinematics of the cluster to fit for the Galactic orbit of the stream, then use the best-fit orbit to generate simulated tidal debris.

We compile measurements of the Pal 13 mean sky position, proper motion, and line-of-sight velocity—(${\alpha }_{2000}$, ${\delta }_{2000}$) = (346685, 12772) (Vasiliev 2019), $({\mu }_{{\alpha }^{* }},{\mu }_{\delta })=(1.615,0.142)\,\mathrm{mas}\,{\mathrm{yr}}^{-1}$ (Vasiliev 2019), and v_los = 25.9 $\mathrm{km}\,{{\rm{s}}}^{-1}$ (Baumgardt et al. 2019)—and their associated uncertainties. We adopt a distance to Pal 13 based on the three cluster RRLs with PS1 measurements (Section 3.2), d = 23.6 ± 0.2 kpc. We adopt recently compiled Galactocentric solar position and velocity measurements (Drimmel & Poggio 2018) to transform between heliocentric and Galactocentric quantities. We use a standard, three-component mass model (implemented in gala; Price-Whelan 2017) to represent the Milky Way, with all disk and bulge parameters fixed to fiducial values (Bovy 2015; Price-Whelan 2017), but we fit for the mass and scale radius of the Navarro, Frenk, and White (NFW) halo component, and treat the cluster distance, proper motion, and radial velocity as free parameters. We include information about the sky track of the stream in our likelihood by adding control points spaced by 1° along ϕ₂ = 0 between ϕ₁ ∈ [−5, 5]°; we evaluate the orbit fit at each of these control points, assuming a constant Gaussian width of the stream as given in Table 1. We note that the spacing of these stream–track control points is arbitrary, but we did not find any major qualitative changes to the orbital fits after varying the number of control points by a factor of two in either direction. We use BFGS optimization (Jones et al. 2001) to maximize the likelihood of the orbit given the data described above.

Table 1.  Measured Properties of Pal 13 Tidal Tails

Parameter	Value
End point (${\alpha }_{2000}$, ${\delta }_{2000}$)	(−9.8, 18.2)°
End point (${\alpha }_{2000}$, ${\delta }_{2000}$)	(−15.7, 8.9)°
Length	109
Width	025

Download table as:  ASCIITypeset image

Our best-fit orbit for Pal 13 has a pericenter of ∼9 kpc and an apocenter of ∼69 kpc and suggests that the cluster passed through pericenter ≈75 Myr ago. We use this best-fit orbit to run an approximate N-body simulation of the stream formation (we do not resolve the internal dynamics or disruption process, but include a mass model for the progenitor cluster and generate tidally stripped stars following Fardal et al. 2015). Figure 4 (top) shows the sky positions (in stream-aligned coordinates) of star particles from this simulation, colored by relative distance to the cluster. The candidate RRL stars associated with the Pal 13 tails—identified from proper motions alone—are generally consistent with the distance trend and on-sky distribution of star particles in the stream model. At most earlier simulation snapshots (i.e., at different orbital phases), the tidal debris tends to be much more diffuse, so the coherence of the tidal tails may be another indicator that the cluster is near pericenter. This also highlights the importance of present-day orbital phase in determining whether tidal tails are observable.

We have applied an isochrone matched-filter technique to stars in DECaLS to detect evidence for tidal tails coincident with the Pal 13 globular cluster. This observation complements previous studies, which have shown evidence for tidal disruption of Pal 13, including the recent discovery of an extended low-density halo beyond the Jacobi radius (Piatti & Fernández-Trincado 2020).

The detected tidal tails extend ∼5° in either direction from the cluster and are well aligned with the proper motion of Pal 13 (Vasiliev 2019). We identify a color–magnitude signal in the on-stream region that is consistent with the stellar population of the Pal 13 cluster. In addition, we find three new RRL stars likely associated with the tidal tails, along with four RRLs known to be associated with the cluster. We generate a model of the tails using an orbit fit to the cluster kinematics and find that the model is consistent with the orientation and spatial distribution of the observed tidal tails, and the distance gradient of the RRLs.

The detection of the Pal 13 tidal tails is complicated by nearby structures in the interstellar dust maps, highlighting the importance of deep, precise photometry and the need for additional phase-space information to confirm and model stellar streams. Future observations, such as radial velocity and metallicity measurements by ${S}^{5}$ (Li et al. 2019), and deeper photometric and proper-motion measurements by the Rubin Observatory (Ivezić et al. 2019), will allow for higher-precision characterization and modeling of the tidal features of Pal 13.

If confirmed, Pal 13 will be one of only a handful of thin, extended stellar streams that has been confidently associated with a bound progenitor, and, conversely, one of the few globular clusters with detectable long tidal tails. Future observations by Gaia, the Rubin Observatory, and WFIRST will allow for the discovery of many more such systems in the coming years. Increasing the population of streams with known progenitors will provide crucial insight into the tidal disruption of globular clusters, the formation of the stellar halo, and the gravitational field of our Galaxy.

This paper is based upon work that is supported by the Visiting Scholars Award Program of the Universities Research Association. N.S. thanks the LSSTC Data Science Fellowship Program; her time as a Fellow has benefited this work. C.M. acknowledges support from the DGAPA/UNAM PAPIIT program grant IG100319.

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.

Software: astropy (Astropy Collaboration et al. 2018), gala (Price-Whelan 2017), matplotlib (Hunter 2007), numpy (Walt et al. 2011), scipy (Jones et al. 2001).