Candidate Members of the VMP/EMP Disk System of the Galaxy from the SkyMapper and SAGES Surveys

Photometric stellar surveys now cover a large fraction of the sky, probe to fainter magnitudes than large-scale spectroscopic surveys, and are relatively free from the target-selection biases often associated with such studies. Photometric-metallicity estimates that include narrow/medium-band filters can achieve comparable accuracy and precision to existing low-resolution spectroscopic surveys such as SDSS/SEGUE and LAMOST. Here we report on an effort to identify likely members of the Galactic disk system among the very metal-poor (VMP; [Fe/H] $\leq$ --2) and extremely metal-poor (EMP; [Fe/H] $\leq$ --3) stars. Our analysis is based on an initial sample of $\sim11.5$ million stars with full space motions selected from the SkyMapper Southern Survey (SMSS) and Stellar Abundance and Galactic Evolution Survey (SAGES). After applying a number of quality cuts to obtain the best available metallicity and dynamical estimates, we analyze a total of $\sim$5.86 million stars in the combined SMSS/SAGES sample. We employ two techniques that, depending on the method, identify between 876 and 1,476 VMP stars (6.9%-11.7% of all VMP stars) and between 40 and 59 EMP stars (12.4%-18.3% of all EMP stars) that appear to be members of the Galactic disk system on highly prograde orbits (v$_{\phi}>150$ km/s). The total number of candidate VMP/EMP disk-like stars is 1,496, the majority of which have low orbital eccentricities, ecc $\le 0.4$; many have ecc $\le 0.2$. The large fractions of VMP/EMP stars associated with the Milky Way disk system strongly suggest the presence of an early forming ``primordial"disk.

In a series of recent papers, An & Beers (2020, 2021a,b) and An et al. (2023) have constructed "blueprints" of the MW's stellar populations from analyses of the orbital rotation (inferred from proper motions and distance estimates alone) as a function of carefully calibrated photometricmetallicity estimates for stars with available broadband ugriz from SDSS/SEGUE and other surveys.This approach has proven quite powerful.Among other results, these authors not only identified previously known substructures and confirmed the presence of the inner-and outer-halo populations but also demonstrated that the metal weak-thick disk (MWTD; Norris et al. 1985;Carollo et al. 2007Carollo et al. , 2010;;Beers et al. 2014) is a separable population with lower metallicity and rotation that lags the canonical thick disk, as shown in Carollo et al. (2019).In addition, they identified a continuous sequence of stars in the rotational velocity vs. metallicity space that may be associated with a starburst event when the earlier disk system encountered Gaia-Sausage-Enceladus (GSE; Belokurov et al. 2018;Haywood et al. 2018;Helmi et al. 2018).Evidence for this starburst event is also reported in Lee et al. (2023).
Whether surveys to identify likely metal-poor (MP) stars are performed with fiber-fed spectrographs such as SDSS or LAMOST or broadband photometric efforts such as SEGUE or Pan-STARRS (Chambers et al. 2016), it is challenging to avoid target selection biases that can confound the relative contributions of stars with different metallicity to the recognized Galactic components.In addition, the first step in surveys dedicated to finding low-metallicity stars is often to limit the regions of the sky under consideration to higher Galactic latitude (e.g., |b| > 30 • ), precluding identification of substantial numbers of VMP/EMP stars in the disk system of the MW.
Over the past few years, photometric surveys based on combinations of narrow-band and medium-to broad-band filters have been (or are being) executed (e.g., SkyMapper; Keller et al. 2007, the Pristine Survey;Starkenburg et al. 2017, Stellar Abundance and Galactic Evolution Survey (SAGES); Zheng et al. 2018, J-PLUS;Cenarro et al. 2019, and S-PLUS;Mendes de Oliveira et al. 2019).Typically, such surveys do not avoid regions of the MW at lower Galactic latitudes, other than those limited by very high interstellar extinction and reddening or crowding.As a result, VMP/EMP stars in the MW's disk system have been increasingly discovered, though their numbers are still relatively small.
The SkyMapper Southern Survey Data Release 2 (SMSS DR2; Onken et al. 2019) was carefully recalibrated by Huang et al. (2021), and used by Huang et al. (2022) to derive stellar parameters, luminosity classifications, and metallicity estimates for over 24 million stars in the Southern Hemisphere.These authors derived effective temperatures (T eff ) by adopting metallicity-dependent T eff -color relations constructed from Gaia (G BP − G RP ) 0 , LAMOST T eff , and [Fe/H].The effective temperature scale of LAMOST has been shown to agree with that of direct measurements (Huang et al. 2015).They adopted Bayesian distance estimates (Bailer-Jones et al. 2021) and ages from the PAdova and tRieste Stellar Evolutionary Code (PARSEC; Bressan et al. 2012;Marigo et al. 2017) isochrones.In addition, an empirical metallicitydependent stellar-locus method (Yuan et al. 2015) was used to estimate the photometric metallicity, with combinations of the SMSS narrow/medium uand v-band filter magnitudes, the G BP magnitude from the Gaia ultra wide band prism spectra, and a maximum likelihood approach (Huang et al. 2022).
The recently completed SAGES (Fan et al. 2023), which employs similar, but not identical, filters to SMSS, has been employed by Huang et al. (2023) to obtain stellar parameters, luminosity classifications, and metallicity estimates for nearly 26 million stars in the Northern Hemisphere.
Here we identify 1,496 VMP and 61 EMP candidate stars with disk-like orbits populating the rapidly rotating disk system of the MW (v ϕ > 150 km s −1 ), selected from a subset of roughly 11.5 million stars from the SMSS and SAGES photometric surveys with available radial velocities (RVs), proper motions, and other astrometric data from which full space motions are derived.We approximately separate stars with disk-like orbits from stars with halo-like orbits by two criteria that have been commonly used in the literature (Haywood et al. 2018;Di Matteo et al. 2020;Mardini et al. 2022a;Bellazzini et al. 2024), and then consider their relative fractions at low metallicities.
This paper is organized as follows: In Section 2, we describe the data sets we employ and the choices made for the adopted metallicity estimates, as well as the derivation of dynamical parameters.In this section we also describe two methods that have been commonly used to identify stars with potential disk-like orbits.In Section 3, we present maps of the orbital rotational velocities of the stars as a function of [Fe/H], where potential VMP/EMP candidates with disk-like orbits can already be seen, and compare their relative fractions as a function of [Fe/H].In Section 4, we present a discussion, along with conjectures on the origins of VMP/EMP disk-like stars based on interpretations from numerical sim-ulations.We conclude with a summary and future prospects in Section 5. et al. (2022) derived stellar parameters, including metallicity estimates, for more than 19 million dwarfs and 5 million giants over essentially the entire Southern Hemisphere from SMSS DR2, including about 731,000 VMP and 27,000 EMP stars.If we restrict their sample to stars with available RVs from Gaia DR3 (Gaia Collaboration et al. 2023) and other sources, the number of stars is about 7.4 million, including roughly 56,000 VMP and 2,300 EMP stars.

Huang
SAGES observed slightly less than half of the Northern Hemisphere.Notably, SAGES did not cover a large fraction of the north galactic pole (NGP), while SMSS covered the entire south galactic pole (SGP).Another crucial difference between SAGES and SMSS is that the central wavelength of the SAGES v-band filter is shifted redward relative to the SMSS v-band filter by about 110 Å, so it fully includes the region of the Ca II H & K lines (Zheng et al. 2018).Huang et al. (2023) used a similar approach to Huang et al. (2022), and obtained effective temperatures, luminosity classifications based on surface gravity, and metallicity estimates for over 26 million stars, including some 874,000 VMP and 13,000 EMP stars from SAGES DR1 (Fan et al. 2023).About 4.1 million stars in this catalog have available RVs, including roughly 41,000 VMP and 1,900 EMP stars.
For this study, we begin with a sample of about 7.4 million stars from SMSS and 4.1 million stars from SAGES with available RVs, proper motions, and distance estimates, as provided in the catalogs from Huang et al. (2022) and Huang et al. (2023), respectively.After combining these data sets, binary stars photometrically classified by Huang et al. (2022Huang et al. ( , 2023) ) and cool dwarfs (T eff < 4,500 K) have been removed.We have additionally applied a more restrictive cut on the bp rp excess factor, < 0.12 × (BP-RP) 0 + 1.13/1.14cuts for dwarfs/giants respectively, as in Xu et al. (2022), and on the renormalized unit weight error (RUWE) < 1.1, in order to exclude possible binary stars.Cuts based on an empirical isochrone, similar to the PARSEC (Bressan et al. 2012;Marigo et al. 2017) isochrone with [Fe/H] = −2 at age = 12 Gyr, were also applied to eliminate the significant contamination from metal-rich stars that could masquerade as VMP/EMP stars.These restrictions removed a total of about 3.3 million stars.
Moreover, we removed about 940 likely stellar globular cluster members based on the catalogs of Harris (2010) 1 and  Baumgardt & Vasiliev (2021).Finally, in order to diminish the effect of reddening on the derived metallicities (of par-1 https://physics.mcmaster.ca/∼ harris/mwgc.datticular importance for stars near the disk), we only included stars with E(B −V ) ≤ 0.3, excluding a total of about 17,000 stars.More discussion about the extinction cut is provided in the Appendix.
Metallicity estimates for the stars in our sample are based on calibrated u−G BP colors and v −G BP colors, a combination of the u/vbands from SMSS/SAGES and the ultra wideband Gaia G BP prism spectra (Huang et al. 2022(Huang et al. , 2023)).As has been noted previously, the colors involving the u band have a greater sensitivity to the presence of enhanced carbon in a star than those involving the v band.For this reason, and in order to provide the best available metallicities, we do not include stars for which only u band metallicity estimates are available, those that have a difference between the u-bandand v-band-based abundances greater than 0.5 dex, and stars with estimated metallicity errors greater than 0.5 dex.See the Appendix for a justification and full discussion of the cuts that are applied to the sample prior to assigning final adopted metallicities.Note that we refer to the photometricmetallicity estimates as [Fe/H] in this study, unless otherwise indicated.

Dynamical Parameters
Orbital parameters for the stars in our combined sample are determined using their 6D astrometric parameters (positions, RVs, proper motions, and distance estimates from Huang et al. 2022Huang et al. , 2023)), as well as their corresponding errors, as inputs to the Action-based GAlaxy Modelling Architecture2 (AGAMA) package (Vasiliev 2019), adopting the solar positions and peculiar motions described in Shank et al. (2022) 3 , and the gravitational potential MW2017 (McMillan 2017).
Similar to Shank et al. (2022), we input quantities through the orbital integration process in AGAMA to calculate the cylindrical velocities (v r , v ϕ , v z ), cylindrical actions (J r , J ϕ , J z ), orbital specific energy (E), r apo , r peri , eccentricity (ecc), Z max (the maximum orbital distances reached by stars from the Galactic plane), and R max (the maximum apocentric distance projected onto the Galactic plane), along with their associated errors 4 .Stars that are possibly unbound (E > 0 km 2 s −2 ) were identified and removed.This resulted in a total of 10.8 million from the initial 11.5 million stars that are suitable for our kinematic analysis.For our present purpose, we only included stars having derived errors less than 25 km s −1 in their orbital rotation velocities and relative errors ≤ 30% in Z max and R max , which removed about 165,000 stars from the combined sample.This produced a final sample of approximately 5.86 million stars, including 4.07 million SMSS and 1.79 million SAGES stars, which we refer to as the SMSS/SAGES sample hereafter.
Figure 1 shows the sky distribution in equatorial and Galactic coordinates for the final 5.86 million stars of the combined SMSS/SAGES sample.The gray filled circles indicate the stars with −2 < [Fe/H] ≤ +0.5, the light-blue filled circles are stars with −3 < [Fe/H] ≤ −2, and the black filled circles represent stars with −4 < [Fe/H] ≤ −3.The ∼ 143, 000 stars in common between the two surveys have differences in [Fe/H] from the v band with a median value of only 0.02 dex; we have adopted an average of these determinations for these stars.
Figure 2 shows the distribution of the errors in the photometric-metallicity estimates (δ [Fe/H]) for the combined SMSS/SAGES sample.The left, middle, and right panels provide the results for the subsamples of stars with [Fe/H] ≤ −1, ≤ −2, and ≤ −3, respectively.The legends in each panel indicate the median errors for all stars in the listed metallicity range, and the errors for stars classified as dwarfs and giants.As can be seen, the errors increase with decreasing metallicity, as expected, but still remain reasonably low (median errors on the order of 0.1 to 0.2 dex).Note that the external errors are somewhat larger, on the order of 0.25-0.35dex (see the discussion in the Appendix).

Separation of Disk and Halo Stars
Previous analyses of the nature of stellar orbits in the MW have used a variety of techniques to separate stars on disklike orbits from stars on halo-like orbits.Two simple approaches are described below.

Maximum Height of Orbits
This approach, employed by Beers et al. (2014), Sestito et al. (2020, 2021), Limberg et al. (2021), Mardini et al. (2022a), andBellazzini et al. (2024), identifies stars in the disk-like and halo-like dynamical populations by assigning stars with Z max ≤ 3 kpc to disk-like orbits and those with Z max > 3 kpc to halo-like orbits.Often, an additional criterion is adopted to identify the stars in the disk system by demanding that they be on highly prograde orbits.We follow a similar approach to that described below, with a further division of the stars on disk-like orbits into those with Z max ≤ 1 kpc, in an attempt to identify possible VMP/EMP thindisk stars.

"Wedges" in the Haywood Diagram
Following Haywood et al. (2018), we have also used plots of Z max vs. arctan (Z max /R max ), which redistributes our sample stars into discrete wedges, corresponding to different dynamical populations, a method also employed by Schuster et al. (2012), Di Matteo et al. (2020), Kim et al. (2021), and Koppelman et al. (2021).
Here R max is defined as the projection of r apo onto the Galactic plane, via the simple geometric relationship R max = r 2 apo − Z 2 max .Note that, for simplicity of notation, below we define an "inclination angle" (IA) to represent arctan (Z max /R max ).It should be kept in mind that Z max and R max are derived from the full ensemble of orbits traced by a given star, so their IA is representative of that complete set, not a single orbit or an average of the orbits.

RESULTS
In this section, we identify about 12,700 VMP/EMP stars over the full range of the rotational velocities of the final 5.86 million stars in the combined SMSS/SAGES sample.Among these metal-deficient stars, we closely examine the 2,150 rapidly rotating VMP/EMP stars, in order to classify them as on halo-like or disk-like orbits.

v ϕ vs. [Fe/H]
Figure 3 shows plots of stellar number density for our sample in the rotational velocity vs. photometricmetallicity plane.The top panel indicates the total combined SMSS/SAGES sample color-coded on a logarithmic scale.The rapidly rotating canonical disk system (comprising both the thin and thick disk) is most visible for [Fe/H] > −1.In addition, as reported in the series of papers by An & Beers (2020, 2021a,b), the MWTD, the Splashed Disk (SD), and a hint of the GSE substructure can be seen in the black dashed ellipses.However, the VMP/EMP stars in the rapidly rotating disk region (v ϕ > 150 km s −1 ; the average value of rotational velocity for the MWTD from Carollo et al. 2010), are less visible in the number density map than the other components.Thus, we represent the ∼2,150 highly prograde VMP/EMP candidates with white circles (out of a total number of about 12,700 VMP/EMP stars).The middle and bottom panels are the same v ϕ vs. [Fe/H] plane, but for dwarfs and giants, respectively.There are about 2,800 total VMP/EMP dwarfs and 9,900 total VMP/EMP giants.The subsets of these stars with v ϕ > 150 km s −1 are roughly 650 VMP/EMP dwarfs and 1,500 VMP/EMP giants, respectively.We point out that, at this stage, we have not separated the VMP/EMP stars with disk-like orbits from those with halo-like orbits.However, from inspection of the middle and bottom panels, it is clear that the distribution of rotational velocity for the VMP/EMP giants (which are expected to contain a greater fraction of halo-like orbits) in the bottom For the purpose of our analysis, we exclude stars identified as likely binaries by Huang et al. (2022Huang et al. ( , 2023)), cool dwarfs (Teff < 4,500 K), metal-rich stars masquerading as VMP/EMP stars (see text), and likely members of recognized globular clusters.For stars in common between the two surveys, we have used the average value of the photometric-metallicity estimates.panel stands in contrast to that of the VMP/EMP dwarfs, seen in the middle panel.

Based on the Zmax Criterion
We first consider MP (MP; [Fe/H] ≤ −1) stars in three regions of Z max : Z max > 3 kpc, Z max ≤ 3 kpc, and Z max ≤ 1 kpc .We assign the stars with Z max > 3 kpc to the halo populations, while those with Z max ≤ 3 kpc and Z max ≤ 1 kpc are candidate members of the MP thick-and thin-disk systems, respectively.The left panel of Figure 4   high eccentricity, for stars kinematically assigned to the halo population is clear, as is the presence of low-eccentricity stars among the stars assigned to the disk system.Note that at this point we have not applied any cuts on v ϕ , only on Z max , so we expect that the subsamples of stars with Z max ≤ 3 kpc and Z max ≤ 1 kpc have some level of contamination from halo stars.
Figure 5 shows the distribution of v ϕ for these three subsamples, but only for stars with [Fe/H] ≤ −2, [Fe/H] ≤ −2.5, and [Fe/H] ≤ −3, from the left to right panels, respectively.The vertical dashed line corresponds to a cut on v ϕ = 150 km s −1 , the average orbital rotation value for the MWTD from Carollo et al. (2010).We note that the adopted limit for the MWTD stars with the lowest v ϕ from Carollo et al. (as well as from An & Beers 2021b) is v ϕ ∼ 100 km s −1 .From inspection, there remains considerable contamination of the VMP/EMP stars with prograde disk-like orbits by stars with prograde halo-like orbits, even with the higher cut at v ϕ > 150 km s −1 (although it is substantially less for the EMP stars shown in the right panel), indicating that a more sophisticated separation methodology is desirable.
Figure 6 shows histograms of the eccentricity distribution for stars with [Fe/H] ≤ −2, ≤ −2.5, and ≤ −3, respectively.We now subdivide the stars into three regions: Z max > 3 kpc, 1 kpc < Z max ≤ 3 kpc, and Z max ≤ 1 kpc, in an attempt to better isolate stars with thick-disk orbits from those with thin-disk orbits.We note that these divisions are imperfect, in that we expect there to be contamination from halo stars at all Z max .Within the 1 kpc < Z max ≤ 3 kpc region there should be few thin-disk stars.Within the cut Z max ≤ 1 kpc there will also remain some contamination from thick-disk stars.
From inspection of the top row of panels in Figure 6, which includes stars on both retrograde and prograde orbits, the VMP/EMP stars do not exhibit prominent low-eccentricity orbits in all three ranges of [Fe/H].However, in the bottom row of panels, the introduction of the v ϕ > 150 km s −1 cut greatly increases the relative dominance of VMP/EMP stars with disk-like orbits, including for stars with ecc ≤ 0.4.We note that similar results for VMP stars have been found by Bellazzini et al. (2024), based on a sample of some 700,000 stars with photometric-metallicity estimates obtained with synthetic Strömgren photometry from Gaia DR3 by Bellazzini et al. (2023).

Based on the Haywood Criterion
There is also evidence for the existence of a VMP/EMP disk system from the Haywood diagram.If we now specialize to the highly prograde stars with orbital velocities v ϕ > 150 km s −1 (indicated by the blue shaded region in Figure 7), the relative dominance of the stars in the disk-like system for VMP/EMP stars becomes even clearer.  of panels, which includes stars on both retrograde and prograde orbits, the candidate VMP/EMP thin disk-like stars are broadly distributed over all eccentricities, while the candidate VMP/EMP thick disk-like stars exhibit similar patterns to halo-like stars in all three ranges of [Fe/H].However, in the bottom row of panels, the v ϕ > 150 km s −1 cut increases the relative dominance of VMP/EMP stars on disk-like orbits, including for stars with ecc ≤ 0.4.
In summary, a total of 1,496 candidate VMP/EMP stars with v ϕ > 150 km s −1 are identified.The total numbers of highly prograde disk-like candidates selected by the Z max and Haywood criteria are 876 and 1,476, respectively.There are 856 stars selected by both methods.

DISCUSSION
Separation based on the Z max criterion, although capable of identifying a relatively pure sample of stars with halo-like orbits, has considerable potential contamination of stars with disk-like orbits by halo-like stars.The Haywood criterion, based on the separation of stars in different dynamical populations, produces a purer sample of disk-like stars, especially when the dominance of stars with thick and thin disk-like or-bits over stars with halo-like orbits at ecc < 0.4 is considered in conjunction.
As we have demonstrated in this paper, the large numbers of stars now available with photometric-metallicity estimates from SMSS and SAGES have increased the numbers of candidate disk system VMP/EMP stars dramatically.
We refer to these stars as candidates for two reasons.First, although the photometric-metallicity estimates have a precision of ∼ 0.1 dex for [Fe/H] ≥ −1 and approximately 0.3 dex at [Fe/H] ∼ −3.5 to −4.0, comparable to those obtained from low to medium resolution (R = λ/∆λ ∼ 1, 800) with a signal-to-noise ratio greater than 20 (Luo et al. 2015;Yanny et al. 2009;Rockosi et al. 2022), they may be influenced by the presence of strong molecular carbon bands, in particular for the most metal-poor stars.We have taken steps to mitigate this behavior, as described in the Appendix, but they should be confirmed by follow-up spectroscopy.Secondly, the question remains whether at least some of the apparent disk system VMP/EMP stars represent members of an early forming in-situ primordial disk system, prior to additional stars being added from accreted dwarf satellites, or are possibly a very/extremely low metallicity tail of the long-recognized MWTD component of the MW.These alternatives may prove difficult to differentiate between based on kinematics alone,  as mergers with dwarf galaxies could readily perturb the orbits of stars that were born in a primordial thin or thick disk.
The best way to distinguish between these two possibilities may be to conduct a thorough study of their elemental abundances and look for differences as a function of declining metallicity.Feltzing & Feuillet (2023) have recently used elemental abundance information from APOGEE, in combination with kinematics, in order to identify the likely presence of an early disk structure in the inner disk of the MW including VMP stars (although they are limited by the lack of lower-metallicity stars in APOGEE, precluding verifica-tion that EMP stars are present as well).Detailed chemical abundances for our candidate VMP/EMP stars would clearly be useful.
Additional information should soon be available from the J-PLUS and S-PLUS photometric surveys, which can obtain estimates for C and Mg (as well as N and Ca, once ongoing calibrations are completed), in addition to [Fe/H], thanks to their narrow/medium-bandpass filters.More complete information will require high-resolution spectroscopic follow-up for at least a subset of the VMP/EMP candidates.Determination of more accurate age estimates than we have at present for candidate VMP/EMP stars on disk-like orbits may also prove illuminating.

Comparison with Simulations
Beyond the identification of the VMP/EMP disk system candidates, we can speculate on their origins by considering numerical simulations of an MW-like galaxies.We analyzed the data from a high-resolution cosmological zoom-in simulation of a Milky Way-like galaxy with a halo mass of 1.2 × 10 12 M ⊙ presented in Hirai et al. (2022).These authors defined the in-situ component as stars formed in the main halo of the central galaxy, whereas the accreted component was defined as stars coming from dwarf galaxy satellites.Data within the Galactocentric distance r GC between 3 and 20 kpc were considered; this region roughly corresponds to the observed region by SMSS and SAGES.
From this simulation, we found that 8% and 92% of stars with v ϕ > 150 km s −1 and [Fe/H] ≤ −2 are formed in the in-situ and accreted components, respectively.We also found that 96% of VMP/EMP stars with v ϕ > 150 km s −1 have ages > 10 Gyr.
Similar results have been shown in the analysis of IllustrisTNG50 simulations by Mardini et al. (2022a) and Carollo et al. (2023).Most recently, Sotillo-Ramos et al. ( 2023) considered a large sample of 138 MW analogs from the TNG50 cosmological simulations and found that, across all of these analogs, about 20% of the VMP/EMP stars have disk-like orbits, with some analogs reaching as high as 30%.Roughly half of their disk-like stars have average ages exceeding 12.5 Gyr, with 70% coming from accreted dwarf galaxies.Taken as a whole, the simulation results suggest that VMP/EMP stars with disk-like orbits comprise stars coming primarily from accreted dwarf galaxies and in-situ stars formed in an early primordial disk, or are associated with the MWTD.
Fractions of stars with prograde orbits can also inform the origin of VMP/EMP stars.Figure 10 compares the prograde fractions as a function of [Fe/H] for the combined SMSS/SAGES sample and the simulation results of Hirai et al. (2022).For [Fe/H] > −2.0, both our SMSS/SAGES sample and the simulation show an increasing trend toward higher metallicity, attributable to disk formation.On the other hand, the prograde fraction in our sample is roughly constant as a function of [Fe/H] for VMP/EMP stars.In the lowest-metallicity regime, the fraction rises to 0.68, reflecting the higher fraction of disk-like orbits among EMP stars (Figures 6 and 9).This tendency is not clearly seen in the simulation.However, note that our sample's prograde fraction for VMP/EMP stars is significantly larger than that of the simulation.It should be kept in mind that the Hirai et al. (2022) simulation is for a single realization, when in fact a variety of galaxy assembly histories are likely to have different outcomes, as demonstrated by other recent simulation studies (e.g., Santistevan et al. 2021).
Prograde orbit fractions higher than 0.5 for [Fe/H] < −2 suggest the accretion of satellites preferentially on prograde orbits or early disk formation at low metallicity.Carter et al. (2021) also reported a high prograde fraction, between 0.7 and 0.8, but with significantly larger error bars owing to the smaller sample they considered (see the top left panel of their Figure 3).They have also shown that the prograde fraction converges to 0.5 with a model assuming an isotropic distribution of orbits in the stellar halo (Rybizki et al. 2018).We confirm their results with a larger sample.Recently, Li et al. (2022b) have shown that 10 out of 12 MW stellar streams with an average [Fe/H] ≈ −2 are on prograde orbits.These enhanced prograde fractions mean that the MW's VMP/EMP stars tend to be formed in accreted components with prograde orbits or in an ancient disk.As discussed above, spectroscopic follow-up of our candidate disk-like VMP/EMP stars may help improve estimates of the relative fractions associated with these differing origins.
For convenience of future comparisons of our observations with those of others and with numerical simulations, Table 1 provides a summary of the numbers, fractions, and orbital characteristics of the SMSS/SAGES sample for different cuts on [Fe/H], Z max , and IA.Note that, except for the first line in each subsection of the table (indicated as "All"), the fractions refer to the total numbers of stars listed on the first line at the top of each column in the subsection (shown in bold).

SUMMARY
We have identified 1,496 candidate VMP/EMP disk system stars in the MW from a subset of the ∼ 50 million stars from SMSS and SAGES with available photometric-metallicity estimates, based on calibrated u−G BP colors and v−G BP colors, a combination of the u/v-bands from SMSS/SAGES and the ultra wide band Gaia G BP prism spectra (Huang et al. 2022(Huang et al. , 2023)).We then trimmed the combined sample, eliminating photometrically identified binaries, cool dwarfs, and likely members of globular clusters.We then obtain the subset of 7.19 million stars in the combined sample with available RVs, proper motions, and distance estimates.
After the determination of dynamical parameters, we remove likely unbound stars and excise stars with errors in their orbital rotation velocities v ϕ > 25 km s −1 and relative errors in Z max (maximum orbital distance from the Galactic plane) and in R max (maximum orbital apocentric distance projected on to the plane) > 30%, leaving a total sample of about 5.86 million stars.
We then apply two methods to separate stars with halo-like and disk-like orbits.The first approach considered stars with Z max > 3 kpc to have halo-like orbits and those with Z max ≤ 3 kpc to have disk-like orbits.Our analysis indicates that there exists a significant population of candidate VMP/EMP disk system stars, moving on rapid prograde orbits (v ϕ > 150 km s −1 ), increasing their relative populations with declining metallicity.We also split the stars with disk-like orbits into the regions 1 < Z max ≤ 3 kpc, and Z max ≤ 1 kpc, in an attempt to better isolate stars with thick-disk orbits from those with thin-disk orbits.Based on this criterion, we find that 28.7% of the VMP stars with Z max ≤ 3 kpc have highly prograde disk-like orbits (707 stars on thick-disk orbits, 169 stars on thin-disk orbits), while 44.9% of the EMP stars have highly prograde disk-like orbits (25 on thick-disk orbits, 15 on thin-disk orbits).These fractions increase further if one also takes the eccentricity of the orbits into account.
The second approach considered the stars populating wedges in the diagram of Z max vs. IA, which redistributes corresponding to different dynamical populations of stars with halo-like and disk-like orbits.Our analysis indicates that there exists a significant population of candidate VMP/EMP disk system stars moving on rapid prograde orbits (v ϕ > 150 km s −1 ), increasing their relative populations with declining metallicity.Based on the Haywood criterion, we find that 21.9% of the VMP stars have highly prograde disk-like orbits (817 stars on thick-disk orbits, 659 stars on thin-disk orbits), while 35.3% of the EMP stars have highly prograde disk-like orbits (28 stars on thick-disk orbits, 31 stars on thin-disk orbits).These fractions increase further if one also takes the eccentricity of the orbits into account.
In the near future, the astrophysical properties and origin of these stars will be examined further with data from the large-scale Javalambre/Southern Photometric Local Universe Surveys (J/S-PLUS).These surveys include additional narrow/medium-band filters that allow for photometric estimates of C, N, Mg, and Ca abundances, once ongoing calibration efforts are completed.Of importance, it will then be possible to greatly reduce the influence of carbon on the metallicity estimates, which affect our current SMSS/SAGES sample, as [C/Fe] can be estimated separately from the [Fe/H].The accuracy and precision of the derived metallicities will be improved as well.
In order to confirm the metallicities and elemental abundance estimates (such as the α-elements or carbonicity, [C/Fe]) for the VMP/EMP stars with disk-like orbits, we require medium-resolution spectroscopic follow-up for the catalog of ∼ 1,500 VMP/EMP stars in our sample.Highresolution spectroscopic follow-up of at least the most interesting subset of these would also be useful.The full catalog of VMP/EMP stars with disk-like orbits is listed in the Appendix and will be made available in the online material.Determinations of age estimates for our candidate VMP/EMP stars would also help to place constraints on their origins.Nevertheless, our present finding that large fractions of VMP/EMP stars are kinematically associated with the rapidly rotating MW disk system (in particular those at low eccentricity) strongly suggests the presence of an early forming "primordial" disk.
The authors express thanks to Deokkeun An for his thoughts and comments on an early version of this paper.We are also grateful to Evan Kirby and Borja Anguiano for their input as well.The authors thank an anonymous referee for their helpful comments that greatly improved this paper.The calibration-corrected G magnitude uncertainty by Huang et al. (2022Huang et al. ( , 2023) ) for the Gaia DR3 [err G C] −

E
The orbital energy as given by AGAMA km 2 s −2 e E The orbital energy uncertainty as given by AGAMA km 2 s −2 Jr,Jphi,Jz The cylindrical actions as given by AGAMA kpc km s −1 e Jr,e Jphi,e Jz The cylindrical actions uncertainty as given by AGAMA kpc km s −1 rperi The Galactic pericentric distance as given by AGAMA kpc e rperi The Galactic pericentric distance uncertainty as given by AGAMA kpc rapo The Galactic apocentric distance as given by AGAMA kpc e rapo The Galactic apocentric distance uncertainty as given by Monte Carlo sampling through AGAMA kpc

Zmax
The maximum height above the Galactic plane as given by AGAMA kpc e Zmax The maximum height uncertainty above the Galactic plane as given by Monte Carlo sampling through AGAMA kpc relerr-Zmax The relative uncertainty of the maximum height above the Galactic plane −

Rmax
The projection of the Galactic apocentric distance onto the Galactic plane as given by r 2 apo − Z 2 max kpc relerr-Rmax The relative uncertainty of the projection of the Galactic apocentric distance onto the Galactic plane.

−
The uncertainty is as given by Monte Carlo sampling through AGAMA −

IA
The inclination angle defined as the arctangent ratio of (Zmax / Rmax) rad ecc The eccentricity as given by (rapo -rperi) / (rapo + rperi) through AGAMA − criterion Flag with "Zmax", "Haywood", or "Both" to indicate the criterion used for the separation of stars on halo-like orbits from − those on thick-disk like and thin-disk like orbits −

Figure 1 .
Figure 1.Mollweide projection of the positions for 5.86 million stars selected from the Southern Hemisphere (SMSS) and Northern Hemisphere (SAGES) photometric surveys, in equatorial (left panel) and Galactic (right panel) coordinates.The metallicity of this final sample is based on calibrated v − GBP colors and a combination of the u/v bands from SMSS/SAGES (Huang et al. 2022, 2023).See the Appendix for a discussion of the cuts that are applied to the sample prior to assigning the final metallicities.The gray filled circles indicate the stars with derived metallicities in the range −2 < [Fe/H] ≤ +0.5, the light-blue filled circles are stars with −3 < [Fe/H] ≤ −2, and the black filled circles represent stars with −4 < [Fe/H] ≤ −3.The stars shown all have available RVs and astrometric information.For the purpose of our analysis, we exclude stars identified as likely binaries byHuang et al. (2022Huang et al. ( , 2023)), cool dwarfs (Teff < 4,500 K), metal-rich stars masquerading as VMP/EMP stars (see text), and likely members of recognized globular clusters.For stars in common between the two surveys, we have used the average value of the photometric-metallicity estimates.

Figure 2 .
Figure 2. Histogram of the errors in photometric-metallicity estimates (δ [Fe/H]) for the final SMSS/SAGES sample.From left to right, the panels correspond to MP ([Fe/H] ≤ −1), VMP ( [Fe/H] ≤ −2), and EMP ([Fe/H] ≤ −3) sub-samples, respectively.The black, blue, and red histograms represent the total, dwarf, and giant stars in each metallicity range, respectively.The number of stars and the median values of δ [Fe/H] are indicated in the legend of each panel.
shows the cumulative numbers of each population.At [Fe/H] ≤ −2, approximately 9,600 halo stars with Z max > 3 kpc (black line) were found, along with about 3,000 stars with Z max ≤ 3 kpc (red line) and 430 stars with Z max ≤ 1 kpc (blue line).The middle panel shows the cumulative distribution function (CDF) for each population, normalized by the number of MP stars.Roughly 11.5% of the MP stars assigned to the halo system are VMP stars, and 0.24% are EMP stars; about 3.5% of the MP stars assigned to the disk system are VMP stars, and 0.08% are EMP stars.It is interesting to note that the CDFs of the stars with Z max ≤ 3 kpc and Z max ≤ 1 kpc are almost identical.The right panel of Figure 4 shows the distribution of orbital eccentricity for [Fe/H] ≤ −1 of these subsamples split on Z max .The broad distribution of eccentricity, peaking at

Figure 3 .
Figure 3. Top panel: rotational velocity distribution (v ϕ ) of the SMSS/SAGES sample as a function of photometric metallicity ([Fe/H]).The number density is color-coded on a logarithmic scale.The MWTD, the SD, and the GSE substructure are marked with black dashed ellipses.The highly prograde VMP/EMP candidates are shown by white circles.Middle panel: same as the top panel, but for dwarfs.Bottom panel: same as the top panel, but for giants.The number of VMP/EMP stars with v ϕ > 150 km s −1 is about 2,150, including 650 dwarfs and 1,500 giants.
Figure 7 shows the distribution of the arctangent of the Z max /R max values (defined here as the inclination angle, IA) for VMP/EMP stars, following Haywood et al. (2018) and Di Matteo et al. (2020).The panels show this distribution for stars with the total numbers of VMP/EMP stars (black line), for stars with prograde orbits (purple line), and for stars with retrograde orbits (orange line), respectively, for stars in the regions with [Fe/H] ≤ −2, ≤ −2.5, and ≤ −3.The vertical dashed line and

Figure 4 .
Figure 4. Left panel: cumulative number distribution of the MP ([Fe/H] ≤ −1) stars as a function of [Fe/H], for stars with Zmax > 3 kpc (black line), Zmax ≤ 3 kpc (red line), and Zmax ≤ 1 kpc (blue line).Middle panel: cumulative distribution functions of [Fe/H] for these subsamples.Each population is normalized on the basis of the number of stars at [Fe/H] = −1.Right panel: eccentricity distribution for these subsamples with [Fe/H] ≤ −1.

Figure 5 .
Figure 5. Number distributions of VMP/EMP stars with prograde orbits as a function of rotational velocity (v ϕ ), for the stars with [Fe/H] ≤ −2 (left panel), ≤ −2.5 (middle panel), and ≤ −3 (right panel), respectively.The vertical dashed line is at v ϕ = 150 km s −1 , which is used to select disk-like stars on highly prograde orbits.Note that the bins on v ϕ used for the stars with [Fe/H] ≤ −3 are twice the size for the more metal-rich stars, due to their lower numbers.dotted-dashedline show the approximate "troughs" in the distributions at IA = 0.25 and 0.65 rad, respectively, which can be used to roughly separate likely halo stars, thick-disk stars, and thin-disk stars.The numbers shown in each region listed in the figure reveal that VMP/EMP stars with prograde orbits dominate over those with retrograde orbits for IA ≤ 0.25 rad and those with 0.25 rad < IA ≤ 0.65 rad, and much less so for IA > 0.65 rad.One can reasonably associate the prograde stars with IA ≤ 0.25 rad with thin-disk orbits, those with 0.25 rad < IA ≤ 0.65 rad with thick-disk orbits, and those with IA > 0.65 rad with halo-like orbits.If we now specialize to the highly prograde stars with orbital velocities v ϕ > 150 km s −1 (indicated by the blue shaded region in Figure7), the relative dominance of the stars in the disk-like system for VMP/EMP stars becomes even clearer.

Figure 8
Figure 8 is a plot of Z max vs. R max for the stars with [Fe/H] ≤ −2, ≤ −2.5, and ≤ −3, in the left, middle, and right panels, respectively.The top panels of this figure show plots of the Z max distribution as a function of R max for the full sample of prograde stars (v ϕ > 0 km s −1 .)The dashed and dotted-dashed lines correspond to the troughs shown in Figure 7 at IA = 0.25 and 0.65 rad, respectively.The number of stars is provided in the legend at the top of each panel.The bottom panels apply to the stars on highly prograde orbits (v ϕ > 150 km s −1 ).Figure 9 shows histograms of the eccentricity distribution for stars with [Fe/H] ≤ −2, ≤ −2.5, and ≤ −3, respectively.The colors represent the same cuts on IA as in Figure 7: halo-like orbits with IA > 0.65 rad are shown with black lines, thick-disk orbits with 0.25 rad < IA ≤ 0.65 rad are shown with red lines, and thin-disk orbits with IA ≤ 0.25 rad are shown with blue lines.From inspection of the top row

Figure 6 .
Figure 6.Top panels: number distributions of VMP/EMP stars on both retrograde and prograde orbits as a function of eccentricity, from the left to right panel, for [Fe/H] ≤ −2, ≤ −2.5, and ≤ −3, respectively.The black, red, and blue solid lines indicate the stars with Zmax > 3 kpc, 1 kpc < Zmax ≤ 3 kpc, and Zmax ≤ 1 kpc, respectively.The number of stars in each region is indicated in the legend in the upper right corner of the panels.The dashed and dotted-dashed lines are shown at ecc = 0.2 and 0.4, respectively.Bottom panels: same as the top panels, but for the highly prograde stars with v ϕ > 150 km s −1 .

Figure 7 .
Figure 7. Number distribution of VMP/EMP stars as a function of IA, for [Fe/H] ≤ −2, ≤ −2.5, and ≤ −3, from the left to right, respectively.The total, prograde, and retrograde orbiting stars are shown with black, purple, and orange lines, respectively.The prograde stars with v ϕ > 150 km s −1 are shaded in blue.The number of stars in each region is indicated in the legend on the top of the panels.The dashed and dotted-dashed lines indicate IA = 0.25 and 0.65 rad, respectively.

Figure 8 .
Figure 8. Top panels: the Zmax distribution as a function of Rmax for the full sample of prograde stars.Dashed and dotted-dashed lines represent IA = 0.25 and 0.65 rad, respectively.The number of stars is shown in the legend at the top right of each panel.Bottom panels: same as the top panels, but for the stars with v ϕ > 150 km s −1 .

Figure 9 .
Figure 9. Top panels: Number distributions of VMP/EMP stars as a function of eccentricity, from the left to the right panel, for [Fe/H] ≤ −2, ≤ −2.5, and ≤ −3, respectively.The black, red, and blue solid lines indicate the stars with IA > 0.65, 0.25 < IA ≤ 0.65, and IA ≤ 0.25 radians, respectively.The dashed and dot-dashed lines are shown at eccentricity = 0.2 and 0.4, respectively.The number of stars in each region is indicated in the legend at the top right of the panels.Bottom panels: Same as in the top panels, but for the highly prograde stars with v ϕ > 150 km s −1 .

Figure 10 .
Figure 10.Fraction of prograde stars as a function of [Fe/H].The black line shows the combined SMSS/SAGES sample.The blue line shows the simulation results of Hirai et al. (2022).The Galactocentric distance of the data considered is confined to between 3 and 20 kpc, roughly corresponding to the combined SMSS/SAGES sample.The error bars shown are calculated using the normal approximation for the binomial proportions; they are quite small owing to the large number of stars in our data set.
and Galactic Evolution Survey (SAGES) is a multiband photometric project built and managed by the Research Group of the Stellar Abundance and Galactic Evolution of the National Astronomical Observatories, Chinese Academy of Sciences (NAOC).
supported in part by grant PHY 14-30152, Physics Frontier Center/JINA Center for the Evolution of the Elements (JINA-CEE), and by OISE-1927130: The International Research Network for Nuclear Astrophysics (IReNA), awarded by the US National Science Foundation.Y.S.L. acknowledges support from the National Research Foundation (NRF) of Korea grant funded by the Ministry of Science and ICT (NRF-2021R1A2C1008679 and RS-2024-00333766).Y.S.L. also gratefully acknowledges partial support for his visit to the University of Notre Dame from OISE-1927130: The International Research Network for Nuclear Astrophysics (IReNA), awarded by the US National Science Foundation.Y. Huang acknowledges support from the National Key R&D Program of China No. 2019YFA0405500 and National Natural Science Foundation of China grants 11903027 and 11833006.Y. Hirai was supported by JSPS KAKENHI Grant Numbers JP22KJ0157, JP21H04499, JP21K03614, and JP22H01259.Numerical computations were in part carried out on Cray XC50 and computers at the Center for Computational Astrophysics, National Astronomical Observatory of Japan.S.X.acknowledges support from the National Natural Science Foundation of China through project NSFC 12222301.K.T. was supported by the National Natural Science Foundation of China under grant Nos.12261141689, 12090044 and 12090040.G.Z. acknowledges support by the National Natural Science Foundation of China under grant nos.11988101 and 11890694.Description of the Candidate VMP/EMP Stars in the Combined SMSS/corrected G magnitude by Huang et al. (2022, 2023) for the Gaia DR3 [G C] − e GCmag in the Right Ascension from Gaia DR3 [pmra] mas yr −1 e pmRA The proper motion uncertainty in the Right Ascension from Gaia DR3 [pmra error] mas yr −1 pmDE The proper motion in the Declination from Gaia DR3 [pmdec] mas yr −1 e pmDE The proper motion uncertainty in the Declination from Gaia DR3 [pmdec error] mas yr −1 pmRApmDEcor The correlation coefficient between the proper motion in Right Ascension and in Declination from Gaia DR3 − Type Flag with "Dwarf" and "Giant" from Huang et al. (2022, 2023) − SubType Flag with "TO" for turn-off stars and "MS" for main-sequence stars from Huang et al. (2022, 2023) − vPHI The rotational velocity as given by AGAMA km s −1 e vPHI

Figure A1 .
Figure A1.Left column: comparison of the combined SMSS/SAGES photometric-metallicity estimates and spectroscopic metallicities for VMP/EMP stars, based on cross-matches to medium/high-resolution spectroscopic samples with available [Fe/H] and [C/Fe].The stars from Yoon et al. (2016), Li et al. (2022a), Placco et al. (2022), Zepeda et al. (2023), and Viswanathan et al. (2024) are shown as red triangles, blue squares, green diamonds, purple pentagons, and gray stars, respectively.The black circles indicate carbon-enhanced metal-poor (CEMP) stars that satisfy [C/Fe] > +0.7.Note that for our present purpose we employ the "as observed" [C/Fe], without applying evolutionary corrections (e.g., from Placco et al. 2014).The black solid line is a linear regression line for all metallicity regions except for the range [Fe/H]Literature ≤ −3, a light-blue shaded region, and the dashed line represents a one-to-one line.The legends in each panel indicate the number of matching stars (N), the p-value and the r 2 value in the Pearson correlation analysis, as well as the biweight location (µ) and scale (σ) of the metallicity residuals.Right column: comparison with the same spectroscopic catalogs, but excluding (recognized) CEMP stars.From top to bottom, the panels indicate matches for stars with metallicities based on the calibrated u − GBP colors, v − GBP colors, and when available, an average of both of these colors (Huang et al. 2022, 2023).

Table A2 .
1,496 VMP/EMP Disk-like Candidate Stars in the Combined SMSS/SAGES Sample NOTE-This table is a stub; the full table is available in the electronic edition.