Millimeter-wave Point Sources from the 2500 Square Degree SPT-SZ Survey: Catalog and Population Statistics

The response of each detector is recorded at 100 Hz as time-ordered data (TOD) as the telescope scans across the sky. We apply a set of filters to the TOD to suppress noise above the temporal frequency corresponding to the map pixel size and low-frequency noise due to atmosphere. The filtering we apply in this work is very similar to that applied in M13. The data are low-pass filtered above a temporal frequency, corresponding to ℓ = 37,500 in the scan direction to remove noise on scales smaller than the chosen map pixel size, 0.25'. To mitigate atmospheric noise, we apply a first-order polynomial subtraction and a high-pass filter below ℓ = 246. Since atmospheric noise will be spatially coherent on the size scale of the detector wedges, we also remove a mean across each wedge of the receiver from all well-performing bolometers at each time step.

The filtered TOD for each bolometer are then coadded into 0.25 by 0.25' pixels by inverse-variance weighting, adding contributions from bolometers to each pixel to form a single map per observation. The weights for each bolometer are calculated from the power spectral density (PSD) of each detector's TOD in the range from 1–3 Hz. We pixelize each field using an oblique Lambert equal-area projection. This choice of projection is important for source finding because it preserves the source shape across the full area of the map. However, it also produces complications in the analysis, since the scan direction rotates with pixel location in the map. The ramifications of this are discussed in Section 3.3.

To make final coadded maps from all the observations, we apply several cuts (which have been previously shown to be useful for SPT data—see, e.g., Schaffer et al. 2011), based on the mean weights and mean rms of the uniform-weight region of each single-observation map. We cut on excessively high median weights, which occurs when a bolometer's TOD has anomalously low noise. In the past, this has been shown to correspond to poor bolometer behavior, such as when a detector changes operating point due to shifts in the amount of loading (Schaffer et al. 2011). For maps with reasonably good weather conditions, the weights scale well with the rms of the map; however, for poor weather days, 1/f noise dominates the rms and the 1–3 Hz range is no longer a good estimate of the weight that should be assigned to that bolometer. Therefore, we also perform a cut on observations where the map rms does not scale properly with the median weight in the map. The single-observation maps that survive the cuts are coadded by inverse-variance weighting each pixel to form a single coadded map per field and per band.

To calibrate the maps, we use both a relative and absolute calibration. The relative calibration of the TOD from one observation to the next is done through repeated observations of the galactic H ii region RCW38 and reference from a thermal calibrator source installed in the bolometer optical path (Schaffer et al. 2011). The absolute calibration is determined from comparisons of the SPT power spectrum to that from the Planck Collaboration et al. (2016b) in the ℓ range from 682–1178. This results in fractional errors in temperature of 1.05%, 1.15%, and 2.24% in 95, 150, and 220 GHz, respectively.

The pointing model used for constructing maps is based on regular measurements of galactic H ii regions in addition to data recorded by sensors on the telescope measuring temperature, linear displacement, and tilt. To check the absolute astrometry, we correct the global pointing of each field by cross-matching the positions of the brightest 40 SPT-SZ sources in each of the three bands per field to source locations in the Australia Telescope 20 GHz (AT20G) Survey catalog (Murphy et al. 2010), which has an rms positional accuracy of 1''. We then fit for a global pointing correction using the cross-matched locations. We iterate on this process of cross-matching and fitting for a correction until the calculated offset is smaller than the residual scatter, and we find the rms residual pointing scatter for the on-average 26 brightest sources with cross-matches in each field to be 4.3'' in decl. and 4.6'' in R.A. · cos(decl.).

Each field is analyzed separately in our pipeline for extracting sources. We then cross-match the single-field catalogs at the end, accounting for places where fields overlap, to form a single catalog for the survey. Field masking is needed to exclude low S/N edges due to turnaround regions of the scan strategy and nonuniform array coverage between bands on the focal plane. However, we also want to define masks such that we have continuous coverage of the full 2500 square degree survey area. This requires that we define separate masks per band for each field, because different bands occupy physically offset locations in different wedges on the focal plane and therefore observe slightly offset regions on the sky. In principle, this choice only adds slightly more complicated bookkeeping for cataloging sources, since now a source detected on the edge of one field in one band could be detected in an adjacent field in a different band. To achieve continuous coverage, we also need to extend the field masks to lower S/N regions compared with M13, making the noise level within each field slightly less uniform.

As the sources we detect are expected to be unresolved by the telescope (except for nearby sources), a source in our maps should manifest in the maps as an SPT beam with the timestream filtering applied. We can improve our S/N for detecting objects with an expected source shape using an appropriate optimal filter. The filter takes advantage of knowledge of the source shape and the noise in the region of the map where the source is located, which includes residual atmosphere, instrument noise, and the primary anisotropies of the CMB, which acts as a source of noise for the detection of compact, extragalactic sources. The first component needed in constructing the source profile is the beam, which is measured using a combination of observations of Jupiter and Venus, as well as the brightest point sources in the fields. The main lobes of the beams are measured to be well described by Gaussian functions with FWHM of 1.7, 1.2, and 1.0' for 95, 150, and 220 GHz, respectively. The sidelobes of the beams are downweighted in the filter, and therefore are unimportant for the point-source analysis pipeline. To model the source profile, we insert a beam into a noiseless map, and then reobserve the source once for each single-observation map using the characteristics of the telescope's performance for that particular observation and the timestream filtering. The Fourier-domain version of this source profile is used as the transfer function for our maps. All the single-observation transfer functions are then coadded into a single transfer function for the coadded map.

Following the formalism in Tegmark & de Oliveira-Costa (1998) and Haehnelt & Tegmark (1996), to maximize the S/N of sources in the map, we filter the map using an appropriately normalized version of the S/N of the source. We apply the optimal filter, ψ, in the Fourier domain given by

$$\begin{eqnarray}&&\psi =\displaystyle \frac{{\tau }^{T}{N}^{-1}}{{\tau }^{T}{N}^{-1}\tau },\end{eqnarray} \tag{ 1 }$$

where τ is the transfer function and N is the 2D noise PSD, resulting in a filtered map still in units of temperature. In addition to the source profile, we also need to characterize the noise around each source. To do this, we find the PSD of the noise of the coadded map by averaging 100 versions of difference maps. Each difference map is constructed by multiplying a randomly chosen half of the individual observation maps by −1 and adding them. The 2D power spectrum of the Fourier transform of each difference map are then averaged to generate a single 2D noise PSD for the coadded map. Because differencing two individual observation maps cancels out the contribution to the noise from the CMB anisotropies, we add back in a Gaussian realization of the best-fit CMB power spectrum from Keisler et al. (2011). Smaller contributions to the noise, such as from secondary CMB anisotropies, and thermal and kinetic SZ effects, are neglected in the filter construction.

The construction of the optimal filter is complicated by two characteristics of the SPT data. Because of the telescope's location at the South Pole, the scan direction is always along constant decl. This means that the effect of timestream filtering is anisotropic in the maps, and we essentially have an anisotropic beam. We account for the smearing of the beam in the scan direction by calculating the transfer function and applying the transfer function during source extraction. But, for point-source work, we use an area-preserving projection, which causes the scan direction to rotate with respect to the axes of the pixel orientation of the maps. Therefore, our anisotropic beam in the map rotates with respect to the pixel x–y location. The second characteristic of the SPT data is that the noise in the maps varies with decl. Because the telescope scans the same distance in azimuth in the same amount of time regardless of elevation, but this distance corresponds to less physical distance farther from the equator, the result is a noise level with a gradient in decl. through our maps, with slightly less noise at higher decl. To account for these two position-dependent complications, we divide up each field in a number of sectors which are small enough that an assumption of zero source rotation and noise uniformity is reasonable. We then calculate separate transfer functions and noise PSDs for each sector independently, construct a single optimal filter for each sector, and extract sources separately per sector. Further description of the process for creating these data products and their salient features can be found in the SPT 2008 data release paper (Schaffer et al. 2011). The number of sectors per field is shown in Table 1.

Essentially, splitting up each field into sectors is a compromise between computation time and accuracy. We test that the sizes of our sectors are appropriate, i.e., that the measured flux density of sources is unaffected by the size of the sector we choose, by applying the transfer function and noise PSD from adjacent sectors to a sector where the effects of noise and scan rotation angle are the most severe, and check that the resultant change in the flux densities of the sources in that sector are below the noise level of the sector. We also test that the noise in different sectors does not differ by more than 5%.

We found in creating and testing our optimal filter that there was residual noise due to incomplete averaging in the creation of the PSD. This resulted in excess noise in the source extraction template, resulting in excess noise in the optimally filtered maps. To mitigate this effect, we apply a smoothing kernel to the optimal filter in the Fourier domain. To test that the strength of the filtering is optimal, we sweep through a range of kernel size while monitoring the noise. As we apply a stronger and stronger smoothing to the filter, we see that the noise level in the optimally filtered map is reduced, indicating that the excess noise being introduced by the filter is being diminished. But, applying stronger smoothing past a certain point eventually causes the noise level to once again rise, as real noise information in the PSD will begin to be cut, and the filter becomes a less realistic description of the actual S/N in the map and therefore less optimal. We take the minimum noise level as our optimized smoothing kernel size.

After optimally filtering the map, we locate and extract source flux densities using a CLEAN algorithm (Högbom 1974). CLEANing was developed originally for radio interferometry, where uneven baseline sampling and a finite number of antennae produce incomplete sampling of the Fourier domain. In turn, this effect produces sidelobes on the beam (a so-called "dirty beam"), which is analogous to the wings on the SPT beam due to the total applied optimal filter. The CLEAN algorithm detects and removes sources iteratively using a template source profile, which allows for the detection of fainter sources hidden underneath the dirty-beam wings of brighter sources. The source template we employ for CLEANing takes into account that we have optimally filtered the map, however, technically this optimal filter is only optimal for a source located at the center of a sector (which is where the simulated beam was placed when calculating the transfer functions for each sector). For sources off center in the sector, this optimal filter is at a slightly incorrect rotation angle. In order to form a template for each source, τ', we rotate the source profile (which is the map space version of the transfer function, τ) to the correct rotation angle for the x–y pixel location in the sector, and then convolve it with the optimal filter for that sector (which is not rotated). Effectively, our source template (in the Fourier domain) is given by

$$\begin{eqnarray}&&\tau ^{\prime} =\psi \tau ,\end{eqnarray} \tag{ 2 }$$

where ψ is the optimal filter function. Each sector of a field has been filtered separately, so we also perform the cleaning separately per sector and then unite the catalogs of detected sources from all fields. Since sources have long wings in the scan direction due to timestream filtering, we need to account for the possibility of false detections from the wings of sources bleeding into a sector from sources just outside the sector. We do this by defining a sector pixel mask to outline the source-finding area for each sector and a second mask that covers a larger area than this sector mask. We define the larger masks such that the extra space on the left and right sides relative to the sector mask edges will be wider than the wings on all but the very obviously brightest sources, which we check by hand if they occur at the edge of a sector. The CLEAN algorithm is applied to the area of the larger mask for each sector, but only the sources that are within the smaller sector pixel mask are saved into the catalog.

To better account for nonuniformity in the noise level across each sector, we construct a scaled-noise map using the weight map for each field's coadded map. We apply the optimal filter for each sector to the inverse of the weight map, and then scale each sector's rms noise by the square root of the ratio of each sector's median weight to its filtered weight map. In essence, we construct a local scaled-noise map, which can be used to construct a local S/N map when combined with the optimally filtered map. Thus, rather than assuming a single noise value per sector when CLEANing, we take into account any local noise nonuniformity and CLEAN down to a locally determined S/N threshold. The most noticeable differences resulting from the implementation of this method arise along the edges of the map, which are noisier than the rms noise of the sectors that include these regions, and fields that were observed with a lead-trail observing strategy and have low-noise strips where the lead and trail observations overlap.

The steps of the CLEANing are as follows:

• 1.  
  Find the location of the brightest pixel in a given sector in the optimally filtered map.
• 2.  
  Rotate the source profile for that sector by the appropriate rotation angle for that x–y pixel location, and convolve it with the optimal filter for that sector. This is the source template.
• 3.  
  Subtract the source template, scaled to the flux of the pixel and multiplied by a loop gain coefficient. The loop gain is a multiplicative factor between 0 and 1 to account for nonideal characteristics of the CLEANing pipeline, such as imperfections in the source model, extended sources, and finite pixelization in the map. We choose a loop gain of 0.1.
• 4.  
  Find the next brightest pixel in the map and repeat the process until all pixels in the map have significance below the chosen S/N detection threshold, in this case 4.5 times the scaled rms noise of that pixel location.

We extract negative sources as well as positive sources during the CLEANing process. Because the CLEANing is performed with a loop gain, bright sources will be broken up into multiple brightest pixels during the CLEANing. Once the CLEANing is finished (i.e., no pixels in the map remain above the chosen significance threshold), the pixels found by the CLEAN are associated into sources using a radius of association that is brightness dependent, scaling from roughly 38'' for detections of 4.5σ up to 2' for detections of 200σ or larger. All of the pixels associated with a single source are used in a centroiding process to find the source's position. The post-CLEAN map, with all sources removed, is called a residual map, and will be used in later steps of the analysis.

Figure 2 gives an example of the CLEAN method, showing cutouts from the raw, coadded map, filtered map, and residual map around a relatively bright source at 150 GHz.

Zoom In Zoom Out Reset image size

After locating sources in the map, we convert from units of CMB temperature fluctuations to units of flux density. Optimally filtering the map is equivalent to fitting the map with a source shape, and the value of the brightest pixel of each source can be used to calculate the integrated source flux. Specifically, we calculate the flux of each source by stacking all of the CLEAN components removed for a given source onto the residual map and taking the maximum in a cutout region. The maps are calibrated in units of CMB temperature fluctuations, so we convert to flux density units using

$$\begin{eqnarray}\begin{array}{rcl}S(\mathrm{Jy}) & = & {T}_{\mathrm{peak}}\cdot {\rm{\Delta }}{{\rm{\Omega }}}_{{\rm{f}}}\cdot {10}^{26}\cdot \displaystyle \frac{2{k}_{{\rm{B}}}}{{c}^{2}}{\left(\displaystyle \frac{{k}_{{\rm{B}}}{T}_{\mathrm{CMB}}}{h}\right)}^{2}\\ & & \times \ \displaystyle \frac{{x}^{4}{e}^{x}}{{\left({e}^{x}-1\right)}^{2}},\end{array}\end{eqnarray} \tag{ 3 }$$

where x = hν/(k_BT_CMB), and ΔΩ_f is the effective solid angle under a filtered source template, given by

$$\begin{eqnarray}&&{\rm{\Delta }}{{\rm{\Omega }}}_{{\rm{f}}}={\left[\int {d}^{2}k\psi ({k}_{x},{k}_{y})\tau ({k}_{x},{k}_{y})\right]}^{-1}.\end{eqnarray} \tag{ 4 }$$

We inspect the detected sources for obvious spurious detections, such as false sources created by the effect of the timestream filtering on bright galaxy clusters and spurious detections very close to extremely bright sources. We also inspect for extended sources, discussed in more detail in Section 4.7. Obvious false detections are trimmed; for the sake of completeness in the catalog, information on extendedness is not used to remove any sources, but is retained as a flag in the catalog.

The raw fluxes in our catalogs are subject to several biases, which must be carefully considered before the fluxes can be used for population statistics. The first is due to the fact that the underlying source number count populations are steep functions of flux. We expect the noise in the map to be Gaussian, but since there are many more dim sources than bright sources, it is much more likely that a detection at a given significance is a dim source on top of a positive noise fluctuation than a bright source on top of a negative noise fluctuation. Therefore, more sources below a significance threshold will be bumped above the threshold and detected as sources due to noise than will be bumped below, resulting in a positive flux bias which most strongly affects low S/N sources. This bias is closely related to Eddington bias, although that term is generally applied to counts as a function of flux rather than the fluxes of individual objects (M13). When applied to individual sources, we refer to this bias as "flux boosting" and its correction as "flux deboosting."

A second bias is due to the fact that we estimate source flux based on peak pixel brightness. A positive noise fluctuation near a source will pull the detected peak position away from the true position and also return a higher flux, whereas a nearby negative noise fluctuation will not have nearly as as strong a corresponding opposite effect on either the returned position or flux (Austermann et al. 2010). For a significance threshold of S/N_meas = 4.5, this is a roughly 5% effect and will be less important for all higher-significance detections (see e.g., Vanderlinde et al. 2010). We therefore neglect this bias in this work.

Finally, a third bias arises from the fact that for sources that we detect only in one or two bands but not all three, the flux(es) for the source in the non-detected band(s) will be subject to a slight negative bias. This is due to the fact that we measure source flux in the non-detection band(s) using a source position determined from a band where the source is detected, and positional uncertainty biases the flux low. This bias is expected to be small given the small positional uncertainty for a 4.5σ detection. We calculate that a 1σ positional offset would result in a flux underestimate of 5%, and therefore neglect this bias.

3.5.1. Bayesian Flux Deboosting

One standard method for dealing with flux deboosting in millimeter and submillimeter surveys is the application of a Bayesian approach, where a posterior probability distribution is calculated given prior knowledge about the underlying source populations (Coppin et al. 2005). The usual Bayesian posterior distribution can be expressed as

$$\begin{eqnarray}&&P({S}_{\mathrm{true}}| {S}_{\mathrm{meas}})\propto P({S}_{\mathrm{meas}}| {S}_{\mathrm{true}})\ P({S}_{\mathrm{true}}),\end{eqnarray} \tag{ 5 }$$

where $P({S}_{\mathrm{true}}| {S}_{\mathrm{meas}})$ is the posterior probability, expressing the probability of the true source flux S_true given the measured value S_meas. $P({S}_{\mathrm{meas}}| {S}_{\mathrm{true}})$ is the likelihood, expressing the probability of measuring a flux S_meas given that the true flux of the source is S_true. Most simply, the likelihood is taken to be a Gaussian with width given by the map noise. P(S_true) is the prior, which expresses previous knowledge about the population of sources being detected, which in our case is proportional to the differential number counts as a function of flux, dN/dS.

Crawford et al. (2010) present an argument for slightly altering the expressions above to account for the fact that we expect the number of sources to rise steeply with decreasing flux and the reality that the telescope observes the sky with some finite resolution (which we further pixelate when creating a map). Therefore, there is a confusion limit due to coexistent faint sources in a single pixel, which contributes to the noise of each detection. The standard Bayesian approach can be modified slightly to account for this:

$$\begin{eqnarray}&&P({S}_{\max }| {S}_{\mathrm{meas}})\propto P({S}_{\mathrm{meas}}| {S}_{\max })\ P({S}_{\max }),\end{eqnarray} \tag{ 6 }$$

where now the posterior, $P({S}_{\max }| {S}_{\mathrm{meas}})$, gives the probability that the highest-flux source contributing to the pixel brightness is S_max given the measured flux of S_meas in that pixel. Similarly, $P({S}_{\mathrm{meas}}| {S}_{\max })$ expresses the likelihood that S_meas will be measured given that the brightest source contributing to that pixel brightness has flux S_max. The likelihood includes the uncertainty in the flux due to the presence of fainter sources. The prior, P(S_max), is still expressed by the differential number counts, dN/dS, but now multiplied by an exponential suppression at low flux representing the probability that no other sources brighter than S_max exist in that pixel.

3.5.2. Simultaneous Three-band Deboosting

Also presented in Crawford et al. (2010) is the framework for expanding the single-band deboosting presented above to a deboosting of fluxes for sources detected in multiple bands simultaneously. Crawford et al. (2010) expands the analysis from one to two bands, and M13 presents the extension to three bands. We use the same method for deboosting as in M13 and present an overview of the methodology below, see M13 for further details.

The goal of multiband deboosting is to estimate the posterior probability for the flux of the source in multiple bands using its measured flux in one or more bands and any prior information known. The simplest way to write this would be

$$\begin{eqnarray}&&\begin{array}{l}P\left({S}_{95}^{\max },{S}_{150}^{\max },{S}_{220}^{\max }| {S}_{95}^{\mathrm{meas}},{S}_{150}^{\mathrm{meas}},{S}_{220}^{\mathrm{meas}}\right)\\ \quad \ \propto \ P\left({S}_{95}^{\mathrm{meas}},{S}_{150}^{\mathrm{meas}},{S}_{220}^{\mathrm{meas}}| {S}_{95}^{\max },{S}_{150}^{\max },{S}_{220}^{\max }\right)\\ \qquad \ \cdot \ P\left({S}_{95}^{\max },{S}_{150}^{\max },{S}_{220}^{\max }\right),\end{array}\end{eqnarray} \tag{ 7 }$$

which would express the three-dimensional posterior probability distribution for the true flux for the detected source in the three bands, given the measured fluxes for that source in three bands. For the multiband prior, one could assume that the priors for each band are independent and therefore could be separated as

$$\begin{eqnarray}&&P\left({S}_{95}^{\max },{S}_{150}^{\max },{S}_{220}^{\max }\right)=P\left({S}_{95}^{\max }\right)\ P\left({S}_{150}^{\max }\right)\ P\left({S}_{220}^{\max }\right).\end{eqnarray} \tag{ 8 }$$

However, in general this assumption would only be accurate if the three bands probed completely separate populations of sources with no overlap. In general, while more synchrotron sources are detected in 95 GHz, and 220 GHz is a stronger probe of dusty sources, there is certainly population overlap between the bands.

To accommodate this issue, we can express the prior as the combination of a prior on flux for one band (for example, 150 GHz), and two priors describing the power law behavior connecting two fluxes as a function of frequency:

$$\begin{eqnarray}\begin{array}{rcl}{S}_{95} & = & {S}_{150}{\left(\displaystyle \frac{{\nu }_{95}}{{\nu }_{150}}\right)}^{{\alpha }_{95\mbox{--}150}}\\ {S}_{220} & = & {S}_{150}{\left(\displaystyle \frac{{\nu }_{220}}{{\nu }_{150}}\right)}^{{\alpha }_{150\mbox{--}220}}.\end{array}\end{eqnarray} \tag{ 9 }$$

Note that the effective band centers for the SPT depend slightly on the assumed spectral index of the source. We assume a flat spectral index of zero, which gives effective band centers of 97.6, 152.9, and 218.1 GHz. M13 found that source fluxes are not affected significantly by making this assumption. Through a change of variables, then, we can express the three-flux prior in terms of one flux and two spectral indices (α):

$$\begin{eqnarray}&&\begin{array}{l}P\left({S}_{95}^{\max },{S}_{150}^{\max },{S}_{220}^{\max }\right)\\ \ =\ P\left({S}_{150}^{\max },{\alpha }_{95\mbox{--}150}^{\max },{\alpha }_{150\mbox{--}220}^{\max }\right)\ \displaystyle \frac{d{\alpha }_{95\mbox{--}150}^{\max }}{{{dS}}_{95}^{\max }}\ \displaystyle \frac{d{\alpha }_{150\mbox{--}220}^{\max }}{{{dS}}_{220}^{\max }},\end{array}\end{eqnarray} \tag{ 10 }$$

where the $\tfrac{d\alpha }{{{dS}}^{\max }}$ can be found from Equation (9).

We then make the assumption that the prior written in this way is made up of three independent components:

$$\begin{eqnarray}&&\begin{array}{l}P\left({S}_{150}^{\max },{\alpha }_{95\mbox{--}150}^{\max },{\alpha }_{150\mbox{--}220}^{\max }\right)\\ \ =\ P\left({S}_{150}^{\max }\right)P\left({\alpha }_{95\mbox{--}150}^{\max }\right)P\left({\alpha }_{150\mbox{--}220}^{\max }\right).\end{array}\end{eqnarray} \tag{ 11 }$$

By separating them, we are assuming that the spectral indices are independent of flux and the two spectral indices are not correlated with each other. Strictly speaking, we know that this assumption of independence is also incorrect—fainter sources tend to have more dust-like spectral indices. More fundamentally, simply changing variables does not change the issue of the priors being correlated, since the amount of information contained in the priors has stayed the same. However, since we are interested in measuring α in this analysis, and allow for the possibility of sources with nontypical spectral indices, expressing the priors in this way allows us to place weak flat priors on both spectral indices between the physically motivated range of −3 ≤ α ≤ 5, allowing the intrinsic population characteristics to emerge.

We now have for the three-dimensional posterior on fluxes:

$$\begin{eqnarray}&&\begin{array}{l}P\left({S}_{95}^{\max },{S}_{150}^{\max },{S}_{220}^{\max }| {S}_{95}^{\mathrm{meas}},{S}_{150}^{\mathrm{meas}},{S}_{220}^{\mathrm{meas}}\right)\\ \quad \ \propto \ P\left({S}_{95}^{\mathrm{meas}},{S}_{150}^{\mathrm{meas}},{S}_{220}^{\mathrm{meas}}| {S}_{95}^{\max },{S}_{150}^{\max },{S}_{220}^{\max }\right)\\ \qquad \ \cdot \ P\left({S}_{150}^{\max }\right)P\left({\alpha }_{95\mbox{--}150}^{\max }\right)P\left({\alpha }_{150\mbox{--}220}^{\max }\right)\ \displaystyle \frac{d{\alpha }_{95\mbox{--}150}^{\max }}{{{dS}}_{95}^{\max }}\ \displaystyle \frac{d{\alpha }_{150\mbox{--}220}^{\max }}{{{dS}}_{220}^{\max }}.\end{array}\end{eqnarray} \tag{ 12 }$$

The likelihood, $P\left({S}_{95}^{\mathrm{meas}},{S}_{150}^{\mathrm{meas}},{S}_{220}^{\mathrm{meas}}| {S}_{95}^{\max },{S}_{150}^{\max },{S}_{220}^{\max }\right)$, is given by a multivariate Gaussian

$$\begin{eqnarray}&&\begin{array}{l}P\left({S}_{95}^{\mathrm{meas}},{S}_{150}^{\mathrm{meas}},{S}_{220}^{\mathrm{meas}}| {S}_{150}^{\max },{S}_{95}^{\max },{S}_{220}^{\max }\right)\\ \quad \ =\ \displaystyle \frac{1}{\sqrt{{\left(2\pi \right)}^{3}\det ({\boldsymbol{C}})}}\exp (-\displaystyle \frac{1}{2}{{\boldsymbol{r}}}^{T}{{\boldsymbol{C}}}^{-1}{\boldsymbol{r}}).\end{array}\end{eqnarray} \tag{ 13 }$$

The noise covariance ${\boldsymbol{C}}$ represents the flux uncertainty due to instrument noise, atmosphere, and uncertainties in the beam and absolute calibration.

The residual vector, ${\boldsymbol{r}}$, is given by

$$\begin{eqnarray}&&{\boldsymbol{r}}=\left[{S}_{95}^{\mathrm{meas}}-{S}_{95}^{\max },{S}_{150}^{\mathrm{meas}}-{S}_{150}^{\max },{S}_{220}^{\mathrm{meas}}-{S}_{220}^{\max }\right].\end{eqnarray} \tag{ 14 }$$

For the flux prior, we estimate the number counts dN/dS based on a sum of synchrotron and dusty population models. Synchrotron populations are calculated using the De Zotti et al. (2005) model at 150 GHz and extrapolated to the other two bands. Dusty populations at 150 and 220 GHz are estimated by use of updated Negrello et al. (2007) models (M. Negrello, private communication). The population at 95 GHz is estimated using an extrapolation of the Negrello et al. (2007) prediction at 850 μm using spectral indices of 3.1 for high-redshift sources (calculated from the spectral energy distribution of the ULIRG Arp 220 shifted to z ∼ 3) and 2.0 for low-redshift sources (from IRAS observations). This is the same method as employed in M13.

There is an asymmetry introduced in our current deboosting algorithm, namely, that one band is chosen to have much stricter prior information applied to it through the flux prior, and the other two bands have much less restrictive priors applied through loose α priors. Therefore, for any given source, with flux information in three bands, the amount of deboosting each band's flux receives depends on the choice made in selecting which band the flux prior is applied to. In Crawford et al. (2010) and M13, this band is termed the "detection band" but this is slightly confusing terminology, since a given source could in fact be detected simultaneously in all three bands or some combination of bands. To avoid this confusion, here we employ the term "flux-prior band" to refer to the band that has the flux prior applied as opposed to a prior on α. In practice, the deboosted fluxes reported in the catalog are calculated using the band with the highest-significance detection in raw flux as the flux-prior band. For number counts, we use the band for which we are calculating number counts as the flux-prior band and then restrict to only sources with a detection in that band.

Since we are interested in calculating posterior distributions for spectral indices in addition to fluxes, we calculate in parallel the posteriors for one flux and two α's:

$$\begin{eqnarray}&&\begin{array}{l}P\left({S}_{150}^{\max },{\alpha }_{95\mbox{--}150}^{\max },{\alpha }_{150\mbox{--}220}^{\max }| {S}_{95}^{\mathrm{meas}},{S}_{150}^{\mathrm{meas}},{S}_{220}^{\mathrm{meas}}\right)\\ \ \ \propto \ P\left({S}_{95}^{\mathrm{meas}},{S}_{150}^{\mathrm{meas}},{S}_{220}^{\mathrm{meas}}| {S}_{150}^{\max },{\alpha }_{95\mbox{--}150}^{\max },{\alpha }_{150\mbox{--}220}^{\max }\right)\\ \quad \ \cdot \ P\left({S}_{150}^{\max },{\alpha }_{95\mbox{--}150}^{\max },{\alpha }_{150\mbox{--}220}^{\max }\right).\end{array}\end{eqnarray} \tag{ 15 }$$

The prior is identical to that used for three fluxes, and the likelihood is very similar:

$$\begin{eqnarray}&&\begin{array}{l}P\left({S}_{95}^{\mathrm{meas}},{S}_{150}^{\mathrm{meas}},{S}_{220}^{\mathrm{meas}}| {S}_{150}^{\max },{\alpha }_{95\mbox{--}150}^{\max },{\alpha }_{150\mbox{--}220}^{\max }\right)\\ \quad \ =\ \displaystyle \frac{1}{\sqrt{{\left(2\pi \right)}^{3}\det ({\boldsymbol{C}})}}\exp (-\displaystyle \frac{1}{2}{{\boldsymbol{r}}}^{T}{{\boldsymbol{C}}}^{-1}{\boldsymbol{r}}),\end{array}\end{eqnarray} \tag{ 16 }$$

the same as before, but where the residual vector is now

$$\begin{eqnarray}\begin{array}{rcl}{\boldsymbol{r}} & = & \left[{S}_{95}^{\mathrm{meas}}-{S}_{150}^{\max }{\left(\displaystyle \frac{{\nu }_{95}}{{\nu }_{150}}\right)}^{{\alpha }_{95\mbox{--}150}^{\max }},{S}_{150}^{\mathrm{meas}}-{S}_{150}^{\max }\right.,\\ & & \left.{S}_{220}^{\mathrm{meas}}-{S}_{150}^{\max }{\left(\displaystyle \frac{{\nu }_{220}}{{\nu }_{150}}\right)}^{{\alpha }_{150\mbox{--}220}^{\max }}\right].\end{array}\end{eqnarray} \tag{ 17 }$$

The likelihood values are identical for the corresponding locations in the different parameter spaces.

From our three-dimensional posterior probability distributions, we marginalize over two of the three parameters in the posterior to find the corresponding one-dimensional posteriors for a parameter of interest. We then integrate the PDFs to the 16%, 50%, and 84% levels in the cumulative distribution to calculate the best-fit values and 1σ error bars.

There are several instances in the analysis pipeline where a cross-match method is employed: cross-matching between the 19 SPT fields within a given band, cross-matching between SPT bands, and cross-matching between the SPT catalog and external catalogs. The same general principle is applied in each case, while the details that differ will be discussed topically in the following sections. A cross-match criterion involving only a radial offset is appropriate when the source densities of the two groups of sources under comparison are comparable, or, in the case of cross-matching with external information, the source density of the external catalog is similar or lower than the SPT catalog and the positional uncertainty is small relative to the typical distance between sources. An appropriate cross-matching radius can then be chosen either analytically or empirically using the measured source density. Depending on the application, either all of the sources within the radial distance are considered associated (in the case of cross-matching between SPT fields for the same band), or the closest candidate within the radial criterion, if one exists, is considered associated (in the case of cross-matching between SPT bands and between the SPT catalog and external catalogs). Selecting a radial threshold that is excessively large will result in falsely associating physically unrelated objects, whereas a radial threshold that is too small risks missing true associations that are shifted in position due to map noise or residual pointing error. Further details of the cross-matching between SPT detections to form a single three-band full-survey-area catalog can be found in the following subsection, and details of cross-matches with external catalogs and redshift information are further detailed in Sections 4.6 and 4.8.

To be included in the source catalog, we require a source to exceed the detection threshold of 4.5σ in raw flux S/N in at least one band. The threshold of 4.5σ was chosen to align with V10 and M13, who calculated purity levels of roughly 90% at 150 GHz and a 4.5σ threshold. We note that the purity level in V10 for 220 GHz was also roughly 90%, whereas it is somewhat lower in the current work (see Section 4.5 below), which is due to the 2008-observed fields having deeper 220 GHz data than the rest of the survey. To form a united single SPT catalog including all fields and bands, we first cross-match across all 19 fields for detections in a single band. About 10% of the full-survey area falls in overlap regions covered by multiple fields, and sources that lie in overlap regions will have repeat detections in different fields. We remove repeat detections by concatenating all fields' detections in each band and employing a radial cross-match as discussed in Section 3.6. We then keep the detection that comes from the map with the lowest noise at that location and throw out the others. To determine the cross-match radius, we use the analytic formalism in Appendix B of Ivison et al. (2007), which takes into account the measured beam FWHM for that band and the source S/N to yield a positional error calculated analytically. Assuming there is equal and uncorrelated error in both positional directions, we use a 3σ positional error for a 4.5σ detection to cross-match, corresponding to 57.8, 40.3, and 34.0'' for 95, 150, and 220 GHz, respectively. (Note: in reality the positional error is correlated in the two orthogonal directions in the map when cross-matching source positions, so assuming errors are uncorrelated is not technically correct.) We test that this is an appropriate radius by comparing to the density of source detections within each single field and also check that we do not associate (and therefore remove) sources detected within the same field. We note that the radial criterion used to associate CLEAN components into sources (as discussed in Section 3.4) corresponds to just under 2σ for a detection of 4.5σ at 95 GHz, the band with the largest beam; therefore, sources within this radius of another source within the same field and band very likely would have been considered a component detection of that source. We remove six sources flagged within the cross-match radius of a source in the same field, five of these are sources detected at 95 GHz, one at 150 GHz, and none at 220 GHz. These six appear to either be multiple detections of the same source or component detections of extended sources.

We additionally remove all sources that lie in regions with overlapping coverage from multiple fields where the source is detected in a field with a higher noise level but not detected in an overlapping field with lower noise, as it is expected that these detections would be false. This step removes 51 sources at 95 GHz, 44 sources at 150 GHz, and 47 sources at 220 GHz, which is roughly 3% of sources in 220 GHz, and a smaller percentage for 95 and 150 GHz. We check that the distribution in S/N of sources trimmed is sensible, i.e., that almost all trimmed sources are near the detection threshold of 4.5σ, and therefore likely to be false detections due to map noise. The one notable exception is a 12.3σ detection in 220 GHz that is removed from field ra1hdec-60 due to overlapping coverage by ra3h30dec-60, which has lower noise at that source location but a non-detection of the source in that band. This source appears to be a flaring radio source that became brighter over the course of observing the 2500 square degree area of the SPT-SZ survey, such that it was brighter in ra1hdec-60, observed in 2010, compared with ra3h30dec-60, observed in 2009. We note that due to detections of this source above 4.5σ in 95 and 150 GHz, this source does survive to the final catalog as SPT-S J015917-6055.9, but with recorded fluxes in the three bands that were not measured contemporaneously.

The next step in creating a multiband catalog is to cross-match across the SPT bands. We employ a radial cross-match method and use a 30'' radius of association, which is chosen similarly to above using the analytical positional uncertainty of SPT sources calculated from the formalism in Ivison et al. (2007). For the band with the widest beam (1.7' at 95 GHz), 30'' is roughly a 1.5σ positional error for a 4.5σ detection. Since the source densities in all three SPT bands are quite low relative to the 30'' association radius, the expected rate of random association of two unrelated sources between bands is also very low. Using just the full-survey average source density, the probability of random association is 0.024%, 0.034%, and 0.012% for 95, 150, and 220 GHz, respectively.

At this step, we also remove any sources with incomplete coverage across the three bands. Because the masks of usable area for each field cover physically offset regions of the sky for the three bands, some area of the survey at the edges will be covered only by one or two bands, and for the sake of consistency, we trim any sources detected in these areas. This removes 115 sources from the final catalog, or about 2% of the catalog.

Our three-band integrated catalog for the full 2530 square degrees of survey area contains 2774 sources detected above 4.5σ at 95 GHz, 3909 at 150 GHz, and 1435 at 220 GHz. An overview of the number of sources above 4.5 and 5.0σ are shown in Table 2. The median purity at this detection threshold across all fields in the survey is 94.4%, 94.8%, and 83.4% at 95, 150, and 220 GHz, respectively; see Section 4.5 for further detail. Cross-matching across SPT bands, this yields a multiband catalog with 4845 total sources detected at a minimum of 4.5σ in at least one band. The noise levels for individual matched-filtered maps are shown in Table 3; taking the median noise level across all fields, 4.5σ corresponds to detections above 9.8, 5.8, and 20.4 mJy in 95, 150, and 220 GHz, respectively. Of the 4845 sources in the catalog, 722 sources are detected at ≥4.5σ in all three bands. 1662 are detected only in 95 and 150 GHz, and 167 are detected only in 150 and 220 GHz. 390 are detected only in 95 GHz, 1358 are detected only in 150 GHz, and 546 are detected only in 220 GHz. Of all the detections in the catalog, roughly 8% have fluxes in different bands drawn from multiple different fields, which is consistent with about 10% of the area of the survey falling in overlap regions covered by multiple fields. Similarly, of all the sources detected above 4.5σ in all three bands, about 9% have fluxes drawn from multiple fields. We compare raw fluxes and deboosted fluxes in the combined catalog in Figure 3. Overplotted are the expected values for spectral indices between the bands, and we see that for the most part, sources follow the characteristic lines for dusty and synchrotron sources. Similarly, we plot α_95–150 versus α_150–220 for both raw spectral indices and deboosted values in Figure 4. We note in these plots, that spectral index does seem to correlate with source brightness, as expected, where the brightest sources are synchrotron dominated. We also note that while there are sources where α_95–150 correlates with α_150–220, there are also numerous sources with spectral indices that are not correlated, indicating sources with a spectral break, which will be discussed further in Section 6. To show the effect of the deboosting, Figure 5 plots deboosted flux as a function of raw flux for each of the three SPT-SZ bands.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

To explore the distributions in spectral indices that we find from deboosting and to separate sources into populations based on spectral index, we normalize each source's posterior probability distribution for α, such that the integral of the marginalized posterior over all possible values of α is unity, and then sum all the posteriors from different sources. In Figure 6, we show these distributions for sources with S/N greater than or equal to 5.0 in both of the bands that a particular spectral index spans. We restrict to higher S/N sources for this part of the analysis to provide a cleaner population separation.

Zoom In Zoom Out Reset image size

From Figure 6, we see that the posteriors for ${\alpha }_{95\mbox{--}150}^{\max }$ show only the presence of a synchrotron population peaking at ${\alpha }_{95\mbox{--}150}^{\max }\sim -0.7\pm 0.6$. As shown in Figure 4, synchrotron sources do dominate the high S/N sources in general, and dusty sources, with a positive spectral index, are much more likely to be below the detection threshold at 95 GHz. In contrast, the posteriors for ${\alpha }_{150\mbox{--}220}^{\max }$ show two peaks in the distribution, representing contributions from both synchrotron and dusty populations, peaking at ${\alpha }_{150\mbox{--}220}^{\max }\sim -0.6\pm 0.6$ and ${\alpha }_{150\mbox{--}220}^{\max }\sim 3.4\pm 0.8$, respectively. Once again, the synchrotron peak is stronger since we are restricting to relatively high S/N detections, which are synchrotron dominated.

We take the minimum of our summed posterior distribution on ${\alpha }_{150\mbox{--}220}^{\max }$ as the dividing criterion to produce separate catalogs of synchrotron and dusty sources. From Figure 6, this produces a population separation at ${\alpha }_{150\mbox{--}220}^{\max }=1.51$. To classify each source as either dusty or synchrotron, we find the probability for each source that ${\alpha }_{150\mbox{--}220}^{\max }\gt 1.51$ from each source's marginalized posterior. If the probability that a source has ${\alpha }_{150\mbox{--}220}^{\max }\gt 1.51$ is less than 50%, we classify the source as synchrotron, and conversely, if the the probability that a source has ${\alpha }_{150\mbox{--}220}^{\max }\gt 1.51$ is greater than or equal to 50%, the source is classified as dusty.

The columns in our catalog are described in Table 4; sources in the catalog are listed in order of detection significance, using the highest-significance detection across all bands. The catalog for the full survey is available online.⁴²

Table 4.  Catalog Column Descriptions

Column	Label	Description
No.
1	Source I.D.	Source IAU identification
2	R.A.	Right ascension (J2000) in degrees
3	Decl.	Declination (J2000) in degrees
4	${S}_{95}^{\mathrm{meas}}/{N}_{95}$	Raw S/N in 95 GHz
5	${S}_{95}^{\mathrm{meas}}$	Raw flux in 95 GHz, (mJy)
6–8	${S}_{95}^{\max }$	Deboosted flux in 95 GHz from integrating 50% of the posterior PDF, with 16% and 84%
		taken as 1σ error bars, (mJy)
9	${S}_{150}^{\mathrm{meas}}/{N}_{150}$	Raw S/N in 150 GHz
10	${S}_{150}^{\mathrm{meas}}$	Raw flux in 150 GHz, (mJy)
11–13	${S}_{150}^{\max }$	Deboosted flux in 150 GHz from integrating 50% of the posterior PDF, with 16% and 84%
		taken as 1σ error bars, (mJy)
14	${S}_{220}^{\mathrm{meas}}/{N}_{220}$	Raw S/N in 220 GHz
15	${S}_{220}^{\mathrm{meas}}$	Raw flux in 220 GHz, (mJy)
16–18	${S}_{220}^{\max }$	Deboosted flux in 220 GHz from integrating 50% of the posterior PDF, with 16% and 84%
		taken as 1σ error bars, (mJy)
19	${\alpha }_{95\mbox{--}150}^{\mathrm{meas}}$	Spectral index between 95 and 150 GHz calculated from the raw 95 and 150 GHz fluxes
20–22	${\alpha }_{95\mbox{--}150}^{\max }$	Deboosted spectral index between 95 and 150 GHz from integrating 50% of the posterior
		PDF, with 16% and 84% taken as 1σ error bars
23	${\alpha }_{150\mbox{--}220}^{\mathrm{meas}}$	Spectral index between 150 and 220 GHz calculated from the raw 150 and 220 GHz fluxes
24–26	${\alpha }_{150\mbox{--}220}^{\max }$	Deboosted spectral index between 150 and 220 GHz from integrating 50% of the posterior
		PDF, with 16% and 84% taken as 1σ error bars
27	$P\left({\alpha }_{150\mbox{--}220}^{\max }\gt 1.51\right)$	Fraction of the spectral index posterior probability distribution above the threshold value
		of 1.51, where a higher value of P means the source is more likely to be dust dominated
28	Type	Classification of a source as either (synch)rotron or (dust) dominated depending on the
		fraction of the integrated 150–220 GHz spectral index posterior above the threshold of 1.51.
		For $P\left({\alpha }_{150\mbox{--}220}^{\max }\gt 1.51\right)\geqslant 0.5$, the source is classified as dusty, for
		$P\left({\alpha }_{150\mbox{--}220}^{\max }\gt 1.51\right)\lt 0.5$, the source is classified as synchrotron.
29	External	Flag on sources with an associated detection in one of the external catalogs we
	counterparts	cross-match: (1): SUMSS, (2): IRAS, (3): AT20G, (4): WISE, (5): PMN,
		(6): AKARI/FIS, (7): AKARI/IRC, (8): RASS. See Section 4.6.
30	Extendedness	Flag on sources that appear to be extended or are multiple members of the same source at
		physically offset locations due to being extended. See Section 4.7.
31	Redshift	Cross-matched redshift from NED, if available. See Section 4.8.
32	Cut classification	Flag indicating a source is a member of the "ext cut" (1), "z cut" (2), SMG list (3), or
		stars (4). See Section 4.10.

Note. The catalog is available for download in the online version and at https://pole.uchicago.edu/public/data/everett20/.

Only a portion of this table is shown here to demonstrate its form and content. A machine-readable version of the full table is available.

Download table as:  DataTypeset image

The completeness of the catalog for a given band is defined as the ratio of the number of sources we detect using the source-finding algorithm compared with the true number of sources in the map for a given flux. Due to the presence of noise in the maps, sources near the detection threshold may be missed by the source finder if they happen to be coincident with a negative noise fluctuation, which pulls their flux below the detection threshold. Completeness is important not only for the robustness of the catalog, but also for calculating number counts, discussed in the following section. The completeness is calculated in practice by performing the source finding on a known population of sources at fixed flux values. At each chosen flux level, we add a set of 100 simulated sources to random locations in the residual map (the optimally filtered map post-CLEANing, which is a good approximation to noise plus a background of sources below the detection threshold of the CLEANing). The source profile used is the real-space version of the transfer function (i.e., a beam with the timestream filtering applied) for the sector that contains the coordinates randomly chosen for the source, rotated to the proper angle. We then run the source finder and cross-match the returned detections with the known inputs. We repeat this process for a broad range of flux levels. The completeness as a function of flux is then given by f_compl (S) = N_recovered/N_input. Since the noise in our maps is to a good approximation Gaussian and sources are rare enough that the noise dominates the distribution of flux in the map, we would expect the completeness to follow an error function of the form

$$\begin{eqnarray}&&{f}_{\mathrm{compl}}(S)=\displaystyle \frac{1}{\sqrt{2\pi {\sigma }^{2}}}{\int }_{S}^{\infty }{e}^{-{\left(S^{\prime} -{S}_{0}\right)}^{2}/2{\sigma }^{2}}{dS}^{\prime} \end{eqnarray} \tag{ 18 }$$

where S₀ is the detection threshold, in this case 4.5 times the mean rms noise in the map for each band. Since this process is computationally expensive, we evaluate the completeness at a few discrete flux levels and fit the error function to those results and use it as a model of our completeness, and we estimate the errors on our completeness estimate using binomial statistics. We repeat this process for each band separately.

Galaxy clusters appear as compact negative signals at 95 and 150 GHz via the thermal SZ effect (see Bleem et al. 2015 for a recent review and catalog release of clusters detected in 2500 square degrees of the SPT-SZ survey; and Sunyaev & Zel'dovich 1972 for background on the SZ effect). Compact clusters with high significance can overlap and therefore cancel out emissive sources, which we do not account for in the completeness calculation. Using an assumed cosmological model and cluster mass function as well as SPT cluster selection functions, M13 calculated an expectation of one cluster large enough to cancel a 4.5σ emissive source per 10 square degrees of the SPT-SZ survey, which corresponds to roughly 10–20 clusters per field or roughly 250 total in the full SPT-SZ survey. Given the relatively low point-source density in the SPT maps above the detection threshold, the likelihood of purely random overlap and cancellation is less than 1% per field for 150 GHz, the band with the highest source density, and even though it is known that point sources and clusters have some preference for clustering, the effect on the completeness due to cluster overlap is expected to be relatively small.

Flux levels averaged over all sectors per each field and band for 50% and 95% completeness are shown in Table 3. The median 95% completeness across all fields is 12.89, 7.60, and 26.83 mJy at 95, 150, and 220 GHz, respectively.

The purity of the catalog as a function of the source S/N is defined as one minus the fraction of sources at that S/N or higher that are expected to be false detections due to noise in the map. To quantify the purity of the catalog, we estimate the number of detections above a given threshold in a simulated noise-only map and compare those with the number detected above the same significance in the real maps. We generate simulated noise maps from difference maps, which contain instrument noise and residual atmosphere. The method for generating difference maps is discussed in Section 3.3. To the noise realizations, we add contributions from the power spectrum of primary anisotropies in the CMB, which is also a source of noise for our source detections. These noise fluctuations have a power spectrum determined from the best-fit ΛCDM model to combined WMAP7 and SPT data (Keisler et al. 2011). We also include an estimate of the thermal SZ effect, as well as contributions from the CIB in terms of a Poisson and clustered component. The component of the noise that we add to our simulations to account for the SZ effect is a Gaussian random field with power spectrum given by fitting measurements in Shirokoff et al. (2011).

Running the source finder on these simulated maps, we calculate the purity as a function of S/N to be

$$\begin{eqnarray}&&{f}_{\mathrm{pure}}=1-{N}_{\mathrm{false}}/{N}_{\mathrm{total}}.\end{eqnarray} \tag{ 19 }$$

Massive clusters in the real maps will contribute to impurity in the source finding because the timestream filtering causes these objects to have positive wings, which can be detected as false sources. However, these false detections are easy to identify and quite rare in the real maps. We remove them from the catalog by hand, and a total of six sources are removed. Thus, there is no need to include them in the purity simulations.

Table 3 shows purity values averaged over all sectors per field and per band for detections ≥4.5σ. The median purity for sources detected at ≥4.5σ across all fields for the full survey is 94.4%, 94.8%, and 83.4% at 95, 150, and 220 GHz, respectively. For sources detected at ≥5.0σ, the median purity across all fields is 98.9%, 97.6%, and 95.1% at 95, 150, and 220 GHz, respectively.

To further characterize the nature of sources in the SPT catalog, we cross-match with seven external catalogs, ranging in wavelength from radio to X-ray. These include:

• 1.  
  The Sydney University Molonglo Sky Survey (SUMSS; Mauch et al. 2003) at 843 MHz.
• 2.  
  The Parkes-MIT-NRAO (PMN) Southern Survey (Wright et al. 1994) at 4850 MHz.
• 3.  
  The Australia Telescope 20 GHz Survey (AT20G, Murphy et al. 2010).
• 4.  
  The Infrared Astronomical Satellite Faint Source Catalog (IRAS-FSC; Moshir et al. 1992) at 12, 25, 60, and 100 μm.
• 5.  
  The Infrared Astronomical Satellite AKARI, IRC Point Source Catalog (Yamamura et al. 2010) at 9 and 18 μm, and the FIS Bright Source Catalog (Ishihara et al. 2010) at 65, 90, 140, and 160 μm.
• 6.  
  The Wide-field Infrared Survey Explorer (WISE) AllWISE Source Catalog at 3.4, 4.6, 12, and 22 μm.
• 7.  
  The ROSAT All-Sky Survey (RASS) Bright Source Catalog (Voges et al. 1999) and Faint Source Catalog (Voges et al. 2000) at X-ray energies 0.1–2.4 keV.

Each external catalog is cross-matched with positions of SPT point sources in the catalog using a radial association criterion, as overviewed in Section 3.6. An appropriate radius for association is determined for each external catalog by looking at the distributions of source separations, selecting a radius such that the probability of a random, false association is approximately 1% and no greater than 2%. The chosen radius for each catalog can be found in Table 5. For most of the external catalogs, the density of sources is low enough that confusion within the SPT beam size is not an issue. For WISE, which has the highest source density, confusion becomes a problem for cross-matching with the detections in the shorter-wavelength WISE bands. Therefore, we restrict the source density in the WISE sources we cross-match with by applying a cut on the WISE catalog using the W4 22 μm band, and restricting to only cross-matching with WISE sources that have W4 flux greater than 5 mJy. We experimented with a more complex cross-matching scheme, incorporating source flux and number density, but found that a simple radial cross-match achieved comparable results.

Figure 7 shows an overview of the wavelengths and detection thresholds of the external catalogs with which we cross-match the SPT-SZ catalog. The figure also shows reference spectral energy distribution (SED) curves for DSFGs, modeled by an Arp 220 profile, and reference SEDs for two examples of flat-spectrum synchrotron sources. Shifting the reference DSFG SED in redshift demonstrates how negative K-correction enables detection of high-redshift dusty sources in millimeter wavelengths. Surveys in the IR observe dusty sources on the Wien side of their SED and therefore shift to a dimmer portion of the spectrum with increasing redshift, in addition to dimming from increasing source distance. In contrast, millimeter-/submillimeter-wavelength surveys observe dusty sources on the Rayleigh–Jeans side of the spectrum and therefore shift to an intrinsically brighter part of the spectrum with increasing redshift, canceling the effect of dimming from increased distance. Table 5 gives an overview of each survey and the number of cross-matches with the SPT catalog. A comparison of cross-matches per catalog, including cross-match overlap between surveys for the total SPT catalog as well as dusty and synchrotron subpopulations is illustrated in Figure 8. The most ubiquitous cross-match for the SPT catalog is with the SUMSS survey in the radio, where 71% of SPT sources have cross-matches in SUMSS. SUMSS is especially useful for cross-matches with synchrotron-dominated sources in the SPT-SZ survey since the wide-field radio survey has full coverage of the SPT-SZ area and is complete to a depth of 6 mJy/beam at 5σ. The SUMSS beam is relatively large, making it unsuitable for cross-matching with high-resolution optical and IR catalogs, but confusion is not a significant issue when compared with SPT, which has a similarly large beam and low source density. For dusty sources, IRAS in the IR is particularly useful for identifying low-redshift dusty galaxies, but both WISE and AKARI overlap IRAS cross-matches considerably, as shown in Figure 8.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Of the 4845 sources in the catalog, 1109 have no cross-matches with external catalogs. 84% of these are detections in only one band, mostly in 150 GHz only or 220 GHz only; 10 sources with no cross-matches in external catalogs have detections in all three bands.

We expect that all extragalactic sources with redshifts greater than z ∼ 0.05 will be unresolved in the maps, given the instrumental beam size of roughly 1'. There is a chance that very nearby sources or bright AGN with extended radio lobes may be resolved in the maps. We take a two-pronged approach to flagging extended sources: first, we fit a cutout around each detected source to a model constructed from the beam profile convolved with a nonsymmetric 2D Gaussian and compare the Δχ² of the fit to a model containing only the beam. Based on looking at the fields with the most obvious extended sources, we use a threshold of Δχ² ≥ 7 to flag sources as extended in the catalog. Second, to ensure that we are catching all sources that are detected as multiple detections of the same source in the CLEANing, we run a by-eye check of all sources within close proximity to other detections and flag sources that appear to be multiple detections at physically offset locations of the same, extended source. Each source flagged as possibly being a multiple detection of the same, extended object, is cross-checked with external catalogs to determine if the detections are indeed from the same object or from distinct objects that appear in our maps with close proximity. For the sake of completeness, we leave all detections in the catalog, but indicate the likelihood that a source is extended. In calculating the number counts, we calculate multiple versions of the counts, including using the extendedness information from both flagging methods. Fluxes for extended sources will be lower limits on the true flux, since the CLEANing is unable to accurately return flux for sources that do not look like our chosen source profile. Using the two methods discussed above, a total of 131 sources from the catalog are flagged as extended.

Redshifts for SPT catalog sources are obtained from a combination of follow-up observations (e.g., Weiß et al. 2013; Strandet et al. 2016) and the literature. We obtain literature redshifts by querying the NASA/IPAC Extragalactic Database (NED) and using an association radius of 0.6'. 743 sources in the catalog have identified redshifts; available redshift information is listed per source in the catalog and shown in Figure 9.

Zoom In Zoom Out Reset image size

A small but interesting subpopulation in the catalog are 10 stars, identified primarily using their cross-matched IRAS flux at 12 μm. In this section, we overview the method for separating stars from other objects in the SPT-SZ catalog; Section 6.1.4 discusses characteristics of the stars we see in the SPT-SZ catalog. As shown in Figures 4 and 10, the stars are not clearly identifiable using SPT data alone: their flux in SPT bands does not set them apart from other SPT sources and their spectral indices in SPT wavelengths span both synchrotron and dusty populations. However, looking at cross-matched flux in IRAS at 12 μm, these sources have considerably higher flux than other sources in the SPT catalog.

The primary selection effect for detecting stars at millimeter wavelengths is a bias toward either large-surface-area and highly luminous stars or stars with excess emission from dust such that they will have sufficient flux at millimeter wavelengths to be detectable. Nine of the 10 stars identified in the SPT catalog are red giants on the asymptotic giant branch (AGB), most of which are late-type M stars. The remaining star is * β Pic, which has a well-known dusty circumstellar disk (Sheret et al. 2004; Riviere-Marichalar et al. 2014). Red giants are luminous, large, and relatively cool, with surface temperatures of order a few thousand Kelvin. Therefore, their stellar flux follows a blackbody distribution peaking around 1 to a few microns (Bedding et al. 1997; Whitelock et al. 1997) with a spectral index of 2 at longer wavelengths, and they are large and bright enough to be detectable at millimeter wavelengths, far into the Rayleigh–Jeans tail of their blackbody spectrum. Dusty galaxies as well as flat- and steep-spectrum synchrotron sources have spectra that rise as a function of wavelength for wavelengths shorter than ∼100 μm, as shown in Figure 7, whereas stars have spectra that are falling between λ = a few μm to ∼100 μm. Therefore, the ratio of IRAS 60 μm/100 μm flux should be less than one for nonstellar objects and greater than one for stars. This ratio can be used with relative success to identify stars, as shown in the right panel of Figure 10, but we find that high IRAS 12 μm flux on its own is a more effective criterion. In theory, it should be possible to use cross-matches with WISE at wavelengths shorter than the IRAS bands to more clearly identify stars, since stars should be even brighter at shorter IR wavelengths than IRAS, closer to the peak of the stellar SED; however, most of the stars observed in the SPT sample are so bright that the WISE flux measurements are saturated and unreliable.

Zoom In Zoom Out Reset image size

To assist with comparing number counts with models and to further characterize source populations within the catalog, we develop three source cuts using extendedness and external cross-match information. Because source fluxes are measured in maps that have been optimally filtered assuming sources are unresolved by the SPT beam, sources that are flagged as extended or measured in the SPT maps as multiple detections will have fluxes that are systematically underestimated and therefore may bias the number counts. We therefore develop two cuts to flag them for removal when calculating the number counts.

First, in the extended cut, or "ext cut," we flag all objects flagged as extended or detected as multiple detections but confirmed to be a single object, using the methods described in Section 4.7. Sources identified as stars are also removed in the counts for this cut, since they are not included in the models with which we compare the counts. The extended cut removes 131 sources from the catalog as a whole, 36 of these are classified as synchrotron dominated and 95 as dust dominated.

Second, we develop a cut to flag all low-redshift objects, using the redshift cross-match information discussed in Section 4.8. Because the extended source flag used in the ext cut involves in part a by-eye inspection of individual sources, a method was sought to remove extended objects more systematically. All extended sources appear large enough in the SPT maps to be resolved by the SPT beam, and therefore should all be relatively local and removable by cutting all objects with low measured redshift. However, cutting below a redshift threshold will remove additional sources as well. The "z cut" trims all sources flagged as stars and all sources with cross-matched redshifts z ≤ 0.1, resulting in flagging 461 sources from the full catalog, of which 248 have a synchrotron classification and 213 have a dusty classification. Looking at the distributions of source angular sizes for SPT sources with NED identifications, we expect that the cut threshold of z < 0.1 will correspond roughly to cutting objects with angular sizes ≳1', roughly the size of the SPT beam. We verify that all sources flagged by ext cut are included in those sources flagged by the z cut.

To more cleanly select sources in the catalog that are likely to be high-redshift SMGs, which in the relatively high-flux range probed by the catalog are likely to be gravitationally lensed, we develop a list of "SPT SMGs" using more strict criteria than the z cut and ext cut source lists. For this cut, we include only dust-dominated sources, we apply the same redshift criterion as the z cut, and we also exclude any remaining detections with IRAS cross-matches. Although a few IRAS detections have been confirmed to be at relatively high redshift (e.g., APM 0827 Irwin et al. 1998), these sources are few in the literature. Furthermore, IRAS detections are unlikely to be high-redshift objects, since, as shown in Figure 7, dusty objects will be observed on the Wien side of the spectrum in the IRAS bands, which will shift to an intrinsically dimmer portion of the spectrum with increasing redshift, in addition to reduced flux from greater distance. In contrast, the negative K-correction of dusty sources in submillimeter/millimeter bands enables the detection of the same luminosity source out to high redshifts, as shown in Figure 11. Additionally, a cut on IRAS objects has been used successfully in previous SPT analyses as a proxy for trimming low-redshift objects (V10 ; M13 ). We also trim SPT detections that are measured as "dipping" in the three SPT bands, meaning that they are sources with a dusty spectral index between 150 and 220 GHz, but a synchrotron spectral index between 95 and 150 GHz. In a reanalysis of SPT number counts from M13, Mancuso et al. (2015) identified a set of sources in the SPT-SZ survey that are relatively bright at 95 GHz and were classified as dusty galaxies using the SPT pipeline in M13. When considering fluxes for each source across a wider range of frequencies than just the SPT data, Mancuso et al. (2015) note that these sources do appear to have significant emission in radio bands indicating synchrotron emission and do not have spectra that would clearly indicate that these sources are DSFGs. The presence of synchrotron emission causes the spectrum in SPT bands to appear dipping, and although we note that the exact classification of these sources remains somewhat unclear, it provides evidence that dipping sources in SPT data are not clearly DSFGs. Furthermore, a set of SPT dipping sources have been observed in preliminary follow-up observations with LABOCA at 870 μm. While a few sources had measured fluxes at 870 μm that would be consistent with the presence of dust, most sources had measured fluxes at 870 μm that were too low to be consistent with a DSFG spectrum, indicating that they may be synchrotron sources with complicated spectra or may be blends of unrelated objects. Follow-up spectroscopy with the Very Large Telescope (VLT) to obtain redshifts for a set of SPT dipping sources have measured redshifts in the range of z = 0.85–2.32, also indicating that these sources are less likely to be high-redshift SMGs. We note that there are a couple dipping sources from V10 and M13 with follow-up observations that confirmed that they are high-redshift lensed objects. However, these objects either appear to have dipping spectral behavior in the SPT bands due to superposition of the high-redshift object with its foreground lens or are a superposition of unrelated objects along the line of sight. Therefore, although we note that high-redshift dusty galaxies may possess significant synchrotron emission, and therefore may manifest in the SPT-SZ catalog with a dipping spectrum, follow-up information on known SPT dipping sources so far has indicated they are less likely to be high-redshift dusty galaxies, and therefore to be conservative, we exclude them from the SMG list. We find a total of 506 sources in the SMG list, of which 73 have detections above 4.5σ at both 150 and 220 GHz.

Zoom In Zoom Out Reset image size

Of the 4845 sources in the catalog, 3980 (82.1%) are classified as synchrotron sources, and 865 (17.9%) as dusty sources, based on the probability that their ${\alpha }_{150\mbox{--}220}^{\max }$ from deboosting is less or greater, respectively, than 1.51, the minimum of the summed posterior distribution of ${\alpha }_{150\mbox{--}220}^{\max }$, as discussed in Section 4.2. 1109 sources in the catalog, or about 23%, have no cross-matches in external catalogs, and of those, 597 are classified as synchrotron and 512 as dusty. 937 or 84% of the sources in the catalog with no external cross-matches are detected in only one band by SPT, and 172 (16%) are detected in at least two bands.

Looking at Figure 4, we see that while a majority of sources in the catalog fit into the paradigm of two populations, dusty and synchrotron, with similar spectral indices between 95–150 and 150–220 GHz, we also see some sources with a spectral break. To categorize different types of behavior, we look at the distributions of ${\alpha }_{95\mbox{--}150}^{\max }$ for dusty and synchrotron sources, and see that ${\alpha }_{95\mbox{--}150}^{\max }$ = 0.5 forms a relatively natural population separation, although this is a somewhat soft threshold. Using ${\alpha }_{95\mbox{--}150}^{\max }$ = 0.5 and ${\alpha }_{150\mbox{--}220}^{\max }$ = 1.51 as population thresholds, we divide the plots in Figure 4 into four quadrants: rising, falling, dipping, and peaking, though we stress that since the population break lines do not fall along α = 0, the behavior of a source in one of the quadrants may not be as simple as the name suggests. For example, a source in the peaking quadrant may have flux that rises with band between all three frequencies but with a spectral index shallow enough that the source was characterized as synchrotron.

6.1.1. Synchrotron Sources

6.1.2. Dusty Sources

6.1.3. Redshift Distribution of SPT-detected Sources

6.1.4. Millimeter-wavelength Star Characterization

Millimeter-wavelength observations of cool stars can provide interesting insight into the nature of these objects. The baseline expected flux measured at millimeter wavelengths is due just to observing the Rayleigh–Jeans tail of the stellar blackbody radiation, where the stars detectable by SPT are those that are large and bright enough that despite observing the SED in a wavelength range where the flux is many orders of magnitude below the star's peak output, it is still detectable above the SPT noise level. As mentioned in Section 4.9, nine of the 10 stars in the SPT-SZ catalog are AGB stars, most of which are M type, and many of which are Mira variables or closely related to Miras, which are known to be very large and luminous, bright enough for the baseline flux in the Rayleigh–Jeans tail to be detectable. However, excess millimeter-wavelength emission above the stellar SED or a spectral break from the expected ν² of the stellar blackbody may indicate the presence of dust, including spinning dust, or potentially stellar winds, though the effects of winds are likely to be subdominant to dust because the stellar atmosphere will be optically thin at millimeter wavelengths (O'Gorman et al. 2017; Tram et al. 2020).

Gathering flux measurements from the literature, Figure 13 shows SEDs for each of the 10 stars detected in the SPT-SZ catalog. Other than * β Pic, which is an A-type star and therefore much hotter than the other stars detected in the catalog (but which has a well-known dusty debris disk causing excess emission in longer wavelengths), the stars detected in the SPT-SZ catalog are expected to have effective temperatures of roughly a few thousand Kelvin, and therefore should have blackbody spectra that peak at ∼1—a few microns. Therefore, to explore the possibility of excess flux at millimeter wavelengths, we fit a simple blackbody model where the model is constrained using only data in the wavelength range of the blackbody peak, 1.25–5 μm, to fit for the blackbody effective temperature and the star's angular diameter.

Zoom In Zoom Out Reset image size

Six stars in the catalog show flux in the SPT bands relatively consistent with blackbody fits to just the baseline stellar SED: * P Dor, * bet Gru, * pi.01 Gru, V* R Hor, V* X Pav, and V* NU Pav. The other four stars identified in the SPT data have somewhat different spectral behavior from following the tail of the stellar blackbody. Two stars show excess emission but similar spectral indices: * β Pic and V* RZ Sgr. As mentioned above, * β Pic has a well-documented debris disk with a median dust temperature of 79 K (Riviere-Marichalar et al. 2014). This star shows strong excess flux in the SPT bands relative to the expected stellar SED and a spectral index slightly steeper than the expected ν² of the stellar blackbody, likely due to dust modification of the blackbody spectrum. V* RZ Sgr shows an excess of flux in millimeter wavelengths but with a typical blackbody spectral index. It is known to have an optical nebula (Whitelock 1994) and a circumstellar shell large enough to be resolved by IRAS at 60 μm with a measured radius of 4.3' (Young et al. 1993), and therefore the extra emission observed in the SPT bands is consistent with the significant presence of dust. Correspondingly, V* RZ Sgr is observed to be extended in the SPT-SZ catalog, and is flagged accordingly. Therefore, the measured SPT fluxes are likely to underestimate the true flux, and this is the likely cause of the disagreement between the measured SPT fluxes and those from Planck (Planck Collaboration et al. 2018), as shown in Figure 13, given the larger Planck beam.

Two stars show spectral indices distinctly different from the stellar blackbody: V* RR Tel and del02 Gruis. V* RR Tel shows quite flat spectral indices as well as excess flux in SPT bands, although we also note that the simple blackbody model is a poor fit to the data. V* RR Tel is known to be a symbiotic nova, with a red giant in mutual orbit with a white dwarf (Ivison et al. 1995), and its distinct spectrum may indicate the presence of significant stellar winds (Güdel 2002) or the effect of ionization from the white dwarf. del02 Gru, a red giant, also has a relatively flat measured spectral index in SPT bands. We note that both V* RR Tel and del02 Gru have spectral indices measured between the SPT bands that are relatively consistent with measured fluxes in radio catalogs, where the measured flux in radio frequencies has clearly departed from the blackbody spectrum. Similarly, V* R Hor, V* NU Pav, and * bet Gru have detections in radio catalogs that also show a break from the blackbody spectrum. The flux measurements in SPT bands for these three are relatively consistent with a blackbody spectrum, but show some departure, potentially consistent with their radio fluxes.

As can be seen in Figure 10, three sources in the SPT catalog have cross-matches with sources in IRAS with fluxes at 12 μm comparable to the stars but are not identified as stars. Looking at each of these objects individually by hand and comparing with data in external surveys, two appear to be likely false cross-matches due to blends of multiple objects superimposed in the SPT maps along the line of sight, making accurate cross-match identification difficult. A third SPT object appears in the catalog as a repeat cross-match with the IRAS source identified as V* RZ Sgr, due to either being a blend of unrelated objects along the line of sight near V* RZ Sgr or possibly a multiple detection of V* RZ Sgr itself, which is known to be extended, as noted above.

6.2.1. Synchrotron Source Population

Differential number counts per band for synchrotron-dominated sources are shown in Figure 14 along with comparison to two models: De Zotti et al. (2005) and Tucci et al. (2011), neither of which have been fit to the SPT counts. The SPT counts shown are calculated using the method described in Section 5 on the SPT source population with the z cut flagged sources removed. Plots comparing non-cut and various cut versions of the synchrotron counts are shown in Figure A1 in the Appendix.

Zoom In Zoom Out Reset image size

The De Zotti et al. (2005) cosmological evolution model includes separate components for multiple synchrotron populations, including primarily steep-spectrum radio sources and two populations of flat-spectrum sources (blazars): FSRQs and BL Lac objects. It describes each population with a comoving luminosity function extrapolated to higher frequencies using a simple power law (α = −0.1) for flat-spectrum sources and some spectral steepening for steep-spectrum sources.

Similar to De Zotti et al. (2005), the Tucci et al. (2011) model extrapolates source counts using spectral behavior measured at low radio frequencies (5 GHz). But the extrapolation is developed using characteristics of the physical mechanisms of emission for different populations, focusing specifically on flat-spectrum sources, which dominate the number counts at centimeter to millimeter wavelengths and fluxes brighter than ∼10 mJy. The spectrum of emission from flat-spectrum sources is expected to break at some frequency in the range of 10–1000 GHz and steepen at higher frequencies due to both electron cooling as electrons are injected from the AGN core into the jets and from a reduction of the apparent size of the optically thick core with increasing observing frequency, such that high-frequency observations become dominated by the optically thin jets, which have a steeper spectrum (Tucci et al. 2011). This effect is most prominent for higher-flux sources, because these are more likely to be flat-spectrum sources. The SPT-SZ counts are compared with the C2Ex version of the Tucci et al. (2011) model, which is the version that best fits data at frequencies ≳100 GHz, as confirmed mainly with comparison to Planck ERCSC counts, which has strong constraining power at the highest-flux ranges due to full-sky coverage.

While historically the De Zotti et al. (2005) model has been broadly successful in extrapolating to higher frequencies (De Zotti et al. 2010), because the De Zotti et al. (2005) model does not include a spectral break for flat-spectrum sources, we expect that it will become less of a good fit to the counts relative to the Tucci et al. (2011) model with increasing observing frequency and increasing source flux. Looking at Figure 14, while the De Zotti et al. (2005) model is in moderate agreement with the SPT-SZ counts at 95 GHz, it becomes an increasingly poor fit to the counts at higher frequencies, particularly at high fluxes, where FSRQs will dominate. The Tucci et al. (2011) model is a reasonably good fit to the data at all three SPT-SZ frequency bands across the flux ranges probed by the SPT-SZ catalog.

6.2.2. Dusty Source Population

Differential number counts per band for sources in the the SMG list are shown in Figure 15 with comparisons to the lensed components of three representative models: Béthermin et al. (2012), Cai et al. (2013), and Negrello et al. (2007). Because measured fluxes in the SPT-SZ catalog for low-redshift sources that are extended will underestimate the true source flux, we restrict our comparisons of the number counts for dusty populations with models to the SPT SMG list, for which comparison with models is the most straightforward. The SPT SMG list comprises dusty sources where all low-redshift sources, IRAS cross-matches, and dipping sources, which may be contaminated by blending with source lenses or unrelated objects along the line of sight, have been removed. We expect this population to be dominated by sources that are magnified by gravitational lensing. None of the models considered here have been fit to the SPT number counts. Figure A2 in the Appendix shows a comparison of different cut versions for the dusty number counts, including no cuts (only stars), ext cut, z cut, and the SMG list. All three models we compare with combine forward-physical and backward-phenomenological components to describe different populations of the observable galaxy population, including late-type warm (starburst) and cold (normal) galaxies, as well as lensed and unlensed spheroidals and protospheroidals.

Zoom In Zoom Out Reset image size

The Negrello et al. (2007) model includes counts for protospheroidals as modeled using the physical model by Granato et al. (2004). The Granato et al. (2004) counts have been rescaled at 850 μm to agree with counts from the SCUBA SHADES survey (Coppin et al. 2006). Protospheroidals virialized generally at z ≳ 1.5, and the z ≲ 1.5 contribution is considered to be dominantly from starburst and disk galaxies, which are modeled from local luminosity functions at 60 μm. The strongly lensed component to the Negrello et al. (2007) model was calculated similarly to Perrotta et al. (2002, 2003) but using the Granato et al. (2004) model for protospheroidals. When comparing with the SMG list, as in Figure 15, we have trimmed sources with fluxes greater than 200 mJy at 60 μm from the Negrello et al. (2007) model. Of the dusty source counts models considered in this work, the Negrello et al. (2007) lensed model is the closest to reproducing the SPT SMG counts, but we note that it is the model that is the most fine-tuned to reproduce counts at submillimeter/millimeter wavelengths.

The Béthermin et al. (2012) model includes main-sequence and starburst galaxies as the two main components of the model, using one SED per component from libraries from Herschel. Because phenomenological or hybrid models are limited by lacking physical underpinnings describing the evolution of the luminosity function, instead the Béthermin et al. (2012) model is based on two distinct star formation mechanisms and their evolution, one for each galaxy component, based on the work in Sargent et al. (2012). The contribution from strong gravitational lensing is accounted for by applying a magnification factor to the luminosity function. The lensed component of the Béthermin et al. (2012) model overestimates the SPT SMG number counts for all but the very lowest fluxes at both 150 and 220 GHz.

The Cai et al. (2013) model is a hybrid model, combining a physical, forward model for spheroidal galaxies and backward-evolution model for late-type galaxies, based on observations that early-type galaxies are dominated by older stellar populations, while late-type galaxies have younger stellar populations. They improve on previous models by considering components of the flux for protospheroidal galaxies from star formation and central AGN in a unified way, rather than being considered separately. Protospheroidal galaxies are modeled using Granato et al. (2004), and low-z galaxy populations are considered in two populations: "warm" starburst galaxies and "cold" late-type galaxies. A magnification factor is applied to account for strong lensing of high-redshift protospheroidals. While the Cai et al. (2013) model counts are lower than Béthermin et al. (2012), this model also overpredicts the SPT SMG counts for generally all but the lowest flux bins probed by the SPT data at both 150 and 220 GHz.

Disagreement of the data with models can help place constraints on the maximum magnification factor for the strong gravitational lensing of protospheroidals. Both the Cai et al. (2013) and Béthermin et al. (2012) models assume no upper bound on strong lensing magnification factor (assuming sources were pointlike). According to Bonato et al. (2014) (see Figure 2), applying a maximum magnification factor of 20–30 (corresponding to a physical source size of slightly less than ≃3 kpc Lapi et al. 2012), provides good agreement with SPT source counts from M13 for dusty sources with no IRAS cross-matches. The SMG list reported in the current work is more conservative in defining SMG candidates relative to M13, restricting by IRAS cross-match but also measured redshift and removing dipping sources. The number counts for the SMG list in the current work are consistent with the dusty counts with no IRAS cross-matches in M13 at 220 GHz and are slightly lower at 150 GHz. A maximum magnification factor of 20–30 is in good agreement with follow-up observations with ALMA of a set of roughly 50 lensed SMGs from the 2500 square degree SPT-SZ survey, where the measured median magnification factor is ∼6, with a maximum of 30 (Spilker et al. 2016). The implied source physical sizes from the Lapi et al. (2012) model for a maximum magnification factor of 30 are also consistent with follow-up observations, where the measured size distribution of strongly lensed sources is found to be consistent with that of unlensed sources (Spilker et al. 2016).

As the third and final compact-source data release from SPT-SZ, the full 2530 square degree analysis covers a factor of 3.3 times the area of the previous release, M13, including 1759 square degrees of previously unanalyzed data. Due to alterations of the source-finding pipeline and differences in mask areas to aid full-survey coverage, the five sky fields covered by the previous two analyses, V10 and M13, have been reanalyzed in the current analysis. For ra5h30dec-55 and ra23h30dec-55, which were originally observed in 2008 and then re-observed in 2010 and 2011 to add 95 GHz coverage and greater depth at 150 and 220 GHz, we have incorporated the previously unanalyzed 2010 and 2011 data. Although a large majority of the sources extracted in the five reanalyzed fields are consistent with past reported catalogs, there are slight differences in sources extracted between M13, V10, and the current analysis. These differences are due mainly to the lower noise in 150 GHz for ra5h30dec-55 and ra23h30dec-55, slight differences in the masks used for each field, and the slightly different treatment of map noise used for source detection, as discussed in Section 3.4. We confirm that the sources that differ generally have S/N values very close to the detection threshold or are located at the edges of the survey area, which are affected by slight differences in masking.

Figure 16 shows a comparison of total number counts between the three generations of compact-source catalog releases with data from SPT-SZ: V10, M13, and the current work. As expected, the increase in sky area with generally comparable noise level also reduces the error bars on the calculated number counts and adds a few flux bins of counts that were either upper limits or missing from M13. The error bars on the uncut version of the counts, which are most directly comparable to the M13 counts, reduce by roughly 50%, consistent with the amount of increase of area.

Zoom In Zoom Out Reset image size

The smaller error bars allow the SPT number counts to be more constraining of the parameters of galaxy evolution models. For the synchrotron counts, the Tucci et al. (2011) model more clearly agrees with the number counts than the older De Zotti et al. (2005) model, as shown in Figure 14, whereas the M13 counts showed a weaker preference between models, particularly at 95 and 150 GHz. Since the main difference between the two models is the inclusion of a spectral break for FSRQs, the greater constraint shows a clear preference for the presence of a spectral break, although the counts are not constraining enough to provide much further information on the models, such as the break frequency, which might further constrain the AGN core size. Similarly, the smaller error bars for the dusty counts also provide clearer constraints, particularly on parameters governing lensed sources, such as the expected lensing magnification factor, as discussed in the previous section.

The upper panels of Figure 17 provide a comparison between the total counts from SPT-SZ from the current work (with only stars removed) with number counts from the all-sky survey from Planck at 100, 143, and 217 GHz (Planck Collaboration et al. 2013), and number counts from the Atacama Cosmology Telescope (ACT; Gralla et al. 2020). The total counts from ACT presented here have been drawn from the sum of synchrotron- and dust-dominated subpopulations, as presented in Table 5 of Gralla et al. (2020). The lower panels of Figure 17 compare synchrotron-dominated counts from SPT-SZ (similarly shown with only stars removed) compared with number counts from the Planck Multi-frequency Catalog of Non-thermal Sources (Planck Collaboration et al. 2018), excluding Galactic-plane sources using the GAL070 mask, and counts for AGN sources from ACT (Gralla et al. 2020). The total and synchrotron-dominated number counts from SPT-SZ are broadly consistent with both Planck and ACT across the full range of fluxes overlapping between the different surveys.

Zoom In Zoom Out Reset image size

As an overview of the different source populations present in the SPT-SZ catalog, Figure 18 shows cumulative number counts at 220 GHz, including total, synchrotron-, and dust-dominated source populations. Contributions to dusty counts are considered in three components: low-z LIRGs, SPT SMGs, and unlensed high-z sources, with empirical counts from the SCUBA-2 instrument (Geach et al. 2017), scaled from 850 μm to 220 GHz using an SED for Arp 220 shifted to z ∼ 2.5. SPT number counts for low-z dusty sources have been calculated using sources that are trimmed by the z cut. We know the SPT flux measurements for these sources will be biased low, and this is one of the causes of the slight discrepancy with the Béthermin et al. (2011) model counts for these sources. The counts for this population are shown here primarily for illustration. The Negrello et al. (2007) lensed-only model agrees well with cumulative number counts calculated from our SPT SMG list.

Zoom In Zoom Out Reset image size