Far-infrared lines hidden in archival deep multi-wavelength surveys: Limits on [CII]-158$μ$m at $z \sim 0.3-2.9$

TL;DR

This study develops LinSimStack, a linear, confusion-limited stacking framework to search for aggregate [CII]-158 μm emission in broadband FIR maps by tomographically stacking COSMOS2020 galaxies. By fitting continuum SEDs and analyzing residuals across four redshift bins around the peak of star formation, the authors place 3σ upper limits on the mean [CII] intensity that are significantly lower than prior Planck-based claims, and find results more consistent with star-formation-rate–driven [CII] emission. They carefully quantify uncertainties via multiple error-estimation methods, assess contaminant lines, and apply completeness corrections using the stellar-mass function, concluding that [CII] contributes only a few percent at these epochs. The work demonstrates a scalable approach to constrain line emission in archival broad-band surveys and discusses implications for upcoming missions TIM, EXCLAIM, and Euclid in resolving tensions and mapping [CII] across cosmic noon. The methodology is extendable to other lines (e.g., CO, [OI], [OIII]) with higher spectral resolution data.

Abstract

Singly-ionized carbon is theorized to be the brightest emission line feature in star-forming galaxies, and hence an excellent tracer of the evolution of cosmic star formation. Archival maps from far-infrared and sub-millimeter surveys potentially contain the redshifted [CII]-158$μ$m, hidden in the much brighter continuum emission. We present a search for aggregate [CII]-158$μ$m line emission across the predicted peak of star formation history by tomographically stacking a high-completeness galaxy catalog on broadband deep maps of the COSMOS field and constraining residual excess emission after subtracting the continuum spectral energy distribution (SED). We obtain constraints on the sky-averaged [CII]-158$μ$m signal from the three Herschel/SPIRE maps: $11.8\pm10.2$, $11.0\pm8.7$, $9.6\pm9.8$, and $9.2\pm6.6$ $k$Jy/sr at redshifts $z\sim 0.65$, $\sim1.3$, $\sim2.1$, and $\sim2.6$ respectively, corresponding to $1-1.4σ$ significance in each bin. Our $3σ$ upper limits are in tension with past $z\sim2.6$ results from cross-correlating SDSS-BOSS quasars with high-frequency Planck maps, and indicate a much less dramatic evolution ($\sim\times7.5$) of mean [CII] intensity across the peak of star formation history than collisional excitation models or frameworks calibrated to the tentative PlanckxBOSS measurement. We discuss this tension, particularly in the context of in-development surveys (TIM, EXCLAIM) that will map this [CII] at high redshift resolution. Having demonstrated stacking in broadband deep surveys as a complementary methodology to next-generation spectrometers for line intensity mapping, our novel methods can be extended to upcoming galaxy surveys such as Euclid, as well as to place upper limits on fainter atomic and molecular lines.

Figures (12)

• Figure 1: Overview of data used in our analysis. We include transmission curves (left vertical, bottom horizontal axes) for the eight broadband maps used, from Spitzer, Herschel, and SCUBA-2. We overplot the rest-frame wavelengths of selected FIR emission lines, including [CII]-158$\mu$m; these lines are redshifted into the broadband maps. The top horizontal axis labels translate the observer wavelengths into the rest-frame redshift $z$ of a [CII]-158$\mu$m emitter. We trace the $z$ distribution of the COSMOS2020 photometric catalog (in bins of $\Delta z = 0.05$), with number counts on the right vertical axis; number counts for the selection of star-forming or quiescent emitters are also shown. Redshifted C+ emission will appear in the SPIRE maps for specific $z$ ranges. Finally, at the top, we demarcate the $z$-binning used for stacking with LinSimStack ; bins at $z\sim0.65$, 1.3, 2.1, and 2.6 were chosen to overlap with the SPIRE bands.
• Figure 2: Zoom-in on the stacking outputs for the Herschel/SPIRE 250 $\mu$m COSMOS map, illustrating our forward modeling LinSimStack methodology. Left to right: 2D vectors for the data $d$, model $m$, and the residuals $r=d-m$. The model $m$ is a linear combination of all beam-convolved hit-maps constructed using the COSMOS2020 photometric catalog. We also plot sources from the Herschel/SPIRE Point Source Catalog (HSPSC) with estimated fluxes in the HSPSC above 30 mJy, which correlate with the residual flux in the 2D vector $r$ (see Section \ref{['sec:completeness']}, Appendix \ref{['app:linsimstack']}, the corresponding figures, and the discussion for further details on model incompleteness).
• Figure 3: Errors in stacked mean intensities estimated via four methods (see Section \ref{['subsec:errors']}). Fractional errors in the stacked fluxes, obtained as outputs of the LinSimStack are color-coded by band. The statistical lower bound is given by the closed-form error formulation, obtained by writing the stacking problem as a linear regression. Bootstrapping the catalog or the pixels captures the additional systematic errors, as it samples different realizations of pixel noise and emitter variance; empirically, the bootstrapping over the catalog (lower left subplot) typically yields the more conservative error estimates. Jack-knifing (JK) provides an estimate of the cosmic variance (CV) at the scale of the map. We conservatively add in, in quadrature, the JK estimate with the maximum of the other three estimates. This effectively lower-bounds our bootstrapping error estimate by the analytical linear formulation limit. Past analyses Viero_20222013Viero have only used catalog bootstrapping as the error estimation method. The main subplot above (top-left) compares our final error estimates with those from catalog bootstrapping alone. Our estimates are more conservative by definition and additionally include a measure of CV.
• Figure 4: Stacked intensity estimates obtained from LinSimStack for individual COSMOS2020 bins, fit with a modified black-body emission continuum SED model. The top panels of each sub-figure show the stacked mean intensities (blue crosses) and corresponding $1\sigma$ envelope of the emcee fits (blue) for a subset of the catalog bins; the bottom panels show the signal-to-noise of the fit residuals. The vertical dashed line represents the wavelength cut, $\lambda_\text{cut}$, between the two piecewise components of the SED continuum model in Eqn. \ref{['eqn:sed']} (i.e., the frequency $\nu_0$). For a specific $z$ interval, [CII]-158$\mu$m emission is redshifted into a SPIRE band and manifests as excess residual emission over the continuum; these are marked in orange. The SPIRE 500$\mu$m map contains redshifted [CII]-158$\mu$m from the highest two $z$ bins. Different rows and columns correspond to different stellar mass $\log(M/M_\odot)$ and redshift $z$ cuts; all bins shown here are from the star-forming selection. The labels indicate the number count within a bin in the catalog and as predicted by the stellar mass function (SMF) from Weaver_2023. We conservatively inflate the variance on the [CII]-158$\mu$m residual emission by the reduced chi-squared, $\chi_r^2$, statistic of each SED fit. The bottom row shows residuals obtained by adding over all bins at the same redshift interval, with potential [CII]-158$\mu$m emission again marked in orange. For brevity, we only show stacking within the three highest stellar mass bins (for each redshift bin with potential [CII]-158$\mu$m emission); these are the predominant contributors to the CIB and [CII]-158$\mu$m emission. (See Appendix \ref{['app:allstacking']} for stacking results in all bins.) The bottom row shows the residuals obtained by adding all stellar mass and star-forming plus quiescent bins at a redshift.
• Figure 5: Measurements of [CII]-158$\mu$m at $z\sim 0.3-2.9$ compared with [CII]-158$\mu$m-LF estimates in the local universe 2017Hemmati and $z\sim5$Yan_2020. Also shown are theoretical predictions from C+ evolution models, including 2012Gong (blue hatches), Padmanabhan_2019 (dotted purple: best model, shaded band: uncertainty), Yang_2022 (dash-dotted green), and bethermin+22 (cyan line: delooze+14 version; cyan dashed: high SFRD version at high $z$. Our $3\sigma$ upper limits disfavor high-temperature collisional excitation frameworks and best-fit empirical models calibrated to the Pullen_2018Yang_2019Planck measurement. The $1\sigma$ results are more consistent with SFR-scaling models, which calibrate C+ luminosity to the SFR of sources. Additionally, we show the COBE/FIRAS measurement of the monopole spectrum of the CIB as a function of wavelength, matched to the rest-frame redshift of [CII]-158$\mu$m emission. The [CII]-158$\mu$m contribution is likely no more than a few percent of the total CIB, with the intensity history evolving by less than an order of magnitude across cosmic moon.