Full-sky Ray-tracing Simulation of Weak Lensing Using ELUCID Simulations: Exploring Galaxy Intrinsic Alignment and Cosmic Shear Correlations

In general, for a source galaxy with the observed angular position ${\boldsymbol{\theta }}$ and its real position ${\boldsymbol{\beta }}$, one can characterize the deformation effect of cosmic shear through the distortion matrix (Schneider et al. 1992; Jain et al. 2000)

$$\begin{eqnarray}{ \mathcal A }({\boldsymbol{\theta }})=\displaystyle \frac{\partial {\boldsymbol{\beta }}}{\partial {\boldsymbol{\theta }}}\equiv \left(\begin{array}{ll}1-\kappa -{\gamma }_{1} & -{\gamma }_{2}-\omega \\ -{\gamma }_{2}+\omega & 1-\kappa +{\gamma }_{1}\end{array}\right),\end{eqnarray} \tag{ 1 }$$

where κ is the convergence, ${\boldsymbol{\gamma }}={\gamma }_{1}+i{\gamma }_{2}$ defines the complex shear in lensing, and the additional antisymmetric quantity, ω, describes the overall rotation in the lensed images. In the weak-lensing regime (i.e., $\kappa ,{\boldsymbol{\gamma }}\ll 1$) and to the linear order, the components of the matrix are related to the second derivatives of the gravitational potential as (Bartelmann & Schneider 2001; Hilbert et al. 2009; Kilbinger 2015)

$$\begin{eqnarray}{{ \mathcal A }}_{{ij}}({\boldsymbol{\theta }},\chi ) & = & {\delta }_{{ij}}-\displaystyle \frac{2}{{c}^{2}}{\displaystyle \int }_{0}^{\chi }d\chi ^{\prime} \displaystyle \frac{r(\chi -\chi ^{\prime} )r(\chi ^{\prime} )}{r(\chi )}\\ & & \times {{\rm{\Phi }}}_{,{ij}}(r(\chi ^{\prime} ){\boldsymbol{\theta }},\chi ^{\prime} ),\end{eqnarray} \tag{ 2 }$$

where ${\delta }_{{ij}}$ is the Kronecker delta, c is the speed of light, χ is the comoving distance, and $r(\chi )$ is the comoving angular diameter distance. According to the Poisson equation, the gravitational potential Φ can be related to the density contrast δ. Hence, the convergence κ can be expressed as a weighted integral of the overdensity δ along the line of sight,

$$\begin{eqnarray}&&\kappa ({\boldsymbol{\theta }},\chi )=\displaystyle \frac{3{H}_{0}^{2}{{\rm{\Omega }}}_{{\rm{m}}}}{2{c}^{2}}{\int }_{0}^{\chi }d\chi ^{\prime} \displaystyle \frac{r(\chi -\chi ^{\prime} )r(\chi ^{\prime} )}{r(\chi )}\displaystyle \frac{\delta (r(\chi ^{\prime} ){\boldsymbol{\theta }},\chi ^{\prime} )}{a(\chi ^{\prime} )},\end{eqnarray} \tag{ 3 }$$

where H₀ is the Hubble constant, ${{\rm{\Omega }}}_{{\rm{m}}}$ is the matter density in units of the critical density, and $a(\chi ^{\prime} )$ is the scale factor at $\chi ^{\prime}$.

In the flat-sky limit, the power spectrum of the convergence ${C}^{\kappa \kappa }({\ell })$ on the modulus ℓis known as the two-point correlation in Fourier space,

$$\begin{eqnarray}&&\langle \tilde{\kappa }({\boldsymbol{\ell }}){\tilde{\kappa }}^{* }({\boldsymbol{\ell }}^{\prime} )\rangle ={(2\pi )}^{2}{\delta }_{{\rm{D}}}({\boldsymbol{\ell }}-{\boldsymbol{\ell }}^{\prime} ){C}^{\kappa \kappa }({\ell }),\end{eqnarray} \tag{ 4 }$$

where ${\delta }_{{\rm{D}}}({\boldsymbol{\ell }})$ is the Dirac delta function. Using Equation (3), one can derive the angular power spectrum of the convergence field in the Limber approximation,

$$\begin{eqnarray}&&{C}^{\kappa \kappa }({\ell })={\int }_{0}^{{\chi }_{{\rm{H}}}}d\chi \displaystyle \frac{W{(\chi )}^{2}}{r{(\chi )}^{2}}{P}_{\delta }\left(k=\displaystyle \frac{{\ell }}{r(\chi )},\chi \right),\end{eqnarray} \tag{ 5 }$$

where ${P}_{\delta }(k,\chi )$ is the 3D power spectrum of the matter distribution at a given comoving distance χ, and the integral is calculated along the line of sight to the comoving horizon distance ${\chi }_{{\rm{H}}}$. Here, the weight function $W(\chi )$ is defined as

$$\begin{eqnarray}&&W(\chi )=\displaystyle \frac{3{H}_{0}^{2}{{\rm{\Omega }}}_{{\rm{m}}}}{2{c}^{2}}\displaystyle \frac{r({\chi }_{{\rm{H}}}-\chi )r(\chi )}{r({\chi }_{{\rm{H}}})}\displaystyle \frac{1}{a(\chi )}.\end{eqnarray} \tag{ 6 }$$

From the nonlinear theoretical models, such as the Halofit model (Smith et al. 2003; Takahashi et al. 2012), one can predict the convergence power spectrum ${C}^{\kappa \kappa }({\ell })$ from the nonlinear ${P}_{\delta }(k)$. Therefore, the weak-lensing survey can be used to probe the gravitational growth of the density structure.

While dealing with full-sky measurements, it is useful to note that the angular power spectrum of the weak-lensing convergence κ and complex shear ${\boldsymbol{\gamma }}$ can be derived from the spin-s spherical harmonics ${}_{s}{Y}_{{\ell }}^{m}$ (Stebbins 1996). The relations of power spectra between the convergence and the shear E- and B-modes (Schneider et al. 2002; Bunn 2003; Bunn et al. 2003; Zhao & Baskaran 2010) have been derived by Hu (2000) for an all-sky lensing deformation tensor field. Here we briefly summarize the spin-s spherical harmonic decomposition of the full-sky lensing, referring the reader to Hu (2000) for detailed discussions of the power spectrum in weak lensing.

As reviewed in the appendix of Becker (2013), the convergence, lensing shear, and rotation in the distortion matrix (Equation (1)) can be decomposed by the spherical harmonics (Hu 2000; Becker 2013),

$$\begin{eqnarray}&&\kappa (\hat{{\boldsymbol{n}}})=-\displaystyle \frac{1}{2}\displaystyle \sum _{{\ell }m}{\ell }({\ell }+1){\phi }_{{\ell }m}{Y}_{{\ell }}^{m}(\hat{{\boldsymbol{n}}})\end{eqnarray} \tag{ 7 }$$

$$\begin{eqnarray}&&{\gamma }_{1}{(\hat{{\boldsymbol{n}}})\pm i{\gamma }_{2}(\hat{{\boldsymbol{n}}})=\displaystyle \frac{1}{2}\displaystyle \sum _{{\ell }m}\sqrt{\displaystyle \frac{({\ell }+2)!}{({\ell }-2)!}}({\phi }_{{\ell }m}\pm i{{\rm{\Omega }}}_{{\ell }m})}_{\pm 2}{Y}_{{\ell }}^{m}(\hat{{\boldsymbol{n}}})\end{eqnarray} \tag{ 8 }$$

$$\begin{eqnarray}&&\omega (\hat{{\boldsymbol{n}}})=-\displaystyle \frac{1}{2}\displaystyle \sum _{{\ell }m}{\ell }({\ell }+1){{\rm{\Omega }}}_{{\ell }m}{Y}_{{\ell }}^{m}(\hat{{\boldsymbol{n}}}),\end{eqnarray} \tag{ 9 }$$

where ϕ is the lensing deflection potential, Ω is the pseudo-scalar potential (as described by Stebbins 1996), and $\hat{{\boldsymbol{n}}}$ denotes a given position on the sky. Consequently, the different power spectra can be related as

$$\begin{eqnarray}&&{C}^{\kappa \kappa }({\ell })=\displaystyle \frac{1}{4}{{\ell }}^{2}{({\ell }+1)}^{2}{C}^{\phi \phi }({\ell }),\end{eqnarray} \tag{ 10 }$$

$$\begin{eqnarray}&&{C}^{\omega \omega }({\ell })=\displaystyle \frac{1}{4}{{\ell }}^{2}{({\ell }+1)}^{2}{C}^{{\rm{\Omega }}{\rm{\Omega }}}({\ell }),\end{eqnarray} \tag{ 11 }$$

$$\begin{eqnarray}&&{C}^{\mathrm{EE}}({\ell })=\displaystyle \frac{1}{{{\ell }}^{2}{({\ell }+1)}^{2}}\displaystyle \frac{({\ell }+2)!}{({\ell }-2)!}{C}^{\kappa \kappa }({\ell }),\end{eqnarray} \tag{ 12 }$$

$$\begin{eqnarray}&&{C}^{{BB}}({\ell })=\displaystyle \frac{1}{{{\ell }}^{2}{({\ell }+1)}^{2}}\displaystyle \frac{({\ell }+2)!}{({\ell }-2)!}{C}^{\omega \omega }({\ell }).\end{eqnarray} \tag{ 13 }$$

Thus, in the flat-sky limit, one has $\tfrac{1}{{{\ell }}^{2}{({\ell }+1)}^{2}}\tfrac{({\ell }+2)!}{({\ell }-2)!}\approx 1$, showing that ${C}^{\mathrm{EE}}({\ell })$ equals ${C}^{\kappa \kappa }({\ell })$ at small scales (Kitching et al. 2016; Kilbinger et al. 2017).

Weak lensing will induce an additional coherent deformation in the intrinsic galaxy shape, which means that the measured ellipticity ${{\boldsymbol{\epsilon }}}^{(\mathrm{obs})}$ ⁹ of a galaxy can be expressed as (Bartelmann & Schneider 2001; Meylan et al. 2006)

$$\begin{eqnarray}&&{{\boldsymbol{\epsilon }}}^{(\mathrm{obs})}={\boldsymbol{g}}+{{\boldsymbol{\epsilon }}}^{({\rm{I}})}+{{\boldsymbol{\epsilon }}}^{(\mathrm{rnd})},\end{eqnarray} \tag{ 14 }$$

where ${\boldsymbol{g}}$ is the reduced shear, defined as ${\boldsymbol{g}}={\boldsymbol{\gamma }}/(1-\kappa )$, and ${{\boldsymbol{\epsilon }}}^{(\mathrm{rnd})}$ denotes the noise part in galaxy shape measurements, which is assumed to be uncorrelated with the other components. In the weak-lensing regime, κ is small and the ${\boldsymbol{\gamma }}\simeq {\boldsymbol{g}}$ assumption is often made. The intrinsic shape of a galaxy is described as ${{\boldsymbol{\epsilon }}}^{({\rm{I}})}$. Ideally, if the intrinsic ellipticities of galaxies are isotropic, the lensing shear ${\boldsymbol{g}}$ can be derived by averaging over a population of galaxies. However, it is not the case for real data because of the presence of correlated IAs of observed galaxies.

The observed two-point shear correlation function consists of the following contributions (Troxel & Ishak 2015; Jee et al. 2016; Krause et al. 2016):

$$\begin{eqnarray}&&\langle {{\boldsymbol{\epsilon }}}_{i}^{(\mathrm{obs})}{{\boldsymbol{\epsilon }}}_{j}^{(\mathrm{obs})}\rangle =\langle {{\boldsymbol{g}}}_{i}{{\boldsymbol{g}}}_{j}\rangle +\langle {{\boldsymbol{\epsilon }}}_{i}^{({\rm{I}})}{{\boldsymbol{g}}}_{j}\rangle +\langle {{\boldsymbol{\epsilon }}}_{i}^{({\rm{I}})}{{\boldsymbol{\epsilon }}}_{j}^{({\rm{I}})}\rangle ,\end{eqnarray} \tag{ 15 }$$

where we assume that the two observed galaxies are located at the redshifts z_i and z_j (with ${z}_{i}\leqslant {z}_{j}$). The first term, $\langle {{\boldsymbol{g}}}_{i}{{\boldsymbol{g}}}_{j}\rangle$, represents the shear–shear correlation, GG, which is the weak-lensing signal we want to extract. The correlation $\langle {{\boldsymbol{\epsilon }}}_{i}^{({\rm{I}})}{{\boldsymbol{g}}}_{j}\rangle$, often named as GI, is the cross term between gravitational shear and intrinsic ellipticity. This correlation comes from the fact that the shape of a distant galaxy j is lensed by the foreground gravitational potential, in which galaxy i is intrinsically aligned with the underlying tidal field (Hirata & Seljak 2004). Since nearby galaxies are affected by the same environment, the intrinsic–intrinsic correlation $\langle {{\boldsymbol{\epsilon }}}_{i}^{({\rm{I}})}{{\boldsymbol{\epsilon }}}_{j}^{({\rm{I}})}\rangle$, often referred to as the II term, may be non-zero. Both the II and GI correlations can contaminate our measurements of cosmic shear and are important to quantify, particularly in accurate shear measurements expected from future large lensing surveys (Krause et al. 2016).

In order to model the II and GI parts in the measurements, Hirata & Seljak (2004), Bridle & King (2007), and Joachimi et al. (2011) developed a nonlinear intrinsic alignment model based on the work of Catelan & Porciani (2001). In this model, the power spectra of the II and GI contributions are related to the nonlinear matter power spectrum as ${P}_{\mathrm{II}}(k,z)\,={f}^{2}(z){P}_{\delta }(k,z)$ and ${P}_{\mathrm{GI}}(k,z)=f(z){P}_{\delta }(k,z)$, respectively. Here, the modification factor f(z) is defined as

$$\begin{eqnarray}&&f(z)=-{A}_{\mathrm{IA}}{C}_{1}{\rho }_{{\rm{c}}}\displaystyle \frac{{{\rm{\Omega }}}_{{\rm{m}}}}{D(z)}{\left(\displaystyle \frac{1+z}{1+{z}_{0}}\right)}^{\eta }{\left(\displaystyle \frac{L}{{L}_{0}}\right)}^{\beta },\end{eqnarray} \tag{ 16 }$$

where ${A}_{\mathrm{IA}}$ is a free parameter, ${C}_{1}=5\times {10}^{-14}{h}^{-2}\,{M}_{\odot }^{-1}\,{{\rm{Mpc}}}^{3},$ ${\rho }_{c}$ is the critical density at present time, and D(z) is the linear growth factor (normalized to unity at z = 0). The free parameters η and β account for the dependence on redshift and luminosity around the pivot redshift z₀ and luminosity L₀. Following the discussion of Joudaki et al. (2017) based on the CFHTLenS data, we fix $\eta =0$ and $\beta =0$ in our model fitting. These formulas are used in Section 5.3.2 to fit the measurements of the GI and II terms from our simulation. There we will see that the sign of ${A}_{\mathrm{IA}}$ actually depends on galaxy type.

We use two sets of different N-body cosmological simulations. The first one is part of the ELUCID project (Li et al. 2016; Wang et al. 2014, 2016; Tweed et al. 2017), which is run with 3072³ dark matter particles in a cubic box with ${L}_{\mathrm{box}}=500{h}^{-1}\,\mathrm{Mpc}$ on each side. This simulation is referred to as L500 in the following. The cosmological parameters of L500 are from the WMAP9 cosmology (Hinshaw et al. 2013). The second simulation is the Pangu simulation (PS-I), performed by the Computational Cosmology Consortium of China (Li et al. 2012), which has the same number of particles as L500, but with a box size of $1000\,{h}^{-1}\,\mathrm{Mpc}$ on each side. The cosmological parameters of PS-I are from WMAP7 (Komatsu et al. 2011). In Table 1, we list the parameters of the two N-body simulations. Both simulations were run using the GADGET-2 code (Springel 2005).

Table 1.  Parameters of the Two N-body Simulations

Simulation	${{\rm{\Omega }}}_{{\rm{m}}}$	${{\rm{\Omega }}}_{{\rm{\Lambda }}}$	h	${\sigma }_{8}$	${L}_{\mathrm{box}}/{h}^{-1}\,\mathrm{Mpc}$	${m}_{{\rm{p}}}/({10}^{10}{h}^{-1}{M}_{\odot })$	${l}_{\mathrm{soft}}/{h}^{-1}\mathrm{kpc}$
PS-I	0.260	0.740	0.710	0.80	1000	0.249	7.0
L500	0.282	0.718	0.697	0.82	500	0.034	3.5

Note. Cosmological parameters are given as ${{\rm{\Omega }}}_{{\rm{m}}},{{\rm{\Omega }}}_{{\rm{\Lambda }}}$, h, and ${\sigma }_{8}$. ${L}_{\mathrm{box}}$ is the box size, ${m}_{{\rm{p}}}$ is the particle mass, and ${l}_{\mathrm{soft}}$ is the softening length.

Download table as:  ASCIITypeset image

With its higher mass resolution, the L500 simulation is used to generate galaxies from a semi-analytical model (Luo et al. 2016). This model is based on the L-Galaxies model developed by the Munich group (e.g., Guo et al. 2013; see Section 4.1 for more details). We do not produce mock galaxies using the PS-I simulation due to its lower mass resolution, but use it as a reference to check our calculation of the convergence power spectrum on large scales.

To perform a ray-tracing simulation with full-sky coverage, we follow the multiplane algorithm developed by Das & Bode (2008), Teyssier et al. (2009), and Becker (2013). In order to control the residual in the solution of the lensing potential, we implement an iterative spherical harmonic analysis scheme, which is different from the multigrid method adopted by Becker (2013). In the following, we briefly summarize the main procedures; more details can be found in Appendix A.

To trace the trajectory of a light beam, we first employ the N-body simulations to build a light cone to redshift ${z}_{\max }\sim 2.0$. In practice, the simulation boxes are divided into sets of small cubic boxes with $\sim 100\,{h}^{-1}\,\mathrm{Mpc}$ on each side. These cell boxes are appropriately piled together so as to cover the past light cone from z = 0 to $z={z}_{\max }$. For both L500 and PS-I described above, a full-sky light cone can be constructed in this manner. Note that the size of our simulation box is relatively small compared with the comoving distance to redshift ${z}_{{\rm{s}}}=2.0$, and periodic effects will show up at several specific directions, especially along the box axes, but the effects disappear very quickly apart from these directions.

Then, each light cone is divided into a set of spherical shells with a thickness of $50\,{h}^{-1}\,\mathrm{Mpc}$ centered at the observer, and the dark matter distribution in the corresponding shells is projected onto pixels defined by the HEALPix¹⁰ tessellation (Górski et al. 2005; Calabretta & Roukema 2007). The HEALPix resolution parameter is set to ${N}_{\mathrm{side}}=8192$, which gives an angular resolution of ∼0.43 arcmin.¹¹ The projected surface mass densities are calculated for each shell using the SPH algorithm (Li et al. 2005; Springel 2010). We use the nearest 64 particles to define the kernel size, but keep the smoothing length larger than two HEALPix cells in high-density regions. The lensing potential for the nth shell, ${\phi }_{{\ell }m}^{(n)}$, is then obtained, using the Poisson equation, from the mass density shell after applying iterative spherical harmonic transformations (referred to as HEALPix predefined functions),

$$\begin{eqnarray}&&-{\ell }({\ell }+1){\phi }_{{\ell }m}^{(n)}=2{\kappa }_{{\ell }m}^{(n)}.\end{eqnarray} \tag{ 17 }$$

To perform multisphere ray-tracing simulations, we set the initial positions of the ray beams at the centers of the HEALPix cells and propagate light rays from the observer to a desired redshift applying deflection angle,

$$\begin{eqnarray}&&{\alpha }_{{\ell }m}^{(n)}=-\sqrt{{\ell }({\ell }+1)}{\phi }_{{\ell }m}^{(n)}.\end{eqnarray} \tag{ 18 }$$

From our lensing simulation, we evaluate the distortion matrix ${ \mathcal A }$ on each lensing shell and construct the full-sky map of the convergence and lensing shear. As an illustration, Figure 1 shows one realization of our simulated full-sky convergence map, κ, for sources at redshift ${z}_{s}=1.0$. In Figure 2, we show the power spectra measured from the PS-I (left panels) and L500 simulations (right panels). The top panels show the angular power spectra of the convergence (red solid line), the shear E-mode (blue) and B-mode (magenta), and the rotation mode (cyan). We also show the prediction from the Born approximation (Cooray & Hu 2002) by stacking the density field along the line of sight in our mock light cone as the gray dashed line and the theoretical prediction from the revised nonlinear Halofit (Takahashi et al. 2012; Peacock & Smith 2014) as the black dashed line. The middle panels of Figure 2 show the relative deviations between the convergence powers measured from our ray-tracing simulation and the theoretical predictions; the relative deviation between the ray-tracing simulation and Halofit model is defined by ${{\rm{\Delta }}}_{\mathrm{Halofit}}=[{C}^{\kappa \kappa }({\ell })-{C}^{\mathrm{Halofit}}({\ell })]/{C}^{\mathrm{Halofit}}({\ell })$.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The measured convergence power from our ray-tracing simulation agrees well with the theoretical prediction of Halofit, and the relative error is less than 10% at ${\ell }\lesssim 4000$. At large scales, ${\ell }\lesssim 10$, the power from the PS-I simulation is in better agreement with theoretical predictions than the L500 simulation, as is expected from the fact that PS-I has a larger box to represent the matter power on large scales. The prediction from Born approximation is closer to that of Halofit at small scales than our ray-tracing simulation, because both Born approximation and Halofit are based on first-order approximations. This good agreement, within 10% for ${\ell }\lesssim 6000$, roughly on scales larger than the smoothing scale in our simulation, indicates that the revised Halofit model (Takahashi et al. 2012) provides a good approximation to the nonlinear matter power spectrum. In addition, since the power spectra from our Born approximation and full-sky ray-tracing simulation are both based on the same convergence κ maps, the difference between the two is not caused by the smoothing effect. The gray lines in the middle panels of Figure 2 show that the Born approximation can cause a deviation of more than 10% for ${\ell }\gt 6000$.

The bottom panels show the ratio of shear E-mode and B-mode power spectra relative to the measured convergence and rotation mode spectra, respectively. As shown in Section 2.2, at high-ℓ(small scales), one has ${C}^{\mathrm{EE}}({\ell })={C}^{\kappa \kappa }({\ell })$ from the full-sky weak lensing. Our corresponding measurements on small scales are indeed consistent with this expectation. At large scales (low-ℓ), we must take into account the extra factor of $({\ell }-1)({\ell }+2)/{\ell }/({\ell }+1)$ to explain the difference between the power spectrum of the shear E-mode and convergence spectrum. As discussed in Becker (2013), we also measure the power spectra of the B-mode and the rotation mode from our lensing simulation. We find that the B-mode power is effectively suppressed relative to the E-mode by more than four orders of magnitude. Moreover, the power ratio between the B-mode and the rotation mode shows that the extra numerical B-mode in our simulation is negligible, and the accuracy of our shear map is only limited by the smoothing length at small scales. We refer the reader to Figure 13 in Appendix A for more details. Compared with the Becker (2013) results, the predictions of our simulations for the convergence and shear power spectra are more accurate extending to higher ℓ.

The model galaxies are produced using the semi-analytical model of Luo et al. (2016), which is based on the Guo et al. (2013) model, a version of the Munich Semi-analytical model called L-Galaxies. In a semi-analytical model, a galaxy population is assigned to dark matter halos on the basis of simple assumptions of many physical processes. As a first step in our implementation of L-Galaxies, dark matter halos are identified in the N-body simulation using the standard Friends-of-Friends (FOF) algorithm. Only halos that contain at least 20 particles are used. The subhalos within each FOF halo are identified with the SUBFIND algorithm (Springel et al. 2001, 2005). Merger trees of these dark matter (sub)halos can be constructed by linking the progenitors of a halo in different snapshots. Galaxies are assumed to form at the centers of the dark matter halos according to analytical prescriptions of relevant physical processes, such as gas cooling, star formation, supernova, and black hole feedback. For the details of L-Galaxies, we refer the reader to Guo et al. (2013). Luo et al. (2016) improved the prescription for low-mass galaxies, especially satellite galaxies, by including additional physics of cold gas stripping and an analytical modeling of orphan galaxies. In this model, the stellar mass function and Hi and H₂ mass functions are tuned to match the observational data (Keres et al. 2003; Zwaan et al. 2005; Baldry et al. 2008; Li & White 2009). The fraction and spatial distributions of the central versus satellite galaxies are reproduced roughly correctly by the model, as shown in Luo et al. (2016) and Guo et al. (2013).

To describe the morphology of a galaxy, we use the ratio between the bulge and the total mass (B/T), which can be predicted from the semi-analytical model, for the classification of "early-type" and "late-type" galaxies (Parry et al. 2009). Following Joachimi et al. (2013a), we adopt $B/T=0.6$ to classify the model galaxies into early or late types. Figure 3 shows the cumulative probability distribution of the B/T of our simulated galaxies at three different redshifts. More than 80% of all galaxies are late types, and the fraction is slightly higher at high redshifts. This fraction of early/late types is consistent with that found in other studies (e.g., Guo et al. 2013). However, compared to hydrodynamical simulations, semi-analytical models are less powerful in predicting the shapes of galaxies. To proceed, we have to assign shapes to galaxies and their images with some simplified prescriptions, as described below.

Zoom In Zoom Out Reset image size

One common assumption is that the shape of an elliptical galaxy roughly follows that of its host dark matter halo, while the rotation axis of a spiral galaxy is determined by the spin of its halo (e.g., Kang et al. 2007; Okumura et al. 2009; Agustsson & Brainerd 2010). Joachimi et al. (2013b, hereafter J13) used this assumption and studied the alignment of galaxies from the Millennium Simulation (Springel et al. 2005). Here, we follow J13 to assign shapes to model galaxies.

To assign a shape to a model galaxy using the mass distribution of its dark matter halo, we need to distinguish between central and satellite galaxies. A central galaxy is assumed to be located at the center of a dark matter halo and its shape may be related to that of the host halo. The shape of a dark matter halo is usually defined using the inertia tensor I_ij (Bailin & Steinmetz 2005),

$$\begin{eqnarray}&&{I}_{{ij}}=\displaystyle \sum _{n=1}^{{N}_{{\rm{p}}}}{m}_{p}{x}_{i,n},{x}_{j,n},\end{eqnarray} \tag{ 19 }$$

where ${N}_{{\rm{p}}}$ denotes the particle number of the FOF halo and ${{\boldsymbol{x}}}_{n}$ is the position of the nth particle with respect to the center of the halo. By diagonalizing the inertia tensor ${\boldsymbol{I}}$, one can obtain the eigenvalues ${\lambda }_{1}\leqslant {\lambda }_{2}\leqslant {\lambda }_{3}$ and the corresponding eigenvectors that define a triaxial ellipsoid and its orientation. It has been argued that a minimum number of ${N}_{{\rm{p}}}=300$ is needed to ensure an accurate measurement of the halo shape (Jing 2002; Bett et al. 2007). Constrained by the resolution of our simulation, we have to reduce the number limit to 100. As shown in J13, a minimum number of 100 can lead to a $\sim 10 \%$ deviation in the axis ratio and a 10^◦ deviation in the orientation angle. Because of the magnitude limitation, the corresponding halos always have more than 100 particles in our mock catalogs for DLS and KiDS-450. Once the three-dimensional shapes and orientations of galaxies are obtained, we project them onto the sky to obtain the projected ellipses (Galletta 1983; Binney 1985). Details on how to make the projections can be found in Appendix B.

In addition to the above model, in which perfect alignments are assumed between elliptical galaxies and their host halos, J13 also considered a misalignment model in which the major axis of the central elliptical galaxy is misaligned with that of the halo, with the misalignment angle obeying a Gaussian distribution with zero mean and a dispersion of 35^◦. This is motivated by the finding that such a misalignment is needed to explain the alignment between luminous red galaxies on large scales (e.g., Kang et al. 2007; Faltenbacher et al. 2009; Okumura et al. 2009; Li et al. 2013). We will come back in Section 5.3.2 to discuss the effect of such a misalignment on shear correlation functions.

For a central late-type galaxy, defined by $B/T\lt 0.6$, J13 assigned its shape according to the angular momentum vector of the host halo,

$$\begin{eqnarray}&&{\boldsymbol{L}}=\displaystyle \sum _{n=1}^{{N}_{{\rm{p}}}}{m}_{p}{{\boldsymbol{x}}}_{n}\times {{\boldsymbol{v}}}_{n},\end{eqnarray} \tag{ 20 }$$

where ${{\boldsymbol{v}}}_{n}$ is the velocity of the nth halo particle relative to the halo center. The angular momentum ${\boldsymbol{L}}$ defines a circular disk in the halo, and the complex ellipticity is obtained by projecting the disk along the line of sight. As described in J13, the apparent axis ratio of the projected ellipse is

$$\begin{eqnarray}&&r=\displaystyle \frac{| {L}_{\mathrm{los}}| }{| {\boldsymbol{L}}| }+{r}_{d}\sqrt{1-\displaystyle \frac{{L}_{\mathrm{los}}^{2}}{| {\boldsymbol{L}}{| }^{2}}},\end{eqnarray} \tag{ 21 }$$

where r_d is the ratio between the disk thickness and diameter, and we set ${r}_{d}=0.25$ following J13; ${L}_{\mathrm{los}}$ is the component of ${\boldsymbol{L}}$ along the line of sight. The ellipticity of the mock galaxy is then given by

$$\begin{eqnarray}&&\epsilon =(1-r)/(1+r).\end{eqnarray} \tag{ 22 }$$

For a satellite galaxy, on the other hand, the original dark matter halo associated with it may have suffered huge mass loss after it is accreted into a big halo, depending on the infall time and orbit (Cole et al. 2000). It is thus unclear how the shapes of satellites are connected to those of the dark matter subhalos associated with them. Some earlier investigations (e.g., Pereira et al. 2008) have shown that the tidal torque of the host halo can induce a correlation between the subhalo orientation and its direction to the host center. Therefore, J13 assigned the shape of an early-type satellite galaxy by randomly choosing the three-dimensional axis ratios from a halo sample with more than 300 particles and then made its major axis point to the central galaxy. For a late-type satellite, its spin is assumed to be perpendicular to the line connecting the satellite to the central galaxy, and the ellipticity is obtained by projecting the disk onto the sky, as done for the central spiral galaxies.

The above assumptions for the shapes and orientations of satellite galaxies are clearly too idealistic. In fact, while orbiting in their host halos, satellite galaxies may have their radial alignments scrambled. Observationally, measurements of the shapes of faint satellite galaxies are difficult and sensitive to the methods used to derive galaxy shapes (e.g., Hao et al. 2011). Currently, there is no consensus on the alignments of satellites. Some studies have reported the detection of the radial alignment of satellites (Pereira & Kuhn 2005; Agustsson & Brainerd 2006; Faltenbacher et al. 2007; Chisari et al. 2014; Singh et al. 2015; Huang et al. 2016), while others have not found such an alignment (Siverd et al. 2009; Hao et al. 2011; Schneider et al. 2013; Sifón et al. 2015). Because of this uncertainty, we also consider a simple case in which the orientations of satellites, regardless of their types, are randomly distributed in their host halos. We will show in Section 5.3 that satellite alignment has a stronger effect on the shear correlation on smaller scales.

The intrinsic shape correlation function of galaxies, $\eta (r)$, is defined as

$$\begin{eqnarray}&&\eta (r)=\langle {\epsilon }_{t}({\boldsymbol{x}}){\epsilon }_{t}({\boldsymbol{x}}+{\boldsymbol{r}})+{\epsilon }_{\times }({\boldsymbol{x}}){\epsilon }_{\times }({\boldsymbol{x}}+{\boldsymbol{r}}){\rangle }_{{\boldsymbol{x}}},\end{eqnarray} \tag{ 23 }$$

where r is the three-dimensional comoving separation between two galaxies (e.g., Heymans et al. 2006). The quantities ${\epsilon }_{t}$ and ${\epsilon }_{\times }$ are the tangential and cross components of the galaxy ellipticity:

$$\begin{eqnarray}&&{\epsilon }_{t}+i{\epsilon }_{\times }=-{\boldsymbol{\epsilon }}{e}^{-2i\varphi },\end{eqnarray} \tag{ 24 }$$

where φ is the angle between the separation vector of a given galaxy pair and the horizontal axis (Bartelmann & Schneider 2001).

Following Joachimi et al. (2013b), we first measure the intrinsic shape correlation of galaxies in our simulation by projecting semi-analytic galaxies along the line of sight parallel to the edges of the simulation box. In order to estimate the error bars in our measurements, we divide our simulation box into eight equal-sized cubic sub-boxes of $250\,\mathrm{Mpc}/h$, and use the jackknife method to estimate the errors. Figure 4 shows the redshift dependence of the correlation function $\eta (r)$ for the early-type (upper) and late-type (lower) central galaxies in our simulation. It is seen that early-type galaxies have a strong correlation, and the correlation is stronger at higher redshifts. In contrast, late-type galaxies do not show any significant correlation of their projected ellipticities, although some weak positive correlation signals can be seen at small scales, $r\lesssim 8\,\mathrm{Mpc}/h$. The weak/null correlation for spiral galaxies in our simulation is consistent with the non-detection in both observations (e.g., Mandelbaum et al. 2011) and in simulation results (e.g., Joachimi et al. 2013b).

Zoom In Zoom Out Reset image size

In Figure 5, we further investigate the dependence of the shear correlation $\eta (r)$ on halo mass and luminosity for early-type galaxies at z = 1.0. To separate the two dependencies, we select galaxies in small ranges of luminosity and halo mass, and divide galaxies into two subsamples in halo mass (the top panel) and in galaxy luminosity (the bottom panel). It can be seen that there is no significant luminosity dependence for a fixed halo mass, but a significant dependence on halo mass is seen at a fixed luminosity. This dependence of the intrinsic ellipticity correlation on mass and luminosity in our simulation is similar to the results found in J13.

Zoom In Zoom Out Reset image size

Weak lensing induces small correlated distortions in observed galaxy shapes. This correlation can be quantified using different statistics. Observationally, the most direct measurement of the lensing signal is the two-point shear correlation function. The shear–shear correlation between galaxies at a given separation ϑ is estimated as

$$\begin{eqnarray}&&{\xi }_{{tt}}(\vartheta )=\displaystyle \frac{{\displaystyle \sum }_{i,j}{w}_{i}{w}_{j}{\epsilon }_{t,i}{\epsilon }_{t,j}}{{\displaystyle \sum }_{i,j}{w}_{i}{w}_{j}}\end{eqnarray} \tag{ 25 }$$

and

$$\begin{eqnarray}&&{\xi }_{\times \times }(\vartheta )=\displaystyle \frac{{\displaystyle \sum }_{i,j}{w}_{i}{w}_{j}{\epsilon }_{\times ,i}{\epsilon }_{\times ,j}}{{\displaystyle \sum }_{i,j}{w}_{i}{w}_{j}},\end{eqnarray} \tag{ 26 }$$

where w_i is the ellipticity weight of the ith galaxy and ϑ is the angular separation between the galaxy pair.

The two linear combinations of ${\xi }_{{tt}}$ and ${\xi }_{\times \times }$ that are also frequently used in the lensing analysis are

$$\begin{eqnarray}&&{\xi }_{\pm }={\xi }_{{tt}}\pm {\xi }_{\times \times }.\end{eqnarray} \tag{ 27 }$$

The convergence power spectrum ${C}^{\kappa \kappa }({\ell })$ can be related to the estimator ${\xi }_{\pm }$ as

$$\begin{eqnarray}&&{\xi }_{\pm }(\vartheta )=\displaystyle \frac{1}{2\pi }{\int }_{0}^{\infty }d{\ell }{\ell }{J}_{\mathrm{0,4}}({\ell }\vartheta ){C}^{\kappa \kappa }({\ell }),\end{eqnarray} \tag{ 28 }$$

where ${J}_{\mathrm{0,4}}({\ell }\vartheta )$ denotes the zeroth and fourth Bessel functions for ${\xi }_{+}$ and ${\xi }_{-}$, respectively (Schneider et al. 2002).

Tomographic measurement is capable of utilizing the redshift dependence of cosmic shear signals to reveal both the cosmological structure growth and the redshift-dependent geometry in the universe (King & Schneider 2003), and has been widely used in weak-lensing observations, such as CFHTLens (Heymans et al. 2013), DLS (Jee et al. 2016), and KiDS-450 (Hildebrandt et al. 2017). In order to perform a similar tomographic analysis of cosmic shear in our mock observation and compare the prediction with the KiDS and DLS observations, we employ the tomographic redshift bins as used by Hildebrandt et al. (2017) for KiDS-450 and Jee et al. (2016) for DLS to mimic their measurements of cosmic shears.

Figure 6 shows the redshift distributions of the source galaxies in the two surveys, where ${z}_{{\rm{B}}}$ denotes the Bayesian point estimates of the photo-z (Jee et al. 2016; Hildebrandt et al. 2017). It is seen that the two distributions are quite different. In KiDS-450, the shape of the distribution is not regular, with more overlaps between different redshift bins, while the distribution for DLS is more regular and different bins are more clearly separated. For a consistent comparison between model predictions and observations, we adopt their redshift distributions for the source galaxies, and we truncate the source galaxies at z = 2.0.

Zoom In Zoom Out Reset image size

To compare with the survey results, we first divide the full sky into a set of small patches with sizes $\sim {3.6}^{^\circ }\times {3.6}^{^\circ }$ and then randomly select 35 patches in total to cover a field of ∼450 square degrees to match the sky coverage of KiDS-450. By setting a limiting magnitude of $\sim 24.5$ mag in the r-band, the effective number density in our light cone is $n\sim 8\ {\mathrm{arcmin}}^{-2}$, similar to that in the KiDS-450 observations. The DLS is much deeper, with a magnitude $\sim 27\mathrm{th}$ in the r-band, producing an effective number density of $\sim 11\ {\mathrm{arcmin}}^{-2}$ of the source population in five tomographic bins. Given that our lensing simulation is performed to redshift ${z}_{\max }\sim 2.0$, we discard the fifth redshift bin, using only tomographic bins 1–4 of DLS, which gives $n\sim 8\ {\mathrm{arcmin}}^{-2}$. Thus, a set of mock galaxy catalog is constructed to mimic the sky coverage and galaxy number density for each of KiDS-450 and DLS. For each mock, we produce 100 realizations by randomly sampling the patches at different positions to estimate the uncertainties of tomographic shear correlations due to the cosmic variance and sampling noise.

We measure the auto-correlation and cross-correlation functions ${\xi }_{\pm }^{({ij})}$ using the public code Athena,¹² which estimates the second-order shear correlation functions from Equation (27). The superscript $({ij})$ denotes the different redshift bins used for the calculation of the correlation function. In our case, there are four redshift bins labeled from 1 to 4 with increasing redshift.

5.3.1. Model Predictions and Comparison with Observations

In Figures 7 and 8, we show the tomographic shear correlations ${\xi }_{\pm }$ from our model and compare them with the KiDS-450 results. We note that here the orientation of the central galaxy is assumed to follow that of the dark matter halo. Namely, for an elliptical central, its major axis follows that of the halo, while for the spiral central, its orientation is determined by halo spin. For satellite galaxies, regardless of ellipticals or spirals, their orientations are randomly distributed on the sky. The black circles are the model predictions for the ellipticity correlation, which can be directly compared with the data (blue circles). As mentioned in Section 2.3, the predicted ellipticity correlations are combinations of the GG, GI, and II correlations. For simplicity, here we only show the GG terms (the red circles connected by the red lines). We will show the contributions of the II and GI terms in Section 5.3.2.

Zoom In Zoom Out Reset image size

Figures 7 and 8 show that in general the model predictions (black circles) agree well with the KiDS-450 results. To quantify the difference between the model and the data, we calculate the reduced ${\chi }^{2}$ defined by

$$\begin{eqnarray}&&{\chi }_{\pm }^{2}=\displaystyle \frac{1}{n}\sum \displaystyle \frac{{({\xi }_{\pm }^{{\rm{M}}}-{\xi }_{\pm }^{{\rm{D}}})}^{2}}{{({\sigma }_{\pm }^{{\rm{D}}})}^{2}}.\end{eqnarray} \tag{ 29 }$$

Here, n is the number of data points in the tomographic measurements; ${\xi }_{\pm }^{M}$ and ${\xi }_{\pm }^{{\rm{D}}}$ represent the predicted and observed tomographic correlations, respectively; ${\sigma }_{\pm }^{{\rm{D}}}$ is the error in the data. The error bars we predicted only contain the intrinsic ellipticity dispersion of the galaxy, cosmic variance, and shot noise, and they represent the dispersion between the results of our 100 realizations, i.e., the uncertainties of one realization. The uncertainties of ${\xi }_{\pm }^{{\rm{M}}}$ that we predicted (the mean of 100 realizations) are very low, so we do not take them into account in Equation (29). We then find the reduced ${\chi }_{+(-)}^{2}=1.70\ (1.82)$ between our model prediction and the KiDS-450 results. If the correlation between ${\xi }_{+}$ and ${\xi }_{-}$ is taken into account, we can estimate the reduced ${\chi }^{2}$ from the correlation matrix,

$$\begin{eqnarray}&&{\chi }^{2}=\displaystyle \frac{1}{n}\displaystyle \sum _{i,j}{\rm{\Delta }}{\xi }_{i}{C}_{{ij}}^{-1}{\rm{\Delta }}{\xi }_{j},\end{eqnarray} \tag{ 30 }$$

where C_ij is the covariance matrix of the data (Hildebrandt et al. 2017) and ${\rm{\Delta }}{\xi }_{i}$ is the difference between the model prediction and the data in the ith separation bin. This gives a reduced ${\chi }^{2}=1.36$ for the full data vector of KiDS-450, which is slightly higher than the reduced ${\chi }^{2}=1.33$ in the fiducial analysis of KiDS-450 (Hildebrandt et al. 2017). If we calculate the reduced ${\chi }^{2}$ for ${\xi }_{+}$ and ${\xi }_{-}$ separately, we find ${\chi }^{2}=1.23$ and 1.61, respectively.

Zoom In Zoom Out Reset image size

The good agreement between our model and the KiDS-450 data is encouraging, as this is the first time that observational results using lensed images of mock galaxies in N-body simulations combined with a realistic model of galaxy formation have been reproduced. However, inspecting Figure 7 carefully, one can see that for some bins, such as 12, 24, and 34, the model predictions are slightly higher than the data on large scales. Note that the error bars are also larger in these bins, and the data points are not well-described by the best-fitting model of Hildebrandt et al. (2017). To further quantify the systematic deviation between our model prediction and the data, we define the weighted mean deviation as

$$\begin{eqnarray}&&{{\rm{\Delta }}}_{{\rm{m}}}=\displaystyle \frac{1}{n}\sum \displaystyle \frac{{\rm{\Delta }}{\xi }_{i}}{{\sigma }_{i}}.\end{eqnarray} \tag{ 31 }$$

Here, only errors from the data are used as the weight ${\sigma }_{i}$, and the correlation between the measurements is not accounted for. Note that the correlation between the measurements will slightly reduce the deviation, similar to previous analyses for the reduced ${\chi }^{2}$. The systematic bias between the simulations and the observed data can then be estimated in terms of the standard deviation as

$$\begin{eqnarray}&&S=\displaystyle \frac{{{\rm{\Delta }}}_{{\rm{m}}}}{{\sigma }_{{\rm{m}}}},\end{eqnarray} \tag{ 32 }$$

where ${\sigma }_{{\rm{m}}}=1/\sqrt{n}$ is the scatter of ${{\rm{\Delta }}}_{{\rm{m}}}$. In this way, we find that the systematic bias between our model prediction and the KiDS-450 result is S = 1.80 and 1.92 for ${\xi }_{+}$ and ${\xi }_{-}$, respectively. This positive deviation suggests that our predicted correlations are slightly, but systematically, higher than the KiDS-450 data with a significance of $\sim 1.8\sigma$. As the error bars are often correlated between different redshift bins, our analyses of these derivations might be too simplistic. However, we do not intend to quantify the difference between the model and the data in detail, but would like to point out that such a difference could be due to the cosmological parameter, ${\sigma }_{8}\sqrt{{{\rm{\Omega }}}_{{\rm{m}}}/0.3}=0.79$, adopted in our simulation, which is slightly larger than that derived, 0.745 ± 0.039, from the KiDS-450 data (Hildebrandt et al. 2017).

The results shown in Figures 7 and 8 assume that satellites have random orientations. In Joachimi et al. (2013b), the orientations of satellite galaxies are assumed to be radially aligned with central galaxies. Thus, for an early-type satellite, its major axis is assumed to point toward the central galaxy, while for a late-type satellite, its spin is assumed to be perpendicular to the line connecting the satellite to the central galaxy. In Figure 9, we compare the model results obtained by assuming such radial alignments to the KiDS-450 results. We can see that the shear correlations in diagonal panels at lower redshift are much higher than the data on small scales, except for the highest redshift, 4–4 bin, where the model prediction is close to the data. This indicates that the radial alignment model of satellites produces too strong of a correlation on small scales. In addition, this model also leads to a strong negative GI term, suppressing the measured total signal on small scales in the cross-correlation of shears, as seen in the off-diagonal bins (13, 14, 24). This strong positive correlation at small scales in the auto-correlation bins and the strong negative correlation in the cross bins can be explained by the balance between the contributions of the II term and the GI term. With the J13 assumption for satellites, the II term is positive and the GI term is negative. In the auto-correlation bins, the II term is stronger than the GI term, so as to make a strong positive correlation of the total signal at small scales, as shown in Figure 9. On the other hand, in the cross-correlation bins, the contribution of the II term is reduced relative to the GI term, so that the correlation at small scales shows a strong negative correlation in the cross-correlation bins. The reduced ${\chi }^{2}$ can be calculated from the covariance matrix of the data, and we find ${\chi }^{2}=3.73$ for the full data vector of KiDS. Calculating ${\chi }^{2}$ separately for ${\xi }_{+}$ and ${\xi }_{-}$, we get ${\chi }^{2}=4.63$ and 2.76, respectively. The results in Figure 9, therefore, suggest that the radial alignment model for satellite galaxies can be rejected. In what follows, we will only show model predictions in which satellites are assumed to have random orientations.

Zoom In Zoom Out Reset image size

DLS is another weak-lensing survey completed recently (Jee et al. 2013, 2016). Compared to KiDS-450, DLS has a smaller sky coverage of 20 square degrees. We produce mock DLS catalogs following its sky coverage, galaxy number density, and redshift distribution of source galaxies (see Figure 6). Our model predictions and comparisons with DLS are shown in Figure 10. Considering that the data for ${\xi }_{-}$ are not available for DLS (Jee et al. 2016), here we only present the results of tomographic correlations ${\xi }_{+}$. It is seen that the model results (black circles) are lower than the DLS data (blue circles), especially in the lower redshift bins. Note that the error bars in the model are slightly larger. For a DLS-like survey, the reduced ${\chi }_{+}^{2}=1.56$, and the systematic bias between the simulation and the data is $S=-2.57$ for ${\xi }_{+}$. The lower reduced ${\chi }_{+}^{2}$ seems to indicate that the agreement between our model and DLS is slightly better than the agreement with KiDS-450. However, it is clear that the lower reduced ${\chi }_{+}^{2}$ is also related to the fact that the error bars in DLS data are much larger than those of the KiDS data. The large error bars in the DLS data are partly due to its small sky coverage and the smaller sample of galaxies. In the tomographic analysis, the KiDS-450 sample is more than 10 times larger than the DLS sample in terms of the total number of galaxies. The strong negative systematic bias, $S=-2.57$, indicates that the DLS data are systematically higher than our model predictions. Since the cosmic variance has been taken into account in the error bars of the DLS data, such a large systematic deviation can hardly be explained by cosmic variance only.

Zoom In Zoom Out Reset image size

It is unclear what causes the discrepancy ($\sim 2\sigma$ in S₈) between the observational results of DLS and KiDS-450. One potential cause might be from the estimations of photo-z. For KiDS (Hildebrandt et al. 2017), the DIR method is used to estimate the redshift distribution of galaxies, while for DLS (Jee et al. 2016), the BPZ method is adopted. As briefly discussed in Hildebrandt et al. (2017), the ${\chi }^{2}$ can increase by ∼10 when switching from the DIR redshift distribution to the BPZ distribution. They also argued that the deeper DLS data are harder to calibrate. It is beyond the scope of our study to discuss the discrepancy between the KiDS-450 and DLS results in detail. We refer the reader to the paper cited above for more discussions.

5.3.2. The Contributions of the II and GI Terms

Our previous model results in Figure 7 show that there is a difference of around 10% between the GG term and the total shear correlation. The difference is due to a combination of the II and GI terms. In Figure 11, we show the contributions from the two components separately. As some data points are negative, we plot $\theta {\xi }_{+}$ in linear scales, and for clarity, we do not show the observational data. The top panel shows our fiducial results for the KiDS-450 mock in some tomographic redshift bins (the same as in Figure 7), and the lower panel shows the ratios of the GI and II terms with the GG term. As one can see from the plot, the II term is very weak, consistent with zero. This is expected because most galaxies in our model are late types, and their IA is very weak (Figure 4). Moreover, the GI term is basically positive, and its contribution could be as large as 15% on large scales. Following the procedures usually adopted in observational work to determine the free parameter ${A}_{\mathrm{IA}}$ (e.g., Heymans et al. 2013), we fit the total signals (black circles) from simulation using the nonlinear IA model (Equation (16)). Note that in our calculation, the GG signal is given using the nonlinear theoretical power spectrum with the given cosmological parameters, and the red dashed line in Figure 11 shows that the theoretical prediction agrees well with the measured GG term from our simulation. The best fit to the total signal (GG+II+GI) gives ${A}_{\mathrm{IA}}=-0.972\pm 0.217$. The fit to each component is also shown as the dashed line in Figure 11.

Zoom In Zoom Out Reset image size

Hildebrandt et al. (2017) constrained the amplitude of the IA to be positive,${A}_{\mathrm{IA}}=1.10\pm 0.64$, in their fiducial analysis of KiDS-450, which gives a negative GI term in the measurements. By contrast, a negative ${A}_{\mathrm{IA}}$ here indicates that the contribution from GI term is positive and acts to increase the overall correlation signals. In fact, a positive GI signal is not surprising and has been reported in major weak-lensing surveys. For example, Fu et al. (2008) found that ${A}_{\mathrm{IA}}=-{2.2}_{-3.8}^{+4.6}$ from third-year CFHTLenS data, and Heymans et al. (2013) reported that ${A}_{\mathrm{IA}}=-{1.18}_{-1.17}^{+0.96}$ from the final CFHTLenS data. Joudaki et al. (2017) found that ${A}_{\mathrm{IA}}=-3.6\pm 1.6$ from a re-analysis of CFHTLenS data. Hildebrandt et al. (2017) found a positive GI term with ${A}_{\mathrm{IA}}=-{1.10}_{-0.70}^{+0.96}$ for the KiDS-450 data when they used the BPZ method to estimate the galaxy photometric redshift. Recently, Troxel et al. (2017) found from the DES data that for spiral galaxies, the GI term is also positive with ${A}_{\mathrm{IA}}=-0.8$ at an 84% confidence level.

Since most galaxies in our model are late types, we conclude that the positive GI signal is contributed by spirals. As a further test of the contributions from spiral and elliptical galaxies, we show, in Figure 12, the ratios of the II and GI terms with the GG term for two ideal cases. In the top panel, we assume that all central galaxies in our model are spirals and that their spins follow the spins of their host dark matter halos. In the bottom panel, we assume that all central galaxies are ellipticals, and their shapes follow the shapes of the host halos defined with the inertia tensor, but with a misalignment given by a Gaussian distribution with a dispersion of 35^◦. In the two cases, all satellite galaxies are assumed to have random orientations. Figure 12 shows that the GI term is indeed positive on all scales for spiral galaxies, although their II contribution is close to zero. For elliptical galaxies, their II term is positive and the GI term is negative. Note here that the error bars are different for the two different ratios. This can be simply explained using Equation (15) by considering the noise in the correlation of the GI and II terms. As show by the GG correlation, gravitational shear can be accurately measured in the mock. So, if we consider the noise (${\boldsymbol{N}}$) in the measurement of the intrinsic shapes of galaxies, the GI and II terms can be expressed as $\langle {{\boldsymbol{\epsilon }}}_{i}^{({\rm{I}})}{{\boldsymbol{g}}}_{j}\rangle +\langle {{\boldsymbol{N}}}_{i}{{\boldsymbol{g}}}_{j}\rangle$ and $\langle {{\boldsymbol{\epsilon }}}_{i}^{({\rm{I}})}{{\boldsymbol{\epsilon }}}_{j}^{({\rm{I}})}\rangle +2\langle {{\boldsymbol{\epsilon }}}_{i}^{({\rm{I}})}{{\boldsymbol{N}}}_{j}\rangle +\langle {{\boldsymbol{N}}}_{i}{{\boldsymbol{N}}}_{j}\rangle$, respectively. The effect of shape noise can contribute to the additional correlations as $\langle {{\boldsymbol{N}}}_{i}{{\boldsymbol{g}}}_{j}\rangle$ for the GI term and $2\langle {{\boldsymbol{\epsilon }}}_{i}^{({\rm{I}})}{{\boldsymbol{N}}}_{j}\rangle +\langle {{\boldsymbol{N}}}_{i}{{\boldsymbol{N}}}_{j}\rangle$ for the II term. Combining with the definition of correlation ${\xi }_{+}$, this difference can be used to explain the different error bars for the two different ratios.

Zoom In Zoom Out Reset image size

A positive GI term from spiral galaxies is not expected from the tidal field model. From linear theory (e.g., Hirata & Seljak 2004), the GI term is found to be negative and has been used as a fiducial model in weak-lensing data analyses (e.g., Joachimi et al. 2011; Heymans et al. 2013). One important assumption in the linear model is that the shape of foreground galaxies is radially aligned with the nearby overdense region. This is on average true for elliptical galaxies. But for spiral galaxies, where alignments are mainly determined by the angular momenta of the dark matter halos through a large-scale tidal field (see Schäfer 2009 and references therein), this may not be true. Observationally, it is found that the spins of spiral galaxies tend to align with those of nearby filaments, but the short axes of the ellipticals are perpendicular to those of the filaments (Jones et al. 2010; Tempel & Libeskind 2013). Both hydrodynamical simulations (e.g., Codis et al. 2015) and N-body simulations (e.g., Kang & Wang 2015; Wang & Kang 2017) also confirmed such a dependence on galaxy type. In particular, Chisari et al. (2015) found from hydrodynamical simulations that spiral galaxies have a significant tendency to be tangentially aligned with overdensity regions. Their Figure 10 demonstrates clearly the alignment of spiral galaxies around overdensity regions and the origin of a positive GI term.

Finally, we note that the II and GI terms are closely related to how we model galaxy shapes and orientations. In this paper, we simply assumed that spiral galaxies follow the spins of dark matter halos. However, as shown in, e.g., Bett et al. (2010), galaxy spins have a broad distribution of misalignment with dark matter halos. This misalignment will reduce the positive GI terms. Furthermore, the total GI and II terms in real data also depend on the fraction of spiral and elliptical galaxies as well as on galaxy luminosities and redshifts. More comprehensive analyses of these factors are needed to quantify their impacts on the GI and II terms. This paper, which makes use of both N-body simulations and galaxies from a semi-analytical model, is a step toward this goal. But here, we only focus on a first comparison of the predicted shear correlations with the data. We will present a more comprehensive investigation on the contribution of GI and II terms in a future paper.

After decomposing the light cone into a set of shells with the width of $\sim 50\,{h}^{-1}\,\mathrm{Mpc}$, we can obtain the surface matter overdensity ${\sigma }^{(n)}$ of the nth shell through

$$\begin{eqnarray}&&{\sigma }^{(n)}({{\boldsymbol{\theta }}}^{(n)})={\int }_{{\chi }_{n-1/2}}^{{\chi }_{n+1/2}}d{\chi }^{{\rm{{\prime} }}}\delta (r({\chi }^{{\rm{{\prime} }}}){{\boldsymbol{\theta }}}^{(n)},{\chi }^{{\rm{{\prime} }}}).\end{eqnarray} \tag{ 33 }$$

The convergence field is then given by

$$\begin{eqnarray}&&{\kappa }^{(n)}({{\boldsymbol{\theta }}}^{(n)})={W}^{(n)}{\sigma }^{(n)}({{\boldsymbol{\theta }}}^{(n)}),\end{eqnarray} \tag{ 34 }$$

where the lensing kernel ${W}^{(n)}$ is defined as

$$\begin{eqnarray}&&{W}^{(n)}=\displaystyle \frac{3}{2}{\left(\displaystyle \frac{{H}_{0}}{c}\right)}^{2}{{\rm{\Omega }}}_{{\rm{m}}}\displaystyle \frac{r({\chi }_{n})}{a({\chi }_{n})}.\end{eqnarray} \tag{ 35 }$$

By solving the Poisson equation, one can obtain the lensing potential in harmonic space,

$$\begin{eqnarray}&&{\phi }_{{\ell }m}^{(n)}=-\displaystyle \frac{2}{{\ell }({\ell }+1)}{\kappa }_{{\ell }m}^{(n)}.\end{eqnarray} \tag{ 36 }$$

In the context of gravitational lensing, the deflection field ${\alpha }_{{\ell }m}^{(n)}$ can be derived from the gravitational lensing potential through Equation (18), and light rays can be propagated to the next shells following (Teyssier et al. 2009)

$$\begin{eqnarray}&&{{\boldsymbol{x}}}^{(n+1)}={ \mathcal R }({{\boldsymbol{n}}}^{(n)}\times {{\boldsymbol{\alpha }}}^{(n)},\parallel {{\boldsymbol{\alpha }}}^{(n)}\parallel ){{\boldsymbol{x}}}^{(n)},\end{eqnarray} \tag{ 37 }$$

where the rays are initialized at the center of each HEALPix cell. The rotation matrix ${ \mathcal R }$ defines the propagation direction between different shells. The lensing distortion matrix ${{ \mathcal A }}^{(n)}$ can be evaluated from (Becker 2013)

$$\begin{eqnarray}{{ \mathcal A }}_{{ij}}^{(n+1)} & = & \left(1-\displaystyle \frac{{D}_{0}^{n}}{{D}_{0}^{n+1}}\displaystyle \frac{{D}_{n-1}^{n+1}}{{D}_{n-1}^{n}}\right){{ \mathcal A }}_{{ij}}^{(n-1)}\\ & & +\displaystyle \frac{{D}_{0}^{n}}{{D}_{0}^{n+1}}\displaystyle \frac{{D}_{n-1}^{n+1}}{{D}_{n-1}^{n}}{{ \mathcal A }}_{{ij}}^{(n)}-\displaystyle \frac{{D}_{n}^{n+1}}{{D}_{0}^{n+1}}{{ \mathcal U }}_{{ik}}^{(n)}{{ \mathcal A }}_{{ij}}^{(n)},\end{eqnarray} \tag{ 38 }$$

where the angular diameter distance ${D}_{n}^{n+1}\equiv r({\chi }_{n+1}-{\chi }_{n});$ ${{ \mathcal U }}_{{ij}}^{(n)}$ is the tidal matrix of the nth shell, which can be related to the second derivatives of the lensing potential, ${{ \mathcal U }}_{{ij}}^{(n)}={\phi }_{,{ij}}^{(n)}$.

For more accurate computations, we perform the spherical harmonic analysis with an iterative algorithm, as is performed in the HEALPix subroutine map2alm_iterative, to control the residual in the solution of the lensing potential. The order of iteration in the analysis is chosen to be 3, as a compromise with the computational time. In Figure 13, we show the statistical measurements with or without this iterative method. Here we define the relative deviation as ${\rm{\Delta }}=[{C}_{0}({\ell })-{C}_{3}({\ell })]/{C}_{3}({\ell })$, where the subscript denotes the order of iteration. Clearly, via the iterative algorithm, the numerical errors are controlled effectively, and the power spectra of the shear B-mode are significantly suppressed, especially for scales smaller than $\sim 1\ \mathrm{arcmin}$. It is also more accurate than the multigrid method used in CALCLENS¹³ as discussed in Becker (2013).

Zoom In Zoom Out Reset image size