PROBING THE PRE-REIONIZATION EPOCH WITH MOLECULAR HYDROGEN INTENSITY MAPPING

Existing cosmological observations show that the reionization history of the universe at z > 6 is likely both complex and inhomogeneous (e.g., Haiman 2004; Choudhury & Ferrara 2006; Zaroubi 2012). While the polarization signal in the cosmic microwave background (CMB) anisotropy power spectrum constrains the total optical depth to electron scattering and the existing Wilkinson Microwave Anisotropy Probe measurements suggest reionization happened around z_ri = 11 (Komatsu et al. 2011), it is more likely that the reionization period was extended over a broad range of redshifts from 20 to 6. Moving beyond CMB, observations of the 21 cm spin-flip line of neutral hydrogen are now pursued to study the spatial distribution of the matter content during the epoch of reionization (e.g., Madau et al. 1997; Loeb & Zaldarriaga 2004; Gnedin & Shaver 2004). Unlike CMB, 21 cm data are useful as they provide a tomographic view of the reionization (Furlanetto et al. 2004; Santos et al. 2005). The anisotropy power spectrum of the 21 cm line emission is also a useful cosmological probe (Santos & Cooray 2006; McQuinn et al. 2006; Bowman et al. 2007; Mao et al. 2008).

While the 21 cm signal primarily traces the neutral hydrogen content in the intergalactic medium during reionization, line emission associated with atomic and molecular lines is of interest to study the physical properties within dark matter halos, such as gas cooling, star formation, and the spatial distribution of the first stars and galaxies. Motivated by various experimental possibilities, we have studied the reionization signal associated with the CO (Gong et al. 2011), C ii (Gong et al. 2012), and Lyα (Silva et al. 2013) lines. As the signal is sensitive to the metal abundance, these atomic and molecular probes are more sensitive to the late stages of reionization, perhaps well into the epoch when the universe is close to full reionization and has a low 21 cm signal (Basu et al. 2004; Righi et al. 2008; Visbal & Loeb 2010; Carilli 2011; Lidz et al. 2011).

Although the end of the reionization era can be effectively probed with H i, CO, C ii, and Lyman-α, it would also be useful to have an additional probe of the onset of reionization at z > 10. Here we consider molecular hydrogen and study the signal associated with rotational and vibrational lines in the mid-IR wavelengths. Molecular hydrogen has been invoked as a significant coolant of primordial gas leading to the formation of the first stars and galaxies (e.g., Haiman 1999; Bromm & Larson 2004; Glover 2005; Glover 2012). While molecular hydrogen is easily destroyed in the later stages of reionization, its presence in the earliest epochs of the cosmological history can be probed with line emission experiments.

This paper is organized as follows: in the next section, we outline the calculation related to the cooling rate of H₂ rotational and vibrational lines. In Section 3, we present results on the H₂ luminosity as a function of the halo mass, while in Section 4, we discuss the mean H₂ intensity and clustering auto and cross-power spectra. The cross-power spectra between various lines are proposed as a way to eliminate the low-redshift contamination and increase the overall signal-to-noise ratio (S/N) for detection. We discuss potential detectability in Section 5. We summarize our results and conclude in Section 6. We assume the flat ΛCDM with Ω_M = 0.27, Ω_b = 0.046, σ₈ = 0.81, n_s = 0.96, and h = 0.71 for the calculation throughout the paper (Komatsu et al. 2011).

The radiation emitted by H₂ will be generated by the heating/cooling of the gas as the collapsing process evolves. Therefore, in order to calculate the H₂ luminosity, we first evaluate the cooling rate of the H₂ rotational and vibrational lines for optical-thin and optical-thick media, respectively. For hydrogen density n_H < 10⁹ cm⁻³, the optical depth is thin for H₂ emission lines. Following Hollenbach & McKee (1979), the cooling coefficient for the rotational and vibrational lines can then be expressed as

$$\begin{equation} \Lambda _{\rm H_2}^{\rm r,v}({\rm H,H_2}) = \frac{\Lambda _{\rm LTE}^{\rm r, v}({\rm H,H_2})}{1+n_{\rm cr}^{\rm r,v}({\rm H,H_2})/n_{\rm H,H_2}}, \end{equation} \tag{ 1 }$$

where $\Lambda _{\rm LTE}^{\rm r, v}({\rm H,H_2})$ is the cooling coefficients of rotation or vibration at the local thermodynamic equilibrium (LTE) for the collisions with hydrogen atoms H or H₂. Here, $n_{\rm cr}^{\rm r,v}$ is the critical density of H or H₂ to reach the LTE, and $n_{\rm H,H_2}$ is the local number density of H or H₂.

In the LTE, we have A_ul = C_uln, where A_ul is the Einstein coefficient, C_ul is the collisional de-excitation rate from upper to lower level, and n is the particle number density (Hollenbach & McKee 1979). Then the rotational and vibrational LTE cooling coefficient can be written as

$$\begin{equation} \Lambda _{\rm LTE}^{\rm r,v}({\rm H,H_2}) = \frac{1}{n_{\rm H,H_2}}A_{J}\frac{g_J}{g_{J^{\prime }}}e^{-\frac{\Delta E}{kT}}\Delta E, \end{equation} \tag{ 2 }$$

where g_J = 2J + 1 is the statistical weight, J denotes the total angular momentum quantum number of the rotational energy level, A_J is the Einstein coefficient for the J → J' transition at the same vibrational energy level or between two vibrational energy levels, and ΔE is the energy difference between E_J and $E_{J^{\prime }}$. The wavenumber for each energy level and the calculation for E_J and wavelength are given in the Appendix. The values of A_J are taken from Turner et al. (1977). Note that we only consider two vibrational level transitions v = 0 and 1 in the following calculation. The number of H₂ at the higher vibrational levels is much smaller than that at v = 0 and 1 in our case, so the strength of these lines is much weaker than those from v = 0 and 1. The n_cr/n term in Equation (1) can be approximated by (Hollenbach & McKee 1979)

$$\begin{equation} \frac{n_{\rm cr}^{\rm r,v}({\rm H,H_2})}{n_{\rm H,H_2}} = \frac{\Lambda _{\rm LTE}^{\rm r,v}({\rm H,H_2})}{\Lambda _{\rm n\rightarrow 0}^{\rm r,v}({\rm H,H_2})}, \end{equation} \tag{ 3 }$$

where $\Lambda _{\rm n\rightarrow 0}^{\rm r,v}$ is the low-density limit of the cooling coefficient, which can be obtained by replacing the Einstein coefficient in Equation (2) by $C_J^{\rm H,H_2}n_{\rm H,H_2}$, i.e.,

$$\begin{equation} \Lambda _{\rm n\rightarrow 0}^{\rm r,v}({\rm H,H_2}) = C_J^{\rm H,H_2}\frac{g_J}{g_{J^{\prime }}}e^{-\frac{\Delta E}{kT}}\Delta E. \end{equation} \tag{ 4 }$$

Here $C_J^{\rm H,H_2}$ is the collisional de-excitation coefficients with H or H₂ for the J → J' transition in the same vibrational level or between v = i and v = j, which are estimated by the fitting formula given in Hollenbach & McKee (1979) and Hollenbach & McKee (1989; see the Appendix). Using Equations (1)–(4), we can estimate the cooling coefficient for a given H₂ line. Note that in Hollenbach & McKee (1979), the fitting formulae of total cooling coefficients for v = 0 rotational and v = 0, 1, and 2 vibrational lines are given. We denote those by $\Lambda _{\rm H_2}^{\rm tot, rv}$, $\Lambda _{\rm LTE}^{\rm tot,rv}$, and $\Lambda _{n \rightarrow 0}^{\rm tot,rv}$ in the cases of the total local, LTE, and low-density limit cooling coefficients, respectively. Then n_cr/n can be expressed by $\Lambda _{\rm LTE}^{\rm tot,rv}/\Lambda _{n \rightarrow 0}^{\rm tot,rv}$, where n_cr denotes the critical density when all energy level transitions are in the LTE.

For the hydrogen density n_H > 10⁹ cm⁻³, the optical depth is thick for the H₂ emission, and we have to consider the absorption effect for the H₂ cooling. Following Yoshida et al. (2006), we make use of the cooling efficiency, which is defined by f_ce = Λ_thick/Λ_thin, to evaluate the cooling coefficient Λ_thick for the optical-thick case. This reduction factor is derived from their simulations and is available for n_H ≲ 10¹⁴ cm⁻³, which is well within the density ranges of our calculation. The f_ce is about 0.02 when n_H ∼ 10¹⁴ cm⁻³ and increases to about 0.1 when n_H ∼ 10¹² cm⁻³. At densities below 10¹⁰ cm⁻³, we have f_ce = 1 (see Figure 4 of Yoshida et al. 2006).

In Figure 1, as an example, we show the cooling coefficients for H₂–H and H₂–H₂ collisions as a function of temperature T. We assume the number density of hydrogen atom H and molecular hydrogen H₂ to be 10⁵ cm⁻³ here. The dashed curves show the cooling coefficients of the rotational lines, which are in red (0–0S(0)), magenta (0–0S(1)), green (0–0S(2)), and blue (0–0S(3)). The solid curves are for two vibrational lines 1–0S(1) in red and 1–0Q(1) in blue, respectively. As can be seen, the rotational cooling dominates at the low temperature (T ≲ 10³ K) and the vibrational cooling dominates at high temperature (T ≳ 10³ K). Also, we find the $\Lambda _{\rm H_2}^{\rm r}$ for H₂–H₂ collision is generally greater than that for H₂–H collision at T ≲ 10³ K, while the $\Lambda _{\rm H_2}^{\rm v}$ for H₂–H₂ collision is less than that for H₂–H collision in this temperature range. This indicates that at low temperature, the total $\Lambda _{\rm H_2}^{\rm r}$ is mainly from H₂–H₂ collisions, and the total $\Lambda _{\rm H_2}^{\rm v}$ is from H₂–H collisions. At higher temperature with T ≳ 10³ K, the cooling rates $\Lambda _{\rm H_2}^{\rm r}$ and $\Lambda _{\rm H_2}^{\rm v}$ for both H₂–H₂ and H₂–H collisions are similar.

Zoom In Zoom Out Reset image size

We now explore the H₂ luminosity as a function of halo mass for rotational and vibrational lines. As discussed in the last section, the H₂ cooling coefficient $\Lambda _{\rm H_2}$ is dependent on the local gas temperature and density of hydrogen and molecular hydrogen. To evaluate H₂ luminosity versus halo mass relation for molecular hydrogen cooling within primordial dark matter halos, we first need to know the radial profile of the gas temperature and density within dark matter halos.

Following the results from numerical simulations involving the formation of primordial molecular clouds (e.g., Omukai & Nishi 1998; Abel et al. 2000; Omukai 2001; Yoshida et al. 2006; McGreer & Bryan 2008), we assume the gas density profile as

$$\begin{equation} \rho (r) = \rho _0\left(\frac{r}{r_0}\right)^{-2.2}, \end{equation} \tag{ 5 }$$

where we set r₀ = 1 pc and ρ₀ is the normalization factor which is obtained by

$$\begin{equation} M_{\rm gas}=4\pi \int _0^{r_{\rm vir}} r^2 \rho (r) dr{.} \end{equation} \tag{ 6 }$$

Here, M_gas = (Ω_b/Ω_M)M is the gas mass in the virial radius of the halo with dark matter mass M, and r_vir is the virial radius which is given by

$$\begin{equation} r_{\rm vir} = \left[ \frac{M}{(4/3)\pi \rho _{\rm vir}} \right]^{1/3}{.} \end{equation} \tag{ 7 }$$

Here, ρ_vir(z) = Δ_c(z)ρ_cr(z) is the virial density, ρ_cr(z) = 3H²(z)/(8πG) is the critical density at z, H(z) is the Hubble parameter, and Δ_c(z) = 18π²  +  82x  −  39x² where x = Ω_M(z)  −  1.

We then derive the number density of gas by n(r) = n_H(r) + n_He(r). Here, n_H(r) = f_Hρ(r)/m_H is the number density of hydrogen, where f_H = 0.739 is the hydrogen mass fraction and m_H is the mass of hydrogen atom. Similarly, n_He(r) = (1 − f_H)ρ(r)/m_He is the number density of helium, where m_He is the mass of helium atom. Also, the temperature–density relation T(n) and the H₂ fraction–density relation $f_{\rm H_2}(n)=n_{\rm H_2}/n$ can be derived from existing numerical simulations. Here, we use the results on T(n) and $f_{\rm H_2}(n)$ from Omukai (2001) and Yoshida et al. (2006), which are available for n ≃ 10⁻²–10²³ cm⁻³ as shown in the left panel of Figure 2.³ The uncertainties of the gas temperature and H₂ fraction are shown in blue regions. These uncertainties are evaluated based on the differences in the far-ultraviolet radiation background from the first stars and quasars (Omukai 2001).

Zoom In Zoom Out Reset image size

As can be seen, the gas temperature does not monotonously increase with the gas density. For instance, it drops from T ∼ 2000 to 200 K between n ≃ 1 and 10⁴ cm⁻³ where molecular hydrogen density is rising to $f_{\rm H_2}\,{\sim}\, 10^{-3}$. This indicates that H₂ cooling is starting to become important in this gas density range. At n ≃ 10⁴ cm⁻³, H₂ cooling saturates and turns into the cooling at the LTE. For n = 10¹⁰–10¹¹ cm⁻³, almost all of the gas particles become molecular hydrogen due to the efficient H₂ three-body reaction (Yoshida et al. 2006), and we find $f_{\rm H_2}\simeq 0.5$ by definition. At n ≃ 10²⁰ cm⁻³ with T ≃ 10⁴ K, H₂ begins to dissociate and the fraction drops quickly to $f_{\rm H_2}<10^{-5}$ when n ≃ 10²³ cm⁻³ and T ≃ 10⁵ K.

Next, with the help of Equation (5), we can evaluate the gas temperature and H₂ fraction as a function of the halo radius, i.e., T(r) and $f_{\rm H_2}(r)$. Once these are established, we can derive $n_{\rm H_2}(r)$, n_H(r), and $\Lambda _{\rm H_2}^{\rm r,v}(r)$ which are needed for the H₂ luminosity calculation. In the right panel of Figure 2, we show the density profile of the gas and molecular hydrogen in blue solid and dashed lines for a dark matter halo with M = 10⁶ M_☉ h⁻¹ at z = 15. The blue region shows the uncertainty of the H₂ density profile which is derived by the uncertainty of $f_{\rm H_2}$ in the left panel of Figure 2. The gas density profile is a straight line with a slope of −2.2 as indicated by Equation (5). On the other hand, the density profile of molecular hydrogen has a more complex shape which is dependent on the relation between $f_{\rm H_2}$ and gas density n. For the outer layer of the gas cloud (r > 10⁻² pc), gas density is less than 10⁵ cm⁻³ and $f_{\rm H_2}\lesssim 4\times 10^{-4}$. Here, H₂ density is much smaller than the gas density. For the inner region with 10^−3.5 < r < 10⁻² pc, we find 10⁵ < n < 10¹¹ cm⁻³ and $f_{\rm H_2}$ begins to rise up quickly with $n_{\rm H_2}$ becoming close to the gas density. For the innermost region at r < 10^−3.5 pc, the gas density is greater than 10¹¹ cm⁻³, and $f_{\rm H_2}\simeq 0.5$ so that almost all of the hydrogen ends up forming molecular hydrogen.

The luminosity of the H₂ rotational or vibrational lines can then be estimated by

$$\begin{eqnarray} L_{\rm H_2}^{\rm r,v}(M,z) &=& 4\pi \int _0^{r_{\rm vir}}dr r^2 n_{\rm H_2}^{\rm r,v}(r) \nonumber \\ &&\times \left[n_{\rm H}(r)\Lambda _{\rm H_2}^{\rm r,v}(\rm H)+n_{\rm H_2}(r)\Lambda _{\rm H_2}^{\rm r,v}(\rm H_2)\right], \end{eqnarray} \tag{ 8 }$$

where $n_{\rm H_2}^{\rm r,v}(r)$ is the number density of the molecular hydrogen that can emit at a given rotational or vibrational line at r. We first evaluate the total $n_{\rm H_2}$ at v = 0 and 1 states by condensing all the rotational levels at a given vibrational state to be a single vibrational level, n_i = n_{i − 1} exp[ − ΔE_{i, i − 1}/(kT)] where i = 1. Here, g₀ = g₁ = 1 for v = 0 and 1, respectively, and ΔE₁₀/k ≃ 5860 K (Hollenbach & McKee 1979). Then we estimate $n_{\rm H_2}$ for a given rotational energy level J in a vibrational level i by ${n_J} = n_{J^{\prime }}\ ({g_J}/g_{J^{\prime }})\ {\rm exp}[-\Delta E_{J,J^{\prime }}/(kT)]$. The fractions of the ortho and para states of total $n_{\rm H_2}$ are assumed to be 0.75 and 0.25, respectively, in our calculation.

In Figure 3, we show the H₂ luminosity of several lines as a function of the halo mass at z = 15. In these lines, we find that the rotational line 0–0S(3) at a rest-frame wavelength of 9.66 μm is the most luminous one. Other lines, such as 0–0S(5), 1–0S(1), 1–0Q(1), and 1–0O(3), are also strong for halos with high mass (see also Table 1). As can be seen, for low halo masses with M ≲ 10⁵ M_☉ h⁻¹, the rotational lines are stronger than the vibrational lines. This is caused by the fact that the halos with low masses have lower mean gas temperature than the massive halos, and the rotational cooling is stronger than the vibrational cooling in such halos, as indicated by Figure 1.

Zoom In Zoom Out Reset image size

Table 1. The Wavelength, ΔJ = J − J', Spontaneous Emission Coefficient A_J, Mean Bias, and Mean Intensity for the H₂ Rotational and Vibrational Lines at z = 15

H2 Line	λ	ΔJ	AJ(s−1)	$\bar{b}_{\rm H_2}$	$\bar{I}_{\rm H_2}$ (Jy sr−1)
	(μm)
0–0S(0)	28.2	+2	2.94 × 10−11	$2.6^{+0.4}_{-0.1}$	$0.08^{+0.28}_{-0.06}$
0–0S(1)	17.0	+2	4.76 × 10−10	$2.8^{+0.3}_{-0.3}$	$1.52^{+4.87}_{-0.83}$
0–0S(2)	12.3	+2	2.76 × 10−9	$3.0^{+0.2}_{-0.3}$	$1.32^{+2.20}_{-0.61}$
0–0S(3)	9.66	+2	9.84 × 10−9	$3.1^{+0.2}_{-0.2}$	$5.90^{+3.60}_{-2.64}$
0–0S(4)	8.03	+2	2.64 × 10−8	$3.2^{+0.2}_{-0.1}$	$1.97^{+0.94}_{-0.92}$
0–0S(5)	6.91	+2	5.88 × 10−8	$3.3^{+0.2}_{-0.2}$	$4.26^{+2.15}_{-2.1}$
0–0S(6)	6.11	+2	1.14 × 10−7	$3.4^{+0.2}_{-0.2}$	$0.78^{+0.52}_{-0.40}$
0–0S(7)	5.51	+2	2.00 × 10−7	$3.5^{+0.2}_{-0.3}$	$1.05^{+1.02}_{-0.54}$
0–0S(8)	5.05	+2	3.24 × 10−7	$3.6^{+0.3}_{-0.6}$	$0.13^{+0.22}_{-0.07}$
0–0S(9)	4.69	+2	4.90 × 10−7	$3.8^{+0.2}_{-0.9}$	$0.13^{+0.44}_{-0.07}$
0–0S(10)	4.41	+2	7.03 × 10−7	$4.0^{+0.2}_{-1.3}$	$0.01^{+0.11}_{-0.01}$
0–0S(11)	4.18	+2	9.64 × 10−7	$4.2^{+0.2}_{-1.6}$	$0.01^{+0.23}_{-0.01}$
1–0S(0)	2.22	+2	2.53 × 10−7	$3.4^{+0.2}_{-0.7}$	$0.24^{+0.87}_{-0.12}$
1–0S(1)	2.12	+2	3.47 × 10−7	$3.5^{+0.2}_{-0.7}$	$0.83^{+3.06}_{-0.42}$
1–0Q(1)	2.41	0	4.29 × 10−7	$3.4^{+0.2}_{-0.5}$	$1.00^{+1.96}_{-0.50}$
1–0O(3)	2.80	−2	4.23 × 10−7	$3.4^{+0.2}_{-0.5}$	$0.99^{+1.96}_{-0.50}$

Download table as:  ASCIITypeset image

Given the relation between H₂ luminosity and the dark matter halo mass, the mean intensity of the H₂ lines can be expressed as (Visbal & Loeb 2010; Gong et al. 2011)

$$\begin{equation} \bar{I}_{\rm H_2}(z) = \int _{M_{\rm min}}^{\infty }dM \frac{dn}{dM}(M,z)\frac{L_{\rm H_2}(M,z)}{4\pi D_{\rm L}^2}y(z)D_{\rm A}^2 \,, \end{equation} \tag{ 9 }$$

where we choose M_min = 10 M_☉ h⁻¹, dn/dM is the halo mass function (Sheth & Tormen 1999), $y(z)=d\chi /d\nu =\lambda _{\rm H_2}(1+z)^2/H(z)$ when χ is the comoving distance, and $\lambda _{\rm H_2}$ is the wavelength of H₂ lines in the rest frame. Our results are not strongly sensitive to the exact value of minimum halo mass. If we increase the minimum halo mass to the level of 10⁶ M_☉ h⁻¹, then the mean intensity we present here decrease by a factor of ∼2 for all H₂ lines.

In Figure 4, we show the mean intensity of the eight strongest H₂ lines as a function of redshift z. The uncertainty in the intensity of the 0–0S(3) line is shown with the shaded blue region, which is derived from uncertainties in the gas temperature and $f_{\rm H_2}$ in the left panel of Figure 2. We find that the mean intensity of the 0–0S(3) rotational line is the strongest for 10 ⩽ z ⩽ 30. This is because the 0–0S(3) line is the most luminous line for low-mass halos, which have a higher number density and dominate the halo distribution at z = 15. Also, the slopes of the intensity–redshift relation for the vibrational lines are steeper than that of the rotational lines, since they have steeper slopes for cooling coefficients with temperature as shown in Figure 1. However, we find the difference in slopes to become smaller for the rotational lines as J increases, indicating that the high-J rotational lines have similar slopes with cooling coefficient when compared to that of the vibrational lines.

Zoom In Zoom Out Reset image size

We note here that we do not consider the dissociation effect of the molecular hydrogen by Population II and Population III stars in our calculation. We expect the formation of these stars to be important at z  <  10 and that there would be a significant amount of H₂ that should be dissociated by the UV photons emitting from the first stars. Thus, H₂ emission could be suppressed significantly at z  ⩽  10. At z  ∼  15, there should still be some dissociation but we ignore it to obtain a safe upper limit estimate on the expected H₂ intensity for experimental planning purposes.

Next, we can derive the clustering power spectrum of the H₂ lines, writing the intensity as $I_{\rm H_2}(z)=\bar{I}_{\rm H_2}[1+b_{\rm H_2}\delta ({\rm \bf x})]$. Here, $b_{\rm H_2}$ is the average H₂ clustering bias, which can be estimated from

$$\begin{equation} \bar{b}_{\rm H_2}(z) = \frac{\int _{M_{\rm min}}^{\infty }dM \frac{dn}{dM} L_{\rm H_2}b(M,z)}{\int _{M_{\rm min}}^{\infty }dM \frac{dn}{dM} L_{\rm H_2}} \,, \end{equation} \tag{ 10 }$$

where b(M, z) is the bias factor for dark matter halos with mass M at z (Sheth & Tormen 1999). The H₂ clustering auto power spectrum is then given by

$$\begin{equation} P_{\rm H_2}^{\rm clus}(k,z) = \bar{I}_{\rm H_2}^2 \bar{b}_{\rm H_2}^2 P_{\delta \delta }(k,z), \end{equation} \tag{ 11 }$$

where P_δδ(k, z) is the matter power spectrum, which is obtained from a halo model (Cooray & Sheth 2002). At a high redshift as z = 15, the structure of matter distribution is extremely linear and the two-halo term dominates the power spectrum.

We can also estimate the shot-noise power spectrum for the H₂ lines, which is caused by the discretization of the spacial distribution of the primordial clouds,

$$\begin{equation} P_{\rm H_2}^{\rm shot}(z) = \int _{M_{\rm min}}^{\infty } dM \frac{dn}{dM} \left[ \frac{L_{\rm H_2}}{4\pi D_{\rm L}^2} y(z) D_{\rm A}^2 \right]^2 \,. \end{equation} \tag{ 12 }$$

In Table 1, we tabulate the rest-frame wavelength, ΔJ = J − J', spontaneous emission coefficient A_J, mean bias, and mean intensity for 13 rotational and 4 vibrational lines at z = 15. The uncertainties of the mean bias and intensity are evaluated by the uncertainty of the gas temperature and $f_{\rm H_2}$ from the simulations. We find that the mean intensity of the 0–0S(3) rotational line at a rest wavelength of 9.66 μm is the strongest among these lines, with a value of around 6 Jy sr⁻¹ and a range from 3 to 10 Jy sr⁻¹. The other rotational lines such as 0–0S(5), 0–0S(4), 0–0S(1), and 0–0S(2) are also bright with total intensities of ∼4.3, 2.0, 1.5, and 1.3 Jy sr⁻¹, respectively, at z = 15. The vibrational lines 1–0S(1), 1–0Q(1), and 1–0O(3) have low mean intensities at the level of 0.83, 1.0, and 0.99 Jy sr⁻¹, respectively. The mean bias factor of these lines lies between 2.6 (for 0–0S(0)) and 4.2 (for 0–0S(11)), and the bias factors of the rotational lines at higher rotational energy level are higher than that at lower level. This is because the lines with high J are stronger at higher mass halos where the temperature is larger.

In the left panel of Figure 5, the clustering auto power spectra of eight H₂ lines at z = 15 are shown. We find that the shot-noise power spectrum P_shot is relatively small compared to the clustering power spectrum P_clus, and would not affect P_clus at the scales of interest. This is easy to understand if we notice that the halo mass function is dominated by halos with low masses which are more abundant.

Zoom In Zoom Out Reset image size

We also calculate the cross-correlation between two different H₂ lines. Such a cross-correlation will reduce the astrophysical contamination from the other sources, such as low-redshift emission lines from star-forming galaxies, including 63 μm [O i] and 122 μm [N ii], among others. At z  ∼  15, the dominant rotational line 0–0S(3) would be observed at a wavelength of 155 μm. Such a line would be contaminated by, for example, z ∼ 0.3 galaxies emitting [N ii]. Thus the auto power spectrum would be higher than what we have predicted given that the line intensities of [N ii] are higher than the H₂ lines. To avoid this astrophysical line confusion, we propose a cross-correlation between two rotational or rotational and vibrational lines of the H₂ line emission spectrum.

The cross-clustering and shot-noise power spectrum for two such H₂ lines i and j can be evaluated as

$$\begin{equation} P_{\rm clus}^{ij}=\bar{I}_{\rm H_2}^i \bar{I}_{\rm H_2}^j \bar{b}_{\rm H_2}^i \bar{b}_{\rm H_2}^j P_{\delta \delta } \end{equation} \tag{ 13 }$$

and

$$\begin{equation} P_{\rm shot}^{ij} = \int _{M_{\rm min}}^{\infty } dM \frac{dn}{dM} \frac{L_{\rm H_2}^i}{4\pi D_{\rm L}^2} \frac{L_{\rm H_2}^j}{4\pi D_{\rm L}^2}y^i(z) D_{\rm A}^2 y^j(z) D_{\rm A}^2 \,, \end{equation} \tag{ 14 }$$

respectively. From these equations, we find that the cross-power spectrum should have a similar magnitude to the auto power spectrum. The clustering cross-power spectra $P_{ij}^{\rm clus}$ for several H₂ lines at z = 15 are shown in the right panel of Figure 5. We choose the strongest 0–0S(3) line to cross-correlate with the other five bright lines, i.e., 0–0S(1), 0–0S(2), 0–0S(4), 0–0S(5), and 1–0Q(1). We find that 0–0S(3) × 0–0S(5) is the largest cross-power spectrum since they are brightest two lines. At z ∼ 15, then we would be cross-correlating the wavelength regimes around 110 and 155 μm. A search for mid-IR lines revealed no astrophysical confusions from low redshifts that overlap in these two wavelengths at the same redshift. Thus, while low-redshift lines will easily dominate the auto power spectra of H₂ lines, the cross-power spectrum will be independent of the low-redshift confusions. In addition to reducing the astrophysical confusions, the cross-power spectra also have the advantage that they can minimize instrumental systematics and noise, depending on the exact design of an experiment.

In this section, we investigate the possibility of detecting these lines based on current or future instruments. We assume a SPICA-like⁴ survey with 3.5 m aperture diameter, 0.1 deg² survey area, 10 GHz band width, R = 700 frequency resolution, 100 spectrometers, and 250 hr total integration time and noise per detector $\sigma _{\rm pix}=10^6\,{\rm Jy\,\sqrt{s}\,sr^{-1}}$ at 100 μm. Such an instrument corresponds to the latest design of the mid-IR spectrometer, BLISS, from SPICA (M. Bradford 2010, private communication).

In Figure 6, we show the errors of the auto power spectrum of the 0–0S(3) line and the cross-power spectrum of 0–0S(3) × 0–0S(5) at z = 12 for two cases, a SPICA/BLISS-like and an experiment with 10× better sensitivity than with the current design of SPICA/BLISS with $\sigma _{\rm pix}=10^5\,{\rm Jy\,\sqrt{s}\,sr^{-1}}$. The noise power spectrum from the instrument and shot-noise power spectrum caused by the discrete distribution of the gas clouds are also shown in long-dashed and dotted lines, respectively. We estimate the noise power spectrum and the errors by the same method described in Gong et al. (2012). We find that the S/N is 0.2 and 5.2 for the auto power spectrum of the 0–0S(3) line in the two cases, and S/N = 0.1 and 4.5 for the cross-power spectrum of 0–0S(3) × 0–0S(5). This indicates that the current version of SPICA/BLISS does not have the sensitivity to measure the intensity fluctuation of the H₂ lines at a redshift around z = 12. We find that the noise requirements suggest an instrument that is roughly 10 times better in detector noise than current SPICA/BLISS for a reliable detection.

Zoom In Zoom Out Reset image size

In addition to measuring the intensity fluctuations, we also explore the detection of the H₂ point sources at high redshifts. In Figure 7, we estimate the number of the H₂ sources, for the 0–0S(3) line per deg² with flux greater than a given value for three redshift ranges 10 ⩽ z ⩽ 15, 15 ⩽ z ⩽ 20, and 20 ⩽ z ⩽ 25. The uncertainty for 10 ⩽ z ⩽ 15 is shown as an example which is derived from the uncertainties of the simulations. The flux limits of a pencil-beam survey with a SPICA/BLISS-like instrument and a 10× better SPICA/BLISS survey for 1σ detection with 250 hr of total integration time are also shown in vertical dash-dotted lines. We find it is hard to detect the H₂ sources even for the redshift range 10 ⩽ z ⩽ 15 using the SPICA/BLISS-like experiment. The number count of the H₂ sources at 10 ⩽ z ⩽ 15 is around 10⁻⁵ per deg² for the SPICA/BLISS-like survey. For the instrument that is 10× better in sensitivity than SPICA, we find that in this first estimate of H₂ counts, we can aim to get about 10 sources per deg² at 10 ⩽ z ⩽ 15.

Zoom In Zoom Out Reset image size

In this paper, we propose intensity mapping of H₂ rotational and vibrational lines to detect the primordial gas distribution at large scales during the pre-reionization epochs at z > 10. At such high redshifts, the molecular hydrogen takes the role of main coolant that leads to the formation of the first stars and galaxies and the detection of the H₂ power spectrum can reveal details about the halo mass scales which first form stars and galaxies in the universe.

We first estimate the cooling rates for both H₂ rotational and vibrational lines with the help of fitting results from Hollenbach & McKee (1979) and Hollenbach & McKee (1989). We find the rotational lines are dominant at low gas temperature, while the vibrational lines are stronger at high gas temperature. Also, the slope of the cooling coefficient–temperature relation for the vibrational lines is steeper than that for the rotational lines. We then derive the gas number density, temperature, and H₂ fraction as functions of the halo radius and estimate the relation of the H₂ luminosity and halo mass.

Next, we calculate the mean intensity for several H₂ lines at different redshifts and find 0–0S(3) is the brightest line for 5 ⩽ z ⩽ 30 (≃ 5.9 Jy sr⁻¹ at z = 15). Note that we do not consider the dissociation effects of the H₂ by the Population III and Population II stars at 5 ⩽ z ⩽ 10 that could suppress the H₂ emission significantly in this redshift range. Finally, we evaluate the clustering and shot noise of the auto and cross-power spectra for the H₂ lines at z = 15. We find 0–0S(3) ×0–0S(5) is the strongest cross-power spectrum at z = 15. We propose such a cross-power spectrum for an experimental measurement as it has the advantage that it can minimize astrophysical line confusion from low-redshift galaxies.

In order to consider potential detection of these mid-IR molecular lines, we evaluate the errors of the H₂ auto and cross-power spectrum at z = 12 for a SPICA/BLISS-like and a design that is 10× better than the current instrumental parameters. We find the S/N for the z = 12 cross-power spectrum detection is around 0.1 and 5 for these two experiments. We also estimate the detectability of H₂ point sources over different redshift ranges. We find that a SPICA/BLISS-like instrument is not able to detect H₂ sources for z > 10, but an instrument with 10× better sensitivity than SPICA/BLISS should be able to detect about 10 sources per deg² for 10 < z < 15. We encourage further work on this topic to fully account for dissociation as stars and galaxies form and additional formation mechanisms, such as shock heating, that results in H₂ line emission from low-redshift galaxies.

This work was supported by NSF CAREER AST-0645427. M.G.S. acknowledges support from FCT-Portugal under grant PTDC/FIS/100170/2008. We thank Matt Bradford for helpful discussions and questions that motivated this paper.

Here, we list the fitting formulae of the collisional de-excitation coefficients $C^{\rm H, H_2}_J$ for both H₂–H and H₂–H₂ collisions from Hollenbach & McKee (1979) and Hollenbach & McKee (1989). For the rotational cooling in v = 0, we have

$$\begin{eqnarray} C_J^{\rm H}(T) &=& \left(\frac{10^{-11}T_3^{0.5}}{1+60T_3^{-4}}+10^{-12}T_3\right) \nonumber\\ && \times \left\lbrace 0.33+0.9\ {\rm exp}\left[- \left(\frac{J-3.5}{0.9}\right)^2 \right] \right\rbrace \,{\rm cm^3\,s^{-1}}, \nonumber\\ \end{eqnarray} \tag{ A1 }$$

$$\begin{eqnarray} C_J^{\rm H_2}(T) &=& (3.3\times 10^{-12}+6.6\times 10^{-12}T_3) \nonumber\\ && \times \left\lbrace 0.276J^2\ {\rm exp}\left[ -\left(\frac{J}{3.18}\right)^{1.7}\right] \right\rbrace \,{\rm cm^3\,s^{-1}}, \nonumber\\ \end{eqnarray} \tag{ A2 }$$

where T₃ = T/1000 K and T is the gas temperature. For the vibrational cooling between v = 1 and 0, we have

$$\begin{equation} C_{\rm 10}^{\rm H}(T) = 1.0\times 10^{-12}\ T^{0.5}\ {\rm exp}[-(1000/T)]\,{\rm cm^3\,s^{-1}}, \end{equation} \tag{ A3 }$$

$$\begin{equation} C_{\rm 10}^{\rm H_2}(T) = 1.4\times 10^{-12}\ T^{0.5} {\rm exp}{-[18100/(T+1200)]}\,{\rm cm^3\,s^{-1}}. \vspace*{5pt} \end{equation} \tag{ A4 }$$

The wavenumbers k = 1/λ of the H₂ energy levels for J = 0, ..., 13 at v = 0 and 1 from Dabrowski (1984) are listed in Table 2. We can get the energy for each level by E_J = h_Pck, where h_P is the Planck constant and c is the speed of light. The wavelengths of the H₂ line hence could be derived by $\lambda _{\rm H_2}=h_{\rm P}c/(E_J-E_{J^{\prime }})$.