BGPS XIII. Physical Properties

The Bolocam Galactic Plane Survey. XIII. Physical Properties and Mass Functions of Dense Molecular Cloud Structures

Timothy P. Ellsworth-Bowers11affiliation: CASA, University of Colorado, UCB 389, Boulder, CO 80309, USA 22affiliation: email: , Jason Glenn11affiliation: CASA, University of Colorado, UCB 389, Boulder, CO 80309, USA , Allyssa Riley11affiliation: CASA, University of Colorado, UCB 389, Boulder, CO 80309, USA , Erik Rosolowsky33affiliation: Department of Physics, 4-183 CCIS, University of Alberta, Edmonton, AB T6G 2E1, Canada , Adam Ginsburg44affiliation: European Southern Observatory, Karl-Schwarzschild-Straße 2, 85748, Garching bei München, Germany , Neal J. Evans II55affiliation: Department of Astronomy, University of Texas, 2515 Speedway, Stop C1400, Austin, TX 78712, USA , John Bally11affiliation: CASA, University of Colorado, UCB 389, Boulder, CO 80309, USA , Cara Battersby66affiliation: Harvard-Smithsonian Center for Astrophysics, 60 Garden Street, Cambridge, MA 02138 USA , Yancy L. Shirley77affiliation: Steward Observatory, University of Arizona, 933 North Cherry Avenue, Tucson, AZ 85721, USA 88affiliation: Adjunct Astronomer at the National Radio Astronomy Observatory , and Manuel Merello99affiliation: Istituto di Astrofisica e Planetologia Spaziali-INAF, Via Fosso del Cavaliere 100, I-00133 Roma, Italy

We use the distance probability density function (DPDF) formalism of Ellsworth-Bowers et al. (2013, 2015) to derive physical properties for the collection of 1,710 Bolocam Galactic Plane Survey (BGPS) version 2 sources with well-constrained distance estimates. To account for Malmquist bias, we estimate that the present sample of BGPS sources is 90% complete above and 50% complete above . The mass distributions for the entire sample and astrophysically motivated subsets are generally fitted well by a lognormal function, with approximately power-law distributions at high mass. Power-law behavior emerges more clearly when the sample population is narrowed in heliocentric distance (power-law index for sources nearer than 6.5 kpc and for objects between 2 kpc and 10 kpc). The high-mass power-law indices are generally for various subsamples of sources, intermediate between that of giant molecular clouds and the stellar initial mass function. The fit to the entire sample yields a high-mass power-law . Physical properties of BGPS sources are consistent with large molecular cloud clumps or small molecular clouds, but the fractal nature of the dense interstellar medium makes difficult the mapping of observational categories to the dominant physical processes driving the observed structure. The face-on map of the Galactic disk’s mass surface density based on BGPS dense molecular cloud structures reveals the high-mass star-forming regions W43, W49, and W51 as prominent mass concentrations in the first quadrant. Furthermore, we present a 0.25-kpc resolution map of the dense gas mass fraction across the Galactic disk that peaks around 5%.

Subject headings:
Galaxy: structure Ð- ISM: clouds -Ð stars: formation – submillimeter: ISM
slugcomment: Accepted Version: April 6, 2015

1. Introduction

As the possible precursors to stellar clusters, OB associations, or smaller stellar groups, molecular cloud clumps and cores have become a primary focus for understanding the process of high-mass star formation (e.g.,  McKee & Ostriker, 2007, and references therein). The recent advent of large-scale continuum surveys of the Galactic plane at (sub-)millimeter wavelengths (BGPS, Aguirre et al., 2011, Ginsburg et al., 2013; ATLASGAL, Schuller et al., 2009; Hi-GAL, Molinari et al., 2010b) have detected tens of thousands of these objects in thermal dust emission. A major payoff of blind surveys in this portion of the electromagnetic spectrum is the derivation of physical properties of regions hosting high-mass star formation. A detailed census of these dense molecular cloud structures can help constrain star-formation and galactic-evolution theories (e.g.,  Kennicutt & Evans, 2012). Studies of the physical properties (e.g.,  Peretto & Fuller, 2010; Giannetti et al., 2013) and mass distributions (e.g.,  Netterfield et al., 2009; Olmi et al., 2013; Gómez et al., 2014) of these objects have begun, but despite the richness of the current data sets, a coherent picture has not yet emerged for the evolution of the dense interstellar medium and the origin of the stellar initial mass and stellar cluster mass functions.

While giant molecular clouds (GMCs) have been studied for several decades using the lowest rotational transitions of the isotopologues of CO and other simple molecules (e.g.,  Scoville & Solomon, 1975; Cohen et al., 1980; Dame et al., 1987, 2001, and references therein), and studies of the apparently uniform stellar and cluster initial mass functions (e.g.,  Bastian et al., 2010) go back to Salpeter (1955), the observational technology (detectors and angular resolution) has only recently developed to study the dense substructures intermediate between these two extremes. And while studies of nearby molecular cloud complexes have recently yielded estimates of the core mass function (e.g.,  Motte et al., 2007; Enoch et al., 2007; Swift & Beaumont, 2010) for the possible progenitors of single stellar systems, it has only been with the recent large-scale blind surveys that the study of the clump mass function has been possible.

Theoretical modeling of molecular cloud structure evolution is beginning to place constraints on the clump mass function (e.g.,  Donkov et al., 2012; Veltchev et al., 2013), and two primary functional forms are discussed in the literature: the power-law and lognormal distributions. Both arise from physical processes in molecular clouds; supersonic turbulence within molecular clouds produces a lognormal density distribution (e.g.,  Padoan et al., 1997, although Tassis et al., 2010 find other means to produce a lognormal), while gravitational collapse of dense structures tends to produce a power law (e.g.,  Padoan & Nordlund, 2002). The distinctions between these have implications for competing theories of high-mass star formation (Elmegreen, 1985). Furthermore, the observed mass function for dense molecular cloud structures should resemble a combination of these forms due to the complex interaction of physical processes (e.g.,  Offner et al., 2014); for a discussion of the details, see Hopkins (2013).

In a broader context, the Galactic distribution of dense molecular gas and star formation has implications for using the Milky Way as ground truth for studies of extragalactic star formation and galaxy formation. Cosmological simulations of galaxy formation and evolution have made remarkable progress in the last two decades reproducing some of the observed properties of galaxies (e.g.,  Katz et al., 1996; Kereš et al., 2009; Stinson et al., 2013), but numerical simulations have largely been unable to a priori produce the strong winds and instantaneous inefficiency of star formation observed in galaxies (Hopkins et al., 2014). By incorporating feedback from stars via ‘sub-grid’ physics recipes (e.g.,  Kim et al., 2014), numerical galaxy evolutions models can now reliably produce spiral galaxies. Very recently, Hopkins et al. (2014) introduced simulations that do include realistic feedback to produce galaxies commensurate with locally observed examples. The continued refinement of these models must rely upon the constraints of derived physical properties of dense molecular cloud structures in the Milky Way. Additionally, Galactic structures can provide a comparison for scaling relationships for star formation rates with respect to more active galaxies (Kamenetzky et al., 2014).

As the first completed dust-continuum millimeter-wavelength survey of the northern Galactic plane, the Bolocam Galactic Plane Survey (BGPS) has provided a foundation for studying dense molecular cloud structures in a variety of Galactic environments. The majority of objects detected in the BGPS correspond to molecular cloud clumps (Dunham et al., 2011, hereafter D11), each of which will likely form a cluster of stars.111Molecular cloud clumps are typically  with sizes  pc (Bergin & Tafalla, 2007). Following D11, we generally adopt the [][]Bergin:2007 classifications for dense molecular cloud structures. Derivation of physical properties of these objects requires robust distance estimates; Ellsworth-Bowers et al. (2013, 2015) developed a distance probability density function (DPDF) methodology for utilizing all available information towards this goal. By combining kinematic distance information with prior DPDFs based on ancillary data (such as Spitzer/GLIMPSE mid-infrared images and trigonometric parallax measurements) in a Bayesian framework, Ellsworth-Bowers et al. (2015, hereafter EB15) presented a catalog of 1,710 BGPS objects with well-constrained distance estimates.

In this work, we build on that distance catalog to derive the physical properties of and fit mass functions to this vital population of sources. The probabilistic description of the distance to each object in the survey catalog may be used to directly propagate uncertainties into the derived physical properties via Monte Carlo methods. We present here a catalog of physical properties for this set of sources, and fit mass functions to the entire sample, as well as astrophysically motivated subsets. Furthermore, we explore the extension of fundamental molecular cloud scaling relationships (Larson, 1981) to dense substructures and map the distribution of dense molecular gas throughout the Galactic disk.

2. Data

The Bolocam Galactic Plane Survey version 2 (BGPS V2; Aguirre et al., 2011; Ginsburg et al., 2013, hereafter G13), is a  mm continuum survey covering 192 deg at 33″ resolution. It is one of the first large-scale blind surveys of the Galactic plane in this region of the spectrum, covering with and selected larger cross-cuts to , plus selected regions in the outer Galaxy. For a map of BGPS V2 coverage and details about observation methods and the data reduction pipeline, see G13. From the BGPS V2 images, 8,594 millimeter dust-continuum sources were identified using a custom extraction pipeline. BGPS V2 pipeline products, including image mosaics and the catalog, are publicly available.222Available through IPAC at

The BGPS data pipeline removes atmospheric signal using a principle component analysis technique that discards common-mode time-stream signals correlated among bolometers in the focal plane array. The pipeline attempts to iteratively identify astrophysical signal and prevent that signal from being discarded with the atmospheric signal. This sky-subtraction behaves like a high-pass filter that may be characterized by an angular transfer function passing scales from approximately 33″ to 6′ (see G13 for a full discussion). The effective angular size range of detected BGPS sources therefore corresponds to anything from molecular cloud cores up to entire GMCs depending on heliocentric distance (D11).

Ellsworth-Bowers et al. (2015) presented a distance catalog for BGPS V2 sources using a Bayesian DPDF framework. The DPDFs describe the relative probability of finding an object at points along the line of sight and are constructed from multiple sources of distance-related information via


where is the kinematic distance likelihood function based on a line-of-sight velocity (), and the are prior DPDFs based on ancillary data. Given circular orbits about the Galactic center, there are two points along a line of sight looking interior to the Solar circle that correspond to a single , the kinematic distance ambiguity (KDA). The for such objects is double-peaked, representing the equal probability of finding the object at either kinematic distance; prior DPDFs are used principally to resolve the KDA, placing an object at either the near or far distance (see EB15 for further discussion). The resulting catalog contains 1,710 sources with well-constrained distance estimates (distance uncertainties generally  kpc), hereafter referred to as the Distance Catalog. The posterior DPDFs333Available with the BGPS release at IPAC. are normalized to unit integral probability () to facilitate marginalization over distance, and allow for the Monte-Carlo propagation of uncertainty in an object’s distance to the derived physical properties by randomly sampling the DPDFs many times.

3. Method

3.1. Computing Physical Properties Using the DPDF

It is straightforward to compute physical properties of dense molecular cloud structures using single-value distance estimators. The complete information contained in an object’s posterior DPDF, however, represents the combination of knowledge and ignorance about its distance based on a wide collection of data encapsulated in the suite of prior DPDFs applied. Use of this full information, therefore, requires careful consideration. Below, we discuss mass calculation for BGPS objects in some detail (§3.1.1), and then apply those methods to the computation of other physical properties (§3.1.2).

3.1.1 Mass Derivation

A simple estimate of a dense molecular cloud structure’s mass may be computed from optically thin millimeter dust continuum data via


where is the source-integrated flux density,  is the estimated heliocentric distance, is the opacity per mass of dust, is the gas-to-dust mass ratio, and is the Planck function evaluated at dust temperature . For the specific case of the BGPS,


where is the  mm source-integrated flux density,  cm g of dust (Ossenkopf & Henning, 1994, Table 1, Column 5),444This table of represents dust grains with ice mantles, coagulating in cold, dense molecular regions at  cm for  yr. This is appropriate for molecular cloud cores, and may be an upper limit for less-dense clumps (e.g.,  Martin et al., 2012). The computed masses may, therefore, be biased low. See a recent study of dust opacity for dense molecular gas by Juvela et al. (2015). = 100 (Hildebrand, 1983), and the equation has been normalized to  K. Standard error-propagation methods may be used to determine the mass uncertainty, but the asymmetric nature of many DPDF error bars (EB15) complicates the issue. The uncertainty in the distance may be marginalized over to produce the expectation value of the mass via


where is proportional to the second moment of the DPDF. The bimodal nature of many posterior DPDFs in the inner Galaxy, however, causes this to be a biased estimator of the mass, in much the same way  (the first-moment distance) can be a biased estimator of the distance. For instance, an object with a well-constrained distance estimate may have up to of the integrated DPDF in the non-favored kinematic distance peak; the difference between  and the maximum-likelihood distance () can be greater than 1 kpc, or the uncertainty in  (see Ellsworth-Bowers et al., 2013 for a complete discussion).

This bias, coupled with the measured source flux-density and assumed dust temperature uncertainties, suggests using a Monte Carlo approach for estimating an object’s mass and associated uncertainty. To determine the maximum-likelihood mass for a given catalog object, Equation (2) is computed many times with each realization drawing a distance, flux density, and temperature from the suitable probability density functions. The resulting (properly normalized) mass probability density function (MPDF) may then be used to estimate the maximum-likelihood mass () and the 68.3% error bars in the same fashion as done with DPDFs in EB15.

Being the probability of an object having a mass between and , the MPDF of a given object is created from the histogram of masses produced by many () Monte Carlo realizations of Equation (2). Given the fractional nature of the uncertainties in the derived mass, we chose to use logarithmic binning for computation of MPDFs and the probability density functions for other physical properties. Note that the computation of DPDFs (EB15) is done with respect to a fixed linear distance scale and is not subject to the need to create histograms from Monte Carlo calculations.

3.1.2 Physical Radius and Number Density

In addition to the masses of dense molecular cloud structures, the catalog of heliocentric distances enables the computation of physical sizes and number densities of these objects. Not only can these physical properties help to disentangle the various populations of objects (i.e., cloud, clump, core) detected by the BGPS and others, but they can also test the extension of various empirical relationships for GMCs down to more dense substructures.

The physical radius is computed from


where is the deconvolved angular radius of a detected object computed from geometric mean of the deconvolved major () and minor () axes of the flux density distribution,


(Rosolowsky et al., 2010). For the BGPS, the rms size of the beam is given by , and is a factor relating the rms size of the emission distribution to the angular radius of the object.555The appropriate value of to be used depends on the true emission distribution of the object and its size relative to the beam. Rosolowsky et al. (2010) computed for a large range of models, varying the emissivity distribution, angular size relative to the beam, and signal-to-noise ratio. The chosen value for the BGPS catalog is the median value from these simulations, but the variations span more than a factor of two. Use of a specific source model in conjunction with the BGPS catalog would require the appropriate value of and rescaled catalog values for Objects whose measured minor axis dispersion is smaller than do not have a real solution to Equation (6); actual physical objects should not fall into this category, but effects of the cataloging process lead to 29% of the full BGPS V2.1 catalog and 20% of the Distance Catalog having a complex computed . These objects are excluded from the analyses of physical radius and number density in §4.

As a means for classifying dense molecular cloud structures, the mean number density offers possible insight into the physical processes, such as structure formation via turbulence or gravitational collapse, that may be at play in a given catalog object. The number density computed with assumed cylindrical geometry (i.e., a circle on the sky with radius , and depth ),666This parameterization of the volume results in ; assuming spherical geometry increases the volume by a factor of 4/3, and assuming a cylinder of depth increases it by a factor of 2. As shown in §4.1, the number density is a poorly constrained quantity for these objects, and source-to-source geometry differences play a significant role. , may be parameterized in terms of mass and radius or observable quantities as


The first parameterization relates the physical quantities in units appropriate to dense molecular cloud structures, while the second is specific to the observational quantities of the BGPS and the assumptions used in Equation (3). Note that the number density is not as strong a function of heliocentric distance as is the mass; incorrect KDA resolution will have a smaller effect on the derived value.

3.2. Constructing Mass Distributions and Estimating Completeness

The mass function describing the relative numbers of dense molecular cloud structures is of particular interest for connecting star-formation theory to observation. Functional forms are fitted to the distribution of masses and we explore some of the possible mass distribution constructions. Additionally, since the range of masses over which conclusions may be drawn is restricted by completeness considerations, we present an analysis estimating the mass completeness levels for a heterogenous data set such as the BGPS.

3.2.1 Mass Distributions

Figure 1.— Comparison of mass computation methods. Left: The mass distribution computed by applying Equation (2) to each object, using the appropriate single distance estimate from the posterior DPDF (as discussed in Ellsworth-Bowers et al. (2013), ‘Simple ’, red) plotted atop the mass distribution using the maximum-likelihood masses from the MPDFs (‘Maximum-Likelihood ’, black). Also shown is the Monte Carlo composite mass distribution, as discussed in §4.2.1 (hashed cyan). Middle: Comparison of with for each object in the Distance Catalog. The red dashed line marks the 1:1 relationship. Right: Histogram of the ratio for each object in the Distance Catalog, with Gaussian fit overlaid in red.

In its simplest form, the mass distribution of dense molecular cloud structures may be compiled from the direct application of Equation (2) to each object with a well-constrained distance estimate, where  is either the maximum-likelihood distance or the first-moment distance, as discussed in Ellsworth-Bowers et al. (2013). This form of the distribution for the Distance Catalog sources is shown as the red histogram in Figure 1 (left) and represents a naïve, but effective, realization of the mass distribution.

While the simple mass distribution utilizes the best distance estimate from the DPDF formalism, it discards much of the information contained in the posterior DPDFs. Given the asymmetric nature of many DPDFs, the most straightforward way to capitalize on this pool of information is through Monte Carlo realizations utilizing the MPDFs. Plotted in black in Figure 1 (left) is the mass distribution computed using the  from the Monte-Carlo generated MPDFs. While there are noticeable deviations in the wings of the distributions, a two-sample KS test between these two distributions reveals that they are indistinguishable to better than 99.9%. Despite the negligible differences between the distributions, the power of the MPDF-based mass distributions lies in the ability to create many Monte Carlo realizations of the mass distribution, randomly drawing a mass from each object’s MPDF for each trial.

If the aggregate mass distributions between the two methods are indistinguishable, what of the masses of individual objects? The comparison of the mass derived from a simple application of Equation (2) using the single-value distance estimate from EB15 () with the maximum-likelihood value from the MPDF () is shown in the middle panel of Figure 1, where the dashed line marks the 1:1 relationship. The histogram of the ratio of these masses () is plotted in the right panel along with a Gaussian fit. The mean of the distribution is consistent with unity, and has a standard deviation of . Differences between the masses are more significant for low-mass objects (), with becoming larger than the simple mass by up to 60% in the most extreme case.

Going beyond constructing a mass distribution for the Distance Catalog using a single mass for each object, a more complete Monte Carlo approach is shown as the hashed cyan histogram in Figure 1 (left). Rather than taking the maximum-likelihood value from the MPDF of each object (as does the black histogram in that panel), each object’s MPDF was randomly sampled 20,000 times to produce a mass distribution containing 34.2 million entries. The resulting histogram was then normalized to represent the objects of the Distance Catalog. The resulting distribution is much smoother at both mass extrema because of the lack of small-number statistics. Interestingly, the mode (peak) of the distribution is shifted to slightly higher mass, likely due to the dependence in Equation (2). This Monte Carlo composite distribution is discussed further in §4.2.1.

3.2.2 Mass Completeness Function

Figure 2.— Computed mass completeness levels. Cumulative distribution of the mass corresponding to  rms flux-density noise and the source distance from the posterior DPDF. The cyan curve represents the ‘Simple’ distribution, using the appropriate single distance estimator from the DPDF (see Ellsworth-Bowers et al., 2013). The black curve is the aggregate of Monte Carlo realizations, with each realization randomly drawing a distance from each object’s posterior DPDF and dust temperature from the lognormal distribution discussed in §4.1.1. Dotted horizontal lines identify 50% and 90% completeness. Light vertical dashed and dot-dashed lines identify the 50% completeness levels for both distributions, and the dark vertical lines mark the 90% completeness levels.

(Sub-)Millimeter continuum surveys of thermal dust emission are flux-limited, so Malmquist bias must be addressed when analyzing the mass distribution of detected objects. It is relatively straightforward to identify flux-density completeness levels, but since higher mass (and therefore brighter) objects are detectable through a larger Galactic volume than low-mass objects, it is imperative to define the mass criteria at which a survey is complete to some level. While small surveys of isolated regions at known distances (e.g.,  Ridge et al., 2006; Enoch et al., 2006; Rosolowsky et al., 2008) can directly estimate strict mass completeness levels, heterogenous surveys of the entire Galactic plane must be content to have a more complex criterion.

We cast the issue of completeness levels in terms of a mass completeness function, which describes the fraction of objects in the survey area at any given mass that should be detected. Mass completeness begins with flux-density completeness, usually defined in terms of the rms noise in survey maps. For instance, the BGPS V1 catalog is 99% complete at the rms noise flux density for any given field (Rosolowsky et al., 2010), a value that is unchanged in version 2 (G13). While surveys that do not have uniform noise properties across the Galactic plane pose an additional complication to the process, variable backgrounds throughout the Galactic plane (especially at shorter wavelengths) will also affect the completeness calculation. The flux-density completeness level must, therefore, be estimated on a field-by-field777The BGPS has a variable rms noise from field to field (see G13, their Figures 1 and 2). (or, better, source-by-source) basis. Estimates of rms noise in survey images, however, usually do not take the confusion limit into account, so these completeness levels should be viewed as lower limits.

To begin, we compute the local flux-density completeness level for each catalog object as the mean of the BGPS noise maps within the object’s catalog boundary (see Rosolowsky et al., 2010). The minimum complete flux density for each object implies that any source at this location in the survey map brighter than this value will be detected of the time. Approximately 13% of the objects in the Distance Catalog have flux densities below their own completeness value, meaning that they are detections. These objects represent a tip of the iceberg; similar objects exist in those same fields but were not recovered.

Combining the completeness flux density with the object’s distance information allows an estimate of the “minimum complete mass” for each object, meaning that any object at the specified () with mass higher than this value will be detected of the time. Repeating this process on a source-by-source basis allows for the effective marginalization over the variable noise properties across the Galactic plane. The mass completeness function is therefore the cumulative distribution of the minimum complete mass for all objects in the sample.

Figure 2 illustrates the mass completeness functions computed using the two techniques described in Figure 1. The completeness function for the single distance estimate from EB15 for each object is plotted in cyan, and an aggregation of Monte Carlo realizations pulling distances from each object’s DPDF is shown in black. The curves are similar, but a two-sample KS test reveals they are not drawn from the same population at greater than 99.5% significance. Shown as horizontal dotted lines are the 50% and 90% completeness levels, and the vertical lines mark the intersections of the curves with these levels. The 50% completeness level is for both curves, and the 90% completeness level is . The importance of these values is in limiting conclusions about the mass distributions (or other derived physical properties) to larger masses; extending conclusions to smaller-mass objects is not supported by the available data.

3.3. Fitting a Functional Form to the Mass Distribution

To compare the observed mass distribution for dense molecular cloud structures with theory and other observations, we fit both power law and lognormal functions to the data. The formulae and procedures for this fitting are discussed below, and are based on maximum-likelihood methods to avoid the pitfalls of working with binned data.

3.3.1 Power Law

The use of power-law functions for describing the mass distributions of stars and their precursors goes back to the original studies of the stellar initial mass function by Salpeter (1955). The mass function is defined as the number of objects either per logarithmic mass interval,


or per linear mass interval,


where , and (Chabrier, 2003). Within this framework, the Salpeter power-law index is , . For a detailed description of this parameterization and its consequences, see Swift & Beaumont (2010), Olmi et al. (2013), and references therein. The mathematical nature of a power law requires a finite lower limit on the mass range over which it is a valid descriptor of the data in addition to the observational constraint that the mass distribution turns over at (e.g.,  Kroupa, 2001).

To fit a power-law function to the observed mass distribution, we use the maximum-likelihood method described by Clauset et al. (2009)888 to estimate the power-law index and range over which it is valid. Briefly, assuming and , the probability (or likelihood) that the observed masses () are drawn from a power-law distribution is proportional to


where is defined as in Equation (9). The maximum-likelihood value is computed by taking the derivative of the logarithm of the likelihood ()999The log-likelihood is used for mathematical expediency, as the product in Equation (10) becomes a sum. While the numerical value of the function is changed by taking the logarithm, the location of the maximum is not. with respect to and finding the root (i.e.). The power-law index estimator is


The minimum mass () for which a power law is a good descriptor of the mass distribution is computed from the data. For each value in the mass distribution, Equation (11) is computed using that mass as . The cumulative distribution of is compared with the analytic power law with index , and the maximum difference between them is computed (i.e., the K-S test -statistic); the value of that minimizes is returned by the algorithm (see Clauset et al., 2009 for a complete description of the algorithm). The results presented here were computed using a python translation of the PLFIT routine from Clauset et al.101010

3.3.2 Lognormal

Predicated on the properties of supersonic turbulence and the turnover in the IMF at , the lognormal form of the mass function is also commonly fit to the mass distributions of dense molecular cloud structures (e.g.,  Swift & Beaumont, 2010; Olmi et al., 2013). Analogous to the power-law case, we use a maximum-likelihood method for computing the best fit to the mass distribution. The likelihood of the observed masses being drawn from a lognormal distribution is proportional to


where the normalization factor


and erfc is the complementary error function evaluated at . The parameters [] are the mean and width of the lognormal Gaussian, respectively, and [] are related to the minimum and maximum masses, respectively, over which the lognormal fit is valid. The maximum-likelihood estimators are found by numerically solving (i.e., numerically maximizing ), where


As with the power-law case, the limits over which the lognormal fit best describe the data must be estimated from the data themselves. We employ the same K-S test -statistic analysis as described above, but use a numerical minimizer to explore the 2-dimensional parameter space including . Following the algorithmic structure of PLFIT3.3.1), we wrote a routine in python to maximize Equation (14) using Powell’s Method inside a numerical minimization of the K-S test -statistic using a downhill simplex optimization.111111Both optimization steps utilized routines from the scipy library (

4. Results

4.1. Physical Properties of BGPS Sources

4.1.1 Catalog

BGPS V2.1 Catalog Properties Derived from the PDFsaaSee §3.1.
Catalog bbThe source-integrated  mm flux density, uncertainties in parentheses. ccThe deconvolved radius of the catalog object, computed using Equation (6). ddThe appropriate distance estimate from EB15;  for objects near the tangent point and  otherwise.
Number (°) (°) (Jy) (″) (kpc) () (pc) (log cm)
2235 7.993 74.3
2254 8.187 + 35.9 0
2256 8.207 + 56.7 0
2261 8.249 + 77.6 0
2262 8.263 + 79.8 0
2265 8.281 + 23.8 0
2295 8.545 35.5 0
2296 8.551
2297 8.579 37.3 0
2300 8.669

Note. – This table is available in its entirety in a machine-readable form.

Table 1Physical Properties of BGPS Distance Catalog Sources

We present a physical properties catalog for the 1,710 sources in the Distance Catalog. The mass, physical radius, and number density for each object were computed from the Monte-Carlo probability density functions as described in §3.1, where a dust temperature was drawn for each realization from a lognormal distribution with a mean of 20 K and full-width at half-maximum of 8 K (Battersby et al., 2011). This choice of temperature distribution was recently affirmed by the analysis of analogous dense molecular cloud structures from the ATLASGAL survey by Wienen et al. (2015). We note that the study of Merello et al. (2015) combining  mm and  µm data finds a lower mean temperature for BGPS sources ( K for dust spectral index , and  K for ). Using a mean  K would lead to an underestimate of the mass of each object by a factor of 1.3 to 1.8 for the and , respectively, of Merello et al. (2015). While the computed mass and number density values are affected by this uncertainty in temperature, Monte Carlo simulations of BGPS source mass distributions by Schlingman et al. (2011) showed the power-law slope of the mass function (§4.2.2) is not very sensitive to the choice of underlying dust temperature distribution.

The pdfs for mass, physical radius, and number density were constructed in log space to capture the spread in values over orders of magnitude, and error bars were computed to enclose 68.3% of the total probability, where the endpoints occur at equal probability (see the discussion of DPDF error bars in EB15). The error bars approximate the Gaussian region, but account for the asymmetric nature of the pdfs. The catalog of physical properties is shown in Table 1, along with relevant information from the BGPS V2.1 catalog (G13) and heliocentric distance (EB15).

Figure 3.— Brief analysis of resolved versus unresolved BGPS V2 sources. Top: Scatter plot of deconvolved radius versus integrated  mm flux density. Evenly spaced contours represent the density of points per 0.1 dex  0.1 dex bin. Only resolved sources are shown. Plotted in green are Equation (17) evaluated at  = 1 kpc (solid), 5 kpc (dotted), and 20 kpc (dashed) (see §5.2 for discussion). Bottom: Histograms of integrated  mm flux density. Black represents all resolved sources (i.e., the projection of the points in the top panel), and the gray diagonally hashed histogram represents the unresolved sources. Red horizontally hashed and blue solid histograms represent the resolved and unresolved populations, respectively, of the Distance Catalog (EB15).

Of particular note in the catalog are the 20% of sources in the BGPS Distance Catalog whose angular extents on the sky do not yield real solutions to Equation (6) for . A brief analysis of these objects in the context of the flux-density distribution for all catalog sources reveals a strong trend whereby objects with non-physical are preferentially dim (Figure 3). The gray and blue histograms represent the “unresolved” sources in the full V2.1 catalog and the Distance Catalog, respectively, while the black and red histograms represent the resolved sources in the same groups.

The obvious correlation in the data points in the top panel (Spearman rank correlation coefficient ) is a direct consequence of a universal scaling relationship in the dense interstellar medium first identified by Larson (1981) whereby the mass of a molecular cloud structure is roughly proportional to the square of its radius. Section 5.2 explores this relationship in more detail, culminating in Equation (17), which is plotted as green lines in Figure 3 for three different heliocentric distances.

The Larson relationships in the figure extend to low flux-density, but since the cataloging process truncates source boundaries where the flux density meets the local rms value (Rosolowsky et al., 2010), source masks do not extend to the theoretical zero-flux-density isophot. The resulting catalog source may have an extent in the map smaller than the beam size, leading to a complex solution to Equation (6). To further aggravate the situation, the beam deconvolution step assumes that sources are spherical (or at least have circular projections on the sky), which is invalid for many of the fainter filamentary structures throughout the Galactic plane (e.g.,  Molinari et al., 2010a). Future cataloging routines may be constructed that divide emission between filamentary and compact sources, and may decrease the number of BGPS sources with invalid deconvolved angular extent. For the present, however, unresolved sources are excluded from the following analyses that rely on physical radius or number density, leaving a sample of objects with a tabulated physical size.

4.1.2 Ensemble Physical Property Distributions

Figure 4.— Physical property histograms for the Distance Catalog, computed using the maximum-likelihood values from the Monte Carlo method described in §3.1. The gray filled and magenta hashed histograms represent sources above the 50% and 90% mass completeness levels, respectively (see §3.2.2). Top: Source mass distribution. Vertical lines mark the locations of the 50% (gray) and 90% (magenta) completeness levels. Middle: Physical radius distribution. The vertical dashed lines identify the boundaries between cloud / clump / core as used by Dunham et al. (2011). Bottom: Number density distribution; vertical lines as in the middle panel.

The ensemble physical property distributions from Table 1 are shown in Figure 4 to illustrate the range of BGPS objects in the Distance Catalog. This subsample is generally representative of the entire BGPS catalog, except that detection in a molecular transition line that traces dense gas (e.g., HCO(3-2)) is strongly correlated with continuum flux density; dim sources are largely absent from the Distance Catalog (see EB15 for a complete discussion). In all three panels, sources above the 50% (70 ) and 90% (400 ) mass completeness levels are indicated by gray solid and magenta hashed histograms, respectively.

The mass distribution (top panel) peaks around 100 , but the completeness levels indicated by the filled histograms imply that the location of the distribution peak is observationally biased. The location of the 50% completeness level demonstrates that the predicted turnover at low mass is not constrained by the current observations. Physical radius is illustrated in the middle panel of Figure 4, with vertical dashed lines marking plausible boundaries between molecular clouds (large ), clumps (intermediate ), and individual cores. We place approximate radius boundary markers at 1.25 pc between cloud and clump and 0.125 pc between clump and core, the same boundaries used by D11 which were based on the general guidelines from Bergin & Tafalla (2007). These divisions, while rooted in physical distinctions based on virial ratio and capacity for Jeans fragmentation, are flexible and the boundaries are easily blurred by the continuum of structure in molecular cloud complexes (see §5.2), and by variations in Galactic environment. Note, however, the near complete lack of Distance Catalog objects with radii  pc). From the relationship in Figure 3, smaller sources have lower flux density and are therefore less likely to be included in the Distance Catalog (EB15).

Definitions notwithstanding, the typical size scale of BGPS-identified dense molecular cloud structures is clearly in the vicinity of 1 pc and may encompass multiple true but unresolved molecular cloud clumps. The mass completeness subsets illustrate the effects of Malmquist bias (coupled with the Larson relationship connecting mass and physical radius) whereby the 90% complete sample is nearly all cloud-scale ( pc) objects, and the BGPS does not detect any large ( pc), lower-mass objects (). This places a rough constraint on the minimum number density ( cm, Equation 7) detectable by the BGPS for large-scale objects.

The bottom panel of Figure 4 depicts the number density distribution of resolved sources, which spreads across a wide dynamic range. Approximate number density boundaries at  cm between cloud and clump, and  cm between clump and core (D11) are indicated by vertical dashed lines. Whereas there is nearly complete overlap between high-mass and large-radius objects (‘clouds’; middle panel), there exist a handful of low-mass () objects with low computed number density (‘clouds’; bottom panel). These objects, along with the spread of higher-mass objects () across classifications in the bottom panel insinuate that number density computed via Equation (7) for BGPS sources is not a reliable indicator of physical classification. It is possible that non-uniform geometry of sources is causing the spread of the completeness samples over a large range, or that the algorithm by which the cataloging routine decomposes emission into discrete sources is doing so in a manner inconsistent with the underlying physical structure in molecular clouds (e.g., by assuming spherical clouds when computing (H)).

4.1.3 Physical Properties as a Function of Heliocentric and Galactocentric Position

Figure 5.— Physical properties as a function of heliocentric distance and Galactocentric radius. Grey dots mark individual objects from the BGPS Distance Catalog, and black squares identify the median values in 2-kpc () or 1-kpc () bins. Top Row: Computed mass for each object. The hashed region represents the mean flux-density completeness level as a function of heliocentric distance. For the right panel, this region is computed for  along the line. Horizontal dashed lines mark the 50% and 90% completeness levels (as in Figure 4). Middle Row: Physical radius. The lower hashed region corresponds to , while the upper marks . For the right panel, the lower region is again computed for  along . Horizontal dashed lines mark the approximate cloud / clump / core boundaries (as in Figure 4). Bottom Row: Number density. Horizontal dashed lines as above.

A major challenge in characterizing the ensemble physical property distributions of detected molecular cloud structures is separating observational and data-processing systematic effects from underlying physical relationships, as systematic effects on the measurement of physical properties can lead to significant biases in the derived quantities. Understanding the physics of dense molecular cloud structures requires well-constrained distance estimates. It is instructive, therefore, to analyze the functional dependences.

The left panels of Figure 5 illustrate the dependence on heliocentric distance of the mass, physical radius, and number density for sources in the Distance Catalog. Sources are plotted as gray dots, with the median value per 2-kpc bin marked by black squares. The systematic dependence on heliocentric distance can be clearly seen in the median values as , , and , following Equations (2), (5), and (7), respectively. For the distribution of mass (top left panel), the hashed region identifies the mass associated with the mean flux density completeness level ( Jy, Section 3.2.2). Especially given the variable rms noise from field to field in the BGPS, sources may be detected in this region but completeness is not assured. The horizontal dashed lines mark the 50% and 90% completeness levels (as in Figure 4). The large collection of sources around   kpc has a median mass , whereas the median grows to for objects  kpc.

The physical radius plot (middle left panel) shows similar effects, with the hashed regions marking (typical minimum angular radius) and (approximate maximum recoverable radius; see G13). The horizontal dashed lines represent the  pc and  pc divisions between physical object categories discussed above. Finally, the number density distribution (bottom left panel) displays the expected trend, except that the median values remain somewhat constant between 500 cm cm over a large range of heliocentric distance. Since this flattening is not observed for physical radius (which also depends linearly on ), this range of may represent a typical value for BGPS dense molecular cloud structures, subject to the geometric and cataloging caveats discussed above.

The right panels of Figure 5 illustrate the dependence of these physical properties on Galactocentric radius. The vertical dashed line identifies our adopted  kpc (Reid et al., 2014) for reference. Black squares mark the median values in 1-kpc bins, and roughly show the systematic effects of heliocentric distance (i.e., objects far from  = are at large ). Gaps in the black function represent bins with  sources, where the median value is skewed by outliers. The gap in sources (gray dots) at  =  kpc is pointed out in EB15 as arising from a unique combination of BGPS coverage region and a kinematic avoidance zone (where kinematic distances are unreliable in the face of non-circular motions about the Galactic center), and is not likely a true feature of the underlying Galactic dense molecular gas distribution. Many of the dense molecular cloud structures observed in the Molecular Ring / Scutum-Centarus Arm feature at   kpc,   kpc, correspond the physical properties typical of molecular cloud clumps.

The systematic effects of heliocentric distance will obscure some of the underlying physical trends in properties as a function of Galactocentric radius. For instance, the spike in median number density near   is due to the population of nearby (  kpc) core-scale objects. The observed decrease in mass for dense molecular cloud structures from  = 4 kpc to  = 7 kpc is consistent with the trend observed for GMCs by Roman-Duval et al. (2010), but is likely primarily an observational bias.

4.2. Mass Function Fits

4.2.1 Mass Distribution Trials

The full information available from the MPDFs (and ultimately the DPDFs from EB15) allow for the utilization of a large number of Monte Carlo trials to fit functional forms to the dense molecular cloud structure mass distribution. For each trial, a single value was drawn from the MPDF of each object (§3.1.1) to create an independent realization of the mass distribution of the Distance Catalog. Both a power law (§3.3.1) and lognormal (§3.3.2) function were fit to each realization, and the fit parameters recorded. This process was repeated times, a value chosen to balance appropriately sampling the parameter spaces for the functional fits with computation time. The collected fit parameters from the Monte Carlo process were then plotted as two-dimensional histograms to identify joint confidence intervals. Below, we discuss the results of the power-law (§4.2.2) and lognormal (§4.2.3) fits to the entire data set. Attempts to subdivide the Distance Catalog into astrophysically meaningful subsets are described in Appendix A, but most fits resemble those of the entire sample. Maximum-likelihood fit parameter values for all fits are listed in Table 2.

Figure 6.— Power-law fit to the aggregate mass distribution. Left: The complementary cumulative distribution function (CCDF, solid black curve) and differential mass distribution (gray histogram) for the Distance Catalog. This distribution is the aggregate of 20,000 Monte Carlo realizations of the mass function, normalized to a sum of to create an average mass distribution that marginalizes over the uncertainties in the mass of each object. The solid cyan and dashed black curves represent the maximum-likelihood power-law index (from right panel), and the vertical gray dashed line marks the maximum-likelihood minimum mass for which a power law describes the data. Right: The two-dimensional histogram of the power-law fit parameters for each of the 20,000 Monte Carlo trials, whereby a power law was fit to a single realization of the mass distribution. Contours enclose 68.3% (highest), 95.4%, and 99.7% of the points in the plane. The yellow marks the largest concentration of points, or the joint maximum-likelihood values of and .

4.2.2 Power-Law Fit

The collected power-law fits to the complete mass distribution are shown in Figure 6. The left panel illustrates the mass distribution both as the complementary cumulative distribution function (CCDF, solid black)121212The CCDF is simply 1 - CDF, where CDF is the typical cumulative distribution function. The CCDF measures the probability of finding an object with mass greater than a given value, whereas the CDF describes the probability of finding an object with mass less than a given value. The CCDF is used to illustrate the distribution in a manner robust against fluctuations caused by finite sample size (Clauset et al., 2009). and as the differential mass distribution (gray histogram). Both distributions represent the aggregate of 20,000 Monte Carlo realizations of a single mass distribution, normalized so the sum equals 1,710 objects. To smooth the effects of small sample size and utilize the full amount of information encoded in the collected set of MPDFs, the aggregate is shown rather than any single one of the Monte Carlo trials. The power-law fits (derived from the data in the right panel) are shown for the CCDF (cyan line) and the differential mass distribution (black dashed line).

The power-law fitting algorithm (PLFIT) returns two parameters, the power-law index and the minimum mass over which the data may be described by a power law (). The two-dimensional histogram of the 20,000 parameter pairs from the Monte Carlo trials are shown in the right panel of Figure 6. Color scale intensity indicates the number of points in each pixel of the histogram, and contours enclose 68.3% (highest contour), 95.4%, and 99.7% of the parameter pairs, meant to illustrate joint confidence regions. The high degree of correlation between and is indicative of the continuously steepening nature of the mass distribution shown in the left panel (where a power law is a straight line). The yellow in the right panel marks the peak of the two-dimensional histogram, or the joint maximum-likelihood parameter values (). This minimum mass is above the 90% completeness level described in §3.2.2 (), meaning the power-law fit is not affected by Malmquist bias, although other selection effects may be at play. Note the secondary peak of points near in the right panel; there appears to be a sharp steepening of the mass distribution at the high-mass end (see §5.3.2 for discussion). The distribution truncates around due to a requirement in the fitting algorithm for a minimum number of data points required for a robust fit (see Clauset et al., 2009). This truncation is a general feature of the astrophysical subset fits discussed in Appendix A.

4.2.3 Lognormal Fit

Figure 7.— Lognormal fit to the aggregate mass distribution. Left: The CCDF (solid black curve) and differential mass distribution (gray histogram) are identical to Figure 6. The blue and black dashed curves represent the maximum-likelihood lognormal fit (from the right panel), and the vertical gray dashed lines mark the maximum-likelihood minimum and maximum mass for which a lognormal describes the data. Right: The two-dimensional histogram of the lognormal fit parameters [] for each of the 20,000 Monte Carlo realizations of the mass distribution. Contours enclose 68.3% (highest), 95.4%, and 99.7% of the points in the plane. The yellow marks the largest concentration of points, or the joint maximum-likelihood values. The inset shows the parameter space distribution for the pair [].

The continuously steepening nature of the mass CCDF makes fitting or interpreting a power law function for the entire data set difficult at best. Furthermore, distance inhomogeneities in large-scale surveys such as the BGPS tend to drive the composite observed mass distribution toward a lognormal (see §5.3.3). We applied our lognormal fitting routine to the Monte Carlo realizations of the mass distribution, and the resulting fit parameter distributions are shown in Figure 7. The CCDF and differential mass distributions in the left panel are identical to Figure 6 (left). The functional fits are drawn from the right panel as was done with the power law, where the cyan line fits the CCDF and the dashed black line fits the differential mass distribution.

Whereas the power-law fit returns two parameters, the lognormal fitting optimizes the functional fit over four parameters: the lognormal characteristic [], and also both ends of the mass range [] over which the lognormal fit is valid. For computational expediency, we divided the fitting into a nested structure following the recipe from Olmi et al. (2013, Appendix B). The lognormal characteristic is fit for a given pair of mass range values, then the mass range is adjusted to find the best match between the model CCDF and the CCDF constructed from the source masses. Shown in the right panel of Figure 7 are the distributions of lognormal best-fit parameters paired in this fashion. The more fundamental is the [] plane, shown in the background, since these define the position and shape of the function. As with the power-law fit, the two-dimensional histogram of the 20,000 Monte-Carlo trial fit parameters is shown, with contours representing the same enclosed fractions of points. Both parameters are tightly constrained, with () and , where the uncertainties represent the marginalized 68.3% confidence intervals. The joint maximum-likelihood values (marginalized over [], yellow ) occur at , and are reflected in the solid cyan and dashed black curves in the left panel. Since the peak occurs at 200 , below the 400  90% completeness level, it may be biased and should be viewed as an upper limit.

Inset in the right panel is the two-dimensional histogram of the [] plane, showing the relationship between the bounding parameters. There is scant correlation, except that the spread along the diagonal from top-left to bottom-right represents expansion and contraction about a central value. The optimum values of these parameters (marginalized over []) are indicated by the yellow in the inset (), and are marked in the left panel by vertical gray dashed lines. This minimum mass extends down to the 50% mass completeness level (), and its location may not be entirely reliable, as the peak of the fitted lognormal function always falls below the 90% completeness level ().

Comparison of the lognormal fit in Figure 7 with the power-law fit in Figure 6 yields some interesting notions. As expected, the minimum mass for the power-law fit is well above the peak of the lognormal distribution. In fact, the 99.7% contour in Figure 6 (right) terminates at about this value. Furthermore, the curvature of the mass distribution is well-matched to a lognormal function over nearly 2.5 orders of magnitude. The continual curve of the mass distribution may be an effect of observational and cataloging biases against the highest mass sources (see §5.3.2).

Power-Law Fit Lognormal Fit
Cut SubsetaaSee Appendix A for definitions of these subsets. bbThe lognormal fitting routine tends to find mass bounds very near the initial guess ( ,  ), but parameter correlations show that the exact bounds have little impact on the [] returned. bbThe lognormal fitting routine tends to find mass bounds very near the initial guess ( ,  ), but parameter correlations show that the exact bounds have little impact on the [] returned.
Type Name () () () () ()
Distance Catalog 1710
Nearby 1328
Distant 382
Sct-Cen Arm 866
Sgr / Local Arms 610
Sgr Arms 330
Per / Outer Arms 234
Clump / Core 807
Cloud 562
Clump / Core 993
Cloud 376
Mean  pc 1444
Mean  pc 358
Peak  pc 392
Peak  pc 1410
 Type Dense Gas 1372
CO 422
“Protocluster” Blind () 1513
 pc 1331
Mixed Protocluster 579
Table 2Mass Distribution Functional Fit Parameter Values

4.2.4 Mass Functions of Astrophysically Motivated Subsets

Theories of star formation hold that a strong power-law behavior should be evident for dense molecular cloud structures, particularly leading toward the stellar IMF (Krumholz et al., 2005; Hennebelle & Chabrier, 2011, 2013; Padoan & Nordlund, 2011; Federrath & Klessen, 2012). The lack of a clear single power law fit to the complete Distance Catalog prompted us to fit functional forms to astrophysically motivated subsets of the data to identify homogenous sets of objects predicted by theory. While the fits to many of the subsets, did not reveal any strongly compelling behavior (see Appendix A), there are two whose power-law fits are fairly well constrained. The set of nearby objects (  kpc) may be described by an power law, and set of criteria aimed at isolating molecular cloud clumps (2 kpc    10 kpc and ) is fit by a power-law. The heliocentric distance cut removes cloud-scale objects that are more strongly affected by the angular filtering function of the BGPS data pipeline, leaving principally molecular cloud clumps, similar to the mixed criteria.

The overall lack of clean fits to either the complete Distance Catalog or most of the subsets indicates that the sample is largely homogenous and the interstellar medium is not cleanly, tightly described by simple relations. These effects are discussed in depth in §5.3. The results of the studied astrophysical subsets are presented in Appendix A, along with the accompanying figures similar to Figures 6 and 7. The fit results are all shown in Table 2.

4.3. The Distribution of Dense Molecular Gas in the Galactic Plane

4.3.1 Face-On View of the Milky Way

Figure 8.— Dense gas mass surface density over the Galactic disk for Distance Catalog sources. Left: Color scale and contours indicate the dense gas surface mass density ( kpc) in each 0.25 kpc  0.25 kpc pixel. The Sun and Galactic Center are identified by stars at  kpc and  kpc, respectively. Right: The same mass surface density contours plotted atop an artist’s rendering of the Milky Way based on Spitzer data (R. Hurt: NASA/JPL-Caltech/SSC). Several of the mass concentrations correspond to prominent star-formation complexes originally discovered by Westerhout (1958, marked with red ‘’s). As discussed in EB15, we have excluded from consideration of much of   kpc for kinematic reasons; sources appearing inside this circle have robust distance estimates through association with a trigonometric parallax measurement.

Ellsworth-Bowers et al. (2015, their Figure 15) presented a map of the locations of dense molecular cloud structures in the Galactic disk with respect to a rendering of Galactic structure based on Spitzer data. Since these objects span several orders of magnitude in mass (Figure 4), their locations alone do not necessarily indicate the mass distribution of dense molecular gas. To estimate the Galactic disk mass surface density of the objects in the BGPS Distance Catalog, we constructed a grid consisting of 0.25 kpc  0.25 kpc bins into which masses computed for individual sources were placed. Since typical distance uncertainties are  kpc, using 0.25 kpc bins approximately Nyquist samples the data. Taking advantage of the information contained in source DPDFs, we produced a map (Figure 8, left) by randomly drawing distances for each object, and adding the mass derived using Equation (3), divided by to reflect the probabilistic contribution of that mass to the appropriate grid cell. This use of the full DPDF leads to some of the smearing effects along the line of sight towards isolated regions in the outer Galaxy (e.g.), where peculiar motions have a larger effect on kinematic distances. Furthermore, given that BGPS coverage for is neither blind nor contiguous, lack of mass in Figure 8 in the outer Galaxy should not be construed to imply a void in the Galactic plane.

Mass surface density contours are plotted atop an image of the Galaxy produced using Spitzer near-infrared data (R. Hurt: NASA/JPL-Caltech/SSC) in the right panel of Figure 8. As with the corresponding map from EB15, there are many caveats in the interpretation of this image. First, there are virtually no BGPS sources with well-constrained distances at   kpc owing to non-circular motions induced by the Galactic bar,131313The DPDF formalism relies primarily upon kinematic distances whereby velocity measurements of dense molecular gas tracers are mapped onto Galactic position using a rotation curve assuming circular motion about the Galactic center. The Galactic bar induces strong radial streaming motions, rendering regions under its influence unusable for kinematic distance measurement. with the exception of objects associated with trigonometric parallaxes measurement of maser transitions in regions of high-mass star formation (e.g.,  Reid et al., 2014). Second, given that the typical uncertainty for well-constrained distance estimates from EB15 is  kpc, the map and contours have a  kpc effective resolution, significantly coarser than the underlying image, although commensurate with Herschel maps of nearby galaxies (mean distance  Mpc; Kennicutt et al., 2011).

Third, there is a significant concentration of mass near  kpc that appears to lie in a void in the underlying image. As discussed in EB15, this could represent either molecular gas at that location not identified in the Spitzer data (which is a “best guess” picture of the Galaxy), or sources in the W43 star-forming complex at the end of the Galactic bar erroneously placed at the far kinematic distance. Comparison of Figure 8 with the Galactocentric locations of GMCs in the CO Galactic Ring Survey (GRS; Jackson et al., 2006) reveals a similar large collection of gas at this location (Roman-Duval et al., 2010, their Figure 5). Most of the BGPS sources in this region are associated with one or more H ii regions with robust KDA resolution from the H ii Region Discovery Surveys (HRDS; Bania et al., 2010, 2012) and the distances to GRS clouds were primarily resolved using the same H i absorption techniques as the HRDS (EB15).

The BGPS-GRS correlation implies either the existence of a large collection of mass in an apparently empty region of the Spitzer image or a fundamental limitation in the application of H i absorption techniques for distance resolution. Evidence in support of the former comes from both Russeil et al. (2011), who find interarm clumps for this region in the Hi-GAL SDP field, and a study of the M51 spiral galaxy that reveals a population of GMCs downstream of a major spiral arm (Egusa et al., 2011). Furthermore, the appearance of interarm gas could also be the effect of a spiral feature possessing a systematic noncircular velocity, biasing kinematic distances.

Identified in the right panel of Figure 8 are the high-mass star-forming regions W43, W51 and W49, which stand out in the mass diagram as high surface density areas. Despite the rough tracing of some spiral structure, this map is not a definitive distribution of dense molecular gas in the Milky Way. First of all, objects used to produce this map represent only 20% of the full BGPS catalog (by number, although 33% by flux density). Second, given the coarse resolution of the map, any conclusions about detected spiral structure should be used cautiously as the map serves to show where the detected mass lies, and not a complete census of mass in the Galaxy. The gap observed between , however, is real and corresponds to a similar minimum in CO column seen in the Dame et al. (2001) maps. The expansion of robust distance determination to other dust-continuum surveys, however, will have high utility for tracing the Galactic dense molecular gas distribution across the entire disk.

4.3.2 Vertical and Radial Distributions of Star-Forming Mass

Figure 9.— Vertical and radial mass distributions. The color-scale image in the center represents the azimuthal integration of mass, computed in analogous fashion to Figure 8, left. Since bins at larger Galactocentric radius encompass larger volumes, the absolute value of the color-scale is not physically meaningful, but rather represents the relative concentrations of mass. At top is the mass distribution versus , marginalizing over vertical position, and at right is the distribution versus , marginalizing over Galactocentric radius. The fit to the mass-weighted histogram is consistent with the number-count histogram from EB15.

The azimuthally integrated vertical and radial distributions of the mass of star-forming clouds provide an alternate view of Galactic star-formation distribution, and are shown in Figure 9. The Galactocentric radius and vertical position were computed using the conversion matrix from Appendix C of Ellsworth-Bowers et al. (2013), which accounts for the 25 pc vertical solar offset above the Galactic midplane. As in Figure 8, masses and Galactocentric positions were derived from independent samples from the DPDF for each object. The color image at the center represents the amount of mass in each 0.25 kpc (radial) 10 pc (vertical) bin. Since this is an azimuthal integration, the color scale intensity of the image is physically meaningless other than to identify overdense regions; pixels at larger  represent a significantly larger volume than those nearer the Galactic center. Identified in the image are the W43, W51, and W49 regions that stand out in the face-on image of the Galaxy. These three regions dominate the mass distribution, but there are still clear collections of mass in the 4-5 kpc range, and around 6 kpc. The warping of the disk to large  is visible in the mass distribution, but we caution that BGPS observations in the outer Galaxy () were targeted towards regions of known star formation, and do not represent the same blind sample as the inner Galaxy coverage.

To the right of the image is the marginalization over , or the the vertical histogram of the mass distribution, shown as total mass per bin. A Gaussian fit to the this mass-weighted histogram ( pc, FWHM =  pc) is consistent with the distribution derived from the number-count histogram from EB15. This implies that the mass is roughly evenly distributed amongst catalog sources regardless of vertical position, and that there is no particular bias in the mass of dense molecular cloud structures near the Galactic midplane.

5. Discussion

5.1. What is a BGPS source?

Based on a subset of BGPS objects throughout the Galactic plane observed in the lowest inversion transitions of NH around 24 GHz, Dunham et al. (2010, 2011) began a discussion about the nature of BGPS sources. Using the larger sample of objects with well-constrained distance estimates from EB15, we continue this discussion.

Observationally, dense molecular cloud structures are typically defined in terms of cloud-, clump-, and core-scale objects. The BGPS is generally sensitive to clump-scale objects, but we investigate to what extent clouds and cores are present in the data. Figure 10 illustrates the relative concentrations of non-clump objects at various heliocentric distances, and is modeled after D11 (their Figure 21). Comparison between the two figures yields insight into the effects of using a larger sample. The thin black lines in Figure 10 mark the fraction of sources closer than a given distance with properties consistent with molecular cloud cores based on physical radius (solid;  pc), mass (dotted; ), and number density (dashed;  cm), as defined in §4.1.2. Many () of the objects at   kpc fall under the molecular cloud core definition, and conversely, the majority of core-scale objects are seen at   kpc. This is commensurate with the completeness and angular transfer functions illustrated in Figure 5.

Figure 10.— Fraction of source type as a function of heliocentric distance. Shown in thin black lines is the fraction of sources in the BGPS Distance Catalog nearer than the indicated heliocentric distance that are cores, and thick blue lines indicate the fraction of sources farther that are clouds. Line styles indicate core and cloud divisions as determined by physical radius (solid), mass (dotted), and number density (dashed). Both ends of the distance range suffer from small number statistics. (This figure is styled after D11, their Figure 21.)

Similarly, thick blue lines in Figure 10 mark the fraction of sources located farther than a given distance that are characterized as clouds ( pc, , and  cm). Of particular note is the disparity in the distributions between different definitions for GMC-scale objects. Beyond 6 kpc, of objects in the Distance Catalog have large physical radius. As found by D11 for the NH-observed subset, the density curve never reaches 90% because there exists a substantial fraction of high-density sources in the distance range 6 kpc    kpc (see Figure 5, bottom left). Furthermore, the geometrical uncertainties in both the cataloging process and derivation of physical properties decrease the reliability of this quantity for object classification. The combination of improved cataloging routines with higher fidelity to the underlying structure and future high angular resolution studies of the Galactic plane (e.g., with CCAT) may alleviate this disparity. One principal deviation from the D11 results, however, is the mass fraction of cloud-scale objects as a function of distance. Whereas for the NH subset the cloud fraction by mass closely tracked the radius fraction, in Figure 10 the cloud fraction by mass is significantly smaller. This may be due either to the targeted nature of the D11 observations or differences between the BGPS version 1 and version 2 catalogs (see G13).

With the enlarged set of physical properties presented here, we confirm the previous conclusions of studies of BGPS objects that the majority of sources detected are clump-scale objects in the distance range 2 kpc    kpc, cores nearer by, and cloud-scale objects beyond. Uncertainties in the geometry of dense molecular cloud structures remains one of the most significant issues in terms of understanding the number density distribution and its role in the evolution of dense regions.

5.2. Larson’s Laws & The Continuum of Structure

5.2.1 The Mass-Radius Relationship

The physical properties of molecular clouds and their substructures appear to follow a series of scaling relationships known as Larson’s Laws, following their introduction in a literature study of molecular clouds by Larson (1981). The most frequently used of these scaling relationships is the size-linewidth relationship, which holds that virial and/or turbulent motions across a molecular cloud structure increase with the square root of the size of that structure. This relationship is valid for single molecular line tracers (e.g.,  Solomon et al., 1987; Heyer et al., 2009), but the linewidths of different molecular transitions observed in the same dense molecular cloud structure may vary by an order of magnitude. The heterogenous nature of the spectroscopic data gathered for BGPS distance estimation therefore prevents their use in deriving such a relationship.

Figure 11.— The mass surface density - mass () relationship for molecular cloud structures. Black dots identify BGPS Distance Catalog objects from Table 1, while cyan squares mark the values for GRS GMCs from Roman-Duval et al. (2010, labeled RD10) and orange triangles the GRS structures discussed in Battisti & Heyer (2014, labeled BH14). The black, cyan, and orange lines plot the respective power-law fits to the individual data sets, and the red line corresponds to the power-law fit for the combined data. The original result from Larson (1981) is shown as a thick gray dot-dashed line. Fit parameters for each line are listed in Table 3.
Data Set aaParameter correlation coefficient .
BGPS 0.15 (0.04) -0.96
RD10 0.16 (0.04) -0.98
BH14 0.15 (0.08) -0.99
Combined 0.12 (0.02) -0.95
Larson (1981) 202
Table 3Power Law Fit Values for Larson-like Surface Density - Mass Relationship

The first two Larson relationships (the size-linewidth and mass-linewidth relationships) may be combined to form a mass-size relationship, sometimes called Larson’s 3rd Law.141414This name is also given to the number density-size relationship similarly derived from the first two Larson (1981) relations. Using the physical radius rather than the overall length of the molecular feature for consistency with the remainder of our analysis, we derive


as the form of the Larson mass-radius relationship. This power-law relationship has been reproduced by simulations of isothermal supersonic turbulence (Kritsuk et al., 2013).

While Equation (15) appears to have theoretical backing, it is important to note the effect of observational systematics on the derived relationship. Combining Equations (2) and (5) leads to a systematic relationship, which has been advanced as the source of Larson’s 3rd Law (e.g.,  Ballesteros-Paredes et al., 2012, and references therein). Owing to the similarity between Equation (15) and the systematic relationship, we investigate the deviation from for both BGPS sources and measurements of GMCs from the CO(1-0) data of the GRS.

To remove the effect of heliocentric distance, we compare the mean mass surface density () with the mass of each object. Figure 11 plots the 1,369 BGPS Distance Catalog objects with real as black circles, where the vertical axis is computed via . The wide scatter in points is likely due in part to the non-uniform geometry of objects (as discussed earlier). To evaluate the correspondence of the relationship for the BGPS objects (i.e., clump-scale structures) with larger structures, we also plot in Figure 11 results of two studies of CO clouds from the GRS. GMCs studied by Roman-Duval et al. (2010, hereafter RD10) are shown as cyan squares, and the clouds containing BGPS sources from Battisti & Heyer (2014, hereafter BH14) are depicted with orange triangles. While both of these studies used the GRS CO data, they used different source identification algorithms, resulting in the differences shown in Figure 11. For each data set, we fitted a linear function in logarithmic space that is robust against outliers151515Rather than a “least squares” fit, we utilized a “least absolute deviation” method to minimize the effects of outliers on the resulting fit (see Press et al., 2007, §15.7.3). to the data as , and the resulting power-law fit values are listed in Table 3. The individual fits for each data set are shown in the same color as the data points. Furthermore, we also fit a power law through the entire collection of points (red) to span a further order of magnitude in each dimension than any one data set alone.

Casting Equation (15) in the same form as the data of Figure 11 yields


indicating a nearly-flat mean mass surface density as a function of mass, and is shown as the dot-dashed gray line in Figure 11. The slope is so shallow as to cause only a factor of two difference in over some six orders of magnitude of mass.

Strikingly, the data from both continuum (BGPS) and CO (RD10, BH14) observations consistently show a trend of increasing with object mass, with the fits to the three individual sets returning nearly identical power-law indices. The similarity of these fit values indicates no significant difference between physical properties computed from dust continuum and molecular transition line observations of dense molecular cloud structures. These fit values significantly differ from the value in Equation (16). Although contrasting with the Larson slope, these fits still represent less than an order of magnitude change in over six orders of magnitude in mass.

Caution should be taken when drawing conclusions from this result. While the intrinsic scatter within each data set is large, none of the data sets could be consistent with the predicted Larson relationship () within the fit uncertainties. However, both the Larson relationship and the fits shown in Table 3 deviate only slightly from the systematic . The fact that the properties of dense molecular cloud structures from three separate studies all reveal nearly the same, slight departure from the systematic indicates a real, though slight, trend in mass surface density as a function of mass. In its simplest interpretation, this relation seems to imply that higher mass GMCs are also denser and therefore could possibly form stars more efficiently. This interpretation must be tempered, however, by the inability of current survey data-reduction techniques and cataloging algorithms to accurately represent large-scale low-surface-brightness structures in the interstellar medium.

5.2.2 The Continuum of Structure

The continuum of structure within molecular cloud complexes is a fundamental property whose connection with the intimate processes of star formation requires further study. Much of the literature discusses objects as clouds, clumps, or cores, and assigns specific properties and physical significance to each (e.g.,  McKee & Ostriker, 2007; Bergin & Tafalla, 2007; Kennicutt & Evans, 2012, and references therein). As seen in Figure 11, however, observed dense molecular cloud structures exhibit a continuum of properties across many orders of magnitude in scale; here we discuss the implications of this continuum.

The near uniformity of mean mass surface density of dense molecular cloud structures provides a crucial postulate for the construction of the Galactic distribution of dense molecular gas. Whether an artifact of observational systematics or a consequence of supersonic turbulence (e.g.,  Hopkins, 2012), the observed uniform relationship indicates that the fractal quality of the molecular interstellar medium161616The interstellar medium, and molecular clouds in particular, display fractal dimensions (e.g.,  Beech, 1987), wherein a self-similarity relationship exists over many orders of magnitude and no preferred scale exists between the limits of Galactic disk scale heights and runaway gravitational collapse of molecular cloud cores. is perhaps more fundamental and a stronger driver of molecular cloud complex evolution than the individual classes of objects (i.e., clouds, clumps, cores) working independently. There exist real physical divisions as the fractal nature of molecular clouds is driven by the the dominance of turbulence (Kritsuk et al., 2013); regimes dominated by other physical processes (i.e., gravitational collapse) obey other scaling relationships.171717Following a Jeans collapse scenario, , given uniform composition and temperature.

Aside from questions of the evolution of the ISM, what are the implications of these scaling relations on observational data? As a basic reference point, we start from the Larson result of Equation (15) and convert the physical properties to observable quantities via Equations (3) and (5). The resulting relationship between observed flux density, deconvolved angular radius, and heliocentric distance is given by


This equation implies that for dense molecular cloud structures, the observed flux density is nearly independent of heliocentric distance, and depends only on the angular radius on the sky, which is the same as saying the mass surface density is independent of mass (as seen in Figure 11). This relationship is plotted in Figure 3 (top) as green lines for  = 1 kpc (solid), 5 kpc (dotted), and 20 kpc (dashed), and matches the slope of data quite well. Therefore, in maps from (sub-)millimeter continuum surveys, the physical scales of detected objects are indistinguishable, meaning that a core looks like a clump looks like a cloud. It is only through robust distance estimation or application of the size-linewidth relationship that the scales of molecular cloud structures may be interpreted. This is an important point because the use of ancillary data is essentially the only means for deriving physically meaningful results from these surveys.

In EB15, it was noted that the surface brightness () distributions for the full BGPS catalog and various subsets (including the Distance Catalog) were largely independent of cataloged flux density. Indeed, by combining the surface brightness equation from EB15 (their Equation 3) with Equation (17), we obtain


which indicates that the surface brightness of dense molecular cloud structures, regardless of size scale and heliocentric distance, should be nearly identical. Indeed, Solomon et al. (1987) comment that for the standard size-linewidth relationship, virial equilibrium translates into constant average surface density, and hence constant average surface brightness. The observed surface brightness along any line of sight is then proportional to the number of structures, or total mass, along that line of sight. If one assumed an axisymmetric model and sampled a large fraction of the plane, the radial dependence (and normalization) of the molecular gas along the line of sight could be solved for.

5.3. Mass Distributions and Mass Functions

5.3.1 Aggregate Mass Distribution versus Subsets

The aggregate mass distribution shown in the left panels of Figures 6 and 7 represents an inhomogenous set of dense molecular cloud structures across the Galactic plane. The entire collection of sources appears to be fitted well by a lognormal function over more than two orders of magnitude, and displays one or two power laws at the high-mass end. Prompted by theoretical predictions for a power-law form of the so-called clump mass function (e.g.,  Donkov et al., 2012; Veltchev et al., 2013), we made a series of astrophysically motivated cuts (see Appendix A) in an attempt to isolate a more-homogenous population of sources.

While the classifications of molecular cloud versus clump versus core are intuitively appealing, the mass distributions of BGPS sources divided into core/clump and cloud groups (Table 2) return functional fits very similar to the aggregate set. Three possible explanations exist, likely working in concert: (1) the dividing lines between object classes are somewhat arbitrary, (2) uncertainties in the heliocentric distances mix objects across classification boundaries, and (3) the distinctions themselves are artificial in light of the continuum of structure in molecular cloud complexes. It is likely that a survey with improved angular resolution may be able to disentangle mass function parameters as a function of source classification or physical properties. Furthermore, use of multi-band Hi-GAL data will obviate the need to assume a dust temperature (or to marginalize over the dust temperature, as was done with the Monte Carlo trials), eliminating a source of uncertainty in derived properties.

The best power-law fits were for subsets primarily limited in heliocentric distance, and are characterized by nearly-symmetric [] joint confidence intervals (similar to the right panel of Figure 6). These subsets were the “nearby” heliocentric distance and “mixed protocluster” cuts, which exclude the most distant objects, with the latter additionally excluding nearby, low-mass objects. While both are strongly power-law the “nearby” subset has a steeper power-law index ( versus ), although the two are marginally consistent with each other.

The lack of a strong power-law fit for the “distant” sample or the two cloud-scale subsets (in and ) may be partially caused by the apparent steepening in the mass distribution at high-mass as discussed below. The subsets that do exhibit power-law structure, however, indicate that the presence of a true underlying power law (if present) is not obscured and converted to a lognormal by sample inhomogeneity (see §5.3.3).

5.3.2 Steepness of the High-Mass Power-Law Fit

One curious aspect of the mass function fits is the propensity for the high-mass end to display a fairly steep power-law tail (). While this power-law index is similar to the high-mass end of the Salpeter (1955) IMF, objects detected by the BGPS in this mass range are generally cloud-scale in physical radius. Since Galactic GMCs tend to follow a power-law mass function with (Blitz, 1993), the present result is somewhat contradictory.

A steep power-law function (i.e.) implies that most of the total mass is in low-mass objects and that high-mass objects are rarer than more shallow functions (i.e., there are fewer per mass bin; Equation 9). Could the apparent steepening be caused by “missing” high-mass sources in the present Distance Catalog? Simple simulations show that removal of sources from a shallow power-law distribution does have the effect of imitating a steeper power law, given that a progressively larger fraction of sources are removed with increasing mass.

Between the atmospheric subtraction in the BGPS data-reduction pipeline removing emission on large angular scales (; G13) and the tendency of the cataloging routine to break up large-scale emission into smaller chunks based on intermediate peaks in the flux-density distribution, the net effect is to remove the largest scale objects and convert single large structures into multiple smaller ones. What effect do these have on the derived mass distribution? Converting the mass-size scaling relationship of Equation (15) into observational quantities via Equation (5) yields