Cosmological constraints from the cross-correlation of DESI Luminous Red Galaxies with CMB lensing from Planck PR4 and ACT DR6

The LRG samples are a subset of the Legacy Survey (LS) DR9 [69], which is currently being used for targeting by the DESI spectroscopic redshift survey. The data are comprised of optical ($g,r,z$) imaging from the Beijing–Arizona Sky Survey (BASS [79]), Mayall $z$-band Legacy Survey (MzLS [69]), Dark Energy Camera Legacy Survey (DECaLS [69]) and the Dark Energy Survey (DES [80]), as well as four mid-infrared ($W1-W4$) bands from the Wide-field Infrared Survey Explorer (WISE [81]).

Photometric selection: We use the “Main LRG sample” whose photometric selection is described in detail in [68]. Briefly, the selection is composed of three color cuts on extinction-corrected $g,\,r,\,z$ and $W1$ magnitudes to mitigate stellar contamination, remove galaxies below $z\lesssim 0.4$ and produce a roughly constant number density out to $z\sim 0.8$, in addition to a cut in $z$-band fiber magnitude to produce a tail in the redshift distribution extending just beyond $z=1$. Due to differences in the photometry of the different imaging surveys, the selection in the “Northern” (BASS+MzLS) and “Southern” (DECaLS+DES) regions differ slightly in their implementation (see Eqs. 1 & 2 of [68]) in an effort to homogenize the LRG sample across the full footprint.

Photometric redshifts: The LRG sample has been subdivided into four photometric redshift bins [64] that we label $z_{1}$ through $z_{4}$. We note that the photometric redshifts presented in [64] differ slightly from those released in [82] and used in previous analyses (e.g. [49]). A detailed description of the photo-$z$ algorithm and its latest improvements is presented in Appendix B of [64]. Most notably, the spectroscopic training data have been updated to include redshifts from DESI’s Survey Validation [83] and Early Data Release [84]. We also note that the definition of the photo-$z$ bins in the Northern and Southern regions differ slightly in an effort to homogenize the sample (see Table 2 of [64]).

Image quality cuts: We use the LRG masks publicly available at this URL [64], which have been constructed to mitigate contamination from unwanted imaging artifacts, stars, and other undesired astrophysical contaminants (e.g. large galaxies, star clusters, and planetary nebulae). These include the BRIGHT, GALAXY and CLUSTER bit masks used in the LS DR9, in addition to a “veto” mask constructed from unWISE [85], WISE [81], Gaia/Tycho-2 [86, 87], and visually-inspected LS DR9 data to remove problematic regions not captured by the LS DR9 bit masks. For a more detailed description of the veto mask, see Section 2.4 and Appendix D of [68]. With this set of photometric selections and masks, only 0.3% of the LRGs have been classified as stars by DESI’s current spectroscopic data. In addition the following image quality cuts have been imposed: pixels with $E(B-V)\geq 0.15$ (determined with the extinction map from ref. [88], hereafter SFD) are masked; pixels where the stellar density is larger than 2500 deg^-2 (determined by the Gaia [89] stellar density map) are masked; and only pixels with at least two exposures in the $g$, $r$, and $z$ bands are included.

Imaging weights: Ref. [64] removes observational modulations to the LRG density maps following a standard procedure that we now describe. Weights are constructed under the assumption that variations in the local mean LRG density due to instrumental and galactic effects can be modeled as a random downsampling of the “true” LRG population. Under this assumption these density trends can be removed via division of the galaxy density map by a suitable weight map reflecting this downsampling factor. In practice, these weights are taken to be a linear combination of a set of templates whose coefficients are chosen via linear regression to remove trends in LRG density with the templates. The fiducial weights provided in [64] use seven templates: depth and seeing in the three optical bands ($g,r,z$) and $E(B-V)$ [88]. Ref. [64] also provides a set of weights that do not use an $E(B-V)$ template, which we use in §6 to assess the importance of Galactic extinction and CIB contamination in the SFD map as a potential systematic. On the scales relevant for our analysis the systematic weights impact the galaxy auto-spectra by less than a percent (Fig. 10).

Map making: We use the LRG density contrast maps publicly available at this URL. Here we briefly summarize the map making performed by ref. [64]. The LRGs were binned into HEALPix [90] pixels with nSide = 2048. Likewise a pixelized map of weighted randoms was constructed using the previously described imaging weights. The masks and image quality cuts (discussed above) were then applied to both the LRG and weighted random maps. In an effort to mitigate the presence of many small holes in the LRG mask, pixels whose (weighted) random density was less than $1/5$ of the mean (weighted) random density over the footprint satisfying the image quality cuts were additionally masked. We follow [49, 64] and choose to not apodize the LRG masks. The LRG density contrast is defined as the LRGs divided by the weighted randoms, mean subtracted and normalized to the mean.

Redshift distribution: The redshift distributions of the four photometric LRG samples have been directly calibrated in [64] using 2.3 million LRG redshifts from DESI’s Survey Validation [83] and first year (Y1) of observations. These redshift distributions account for the aforementioned imaging weights, masks and image quality cuts, as well as the spectroscopic systematics weights for DESI Y1 observations. In addition, each spectroscopic galaxy is weighted by its inverse success rate (as defined in section 4.4 of [68]) to mitigate the impact of redshift failures²²2The spectroscopic sample used to calibrate these redshift distributions has a redshift failure rate of 1.3% and catastrophic failure rate (defined as deviating by more than 1000 km/s from the true redshift) of 0.2%.. The normalized redshift distributions of the four LRG samples, averaged over the full imaging footprint (determined by area weighting the redshift distributions calibrated on the North, DECaLS and DES footprints), is shown in the top right panel of Fig. 2. In the bottom three panels we show the variations of the redshift distribution in the three imaging regions. We propagate these variations to $\Delta C^{gg}_{\ell}$ in §6, and examine their impact on cosmological constraints in §7.

Number count slope: We make use of the number count slopes $s_{\mu}$ measured in [64] for the “combined” sample, which are listed in Table 1. The logarithmic change to the total number of LRGs for a small change in magnitude due to magnification is defined as $s_{\mu}=d\log_{10}N/dm$. These are measured via finite difference by shifting the ($g,\,r,\,z,\,W1$) LRG magnitudes by $\delta m=\pm 0.01$, computing the appropriate $z$-band fiber magnitude shift $\delta m_{\text{fiber}}$ for each galaxy (see Appendix C of [64] for details regarding the $\delta m_{\text{fiber}}$ calculation), reapplying the photometric selection, recomputing the photometric redshifts with the shifted magnitudes and galaxy sizes, and finally rebinning the LRGs into photometric redshift bins.

Changes to the LRG sample since White et al. [49]: These are listed in Appendix A of [64] and briefly summarized here. The LRG photo-$z$’s have been improved using more training data, the redshift distributions are estimated with a larger number of spectra, the zero point offset for DEC $<-29.25\mathrm{\SIUnitSymbolDegree}$ in the $r$ and $z$ bands has been corrected using the offset map from [91], and the imaging systematics no longer include a template for WISE $W1$ depth (due to a lack of an observed density trend with $W1$ depth and a discrepancy in the $W1$ depths in the data and randoms). With these changes, the LRG auto-correlation differs by at most 2% from the measurements used in [49], which is primarily driven by the improved photo-$z$’s and hence narrower redshift distributions.

We consider both Planck 2018 [44] (PR3) and PR4 [45] CMB lensing measurements in this work. For both releases we low pass filter the CMB lensing convergence multipoles³³3As in [49] we multiply the reconstructed CMB lensing multipoles $\hat{\kappa}_{\ell m}$ by $\text{Exp}[-(\ell/2500)^{6}]$. This filter impacts the cross-correlation measurements by less than $0.02\%$ for $\ell<600$, thus we neglect this filtering in our forward model. We do however include the filter in the fiducial CMB lensing noise curves used for our covariance. $(\hat{\kappa}_{\ell m})$ and rotate them from Galactic to equatorial coordinates using healpy’s rotate_alm routine before constructing nSide = 2048 HEALPix convergence maps via alm2map. Following [49] we apodize the (binary) CMB lensing masks with a 30 arcmin “C2” filter using the NaMaster \faGithub [92] routine mask_apodization before rotating the masks from galactic to equatorial coordinates with the healpy routine rotate_map_alms.

Public Release 3 (PR3): We reanalyze the cross-correlation with PR3 to make direct contact with a previous analysis [49] (which found $S_{8}=0.725\pm 0.030$) using a similar sample of LRGs. In §7.2 we quantify the impact of the improved LRG sample, new binning scheme, and updated modeling choices on these results. The Planck PR3 lensing measurement [44] reconstructs the lensing convergence from a linear combination of minimum variance quadratic estimators [93] whose inputs are Wiener-filtered and inverse-variance-weighted SMICA [94] temperature and polarization multipoles with $100\leq\ell\leq 2048$. Specifically, we use the minimum variance lensing reconstruction⁴⁴4COM_Lensing_4096_R3.00 from the Planck legacy archive. along with the provided CMB lensing mask and effective reconstruction noise curve $N^{\kappa\kappa}_{\ell}$.

Public Release 4 (PR4): The Planck PR4 lensing map [45] is reconstructed using a more optimal Global Minimum Variance estimator [95] from the latest NPIPE processing pipeline, which uses $\sim 8\%$ additional measurement time and a more optimal (anisotropic) filtering scheme (using the same scales $100\leq\ell\leq 2048$ as PR3). The resulting CMB lensing map is signal dominated (per mode) out to⁵⁵5We use context rather than different symbols (e.g. with $L$ instead of $\ell$, as in the companion paper [66]) to differentiate CMB lensing and primary CMB multipoles, since we use $L$ in §3 to denote bandpowers. From §3 onward $\ell$ always refers to a CMB lensing multipole. $\ell\sim 70$, compared with $\ell\sim 40$ with the previous PR3 map, resulting in an overall 20% improvement to the signal-to-noise ratio and a $42\sigma$ detection significance of the CMB lensing auto-correlation. Specifically, we use the reconstructed convergence multipoles, mask and reconstruction noise curve provided on NERSC⁶⁶6See github.com/carronj/planck_PR4_lensing for more information. Specifically, we use PR4_klm_dat_p.fits, mask.fits.gz and PR4_nlkk_p.dat..

Our headline analysis uses the latest CMB lensing measurement⁷⁷7We use the baseline reconstructed multipoles ($\hat{\kappa}_{\ell m}$, kappa_alm_data_act_dr6_lensing_v1_baseline.fits) and corresponding mask (mask_act_dr6_lensing_v1_healpix_nside_4096_baseline.fits) that are available at this URL. We downgrade the mask from an nSide of 4096 to 2048 using healpy’s ud_grade routine. (which detected the CMB lensing auto-correlation with a 43$\sigma$ significance) from the Atacama Cosmology Telescope (ACT) DR6 [46, 70, 63], which we low pass filter⁸⁸8We set $\hat{\kappa}_{\ell m}=0$ for $\ell>3000$, and leave modes with $\ell\leq 3000$ untouched. before constructing a $\kappa$ map using healpy’s alm2map routine with nSide = 2048. The data used to construct this map consists of nighttime observations made through 2021 at 98 and 150 GHz. The CMB lensing convergence is reconstructed from a (nearly optimal) linear combination of temperature- and polarization-based quadratic estimators, which in turn have been separated into a linear combination of estimators measured with disjoint data splits [96] to mitigate noise biases. All estimators use the scales $600<\ell<3000$ to mitigate contamination from Galactic emission and extragalactic foregrounds. To further mitigate extragalactic contamination in the temperature-based reconstruction, a profile-hardened estimator [97, 98, 99] has been adopted, which has been shown to efficiently suppress tSZ- and CIB-induced biases without requiring a finely-tuned model for the mean tSZ intensity profile [100].

The models discussed in the following subsections are built on the Aemulus $\nu$ simulations [65], which are a set of 150 $N$-body simulations spanning a broad range of cosmologies (in particular, $1.1<10^{9}A_{s}<3.1$ and $0.08<\omega_{c}<0.16$) including both massive neutrinos and variations to the dark energy equation of state. The simulations treat the CDM+baryon ($cb$) and neutrino fields as two sets of $1400^{3}$ collisionless particles that are evolved in a $1.05h^{-1}\text{Gpc}$ box, yielding a $cb$ mass resolution of $3.51\times 10^{10}(\Omega_{cb}/0.3)h^{-1}M_{\odot}$. The initial displacements ($z_{\text{ini}}=12$) of the $cb$ field are calculated using third-order Lagrangian Perturbation Theory (as implemented in monofonic [103, 104]) from an initial Gaussian density field whose power spectrum is taken to be $P_{cb}(k,z_{\text{ini}})=D_{cb}(k,z_{\rm ini})^{2}P_{cb,\text{lin}}(k,z=0)$, where $P_{cb,\text{lin}}(k,z=0)$ is computed with CLASS [105] and the scale-dependent growth factor is computed using first-order Newtonian fluid approximation as implemented in zwindstroom [106], while the initial conditions for the neutrinos are set with fastDF [107]. The $cb$ and neutrino particles are evolved from $z_{\text{ini}}=12$ to $z=0$ using a modified version of Gadget-3 [108] that enables modifications and includes relativistic corrections to the background evolution $H(z)$.

We adopt a Hybrid Effective Field Theory (HEFT) model [109, 67, 110, 111, 112, 113] in our fiducial analysis. Within this framework the “initial” (Lagrangian) galaxy density is modeled as a biased tracer of the Lagrangian CDM+baryon density. Proto-galaxies at Lagrangian position $\bm{q}$ are advected to their later (real-space) positions $\bm{x}(t)=\bm{q}+\bm{\Psi}_{cb}(\bm{q},t)$ following the (non-linear) displacement vector $\bm{\Psi}_{cb}(\bm{q},t)$ of the CDM+baryon field [114, 115, 116, 117, 118, 119, 120, 56] that is measured from simulations. From number conservation, the (real-space) late-time galaxy density contrast is given by

where the weight function $F(\bm{q})$ characterizes the initial galaxy density fluctuations and $\delta^{D}$ is the Dirac delta function. Here we take $F(\bm{q})$ to be a perturbative expansion in scalar combinations of the tidal $(\partial_{i}\partial_{j}\Phi)$ and velocity $(\partial_{i}v_{j})$ tensors to next-to-leading order in derivatives and perturbations [121]

where $s^{2}_{cb}=\sum_{ij}s^{ij}_{cb}s^{ij}_{cb}$ is the square of the shear field $s^{ij}_{cb}=(\partial_{i}\partial_{j}/\partial^{2}-\delta^{K}_{ij}/3)\delta_{cb}$ and $\mathcal{E}(\bm{q})$ is a small-scale stochastic component. Within this approximation the late-time galaxy density contrast takes the form

where $X\in\{1,2,s\}$, $F_{X}(\bm{k})$ corresponds to the respective term (post-advection and Fourier transformed) in the weight function, and we have approximated the post-advection $\nabla^{2}\delta_{cb}$ term as $-k^{2}\delta_{cb}(\bm{k})$. The galaxy power spectrum is then

where $P_{cb}$ is the power spectrum of the CDM+baryon field, $P_{cb,X}$ is the cross-correlation of the CDM+baryon field with $F_{X}$, $P_{XY}$ is the cross-spectrum of $F_{X}$ and $F_{Y}$, SN${}^{\text{3D}}$ is a shot noise contribution, we introduced the “counterterm” $\alpha_{a}$ for the auto-correlation, and we have neglected $\mathcal{O}(k^{2}\delta^{3}_{cb})$ contributions. Likewise the cross-correlation of the galaxies with the matter is approximately

where $P_{cb,m}$ is the cross spectrum of the CDM+baryon field with the total matter field, $P_{m,X}$ is the cross correlation of the total matter field with $F_{X}$, and we have defined the counterterm $\alpha_{x}$ for the cross-correlation.

The power spectra ($P_{cb},P_{cb,X}$, $P_{XY}$, $P_{cb,m}$, and $P_{m,X}$) appearing in Eqs. (4.4) and (4.5) are calculated in the Aemulus $\nu$ simulations described in §4.1. The sample variance of these measurements on perturbative scales has been suppressed using the Zel’dovich density field (implemented using the ZeNBu \faGithub and velocileptors \faGithub [122] codes) as a control variate¹⁴¹⁴14Specifically, [65] estimates a power spectrum (e.g. $P_{X,Y}$) using $\hat{P}^{\text{CV}}_{X,Y}\equiv\hat{P}^{\text{Aemulus}}_{X,Y}-\beta(\hat{P}^{\text{LPT}}_{X,Y}-P^{\text{LPT}}_{X,Y})$, where $\hat{P}^{\text{Aemulus}}_{X,Y}$ is the measured power spectrum from the Aemulus $\nu$ simulation, $\hat{P}^{\text{LPT}}_{X,Y}$ is the measured power spectrum of the (Zel’dovich) LPT evolution of the same initial conditions used in the simulation and $P^{\text{LPT}}_{X,Y}$ is the analytic prediction for the ensemble average (which is known to machine precision). The variance of $\hat{P}^{\text{CV}}_{X,Y}$ is (optimally) suppressed by a factor of $1-\rho^{2}_{\text{Aemulus},\text{LPT}}$ relative to the variance of $\hat{P}^{\text{Aemulus}}_{X,Y}$ if one chooses $\beta=\text{Cov}[\hat{P}^{\text{Aemulus}}_{X,Y},\hat{P}^{\text{LPT}}_{X,Y}]/\text{Var}[\hat{P}^{\text{LPT}}_{X,Y}]$, where $\rho_{\text{Aemulus},\text{LPT}}$ is the correlation coefficient of the Aemulus and LPT measurements.. We use the Aemulus $\nu$ emulator \faGithub in our calculations, which has been shown to be accurate to within $0.25\%$ for the scales $(k\lesssim 0.6\,h\,{\rm Mpc}^{-1})$ and redshifts $(z\lesssim 1)$ relevant for our analysis [65].

While we have introduced two counterterms ($\alpha_{x}$ and $\alpha_{a}$) in Eqs. (4.4) and (4.5), these parameters are not independent in the limit that the bias expansion (Eq. 4.2) accurately describes the field-level galaxy distribution and the simulated $\delta_{m}(\bm{x})$ accurately describes the true matter distribution¹⁵¹⁵15This is analogous to multi-tracer studies where the coefficients of counterterm-like contributions to the cross-correlation of two biased tracers can be expressed in terms of counterterm coefficients appearing in the two individual auto-spectra (shown explicitly in Appendix C of [123]).. The former approximation is valid up to baryonic feedback and the finite mass resolution of the simulations. Both of these effects are small on the scales of interest and have the qualitative impact of damping small-scale power that we approximate by introducing an effective “derivative bias” in the matter distribution: $\delta_{m}(\bm{k})\to\big{(}1-\tilde{\epsilon}k^{2}/4\big{)}\delta_{m}(\bm{k})$. Baryonic feedback scenarios capable of resolving the $S_{8}$-tension in galaxy-lensing measurements (Fig 6 of ref. [124], C-OWLS AGN with $\log_{10}(\Delta T_{\text{heat}}/{\rm K})=8.7$ at $z=0$) require a $\sim 10\%$ suppression to the matter power spectrum at $k\simeq 0.4$ $h\,$Mpc^-1 corresponding to $\tilde{\epsilon}\simeq 1.3$ $h^{-2}$Mpc², while the effective grid scale in the Aemulus $\nu$ simulations is $\approx 1\,h^{-1}$Mpc corresponding to $\tilde{\epsilon}\simeq 1$ $h^{-2}$Mpc². With the inclusion of the effective derivative bias in the matter distribution the $\mathcal{O}(k^{2})$ components of the galaxy-auto and cross-correlation with matter are

where we neglect $\mathcal{O}(k^{2}\delta^{3})$ contributions (in particular, $F_{1}(\bm{k})=\delta_{cb}(\bm{k})$ under this approximation), we introduced the Eulerian bias $b_{1}^{E}\equiv 1+b_{1}^{L}$, the symbol $\supset$ denotes a subset of the terms appearing in the power spectra predictions (in this case, the $k^{2}$ terms), and primed brackets correspond to an ensemble average with the momentum-conserving delta function removed, e.g. $\langle\delta(\bm{k}_{1})\delta(\bm{k}_{2})\rangle\equiv(2\pi)^{3}\delta^{D}(\bm{k}_{1}+\bm{k}_{2})\langle\delta(\bm{k}_{1})\delta(\bm{k}_{2})\rangle^{\prime}$. Comparing Eq. (4.6) to Eqs. (4.4) and (4.5) we see that $\alpha_{a}=b^{E}_{1}b^{L}_{\nabla^{2}}$ and $2\alpha_{x}=b_{1}^{E}\tilde{\epsilon}+b^{L}_{\nabla^{2}}$, which in turn implies

where $\epsilon\equiv b_{1}^{E}\tilde{\epsilon}/2$ ($b_{1}^{E}/2\simeq 1$ for the LRGs considered here) has a “typical” value of $\simeq 1\,h^{-2}$Mpc², which we marginalize over with an informative Gaussian prior centered at zero with width $\simeq 2\,h^{-2}$Mpc² in our fiducial analysis (see §5.2). We note that this prior is sufficiently broad to permit an enhancement of small-scale power (within one prior $\sigma$) rather than a suppression, even when taking into account the suppression already present in the model prediction arising from the effective grid scale in the Aemulus $\nu$ simulations.

Non-linear corrections to the matter power spectrum are $\leq 1.5\%$ for $k\leq 0.1\,h\,{\rm Mpc}^{-1}$ and $z>0.3$. Provided that the length scale associated with the bias expansion (e.g. the Lagrangian radius of LRG host halos) is smaller than that associated with non-linear (gravitational) dynamics¹⁶¹⁶16We note that this approximation is implicitly made when fitting to higher $k$ with HEFT than with “pure” perturbation theory. A caveat to this claim is that some higher-order contributions, such as $\langle\delta^{2}\delta^{2}\rangle$, are non-zero in the $k\to 0$ limit. However, these contributions are scale-independent in this limit and thus degenerate with shot noise., one can accurately model the galaxy distribution using linear theory on these scales, which we consider as an alternative to our fiducial HEFT model in §7.

Rather than using pure linear theory for the galaxy cross- and auto-correlations, we choose to use the HEFT prediction with $b_{2}^{L}=b_{s}^{L}=0$, and adjust our scale cuts such that $\ell_{\rm max}(z)\simeq\chi(z)\times(0.1\,h\,{\rm Mpc}^{-1})$ (see §5.1 and Table 2). This is done purely for convenience (e.g. the emulator is faster than a CLASS evaluation), since swapping the HEFT predictions for linear theory on these scales would have a negligible impact on our constraints. To approximately account for the residual systematic error in $P_{gg}$ and $P_{gm}$ from neglecting higher-order contributions on these scales we choose to marginalize over $\alpha_{a}$ and $\alpha_{x}$ independently with an informative prior chosen such that a 1.5% correction is allowed at $k=0.1\,h\,{\rm Mpc}^{-1}$ at $1\sigma$ (see §5.2 and Table 3). As for our fiducial model, we use the non-linear $P_{mm}$ prediction when calculating magnification contributions (§4.5).

As yet another alternative to HEFT and linear theory, we consider a “model independent” (or “fixed shape”) approach to constrain $\sigma_{8}(z)$. Under this approach we use the same scale cuts as for our linear theory analysis, fix the cosmology, and allow the linear bias for the galaxy auto ($b^{L,a}_{1}$) and cross ($b^{L,x}_{1}$) to differ by a factor

Within linear theory one constrains $\sigma_{8}$ schematically through the ratio $C^{\kappa g}/\sqrt{C^{gg}}\sim\alpha_{8}\sigma_{8,\text{ fid}}$, where $\sigma_{8,\text{ fid}}$ is the value of $\sigma_{8}$ for the fiducial cosmology. In §7.6 we consider fitting to each redshift bin independently, and interpret our constraint on $\alpha_{8}\,\sigma_{8,\text{ fid}}(z_{i})$ as a model independent measurement of $\sigma_{8}(z_{i})$.

We use the Kaiser-Limber approximation [125, 126] (hereafter Limber approximation) to model the angular galaxy auto- and cross-correlation:

where $k=(\ell+1/2)/\chi$ [127] and $z$ is implicitly a function of $\chi$. The projection kernels for the galaxies, magnification contribution, and CMB lensing are [128, 129, 130]

where $\phi(z)\propto dN/dz$ is the normalized ($\int dz\phi(z)=1$) redshift distribution, $\chi_{*}$ ($z_{*}$) is the comoving distance (redshift) to the surface of last scattering, and $s_{\mu}$ is the number count slope. Due to the narrow redshift distributions of the LRG bins, the magnification contributions to the galaxy auto-spectra amount to less than a percent of the total signal, while for the cross-correlation with CMB lensing the magnification contribution is at most $\simeq 10\%$ for the highest redshift bin (and smaller for the remaining bins). We compute background quantities (i.e $\Omega_{m}$, $\chi_{*}$, $H(z)$, and $\chi(z)$) relevant for Limber integration from our sampled cosmological parameters (§5.2 and Table 3) using CLASS [105].

In practice, we neglect the evolution of $P_{gg}$ and $P_{gm}$ in the Limber integrals, and evaluate them at the effective redshift

This approximation, which we expect to be accurate to subpercent precision for our LRG sample (see Fig. 4) due to its narrow redshift distributions, has the advantage of making our analysis largely agnostic to the assumed redshift evolution of the galaxy nuisance parameters within each redshift bin. Instead, we treat e.g. $b_{1}^{L}(z_{i})$ as a single number evaluated at the effective redshift for the $i$’th redshift bin.

Deviations from the Limber approximation are primarily relevant on large scales where linear theory is a good approximation. Within this approximation and neglecting magnification bias, the spherical harmonic coefficients $g_{\ell m}$ of the projected galaxy density contrast can be split in two pieces $g^{\text{real}}_{\ell m}+g^{\text{RSD}}_{\ell m}$ where [131, 132]

$\delta_{g}(\bm{k})$ is the real-space 3D galaxy density contrast, $\beta(z)=f(z)/b_{1}^{E}(z)$ where $f(z)$ is the linear growth rate and $j^{\prime}_{\ell}(x)=\partial_{x}j_{\ell}(x)$ is the derivative of the spherical Bessel function of the first kind. The cross-correlation $C^{XY}_{\ell}\equiv\langle g^{X}_{\ell m}g^{Y*}_{\ell m}\rangle$ is given by the “full integral”

where $\langle\delta_{g}(\bm{k},\chi)\delta_{g}(\bm{k}^{\prime},\chi^{\prime})\rangle\equiv(2\pi)^{3}\delta^{D}(\bm{k}+\bm{k}^{\prime})P_{gg}(k;\chi,\chi^{\prime})$.

In Fig. 4 we show the full integral for both the auto- and cross-correlations assuming linear theory: $P_{gg}(k;\chi,\chi^{\prime})=D(z)b(z)D(z^{\prime})b(z^{\prime})P_{0}(k)$ where $D(z)$ is the linear growth factor normalized to $D(0)=1$ and $P_{0}(k)$ is the $cb$ linear power spectrum at $z=0$ (we are ignoring shot noise). We take $b(z)$ to be a linear interpolation of the best-fit values listed in Table 1 of [49]. We also show the Limber approximations, including the effective redshift approximation. For $\ell\lesssim 75$, the largest beyond Limber effect is due to RSD, while above $\ell\gtrsim 75$ the largest effect is due to the effective redshift approximation, which is still $<1\%$ and neglected in our analysis.

The LRGs considered here have been placed on a HEALPix grid¹⁷¹⁷17See [133] for an alternative approach that bypasses placing galaxies on a discrete grid, and hence the pixel window function and aliasing.. This pixelization has the approximate effect of convolving the galaxies with the pixel window function, but in addition gives rise to aliasing of power [134, 135, 133, 136]. The combination of these effects can be well approximated by taking $C^{\kappa g}_{\ell}\to w_{\ell}C^{\kappa g}_{\ell}$ and $C^{gg}_{\ell}\to w_{\ell}^{2}(C^{gg}_{\ell}-\text{SN}^{\rm 2D})+\text{SN}^{\rm 2D}$, where $w_{\ell}$ is the pixel window function (as determined by healpy’s pixwin for the appropriate nSide). Previous analyses (e.g. [49, 50]) have approximately corrected for the pixel window function at the data level, which requires assuming a fiducial value for the projected shot noise. This is somewhat undesirable (but still reasonable) given that one fits for the shot noise in practice. Here we opt to forward model the impact of pixelization using the aforementioned approximation. For nSide = 2048 the window function differs from one by less than $0.4\%$ for $\ell<600$, and as a result this improved treatment of the pixel window function has a minuscule impact on our final constraints.

Throughout §7 we use (a subset of) the bandpowers discussed in §3 spanning the range $20\leq\ell<600$. The fiducial scale cuts used for each our modeling choices (sections 4.2, 4.3 and 4.4) are summarized in Table 2. As discussed in §4.5, we adopt the Limber approximation and neglect the impact of redshift space distortions in our model, which thus acts as a systematic on large scales. For both the galaxy auto- and cross-correlation with Planck, we adopt an $\ell_{\rm min}$ ($20$ for $\kappa g$ and $79$ for $gg$) where the “beyond Limber” corrections are less than a percent (see Fig. 4). On large scales the ACT DR6 lensing reconstruction is contaminated by a mean field contribution that greatly exceeds the CMB lensing signal for $\ell\lesssim 30$ (e.g. Fig. 9 of [46]), making it difficult accurately correct for in this regime. Following suit with the DR6 auto-correlation measurement [46] (which adopted $\ell_{\rm min}=40$) we choose to adopt a larger $\ell_{\rm min}=44$ for the ACT cross-correlation to mitigate spurious correlations with the mean field.

At high $\ell$ our scale cuts vary with the model being used. Previous work [67] has found that our fiducial HEFT model (§4.2) is accurate to subpercent precision for $k\leq 0.6\,h\,{\rm Mpc}^{-1}$ when fitting to dark matter halos with the same characteristic masses $(10^{12.5}<\log_{10}M_{h}/M_{\odot}<10^{13})$ and redshifts as expected from recent HOD fits [148] for the host halos of our LRG sample. This motivates $\ell_{\rm max}\simeq\chi(z_{\rm min})\times(0.6\,h\,{\rm Mpc}^{-1})$ where $z_{\rm min}$ is the lower edge of a given redshift bin. In particular, for the first redshift bin $z_{\rm min}\sim 0.35$ and thus $\ell_{\rm max}\simeq 600$. In principle one could reliably extend to a higher $\ell_{\rm max}$ for the higher redshift samples, however in practice we expect limited gains from higher $\ell$ due to shot noise, lensing reconstruction noise, and degeneracies with galaxy bias. For simplicity we adopt a redshift-independent $\ell_{\rm max}=600$ in our fiducial HEFT analysis for both the galaxy auto- and cross-correlation.

As discussed in §4.3 we tune our linear theory scale cuts such that $\ell_{\max}\simeq\chi(z_{{\rm eff},i})\times(0.1\,h\,{\rm Mpc}^{-1})$, where the differences between the linear and non-linear matter power spectrum are $\leq 1.5\%$. For the lowest redshift bin applying this scale cut removes all but one bandpower for the galaxy auto-correlation ($79<\ell<124$), which we opt to discard when quoting linear theory constraints. For bins $z_{2}$, $z_{3}$ and $z_{4}$ we adopt $\ell_{\max}=178,\,243,\,243$ corresponding to $k_{\rm max}\simeq 0.11$, $0.13$, $0.11\,h\,{\rm Mpc}^{-1}$. We use the same $\ell_{\rm max}$’s used in our linear theory fits for our model independent constraints. Unlike for linear theory, we consider a model independent $\sigma_{8}(z_{1})$ measurement in the lowest redshift bin, for which we take $\ell_{\max}=178$. In §5.5 we verify that all of our models recover unbiased results on simulations when adopting these scale cuts.

Our data (Fig. 3) are particularly powerful at extracting the relative amplitude between $C^{\kappa g}_{\ell}$ and $\sqrt{C^{gg}_{\ell}}$ on large scales, roughly corresponding to $S_{8}$. Since our goal is to assess the consistency of a low redshift structure growth measurement with that predicted within $\Lambda$CDM conditioned on primary CMB data, we fix the remaining relevant $\Lambda$CDM parameters to their Planck 2018 [149] mean values¹⁹¹⁹19We use the values from the TT,TE,EE+lowE column in Table 2 of [149], which contains no low redshift information modulo lensing and other secondary anisotropies.. Specifically, we fix the baryon abundance $\Omega_{b}h^{2}=0.02236$, spectral index $n_{s}=0.9649$, $\Omega_{m}h^{3}=0.09633$ as a proxy for the angular acoustic scale $\theta_{*}$ [150, 149] and additionally fix the sum of the neutrino masses to $\sum m_{\nu}=0.06$ eV. We directly sample the (log) primordial power-spectrum amplitude $\ln(10^{10}A_{s})$ and dark matter abundance $\Omega_{c}h^{2}$ with uniform $\mathcal{U}(2,4)$ and $\mathcal{U}(0.08,0.16)$ priors respectively. When quoting model-independent constraints these parameters are fixed to $3.045$ and $0.1202$ respectively [149] and we vary $\alpha_{8}$ (Eq. 4.8) with a $\mathcal{U}(0.5,1.5)$ prior. We note that the $\ln(10^{10}A_{s})$ prior extends beyond the range of the Aemulus $\nu$ simulations (§4.1), however, our data are sufficiently constraining such that for nearly all²⁰²⁰20The one exception is when analyzing $z_{1}$ independently without a BAO prior. Restricting the $\ln(10^{10}A_{s})$ prior to $\mathcal{U}(2.4,3.4)$ (corresponding to the Aemulus $\nu$ range) truncates the $\sigma_{8}-\Omega_{m}$ contour shown in Fig. 14, however, the mean $S_{8}$ is largely unaffected ($\Delta S_{8}=0.003$, corresponding to a $<\sigma/10$ shift). of the results presented in §7 the $2\sigma$ credible intervals lie within the Aemulus $\nu$ range.