Spatial matter density mapping of the STAGES Abell A901/2 supercluster field with 3D lensing

We present weak lensing data from the HST/STAGES survey to study the three-dimensional spatial distribution of matter and galaxies in the Abell 901/902 supercluster complex. Our method improves over the existing 3D lensing mapping techniques by calibrating and removing redshift bias and accounting for the effects of the radial elongation of 3D structures. We also include the first detailed noise analysis of a 3D lensing map, showing that even with deep HST quality data, only the most massive structures, for example M200>~10^15 Msun/h at z~0.8, can be resolved in 3D with any reasonable redshift accuracy (\Delta z~0.15). We compare the lensing map to the stellar mass distribution and find luminous counterparts for all mass peaks detected with a peak significance>3\sigma. We see structures in and behind the z=0.165 foreground supercluster, finding structure directly behind the A901b cluster at z~0.6 and also behind the SW group at z~0.7. This 3D structure viewed in projection has no significant impact on recent mass estimates of A901b or the SW group components SWa and SWb.


INTRODUCTION
Our current understanding of cosmology describes an intricate web of dark matter spanning the Universe dictating when and where galaxies form. These large-scale structures of dark matter have long been simulated in large N-body computations (Springel 2005), but due to the nature of dark matter, linking this feature of dark matter theory directly to observations has been a challenging task.
Historically, the presence of dark matter as an essential ingredient of the standard model of cosmology was deduced indirectly from the rotation of galaxies or the stability of galaxy clusters as the largest known cosmologi-cal objects. Nowadays, the dark matter density is inferred from the break-scale in the galaxy power spectrum, and by combining fluctuations in the cosmic microwave background, the luminosity-density relation of supernovae Type Ia, the Lyman-α forest, or the bulk motions of galaxies. For an overview, see for example Dodelson (2003) and Peacock (1999) and references therein. Nevertheless, a direct detection of dark matter particles as such is still an open issue.
Weak gravitational lensing offers a unique way to chart dark matter structures in the Universe as it is sensitive to all forms of matter. Weak lensing has been used to map the dark matter in galaxy clusters (see for example Clowe et al. 2006) with high resolution reconstructions recovered for the most massive strong lensing clusters (see for example Bradač et al. 2006). Several lensing studies have also mapped the projected surface mass density over degree scale fields (Gavazzi & Soucail 2007;Schirmer et al. 2007;Kubo et al. 2009) to identify shear selected groups and clusters. The minimum mass scale that can be identified is limited only by the intrinsic ellipticity noise in the lensing analysis and projection effects. Using a higher number density of galaxies in the shear measurement reduces this noise, and for this reason the Hubble Space Telescope (HST) is truly unique for this area of research. Compared to ground-based observatories, observations from space more than triple the number density of galaxies that can be resolved with a significantly enhanced imaging quality, permitting high resolution reconstructions of dark matter in the field Heymans et al. 2008) and the study of lenses at higher redshift.
The lensing distortion experienced by a source galaxy depends on the total projected surface mass density along the line of sight (l.o.s.) and a geometrical factor that increases with source distance. The highest redshift galaxies therefore experience, on average, the strongest shear distortions. This redshift dependence can be used to recover the full three-dimensional gravitational potential of the matter density (Taylor 2001;Bacon & Taylor 2003) and has been applied to the COMBO-17 survey in Taylor et al. (2004) and the COSMOS survey in Massey et al. (2007). All these 3D mass reconstruction methods however require the use of a prior based on the expected mean growth of matter density fluctuations. Without the inclusion of such a prior, Hu & Keeton (2002) show that a Wiener filter approach is unable to reasonably constrain the radial matter distribution, even for densely sampled space-based quality lensing data.
An advancement of the method of Hu & Keeton (2002) is presented in Simon et al. (2009), hereafter STH09, which details a non-parametric reconstruction method that directly reconstructs the 3D matter density distribution. This is in contrast to previous methods which have mapped the gravitational potential. In broad terms one can view the gravitational potential as a heavily smoothed realisation of the density field. As an example, a density map of a single point mass will have a broad corresponding gravitational potential map that radially decays as 1/r in 3D (as ln r in a 2D on the sky projection) from the position of the source. Directly reconstructing the density field allows one to see the sources that produce the smoother gravitational potential through Poisson's equation and therefore this is ultimately what we wish to do. Historically, the 3D potential field has been mapped because the signal-to-noise was higher. Applying a Laplacian to a noisy potential field appeared unfeasible on cluster scales. However, the work of STH09 showed that the Wiener filter method could also be applied in this regime to recover the matter-density field. Since the matter-density is the quantity we are interested in, and because the mass sheet degeneracy leads to unconstrained quadratic terms in the potential field, we now focus purely on the reconstruction of the 3D matter-density field. Since this is the first application of this method, we shall pay close attention to the changes in signal-to-noise.
Another new, computionally faster singular value decomposition technique has recently been presented in Van-derPlas et al. (2011). This paper proposes to use Eigenmodes in the reconstructed 3D mass density field, rejecting modes above a given noise level, thereby achieving a radial smoothing without a Wiener filter. On a deeper level, this is what a Wiener filter is doing as well, as it suppresses density modes according to their signal-to-noise, for which the prior is used to determine the expected signal in the reconstruction.
The Wiener filter employed by STH09 uses a ΛCDM prior to optimally smooth and reconstruct the 3D map. The smoothing does, however, have an unwanted impact as it (1) spreads out structures in radial (redshift) direction, (2) shifts on average the centre of structures to different redshifts (here-after referred to as z-shift bias) and (3) scales the overall amplitude of matter density fluctuations. We note here that while the radial spread and scaling of the amplitude of recovered perturbations was noted and detected in the 3D potential mapping approach, shifts in the position of structure were first pointed out in STH09 and therefore not accounted for in these studies. The information about how the true matter distribution is smoothed in the resulting map is clearly defined by a radial point spread function (p.s.f). This radial p.s.f needs to be taken into full account for any correct interpretation of the map. In addition, as also outlined in STH09, the redshift resolution of all 3D lensing reconstruction will be low owing to the weak response of lensing to variations in lens redshift encoded in the broad lensing kernel.
The STH09 method has been optimised for observational data as it utilises the lensing information from all galaxies, whether they are with or without a photometric redshift estimate. It also has inbuilt tests for the effects of systematic errors in the lensing measurement, using an E/B mode decomposition. Crucially it also produces a 3D noise map so the significance of structures in the resulting maps can be assessed.
In this paper, we present the first application of the 3D matter density non-parametric reconstruction method, analysing the HST STAGES survey of the Abell 901/902 supercluster field, hereafter A901/902. We choose this field to test this method as it contains four known massive structures at z ∼ 0.165 detected at high significance in a 2D map (Heymans et al. 2008) and a fifth known structure in projection at z ∼ 0.46 (Taylor et al. 2004). In addition, STAGES is one of the largest HST surveys with one of the highest quality lensing data sets for this type of analysis.
The 3D mass map is presented as a raw Wiener 3D reconstruction that reveals the expected low redshift resolution, based on our knowledge of the radial p.s.f and lensing kernel. Assuming the local maxima along each l.o.s. are a good estimate of the position of the peak density, we present 3D visualisations of the dark matter density in this field. Highlighting the most likely positions of structures in this way produces remarkable images that can be compared to other published 3D lensing reconstructions. However we voice strong caution in using any of these visualisations for scientific purposes as they present a very noisy estimate of the density field and contain no error information. Moreover, structure along the same l.o.s. can overlap and thereby merge in the Wiener map to one single structure at a combined average redshift.
For scientific purposes, we propose the use of raw Wiener reconstructions. We identify interesting l.o.s. across the field from a 3D cross-correlation analysis of mass and light and then fully take into account the radial biasing and smoothing of the Wiener reconstruction to determine each lens redshift probability distribution. We then compare this probability distribution to over-densities of stellar mass along the same line l.o.s. and discuss our findings for the 3D matter distribution in the STAGES field.
The paper is set out as follows. In Sect. 2 we review the STAGES data set and the 3D reconstruction method and present our results in Sect. 3. We discuss our findings and conclude in Sect. 4. As this analysis places new requirements on the accuracy of the weak lensing observations, the Appendix details our tests of the observational analysis that ensure the results presented in this paper are robust.

3D MASS RECONSTRUCTION OF STAGES
In this section, we focus on the application of the STH09 non-parametric 3D matter density reconstruction method to the STAGES A901/902 supercluster field, referring the reader to STH09 for details on the full theory that underpins this method. This field has been the subject of several weak lensing analyses (Gray et al. 2002;Taylor et al. 2004;Schirmer et al. 2007;Heymans et al. 2008), investigating different aspects of the dark matter environment of this dense field. The most recent analyses have focused on data from the Space Telescope A901/902 Galaxy Evolution Survey (STAGES, Gray et al. 2009). In this analysis, weak lensing distortions have been measured for over 60, 000 galaxies (∼ 65 galaxies per square arcmin; Heymans et al. 2008) of which ∼ 15% of the sample have accurate photometric redshift and stellar mass estimates from the COMBO-17 survey (Wolf et al. 2004). Lensed galaxies without photometric redshifts estimates are incorporated using magnitude bins for which the redshift distribution is estimated using Schrabback et al. (2007).
The 3D density reconstruction presented in this paper is the most complex lensing analysis of these data and places new demands on the accuracy and reliability of the weak lensing shear measurement and the accuracy of the photometric redshifts. In the Appendix, we therefore document the steps taken to confirm the reliability of the weak lensing shear catalogue that forms the basis of the 3D reconstruction.
At the very end of this section, in subsection 2.9, we give an executive summary of the main results to highlight the new results of this section.

3D mapping in a nutshell
In weak gravitational lensing, gravity-induced shape distortions of galaxy images are linearly related to the projected matter distribution along a line of sight (l.o.s.) direction θ. Let us start with a simple case in which we bin the complex ellipticities of all sources in a full shear catalogue onto the same regular grid. The grid covers the field-of-view of the survey, the mean ellipticity of the ith grid cell is denoted by ǫ(θi); θi is the l.o.s. direction of the ith cell. We arrange the binned ellipticities as one vector (1) In the working regime of weak gravitational lensing, the convergence κ(θ) of light bundles is the projection of the 3D matter density fluctuations δ(r ⊥ , χ), where r ⊥ is a comoving distance vector perpendicular to a fiducial l.o.s., and χ is the comoving radial distance : Here c/H0 is the Hubble radius, f k (χ) the angular diameter distance of χ, Ωm the mass density parameter and a(χ) the scale factor at distance χ. The function W (χ) expresses the lensing kernel averaged for a distribution of sources in χ: pχ(χ) is the probability density of source distances. As Ansatz, we chop the matter distribution δ into a set of N lp matter slices or lens planes. Within each matter slice, δ is assumed to be constant along a l.o.s. θ, so that the jth matter slice is represented by δ (j) (θi), the average δ in direction θi (lens plane), and the width ∆χj = χj+1 − χj of the slice; χj defines the boundaries of adjacent slices. By means of we arrange the grid values of δ of every slice as a vector. Note that we assume a flat sky where every lens plane is a tangential plane to the celestial sphere in the direction of the fiducial l.o.s..
Likewise, we arrange the κ(θ) into a vector: With this Ansatz, Eq. 2 can now be written as linear projection of the δ-slices into the κ grid: where In the last step of Eq. 6, we denoted for convenience the weighed sum of δ (i) 's by Qδ.
In order to relate the convergence on the grid to the binned ellipticities of our survey, we have to transform convergence to shear in the next step. On a flat sky, the transformation is expressed by the well-known Kaiser & Squires (1993) relation, which is written down for a regular grid here (Hu & Keeton 2002): where we have for the convergence-to-shear transformation matrix with θij = θi −θj and A being the solid angle of a grid pixel. Note that for the 2D vector θ we are using the commonly employed complex notation. So far, nothing has been said about the noise, which is necessarily included in ǫ owing to the intrinsic shape of the sources and the non-uniform distribution of sources resulting in occasional cases of grid cells with zero sources (for example in a masked region): Of n only the statistical properties are known: it is assumed to be uncorrelated to δ, has a vanishing statistical average and the covariance N = nn † . In the case of missing sources, the noise n(θi) associated with the ith grid cell is infinite.
In 3D lensing, we have a survey that can be split into various subsamples of sources with distinct and known redshift distributions. Now each source subsample provides a different grid ǫi, which every time is a projection of the same matter cube δ but with a different set of weights Q i and noise grids ni. Therefore, employing a compact notation, we now have The map making algorithm utilises this relation to define an estimator for δ given ǫ and the noise properties N. Although Eq. 13 has no exact solution due to the unknown noise component n and the non-invertible Pγκ (mass sheet degeneracy), one can find an optimised solution, i.e. a linear filter H of ǫ, that minimises the average deviation between the estimated minimum variance δmv = Hǫ and the true δ, i.e. δmv − δ 2 = min.
(By x 2 = x t x we mean the Euclidean norm of x.) The solution to this optimisation problem is in our case (e.g. Zaroubi et al. 1995): where R ≡ PγκQ, and α is an additional constant, usually between α ∈ [0, 1], that allows to change the signal-to-noise in the resulting map by regulating the smoothing of the filter (Saskatoon filter; Tegmark 1997). The signal covariance S ≡ ǫǫ † − N has to be given to the filter as prior and is computed for a particular fiducial cosmological model. The practical challenge is the computation of δmv = Hǫ within a reasonable time, see the Appendix of STH09 for details.
To evade the problem of the mass-sheet degeneracy, we restrict ourselves to solutions for which the average of δmv vanishes for every redshift slice individually.

Reconstruction prior and grid parameters
The prior can be implemented in the form of a radial or a transverse filter; a transverse filter varies as a function of angular wave number ℓ (STH09). For S, we choose a radial Wiener filter where F (x) = 2x −1 J1(x) is the Fourier transform of a tophat window, Jn(x) a Bessel function of first kind and P (i) δ (ℓ) is the power of matter fluctuations, projected along the width of the ith matter slice. We approximate our square grid pixels by circular pixels of identical area A = πΘ 2 s . As a fiducial matter power spectrum we use the approximation of Smith et al. (2003) and Limber's equation for converting 3D to 2D power on the sky. The fiducial cosmology is a ΛCDM model (adiabatic fluctuations) with Ωm = 0.24, Ω b = 0.0416, ΩΛ = 1 − Ωm and H0 = h 100 km s −1 Mpc −1 with h = 0.732. The normalisation of the matter fluctuations within a sphere of radius 8 h −1 Mpc at redshift zero is assumed to be σ8 = 0.76. For the spectral index of the primordial matter power spectrum we use ns = 0.96. These values are consistent with the third-year WMAP results (Spergel et al. 2007).
The smoothing of the Wiener filter can be regulated by rescaling the noise covariance of the data via α: For α = 0 Wiener filtering is switched off, i.e. the prior S cancels out at the expense of the lowest signal-to-noise of the reconstruction but for the benefit of an unbiased estimate of δ. For α = 1 the Wiener filter boosts the signal-to-noise to a maximum at the expense of heavy smoothing in the radial direction along with a strong z-shift bias and loss in redshift resolution. We find α = 10 −3 provides a good balance to maximise signal-to-noise and minimise bias, which we fully account for below. Note, a radial prior does not produce transverse smoothing so an additional transverse smoothing is required for every lens plane. For post-reconstruction smoothing we use a Gaussian kernel of σrms = 40 arcsec, which is the size of several grid pixels.
The number of matter lens planes used here is N lp = 22, covering a range from z = 0 . . . 1.15 with an associated slice width of ∆z = 0.05. There is one exception for the closest plane which has ∆z = 0.1. The width of the slices is finer than the radial resolution we can possibly hope for. Nevertheless, the large number of lens planes was chosen to avoid numerical instabilities in the reconstruction. The lens planes are gridded with an angular resolution of 128 × 128 pixel 2 , one pixel has a square size of 16.34 arcsec. The vector δ hence has 128 2 × N lp ≈ 3.6 × 10 5 degrees of freedom.

The radial point spread function
As outlined in STH09, the reconstruction algorithm does not provide an unbiased representation of the true radial matter distribution. The bias can however be quantified.
Ignoring the noise in a particular reconstruction, on average the map δmv will be a smoothed image of the true matter distribution δ: where . . . is the ensemble average over all source noise realisations while keeping the actual matter distribution δ constant. Hence the matrix B = HPγκQ holds all information on the smoothing involved in the reconstruction process. In the optimal case, one would have B = 1. Here, we do not have B = 1, though, but still B contains all information needed for a quantitative assessment of the resulting map.
Unfortunately, as discussed in STH09, simply applying the inverse of B to δmv, even if B was regular, would undo the Wiener filter, leaving behind a heavily oscillating reconstruction with a significant noise level because low signalto-noise modes in the map would no longer be suppressed.
Here, the full matrix B has of the order of 10 5 × 10 5 = 10 10 elements. It turns out, however, that in practise most of its elements are negligible because the devised Wiener filter smoothes the matter distribution only in the radial direction.
Concretely, consider that we have just one point source overdensity inside the entire 3D density field δ, represented by one unity-valued grid pixel in direction θi on the jth lens plane. To lift the mass-sheet degeneracy, the algorithm assumes that the density contrast averaged over every lens plane vanishes. Hence, to stay inside the space of considered δ-configurations all other pixels of the jth lens plane have to be set to −(127 × 128) −1 for our (128 pixel) 2 grids; pixels of all other lens planes i = j are zero. We denote this special density contrast configuration by e the vector n denotes the shear noise in a particular realisation.
Essentially the original matter point source e (j) i is smeared out in the l.o.s. direction θi, but not in the transverse. In STAGES the radial smearing will be roughly unchanged when varying θi inside the same lens plane i. We checked this by computing Be (j) i for a set of N los = 10 × 10 l.o.s. directions for all lens planes covering the inner 80% of the field. Therefore, instead of keeping the full Be where [A, θi] denotes the linear operation that takes only elements from the 3D grid A that are along the l.o.s. direction θi and arranges them as an N lp element vector in order of ascending lens plane redshifts. Note that Be (j) i is calculated in practise by applying the reconstruction algorithm to PγκQe (j) i so that in total N 2 los × N lp = 2200 individual reconstructions were performed. Fig. 1 shows the resulting radial p.s.f. for a selected set of lens planes. Radially located peaks in the true matter distribution will have a cigar-like appearance in the map due to the evident width of the p.s.f. The spread will be reduced for better signal-to-noise in the data, but will always be present to some degree when a Wiener filter is employed.

Z-shift bias
An obvious application of the Wiener map is to estimate the redshifts of 3D peaks by the local maximum of the cigar density profile or to compare the redshift estimate to other  maps of the galaxy number density or luminosity. However, as can be seen from Fig. 1, generally the maximum of the p.s.f. does not coincide with the true redshift of structures, which we hereafter refer to as z-shift bias. To determine the extent of z-shift bias we use Fig. 1 to determine where, on average a point source mass is recovered in the reconstruction. The result is shown in Fig. 2, which displays the redshift of a p.s.f. local maximum (apparent redshift) and error in relation to the original point source redshift in e (j) i . We use the data points in this figure to construct a new Wiener map free of the z-shift bias in the following way: We take a set of new lens planes with "true redshifts" and fill for every θi these planes with values from the uncorrected Wiener map at the "apparent redshift" given by Fig. 2.

Quantifying noise and systematics
To assess the S/N of features inside the map, we divide map values by the variance observed in 10 3 map noise realisations before z-shift bias correction. Noise maps are realised by randomly rotating the orientations of sources and rerunning the reconstruction algorithm. During this process, we also average the noise covariance N of density pixels along the same l.o.s.. Dividing by the noise variance inside a lens plane changes the shape of the p.s.f. and the z-shift bias, see open squares in Fig. 2. In the case of STAGES, the S/N-map shifts all p.s.f. peaks essentially to the middle of the map (in redshift), which can be explained by larger levels of noise variance at low and high redshifts. Due to the strong impact of the z-shift bias, the S/N-map is in essence thus only a 2D map stretched out over the entire redshift range. For that reason, we will use it in the following only as a 2D projection to quote the maximum S/N expected in a l.o.s. direction.
To test the impact of systematic errors in the lensing shear measurement, B-mode maps are produced. These are obtained by rotating the source ellipticities by 45-degrees followed by a repetition of the reconstruction process. Ideally, these maps should be pure noise maps as a gravitational tidal field is unable to generate a B-mode distortion by shearing galaxy images. Instead, they arise from systematics related to the source shape p.s.f. correction procedure, the intrinsic clustering of sources (Schneider et al. 2002) and intrinsic galaxy alignment correlations of physically close sources or shape-shear-correlations , and references therein). Moreover, as discussed in STH09, matter fluctuations outside the field-of-view but close to the field-of-view may be capable of generating large-scale B-modes within the reconstruction volume. Therefore, an entirely vanishing level of B-modes is never guaranteed per se. To partly suppress this effect, we take a reconstruction field-of-view somewhat larger than the actual coverage of the survey and later remove the area not covered by observations. This allows unknown sources of shear to be outside the field-of-view and to be properly accounted for in the model reconstruction.
The peak and trough statistics in Sect. A3 shows that the B-mode map is consistent with noise realisations, whereas the E-mode map is inconsistent with noise at high confidence as it should be. This leads us to the conclusion that any systematic errors are low and within our statistical errors.

Redshift uncertainties
The main drawback of the Wiener method is that the radial spread of the p.s.f. hampers the ability to attach a redshift to a cigar-shaped profile in the map. This would be true even if a p.s.f. was very narrowly peaked as the uncertainty of the estimate still needs to be determined by accounting for the noise in the map, which is not directly visible in the Wiener estimate for the 3D matter density contrast.
What really determines the redshift accuracy is the ability to distinguish a radial profile p (i) from p (j) , Eq. 19, despite of the noise superimposed onto the average profile in a concrete realisation.
To get a handle on the principle redshift accuracy, let us consider a simple case where the map density profile in the l.o.s. direction is dominated by one single structure on the ith lens plane, i.e. δ (i) = δ0p (i) + n; δ0 is the grid pixel density contrast and n the noise vector in that particular realisation. To pin down the redshift of the structure from the noisy vector δ (i) means whether we can statistically distinguish, given the data, a model (a) in which the mass is presumed to be on lens plane j from the correct model (b) where the mass is on lens plane i. Accepting model (a), we would estimate the density δ ′ 0 = λδ0 at redshift zj by minimising We are assuming a Gaussian noise model with covariance N = nn t here. Since the source noise is roughly Gaussian and the Wiener filter linear this is a valid assumption here. The minimum solution for the previous equation is so that the χ 2 of the best-fit of model (a) has which on average for all realisations n is Therefore, to distinguish model (a) and model (b), we have to decide when χ 2 ij is inconsistent with a χ 2 -distribution with N lp degrees of freedom, i.e. the χ 2 -distribution of model (b). Since model (b) hasλmin = 1 and χ 2 ij = N lp , the difference in χ 2 is simply the above ∆χ 2 ij . To get fiducial values of ∆χ 2 ij for STAGES, we model δ0 as the centre of a singular isothermal sphere (SIS) with mass profile (e.g. ) where χ h is the SIS comoving distance and σ 2 v the velocity dispersion of the SIS. The SIS is smoothed it about the centre with a Gaussian kernel of size σrms = 40 arcsec to the Wiener map transverse smoothing scale: Note that we also have to divide the smoothed δsis (SIS density contrast integrated along the l.o.s.) by the width ∆χ h of the matter slice it resides in, to comply with the definition of the mean density contrast in our mapping algorithm, Eq. 6. Based on the Ansatz outlined above, Fig. 3 displays the expected accuracy with which redshifts of structures can be determined from a Wiener filtered 3D map. In this plot, masses are defined as SIS-M200 masses, i.e. the mass contained inside a radius r200 within which the mean density equals 200 times the critical density: For comparison, the two most massive clusters in the field, A901a and A901b, have approximately a mass of M200 ≈ 2 × 10 14 M⊙h −1 (Heymans et al. 2008). The figure adopts the average radial p.s.f. of STAGES. For a STAGES-like survey we find the redshift of a massive galaxy cluster of M200 ∼ 10 15 M⊙ h −1 could still be measured with accuracy of ∆z ≈ ±0.08 (±0.15) at z ∼ 0.6 (0.8), whereas a cluster of mass with M200 ∼ 3 × 10 14 M⊙ h −1 has ∆z ≈ 0.05 at z ∼ 0.2 but already ∆z ≈ 0.15 at z ∼ 0.6. For masses smaller or equal than M200 ∼ 10 14 M⊙ h −1 , the accuracy quickly deteriorates with redshift, at most imposing upper limits on the redshift. This lack of 3D resolution for lower mass halos poses the natural question of how useful a 3D lensing reconstruction is, which we discuss in Section 4.

Detection limit for mass concentrations
In order to gauge the mass detection limit of the STAGES Wiener map, we have produced Fig. 4. For this plot, we place the central pixel δ0(zi) of a SIS at lens plane redshift zi and determine the maximum S/N ratio, along the SIS l.o.s. in the Wiener map; STAGES survey parameters are assumed. For a given zi, M200 is changed until Rsn reaches a desired S/N limit. This limit is easily found due to the scaling M200 ∝ R 1.5 sn , which follows from Eqs. 28 and 27.
For A901a/b at z = 0.165 we expect a S/N between 6 − 7σ in the Wiener map, which indeed is found. We can also see from this plot how well suited the STAGE field is for this type of analysis as we are primarily sensitive to structure at low and moderate redshift for M200 10 14 M⊙h −1 . To reconstruct mass at higher redshift requires a significant structure to lie in the field.

Product Map Mass and Stellar Mass comparison method
For a further analysis, we will compare the lensing mass map to the distribution of stellar mass inside the data cube. The stellar mass map is derived directly from the COMBO-17 galaxy catalogue, which is binned onto the 3D grid used for the mass reconstruction. We compare the two 3D gridded distributions on a pixel-by-pixel basis with the aim to find local matches and discrepancies between the light and mass distribution.
For the mass-light comparison, we restrict the light map to a selected redshift range: where the indices l, k = 1 . . . N lp , l k choose the redshift range in question; ei denotes a vector with N lp elements being zero everywhere except for the ith element. The range is z ∈ [z k , z l ], where zi is the redshift of the ith lens plane.
We understand the lensing map as a convolved underlying matter distribution, The convolution of the mass map is an unavoidable byproduct of the Wiener filter reconstruction. Therefore, in contrast to the light map, for the mass map we do not have the choice of selecting a certain redshift range. We therefore cross-correlate stellar mass with the full lensing mass map. This is reasonable because distinct lens planes are well separated enough in physical distance to have no relevant astrophysical correlation, i.e. for i = j Our analysis tool is a product map where the product is taken of the stellar mass and total mass density profiles in radial direction: which on average is where we define the weights wi ≡ [p (i) ] t ei. From the latter it becomes obvious that Eq. 33 can be made an estimator for the weighted l.o.s. product of mass and light distribution within the selected redshift range [z k , z l ], namely by the product map (PM) which for statistically independent lens planes has the average We use product maps to obtain candidates for l.o.s. directions for which features in the mass map could be related to stellar mass (PM > 0), directions that are significant in the mass map but apparently uncorrelated to the stellar mass map (PM = 0) or even anti-correlated (PM < 0). The candidates will be subject to further modelling of the radial density profile to determine the probability distribution of mass along the l.o.s., Sect. 3.2.
The product method has a limited radial resolution due to the broad radial p.s.f. p (i) . For illustration, consider a l.o.s. that has exactly one non-zero pixel ∆ l in δ l on the ith lens plane. On the other hand, the matter distribution may have the only non-zero pixel ∆m on the jth lens plane, hence at a different redshift. Despite an obvious mismatch between light and mass, this case will contribute to the product map with ∆ l ∆m[ei] t p (j) (Eq. 33), if p (j) is broad enough to reach the ith lens plane.
In Fig. 5, we show the response, or signal, in the product map to a stellar mass point source at discrete redshifts . This is an effect of the broad radial smoothing kernel the Wiener map is subject to. In an ideal world, the solid lines would be strongly peaked about the triangles, i.e. [p (i) ] t e j ≈ δ K ij . Thicker lines belong to higher stellar mass redshifts. For comparison, the dashed-dotted line shows the response we would expect for a purely 2D (convergence) map; the response curve is independent of the stellar mass redshift in this case. z = 0.05, 0.33, ...1.08 (shown with triangles) from a lensing mass spread over the full redshift range 0 < z < 1. Mass along the full redshift range will contribute at some level to the product map, even though it is not associated with the stellar mass point source. This shows that our interpretation of the product map is non-trivial as there is always considerable confusion present for the STAGES map. However it does show that we have relatively little confusion from low redshift lensing mass features to stellar mass redshift range is z 0.5. This is important as it means one can, in principle, distinguish in the product map the supercluster structure at z ∼ 0.165 from structures in the higher redshift regime. The response could be further improved by using less radial smoothing in the Wiener map by means of devising a smaller tuning parameter. This, however, also decreases the signal-to-noise in the product map.
A cross-correlation of a stellar map with a 2D convergence map, on the other hand, always produces the same response irrespective of the stellar mass redshift. This response is shown in Fig. 5 by the dotted-dashed blue (amplitude was rescaled for plotting purposes). Although the cross-correlation with the 3D Wiener map displays some improvement over a cross-correlation with a convergence map, it is unclear how the still broad response can be properly dealt with. For instance, a massive mass peak at low redshift could generate a strong response with a large stellar mass peak at high redshift.
To advance the 3D product map, it is conceivable to define a radial smoothing of the light map in order to optimise the response in the product map by substituting ei by an appropriate vector. In addition, one could properly account for the z-bias so that the maximum response is at the location of the stellar mass peak or trough, which is not exactly the case here (Fig. 5). We do not pursue such a strategy for this work, though, as we are not expecting a significant improvement for STAGES. In particular, we do not use a 3D product map to draw strong conclusions about radial matches of stellar mass and matter mass. Instead, we solely employ the product map in combination with the Wiener S/N-map to identify interesting l.o.s. to be further analysed by a purely lensing method Sect 3.2.
Note that smoothing the light map with the Wiener map p.s.f. before multiplication, i.e. set ei = p (i) , results in a worse response than the one adopted for this analysis; the response would be the product of pairs of Wiener p.s.f. equalling their overlap, i.e. wi = [p (i) ] t p (i) in this case.
Finally, the signal-to-noise of the product map is assessed by dividing Eq. 35 by the variance seen in noise realisations, which are obtained by randomising the positions of galaxies in the M * map. Note that the signal-to-noise levels in the product map can become considerably larger than in a pure lensing map because for uncorrelated noise of the lensing and M * map the signal-to-noise in the individual maps just multiplies.

Executive Summary
We have reviewed the theory presented in STH09, then presented details of the extensions and new analysis tools. These tools were developed for the first application of this method to real data. Here we highlight the new results presented above for different readers.
For readers optimising future survey designs, the most relevant results are contained in Figs. 3 and 4, which show the best 3D mass resolution that can be attained on the considered angular scales with space-based quality data with a galaxy number density of 65 galaxies per square arcmin matched by 10 galaxies per square arcmin with groundbased photometric redshift estimates with a redshift error of σz/(1 + z) = 0.02. Note these photometric redshifts are limited at z < 1.4 in the absence of near infra-red photometry to extend reliably to higher redshifts. This survey characterisation is similar to planned deep imaging areas of space-based surveys with matching photometric redshifts from shallower multi-colour optical surveys such as PanSTARRS, DES or VST-KIDS. These reveal rather relatively poor expectations for attaining mass resolution in 3D using lensing alone for anything other than the most massive structures (at least several 10 14 M⊙h −1 ). It therefore motivates future surveys to invest in sufficiently deep photometric ground-based surveys such as LSST and matching deep near infra-red photometry.
For readers focused on applying the STH09 method to a new data set, Sect. 2.3 and 2.4 detail how to make the nontrivial measure of Z-shift bias and the radial p.s.f from real data that was not discussed in STH09. Fig. 2 also shows how the application of standard 2D mass mapping noise estimate techniques to a 3D mass map biases the redshift measure with the implication that a 3D S/N map, which has not previously been calculated in the literature, is essentially a 2D map for lensing surveys like STAGES. Section 2.8 therefore outlines a new product map method that enables one to compare significant peaks detected in the essentially 2D S/N mass to the 3D stellar mass distribution and the 3D lensing density distribution.
We show that the product map offers theoretical improvements over a cross-correlation of a 3D stellar mass map with a 2D convergence map, but it remains unclear at this point how to clearly exploit this map. Nevertheless, we employ the product map to select interesting l.o.s. to be analysed further in the following.

RESULTS
In this section, we present the results of our 3D reconstruction of the STAGES A901/902 supercluster field and our analysis of that reconstruction. We start by presenting the raw Wiener reconstruction in Fig. 6, splitting the field into two parts: the northern half of the field, comprising A901a/b (left column), and the southern part, comprising A902 and the SW group (right column). In the upper panels we present the 3D density map projected on the plane of the sky showing with high significance the four structures in the supercluster. This projection shows the same result as the 2D reconstruction in Heymans et al. (2008) and we overlay signal-to-noise contours obtained by comparing the mass map to noise realisations. Dividing by the noise in the 3D mass map, changes the z-shift bias in the new resulting map, rendering it essentially a 2D map for STAGES. For this reason, the S/N map is shown only in 2D as maximum along a l.o.s. direction. In the lower panels, we present the 3D map now projected along the declination. This reconstruction shows the redshift-density distribution convolved with the the radial p.s.f. of the map, shown in Fig. 1, corrected for z-shift bias, shown in Fig. 2. As this radial p.s.f is so broad, structures are radially spread out in the map and overlap in redshift. Nevertheless, in cases where structures do not lie in close projection, once the z-shift bias is corrected local maxima are unbiased estimators of the redshift of mass concentrations. Whilst unbiased however they have limited redshift accuracy, reflected by Fig. 3.
In order to better visualise the 3D matter distribution from the raw Wiener map, we make the general assumption that the local maxima are good estimators of the peaks in the density field. We then create a 3D map retaining only the density pixels that are a local maxima along each l.o.s.. We then smooth the visualisation with a Gaussian kernel (0.2 arcmin) 2 in r.m.s. size on the sky and 0.025 in redshift with the result shown in Fig. 7. This cleaning process produces a striking image of the 3D density field revealing structure at high redshift (marked SWc and b/3) behind the foreground supercluster (marked A901a, A901b, A902, SWa/b). It is this realisation of the noisy 3D density field that should be compared with other 3D mass reconstructions in the literature Taylor et al. 2004). A one-to-one comparison is still not easy, though, since other reconstructions show the estimated gravitational potential. We strongly caution the reader over the direct use of 3D lensing maps, however, as they are noisy estimators of the redshift distribution and contain no error information. In particular, owing to the broad smoothing kernels, structure along the same l.o.s. can merge in this visualisation to produce a density profile with a local maximum at an intermediate redshift. As our detailed analysis of the raw Wiener map will show however, after the radial p.s.f has been correctly accounted for, the higher redshift structures picked out in this visualisation are still seen with a high probability. Figure 6. Different projections of the Wiener mass map (top row : projection onto the sky; bottom row : radial projection of the area in the top row; z-shift bias corrected) The colour scales shows the maximum density contrast in the map along the projection axis. The contours in the top panels correspond the S/N levels 2.5σ, 3σ, 4σ, 5σ and 6σ as determined from noise realisations of the data. Contours in the bottom panels are density contrast levels included to highlight the structure in the radial map. Structures appear radially stretched due to the impact of the Wiener filter smoothing kernel.

3D comparison of mass and light
One of the main science drivers for 3D mass mapping is to investigate how well the total mass traces light or stellar mass in 3D. We use the product map technique outlined in Sect. 2.8 to compare total mass and stellar mass in redshift slices and identify interesting l.o.s. in our 3D map.
We use 17-band COMBO-17 stellar mass M * estimates for the full galaxy sample down to R < 24 mag. The M * -map is produced by summing the stellar masses Mi of galaxies within the same angular bin and lens plane (pixel). Bins containing no galaxies have a stellar mass of zero. This matrix of binned masses is divided by the average pixel mass, M = i Mi/Np where Np is the total number of unmasked pixels. All pixels of the same l.o.s. are flagged unmasked if there is at least one pixel along the l.o.s. containing one or more galaxies. The majority of masked pixels are found at the edges of the field. The effect of masking in the central part of the map is hence small. The normalised M * -map of the lens plane is in the unmasked regions, up to a constant, an estimator for the M * density contrast, M (θ)/M − 1. In a next step, the average M (θ)/M is subtracted from all pixels.
In Fig. 8    Intensity levels in the background show the signal-to-noise in the product map with black regions having the most negative values and red the most positive levels. Signal-to-noise mass density contours are shown in red with levels 2.5σ, 3σ, 4σ, 5σ and 6σ, depicting the maximum significance along a l.o.s. in the lensing mass map. Note the source of the signal may be at any redshift as the S/N maps are almost completely spread out in redshift. White contours highlight regions of high stellar mass M * = 10 10 , 2 × 10 10 , 5 × 10 10 , 10 11 M ⊙ , h −1 in the M⋆ map smoothed with a 1.2 arcmin Gaussian kernel. The map pixel volume increases towards greater redshift. Open diamonds are reference points marking prominent bright cluster galaxies in A901a/b, A902 and the SW group. Numbers on the x-and y-axis denote the right ascension and declination, respectively (J2000.0). mass. The 3D stellar mass distribution is shown with white contours. We find a clear correlation of mass and stellar mass in the z = 0.165 A901/902 supercluster in the lowest redshift panel (upper left) and also the CBI cluster at z = 0.46 (upper right).
The interpretation of this map is however non-trivial as the radial p.s.f. of the map causes overdensities in the mass map to contribute to the product map even when they are at a different redshift to the concentration in the stellar mass map (Sect. 2.8). We can therefore only use this product map to identify interesting l.o.s., which we then follow up with a detailed analysis that incorporates knowledge of the radial p.s.f. of the Wiener map. We identify interesting l.o.s. as follows. First, we only consider pixels in the total mass map that are detected with a S/N > 3 in projection. We then choose all pixels where we find a peak positive or negative correlation of significance greater than 1σ with the stellar mass. The majority of the lines of sight above the 3σ threshold exhibit a positive correlation in the product map. According to the noise peak statistics in Appendix A3, however, we expect at least some of the matches below 3.5σ to be confusions of random lensing map fluctuations with the 3D stellar mass map. We find that the cases of negative correlation in the lensing map with S/N > 3 are all behind foreground structures A901a, A901b, A902 and SW, which are likely to be produced by the leakage effect, outlined in Sect. 2.8. With the exception of A901a, these negative correlations have a product map signal-to-noise of < 3σ and are therefore statistically insignificant. They are all included in the analysis of the selected l.o.s. with positive correlations in the foreground, and therefore automatically included in the following analysis. Our ten chosen l.o.s. are listed in Table  1 with identifiers based on the known foreground z = 0.165 structures in this group, A901a, A901b, A902, SW, and also the higher redshift CBI cluster at z = 0.46. 1 Before continuing with our detailed line of sight modelling of the ten structures marked in the product map, we briefly discuss other strong features. In the two highest redshift slices, we find two stellar mass overdensities at M * = 10 11 M⊙ (with 1.2 arcmin smoothing). In the z ∈ [0.50, 0.70] bin, we see a strong positive correlation P1 of total mass and stellar (shown white) where the total mass is detected at 2.5σ. In the z-bin z ∈ [0.70, 1.00], we find significant anti-correlation P2 (∼ 6σ in the product map), which also is a stellar mass overdensity matched by a lensing mass underdensity. The M * overdensity region is found to be overlapping with gaps produced by three stars in this area, which may be interfering with the lensing map. Moreover, as shown in Fig. 4, at these high redshifts we would only expect to reliably detect haloes of mass M200 ∼ 2 × 10 14 M⊙, a higher mass range than what we would expect for these stellar mass overdensities. The same argument holds for the feature P3. These stellar mass overdensities therefore demonstrate the effect of noise in our mass maps where we are detecting mass with lensing only at the 2.5σ level. We also note that we find Table 1. Compiled catalogue of l.o.s. candidates of positive masslight correlations with projected l.o.s. significance higher than 3σ in the lensing map (shown as S/Nκ); S/N PM lists the corresponding signal-to-noise in the product map. Lines of sight with no product map match are not listed. The equatorial coordinates (J2000.0) are the peak positions in the product map, redshift constraints are estimates with 68% credibility from a thorough l.o.s. modelling in Sect. 3.2. A dagger † denotes redshift estimates from a two component model rather than a one-component mass model. −0.12 no significant ( 3.5σ) lensing mass overdensity matched by a stellar mass underdensity in the entire 3D map ("dark clump").

Line of sight modelling
Using the product maps in Fig. 8, we identify interesting l.o.s. in the map to perform a more detailed analysis of the 3D structure that takes fully into account the radial p.s.f. of the Wiener map. For each l.o.s., listed in Table 1, we model the matter distribution by assuming that it consists of one individual non-zero mass halo ∝ 1 + δ0 on the lens planes. The model profile in the convolved 3D Wiener map is given bŷ where δ0 > −1 is the density contrast of the non-zero pixel on the ith lens plane; all other pixels have δ (i) m = 0 which amounts to average density in the corresponding epoche. The radial matter densityδm in direction θ of the Wiener map with noise covariance N enables us to write down the posterior likelihood of the mass pixel lens plane index i and its contrast δ0: where the Heaviside function prior H(x) asserts a positive density 1 + δ0 0. Marginalising over δ0, we employ the MCMC technique to construct a probability distribution P (i|δm) for each lens plane i, as shown as black lines in Fig. 9. This probability distribution gives the true probability of a density peak positioned along each l.o.s. at a given redshift. This result correctly accounts for the broad radial p.s.f that has hampered the interpretation of the maps presented up to this point. The blue error bars denote a 68% confidence regions about the median redshift of the non-zero pixel. For comparison, the density contrast from the M * map for the same pencil beam is shown as grey filled bars; the beams have a radius of ∼ 40 arcsec. The bars are normalised such that the maximum density reaches the half height of the diagram.
Motivated by the results of the product map in Fig. 8, we also test an alternative two-halo mass model for most of the l.o.s. , which seem to require a more complex modelling when comparing the one-component p.d.f. to the M * distribution in the pencil beam. The posterior likelihood of the two-component model is similar to Eq. 38, but now contains two positions i1 and i2 and two densities. As additional prior, we demand that i1 < i2. The resulting marginalised probability distributions are shown as black solid and dashed lines, the two error bars denote the constraints on redshifts in the two-component model. The orange solid line shows the probability distribution for the one-component model. Note that in a two-component model fit with no signal in the data, the probability densities of both components are not simply flat: Due to the adopted prior z1 < z2, the null signal distribution has a saw-tooth shape as in Fig. 10, which is produced due to the correlation of the probability distributions of z1 and z2.
In Fig. 9 we see good agreement between the peaks in the density redshift probability distribution and the overdensities in stellar mass. We also see clear evidence for l.o.s. projections in the stellar mass distribution which have a significant probability in the 3D density map. In almost all cases, a two-component model is not clearly preferred over a one-component model with the slight exception of A901α that prefers a two-component model (d.o.f.=18) with ∆χ 2 = 1.90 (∼ 68% confidence) as compared to a onecomponent model (d.o.f.=20). The best-fit one-component model of SWc has χ 2 = 17.23 as opposed to χ 2 = 16.70 for the two-component model, which is statistically insignificant. For SWc, however, a model with all mass concentrated at z = 0.165, the distance of SWa/b, is excluded with ∼ 68% confidence (χ 2 = 18.69; d.o.f.=21). Combining this information with our knowledge of the galaxy distribution, we favour the two component model which reveals a new high redshift structure SWc at z ∼ 0.7, lying in projection with the foreground SW group at z = 0.165.

DISCUSSION AND CONCLUSIONS
This paper uses 3D lensing data from the STAGES project to study the spatial distribution of matter density and galaxies in the Abell 901/2 supercluster field. We have presented a non-parametric lensing-based 3D mass density reconstruction. This application of the 3D mapping algorithm to data is the first of its kind as earlier attempts Taylor et al. 2004) to do cosmography with 3D lensing data focused on a reconstruction of the gravitational potential rather than the mass density itself. In addition to directly mapping the 3D density field, the new algorithm presented in this paper accounts for redshift bias, a term previously neglected in the literature, and it optimally includes all lensed sources including those without reliable photometric redshift estimates. 2 We present both our raw Wiener reconstruction of the STAGES 3D density field in Fig. 7 and a cleaned visualisation of the 3D field in Fig. 6. Using the visualisation is a valid tool to describe the structure in the field, in particular because, in this instance, the z-shift bias has been correctly removed. Further it is this cleaned map which should be compared to previous 3D gravitational potential maps in the literature where the displayed isophotal potential levels have been optimised to highlight structures. However we urge caution in using any of these cleaned visualisations for scientific analysis as they are very noisy estimators of the density field and do not show the redshift error in the map or account for the broad radial smoothing of the field, shown in Fig. 1. However we have created a striking visualisation tool of the 3D dark matter Universe whose impact should not be underestimated even if its direct use for scientific analyses is somewhat limited. Note, as discussed in STH09, for deeper observations and higher number densities of sources, the radial smearing effect is reduced but it will always be present as a generic feature of this Wiener filter 3D lensing method, no matter how good a data set is analysed. A future improvement, albeit not explored for this paper, may be to regularise the Wiener filter, e.g. by a prior that forces structure to be confined to one lens plane. Although this resulted in a more realistic estimate for the 3D matter distribution, it would not change the considerable uncertainties in the visualisation.
From the raw 3D Wiener lensing reconstruction we infer that most of the structure in the STAGES field belongs to the Abell supercluster at z ∼ 0.165 since most l.o.s. density profiles in the map are locally peaked at around this redshift. There are a few exceptions however, of which the most prominent is in the SW group region. Here, we find a strong peak at around z ∼ 0.7 apparent in the raw Wiener map.
In order to interpret the Wiener map result we have compared the total mass to the stellar mass in the field using a product map technique outlined in Sect. 2.8. In the resulting product map in Fig. 8, significant matches or mismatches of density fluctuations between total mass and M * from the selected redshift range are identified and followed up with a more detailed analysis in Sect. 3.2. The conclusions we can draw from this analysis are as follows: (i) We find no detection of "dark clumps" in the 3D map, i.e. matches of lensing mass overdensities with stellar mass underdensities, which have been discussed in the literature (von der Linden et al. 2006;Erben et al. 2000;Massey et al. 2007).
(ii) As expected from Fig. 3, A901b is well confined in redshift by the 3D lensing data, and the estimated matter distribution in the pencil beam fits well to the M * distribution.
(iii) A901a, the second massive cluster in the field is less well constrained in redshift, possibly explained by the presence of an extra mass peak at around z ∼ 0.3 as seen in M * and also highlighted in the product map (Fig. 8, top right  panel).
(iv) The SW group pencil beam is equally well fit by a one or two-halo model in projection. However the single-halo model exhibits an additional peak at higher redshift and is slightly inconsistent with all mass concentrated at z = 0.165. We therefore conclude that along this l.o.s. there are indeed two structures; the well documented foreground SW group at z = 0.165 and a second possible mass concentration at z = 0.7 ± 0.3, which is confirmed by a corresponding stellar mass peak. SWc is also the most prominent high-z peak in the raw Wiener map. This finding was missed in the past 2D lensing analyses of Gray et al. (2002);Heymans et al. (2008) and the smoother 3D map of the gravitational potential in Taylor et al. (2004).
(v) We find a similar result for the the A901α X-ray group known to be infalling on A901a at a redshift of z ∼ 0.17. In this case the one-halo model again rules out a lowredshift structure at moderate significance, suggesting that the equally well fit two-halo model best represents the data with widely scattered matter behind z = 0.165 as indicated by the M * map.
(vi) A902 is equally well fit by a low-redshift one-halo model consistent with the known low-redshift of the cluster, or a two-halo projection model with one at z ∼ 0.2 and one at z ∼ 0.8.
(vii) CBI, clearly visible as M * peak at z = 0.46, is weakly detected in the mass map, which indicates a mass feature at a consistent redshift z = 0.68 +0.30 −0.45 . (viii) We also find evidence for mass concentrations behind the A901b cluster termed as A901b/2 and A901b/3. The latter is more prominent in the product map, concentrated at z ∼ 0.7.
Note that none of these conclusions could be drawn from the purely 2D analysis of Heymans et al. (2008). It is the non-parametric, purely lensing-based modelling of selected l.o.s., Fig. 9, that has to be used to investigate these significant features. In contrast to our 3D Wiener map and comparable visualisations (Hu & Keeton 2002; Tay-Table 2. Constraints on the NFW virial masses M 200 and the virial radii r 200 for A901b and SWa,b in a simple model fit to the 3D lensing data and additionally A901b/2-3 and SWc in an updated model with more components. The redshift of the foreground mass peaks is fixed with z = 0.165 and with z = 0.8, 0.75, 0.75 for A901b/2-3 and SWc, respectively. Within the fit an uncertainty of ∆z = 0.05 is factored in. Quoted are 68% credibility regions about the mean of the posterior likelihood. As fiducial cosmology here Ωm = 0.3 and Ω Λ = 0.7 were assumed.  Taylor et al. 2004;Massey et al. 2007;Vander-Plas et al. 2011), information on the redshift uncertainty of lensing mass peaks are directly visible in this Bayesian approach. From here we can see that even for a survey as good as STAGES radial distances of objects at higher redshifts are usually weakly constrained, most strikingly visible for the background objects in a two-component fit. This confirms the theoretical estimate in Fig. 3, which is based on SIS matter haloes. Nevertheless, the influence of background matter on the A901b, A902 and especially SW group lensing signal is visible in Fig. 9; compare the foreground p.d.f. in the two-component model to the p.d.f. of the one-component model. This influence is hidden in a 2D lensing analysis. All in all the modelling and its comparison to the stellar mass distribution confirms that the cleaned visualisation of the 3D lensing map correctly identifies previously known structure (e.g. Heymans et al. 2008) at low redshift, but also seems to identify hitherto unknown structure at larger redshift, most notably A901b/3 and SWc. These overlapping mass concentrations are a possible source of error when it comes to mass determinations by lensing, a method that naturally measures the integrated mass within cylinders along the l.o.s. To gauge this effect, we fit Navarrow-Frenk-White (NFW) matter density halo profiles (Navarro et al. 1997;Wright & Brainerd 2000) with fixed redshifts and positions on the sky to the 3D lensing data. A method similar to the one in Heymans et al. (2008) is employed and hence not detailed here; the redshift p.d.f. of our shear catalogue subsamples are utilised, and constrains from all subsamples are combined. In one fit assuming a simpler mass model, we include SWa/b, A901b and A901α (all at z = 0.165) only, and in a second fit, invoking a more complex mass model based on the previous conclusions, we also include SWc (z = 0.75) and A901b/2-3 (z = 0.8, 0.75) as possible carrier of mass. The fit results can be found in Table  2. Comparing the model fits with and without 3D structure in projection shows a slight drop in the mass estimates for A901b and SWa/b but this decrease is not statistically distinguishable from the simple mass model. This demonstrates that previous 2D mass estimates in the literature may have overestimated mass measures but not significantly within the statistical errors. This result is expected based on Fig. 9 where single components are as equally well fit as a multiple component model along these lines of sight.
There is also the possibility that the observed features in the background are actually belonging to the foreground supercluster and are moved to higher redshift in the Wiener map owing to the noise in the lensing signal; their mass is relatively small (Table 2) in the light of Fig. 3. In this case, the match with stellar mass at higher redshift would be a coincidence.
We have presented a method that can be used to account for the radial smearing of the 3D mass map and produce full redshift probability distributions along individual l.o.s.. In future work an alternative, potentially more fruitful method could be to first identify significant lensing detections in a 2D map and then produce smaller scale 3D reconstructions in these interesting regions. The benefit of this method would be that each 3D reconstruction could be optimised for each source by tuning the α parameter in the Wiener reconstruction (see Sect. 2.2) depending of the mass and redshift of the source. In addition, information from measurements of higher order weak lensing flexion distortions could be included in the analysis to boost the signal to noise in the reconstruction (Velander et al. 2011). Moreover, for deeper surveys, which have reliable photometric redshifts down to z ∼ 2 or more, we also anticipate an moderately improved mass density recovery at higher redshift based on simulations ( Fig. 13 in STH09).
Throughout this paper, we have neglected errors caused by the intrinsic alignment of galaxies, which for this analysis would be dominated by the 'GI' term, the anti-correlation of sheared background sources with the intrinsic shape of foreground galaxies. Based on analysis of simulations and data, the systematic error contributed by this astrophysical effect is expected to be low on average (Joachimi & Bridle 2010) and probably negligible for the results presented in this paper. That said, the intrinsic alignment of galaxies in cluster environments has yet to be studied in data but has been explored using analytical prescriptions (Schneider & Bridle 2010). If cluster galaxies strongly align as they infall this could result in a strong GI term that could effect our results by lowering the lensing signal measured behind the foreground cluster. The opposite however could also occur whereby galaxy merging events in the dense cluster environment randomise the intrinsic ellipticities resulting in a very weak GI term. As discussed in STH09, the contribution of GI and II correlations to the shear field can in principle be incorporated into the minimum variance estimator for the 3D map. At this stage, however, we expect this effect to be small. This uncertainty is something that will be investigated further using the STAGES data in future work.
The driving motivation behind the development of our 3D reconstruction technique was to enable an unbiased 3D comparison of mass and light. Dark haloes for example would only be detected in this manner. However with the detailed analysis of noise and the radial p.s.f in the 3D lensing reconstructions presented in this paper, we have demonstrated how inherently noisy the process is. Even with one of the best HST data sets and a very massive foreground system we have found it difficult to reconstruct the 3D mass distribution at any high degree of spatial and redshift resolution using weak lensing alone, especially at higher redshift. Given the limitations of the method to resolve only the most massive structures in 3D the future direction for the application of this method should be to reconstruct larger scale structures in the 3D density field using more heavily spatially smoothed data. With upcoming wide area surveys such as CFHTLenS, DES and VST-KIDS, which are several orders of magnitude wider than the STAGES survey presented in this paper, we can expect higher quality 3D resolution reconstructions on degree scales because on these scales the significance of modes in a 3D mass density reconstruction are increased (STH09).
p.s.f. model is then applied to the data. CH also chooses to maximise the signal-to-noise in the lensing catalogue fitting the isotropic correction term P γ as a function of galaxy size (see Heymans et al. (2008) for more details).
In the analysis of the dense stellar fields, TS has shown the HST ACS p.s.f. can vary rapidly with time. To account for this, TS employs a fully time-dependent correction for the spatially and temporally varying p.s.f. based on individual exposures, using a p.s.f. template library obtained from stellar fields (Schrabback et al. 2007). In addition, TS applies corrections for charge transfer inefficiency and signal-to-noise dependent shear calibration bias as detailed in Schrabback et al. (2010). In this analysis of the STAGES data TS employs conservative galaxy selection criteria including a high cut on the shape measurement signal-to-noise ratio (Erben et al. 2001) The resulting catalogues differ in the number density of objects and the noise on the shear measurement quantified through the root-mean square variation of the measured galaxy ellipticities. CH, hereafter Sample I, contains 65 galaxies per square arcmin with σǫ = 0.26. TS, hereafter Sample II contains 34 galaxies per square arcmin with σǫ = 0.29.

A2 Photometric redshifts: COMBO-17 sample
The STAGES HST data is complemented by the 17-band optical COMBO-17 survey. Using a combination of broad and narrow-band filters, COMBO-17 uses template-fitting to estimate photometric redshifts (Wolf et al. 2004). For the A901/902 field the redshift error in units of σz/(1 + z) is estimated at ≈ 0.02 to R < 23 with an increase to 0.035 in the interval R = [23, 24] (see Gray et al. 2009). Outlier rates also increase with magnitude: we estimate roughly 1% outliers with |∆z|/(1 + z) > 0.10 at R < 22 and 7% outliers with |∆z|/(1+z) > 0.15 at R = [23,24]. Above a magnitude R = 24 we do not use COMBO-17 redshift estimates due to rapidly increasing errors and outlier rates.
Owing to the width of the lensing kernel in redshift space, the accuracy of the COMBO-17 redshifts are more than sufficient. Redshift estimates, however, only exist for a small fraction of the galaxies (13 galaxies per square arcmin for Sample I and 10 galaxies per square arcmin for Sample II). For the remaining galaxies we assign a magnitude dependent redshift distribution given by Schrabback et al. (2007) (see Fig. A1). All this information is combined for the final mass map.
The stellar mass estimates used in this analysis derive from low resolution 17-band spectra fits to parameterised star formation history models as described in Borch et al. (2006).
A3 Quality control; comparison of E and B mode signals to pure noise maps A reconstruction of signal and noise, based on the shear Sample I and shear Sample II, is carried out following STH09 with the resulting map is shown in Fig. 6 for Sample I and Fig. A2 for Sample II. We find that both samples detect the same significant structures at the positions of the known galaxy clusters A901a, A901b and A902, in addition to a detection of a group of galaxies called the South West group (SW) as seen before in earlier 2D reconstructions of the projected matter density (Heymans et al. 2008;Gray et al. 2002). As expected from the galaxy number counts in the two samples, the signal-to-noise in the maps is higher for Sample I than Sample II. All structures that differ between the maps are at a level consistent with noise.
For an analysis of the systematics, we count the number of peaks and troughs in the E-and B-mode S/N 3D maps given signal-to-noise range and compare this to the expected numbers in pure noise maps. The noise statistics are obtained from the randomised data. In order to have an easily computable and objective statistic for "peaks" and "troughs", we count the number of grid pixels in the Wiener map under the condition that the pixels are local maxima within a cuboid volume centred on the pixel. The size of the cuboid is ±2 pixels in transverse and ±2 lens-planes in radial direction. For S/N-thresholds 3σ this definition fits reasonably well the number of peaks one would select by eye. Troughs are considered negative peaks (local minima with negative signal-to-noise). Fig. A3 shows the results for Sample I and Sample II. Clearly, for thresholds 3.5σ, E-mode peaks (troughs) are above the average noise levels and above 1σ-noise variance for 4.0σ, whereas B-mode counts are consistent (Sample I) with the noise expectation. From the noise realisations, we estimate the covariance of the noise counts which is used for a χ 2 -test, peak and trough counts combined, to check consistency with noise. For Sample I we find χ 2 /d.o.f. = 177.4 (E-mode map) and χ 2 /d.o.f. = 0.40 (B-mode map), showing that the E-mode map is inconsistent with noise, while the B-mode map is consistent with noise. Corresponding values for Sample II are χ 2 /d.o.f. = 220.0 (E-mode) and χ 2 /d.o.f. = 0.29 (B-mode) leading us to a similar conclusion.
We therefore find the B-modes for both Samples to be low and consistent with noise.This is potentially a surprise for Sample I, bearing in mind the rapid variations in p.s.f. ellipticity seen in Schrabback et al. (2007) but unaccounted for in the Sample I analysis. The success of the semi-time dependent p.s.f. subtraction method of Sample I may be explained by the fact that HST was significantly under-focus during the STAGES observations with the result of a relatively stable p.s.f., once averaged over intra-orbit variations.
Owing to the larger Sample I, there are peaks up to 7σ, unseen for Sample II. Note that in the high-S/N regime the number of peaks or troughs becomes rare and statistically Poissonian (one or no peak) which may not be properly reflected by the symmetric error bars in the figure.

A4 Redshift Scaling
Recent lensing analysis of ground-based data using the KSB method found errors in the expected redshift scaling of the measured cosmic shear signal (Kilbinger et al. 2009). The reason for this bias is unknown but is thought to lie in a weak signal-to-noise and galaxy-to-psf size dependent calibration bias inherent to shear measurement methods Bridle et al. 2010;Semboloni et al. 2009).
As a further quality check we test the redshift scaling of the KSB-measured space-based lensing signal for Sample I and II. We select a set of peaks that are known to be at redshift z ≈ 0.165 (A901a/b/α, A902/1 and SW/a/b). For each peak the mean tangential shear within an aperture, of radius between one and two arcmin, is measured from the data: where θc is the l.o.s. direction of the aperture centre, θi the angular position of a source in the common complex notation and ǫi the complex ellipticity of the source. Only sources within the aperture and within some redshift range are considered. The statistical uncertainty of the estimator for the mean shear is the variance of the mean as obtained from the sources included in the sum. The mean redshift,z, of included sources determines the redshift associated with γt. The statistics from all lensing peaks is combined and weighed by their statistical uncertainties. Thus we are stacking the cluster signals to increase the signal-to-noise of this test. Fig. A4 shows aperture statistics for Sample I and Sample II by plotting the real parts ofγt (E-mode) and the imaginary parts (B-mode; radial shear) separately. If the lensing signal within the apertures stems in essence only from structure located at the cluster redshift z ∼ 0.165, we expect a redshift scaling dictated by the lensing efficiency (with fixed lens distance): where f k (w) is the angular diameter distance as function of radial comoving distance w; w h , here constant, andws are the halo and average source distance, respectively. The last step is only approximately correct. This approximation is used for the fit in Fig. A4. This behaviour is independent of the intrinsic density Figure A4. Test of redshift-scaling of the lensing signal in the data. Mean tangential shear as function of source redshift within an aperture of radius 2 arcmin (A901a and A901b) or 1.4 arcmin (A901α, A902/1 and SWa/b) centred on the peak centres; the signal for all four apertures is averaged. The sources are binned in redshift for the samples with photometric galaxy redshifts (z 1.2). The three highest z-bins belong to the magnitude samples for which redshifts correspond to the mean as determined from the redshift distribution of sources. With most matter inside the aperture located at the same redshift, the mean tangential shear (solid points) should scale according to the black solid line in the plot. The shaded region around the line is the 1σ-confidence region of the amplitude. Open crosses denote the average radial shear inside the aperture, which should be vanishing. The shaded red area marks the 1σ-area for the mean B-mode over all redshifts. profiles of the clusters, albeit its amplitude depends on the cluster details. As we are solely interested in verifying the correct redshift scaling behaviour, we fit Eq. A2 with the amplitude as unknown variable to the E-mode data points which for χ 2 /d.o.f. = 0.46(0.50) is an excellent fit to Sample I (Sample II). B-modes should be consistent with zero for all redshifts. This we test by fitting a straight horizontal line, with unknown amplitude, to the average radial shear. We find that the fit is both consistent with zero on a 1σ-level and well fitted by a model with redshift-independent signal, yielding χ 2 /d.o.f. = 0.97(0.87).

A5 Conclusion
In conclusion, we find both shear catalogues to have a low level of systematics in the B-mode channel consistent with noise and a redshift-scaling of the shear signal that is consistent with ΛCDM cosmology. Moreover, we find no significant differences between the E-mode channel maps of Sample I and Sample II. As Sample I yields an overall more significant map owing to the higher number density of sources, we present results in the main section of the paper only from the data of Sample I (CH).