THE IMPORTANCE OF WIDE-FIELD FOREGROUND REMOVAL FOR 21 cm COSMOLOGY: A DEMONSTRATION WITH EARLY MWA EPOCH OF REIONIZATION OBSERVATIONS

A major goal of modern experimental cosmology is the detection of 21 cm emission from neutral hydrogen at high redshifts. Depending on the redshifts studied, these observations can probe a wide range of physical and astrophysical phenomena. Observations at ∼100–200 MHz ($z\sim 6\mbox{--}13$ in the 21 cm line) probe the Epoch of Reionization (EoR)—the reionization of the intergalactic medium (IGM) by ultraviolet photons emitted by the first stars and galaxies. Observations at higher frequencies (lower redshifts) trace the neutral hydrogen that remains in galactic halos, and provide a low-resolution "intensity map" of large-scale structure and, potentially, the baryon acoustic oscillation (BAO) features in the power spectrum. At lower frequencies (higher redshifts), one begins to trace the birth of the first stars during "Cosmic Dawn" and even the preceding Dark Ages. For reviews of the 21 cm cosmology technique and the associated science drivers, see Furlanetto et al. (2006), Morales & Wyithe (2010), Pritchard & Loeb (2012), and Zaroubi (2013).

A large number of experiments seeking to detect the power spectra of 21 cm fluctuations are already operational or being commissioned, including the LOw Frequency ARray (LOFAR; Yatawatta et al. 2013; Wise et al. 2013)²⁶ , 21 CentiMeter Array (21CMA; Zheng et al. 2012)²⁷ , the Giant Metrewave Radio Telescope EoR Experiment (GMRT; Paciga et al. 2013)²⁸ , the MIT Epoch of Reionization Experiment (MITEoR; Zheng et al. 2014), the Donald C. Backer Precision Array for Probing the Epoch of Reionization (PAPER; Parsons et al. 2010)²⁹ , and the Murchison Widefield Array (MWA; Lonsdale et al. 2009; Bowman et al. 2013; Tingay et al. 2013)³⁰ , all of which are targeting the signal from the EoR. A number of additional experiments are also under construction or planned, such as the low-frequency Square Kilometre Array (SKA-low; Mellema et al. 2013)³¹ and the Hydrogen Epoch of Reionization Array (HERA; Pober et al. 2014)³² at EoR and Cosmic Dawn redshifts, and BAOs from Integrated Neutral Gas Observations (BINGO; Battye et al. 2013), TianLai³³ , BAORadio (Ansari et al. 2012a, 2012b), the Canadian Hydrogen Intensity Mapping Experiment (CHIME; Shaw et al. 2014)³⁴ , and the BAO Broadband and Broad-beam experiment (BAOBAB; Pober et al. 2013b) at lower redshifts.

At all redshifts, however, 21 cm experiments are limited by both the inherent faintness of the cosmological signal and the presence of foregrounds, which can exceed the 21 cm emission by as much as five orders of magnitude in brightness temperature (Santos et al. 2005; Bernardi et al. 2013; Pober et al. 2013a; Yatawatta et al. 2013). As such, the only current detection of H i at cosmological distances comes from cross-correlation studies using maps from the Green Bank Telescope and optical galaxy surveys (Chang et al. 2010; Masui et al. 2013; Switzer et al. 2013). Analysis techniques for recovering the signal focus on the relative spectral smoothness of the foreground emission as an axis for distinguishing these contaminants from the 21 cm emission. Over the past decade, a large body of literature has worked to develop pipelines that can subtract foreground sources from 21 cm data sets (e.g., Morales et al. 2006; Bowman et al. 2009; Liu et al. 2009; Liu & Tegmark 2011; Chapman et al. 2012, 2013; Dillon et al. 2013; Wang et al. 2013). More recently, however, studies of the chromatic interaction of an interferometer with foreground emission have demonstrated that smooth-spectrum foregrounds will occupy an anisotropic wedge-like region of cylindrical (${k}_{\perp },{k}_{\parallel }$) Fourier space, leaving an "EoR window" above the wedge where the 21 cm signal can be cleanly observed (Datta et al. 2010; Morales et al. 2012; Parsons et al. 2012b; Trott et al. 2012; Vedantham et al. 2012; Thyagarajan et al. 2013; Liu et al. 2014a, 2014b). These predictions have since been confirmed in data sets from PAPER and the MWA (Pober et al. 2013a; Dillon et al. 2014; Parsons et al. 2014; Ali et al. 2015; Jacobs et al. 2015; Thyagarajan et al. 2015a), although significantly more sensitive observations will be necessary to see if the window remains uncontaminated down to the level of the 21 cm signal. Pober et al. (2014) demonstrate that while current EoR observatories (PAPER, the MWA, and LOFAR) do not possess the sensitivity to detect the 21 cm signal with this pure "foreground avoidance" technique, next-generation experiments such as HERA and the SKA-low can yield high-fidelity power spectrum measurements using this approach, and begin to place constraints on the physics of reionization.³⁵ However, the cosmological signal strength peaks on large scales, so that k modes within the wedge can have significantly more 21 cm power than modes within the window. Pober et al. (2014) show that if foregrounds can be subtracted from 21 cm data sets, allowing the recovery of k modes from within the wedge, then the significance of any power spectrum measurement can be substantially boosted—enabling the current generation of 21 cm experiments to make a detection.

Continued research into foreground subtraction algorithms is therefore clearly well motivated. As yet, no technique—whether subtracting a model of the sky or using a parameterized fit in frequency—has demonstrated that foreground emission in actual observations can be removed to the thermal noise level of current instruments (although the EoR window has, to date, proven relatively free of foregrounds when care is taken to limit leakage from the wedge; Pober et al. 2013b; Parsons et al. 2014; Ali et al. 2015; Jacobs et al. 2015). The purpose of this work is to investigate some of the wide-field effects that complicate the removal of foreground emission using data from the MWA. In particular, we focus on the contribution of sources outside the main lobe of the instrument's primary beam (in this work, we use the term primary beam to refer to the all-sky power pattern of the antenna or tile element, including sidelobes). Far from the pointing center, chromatic effects in the interferometer response become stronger; sources out in the sidelobes of the primary beam therefore create foreground contamination in higher k_∥ modes than sources near the pointing center. Here, we explore this effect in more detail.

This paper is structured as follows. In Section 2, we lay out a heuristic derivation of how the instrument's primary beam enters in measurements of the 21 cm power spectrum and how foregrounds are distributed throughout the $({k}_{\perp },{k}_{\parallel })$ plane. In Section 3, we briefly describe the MWA and the data analyzed in this study. In Section 4, we build on the pedagogical nature of the previous analysis through simulated MWA power spectra using a sky model containing a single point source of emission. By changing the location of this source, we demonstrate these primary-beam effects in a realistic but controlled fashion. In Section 5, we describe the calibration, preprocessing, and foreground subtraction applied to the observed data before making a power spectrum. The main result is presented in Section 6, where we compare two distinct power spectra made from our data: one where we have subtracted a foreground model that includes sources in the beam sidelobes, and one where only sources in the main field of view are removed. We discuss the implications of these results for future foreground subtraction efforts in Section 7 and conclude in Section 8.

Although many 21 cm experiments have wide fields of view, only recently have studies focused on how wide-field effects might complicate measurements of the 21 cm power spectrum. Theoretical work has identified the foreground wedge described above and provided a formalism for mapping the position of foreground emission on the sky to k modes of the 21 cm power spectrum (Morales et al. 2012; Parsons et al. 2012b; Trott et al. 2012; Vedantham et al. 2012; Thyagarajan et al. 2013; Liu et al. 2014a, 2014b). Broadly speaking, there are two flavors of 21 cm power spectrum analysis: a "delay spectrum" approach, where the line-of-sight Fourier transform is done on individual visibilities—and is therefore not strictly orthogonal to the transverse directions because of the frequency dependence of an individual visibility—and an "imaging" approach, where the line-of-sight Fourier transform spans multiple visibilities and is truly orthogonal to the transverse directions on the sky.³⁶ A full discussion of the differences between these two approaches is outside the scope of this work, but previous analyses have shown that the wedge and the mapping from foreground sky position to k_∥ modes of the power spectrum remain valid for both frameworks. In the delay spectrum approach, the chromatic dependence of an individual baseline is completely preserved, so that all foreground emission at a given location maps to a given k_∥ mode. An imaging approach, however, removes the mapping between delay and sky position by projecting out the frequency sine wave for a known geometric delay. In an imaging power spectrum, frequency structure is dominated by the intrinsic spectra of the sources, so that a significant amount of foreground emission maps to low k_∥ modes, reflective of their inherent (smooth) frequency spectrum. However, the chromatic response of the interferometer still affects the observed emission, leading to a wedge feature analogous to that of the delay spectrum approach, but with more of the emission concentrated at low k_∥ (Morales et al. 2012; Dillon et al. 2015). Because of the brightness of foreground emission, this wedge still dominates any 21 cm signal in the modes it occupies.

Explorations of these wide-field effects in actual data have been more limited. Thyagarajan et al. (2015a, 2015b) studied both simulated and actual MWA observations using the delay spectrum technique and found an excellent match between the two, demonstrating a good understanding of both foreground emission and the primary beam of the MWA. They also found that the foreshortening of baseline lengths when projected toward the horizon creates sensitivity to diffuse emission normally resolved out on longer baselines. Diffuse foregrounds are bright enough that they can be detected despite the small (but non-zero) response of the MWA element toward the horizon. This led to what they dubbed the "pitchfork" effect, a foreground signature in delay space where bright emission from within the main field of view appeared at low delays, and emission from the horizon at high delays.

This work studies similar effects using an imaging power spectrum approach and will confirm that the sky-position to k_∥ mapping still holds. We will also focus on the ability to subtract foreground emission away from the main field of view to lower the contamination in high k_∥ modes. In this section, however, we use the delay-spectrum formalism (Parsons & Backer 2009; Parsons et al. 2012b) to provide a general framework for understanding these effects. We stress that the delay spectrum provides a straightforward, pedagogical way to interpret power spectrum results, since the wide-field chromatic effects appear at first order. As argued in Morales et al. (2012), Trott et al. (2012), and Liu et al. (2014a, 2014b), and as will be confirmed with data below, these wide-field effects are generic to all interferometric 21 cm experiments.

The basic premise of the delay-spectrum technique presented in Parsons et al. (2012b) is that the square of the frequency Fourier transform of a single baseline's visibility spectrum (i.e., the delay spectrum, ${\tilde{V}}_{b}(\tau )$) approximates a measurement of the cosmological power spectrum (to within a proportionality factor):

$$\begin{eqnarray}&&| {\tilde{V}}_{b}(\tau ){| }^{2}\propto P({k}_{\perp },{k}_{\parallel }),\end{eqnarray} \tag{ 1 }$$

where

$$\begin{eqnarray}&&{\tilde{V}}_{b}(\tau )=\int d\nu \;{V}_{b}(\nu )\;{e}^{2\pi i\nu \tau }\end{eqnarray} \tag{ 2 }$$

is the delay spectrum, τ is delay, ν is frequency, V is a visibility, and the subscript b indicates that the visibilities are from a single baseline.

Intuitively, this relation is well motivated. To a good approximation, a single baseline b probes a single transverse scale, and thus a single k_⊥ mode. And, since cosmological redshifting of the 21 cm line maps observed frequencies into line-of-sight distances, the Fourier transform of the frequency spectrum approximates a range of k_∥ modes. Put more succinctly, for an interferometer, baseline length b maps to cosmological k_⊥ and delay τ maps to k_∥.

The power of this simple formalism is that we can now map the effects of the primary beam, which enter into a visibility measurement in a well-known way, to cosmological Fourier space and the power spectrum $P({k}_{\perp },{k}_{\parallel })$. We begin with the form of a visibility in the flat-sky approximation (Thompson et al. 2001)³⁷ :

$$\begin{eqnarray}&&{V}_{b}(\nu )=\int \;{dl}\;{dm}\;A(l,m,\nu )\;I(l,m,\nu )\;{e}^{-2\pi i({ul}+{vm})},\end{eqnarray} \tag{ 3 }$$

where A is the primary beam, I is the sky brightness distribution, l and m are direction cosines on the sky, ν is frequency, and u and v are the projected baseline lengths on the ground plane measured in wavelengths. We can rewrite this expression in terms of the geometric delay ${\tau }_{g}$ (Parsons & Backer 2009):

$$\begin{eqnarray}&&{V}_{b}(\nu )=\int \;{dl}\;{dm}\;A(l,m,\nu )\;I(l,m,\nu )\;{e}^{-2\pi i\nu {\tau }_{g}},\end{eqnarray} \tag{ 4 }$$

where

$$\begin{eqnarray}&&{\tau }_{g}=\displaystyle \frac{{\boldsymbol{b}}\cdot \hat{s}}{c}=\displaystyle \frac{1}{c}({b}_{x}l+{b}_{y}m),\end{eqnarray} \tag{ 5 }$$

${\boldsymbol{b}}\equiv ({b}_{x},{b}_{y})$ is the baseline vector measured in meters (i.e., ${\boldsymbol{u}}\equiv (u,v)=\nu {\boldsymbol{b}}/c)$, and $\hat{s}\equiv (l,m)$. Performing the delay transform given by Equation (2) gives us a delay spectrum:

$$\begin{eqnarray}&&{\tilde{V}}_{b}(\tau )=\int \;{dl}\;{dm}\;d\nu \;A(l,m,\nu )\;I(l,m,\nu )\;{e}^{-2\pi i\nu ({\tau }_{g}-\tau )}.\end{eqnarray} \tag{ 6 }$$

If we make the pedagogical assumption that both A and I are independent of frequency, we can straightforwardly perform the delay transform integral³⁸ :

$$\begin{eqnarray}&&{\tilde{V}}_{b}(\tau )=\int \;{dl}\;{dm}\;A(l,m)\;I(l,m)\;\delta ({\tau }_{g}-\tau ).\end{eqnarray} \tag{ 7 }$$

Since Equation (5) relates the geometric delay ${\tau }_{g}$ to a specific set of sky direction cosines (l, m), the delta function selects a subset of sky positions that contribute to each τ mode in the delay spectrum, albeit with a baseline-dependent non-trivial mapping between sky position and τ. It is always true, however, that sources that appear at high delays are those that are far from the pointing center of the instrument (hence the name "horizon limit" given to the maximum delay a source can appear at in Parsons et al. 2012b). Following Equation (1), we can say

$$\begin{eqnarray}&&P({k}_{\perp },{k}_{\parallel })\propto {\left[\int {dl}{dm}A(l,m)I(l,m)\delta ({\tau }_{g}-\tau )\right]}^{2},\end{eqnarray} \tag{ 8 }$$

where the length of baseline b sets k_⊥ and $\tau \propto {k}_{\parallel }$. This analysis therefore implies that sources at large delays (i.e., sources near the edges of the field of view, by Equation (5)) contaminate the highest k_∥ modes of the wedge $({k}_{\perp },{k}_{\parallel })$ space. Although not always stated as directly, this result was also found in Morales et al. (2012), Vedantham et al. (2012), Thyagarajan et al. (2013), and Liu et al. (2014a) using entirely independent formalisms.

Equation (8) also shows the main result we wished to derive in this section: the (smooth-spectrum) sky emission $I(l,m)$ that appears in each delay mode is multiplicatively attenuated by the primary beam of the instrument. Therefore, the foreground emission that contaminates those k_∥ modes measured by a single baseline will itself be attenuated by a (distorted) slice through the square of the primary beam of the instrument. This result is illustrated schematically in Figure 1. Note that, because delay space is a one-dimensional projection of the sky coordinates (see Equation (5)), the attenuating beam in a k_∥ spectrum will vary depending on the orientation of the baseline. On an east/west baseline, for example, the delay axis probes the relative east/west position of the source and is insensitive to north/south translations in source positions. Such a baseline will therefore clearly show the effects of the eastern and western sidelobes of the primary beam in its k_∥ spectrum. Similar logic applies to a north/south baseline and the northern and southern sidelobes of the primary beam. Delays on a northeast/southwest baseline, however, probe northeast/southwest sky position, and thus the east/west translation of a source through the eastern and western sidelobes does not cause as rapid a change in k_∥. The net effect is that when all baselines of the same magnitude are averaged into a k_⊥ bin, these different k_∥ sidelobe patterns add up and smear out the location of the sidelobes.

Zoom In Zoom Out Reset image size

An important but subtle point is that the above derivation for mapping sky coordinates into k space was strictly for flat-spectrum emission. As shown in Parsons et al. (2012b), any spectral structure—whether intrinsic to the source or the instrumental response—introduces a convolving kernel that broadens the footprint of each k_∥ mode in cosmological Fourier space. While this kernel is narrow for smooth-spectrum foregrounds, spectral structure in the 21 cm signal spreads the 21 cm power across a wide range of k_∥ modes. This is equivalent to saying that the 21 cm signal intrinsically has power on these cosmological scales. The situation for foreground emission is different, however. Although power spectrum plots are labeled with axes of $({k}_{\perp },{k}_{\parallel })$ with units of $h{{\rm{Mpc}}}^{-1}$, these cosmological scalings apply only to the 21 cm signal. The analysis above shows how foregrounds map into this space, and how the primary beam affects this mapping. The primary beam of the instrument does still act as a window function and can affect high k_∥ modes of the cosmological signal; however, the cosmological signal has been shown to be relatively featureless on the scale of this kernel (see Parsons et al. 2012b), rendering this effect very small. Nevertheless, the 21 cm signal is an all-sky signal with real intrinsic k_∥ structure. There is therefore always a 21 cm signal at the peak beam response, so there will always be power at all k_∥ modes truly intrinsic to the cosmological signal. This point will be discussed further in Section 7, where we consider the possibility of detecting 21 cm emission at k_∥ modes where the foregrounds fall in the nulls of the primary beam.

The very wide and relatively smooth primary beam of the PAPER instrument makes the predicted foreground attenuation difficult to see in the analysis of Pober et al. (2013a). However, for instruments such as the MWA and LOFAR, which use tiles of dipoles to increase the system gain and narrow the size of the primary beam, there should be two clear effects visible in the power spectra. First, there should be significant attenuation of the wedge foreground emission before the horizon limit, since the instrument's field of view is significantly smaller than $2\pi$ steradians, as is seen in Dillon et al. (2014). Second, at higher k_∥ values than those corresponding to the main beam of the instrument, foreground emission should appear coming from the sidelobes of the primary beam. These two effects can be seen in the delay-space simulations of different antenna elements presented in Thyagarajan et al. (2015b). For an imaging power spectrum technique that averages baselines together, the second effect will be less clear for an instrument such as the MWA, in which all the dipoles and tiles are oriented in the same direction. In this case, the sidelobes are always oriented north/south and east/west; as explained above, however, the beam footprint in k_∥ will differ from baseline to baseline, depending on each baseline's orientation relative to the sidelobe pattern. This will have the effect of smearing out the sidelobe across a wider range of k_∥ modes than would be seen in an instrument with circularly symmetric sidelobes, but as we will show, the feature is still quite visible in the power spectrum.

The structure of the remainder of this paper is as follows. First, in Section 3, we describe the MWA instrument and observations in more detail. With this context provided, we provide the results of two principal analyses. In Section 4, we present simulated MWA power spectra made from a sky consisting of a single point source. By moving the position of this source from simulation to simulation, we can see the primary-beam effects described above in a controlled fashion. In Section 5, we use observations from the MWA to analyze these primary-beam features and present the power spectra of these data in Section 6. In particular, we focus on the effect of subtracting sources from sidelobes outside the primary field of view.

The MWA in Australia consists of 128 tile elements, and each tile is composed of 16 dual-polarization dipole antennas; the array configuration is shown in Figure 2. The tile element has the effect of significantly narrowing the MWA's field of view over that from a single dipole, but also introduces significant regular sidelobes in the primary beam. Figure 3 shows three MWA tiles; every tile is aligned north/south, so the sidelobes from each tile appear with nearly the same orientation.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The data used in this work were taken with the MWA on 2013 August 23 (Julian Date 2456528) over the course of approximately three hours from 16:47:27 to 19:58:24 UTC. The observations were taken over a frequency band centered on 182.415 MHz, with a total bandwidth of 30.72 MHz divided into 24 1.28 MHz coarse channels, which are each further divided into 768 40 kHz fine channels.

The data used in this analysis span a total of six pointings each 30 minutes long, where an analog beamformer steers the main lobe of the primary beam to nearly the same sky coordinates for each pointing. The sky is then allowed to drift overhead for 30 minutes before repointing. The data within each pointing are saved as individual "snapshot" observations, each lasting 112 s, with individual integrations of 0.5 s. Figure 4 shows the tile primary beam at three different beamformer pointings: the beginning of the observation, a zenith-phased pointing, and the end of the observation. Since each pointing changes the overall primary-beam response of the instrument, the sidelobe patterns in the final integrated power spectrum will be smeared. As will be shown below, however, the effects of the sidelobes are still quite visible despite the changing shape of the primary beam.

Zoom In Zoom Out Reset image size

Before presenting the full analysis of this data set, we will first investigate the effects of the location of celestial emission on the cosmological power spectrum and the wedge in particular. In this section, we will simulate visibilities for a single point source and calculate the dependence of the power spectrum on the source's location. Visibilities are simulated using the Fast Holographic Deconvolution (FHD) software package.³⁹ Visibility simulation is one of several functions in FHD; as described below, FHD also performs calibration and source subtraction on our actual data. As a simulator, FHD constructs a model of the sky in uv space and integrates small regions of the uv plane using the holographic beam kernel (Morales & Matejek 2009) to create model visibilities.

For this analysis, we simulate visibilities for all the baselines in the MWA in 768 fine frequency channels spanning the observed 30.72 MHz frequency band. We simulate only one 112 s snapshot when the primary beam is pointed at zenith (i.e., the snapshot shown in the middle panel of Figure 4). In addition to reducing the computational demand of the simulations, using only one snapshot allows us to see the sidelobe patterns most clearly, since integrating over a longer time means including data when the array had a different pointing and primary beam.

We conduct four simulations, each consisting of one radio point source at a different location on the sky; the locations simulated are shown in Figure 5. For each simulation, the inherent flux density of the source is increased relative to source D (located at zenith) by the inverse of the primary beam response at its location. In other words, each source simulated has the same apparent flux density. This choice places all the final power spectra on the same scale, allowing for easier comparison.

Zoom In Zoom Out Reset image size

In Figure 6, we show the 2D $({k}_{\perp },{k}_{\parallel }$) power spectra for each of the four simulations described above. We show the power spectrum from only one of the two linearly polarized dipoles of the MWA; the power spectrum for the other polarization is quite similar. Letters correspond to the source labels in Figure 5. To make the power spectrum, the simulated visibilities are imaged by FHD and then analyzed by the pipeline described in B. J. Hazelton et al. (2015, in preparation).⁴⁰ For more information on the data products transferred between FHD and , see D. C. Jacobs et al. (2015, in preparation).

Zoom In Zoom Out Reset image size

The effect of source position on the power spectrum is clear and agrees with the intuition developed in Section 2. Source D is located directly at zenith, with the subsequent sources offset to higher declination (with right ascension held fixed). In Figure 6, source D exhibits no wedge feature. (The power at high k_⊥ values is due to poor uv coverage on these scales and is described in more detail below.) Sources C and B show a clear wedge feature arise as the source is moved away from zenith, and the power spectrum of source A—where the source is located in the sidelobe of the primary beam—shows a concentration of power outside the main field of view (indicated by the dashed black line) but inside the horizon limit (solid black line). This feature is in exact accord with our predictions. Simulations using sources offset in right ascension (instead of declination) show the same effect, as do sources with offsets in both right ascension and declination: power moves to higher k_∥ as the source moves further from the field center.

In this section we present the full analysis of the three hours of MWA data described in Section 3. The data are processed through the same imaging and power spectrum analyses (done by FHD and , respectively) applied to the simulations. However, there are initial preprocessing, calibration, and foreground subtraction steps applied to the data, which we describe here.

5.1. Preprocessing

5.2. Calibration and Imaging

5.3. Foreground Subtraction

It is through FHD that model visibilities are also subtracted from the data. We use two sets of model visibilities generated from a custom-made point source catalog. In the main field of view, the catalog contains sources generated from FHD deconvolution outputs and an advanced machine-learning source identifier designed to reject spurious sources (P. A. Carroll et al. 2015, in preparation). Outside the main field of view, our catalog combines sources from the MWA Commissioning Survey (MWACS; Hurley-Walker et al. 2014), the Culgoora catalog (Slee 1995), and the Molongolo Reference Catalog (Large et al. 1981). In one model we include only ∼4600 sources that fall within the primary lobe of the MWA beam; in the other, we include all sources up to and including the first sidelobe (∼8500 sources). An image of all the sources included during the zenith-phased snapshot is shown in Figure 7. There are two effects that serve to limit the number of sources included in our model. First, we use a primary-beam threshold cut: any sources that fall where the beam response is less than 1% of the peak response are not included in the model. Second, because it is a composite of several surveys, the completeness of our catalog is not uniform over the sky. In particular, MWACS does not cover the full declination range of the observations here; the effect is that fewer sources are removed from the lower declinations of the southern sidelobe, and very few are in the northern sidelobe. This has the effect of introducing a small time dependence in the number of sources included in our model, since the declination coverage of the primary beam does change with pointing (see Figure 4). MWACS also avoids the Galactic plane, which reduces the number of sources in the model at the early and late pointings to ∼7000.

Zoom In Zoom Out Reset image size

It is also important to note that our sky model assumes a fixed spectral index of −0.8 for each source. Although the actual sources on the sky will have some spectral structure, the fact that we include minimal frequency dependence in the model serves to strengthen the arguments below: subtracting a nearly achromatic foreground model removes power from chromatic (i.e., high k_∥) modes of the power spectrum. This is a clear demonstration of the inherent chromaticity of the interferometer response pattern.

We now present the power spectra of these data generated by the code. With observational data, empirically calculates the noise level in the visibilities and fully propagates errors in the visibilities through to the 3D power spectrum. The important results here are the cylindrically averaged 2D power spectra, shown in Figure 8. In each row in this figure, the left-hand panel shows the power spectrum with only sources in the primary lobe removed, while the center panel shows the power spectrum where sources are also subtracted from the sidelobes. In order to enhance the subtle difference between the two panels, we subtract the power spectrum where sidelobe sources were removed from the power spectrum where only main-lobe sources were removed (i.e., we subtract the center panel from the left-hand panel). Note that we perform the subtraction in full 3D $({k}_{x},{k}_{y},{k}_{z})$ space before binning into 2D $({k}_{\perp },{k}_{\parallel })$ space. We plot the result of this subtraction in the right-hand panel in each row of Figure 8. Most of the difference fluctuates randomly between positive (blue) and negative (red) values, showing no systematic change of the power spectrum in these regions. However, the consistently blue region shows that subtracting sources from the sidelobes removes a non-trivial amount of power (as much as 10% compared to the power spectrum with no subtraction of sidelobe sources, although typical values are $\sim 1\%$) from the region where the sidelobe is expected—outside the main lobe (dashed black line) but within the horizon (solid line). Since the size of the main lobe is dependent on frequency and pointing, the dashed black line is only an approximate marker; the power that is removed from k_∥ modes below this line is consistent with being sidelobe power from a range of frequencies and pointing centers.

Zoom In Zoom Out Reset image size

Although not the primary goal of this paper, there are a few additional features in the power spectra that warrant explanation.

• 1.  
  The horizontal lines running across the EoR window are the effect of the coarse channelization used by the MWA. Between any two 1.28 MHz coarse channels are two 80 kHz channels that are flagged due to low signal response and potential aliasing concerns. This flagging in frequency has the effect of introducing covariance into the line-of-sight k_∥ modes, which are effectively produced by a Fourier transform of the frequency axis. This additional covariance has the effect of coupling power from the wedge into higher k_∥ modes. Because the flagging is at regular intervals, this additional power also appears at regular intervals in k_∥ (the appearance of irregular spacing comes from the logarithmic scale on the y axis). Work is underway on algorithms that can reduce this covariance using priors on the fact that the power comes from k_∥ modes within the wedge.
• 2.  
  The vertical lines, which are especially prevalent at high k_⊥ modes, come from the uv coverage of the MWA. The MWA has exceptionally dense coverage at low k_⊥ due to its large number of short baselines. However, at higher ${k}_{\perp }\;(\gtrsim {10}^{-1}\;h{{\rm{Mpc}}}^{-1})$ there are gaps in the coverage, which results in particularly noisy measurements of certain modes. Therefore, while these modes appear to have very high power, they also have very large associated error bars. A plot of the errors calculated by for the XX polarization is shown in Figure 9.
• 3.  
  There are blue/purple regions outside the wedge that are negative. This is because cross-multiplies the even time samples in the data set with the odd time samples (with the samples interleaved on a timescale of 2 s); this has the effect of removing the positive-definite noise bias that would result from squaring the entire data set. Alternating positive and negative values correspond to noise-dominated regions.
• 4.  
  Most obviously, a large amount of foreground power remains in the power spectra. This is not surprising, because our analysis subtracted only a few thousand point sources, ignoring diffuse emission from both the Galaxy and unresolved point sources. Subtracting models of this emission will clearly be necessary for any possibility of recovering the 21 cm signal from inside the wedge. The effects concerning sidelobes presented in this work, however, are still quite important: the additional fraction of emission removed when including sidelobe sources is more than enough to swamp the EoR signal, which might have a peak power spectrum brightness of the order of ${10}^{6}\;{{\rm{mK}}}^{2}\;{h}^{-3}{{\rm{Mpc}}}^{3}$.

Zoom In Zoom Out Reset image size

Through the advances in our understanding of EoR foregrounds (i.e., the "wedge" and "EoR window" paradigm), we now have a model for the detailed impact of sources far from pointing center on the recovery of the 21 cm power spectrum. This work demonstrates that sources outside the main field of view are a significant contaminant of the modes of interest in the 21 cm power spectrum, even for an "imaging" power spectrum analysis. It is therefore worthwhile to heuristically consider the detailed pattern that sources far from pointing center leave in cylindrically averaged (${k}_{\perp },{k}_{\parallel }$) space. While not all of the conclusions below directly follow from the empirical power spectra analyzed here, the formalism presented in Section 2, the delay spectrum analyses in Thyagarajan et al. (2015a, 2015b), combined with the results of our sidelobe source subtraction from MWA observations suggest several interesting lines of reasoning. To guide this discussion, we divide the sky into three rough categories—the primary field of view, the sidelobes, and the nulls—and discuss the effects of sources appearing in each regime.

• 1.  
  Primary Field of View: These are the sources that are traditionally considered when treating foreground removal from 21 cm experiments. Although the primary field of view may not be at zenith for phased array telescopes (e.g., MWA and LOFAR), phase rotation still places these sources at low delays, and therefore low k_∥ in power spectrum measurements. The dashed line in Figure 8 roughly indicates the edge of the main lobe of the MWA, and the brightness of emission can be seen clearly to fall as one moves to higher k_∥ modes. Since this emission is located where the beam response is at a maximum, it appears as the brightest contaminant in the 21 cm power spectrum. Because they are detected at high signal-to-noise ratio, point sources in the primary field of view are often used to simultaneously calibrate the instrument response while they are subtracted from a visibility model (e.g., Yatawatta et al. 2013). Diffuse emission and unresolved point sources generally dominate the total foreground power, requiring additional models or parametric methods for removal (e.g., Chapman et al. 2012).
• 2.  
  Sidelobes: As seen in the present work, emission in the sidelobes appears at higher k_∥ values than emission from inside the primary field of view. Therefore, a model of the sidelobes and the emission that falls within them must be subtracted from the data in order to recover these k_∥ modes closer to the EoR window. The primary beam attenuation has the effect of reducing the fractional accuracy required in modeling these sources, since only their residual apparent flux density contaminates the power spectrum. Since emission in sidelobes contaminates higher k_∥ modes than emission in the main beam, however, removing wide-field emission may be more effective at reducing leakage into the EoR window than removing emission from the primary field of view.
• 3.  
  Nulls: Lastly, one might assume that the nulls in the primary beam might serve as good "EoR windows", just like the regions of k space outside the horizon. (Recall the discussion in Section 2: while nulls on the sky clearly attenuate all emission from that sky position, when we refer to beam nulls in k_∥, however, these nulls are confined to the foreground emission. The EoR signal from other positions on the sky still produces unattenuated power at these k_∥ modes.) And while the present analysis does indeed suggest that areas of k space corresponding to sky positions of exceedingly low beam response will be free from foreground contamination, some caveats must be issued. First, the nulls between sidelobes are likely not as deep as analytic models suggest, due to effects such as mutual coupling between the dipoles in a station and group delay errors (Neben et al. 2015). Second, as derived in Parsons et al. (2012b), there is a non-negligible k-space point-spread function convolving each source of emission. Therefore, while there may be narrow nulls between sidelobes, emission within the sidelobes can contaminate these nulls due to this spillover effect. Finally, it is worth remembering that the mapping from zenith angle to delay (which, recall, maps to k_∥) is nonlinear. A delay bin near the horizon corresponds to much lower elevation than a delay bin near zenith. This means that while emission in bins far from the main lobe of the beam is strongly attenuated by the beam response, the total aggregate sum of emission in that k_∥ bin can still be large, since it corresponds to a large area of sky (Thyagarajan et al. 2015b). It may still be that even with instruments such as SKA and HERA, which have narrower fields of view and less sensitivity to foreground emission away from pointing center (Thyagarajan et al. 2015b estimate SKA and HERA will suppress emission near the horizon 40 dB more than MWA), the only safe place for foreground "avoidance" is beyond the horizon.

These arguments have important ramifications for experiments looking to subtract foreground emission and recover k modes from inside the wedge. In particular, they suggest that experiments looking to "enlarge" the EoR window and remove foreground emission from modes near the horizon will benefit most from subtracting emission outside the main lobe. For this subtraction to be effective, accurate wide-field primary beam calibration is necessary to properly characterize the sidelobe patterns as a function of frequency. Such wide-field calibration may require new techniques, e.g., Pober et al. (2012), Yatawatta et al. (2013), and Neben et al. (2015). Additionally, these arguments motivate the need for low-frequency, wide-field sky surveys, especially in the Southern Hemisphere, where the vast majority of EoR-frequency 21 cm experiments are being constructed. Experiments such as HERA do not have a steerable beam, and thus will have difficulty measuring the flux densities of sources in their sidelobes. An accurate catalog produced by another survey covering a larger area (e.g., Jacobs et al. 2011; Williams et al. 2012; Hurley-Walker et al. 2014; Wayth et al. 2015) will be highly valuable for subtracting sources outside the main field of view. Northern Hemisphere experiments such as LOFAR and GMRT may also be valuable for characterizing the foregrounds at higher declinations.

Another important conclusion of this work is the implication that foreground subtraction cannot simply target the removal of some total amount of flux density independently of the position of that emission on the sky. Even if all emission from inside the primary field of view could be perfectly removed, sources in the sidelobes would continue to dominate higher k_∥ modes. To recover all modes of the 21 cm power spectrum, foreground models must extend into any primary-beam sidelobes where the level of beam attenuation does not reduce the foreground power below that of the EoR signal. While attention has been paid to the removal of bright off-axis point sources (Offringa et al. 2012), the remaining diffuse emission and confused source background will still have significant spectral structure from the instrumental effects we have described. Given the extremely high foreground-to-signal ratio in 21 cm experiments, emission far from the pointing center cannot be neglected even if it is largely attenuated by the instrument's primary beam.

In practice, wide-field source subtraction at the level needed to recover the EoR signal may require more than an accurate foreground model, especially for experiments with very wide fields of view such as PAPER. First, curved-sky effects become important near the horizon; an imaging-based analysis that does not correctly handle the curved sky (e.g., with w-projection; Cornwell et al. 2008) could wash out the input source model and reduce the amount of power subtracted. Second, ionospheric effects may become important for the level of accuracy needed in subtraction (Mitchell et al. 2008; Bernardi et al. 2009; Intema et al. 2009). The ionosphere may introduce frequency-dependent, time-dependent, and direction-dependent gains, all of which could lead to errors in model-based source subtraction if not corrected for. None of these effects alleviates the need for wide-field foreground subtraction, however; rather, they increase the difficulty of implementing a scheme that potentially removes foregrounds to a level below the EoR signal.

In this work, we have presented a heuristic description of the imprint of the primary beam in power spectrum measurements from 21 cm interferometers. In particular, we find that wide-field effects—especially primary-beam sidelobes—leave a highly chromatic imprint, so that even smooth-spectrum emission that falls within the sidelobes corrupts high k_∥ modes of the power spectrum. We further demonstrate this effect both with pedagogical simulations using single point sources and by removing a source model from MWA observations. When the model includes sources out to the first primary-beam sidelobes, it produces a significant (percent-level) reduction in power at high k_∥ values.

This result has significant implications for experiments looking to measure the power spectrum of 21 cm emission from the EoR or any other epoch. In particular, it shows that foregrounds must be considered as a wide-field contaminant. Removing foregrounds from just the primary field of view will not reduce power in high k_∥ modes corresponding to the primary-beam sidelobes. As a corollary, pipeline simulations that include only foregrounds within the primary field of view are missing a major contaminant of the EoR signal. Experiments looking to use foreground subtraction to enlarge the EoR window must also pay particular attention to emission in the sidelobes, since it is this emission that corrupts modes closest to the EoR window.

The authors wish to thank Adrian Liu for helpful conversations and our referee for a number of helpful suggestions that improved the paper. J.C.P. is supported by an NSF Astronomy and Astrophysics Fellowship under award AST-1302774. This scientific work makes use of the Murchison Radio-astronomy Observatory, operated by CSIRO. We acknowledge the Wajarri Yamatji people as the traditional owners of the Observatory site. Support for the MWA comes from the U.S. National Science Foundation (grants AST-0457585, PHY-0835713, CAREER-0847753, and AST-0908884), the Australian Research Council (LIEF grants LE0775621 and LE0882938), the U.S. Air Force Office of Scientific Research (grant FA9550-0510247), and the Centre for All-sky Astrophysics (an Australian Research Council Centre of Excellence funded by grant CE110001020). Support is also provided by the Smithsonian Astrophysical Observatory, the MIT School of Science, the Raman Research Institute, the Australian National University, and the Victoria University of Wellington (via grant MED-E1799 from the New Zealand Ministry of Economic Development and an IBM Shared University Research Grant). The Australian Federal government provides additional support via the Commonwealth Scientific and Industrial Research Organization (CSIRO), National Collaborative Research Infrastructure Strategy, Education Investment Fund, and the Australia India Strategic Research Fund, and Astronomy Australia Limited, under contract to Curtin University. We acknowledge the iVEC Petabyte Data Store, the Initiative in Innovative Computing and the CUDA Center for Excellence sponsored by NVIDIA at Harvard University, and the International Centre for Radio Astronomy Research (ICRAR), a Joint Venture of Curtin University and The University of Western Australia, funded by the Western Australian State government.