{'fraction_non_alphanumeric': 0.02828631594923261, 'fraction_numerical': 0.020375253742103167, 'mean_word_length': 4.326596604688763, 'pattern_counts': {'":': 0, '<': 13, '<?xml version=': 0, '>': 3, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
HI as a Probe of the Large Scale Structure in the Post-Reionization Universe

We model the distribution of neutral Hydrogen (HI hereafter) in the post-reionization universe. This model uses gravity only N-Body simulations and an ansatz to assign HI to dark matter haloes that is consistent with observational constraints and theoretical models. We resolve the smallest haloes that are likely to host HI in the simulations, care is also taken to ensure that any errors due to the finite size of the simulation box are small. We then compute the smoothed one point probability distribution function and the power spectrum of fluctuations in HI. This is compared with other predictions that have been made using different techniques. We highlight the significantly high bias for the HI distribution at small scales. This aspect has not been discussed before. We then discuss the prospects for detection with the MWA, GMRT and the hypothetical MWA5000. The MWA5000 can detect visibility correlations at large angular scales at all redshifts in the post-reionization era. The GMRT can detect visibility correlations at lower redshifts, specifically there is a strong case for a survey at z=1.3. We also discuss prospects for direct detection of rare peaks in the HI distribution using the GMRT. We show that direct detection should be possible with an integration time that is comparable to, or even less than, the time required for a statistical detection. Specifically, it is possible to make a statistical detection of the HI distribution by measuring the visibility correlation, and, direct detection of rare peaks in the HI distribution at z = 1.3 with the GMRT in less than 1000 hours of observations.


INTRODUCTION
Large scale structures in the universe are believed to have formed by gravitational amplification of small perturbations (Peebles 1980;Shandarin & Zeldovich 1989;Peacock 1999;Padmanabhan 2002;Bernardeau et al. 2002). Much of the matter in galaxies and clusters of galaxies is the so called dark matter that is believed to be weakly interacting and non-relativistic (Trimble 1987; Komatsu et al. 2009). Dark matter responds mainly to gravitational forces, and by virtue of larger density than baryonic matter, assembly of matter into haloes and large scale structure is driven by gravitational instability of initial perturbations.
Galaxies are believed to form when gas in highly over-dense haloes cools and collapses to form stars in significant numbers (Hoyle 1953;Rees & Ostriker 1977;Silk 1977;Binney 1977 Bromm & Larson 2004) in turn leads to emission of UV radiation that starts to ionize the inter-galactic medium (IGM). The period of transition of the IGM from a completely neutral to a completely ionized state is known as the epoch of reionization (EoR), e.g., see Loeb & Barkana (2001). Observations indicate that the process of reionization was completed before z ∼ 6 (Fan, Carilli, & Keating 2006;Becker et al. 2001;Fan et al. 2006). Many possible sources of ionizing radiation have been considered, although stellar sources are believed to be the most plausible candidates, see, e.g., Bagla, Kulkarni, & Padmanabhan (2009).
Prior to the EoR, almost all the Hydrogen in the universe is in atomic form. Through the EoR Hydrogen is ionized till we are left with almost no HI in the inter-galactic medium (IGM) and almost all the HI resides in the inter stellar medium (ISM) of galaxies. We focus on the post-reionization era in this work and our aim is to make predictions about the distribution of HI in this regime. It had been proposed by Sunyaev & Zeldovich (1972, 1975) that the hyperfine transition of Hydrogen may be used to probe primordial galaxies. This problem has been approached in past from the perspective of making specific predictions for existing instruments arXiv:0908.3796v3 [astro-ph.CO] 28 Apr 2010 like the Giant Meterwave Radio Telescope (GMRT) 1 (Subramanian & Padmanabhan 1993;Kumar, Padmanabhan, & Subramanian 1995;Bagla, Nath, & Padmanabhan 1997;Bagla 1999;Bharadwaj, Nath, & Sethi 2001;Bagla & White 2003;Bharadwaj & Srikant 2004;Bharadwaj & Ali 2005), or upcoming instruments like the the MWA 2 , ASKAP 3 MeerKAT 4 and SKA 5 , e.g., see Wyithe, Loeb, & Geil (2008). Much work in the last decade has focused on making predictions for the power spectrum of fluctuations in HI, with an implicit assumption that it is easier to make a statistical detection than a direct detection. Further, it has been argued that HI can be used as a tracer of the large scale structure and observations can be used to constrain cosmological parameters with a special emphasis on observations of the baryon acoustic oscillations (BAO) in the power spectrum (Visbal, Loeb, & Wyithe 2008;Chang et al. 2008;Bharadwaj, Sethi, & Saini 2009). Most of the upcoming instruments are sensitive to the power spectrum of fluctuations at large scales, certainly larger than the scale of nonlinearity at the relevant epoch, and hence a significant fraction of the work done in terms of making predictions is based on linear theory or approximations that work well in the linear and quasilinear regime. The halo model has also been used for predicting the spectrum of fluctuations in the post-reionization universe .
In this study we revisit the issue of predicting fluctuations and use high resolution N-Body simulations. This allows us to study fluctuations at scales comparable to, or even smaller than the scale of non-linearity. We use dark matter simulations along with simple ansatz for assigning neutral Hydrogen in order to make predictions of fluctuations in surface brightness temperature. Our model is discussed in detail in §2 where we summarize our knowledge of the HI distribution at high redshifts and motivate the assignment schemes we use in this work . The N-Body simulations used by us are described in §2.1. Results of the HI signal are presented in §3. We then move on to review the relation between the source flux and observables in radio interferometers, i.e. visibilities, in §4 and also discuss the sensitivity of interferometers. Finally we look at the prospects of detection, both statistical and of rare peaks in §5. We conclude with a discussion in §6.

MODELING THE HI DISTRIBUTION
In this section we describe our model of the HI distribution at high redshifts. Our knowledge of the HI distribution in the universe is derived mainly from QSO absorption spectra, where the gas absorbs in the Lyman-α transition of the Hydrogen atom. We know from observations of these absorption spectra that much of the inter-galactic medium (IGM) is highly ionized and does not contain a significant amount of neutral Hydrogen. Most of the neutral Hydrogen resides in relatively rare damped Lyman-α systems (DLAS) (Wolfe, Gawiser, & Prochaska 2005). DLAS and other high column density absorption features are believed to arise due to gas within galaxies (Haehnelt, Steinmetz, & Rauch 2000;Gardner et al. 2001). It is possible to make a quantitative estimate of the total neutral Hydrogen content in DLAS and study the evolution of the total neutral Hydrogen content of the universe (Storrie-Lombardi, McMahon, & Irwin 1996;Rao & Turnshek 2000;Péroux et al. 2005). These observations indicate that at 1 ≤ z ≤ 5, the neutral Hydrogen content of the universe is almost constant with a density parameter of ΩHI 0.001 6 .
At low redshifts, the HI content can be estimated more directly through emission in the Hyperfine transition. Observations in the local universe indicate a much lower neutral Hydrogen content than seen at z ≥ 1 (Zwaan et al. 2005). The neutral gas fraction in galaxies at intermediate redshifts appears to be much higher than in galaxies in the local universe (Lah et al. 2007(Lah et al. , 2009, i.e., the neutral gas fraction in galaxies appears to increase rapidly with redshift z. The spin temperature couples to the gas temperature through collisions of atoms with other atoms, electrons, ions and also the Wouthuysen-Field effect (Purcell & Field 1956;Field 1958Field , 1959Wouthuysen 1952a,b;Zygelman 2005;Furlanetto, Oh, & Briggs 2006;Furlanetto & Furlanetto 2007a,b). Observations of 21 cm absorption by DLAS indicate that the spin temperature is orders of magnitude higher than the temperature of the cosmic microwave background radiation (CMBR) at corresponding redshifts (Chengalur & Kanekar 2000;Kanekar et al. 2009). This implies that the emission in the 21 cm hyperfine transition can be safely assumed to be proportional to the density of neutral Hydrogen, see, e.g., Furlanetto, Oh, & Briggs (2006).
where xHI is the fraction of Hydrogen in neutral form, δ is the density contrast of the gas distribution, T cmb is the temperature of the CMBR Ts is the spin temperature defined using the relative occupation of the two levels for the hyperfine transition: where subscripts 1 and 0 correspond to the excited and ground state levels of the hyperfine transition, T = hν/kB = 68 mK is the temperature corresponding to the transition energy, and (g1/g0) = 3 is the ratio of the spin degeneracy factors of the levels.
Observations indicate that neutral gas is found only in galaxies in the post reionization universe. We also know that at very low redshifts galaxies in groups and clusters do not contain much neutral gas. As cold gas is associated with galaxies, we may assume that it exists only in haloes that are more massive than the Jeans mass. Further, we may assume that neutral gas is predominantly found in galaxies and not in larger haloes that may contain several large galaxies. Jeans mass for haloes in a photo-ionizing UV background depends on the shape of the spectrum of the ionizing radiation, typically we expect gas in haloes with a circular velocity in excess of 60 km/s to cool, fragment and form stars. The mass within the virial radius is related to the circular velocity and the collapse redshift as: Simulation studies show that DLAS can reside in haloes with even lower circular velocities (Pontzen et al. 2008), i.e., haloes with mass lower than the expected Jeans mass can contain significant amount of neutral Hydrogen. The gas in these haloes is able to self shield and maintain a significant amount of HI even though the amount of gas is insufficient for sustaining star formation 7 . We use this input and impose a lower cutoff in circular velocity of 30 km/s, i.e., haloes with a lower circular velocity are not assigned any HI. This lower limit is appropriate for z ∼ 3 and it is likely that the limit shifts to higher values at lower redshifts. We keep the value fixed at 30 km/s in order to make a conservative estimate of clustering of such haloes. Pontzen et al. (2008) also found that haloes much more massive than a few times 10 11 M do not host significant amount of HI. This is consistent with observations in the local universe where galaxies in groups and clusters of galaxies contain very little HI. Some observations suggest that the neutral fraction in galaxies in outer parts of clusters of galaxies at intermediate redshifts may be high compared to corresponding galaxies in the local universe (Lah et al. 2009). Assignment for HI in more massive haloes is done in such a way as to ensure that very massive haloes have a zero or negligible fraction of the total gas in neutral form. In case larger haloes have more HI than we assume, our results under estimate the true signal and may be considered as a lower bound. The transition scale for the higher masses is chosen to coincide with 200 km/s. Given that it is easier to estimate mass than circular velocity for haloes identified in simulations, in particular for haloes with a small number of particles, we choose to translate the threshold in circular velocity to a threshold in mass assuming the collapse redshift and the redshift at which the system is being observed to be the same, i.e., zc = z. While this may lead to inaccuracies regarding inclusion of haloes close to the low mass end, we do not expect it to influence the results in a significant manner. The minimum and the maximum masses for the redshifts considered here are listed in Table 1. Given a halo of mass M , all particles in it are assigned an equal HI mass which is a fraction of their total mass. We describe three kinds of mass assignments here. Here F k (M ) is the mass fraction of HI in a given halo. The fraction of baryons that is in the form of neutral Hydrogen is given by where Ω b is the density parameter of Baryons and Ωnr is the density parameter for non-relativistic matter. For the best fit model for WMAP-5 (Komatsu et al. 2009) the ratio has the value Ωnr/Ω b 1/0.17 6. In all three cases the constant f k is determined by normalizing the HI mass in the simulation volume to ΩHI = 0.001, the density parameter of HI indicated by observations at relevant redshifts (Péroux et al. 2005). It is noteworthy that this normalization requires a significant fraction of gas in haloes in the mass range between Mmin and Mmax be neutral at high redshifts. In the assignment scheme with the sharp cutoff (scheme 1, Eqn. 4), the HI fraction is the same in all haloes with mass between Mmin and Mmax. In the second assignment scheme the HI fraction in a halo decreases monotonically with increasing halo mass, and for large halo masses the HI mass goes to zero. This scheme applies to the physical situation where very large haloes do not have any neutral Hydrogen. In scheme 3, the HI fraction decreases monotonically with increasing halo mass but the HI mass goes to a constant value. The motivation of this scheme is to allow for some neutral Hydrogen surviving in massive haloes in galaxies. In the following section, we compare the resulting distribution of HI with the three assignment schemes described above. The present work is the first one where mass resolution of simulations is adequate for resolving the smallest haloes that may host HI and the simulation volume is also sufficient to make statements about the large scale distribution. We also ensure that the finite size of the simulation volume does not affect the results presented here.

N -Body Simulations
We use gravity only simulations run with the TreePM code (Bagla 2002;Bagla & Ray 2003;). The suite of simulations used here is described in the Table 2. The cosmological model and the power spectrum of fluctuations corresponds to the best fit model for WMAP-5: Ωnr = 0.26, ΩΛ = 0.74, ns = 0.96, We use the Friends-of-Friends (FOF) (Davis et al. 1985) algorithm with a linking length l = 0.2 to identify haloes and construct Table 2. Columns 1 and 2 list the size of the box and the number of particles used in the simulations. Columns 3 and 4 give the mass and force resolution of the simulations, while columns 5 and 6 tell us the redshift at which the simulations were terminated and the redshifts for which the analysis were done a halo catalog. The HI assignment schemes are then used to obtain the distribution of neutral Hydrogen. Several existing and upcoming instruments can probe the postreionization universe using redshifted HI emission. The GMRT can observe redshifted HI emission from a few selected redshift windows whereas most other instruments have continuous coverage in redshift bounded on two sides. We choose to focus on the GMRT windows, as these are representative of the range of redshifts in the post-reionization universe. In particular we will focus on the following redshift windows of GMRT: zout = 5.04, 3.34, 1.33. Finite box effects can lead to significant errors in the distribution of haloes that host galaxies, apart from errors in the abundance of haloes of different masses (Bagla & Prasad 2006;Bagla, Prasad, & Khandai 2009). The choice of simulations used in this work ensures that such effects do not contribute significantly. Previous studies have indicated that at z 0, we need a simulation box with L box ≥ 140 h −1 Mpc for the finite size effects to be negligible (Bagla & Ray 2005;Bagla & Prasad 2006). On the other hand the mass resolution of particles decreases as the cube of simulation volume. We balance the requirements of high mass resolution and a sufficiently large box size by using different simulations for studying the HI distribution at different redshifts. Details of the simulations are given in the table 2.

RESULTS
Given an HI mass assigment for particles in a simulation, we can proceed to compute the expected signal from the HI distribution by making mock radio maps and spectra. We also compute the power spectrum in both the real and the redshift space. For redshift space calculations, we use the peculiar velocity information of particles in haloes. Unlike earlier studies, we resolve haloes of individual galaxies and hence the internal velocity dispersion is naturally accounted for and there is no need to add it by hand (Kumar, Padmanabhan, & Subramanian 1995;Bagla, Nath, & Padmanabhan 1997;Bagla & White 2003).
We use the clouds-in-cell (CIC) smoothing to interpolate densities from particle locations to grid points for the purpose of computing densities and the power spectrum. It is convenient to express the HI power spectrum in terms of the brightness temperature, though in comparisons with the dark matter power spectrum we revert to the usual dimensionless form.
We define the real and redshift space scale dependant bias by Figure 1. Effect of HI mass assignment type on bias. Solid,dashed and dotdashed lines are for HI assignment types 1, 2 and 3 as described in eqs. 4-6.
Bias is computed at z = 3.34 for the 51.2h −1 Mpc box with Mmax = 10 11.5 M . . the ratio of the corresponding dark matter and HI power spectra.
We start by checking the effect of HI assignment scheme on the resulting distribution. We have computed the bias b(k) at z = 3.34 for the L box = 51.2h −1 Mpc simulation using the three schemes introduced above. The results are shown in Figure 1 for the three assignment types described in Eqs. (4-6). We see that at large scales the three assignment schemes give very similar results, whereas there is some disagreement between the assigment scheme Eqn.(4) and the other two. The difference can be attributed to the fact that the scheme described in Eqn.(4) puts more HI mass in haloes with masses near Mmax. The differences between the HI assignment schemes are relatively minor. We choose to work with the scheme described in Eqn.(6) due to a better physical justification, as discussed in §2. As mentioned before, if large mass haloes happen to have a larger HI fraction than we assume then the signal will be larger than what we predict here.
It is noteworthy that while the bias is strongly scale dependent at large k (small scales), it flattens out to a constant value at small k (large scales).
In order to check for the effects of a finite box size, we carried out a test for the 23 h −1 Mpc simulation box. Instead of using the range of halo masses for HI assignment that is appropriate for z = 5.1, we work with a slightly smaller range so that haloes of these masses can be found in the 51.2 h −1 Mpc simulation as well. Figure 2 shows the dark matter power spectrum (Top panel) and the HI power spectrum (Lower panel) with the HI assignment restricted to a smaller range of masses, as described above. We see that the dark matter power spectra from the two simulations agree through the range of scales where there is an overlap. The HI power spec- The agreement in the two curves indicates that the finite box size effects do not lead to errors in the smaller box at the redshift of interest. Bottom: HI power spectrum in the same simulations, see text for details. Show the effect of finite box size on the simulated HI power spectrum, good agreement between the two curves imples that the error is small. tra also agree, though not as well as the dark matter power spectra. These differences are so small that we do not expect these to affect the final results in a significant manner. The difference can be attributed to the fact that clustering of haloes is affected more strongly by the box size effects in simulations (Gelb & Bertschinger 1994;Bagla & Ray 2005).
These tests validate our approach for assignment of HI to haloes, and also show that the effects of a finite box-size are not significant at the level of the power spectrum or the mass function.   Table 1. The dashed line shows b(k) when we use the reference value for Mmax but increase M min to twice the reference value. The dot-dashed curve shows b(k) when Mmax is chosen to be higher: 10 11.9 M , while M min is kept at the reference alue. We see that as the characteristic mass of haloes with HI increases, b(k) increases. All curves are for z = 3.34 and have been computed with the simulation with L box = 51.2 h −1 Mpc.

Bias
This is amongst the first studies of the HI distribution at high redshifts where we resolve the smallest haloes that can host significant amount of HI while ensuring that the finite box size effects do not lead to an erronous distribution of haloes. One of the points that we can address here is the effect of non-linear clustering on the HI distribution and the scale dependance of bias. Some of these effects are illustrated in Figure 3 where we have plotted the power spectrum of fluctuations: the solid line shows the linearly extrapolated power spectrum, the dashed line shows the non-linear dark matter power spectrum and the dot-dashed line shows the HI power spectrum. All power spectra are for z = 3.34. The dark matter and the HI power spectra have been computed with the simulation with It is apparent that at k > 0.5 h Mpc −1 the effects of non-linear clustering significantly enhance the dark matter power spectrum. At k ∼ 10 h Mpc −1 , the enhancement is close to an order of magnitude. This is of utmost interest for upcoming instruments that can resolve small angular scales.
We find that bias b(k) for the HI distribution is much greater than unity at high redshifts. This is to be expected of galaxies at high redshifts (Fry 1996;Mo & White 1996;Bagla 1998a,b;Mo, Mao, & White 1999;Baugh et al. 1999;Magliocchetti et al. 2000;Benson et al. 2000;Roukema & Valls-Gabaud 2000;Sheth, Mo, & Tormen 2001;). We note that the bias is scale dependent, and leads to a larger enhancement in the power spectrum at very small scales.
The value of bias depends strongly on the choice of the characteristic mass of haloes with HI. This is shown in Figure 4. All curves are for z = 3.34 and have been computed with the simulation with L box = 51.2 h −1 Mpc. This figure illustrates that the bias at all scales varies monotonically with the charactersitic mass of haloes with HI. Variation is gentle at large scales but fairly strong at small scales. Therefore it is extremely important to have an accu-rate estimate of the charactersitic masses of such haloes. Indeed, it has been pointed out that observations of the amplitude of clustering in the HI distribution can be used to constrain masses of haloes that host DLAS (Wyithe 2008).
While the preceeding figures describe the statistical bias computed from the ratio of power spectra, Figure 5 shows the stochasticity of bias (Dekel & Lahav 1999) in the HI distribution. This figure shows a scatter plot of δHI smoothed at a scale of 3 h −1 Mpc with a spherical top hat window, plotted as a function of δDM smoothed at the same scale. The figure shows a random subset of points from the simulation with L box = 51.2 h −1 Mpc at z = 3.34. The scatter about the average trend in the δHI −δDM is significant, and increases as we go towards large overdensities in dark matter. The scatter becomes small if we smooth the density distribution at much larger scales.
Similar results concerning the effect of non-linear gravitational clustering and bias are obtained for other redshifts. We concentrate on the evolution of bias and the power spectrum rather than go into a very detailed discussion of the HI distribution at each redshift. We discuss simulated maps in the following subsection. Figure 6 (left panel) shows the evolution of the redshift space power spectrum for the HI distribution. Curves show the ∆ 2 HI (k) as a function of k for z = 5.1 (solid line), z = 3.34 (dashed line) and z = 1.3 (dot-dashed line). Unlike the power spectrum in real space, the enhancement at very small scales is less strong and this is due to velocity dispersion within haloes at smaller scales. Power at larger scales is enhanced due to the Kaiser effect (Kaiser 1987), leading to a relatively stronger enhancement at larger scales, hence the difference in shapes of the HI power spectrum in the real space and the redshift space. Even though the amplitude of density perturbations increases with time, we see that the amplitude of brightness temperature fluctuations decreases instead. This is a manifestation of the decreasing correlation function for galaxies, see, e.g., Bagla (1998b); Roukema & Valls-Gabaud (1998).
The evolution of bias at large scales, denoted as b lin (k) to emphasise that it refers to scales where clustering is still in the linear regime, is shown in the right panel of Figure 6. We see that this decreases from around 2.5 at z = 5.1 to 1.2 at z = 1.3, the variation is slightly steeper than (1 + z) 1/2 and hence the gradual lowering of the amplitude of the brightness temperature power spectrum. Figure 7 shows the simulated radio maps and spectra. The usual conversion from HI density to signal has been used for this (Furlanetto, Oh, & Briggs 2006). The top row shows simulated radio maps, the top-left panel is for z = 3.34 and the top-right panel is for z = 1.3. The pixel size is chosen to be the same as the resolution of the central square of the GMRT. The bandwidth of the map is 0.5MHz and 1.0MHz for z = 3.34, 1.33 respectively. The corresponding contours mark regions with signal (9.2, 5.75, 1.15) µJy for z = 3.34 and (16.8, 10.5, 2.1) µJy for z = 1.33. As we can see in these maps, bright regions have a large angular extent. This fact may be used to enhance prospects of detection §5.

Radio Maps
The lower panel shows spectra from the mock radio maps. These spectra have been plotted for regions that are the size of the resolution of the central square of the GMRT. We have chosen to plot the spectrum through the brightest region in the simulation. The lower-left panel corresponds to z = 3.34 and the lower right panel shows the spectrum for z = 1.3. The peak flux is around 0.025 mJy and the FWHM of the line is close to 300 kHz. The strength of the signal for z = 3.34 is comparable to what we found in an earlier study, although the line width is much smaller. The reason for a smaller line width is likely to be the supression of the HI fraction in very massive haloes. Considering the peak for which the spectrum has been plotted here, it may be possible to make a 2 σ detection with about 2 × 10 3 hours of observations with the central square of the GMRT ( §5). Given that the GMRT observes a much larger volume at these redshifts that the simulation volume, it is highly likely that even rarer peaks in the HI distribution will be observed in a generic pointing (Subramanian & Padmanabhan 1993).
The strentgh of the predicted signal at z = 1.3 is of the order of 0.05 mJy, with an FWHM of around 1 MHz.
Detection of rare peaks in the HI distribution is an exciting possibility. The size of the region represented in these rare peaks is fairly large, and it should be possible to establish the total mass contained in these regions using observations in other wavebands. This will allow us to estimate ΩHI in emission, and hence provide an independent measurement of the cold gas fraction.
The time required for detection of HI at high redshifts with the GMRT increases rapidly as we go to higher redshifts. This trend continues as we move to z 5.1 and therefore we do not discuss those results in detail here.
Instruments like the GMRT, MWA can in principle detect signal from the HI distribution at high redshifts. The angular resolution of the MWA is fairly poor and hence the effects of non-linear clustering and scale dependant bias do not make a significant impact on predictions. We shall discuss the prospects for detection with the GMRT and the MWA in detail in the §5.

SIGNAL AND NOISE IN INTERFEROMETERS
In this section we briefly review the relation between the flux from sources and the observable quantities for radio interferometers. We then proceed to a discussion of the sensitivity of interferometers and the corresponding limitations arising from that.

Visibility Correlation
The quantity measured in radio-interferometric observations is the visibility V (U, ν) which is measured in a number of frequency channels ν − (ν + ∆ν) across a frequency bandwidth B for every pair of antennas in the array. The visibility is related to the sky specific intensity pattern Iν (Θ) as The baseline U is d/λ where d is the antenna separation projected in the plane perpendicular to the line of sight. Θ is a two dimensional vector in the plane of the sky with origin at the center of the field of view, and A(Θ) is the beam pattern of the individual antenna. For the GMRT this can be well approximated by Gaussian A(Θ) = e −θ 2 /θ 0 2 , where θ = |Θ|. We use θ0 = 0.56 • at 610 MHz for the GMRT, and it scales as the inverse of frequency.
We next consider the visibility-visibility correlation (hereafter only visibility correlation) measured at two different baselines U and U + ∆U and at two frequencies ν and ν + ∆ν. As argued in Bharadwaj & Ali (2005) the visibilities at baselines U and U + ∆U will be correlated only if |∆U| < 1/πθ0. Visibilities at different frequencies are expected to be correlated only if the flux from sources is correlated. This is certainly true of continuum radiation. In case of spectral lines, like the redshifted HI 21 cm line being considered here, visibilities are expected to be correlated over a range of frequencies comparable to width of spectral lines. For generic applications, we can take this to be about 1 MHz though we can use the simulated radio maps for a more refined model. As a first approximation, the visibilities at frequencies ν and ν +∆ν can be assumed to be uncorrelated for ∆ν > 1 MHz (Datta, Choudhury, & Bharadwaj 2007) for range of baselines of our interest. We use the power spectra to compute visibility correlations using the relation enunciated by . To calculate the strength and nature of visibility correlation of the HI fluctuations we consider the visibility correlation at same baselines and Figure 7. The top row shows mock radio maps, the top-left panel is for z = 3.34 and the top-right panel is for z = 1.3. The pixel size is chosen to be the resolution of the central square of the GMRT. The bandwidth of the map is 0.5MHz and 1.0MHz for z = 3.34, 1.33 respectively. The corresponding contours mark regions with signal (9.2, 5.75 ,1.15) µJy and (16.8, 10.5, 2.1) µJy for z = 3.34 and 1.33 respectively. The mock maps show the brightest regions from the simulation. As we can see, these regions have a large angular extent this fact may be used to enhance prospects of detection. The lower panel shows spectra from the mock radio maps. These spectra have been plotted for regions that are the size of the resolution of the central square of the GMRT for z = 3.34 and twice the resolution of the central square of GMRT for z = 1.33. We chose to plot the spectrum through the brightest region in the simulation. The lower-left panel corresponds to z = 3.34 and the lower right panel shows the spectrum for z = 1.3. The strength of the signal for z = 3.34 is comparable to what we found in an earlier study. The strentgh of the predicted signal at z = 1.3 is or the order of 0.05 mJy, with an FWHM of around 1 MHz. The best rms sensitivity achieved in this band at the GMRT is 0.02 mJy. This has been achieved with a bandwidth 32 MHz.
frequencies which can be written as (for details see Bharadwaj & Ali (2005) whereĪν is mean HI specific intensity, rν is comoving distance corresponds to frequency ν = 1420/ (1 + z). k = k m + 2πU/rν , where m is unit vector along the line of sight. Here the quantityĪ is proportional to ΩHI , hence the visibility correlation scales as the square of ΩHI .
We next shift our focus on the redshift space HI power spectrum PHI (k, z) (Bharadwaj & Ali 2004). We calculate PHI (k, z) from N-Body simulations and put this in Eqn. 10 to calculate the visibility correlation. As simulations have a finite box size, we Table 3. Characteristics of the GMRT, MWA and MWA 5000 that have been used for calculation of the visibility correlation and noise in the visibility correlation as well as noise in image are listed here. The values of the system parameters for GMRT and MWA were taken from http://gmrt.ncra.tifr.res.in/ and Bowman (2007)   patch it with a linearly extrapolated power spectrum with a constant bias at small k. This is especially appropriate for high redshifts where the field of view of the GMRT is larger than the simulation box. MWA, of course, has a much larger field of view and the statement applies to that as well.
At sufficiently large scales one can write PHI (k, z) as (Kaiser 1987) Here b is the scale independent bias at large scales and β(k) = f /b where f = Ω 0.6 nr is the redshift distortion parameter. P (k, z) is matter power spectrum at redshift z and is dominated by the spatial distribution of dark matter.

System Noise in the Visibility Correlation
The noise rms. in real part in each visibility for single polarization is where Tsys is the total system temperature, kB is the Boltzmann constant, A ef f is the effective collecting area of each antenna, ∆ν is the channel width and ∆t is correlator integration time. This can be written in terms of the antenna sensitivity K = A ef f /2 kB as The noise variance in the visibility correlation is written as (Ali et al. 2008) where Np is the number of visibility pairs in a particular U bin and frequency band B. The HI signal is correlated but the system noise is taken to be uncorrelated in visibilities. We next consider how Np is calculated. We consider an elementary grid with cell of area ∆U 2 in the u − v plane where the signal remains correlated. Total observation time we consider is T . We further define ρ(U, ν) (circularly symmetric) to be the baseline fraction per unit area per in the baseline range U and U + ∆U and is normalised: d 2 U ρ(U, ν) = 1. Number of visibility pairs in a given U , ∆ν bin in the grid is Here N is the total number of antennae Note that the baseline fraction function ρ(U, ν) is different for different interferometric arrays and plays an important role in determining the sensitivity of an array. This function is also frequency dependent. The signal and the system noise is expected to be isotropic and only depends on the magnitude of U . We take average over all cells of area ∆U 2 in the circular annulus between U and U + ∆U . Number of independent elementary bins in a circular annulus between U and U + ∆U is 2πU ∆U/∆U 2 . Moreover there are B/∆ν independent channels for the total bandwidth B. Combining all these we have The noise rms in the visibility correlation can then be written as The bin size ∆U is determined by the field of view of the antenna and for a Gaussian antenna beam pattern A(θ) = e −θ 2 /θ 2 0 , one can show that ∆U = 1/πθ0. The equation 17 gives the system noise in the visibility correlation for any radio experiment. Note that, the system noise rms in the visibility correlation scales as ∼ 1/N (N − 1)T unlike the noise rms in the image which scales as ∼ 1/ N (N − 1)T .
There are two effects that we ignore in the analysis here.
• At small U , comparable with the inverse of the field of view, cosmic variance also limits our ability to measure the visibility correlations. As the field of view of the GMRT is (much) smaller than the MWA, this concern is more relevant in that case. Cosmic variance is subdominant to the system noise as long as we work with scales that are much smaller than the field of view, or work with U that are much larger than the one that corresponds to the field of view.
• For MWA, the field of view is very large and one needs to take the curvature of sky into account (McEwen & Scaife 2008;Ng 2001Ng , 2005. We ignore this as we are not estimating the signal at very large angles. We consider three instruments: the GMRT, MWA and a hypothetical instrument MWA5000. The GMRT has a hybrid distribution with 14 antennae are distributed randomly in the central core of area 1 km ×1 km and 16 antennae are placed along three arms (Y shaped). The antennae are parabolic dishes of diameter 45 m. The MWA will have 500 antennae distributed within a circular region of radius 750 m. The antennae are expected to follow ρant(r) ∝ 1/r 2 distribution with a 20 m flat core. The hypothetical instrument MWA5000 is almost similar to the MWA but with 5000 antennae. We assume antennae are distributed as ρant(r) ∝ 1/r 2 within 1000 m radius and with a flat core of radius 80 m (McQuinn et al. 2006). The normalized baseline distribution ρ(U, ν) for frequency ν = 330 MHz is presented in the Figure 8 for the three instruments (see  for details). Table 3 lists the wavebands and other instrument parameters for which analysis is done in this paper. We assume the effective antenna collecting area (A ef f ) to be same as the antenna physical area. We also list parameters used for calculation of noise in visibility correlations.

Noise in Images
Observed visibilities are used to construct an image of the sky for each frequency channel. The size of the image is of the order of the primary beam θ0 and the resolution depends on the largest baselines used. Both the size of the image and resolution, or the pixel size depend on the frequency of observation in the same manner.
Unlike visibilities, that are uncorrelated outside of the small bins, pixels in an image are not completely uncorrelated. The correlation of pixels in a raw image is similar to the dirty beam. Therefore the analysis of noise in images as a function of scale is a fairly complex problem. Scaling of noise with bandwidth is easy as the noise for different frequency channels is uncorrelated. rms noise for a given pixel can be computed using (assuming two polarizations) (Thompson, Moran, & Swenson 2001) The symbols in the above equation have the same meaning as described in §4.2. We find that for the z = 3.3 window of the GMRT σimage 3.6 µJy for an integration time of 2000 hours and a bandwidth of 1 MHz. This has been computed for the central square of the GMRT with N = 14. The corresponding number of z = 1.3 is σimage 7.8 µJy for an integration time of 400 hours and a bandwidth of 1 MHz.
As a first approximation, we assume that noise in images is uncorrelated and that it scales as the square root of the number of pixels over which signal is smoothed. This allows us to estimate signal to noise ratio for extended structures, but it must be noted that this analysis is approximate.

PROSPECTS FOR DETECTION: SIMULATED MAPS AND SIGNAL
In this section we describe results from simulated HI maps as described in the previous sections.

Visibility Correlation
We use the power spectrum of the HI distribution derived from our N-Body simulations to compute the visibility correlations. For large scales, larger than the simulation size, we patch this with the linearly extrapolated power spectrum with a constant bias. This is the only input required in Eqn. (10) for calculation of visibility correlation. Visibility correlation for various redshifts is shown in Figure 9. We have also plotted the expected system noise in each of the panels. We use the the functional form of ρ(U, ν) from  (also see Figure 8). The expected visibility correlation is shown as a function of U . The expected signal is shown with a solid curve and the expected noise in the visibility correlation is shows by a dashed curve. This has been shown for three instruments: GMRT, MWA and a hypothetical instrument MWA5000. The noise in each panel has been computed assuming 10 3 hours of integration time and visibilities correlated over 1 MHz. We find that MWA5000 should detect visibility correlations for baselines U < 500 at z 5 and at U < 400 at z 3.7. Higher observation time is needed for the MWA to detect the visibility correlations. The detection at smaller U should be possible at high significance level. Whereas the prospects for detection with the GMRT are encouraging for z = 1.3 where it should be possible to make a detection for baselines U < 600.
It is possible to enhance the signal to noise ratio by combining data from nearby bins in U , thus detection may be possible in a shorter time scale as well. The GMRT has been in operation for more than a decade, and its characteristics are well understood. Therefore it is important to make an effort to observe HI fluctuations at z 1.3 with the instrument. The observations need not pertain to the same field as the region observed by the GMRT in a single observation is fairly large and we do not expect significant fluctuations in the HI power spectrum and hence the visibility correlation from one field to another. Thus good quality archival data for a few fields can be combined, in principle, to detect visibility correlations.
Detection of visibility correlations in turn gives us a measurement of the power spectrum for the HI distribution . This can be used to constrain various cosmological parameters (Bharadwaj, Sethi, & Saini 2009;Visbal, Loeb, & Wyithe 2008;. Recently, it has also been pointed out that the power spectrum can also be used to estimate the characteristic mass of damped lyman alpha systems (DLAS) (Wyithe 2008).

Rare Peaks
We now discuss the possibility of directly detecting rare peaks in the HI distribution. This is an interesting possibility for instruments like the GMRT and ASKAP that can resolve small angular scales. At such scales, non-linear gravitational clustering and the large bias for the HI distribution combine to enhance the amplitude of fluctuations significantly ( §3). It is also of interest to see whether this enhancement can make detection of rare peaks as easy as statistical detection via visibility correlations. Figure 10 shows the expected signal to noise ratio (SNR) for Figure 9. This figure shows the expected visibility correlation as a function of U . The expected signal is shown with a solid curve and the expected noise in the visibility correlation is shows by a dashed curve. This has been shown for three instruments: GMRT, MWA and a hypothetical instrument MWA5000. The noise in each panel has been computed assuming 10 3 hours of integration time and that visibilities are correlated over 1 MHz. We note that the MWA5000 should detect visibility correlations at U < 500 at z 5 and at U < 400 at z 3.7. The detection at smaller U should be possible at a high significance level. Prospects of detection with the GMRT are encouraging for z = 1.3 where it should be possible to make a detection for U < 600. It is possible to enhance the signal to noise ratio by combining data from nearby bins in U , thus detection may be possible in a shorter time scale as well. Figure 10. This figure shows the expected signal to noise ratio (SNR) for direct detection of rare peaks in the HI distribution at high redshifts. All the plots pertain to the GMRT and we use the expected noise in the image for the central square. The left column is for z 3.3 and the right column is for z 1.3. The top row is for variation of the SNR with bandwidth for one pixel, whereas the bottom row shows the variation of SNR with the square root of the solid angle over which signal has been smoothed for a given bandwidth. These figures show that a 4 − 5σ detection of a rare peak is possible in 400 hours for z 1.3 and in 2000 hours for z 3.3. the brightest region in the simulated maps if it is observed using the central square of the GMRT. The left column is for z 3.3 and the right column is for z 1.3. The top row shows variation of the SNR with bandwidth for one pixel, whereas the bottom panel shows the variation of SNR with the square root of the solid angle over which signal has been smoothed for a given bandwidth. In plots of variation of SNR with the bandwidth, the peak is very close to the FWHM of spectral lines for large structures (see spectra of rare peaks in §7). Although SNR falls off towards smaller bandwidths, it does not decline significantly towards larger bandwidths and this can be attributed to the presence of nearby lines and is a reflection of clustering at small scales.
These figures show that a 4 − 5σ detection of a rare peak is possible in 400 hours for z 1.3 and in 2000 hours for z 3.3. The time for detection here is comparable to, indeed slightly less than the time required for a statistical detection. This is a very exciting possibility as it allows us to measure the HI fraction of a fairly large region, potentially leading to a measurement of ΩHI . However, unlike statistical detection of the HI distribution, the in-tegration time cannot be divided across different fields for a direct detection.
The volume sampled in observation of one field with the GMRT at z 3.3 is much larger than the volume of the simulation used here. Thus there is a strong possibility of finding even brighter region in a random pointing at that redshift. On the other hand, the volume sampled by the GMRT at z 1.3 is smaller by a factor of more than 5 as compared to the volume of the simulation used for making mock maps. We have studied the distribution function for signal in pixels and find that there are a sufficiently large number of bright pixels in the simulated map and a random pointing should have a rare peak comparable to the one discussed here within the field of view. This is shown in Figure 11 where we have shown the number density of pixels in the simulated radio map with signal above a given threshold. Here, each pixel corresponds to the angular resolution of the central square of the GMRT and a bandwidth of 1 MHz. The volume covered by one GMRT field with a bandwidth of 32 MHz should contain more than ten pixels brighter Figure 11. This figure shows the number density of pixels in the simulated radio map with signal above a given threshold. Here, each pixel corresponds to the angular resolution of the central square of the GMRT and a bandwidth of 1 MHz. than 10 µJy. For comparison, the brightest pixel in the simulated map gives out twice this flux.
We do not discuss the possibility of detection of rare peaks with the MWA as it has a poorer sensitivity at small scales.
Note that in our results we have ignored the foregrounds. The foregrounds at relevant frequencies will dominate over the cosmological HI 21-cm signal and needs to be subtracted. Line of sight modes k 0.065 4.5 1+z 8 MHz B Mpc −1 can not be measured because of foreground subtraction, B is the bandwidth over which foregrounds are subtracted (McQuinn et al. 2006). This would effect the baselines U 80 for a 8 MHz subtraction bandwidth. Rare peaks are much smaller objects (< 10 arcmin) than the scales where the foreground subtraction would affect and hence we expect that our results on rare peaks should not be affected by foregrounds.

DISCUSSION
Several attempts have been made in recent years to model the HI distribution in the post-reionization universe (Scott & Rees 1990;Subramanian & Padmanabhan 1993;Kumar, Padmanabhan, & Subramanian 1995;Bagla, Nath, & Padmanabhan 1997  , with the recent resurgence in interest being due to upcoming radio telescopes and arrays. In this paper we have revisited the issue using dark matter simulations and an ansatz for assigning HI to dark matter haloes. This is amongst the few attempts made using simulations where the smallest haloes that may contain significant amount of HI are resolved, and, the simulations are large enough to limit errors due to a finite size of the simulation box. The prospects for detection of the HI distribution at high red-shifts are very promising. The GMRT can be used to measure visibility correlations at z 1.3. We can expect the hypothetical MWA5000 to measure visibility correlations at z > 3.7 and provide an estimate of the power spectrum. The integration time required is not very large. This is very encouraging as no other existing or upcoming instrument can probe this redshift, though ASKAP and MeerKAT will be able to probe the HI distribution at z < 1. Measurement of visibility correlations, and hence the HI power spectrum can be used to constrain several cosmological parameters. If we are able to make measurements at several redshifts then it becomes possible to constrain models of dark energy. Observations of visibility correlations can also be used to put constraints on the characteristic mass of damped Lyman alpha systems (Wyithe 2008). The quadratic dependence of the visibility correlations on ΩHI can be used to determine the density parameter for HI at high redshifts.
An interesting potential application of the method used here is to make predictions for z < 1 where statistical and individual detection with the GMRT and ASKAP is likely to be much easier. It is obvious that detection of rare peaks of HI with the GMRT may be even easier at lower redshifts (z ≤ 0.4). The main issue to be addressed at these redshifts is that the observed volume becomes small and hence we can expect significant scatter in the observed visibility correlations from one field to another, and across frequency channels. This is true even for rare peaks and one may need to observe several fields before finding a very bright object. Or one may need to observe for longer periods of time in order to observe a not so rare peak. This is particularly relevant for the GMRT where the primary beam covers a solid angle that is nearly 14 times smaller than that covered by ASKAP. We propose to study this issue in a later publication, where we will also discuss prospects of observing HI using ASKAP and MeerKAT. We are also working on using wavelet based methods for detection of rare peaks in the HI distribution.
In this work we have used a fairly simple HI assignment scheme while ignoring the evolution of gas content of haloes. We are working on rectifying this shortcoming and expect to use a semi-analytical model to improve this aspect of modeling.
Key conclusions of this paper may be summarized as follows: • We find that non-linear gravitational clustering enhances the amplitude of perturbations by a significant amount at small scales.
• HI distribution is strongly biased at high redshifts and this enhances the HI power spectrum significantly as compared to the dark matter power spectrum. Bias b(k) is scale independent at large scales k k nl , where k nl is the scale of non-linearity. • Bias decreases sharply from high redshifts towards low redshifts. This leads to a gradual decrease in the brightness temperature power spectrum. The change in the amplitude of the brightness temperature power spectrum is slow, being less than a factor two between z = 5.1 and z = 1.3.
• Brightest regions in the simulated radio maps are found to be extended, with an angular extent much larger than the resolution of the GMRT. This can be used to enhance prospects of detection.
• Spectra in the simulated maps appear to have an FWHM of around 1 MHz for z = 1.3 and about half of this for z = 3.34. A comparison with earlier work for z = 3.34 shows that the FWHM is found to be smaller in the present work. We attribute this to the low HI fraction assigned to the most massive haloes, whereas such haloes were not excluded in earlier work.
• The MWA5000 should detect visibility correlations at U < 300 at z 5 and at U < 700 at z 3.7. The detection at smaller U should be possible at a high significance level.
• Prospects of detection with the GMRT are encouraging for z = 1.3 where it should be possible to make a detection for U < 600.
• Signal for the brightest source of redshifted 21 cm radiation from z = 1.3 appears to be significant, and within reach of an instrument like the GMRT.
• Detection of rare peaks in the HI distribution offers exciting possibilities since it does not require a much longer integration time as compared to statistical detection. The size of the region represented in these rare peaks is fairly large, and it should possible to establish the mass contained in these regions using observations in other wavebands. The direct detection at z 1.3 requires only a few hundred hours of observing time with an existing instrument. This may allow us to estimate ΩHI in emission, and hence provide and independent measurement of the amount of cold gas.