High Resolution Imaging of PHIBSS $z\sim 2$ Main Sequence Galaxies in CO $J=1\rightarrow 0$

MD94 is clearly resolved in the observations (c.f., Table 4). The simple Gaussian fit (CASA imfit) yields an elongated source in the N-S direction, with deconvolved FWHM diameter (after averaging the C-only and C+D robust results) of $\sim 1\farcs 3\times 1\farcs 0$ at high significance ($\geq 5\sigma$). With an image scale of $\sim 8.4$ kpc per arcsecond at this redshift, this corresponds to $\sim 11\times 8$ kpc. The residual RMS in the central region after removing the Gaussian fit is very similar to that measured in empty regions of the original image, showing that a 2D Gaussian is a good model for the source at this signal-to-noise. BX610 is clearly resolved along its major axis, with a deconvolved CO FWHM diameter of $1\farcs 6\times 0\farcs 6$ corresponding to $\sim 13\times 5$ kpc (for this source we adopt the results of the C+D robust combination, which despite the larger beam is better behaved than the C-only observation). We also did UV-plane model fitting using the uvmodelfit task in CASA, and the results are essentially identical to those from image-plane fitting.

In $z=0$ galaxies there is a good relation between the sizes of molecular and stellar disks (e.g., Regan et al., 2001). The CO emission is found in an approximately exponential disk with a $R_{\rm CO}$ scale-length that is $R_{\rm CO}\simeq 0.2\,R_{25}$, the isophotal radius of the B-band light at 25^th magnitude (Young & Scoville, 1991; Young et al., 1995; Schruba et al., 2011). Consequently half of the CO emission in local disks is contained within a radius $\mbox{$R_{\rm 1/2,CO}$}=0.34\,R_{25}$. For optical light it is usually assumed that the stellar disk scale length is $R_{*}\sim 0.25\,R_{25}$ (implying an effective radius containing half the optical light of $\mbox{$R_{\rm 1/2,Opt}$}\sim 0.4R_{25}$), so the stellar and molecular disks track each other very well in their exponential part.

How do the CO and optical sizes relate in our high-$z$ SFGs? For a Gaussian source 50% of the flux is contained within the FWHM, therefore the measured $R_{\rm 1/2}$ CO radius for MD94 and BX610 is respectively $\mbox{$R_{\rm 1/2,CO}$}\sim 4.8$ and $4.1$ kpc. We can compare directly these numbers with near-IR NICMOS observations of BX610 (Figure 2) that sample the optical disk of that galaxy (Förster Schreiber et al., 2011). Rest-frame 5000 Å observations of BX610 yield an effective radius of $\mbox{$R_{\rm 1/2,Opt}$}\approx 4.6$ kpc, which within the errors is very comparable to the $J=1-0$ $R_{\rm 1/2,CO}$ we obtain. To put disk optical sizes in the $z\sim 0$ context, the typical half-light optical radius of a $M_{*}\sim 10^{11}$ M_⊙ in the SDSS sample analyzed by Shen et al. (2003) is $\mbox{$R_{\rm 1/2,Opt}$}\sim 4$ kpc. The existing HST observations for MD94 (Figure 2) are rest-frame 2500 Å (Peter et al., 2007), which combined with the AGN presence makes it difficult to obtain a reliable disk size (the 4.8 kpc optical radius reported by Tacconi et al., 2013, is from H$\alpha$). Surface photometry shows that 80% of the 2500 Å light is in an exponential component with $\mbox{$R_{\rm 1/2,Opt}$}=0.37\pm 0.08$ arcseconds ($3.1\pm 0.7$ kpc) using the Fisher & Drory (2008) methodology. Converting these measurements to an optical size is, however, not straightforward even in local galaxies (which display a range of scalings, Taylor et al., 2005; Muñoz-Mateos et al., 2007), much less in the poorly characterized high-$z$ population.

It is interesting to compare our results to some other resolved measurements of high-$z$ SMGs. Bothwell et al. (2010) analyze resolved observations of intermediate excitation CO $J=3-2$ and $J=4-3$ in several such galaxies at $z\sim 1.2$ to $z\sim 2$. Their typical CO 50% emission radii (computed as the harmonic mean of the major and minor axis) are $\mbox{$R_{\rm 1/2,CO}$}\sim 3$ kpc, although the observations span a range of $\mbox{$R_{\rm 1/2,CO}$}\sim 1-5$ kpc. Similar measurements are reported by Tacconi et al. (2008) and Engel et al. (2010), who find typical CO sizes $\mbox{$R_{\rm 1/2,CO}$}\lesssim 4$ kpc. Sharon et al. (2015) find lensing-corrected source sizes of $\mbox{$R_{\rm 1/2,CO}$}\sim 1.7$ kpc in $J=3-2$ and $\mbox{$R_{\rm 1/2,CO}$}\sim 2.0$ kpc in $J=1-0$ for the two components of the $z\sim 2.7$ SMG they study. The $\mbox{$R_{\rm 1/2,CO}$}\sim 4.1$ and 4.8 kpc sizes measured here suggest that the CO $J=1-0$ disks of main sequence galaxies tend to be somewhat more extended than the region producing mid-J CO emission in SMGs at comparable redshift. On the other hand, resolved CO $J=1-0$ measurements of five $z\sim 2$ SMGs by Ivison et al. (2011) find complex morphologies with radii of $\gtrsim 6$ kpc for four of the objects (the fifth is unresolved) and a median FWHM linewidth of 540 km s^-1. These SMG CO $J=1-0$ sizes and linewidths are larger than those measured for the same objects in the mid-J transitions, and also larger than the disk sizes we measure here.

The similarity between the molecular and optical sizes in main sequence galaxies is also seen in the CO $J=3-2$ measurements reported by Tacconi et al. (2013), mostly obtained for $z\sim 1$ sources (Fig. 3). These measurements suggest that $z\sim 2.3$ main sequence galaxy disks follow a scaling between optical light and CO $J=1-0$ similar to that in local disks. The CO is not measurably more extended than the stars, as may occur if a significant amount of pre-enriched newly-arrived material is present in the disk outskirts. Such pre-enriched material is unlikely to arrive through smooth accretion; it seems much more likely that occurs through mergers, and in particular minor mergers given the placement of these galaxies on the main sequence. The fact that the scaling between CO-emitting molecular gas and the stars is similar to that observed in local disks suggests that we are not catching these $z\sim 2$ main-sequence disks at an out-of-equilibrium stage. Rather, we observe them in an equilibrium situation where the established stellar population and the enriched CO-bright molecular gas that fuels the current star formation have similar radial distributions.

The excellent signal-to-noise of the observations for MD94 and BX610 permits, in principle, a precise determination of their CO $J=1-0$ luminosity. Unfortunately there is a distressing level of inconsistency between the fluxes measured in the C and D configurations, particularly for MD94, with the former being higher than the latter. After detailed inspection of the data and calibrations it is not clear what to attribute these differences to. In particular, the changes in the flux measured for the sources between D and C configurations are not mirrored by changes in the flux of their respective gain calibrators which could suggest systematic problems with the flux calibration (for MD94 the change is in the opposite direction), and are not due to the phase noise as measured on the gain calibrator (which shows RMS$\sim 10-20$ degrees in either configuration). After application of the calibrations we also recover the expected flux for the secondary calibrator used in D-array, as noted in §2. The flux inconsistency in MD94 is at the $\sim 2\sigma$ level, which suggests that it could be simply due to noise, much higher for the C-array data. The fact that the C-array observations were obtained during a transition period for the instrument may also have had an effect on the accuracy of those fluxes. Because of these considerations we base our discussion on the fluxes from the naturally weighted combination of C+D observations, which gives preferential weighting to the D-array data and has the best signal-to-noise.

To compute the molecular masses we use Eq. 3 in Bolatto et al. (2013), which is applicable to CO $J=1-0$, with a Galactic CO-to-H₂ conversion factor. This conversion relies fundamentaly on two assumptions: 1) the molecular gas is primarily self-gravitating, and 2) the ratio $\sqrt{\rho}/T$ between gas density, $\rho$, and temperature, $T$, is roughly similar to that in local galaxy disks (Bolatto et al., 2013, and references therein). We use $\mbox{$\alpha_{\rm CO}$}=4.36$ M_⊙ (K km s^-1 pc²)^-1 which includes the Helium contribution to the mass as is customary ($\mbox{$\alpha_{\rm CO}$}=3.2$ M_⊙ (K km s^-1 pc²)^-1 without the Helium contribution). A number of arguments point to an approximately Galactic value of the conversion factor for main sequence high-$z$ galaxies (Tacconi et al., 2010; Daddi et al., 2010a; Genzel et al., 2012), including the similarity between the gas masses inferred from CO and from dust (Genzel et al., 2015). This is supported by detailed numerical modeling of high-$z$ galaxies, which suggests that the average $\alpha_{\rm CO}$ in these systems is not significantly different from $\alpha_{\rm CO}$ in local disks (Bournaud et al., 2014a). We will see in §3.3 that there are indications that the molecular gas temperature is higher in these sources than the average in GMCs in the Galaxy. It is important to remark that this does not automatically translate into a lower-than-Galactic $\alpha_{\rm CO}$. In a self-gravitating cloud what matters is the $\sqrt{\rho}/T$ ratio, not just the temperature, as stated in assumption number two above. Indeed, there are indications that densities (certainly column densities) are much higher in these objects than in Galactic GMCs.

Under these assumptions, the molecular masses inferred from CO $J=1-0$ for MD94 and BX610 are $M_{\rm mol}\approx 1.7\times 10^{11}$ and $1.1\times 10^{11}$ M_⊙ respectively. These masses are similar to the stellar masses determined for these galaxies (Tacconi et al., 2013), resulting in molecular gas fractions of $f_{gas}\sim 50\%$. Bolatto et al. (2013) discuss a surface density correction that is applicable if the CO-emitting gas is not self-gravitating but bound to the overall potential of the galaxy as is applicable in a merger and some galaxy centers (c.f., their Eq. 31). For our galaxies, application of this correction would amount to reducing their molecular mass by a factor of $\sim 3$. This is consistent with what detailed modeling recovers for very turbulent starbursting mergers, but it is likely too large a correction for the galaxies presented here (by contrast modeling of very actively star-forming clumpy main sequence galaxies finds $\alpha_{\rm CO}$ within 20% of the Galactic value Bournaud et al., 2014a).

Given the sizes we measure for the molecular emission, the molecular masses inferred for MD94 and BX610 imply corresponding average deprojected surface densities of $\Sigma_{\rm mol}\sim 900$ M_⊙ pc^-2 ($i\sim 40^{\circ}$) and $400$ M_⊙ pc^-2 ($i\sim 65^{\circ}$) assuming the disk size is that of the deconvolved major axis and the measured elongation is due to projection effects ($\mbox{$\Sigma_{\rm mol}$}\sim 1200-1100$ M_⊙ pc^-2 in projection). The surface density found for BX610 is approximately the typical $\mbox{$\Sigma_{\rm mol}$}\sim 400$ M_⊙ pc^-2 inferred by Freundlich et al. (2013) using position-velocity information in four $z\sim 1.2$ galaxies, as well as spatially-resolved measurements in $z\sim 1.5$ galaxies (Tacconi et al., 2010; Genzel et al., 2013). MD94, on the other hand, exhibits higher deprojected disk surface densities than BX610 mostly as a consequence of its smaller inferred inclination, that place it at the upper end of the Freundlich et al. measurements. We note that these “deprojected” estimates are very tentative and, although more physically interesting than projected quantities, also considerably more uncertain because of the poor knowledge of the geometry. In particular, using the H$\alpha$-derived inclination of BX610 (see next paragraph) would result in a deprojected $\Sigma_{\rm mol}\sim 900$  M_⊙ pc^-2.

We can also extract size and velocity information from the data cubes themselves, and obtain dynamical mass estimates. We use GalPak3D, a Bayesian multiparameter Markov Chain Monte Carlo fitter for 3D galaxy data that takes into account the effects of instrumental resolution (Bouché et al., 2015), to fit the robust C+D cubes for MD94 and BX610. The process requires a starting set of assumptions (source center, size, inclination, velocity dispersion, etc), that are then used to produce a 3D model that is compared to the data to compute a reduced $\chi^{2}$ figure-of-merit which is then minimized. All resulting estimates are extremely tentative given the resolution and signal-to-noise of the data (Figure 4). The reduced $\chi^{2}$ resulting from the minimization are 1.15 and 1.05 for MD94 and BX610 respectively, showing that the model is a good fit to the data. The minimum, however, is shallow, showing that the models are not unique. BX610 is fit by a (thick) rotating disk with a maximum velocity of $v_{rot}\sim 220$ km s^-1 and intrinsic dispersion $v_{disp}\sim 40$ km s^-1, inclination of $i\sim 50^{\circ}$, and radius $\mbox{$R_{\rm 1/2,CO}$}\sim 0.4\arcsec$. These parameters are comparable to those obtained from the analysis of its H$\alpha$ kinematics (Cresci et al., 2009, Table 2), although both the rotational velocity and radius are lower (the latter is also lower than the radius we obtain from the integrated intensity map), and the inclination is higher. The resulting dynamical mass computed as $M_{dyn}=2\,\mbox{$R_{\rm 1/2,CO}$}v_{rot}^{2}/G$ using the more robust $\mbox{$R_{\rm 1/2,CO}$}=4.1$ kpc estimate from Table 4 is $M_{dyn}\sim 1.1\times 10^{11}$ M_⊙, lower than found by Cresci et al. (2009) and too low compared to the inferred baryonic mass $M_{*}+M_{\rm mol}\simeq 2.1\times 10^{11}$ M_⊙. This is, however, entirely attributable to the inclination found in the CO fit. If we instead use the H$\alpha$- derived inclination ($i\sim 33^{\circ}$) we correct $v_{rot}$ to $\sim 310$ km s^-1, resulting in $M_{dyn}\sim 2.0\times 10^{11}$ M_⊙ in excellent agreement with the estimate of the baryonic mass (we do not expect dark matter to make a significant contribution inside $R_{1/2}$). MD94, on the other hand, is fit as a dispersion-dominated system ($v_{disp}\sim 180$ km s^-1 and $v_{rot}\sim 60$ km s^-1 with $i\sim 45^{\circ}$). The corresponding dynamical mass is $M_{dyn}\sim 0.9\times 10^{11}$ M_⊙, too low to accommodate even the stellar mass. Unfortunately there are no H$\alpha$ kinematics to compare with, and the example of BX610 highlights the uncertainties of the analysis.

An interesting local comparison may be provided by the starbursting region of the nearby galaxy NGC253, which has 50% of its CO emission within an area of $150\times 50$ pc in projected size with a luminosity L${}_{\rm CO}\approx 3.3\times 10^{8}$ K km s^-1 pc² (Leroy et al., 2014, Table 5). That region of NGC253 has a star formation rate SFR$\sim 2-3$ M_⊙ yr^-1, while MD94 and BX610 have inferred integrated SFR$\sim 271$ and $212$ M_⊙ yr^-1 respectively (Tacconi et al., 2013). About $\sim 100$ NGC253-like complexes filling $\sim 2\%$ of the area of these main sequence disks would be compatible with our CO luminosity and surface brightness, while simultaneously fueling the observed star formation activity. This conceptual “picture” of an agglomeration of starburst-like regions also fits the $r_{31}$ excitation as we discuss in the next section. Note, however, that there is no deep physical significance to the inferred beam filling area other than the fact that it reproduces the observed integrated surface brightness.

The Cosmic Microwave Background (CMB) has a temperature of $T_{CMB}\approx 9$ K at $z\approx 2.3$, which sets an effective floor for the physical temperature of the cold molecular gas (see e.g., da Cunha et al., 2013). In terms of the Rayleigh-Jeans brightness temperature (radiation temperature) measured by the interferometer, the CMB presents a background that lowers the contrast. Thus even for optically thick, thermalized emission the observed Rayleigh-Jeans brightness temperature is systematically lower than the excitation temperature, $T_{ex}$, and vanishes when $T_{ex}=T_{CMB}$ in the rest frame (see Eq. 6 in Bolatto et al., 2013, which is applicable to the rest frame). The fact that we see relatively bright CO $J=1-0$ and $J=3-2$ emission already says that $T_{ex}\gg 9$ K.

The effect of the CMB on the ratio of two transitions is small unless the excitation temperature of the emission is very low. The $S_{\rm CO}\Delta v$ integrated flux ratio between the $J=3-2$ and $J=1-0$ transitions is $10.6\pm 1.5$ for MD94 and $8.3\pm 1$ for BX610. Converting to ratios in brightness temperature units requires dividing by the ratio of upper level quantum numbers, $J^{2}$, resulting in $r_{31}\approx 1.17\pm 0.17$ and $r_{31}\approx 0.92\pm 0.11$ for MD94 and BX610 respectively. These data have much better signal-to-noise in CO $J=1-0$ than observations used by Aravena et al. (2014) to measure $r_{31}=0.58^{+0.21}_{-0.13}$ in BX610, and the $r_{31}$ we obtain is at the $1.5\sigma$ upper boundary of the previous measurement.

Under the assumption of coextensive emission, these integrated numbers suggest that the $T_{ex}$ of the $J=1-0$ and $J=3-2$ transitions are very similar, and both transitions arise from optically thick warm gas (that is, gas with excitation and kinetic temperature significantly higher than the $E_{32}/k\approx 16.6$ K energy necessary to collisionally jump from $J=2$ to the $J=3$ level, where $k$ is Boltzmann’s constant).

The data strongly suggest that the emission in CO $J=1-0$ and $J=3-2$ is coextensive in these objects. The source sizes are very similar in both CO transitions in BX610, where the measured sizes of the CO emission for the $J=1-0$ and $J=3-2$ transitions are the same within the uncertainties ($R_{1/2}=4.1\pm 1.1$ kpc for $J=1-0$, $R_{1/2}=3.2\pm 0.6$ kpc for $J=3-2$). Moreover, the spectral shape of the emission is very similar for both transitions in MD94 and BX610, within the uncertainties associated with the noise (Figure 5). The respective linewidths in CO $J=1-0$(CO $J=3-2$) from Single Gaussian fitting are $296\pm 66$ ($430\pm 61$) km s^-1 and $294\pm 49$ ($266\pm 35$) km s^-1 respectively (the difference in the case of MD94 is driven by the blue side of the spectrum, where the fit to $J=3-2$ incorporates as line a wider range of velocities resulting in a shift of the centroid). Averaging these results yields a typical $\Delta v(\mbox{$J=1-0$})/\Delta v(\mbox{$J=3-2$})=0.90\pm 0.15$ on average for these two sources. This suggests that, within our ability to quantify them with the present data, no dramatic excitation gradients are present, as both transitions sample the same gas kinematics. This is in contrast with what may be the situation in SMGs, where Ivison et al. (2011) find a measurably larger typical linewidth for the CO $J=1-0$ in comparison to $J=3-2$, suggesting the colder gas is more extended. Note, however, that this is a subtle $\sim 15\%$ effect in the SMG sample, and so it could go undetected in our data.

These results point to a molecular gas excitation that is high and uniform in MD94 and BX610, up to the $J=3$ level. Given the widespread molecular emission, the large average gas surface densities, and the inferred SFRs, this should not be surprising. The excitation requirements for the $J=3-2$ transition are not particularly stringent: $T\gtrsim E_{32}/k\approx 16.6$ K, and an effective critical density $n_{crit}\approx 2\times 10^{4}/\tau_{\rm CO\,\,3-2}\sim{\rm few}\times 10^{3}$ cm^-3, making the $r_{31}$ ratio a rather blunt indicator of excitation. By comparison, the average relation for dust temperature along the main sequence predicts $T_{dust}\sim 31$ K for these objects (Genzel et al., 2015, Table 5), which can be taken as a proxy for the gas temperature in dense environments and certainly meets the requirement for excitation to $J=3$. The measured $r_{31}$ values show that there is not a large reservoir of cold and/or diffuse low excitation molecular gas, bright in CO $J=1-0$ but not apparent in $J=3-2$, present in these sources.

In a study using a library of galaxy simulations, Narayanan & Krumholz (2014) make the case for the link between $\Sigma_{\rm SFR}$ and CO line ratios in galaxies. They show that in moderately active star-forming galaxies ($\mbox{$\Sigma_{\rm SFR}$}\gtrsim 0.2$ M_⊙ yr^-1 kpc^-2) the molecular gas has average temperatures larger than the energy of the $J=3$ level, and given its large typical optical depth it is very easy to excite this transition through collisions and radiative trapping. In particular, they parametrize the expected line ratios as a function of $\Sigma_{\rm SFR}$. Using their Eq. 19 and the $\mbox{$\Sigma_{\rm SFR}$}\gtrsim 2-3$ M_⊙ yr^-1 kpc^-2 measured for our galaxies, we would expect $r_{31}\sim 0.9$, very similar to what we measure.

Interestingly, our results are also consistent with the $r_{31}$ measured for CO in the central regions of NGC253 (Bradford et al., 2003). Another point of comparison is provided by the nearby Seyfert NGC1068. There the line ratio averages to $r_{31}\sim 1.2$ in the starburst ring, but with a large range of values correlated with the local star formation activity measured in Paschen $\alpha$: $r_{31}\sim 0.7-1$ in regions with low SFR, and $r_{31}\sim 2-3$ in regions with high SFR (Garcia-Burillo et al., 2014). The compact circumnuclear disk region, heavily influenced by the AGN, shows $r_{31}\sim 2.7$.

Actively star-forming $z=0$ galaxies span a range of $r_{31}$, and although a typical value is $r_{31}\sim 0.66$ there exist several local examples with $r_{31}\geq 0.9$ (Figure 6; Mauersberger et al., 1999; Yao et al., 2003; Papadopoulos et al., 2011). There has been a tendency to attribute the high $r_{31}$ measured in some high-$z$ sources to the presence of a strong AGN (e.g., Riechers et al., 2011). Nonetheless, it is clear now that at least some high-$z$ sources display high $r_{31}$ ratios with no clear AGN activity (Sharon et al., 2013, 2015). That is certainly the case also in the NGC253 starburst, where no strong AGN is detected, and in the starburst ring of NGC1068 where the influence of the AGN is negligible. It is more likely that $r_{31}$ is not so much acting as an indicator of overall density and temperature, but rather showing how well mixed the star formation activity (which provides the bulk of the gas heating) and the molecular reservoir are. Galaxies where part of the molecular reservoir is quiescent (thus cold and/or diffuse) will show lower $r_{31}$, while galaxies where most of the molecular gas in experiencing star formation will exhibit higher $r_{31}$ ratios.