© ESO, 2015

1. Introduction

Massive star clusters, those containing > 10⁴M_⊙ of stars, are among the most visually outstanding features in the night sky (see review by Longmore et al. 2014). In other galaxies, they are useful probes of the star formation history and can be individually identified and measured (Bastian 2008). Locally, they are the essential laboratories in which we can study the formation of massive stars (Davies 2012).

In order to utilize these clusters as laboratories, we need to understand their formation in detail. Clusters are often assumed and measured to be coeval to within a narrow range (e.g. 10⁵ yr; Kudryavtseva et al. 2012), but uncertainties remain (Beccari et al. 2010). In the most massive clusters, there are predictions that multiple generations or an extended generation of stars should form prior to gas expulsion because the gas will remain gravitationally bound (Bressert et al. 2012). Feedback from and within young massive clusters is an active field of numerical study (Rogers & Pittard 2013; Dale et al. 2005, 2013, 2012; Dale & Bonnell 2008; Parker & Dale 2013; Myers et al. 2014; Krumholz et al. 2014). While only 5−35% of all stars form in bound clusters¹ (Kruijssen 2012), these clusters form the basis of our understanding of stars and stellar evolution (Kalirai & Richer 2010), and understanding their formation is therefore crucial.

The results of cluster formation may be decided before the first stars are formed. The starless initial conditions of massive clusters have not yet been definitively observed (Ginsburg et al. 2012) though there are viable candidates such as G0.253+0.016 (Longmore et al. 2012). The initial conditions for star formation on any scale are clearly turbulent. However, there is no evidence whether these initial conditions differ in any qualitative way from turbulence in local, low-mass star-forming regions.

Star formation appears to occur efficiently only above a volume density threshold in molecular gas, specifically n(H₂) ≳ 10⁴cm^-3, in the Galactic disk (Lada et al. 2010, who advocate a column density threshold corresponding to this density); note however that the threshold is smooth, not a step function as implied by the term (Padoan et al. 2014). The dense gas mass fraction has been measured for various definitions of dense (Battisti & Heyer 2014; Wu et al. 2005), but these definitions don’t always correspond to the associated threshold (Kauffmann & Pillai 2010; Parmentier et al. 2011; Parmentier 2011). Often, the dense gas threshold is observationally defined as the threshold to see a given molecule, which is sometimes incorrectly assumed to correspond to a fixed critical density for the molecule. The existence of a universal density threshold for star formation is contentious (Burkert & Hartmann 2013; Clark & Glover 2013); the star formation threshold likely varies with the turbulent properties of clouds (Longmore et al. 2013; Kruijssen 2014; Hennebelle & Chabrier 2011, 2013; Padoan & Nordlund 2011; Krumholz et al. 2005). Star formation thresholds can only be evaluated when measurements of density and column density are simultaneously available.

1.1. W51

The W51 cloud complex (Fig. 1), containing the two massive protocluster candidates W51 Main and W51 IRS 2 from Ginsburg et al. (2012), is located at ℓ ~ 49,b ~ − 0.3, very near the Galactic midplane² at a distance of 5.1 kpc (Sato et al. 2010). It is a well-known and thoroughly studied collection of clouds massing M> 10⁶M_⊙ (Carpenter & Sanders 1998; Bieging et al. 2010; Kang et al. 2010; Parsons et al. 2012). The radio-bright regions are generally known as W51 A to the east, W51 B to the west, and W51 C for the southern component, known to trace a supernova remnant (Koo et al. 1995; Brogan et al. 2000, 2013).

1.2. Formaldehyde

Formaldehyde (H₂CO) has been recognized as a useful probe of physical conditions in the molecular interstellar medium for decades (Mangum & Wootten 1993). The centimeter lines, H₂CO1₁₀ − 1₁₁ (6.2 cm, 4.82966 GHz) and 2₁₁ − 2₁₂ (2.1 cm, 14.48848 GHz), have a peculiar excitation process in which collisions overpopulate the lower of the two K_c rotational states, where K_c is the quantum number representing the angular momentum projected onto the long axis of the molecule. The overpopulated lower-energy state leads to amplified absorption, or effective excitation temperatures less than the background temperature, allowing H₂CO centimeter absorption to be seen against the cosmic microwave background (CMB). Because the 1₁₀ − 1₁₁ and 2₁₁ − 2₁₂ level pairs populate differently depending on the volume density of the colliding partner (a mix of p-H₂, He, and o-H₂), their ratio is sensitive to the local gas volume density³.

H₂CO1₁₀ − 1₁₁ has been observed in the W51 Main region with the VLA (Martin-Pintado et al. 1985a) and Westerbork (Arnal & Goss 1985), and this data was used to gain some early constraints on the geometry of the region (e.g. Carpenter & Sanders 1998). Henkel et al. (1980) presented observations of the 1₁₀ − 1₁₁ and 2₁₁ − 2₁₂ lines, and Martin-Pintado et al. (1985b) presented single-dish mapping observations of the H₂CO2₁₁ − 2₁₂ line toward the W51 Main region, but both treated the region as a single-density structure.

1.3. Paper overview

We present a detailed examination of the dense gas in the W51 cloud complex. Section 2 presents the observations and data reduction. Section 3 describes the analysis techniques, including measurement of the dense gas mass fraction and derivation of the cloud geometry. Section 4 discusses the implications of the density measurements for Galactic clouds and extragalactic interstellar media. Section 5 concludes. There are two appendices: Appendix A describes the radio recombination line and continuum data. Appendix B describes details of individual regions, with a focus on the cloud and H ii region line-of-sight geometry.

2. Observations and data reduction

The W51 survey was performed in September 2011 on the Green Bank Telescope (GBT) and in 2012 using the Arecibo Observatory. The data was reduced using custom-made scripts based off of both GBTIDL’s mapping routines by Glen Langston⁴ and Phil Perillat’s AOIDL routines. The code is available at https://github.com/keflavich/sdpy. The data reduction code and workflow are included in a corresponding git repository⁵.

2.1. Arecibo 6 cm

The Arecibo data were taken as part of project A2705 over the course of 4 nights, September 10, 11, 12, and 15 2012. The Mock spectrometer was used to cover the range 4.6 to 5.4 GHz with a spectral resolution ~1 km s^-1, including the o - H₂CO 4.82966 GHz line and the H107-H112α recombination lines. On the first night, September 10 2012, a significant fraction of the data was lost due to an internal instrument error within the Mock spectrometer, which resulted in a loss of the high spectral resolution component of the H₂CO data for that night. As a result, we have focused our study on the lower-resolution (~ 1 km s^-1) data.

The fields were observed with east-west maps using the C-band receiver. No crosshatching was performed with Arecibo.

The Arecibo data reduction process for W51 presented unique challenges: at C-band, the entire region surveyed contains continuum emission, so no truly suitable off position was found within the survey data. Similarly, H₂CO is ubiquitous across the region, so it was necessary to mask out the absorption lines when building an off position. This was done by interpolating across the line-containing region of the spectrum with a polynomial fit. The fits were inspected interactively and tuned to avoid over-predicting the background.

The Arecibo data were corrected to main beam brightness temperature T_MB using a main-beam efficiency as a function of zenith angle in degrees (za):

${\mathit{\eta }}_{\mathrm{MB}}\mathrm{\left(}\mathit{za}\mathrm{\right)}\mathrm{=}\mathrm{0.491544}\mathrm{+}\mathrm{0.00580397}\mathit{za}\mathrm{-}\mathrm{0.000341992}\mathit{z}{\mathit{a}}^{\mathrm{2}}\mathit{.}$This is a fit to 5 years worth of calibration data acquired at Arecibo and assembled by Phil Perrilat⁶.

The maps were made by computing an output grid in Galactic coordinates with 15′′ pixels and adding each spectrum to the appropriate pixel⁷. In order to avoid empty pixels and maximize the signal-to-noise (S/N), the spectra were added to the grid with a weight set from a Gaussian with FWHM = 20″, which effectively smooths the output images from FWHM ≈ 50″ to ≈54″. See Mangum et al. (2007) for more detail on the on-the-fly mapping technique used here.

The Arecibo data were taken at a spectral resolution of 0.68 km s^-1. The spectra were regridded onto a velocity grid from − 50 to 150 km s^-1 with 1 km s^-1 resolution. To achieve this, they were first Gaussian-smoothed to FWHM = 1 km s^-1 then downsampled appropriately.

The position-position-velocity (PPV) cubes were created with units of brightness temperature. The Arecibo cubes have contributions from 15−20 independent spectra in each pixel, though this hit rate varies in a systematic striped pattern parallel to the Galactic plane. The small overlap regions between different maps have a significantly higher number of samples; these regions constitute a small portion of the map. The resulting noise level is rms ~ 50 − 60 mK except toward the H ii regions, where it peaks at about 400 mK. The continuum is derived by averaging line-free channels; its S/N peaks at ~ 900.

The Arecibo data have smaller systematic continuum offsets than the GBT data (Fig. 2), but they are more visually pronounced because there is much more diffuse emission at 6 cm. The continuum zero-point level in the Arecibo data was set to be the 10th percentile value of each scan, which is effectively the minimum value across each scan but with added robustness against noise-generated false minima. In the eastmost and westmost blocks, this strategy was very effective, as there are clearly areas in each scan that see no continuum. However, in the central block, this approach resulted in a vertical negative filament that almost certainly represents a local minimum that should be positive. This negative filament has values ≳− 0.08 K. Given the excellent agreement between the three independently observed regions in the areas that they overlap, it is clear that the continuum is reliable above ≳0.5 K, which is the entire regime in which it is a significant contributor to the total background emission (at lower levels, the CMB is dominant).

2.2. GBT 2 cm

The GBT data were taken as part of program AGBT10B/019. We used the GBT Spectrometer with 4 windows covering the o - H₂CO 14.488479 GHz line, H77α (14.12861 GHz), and two others that were not reduced targeting H₂CN (14.82579 GHz) and NaCl and SO (13.03606 GHz); online examination of the latter windows suggested that we did not detect any emission or absorption. The data presented in this paper include sessions 10, 11, 14, 16, 17, 20, 21, and 22; the other sessions from this project include maps of outer galaxy regions and a single-pointing survey of Galactic plane sources that will be presented in another paper.

Data were taken in on-the-fly mode with the GBT Ku-band dual-beam system. Cross-hatched north-south and east-west maps were created in Galactic coordinates. Each spectrum was calibrated using the first and last scans of each observation as off positions. The background level to be subtracted off of the continuum was determined by linearly interpolating between these scans.

The Green Bank data have a main beam efficiency η_MB = 0.886, or a gain of 1.98 K/Jy assuming a 51′′ beam (see Mangum et al. 2013, for additional discussion). The GBT data were also corrected for atmospheric opacity using Ron Maddalena’s getForecastValues⁸ with a typical zenith optical depth τ_z ≈ 0.02; this correction was never more than ~5%. The GBT data were taken with a spectral resolution Δv = 0.25 km s^-1

Typical noise levels were ~10−20 mK per 1 km s^-1 channel; the levels vary across the map.

The GBT data were mapped in an orthogonal grid pattern, so the hit coverage is more uniform on small scales than the Arecibo data, but because of the dual-beam Ku-band system, the overall noise levels are much more patchy. Additionally, the nights with better weather yielded a lower noise level. The noise ranges from ~7 mK in the W51 Main region to ~20 mK in the westmost region. As with the Arecibo data, the H ii region adds noise, but the peak noise towards an H ii region is only ~20 mK. This difference is because the diffuse H ii region is fainter at 2 cm. The S/N in the continuum peaks at ~400.

While the noise in the continuum is nominally quite low, there are significant systematic effects visible in the continuum maps. The continuum zero-point of each GBT map was determined by assuming that the first and last scan both observed zero continuum and that the sky background can be linearly interpolated between the start and end of the observations. In general, these are good assumptions, but they leave in systematic offsets of up to ≲− 0.15 K in the maps, most likely because there is a ~ 0.15 K variation in the diffuse background emission.

2.3. Optical depth data cubes

The data cubes were converted into optical depth data cubes by dividing the integrated H₂CO absorption signature by the measured continuum level. We added a fixed background of 2.73 K to the reduced images to account for the CMB, which is absent from the images due to the background-subtraction process. We define an observer’s optical depth ${\mathit{\tau }}_{\mathrm{obs}}\mathrm{=}\mathrm{-}\ln \left[\frac{{\mathit{T}}_{\mathrm{MB}}}{{\mathit{T}}_{\mathrm{bg}}}\right]$(1)as opposed to the true optical depth, which is modeled in radiative transfer calculations $\mathit{\tau }\mathrm{=}\mathrm{-}\ln \left[\frac{{\mathit{T}}_{\mathrm{MB}}\mathrm{-}{\mathit{T}}_{\mathrm{ex}}}{{\mathit{T}}_{\mathrm{bg}}\mathrm{-}{\mathit{T}}_{\mathrm{ex}}}\right]\mathrm{·}$(2)The approximation τ_obs = τ is valid for T_ex ≪ T_bg, which is true when an H ii region is the backlight but generally not when the CMB is. Displaying the data on this scale makes regions of similar gas surface density appear the same, rather than being enhanced where there are backlights. The noise is correspondingly suppressed where backlighting sources are present.

H₂CO absorption is ubiquitous across the map. In the Arecibo data, 8044 of 17 800 spatial pixels have peak optical depths >5σ, and 14 547 have peaks >3σ, so H₂CO absorption is detected at ~ 80% of the observed positions.

The GBT H₂CO2₁₁ − 2₁₂ data have lower peak signal-to-noise because the continuum background is lower. Additionally, the 2₁₁ − 2₁₂ line is expected to trace denser gas and therefore be detected along fewer lines of sight. The 2₁₁ − 2₁₂ line is detected with a peak at >5σ in 3497 pixels (20%) and > 3σ in 12 254 pixels (69%). The high detection rate validates H₂CO as an efficient dense-gas tracer.

There were no detections of H₂CO1₁₀ − 1₁₁ or 2₁₁ − 2₁₂ emission. The significance of the nondetection of emission is discussed in Sect. 4.3.

2.4. A note on nondetections

H₂¹³CO was not detected anywhere in the W51 complex in either the 1₁₀ − 1₁₁ or 2₁₁ − 2₁₂ lines. The peak signal-to-noise in the H₂CO1₁₀ − 1₁₁ cube was 180 at 1 km s^-1 resolution (corresponding to an optical depth τ_obs ~ 1/5 at the peak continuum detection point), so we report a 3σ upper limit on the H₂CO/H₂¹³CO ratio R> 60, which is consistent with solar values of the ¹²C/¹³C ratio.

3. Analysis

3.1. H₂CO modeling

The H₂CO line ratio can be transformed into a volume density of hydrogen n(H₂) using large velocity gradient model grids. The column density of o - H₂CO, the ortho-to-para-ratio of H₂, and the gas temperature are the three main nuisance parameters that can be marginalized over.

The o - H₂CO column density is degenerate with the velocity gradient in LVG models. The H₂ column density is degenerate with this gradient and the abundance of o - H₂CO. Precise measurements of the H₂CO abundance are not generally available, but typical values of X_{o - H2CO} = 10^-10 − 10^-8 relative to H₂ are generally assumed (Mangum & Wootten 1993; Ginsburg et al. 2011; this paper; Ao et al. 2013) and found to be consistent with the observations. Nonetheless, little is known about local variations in o - H₂CO abundance, except that it freezes out in cold, dense cores (Young et al. 2004). The abundance was left as a free parameter in the model fitting.

The model grids were generated using RADEX LVG models (python wrapper⁹; original code van der Tak et al. 2007) and the grids were fit using https://github.com/keflavich/h2co_modeling. The RADEX models assume a velocity gradient of 1 km s^-1pc^-1. Ginsburg et al. (2011, this paper) discussed the effect of a local gas density distribution on the molecular excitation, but due to the complexity involved in accounting for these effects, we ignore them here. The derived physical parameters are moderately sensitive to the input collision rates, with a 50% error in collision rates yielding a factor of 2 error in derived density (Zeiger & Darling 2010), but the error on the collision rates in the low temperature regime we are modeling have recently been improved from the previously used Green (1991) rates and should not be a dominant factor in our calculations (Troscompt et al. 2009; Wiesenfeld & Faure 2013).

3.2. H₂CO observables

Figure 3 shows the most important observed properties of the H₂CO lines. The figures show the peak observed optical depth τ_obs = −ln(T_MB/T_continuum) in each line along with the ratio of the 1₁₀ − 1₁₁ to the 2₁₁ − 2₁₂ optical depth. The noise is computed by measuring the standard deviation over a signal free region (− 50 to 0 km s^-1) along each spatial pixel. The cubes were masked to show significant pixels determined by:

• 1.

  Selecting all voxels withS/N> 2 in both images or S/N> 4 in either and with at least 7 (of 26 possible) neighbors also having S/N> 2

• 2.

  Selecting all voxels with ≥10 neighbors having S/N> 2

• Growing (dilating) the included mask region by 1 pixel in all directions

• 4.

  Selecting all voxels with ≥5 neighbors marked as included by the above steps (this is a closing operation)

• 5.

  (2D only) When used to mask 2D images, the selection mask is then collapsed along the spectral axis such that any pixel containing at least one voxel along the spectral axis is included.

This approach effectively includes all significant pixels and all reliably detected regions within the data cube, though the number of neighbors used at each step and the selected growth size are somewhat arbitrary and could be modified with little effect.

Figures 3−5 each contain peak optical depth maps and two ratio maps. The first ratio map shows the observed optical depth ratio, while the second shows the true optical depth ratio assuming an excitation temperature for each line, T_ex(1₁₀ − 1₁₁) = 1.0 K and T_ex(2₁₁ − 2₁₂) = 1.5 K. These excitation temperatures are representative of those expected for most of the modeled parameter space in which absorption is expected. Fitting of individual lines-of-sight confirm that good fits can be achieved using these assumed temperatures.

However, there are some clear outliers within the map: the clouds at G48.9-0.3 and G49.4-0.2 both show very low 1₁₀ − 1₁₁/2₁₁ − 2₁₂ ratios over a broad area. As discussed in Sects. B.2 and B.8, these two regions have H ii regions in the foreground of the molecular gas. The ratios displayed in Fig. 3 are therefore computed with an incorrect background assumed; we correct for the different background in the next section.

3.3. Density maps

We computed the density in each voxel using the χ² minimization technique from Ginsburg et al. (2011). We measure χ² over the full 4D parameter space (density, column density, gas temperature, and ortho-to-para ratio [OPR]) $\begin{array}{ccc}{\mathit{\chi }}^{\mathrm{2}}& \mathrm{=}& {\left(\frac{{\mathit{T}}_{\mathrm{B}}\mathrm{\left(}{\mathrm{1}}_{\mathrm{10}}\mathrm{-}{\mathrm{1}}_{\mathrm{11}}\mathrm{\right)}\mathrm{-}{\mathit{T}}_{\mathrm{model}}\mathrm{\left(}{\mathrm{1}}_{\mathrm{10}}\mathrm{-}{\mathrm{1}}_{\mathrm{11}}\mathrm{\right)}}{\mathit{\sigma }\mathrm{\left(}{\mathit{T}}_{\mathrm{B}}\mathit{,}{\mathrm{1}}_{\mathrm{10}}\mathrm{-}{\mathrm{1}}_{\mathrm{11}}\mathrm{\right)}}\right)}^{\mathrm{2}}\\ & & \mathrm{+}{\left(\frac{{\mathit{T}}_{\mathrm{B}}\mathrm{\left(}{\mathrm{2}}_{\mathrm{11}}\mathrm{-}{\mathrm{2}}_{\mathrm{12}}\mathrm{\right)}\mathrm{-}{\mathit{T}}_{\mathrm{model}}\mathrm{\left(}{\mathrm{2}}_{\mathrm{11}}\mathrm{-}{\mathrm{2}}_{\mathrm{12}}\mathrm{\right)}}{\mathit{\sigma }\mathrm{\left(}{\mathit{T}}_{\mathrm{B}}\mathit{,}{\mathrm{2}}_{\mathrm{11}}\mathrm{-}{\mathrm{2}}_{\mathrm{12}}\mathrm{\right)}}\right)}^{\mathrm{2}}\mathrm{·}\end{array}$(3)The modeled brightness temperature T_model is different for each spatial pixel in order to account for the varying continuum background, although for some pixels with significant continuum detected, we still use the CMB as the background continuum because the molecular gas is behind the other continuum emission; see Sect. 3.4.

We have not enforced any constraints on the column density, temperature, or ortho-to-para ratio when fitting. The best-fit value of each of these parameters is taken to be the mean of those parameters over the range $\mathrm{\Delta }{\mathit{\chi }}^{\mathrm{2}}\mathrm{=}{\mathit{\chi }}^{\mathrm{2}}\mathrm{-}{\mathit{\chi }}_{\min }^{\mathrm{2}}\mathit{<}\mathrm{1}$, where ${\mathit{\chi }}_{\min }^{\mathrm{2}}$ is the minimum χ² value.

The gas temperatures returned from the fitting process are, as expected, purely noise: the H₂CO1₁₀ − 1₁₁/2₁₁ − 2₁₂ ratio provides virtually no constraint on the gas temperature and therefore leaving it as a free parameter has no effect on the fitted density. Similarly, the ortho-to-para ratio of H₂ is unconstrained in our data. In principle, the H₂ OPR has some effect on H₂CO excitation (Troscompt et al. 2009), but in the regime we have modeled and observed, no effect is apparent.

The o - H₂CO-column-weighted volume-density along each line of sight is shown in Figs. 6 and 7. The former shows the weighted density over all voxels and the latter shows the weighted density over the two velocity ranges previously discussed. These projections include no information about the errors in the individual fits, which are available from Eq. (3), but by weighting by column density, we have effectively selected the highest S/N points; the statistical errors are therefore negligible relative to the systematic (i.e., those caused by invalid assumptions about the single-zone nature of the models) in these maps. The maps are shown split into two velocity components, v< 62 km s^-1 and v> 62 km s^-1, which approximately separates out the 68 km s^-1 cloud from other components, though because the lines are quite broad the separation is imperfect. Figures B.3 and B.4 show the velocity separation in more detail.

The overall picture is of a central protocluster region (W51 Main and IRS2) with most of the gas mass at a density n ~ 10^5.5cm^-3 within a diameter of ~ 3 pc, surrounded by a rich cloud in which most of the mass is at a density ~ 10⁴cm^-3 out to a diameter d ~ 14 pc.

3.4. Model fitting and geometry

Both H₂CO lines are seen only in absorption. However, in some cases the absorption is against a continuum background, while in others the absorption may be only against the CMB.

We have fit the H₂CO lines constrained by the LVG models (Sect. 3.1) to spectra averaged over apertures with coherent molecular absorption signatures (e.g., Fig. 8). We compared the χ² values for fits with the observed continuum as the background to those with the background fixed to T_BG = T_CMB. We then selected the better of the two fits as representative of the physical conditions. Figures 9 and 10 show the qualitative version of this analysis, highlighting CO-bright regions that lack the expected H₂CO absorption signatures and therefore have a different geometry. Regions with the continuum emission in front of the molecular absorption were converted into masks that were then used to decide which models to use for the per-pixel fitting process described in the previous section. The geometric analysis of each region is discussed in detail in Appendix B.

It is possible that there are multiple continuum emitters along the line of sight in many cases, with the absorbing molecular gas somewhere in the middle. While this possibility adds uncertainty to the measurements, there are some cases in which the dominant continuum can unambiguously be assigned a foreground or background position.

Summary figures of our geometric analysis are shown in the cartoon Fig. 11, with an accompanying labeled on-sky map in Fig. 12.

3.5. Dense gas mass fractions

Because the H₂CO densitometer yields a mass-weighted¹⁰ measurement of the gas volume density, it is difficult to connect directly to the total gas mass, which is the quantity of interest when determining bulk properties like star forming efficiency. However, because the H₂CO and CO chemically trace the same gas, the H₂CO-derived density can be applied to the total mass measured by CO. We use the the FCRAO Galactic Ring Survey of ¹³CO (Jackson et al. 2006) and assume that each ¹³CO PPV voxel has a mass proportional to its integrated intensity and a density given by the n(H₂) delivered from the H₂CO densitometer.

The dense gas mass fraction (DGMF) is an oft-quoted measurement used to argue about the speed of the star formation process, the existence of density thresholds, and turbulent properties of the ISM (e.g. Fig. 5 of Krumholz et al. 2007; Battisti & Heyer 2014; Kainulainen et al. 2013; Juneau et al. 2009; Muraoka et al. 2009; Hopkins et al. 2013). However, these fractions are most often quoted as mass of gas at a single density divided by the total mass. We present an improvement on these measurements, showing the continuous distribution of the dense gas mass fraction.

Figure 13 shows the result of using our H₂CO PPV density cubes to measure the DGMF from ¹³CO. We use a range of density thresholds from ~ 10³ to ~ 10⁶cm^-3. At each density, we identify all voxels in the ¹³CO PPV cube above that density and integrate those. We then divide by the total integrated ¹³CO brightness to get the mass fraction: $\begin{array}{ccc}& & \mathit{N}\mathrm{(}\mathit{n}\mathit{>}{\mathit{n}}_{\mathrm{0}}\mathrm{)}\mathrm{=}\mathit{C}\mathrm{\int }{\mathit{T}}_{\mathrm{13}\mathrm{CO}}\mathrm{\left(}\mathit{v,\ell ,b}\mathrm{\right)}\mathit{G}\left(\mathit{n}\mathrm{\left(}\mathit{v,\ell ,b}\mathrm{\right)}\right)\mathrm{d}\mathit{v}\mathrm{d}\mathit{\ell }\mathrm{d}\mathit{b}\\ & & \mathit{G}\mathrm{\left(}\mathit{n}\mathrm{\right)}\mathrm{=}\{\begin{array}{c}\\ & \\ & \end{array}\\ & & \mathit{N}\mathrm{\left(}\mathrm{total}\mathrm{\right)}\mathrm{=}\mathit{C}\mathrm{\int }{\mathit{T}}_{\mathrm{13}\mathrm{CO}}\mathrm{\left(}\mathit{v,\ell ,b}\mathrm{\right)}\mathrm{d}\mathit{v}\mathrm{d}\mathit{\ell }\mathrm{d}\mathit{b}\\ & & \mathrm{DGMF}\mathrm{=}\mathit{N}\mathrm{(}\mathit{n}\mathit{>}{\mathit{n}}_{\mathrm{0}}\mathrm{)}\mathit{/}\mathit{N}\mathrm{\left(}\mathrm{total}\mathrm{\right)}\end{array}$where N is the column of H₂ and C = 8.1 × 10²⁰ cm^-2/ (K km s^-1) is a constant used to convert ¹³CO brightness to H₂ column; it cancels in the DGMF equation.

Figure 14 shows the same results, but for two individual regions: the W51 Main/IRS2 protoclusters and the W51 B region. Within about 10 pc of W51 Main, around half of the mass is at density n> 10⁴cm^-3. By contrast, the rest of the molecular cloud shows a consistent fraction f(n> 10⁴ cm^-3) ~ 10%.

These DGMF plots are similar to cumulative distribution functions of the column density (e.g. Battersby et al. 2014, Fig. 6), but with the added advantage of assigning a density to each resolution element in both velocity and position.

The multi-density DGMF presented here can be compared to models of global collapse in which progressively more gas should be observed in denser structures over time. They are effectively a gas density cumulative distribution function. However, to understand the systematic effects of line-of-sight stacking of different velocity components (and corresponding radiative transfer issues), similar analysis should be performed on hydrodynamic simulations.

3.5.1. Dense gas fraction assumptions and caveats

The DGMF analysis relies on the ¹³CO being optically thin and thermally excited, both of which are generally good assumptions for the majority of the mass. The molecular cloud probably includes no more than ~ 20% of its mass in ¹²CO-dark gas (Pineda & Teixeira 2013; Langer et al. 2014; Smith et al. 2014), which adds little to the overall uncertainty. Similarly, we expect that there is little ¹³CO-dark molecular gas (the ratio of ¹²CO to ¹³CO should not vary significantly in molecular clouds; Visser et al. 2009).

In Sect. 3.3, we discussed the various caveats and issues related to H₂CO density fitting. To account for the full range of errors in that analysis, we have plotted the DGMF calculated using the minimum and maximum values of the H₂CO-derived density consistent with the data at the 1σ level in Figs. 13 and 14.

We have assigned all of the mass associated with a given PPV voxel with a single, fixed density in this analysis. There is certainly some mass at a lower density in each PPV pixel associated with each voxel along that line of sight. This additional mass biases the measured DGMF upwards, but probably only by a small amount at each threshold. This systematic bias can be better characterized by performing a similar analysis on molecular cloud simulations projected into PPV space (as demonstrated for other analysis techniques by Beaumont et al. 2013).

4. Discussion

The W51 cloud complex includes a full range of star forming conditions. In the west, W51 B, there is an older generation of stars including at least one supernova remnant. In the east, there is a pair of forming, still-embedded massive clusters. We have described the geometry of these regions and features of the cloud structures, now we speculate on the broader implications of these observations.

The gas in the W51 B region, while clearly impacted by the expanding W51 C supernova, is less dense than most of the gas in the W51 A region. The supernova feedback is, if anything, destructive; a collect-and-collapse scenario does not fit the observed gas structure since there is less dense gas in the vicinity of the SNR. Given the large radius of the W51 C SNR, ≳100 pc, collection and the early stages of collapse should have happened by now if they are to happen at all.

The proximity of the 68 km s^-1 filamentary high velocity stream and the W51 Main protocluster and their relative line-of-sight positions have been presented as evidence for a cloud-cloud collision (Kang et al. 2010). Examination of the H₂CO line ratios has shown that the protocluster is embedded in the ~55 km s^-1 molecular cloud (e.g. Fig. 10). The velocity difference and their relative positions along the line of sight unambiguously indicate that these components are approaching each other, which is consistent with the cloud-cloud collision hypothesis, though the distinction of these velocity components as individual clouds is somewhat arbitrary since they are components of the same hierarchical medium.

It is possible that the 68 km s^-1 cloud is streaming in to a spiral arm, while the 50 km s^-1 cloud is slowed down as it exits the spiral arm on the far side. In this scenario, W51 Main is in the deepest part of the spiral arm potential. Gas is accumulating at W51 Main, becoming compressed and undergoing a mini-starburst. The W51 B/C region, with its mature HII regions and SNR, represents a slightly older generation than W51 Main.

The line-of-sight length of the W51 complex is still uncertain, despite our constraints on the relative geometry of different regions. The best prospect for resolving the line-of-sight structure of the region is via precise constraints on distances to the individual regions. Spectrophotometric surveys of the individual stellar sub-clusters may be able to provide this and should be undertaken. Maser parallax observations of different zones may also provide differential distance estimates.

4.1. Gas density and its relation to star formation

To assess the relation between gas density and star formation rate or efficiency, one can multiply the DGMF by the total mass to get a naïve measurement of the total mass directly involved in star formation (Fig. 15 transformed to Fig. 16). We have compared this map to the distribution of Class I and flat-spectrum Spitzer YSO candidates from Kang et al. (2009) in Fig. 16. While there are some regions of decent agreement between the YSO density and the star-forming gas mass, e.g. in the W51 B cloud and some of the more diffuse filaments, the densest pockets of star-forming gas contain few or no YSOs. Spitzer mid-infrared sources are either too confused or obscured to be detected in these regions. The Herschel images are also too crowded to identify a full sample of individual YSOs.

We also compare the star-forming gas surface density to a 21+24 μm-derived star formation rate surface density in Fig. 16 using standard extragalactic SFR calibrations (Price et al. 2001; Carey et al. 2009; Rieke et al. 2009; Kennicutt & Evans 2012)¹¹. The map was created by filling the saturated regions of the 24 μm MIPS map with empirically scaled 21 μm values. The star formation rate density is only well-correlated with the dense gas surface density in the central portion of W51 A. At lower gas surface densities, the SFR and gas are nearly anticorrelated. This offset is due to an age difference: the gas is tracing star formation yet to begin up to ~ 1 Myr, while the infrared emission traces older (0−100 Myr, with a peak at 5 Myr; Kennicutt & Evans 2012) star formation. The small-scale anticorrelation is consistent with the 5−20′ scales required to recover a star formation law from Galactic plane survey observations (Vutisalchavakul et al. 2014) because of the Kruijssen & Longmore (2014) uncertainty principle for star formation. The improved correlation at the highest densities also implies that the massive star formation threshold is closer to 10⁵cm^-3 than the 10⁴cm^-3 often measured for nearby, low-mass star-forming regions.

In order to evaluate star formation efficiency as a function of dense gas fraction, a more complete assessment of the present-day or very recent, t< 5 Myr, star formation in the W51 clouds is needed.

4.2. The future evolution of W51

The low DGMF associated with the W51 B cloud indicates that it has a low star formation potential despite containing significant mass (M ≳ 1 × 10⁵M_⊙). The presence of a supernova remnant and old, diffuse H ii regions indicates that the cloud did previously (and recently) form stars, but is now being destroyed. The fact that this cloud contains 10−30% of the total mass of W51, but has only a tiny fraction of its total mass above the purported star forming thresholds supports this story.

This cloud is therefore a good region to examine the effects of different feedback mechanisms (radiation, ionization, and supernovae) in parallel. It may also be a good location to examine how star formation comes to an end in the presence of massive star feedback.

By contrast, the W51 Main region has a dense gas fraction ~ 1 in its center, or at least ~ 50% out to nearly 5 pc. It has presumably only formed a small fraction of its total potential. It exceeds all of the various star formation and massive star formation thresholds (e.g. Lada et al. 2010; Krumholz & McKee 2008; Kauffmann & Pillai 2010), and therefore is expected to form additional stars efficiently, up to ~ 10⁴ − 10⁶M_⊙.

4.3. Implications for extragalactic observations of H₂CO

Although W51 is one of the most massive and active GMCs in the galaxy, containing 7% of the present-day massive star formation galaxy-wide (Urquhart et al. 2014), its star-forming gas mass is predominantly at a moderate density, n ~ 5 × 10⁴cm^-3; there is very little gas above 10⁶cm^-3 even in W51 Main. There were no detections of H₂CO emission on the ~ 1.25 pc (50′′) scales observed.

The nondetection is not just a geometric effect, as the only gas within W51 that is capable of emitting is tightly associated with the deeply embedded UCH ii regions within W51 Main. The dense gas is spatially compact and very likely to reside along the line of sight toward the continuum sources independent of viewing angle. This geometry is confirmed by existing interferometric observations that also show only absorption (Martin-Pintado et al. 1985a; Arnal & Goss 1985) and observations in other tracers showing that the dense gas is associated on ~ 0.1 pc scales with the UCH ii regions (Zhang & Ho 1997; Zhang et al. 1998).

By contrast, in extragalactic observations of starburst galaxies, there have been detections of H₂CO emission on ~ 100 pc scales. Mangum et al. (2013) report detections of H₂CO1₁₀−1₁₁ emission in NGC 3079, IC 860, IR 15107+0724, and Arp 220 on ~ 10 kpc scales. The implied local column densities from their analysis are modest, but the densities are extreme: their observations imply that the local-scale chemical conditions are comparable to W51, but the densities are different.

The only location in which densities n ≳ 5 × 10⁵cm^-3 (comparable to n(Arp220), etc.) are observed in W51 are in the central W51 Main region. We do not observe emission because of the bright continuum background source. It is therefore not possible to explain a H₂CO-emission galaxy by constructing it from collections of UCH ii regions; such a galaxy would be continuum-bright and show only H₂CO absorption. Instead, they must be assembled from huge quantities of high-density, non-star-forming gas. The presence of this gas in turn implies that the density threshold for star formation must be higher in starburst galaxies. This result is in contradiction to the idea that Giant H ii Regions are the building blocks of starburst galaxies (e.g. Miura et al. 2014).

One interpretation of the difference between the Galactic and extragalactic H₂CO is that the H ii regions that form in starburst galaxies have their radio luminosity significantly suppressed. The most straightforward explanation for a lower luminosity is that they are smaller and optically thick, compressed by a much more massive and overpressured ISM. Tightly squeezed H ii regions such as these are locations where radiation pressure is likely the dominant form of feedback, with multiple photon scatterings transferring additional momentum to the surrounding medium as described in Murray et al. (2010); this scenario does not occur in the optically thin H ii regions in the Galaxy where ionized gas pressure dominates.

5. Conclusion

We have presented maps of the H₂CO1₁₀ − 1₁₁ and 2₁₁ − 2₁₂ and H77α and H110α lines covering the W51 star forming complex and used these maps to examine the geometry and density structure of the complex. For the recombination lines, see Appendix A.

The H₂CO1₁₀ − 1₁₁/2₁₁ − 2₁₂ line ratio was used to measure gas volume densities and dense gas mass fractions. The W51 protoclusters have the majority (> 70%) of their mass in gas with density n(H₂) > 5 × 10⁴cm^-3. The rest of the cloud has a small dense gas fraction, with f(n> 1 × 10⁴ cm^-3) ~ 10%. The W51 B cloud therefore appears to be at the end of its star-forming lifetime, while W51 A will continue to form stars efficiently in the future. The relative weakness of present and future star formation in the W51 B/C region suggests that the collect-and-collapse mechanism is operating inefficiently or not at all.

Present-day high-mass star formation is associated only with W51 A, in which most of the molecular gas has density n ≳ 1 × 10⁵cm^-3. The highest gas density is closely associated with bright mid-infrared emission. W51 B and the outskirts of the W51 cloud have some gas with n> 10⁴cm^-3, but these regions exhibit limited star formation and show infrared emission anticorrelated with the dense gas. The tighter correlation with massive-star driven star formation indicators at high densities suggests that the density threshold for high-mass star formation is higher than that for low-mass star formation.

The H₂CO lines and their ratios have also been used to constrain the geometry of the W51 GMC and the associated H ii regions. The Galactic-face-on view of W51 is presented in more detail than has previously been possible. Analysis of our H₂CO data lead to the following conclusions about the structure of the GMC:

• The high velocity68 km s^-1 cloud is infront of the51 km s^-1 cloud andthe rest of the W51 GMC complex.

• The most luminous clusters and associated HII regions are in between the 51 km s^-1 and 68 km s^-1 clouds.

• There is molecular absorption associated with W51 B both in front of and behind the W51 C supernova remnant; W51 C is therefore within the W51 B cloud.

It is possible that the foreground 68 km s^-1 cloud is falling into the spiral potential from the near side, interacting with the 51 km s^-1 cloud exiting the potential on the far side to produce the W51 protoclusters.

Finally, we did not detect any H₂CO emission throughout the entire W51 GMC. The nondetection on scales from ~ 1 − 100 pc implies that detections of H₂CO1₁₀ − 1₁₁ emission in other galaxies comes from gas that does not surround bright H ii regions. These galaxies may therefore have interstellar media dominated by very high-density gas (n(H₂) > 10^5.5cm^-3) that is not presently forming stars. The density threshold for star formation in these galaxies must therefore be larger than in the Galactic disk, confirming earlier empirical (Longmore et al. 2013) and theoretical (Krumholz et al. 2005; Hennebelle & Chabrier 2013; Padoan & Nordlund 2011) results that such a threshold cannot be universal. These results are inconsistent with a universal density threshold for star formation observed in studies of nearby clouds (Lada et al. 2010, 2012; André et al. 2013).

The data are made available in FITS cubes and images hosted at the CfA dataverse¹². The entire reduction and analysis process and all associated code and scripts are made available via a git repository hosted on github¹³, with a snapshot of the publication version available from zenodo¹⁴.

Code Bibliography:

• The GBT KFPA Pipeline https://safe.nrao.edu/wiki/bin/view/Kbandfpa/ObserverGuide

• aoIDL http://www.naic.edu/~phil/download/aoIdl.tar.gz

• gbtidl http://gbtidl.nrao.edu/

• astropy www.astropy.org (Astropy Collaboration et al. 2013)

• astroquery astroquery.readthedocs.org (http://dx.doi.org/10.5281/zenodo.11656)

• sdpy https://github.com/keflavich/sdpy

• FITS_tools https://github.com/keflavich/FITS_tools

• aplpy http://aplpy.github.io

• image-registration http://image-registration.rtfd.org

• wcsaxes wcsaxes.rtfd.org (>=0.3.dev409)

• pvextractor pvextractor.rtfd.org

• agpy https://code.google.com/p/agpy/

• pyspeckit pyspeckit.bitbucket.org (Ginsburg & Mirocha 2011)

• ipython http://ipython.org/ (Pérez & Granger 2007).

Online material

Appendix A: Continuum and RRL data

We compare various data sets to assess calibration uncertainties and provide details of reduction for archival purposes. The data presented in this section are suitable for comparisons of Galactic to extragalactic star formation rate measurements, for example, since they are among the largest angular scale maps of radio recombination lines available with sub-parsec resolution.

Appendix A.1: Comparison between GBT and GPA data

The Galactic Plane A survey (Langston et al. 2000) covered the Galactic plane at 14.35 GHz using the Green Bank Earth Station (GBES) 13.7 m telescope, with a reported FWHM beam size of 6.6′. The published images were released with a FWHM resolution of 8′. We compared our GBT continuum observations to theirs in order to determine whether a significant DC component is missing from our data. Because the GPA used 10deg long scans in Galactic latitude, it should fully recover all diffuse Galactic Plane emission. In the released brightness temperature maps, brightness down to a scale of 1.5deg is recovered. However, because the GPA data undersampled the sky (its 5′ steps between scans were larger than the Nyquist sampling scale of the 14.35 GHz beam), point source fluxes in the GPA are underestimated by 19% and flux on small angular scales may be unreliable.

We resampled the GPA image onto the GBT grid using cubic spline interpolation, then smoothed both data sets to 9.5′. There are image artifacts (particularly vertical streaking) in the GPA data that are diminished by this large smoothing kernel.

We compared the surface brightness in the GPA and GBT data, and found that the GPA data was ~ 0.2 K brighter than the GBT in the diffuse portion of the W51 Main region; the offset

is not consistent with a purely multiplicative offset (Fig. A.1). The GBT observed the W51 Main peak to be moderately brighter, which is likely a result of the sparse sampling in the GPA. The morphological agreement between the maps is imperfect, perhaps in part because of the small area mapped in our GBT data, though there also appears to be vertical (along a line of constant longitude) stretching of the W51 Main source in the unsmoothed GPA data that is not consistent with the GBT observations.

Appendix A.1.1: Comparison between Arecibo and Urumqi data

We compare the 6 cm continuum to the Urumqi 25 m data from Sun et al. (2007) and Sun et al. (2011a). Figure A.2 shows the comparison of the Urumqi data and the Arecibo continuum data smoothed to 9.5′ resolution. The Arecibo and Urumqi data agree well as long as the main beam efficiencies of the respective telescopes (0.5 and 0.67) are accounted for.

Appendix A.1.2: Comparison of GBT and Arecibo data

In order to compare the Green Bank and Arecibo continuum data, we converted the brightness temperature maps to Janskys assuming a beam FWHM of 50′′ for both surveys and central frequencies of 4.8 and 14.5 GHz for Arecibo and Green Bank respectively. Measured beam widths for both telescopes were ~ 49 − 54″, so the relative error from assuming the same beam size should be ≲10%. In this section, the target frequencies are referred to as S_{5 GHz} and S_{15 GHz} for brevity.

The data are well-correlated, with S_{5 GHz} ~ 1.4S_{15 GHz} (S_{15 GHz} ~ 0.7S_{5 GHz}; Fig. A.3), consistent with a spectral index α_ν = −0.3 slightly steeper than usually observed for optically

thin brehmsstrahlung and consistent with there being some contribution from synchrotron emission. The lower-brightness regions have a lower S_{15 GHz}/S_{5 GHz}, indicating that these regions are more affected by synchrotron. In Fig. A.4a, a great deal of structure in the S_{15 GHz}/S_{5 GHz} ratio is evident in the vicinity of W51 Main: the ratio is higher towards the continuum peaks, indicating that the peaks have higher free-free optical depths, or lower relative contributions from synchrotron emission, than their envelopes.

We additionally compare the radio recombination lines observed simultaneously with the continuum and H₂CO. Hydrogen RRLs are often extremely well-correlated with the continuum and are therefore good indicators of the calibration quality.

In Fig. A.4, we show the ratios between the two frequencies in RRLs and continuum and the line-to-continuum ratios at both frequencies. The line values are the integrated flux densities over the range 20 to 100 km s^-1, which includes all Hα emission but no Heα.

The ratios between the x and y axis in each plot in Fig. A.3 are fitted using a total least squares approach with uniform errors for each data point. The line-to-continuum ratio is L/C(H77α) ~ 0.15 and L/C(H112α) ~ 0.04; in both cases there is little evidence for deviation from a linear relationship.

Appendix A.1.3: Comparison of the RRL and continuum data

Radio recombination lines are generally observed to be well-correlated with the corresponding radio continuum, particularly at low frequencies. At 5 and 15 GHz, the population level departure coefficients are close to 1, b_n> 0.95 (Wilson et al. 2009; Walmsley 1990).

While radio recombination lines are purely thermal in nature, the large-scale continuum may include a contribution from synchrotron emission. The morphological similarity between the 90 cm and 4 m (meter – i.e., 74 MHz) images presented by Brogan et al. (2013) and our 6 and 2 cm data hint that synchrotron emission could be significant. However, the high degree of correlation between the 2 and 6 cm described below suggest that synchrotron contamination is minor at both wavelengths.

Figure A.3 shows a comparison between the integrated RRL surface brightness and radio continuum at both 2 and 6 cm¹⁵. The figure shows the total least squares best-fit slopes to the data assuming uniform error, which yield a measurement of the line-to-continuum ratio.

We use the line-to-continuum ratio in both bands to measure the electron temperature using Eq. (14.58) of Wilson et al. (2009), which assumes a plane-parallel, optically-thin emission region with lines formed in local thermodynamic equilibrium (the ∗ in ${\mathit{T}}_{\mathrm{e}}^{\mathrm{\ast }}$ is meant to indicate these three assumptions are made). The two lines yield consistent measurements, with mean ${\mathit{T}}_{\mathrm{e}}^{\mathrm{\ast }}\mathrm{~}\mathrm{7000}\mathrm{-}\mathrm{8000}$ K; these measurements are consistent with smaller-scale measurements using the VLA with H92α (Mehringer 1994). There is little structure in the ${\mathit{T}}_{\mathrm{e}}^{\mathrm{\ast }}$ maps, with a hint of higher temperatures around G49.1-0.4, coincident with the W51C supernova remnant. Other structures are most likely due to the limited S/N.

Finally, we fit a single-component Gaussian to each pixel to produce velocity maps. These are discussed in Sect. B.11.

Appendix A.2: Carbon and helium RRLs

Helium RRLs were prevalent and reasonably well-correlated with the hydrogen RRLs, but we did not examine them in detail. He77α is detected at much higher S/N than than He107-112α. There were no clear detections of C77α or C107-112α, though there is a possible C77α signal at G49.366-0.304 with v_lsr ≈ 55 km s^-1 and a possible detection toward W51 Main along the wing of the He77α line. The He77α line detections are associated with regions of high H_nα but not regions of different ${\mathit{T}}_{\mathrm{e}}^{\mathrm{\ast }}$.

Appendix B: Geometry of individual regions

Appendix B.1: W51 main and W51 IRS 2

The W51 main and IRS 2 spectra show that both have ionized gas components at v_LSR ~ 55 km s^-1. This velocity approximately coincides with the peak of the ¹³CO emission.

The H₂CO1₁₀ − 1₁₁ spectra are deepest at ~ 68 km s^-1, while the 2₁₁ − 2₁₂ have depths approximately equal between the ~ 58 km s^-1 and ~ 68 km s^-1 components. The 55−60 km s^-1 components are too deep to be entirely behind the H ii regions. This indicates that the 55 km s^-1 ionized gas must be embedded within the molecular cloud, with molecular gas on both sides of the ionized gas along the line of sight.

Because these are well-studied regions, the low spatial resolution H₂CO spectra we present here add little new information about the gas kinematics. However, all of the velocity components observed in the W51 region are apparently kinematically connected to the W51 clusters.

Appendix B.2: The W51 B filament

The W51 B filament (right side of Fig. 11 at ~ 68 km s^-1), exhibits bright CO emission (${\mathit{T}}_{\mathrm{A}}^{\mathrm{\ast }}\mathrm{~}\mathrm{30}\mathrm{-}\mathrm{50}$ K in the Parsons et al. 2012 CO 3−2 data) but has relatively weak H₂CO absorption. The absorption models are inconsistent with the molecular gas being in front of the continuum emission, so Fig. 11 shows the continuum sources in front of the cloud at lower ℓ. Figure 8 shows an example model fit with the continuum assumed to be in front and in back, illustrating that the best-fit model parameters with continuum in the back do not reproduce the data. The relative positioning of the molecular gas behind the H ii regions suggests that the molecular gas is also behind the W51 C supernova remnant.

Appendix B.3: The edge of W51 C

W51 C is a supernova remnant that spatially overlaps with the W51 B star forming region. Brogan et al. (2013) argue that the supernova remnant must be in front of the H ii region G49.20-0.35 because the H ii-region has not absorbed all of the 4 m (74 MHz) nonthermal emission. The G49.1-0.4, G49.0-0.3, and G48.9-0.3 regions, however, show 4 m absorption signatures and may be in the foreground. There are clumps aligned along the 68 km s^-1 filamentary cloud with very high CO and H i velocities (Koo & Moon 1997b,a; Brogan et al. 2013), indicating that the SNR is interacting with the molecular gas.

The clumps at G49.1-0.3, ~ 68 km s^-1 are either lower density (n< 1.5 × 10⁴cm^-3) and in the background of the H ii region or high density (n> 1.5 × 10⁵cm^-3), low-column density and in the foreground. The 62 km s^-1 clumps have densities a few times higher, n ~ 4 × 10⁴cm^-3, and are clearly in the foreground of the continuum emission because their absorption depths are ~ 2.5 K, which cannot occur for absorption against the CMB. Figure B.1 shows a model spectrum fitted assuming the continuum lies between the two molecular velocity components. The relative strength of the ¹³CO and the H₂CO also suggests that the 68 km s^-1 component is behind the continuum.

We are seeing molecular gas both in front of and behind the supernova. This geometry can be readily confirmed by looking for molecular absorption at much lower frequencies where the SN synchrotron emission dominates over the H ii region free-free emission, i.e. the 335 and 71 MHz p - H₂CO lines.

Appendix B.4: G49.20-0.35 and G49.1-0.4

Tian & Leahy (2013) focus on the H ii regions G49.20-0.35 and G49.10-0.40 (called G49.10-0.38 in their work) to determine the relative geometry of the W51 C SNR and the W51 B H ii/star-forming region. They observe that the high-velocity H iis not detected toward either of these sources, indicating that the H ii regions must be behind the high-velocity H i features.

We detect H₂CO1₁₀ − 1₁₁ at ~ 58 and ~ 63 km s^-1 toward G49.10-0.40, with line ratios that are consistent with the H ii region being behind the molecular cloud complex. It also has an extreme RRL velocity, v_110α ≈ 72 km s^-1, the most redshifted seen in the entire W51 region (see Figs. B.3 and B.2).

G49.20-0.35 is also clearly behind the molecular cloud, as evidenced both by H₂CO absorption depth and the IRDC absorption in the foreground. It has an RRL velocity v_110α ≈ 70 km s^-1.

Because both H ii regions are extremely redshifted, they are most likely associated with the W51 B cloud complex, contrary to the interpretation by Tian & Leahy (2013) in which they are unrelated background clouds. The Galactic rotation curve doesn’t allow for velocities red of ~ 60 km s^-1, and almost none of the molecular gas exceeds ~ 70 km s^-1 even on the wings. The H ii regions are therefore probably shooting out the back side of the molecular cloud in a champagne flow, perhaps accelerating ionized gas from the ~ 68 km s^-1 component further to the red.

Appendix B.5: The 66 km s^-1 8 μm dark cloud

Between W51 A and W51 B, there is a component of the 68 km s^-1 cloud that is filamentary and in the foreground of all of the free-free emission. This cloud component is evident as an IRDC in the Spitzer GLIMPSE images from ℓ = 49.393, b = −0.357 to ℓ = 49.207, b = −0.338; it is labeled as IRDC G49.37-0.35 in Figs. 11 and 12.

The H ii region G49.20-0.35 is clearly behind the IRDC, though there are strong morphological hints that it is interacting with and truncated by the cloud.

Appendix B.6: G49.27-0.34

The UCH ii region G49.27-0.34, which was considered a candidate extended green object (EGO) and subsequently rejected for lack of H₂ emission (De Buizer & Vacca 2010; Lee et al. 2013), exhibits a second velocity component at ~ 68 km s^-1, slightly but clearly redshifted of the rest of the IRDC. The dust component contains a gas mass ~ 2 × 10³M_⊙ based on the BGPS flux and using the assumptions outlined in Aguirre et al. (2011), suggesting that the high velocity could be due to infall or virialized gas within a deep potential. The virialized velocity width, given the radius and mass from the BGPS data, is σ_vir = 8.8 km s^-1, while the measured H₂CO linewidth is FWHM(H₂CO) = 7.2 km s^-1, wider than in any other part of the cloud except W51 Main.

Both radio continuum and RRLs are detected toward this source. The H77α RRL velocity is ~ 58 km s^-1, significantly blueshifted from the molecular gas. The H₂CO lines do not independently distinguish between the continuum source being in the front or back of the cloud, but the mean density from the BGPS mass and radius n ~ 2.5 × 10⁴cm^-3 is within a factor of 2 of the H₂CO-derived density, n ~ 1.4 × 10⁴cm^-3, if the continuum source is behind the gas, while the H₂CO-derived density is too low, n ~ 2 × 10³cm^-3 if the continuum source is in front.

The implied geometry therefore has the H ii region behind the molecular gas, plowing toward it at a velocity difference Δv ~ 10 km s^-1. Such a high velocity difference may indicate that the H ii region is confined by the molecular gas and on a plunging orbit into the cloud.

Appendix B.7: G49.4-0.3f, aka G49.34-0.34, aka IRAS 19209+1418

The H ii region centered at 49.34-0.34 was identified by Mehringer (1994) as part of the G49.4-0.3 complex. There are 3 distinct H₂CO line components at 51, 63.70, and 68.47 km s^-1. The 51 km s^-1 component is behind the H ii region; the ¹³CO line is detected at comparable brightness at 51 km s^-1 and 63 km s^-1, while the H₂CO1₁₀ − 1₁₁ line is ~ 10 × deeper at 63 km s^-1. The RRLs associated with this source are at v_LSR = 58 ± 1 km s^-1.

The H₂CO lines are moderately well-fit by the two-velocity-component model, but there is a relative excess of 2₁₁ − 2₁₂ absorption at 66 km s^-1 (associated with the 68 km s^-1 component). The extra absorption may indicate that there is a high-density, low-column component at this velocity.

The 8 μm GLIMPSE image shows that the 68 km s^-1 IRDC crosses in front of this source. Herschel Hi-Gal 70 μm images reveal a ring structure that is hinted at in the 8 μm image. There is no evidence for interaction between the ring feature and the IRDC. This intriguing feature will likely be difficult to study in detail because the dusty, molecular gas feature lies in front of it.

Appendix B.8: G49.4-0.3

The collection of H ii regions around G49.4-0.3 vaguely resembles a cartoon mouse. As noted in Carpenter & Sanders (1998), the molecular gas in this region is separated into two distinct components, one at 51 km s^-1 and the other at 64 km s^-1. The 64 km s^-1 component is in the foreground, while the 51 km s^-1 is in the background of most of the H ii regions.

Both cloud components are in the foreground of the central H ii regions at G49.36-0.31, the eyes of the mouse. The density of the 51 km s^-1 component is an order of magnitude higher than that in the 64 km s^-1 component in this region, suggesting that the gas is being compressed by the H ii region. The clean separation between the 64 and 51 km s^-1 cloud components suggests that they are not interacting at this location.

Based on the absorption line depths, the G49.38-0.30, IRAS 19207+1422, and G49.37-0.30 H ii regions are behind the 51 km s^-1 cloud. The 8 μm absorption features are associated with the 64 km s^-1 cloud and are in front of all of the H ii regions.

The 8 μm morphology of G49.42-0.31 is bubble-like, so it is plausible that the H ii region is neither in front nor behind the 51 km s^-1 cloud but embedded within it, blowing a hole in the cloud.

Appendix B.9: 8 μm dark cloud G49.47-0.27

The cloud to the north of W51 Main/IRS2 appears as a dark feature in Spitzer GLIMPSE 8 μm maps. It is detected in H₂CO from 54 to 64 km s^-1. Throughout, it has a high 1₁₀ − 1₁₁/2₁₁ − 2₁₂ ratio, ≳7 in most voxels, indicating a low density n ≲ 10³cm^-3.

Centered at 60.6 km s^-1, the region has a line FWHM 5−7 km s^-1, indicating that it is quite turbulent, with 3D Mach number in the range 10 < ℳ < 20 for an assumed 10 <T< 20 K. At its centroid velocity, it is connected to the W51 main cloud.

There is a previously unreported bubble HII region in the north part of this cloud, which we designate G49.47-0.26, with radius ~ 70″ (1.7 pc). The H ii has RRL velocities v_lsr ≈ 50 km s^-1. Because it is not detected in Brackett γ emission (from the UWISH2 survey: Froebrich et al. 2011), it is most likely behind the cloud.

The Kang et al. (2009)Spitzer survey of YSOs in the region indicates that there are no YSOs within the boundaries of this cloud; it is very likely non-star-forming at present (see also Fig. 16).

Because the cloud is continuous with the W51 Main region in velocity and is 8 μm-dark, it is most likely at the same distance as W51 Main and associated indirectly with the massive cluster forming region.

Appendix B.10: The 40 km s^-1 clouds

There are clouds observed at 40 km s^-1 that show only weak H₂CO absorption spread across nearly the entire region. These molecular clouds are behind nearly all of the H ii regions in the W51 complex. There are additional 40 km s^-1 clouds clearly seen in H i absorption (Stil et al. 2006) that are not associated with these molecular clouds, but instead represent a foreground population of neutral atomic medium clouds.

Appendix B.11: Kinematic maps

Maps showing the overall kinematics of the region are shown in Figs. B.3 and B.4. Figure B.3 shows the velocity at peak absorption of the H₂CO1₁₀ − 1₁₁ line and the fitted radio recombination line centroid velocity. Figure B.4 shows the best simultaneous fit to the H₂CO1₁₀ − 1₁₁ and 2₁₁ − 2₁₂ absorption features over two different velocity ranges. The 1₁₀ − 1₁₁ absorption velocity in Fig. B.3a approximately shows the velocity of the front-most molecular clouds along the line of sight at each position.

Appendix B.12: Abundances

The LVG modeling also yields measurements of abundance that are degenerate with the assumed velocity gradient. The abundance within the LVG model is defined as $\mathit{X}\mathrm{=}\frac{\mathit{N}\mathrm{\left(}{\mathrm{o}\mathrm{-}\mathrm{H}}_{\mathrm{2}}\mathrm{CO}\mathrm{\right)}\mathrm{[}{\mathrm{cm}}^{-2}{\mathrm{􏼮}{\mathrm{(}{\mathrm{km}\mathrm{s}}^{-1}\hspace{0.17em}{\mathrm{pc}}^{-1}\mathrm{)}}^{\mathrm{]}}}^{}}{\mathit{n}\mathrm{\left(}{\mathrm{H}}_{\mathrm{2}}\mathrm{\right)}\left[{\mathrm{cm}}^{-3}\right]}\left(\frac{\mathrm{1}{\mathrm{km}\mathrm{s}}^{-1}\hspace{0.17em}{\mathrm{pc}}^{-1}}{\mathrm{pc}}\right)\mathrm{·}$(B.1)\newpage\noindentWe show the distribution of fitted abundances under this definition in Fig. B.5. The histogram shows the average abundance along each line of sight derived from the likelihood-weighted density divided by the likelihood-weighted column. We also show a two-gaussian fit to the abundance distribution: 30% of the area is consistent with an abundance X(o - H₂CO) = 10^{− 8.4 ± 0.5} and 70% with X(o - H₂CO) = 10^{− 9.9 ± 0.6}. We caution that these abundance measurements are highly uncertain and are contingent on both the backlighting source brightness and the assumed velocity gradient. The individual abundances going in to the histogram are generally not well-constrained. A more accurate abundance measurement could be obtained by measuring the millimeter lines of o - H₂CO at 140 and 150 GHz.

Appendix B.13: Supplementary picture

We include a large-scale WISE + BGPS composite image of the entire W51 cloud and its surroundings. The figure highlights the W51 C supernova remnant.

---

Acknowledgments

We thank Xiaohui Sun for providing the Urumqi 6 cm Stokes I image prior to its availability on the survey website. The paper benefitted from discussions with Jonathan Tan, Neal Evans, and Diederik Kruijssen. We thank our referee Jeff Mangum for a helpful and rapid referee report. We are grateful to the editor, Malcolm Walmsley, for additional comments and for catching many small but significant errors. E.R. is supported by a Discovery Grant from NSERC of Canada. This research has made use of the VizieR catalogue access tool, CDS, Strasbourg, France. The original description of the VizieR service was published in (Ochsenbein et al. 2000). This publication makes use of molecular line data from the Boston University-FCRAO Galactic Ring Survey (GRS). The GRS is a joint project of Boston University and Five College Radio Astronomy Observatory, funded by the National Science Foundation under grants AST-9800334, AST-0098562, AST-0100793, AST-0228993, & AST-0507657. This research has made use of the NASA/ IPAC Infrared Science Archive, which is operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with the National Aeronautics and Space Administration.

References

All Figures