© ESO 2018

1. Introduction

Over the past decades we have obtained an increasingly clear picture of the formation and evolution of galaxies in the Universe. The galaxy population can be broadly divided in two distinct types. Star-forming galaxies have blue colours, typically disk-like morphologies and a relatively high star formation rate, whereas quiescent galaxies have redder colours, more spheroidal morphologies, and a (near) absence of star formation. Generally speaking, star-forming galaxies are found to dominate in abundance at relatively early times and at low stellar masses (Kauffmann et al. 2004; Baldry et al. 2006; Peng et al. 2010; Muzzin et al. 2013a). A central question is how galaxies transform from actively star-forming systems to passive quiescent galaxies. An important factor in this quest is to understand how the environment of a galaxy affects this transformation process, since at fixed stellar mass, galaxies in groups and clusters are found to have a higher probability of being quiescent than galaxies in the field (Dressler 1980; Blanton et al. 2005; Woo et al. 2013).

In recent years it has been shown that the effect of the environment is largely separable from the quenching processes that act internally, at least in the local universe (Baldry et al. 2006; Peng et al. 2010; Kovač et al. 2014). Specifically, Peng et al. (2010) introduced terms of “mass quenching”, which would be operating independently of the environment and most effectively quench massive galaxies, and “environmental quenching”, which would be operating independently of stellar mass and quench galaxies preferentially in overdense regions (also see e.g. Davies et al. 2016; Kawinwanichakij et al. 2017). The total quenching effect would simply be the product of these two contributions.

The physics behind the two separate quenching routes is a matter of debate. “Mass quenching” is often associated with feedback, either from supernovae and galactic winds (Oppenheimer et al. 2010), or from active galactic nuclei (Bower et al. 2006; Croton et al. 2006; Tremonti et al. 2007). Many mechanisms have been proposed to be responsible for “environmental quenching”. In the overdensities of groups and clusters, stripping of cold gas (through ram-pressure: Gunn & Gott 1972), or hot gas (“strangulation”/“starvation”: Larson et al. 1980), as well as galaxy harassment or (dry) mergers (Moore et al. 1996), could be responsible for the observed trends. To make progress in understanding these environmental processes, it is essential to understand in which environments and on which timescales they operate (Balogh et al. 2004).

Most studies that aim to understand the process of environmental quenching are focussed on large cosmological volumes in, for example, the COSMOS field (Peng et al. 2010; Darvish et al. 2016), VIPERS (Davidzon et al. 2016), 3D-HST (Fossati et al. 2017), the Subaru Hyper Suprime-Cam survey (Jian et al. 2018), and ZFOURGE (Kawinwanichakij et al. 2017; Papovich et al. 2018). Most of these studies separate the galaxy population into four density quartiles, so that the environmental effects can be studied between the different quartiles. At least at high redshift (z  ≳  1) it is found that the environmental quenching process is not working completely independently of stellar mass (Papovich et al. 2018). A starvation/strangulation scenario in which the supply of hot gas is cut off from a galaxy would be highly effective at these redshifts, where star formation rates and outflows deplete the cold gas supply (the “overconsumption model”, McGee et al. 2014). This would likely introduce a mass-dependent effect, since the gas depletion time of higher-mass galaxies ought to be shorter, and this could explain the measured trends. Also for lower mass galaxies (M_⋆  ≲  10^9.5 M_⊙), measured quenching timescales are comparable to the gas depletion timescale (Davies et al. 2016).

Whereas such a division in density quartiles allows for a study of environmental quenching in reasonably overdense regions, the most extreme environments are either not probed, or are washed out by more moderate overdensities. And yet it is in these regions where the physics of the quenching may be notably different (Balogh et al. 2016; Kawinwanichakij et al. 2017). Quenching processes may be more violent, leading to very high quenched fractions of cluster galaxies compared to the field at the same redshift (De Lucia et al. 2004; van der Burg et al. 2013; Annunziatella et al. 2014; Balogh et al. 2016; Nantais et al. 2016, 2017; Wagner et al. 2017). Dynamical processes in clusters may become more important at later time to quench star formation, and these may act in a mass-independent fashion. A notable example of this is stripping by ram pressure, which would directly remove the cold gas supply from a galaxy and quench its star formation on a very short timescale (Zinger et al. 2018; Fossati et al. 2016; Bellhouse et al. 2017; Jaffé et al. 2018).

To understand where and when the most extreme environmental quenching is taking place in clusters, several studies have focussed on a clustercentric-distance dependent study, some even focussing on projected phase space of different galaxy populations (such as those of galaxies in the “transition” phase; Oman et al. 2013; Muzzin et al. 2014; Poggianti et al. 2016; Jaffé et al. 2018). To study the relative excess of already-quenched galaxies in clusters compared to the field, single values of the environmental quenching efficiency are typically reported for satellite galaxies (cf. Fig. 7 in Nantais et al. 2016, and references therein), even though there is likely to be a substantial trend with radius. Such a study as a function of radius can probe the “backsplash” of already-quenched ejected cluster satellites (Wetzel et al. 2014), and pre-processing/quenching of future cluster satellites in the surrounding large-scale overdensity (Fujita 2004).

While new surveys push the frontiers of these studies to higher redshifts (e.g. Balogh et al. 2017), their samples are thus-far limited in size, and clusters are of moderate mass and overdensity. This renders a radial-dependent study of environmental quenching difficult. In this paper we focus instead on highly massive and thus over-dense clusters at intermediate redshifts. This allows us to specifically study the environmental quenching excess as a function of radial distance from the cluster centres. The sample we study is composed of 21 massive clusters detected with the Planck SZ survey at redshifts 0.5 <  z <  0.7. Since Planck is an all-sky survey (even though we only consider the northernmost two thirds here), we probe the highest-density environments at these intermediate redshifts, where environmental effects are expected to be substantial. We concentrate on a photometric data set spanning from the u- to the K_s-band for each cluster, which we use to estimate stellar masses for individual galaxies, and separate them by type based on their best-fit SED.

Our starting point is a measurement of the galaxy stellar mass function (SMF), which describes the number density of galaxies as a function of their stellar mass, and which is a key observable to study the formation and evolution of galaxies. By comparing the SMF to the underlying halo mass function, the efficiency with which galaxies form can be measured, and this is an essential test and diagnostic tool for large hydrodynamical simulations such as Illustris (Genel et al. 2014) and EAGLE (Schaye et al. 2015). Measuring the SMF in these massive clusters provides further constraints for the next generation of hydrodynamical simulations, in which large overdensities can be specifically focussed on to study the influence of such environments on the evolution of galaxies (e.g. Bahé et al. 2017). These SMFs are the main ingredients to estimate environmental quenching locations and timescales, which we describe and discuss in the remainder of this work.

The structure of this paper is as follows. Section 2 describes our cluster sample and photometric data set. Section 3 lays out the main analysis, ranging from photometric redshift measurements to a statistical accounting of fore- and background galaxies in our cluster galaxy sample. The main results, measurements of the SMF and environmental quenching efficiency, are presented in Sects. 4 and 5, respectively. We discuss our findings in Sect. 6, and conclude in Sect. 7.

All magnitudes we quote are in the AB magnitudes system, and we adopt ΛCDM cosmology with Ω_m  =  0.3, Ω_Λ  =  0.7 and H₀  =  70 km s⁻¹ Mpc⁻¹. Uncertainties are given at the 1−σ level, unless explicitly stated otherwise.

2. Cluster sample and data

The clusters we study are drawn from a sample of 33 clusters that were detected with Planck, and confirmed by autumn 2011 to be at z >  0.5. This sample was the target of an XMM-Newton Large Programme “Unveiling the most massive galaxy clusters at z >  0.5 with Planck and XMM-Newton” (PI M. Arnaud) in AO-11. The properties of the sample, in particular regarding their morphological properties, Intra Cluster Medium and the cluster scaling relations, are outlined in Arnaud et al. (in prep.).

In this paper, we study the galaxy content in a sub-sample of 21 clusters that make up the northernmost (Dec >  −25^∘) part of this parent sample. Several of these clusters have already been studied, particularly as part of X-ray selected samples of clusters in this redshift range (Böhringer et al. 2000; Ebeling et al. 2007; Piffaretti et al. 2011). Several are in optical catalogues constructed using SDSS data (Wen et al. 2012; Rykoff et al. 2014).

Table 1 presents the main characteristics of the sample, making a comparison of the mass estimated from the Planck SZ signal and mass based on the deep X-ray maps (M − Y_X relation). Even though the SZ mass proxy is blindly extracted (i.e. without prior knowledge on the location of the cluster; cf. Planck Collaboration Int. XXXVI. 2016), both proxies are consistent at the massive end (within ∼10% in mass); they slightly diverge at the low-mass end due to Eddington bias in the SZ mass proxy (cf. e.g. van der Burg et al. 2016). The differences between Y_X and Y_SZ will be discussed in more detail in Arnaud et al. (in prep.). The current paper refers to cluster masses as M₅₀₀ ¹, based on the Y_X scaling relation.

To support our analysis, we combine different sources of spectroscopic information in the 21 fields we study. Fourteen of the clusters are covered in DR14 of SDSS (Abolfathi et al. 2018). For nine clusters, we have obtained redshifts with the Nordic Optical Telescope (Planck Collaboration XXVII 2016, Dahle et al., in prep.) or Gemini (Planck Collaboration XXIX 2014). Ebeling et al. (2014) publish a catalogue with hundreds of spectroscopic redshifts in two of the fields we study (PSZ2 G180.25+21.03 and PSZ2 G228.16+75.20). For PSZ2 G046.13+30.72 and PSZ2 G155.27-68.42, we have obtained spectroscopic redshifts from a program undertaken with the Canary Islands observatories (Planck Collaboration Int. IV 2013). We searched the NED database² for any spec-zs we may have missed in the literature. The only two clusters that remain without a previously measured spectroscopic redshift (PSZ2 G193.31-46.13 and PSZ2 G219.89-34.39) have been observed using VLT/FORS2 multi-object spectroscopy (PID = 090.A-0925, PI=Bohringer). Using our own custom pipeline, we reduced these spectra and measured the redshifts of several member galaxies. The number of (unique) spectroscopic targets and cluster members for all clusters are listed in Table 1.

2.1. Cluster photometry

Deep archival follow-up imaging data are available for many of these clusters (in particular data retrieved from the SMOKA science archive; Baba et al. 2002), but these tend to be drawn from X-ray selected samples. Our analysis benefits from a homogeneous wavelength coverage of the full sample of 21 clusters. Studying the full sample not only statistically enhances the results presented in this study, but also ensures that we sample the full range of cluster dynamical properties, approximating a mass-selected sample. Since the spectroscopic data are of varying quality and completeness, we base our analysis almost entirely on the photometry. To be able to measure accurate and precise photometric redshifts of galaxies in the cluster fields, we require photometric coverage with at least seven filters per cluster, ranging from the u-band (∼3000 Å) to the K_s-band (∼22 000 Å). Accounting for the deep archival data for a subset of the sample, we obtained the remaining photometry through different time allocations on the wide-field imagers at CFHT (PI vdBurg, PIDs: 15AF006, 15BF005, 16BF013) and Subaru (PI Dahle, PIDs: S12B-164S, S13A-120, S15A-118).

An overview of all imaging data is given in Table 2. The photometric data we use have been taken over more than 10 yr, and some of the instruments and filter sets have been upgraded over this time span. These differences are indicated in the footnotes of Table 2, and the different wavelength responses of each filter and instrument are taken into account in our analysis.

We perform basic steps to reduce the optical imaging data, such as bias, flat-field correction and cosmic-ray removal. As an additional step, we remove background patterns, particularly fringe residuals, by using the dithered pattern of observations to differentiate signals that are fixed in position on the ccd array from sky-bound signals. This is explained and illustrated in Fig. 1 of van der Burg et al. (2016). For the near-infrared (NIR) data from our own WIRCam program (PI vdBurg) we have followed an extended dither pattern, where the cluster centre is dithered from chip to chip. This ensures a very clean background subtraction, even on scales of the intra-cluster light (ICL).

Astrometric registering has been performed with SCAMP (Bertin 2006) using the USNO-B1 catalogue as reference, or with external catalogues from Pan-STARRS (Chambers et al. 2016) for clusters that have been observed after the public release of their data. The astrometric precision between filters is better than 0.10″, ensuring that colour measurements can be performed accurately.

We automatically place masks on bright stars based on their locations in the guide star catalogue 2.3 (Lasker et al. 2008). After this, the images are inspected manually, and masks are placed on additional diffraction spikes, reflective haloes, and image artefacts. During our analysis we take account of the reduced effective areas after masking.

2.2. Object detection and colours

We perform object detection in the original K_s-band. Since the range in M/L between the different galaxy types is smallest in this band, this ensures a catalogue that is close to being stellar-mass selected. We use SExtractor to detect objects, following the criterion that at least five adjacent pixels have a flux density that is > 1.5σ relative to the local background RMS.

To be able to measure colours of the same intrinsic part of the galaxies we study, differences in the PSF between clusters and filters need to be accounted for. We use PSFEx (Bertin 2011) to determine convolution kernels that homogenise the PSF for each cluster between different filters. We then measure colours for the K_s-band-selected sources by performing photometry in circular apertures with a 2″ diameter on these PSF-matched images.

A benefit of these wide-field images is that they contain a large population of Galactic stars that can be used to calibrate the flux scale (e.g. High et al. 2009; Gilbank et al. 2011; van der Burg et al. 2013; Kelly et al. 2014). We use the effective wavelength-response curves of each detector, filter³ and atmosphere model in combination with the stellar spectral library of Pickles (1998) to create a reference stellar locus. The stellar spectral library we use⁴ is updated with additional information in the NIR from Ivanov et al. (2004), in addition to the original library from Pickles (1998). After applying magnitude offsets to our instrumental magnitudes to match this reference locus, we reach a photometric calibration with an uncertainty that depends on the filter, but is generally ≲0.05 magnitudes. While this calibrates the colours of sources in our catalogues, we perform absolute flux calibration in the K_s-band with respect to the 2MASS all-sky reference catalogue (Cutri et al. 2003).

Uncertainties on aperture flux measurements of faint galaxies are dominated by fluctuations in the sky background. We estimate this noise component by randomly placing apertures on sky positions that do not overlap with sources detected in the K_s-band. This procedure takes into account the correlated noise properties between adjacent pixels that originated from the convolution and the re-binning of the data to a common grid and PSF.

2.3. UltraVISTA reference field

There are two reasons to study a reference field, that is, a field without a massive cluster, in our work. First, to study the impact on the evolution of galaxies by their massive host haloes, properties of galaxies at the same redshift between cluster and field are compared. Second, to study the properties of galaxies that are part of the massive clusters in our sample, we have to consider projection effects; i.e. fore- and background galaxies that enter our sample of cluster galaxies. By performing the same selection criteria on galaxies in a parallel field, we can take these projection effects into account statistically. We refer to this as “background subtraction”. The more information (photometric redshift, distance from cluster centre) one can use, the more cleanly the background can be accounted for.

We make use of the COSMOS/UltraVISTA field, which has been extensively studied and for which a multi-band photometric catalogue is publicly available (Muzzin et al. 2013b). How we use this catalogue for the background correction, taking cosmic variance uncertainties due to this single field into account, is described in the Sect. 3.3.

3. Analysis

3.1. Photometric redshifts

We use the template-fitting code EAZY (Brammer et al. 2008) to estimate photometric redshifts (photo-zs) for each object. The photo-zs correspond to the peak (z_peak) of the posterior probability distribution P(z) given by EAZY. Figure 1a shows a comparison between the spectroscopic redshift and photometric redshift. We define a relative scatter $\mathrm{\Delta }z=\frac{z_{\mathrm{phot}}-z_{\mathrm{spec}}}{1+z_{\mathrm{spec}}}$ for each object with a reliable spectroscopic redshift z_spec. There are ∼3.5% outliers, defined as objects for which |Δz|> 0.15. For the remaining galaxies, we measure the mean of Δz and the scatter around this mean, σ_z, finding a bias of |Δz|  =  0.008 and scatter of σ_z = 0.029.

Since the galaxies for which spectroscopic redshifts have been measured are generally bright, and preferentially have emission lines, it is not immediately clear if the reported performance is representative of the entire galaxy sample down to the detection limit. In Fig. 1b we separate the galaxies by class (star-forming vs. quiescent; cf. Sect. 3.5), and plot the differences as a function of K_s-band magnitude, restricted to the redshift range where our clusters are located, 0.5 <  z_spec <  0.7. Within this spectroscopic redshift range, the photo-z scatter is σ_z = 0.028 for star-forming galaxies, and σ_z = 0.023 for quiescent galaxies.

This figure shows that the success rate of measuring spectroscopic redshifts of faint galaxies is higher for star-forming than for quiescent galaxies, since the former have typically strong emission lines. The scatter, shown by the solid curve, based on a running bin width of 1 magnitude, does not significantly depend on magnitude. This suggests that the flux measurements that define the SED are precise, and the photo-z determination is limited by the representativity of the templates, the accuracy of the filter curves, and the accuracy of the flux calibration. Since photo-zs for quiescent galaxies are slightly more precise than for star-forming galaxies (due to a stronger spectral feature around the characteristic 4000 Å break), we expect the photo-z performance of quiescent galaxies to be at least similar to that of star-forming galaxies at faint magnitudes.

We use the broad-band colours to identify and flag stars in our catalogue (e.g. Whitaker et al. 2011; van der Burg et al. 2013), without having to make a selection based on size or morphology. That is because galaxies have very different spectral energy distributions (SEDs) from stars, particularly towards NIR wavelengths.

We use the following colour criterion, which is similar to the ones used in aforementioned studies, to select the sample of galaxies:$$\begin{array}{c}J-K_{\text{s}\hspace{0.17em}}>0.18\times (u-J)-0.60\cup J-K_{\text{s}\hspace{0.17em}}>0.08\times (u-J)-0.30\end{array}$$(1)

Fig. 1. Left panel: spectroscopic vs. photometric redshifts for the 21 cluster fields. Outliers, objects for which |Δz|> 0.15, are marked in red. The outlier fraction is 3.5%, the scatter of the remaining objects is σz  =  0.029. Right panel: scatter in Δz as a function of Ks-band magnitude, for sources with 0.5 <  zspec <  0.7. Quiescent and star-forming galaxies are marked in red and blue, respectively. The source detection limits of the 21 Ks-band stacks are indicated.

3.2. Identifying brightest cluster galaxies

We select the brightest cluster galaxy (BCG) in the K_s-band that is located within 1′ of the X-ray peak, and that has a photometric redshift within 0.10 of the cluster redshift. In all but one case this automatic identification corresponds to what we would have selected by hand based on the colour images in Appendix B. For PSZ2 G219.89-34.39, there is a misidentification since the apparent BCG has a blue core, and is likely contaminated by blue light from a nearby source (cf. Fig. B.1, as is also apparent from the VLT/FORS2 spectrum of this galaxy). The photo-z is 0.56, which is 0.17 lower than the cluster redshift. We select only this BCG by hand, and all the others following the criteria above. The BCGs are marked in Fig. B.1, and their positions are reported in Table 1.

3.3. Stellar masses and background correction

We measure stellar masses for each galaxy using the SED-fitting code FAST (Kriek et al. 2009). We use stellar population synthesis models from Bruzual & Charlot (2003), and assume a Chabrier (2003) IMF, solar metallicity, and the Calzetti et al. (2000) dust law. The star formation history is parametrised as SFR  ∝  e^−t/τ, swhere the timescale τ is allowed to range between 10 Myr and 10 Gyr. These settings are identical to those used to measure stellar masses of the UltraVISTA sample, which we use to provide a field comparison. For an appropriate analysis that relies on a statistical subtraction of galaxies in the clusters’ fore- and background, as described immediately below, we fix the redshift of each individual galaxy to the clusters’ mean redshift (cf. Table 1).

We also construct tailored catalogues from the UltraVISTA main catalogue for each cluster, where we select only eight filters: uBVrizJK_s, and add artificial noise to the aperture flux measurements to match the depth of the data in the cluster fields. We verify that the performance of EAZY (scatter in photo-zs vs. spec-zs) is then similar to that of the cluster fields for similarly bright galaxies. To perform the statistical field subtraction of fore- and background galaxies, we run FAST on the tailored catalogues to estimate stellar masses, while also fixing the redshift of each galaxy to the cluster’s mean redshift. Identical settings are used when running FAST on the cluster fields, as on the reference field.

Cluster galaxies are initially selected to have a photometric redshift within 0.07 from the cluster mean spectroscopic redshift, which is several times larger than the photo-z scatter. In Appendix A.1, we test the effect of this choice, and show that the reported results are robust.

The uncertainty due to cosmic variance in the reference field is taken into account following the prescription of Moster et al. (2011), based on the volume subtended by the UltraVISTA area in the redshift range 0.5 <  z <  0.7. We find that this has no significant effect on the cluster SMF, except in the outskirts, where the overdensity of the cluster field with respect to the background decreases. The estimated uncertainties are shown in the main figures presented in this work.

3.4. Completeness correction

We identify and correct for two observational effects that affect galaxies around the detection limit in the K_s-band stack. First, the detection rate of objects drops towards fainter magnitudes due to noise fluctuations. Second, the objects that are detected may have a flux measurement that is biased compared to their intrinsic brightness.

To measure the influence of these effects, we study the recovered fraction and fluxes of simulated sources in the detection band. We inject sources with an exponential (i.e. Sérsic-n = 1) profile and half-light radii in the range 1–3 kpc (uniform distribution), ellipticities uniformly drawn between 0.0 and 0.2, and a uniform magnitude distribution. These values are appropriate for sources around our detection limit. Since the depth of the detection image is not uniform, we consider the region within 6′ radial distance of the cluster centres. This corresponds to 2.2 Mpc (2.6 Mpc) at z  =  0.50 (z  =  0.70) and covers the parts of the clusters that are relevant for this study. On one occasion we probe the properties of galaxies up to larger cluster-centric distances; in Appendix A.3 we study the effect of the reduced depth in the cluster outskirts on this result.

We run exactly the same source detection algorithm as for the main analysis (cf. Sect. 2.2) on the K_s-band stacks that include the simulated sources. The K_s-band magnitude limits at which 80% of the simulated sources are recovered, are listed in Table 2 and also indicated in the right hand panel of Fig. 1. We find that around this limit, sources are measured to be 0.12 magnitudes fainter than they are intrinsically, while for the brightest sources this difference is negligible (∼0.01). We correct the measured fluxes by these magnitude-dependent corrections.

Stellar mass limits that correspond to the 80% completeness limit in the K_s-band are also listed in Table 2. We base these on a Bruzual & Charlot (2003) template with a formation redshift of z = 3.0. Since such an old stellar population is relatively faint for their stellar mass, this corresponds to a conservative limit when galaxies with more recent star formation are also considered.

3.5. Star-forming vs. quiescent galaxies

We measure rest-frame U−V and V−J colours of the best-fit SEDs from EAZY, while fixing the redshift to the cluster mean redshift. These colours have been shown to be effective to separate star-forming from quiescent galaxies, even in the presence of dust reddening (e.g. Wuyts et al. 2007; Williams et al. 2009; Patel et al. 2012).

There are small offsets in the UVJ colour distribution between the different clusters, and with respect to the UltraVISTA field sample. Such differences are not uncommon, and several studies have applied corrections to the selection criteria to rectify this (Whitaker et al. 2011; Muzzin et al. 2013a; Skelton et al. 2014, also see the discussion in Appendix A of Lee-Brown et al. (2017)). Rather than adapting the selection criteria of quiescent galaxies for each cluster, we applied offsets to the rest-frame colours to shift them back to the UltraVISTA reference in redshift range 0.5 <  z <  0.7. The mean absolute shift applied is 0.06 magnitude to the U−V colour, and 0.04 to the V−J colour. This brings the colour distribution of cluster galaxies into good agreement with those of field galaxies in UltraVISTA; see Fig. 2. In Appendix A.2 we test the effect of a possible residual shift in colour between the cluster fields and the reference field to the main results presented in this paper.

Inspecting the bimodal galaxy distribution by eye, we select a sample of quiescent galaxies following the following criteria:$$U-V>1.3\cap V-J<1.6\cap U-V>0.55+(V-J).$$(2)

We note that Fig. 2 shows the UVJ colour distribution of the UltraVISTA galaxies based on all the photometric information. For the purpose of a statistical background subtraction we run a similar analysis on a subset of eight UltraVISTA filters, with noise added to resemble the photometric quality of the cluster fields (cf. Sect. 3.3).

Fig. 2. Rest-frame U−V vs. V−J diagram for galaxies with stellar masses M⋆  ≥  1010 M⊙. Grey distribution: UltraVISTA field galaxies with redshifts 0.5 <  z <  0.7. Red points: galaxies from the cluster sample studied here, within redshift |Δz|  ≤  0.05 from the cluster mean redshift, and within projected R  ≤  R500 from the cluster centres.

4. The stellar mass function

The galaxy SMF is a fundamental measure of any population of galaxies, and a critical measurement against which galaxy-formation models are tested. The SMF is measured for the cluster galaxies up to a cluster-centric radius of 2  ×  R₅₀₀ and shown in the left panel of Fig. 3. For this measurement all 21 clusters are stacked, excluding the BCGs. The background is subtracted and an incompleteness correction is applied. For each bin we only use clusters down to their 80% stellar mass completeness limit (cf. Table 2). To compensate for clusters missing in the lowest-mass bins, we increase the weight of galaxies in the remaining clusters. To weigh each cluster properly, we base these on the richnesses measured for each cluster. Richnesses are measured following the definition given in Rykoff et al. (2014), and we discuss this mass proxy in van der Burg et al. (in prep.).

We measure the “average” cluster all the way down to the stellar mass limit of 10^9.5 M_⊙. Error bars denote uncertainties estimated from 100 bootstrap re-samplings of all galaxies in which we draw galaxies with replacement. The number of galaxies we draw in each realisation follow a Poisson distribution with mean equal to the number of galaxies in the stack. To make sure that the uncertainties are not dominated by the (perhaps too low) number of clusters in our sample, we perform 25 additional re-samplings of the clusters themselves. We find that both bootstrap procedures result in comparable uncertainties; the sample of clusters is therefore large enough that we would have obtained the same results as presented here, if we would have observed 21 different clusters taken from a similar parent sample.

From the separation between quiescent and star-forming galaxies in the left panel of Fig. 3, we find that the galaxy population in these massive clusters is completely dominated by quiescent galaxies, all the way down to the stellar mass limit (10^9.5 M_⊙). To further interpret our findings of the cluster SMF, we make a comparison with the SMF of field galaxies at the same redshift as the clusters, in the right-hand panel of Fig. 3. From a comparison of cluster and field, it is clear that the fraction of quenched galaxies in the clusters is substantially higher than in the field, at each stellar mass. This is quantified further in the lower panels, where the relative fractions of star-forming and quiescent galaxies are presented as a function of stellar mass.

We model the SMF by fitting a Schechter (Schechter 1976) function to the data. This function is parameterized as$$\begin{array}{c}\hfill \mathrm{\Phi }(M)=ln(10){\mathrm{\Phi }}^{\ast }\left[{10}^{(M-M^{\ast })(1+\alpha )}\right]exp[-{10}^{(M-M^{\ast })}],\end{array}$$(3)

with M^* being the characteristic mass, α the low-mass slope, and Φ^* the overall normalisation. We follow a maximum-likelihood approach to estimate the parameters that define the shape of the Schechter functions, M^* and α, along with their uncertainties. For this we use the un-binned data points, and include the completeness correction to the individual galaxies, by setting their weights in the likelihood maximisation. For low stellar masses, these weights also compensate for clusters that are not complete down to these limits. For this purpose, each cluster is scaled by its richness. The background galaxies, from the reference field, are included in the same likelihood and have a negative weight.

The normalisation of the Schechter function, Φ^*, is evaluated by requiring that the integral over the considered stellar mass range (i.e. stellar masses larger than 10^9.5 M_⊙) equals the number of all cluster galaxies (or more specifically, the sum of all weights). The best-fitting Schechter parameters are listed in Table 4, and the corresponding functions are shown in the top panels of Fig. 3. The reported goodness-of-fit values indicate that the Schechter functions provide reasonable fits to the data points. However, there seem to be some systematic residuals, especially towards the low-mass end of the SMF. Indeed, some literature studies fit double Schechter functions, but since we primarily work with the data points from now on, this paper does not discuss whether a fit can be improved with more degrees of freedom. The data points themselves are listed in Table 3.

Fig. 3. Top row, left panel: SMF of cluster galaxies within R  ≤  2  ×  R500 from the cluster centres. Blue and red data points: the population of star-forming and quiescent galaxies, respectively. Black points: Full galaxy population. Top row, right panel: SMF of field galaxies from the COSMOS/UltraVISTA survey at the same redshifts (0.5 <  z <  0.7). Lower panels: relative fraction of quiescent galaxies and star-forming galaxies as a function of stellar mass.

4.1. Radial-dependence of SMF

A study of the SMF of galaxies at different radial distances from the cluster centres would allow a more detailed understanding of what happens to the galaxies as they are accreted by the clusters. With the current sample, we have the statistics to do this, and Fig. 4 shows the former cluster SMF, split into four radial bins. The best-fitting Schechter functions are overplotted. Qualitatively, a strong trend is immediately visible; the quenched fraction of galaxies drops with radial distance from the cluster centre, at each stellar mass.

Realising that the uncertainties on the best-fitting Schechter parameter are not independent of each other, we plot the two-dimensional (2D) 68% and 95% confidence regions on α and M^* in Fig. 5. We note that for star-forming galaxies the overdensity compared to the reference field is low, especially in the outskirts. This results in large uncertainties on their SMF, and on the fitted Schechter parameters. The cosmic variance uncertainty is visualised by the lines that are superimposed on the ellipses; these connect the ±1σ systematic uncertainties due to cosmic variance with the nominal measurement. These systematics are so large that the shape of the SMF of star-forming galaxies is consistent with being independent of environment.

The quiescent galaxies have a much higher overdensity compared to the reference field, and their SMF can therefore be measured more accurately since cosmic variance plays a negligible role. The shape of the SMF of quiescent galaxies does not vary significantly with the radial distance from the cluster centre. Also the shape of the total SMF, which is always completely dominated by quiescent galaxies, does not vary significantly between radial bins. There is a significant difference, however, between the SMF of quiescent galaxies in the clusters and in the average field; there are relatively more low-mass quiescent galaxies in the clusters. Some quenching models (e.g. Peng et al. 2010) interpret this as indicative of another quenching mechanism of galaxies in the field as well as those in clusters. If galaxies are “environmentally” quenched, one would expect a steeper SMF of quenched galaxies at low masses. That is because the SMF of star-forming galaxies is steep, and because quenching due to environment is supposed to be largely mass-independent. Our findings are broadly consistent with that picture.

Fig. 4. Top panels: SMF of cluster galaxies at different radial distances from the cluster centres. Blue and red data points: the population of star-forming and quiescent galaxies, respectively. Black points: full galaxy population. Lower panels: relative fraction of quiescent galaxies and star-forming galaxies as a function of stellar mass.

Fig. 5. 68% and 95% likelihood contours for the Schechter parameters α and M*, for the different galaxy types in different panels. Black: schechter parameters of field galaxy population. Other colours: schechter parameters of the cluster galaxies, in different radial bins, as indicated. The lines within the contours indicate how cosmic variance (+/−1σ) affects the best-fit Schechter parameters. That effect is strongest for the star-forming galaxies at large clustercentric distances, where the overdensity compared to the background is lowest.

4.2. Normalisation of the SMF

The cluster SMFs presented in Figs. 3 and 4 are normalised in number of galaxies per cluster. To understand and compare the efficiency of galaxy formation in clusters and in the field, we make a more direct comparison in the normalisations between the different environments. Since clusters have, by definition, a very high volume density of galaxies compared to the field, normalising the SMF over the total volume may not be insightful. Instead, following van der Burg et al. (2013), we normalise the SMF of field and clusters to the total amount of matter associated with each galaxy population in Fig. 6.

For the field sample, we take the total comoving volume in the redshift range 0.5 <  z <  0.7 covered by the UltraVISTA unmasked survey area of 1.62 square degrees (Muzzin et al. 2013b), amounting to 1.5  ×  10⁶ Mpc³. Multiplying this by the average matter density of the Universe in our cosmology, we find a total amount of matter (i.e. dark matter + baryonic) of 6.1  ×  10¹⁶ M_⊙.

For the clusters, we take the total mass within a projected radius of R₅₀₀, but integrated along the line-of-sight. For this we assume that the galaxies follow the total mass, approximated by an NFW profile with concentration c₅₀₀  ≈  2 − 3 (e.g. Dutton & Macciò 2014; Klypin et al. 2016, van der Burg et al., in prep.). To go from the mass within a sphere of radius R₅₀₀ (which equals M₅₀₀, by definition) to the mass within this cylinder, one multiplies M₅₀₀ with a constant factor of ∼1.43. The total mass associated with the 21 clusters, within a projected radius of R₅₀₀ is then 2.0  ×  10¹⁶ M_⊙.

Figure 6 shows the resulting SMF, normalised by the total mass. Per unit total mass there is a clear overdensity of galaxies in the clusters compared to the field (qualitatively similar to what was found by van der Burg et al. 2013). This shows that it is impossible to create clusters by simply accreting an average field population of galaxies, since the latter includes low-density regions such as voids, where the star-formation efficiency is very low. Since groups are typically found in the vicinity of clusters, it is likely that the accretion of these systems caused a galaxy excess in clusters compared to the field (cf. discussions in Hennig et al. 2017; Chiu et al. 2018).

Fig. 6. SMF of galaxies in clusters and field, normalised by the total mass associated with the cluster population (black) and the total mass in the volume from which the field SMF is measured (purple). Squares: SMF including BCGs, to provide a full accounting of stellar mass in the cluster galaxies.

4.3. Literature comparison

The general trends we observe are in line with previous measurements. Vulcani et al. (2013) study the SMF of cluster galaxies in the same redshift range as we do. For each galaxy type, they find no significant difference in the shape of their SMF between cluster and field. Contrary to their work, however, we do find a significant difference between the SMF of quiescent cluster and field galaxies. It is possible that this is owing to our data being substantially deeper (∼1 dex), allowing us to probe the low-mass regions where the differences are most pronounced.

Papovich et al. (2018) also study the SMF over a redshift range that overlaps with ours, in different environments probed with ZFOURGE and in the NMBS. They have data with similar depth, or even slightly deeper than ours, and base their density estimates on a nearest-neighbour approach. Similar to this work, they find a steepening in the low-mass slope of the SMF of quiescent galaxies in overdensities compared to the field. However, Papovich et al. (2018) also find an increase in M^* towards higher densities (also see Davidzon et al. 2016, who base their study on the VIPERS data set), which we do not find in the cluster environments. A possible explanation for the apparent discrepancy with these studies is that they probe more moderate overdensities, and include galaxies that are central to their own haloes, while we have deliberately not taken the BCGs into account. Central galaxies are expected to grow from the merging of in-falling satellites, and this may explain the increase in M^*. Also, in the more moderate environments probed by ZFOURGE and VIPERS, mergers between satellites may be more frequent than in the cluster environment, where relative velocities are expected to be too high for mergers to occur.

Annunziatella et al. (2014) study the SMF of a massive CLASH cluster at z  =  0.44. They find that the galaxy population is dominated by quiescent galaxies, but only for stellar masses M_⋆  ≳  10¹⁰ M_⊙. This apparent lack of cluster quiescent galaxies may be the result of a different way of separating star-forming from quiescent galaxies, compared to our method. Also cluster-to-cluster variations may play a role here, as suggested by a similar study of the more local cluster Abell 209 (Annunziatella et al. 2016).

At higher redshift, z  ∼  1, van der Burg et al. (2013) measure the SMF of ten clusters, and find that these systems are already dominated by quiescent galaxies down to stellar masses of ∼10¹⁰ M_⊙. Contrarily to more local studies, they find that the shape of the SMF of quiescent galaxies is similar between clusters and the field. A likely explanation, apart from a possible evolution with redshift, is that a stellar mass limit of 10¹⁰ M_⊙ may not be low enough to probe any differences in the SMF where they are expected (cf. Nantais et al. 2016, for a study at even higher redshift).

5. Environmental quenching efficiency

A significant result of this work is that the fraction of quiescent galaxies is much higher in the cluster environment than in the field, at the same redshift, and at a given stellar mass. Here we quantify this using the environmental quenching efficiency (f_EQ), which can be thought of as the fraction of galaxies that would normally be star-forming in the field, but are quenched by their environment. Specifically,$$\begin{array}{c}\hfill f_{\mathrm{EQ}}=\frac{f_{\text{q,cluster}}-f_{\text{q,field}}}{1-f_{\text{q,field}}},\end{array}$$(4)

where f_q, cluster and f_q, field are the quiescent fraction of galaxies in the cluster and field environment, respectively. The quenched fractions are a function of both stellar mass and environment, and therefore also f_EQ may depend on stellar mass and environment. Such a term was already introduced by van den Bosch et al. (2008), and sometimes the term “conversion fraction” is used for a similarly defined quantity (Balogh et al. 2016; Fossati et al. 2017).

A compilation of f_EQ in groups and clusters is shown in Fig. 7 of Nantais et al. (2016). The general trend is that f_EQ increases with increasing halo mass (broadly speaking, from group to cluster environments). There is also a hint that f_EQ increases towards the local universe, at a given halo mass. This is in line with results from Darvish et al. (2016), for example, who find no evidence for environmental quenching in more moderate overdensities in the COSMOS field at redshift z  ≳  1.

5.1. Environmental quenching efficiency versus stellar mass

Figure 7 shows the f_EQ as a function of stellar mass, in four radial bins from the cluster centres. Error bars are the 68% confidence regions from 100 bootstrap resamplings, where galaxies within the clusters are drawn with replacement. In addition, we perform 25 bootstrap resamplings in which the clusters themselves are drawn with replacement. The two bootstrap runs lead to similar uncertainties. We also show that cosmic variance uncertainty (indicated by the shaded regions) increases towards larger clustercentric radii, where the cluster over-density is lower.

In none of the radial bins does the f_EQ show a systematic trend with stellar mass (either increasing or decreasing). However, there are significant wiggles around the mean values (indicated by the dotted line); in Sect. 6 we quantify these as signatures of merging of cluster galaxies compared to the field, which likely happened in their pre-processing environment.

A lack of clear stellar-mass dependence of the f_EQ in any of the radial bins suggests that the main quenching process happening in clusters is mass-independent (also see e.g. Kawinwanichakij et al. 2017). One candidate would be a stripping process of the star-forming gas that is so rapid that it is essentially mass-independent, namely ram pressure stripping.

Fig. 7. Environmental efficiency as a function of stellar mass, for different radial bins. While there are wiggles around a flat relation, there is no significant monotonic mass-dependent trend.

5.2. Environmental quenching efficiency versus radius

The f_EQ being independent of stellar mass, we study its dependence on clustercentric radius. We combine the f_EQ of all galaxies with stellar masses in the range 10^9.5  ≤  M_⋆/M_⊙  ≤  10¹¹ in Fig. 8. Using logarithmic bins, we study the f_EQ from 0.01  ×  R₅₀₀ to 4  ×  R₅₀₀. There is a clear signature of environmental quenching that depends on radial distance (which scales with local density). The f_EQ in the cluster centres (where R/R₅₀₀ <  0.1) are extremely high, ∼90%, even in projection where galaxies on the cluster periphery are mixed along the line-of-sight. There is a steep drop in the f_EQ with radial distance, especially in the range 0.2 <  R/R₅₀₀ <  2.0.

Interestingly, the f_EQ does not drop all the way to zero towards the cluster periphery, but rather converges to a value of ∼0.35 at the largest clustercentric distances that we probe. In Appendix A.3 we demonstrate that this measurement is robust, even though we study galaxies detected in a slightly shallower region of the K_s-band stack. Using a cosmological N-body simulation, Wetzel et al. (2014) show that this excess quenching at large radii may be the result of ejected cluster satellites, which orbit even beyond the clusters’ virial radii. Another possible contributor to this observed excess quenching on larger scales comes from galaxies that have been pre-processed in the rich group environment that surrounds galaxy clusters (e.g. Haines et al. 2015; Bianconi et al. 2018). If we re-define the environmental quenching efficiency in Eq. (4) with respect to the quenched fraction in the cluster periphery, i.e. substitute f_q, field with f_{q, periphery}, we obtain the dashed curve shown in Fig. 8. This curve is based on a pre-processed value of 0.35, and illustrates the effect of the main quenching mechanism in the cluster.

At first glance, the measured strong dependence of f_EQ on radius suggest that, whatever physical process is responsible, quenching must happen on a reasonably rapid timescale, at least when galaxies approach the cluster centres. If quenching were a slow process, freshly accreted star-forming galaxies would have time to migrate to the cluster centres while still forming stars, and this would lower the observed f_EQ in the cluster centres. We quantify these statements in the following subsection, in which we employ a simple quenching model to put the observations into context.

Fig. 8. Environmental quenching efficiency as a function of radial distance from the cluster centres, measured in a single broad mass bin. We can do this because environmental quenching, at least in this regime, seems to be a process that is largely mass-independent (cf. Sect. 5.1). Dashed curve: same but increasing the quenched fraction in the field (cf. Eq. (4)) so that fEQ is consistent with zero in the outermost bin.

5.3. A simple quenching model

We consider a model to identify the approximate timescale over which a galaxy is environmentally quenched in the cluster, and the location where this quenching process is triggered. Our basis is a set of N-body simulations of four galaxy clusters, introduced in Taranu et al. (2014). The four most massive clusters were identified from a large cosmological N-body simulation with 256³ particles in a cosmological box of side length 512 h⁻¹ Mpc. Particles in the re-simulation have masses of 6.16  ×  10⁸ M_⊙, meaning that sub-haloes down to relatively low halo mass can be resolved and traced in time from z  =  3 to z  =  0.

Using this simulation, we investigate at which distances from the cluster centres a quenching transformation process is likely to start, and how long it would take for a galaxy to show a signature of quenched star formation. Following a similar approach to that in Muzzin et al. (2014), in which phase-space distribution of specific subhaloes in these simulations were tracked, we now only consider the projected clustercentric distances of a population of subhaloes in the simulation. Subhaloes are marked that have passed, for the first time, a clustercentric distance r_3D, quench/R₅₀₀ at least T_quench Gyr ago. Projecting each cluster in three directions (x, y, z), we mark the fraction of subhaloes that satisfy these criteria, as a function of projected clustercentric radius.

The results are in Figs. 9 and 10, where one parameter in the model is kept constant, while the other is varied. We note that T_quench has to be interpreted as a delay time + quenching time, and that the quenching time itself is supposed to be a rapid process due to the absence of a significant fraction of green valley (transition) galaxies (Peng et al. 2010; Wetzel et al. 2013). The similarity between Figs. 9 and 10 indicates a degeneracy between quenching radius and timescale.

While any ejected satellites that have been quenched by the same (cluster-specific) mechanism would show up in the projected distributions, we note that any possible pre-processing of galaxies in the large-scale overdensity surrounding the clusters is by definition not shown in this simplistic model. We perform a maximum likelihood comparison between the results from the simulation and the actual data, as a function of quenching timescale and quenching radius, while we include the pre-processed fraction of satellites separately as an extra free parameter. For each value of r_3D, quench/R₅₀₀ and T_quench in the model, we therefore marginalise over this pre-processed fraction of galaxies. As an example, the data points for a pre-processed fraction of 0.35 are shown in Fig. 8. The resulting confidence regions on the parameters that represent the best model are shown in Fig. 11 and these confirm the strong degeneracy between r_3D, quench/R₅₀₀ and T_quench.

The contours shown in Fig. 11 follow an intriguing degeneracy; quenching appears to happen either on a short timescale at small clustercentric radii, or over a longer timescale at larger clustercentric radii. Both scenarios reproduce the data similarly well. We note that the bimodality in the 1−σ contour is likely due to the limited number of four clusters that we used in the simulation (albeit as studied from three orthogonal sight lines). The pre-processed fractions of galaxies are not shown in Fig. 11, and this fraction gradually decreases when going from small radii (and corresponding short times), to larger radii (and longer times); for example, the pre-processing level at r_3D, quench  =  1  ×  R₅₀₀ and T_quench = 1 Gyr is $31.7_{-2.4}^{+2.1}$, while it decreases to $19.8_{-2.9}^{+3.2}$ at r_3D, quench  =  7  ×  R₅₀₀ and T_quench = 3 Gyr.

Even though there is a strong degeneracy between the parameters, we can put a firm lower limit to the radial distance from the cluster centre where quenching is triggered; r_quench >  0.67 R₅₀₀ (95%CL). Using our deep XMM-Newton observations we can link this to an ICM density where quenching occurs. Leaving a more detailed study of the radial gas density distribution of this cluster sample to van der Burg et al. (in prep.), the gas density at R  =  0.7 R₅₀₀ is around ρ_ICM  ≈  2  ×  10⁴ M_⊙ kpc⁻³.

Following their Eq. (62), Gunn & Gott (1972) estimate at which ICM density ram pressure stripping is expected to become effective in the stripping of the interstellar material of a typical infalling galaxy. They find this to happen at an ICM density of ∼5  ×  10⁻⁴ atoms cm⁻³. The density at R  =  0.7 R₅₀₀ is already about twice this value, making ram pressure stripping a likely contributor given our constraints on the location and timescale of the main stripping process. Similarly, based on a cosmological simulation, Zinger et al. (2018) find that a significant removal of star-forming gas happens at r  ≲  0.5 R_vir, which is a similar fraction of R₅₀₀.

Another likely contributor to the quenching process is strangulation/starvation, which is a cut-off from the cosmological accretion of hot gas after a galaxy is accreted by the cluster’s main halo. Star formation is then expected to quench after the reservoir of molecular gas is depleted. This timescale is on the order of 1 Gyr for local low-mass galaxies (Tacconi et al. 2018), and likely shorter at higher redshift. This is a slower process than ram pressure stripping, but its timescale is also consistent with our quenching constraints in Fig. 11. Since this process will be triggered at larger distance (or equivalently, at earlier times) than ram pressure stripping, both processes may be contributing to the observed elevated fraction of quenched galaxies.

Fig. 9. Radial dependence of fEQ from the model, where we vary the quenching time (as indicated) while fixing the quenching location to R500.

Fig. 10. Radial dependence of fEQ from the model, where we vary the quenching location (as indicated) while fixing the quenching time to 1 Gyr.

Fig. 11. 1-, 2-, and 3-σ contours surrounding the best-fitting location and timescale of quenching, so that our quenching model reproduces the observed radial-dependence of the environmental quenching efficiency. Pre-processing is left as a free parameter in the likelihood maximisation.

5.4. Caveats

Our modelling has shown that the radial distribution of f_EQ combined with a realistic simulation of the orbits of cluster galaxies can provide meaningful constraints on the r_3D, quench/R₅₀₀ and T_quench of cluster galaxies. There are several caveats in our analysis that could play a role in the interpretation, and these are highlighted below.

In this paper we measure the f_EQ as the excess quenching with respect to the field, i.e. a representative section of the universe at the same redshift. The field sample therefore already includes regions of higher densities, such as groups. Several studies have quantified the f_EQ (or conversion fraction) as the excess quenching with respect to the lowest mass density found in the field (typically the lowest-density quartile, cf. Papovich et al. 2018). The quenched fractions we report for the field include some “pre-quenching” in moderate overdensities, and the reported f_EQ values would have been even higher if defined with respect to the lowest-density regions. However, since our modelling marginalises over this pre-processing component, this assumption does not affect our estimated parameters that describe the main quenching process in the clusters.

We have defined f_EQ with respect to the field at the same epoch as the clusters are observed at. An alternative approach is to consider the field quenched fraction at the time of galaxy accretion (see the discussion in Balogh et al. 2016), as for example Foltz et al. (2018) have done in their modelling of the quenching time scales. Both approaches have their uses; we have chosen the former so that we can isolate what happens in clusters, separately from and in addition to what would have happened to the star-forming properties of the galaxies if they had remained central in their own haloes.

Furthermore, our definition of f_EQ is interpreted in the absence of mergers between galaxies. We know mergers are happening, especially for galaxies in group scales, and mergers may even contribute to the quenching of satellites (e.g. Peng et al. 2010; Darvish et al. 2016). Mergers may also have affected the observed f_EQ dependence of clustercentric radius. Moreover, mergers lead eventually to growth of the central BCGs. In the following section we discuss these limitations and present a “transformation function” that describes how additional processes such as mergers are affecting the galaxy population to lead to the SMF we observe for cluster galaxies.

6. Discussion

The previous section explored the stellar-mass dependence, and the radial dependence of f_EQ to infer a simple quenching scenario for which we have constrained the approximate location and time scale. In this section we study the impact of a basic environmental-quenching scenario on the SMF, which we measured and studied in Sect. 4.

6.1. A simple quenching model

The simplest environmental quenching model that works in a mass-independent fashion (as suggested by Fig. 7) leaves the shape of the SMF of star-forming galaxies independent of environment (Peng et al. 2010). This is consistent with our findings, considering the measurement uncertainties associated with measuring the SMF of star-forming galaxies in the cluster periphery (Fig. 5). The population of quenched cluster galaxies is then a combination of galaxies that would also have quenched outside of the cluster (mass-quenched galaxies), and the environmentally quenched galaxies, which in principle follow a similar mass distribution as the star-forming galaxies. We start by employing this quenching model in its basic form, following van der Burg et al. (2013).

Figure 12 shows the SMF of cluster galaxies within R  ≤  R₅₀₀, that is, a combination of the first two panels of Fig. 4, or a differently normalised version of Fig. 6. The plotted curves are not fits to the plotted data, but rather adapted versions (only in normalisation) of the best-fitting Schechter functions to the field data from the UltraVISTA survey. The total normalisation of the (red+blue) curve is set by the total stellar mass in cluster galaxies (all data points, including the BCGs). Since we assume that the environmental quenching process affects the star-forming population in a mass-independent manner, we set this by parameter f_EQ, so that the relative normalisation of the star-forming SMF is (1 − f_EQ) compared to the one in the field. Again requiring that the normalisation (total stellar mass) of the blue Schechter function is the same as that of the blue cluster galaxies (data points), we find that f_EQ  =  0.80. The quenched part of the star-forming population gives the “environmentally quenched galaxies”, and when these are combined with the “mass quenched galaxies”, we arrive at the total population of quenched cluster members, as shown by the solid red curve in Fig. 12.

In its basic form, this quenching model over-predicts the abundance of quiescent low-mass galaxies with M_⋆  ≲  10^10.2 M_⊙. Interestingly, a similar trend is revealed in Fig. 10 of van der Burg et al. (2013), in which clusters at slightly higher redshift were studied. However, contrary to this earlier work, we now have the statistics to explore this regime in more detail. The fact that low-mass galaxies show a deficit in clusters compared to this simple quenching model may be related to their destruction, potentially leading to a build-up of the intra-cluster light (ICL). Given the negative colour gradients that are observed in the ICL of massive clusters, dwarf galaxies are likely contributors to the ICL at large clustercentric distances (DeMaio et al. 2015).

Furthermore, milder interactions and mergers between galaxies (likely also in their pre-processing environment, cf. Tomczak et al. 2017) may also affect the SMF of galaxies, and cause it to diverge from the field. For instance, Rudnick et al. (2012) invoke a model that includes mergers between galaxies to reproduce the luminosity function of clusters in the local universe (starting from a distant cluster at z  =  1.62). We attempt to encompass all these processes in a “transformation function”, shown in the lower panel of Fig. 12. Plotted is the difference between the data points for the quiescent galaxies, and the simple quenching model (solid red line). The plot is normalised in stellar mass per bin and per cluster, and highlights again the relative destruction of low-mass galaxies, in favour of the growth of more massive galaxies such as the central BCGs. Since this plot adds up to zero, by construction, it does not include the build-up of ICL.

Fig. 12. Top panel: cluster data points for the different galaxy types within R500. The curves show a basic quenching model inspired by the formalism presented in Peng et al. (2010), with the UltraVISTA SMF in the redshift range 0.5 <  z <  0.7 g point. Bottom panel: difference between the data (total distribution of quenched galaxies) in the cluster compared to the prediction from the model. This suggests a relative removal/destruction of low-mass galaxies and the expected build-up of BCGs.

6.2. How to build a massive cluster of galaxies?

Given our study of 21 massive clusters at 0.5 <  z <  0.7, we summarise some of the steps required to assemble the galaxy population observed within these systems, as opposed to the general field:

• As shown in Fig. 6, where we have normalised the SMFs of clusters and field with respect to the total amount of mass associated with the respective galaxy populations, galaxies form relatively efficiently in (future) cluster environments compared to the average Universe. The galaxy abundance in the clusters we study, per unit total mass, is about twice average.

• At the same redshift, and for each stellar mass we probe, the quenched fraction is much higher in these clusters than in the field; it is also elevated compared to the pre-processing we find happening in the cluster surroundings. This quenching process happens in a largely mass-independent fashion (cf. Fig. 7).

• There is a significant and strong radial trend in the quenched fraction of cluster galaxies, which we describe as the environmental quenching efficiency f_EQ (cf. Fig. 8). A comparison with a model that is based on orbits taken from an N-body simulation suggests that the quenching process likely involves strangulation/starvation after cut-off from cosmological accretion, and ram pressure stripping at smaller clustercentric distances to “finish the job”. Each of these processes ought to happen on timescales that are roughly consistent with what we find in the model.

• Additional transformations are required to reproduce the observed cluster SMF. These are likely largely caused by merger events, the possible destruction of low-mass galaxies in the clusters, and an effective build-up of the BCGs. Mergers are happening likely before galaxies are being accreted into the clusters, since relative velocities in the cluster environments are too high for galaxies to merge there. We have quantified the combination of these effects in the lower panel of Fig. 12.

7. Summary and conclusions

We have studied the galaxy population in a sample of 21 high-mass clusters at 0.5 <  z <  0.7, found in the Planck SZ survey. Using multi-band photometry spanning from the u- to the K_s-band for each cluster, we have defined a sample of cluster galaxies that are highly overdense compared to the back- and foreground.

The data allow for a precise measurement of the galaxy SMF in clusters of intermediate redshift. We have identified differences in the SMF between the cluster population, and galaxies in the field at the same redshift. Normalising the SMF to the total amount of matter associated with each galaxy population, we find that clusters have a higher galaxy content, per unit total mass, than the average field.

The most significant differences between the galaxy population in cluster and field arise when we separate the galaxy population into star-forming and quiescent galaxies by means of their rest-frame U−V and V−J colour distributions. The shape of the SMF of star-forming galaxies does not depend on environment. On the contrary, the SMF of quiescent galaxies is significantly different between the cluster and field; there is a relatively higher fraction of low-mass quiescent galaxies in the clusters. Moreover, the fraction of passive galaxies is much higher in the clusters than in the field, and we quantify how this fraction rises steeply with decreasing cluster-centric radius.

We measured the environmental quenching efficiency (f_EQ), which describes the fraction of galaxies that would be forming stars in the field, but are quenched solely due to their environment. At fixed radial distance from the cluster centre, the f_EQ does not depend on stellar mass. Contrarily, the f_EQ shows a strong radial dependence within the cluster environment.

We interpret the observed radial-dependence of the f_EQ with a simple quenching model based on an N-body simulation. We use this model to constrain the characteristic location and timescale of the main environmental quenching process and find a strong degeneracy between those two parameters. According to the model, quenching may already be triggered at r_3D, quench  ≈  7  ×  R₅₀₀, and would then happen on a long timescale of T_quench  ≈  3 Gyr. If quenching is triggered at shorter radial distances, r_3D, quench  ≈  1  ×  R₅₀₀, it happens on roughly the molecular gas depletion timescale, T_quench  ≈  1 Gyr. Part of the observed quenching excess is therefore likely due to “starvation”/“strangulation” of galaxies using up their cold gas supply. Interestingly, the model rules out a quenching location r_quench <  0.67 R₅₀₀ at 95% confidence. Our XMM-Newton data show that the gas density at this clustercentric distance is so large that ram pressure stripping should be effective, and is likely responsible for the satellite quenching there. This process may therefore “finish the job” whenever the starvation mechanism does not operate rapidly enough.

Acknowledgments

We thank Andrea Biviano and Gabriella De Lucia for insightful discussions, and Michael Balogh for constructive feedback on the manuscript. The research leading to these results has received funding from the European Research Council under the European Union’s Seventh Framework Programme (FP7/2007-2013) / ERC grant agreement n^∘ 340519. HD acknowledges financial support from the Research Council of Norway. Based on observations obtained with MegaPrime/MegaCam, a joint project of CFHT and CEA/DAPNIA, at the Canada-France-Hawaii Telescope (CFHT) which is operated by the National Research Council (NRC) of Canada, the Institute National des Sciences de l’Univers of the Centre National de la Recherche Scientifique of France, and the University of Hawaii. This work is also based on data collected at Subaru Telescope, which is operated by the National Astronomical Observatory of Japan. Based in part on data collected at Subaru Telescope and obtained from the SMOKA, which is operated by the Astronomy Data Center, National Astronomical Observatory of Japan. Based on observations obtained with XMM-Newton, an ESA science mission with instruments and contributions directly funded by ESA Member States and NASA. The development of Planck has been supported by: ESA; CNES and CNRS/INSU-IN2P3-INP (France); ASI, CNR, and INAF (Italy); NASA and DoE (USA); STFC and UKSA (UK); CSIC, MICINN and JA (Spain); Tekes, AoF and CSC (Finland); DLR and MPG (Germany); CSA (Canada); DTU Space (Denmark); SER/SSO (Switzerland); RCN (Norway); SFI (Ireland); FCT/MCTES (Portugal); and PRACE (EU). This research made use of the following databases: the NED and IRSA databases, operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with the NASA; SIMBAD, operated at CDS, Strasbourg, France; SZ-cluster database operated by IDOC at IAS under contract with CNES and CNRS.

References

Appendix A:

Robustness tests

In this section we study the effect of several assumptions we had to make in our analysis on the final results.

A.1. Photo-z selection

One of the choices in our analysis is the initial photo-z selection window of cluster members, before the background subtraction is performed. If such a window is chosen to be too small, uncertainties in photometric redshifts of cluster members may cause them to scatter out of the selection window. On the other hand, if the photo-z selection window is chosen to be larger than necessary, we introduce additional noise while performing the statistical background subtraction. Since we use a single reference field to do the statistical background subtraction, especially the contribution from cosmic variance in that field would increase with a larger photo-z cut (since the overdensity will go down).

From all 1527 spectroscopically confirmed cluster members in our sample, we find that 89.6% satisfy the photo-z cut of 0.07 (i.e. only 10.4% have scattered out). Therefore, to first order, the normalisation of the cluster SMF is higher by 1/0.896 compared to our measurement. The scattering may also bias the shape of the SMF, since the photo-z scatter depends (in principle, though only slightly in practise) on the galaxy stellar mass (cf. Fig. 1). To test the effect of this, we have performed the analysis with photo-z selections of |Δz|< 0.10 and 0.13. This leads to percentages of galaxies that scatter out, of 5.7% and 4.5% respectively, approaching the percentage of catastrophic outliers. We checked that the main results in this paper, i.e. the SMF and the behaviour of the f_EQ , do not change significantly (i.e. by more than the reported uncertainties) when a broader photo-z selection is chosen. The only exception is the SMF of star-forming galaxies in the outskirts (R >  R₅₀₀) of the clusters, where the overdensity is very low. We note, however, that since the overdensity of cluster galaxies with respect to the background drops for a broader photo-z selection, uncertainties of all measurements grow substantially. We have therefore chosen a cut of 0.07 in |Δz|.

A.2. UVJ division

The division between star-forming and quiescent galaxies is a critical part of the analysis. Our analysis makes use of the U−V and V−J rest-frame colours, which ensures that we separate the effect of dust reddening from the reddening due to lack of star formation. In Sect. 3.5 we described the small corrections (on average 0.04–0.06) we made to the rest-frame U−V and V−J colours of the cluster galaxies, to match them to the colour distribution of the reference field. Here we perform a test to check how a residual systematic colour offset between field and cluster data would impact our results. We increase all U−V rest-frame colours of cluster galaxies, and reference background galaxies, by ±0.05, and re-measure the quenched fractions of galaxies in the cluster. Even such a large offset, compared to the residuals we expect, changes quenched fractions of cluster galaxies by at most 10%. When we measure f_EQ , this has a larger effect, since we have not changed the rest-frame colours of the field galaxies. The result on the radial dependence of f_EQ is shown in Fig. A.1. Even though a significant change is notable, we note that this is based on a rather extreme systematic error between rest-frame colours of clusters and field galaxies.

Fig. A.1. Environmental quenching efficiency as a function of radial distance from the cluster centres, measured in a single broad mass bin. This is a similar plot as Fig. 8, but now we explore the effect of a rather extreme residual shift in rest-frame colour between the field and cluster data.

A.3. Radial dependence of f_EQ

The main analysis of this paper studies the properties of galaxies within 2  ×  R₅₀₀ from the cluster centres. These regions fall on the deep part of the K_s-band stacks, and we have characterised the completeness in Sect. 3.4. Figure 8, however, shows the environmental quenching efficiency up to 4  ×  R₅₀₀, which covers part of the shallower regions. Due to the dither strategy we chose for the K_s-band imaging, the depth in the K_s-band drops by a maximum amount of 0.7 magnitudes towards a distance of 4  ×  R₅₀₀ from the cluster centres, corresponding to 0.3 dex in stellar mass. We note that the optical data, which are essential for precise estimates of photometric redshifts, extend to a larger region around the clusters at uniform depth. We studied the impact of a slight decrease in depth of the detection band on the results plotted in Fig. 8. The measurements move within the plotted uncertainties, when galaxies in the range 10^9.8  ≤  M_⋆/M_⊙  ≤  10¹¹ are considered (instead of 10^9.5  ≤  M_⋆/M_⊙  ≤  10¹¹). In particular, we verified that the “plateau” in f_EQ at radii R  ≳  2  ×  R₅₀₀ is robust, and not an effect of this decrease in depth.

Appendix B:

Cluster gallery

Colour-composite images of all 21 clusters are shown in Fig. B.1. They are composed of g- or B-, i- or I_c-, and K_s-band imaging. Regions around bright stars, and their diffraction spikes, are clearly visible here, but these are all masked and not considered in our analysis.

Fig. B.1. Colour composite images of the clusters in our sample, based on g- or B-, i- or Ic-, and Ks-band imaging.

Fig. B.1. continued.

Fig. B.1. continued.

Fig. B.1. continued.

Overplotted are X-ray surface brightness contours from the deep XMM-Newton observations, which we have available for all clusters. They are based on (adaptively-) smoothed surface brightness maps, which are background subtracted, exposure corrected, and from which point sources are excised. Contours are logarithmically spaced with 0.2 dex increments. These data form the basis of the X-ray morphological analysis, which is presented in Arnaud et al. (in prep.).

All Tables

All Figures

Fig. 1. Left panel: spectroscopic vs. photometric redshifts for the 21 cluster fields. Outliers, objects for which |Δz|> 0.15, are marked in red. The outlier fraction is 3.5%, the scatter of the remaining objects is σz  =  0.029. Right panel: scatter in Δz as a function of Ks-band magnitude, for sources with 0.5 <  zspec <  0.7. Quiescent and star-forming galaxies are marked in red and blue, respectively. The source detection limits of the 21 Ks-band stacks are indicated.

In the text

	Fig. 2. Rest-frame U−V vs. V−J diagram for galaxies with stellar masses M⋆  ≥  1010 M⊙. Grey distribution: UltraVISTA field galaxies with redshifts 0.5 <  z <  0.7. Red points: galaxies from the cluster sample studied here, within redshift |Δz|  ≤  0.05 from the cluster mean redshift, and within projected R  ≤  R500 from the cluster centres.
In the text

Fig. 3. Top row, left panel: SMF of cluster galaxies within R  ≤  2  ×  R500 from the cluster centres. Blue and red data points: the population of star-forming and quiescent galaxies, respectively. Black points: Full galaxy population. Top row, right panel: SMF of field galaxies from the COSMOS/UltraVISTA survey at the same redshifts (0.5 <  z <  0.7). Lower panels: relative fraction of quiescent galaxies and star-forming galaxies as a function of stellar mass.

In the text

	Fig. 4. Top panels: SMF of cluster galaxies at different radial distances from the cluster centres. Blue and red data points: the population of star-forming and quiescent galaxies, respectively. Black points: full galaxy population. Lower panels: relative fraction of quiescent galaxies and star-forming galaxies as a function of stellar mass.
In the text

Fig. 5. 68% and 95% likelihood contours for the Schechter parameters α and M*, for the different galaxy types in different panels. Black: schechter parameters of field galaxy population. Other colours: schechter parameters of the cluster galaxies, in different radial bins, as indicated. The lines within the contours indicate how cosmic variance (+/−1σ) affects the best-fit Schechter parameters. That effect is strongest for the star-forming galaxies at large clustercentric distances, where the overdensity compared to the background is lowest.

In the text

	Fig. 6. SMF of galaxies in clusters and field, normalised by the total mass associated with the cluster population (black) and the total mass in the volume from which the field SMF is measured (purple). Squares: SMF including BCGs, to provide a full accounting of stellar mass in the cluster galaxies.
In the text

	Fig. 7. Environmental efficiency as a function of stellar mass, for different radial bins. While there are wiggles around a flat relation, there is no significant monotonic mass-dependent trend.
In the text

Fig. 8. Environmental quenching efficiency as a function of radial distance from the cluster centres, measured in a single broad mass bin. We can do this because environmental quenching, at least in this regime, seems to be a process that is largely mass-independent (cf. Sect. 5.1). Dashed curve: same but increasing the quenched fraction in the field (cf. Eq. (4)) so that fEQ is consistent with zero in the outermost bin.

In the text

	Fig. 9. Radial dependence of fEQ from the model, where we vary the quenching time (as indicated) while fixing the quenching location to R500.
In the text

	Fig. 10. Radial dependence of fEQ from the model, where we vary the quenching location (as indicated) while fixing the quenching time to 1 Gyr.
In the text

	Fig. 11. 1-, 2-, and 3-σ contours surrounding the best-fitting location and timescale of quenching, so that our quenching model reproduces the observed radial-dependence of the environmental quenching efficiency. Pre-processing is left as a free parameter in the likelihood maximisation.
In the text

Fig. 12. Top panel: cluster data points for the different galaxy types within R500. The curves show a basic quenching model inspired by the formalism presented in Peng et al. (2010), with the UltraVISTA SMF in the redshift range 0.5 <  z <  0.7 g point. Bottom panel: difference between the data (total distribution of quenched galaxies) in the cluster compared to the prediction from the model. This suggests a relative removal/destruction of low-mass galaxies and the expected build-up of BCGs.

In the text

	Fig. A.1. Environmental quenching efficiency as a function of radial distance from the cluster centres, measured in a single broad mass bin. This is a similar plot as Fig. 8, but now we explore the effect of a rather extreme residual shift in rest-frame colour between the field and cluster data.
In the text

	Fig. B.1. Colour composite images of the clusters in our sample, based on g- or B-, i- or Ic-, and Ks-band imaging.
In the text

	Fig. B.1. continued.
In the text

	Fig. B.1. continued.
In the text

	Fig. B.1. continued.
In the text