The Atacama Cosmology Telescope: measurement and analysis of 1D beams for DR4

We describe the measurement and treatment of the telescope beams for the Atacama Cosmology Telescope's fourth data release, DR4. Observations of Uranus are used to measure the central portion (<12') of the beams to roughly -40 dB of the peak. Such planet maps in intensity are used to construct azimuthally averaged beam profiles, which are fit with a physically motivated model before being transformed into Fourier space. We investigate and quantify a number of percent-level corrections to the beams, all of which are important for precision cosmology. Uranus maps in polarization are used to measure the temperature-to-polarization leakage in the main part of the beams, which is ≲ 1% (2.5%) at 150 GHz (98 GHz). The beams also have polarized sidelobes, which are measured with observations of Saturn and deprojected from the ACT time-ordered data. Notable changes relative to past ACT beam analyses include an improved subtraction of the atmospheric effects from Uranus calibration maps, incorporation of a scattering term in the beam profile model, and refinements to the beam model uncertainties and the main temperature-to-polarization leakage terms in the ACT power spectrum analysis.


INTRODUCTION
The Atacama Cosmology Telescope (ACT) is a 6 m off-axis Gregorian telescope located at an altitude of 5190 m in the Atacama Desert of northern Chile.It is designed for millimeter wavelength observations of the cosmic microwave background (CMB) at arcminute resolution.The telescope and receiver are described in Fowler et al. (2007) and Thornton et al. (2016) respectively.This paper describes the treatment of the ACT beams for its fourth data release, referred to as DR4.This data release includes temperature and polarization data collected by ACT between 2013 and 2016, covering seven regions of the sky spanning roughly 18,000 deg 2 , in frequency bands centered on 98 and 150 GHz (Thornton et al. 2016).The DR4 data were obtained using two 150 GHz detector arrays, PA1 and PA2, and one dichroic detector array, PA3, operating at 98 and 150 GHz.The three array positions are not identical optically.Each year of data is referred to as a season (S13-S16) and was analyzed separately.The DR4 data products and some of their analyses are presented in Choi et al. (2020), Aiola et al. (2020), Darwish et al. (2020), Han et al. (2021), Madhavacheril et al. (2020), Namikawa et al. (2020), andMallaby-Kay et al. (2021).
The beams of the telescope determine the instrument's response to different scales on the sky.Quantifying these beams and the uncertainties on the beam measurements is one of the most challenging aspects of characterizing the instrument.An incorrect calibration of the beams would bias virtually all of the ACT science analyses, including the precision measurement of the CMB power spectrum.
The previous paper that focused on ACT beams, Hasselfield et al. (2013), dates from DR2, which included data through 2010.The DR3 analysis papers (Louis et al. 2017;Naess et al. 2014) referenced Hasselfield et al. (2013) and described the small changes in methodology since then.In recent years, there have been several modifications to our beam pipeline that merit further discussion.This paper provides a comprehensive guide to the analysis of the ACT beams for DR4.
The paper is organized as follows.In §2 we describe the planet observations used to measure the main portion of the ACT beams.In §3 we explain how these observations are made into maps.In §4 we describe the steps of the beam pipeline from planet maps to a model of the ACT beams and their covariance.§5 accounts for the polarization component of the beams.Finally, in §6 we discuss the publicly released beam products, assumptions made in the analysis, and future directions for ACT beam analyses.Part of the information in this paper overlaps with that in previous ACT papers, including Aiola et al. (2020), Choi et al. (2020), and Hasselfield et al. (2013), and is included here for clarity and completeness.

OBSERVATIONS
Observations of sufficiently bright point-like objects are required in order to properly characterize the shape of the beams out to large angles.For ACT, planets are the best candidates for this purpose, though not all equally so.The proximity of Mercury and Venus to the Sun make them difficult to observe at night; Mars' temperature is not sufficiently constant (see, for example, Wright 1976), so it intermittently becomes too bright; Jupiter's apparent size is too large, compared to the instrument's angular resolution, to be considered a point source (and it is too bright, as is Saturn); Saturn, though suitable for sidelobe analysis, is bright enough to induce a non-linear detector response near peak amplitude, resulting in signal suppression within the main lobe of the beam;1 Neptune's small angular diameter relative to the beamwidth dilutes its brightness, reducing signal-to-noise.We have thus chosen to base the majority of the ACT beam analysis on observations of Uranus, which achieves adequate signal-to-noise without exceeding the dynamic range of the instrument and can be approximately treated as a point source. 2As described in §5.2, Saturn is used to study the beam's polarized sidelobes.
Planet observations are done by briefly interrupting the CMB survey scans, changing the azimuth at which the telescope is pointing, and scanning again at a fixed elevation.Uranus is observed roughly once per night during the observing season, but only the higher quality observations are used to characterize the beam.Before making maps of the observations, we apply the same data selection criteria as for the main maps used for cosmological analysis.For example, we discard data taken during the daytime, data taken when the weather was particularly bad, or data from detectors that do not meet criteria summarized in Aiola et al. (2020) and detailed in Dünner et al. (2013).Once maps of the observations are made, further selection criteria are applied, as described in §4.1.

PLANET MAP-MAKING
3.1.Filtering And Noise Uranus maps were made using a special filter-and-bin map-maker built on the same framework used for building our normal maximum-likelihood maps (enki). 3The reason we do not use maximum-likelihood maps for planets is that these are not robust to "model errors", as described in e.g.Naess (2019).The difference between the true continuous sky and the pixelated model is interpreted as noise, leading to a bias which propagates along the scan direction for a noise correlation length.This is a significant effect for a bright, localized signal such as a planet.
It is possible to make maximum-likelihood models robust to these effects at significant computational and implementational cost (Piazzo 2017), but we here choose a simpler method that exploits the compactness of the planet signal relative to the atmospheric noise correlation length to separate the two.The idea is to split the data into two regions: a target region (including the planet) which we want to map, and a reference region (the rest of the data) which is used to estimate the behavior of the atmospheric noise. 4The atmosphere in the target region is then estimated by interpolating the atmospheric modes from the reference region, and it is then subtracted from the target region of the map. 5 To the extent that the signal in the target region is uncorrelated with any signal in the reference region, then this subtraction will not introduce bias on any scale in the target region.A point source like Uranus is very close to fulfilling these criteria, but in practice the outermost parts of the beam extend into the reference region, introducing a slight bias that we address in §3.2.
The target region must be large enough to capture all of the main beam, but small enough that the atmosphere's behavior in the target region can be adequately predicted using only data from the reference region.We choose a 12 radius area centered on the planet as the target region, and estimate the atmosphere's behavior here by exploiting its spatial correlation structure.Nearby detectors see similar parts of the atmosphere, so if at one moment in time some detectors are seeing the target region and some are seeing the reference region, the data from the latter are used to help determine the atmospheric modes in the former. 6This new interpolation technique significantly improves the noise modelling compared to the technique used in Louis et al. (2017), which used a simpler variant of filter-and-bin where the detector common mode was subtracted.
The resulting planet maps are cleaner and dominated mostly by white noise.In addition, we now inverse variance weight the detectors when making the planet maps, which was not done previously.This improved mapping has contributed to a more detailed understanding of the temperature-to-polarization leakage, which will be described in §5.1.

Bias From Map-Making
While it is approximately true that all the planet signal is contained inside the target region, a small part of it extends into the reference region due to the faint wings of the beams.This introduces a small but nonnegligible bias in the target region when the atmosphere is subtracted.
To study this bias, a set of simulated beam-convolved Uranus observations are processed through the mapmaker.The beam-convolved planet signal is simulated using a 2D beam model7 which uses the inverse Fourier transforms of two-dimensional Zernike polynomials as its basis set, similar to the 1D model of the beam described in §4.2.3.This simulated planet signal is then injected into several subsets of real time-domain data taken from CMB observations.The resulting data are then processed through the planet map-making pipeline, which includes estimating and subtracting the atmospheric noise.Both the CMB observations alone and the CMB observations with the planet signal injected are mapped separately, and then the former are subtracted from the latter.The result is a set of noise-free simulated maps of Uranus observations. 8 An example of a map of a simulated Uranus observation is shown in Figure 1.Radial profiles are constructed by taking the azimuthal average of the input beam model and the mean output of the simulated observations.The difference between these two radial profiles represents the bias due to the atmospheric subtraction.As seen in Figure 2, other than small residual fluctuations near the center of the beam, the bias can be well approximated by a constant offset. 9In these simulations the target region for atmosphere subtraction was chosen to have a 16 radius, but the qualitative conclusions we draw are applicable to the final planet maps used for the beam analysis which have a target region with a 12 radius. 10As shown in Figure 2, the effect of the atmosphere removal is almost exactly as expected: the value of the input beam at the target radius (16 in this case) is uniformly subtracted within that radius to produce the output beam.
So the planet map-making bias is well approximated by a constant offset in each map, determined by the value of the Uranus signal at the edge of the target radius used for noise estimation.This offset is estimated and subtracted prior to analysis, as described in §4.2.1.

Map Selection And Processing
A map of Uranus is made for each detector array, frequency, and observation. 11Only a fraction of these maps are used in the final beam analysis.We select maps made from observations done at night (23-11 UTC) taken with the secondary mirror in its final focus position for the season in question.We pick the maps where there are not too many zero-weight pixels (which are indicative of poor scan coverage), where the signal-to-noise is sufficiently high (> 10), and where there are not too many horizontal stripes (determined by visual inspection of the maps).This striping in the maps is caused by large-scale variations due to residual atmosphere.The rejection of maps with too much striping is new, and is done to simplify the analysis.It results in fewer maps being selected 8 In the end, it was discovered that including the atmospheric noise in the simulations does not affect the results, since the atmosphere removal is approximately linear.
9 These small residual features only appear in the T maps when solving for the maps in T , Q, U simultaneously.When we solve only for T and set Q = U = 0, the fluctuations disappear.This seems to be due to an issue with the conditioning of the covariance matrix in the simulations, and is not a physical effect.
10 During preliminary beam studies, maps of both Uranus and Saturn were made, with target regions both of radius 12 and 16 .Saturn was found to be too bright.In the Uranus maps with a target region of radius 16 , the signal-to-noise in the radial profiles past 12 was insufficient to be of use in the analysis.Hence, the Uranus maps with a target region of radius 12 were used, in order to make the best possible prediction of the atmosphere in the region of interest.
11 So we do not make any per-detector maps.For each detector array, frequency, and observation of Uranus, we combine the data from all the detectors that meet the quality criteria mentioned in §2 into one map.The difference between the "input" beam model (left) and the map of the "output" simulated observation (center).Note that in each case the color scale is a symmetrical log scale a in dB (with a linear threshold of −50 dB for the input and output maps and −45 dB for the difference map) and negative values are enclosed by parentheses.Here the target region where the atmosphere was estimated and removed during the map-making process was chosen to have a radius of 16 .The large-scale residuals seen are mostly constant within this target region of the map.The slight asymmetry in the beam (which may be safely ignored for DR4, as explained in §6.2) is due primarily to the positions of the telescope's optics tubes in the focal plane.
a A symmetrical log scale is logarithmic in both the positive and negative directions, with a small linear range around zero to avoid having values tend to infinity (https://matplotlib.org/stable/api/scale_api.html#matplotlib.scale.SymmetricalLogScale).
Fig. 2.-The average radial profile of the simulated Uranus observations for S16 PA2 at 150 GHz, comparing the azimuthal averages of the "input" 2D beam model and the mean map of the "output" simulated observations.The vertical dashed line indicates the 16 radius target region of the map within which the atmosphere was estimated and subtracted.(For the Uranus maps used to characterize the beam, a 12 radius target region is used, but the qualitative conclusions we draw here still apply.)Other than the small (roughly −40 to −30 dB) variations near the beam center, the difference is well approximated by a constant offset in the region from 3.5 to 10.0 where we fit for it.As described in §4.6, we adjust our beam covariance matrices to account for possible uncertainty due to variations in the range over which we fit the offsets, exploring the three independent ranges of 3.5 -5.0 , 5.0 -7.0 , and 7.0 -10.0 .
than in Louis et al. (2017).For each season, array, and frequency, the number of maps discarded due to striping ranges from 3 to 42 (on average 19), which represents between 4% and 36% of the maps (on average 16%).
The number of Uranus maps used for the beam analysis versus the total number of observations is shown in Figure 3 and Table 1.For example, in 2014 (S14) we made 129 observations of Uranus.In the case of PA2, we discarded 34 observations because at the time the telescope pointing was optimized for PA1 (so Uranus was not well Fig. 3.-Distribution of the total number of Uranus observations that were made (faded colors) and ultimately became part of the final beam analysis (solid colors) for all arrays combined, shown by observing seasons from 2013-16.The shaded regions between 0 UTC and 11 UTC, as well as 23 and 24 UTC, demarcate the ACT nighttime dataset (note that local time at the observing site fluctuates between UTC-3 and UTC-4).See Table 1 for a summary of the number of observations used vs total for each detector array and season.measured by PA2), 36 observations because the resulting maps had too much striping, and 8 because the signal-tonoise was too low.This leaves 51 Uranus observations to measure the beam.An example "good" map which was selected for the beam analysis and a "bad" map which was discarded due to too much striping are shown in Figure 4.In S13, a large fraction of the Uranus maps were thrown out because the observations were made in the early commissioning phase of the telescope, before it had achieved its final focus.
For the purposes of the beam measurements, we recenter each Uranus map by fitting a 2D Gaussian in the Fig. 4.-Two maps of individual Uranus observations for S15 PA2 at 150 GHz.(Top) A "good" map which was selected for the beam analysis.(Bottom) A "bad" map which was thrown out due to having too many horizontal stripes (caused by largescale variations due to residual atmosphere remaining after the atmosphere subtraction described in §3.1).In each case, the color scale is a symmetrical log scale in dB, with a linear threshold of −30 dB and negative values enclosed by parentheses.Any slight asymmetry visible here may be safely ignored for DR4, as explained in §6.2, and is due primarily to the positions of the telescope's optics tubes in the focal plane.
vicinity of the planet.The amplitudes from these fits are used to normalize the beam profiles to have peak values of unity.We also estimate the noise level in each map by computing the standard deviation of the pixel values outside the target region (the 12 radius area centered on the planet, described in §3.1).

Radial Profiles
The instantaneous ACT beams are slightly elliptical.However, when fitting the beam we work with azimuthally averaged ("symmetrized") radial beam profiles, treating the beams as if they were circular.In our case, the cross-linking of the scans (visible in Figure 2 of Choi et al. 2020) means that the telescope beams contribute to the maps with roughly equal weight at two different orientations that are approximately 90 degrees apart for a large part of the ACT sky coverage.The effective beams are thus circularized. 12his symmetrizing effect works well in temperature, but it does not help with the temperature-to-polarization leakage caused by beam asymmetry.However, this effect has been quantified and may be safely ignored for DR4, as explained in §6.2.
The ACT beams are small enough that we use a flat sky approximation for modeling them.We denote each beam map by B(θ, φ), where θ is the radial distance from the beam center and φ is the polar angle.We set B(0, φ) = 1.The symmetrized radial beam profile is then Each map, with a target region of radius 12 , is binned into a symmetrized radial profile with bins of width 11 out to a radius of 10 , for a total of 55 bins.The dominant source of noise in the Uranus maps comes from variations in the atmosphere which occur at relatively large angular scales, so there is significant correlation between the bins in each radial profile.
Before combining the radial profiles into one mean profile for each season, array, and frequency, each map's profile must be corrected for the offset due to the bias induced by the planet map-maker that was described in §3.2.This is done by fitting the model α/θ 3 + β (motivated in §4.2.4) plus a scattering term (described in §4.2.5) to each profile in the range 3.5 -10.0 and then subtracting the measured offset, β, from each.An example of such a fit is shown in Figure 5.
We then compute an average radial profile independently for each season, detector array, and frequency.When averaging radial profiles, each profile is weighted by 1/N where N is the white noise variance estimated from the corresponding map.The averaged radial profiles extend to roughly −40 dB of the peak, or 35 dBi (40 dBi) at 150 GHz (98 GHz).
An advantage of computing each radial profile then taking the average, compared to making one average Uranus map then computing its radial profile (as done in Hincks et al. (2010), for example), is that it more easily allows for the estimation of the full covariance matrix for the averaged radial profile.The off-diagonal elements of this matrix are important since there is significant noise correlation between the radial bins, which propagates into the low-uncertainty of the beam window function described in §4.3.

Radial Profile Covariance Matrix
As the averaged radial profile is computed for each season, array, and frequency, the covariance matrix between the amplitudes of the 55 radial bins is also computed in order to account for the covariant uncertainty on large angular scales.Given the small sample size, with only between 6 and 51 profile measurements used to estimate each matrix, a shrinking algorithm (Schäfer & Strimmer 2005) is used to down-weight the off-diagonal terms of the matrix.This method is described in detail in Appendix A. In our case this shrinking procedure is necessary as a standard estimate of the covariance matrix is often so noisy it becomes ill-conditioned and cannot be inverted.
The shrinkage technique works by combining an empirical estimate of the covariance matrix S (a highdimensional estimate of the underlying covariance with little or no bias) with a model T (a low-dimensional estimate which may be biased but has much smaller variance) to minimize the total mean squared error (sum of bias squared and variance) with respect to the true underlying covariance: (2) Here λ * is the parameter (often referred to as shrinkage intensity) that determines the contribution of each matrix.In our case, the covariance matrix S is an unbiased empirical estimate of the covariance, the sample covariance matrix, and the model matrix T is given by the diagonal of S: a common choice.We analytically calculate the optimal combination of the low and high dimensional estimates by determining λ * from S: where Var(S ij ) is an estimate of the variance of each covariance matrix element.The derivation of λ * can be found in Appendix A. For the analysis presented here, λ * ranges from 0.24 to 1.

Core Model
Since ACT's primary mirror has a diameter of D = 6 m, at a given wavelength λ the diffraction limited optical response of the telescope is restricted to spatial frequencies below max 2πD/λ.The Fourier transform of each beam is therefore compact on a disk.Thus, a natural choice of basis functions to model the beams in harmonic space is the Zernike polynomials, an orthonormal set of basis functions that is complete on the unit disk. 13 Expressed in polar coordinates ρ and φ in the aperture plane (where ρ is the radial distance 0 ≤ ρ ≤ 1, and φ is the azimuthal angle 0 ≤ φ ≤ 2π), the Zernike polynomials may be written as where m and n are integers such that n ≥ |m| and R m n (ρ) are a set of radial functions defined as (6) For an azimuthally symmetric beam, it is only necessary to consider the polynomials for which m = 0, which may be expressed as where P n (x) are Legendre polynomials (Born & Wolf (1980)).The Fourier transform of these polynomials is where J 2n+1 is a Bessel function of the first kind.Motivated by this, and following previous analyses (Hasselfield et al. 2013), we adopt the basis functions to fit the core of the beams in real space, where the parameter max varies the scale of the function's argument and n is a non-negative integer.

Asymptotic Behavior
On scales smaller than a few arcminutes, the noise is subdominant to the planetary signal, and so the beam profiles are well measured even by a single observation of Uranus.However, at larger angular scales, the noise 13 The Zernike polynomials are usually used to fit effects in the electric field, not in intensity, but since they are a basis set, we can use them to fit the intensity beams.So the asymptotic behavior of fn in Equation 9 is 1/θ 1.5 , not 1/θ 3 as we expect for the beams (as described in §4.2.4).
quickly becomes non-negligible.Thus, the asymptotic behavior of the beams cannot be separated from the background without making some assumptions about the shape of the beams far from the core.
The illumination of the ACT optics is controlled by a cryogenic Lyot stop.The beam's Fraunhofer diffraction pattern for monochromatic radiation is described by the squared modulus of the Fourier transform of this circular aperture (Born & Wolf 1980).The resulting intensity response with angle, or Airy pattern, is where k = 2π/λ is the wavenumber.14For ACT, k(D/2) = πD/λ is large, since the primary is several thousand wavelengths across.When the argument of the Bessel function is 1 and real, as is the case for θ > 2 , we can make the approximation (Abramowitz & Stegun 1972).In this regime, the Airy pattern asymptotically approaches such that the envelope of A(θ) falls as 1/ sin 3 θ.This implies that at angles larger than the "near sidelobes" (described in §5.2) and neglecting the effects of scattering, the beams should fall as 1/ sin 3 θ 1/θ 3 , since θ is small in the region considered for the beam models.
A different, more general, way to derive this asymptotic behavior of the beams comes from the geometric theory of diffraction from Keller (Keller 1956(Keller , 1959(Keller , 1962)).As detailed in Page et al. (2003), for an illuminated disk, the diffraction pattern from the rim, regardless of the interior, leads to a gain (response) of the shape for an electric field parallel (+) or perpendicular (−) to the edge, in the far field-limit.For unpolarized light, one can simply average over both polarizations.
The beam behavior derived here assumes perfectly smooth and uniform surfaces within the telescope (i.e.perfect coherence).Nonuniformity, including mirror surface roughness, leads to diffuse scattering (i.e.decoherence of the fields), and thus a reduction in the telescope gain (Ruze 1966).As described in the next section, an additional term is included in our beam model to account for this effect.

Scattering
The primary telescope surface can be characterized by an rms surface roughness, , and a correlation length, c.We measured these properties for ACT using photogrammetry,15 placing 284 targets on the corners of the primary's 71 panels where the panels join.The rms of the raw measurements was found to be approximately 30 µm.However, these measurements are not a representative sample of the rms of the surface, since the targets are at the most extreme points on the corners of the panels.By fitting a polynomial surface to the measurements and looking at the residuals over the entire surface, we find an rms of = 20 µm.Using the photogrammetry software, the uncertainty on the rms measurements is estimated to be 10 µm.This does not include any misalignments due to macroscopic deformations of the telescope.
From these measurements, we can also compute the correlation length, c.The correlation function is where r is the position on the surface, z is the meansubtracted departure from the best-fit designed surface (described in Section 3 of Fowler et al. 2007) is the separation, and the sum is over all measurement pairs with a separation d.We find that the shape of the correlation function follows the above form reasonably well.By fitting the photogrammetry measurements with Equation 14, for the ACT primary, we find c = 280 mm.This scattering leads to another angular scale (or shoulder) in the beam response, not governed by the diameter of the Lyot stop.
For a beam B(θ) normalized to unity at θ = 0, the corresponding solid angle Ω is The expression for the gain due to scattering off a rough surface, G(θ), given by Equation 8in Ruze (1966) 16 can be re-written in terms of the unit-normalized beam for an ideal reflector, B 0 (θ), its corresponding solid angle, Ω 0 , the total beam, B(θ), and its corresponding solid angle, Ω, by using the relation G(θ) = (4π/Ω)B(θ).We then obtain the equation with the S(θ) term given by where δ 2 = (4π /λ) 2 is the variance in the surface phase (and and c are the surface rms and correlation length, as described earlier).
The sum converges quickly, with four terms sufficient for our purposes.We have simulated the ACT beam by taking the Fourier transform of a perturbed surface with deformations of a Gaussian shape and verified that we indeed recover the Ruze equation as it is written above.17 The first term in Equation 16above is expected to decay as 1/θ 3 , but the same does not hold true for S(θ), so the inclusion of this term in our fitting procedure is important as it affects how the beam fits are extrapolated.This scattering term S(θ), which we refer to as the Ruze beam, contains roughly 1.5% of the solid angle of the main beams at 150 GHz.The addition of this term to the beam model is new to the DR4 analysis.

Beam Model Fitting
We separate the main beam for each season, array, and frequency into two domains: the core, where we fit the set of basis functions described in §4.2.3, and the wing, where the model is composed of a 1/θ 3 term in addition to the Ruze beam (where we use the measured values for and c described in the previous section).As explained in §4.2.1, a constant offset has already been subtracted from each measured beam profile before combining them into an average radial profile for each season, detector array, and frequency.The full model for each averaged beam profile can be written as The initial linear least squares fits are performed for the a n and α for a range of values of max , θ 1 , and n mode (allowing n mode to vary up to a maximum value of half the number of data points in the core).We first sample a range of values for the scaling parameter max and then allow this non-linear parameter to vary later on.Based on the optics of the system, at 150 GHz (98 GHz) we expect to have max 2πD/λ 19,000 (12,500).We sample max in integer steps from 1 to 30,000.
We fit the model out to 10 , a value which is chosen because this is the radius out to which the beam is measured to sufficient signal-to-noise.Despite the atmosphere subtraction in the target region that extends out to 12 , residual modes remain, so the signal-to-noise out to 12 is not sufficient.By only fitting out to 10 we ensure that the fit remains in the region where the measurement is dominated by the beam, not atmosphere.
The model is fit to each averaged radial profile, with its associated covariance matrix.We do not include continuity conditions at θ 1 , but any small discontinuities are well below ACT's resolution. 18n order to select the optimal number n mode for each beam profile, we must strike the right balance between minimizing the χ 2 and not adding too many parameters.To do this, we use the corrected Akaike information criterion (AICc), which estimates the relative quality of models based on both their goodness of fit and their simplicity.For each choice of max , θ 1 and n mode , the AICc is computed for the best-fitting model.
The uncorrected AIC is given by where k is the number of estimated parameters in the model and L is the maximum value of the likelihood function for the model (Akaike 1973(Akaike , 1974)).For small sample sizes, as is the case here, the AIC can exhibit a large bias.To account for this, we use the corrected criterion AICc, where and n is the sample size (Hurvich & Tsai 1989).
We select the values for max , θ 1 , and n mode corresponding to the lowest AICc. 19,20For the DR4 beams, the best-fit values for max range from 15,604 to 17,936 (from 10,960 to 11,285) at 150 GHz (98 GHz), the bestfit values for θ 1 range from 3.58 to 7.79 , and the best-fit values for n mode range from 9 to 13.We then use MCMC to sample the posterior distribution of the parameters max , the a n , and α, assuming uniform priors (Metropolis et al. 1953).This method produces an estimate of the parameter means and the full covariance between all the basis functions and wing parameters, including the non-linear scaling parameter max .
The term S(θ) in Equation 18depends on the beam's solid angle, Ω, but the calculation of Ω depends on the beam model.To get around this issue, the beam fit is performed iteratively, first using a rough estimate for Ω, and then re-computing Ω once a beam model is obtained, and then performing the fit again.In total, we fit the beam four times, re-estimating Ω each time.By the last iteration, the change in Ω becomes undetectable, and so we have converged to a final value for the beam solid angle.
We did test whether the data prefer an asymptotic wing fit term that differs from α/θ 3 , by fitting for α/θ b instead.While the best-fitting value for b often differs (by 10-20%) from three, it is not strongly constrained by the data.The AICc indicates that the addition of this new parameter is not warranted.
An example of the beam model fit to radial profile data is shown in Figure 6, for the S15 PA2 Uranus maps at 150 GHz, together with the residuals.The core funcbeam profiles shown in Figure 8.Any small discontinuity that may have been present in the initial fit disappears due to the resolution with which we do the transform.tions, α/θ 3 term, and the Ruze beam are indicated in the lower panel, with the core functions used out to a radius of θ 1 = 5.59 , then transitioning to the α/θ 3 plus Ruze beam model at larger radii, and using this model to extrapolate past 10 .We measure the beam profiles down to −40 dB from the peak, which would leave a few percent of the beam's solid angle unaccounted for if we did not use our fit to extrapolate past 10 .
Figure 7 shows how the symmetrized model compares to the average Uranus map data for this case.The full set of radial profile fits for the S15 data is shown in Figure 8 (with some small additional corrections, as described in §4.4).We also obtained model fits for the S13, S14 and S16 data that make up DR4.The solid angles, gains, and FWHMs for all seasons are reported in Table 2 (again, with the small corrections from §4.4).For reference, 1 nsr 0.0118 arcmin 2 .Even though we make a high signal-to-noise measurement of each beam, the uncertainty on the solid angle is limited to ∼ 2% because of the uncertainty in extrapolation, which depends on the model.The fractional uncer- tainties in the solid angles at 98 GHz may be larger than those at 150 GHz in part because at 98 GHz the beam is broader, so a larger fraction of the beam is affected by the uncertainty in our extrapolation.When using the beam solid angle to calculate the flux density of a point source, the uncertainty on the measurement may be reduced by first applying a matched filter to both the beam and the source in order to remove the scales associated with the extrapolation uncertainty.
While much of this fitting method follows the approach used in Hasselfield et al. (2013) and Louis et al. (2017), notable improvements are the addition of the scattering beam term in the model, the use of MCMC sampling that includes estimating the non-linear parameter max , the exploration of different radii (θ 1 ) out to which the core model is fit, and the use of the AICc to choose the final best-fit model.

Beam Window Functions
In spherical harmonic space, the beam information is encoded in the harmonic transform b and the window function w = b 2 , which describes the instrument's response to different multipoles, .This window function is an essential component of the DR4 power spectrum analysis in Choi et al. (2020).
The harmonic transform b is the Legendre transform, or more accurately the Legendre polynomial transform, of the beam radial profile: For small beams such as that of ACT, this is effectively a Fourier transform.The derivation of the Legendre transform and details about how the transform is computed are presented in Appendix B.
We use b instead of B to indicate the division by Ω, which normalizes b to unity at = 0 (since P 0 = 1).B = Ωb has units of sr, whereas b is dimensionless.We extrapolate the model beyond the fit radius of 10 when computing the transform.This is necessary to capture the low-part of the window function, and to account for the part of the beam solid angle that is beyond the range we fit.
A subset of the beam transforms from DR4 is shown in Figure 9.A similar figure in Aiola et al. (2020) shows window functions, b 2 , which are used to correct the power The difference between the measured beam (left) and the beam model (center).Note that the color scale for the residuals is different from the other two maps.Both PA2 and PA3 are shown here, since PA2 serves as an example of a less elliptical beam, and PA3 serves as an example of a more elliptical one.The corresponding maps for S15 PA1 at 150 GHz resemble the ones shown for PA2 here, and the maps for S15 PA3 at 98 GHz are similar to the ones shown here for PA3 at 150 GHz, but broader.For each individual detector array, the residuals are fairly constant from one season to the next.As shown in Figure 6 for PA2, the azimuthal average of the residuals is consistent with zero, which is why the fit is successful, despite the residuals visible in the maps.These residuals are expected to have a quadrupole-like shape, since we know our beams are slightly elliptical, and the quadrupole is the dominant asymmetric azimuthal mode for an elliptical beam.
spectra.For a given array and season, if the beam transforms for PA1 and PA2 are, respectively, b PA1 and b PA2 , then for the auto-power spectrum of the PA1 or PA2 maps the window functions are (b PA1 ) 2 and (b PA2 ) 2 , and for the cross-spectrum of the PA1 and PA2 maps the window function is b PA1 b PA2 .

Additional Corrections
The resulting beam models and covariance matrices are an accurate description of the binned radial beam profiles, but they must be corrected for some systematic effects.Following the same approach as in Hasselfield et al. (2013), corrections are applied to account for the pixelization of the planet maps, the binning of the maps into radial annuli, Uranus' angular diameter, and the planet's effective frequency.The effect of each of these corrections is shown in Figure 10 for S15 PA2 at 150 GHz, as a typical example.
To correct for the pixelization of the planet maps, we divide the beam transform by the azimuthal average of the pixel window function, sinc(pk x /2π) sinc(pk y /2π), where {k x , k y } are the spatial frequencies and p is the pixel size in radians.This is a 0.1% effect for < 10, 000.Then, we estimate and correct for the transfer function induced due to binning the planet maps into radial annuli.This is done by simulating planet maps, using the best-fitting beam profile model as the input, and then estimating their radial profiles following the same procedure as with the data.Comparing the input model with the output radial profile gives an estimate for the transfer function resulting from the radial binning.This is a 1% effect for < 10, 000.We also correct for Uranus' angular diameter, since Uranus is large enough that it cannot be treated as a point source given the precision to which we measure the beam.For each season, array, and frequency, we assume Uranus is a disk with radius equal to the weighted mean of the radii for all the Uranus observations contributing to the beam measurement.We then deconvolve this shape in harmonic space using the function 2J 1 ( r)/ r, where r is the radius of Uranus' disk in radians.The factor of 2 normalizes the function to unity as → 0. This is a 0.1% effect for < 10, 000.
The beam at this point describes the telescope's response to a point-like source with an approximately Rayleigh-Jeans (RJ) spectrum.Near our frequencies, the temperature spectrum of Uranus goes roughly as ν −0.25 (Planck Collaboration Int.LII et al. 2017).It is sufficiently close to the RJ limit (ν 0 ) for our purposes.We apply a simple, first-order correction to obtain the relevant beam for the CMB blackbody spectrum using the band effective frequencies from Thornton et al. (2016). 21For this radiation with a band effective frequency ν CMB , the beam is taken to be B ( ) = B( ν RJ /ν CMB ), where ν RJ is the effective frequency for radiation with a Rayleigh-Jeans spectrum.This is a 1.5% effect for < 10, 000.
The resulting beam models are referred to as the "instantaneous" beams, and it is in fact these corrected beams that are shown in Figures 8 and 9.

Jitter Beams
In practice, the effective beam for a given sky region, season, detector array, and frequency is broader than the instantaneous beam.This is due to combining observations taken on multiple different nights throughout each season, so the resulting effective beam for each map is not as sharp as the beam inferred from planet observations which have been carefully recentered before co-adding.Broadening can be caused by variations in the pointing and global alignment, as well as possible small changes in the beam over the course of a season.
As described in Aiola et al. (2020), ACT's blind pointing accuracy for DR4 is comparable to the average beam FWHM.If left uncorrected, this would significantly broaden the effective beams.Instead, as was done in Louis et al. (2017), we correct the pointing by comparing the observed positions of bright point sources to their known catalog positions.This is done for each 10-minute section of the time-ordered data from each detector array, as described in Aiola et al. (2020).Instead of performing the fit in the time domain as in Louis et al. (2017), due to the larger data volume the fit is now performed in mapspace.The resulting fit is obtained more quickly, but is slightly less accurate, leaving a larger residual variation in the beam due to pointing uncertainty.
This residual variation is captured by the effective beam, B eff , which is parametrized in terms of the instantaneous beam, B , and a correction term that is Gaussian in , as If we were to interpret the Gaussian correction term as arising purely from residual pointing errors, then V would be the residual pointing variance in square radians, which is why it is also referred to as "pointing jitter."This variance, which in practice also includes the additional errors due to alignment and seasonal changes in the beam, is expected to be different for each sky region used in the power spectrum analysis.These different regions (Deep1, Deep5, Deep6, Deep56, Deep8, BOSS, and AdvACT) are shown in Figure 2 of Choi et al. (2020).
To estimate V , in each region we create a catalog of the brightest point sources that are found in the maps using a matched filter, that have a signal-to-noise of at least 10 in each map, and that are matched to catalogs of known sources.For each of these sources, we then obtain an estimate of the variance V and its associated  uncertainty for each season, array, and frequency.To do this, we crop a section (10 × 10 ) of each map around the source and remove large-scale variations by high-pass filtering the data (by multiplying the data by 1−G, where G is a Gaussian filter with a FWHM of 8 ).We then compute the effective beam using Equation 22, multiply this model by the local pixel-window of the cropped map, transform it into map space, high-pass filter it, and then compute the χ 2 using an estimate of the local white noise amplitude.An example of such an individual source fit is shown in Figure 11.
This procedure was validated by verifying that when simulated point sources convolved with Equation 22(with a known value for V ) are injected into real CMB data and run through the pipeline, the input value for V is recovered.
For convenience, we exclude sources identified as galaxies, pairs of galaxies, and planetary nebulae in the known source catalogs from this fitting procedure, as these types of sources sometimes appear to ACT to be extended.This leaves approximately 20 sources (for the smaller, deep regions) to 520 (for the largest region) that are included in the fits.They are primarily quasars and radio sources.
Once all N sources have been fit using this method, we use the resulting set of pointing jitter estimates d and their uncertainties σ to obtain an estimate for the mean and intrinsic scatter, V and σ V , of the effective pointing jitter in each map.The likelihood for this is written as where is a normalization factor.We estimate the posterior distribution for V and σ V using MCMC with the Metropolis-Hastings algorithm (Hastings 1970), assuming a uniform prior on V .An example of these posterior distributions for one of the maps is shown in Figure 12.
This likelihood, which was not used for previous ACT analyses, better takes into account the intrinsic scatter in the residual pointing variance, resulting in an improved estimate.
After computing the pointing variance independently in each map, we use Equation 22to obtain an effective beam, and associated covariance matrix, for each season, region, detector array, and frequency.An example of the effect of the pointing jitter correction on the beam transform is shown in Figure 10.The jitter values and uncertainties are given in Table 3 and the resulting solid angles and uncertainties for the effective beams are given in Table 4.
The jitter values for PA3 at 98 GHz are significantly greater than those at 150 GHz.This empirical finding serves as a useful reminder that this "jitter" encompasses not only pointing errors, but also any changes in the beam throughout each season.These changes may be greater at 98 GHz since that is where the beam is broader.
Some of the best-fit jitter values for the Deep8 region are negative, but considering their uncertainties, they are consistent with positive values.We allow for negative jitter values in the fits in order to account for all possible effects on the beams, not just those due to pointing.While we have not investigated whether these negative values are related to a data quality issue, it should be noted that the data from Deep8 were not used in the cosmological analysis due to the poor cross-linking in the region.
In Appendix C we provide the effective beam solid an- gles for beams at different frequencies.We applied the first-order spectral correction from §4.4 with effective central frequencies corresponding to synchrotron emission, dust, and the thermal Sunyaev-Zel'dovich effect.These beam solid angles differ from the CMB estimates by 0.3 to 5.3%.

Beam Transform Covariance Matrix
There have been several modifications since Hasselfield et al. (2013) in how we estimate uncertainties in the beam model.Previously, the non-linear scaling parameter max was varied in the fits, but the uncertainty in this parameter was not propagated to the covariance matrices.We now estimate the covariance of each beam transform, b , by computing the Legendre transform of the sampled posterior distribution for the parameters describing the radial profile, as in Equation 21, and then applying the small-order corrections from §4.4.We then estimate the covariance matrix from this suite of -space beam transforms and the added uncertainty associated with the jitter correction.These matrices are large, since the beam transforms are computed at each integer from 0 to 30,000, so we do not store them in their entirety.Rather, we decompose each matrix into independent modes (via singular value decomposition) and store the largest 10 modes.We find that this is sufficient to capture the majority of the covariance (we do not discard any modes with a singular value larger than 10 −3 of the maximum value).
We include additional modes to account for possible uncertainty due to variations in the surface rms (which enters into the calculation of the Ruze beam, as described in §4.2.5, and was fixed to our estimate of 20 µm),22 as well as the range over which the model for each beam offset is fit (which initially was 3.5 -10.0 ).For the surface rms, the values explored are 20 µm and 30 µm.For the region over which the offsets are fit, the three independent ranges explored are 3.5 -5.0 , 5.0 -7.0 , and 7.0 -10.0 .We store the top 3 modes associated with these model variations.The final uncertainty for each beam is thus composed of a total of 13 modes.These beam uncertainties are later added to the data covariance matrix as part of the power spectrum analysis pipeline.
An example of the effect of these additional modes on the beam transform uncertainties is shown in Figure 13.The significant increase in the uncertainty at low multipoles is mainly due to the inclusion of the different fit ranges for the offsets.Previously, as described in Hasselfield et al. ( 2013), the uncertainties were simply doubled from their formal values to account for potential systematic variations due to different fitting ranges.Even though this earlier method lead to a smaller estimate of the beam uncertainties at low multipoles, this did not have a significant effect on our results for DR3, since the uncertainty on the beams is subdominant in the power spectrum analysis, and a low-cutoff of 500 (350) was applied to the T T (T E and EE) data.
Another difference compared to the DR3 analysis in Louis et al. (2017) is the treatment of calibration uncertainty.As described in Choi et al. (2020), for each season, region, array, and frequency the angular power spectra from ACT are calibrated to the Planck tempersured parameter for the Ruze beam, the correlation length c, since it is well constrained by our measurements.ature maps in the range 600 < < 1800.The Louis et al. (2017) analysis factored out the beam amplitude and uncertainty at an "effective" calibration scale (chosen in that case to be = 1400).23However, we have changed this procedure for DR4 to reflect the fact that we are not calibrating the data at a single , but over a range of values.We now simply normalize the beam to unity at θ = 0 and treat the calibration, and its uncertainty, separately in the ACT analysis.This is why the fractional beam uncertainties in Figure 9 are no longer smallest at = 1400, differing from those reported for ACT DR3.

Main Beam
Measurements of Uranus in polarization at visible and near-infrared wavelengths have shown that its diskintegrated polarization is less than 0.05% (Schmid et al. 2006).While we are not aware of any similarly precise data at millimeter wavelengths, it is expected that the relevant scattering effects in the planetary atmosphere would be much weaker, resulting in net polarization levels considerably lower than the measurement cited above.
Since we do not expect Uranus to be significantly polarized in the bands observed by ACT, we interpret any polarization response measured from Uranus as being due to temperature-to-polarization (T-to-P) leakage.Although this leakage is relatively small in magnitude, the ACT DR4 data are now sensitive enough that we must account for it in our analysis.To do this, we use observations of Uranus to build an -space T-to-P leakage function for each season, array, and frequency.This measurement and correction for the leakage in the main part of our beams is new to DR4 (see also Aiola et al. 2020;Choi et al. 2020).
We begin by making maps of the Q and U Stokes parameters for the same set of Uranus observations chosen in §4.1 to fit the main beam.We then convert each set of {Q, U } polarization maps to the radial Stokes parameters {Q r , U r } (following the third definition of cross polarization in Ludwig 1973) using the flat-sky approximation: where θ ≡ (θ, φ θ ) are standard polar coordinates with the beam centroid as their origin and φ θ increases clockwise from the positive y-axis (assuming that one uses the convention in which +x points to the right and +y points upward).Here Q and U follow the COSMO convention (Górski et al. 2005a), whereas for the ACT maps released as part of DR4, the polarization components are defined by the IAU convention (Hamaker & Bregman 1996).The polarization convention initially used for {Q, U } does not matter once we have transformed to {Q r , U r }.While it is not obvious how the polarized E or B beams ought to behave at larger radii, we do have some intuition about Q and U .For an unpolarized point source such as Uranus, any polarized signal will be the result of beam differences between the two axes of a polarimeter.Although the beams may be slightly different, due to e.g., differential ellipticity, we still expect each of them to decay radially as 1/θ 3 as they are part of the same diffraction-limited optical system.The Q and U maps are essentially just radially independent linear combinations of the various detector axes in an array, appropriately weighted by the map-maker, so the same asymptotic behavior should apply.This asymptotic behavior would then hold true for the Q r and U r maps as well.This was confirmed with simulations where the beams for individual polarimeter axes were constructed using an "elliptified" version of the azimuthally averaged intensity beam for a given season and array.We thus use the same basis functions in the core and the α/θ 3 term in the wing to fit the beam in polarization as we did for the main beam in §4.2.6.
Since we are ultimately interested in how leakage manifests itself in the angular power spectra, we need to translate any polarized beam models of the azimuthally averaged radial profiles of Q r and U r , Qr and Ũr , to an -space representation of E and B. Conveniently, Qr and Ũr have a direct correspondence to the azimuthally averaged -space E and B transforms, Ẽ( ) and B( ).As shown in Appendix D, there exists a simple relation between these components, which in the flat-sky approximation takes the form of a second-order Hankel transform: From the set of Uranus maps in the {Q r , U r } basis, examples of which are shown in Figure 14, we thus construct average radial profiles for each season, array, and frequency, and we fit them in a similar way as the main beam is fit in §4.While for the usual temperature beam fitting pipeline we fit an offset to each individual Uranus profile before taking an average, we do not do this in polarization since we expect the mapping transfer function in that case to be zero. 24The only difference in the fitting procedure here is that we do not include a scattering term in the polarized beam model, since it is not expected to matter to first order, and it is unclear how it would vary for the different detector polarizations. 25xamples of the radial profile model fits in polarization are shown in Figure 15.
As can be seen in Figures 14 and 15, while our model is a good fit to the Q r and U r radial profiles, there are significant, quadrupole-like residuals in the maps when our model is subtracted.It turns out that the leakage beams are far less azimuthally symmetric than the temperature beams (as can be seen by comparing Figures 7  and 14).The treatment of this leakage will be revisited for upcoming ACT beam analyses and may be improved upon, for example, by fitting a 2D model to the polarized beam profiles in order to properly capture the asymmetry. 26In the meantime, the fitting done here is sufficient.The level of residuals seen in Figure 14 has an insignificant effect on the results of the power spectrum analysis for DR4 (see §6.2).The features visible in the residual maps can be expected to occur due to physical effects.For example, we simulated the T-to-P leakage for a point source due to both differential beam ellipticity between the two axes of a polarimeter and a polarization angle off-set.The resulting maps of the simulated leakage beam in Q r and U r had strong quadrupole-like features.As shown in Figure 15, the azimuthal averages of the residuals in Q r and U r are close to zero, which is why the fits works well, despite the appearance of the residuals in the maps.
The Q r and U r radial profiles fits and their uncertainties are transformed to -space using a similar approach as for the temperature data, except using Equation 25 instead of the usual Legendre transform.The transforms for our example case are shown in Figure 16.These transforms do not include corrections for small systematic effects, as was done for the main beam fitting pipeline in §4.4.This is because we are ultimately interested in the ratios of transfer functions for the leakage beams, in which these -space corrections cancel out.To determine the polarization leakage beams, an example of which is shown in Figure 17, we divide the E and B transforms by their corresponding (uncorrected) T transform.For PA1 and PA2 the leakage values are within 1.5% (comparable to the leakage values for SPT-3G, Dutcher et al. 2021), whereas for PA3, the leakage increases around 4000-6000, and reaches over 6% (4%) at ∼10,000 at 150 GHz (98 GHz).We attribute the higher leakage for PA3 to an imperfectly optimized horn design in this first generation of multichroic polarimeters.
The top 10 modes from the leakage beam covariance matrices are stored.In addition, similar to the main beam analysis, two modes are added to account for variations in the model for the temperature beam. 27The final leakage beam uncertainties are thus comprised of 12 modes.
The leakage beams were used in the DR4 power spectrum likelihood, as explained in §5.3.

Polarized Sidelobes
As described in Louis et al. (2017) and Aiola et al. (2020), we detect polarized sidelobes of the main ACT beams.Although weak in amplitude, these sidelobes cause noticeable T-to-P leakage.The sidelobes for PA1 and PA2 were shown in Louis et al. (2017).While PA3 was not used in that analysis, it was mentioned at the time that polarized sidelobes were not detected for PA3.However, using additional observations of Saturn to conduct a more thorough analysis, we have detected sidelobes in PA3 at 150 GHz, with an amplitude roughly 10% that of the sidelobes in PA1 and PA2.The general features of the sidelobes are common to all detector arrays; they consist of a group of compact lobes, each resembling a slightly elongated image of the main beam, with approximate four-fold symmetry, strongly polarized perpendicular to the radius from the beam center (which corresponds to −Q r ).This results in most of the leakage due to the sidelobes being from temperature into E-mode polarization rather than into B-mode polarization.These sidelobes are stable in time.We also observe that sidelobes from Saturn are only seen if Saturn lies in a focal plane's field of view.As explained in Section 3.8 of Aiola et al. (2020), this is consistent with the sidelobes being due to an optical effect inside the receiver.
The sidelobes for PA1 and PA2 are at a distance of roughly 15 from the beam centroid, with an additional set visible at 30 for PA1.The sidelobes for PA3 are located approximately 30 to 40 from the beam centroid.This larger angular separation means the strongest Tto-P leakage occurs at 300 for PA3 compared to 500 for PA1 and PA2.Also, since the sidelobes only map to the sky when the main beam is also in the field of view, fewer detectors are affected by each sidelobe for PA3.We do not detect any sidelobes in PA3 at 98 GHz, which implies they must either be of significantly lower amplitude than those at 150 GHz, or in a different position.For PA1, PA2, and PA3, the amount of solid angle contained in these sidelobes is roughly 1.2%, 3.4%, and 0.15%, respectively, of that in the main beam.
Preliminary studies suggest that these sidelobes are due to diffraction caused by the arrays of metal elements that make up ACT's optical filters.We note that similar filters (Ade et al. 2006) are also used for the South Pole Telescope (SPT), POLARBEAR, and BICEP/Keck (Padin et al. 2008;Arnold et al. 2010;Keating et al. 2003) at varying locations in the optical path. 28f the sidelobes were indeed due to diffraction from the filters, we would expect them to occur at a larger radius for lower frequencies.Based on our calculations, at 98 GHz, the sidelobes would appear starting at a distance of approximately 47 from the beam centroid.Since this is roughly the size of our field of view, the sidelobes would then map to the sky for only a small fraction of the detectors, the ones at the edges of the focal plane.This is consistent with the lack of detectable sidelobes at 98 GHz.We do not apply any sidelobe corrections to the PA3 data at 98 GHz.
As mentioned in Aiola et al. (2020), to study the sidelobes we use observations of Saturn.While Saturn's brightness appears to induce a non-linear response near peak amplitude, it is useful for studying the relatively weak sidelobes.There are no issues due to non-linearity when observing the sidelobes.This is confirmed by the fact that the amplitude of the sidelobes seen by a detector is consistent, whether they are seen before or after the detector sees Saturn's peak.
We apply the same treatment to the sidelobes for all detectors at 150 GHz.In short, we model the sidelobes as a sum of polarized, spatially shifted copies of the main beam and fit the amplitudes of the T , Q, and U components of each beam instance using maps of Saturn.We then use this model to deproject the sidelobes from the time-ordered data prior to map-making.The idea is to subtract the total flux in the sidelobes, even if their shapes are not exactly zeroed out in the maps.It can be shown that this removes the low-T-to-P leakage.Maps of the sidelobes are shown in Figure 18, along with the T-to-P leakage functions.
For PA1 and PA2, we re-use the sidelobe models from DR3, constructed using observations of Saturn from 2014.For PA3, which was not part of DR3, we use Saturn observations from 2015 to characterize the sidelobes The difference between the measured beam (left) and the beam model (center).Features in the residuals could arise due to physical effects, such as differential beam ellipticity between the two axes of a polarimeter or a polarization angle offset (Hu et al. 2003).The level of residuals seen in the last column has an insignificant effect on the results of the power spectrum analysis.As shown in Figure 15, the azimuthal averages of the residuals are close to zero, which is why the fits are successful.
for DR4 in the same way as was done for DR3.
Looking more closely at how these sidelobe models are constructed, we begin with a series of observations of Saturn to which we apply the same data selection criteria as we did to Uranus, as described in §2.We then map the chosen observations with the moby229 map-maker (the same map-maker that was used for the beam analysis in Louis et al. 2017) and coadd the maps together to produce one map of Saturn for each detector array (PA1, PA2, PA3).The maps are coadded with weights based on an estimate of their white noise level determined outside the planet region.
We then visually identify regions in the maps containing compact polarized sidelobes.As can be seen in Figure 18, the pattern of the sidelobes in each map has approximate four-fold symmetry, four groups of sidelobes appearing at roughly equal distances from the beam centroid.For each sidelobe we choose how many copies of the main beam should be used to model it.Stronger, elongated sidelobes are modelled by two copies of the main beam, to capture the elongation; weaker sidelobes are modelled with a single copy of the beam.Then for each of the four groups of sidelobes we fit the position and amplitude of the copies of the main beam.In theory this fit could be done with any of the signals, but currently we perform this fit to a map of P = Q 2 + U 2 because that is a bright, clear signal.With the model positions fixed, we then fit the amplitudes for each signal (T , Q r , U r ) independently.The result is a base model for the sidelobes, but this model does not yet contain per-detector detail.
We next account for the fact that the sidelobes do not always appear in all the detectors.This is relevant to the extent that the per-detector weights in the planet maps are different from the per-detector weights in the survey maps (which are used for the CMB analysis, for example).This could possibly be a significant (up to tens of percent) effect, since the moby2 map-maker used to make the planet maps doesn't weight the detectors by their noise, whereas the enki map-maker used for the survey maps does.The sidelobes occur primarily in the detectors at the periphery of the array, which is also where the noise tends to be higher, so this is a correlated effect.
To deproject the sidelobes from the per-detector timeordered data, we need to estimate whether or not the sidelobes are seen by individual detectors, based on their position on the focal plane.The parameter we want to estimate is the radius of the "aperture", or circle, around a detector such that if Saturn falls within the circle, the beam sidelobes are visible to that detector.To estimate this radius, we first divide the detectors into subsets.For PA3, the subsets are the three hexagonal wafers of detectors, described in Thornton et al. (2016).For each detector subset, we re-make maps of Saturn, and measure the total sidelobe flux in the resulting coadded map.We then compare these measurements to a model for the sidelobe flux as a function of radius, which scales with the fraction of detectors in each subset that would see each sidelobe.
As for the main beam analysis in §4.4,small corrections are made to the model to account for systematic effects.
The sidelobe removal described here is not perfect, so there is still residual T-to-P leakage in the ACT data due to the sidelobes.The residual T E and T B leakages, shown in Figure 18, are estimated by making new maps of Saturn after our sidelobe removal.Since the sidelobes for PA3 are already significantly weaker than PA1 and PA2, we do not compute residuals for PA3.We add the residuals for PA1 and PA2 to the main beam leakage for use in the power spectrum analysis, as described in the next section.
The residual sidelobe leakage in Figure 18 can be compared to the leakage in the main beam shown in Figure 17.The effect of each of these components on the final spectra (if there is no leakage correction as in §5.3) is shown in Figure 19.At low , the residual sidelobe leakage dominates, whereas the main beam leakage grows larger by ∼ 2000, and dominates at high .The main beam leakage and the residual leakage from the sidelobes are included in the power spectrum likelihood (see Section 12 of Choi et al. 2020) by making use of a leakage-corrected model for the T E and EE theory spectra, computed each time the likelihood is estimated.The corrected model spectra, T i E j and E i E j , are related to the input theory spectra, T i T j , T i E j and E i E j , via

Leakage Correction
Here the i and j subscripts denote different spectra and the γ factors encode the -dependent leakage (the sum of the main beam leakage from §5.1 and the residual sidelobe leakage from §5.2) and are shown in Figure 20.At 150 (98) GHz, the amplitude of γ is never greater than 0.0035 (0.038).Around ∼ 1000-3000, this leakage correction is roughly a few-percent (∼ 1-4%) effect for T E and a sub-percent (∼ 0-0.4%) effect for EE.
As explained in Choi et al. (2020), in the power spectrum analysis we first compute a power spectrum for each season, sky region, detector array, and frequency.We then coadd over seasons and arrays to obtain one spectrum per region and frequency.Finally, the regions are  -Measured temperature-to-polarization leakage in the main beams for S15 for all arrays, showing the temperature leaking into E-mode and B-mode polarization.Note the different y-axis ranges.For PA1 and PA2 the leakage values are within 1.5%, whereas for PA3, the leakage increases around 4000-6000, and reaches over 6% (4%) at ∼10,000 at 150 GHz (98 GHz).We attribute the higher leakage for PA3 to an imperfectly optimized horn design in this first generation of multichroic polarimeters.The uncertainties (lower panels) include the adjustments for model variability.
divided into two groups based on the detection thresholds for point sources: deep (Deep1, Deep5, Deep6, and Deep8) and wide (BOSS and AdvACT).The spectra for the regions in these two groups are coadded, resulting in a single deep and a single wide power spectrum at each frequency.
To obtain the γ factors, we coadd the leakage beams for individual seasons and detector arrays using the same weights used to coadd the spectra, giving an effective leakage beam for both the deep and wide coadded spectra at each frequency.The uncertainties in the γ factors are incorporated in the data covariance matrix by using the main leakage beam uncertainties and the sidelobe residuals to compute another covariance matrix (similar to the main temperature beam covariance matrix) which is then added (in quadrature) to the data covariance matrix.
As described in Aiola et al. (2020); Choi et al. (2020), including this leakage correction in the likelihood reduces the residuals compared to the best-fitting ΛCDM model but does not have a significant effect on the inferred cosmology.Choi et al. (2020) (in Section 12.3) also carried out a test by fitting for two scaling factors (one at 150 GHz and one at 98 GHz), that multiply the nominal values of the γ factors.The data support the baseline model, where the scaling factors are unity.
In the DR3 analysis in Louis et al. (2017) the residual leakage from the sidelobes was similarly treated as a systematic uncertainty in the cosmological power spectrum analysis, but the correction for the main beam polarized leakage is new to DR4. 30 6. DISCUSSION 6.1.Beam Products As part of DR4, several ACT data products were made publicly available on the NASA Legacy Archive for Mi-30 As described in Choi et al. (2020), for DR4 the ACT team adopted a blinding strategy in an attempt to prevent confirmation bias on the cosmological parameters.The analysis pre-unblinding used maps that deprojected the polarized sidelobes, but did not include the additional polarized leakage correction from the main beam or the sidelobe residuals.The post-unblinding analysis revealed some features in the T E residuals to ΛCDM that led us to a more thorough search for possible sources of systematic uncertainty in the T E spectrum.As a result, it was at this point that we added the correction in the likelihood for the T-to-P leakage from both the polarized sidelobe residuals, and the main beam.a Saturn's brightness appears to compress the detector gain by a few percent.The effect can be ignored here because it is a small correction on the sidelobe leakage, which itself is a small effect.
crowave Background Data Analysis 31 (LAMBDA) and at the National Energy Research Scientific Computing Center (NERSC).Details of all these data products are given in Mallaby-Kay et al. (2021).
The beams described in this paper and used for the analysis in Choi et al. (2020) and Aiola et al. (2020) are included in the "ancillary products" section of the data release.For each season/region/array/frequency combination, there are multiple beam files.Both the real-space radial beam profiles and their harmonic-space transforms are provided, both for the instantaneous beams and for the beams with the jitter corrections included, and 31 lambda.gsfc.nasa.gov/product/act/actpol_prod_table.cfm  in each case are available both with and without the Rayleigh-Jeans-to-CMB spectral correction.
Note that the instantaneous beams are not regiondependent, so only the jitter-corrected beams contain region flags in their filenames.The instantaneous beams are suitable for time-domain analysis, but they should not be used when analyzing maps.The maps released as part of DR4 have not been corrected for the instrument beam.When working with the maps, one should use the jitter-corrected beams where the pointing variance and its uncertainty have been accounted for.The beam harmonic-space transforms can be used to correct for the beam effects in harmonic space (Bond & Efstathiou 1987).
The T E and T B harmonic-space leakage beams for the main beam and the sidelobes and their residuals are also included in the ancillary products for DR4.
Since DR4 includes the data from DR3 as a subset, it is possible to compare some of the DR4 beams with the beams for the same season, sky region, detector array, and frequency released as part of DR3.As shown in Appendix E, despite the differences in the analyses, the beams from DR3 and DR4 are consistent.

Beam Asymmetry
The asymmetry of the ACT main beams is relatively well described as elliptical, with aspect ratios varying from 1% to 20% across the different arrays, as shown in Table 1 of Choi et al. (2020).Here we give a brief summary of an investigation into this azimuthal asymmetry and the associated spurious signal.
Beam asymmetry is, to leading order, responsible for two effects.The asymmetric convolution creates anisotropy in the sky maps that distorts the shape of point sources and introduces statistical anisotropy in the inferred CMB.The magnitude of the effect is reduced by increased cross-linking: as the telescope observes a position on the sky using approximately orthogonal scan directions, the spurious signal from one scan roughly cancels with that of the other scan.The cross-linking in the ACT maps is sufficient to reduce the contamination to the power spectrum from this effect to an insignificant amount.
The second effect of beam asymmetry is to introduce T-to-P leakage.The beam asymmetry introduces a dependency in the time-ordered data on the position angle of the instrument: observations with different scan directions yield systematically different data.For asymmetric beams with a quadrupole shape, which is the dominant asymmetric azimuthal mode of our approximately elliptical beams, the dependence on the position angle is interpreted by the map-maker as a linearly polarized sky component (Hu et al. 2003).In contrast to the first effect, the T-to-P leakage is not averaged down by cross-linking.Instead, the leakage is reduced by the instrument's orthogonally polarized co-pointing detectors: the leakage picked up by one detector approximately cancels with the leakage picked up by its partner.32As mentioned in Choi et al. (2020) the residual leakage causes an additive bias to the T E and T B power spectra that is roughly constant with multipole with an amplitude that is less than 0.2σ away from zero.No attempt to remove this leakage has been made.

Conclusion
In this paper, we have presented the analysis of the ACT beams for DR4, which includes data from 2013-16.Improvements to the beam pipeline include: better atmosphere subtraction for the Uranus maps, a new scattering term in the model that is fitted to the main beams, a better estimate of the uncertainty in these fits, and residual T-to-P leakage terms that are included in the ACT power spectrum likelihood.Considerable effort was spent developing a realistic model of the beams (including optical effects) and the mapping process, in order to study all elements of the analysis, including the propagation of systematic uncertainties.
The DR4 beams presented here were also used for the 2013-16 data that were part of the ACT DR5 maps (Naess et al. 2020).Finally, the pipeline presented here was used to obtain preliminary beams for the 2017-18 Advanced ACT data included in DR5.Looking forward, as we collect and analyze more ACT data and approach the cosmic variance limit, we will become increasingly sensitive to details of the instrument beams.While some details may change 33 , we expect to adopt similar methods as described here for analysis of the post-2016 ACT data.profile, n p, which is precisely when shrinkage is the most important.The term Cov(T ij , S ij ) in the numerator accounts for the fact that S and T are estimated from the same data, so if some elements of T are equal to elements of S then the two terms in the numerator cancel for those elements, and they do not affect the estimate of λ * .The denominator in Equation A5 ensures that if our choice of T is very different from S, then λ * will be small and C will be close to S. A poor choice of target matrix should not therefore negatively affect C.
If Equation A5 leads to a value of λ * that is greater than one, then we set λ * = 1, which means our final estimate of the covariance matrix is composed of only the target matrix.On the other hand, if the equation leads to a value of λ * that is less than zero, we set λ * = 0, and our final estimate of the covariance matrix is equal to the empirical covariance matrix.
Several common choices for the target matrix are listed in Table 2 of Schäfer & Strimmer (2005), along with simplified expressions for the associated shrinkage intensities λ * .In this analysis we use a "diagonal, unequal variance" target matrix, that is, a matrix whose diagonal elements are equal to the diagonal elements of our empirical estimate S of the covariance matrix and whose off-diagonal elements are zero: which is a Hankel transform of order zero, and where the last line follows from the identity J n (−z) = J n (z) for integer n.
In order to compute the harmonic transform of our beam profile, we evaluate the expression above separately for the three main terms in our beam profile fit: the core term (composed of the sum of basis functions), the scattering term, and the 1/θ 3 asymptotic term.The integrals for the core and scattering terms are computed numerically, but we derive an analytic expression for the integral of the 1/θ 3 term, shown below.
Given a fit amplitude α, the Hankel transform for the 1/θ 3 term may be written as The analytic expression we use for this integral is J 0 (θ ) θ 2 dθ = J 1 (θ ) − J 0 (θ ) where H n (x) is the Struve function.

C. BEAM SOLID ANGLES FOR DIFFERENT EFFECTIVE FREQUENCIES
To indicate the effect of the passbands on the beams, we tabulate the solid angles for select effective frequencies.The effective frequencies of the band centers are currently uncertain to approximately 2.4 GHz.However, since this is in part due to systematic errors in the measurements, the relative uncertainties are smaller.
Similar to the correction we make to the main beams for use with the CMB, in each case we take the beam to be B ( ) = B( ν RJ /ν S ), where ν RJ is the effective frequency for radiation with a Rayleigh-Jeans spectrum and ν S is the effective frequency for the source of interest (either CMB, synchrotron emission, dust, or the thermal Sunyaev-Zel'dovich (tSZ) effect).
As described in §4.4,for the beam analysis for DR4, the effective frequencies from Thornton et al. (2016) were used for the RJ-to-CMB beam spectral correction.Subsequently, the effective frequencies were re-computed for the foreground modeling for DR4, using improved passband data and upgraded code, as detailed in Appendix D of Choi et al. (2020).These updated frequencies are shown here in Table 5 and the corresponding beam solid angles are shown in the following tables.Considering the uncertainties on the passbands, the frequencies from Thornton et al. (2016) and Choi et al. (2020) are consistent.In addition, given that the uncertainty on the beams is subdominant in the power spectrum analysis, which of these effective frequencies one uses for the beam spectral correction does not have a significant effect on the results.Still, the solid angles for the CMB in Table 6 below are slightly different from those in Table 4.
For any particular season/region/detector array/frequency, the uncertainty on the solid angles is 2-5%.This includes both the instantaneous beam uncertainty and the uncertainty due to the region-dependent jitter correction.As is apparent here, the derived solid angle and passband are intimately connected.
In addition, the beam and the passband are coupled, an effect which we considered for the first time in detail in Madhavacheril et al. (2020), and described in Appendix A of said paper.In short, the beam shape evolves as a function of frequency across the passband, and so our customary assumption of separability of these two components does not hold to high precision.This beam-bandpass coupling was modeled as part of a systematic check of the DR4 results in Choi et al. (2020) and was found to not have a significant effect on the inferred cosmological parameters.

E. DR3 VS DR4 BEAMS
The beam transforms made publicly available for the DR3 and DR4 releases have been compared for each season, sky region, detector array, and frequency in common, as shown in Figure 21.Despite the changes in the analyses, the beam transforms are consistent.

Fig. 1 .
Fig. 1.-Example of a simulated planet observation used to characterize the mapping transfer function.This example is for S16 PA2 at 150 GHz.(Left) The input 2D beam model.(Center) One simulated observation.(Right)The difference between the "input" beam model (left) and the map of the "output" simulated observation (center).Note that in each case the color scale is a symmetrical log scale a in dB (with a linear threshold of −50 dB for the input and output maps and −45 dB for the difference map) and negative values are enclosed by parentheses.Here the target region where the atmosphere was estimated and removed during the map-making process was chosen to have a radius of 16 .The large-scale residuals seen are mostly constant within this target region of the map.The slight asymmetry in the beam (which may be safely ignored for DR4, as explained in §6.2) is due primarily to the positions of the telescope's optics tubes in the focal plane.aA symmetrical log scale is logarithmic in both the positive and negative directions, with a small linear range around zero to avoid having values tend to infinity (https://matplotlib.org/stable/api/scale_api.html#matplotlib.scale.SymmetricalLogScale).

Fig. 5 .
Fig. 5.-Example of the offset fit for the radial profile of one Uranus observation for S14 PA2 at 150 GHz.The points are found by taking the azimuthal average of the map for each radial bin.Uncertainties on these points are not computed.The curve is the best-fit model α/θ 3 + β, where in this case α = 0.08934 and β = 0.00012.

Fig. 6 .
Fig. 6.-Top: Measured radial beam profile (blue) for S15 PA2 at 150 GHz and the best-fit model (red), with the (narrow) red shaded region indicating the 1σ model uncertainty bounds, and residuals shown in the middle panel.Note that the bins are correlated.Bottom: The terms in the model for the beam.Here the final best-fit value of θ 1 , where the model for the beam changes, is 5.59 .The number of modes in the core fit (n mode ) is 13.

Fig
Fig. 7.-Maps of the main beam for S15 PA2 (top) and PA3 (bottom) at 150 GHz.In each case, the color scale is a symmetrical log scale in dB, with a linear threshold of −30 dB and negative values enclosed by parentheses.(Left) The weighted average of the Uranus observations used to characterize the beam.(Center) The model which is fit to the radial profile of the measured beam.(Right)The difference between the measured beam (left) and the beam model (center).Note that the color scale for the residuals is different from the other two maps.Both PA2 and PA3 are shown here, since PA2 serves as an example of a less elliptical beam, and PA3 serves as an example of a more elliptical one.The corresponding maps for S15 PA1 at 150 GHz resemble the ones shown for PA2 here, and the maps for S15 PA3 at 98 GHz are similar to the ones shown here for PA3 at 150 GHz, but broader.For each individual detector array, the residuals are fairly constant from one season to the next.As shown in Figure6for PA2, the azimuthal average of the residuals is consistent with zero, which is why the fit is successful, despite the residuals visible in the maps.These residuals are expected to have a quadrupole-like shape, since we know our beams are slightly elliptical, and the quadrupole is the dominant asymmetric azimuthal mode for an elliptical beam.

Fig. 8 .
Fig. 8.-The estimated instantaneous beam profiles for S15, for the three arrays (PA1, PA2 and PA3), and for both 98 GHz and 150 GHz for PA3.The shaded bands indicate the 1σ uncertainty bounds.The uncertainties are strongly correlated between radial bins.The S13, S14 and S16 beam profiles are similar, and are included in the public data release.At 150 GHz, the beam for PA3 is quite different from those for PA1 and PA2.This is because PA3 was not focused as well, due to its position, further off-axis, compared to PA1 and PA2.

Fig. 9 .
Fig. 9.-Instantaneous beam transforms and their uncertainties for the S15 data, for all the detector arrays.The uncertainties are strongly correlated between multipoles.For context when looking at this figure along with Figure 8, in the mapping from angle to multipole ( π/θ), 1 corresponds to 10800, 1.8 corresponds to 6000, and 10.8 corresponds to 1000.

Fig. 10 .
Fig.10.-Change in the beam transform (in %) due to applying different corrections for S15 PA2 at 150 GHz.For this plot, the beam transforms have all been normalized at = 1400.This corresponds roughly to our effective calibration scale, so any change in our beam at this multipole would be corrected by our subsequent calibration to Planck.

Fig. 11
Fig. 11.-Example of the residual pointing variance, estimated by fitting an effective beam to a point source.The top panels compare the flux from the best-fitting model (right) to the data (left), with the residuals (bottom left), and the reduced χ 2 as a function of pointing variance (bottom right).In each case, the data is high-pass filtered, and the colorbars are in units of µK.The best-fit variance for this source is 41 ± 93 arcsec 2 .

Fig. 12
Fig. 12.-Distribution of the residual pointing variance, V , and pointing variability, σ V , for the S15 PA2 beam at 150 GHz in the Deep56 region, resulting from fitting a pointing variance model to point sources.Notice here that the intrinsic scatter (measured with σ V ) is significantly larger than the measurement uncertainty on V .

Fig. 13 .
Fig. 13.-Uncertainties on the beam transform (in %) for S15 PA2 at 150 GHz, before (blue) and after (red) adding the additional modes to account for different surface rms values and different fit ranges for the beam offset levels.These uncertainties are included in Figure 9.

Fig. 14
Fig. 14.-Maps of the leakage beam in Qr (top) and Ur (bottom) for S15 PA2 at 150 GHz.In each case, the color scale is a symmetrical log scale in dB, with a linear threshold of −30 dB and negative values enclosed by parentheses.(Left) The weighted average of the Uranus observations used to characterize the beam.(Center) The model that is fit to the radial profile of the measured beam.(Right)The difference between the measured beam (left) and the beam model (center).Features in the residuals could arise due to physical effects, such as differential beam ellipticity between the two axes of a polarimeter or a polarization angle offset(Hu et al. 2003).The level of residuals seen in the last column has an insignificant effect on the results of the power spectrum analysis.As shown in Figure15, the azimuthal averages of the residuals are close to zero, which is why the fits are successful.

Fig. 15 .
Fig. 15.-Measured radial polarization profiles (blue) in Qr (top) and Ur (bottom) for S15 PA2 at 150 GHz and the model we fit to each profile (red), with the red shaded region indicating the 1σ model uncertainty bounds.The profiles have been normalized by the peak amplitude of the corresponding T beam (at the beam centroid, θ = 0).Note that the bins are correlated.

Fig. 16 .
Fig. 16.-The transform of the radial beam profiles in E (top) and B (bottom) for S15 PA2 at 150 GHz.These transforms are normalized by the amplitude of the corresponding T beam transform at = 0.The lower panels show the dominant independent modes of each transform's covariance matrix (in color) and the magnitude of the diagonals (the gray shaded regions).

Fig. 17
Fig. 17.-Measured temperature-to-polarization leakage in the main beams for S15 for all arrays, showing the temperature leaking into E-mode and B-mode polarization.Note the different y-axis ranges.For PA1 and PA2 the leakage values are within 1.5%, whereas for PA3, the leakage increases around 4000-6000, and reaches over 6% (4%) at ∼10,000 at 150 GHz (98 GHz).We attribute the higher leakage for PA3 to an imperfectly optimized horn design in this first generation of multichroic polarimeters.The uncertainties (lower panels) include the adjustments for model variability.
Fig. 17.-Measured temperature-to-polarization leakage in the main beams for S15 for all arrays, showing the temperature leaking into E-mode and B-mode polarization.Note the different y-axis ranges.For PA1 and PA2 the leakage values are within 1.5%, whereas for PA3, the leakage increases around 4000-6000, and reaches over 6% (4%) at ∼10,000 at 150 GHz (98 GHz).We attribute the higher leakage for PA3 to an imperfectly optimized horn design in this first generation of multichroic polarimeters.The uncertainties (lower panels) include the adjustments for model variability.

Fig. 18 .
Fig. 18.-Sidelobes in the PA1, PA2, and PA3 detectors at 150 GHz.(Left) Maps of the polarized sidelobes, obtained by stacking observations of Saturn.The main part of the beam at the center of each image is masked.In each the grayscale is linear, indicating the sidelobe amplitude in the range from −0.002 (black) to +0.001 (white) relative to the main beam peak, a with positive (negative) numbers corresponding to polarization parallel (perpendicular) to the radial direction.The complementary polarization leakage (corresponding to T B leakage) is smaller and not shown in the maps, but it is also estimated.(Right) The effect of the sidelobes on our measurements, expressed as beam transfer functions B T →E and B T →B .Note that the scales on the right differ for all three arrays.The sidelobes are projected out of the time-ordered data prior to mapping, so the effective (residual) leakage is what affects our spectra.Since the sidelobes for PA3 are significantly weaker than PA1 and PA2, we do not compute residuals for PA3.No sidelobes are seen for PA3 at 98 GHz.a Saturn's brightness appears to compress the detector gain by a few percent.The effect can be ignored here because it is a small correction on the sidelobe leakage, which itself is a small effect.

Fig. 19 .
Fig. 19.-An example of the size of various effects on the inferred spectra, compared to our final uncertainties.Here we used the leakage in the main beam for S15 PA2 at 150 GHz, the residual sidelobe leakage for PA2 at 150 GHz, simulated leakage due to asymmetry for S15 PA2 at 150 GHz, and the uncertainty in the main beam for S15 PA2 deep56 at 150 GHz.The best fit ΛCDM plus foreground model for the deep 150 × 150 GHz spectra for ACT only from Choi et al. (2020) was used to compute the resulting change in D T T , D T E , and D EE due to these elements.These changes are then divided by the corresponding uncertainties σ T T , σ T E , and σ EE on the final binned deep spectra.At larger multipoles, eventually the uncertainty on the main beam would come to dominate.Note that we correct for the main beam leakage and the residual sidelobe leakage in the power spectrum likelihood, as explained in §5.3.The asymmetry leakage (simulated using an allsky beam convolution code, as mentioned in §6.2) is negligible for all detector arrays.Sudden changes around = 2000, as can be seen in the top panel, for example, are due to a change in the bin width at this scale.

Fig. 20 .
Fig. 20.-The γ factors encoding the -dependent T-to-P leakage, used in the power spectrum likelihood for both the deep and wide spectra at 150 GHz (top) and 98 GHz (bottom).The γ factors for 150 GHz include both the main beam leakage from §5.1 and the residual sidelobe leakage from §5.2.Since no sidelobes are seen at 98 GHz, those γ factors consist of only the main beam leakage.

Fig. 21 .
Fig. 21.-Ratio of the beam transforms for DR3 and DR4 at 150 GHz, for S13 (top) and S14 (bottom).The shaded bands indicate the 1σ uncertainty bounds, determined using the maximum uncertainty of the two transforms being compared.For these plots, the transforms have been normalized at = 1400, which corresponds roughly to the effective calibration scale.

TABLE 1
Summary of Uranus observations -number used/total.

TABLE 2
Instantaneous beam solid angles, gains, and FWHMs