Lensing and X-ray mass estimates of clusters (SIMULATION)

[Abridged] We present a comparison between weak-lensing (WL) and X-ray mass estimates of a sample of numerically simulated clusters. The sample consists on the 20 most massive objects at redshift z=0.25 and Mvir>5 x 10^{14} Msun h^{-1}. They were found in a cosmological simulation of volume 1 h^{-3} Gpc^3, evolved in the framework of a WMAP-7 normalized cosmology. Each cluster has been resimulated at higher resolution and with more complex gas physics. We processed it thought Skylens and X-MAS to generate optical and X-ray mock observations along three orthogonal projections. The optical simulations include lensing effects on background sources. Standard observational tools and methods of analysis are used to recover the mass profiles of each cluster projection from the mock catalogues. Given the size of our sample, we could also investigate the dependence of the results on cluster morphology, environment, temperature inhomogeneity, and mass. We confirm previous results showing that WL masses obtained from the fit of the cluster tangential shear profiles with NFW functionals are biased low by ~ 5-10% with a large scatter (~10-25%). We show that scatter could be reduced by optimally selecting clusters either having regular morphology or living in substructure-poor environment. The X-ray masses are biased low by a large amount (~25-35%), evidencing the presence of both non-thermal sources of pressure in the ICM and temperature inhomogeneity, but they show a significantly lower scatter than weak-lensing-derived masses. The X-ray mass bias grows from the inner to the outer regions of the clusters. We find that both biases are weakly correlated with the third-order power ratio, while a stronger correlation exists with the centroid shift. Finally, the X-ray bias is strongly connected with temperature inhomogeneities.

effects of masses along the line of sight to the clusters. However, it is still not well established how safely the hydrostatic equilibrium approximation can be made (Rasia et al. , 2006. As the highest mass concentrations in the universe, galaxy clusters are the most efficient gravitational lenses on the sky. Their matter distorts background-galaxy images with an intensity that increases from the outskirts to the inner regions. Strong distortions occur in the cores of some massive galaxy clusters, leading to the formation of "gravitational arcs" and/or to the formation of systems of multiple images of the same source. Weak distortions, which can be only measured statistically, are impressed on the shape of distant galaxies that lie on the sky at large angular distances from the cluster centers (e.g. Bartelmann & Schneider 2001). Both these lensing regimes can be used to map the mass distribution in galaxy clusters. Gravitational lensing can directly probe the cluster total mass without any strong assumptions on the equilibrium state of the lens. Further, mass profiles can be measured over a wide range of scales, from 100 kpc out to the virial radius. However, lensing measures the projected mass instead of the 3D mass and it is sensitive to projection effects, such as triaxiality and additional concentrations of mass along the line of sight. Given the pros and cons of each method, we can conclude that lensing and X-ray are complementary in many ways. In particular, the comparison of these two mass estimates can greatly help in improving the accuracy of the measurements and understanding the systematic errors.
Numerical simulations provide a unique way to investigate the performance of the lensing and X-ray techniques for measuring the mass profiles of galaxy clusters. Several studies were performed in the past which made use of relatively simple descriptions of galaxy clusters and simulation set-ups, but still were able to assess some fundamental limits of these techniques and possibly suggest improvements (see e.g. Metzler et al. 2001;Clowe et al. 2004;Becker & Kravtsov 2011;Rasia et al. 2004;Piffaretti et al. 2003). In the last years, using the increasing number of observational constraints and profiting of the huge increment of computational efficiency, the simulations have become even more sophisticated and can now include a large number of realistic and important features. These improvements regard both the description of the physical processes determining the evolution of the cosmic structures (see Borgani & Kravtsov 2009, for a review) and the interface between simulations and observations. In particular, few pipelines have been developed which produce simulated observations of the numerically modeled clusters at different wavelengths (Gardini et al. 2004;Nagai et al. 2007;Rasia et al. 2008;Meneghetti et al. 2008;Heinz & Brüggen 2009). These pipelines can simulate observations with a variety of existing and future instruments and include most observational noises which typically affect and limit real measurements. Thus, they are ideal for testing data reduction pipelines (Rasia et al. 2006;Mazzotta et al. 2004;Nagai et al. 2007).
In  (M10, hereafter) we combined our optical simulator, SkyLens (Meneghetti et al. 2008), with our Xray one, XMAS (Rasia et al. 2008;Gardini et al. 2004), to study the systematic effects in mass measurements encountered following standard lensing and X-ray analysis. In that work, we used three simulated clusters and study them along three independent lines-ofsight. In this paper, we extend that work to a much larger sample. We consider 60 mock optical and X-ray images (20 independent clusters for three orthogonal lines-of-sight). Throughout the paper the quoted errors correspond to 1σ level.
The paper is structured as follow. Section 2 contains a short review of the results obtained in previous numerical studies, especially in M10. Section 3 presents a description of the simulated clusters. Section 4 and Section 5 describe the lensing and the X-ray simulation pipelines and the methods of analyses. We present the results in Section 6, where we first discuss the lensing and X-ray mass estimates individually. We show how the bias and the scatter of the mass measurements depend on the cluster morphology and environment in Section 7. Finally, we discuss our results in Section 8.

Strong lensing
In M10, we used the parametric code Lenstool to construct mass models from the multiple image systems detected in synthetic Hubble-Space-Telescope (HST) observations. In the region where strong lensing constraints were found (within the Einstein radius), the mass profiles recovered agree with the input mass distributions with an accuracy of a few percent. Similar results were obtained by Jullo et al. (2007) testing the performances of Lenstool with lens models produced using semi-analytical methods. The stronglensing models are constructed by combining a main halo component and additional massive clumps associated to star clumps (the galaxies of the cluster). Fundamental, in this process, is the modeling of the central galaxy, BCG (Donnarumma et al. 2011(Donnarumma et al. , 2009Comerford et al. 2006). M10 demonstrated that a wrong parametrization of the BCG leads to a severe under-or over-estimate of the strong-lensing masses extrapolated at large radii. Indeed, when the central galaxy was excluded during the creation and the analysis of the synthetic optical images, both the bias and scatter were largely reduced. Discrepancies were seen already at R 2500 , a radius which is typically 2-3 times larger than the Einstein radius. The parametrization of the BCG is also important for a more realistic estimate of the lensing cross-section: its presence increases the strong lensing signal up to a twenty percent in cluster size haloes (Giocoli et al. 2011) .
The necessity of having an accurate model for the BCG in order to extrapolate the strong lensing mass at large radii makes the comparison between strong lensing and X-ray mass estimates highly uncertain. Indeed, X-ray emission from the central region (∼ 70 − 100 kpc) is often excluded from the X-ray analysis because more difficult to model (see e.g. Vikhlinin et al. 2006).

Weak lensing
The weak-lensing analysis is based on the measurement of the shape of galaxies in the background of the clusters, whose ellipticity can be used to estimate the shear produced by the lens. Details on the weak-lensing analysis can be found in Section 4. In M10, we found that fitting the reduced tangential shear profile with a Navarro-Frenk-White (NFW, Navarro et al. 1997) functional or using the aperture mass densitometry produce quite similar results. The measured projected mass is accurate at level of ∼ 10% for those clusters that do not show any massive substructures nearby. Two lens planes presented massive clumps just outside the virial radius of the cluster. This dilutes the shear tangential to the main cluster clump even at smaller radii. As a consequence, the mass profiles of these two lenses resulted to be severely under-estimated. Such problem affects the methods where the shear is measured tangentially. Instead, we tested that techniques which combine strong-and weak-lensing , such as those by Cacciato et al. (2006) and by Merten et al. (2009) are not influenced. To reconstruct the lensing potential the latter method uses an adaptive grid (see also Bradač et al. 2005;Diego et al. 2007;Merten et al. 2011) and naturally incorporates the effects of substructure. As result, the scatter in the projected mass measurement is reduced and limited to 10%.
The main causes of substantial scatter in the deprojected masses are triaxiality and presence of substructure. Under the standard assumption of spherical symmetry, three-dimensional masses are over-or under-estimated, depending on the orientation of the cluster major axis with respect to the line-of-sight. For systems whose major axis points toward the observer, masses are typically over-estimated. The opposite occurs for clusters elongated in the plane of the sky. In M10, the resulting scatter of our sample is of order ∼ 17%. For clusters with substructure, the unknown location of the substructure along the line of sight also makes the three-dimensional mass estimate highly unsure.
More recently, Becker & Kravtsov (2011) used a large number of simulated halos extracted from a large cosmological box to discuss the accuracy of weak lensing masses measured by fitting the cluster shear profiles with NFW models. Their results are consistent with ours. Given the large size of their sample, they significantly probe that weak-lensing masses measured by fitting the tangential shear profiles are biased low, concluding that the NFW model is actually a poor description of the actual shear profiles of clusters at the radii used in the fitting. At the radius enclosing an overdensity of 500 times the critical density of the universe, R 500 , they found that the bias amounts to ∼ 10% for both clusters at z = 0.25 and z = 0.5. They varied the integration length to see the dependence on large-scale structure on the deflection field. The scatter found using our integration length (i.e. 20 h −1 Mpc) is comparable to ours.
Within the integration depth we chose, the large-scale structure can be considered correlated. If the integration length is larger, we will include also uncorrelated structures. Their effects on the weak lensing mass estimates have been discussed in detail in several papers: uncorrelated structures introduce a noise in the mass estimates and their contribution to the total error budget is comparable to the statistical errors (Cen 1997;Metzler et al. 1999;Hoekstra 2001Hoekstra , 2003White & Vale 2004;Hoekstra et al. 2011). Becker & Kravtsov (2011), specifically, showed that the scatter in the weak-lensing masses of low-mass halos increases more than for highmass halos as a function of line-of-sight integration length because the high-mass halos generate more shear than the low-mass halos. The large scale structure has different impact depending on the redshift of the lenses and on the depth of the observations (Hoekstra 2003).
Finally, Becker & Kravtsov (2011) also discussed how the scatter and the bias changes under varying number density of background sources, n g . They show that as the number density increases, the shape noise contribution (due to the intrinsic ellipticity of the sources) to the scatter decreases and eventually becomes subdominant with respect to the intrinsic scatter in weak-lensing mass measurements. For clusters at z = 0.25, the total scatter on M 500 changes from ∼ 37% for n g = 10 gals arcmin −2 to ∼ 25% for n g = 40 gals arcmin − 2. As for the bias, they found that fitting the NFW functional form within R 500 can reduce the bias by ∼ 5%.

X-ray
Regarding the X-ray analysis, we tested two different approaches in M10 that we dubbed the backward and forward methods. The backward procedure assumed a priori a functional form for the mass (such as NFW), spherical symmetry and hydrostatic equilibrium (Eq. 26 in M10): − Gµm p n gas M tot (< r)/r 2 = dP/dr = d(n gas × T )/dr (1) where G is the gravitational constant, µ = 0.59 is mean molecular weight in a.m.u., m p is the proton mass, k is the Boltzmann constant, and n gas and T are the gas density and temperature profiles. These are estimated at once by geometrically de-projecting the measured X-ray surface brightness and temperature data. The 3D temperature is computed following the recipe by Mazzotta et al. (2004). More details can be found in (Ettori et al. 2002;Morandi et al. 2007, R06, M10). The forward method, instead, uses a complex parametric formulae to fit the projected surface brightness and temperature profiles. Subsequently, the analytic 2D expressions are de-projected assuming sphericity and finally the total mass is computed through Eq.1. The two methods are based on the same basic hypothesis: spherical symmetry and hydrostatic equilibrium. They differ for the quantity they analytically parametrized. The first method imposes a fixed mass profile (usually NFW) while the second uses parametric formulas for the surface brightness and the temperature distributions (see Sec.5.2). In this way, the forward approach has smoother radial profiles to be derived, but also more parameters. The two procedures consistently reconstruct both the total and the gas masses, as we demonstrated. For this reason, here, we limit our X-ray analysis to the forward method (presented with more detail later on the paper). The X-ray masses were shown to systematically underestimate the true mass of the simulated clusters by 5-20% with an average bias of 10% between R 2500 and R 200 and a scatter of 6%. The gas masses reconstructed were usually 5% higher than the true ones within the region with sufficient signal. Thanks to the high exposure time used (500 ks) and the field of view of the images, we compared the mass profiles up to R 200 . M10 results were similar to the findings of Nagai et al. (2007). In the same fashion, the authors created mock X-ray images of 16 objects simulated with an adaptive mesh code. The exposure time was 1 Msec and the field of view selected extended well beyond R 200 . Processing the images, they followed the forward method. In their whole sample, the average difference between the total mass and the X-ray derived mass was 16% at R 500 with a scatter of 9%, reducing to 13% ± 10% for regular systems. Rasia et al. (2006) studied a smaller sample of five objects. We follow the backward method assuming different parametrization for the total mass: β model either isothermal or with polytropic temperature profile, NFW and the model presented by Rasia et al. (2004). Under the condition of perfect background subtraction, we found an averaged bias of 23% at R 500 and 20.6% at R 2500 . The causes of the bias were double: the neglect contribution of the gas bulk motions to the total energy budget and the temperature bias towards lower values of the X-ray temperatures. The contribution of the last factor was confirmed by Piffaretti & Valdarnini (2008) who analyzed more than 150 SPH-simulated clusters. Both papers found a temperature bias of 10-15% (see also Rasia et al. 2005). Ameglio et al. (2009) pointed out the direct correlation between this bias and the cluster mass (or temperature): the bias is higher in most massive systems because they have a larger spread in temperature.

SIMULATIONS
Our analysis is based on 20 simulated clusters identified at z = 0.25, all having virial mass M vir > 5 × 10 14 h −1 M at that redshift, and each observed along three orthogonal projection directions. These clusters belong to the set of radiative simulations presented by Fabjan et al. (2011), whose initial conditions have been described in details by Bonafede et al. (2011). The Lagrangian regions around each of these clusters have been identified within a low-resolution N-body cosmological simulation, that followed 1024 3 DM particles within a box having a comoving side of 1 h −1 Gpc. The cosmological model assumed is a flat ΛCDM one, with Ω m = 0.24 for the matter density parameter, Ω bar = 0.04 for the contribution of baryons, H 0 = 72 km s −1 Mpc −1 for the present-day Hubble constant, n s = 0.96 for the primordial spectral index and σ 8 = 0.8 for the normalization of the power spectrum, thus consistent with the CMB WMAP7 constrains (Komatsu et al. 2011). Within each Lagrangian region mass resolution is increased following the Zoomed Initial Condition (ZIC) technique (Tormen et al. 1997). Resolution is progressively degraded outside such regions, so as to save computational time, while preserving a correct description of the large-scale tidal field. Within the high-resolution region, it is m DM = 8.47 × 10 8 M h −1 and m DM = 1.53 × 10 8 M h −1 for the masses of the DM and gas particles, respectively.
Simulations have been carried out using the TreePM/SPH GADGET-3 code, a newer and more efficient version of the GADGET-2 code originally presented by Springel (2005). A Plummer-equivalent softening length for the computation of the gravitational force in the high-resolution region was fixed to = 5 h −1 kpc in physical units at redshift z < 2, while being kept fixed in comoving units at higher redshift. As for the computation of hydrodynamic forces, we assume the SPH smoothing length to reach a minimum allowed value of 0.5 . Our simulations include metal-dependent radiative cooling and cooling/heating from a spatially uniform and evolving UV background, according to the prescription presented by Wiersma et al. (2009). Following the star-formation model by Springel & Hernquist (2003) gas particles whose density exceeds a given threshold value are treated as multi-phase particles, where a hot ionized phase coexists in pressure equilibrium with a cold phase, which is the reservoir for star formation. We also include a detailed description of metal enrichment from different stellar populations, using the model originally described by Tornatore et al. (2007). The effect of SN feedback is included through the effect of galactic winds having a velocity of 500 km s −1 .
The cluster significant radii (R 2500 , R 1000 , R 500 , R 200 , R vir ) 16 and the corresponding masses are listed in Table 1. To compute these quantities we center on the minimum of the potential well as done in Rasia et al. (2011). In the following, we will refer to these numbers as the true or the intrinsic values. To simulate their lensing effects on a population of background sources, we process the halos using our well tested optical simulation pipeline SkyLens (e.g. Meneghetti et al. 2008;Meneghetti et al. 2011, and M10). Here, we briefly summarize the basic steps toward the realization of the simulated images and refer the reader to those papers for further details.
We begin selecting particles falling into a cylinder centered on the cluster and having its width and depth set equal to 10 and 20 h −1 Mpc, respectively. This ensures to include in the simulation the effects of filaments apart from the cluster and of additional mass clumps that could produce addition shear signal. Since we are focusing on high resolution re-simulated clusters, we do not include the effects of un-correlated large-scale-structures. As matter of fact, the importance of matter along the line of sight is fairly small for rich clusters at intermediate redshifts, like those in our sample, provided that the bulk of the sources are at high redshift compared to the cluster (see Section 1).
We project the mass distribution (i.e. the selected particles) on a lens plane at the redshift of the cluster, z L = 0.25. For each cluster in the sample we derive three lens planes, corresponding to the projections (named 1, 2, and 3) along the three axes of the simulation box. The final number of lens planes used in this study is 60. This is a factor of ∼ 7 larger than the sample previously investigated in M10.
The deflection field of each cluster is determined by tracing a bundle of 4096 × 4096 light-rays from the observer position through the lens plane (see M10 for the description of the tree-code). The final deflection matrix is used to further trace the light rays toward the background sources, allowing us to reconstruct their distorted images. In short, the code uses a set of real galaxies decomposed into shapelets (Refregier 2003) to model the source morphologies on a synthetic sky. In the current version of the simulator, the shapelet database contains ∼ 3000 galaxies in the z-band from the GOODS/ACS archive (Giavalisco et al. 2004) and ∼ 10000 galaxies in the B,V, i, z bands from the Hubble-Ultra-Deep-Field ( HUDF) archive (Beckwith et al. 2006). Most galaxies have spectral classifications and photometric redshifts available (Benitez 2000;Coe et al. 2006), which are used to generate a population of sources whose luminosity and redshift distributions resemble those of the HUDF.
SkyLens allows us to mimic observations with a variety of telescopes, both from space and from the ground. For this work, we simulate wide field observations, on which we carry out a weak lensing analysis, using the SUBARU Suprime-Cam. All simulations include realistic background and instrumental noise. The galaxy colors are realistically reproduced by adopting 22 SEDs to model the background galaxies, following the spectral classifications published by Coe et al. (2006).
Compared to M10, we use here a different setup. First, we assume an exposure time of 2000s in the I-band, which is a factor of three shallower than in M10. This is aimed at testing the weak lensing analysis under more realistic conditions. Second, we use real stars observed with SUBARU to model the PSF. The PSF model is characterized by a FWHM of ∼ 0.6 . M10 used an isotropic gaussian PSF instead. For all lens planes, we produce wide-field images covering a region of 2400 × 2400 around the cluster center. This allows us to measure the shear signal up to a distance of ∼ 3.5 h −1 Mpc at z = z L , well beyond the virial radius of any cluster in the sample.

Weak-lensing analysis
The weak lensing measurements are done using the standard Kaiser-Squires-Broadhurst (KSB) method, proposed by Kaiser et al. (1995) and subsequently extended by Luppino & Kaiser (1997) and by Hoekstra et al. (1998). The galaxy ellipticities are measured from the quadrupole moments of their surface brightness distributions, corrected for the PSF, and used to estimate the reduced shear under the assumption that the expectation value of the intrinsic source ellipticity vanishes.
In this study we use the publicly available pipeline KSBf90 by C. Heymans 17 to process our images and measure the shear fields. The final galaxy catalogs are constructed by selecting only the galaxies with signal-to-noise ratio S/N > 10 (as provided by SExtractor, Bertin & Arnouts 1996), half-light radius larger than 1.15 times the PSF size, and reduced shear |g| < 1. Given the above mentioned exposure time and seeing conditions, the effective number density of galaxies in the final shear catalogs is ∼ 17 arcmin −2 . In observations exploiting a depth similar to our simulated images, the number of sources available for the lensing analysis may be smaller because the light emission from cluster galaxies (not included in these simulations) are potential contaminants. These non-lensed galaxies bias low the lensing signal if they are accidentally included in the shear catalogues. Color based techniques (e.g. Medezinski et al. 2007Medezinski et al. , 2010 allow to separate the foreground and the background galaxy populations efficiently, but, to be conservative, several sources which may have a dubious classification are usually excluded from the lensing catalogs. Among them several background galaxies. In M10, we considered several methods to measure the total mass using the observed shear field. As discussed above, we found that the most precise mass measurements are obtained by combining weak and strong lensing non-parametrically (see e.g. Merten et al. 2009). The disadvantage of this approach is that it is very expensive both in terms of time needed to carry out the analysis and in terms of data requirements. The identification of the strong-lensing features, used to constrain the model in the inner region, usually requires deep and high-resolution Hubble-Space-Telescope imaging. Moreover, strong-lensing clusters are relatively rare and known to be affected by many biases Hennawi et al. 2007). Fitting the tangential shear profiles with functionals describing the cluster density profiles is a very common and easy alternative to measure the mass (e.g. Clowe & Schneider 2002;Hoekstra et al. 2000;Jee et al. 2005;Dahle 2006;Paulin-Henriksson et al. 2007; Bardeau et al. 2007;Oguri et al. 2009;Kubo et al. 2007;Romano et al. 2010;Umetsu et al. 2011;Zitrin et al. 2011). Further, this method can be applied to clusters down to relatively small mass limits and in absence of strong lensing features.
Here, we assume that the density profiles of clusters are well described by the Navarro-Frenk-White profile (Navarro et al. 1997), where ρ s and r s are the characteristic density and the scale radius, respectively. The characteristic density is often written in terms of 17 http://www.roe.ac.uk/∼heymans/KSBf90/Home.html the concentration parameter, c 200 = r 200 /r s , as We derive the mass by fitting the one-dimensional reduced tangential shear profile with the corresponding NFW functional (Bartelmann 1996; Wright & Brainerd 2000;Meneghetti et al. 2003). The tangential shear profile is derived from the data by radially binning the galaxies and averaging the tangential component of their ellipticity within each bin. The tangential and cross components are, respectively, defined as (4) The angle φ specifies the direction from the galaxy centroid towards the center of the cluster, which we identify with most bound particle in the simulation. When averaging over many galaxies, the expectation value of the intrinsic source ellipticity vanishes, and the reduced tangential shear is given by g + = + . On the contrary, in absence of systematics the averaged cross component of the ellipticity should be zero.

X-MAS simulations
Before producing the X-ray synthetic catalogue, we have applied the technique described in Appendix A to remove over-cooled particles. Subsequently, our clusters are processed through X-MAS to obtain Chandra mock images. The characteristics of this software package are described in detail in other works (Gardini et al. 2004;Rasia et al. 2008). To create the photon event file, we assumed the Ancillary Response Function (ARF) and Redistribution Matrix Function (RMF) typical of the ACIS-S3 detector aimpoint. We consider the redshift as that of the simulated time frame (z = 0.25) and the metallicity constant and equal to 0.3 solar in respect to the tables by Anders & Grevesse (1989) 18 . The field of view of our images has a side of 16 arcmin. For our cosmology and redshift, this corresponds to 2561 h −1 kpc. All the clusters have their R 500 regions within the field of view (see Table 1), even if some of them at that radius do not emit a sufficient number of photons to allow a precise spatial and spectral analysis (see more in the next session). We account for the emission by all the particles within a depth of 10 Mpc along the line of sight direction and centered on the cluster. The exposure time chosen is 100 ksec. This setting differs from what adopted in M10: the reduction of the exposure time allows a more realistic comparison with observed data. In the final event files, we add a contribution for the galactic absorption by a WABS model with N H = 5 × 10 20 cm −2 . As in M10, we do not include the influence of the background since R06 proved that its net effect is to enlarge the error on the mass estimates without introducing an extra bias. Furthermore, new background models are capable to predict the spatial variation of the Chandra background with an accuracy better than 1% (Bartalucci et al. in preparation).
We note that tools as X-MAS are not suitable to address calibration problematics since the same response files are used both to create and to analyze the data. In this sense, in our analysis we assume a perfect knowledge of the instrument calibration and the results do not depend on the instrument reproduced. In the analysis of real observational data, systematic instrumental uncertainties are highly important, in particular in situation of high statistics (high number of counts). To treat them correctly one needs to include them in the analysis. Lee et al. (2011) have recently provided a bayesian statistical method to tackle this problem.

X-ray analysis
Using CIAO tool (Fruscione et al. 2006) we extract soft band images in the [0.7 − 2] keV band. We apply the wavelet algorithm of Vikhlinin et al. (1998) to identify clumps. These and any major substructure have been masked and excluded from the following analysis. The surface brightness profiles are centered in the X-ray centroid ) and account 15-30 linearly spaced annuli with at least 100 counts. The innermost annulus is selected outside the central 10% of R 500 , the outermost one is always beyond R 1000 and reaches R 500 in the majority of the cases (see Table 3). The radial coverage is comparable to recent observations, some of which extend beyond R 500 (Neumann 2005;Ettori & Balestra 2009;Leccardi & Molendi 2008;Vikhlinin et al. 2009). The temperature profile is calculated in 6-10 annuli spanning over the same radial range of the surface brightness profile. The minimum number of photon per temperature annulus is 1000. The spectra are grouped and fitted by a single-temperature MEKAL model in the XSPEC package (Arnaud 1996). The statistics used is χ 2 and the energy band considered is [0.8-7] keV. In the pipeline the values of galactic absorption, redshift, hydrogen column density and metallicity are fixed equal to the input ones.
To compute the total mass from the X-ray analysis, we follow the "forward" method of M10 (see also, Vikhlinin et al. 2006). The surface brightness and the temperature profiles are fitted by the analytic formulae: Since we exclude the cluster central part from our analysis, we do not model the cooling core region as done in Vikhlinin et al. (2006). This excision is common in both simulations (e.g. M10, Nagai et al. 2007) and observations (e.g. Vikhlinin et al. 2009;Ettori et al. 2010) to avoid, in the former case, the influence of the overcooled central region and, in the latter case, the presence of central active galaxy, cool-core regions, gas sloshing. The 2D analytic formulae are deprojected and the total mass is recovered by assuming hydrostatic equilibrium, that from Eq. 1 it can be written as where T and ρ are the deprojected 3D analytic profiles. Following this procedure, we obtain the X-ray mass that we compare with the true mass of the simulated cluster. The uncertainties in the estimate of this mass, eM X , were obtained through Monte Carlo simulations. In each Monte Carlo realization, surface brightness and temperature profiles were varied within their measured errors. Each time a new mass was then derived. The resulting uncertainty were defined as the standard deviation computed over 100 such realizations.

Weak-lensing mass estimates
Weak-lensing allows us to measure the mass of the cluster projected on the plane of the sky. The NFW analytic formula of the integral along the line of sight of the mass contained in a cylinder is M( The profile parameters r s and ρ s are obtained from fitting the tangential shear profile, as discussed above. Unfortunately, in most cosmological applications the projected mass is not the interesting quantity. Rather, we need to measure the three-dimensional mass . To derive it by de-projecting the two-dimensional mass, one needs to make strong assumptions about the shape of the cluster, which is usually assumed spherical. In this case, the NFW model is given by where y = r/r s . To both the 2D and the 3D masses we associate errors, eM W L,2D and eM W L,3D , computed by propagating the errors on r s and ρ s , as obtained from fitting the tangential shear profiles.
In the following, we show how well we measure projected and de-projected mass profiles of the clusters in our sample. In both cases, the quality of the mass measurement is assessed by means of the ratio Q WL between measured and true mass: Q WL = M WL /M true . The uncertainty on this ratio, eM W L /M true , accounts only for the errors in the weak-lensing mass since the true mass is perfectly known from simulations. The weighted-mean of both the 2D and the 3D weak-lensing bias radial profiles, < Q WL >, is shown by the solid red line in Fig. 1. Its scatter is quantified by the standard deviation of its distribution, and is represented by the shaded yellow region. In formulae: The profiles are plotted in units of R vir . The over-imposed crosses refer to the weighted averaged Q WL computed at the significant radius of each object, with average computed at a radius corresponding to a fixed overdensity ∆, (R ∆,i /R vir,i ), and not over the whole radial profile (R/R vir,i ). They are located at the respective averaged over-density radii. The cross horizontal bars show the dispersion around the radii in units of R vir . The quantitive version of Fig. 1 is reported on Table 2 where we resume all our results. Each value of Q 3D,WL and its corresponding uncertainty is listed in Table B 6 of Appendix B.
Two important conclusions emerge from this analysis. First, the mass measured fitting the reduced tangential shear profile with an NFW functional is biased low. The bias amounts to ∼ 7 − 10% between R 2500 and R 500 and reaches ∼ 13% for larger distances from the cluster center. Second, the scatter in both two-and three-dimensional masses ranges between ∼ 10% at small radii and ∼ 25% at larger radii, being 20% at R 500 .
These results agree with the findings of M10 and Becker & Kravtsov (2011), confirming that the weak-lensing analysis via the KSB pipeline does not introduce significant systematics 6.2. X-ray mass estimates Contrarily to the optical mass measurements, the X-ray mass derivation gives directly the 3D mass profile. Therefore, we can straightly define the ratio between the X-ray mass and the true ones: Q X = M X /M true . Similarly to weak lensing analysis, we compute the weighted average of Q X over the whole sample and the standard deviation of the distribution (Eq. 9). In Fig. 2, we plot the values only within R 500 , which is the radius reached by most of the surface brightness and temperature profiles ( Table 3).The cross shown at R 200 is the result of the extrapolation of the analytic formulae. On the third part of Table 2, we report the weighted-averaged Q X and its scatter. Each Q X value and its corresponding uncertainty is listed in Table B 7 of Appendix B. Fig. 2 confirms some previous findings that we will synthesize here postponing to Section 8 a more profound discussion. The average X-ray mass is consistently underestimating the true mass. The X-ray bias is around 25% at the center and 30-35% at R 500 . The decline in the most external regions is expected since the cluster outskirts present a more dramatic lack of hydrostatic equilibrium (Lau et al. 2009) and a stronger influence of gas clumpiness . The presence of gas clumps affects the X-ray mass determination in two ways: it shallows the surface brightness profile and it cools X-ray temperature (see more on this in the Section 8.) Massive systems, as the ones studied in this paper, are expected to be still accreating and therefore far for an equilibrium state. Moreover, the temperature bins in the external regions, where the temperature profile declines more steeply, are usually larger containing more temperature structure. Finally, the large bias on the most external region has to be taken with caution since it is not the result of a measurement but of an extrapolation. The dispersion around the average is quite small. The standard deviation is less than 10% at all radii apart from R 2500 where it is 12%. These numbers are twice or three times smaller than the gravitational lensing dispersion.  We investigate in this Section the efficiency in reducing bias and scatter on both X-ray and gravitational lensing masses of two selecting criteria. We create different sub-samples determined by the morphology of the X-ray images or by the presence of substructures on their environment.

Masses and X-ray morphology
To limit the impact of the non-thermal processes on the X-ray mass estimates, clusters are often selected on the basis of their appearance. The literature is rich of studies where clusters have been classified into relaxed, or regular, and unrelaxed, or disturbed, because of their X-ray morphology (e.g. Zhang et al. 2008;Vikhlinin et al. 2009). Most of the time, the classification is done "visually", i.e. simply quantifying the regularity of an object from the X-ray image in the soft band. More objective criteria, proposed in the past, are the power ratios, centroid-shift, asymmetry and fluctuation parameters, and hardness ratio. We test all of them and present here our result.
Third order power ratio and centroid shift. Buote & Tsai (1995) suggested to decompose the surface brightness distribution in multipoles. The high order multipoles, usually normalized by the monopole and called power ratios, are used to quantify the contribution of different scale components (asymmetries and substructures) to the surface-brightness power spectrum relative to the large-scale smooth cluster emission. Most information in the power spectrum is contained in the first four multipoles. P 0 is the monopole. The power ratio P 1 /P 0 measures the dipole of the X-ray emission, which is zero if measured with respect to the X-ray centroid. The power ratio P 2 /P 0 measures the ellipticity (quadrupole). The third order power ratio P 3 /P 0 can be used to quantify asymmetries and is the best indicator of clusters with multimodal distributions. Substructures on smaller scales contribute to higher order multipoles.
Another indicator of the dynamical state and of the asymmetry of the X-ray emission is the centroid-shift, i.e. the shift of the surface brightness centroid in apertures of increasing size. This parameter points out the dynamical state of the cluster as well as the asymmetry. Following Poole et al. (2006) and Maughan et al. (2008), we define the centroid-shift as where R max is the radius of the largest aperture, and ∆ i = R c,i − R c,max is the shift of the centroid in the i−th aperture with respect to the centroid in the largest aperture, R c,max . ∆ is the mean value of the various ∆ i and the sum is done over all the N apertures with radii up to R max . In this work we assumed N = 17 apertures with radii ranging between R min = 0.15 × R 500 and R max = R 500 . The third-order power ratio and the centroid shift were shown to be effective in classifying clusters by two recent works by Cassano et al. (2010) and Böhringer et al. (2010). Clusters are located in a rather well defined region in the P 3 /P 0 − w plane: objects with small centroid shift and small P 3 /P 0 are classified as "regular". The majority of them are cool core systems, not very dynamically active and showing absence or very little radio emission. For all these reasons, often, these objects are referred as "relaxed".
In this work, we compute the power ratio, P 3 /P 0 , and the centroid-shift, from the signal of the region within R 500 of the masked images. With this attention, we aim to evaluate the "irregularity" of the actual portion of the image that we use to retrieve the mass.
Both the values and their uncertainties are derived from Monte Carlo simulations. We create 100 new images where the photons are re-distribuited accordingly to a Poisson statistics. We evaluate the estimators in each image. Finally, we extract the medians and the 16 th and 84 th percentile of the Monte Carlo distributions to represent the final values of the morphological estimators and their uncertainties. FIG. 2.-Comparison between X-ray and true masses using the whole sample. The meaning of lines, crosses and shaded regions is the same of Fig.1.
The third order power ratios and the centroid shifts of our sample with their uncertainties are shown in the left panel of Fig. 3. We recognize the region with regular systems by slightly relaxing the criteria adopted by Cassano et al. (2010). In our work, we define a cluster to be regular when w < 0.03 and P 3 /P 0 < 2 × 10 −7 . Our choice is motivated by the fact that our aperture is equal to R 500 , thus larger than the 500 kpc aperture radius analyzed by Cassano et al. (2010). Reducing this radius, they naturally measured lower values of the morphological estimators and, in particular, of the power ratios (Böhringer et al. 2010).
The 17 regular objects are denoted by asterisks in the figure. In most cases, these are systems with small companions or some minor irregularity in the surface brightness map. The full classification is listed in Table 3 below the column "P 3 /P 0 , w". Two extreme examples are represented in the right panel of Fig. 3. The X-ray images of the most disturbed system of our sample (CL 13) are shown on the top panels, while the bottom panels refer to a relaxed cluster (CL 9).
The uncertainties on our parameters are smaller than those of Böhringer et al. (2010) because of the better spatial resolution of Chandra with respect to XMM-Newton (see Böhringer et al. 2010, for a detailed discussion about the influence on spatial resolution or point-spread function). The comparison between their work and that by Cassano et al. (2010), based on Chandra data, confirms this statement. The large exposure time assumed in our mock observations ensures a high counts statistics and, therefore, a further reduction on the uncertainties. If future missions will reproduce the great spatial resolution of Chandra, both power-ratios and centroid shifts will be available with sufficient accuracy for a large number of objects, thus allowing highy detailed studies of cluster morphologies.
We now proceed by checking whether the lensing and true masses biases improve when selecting only the X-ray regular systems. Similarly to Fig. 1, we show in Fig. 4, Q WL as a function of the distance from the cluster center for the sub-samples of regular systems. Quantitative results are listed in Table 2 including those for Q X . The scatter on lensing bias is reduced by 20%-40%. However, the bias itself worsen with respect to the whole sample. Clearly, a selection based on the X-ray morphology is not optimal for lensing purposes. The reason of this behavior can be explained by comparing Table 3 and Table B 6. Among the X-ray regular clusters, there are three images (projection 2 of CL1, CL9, and CL20) whose lensing measurements are severely under-estimated. All of them present in the outskirts of the optical images filaments or falling substructures, that do not have any obvious counterpart in the X-ray images or are lying outside the Chandra field of view. Jeltema et al. (2008) and Piffaretti & Valdarnini (2008) claimed to find a significant trend of the X-ray masses bias on the morphological estimators. We investigate this aspect further by including an analysis of the biases in the weak-lensing mass reconstruction. To this purpose, we considered the absolute values of |(1 − Q 3D,W L )| to evaluate the dependence for any deviation. A linear fit between the mass biases and the morphological estimators has been computed accounting for measurement errors in both variables. The results are reported in Table 4. The centroid shift performs better that the third order power ratio. The slopes of the Q − w relations are always significantly different from zero with the correct sign (negative for the X-ray bias and positive for the weak-lensing deviations). The best fit values are similar to those found by Jeltema et al. (2008) and Piffaretti & Valdarnini (2008). We further quantify the correlation between the mass biases and P 3 /P 0 or w by means of the Pearson correlation coefficient. We find always a negative correlation, being the values between -0.3 and -0.4 for w and around -0.2 and -0.3 for P 3 /P 0 .
In Fig. 5 we present the best combinations: Q X − w for R 2500 , and |1 − Q 3D,W L | − w for R 1000 and R 500 . The top panels, are similar to Fig. 3 where the different colors and symbols refer to different values of the mass bias. The red triangles refer to clusters whose X-ray mass biases, Q X , are within 20% or whose weak-lensing-masses deviations, |(1 − Q 3D,W L )|, are within 10%. In all three top panels, we distinguish no segregation of colors. This implies that a better estimate of the total mass (red triangles) does not necessarily come from regular clusters defined on the basis of P 3 /P 0 and w values. However, for the centroid shift, even if this condition is not be necessary, it is sufficient at all radii: the weighted-average bias for clusters whose centroid shift is lower than 0.3 is 15-20% lower than those with w > 0.06 (see difference on horizontal lines in the central panels). As confirmation, the third order power ratio weakly discriminate between good and bad estimates.  projection 1 projection 2 projection 3 cluster R 500,X P 3 /P 0 , w I.p. O.p. R 500,X P 3 /P 0 , w I.p. O.p. R 500,X P 3 /P 0 , w I.p. O.p. CL1 √ Finally, we consider the distance between the centers, ∆C used in the our X-ray and weak-lensing analysis. In the former case, we considered the X-ray centroid, while in the latter we used the center of the BCG, which is also coincident with the minimum of the DM-potential well. A shift between the two centers might testimony a recent merger able to separate the two components, as for the Bullet cluster (Markevitch et al. 2004).
All these parameters have been compared with the X-ray and weak lensing bias measures. For the weak lensing case we consider both the values of Q 3D,W L and of |Q 3D,W L − 1| to evaluate general deviation from the true mass. As done for the centroid-shift and the third order power ratio, we measured the correlation between all the parameters and the biases. We found that none of the new parameters is more strongly related to the bias than the centroid shift and the third order power ratio. On the opposite, their Pearson correlation coefficient is always smaller than 0.1 in absolute values.

Masses and cluster environment
As shown in M10, the scatter in WL mass measurements is due to substructures and triaxiality. The combination of X-ray and lensing data may help to further identify the most spherical systems or to correct the mass estimates for triaxiality effects (Morandi et al. 2011;Sereno et al. 2010). However, these techniques are still model-dependent and subject to strong assumptions. Furthermore, they require a certain amount of handling of the data which cannot be applied to a large sample of objects. Substructures, instead, may be more easily identified.
To classify clusters on the basis of the level of substructures in their surroundings, we visually inspect the projected mass maps of each cluster in our sample. We identify the objects whose environment is poor of substructures within a region of 5 h −1 Mpc around their centers. At first, we look directly at the intrinsic density map from the simulations. We find 10 cluster projections that match this criterium. We named this sample "I. poor" (Intrinsically poor). Often a regular cluster is not part of the poor environment class FIG. 5.-Distribution of the clusters power ratio P 3 /P 0 and centroid shift in dependance of the mass biases. In the first column the bias shown refer to X-ray mass estimates computed within R 2500 . In the second and third column, we plot the absolute deviation of the 3D weak-lensing mass bias computed at R 1000 and R 500 , respectively. With red triangles and blue asterisks we show the weakest and strongest mass biases. The intermediate situations are shown with yellow rhombi and cyan squares. In the central and lower panel, we separate the dependance of the bias by each morphological estimator. The horizontal lines represent the weighted average of the bias for particular values of the parameters (w below 0.03 and above 0.06, P 3 /P 0 below 2. × 10 −7 and above 10 −6 .) 0.79 ± 0.01 -0.88 ± 0.18 -0.73 ± 0.005 -0.20 0.10 -0.68 ± 0.0049 -0.27 ± 0.091 -Q X , P 3 /P 0 0.51 ± 0.03 -0.04 ± 0.005 -0.67 ± 0.02 -0.007 0.003 -0.71 ± 0.020 0.006 ± 0.003 + |1 − Q 3D,W L |, w 0.10 ± 0.02 1.30 ± 0.47 + 0.11 ± 0.02 1.54 ± 0.50 + 0.12 ± 0.03 1.68 ± 0.56 + |1 − Q 3D,W L |, P 3 /P 0 0.32 ± 0.08 0.03 ± 0.01 + 0.32 ± 0.10 0.02 ± 0.02 0 0.28 ± 0.12 0.01 ± 0.02 0 Q X , Q T -0.49 ± 0.15 1.48 ± 0.17 + -0.13 ± 0.10 1.04 ± 0.11 + -0.21 ± 0.10 1.10 ± 0.12 + (see characterizations on Table 3). This is easily explained by presence of substructures outside the X-ray field of view which is limited to ∼ 2.5h −1 Mpc. One of such clusters is the already-mentioned projection 2 of CL9. Despite being X-ray regular (Fig. 3) it shows evidence of filaments in its surroundings causing a strong underestimate of the weak-lensing masses in that projection (see Table B 6). This is the main reason for which some X-ray regular objects show a large weak-lensing bias: they are lying in a rich environment that cannot be detected in the X-ray images. The mass bias of the clusters classified 'I. poor' is reported on Table 2 and showed in the middle panels of Fig. 4 for the weaklensing masses. The exciting result is that the projected true masses are almost exactly recovered. For these systems Q W L deviates from unity by only a few percent for the 2D mass and by less than 5% for the 3D mass. The scatter is strongly reduced especially among the two-dimensional masses, being only of order ∼ 5% over a wide range of radii. This is smaller than that found in M10 combining SL and WL non-parametrically. For the three-dimensional masses, the scatter at the most external radii (R 500 and R 200 ) decreases by ∼ 50−60% with respect to the whole sample and by 20-40% with respect to the systems with regular X-ray morphology. As shown in M10, this residual bias is caused by triaxiality. It may be alleviated by means of introducing a parameter describing the elongation along the line of sight in the fitting model. This requires a combination of different probes, as proposed by Morandi et al. (2011) andSereno et al. (2010), who combine lensing and X-ray data. However, a large uncertainty still remains, due to possible non-thermal pressure in the ICM, which is degenerate with the cluster triaxiality.
As second step, we visually inspect the projected mass maps reconstructed from the synthetic weak lensing observations. The method used for the reconstruction is described in Cacciato et al. (2006) and in Merten et al. (2009), although we do not make use of the SL systems in this test. This approach, even if more "observationally-oriented", is still subjective. Furthermore, the visual classification is more challenging because the resolution smears out possible features and the noise in the optical images reduces the detectability of clumps. Clusters that appear isolated in the reconstructed mass maps are called "O. poor". This classification is also presented in Table 3. The results on the weak-lensing and X-ray mass biases are reported in Table 2 and shown in the bottom panels of Fig. 4. The net improvement in terms of bias and scatter of both the 2D and 3D masses is now much less evident, but still masses are better recovered than in the sub-sample of X-ray regular clusters. In particular, the bias is reduced by about a factor two at R 500 and R 200 , for both 2D and 3D masses.
Notice from Table 2 that the X-ray bias and scatter does not vary substantially in the three samples ( regular, I. poor, and O. poor). This outcome is expected because, in general, X-ray masses have small scatter. The intra-cluster medium is generally more spherical than the DM or the galaxy distribution, especially outside the core . As a consequence, the X-ray method is less prone to triaxiality. This implies that removing/adding a few objects, as long as they are not very disturbed, does not change significantly the result.

DISCUSSION AND CONCLUSION
This paper is an extension of the work of M10. We compared galaxy cluster masses derived from gravitational lensing and X-ray using 20 new massive halos simulated at high resolution including radiative gas physics. Each halo was observed along 3 different lines of sight and located at redshift 0.25. We used an optical and an X-ray simulators, namely Skylens and X-MAS, to build both optical and X-ray mock images mimicking, Subaru and Chandra observations, respectively. To perform the weak lensing analysis, we measured the galaxy shapes using the KSB method and we derived the masses by fitting the tangential shear profiles using NFW functionals. For the X-ray, instead, we used the forward approach described in M10 and derived the mass under the hypothesis of hydrostatic equilibrium. Then, we selected a subsample of regular clusters on the basis of the X-ray morphological estimators (the third moment of power ratios and the centroid shift). We further classified the objects in our sample based on the presence of substructures in the cluster environment. This classification was based on the visual inspection of both the true and the lensing-reconstructed projected mass maps of the systems under investigation.
In the following, we discuss our main results: Q W L and Q X for the whole sample. The weak lensing mass bias is less than 10% within R 500 , and it grows to 13% in the most external region. The X-ray bias is around 25% in the central region and increases to 33% at R 500 . The scatter of the bias is always higher by at least a factor of two for weak lensing than for X-ray mass measurements. The weak lensing bias and its large scatter is caused by presence of substructures in the cluster surroundings and by the triaxiality of the systems. The X-ray bias, instead, is mainly due to the lack of hydrostatic equilibrium, presence of clumps (in the external regions) and temperature dis-homogeneity (see further discussion).
Q and morphological parameters. We evaluate the effectiveness of some morphological estimators for reducing the mass bias. We found that a selection based on centroid shift (w < 0.03) and third order power ratio (P 3 /P 0 < 2 × 10 −7 ) reduces the X-ray bias, especially in the central regions. This selection has no effect on the weak lensing bias itself but it decreases the scatter by 20-40% (see Table 2). Among the different morphological parameters used to identify disturbed morphologies (including asymmetry and fluctuation parameters, two hardness ratios, and the optical-X-ray center offset), the only one that shows a mild correlation with the bias is the centroid-shift. This is true also for the weak-lensing masses bias. In terms of future X-ray missions, an optimal use the centroid-shift to identify "ideal" clusters to maximizing the efficiency of the mass measurements would require an imaging quality comparable to that of Chandra, but over a larger field of view so that the area probed by X-ray and optical observations become comparable.
Q W L substructures. We established already that weak lensing methods based on single model fitting can severely fail to measure the mass of clusters in presence of massive substructures (M10). Working on single objects, the effect of substructures could be minimized by adopting multi-halo fitting techniques , provided that substructures can be clearly identified as peaks in the weak lensing maps. Filtering techniques might also offer a possibility to mitigate the effect of structures perturbing the cluster shear profiles (e.g. Gruen et al. 2011). We verified that removing from the sample those clusters which live in environments rich of substructures allows to minimize both the bias and the scatter in the weak lensing mass estimates. Unfortunately, the identification and the characterization of the substructures is a very difficult task. We will dedicate to another work the study of substructures detectability. Several galaxy surveys are planned for the next years which will scan large portions of the sky (see e.g. The Dark Energy Survey Collaboration 2005; LSST Science Collaborations et al. 2009;Refregier et al. 2010). The data are expected to provide galaxy number densities in the range 15 − 40 gals arcmin −2 , allowing to measure the shear signal of several thousands of clusters. Having deep and sharp observations over a large field of view and with good spatial resolution would make it possible to detect substructure in a more efficient way than what presented in this analysis, thus enabling to virtually identify all the relevant substructures. Detection and mass measurements of sub-structures are already possible in the Coma cluster . Methods based on higher order lensing distortion measurements (lensing flexion) also seem very promising (e.g. Okura et al. 2007;Velander et al. 2011).
Q W L triaxiality. Triaxiality introduces a further scatter and bias in the three-dimensional lensing mass estimates. Even minimizing the impact of substructures, by restricting the analysis to the poor environment clusters, we notice a tendency to underestimate the total mass on average. This is due to the fact that a large fraction of the systems in the sample is mostly elongated on the plane of the sky. Under this circumstance, we expect to under-estimate the mass in the de-projection phase (see e.g. Feroz & Hobson 2011). Conversely, the mass is over-estimated in clusters seen along their major axis. Studies based on simulations showed that clusters forming in a CDM framework are generally prolate systems (Shaw et al. 2006;Allgood et al. 2006;Lau et al. 2011). For such mass distributions there is a larger chance to infer a smaller mass within a given radius from the projected density field and explains both the presence of the 2D bias and its increase after the de-projection. This result depends on the selection of our sample which is mass limited. If our clusters were chosen for their lensing signals, and in particular for strong-gravitational lensing, we would have had more objects strongly aligned along the line of sight (e.g. . Then, the measured masses would have been over-estimated on average. Q X and Q W L and dependence on cluster mass. The mass range in our sample is too narrow to make a robust statistical analysis on the dependence of Q X and Q W L on the cluster masse. However, we can attempt to evaluate this effect by averaging these values on the three most (CL2, CL10, and CL11, all with M 500 > 7.5 × 10 14 h −1 M ) and least massive systems (CL 4, CL18, and CL19, all with M 500 < 3.5 × 10 14 h −1 M ). The main result is that the bias Q 3D,W L of the smallest clusters drops by almost a factor of 5 going from R 2500 ( < Q 3D,W L >= 0.94 ± 0.20) to R 200 (< Q 3D,W L >= 0.83 ± 0.39), while for massive clusters there is no change in the bias, always equal to 0.95. This indicates that at large distances from the cluster center, the mass estimates of the least massive clusters will be affected by a stronger bias compared to those of the most massive systems. This is not surprising, because M10 already showed that, when the tangential shear is used to measure the mass at a given radius, additional sources of shear, like massive substructures, located outside that radius lead to under-estimate the mass. The impact of such perturbers is obviously more significant in clusters of smaller mass. Repeating the same test for the X-ray bias we do find the opposite trend. The bias of larger systems goes from < Q X >= 0.77 at R 2500 to < Q X >= 0.70 at R 500 (with error of order of 0.05), while for smaller objects the bias is constantly equal to 0.74. This can be ascribed to two main reasons: the massive objects are the most disturbed ones and they have a complex temperature structure (Ameglio et al. 2009). Its importance on the X-ray mass determination will be investigated more on the next paragraph.
Q X and temperature distribution. The values of Q X that we find in this work are consistent with what previously found in R06. In M10 we find smaller deviations being the average bias around 10%. There are some differences between this paper and M10 which could affect the analysis, such as five times smaller exposure time and a narrower field of view. However, we believe that most of the difference is not due to the X-ray preparation and analysis but to the physics adopted in the simulation. The simulations analysed in M10 included a description of thermal conduction, which is instead set to zero in the simulations presented here and in those analysed by R06. Before elaborating more on the effect of changing the physical description of the ICM, we remind that the X-ray mass is derived from the hydrostatic equilibrium equation (Eq. 1) where three terms are present: the derivative of the gas density, the derivative of the temperature, and the temperature itself at the radius considered. The over-or under-estimate of the temperature leads to an overor under-estimate on the X-ray mass of identical amplitude. If the temperature structure is spherically homogenous, the measured temperature will be independent on the method used to derived it. However, when the plasma presents temperature structures in the annulus, the X-ray measurements are biased low because the X-ray detectors of Chandra and XMM-Newton have an higher efficiency on the soft band and, thus, weight more colder gas . In presence of inhomogeneity, the X-ray temperature is, therefore, lower than the mass-weighted one. As consequence, the hydrostatic masses computed directly using the intrinsic gas density and mass-weighted temperatures of the simulated clusters are higher than those obtained following entirely the X-ray procedure. This is illustrated in Fig. 6, where we plot the values of Q X reported in Table. 2 in black and a similar ratio related to intrinsic calculation in red. For all clusters, the intrinsic bias is 12.6% ± 5.9%, 18.3% ± 4.5%, and 22.6% ± 5.1% at the three radii, respectively. In the case of regular and poor systems, these values decrease to ∼ 10%, ∼ 15% and ∼ 17% and are in agreement with previous works based only on intrinsic evaluation of the hydrostatic equilibrium mass using the mass-weighted temperature (e.g. Piffaretti & Valdarnini 2008;Jeltema et al. 2008;Ameglio et al. 2009;Rasia et al. 2004;Lau et al. 2009).
To stress more the dependence of the X-ray mass bias on the temperature difference, we plot on the right panels of the Fig. 6 Q X versus Q T = T X /T MW for the three X-ray-significant over-densities. The uncertainty associated with the temperature bias are equal to 1 σ error obtained from the spectroscopic analysis. The over-plotted line is the best fit to the relation: Q X = A + B × Q T calculated excluding very disturbed objects (blue asterisks) and considering the errors in both coordinates. The values of B, for R 2500 , R 1000 and R 500 , are largely different from zero (see Table 4) and very close to 1. The biases are strictly correlated one to the other showing Pearson coefficients equal to −0.7, −0.8, −0.7.
These results also suggest an explanation for the difference between our results and those Nagai et al. (2007): their average bias is consistent with M10 and lower than the result of this paper. Their analysis and procedure are almost identical to ours, therefore, our results can be compared straightforwardly, despite the instrument-setting is slightly different (they consider reproduce synthetic observation of ACIS-I and not ACIS-S3). Indeed in both our works, the choice of the instrument is not as important as in real observation. Our mock images are, by construction, not affected by any calibration issues since we assume the same response files when generating and analyzing the images. In this condition of perfect calibration, we could have some differences only if the shape of the responses will be extremely different in shape (not in normalization), e.g. if an instruments weight a lot more plasma at 5 keV in respect to the plasma at 8keV. However, as demonstrated by Mazzotta et al. (2004) in their Fig. 8 the detectors on board of Chandra and XMM-Newton give consistent answers in this respect. Furthermore, Nevalainen et al. (2010) in a study focused on real data showed that both the Coma Cluster and A1795 have broad band temperatures consistent within the 3% statistical uncertainties and broad band fluxes may differ by up to 2-3%. We, therefore, exclude the possibility that the difference between our results and Nagai et al. (2007) can be due to the instrument chosen.
As for the difference with respect to M10, as already mentioned the analysis presented in that paper was based on simulations which included the effect of thermal conduction, with conduction coefficient set to one-third of the Spitzer value. As discussed by Dolag et al. (2004), thermal conduction in hot clusters is quite effective not only in removing cold blobs, but also in making the FIG. 6.-Left panel: Q X (black), intrinsic Q X assuming mass-weighted temperature (red), and Q 3D,W L (green) for all clusters, only regular, and only in poor environment. The values and scatter of Q X and Q 3D,W L are reported in Table 2, while those of the intrinsic hydrostatic-mass bias are listed on the text (Section 8). Right panel: relation between Q X and Q T = T X /T MW at the three significant radii: R 2500 (top panel), R 1000 (central panel), and R 500 (bottom panels). The lines represent the best-fit in the form of Q X = A + B × Q T excluding the very disturbed clusters (blue asterisks). thermal structure of the ICM more homogeneous. This leads to an increase of the spectroscopic temperature and, therefore, of the hydrostatic mass.
As for the comparison, with Nagai et al. (2007), while their simulations do not include thermal conduction, they are based by an Eulerian Adaptive Mesh Refinement hydrodynamical scheme. As discussed by several authors (e.g. Agertz et al. 2007;Mitchell et al. 2009;Vazza et al. 2011), Eulerian hydrodynamics leads to a more efficient mixing of gas entropy. Therefore, one expects in Eulerian simulations that low-entropy gas residing in high-density clumps becomes more efficiently mixed, than in SPH simulations, to the high-entropy ICM. Again, this should result in a more homogenous temperature distribution, with a smaller bias introduced in the estimate of X-ray temperature.
To summarize our X-ray mass bias derived intrinsically assuming hydrostatic equilibrium and mass-weighted temperature are comparable to previous results. Following the X-ray approach, instead, we confirm the findings of R06 but we find a stronger bias in comparison to M10 and Nagai et al. 2007. The further ∼ 10-15% is caused by temperature inhomogeneity (see also Piffaretti & Valdarnini 2008;Ameglio et al. 2009). This result highlights that, while hydrodynamical simulations are powerful tools to understand biases in observational mass estimates, a detailed assessment of this bias (e.g. its precise value) is still uncertain depending on the physical processes included in the simulations and on their numerical description. In this respect, we remind that none of the recent theoretical studies on X-ray mass biases includes the effect of feedback from AGNs.
Lensing and X-ray masses comparison. We compare the gravitational lensing masses with the X-ray masses following the same fitting procedure described in Mahdavi et al. (2008) and M10 to whom we address the curious reader for a detailed description. In brief, for each over-density ∆, we define a parameter, a ∆ , as M X = a ∆ × M 3D,W L . The error associated correspond to 68% confidence level. In Table. 5, we report our results for all clusters and for the relaxed sample. We compared the lensing masses (green crosses in Fig. 6) to both the X-ray-strictly derived masses (black crosses in Fig. 6) and the intrinsic ones (red crosses in Fig. 6). For reference, we include the values found in M10 and in two observational papers: Zhang et al. (2008) for the Locuss sample and Mahdavi et al. (2008) for the CCCP sample.
Our intrinsic results are consistent within the errors with the observational data, especially for regular clusters. Our "observed" X-ray to weak lensing mass ratios are, instead, consistent only with Mahdavi et al. (2008) for R 1000 and R 500 . In all the other cases, our ratios are lower than the observed ones. This is due, once again, to the temperature in-homogeneity detected in our simulated clusters.
This last comparison has three consequences. On one hand it could be that SPH codes generate more temperature structures, i.e. deviation from a spherically symmetric temperature distribution, than present in real clusters (Sijacki et al. 2011). Unfortunately, current X-ray telescopes cannot provide detailed temperature maps for large sample of clusters with enough spatial resolution to confirm or dismiss this hypothesis. Indeed, previous observational techniques to evaluate temperature structures Zhang et al. 2009) have been applied to a limited number of nearby objects. Increasing the size of the samples for which detailed observational studies are carried out would provide a unique opportunity to test the reliability of simulations in describing the complexity of the ICM thermal structure. Despite this precise match cannot be done at present, some observations have already evidenced a significant level of gas clumpiness around the virial radius of some clusters (Simionescu et al. 2011, Vikhlinin et al. in prep.), leading to the conjecture that the plasma may be not completely thermalized at R 500 . Another possible way to explain the mis-match between our results and the observations is that weak-lensing masses of real clusters suffer of some additional bias. For example, as mentioned above, the contamination of shear catalogues by foreground galaxies may bias low the masses. The impact of this source of contamination and its possible correction using color-selection techniques is currently under investigation. Finally, it is likely that the inclusion of AGN should significantly reduce the temperature in-homogeneity. In that case, the X-ray temperatures will be closer to the mass-weighted ones, leading to 5%-10% difference between X-ray and weak-lensing masses (see Table5). A common characteristics of hydrodynamical simulations of galaxy clusters, is that gas associated to merging galaxies or small groups keep their identity for longer time within the hot ICM atmosphere. This is especially true for simulations including radiative cooling, without an efficient feedback mechanism, and for simulations based on SPH, which provide a rather inefficient mixing between high-and low-entropy gas phases. These structures are revealed in the X-ray images as compact sources of dense gas, which are characterized by a strong emission. The majority of them is located in the vicinity of the cluster center. A consequence of the overcooling problem is also that the central galaxy shows an extremely powerful X-ray peak. Since the presence of these features is mostly due to unaccounted physical processes, the standard procedure is to exclude them after their identification through a wavelet-detection algorithm (Vikhlinin et al. 1998) and to excise the central 15% of R 500 (R06, Rasia et al. 2008;Nagai et al. 2007).
This procedure to remove high-density cold gas clumps is rather time-consuming, in particular for bright clusters, as those we are analyzing here, since they host a large number of sub-clumps. As a more efficient approach, we explore a method to exclude a priori particles which have short cooling time. These particles, which have high density and low temperature, can be identified in a well-determined region of the phase diagram of temperature, T p , and gas densitie, ρ p . Empirically, we found that all clusters in our sample have a separated phase of cooling particles which satisfy the following condition: (A1) The normalization factor, N, depends on the mass of the cluster being higher for the more massive systems. In our sample it does not vary substantially since the mass range considered is relatively small. Therefore, we consider a fixed value of N equal to 3 × 10 6 with density expressed in units of (g cm −3 ) and temperature in keV. This normalization is conservative, meaning that all the particles belonging to a genuine hot phase of all our cluster lye above the relation set by this limit. The exponential factor, 0.25, depends on the polytropic index. If we assume that the pressure, P ≡ constant × T × ρ gas , is related to the gas density through the polytropic equation: P ∝ ρ γ gas , we obtain T ∝ ρ γ−1 gas . The polytropic index, γ, physically can vary between γ = 1 (isothermal plasma) and γ = 5/3 (adiabatic gas) (see Sarazin 1988, for a review). For simulated clusters, the polytropic index lies between 1.15 and 1.25 (Ascasibar et al. 2003;Rasia et al. 2004;Ostriker et al. 2005;Bode et al. 2009), constraining the exponential factor in the equation A1 between 0.15 and 0.25. This prescription to remove gas particles belonging to a cold phase is visually illustrated in the top panels of Fig. A7. The right panel shows all the particles centered on CL1, in a field of view of 16 × 16 arcmin 2 and within 10 Mpc h −1 along the line of sight (projection 1). The corresponding X-ray synthetic soft energy image is in the central panel. The black line in the right panel corresponds to the cut applied in our sample (see Eq. A1). Only the particles above that line contribute to the X-ray emission shown in second soft map (right panel), which is the image used for the analysis presented in this work. The small green circles are the regions identified and removed with the wavelet algorithm. The annuli in the two images have radii equal to 0.15×R 500 and R 500 . The difference between the two emission maps is evident. Our cutting technique allow us to clean the image from about 30 small blobs while does not affect the emission of the clusters or of the large sub-clumps. To test this last statement, we analyzed both synthetic images (considering and excluding the overcooled particles) for three clusters in our sample. In the six resulting images, we run the wavelet algorithm to identify the X-ray peaks and exclude them. Subsequently, we extract both surface brightness and temperature profiles and confirm that the cutting technique does not introducing any bias. In the bottom panels of Fig.A7 we plot the comparison for CL1, projection 0. Both the surface brightness profiles and the temperature profiles are consistent with each other.

B. MEASURED MASSES
In Tables B6 and B7 we report the values of Q 3D,W L and Q X computed at different radii: R 2500 , R 1000 , R 500 , and R 200 . Uncertainties are computed following the procedure described in Sections 6.1 and 5.2, respectively. The condition expressed in eq.A1 is represented by the black line. In the bottom panels we compare the surface brightness profiles of the two soft X-ray images once excluded the regions marked (on the left) and the temperature profiles on the right.