---
abstract: 'We present new observations from Z-Spec, a broadband 185 - 305 GHz spectrometer, of sub-millimeter bright lensed sources recently detected by the Herschel Astrophysical Terahertz Large Area Survey (H-ATLAS). Four out of five sources observed were detected in CO, and their redshifts measured using a new redshift finding algorithm that uses combinations of the signal-to-noise of all the lines falling in the Z-Spec bandpass to determine redshifts with high confidence, even in cases where the signal-to-noise in individual lines is low. Lower limits for the dust masses ($\sim$a few 10$^{8}$ M$_{\odot}$) and spatial extents ($\sim$1 kpc equivalent radius) are derived from the continuum spectral energy distributions, corresponding to dust temperatures between 54 and 69 K. The dust and gas properties, as determined by the CO line luminosities, are characteristic of dusty starburst galaxies, with star formation rates of 10$^{2-3}$ M$_{\odot}$ yr$^{-1}$. In the LTE approximation, we derive relatively low CO excitation temperatures ($\lesssim 100$ K) and optical depths ($\tau\lesssim1$). Using a maximum likelihood technique, we perform a non-LTE excitation analysis of the detected CO lines in each galaxy to further constrain the bulk molecular gas properties. We find that the mid-$J$ CO lines measured by Z-Spec localize the best solutions to either a high-temperature / low-density region, or a low-temperature / high-density region near the LTE solution, with the optical depth varying accordingly.'
author:
- |
    R. E. Lupu,$^{1\ast}$ K. S. Scott,$^{1}$ J. E. Aguirre,$^{1}$ I. Aretxaga,$^{2}$ R. Auld,$^{3}$ E. Barton,$^{4}$ A. Beelen,$^{5}$ F. Bertoldi,$^{6}$ J. J. Bock,$^{7,8}$ D. Bonfield,$^{9}$ C. M. Bradford,$^{7,8}$ S. Buttiglione,$^{10}$ A. Cava,$^{11,12}$ D. L. Clements,$^{13}$ J. Cooke,$^{4,8}$ A. Cooray,$^{4}$ H. Dannerbauer,$^{14}$ A. Dariush,$^{3}$ G. De Zotti,$^{10,15}$ L. Dunne,$^{16}$ S. Dye,$^{3}$ S. Eales,$^{3}$ D. Frayer,$^{17}$ J. Fritz,$^{18}$ J. Glenn,$^{19}$ D. H. Hughes,$^{2}$ E. Ibar,$^{20}$ R. J. Ivison,$^{20,21}$ M. J. Jarvis,$^{9}$ J. Kamenetzky,$^{19}$ S. Kim,$^{4}$ G. Lagache,$^{22,23}$ L. Leeuw,$^{24,25}$ S. Maddox,$^{16}$ P. R. Maloney,$^{19}$ H. Matsuhara,$^{26}$ E. J. Murphy,$^{27}$ B. J. Naylor,$^{7}$ M. Negrello,$^{28}$ H. Nguyen,$^{7}$ A. Omont,$^{29}$ E. Pascale,$^{3}$ M. Pohlen,$^{3}$ E. Rigby,$^{16}$ G. Rodighiero,$^{30}$ S. Serjeant,$^{28}$ D. Smith,$^{16}$ P. Temi,$^{31}$ M. Thompson,$^{9}$ I. Valtchanov,$^{32}$ A. Verma,$^{33}$ J. D. Vieira,$^{8}$ J. Zmuidzinas$^{7,8}$\
    \
    \
    \
    \
    \
    \
    \
    \
    \
    \
    \
    \
    \
    \
    \
    \
    \
    \
    \
    \
    \
    \
    \
    \
    \
    \
    \
    \
    \
    \
    \
    \
title: 'MEASUREMENTS OF CO REDSHIFTS WITH Z-SPEC FOR LENSED SUBMILLIMETER GALAXIES DISCOVERED IN THE H-ATLAS SURVEY'
---

INTRODUCTION
============

Galaxies detected by their thermal dust emission at submillimeter (submm) and millimeter (mm) wavelengths ($\lambda \approx 250-2000\,\mu$m) comprise an important population of massive systems in the early Universe that are thought to be undergoing an early phase of intense star formation in their evolution [@blain02]. Dust grains within star-forming regions in these galaxies are heated by incident optical and ultraviolet (UV) radiation from young stars and thermally re-radiate this energy at far-infrared (far-IR) to mm wavelengths, with the peak of dust emission occurring at $\sim60-200\,\mu$m in the rest-frame [@soifer91]. It is estimated that about half of all star-formation in the Universe is heavily obscured by dust and therefore difficult to identify in even the deepest surveys at optical/ultraviolet wavelengths [@Puget1996]. Observations at submm/mm wavelengths sample the Rayleigh-Jeans tail of the thermal dust spectrum, which rises steeply with frequency $\sim\nu^{3.5}$ [@dunne00]. For observations at $\lambda > 500\,\mu$m, the climb up this steep spectrum with increasing redshift roughly cancels the effect of cosmological dimming with increasing distance [e.g., @blain02], making galaxies with a fixed luminosity have roughly the same observed flux density at submm/mm wavelengths for redshifts between $1 < z < 10$. This allows a distance-independent study of dust-obscured star-formation and galaxy evolution spanning the epoch of peak star formation activity in the Universe [$z\sim2-3$, e.g., @chapman05]. Although first predicted by @Low68, the population of high-redshift and heavily dust-obscured galaxies (submillimeter galaxies, SMGs) was first revealed a decade ago [@smail97], and several wide-area surveys at 850$\mu$m – 1.2mm have been carried out since then [e.g., @weiss09b; @austermann10; @Coppin2006; @Bertoldi2007; @scott08], mapping a total of $\sim4$deg$^2$ of sky. More recently, much larger area surveys have been undertaken with the South Pole Telescope [SPT, @vieira09] at $\lambda=1.4-2$mm, the Balloon-borne Large Aperture Submillimeter Telescope [BLAST, @pascale08] at $\lambda=250-500\,\mu$m, and the [*Herschel Space Observatory*]{} [@Pilbratt10] at $\lambda=55-670\,\mu$m. Mapping a total area of $\sim200$deg$^2$ to date [@pascale08; @vieira09; @Eales10], these surveys have uncovered a population of rare, and unusually bright, distant galaxies. Their inferred IR luminosities and high redshifts are consistent with a significant fraction of these extremely bright submm/mm galaxies being gravitationally lensed [@negrello07], but proof requires extensive multi-wavelength follow-up campaigns. Their observed flux densities can be magnified by factors $> 10$ due to lensing by intervening foreground galaxies or clusters, as observed in similarly bright systems [e.g., @Swinbank2010; @Solomon2005]. By targeting lensed objects, we can study the typical properties of the star forming galaxies in the early Universe that would otherwise be inaccessible in a blank survey due to sensitivity limitations and source confusion. The ongoing Herschel-Astrophysical Terahertz Large Area Survey [H-ATLAS, @Eales10] in the Science Demonstration Phase (SDP) has already covered 14.4 deg$^{2}$ out of the $\sim$550 deg$^{2}$ planned, resulting in $\sim$6600 sources [@Clements10; @Rigby10] with fluxes measured at 250, 350, and 500 $\mu$m using the Spectral and Photometric Imaging Receiver [SPIRE, @Griffin10; @Pascale10], and fluxes at 100 and 160 $\mu$m obtained with the Photodetector Array Camera and Spectrometer [PACS, @Poglitsch10; @Ibar2010]. Given the large areal coverage, H-ATLAS can detect the brightest (i.e. rarest) distant submm galaxies and is the first example where the efficient selection of lensed galaxies at submm wavelengths has been demonstrated [@Negrello10]. To understand the nature of these galaxies, in particular whether they represent a previously undiscovered population of intrinsically bright sources [e.g., @devriendt10] or are relatively normal starburst galaxies lensed by foreground structures [e.g., @negrello07], requires both complementary data at other wavelengths and measurements of their redshifts. However, measuring spectroscopic redshifts for these sources is challenging: their positional accuracy from submm/mm imaging is often poor due to diffraction limitations at these long wavelengths, and they tend to be highly extincted by dust, making spectroscopic measurements from optical ground-based telescopes difficult [e.g., @chapman05]. Photometric redshifts obtained using submm bands are very useful for estimating the high redshift nature of such sources, but suffer from errors due to the degeneracy between the dust temperature and the redshift, which limit their precision to $\Delta z  \approx  0.3$ [@Aretxaga07; @Hughes02]. When the photometric redshift estimates involve SED template fitting, errors can also arise from our limited knowledge of the intrinsic SMG SED from FIR to radio, and its evolution with redshift. Combined with multi-wavelength data, training sets of spectroscopic redshifts may prove useful for reducing these errors for application to the large photometric datasets from ongoing and future surveys.

SMGs contain large reservoirs of molecular gas [$10^{10-11}$M$_{\odot}$, @tacconi08], whose cooling is dominated by the rotational lines of CO, almost equally spaced by $\sim$115 GHz in the rest frame. Thus, the CO line detections at wavelengths beween 1 cm and 1 mm (30-300 GHz) offer the most direct measurement of their redshifts. However, with the exception of only three other CO redshifts measured to date [@daddi09a; @weiss09a; @Swinbank2010], CO detections have largely been limited to SMGs whose redshifts were already known from optical spectroscopy [e.g., @Frayer98], as a consequence of the small instantaneous bandwidth of typical mm-wavelength receivers. Z-Spec overcomes this limitation due to its large bandwidth, covering the entire 1-1.5 mm atmospheric window, allowing simultaneous observations of multiple CO lines for galaxies at redshifts $z > 0.5$. Although the potential of using the CO ladder for redshift determination is well known [e.g., @silk1997], due to sensitivity limitations of current instruments, only large area submm surveys can provide a significant number of sources bright enough for such measurements.

These spectra can be used not only for an efficient redshift determination, but also to constrain the physical properties of the gas and dust (e.g., mass, density, temperature) in these galaxies [e.g., @bradford09], by measuring the CO line strengths and the continuum slope. The analysis of the CO properties requires measurements of multiple CO lines, often involving the use of multiple instruments. To date, several spectral line energy distributions (SLEDs) for the CO molecule have been constructed for small mixed samples of galaxies and quasars [@Papadopoulos2010b; @Wang2010; @Bayet2009], or individual objects. Relatively well sampled CO SLEDs have been constructed from the ground for bright quasars [@Weiss2007; @bradford09; @Danielson2010], while complete CO SLEDs have been measured by the $Herschel$ $Space$ $Observatory$ in low redshift galaxies [@Panuzzo2010; @vanderWerf2010]. Most SMGs have been observed in only one or two CO lines [see e.g., @Harris10; @Ivison2010; @Aravena2010; @tacconi08; @greve05; @Solomon2005], and their physical properties remain largely unknown. This situation has improved in recent years, with observations of multiple CO lines in individual SMGs [@Ao2008; @Carilli2010; @Lestrade2010; @Riechers2010; @Danielson2010]. The best sampled CO SLEDs show that multiple CO components are required to explain the full line luminosity distribution, where most of the mid-$J$ CO emission can generally be fit by a warm component, with kinetic temperatures of 40-60 K and gas volume densities of 10$^{3}$-10$^{4}$ cm$^{-3}$. However, solutions with kinetic temperatures of a few$\times$100 K and lower densities are also allowed by the data [@Ao2008; @Weiss2007; @Bayet2009], and this region of the parameter space has been insufficiently explored. Z-Spec has an important advantage over other instruments, since it covers a large portion of the CO SLED in a single observation (depending on the redshift) and with a common calibration for the entire bandpass.

This paper describes observations of five H-ATLAS sources undertaken with Z-Spec. Based on the CO emission detected by Z-Spec, we successfully determined the redshifts of four out of five lensed galaxies. The Z-Spec observations are described in Section \[sobs\], followed by the description of the algorithm for redshift determination in Section \[sredshift\]. We use the measured redshift to constrain the spectral energy distribution (SED) of these galaxies, estimating the dust temperature and emissivity index, as well as the total infrared luminosity. We perform an analysis of the partial CO SLEDs, constructed from the lines observed by Z-Spec, to constrain the physical conditions of the molecular gas. The analysis of the galaxy SEDs and CO emission lines is presented in Section \[sgalaxy\], and a summary of our results can be found in Section \[sconcl\]. Throughout the paper we assume a standard $\Lambda$CDM cosmology, with H$_{0}$=71 km s$^{-1}$ Mpc$^{-1}$, $\Omega_{\mathrm{M}}$=0.27, $\Omega_{\Lambda}$=0.73 [@Spergel07].

OBSERVATIONS AND DATA REDUCTION\[sobs\]
=======================================

We selected five high-z candidates among submm-bright galaxies with F(500$\mu$m)$>$100 mJy (Table \[tab:obs\]) from the H-ATLAS survey for follow-up observations with Z-Spec on the 10-m Caltech Submillimeter Observatory (CSO). Z-Spec is a single spatial pixel grating spectrometer with 160 silicon-nitride micro-mesh bolometer detectors (i.e. channels) operating from $190-308$GHz [@naylor03; @earle06; @bradford09]. The frequency response of the Z-Spec channels is approximately gaussian, with a variable FWHM from 720 to 1290 km s$^{-1}$ over the bandpass, that is roughly equal to the channel separation [@earle06]. This follow-up was undertaken to confirm that these sources are lensed by lower redshift galaxy clusters or individual galaxies by measuring their redshifts directly at submm wavelengths. The redshifts of the foreground objects have been measured in the optical, and found to be much lower than the redshifts of the submm galaxies, consistent with the lensing scenario [@Negrello10]. For convenience, throughout the paper we identify our targets by their names used in the SDP H-ATLAS catalogue (SDP.9, SDP.11, SDP.17, SDP.81, and SDP.130).

We carried out the Z-Spec observations of H-ATLAS sources at the CSO from 2010 March 07 - May 14 under generally good to excellent observing conditions, accumulating from 6.8 to 22.5 hour integrations on each target. The zenith opacity at 225GHz (monitored by the CSO tau meter) was $\tau_{225\mathrm{GHz}} = 0.06$ on average, and $\tau_{225\mathrm{GHz}} \le 0.07$ for 75% of the observations. A summary of the observations, including the total integration time on each source is given in Table \[tab:obs\].

The Z-Spec data were taken using the standard “chop-and-nod” mode in order to estimate and subtract the atmospheric signal from the raw data. The secondary mirror was chopped on- and off-source at a rate of 1.6Hz with a chop throw of 90arc-sec while stepping through a 4-part nod cycle which position switches the primary mirror, integrating for 20sec at each nod position. The chopping removes atmosphere fluctuations and the nodding removes instrumental offsets due to imperfect match between the two chopped positions. We checked the pointing every $2-4$ hrs by observing quasars and other bright targets located close in elevation to the H-ATLAS targets, making small (typically $<$ 10arc-sec) adjustments to the telescope pointing model in real-time. We analyze the data using customized software in the same manner as described in @bradford09. For each channel, the nods are calibrated and averaged together, weighting by the inverse variance of the detector noise. Absolute calibration is determined by observations of Mars once per night, which we use to build a model of the flux conversion factor (from instrument Volts to Jy) as a function of each detector’s mean operating (“DC”) voltage [@bradford09]. Since the DC voltage depends on the combination of the bath temperature and the total optical loading on the detectors, we use these curves to determine appropriate calibration factors to apply to each nod individually. Based on the root mean square (rms) deviations of the Mars measurements from the best-fit curves, the channel calibration uncertainties are $3-8$%, excluding the lowest frequencies for which a clean subtraction of the atmosphere is hindered by the pressure-broadened 183GHz atmospheric water line. These uncertainties are propagated through the data reduction. The median rms uncertainties on the final co-added spectra for the H-ATLAS galaxies are listed in Table \[tab:obs\]. These errors do not include the $\sim5$% uncertainty on the brightness temperature of Mars [@wright07]. The calibrated Z-Spec spectra of the five ATLAS galaxies are shown in Figures \[spectra\] and \[spectra17\]. The redshifts of these sources are determined using a custom algorithm, tailored specifically for multiple lines observed simultaneously in the same bandpass. This algorithm is presented in the next section.

REDSHIFT DETERMINATION\[sredshift\]
===================================

Algorithm Description
---------------------

The redshift determination relies on the simultaneous detection of multiple CO and atomic lines. Using the multiple CO lines falling in the Z-Spec bandpass, we developed a redshift-finding algorithm that is capable of handling cases where the signal-to-noise in individual lines is low. The number of CO lines redshifted in the Z-Spec bandpass grows from 2 at $z=0.51$ (CO(3-2) and CO(4-3)) to 4 or more at $z>2$ (starting at CO(5-4) through CO(8-7)). The measurement of redshifts greater than $\sim3$ requires the presence of high-excitation, warm CO gas, where transitions between higher rotational levels become important.

As the width of the Z-Spec channels varies from 720 to 1290 km s$^{-1}$ over the bandpass, larger than most observed line widths, most of the signal from one line will be concentrated in a single channel. The redshift finding algorithm uses two test statistics, $E_{1}$ and $E_{2}$ (Eqs. \[e1eq\] and \[e2eq\]), constructed from combinations of the detection significance in those channels in which a reference line would be observed by Z-Spec at a given redshift $z$. The reference line list contains all CO rotational lines up to the CO(17-16) transition, the \[\] 492.16 GHz line, the \[\] 1458.8 GHz line, and the \[\] 1900.569 GHz line. For redshifts $z>1$ the \[\] 809.342 GHz and CO(7-6) 806.651 GHz lines can fall in the same Z-Spec channel, and are therefore degenerate for the purpose of this procedure. The values of the test statistics are related to the probability that the lines from the reference list, redshifted by a factor of $(1+z)$, are present in the spectrum. The null hypothesis is that there is no signal detected in the spectrum for any value of the redshift, and the distribution of the test statistics for the null hypothesis can be obtained using blank sky spectra, flat continuum sources, or noise simulations. Both statistics are maximized when the redshifted frequencies of the reference lines match the frequencies where the largest signal is present in the spectra, after continuum subtraction. The continuum subtraction uses a fourth degree polynomial to better account for local deviations from a power-law. The use of two statistics instead of one helps reduce the number of false redshifts that can be due to random noise fluctuation in the spectra.

We search a redshift range between 0.6 and 5.6 in steps of 0.001. Let $n(z)$ be the number of reference lines that would fall in the Z-Spec bandpass at redshift $z$. The algorithm loops through all the $z$ values, redshifting all the lines in the line list, and finding the set of $n(z)$ Z-Spec channels corresponding to the lines in the bandpass at that redshift. The two test statistics, $E_{1}$ and $E_{2}$, are evaluated for each redshift using the signal and noise in the set of $n(z)$ channels determined in the previous step. The first test statistic, $E_{1}(z)$, is defined as the ratio of the signal to the noise summed only over the Z-Spec channels that correspond to a line in our list when redshifted to redshift $z$,

$$\label{e1eq}
 E_{1}(z)=\frac{\sum_{i} S_{i}}{\sqrt{\sum_{i} \sigma_{i}^{2}}},$$

where the sum is taken from 1 to $n(z)$, and $S_{i}$ and $\sigma_{i}$ are the continuum subtracted signal and noise, respectively, for the channel corresponding to line $i$. The second test statistic, $E_{2}(z)$, is defined as

$$\label{e2eq}
 E_{2}(z)=\mathrm{median}\{ f_{ij} | f_{ij}=0.5(S_{i}/\sigma_{i}+S_{j}/\sigma_{j}), 1\leq i,j\leq n(z), i < j\},$$

where the set contains all possible pairs of lines in the Z-Spec bandpass at the corresponding redshift. The $E_{2}$ statistic can be easily understood as a rejection criterion. Setting $E_{2}=\epsilon$ as a lower limit is equivalent to discarding all redshifts for which the median value of the average S/N in any pair of channels is less than $\epsilon$. Note that $E_{1}$ would be a reasonable statistic for single line detections, but $E_{2}$ relies on multiple lines being present in the Z-Spec bandpass at a given redshift. For redshifts $<0.6$, $E_{1}$ is still a useful statistic, though the presence of only one line does not make for an unambiguous redshift determination. No a priori knowledge of the relative line strengths is assumed, and therefore the algorithm gives equal weight to all the lines in the search list. As such, the significance of the determined redshift is dependent on which lines are detected in the spectrum relative to the lines that were expected based on the reference list. This is the case for SDP.17b (see Section \[iz\]), where the redshift significance increases greatly if we include the water line on our line list. However, such an extension of the reference line list is not always justified, since the significance for the redshift of another galaxy where the water line is not detected would be unnecessarily diminished. Care must also be taken in using the current line list at higher redshifts, where the high $J$ CO lines are likely to have much lower significance relative to the \[\] and \[\] lines. For example, the strength of the CO lines drops with increasing $J$ for starburst galaxies [e.g., @Danielson2010], but remains relatively constant for AGN-dominated galaxies and quasars [e.g., @bradford09; @vanderWerf2010]. Current knowledge of line strengths and the nature of the sources at high redshift does not warrant including this information in the redshift search.

Testing the Method: Simulations on Blank Sky
--------------------------------------------

The distributions of $E_{1}$ and $E_{2}$, as well as the effect of noise on the positions of their maxima, have been confirmed by simulating multiple realizations of [*blank sky*]{} spectra. The simulations take an input spectrum and generate new realizations by choosing random values of the signal in each channel, normally distributed around the value of the input signal in that channel, with a standard deviation equal to the measurement error in the same channel. The error array associated with the new realization is scaled accordingly, such that the signal-to-noise in each channel is equal to that of the input spectrum. As the input spectrum we used a 6.5 hour integration on blank sky, recorded on May 11 and 12, 2010. Other featureless spectra can also be used for this test. For each realization, we calculate the values of both test statistics and determine the redshifts corresponding to their maxima. The distribution of these redshifts is shown in Figure \[simskyzearch\], derived from a simulation containing 1000 realizations of the sky spectrum. Ideally, this distribution should be uniform, but the input sky spectrum is not perfectly smooth, possibly due to small instrumental artefacts, which are reflected by the probability variations. The noise cross-correlation between channels are present in the initial sky spectrum, but not in the simulated spectra. We find that although the $values$ of $E_{1}$ and $E_{2}$ for each redshift are mostly correlated, the [*locations of their maxima*]{} are not (Pearson correlation coefficient 0.15). Figure \[simskyzearch\] conveys this fact by comparing the maximum joint probability that $E_{1}$ and $E_{2}$ select the same redshift, relative to the maximum probability of a redshift found by $E_{1}$ independently ($\sim$3% vs. $\sim$17%). This shows that the combination of two test statistics is more robust against random fluctuations, considerably reducing the noise floor across the redshift range.

A realization of the two test statistics using a blank sky spectrum is shown in Figure \[skyzearch\]. In this case, the maxima of the two statistics occur for different redshifts, and no value of $E_{2}$ is above the threshold of 2. For noise-dominated spectra, such as blank sky, $E_{1}(z)$ is well fit by a normal distribution with $\sigma_{E1}$=0.97, and $E_{2}(z)$ by a t distribution with 10 degrees of freedom (Figure \[estdist\]). For our redshift determination, we first ask that the maximum value of $E_{1}$ have a significance $\geq3\sigma$. Equivalently, the null is rejected at 99.73% confidence level. However, as noise fluctuations can still lead to random, albeit significant, peaks in these statistics even for blank sky spectra, the redshift confirmation also requires that the maxima of the two test statistics occur at the same redshift. To summarize, a redshift $z_{0}$ is accepted when the following conditions hold:

$$\label{acceq}
 \begin{aligned}
 E_{1}(z_{0})\geq 3\sigma_{E1},\\
 E_{1}(z_{0})=max(E_{1}),\\
 E_{2}(z_{0})=max(E_{2}).
 \end{aligned}$$

For a line pair, the $E_{1}$ 3$\sigma$ threshold corresponds to an average S/N per channel of 2.12. Z-Spec reaches a sensitivity of about 0.5 Jy s$^{1/2}$ per channel for an atmospheric optical depth $\tau_{225}$=0.068 [@inami08]. Combining the 3$\sigma$ threshold criterion with the measured sensitivity of Z-Spec, we estimate that a redshift can be determined in less than 1.5 hours of integration time if the line flux densities per channel are on the order of 15 mJy, but can require more than 12.5 hours if the flux density is less than 5 mJy. For our galaxy sample, the mean integrated CO line flux (Table \[tab:lco\]) is $\sim$18 Jy km s$^{-1}$, while the average width of the channels is 950 km s$^{-1}$. However, the flux density per channel could be only $\sim$10 mJy if the line flux happens to be split between two adjacent channels. In this case, the typical integration time for obtaining a redshift with Z-Spec would be at least 3.5 hours.

Comments on Individual Redshifts \[iz\]
---------------------------------------

The results of applying this algorithm to our galaxy sample are shown in Figure \[A981zearch\]. The redshift value and its uncertainty (Table \[tab:cont\]) are determined from the position and width of the peak of the $E_{1}$ test statistic, and the significance of the redshift peak is tested against the null hypothesis (Figure \[e1hist\]). We fit a gaussian to the peak of the $E_{1}(z)$ statistic, and define the redshift error as the upper limit for the standard deviation of this gaussian. When lines are present in the spectrum, aside from the main peak due to the true redshift, secondary peaks will arise in the $E_{1}$ and $E_{2}$ distributions, corresponding to redshifts where some of the lines in the line list fall on the same channels as the observed lines. The secondary peaks are marked by blue asterisks for each source in Figure \[A981zearch\]. The real redshift will have higher significance than these secondary peaks, since the largest number of lines add their contribution to the total signal in this case. Excluding the true redshift and all the redshifts corresponding to secondary peaks from the test statistics, we are left with the null distribution, as shown in Figure \[e1hist\]. This figure also shows the significance of the redshift determination for each source, by marking the position of the peak value of the test statistic relative to the 3$\sigma$ threshold.

[*SDP.81*]{} The redshift for SDP.81, z=3.037$\pm$0.01, obtained by this method at 3.8$\sigma$ significance on 19 March 2010, was confirmed (z=3.042$\pm$0.001) with follow-up observations with the IRAM Plateau de Bure Interferometer on 23 March 2010 [IRAM/PdBI, @Negrello10; @Neri10] and with an independent blind search on 25 March 2010 by the Zpectrometer instrument at the Green Bank Telescope [GBT/Zpectrometer, @Negrello10; @Frayer10], informed by a concurrent photometric redshift estimate (2.9$^{+0.2}_{-0.3}$). The significance of the test statistics for our redshift determinations shows that the integration time necessary to secure a redshift can be much shorter than the values listed in Table \[tab:obs\], and future submm instruments with better sensitivity will be able to obtain the redshifts of such galaxies even faster.

[*SDP.17*]{} Given the size of the Z-Spec beam (FWHM$\approx$30) and the possible presence of lensing or other foreground structures in the same beam, the observed spectrum could be a combination of features from multiple objects. We choose this interpretation for the spectrum of SDP.17, best described by two components at different redshifts (both listed in Table \[tab:cont\]). The first redshift found by our algorithm is 2.308 (SDP.17b). After fitting the CO lines at this redshift and subtracting them from the spectrum, we perform a second redshift determination, identifying a second component with a redshift of 0.942 (SDP.17a). This combination explains all the features present in the spectrum (see Figure \[spectra17\]), and is consistent with the interpretation of the 299 GHz feature as the restframe 987 GHz water line at a redshift of 2.308. This water line has been seen to be very strong in other AGN and star-forming galaxies at low redshift, such as Mrk231 and Arp 220 [@Gonzales2010], and it has been tentatively detected in the Cloverleaf quasar at z=2.56 by @bradford09. The redshifts of SDP.9 and SDP.17b have recently been confirmed by follow-up observations of the CO(2-1) and (3-2) lines, respectively (L. Leeuw, private communication), with the Combined Array for Research in Millimeter-wave Astronomy (CARMA). The second redshift (SDP.17a) has a much lower significance, but it is in agreement with the photometric and spectroscopic optical redshifts [0.77$\pm$0.13 and 0.9435$\pm$0.0009, respectively @Negrello10]. Alternatively, the peak now identified with the CO(5-4) line at z=0.94 could be a wing of this water line, due to a $\sim$1600 km s$^{-1}$ outflow. The presence of multiple ULIRGs in a single line of sight is intriguing, and is an example of discoveries that can be made possible by Z-Spec’s broad bandwidth.

[*SDP.130*]{} SDP.130 has a redshift of 2.6260$\pm$0.0003, measured by GBT/Zpectrometer [z=6.625$\pm$0.001, @Frayer10], and made more precise with PdBI/IRAM [@Negrello10; @Neri10]. So far, three CO lines have been measured in this galaxy at this redshift, namely the CO(1-0) line observed with the Zpectrometer, and the CO(3-2) and CO(5-4) lines observed with PdBI [@Negrello10], on a tuning that was successfully guided by the sub-mm photometric redshift of $z=2.6^{+0.4}_{-0.2}$ [@Negrello10]. However, we do not detect any of the higher $J$ transitions ($J_{u}>$6) that would fall in the Z-Spec bandpass at this redshift. This non-detection, which places upper limits on the integrated fluxes of the CO (6-5) through (9-8) lines at $<$ 12.5 Jy km s$^{-1}$, suggests a low ($<$ 50 K) gas temperature in the z=2.626 galaxy. We attempted to identify the line at 277 GHz, marked in red in Figure \[spectra\], with the CO (3-2) transition at z=0.25, but that would be inconsistent with the optical spectroscopic redshift of the lensing galaxy [0.220$\pm$0.002, @Negrello10] by more than 7000 km s$^{-1}$, as well as inconsistent with the other observed properties of the foreground galaxy. Similarly, identifying this feature with the 987 GHz water line at z=2.626 would require a velocity offset of $\sim$4200 km s$^{-1}$, and usually the presence of highly excited CO gas, which is not observed. This feature remains unidentified.

[*SDP.9 and SDP.11*]{} The redshifts for these galaxies have been measured at 6.5$\sigma$ and 5.3$\sigma$ significance, respectively. The redshift of SDP.9 has been confirmed by CARMA observations, and more follow-up observations are currently planned for both SDP.9 and SDP.11.

GAS AND DUST PROPERTIES\[sgalaxy\]
==================================

A model including the lines and power-law continuum is fit to each spectrum in Figures \[spectra\] and \[spectra17\], allowing the line intensities, redshift, and continuum slope to vary. The best fit power-law index $\alpha$ for each galaxy is listed in Table \[tab:cont\]. The initial estimate for the redshift is provided by the algorithm described above, and the fit is constrained by the requirement that all the lines be at the same redshift. In cases where some of the lines are blended, we first fit only the unblended lines to obtain a more precise value for the redshift, and then we fit all the lines simultaneously, with the redshift kept fixed, to get the integrated line strengths, listed in Table \[tab:lco\]. Although the lines are not resolved, the signal from one line can be spread among adjacent channels due to the overlap of their frequency responses. We measure only the integrated line strengths, taking into account the frequency response of each Z-Spec channel, weighted according to the line width. On average, line widths below $\sim$1000 km s$^{-1}$ are not resolved by Z-Spec, and we choose a value of 300 km s$^{-1}$ in fitting the integrated line strengths. This value closely matches the width of the lines for SDP.81 and SDP.130 at PdBI [@Neri10], but is relatively low compared to the range found by interferometric measurements of other lensed high-redshift galaxies [@greve05; @Knudsen2009]. However, the determination of the integrated line fluxes is not sensitive to the choice of the line width up to values of the order of the channel width. The largest uncertainties in the integrated line strengths arise in the case of line blending, such as the CO(7-6) and \[\] $^{3}P_{2}\rightarrow^{3}$$P_{1}$ lines, or the overlapping lines at different redshifts in SDP.17 (blended lines are indicated in Table \[tab:lco\]).

Continuum Spectral Energy Distributions\[sseds\]
------------------------------------------------

The continuum data for all 5 galaxies is shown in Figure \[seds\]. The measured continuum flux from Z-Spec is found to be in good agreement with the MAMBO 1.2 mm photometry [@Negrello10], except for SDP.9. Estimates of the total amount of dust and star formation rates in each galaxy can be obtained by fitting their far-infrared (far-IR) to submillimeter spectral energy distribution (SED). In this fit we include the Z-Spec data along with the Herschel-SPIRE and Herschel-PACS photometric points, as well as the Submillimeter Array (SMA) measurements at 880 $\mu$m for SDP.81 and SDP.130 [@Negrello10].

The far-IR rest frame SED can be described by a modified blackbody function, defined as

$$\label{bbody}
F_{\nu}=Q_{\nu}(\beta)B_{\nu}(T_{d})\Omega_{d}=(1-e^{-\tau{(\nu_{0})}(\nu/\nu_{0})^{\beta}})\frac{2h\nu^{3}}{c^{2}}\frac{1}{e^{h\nu/kT_{d}}-1}\Omega_{d}=\frac{L_{IR}}{4\pi d^{2}}\frac{Q_{\nu}(\beta)B_{\nu}(T_{d})}{\int Q_{\nu}(\beta)B_{\nu}(T_{d})d\nu},$$

where $B_{\nu}$($T_{d}$) is the Planck function, $\tau(\nu_{0})=1$ is the optical depth at $\nu_{0}$, $\Omega_{d} $ represents the observed solid angle of the dust emitting region, $d$ is the (known) distance to the source, and $h$ and $k$ denote the Planck and Boltzmann constant, respectively. The fit can be performed with three parameters: $T_{d}$, $\beta$, and a scale factor. The overall scale of the SED can be parametrized either in terms of the solid angle $\Omega_{d}$, or the total infrared luminosity ($L_{IR}$), defined as the integral of the SED from 8 to 1000 $\mu$m (rest frame). The $L_{IR}$ derived in this manner underestimates the true total infrared luminosity, due to the likely presence of warmer dust components that contribute at shorter wavelengths. Including $\nu_{0}$ as a fourth parameter in the fit leads to a value of 1251$\pm$130 GHz for SDP.9, but no strong constraints are found for the rest of the sample. The low value found for $\nu_{0}$, and the observed flattening of the peak of the SEDs in the far-IR suggest that the SEDs of the galaxies in our sample can be modeled either as combinations of multiple graybodies with different temperatures, or as a single graybody with a large optical depth at far-IR wavelengths.

The simplest model that can reproduce the data for the entire sample has fixed $\beta=2$ [@Priddey2001] and $\nu_{0}$=1300 GHz, in agreement to the value found for SDP.9. The dust emissivity index $\beta=2$ is also consistent with the error bars of the Z-Spec spectra. The best-fit models are shown in Figure \[seds\], and the corresponding values for $T_{d}$ are listed in Table \[tab:cont\]. With dust temperatures between 54 and 69 K, the peak of the dust SED is found in a narrow range of wavelengths (73 to 92 $\mu$m) for all lensed galaxies in the sample. It is important to bear in mind that this fitting function for the SED is largely empirical, and the degree to which $T_{d}$ and $\beta$ represent physical quantities is complicated by the spatial averaging over the entire galaxy and the degeneracy between a distribution in dust temperature and a distribution of dust types (represented by $\beta$). The formal errors for the fitted parameters ($T_{d}$) should not be interpreted as errors on physical quantities, due to these caveats. Even though such SED fits could be obtained using just the photometric points, the addition of Z-Spec data not only strongly constrains the continuum slope, but also breaks the degeneracy between $T_{d}$ and redshift [@Blain:1999], by independently determining the latter. Since Z-Spec has determined the redshift, we are able to obtain $T_{d}$ from the continuum fit, which otherwise would constrain the quantity $T_{d}/(1+z)$ [e.g., @Amblard10]. This degeneracy can lead to significant variations in the derived $T_{d}$ if the redshift is not measured independently. The implications of the continuum slope measured by Z-Spec for the dust composition will be analyzed in a future work. Using Eq. \[bbody\] corrected for redshift and the derived $T_{d}$, we can estimate the observed size of the dust emitting region. This solid angle will be affected by the lensing magnification factor. If the dust optical depth at submm wavelengths is low, as is often the case, $\Omega_{d}$ will be correlated with $\tau$, and therefore with $\beta$ [@Hughes1993]. However, we break this degeneracy by fixing $\beta$. The resulting $\Omega_{d}$ ranges between 0.30 arcsec$^{2}$ for SDP.17b and 1.44 arcsec$^{2}$ for SDP.17a. These values may underestimate source sizes that are resolved in the SMA images with a resolution of $\sim$0.8 arcsec at 340 GHz. Using the magnification factors from Table \[tab:cont\], the intrinsic size of the dust emitting region will have an equivalent radius of 0.7 kpc for SDP.81 and 1.3 kpc for SDP.130. Note however, that $\Omega_{d}$ corresponds to the effective solid angle of the dust emitting region, such as the total area of small clumps spread over a larger region. In an image where these clumps are unresolved, the total observed solid angle can appear to be larger. Having estimated the source size, the total dust mass follows from the relationship $\tau(\nu)=\kappa(\nu)M_{d}/D_{A}^{2}\Omega_{d}$, with $\kappa(\nu)=0.4(\nu/250GHz)^{\beta}$ [e.g., @Weiss2007], where D$_{A}$ is the angular diameter distance. The dust mass can also be estimated in the optically thin limit ($1-e^{-\tau}\simeq\tau$) without the additional step of deriving $\Omega_{d}$, by substituting $\tau(\nu)$ directly in Equation \[bbody\]. This is a good approximation at 250 GHz (1.2 mm), in the middle of the Z-Spec bandpass. Calculated in the optically thin limit, the dust mass is a robust estimate of the lower limit for the total dust mass in the galaxy, $M_{d,lim}$. Using the optically thin approximation and the 250 GHz flux density measured by Z-Spec, we derive values for the magnified $M_{d,lim}$ of a few$\times$10$^{9}$ M$_{\odot}$, as listed in Table \[tab:cont\]. If the dust is optically thick, as suggested by $\nu_{0}$=1300 GHz, the calculated $M_{d,lim}$ will underestimate the true dust mass for our galaxy sample by at most 30%. The dust mass is also inversely correlated with the assumed temperature, and will be underestimated when using the dust temperature corresponding to the peak of the SED. This temperature is likely too large to represent the bulk of the dust. Assuming that the 250 GHz flux is partially due to a dust component with a temperature as low as 20 K, and taking into account the optical depth corrections, we estimate that the total dust mass could be larger than $M_{d,lim}$ by up to a factor of $\sim$4. To summarize, with good approximation, the true dust masses for these galaxies will be found in the interval $[1,4]\times M_{d,lim}$. The uncertainties in $M_{d,lim}$ are mostly due to uncertainties in the expression for $\kappa(\nu)$. Note that the quantity $M_{d}/\Omega_{d}$ is proportional to $\tau$ and independent of temperature; for a given $\tau$, a lower limit for $M_{d}$ implies a lower limit for $\Omega_{d}$, but this limit for $\Omega_{d}$ will decrease with increasing $\tau$.

A more realistic approach is to fit the photometric and continuum points with a library of SEDs, taking into account the transmission curve of each instrument. We apply this method to our galaxy sample, using the SED libraries of @chary01 (CE01), and @dale02 (DH02). We find that the IR luminosities derived from the modified blackbody fitting are at most factor of $\sim$2 lower than those when we use the SED libraries, and within 20% from the $L_{IR}$ obtained assuming the models of @dacunha08, calibrated for ULIRGs, with $A_{V}>2$ [@Negrello10]. The variations between the values of $L_{IR}$ obtained by different methods reflect the systematic uncertainties in deriving this quantity. Similar underestimates have been found by others, and are due to the fact that the submm photometry does not measure the warm dust component of the SED, if one is present [@Swinbank2010; @Ivison10]. This is emphasized by the poor fit of the modified blackbody curve to the Herschel-PACS data points, and the significant flux at shorter wavelengths predicted by the CE01 and DH02 models (Figure \[seds\]). Any derivation of $L_{IR}$ is model dependent, with the largest differences arising from the presence of a warm dust component in the SED libraries. In Table \[tab:cont\] we list the $L_{IR}$ values derived from the SED template fitting method, as they represent a more accurate description of the total infrared energy output than the modified blackbody. On average, our $L_{IR}$ values are about 25% lower than those found by @Negrello10, but these differences are difficult to judge without data shortward of 100 $\mu$m. Note that these values are rather smaller than typical IR luminosities of classical SMGs and unbiased sources detected by Herschel surveys.

Except for SDP.17, we attribute all the submm flux density to the high redshift galaxy. The foreground lenses for the other galaxies in our sample have optical properties consistent with being quiescent elliptical galaxies, and are therefore unlikely to have a significant submm emission. We have attempted a decomposition of the SDP.17 SED, using the two measured redshifts and a wavelength-independent scaling factor for each of the two components. The $\chi^{2}$ value for the SED template fits is minimized when the observed flux density is split in half between the two components. This factor has been taken into account in Figure \[seds\] and in deriving the $L_{IR}$ for SDP.17a and SDP.17b, as listed in Table \[tab:cont\]. However, the large dust mass inferred for SDP.17a could be an indication that the SED decomposition between SDP.17a and SDP.17b overestimates the contribution of SDP.17a.

We estimate the star formation rates for our galaxy sample using the conversion factor $SFR $($M_{\odot}$ yr$^{-1}$)=1.5$\times$10$^{-10}$($L_{IR}$/$L_{\odot}$) [@Solomon1997]. Since the selected galaxies are lensed by foreground objects with magnification factors $\sim$10 [@Negrello10], the intrinsic IR and CO line luminosities will be $\sim$10 times lower than the direct conversion from the measured fluxes. SDP.81 and SDP.130 have magnification factors of 25 and 6, respectively, as derived from the best-fit lens model to the high-resolution sub-mm images available for these two objects [@Negrello10]. In Tables \[tab:cont\] and \[tab:sfr\], we left the quantities affected by gravitational lensing magnification unmodified, for reference, but the presence of this contribution is indicated by the letter $\mu$ in front. Based on model predictions [@negrello07], a typical amplification factor of 10 can be applied to these values. Once corrected for magnification, the infrared luminosities and corresponding SFRs are those typical of Ultra Luminous Infrared Galaxies (ULIRGs).

CO Line Luminosities and Spectral Energy Distributions\[sco\]
-------------------------------------------------------------

The measurements of CO lines reveal important information about the physical properties and excitation conditions of the molecular gas, as well as the total gas budget in these galaxies. These parameters can be used to investigate the link between star formation and gas properties. Higher gas temperatures and lower densities would suggest that the star formation has been quenched, while lower temperatures associated with higher densities would show that most of the gas in the galaxy is still collapsing and forming stars. This follows from the fact that the Jeans mass increases rapidly at higher temperatures and low densities.

A useful quantity describing the CO lines is the velocity-integrated brightness temperature scaled by the area of the source, $L_{CO}^{'}$, in units of K km s$^{-1}$ pc$^{2}$. In what follows, the brightness temperatures are computed in the Rayleigh-Jeans limit. If the CO is thermalized and the lines are optically thick, $L_{CO}^{'}$ will be the same for all rotational transitions for which the Rayleigh-Jeans approximation holds. This quantity is traditionally derived from the CO(1-0) transition and related to the total molecular mass via the empirical relation $M_{gas}=\alpha L_{CO}^{'}$, where $\alpha$ is 4.6 M$_{\odot}$ (K km s$^{-1}$ pc$^{2}$)$^{-1}$ for the Galaxy [@Solomon1997], and 0.8 M$_{\odot}$ (K km s$^{-1}$ pc$^{2}$)$^{-1}$ for ULIRGs [@Solomon2005; @tacconi08]. Following @Solomon2005, we can use the latter value for $\alpha$ and the $L_{CO}^{'}$ for the lowest observed CO transition to determine the gas masses. This procedure assumes that all transitions from CO(1-0) up to the lowest observed are thermalized, which might not necessarily be the case. A recent comparison of the CO(3-2) and (1-0) lines [@Harris10] shows that the ratio of the brightness temperatures for these two lines averages to 0.6 rather than 1, due to the presence of multi-phase CO gas. Moreover, the mid-$J$ CO transitions do not account for the possible presence of a colder gas component, making the $M_{gas}$ derived in this manner a lower limit for the total gas mass in the galaxy. Assuming that the lines with $J_{u}>$3 are thermalized, corresponding to a warmer gas component, we apply this correction factor to the lowest CO transition measured, and obtain the gas masses listed in Table \[tab:sfr\]. These values result in an average molecular gas-to-dust ratio for the lensed galaxies of 34$\pm$10, independent of magnification. The mean does not include the foreground SDP.17a, and is in agreement with the values found for other SMG samples [@Kovacs2006; @Michalow2010; @Santini2010].

Using the SFRs derived from the IR luminosities, the gas reservoir probed by CO implies a gas depletion time \[$M_{gas}$/$L_{IR}$\] in these objects of $\sim10^{7}$ years, similar to other known SMGs [@Solomon2005; @greve05]. This can be interpreted as the starburst lifetime under the assumptions of constant SFR and no gas inflow. Note that this estimate of the gas depletion time is independent of lensing magnification. The star formation efficiency can be expressed directly in terms of $L_{IR}/L_{CO}^{'}$, without the need for a gas mass conversion factor. After accounting for the lensing magnification factor, $L_{IR}$ and $L_{CO}^{'}$ for our sample follow the same relationship as other SMGs and ULIRGs [@greve05; @Wang2010], within the scatter. In order to derive the physical characteristics of the gas in these galaxies, including the gas temperature, density, pressure, and CO column density, we need measurements of multiple CO transitions, sampling the rotational ladder as well as possible. The spectral line energy distribution (SLED) for the CO molecule has been constructed in a few cases for nearby and low redshift galaxies [e.g., @Panuzzo2010]. In Figure \[sleds\] we show the partial SLEDs for our galaxy sample, constructed from the lines detected in the Z-Spec bandpass. This plot favors a distribution with the brightest lines between CO(5-4) and (7-6), similar to the distribution observed for other SMGs and starburst galaxies [@Weiss07; @Danielson2010]. The shape of the line luminosity distribution does not reflect only the CO excitation temperature, but also the gas density, and the effects of the optical depth at the line frequency [@Goldsmith1999; @Papadopoulos2010b]. In the optically thin limit, the CO column density scales with the absolute value of the line intensity, assuming that the source size and the magnification factor are known. Under the assumption of local thermodynamic equilibrium (LTE), all CO transitions have the same excitation temperature, $T_{ex}$, also equal to the gas kinetic temperature $T_{kin}$, signifying that all rotational levels are populated according to the Maxwell-Boltzmann distribution at temperature $T_{ex}$. In Section \[slte\] we estimate these parameters by fitting the partial SLEDs, using the relationship between the integrated line brightness temperature, column density and excitation temperature, under LTE. Although this case is limiting due to the assumption of constant $T_{ex}$ for all levels, it is interesting to compare the predictions of this model to the more general non-LTE models, given its simple physical interpretation. In the non-LTE case, presented in Section \[snlte\], the models involve a larger number of parameters, and are less well constrained. We use RADEX [@vanderTak07] to compute the brightness temperatures of the CO lines and estimate the likelihood distribution over the parameter space. These distributions allow us to assess if the available data are able to distinguish between the LTE and non-LTE models.

### LTE models\[slte\]

The integrated line flux $S_{\nu}\Delta$v (in Jy km s$^{-1}$) in the observer’s frame is related to the velocity-integrated Rayleigh-Jeans source brightness $W(J)$ by [e.g., @Solomon1997]

$$\label{ttojy}
W(J)=\frac{\lambda^{2}_{J,J-1,rest}(1+z)^{3}}{2k \Omega_{a}}S_{\nu}\Delta\mathrm{v}\frac{\Omega_{a}}{\Omega_{s}}$$

where $\Omega_{s}$ and $\Omega_{a}$ are the solid angles of the source and the antenna, respectively. $W(J)$ is in units of K km s$^{-1}$. The last fraction represents the inverse of the beam filling fraction. The contribution of the gravitational lensing magnification should cancel out in this expression, as it contributes to both $S_{\nu}$ and $\Omega_{s}$, but the true $\Omega_{s}$ is not known. In principle, the same approach taken for the continuum (Section \[sseds\]) could be used to determine the source size. However, such a fit requires a minimum of three parameters, and will not be well constrained by the number of CO lines in our SLEDs. In addition, the optical depth depends directly on the column density, and cannot be estimated independently, in the same way that the dust optical depth was determined by the continuum slope. We assume an intrinsic source size of $\sim$2 kpc, consistent with the angular diameter of 0.2 found for an SMG at a redshift of 2.3259 [@Swinbank2010; @Negrello10], and similar to the size of the dust emitting region found in Section \[sseds\]. The corresponding beam filling fractions are listed in Table \[tab:sfr\]. As this source solid angle now represents the intrinsic size, and not the magnified one, we must correct the observed flux densities by the lensing magnification factors. We use the values listed in Table \[tab:cont\], when available, and assume a value of 10 in all other cases. The case of SDP.17a is treated differently, as it is assumed to be a foreground galaxy, not affected by gravitational lensing. For the intrinsic size of SDP.17a, we use a value of 1.54 arcsec$^{2}$, which approximates the size of the optical image. @Negrello10 identify two galaxies in the i-band image of SDP.17 and fit both light distributions with the GALFIT software. As the presence of two galaxies could indicate a possible merger, we choose the source size of SDP.17a to be the sum of the areas of these two galaxies. The distribution of the velocity-integrated brightness temperatures for the CO lines can be constructed starting from the CO column density and gas temperature, under the assumption of LTE. Following @Goldsmith1999, the velocity-integrated Rayleigh-Jeans source brightness is given by

$$\label{btemp}
W(J)=N_{J}\frac{hc^{3}A_{J,J-1}}{8\pi k\nu^{2}}\frac{1-e^{-\tau_{J,J-1}}}{\tau_{J,J-1}},$$

where $\tau_{J,J-1}$ is the line center optical depth, and $A_{J,J-1}$ is the Einstein $A$ coefficient for the transition. In LTE, the column density of molecules in the upper level, $N_{J}$, is related to the total column density $N$, by $$\label{nupper}
N_{J}=\frac{N}{Z}g_{J}e^{-E_{J}/kT_{ex}},$$

where $Z$ is the partition function, $E_{J}$ is the energy of level $J$, and $g_{J}=2J+1$ is the degeneracy of level $J$. The line center optical depth can be expressed as a function of column density, temperature, and line width $\Delta$v as

$$\label{tauline}
\tau_{J,J-1}=A_{J,J-1}\frac{c^{3}}{8\pi\nu^{3}\Delta\mathrm{v}}N_{J}(e^{h\nu/kT_{ex}}-1).$$

We fit Equation \[btemp\] to the measured $W(J)$ distribution, with the column density and gas temperature as free parameters, and $\Delta$v = 300 km s$^{-1}$. We find that the best fit models have relatively low optical depths ($\lesssim 1$) such that the choice of the line width has only a small effect on the fitted parameters. For the lensed galaxies, the measured CO SLEDs and the range of SLEDs allowed by the formal 1$\sigma$ interval for the gas temperature are are shown in Figure \[sleds\]. The SLEDs can be characterized by an overall scale and line ratios. The scale of the observed SLEDs is mainly a result of the CO column density and the beam filling fraction, while the line ratios depend on the CO temperature and gas (H$_{2}$) density. The parameters in each pair are therefore largely degenerate and anti-correlated. This degeneracy is characteristic to CO and other molecular SLEDs, regardless of galaxy type. The last correlation (between temperature and gas density) only exists until LTE is reached, and the temperature becomes fixed. By making assumptions on the beam filling fraction and gas density, we can place limits on the remaining parameters. The error bars on the column densities derived in this manner are correlated with the errors in the beam filling fraction, which are not known. Similarly, by making the assumption of LTE for all transitions up to CO (7-6), we are constraining the gas density to be greater than the critical density for this transition ($n$\[H$_{2}$\] $\gtrsim$ 3$\times$10$^{5}$ cm$^{-3}$). At densities $n$\[H$_{2}$\] $\gtrsim10^{6}$ cm$^{-3}$, considerably larger than the average value observed in Galactic molecular clouds, all observed lines should be in LTE. Values of the gas density more typical for Galactic molecular clouds (10$^{3}$-10$^{4}$ cm$^{-3}$) correlate with higher gas temperatures, of a few hundred degrees, in order to reproduce the observed line ratios.

The best-fit LTE CO column densities are $\sim$few$\times$10$^{18}$ cm$^{-2}$, and the gas temperature ranges between 41 and 115 K, as listed in Table \[tab:sfr\]. Taking into account the assumed source size, we estimate total CO masses of a few$\times$10$^{6}$ M$_{\odot}$, only $\sim$10$^{-4}$ of the total gas mass. Together with the large pressure associated with this model ($\sim$10$8$ K cm$^{-3}$), this suggests that the LTE scenario can only describe a small portion of the gas. Other regions of the parameter space are associated with non-LTE gas excitation, explored with the RADEX modeling in the next section.

### Non-LTE radiative transfer models of CO line excitation\[snlte\]

In general, the rotational levels of the CO molecule might not be populated according to a single temperature, and the gas is not necessarily in equilibrium with the radiation field. By dropping the LTE assumption, we allow the excitation temperature to be a function of transition, being determined by the level populations for each line, while the kinetic temperature will be the global quantity describing the thermal energy of the gas. The level populations are found by solving the detailed balance equations including both radiative and collisional rates, and the output intensities are calculated by solving the radiative transfer equations. Usually, these equations are strongly coupled, involving large spatial and frequency grids, and further complicated by the number of molecules and transitions involved. Simplifying assumptions are usually made to reduce the computing time, depending on the problem at hand.

We use RADEX to estimate the range of physical parameters consistent with the measured line strengths when dropping the LTE assumption. RADEX is a one dimensional, non-LTE radiative transfer code, that solves for the level populations iteratively, employing the escape probability approximation for the radiative transfer [@vanderTak07]. The medium is assumed homogeneous and isothermal, and the number, type, and abundance of the participating molecules is selectable by the user. The input parameters are the kinetic temperature, $T_{kin}$, the number density of molecular hydrogen, $n$\[H$_{2}$\], as the collisional partner, and the column densities per unit line width of the participating molecules, only CO in our case. The background radiation field is the CMB, redshifted according to the redshift of each galaxy. The output contains the predicted line excitation temperatures, optical depths, and line intensities. The output line fluxes are scaled by an additional factor $\phi$, that represents fractional corrections to the size of the emitting region and to the gravitational lensing magnification factor. It would correspond to the area filling fraction of the emitting region, if the size and lensing magnification factor of the source were known precisely. A value $\phi>1$ would suggest that the assumed source size was underestimated. We compare the measured flux densities with the line intensities output by RADEX using the values for source sizes, line widths, and lensing magnification factors assumed in Section \[slte\] for the LTE model. For the case $\phi=1$ and $n$\[H$_{2}$\]$\gg n_{crit}$ for all transitions, RADEX should recover the LTE SLED as determined from $T_{ex}$ and N\[CO\] in Section \[slte\]. This set of parameters represents the model m01, listed in Table \[tab:radex\] for each galaxy. The high $n$\[H$_{2}$\] chosen for this model (10$^{6}$ cm$^{-3}$) insures that all measured lines will be in LTE. The SLED computed by RADEX for m01 is shown by the blue line in Figure \[sleds\], showing a good agreement with the LTE method described in the previous section. We run RADEX for a range of input models, parametrized by $T_{kin}$, $N$\[CO\], $n$\[$H_{2}$\] and $\phi$, and compute the likelihood density function for all models following the method described in @Ward03. Weak priors are set to rule out unphysical solutions, keeping the total molecular mass smaller than the dynamical mass, and the length of the CO column smaller than the physical size of the galaxy [@Ward03; @Panuzzo2010]. The dynamical mass cut-off is estimated choosing the line width of 300 km s$^{-1}$. We also impose a limit for the kinetic temperature at 3000 K, where collisional dissociation of CO starts to dominate, weakly dependent on the gas density.

We map the surface of the likelihood distribution and determine the location of its maximum by running a Markov chain Monte Carlo (MCMC) algorithm, described in detail in @Scott10. The 2D marginal probability contours obtained from the MCMC algorithm are shown in Figures \[tnmcmc\] and \[nnmcmc\], with the position of the 4D maximum likelihood indicated by the dotted line, and the parameters of the m01 model, equivalent to the LTE solution, shown by the dashed lines. Note that due to projection effects the parameters corresponding to the maximum of the likelihood in the 4D parameter space are not necessarily the same as the coordinates of the maximum for the 2D marginal probability distributions. The set of parameters that maximizes the 4D likelihood for each galaxy is listed in Table \[tab:radex\] as model m02, and the line luminosities predicted by this model are shown in blue in Figure \[sleds\]. The 68% credible regions are calculated as the smallest intervals containing 68% of the 1D marginal probability for each parameter, around the value corresponding to the 4D maximum likelihood. The 4D probability distributions are highly non-gaussian, and the maxima of the marginalized distributions are not a good representation of the true maximum of the 4D distribution.

Due to the aforementioned degeneracies (see Section \[slte\]), the product between the kinetic temperature and gas density on one hand, and CO column density and $\phi$ on the other hand, are better constrained than individual parameters. These products lead to the values for gas pressure and total gas mass listed in Table \[tab:radex\]. Taking into account the lensing magnification factor, the inferred gas mass for non-LTE models is in good agreement with the mass derived from $L_{CO}^{'}$ (Table \[tab:sfr\]). The credible regions for $\phi$ suggest that the size of the emitting region could be larger than assumed for the LTE models, which would also correspond to a larger area characterized by lower gas density and pressure than the LTE case.

The properties of these models can be compared by calculating the $total$ CO luminosity, summed over all transitions in the model, as a proxy for the indicators of the star formation rate and efficiency. The ratios $L_{IR}$/$L_{CO}^{'}$ for each model listed in Table \[tab:radex\] suggest that the star formation efficiency is higher for the LTE models. Using the correlation between total $L_{CO}$ and $L_{IR}$ derived for a mixed sample of nearby and high redshift galaxies by @Bayet2009, we find that the total $L_{CO}$ obtained by integrating the models will overpredict the measured $L_{IR}$ in both cases, but are consistent within the scatter. The deviation between the predicted and measured $L_{IR}$ is larger for model m02, slightly disfavoring this solution. While the LTE models tend to underpredict the total gas masses, the non-LTE models potentially underpredict the total infrared luminosity.

The region of the parameter space that is most consistent with the observed line strengths is enclosed by the likelihood contours in Figures \[tnmcmc\] and \[nnmcmc\]. We see that the LTE solution implies large gas densities ($n$\[H$_{2}$\]), as expected, but a class of non-LTE solutions can be found in the region defined by the likelihood contours in parameter space. The likelihood space roughly splits into high density/low temperature, and low density/high temperature solutions. One possible additional complication to the interpretation can arise from high dust optical depths @Papadopoulos2010b, not included in our models, that can lead to the suppression of mid-$J$ CO lines and to an underestimate of the excitation temperature. However, the likelihood distribution is relatively shallow over the whole region, reflecting the insufficient amount of information in our data. Comparatively, other studies of high redshift SMGs find a warm CO component with $n$\[H$_{2}$\] around 10$^{4}$ cm$^{-3}$ and temperatures between $\sim$40 and 60 K [@Riechers2010; @Carilli2010; @Danielson2010], a region marginally allowed by our contours.

The constraints on the parameter space for the non-LTE models are weak, as expected given the average sampling of the SLED, and cannot distinguish at this point between the LTE and non-LTE scenarios. To emphasize the insight gained by including additional lines in the fit, we add to the SLED of SDP.81 the CO(1-0) integrated flux from @Frayer10. A likelihood analysis for the new set of lines results in the best fit parameters listed in Table \[tab:radex\] as m03. The 2D marginalized likelihoods for this case are shown by the light grey contours in Figures \[tnmcmc\] and \[nnmcmc\]. The tightening of the likelihood contours is substantial with just one line added to the data, and the LTE region of the parameter space becomes less favored. However, the limitation of this model is that it assumes a single gas component, while most of the emission in the CO(1-0) line could be originating from cold molecular gas.

The brightness temperatures predicted by the RADEX models m01 and m02 are shown in Figure \[tempsone\], to emphasize the large deviations beween the predictions of the two models, especially for lower $J$ transitions. The measured data points have been scaled by the lensing magnification factors listed in Table \[tab:cont\] when available, and by a factor of 10 in all other cases. This figure shows that the constraints on the model parameters can be tightened by measurements of lower $J$ transitions, especially the CO(1-0) line. Even if most of the CO(1-0) emission comes from a colder gas component, using this value as an upper limit will help rule out some regions of the parameter space, as in our example for SDP.81. Distinguishing between the different regions in parameter space will clarify the state of the ISM in these galaxies, and thus their star formation histories. Specifically, hot/low-density gas may signal the action of a feedback process quenching further star formation, by increasing the Jeans mass. This high temperature/low density solution has not been fully investigated, but recent studies show that other CO SLEDs can be consistent with it [@Panuzzo2010; @Weiss2007; @Ao2008; @Bayet2009]. The CO SLED in M82 is fit by a small CO component with a kinetic temperature of almost 600 K [@Panuzzo2010], while solutions with $T_{kin}$ of a few$\times$100 K are found by @Bayet2009, and can be allowed by the LVG models for IRAS F10212+4724 [@Ao2008] and APM 08279+5255 [@Weiss2007]. Such temperatures suggest energy input from outflows or AGN activity. The presence of an AGN component in SDP.17b is supported by the relatively flat SLED from CO(6-5) to CO(8-7), similar to the Cloverleaf quasar or Mrk231 [@bradford09; @vanderWerf2010], and the emission line of water, also observed in galaxies with an AGN component, such as Mrk231 [@Gonzales2010].

CONCLUSIONS\[sconcl\]
=====================

Far-IR / submillimeter-wave surveys are revealing submillimeter-bright galaxies from the first half of the Universe by the tens of thousands, but their detailed study requires spectroscopic redshift measurements. We have studied a sample of the brightest sources and have demonstrated a new redshift-measurement technique with our broadband millimeter-wave grating spectrometer, Z-Spec. Z-Spec measures multiple rotational transitions of carbon monoxide, the dominant coolant of molecular gas in galaxies, and thus is not dependent on optical counterparts which are often absent or hard to identify, as is the case for these galaxies. We find redshifts ranging roughly between 1 and 3, in line with previous determinations [@chapman05], reaching back to an era when the Universe was 15% is present age. Their fluxes are proven to be amplified by gravitational lensing [@Negrello10], making them ideal targets for spectroscopic follow-ups. From the observed CO line luminosities and integrated $L_{IR}$, typical conversion factors reveal that these galaxies each house roughly 10$^{10}$ M$_{\odot}$ of molecular gas, and have SFRs between 10$^{2}$ and 10$^{3}$ M$_{\odot}$yr$^{-1}$, after correcting for lensing magnification. Regardless of the magnification details, we are clearly witnessing a rare episode of rapid star formation in these galaxies, since the timescale over which the observed luminosity can be generated by converting the inferred mass of gas into stars is only a few tens of millions of years (depending on the details of the star formation and the accretion of more gas), which is a small fraction of the Universe’s age even at this early epoch. We estimate that the dust masses in our sample of lensed galaxies are around a few$\times$10$^{8}$ M$_{\odot}$, and the wavelengths corresponding to the peaks of their dust SEDs fall within a narrow range, between 73 and 92 $\mu$m in the rest frame. For this initial set of lensed submm galaxies both the dust properties derived from the IR SED, and the physical conditions of the molecular gas probed by the CO lines, are broadly comparable to those in known SMGs [@greve05; @Solomon2005; @Casey09], with excitation temperatures in the 30-120 K range, and $L_{CO}^{'}$/$L_{IR}$ between 1 and 3$\times$10$^{-3}$ K km s$^{-1}$ pc$^{2}$/L$_{\odot}$, as measured from the mid-$J$ CO lines.

The partial SLEDs for the CO molecule constructed from the lines observed by Z-Spec cannot distinguish between different models of CO excitation. The simplest assumption is that of local thermodynamic equilibrium (LTE), under which we can derive the gas column density and excitation temperature. We find that the relative line strengths can be reproduced by relatively low excitation temperatures ($<$ 100 K), and optical depths ($<$ 1). In the non-LTE case, other parts of the parameter space are allowed, including higher optical depths, while measurements of the lower rotational transitions are essential in confirming such models.

By being able to characterize galaxies that can be inaccessible at other wavelengths, the combination of large-area submm surveys and spectroscopic follow-ups of the CO emission lines will lead to substantial progress in our understanding of high redshift galaxies and their evolution. These results suggest the possibility of a rapid growth in our understanding of SMGs, independent of their optical or radio counterparts, but enabled by strong gravitational lensing magnification.

We are indebted to the staff of the Caltech Submillimeter Observatory for their unflagging support. This work was supported by NSF grant AST-0807990 to J. Aguirre and by the CSO NSF Cooperative Agreement AST-0838261. Support was provided to J. Kamenetzky by an NSF Graduate Research Fellowship. Z-spec was constructed under NASA SARA grants NAGS-11911 and NAGS-12788 and an NSF Career grant (AST-0239270) and a Research Corporation Award (RI0928) to J. Glenn, in collaboration with the Jet Propulsion Laboratory, California Institute of Technology, under a contract with the National Aeronautics and Space Administration. We acknowledge Peter Ade and his group for their filters and Lionel Duband for the 3He / 4He refrigerator in Z-Spec, and are grateful for their help in the early integration of the instrument. R. L. wishes to thank Tom Loredo for useful discussions regarding the significance of the redshift determination. We appreciate the help of Robert Hanni and Jon Rodriguez with observing.

[cccccc]{} H-ATLAS J090740.0-004200 &SDP.9 & Apr 27 - May 14 & $0.05-0.21$ & 10.6 & 4.0\
H-ATLAS J091043.1-000322 &SDP.11 & Apr 28 - May 4 & $0.05-0.18$ & 6.8 & 5.5\
H-ATLAS J090302.9-014128 &SDP.17 & Mar 28 - Apr 1 & $0.04-0.08$ & 18.2 & 2.9\
H-ATLAS J090311.6+003905 &SDP.81 & Mar 7 - Mar 12 & $0.02-0.05$ & 22.5 & 2.3\
H-ATLAS J091304.9-005344 &SDP.130 & Mar 21 - Mar 22 & $0.04-0.08$ & 8.6 & 4.4\

[cccccccccc]{} SDP.9 & ... & 1.577$\pm$0.008 & 1.1e-11 & 4.4$\pm$0.6 & 57$\pm$1 & 3.8$\pm$0.2 & 2.5& 0.65 & 6.6$\pm$0.9\
SDP.11 & ... & 1.786$\pm$0.005 & 3.1e-8 & 7.8$\pm$0.9 & 69$\pm$1 & 5.7$\pm$0.4 & 1.7 & 0.43 & 11.7$\pm$1.4\
SDP.17a & ... & 0.942$\pm$0.004 & 1.0e-4 & 0.4$\pm$0.1 & 27$\pm$1 & 2.9$\pm$0.1 & 4.9 & 1.44 & 0.6$\pm$0.2\
SDP.17b & ... & 2.308$\pm$0.011 & 1.9e-7 & 3.9$\pm$1.0 & 66$\pm$1 & 2.9$\pm$0.1 & 1.1 & 0.30 & 5.8$\pm$1.5\
SDP.81 & 18-31 & 3.037$\pm$0.010 & 6.7e-5 & 6.4$\pm$0.3 & 58$\pm$1 & 3.2$\pm$0.1 & 2.2 & 0.69 & 9.6$\pm$0.5\
SDP.130 & 5-7 & 2.626$\pm$0.0003 & N/A & 4.3$\pm$0.2 & 55$\pm$1 & 2.7$\pm$0.3 & 1.6 & 0.47 & 6.5$\pm$0.3\

[ccccccc]{} CO (4-3) & 461.041 & ... & ... & 15$\pm$9 &... & ...\
CO (5-4) & 576.268 & 25$\pm$5 & 23$\pm$8 & 29$\pm$9 &... & ...\
CO (6-5) & 691.473 & 33$\pm$7 & 29$\pm$10 & ... &17$\pm$5 & ...\
CO (7-6) & 806.652 & ... & 18$\pm$14 & ... &11$\pm$7 & 12$\pm$4\
CO (8-7) & 921.800 & ... & ... & ... &16$\pm$6 & 5$\pm$3\
CO (9-8) & 1036.91 & ... & ... & ... &... & 6$\pm$3\
CO (10-9) & 1151.99 & ... & ... & ... &... & $<$ 6.5\
\[\] $^{3}P_{1}\rightarrow^{3}$$P_{0}$ & 492.160 & $<$ 11 & ... & $<$ 6 &... & ...\
\[\] $^{3}P_{2}\rightarrow^{3}$$P_{1}$ & 809.342 & ... & 31$\pm$14 & ... &13$\pm$7 & $<$ 6.5\
H$_{2}$O $2_{0,2}-1_{1,1}$& 987.914 & ... & ... & ... &19$\pm$7 & ...\

[ccccccccc]{} SDP.9 & 13$\pm$3 & 1.7 & 2.6 & 19$\pm$4 & 76$\pm$24 & 0.3 & 0.251 & 1.4$\pm$0.3\
SDP.11 & 15$\pm$5 & 2.0 & 1.7 & 47$\pm$14 & 41$\pm$6 & 0.3 & 0.983 & 3.5$\pm$1\
SDP.17a & 4$\pm$3 & 0.5 & 8.3 & 2.3$\pm$0.7 & 114$\pm$99 & 1.3 & 0.015 & 0.17$\pm$0.05\
SDP.17b & 12$\pm$3 & 1.6 & 2.8 & 19$\pm$6 & 65$\pm$10 & 0.3 & 0.298 & 1.4$\pm$0.4\
SDP.81 & 10$\pm$3 & 1.3 & 1.4 & 11$\pm$3 & 53$\pm$5 & 0.6 & 0.149 & 0.8$\pm$0.2\

[cccccccc]{} SDP.9 &&&&&&&\
m01 & 76. & 18.27 & 6.20 & 1.00 & 8.1 & 8.8 & 55\
m02 & 2628. & 20.50 & 3.37 & 0.05 & 6.8 & 9.7 & 29\
68% credible region & 149-2631 & 19.03-21.10 & 3.37-6.06 & 0.03-0.58 & 6.29-8.64 & 8.71-10.14 & 12-114\
SDP.11 &&&&&&&\
m01 & 41. & 18.67 & 6.20 & 1.00 & 7.8 & 9.2 & 78\
m02 & 21. & 19.74 & 6.50 & 1.66 & 7.8 & 10.5 & 67\
68% credible region & 21-611 & 18.42-20.17 & 3.61-6.55 & 0.19-1.82 & 5.85-7.83 & 9.14-10.50 & 27-520\
SDP.17b &&&&&&&\
m01 & 63. & 18.32 & 6.20 & 1.00 & 8.0 & 8.8 & 55\
m02 & 1988. & 20.03 & 1.99 & 0.85 & 5.3 & 10.5 & 21\
68% credible region & 136-1992 & 19.11-21.03 & 1.99-5.14 & 0.04-1.02 & 5.29-7.59 & 8.93-10.51 & 0-626\
SDP.81 &&&&&&&\
m01 & 53. & 18.04 & 6.20 & 1.00 & 7.9 & 8.6 & 71\
m02 & 213. & 20.33 & 2.63 & 0.38 & 5.0 & 10.4 & 25\
68% credible region & 63-855 & 18.53-20.43 & 2.55-5.17 & 0.05-1.12 & 4.95-7.24 & 8.44-10.44 & 8-640\
m03 & 481. & 16.71 & 3.20 & 130. & 5.89 & 9.35 & 18\
68% credible region & 368-2749 & 14.85-17.08 & 2.56-3.42 & 35-6600 & 5.80-6.04 & 9.22-9.43 & 12-55\

![The Z-Spec spectra of four submillimeter bright H-ATLAS galaxies. The fit to the continuum and CO lines at the measured redshift is overplotted in red, and the positions of the strongest lines falling in the Z-Spec bandpass are indicated by the vertical blue lines. The line indicated in red in the spectrum of SDP.130 is unidentified. \[spectra\]](f1.eps)

![The Z-Spec spectrum of the H-ATLAS source SDP.17. The fit to the continuum and CO lines at $z=0.94$ is overplotted in red in the upper panel, and the rotational CO lines are indicated by the vertical blue lines. These lines have been subtracted from the spectrum shown in the lower panel. The red line in the lower panel shows the fit including the lines identified at $z=2.308$. \[spectra17\]](f2.eps)

![ Blank sky probability distributions for the positions of the $E_{1}$ and $E_{2}$ maxima as a function of redshift. The redshift bin size is 0.003, corresponding to the average width of the Z-Spec channels. These probabilities have been obtained by running the redshift finding algorithm on 1000 realizations of the blank sky spectrum. The highest probability corresponds to the redshift most likely to be determined by each statistic. The last panel shows the probability distribution when the same redshift is found jointly by $E_{1}$ and $E_{2}$. The maximum probability that a false positive is obtained simultaneously by the two test statistics is $\sim$3%. \[simskyzearch\]](f3.eps)

![Null test for the redshift finding algorithm, using a blank sky spectrum. Note that the values of the $E_{2}$ test statistic are systematically below 1.5, and the positions of the maxima of the two test statistics do not coincide. \[skyzearch\]](f4.eps)

![Distributions of the two test statistics derived from blank sky spectra. The histograms for $E_{1}$ and $E_{2}$ are shown with a dashed and continuous line, respectively. The fits with a normal distribution for $E_{1}$ and a t distribution for $E_{2}$ are overplotted as the smooth red curves. \[estdist\]](f5.eps)

![Results of running the redshift-finding algorithm for all the H-ATLAS sources in our sample. The $E_{2}$ test statistic has been offset vertically by 7 units, for clarity. The blue asterisks show the positions of the largest secondary peaks arising from coincidences with the lines from the actual redshift (see text). These peaks contain the same information as the main peak. In the SDP.17 panel, we note the extra peaks that do not match the secondary peaks corresponding to the first selected redshift. The SDP.17a panel shows the determination of the second redshift from the same spectrum, after subtracting the high-redshift component. No redshift is determined for SDP.130. \[A981zearch\]](f6.eps)

![$E_{1}$ distribution for each galaxy for which we determined the redshift. The original $E_{1}$ distribution is shown in black, the red distribution being obtained after subtracting the main peak of $E_{1}$ and all its secondary peaks. The normal distribution fit to the null $E_{1}$ distribution is shown in blue. The red dashed line shows the 3$\sigma$ position for the null gaussian, and the black dashed lines show the positions of the maximum of the test statistic for the measured redshifts, to emphasize the significance of our detections. \[e1hist\]](f7.eps)

![The best-fit SED models for the five H-ATLAS galaxies in our sample. The continuous line shows the modified blackbody spectrum with $\nu_{0}=1300$ GHz and $\beta=2.0$, while the dotted and dashed lines show the SEDs obtained from the SED libraries of CE01 and DH02, respectively. The total infrared luminosities are calculated as the average between the CE01 and DH02 SED template fits, to account for emission above the blackbody spectrum at higher frequencies. The parameters for the modified blackbody fits are also listed in Table \[tab:cont\].\[seds\]](f8.eps)

![Spectral line energy distributions, uncorrected for gravitational lensing magnification. The Z-Spec measurements are shown connected by the black histogram. The data point for the CO(1-0) line measured by @Frayer10 in SDP.81 falls at the bottom of the panel, and is better seen in Figure \[tempsone\]. The red lines show the SLEDs predicted by the best-fit LTE model (continuous), and the LTE models corresponding to the limits of the 1$\sigma$ standard confidence interval for $T_{ex}$ determined from the fit (dashed). The blue line corresponds to the SLED predicted by RADEX using the same temperature and column density as the LTE model, while the green line shows the SLED predicted by RADEX with the parameters given by the 4D maximum likelihood solution. The parameters for these two RADEX models, m01 and m02, respectively, are listed in Table \[tab:radex\]. \[sleds\]](f9.eps)

![Contour plots of the ($T_{kin}$, $n$\[H$_{2}$\]) 2D marginal likelihood distributions, generated by a MCMC sampling of the parameter space for RADEX models. The contours are in n$\sigma$-equivalent steps, enclosing 68.3%, 95.4%, 99.7%, and 99.99% of the probability, respectively. The dashed lines correspond to parameters that reproduce the LTE solution (model m01), and the dotted lines indicate the parameters corresponding to the RADEX 4D maximum likelihood solution (model m02). Note that the 2D marginal distributions will not necessarily have the same maximum as the 4D distribution. The kinetic temperature is limited to 3000 K, where collisional dissociation of CO becomes important. In the SDP.81 panel, the lighter contours show the probability levels for a model including the CO(1-0) from @Frayer10. The parameters for this model are listed as model m03 in Table \[tab:radex\]. \[tnmcmc\]](f10.eps)

![Same as Figure \[tnmcmc\] for the ($N[CO]$, $n$\[H$_{2}$\]) 2D marginal likelihood distributions. \[nnmcmc\]](f11.eps)

![$W(J)$ as a function of transition for four of the galaxies in our sample. For clarity, the values for the same transition in different galaxies have been slightly offset around the position of the upper $J$ level. The triangle point represents the intensity of the CO(1-0) line for SDP.81 measured by @Frayer10. The $W(J)$ distributions predicted by the RADEX models m01 and m02 are shown with a dashed and continuous line, respectively. These lines emphasize the constraints on the allowed parameter space that can be gained by having measurements of both higher and lower $J$ transitions. \[tempsone\]](f12.eps)
