Mapping the Milky Way with SDSS-V: Galactic Genesis Selection Function and Mock Catalog

Introduction

The Milky Way (MW) offers a uniquely detailed case study for understanding galaxy formation and evolution. With billions of individually resolved stars whose chemical and kinematic properties encode the Galaxy’s formation history, it remains an unparalleled laboratory for galactic archaeology (Grebel 1999; Feast & Whitelock 2000; Freeman & Bland-Hawthorn 2002; Sandage & Tammann 2006; Catchpole et al. 2016).

The rise of large-scale surveys have enabled a revolution in chemo-dynamic studies of the MW. The Gaia survey has provided parallaxes and proper motions for billions of stars [48, 48, 50]. Large spectroscopic surveys, e.g., RAVE (Steinmetz et al. 2006, 2020), LAMOST (Deng et al. 2012; Zhao et al. 2012), SDSS-III/IV (Eisenstein et al. 2011; Blanton et al. 2017; Kollmeier et al. 2026), and GALAH [35] (Martell et al. 2017), provide complementary radial velocities and chemical abundances.

These surveys have given rise to numerous studies mapping the энд correlations between positions, velocities and chemistry of stars in the MW. One of the most fundamental and informative correlations is the distribution of [$\alpha$/Fe]-[Fe/H] as a function of galactocentric cylindrical radius ($R$) and height ($Z$), which constrains the role of various processes governing the galactic disk, including chemical evolution and radial migration ([51] (Anders et al. 2014; Hayden et al. 2015)). More recently azimuthal variations in mean chemical abundances have been characterized and theorized to trace dynamical process, like bar-driven spiral density waves or moving groups ((Eilers et al. 2022; Poggio et al. 2022; Hawkins 2023; Filion et al. 2023; Hackshaw et al. 2024)). Moving beyond the Galactic disk, numerous works have focused on the halo in order to trace the MW’s accretion history (e.g., (Carollo et al. 2007; Nissen & Schuster 2010; Lee et al. 2017; Han et al. 2022; Chandra et al. 2025; Lucey et al. 2026b)).

The inferences drawn from these studies depend on how well the underlying samples represent the true stellar population, so it is essential to account for their incompleteness when mapping the MW. Each survey has their own targeting strategy, observational constraints, data processing pipeline, etc. which can all impact exactly which stars are included in the final released catalogs. Furthermore, depending on the specific science case, samples are further down-selected from the catalogs, typically including a quality cut, or identifying a certain stellar type. If not modeled carefully, these processes can introduce biases into the analysis. By defining the selection function, which gives the probability a given star is included in the sample as a function of relevant observational properties, these biases can be properly accounted for ((Bovy et al. 2016; Mints & Hekker 2019; Everall & Das 2020; Rix et al. 2021)).

Since selection functions are crucial for Galactic chemical cartography, they are already modeled for at least some part of the data for many of the large spectroscopic surveys, including SDSS-III/IV/V ((Bovy et al. 2012, 2014; Nidever et al. 2014; Cantat-Gaudin et al. 2024; Medan et al. 2026)), RAVE ((Wojno et al. 2017)), and LAMOST ((Chen et al. 2018)). For APOGEE, the main survey targets were drawn from 2MASS ((Zasowski et al. 2013, 2017)), so the empirical selection function can be estimated as the ratio of observed stars to available 2MASS sources in each field and color-magnitude bin ((Cantat-Gaudin et al. 2024)). Similarly, SDSS-V BOSS data can mostly be treated as a subset of Gaia ((Medan et al. 2026)). These empirical selection function models are projections on the sky, and can be converted into effective selection functions as a function of distance using a 3D dust map and a stellar population model drawn from isochrones ((Bovy et al. 2016; Cantat-Gaudin et al. 2024; Imig et al. 2025; Medan et al. 2026)). As the Gaia survey has provided pivotal data for Galactic structure modeling, there has been huge effort to quantify its selection function ((Bovy 2017; Boubert et al. 2020; Boubert & Everall 2020; Boubert et al. 2021; Everall et al. 2021; Rybizki et al. 2020, 2021; Everall & Boubert 2022; Cantat-Gaudin et al. 2023; Castro-Ginard et al. 2023)). For example, Rybizki et al. (2020) provided all-sky magnitude limit maps for Gaia DR2 subsets, Rybizki et al. (2021) derived the selection function of the radial velocity sample from its completeness relative to the full Gaia catalog, and Cantat-Gaudin et al. (2023) built an empirical model of the Gaia DR3 source catalog selection function in terms of $G$ magnitude and sky position. Selection functions have also been determined for specific stellar types like RR Lyrae, where completeness is estimated empirically by comparing independent catalogs as a function of sky position and distance ((Mateu et al. 2020)). This was updated for Gaia DR3 and a combined Gaia+PS1+ASAS-SN catalog by Mateu (2024).

One method commonly used to vet, and utilize selection functions is to create mock catalogs. There are generally two approaches to this. One samples empirically derived density distributions and typically assumes dynamical equilibrium ((Robin et al. 2003; Girardi et al. 2005; Jurić et al. 2008; Rybizki et al. 2018, 2020)). However, in order to include dynamical disequilibria beyond analytical models, the other approach samples star particles from MW-mass cosmological simulations ((Sharma et al. 2011; Lowing et al. 2015)). This method had been used to make mock Gaia catalogs ([61] (Sanderson et al. 2020; Nguyen et al. 2024)), as well as a mock SDSS-IV catalog ((Nikakhtar et al. 2021)).

As one of the main goals of the SDSS-V survey is to map the MW ((Kollmeier et al. 2026)), determining the selection function of the survey is crucial. We specifically focus on the Galactic Genesis (GG) survey (Bird et al. in prep), which uses the 2MASS ((Skrutskie et al. 2006)) and Gaia DR2 ([48]) surveys for targeting as described in § 2. The FIRE-2 MW-mass simulation, which are used to make a mock GG catalog, is presented in § 3. In § 4, we derive the selection function empirically, and in § 5 create and verify the mock catalog. Last, § 6 summarizes this work and the conclusions.

SDSS-V Milky Way Mapper

The SDSS-V survey is divided between three mappers: The Black Hole Mapper (BHM), The Local Volume Mapper (LVM) and the Milky Way Mapper (MWM). One of the main goals of the MWM is to study the Galaxy’s formation and evolutionary history using chemical cartography and galactic dynamics. This is achieved by collecting millions of optical and infrared spectra using the SDSS telescope ((Gunn et al. 2006)) in the northern hemisphere and the du Pont telescope ((Bowen & Vaughan 1973)) in the southern hemisphere. In this work, we focus on the infrared data from the APOGEE instruments, which have wavelength coverage of 1.51 – 1.70 $\mu$m and $R \approx22,500$ ((Wilson et al. 2019)). These spectra are the primary data for Galactic chemodynamic studies in the MWM survey.

The targeting of APOGEE data in SDSS-V is complex due to a wide variety of science goals in the survey. These goals are divided between programs, each with its own target selection procedure. For more details of all the programs and their targeting strategy we refer the reader to Kollmeier et al. (2026) and the SDSS website¹. Furthermore, a number of programs have released papers describing their target selection in more details, including the Magellanic Genesis program ((Nidever et al. 2026)) and the Young Stellar Objects program ((Kounkel et al. 2023)). However, the largest program in MWM is the GG survey, which we focus on in this work as it is specifically designed for population and MW structure studies.

Once targets are selected by a science program, there are still a number of steps that reduce the number of stars to what actually makes it into the final data release catalog. Firstly, there is a complex procedure for which targets actually get assigned to a fiber in the observing program, taking into account priority and fiber positioning constraints ((Medan et al. 2025; Blanton et al. 2026)). As all of the observations for BHM and MWM are planned at once, there are many competing programs which impact the fiber assignment in complicated ways. Weather and technical issues may cause deviations from this plan that can also impact which stars get observed and included in the final catalog. Finally, there may be steps in the data analysis pipeline that result in certain stars not being included in specific catalogs for a variety of reasons, including data artifacts, or poor calibrations. These can all be difficult to track which motivates an empirical estimate of the selection function, which is informed by the targeting strategy, rather than trying to model directly each of these effects.

The Galactic Genesis Survey

The main goal of the GG survey is to obtain millions of APOGEE spectra of red giant stars across the sky with minimal selection biases, in order to study the MW’s history through the chemodynamics of its stars (Bird et al. in prep). Therefore, the targeting strategy is intentionally simple, using 2MASS and Gaia DR2 photometry. However, SDSS-V includes data from SDSS-IV which used different targeting strategy and fiber plug plates (Zasowski et al. 2017; Beaton et al. 2021; Santana et al. 2021). Furthermore, small changes in targeting selection criteria occurred after SDSS-V DR19. The targeting flags² in the final data products delivered by the astra pipeline (Casey et al. in prep) label any data that pass the targeting selection criteria for a given program agnostic to the targeting strategy used when it was observed. Therefore, we include SDSS-IV data that meets the SDSS-V DR19 targeting criteria:

[139] (Almeida et al. 2023). For data releases beyond DR20, this criteria will expand to:

Therefore, in data releases beyond DR20, the SDSS-V GG survey may include more data from SDSS-IV because of this additional or statement. In this work, we use the latest public APOGEE release which is DR19, since there was no new APOGEE data released in DR20 [138].

APOGEE also has a nominal bright limit of $H > 7$ but, using fiber offsets, a handful of stars with $4.5 < H \le7$ have been observed. This selection targets very red stars, including stars too red to be observed by Gaia, as GG is primarily interested in cool giants embedded in the dusty Galactic disk, particularly in the highly extincted inner Galaxy. The parallax cut is intended to select giant stars and corresponds to $H_{\mathrm{abs}} < -2.5$. It is important to note that the “or” statements mean that not every star in the sample meets all three criteria within the bracket, but at least meets one of them. For example, stars with $G-H < 5$ may be included if they satisfy the parallax-$H$ criterion.

In this work, we focus primarily on results from the ASPCAP pipeline [53] (Holtzman et al. 2018; Jönsson et al. 2020), which provides stellar parameters and abundances for APOGEE spectra and is the most widely used pipeline for Galactic studies (e.g., [129] (Hayden et al. 2015; Feuillet et al. 2019; Imig et al. 2023; Griffith et al. 2024)). We use this final catalog (astraAllStarASPCAP-0.6.0.fits) rather than an earlier-stage product, such as the list of all observed stars, because we aim to capture the full impact of the data analysis pipeline on the selection function. These are the products actually used in scientific analyses.

FIRE-2 Cosmological Zoom-in Simulations

We create the mock SDSS-V catalog by employing the publicly available³ Wetzel et al. 2023 Latte suite of seven isolated MW-mass galaxies (Wetzel et al. 2016), together with the ELVIS suite of three Local Group-like pairs of comparable mass (Garrison-Kimmel et al. 2019a). Each simulation adopts the FIRE-2 physics model (Hopkins et al. 2018b) and is evolved with the GIZMO⁴ gravity and hydrodynamics code, run in meshless finite-mass (MFM) mode (Hopkins 2015). Full implementation details are given in these references; here we summarize the properties most relevant to this work. All runs assume a flat $\Lambda$-CDM cosmology consistent with Planck Collaboration et al. [125]. Except for m12w, the Latte suite adopts $\Omega_{\mathrm{m}}=0.272$, $\Omega_{\mathrm{b}}=0.0455$, $\sigma_8=0.807$, $n_{\mathrm{s}}=0.961$, $h=0.702$. The ELVIS pairs Thelma & Louise and Romulus & Remus retain the cosmology of the original dark-matter-only ELVIS suite: $\Omega_{\mathrm{m}}=0.266$, $\Omega_{\mathrm{b}}=0.0449$, $\sigma_8=0.801$, $n_{\mathrm{s}}=0.963$, $h=0.71$. The remaining pair, Romeo & Juliet, along with m12w, instead follow the updated parameters Planck Collaboration et al. [126]: $\Omega_{\mathrm{m}}=0.31$, $\Omega_{\mathrm{b}}=0.048$, $\sigma_8=0.82$, $n_{\mathrm{s}}=0.97$, $h=0.68$.

The mass distribution of simulated galaxies is known to depend sensitively on how stellar feedback is treated (Pontzen & Governato 2014; Lazar et al. 2020). In FIRE-2, star formation is self-regulated through a combination of stellar winds, radiation pressure from young stars, Type II and Type Ia supernovae, photoelectric heating, and photoionization. Star formation occurs above a gas density threshold of $n_{\mathrm{SF}} > 1000\ \mathrm{cm}^{-3}$, with feedback rates, luminosities, energies, and mass-loss rates taken directly from stellar evolution models (STARBURST99; Leitherer et al. 1999). Across the 13 simulated galaxies, $z=0$ dark matter halo masses span $M_{200\mathrm{m}}=1-2.1\times10^{12}\ M_{\odot}$ (Samuel et al. 2020; Barry et al. 2023). Baryonic resolution differs by suite: initial stellar particle masses are $7070\ M_{\odot}$ for Latte and $3500\ M_{\odot}$ for the higher-resolution ELVIS runs, with star particle softening lengths of $\approx4$ pc and dark matter softening of $\approx40$ pc. Corresponding dark matter particle masses are $35{,}000\ M_{\odot}$ for Latte and $19{,}000\ M_{\odot}$ for ELVIS. These choices give numerical convergence at the $\approx100$ pc scale (Wetzel et al. 2016), well below the $\lesssim950$ pc disk scale heights these simulations resolve (Yu et al. 2021). Numerical effects should therefore not strongly impact our analysis of the stellar halo, which focuses on scales $>1$ kpc.

FIRE-2 galaxies provide useful mock catalogs of the MW because they broadly match several key MW properties, including the stellar halo mass fraction (Sanderson et al. 2018), the presence of a metal-rich in-situ halo component (Bonaca et al. 2017), and the radial and vertical structure of the stellar disk (Ma et al. 2017; Sanderson et al. 2020; Bellardini et al. 2021; McCluskey et al. 2024). Simulated satellite populations around both the MW and M31 analogs are likewise consistent with observations (Wetzel et al. 2016; Samuel et al. 2020; Garrison-Kimmel et al. 2019b; Panithanpaisal et al. 2021; Cunningham et al. 2022). Some tensions with MW data remain, however: although FIRE-2 reproduces the observed stellar masses and number of tidal streams, Shipp et al. (2023) found that stream orbital properties are not well reproduced. Latte disks also tend to assemble later than observational estimates of the MW disk age suggest (McCluskey et al. 2024). While these differences should be kept in mind, they do not invalidate the utility of the FIRE-2 simulations for constructing mock catalogs.

Selection Function

As described in § 2, the selection of targets in SDSS-V is a complex process. The data collection, reduction and analysis process further impacts the final sample available for scientific inferences. As these effects are difficult to track, it is best to instead empirically model the selection function. To achieve this we follow a procedure similar to Cantat-Gaudin et al. (2024) and Medan et al. (2026), where we estimate the completeness as a function of sky position, and apparent magnitude ($H$) with an optional additional color axis ($G-H$). Therefore, our selection function assumes that stars are exchangeable within each cell of sky position, $H$, and $G-H$, so each cell value is the probability, averaged over the cell, that a 2MASS source reaches the final ASPCAP catalog, even though individual probabilities vary with factors such as fiber collisions and competing targets. Our spatial cells range from $\sim1^\circ$ across at $N_{\mathrm{side}}=64$ to $\sim7^\circ$ across at $N_{\mathrm{side}}=8$ where counts are sparse, comparable to the $3^\circ$ diameter of an SDSS-V field, so each cell averages over the small-scale effects of fiber collisions. This within-cell variation does not bias population-level inference provided it is independent of the stellar properties being inferred, which holds for fiber assignment because it depends on sky position, carton priority, and robot reach rather than on stellar chemistry or kinematics. However, this is only valid within the possible parameter range of the ASPCAP catalog given by the model grid (Jönsson et al. 2020). For example, ASPCAP fails for stars with $[\mathrm{M}/\mathrm{H}]<-2.5$, which are therefore absent from the catalog regardless of their position, magnitude, or color.

As SDSS-V primarily uses the 2MASS survey (Cutri et al. 2003; Skrutskie et al. 2006) for APOGEE targeting, we can treat it as a subsample of the 2MASS catalog and use the following equation:

where $S_C(\mathbf{q}\mid\mathbf{q}\ \text{in parent})$ is the probability a source with attributes $\mathbf{q}=\{l,b,H,\ldots\}$ will be in the subsample given that it is in the parent catalog. $S_C^{\mathrm{parent}}(\mathbf{q})$ is the corresponding probability that a given source will be in the parent catalog.

The 2MASS point source catalog is 99% complete at H=15 in high-latitude regions where source confusion is negligible⁵. Towards the Galactic center, the completeness limit occurs at significantly brighter magnitudes. While the completeness here is not measured as robustly, the analysis available indicates that it is safe to assume the 2MASS point source catalog is complete up to H=12⁶. If this assumption over estimates the completeness of 2MASS the SDSS-V selection function will also be over-estimated. However, there are limited sources in SDSS-V with H>12 and the GG survey specifically only targets stars with H<11. Therefore, the impact on the overall calculation of the selection function should be limited and we assume $S_C^{\mathrm{parent}}(\mathbf{q})=1$ within our $\mathbf{q}$ bounds.

For sources with attributes $\mathbf{q}$, we define $p(\mathbf{q})$ as the probability that a single source is selected into the sub-sample. Given $n$ sources in the parent and $k$ sources in the subsample catalog sharing attributes $\mathbf{q}$, the likelihood of $p(\mathbf{q})$

Using an uninformative uniform prior, the posterior can be described as a Beta$(k+1,n-k+1)$. We use two methods of deriving a point estimate for $S_C^{\mathrm{parent}}(\mathbf{q})=p(\mathbf{q})$ from this posterior. One method is the mean or expected value of the posterior which is given by:

We also use the maximum a posteriori (MAP) which is given by:

The expectation value is typically used as the estimate of $S_C^{\mathrm{parent}}(\mathbf{q})$ (e.g., Castro-Ginard et al. 2023; Medan et al. 2026). However, the expectation value is never zero even in bins with zero observed targets. In order to create a mock catalog which accurately recreates the sky and magnitude distribution of the observed DR19 GG sample, we require bins with zero targets to have a selection function value of zero, so mock stars in regions of the sky and magnitude space that SDSS did not observe do not leak into the catalog. The MAP value achieves this (see § 5 and Figure 4 and Figure 5). On the other hand, however, the MAP value does not distinguish a cell with a genuinely zero selection probability from a sparsely sampled cell in which no stars happened to be observed.

We calculate the expectation values and MAP in bins of $\mathbf{q}$. Given that the targeting strategy for the GG survey contains the color term $G-H$ (see Equation (1)), we also include it as an optional additional dimension of $\mathbf{q}$. However, there are sources in 2MASS that do not have a reported Gaia DR2 $G$ magnitude, likely mostly because of dust extinction. For those sources, we instead calculate their minimum G-H, given their H magnitude and the Gaia DR2 faint magnitude limit of G=21 (Gaia Collaboration et al. [48]). We note that only 2.5% of 2MASS sources ($\approx350,000$ stars) with $4.5<H<11$ do not have a match in Gaia. Therefore, this $G-H$ approximation will not have a large effect on our selection function calculation, but it is still important to include in the calculation of $n$ so that the completeness is not overestimated.

The GG targeting strategy also includes stars that do not have a reported Gaia DR2 $G$ magnitude. When the $G-H$ dimension is included in $\mathbf{q}$, we exclude these stars from the numerator ($k$), so in this case we are specifically calculating the selection function of the GG sample with Gaia photometry. These stars make up a small fraction of the GG sample, only 5,967 out of 235,219 stars (2.5%) in DR19. Their missing Gaia photometry in SDSS-V catalogs may arise in part from differences between cross-match chains or from bookkeeping issues, rather than a true non-detection. Therefore, assigning these stars $G=21$ could spuriously inflate the selection function at red $G-H$, so we exclude them from the color selection function.

When the $G-H$ term is not included in $\mathbf{q}$, we simply do not use Gaia data at all and the selection function is computed only as a function of sky positions and $H$ magnitude, rather than marginalizing over $G-H$. Therefore, stars missing Gaia $G$ magnitudes can easily be included in both numerator ($k$) and denominator ($n$).

We nominally use bins of 0.5 mag in both H and G-H. We note that in the accompanying Python package sdssv_apogee_sf (Lucey 2026) (see Appendix E), users can specify their own binning. The binning in sky coordinates ($l,b$) is performed using HEALPix (Górski et al. 2005; Zonca et al. 2019). We use adaptive binning to decrease the resolution until a minimum of 20 2MASS objects per bin ($n\geq20$). With this limit, the statistical uncertainty on the selection probability in any cell is $\sim0.1$, and a cell with no observed stars still constrains the selection probability to $\lesssim0.14$ at 95% confidence. In sparsely populated ($H,G-H$) bins, the adaptive scheme must coarsen the HEALPix resolution until each cell contains at least 20 2MASS stars, which can wash out real small-scale structure in the survey footprint. To recover this structure, we assume that the small-scale angular pattern of the completeness varies only weakly with $H$ over broad magnitude ranges, so that it can be measured in wider $H$ bins and shared among the native bins they contain. We therefore merge the native $H$ bins into two such bins ($H < 7$ and $7 \le H < 11$), whose larger star counts allow the completeness to be estimated at finer angular resolution. Within each native cell, we then redistribute the original completeness across the finer sub-pixels in proportion to the completeness of the corresponding merged bin, normalized so that the expected number of observed stars in the cell is exactly preserved. The hybrid therefore leaves the selection function unchanged on the original coarse scale and only adds angular detail within each cell. In high stellar density regions, we achieve fine HEALPix resolution with $n_{\mathrm{side}}=64$, which corresponds to angular resolution of

$0.9^\circ$ and a physical resolution of $\approx130$ pc at the Galactic center. In the lowest density regions, we set a minimum bin size of $nside=8$ with angular resolution of $7.33^\circ$ and a physical resolution of $\approx1$ kpc at the Galactic center. In Figure 1, we show the observed DR19 GG survey sky distribution (top panel) compared to the empirical selection function. The middle panel shows the mean of the posterior version of the selection function, while the bottom panel shows the MAP version. In general, the selection function has the highest completeness in regions of the sky where the GG sample has its highest density. However, this is not always the case, especially in the inner Galaxy, where the GG sample has many targets, but the total underlying number of stars there is very high so that the completeness is still quite low. Certain areas of the sky have more complete coverage because of the observing strategy, including the Magellanic clouds which were prioritized in SDSS-IV. We also see the remnants of the SDSS-IV pencil-beam footprint, where SDSS-V is now filling in with full-sky coverage. The mean and MAP versions of the selection function are generally very similar except for one important difference. In the regions of sky where there are no observed GG targets, the mean selection function value is small but nonzero (because of the +1 term in the numerator), while the MAP version report exactly zero. This has important consequences when applying the selection function for creating a mock catalog as shown in § 5.

We further investigate the color and magnitude dependence of the selection function in Figure 2. The left column shows the mean version of the selection function, while the right column shows the MAP version. The top panels shows the bright end of APOGEE observations with $7 \le H < 7.5$, while the second row shows the faint end with $10 \le H < 11$. The bottom rows show the selection function evaluated at the faint end with two different $G-H$ bins. Specifically, we show one of the most densely populated bins that still has full sky coverage in the third row from the top ($2.5 \le G-H < 3$). The $H$, $G-H$ bin with the highest density of targets is $10.5 \le H < 11$ and $3.5 \le G-H < 4$ and is shown in the bottom row. The selection function is only calculated here along the disk and the Magellanic clouds, because of the lack of 2MASS targets elsewhere in the galaxy that are sufficiently red. These stars only exist in areas with significant dust extinction.

Figure 3 shows the observed distribution of $H$, $G-H$ (leftmost panel) compared to the empirical selection function. The middle and right panels show the selection function as a function of $H$, $G-H$ marginalized across the sky with the mean version in the middle and the MAP version in the right panel. The bright and faint limits ($7 < H < 11$) of the APOGEE GG survey are clearly visible here. There are a few stars with $H < 7$, observed with fiber offsets, but they are a small fraction of the total sample. The grey bins do not have a selection function value because they do not meet the criteria of having at least 20 stars in the 2MASS catalog. Again the difference between the mean and MAP versions of the catalog is visible here where we see bins with no GG observations have selection function value of exactly zero in the MAP version, but have small nonzero values in the mean version.

Mock Catalogs

Using this empirical selection function, we build mock catalogs designed to reproduce the SDSS-V APOGEE survey. Mock catalogs of this kind are valuable tools for comparing simulations or theoretical models directly to MW data ([61] (Robin et al. 2003; Sharma et al. 2011; Sanderson et al. 2020)). By imposing the same selection effects present in the real survey, they let us test which physical quantities are actually measurable given the data in hand, and identify observational selection effects that could bias the resulting scientific conclusions (Nandakumar et al. 2017; Mints & Hekker 2019; Everall & Das 2020).

As we are interested in the SDSS-V APOGEE sample that is in Gaia DR3, we use the mock Gaia DR3 catalogs from Nguyen et al. (2024) as the starting point for our SDSS-V mocks. We refer the reader to Nguyen et al. (2024), for the full details of the process to create these catalogs. In short, these catalogs are built from the Latte suite of FIRE-2 cosmological zoom-in baryonic simulations of MW-mass galaxies (Wetzel et al. 2016), part of the Feedback In Realistic Environments (FIRE) simulation project, which were run using the Gizmo gravity plus hydrodynamics code in meshless finitemass (MFM) mode (Hopkins 2015) and the FIRE-2 physics model (Hopkins et al. 2018a). They use the Ananke framework to transform star particles from the simulations into synthetic stars with Gaia DR3-like parallaxes and photometry (Sanderson et al. 2020; Thob et al. 2024). Specifically, Ananke uses the Galaxia code (Sharma et al. 2011) with the Kroupa (2001) initial mass function (IMF). Nguyen et al. (2024) specifically uses the updated Gaia DR3 Padova CMD v3.6⁷ isochrones. Since these simulations do not directly model the formation of dust grains, the extinction is assumed to be proportional to the metal-weighted column density of hydrogen along the line of sight (Rachford et al. 2009; Sanderson et al. 2020). The detailed Gaia DR3 empirical selection function (Cantat-Gaudin et al. 2023) is not applied and instead a simple G cut ($3<G<21$) based on the Gaia DR3

limiting magnitudes is employed. Therefore, this mock overselects stars, as it does not account for the lower completeness in regions affected by crowding and the Gaia scanning law. However, this does not impact its usefulness as a starting point for our mock catalog where we will apply the SDSS-V empirical selection function.

The first step to build the SDSS-V APOGEE catalogs is to infer 2MASS magnitudes for the mock Gaia DR3 catalogs. We use the Galaxia code to calculate the 2MASS magnitudes using Padova isochrones (Girardi et al. 2002; Cohen et al. 2003; Marigo et al. 2008). We then calculate the dust extinction using $A_J/A_V$, $A_H/A_V$, $A_{K_s}/A_V$ from Wang & Chen (2019). We calculate two versions of the true apparent magnitudes. The first uses the E(B-V) values reported in Nguyen et al. (2024), which are based on the distribution of metal-rich gas in the simulations, so that the star and dust distributions are self-consistent. This self-consistency is most important for regions with ongoing star formation, where the dense gas and young stars are most strongly correlated (Sanderson et al. 2020). The second version instead calculates the E(B-V) from MW

dust maps, specifically using the Bayestar package (Green 2018) using the DECaPS2-based map for the southern hemisphere (Zucker et al. 2025) and the map from Green et al. (2019) for the northern hemisphere with Pan-STARRS coverage (Chambers et al. 2016). In the MW dust case we also calculate the extinction in the Gaia G band using the dustapprox package (Fouesneau et al. 2025) which calculates the $k_G = A_G/A_0$ as a function of $A_0$ (the monochromatic extinction at 550 nm) and effective temperature ($T_{\mathrm{eff}}$), consistent with Nguyen et al. (2024). Since the dust in both the MW dust map and the self-consistent extinction estimation primarily lies on the Galactic plane, either version will give similar results for older stellar populations.

To inject realistic photometric noise into synthetic 2MASS/APOGEE catalogs, we characterize the magnitude-dependent uncertainty reported by the 2MASS Point Source Catalog (PSC) empirically, rather than assuming an analytic error curve. For each of the three 2MASS bands ($J$, $H$, $K_s$), we query the full PSC and bin every source by its catalog magnitude in that band using 0.1 mag bins. spanning $-5 \le m \le17$, wide enough to cover the saturated-bright to noise-floor-faint ends of all three bands. Within each bin we compute the source count and the mean $\mu(m)$ and standard deviation $s(m)$ of that band’s combined photometric-uncertainty, reducing the empirical uncertainty distribution in each bin to a Gaussian approximation, $\sigma\sim\mathcal{N}(\mu(m),s(m))$. Bins containing fewer than 20 sources are flagged as unreliable and discarded, so small-number statistics do not drive the model at the sparsely populated bright and faint extremes. See Appendix § 9 for a plot of the noise models. For the recalculated Gaia G magnitude we also recalculate the Gaia DR3 G mag uncertainty using the pygaia uncertainty model (Riello et al. 2021). To simulate photometry for the mock catalogs, the uncertainty ($\sigma$) then perturbs the true magnitude via $m_{\mathrm{obs}} = m_{\mathrm{true}} + \mathcal{N}(0,\sigma)$.

Each mock star’s position and (noisy) $H$, $G-H$ are passed to the selection function query, returning the probability that a real 2MASS star at that position and magnitude would have been targeted and observed by APOGEE. This probability is then applied per star as a Bernoulli accept/reject draw (retained with probability $S$, discarded otherwise), producing a synthetic "APOGEE-observed" mock catalog that carries the same selection effects, including its stochastic incompleteness pattern, as the real data. The resulting mock catalog can then be compared directly to the observed APOGEE sample without needing to reproduce the underlying target-selection algorithm itself.

In addition to the stellar parameters ($T_{\mathrm{eff}}$ and $\log g$) the FIRE-2 simulations have a number of elemental abundances that are in common with APOGEE. These elements are Fe, O, Mg, Si, Ca, S, C and N. We include these elemental abundances for our mock stars with injected APOGEE/ASPCAP-like noise. Similar to the procedure for 2MASS magnitudes, we approximate the distribution of APOGEE signal-to-noise (SNR) as a function of H band magnitude with a Gaussian mean and standard deviation and use this to probabilistically assign an SNR to each mock star. We again follow a similar procedure for deriving an empirical relation between APOGEE SNR to ASPCAP stellar parameters ($T_{\mathrm{eff}}$ and $\log g$) as well as the chemical abundances, [Fe/H], [O/H], [Mg/H], [Si/H], [Ca/H], [S/H], [C/H] and [N/H]. While this procedure reliably recreates ASPCAP-like uncertainties, it is well-known that DR19 ASPCAP uncertainties are underestimated and do not take into account systematics (Mészáros et al. 2025).

In Figure 4 we show the sky distribution colored by mean metallicity of the DR19 ASPCAP GG sample (top panel) compared to the mock catalogs calculated with the mean (middle panel) and the MAP (bottom panel) versions of the selection function. Although some fields are more complete, much of the DR19 sky distribution is dominated by the pencil beam design of SDSS-IV. SDSS-V will continue to fill this in more in future releases. The mock catalog from the MAP version of the selection function better matches the sky distribution of the real data. This is because the MAP selection function is exactly zero in pixels with zero observed stars, but for the mean version the minimum value is $3.95 \times10^{-5}$ because the numerator is never zero (see Equation (5)). The MAP version of the mock is still more dispersed across the sky than the observed sample because of the spatial resolution in the selection function estimate. Higher spatial resolution could be achieved at the cost of coarsening the magnitude binning or decreasing the minimum number of 2MASS observations per bin and therefore increasing the noise in the completeness estimate.

Overall, the mock contains more stars out of the plane than the DR19 data. Specifically, only 5.7% of the DR19 stars have $|b| > 15^\circ$, while this represents 23.5% of the MAP mock. This is partially driven by the lower spatial resolution that can be achieved in the selection function at these higher latitudes where 2MASS stellar counts are more sparse. Furthermore, there is evidence that the FIRE-2 cosmological simulations have more massive stellar halos than the MW (Lucey et al. 2026a). This may also contribute to why the average [Fe/H] at high galactic latitudes ($|b| > 15^\circ$) is +0.35 dex more metal-rich in the MAP mock than the DR19 data. The average [Fe/H] within the galactic plane ($|b| < 5^\circ$) for the MAP mock and DR19 data agree to within 0.06 dex.

SDSS-IV completed an extensive mapping of the LMC and SMC with APOGEE (Nidever et al. 2020; Hasselquist et al. 2021). A lot of this data also satisfies the target selection of the GG program, leading to the MC coverage as a significant feature in the DR19 ASPCAP GG sky distribution (see the top panel of Figure 4). While we recreate this feature in the sky distribution, it is important to note that we do not actually ensure the presence of an MC-like satellites at those sky positions in the FIRE-2 simulations or resulting mocks. Therefore, it is essentially just a random halo field which is significantly more metal-rich than the MCs.

Figure 5 compares the apparent magnitude and color distribution of the ASPCAP GG DR19 sample (top left panel) to the mock catalogs calculated from the mean (top middle panel) and MAP (top right panel) versions of the selection function. The bottom panels show the corresponding Kiel diagram. There are clear differences between the observed and the mock catalogs’ color magnitude distributions, although the MAP version is more similar to the observed than the mean version. These differences are partially driven by the coarseness of magnitude binning in the selection function estimate (see Figure 3). The mock catalog contains more targets on the edges of the observed distribution where the selection function estimate is small but not zero. This effect is heightened in the mean version where bins with zero observed targets have a non-zero, but small selection function estimate.

Overall the $T_{\mathrm{eff}}$ distribution of stars is much hotter in the mock than the DR19 observed data. For example, 97.7% of the DR19 data has $T_{\mathrm{eff}} < 5000$ K, while the same is true for only 82.9% of the MAP mock. Furthermore, there are 4,512 stars in the MAP mock with $T_{\mathrm{eff}} > 7000$ K, while the DR19 data has only 42. The mocks also contain a higher percentage of turn-off/sub-giant stars, relative to the DR19 data. 17% of the mock has $3 \le\log g \le4$, while the same is true for only 1% of the DR19 data.

Some of the differences in the Kiel diagram between the observed ASPCAP DR19 GG sample and the mocks are driven by the difference in the present-day star formation rate between the MW and the simulated galaxies; the simulations have a local SFR several times that of the MW (Sanderson et al. 2020) resulting in a larger proportion of young stars in their mock catalogs. Additionally, while the targeting is based on H and G-H, the ASPCAP analysis pipeline further limits the $T_{\mathrm{eff}}$ and $\log g$ range. ASPCAP cannot measure $\log g > 5$ or $\log g < -0.5$. While ASPCAP can handle stars with 2500 K $\le T_{\mathrm{eff}} \le20,000$ K, but only for stars with $\log g < 3$. The hot star range is clipped to $T_{\mathrm{eff}} \le8,000$ K for stars with $3 \le\log g \le1$ and further clipped to $T_{\mathrm{eff}} \le6,000$ K for $1 \le\log g \le0$ and $T_{\mathrm{eff}} \le4,000$ K for $0 \le\log g \le-0.5$. Furthermore, ASPCAP cannot measure parameters for stars with metallicity $< -2.5$, which may contribute to the lack of blue horizontal branch stars in the observed sample. While ASPCAP may be capable of deriving parameters within the range described above, it can fail more frequently in certain regimes. For example, 63% of stars with G-H$<2$ are flagged as bad by ASPCAP, which roughly corresponds to $T_{\mathrm{eff}} > 6000$ K. The mocks also have sharper features (e.g. red clump) in the Kiel diagrams because the ASPCAP uncertainties are underestimated, including only the statistical contribution and excluding the larger systematic uncertainties (see Appendix § 9; Mészáros et al. 2025).

The evolved stellar populations are also sensitive to the choice of stellar models used to create the isochrones that sample the stellar populations for the mocks (e.g. Marigo et al. 2013; Rosenfield et al. 2016; Pastorelli et al. 2019, 2020), and to the choice of models of the most massive stars (Tang et al. 2014). These models have been calibrated against observations but only in a limited range of metallicities compared to the range needed to sample all the simulated stellar populations. Any of these can lead to differences with the observed ASPCAP catalog in the regions where the discrepancies are most apparent in Figure 5. This specifically drives the strong diagonal feature near the top of the giant branch which are mostly stars with age $\ge10$ Gyr.

Azimuthal Variations Example

One of the most prominent observed chemical features in the MW disk is the negative radial metallicity gradient (an incomprehensive list includes: Mayor 1976; Andrievsky et al. 2002; Magrini et al. 2009; Luck & Lambert 2011; Bergemann et al. 2014; Xiang et al. 2015; Yan et al. 2019; Hawkins 2022; Imig et al. 2023). This gradient is theorized to be the result of the ‘inside-out’ formation of the Galaxy, by which the inner Galaxy forms first and fast, igniting rapid star formation and increasing the metal content of stars located towards the Galactic center (e.g., Matteucci & Francois 1989; Chiappini et al. 1997).

On top of this gradient, azimuthal substructure in metallicity has been reported on the $\sim0.1$ dex-level. Different tracers have exhibited this azimuthal structure, including hot OBA stars (Hawkins 2022) in LAMOST, hot giant stars in Gaia (Poggio et al. 2022), and open clusters (Otto et al. 2026). While the origin of these azimuthal metallicity variations likely have dynamical origins (Jurado et al. 2026), it is apparent that they are a key physical feature of the disk.

To identify if this azimuthal substructure is present in the presented mock catalogs, we follow the methodology put forth in Hackshaw et al. (2024) (hereafter H24). H24 used APOGEE DR17 stars to identify a negative radial metallicity gradient of $\Delta[\mathrm{Fe}/\mathrm{H}]/\Delta R \sim-0.0678 \pm0.0004$ dex/kpc. They identified azimuthal variations on top of this gradient, with the magnitude of these azimuthal variations correlating with various dynamical parameters. Similarly, we identified cool ($3500 < T_{\mathrm{eff}} < 5000$ K) giant (log $g < 3.6$ dex) stars in the DR19 astroNN GG catalog. We removed any stars with errors greater than 0.1 dex on [Fe/H] and [Mg/H], as well as any stars errors greater than 150 K and 0.3 dex on $T_{\mathrm{eff}}$ and log $g$, respectively. We then selected stars within 0.3 kpc of the Galactic midplane and $3.5 < R_{\mathrm{galactocentric}} < 15.5$ kpc to select stars confined to the plane and minimize the effects of the bar in the inner Galaxy (i.e., Hammersley et al. 1994; Wegg et al. 2015; Lucey et al. 2023).

To further investigate the impact of the selection function on the derived azimuthal substructure, we create a mock version that uses a uniform selection which sets the ground truth. Specifically, this selection randomly down samples the Gaia DR3 mock with a 5% acceptance rate. As the Gaia mock merely applies a selection based on G band magnitude that is not a function of the sky position or color, this sample gives an unbiased representation of the true underlying disk population.

Following H24, we identify negative radial metallicity gradients in the resulting DR19 astroNN GG sample and both the m12i mock samples. Azimuthal substructure in $\delta[\mathrm{Fe}/\mathrm{H}]$ is identified in all samples after subtracting the radial gradient off of the sample metallicities. Figure 6 shows the gradients, model metallicities, and resulting azimuthal variations in both samples. It is already well established that FIRE-2 galaxies do not reproduce the strength of the MW’s radial metallicity gradient (Graf et al. 2025). Here, we confirm that accounting for selection effects does not resolve this discrepancy. However, the azimuthal substructure in the m12i MAP mock differs substantially from that of the uniform m12i mock. This suggests that the selection function also shapes the observed azimuthal structure of the MW and must be accounted for to recover its substructure.

Summary and Conclusions

The MW remains an unparalleled laboratory for detailed studies of galaxy formation, and large-scale surveys like Gaia and SDSS-V are essential for mapping its chemodynamic structure. Robust mapping requires accounting for each survey’s selection function, the probability that a given star is included in the sample as a function of observational properties. While selection functions have been derived for many major surveys (e.g., Wojno et al. 2017; Chen et al. 2018; Cantat-Gaudin et al. 2023, 2024), this had not yet been done for the SDSS-V GG survey, but is a key ingredient to achieving many of the survey’s science goals.

In this work, we calculate the SDSS-V GG selection function as a subsample to its 2MASS parent catalog, which we assume to be complete within the relevant magnitude range of $4<$H$<11$. Based on the targeting strategy, we calculate the selection function as a function of sky position, H and G-H. We calculate two versions, one using the mean of the posterior consistent with previous works (e.g., Castro-Ginard et al. 2023; Medan et al. 2026) and one with the MAP. The main difference is that the MAP version has the selection function equal to zero in bins where zero stars are observed while the mean version instead reports a small but nonzero number.

With this selection function, we also utilize FIRE-2 MW-mass galaxy simulations to create mock DR19 GG ASPCAP catalogs, building upon the Gaia DR3 mocks from Nguyen et al. (2024). We create two versions of the model photometry one with MW dust maps (Green et al. 2019; Zucker et al. 2025) and another with the dust distribution consistent with the simulation. We inject noise by deriving a 2MASS and ASPCAP noise model based on apparent magnitude. We also provide further mock catalog variations based on whether the mean or MAP version of the selection function is used, as well as using different FIRE-2 simulations.

We work a short example use case of the mock catalogs by comparing the azimuthal variations between the observed and mock data following the procedure presented in (Hackshaw et al. 2024). This comparison shows that the selection function substantially alters the observed azimuthal substructure, highlighting the importance of modeling it when interpreting the MW’s chemical structure. In general, the selection function and corresponding mock catalog provide crucial tools for reaching the full scientific capability of the SDSS-V GG survey. For future data releases (DR21 and DR22), the mock catalogs will be provided as a Value Added Catalog (VAC) and the selection function will be updated and available on the corresponding Github package⁸ (see Appendix E).

SOFTWARE

This research made use of Astropy (Astropy Collaboration et al. [4, 5, 6]), NumPy (Harris et al. 2020), Matplotlib (Hunter

2007), pandas (Wes McKinney 2010; pandas development team 2020), astroquery (Ginsburg et al. 2019), healpy (Zonca et al. 2019; Górski et al. 2005), and h5py (Collette et al. 2022). Generative AI (Anthropic’s Claude, via Claude Code) was used to help write and debug analysis code and the sdssv_apogee_sf package, and edit portions of the text. All AI-assisted output was reviewed, tested, and revised by the authors, who takes full responsibility for the content of this work.

Acknowledgements

This material is based upon work supported by the National Science Foundation under Award No. 2303831. RES gratefully acknowledges travel support from the Simons Foundation under award 1018462. APJ acknowledges support from NSF grant AST-2510795, the Alfred P. Sloan Research Fellowship, and the University of Chicago’s Research Computing Center. J.I. gratefully acknowledges support from STScI DDRF project D0001.82536

Funding for the Sloan Digital Sky Survey V has been provided by the Alfred P. Sloan Foundation, the Heising-Simons Foundation, the National Science Foundation, and the Participating Institutions. SDSS acknowledges support and resources from the Center for High-Performance Computing at the University of Utah. SDSS telescopes are located at Apache Point Observatory, funded by the Astrophysical Research Consortium and operated by New Mexico State University, and at Las Campanas Observatory, operated by the Carnegie Institution for Science. The SDSS web site is www.sdss.org.

SDSS is managed by the Astrophysical Research Consortium for the Participating Institutions of the SDSS Collaboration, including the Carnegie Institution for Science, Chilean National Time Allocation Committee (CNTAC) ratified researchers, Caltech, the Gotham Participation Group, Harvard University, Heidelberg University, The Flatiron Institute, The Johns Hopkins University, L’Ecole polytechnique fédérale de Lausanne (EPFL), Leibniz-Institut für Astrophysik Potsdam (AIP), Max-Planck-Institut für Astronomie (MPIA Heidelberg), Max-Planck-Institut für Extraterrestrische Physik (MPE), Nanjing University, National Astronomical Observatories of China (NAOC), New Mexico State University, The Ohio State University, Pennsylvania State University, Smithsonian Astrophysical Observatory, Space Telescope Science Institute (STScI), the Stellar Astrophysics Participation Group, Universidad Nacional Autónoma de México, University of Arizona, University of Colorado Boulder, University of Illinois at Urbana-Champaign, University of Toronto, University of Utah, University of Virginia, Yale University, and Yunnan University.

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.

Data Availability

The many ways to access SDSS-V data are described here: https://www.sdss.org/dr20/data_access/. The FIRE-2 public data releases (Wetzel et al. 2023) are available on the Globus platform.⁹

2MASS and ASPCAP Noise Model

In order to make mock 2MASS magnitudes and APOGEE ASPCAP outputs, we add noise to the corresponding simulated star properties and report a accompanying uncertainty. This procedure is described in § 5, but here we provide further details on the derived empirical noise models.

In Figure A1, we show the mean and standard deviation (shaded region) of the 2MASS catalog magnitude uncertainties as a function of magnitude for the $J$ (blue), $H$ (green), and $K_{\mathrm{s}}$ (pink) bands. These are calculated using 0.1 mag bins. In the bottom panel, we show the number of sources per bin as a function of magnitude. In general, we can see a jump in the average uncertainty at magnitudes $<5$ due to the 2MASS saturation limit. Between magnitudes of 5 to 10 we see the average uncertainty stays roughly constant at values between 0.02–0.06 mag. At magnitudes $>10$, we see the average uncertainty increase with increasing magnitude as expected.

Figure A2 we show the average and standard deviations of the APOGEE SNR for the DR19 GG data as a function $H$ magnitude. Similar to Figure A1, we use 0.1 mag bins. At $H<7$, the average SNR stays roughly constant at 250–300. These stars are brighter than the nominal APOGEE bright limit and are therefore targeted with a different strategy than the rest of the survey. At $H>7$, we see the average SNR decreases sharply with $H$ magnitude. This is because all targets receive the same 15-minute exposure, and the observing strategy calls for a single epoch for the majority of them. Only a subset of stars have repeat observations.

The average and standard deviations of the ASPCAP stellar parameters and abundances as a function of APOGEE SNR for the DR19 GG data are shown in Figure A3. For all the parameters and abundances we see the uncertainty decreases as the SNR increases. The $T_{\mathrm{eff}}$, $\log g$, [M/H], [C/H] and [N/H] are calculated simultaneously through $\chi^2$ spectral template fitting. This is why [M/H], [C/H], and [N/H] have weaker relationships with SNR and larger standard deviations than the other elemental abundances. The average $T_{\mathrm{eff}}$ uncertainty ranges from 20 K at the low SNR end to 1 K at high SNR. The average $\log g$ uncertainty ranges from 0.07 dex to 0.004 dex. [M/H] average uncertainty ranges from 0.001–0.02 dex, while

the [N/H] average uncertainty is much larger ranging from 0.003–0.4 dex. The [C/H] average uncertainty is in the middle with values between 0.0025–0.07 dex. Of the individually measured elemental abundances, Fe is the most precise, with uncertainties of 0.001 to 0.01 dex. This is expected, as the Fe has the most lines and the fit uses almost an order of magnitude more spectral windows than any other atomic species. The ASPCAP abundance uncertainty floor is 0.001 dex, which causes the flattening in the [Fe/H] uncertainty at SNR$>150$. For all the other elements (Mg, Si, S, O and Ca) the average uncertainty ranges from 0.025–0.002 dex. We note that these values encompass the noise driven statistical uncertainties but underestimate the total uncertainty because the systematic uncertainties are not taken into account ((Mészáros et al. 2025)).

Chemical Abundances

The abundances included in the mock catalog are the seven elements (C, N, O, Mg, Si, S, and Ca) that the FIRE-2 simulations track that are also included in the APOGEE ASPCAP catalogs. The [X/Fe] ratios for these elements are shown as a function of metallicity in Figure B1 for the nominal mock catalog that uses the MAP selection function. The larger uncertainties for [C/Fe] and [N/Fe] (see Figure A3) cause significant spread in the abundances for these elements. The reasons that the FIRE-2 abundances do not match the MW are complex, from limitations to the input yields and delay-time distributions for different channels of enrichment ([52] (Hopkins et al. 2018b; Hopkins et al. 2020; Muley et al. 2021; Hopkins et al. 2023), Gandhi et al. in prep) to the method of implementation of metal transport through diffusion (Escala et al. 2018) to the choice to initialize gas metallicity at a set value rather than model Pop III stars directly (Hopkins et al. 2018b, 2023, Meizani et al. in prep). To bring the mock abundances closer to the MW, we rescale the [X/H] values so that the median [X/Fe]=0 between -0.5 < [Fe/H] <+0.5. This is not necessarily true in the MW, but it is true for the solar neighborhood, and this change brings the abundances closer to MW trends in general. The rescaling is unlikely to be appropriate for stars with abundances very far from solar given the complexity of the modeled enrichment processes, but at least allows for a more informative comparison near the Solar circle.

Figure B2 is the same as Figure B1, but with the rescaled abundances. The value that each are rescaled by is shown in the lower left corner. Ca is the only element that gets scaled with a positive value in that the abundances needed to increase. This is also the smallest shift with a value of 0.047 dex. The largest shift is -0.476 dex applied to the O abundances. The average shift is -0.27 dex.

Milky Way Dust Comparison

For the nominal version of the catalog we prioritize self-consistency and calculate the dust extinction based on the metal-weighted hydrogen column density along the line of sight to each mock star (see § 5). However, we also provide a version where we utilize the MW 3D dust map (Green et al. 2019; Zucker et al. 2025), as described in § 5. Since the H and G-H magnitudes change with the dust treatment, applying the selection function yields slightly different stellar samples in each mock. We see the MW dust extinction keeps more stars within the survey observing strategy with a total of 118,256 stars in the final catalog, while the nominal catalog using the simulation-based dust estimate has 118,088 stars.

We compare the mocks with different dust treatment in Figure C1. In the left panel, we show the mean metallicity difference as a function of sky position between the MW dust and simulation dust catalogs. While some pixels change as much as 0.3 dex, there are generally no significant trends. The right panel shows the difference in number of stars between the catalog as a function of T$_{\rm eff}$ and log $g$. We can see the largest differences are mostly for cool giants. Because of the APOGEE magnitude limits, these stars are generally more distant than the dwarfs, including in the inner Galaxy where dust extinction is highest. Therefore, they are the most likely to be effected by changes in the dust extinction treatment. We also confirm there are no significant differences in the azimuthal structure between the mocks with different dust treatment.

!Difference in mean [Fe/H] and star counts between dust prescriptions bbox=0.13,0.08,0.86,0.49

FIRE-2 Mock Variations

In addition to the nominal m12i mock catalog, we calculate mock catalogs for two other FIRE-2 galaxies which have Gaia DR3 mocks presented in Nguyen et al. (2024). These are all calculated using the MAP estimate of the selection function and using the dust model from the simulations, rather than using a MW dust map.

One of the largest differences between the catalogs are the total number of stars. M12m has the highest stellar mass ($M_* = 10 \times10^{10} M_\odot$), then m12f ($M_* = 6.9 \times10^{10} M_\odot$) and m12i has the lowest ($M_* = 5.3 \times10^{10} M_\odot$) (McCluskey et al. 2024). As the selection function draws a fraction of stars from the simulated galaxy, the mock catalog size scales with the stellar mass of the simulated galaxy. The m12m mock has a total of 337,763 stars, while m12f has 199,284 stars and m12i has 118,088 stars.

In order to briefly explore the differences between these catalogs, we remake Figure 6, but with the other mock galaxies. The top row of Figure D1 is identical to the bottom row of Figure 6, with the middle and bottom rows showing the equivalent for m12f and m12m, respectively. Consistent with expectations for the highest stellar mass, m12m is also the most metal-rich and the average metallicity of stars in the mock catalog is -0.03 dex. Although, m12f has slightly more stellar mass, it is a little bit more metal-poor than m12i where the average metallicity of stars in the mock for m12f is -0.14 dex while the same for m12i is -0.13 dex. There are also a slight difference in the radial metallicity gradient as measured from the mocks. M12i has the strongest radial metallicity gradient (-0.0156 dex/kpc) which is most similar to the estimate from the SDSS-V DR19 GG observed sample (-0.0572 dex/kpc). Therefore, we adopt it as our nominal mock catalog in this work. M12m has the weakest radial metallicity gradient at -0.0098 dex/kpc.

In the right panels of Figure D1, we can see that the azimuthal variations in the mean metallicities varies between the FIRE-2 galaxies. M12m especially appears to show spiral structure-like variations, but further investigation that is beyond the scope of this work is needed to confirm this. Although each galaxy shows unique variations, in all cases the solar neighborhood is more metal-poor than expected from the radial metallicity gradient, while other azimuths at the same radius are more metal-rich. This is likely induced by the selection function, as supported by the absence of this feature in the uniform mock (Figure 6). In contrast, the MW data show a relatively metal-rich solar neighborhood. However, we cannot assume the selection function has the same effect on the MW data, since it depends on the underlying stellar populations. Further work modeling the MW’s azimuthal substructure while accounting for the selection function is required to fully understand this signal and its effect on MW measurements.

The sdssv_apogee_sf Package

The method of § 4 is implemented in the Python package sdssv_apogee_sf¹⁰ (Lucey 2026). In this paper we use it in a single configuration: the selection function of the DR19 GG sample as a function of sky position, $H$, and $G-H$, reported as the map and mean, and with a minimum of 20 2MASS sources per cell, which we then apply to the FIRE mock catalogs in § 5. The package itself is more general. In this appendix we summarize its capabilities (Appendix § 13.1) and give a short example of its use (Appendix E2).

Capabilities

The only input the user provides is a table of observed stars with sky positions, $H$ magnitudes, and, optionally, Gaia $G$ magnitudes. The 2MASS denominator, which is the expensive part of the calculation, ships with the package and contains $5.46 \times10^{8}$ 2MASS sources counted on a HEALPix $N_{\mathrm{side}} = 64$ grid in 0.5 mag bins over $-3 \le H < 18$ and $-2.5 \le G-H < 14.5$. Because the observed sample is supplied by the user, the selection function can be computed for any SDSS-V APOGEE program, data release, or subset defined by quality or stellar-parameter cuts, not only for the sample considered here. However, it is important to note that the underlying assumption that 2MASS is complete begins to break down at $H > 12$ towards the inner Galaxy so selection functions calculated with data beyond that cut are likely overestimating the completeness. The user can also choose (i) a selection function in $H$ alone or in both $H$ and $G-H$; (ii) the magnitude range and bin widths, in multiples of 0.5 mag; (iii) the minimum number of 2MASS sources per cell before pixels are merged, and the coarsest HEALPix resolution allowed; and (iv) whether each cell reports the posterior mean, $(N_{\mathrm{obs}}+1)/(N_{\mathrm{2MASS}}+2)$, or the posterior mode, $N_{\mathrm{obs}}/N_{\mathrm{2MASS}}$. Building a selection function for the full DR19 GG sample ($\approx2.4 \times10^5$ stars) takes seconds on a laptop, so it is practical to rebuild it for each subsample in an analysis. The result can be saved to and reloaded from disk, and its query method returns the selection probability for any set of sky positions and magnitudes. This is what is needed to forward model a survey, whether to draw an observed sample from a model or simulation, as we do for the mock catalogs, or to include the selection in the likelihood when fitting a model to the observed stars.

Example Usage

Below we compute the selection function for the DR19 GG sample from public data. A longer version of this example, including maps of the selection function, is provided as a tutorial notebook in the repository. The package is installed with

pip install git+https://github.com/mlucey/
    sdssv_apogee_sf

and the 2MASS denominator is included in the installation. The observed sample is taken from the DR19 ASPCAP summary file¹¹.

We select GG targets with the SDSS-V targeting flags and keep one entry per 2MASS source, matching how the denominator is counted.

import numpy as np
from astropy.table import Table, unique
from sdss_semaphore.targeting import TargetingFlags

dr19 = Table.read('astraAllStarASPCAP-0.6.0.fits.gz',
                  hdu=2)
flags = TargetingFlags(dr19['sdss5_target_flags'])
gg = flags.in_program('mwm_galactic') | flags.in_program
    ('mwm_gg')
obs = unique(dr19[gg], keys='obj', keep='first')    # one
    # row per 2MASS source

h = np.asarray(obs['h_mag'], dtype=float)
g = np.asarray(obs['g_mag'], dtype=float)
h[(h < -100) | (h > 20)] = np.nan                # remove sentinel values
g[(g <= 0) | (g > 22)] = np.nan
observed = Table({'ra': obs['ra'], 'dec': obs['dec'],
                  'h_mag': h, 'g_mag': g})

The selection function is then built from this table. The settings below match those used for the mock catalogs in § 5; omitting g_col would instead give a selection function in $H$ only.

from sdssv_apogee_sf import APOGESelectionFunction

sf = APOGESelectionFunction.from_observed(
    observed, ra_col='ra', dec_col='dec',

h_col='h_mag', g_col='g_mag',     # selection in H and
                                  G-H
estimator='map',                  # posterior mode,
  N_obs / N_2MASS
min_count=20, nside_min=8)        # pixel-merging
settings

sf.write('gg_sf_dr19.npz')        # reload with
    APOGEESelectionFunction.read

Finally, query evaluates the selection function at any sky positions and magnitudes (scalars or arrays), returning NaN outside its range:

from astropy.coordinates import SkyCoord

c = SkyCoord(l=30.0, b=5.0, unit='deg', frame='galactic'
             )
S = sf.query(c, H=10.0, GH=3.0)    # probability that
    such a star was observed

This paper has been typeset from a TeX/LaTeX file prepared by the author.