Precision calculations of the cosmic shear power spectrum projection

TL;DR

This paper evaluates the accuracy of various approximations in calculating the cosmic shear power spectrum, demonstrating their negligible impact on current surveys and sufficiency for future high-precision measurements, with publicly available code.

Contribution

It provides a comprehensive comparison of projection approximations for the cosmic shear power spectrum and offers validated methods for future surveys.

Findings

Approximations have negligible impact on current survey constraints.

Second-order Limber approximation is sufficient for future high-precision surveys.

Flat-sky Hankel transformation errors are below 2% compared to full spherical calculations.

Abstract

We compute the spherical-sky weak-lensing power spectrum of the shear and convergence. We discuss various approximations, such as flat-sky, and first- and second- order Limber equations for the projection. We find that the impact of adopting these approximations is negligible when constraining cosmological parameters from current weak lensing surveys. This is demonstrated using data from the Canada-France-Hawaii Telescope Lensing Survey (CFHTLenS). We find that the reported tension with Planck Cosmic Microwave Background (CMB) temperature anisotropy results cannot be alleviated. For future large-scale surveys with unprecedented precision, we show that the spherical second-order Limber approximation will provide sufficient accuracy. In this case, the cosmic-shear power spectrum is shown to be in agreement with the full projection at the sub-percent level for l > 3, with the corresponding errors an order of magnitude below cosmic variance for all l. When computing the two-point shear correlation function, we show that the flat-sky fast Hankel transformation results in errors below two percent compared to the full spherical transformation. In the spirit of reproducible research, our numerical implementation of all approximations and the full projection are publicly available within the package nicaea at http://www.cosmostat.org/software/nicaea.

Tables

Table 1. Table 1: The shear power spectrum approximations studied in this paper. ‘ID’ is the label used in the text and figures. The sixth (seventh) column indicates the mode ℓ x subscript ℓ 𝑥 \ell_{x} such that for ℓ ≤ ℓ x ℓ subscript ℓ 𝑥 \ell\leq\ell_{x} , the approximated power spectrum is more accurate than x 𝑥 x , with x = 0.1 𝑥 0.1 x=0.1 ( 0.01 0.01 0.01 ).

Case	ID	equation	$p\left(\ell \right)$	$\nu \left(\ell \right)$	${\ell }_{0.1}$	${\ell }_{0.01}$	Comment
$1^{\mathrm{st}}$-order standard Limber, flat	L1Fl	(43)+(i)	$1$	$\ell$	$4$	$60$	Pre-2014 CFHTLenS and DLS
$1^{\mathrm{st}}$-order extended Limber, flat	ExtL1Fl	(43)+(ii)	$\frac{{\ell }^4}{{\left(\ell +1/2\right)}^4}$	$\ell +\frac{1}{2}$	$20$	$200$	Converges only with $𝒪\left({\ell }^{-1}\right)$
$1^{\mathrm{st}}$-order extended Limber, hybrid, flat	ExtL1FlHyb	(43)+(iii)	$1$	$\ell +\frac{1}{2}$	$4$	$10$	Post-2014 CFHTLenS, DES-SV and KiDS
$1^{\mathrm{st}}$-order extended Limber, spherical	ExtL1Sph	(38)	$\frac{{\text{ł}}^2(\ell ,2)}{{\left(\ell +1/2\right)}^4}$	$\ell +\frac{1}{2}$	$3$	$12$
$2^{\mathrm{st}}$-order extended Limber, flat	ExtL2Fl	(43)+(44)+(iii)	$\frac{{\ell }^4}{{\left(\ell +1/2\right)}^4}$	$\ell +\frac{1}{2}$	$19$	$200$	Converges only with $𝒪\left({\ell }^{-1}\right)$
$2^{\mathrm{st}}$-order extended Limber, hybrid, flat	ExtL2FlHyb	(43)+(44)+(iii)	$1$	$\ell +\frac{1}{2}$	$5$	$16$	Best flat-sky approximation
$2^{\mathrm{st}}$-order extended Limber, spherical	ExtL2Sph	(38)+(42)	$\frac{{\text{ł}}^2(\ell ,2)}{{\left(\ell +1/2\right)}^4}$	$\ell +\frac{1}{2}$	$2$	$5$	Best approximation
Full spherical	FullSph	(21)	-	-	-	-	Correct projection

Table 2. Table 2: Mean and 68% credible interval for σ 8 ​ ( Ω m / 0.27 ) 0.6 subscript 𝜎 8 superscript subscript Ω m 0.27 0.6 \sigma_{8}(\Omega_{\rm m}/0.27)^{0.6} and σ 8 ​ ( Ω m / 0.3 ) 0.6 subscript 𝜎 8 superscript subscript Ω m 0.3 0.6 \sigma_{8}(\Omega_{\rm m}/0.3)^{0.6} for various approximations to the lensing power spectrum projections listed in Table  1 .

ID	${\sigma }_8{\left({\Omega }_m/0.27\right)}^{0.6}$	${\sigma }_8{\left({\Omega }_m/0.3\right)}^{0.6}$
L1Fl	${0.787}_{-0.033}^{+0.031}$	${0.739}_{-0.031}^{+0.029}$
ExtL1Fl	$0.792\pm 0.032$	$0.744\pm 0.030$
ExtL1FlHyb	${0.788}_{-0.033}^{+0.031}$	${0.740}_{-0.031}^{+0.029}$
ExtL2FlHyb	${0.788}_{-0.033}^{+0.031}$	${0.740}_{-0.031}^{+0.029}$
ExtL2Sph(Hankel)	${0.789}_{-0.032}^{+0.031}$	${0.740}_{-0.030}^{+0.029}$

Equations

$$\psi (\theta ,\phi )=\frac{2}{{\mathrm{c}}^2}\underset{0}{\overset{\infty }{\int }}\frac{\mathrm{d}\chi }{\chi }\Phi (\chi ,\chi \theta ,\chi \phi ;\chi )q(\chi ),$$

$$q(\chi )=\underset{\chi }{\overset{{\chi }_{\mathrm{h}}}{\int }}\mathrm{d}{\chi }^{\prime }n({\chi }^{\prime })\frac{{\chi }^{\prime }-\chi }{{\chi }^{\prime }},$$

$$\hat{\Phi }(\mathbf{k};\chi )=$$

$$\Phi (\mathbf{r};\chi )=$$

$$\psi (\theta ,\phi )=\sum _{\ell =0}^{\infty }\sum _{m=-\ell }^{\ell }{\psi }_{\ell m}{\mathrm{Y}}_{\ell m}(\theta ,\phi );$$

$${\psi }_{\ell m}={\int }_{{\mathbb{S}}^2}\mathrm{d}\Omega \psi (\theta ,\phi ){\mathrm{Y}}_{\ell m}^{\ast }(\theta ,\phi ).$$

$$⟨{\psi }_{\ell m,i}{\psi }_{{\ell }^{\prime }{m}^{\prime },j}^{\ast }⟩={\delta }_{\ell {\ell }^{\prime }}{\delta }_{m{m}^{\prime }}{C}_{ij}^{\psi }(\ell ).$$

$$C_{ij}^{\psi }(\ell )=$$

$$=$$

$$\mathcal{A}=\frac{3}{2}{\Omega }_{\mathrm{m}}{\left(\frac{H_0}{c}\right)}^2.$$

$$⟨\hat{\Phi }(\mathbf{k};\chi ){\hat{\Phi }}^{\ast }({\mathbf{k}}^{\prime };{\chi }^{\prime })⟩=$$

$$⟨\hat{\delta }(\mathbf{k};\chi ){\hat{\delta }}^{\ast }({\mathbf{k}}^{\prime };{\chi }^{\prime })⟩=$$

$$\gamma (\theta ,\phi )=$$

$${\gamma }^{\ast }(\theta ,\phi )=$$

$$({\gamma }_1\pm \mathrm{i}{\gamma }_2)(\theta ,\phi )=$$

$${}_2{\gamma }_{\ell m}=$$

$${}_{-2}{\gamma }_{\ell m}=$$

$$\text{l}(\ell ,s){}_s{\mathrm{Y}}_{\ell m}(\theta ,\phi )=$$

$$\text{l}(\ell ,s){}_{-s}{\mathrm{Y}}_{\ell m}(\theta ,\phi )=$$

$$\text{l}(\ell ,s)=\sqrt{\frac{(\ell +s)!}{(\ell -s)!}}.$$

$${}_{\pm 2}{\gamma }_{\ell m}=\frac{1}{2}\text{l}(\ell ,2){\psi }_{\ell m}.$$

$$⟨{}_2{\gamma }_{\ell m,i}{}_2{\gamma }_{{\ell }^{\prime }{m}^{\prime },j}^{\ast }⟩={\delta }_{\ell {\ell }^{\prime }}{\delta }_{m{m}^{\prime }}{C}_{ij}^{\gamma }(\ell ).$$

$$C_{ij}^{\gamma }(\ell )=\frac{1}{4}\text{l}{}^2(\ell ,2)C_{ij}^{\psi }(\ell )=$$

$$\kappa (\theta ,\phi )=\frac{1}{2}ðð{}^{\ast }\psi (\theta ,\phi )=\frac{1}{2}{\nabla }^2\psi (\theta ,\phi ).$$

$${\nabla }^2{\mathrm{Y}}_{\ell m}(\theta ,\phi )=-\ell (\ell +1){\mathrm{Y}}_{\ell m}(\theta ,\phi )=-\text{l}{}^2(\ell ,1){\mathrm{Y}}_{\ell m}(\theta ,\phi ).$$

$$C_{ij}^{\kappa }(\ell )=\frac{1}{4}\text{l}{}^4(\ell ,1)C_{ij}^{\psi }(\ell )=\frac{\ell (\ell +1)}{(\ell -1)(\ell +2)}{C}_{ij}^{\gamma }(\ell ).$$

$$\psi (\mathbf{\vartheta })=$$

$$\hat{\psi }(\mathbf{\ell })=$$

$$⟨{\hat{\psi }}_i(\mathbf{\ell }){\hat{\psi }}_j^{\ast }({\mathbf{\ell }}^{\prime })⟩=(2\pi {)}^2{\delta }_{\mathrm{D}}(\mathbf{\ell }-{\mathbf{\ell }}^{\prime })P_{ij}^{\psi }(\ell ).$$

$${\ell }^2{}_{\pm 2}{\mathrm{Y}}_{\ell m}(\theta ,\phi )\approx {\mathrm{e}}^{\mp 2\mathrm{i}{\varphi }_{\ell }}({\partial }_1\pm \mathrm{i}{\partial }_2{)}^2{\mathrm{Y}}_{\ell m}(\theta ,\phi ).$$

$$P_{ij}^{\gamma }(\ell )=$$

$${\mathrm{j}}_{\ell }(x)=\sqrt{\frac{\pi }{2x}}{\mathrm{J}}_{\ell +1/2}(x)$$

$${\int }_0^{\infty }\mathrm{d}\chi f(\chi ){\mathrm{J}}_{\nu }(k\chi )={\int }_0^{\infty }\mathrm{d}x{k}^{-1}f(x/k){\mathrm{J}}_{\nu }(x)$$

Full text

††pagerange: Precision calculations of the cosmic shear power spectrum projection–C††pubyear: 2017

1 Introduction

The measurement of weak gravitational lensing by large-scale structures provides a powerful cosmological probe of dark matter, dark energy, and modifications to gravity. As such it is the primary science goal of several current (KiDS, DES, HSC111KiDS: , DES: , HSC: missing) and future (Euclid, LSST, WFIRST222Euclid: , LSST: , WFIRST: missing) surveys. Interest in the results from these surveys is high as statistically significant deviations have been found between the cosmological parameter constraints from the CMB Planck experiment (Planck Collaboration et al. 2016) in comparison to weak lensing constraints from both the Kilo-Degree Survey (KiDS; Hildebrandt et al. 2017) and the Canada-France-Hawaii Telescope Lensing Survey (CFHTLenS; Joudaki et al. 2017). If the source of this tension is not a result of so-far unconsidered sources of systematic errors in one or all experiments, extensions to the standard flat $\Lambda$CDM cosmological models need to be considered. Joudaki et al. (2016) have shown, for example, that the tension can be resolved with an evolving dark energy model.

In the era of the upcoming large-scale surveys that will provide measurements of cosmic shear with unprecedented precision, one needs to revisit the theoretical predictions of the observables to ensure that the accuracy of the models meet the accuracy of the observations. In this paper, we examine the widely used Limber approximation for the projected weak-lensing power spectrum. We consider spherical coordinates and the flat-sky approximation, and compute the full projection of the lensing power spectrum. The first-order extended Limber approximation provides sub-percent accuracy for $\ell>10$ and is sufficient for present surveys. The associated errors are sub-dominant even for future large surveys.

We further show that the second-order Limber approximation is an accurate representation of the full projection, with better than percent level precision for scales $\ell>3$. Since this approximation involves only 1D integrals over the matter power spectrum, it is very fast to calculate numerically and can readily be employed in Monte-Carlo sampling methods to obtain precision constraints on cosmological parameters. We also compute the shear correlation function using a spherical transformation, and compare this to the flat-sky approximated, commonly used fast Hankel transformation.

This paper is organised as follows. In section 2 we provide a pedagogical introduction to weak gravitational lensing theory, projections and power spectra for the flat-sky and spherical cases, following the seminal work by Hu (2000) (see also Castro et al. 2005). In section 3 we derive weak-lensing observables using a second-order Limber approximation first introduced by Loverde & Afshordi (2008). We compare the shear power spectrum and the commonly-used two-point shear correlation function for the full solution to a range of different approximations in section 4, providing cosmological constraints for each case using CFHTLenS data from Kilbinger et al. (2013). This paper draws from several sources of literature which have previously discussed the full projection, or reviewed the accuracy of the Limber approximation for weak lensing studies, namely Schmidt (2008); Bernardeau et al. (2012); Giannantonio et al. (2012); Kitching et al. (2017), see also Lemos et al. (2017) for a more recent work. We present a discussion and comparison of our results to these papers in Appendix B.

2 Weak-lensing projections and power spectra

In this section we review the basic weak-lensing projection expressions, and compute lensing power spectra for a spherical case, and in the flat-sky approximation. For completeness we provide a derivation of the weak lensing power spectra in Appendix A.

2.1 The lensing potential

The lensing potential $\psi$ at a position on the sky $(\theta,\varphi)$ in the Born approximation is defined as the projected 3D metric potential $\Phi$ along the line of sight of a flat Universe (Kaiser 1998; Bartelmann & Schneider 2001),

[TABLE]

where the lensing efficiency $q$ is given by

[TABLE]

corresponding to a population of lensed galaxies with a normalised source redshift distribution $n_{z}(z){\rm d}z=n(\chi){\rm d}\chi$, with the limit being the comoving distance to the horizon $\chi_{\rm h}$333In (1) we have replaced this limit without loss of generality with $\infty$, since $q(\chi>\chi_{\rm h})=0$.. Here, ${\rm c}$ is the speed of light, and the projection is carried out over comoving distances $\chi$. The last argument of the potential $\Phi$ is not to be understood as coordinate, but rather as a substitute of cosmic time, $t(\chi)$, to express the time at which the potential is evaluated. This is true in the following for all fields and functions thereof that dynamically change with cosmic epoch.

The form of the lensing efficiency $q$ in equation (2) assumes a homogeneous galaxy distribution without clustering, so that the redshift distribution in this approximation does not depend on the direction on the sky. Accounting for this position-dependence leads to corrections of weak-lensing quantities due to clustering of source galaxies with other sources (Schneider et al. 2002), and with galaxies associated to lens structures (Bernardeau 1998; Hamana et al. 2002).

The 3D potential is related to the density contrast $\delta$ via the Poisson equation. Assuming General Relativity, this relation is written in Fourier space as

[TABLE]

where $\Omega_{\rm m}$ is the matter density parameter, $H_{0}$ the Hubble constant, $\bm{k}$ a 3D Fourier wave vector with modulus $k$ being the comoving wave number, and $a$ the scale factor with $a=1$ today. The Fourier transform of the potential and its inverse are defined as

[TABLE]

where the integration range for both integrals is $\mathbb{R}^{3}$.

2.2 Lensing power spectra in the spherical case

2.2.1 Lensing potential power spectrum

Following Hu (2000) we decompose the lensing potential (equation 1) into spherical harmonics, in analogy to the cosmic microwave background (CMB) temperature, both of which are scalar functions on the sphere. This decomposition and its inverse are

[TABLE]

Complex conjugation is denoted with the superscript ∗. To specify tomographic redshift bins $i=1\ldots N_{z}$, we introduce a family of lensing efficiency functions $q_{i}$ defined by a corresponding family of redshift distributions $n_{i}$ via equation (2). The resulting lensing potential is denoted by $\psi_{\ell m,i}$. The tomographic power spectrum of the lensing potential between two redshift bins $i$ and $j$, $C_{ij}^{\psi}$ (Peebles 1980) is then defined as

[TABLE]

Note that the argument $\ell$ is a discrete integer variable, and is often written as index, $C_{\ell}$. Using the properties of the spherical harmonics (see App. A for details) the power spectrum can be written as

[TABLE]

where ${\rm j}_{\ell}$ is the spherical Bessel function of order $\ell$. For convenience we introduce the normalisation constant ${\cal A}$ as

[TABLE]

The first line expresses $C_{ij}^{\psi}$ in terms of the 3D potential power spectrum $P_{\Phi}$, defined as

[TABLE]

The second line uses the 3D matter density power spectrum $P_{\rm m}$, which is defined analogously as

[TABLE]

and is related to $P_{\Phi}$ via the absolute square of the Poisson equation (3).

This type of cross-power spectrum between different cosmological epochs $\chi$ and $\chi^{\prime}$ was introduced in Castro et al. (2005). In Sects. 3.1.1 and 3.1.2 this unequal-time cross-spectrum (Kitching & Heavens 2017) will be further evaluated and simplified in the context of the Limber approximation. The oscillating Bessel functions in equation (10) ensure that only relatively close epochs contribute to the lensing potential correlation. This make sense since observed light rays from two galaxies at different positions on the sky that necessarily converge at the observer today, pick up the density fluctuations at similar times while propagating through the large-scale structure. A similar argument has been made in Bartelmann & Schneider (2001): since the matter power spectrum scales with $k$ for $k\rightarrow 0$, there is decreasing power towards larger and larger scales. In particular, the correlation of cosmic fields decreases strongly above a coherence scale |\chi-\chi^{\prime}|\begin{array}[]{c}>\\ \sim\end{array}L_{\rm coh} which is significantly smaller than the horizon scale $\chi_{\rm h}$.

In the following section we will discuss the relations between shear and convergence to the lensing potential on the sphere, and derive the power spectrum of the former two fields.

2.2.2 Shear power spectrum

The shear $\gamma=\gamma_{1}+{\rm i}\gamma_{2}$ is related to the potential at linear order by the trace-free part of the Jacobi matrix. The involved differential operator on the sphere is called edth derivative, $\,\eth\,$, see Castro et al. (2005) for an in-depth mathematical discussion of this concept. The edth operator $\,\eth\,$ ($\,\eth\,^{\ast}$) raises (lowers) the spin $s$ of an object. Twice applying this operator to the scalar (spin-0) potential creates the spin-2 shear:

[TABLE]

To write the shear on the sphere in terms of the lensing potential $\psi$, we insert the harmonics expansion of the potential (6). This requires the calculation of second derivatives of the spherical harmonic functions. This operation defines a new object, the spin-weighted spherical harmonic ${}_{s}{\rm Y}_{\ell m}$. The shear can be written on the sphere in terms of these functions as a spherical harmonics multipole expansion with coefficients ${}_{\pm 2}\gamma_{\ell m}$. This expansion together with its inverse is

[TABLE]

The spin-weighted spherical harmonics ${}_{s}\!{\rm Y}_{\ell m}$ that are the basis function in the expansion of the shear (equation 15) can be calculated via the relations

[TABLE]

for $0\leq s\leq\ell$, with the spin pre-factor (Bernardeau et al. 2012)

[TABLE]

Inserting the lensing potential expansion (equation 6) into the expression for the shear (equation 14), and using equation (17) to compute the derivatives, we find for the shear expansion coefficients (Hu 2000; Taylor 2001)

[TABLE]

The two coefficients ${}_{+2}\gamma_{\ell m}$ and ${}_{-2}\gamma_{\ell m}$ are identical since the potential $\psi$ is a real function.

The tomographic shear power spectrum, in analogy to equation (10), is defined by

[TABLE]

This is given by

[TABLE]

where we use the label ‘FullSph’, see Table 1 for a list of cases discussed in this work. The spherical spin pre-factor for the full projection shear power spectrum is $\hbox{\it\l}\,^{2}(\ell,2)$, which will be modified under flat-sky and Limber approximations below.

2.2.3 Convergence power spectrum

The convergence is related to the lensing potential on the sphere via the product of spin-raising and spin-lowering edth operators, which are identical to the spherical Laplacian differential operator.

[TABLE]

The spherical harmonics are eigenfunctions of the Laplacian,

[TABLE]

The convergence power spectrum is then similar to the shear power spectrum (equation 20) with a different spherical pre-factor, (Hu 2000; Joudaki & Kaplinghat 2012)

[TABLE]

The convergence power spectrum is thus larger than the shear power spectrum, by 10% for $\ell=4$, 1% for $\ell=14$, and less than 0.1% for $\ell>45$.

2.3 Flat-sky approximation

The majority of cosmic shear analyses have used the predicted lensing power spectrum approximated on a flat sky, neglecting the sky curvature. This is a valid approach for past and current survey areas with an extent less than 10 degrees. To account for the sky curvature of the observed data, the shear correlation functions from observed galaxy ellipticities are now routinely computed using spherical coordinates, since projecting to a Cartesian plane has been shown to cause significant biases of the two-point funtion on large scales (Fu et al. 2008), and lead to spurious B-modes (Asgari et al. 2017). Here we examine the effect of sky curvature on the theoretical models of the power spectrum, and the effect on cosmological parameter inference (see Sect. 4.3).

For a flat-sky, the spherical harmonic expansions are approximated by Fourier transforms. The flat-sky equivalent of equations (6) and (7) are

[TABLE]

where $\bm{\vartheta}=(\theta,\varphi)$ is the vector describing a 2D angle on the sky. Instead of a harmonics coefficients $\psi_{\ell m}$, the Fourier representation of the potential $\hat{\psi}$ now depends on the vector $\bm{\ell}\in\mathbb{R}^{2}$.

The flat-sky power spectrum, i.e. the flat-sky analogue of equation (8), is defined by

[TABLE]

Hu (2000) has shown that for small angles the harmonics expansion (equation 6) can be approximated by the Fourier representation in equation (25). They also demonstrated that the power spectra are approximately equal, $P^{\psi}\approx C^{\psi}$.

For a spin-2 field, Hu (2000) approximates the edth operator by Cartesian derivatives, and approximates equation (17) as

[TABLE]

The spin pre-factor $\hbox{\it\l}\,(\ell,2)=\sqrt{(\ell-1)\ell(\ell+1)(\ell+2)}$ is replaced by $\ell^{2}$ in flat co-ordinates, an approximation that holds for large $\ell$, since sky curvature can be neglected for small angular scales. We find for the flat-sky shear power spectrum

[TABLE]

with flat-sky pre-factor $\ell^{4}$. See App. B.4 and B.5 for discussions of alternative expressions for the flat-sky power spectrum.

3 Second-order Limber approximation for weak lensing

3.1 Spherical case

We follow Loverde & Afshordi (2008) who derive the second-order Limber expansion for general projections from 3D to 2D scalar fields in the spherical, all-sky case. We apply their general derivation to the case of weak lensing, and contrary to Loverde & Afshordi (2008) account for a time-dependent power spectrum using two approaches presented in Sects. 3.1.1 and 3.1.2.

First, we use the identity of Bessel functions

[TABLE]

in equation (21), where ${\rm J}_{\nu}$ is the Bessel function of the first kind and order $\nu$. Next, Loverde & Afshordi (2008) solve integrals of the form

[TABLE]

by performing a Taylor expansion of an arbitrary differentiable function $f$ around $x=k\chi=\nu=\ell+1/2$, where the Bessel function has its approximate maximum.

3.1.1 Geometric mean cross-correlation power spectrum

To separate the $k$- and $\chi,\chi^{\prime}$-terms in equation (21), we first approximate the matter power cross-spectrum between two distances by the geometric mean of the two involved distances (Castro et al. 2005; Kitching & Heavens 2017),

[TABLE]

This form is justified when considering the linear power spectrum, and follows when inserting the linearly evolving density contrast $\hat{\delta}(\bm{k};\chi)=D_{+}(\chi)\hat{\delta}_{0}(\bm{k})$ into equation (13), where $\delta_{0}$ is the present-day linearly extrapolated density contrast, and $D_{+}$ the linear growth factor with $D_{+}(0)=1$. This is a good approximation also in the non-linear case as shown in Kitching & Heavens (2017).

With this, equation (21) is written as

[TABLE]

Note that this equation has a pre-factor $\hbox{\it\l}\,^{2}(\ell,2)$ corresponding to a spin-2 field, in contrast to Loverde & Afshordi (2008) who show calculations for a scalar field.

Following Loverde & Afshordi (2008) we expand to third order

[TABLE]

with $g(x)=k^{-1}f(k,\chi)$, $\chi=x/k$, and its derivatives $g^{(n)}(x)=k^{-1-n}f^{(n)}(k,\chi)$ for a given $k$, where the derivatives of $f$ are with respect to the second argument $\chi$. In this series expansion, the replacement of the integral with the evaluation of $g$ and its derivatives at the maximum of the Bessel function is a good approximation of the integral if $g$ is varying more slowly than the oscillating Bessel function. As we will show below, $f$ is a slowly varying function of the comoving distance. In our case the projection kernel is

[TABLE]

In the tomographic case, the index $i$ is added to $q$ and $f$. Replacing both distance integrals in equation (33) by their Taylor-expansions around the maxima $\nu(\ell)=\ell+1/2$ of the two Bessel functions, which are $k\chi$ and $k\chi^{\prime}$, respectively, yields

[TABLE]

Changing the integration to $\chi=\nu/k$ and collecting terms according to their $\nu$-dependence:

[TABLE]

The first term corresponds to the well-known first-order Limber approximation (Limber 1953; Kaiser 1992), which is widely used in weak gravitational lensing. We retrieve the (spherical) standard expression by inserting back the projection kernel (equation 35),

[TABLE]

In the Limber approximation, modes between structures at different epochs do not contribute to the single line-of-sight integration.

The second-order Limber term in equation (37) has an additional $\nu^{-2}$-dependence, and is therefore strongly suppressed for large $\ell$,

[TABLE]

The higher-order derivatives of the filter functions have to be computed numerically in the general case. These suffer from numerical noise and are sensitive to the set-up, for example the step size. The tabulation and interpolation of those derivatives is time-consuming since they depend on two arguments, $\nu$ and $\chi$. In the following section, we will separate the $k$- and $\chi$-dependent parts of the power spectrum to make the numerical derivatives faster and more smooth.

3.1.2 Approximation for small $\ell$

To further develop equation (32), we divide out the growth factor of the power spectrum,

[TABLE]

The function $\tilde{P}_{\rm m}$ is in general not time-independent, except in the case of a linear matter power spectrum in the absence of massive neutrinos, and in General Relativity. However, the second-order Limber terms are expected to be important only for small $\ell$, since for large $\ell$ the first-order Limber equation (38) is dominant. Eq. (40) becomes a good approximation for small $\ell$, since that means either small $k$ where the evolution is linear, or small $\chi$, where the lensing efficiency is small. The accuracy of a very nearby tomographic bin with low mean redshift should be further examined, but this is not the case for CFHTLenS.

With the definition (40), we factor out the function $\tilde{P}_{\rm m}(k,\chi)$ from equation (35), and we define the separated kernel function $f_{\rm s}$ as

[TABLE]

The Limber equation up to second order of the shear power spectrum can then be approximated as

[TABLE]

We compute the numerical higher derivatives as follows. The functions $f_{\rm s}(\chi)$ are fitted as power laws with index $\approx-1.5$, which is expected if the growth suppression factor $D_{+}(a)/a$ varies slowly with $\chi$, and the lensing efficiency $q\approx 1$ for small and medium $\chi$, given the CFHTLenS redshift range. The fit is carried out between $\chi_{\rm min}=0.001\,{\rm Mpc}/h$ and $\chi_{\rm max}=500\,\mbox{Mpc}/h$. At larger comoving distances the kernel decreases faster than a power law, so we exclude this range from the fit. Even though on those scales the derivatives are larger due to the steeper decline, the associated errors are very small as these scales are down-weighed by the kernel function $f_{\rm s}$ itself. At $\chi=500\,\mbox{Mpc}/h$ the filter function is less than $10^{-4}$ of its value at $1$ Mpc$/h$.

3.2 Flat-sky

The extended flat-sky Limber approximation is readily derived from the spherical case equations (38, 42), by replacing the pre-factor $\hbox{\it\l}\,^{2}(\ell,2)$ with $\ell^{4}$,

[TABLE]

Further approximations can be made for the pre-factor $p(\ell)=\ell^{4}/\nu^{4}(\ell)$ and $\nu(\ell)$:

$p(\ell)=1$, this corresponds to $\nu(\ell)=\ell$, which is the standard Limber approximation (L1Fl) Until recently, i.e. for all pre-2014 CFHTLenS results and DLS (Jee et al. 2013) analyses, this was the approximation of choice. Note that we do not discuss the second-order Limber approximation with $p(\ell)=1$. 2. 2.

$p(\ell)=\ell^{4}/(\ell+1/2)^{4}$. This corresponds to the extended Limber approximation (ExtL1Fl, ExtL2Fl) with $\nu(\ell)=\ell+1/2$; however the following case is typically employed: 3. 3.

$p(\ell)=1$, but keeping as argument of the power spectrum $\nu(\ell)=\ell+1/2$. This is a hybrid between standard and extended Limber approximation (ExtL1FlHyb, ExtL2FlHyb), and the first-order case was used in Hildebrandt et al. (2017); Joudaki et al. (2017, 2016); Abbott et al. (2016). As is shown below, this is a better approximation to the full projection than case (ii). In equation (44) the second-order suppression factor is also left to be $\nu^{-2}(\ell)=(\ell+1/2)^{-2}$, providing a slightly more accurate approximation compared to $\nu^{-2}(\ell)=\ell^{-2}$.

4 Results

4.1 Comparison of the approximations for the lensing power spectrum

In Fig. 1 we present the full spherical projection of the shear power spectrum in comparison to shear power spectra derived assuming the range of different approximations listed in Table 1. The adopted redshift distribution corresponds to CFHTLenS (Kilbinger et al. 2013) and we assume their best-fit flat $\Lambda$CDM cosmology with $\Omega_{\rm m}=0.279$, $\Omega_{\rm b}=0.046$, $\sigma_{8}=0.79$, $h=0.701$, $n_{\rm s}=0.96$. For $\ell>100$ we find that all shear power spectra predictions agree with the full spherical solution to better than one percent, with the majority of the approximations tested accurate to better than 0.1 per cent.

Considering first the flat-sky cases, the standard first-order Limber approximation, (L1Fl), that was adopted for all pre-2014 CFHTLenS analyses and DLS analyses, we find it to be accurate to better than 10% for $\ell>3$, converging slowly to the true projection with percent level precison at $\ell>100$. For the extended Limber approximations ‘hybrid’ cases (ExtL1FlHyb and ExtL2FlHyb), despite decreased accuracy for $\ell<8$ in comparison to the standard first-order Limber case, the errors with respect to the true power spectrum decrease much faster, as $\ell^{-2}$, such that percent-level precision is reached at $\ell>15$. The first-order extended Limber approximation ‘hybrid’ cases (ExtL1FlHyb) was adopted by Joudaki et al. (2017), Joudaki et al. (2016), DES-SV (Abbott et al. 2016) and Hildebrandt et al. (2017)444We confirm that there is a typographical error in equation 4 of Hildebrandt et al. (2017) and in equation 2 of Abbott et al. (2016) (private communication with Joe Zuntz) which should include the extra term of ‘+0.5’ in $\nu(\ell)$ that was incorporated in both cosmological analyses..

The outlier in the flat-sky cases is the extended Limber approximation (ExtL1Fl) which performs relatively poorly, and reaches 10% precision only at $\ell>100$. The same slow convergence can be observed for the corresponding second-order flat case, ExtL2Fl. To our knowledge this form of the flat-sky approximation has not been used in any cosmic shear studies to date, and should not be used in any future studies given these results. The poor behaviour of this case, in comparison with the ‘hybrid’ case, (e.g. ExtL1FlHyb) can be understood by considering Taylor expansions of the different pre-factors. The spherical pre-factor,

[TABLE]

can be compared to the flat-sky extended Limber pre-factor

[TABLE]

showing it to be more deviant from the full spherical solution, than the ‘hybrid’ $p(\ell)=1$ case.

Considering now the spherical-sky cases, we find that including the spherical pre-factor decreases the difference between the extended Limber approximated cases (ExtL1Sph and ExtL2Sph) and the full spherical solution by a factor of a few for $\ell<5$. We find that using the spherical-sky second-order extended Limber approximation (ExtL2Sph) yields percent accuracy down to $\ell>3$. The numerical calculation of the second-order extended Limber approximation is a factor of 15 times faster than the calculation of the full spherical solution (averaged over the first $18$ $\ell$-modes). We note that the sub-$0.1\%$-fluctuations seen in the right panel of Fig. 1 is due to numerical noise arising from numerical integration errors in the calculation of the full spherical solution when $\ell>20$.

We note that in all cases the second-order Limber expansion adds power to the first-order term. In the flat hybrid case, which over-estimates the full spherical solution, this results in the second-order expansion being less precise compared to first-order.

Compared to the statistical power of current and future surveys, all approximations discussed here are subdominant to the cosmic variance, $\Delta C(\ell)/C(\ell)=[f_{\rm sky}(2\ell+1)]^{-1/2}$ (Kaiser 1992), where $f_{\rm sky}$ is the fraction of sky observed by the survey. The uncertainties from the Limber approximation in the case of ExtL2Sph is an order of magnitude below the cosmic variance of a Euclid-like survey (sky area $15,000$ square degrees) for all $\ell$.

4.2 Effects on the shear correlation function

The majority of cosmic shear analyses to date have adopted real-space correlation statistics, since these can be measured directly from an observed galaxy shape catalogue. The baseline quantity is the two-point correlation function (Miralda-Escude 1991; Kaiser 1992; Bartelmann & Schneider 2001), given in the flat-sky approximation by

[TABLE]

The flat-sky shear power spectrum $P^{\gamma}$ can be related to the underlying matter power spectrum through equation (43) when adopting a first-order extended Limber approximation, or equation (44) when adopting a second-order extended Limber approximation.

On a sphere, correlation functions formally cannot be related to the power spectrum by the Hankel transform in equation (47), and should be replaced by the spherical transform (Ng & Liu 1999; Chon et al. 2004)

[TABLE]

where ${\rm d}^{\ell}_{m\,n}$ are the reduced Wigner ${\rm D}$-matrices, see App. C for details on their numerical calculation.

The spherical power spectrum is formally not defined for non-integer $\ell$ (see Castro et al. 2005, for alternative spherical-sky formulae for the two-point correlation function), as functions defined on the sphere are necessarily periodic. As we have shown in section 4.1, however, the spherical second-order extended Limber approximation provides a percent-level precision representation of the full spherical projection for $\ell>3$. The spherical pre-factor $\hbox{\it\l}\,(\ell,2)$ (equation 18) can be generalised to non-integer arguments and is positive for $\ell\geq 2$. We can thus use the spherical power spectrum with the Hankel transformation in equation (47) to compute the two-point correlation functions. This has the advantage being able to employ fast FFT numerical implementations of the Hankel tranforms (Hamilton 2000), when Monte-Carlo sampling.

We compare predictions for the two-point shear correlation function using the Hankel transformation and the full spherical transformation in Fig. 2. We show the full projection and the best approximation, ExtL2Sph. Note that for the ‘FullSph’ case we replace the full projection with ExtL2Sph for $\ell>200$ to reduce computation time. We find that the Hankel transform (equation 47) is accurate to better than $5$ ($0.2$) percent for $\xi_{+}$ ($\xi_{-}$). The difference between the second-order Limber and full projection solution using the spherical transform (equation 48) is well below one ($0.03$) percent on scales of $\vartheta<300$ arcminutes for $\xi_{+}$ ($\xi_{-}$). The red lines in Fig. 2 present the limit of precision that we can achieve with the current fast Hankel transform implementation for the correlation function. The green lines show the limit of the second-order Limber approximation on the correlation function.

Figure 3 shows the two-point correlation functions $\xi_{+}$ (left) and $\xi_{-}$ (right) using the different cases for the shear power spectrum listed in Table 1. The component $\xi_{+}$ is calculated with the appropriate transformation, i.e. Hankel for the flat cases, and spherical involving Legendre polynomials for the spherical cases. For $\xi_{-}$ which is less dominated by large scales and thus Limber and flat-sky approximations, in comparison to $\xi_{+}$, we use the Hankel transform in all cases. In addition, for $\xi_{-}$ the approximation ‘ExtL2Sph’ is our reference case. The adopted CFHTLenS redshift distribution and fiducial cosmological model are the same as in Figure 1. As is clear by the red dotted curve, using the Planck cosmology (Planck Collaboration et al. 2016) induces a significant change in the amplitude of the shear correlation function in comparison to the different projection methods.

4.3 Application to CFHTLenS data

The Canada-France-Hawaii Telescope Lensing Survey (CFHTLenS) represented a major step forward for the field of weak gravitational lensing, in terms of improved accuracy in data reduction (Erben et al. 2013), the implementation of PSF-Gaussianised matched multi-band photometry (Hildebrandt et al. 2012), cross-correlation clustering analysis between photometric redshift slices to verify tomographic redshift distributions (Benjamin et al. 2013), accurate calibrated shape measurements (Miller et al. 2013) and a full suite of informative systematic tests to select a clean data set (Heymans et al. 2012). Since the public release of this survey in 2013, the community has continued to scrutinise and advance our understanding of CFHTLenS by identifying a number of areas where analyses could improve:

• •

Choi et al. (2016) identified significant biases in the tomographic photometric redshift distributions using a more effective clustering analysis, in comparison to Benjamin et al. (2013), by incorporating newly overlapping spectroscopic data from the Sloan Digital Sky Survey. The conclusion of this work was that any re-analysis of CFHTLenS should include systematic error terms to account for bias and scatter, with a prediction that accounting for these biases would reduce the recovered amplitude of $\sigma_{8}$ by $\sim 4$%. Additional new techniques to calibrate the redshift distribution of tomographic bins was introduced recently in Hildebrandt et al. (2017).

• •

The CFHTLenS tomographic cosmological analysis was then revisited by Joudaki et al. (2017) in order to include a full redshift error analysis based on the results from Choi et al. (2016). The impact of correcting for these biases, including their associated errors, served to reduce the overall constraining power of the survey and hence also the tension between CFHTLenS and CMB constraints.

• •

Asgari et al. (2017) used the stringent COSEBI statistic (Schneider et al. 2010) to identify significant non-lensing B-mode distortions when the CFHTLenS data was split into tomographic slices.

• •

Kuijken et al. (2015) showed that the CFHTLenS shear calibration corrections derived in Miller et al. (2013) were underestimated as a result of an imperfect match between the galaxy populations in the data and image simulations.

• •

Fenech Conti et al. (2017) demonstrated that the CFHTLenS data would have been subject to a weight bias that favours galaxies that are more intrinsically oriented with the point-spread function. They also showed that the impact of calibration selection biases, that were not considered in Miller et al. (2013), would have lead to the over-correction of multiplicative shear bias in the CFHTLenS analyses, by a few percent.

• •

Joudaki et al. (2017) updated the CFHTLenS covariance matrices using larger-box numerical simulations that were less subject to the lack of power on large scales. A complementary accurate estimate of the covariance matrix using analytical methods will be published soon (Joachimi et al. in prep.)

• •

Takahashi et al. (2012) provided a more accurate non-linear power spectrum correction than that used in the original CFHTLenS analyses, and the halo model from Mead et al. (2015) allowed for simultaneous modelling of baryonic modifications to the non-linear power spectrum.

All these advances in our understanding were incorporated and accounted for in the recent KiDS cosmic shear analysis (Hildebrandt et al. 2017) which reports a $2.3\sigma$ tension with Planck. Efforts are now underway to fully re-analyse CFHTLenS using the advanced KiDS analysis pipeline with revised shape measurements and calibrations for the shear and photometric redshifts. Until this analysis is complete we note that these known shortcomings with the original CFHTLenS results impact in different ways the cosmological conclusions that one can draw. As CFHTLenS has similar statistical power to current weak lensing surveys, however, it nevertheless provides a very useful testbed with which to demonstrate the impact of adopting different approximations when constraining cosmological parameters.

In this work, we focus on the weak-lensing power spectrum projection, and assess the impact of various approximations on cosmological constraints from CFHTLenS. For consistency with the original analysis (Kilbinger et al. 2013), we adopt the same priors and non-linear power spectrum corrections from Smith et al. (2003).

We re-analyse the 2D CFHTLenS measurement of the two-point shear correlation function $\xi_{\pm}(\theta)$ from Kilbinger et al. (2013), defined in equation (47). As in Kilbinger et al. (2013) we fit both components $\xi_{+}$ and $\xi_{-}$ between angular scales $\theta=0.8$ and $350$ arc minutes, and use a $N$-body simulation estimate of the non-Gaussian covariance including the cross-covariance between both components. Bayesian Population Monte-Carlo parameter sampling is performed using the publicly available software CosmoPMC555http://www.cosmostat.org/software/cosmopmc (Wraith et al. 2009; Kilbinger et al. 2010). The cosmological modelling part includes the various lensing projections, calculated using the software library nicaea666http://www.cosmostat.org/software/nicaea.

For a first-order standard Limber flat-sky approximation (L1Fl) we find $\sigma_{8}(\Omega_{\rm m}/0.27)^{0.6}=0.787^{+0.031}_{-0.033}$, the same result that was published in Kilbinger et al. (2013). Using the second-order extended Limber flat-sky hybrid approximation (ExtL2FlHyb) results in $\sigma_{8}(\Omega_{\rm m}/0.27)^{0.6}=0.788\pm 0.032$, a negligible change of the amplitude that is well within the Monte-Carlo sampling noise. The largest difference is measured with the deprecated case ExtL1Fl, for which the recovered amplitude is larger by $16\%$ of the statistical error. These negligible changes to the error bars were to be expected owing to the high level of statistical noise and cosmic variance in comparison to the low-level impact of the various approximations shown in Fig. 1.

Table 2 lists the mean and 68% credible interval for $\sigma_{8}\Omega_{\rm m}^{0.6}$ for the various approximations to the lensing power spectrum projections listed in Table 1. Note again that these values do not represent the state-of-the-art cosmological results, since many of the above listed analysis advancements made since 2013 have not been taken into account. As an example of a significant effect, when using the revised non-linear power spectrum of Takahashi et al. (2012) in place of Smith et al. (2003), there is a decrease of $0.6\sigma$ with $\sigma_{8}(\Omega_{\rm m}/0.27)^{0.6}=0.768^{+0.029}_{-0.031}$, using L1Fl.

Considering the cosmological constraints from tomographic Kilo-Degree Survey (KiDS), we conclude that these are robust to flat-sky and Limber approximations. The case ExtL1FlHyb that was used for the analysis of KiDS data in Hildebrandt et al. (2017) and Joudaki et al. (2016) introduces errors that are more than an order of magnitude lower than the cosmic variance for that survey, and thus this approximation has a negligible impact on the cosmological parameters.

4.4 Alternative two-point shear statistics; the mass aperture statistic

and COSEBIs

The two-point shear correlation function $\xi_{\pm}$ represents the current baseline observable for cosmic shear measurements. As shown in Figure 3, however, using the standard first-order extended Limber flat-sky approximation (ExtL1FlHyb) can result in errors exceeding 10 percent, on angular scales $\theta>300$ arcmin. This is a result of the weight given to low $\ell$ modes in the $\xi_{+}$ statistic, as illustrated in Figure 4 which shows the integrand of $\xi_{+}$ and $\xi_{-}$ (upper two panels) for two cases ( $\theta=100$ and $\theta=350$ arcmin), normalised to their maximum value. This error does not impact CFHTLenS analyses, given the low signal-to-noise of the measurements on these scales. It will however become increasingly important for upcoming wider-field surveys that will accurately probe these scales.

In this paper we provide a solution in the form of the second-order extended Limber approximation, but another option to consider is the use of alternative two-point shear statistics that are less sensitive to accuracy in shear power spectrum measurement at low $\ell$. Both the aperture-mass dispersion, $\langle M_{\rm ap}^{2}\rangle$ (Schneider et al. 1998), and the Complete Orthogonal Sets of E/B-mode Integrals (COSEBIs), $E_{n}$ (Schneider et al. 2010) statistics satisfy this requirement and are linearly related to the shear power spectrum in the flat-sky approximation via integrals of the form

[TABLE]

where the Fourier-space filter functions $\hat{U}$ and $W_{n}$ are defined in Schneider et al. (1998) and Schneider et al. (2010), respectively. Figure 4 shows the integrands of these statistics, again normalised to their maximum value, where the integrands are of the form $\ell F(\ell)P^{\gamma}(\ell)$. The lower middle panel in Figure 4 shows the COSEBIs integrands for two angular ranges, $[1^{\prime},100^{\prime}]$ and $[0.8^{\prime},350^{\prime}]$, where we only show the integrands for the lowest COSEBIs mode, $E_{1}$, as the higher modes generally probe larger $\ell$-modes. The lowest panel shows the integrands of aperture mass dispersion statistics, for the same two maximum angular ranges.

Note that the development of the aperture-mass dispersion statistic, $\langle M_{\rm ap}^{2}\rangle$ was initially motivated to enable the separation of the measured signal into an E-mode (cosmological signal) and B-mode (systematics). This statistic is, however, a lossy conversion and is biased by small angular separations, where blending of galaxies makes shear measurement challenging (Kilbinger et al. 2006). The COSEBIs statistic tackles both these shortcomings. Kilbinger et al. (2013) present a detailed comparison of cosmological constraints obtained from this range of different two-point shear statistics finding consistent results.

5 Conclusions

In this paper we evaluate precision theoretical calculations for cosmic shear observables, bringing together different sources from the literature to provide a pedagogical review of the impact of adopting flat-sky and Limber approximations. We demonstrate that for current surveys, such as CFHTLenS and KiDS, these approximations have a negligible impact on cosmological parameter constraints.

For future surveys, the decrease in statistical errors places higher requirements on the accuracy of the theoretical modelling. There is also, however, the need to be able to rapidly sample the multi-dimensional cosmological parameter likelihoods. This requirement for computational speed is incompatible with a theoretical analysis that calculates a full spherical solution for the shear power spectrum, without adopting any approximation. We therefore present alternative solutions, revisiting the work of Bernardeau et al. (2012) who showed that adopting the second-order extended Limber approximation of Loverde & Afshordi (2008) provides a representation of the full spherical solution for the shear power spectrum that is accurate at the sub-percent level for $\ell>3$. We have verified this result and provide to the community our fast numerical implementation of all the approximations studied in this analysis, and the slow calculation of the full projection within the publicly available package nicaea at http://www.cosmostat.org/software/nicaea.

Finally we propose that future surveys seek to optimise the statistical analyses of their cosmic shear data. For example moving from the standard two-point shear correlation function statistic to the more stringent ‘COSEBI’ statistic (Schneider et al. 2010) renders the cosmic shear measurement insensitive to the low-$\ell$ scales where the Limber and flat-sky approximations have an impact on the precision of the theoretical modelling.

We have considered a flat Universe throughout this paper. To generalise the calculations to non-flat models, one needs to modify the comoving angular diameter distance to account for the spatial curvature $K\neq 0$. In addition, the spherical Bessel functions are replaced with hyperspherical Bessel functions (Abbott & Schaefer 1986). In a universe with positive curvature $K>0$ the 3D wave modes become discrete integer variables. For non-flat models, we do however not expect qualitative differences from our results.

Acknowledgments

The authors thank Tommaso Giannantonio, Sarah Bridle, Ami Choi, Donnacha Kirk, Lance Miller, Benjamin Joachimi, Chris Blake, Joe Zuntz, Joanne Cohn, Alex Hall and Adam Amara for very helpful discussions. This work was financially supported by the DFG (Emmy Noether grant Hi 1495/2-1; SI 1769/1-1; TR33 ‘The Dark Universe’), the Alexander von Humboldt Foundation, the ERC (grants 279396, 647112), the Seventh Framework Programme of the European Commission (Marie Sklodwoska Curie Fellowship grant 656869), the Netherlands Organisation for Scientific Research (NWO) (grants 614.001.103), NSERC, CIfAR, and the World Premier International Research Center Initiative (WPI), MEXT, Japan. Parts of this research were conducted by the Australian Research Council Centre of Excellence for All-sky Astrophysics (CAASTRO), through project number CE110001020i, and by the German Federal Ministry for Economic Affairs and Energy (BMWi) under project no. 50QE1103.

Appendix A Derivation of the weak-lensing power spectra

The following derivations are detailed in Hu (2000) and Castro et al. (2005), and are provided here for completeness.

A.1 Spherical case

A.1.1 Lensing potential power spectrum

To obtain the power spectrum of the lensing potential, we insert the lensing projection (equation 1) into the inverse harmonics expansion (equation 7) and write the 3D potential as its Fourier transform (equation 5) to get

[TABLE]

The 3D position vector $\bm{r}$ is a 3D position vector with polar coordinate $r=\chi$ and polar angles $(\theta,\varphi)$. Similarly we denote with $\theta_{k},\varphi_{k}$ the polar angles of the 3D Fourier vector $\bm{k}$. We insert the expansion of a plane wave into spherical harmonics,

[TABLE]

Making use of the orthogonality of the spherical harmonics

[TABLE]

the expression for $\psi_{\ell m}$ simplifies to

[TABLE]

To obtain the potential power spectrum, we take the absolute square of the last equation and use the definition of the 3D potential power spectrum (equation 12). The delta-function resolves one 3D Fourier integral. We split the second integration into radial and spherical coordinates, ${\rm d}^{3}k={\rm d}kk^{2}{\rm d}\Omega_{k}$ and use once again the orthogonality of the spherical harmonics to resolve the spherical integral. This leads to the potential power spectrum in equation (9).

Appendix B Discussion and comparison to previously published work

In this section we briefly discuss previously published work on the full projection and second-order Limber equation of the lensing power spectrum, cross-checking and comparing their results with our independent findings.

B.1 Kitching et al. 2016 (version 1)

Kitching et al. (2017) compute the full projection of the weak-lensing power spectrum, which they present as spherical-radial representation of the 3D shear field. Our results in equation (21) corresponds to their Equations (7) and (8) assuming a flat Universe and the case of perfect photometric redshifts, with $p(z|z_{\rm p})=\delta_{\rm D}(z-z_{\rm p})$, and for a bin function that selects the redshift bin of $z$, $W^{\rm SR}(z,z_{\rm p})$ is unity if $z_{\rm p}$ is in the redshift bin denoted by $z$, and zero otherwise. We find that Equation (7) in Kitching et al. (2017) is missing a factor $2/\pi$.

Kitching et al. (2017) derive the spherical and extended Limber approximation starting from the full spherical projection in their Appendix. A. We find that the filter function $q$ defined in their Equation (31) has an additional factor of comoving distance $r$, and an additional factor of $\pi/2$.

As shown in this paper, we are unable to reproduce the differences that Kitching et al. (2017) report, between the full spherical solution and the different approximations, neither for the power spectrum nor for the shear correlation function.

B.2 Bernardeau, Bonvin, Van de Rijt, Vernizzi (2012); Van de Rijt (2012)

Bernardeau et al. (2012) present the non-tomographic full projection $C(\ell)$ in the approximation of the 3D potential power spectrum $P_{\Phi}$ separating into $k$- and $\chi$-dependent functions in their Equation (44). Their expression holds for a single source redshift.

The PhD thesis of Van De Rijt (2012) presents an explicit calculation of the second-order Limber approximation. They carry out the derivatives of the kernels $f$ under the assumption of a constant growth-suppression factor $D_{+}(a)/a$. Then, for a constant source comoving distance $\chi_{\rm S}$, the lensing efficiency in equation (2) is $q(\chi)=(\chi_{\rm S}-\chi)\chi_{\rm S}^{-1}$, and the derivative of the separated kernel function in equation (41) can be calculated analytically,

[TABLE]

This result confirms the power-law behaviour for $\chi\ll\chi_{\rm S}$ of the filter function, which we exploited earlier to fit this function.

Inserting equation (55) into the second-order Limber power spectrum in equation (42) and using the inverse Poisson equation to replace the matter with the potential power spectrum we obtain the same expression as Van De Rijt (2012) (their Equation 7.19).

In Fig. 5 we reproduce Fig. 7.3 from Van De Rijt (2012) using a similar set up, a flat $\Lambda$ Universe with $\Omega_{\rm m}=0.3,h=0.65,\Omega_{\rm b}=0.0461,\sigma_{8}=0.8,n_{\rm s}=0.96$. All source galaxies are at redshift $z_{\rm S}=1$. The non-linear 3D matter power spectrum from Takahashi et al. (2012) is used. The ratio of the first-and second-order Limber approximated power spectra to the full projection shows excellent agreement at the sub-percent level.

B.3 LoVerde & Afshordi (2009)

This paper introduces the extended Limber approximation to second order that we apply in this work. Although they present the specific case of 2D galaxy clustering, their calculations are general enough to apply to a weak-lensing context. Their Equation (5) is a spherical cross-power spectrum of two scalar fields $A$ and $B$, projected from 3D to 2D via projection kernels $F_{A}$ and $F_{B}$ defined in their Equation (4). Comparing their expressions with the weak-lensing potential (1), we set $F_{A}(\chi)=F_{B}(\chi)=2c^{-2}D_{+}(\chi)q(\chi)\chi^{-1}$. With that, their Equation (5) is identical to equation (9).

The second-order Limber approximation in Loverde & Afshordi (2008) is presented in eq. (12). This is consistent with our first-order (equation 38) and second-order (equation 39) Limber approximation terms, when accounting for the difference between lensing potential and shear 2D power spectrum, and 3D potential and matter power spectrum.

B.4 Schmidt (2009)

Schmidt (2008) derive the lensing power spectrum in the flat-sky limit, see their Equation (9). Inserting the Poisson equation and the growth function, and writing the redshift filter function $W_{\kappa}$ (eq. 10) in terms of the lensing efficiency and comoving distances, $W_{\kappa}[z(\chi)]=H(z)^{-1}\chi q(\chi)$, their expression, using our notations, reads

[TABLE]

This is consistent with our equation (29) under the additional assumption that mainly modes with $k\chi\approx k\chi^{\prime}\approx\ell$ contribute to the integral; that is modes around the maxima of the Bessel functions. Then, we can draw out of the integral the factor $\ell^{4}\approx k^{4}\chi^{2}{\chi^{\prime}}^{2}$, to recover equation (29). In fact, using the approximation $k\chi=\ell$ seems to go too far, this is already halfway the Limber approximation.

B.5 Giannantonio et al. (2012)

Giannantonio et al. (2012) derive the flat-sky lensing power spectrum in their Equations (25) and (26). Their window function $W^{\varepsilon_{i}}$ defined in their Equation (25) for a flat Universe equals ${\cal A}q(\chi)/a(\chi)$, since ${\rm d}z\,({\rm d}N(z)/{\rm d}z)={\rm d}\chi n(\chi)$; however due to a typo there is a factor $r_{K}[r(z)]$ missing in the window function, which translates into a missing $r_{K}[r]r_{K}[r^{\prime}]$ in the full projection integral777Typo confirmed by T. G., priv. comm.. This also leads to an errornous $r_{K}^{-2}(r)$ in their Limber equation (27). With these factors accounted for, and making the additional approximation $k\chi\approx k\chi^{\prime}\approx\ell$ (see their Appendix B.4) we reproduce the expressions of Giannantonio et al. (2012).

Appendix C Fast evaluation of the shear correlation function on the sphere

The calculation of the shear correlation function on the sphere requires the estimation of the reduced Wigner ${\rm D}$-matrices ${\rm d}^{\ell}_{2\,2}$ and ${\rm d}^{\ell}_{2\,-2}$. The general calculation of ${\rm d}^{\ell}_{m\,n}$ is cumbersome, but there are quick and numerically stable recurrence relations if we are interested in only a subset of these matrices. In particular, following Blanco et al. (1997) we can show that

[TABLE]

which allows to calculate all required reduced Wigner ${\rm D}$-matrices from the first two elements.