ABSTRACT We present new measurements of dark matter distributions in 25 X-ray luminous clusters by making a full use of the two-dimensional (2D) weak lensing signals obtained from high-quality Subaru/Suprime-Cam imaging data. Our approach to directly compare the measured lensing shear pattern with elliptical model predictions allows us to extract new information on the mass distributions of individual clusters, such as the halo ellipticity and mass centroid. We find that these parameters on the cluster shape are little degenerate with cluster mass and concentration parameters. By combining the 2D fitting results for a subsample of 18 clusters, the elliptical shape of dark matter haloes is detected at 7σ significance level. The mean ellipticity is found to be 〈e〉 = 〈1 − b/a〉 = 0.46 ± 0.04 (1σ), which is in excellent agreement with a theoretical prediction based on the standard collisionless cold dark matter model. The mass centroid can be constrained with a typical accuracy of ∼ 20 arcseconds (∼ 50h−1 kpc) in radius for each cluster. The mass centroid position fairly well matches the position of the brightest cluster galaxy, with some clusters showing significant offsets. Thus the 2D shear fitting method enables to assess one of the most important systematic errors inherent in the stacked cluster weak lensing technique, the mass centroid uncertainty. In addition, the shape of the dark mass distribution is found to be only weakly correlated with that of the member galaxy distribution or the brightest cluster galaxy. We carefully examine possible sources of systematic errors in our measurements including the effect of substructures, the cosmic shear contamination, fitting regions, and the dilution effect, and find none of them to be significant. Our results demonstrate the power of high-quality imaging data for exploring the detailed spatial distribution of dark matter, which should improve the ability of future surveys to conduct cluster cosmology experiments.

Gravitational lensing is one of the most powerful methods to constrain the cluster mass distribution. This is because gravitational lensing probes the matter distribution directly, regardless of whether it is luminous or dark. This is in marked contrast with other probes such as X-ray and Sunyaev-Zel’dovich effect for which assumptions on gas state and other pressure components are crucial in extracting the mass distribution. In particular the cluster mass distribution has extensively been studied using weak lensing technique which makes use of small distortions of background galaxies produced by a foreground cluster (see Bartelmann & Schneider 2001, for a review). Indeed the weak lens technique has been quite successful in constraining dark matter distributions in clusters independently of the distribution of the hot gas component (Kneib et al. 2003; Clowe et al. 2006; Jee et al. 2007; Mahdavi et al. 2007, 2008; Okabe & Umetsu 2008; Bradac et al. 2008).

While the weak lensing technique allows the non-parametric reconstruction of the mass map (Kaiser & Squires 1993), azimuthally-averaged onedimensional (1D) tangential shear profiles have been used in most of quantitative studies so far (e.g., Bardeau et al. 2007; Broadhurst et al. 2008; Hamana et al. 2009; Oguri et al. 2009; Okabe et al. 2010). The main reason for this is that the 1D profile fitting provides an easy and efficient way to extract several important cluster properties, such as mass, concentration, and radial density profiles of clusters. However, the highly non-spherical nature of CDM haloes predicted by N-body simulations indicates that the full two-dimensional (2D) study of weak lensing data is crucial for further tests of the cluster mass distribution. Moreover, an important disadvantage of the 1D profile fitting is that we have to assume the centre of the cluster a priori. It has often been assumed that the mass centroid coincides with the location of the brightest cluster galaxy (BCG), despite the fact that offsets of BCGs from the mass centroids appear to be common (e.g., Katayama et al. 2003; Lin & Mohr 2004; Koester et al. 2007; Ho et al. 2009; Sanderson et al. 2009). The 2D weak lensing studies with non-spherical mass models have been attempted, but only for a handful of lensing clusters (e.g., Cypriano et al. 2004; Deb et al. 2010). Fitting of weak lensing data by taking fully into account a 3D (triaxial) cluster shape and the projection along arbitrary line-of-sight has been performed by Oguri et al. (2005), and later by Corless et al. (2009), but their interests lied in the determination of masses and concentration parameters rather than the shape of the projected mass distribution. Given the strong selection bias of the apparent 2D cluster ellipticity for lensing clusters (Oguri & Blandford 2009), more systematic 2D weak lensing studies of massive clusters should be conducted using well-defined statistical cluster samples.

Cluster mass distributions have also been studied in details using stacked weak lensing data (Sheldon et al. 2007a,b; Johnston et al. 2007; Mandelbaum et al. 2006b, 2008). Although the stacked weak lensing technique is quite successful in extracting average properties of clusters, its main disadvantage comes from the centroid problem (e.g., Johnston et al. 2007; Mandelbaum et al. 2010); in stacking clusters one has to assume a priori the mass centroid for each cluster. The misidentification of the mass centroid results in the suppression of weak lensing signals, particularly near the cluster centre. Recently, Evans & Bridle (2009) measured the ellipticity of isolated clusters with stacked weak lensing technique to obtain the average cluster ellipticity of e ∼ 0.5. However, their result is built on the assumption that the cluster mass distributions is aligned with the spatial distribution of member galaxies, and also on the assumption about the cluster centre, which makes the interpretation somewhat difficult.

In this paper, we present a systematic 2D study of weak lensing maps. We analyse 25 clusters at 0.15 < z < 0.3 presented by Okabe et al. (2010), who examined azimuthallyaveraged 1D tangential profiles to study the radial mass profiles of the clusters. The cluster sample is essentially X-ray flux-limited statistical sample from the Local Cluster Substructure Survey (LoCuSS; G. Smith et al., in preparation). In this paper we extend the lensing analysis of Okabe et al. (2010) and perform a full 2D fitting of the shear fields with particular emphasis on the projected 2D ellipticities of the clusters. All the weak lensing data are based on optical images taken with Subaru/Suprime-cam (Miyazaki et al. 2002) which is known to be the current best ground facility for weak lensing studies because of its exquisite image quality and wide field-of-view (e.g., Kasliwal et al. 2008). In our analysis we leave the centre of the mass profile as a free parameter and fit it simultaneously in order to see how the 2D analysis of shear maps constrains the cluster mass centroid which has often been assumed to coincide with the location of the BCG. Our approach to fit individual clusters with the centroid and orientation of the cluster as free parameters not only overcomes the problems inherent to stacked weak lensing technique but also allows us to even test the assumptions on the centroid and orientation from weak lensing data alone.

The structure of this paper is as follows. We describe our fitting procedure in §2, and present results in §3. In §4 we discuss possible systematic effects, and in §5 we draw our conclusions. Throughout the paper we assume a flat universe with the matter density ΩM = 1 − ΩΛ = 0.274, the dimensionless Hubble constant h = 0.704, the baryon the normalisation of the matter power spectrum σ8 = 0.812 (Komatsu et al. 2009).

2 DATA ANALYSIS

2.1 Cluster sample

In this paper, we study the same cluster sample as in Okabe et al. (2010). The sample comprises 30 massive clusters observed with Subaru/Suprime-cam, and represents an unbiased subsample of the LoCuSS clusters, an all-sky sample of ∼ 100 X-ray luminous (LX & 2 × 1044ergs−1) galaxy clusters located at 0.15 < z < 0.3 selected from the ROSAT All Sky Survey catalogues (Ebeling et al. 1998, 2000; Bohringer et al. 2004). Since X-ray properties of clusters are not very much affected by the degree of the elongation along the line-of-sight (e.g., Gavazzi 2005), the orientation bias, which is important for lensing-selected clusters

The properties of the clusters and the basic parameters of the Subaru/Suprime-cam observations are summarised in Table 1 of Okabe et al. (2010). All the 30 clusters were observed with sub-arcsecond seeing conditions. Among these 25 clusters have two-filter images (mostly V - and i-bands, with typical exposure time of ∼ 20 − 40 minutes for each filter), whereas for the other 5 clusters only one-filter images are available. The colour information is very important to minimise the dilution effect by cluster member galaxies (see also Broadhurst et al. 2005). In addition, the colour information allows to define a sample of member galaxies from galaxies around the red sequence locus in the colourmagnitude diagram, and to compare the member galaxy distribution with the mass distribution constrained from lensing observables. Thus, in this paper we restrict our analysis to the 25 clusters with colour information. As a fiducial background galaxy population, we adopt “red+blue” galaxy sample defined as faint galaxies with colours redder or bluer than the red-sequence by properly chosen offsets. As shown in Okabe et al. (2010), the dilution effect by member galaxies appears to be quite small for this source galaxy sample. The effective source redshift, which is defined such that the

lensing depth Dls/Ds (Dls and Ds being the angular diameter distances from the lens to the source and from the ob- server to the source, respectively) at that redshift becomes equal to the mean lensing depth averaged over the background galaxy population, is estimated by matching magnitudes and colours of the background galaxies to COSMOS photometric redshift catalogue of Ilbert et al. (2009). See Okabe et al. (2010) for more detailed descriptions.

2.2 Weak lensing data

Weak lensing distortion measurements enable a direct reconstruction of the projected 2D mass distribution (Kaiser & Squires 1993). However, the non-local nature of the mass reconstruction indicates that reconstructed mass densities between different pixels in the mass map are highly correlated. Thus one has to take account of the full covariance matrix of the mass map in order to extract proper information from the mass map, as done in Oguri et al. (2005) and Umetsu & Broadhurst (2008). In this paper, we avoid this complication by working directly on the 2D distortion maps.

We use the weak lensing distortion measurements presented by Okabe et al. (2010). For each cluster field, the reduced distortion, gα = γα/(1 − κ), was estimated from the shape of each source galaxy by analysing the

Subaru/Suprime-cam images based on the algorithm described in Kaiser et al. (1995). Throughout this paper we employ the convention that the Greek subscripts denote the two components of distortion, e.g., α = 1 or 2. The lensing distortion is not measurable from individual galaxy image due to the dominant intrinsic shape noise. Instead the signal is measurable in a statistical sense, i.e., by averaging the galaxy shapes over a sufficient number of galaxies. Given the deep Subaru images and the distortion strength in a cluster region, the angular resolution of weak lensing distortion is an arcminute scale. In order to conduct the full two-dimensional analysis of distortion fields, we employ the pixelised data of background galaxy shapes. Specifically, the reduced shear in the l-th pixel (its angular position θl) is estimated as where θi denotes the angular position of the i-th source galaxy and the summation runs over source galaxies con- tained within the l-th pixel. Following Okabe et al. (2010) we use the weighting of each galaxy shape such that a galaxy whose shape is more reliably measured is assigned a larger weight. The weight wi for the i-th galaxy is given by with a = 0.4 and σg(i) being the uncertainty of shape measurement for each galaxy (see Okabe et al. 2010). As stated above, the dominant source of distortion measurement error is the intrinsic galaxy shape. The shape noise in each pixel is estimated as

For all the clusters we adopt a grid size of 1′ × 1′ (we tried smaller grid sizes and found that the results are almost unchanged; see also Appendix A). We do not use four innermost grids (2′ × 2′ box) for the fitting, because source galaxies are obscured by dense distribution of cluster member galaxies especially in the central region, and also because our assumption of single source redshift may become inaccurate near the cluster centre due to the fewer sampling of source galaxies. Moreover, the weak lensing approximation breaks down near the cluster centre. Although the field-ofview of Subaru/Suprime-cam is 34′ × 27′, we conduct our fitting only in a 20′ ×20′ region (20′ corresponds to physical transverse sizes of 2.2h−1 Mpc for z = 0.15 and 3.8h−1 Mpc for z = 0.3, respectively) centred at the BCG, which roughly corresponds to virial radii of clusters in our sample, in order to reduce the projection effect, i.e. the effect of different structures along the same line of sight, which is more prominent in the boundary region where the cluster lensing signal is very weak. However, it should be noted that this restriction of fitting region little affects the final results, as we will discuss later in more detail. Thus the total number of girds used for the fitting is Npixel = 396. An example of our weak lensing shear map is given in the left panel of Figure 1.

For the range of angular scales we use for the fitting, the measured distortion field is nearly the shear field, gα ' γα, which we will simply assume in the following analysis.

2.3 Cluster mass model

We assume that the cluster mass distribution can be described by a single halo component with its radial profile being so-called NFW density profile (Navarro et al. 1996, 1997). The spherical NFW model is fully specified by two parameters, the halo concentration and mass parameters. Although the density profile was obtained from N-body simulation by spherically averaging the halo mass distribution,

4 M. Oguri et al.

Figure 1. Left panel: An example of our weak lensing measurement for A2390. The size of Each panel is 12′×12′. The stick in each 1′×′1 pixel shows the distortion field estimated from background galaxy images contained within the pixel, where a background galaxy image is deformed along the stick direction, and the length is proportional to the shear amplitude. The shear field in this panel is smoothed with a Gaussian with the full width at half maximum of ' 1.6′ for illustrative purpose. Overplotted is the surface mass density map reconstructed from the weak lensing shear measurement (see Okabe et al. 2010). North is up and East is left. Right panel: The shear field predicted by our best-fit elliptical NFW model (see also Figure 2 and Table 1), while the contours are the isodensity map. The best-fit ellipticity of the projected mass density is e ≡ 1 − b/a = 0.598.

we can construct an elliptical lens model simply by introducing an ellipticity in the isodensity contour. Specifically, we adopt the following mass model in our analysis:

where κsph(r) is the radial convergence profile for the spherical NFW profile (e.g., Bartelmann 1996). Here the coordi- nate origin is taken as the halo centre (xc, yc). The halo ellipticity e is related with the major (a) and minor (b) axis lengths of the isodensity contour as e = 1 − b/a. Throughout the paper we adopt the coordinate system with the x- and y-axes being aligned with West and North respec-

tively. With this coordinate system the position angle θe is measured East of North. The lensing shear is computed by solving the two-dimensional Poisson equation whose source term is given by the convergence κ(x,y), as described in Schramm (1990). We note that this elliptical model includes a triaxial halo model which better describes haloes in N- body simulations than the spherical model (Jing & Suto 2002; Kasun & Evrard 2005; Allgood et al. 2006), because the convergence map of a triaxial halo has elliptical isodensity contours when projected along arbitrary directions (Oguri et al. 2003; Oguri & Keeton 2004).

In summary, the elliptical NFW model is specified by 6 model parameters:

Unless otherwise specified, we adopt the virial overdensity

∆vir ' 110 (with respect to the critical density of the universe) computed using the spherical collapse model to define the mass Mvir and the concentration parameter cvir.

2.4 2D weak lensing fitting

The two-dimensional pixelised distortion field described in § 2.2 can be compared with the two-dimensional mass model in § 2.3 in order to constrain properties of halo mass distribution including the halo ellipticity. In this paper we employ the χ2 fitting given as αβ,kl where the indices α and β run over the two components of distortion (α,β = 1, 2) and the indices k and l denote

the pixel position (k,l = 1,

, Npixel). The matrix C de-

notes the error covariance matrix and C−1 is the inverse matrix (see below).The best-fit model parameters are found by minimizing the χ2 value given the distortion data.

For the covariance matrix we include two contributions: the intrinsic ellipticity noise, which is the primary source, and the cosmic shear contamination arising from large-scale structures at different redshifts, but along the same line of sight. The covariance matrix is expressed as

Direct measurement of dark matter halo ellipticity 5

Table 1. Results of weak lensing analysis for 25 clusters, obtained by fitting the 2D shear maps with predictions of the elliptical NFW model. Note that the degree of freedom of the χ2 fitting is 786. The coordinate origin in the fitting is taken as the BCG position for each cluster. The position angle θe is defined East of North. Each error indicates the 1σ error marginalised over other parameter uncertainties.

Name z χ2 xc yc Mvir cvir e θe (arcsec) (arcsec) (1014h−1M ) (deg)

a Removed for the analysis of cluster ellipticities, because the elliptical NFW model does not give an acceptable fit to the measurement (see text for more details).

The intrinsic ellipticity noise is expected to be uncorrelated between different galaxies. Therefore, for the pixelised map, the shape noise covariance matrix has only diagonal terms:

shape ] αβ,kl where δK αβ or δK kl denotes the Kronecker delta function, and σg is given by equation (3). On the other hand, the covariance due to large-scale structure is given in Dodelson (2004) as (see also Hoekstra lss ] where ξ is the cosmic shear correlation functions, and we have assumed that ξ is given as a function of the length of vector connecting the two points θk and θl due to the statistical isotropy of the Universe. Specifically, the shear correlation functions constructed from combinations of two shear components are expresses as where r = |θk − θl|, φ is the position angle between the coordinate x-axis and the vector θk − θl, and ξ++ and ξ× denote the tangential and cross component shear correlation functions (e.g., Bartelmann & Schneider 2001). For a given cosmological model, ξ++ and ξ× are computed once the nonlinear mass power spectrum Pδ(k) and the mean source redshift are given. We assumed the concordance ΛCDM model (Komatsu et al. 2009) and used the fitting formula in Smith et al. (2003) to compute the nonlinear mass power spectrum.

Note that we ignored the error contribution arising from structures surrounding a cluster of interest. The surrounding structures cause the so-called two-halo contribution to the shearing effect on background galaxies (Sheldon et al. 2007a). However, we have checked that the two-halo term contribution is negligible over angular scales we have considered.

The dimension of the covariance matrix is given by twice of the number of pixels in the distortion map. For our case,