© ESO, 2012

1. Introduction

Measuring the structure and dynamics of the dense, molecular interstellar medium (ISM) is of crucial importance for our understanding of star formation from low- to high-mass stars and star clusters. There are a number of outstanding open questions in star formation that come down to the coupling between star formation and the fundamental physical processes acting in molecular clouds. What determines the star-formation rates and efficiencies in molecular clouds and why are those parameters generally low? What is the relative importance of gravity, turbulence, and magnetic fields during molecular cloud evolution and does this vary for low- and high-mass star formation? These questions can be approached by studying the physical characteristics of molecular clouds at different evolutionary stages. Of especial importance are studies of clouds that have not yet started to form stars or are in a very early stage of doing so. This is because such clouds still bear the imprints of the processes that initially formed them and are responsible for driving them toward star formation. For reviews of these topics, we refer to McKee & Ostriker (2007).

Our knowledge on the maternal structure of clouds potentially forming high-mass stars (M ≳ 8 M_⊙) and star-clusters is relatively scarce. This is, in part, because of the rapid destruction of the maternal cloud structure by stellar winds, outflows, ionization, and radiation pressure from the first high-mass stars to form in the cloud. Once the presence of high-mass stars in a cloud can be observationally confirmed, the maternal cloud structure has already been heavily affected by them. In part, the difficulties also originate from trivial observational restrictions: most sites that are thought to be favorable for future high-mass star and star cluster formation are located at distances of D ≳ 1.5 kpc, and they are concentrated in the Galactic plane. This makes attaining sensitive, high-spatial resolution data on their structure a very challenging task.

Only relatively recently, a class of dense molecular cloud clumps now commonly referred to as infrared dark clouds (IRDCs, hereafter) was discovered (Perault et al. 1996; Egan et al. 1998). These objects are currently regarded as the best candidates of future high-mass star formation (e.g., Carey et al. 2000; Garay et al. 2004; Beuther et al. 2005; Rathborne et al. 2005, 2006), although how commonly they might form massive stars remains still unclear (e.g., Kauffmann & Pillai 2010). In the most recent years, the study of IRDCs has seen the coming of large-area Galactic plane surveys that examine the Galactic-scale distribution and properties of IRDCs through comprehensive, statistical samples (e.g., Simon et al. 2006; Jackson et al. 2008; Peretto & Fuller 2010; Aguirre et al. 2011).

One persistent obstacle in defining the physical conditions in the maternal environment of IRDCs is our inability to measure, in a uniformly defined manner, the mass distribution over the wide dynamic range of (column-/mass surface-) densities present in IRDCs. Various observational techniques are commonly used in probing the mass distributions: CO emission line observations (e.g., Roman-Duval et al. 2010), thermal dust emission at (sub)mm wavelengths (e.g., Schuller et al. 2009; Molinari et al. 2010), dust extinction mapping using background stars in near-infrared (NIR, e.g., Kainulainen et al. 2011a), and dust extinction mapping in mid-infrared (MIR) using either point source data (e.g., Rygl et al. 2010) or surface brightness data (e.g., Butler & Tan 2009; Ragan et al. 2009; Vasyunina et al. 2009). All these mass-tracing techniques have their specific properties that restrict the parameter space (usually, column density range) within which they can provide accurate mass measurements. For example, CO line emission is known to be linearly correlated with column density only in a very narrow column density¹ range (A_V ≈ 3−6 mag) (e.g., Goodman et al. 2009). Dust emission measurements depend on the temperature of the dust grains that can vary by a factor of  ~2 in the clouds (e.g., Bontemps et al. 2010; Stutz et al. 2010; Arzoumanian et al. 2011) and of their emissivity index. Dust extinction measurements in the near-infrared are limited to relatively low column densities, typically to A_V ≈ 1−25 mag in (e.g., Lombardi & Alves 2001; Kainulainen et al. 2011a). For MIR surface brightness extinction mapping, higher column densities can be probed up to A_V ≈ 100 mag, but, as discussed below, systematic uncertainties limit the accuracy below A_V ≈ 10 mag. As a result, it is not straightforward to relate the mass distribution at small, dense, star-forming scales to the surrounding larger-scale, more diffuse environment.

In our previous works, we have developed techniques to probe the mass distributions of IRDCs by measuring the NIR (Kainulainen et al. 2011a) and MIR (Butler & Tan 2009, 2012, BT09 and BT12, hereafter) dust extinction through them. While the NIR technique has proven to be relatively sensitive, reaching column densities of A_V ≈ 25 mag at a resolution of  ~30″, the MIR data probe reliably column densities between A_V ≃ 10−100 mag at a resolution of  ~2″. Both of the techniques are essentially based on a model for dust optical properties in infrared which are relatively well-known, and the dynamical ranges probed by them clearly overlap. This provides a good basis for a conception of combining the two techniques with aim of achieving both extended dynamical range and high spatial resolution. This goal is further motivated by the fact that the MIR technique is known to have caveats affecting especially low-column densities. The technique effectively filters out the large-scale, diffuse (A_V ≲ 5−10 mag) extinction component, because of the estimation of the background intensity (for details, see Butler & Tan 2009, 2012, BT09 and BT12 hereafter, respectively). The details of the adopted filtering scheme can significantly affect the performance of the technique (e.g., BT09 Wilcock et al. 2012). However, the diffuse dust component that is filtered away is well recovered by the NIR data, and it would seem feasible to use the NIR data to “re-calibrate” the MIR data, thus alleviating the uncertainty of the MIR technique arising from the large-scale filtering needed to estimate the MIR background intensity.

We develop in this paper a new technique for probing the mass distributions of IRDCs that takes the advantage of both NIR point source data and MIR surface brightness data. The technique is based on data that are readily available for a large portion of the Galactic plane, i.e., NIR data from the UKIDSS/Galactic Plane Survey (Lucas et al. 2008) and MIR data from the Spitzer/GLIMPSE survey (Benjamin et al. 2003). We will demonstrate the use of the technique by deriving column density data for ten IRDCs and their surroundings.

We also demonstrate in this paper the scientific potential of high-dynamic-range column density data in examining the probability density functions (PDFs) of column densities in the clouds. The column density PDFs can be used to study the roles of gravity and turbulence in molecular clouds (e.g., Froebrich et al. 2007; Goodman et al. 2009; Kainulainen et al. 2009; Froebrich & Rowles 2010), generally building on the fact that simulations of supersonic turbulent media predict a log-normal shape for the density PDF (e.g., Vazquez-Semadeni 1994; Padoan et al. 1997a; Ballesteros-Paredes et al. 2011; Kritsuk et al. 2011) (1)where x = ρ/ρ₀ is the mean-normalized volume density and μ and σ_lnx are the mean and standard deviation in logarithmic units. Deviations from the log-normal shape are predicted for strongly gravitating systems, inducing higher relative fraction of high-column density material (e.g., Klessen 2000; Federrath et al. 2008a). In K09 we presented an analysis of the PDFs of a sample of nearby molecular clouds and concluded that the PDFs of star-forming clouds are consistent with a log-normal (or normal) shaped function at lower column densities (A_V ≲ 2−5 mag) and a power-law like extension at higher columns (A_V ≳ 2−5 mag). In contrast, the PDFs of non-star-forming clouds were usually very well described with log-normals. This led us to suggest that quiescent clouds are primarily shaped by turbulence, while the effect of self-gravity is evident in the PDFs of star-forming clouds. However, in Kainulainen et al. (2011b) we analyzed the data further and suggested that the PDF tail could be related to clumps that are not self-gravitating but rather supported significantly by external pressure from the surrounding, gravitating cloud material.

In this paper, we also use the derived column density data together with CO line emission data to examine the relationship between density fluctuations, measured by the column density dispersion σ_N/⟨N⟩, and turbulence energy, measured by the sonic Mach number ℳs, in the clouds. These variables are expected to be correlated in turbulent media and a correlation coefficient between them is assumed by numerous analytic models of star formation (e.g., Padoan et al. 1997a; Padoan & Nordlund 2002; Krumholz & McKee 2005; Hennebelle & Chabrier 2008; Elmegreen 2011). We measure the relationship between ℳ_s and σ_N/⟨N⟩ in the IRDCs and also in a larger sample that is supplemented by including a set of nearby molecular clouds we have previously derived column density data for (Kainulainen et al. 2009, hereafter K09). Our results demonstrate that high-dynamic-range data can, and are in fact required to, probe the relation.

2. Data

2.1. The IRDC sample

We use as a sample the ten IRDCs selected by BT09. The coordinates of the clouds are listed in Table 1. The clouds were originally chosen from the sample of Rathborne et al. (2006) as massive, relatively nearby clouds with high contrast against the Galactic background, being thereby well-suited for the 8 μm extinction mapping. The clouds are distributed between the Galactic longitudes of l =  [18°,53°]  and kinematic distances of D ≈ 1.8−5.7 kpc (derived by Simon et al. 2006). The most massive cloud (D) has a MIR-derived mass of M_MIR = 5.34  ×  10⁴ M_⊙. The masses are listed in Table 1. The masses were derived by BT12 and they comprise the material inside the ellipsoids defined by Simon et al. (2006), which is likely to be  ~10% of the total mass of the clouds estimated over larger scales using NIR extinction mapping (Kainulainen et al. 2011a). Thus, the most massive clouds in the sample are somewhat more massive than the nearby Orion A molecular cloud (M ≈ 10⁵ M_⊙, e.g., K09). All clouds contain mm-cores at 11″ resolution (Rathborne et al. 2006), most of which also appear as single cores at the 2″ resolution of the MIR dust extinction data (BT12). The masses of the cores are typically some tens of solar masses, but depending on how much of the surrounding envelope is included, the masses can also be higher by a factor of a few (BT12). While it is currently not well-constrained how massive a core is needed for high-mass star formation, the cores in principle contain enough material to form such stars.

2.2. NIR dust extinction data

We use in this work NIR dust extinction data as a tracer of low-extinction material (A_V < 10 mag) in the IRDCs. The data were derived using the NIR color-excess mapping scheme presented by Kainulainen et al. (2011a), which is based on the so-called NICER technique (Lombardi & Alves 2001), but refined to deal with the diffuse extinction component in the Galactic plane and with the large number of foreground sources towards IRDCs. For the detailed description of the technique, we refer to Kainulainen et al. (2011a) and Lombardi & Alves (2001). Kainulainen et al. (2011a, see their Appendix A) presented the NIR extinction maps for the clouds A, B, G, H, and I. The maps for the remaining clouds were computed for this work.

The NIR extinction mapping technique was used in conjunction with the JHK_S band data from the UKIDSS/Galactic Plane Survey/Data Release 6 Plus. The UKIDSS project is defined in Lawrence et al. (2007). UKIDSS uses the UKIRT Wide Field Camera (WFCAM; Casali et al. 2007) The photometric system is described in Hewett et al. (2006), and the calibration is described in Hodgkin et al. (2009). The science archive is described in Hambly et al. (2008). The derived NIR extinction maps have the sensitivity (statistical uncertainty) of σ(A_V) ≈ 0.4−0.7 mag, depending somewhat on the Galactic coordinates and the distance of the clouds. The zero-point of the data is rather uncertain because of the uncertainties in the foreground source removal and in estimating the extinction because of the diffuse dust in the Galactic plane; Kainulainen et al. (2011a) give the value of σ₀(A_V) ≈ 1 mag for the zero-point uncertainty. The spatial resolution of the data is FWHM = 30′′ and at that resolution they typically span extinction values between A_V ≈ 1−30 mag. We note that in derivation of the NIR extinction maps a wavelength dependency of optical depth, i.e., the reddening law, is assumed. Kainulainen et al. (2011a) used the reddening law as derived by Cardelli et al. (1989)(2)It should be emphasized that the NIR technique determines the extinction in NIR, e.g., A_K. In order to estimate column densities, we convert the NIR extinction to visual extinction using the coefficient from Cardelli et al. (1989)(3)We note that we combine our results later in this paper with K09 who used the same relative extinction law. Then, the calculated visual extinction can be converted to gas column density (Savage et al. 1977; Bohlin et al. 1978; Vuong et al. 2003) (4)Column density can further be expressed in terms of the total mass surface density by adding 40% (implying the H:He number ratio of 10:1) of helium (5)We can also express the empirical conversion from extinction to mass surface density (i.e., Eqs. (3)–(5)) with the help of an “effective” opacity, , through the definition (6)Inputing Eqs. (3)–(5) to this results in the effective opacity (7)Corresponding effective opacities for H and J bands result directly from Eq. (2)

2.3. MIR dust extinction data

We use the MIR dust extinction data derived using the Spitzer/GLIMPSE survey 8 μm images by BT12 (see also BT09) as another column-density tracer for the IRDC complexes. The MIR mapping technique is based on measuring the decrease of surface brightness towards the IRDCs, which can be interpreted in terms of 8 μm optical depth. For the details of the mapping technique we refer to BT12. The derived optical depths are used to estimate mass surface densities, Σ₈, with an equation analogous to Eq. (6) and adopting a dust opacity at 8 μm (10)BT09 and BT12 chose the value κ₈ = 7.5 cm² g^-1, which was taken as an approximate mean of five prospective dust models from Ossenkopf & Henning (1994) and Weingartner & Draine (2001), but closest to the Ossenkopf & Henning (1994) moderately coagulated thin ice mantle dust model (see BT09 for details). BT09 also used weighting by the MIR background spectrum and the filter response function in calculating the opacity values. Figure 1 illustrates the opacities κ₈ and  together with the opacity laws from Weingartner & Draine (2001) and Ossenkopf & Henning (1994). The figure shows that both in NIR and MIR the adopted opacities represent an average of the shown opacity laws. It is not possible to estimate the accuracy of the adopted 8 μm opacity, but the standard deviation of the opacities of the five models is 22%. Also, it is not straightforward to estimate the uncertainty due to variations in the MIR background spectrum, but if we re-calculate the mean opacity using the case of flat background spectrum, a value κ₈ = 9.0 cm² g^-1 follows, indicating 20% larger value than the mean with the chosen background spectrum. While the flat background spectrum is a rather extreme case, it seems reasonable to use the total uncertainty of  ~30% for the absolute value of the MIR opacity. However, one should note that unless average dust grain properties significantly alter from cloud-to-cloud, the relative opacity (or extinction/mass surface density) measurements between clouds is significantly more robust than this value.

The Spitzer 8 μm images have the resolution of FWHM = 2″ ( pixel size) and that is also the resolution of the MIR extinction data. At this resolution, the technique probes the mass surface density values up to Σ ≲ 0.4−0.5 g cm^-2, or equivalently, A_V ≲ 92−115 mag depending on the cloud. These upper limits, or “saturation” values, are estimated empirically as described in BT12 and listed in Table 1. The sensitivity (the absolute uncertainty) of the technique is not straightforward to determine, because the main sources of uncertainty are the variations of the background and foreground surface brightness. BT09 discusses that the technique becomes unreliable at Σ ≲ 0.013 g cm^-2, or A_V ≲ 3 mag. However, at column densities lower than Σ ≲ 0.04 g cm^-2, or A_V ≲ 9 mag, the MIR extinction maps show systematic variations that are not present in the NIR data, likely resulting from variations in the MIR background, and hence, from the zero-point uncertainty of the MIR data. The MIR mapping technique is also affected by bright point sources and extended emission (nebulae) in the 8 μm data. These were masked out during the mapping procedure as detailed in BT09, i.e., we do not estimate extinction in the regions that suffer from significant nebulosity. We also note that the GLIMPSE data contains some stripe-like artifacts resulting from non-uniform calibration.

The properties above reflect the complimentary of the NIR and MIR dust extinction techniques. While the NIR data are more sensitive and less affected by systematic uncertainties at A_V ≲ 10 mag, the MIR data can probe significantly higher column densities. There is also an overlap range between the dynamic ranges of the techniques at about A_V ≈ 5−15 mag.

2.4. CO line emission data

We use the publicly available ¹³CO (J = 1–0) line emission data from the FCRAO/Galactic Ring Survey (Jackson et al. 2006, hereafter GRS) to examine the velocity structure of the IRDC complexes. The GRS data have a spatial resolution of FWHM = 46″, spectral resolution of 0.2 km s^-1, and typical sensitivity of σ(T_A,rms) = 0.4 K.

Since the IRDCs are located in the Galactic plane, the lines of sight toward them can show several emission components. We first determined which component corresponds to the IRDC by comparing their morphologies to the extinction data. This was done subjectively by eye, but we note that in almost all cases the choice was trivial; the IRDCs were originally chosen to be in somewhat isolated positions. In nine out of ten IRDCs the strongest ¹³CO emission peak matched well the morphology of the IRDC as seen by dust extinction. Toward the Cloud E the velocity structure of ¹³CO was more complicated (see Fig. A.1). Because the dust extinction techniques can probe the total column density up to D ≲ 8 kpc, it appears likely that the extinction data for cloud E is affected by several, physically unrelated clouds. Therefore, we excluded the cloud E from the analyses performed in Sect. 4.

We used the ¹³CO line emission to estimate the velocity dispersions in the clouds. The estimation was performed inside two regions for each cloud: an ellipsoid as defined by Simon et al. (2006), and a box-like region that (often, but not always) contains a closed A_V = 7 mag contour. We describe these regions more closely in Sect. 3.2 and show the average spectra in Appendix A. For both regions, we estimated the one-dimensional velocity dispersion, σ_v, using the following two approaches. The first approach was based on fitting a single Gaussian to the mean spectrum over the velocity interval of the component. The velocity interval was defined using an average spectrum over the region. We calculated the first derivative of the average spectra, and determined the first zero-points of the derivative from both left and right side of the velocity of the maximum intensity. These zero points defined the interval limits. We also required the interval limit to be at most 15 km s^-1 offset from velocity of the peak. The intervals are indicated in Figs. A.1 and A.2 with vertical dotted lines. Then we fitted a Gaussian to the spectrum within this interval. To get a handle on how much the velocity dispersion depends on the chosen velocity interval, we also fitted a Gaussian to the interval v =  [v_peak − 1.5σ_v,v_peak + 1.5σ_v] , where σ_v is the dispersion from the first Gaussian fit. These two velocity dispersions (and the corresponding fits) are also shown in Figs. A.1 and A.2 and we give in Table 1 their mean value. The mean value is also used in all further analyses.

The second approach for calculating the velocity dispersions was calculating the standard deviation of the spectrum in the fixed velocity interval whose definition is described above. As the mean in this calculation, we used the intensity-weighted mean, which represents the center-of-mass in the case where intensity is linearly proportional to the column density along the line of sight. The linear proportionality assumes that the ¹³CO emission is optically thin and does not suffer from depletion (see Hernandez & Tan 2011; Hernandez et al. 2011, for a study of these effects in the clouds F and H). Only channels with intensity higher than three times the root-mean-square noise level were included in the calculation. The resulting velocity dispersions are given in Table 1.

In addition to the GRS data, we use in Sect. 5.2.2¹²CO data from the COMPLETE survey (Ridge et al. 2006) and the survey of Dame et al. (2001) to examine the velocity dispersions of nearby molecular clouds. The COMPLETE data was available for Ophiuchus and Perseus. The Dame et al. (2001) data was adopted for the clouds in the sample of K09 (except for the Lupus clouds and LDN134 that are not covered by Dame et al.). The CO data for these nearby clouds do not suffer from the line-of-sight confusion (they are outside the Galactic plane) and the linewidths toward them were computed from a Gaussian fit to the mean spectrum. In addition, we adopt the mean linewidth for Taurus from Goldsmith et al. (2008), i.e., σ(¹³CO)  = 1.0 km s^-1.

3. Methods

3.1. Combination of NIR and MIR data

We combine the NIR and MIR dust extinction maps with the aim of achieving high-dynamic-range extinction data. We first examined the correlation between the data sets by convolving the MIR data (FWHM = 2″) to the significantly lower resolution of the NIR data (FWHM = 30″) and making pixel-to-pixel comparisons. As an example, Fig. 2 shows the comparison for the cloud B. Despite the rather high scatter, the data sets correlate at column densities below A_V ≲ 10−15 mag. In this range, the column densities measured by the MIR data are systematically about  ~4 mag lower than the column densities measured from NIR data. The offset is present in all clouds and varies between A_V ≈ 2−8 mag. The offset is most likely a result of the spatial filtering of the MIR mapping technique which effectively subtracts the diffuse cloud component from the MIR maps. At column densities higher than A_V ≳ 1−15 mag the NIR data begins to underestimate the column density, which results from the fact that the NIR data does not efficiently sample high (A_V ≳ 30 mag) extinction lines-of-sight (none, or very few, background stars are visible through higher-extinction clumps/cores).

Since the NIR and MIR data showed correlation in the column density range A_V ≲ 10−15 mag, we conclude that it is reasonable to assume that the difference between the MIR and NIR extinction data below this range is a good representation of the background component that the MIR technique filters out. With this in mind, we combined the NIR and MIR data with the following procedure. First, we re-binned and smoothed the MIR data to the resolution of NIR data (FWHM = 30″, pixel size 15″) and subtracted it from the NIR image to generate “a correction image” (11)When built in this way, the correction image, , approximates the large-scale component filtered out by the MIR technique at spatial scales equal or larger than the background filter size used in the MIR technique (in the case of BT12, a few arcmin). Thus, the correction image clearly resolves the filtering scale of the MIR technique. Since the NIR data underestimates the column density above A_V ≳ 10−15 mag, we blanked the pixels in the correction image in which the extinction is higher than A_V ≥ 10 mag. Values for these pixels were estimated as a mean within an aperture of 75″ (5 pixels of the NIR data). A Gaussian weighting function with FWHM = 30″ was used to weight the values when calculating the mean. Then, the correction image was re-binned to the pixel scale of the MIR data () and added to the original MIR column density map (12)This resulted in column density maps, , that take advantage of the signals both in the NIR and MIR data. The pixels that lack either NIR or MIR data (e.g., because of the incomplete coverage of the NIR data, bright MIR point sources or extended MIR nebulosity) were flagged to have no data in the maps. As a result, extinction is not determined for every pixel but only for regions in which we consider both the NIR and MIR techniques to be unaffected by any artifacts. We describe the results further in Sect. 4.1. In the remainder of this paper, when we refer to dust extinction/column density/mass surface density maps, or data, we refer to these combined maps.

We finally note that the combination procedure is not sensitive to the absolute values of the adopted dust opacities (see Fig. 1), only to the ratio of dust opacities in NIR and MIR. Only the scaling of the observed extinctions to the units of mass surface density requires knowledge of the absolute dust opacity. In our approach the NIR-MIR dust opacity ratio equals to (13)As shown in Fig. 1, this ratio represents an average behavior of the displayed models. We quote as uncertainty of the ratio the standard deviation of the models used in Fig. 1, which is 30% for the 8 μm (see Sect. 2.3) band and 15% for the K band. Combining these values quadratically results in uncertainty of 36%. However, we again point out that this measure of uncertainty reflects a selection of dust models. If the dust grain properties are relatively similar from cloud-to-cloud, the ratio in considerably less uncertain. Currently, there is no detailed information available on this ratio in the particular environment of IRDCs, and hence, aligning the relative reddening law more accurately based on actual measurements is not possible.

3.2. Definition of the cloud area

It is not straightforward to define an object such as “a cloud” in the relatively crowded fields toward the Galactic plane in which the IRDCs reside. In Sects. 5.2.1, 5.2.2 we examine statistics derived from the column density PDFs of the IRDCs. Despite the seemingly straightforward definition of the PDF, it is not trivial to measure its shape from observational data in a robust way. An ideal target would be a totally isolated cloud in which the lowest column density included in the PDF forms a closed contour in the column density map. Such a closed contour can be regarded as the completeness limit of the PDF; the shape of the PDF below the level of the lowest closed contour can clearly suffer from incompleteness. For example, in K09 we used as a definition of a cloud the level A_V = 1 mag, which in most cases was below the mode of the PDF. In the case of the IRDCs that are located in the Galactic plane the confusion from physically unrelated clouds that are nearby in projection makes it more difficult to define the region which should be included in the PDF.

With this caveat in mind, we adopted two definitions of the cloud area for our analysis. The first definition is to use the ellipsoids defined by Simon et al. (2006). We chose from the catalog the ellipsoids with largest area within the mapped area, which in practice roughly encircles the high-column density (A_V ≳ 10−15 mag) regions in the IRDCs. The ellipses are shown in the figures of Appendix B. The second definition of the cloud area was constructed with the aim in being as complete as possible down to relatively low column densities. In the second definition, we defined a rectangular area that surrounds the Simon et al. (2006) ellipsoid, and makes an effort to completely surround a closed A_V = 7 mag contour. However because of blending of neighboring clouds in two-dimensional data, this contour is not necessarily closed in all clouds. In cases where the blending was severe, we used the velocity structure as traced by the ¹³CO data to help disentangling physically unrelated components. The final rectangular regions chosen are shown in the maps in Appendix B. We conclude that the data within these areas should be close to complete above A_V > 7 mag. However, we also note that in some, more isolated clouds slightly off the Galactic plane, such as clouds C and J, the lowest closed contour can be as low as A_V ≈ 2 mag. For brevity, we refer to this second cloud definition as a “A_V = 7 mag box” in the remainder of the paper.

4. Results

4.1. High-dynamic range column density maps of IRDCs

We used the method described in Sect. 3.1 to derive column density maps for the ten IRDC complexes of our sample. We present the maps in Appendix B (online only). The pixels in the maps for which there was either no MIR or NIR data available are blanked to black. These blanked areas can result from strong MIR point sources and nebulosity, or missing UKIDSS data (due to the incomplete coverage of the survey).

Figure 3 shows as an example the column density map of the IRDC H and, for comparison, the separate NIR- and MIR-derived column density maps of the complex (Kainulainen et al. 2011a, BT12, respectively) in the same color-scale. The combined map clearly reflects the characteristics of both data sets: it recovers the large-scale structures in the cloud with better sensitivity than the MIR data alone, and it retains the high angular resolution of high-column density structures of the original MIR-derived maps. The combined maps probe column densities from the sensitivity limit (σ(A_V) ≈  1 mag) to the saturation limit of the MIR data (A_V,sat ≃ 100 mag).

We list in Table 1 the new masses of the IRDCs calculated from the combined column density data. The masses from the combined data are on average 1.4 times larger than the masses calculated from the MIR data alone. We also show in Table 1 also the virial parameters of the clouds, defined as (Bertoldi & McKee 1992) (14)In this σ_v is the velocity dispersion in km s^-1 for which we use the dispersions derived using Gaussian fits and ellipsoidal cloud definition (i.e., σ_G,ell in Table 1). R is the radius of the ellipsoid in pc for which we use the expression , and M is the mass in solar masses. G = 1/232 is the gravitational constant in corresponding units. Five out of ten clouds have their virial parameters within the factor of 2 of unity. The cloud G has the largest virial parameter, α_vir = 4.7. None of the clouds have their virial parameters more than factor of 2 below unity. Modulo the effects of surface pressure terms, internal density structure, and cloud elongation, these results suggest that IRDCs (at least those typical of this sample) tend to be self-gravitating, gravitationally bound objects (see also Hernandez et al. 2012).

4.2. Probability distributions of column densities

We used the derived maps to examine the column density PDFs in the IRDCs. The PDFs describe the probability of observing a column density value A_V in the clouds. Because of the completeness issues in defining the clouds with the Simon et al. (2006) ellipsoids, we only examine the PDFs within the A_V = 7 mag boxes (see Sect. 3.2), and analyze them only above A_V ≥ 7 mag. We also note that we do not use all the dynamic range present in the column density maps which extend typically to A_V = 140−150 mag. As discussed by BT12 the saturation of the mid-infrared technique starts to affect the column density values clearly earlier. The saturation limits, Σ_sat, below which BT12 conclude the maps to be relatively unaffected by saturation, are given in Table 1. The error in column density due to saturation is for most clou  ≲ 10% when Σ ≲ 0.8Σ_sat. We only include the column densities up to the saturation limits Σ_sat to the PDFs.

Figure 4 shows the PDFs of eight IRDCs of our sample. The PDFs were normalized for the plot by first fitting a line to the range lnA_V = 2−3 (A_V = 7.4−20 mag) and then setting a boundary condition p(lnA_V = 2.0) = 1. In general, the PDFs of all IRDCs are very similar to each other below A_V ≲ 30−40 mag. Above that, they differ from each other, with the PDFs of some clouds decreasing more rapidly than others (i.e., they contain a lower fraction of high column densities). For reference, we show in Fig. 4 a log-normal function with μ = 1.2, σ_lnN = 0.9 (see Eq. (1)). The PDFs of all but two IRDCs (C and D) are within a factor of three of this curve throughout their dynamic range. The individual clouds show varying behaviors, e.g., the cloud C has almost a power-law like PDF until A_V ≈ 70 mag, after which it drops abruptly. It is not easy to assess whether the PDFs are better fit with a power-law between lnA_V ≈ 2−3.5 and a drop afterwards, or by a log-normal. This is because of the similarity of the power-law and a wide log-normal function whose peak (mode) is not covered. We analyze the properties of the PDFs further in Sects. 5.2.1, 5.2.2.

We also note as a caveat that the pixels affected by bright MIR sources (YSOs) and nebulae are blanked in the extinction data and can affect the PDFs (see Sect. 3.1). In principle, it is possible that especially the high-extinction portions of the PDFs are affected by the missing data, because high-column density regions could presumably harbor more active star formation than low-column density regions. There are bright MIR point sources, likely local to the cloud, present inside the Simon et al. (2006) ellipsoid areas of the clouds A, C, D, E, and H. However, the area covered by these blanked regions is only a fraction of, e.g., the area above A_V > 40 mag, which gives some suggestion that their effect should not be very large. Similarly, while almost every cloud has some MIR-bright nebulosity in their area, only in the case of Cloud I is the fraction of blanked pixels higher than a few percent (for the cloud I about 10% of the “7 mag box” is affected, see Sect. 3.2). Therefore, we note that some caution is in place when regarding the PDF of the Cloud I, and we also exclude it from the analyses of the PDFs performed later in this paper.

4.3. Fractional amount of high-column density material

We examined the fractional amount of high-column density material in the IRDCs by calculating the cumulative distribution functions (CDFs) of their column densities. Here we define the CDF as the fraction of mass above the reference value (15)where M(>A_V) is the mass above the threshold A_V and M_TOT is the total mass above A_V = 7 mag. The CDFs are constructed from the same pixels as the PDFs in Sect. 4.2, i.e., for the A_V = 7 mag box. The resulting CDFs are shown in Fig. 5 (left panel). The CDFs show a relatively linear shape in the log-linear scale, i.e., their shapes are consistent with exponential functions. We find that on average, the CDFs of the clouds follow the form (16)This relation was obtained through a fit to the mean CDF of all IRDCs. However, we note that the CDF, in general, depends on the dynamic range and resolution covered by the observations. Remembering this, we compared the CDFs derived for IRDCs with those of nearby molecular clouds from K09. To facilitate the comparison, the dynamic ranges of the IRDC data and the data from K09 were truncated to the range covered by both data sets. The lower limit for this range is set by the IRDC data, which is complete only above A_V ≳ 7 mag. The upper limit is set by the data of K09 that probe only values up to about A_V ≲ 30 mag. Therefore, we constructed the CDFs for this comparison using only the range A_V = 7−30 mag.

Figure 5 (right panel) shows the resulting CDFs for the IRDCs and three nearby giant molecular clouds (Orion A, Orion B, and the California Cloud). These nearby clouds have similar masses (M ≈ 10⁵ M_⊙), but their star-forming activities differ significantly. Orion A is the most active cloud with 2862 member sources, Orion B has 635 members, and the California Cloud only 279 members (adopting the member counts from Lada et al. 2009). Fits of Eq. (16) to the CDFs of these clouds result in slopes −0.18, −0.22, and −0.49 for the Orion A, Orion B, and California Cloud, respectively. The average CDF of the IRDCs has an exponent of −0.12, which is even flatter than that of Orion A. The main source of uncertainty in the CDF of the IRDCs originates from the chosen absolute value of κ₈ (see Sect. 3.1). To first-order, the uncertainty in κ₈ reflects as a multiplier to the CDF: for example, underestimation of κ₈ by a factor of 2 results in over-estimation of column densities by a factor of roughly two, and hence, the exponent of the CDF will be under-estimated by the same factor. In Sect. 2.3 we discussed that the absolute uncertainty of κ₈ is  ~30%. This translates to uncertainty of  ~30% for the exponent of the CDF. Hence, the exponent −0.12 we derive from the IRDCs is roughly 2-σ higher than that of Orion A. We conclude that this implies that the fraction of high-column density material in IRDCs is at least comparable, and possibly higher than in nearby active star-forming clouds. This result is similar to what we found out in Kainulainen et al. (2011a) using only NIR data. With the data presented in this paper we quantify the form of the CDF above A_V > 7 mag and establish that it on average has an exponential shape that continues up to A_V ≲ 100 mag.

5. Discussion

5.1. Column density probability distributions of IRDCs

How do the PDFs we derive for IRDCs in this paper compare with earlier findings of the column density PDFs in molecular clouds? We show in Fig. 4 the PDFs of the IRDCs and that of the Ophiuchus cloud from K09, which is representative of most star-forming clouds. The PDF of Ophiuchus is very similar to those of the IRDCs, and consistent with either a power-law or a log-normal (μ = 1.2, σ_lnN = 0.9). Thus, the PDFs of IRDCs seem to resemble the PDFs of nearby star-forming clouds, although they are slightly shallower, i.e., they contain relatively more high-column density material.

We note that in K09 we refer to the high-column density tails (at A_V ≳ 2−5 mag) of the PDFs as “power-law” tails. We elucidate the nomenclature and the dynamic ranges of the K09 data and this work in Fig. 6. The power-law-like behavior of the PDF of Ophiuchus is obvious in Fig. 4. However, the PDF of Ophiuchus is, as well as the PDFs of the IRDCs are, consistent with the log-normal function shown in the figure. This log-normal function is significantly wider than those used by K09 to describe the PDFs of low-column density (A_V ≲ 2−5 mag) material. Thus, the data we present in this paper suggests that also the high column densities in molecular clouds can be described by a log-normal PDF, however, one clearly wider than the log-normal that we used in K09 to describe low column densities. It remains to be explored in future studies whether the column density PDFs of IRDCs are well-described by log-normals down to low column densities (A_V ≈ 1 mag), or whether they show a similar change in shape at around A_V = 2−5 mag as the nearby star-forming clouds do.

In K09 we showed that the non-star-forming clouds do not show similar high-column density tails as the star-forming clouds do (or, considering the discussion above, their PDFs are clearly narrower than those of the star-forming clouds). This could be interpreted as gravity having a significant effect in shaping the cloud structure as has been established in numerical studies (e.g., Klessen 2000; Ballesteros-Paredes et al. 2011), but also possibly as confinement of non-self-gravitating clumps by the pressure imposed to them by their surrounding medium (Kainulainen et al. 2011b). What is common in both of these frameworks is that the shape of the PDFs, and especially the presence of a wider, log-normal (or power-law-like) component, is intimately linked to the onset of star formation. Further, Rowles & Froebrich (2011) analyzed the PDF shapes and embedded populations of nearby clouds and suggested that the PDF shape is also linked to the star-formation mode (isolated/dispersed vs. clustered) of the clouds. The PDFs of the IRDCs resemble more closely the PDFs of the star-forming clouds than non-star-forming clouds (cf., Figs. 4, 5) in a sense that the fraction of high-column densities is high (regardless whether the PDF shape is log-normal or power-law). This alone suggests that substantial (low-mass) populations may be present in the IRDCs of our sample.

We established in Sect. 4.3 that the fractional amount of high-column density material in IRDCs in the range A_V = 7−30 mag at least comparable to, and possibly even higher than in active nearby clouds. While there is a well-known correlation between the gas surface density and star-forming rate (SFR; Schmidt 1959; Kennicutt 1998), it was recently proposed by Lada et al. (2012) that there exists a more fundamental scaling between the surface density of dense gas and SFR. In particular, Lada et al. (2012) suggested that the surface density of SFR depends linearly on the mass surface density of H₂(17)where f_DG is the fraction of dense gas in the clouds, defined as (18)Here, A_K is the extinction in the K band (A_K = 1 mag is equivalent to A_V ≈ 8.8 mag, Eq. (2)). The fraction f_DG was measured and analyzed in a number of nearby clouds by Lada et al. (2010), who derived the mean value ⟨f_DG⟩ = 0.1 ± 0.06, with maximum fractions of f_DG ≈ 0.2 in the Orion A and Corona Australis clouds. Because the PDFs of the IRDCs we present in this paper are complete only above A_V > 7 mag, we cannot directly calculate the values of f_DG from our data. However, we can present an estimate for the f_DG in the following manner. We re-examined the CDF data presented by K09. We fitted exponential functions to the high-column density slopes of the CDFs of all star-forming clouds (13 clouds) and extrapolated the clouds’ masses above A_K > 0.1 mag using the fitted functions. This test showed that the cloud masses were underestimated on average by a factor of  using this procedure. The scatter in the factor had the 3-σ interval of c =  [0.7,5.2] .

Since the cloud masses are systematically underestimated, we applied the factor ⟨c⟩ = 1.9 to the masses of IRDCs derived by extrapolating their exponential CDFs to lower column densities. This resulted to the average dense gas mass fraction of f_DG ≈ 0.2 of our IRDC sample, which is about two times larger than the average fraction derived for nearby clouds by Lada et al. (2010). Given the uncertainty of our analysis, we find that this implies that the f_DG in IRDCs is of similar level, and possibly even higher, than in nearby clouds. This, in turn, suggests that the SFR of the IRDCs should be at least comparable to equally massive nearby molecular clouds on average. We recognize that this prediction remains suggestive, because we could not directly calculate the dense gas mass fraction, but had to rely on a correction for the missing mass below A_V ≲ 7 mag.

5.2. Relationship between the sonic Mach number and column density variance in molecular clouds

We demonstrate in this section the potential of high-dynamic-range column density data to examine the relationship between the column density dispersion and the sonic Mach number in molecular clouds. This relation is interesting, because it is directly linked to the fundamental relation between density fluctuations and turbulent energy in supersonically turbulent media: schematically, higher turbulent energy (higher Mach number) is expected to induce stronger density fluctuations, and hence, wider density PDF. In the following, we describe the formalism of this connection in more detail.

Numerical simulations predict a log-normal form for the PDF of volume densities in non-gravitating, isothermal media (e.g. Vazquez-Semadeni 1994; Padoan et al. 1997a). The dispersion of the logarithmic standard deviation σ_lnx (see Eq. (1)) is linked to the turbulent energy, i.e., the (non-thermal) sonic Mach number, ℳ_s, via the equation (e.g., Nordlund & Padoan 1999; Ostriker et al. 1999) (19)where b is a proportionality constant. This equation translates to the dispersion of non-logarithmic, mean-normalized volume densities (20)Thus, the proportionality constant b concretely describes how strong is the coupling between the turbulent energy and density fluctuations. Exploring the value of b in numerical simulations (e.g., Padoan et al. 1997a; Kritsuk et al. 2007; Federrath et al. 2008b; Lemaster & Stone 2008; Price et al. 2011) and analytically (e.g., Molina et al. 2012) has received ample attention recently, owing to its intimate connection to the analytical formulations of the dense core formation rate and the possible origin of the IMF (e.g., Padoan & Nordlund 2002; Krumholz & McKee 2005; Hennebelle & Chabrier 2008; Elmegreen 2011).

While observations of column density only probe the 2D density statistics, it is established that the log-normal nature of the volume density PDF remains in 3D to 2D transformation (e.g., Vázquez-Semadeni & García 2001; Fischera & Dopita 2004; Federrath et al. 2010). Further, Fischera & Dopita (2004) and Brunt et al. (2010) examined the relation between the 3D and 2D density variances and presented formalisms to connect the observed column density variance to the underlying volume density variance. In particular, the method of Brunt et al. (2010) is based on the possibility to estimate the 3D-2D variance ratio (21)using solely observations of 2D density field. In the equation σ_N/⟨N⟩ is the standard deviation of mean-normalized column densities. We note that Brunt et al. (2010) examined the values of R in a set of isothermal simulations of supersonic turbulence and found that in all their models R was between R =  [0.03,0.15] . The results of Fischera & Dopita (2004) indicate a similar value of R = 0.13 (for the case in which the thickness of the turbulent medium equals to its maximum turbulent length-scale). If we now combine Eqs. (21) with (20), we have an expression for b which depends only on 2D statistics (22)In this equation, a₁ is the ratio between the sonic Mach number and column density variance, i.e., the slope potentially measurable (however, see the discussion below) through a linear fit to the  [ℳ_s,σ_N/⟨N⟩]  data (23)It remains to be explored in future work whether the quality of our combined NIR+MIR data is suitable to directly estimate the variance ratio R, and from therein, the 3D variance in IRDCs. Nevertheless even in the absence of the transformation, relationship between the 2D density variance with the sonic Mach number is a proxy of the relationship in 3D. Supporting this approach, recently Burkhart & Lazarian (2012) showed that in MHD simulations the ℳ_s − σ_lnx relation is preserved in transformation from 3D to 2D variables, with the details of the correlation depending somewhat on the simulation parameters.

One should emphasize that interpreting observational (ℳ_s,σ_N/⟨N⟩) data in terms of the formalism presented above assumes that the measurements of the column density variance and the mean column density are performed over such dynamic range that they represent the true underlying total variance and mean. This is clearly not the case if the column densities are measured only over a relatively narrow dynamic range, as would be the case with, e.g., ¹³CO line emission observations (A_V ≈ 3−6 mag, e.g., Goodman et al. 2009).

With this background in mind, we present in the following two approaches for analyzing the ℳ_s − σ_N/⟨N⟩ relation. First, we describe a calculation in which we simply correlate the column density dispersions measured from the A_V data above A_V ≳ 7 mag with ℳ_s. This experiment can serve as a point of comparison for other numerical/observational works, but we cannot directly connect it with the formalism presented earlier in this section. Second, we make an effort to estimate the mean and variance of the observed column density fields in a dynamic-range-independent manner by assuming that the underlying column density PDF is log-normal. This assumption allows us to connect the relation directly to Eqs. (20)–(23) and to estimate how strongly the density fluctuations couple with the turbulent energy, i.e., estimate the coefficient b (Eq. (22)).

5.2.1. Model-independent approach

We first examine the relationship between the sonic Mach number, ℳ_s, calculated from the ¹³CO line emission and the standard deviation of the mean-normalized column density data, , defined as²(24)where n is the number of column density values N that enter the calculation and ⟨N⟩ is the mean column density in the field. We examined ℳ_s and only in the case of the A_V = 7 mag box, because of potential completeness issues in with ellipsoidal cloud definition. Only data above A_V ≥ 7 mag was considered. The Mach numbers were calculated from the one-dimensional velocity dispersions derived in Sect. 2.4 (given in Table 1) with the expression (25)where c_s = (kT/μm_H)^1/2 is the isothermal speed of sound and the velocity dispersion of CO is assumed to represent turbulent velocity. In the lack of measurements of the gas kinetic temperatures, we chose T = 15 K, consistent with the observations of Carey et al. (1998) and Pillai et al. (2006), and μ = 2.33 for the calculation of the speed of sound for all IRDCs. Equation (25) assumes that the velocity field is isotropic and the 3D velocity dispersion is related to the line-of-sight velocity dispersion with .

Figure 7 shows the relation between ℳ_s and in the IRDCs. The figure shows the data separately for the two different methods to calculate ℳ_s (see Sect. 2.4). The Mach numbers of the clouds span the range ℳ_s = 8−35, while the column density dispersions extend over a narrower range of . This range is in agreement with BT09 who concluded that there are differences in the widths of the column density PDFs of these clouds (BT09 used the “ellipsoid” area in determining their cloud region; see Sect. 3.2). The small number of clouds in our sample limits the ability to establish whether the two variables or correlated or not. Also, the uncertainties in both variables are likely dominated by the choices made in defining the cloud region instead of statistical uncertainty. The Pearson’s correlation coefficients for both choices of velocity dispersions are below 0.8 (which corresponds to the p-value of p ≈ 0.02 when degrees of freedom is 6). We also fitted the data with the linfitex routine from the mpfitfun library (Markwardt 2009). In the fitting, we used as uncertainties the relative uncertainty of 20% for (originating mostly from the zero-point error of the data, σ(A_V) ≈ 1 mag) and 30% for ℳ_s (originating mostly from the uncertainty in temperature, which is typically 10–20 K in IRDCs). The resulting slopes (and their errors) imply that the correlation of variables is not significant. However, we emphasize here that, as explained in Sect. 5.2, this relation cannot be straightforwardly linked to the theory describing the ℳ_s − σ_ρ/⟨ρ⟩ relation in turbulent media. We also note that uncertainty in opacity-law somewhat increases the uncertainty of the derived relationship (see discussion in Sect. 5.2.2).

5.2.2. Approach based on a log-normal column density PDF

We can connect the observed column density PDFs and ℳ_s values intimately into the theory of turbulent media (Sect. 5.2) by interpreting the PDFs shown in Fig. 4 as the high-column density sides of log-normal-shaped PDFs and estimating the shape parameters of the underlying PDFs through fits of log-normal functions. The log-normal function is given by Eq. (1), and we use in it the variable x = N, which is the (non-normalized) column density. With this definition, μ and σ_lnN are then the mean and standard deviation (in logarithmic units) of the log-normal distribution of column densities. Note that this formulation describes the PDF of non-normalized column densities: we do not know a priori the underlying mean of the column density field which should be used for normalization. From the fit, we determine the parameters (μ,σ_lnN), which are used to calculate the mean, ⟨N⟩, and standard deviation, σ_N, in non-logarithmic units \arraycolsep1.75ptThe standard deviation of the mean-normalized, non-logarithmic column densities is finally the standard deviation of non-normalized data divided by the mean m(28)The parameters resulting from the fit are listed in Table 2. Figure 8 shows the resulting ℳ_s − σ_N/⟨N⟩ relation. The Mach numbers were calculated as described earlier in Sect. 5.2.1. We show the relation only for the case when the cloud is defined with the A_V = 7 mag box, because the PDFs of the ellipsoidal definition can suffer from incompleteness (see Sect. 3.2).

The advantage of using a log-normal model to estimate the shape parameters of the PDF is that allows estimation of the total standard deviation of column densities, and from therein, directly connect the shape parameters to the theory of turbulent density fluctuations (see Sect. 5.2). Under this scheme, we can also combine the ℳ_s − σ_N/⟨N⟩ measurements of the IRDCs with additional data points from K09 (listed in Table 3, calculated from the data presented in Table 1³ of K09), who measured the standard deviation of column densities in a sample of nearby (D < 200 pc) clouds. K09 data allow measurements of the standard deviation of the column densities in the range 1  mag < A_V ≲ 25  mag. Importantly, the standard deviation measured at this range measures reasonably accurately the standard deviation measured from any wider range, assuming that the shape of the PDFs is even approximately correct. This is because the probability of having column densities above A_V > 25 mag is, in fact, so small that the contribution of those column densities to the standard deviation becomes negligible. For example, assuming that the true underlying PDF is similar to that of Ophiuchus (see Fig. 4), but extends onwards from A_V > 25 mag in a power-law (or log-normal) fashion, the standard deviation calculated from the range 1  mag < A_V ≲ 25  mag is only 5% smaller than the standard deviation calculated from the all data above A_V > 1 mag. In summary, we can consider the K09 data compatible with the measurements of σ_N/⟨N⟩, because the total standard deviation seems to be reasonably well recovered by them. However, we acknowledge that the measurement is based on different approach than the one we use in the case of IRDCs. It remains to be explored with a larger sample of IRDCs which exactly are the conditions under which the assumption of the compatibility is valid.

We show in Fig. 8 the ℳ_s − σ_N/⟨N⟩ values for 14 nearby clouds with red triangles, calculated using the K09 PDF data and Dame et al. (2001) CO data. We list in Table 3 the Mach numbers derived for nearby clouds that span the range ℳ_s ≈ 7−22. The temperature of T = 10 K was adopted in the derivation. We also show in the plot the values derived by Brunt (2010) and Padoan et al. (1997b) that fall well into the range of our measurements. We note that the velocity dispersions for nearby clouds are derived from ¹²CO data that generally probes slightly larger scales than ¹³CO data. However the effect of this is relatively small; we calculated that in Perseus, Ophiuchus and Taurus the mean difference between the ¹²CO and ¹³CO linewidths is a factor of 1.2. We made an experiment in which we multiplied the linewidths of IRDCs by this factor and it does not significantly affect the relations derived in this section.

Again, judging whether the ℳ_s − σ_N/⟨N⟩ data of our IRDC sample are correlated or not suffers from the low number of clouds. If we treat the measurements for IRDCs and nearby clouds as equals (i.e., assume that both approaches trace reasonably well the total column density fluctuations), the Pearson’s correlation coefficients for the data are 0.61 and 0.66 for the Gaussian and fixed range velocity dispersions. This coefficient indicates correlation of variables with high significance (degrees of freedom = 20, p ≲ 0.001). We also fitted a linear model to the data again using the linfitex procedure. In this calculation, we used the relative errors of 30% for the ℳ_s values, and 30% for the σ_N/⟨N⟩ values (the error of σ_N/⟨N⟩ originates almost entirely from the uncertainty in the mean μ. We use the typical value given by the linfitex procedure). The slopes resulting from this were a₁ = 0.051 ± 0.016 and a₁ = 0.042 ± 0.013 for the Gaussian and fixed range velocity dispersions, respectively. Thus, the derived slopes are approximately 3-σ offset from zero. The error-weighted mean of the slopes is a₁ = 0.047. We note that choosing the temperature of T = 15 K for nearby clouds (uniformly with the value chosen for IRDCs) instead of T = 10 K does not affect the significance of the correlation and results to similar slopes (a₁ = 0.050 ± 0.015 and a₁ = 0.039 ± 0.012 for the Gaussian and fixed cases).

We note that the uncertainty of our adopted  ratio is expected to play only relatively minor role in the uncertainty of the ℳ_s − σ_N/⟨N⟩ relation. This is because we correlate specifically the mean-normalized column density variances with ℳ_s. Even though possible use of an erroneous opacity-law between NIR and MIR reflects directly to the column densities, their effect to mean-normalized column densities should be a second-order effect. This is because, to first order, the column densities used in the PDFs (i.e., A_V ≥ 7 mag) depend linearly on the opacity-law: for example, increase of κ₈ by a factor of 1.5 would lower the maximum column densities approximately by that factor. And the PDF of a mean-normalized variable is invariant to multiplication of the variable by a constant factor. However, the uncertainty of the opacity-law does introduce a second-order effect to the shapes of the PDFs, because the response of column densities to the opacity-law is not exactly linear. Below A_V ≲ 10−20 mag the large-scale component of the column density field becomes significant compared to the smaller-scale component and the signal of the combined column density mapping technique (by construction, see Sect. 3.1) switches from being MIR-dominated to NIR-dominated. We illustrate the response of the combined column densities to the adopted ratio in Fig. C.1 (online only). In this test, we kept the NIR-opacity constant and varied the MIR opacity by 15%, 30%, and 60%. The experiment shows that in these cases the PDF shape of even a mean-normalized column densities below A_V ≲ 10−20 mag could be altered by the combination procedure. However, we include in the PDFs only column densities above A_V ≥ 7 mag, and therefore it seems justified to assume that this kind of bias is clearly a second-order effect. For example, in the test above (which is for Cloud A), the σ_N/⟨N⟩ changes about 5% and 20% for the cases where MIR opacity is underestimated by 30% and 60%, respectively.

Based on the experiments described above, we conclude that the observational data tentatively suggest correlation between the sonic Mach number and column density dispersion in molecular clouds. These measurements are, to our knowledge, the first direct measurement of the relationship of these variables and also the first direct detection of the correlation. However as also described above, we recognize the uncertainty of our measurements, rising mainly from the very small number of clouds in our sample and the uncertainties in defining “the cloud area” in our sample that is not tailored for the measurements of this particular quantity. It would require a larger sample of IRDCs to strengthen this conclusion. In practice, such a study would need to target clouds in locations that are somewhat offset from the Galactic plane and as isolated as possible.

Finally, we briefly speculate about the interpretation of the slope value we derive from the ℳ_s − σ_N/⟨N⟩ data. As noted earlier, Brunt et al. (2010) examined the 3D-2D variance ratio in a sample of isothermal turbulence simulations and measured the ratios to be between R =  [0.03,0.15] . If we simply assume that the ratio is in this range, the mean slope gives the values  with the mean at b = 0.20. We report the 3-σ confidence range for this mean as the 3-σ deviation from the highest and lowest possible values: upper limit is given by  and the lower limit by . Thus, we report as our measurement . This result points toward a somewhat lower value of b than derived by Brunt (2010) (b = 0.49⁴) and Padoan et al. (1997b) (b ≈ 0.5). Low values of b are generally predicted by (non-magnetized) isothermal simulations in which the turbulence driving is purely solenoidal (b = 1/3, e.g., Federrath et al. 2008b, 2010), while increasing the fraction of compressive modes increases the value of b toward unity (e.g., Federrath et al. 2010). The presence of magnetic fields also modifies the ℳ_s − σ_ρ/⟨ρ⟩ correlation, and in general, marginally flattens the slope (lowers the b) compared to non-magnetized cases (e.g., Price et al. 2011). Also, strong density dependence of the magnetic field makes the density fluctuations more weakly dependent on ℳ_s (e.g., Ostriker et al. 2001; Price et al. 2011; Molina et al. 2012). Even further, when the density dependence of the magnetic field strength approaches B ∝ ρ, the ℳ_s − σ_N/⟨N⟩ is predicted to show no correlation at Mach numbers greater than ℳ_s ≳ 7 (Molina et al. 2012). In this context, the fact that we detect a ℳ_s − σ_N/⟨N⟩ correlation up to high Mach numbers supports the picture in which the magnetic field depends relatively weakly on density. This has also been suggested through Zeeman effect observations (B ∝ ρ^0.65 at n > 300 cm^-3, Crutcher et al. 2010).

5.3. Comparison with earlier work

How do our results regarding the ℳ_s − σ_ρ/⟨ρ⟩ relation compare with the previous works on the topic? The ℳ_s − σ_ρ/⟨ρ⟩ relation (or ℳ_s − σ_N/⟨N⟩ relation) has not been widely examined by observations. The most comprehensive work to our knowledge is that of Brunt (2010) who exploit CO line emission data and near-infrared dust extinction data to estimate ℳ_s and σ_N/⟨N⟩, respectively, and use the method of Brunt et al. (2010) to infer σ_ρ/⟨ρ⟩ and the proportionality constant b. Their observations indicate the value of b = 0.49 (see Footnote ⁴). A similar technique was applied earlier to the dust extinction data of IC5146 by Padoan et al. (1997b) who estimated the value b ≈ 0.5. Both of these works, however, only measure one ℳ_s − σ_N/⟨N⟩ ratio pair, and not explicitly the correlation between the two variables.

Goodman et al. (2009) examined the PDFs and variances of the column density field in six sub-regions of the Perseus molecular cloud. The data (and the method) used by their work were very similar to what we employ in this paper. They tested if ℳ_s increases with increasing 2D variance (and mean) and they conclude that it does not. Thus, the result of Goodman et al. (2009) suggests no correlation between σ_ρ/⟨ρ⟩ and ℳ_s. In the scope of the data we show in Fig. 7, it seems clear that the contradiction between our and Goodman et al. (2009) results originate from the fact that the velocity dispersions in the sub-regions within Perseus span only a narrow range (Goodman et al. do not quote the linewidths of each sub-region, but the integrated linewidth over the entire Perseus is  ≈ 1 km s^-1, inevitably restricting the values to a rather narrow range). This leaves little leverage for the correlation to be detected, given the large scatter in the data shown in Fig. 8.

6. Conclusions

We present in this paper a new technique for probing the column density structure of IRDCs. The technique is based on combining column density maps derived from NIR observations of background stars shining through the IRDCs (UKIDSS/Galactic Plane Survey data) and from MIR (8 μm) shadowing features against the Galactic background (Spitzer/GLIMPSE survey). We demonstrate the use of the method by deriving column density data for ten IRDC complexes. We then use the resulting maps to examine the relation between the sonic Mach number and column density variance in the clouds. The main conclusions of our work are as follows.

• 1.

  The technique we present provides temperature-independentcolumn density data over a unique combination of spatialresolution (~2″) and dynamic range (A_V ≈ 1−100 mag). Unlike the mid-infrared technique alone, the new method is well-calibrated at column densities (1 mag  < A_V ≲ 10 mag). The new information on the low-column density surroundings of the IRDCs can significantly improve our understanding of the impact of the envelope material to the physics of the IRDCs (e.g. Hernandez et al. 2012). The technique is applicable to clouds up to about D ≲ 8 kpc, thereby potentially reaching a great number of IRDCs.

• 2.

  The column density PDFs of the IRDCs in the range A_V = 7−100 mag are, on average, consistent with a wide log-normal shape that has a width of σ_lnN ≈ 0.9. This shape is quite similar, although slightly shallower, to the PDFs of nearby active star-forming clouds reported by K09 in the column density range in which the data overlap, i.e., A_V = 7−25 mag.

• 3.

  The IRDCs contain relatively high amount of high-column density material, comparable to the fraction observed in nearby active star-forming clouds such as Orion A, and possibly even exceeding it. The average cumulative distribution function of IRDCs between 7  mag < A_V < 100 mag behaves exponentially, CDF ∝ e^{−0.07 × N}. In the framework of the recent work on the star-forming rates of molecular clouds (Lada et al. 2012), this suggests that the star-forming rates (per unit cloud mass) of IRDCs are at least comparable to active nearby clouds.

• 4.

  We present the first direct observational measurement of the relation between the turbulent energy (sonic Mach number, ℳ_s), and density fluctuations (column density dispersion, σ_N/⟨N⟩), in molecular clouds. We report a tentative detection of correlation with about 3-σ confidence. Linear model fits to the (ℳ_s,σ_N/⟨N⟩) data yields the slopes of σ_N/⟨N⟩ = ({ 0.051,0.042 }  ±  { 0.016,0.012 } ) × ℳ_s for the two methods we use for measuring ℳ_s.

• 5.

  Our results suggest the correlation coefficient between the sonic Mach number and volume density variance, σ_ρ/⟨ρ⟩, in molecular clouds. In this expression, the quoted uncertainties indicate the 3-σ confidence interval assuming that the 2D-3D variance ratio R is between R =  [0.03,0.15]  (as measured by Brunt et al. 2010). The value of b we derive is in agreement with earlier derivations by Padoan et al. (1997b) and Brunt (2010), that are, however, based on measurements of only single (ℳ_s,σ_N/⟨N⟩) ratio pairs.

• 6.

  When combined with recent numerical and analytical predictions (Price et al. 2011; Molina et al. 2012), the detection of correlation between ℳ_s and σ_N/⟨N⟩ is suggestive of relatively weak density dependence of the magnetic field strength in the clouds. This is in agreement with the field strength measurements of Crutcher et al. (2010).

While the conclusions given above remain tentative, the technique we present provides a good basis for measuring the relationship between ℳ_s − σ_N/⟨N⟩ in a statistically significant sample of molecular clouds in the future.

---

Acknowledgments

The authors are grateful to Henrik Beuther, Blakesley Burkhart, Christoph Federrath, Simon Glover, Thomas Henning, Ralf Klessen, and Lukas Konstandin for fruitful discussions and comments regarding the manuscript. We thank the referee whose comments improved the manuscript significantly. The work of J.K. was supported by the Deutsche Forschungsgemeinschaft priority program 1573 (“Physics of the Interstellar Medium”). J.K. gratefully acknowledges support from The Finnish Academy of Science and Letters, Väisälä Foundation for this work. J.C.T. acknowledges support from NSF CAREER grant AST-0645412; NASA Astrophysics Theory and Fundamental Physics grant ATP09-0094 and NASA Astrophysics Data Analysis Program ADAP10-0110. This work was assisted by the ASTROWIN workshop program supported by the U. Florida Dept. of Astronomy. 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.

References

Online material

Appendix A: CO spectra towards the IRDCs

Appendix C: Response of the technique to changes of the opacity-law

All Tables

All Figures