Abstract

We present a new upper limit to the 21-cm power spectrum during the Epoch of Reionization (EoR) which constrains reionization models with an unheated IGM. The GMRT-EoR experiment is an ongoing effort to make a statistical detection of the power spectrum of 21-cm neutral hydrogen emission at redshift z∼ 9. Data from this redshift constrain models of the EoR, the end of the Dark Ages arising from the formation of the first bright UV sources, probably stars or mini-quasars. We present results from approximately 50 h of observations with the Giant Metrewave Radio Telescope in India from 2007 December. We describe radio-frequency interference (RFI) localization schemes which allow bright sources on the ground to be identified and physically removed in addition to automated flagging. Singular-value decomposition is used to remove the remaining broad-band RFI by identifying ground sources with large eigenvalues. Foregrounds are modelled using a piecewise linear filter and the power spectrum is measured using cross-correlations of foreground-subtracted images.

1 INTRODUCTION

Between the recombination epoch at z∼ 1100 and the first round of star formation at z∼ 10 the Universe was filled with neutral hydrogen. This neutral gas is thought to have produced 21-cm hyperfine emission with an effective continuum brightness temperature between −500 and 30 mK (Furlanetto, Oh & Briggs 2006) using z = 8.6 and WMAP 7-yr cosmological parameters (Komatsu et al. 2011), which have formal error bars of roughly 10 per cent. As the first stars formed, beginning at local peaks in the matter density, the hydrogen gas was locally ionized. This is the start of a period of cosmic evolution called the Epoch (or Era) of Reionization (EoR). Over time the ionized cells grew and overlapped, creating a patchwork of ionized and neutral cells. This general patchy topology is well motivated by theory and simulations (e.g. Furlanetto, Zaldarriaga & Hernquist 2004; Iliev et al. 2006; McQuinn et al. 2007; Zahn et al. 2007; see Trac & Gnedin 2009 for a review) though the exact properties are poorly constrained. Eventually, by z∼ 6, the ionization of the broadly distributed material was complete, leaving only rare pockets of neutral gas.

In addition to directly constraining the redshift of last scattering zls∼ 1100, cosmic microwave background (CMB) polarization data have been used to constrain the redshift to the range 8.0 < z < 12.8 (WMAP7; Komatsu et al. 2011) under the assumption that it was instantaneous, although evidence from the CMB indicates that reionization was an extended process (Dunkley et al. 2009). There is also substantial information on the neutral content of the intergalactic medium (IGM) since z∼ 6 via observations of the Gunn–Peterson trough (Gunn & Peterson 1965) in quasar spectra (Becker et al. 2001; Djorgovski et al. 2001; Fan et al. 2006; Willott et al. 2007). Lyman alpha (Lyα) absorption provides a sensitive probe of neutral gas density and can be used to constrain reionization (e.g. Fan et al. 2002). With these two techniques, we can observe the start and the end of the neutral era. However, neither CMB nor quasar observations allow detailed examination of the reionization era itself. For this the use of 21-cm fluctuations in the brightness temperature of neutral hydrogen has long been proposed (Sunyaev & Zeldovich 1975; Hogan & Rees 1979). The 21-cm power spectrum, resulting from a combination of the patchy ionized and neutral medium and the underlying mass power spectrum, is generally considered one of the most promising signals (Scott & Rees 1990; Zaldarriaga, Furlanetto & Hernquist 2004) and much attention has been paid in the literature towards designing suitable experiments to detect it (e.g. Morales & Hewitt 2004; Morales 2005; Bowman, Morales & Hewitt 2006; Harker et al. 2010).

Here the observing frequency is 1420 MHz(1 +z)−1, in the high frequency (VHF) range of the radio spectrum near 150 MHz. Several programmes are underway to study the EoR in the VHF band including LOFAR1 (Kassim et al. 2004; Röttgering et al. 2006; for the EoR case, see e.g. Zaroubi & Silk 2005; Harker et al. 2010), MWA2 (Lonsdale et al. 2009; for the EoR case, see e.g. Bowman et al. 2006; Morales et al. 2006a; Lidz et al. 2008), PAPER3 (Parsons et al. 2010), 21CMA4 (also known as PaST; Peterson, Pen & Wu 2004), in addition to the current GMRT5 programme, and the first goal of such efforts is to determine the redshift at which roughly half the volume of the Universe is ionized. Because the 21-cm emission is an isolated single line, these observations can provide three-dimensional information over a large redshift range (Loeb & Zaldarriaga 2004; Furlanetto et al. 2006), allowing a broad search for the half-ionization redshift. In the future, 21-cm tomography has the potential to add strong constraints on cosmological parameters (e.g. McQuinn et al. 2006; Cooray, Li & Melchiorri 2008; Mao et al. 2008; Furlanetto et al. 2009; Masui, McDonald & Pen 2010).

The brightness temperature of the 21-cm line relative to the CMB is determined by the ionization fraction of hydrogen and the spin temperature of the neutral population, which is in turn governed by the background radiation and the kinetic temperature of the gas (Purcell & Field 1956; Field 1959; Furlanetto et al. 2006). Reionization requires a minimum expenditure of 13.6 eV of energy per hydrogen atom. In contrast, Lyα photons, when absorbed by a neutral atom, are quickly re-emitted. Lyα photons are said to undergo ‘resonant scattering’ and each (multiply scattered) photon can affect many atoms before cosmic redshifting makes them ineffective (Wouthuysen 1952; Field 1959; Chuzhoy & Shapiro 2006). This means Lyα pumping of the hyperfine transition requires only about 1 per cent of the UV flux required for ionization. Assuming a gradual increase in UV flux with time, well before flux levels for reionization are reached, Lyα pumping will couple the spin temperature to the kinetic temperature of the gas. If there was no source of heat in that era other than a weak UV flux, the gas kinetic temperature must have been at its adiabatic expansion value of 1.7 K at z = 8.5, and neutral gas will produce a signal due to the absorption of CMB photons in excess of stimulated emission, and ionized structures would be seen as low-brightness regions in the sky (Chen & Miralda-Escudé 2004). In such a cold-gas model the brightness temperature of the neutral gas against the CMB can be as low as −500 mK.

A small fraction of the mass will have collapsed into minihaloes, which have a temperature higher than this adiabatic temperature (Shapiro et al. 2006). This fraction may even be substantially smaller due to non-perturbative velocity flows (Dalal, Pen & Seljak 2010; Tseliakhovich & Hirata 2010). No study was found on the impact of these non-linear effects on the Lyα pumped IGM temperature. In linear order, the positive and negative over- and underdensities cancel, and the non-linear collapse fraction is still small, so we expect the realistic value to be similar to the adiabatic prediction.

Alternatively, X-rays from supernovae or QSOs might have heated the IGM above the CMB temperature of ≈30 K before reionization (Madau, Meiksin & Rees 1997; Chen & Miralda-Escudé 2004). X-ray heating between 30 and 10 000 K will result in a largely neutral but warm IGM, which would be seen in emission. In the limit that Ts≫ 30 K, the volume emissivity becomes independent of temperature. Patchy X-ray heating can also result in large angular scale structure (Alvarez, Pen & Chang 2010). The sky brightness temperature of the neutral gas has an asymptote at ∼30 mK (Furlanetto et al. 2006).

The cosmic luminosity of X-rays at z∼ 9 is not known and difficult to estimate (see e.g. Dijkstra, Haiman & Loeb 2004; Salvaterra, Haardt & Volonteri 2007; Guo et al. 2009). If the rate of core collapse SNe at high redshift matches that of today, the X-ray output at z∼ 9 would have been sufficient to raise the IGM temperature above the CMB temperature at the onset of reionization. However, the mechanism of and factors affecting core-collapse are not known and most numerical models appear to generically not result in SNe at all (Mezzacappa 2005). It is possible that at high z, core collapse SNe were not as abundant as today, and the IGM was still in absorption during the EoR. We therefore consider two limits: one where the IGM is still cold and Tb=−500 mK, and the other where it is heated above CMB and Tb = 30 mK. A general parametrization of these scenarios was recently proposed by Pritchard & Loeb (2010).

The GMRT-EoR project uses the Giant Metrewave Radio Telescope (GMRT; Ananthakrishnan 1995) located near Pune, India, to make measurements of the power spectrum of the neutral hydrogen signal with the hope of characterizing the structure in the range of 8.1 < z < 9.2. The FWHM of the GMRT primary beam at 150 MHz is forumla which provides a cylindrical comoving survey volume of (280 h−1 Mpc)3, with about equal dimensions in three directions. The primary sensitivity comes from the compact central core which is contained within ≈1 km, or 30 dish diametres, which gives images with ∼302 resolution elements. As recorded, the data have high spectral resolution along the line of sight, but we bin the data for comparable transverse and radial resolution. Each data cube has ∼303≈ 30 000 resolution elements so even though the signal-to-noise ratio of each element is less than 1, the spatial power can be measured to 1 per cent. A statistical measurement of the power spectrum provides the best hope of pinning down the EoR half-ionization redshift.

In Section 2 we describe the data and analysis, including radio-frequency interference (RFI) removal and foreground filters used. In Section 3 we extend the analysis to cross-correlations of multiple nights, and a measurement of the power spectrum. We then conclude in Section 4. All distances are in comoving coordinates.

2 OBSERVATIONS

2.1 GMRT and data description

The GMRT is a radio interferometer consisting of 30 antennas, each with a diametre of 45 m. 14 of these are arranged in a dense central core within 1 km which allows the high brightness sensitivity required to search for the dim EoR signal (Pen et al. 2009). The longest separation between antennas is about 25 km. For this experiment the telescope was operated at 150 MHz with a 16.7 MHz bandwidth. For 21-cm (1421 MHz) emission, this frequency range probes a redshift range of 8.1 < z < 9.2. The highest angular resolution at this observing frequency at GMRT is about 20 arcsec.

Each antenna provides a pair of left- and right-circularly polarized signals, which are passed through a new signal processing system, built in part for this project, called the GMRT Software Backend (GSB; Roy et al. 2010). The final recorded visibilities have a resolution of 7.8 kHz and 1/4 s in frequency and time, respectively.

The field is centred on the pulsar B0823+26. This pulsar has a period of 0.53 s and an average flux of 350 mJy at 150 MHz (Hobbs et al. 2004). It is situated in a relatively cold part of the sky at galactic latitude ≈30° with few nearby bright sources. The on-pulse flux is about 6 Jy, brighter than all other sources in the field, making it a good calibrator. The calibration of these data was described in Pen et al. (2009) and is summarized below. The positions of the brightest sources in the field are shown in Fig. 1.

Figure 1

The 12 brightest sources in the field used for these observations. Pulsar PSR B0823+26 is at the centre of the field. The circle denotes the half-power diametre of the main beam.

Figure 1

The 12 brightest sources in the field used for these observations. Pulsar PSR B0823+26 is at the centre of the field. The circle denotes the half-power diametre of the main beam.

The data are folded in time into 16 ‘gates’ such that the pulse from the pulsar is contained within one gate. The period of the pulsar is much shorter than the time-scale at which ionospheric fluctuations dominate, so by comparing the ‘on’ gate with the neighbouring ‘off’ gates, everything in the field that is constant over the period of the pulsar can be removed. This includes sky sources and most RFI. By comparing the pulsar signal received by each antenna, the relative system gain of each antenna can be calibrated. This technique allows calibration of both phase and polarization in real time, with 100 per cent observing efficiency, since no additional time is needed for phase or flux calibration. Errors in clock synchronization sometimes cause the pulsar pulse to straddle two gates. In these cases it is necessary to include the flanking gates when identifying the pulsar signal. This is not ideal since the pulsar signal is diluted over multiple gates, but the correction is limited to only the portion of observations that require it.

Since the pulsar amplitude varies from pulse to pulse, the absolute system gain needs to be measured separately. This is done using a noise injection system, which the GSB decodes to calibrate the absolute gain of the system, while a sky radio source is used to transfer the noise source calibration on to the sky. The primary flux calibrator used was 5C 7.245, a radio galaxy at z≈ 1.6 (Willott, Rawlings & Blundell 2001), located within the field of view (see Fig. 1). However, since this radio galaxy is an extended object, it can only be used as a calibrator for short baselines where its structure is not resolved. For antennas in longer baselines we determine the relative gain calibration using the pulsar.

2.2 Data analysis

Data from 2007 December 10, 11, 14, 16, 17 and 18 are included in the current analysis. A software pipeline has been developed to automate the calibration and interference removal steps described below. The current configuration takes approximately 11 h to process 1 h of data on the CITA Sunnyvale computing cluster, or 20 h on the National Centre for Radio Astronomy (NCRA) HP computing cluster in Pune, India. The NCRA facility can process 4 h in parallel, completing one night of observations in slightly less than 2 d. While the Sunnyvale cluster has the capacity to process as much as 12 h of data in parallel, storage capacity, shared usage, and other bottlenecks reduce this significantly.

One of the limiting factors for measurement of the EoR signal is broad-band RFI, which dominates the signal at 150 MHz. In addition to the standard procedure of flagging bad data, we also use a singular value decomposition (SVD) to remove broad-band RFI from the data, and also to identify and remove interference sources. We believe this approach is unique among other RFI mitigation strategies in the radio community.

Narrow-line interference was removed by masking points in each frequency bin with an intensity above some threshold. This was done twice in the data reduction pipeline. Input data are initially flagged with a threshold of 8σ on a Gaussian scale, and then again at 3σ after broad-band interference was removed. The first mask cannot be too aggressive or the techniques to remove broad-band interference, discussed below, will fail.

SVD is used to separate broad-band radio sources on the ground from those in the sky. Ground-based sources contribute most to the largest eigenvalues since they do not move as a function of time with respect to the array, while sky sources rotate. 1 h of data has 14 396 time records, and each record has about 7.5 million entries corresponding to the number of frequency channels and baselines between the 60 antennas. This is treated as a matrix with each time record as one row. The 50 largest eigenvalues are identified through an SVD and flagged as noise to be removed. A sample of the data at a few intermediate stages showing the successful removal of both line and broad-band RFI can be seen in Fig. 2. The RFI patterns in (u, v) space both before and after the SVD are illustrated in Fig. 3.

Figure 2

Representative stages in the data reduction pipeline. In each panel, the horizontal axis is frequency, covering approximately 1 MHz, and the vertical is time, increasing downward and covering approximately 1 h. The grey-scale is the cross-power spectral density. The top row is the C0–C8 baseline and the bottom row the C0–W4 baseline, which are approximately 560 and 9400 m, respectively. The large bright patches indicate broad-band RFI, and the vertical lines indicate line RFI. The first panel in each row is the initial input data. The second is after the initial 8σ mask, with most interference still visible. The third panel is after removing the largest eigenvalues in the SVD. The broad-band interference is no longer visible. The last panel shows the final 3σ mask removing the line RFI, leaving a nearly uniform image.

Figure 2

Representative stages in the data reduction pipeline. In each panel, the horizontal axis is frequency, covering approximately 1 MHz, and the vertical is time, increasing downward and covering approximately 1 h. The grey-scale is the cross-power spectral density. The top row is the C0–C8 baseline and the bottom row the C0–W4 baseline, which are approximately 560 and 9400 m, respectively. The large bright patches indicate broad-band RFI, and the vertical lines indicate line RFI. The first panel in each row is the initial input data. The second is after the initial 8σ mask, with most interference still visible. The third panel is after removing the largest eigenvalues in the SVD. The broad-band interference is no longer visible. The last panel shows the final 3σ mask removing the line RFI, leaving a nearly uniform image.

Figure 3

Raw visibilities with (u, v) distance |b| < 200 before the SVD RFI removal step (left) and after the RFI removal (right). The bright patterns in the central vertical strip (low |u|) are caused by RFI, and are significantly reduced by the SVD procedure. The effect of this process on the final power spectrum is considered in Section 3.2.

Figure 3

Raw visibilities with (u, v) distance |b| < 200 before the SVD RFI removal step (left) and after the RFI removal (right). The bright patterns in the central vertical strip (low |u|) are caused by RFI, and are significantly reduced by the SVD procedure. The effect of this process on the final power spectrum is considered in Section 3.2.

2.3 Physical removal of radio interference sources

While some RFIs can be removed in the data analysis, this risks removing sky signal as well. Ideally, one would prevent RFIs from occurring at all by identifying and correcting the physical sources. We can take advantage of the fact that as an interferometer GMRT is able to make images of the near-field, and create maps of bright RFI sources near the GMRT itself.

Candidate sources detected by a single baseline appear as a hyperbola of equal light arrival time in near-field images of a single SVD mode. When a source is detected by many baselines, the corresponding hyperbolas intersect at a single point in the image. These near-field images, an example of which is shown in Fig. 4, become the ‘RFI maps’ which are used to isolate interference sources. For sources with a high enough duty cycle, the GPS position given by the RFI maps are accurate enough that a handheld yagi antenna and portable radio receiver can be used to find the precise source. With the handheld antenna, the direction to a source can be determined by ear to better than 30°, and thereafter localized by successive triangulation.

Figure 4

(a) A near-field image, measuring about 40 km on each side and centred on GMRT. An RFI source identified through large SVD modes is clearly visible as a bright spot in the upper right. (b) A map of GMRT covering the same area as in (a). Each black circle indicates an antenna location, and the open circle corresponds to the RFI source, #14 in Table 1, identified in (a).

Figure 4

(a) A near-field image, measuring about 40 km on each side and centred on GMRT. An RFI source identified through large SVD modes is clearly visible as a bright spot in the upper right. (b) A map of GMRT covering the same area as in (a). Each black circle indicates an antenna location, and the open circle corresponds to the RFI source, #14 in Table 1, identified in (a).

By comparing the position of a noise transmitter relative to the candidate interference source both physically and in the RFI map, we can confirm whether we have identified the correct source. Occasional misidentifications are likely due to the intermittent nature of the sources; there is no guarantee that the radiating source in the images will still be active when an attempt at identification is made. Since the beginning of this effort, calibration has been improved to locate sources within about 100 m. Sources could in principle be located using only the array with a precision 10 m, but the accuracy of available maps, GPS equipment and other factors limit the real-world precision.

Candidate RFI sources identified so far include transformers, power line junctions and loose wires in contact with power lines. Table 1 lists some candidate sources recorded on 2008 December 6, and descriptions. In 2009 February, GMRT began a collaboration effort with MSEDCL, which controls the power transmission lines in the area around GMRT, to assist us in removing local interference sources. A wire hanging over 132-kV high-tension power lines (#1 in Table 1) which was by far the brightest interference source within 20 km, appearing in 30 of the 50 largest eigenmodes, was removed in 2009 February. Other bright sources have also been removed since then. Although these sources are still present in the raw data from 2007, this procedure has shown that the SVD algorithm applied in the present work is successful in identifying real RFI sources during analysis.

Table 1

List of candidate RFI sources based on observations on 2008 December 6, identified by the strongest SVD mode number from which the GPS coordinates are derived. The two brightest sources, modes number 1 and 9, were removed in early 2009.

Mode Longitude Latitude Description 
74.034 653 19.147 676 Conclusively identified as a wire hanging over 132-kV high-tension power lines, clearly visible using the handheld yagi antenna several kilometres away. Wire removed with assistance from MSEDCL on 2009 February 26. 
74.085 335 19.171 797 Identified with wires hanging from an unused telephone pole, positioned directly under high-tension lines. 
13 74.034 271 19.037 518 Source in the area transmitted very intermittently, making identification difficult. Potentially a transformer 20-m north of coordinate. 
14 74.096 581 19.209 597 Identified as a small wire hanging on a 500 kV power line, removed in February 2010. 
15 74.120 018 19.181 877 Two possible sources, both transmission towers, separated by about 300 m. The one closest to the coordinate is a T-junction of two high-voltage lines. 
17 74.072 952 19.088 457 Initially identified to be a small pump approximately 150-m west of this coordinate, though other candidates include two nearby transformers. 
33 74.075 623 19.108 976 Transformer approximately 100-m northwest, unambiguously radiating at low levels. 
Mode Longitude Latitude Description 
74.034 653 19.147 676 Conclusively identified as a wire hanging over 132-kV high-tension power lines, clearly visible using the handheld yagi antenna several kilometres away. Wire removed with assistance from MSEDCL on 2009 February 26. 
74.085 335 19.171 797 Identified with wires hanging from an unused telephone pole, positioned directly under high-tension lines. 
13 74.034 271 19.037 518 Source in the area transmitted very intermittently, making identification difficult. Potentially a transformer 20-m north of coordinate. 
14 74.096 581 19.209 597 Identified as a small wire hanging on a 500 kV power line, removed in February 2010. 
15 74.120 018 19.181 877 Two possible sources, both transmission towers, separated by about 300 m. The one closest to the coordinate is a T-junction of two high-voltage lines. 
17 74.072 952 19.088 457 Initially identified to be a small pump approximately 150-m west of this coordinate, though other candidates include two nearby transformers. 
33 74.075 623 19.108 976 Transformer approximately 100-m northwest, unambiguously radiating at low levels. 

2.4 Single-night images

After RFI removal has been completed in each 1-h scan, polarization calibration is used to combine these into images of an entire night, typically of about 8 h, as can be seen in Fig. 5. This step accounts for leakage between the left and right polarization signals, as well as effects caused by the relative rotation of the array and sky source by taking into account the changing parallactic angle with time.

Figure 5

Dirty sky image made with 8 h of data from 2007 December 10. The maximum (u, v) distance is 2000 with grid size of approximately Δu = 8 giving a field of view of 7forumla3 on each side. This same (u, v) distance cut is used to calibrate the peak flux of 5C 7.245 to 1.6 Jy for each day. The rms of this image is 44 mJy within the primary beam.

Figure 5

Dirty sky image made with 8 h of data from 2007 December 10. The maximum (u, v) distance is 2000 with grid size of approximately Δu = 8 giving a field of view of 7forumla3 on each side. This same (u, v) distance cut is used to calibrate the peak flux of 5C 7.245 to 1.6 Jy for each day. The rms of this image is 44 mJy within the primary beam.

Although the array elements very nearly lie in a single plane, the rotation of the Earth turns this plane over time with respect to the line of sight. To correctly stack many hours together, one must consider the frequency and directional dependence of the measurement. One effect to consider is the change of the primary beam with frequency. For a given frequency, the primary beam is well defined and independent of (u, v). Deconvolution is difficult, but one can pick the frequency with the smallest beam and restrict the field of view at other frequencies to match. Computationally, this can be done by convolving with the ratio of the beams, which is possible because of the commutativity of the primary beam operators.

However, this is a small effect compared to the effects of a w term. Typically, relative antenna positions are considered to be on a two-dimensional plane with coordinates parametrized by (u, v) in units of the observing wavelength. Since baselines become non-coplanar over the course of the night, a third dimension, w, must be included, without which the image becomes blurred.

For a thorough explanation of w term issues, see for example Cornwell, Golap & Bhatnagar (2008) and, for a more general overview, Thompson, Moran & Swenson (2001). As the observing frequency changes across the band, a single baseline samples different points in (u, v) space. At any given frequency, the (u, v) plane is typically sparsely sampled, but with many frequency channels there will often be data at a single (u, v) point for at least a few channels. One can subtract the average from these, and under the assumption that foreground sources have the same spectral index as the calibrator, this removes the foregrounds very well. Unfortunately, when considering the full (u, v, w) cube, the data are too sparse for this to work. Since the w coordinate is dependent on the position in (u, v) space, and can become very large, the effects of the w term cannot be corrected at different frequencies in the same way as the beam.

Strategies for correcting this effect for GMRT are being developed, which build on those used in CMB studies (Hobson & Maisinger 2002; Myers et al. 2003), but we have not yet applied these to our data. In the interim we restrict ourselves to the short baselines for which w is small and such corrections are not needed. The short baselines are the ones most sensitive to the EoR signal, so we retain most of the EoR sensitivity of the array.

The flux scale of the field is set by calibrating relative to the pulsar at the centre of the field for every baseline, receiver and frequency. Since the pulsar flux is known to be variable, to set the absolute flux scale we identify the source with the highest flux in the sky image, and set this peak value equal to the known flux of that source. The exact value found depends on the angular resolution, determined by the maximum baseline length used. For all days analysed in this work, the brightest apparent source is the radio galaxy 5C 7.245. The flux of this source was measured as 215.1 ± 12.6 mJy at 1415 MHz by Willis, Oosterbaan & de Ruiter (1976) at the Westerbork Radio Telescope, and as 685 ± 47 mJy at 408 MHz in the 5C catalogue (Pearson & Kus 1978) for a spectral index of 0.93 ± 0.07. Following the same argument as Pen et al. (2009), to achieve the best match with other bright sources in the field, we adopt a value of 1.6 Jy for 5C 7.245. Although the flux of this source changes by as much as 10 per cent across the band, this change is less than that due to the uncertainty in the spectral index. Additionally, this galaxy has two components separated by 12 arcsec (Willott et al. 2001). GMRT is capable of resolving this structure with |b| ≳ 2600, where b = (u, v) is the distance between the two antennas at the baseline in units of the observing wavelength. To calibrate the flux we require a point source so that the flux is not spread over multiple image pixels. This is achieved by limiting the maximum baseline length used while calibrating to |b| ≤ 2000. At much smaller |b|, we become similarly limited by confusion.

Removing the foregrounds adequately is essential for detecting the EoR signal (see e.g. Wang et al. 2006) and much work has been done on simulating foregrounds (e.g. Bowman, Morales & Hewitt 2009; Jelić et al. 2010) and designing removal strategies (Morales, Bowman & Hewitt 2006b; Harker et al. 2009; Liu et al. 2009a, and others). For a review, see Morales & Wyithe (2010) and references therein. To remove foreground sources, we take advantage of the fact that the flux of such sources do not vary greatly with frequency. A signal originating from the EoR will appear as additional variation on top of the foreground signal. To model spectrally smooth sources, for every time-stamp and baseline we subtract a piecewise linear fit between the median fluxes in frequency bins of 8, 2 or 0.5 MHz. This is simple to implement and has the advantage over a polynomial fit of a having straightforward and local window function which is simple to interpret. Foregrounds which do not vary over this range are removed, while features with variability at higher frequencies remain. This method results in an upper limit to the EoR signal since the measured power may still include residual foreground variation. By subtracting the mean flux over 2 MHz to remove galactic foregrounds, noise levels can be lowered to about 2 mJy for maps of approximately 10°. Bright point sources may play a significant role (Di Matteo, Ciardi & Miniati 2004; Liu, Tegmark & Zaldarriaga 2009b; Datta, Bowman & Carilli 2010), but are not treated separately here.

To first order, this is the same as a boxcar average6 which we can write as  

1
formula
where D(ν) is the input data, forumla is the data after filtering and the window function is  
2
formula
Here, δ(x) is the Dirac delta function and H(x, Δν) is a step function centred at x with width Δν corresponding to our chosen filter, and of unit area.

In Fourier space, this window function is  

3
formula
shown in Fig. 6, where we have denoted the wave number by kν to make explicit that it is in the units of inverse frequency. This can be rewritten in terms of k using the fact that in our redshift range 1 MHz ≈ 11.6 h−1 Mpc. When converted to units of h Mpc−1, we write the wave number as k to emphasize that the window function acts on structure along the line of sight only.

Figure 6

Fourier transform of our chosen window function for Δν = 0.5 MHz (solid line). Dashed lines indicate where w(k) = 0.5, which sets the minimum k to which we are sensitive. This value scales inversely with Δν. The peak of w(k)2/k3 (the dot–dashed line) indicates the k to which we are most sensitive under the assumption that the power spectrum is proportional to k−3.

Figure 6

Fourier transform of our chosen window function for Δν = 0.5 MHz (solid line). Dashed lines indicate where w(k) = 0.5, which sets the minimum k to which we are sensitive. This value scales inversely with Δν. The peak of w(k)2/k3 (the dot–dashed line) indicates the k to which we are most sensitive under the assumption that the power spectrum is proportional to k−3.

After filtering in this way, we describe the kν for which w(kν) = 0.5 as corresponding to the minimum k along the line of sight to which we are sensitive, while smaller k are removed by the filter and will not contribute as strongly to the power spectrum. The three filters of 8, 2 and 0.5 MHz correspond to a minimum k of 0.04, 0.16 and 0.65 h Mpc−1, respectively. Since the power spectrum is thought to be proportional to k−3 (Iliev et al. 2008), the line w(k)2/k3 indicates the k to which we are most sensitive.

The 0.5 MHz foreground filter reduces the peak and rms values by a factor of 50, as can be seen in Fig. 7. The filter is most effective within the primary beam. Fig. 8 demonstrates the dominance of RFI in the filtered dirty maps, and the improvement that the SVD step provides. Without the SVD, the maps in Fig. 7 would have been a factor of 4 noisier.

Figure 7

Data from December 10. The top row is the sky image with a 11.4° field of view, with |b| < 200 and |u| > 60 binned in the (u, v) plane by Δu = 5. The bottom row is the visibilities in the same range with Δu = 0.4 to show structure. The left-hand column is before any foreground subtraction. In this image the dominant source is B2 0825+24 (or 4C 24.17) just south of the FWHM of the primary beam, with a peak value of 2.2 Jy. rms within the beam is 343 mJy. The right-hand column is after a 0.5 MHz subtraction. The peak value of this image is 47 mJy with an rms of 6.2 mJy, lower by a factor of about 50. If put on the same grey-scale as the image without foreground subtraction no features would be visible. The dominant source after this filter is 3C 200, well outside the beam.

Figure 7

Data from December 10. The top row is the sky image with a 11.4° field of view, with |b| < 200 and |u| > 60 binned in the (u, v) plane by Δu = 5. The bottom row is the visibilities in the same range with Δu = 0.4 to show structure. The left-hand column is before any foreground subtraction. In this image the dominant source is B2 0825+24 (or 4C 24.17) just south of the FWHM of the primary beam, with a peak value of 2.2 Jy. rms within the beam is 343 mJy. The right-hand column is after a 0.5 MHz subtraction. The peak value of this image is 47 mJy with an rms of 6.2 mJy, lower by a factor of about 50. If put on the same grey-scale as the image without foreground subtraction no features would be visible. The dominant source after this filter is 3C 200, well outside the beam.

Figure 8

Dirty sky images before the SVD RFI removal step (left) and after the RFI removal (right) with |b| < 200, after the 0.5 MHz foreground subtraction. The field of view is 11forumla4. RFI clearly dominates the map on the left, with peak value of 73 mJy and an rms of 18 mJy. After RFI removal, the peak drops to 32 mJy with an rms of 4.2 mJy.

Figure 8

Dirty sky images before the SVD RFI removal step (left) and after the RFI removal (right) with |b| < 200, after the 0.5 MHz foreground subtraction. The field of view is 11forumla4. RFI clearly dominates the map on the left, with peak value of 73 mJy and an rms of 18 mJy. After RFI removal, the peak drops to 32 mJy with an rms of 4.2 mJy.

3 MULTIDAY ANALYSIS

3.1 Differences between nights

Since RFI and foreground signals are so much larger than the EoR signal, small errors in the subtraction of foregrounds could easily result in a spurious signal. To avoid this we use cross-correlations between multiple days to make a statistical measurement.

To gauge how successful a cross-correlation measurement might be, we are interested in the relative similarity of the different nights, which we can gauge by taking the difference of visibilities. Days which subtract well will show mostly noise, while days which subtract poorly will still show evidence of bright sources. Visibilities are only used when data exist at the same position in both days. This serves as a test of the importance of ionospheric fluctuations. Though the ionosphere is calibrated along the line of sight to the pulsar at the centre of the field, the ionosphere could change across the field of view. Pairs of days for which this change away from the field of view is different will not subtract well. The results of the subtractions are shown in Table 2.

Table 2

Peak flux and rms of the difference of each 2-d pair available, in mJy, with a maximum (u, v) distance of 600, which corresponds to a maximum baseline of 1.2 km. All values are in mJy. To remove foregrounds, a 2-MHz linear filter was applied.

  Unfiltered With 2-MHz filter 
Subtracted pair Peak flux rms Peak flux rms 
December 10 December 11 1856.1 227.3 137.9 19.9 
December 10 December 14 888.4 185.6 62.5 11.6 
December 10 December 16 278.4 49.8 27.1 3.2 
December 10 December 17 447.9 64.9 26.2 3.6 
December 10 December 18 548.3 70.4 32.5 4.2 
December 11 December 14 1037.8 145.1 256.1 56.5 
December 11 December 16 2148.3 256.7 158.1 23.9 
December 11 December 17 2528.2 318.5 106.1 21.1 
December 11 December 18 2997.5 356.1 159.1 12.4 
December 14 December 16 1221.5 193.9 50.2 7.0 
December 14 December 17 1311.1 257.2 50.4 7.6 
December 14 December 18 1433.5 233.0 67.1 4.9 
December 16 December 17 202.0 78.9 23.8 2.8 
December 16 December 18 224.2 31.2 24.4 3.4 
December 17 December 18 206.9 60.6 24.6 3.2 
  Unfiltered With 2-MHz filter 
Subtracted pair Peak flux rms Peak flux rms 
December 10 December 11 1856.1 227.3 137.9 19.9 
December 10 December 14 888.4 185.6 62.5 11.6 
December 10 December 16 278.4 49.8 27.1 3.2 
December 10 December 17 447.9 64.9 26.2 3.6 
December 10 December 18 548.3 70.4 32.5 4.2 
December 11 December 14 1037.8 145.1 256.1 56.5 
December 11 December 16 2148.3 256.7 158.1 23.9 
December 11 December 17 2528.2 318.5 106.1 21.1 
December 11 December 18 2997.5 356.1 159.1 12.4 
December 14 December 16 1221.5 193.9 50.2 7.0 
December 14 December 17 1311.1 257.2 50.4 7.6 
December 14 December 18 1433.5 233.0 67.1 4.9 
December 16 December 17 202.0 78.9 23.8 2.8 
December 16 December 18 224.2 31.2 24.4 3.4 
December 17 December 18 206.9 60.6 24.6 3.2 

It can be seen in Table 2 that December 11 gives consistently poor results in both the unfiltered and filtered images, and was thus excluded from all subsequent analyses. From these maps we can also conclude that RFI is not dominating the differences between one day and the other. If this had been the case, since RFI is typically isolated in (u, v) space, it would be visible across the whole sky image, including far outside the primary beam. However, it can be seen in Fig. 7 that structure decreases rapidly away from the centre of the field, indicating an astronomical source. This is true of all subtraction pairs.

3.2 Cross-correlations and power spectra

The power spectrum of sky structure can be determined directly from the visibilities (Zaldarriaga et al. 2004). To find the cross-power of 2 d, we take the product of the visibilities after gridding with a cell spacing equal to the size of the beam (Δu = 20). Calculating the power spectrum from the visibilities instead of from the correlations in the sky image takes advantage of the fact that the visibilities have a nearly diagonal correlation matrix in the noise (White et al. 1999). We then find the weighted average of visibilities in the annuli of (u, v) space which gives the power in units of Jy2. Since the amount of data at large |b| decreases, we increase the width of each successive annulus by 60 per cent with increasing |b|. The smallest binwidth is equal to the size of the beam. This prevents a large artificial variability in power due to sparse sampling. Visibilities are weighted ‘naturally’, by the inverse of the noise. Since we expect our sensitivity to the EoR signal to diminish rapidly with increasing baseline length, we look only at the first few points, averaged over all possible cross-correlations. Additionally, it is known that the SVD will introduce a loss of power at low |u|. To avoid this, we impose a limit of |u| > 60 when taking the cross-correlations, determined by requiring that the power spectrum before and after the SVD differ by less than 1σ as shown in Fig. 9. The part of the (u, v) plane that is lost with this cut can be seen in Fig. 3.

Figure 9

Comparison of the cross-power of December 10 with all other days under four different conditions. The dot–dashed and double-dot–dashed lines are before and after the SVD RFI removal, respectively, including all u. It can be seen that power is lost in the SVD. Similarly, the solid and dashed lines are before and after the SVD step, respectively, this time with a |u| > 60 limit imposed. In this case the two lines are almost identical, meaning the SVD had little effect on the total power.

Figure 9

Comparison of the cross-power of December 10 with all other days under four different conditions. The dot–dashed and double-dot–dashed lines are before and after the SVD RFI removal, respectively, including all u. It can be seen that power is lost in the SVD. Similarly, the solid and dashed lines are before and after the SVD step, respectively, this time with a |u| > 60 limit imposed. In this case the two lines are almost identical, meaning the SVD had little effect on the total power.

The power spectrum of the cross-correlation can be converted to units of K2 using  

formula
where l = 2 π |b|, b = (u, v) is the visibility coordinate in the units of wavelength, Cl is the power measured in K2, θb is the primary beam size and ν is the wavelength (Pen et al. 2009). The quantity in the angle-brackets on the right is equal to the power in Jy2 found above. This conversion is written in terms of the GMRT observations with a primary beam of forumla and ν = 150 MHz. If gridding in the (u, v) plane is too fine, the data become noisy, while very coarse gridding requires the assumption that the data are constant across the whole cell. As mentioned, we use a gridding equal to the size of the beam. A fully optimized estimate would require a maximum likelihood code which is being worked on in the context of the Cosmic Background Imager (CBI) gridder (Myers et al. 2003). Fig. 10 shows the weighted average of all cross-correlation pairs, excluding December 11. Bernardi et al. (2009) reported a power spectrum of foregrounds without subtraction in the galactic plane at a level comparable to our measurement both here and in Pen et al. (2009).

Figure 10

Average power spectrum in units of K2 of all combinations of days, excluding December 11, as a function of the multipole moment l. Each point is shown with a 2σ upper limit derived from a bootstrap error analysis, which is in most cases smaller than the size of the point. The points are logarithmically spaced as described in the text, from left to right covering the ranges 377 < l < 578, 578 < l < 899 and 899 < l < 1414. Triangles are the power before subtracting foregrounds, diamonds are after 8 MHz mean subtraction, squares are after 2 MHz mean subtraction and circles are after 0.5 MHz subtraction. The curved solid line is the theoretical EoR signal from Jelić et al. (2008), and the dashed line is the theoretical EoR signal with a cold absorbing IGM as described in the text.

Figure 10

Average power spectrum in units of K2 of all combinations of days, excluding December 11, as a function of the multipole moment l. Each point is shown with a 2σ upper limit derived from a bootstrap error analysis, which is in most cases smaller than the size of the point. The points are logarithmically spaced as described in the text, from left to right covering the ranges 377 < l < 578, 578 < l < 899 and 899 < l < 1414. Triangles are the power before subtracting foregrounds, diamonds are after 8 MHz mean subtraction, squares are after 2 MHz mean subtraction and circles are after 0.5 MHz subtraction. The curved solid line is the theoretical EoR signal from Jelić et al. (2008), and the dashed line is the theoretical EoR signal with a cold absorbing IGM as described in the text.

We have plotted the power spectra with 2σ bootstrap errors (Efron 1979) which were derived as follows. Using 5 d of data, there are 10 possible cross-correlations. From these, 10 are randomly sampled, with replacement, resulting in a slightly different power spectrum. This is repeated 104 times, and the variance on this set of power spectra is calculated to give the error on the original. Formally this quantifies the error when taking independent samples of a statistical distribution. In our case this is not a rigorous error, but a suitable straightforward estimate given the main complications discussed earlier.

3.3 Comparison to models

Fig. 10 can be compared to simulated results from the Low Frequency Array (LOFAR) EoR project in Jelić et al. (2008), which assumes TsTCMB. At low l, their simulated EoR signal is approximately (10 mK)2, while our lowest point with a similar 0.5 MHz bandwidth filter is (50 mK)2 with a 2σ upper limit of (70 mK)2. These results are comparable to the sensitivities LOFAR expects after 400 h of the EoR project. We have also considered the case where reheating of the IGM does not occur, so the spin temperature remains coupled to the kinetic temperature of the gas (Ciardi & Madau 2003). In this case, the IGM cools adiabatically after decoupling from the CMB at z≈ 150. The temperature fluctuations scale with (1 +TCMB/Ts), and the power scales with the same factor squared. Using Tk=TCMB(1 +z)/150 at z = 8.6, the power becomes approximately 275 times larger. This line is shown in Fig. 10, and is comparable to the data. Scaling up the warm IGM power spectrum from Iliev et al. (2008) or Jelić et al. (2008) in this way is a reasonable approximation to the expected signal in a cold IGM. For a more detailed study of the signal in such an absorption regime, we refer the reader to Baek et al. (2009, 2010).

The power spectrum of reionization is intrinsically three-dimensional (Morales & Hewitt 2004). The strongest constraints on the 3D power spectrum for the Δν = 2 and 0.5 MHz foreground filter case are shown in Fig. 11. This uses k2=k2+k2, where k is given by the windowing function of the filter and k=l/6 h−1 Gpc. When comparing this to the prediction from Iliev et al. (2008), one should note that our windowing function will also reduce the predicted signal by at most a factor of 2.

Figure 11

3D power spectrum for the data shown in Fig. 10, using k2=k2+k2. This is dominated by k, so the binwidth in k does not influence the horizontal position of the limits. The strongest constraints from the 2 and 0.5 MHz filters are shown (square and circle, respectively). Upper limits are 2σ bootstrap errors. Three possible signals are shown. The dashed line is the prediction from Iliev et al. (2008) and the double-dot–dashed line is the same for a cold IGM. The solid line comes from the single-scale bubble model as described in the text for a cold IGM, using k = 2.5/R to show the maximum power at all k. For the two points shown, the bubble diametres which achieve this maximum power are 27 and 7.4 h−1 Mpc, respectively. Only the 0.5-MHz point imposes a limit on the diametre. For a warm IGM case, this signal would be reduced by the same factor as in the two dashed lines.

Figure 11

3D power spectrum for the data shown in Fig. 10, using k2=k2+k2. This is dominated by k, so the binwidth in k does not influence the horizontal position of the limits. The strongest constraints from the 2 and 0.5 MHz filters are shown (square and circle, respectively). Upper limits are 2σ bootstrap errors. Three possible signals are shown. The dashed line is the prediction from Iliev et al. (2008) and the double-dot–dashed line is the same for a cold IGM. The solid line comes from the single-scale bubble model as described in the text for a cold IGM, using k = 2.5/R to show the maximum power at all k. For the two points shown, the bubble diametres which achieve this maximum power are 27 and 7.4 h−1 Mpc, respectively. Only the 0.5-MHz point imposes a limit on the diametre. For a warm IGM case, this signal would be reduced by the same factor as in the two dashed lines.

We also consider an idealized case in which the ionized bubbles during reionization are of uniform scale and non-overlapping. Then for a given k there will be a characteristic bubble radius R at which the power is maximized. By taking the 3D Fourier transform of a perfectly ionized bubble, and requiring that the universe is 50 per cent ionized, one can show that the power is given by  

4
formula
and is maximized when kR≈ 2.5. In this case, k3/(2π2) P(k) ≈T2b/5, where Tb is the brightness temperature ≈30 mK in an X-ray heated IGM or almost -500 mK in a cold absorbing IGM. This signal would be more than an order of magnitude larger than the predictions by e.g. Iliev et al. (2008) or Jelić et al. (2008). In Fig. 11 we have included the power spectrum from this model with k chosen to maximize the power in the range of interest. The data currently impose a limit on the size of bubbles in this single-scale model. Our upper limits with a 0.5 MHz foreground subtraction rule out bubbles with diametres from 2.2 to 12.4 h−1Mpc in the redshift range of 8.1 < z < 9.2. The cold IGM constraint is applicable even in the case of simulations, since only UV photons are included which themselves do not heat the IGM.

4 CONCLUSION

The data analysis has been completed in 6 d from December 2007 with a noise level of approximately 2 mJy on most nights. The SVD removal strategies for broad-band RFI used lower noise by a factor of 4 in temperature, or 16 in power, which flagging alone cannot achieve. We have also tested for ionospheric variations and found that our pulsar calibration is sufficient for dealing with these effects. After RFI removal and foreground subtraction, we have measured a power spectrum which represents a new upper limit on the 21-cm brightness temperature fluctuations during the epoch of reionization. These results can be used to constrain assumptions about the state of the IGM at these times, particularly in the case of a Lyα-pumped, but cold, IGM.

The previous best limit on 21-cm signal at a comparable redshift was obtained by Bebbington (1986), who reported no features down to 5 K at z = 8.4. Parsons et al. (2010) have reported a similar limit of about 5 K using PAPER. The upper limit we present here is approximately 70 mK on the variance in 21-cm brightness temperature at z = 8.6, almost 2 orders of magnitude better than these previous limits. Residual foreground contamination and RFI may still be contributing to this power, but the EoR signal cannot be larger than this.

The analysis of additional observations of the B0823+26 field made since 2007, containing approximately two times more data than has been treated in this work, is going on. We expect to continue to improve these results, and the planned addition of a second field will strengthen this statistical measurement. The dominant uncertainty source remains RFI and foreground modelling.

These results provide a first-look at the progress made at GMRT in detecting EoR. The GMRT-EoR team has been continuing observations, most recently with an additional 150 h allocated in the summer of 2010. We continue to improve both the system temperature of GMRT antennas and the RFI monitoring and removal strategies and expect to improve upon these results as the analysis of the newer data progresses.

6
More precisely, it is a boxcar average with a variable width, and therefore a somewhat different curvature.

We thank Chris Hirata for his contributions to the analysis pipeline and a reading of the paper, and the anonymous reviewer for many useful comments. We also thank the staff of the GMRT that made these observations possible. GMRT is run by the National Centre for Radio Astrophysics of the Tata Institute of Fundamental Research. The computations were performed on CITA's Sunnyvale clusters which are funded by the Canada Foundation for Innovation, the Ontario Innovation Trust and the Ontario Research Fund. The work of U-LP is supported by the National Science and Engineering Research Council of Canada. The map image in Fig. 4 is copyright by OpenStreetMap7 contributors and is used under license.8

REFERENCES

Alvarez
M. A.
Pen
U.
Chang
T.
,
2010
,
ApJ
 ,
723
,
L17
Ananthakrishnan
S.
,
1995
,
JA&A
 ,
16
,
427
Baek
S.
Di Matteo
P.
Semelin
B.
Combes
F.
Revaz
Y.
,
2009
,
A&A
 ,
495
,
389
Baek
S.
Semelin
B.
Di Matteo
P.
Revaz
Y.
Combes
F.
,
2010
,
A&A
 ,
523
,
A4
Bebbington
D. H. O.
,
1986
,
MNRAS
 ,
218
,
577
Becker
R. H.
et al.,
2001
,
AJ
 ,
122
,
2850
Bernardi
G.
et al.,
2009
,
A&A
 ,
500
,
965
Bowman
J. D.
Morales
M. F.
Hewitt
J. N.
,
2006
,
ApJ
 ,
638
,
20
Bowman
J. D.
Morales
M. F.
Hewitt
J. N.
,
2009
,
ApJ
 ,
695
,
183
Chen
X.
Miralda-Escudé
J.
,
2004
,
ApJ
 ,
602
,
1
Chuzhoy
L.
Shapiro
P. R.
,
2006
,
ApJ
 ,
651
,
1
Ciardi
B.
Madau
P.
,
2003
,
ApJ
 ,
596
,
1
Cooray
A.
Li
C.
Melchiorri
A.
,
2008
,
Phys. Rev. D
 ,
77
,
103506
Cornwell
T. J.
Golap
K.
Bhatnagar
S.
,
2008
,
IEEE J. Selected Topics Signal Processing, 2
 ,
647
Dalal
N.
Pen
U.
Seljak
U.
,
2010
,
J. Cosmol. Astropart. Phys.
 ,
11
,
7
Datta
A.
Bowman
J. D.
Carilli
C. L.
,
2010
,
ApJ
 ,
724
,
526
Di Matteo
T.
Ciardi
B.
Miniati
F.
,
2004
,
MNRAS
 ,
355
,
1053
Dijkstra
M.
Haiman
Z.
Loeb
A.
,
2004
,
ApJ
 ,
613
,
646
Djorgovski
S. G.
Castro
S.
Stern
D.
Mahabal
A. A.
,
2001
,
ApJ
 ,
560
,
L5
Dunkley
J.
et al.,
2009
,
ApJS
 ,
180
,
306
Efron
B.
,
1979
,
Ann. Statist.
 ,
7
,
1
Fan
X.
Narayanan
V. K.
Strauss
M. A.
White
R. L.
Becker
R. H.
Pentericci
L.
Rix
H.
,
2002
,
AJ
 ,
123
,
1247
Fan
X.
et al.,
2006
,
AJ
 ,
132
,
117
Field
G. B.
,
1959
,
ApJ
 ,
129
,
536
Furlanetto
S. R.
Zaldarriaga
M.
Hernquist
L.
,
2004
,
ApJ
 ,
613
,
1
Furlanetto
S. R.
Oh
S. P.
Briggs
F. H.
,
2006
,
Phys. Rep.
 ,
433
,
181
Furlanetto
S. R.
et al.,
2009
, in
astro 2010: The Astronomy and Astrophysics Decadal Survey, Vol. 2010
 . p
82
Gunn
J. E.
Peterson
B. A.
,
1965
,
ApJ
 ,
142
,
1633
Guo
Q.
Wu
X.
Xu
H. G.
Gu
J. H.
,
2009
,
ApJ
 ,
693
,
1000
Harker
G.
et al.,
2009
,
MNRAS
 ,
397
,
1138
Harker
G.
et al.,
2010
,
MNRAS
 ,
405
,
2492
Hobbs
G.
Lyne
A. G.
Kramer
M.
Martin
C. E.
Jordan
C.
,
2004
,
MNRAS
 ,
353
,
1311
Hobson
M. P.
Maisinger
K.
,
2002
,
MNRAS
 ,
334
,
569
Hogan
C. J.
Rees
M. J.
,
1979
,
MNRAS
 ,
188
,
791
Iliev
I. T.
Mellema
G.
Pen
U.
Merz
H.
Shapiro
P. R.
Alvarez
M. A.
,
2006
,
MNRAS
 ,
369
,
1625
Iliev
I. T.
Mellema
G.
Pen
U.-L.
Bond
J. R.
Shapiro
P. R.
,
2008
,
MNRAS
 ,
384
,
863
Jelić
V.
et al.,
2008
,
MNRAS
 ,
389
,
1319
Jelić
V.
Zaroubi
S.
Labropoulos
P.
Bernardi
G.
de Bruyn
A. G.
Koopmans
L. V. E.
,
2010
,
MNRAS
 , p.
1369
Kassim
N. E.
Lazio
T. J. W.
Ray
P. S.
Crane
P. C.
Hicks
B. C.
Stewart
K. P.
Cohen
A. S.
Lane
W. M.
,
2004
,
Planet. Space Sci.
 ,
52
,
1343
Komatsu
E.
et al.,
2011
,
ApJS
 ,
192
,
18
Lidz
A.
Zahn
O.
McQuinn
M.
Zaldarriaga
M.
Hernquist
L.
,
2008
,
ApJ
 ,
680
,
962
Liu
A.
Tegmark
M.
Bowman
J.
Hewitt
J.
Zaldarriaga
M.
,
2009a
,
MNRAS
 ,
398
,
401
Liu
A.
Tegmark
M.
Zaldarriaga
M.
,
2009b
,
MNRAS
 ,
394
,
1575
Loeb
A.
Zaldarriaga
M.
,
2004
,
Phys. Rev. Lett.
 ,
92
,
211301
Lonsdale
C. J.
et al.,
2009
,
IEEE Proc.
 ,
97
,
1497
Madau
P.
Meiksin
A.
Rees
M. J.
,
1997
,
ApJ
 ,
475
,
429
Mao
Y.
Tegmark
M.
McQuinn
M.
Zaldarriaga
M.
Zahn
O.
,
2008
,
Phys. Rev. D
 ,
78
,
023529
Masui
K. W.
McDonald
P.
Pen
U.
,
2010
,
Phys. Rev. D
 ,
81
,
103527
McQuinn
M.
Zahn
O.
Zaldarriaga
M.
Hernquist
L.
Furlanetto
S. R.
,
2006
,
ApJ
 ,
653
,
815
McQuinn
M.
Lidz
A.
Zahn
O.
Dutta
S.
Hernquist
L.
Zaldarriaga
M.
,
2007
,
MNRAS
 ,
377
,
1043
Mezzacappa
A.
,
2005
,
Ann. Rev. Nuclear Part. Sci.
 ,
55
,
467
Morales
M. F.
,
2005
,
ApJ
 ,
619
,
678
Morales
M. F.
Hewitt
J.
,
2004
,
ApJ
 ,
615
,
7
Morales
M. F.
Wyithe
J. S. B.
,
2010
,
ARA&A
 ,
48
,
127
Morales
M. F.
Bowman
J. D.
Cappallo
R.
Hewitt
J. N.
Lonsdale
C. J.
,
2006a
,
New Astron. Rev.
 ,
50
,
173
Morales
M. F.
Bowman
J. D.
Hewitt
J. N.
,
2006b
,
ApJ
 ,
648
,
767
Myers
S. T.
et al.,
2003
,
ApJ
 ,
591
,
575
Pearson
T. J.
Kus
A. J.
,
1978
,
MNRAS
 ,
182
,
273
Parsons
A. R.
et al.,
2010
,
AJ
 ,
139
,
1468
Pen
U.
Chang
T.
Hirata
C. M.
Peterson
J. B.
Roy
J.
Gupta
Y.
Odegova
J.
Sigurdson
K.
,
2009
,
MNRAS
 ,
399
,
181
Peterson
J.
Pen
U.
Wu
X.
,
2004
, preprint (0404083)
Pritchard
J. R.
Loeb
A.
,
2010
,
Phys. Rev. D
 ,
82
,
023006
Purcell
E. M.
Field
G. B.
,
1956
,
ApJ
 ,
124
,
542
Röttgering
H. J. A.
et al.,
2006
, preprint (0610596)
Roy
J.
Gupta
Y.
Pen
U.
Peterson
J. B.
Kudale
S.
Kodilkar
J.
,
2010
,
Exp. Astron.
 ,
28
,
25
Salvaterra
R.
Haardt
F.
Volonteri
M.
,
2007
,
MNRAS
 ,
374
,
761
Scott
D.
Rees
M. J.
,
1990
,
MNRAS
 ,
247
,
510
Shapiro
P. R.
Ahn
K.
Alvarez
M. A.
Iliev
I. T.
Martel
H.
Ryu
D.
,
2006
,
ApJ
 ,
646
,
681
Sunyaev
R. A.
Zeldovich
I. B.
,
1975
,
MNRAS
 ,
171
,
375
Thompson
A. R.
Moran
J. M.
Swenson
G. W.
Jr
,
2001
,
Interferometry and Synthesis in Radio Astronomy, 2nd edn
 .
Wiley
, New York
Trac
H.
Gnedin
N. Y.
,
2009
, preprint (0906.4348)
Tseliakhovich
D.
Hirata
C.
,
2010
,
Phys. Rev. D
 ,
82
,
083520
Wang
X.
Tegmark
M.
Santos
M. G.
Knox
L.
,
2006
,
ApJ
 ,
650
,
529
White
M.
Carlstrom
J. E.
Dragovan
M.
Holzapfel
W. L.
,
1999
,
ApJ
 ,
514
,
12
Willis
A. G.
Oosterbaan
C. E.
de Ruiter
H. R.
,
1976
,
A&AS
 ,
25
,
453
Willott
C. J.
Rawlings
S.
Blundell
K. M.
,
2001
,
MNRAS
 ,
324
,
1
Willott
C. J.
et al.,
2007
,
AJ
 ,
134
,
2435
Wouthuysen
S. A.
,
1952
,
AJ
 ,
57
,
31
Zahn
O.
Lidz
A.
McQuinn
M.
Dutta
S. A.
,
2007
,
ApJ
 ,
654
,
12
Zaldarriaga
M.
Furlanetto
S. R.
Hernquist
L.
,
2004
,
ApJ
 ,
608
,
622
Zaroubi
S.
Silk
J.
,
2005
,
MNRAS
 ,
360
,
L64