© ESO 2018

1 Introduction

Matter components in galaxies and galaxy systems contribute to the common gravitational potential with roughly concentric distributions. This matter distribution hinder the analysis of dark- and baryonic-mass distributions independently. Interacting galaxy clusters allow us to study the different mass components separately, given that, in the standard “cold dark matter” (CDM) model, the intracluster gas and the dark matter component can be spatially segregated in the merging process (Furlanetto & Loeb 2002). The galaxies and the dark matter of the merging systems behave as collisionless particles, while the gas component is fluid-like and experiences ram pressure. Therefore, the analysis of merging clusters provides independent information regarding the baryonic component, mostly constituted by the intracluster plasma, and the gravitational potential traced by the dark-matter distribution.

Merging galaxy clusters are the largest and most energetic events in the Universe since the Big Bang (Sarazin 2002). Also, they provide evidence that supports the existence of dark matter (Dawson 2013; Dahle et al. 2013; Harvey et al. 2015) since, when clusters collide, galaxies behave as collisionless particles while the X-ray emission plasma is slowed by ram pressure. This scenario allows us to test the existence of dark matter. Without dark matter, the gravitational potential would be traced by the dominant baryonic component, the X-ray plasma. On the other hand, if the principal component in clusters was collisionless dark matter then the gravitational potential would be traced by this component which is expected to follow the galaxy distribution. Thus, deriving the gravitational potential distribution allows us to distinguish between these possibilities.

One of the hottest and most X-ray-luminous clusters observed is 1E0657-558, the Bullet cluster. This cluster, at z = 0.296, was discovered by Tucker et al. (1995) using Chandra X-ray observations. The most probable scenario would be that the Bullet cluster is the result of the merging process of two systems, taking place in the sky plane (Barrena et al. 2002; Markevitch 1998). X-ray observations revealed the presence of a fairly cold (T ~ 6 keV) bullet-like structure, just leaving the core of the dominant system (T ~ 14 keV) with a velocity of 4500 km s⁻¹, producing a prominent bow shock (Markevitch 1998). Taking into account that the offset between both systems is 0.66 Mpc, the closest approach would have occurred 0.1–0.2 Gyr ago. Weak lensing analysis of this system showed that the gravitational potential does not follow the plasma distribution, but rather is consistent with the galaxy distribution (Clowe et al. 2006). The discrepancy between the total mass center and the baryonic center provides evidence for the existence of dark matter with an 8σ confidence level.

Several works have analyzed Bullet-like clusters obtaining direct empirical proof of the existence of dark matter (Jee et al. 2014; Bradač et al. 2008; Dawson et al. 2012; Andrade-Santos et al. 2015; van Weeren et al. 2017). Even though the standard cosmological model has proved to be successful at large cosmological scales, there have been differences between the predictions of this model and the observations at smaller scales. These disagreements could be overcome if dark-matter particles were self-interacting (Spergel & Steinhardt 2000; Rocha et al. 2013; Zavala et al. 2013). The analysis of merging systems could provide an upper limit for the ratio between the self-interacting cross-section and the mass of the dark-matter particle, which could be obtained by taking into account the separations of the matter components after the interaction between the galaxy systems. Harvey et al. (2015) obtained ${\sigma }_{\text{DM}}/m=-{0.25}_{-0.43}^{+0.42}$ cm² g⁻¹ with a confidence level of 68%, or σ_DM ∕m < 0.47 cm² g⁻¹ with a confidence level of 95% based on the analysis of 30 merging clusters. This result rules out parts of model space of hidden sector dark-matter models that predicts σ_DM ~ 0.6 cm² g⁻¹ on cluster scales through a long-range force.

With the aim of testing the standard cosmological model ΛCDM, different studies used numerical simulations to analyze the probability of observing a system with similar properties to those observed in the Bullet Cluster. There is extensive debate regarding this issue. Several works conclude that the velocity observed for the less massive system would not be physically possible in a ΛCDM scenario (Markevitch 1998; Springel & Farrar 2007). Even taking into account that the matter velocity (~ 2600 km s⁻¹, Milosavljević et al. 2007) would be lower than the plasma velocity, ΛCDM models would not be able to produce Bullet-like systems (Lee & Komatsu 2010; Thompson & Nagamine 2012). Nevertheless, results based on larger cosmological simulations demonstrate that these systems are not extreme, thus, there would be no tension between the observationsand the standard cosmological model. Interestingly, Thompson & Nagamine (2012) conclude that a volume of ${(4.5h^{-1}\text{Gpc})}^3$ would be required in order to observe Bullet-like clusters.

Taking into account the importance of the analysis of merging systems, it is necessary to study uniformly selected samples to derive robust constraints. de los Rios et al. (2016) have recently presented a method to identify galaxy merging systems based on redshift and photometric galaxy catalogs. This methodology provides a highly reliable sample of candidates for merging systems with low contamination and some estimated properties. This paper is the first in a series of works that analyzes merging systems cataloged by de los Rios et al. (2016). The analysis includes photometric, lensing, X-ray, radio, Sunyaev–Zeldovich, and dynamical studies. Our aim is to characterize the cataloged merging candidate clusters.

In this paper, the first of a series devoted to characterizing galaxy clusters in merger processes, we present a dynamical and a lensing analysis of the galaxy clusters A1204 and A2029/A2033. We derive the lensing masses for each identified interacting structure, then we characterize the dynamical state of these clusters according to a two-body model. We also perform a qualitative description of the gas and the mass distribution in the field of the systems based on a lensing and an X-ray analysis. These systems were randomly selected from the merging system catalog, taking into account the publically available data for the lensing and X-ray analysis. Details regarding the redshift catalogs for the dynamical analysis, the observations used for the lensing, and the X-ray studies are presented in Sect. 2. In Sect. 4 we characterize the dynamic state of these systems. In Sect. 3 we describe the lensing methodology. We present and discuss our results in Sect. 5 and finally our conclusions are presented in Sect. 6. When necessary, we adopt a standard cosmological model: H₀ = 70 km s⁻¹ Mpc, Ω_m = 0.3 and Ω_Λ = 0.7.

2 Data acquisition

2.1 Photometric observations

Details of the studied objects and images for the weak lensing analysis are presented in Table 1. There we show the location of the primary components of the merging systems and their corresponding images used for the analysis.

For the lensing analysis, in the case of A1204, we use optical images retrieved from the SMOKA (Subaru Mitaka Okayama Kiso Archive) database. Observations were obtained using the Suprime-Cam (Miyazaki et al. 2002)with the wide-field prime focus of the 8.2m Subaru telescope located at the Mauna Kea Observatory, under the observation program COSMOS (Cosmic Evolution Survey). The Suprime-Cam is an 80-megapixel (10 240 × 8192) mosaic CCD camera, covering 34^′× 27^′ with a resolution of 0.202^″∕pixel. Retrieved images are in the i+ and V Subaru filters. Reduction and astrometric calibration is performed using SDFRED2 in a standard mode, by combining 90 and 70 frames for filters i+ and V , respectively.

In the case of A2029/A2033, observations for the lensing analysis were retrieved from the Canadian Astronomy Data Centre (CADC), taking into account the optical requirements for the lensing measurements: long exposure times, red optical filters, and a large celestial coverage. Therefore we use reduced and calibrated data observed with the MEGACAM camera mounted on the 3.6m CFHT telescope in r^′ and i^′ bands. MEGACAM consists of 36 CCDs with 2048 × 4612 pixels with aresolution of 0.187^″, covering a full 1 × 1 square degree field of view. Observations were carried out through queued service observing (QSO). The data were pre-processed (removal of the instrumental signature) and calibrated (photometry and astrometry) using Elixir version 2.0. In order to perform the lensing analysis, we use a parallelized pipeline, to carry out the source detection, photometry, classification and the shape measurements in the 36 frames in parallel.

2.2 X-ray observations

In this workwe intend to describe the gas distribution in the field of the analyzed systems from the X-ray brightness contours. Hence, our present work does not focus on the detailed spatial distribution and spectral characteristics of the hot intra-cluster gas. Taking this into account, for our X-ray analysis we choose to make use of the ROSAT data as a compromise between sensitivity, spatial resolution, and field of view.

ArchivalROSAT¹ imaging data for our sample clusters were downloaded from the High Energy Astrophysics Science Archive Research Center (HEASARC) on-line repository. Data acquisition details are showed in Table 1. The ROSAT X-ray data are in the form of event files, which we use to generate images. The ROSAT Position Sensitive Proportional Counter (PSPC) observations in the soft and hard bands were used to construct a broad band image. These data sets were prepared for analysis by first eliminating sources of contamination and selecting good time intervals. We calculate the background, assuming a Poisson statistics with the tool ximage from the HEASARC software, to subtract it from the surface brightness. For the broad band image of A1204 we obtain a root mean square (rms) noise of 0.32 (cts/im pix) and from the data of A2029 the rms noise obtained is of 0.72 (cts/im pix). The X-ray images in the 0.1–2.4 keV band are used to measure the X-ray surface brightness profiles of the analyzed clusters.

2.3 Redshift catalogs

Here we describe the redshift catalogs used for dynamical analysis. For the A1204 cluster we use the spectroscopic data from the Hectospec Cluster Survey (HeCS). This catalog was constructed by Rines et al. (2013) using the Hectospec instrument (Fabricant et al. 2005) mounted on the MMT 6.5m telescope. The HeCS clusters were selected from X-ray catalogs based on the ROSAT All Sky Survey (RASS, Voges et al. 1999), restricting to systems with 0.1 ≤z ≤ 0.3. The final HeCS sample consists of 53 clusters from the Bright Cluster Survey (Ebeling et al. 1998) and ROSAT-ESO FLux-LimitEd X-ray cluster survey (Böhringer et al. 2004) catalogs within the SDSS DR6 photometric footprint and with f_X ≥ 5 × 10⁻¹² erg s⁻¹ cm⁻² and four clusters with f_X ≥ 3 × 10⁻¹² erg s⁻¹ cm⁻².

In the case of the A2029/A2033 cluster, we use the spectroscopy data from the SDSS DR7 catalog (York et al. 2000; Abazajian et al. 2009) restricting to galaxies with apparent r magnitude between 14.5 and 17.77. In order to identify the galaxy clusters we perform a friends-of-friends (FoF) algorithm as described in Merchán & Zandivarez (2002, 2005). We use a transversal linking length corresponding to an overdensity of δρ∕ρ = 80, a line-of-sight linking length of V ₀ = 200 km s⁻¹, and a fiducial distance of D_f = 10 h⁻¹ Mpc. We also compute a scaling factor R using a galaxy luminosity function fitted with a Schechter function with parameters (α = −1.05 ± 0.01, M_* − 5logh = −20.44 ± 0.01) given by Blanton et al. (2003). It is worth noting that this algorithm joined A2029 with A2033, indicating that they may be physically connected.

3 Weak-lensing analysis

In order to obtain the projected mass density distribution and the total mass of each identified structure, we perform a similar analysis as the one presented in Gonzalez et al. (2015). In the following subsections, we give the details of the criteria used for the source detection and classification, and we describe how we obtain the shape parameters of galaxies classified as background galaxies. Finally, we obtain the two-dimensional (2D) projected density distribution and estimate the total mass of each subcluster, considering the measured ellipticities of these galaxies.

3.1 Source identification and classification

Source detection and photometry is carried out using SExtractor (Bertin & Arnouts 1996) in a two-pass mode as described in Gonzalez et al. (2015). We detect and compute the photometry for the brightest objects by using a threshold of five times the standard deviation of the background. This first run is made in order to compute the seeing and the saturation level of each frame. The seeing is obtained by averaging the full width at half maximum (FWHM) of the point-like objects, according to their position in the FWHM/magnitude diagram. The saturation level is estimated as 0.8 times the maximum value of the FLUX_MAX parameter, given by SExtractor and defined as the peak flux above background. Then, we run SExtractor considering the obtained seeing and saturation level in the configuration file and with a detection threshold of 1.5.

We classify the detected sources in galaxies and stars, and discard the false detections taking into account their position in the magnitude/central flux diagram, the FWHM respect to the seeing, and the stellarity index given by SExtractor. In Fig. 1 we show the diagrams used for the source classification. Visual inspection of the discarded sources revealed that these objects are mainly hot pixels, sources at the edges of the CCD, spikes of saturated stars, cosmic rays and overly dim sources.

In orderto perform the lensing analysis we need to identify the background galaxies which would be indeed distorted by the lensing effect. Since we do not have redshift estimations for all the identified sources, we use a photometric criterion to classify background galaxies. Hence, we consider a galaxy as a background galaxy if its apparent magnitude satisfies m_p < mag < m_MAX + 0.5, where mag is the measured magnitude in r^′ and i+ bands for CFHT and Subaru frames, respectively, m_p is defined as the magnitude for which the probability that the galaxy is behind the cluster is larger than 0.7, and m_MAX corresponds to the peak of the mag distribution of all the identified galaxies. The upper cut ensures that we are not taking into account overly faint galaxies with higher uncertainties in the shape measurements. We obtain m_MAX = 25.7 and m_MAX = 22.7 for the Subaru i+ and the CFHT r^′ images, respectively. A discrimination by color is taken into account to discard blue galaxies which are mainly foreground galaxies that dilute the lensing signal. Therefore, we consider for the analysis only galaxies with r^′ − i^′ > −0.5 and V − i+ > −1.0 identified in the CFHT and Subaru frames, respectively.

The lensing efficiency depends on the geometrical factor defined as β := D_LS∕D_S, where D_LS and D_S are the angular diameter distances from the lens to the source and from the observer to the source, respectively. To estimate m_p and ⟨β⟩ (where ⟨… ⟩ express the average over the considered background galaxies for the lensing estimator) we use catalogs of photometric redshifts. In the case of the lensing analysis performed with CFHT frames we use the Coupon et al. (2009) photometric catalog, based on the public release Deep Field 1 of the Canada-France-Hawaii Telescope Legacy Survey (CFHTLS Deep1), which has an 80% completeness limit of $m_{r^{\prime }}\approx 26.$ and covers a sky region of roughly 1 deg². On the other hand, for the analysis with Subaru frames we use the catalog of photometric redshifts given by Laigle et al. (2016) which contains precise photometric redshifts over the 2 deg² COSMOS field and with a limiting magnitude m_i+ = 26.2. We compute the fraction of galaxies with z > z_cluster in magnitude bins of 0.25 mag and then we chose m_P as the lowest magnitude for which the fraction of galaxies was greater than 0.7, obtaining m_P = 18.2 and m_P = 18.5 for i+ and r^′ images, respectively. Then we apply the photometric selection criteria (m_P < m_mag < m_MAX + 0.5) to the catalog and we compute β for the whole distribution of galaxies. To take into account the contamination by foreground galaxies given our selection criteria, we set β(z_phot < z_cluster) = 0 which outbalances the dilution of the shear signal by these unlensed galaxies.

To estimate the error in ⟨β⟩ regarding the cosmic variance, we divide CFHTLS Deep1 and COSMOS fields into 25 and 64 non-overlapping areas of ~144 and ~160 arcmin², respectively. Then we compute ⟨β⟩ for each area considering the redshift of the analyzed clusters according to the used images. The uncertainties in ⟨β⟩ due to cosmic variance are estimated according to the scatter among the values for each area, obtaining ~ 0.012 for CFHTLS DEEP1 field and ~0.010 for COSMOS field. These uncertainties were taken into account in the error estimation of the fitted parameters, and propagated to the resulting system masses.

In orderto take into account the contamination of foreground galaxies in the catalog, we compute for each galaxy the probability that it is behind the cluster. This probability is computed using the described photometric catalogs considering the fraction of galaxies with z > z_cluster for each bin in magnitude, mag, and color (V − i+ and r^′− i^′). Hence, given the magnitude and the color of each galaxy, we assign to it a weight, w, as the fraction of galaxies with z > z_cluster in that bin. This weight is applied to compute the mass and the 2D density distribution.

Fig. 1 Classification of the detected sources in the Subaru frame (upper panel) and in one of the frames of the CFHT mosaic (lower panel). Here stars are represented with blue points, galaxies by green points and false detections by black crosses. Left panels show the central surface brightness of the objects (μMAX ) as a function of the magnitude. Stars are confined in the region marked by the red lines. Right panels show the source distribution in the diagram FWHM as a function of the magnitude. As it can be seen, the FWHM of the sources classified as stars remains roughly invariant for the total range in magnitude.

3.2 Shape measurements

Gravitational lensing effects are characterized by an isotropic stretching called convergence, κ, and an anisotropic distortion given by the complex-value lensing shear, γ = γ₁ + iγ₂. Using the second derivative of the projected gravitational potential to express the shear and convergence, one can show that for a lens with a circular-symmetric projected mass distribution, the tangential component of the shear, γ_T , is related to the convergence through Bartelmann (1995): $${\gamma }_{\text{T}}(r)=\overline{\kappa }(<r)-\langle \kappa \rangle (r),$$(1)

where $\overline{\kappa }(<r)$ and ⟨κ⟩(r) are the convergence averaged over the disc and circle of radius r, respectively. On the other hand, the cross component of the shear, γ_× defined as the component tilted at π∕4 relative to the tangential component, should be exactly zero.

If lensing is weak (κ ≪ 1), the image of a circular source appears elliptical with a and b as major and minor semi axes, respectively. The induced ellipticity could be directly related to the shear, γ ≈ e. Here we define the ellipticity as a complex number, e = e₁ + ie₂, with magnitude |e| = (a − b)∕(a + b) and orientation angle determined by the direction of the major elliptical axis. If the source has an intrinsic ellipticity e_s the observed ellipticity would be e ~ e_s + γ (Bartelmann & Schneider 2001). If we consider many sources with intrinsic ellipticities randomly orientated so that ⟨e_s ⟩ = 0, the ensemble average ellipticity after lensing gives an unbiased estimate of the shear: ⟨e⟩ ≈ γ. Therefore, we can estimate the shear by averaging the shape of background sources that would be affected by the lensing effect.

In order to perform the lensing analysis and, in turn, to estimate the projected masses of the galaxy systems, it is necessary to obtain the shape parameters of background galaxies. It is important to take into account the roundness effects due to the atmosphere presence, as well as the distortion caused by the telescope effects. All of these are considered by the point spread function (PSF), which is convolved with the true galaxy intensity light distribution.

To obtain the shape parameters we use im2shape (Bridle et al. 2002), which models the galaxies as a sum of Gaussians convolved with a PSF, which is, in turn, a sum of Gaussians. For simplicity both the PSF and the object are modeled with a single elliptical Gaussian profile. To estimate the PSF at the position of each galaxy we average the shape parameters of the five closest sources classified as stars, since they are considered point-like sources (this procedure is described in Gonzalez et al. 2015). To check our PSF correction we obtain the shape parameters of the stars with and without considering the PSF. In Fig. 2 we show the semi-major axis distribution of the stars before and after considering the computed PSF for obtaining the shape parameters. As can be seen from the figure, when we consider the estimated PSF at the position of each star, the semi-major axes are more randomly distributed and their lengths are substantially smaller, according to the distribution that would be presented by point-like sources.

To compute the lensing masses and the 2D density distributions we only kept the galaxies classified as background galaxies and with FWHM > 6 pix, to ensure galaxies with a good pixel sampling and, therefore, with better-shaped parameter estimations. Also, we discard galaxies with σ_e > 0.2, where σ_e is defined as the error in the measured ellipticity and is computed as ${\sigma }_e=\sqrt{{\sigma }_{e_1}^2+{\sigma }_{e_2}^2}$, where ${\sigma }_{e_1}$ and ${\sigma }_{e_2}$ are the errors in the ellipticity components, e₁ and e₂, respectively,provided by im2shape. With these criteria we obtain a density of background galaxies of ~17.4 galaxies arcmin⁻² for Subaru image and ~ 3.5 galaxies arcmin⁻² for CFHT mosaic image.

Fig. 2 PSF treatment applied to stars in the Subaru frame (upper panel) and in one of the frames of the CFHT mosaic (lower panel). Semi-major axis distribution in the CCD, of the sources classified as stars, before (left) and after (right) the PSF deconvolution. For a better visualization of the plot, semi-major axes have a different scale than the CCD axes. Semi-major axis scale is given by the red thicker segment in the upper left-hand corner and corresponds to 3 pixels. We highlight that the semi-major axis is more randomly distributed and is significantly smaller after taking into account the PSF.

3.3 Individual masses

Once we obtain the shape parameters of the background galaxies we can compute the shear profile which couldbe fitted by assuming a density distribution for each structure identified in the analyzed systems. Spherical symmetry implies that the average in annular bins of the tangential component ellipticity, e_T , of the lensed galaxies traces the reduced shear. On the other hand, the average in annular bins of the component tilted at π∕4 relative to the tangential component, e_×, should be exactly zero for the case of perfect symmetry (e.g., Bartelmann & Schneider 2001, Sect. 4). Since galaxies have an intrinsic ellipticity, the error in the shear estimator ⟨γ(r)⟩, obtained by averaging the tangential ellipticity component of the N galaxies at adistance r ± δr from the center considered, ⟨e_T(r)⟩, would be (Schneider et al. 2000): $${\sigma }_{\gamma }\approx \frac{{\sigma }_{ϵ}}{\sqrt{N}},$$(2)

where σ_ϵ is the dispersion of the intrinsic ellipticity distribution (σ_ϵ ≈ 0.3, Leauthaud et al. 2007, which is fixed at this value for the rest of the analysis). The brightest galaxy member identified from each structure, according to the membership presented in the catalog of de los Rios et al. (2016), is considered as the center to build the profiles. Shear profiles are computed using non-overlapping logarithmic annuli, in order to have a similar signal-to-noise ratio in each annuli. To obtain the shear estimator, ⟨γ(r)⟩, we compute a weighted average of tangential ellipticity components of the N galaxies in the bin center at a radius r. The weight considered is the probability computed for each galaxy to be behind the cluster, as was described in Sect. 3.1. We test different annuli sizes but the final mass results do not show a strong dependence on this parameter. We fix the size for the one for which we obtain the lowest χ² value for the fitted parameters. Profiles are fitted from the inner part where the signal becomes significantly positive, to reduce the impact of miss-centering, up to ~1.5 Mpc to avoid the region where the signal relating the companion structure becomes significant.

Two mass models are used to fit the resultant shear profiles: a singular isothermal sphere (SIS) and a NFWprofile (Navarro et al. 1997). The SIS profile is the simplest density model for describing a relaxed massive sphere with a constant value for the isotropic one dimensional velocity dispersion, σ_V . This is mainly described by the density distribution: $$\rho (r)=\frac{{\sigma }_{\text{V}}^2}{2\pi G{r}^2},$$(3)

where G is the gravitational constant. This model corresponds to a distribution of self-gravitating particles where the velocity distribution at all radii is a Maxwellian with one dimensional velocity dispersion, σ_V . From this equation, we can get the critical Einstein radius for the source sample as: $${\theta }_{\text{E}}=\frac{4\pi {\sigma }_{\text{V}}^2}{c^2}\langle \beta \rangle ,$$(4)

in terms of which one obtains: $$\kappa (\theta )=\gamma (\theta )=\frac{{\theta }_{\text{E}}}{2\theta },$$(5)

where θ is the angular distance to the cluster center. Hence, fitting the shear for a different radii, we can estimate the Einstein radius, and from that, we can obtain an estimation of the mass M₂₀₀ as (Leonard & King 2010): $$M_{200}=\frac{2{\sigma }_{\text{V}}^3}{\sqrt{50}GH(z)},$$(6)

where H(z) is the redshift dependent Hubble parameter.

The NFW profile is derived from fitting the density profile of numerical simulations of CDM halos (Navarro et al. 1997). This profile depends on two parameters, the radius, R₂₀₀, that encloses a mean density equal to 200 times the critical density (ρ_crit ≡ 3H²(z)∕8πG), and a dimensionless concentration parameter, c₂₀₀: $$\rho (r)=\frac{{\rho }_{\text{c}}{\delta }_{\text{c}}}{(r/r_{\text{s}}){(1+r/r_{\text{s}})}^2},$$(7)

where r_s is the scale radius, r_s = R₂₀₀∕c₂₀₀ and δ_c is the characteristic overdensity of the halo, $${\delta }_{\text{c}}=\frac{200}{3}\frac{c_{200}^3}{\ln (1+c_{200})-c_{200}/(1+c_{200})}.$$(8)

We use the lensing equation for the spherical NFW density profile from Wright & Brainerd (2000). If we fit the shear for different radius we can obtain an estimation of the parameters c₂₀₀ and R₂₀₀. Once we obtain R₂₀₀ we can compute the M₂₀₀ mass. Nevertheless, there is a well-known degeneracy between the parameters R₂₀₀ and c₂₀₀ when fitting the shear profile in the weak lensing regime. This is due to the lack of information on the mass distribution near the cluster center and only a combination of strong and weak lensing can raise it and provide useful constraints on the concentration parameter. Due to the lack of strong lensing modeling in our sample, we follow van Uitert et al. (2012), Kettula et al. (2015) and Pereira et al. (2018), by fixing the concentration parameter according to the relation c₂₀₀ (M₂₀₀, z) given by Duffy et al. (2011): $$c_{200}=5.71{(M/2\times {10}^{12}{h}^{-1})}^{-0.084}{(1+z)}^{-0.47}.$$(9)

The particular choice of this relation does not have a significant impact on the final mass values, with uncertainties dominated by the noise of the shear profiles. Thus we fit the profile with only one free parameter: R₂₀₀.

To derive the parameters of each mass model profile we perform a standard χ² minimization: $${\chi }^2={\displaystyle \sum _i^N\frac{{(\langle \gamma (r_i)\rangle -\gamma (r_i,p))}^2}{{\sigma }_{\gamma }^2(r_i)}},$$(10)

where the sum runs over the N radial bins of the profile and the model prediction. p refers to either σ_V for the SIS profile, or R₂₀₀ in the case of the NFW model. Errors in the best-fitting parameters are computed according to the variance of the parameter estimate.

3.4 Two-dimensional density distribution

We obtain the 2D projected density distribution for the two analyzed galaxy systems by using the LensEnt2 code (Bridle et al. 1998; Marshall et al. 2002). This code applies maximum-entropy method algorithm for reconstructing the 2D density distribution from the two ellipticity components, e₁ and e₂, of the background galaxies and their respective uncertainties. We consider σ_ϵ ∕w as the uncertainty for both components, where w is the assigned weight to each galaxy, as was explained in Sect. 3.1. It is expected that galaxy clusters have smooth, extended projected mass distributions, therefore, an intrinsic correlation function (ICF) is included in the analysis, so the physical projected density distribution, Σ, is expressed as the convolution with a broad kernel, given by the ICF. In other words, the ICF smooths the resultant distribution. We adopt a circularly symmetric Gaussian as the ICF with a width ω. This width could be selected to maximize the Bayesian evidence, nevertheless this could lead to a sub-estimation of the substructure in the field. Hence, a visual inspection of the obtained maps and a comparison with other constraints (X-ray contours, galaxy member distribution and the position of the BCG) could be important in order to select an adequate value. Taking this into account we obtain the density distributions for different widths, ω, and then we fix this value considering the resultant Bayesian evidence and according to a visual inspection of the reconstructed maps. We fix ω at 550 and 410 kpc for A1204 and A2029/A2033, respectively.

Finally to obtain the projected density distribution, Σ, it is necessary tocompute the critical density defined as: $${\text{UpSigma}}_{\text{c}}=\frac{c^2}{4\pi G}\frac{1}{D_{\text{L}}\langle \beta \rangle },$$(11)

where we adopt the parameters D_L and ⟨β⟩ for the primary component.

The code also provides an error map according to the covariance matrix of the density reconstruction. We consider as the error in our 2D density distribution the mean of the values in the error map. Contours in the 2D distribution are built from a 3σ confidence level, while the distribution is shown from a 2σ confidence level.

4 Dynamical analysis

In a previous work de los Rios et al. (2016) studied the dynamical properties of the clusters A2029/A2033 and A1204 using photometric and spectroscopic information. They estimated some parameters that correlate with the dynamical status of the clusters and then, by using machine learning algorithms, classified them as candidates for merging systems (see Table 1 in de los Rios et al. 2016). In this section, using the publicly available Merging Systems Identification Algorithm² , we extend their analysis and describe the dynamical state of A1204 and A2029/A2033.

For the A1204 cluster we study its dynamical status with the merging system identification algorithm described in de los Rios et al. (2016) and we find that this is a merging cluster candidate. In order to study the stability of our result, we apply the same technique but deleting some randomly selected galaxies and then we analyze how the classification changes with the percentage of dropped galaxies. In particular, we perform 100 realizations for each value of completeness (99%, 95% and 90%). In all cases, we classify the system as a merging cluster in ~50% of the realizations, suggesting that the classification is not stable, therefore if there is an interaction between the components, it will not be strong. On the other hand, we also apply the identification algorithm separately to the two components identified by the algorithm in the first study. We found that each individual component is classified as not being in the process of merging, indicating that this is not a case of a multiple merger.

In the case of the A2029/A2033 cluster, we analyze the cluster identified by the FoF algorithm (i.e., Abell 2029 and Abell 2033 as a single cluster) with the Merging Systems Identification Algorithm. We find this system to be a merging cluster candidate and that A2029 and A2033 are the two substructures that are interacting. As in the case of A1204, we study the stability of the classification by reclassifying the dynamical status of the FoF cluster in random realizations dropping a percentage of member galaxies. In all the realizations, the algorithm classified the cluster as a merging system, hence, the classification is stable. We also analyze each system (A2029 and A2033) separately applying the Merging Systems Identification Algorithm and we obtain that each cluster is classified as not being in the process of merging. This is in agreement with previous classifications of the individual clusters as relaxed systems (Wen & Han 2013; Mantz et al. 2016).

5 Results and discussion

In this section we present the results of the lensing analysis, together with the X-ray and dynamical characterization of the systems. In Table 2 we show the obtained lensing masses for each structure according to the fitted shear profiles (Figs. 3 and 4). In order to characterize the gas distribution, we compute brightness contour maps at (3, 5, 7, 9, 12) times the background and smoothed over 4 pixels. The derived 2D projected density distributions together with X-ray contours are shown in Fig. 5.

Fig. 3 Shear radial profiles as a function of cluster-centric projected distance for the primary (upper panel) and the secondary component (bottom panel) of A1204, identified by de los Rios et al. (2016). The solid and the dashed lines represent the best fit of SIS and NFW profiles, respectively, with the fitted parameters given in the right-upper corner. The points and crossings show the tangential and cross ellipticity components of the selected background galaxies, averaged in annular bins, respectively. Error bars are computed according to Eq. (2).

Fig. 4 As in Fig. 3 but for A2029 (upper panel) and A2033 (bottom panel).

Fig. 5 Projected density distribution in the field A1204 (upper panel) and in the field of A2029 and A2033 (lower panel) obtained with the weak lensing analysis. The scale, marked at bottom left, is given in h70   M⊙   pc−2. Red contours corresponds to a projected density above 3σ significance level (250 h70 M⊙ pc−2 for A1204 and 400 h70 M⊙ pc−2 for A2029/A2033) with steps of 50 M⊙ pc−2. X-ray contours are plotted in blue, the contour levels are (3, 5, 7, 9 and 12) times the rms noise. Red dots and green squares are the BCGs positions and the dynamical centres, respectively. Pink and light-blue points are the positions of the galaxies classified as members, for the primary and secondary component of A1204 and for A2029 and A2033, respectively.

5.1 A1204

For this system we obtain the lensing masses of both identified substructures and the ratio between derived masses is ~ 2. Babyk et al. (2012) computed the mass for A1204 by using X-ray Chandra observations and assuming a NFW profile, obtaining $M_{200}={3.18}_{-0.24}^{+0.34}\times {10}^{14}{h}_{70}^{-1}{M}_{\odot }$. This result is in good agreement with the estimated mass in this work ($M_{200}=4.0\pm 1.8\times {10}^{14}{h}_{70}^{-1}{M}_{\odot }$).

In the 2D density distribution (Fig. 5, upper panel) we can distinguish only the primary structure. Taking into account the measured mass according to the fitted shear profile for this component, we compute the corresponding projected density profile, Σ(r), assuming a NFW distribution according to the equations given by Wright & Brainerd (2000). All of the obtained values are below the 3σ detection level (250 h₇₀ M_⊙ pc⁻²) and only the inner region, corresponding to a radius of 60 kpc would present a density larger than the 2σ detection level (180 h₇₀ M_⊙ pc⁻²). Thus, we consider that the mass of this structure is near the detection threshold of the 2D density contrast distribution. Also, we do not detect significant X-ray emission in the secondary structure region above the threshold adopted to build the brightness contours. In order to establish if this is due to the observing detection limit, we consider the lowest detected flux of RASS-based catalogs according to Piffaretti et al. (2011), which contains the lowest X-ray emission clusters identified using this data. This corresponds to a flux F_lim = 1.5 × 10⁻¹² erg s⁻¹ cm⁻², in the 0.1–2.4 keV band of the SGP catalog (Cruddace et al. 2002). Considering this flux and the redshift of the secondary structure, we obtain a limiting luminosity of 1.5 × 10⁴² $h_{70}^{-2}$ erg s⁻¹. Taking into account the obtained lensing mass, we compute the expected X-ray luminosity for the system according to Leauthaud et al. (2007) M_ L_X relation, obtaining (3.7 ± 2.4) × 10⁴³ $h_{70}^{-2}$ erg s⁻¹, only 1.5σ above the detection limit. For this calculation we use the derived SIS mass since it has the lowest χ² value and errors correspond to the propagation of the M_L_X parameter errors. Therefore, according to the lensing estimated mass the secondary structure would be roughly at the detection limit of the X-ray observation. We considered using more sensitive X-ray data, therefore we checked in the available databases of Chandra and XMM; however this component is not present in the field-of-view of these surveys. Considering the X-ray luminosity emission upper limit (6.1 × 10⁴³ $h_{70}^{-2}$ erg s⁻¹) the emission is 2.5σ above the detection limit given by F_lim (Piffaretti et al. 2011). Therefore, it is important to highlight that a low X-ray emissivity could be explained by a low density of the intracluster gas, which might be produced by a past interaction between the structures.

The cluster A1204 has recently been classified as a system with a strong cool core (Zhang et al. 2016) which is in good correspondence with the X-ray contours (Fig. 5, upper panel). Thus, the cluster does not show any evidence of having suffered a recent merger event. Also, for this cluster the classification as a merging system is unstable as was stated in Sect. 4. Nevertheless, we detect lensing signal for the secondary structure through the lensing profile (Fig. 3) from which we obtain the total mass, and, also, the density contours exhibited in Fig. 5 are elongated in the direction of the secondary component. Furthermore, assuming a relative velocity for these systems similar to the velocities measured in other merging clusters (1000 km s⁻¹, Thompson & Nagamine 2012), the time necessary for the structures to reach the observed distance between the centers of the identified substructures (~ 2.4 Mpc) is ≤ 3 Gyr. Hence, taking into account that there are no other signatures of a collision, it is not likely that this could be a past merger event. Nevertheless this scenario could not be discarded if we were to consider lower velocities for the components; we also do not discard the possibility that the interaction between these structures could be in process.

5.2 A2029/A2033

For this system we obtain the total mass for both structures which can be identified in the 2D projected density distribution. Both clusters show X-ray emission in good correspondence with the density distribution (Fig. 5). The obtained total mass for A2029 ($M_{200}=9\pm 6\times {10}^{14}{h}_{70}^{-1}{M}_{\odot }$) is in good agreement with that estimated by Walker et al. (2012) using X-ray observations ($M_{200}=10.1\pm 0.6\times {10}^{14}{h}_{70}^{-1}{M}_{\odot }$).

A2029 is a relaxed galaxy cluster which has been studied extensively in X-rays (e.g., Sarazin et al. 1998; Lewis et al. 2002; Clarke et al. 2004; Vikhlinin et al. 2005; Snowden et al. 2008). It has a large cD galaxy (Uson et al. 1991) whose major axis is aligned in the NE to SW direction, in approximately the same direction as that joining it to nearby A2033. These two systems, together with A2028 and A2066 form a small supercluster (Einasto et al. 2001). Studies examining whether or not these systems are connected by a filament structure show that this is not the case (Walker et al. 2012; Planck Collaboration VIII 2013). Therefore, there is no conclusive evidence showing that these systems have interacted in the past or that they are now interacting, given the observed relaxed state of both structures. Nevertheless, if we take into account that they belong to the same supercluster and that the dynamical classification as interacting system is stable (see Sect. 4), it could be expected that these clusters would interact.

5.3 Two-body model

In orderto obtain information regarding the state of evolution of the studied systems, we apply a Newtonian gravitational binding criterion that the two-body system is bound if the potential energy of the system is equal to or greater than the kinetic energy. This dynamical model was described in detail by Beers et al. (1982) and Gregory & Thompson (1984) and also applied to the analysis of several bimodal galaxy systems (e.g., Cortese et al. 2004; Hwang & Lee 2009; Yan et al. 2014; Andrade-Santos et al. 2015; Caglar & Hudaverdi 2017). This model assumes radial orbits for the identified structures, which start their evolution at time t₀ = 0 with separation R₀ = 0, and are moving apart or coming together for the first time in their history. With the obtained lensing masses, this method allows us to estimate the probability that (1) the system is bound but still expanding, (2) the system is collapsing, or (3) the two structures are not bound to each other but are merely close together on the sky. It is important to highlight that this model does not considerthe angular moment of the system, the distribution of matter inside each cluster, or the gravitational interaction of the infalling matter outside the cluster pair, since it is assumed that the masses are constant since their formation time (Nascimento et al. 2016). Nevertheless, this model puts another constraint on the dynamical state and the evolution of the systems.

The Newtonian criterion for gravitational binding can be stated in terms of the projected separation, R_p , the radial velocity difference, V _r and the total mass of the combined system, M, as: $$V_{\text{r}}^2{R}_{\text{p}}\le 2GM{\sin }^2(\alpha )\cos (\alpha ),$$(12)

where α is the projection angle between the plane of the sky and the line connecting the components of the system. We can express the true (3D) velocity, V , and separation, R, as: $$R_{\text{p}}=R\cos (\alpha ),V_{\text{r}}=V\sin (\alpha ).$$(13)

We compute V _r considering the redshifts of the galaxy members of each component within the R₂₀₀ radius, computed according to SIS masses. R_p is computed as the distance between each BCG (which is the adopted center for the lensing analysis) and the combined mass of the system is obtained by adding the SIS masses of each structure, since these have the lower errors.

Using Beers et al. (1982) equations of motion for unbound and bound systems, we can obtain the projected angle α as a function of the radial velocity difference V _r. The equation relates the time (which is assumed for each system as the age of the Universe at the mean redshift of the considered structures), the velocity, V , and the separation, R, with the developmental angle, χ, and the maximum separation, R_m, of the system components for the bound solution and the asymptotic expansion velocity for the unbound solution. The obtained V _r − α relations areplotted in Fig. 6. The solid black curve separates the bound and unbound solutions according to the Newtonian criterion (Eq. (12)), the blue curves are the solutions of the equations of motions for bound (solid line) and unbound systems (dashed line) and the solutions for each system are marked with open circles, according to the observed V _r . For both systems the solutions are bound, defined by the intersections between the solutions of the equations of motion and the V _r . For each solution, i, we compute the probabilities given by: $$p_i={\displaystyle {\int }_{{\alpha }_{\text{inf},i}}^{{\alpha }_{\text{sup},i}}\cos }(\alpha )\text{d}\alpha ,$$(14)

then we normalize the obtained probability to obtain P_i = p_i∕(∑ _ip_i).

Solutions are presented in Table 3. For A2029 we obtain one bound-outgoing solution, close to the unbound solution considering the errors inV _r. According to this, the components are expanding separated by a distance of $13.6h_{70}^{-1}$ Mpc. On the other hand, for A1204 we obtain three bound solutions, two incoming and one outgoing. The most probable is the incoming solution for which the components have reached their maximum distance of $4.4h_{70}^{-1}$ Mpc and are now separated by $2.4h_{70}^{-1}$ Mpc collapsing with a velocity of 1000 km s⁻¹.

Fig. 6 Projection angle (α) as a functionof the relative radial velocity difference (V r) computed according to the equations of motion for A1204 (left panel) and A2029/A2033 (right panel). Solid black curves separate the bound and unbound solutions according to the Newtonian criterion. Blue curves are the solutions of the equations of motions for bound (solid line) and unbound systems (dashed line) and the solutions for each system are marked with open circles, according to the observed V r marked with the red line where the gray region corresponds to the uncertainty in these values delimited by the dashed black line. Point curves express the uncertainties in the computed curves considering the errors in the adopted lensing masses.

6 Summary and conclusions

In this work we present the analysis of two candidates for interacting galaxy systems, A1204 and A2029/A2033. For this analysis we include spectroscopic information of the galaxies identified as members, X-ray observations, and a weak-lensing study.

The dynamical analysis suggests that A1204 classification as a merging system is not stable and that each detected substructure is dynamically identified as not being in the process of merging. From the lensing analysis we were able to obtain the total mass of each component and the 2D projected density map. In the latter we observe a concentric distribution of the mass centered at the BCG of A1204, compatible with the gas distribution, according to the X-ray brightness contours. Therefore, this component is relaxed, hence there is no evidence that A1204 has passed through a recent merger event. Nevertheless, according to the two-body model, the identified structures are bound and collapsing. Therefore the identified components could interact in the future.

On the other hand, A2029/A2033 classification as a merging system is stable according to the dynamical analysis and each cluster is classified separately as not being in the process of merging. By using our lensing analysis we derived the masses for each cluster and both components are visible in the 2D projected density distribution. The gas distribution given by the X-ray contours are in good agreement with the mass distribution and both systems seem to be relaxed. Hence, we conclude that for this system there is no evidence that these clusters have passed through a recent merger event. Furthermore considering the result of the two-body model, these clusters are not in the process of merging, since the solution describes an outgoing motion for the components.

Although the analyzed systems were classified as candidates for interacting galaxy clusters, they do not show evidence of having passed through a recent merger event. However, A1204 seems to be in process of interaction. Thus in order to accurately classify the dynamical state of the systems, more constraints need to be taken into account. With this purpose we will continue the analysis of the identified systems in order to improve the algorithms of detection and also to provide a better understanding of the classified systems. Nevertheless, we highlight the performance of the merging system identification algorithm and the forthcoming works that would be based on its resultantcatalog, since it would provide a uniform selected and analyzed sample of these interesting galaxy systems.

Acknowledgements

We thank the anonymous referee for the useful comments and suggestions which helped to improve the manuscript. This research is based on data collected at Subaru Telescope, which is operated by the National Astronomical Observatoryof Japan and on observations obtained with MegaPrime/MegaCam, a joint project of CFHT and CEA/IRFU, at the Canada-France-Hawaii Telescope (CFHT) which is operated by the National Research Council (NRC) of Canada, the Institut National des Science de l’Univers of the Centre National de la Recherche Scientifique (CNRS) of France, and the University of Hawaii. This work is based in part on data products produced at Terapix available at the Canadian Astronomy Data Centre as part of the Canada-France-Hawaii Telescope Legacy Survey, a collaborative project of NRC and CNRS. Ww have made use of data and/or software provided by the High Energy Astrophysics Science Archive Research Center (HEASARC), which is a service of the Astrophysics Science Division at NASA/GSFC and the High Energy Astrophysics Division of the Smithsonian Astrophysical Observatory. We have made use of the ROSAT Data Archive of the Max-Planck-Institut für extraterrestrische Physik (MPE) at Garching, Germany. We have also use SDSS-III data. Funding for SDSS-III has been provided by the Alfred P. Sloan Foundation, the Participating Institutions, the National Science Foundation, and the U.S. Department of Energy Office of Science. The SDSS-III web site is http://www.sdss3.org/. SDSS-III is managed by the Astrophysical Research Consortium for the Participating Institutions of the SDSS-III Collaboration including the University of Arizona, the Brazilian Participation Group, Brookhaven National Laboratory, Carnegie Mellon University, Universityof Florida, the French Participation Group, the German Participation Group, Harvard University, the Instituto de Astrofísica de Canarias, the Michigan State/Notre Dame/JINA Participation Group, Johns Hopkins University, Lawrence Berkeley National Laboratory, Max Planck Institute for Astrophysics, Max Planck Institute for Extraterrestrial Physics, New Mexico State University, New York University, Ohio State University, Pennsylvania State University, University of Portsmouth, Princeton University, the Spanish Participation Group, University of Tokyo, University of Utah, Vanderbilt University, University of Virginia, University of Washington, and Yale University. We made an extensive use of the following python libraries: http://www.numpy.org/, http://www.scipy.org/, http://roban.github.com/CosmoloPy/ and http://www.matplotlib.org/. This work was partially supported by the Consejo Nacional de Investigaciones Científicas y Técnicas (CONICET, Argentina) and the Secretaría de Ciencia y Tecnología de la Universidad Nacional de Córdoba (SeCyT-UNC, Argentina). J.L.N.C. is grateful of the Southern Office of Aerospace Research and Development (SOARD), a branch of the Air Force Office of the Scientific Researchs International Office of the United States (AFOSR/IO), for the financial support obtained through grant FA9550-15-1-0167. J.L.N.C. also acknowledges financial support from the Dirección de Investigación de la Universidad de La Serena through the Programa de Incentivo a la Investigación Académica (PIA-DIULS).

References

All Tables

All Figures

Fig. 1 Classification of the detected sources in the Subaru frame (upper panel) and in one of the frames of the CFHT mosaic (lower panel). Here stars are represented with blue points, galaxies by green points and false detections by black crosses. Left panels show the central surface brightness of the objects (μMAX ) as a function of the magnitude. Stars are confined in the region marked by the red lines. Right panels show the source distribution in the diagram FWHM as a function of the magnitude. As it can be seen, the FWHM of the sources classified as stars remains roughly invariant for the total range in magnitude.

In the text

Fig. 2 PSF treatment applied to stars in the Subaru frame (upper panel) and in one of the frames of the CFHT mosaic (lower panel). Semi-major axis distribution in the CCD, of the sources classified as stars, before (left) and after (right) the PSF deconvolution. For a better visualization of the plot, semi-major axes have a different scale than the CCD axes. Semi-major axis scale is given by the red thicker segment in the upper left-hand corner and corresponds to 3 pixels. We highlight that the semi-major axis is more randomly distributed and is significantly smaller after taking into account the PSF.

In the text

Fig. 3 Shear radial profiles as a function of cluster-centric projected distance for the primary (upper panel) and the secondary component (bottom panel) of A1204, identified by de los Rios et al. (2016). The solid and the dashed lines represent the best fit of SIS and NFW profiles, respectively, with the fitted parameters given in the right-upper corner. The points and crossings show the tangential and cross ellipticity components of the selected background galaxies, averaged in annular bins, respectively. Error bars are computed according to Eq. (2).

In the text

	Fig. 4 As in Fig. 3 but for A2029 (upper panel) and A2033 (bottom panel).
In the text

Fig. 5 Projected density distribution in the field A1204 (upper panel) and in the field of A2029 and A2033 (lower panel) obtained with the weak lensing analysis. The scale, marked at bottom left, is given in h70   M⊙   pc−2. Red contours corresponds to a projected density above 3σ significance level (250 h70 M⊙ pc−2 for A1204 and 400 h70 M⊙ pc−2 for A2029/A2033) with steps of 50 M⊙ pc−2. X-ray contours are plotted in blue, the contour levels are (3, 5, 7, 9 and 12) times the rms noise. Red dots and green squares are the BCGs positions and the dynamical centres, respectively. Pink and light-blue points are the positions of the galaxies classified as members, for the primary and secondary component of A1204 and for A2029 and A2033, respectively.

In the text

Fig. 6 Projection angle (α) as a functionof the relative radial velocity difference (V r) computed according to the equations of motion for A1204 (left panel) and A2029/A2033 (right panel). Solid black curves separate the bound and unbound solutions according to the Newtonian criterion. Blue curves are the solutions of the equations of motions for bound (solid line) and unbound systems (dashed line) and the solutions for each system are marked with open circles, according to the observed V r marked with the red line where the gray region corresponds to the uncertainty in these values delimited by the dashed black line. Point curves express the uncertainties in the computed curves considering the errors in the adopted lensing masses.

In the text