A&A 417, 115-133 (2004)
DOI: 10.1051/0004-6361:20031733

The circumstellar environments of high-mass protostellar objects

I. Submillimetre continuum emission^,

S. J. Williams¹ - G. A. Fuller¹ - T. K. Sridharan²

1 - Department of Physics, UMIST, PO Box 88, Manchester, M60 1QD, UK
2 - Harvard-Smithsonian Center for Astrophysics, 60 Garden Street, MS 78, Cambridge, MA 02138, USA

Received 27 March 2003 / Accepted 3 December 2003

Abstract
We present maps of the 850 $\mu$m and 450 $\mu$m continuum emission seen towards a sample of 68 high-mass protostellar candidates with luminosities ranging from 10 $^{2.5}~L_{\odot}$ to ${\sim} 10^{5}$ $L_{\odot}$. Most of these candidate high-mass stars are in the earliest stages of evolution, and have not yet developed an ultra-compact HII region. We observe a variety of continuum emission morphologies, from compact symmetric sources through to multiple cores embedded in long filaments of emission. We find on average there is a 65% probability of an IRAS point-source having a companion detection at submillimetre wavelengths. The ratio of integrated flux to peak flux for our detections shows no strong dependence on distance, suggesting the emission we have observed is primarily from scale-free envelopes with power-law density structures. Assuming a near kinematic distance projection, the clumps we detect vary in mass from ${\sim}1~M_{\odot}$ to over 1000 $M_{\odot}$, with a mean clump mass of 330 $M_{\odot}$, column density of $9\times10^{23}$ cm^-2 and diameter of ${\sim}0.6$ pc. The high luminosity and low mass of the smallest clumps suggests they are accompanied by a minimal number of stellar companions, while the most massive clumps may be examples of young protogroups and protoclusters. We measure the spectral index of the dust emission ($\alpha$) and the spectral index of the dust grain opacity (${\beta }$) towards each object, finding clumps with morphologies suggestive of strong temperature gradients, and of grain growth in their dense inner regions. We find a mean value for ${\beta }$ of 0.9, significantly smaller than observed towards UCHII regions.

Key words: stars: formation - stars: circumstellar matter - ISM: clouds - ISM: dust, extinction

1 Introduction

Our understanding of the processes involved in low-mass star formation has matured steadily over the last forty years, and the concepts of gravity-driven collapse and accretion-driven evolution appear to consistently explain how a low-mass pre-main-sequence star can form from a cloud core. However, our knowledge of how high-mass stars form has remained limited, primarily due to a lack of candidate high-mass protostars to study. At the upper reaches of the initial mass function, high-mass stars are statistically rare, and coupled with a characteristically short evolutionary timescale, there have been few chances to observe massive stars at the instance of formation.

Until recently, most young high-mass stars were first identified through the detection of a radio-bright ultra-compact HII (UCHII) region, considered a beacon pointing to the presence of a young high-mass star. As high-mass protostars increase in mass and luminosity, they emit an ever larger number of high energy UV photons which ionize the protostar's immediate surroundings, hence the small, compact nature of a UCHII region is usually considered evidence of the youthful status of the driving source (although debate continues about the exact timescale of the UCHII stage; for a review, see Kurtz et al. 2000).

Unfortunately, a powerful protostar and UCHII region soon act to disrupt and confuse their surroundings, so the initial conditions of the natal cloud and the mechanisms that led to the formation of the massive protostar cannot be unambiguously reconstructed. As a result, many questions about high-mass protostars remain - in particular, do they form via processes similar to their low-mass counterparts? To address the mechanisms that create and shape high-mass stars, we must observe before they have formed a UCHII region, during the initial collapse of the star-forming core.

1.1 The search for precursors of UCHII regions

Clumps bearing the youngest high-mass protostars have proved particularly difficult to find, even though their identifying characteristics have been known for almost 25 years. For example, Habing & Israel (1979) predicted that candidate high-mass protostars should be founded embedded in dense environments, and although highly luminous, they should not at this stage be associated with HII regions. However, despite knowledge of these distinguishing features, it was not until the last decade that samples of candidate high-mass protostellar objects (HMPOs) were finally compiled.

Observations of these preliminary samples of HMPOs have allowed the first glimpses of high-mass protostars in their earliest evolutionary states: a typical core not yet associated with an ionised region is found to be larger, more massive, and more turbulent than a UCHII-class protostar, with a typical diameter of around 0.5-1.0 pc and a mass that may range from a few tens to a few thousand solar masses (Brand et al. 2001; Beuther et al. 2002a). They are cooler, with typical dust temperature averaging around 30-40 K (Sridharan et al. 2002; Molinari et al. 2000), while the dust opacity usually has a spectral index of around 2, suggestive of silicate dust grains (Molinari et al. 2000). Self-absorption profiles towards a number of candidate HMPOs suggest infall may be an important part of the formation mechanism (e.g. Brand et al. 2001; Fuller et al. 2004), while outflow observations suggest that accretion is a significant process (Zhang et al. 2001; Beuther et al. 2002b; Molinari et al. 2002). Water maser emission has also been detected towards candidate protostars, a feature thought to be missing from more evolved sources (e.g. Palla et al. 1993; Sridharan et al. 2002). Despite these advances, there is still much to be learned about the pre-UCHII stage of high-mass star formation, and there remains a need for additional candidates and further observations.

1.1.1 A new sample of high-mass protostars

Recently, Sridharan et al. (2002; SBSMW hereinafter) identified a new sample of HMPOs. The SBSMW sample is a flux-limited sample, constructed through an analysis of the IRAS point-source catalogue: as young high-mass stars are usually associated with UCHII regions, SBSMW began by initially selecting bright IRAS detections (S₆₀>90 Jy and S₁₀₀>500 Jy) with colour characteristics similar to known UCHII regions (they conform to the Wood & Churchwell (1989) FIR colour criteria that selects UCHII regions, and they also satisfy the additional Ramesh & Sridharan (1997) criteria). Candidate sources detected in Galaxy-wide 5GHz continuum surveys were removed, thereby rejecting sources already sufficiently evolved to have ionized their surroundings. As a final requirement, successful candidates must also be associated with CS(2-1) emission, an indicator of dense molecular gas (Bronfman et al. 1996).

In total, sixty-nine IRAS point sources satisfied these cumulative criteria, identifying these sources as potentially among the most massive and deeply embedded pre-UCHII protostars in our Galaxy. The SBSMW sample has been studied in detail over the last few years, and their status as high-mass candidate protostars has been supported through observations of 1.2 mm and 3.6 cm continuum emission (Beuther et al. 2002a; SBSMW), molecular line emission (CS, CO and NH₃), and H₂O and CH₃OH maser transitions (SBSMW) towards the sources.

This paper presents the results of a new set of submillimetre (submm) observations of the SBSMW sample of candidate high-mass protostars. All but one (IRAS 18517+0437) of the SBSMW candidate HMPOs were observed. An additional source, IRAS 18449-0158, was observed but this source does not satisfy the SBSMW criteria and is not included in any analysis. Our observations are detailed in Sect. 2, with maps of the reduced data found in Sect. 3. We measure the multiplicity of the detections in Sect. 3.1, commenting on the position and morphologies of the sample in Sect. 3.4. We analyse the dust optical depth in Sect. 4, and use the spectral index of the emission to investigate the nature of the dust in Sect. 4.1. We calculate the mass characteristics of our sample in Sect. 5, and consider the implications of the cumulative mass spectrum in Sect. 5.1. After a brief discussion and comparison of our results with the IRAM 1.2 mm continuum observations of Beuther et al. (2002a) in Sect. 6.1, we conclude in Sect. 7 with a summary of our results.

This paper presents the first half of our study and analysis of the dust emission; the companion to this paper presents the results of radiative transfer modelling of the clumps (Williams et al. 2004).

Figure 1: Jiggle-maps of Uranus measured at 850 $\mu$m (left-hand image) and at 450 $\mu$m (right-hand image). Contours were chosen to highlight the JCMT beam structure, and are drawn every 0.8 mag below 65 Jy for the 850 $\mu$m image and every 0.8 mag below 150 Jy for the 450 $\mu$m map. The primary error beam can clearly be seen as a ring encircling the main beam, containing $\sim$9% and $\sim$25% of the total flux at 850 $\mu$m and 450 $\mu$m respectively.

Open with DEXTER

   
2 Observations and data reduction

The sample of HMPOs was observed at 850 $\mu$m and 450 $\mu$m between March 2000 and June 2000 using the Submillimetre Common-User Bolometer Array (SCUBA) on the James Clerk Maxwell Telescope (Holland et al. 1999). The SBSMW sample target co-ordinates and date(s) of observation are listed in Table 1. The SCUBA array covers a hexagonal 2.5' field of view with 97 and 37 pixels at 450 $\mu$m and 850 $\mu$m respectively. Maps were formed simultaneously at both frequencies using the "jiggle'' mode, in which the telescope beam is moved around a 64-point pattern by the secondary mirror in order to fully sample the sky.

The data were reduced using the SCUBA User Reduction Facility (SURF; Jenness & Lightfoot 1998). Correlated sky noise was removed using the REMSKY routine, based on the signal from a hand-picked sample of bolometers considered free from source emission. Maps were extinction calibrated from skydips and flux calibrated in terms of Jy beam^-1 from maps of Uranus, IRAS 16293-2422, CRL 618, and CRL 2688, following the procedures defined by Sandell et al. (2001).

Zenith opacities at 225 GHz ranged from 0.05-0.12 during the observations, but usually averaged around 0.10. Telescope pointing was calibrated many times during each observing run, and telescope drift was minimal, requiring very small ( $\sigma=1.7''$) corrections overall. We measured the JCMT beam size from observations of Uranus, finding a full-width half-maximum of $\theta_{\rm beam}=8.0''$ at 450 $\mu$m and $\theta_{\rm beam}=14.4''$ at 850 $\mu$m. An average 1-$\sigma$ rms noise level of 0.03 Jy beam^-1 and 0.69 Jy beam^-1 was found at 850 $\mu$m and 450 $\mu$m respectively. The rms noise level measured in each jigglemap is listed in Table 2 as the uncertainty in the peak flux.

Clumps were identified using object detection routines in the software package GAIA (Chipperfield & Draper 2001). We define a positive detection as a group of pixels subtending at least the area of the JCMT beam with emission above a 3-$\sigma$ level, where $\sigma$ is the rms noise level of the jigglemap. The validity of each detection was also confirmed manually. Clumps not quite bright enough to be automatically detected were examined, and if deemed worthy of inclusion, added to the list of detections. These lower sigma detections are labelled by a note in Table 2.

We list the peak flux per beam and the integrated flux for each detection. The peak flux per beam gives the peak flux level averaged in a 14.4'' beam for 850 $\mu$m maps and in an 8.0'' beam for 450 $\mu$m maps, while the integrated flux of a detection measures the total flux inside an isophote tracing the 3-$\sigma$ rms noise level around the detection. We quote the position of each detection as the location of peak emission, not as the centroid of the 3-$\sigma$ isophote.

Calibrating the integrated flux of a detection required additional consideration, as the JCMT beam structure is complex (Fig. 1), so the number of detector counts recovered within an aperture is also a function of aperture size. We quantified the extent of this relationship using maps of Uranus (which we consider a point source), calibrating detector counts recovered inside circular apertures of increasing radius. We did not include the small number of non-planetary flux calibrators in the calibration of integrated flux. The resulting function measures increasing counts with aperture size, asymptotically reaching maximum counts once the aperture has expanded to encompass the JCMT beam and its primary error beam. For each detection, we then converted n counts recovered inside an isophote of area A to Janskys by multiplying n by the counts-to-Jy conversion factor derived from a circular aperture of equivalent area.

Our observations were performed on seven nights over a period of three months. Despite the protracted nature of our observations, a comparison of the counts-to-Jy conversion factor calculated for each night showed it usually remained consistent with the published JCMT response. Where the conversion factor appeared inconsistent and no other recent flux calibrator maps were available, we assumed a conversion factor equal to the mean value for our run. A comparison with the independent 1.2 mm continuum observations of Beuther et al. (2002a) shows the data to be consistently calibrated (Sect. 6.1), and we estimate the absolute flux uncertainty to be $\pm$10% at 850 $\mu$m and $\pm$30% at 450 $\mu$m.

2.1 Telescope and reduction artefacts

We masked a small number of consistently noisy 850 $\mu$m bolometers during data reduction; these bad bolometers usually fell in a region of background sky or faint extended emission. Unfortunately, the secondary source seen towards IRAS 05358+3543 fell on a noisy 850 $\mu$m bolometer, but we consider the companion source an important feature of the jigglemap and so leave the bolometer unmasked. In addition, the map of IRAS 18553+0414 forms an isolated case that suffers from an unusually large number of bad 850 $\mu$m bolometers; we still include this data as the 450 $\mu$m map reveals the majority of flux has been recovered by good bolometers.

Our observations used a 120'' chop to sky to measure and remove the background emission. However, in crowded regions, the 120'' chop-throw sometimes points the telescope towards an occupied region of sky rather than an empty field. When this occurs, emission from objects in the sky reference beam is subtracted from the target field emission, resulting in negative images of clumps seen towards the reference position superimposed onto the final map. Some of our maps contain these artefacts, which are usually seen away from regions of interest (e.g. IRAS 18151-1208, IRAS 18431-0312), but chopping onto emission altered the map of IRAS 18454-0158 to such an extent that no reliable measurement was possible, and this source was removed from our analysis.

We occasionally observed additional jigglemaps offset from the target position to map fields with emission continuing outside the ${\sim} 120''$ SCUBA field of view. These additional maps were calibrated as individual jigglemaps before they were combined into a mosaic, weighting the contribution of each map to intersecting areas by 1/ $\sigma^{2}$, where $\sigma$ is the rms noise level in the map. Detections within the mosaic are still defined as clumps with emission above a 3-$\sigma$ limit over an area the size of the JCMT beam, but using the rms local to the section of mosaic being measured.

   
3 Results

Submillimetre emission at 850 $\mu$m was detected towards all the IRAS sources in our sample, although not all sources were bright enough to be detected above the increased background emission at 450 $\mu$m. When sources were bright enough to be detected at 450 $\mu$m, the increased resolution of the 450 $\mu$m observations sometimes resolved additional peaks within the area of a single 850 $\mu$m detection. A complete presentation of the reduced jigglemaps, calibrated in Jy beam^-1, alongside maps of $\alpha$, the spectral index of the dust emission (detailed in Sect. 4.1), can be found in Fig. 2. The position and flux measured for each detection is listed in Table 2.

Figure 2: Maps of S850, S450, and the $\alpha$ distribution in the left-hand, center, and right-hand panels respectively. Greyscale limits are chosen to emphasise flux levels between -2 $\sigma$ and +7 $\sigma$, where $\sigma$ is the rms noise level in the map. Contours trace the intensity in units of Jy beam-1, using the step-size listed in the bottom left-hand corner of each map. The first contour is drawn at the first step above zero Jy beam-1 unless the index is marked with an asterix: this signifies that an additional contour is plotted at half a step above zero Jy beam-1. Triangular symbols plot the location of MSX point sources in each field of view. The spectral index of the dust emission ($\alpha$) is plotted in the right-hand greyscale maps: $\alpha$ is masked outside the boundary of the first 850 $\mu$m and 450 $\mu$m contours. Contours on the $\alpha$ map directly mirror those drawn on the 850 $\mu$m submillimetre emission map.

Open with DEXTER

The target sources IRAS 19266+1745 and IRAS 18553+0414 displayed a gas+dust mass incompatible with the luminosity of the driving protostar, unless these sources are projected to the far kinematic distance (Sect. 6.2). Therefore, we reject the near kinematic distance for these objects and consider them resolved to the far kinematic distance for all subsequent analysis.

   
3.1 Companion clump fraction

The majority of candidate HMPOs were detected as companionless, compact and approximately spherically symmetric submm clumps (e.g. IRAS 05553+1631), although a significant number exhibited submm nebulosity (e.g. IRAS 18566+0408), appeared in filaments of emission (e.g. IRAS 18437-0216) and/or existed with multiple additional detections within the field of view (e.g. IRAS 23545+6508). While our sample was explicitly constructed to consist of isolated, companionless candidates away from sources of confusion, we found only 38 of the 68 target IRAS fields contained a single, companionless clump. The remaining IRAS fields contained more than one submm clump, usually two detections, with the mosaic map towards IRAS 18089-1732 containing the most companions, where 5 separate clumps were resolved. This demonstrates the difficulty in locating truly isolated candidate HMPOs. We can characterise the multiplicity of our detections by calculating the companion clump fraction (CCF), expressed by the formula

$$CCF = \frac{B+2T+3Q+4P}{S+B+T+Q+P},$$ (1)

where S, B, T, Q and P are the number of single, binary, triple, quadruple and quintuplet clumps in our sample. If all clumps in our sample were solitary clumps, the CCF would be 0.0, while if all clumps had one companion the CCF would be 1.0.

The CCF for our sample is $0.65\pm0.1$, where the uncertainty comes from $\sqrt{N}$ counting statistics. In reality, the absolute value of the CCF and the quoted uncertainty are both lower limits, as they are calculated assuming we have detected, and are uniformly sensitive, to all companions. This is not the case, as our limited angular resolution precludes the detection of companions closer than around a beamsize, plus the finite field of view means companions of greater than ${\sim}60''$ separation (assuming a clump central in the jigglemap) will not be detected. Additionally, when coupled with the large difference in projected distance (the most distant sources being more than 15 times further away than the closest sources), our companion mass sensitivity also bears a dependence on distance.  

We examined the effect of different distance projections by sorting our candidate HMPOs into four bins, containing sources <2 kpc, 2 to 4 kpc, 4 to 8 kpc and >8 kpc distant, respectively. The CCF of these subsamples remains remarkably consistent, each section in agreement with the full sample CCF within the uncertainty limits. This is true regardless of whether distance-unresolved sources are projected to the near or far kinematic distance, with the exception of the <2 kpc bin projected to the far kinematic distance, and suggests clumps have a similar number of companions over a wide range of distance scales.

The expression of multiplicity given in Eq. (1) is usually used as a diagnostic of more evolved stars, in particular to quantify the number of companions a low-mass star is born with (e.g. Beck et al. 2003; Patience et al. 2002), whereas in this study the CCF can be interpreted as the likelihood of finding additional potentially star-forming clumps when observing Galactic HMPOs identified by a similar flux-limited criteria. While it remains difficult to constrain the statistics of such a disparate sample, the CCF does emphasise that most clumps do not form in isolation, and that a single IRAS detection is usually resolved into several submm clumps. New samples of protostars comparable to our sample are hard to compile, but the strong likelihood of detecting additional clumps in the locality of our sample suggests that wide-field surveys towards existing high-mass protostars may also be a productive way of locating new protostellar candidates.

3.2 Clump morphology statistics

Forming an unbiased statistical analysis of the morphological features is difficult, as distinguishing features are mainly found in the appearance of low-level extended emission. As such, measuring the FWHM of detections is of limited use, as it is not sensitive to the faint emission features we wish to characterise. Instead, we formed a simple statistic that indicates how much mass lies outside the central beam by measuring Y, the ratio of integrated flux ( $S_{\rm int}$) to peak flux ( $S_{\rm peak}$) at 850 $\mu$m. We measure at 850 $\mu$m because of higher signal-to-noise than in the corresponding 450 $\mu$m maps. For a point source, Y equals unity. The value of Y for each detection is listed in Col. 9 of Table 2.

In Fig. 3 we plot the Y distribution for our submm detections, dividing our sample into two groups: subsample A, containing detections with a high confidence of being solitary cores (having just one detection within an IRAS, MSX, SCUBA and IRAM field of view) and subsample B, the remainder. While the solitary detections in sample A do not display the extended distribution tail seen in subsample B, we see both groups peak at an intensity ratio of ${\simeq} 3$, which in light of the factor 15 range in distance suggests that the envelope structures may be scale-free.

	Figure 3: A histogram of the $S_{\rm int}/S_{\rm peak}$ flux ratio for "confirmed'' isolated sources (sample A: shaded bars) and sources potentially with companions (sample B: outline).
Open with DEXTER

The large scatter in Fig. 4 means we do not find any strong correlation with distance at the near distance projection (which is the most likely projection for our sample), although considering the distance-resolved sources alone does reveal a trend beyond d=4 kpc, where $Y \propto d^{-0.8}$. This fall-off does not reveal a physical change, but reflects the diminishing level of integrated emission as the 3-$\sigma$isophote encloses less of the envelope for more distant sources.

We must qualify a number of uncertainties that could affect the distribution of Y, not least our variable sensitivity to additional embedded sources. While the large-scale envelope structures appear to be scale-free, on the small scale there are indications that we are still undersampling the number of companions separated by less than a beam width. A number of apparently single detections at 850 $\mu$m are barely resolved as multiple sources at 450 $\mu$m (e.g. IRAS 05490+2658); a reminder that further clustering on size scales less than a 450 $\mu$m beam width may also be present. Overall, it is inevitable that with limited resolution we misclassify some multiple cores as solitary detections, a point demonstrated in Fig. 5, where all sources (bar one) with 850 $\mu$m flux ratios larger than 7 are resolved as multiple detections at 450 $\mu$m. IRAS 22551+6221 provides the most visible demonstration of this effect, where the high 850 $\mu$m flux ratio arises from the inclusion of flux from a bright neighbouring source that is only fully resolved at 450 $\mu$m.

We also tend to overestimate the flux of multiple detections, as the elliptical apertures used for photometry could include emission shared with a companion source. Although the intersection of apertures around adjacent components was minimised where possible, it remains a potential cause of uncertainty. Finally, the ratio for extended sources is likely to be a lower limit, as emission from a large, extended envelopes is more likely to project emission onto a noisy bolometer, and flux incident on these noisy bolometers is masked during jigglemap reduction. As a result, the quoted integrated emission is a lower limit, and the flux ratio is underestimated.

To conclude, while these concerns affect the quantative results, qualitatively we still observe that a significant fraction of the total mass lies outside the central "core'' at this stage of evolution.

3.3 Surface density

Stars generally form in association with other stars, and the spatial distribution of these groups of stars can provide information on how the natal molecular cloud fragmented. One way to probe the distribution of sources is to use the mean surface density of companion sources (MSDC), a method which has been successfully used to probe the transition from the formation of binary stars to star clusters in low mass star forming regions (e.g. Gomez et al. 1993; Larson 1995; Nakajima et al. 1998).

Figure 4: A plot of $Y=S_{\rm int}/S_{\rm peak}$, the ratio of 850 $\mu$m integrated flux to peak 850 $\mu$m flux measured in a 14.4'' beam, against the kinematic distance of each detection. Distance resolved sources are plotted by filled circles, while distance unresolved sources are projected to near and far kinematic distances and plotted with open rectangles and triangles. The curve in the distance-resolved plot displays a power law of the form $S_{\rm int}/S_{\rm peak} \propto d^{-0.8}$.

Open with DEXTER

The MSDC of our detections is shown in Fig. 6. The MSDC was calculated by measuring the linear separation r of each detection to its companions. The separation of each companion pair was binned into annuli of separation r to $r+\delta r$. The number of pairs N within each annulus was then divided by the area of the annulus and the total number of sources N_* to give the MSDC, $\Sigma(r)$, as

$$\Sigma(r) = N / (2\pi r \delta r N_{*}).$$ (2)

Above 0.4 pc, the MSDC of our detections has a measured gradient $\gamma$ of -1.7, which is roughly halfway between the power-law indices of binary pre-main-sequence populations (where $\gamma \sim -0.5$; Nakajima et al. 1998) and that of more distant companions (where $\gamma\sim-2.2$; Nakajima et al. 1998). The MSDC we observe appears to turn over below 0.4 pc, but the validity of this turnover is questionable as we have few measurements in this region. The clump MSDC break is also quite distinct from the stellar MSDC power-law break found at ${\sim} 0.04$ pc (Gomez et al. 1993; Larson 1995; Nakajima et al. 1998), below which the MSDC steepens from the inclusion of close binaries, as the stellar MSDC breakpoint occurs on much smaller scales than are detectable by our survey.

	Figure 5: A comparison of the peak flux and integrated flux of each 850 $\mu$m detection. The solid line traces the $Y=S_{\rm int}/S_{\rm peak}=3$ distribution peak found in Fig. 3, while the dotted line traces the ratio Y=7.
Open with DEXTER

To identify the space density distribution consistent with the clump MSDC, we modelled a number of systems, each containing a collection of 10³ sources randomly distributed according to a number of power-law space density distributions. The projected MSDC of these simulations are also shown in Fig. 6. For separations above 0.4 pc, the MSDC power-law slope of -1.7 most closely corresponds to a space density distribution with the number density per unit volume $N_{\rm vol}(r) \propto r^{-0.75}$.

Figure 6: MSDC for single power-law distribution models (light dotted lines), with $p(N<r)\propto r^{\gamma }$ from $\gamma =-2$ (uppermost dotted line) to $\gamma =3.5$ (bottommost dotted line) with $\delta \gamma =0.5$. The square symbols represent the observed MCSD (multiplied by 100), with a thick dashed line plotting the line of best fit for the observed MSDC above the break-point. The best fit power-law distribution ( $\gamma =-0.75$) is plotted by a thick black line.

Open with DEXTER

We must consider the MSDC statistic cautiously, for the MSDC is constructed from observations of widely different companion seperation sensitivities and companion flux sensitivities, in a similar vein to the CCF (Sect. 3.1). Complications arise from the wide range of projected distances to our HMPO candidates: specifically, our observations are sensitive to angular separations from around a beamwidth up to the upper limit of a 120'' field of view, but as most HMPO candidates in our sample are less than 4 kpc away, there are only a very small range of uniformly sampled linear separations for our candidates. Towards more distant IRAS fields, we are uniformly sensitive to a larger range of linear separations, but we then suffer from fewer measurements and from a reduced sensitivity to close companions.

We examined the significance of variable sensitivity using a procedure similar to that used for the CCF (Sect. 3.1), comparing the MSDC slope of groups of sources with similar distances, finding the slope of each MSDC segment agrees with the MSDC of the whole sample within the uncertainty limits. We suggest the MSDC as calculated provides at least a basic estimate of the clustering properties of these clumps. Ultimately, the scarcity of high-mass protostars means there will always be a large range of distances in samples of HMPOs, and we may never be able to construct a set of uniformly sampled observations to the extent possible with low-mass protostars.

   
3.4 Clump positions

For each candidate HMPO, the telescope was pointed so that the IRAS point source was central in the field of view, hence the location of the IRAS source in each map is given implicitly by the map center. We indicate the position of neighbouring MSX sources in the field of view with triangular symbols in Fig. 2, finding that many IRAS sources and MSX detections are roughly coincident with the submm clumps, but in agreement with SBSMW we also find that some submm detections and IR detections cannot be coincident within the positional uncertainties of our survey and the IRAS/MSX surveys. For example, towards IRAS 23545+6508 there are two MSX point sources found within $\sim$30'' of the two submm clumps, an offset greater than the expected absolute uncertainty, plus the MSX point sources also have smaller relative separation between components than the corresponding SCUBA detections. Sources with large SCUBA/MSX offsets do not appear to have ${\beta }$ or S₈₅₀ characteristics different to more coincident detections.

Considering each clump as a protocluster may explain displaced IRAS/MSX and submm detections, as additional stars embedded in the less dense, more transparent outer reaches of the envelope will not encounter the same degree of opacity, providing a mechanism for shorter wavelength photons to pass. To examine this possiblity and resolve whether MSX and SCUBA detections trace the same body of material will require further high-resolution IR observations.

   
4 Dust optical depth

The dust optical depth $\tau$ can be determined from the flux density using the expression

$$\tau_{\nu} = -\ln \left[1-\frac{S_{\nu}}{B_{\nu}(T)\Omega}\right],$$ (3)

where $B_{\nu}(T)$ is the Planck function at temperature T, and $\Omega$ is the solid angle subtended by the telescope beam, all measured at frequency $\nu$. This relationship assumes the emission fills the telescope beam, which may not be the case, so the derived value of $\tau_{\nu}$ must be considered the beam-averaged optical depth.

To calculate the dust optical depth for our detections using Eq. (3) we assumed dust temperatures equal to the SBSMW cold-component dust temperatures. In the SBSMW study, the spectral energy distribution (SED) of each IRAS source was successfully modelled as a composite of two greybodies: one greybody representing a cold dust component, accounting for the >60 $\mu$m flux, while a separate hot dust component contributes the majority of near-IR flux. As SCUBA is only sensitive to emission from the cold dust greybody, we set T equal to the temperature of the cold component ( $T_{\rm cd}$) as given in Table 1 of SBSMW.  

The beam-averaged optical depth and other parameters derived from the flux density are listed in Table 3. While our sources consist of very dense clumps, the beam-averaged 450 $\mu$m and 850 $\mu$m optical depths show that they are usually optically thin at submm wavelengths. The optical depth at 850 $\mu$m spans almost two orders of magnitude, from $10^{-3} \le \tau_{850} \le 10^{-1}$. At 450 $\mu$m, $\tau_{450}$ is found within the range $2\times10^{-1} \le \tau_{450} \le 10^{-2}$ with three exceptions: detection #8 (IRAS 18089-1732) appears optically thick at 450 $\mu$m with $\tau_{450}=1.1$ but this is by far the brightest detection of our survey. Two other detections (#16: IRAS 18182-1433 and #19: IRAS 18264-1152) have high S₄₅₀, leading to higher  $\tau_{450}$ than the majority of the detections, but they remain with $\tau_{450}<1$.

Table 3: Derived parameters of the submm detections resolved by this survey. The mass of each clump is calculated from the 850 $\mu$m integrated flux value, assuming the dust grains have thin ice mantles, and using a 100:1 gas-to-dust ratio. The column density refers to the total gas column density (i.e. n(H+H₂)).

We have modelled the submm emission seen towards our sample using a one-dimensional radiative-transfer code (details and results can be found in the companion to this paper: Williams et al. 2004). These models assume a fixed dust grain chemical composition, with a silicon-to-graphite ratio half that of the interstellar medium (Mathis et al. 1977) and a standard Draine & Lee (1984) dust grain size distribution. We used these model dust grains to predict the optical depth of the cores as a function of wavelength, which when scaled to match the observed 850 $\mu$m and 450 $\mu$m optical depths gives an estimate of the wavelength at which the submm detections become optically thick (listed in Col. 7 of Table 3). We find that for the average clump, $\tau>1$ for wavelengths shorter than $\sim$90 $\mu$m.

   
4.1 Spectral index of the dust emission

As noted in Sect. 4, the dust emission at submm wavelengths can be well represented by a greybody, with intensity varying smoothly as a function of frequency. Having measured the intensity of dust emission at two submm frequencies, we can characterise the SED using $\alpha$, the spectral index of the dust emission. This is defined

$$\alpha = \frac{\ln(S_{1}/S_{2})}{\ln(\nu_{1}/\nu_{2})},$$ (4)

where S₁ and S₂ are flux densities at frequencies $\nu_{1}$ and $\nu_{2}$ respectively. By dividing calibrated images measured at frequencies $\nu_{1}$ and $\nu_{2}$, we can determine the spatial distribution of $\alpha$and examine its relationship with the intensity of emission and density of gas and dust.

However, first we must take into account the different JCMT response and beam patterns at 850 $\mu$m and 450 $\mu$m, for as seen in Fig. 1, more flux lies outside the main beam at 450 $\mu$m compared to the response at 850 $\mu$m. We accounted for these differences by following the procedures defined in Hogerheijde & Sandell (2000), normalizing the images to a common response before finally determining $\alpha$ as the ratio of images. In detail, we described the JCMT beam at each wavelength as a superposition of three Gaussians, the parameters of which were found by a fit to the azimuthal average response to Uranus. The amplitude and FWHM of these components are listed in Table 4. We then deconvolved the 850 $\mu$m and 450 $\mu$m SCUBA images with the corresponding beam model, smoothing the deconvolved images with a single Gaussian to achieve a final, uniform, spatial resolution of 15.0'' before forming $\alpha$ as given in Eq. (4).

	Figure 7: A plot of the typical cross-sections seen for peaked $\alpha$ morphologies (IRAS 05358+3543; left-hand plots) and negative-dip morphologies (IRAS 05490+2658; right-hand plot). The upper row displays the orientation of the cross-sections, while the bottom row displays the $\alpha$ index measured along the labelled cut, using a pixel scale of 3'' per pixel.
Open with DEXTER

The spatial distribution of $\alpha$ can be seen to the right of the submm emission maps in Fig. 2; $\alpha$ is blanked for emission lying outside the first 850 $\mu$m and 450 $\mu$m contours, so only the intersecting area of the maps with good signal to noise is both displayed and analysed. In Table 3 we list the spectral index both in terms of the index at the position of peak 850 $\mu$m flux ( $\alpha _{\rm peak}$; Col. 13) and the mean value of the $\alpha$ map ( $\alpha_{\rm mean}$; Col. 12). Unless we state otherwise, the spectral index is discussed in terms of $\alpha _{\rm peak}$, the index at the location of peak submm flux. Also, we do not include detection #13 (IRAS 18102-1800) in the analysis, as it appears an outlier with significantly lower $\alpha$ than any other detection; we suspect this is caused by the suspiciously weak 450 $\mu$m emission seen towards this source, most likely due to a transient telescope or calibration error as the 450 $\mu$m flux appears inconsistent with the 100 $\mu$m, 850 $\mu$m and 1.2 mm flux constraints.

The spectral index measured at the position of peak 850 $\mu$m emission varies from $\alpha_{\rm peak}=1.1$ (IRAS 19411+2306) to $\alpha_{\rm peak}=3.7$ (IRAS 20293+3952), though overall this index is fairly uniform with a sample mean of $\overline{\alpha_{\rm peak}} = 2.6 \pm 0.4$. Averaging all spectral index data around a detection slightly reduces the statistical variability, so that $\alpha_{\rm mean}$ ranges from 1.1 to 3.5, and the sample mean falls to $\overline{\alpha_{\rm mean}} = 2.4 \pm 0.3$. As implied by these statistics, the majority of detections display an $\alpha$ distribution that peaks towards the location of maximum 850 $\mu$m emission, although some sources display an anticorrelation with intensity.

Cross-sections of the $\alpha$ distributions towards IRAS 05358+3543 and IRAS 05490+2658 are presented in Fig. 7. Although these sources display very different $\alpha$ morphologies, we note that the positive and negative features seen towards the location of peak submm emission are roughly comparable in depth and width. As these sources are at roughly the same kinematic distance, the cause of these features could potentially occur on a similar spatial scale.

Table 4: Parameters used for a three Gaussian component description of the JCMT beam.

4.1.1 Peaked $\alpha$ distributions

The spectral index of the dust emission depends on a combination of the beam-averaged values of dust temperature, opacity, and spectral index of the dust opacity (${\beta }$, defined in the sense  $\tau \propto \nu^{\beta}$). As a result, there are three mechanisms which may explain the spatial distributions we observed:

1.

temperature variations through the dust envelope;

2.

emission originating from optically thick regions;

3.

changes in the composition of the dust grains themselves.

The only currently available estimate of dust temperature was derived by SBSMW using greybody fits to the SED at long wavelengths. The majority of dust cores have a temperature similar to the mean of the sample ( $\overline{T_{\rm dust}}=44$ K), although a number of cores are associated with higher temperatures: IRAS 18440-0148 with $T_{\rm dust}=97$ K, IRAS 20319+3958 with $T_{\rm dust}=73$ K, and IRAS 23151+5912 with $T_{\rm dust}=68$ K are the most prominent higher temperature cores. However, with only a single dust temperature estimate for each source, the magnitude of any dust temperature gradient across the protostellar envelopes remains unknown, which leaves the contribution of any temperature-dependent mechanism to variations in $\alpha$ unclear.

Single-dish NH₃ observations have also been conducted towards our sample, tracing gas within the cooler, extended envelope (SBSMW). From these observations, SBSMW found a mean temperature of $T_{\rm NH_{3}}=19$ K, around 25 K lower than the dust temperature in an average core. However, a core containing warm dust and an extended envelope characterised by cool gas does not prove the existence of a temperature gradient, for it is very difficult to make the gas temperature close to that of warm dust, even with the high densities ( $N_{\rm H}\simeq 10^{6}$ cm^-3) seen towards typical protostellar candidates (Goldsmith et al. 1997). However, CH₃OH and CH₃CN molecular tracers, pointing to high temperature, high density regions and indicating the presence of a hot core, have been detected towards a number of our candidate HMPOs (SBSMW). These detections imply that there are regions within the clumps of much higher temperature than the beam-averaged dust temperature alone would suggest, so we expect a strong temperature gradient must be present towards at least some of our sources.

We examined the significance of a temperature gradient by forming simulated $\alpha$ maps using the 450 $\mu$m and 850 $\mu$m continuum images created by our best-fit radiative transfer models (Williams et al. 2004). These models assume a single luminous protostar embedded in a dense, dusty envelope, and form excellent fits to the observed emission while maintaining constant dust grain characteristics (i.e. opacity and variable grain composition are not a factor in the simulated $\alpha$ map). In general, our best fit models suggest the presence of dust envelopes with temperatures around 300-500 K at the inner boundary, falling to around 10-15 K at the outer boundary. With a temperature gradient as the only factor, the simulated $\alpha$ distribution peaks towards the hottest, densest, most central regions. These centrally peaked $\alpha$ morphologies are similar to those seen towards the majority of our sample (e.g. IRAS 18247-1147; IRAS 18306-0835), suggesting the observed $\alpha$ features are dominated by temperature gradients across the envelope. This result emphasises that we must know the spatial temperature distribution of the clumps if we are to refine our investigation and accurately quantify the contribution of other factors towards these sources.

4.1.2 $\alpha$-dip distributions

On the other hand, peaked $\alpha$ distributions are not the only morphology observed: the IRAS sources 05490+2658, 18290-0924, 18530+0215, 19413+2332 and 20051+3435 form notable exceptions where $\alpha$ falls towards the location of maximum 850 $\mu$m emission. In terms of temperature gradients, these $\alpha$ morphologies run counter to the $\alpha$ distribution expected for an internally heated core: if no other factors are involved, they imply that the inner core must be cooler than the surrounding envelope. But is the formation of a hot envelope and cool inner core a realistic possibility? To form a typical negative dip morphology with an $\alpha$ valley depth of $\Delta \alpha=-0.3$, while maintaining an average dust temperature of 44 K, would require an inner core temperature of around 26 K with a surrounding envelope of ${\sim}60$ K. It would be hard to explain such low inner temperatures in the presence of large, luminous protostars when "hot-cores'' associated with typical pre-UCHII protostars have temperatures of >100 K (Kurtz et al. 2000). Alternatively, external heating could warm the outer layers of the clump relative to the inner core, but the interstellar radiation field alone is not capable of heating such dense dust to such high temperatures, and while nearby luminous stars could conceivably heat the exterior to higher temperatures our radiative transfer modelling shows the submm emission profiles are well matched by low temperatures ( ${\sim}10{-}15$ K) at the external envelope boundary (Williams et al. 2004). Furthermore, low core temperatures may preclude the very formation of a massive protostar, as the Jeans mass becomes much lower within cooler cores, suggesting that a series of lower mass protostars would form instead. This does not exclude the possiblity that a massive star could form through the coalescence of discrete low-mass protostars, but considering the weight of evidence it is hard to envisage how "cool cores'' may cause the observed variations.

An alternative explanation is that these $\alpha$-dip cores are optically thick. This possibility is unlikely, as Sect. 4 shows that all cores (with the exception of Source #8, the main component towards IRAS 18089-1732) are optically thin even at 450 $\mu$m. Even so, we recognise that this statement is based on the beam-averaged values, and there may be much denser, optically thick regions present on scales smaller than our observations can probe. For example, a circumstellar disk would lead to a large density concentration in the very inner envelope, but would remain unresolved by our observations. Then again, the presence of circumstellar disks has been confirmed towards IRAS 20126+4104 (Cesaroni et al. 1997) and IRAS 05553+1631 (Shepherd & Kurtz 1999), and the $\alpha$ distribution towards these sources remains strongly peaked. Overall, we conclude that optically thick regions do not significantly affect the $\alpha$ distribution at the spatial resolution of our measurements.

Finally, variations in the properties of the dust-grains themselves could help explain the trends. The optical properties of dust grains can be quantified by ${\beta }$, the spectral index of the dust opacity. This is often a quantity of interest as it may give information on the composition and evolutionary history of dust grains within the envelope. There are many models that predict ${\beta }$ for different grain characteristics, and the majority of grain compositions result in a spectral index of $\beta\sim1.5{-}2$ (e.g. Gezari et al. 1973; Draine & Lee 1984; Krügel & Siebenmorgen 1994), although it may range from $\beta\sim1$ (Mathis & Whiffen 1989) up to $\beta=3$ (Aannestad 1975). To observe an $\alpha$-dip morphology, dust grains within the central core must be of lower ${\beta }$ than grains in the surrounding envelope. Low ${\beta }$ and ${\beta }$ distributions that fall towards regions of high density are usually attributed to grain growth in these dense, innermost regions (e.g. Mannings & Emerson 1994; Beckwith & Sargent 1991; Goldsmith et al. 1997), and our observed $\alpha$-dip morphologies are generally consistent with this grain growth interpretation, as $\alpha$ (and thus ${\beta }$) fall preferentially towards the center of the dense cores we have observed.

However, the models of Ossenkopf & Henning (1994) predict that ${\beta }$ will only change if the dust grains do not have ice mantles, which would require inner cores with dust temperatures $\geq$100 K. The detection of CH₃OH and CH₃CN towards the $\alpha$-dip detections IRAS 19413+2332 and IRAS 18530+0215 signifies the presence of a hot core of sufficient temperature to melt ice mantles, thus permitting grain growth. However, the detection of CH₃OH and CH₃CN towards a large number of candidate protostars with centrally peaked $\alpha$ morphologies raises an interesting question: these molecular tracers imply a high central temperature - certainly high enough to melt ice mantles and permit changes in ${\beta }$, so why are there no signs of grain growth? Why do the majority of these sources have positively peaked $\alpha$ morphologies? The strongly peaked $\alpha$distributions we observe suggests that the temperature gradient effect outweighs any contribution from grain evolution. On the other hand, as it takes time for ice mantles to melt, perhaps these icy grains remain towards cores only recently heated, and perhaps these cores are younger than those associated with $\alpha$-dip distributions.

Also, CH₃OH and CH₃CN have not been detected towards IRAS 05490+2658 and IRAS 18290-0924, suggesting a hot core has not formed, yet these objects are still found with $\alpha$-dip morphologies. Clearly, neither temperature gradients or variable dust grain composition taken alone can fully explain the observed $\alpha$ morphologies. It is not clear that hot cores are strongly correlated with grain growth nor with $\alpha$-dip morphologies, and accurate high-resolution measurements of the temperature of the clumps are vital if we are to determine the magnitude of grain growth towards our sample.

   
4.2 ${\beta }$ and signs of grain growth

As we have an estimate of the temperature of dust grains within the cores, we can calculate the spectral index of the dust opacity using the equation

$$\beta = (\alpha + \Delta\alpha) - 2,$$ (5)

where $\Delta\alpha$ is a Rayleigh-Jeans correction factor, necessary because the dust body temperatures are closer to the equivalent temperature $T_{\nu}$of our observations than a Rayleigh-Jeans approximation would permit (where temperature $T_{\nu}$ at frequency $\nu$ is $T_{\nu}=h\nu/k$, giving $T_{850} \approx 17$ K and $T_{450} \approx 32$ K). Assuming the equivalent temperatures $T_{\nu_{1}}$ and $T_{\nu_{2}}$ for the observed frequencies $\nu_{1}$ and $\nu_{2}$, the required correction can be expressed

$$\Delta\alpha = \frac{ \ln \left[ \frac{\exp(T_{\nu_{2}} / T_... ..._{\rm dust})-1}\right] }{\ln (T_{\nu_{2}}/T_{\nu_{1}})} - 1.$$ (6)

For example, with a median cold-component dust temperature of 43 K, the average correction is around +0.3. By taking our $\alpha$ maps and adding $(\Delta \alpha - 2)$ to each image, the $\alpha$ maps are transformed to display the spatial distribution of ${\beta }$. However, as we only have an estimate of the dust temperature in the inner core (setting $T_{\rm dust} = T_{\rm cd}$, the cold-component dust temperature from SBSMW), likewise the estimate of ${\beta }$ obtained from these maps is only valid towards the innermost regions.

The value of ${\beta }$ at the location of peak 850 $\mu$m submm emission is listed individually for each detection in Table 3. The average grain opacity index for our sample is $\overline{\beta} = 0.9 \pm 0.4$. This index is smaller than for that seen towards high-mass stars associated with UCHII regions (Hunter 1997), and the ${\beta }$distribution is substantially shifted to lower indices compared to the more evolved objects (Fig. 8). Low ${\beta }$ is often associated with young, less evolved sources, which would point to further evidence that our sample of young high-mass stars are at an earlier stage of evolution than their UCHII counterparts. An inadequate Rayleigh-Jeans correction could raise our estimate of ${\beta }$, but even if we have globally overestimated core temperatures by 20 K, $\overline{\beta}$ would only to rise to ${\sim} 1.2$, still lower than ${\beta }$ towards more evolved UCHII sources. No relationship is found between $\alpha$, ${\beta }$ and distance (Fig. 9), suggesting resolution is not an issue.

	Figure 8: A histogram of the $\beta _{\rm peak}$ distributions for our sample of HMPOs (filled histogram) and for young high-mass stars associated UCHII regions (clear bars). The latter distribution is derived from the cold-component ${\beta }$ values listed in Table 3.5 in Hunter (1997).
Open with DEXTER

	Figure 9: A plot of $\alpha _{\rm peak}$ (left-hand plot) and $\beta _{\rm peak}$ (center and right-hand plot) against distance. Distance resolved sources are plotted by filled circles, while the square and triangular symbols represent distance unresolved sources projected to the near and far kinematic distancerespectively.
Open with DEXTER

Figure 10 shows that despite the absence of any correlation between $\alpha$ and $\tau _{850}$ or $\alpha$, ${\beta }$, $\tau _{850}$ or $T_{\rm dust}$ and source distance, there is significant correlation between ${\beta }$ and $\tau _{850}$, albeit with a large scatter. Our estimate of ${\beta }$ is dependent on the SBSMW dust temperature, which was calculated assuming $\beta=2$. Had SBSMW used lower values of ${\beta }$, they would have derived higher dust temperatures, which in turn both reduces the Rayleigh-Jeans correction and decreases the implied dust optical depth. As a result, low-${\beta }$ points in Fig. 10 would move down (due to the smaller Rayleigh-Jeans correction) and to the left (due to decreased optical depth), a shift greatest for points towards the bottom-left quadrant as these objects would be associated with the largest dust temperature increase, hence strengthening the correlation between ${\beta }$ and $\tau _{850}$.

	Figure 10: A plot of the relationship of $\alpha _{\rm peak}$ (left-hand plot) and $\beta _{\rm peak}$ (right-hand plot) to $\tau _{850}$, the optical depth at 850 $\mu$m.
Open with DEXTER

Low values of ${\beta }$ are typically attributed to growth and evolution of the dust grains within dense, dusty regions, so the observed trend of low values of ${\beta }$ at low optical depths is initially surprising. We would expect the densest, most massive clumps to undergo the most significant grain growth, resulting in lowest ${\beta }$ for these cores, whereas actually we seem to observe the opposite trend. In addition, Ossenkopf & Henning (1994) found that dust grains within a protostellar core remain below the Rayleigh scattering size limit after 10⁵ years (the typical age expected for our sample; Behrend & Maeder 2001), with the dust grain opacity changing only by a factor of ${\sim} 2$ at 850 $\mu$m. However, the results in Fig. 10 suggest a change of close to an order of magnitude in opacity: for values of ${\beta }$ around 0.5, $\tau _{850}$ has a value of around $4 \times 10^{-3}$, increasing to $\tau_{850} \sim 3 \times 10^{-2}$ for dust grains with $\beta = 1.5$. Only in the extremely dense and cold regions within circumstellar disks is it believed that a significant number of grains can grow beyond the Rayleigh limit, allowing the large shift in opacity (e.g. Schmitt et al. 1997).

These inconsistencies can be understood by considering the inhomogeneous nature of our sample, and how the clumps we have observed vary in mass (Sect. 5) and most probably evolutionary status (SBSMW). Sources within the high $\tau _{850}$, high ${\beta }$ quadrant of Fig. 10 show considerable "excess'' mass compared to the mass of an equivalent luminosity main-sequence star (cf. Fig. 11), suggesting these clumps could easily be forming protogroups or protoclusters (Sect. 5.1). In any case, the majority of dust grains within these high-mass clumps will not be intimately associated with the high-mass protostar, and will most likely remain outside the T > 100 K boundary necessary to melt ice mantles and allow large variations in ${\beta }$.

	Figure 11: A comparison of detection mass and ${\beta }$ at the location of peak 850 $\mu$m emission ( $\beta _{\rm peak}$). Sources with known distance are plotted by filled circles, while sources whose distance remains unresolved are projected to the near kinematic distance (square symbols) and the far kinematic distance (triangular symbols).
Open with DEXTER

In contrast, the much lower mass of the low $\tau _{850}$, low ${\beta }$ detections means these clumps may be forming solitary high-mass stars (solitary for high-mass stars being a relative term, which we interpret as existing with only a small number of lower-mass stellar companions), with fractionally much more dust lying inside the massive protostar's sphere of influence. This distinction becomes important when we consider that the  $\tau _{850}$ and ${\beta }$ we observe are actually the optical depth and dust opacity index averaged along the line of sight. Within the largest clumps, low ${\beta }$ grains in the vicinity of the high-mass protostar will be rendered less detectable, swamped by the higher ${\beta }$ dust grains lying within the envelope of the proto-group/cluster, whereas fractionally there will be many more low-${\beta }$ grains along the line of sight towards low mass, low  $\tau _{850}$ cores, making these evolved grains appear more prominently towards less massive cores.

Secondly, grain growth will occur in parallel with evolution of the high-mass protostar. The mass of material in the core will decrease as material either collapses to form stars or has been removed from the core by the action of stellar jets and winds. These mechanisms will act to reduce the opacity towards more evolved cores by removing dust and gas. For example, the outflows of high-mass stars may disperse up to eight times the mass of material that falls onto the star (Churchwell 1997), allowing the dust opacity to change by a larger degree than that possible via grain growth alone.

   
5 Mass

The majority of photons we detect with SCUBA originate from optically thin dust emission, so flux contours also trace the mass and column density within each field of view. The combined mass of the gas and dust can be calculated from the expression

$$M = g S_{\nu}d^{2}/\kappa_{\nu}B_{\nu}(T_{\rm dust}),$$ (7)

where $S_{\nu}$ is the flux density, d is the distance to the source, $\kappa_{\nu}$ is the absorption coefficient per unit mass of dust, g the gas-to-dust ratio and $B_{\nu}(T_{\rm dust})$ represents the Planck function for a blackbody of temperature $T_{\rm dust}$, all measured at frequency $\nu$. We adopt SBSMW kinematic distances and cold-component dust temperatures, but as companion clumps do not have individual temperature measurements we have to assume all multiple detections in a field of view have the same dust temperature. We derive the absorption coefficient at 850 $\mu$m ( $\kappa_{850}$) from the evolutionary opacity models of Ossenkopf & Henning (1994), assuming an initial hydrogen number density of $n_{\rm H}=10^{6}$ cm^-3, thin ice mantles and a formation timescale of 10⁵ years to get $\kappa_{850}=1.54 \times 10^{-2}$ cm² g^-1 and $\kappa_{450}=5.23 \times 10^{-2}$ cm² g^-1. For parity with other continuum observations, we assume a gas-to-dust ratio of 100:1, and we list the derived clump masses in Table 3. However, note that the mass may be uncertain by a factor of 2.5 or more, due to large variation in the gas-to-dust ratio (e.g. Hildebrand 1983; McCutcheon et al. 1995).

Figure 12 shows the histogram of clump mass, where we see the majority of clumps have a mass of less than 500 $M_{\odot}$, regardless of whether distance-unresolved cores are projected to the near kinematic distance or the far kinematic distance. The distance-resolved detections in our sample have a mean clump mass of around 350 $M_{\odot}$, although the median mass is less at around 100 $M_{\odot}$. Assuming the near kinematic distance for distance-unresolved sources results in a mean clump mass of $\sim$330 $M_{\odot}$, with a median roughly half this value at 143 $M_{\odot}$, while projecting to the far kinematic distance results in a mean clump mass of 1120 $M_{\odot}$ and a median mass of 460 $M_{\odot}$.

The mass of the distance-resolved and near distance-projected clumps are comparable to other continuum-derived mass estimates of similar massive protostars. For example, Mueller et al. (2002) find an average clump mass of 209 $M_{\odot}$, while Molinari et al. (2002) measure $\overline{M}=235$ $M_{\odot}$. These values support claims that natal clumps bearing massive stars may contain up to 100 times the mass of the most massive adult star that emerges (Churchwell 1997). This is a very different scenario to that seen in the primarily low-mass star-forming region of Rho Ophiucus (Motte et al. 1998), where a high fraction of the initial clump mass is seen to transfer onto the resulting low-mass protostars.

Using the mass and column density relationships defined in Hildebrand (1983), we can also translate the mass of each detection to a beam-averaged gas column density. The values we derive are presented in Table 3, where we find an average H+H₂ column density of $9\times10^{23}$ cm^-2. With an average clump diameter of 0.6 pc and projected distance of 4 kpc, this translates to a mean hydrogen number density of $3\times10^{5}$ cm^-3 through the clump.

Figure 12: A histogram of the object mass, assuming the near kinematic distance (left-hand and center plot) and far kinematic distance (right-hand plot) for distance-unresolved IRAS fields. The contribution of distance resolved detections to each histogram is plotted by an outlined histogram. The mass of each detection is calculated from 850 $\mu$m continuum emission using a 100:1 gas-to-dust ratio. The center panel display an expanded view of the near kinematic distance projection for clumps with $M<300~M_{\odot}$.

Open with DEXTER

   
5.1 Cumulative mass spectrum

Figure 13 presents the cumulative mass spectrum of our 850 $\mu$m detections, distance-limited to IRAS fields less than 5 kpc distant to increase the region of complete sampling. We estimate a completeness limit of 10 $M_{\odot}$, calculated by determing the mass of a 3-$\sigma$ detection at the upper distance limit of 5 kpc, assuming the dust temperature of the detection equals that of the sample average dust temperature, with $\overline{T_{\rm dust}}=44$ K. The best fit power-law to the mass spectrum is fairly flat below ${\sim} 80$ $M_{\odot}$, with $N({>}m) \propto m^{-0.14}$. A break in the spectrum is seen around 100 $M_{\odot}$, above which the mass spectrum steepens to $N({>}m) \propto m^{-1.32}$. The mass distribution breakpoint is found comfortably above the completeness limit, suggesting this is not an observational artefact.

Figure 13: Cumulative mass distribution of all detections found at a kinematic distance of <5 kpc, incorporating 65 clumps when assuming the near kinematic distance (upper curve) and 32 clumps for the far kinematic distance (lower curve). The mass of each object is calculated from the 850 $\mu$m emission using a 100:1 gas-to-dust ratio. The error bars correspond to $\sqrt []{N}$ counting statistics. The thick dotted and thick dashed lines plot the lines of best fit for the first seven bins and the last four bins respectively.

Open with DEXTER

We compared the index of these power-law fits to other observations, finding the distribution below 100 $M_{\odot}$ to be significantly flatter than that seen in other studies. This points to either an absence or accelerated evolution of the lower-mass clumps, both of which would reduce the submm emission observed in the lower mass range. Most likely, this reflects an absence of low-mass clumps as we have only observed high-mass candidates, preferentially sampling only the high-mass tail of the initial mass distribution, and additionally many close low-mass clumps undoubtedly lie unresolved. Correcting for these factors could raise the power-law index to a more typical IMF-like power-law slope, whereas at the moment the most similar power-law index for this region is that of clumps within molecular clouds, where the index averages around -0.60 (Kramer et al. 1998; Kramer et al. 1996). This is still a factor of four greater than our study, and without knowing how far we can minimize this difference we must question their true degree of similarity.

Above the 100 $M_{\odot}$ breakpoint, the slope of our sample becomes very similar to that of the field star IMF ( $N({>}m) \propto m^{-1.3}$ for $m/M_{\odot} \ge 0.5$; Kroupa 2001). If the IMF-like distribution is valid, the apparently similar power-law index of main-sequence and massive pre-stellar clumps would suggest the star formation efficiency within these protoclusters is relatively mass invariant; it would appear that just the breakpoint shifts to lower mass as the core fragments and additional stars form.

Assuming the multiple power-law IMF of Kroupa 2001, we may estimate the number of sources within an average 350 $M_{\odot}$ clump. Even if the star-formation efficiency $\epsilon$ for a 10 $M_{\odot}$ star is as low as 5%, this still leaves sufficient mass to form an additional ${\sim} 100$lower-mass objects created with efficiency $\epsilon=30\%$, distributed via the IMF number ratios given in Table 5. With a typical clump radius of ${\sim}0.5'$, this would result in a typical stellar volume density of ${\sim}150$ stars pc^-3.

Table 5: Stellar IMF, taken from Kroupa (2001). The third and fourth columns tabulate the contribution of each object type to the population and mass total.

With a median mass of $\sim$143 $M_{\odot}$, the typical median clump could easily form at least one >8 $M_{\odot}$ star assuming an average 50% star-formation efficiency and an IMF-like mass distribution, assuming a high-mass star accounts for 17% of the total cluster mass (Table 5). However, the limited mass reservoir of the lowest-mass companionless clumps means the star-formation efficiency $\epsilon$ within these objects must be fairly high, as the total mass of gas and dust is close to that of a high-luminosity protostar. The limited residual mass reservoir would also mean that not many lower-mass protostars can co-exist within the clumps. As a result, the stellar mass spectrum within these clumps must be skewed compared to the IMF of field stars or that within higher mass clumps. Although an isolated high-mass protostar has yet to be found, we suggest these clumps form ideal candidates of high-mass protostars with a minimal number of stellar companions.

	Figure 14: A plot displaying the correlation between mass (upper plot) and peak flux (lower plot) of coincident IRAM 1.2 mm detections and JCMT 850 $\mu$m detections. The dashed line in the lower plot traces the 1:1 mass ratio.
Open with DEXTER

6 Discussion

    
6.1 Comparison to 1.2 mm observations

The 850 $\mu$m and 450 $\mu$m emission we observe with the JCMT is found to be similar to the 1.2 mm IRAM continuum observations of Beuther et al. (2002a). The detections are comparable in size and morphology, although the extended IRAM field of view reveals more sub-clumps per target field, strengthening the observation that high-mass star-forming clumps do not exist in isolation.

The mass of coincident 850 $\mu$m detections and 1.2 mm detections correlate well when the same gas and dust characteristics are assumed, as can be seen in the upper plot of Fig. 14. Points lying away from the main trend generally represent sources with companions lying in extended emission, suggesting the difference originates in the different techniques used to measure the integrated emission around extended and multiple component sources. Admittedly, the choice of where once source ends and another begins can be subjective, and we do not consider this a cause for concern.

There is good agreement in the peak flux of 850 $\mu$m and 1.2 mm detections, which implies there is no large optical depth gradient between these wavelengths, and additionally that the studies are calibrated consistently with respect to one another, measuring approximately 7 Jy/14.4'' beam at 850 $\mu$m per 1 Jy/11'' beam measured at 1.2 mm (lower plot of Fig. 14). This is another indication that the 850 $\mu$m and 1.2 mm observations are detecting the same material and thus trace the same amount of mass. By confirming the mass of the clumps, we may place more confidence in the mass-luminosity relationship derived by Sridharan et al. (2002).

   
6.2 Distance ambiguities

We remain conscious that the distance ambiguities remaining towards a number of sources could alter the conclusions we are able to draw. Attempts to resolve the uncertainties using the C³⁴S(2-1) and CS(2-1) linewidths of Beuther et al. (2002a), comparing the continuum-derived mass estimate with the virial mass at both kinematic distances, proved fruitless. However, the simulations of Walsh et al. (2001) require the majority of our distance unresolved sources to be placed at the near kinematic distance if we are to obtain M/L ratios similiar to those expected for the majority of $10^{2}{-}10^{3}~ M_{\odot}$ clusters.

The greybody analysis of SBSMW also gives an estimate of the bolometric luminosity of each IRAS source. Assuming $M_{*} \propto L_{\rm bol}^{0.25}$ (where M_* is the mass of the protostar and $L_{\rm bol}$ is the SBSMW luminosity), similar to the mass-luminosity relationship for massive stars on the main sequence, it is simple to calculate the mass of the illuminating source. In comparing M_* to the mass of the clump at both the near and far kinematic distance projections, we found the luminosity of detections #64 (IRAS 18553+0414) and #74 (IRAS 19266+1745) to be incompatible with their near-distance clump mass, suggesting that these sources are actually located at the far kinematic distance.

   
7 Conclusions

We observed a sample of candidate high-mass protostars with the JCMT. Dust continuum emission was detected towards all sources, and from analysis of the 850 $\mu$m and 450 $\mu$m maps we reach the following conclusions:

• The average clump has a mass of ${\sim}330$ $M_{\odot}$, a total hydrogen column density of ${\sim}9\times10^{23}$ cm^-2, and an angular diameter of 30'', which projected to the mean projected distance of 4 kpc equals a linear diameter of 0.6 pc. Assuming spherical symmetry, these values translate to an average hydrogen number density of $3\times10^{5}$ cm^-3 in a typical clump.

• The clumps we detect vary in mass from around 1 $M_{\odot}$ to over 1000 $M_{\odot}$. The high luminosity and low mass of the smallest clumps suggests they will only form a minimal number of stellar companions, and that we are close to probing the final mass of the most luminous protostar. These small, isolated clumps may be analogues of the low-mass NH₃ cores (Myers & Benson 1983), from which we have learned a great deal about low mass star formation. As such, these small high luminosity cores may be the most promising sites for further high-resolution observations. The extremely large mass of the largest clumps suggest they may be forming proto-groups or proto-clusters.

• Above 100 $M_{\odot}$, the mass spectrum of the submm detections displays a distribution very similar to that of the field star IMF, suggesting the clump-mass to stellar-mass transfer efficiency is relatively mass invariant within these protoclusters.

• A large fraction of the submm clumps we detect are roughly coincident with IRAS and MSX detections in the region. This suggests the short wavelength photons and submm photons have the same origin, despite the typical envelopes being optically thick at wavelengths below 90 $\mu$m. The coincident detections could be conceivably be reconciled by non-uniformity of the envelope, e.g. where outflow cavities provide a low-opacity escape route for short wavelength photons.

• We also find some submm detections are significantly offset from their associated MSX and IRAS counterparts, in agreement with Sridharan et al. (2002). This raises the question of whether these objects are related, or whether they are actually unique objects with different submm and IR properties. If the clumps we detect are potential proto-groups or proto-clusters, it is possible that young stars in the less dense, more transparent extremeties of the envelope may be more visible at shorter wavelengths. This possibility could explain both coincident and non-coincident long/short wavelength detections.

  The degree of coincidence may point to evidence of evolution, from the oldest, most evolved sources with MSX detections and no coincident submm flux, to MSX detections with some degree of submm flux, to the youngest, most embedded sources with potentially no MSX detection and high IR optical depth. Further investigation of these possibilities will require high resolution observations at submm and far infra-red wavelengths.

• We find a companion clump fraction of ${\sim}0.6$, emphasising that clumps bearing high-mass stars do not form in isolation, and that they may lie in close proximity to other potentially star-forming clumps.

• We measure the mean surface density of companions, and find the clump spatial density distribution peaks at a separation of $r \sim 0.4$ pc. The projected distribution ( $N_{\rm proj}(r)\propto r^{-1.7}$) corresponds to a volume density distribution such that $N_{\rm vol}(r) \propto r^{-0.75}$.

• The mean spectral index of the dust emission, $\alpha$, at the position of peak submm emission is 2.6 $\pm$ 0.4, and we observe both $\alpha$ morphologies both correlated and anti-correlated with submm intensity.

  Peaked $\alpha$ morphologies can be reproduced with simple internal heating of a dusty envelope, with the observed $\alpha$ gradient resulting from the corresponding temperature gradient through the envelope.

  An $\alpha$-dip morphology could occur if the inner region of the clumps are cool with respect to their surroundings or if substantial grain growth occurs in these dense central regions. Although we suggest grain growth is the most likely factor, identifying the dominant mechanism will require high-resolution observations of temperature tracers towards the regions.

• The mean spectral index of the optical depth, ${\beta }$, at the position of peak submm emission is 0.9 $\pm$ 0.4. This is lower than ${\beta }$ towards in most dusty regions, but could be consistent with grain growth occuring in the dense clumps we observe.

Acknowledgements

We would like to thank Claire Chandler for her insight and invaluable comments that helped improve this paper. The James Clerk Maxwell Telescope is operated on a joint basis between the United Kingdom Particle Physics and Astronomy Research Council (PPARC), The Netherlands Organization for the Advancement of Pure Research (ZWO), the Canadian National Research Council (NRC), and the University of Hawaii (UH). This research has made use of NASA's Astrophysics Data System Bibliographic Services and the SIMBAD database, operated at CDS, Strasbourg, France.

References

• Aannestad, P. A. 1975, ApJ, 200, 30 [CrossRef] (In the text)
• Beck, T. L., Simon, M., & Close, L. M. 2003, ApJ, 583, 358 [CrossRef] (In the text)
• Beckwith, S. V. W., & Sargent, A. I. 1991, ApJ, 381, 250 [CrossRef] (In the text)
• Behrend, R., & Maeder, A. 2001, A&A, 373, 190 [EDP Sciences] [CrossRef] (In the text)
• Beuther, H., Schilke, P., Menten, K. M., et al. 2002a, ApJ, 566, 945 [CrossRef] (In the text)
• Beuther, H., Schilke, P., Sridharan, T. K., et al. 2002b, A&A, 383, 892 [EDP Sciences] [CrossRef] (In the text)
• Brand, J., Cesaroni, R., Palla, F., & Molinari, S. 2001, A&A, 370, 230 [EDP Sciences] [CrossRef] (In the text)
• Bronfman, L., Nyman, L. A., & May, J. 1996, A&AS, 115, 81 (In the text)
• Cesaroni, R., Felli, M., Testi, L., Walmsley, C. M., & Olmi, L. 1997, A&A, 325, 725 (In the text)
• Chipperfield, A. J., & Draper, P. W. 2001, Starlink User Note 226 CCLRC, Rutherford Appleton Laboratory, Particle Physics & Astronomy Research Council (In the text)
• Churchwell, E. 1997, ApJ, 479, L59 [CrossRef] (In the text)
• Draine, B. T., & Lee, H. M. 1984, ApJ, 285, 89 [CrossRef] (In the text)
• Fuller, G. A., Williams, S. J., & Sridharan, T. K. 2004, in preparation (In the text)
• Galli, D., & Shu, F. H. 1993, ApJ, 417, 220 [CrossRef]
• Gerazi, D., Joyce, R., & Simon, M. 1973, ApJ, 179, L67 [CrossRef] (In the text)
• Goldsmith, P. F., Bergin, E. A., & Lis, D. C. 1997, ApJ, 491, 615 [CrossRef] (In the text)
• Gomez, M., Hartmann, L., Kenyon, S., & Hewett, R. 1993, AJ, 105, 1927 [CrossRef] (In the text)
• Habing, H. J., & Israel, F. P. 1979, ARA&A, 17, 345 (In the text)
• Hildebrand, R. H. 1983, QJRAS, 24, 267 (In the text)
• Hogerheijde, M. R., & Sandell, G. 2000, ApJ, 534, 880 [CrossRef] (In the text)
• Holland, W. S., Robson, E. I., Gear, W. K., et al. 1999, MNRAS, 303, 659 [CrossRef] (In the text)
• Hunter, T. 1997, Ph.D. Thesis, Caltech
• Jenness, T., & Lightfoot, J. F. 1998, in Astronomical Data Analysis Software and Systems VII, ed. R. Albrecht, R. N. Hook, H. A. Bushouse, ASP Conf. Ser., 145 (In the text)
• Kramer, C., Stutzki, J., & Winnewisser, G. 1996, A&A, 307, 915 (In the text)
• Kramer, C., Stutzki, J., Röhrig, R., & Corneliussen, U. 1998, A&A, 329, 249 (In the text)
• Krügel, E., & Siebenmorgen, R. 1994, A&A, 288, 929 (In the text)
• Kroupa, P. 2001, MNRAS, 322, 231 [CrossRef] (In the text)
• Kurtz, S., Cesaroni, R., Churchwell, E., Hofner, P., & Walmsley, C. M. 2000, in Protostars and Planets IV, ed. V. Mannings, A. P. Boss, & S. S. Russell (Tucson: Univ. Arizona Press), 299 (In the text)
• Larson, R. B. 1995, MNRAS, 272, 213 (In the text)
• Mannings, V., & Emerson, J. P. 1994, MNRAS, 267, 361 (In the text)
• Mathis, J. S., Rumpl, W., & Nordsieck, K. H. 1977, ApJ, 217, 425 [CrossRef] (In the text)
• Mathis, J. S., & Whiffen, G. 1989, ApJ, 341, 808 [CrossRef] (In the text)
• McCutcheon, W. H., Sato, T., Purton, C. R., Matthews, H. E., & Dewdney, P. E. 1995, AJ, 110, 1762 [CrossRef] (In the text)
• Molinari, S., Brand, J., Cesaroni, R., & Palla, F. 2000, A&A, 355, 617 (In the text)
• Molinari, S., Testi, L., Rodríguez, L. F., & Zhang, Q. 2002, AJ, 570, 758 [CrossRef] (In the text)
• Motte, F., André, P., & Neri, R. 1998, A&A, 336, 150 (In the text)
• Motte, F., & André, P. 2001, A&A, 365, 440 [EDP Sciences] [CrossRef]
• Mueller, K. E., Shirley, Y. L., Evans N. J. II, & Jacobson, H. R., in Hot Star Workshop III: The Earliest Phases of Massive Star Birth, ed. P. A. Crowther, ASP Conf. Ser. (In the text)
• Myers, P. C., & Benson, P. J. 1983, ApJ, 266, 309 [CrossRef] (In the text)
• Nakajima, Y., Tachihara, K., Hanawa, T., & Nakano, M. 1998, ApJ (In the text)
• Ossenkopf, V., & Henning, Th. 1994, A&A, 291, 943 (In the text)
• Palla, F., Cesaroni, R., Brand, J., et al. 1993, A&A, 280, 599 (In the text)
• Patience, J., Ghez, A. M., Reid, I. N., & Matthews, K. 2002, AJ, 123, 1570 [CrossRef] (In the text)
• Ramesh, B., & Sridharan, T. K. 1997, MNRAS, 284, 1001 (In the text)
• Salpeter, E. E. 1955, ApJ, 121, 161 [CrossRef]
• Sandell, G., Jessop, N., & Jenness, T. 2001, The SCUBA map reduction cookbook, CCLRC, Rutherford Appleton Laboratory, Particle Physics & Astronomy Research Council (In the text)
• Schmitt, W., Henning, T., & Mucha, R., A&A, 325, 569 (In the text)
• Shepherd, D. S., & Kurtz, S. E. 1999, ApJ, 523, 690 [CrossRef] (In the text)
• Sridharan, T. K., Beuther, H., Schilke, P., Menten, K. M., & Wyrowski, F. 2002, ApJ, 566, 933 (In the text)
• Williams, S. J., Fuller, G. A., & Sridharan, T. K. 2004, in preparation (In the text)
• Wood, D. O. S., & Churchwell, E. 1989, ApJ, 340, 265 [CrossRef] (In the text)
• Walsh, A. J., Bertoldi, F., Burton, M. G., & Nikola, T. 2001, MNRAS, 326, 36 [CrossRef] (In the text)
• Zhang, Q., Hunter, T. R., Brand, J., et al. 2001, ApJ, 552, L167 [CrossRef] (In the text)

  
Online Material

Table 1: Positions of the IRAS point sources satisifying the Sridharan et al. (2002) criteria, precessed to J2000 co-ordinates, alongside the kinematic distance of the IRAS source and date(s) of observation. All sources were observed during the spring and summer of 2000. All distances are taken from Sridharan et al. 2002 with the exception of IRAS 18449-0158, for which we derive the kinematic distance using $V_{\rm LSR}$ from the CS(2-1) observations of Bronfman et al. (1996). Candidates with only the far kinematic distance listed have had their distance uncertainty resolved.

	Figure 2: continued. Maps of the 850 $\mu$m and 450 $\mu$m continuum emission ( left and center respectively) alongside the spatial distribution of $\alpha$ ( right).
Open with DEXTER

	Figure 2: continued.
Open with DEXTER

	Figure 2: continued.
Open with DEXTER

	Figure 2: continued.
Open with DEXTER

	Figure 2: continued.
Open with DEXTER

	Figure 2: continued.
Open with DEXTER

	Figure 2: continued.
Open with DEXTER

	Figure 2: continued.
Open with DEXTER

	Figure 2: continued.
Open with DEXTER

	Figure 2: continued.
Open with DEXTER

	Figure 2: continued.
Open with DEXTER

	Figure 2: continued.
Open with DEXTER

	Figure 2: continued.
Open with DEXTER

	Figure 2: continued.
Open with DEXTER

	Figure 2: continued.
Open with DEXTER

	Figure 2: continued.
Open with DEXTER