CANDIDATE CLUSTERS OF GALAXIES AT z > 1.3 IDENTIFIED IN THE SPITZER SOUTH POLE TELESCOPE DEEP FIELD SURVEY

The SSDF, a post-cryogenic Spitzer Exploration Science program, obtained 120 s depth observations in the 3.6 μm and 4.5 μm IRAC bandpasses (hereafter, [3.6], [4.5]). Ashby et al. (2013) provide detailed information on the survey design, observations, processing, source extraction, and publicly available data products. Our study is based on the 4.5 μm selected Spitzer/IRAC band-merged ([3.6], [4.5]) catalog, which contains ∼3.7 million distinct sources down to the SSDF 5σ sensitivity limit of 21.46 AB mag (9.4 μJy) at 4.5 μm; the corresponding 5σ sensitivity of the 3.6 μm bandpass is 21.79 AB mag (7.0 μJy). Throughout this paper, we use aperture-corrected, 4'' diameter aperture magnitudes for the IRAC data.

The SuperCOSMOS survey (Hambly et al. 2001) provides all-sky optical photometry based on scans of photographic Schmidt survey plates from the UK Schmidt Telescope (UKST) Southern Surveys (Hartley & Dawe 1981; Cannon 1984) and Palomar Oschin Schmidt Telescope Surveys (POSS; Minkowski & Abell 1963; Reid et al. 1991). These shallow data provide I-band magnitudes down to I ∼ 20.45 mag (AB) in the SSDF (Hambly et al. 2001). As discussed above, portions of the SSDF have deeper optical data from the BCS and, eventually, the entire SSDF will be covered by the DES. However, in the interest of uniformity in our cluster search over the widest possible area, we only consider the relatively shallow SuperCOSMOS optical data in the following analysis.

We select candidate distant galaxy clusters based on their [3.6]–[4.5] galaxy color, following the approach of Papovich (2008); that methodology was proven effective by discovering a z = 1.62 galaxy cluster in the SWIRE-XMM field (Lonsdale et al. 2003) using data of a depth similar to what is available in the SSDF (Papovich et al. 2010). The method takes advantage of the fact that the [3.6]–[4.5] color is a linear function of redshift between 0.7 ≲ z ≲ 1.5 (Figure 1), and thus can be used as an effective redshift indicator. At z ≳ 1.5, the color reaches a plateau out to z ∼ 3. The basis of the selection is that galaxy stellar populations with ages >10 Myr have a prominent bump at ∼1.6 μm due to a minimum in the opacity of the H⁻ ion present in the atmospheres of cool stars. This feature is seen in the spectral energy distribution (SED) of essentially all galaxies, largely independent of star formation history or age (John 1988; Simpson & Eisenhardt 1999; Sawicki 2002; Sorba & Sawicki 2010). However, above redshift z ∼ 1.5 the color is not precise enough to be useful in estimating redshifts other than to constrain the redshift to be larger than z ≳ 1.5.

Zoom In Zoom Out Reset image size

As shown in the middle panel of Figure 1, while an IRAC color cut [3.6]–[4.5]  >   − 0.1 is effective at distinguishing galaxies at z > 1.3, having at least one relatively shallow optical band is very useful for alleviating contamination from foreground interlopers at z ∼ 0.3 (see discussion in Muzzin et al. 2013b). With this aim, we apply magnitude cuts in the I and [4.5] bands to remove most z < 0.4 galaxies (see bottom and top panel of Figure 1). Combined, these cuts effectively remove the bulk of the foreground galaxy population at z < 1.3. The remaining sources consist predominantly of high-redshift galaxies. We note that other sources of contamination are cool brown dwarfs (Stern et al. 2007) and powerful active galactic nuclei at all redshifts (Stern et al. 2005); however, these are expected not to be contaminants (see discussion in Galametz et al. 2012).

Therefore, to discover distant clusters of galaxies at z > 1.3, we have implemented a simple three-filter algorithm to search for overdensities of galaxies based on their Spitzer IRAC color ([3.6]–[4.5] > −0.1), their 4.5μm magnitude ([4.5] > 19.5), and requiring non-detection in the SuperCOSMOS I-band data (I > 20.45). Similar algorithms have been demonstrated to be effective using several programs (Papovich et al. 2010; Galametz et al. 2010, 2013; Gettings et al. 2012; Muzzin et al. 2013b).

To identify candidate galaxy clusters, for each galaxy in the SSDF catalog that meets the aforementioned criteria we count the number of similarly selected companions within 10 radius cells. We correct these counts in cells for completeness using the values published in Table 5 of Ashby et al. (2013), who used a Monte Carlo approach to estimate the survey completeness and sensitivity by placing hundreds of simulated sources with randomly assigned Vega magnitudes (between 10 and 21 mag) in the SSDF mosaics at random locations (accounting realistically for the effects of source confusion). Ashby et al. (2013) then photometered the artificial sources in a manner identical to that used for the original mosaics. The process was iterated, so that a total population of 20,000 simulated sources was ultimately analyzed in each of the 0.5 mag wide bins. This procedure was repeated in both bands. The results are given in their Table 5 and shown in their Figure 5. We also note that at the depths considered in our analysis, the completeness corrections vary, as a function of magnitude, between 1.3 and 11.5%.

When creating an IRAC mosaic, the mapping strategy adopted results in some variations in depth across a large field. Therefore, we search for cells with the highest significance overdensity of IRAC-selected sources above the local background. With this aim, for each selected object in the catalog, we also count the completeness-corrected number of similarly selected companions in an outer annulus defined by radii of 5' and 8'. After taking into account the difference in surface areas, we then derive the excess number of selected companions with respect to the local background, ΔN.

Figure 2 shows the distribution of the completeness-corrected excess number of companions, ΔN, within 10 radius. The mean number of excess galaxies is 〈ΔN〉 = 0.4. The distribution is skewed toward objects with higher-than-average numbers of companions, suggestive of strong clustering. Similar to Papovich (2008), we fit a Gaussian to the distribution, iteratively clipping at 2σ. The best fit is shown in Figure 2 and it matches well with the low-excess half of the distribution. We find the width of the fitted Gaussian to be σ_ΔN = 3.7.

Zoom In Zoom Out Reset image size

The minimum detection significance to adopt, $X_{f_{\rm min}}$, is somewhat arbitrary, trading between completeness and reliability in the derived candidate cluster catalog. We based our conservative choice on extensive tests performed on similar-depth IRAC and SuperCOSMOS data in the Boötes field, where ancillary spectroscopic and much deeper multi-wavelength photometric coverage is also present. We describe these tests in Section 3.2, where we estimate that the threshold adopted here ensures a purity of ∼80%. Thus, we define our sample of candidate clusters of galaxies to be those with $\Delta N \ge \langle \Delta N \rangle + X_{f_{\rm min}} \cdot \sigma _{\Delta N}$, with $X_{f_{\rm min}} = 5.2$. For the SSDF field, this corresponds to objects with ΔN ⩾ 19.5 excess companions with respect to the local background (indicated by the vertical dot–dashed red line of Figure 2).

By design, many of the objects in overdense regions will be identified as companions to multiple sources by our algorithm. To address this, we adopt an approach similar to that in Papovich (2008). We merge cluster candidates by applying a friends-of-friends algorithm with a linking length of 10. This results in a final catalog containing N_obs = 279 candidate clusters; their spatial distribution and color-coded detection significance are shown in Figure 3. We note that none of our clusters is found in the list of SPT SZ-detected clusters reported in Reichardt et al. (2013). This is not unexpected, since the two samples do not overlap in mass–redshift space: the highest-redshift cluster in Reichardt et al. (2013) is at z = 1.075, below the range to which our cluster-finding algorithm is sensitive, while the typical mass of the clusters we find is below the mass threshold of the Reichardt et al. (2013) catalog.

Zoom In Zoom Out Reset image size

In Table 1, we present the list of the top 10 most significant cluster candidates in the SSDF field. Figures 4 and 5 present their images and color–magnitude diagrams. The coordinates listed in Table 1 refer to the galaxy at the center of the cell containing the highest number of cluster member candidates.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

For a large sample of candidate clusters to be useful for galaxy evolution and cosmological studies, it is important to determine its purity f_pure as a function of the detection significance X_f, defined as

$$\begin{equation} f_{{\rm pure}} (X_{f})= \frac{N_{{\rm real}}}{N_{{\rm tot}}} = 1 - \frac{N_{{\rm false}}}{N_{{\rm tot}}}, \end{equation} \tag{ 1 }$$

where N_tot is the total number of cluster candidates above the detection threshold X_f, N_real is the number of candidates corresponding to real clusters, and N_false is the number of false detections. With this aim, we run our cluster-finding algorithm on comparable-depth observations from the IRAC Shallow Cluster Survey (ISCS; Eisenhardt et al. 2004).

The ISCS is a wide-field IR-selected galaxy cluster survey carried out using 90 s Spitzer/IRAC imaging of the 8.5 deg² Boötes field of the NOAO Deep, Wide-Field Survey (NDWFS; Jannuzi & Dey 1999). SuperCOSMOS data are available in Boötes of a depth comparable to the data available in the SSDF field. Over the last decade, the Boötes cluster candidates have been the target of extensive ground- and space-based spectroscopic and photometric campaigns (e.g., Eisenhardt et al. 2008; Ashby et al. 2009; Stern et al. 2010; Brodwin et al. 2011; Stanford et al. 2012; Zeimann et al. 2012; Gonzalez et al. 2012; Brodwin et al. 2013). Using a wavelet search algorithm operating on photometric-redshift probability distribution functions, as described in Brodwin et al. (2006), more than 100 rich cluster candidates at z > 1 were identified. To date, 18 of these have been spectroscopically confirmed out to z = 1.9. Note that very accurate photometric redshift measurements are based on data from the Spitzer Deep, Wide-Field Survey (SDWFS; Ashby et al. 2009, 4× deeper than the IRAC data used here), deep optical imaging in B_W, R, I bands from the NDWFS, and NIR photometry from the FLAMEX survey (Elston et al. 2006) in the J and K_s bands. Note, however, that the ISCS is ongoing, with many cluster candidates still awaiting confirmation. In the absence of a complete, spectroscopically confirmed catalog of high-redshift clusters, we deem the ISCS the best available survey to test our algorithm, as it contains the largest sample (>10) of spectroscopically confirmed clusters at z > 1.3.

We find the purity of our sample, f_pure, to be a monotonic function of X_f reaching f_pure = 0.8 for $X_{f_{\rm min}}=5.2$. At this very high detection significance, our algorithm identifies 14 candidate clusters at z > 1.3 in the Boötes field, of which 3 are spectroscopically confirmed at redshifts of z = 1.37 (Brodwin et al. 2013; Zeimann et al. 2013), z = 1.41 (Stanford et al. 2005; Brodwin et al. 2011), and z = 1.75 (Stanford et al. 2012), and an additional eight have accurate photometric redshifts 1.3 < z_phot < 2.3. In the left panel of Figure 6, we show the [3.6]–[4.5] color versus [3.6] magnitude for the cluster IDCS J1426.5+3508 at z = 1.75 as found by our algorithm with X_f = 5.3. This cluster, spectroscopically confirmed with HST/WFC3 grism observations (Stanford et al. 2012), was detected in both archival 8.3 ks Chandra imaging of the field (Stanford et al. 2012) as well as follow-up SZ observations with the CARMA array (Brodwin et al. 2012). The cluster also has a giant arc in HST imaging, implying it is a lensing cluster; this is particularly surprising given the cluster redshift and the small area of sky surveyed (Gonzalez et al. 2012). At the time of its discovery, IDCS J1426.5+3508 was the highest redshift cluster for which the SZ effect had been measured. One of the other confirmed clusters, ISCS J1438.1+3414 at z = 1.41, was the most distant cluster known at the time of its discovery (Stanford et al. 2005). Follow-up Chandra observations detected the cluster in the X-ray (Andreon et al. 2011; Brodwin et al. 2011). This cluster is found by our algorithm with the highest significance (X_f = 8.64; see the middle panel of Figure 6) in the Boötes field. The only other confirmed high-redshift Boötes cluster with X-ray emission is ISCS J1432.4+3250 at z = 1.49, which is also identified as an overdense region by our algorithm, though at a level slightly below the very conservative threshold adopted here (X_f = 4.2; see the right panel of Figure 6). The red-sequence color and scatter for the known spectroscopically confirmed clusters in Boötes is very similar to the one shown by the SSDF candidates throughout the entire range of detection significance, down to $X_{f_{\rm min}}=5.2$ (as illustrated in Figure 7). To summarize, the conservative cut in detection significance we have derived from this analysis yields many previously identified clusters in the Boötes survey out to z = 1.8, lending confidence in the effectiveness of our cluster-finding algorithm when applied to similar-depth data in the SSDF.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

In order to derive the comoving number density of our cluster sample, we first need to infer the expected redshift distribution ϕ(z) of our sample. With this aim, we employ the COSMOS/UltraVista photometric redshift catalog from Muzzin et al. (2013a). This catalog includes 3.6 μm and 4.5 μm photometry, reaching 2 magnitudes deeper than the SSDF. The central panel in Figure 8 shows the smoothed [3.6]–[4.5] color versus redshift distribution of the COSMOS catalog; details on this filtering procedure can be found in Martinez-Manso et al. (2015). The white line marks the [3.6]–[4.5] > −0.1 selection, and the magenta curve follows the peak of the full galaxy distribution. We parameterize this curve with redshift and denote it as $\mathcal {S}(z)$. At a given redshift, we expect the color distributions of cluster and field galaxies to have the same centroid (given by $\mathcal {S}(z)$), but with a smaller scatter for the clusters given their more homogeneous star formation histories and faster evolution (Rettura et al. 2010). As pointed out in Section 3.1, a significant contamination by z ∼ 0.3 occurs with a simple IRAC color cut. As we remove this contamination with optical and [4.5] magnitude cuts, in the following analysis we assume that all our clusters are at z > 1.

Zoom In Zoom Out Reset image size

To analyze the colors of our cluster candidates, we consider the 10 brightest galaxies in [4.5] from each of them (the total number of member candidates is ∼40 per cluster). The goal of this selection is twofold: it removes faint galaxies, which are more likely to be contaminants due to photometric errors, and keeps the most massive galaxies, which are more likely to share a similar star formation history. The typical scatter in color in clusters after applying this selection is ∼0.12 mag. In comparison, field galaxies at a given redshift have a much larger color scatter of ∼0.25 mag, confirming our expectation that the cluster members are a more homogeneous population.

At this point, we have established that our clusters are likely to follow the redshift–color relation given by $\mathcal {S}(z)$. In order to derive the cluster redshift distribution and its errors, we apply a bootstrap procedure on the cluster colors. In each iteration, a new color for each cluster is drawn from a normal distribution, with the center and standard deviation equal to the mean and standard error of the colors of the 10 brightest member candidates, respectively. These cluster colors are mapped to redshift using $\mathcal {S}(z)$ and binned onto a fixed redshift grid. After 100 iterations, we compute the mean and standard deviation of all values in each redshift bin, yielding the cluster redshift distribution shown in the top panel of Figure 8.

We can calculate the spatial number density of the cluster sample at the pivot redshift z_p ≡ 1.5 by combining the redshift distribution ϕ(z), the SSDF survey area, and the number of cluster candidates, N_obs. The number of clusters within z_p ± δz/2 can be written as

$$\begin{equation} \Delta N = N_\mathrm{obs} \frac{\phi (z_p)}{\int \phi (z^\prime)dz^\prime } \delta z. \end{equation} \tag{ 2 }$$

The sampled volume is

$$\begin{equation} \Delta V = \frac{dV(z_p)}{dz}\delta z = \frac{c \, \Omega \, \chi ^2(z_p)}{H(z_p)}\delta z, \end{equation} \tag{ 3 }$$

where χ(z) is the comoving radial distance, H(z) is the Hubble function, c is the speed of light, and Ω = 0.0271 sr is the solid angle subtended by the SSDF survey. Hence, the number density for our SSDF cluster sample, n_c, at z = 1.5 is

$$\begin{equation} n_{\mathrm{c}} = \frac{\Delta N}{\Delta V} = (7.6 \pm 0.6) \times 10^{-7} \; h^{3} \; \mathrm{Mpc}^{-3}. \end{equation} \tag{ 4 }$$

Our infrared selection of galaxies and the ranking system we adopt for our cluster search algorithm is expected to produce a nearly mass-limited cluster sample at z > 1.3. Thus, we assume that our cluster sample is comprised of mostly dark matter halos above some characteristic minimum mass M_min, following the halo occupation framework (Ma & Fry 2000; Seljak 2000; Peacock & Smith 2000; Scoccimarro et al. 2001; Cooray & Sheth 2002; Berlind & Weinberg 2002; Berlind et al. 2003; Kravtsov et al. 2004; Zheng et al. 2005). We model the probability of a halo of mass M to be part of the cluster sample as

$$\begin{equation} N_c(M) = \left\lbrace \begin{array}{@{}ll@{}} 1 & \quad {\rm if} \; M \ge M_\mathrm{min}\\ 0 & \quad {\rm if} \; M < M_\mathrm{min} \end{array} \right\rbrace . \end{equation} \tag{ 5 }$$

Then, the comoving number density of clusters is determined to be

$$\begin{equation} n_c = \int _{M_\mathrm{low}}^{M_\mathrm{high}}dM \frac{dn}{dM}(M) N_c(M|M_\mathrm{min}), \end{equation} \tag{ 6 }$$

where (dn/dM)(M) is the halo mass function from Tinker et al. (2010) and the integration limits are M_low = 10¹¹M_☉ and M_high = 10¹⁶M_☉. In addition, we can calculate the mean mass of the cluster sample as an average over the occupied halos:

$$\begin{equation} M_\mathrm{mean} =\frac{1}{n_c}\int dM \frac{dn(M)}{dM}N_c(M|M_\mathrm{min})\,M. \end{equation} \tag{ 7 }$$

An additional observable that links the distribution of dark matter halos and our clusters is the measurement of their clustering. In particular, we focus on the modeling of the two-point spatial correlation function (SCF, formally ξ(r); Peebles 1980), which represents the excess probability of finding two objects at a separation r with respect to the case of a randomly distributed sample. The formation of dark matter halos of a given mass M is directly correlated with the dark matter overdensity at their location (Kaiser 1984; Fry & Gaztanaga 1993; Mo & White 1996). This translates into a scaling between the dark matter and halo SCFs, which is called the halo bias b_h:

$$\begin{equation} \xi _h(M,r,z) = \xi _{dm}(r,z) \, b_h^2(M,z). \end{equation} \tag{ 8 }$$

Our definition of halo mass is that enclosed by a sphere with a density 200 times larger than the critical density of the universe. We obtain the dark matter correlation function from the CAMB software (Lewis & Bridle 2002) and the halo bias from Sheth et al. (2001), using the updated parameters from Tinker et al. (2005). The bias of the clusters with respect to dark matter is then

$$\begin{equation} b_c = \frac{1}{n_c}\int dM \frac{dn(M)}{dM}N_c(M|M_\mathrm{\rm min})b_h(M). \end{equation} \tag{ 9 }$$

Since the observed configuration space of our cluster sample is the celestial sphere, we measure its clustering via the angular correlation function (ACF, formally ω(θ)). The ACF can be considered the radial projection of the SCF, which can be computed with the knowledge of the sample's redshift distribution (Limber 1953; Phillipps et al. 1978):

$$\begin{equation} \omega (\theta) = \frac{2}{c} \int _0^\infty dz H(z) \phi ^2(z) \int _0^\infty dy \,\, \xi _c (r=\sqrt{y^2+D_c^2(z)\theta ^2}), \end{equation} \tag{ 10 }$$

where ϕ(z) is the normalized redshift distribution, D_c(z) is the radial comoving distance, c is the speed of light, and θ is the angular separation given in radians. We measure ω(θ) with the estimator presented in Hamilton (1993), which counts the number of galaxy pairs with respect to those of a random sample distributed in the same geometry:

$$\begin{equation} \hat{\omega }(\theta)=\frac{\mathrm{RR}(\theta) \mathrm{GG}(\theta) }{\mathrm{(GR)}^2(\theta)} - 1, \end{equation} \tag{ 11 }$$

where GG, GR, and RR are total number of galaxy–galaxy, galaxy–random and random–random pairs separated by an angle θ. There is no need to include a correction for the integral constraint (Peebles 1980), since it was shown in Martinez-Manso et al. (2015) that it is negligible for this very wide field survey geometry. We estimate ω(θ) errors using jackknife resampling. For this, the entire sample is divided into N_jack = 32 spatial regions of equal size. Then, the correlation is run N_jack times, each time excluding one of those regions from the sample. The value of the estimator is the average $\bar{\omega }(\theta)$ of those iterations and the covariance between angular bins is given by

$$\begin{equation} C_{jk} = \frac{N-1}{N}\sum _{i=0}^N \, [ \hat{\omega }_i(\theta _j)-\bar{\omega }(\theta _j) ] \, [ \hat{\omega }_i(\theta _k)-\bar{\omega }(\theta _k)]. \end{equation} \tag{ 12 }$$

The observed ω(θ) function is shown in Figure 9, where the error bars are derived from the diagonal elements of the covariance matrix.

Zoom In Zoom Out Reset image size

4.3.1. Results

The halo model of the clustering described so far can be fully specified either by fixing n_c through the observed cluster counts (Equation (4)), or by fixing b_c through the amplitude of the observed clustering. We refer to these two approaches as the density and clustering fits, respectively. Thus, comparing these sets of results is an excellent way to test the consistency of the model and observations. Our fits are based on the maximization of the likelihood of the model given the data, $\mathcal {L}($mod|data$)=e^{-\chi ^2/2}$. For the clustering fit, this is specified by

$$\begin{equation} \chi ^2 = \sum _{i=0}^K \sum _{j=0}^K [ \omega _m(\theta _j)-\bar{\omega }(\theta _j) ] C_{ij}^{-1} [ \omega _m (\theta _k)-\bar{\omega }(\theta _k)]. \end{equation} \tag{ 13 }$$

Here, ω_m and $\bar{\omega }$ are the predicted and observed ACFs (Equations (10) and (11), respectively), C_ij is the covariance matrix from Equation (12), and K = 4 is the number of angular bins.

For the density fit,

$$\begin{equation} \chi ^2 = \left(n_c^\mathrm{mod}-n_c^\mathrm{data} \right)^2/\sigma ^2_{n_c}, \end{equation} \tag{ 14 }$$

where the variance $\sigma ^2_{n_c}$ is derived from Poisson statistics in the cluster counts. We have neglected the errors from ϕ(z) after checking that their contribution to the final error budget is almost an order of magnitude smaller than the other sources being accounted for in either fit. Results are shown in Table 2. We have included the calculation of r₀, which marks the distance scale length where ξ_c(r₀) = 1, with $\xi _c = \xi _{dm}(z=z_p) b_c^2$ for both methods.

Table 2. Best-fit Results for the Number Density, the Characteristic Minumum Mass, the Mean Mass, the Bias, and the Angular Correlation Length of the SSDF Cluster Sample as Found by the Density and the Clustering Fits

	Density	Clustering
nc (10−7 h3 Mpc−3)	7.6 ± 0.6	$0.7^{+6.3}_{-0.6}$
Mmin (1014 h−1 M☉)	0.78 ± 0.02	$1.5^{+0.9}_{-0.7}$
Mmean (1014 h−1 M☉)	1.08 ± 0.02	$1.9^{+1.0}_{-0.8}$
bc	8.2 ± 0.1	10.8 ± 2.5
r0 (h−1 Mpc)	25.9 ± 0.4	32 ± 7

Download table as:  ASCIITypeset image

Overall, there is a reasonable agreement between the density and clustering sets, differing by less than 2σ. Throughout the paper, we make the conservative choice of adopting the clustering sets as fiducial values since they have the larger errors. The two angular correlation functions, ω(θ), one corresponding to the prediction based on the number density (solid curve) and one based on the fit of the measured ACF points (dashed curve), are displayed in Figure 9 and found to be consistent within the errors.

The values of n_c and r₀ found for our sample of candidate clusters are consistent with those found by previous observational studies in the literature (Abadi et al. 1998; Collins et al. 2000; Bahcall et al. 2003; Brodwin et al. 2007, see also Figure 9 of Papovich 2008) and are well in agreement with predictions based on ΛCDM cosmological models (Springel et al. 2005). This analysis lends further confidence in the effectiveness of our cluster-finding algorithm.