Mon. Not. R. Astron. Soc. 356, 829–838 (2005)

doi:10.1111/j.1365-2966.2004.08512.x

Statistical distribution of gravitational-lensing excursion angles: winding
ways to us from the deep Universe
Takashi Hamana,1  Matthias Bartelmann,2 Naoki Yoshida3 and Christoph Pfrommer4
1 National Astronomical Observatory of Japan, Mitaka, Tokyo 181-8588, Japan
2 ITA, Universität Heidelberg, Tiergartenstr. 15, D–69121 Heidelberg, Germany
3 Department of Physics and Astrophysics, Nagoya University, Nagoya 464-8602, Japan
4 Max-Planck-Institut für Astrophysik, Karl-Schwarzschild-Strasse 1, D–85748 Garching, Germany

Accepted 2004 October 14. Received 2004 July 29

ABSTRACT

We investigate statistical distributions of differences in gravitational-lensing deflections between two light rays, the so-called lensing excursion angles. A probability distribution functionof the lensing excursion angles, which plays a key role in estimates of lensing effects on angular clustering of objects (such as galaxies, quasi-stellar objects and also the cosmic microwave
background temperature map), is known to consist of two components: a Gaussian core andan exponential tail. We use numerical gravitational-lensing experiments in a CDM cosmology for quantifying these two components. We especially focus on the physical processes
responsible for generating these two components. We develop a simple empirical model forthe exponential tail which allows us to explore its origin. We find that the tail is generated
by the coherent lensing scatter by massive haloes with M > 1014 h −1 M at z < 1 and thatits exponential shape arises due to the exponential cut-off of the halo mass function at thatmass range. On scales larger than 1 arcmin, the tail does not have a practical influence on thelensing effects on the angular clustering. Our model predicts that the coherent scatter may have
non-negligible effects on angular clustering at subarcminute scales.
Key words: gravitational lensing – cosmology: theory – dark matter – large-scale structureof Universe.

1 INTRODUCTION
Light rays are deflected when they propagate through an inhomogeneous gravitational field, such as the real Universe we live in.
The lensing deflection angle varies from one direction to another,and thus the difference in deflection angles between two light rays,
which we call the ‘lensing excursion angle’, does as well. Consequently, it is not easy to infer the transverse distance between twocelestial objects at a cosmological distance from their angular separation in the sky. Strictly speaking, lacking complete knowledge ofthe matter distribution in the Universe, this is impossible to do.
Because distance is one of most fundamental physical quantities,
the lack of a precise distance measure to the far Universe may prevent
us from a detailed understanding of the Universe. A well-known caseis that angular correlations of distant galaxies and of the temperature
map of the cosmic microwave background are altered by lensing
deflections (see section 9 of Bartelmann & Schneider 2001 for areview and references therein).
Although we cannot know the lensing excursion angle for an individual pair of light rays, knowledge of their statistical distribution
 E-mail: hamanatk@cc.nao.ac.jp
C 2004 RAS

greatly helps us in estimating the order of magnitude of lensing effects. In addition, it allows the intrinsic angular correlation functions
to be deconvolved from measured correlation functions (Bartelmann& Schneider 2001). It is thus of fundamental importance to understand in detail the statistical distribution of the lensing excursion
angles.
The analytic model for computing the variance of lensing excursion angles in the framework of modern cosmological models
was developed by Seljak (1994, 1996), based on linear perturbation theory (the so-called power spectrum approach). Hamana &Mellier (2001) performed numerical experiments of the gravitational lensing deflections in cold dark matter models and examined
the statistical properties of the lensing excursion angles. They foundthat the probability distribution function (PDF) of the excursion angles consists of two components, a Gaussian core and an exponential
tail, and that the variance of the Gaussian core component agrees
well with the prediction by the power spectrum approach. They argued that the exponential tail may be generated by coherent lensing
scattering by massive haloes which is not taken into account in the
power spectrum approach.
The purpose of this paper is twofold. The first is to explore theorigin of the exponential tail of the lensing excursion-angle PDF.
The second is to develop an empirical model for the exponential

830

T. Hamana et al.

tail. To pursue these purposes, we first examine in detail propertiesof the exponential tail of the excursion-angle PDF using numericalexperiments in Section 2. Then in Section 3, we describe a model forthe exponential tail which is based on the assumption that the tailoriginates from coherent lensing scattering by individual massive
dark matter haloes, and compare the model predictions with numerical results. We also discuss a general picture of light propagation in
the Universe paying special attention to the role of secular random
deflections by either large- or small-scale structures and coherent
scatter by a massive halo. Finally, we give a summary and discussion
in Section 4.
2 R AY- T R AC I N G S I M U L AT I O N
2.1 VLS N-body simulation
We performed weak lensing ray-tracing experiments in a Very Large
N-body Simulation (VLS) carried out by the Virgo Consortium(Jenkins et al. 2001, and see also Yoshida, Sheth & Diaferio 2001 forsimulation details). The simulation was carried out using a parallel
P3 M code (Macfarland et al. 1998) with a force softening length ofl soft ∼ 30 h −1 kpc. The simulation employed 5123 CDM particles
in a cubic box of 479 h −1 Mpc side-length, which gives a particle
mass of m part = 6.86 × 1010 h −1 M . It uses a flat cosmologicalmodel with a matter density  0 = 0.3, a cosmological constant  = 0.7, and a Hubble constant H 0 = 100h with h = 0.7. The initial matter power spectrum was computed using CMBFAST (Seljak &Zaldarriaga 1996) assuming a baryonic matter density of  b = 0.04.
The normalization of the power spectrum is taken as σ 8 = 0.9.

scheme (Hockney & Eastwood 1988). The grid size is 0.23 h −1 Mpc
which is chosen to maintain the resolution provided by the N-bodysimulation and removing at the same time the shot noise due todiscreteness in the N-body simulation [this choice is equivalent tothe ‘large-scale smoothing’ in Ménard et al. (2003), and we refer
the reader to this reference for further examination of the effective
resolution of the ray-tracing simulation]. Its computation follows
the procedure described in Hamana & Mellier (2001).
Having produced surface density fields on all lens planes, 10242
rays are traced backwards from the observer’s point using themultiple-lens plane algorithm (e.g. Schneider & Ehlers 1992). The
initial ray directions are set on 10242 grids with a grid size of
0.25 arcmin, thus the total area covered by rays is 4.272 deg2 . We
produced 36 realizations of the underlying density field by randomly
shifting the simulation boxes in the direction perpendicular to theline of sight using the periodic boundary conditions of the N-bodyboxes.
The 36 realizations are not perfectly independent because they are
generated from the same N-body outputs (but using different combinations of random lines of sight) which come from two runs of
N-body simulation. Therefore the generated lensing data (the lensing deflection field, lensing convergence and shear map) are subject
to sample variance. In order to test its magnitude, we compare theconvergence two-point correlation function with its theoretical prediction (Jain & Seljak 1997) in Fig. 1. The measurements from the
ray-tracing experiment are plotted by symbols with error bars which
represent the mean and root mean square among the 36 realizations,while the solid lines show the prediction. Note that the measurement
for z s = 1100 should be compared with the dotted line which shows
the theoretical prediction for z s = 1100 but the contribution from

2.2 Weak lensing ray-tracing simulation
The multiple-lens plane ray-tracing algorithm we used is detailedin Hamana & Mellier (2001); see also Bartelmann & Schneider(1992); Hamana et al. (2000); Jain et al. (2000); Vale & White (2003)
and Hamana et al. (2004) for the theoretical basics and technicalissues; thus in the following we describe only aspects specific to the
VLS simulation data and to ray-tracing experiments in this study.
We use 13 snapshot outputs from two runs of the N-body simulation which differ only in the realization of the initial fluctuation
field. A stack of these outputs provides the density field from z = 0
to z = 6.8. We do not use further higher redshift outputs for two reasons: (1) discreteness effects of particles (Hamana, Yoshida & Suto
2002), and (2) an artificial power excess in the density power spectrum due to the ‘glass’ initial condition (White 1996) at around themean separation length of particles; both of them are significant atsuch high redshifts. For higher redshifts up to the last scattering surface (z ≈ 1100), we simply consider a homogeneous density field.
Thus within 6.8 < z < 1100 rays propagate as in a perfectly homogeneous Universe. This treatment misses lensing contributions fromstructures at that redshift range. It has turned out that this approximation causes only a minor effect. We will discuss its influences on
our analyses later.
Each N-body box is divided into four subboxes with an equal
thickness of 119.75 h −1 Mpc. The N-body particles in each subbox
are projected on to lens planes. In this way, the particle distribution
between an observer and z = 6.8 is projected on to 50 lens planes.
Note that, in order to minimize the difference in redshift between alens plane and an output of N-body data, only one half of the data
(i.e. two subboxes) of z = 0 output is used. The particle distribution
on each plane is converted into the surface density field on a 20482
regular grid using the triangular shaped cloud (TSC) assignment

Figure 1. The convergence two-point correlation functions. Measurementsfrom ray-tracing experiments are shown by symbols with error bars which
represent the mean and root mean square among the 36 realizations. Crosses,
filled triangles and open circle are for z s = 1, 3 and 1100, respectively. The
error bars of the cases for z s = 3 and 1100 are not displayed for clarity
but are of similar magnitude to the z s = 1 case. The solid lines show thetheoretical prediction (e.g. Jain & Seljak 1997) in which the fitting functionof the non-linear power spectrum by Peacock & Dodds (1996) is used toinclude the effect of the non-linear growth of the density field. The dottedline shows the theoretical prediction for z s = 1100 but the contribution fromdensity fluctuation at 6.8 < z < 1100 is not integrated.

C 2004 RAS, MNRAS 356, 829–838

Statistics of lensing excursion angles
the density fluctuations between z = 6.8 and 1100 is ignored. The
measurements are in good agreement with the theoretical prediction in shape, but are slightly higher in amplitude. This excess may
be mostly attributed to the sample variance and implies that there
exists an excess power in the lensing potential field. Because lensing deflections result from the same potential field, it is expectedthat there exists, to a similar extent, an excess in the deflection angle
statistics. On scales smaller than 1 arcmin, the slope of the measured
correlation function becomes flatter than predicted; this is due to thelimited resolution of the N-body simulation. The effective angularresolution of the convergence field is about 1 arcmin for lower redshift (z s < 3) and is slightly better for higher redshifts (see Ménard
et al. 2003 for further discussion on the resolution issue).

2.3 PDF of the lensing excursion angles
Using weak-lensing experiments, we study the statistics of differences in deflection angles between two light rays, which we refer toas the ‘lensing excursion angle’. The deflection angle of a light ray,which is computed by the lens equation, is simply the difference
between its positions θ I and θ S on the image and source planes,
respectively. Denoting the deflection angle of two rays by α1 and
α2 , respectively, we write the lensing excursion angle between theserays as δα = α1 − α2 . Similarly we denote their intrinsic separationby θ 12 = |θ 1S − θ 2S |.
Let us first look into the probability distribution function (PDF)of the lensing excursion angles which is one of most fundamental
statistics. Fig. 2 shows the PDFs of the lensing excursion angles
normalized by its intrinsic separation (i.e. δα/θ 12 ). Because thevector field δα has no preferred direction, we use both components,
δα 1 and δα 2 , to compute the PDFs. The dotted lines in each plotof Fig. 2 show the Gaussian PDF with its standard dispersion (σ 2 )831

computed from the PDF itself, i.e. σ 2 = dx x 2 PDF(x). As was
first pointed out by Hamana & Mellier (2001), the PDFs consist of
two components, a Gaussian core and the exponential tail, which are
generated by different physical processes as we will discuss below.
The origin of the Gaussian core is explained as follows. Light
rays from a cosmological distance undergo many (either strong orweak) gravitational lensing deflections. Because the spatial distribution of lenses at large separation are uncorrelated, rays basicallyundergo many uncorrelated deflections. Provided the separation between two rays is so large that effects of coherent scattering can beignored, two light rays undergo independent deflections. According
to the central limit theorem, the statistical distribution of the lensing excursion angles of such light-ray pairs is given by a Gaussian.
A necessary condition for the central limit theorem to hold is thatthe parent distribution of the individual events which are being superposed has finite variance. The deflection angle calculated in the
weak-lensing regime using the power spectrum approach has finite
variance (Seljak 1994, 1996), but it is based on linearized gravity
and ignores lensing by individual haloes. On the other hand, numerical gravitational lensing experiments show that the PDF consists of
the Gaussian core and the exponential tail. In order to understandwhether the Gaussian core can indeed be caused by the superposition of many deflections, we need to investigate the variance of the
excursion angle. In particular, the numerical experiments miss theinfluence of numerous distant lenses, because of their necessarily
finite volume. We will now show, using a simple approach, that accumulated contributions from very distant lenses do not significantly
affect the excursion-angle variance, thus it remains finite.
Consider first a single light ray passing the lens plane at the origin.
There is a finite number of lenses close to the ray, thus we canrestrict ourselves to distant lenses as we are investigating whether the
deflection-angle variance is finite or not. Axially symmetric lenses
more distant than their (e.g. virial) radii act as point lenses, thus wecan approximate their deflection angles by α i = θ i /θ i2 . Assuming
the lenses have a number density n and are randomly distributed,
the variance of the total deflection angle contributed by lenses in aring around the origin with radius θ and width dθ isα2  =

 N

 1
i=1

θ2

=

2πnθ dθ= 2πn d ln θ.
θ2(1)

Integrating over θ shows that the variance diverges logarithmically.
The situation changes for the excursion angle. Consider two lightrays piercing the lens plane at positions θ 1,2 = (∓d/2, 0). Specializing again to distant lenses, we can approximate the individual
deflection angles by those of point lenses. The excursion angle ofthe lenses in a ring of radius θ and width dθ around the origin canthen be expanded to lowest order in d/θ,δα(θ) =

N

d
i=1

θ2

(−cos 2φi , sin 2φi ) ,(2)

where φ i is the polar angle of the ith lens. Again assuming randomly
distributed lenses, the variance of the excursion angle contributedby the lenses in the ring is thusδα (θ) =2

Figure 2. The solid curves show the probability distribution function of
the lensing excursion angles normalized by the intrinsic ray separation. The
source redshift is z s = 1100. The ranges of the intrinsic separation of thelight-ray pairs θ 12 are given in each panel. The dotted curves show Gaussianwith their σ computed from the measured PDFs, i.e. σ 2 =distributions
dx x 2 PDF(x).

C 2004 RAS, MNRAS 356, 829–838

 N πd 2i=1

θ4=

2π2 d 2 nθ dθ
,
θ4

(3)

i.e. the excursion-angle variance converges like θ −2 when integrated
over θ to infinity. Thus, we can apply the central limit theorem tothe excursion angle, which we could not do for the deflection angle
itself.

832

T. Hamana et al.

Figure 3. The standard dispersion of the lensing excursion angles plotted as
the function of the separation angle. Symbols and error bars show the mean
and rms among 36 realizations of the ray-tracing numerical experiments.
Solid lines represent the theoretical predictions from the power spectrum
approach (Seljak 1996), in which the fitting function of the non-linear power
spectrum by Peacock & Dodds (1996) was used. The source redshifts are,from lower to higher, z s = 1, 3, 6.8 and 1100, respectively. The dotted line
shows the theoretical prediction for z s = 1100 but the contribution fromdensity fluctuation at 6.8 < z < 1100 is ignored.

Now we test the theoretical model prediction for the variance ofthe excursion angles developed by Seljak (1994, 1996) against ournumerical results. Fig. 3 compares the standard dispersion measured
form the numerical experiments with the theoretical prediction. The
dispersion becomes larger as the light rays travel a longer distance,
because the rays can undergo more deflections. It is found in theplot that the measurements are slightly larger than the prediction.
However, a similar excess is seen in the convergence correlation
function (Fig. 1), thus this is mostly due to the sample variance. We
may therefore conclude that the power spectrum approach provides
a good prediction even for z s = 1100, and the non-Gaussian tail has
no strong contribution to the variance. It is important to notice that
coherent scattering by lensing due to massive haloes that generate
the exponential tail contribute only very little to the excess in themeasured dispersion over the prediction.
Let us now turn to the tail of the lensing excursion-angle PDF.
Fig. 4 compares the PDF obtained from the ray-tracing numericalexperiments for three source redshifts, z s = 1, 3 and 1100, andfor various ranges of ray separations. This figure represents major
characteristics of the non-Gaussian tail: (a) it has an approximately
exponential slope; (b) it changes little with the source redshift, but
its amplitude increases with the source redshift, at least within theredshift range we consider (z s > 1). We fit the tail of the PDFs tothe exponential distribution:
E(x) = p exp(−q x).

(4)

To do this, we take two points x 1 and x 2 such that PDF(x 1 ) = 1.0 ×
10−2 and PDF(x 2 ) = 1.0 × 10−3 . Comparison between the results ofz s = 3 and z s = 1100 plotted in Fig. 5 confirms the above point (b)Figure 4. The PDF of lensing excursion angles normalized by its intrinsic separation (plotted positive side only). The source redshifts are, fromnarrower to broader PDF, z s = 1, 3 and 1100, respectively. As these plotsshow, for a given θ range, the amplitude of the PDF tail becomes higher
as z s increases, while their slope is almost unchanged for plotted redshifts
z s > 1.

Figure 5. Parameters in the exponential distribution (equation 4) computedby fitting to the non-Gaussian tail of PDFs from ray-tracing experiments.
The filled circles and crosses show the results from the z s = 1100 and z s =3 cases, respectively. The error bars represent the ranges of ray separations
taken to compute the PDFs.

in a quantitative manner. Note that fitting the exponential function
to the PDF tails becomes poor for large ray separations because the
non-Gaussian tail does not appear prominently due to the limited
statistics. This accounts for the steep rise of both p and q at large
θ 12 which rather reflects the slope of the Gaussian core.
C 2004 RAS, MNRAS 356, 829–838

Statistics of lensing excursion angles

833

multiple strong lensing events is negligible. The former is reasonable because a process that is not taken into account in the power
spectrum approach could generate non-Gaussian features. Also thevisual impression from Fig. 6 could be support for that idea. The
latter is validated by the observational fact of the small cross-section
for strong lensing events by a single lens (either a galaxy or clusterof galaxies) such as multiply imaged quasi-stellar objects (QSOs)
and strongly lensed arc-like images of distant galaxies. Thus, it iscertain that multiple scattering by more than one massive halo isvery rare.
We consider the same CDM cosmology as the one adopted
for the numerical experiments in Section 2. We denote the PDF
for finding a ray pair with θ 12 having the excursion angle δα by
PDF(δα|θ 12 ).
3.1 Lensing deflection by a universal density profile halo

Figure 6. The celestial distributions of ray pairs having a large lensing
excursion angle and of massive haloes in one realization of the numericalexperiment. See Hamana et al. (2004) for a detailed description of construction of the halo catalogue on a light-cone. Black symbols represent
haloes; the large filled circles, small filled circles and dots are for haloes with
M halo > 4 × 1014 , 4 × 1014 > M halo > 4 × 1013 and 4 × 1013 > M halo
> 1 × 1013 h −1 M , respectively. Only the haloes within the redshift interval between 0 and 1 are displayed. Red dots represent the ray pairs whichobey the following criteria: an unlensed ray separation of θ 12 > 1 arcmin
and an excursion angle of |δα|/θ 12 > 0.7. The mid-points of ray pairs are
displayed.

In the next section, we develop a model of the exponential tailand explore its origin. To do this, a visual impression from Fig. 6
could be informative. In this figure, the celestial distributions of
ray pairs having a large lensing excursion angle and of massive
haloes (M halo > 1013 h −1 M ) in one realization of the numericalexperiment are displayed (see Hamana et al. 2004 for a detailed
description of the construction of a halo catalogue on a light-cone).
Only the haloes within the redshift interval between 0 and 1 areplotted using black symbols. Red dots represent the ray pairs whichobey the following criteria: an unlensed ray separation of θ 12 > 1
arcmin and an excursion angle of |δα|/θ 12 > 0.7. Clearly, most of
the large excursion-angle ray pairs pass very close to a massive halo.
We argue in the following section that the exponential tail resultsfrom coherent strong deflections of two nearby rays by a massive
halo, and explain the origin of the above characteristics using simple
models.
3 O R I G I N O F T H E E X P O N E N T I A L TA I L O F
THE EXCURSION-ANGLE PDF
In this section we explore the origin of the exponential tail of the
excursion-angle PDF found in the ray-tracing numerical experiments. For this purpose, we focus on the tail part and do not consider
the Gaussian core whose origin has been investigated in the literature(Seljak 1994, 1996; see also chapter 9 of Bartelmann & Schneider)and also in the last section. We develop a theoretical model from two
assumptions: large excursion angles are mainly caused by the strong
lensing of a massive halo, and the probability for a ray to undergo
C 2004 RAS, MNRAS 356, 829–838

Navarro, Frenk & White (1996, 1997, NFW hereafter) found from
N-body simulations that the density profile of dark matter haloes canbe fitted by a universal form regardless of their mass and redshift.
We adopt a truncated universal profile:ρsr
ρ(x) = s, x= ,(5)x (1 + x)3−srs
for r < r vir and 0 otherwise, where r s and r vir are the scaleradius and virial radius, respectively. It is convenient to introduce
the concentration parameter c vir = r vir /r s . Navarro et al. (1996)
proposed the universal inner slope of s = 1, while a steeper slope
was claimed by later studies using higher resolution N-body simulations; Moore et al. (1998, 1999), Ghigna et al. (2000) andFukushige & Makino (2001, 2003) found larger values such as s
= 1.5, while Jing (2000) and Jing & Suto (2000) pointed out that
it varies from 1.1 to 1.5 and argued a possible weak dependence
on the halo mass. In this paper, we consider the two cases s= 1
and 1.5. Navarro et al. (1997) and Bullock et al. (2001) have extensively examined a relation between the concentration parameter and the halo mass and its redshift evolution adopting a fixedvalue of s = 1. We adopt a generalized mass–concentration relation proposed by Oguri, Taruya & Suto (2001); see also Keeton &Madau (2001):c∗cvir (M, z) = (2 − s)
1+z



M
1014 h −1 M−0.13.

(6)

Bullock et al. (2001) suggested c ∗ ∼ 8 for the CDM model, which
we adopt as a fiducial choice. We note that there is a relatively
large scatter in this relation (Bullock et al. 2001; Jing 2000). The
virial mass (defined by the mass within the virial radius r vir ) of theuniversal halo is given by3
Mvir = 4πρs rvirm(cvir , s)
,
3
cvir

(7)

with

 cvirm(cvir , s) =

dx
0

x 2−s
.
(1 + x)3−s(8)Because the spherical collapse model indicates that Mvir =34πrvirδvir (z)ρ̄0 /3, where δ vir is the overdensity of collapse (seeNakamura & Suto 1997 and Henry 2000 for useful fitting functions), one can express ρ s in terms of δ vir (z), c vir and s:
ρs =

3
cvir
δvir ρ̄0.
3 m(cvir , s)

(9)

834

T. Hamana et al.

Let us summarize the basic equations for gravitational lensing
properties of the truncated universal profile halo (Takada & Jain
2003; see Bartelmann 1996; Wright & Brainerd 2000 and Oguri et al.
2001 for lensing properties of the non-truncated universal profile
lens model). The surface mass density of the truncated universal
profile halo is given by


(y) =



2 −y 2
cvir



−

with


f (y) =

dz ρ(y, z) = 2ρs rs f (y),y=2 −y 2cvir

2 −y 2
cvir

dz
0

1
(y 2 + z 2 )s/2 (1 +



r
,
rs

y 2 + z 2 )(3−s)

(10)

,

(11)

for y  c vir and f (y) = 0 otherwise. The projected mass within aradius b is

 bM(< b) = 2π

dy y
0

(y )

 h

= 4πrs3 ρs(12)dy y f (y),0where h = b/r s . We perform the above integration numerically. The
thin lens equation is written byθS = θI − α(Dl θI ),Figure 7. Deflection angle profiles of the truncated universal density profile
lens, g(θ ), defined in equation (16), as a function of the impact parameter(normalized by the angular virial radius) for various concentration parameters denoted in the plot. The left panel is for s = 1 and right panel for s= 1.5. Note that g(θ ) has an asymptotic inner slope of ∝ θ 0.82 (∝ θ 0.48 ) fors = 1(s = 1.5).

(13)
3.2 Lensing excursion angles

withα(b) =

4G M(< b) Dls
.
c2 b
Ds

(14)

In the last expression, the origin of the coordinates is taken at thelens centre, b = D l θ I is the impact parameter, and D l , D ls and D s
are the angular diameter distances from observer to lens, from lensto source, and from observer to source, respectively. The deflection
angle of the truncated universal profile halo is given byα(θ ) = α∗ g(θ ),
with
cvir
g(θ ) =m(cvir , s)andDls
α∗ = 2m
Ds



4



(15)

x

dy y f (y)

0

x
H0
c

,

x=

Dl θ
rs

(16)

2
2
rvir
δvir

Mvir
1014 h −1 M



rvir
1 h −1 Mpc

−1

Dls
.
Ds

(17)

Note that α ∗ ∝ M . It is important to notice that the dependenceof the halo profile parameters on the deflection angle enters only
through the function g(θ ). Note that for θ  θ vir it reduces to g(θ) =θ vir /θ (where θ vir is the angular virial radius defined by θ vir =r vir /D l ). In Fig. 7, the function g(θ ) is plotted for various valueof c ∗ . As one may see in the figure, the deflection angle profileg(θ ) peaks at θ ∼ θ vir /c ∗ = θ s (θ s = r s /D l ) and the peak value
does not strongly depend on the inner slope s. It is also found thatthe peak value relates to the concentration parameter roughly by
g max ∼ 0.1 c ∗ + 1. Therefore, in a reasonable range of c ∗ the maximum deflection angle by a single universal halo lens is α max = (1 −3) α ∗ . One may also find that at the inner part, the deflection angle
is larger for a larger c ∗ or for a steeper inner slope (thus for more
centrally concentrated haloes). We find that g(θ) has an asymptotic
inner slope of ∝ θ 0.82 (∝ θ 0.48 ) for s = 1 (s = 1.5).
2/3

Because the deflection angle of a universal density profile halo is
finite, the excursion angle is finite as well. Clearly, the largest excursion angle is 2α max which happens when one ray passes at the
distance ∼θ s from the lens and the other ray passes at the same
distance in the opposite side of the lens. Thus this happens only if
θ 12 = 2θ s , and is very rare. Let us consider the maximum excursion
angle produced for other ray separations. If θ 12 > θ s , the largest
excursion angle is in the range α max < δα < 2 α max . This happens
when one ray passes at ∼θ s and the another ray passes at the opposite side of the lens in the direction connecting the lens centre andthe first ray, while if θ 12 < θ s , the largest excursion angle is smaller
than α max . An important consequence of this is that for ray pairs
with separation angle larger than θ s , the maximum excursion angle
does not strongly depend on the ray separation but lies in a small
range of (1 − 2)α max , and it scales with the halo mass as ∝ M 2/3 .

3.3 PDF of δα
Let us first consider the probability distribution induced by one halo,
which we denote by PDF 1 (δα|θ 12 ). Let p(θ) be the probability of
a ray passing at small area θ → θ + δθ from a lens centre, which is
given by the cross-sectional area (denoted by A) normalized by theunit solid angle (d):p(θ) =

A(θ → θ + δθ)
.
d

(18)

Then PDF 1 (δα|θ 12 ) is given by the joint probability



PDF1 (δα|θ12 ) =

d2 θ1
d



dφ12
p(θ 1 ) p(θ 1 + θ 12 ),
2π

(19)

where θ 12 = {θ 12 cos(φ 12 ), θ 12 sin(φ 12 )}. Note that as we are considering lensing by a single halo having a certain density profileand mass, given a configuration of a light ray pair, its excursion
angle is uniquely determined. The total PDF(δα|θ 12 ) is obtained by
C 2004 RAS, MNRAS 356, 829–838

Statistics of lensing excursion angles

835

summing PDF 1 (δα|θ 12 ) over haloes within a light-cone volume,

PDF(δα|θ12 ) =dVdM n halo (M, z)PDF1 (δα|θ12 ) r (zs )=dr r 2dM

0× n halo (M, r [z])PDF1 (δα|θ12 ),(20)

where r is the comoving radial distance, dV = r 2 dr (this expression is valid only for a flat cosmological model) is the unit volume
element and n halo (M, z) is the halo mass function. We adopted themass function by Sheth & Tormen (1999). Note that this approach
breaks down for small excursion angles where secular scattering by
distant and/or small haloes is important; this generates the Gaussiancore of the excursion-angle distribution.

3.4 Results
Let us start with the excursion-angle PDFs from one halo plottedin Fig. 8 which helps us to understand the origin of the exponential
tail. The most important point which should be noticed is the sharp
cut-off in a large excursion angle. This is a natural consequence ofthe fact that the deflection angle of the universal density profile halo
is finite (see Section 3.1). The maximum excursion angle scales withthe halo mass roughly as ∝ M 2/3 (with a small correction by the massdependence of the concentration parameter) as long as the separation
angle is larger than the angular scale radius (θ s ) of a lensing halo, as
explained in the last subsection. The other important point is that forthe mass-independent double power-law PDF, its power-law slope
is ∝ δα −2 for smaller excursion angles and ∝ δα −2.64 for larger
angles. The former is generated by ray pairs in which both rays passoutside of the virial radius, while the latter is generated by pairs ofrays of which one passes outside of the halo, and the other inside.

Figure 8. PDF of the lensing excursion angles from one halo, PDF 1 (δα),
for three haloes masses M halo = 1013 , 1014 and 1015 h −1 M . The normalization is arbitrary. The lens and source redshifts are z l = 0.3 and z s = 1100,
respectively. The solid lines are for s = 1 and dashed lines are for s = 1.5.
The concentration parameter taken is c ∗ = 8 for all cases.
C 2004 RAS, MNRAS 356, 829–838

Figure 9. Model prediction of the lensing excursion-angle PDF for z s =1100, c ∗ = 8 and s = 1. The solid line shows the full PDF, while brokenlines show contribution from a limited range of the halo mass denoted in thepanel.

Under the assumptions stated in the last subsection, the excursionangle PDF is obtained by summing up contributions from single
haloes over a wide range of the halo mass and integrating over theredshift of haloes as defined by equation (20). We plot the PDF
computed from such a model in Fig. 9. The solid line shows thetotal PDF, while broken lines represent contributions from narrow
limited mass ranges. This figure clearly illustrates the origin of the
exponential tail. There are two key points; a large excursion angle
can only be generated by massive haloes with mass typically larger
than 1014 h −1 M , and at such a halo mass range, the mass function
decreases exponentially. Accordingly, the number of more massive
haloes that can contribute to a larger excursion angle decreases exponentially, and as a result, the exponential slope of the PDF arises.
In order to examine what redshift range of haloes makes a major
contribution to the exponential tail, we plot in Fig. 10 the excursionangle PDFs computed for limited ranges of the lens redshifts (upper
panel) and their percentage of the contribution (lower panel). It isseen in the lower panel that contributions from haloes at redshifts
below 1 account for almost the full amplitude of the PDF tail. Thisis explained by a rapid evolution of the halo mass function at thehigh-mass end. In fact, the number density of massive haloes with
M halo  2 × 1014 h −1 M decreases by more than one order ofmagnitude from z = 0 to 1.
The above results lead to the following explanation for the origin of the lensing excursion-angle PDF which has two components,
a Gaussian core and an exponential tail. The light rays emitted athigh redshifts undergo many gravitational deflections by (large- or
small-scale) structures on their way to us, which after many uncorrelated deflections produces the Gaussian core. Some small part of
the rays are strongly lensed by a massive halo with mass larger than1014 h −1 M at a low redshift of z < 1, and the coherent deflections
caused by the strong lensing produce the exponential tail. Therefore
even if a ray pair encounters a strong coherent deflection by a single

836

T. Hamana et al.

Figure 10. Model prediction of the lensing excursion-angle PDF for z s =1100, c ∗ = 8 and s = 1. The solid line shows the full PDF, while the brokenlines show contribution from a limited range of the lens redshift denoted in
the panel.

Figure 11. The model predictions of the lensing excursion-angle PDF after
the correction to the accumulative lensing effect being made (using the approximate way of equation 22). The halo parameters are c ∗ = 8 and s = 1.
The ray separations are denoted in each plot. The source redshifts are, frominner to outer, z s = 1, 3 and 1100, respectively.

massive halo, its excursion angle is not solely determined by thestrong halo lensing but the random scattering also contributes to it
to a smaller extent. This effect on the PDF is taken into account by
the convolution

1.0 × 10−3 to the exponential function. Note that after this correction
the amplitude of the tails can be increased more than one order ofmagnitude because of the steep slope of the exponential tail (thusfor a large q).
Four panels of Fig. 12 show the corrected model PDFs for various
values of c ∗ and s, which help us to understand the dependencies ofthe shape of the PDF tail on the halo parameters. Clearly, the broader
tail appears for models with a larger c ∗ or a larger s, because suchmodels generate a larger maximum deflection angle. It is important
to notice that a small change of c ∗ or s causes a very large, nonlinear change in the shape and amplitude of the PDF. Therefore,choosing mean values of c ∗ and s does not provide a mean PDF
but gives a lower amplitude. It should be noted that the halo model
parameters indeed have a large scatter (Jing 2000; Jing & Suto 2000;
see also fig. 9 of Hamana et al. (2004) which clearly shows that thecompactness of halo mass distributions has a large scatter).
Finally, we compare in Fig. 13 the parameters in the exponentialfunction (22) measured from the model PDFs (symbols with lines),
with the results from the ray-tracing experiments. We evaluate p
and q by fitting the corrected model PDF, i.e. after correction byequation (22) being made. As shown in Fig. 13, the parameter q
correlates with the slope parameter q as expected. The measured
exponential slope parameter q plotted in the lower panel is larger
than the results from the numerical experiments, though the slopewith the separation angle is very similar. The discrepancy is smaller
for models with a larger c ∗ or a larger s. Therefore the model prediction may be improved if one takes into account the scatter in thehalo parameters c ∗ or a larger s.
Quantitatively, none of the four models plotted in Figs 12 and 13
are in very good agreement with the simulation results. The discrepancy is partly due to the scatter in the halo model parameters as has
been discussed above. Also, a deviation in the halo mass distribution
from spherical symmetry could partly account for it. Actually, themass distribution of most of the haloes significantly deviates fromE(x) =dy G(y)E (x − y),(21)

where G(y) denotes the Gaussian distribution, and E (x) denotes
the tail part produced by the halo lensing without considering thecontribution from the random scattering (which could be computedby the model described in this section). Because, as shown in Fig. 9,
the tail is well approximated by the exponential shape with a constant index over a wide range of δα, it is reasonable to approximate
E (x) p exp(−q x) and
E(x)

p



√
2πσ



y2
dy exp − 22σexp[−q(x − y)](22)= p exp(−q x),
with the boosted amplitude
p = p exp

(σ q)2.
2(23)

Thus under the above approximation, the slope of the exponential
tail is unchanged but the amplitude is increased. This combinedwith the fact that almost all contributions to the PDF tail come fromcoherent scattering by massive haloes at z < 1 (Fig. 10) accounts forthe trend observed in the excursion-angle PDFs obtained from the
ray-tracing simulation that its slope does not depend strongly on thesource redshift but its amplitude increases with the source redshift.
Actually, a very similar trend is observed in the model PDFs plottedin Fig. 11 which shows the corrected PDF tails for z s = 1, 3 and1100 and for four ray separations. Here, in order to compute theboost factor of equation (23), we compute p and q by fitting themodel PDFs at two points, PDF(x 1 ) = 1.0 × 10−2 and PDF(x 2 ) =
C 2004 RAS, MNRAS 356, 829–838

Statistics of lensing excursion angles

837

spherical symmetry (Hamana et al. 2004). In addition, the spatial
correlation of haloes may have an influence on the excursion-angle
PDF, because massive haloes are strongly clustered. A close lookat the sky distributions of ray pairs having a large excursions angle
and of massive haloes shown in Fig. 6 reveals that those deviations
from our simple model should indeed play a role; namely, it isseen in the figure that a small part of most massive haloes does notproduce a large lensing excursion event, and that a small part oflarge excursion-angle ray pairs does not intersect a very massivehalo.
We may conclude, from what has been seen above, that our simple model succeeds in getting the essential mechanism of generating the exponential tail and in explaining the origin of the major
characteristics of the tail. The model predictions are in reasonableagreement with the simulation results. Further modifications of themodel taking into account details of halo properties, such as scatterin the halo model parameters, deviations of the halo mass distribution from spherical symmetry and clustering of haloes, are needed
to improve the accuracy of the model prediction.
4 S U M M A RY A N D D I S C U S S I O N
Figure 12. The model predictions of the lensing excursion-angle PDFs. We
have applied the corrections given by equations (22) and (23). Ray separations are θ = 2, 8 and 32 arcmin (solid, dotted and dashed, respectively). The
source redshift is taken by z s = 1100. Halo model parameters are denotedin each panel. The thin lines in the top right and two bottom panels show,
for comparison, the predictions of the fiducial model plotted in the top left
panel.

Figure 13. Parameters in the exponential distribution (equation 22) computed from the fit to the model predictions of excursion-angle PDFs forz s = 1100. Model parameters are denoted in the plot. Dotted lines show theresults from ray-tracing simulation.

C 2004 RAS, MNRAS 356, 829–838

We have investigated the statistical distribution of lensing excursion
angles, paying special attention to the physical processes that are responsible for generating two components of the PDFs: the Gaussiancore and the exponential tail. We have used numerical gravitational
lensing experiments in a CDM cosmology to quantify these twocomponents.
The origin of the Gaussian core is explained by the random lensing
deflections by either linear or non-linear structures. The variancesof the Gaussian core measured from the results of the numericalexperiments are found to be in a good agreement with the prediction
by the power spectrum approach (Seljak 1994, 1996).
The presence of the exponential tail was first found by Hamana &Mellier (2001) but its origin remains unrevealed. The tail is characterized by two parameters: the slope and amplitude. We have foundfrom the numerical experiments that the slope changes little with thesource redshift while the amplitude becomes greater as the source
redshift increases, at least within the redshift range we consider(1 < z s < 1100). Because the random lensing deflections resultin the Gaussian core, the exponential tail is most likely to result
from coherent deflections. In addition, in order to generate a large
excursion angle, massive virialized objects should be responsible
for the exponential tail. Therefore, we supposed that the exponential tail originates from coherent lensing scatter by single massive
haloes.
We have developed a simple empirical model for the exponential
tail of the lensing excursion-angles PDF. We used the analytic models of the dark matter haloes, namely the modified Press–Schechter(Press & Schechter 1974) mass function (Sheth & Tormen 1999)
and the universal density profile first proposed by Navarro et al.
(1996). Although we only consider coherent lensing scatter by asingle massive halo and did not take into account the scatter in thehalo parameters (the concentration parameter c ∗ and the inner slopes), our model reasonably reproduces the exponential tails computedfrom the numerical experiments. It is found that massive haloes with
M > 1014 h −1 M are responsible for the tail and that the exponential shape arises as a consequence of the exponential cut-off of thehalo mass function at this mass range. Almost all contributions to
the tail come from the haloes at redshifts below 1. Therefore, theslope of the tail is formed by the haloes at z < 1. On the other hand,the amplitude of the tail is determined by the convolution of two

838

T. Hamana et al.
ray pairs should be investigated in future work with a gravitational
numerical experiment having a higher resolution.
AC K N OW L E D G M E N T S
The N-body simulations used in this work were carried out by the
Virgo Consortium at the computer centre at Max-Planck-Institut,
Garching (http://www.mpa-garching.mpg.de/NumCos). The numerical computations presented in this paper were partly carried out
at ADAC (the Astronomical Data Analysis Centre) of the National
Astronomical Observatory Japan, at the Yukawa Institute Computer
Facility and at the computing centre of the Max-Planck Society in
Garching.
REFERENCES

Figure 14. The empirical model predictions of the tail part of the excursionangle PDF for subarcmin ray separations. The source redshift is z s = 1100,
the concentration parameter is c ∗ = 8 and the inner slope is s = 1 (left) ands = 1.5 (right).

contributions, the coherent scatter and the random deflections. Because the contribution from the random deflections becomes greater
as the source redshift increases, the amplitude of the tail becomes
greater for a higher source redshift. These explain the redshiftindependent slope and the redshift-dependent amplitude found fromthe numerical experiments.
Does the exponential tail have an influence on the angular power
spectrum of the temperature map of the cosmic microwave background (C )? As far as angular scales larger than 1 arcmin are concerned, the answer is no. When one computes the lensed C one
can safely use the approximate convolution equation given by Seljak (1996), because the key assumption in the approximation madefor deriving the convolution equation is not the Gaussianity of the
excursion-angle PDF, but that its variance is small (cf. section 9of Bartelmann & Schneider 2001). Because the smallness of thevariance is also the case for lower redshifts, the same convolution
technique can be applied to other angular correlation functions suchas that of galaxies and QSOs.
Before closing this paper, we present predictions for tails of the
excursion-angle PDF for subarcminute ray separations. We plot in
Fig. 14 our empirical model predictions. It is found that for ray separations θ 12 > 10 arcsec, the amplitude of the tail keeps increasing
with the slope becoming flatter in a similar rate of larger separations.
However, below that separation, the growth rate becomes graduallysmaller. The standard dispersions of these distributions are smaller
than unity but can be of order O(0.1) which is comparable to that of
the Gaussian core, though our simple model adopting average haloparameters tends to predict greater amplitude than the results fromnumerical experiments (see Section 3.4). Therefore, it is possiblethat on arcsecond scales coherent lensing deflections have a nonnegligible influence on the angular clustering of objects in the distant
Universe. Note that even if the exponential tail is taken into account,the standard dispersion of the excursion angles is less than unity,
thus the approximate convolution equation can still be valid but thecontribution from the tail to the dispersion should be included. We
notice that it is, however, not clear whether the assumptions in ourmodel are still valid for such small ray separations. The statistical
distribution of lensing excursion angles for arcsecond-separation

Bartelmann M., 1996, A&A, 313, 697
Bartelmann M., Schneider P., 1992, A&A, 259, 413
Bartelmann M., Schneider P., 2001, Phys. Rep., 340, 291
Bullock J. S., Kolatt T. S., Sigad Y., Somerville R. S., Kravtsov A. V., Klypin
A. A., Primack J. R., Dekel A., 2001, MNRAS, 321, 559
Fukushige T., Makino J., 2001, ApJ, 557, 533
Fukushige T., Makino J., 2003, ApJ, 588, 674
Ghigna S., Moore B., Governato F., Lake G., Quinn T., Stadel J., 2000, ApJ,
544, 616
Hamana T., Mellier Y., 2001, MNRAS, 169, 176
Hamana T., Martel H., Futamase T., 2000, ApJ 529, 56
Hamana T., Yoshida N., Suto Y., 2002, ApJ, 568, 455
Hamana T., Takada M., Yoshida N., 2004, MNRAS, 2004, 350, 893
Henry J. P., 2000, ApJ, 534, 565
Hockney R. W., Eastwood J. W., 1988, Computer Simulation Using Particles.
Adam Hilger, Bristol
Jain B., Seljak U., 1997, ApJ, 484, 560
Jain B., Seljak U., White S. D. M., 2000, ApJ, 530, 547
Jenkins A., Frenk C. S., White S. D. M., Colberg J. M., Cole S., Evrard
A. E., Couchman H. M. P., Yoshida N., 2001, MNRAS, 321, 372
Jing Y. P., 2000, ApJ, 535, 30
Jing Y. P., Suto Y., 2000, ApJ, 529, L69
Keeton C. R., Madau P., 2001, ApJ, 549, L25
Macfarland T., Couchman H. M. P., Pearce F. R., Pichlmeier J., 1998, New
Astron., 3, 687
Ménard B., Hamana T., Bartelmann M., Yoshida N., 2003, A&A, 403,
817
Moore B., Governato F., Quinn T., Statal J., Lake G., 1998, ApJ, 499,
L5
Moore B., Quinn T., Governato F., Statal J., Lake G. 1999, MNRAS, 310,
1147
Nakamura T. T., Suto Y. 1997, Prog. Theor. Phys., 97, 49
Navarro J., Frenk C., White S. D. M., 1996, ApJ, 462, 564
Navarro J., Frenk C., White S. D. M., 1997, ApJ, 490, 493
Oguri M., Taruya A., Suto Y., 2001, ApJ, 559, 572
Peacock J. A., Dodds S. J., 1996, MNRAS, 280, L19
Press W. H., Schechter P., 1974, ApJ, 187, 425
Schneider P., Ehlers J., Falco, C. C., 1992, Gravitational Lenses. SpringerVerlag, New York
Seljak U., 1994, ApJ, 436, 509
Seljak U., 1996, ApJ, 463, 1
Seljak U., Zaldarriaga M., 1996, ApJ, 469, 437
Sheth R. K., Tormen G., 1999, MNRAS, 308, 119
Takada M., Jain B., 2003, MNRAS, 340, 580
Vale C., White M., 2003, ApJ, 592, 699
White S. D. M. 1996, in Schaefer R., Silk J., Spiro M., Zinn-Justin J., eds,
Cosmology and Large-Scale Structure. Elsevier, Dordrecht
Wright C. O., Brainerd T. G., 2000, ApJ, 534, 34
Yoshida N., Sheth R. K., Diaferio A., 2001, MNRAS, 328, 669
This paper has been typeset from a TEX/LATEX file prepared by the author.
C 2004 RAS, MNRAS 356, 829–838