text:
# Detecting Clusters of Galaxies in the Sloan Digital Sky Survey I : Monte Carlo Comparison of Cluster Detection Algorithms

## 1 Introduction

Over the past few decades, clusters of galaxies have been used as valuable tools for cosmological studies: they are tracers of large scale structure (Bahcall 1988, Huchra et al. 1990, Nichol et al. 1992, Nichol et al. 1994, Peacock & Dodds 1994, Collins et al. 2000), their number density is a constraint on cosmological models (Gunn & Oke 1975, Hoessel, Gunn & Thuan 1980, Evrard 1989, Bahcall et al. 1997, Carlberg et al. 1997, Oukbir & Blanchard 1997, Reichart et al. 1999), and they act as laboratories for probing the formation and evolution of galaxies and their morphologies (Dressler 1980, Butcher & Oemler 1984, Gunn & Dressler 1988, Dressler et al. 1997). While it is sometimes sufficient simply to have a large sample of clusters, most cosmological studies require a homogeneous sample with accurate understanding of the selection biases and the completeness of the catalog.

The Abell cluster catalog is by far the most widely used catalog to date (Abell 1958; Abell et al. 1989, hereafter ACO). Like some of the other cluster catalogs that are available (e.g., Zwicky et al. 1968; Gunn et al. 1986), the Abell catalog was constructed entirely by visual inspection of photographic plates. Although the human eye is a sophisticated and efficient detector for galaxy clusters, it suffers from subjectivity and incompleteness, and visual inspection is extremely time consuming. For cosmological studies, the major disadvantage of such visually constructed catalogs is that it is difficult to quantify selection biases and the selection function.

The main motivation for automated cluster finding schemes is thus to overcome these major drawbacks of visual catalogs, namely speed, objectivity and reproducibility. The first automated cluster finding in optical surveys was attempted by Shectman (1985), and was followed by several variants of this peak finding method; Lumsden et al. (1992) with the EDSGC (Nichol, Collins, & Lumsden 2001) and Dalton et al. (1992) with the APM survey (Maddox et al. 1990). Lahav & Gull (1989) introduced the concept of fitting the observed apparent diameter distribution to the absolute distribution to obtain an estimate of the distance to each cluster, which led to the development of the Matched Filter Algorithm (Postman et al. 1996, hereafter P96). The matched filter technique has been widely used, and several variants have been further developed (Kawasaki et al. 1998; Schuecker & Böhringer 1998; Lobo et al. 2000), including the Adaptive Matched Filter (Kepner et al. 1999; hereafter K99, Kepner & Kim 2000). In addition, thanks to the availability of multi-band CCD imaging data, several cluster finding methods have been put forward that utilize the color characteristics of galaxy clusters (Annis et al. 2001, Gladders & Yee 2000, Goto et al. 2001, Nichol et al. 2000b). While these methods can be efficient, each inevitably has its own biases, depending on the nature and extent of constraints that have been imposed for the selection. Little effort has been made to date to compare the performances of these different methods.

In this paper, we present a comparison of three different cluster finding algorithms, using a Monte Carlo experiment with simulated clusters. The algorithms we analyze are the Matched Filter algorithm (P96), the Adaptive Matched Filter algorithm (K99), and the Voronoi Tessellation Technique, which is introduced in detail in §2.2. The last technique, which uses color information, is in part based on previous cluster finding efforts that use the classical Voronoi Tessellation as a peak finder (Ramella et al. 1998, Ebeling 1993, Kim et al. 2000), but is introduced here for the first time in its current form. These three methods constitute the basis of a cluster catalog derived from 150 $`\text{deg}^2`$ of Sloan Digital Sky Survey (York et al. 2000; hereafter SDSS) commissioning data, which we present in Paper II (Kim et al. 2001). Work in progress (Bahcall et al. 2001) will present cluster catalogs using a wide range of finding techniques. This paper is thus geared towards understanding the behaviour of the cluster finding algorithms in the SDSS data.

In §2 we describe the various cluster detection algorithms, and in §3 we present the Monte Carlo experiment in which these cluster finding algorithms were run. The results and comparison between methods are presented in §4. We summarize in §5. Throughout this paper we assume a cosmology in which $`\mathrm{\Omega }_m=0.3`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$ and $`H_0=70\text{km}\text{s}^1\text{Mpc}^1`$ unless noted otherwise.

## 2 Cluster Selection Algorithms

The Matched Filter algorithm, which was first presented as a fully automated cluster finding scheme in P96, has been widely adapted for a variety of cluster detection efforts (e.g., Olsen et al. 1999, Scodeggio et al. 1999, Bramel et al. 2000; hereafter BNP00, Postman et al. 2001, Willick et al. 2001). The Adaptive Matched Filter algorithm (K99, Kepner & Kim 2000), which is described below, should be a substantial improvement both in content and efficiency over the original Matched Filter algorithm; the major changes being the adoption of a full likelihood function and the incorporation of 3-D (redshift) information. Nevertheless, unlike the P96 Matched Filter, the Adaptive Matched Filter algorithm has not yet been applied to real data, and therefore lacks the optimizations and adaptations that exposure to real data would give. Thus we have chosen to apply both the original Matched Filter (hereafter MF) and the Adaptive Matched Filter (hereafter AMF) for comparison and cross-checking purposes. Although the AMF should in theory converge to the MF results in the 2-D case, there are various differences in the details of the two codes (e.g., peak selection criteria, final parameter evaluations, likelihood function detailed further below) that can cause the final results to differ somewhat. The fact that they are two completely independent codes written in two different languages (C and IDL) also makes the cross check particularly useful.

The Matched Filter technique is an efficient likelihood method for finding clusters in two dimensional imaging data. A model cluster radial profile and galaxy luminosity function are used to construct a matched filter in position and magnitude space from which a cluster likelihood map is generated. Using the magnitudes rather than simply searching for density enhancements suppresses false detections that occur by chance projection. Of course, the results are dependent on the assumed filter shape and extent; they are thus affected by the assumptions made for the universal radial profile and cluster luminosity function used in the algorithm. In other words, the cluster parameters that are derived from the algorithm are somewhat dependent on the model that was assumed for the cluster, and the probability for detecting a certain cluster may differ as the cluster shape varies (e.g., spherical vs. elongated), or as the cluster parameters deviate from the assumed cluster model. In fact most clusters are elongated, and are known to have a variety of shapes (e.g., Bautz-Morgan (1970) type, Rood-Sastry (1971) type). The third method, the Voronoi Tessellation Technique (hereafter VTT), was in part motivated by this model-dependent aspect of the Matched Filter algorithms (MF and AMF), to examine if any severe biases occur in the selection due to model assumptions. In addition to the fact that the VTT assumes no intrinsic cluster properties (except for very mild constraints in color-magnitude space, see below), it is simple and fast (for $`n=10^6`$ galaxies, it takes seconds on a 400 MHz CPU to evaluate Voronoi Tessellation).

The details of the MF and the AMF are given in the respective references, so here we only describe them briefly (§2.1), mainly comparing them carefully and emphasizing their differences. This is followed by a detailed recipe for the VTT in §2.2.

### 2.1 The Matched Filter Techniques

The foundation for both matched filter techniques is the model for the total number of galaxies per unit area per unit observed flux $`l`$:

$$n_{\text{model}}(\theta ,l,z_c)dAdl=[n_f(l)+\mathrm{\Lambda }_{cl}n_c(\theta ,l,z_c)]dAdl,$$
(1)
which consists of contributions from background field galaxies ($`n_f`$) and the cluster galaxies ($`\mathrm{\Lambda }_{cl}n_c`$) at redshift $`z_c`$, where $`\mathrm{\Lambda }_{cl}`$ is the cluster richness measure (see below). Here, $`\theta `$ is the angular distance from the cluster center, and $`dA=2\pi \theta d\theta `$. The background number density ($`n_f`$) is simply taken from the global number counts of the survey. The essential ingredient of the matched filter is the model for the cluster number density, $`\mathrm{\Lambda }_{cl}n_c`$, which is a product of a projected cluster density profile and a luminosity profile for a cluster at redshift $`z_c`$:

$$\mathrm{\Lambda }_{cl}n_c(\theta ,l;z_c)=\mathrm{\Lambda }_{cl}\mathrm{\Sigma }_c(r)\left(\frac{dr}{d\theta }\right)^2\varphi _c(L)\left(\frac{dL}{dl}\right),$$
(2)
where $`r`$ is the projected comoving radius and $`L`$ is the absolute luminosity. The conversion from physical units to apparent units includes proper treatment for cosmology and the K-correction. For the projected density profile we use the modified Plummer law model (see K99 Appendix A),

$$\mathrm{\Sigma }_c(r)\{\begin{array}{c}[1+(r/r_c)^2]^{n/2}[1+(r_{\text{max}}/r_c)^2]^{n/2}\text{for}(r<r_{\text{max}})\hfill \\ 0\text{for}(rr_{\text{max}})\hfill \end{array}$$
(3)
where $`n2`$, the slope of the profile, $`r_c`$ is the core radius and $`r_{\text{max}}`$ is a cutoff radius which approximates the extent of a cluster. This cutoff radius naturally constrains the radius of the search; forcing $`\mathrm{\Sigma }_c(r_{\text{max}})=0`$ is equivalent to putting a “taper” on the filter which reduces the “sidelobes” and narrows the “beam”. This is standard practice in detection theory, which improves the spatial accuracy and makes the process more robust (at the cost of slightly reducing the sensitivity of the filter), as well as reducing the contamination that arises from other nearby clusters.

Any method of smoothing the data consists of choosing a filter shape and a filter bandwidth. Numerical experiments show that the efficiency of the estimator is much more sensitive to the filter bandwidth than is the filter shape itself (Silverman 1986). Hence determining the appropriate values of $`r_{\text{max}}`$ and $`r_c`$ is more important than our particular choice of the cluster profile functional form. P96 (see §4 in their paper) discuss the effect that the cutoff radius has on their detection efficiency and conclude $`r_{\text{max}}=1`$$`h^1`$ Mpc as an optimal choice, which we adopt as an appropriate value. Increasing the value of $`r_{\text{max}}`$ further will significantly degrade the signal-to-noise ratio; the cluster signal will go down since less weight is given to the core, and noise from the background and nearby clusters will increase as more weight is given to galaxies at larger radii.

For the luminosity profile we adopt a standard Schechter luminosity function (Schechter 1976),

$$\varphi _c(L)dL(L/L^{})^\alpha e^{L/L^{}}d(L/L^{}).$$
(4)
The overall normalizations of $`\mathrm{\Sigma }_c`$ and $`\varphi _c`$ are chosen such that the cluster richness measure $`\mathrm{\Lambda }_{cl}`$ is the total cluster luminosity within $`r_{\text{max}}`$ in units of $`L^{}`$, i.e., $`L_{cl}(r_{\text{max}})=\mathrm{\Lambda }_{cl}L^{}`$ (see P96 and K99 for details). The Abell richness $`N_A`$ is defined by the number of cluster galaxies (within $`r<1.5`$$`h^1`$ Mpc) with magnitudes between $`m_3`$ and $`m_3+2`$, where $`m_3`$ is the magnitude of the third brightest galaxy in the cluster. The Abell Richness Class (RC; Abell 1958) is determined by this quantity $`N_A`$; $`30N_A49`$ for RC = 0, $`50N_A79`$ for RC = 1, $`80N_A129`$ for RC = 2, and $`130N_A199`$ for RC = 3. The relation between $`N_A`$ and $`\mathrm{\Lambda }_{cl}`$ (with $`r_{max}=1`$$`h^1`$ Mpc) is found to be approximately $`N_A(2/3)\mathrm{\Lambda }_{cl}`$ (Bahcall & Cen 1993, P96) but with large scatter (P96). This relation is addressed further in Paper II with the clusters detected from the Sloan Digital Sky Survey imaging data.

The difference between the MF and the AMF starts from the definition of the likelihood function. The MF adopts a Gaussian likelihood function, which is based on the assumption that there are enough galaxies in each virtual bin ($`j`$) in (position, magnitude) space that the Poisson probability distribution can be approximated by a Gaussian. Furthermore, it assumes that the background galaxy distribution ($`n_f`$) is uniform and large enough to dominate the noise, therefore the likelihood can be written as,

$$(\theta )=\underset{j}{}2\mathrm{ln}P_j=2\underset{j}{}\mathrm{ln}\left[\frac{1}{\sqrt{2\pi n_f^{(j)}}}\mathrm{exp}\left(\frac{(n_{\text{data}}^{(j)}n_{\text{model}}^{(j)}(\theta ))^2}{2n_f^{(j)}}\right)\right].$$
(5)
By approximating the summation with an integral, using Eq. (1) for $`n_{\text{model}}`$, setting $`\delta n_c/n_f`$, and dropping all terms irrelevant to $`\theta `$, one obtains a likelihood function that is linear in the data and easy to calculate:

$$=\mathrm{\Lambda }_{cl}\underset{i}{}\delta _i,\text{where}\mathrm{\Lambda }_{cl}=\frac{\underset{i}{}\delta _i}{\delta (\theta ,l)n_c(\theta ,l)𝑑A𝑑l},$$
(6)
where $`i`$ stands for each galaxy, and the sum is over every galaxy within $`r_{\text{max}}`$. The richness measure $`\mathrm{\Lambda }_{cl}`$ is obtained by first solving the equation $`/\mathrm{\Lambda }_{cl}=0`$ (see K99 Appendix C or P96 for details).

In the AMF, this likelihood function is referred to as the “coarse filter”. Since it is simple and easy to calculate, AMF uses Eq. (6) to construct a “coarse” likelihood map in order to select clusters, by identifying peaks from this map. However, the Gaussian approximation breaks down when there are not enough galaxies, i.e., especially for poor clusters or those at high redshift, and in general the parameters that are evaluated by this “coarse filter” are found to be somewhat biased (P96 discuss empirical correction factors for this bias, which we will not discuss here in detail). Hence, the AMF defines a second likelihood function that assumes a Poisson likelihood instead of a Gaussian (Dalton et al. 1994). This is called the “fine filter” and reduces to:

$$_{\text{fine}}=\mathrm{\Lambda }_{cl}N_c+\underset{i}{}\mathrm{ln}(1+\mathrm{\Lambda }_{cl}\delta _i),$$
(7)
where $`\mathrm{\Lambda }_{cl}`$ is obtained by solving

$$N_c=\underset{i}{}\frac{\delta _i}{1+\mathrm{\Lambda }_{cl}\delta _i}.$$
(8)
$`N_c`$ is the number of galaxies expected in a $`\mathrm{\Lambda }_{cl}=1`$ cluster (see K99 Appendix C for the derivation). The resulting fine filter likelihood function is nonlinear and requires more computations to evaluate. However, as it is based on correct statistics, the evaluated cluster parameters (redshift and richness) are expected to be more accurate. The AMF thus adopts a two-layered approach, first to identify clusters by peaks in the coarse filter likelihood map, and then to evaluate proper parameters by the fine filter on those selected cluster positions. The computing time difference between the coarse and the fine filter is due to solving Eq. (8), which must be done at every grid position; using only the fine filter would take approximately 10 times longer. In addition, the AMF approach allows an internal cross check for the evaluation of $`z`$ and $`\mathrm{\Lambda }_{cl}`$.

The MF and AMF have further differences that turn out, as we will see, to be quite important. The MF uses a uniform grid in 2D space on which the likelihood function is evaluated for a series of redshift values that span the desired redshift range for the cluster search. The grid size is a function of redshift assumed; we use a grid size of 1/2 the core radius at each redshift to ensure proper sampling. Hence the result is a series of likelihood maps evaluated for each assumed redshift, i.e., a map indicating the likelihood for the existence of a cluster at that redshift. The AMF instead uses a so-called “naturally adaptive grid”, evaluating the likelihood function at the galaxy positions themselves, meaning that as we go to higher redshift, the effective resolution becomes finer as necessary. Such a non-uniform grid is slightly more complicated to handle, but ensures that a unique galaxy lies at the center of each cluster, and no computation is wasted on unnecessary positions. As a result, instead of producing series of likelihood maps to be stored at each assumed redshift value, the AMF calculates the likelihood values on a redshift grid for each galaxy separately, and records the redshift and likelihood at which the likelihood is maximized. In other words, the outcome is two quantities at each galaxy: the peak likelihood that a cluster lies at this point, and the corresponding estimated redshift. This saves a significant amount of disk space relative to the MF, and is cleverly structured to reduce an intermediate step in the analysis (step 5 in Table 1), but as we shall see below, this approach has drawbacks.

Once the coarse likelihood map is generated, the two algorithms differ in the final steps of the cluster selection process. The MF, having stored likelihood maps for each redshift bin, finds local maxima in each map, and registers them as candidate clusters when they lie above a prescribed threshold for each map (e.g., 95 percentile within the map; approximately $`2\sigma `$ level). This is repeated for every redshift bin, then all the cluster candidates from all redshift bins are combined together to filter out overlaps and to find the most likely redshift for each cluster – the redshift where the peak signal is maximized. Each cluster also has a significance of detection $`\sigma _{\text{det}}`$ which is translated from the pixel distribution of the likelihood signals at the final redshift assuming Gaussian statistics. We perform a final cut on this quantity ($`\sigma _{\text{det}}\sigma _{\text{cut}}`$), to select significant detections. On the other hand, the AMF has already stored the most likely redshift of each galaxy point, and therefore simply locates the position of the highest likelihood signal, registers it as a cluster, eliminates all galaxies around this point within a sphere of a given cluster size ($`r_{\text{max}}`$), and then looks for the next highest peak, and so on until the likelihood signal ($`_{\text{coarse}}`$) drops below a prescribed threshold ($`_{\text{cut}}`$). The difference between the two algorithms is not merely the order in which the procedures are executed (the MF locates clusters in angular space first and then determines the redshift, while the AMF determines the redshift first for every point in space and then filters out the clusters in angular space), but that the AMF uses a threshold in likelihood that is constant, thus redshift independent, while the MF uses a threshold that differs for each redshift. The effect that is caused by these differences is discussed in §5.

Finally, the AMF is completed by evaluating the fine filter on the cluster positions that have been selected, then determining the redshift at which the peak of the fine likelihood occurred, and calculating the richness $`\mathrm{\Lambda }_{cl}`$ for that redshift. Table 1 summarizes the procedures and the parameters outlined above for the two Matched Filter algorithms.

### 2.2 The Voronoi Tessellation Technique

The Voronoi Tessellation made its debut in astrophysics as a convenient way of modeling the large-scale structure of the universe (Icke & van de Weygaert 1987, Ling 1987). With a distribution of seeds (nuclei), Voronoi Tessellation creates polyhedral cells that contain one seed each, enclosing all the volume that is closest to its seed. This is the definition of a Voronoi cell. This natural partitioning of space by Voronoi Tessellation has been used to model the large scale distribution of galaxies. This is achieved by envisioning the seeds to be the expansion centers of “voids”, the planes that intersect two adjacent cells as “walls”, the ridges where three walls intersect as “filaments”, and the vertices where four filaments come together as galaxy clusters (van de Weygaert & Icke 1989).

A slightly different application of the Voronoi Tessellation is to identify X-ray sources by locating the overdensity in X-ray photon counts (Ebeling 1993, Ebeling & Wiedenmann 1993); this is directly related to our application for cluster finding, as we now describe. The galaxy positions are input as the seeds for the Voronoi Tessellation, and the Voronoi cell around each galaxy is interpreted as the effective area that each galaxy occupies in space. Taking the inverse of these areas gives a local density at each galaxy in two dimensions. This information is then used to threshold and select galaxy members that live in highly over-dense regions, which we identify as clusters. We do so by calculating the density contrast at each galaxy position $`\delta (\rho \overline{\rho })/\overline{\rho }=(\overline{A}A)/A`$, where $`A`$ is the area of the Voronoi cells, and $`\overline{A}`$ is the mean area of all cells. We then impose a cut in the density contrast $`\delta >\delta _c`$ to select galaxies in high density environments. One can, in fact, use a more rigorous statistical approach for the detection criteria, using statistics of Voronoi Tessellation for a random distribution of seeds (Kiang 1966; see Ramella et al. 1999, 2001 for details). However, as our approach described below is empirical, we adhere to a simple cut in constant density contrast whose value is tested by a Monte Carlo method using simulated clusters, as described in section 3.3.

The SDSS is currently working on determining photometric redshifts for galaxies in the imaging data, but until they are available, we are confined to working with the two dimensional projected distribution; therefore we need to divide the galaxy sample to group them into comparable redshifts in order to enhance the cluster detectability. Therefore, our recipe for the Voronoi Tessellation Technique utilizes a priori knowledge of characteristics of cluster member galaxies, namely, the color-magnitude relation.

Galaxies within a cluster usually exhibit a tight correlation between their colors and magnitudes. It is well known that the core of a typical rich cluster consists mainly of early type galaxies (i.e., Hubble Types E and S0; Hubble 1936, Oemler 1974, Postman & Geller 1984, Dressler 1980, 1984) which all have very similar red colors. This includes the Brightest Cluster Galaxy (hereafter BCG), whose properties have been well studied (Schneider, Gunn & Hoessel 1983, Postman & Lauer 1995). Figures 1a and b show the color-magnitude relation (C-M diagram) of Abell Clusters 168 and 295 respectively, observed with the SDSS camera; the diagram shows only those galaxies whose cluster membership has been confirmed spectroscopically by the ENACS survey (Katgert et al. 1998). The BCG is marked with a cross, and the narrow horizontal line of galaxies at nearly constant color is referred to as the E/S0 ridgeline (Visvanathan & Sandage 1977, Annis et al. 1999). The color-magnitude relation for E/S0 galaxies has been well known since Baum (1959), who first noted that fainter early type galaxies tend to be bluer, showing a negative slope of the E/S0 ridgeline in C-M diagrams (see also Visvanathan & Sandage 1977, Lugger 1984, Metcalfe et al. 1994). It has been suggested that this slope evolves with redshift (Gladders et al. 1998) and even differs with richness (Stoughton et al. 1998). However, as we see in the C-M diagrams, the slope is very shallow in $`(g^{}r^{})`$ vs. $`r^{}`$ space, and is a negligible effect for the recipe that we now describe.

We will use this characteristic shape in the C-M diagram to select against foreground and background galaxies. Thus we carry out the following approach: first select a redshift, then apply a Voronoi Tessellation on all galaxies in a restricted region of the C-M diagram. This region is shown as solid lines in Figures 1a and b, for a redshift of $`z=0.045`$, whose limits enclose most of the galaxies that are confirmed members of the cluster. Figures 1c and d show similar C-M diagrams of cluster candidates at higher redshifts, found by both the MF and the AMF from the SDSS. These clusters were visually confirmed with SDSS images, and the redshifts of their BCGs were obtained by the SDSS spectroscopic survey (York et al. 2000). Since membership information is not available, we simply plot all the galaxies that are within 1$`h^1`$ Mpc radius from the cluster center as circles. For comparison, the C-M distribution of all the galaxies in a survey region of 25 $`\text{deg}^2`$ is shown in contours and small dots. The difference in these distributions illustrates the efficiency of using C-M information in discriminating cluster members against the background population. Here again, the areas enclosed by the solid lines are those that are selected for detecting a cluster at $`z=0.22`$ and $`0.35`$ respectively, which includes most of the galaxies within $`r=1`$$`h^1`$ Mpc around the cluster center.

The empirical C-M limits we use are as follows:

$`r_{bcg}^{}1<`$ $`r^{}`$ $`<r_{bcg}^{}+5`$ (9)
$`ar^{}+b<`$ $`(g^{}r^{})`$ $`<(g^{}r^{})_{bcg}+0.3`$ (10)
where $`a`$ and $`b`$ are given by a simple linear relation with redshift,

$$a=1.295z0.084,b=30.13z+0.88$$
(11)
where $`r_{bcg}^{}`$ and $`(g^{}r^{})_{bcg}`$ are the Petrosian $`r^{}`$ magnitude and the model $`g^{}r^{}`$ color of the BCG, respectively (see §3.1 for description of magnitudes and colors). These relations are established empirically, by examining the C-M relations of known clusters (e.g., Fig. 1 a,b) and rich clusters of galaxies found by the Matched Filter algorithms and the MaxBCG technique (Annis et al. 2001) from the SDSS data itself. The limits are chosen to generously include most of the galaxies within 1$`h^1`$ Mpc of the cluster center. The lower limit in magnitude, 5 magnitudes fainter than that of the BCG, is chosen to cover the magnitude range of spectroscopically confirmed cluster members from the ENACS (Katgert et al. 1998) for low redshift clusters ($`z<0.1`$). This is a moderate coverage of cluster galaxies for a typical cluster luminosity function (corresponding to $`M_{lim}L^{}+3`$); also note the recent findings of fossil groups whose difference between the BCG and the second luminous galaxy can be as large as 4-5 magnitudes (Zabludoff & Mulchaey 1998).

These C-M limits are shown in summary for a range of redshifts, from $`z=0.04`$ to $`z=0.5`$, in Figure 2. Here, the track of BCGs for the above redshifts are also shown as large dots (filled and open), which provides the basis for determining these limits for each redshift (Eqs. (9) and (11)). Finally, the C-M distribution of all galaxies in the survey region are shown as contours and small dots for comparison. The BCG track is computed assuming a constant absolute magnitude $`M(r^{})_{bcg}=23`$ with no evolution (see Eisenstein et al. 2001). We use the PEGASE evolutionary synthesis model (Fioc & Rocca-Volmerange 1997) to generate the spectral energy distribution of the BCG; this is then redshifted and convolved with the SDSS filter responses to obtain color K-corrections.

Once we have applied these cuts for a given redshift, we apply the Voronoi Tessellation on the resulting distribution of galaxies. We then select all galaxies that satisfy $`\delta >\delta _c`$, where $`\delta _c`$ is a constant overdensity threshold defined above. We will set this threshold to $`\delta _c=3`$. Figures 3 and 4 illustrate this procedure: Figure 3 shows Voronoi Tessellation executed on all galaxies in the region of Abell 957, whereas Figure 4 shows only those galaxies that satisfy the C-M cut for the cluster redshift ($`z=0.045`$). The large circle indicates a 1$`h^1`$ Mpc radius around the cluster center. In both cases, the dots highlight the overdense galaxies with $`\delta >3`$. Figure 4 shows remarkable enhancement of the cluster while Figure 3 does not. As we can see, although the C-M limits are chosen generously in order to ensure proper coverage of the observed cluster galaxy characteristics, using these limits provides an important enhancement in the efficiency of cluster selection due to the elimination of a significant background population.

Once we have selected galaxies that highlight densely populated regions, we have to identify regions with a concentration of these “high-density” galaxies. This is done by selecting regions in which the number of such galaxies, $`N_{hdg}`$, within a radius of 0.7$`h^1`$ Mpc at the assumed redshift, exceeds a certain threshold, such that $`N_{hdg}N_{cut}`$. This is executed around each “high-density” galaxy, and the center of a cluster is determined by the one that gives the largest $`N_{hdg}`$. We repeat this process as a function of redshift. Once we obtain all the cluster candidates for each redshift bin, we go through a process to filter out significant overlaps to finalize our selection of clusters. For every significant overlap we choose the final cluster and its redshift to be the one that yields the largest value of $`N_{hdg}`$. The distribution of $`N_{hdg}`$ with respect to redshift is generally highly peaked and therefore justifies this simple method of determining the redshifts.

Here, we take the high-density galaxy number cut $`N_{cut}`$ to be a constant value, independent of redshift. This may introduce a bias with redshift; the same cluster at higher redshift will contain fewer galaxies, since the faint magnitude limit in the C-M cuts (5 magnitudes fainter than the BCG) soon exceeds the survey magnitude limit as we go to higher redshift. In order to find the optimal threshold that changes with redshift, the algorithm needs to be tested carefully to assure ourselves that the contamination level (false-positives) does not increase too much. We do not carry out such an analysis here, as it requires quantitative assessment of false-positives which is only possible with full N-body simulations, with complete knowledge of cluster identities. This would also require proper assignments of colors for the background and cluster galaxies. Another way to investigate this matter is to use the real data itself; although properly identifying false-positives can be a slightly tricky business, we do address this issue of variable $`N_{cut}`$ using visual inspection of cluster candidates in Paper II. However for the current paper, we keep our VTT threshold $`N_{cut}`$ constant. Therefore, as with every method, having such potential biases motivates us to evaluate the selection function to assess the fraction of clusters selected at each redshift and richness (§4.2). In addition, $`N_{hdg}`$ is not intended to measure the richness of the cluster, but is rather a measure of the significance of detection. In Paper II, the final VTT selected clusters will be run through the AMF fine filter for consistent estimation of the cluster richness and redshift.

## 3 Testing the Algorithms Using Simulated Clusters

The MF was originally developed for the Palomar Distant Cluster Survey (PDCS; P96) which is deep ($`V23.8`$) and narrow (5.1 $`\text{deg}^2`$), and has been applied to similar data by others (e.g., Scodeggio et al. 1999 : 12 $`\text{deg}^2`$, $`I<23`$; Postman et al. 2001 : 16 $`\text{deg}^2`$, $`I_{AB}<24`$). Naturally, the algorithm is optimized for this type of survey, whereas the SDSS is shallower ($`z<0.5`$) and much wider (Paper II will present results for data over 150 $`\text{deg}^2`$, which is less than 2% of the complete SDSS survey). The shallower depth makes the large-scale structure variations much more pronounced in the two dimensional distribution, which can affect the matched filter algorithm’s performance, as it assumes a uniform and homogeneous background (BNP00). Second, covering a large area increases the probability of intersecting very nearby clusters ($`z\stackrel{<}{}0.05`$) which have angular extents as large as a few degrees. This increases the rate of cluster overlaps, especially because we go to $`z0.5`$ for the SDSS. In a 2D projection, discriminating between two different overlapping clusters at two different redshifts using the algorithms outlined above can be difficult unless the redshift gap is sufficiently large (usually $`\mathrm{\Delta }z\stackrel{>}{}0.4`$); this can affect the completeness of intermediate redshift clusters. The narrow pencil beam surveys on the other hand, are carried out in regions known not to have foreground clusters so they suffer less from this effect. Wide angle, shallow surveys such as the APM survey or the EDSGC have less range in redshift, and thus also suffer less from this effect. Combining the local space density of Abell clusters and P96 results, and assuming an unclustered distribution of clusters for simplicity yields a $`15\%`$ rate of overlap for Abell richness class $`0`$ clusters (using a 1$`h^1`$ Mpc radius for each cluster) for the redshift range of $`0.05z0.5`$. This rate of overlap will be further enhanced by taking the non-zero two point correlation function into account, and shall be addressed further in future work with the observed space density of clusters from the SDSS.

Each method for detecting clusters has its own biases; moreover, the sample selected depends sensitively on the detection threshold. For example, a 20% change in $`\sigma _{cut}`$ can result in doubling the final number of clusters in a given field. Therefore it is crucial for any cluster identification study to understand the behavior of the results with respect to the selected thresholds. Hence in this section we attempt to understand the parameters and their limits that play a crucial role in determining the final cluster selection, for the three algorithms.

Tests of the algorithms are performed on a set of artificial clusters embedded in two different versions of a background galaxy distribution: a uniform background and a clustered one. The latter is taken to be the real galaxy distribution itself, from the SDSS. The major objectives for these tests are to determine the best detection thresholds for the final cluster catalog, to evaluate a realistic selection function for these thresholds and to examine biases with respect to the local background density.

In the following we describe these tests, starting with a brief description of the SDSS data that we used for the background.

### 3.1 Sloan Digital Sky Survey Imaging Data

The SDSS imaging data is taken with an imaging camera (Gunn et al. 1998) on a wide-field 2.5 telescope, in five broad bands ($`u,g,r,i,z`$) centered on \[3551Å, 4686Å, 6166Å, 7480Å, and 8932Å\] respectively (Fukugita et al. 1996, Stoughton et al. 2001), in drift-scan mode at sidereal rate. This results in an effective exposure time of $`54.1s`$, which yields a point source magnitude limit of $`r^{}22.5`$ (at $`1.5^{\prime \prime }`$ seeing). More details of the observations and data are covered in York et al. (2000), Stoughton et al. (2001) and Paper II (and references therein).

The data used here and in Paper II are Equatorial scan data taken in September 1998 during the early part of the SDSS commissioning phase, and are part of the SDSS Early Data Release (runs 94 and 125; Stoughton et al. 2001). A contiguous area of about 150 $`\text{deg}^2`$ was obtained during two nights, where the seeing varied from $`1.1^{\prime \prime }`$ to $`2.5^{\prime \prime }`$ (85% of the data was below $`1.8^{\prime \prime }`$), and coordinates ranging from $`5^{}<\text{RA}<55^{}`$ and $`1.25^{}<\text{Dec}<1.25^{}`$. We include galaxies to $`r^{}=21^m`$ (Petrosian magnitude; see below), a conservative limit at which star-galaxy separation is reliable (see Paper II for details). For the present paper, a subset of 25 $`\text{deg}^2`$ from this data was taken for the test region described in §3.3. The coordinates of this region are $`10^{}<\text{RA}<20^{}`$ and $`1.25^{}<\text{Dec}<1.25^{}`$, which was chosen because it exhibits prominent large scale structure and clumpiness.

The magnitude of all objects quoted here are measured in Petrosian quantities (Petrosian 1976) through the SDSS photometric pipeline (Lupton et al. 2001). However, the colors of each object quoted are computed from “model magnitudes”. Each galaxy is fit to two profiles in the $`r`$ band : a de Vaucouleurs law and an exponential law. The model magnitudes in all five bands are computed from the better of the two $`r`$ band fits; the colors obtained from model magnitudes are thus a meaningful quantity as it uses the same profile in all bands.

### 3.2 The Uniform Background Case and the Hybrid Matched Filter

In order to produce a uniform galaxy background distribution, we took the 25 $`\text{deg}^2`$ region from the SDSS data described above, and randomly repositioned the galaxies while keeping their photometric properties fixed. This creates a uniform galaxy background distribution while ensuring that the galaxies otherwise have SDSS-like properties (luminosity function and colors). The number-magnitude relation $`N(r^{})`$ to the limiting magnitude $`r^{}=21`$ is presented in Yasuda et al. (2001): it shows power-law behaviour to $`r^{}=17`$, curving below this model at fainter magnitudes, as expected from K-corrections and relativistic corrections.

Artificial clusters with 6 different richnesses were embedded at 8 different redshifts each, giving a total of 48 clusters with different properties. The clusters were generated with a Schechter luminosity function ($`\alpha =1.1,M_{}(r^{})=21.5`$; Blanton et al. 2001) and a modified Plummer law radial profile ($`r_{\text{max}}=1h^1\text{Mpc},r_{\text{core}}=0.1r_{\text{max}}`$; see Eq. (3)). Each cluster, properly normalized according to its richness $`\mathrm{\Lambda }_{cl}`$, was placed at the corresponding redshift and then trimmed to the survey magnitude limit ($`r^{}=21`$). The ($`g^{}r^{}`$) color for each galaxy was assigned characteristic for a cluster at each redshift (as shown in Fig. 2). The insertion of simulated cluster galaxies increased the total number of galaxies in the 25 $`\text{deg}^2`$ region by 7%, to $`105,600`$ galaxies in total. The upper panel of Figure 5 shows the distribution of the 48 artificial clusters over a $`10\times 2.5`$ deg area (the parameters of the clusters themselves are given in the figure caption), and the resulting distribution when inserted into our uniform background is shown in the middle panel.

All three cluster finding algorithms (MF; AMF; VTT) were first run on this distribution of clusters in a uniform background of galaxies. The goal for this was to test each algorithm in the simplest case, and to find a reasonable detection threshold for each of them that maximizes the number of successful detections, while keeping the false detection rate minimal. Although this uniform background case is far from realistic, it has the advantage of unambiguously recognizing false detections. However, note that the number of clusters inserted and their distribution of $`z`$ and $`\mathrm{\Lambda }_{cl}`$ are arbitrary. Thus, neither the fraction of false detections nor the absolute value of the recovery fraction that we quote below have physical significance; it is the relative values for different algorithms that should be noted.

Figure 6 shows the detection efficiencies for all three algorithms; each panel shows the number of successfully recovered clusters (solid curve) and the number of false detections (dotted curve) for one of the algorithms as a function of the detection threshold for the 48 clusters inserted into the uniform background. As the detection thresholds are decreased (going rightward in each panel of Figure 6), the number of successful detections increases, but naturally the rate of false detections due to Poisson statistics increases as well. Beyond a certain threshold value, the number of false detections start to increase extremely rapidly, while the number of successful detections only increases slowly; this implies a large drop in efficiency. Figure 7 shows this in another way, plotting the rate at which the number of successful detection increases as a function of the number of false detections for all three algorithms. We thus find a cut at which the success rate starts to flatten out with respect to the false detection rate; the vertical dotted line in Figure 7 shows an appropriate choice drawn just after the steepest part of the efficiency curve. This cut gives 14 false detections over 25 $`\text{deg}^2`$ for all three methods, $`0.56\text{deg}^2`$, an acceptable level given the expected surface density of real systems ($`>5\text{deg}^2`$; see Paper II). The vertical dashed lines in each panel of Figure 6 also show this cut, corresponding to $`\sigma _{cut}=5.5`$ for the MF and $`_{cut}=210`$ for the AMF. For the VTT the line indicates $`N_{cut}=9`$ (with a constant density contrast cut of $`\delta _c=3`$). We will refer to this cut in subsequent sections as the detection limit chosen for the uniform background case (note that $`\sigma _{cut}`$ will be slightly lower for the final results; §4.1). As Figure 7 demonstrates, of all methods, the MF is the most efficient in recovering clusters for a given number of false detections.

The additional clusters detected by the MF are clusters with a weaker signal (low $`\mathrm{\Lambda }_{cl}`$, high $`z`$), as seen in Figures 8 a and b. These results are obtained using the cuts determined above. Although the AMF (crosses in Figure 8a) should in principle converge to the MF (squares in Figure 8a) in the 2-D case, they differ for weak clusters due to the fact that their peak selection method in the final step is done differently (refer to §2.1; Table 1). The parameter evaluation, on the other hand, is most accurate with the AMF fine filter. Figures 8 c and d show the input parameters for the clusters, $`z`$ and $`\mathrm{\Lambda }_{cl}`$, plotted against the parameters evaluated by the AMF fine filter. They include the additional clusters found by the MF. The MF $`z`$ estimates show a factor of $`1.5`$ increase in standard deviation for higher redshift clusters ($`z>0.25`$), which translates into a larger deviation in $`\mathrm{\Lambda }_{cl}`$ as well, by a similar amount. This result thus calls for a marriage of the two algorithms in order to maximize the efficiency of the final result. We have therefore adopted the following hybrid method for the final Matched Filter based cluster finder (HMF: Hybrid Matched Filter),
$``$ First the MF creates likelihood maps.
$``$ We threshold at $`\sigma _{det}=\sigma _{cut}`$ according to the MF to choose final cluster candidates.
$``$ Finally, we evaluate the AMF fine filter on these cluster positions to determine $`\mathrm{\Lambda }_{cl}`$ and $`z_{est}`$.
We use this recipe as the standard Matched Filter method throughout the rest of this paper.

### 3.3 Monte Carlo with a “Realistic” Background

The previous test, using a uniform background, allowed us to determine the detection threshold necessary to minimize false detections due to projections in random fields. We need to apply a realistic background in order to get an accurate determination of the selection function, including the effect of large-scale structure and cluster projection along the line-of-sight. In relatively shallow surveys like the APM or the EDSGC the cluster detection efficiency has been shown to depend significantly on the local background density (see BNP00). The SDSS is considerably deeper ($`z0.5`$ compared to $`z0.2`$) and therefore this effect should be less dramatic, due to projection in 2D of the real large-scale structure. Nevertheless as the lower panel of Figure 5 shows, the SDSS galaxy distribution is still far from uniform, and the effect of large-scale structure on the selection needs to be quantified.

We now place the same set of simulated clusters on the real SDSS distribution as shown in the lower panel of Figure 5. The effect of the real background is immediately noticeable even by eye; the low richness clusters ($`\mathrm{\Lambda }_{cl}40`$) are nearly washed out and the intermediate richness clusters ($`\mathrm{\Lambda }_{cl}=70,110`$) start to blend in with the clumpy background.

The 48 simulated clusters with the properties given in Figure 5 were inserted into the data at random positions (unlike the grid distribution in the lower panel of Figure 5), but avoiding overlap with one another. Both the HMF and the VTT were run on this catalog. This was repeated 100 times. We use the results to test both the detection and recovery of clusters, and the evaluation of their parameters, $`z`$ and $`\mathrm{\Lambda }_{cl}`$. In the next section we discuss the results from these tests in three parts. First, we study the effect of the imposed detection limits; second, we evaluate the selection function; and finally we examine the dependence of the detection efficiency and the recovered parameters on the local background.

## 4 Results

### 4.1 Detection Limits

A set of clusters embedded in a uniform background distribution is an ideal case. In the real universe, non-uniformity comes in many forms: from large scale modulations (great wall, voids etc.) down to small scale fluctuations (e.g., compact groups, close pairs and triplets etc.), in addition to actual clusters of galaxies. It is mainly the small scale features aided by projection that will cause the false-positive rate to increase. Therefore, the false detection rate that was determined with the uniform background ($`0.56/\text{deg}^2`$ for $`\sigma _{cut}=5.5`$ and $`N_{cut}=9`$) can be regarded as a lower limit for those thresholds in a realistic background. In addition, the structure exhibited on a range of scales diminishes the contrast of real clusters, making them more difficult to recover. This implies that these detection thresholds are upper limits on what they should be in the real universe. As we will see below, the detection threshold has to be set lower with a realistic background in order to recover a similar range of clusters.

For each cluster finder, we determine how many of the clusters that were detected in the uniform background case (with the corresponding detection thresholds $`\sigma _{cut}=5.5`$ and $`N_c=9`$; see Figure 8), were also recovered in each Monte Carlo realization with the SDSS galaxy background. This recovery ratio was evaluated for different detection thresholds, and averaged over all realizations. Figure 9 shows the results for both the HMF and the VTT, for a range of detection thresholds that are equivalent to those used in Figure 7. This recovery ratio (relative to the uniform background case) changes quite rapidly as a function of $`\sigma _{cut}`$ for the HMF, but stays rather robust as a function of $`N_{cut}`$ for the VTT. This comparison can be misleading since there is no formal correspondence between $`N_{cut}`$ and $`\sigma _{cut}`$. However, as noted above, these ranges are calibrated to yield approximately the same number of false detections in a uniform background distribution (c.f., Fig. 6). Given this relation as a yardstick, the noticeable difference in the slope of the two curves suggests that the Matched Filter algorithm results will be more sensitive to the detection threshold that is chosen, while a sample based on the VTT is robust to the exact value of $`N_{cut}`$. This implies in addition, that the Matched Filter is more subject to the effects of a non-uniform background, i.e., the detection efficiency is affected by the background distribution.

This argument is further supported by the lower recovery rate of the HMF shown in Figure 9 at the thresholds of $`\sigma _{cut}=5.5`$ for the HMF and $`N_{cut}=9`$ for the VTT (dotted lines; where both algorithms yielded 14 false detections in the uniform background case). Thus the detection efficiency of the HMF appears to be more affected by a non-uniform background than is the VTT. Indeed, this is not too surprising as the Matched Filter algorithm explicitly assumes a uniform background in its model. In order for the HMF to achieve the same recovery rate as the VTT, the detection threshold would need to be lowered to a value of $`\sigma _{cut}4.7`$ (which greatly increases the number of false detections, this would yield $`60`$ false detections in a uniform background; Fig. 6). These fractional recovery rates refer to the sample of 48 clusters over their entire range of $`z`$ and $`\mathrm{\Lambda }_{cl}`$; hence it does not apply to an observed sample of clusters with a true richness function distribution. The clusters that are missed in the selection are those at the weaker end of the distribution of signal, which are the most abundant in the universe (i.e., poor, distant), therefore the gap between the recovery fractions shown above could be even larger for a more realistic configuration.

These result do not immediately mean, however, that the VTT does better overall. Recall from Figure 6 that the HMF was more efficient in recovering clusters in the uniform background to begin with. Given the right choice of detection thresholds we show that their performances are in fact similar in the end. Figure 10 shows the absolute recovery fraction of inserted clusters, rather than that relative to the uniform background (Fig. 9), for both the HMF (solid) and the VTT (dashed). This is shown in four different subgroups in cluster parameter space: low $`z`$, low $`\mathrm{\Lambda }_{cl}`$ (lower left); high $`z`$, low $`\mathrm{\Lambda }_{cl}`$ (upper left); low $`z`$, high $`\mathrm{\Lambda }_{cl}`$ (lower right); and high $`z`$, high $`\mathrm{\Lambda }_{cl}`$ (upper right). For clusters with the strongest signals (rich and nearby), both algorithms agree well with very high efficiencies, constant with respect to the detection thresholds. For clusters with weakest signals (poor, high $`z`$) both methods have very low recovery rate. Apart from the rich nearby regime, the HMF indeed shows slightly lower efficiencies than the VTT for the thresholds determined above ($`\sigma _{cut}=5.5`$, $`N_{cut}=9`$) suggesting a drop in efficiency due to the presence of a non-uniform background. However, the slopes of the HMF efficiency curves are much steeper, and allow us to bring up the detection efficiency easily to the VTT performance level, by lowering the detection threshold slightly to $`\sigma _{cut}=5.2`$. This naturally increases the false detection rate of the HMF in the uniform background case (refer to Fig. 6) by approximately 80% to $`1\text{deg}^2`$, but as it is less than 20% of the expected surface density of real clusters (Paper II) we adopt this new value for $`\sigma _{cut}`$ as appropriate for selecting clusters from the SDSS imaging data. The final detection thresholds determined for both algorithms, $`\sigma _{cut}=5.2`$ and $`N_{cut}=9`$, are shown as dotted lines in Figure 10, indicating that the performance of the two algorithms is very similar.

### 4.2 Cluster Selection Function

With these detection thresholds, we now assess the overall performances of both cluster finders for these values in more detail, an important task for any study that requires a complete sample. We present selection functions evaluated from the fraction of clusters in each redshift and richness class that are recovered in the Monte Carlo simulations (using a realistic background). These selection functions for the HMF and the VTT are shown in Figures 11 and 12 with $`\sigma _{cut}=5.2`$ and $`N_{cut}=9`$, respectively. The two methods are qualitatively similar in terms of their overall performances; for clusters with $`\mathrm{\Lambda }_{cl}70`$, over 80% are recovered to redshift $`z0.45`$, falling off rapidly as one goes to $`z0.5`$. There are still slight differences; the HMF does better in the high $`\mathrm{\Lambda }_{cl}`$ intermediate $`z`$ ($`0.2\stackrel{<}{}z\stackrel{<}{}0.4`$) domain, while the VTT seems to perform slightly better for clusters with low $`\mathrm{\Lambda }_{cl}`$ and low $`z`$ ($`z\stackrel{<}{}0.2`$). This can be explained in terms of the generous C-M cuts adopted for the VTT, intended to account for possible fluctuations in the cluster color-magnitude properties. As a result, as one goes to higher redshift, where the C-M limits start to move in towards the core of the C-M distribution of normal galaxies (see Fig. 2), the population of galaxies that are rejected by this procedure reduces significantly, making this filtering less effective and therefore reducing the efficiency of recovering the clusters. At low redshifts $`z\stackrel{<}{}0.2`$, on the other hand, it rejects most of the galaxy background, making the search most efficient. If we were to narrow these C-M limits, we might bias ourselves against clusters with unusual properties. This will be tested with real SDSS clusters by investigating their color-magnitude properties, to determine how tightly we can impose the limits without biasing the cluster selection.

We compare our selection functions to those of BNP00 (see their Fig. 2), evaluated for cluster detections in the EDSGC data. The redshift range probed with the EDSGC is much shallower ($`z\stackrel{<}{}0.15`$; $`b_j<20.5`$, which corresponds roughly to $`r<19`$ for a typical elliptical galaxy) than that of the SDSS, however, they use the same Monte Carlo technique by inserting simulated clusters in the real data itself, making the comparison meaningful. Recall that the cluster detection completeness depends significantly on the use of different background distributions, as we have shown in Figures 9 and 10.

Figure 2 of BNP00 shows the selection function for four different richnesses as a function of redshift. Their richnesses correspond to $`\mathrm{\Lambda }_{cl}=10.15,20.3,40.6,\text{and}81.2`$, allowing a straightforward comparison to Figures 11 and 12. For $`\mathrm{\Lambda }_{cl}=20`$ clusters, their efficiency drops from 70% to 40% over the redshift range of $`0.05z0.15`$. Neither the HMF nor the VTT does significantly better for the same range of parameters, and beyond $`z0.2`$, clusters with $`\mathrm{\Lambda }_{cl}=20`$ are virtually invisible even in the SDSS. However, for richer clusters, the SDSS clearly does better. For $`\mathrm{\Lambda }_{cl}=40`$ clusters, BNP00 shows a steep drop in efficiency to $`60\%`$ at $`z=0.15`$, whereas the SDSS efficiencies stay at $`80\%`$ out to $`z=0.2`$. Similarly, the recovery of $`\mathrm{\Lambda }_{cl}=70`$ clusters in the SDSS stays highly efficient ($`\stackrel{>}{}90\%`$) and constant out to $`z0.4`$, while BNP00 efficiency for $`\mathrm{\Lambda }_{cl}=80`$ clusters shows a clear drop from 100% to 80% already at $`z=0.15`$. The above comparison conforms to what we would expect: While going deeper in magnitude helps the recovery of clusters at higher redshifts, there is a limit in the cluster richness for which this is true. In other words, poor clusters are hard to find at high redshifts no matter what the depth of the survey is. Apart from those differences in the SDSS and EDSGC, the overall performances are qualitatively similar, which is reassuring given that our HMF is quite similar to their method.

### 4.3 Dependence on Local Background

The final issue we examine from this Monte Carlo experiment is the extent to which the results depend on the local background density. First, we look at the behavior of the detection efficiency at different locations. We divide our $`2.5\times 10`$ degree Monte Carlo region into $`5\times 20`$ subregions, $`0.5`$ deg on a side. This bin size is chosen to be large enough to get an appreciable number of clusters in each bin, for reliable statistics. We then count the number of background galaxies within each bin $`N_g`$ (excluding cluster galaxies), and evaluate the detection efficiency of clusters that fall into each bin. Figure 13 shows a scatter plot of the cluster detection efficiency as a function of the background density contrast $`\delta _g=(N_g\overline{N})/\overline{N}`$. We have performed Spearman Rank Correlation tests (Press et al. 1990) and find that $`r_s=0.17,t=1.72`$ for the HMF and $`r_s=0.14,t=1.39`$ for the VTT, where $`r_s`$ is the Spearman rank-order correlation coefficient, and $`t=r_s\sqrt{(N2)/(1r_s^2)}`$ which is distributed like a Student’s t distribution with $`N2`$ degrees of freedom. This translates into correlations at only 5% and 8% confidence for the HMF and VTT respectively, confirming that neither method shows any significant correlation between their detection efficiency and the local background density. This is very encouraging, considering the unambiguous dependence that BNP00 have shown for their shallower sample (see Figure 3 in their paper). This suggests that our SDSS sample is deep enough to subdue the background fluctuations by projection effects, to a level that it does not affect the detection efficiency significantly. However, we have seen in §4.1 that the overall performance is still highly influenced by the presence of a non-uniform background, especially for the Matched Filter.

We next investigate how HMF parameter evaluations depend upon the local background density. First, we show in Figure 14 the input and output values of $`z`$ and $`\mathrm{\Lambda }_{cl}`$ for all detections in the Monte Carlo simulation. The distribution of $`\mathrm{\Delta }z`$ is somewhat positively skewed, implying an overestimation of redshifts, which is pronounced as we go to higher $`z`$: the median values of $`\mathrm{\Delta }z`$ are 0.007, 0.02 and 0.05 for input redshift ranges $`0<z0.2,0.2<z<0.35`$ and $`z>0.35`$, respectively. Also, half of the clusters with input redshifts $`z\stackrel{>}{}0.4`$ are recovered with $`z`$ estimates at the upper $`z`$ limit ($`z=0.5`$). This trend was noted in the Matched Filter algorithm by P96, which usually happens with weak signals (poor and/or high redshift clusters), or in the field where there are no clusters, where the number-magnitude distribution of the background is such that the cluster likelihood calculated yields a monotonically increasing function with redshift. This could possibly be remedied by tweaking the input luminosity function to avoid such a trend in the likelihood function (e.g., flatter faint end slope), but this could bias parameter estimation in other ways, and we will not pursue this here.

We also find that the overestimation of redshift is amplified by the use of a non-uniform background. We have carried out the experiment of running the HMF on a given set of simulated clusters in both a uniform and a clustered background. We find that using the real background indeed increases the median redshift by about $`1015\%`$, but also creates a long tail of under-estimated redshifts. We suggest that this is partly due to real clusters in the data that intercept the inserted clusters, thus affecting the output redshifts. The other concern pertains to the centroid of recovered clusters. While we may recover unbiased values exactly at the input center, the recovered clusters are always somewhat “offset” from the real center, and this could lead to some systematic effects. We thus repeated the Monte Carlo simulation using the true centers and the recovered centers, and found no sigificant systematic bias in the richness between the two runs. However, we do find that the number of clusters falsely estimated to be at $`z=0.5`$ was reduced by $`15\%`$ with the true centers.

The $`\mathrm{\Lambda }_{cl}`$ input and output values are tightly correlated and do not show much systematic effects, but there is an extended tail of clusters with overestimated $`\mathrm{\Lambda }_{cl}`$. These are the clusters that were falsely assigned to $`z=0.5`$, which boosts the estimated $`\mathrm{\Lambda }_{cl}`$ to make up for the loss of galaxies beyond the magnitude limit. In what follows, we have left out all detections with $`z_{est}=0.5`$ to avoid such complications.

For each cluster recovered by the HMF, Figure 15 shows the difference between the input value and recovered value of $`z`$ and the ratio of the input and output values of $`\mathrm{\Lambda }_{cl}`$, plotted against the density contrast of the local background, $`\delta _g`$. The value of $`\delta _g`$ at each cluster position is evaluated from counts in cell statistics on the background with a radius of 5 arcmin. This radius is chosen to represent the immediate local background of a cluster ($``$ 1 Abell radius at the median redshift of 0.2), much smaller than what was used for Figure 13, $`(0.5\text{deg})^2`$, where the bin was kept large in order to contain enough clusters to do statistics. The figure shows a slight correlation between the parameters and the local background density contrast. This is to some degree expected, since the model assumes a uniform background, therefore, a higher background density results in a higher estimation of richness $`\mathrm{\Lambda }_{cl}`$, and vice versa. The redshift dependence shown in the lower panel of Figure 15 is due to the way in which the evaluation of the two parameters $`z`$, $`\mathrm{\Lambda }_{cl}`$ are correlated. There are two effects operating. If the redshift is overestimated, then the angular extent of the cluster ($`r_{\text{max}}`$), within which $`\mathrm{\Lambda }_{cl}`$ is calculated, is underestimated, resulting in a underestimated value of $`\mathrm{\Lambda }_{cl}`$. However, if the redshift is overestimated, the value of $`\mathrm{\Lambda }_{cl}`$ is over-compensated for the loss of galaxies fainter than the survey limit, but this effect is only significant near the faint limit of the survey (at the high-$`z`$ limit). Since the majority of recovered clusters are those at lower redshifts, the results exhibit an anti-correlation between the evaluated $`z`$ and $`\mathrm{\Lambda }_{cl}`$. Figure 15 also demonstrates that most of the outliers in the parameter evaluations are high redshift clusters ($`z>0.3`$), shown as circles in the scatter plot.

## 5 Discussion and Summary

We have presented a comparison of three cluster finding algorithms which are being used to define a cluster catalog from commissioning data of the Sloan Digital Sky Survey. The three algorithms are two Matched Filter Algorithms (MF and AMF), and the Voronoi Tessellation Technique (VTT) which is introduced in its current form in this paper. By applying both Matched Filters on the same galaxy distribution, we have found that the MF is more efficient in locating the clusters, whereas the AMF evaluates the cluster parameters more accurately. This has motivated us to put forward a hybrid method (HMF) which uses the MF to select clusters and the AMF to evaluate the redshifts and the richnesses for those clusters.

The MF is more efficient in selecting clusters than the AMF is because its thresholding method is redshift dependent. The AMF locates peaks in redshift space first, and selects candidates with a signal above a prescribed threshold regardless of the redshift. On the other hand, the MF selects cluster candidates from peaks in the likelihood map of each assumed redshift; those lying above a threshold that is newly determined from each map. As similar clusters at different redshifts will have very different signals (weakening as one goes to higher $`z`$), it is not surprising to find that a redshift-dependent cut results in a better performance. Hence in future work, the AMF should be modified to adopt the peak selection procedure of the MF for better performance, in order to make use of the further advantages in the AMF (e.g., using three dimensional information and a parallelization scheme; Kepner & Kim 2000).

In the VTT method we have applied a filter in color-magnitude space to select galaxies that are most likely members of clusters at a certain redshift. This greatly enhances their contrast relative to the background. This idea is in principle yet another example of a matched filter; this time in color-magnitude space, although it is not a maximum likelihood method like the MF or the AMF. There are a number of existing algorithms that efficiently use this color-magnitude relation as a filter, but are more restrictive: the red-sequence method of Yee & Gladders (1999; Gladders & Yee 2000) and the maxBCG technique (also applied to the SDSS) of Annis et al. (2001). While using a restrictive color-magnitude relation could enhance the efficiency quite a bit, it is also likely to suffer from selection biases, such as missing clusters with significant blue populations of galaxies, namely, the Butcher-Oemler (1984) clusters.

Our current C-M filter for the VTT is generous, so that we can focus on the differences from the Matched Filter algorithms, where a specific cluster model is used. These differences will be further discussed in Paper II with real clusters; in the current paper we investigate the performance of this new method only with simulated clusters that exactly follow the spatial and luminosity profile assumed by the Matched Filter. Despite this advantange for the MF, the VTT has a similar selection function and better false positive rate compared to the MF, which points to the power of the technique, and suggests that photometric redshifts (which use color information) will significantly improve the performance of the AMF.

A Monte Carlo test for the HMF and the VTT was carried out with simulated clusters inserted in $`25\text{deg}^2`$ of SDSS background. We found that the HMF shows a larger drop in detection efficiency in the presence of a non-uniform background than does the VTT. This effect may be due not only to the non-uniformity, but also to overlapping foreground and background clusters in the data, which can cause some of the inserted clusters to be overlooked. This effect is roughly 15% for Abell richness class $`0`$ clusters in the redshift range of our interest, if we assume a random distribution of clusters. The VTT uses color information (although generous) that provides stability against such contamination, while the HMF assumes a uniform background and hence is affected more by this non-uniformity. Thus using a proper model of the background as a function of position can further improve the efficiency of the HMF (Lobo et al. 2000).

We have determined appropriate detection thresholds for the final cluster catalog that is to be drawn from the SDSS data itself. These thresholds are: $`\sigma _{cut}=5.2`$ and $`N_{cut}=9`$ for the HMF and the VTT respectively. These values give a reasonable recovery fraction; only 15% of clusters detected in the uniform background case are not recovered with the realistic background, while the lower limit on false detections, being a small fraction of the expected surface density of real systems, is acceptable ($`1\text{deg}^2`$ for the HMF and $`0.56\text{deg}^2`$ for the VTT).

The selection functions for both algorithms were evaluated using these detection thresholds. The performance of both methods are very similar, although the VTT efficiency tends to drop in the intermediate redshift range compared to the HMF, and is slightly better for lower redshift. Both methods are complete for rich clusters ($`\mathrm{\Lambda }_{cl}70`$) up to $`z\stackrel{>}{}0.4`$. We compare our selections functions to those of BNP00, where a similar cluster finding technique was used to find clusters in shallower data. We find that the performances are similar for the very low richness clusters ($`\mathrm{\Lambda }20`$), while the SDSS outperforms BNP00 by going deeper for the richer clusters ($`\mathrm{\Lambda }40`$).

Finally, we have shown that the detection efficiencies of both the HMF and the VTT are nearly independent of the local density of the background, while the estimated redshift and richness from the HMF are only slightly biased as a function of the local background density.

In Paper II, we present the cluster catalog compiled from the SDSS data using these two methods, the HMF and the VTT. With this we will be able to test various properties of the algorithms using real clusters. Photometric redshifts for the SDSS, soon available, will significantly improve the AMF in particular. In future work, it should be very interesting to compare the above methods with other existing techniques, e.g., the maxBCG technique (Annis et al. 2001) that uses far more restrictive color-magnitude information than the VTT, Cut and Enhance (Goto et al. 2001) that uses the proximity in color-space as an enhancement method, or the Expectation Maximization algorithm (Nichol et al. 2000), which does not impose any model constraints at all. Different algorithms show different results mainly for the fainter clusters – poor or distant clusters – the regime where clusters are abundant but our understanding is poor. Their detection will inevitably depend on different aspects of the method that is used. If one were to explore the properties of cluster parameters such as the luminosity function or the density profile, using the Matched Filter that constrains these models a priori will be highly inappropriate; when exploring the density-morphology relation of clusters, one should avoid using any color constraints in their selection. Therefore, there will be no one technique ideal for all aspects of cluster science; each cluster catalog should be accompanied with a proper understanding of the nature of the detection method, so that it can be used for appropriate cluster studies.

The Sloan Digital Sky Survey (SDSS) is a joint project of The University of Chicago, Fermilab, the Institute for Advanced Study, the Japan Participation Group, The Johns Hopkins University, the Max-Planck-Institute for Astronomy (MPIA), the Max-Planck-Institute for Astrophysics (MPA), New Mexico State University, Princeton University, the United States Naval Observatory, and the University of Washington. Apache Point Observatory, site of the SDSS telescopes, is operated by the Astrophysical Research Consortium (ARC). Funding for the project has been provided by the Alfred P. Sloan Foundation, the SDSS member institutions, the National Aeronautics and Space Administration, the National Science Foundation, the U.S. Department of Energy, the Japanese Monbukagakusho, and the Max Planck Society. The SDSS Web site is http://www.sdss.org/. We thank the referee, Chris Collins, as well as Donald P. Schneider and Martin Kerscher for very useful comments on the manuscript. We also thank David Spergel for discussion and help on RSJK’s thesis, of which this paper is a part. RSJK and MAS acknowledge the support of NSF grant AST96-16901 and AST-0071091 and the Princeton University Research Board, RSJK acknowledges the support of grants NSF AST98-02980 and NASA LTSA NAG5-3503.


format: text

text:
# A New Approach to 𝜂' on the Lattice

## 1 Introduction

The masses of light non-singlet mesons have been determined quite accurately in full QCD simulations. The status of the $`\eta ^{}`$ meson, however, is much less satisfactory, since it is affected for various reasons by severe noise problems: (a) the $`\eta ^{}`$-propagator requires the evaluation of fermion loop operators, which so far has been achieved only with stochastic estimator methods; (b) the Zweig-rule forbidden piece, which is a loop-loop correlator, suffers from gauge field fluctuations; (c) from Wick contraction the $`\eta ^{}`$ propagator $`C_\eta ^{}(\mathrm{\Delta }t)`$ is computed as the numerical difference between its connected and disconnected pieces, leading to loss of signal at large time separations $`\mathrm{\Delta }t`$.

The challenge is to fight noise by exploiting the information from the $`\eta ^{}`$ correlator at small time separations. In this contribution, we intend to revisit the issue of the $`\eta ^{}`$ signal by testing the potential of a spectral representation. The question is whether in the current regime of full QCD calculations (we use SESAM and T$`\chi `$L data), low eigenmodes are actually saturating the fermion loop expressions. Another concern is to improve the groundstate projection in the $`\eta ^{}`$-channel, aiming at very early mass plateauing.

## 2 The spectral approach

Our truncated eigenmode approach (TEA)111For details see  applies the Lanczos method on the Hermitian form of Wilson-Dirac matrix, $`M`$: $`Q\gamma _5M`$ and computes the lowest (in terms of their moduli) eigenmodes, $`|\psi _i`$: $`Q|\psi _i=\lambda _i|\psi _i`$.

In general we need the spectral representation of the quark propagator, $`P`$, which implies the eigenfunctions

$$P(x,y)=\underset{i=1}{\overset{V}{}}\frac{1}{\lambda _i}\frac{\gamma _5|\psi _i(x)\psi _i(y)|}{\psi _i|\psi _i}.$$
(1)
We use the 300 lowest-lying eigenmodes, which at our lightest quark mass amounts to about an equal computational effort as the use of 400 stochastic sources. For details of our analysis, see the parameters given in Table (1).

The $`\eta ^{}`$ propagator has the form

$$C_\eta ^{}(\mathrm{\Delta }t)=C_\pi (\mathrm{\Delta }t)2T(\mathrm{\Delta }t),$$
(2)
where the connected piece, $`C_\pi `$, is readily computed by use of linear solvers. The two-loop spectral expression, $`T`$, reads:

$`T(\mathrm{\Delta }t)=`$
$`{\displaystyle \underset{t}{}}{\displaystyle \underset{i,j}{}}{\displaystyle \frac{1}{\lambda _i}}{\displaystyle \frac{\psi _i(t)|\psi _i(t)}{\psi _i|\psi _i}}{\displaystyle \frac{1}{\lambda _j}}{\displaystyle \frac{\psi _j(t^{})|\psi _j(t^{})}{\psi _j|\psi _j}},`$ (3)
with $`t^{}=t+\mathrm{\Delta }t`$.

We propose to project out excited states from the connected piece by the replacement $`C_\pi C_\pi ^g`$, where $`C_\pi ^g`$ stands for the fitted ground state correlation function , and to extract local $`\eta ^{}`$ masses from the combination of ‘pseudodata’

$$\stackrel{~}{C}_\eta ^{}(\mathrm{\Delta }t)=C_\pi ^g(\mathrm{\Delta }t)2T(\mathrm{\Delta }t).$$
(4)
Using Eq. (4) we obtain local mass plateaus that start at very small values of $`\mathrm{\Delta }t`$, as illustrated in Fig. (1), where both TEA and stochastic estimator results  are shown. Note that the statistical accuracy is sufficient to discriminate the flavour non-singlet pseudoscalar mass very well from the singlet one. We find the data from the TEA analysis to be less noisy than the SET one.

Given the very early onset of the $`\eta ^{}`$ mass plateaus, we can perform single cosh fits to the flavour singlet correlator with safe control on the underlying $`t`$-range. The fit ranges in $`t`$ as well as the $`\eta ^{}`$-masses in lattice units for SET and TEA are listed in Table (1). We find perfect agreement between the results from the SET and TEA analyses.

## 3 Chiral extrapolation

We translate the lattice mass values to physical units through the lattice spacing $`a_\rho ^1(\kappa _{light})=2.302(64)\text{GeV}`$ derived from our previous light spectrum analysis . The critical and physical light quark hopping parameters are $`\kappa _c=0.158507(44)`$ and $`\kappa _{light}=0.158462(42)`$, respectively.

We perform the chiral extrapolation along $`\kappa _{sea}=\kappa _{val}`$, assuming either $`m_\eta ^{}`$ or $`m_\eta ^{}^2`$ to be linear in the quark mass. Our fits favour the dependence $`m^2=c+c^{}m_q`$ from the $`\chi ^2/d.o.f.`$ in Table (2).

To compare with the physical $`\eta ^{}`$ mass, we have to transform the mass from the real world of three flavours into our $`N_f=2`$ setting. This can be done by making use of the experimental mass splitting between flavour singlet and non-singlet states

$$M_{0;N_f=3}^2=M_{\eta ^{};N_f=3}^2\left(2M_K^2M_\eta ^2\right)$$
(5)
and the Witten-Veneziano formula222for a recent discussion on the lattice version of this eq., see

$$M_0^2=2N_f\chi /F_\pi ^2.$$
(6)
For $`N_f=2`$ all non-singlet masses are degenerate. Hence the above equations translate to

$$M_{\eta ^{};N_f=2}^2=2/3M_{0;N_f=3}^2+M_\pi ^2.$$
(7)
Inserting the physical mass values on the right hand side of Eq. (7) leads to the two flavour $`\eta ^{}`$ mass

$$M_{\eta ^{};N_f=2}=715\text{MeV}.$$
(8)
By inspection of Table (2) we find that (at our lattice spacing!) the linear (in quark mass) chiral extrapolations of both $`m_\eta ^{}`$ and $`m_\eta ^{}^2`$ are definitely above the pion mass, yet significantly below the two-flavour pseudoexperimental value estimated in Eq. (8).

## 4 Summary and Conclusion

We have presented a new approach to the computation of flavour singlet masses, where the two-loop correlators are estimated by a spectral representation. We have shown that $`𝒪(300)`$ low eigenmodes suffice to saturate the fermion loop contribution to the $`\eta ^{}`$ propagator in the regime of quark masses of the SESAM and T$`\chi `$L settings.

We found that, in both the TEA and SET analyses, most of the excited state contamination to the $`\eta ^{}`$ propagator is related to its connected contribution: by a conventional ground state projection of the latter, one can achieve a strikingly early onset of plateau behaviour in the effective $`\eta ^{}`$ mass. This then leads to 14 % statistical accuracy of $`M_\eta ^{}`$ after chiral extrapolation.

In our setting, TEA and SET are equally cost effective. It is obvious that TEA is bound to win as we shall go to lighter quark masses in the future; there iterative solvers – unlike the Lanczos procedure – will suffer severe convergence problems due to large condition numbers, while eigenmode expansions will converge even better.

It appears worthwhile to investigate the potential of TEA for disconnected matrix elements. Research in this direction is under way.