MASSIVE CLUMPS IN NGC 6334

# Massive Clumps in the Ngc 6334 Star Forming Region

Diego J. Muñoz11affiliation: dmunoz@cfa.harvard.edu , Diego Mardones, Guido Garay and David Rebolledo Departamento de Astronomía, Universidad de Chile, Casilla 36-D, Santiago, Chile Kate Brooks Australia Telescope National Facility,P.O. Box 76, Epping NSW 1710 Australia Sylvain Bontemps Observatoire de Bordeaux, 2 rue de l’Observatoire,33270 Floirac, France
###### Abstract

We report observations of dust continuum emission at 1.2 mm towards the star forming region NGC 6334 made with the SEST SIMBA bolometer array. The observations cover an area of square degrees with approximately uniform noise. We detected 181 clumps spanning almost three orders of magnitude in mass (3 ) and with sizes in the range 0.1–1.0 pc. We find that the clump mass function is well fit with a power law of the mass with exponent (or equivalently ). The derived exponent is similar to those obtained from molecular line emission surveys and is significantly different from that of the stellar initial mass function. We investigated changes in the mass spectrum by changing the assumptions on the temperature distribution of the clumps and on the contribution of free-free emission to the 1.2 mm emission, and found little changes on the exponent. The Cumulative Mass Distribution Function is also analyzed giving consistent results in a mass range excluding the high-mass end where a power-law fit is no longer valid. The masses and sizes of the clumps observed in NGC 6334 indicate that they are not direct progenitors of stars and that the process of fragmentation determines the distribution of masses later on or occurs at smaller spatial scales.

The spatial distribution of the clumps in NGC 6334 reveals clustering which is strikingly similar to that exhibited by young stars in other star forming regions. A power law fit to the surface density of companions gives .

stars: formation — molecular clouds: individual( (catalog NGC 6334))—stars: initial mass function
slugcomment:

## 1. Introduction

### 1.1. Massive Star Formation and Molecular Cloud Structure

Half of the mass in the interstellar medium is in the form of molecular gas exhibiting a broad range of structures, ranging from small isolated clouds, with masses of a few and subparsec sizes, to Giant Molecular Clouds (GMCs), with masses of several times and sizes of 100 pc (Blitz, 1993). We adopt the nomenclature used in the review of Williams, Blitz & McKee (2000) to refer to the different observed molecular structures. GMCs are the sites of most of star formation activity in the Milky Way, and in particular of high mass stars, which are usually born in clusters within massive cores. In order to understand how massive cores form from the GMC complexes we must understand how fragmentation and condensation proceed within them. It is essential for this purpose to determine the physical properties of complete samples of massive cores within GMCs.

Molecular line surveys, at millimeter and sub-millimeter wavelengths, have revealed the structure of GMCs to be highly inhomogeneous and clumpy (Blitz, 1993; Evans, 1999). These surveys have shown that the mass spectra of clouds (Sanders et al., 1985; Solomon et al., 1987), clumps (Blitz, 1993; Kramer et al., 1998; Lada, 1999; Williams, Blitz & McKee, 2000) and total mass of embedded clusters (Lada & Lada, 2003) are similar to one another. These mass spectra are notably different than the stellar mass spectrum: the Initial Mass Function (IMF; Salpeter, 1955. In particular, several works on molecular line mapping of GMCs, show that their mass spectra follow a power law with nearly the same exponent, , where (e.g. Blitz, 1993; Williams, Blitz & McKee, 2000).

The development of large bolometer arrays during the last ten years has permitted to carry out extended millimeter and sub-millimeter continuum surveys, allowing the direct dust mass determination of clumps. Motte et al (1998) mapped the Ophiuchus cloud with the IRAM 30-m telescope and found a distribution of clump masses similar to Salpeter’s IMF. Several works (Johnstone et al., 2001; Beuther & Schilke, 2004; Mookerjea et al., 2005; Reid & Wilson, 2005; Johnstone et al., 2006) have reported mass spectra of dust cores with indexes similar to Salpeter’s IMF, and different from those derived from molecular line studies. If the mass function of cores is similar to the Salpeter IMF, independent of the range of masses involved, then the star formation process within GMCs would be defined in the earliest stages as a result of cloud fragmentation.

In this paper we present a large scale 1.2 mm continuum study of the NGC 6334 GMC aimed to find and study a complete sample of massive cores.

### 1.2. NGC 6334 and NGC 6357

NGC 6334 is one of the nearest and most prominent sites of massive star formation, at a distance of only 1.7 kpc (Neckel, 1978). The central region of NGC 6334 consists of a  pc long filament with seven sites of massive star formation. Within them there is a wide variety of activity associated with star formation, such as water masers, HII regions (Rodríguez, Cantó & Moran, 1982; Carral et al., 2002), and molecular outflows.

Early FIR (McBreen et al., 1979) and radio (Rodríguez, Cantó & Moran, 1982) surveys characterized the overall properties of the main sites of massive star formation, whereas NIR studies (Straw et al., 1989; Straw & Hyland, 1989) revealed the cluster forming nature of many of these sites. Kraemer & Jackson (1999) and Burton et al. (2000) review the different notations used by previous authors to identify bright sources. Sub-millimetric and millimetric results on NGC 6334 have been previously reported by Gezari (1982), Sandell (2000) and McCutcheon et al. (2000), all of which focused in the northern portion of the main filament: sources I and I(N) which have been recently resolved into smaller cores (Hunter et al., 2006). In the case of Sandell (2000), the sources I and I(N) were redefined as the central peaks of each respective source. This after observations at 350 m, 450 m, 380 m, 1.1 mm and 1.3 mm showed inner structure for Gezari’s cores. In Sandell (2000) sources I and I(N) have angular sizes of and respectively. These sizes are times smaller than our massive clumps cl1 and cl2 identified as I and I(N) respectively. For the present work, cl1 encloses sources I and I(NW) of Sandell’s nomenclature, while cl2 encloses I(N), SM1,SM2,SM4 and SM5 as well as considerable extended emission in both cases. Hence, the comparison between that work and ours is not straighforward.

Following the evolutionary sequence proposed by Beuther et al. (2006), NGC 6334 and NGC 6357 appear to be in an intermediate phase of Massive Starless Clumps and Protoclusters. The latter region seems more evolved that the former since compact HII regions are still prominent in NGC 6334 at radio wavelengths while the structure of NGC 6357 is more disrupted, suggesting that massive stars have already shaped the mother clouds. This is also supported by the presence of more infrared sources in NGC 6357. From this evidence, and the lack of significant amounts of cold material in between the two regions, we will consider NGC 6334 and NGC 6357 as two independent regions, defining NGC 6334 as the southeastern portion of the map as shown in Figure 1.

Even though NGC 6334 is one of the closest GMCs, it is still far compared to low-mass star-forming regions (e.g. Ophiucus, Orion B, Taurus). As a consequence, the detected clumps in NGC 6334 are larger and considerably more massive, and can be considered likely cluster-forming cores (Motte et al., 2003; Ward-Thompson et al., 2006; Beuther et al., 2006) and we are unable to resolve their inner structure.

## 2. Observations and Data Reduction

The regions NGC 6334 and NGC 6357 were mapped using the 37 channel SEST Imaging Bolometer Array (SIMBA) in the fast mapping mode in three different epochs: July 2002, September 2002 and May 2003. The passband of the bolometer has an equivalent width of 90 GHz and is centered at 250 GHz (1.2 mm). The half-power beamwidth of the instrument is 24″ giving a spatial resolution of  0.2 pc. Ninety five observing blocks were taken towards the NGC 6334 and NGC 6357 regions, with typical extension of , to sample a total area of square degrees between and . Skydip observations were done approximately every two hours to determine the zenith opacity at 250 GHz. Typical opacities were with values ranging from 0.17 to 0.4 in a few cases. We also checked pointing on  Carinae every two hours and found a typical rms deviation of 3-5” in azimuth and elevation. Every night we observed Uranus for flux calibration.

The SIMBA data were reduced using the MOPSI program written by Robert Zylka after conversion by the simbaread program written at ESO. The SIMBA raw data consists of a time series for each of the 37 bolometers (channels) in the array. The time series includes the counts per channel and sky position. The reduction procedure first removes the brightest data spikes. Next a low order baseline in time is fit to the full observation file for each channel, and a zero order baseline is fit in azimuth for each channel. The data is then deconvolved by the time response function of each channel, as measured by the SEST staff. Gain elevation and extinction corrections are applied next. An iterative sky noise reduction procedure is then applied, where the counts of each channel are correlated with those of the other 36 channels, yielding a so-called flat field correction to calibrate the relative sensitivity of each pixel. A source image is finally produced by averaging the flux of all channels as they pass through the same position on the sky.

The sky noise reduction algorithm includes the flux coming from both the source and the sky simultaneously. If the source extends over several channels, this introduces spurious correlated flux which hampers the sky noise reduction procedure. To avoid this, a smoothed model of the source flux distribution on the sky is subtracted from the raw time series data. Thus, the sky noise reduction procedure can be repeated, finding a better source model with each iteration.

The calibration was derived from maps of Uranus. The resulting multiplicative factor varied between 0.06-0.09 Jy count beam. Finally, reduced images were combined with the MOPSI software to produce the final map (Figure 1). In order to reduce noise further, a Gaussian smoothing of 30” was applied to the final image and the map edges were removed. The typical rms noise of the final map is 25 mJy beam. SIMBA observations usually have an absolute flux uncertainty of 20% (Faúndez et al., 2004).

## 3. Results

### 3.1. Clump-finding Algorithms

Different clump-finding algorithms have been used to study the substructure in molecular clouds. Among the most used are Clumpfind (Williams, de Geus & Blitz, 1994) and Gaussclumps (Stutzki & Gusten, 1990). Both algorithms have different biases, but find similar clump distributions (e.g. Schneider & Brooks, 2004). We use the clumpfind algorithm because it makes no assumptions about inherent clump shapes. Clumpfind first finds the brightest emission peak in the image, then it descends to a lower contour level and finds all the image pixels above this level, associating them to the first peak (clump) if contiguous, or else defines one or more additional clumps. The spatial separation of clumps is defined along saddle points. We used a conservative lower detection threshold of 75 mJy () per beam. We found 347 clumps in the whole image, 182 of which are in the NGC 6334 region. Only two clumps are likely to be fictitious based on their small effective radii and location close to the borders, and only one of them is in the NGC6334 region. Thus, we will use the remaining 181 clumps in our analysis of the region NGC 6334 in this paper. Table 1 lists the clumps in our sample associated with previously known radio and IR sources. These sources are all located in the brightest part of the filament as seen in Figure 2.

### 3.2. Clump Size Distribution

Figure 3 shows a histogram of the clump size distribution. The effective radius, or size, is determined from the angular area encompassed by each clump (from the clfind2d output) assuming a distance to NGC 6334 of 1.7 kpc. The sizes range from 0.1 to 1.0 pc, with a median of 0.36 pc. These values are similar to those derived for clumps within GMCs (Blitz, 1993; Williams, Blitz & McKee, 2000; Pudritz, 2002; Faúndez et al., 2004; Beuther et al., 2006).

### 3.3. Mass Estimates

Since the dust emission at 1.2 mm is most likely to be optically thin (e.g. Garay et al. 2002), the mass of each clump can be estimated from the observed flux density. For an isothermal dust source, the total gas mass is related to the observed flux density at an optically thin frequency as (c.f. Chini, Krügel & Wargau, 1987)

 Mg=SνD2RdgκνBν(Td), (1)

where is the dust mass absorption coefficient, is the source distance, is the dust to gas mass ratio, and is the Planck function. In more convenient units the gas mass can be written as

 Mg=20.4(S250GHzJy)(Dkpc)2(0.01Rdg) ×(1cm2g−1κ250GHz)[e12KTd−1]M⊙. (2)

Using a dust mass absorption coefficient of (Ossenkopf & Henning, 1994); a dust-to-gas ratio of 0.01; a distance of 1.7 kpc, and a dust temperature  K 111 We choose a value of 17 K because it lies in the typical range for cold dark clouds (Pudritz, 2002). It is also also a factor of two smaller 34 K Ð the average temperature of clouds with infrared counterparts (Faúndez et al., 2004)Ð making the comparison between both temperatures easier. , we computed the mass of each clump found by clfind2d. We find that the clump masses in NGC 6334 range from 3 to , with a mean value of and a median value of . The total mass of the clumps in NGC 6334 is .

### 3.4. Clump Mass Spectrum

The Clump Mass Function (CMF), , is defined as the number of clumps per unit mass,

 ξα(M)=dNdM≈ΔNΔM, (3)

where and are used to indicate observational estimates. Some workers prefer to use a logarithmic CMF, , defined as the number of clumps per unit logarithmic mass,

 ξx(M)=dNd(logM)≈ΔNΔ(logM). (4)

These two functions are related by the expression . If has a power law dependence with mass, , then the logarithmic CMF should also have a power law dependence with mass, , where the power law exponents are related by . A Salpeter slope corresponds to or .

The CMF is usually estimated from a histogram of the derived clump masses assuming it has a power law form. The distinction between the observationally derived and is made by Scalo (1998), Kroupa (2001) and Larson (2003). Its implications in the type of binning used in the histogram is also mentioned in Klessen & Burkert (2000). Figure 4 shows an histogram of the mass of the clumps within NGC 6334. The bin size has a constant value of 0.5. The completeness limit is estimated to be . A least squares linear fit to the versus relationship gives a slope of .

##### Variations of the CMF within NGC 6334.

Given the large number of clumps in our survey, we can assess possible changes in the CMF as a function of position within the cloud. We consider the three subregions within NGC 6334 defined in Figure 1 and construct a mass spectrum for each of them (see Figure 5). The CMF in NGC 6334a, which encompasses the central filament containing the most massive clumps and where star formation is clearly taking place, is well fitted with a power law with an exponent of , significantly shallower than the value determined for the whole sample. The exponent steepens for NGC 6334b (), which covers a much larger region than the central filament. For NGC 6334c, which excludes the main filament, the slope is even steeper, . Thus, the slope of the CMF in NGC 6334 depends on the location of the clumps within the cloud.

Out of the total mass of in 181 dense massive clumps in NGC 6334, are contained within the NGC 6334a region and are outside. The slope of the CMF of the whole region is dominated by the relatively low mass clumps in the outer region since they dominate by number. The bulk of the cloud mass is located in a few inner clumps which do not affect significantly the exponent of the derived CMF. This could shed light on the effects of the star formation activity, mass segregation, or coalescence in the clump mass spectrum within GMCs. We know that it is more likely to miss low mass clumps due to confusion within the region NGC 6334a than in the outer regions under study. This bias naturally flattens the slope of the mass spectrum in the inner region. However, we find in §3.6.1 that clumps are not only concentrated by mass towards the center, but also by number, hinting that this result is partly real. Higher angular resolution observations will be needed to settle this issue.

### 3.5. Possible Uncertainties in the CMF

In the above analysis we assumed a single temperature for the whole ensemble of clumps. This is clearly an approximation; the clumps themselves are not isothermal and the temperature is likely to be different from clump to clump. In addition, we assumed that all the detected 1.2-mm emission is due to dust thermal emission. It is possible, however, that some of the 1.2-mm emission is due to free-free emission from ionized gas. In what follows we assess these two assumptions and quantify their effects on the derived CMF.

#### 3.5.1 Temperature

Assuming that all clumps within a GMC are isothermal and have the same temperature is clearly a rough approximation. In particular, clumps with already formed stars are expected to be warmer than clumps with no signs of embedded objects. Beuther & Schilke (2004) have argued that this assumption introduces an uncertainty in the derived slope of the mass spectrum. If higher temperatures are adopted for the more massive clumps their derived masses will decrease, whereas if lower temperatures are adopted for the less massive clumps their derived masses will increase: the change would steepen the slope of the spectrum.

To study the possible effects of temperature differences, we use MSX mid-infrared and/or ATCA cm-continuum observations to determine the presence of embedded heating sources which are likely to be responsible for temperature differences between clumps. Figures 6 and  7 show, respectively, images of the MSX and ATCA emission overlayed with contours of the SIMBA 1.2-mm continuum emission. We identified by visual inspection clumps associated with extended or point-like infrared counterparts (Figure 6), or clumps associated with significant radio continuum (free-free) emission (Figure 7). These clumps are likely to have embedded sources and are indicated with green contours in both images. Clumps in red contours appear to be free of embedded infrared or radio sources.

For the purposes of constructing the CMF, the power-law index does not depend on the temperature chosen for the whole ensemble of clumps, and the choice of higher temperatures only displaces the histogram to the left. Any slight variation of is due to the intrinsic problem of binning the data. Nevertheless, the choice of two different temperature can change the shape of the histogram, but essentially depending on the ratio between the temperatures chosen rather than the values themselves.

We identified 16 out of 181 cores that appear to be warmer than the rest. For these clumps we adopt a temperature of 34 K, the average temperature of massive clumps with embedded IRAS sources as determined by Faúndez et al. (2004). In particular, the temperatures assigned to the cores associated with NGC 6334 I and NGC 6334 I(N), of 34 K and 17 K respectively, are in accordance with the values given by Gezari (1982). However,  Sandell (2000) argues that source I is much hotter ( K) while his estimate for I(N) is 30 K, just a factor of two larger than our estimate. His estimate of the mass of NGC 6334 I of 200 Ð in contrast to the in the present workÐ is not solely explained by the temperature difference but also by the integrated flux. With a resolution of 6” at 800 m, his estimate for the size of source I is , implying an effective radius times smaller than our estimate. Similarly, the definition of Sandell (2000) of NGC 6334 I(N) corresponds to the central peak of the larger clump detected by Gezari (1982). Sandell (2000) obtains a mass of for source I(N) after assuming a temperature of 30 K (see Sandell, 2000 and references therein). Adding up the contributions from the different cores resolved within I(N) and including the surrounding cloud, Sandell (2000) finds that the mass of the whole I(N) region is This value is in agreement with other results found in the literature and is consistent with our result after taking into account that our choice for the temperature is 17 K. We remark that estimating an exact value of the temperature for each of the clumps of the cloud is not relevant for the statistical analysis we carry out in the present work.

The newly adopted temperatures imply a repositioning of 10% of the clumps in the mass histogram. The new mass histogram, made assuming a two temperature cloud ensemble, is shown in Figure 8. A linear regression yields a best fit value of (dotted line), but a fit with Poissonian error bars yields (solid line). Thus, even though the warmer clumps are also preferentially the most massive, they do not concentrate solely on the most massive bin and do not affect the derived slope of the mass spectrum significantly.

The repositioning of clumps in the histogram could have a more dramatic effect when the number of clumps is considerably smaller. For example, the area NGC 6334a includes 14 of the 16 warmer clumps in NGC 6334 and a total of only 40 clumps (figure 5). When using two temperatures in NGC 6334a, the shape of the histogram indeed steepens, but the large relative errors yields a slope . Excluding the last bin, containing only one object, we obtain a slope of . In both cases the slope remains consistent with the value from a single temperature clump mass distribution.

Free-free emission is the main continuum contributor at radio frequencies. At frequencies of hundreds of GHz, the contribution from bremsstrahlung emission is usually neglected and all the detected emission is assumed to be due to dust. Since massive stars are being formed within NGC 6334, it is plausible that the 1.2-mm emission is not completely due to cool dust emission and that the ionized gas excited by these stars contribute an important amount. Brooks et al. (2005) studied the Keyhole nebula and found that the 1.2 mm emission towards the HII region Car-II is strongly correlated with the 4.8 GHz continuum emission and that there is a lack of molecular-line emission. They concluded that the 1.2-mm flux from the components of Car-II arise from free-free emission associated with ionized gas and not from cool dust emission associated with molecular gas.

Here we make a correction to the observed emission at 250 GHz by free-free contamination by estimating the expected ionized gas flux density at 250 GHz from the observed flux density at 1.6 GHz. Assuming that the free-free emission is optically thin at both frequencies, the ratio of the emissivities is proportional to the ratio of the factors (Rybicki & Lightman, 1979), where is the Gaunt factor. The exponential is essentially 1 at both frequencies. We did not use the usual radio approximation for the Gaunt factors given by Altenhoff et al. (1961), but computed them more precisely using quantum mechanical calculations following the work of Menzel & Pekeris (Menzel & Pekeris, 1935; Sommerfeld, 1953). Table 2 lists the calculated gaunt factors averaged over a Maxwell-Boltzmann distribution of velocities with temperatures between and  K, typical for HII regions in massive star forming regions (Beckert et al., 2000). At 1.6 GHz the computed Gaunt factors are only % lower than the usual Altenhoff et al. approximation. On the other hand, the calculated Gaunt factors at 250 GHz are typically 15% smaller than the values obtained from the Altenhoff et al. approximation.

Table 3 summarizes the values of the estimated free-free emission at 250 GHz from selected clumps associated with regions (see Figure 2). We subtracted this flux from the measured 1.2-mm flux density to estimate the actual contribution from dust, rederiving the mass of each clump. The high mass bins are the most affected, but the best fit exponent () is consistent with the previous value.

#### 3.5.3 Binning

Given a distribution of clump masses, the binning process invariably loses information (Rosolowsky, 2005). The mass histogram can be interpreted as the derivative of a cumulative number function which counts clumps with mass greater than ,

 N(M)=∫∞Mξα(m)dm∫∞0ξα(m)dm . (5)

Some authors argue in favor of using the cumulative number function (Johnstone et al., 2000b, 2001; Kerton et al., 2001; Tothill et al., 2002) to avoid the loss of information. If the cumulative number function has a power-law dependance with mass, , then the CMF is also a power law . From the definitions in section 3.4, it is clear that

However, if there is an upper mass limit in the distribution of cores, then the cumulative mass function is not a power-law , showing considerable curvature at the high mass end. can be approximated by a power-law with index only at masses (or when ). For finite upper mass limit,

 N(M)=C1M−x+C2 , (6)

where becomes unimportant for small masses (see Rosolowsky, 2005; Reid & Wilson, 2006b; Li et al., 2006). A power-law form, , approximates the true CMF asymptotically towards low masses. Thus, a power law fit to the CMF must be applied within a range of masses that avoids both the incomplete low-mass range and the cutoff high-mass range. When fitting a power-law function to an observed cumulative mass function in a mass range the slope invariably increases, explaining the apparent Salpeter-like slopes found in previous studies. Furthermore, the slope thus obtained is strongly dependent on the break-point chosen to fit the high mass end of the cumulative mass function and does not reflect the underlying differential clump mass distribution ( or ).

For sample sizes of clumps or fewer, binning becomes an important factor in fitting a power law to differential mass functions and the use of a cumulative mass function is preferred. Johnstone et al. (2000b, 2001, 2006) and Reid & Wilson (2005, 2006a) use cumulative mass distribution functions to avoid this problem. Our sample (181 clumps) is large enough to analyse the data using either the cumulative or the differential mass functions. Figure 9 plots the normalized cumulative number function of clumps in NGC6334. A single power-law fit, which we have shown does not represent the true underlying clump mass function, gives . The dashed line in Figure 9 shows that the slope derived from the histogram (Figure 4) represents a good asymptote to the CMF, as predicted by the theory (Eq. 6). We also show in Figure 9 a fit to the top 10% of the cloud mass (, indicated by an arrow) which yields a slope (dotted line) much closer to Salpeter’s value. The best fit exponent changes to for , and to for . Even though these slopes are consistent with Salpeter’s IMF, they are an artifact of having an upper clump mass limit in the sample and do not reflect the true clump mass function

In summary, one must be aware of introducing biases when using the high mass range to fit power laws or when fitting broken power laws to the cumulative mass function. These problems are minimized (but still present) with samples larger than 100 clumps.

### 3.6. Spatial Distribution of Clumps

In order to understand the process of fragmentation, we need to explain how masses are distributed in clumps and how they are positioned in space. A complete theory of star formation must not only reproduce the mass function, but it must explain it in all its physical implications including how clumps, cores and stars are distributed spatially during the evolution of the GMC (Bonnel et al. 2006). The exhaustive study of the spatial distribution of young stars in the Taurus region by Gómez et al. (1993) was extended by other authors and compiled by Larson (1995). Here we undertake a study the degree of clustering of clumps in NGC 6334 using a similar approach.

#### 3.6.1 Clustering and Segregation of Clumps

The number density of clumps can give us insight about the actual state of fragmentation and how the clumps are distributed spatially independent of their mass. We study the number density of clumps by means of the Simple Grid222The Grid Method consists in binning the two-dimensional space with squares of side and then dividing the numbers of sources lying within each square by the area of it , to obtain number density in units of length and the Kernel Methods333The Kernel Method (Silverman, 1986) uses a kernel function offering the advantage that the density distribution is smoothed. In each point of , the kernel density estimator determines the density due to the contributions of all data points.. Both methods require a free parameter which determines the “resolution” of the number density estimator: the binning length in the grid technique, and the smoothing length in the Kernel Method, where a kernel is defined at each pixel of the map in Figure 1 by

 K(x,xi,y,yi)=12πe−(x2+y2)/2h2 , (7)

with the kernel density estimator defined by

 D(x,y)=1h2n∑i=1K(x,xi,y,yi) , (8)

where and are measured in pc ignoring the sky curvature. The grid bin was taken to be 2 pc while was chosen in order to smooth the distribution over an area similar in size to the area over which the simple grid technique smoothed the data points ().

We applied this analysis to the subregion NGC 6334b which covers and area of  pc. We find that the probability that the clumps are distributed at random within this subregion is , which argues in favor of clustered fragmentation at scales between 0.1 and 10 pc. Thus, we conclude that there is spatial segregation in clump number in addition to clump mass.

#### 3.6.2 The Nearest-Neighbor Distribution in NGC 6334

We calculate the nearest-neighbor distribution (i.e. the frequency distribution of the linear distance to the nearest neighbor of each clump) for the clumps in NGC 6334. We neglect those clumps located too close to the edge of the mapping area but do not exclude them from the total sample since they can be the nearest neighbor for an inner clump. We binned the nearest neighbor distances in intervals of 0.2 pc to construct the histograms shown in Figure 10. The nearest neighbor distribution is strongly skewed to small separations and very different from the distribution expected from random positions at the same mean density. The differential probability of observing at least one event in the interval which defines a ring surrounding a central source is (Gómez et al., 1993)

 ΔP(r1,r2,η)=e−πηr21−e−πηr22 . (9)

where is the average density, in this case given by  clumps pc. Figure 10 compares the Poisson distribution histogram (dashed line) with the normalized nearest-neighbor distribution for clumps. The median for the random distribution (1.7 pc) is a factor of greater than that derived from the actual clump distribution (dotted vertical lines in Figure 10). We would need to increase the average density by a factor of 5 to produce a random distribution with the same median separation as the actual clump distribution. Thus, we conclude that the median separation of the clump distribution is significantly smaller than the expected median separation of a random distribution.

### 3.7. The Average Angular Surface Density of Clumps in NGC 6334

Most of the newly formed stars in nearby regions of star formation are located in groups or clusters. The degree of clustering of pre-main sequence stars can be obtained by measuring their surface density as a function of angular distance, , from each star. This surface density can be fit by power-laws (Gomez et al. 1993); Larson (1995) found that there is a characteristic spatial scale ( pc) where the surface density changes slope. This scale could mark a transition between the regime of cores within molecular clouds and protostars within cores.

We applied the approach followed by Larson (1995) to the clumps in NGC 6334, which have typical sizes  pc and separations of  pc. We compute by taking each clump and dividing the surrounding area of the sky into a set of annuli of radius (with ), and counting the number of companion clumps in each annulus (Kitsionas et al., 1998). Then,

 Σ(¯¯¯¯θi)=1N∑Nk=1Nk(θi)π(θ2i−θ2i−1) ,¯¯¯¯θi=(θi+θi−1)2 . (10)

The results are plotted in Figure 11. The surface density is poorly fitted by a single power-law, but the fit is much improved when using a broken-power law. However, this reflects a sampling problem at large separations, since the area observed around NGC 6334 is not square, Its narrowest part has a width of  pc, just where the power-law apparently breaks. The reliable portion of the power-law has a slope of

 Σc∝θ−0.62 0.6pc≤θD≤10pc , (11)

which is remarkably similar to the one found by Larson (1995) for separations larger than 0.04 pc calculated from young stars in the Taurus-Auriga region:

 Σc∝θ−0.62 0.04pc≤θD≤2.5pc. (12)

These two results cover two different separation ranges in regions with widely different physical properties, yet the surface density distribution is the same. If clustering above the characteristic length of 0.04 pc maintains a self-similar behavior up to 10 pc, then the grouping of stars at large separations (between subclusters for example) could be determined from early stages of cloud fragmentationÐ before stars are formed.

## 4. Discussion

##### The CMF and its relation to the IMF.

The clump mass function is, on one hand tied to the formation and evolution of their parent molecular clouds, and on the other to the formation of their daughter stars. Provided that the mechanism of clump fragmentation and collapse to form stars is universal, then the IMF should be a direct consequence of the CMF. However, since a single clump generally forms multiple stars, we can’t expect the distribution of clump masses and stellar masses to have identical functional form.

We find that the mass spectrum of clumps in NGC 6334 has a power-law dependence with mass with an index (). This value is similar to values derived from dust continuum observations for other massive star forming regions (Kerton et al., 2001; Tothill et al., 2002; Mookerjea et al., 2005). It is also similar to values derived from molecular line observations, usually isotopomeric lines of CO, for clouds with similar total mass. For instance, Nozawa et al. (1991) finds in Ophiucus North (mass range of 4–250 ), Stutzki & Gusten (1990) finds the same exponent in M17SW (mass range of 10–3200 ), and Kramer et al. (1998) reports for NGC 7538 (mass range of 50–3.9 ). These results indicate that the clumps in all these massive star forming regions are not the direct progenitors of individual stars. The fact that the massive clumps follow a power-law mass function over a wide range in mass is remarkable. The similarity of the ”high mass” CMF slope to the Mass Function of GMC’s in our Galaxy (Sanders et al., 1985; Solomon et al., 1987; Pudritz, 2002) seems to support a hierarchical or fractal model of the distribution of gas in the Milky Way, where fragmentation and mass distribution can be interpreted as scale free.

In order to build a complete theory of star formation we must understand the process of clustered star formation in clumps satisfying the observed CMF () and leading to the observed IMF (). This has stimulated studies about the process of fragmentation of GMCs and clumps both theoretically and observationally. Hydrodynamic and magnetohydrodynamic numerical calculations suggest that turbulence might play a major role in the clump fragmentation process. For example Bonnell & Davies (1998) suggest that stars are formed through the competitive accretion of gas onto proto-cores within molecular clouds.

Recent sub-millimeter and millimeter continuum observations of low-mass star forming regions (e.g. Ophiucus, Serpens, Orion B) (Motte et al., 1998; Testi & Sargent, 1998; Johnstone et al., 2000b, 2001; Beuther & Schilke, 2004) have revealed a CMF with a slope similar to that to the IMF. The different values obtained from the molecular line observations described above and these continuum observations probably reflects the different mass ranges sampled. In fact, dust continuum emission observations toward more distant massive star forming regions have revealed CMF slopes consistent with a value of -0.6 (Kerton et al., 2001; Tothill et al., 2002; Mookerjea et al., 2005; this work), similar to that obtained from molecular line observations. From observations of a sample of massive clumps toward the massive star forming region NGC 7538, Reid & Wilson (2005) concluded that a Salpeter-like mass function is already established at the earliest stages of star formation. This is a surprising result since many of their clumps, with masses of , are still likely to be undergoing the process of fragmentation and can not be direct progenitors of individual stars. Even if the structure in ISM were fractal, self-similarity must break on small scales, where star formation is taking place. We note, however, that the result of Reid & Wilson (2005) comes from a fit to the high mass end of the cumulative mass function. We showed in section 3.5.3 that these fits are biased towards larger exponents. In addition to molecular line and dust continuum surveys, extinction maps can be used to map the dense molecular cores in star forming regions. Using this technique, Alves et al. (2007) have found a mass spectrum for cores in the Pipe Nubula that is surprisingly similar to the IMF. Their CMF displaced to higher masses with respect to the stellar IMF only by a factor of 3, suggesting a one-to-one mapping from cores to stars with a star formation efficiency of %. In their case, the cores detected span a range of masses of 1-10  approximately. This clearly indicates a much smaller scale than the clumps in our work whether determined by the tracer or just the spatial resolution of the observations. It is not clear if the scale of the observations determines the observed fragmentation conditions Ð i.e. the distribution of masses Ð or if the IMF is a result of the complex evolution of the accreting cores and their interplay with the harboring molecular cloud and other companions. The low star forming activity in the Pipe Nebula has led al07 to suggest that this CMF can be determined at early evolutionary stages.

Assuming that at scales of 0.4 pc and masses from to the Blitz slope () is valid, the question that arises is: What must happen to change the slope from Blitz-like to Salpeter-like? In the context of gravitational opacity limited fragmentation, gravitational collapse starts from density inhomogeneities and proceeds with cooling, which in turn produces smaller Jeans masses in the colder regions, favoring gravitational fragmentation on small scales. In a strictly self-similar regime, fragmentation should occur maintaining the same slope in the mass spectrum at any scale. But at some point the gas cores will not be supported by thermal or non-thermal motions, collapse will occur and stars will form, halting the clump fragmentation, while larger and less dense clumps will continue fragmenting. Eventually large clump masses will be depleted, and the number of small mass clumps will increase steepening the slope to eventually reach a Salpeter value. The exact form in which this happens is probably a combination of many of the mechanisms proposed.

##### Will these clumps form stars?

The discussion above assumes that all clumps fragment to form stars in order to obtain the IMF from a CMF. However, this might not be the case, in particular for some of the least massive clumps further away from the cloud center. Are these clumps gravitationally bound? Are they hence likely to collapse or are they only transient structures or overdensities triggered by turbulent compressive shocks? These questions have yet to be answered,either observationally or theoretically.

Numerical simulations commonly report that many of the lower-mass cores formed are not gravitationally bound (Klessen, 2001; Clark & Bonnell, 2005; Tilley & Pudritz, 2005). Furthermore, the mass spectra can be understood as due to purely hydrodynamical effects without gravity (e.g. Clark & Bonnell, 2006). As Padoan & Nordlund (2002) remark, the mass spectra resulting from turbulent fragmentation is different from the one that considers collapsing or unstable cores. Only the latter form stars. These calculations show that many of the clumps could be transient structures; indeed a significant fraction of the cores end up re-expanding rather than collapsing (Vázquez-Semadeni et al., 2005; Nakamura & Li, 2005). This implies that fragmentation is not sufficient to trigger star formation. If supersonic turbulence generates the initial density enhancements from which cores develop, then these cores might not necessarily approach hydrostatic equilibrium at any point in their evolution (Ballesteros-Paredes et al., 1999).

Can we distinguish observationally between bound and transient clumps? To estimate if clumps are likely or not to collapse, we need to know more than its mass. A clump’s column density is important as a diagnostic of whether the clump is likely to collapse (Pudritz, 2002). Values higher than  cm are observed for cores with embedded sources. The virial parameter, , where is a mass derived from column density and is the virial mass derived from the cloud radius and velocity dispersion(Bertoldi & McKee, 1992), is known to have values close to 1.0 in star forming clouds (Onishi et al., 1996; Yonekura et al., 2005). However, given the observational uncertainties it is not yet clear if this parameter is a good diagnostic for star forming vs. transient clumps.

##### Preferred spatial scale.

As Larson (1995) found in Taurus, there might exist a preferred scale of star formation at which clustering changes. This length can be related to the Jeans length and the Jeans mass. The surface number density of the clumps as a function of separation can reveal a characteristic spatial scale, marking a transition between clumps and cores. At separations larger than 0.6 pc and smaller than 10 pc, we find the same slope for the power-law fit to the surface density of companions as the protostars in Taurus do at large separations: . This suggests that large separations in stellar systems are determined by the position of their progenitor clumps.

## 5. Conclusions

We made observations of the 1.2-mm dust emission toward the Giant Molecular Cloud NGC 6334 using the SIMBA bolometer at the SEST. The main results and conclusions presented in this paper are summarized as follows.

We find 181 clumps, distributed in an area of square degrees centered on the main filament, which harbors most of the star-forming activity. The clumps range in size from 0.1 to 0.9 pc, with a median of 0.35 pc. This range is similar to that found by Faúndez et al. (2004) and Plume et al. (1997) for clumps in different high mass star forming regions. The clump masses, assuming they are isothermal, range from 3 M to , with a completeness limit of  M (assuming  K).

The Clump Mass Function (CMF) is well fit with a power-law dependence with mass with an index () in the mass range between  M to  M. The slope differs from the stellar IMF slope, indicating that clumps are not direct progenitors of stars. Therefore other processesÐ besides fragmentationÐ must be important in setting up the IMF from the CMF.

We assessed possible effects on the derived slope of the CMF due to changes on the temperature assumptions and due to the contribution of free-free emission from ionized gas to the 1.2-mm emission. Although of the clumps are likely to be significantly warmer than 17 K and are associated with regions of ionized gas, the correction for temperature and free-free emission has little effects on the derived slope.

We investigated possible differences in the value of among different sub-regions of NGC 6334. We find that the slope is significantly shallower toward the central filament (), which contains the most massive clumps and represents the minimum of the gravitational potential in the GMC. As we cover more extended regions, with clumps not actively forming stars, the slope steepens (), revealing that the bulk of the clumps are located in the outer areas of the molecular cloud and that these low-mass clumps predominantly determine the shape of the mass function.

We caution about the power-law fitting procedures to the mass function. The differential CMF is sensitive to bin size and to low-number statistics in the last bin (high-mass end) as well as to completeness limits in the low-mass end. On the other hand, fitting a power-law to the high mass end of the cumulative CMF is incorrect due to its high mass cutoff. Both the low-mass and high-mass ends of the cumulative CMF must be avoided in fitting power-laws.

The spatial analysis performed on the two-dimensional distribution of clumps reveals that they are not distributed randomly. They are concentrated toward the center of the filament, indicating not only a segregation in mass but also a segregation in number which could suggest a possible coalescence of massive clumps towards the gravitational potential minimum. In addition, we study the surface density of companions as a function of separation. This is well fit by a power-law with a similar exponent to the one found for proto-stars in Taurus at large angular separations. This suggests that the position of stars in clusters is determined in the fragmentation and star formation stage rather than after dynamical relaxation.

We thank Simón Casassus and Dieter Nurnberger for helpful comments. D.J.M., D.M., G.G., and D.R. gratefully acknowledge support from the Chilean Centro de Astrofísica FONDAP No. 15010003.

## References

• Altenhoff et al. (1961) Altenhoff, W., Mezger, P.G., Wendker, H., & Westerhout, G., 1961, Veroffentl. Sternwarte, 59, 48
• Alves et al. (2007) Alves, J., Lombardi, M., & Lada, C. J. 2007, A&A, 462, L17
• Ballesteros-Paredes et al. (1999) Ballesteros-Paredes, J., Vázquez-Semadeni, E., & Scalo, J. 1999, ApJ, 515, 286
• Beckert et al. (2000) Beckert, T., Duschl, W. J., & Mezger, P. G. 2000, A&A, 356, 1149
• Bertoldi & McKee (1992) Bertoldi, F., & McKee, C. F. 1992, ApJ, 395, 140
• Beuther & Schilke (2004) Beuther, H. & Schilke, P. 2004, Science, 303, 1167
• Beuther et al. (2006) Beuther, H., Churchwell, E. B., McKee, C. F., & Tan, J. C. 2006, ArXiv Astrophysics e-prints, arXiv:astro-ph/0602012
• Bonnell & Davies (1998) Bonnell, I. A., & Davies, M. B. 1998, MNRAS, 295, 691
• Blitz (1993) Blitz, L. 1993, Giant Molecular Clouds, in Protostars and Planets III, eds. E.H. Levy & J.I. Lunine (Tucson: Univ. of Arizona Press)
• Brooks et al. (2005) Brooks, K. J., Garay, G., Nielbock, M., Smith, N., & Cox, P. 2005, ApJ, 634, 436
• Burton et al. (2000) Burton, M.G., Ashley, M.C.B., Marks, R.D., Schinckel, A.E., Storey, J.W.V., Fowler, A., Merril, M., Sharp, N., Gatley, I., Harper, D.A., Loewenstein, R.F., Mrozek, F., Jackson, J.M., & Kraemer, K.E. 2000 ApJ,542, 359
• Carral et al. (2002) Carral, P., Kurtz, S.E., Rodríguez, L.F., Menten, K., Cantó, J., & Arceo, R. 1978, AJ, 123, 2574
• Chini, Krügel & Wargau (1987) Chini, R., Krügel, E. & Wargau, W. 1987, A&A, 181, 378
• Clark & Bonnell (2006) Clark, P. C., & Bonnell, I. A. 2006, MNRAS, 368, 1787
• Clark & Bonnell (2005) Clark, P. C., & Bonnell, I. A. 2005, MNRAS, 361, 2
• Evans (1999) Evans, N. J., II 1999, ARA&A, 37, 311
• Faúndez et al. (2004) Faúndez, S., Bronfman, L., Garay, G., Chini, R., Nyman, L.-Å., & May, J. 2004, A&A, 426, 97
• Garay et al. (2002) Garay, G., Brooks, K.J., Mardones, D., Norris,R.P., & Burton, M.G. 2002, ApJ, 579, 678
• Gezari (1982) Gezari, D.Y. 1982, ApJ, 259, L29
• Gómez et al. (1993) Gómez, M.,Hartmann, S., Kenyon, S.J.,& Hewett,R. 1993, AJ, 105, 5
• Hunter et al. (2006) Hunter, T. R., Brogan, C. L., Megeath, S. T., Menten, K. M., Beuther, H., & Thorwirth, S. 2006, ApJ, 649, 888
• Jackson & Kraemer (1999) Jackson, J.M. & Kraemer, K. E., 1999, ApJ, 512, 260
• Johnstone et al. (2000a) Johnstone, D., Wilson, C. D., Moriarty-Schieven, G., Giannakopoulou-Creighton, J., & Gregersen, E. 2000a, ApJS, 131, 505
• Johnstone et al. (2000b) Johnstone, D., Wilson, C.D., Moriarty-Schieven, G., Joncas, G., Smith, G., Gregersen, E. & Fich, M. 2000b, ApJ, 545, 327
• Johnstone et al. (2001) Johnstone, D.,Fich, M. Mitchell, G.F., & Moriarty-Schieven, 2001, ApJ, 559, 307
• Johnstone et al. (2006) Johnstone, D., Matthews, H., & Mitchell, G. F. 2006, ApJ, 639, 259
• Kerton et al. (2001) Kerton, C.R., Martin, P.G.,, Johnstone, D., & Ballantyne, D.R. 2001, ApJ, 552, 601
• Kitsionas et al. (1998) Kitsionas, S., Gladwin, P. P., & Whitworth, A. P. 1998, ASP Conf. Ser. 132: Star Formation with the Infrared Space Observatory, 132, 434
• Klessen (2001) Klessen, R. S. 2001, ApJ,556, 837
• Klessen & Burkert (2000) Klessen, R. S., & Burkert, A. 2000, ApJS, 128, 287
• Kramer et al. (1998) Kramer, C., Stutzki, J., Rörig, R.,& Corneliussen, U. 1998, A&A, 329, 249
• Kraemer & Jackson (1999) Kraemer, K. E., & Jackson, J. M. 1999, ApJS, 124, 439
• Kraemer et al. (2000) Kraemer, K.E., Jackson, J.M. 1999, Lane, A.P. & Paglione, T.A.D. ApJ, 542, 946
• Kroupa (2001) Kroupa, P. 2001, MNRAS, 322, 231
• Lada (1999) Lada, E. A. 1999, NATO ASIC Proc. 540: The Origin of Stars and Planetary Systems, C.J. Lada, N.D. Kylafis, Eds. (Kluwer, Dotrecht, Netherlands, 200), p 441
• Larson (1995) Larson, R. 1995, MNRAS, 272, 213
• Larson (2003) Larson, R. B. 2003, ASP Conf. Ser. 287: Galactic Star Formation Across the Stellar Mass Spectrum, J.M. De Buizer abd N.S. van der Blick, (eds), 287, 65
• Li et al. (2006) Li, D., Velusamy, T., Goldsmith, P. F., & Langer, W. D. 2006, ArXiv Astrophysics e-prints, arXiv:astro-ph/0610634
• McBreen et al. (1979) McBreen,B., Fazio, G.G., Stier, M., & Wright, E.L. 1979, ApJ, 232, L183
• McCutcheon et al. (2000) McCutcheon, W.H., Sandell, G. Matthews, H.E., Kuiper, T.B.H., Sutton, E.C., Danchi, W.C. & Sato, T. 2000, MNRAS, 316, 152
• Menzel & Pekeris (1935) Menzel, D. H., & Pekeris, C. L. 1935, MNRAS, 96, 77
• Mookerjea et al. (2005) Mookerjea, B., Kramer, C., Nielbock, M., & Nyman, L.-Å. 2005, A&A, 426, 119
• Moran et al. (1990) Moran, J. M., Greene, B., Rodriguez, L. F., & Backer, D. C. 1990, ApJ, 348, 147
• Motte et al. (1998) Motte, F., Andre, P., & Neri, R. 1998, A&A, 336, 150
• Motte et al. (2003) Motte, F., Schilke, P., & Lis, D.C., 2003, ApJ, 582, 277
• Muñoz (2006) Muñoz, D.J.  2006, M.Sc. Thesis, Universidad de Chile
• Nakamura & Li (2005) Nakamura, F., & Li, Z.-Y. 2005, ApJ, 631, 411
• Neckel (1978) Neckel, T. 1978, A&A, 69, 51
• Nozawa et al. (1991) Nozawa, S., Mizuno, A., Teshima, Y., Ogawa, H. & Fukui, Y. 1991, ApJS, 77,647
• Onishi et al. (1996) Onishi, t., Mizuno, A., Kawamura, A., Ogawa, H., & Fukui, Y. 1996, ApJ, 465, 815
• Ossenkopf & Henning (1994) Ossenkopf, V., & Henning, T. 1994, A&A, 291, 943
• Padoan & Nordlund (2002) Padoan, P., & Nordlund, Å 2002, ApJ, 576, 870
• Plume et al. (1997) Plume, R., Jaffe, D. T., Evans, N. J., II, Martin-Pintado, J., & Gomez-Gonzalez, J. 1997, ApJ, 476, 730
• Pudritz (2002) Pudritz, R. E. 2002, Science, 295, 68
• Reid & Wilson (2005) Reid, M. A., & Wilson, C. D. 2005, ApJ, 625, 891
• Reid & Wilson (2006a) Reid, M. A., & Wilson, C. D. 2006a, ApJ, 644, 990
• Reid & Wilson (2006b) Reid, M. A., & Wilson, C. D. 2006b, ApJ, 650, 970
• Rodríguez, Cantó & Moran (1982) Rodríguez, L.F., Cantó, J., & Moran, J.M. 1982, ApJ, 255, 103
• Rosolowsky (2005) Rosolowsky, E. 2005, PASP, 117, 1403
• Rybicki & Lightman (1979) Rybicki, G. B., & Lightman, A. P. 1979, Radiative Processes in Astrophysics, New York, Wiley-Interscience, 1979. 393 p.,
• Salpeter (1955) Salpeter, E. E. 1955, ApJ, 121, 161
• Sandell (2000) Sandell, G. 2000, A&A, 358,242
• Sanders et al. (1985) Sanders, D. B., Scoville, N. Z., & Solomon, P. M. 1985, ApJ, 289, 373
• Scalo (1998) Scalo, J. 1998, ASP Conf. Ser. 142: The Stellar Initial Mass Function (38th Herstmonceux Conference), G. Gilmore and D. Howell (eds),142, 201
• Schneider & Brooks (2004) Schneider, N., & Brooks, K. 2004, Publications of the Astronomical Society of Australia, 21, 290
• Silverman (1986) Silverman, B. W. 1986, Monographs on Statistics and Applied Probability, London: Chapman and Hall, 1986
• Solomon et al. (1987) Solomon, P. M., Rivolo, A. R., Barrett, J., & Yahil, A. 1987, ApJ, 319, 730
• Sommerfeld (1953) Sommerfeld, A.J.F., Atombau und Spektrallinien, Vol.2. Ungar, New York, 1953.
• Straw & Hyland (1989) Straw, S. M., & Hyland, A. R. 1989, ApJ, 340, 318
• Straw et al. (1989) Straw, S. M., Hyland, A. R., & McGregor, P. J. 1989, ApJS, 69, 99
• Stutzki & Gusten (1990) Stutzki, J., & Güsten, R. 1990, ApJ, 356, 513
• Testi & Sargent (1998) Testi, L., & Sargent, A. I. 1998, ApJ, 508, L91
• Tilley & Pudritz (2005) Tilley, D. A., & Pudritz, R. E. 2005, ArXiv Astrophysics e-prints, arXiv:astro-ph/0508562
• Tothill et al. (2002) Tothill,N.F.H., White, G.J., Matthews, H.E., McCutcheon, W.H., McCaughrean, M.J., & Kenworthy, M.A. 2002, ApJ, 580, 285
• Vázquez-Semadeni et al. (2005) Vázquez-Semadeni, E., Kim, J., Shadmehri, M., & Ballesteros-Paredes, J. 2005, ApJ, 618, 344
• Ward-Thompson et al. (2006) Ward-Thompson, D., Andre, P., Crutcher, R., Johnstone, D., Onishi, T., & Wilson, C. 2006, ArXiv Astrophysics e-prints, arXiv:astro-ph/0603474
• Williams, Blitz & McKee (2000) Williams, J.P., Blitz,L. & McKee, C.F. 2000, in Protostars and Planets IV, ed. V. Mannings, A.P. Boss, & S.S. Russell (Tucson: Univ. of Arizona), 97
• Williams, de Geus & Blitz (1994) Williams, J.P., de Geus ,E.J., & Blitz, 1994, ApJ, 428, 693
• Yonekura et al. (2005) Yonekura, Y., Asayama, S., Kimura, K., Ogawa, H., Kanai, Y., Yamaguchi, N., Barnes, P. J., & Fukui, Y. 2005, ApJ, 634, 476