FIVE-YEAR WILKINSON MICROWAVE ANISOTROPY PROBE^* OBSERVATIONS: BEAM MAPS AND WINDOW FUNCTIONS

The fundamental beam data are Jupiter measurements extracted from the mission's TOD archive. These measurements are differential like all WMAP data, but the presence of Jupiter in only one of the two beams of a given DA means that they are effectively single-dish, after subtraction of a differential sky background. Five-year WMAP full-sky maps, which omit Jupiter, are used to estimate this background. Also, measurements affected by a bright source other than Jupiter in either set of optics are omitted. The apparent measured Jupiter temperatures are scaled to a standardized distance of 5.2 AU and are binned by the position of Jupiter in a coordinate system attached to the spacecraft. Planetary coordinates are obtained from the Jet Propulsion Laboratory ephemeris DE200¹⁴ (Standish 1990).

2.2.1. Software Structure

2.2.2. Parameters of Fit

The coordinates and Euler angles of the WMAP optical components are furnished to DADRA as constants, as are the measured mirror shapes and the beam patterns of the feed horns. These constants are results of the preflight structural thermal optical (STOP) performance analysis (Page et al. 2003b). The fit is done by varying a set of mirror distortion parameters, which are defined as small-amplitude Fourier or Bessel modes added to the input shape. The emphasis in the fitting is on the primary mirror, because in preflight laboratory measurements, the shape of the secondary was less susceptible to temperature changes. Also, the illumination pattern occupies a much larger proportion of the primary mirror than of the secondary (Figure 2 of Page et al. 2003b). For these reasons, the perturbation space chosen to characterize the primary mirror encompasses more modes than that for the secondary. Moreover, Fourier modes in x and y are a natural choice because the morphology of the distortions is dominated by a rectangular frame that is part of the backing structure. On the other hand, the mount of the secondary is cylindrically symmetric with a more extensive backing structure, consistent with fitting to a restricted set of Bessel modes.

The primary distortions have the form

$$\begin{eqnarray} \delta Z & = & \sum C(k_x,k_y) \cos \left({2 \pi k_x x\over L} + {2 \pi k_y y \over L} \right)\nonumber \\ & & + \sum S(k_x,k_y) \sin \left({2 \pi k_x x\over L} + {2 \pi k_y y \over L} \right)\nonumber \\ & & + \sum _{\mathrm{special}} A_{\mathrm{sp}} F_{\mathrm{sp}}(x,y).\nonumber \end{eqnarray}$$

The distortion modes are characterized by the spatial frequency indexes k_x and k_y, where k = 1 corresponds to L = 280 cm, which is twice the width of the mirror. The reason for setting L = 280 cm rather than 140 cm is to remove the requirement that the solution be periodic on the circumscribed square; in particular, an approximate tilt of the 140 cm mirror is expressed naturally as half a sine wave with a 280 cm period. Some of the Fourier modes specified in this way are redundant, so the set that is used in the fit is culled to contain only even values of k_x and k_y, with k_x > 0. A few additional modes with k_x and k_y of ±1 are the ones that represent tilt. Special distortion modes, F_sp(x, y), are also allowed, with amplitudes A_sp. One special mode is a scalar offset of the whole mirror surface. A second special mode is simply a map of the mirror surface as measured preflight; in practice, this mode plays no significant role in the fit. The phase and strength of the modes are specified by the sine and cosine amplitudes, S(k_x, k_y) and C(k_x, k_y), and the amplitudes of the special modes, A_sp. Below, we make use of the power per mode, P(f) = C(k_x, k_y)² + S(k_x, k_y)², where f is the spatial frequency in cm⁻¹ and f² = (k_x/280)² + (k_y/280)².

The three-year fit, carried out for the A side only, reached a maximum spatial frequency index of k_max = 12, corresponding to a distortion wavelength of 280/12 ∼ 23 cm. By contrast, the five-year fit reaches k_max = 24, or ∼12 cm, for both the A and the B sides. The number of Fourier modes goes as k²_max, so the five-year fits include ∼400 modes on each primary as compared to ∼100 modes in previous analyses. This extension in k space means that the primary mirror distortions are fitted nearly to the surface correlation length measured in the laboratory under cold conditions, i.e., 9.3 cm on the A side and 11.4 cm on the B side (Page et al. 2003b).

The secondary mirror distortions are described by Bessel functions:

$$\begin{equation} \delta Z = \sum _{n,k} J_n(\nu _{n,k}\rho /L)[ C_{n,k} \cos (n \phi) + S_{n,k} \sin (n \phi) ], \end{equation} \tag{ 1 }$$

where ρ and ϕ are cylindrical coordinates and L is the radius of the mirror. The kth alternating zero of J_n, J'_n is denoted ν_n,k; zeroes of J_n have k odd, and zeroes of J'_n have k even. Inclusion of the zeroes of J'_n removes any constraint on the edge of the mirror. The fitted parameters are C_n,k for n ⩾ 0 and S_n,k for n ⩾ 1. Two different pairs of (n_max, k_max) are used: (1, 3), resulting in nine parameters, or (2, 6), resulting in 30 parameters.

2.2.3. Fitting Method

The optimization is done by a modified conjugate gradient method, a deterministic descent into a χ² valley. Avoidance of local χ² minima is attempted by exploiting the Fourier description of the primary mirror distortions. The largest scale distortions are fitted first, and each result is used as a starting point for the next fit, in which finer scale modes are included. At each stage, distortions at the large scales that have already been fitted are not held constant, but rather, they are refitted together with the small-scale distortions that are newly included.

As defined above, the primary mirror modes do not compose an orthogonal basis. However, orthogonality is desirable in order to make the fit as efficient as possible. Consequently, the primary mirror modes are not fitted directly, but first are orthogonalized with respect to the area inside the circular boundary of the mirror, using a modified Gram–Schmidt method. When the k_max of the included Fourier modes is increased, a new orthogonalization is performed that generates a completely new set of linear combinations of the Fourier and special modes. The conjugate gradient algorithm therefore navigates in a space consisting of two groups of parameters: (1) the amplitudes of the orthogonalized primary mirror modes and (2) the amplitudes of the Bessel modes on the secondary.

Each time χ² is calculated, two types of adjustment are made to the model beams. First, the pointing in the coordinate frame attached to the spacecraft is matched to that implied by the Jupiter observations. Second, the peak sensitivity of the model is scaled to match the peak observed Jupiter temperature. The pointing adjustment may be done for the 10 beams as a group, without altering their mutual displacements, or it may be done for each beam separately. However, the peak scaling is always done separately for each beam.

These adjustments are not parameters of the fit, because corresponding dimensions of χ² space do not exist. Rather, their purpose is to absorb errors in the input coordinates and angles of optical components and prevent them from being projected into the mirror distortions. Ideally, this problem would be avoided by solving a more complicated problem, i.e., by directly fitting the mechanical parameters of the WMAP components. However, a full set of mechanical parameters would include many degeneracies with respect to the beam morphology. By limiting the parameter set to ∼100–400 mirror perturbation modes and ignoring "nuisance" information, we converge on acceptable values of χ².

For a given fitting run, either the primary or secondary mirror parameters can be held constant at the starting values. The DA microwave frequencies can also be fitted as parameters, but are held constant in practice, since they are accurately determined.

The fits for side A begin with the inherited three-year solution (Jarosik et al. 2007). The final run for each value of k_max is listed in Table 1. An indication of the quality of fit after each step is given by the χ²_ν column in the table.

Table 1. Abbreviated Fit History for Side A

Primary kmax	Primary Modes	Secondary Modes	χ2νa
12	122	30	1.15
14	158	30b	1.17
16	206	30b	1.11
18	262	30b	1.11
20	326	30b	1.12
22	386	30b	1.11
24	450	30b	1.05

Notes. ^aχ²_ν is approximate and indicates the progress of the fit. Residuals are shown in Figure 2. ^bNot refitted.

Download table as:  ASCIITypeset image

The B side of the instrument is characterized by an overall shift of the 10 beam pointings in relation to their preflight positions, by ∼01. This shift complicates the fitting strategy. The overall fitting history for side B comprises several different sequences of fits. The most important sequence, leading to the adopted beams, is similar to the A-side fitting sequence, in that k_max is increased in stages, with the secondary mirror distortions held constant after being fitted early in the sequence (Table 2). In other sequences, a different form of the secondary was tried, the floating shift in elevation and azimuth was disabled, or the secondary alone was fitted from various initial conditions. None of these variations improved χ² for the resulting beams, as compared to the adopted fitting sequence.

Table 2. Abbreviated Fit History for Side B

Primary kmax	Primary Modes	Secondary Modes	χ2νa
4	22	0	2.50
6	38	9	1.26
6	38b	30	1.27
8	58	30b	1.22
10	90	30b	1.17
12	122	30b	1.17
14	158	30b	1.16
16	205	30b	1.15
18	262	30b	1.15
20	326	30b	1.15
22	386	30b	1.11
24	450	30b	1.07

Notes. ^aχ²_ν is approximate and indicates the progress of the fit. Residuals are shown in Figure 2. ^bNot refitted.

Download table as:  ASCIITypeset image

The fitted beams and residuals for sides A and B are shown in Figure 2, which can be compared to Figure 9 of Jarosik et al. (2007).

Zoom In Zoom Out Reset image size

2.2.4. Instrument Parameter Results

False-color renditions of the final A and B side mirror surface fits are shown in Figure 3. Two natural length scales for the surface of the primary mirrors are 0.5 cm for a hexagonal mesh that composes one layer of each mirror, and 30 μm for the correlation length of the reflector surface roughening, which was done to diffuse visible solar radiation (Page et al. 2003b). Both of these length scales are too small to be probed either by the direct fitting of the main beam or by the sidelobe observations of the Moon. The main feature of each fitted primary mirror figure is the backing structure, dominated by members that form a rectangular frame near the center of the mirror. In the center part of this rectangle, the primary mirror appears to be depressed by ∼0.5–1 mm. Also seen are hints of the stiffening lugs near the edge of each backing structure. The rms distortions of each primary mirror model are σ_z = 0.023 cm and 0.022 cm for the A and B sides, respectively. preflight cold-measured values on the real mirrors, as extrapolated to the flight temperature of 70 K, were σ_z = 0.023 cm and 0.024 cm, respectively (Page et al. 2003b).

Zoom In Zoom Out Reset image size

For the A side, the measured centroids of Jupiter beam data are displaced by ∼003 ±  002 from the nominal preflight beam positions on the sky, where the error term is the 1σ scatter among beams. However, for the B side, the corresponding displacement is 013 ±  003. The χ² computation in the beam-fitting algorithm includes a floating elevation-azimuth adjustment that is intended to soak up such discrepancies without converting them into parameters of the fit. For the final adopted beam models, the floating displacement amounts to 009 in combined elevation and azimuth for side A, and 021 for side B. The difference between these two values agrees with the raw pointing difference between the beams on each of the two sides.

We emphasize that these displacements are unrelated to the estimated pointing errors of <10'' in the WMAP TOD (Jarosik et al. 2007). The A-side and B-side boresight vectors are accurately determined from flight data as part of the TOD processing, and are not influenced either by the beam fitting or by preflight predictions. The WMAP pointing model is described in Limon et al. (2008).

The mirrors are constrained only where they are substantially illuminated by the feed horns; see Figure 2 of Page et al. (2003b). Thus, for example, the secondary mirror for the B side appears as a bullseye partly because the fit is only constrained in the center (Figure 3). In actuality, the fitted shape consists mostly of a tilt of ∼025; however, apparent mirror tilts reflected in the fitted parameters are difficult to interpret because they are coupled to the floating elevation-azimuth offsets.

The polarization characteristics of the main beam models are consistent with previous results. The morphology of the co- and cross-polar components of both the A- and the B-side models is similar to that found for the A side in the three-year analysis. Similarly, cross-polar suppression, as calculated from peak model values, is within ∼0–2 dB of previously reported A-side values, depending on the DA. However, for Q, V, and W bands, the polarization isolation of the orthomode transducer (OMT) dominates the end-to-end cross-polar response, which was measured preflight (Jarosik et al. 2007).

The modeling of primary mirror distortions with k as high as 24, which affect the beam at relatively wide angles, raises the possibility of comparing DADRA-computed main beams to the innermost parts of the far sidelobe patterns, which are obtained from observations of the Moon (Barnes et al. 2003).

The beams at angles greater than 5°–10° from each boresight may be affected by unmodeled primary mirror distortions with 24 < k ≲ 250. However, extending the models to fit these distortions directly is computationally unmanageable, because the required number of Fourier modes is of order 10⁴, which is ∼ 100 times the number of modes in our normal fits. Nevertheless, the observed sidelobe data provide a constraint on the contribution of such modes to the primary mirror surface shapes, and hence to the main beams. To apply this constraint, we need an appropriate choice of sidelobe data together with a method for extending each main beam model to the low gain levels outside the main beam radius.

The choice of sidelobe data is important, because the full-sky sidelobe patterns are dominated by features that are not captured in the beam models. The WMAP sidelobe patterns are depicted in Figure 2 of Barnes et al. (2003). The most vivid features are formed by reflections from the radiator panels, by diffraction around the edges of the primary mirrors, and by reflection from the focal plane assembly. Despite their striking appearance in the referenced figure, these features are at least 40–60 dB below the peak gain of each main beam; however, they must still be excluded from the comparison.

Only a small region of each sidelobe map, a smooth region in the "shadow" of the primary mirror, is suitable for comparison with the main beams. Although we choose the boundary of this region as conservatively as possible, the choice is subjective and a source of systematic error. Essentially, we draw the boundary according to a combination of sidelobe morphology and angular proximity to the main beam. However, the morphological criterion by nature cannot exclude an extraneous component of sidelobe response that happens to be smooth, such as might arise from a diffuse reflection off the top of the structure that holds the feeds. For this reason, the sidelobe patches chosen for the comparison are regarded as upper limits.

Additionally, these radiometric observations cannot be calibrated as well as the CMB data, primarily because they were taken when the spacecraft was thermally unsettled, during the phasing loops between the Earth and the Moon. The instrumental gain in this part of the mission is estimated to be known to ≲10%. A total of 4.3 days of Moon data were obtained covering ∼1.2π sr of sidelobe area. The calibration standard for the Moon observations is the COBE DMR model of lunar microwave emission as a function of the phase angle (Bennett et al. 1992).

In order to achieve a comparison with the Moon observations, the spatial frequency of the modeled primary mirror distortions is pushed to as high a value as possible. Figure 4 (lower left panel) shows radial profiles of the symmetrized main beam models in comparison to the sidelobe sensitivity pattern for one example beam, V2 on the A side. The sensitivity profile at intermediate angles of 15–25 is directly related to k_max of the fit, as seen also in the grayscale images of the model beams (top row). A natural hypothesis is that an even greater increase in k_max might give a model that joins smoothly to the Moon data. Unfortunately, the direct fitting algorithm cannot accommodate an indefinite increase in k_max, because the number of Fourier modes goes as k²_max.

Zoom In Zoom Out Reset image size

To cope with this difficulty, an approximation is used for Fourier modes with 25 ⩽ k ⩽ 250. Power spectra, P, of the fitted A-side and B-side primary mirror distortions are shown in Figure 5 as a function of spatial frequency, f = k/280 cm⁻¹. The form of the power spectrum expected from ground-based measurements of the mirrors is also shown, under the assumption of a Gaussian two-point correlation function (Page et al. 2003b). As the spatial frequency increases, P decreases. Our approach is to extrapolate the power spectrum P(k) to smaller scales assuming a power-law form, P ∝ (1/k)^α, with 3 ≲ α ≲ 6. Random phases are used to convert the extrapolated spectrum to sine and cosine amplitudes.

Zoom In Zoom Out Reset image size

A grid of beam models is assembled as a function of two variables: the slope α, and the random number seed s that selects the phases used for the extrapolated Fourier modes. Separately for each s, a χ² minimization is used to fit the beam models to the Moon sidelobe data as a function of the slope α. The sidelobe maps are in HEALPix format (Gorski et al. 2005) with $\texttt {nside}=512$, and the beam models are resampled onto the HEALPix grid. The value of χ² is computed from measured and predicted gains, g, as ∑(g_pred − g_Moon)²/σ²_Moon, where the sum is over pixels in the region of overlap between the beam model and the Moon data. The fit uses all of the Q, V, and W DAs together. For each side, the fitted values of α from five values of the seed s are averaged to get the adopted slope. Using this slope value, five new beam models are computed with the original phases, and these models are averaged. The result is termed the augmented beam model.

Figure 6 shows radial profiles of the A-side augmented beam models compared to the fitted subset of Moon data. The K1 and Ka1 beams are not used in the fit because of the relatively strong diffuse light that is seen as the bright profile in the top two panels of Figures 6 and 7. This component, which may result from reflection off the focal plane assembly, is seen in Figure 2 of Barnes et al. (2003), where it appears as a haze in the region above the main beam. In addition, the Q, V, and W bands appear subject to systematic errors depending on the individual DA.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

For CMB analysis, the main implication of the augmented models is an increase in main beam solid angle as compared to the ordinary fitted models with k_max = 24. The effect on the "tail" part of the main beam model is illustrated in the bottom right panel of Figure 4. To compute accurately the effect that the augmented models would have on the symmetrized beam profiles used to compute b_ℓ, a hybridization with Jupiter data is required, as described below (Section 3.1). If augmented rather than ordinary k_max = 24 beam models are used in the hybridization, the resulting increase in main beam solid angle is just ∼0.1%–0.3% depending on DA.

There are two arguments for treating this main-beam solid angle increase as an upper limit. One argument is the one already made above, namely, that the Moon data in all bands may include a diffuse reflected component in addition to the extended main beam response, as is apparent for K1 and Ka1, and which our procedure cannot exclude. The other argument invokes the thermal nature of the CMB power spectrum, which requires that the power spectrum be the same in all microwave frequency bands. Section 3.4 below shows that the ordinary models with k_max = 24 maintain a tighter consistency of the CMB C_ℓ across Q, V, and W than do the augmented models. Consequently, the adopted five-year beam transfer functions incorporate only the ordinary models, whereas the augmented models are used to characterize the innermost part of the sidelobe response (Section 3.2).

To mitigate sensitivity of the window functions to observational noise, we use a beam hybridization technique similar to that employed in the three-year data analysis (Jarosik et al. 2007). In this method, a hybrid beam is constructed for each DA on each side by combining Jupiter observations with the physical optics models. Jupiter observations are used in the central portions of the beam where S/N is high. Model points are substituted for the data in the outlying regions of low signal, called the tail.

The five-year analysis differs from that of Jarosik et al. (2007) in the method for choosing the hybridization threshold, B_thresh, which defines the tail region. In the three-year analysis, the threshold for each DA was chosen to replace noise-dominated parts of the beam with model values, in order to facilitate fitting the beam profile with a smooth function. However, the improvements in the five-year data and modeling open the possibility of extrapolating the main beam to wider angles and subsuming more of the full-sky beam solid angle into the main-beam treatment, rather than the sidelobe pattern. As a result, lower-signal parts of the beam are included in the beam transfer functions, and we require an explicit optimization of S/N in the hybrid beams.

The effect of the beam tail on science occurs through the normalization of the CMB power spectrum, C_ℓ. An increase in main beam solid angle raises the beam-corrected C_ℓ by a constant factor for ℓ ≳ 100. Hence, the error in the solid angle is a convenient indicator of the error induced in the high-ℓ part of the CMB power spectrum via the hybridization. The solid angle error is therefore a natural figure of merit for optimizing B_thresh.

A grid of simulations is run that evaluates solid angle error in the hybrid beam as a function of threshold level. Let the true solid angle of a given beam be Ω, and the solid angle of the hybrid beam be Ω_h(t), where for conciseness we use t to stand for B_thresh. Then, Ω_h(t) = Ω_d(t) +  Ω_m(t), where Ω_d is the portion taken from data, and Ω_m is the portion taken from the model. Another way of decomposing Ω_h is into a true solid angle and error terms, i.e., Ω_h = Ω + e_d(t) + e_m(t), where e_d(t) is the error in the data portion for threshold t, and e_m(t) is the error in the model portion. The model error can be parametrized as a fraction of the model solid angle, such that e_m(t) = aΩ_m(t). If e_m is uncorrelated with e_d, then the fractional variance in Ω_h is

$$\begin{equation} \mathrm{Var}\lbrace \Omega _h(t)/\Omega \rbrace = \mathrm{Var}\lbrace e_d(t)/\Omega \rbrace + \mathrm{Var}\lbrace a\rbrace (\Omega _m(t)/\Omega)^2, \end{equation} \tag{ 2 }$$

and the hybridization threshold is chosen to be the value minimizing this variance.

Essentially, the variable a in the above discussion is a scaling error that is common to all of the model points incorporated in the hybrid beam. A conservative method of estimating systematic error is to assume that it is of the same order as the quantity estimated. In the above formulation, we represent this estimate by setting Var{a} = 1.

Figure 7 shows the fractional error in the hybrid beam solid angle for the V2 beam on the A side, as a function of B_thresh. To avoid a selection bias, B_thresh is referred to the model rather than the data. The errors contributed by the data portion and the model portion are shown along with the total error, which is computed using Equation (2) with Var{a} = 1. For the data, the solid-angle error is obtained as a function of B_thresh from 100 Monte Carlo simulations in which model input beams are combined with white noise appropriate to the Jupiter data for each DA. The contribution of the data portion increases with lowered threshold as more of the noisy data are included in forming the hybrid beam. Conversely, the contribution of the tail portion decreases with lowered threshold as less of the model is included. The adopted B_thresh values for the five-year analysis are obtained from the locations of minimum total error in similar plots made for all the A- and B-side beams. These values are shown in Table 3 together with the three-year equivalents. The five-year thresholds are lower than the three-year thresholds by some 5–10 dB, depending on DA. Thus, we use significantly more of the data than we have in the past.

For the five-year analysis, changes have been made in the sidelobe sensitivity patterns that are distributed as part of the data release. These patterns are in linear units of gain relative to isotropic, and they are assembled from several types of data, including the Moon observations described above. In previous analyses, the conversion of observed Moon brightness to gain assumed a temperature of 175 K for the Moon. However, in the K band the resulting gains were divided by 1.3 to match adjacent pieces of the sidelobe pattern that were derived from other data (Barnes et al. 2003).

For the five-year analysis, the calibration has been improved by integrating the COBE DMR model of lunar microwave emission (Bennett et al. 1992) over the spacecraft and Moon ephemerides. This procedure lowers the gain values in the Moon-derived parts of the sidelobe pattern by factors of 1.35, 1.39, 1.45, and 1.51 for the Ka, Q, V, and W bands, respectively, while confirming the earlier adopted calibration of the K band.

Additionally, some ground-based beam measurements, used within ∼10° of each boresight, have been replaced with the augmented main beam models described above (Section 2.3). The ground-based data were taken under ambient conditions in the Goddard Electromagnetic Anechoic Chamber (GEMAC), where the primary mirror distortions are different from those under flight conditions. This replacement affects at most ∼0.5% of the sky for any given DA.

In the sidelobe response patterns, the area inside the main beams is set to zero. This area is expanded for the five-year analysis (Section 3.3).

The differential signals tabulated in the WMAP TOD archive are corrected for sidelobe contamination. The overall effect can be summarized in one number for each DA, called the sidelobe recalibration factor, which is the factor by which the correction changes the instrumental gain (since the dipole is detected in the sidelobes along with other sources). A sidelobe recalibration factor of unity means that the sidelobe response is zero. In the three-year analysis, these factors differed from unity by 0.3%–1.5% (Jarosik et al. 2007), and for the five-year analysis, they differ from unity by 0.05%–1.4% (Hinshaw et al. 2009). The decrease in the sidelobe correction is caused by the increased main beam area together with the revised calibration of the Moon data.

If the WMAP beam patterns could be well sampled in flight over 4π steradians, then the distinction between main beams and sidelobes would be arbitrary. However, the two regimes are measured by different methods, they are treated differently in the beam analysis, and they are applied differently in the WMAP data reduction, so that some reasonable boundary needs to be drawn. We do so by using the beam models for k_max =  24 to define a transition radius centered on each boresight. With the fitted mirror distortions, a separate DADRA computation is done to extend each beam model into a wide angular field, 11°–13° on a side. Cumulative beam solid angle is computed as a function of radius, and the radius containing 99.9% of the solid angle in the model is determined. The transition radius is then fixed at a round number encompassing the computed radii for both the A and the B sides. The adopted values are 70, 55, 50, 40, and 35 for the bands K, Ka, Q, V, and W, respectively. Compared to the one-year and three-year analyses, the transition radius is increased, as shown in Table 3.

This expansion of the main beam region has the useful consequence of mitigating the sidelobe correction. However, the main beam now includes lower S/N observations, to which the main beam solid angle is sensitive. Similarly, the profile-fitting algorithms of the first- and third-year analyses can no longer be used as previously implemented, because the fitting function is difficult to constrain over the entirety of the new main-beam radius.

In the one- and three-year analyses of azimuthally symmetrized beams, the radial profiles were modeled with basis functions of the form

$$\begin{equation} H_{2n}\left(\frac{\theta }{\sigma _h}\right)\! \exp\! \left(-\frac{\theta ^2}{2\sigma _h^2}\right), \end{equation} \tag{ 3 }$$

where θ is the angle from the beam center, H_2n(x) is a Hermite polynomial of even order, and σ_h determines the width of the Gaussian. The first of these basis functions is a pure Gaussian, which is a good fit to the main lobe of the beam, both theoretically and in reality. The other basis functions parameterize the deviations from Gaussianity. The Hermite fit is limited to the well characterized part of the beam, within a given cutoff angle θ_c of the beam peak. In the five-year analysis, θ_c is increased as compared to the three-year value because of improvements in the data and analysis described above. However, this presents two problems with the Hermite polynomials.

First, the basis functions do not extend far enough in θ. Since the Hermite polynomial H_2n(x) is a polynomial of order 2n, and the function x²ⁿexp(−x²/2) has peaks at $x=\pm \sqrt{2n}$, the basis functions of order 2n extend to $\theta \approx \sigma _h \sqrt{2n}$, beyond which the function is exponentially suppressed. Because this angle increases only with the square root of the order, many basis functions are required to cover the required domain, e.g., θ ≲ 40σ_h in the W band. Second, numerical problems arise in computing the Hermite polynomials of higher order than ∼150. The combination of these two problems rules out the use of Hermite functions in the five-year analysis.

The use of a fitting basis provides a smooth fit through noisy portions of data, and also provides a convenient mechanism for the derivation of a beam covariance matrix via the formal statistical errors in the fit. Because of these benefits, a number of possible sets of basis functions have been explored for the five-year beam data using simulations.

Beam profile simulations test the accuracy to which various sets of basis functions reproduce the known input beams and window functions. A variety of noisy simulated beams is constructed, then fitted. The simulations include pure DADRA models as well as hybrids of two DADRA models. In the case of hybrids, one beam model with noise added is used to represent the Jupiter data, and another model without noise is used for the beam tails. The hybridization thresholds (Section 3.1) are also varied, as is the overall scaling of the beam tails. The result of this testing is that ultimately, no one set of basis functions recovers the input beam solid angle and window function. One of the impediments seems to be the nature of the five-year hybrid itself, which is noisy at intermediate angular scales within the transition radius. Functional fit residuals in that region typically cause a bias of ∼0.5% in the recovered solid angle.

Owing to this difficulty, the method of basis function fitting is not used in the five-year analysis. The adopted hybrid beam profiles are left in the radially binned form, in spite of the noise that remains at low gain levels. Simulations show better recovery of solid angle from the resulting beams than from any of the attempted basis function fits.

The applicability of symmetrized beam profiles depends on the degree to which the assumption of azimuthal beam symmetry is justified. The WMAP scan strategy mitigates the effect of noncircularity in the beams by sampling most sky pixels over a wide range of azimuth angles. The effects of residual noncircularity are of potential importance for CMB power spectrum analysis primarily in the Q band at ℓ ≳ 500, where the effect can reach several percent in C_ℓ; however, Q-band data have low statistical weight in this ℓ range and are not used in the TT power spectrum analysis (Nolta et al. 2009). Mathematical details, together with plots of the relevant correction to w_ℓ, are given in Appendix B of Hinshaw et al. (2007).

To compute beam transfer functions, the radial distance, θ, from the beam centroid is computed for each of the data points in the two-dimensional A- and B-side hybrid beams. Radial profiles b^S(θ) are constructed by sorting all hybrid data into equally spaced bins of Δθ = 025 in width, and taking the mean of each bin. The radial profile only extends out to the transition radius. The beam transfer functions are evaluated using the Legendre transform:

$$\begin{equation} B_\ell = \Omega _B b_\ell = 2 \pi \int b^S(\theta) P_\ell (cos\,\theta) d \cos \theta. \end{equation} \tag{ 4 }$$

Numerically, the integration is performed by summing over rectangular bins of Δθ = 025.

As described above (Section 2.3), the sidelobe data of the Moon motivate an attempt to augment the fitted distortions of the primary mirror (k_max = 24) with random distortions that are extrapolated to finer spatial scales, i.e., k_max = 250. These added distortions affect the hybridized beam through their effect on the outermost, low-gain part of the main beam model. One way of testing the effect of systematic error in the modeled beam tail is to rescale the extrapolated distortion amplitudes up or down as a group, with 100% correlation. The resulting distortion amplitudes are used to compute new model beams, which are processed through hybridization with flight data and transformation to b_ℓ.

The effect of this type of distortion rescaling on the beam transfer functions is shown in the bottom panel of Figure 8 for flight data, and Figure 9 for a noiseless simulation. Rescaling the added mirror distortions changes the slope of b_ℓ between ℓ = 0 and ℓ ∼ 100, while shifting b_ℓ up and down for ℓ ≳ 100. Scale factors in the range 0–2 result in a ∼0.3% total range of variation in the high-ℓ value of b_ℓ.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

These scalings of the added distortions have been tested for their ultimate effect on CMB power spectra. Figure 10 shows the results of this test for scale factors of 2, 1, 0.5, 0.1, and 0, respectively. Each panel shows a mean of year-by-year CMB cross-power spectra computed from the five-year data set for each of the 8 WMAP DAs Q1–W4. The spectra are all computed using the MASTER estimator, and they are corrected for b_ℓ derived from augmented beams, characterized by distortion scale factors as indicated. For plotting, each such power spectrum is divided by the final MASTER power spectrum from the three-year WMAP analysis.¹⁶ In each case, a contribution from unresolved point sources is fitted and removed. In general, the result is that lower values of the scale factor give better consistency between microwave frequency bands for the CMB. Indeed, on this criterion, there is no clear reason to prefer a scale factor greater than zero.

Zoom In Zoom Out Reset image size

As a result, the way the extrapolated random-phase mirror distortions are handled is by omitting them from the adopted beams and b_ℓ, while the actual fitted mirror distortions with k_max ⩽ 24 are retained, via the model part of the hybrid beam. However, we incorporate into the error analysis an estimate of the systematic error in the faint part of the model, by assuming that this error is of the same order as the adopted model, just as we do for optimizing B_thresh (Section 3.1). Monte Carlo experiments done on the primary mirror distortions suggest that this 100% scaling error is conservative.

Combined errors in b_ℓ, which arise both from observational scatter in the Jupiter measurements and from the scaling error in the model, are estimated using Monte Carlo simulations of the hybridization. The DADRA flight models are used to represent the true input beams. These models are sampled to match the observed beam positions in the five-year flight archive. Based on the chosen hybridization threshold, white noise is added to the model for the points that would be taken from Jupiter observations in the actual analysis. The model points that are substituted for the low-gain tail are multiplied by a common, normally distributed scaling factor with μ = σ = 1. Because the beam modeling is expected to be common mode across frequencies, the same scaling factor is used for every DA in a given Monte Carlo realization. A total of 5000 beam realizations is computed, each comprising all ten DAs. The beam transfer function is computed for each beam realization as described above for flight data, and the standard deviation of the realizations at each ℓ is used as the diagonal of the covariance matrix. Figure 11 shows the resulting error bands for each DA, compared to the adopted errors from the three-year analysis.

Zoom In Zoom Out Reset image size

The chosen radial bin size of 025 is the smallest width permitting all the bins to be populated. However, the bin sizes 05 and 075 have also been tested, in order to ensure that window function results are not affected. For the ℓ range of interest, the bin size contributes negligibly to the smearing of the beam, even in the W band, and has no effect on the estimated error bars.

Monte Carlo realizations of b_ℓ for each DA are used in estimating the full covariance matrix of the coadded TT power spectrum in V and W. For CMB analysis, the error in the error is important. From several independent Monte Carlo runs of 5000 realizations apiece, the combined VW window function covariance has an error of ∼±3% in the diagonal elements (Hinshaw et al. 2009). The beam-related errors in the coadded TT power spectrum are shown as a function of ℓ in Figure 12.

Zoom In Zoom Out Reset image size

The comparison of three-year and five-year beam transfer functions is shown in Figure 13. Plotted is the relative change in b_ℓ, i.e., (b^3yr_ℓ − b^5yr_ℓ)/b^3yr_ℓ. Each beam is integrated out to the five-year transition radius; for this plot, the three-year beams are extended using the three-year sidelobe response patterns. Also plotted are ±1σ errors of b^3yr_ℓ. The three-year and five-year b_ℓ are consistent, with ∼1σ changes in V2 and W1–W4. For these DAs, the change has the form of a plateau for ℓ ≳ 200, reflecting an increase in the main beam solid angle for five years. This increase raises B₀ relative to B₂₀₀. The differences plotted in Figure 13 are taken with a sign convention reflecting the difference in the final power spectrum. Thus Δb/b_ℓ = 1% implies a 2% change in w_ℓ, in the sense that the high-ℓ power spectrum increases by ∼2%.

Zoom In Zoom Out Reset image size

Solid angle changes are necessarily attributable to changes in the symmetrized beam profiles. Selected beam profiles are compared in the left column of Figure 14. The right column partitions the solid angle difference between five years (Ω^5yr_B) and three years (Ω^3yr_B) into 1° radial bins. The contribution of each bin is plotted as a percentage of Ω^5yr_B. In this case, the beam profiles are extended to a radius of 10° using the far sidelobe patterns. Also, solid angles are normalized in such a way as to equalize b₂₀₀ between the resulting three-year and five-year beam transfer functions. Most of the beam profile change, and therefore most of the solid angle change, are just inside the three-year transition radius (dashed line). Compare with Figure 4, which shows, for the pure model case, how the increase in the fitted k_max of the mirror distortions for five years of data increases the solid angle inside the three-year radius.

Zoom In Zoom Out Reset image size

In summary, the solid angle increase appears to result primarily from the improved beam modeling, together with the extension of the main beam treatment to larger radii, both resulting from the increased S/N of the Jupiter data. Optimally hybridized two-dimensional beams are symmetrized and reduced to radial profiles in an unbiased way by averaging in annuli, and the resulting profiles are transformed directly to b_ℓ.

We adopt the analysis approach described in Page et al. (2003a) for reduction of the planet observations. The WMAP full-sky maps exclude observations made with Mars, Jupiter, Saturn, Uranus, or Neptune near any main beam boresight, with an exclusion radius of 15. In turn, the sky maps are used to remove the background sky signal from planet data. Since the solid angle of each planet is much less than that of the WMAP beams, a beam map is built up by binning observed antenna temperatures for a planet in a focal plane coordinate system. Rather than being normalized, this beam map may be left in antenna temperature and Legendre transformed (Equation (4)) to produce an unnormalized beam transfer function, T^m_JB_ℓ, where T^m_J is the peak antenna temperature of Jupiter. But B₀ is the beam solid angle Ω_B, so that T^m_JB₀ = T^m_JΩ_B (Page et al. 2003a), and T^m_JB₀/Ω^ref_J = T_J, where T_J is the brightness temperature of the Jupiter disk, and Ω^ref_J is the fiducial solid angle 2.481 × 10⁻⁸ sr for a Jupiter–WMAP distance of 5.2 AU (Griffin et al. 1986). The error in T_J is then the sum in quadrature of the error in Ω_B with the estimated WMAP gain calibration error of 0.2% (Hinshaw et al. 2009). The results of this procedure are given in Table 5 for five years of Jupiter data in each DA.

The main difference between the Legendre transform method and a direct integration of the two-dimensional beam map is that the Legendre transform uses a symmetrized beam profile. Integration of the beam map yields solid angles within the errors of the above approach. The results in Table 5 are consistent with the band averaged ones previously reported by Page et al. (2003a, Section 2.4). Currently, the error in T_J is limited by the 0.5% error in beam solid angle and the 0.2% gain uncertainty.

Season-by-season radiometry of Jupiter is given in Table 6. The values are computed using a template-fitting technique. Radial profiles are produced for each DA for each Jupiter season, then fitted linearly against the mean five-year Jupiter radial profile. Season 2 is omitted because Jupiter is approaching the Galactic anticenter, making background subtraction problematic. Our data place an upper limit on the time variability of T_J as a function of orbital phase of 0.3% ± 0.5%. We conclude that our radiometric observations are consistent with the absence of variability in the Jupiter brightness temperature at this level.

Table 6. Jupiter Temperature Changes by Season

			ΔT/T (%)
Seasona	Start	End	Meanb	Scatterc	$r/\overline{r}$ d
1	2001 Oct 8	2001 Nov 22	0.16	0.12	0.96
3	2002 Nov 10	2002 Dec 24	−0.14	0.23	0.97
4	2003 Mar 15	2003 Apr 29	−0.25	0.49	0.99
5	2003 Dec 11	2004 Jan 23	−0.02	0.24	1.00
6	2004 Apr 15	2004 May 30	−0.01	0.30	1.01
7	2005 Jan 9	2005 Feb 21	0.03	0.19	1.02
8	2005 May 16	2005 Jul 1	0.01	0.27	1.02
9	2006 Feb 7	2006 Mar 24	0.05	0.19	1.02
10	2006 Jun 16	2006 Aug 2	0.08	0.35	1.02

Notes. ^aSeason 2 omitted because Jupiter is in the Galactic plane. ^bMean of the percentage brightness temperature change among the DAs for each season, relative to the 5-year mean. ^c1σ scatter in the percentage temperature change among the DAs for each season. ^dMean Jupiter–WMAP distance for each season, relative to the 10-season mean = 5.34 AU.

Download table as:  ASCIITypeset image

In view of the stability and low errors of these measurements, Jupiter radiometry is the preferred method of transferring the WMAP dipole calibration to another microwave instrument. The key requirement for such an effort is knowledge of the beams. The error values given in Table 5 include WMAP beam errors via error in solid angle, as well as the fundamental gain uncertainty relative to the dipole.

Millimeter-wave brightnesses of other planets are also of potential interest as calibrators. For example, for the WMAP W band beams (3.2 mm), peak antenna temperatures of ∼ 200 mK, ∼35 mK, and ∼6 mK are produced by Jupiter, Saturn, and Mars, respectively. A preliminary analysis of the WMAP five-year Mars and Saturn observations has been undertaken.

Mars is attractive as a calibration source because it is relatively bright. However, significant variations in the observed brightness temperature can occur because of the viewing geometry. Moreover, the radiating properties of the inhomogeneous, pitted planetary surface complicate the determination of an appropriate reference brightness. Figure 15 shows a thermal model developed for the infrared by Wright (1976, 2007). The model is evaluated at 3.2 mm (W band) as a function of the time within the WMAP five-year timeline, which includes five Mars observing seasons (fewer than for Jupiter because of the relative orbital velocity of Mars). The predicted variation in brightness temperature over an observing season can be as much as ∼20 K. The mean 3.2 mm temperature and the scatter among the four W DAs are also shown for WMAP data binned by Mars observing season (Table 7). The model is higher than the WMAP measurements by ∼10%, so that a renormalization factor of 0.9 is applied to the model in the plot. We use this Mars model partly because of its previous use for Earth-based infrared calibration and the convenient availability of the code (Wright 2007); for a model including the effects of a dusty atmosphere and polar caps, see Simpson et al. (1981). The Mars data are referenced to a fiducial distance of 1.5 AU and a solid angle of Ω^ref_Mars = 7.156 × 10⁻¹⁰ sr (Hildebrand et al. 1985).

Zoom In Zoom Out Reset image size

Saturn's apparent brightness is even greater than that of Mars, but the theoretical understanding of the radiometry is less developed (Ulich 1981; Epstein et al. 1980; Hildebrand et al. 1985). A special consideration is Saturn's ring system, of which the viewing aspect from the Earth changes over the course of Saturn's 29 year orbital period. In Figure 16, mean seasonal W-band brightness temperatures as measured by WMAP (Table 7) are shown as black diamonds. These data show a clear decrease in observed temperature with time, a trend which correlates extremely well with the decreased viewing cross-section of Saturn's rings in the same time interval. A simple model of the form T_Sat = T₀ + αsin i, where i is the inclination of the ring plane from our line of sight, is fitted to the data and plotted in red. The fit results are α = −132 ± 16 and T₀ = 102 ± 7. Possible physical causes for the temperature decrease include the decreasing projected radiating area of the rings, a less favorable viewing angle for the "hot spot" at the south pole of Saturn, and Saturn's oblateness. These causes will be the subject of future investigation. The Saturn data are referenced to a fiducial distance of 9.5 AU, corresponding to a Saturn solid angle of Ω^ref_Sat = 5.101 × 10⁻⁹ sr (Hildebrand et al. 1985).

Zoom In Zoom Out Reset image size

Clearly, Jupiter remains the only WMAP source that can be recommended as an instrument calibrator at the 1% level. However, our preliminary results for Mars and Saturn suggest that, with additional analysis and observations, both of these sources may be similarly useful in the future.