Clustering of Sloan Digital Sky Survey III Photometric Luminous Galaxies: The Measurement, Systematics and Cosmological Implications

The Sloan Digital Sky Survey (SDSS) surveyed 14,555 square degrees, and delivered over a trillion pixels of imaging data. We present a study of galaxy clustering using 900,000 luminous galaxies with photometric redshifts, spanning between $z=0.45$ and $z=0.65$, constructed from the SDSS using methods described in Ross et al. (2011). This data-set spans 11,000 square degrees and probes a volume of $3h^{-3} \rm{Gpc}^3$, making it the largest volume ever used for galaxy clustering measurements. We present a novel treatment of the observational systematics and its applications to the clustering signals from the data set. In this paper, we measure the angular clustering using an optimal quadratic estimator at 4 redshift slices with an accuracy of ~15% with bin size of delta_l = 10 on scales of the Baryon Acoustic Oscillations (BAO) (at l~40-400). We derive cosmological constraints using the full-shape of the power-spectra. For a flat Lambda CDM model, when combined with Cosmic Microwave Background Wilkinson Microwave Anisotropy Probe 7 (WMAP7) and H_0 constraints from 600 Cepheids observed by HST, we find \Omega_\Lambda = 0.73 +/- 0.019 and H_0 to be 70.5 +/- 1.6 km/s/Mpc. For an open Lambda CDM model, when combined with WMAP7 + HST, we find $\Omega_K = 0.0035 +/- 0.0054, improved over WMAP7+HST alone by 40%. For a wCDM model, when combined with WMAP7+HST+SN, we find w = -1.071 +/- 0.078, and H_0 to be 71.3 +/- 1.7 km/s/Mpc, which is competitive with the latest large scale structure constraints from large spectroscopic surveys such as SDSS Data Release 7 (DR7) (Reid et al. 2010, Percival et al. 2010, Montesano et al. 2011) and WiggleZ (Blake et al. 2011). The SDSS-III Data Release 8 (SDSS-III DR8) Angular Clustering Data allows a wide range of investigations into the cosmological model, cosmic expansion (via BAO), Gaussianity of initial conditions and neutrino masses. (abridged)


INTRODUCTION
The distribution of light in the universe has long been used as a probe into the structure of the universe. Einstein wrote of the distribution of stars as possibly being uniform on average over large enough distances in 1917 when he discussed the structure of the universe. Hubble tested the uniformity of distribution of faint nebulae in 1926. As the structure of the universe unfolds, the distribution of light from objects such as galaxies has remained a powerful cosmological probe (Peebles 1973;Groth 1973;Wang et al. 1999;Hu 1999;Eisenstein et al. 1999).
Smoothed over large scales, we expect galaxy density to have a simple relationship to the underlying matter density; this implies that the clustering of galaxies at large scales is directly related to the clustering of the underlying matter and is thus a sensitive probe of both the initial conditions of the universe and its subsequent evolution. It is therefore not surprising that a large fraction of the effort in observational cosmology has been devoted to measuring the spatial distribution of galaxies, as in the CfA Redshift Survey (Huchra et al. 1983), the APM Galaxy Survey (Maddox et al. 1990), the DEEP survey 30 (Koo 1998), VIMOS-VLT Deep Survey (Le Fèvre et al. 2005), Two-Degree Field Galaxy Redshift Survey (2dFGRS; Cole et al. 2005), the Two Micron All sky Survey (Skrutskie et al. 2006), COSMOS (Scoville et al. 2007), the Canada-France-Hawaii Telescope Legacy Survey 31 (Ilbert et al. 2006), the Galaxy And Mass Assembly survey (Driver et al. 2009), the WiggleZ Survey (Blake et al. 2010), and the Sloan Digital Sky Survey (SDSS; York et al. 2000). By 2008, the SDSS 32 had probed ∼1.5 Gpc 3 with galaxies, while the current SDSS-III (Eisenstein et al. 2011) will finish surveying ∼15 Gpc 3 in 2014. The planned Large Synoptic Survey Telescope (LSST) 33 will observe ∼1000 Gpc 3 of the universe.
Hidden in the ever-increasing volume of the surveyed universe is a wealth of cosmological information that has not been fully exploited. In particular, the large-scale clustering of any mass tracer, usually characterized by its power spectrum, in the universe contains three features that are of significant interest to contemporary cosmologists. The first distinguishing feature is oscillations in the power spectrum caused by acoustic waves in the baryon-photon plasma before hydrogen recombination at z ∼ 1000, called baryon acoustic oscillations (BAOs; Peebles & Yu 1970;Sunyaev & Zeldovich 1970;Bond & Efstathiou 1984;Holtzman 1989;Hu & White 1996;Eisenstein & Hu 1998). The BAO technique has emerged as the new precision cosmology probe, especially in discerning the properties of this unknown dark component of the universe, "dark energy." The BAO was first observed in early 2005 in the SDSS luminous red galaxy (LRG) sample (Eisenstein et al. 2005), the 2dFGRS data (Cole et al. 2005), and in 2006 by using photometric LRGs in 3500 deg 2 of SDSS (Padmanabhan et al. 2007).
Second, the largest scales of the power spectrum can be used to constrain the primordial potential of the universe, thus testing inflation. In particular, Dalal et al. (2008) have pointed out the relationship between the non-Gaussianity of the potential in the early universe (due to various possible inflationary scenarios) and the large-scale power of mass tracer in the universe. Finally, at k ∼ 0.01 h Mpc −1 , the power spectrum turns over from a k 1 slope (for a scale invariant spectrum of initial fluctuations) to a k −3 spectrum, caused by modes that entered the horizon during a radiation-dominated era and were therefore suppressed. The precise position of this turnover is thus determined by the size of the horizon at matter-radiation equality. It corresponds to a physical scale determined by the total matter (Ω M h 2 ) densities and radiation densities (Ω γ h 2 ). In particular, with a large survey such as SDSS, various groups have used the large-scale power spectrum to put stringent constraints on cosmological parameters, most notably Zehavi et al. (2002), Tegmark et al. (2004), Eisenstein et al. (2005), Padmanabhan et al. (2007), Percival et al. (2010), and Reid et al. (2010).
The SDSS has now surveyed 14,555 deg 2 , and with appropriate photometric selection, we can construct a large uniform sample of the photometric LRGs (Ross et al. 2011) and their photometric redshifts, which can be easily calibrated using the acquired spectroscopic redshifts of a uniform sub-sample (∼10%) of the photometric galaxies. This approach allows for the possibility of using both standard rulers (from the turnover scale of power spectrum and also the BAOs) to acquire cosmological constraints.
We make use of this opportunity to derive one of the most accurate measurements of the galaxy angular power spectra achieved to date. We start with the five-band imaging of the SDSS-III DR8 (Aihara et al. 2011;Eisenstein et al. 2011) and photometrically select a sample of LRGs, following the CMASS galaxy selection detailed in White et al. (2011); the details of the construction of the sample and the redshift distribution are described in Ross et al. (2011). We then measure the angular clustering power spectra as a function of redshift with an optimal quadratic estimator, which has been proved to provide the best statistical error bar when the field is Gaussian. The galaxy density field is not Gaussian on small scales, due to nonlinear evolution; however, at relatively large scales, which are the scales we are concerned with here, the field is close to Gaussian. We will discuss this issue in detail in the paper. With such a large volume of data, we realize that the effects of large-scale systematics are not negligible. To gauge and correct for the effects of large-scale systematics, we develop a novel method in correcting the large-scale systematics given that we know the list of possible systematics. We construct maps of various systematics and calculate their cross-correlation with the galaxy density, the systematic auto-correlations, and cross-correlations. We can then correct for these systematics by applying this new method.
The paper is organized as follows: Section 2 describes the construction of the sample; Section 3 presents the theory and measurement of the angular power spectra; Section 4 discusses the various potential systematics and the novel method applied in correcting for the observational systematics. Section 6.1 describes the validation of the cosmological parameter fitting method, and Section 7 summarizes the cosmological constraints themselves. We conclude in Section 8.

Observations
The SDSS (York et al. 2000;Eisenstein et al. 2011) mapped over a quarter of the sky using the dedicated Sloan Foundation 2.5 m telescope located at Apache Point Observatory in New Mexico (Gunn et al. 2006). A drift-scanning mosaic CCD camera (Gunn et al. 1998(Gunn et al. , 2006 imaged the sky in five photometric bandpasses (Fukugita et al. 1996;Smith et al. 2002) to a limiting magnitude of r 22.5. The imaging data were processed through a series of pipelines that perform astrometric calibration (Pier et al. 2003), photometric reduction (Lupton et al. 2001), and photometric calibration . In particular, Baryon Oscillations Spectroscopic Survey (hereafter BOSS), which is a part of SDSS III (Eisenstein et al. 2011;Aihara et al. 2011), has completed an additional 3000 deg 2 of imaging and is now obtaining spectra of a selected subset of 1.5 million galaxies. The targets are assigned to spectroscopic plates (tiles) using an adaptive tiling algorithm based on Blanton et al. (2003), and observed with a pair of fiber-fed spectrographs.
The availability of a large uniform photometric data set prompted the start of this project, and thus a series of papers, starting with the generation of the photometric redshift catalog by Ross et al. (2011), which uses 112,778 of the BOSS spectra as a training sample for the photometric catalog. The photometric redshift catalog contains over 1.6 million objects, and 900,000 of these objects lie within our imaging mask and the selected redshift range (0.45 < z < 0.65). The redshift range is selected so that it is nearly completely independent from the DR7 analysis of LRG clustering using spectroscopy which stops at z ∼ 0.4 Percival et al. 2010); this allows for the possibility of trivial combination of likelihoods. These galaxies are among the most luminous galaxies in the universe and trace a large cosmological volume while having a high enough number density to ensure that shot noise is not a dominant contributor to the clustering variance. The majority of the galaxies have spectral energy distributions (∼85%, see Masters et al. 2011 andC. Maraston 2011, private communication with the BOSS galaxy-evolution group) that are distinctive of old stellar populations.

Defining Luminous Red Galaxies
We make use of the CMASS sample from BOSS, which is defined in White et al. (2011) and Ross et al. (2011), and we write down the criteria here again for convenience: The magnitudes denoted by "cmod" are "cmodel magnitudes" (see White et al. 2011 for more discussions), and the colors are defined with model magnitudes, except for i fiber2 , which is the magnitude in the 2 spectroscopic fiber (Stoughton et al. 2002;Abazajian et al. 2004). Note that we applied i fiber2 < 21.7, although the current BOSS target selection has moved the limit from 21.7 to 21.5. All magnitudes are extinction corrected using the maps of Schlegel et al. (1998).
In addition to constructing galaxy density maps, we created several additional maps that we use to reject regions heavily affected by sample systematics such as poor sky or stellar density, and to make sure our final power spectra are free of systematics. These include (1) a map of the FWHM of the point-spread function in the r band; (2) a map of stellar density (18.0 < r < 18.5 stars); (3) a map of sky brightness in the i-band in nanomaggies 35 /arcsec 2 ; (4) three maps of the color offsets in u−g, g−r, and r−i from Schlafly et al. (2010); and (5) a map of Galactic extinction simply rescaled from the extinction maps from Schlegel et al. (1998).

Angular and Redshift Distributions
To interpret the clustering of any sample, one must characterize the expected distribution of the sample as if it is completely random. This involves understanding both the angular and radial selection function in addition to the expected galaxy density, which is characterized by its mean density.
To characterize the angular window function, we generate the complete angular mask of the survey following the procedures described below. The observed sky is defined as a union of all fields. Determining the window function requires identifying the fields that cover each position on the sky and deciding which of those fields should be considered as primary at that position. There is a unique set of disjoint polygons on the sky defined by all of the field boundaries, which are calculated using the MANGLE package 36 (Hamilton 1993;Hamilton & Tegmark 2004;Swanson et al. 2008), and each field can be divided into multiple polygons. We now must decide which fields are primary for each polygon in the sky; the process is described in Aihara et al. (2011) in detail. Once we determine which fields are primary for all the polygons in the sky, we make a cut on the field observing conditions (SCORE >= 0.6; for more details on SCORE, see Aihara et al. 2011 or the SDSS-III Web site 37 ). We now have a unified MANGLE polygon file that includes all of the fields that are imaged in the entire SDSS footprint, with the correctly assigned primary fields with good observing conditions. We call this the "full imaging mask," as plotted in Figure 1. The color in Figure 1 represents the date at which the imaging data were taken. The striped pattern perpendicular to the scanning direction is easily visible, and we can also see that the north and south Galactic caps are observed at significantly different epochs of the survey. This provides a hint as to what potential observational systematic effects would look like. To create a more restrictive mask that caters toward photometric red galaxies, we proceed to exclude regions where E(B −V ) > 0.08 35 http://data.sdss3.org/datamodel/glossary.html#nanomaggies 36 http://space.mit.edu/∼molly/mangle 37 sdss3.org  . Preliminary imaging mask after applying primary selection cuts such as cuts on seeing and the bright star mask on the full imaging angular mask. (Scranton et al. 2002;Ross et al. 2006;Padmanabhan et al. 2007;Ho et al. 2008), which is almost identical to A r > 0.2, when seeing in the i band exceeds 2. 0 in FWHM, and masking regions around stars in the Tycho astrometric catalog (Høg et al. 2000). We vary the size of the bright star mask in the Tycho catalog. For brighter stars, we have a bigger mask. We include Tycho catalog stars brighter than B = 13, and vary the radius of the mask between B = 6 and B = 11.5. The final angular selection function covers a solid angle of ∼11,000 deg 2 and is shown in Figure 2.
Applying the selection criteria in Section 2.2 to the 14,555 deg 2 of photometric SDSS imaging considered in this paper yields a catalog of approximately 1,500,000 galaxies. Applying the angular selection function as shown in Figure 2 to the Ross et al. (2011) photometric redshift catalog yields a sample of 872,921 objects, 96% of which are believed to be galaxies (3% are stars, and 1% are quasars, according to statistics gathered in the spectroscopic sub-sample; Ross et al. 2011). For every object, photometric redshifts and the probability of being a galaxy were determined using the ANNz Neural Network (Collister & Lahav 2004;Firth et al. 2003). The calibration and accuracy of these data are discussed in detail in Ross et al. (2011). In the range considered in this paper, the redshifts have calibrated errors ∼0.04 at z ∼ 0.45 and ∼0.06 at z ∼ 0.65. We pixelize these galaxies as a weighted (with the probabilities of being a galaxy) number overdensity, δ g = δn/n, onto a HEALPix pixelization (Gorski et al. 1999) of the sphere, with 12,582,912 pixels over the whole sphere (HEALPix resolution 10, nside = 1024), where each pixel covers a solid angle of 11.8 arcmin 2 . These pixelized maps are used directly to compute the angular power spectra  Notes. Bias parameters are deduced from marginalizing over all the other cosmological parameters (and a free shot-noise term) from combining WMAP 7 + HST + DR8 angular power spectra likelihood using only 30 < < 150 multipoles. The first and last bins are dropped from here on due to the small number of galaxies in those bins.
using the optimal quadratic estimator. The optimal quadratic estimator does not downsample input pixelized maps, rather, it computes the covariance matrix directly from these pixelized maps, and this will be discussed further in Section 3.4. The sample is divided into six photometric redshift slices of thickness Δz = 0.05 starting at z = 0.4 for the CMASS sample (CMASS 0 through CMASS 5; see Table 1 for details), and the underlying redshift distributions for each slice are calculated using the BOSS spectroscopic redshift of the same sample. The redshift distribution of the sample is plotted in Figure 3. We can see that although the majority of the objects in one photometric redshift bin is in their corresponding true redshift bin, a significant fraction of them fall into neighboring bins. The comparisons of these photometric redshifts to the spectroscopic redshifts (obtained via SDSS III spectra) are plotted in Figure 4, while properties of the different slices are summarized in Table 1. We see that the numbers of galaxies in both the first and the last bins are significantly smaller than the others, and therefore we decide to drop these two bins from our analysis. This decision is also facilitated by the fact that we wish to have a nearly independent sample from the Reid et al. (2010) and Percival et al. (2010) LRG clustering analysis, thus allowing for a simple combination of their likelihoods in the cosmological parameter analysis.

Sample Systematics
There are a number of potential systematic effects in photometric samples that contaminate clustering: stellar contamination and obscuration, seeing variations, sky brightness variations, extinction, and color offsets (such as those described in Schlafly et al. 2010). Ross et al. (2011) had extensive discussions concerning these potential systematics; we will concentrate on the particular effects from various systematics on the angular power spectra in the range of scales that affects our science analysis.
The above cuts remove only the parts of the sky that are significantly affected by extinction and seeing variations. With such a large sky coverage, an accurate determination of the angular power spectra of the large-scale tracer is only possible through a thorough understanding of the systematics. However, if we only retain parts of the sky that have the minimum systematic effects, then we must remove most of the coverage, as we have demonstrated in Ross et al. (2011). Therefore, we developed a novel way of dealing with residual sample systematics which we will discuss in detail in Section 5.

THE ANGULAR POWER SPECTRUM
As was noted in the Introduction, the angular power spectrum contains information about both the growth and the expansion of the universe, and thus the shape of the power spectrum, through two standard rulers of the universe: the BAOs and the matter-radiation equality turnover scale. In this section, we will describe both the theory and the computation of the angular power spectrum.

From Galaxy Distributions to Angular Power Spectrum
The intrinsic angular galaxy fluctuations are given by where b(z) is an assumed scale-independent bias factor relating the galaxy overdensity to the mass overdensity, i.e., δ g = b δ, N (χ (z)) is the normalized selection function, and χ (z) is the comoving distance to redshift z. We focus on the auto-power spectrum of the galaxies: where P (k) = P (k, z = 0) is the matter power spectrum today as a function of the wave number k, and the function [g] is The Limber approximation, which is quite accurate when is not too small ( ∼ > 10), can be obtained from Equation (4) by setting P (k) = P (k = ( + 1/2)/χ (z)) and using the asymptotic formula (2/π ) k 2 dkj (kχ)j (kχ ) = (1/χ 2 )δ(χ − χ ) (when 1). We find that the substitution k = ( + 1/2)/χ (z) is a better approximation to the exact expressions than k = /χ (z). Note that j l (x) is the lth-order spherical Bessel function. On large scales where the mass fluctuation δ 1, the perturbations grow according to linear theory δ(k, z) = δ(k, 0)D(z)/D(0).
For auto-correlation, applying the Limber approximation will change Equation (4) to the following: For cross-correlation between two different large-scale structure samples (be it different selection functions, redshift distributions, different biases), we can write the cross-correlation as follows: where g can have different biases, redshift dependence, etc. We have not yet distinguished between the galaxy and the matter angular power spectrum. Throughout this paper, we simply assume where C g ( ) and C( ) are the galaxy and matter angular power spectra, b g is the linear galaxy bias, N shot is a constant shot-noise term which is estimated by the optimal quadratic estimator, and a is a constant term that is fitted as a freely floating parameter. This is a good approximation on large scales, but breaks down on smaller scales; we defer a discussion of its regime of validity, as well as the nonlinear evolution of the power spectrum, to a later section of this paper (Section 3.3). Throughout the paper, we adopt this linear redshift independent (within our redshift slice) bias model with a constant shotnoise term. The bias and the shot-noise term of the galaxy sample for the various redshift slices are fit as extra parameters in Cosmological Monte Carlo (COSMOMC; Lewis & Bridle 2002) chains to ensure that we do not bias our cosmological models via fixing of any particular pre-computed bias.

Redshift-space Distortions
The position of observed galaxies can be inferred from their redshift, and hence the peculiar velocity along the line of sight can, in principle, affect our angular power spectrum. So far, we have neglected the effect of the peculiar velocity, i.e., the redshift-space distortion (RSD) effect on the angular power spectrum. In the three-dimensional (3D) redshift-space power spectrum measured with spectroscopic surveys, the modeling of RSD is still challenging due to the fact that the mapping process from real to redshift space is nonlinear in terms of peculiar velocity. For recent efforts, see, for example, Scoccimarro (2004), Taruya et al. (2010), Seljak &McDonald (2011). It is comparatively easy to model the RSD effect on the angular power spectrum because the RSD information along the line of sight is projected out in the angular clustering. Padmanabhan et al. (2007) formulated the RSD for the angular power spectrum at the linear level, and showed that the linear RSD effect can be seen only at large scales ( < 20). However, we could imagine that if we select thin redshift slices, the nonlinear RSD effect may not be projected out and may become non-negligible at small scales. S. Saito et al. (2012, in preparation) show that such nonlinearities become important only in the case when σ z < 0.01 at > 500, but this is not the case here.
Here, we include the linear RSD effect following Padmanabhan et al. (2007). To be complete, let us review some of the important details from Padmanabhan et al. (2007): where we have now written the normalized selection function as a function of redshift-space distance, s = χ + v ·θ with the peculiar velocity component, v. Assuming the peculiar velocities are small compared with the thickness of the redshift slice, we Taylor expand the selection function to linear order, Substituting this expression into Equation (9), we express separately the two-dimensional (2D) galaxy density field in two terms, g = g 0 + g r , where g 0 is the term discussed in the previous section, while g r is the linear RSD correction. With the help of the linear continuity equation, we have the Legendre coefficient as The component is given by where β is the growth parameter defined by β ≡ d ln D/d ln a/b g , and j is the derivative of the spherical Bessel function with respect to its argument. We can then apply the fact that C ≡ g g * , and calculate the redshift-space distorted angular power spectra.

Nonlinearities
Nonlinearities in the power spectrum are caused by the nonlinear evolution of components of the universe, especially the late time evolution of matter and baryons. To capture the full extent of the nonlinearities, with a lack of a fullfledged nonlinear evolution theory, one will need to simulate the evolution of most if not all of the components of the universe. Extensive research and discussion have been carried out on multiple fronts (Sánchez et al. 2008(Sánchez et al. , 2009), whether it is by perturbation theory (Carlson et al. 2009), dark matter simulations (Hamaus et al. 2010;Heitmann et al. 2009), or fitting functions suggested by dark matter simulations (Smith et al. 2003). Historically, there are a few ways to deal with nonlinearities in utilizing the power spectrum to constrain cosmology, such as comparing the nonlinear power spectrum to the linear power spectrum (usually for a specific cosmological model) and keeping only scales that are believed to be linear Padmanabhan et al. 2007); utilizing the halo occupation model to convert a galaxy power spectrum into a halo power spectrum, which can be easily compared to a halo power spectra from dark matter simulations ); or using a variety of fitting functions developed by Carlson et al. (2009) to fit its observed galaxy power spectra (Blake et al. 2010). Our project both benefits and suffers from the fact that it is a photometric survey. On the one hand, its BAO signal is smeared, as we do not have accurate redshifts; on the other hand, the integration along lines of sight ameliorates the nonlinearities that would have been considerably stronger. Therefore, traditionally, angular power spectra analysis usually only applies a simple cut on the angular scale that roughly corresponds to k = 0.1k Mpc −1 (Padmanabhan et al. 2007). In this paper, we take a small step forward in terms of the nonlinearity treatment of the overall shape of the angular power spectrum, and also adopt a similar treatment as in Eisenstein et al. (2007) and Blake et al. (2010) for the nonlinear treatment on the BAO scales.

Nonlinear Effects on the Overall Shape of the Power Spectrum
There is an extensive body of literature discussing how one can model the linearities of a 3D power spectrum over a large range of scales (Sánchez et al. 2008;Carlson et al. 2009;Hamaus et al. 2010). This paper does not intend to address the issue of fully modeling the nonlinearities in a 3D power spectrum; we do, however, take a simple model that happens to perform quite satisfactorily for the 2D angular power spectrum. We adopt the simple linear redshift-independent biasing model (with shot noise subtracted for every single angular power spectra). Therefore, in addition to the cosmological parameters that are of interest for each model, we include three extra parameters for each redshift slice (b, N shot , and a) as shown in Equation (8).
We test the sufficiency of this model in multiple ways. We test this model by fitting only 2 < < 150 and 2 < < 200 using simulated CMASS mocks (as is discussed in Section 3.5). We compute optimally quadratic estimated power spectra of simulated data (a total of 160 realizations from 20 independent simulation boxes, with 8 lines of sight each), and then we compute 8 averaged (over 20 independent simulations) power spectra, and combine it with a pseudo-WMAP7 likelihood (which has the covariances of WMAP7 likelihood, but with cosmological parameters centered on the input parameters of the CMASS mocks). We find that when using the above-mentioned model for the averaged power spectra, in combination with pseudo-WMAP7, we recover all input cosmological parameters of the CMASS mocks for all eight averaged power spectra to within 1.5σ . We conclude that a spread over 1.5σ is reasonable. The bias parameters recovered are also similar to the input bias of the CMASS mocks as described in White et al. (2011). We therefore conclude that this model is accurate in recovering cosmological parameters when used in the range of angular scales as specified above. To justify our choice of scale for fitting our cosmological parameters, and the model we adopted, we show how well a simple model such as b 2 P halofit (k) + 1/n + a can fit fairly well up to k = 0.2 h Mpc −1 . Top panel shows the nonlinear power spectrum of halos in cosmological simulations (dots with error bars), the model being considered here is b 2 P halofit (k) + 1/n + a (solid line), and the dashed line shows what happens if we only use b 2 P NL (k) instead. Bottom panel shows a/P hh (k) where a = P halofit (k) + 1/n, we can see that the ratio is fairly consistent with 0.0 for a large range of mass, and it starts to deviate from 0.0 starting at k = 0.1 Mpc h −1 . Therefore, we find it prudent to include an extra parameter in our formalism, as we do with some modes where k is larger than k = 0.1 Mpc h −1 . The different color lines correspond to different halo mass ranges. The largest halos are those with highest bias, which also gives the largest deviations from the model. Note that each mass bin has a number density of n = 3.7 × 10 −5 (H/Mpc) 3 .
We further test this model by comparing this model with Hamaus et al. (2010); we found that our simple method fits the nonlinear power spectrum derived from cosmological simulations quite well even up to k = 0.2 h Mpc −1 . In Figure 5, we plot the nonlinear power spectrum from numerical simulations of halos (points with error bars), and while the solid lines are the power spectrum of various halo mass bins calculated using our simple model b 2 P non−lin (k)+1/n, the model fits the nonlinear power spectrum quite well over a significant range in k even when we have not yet added the additional constant term a. The dashed lines show the results without the shot-noise term for various halo mass bins. Our model of the nonlinear power spectrum is based on HALOFIT (Smith et al. 2003), so in order to not confuse the reader, we will call P non−lin (k) by P halofit (k). The lower panel shows the ratios between a ≡ P hh (k) − (b 2 P halofit (k) + 1/n) and P hh (k) plotted as lines. The nonlinear bias is fairly well fit by our simple model even if we do not include the extra constant bias term. We decided to include the extra constant term to help remove the residual difference between P hh (k) and b 2 P halofit (k) + 1/n.

Nonlinear Effects on the BAO
We test the effect on our results of nonlinear evolution on the smearing of the BAO feature by assuming that the nonlinear matter power spectrum follows the expression in Eisenstein et al. (2007): where Σ nl = G(z)/G(0) × 7.527 h −1 Mpc, P wiggle (k) is the linear theory power spectrum (which includes the BAO), and P no−wiggle (k) is a smooth power spectrum, with the same shape Figure 6. We applied two different methods in calculating the power spectrum (including the BAO) with nonlinear effects taken into account, and find that it makes very little difference in the cosmological parameter constraints.
as P wiggle (k) but without any baryonic oscillations, which is computed using the approximation described in Eisenstein & Hu (1998). Both the wiggle and the no-wiggle part have been computed in linear theory; we then added the corresponding nonlinear ratios as a function of the scale to both of them. The nonlinear ratios as a function of the scale follow Blake et al. (2010). In short, we calculate the nonlinear enhancement of power using the input no-wiggles reference spectrum rather than the full linear model including the baryon oscillations. For more details, please see Blake et al. (2010). This approach significantly enhances the power on small scales. We find that the results are not very sensitive to the exact value of Σ nl , provided that it is in the range of 5.527-9.527 h −1 Mpc (Eisenstein et al. 2007). In principle, Σ nl is cosmology dependent, and thus can change our cosmological constraint if it is kept as a free parameter. We have therefore examined our constraints on cosmological parameters using different Σ nl . We test this issue by fitting the full set of cosmological parameters using the MCMC fitting method with COSMOMC, with Σ nl set to 2 h Mpc −1 higher and lower than its currently chosen value (7.527), and find that when we fit for a ΛCDM model in combination with WMAP7, there is less than a 5% change in the error for any of the parameters. The addition of the nonlinear ratios is quite important, not only because the power on small scales in the angular power spectrum at high multipoles is not expected to account for the shot noise due to a finite number of galaxies, but also because the small shift in the BAO wiggles can slightly modify the best-fit shape of the power spectrum, and hence return a different value of Γ ≡ Ω m h. We applied two different methods in calculating the power spectrum (including the BAO) with nonlinear effects taken into account. The first method follows Blake et al. (2010) and the second method is a naive application of HALOFIT (Smith et al. 2003) on the power spectrum computed following Eisenstein et al. (2007). We find that it makes essentially no difference (see Figure 6).

Optimal Estimation of Angular Power Spectrum
The theory behind optimal power spectrum estimation is now well established, so we limit ourselves to details specific to this discussion and refer the reader to the numerous references on the subject (Hamilton 1997;Seljak 1998 and references therein).
We also refer the reader to the Appendix for more specific details that relate to our paper directly.
We start by parameterizing the power spectrum with 20 step functions in l, where the p i are the parameters that determine the power spectrum. We form quadratic combinations of the data, where x is a vector of pixelized galaxy overdensities, C is the covariance matrix of the data, and C i is the derivative of the covariance matrix with respect to p i . The covariance matrix requires a prior power spectrum to account for cosmic variance; we estimate the prior by computing an estimate of the power spectrum with a flat prior and then iterating once. We also construct the Fisher matrix, The power spectrum can then be estimated,p = F −1 q, with covariance matrix F −1 .
We also refer the reader to the Appendix for details more specific to our project.

Tests with Simulations
To test whether the errors estimated by the quadratic estimator employed here are accurate or not and to test the results of our pipeline, we must compute the errors obtained via a series of simulations.
One way to do this is to generate a Gaussian random field using the prior power spectra for each redshift slice to simulate over the entire sphere. We can Poisson distribute galaxies with a probability of (1 + δ)/2 over the survey region, trimmed with the angular selection function. Padmanabhan et al. (2007) have tested this pipeline with the Gaussian random fields simulations; thus, what we need to test here is whether the errors estimated by the quadratic estimator are appropriate, considering that the power spectrum measurement is only a minimum variance measurement when the field is Gaussian, which is not the case here. Given the non-Gaussianity of the field, we need to determine how close we are to minimum variance measurement.
As we would like to simulate our galaxy sample as closely as possible, we employed CMASS mock catalogs from White et al. (2011) to test the accuracy of the optimal quadratic estimator. White et al. (2011) have produced a series of mock catalogs that use the best-fit Halo Occupation Distribution (HOD) models from White et al. (2011) and populate a series of N-body simulations . The majority of the galaxies are central galaxies living in halos of mass 10 13 h −1 M . We generate eight lines of sight from each corner of each of the 20 independent CMASS simulations from White et al. (2011). These mock catalogs are then processed in the same manner as the real data through the quadratic estimator code, and analyzed in the same manner as the real data set. The mock angular power spectra are thus optimally estimated angular power spectra.
We plotted the distribution of the power spectrum from each simulation that is estimated by the quadratic estimator code, and compare these results to the averaged error bar of the simulation (see Figure 7). When comparing the expected error to the distribution of the estimated power spectrum from each simulation and the averaged measured error from each simulation, we conclude that the averaged measured error is a good measure of the expected error. We have plotted the estimated error (red crosses of the middle panel of Figure 7) by examining the variance of the estimated power spectrum from each simulation at each -mode. We have also plotted the averaged measured error (green points), and it is a bit higher than the estimated error from the variance of the simulations (bottom panel). This is probably due to the fact that there are only 20 simulation boxes, with 8 lines of sights overlapping slightly within each box. Therefore, the variance of simulated C l is probably slightly smaller than it should be at all scales, due to the correlations between lines of sight. Regardless, we show that the quadratic estimator code can estimate the errors of the power spectrum in the scale of interest here with reasonably high accuracy. It is important to note that at all scales of interest (to the current paper), the estimated error from the quadratic estimator code is not underestimated.

The Optimally Estimated Angular Power Spectrum
The angular power spectra estimated using the methodology described in Section 3.4 are displayed in Figure 8. In particular, we plot separately the north (Galactic), south, and full angular power spectra of these four redshift bins (CMASS 1-4, from z = 0.45-0.65). We plotted the north and south separately to investigate possible systematic differences due to the long separation of observation time between north and south galactic caps. For the scales of interest (30 < < 150), the north and south are not different enough to prompt separate analyses. Nevertheless, this does not preclude the possibility of systematic differences at the largest scales (at low multipole) in the angular power spectrum. This is only possible because the estimated power in each -bin is not correlated, therefore a systematic difference in one -bin does not affect another.
To test the similarity of north and south region on scales of interest (30 < < 150), we find that all best-fit cosmology parameters (with combined with WMAP7 via MCMC chain using COSMOMC) found by north and south alone, respectively, are consistent with each other. It is interesting to note that the south has a smaller area than the north, and consequently there is less information per -bin, thus the error bars in the south are significantly larger than the north. It will also be discussed later in Section 4 as the systematic treatment presented in this paper will, in principle, correct systematic variations even when the full survey is analyzed in one piece. We can also see the evolution of the angular power spectra over different redshift slices, as it is expected.
By looking at the cross-correlations between two different photometric redshift bins, which has overlap in spectroscopic redshifts would provide additional information on cosmological parameters; while examining cross-correlations between two photometric redshift bins which has no overlap in their actual objects would give us hints on the additional spurious power that are caused by systematics which would be common across different redshift bins even though there were no overlap in volume. As shown in Figure 3, we need to investigate the potential effects of overlapping redshift distributions. We calculate the cross-power of various redshift combinations, and they are shown in Figure 9.
When we examine cross-power across various redshift bins, any difference between the measured power and the expected power (from galaxy auto-correlations in the same redshift range) Figure 7. Top: the estimated (from optimal quadratic estimator) power spectrum from 160 simulated realizations (red crosses) out of 20 independent dark matter only simulations ). The green points show the averaged recovered power spectrum from each of the simulations, while the blue points show the averaged measured error as estimated from the optimal quadratic estimator. We conclude that the averaged measured error from the optimal quadratic estimator is a good measure of the expected error. Middle: the estimated error (red crosses) is based upon the variance of the estimated power spectrum from each simulation at each ell-mode. We have also plotted the averaged measured error (green points), and it is a bit higher than the estimated error from the variance of the simulations (bottom panel). Nevertheless, we show that the quadratic estimator code can estimate the errors of the power spectrum in the scale of interest here with reasonably high accuracy.
can be used as a measure of the effects of systematics. In the top panel, there is significant extra power at the large scale, and also negative correlations (which cannot come from galaxy auto-correlations); therefore, we know that there are significant systematics within CMASS 1. The bottom panel shows that the high-redshift slice CMASS 4 also has substantial effects from systematics.
Finally, to estimate whether it is worth including the crosspower of various redshift slices in the cosmological analysis, we performed a simple Fisher analysis. We calculated Fisher matrices using angular spectra from the four redshift bins (CMASS 1-4), with the redshift distributions given in Figure 3. A standard ΛCDM cosmology is employed to calculate the fiducial spectra. We used the Limber approximation (where the input power spectrum was given by CAMB 38 linear power spectrum and HALOFIT) and ignored RSDs. We employ the standard Gaussian expression for the covariance matrix of the spectra. The shot-noise term was calculated assuming N l = 1/n (withn being the number of galaxies per steradian of the individual bin). Finally, to construct the Fisher matrix, we used the range l = 30-300. The parameter space is given by Ω b , Ω c , Ω ν , Ω Λ , σ 8 , n s , b 1 , b 2 , b 3 , and b 4 (b N refers to biases of galaxy 38 http://camb.info/ sample at redshift slice N). The Fisher matrix is then added to the WMAP7 Fisher matrix and inverted to find the covariance matrix for the parameters. We then consider two cases: (1) using only the auto-spectra as observables and (2) using both autoand cross-spectra as observables. The errors on all parameters improve by less than 5% in going from (1) to (2). We also found that ignoring covariances between different auto-spectra (we do include the covariance between auto-power spectra in the analysis) makes less than 5% difference. This suggests that when we include these covariances in the MCMC, the errors will not change significantly. We therefore adopt a conservative approach where we do not include the cross-power as an extra signal, but we include the bin-to-bin covariance which can, in principle, be double counted due to the overlap of redshift slices.

POTENTIAL SAMPLE SYSTEMATICS
Without accurately addressing known potential systematics on the observed number density of objects in our sample, we cannot claim to understand its expected angular power spectra, nor can we extract cosmological information from it. The treatment of systematics is especially crucial for the overall shape of the power spectrum, since the shape does not deviate much from power laws and has no specific features such as Figure 8. Measured angular power spectrum for the four redshift bins using the methodology described in Section 3.4. We have plotted the full angular power spectra, which takes in the whole sky in the top panel, the north (Galactic) angular power spectra, and the south (Galactic) angular power spectra. Within the range of interest, the north and south angular power spectra are consistent, suggesting that the systematics are at a relatively low level in the scales of interest if they affect the north and south differently. those in BAO. The oscillatory nature of the BAO signal helps it from being contaminated by any systematic signal that does not have oscillatory features. Moreover, most BAO detection methods attempt to minimize any influence from the shape directly (Eisenstein et al. 2007), thus further shielding the BAO technique from any systematic effects. We will propose a novel method of systematic corrections in Section 5 which helps mitigate the effects caused by any systematics in the power spectrum.

Description of Systematics
Here, we consider not only sample systematics, but, in particular, the systematics that may contribute to extra (or deficit) power in the angular scale under consideration. It is especially important to consider the potential systematics in photometric surveys, since we do not have the radial modes and spectral identifications that are available to spectroscopic surveys.

Stellar Contamination and Obscuration
Stars can, in principle, mimic galaxies given the right colors, or give rise to obscuration due to possible foreground subtraction issues due to the presence of a star. As it was pointed out in Ross et al. (2011) that the magnitude range of stars does not change its effect on the galaxy number density, we pick stars of magnitude 18 < r < 18.5 and investigate its influence on the galaxy auto-power spectrum. The auto-power spectrum of the stellar density map is plotted in Figure 10 while the stellar density map is plotted in Fiugre 11. We use the same mask as the CMASS samples, since the stars can only affect the galaxy power spectrum where the two overlap. We calculate the cross-power spectra between the stars and the various CMASS redshift slices and find that there is a significant correlation at several angular scales (see Figure 12), mostly at large scales. In particular, there are strong angular correlations ( < 10 for CMASS 1, < 20 for CMASS 2) between stars and the galaxies at large scales, while we observe that the number density of galaxies is lower when it is closer to a star (as also discussed in Ross et al. 2011). In this paper, we do not include scales that are smaller than < 8 since it is much larger than the scales we are interested in. However, we will discuss the larger scales further in a future publication on primordial non-Gaussianities, as the largest angular scales contain more information concerning primordial non-Gaussianities.
There is an extensive discussion on the stellar contamination in the CMASS catalog in Ross et al. (2011). The fundamental conclusions are that there are two separate effects: (1) stars can be confused as galaxies, thereby contaminating the sample and inducing a positive correlation between the densities of stars and our sample and (2)  . Measured angular cross-correlations for the first three redshift bins with other slices. We do not show repeats of the cross-correlated pairs. When we examine cross-power across various redshift bins, any difference between the measured power and the expected power (from galaxy cross-correlations) can also be used as a measure of the effects of systematics. In the top panel, there is significant extra power at the large scale, and also negative correlations (which cannot come from galaxy auto-correlations), therefore, there are significant systematics within CMASS 1. In the bottom panel, we observe that the high-redshift slice CMASS 4 also has some substantial effects from systematics at large scales. The CMASS 2 and CMASS 3 samples are fairly clean from systematics in the scales of interest. The fitted line is the best-fit model from our best-fit flat ΛCDM model. are highly correlated across bins separating distances, the two effects together impart a slightly negative correlation between the number density of stars and our sample. In our analysis detailed in this paper, the estimated angular power spectra are designed to have minimal correlation across bins; therefore, we can see both positive and negative correlations over different scales ( -bin), as seen in Figure 12. Given that we know stars are likely to contribute to the observed number densities, we can take into account the amount of contamination by using the above discussed technique. Our results are consistent with Ross et al. (2011) even though we do not detect smaller scale correlations between the stars and the galaxies, since the estimator employed in this paper produces estimates of angular band powers that are  minimally correlated with other bins of band powers, while estimates are highly correlated across bins in the analysis of Ross et al. (2011). Therefore, although the correlations between stars and galaxies are concentrated in the largest scales, they appear in smaller angular scales such as those seen in Ross et al. (2011).

Sky Brightness
The sky brightness from SDSS is presented in Figure 13. There are at least two scales on which the sky brightness would affect the expected number of galaxies. The first is the width of the scan of the SDSS camera, the second comes from the fact that the southern cap has a brighter sky, since we observe the south at higher zenith angle, and thus more re-emission from the optically thin lines that were pumped by the sun originally.
We use the auto-power spectrum of the sky brightness (shown in Figure 10), and its cross-correlations with galaxy densities at various redshifts (shown in Figure 14), to estimate the amount of contamination that can come from the sky. We discuss the corrections applied arising from sky brightness in Section 5.

Seeing Variation
Since SDSS uses a ground-based telescope at the Apache Point Observatory, it is expected that the image quality, primarily due to atmosphere seeing, will affect the number of galaxies detected in any part of the sky. To quantify this, we plotted the seeing variations in the sky in Figure 15. There is a striped pattern as different parts of the sky are observed in different nights, which have different atmosphere seeing. We use the autopower spectrum of the seeing variations and its cross-correlation (as shown in Figure 16) with the galaxy density to determine the effects of seeing on the galaxy overdensity clustering power. Since we can see that there are statistically significant but mild cross-correlations between the galaxies and seeing in several of the angular band powers, we correct for the seeing variations as discussed in Section 5.

Extinction
We check for any residual effects on the observed galaxy overdensities due to Galactic extinction by computing the crosscorrelations between the galaxy overdensities and the extinction map (Schlegel et al. 1998; see Figure 17). Since SDSS avoids most heavily extincted areas, we only have a small overlapping area where there is significant extinction and galaxy data. We do not see a statistically significant cross-correlation between the galaxy (except for scale of < 20 of CMASS 1) as shown in Figure 18 and the extinction field, therefore, we conclude that we will drop galactic extinction from the list of possible contributing systematic effects as long as the range of interest in this analysis remains smaller than > 20. Schlafly et al. (2010) reported various color offsets for the SDSS footprint, in particular a north/south offset. As discussed in Schlafly et al. (2010), the photometric offsets can be estimated via two different ways: (1) using the color of stars in the imaging data; and (2) using the stellar spectra to determine spectral classes, and then calculate differences between the observed and expected colors of stars. We adopted the latter method, since it will not be sensitive to the intrinsic variations of stellar properties. However, this approach requires spectroscopy of stars, which is lacking in significant parts of the SDSS southern sky. We still pursue it though, and find that there are no significant detections between the galaxy density map that overlaps with the offset map (which is lacking in southern sky coverage as shown in Figure 19) in Figure 20. We therefore conclude that this is not an important systematic in our sample. However, this is only a statement that at the sky which is overlapped between the galaxy density map and the offset map, there are not significant correlations. We will need more data in terms of the southern sky offsets before we can conclude on the effects of color offsets and galaxy densities.

Magnitude Errors
As mentioned in Ross et al. (2011), the model magnitude errors in the southern cap are larger than the northern cap by ∼10%, which may introduce a possible power excess (or deficit) at the lowest multipole. This issue, however, should not affect any of the other multipoles which are the focus of our paper.

NOVEL TREATMENT OF SYSTEMATICS
Assuming that residual systematics will exist in the best catalog we can construct without a serious loss of sky, how should we handle the remaining systematics? Without any evidence of the possible nonlinear effects of systematics on the observed density fields, we adopt the simplest approach: a linear relationship between the systematics and the observed galaxy density fields.
We start from the following. We transform fields from real space into spherical harmonic space (or -space in particular), so that δ g δ g = C l : where δ t g is the true galaxy density (in the -space), and each δ i is the fluctuation of the map of the ith systematic, while i characterizes how much the ith systematic contributes. With the lack of a better model, we assume a linear relationship between the systematics and the observed galaxy number overdensity, but in principle this model can be modified to include higher order contamination due to the systematics.
For a simple demonstration, we assume that we have only two systematics contributing to the observed galaxy density, so that i = 2. Assuming that the true galaxy density is unrelated to any of our systematics, we have the following: δ o g δ o g = δ t g δ t g + 2 1 δ 1 δ 1 + 2 2 δ 2 δ 2 + 2 1 2 δ 1 δ 2 . (18) Furthermore, we have all the cross-correlations between the systematic and the observed galaxy density map: Since we calculate all the auto-and cross-correlations of all the systematics (on top of the cross-correlations between the systematics and the observed galaxy density), we can solve for 1 and 2 (and they will be functions of ).
In Figure 21, we show the cross-correlations of all the systematics (which contaminate the observed galaxy field), we find that the correlations across different systematics are far from zero, and we must include the cross-correlations among systematics in our model.
For simplicity of demonstration, we show the result of applying the correction from only one systematic (stars) in Figure 22. Since we only include of one systematic, the correction depends only on the auto-power of the stars and the cross-correlation between the stars and the observed density field. Although the star is one of our dominant systematics, its effects on scales of interest ( > 30) are quite minimal. This implies that with the appropriate estimator (which does not correlate powers in various scales), the effects from systematics can be corrected relatively easily. This result further encourages us in terms of the cosmological constraining power that can be harnessed from future imaging surveys that will go deeper and wider.
Finally, we calculate the final systematic-corrected power by including all three systematics that are found to have significant correlations with the observed galaxy density field. We include the auto-power (of both systematics and galaxy density fields) and all cross-powers (among systematics and galaxy fields) at each angular scale. The final corrected power is shown in Figure 23.
Nevertheless, the optimal quadratic estimator produces optimal errors and unbiased measurement only when the field is Gaussian. In the case of highly non-Gaussian fields, the estimate is still unbiased, but the error is not optimal and may not be accurate. We can only test the validity of the estimated error by the quadratic estimator if we simulate a large number of systematic mocks, and carry out variance tests such as the ones carried out for the galaxy density fields.
While we understand the construction of mock galaxy catalogs, we do not fully understand how to construct mock systematic fields. We lack the "theory" of systematics fields (except probably the stellar density map). Therefore, it is unlikely that we can achieve optimal error on the systematic corrections within the scope of this analysis. The estimated values are unbiased, but the error can be overestimated or maybe incorrect (Hamilton 1997). Nonetheless, when the systematic corrections are small, the uncertainty related to the correction cannot be larger than the correction itself. Therefore, we conclude that the most conservative way would be to include only power from multipoles that have relatively small corrections. Lacking a better model for the systematics, we adopt the following simplistic model of estimating the covariance of the systematic-corrected power spectra for multipole bins which require small corrections. We assume Gaussianity for the fields involved, and thus use the following relationship: We modify the above equation by adding the "correctional power" due to systematics: This is for each δ = 1, so we take into account the fact that the δ is not 1 in all of our bins. The quantity ΔC j ( ) is the correctional power contributed by systematics. This method assumes the Gaussianity of the fields, which is not a satisfactory  Schlegel et al. (1998). By comparing with the full mask of the sky, we can see that for regions of maximum extinction, the Galactic plane is completely avoided. assumption; the optimal quadratic estimator can in principle project out powers that are understood, such as those timedependent systematics that can be projected out in the cosmic microwave background (CMB) map making. Nonetheless, fully modeling the systematics and then projecting them out using the optimal quadratic estimator is a much larger undertaking, which will be left to future work.
We also show for completeness purposes the power spectra of various redshift slices before and after the corrections in Figure 24. The reduced chi-square of the fit shown in Figure 24 is 1.14. The fit is computed from the best-fit model using CMB+HST+DR8.
We also compute the correlation function of our systematiccorrected power spectra. One of the systematic correlation correction methods employed in Ross et al. (2011) follows our paper, 39 and thus it is not surprising that we have achieved the same systematically corrected correlation function, even though the computation of the correlation function is independent. Figure 25 shows how our computed correlation function from our systematically corrected optimally estimated angular power spectra is completely consistent (to within 1.5σ ) with the measurement of Ross et al. (2011). We would also note that our correlation function, shown in black lines in Figure 25 (our w(θ ) 39 Even though this analysis is submitted at a later date, the Ross et al. (2011) paper is part of the DR8 clustering project in SDSS III, and thus Ross et al. (2011) have applied our method described here. Figure 18. Cross-correlations between galaxies of various redshift slices with the Galactic extinction map. Since there are not significant correlations between the galaxies overdensities and the Galactic extinction map, we can drop the extinction from the list of potential systematic effects for CMASS 1, CMASS 3, and CMASS 4. The CMASS 2 sample has a significant contribution at 1 multipole, but for reasons that we will discuss in later sections we will not include multipoles at < 30, so this one multipole at CMASS 2 would not affect our cosmological analysis. Figure 19. Color offsets in r−i calculated from sources provided by Schlafly et al. (2010). These are derived from the spectra of the stars, and thus do not include an intrinsic variation due to the metallicity gradient of stars within our Galaxy. calculated using the angular power spectra), has no significant large-scale power, which suggests no exotic inflation scenario, or residual systematics.

Method
As described in Section 3.3, we adopt the simple linear redshift-independent biasing model (with shot noise subtracted for every single angular power spectra). Therefore, in addition to the cosmological parameters that are of interest for each model, we include three extra parameters for each redshift slice (b, N shot , and a) as shown in Equation (8).
We have estimated that the nonlinear RSD effects are minimal in our case (S. Saito et al. 2012, in preparation), therefore we include the full linear RSD following Padmanabhan et al. (2007) as discussed in Section 3.2. However, calculating the full linear RSDs requires significant time, and it is different from the Limber approximation at l < 30; therefore, we made a decision to employ the Limber approximation for multipole ranges at l 30, and employ the full linear RSD calculation only at l < 30.
The measured band powers from the quadratic estimator have contributions from a range of wave numbers, even though they are highly concentrated in their own -bin. The quadratic estimator is designed to compute nearly anti-correlated power spectra (Padmanabhan et al. 2003) across different multipole bins, but it still has a very small (<∼5%) contribution from other multipole bins. We take this effect into account by convolving the theory power spectra with the window function before calculating the likelihood by (d − t) T C −1 (d − t), where d represents the measured power spectra, t represents the theory power spectra convolved with the window function, and C is the covariance across different bands and redshifts as output from the quadratic estimator.
The fitting of all of the cosmological parameters is done through MCMC with COSMOMC (Lewis & Bridle 2002). As was discussed earlier, we do not know the accuracy of the error estimated using a highly non-Gaussian field with a quadratic estimator. We know that it would not be biased, but the error bar can be significantly misestimated. When the systematic corrections are small, we can assume that the error involved with the correction cannot be significant. However, at the largest angular scales, some of the systematic corrections are large enough that only proper error propagation that involves a significant undertaking in the full modeling of systematics would provide sufficient accuracy on the error of the correction. There are also concerns about the validity of the assumption of linear effects of systematics on the galaxy density field, especially when the corrections are large. We do not have reason to believe that the effects of systematics on the galaxy density field are linear or nonlinear. Therefore, we decided to avoid the scales that involve large systematic corrections, which are concentrated at the large scales, thus setting the start of the multipole range from l = 30.
Since our data set is a photometric sample, the nonlinear effects are mitigated by the fact that when we examine the data set, it is already integrated along all of the lines of sight, thus decreasing its nonlinearities. We therefore apply the HALOFIT routine for computing the nonlinear power spectrum, and limit ourselves to multipoles that are at relatively large scales, while nonlinearities remain a small effect. We choose the multipole range by simply testing a variety of ranges using our Figure 20. Cross-correlations between the galaxy overdensities and the color offsets in g−r (top two panels) and r−i (bottom two panels). We can see that there is no significant correlation between the galaxies and the color offsets, removing these offsets from our list of potential systematics.

Notes.
We also show the k max corresponding to the max considered for each redshift slices. z mid is the midpoint of the redshift interval.
simulations described in Section 3.5. In order for the MCMC chain of the simulation to converge quickly, we average the simulations across different simulation boxes (so that they are not correlated). We tested a large range of multipole ranges: 1 < l < 600, 20 < l < 150, 20 < l < 300, and 20 < l < 200, for example, and found 20 < l < 200 and 20 < l < 150 to return results within 1σ of the input parameters. It is also of interest to show the corresponding k-limit for the various ranges, since the range in k may provide an easier reference for the nonlinearities of interest. We list the corresponding k range for the range we use for the redshift bins considered in our analysis in Table 2.
Combining both the low and the high l limit, we conclude that 30 < l < 150 and 30 < l < 200 are both conservative choices for the fitting of cosmological parameters.

Cosmological Parameter Fitting Method Validation
In this section, we present some of the tests used to check the C l -likelihood routine for COSMOMC. To perform such a test, we need an angular power spectrum whose input cosmology (i.e., the cosmological model with which it was generated, modulo cosmic variance) is completely known. Therefore, the mock angular power spectrum described in Section 3.4, being derived from a cosmological simulation with initial conditions given by a known set of values for the cosmological parameters, provides an excellent testbed for our fitting routine.
The most straightforward test is to fit each individual angular power spectrum from the mocks and check that every one of them (out of the 160 available) returns the input cosmology. However, running 160 MCMC chains with only one mock power spectra each is computationally intensive, especially since each angular power spectrum has the power of ∼0.5-1 actual redshift slice from the data, thus it will take significant time for the chains to converge if they converge at all.
We therefore need to combine these mock C l with an additional data set, and the most obvious choice is the CMB data from WMAP7 (Larson et al. 2011). The WMAP7 best-fit parameters, however, are not exactly the same as the simulation input  parameters, and thus we replace the standard CMB likelihood by a much simpler one, in which we compute the value of χ 2 from the actual covariance matrix from WMAP7, but not using the actual parameters themselves. To fully validate the fitting method, we need to use mocks that show a signal-to-noise ratio similar to those observed in the data. We combine individual mocks by using different simulations (and not different lines of sight in the same simulation).

Building Covariance of Mocks
We compare the Gaussian covariance matrix of the power spectrum from OQE (i.e., "OQE covariance matrix") with  For the first plot, we show the systematic-corrected power spectra and their corresponding best-fit model (from flat ΛCDM model) when we include all three dominant systematics (stars, sky, and seeing). For the last plot, we included the measured (as described in Equation (18)) for all the dominant systematics. the dispersions among 160 mocks (i.e., "N-body covariance matrix"). Note that the 160 mocks are not strictly independent from each other, as different lines of sight share a small but nonzero amount of volume. To exclude the artificial covariance between different lines of sight, we derive the covariance matrix of the 20 independent mocks per each line of sight: we then average the 8 covariance matrices for 8 lines of sight. As a comparison, a straightforward dispersion among 160 mocks gives an almost identical result, implying that the different lines of sight share very little volume.
In the upper panel of Figure 26, the red points are the square roots of the diagonal elements of the OQE covariance matrix and the black squares are from the N-body covariance matrix. The diagonal elements of the OQE covariance matrix can be analytically calculated based on the smooth fit to the measured power spectrum and the number of independent modes for each wave number band assuming Gaussianity, if the matrix is diagonal. However, the OQE covariance matrix includes the effect of the window function due to the survey geometry, and the covariance matrix therefore is not strictly diagonal and has a small anti-correlation between neighboring bins. Indeed, we find that there is a small deviation between the OQE covariance matrix and a naive Gaussian error calculation without accounting for the window function. The difference is expected since the naive Gaussian error calculation does not include the effects of the actual survey geometry. The black dashed lines in the figure are the theoretical, expected errors derived based on Gaussianity: we have rescaled it with an empirical boost factor of 1.1 to better match the observed dispersion. The dispersions between mocks are systematically lower than the OQE on large scales but appear to lie between the OQE expectation and the boosted Gaussian approximation.
The lower panel shows the off-diagonal elements of the N-body covariance matrix in comparison to the OQE covariance, for a slice at l = 185. We observe fluctuations up to 20% in the measured off-diagonal terms but find no obvious indication that it disagrees with the OQE covariance matrix. Therefore, we conclude that the OQE covariance matrix based on the Gaussian assumption does not underestimate the true error of the 2D projection of the nonlinear galaxy field.
For the real data, we use the covariance matrix from the OQE for the auto-power for each redshift bin. The upper panel of Figure 27 shows the square roots of the diagonal elements of the OQE covariance matrix of the real data (open circles) in comparison to the prediction based on the Gaussian prediction (after boosted by 1.1: solid squares) for all four redshift bins. The agreement is even better than the mock case, and it is probably due to the larger survey area of the real data that decreases the cross-correlation between different bins that our simple Gaussian approximation cannot access. The lower Figure 25. Correlation function computed from systematically corrected angular power spectra (black lines) is compared to the measured and systematically corrected correlation function (red squares with error bars) in Ross et al. (2011). The correlation function in Ross et al. (2011) also is systematically corrected following our method detailed in this paper. Following the black lines, we can see hints of the BAO at nearly all redshift bins, and that there is no significant large-scale power seen by other previous analyses such as Thomas et al. (2011) prior to our DR8 analysis. panel shows the cross-correlation between different bins for four different slices of the OQE covariance matrix. We overplot covariance matrices for CMASS1 (black), CMASS2 (red), CMASS3 (blue), and CMASS4 (magenta). Note that the covariance structure is identical for the four redshift bins, which is reasonable as the four redshift bins are subject to the same mask. The same structure should apply to the covariance between different redshift bins as well. We therefore build the cross-covariance between different redshift bins by combining the diagonal elements from the Gaussian assumption and the covariance structure in the right panel of Figure 27. The diagonal elements are constructed using smooth fits to the measured auto and the cross-power spectra of and between redshift bins, and boosted by 1.1 based on the results of the auto-power spectra: where i and j indicate a redshift slice, f sky is the fraction of the sky, N mode is the number of wave modes within the band, and C ij ( ) is a smooth fit to the auto-or cross-power spectrum. We include the shot-noise contribution to C ij ( ) in the case of the auto-power spectra (i.e., i = j ). Cov G ii,jj is the covariance between the auto-power spectrum C i,i and C j,j . The factor a fac is the empirical factor of 1.1 that we introduce to match the OQE covariance matrix and Equation (23). We use this equation to build the covariance between different redshift slices, while using the OQE covariance matrix for the covariance within the redshift slice.

Mock Test Results
With the combined average of 20 spectra each in combination with the pseudo-WMAP7, we find that the above model recovers all input cosmological parameters of the CMASS mocks for all eight averaged power spectra to within 1.5σ . The recovered bias parameters are also very similar to the input bias of the CMASS mocks as described in White et al. (2011). We therefore conclude that this model accurately recovers cosmological parameters when used in the range of angular scale specified above.

Constraints on Cosmological Models
The angular clustering measurement can be used to constrain the cosmological model in several different ways: through standard rulers such as the matter-radiation turnover scale, the BAOs, or through large-scale power that would constrain primordial non-Gaussianities (Dalal et al. 2008;Slosar et al. 2008). In a companion paper (Seo et al. 2012), we examine only BAOs, and remove any contribution from the overall shape of the power spectrum: In this paper, we include the overall shape of the power spectrum and parts of the BAOs to derive constraints on cosmological models. There is a companion publication on the neutrino mass constraints using the same angular power spectrum (de Putter et al. 2012).  for the four redshift bins. The two agree very well. The bottom panel shows the uniformity of the off-diagonal structure of the four redshift bins for four slices of the covariance matrix. The uniformity arises from the common mask. We therefore use a Gaussian assumption and such a uniform off-diagonal structure to build a cross-covariance between different redshift bins. Here, we choose to consider a variety of cosmological models, although we do not intend to exhaust all possibilities. We include several other data sets to help break cosmological degeneracies, such as WMAP7 (Larson et al. 2011); the "Union 2" supernova data set (hereafter SN), which includes 557 supernovae from Hamuy et al. (1996), Riess et al. (1999Riess et al. ( , 2007 The prior on H 0 is 74.2 ± 3.6 km s −1 Mpc −1 .

Flat CDM Model with a Constant Equation of State
We investigate the flat Cold Dark Matter (CDM) model with a constant equation-of-state parameter (w) to characterize Dark Energy with the combination of our angular power spectra from SDSS-III Data Release 8 (DR8) and other data sets. When we combine our DR8 observed angular power spectra with the WMAP7 +HST + SN data set, we find w = −1.07 ± 0.0775, Ω m = 0.2699 ± 0.0166, and σ 8 = 0.85 ± 0.044 (see Figure 28). We also combine our systematic-corrected DR8 angular power spectra with the WMAP7 +HST +SN data set, and we find w = −1.064 ± 0.0757, Ω m = 0.267 ± 0.0163. The systematically corrected angular power spectra give consistent results compared with the observed angular power spectra.
We compare our results with other large-scale structure data sets, such as the latest large-scale structure constraints from galaxy clustering in Blake et al. (2010), which has detected a BAO at z ∼ 0.6 using the spectroscopic survey WiggleZ which includes 200,000 galaxy spectra over 800 deg 2 . They found a similar constraint on the equation of state of dark energy: w = −1.03 ± 0.08 when combined with WMAP7 + SN. This implies that our data set, even though it is purely imaging data, gives a similar constraining power when compared to the latest spectroscopic surveys such as WiggleZ.
We also compare our results with the BAO constraints from SDSS-DR7. When combined with WMAP7 + HST Percival et al. 2010), they found w = −1.10 ± 0.14, while Montesano et al. (2012) used the full shape of P (k) from SDSS-DR7. When they combined it with CMB + HST, they found w = −1.07±0.11. When we combine with the same data set (CMB+HST), we find w = −1.165 ± 0.12, which implies that our data set gives a similar constraining power as with the full 3D DR7 spectroscopic sample (at z < 0.45), while our purely imaging data set is at a higher redshift range (0.45-0.65).
When we compare our measurement to other largescale structure measurements, such as DR7-BAO constraints, we find very similar constraints; for example, on Ω K WMAP7+HST+DR7 gives Ω K = −0.0023 ± 0.0055.

Flat ΛCDM Model
For a flat ΛCDM model, when combined with WMAP + HST, we find Ω Λ = 0.73 ± 0.019, σ 8 to be 0.817 ± 0.023, and H 0 to be 70.5 ± 1.6 s −1 Mpc −1 km, which are consistent with WMAP+HST only, while improving the accuracy over just WMAP+HST by ∼5% for all parameters. We show the improvement on cosmological constraining power over the combination of WMAP7 and HST in Figure 30. When compared with "WMAP+HST+SN," we improve the accuracy by a relatively small amount, and the result is shown in Table 3. We also show the best-fit values of the galaxy bias and the nuisance parameter a for this model in Table 4.

Companion Results
In this paper, we utilize the full power spectrum over 30 < l < 150, including both the broadband shape and the BAO. As a comparison, our companion Paper II (Seo et al. 2012) derives the angular diameter distance scale using the BAO feature over 30 < l < 300 as a standard ruler, while excluding nearly all the non-BAO information. To summarize their method and result, they use the angular power spectra and the covariance matrix shown in this paper, build a reasonable template power spectra based on the estimate, use true galaxy distribution as a function of redshift and the concordance cosmology, and fit for the distance scale, while marginalizing over many free parameters that account for the shape of the broad band. They derive D A (z)/r s = 9.212 +0.416 −0.404 at z = 0.54 and the result is shown to be robust against assumptions they make during the  fitting process. Figure 31 summarizes the BAO fits they derived before and after the systematics correction. In a second companion paper, de Putter et al. (2012), the angular spectra discussed in the present paper are used to derive a strong new upper bound on the sum of the neutrino masses. As neutrinos suppress growth of structure on scales above the neutrino free streaming length, they leave a characteristic signature in the power spectra. To exploit this signal, de Putter et al. (2012) model the galaxy spectra and their neutrino mass dependence, test the model using mocks, and show that it can be safely applied to the multipole range = 30-200, while also considering the conservative range = 30-150. The angular clustering galaxy data are then combined with priors from WMAP7, the HST Hubble parameter measurement, supernova distances, and the (low redshift) SDSS BAO measurement, and the resulting upper bounds are discussed. We quote here the conservative bound, from DR8+CMB+HST+BAO+SN, of Figure 30. 1D marginalized constraints of Ω Λ and H 0 when compared to using only WMAP7 + HST. Figure 31. BAO-only fits derived in Paper II: a stacked C l /C l,sm data of the four redshift bins before (left) and after systematics correction(right). α wg means the best fit D A (z)/r s /8.584 Mpc each case. The solid red line is the best fit for LRG2 as a comparison after the wavenumber is rescaled to z = 0.54. For more details, see Paper II.
Σm ν < 0.35eV at 95% CL, and refer the reader to the article for more details.

CONCLUSION AND DISCUSSION
We have measured the 2D clustering power spectrum of LRGs using the SDSS-DR8 photometric survey. The principal results of this analysis are summarized and discussed below.
Using photometric redshifts, we constructed a large uniform sample of galaxies between redshifts z = 0.45-0.65. This probes a cosmological volume of ∼3 h −3 Gpc 3 , making this the largest cosmological volume ever used for a galaxy clustering measurement. The large volume allows us to obtain a measurement of power going from the smallest scales to the largest. In particular, we probe all the way from the smaller scales such as the BAOs scale out to the scale of matter-radiation equality with one of the most accurate measurements of angular clustering at achieved to date z = 0.45-0.65.
We applied a novel approach of treating systematics by incorporating both the cross/auto-power of systematics with themselves and cross-power with galaxies. This allows us to not only understand the impact of various systematics on observed galaxy number densities, but also allows the application of small corrections at scales where the corrections are small, and thus the uncertainty related to the corrections is negligible. Since we choose only scales that are minimally affected by systematics, we expect that the final cosmological constraints from both pre-/and post-systematic corrections are consistent with each other, which is indeed the case. This method can be improved drastically by two extra components, which will be left to future work. First, we should be able to project out the appropriate modes that are contributed by the systematics in a similar way as is done in CMB map making (Stompor et al. 2002), which can be done when one is estimating the optimally estimated power spectra. Second, we should be able to model the distribution of the systematics (for example, whether it is Gaussian or not) by investigating multi-epoch data that are available in SDSS-DR8. In particular, this is available in the Stripe 82 area, which even though small (in comparison to DR8 full footprint), it contains ∼250 deg 2 which were multiply scanned (∼15-20 repeats for each field). This will allow us to estimate the uncertainty of our systematic corrections properly. These two components will not only significantly improve our understanding of the systematics, they will allow us to push our analysis to include larger angular scales (which are affected by the systematics more significantly). However, both of them require significant undertakings in data collection, code development, and simulations. The improvements will not dramatically change our cosmological interpretation in this paper, therefore, we leave these two components for future projects.
For a flat ΛCDM model, combining our data with WMAP7 + HST, we find Ω λ = 0.7301 ± 0.019 and H 0 to be 70.5 ± 1.6s −1 Mpc −1 km. For an open ΛCDM model, when combined with WMAP7 + HST, we find Ω K = 0.003476 ± 0.00538, improved over WMAP7+HST alone by 40%. For a wCDM model, when combined with WMAP7+HST+SN, we find w = −1.071 ± 0.0775, and H 0 to be 71.31 ± 1.65 s −1 Mpc −1 km, which is competitive with the latest large-scale structure constraints from spectroscopic surveys, such as those by WiggleZ (Blake et al. 2010) and SDSS DR7 spectroscopic surveys, especially in the analysis led by Reid et al. (2010), Percival et al. (2010), and Montesano et al. (2012). This result implies that our data set, even though it is purely imaging data, possesses a similar constraining power as the spectroscopic surveys such as WiggleZ or SDSS-DR7. What we lack in redshift precision, we compensate for by shear volume. This suggests that future and upcoming imaging surveys such as PanStarrs, 40 DES, 41 and LSST 42 can achieve significant cosmological constraints via large-scale structure clustering even when compared to other spectroscopic surveys. This is Paper I of the project, which mostly describes the construction of the data set, treatment of systematics, estimation of the angular power spectra, and, finally, using the overall shape of the angular power spectra over a large range of angular scales to derive constraints on our cosmological models. We refer readers to Paper II (Seo et al. 2012) of the project, which uses only the BAO feature to fit for various cosmological parameters. We also refer to Paper III (de Putter et al. 2012) of the project, which uses the overall shape of the power spectra to fit for various neutrino models.
We thank Nico Hamaus for testing our treatment for nonlinearities using N-body simulations, Pat McDonald for fruitful discussion on our systematic treatments, and Uros Seljak for useful discussions. S.H. acknowledges Martin White for all the useful discussions and encouragement, even though he is already a coauthor. S.H. also thanks UC Berkeley, Department of Energy and Lawrence Berkeley National Laboratory for their support through the Seaborg and Chamberlain Fellowship. This work was supported by the U.S. Department of Energy under Contract No. DE-AC03-76SF00098 and in part by the facilities and staff of the Yale University Faculty of Arts and Sciences High Performance Computing Center. This research used resources of the National Energy Research Scientific Computing Center, which is supported by the Office of Science of the U.S. Department of Energy under Contract No. DE-AC02-05CH11231. Funding for SDSS-III has been provided by the Alfred P. Sloan Foundation, the Participating Institutions, the National Science Foundation, and the U.S. Department of Energy Office of Science. The SDSS-III Web site is http://www.sdss3.org/.
SDSS-III is managed by the Astrophysical Research Consortium for the Participating Institutions of the SDSS-III Collaboration including the University of Arizona, the Brazilian Participation Group, Brookhaven National Laboratory, University of Cambridge, Carnegie Mellon University, University of Florida, the French Participation Group, the German Participation Group, the Instituto de Astrofisica de Canarias, the Michigan State/Notre Dame/JINA Participation Group, Johns Hopkins University, Lawrence Berkeley National Laboratory, Max Planck Institute for Astrophysics, New Mexico State University, New York University, Ohio State University, Pennsylvania State University, University of Portsmouth, Princeton University, the Spanish Participation Group, University of Tokyo, University of Utah, Vanderbilt University, University of Virginia, University of Washington, and Yale University.

QUADRATIC ESTIMATOR
Consider a Gaussian random field x i with x i = 0 and covariance We wish to form an estimator, p α , of p α which is quadratic in the data where Q is symmetric. Requiring the estimator to be unbiased This problem is easiest if we consider a single parameter at a time, with all other parameters held fixed (and absorbed into C (0) ). Thus, we wish to minimize tr[CQ (α) CQ (α) − 2λC (α) Q (α) ].
Taking derivatives with respect to the components of Q (α) gives If the dependence of C on p α is not linear, then we can use a Newton-Raphson iteration, where C (α) is the derivative of C evaluated at the current best value of p. Iterating by replacing p → p + δp until the best-fit δp = 0 results in a maximum likelihood solution. In practice, it only takes a few iterations to achieve the maximum likelihood solution. This approach also results in another fact that is underappreciated in the literature. The above choice of Q (which can have a slightly different form; see Table 4 of Padmanabhan et al. 2003 for more details) produces error bars that are anti-correlated across different band powers. In this paper, we include the window function (which is mostly affected by the mask) before we compare the observed power and the theoretical power.