A Corona Australis cloud filament seen in NIR scattered light
Key Words.:ISM: clouds – Infrared: ISM – Radiative transfer – Submillimeter: ISM
Context: The dust is an important tracer of dense interstellar clouds but its properties are expected to undergo changes affecting the scattering and emitting properties of the grains. With recent Herschel observations, the northern filament of the Corona Australis cloud has now been mapped in a number of bands from 1.2 m to 870 m. The data set provides a good starting point for the study of the cloud over several orders of magnitude in density.
Aims: We wish to examine the differences of the column density distributions derived from dust extinction, scattering, and emission, and to determine to what extent the observations are consistent with the standard dust models.
Methods: From Herschel data, we calculate the column density distribution that is compared to the corresponding data derived in the near-infrared regime from the reddening of the background stars, and from the surface brightness attributed to light scattering. We construct three-dimensional radiative transfer models to describe the emission and the scattering.
Results: The scattered light traces low column densities of 1 better than the dust emission, remaining useful to . Based on the models, the extinction and the level of dust emission are surprisingly consistent with a sub-millimetre dust emissivity typical of diffuse medium. However, the intensity of the scattered light is very low at the centre of the densest clump and this cannot be explained without a very low grain albedo. Both the scattered light and dust emission indicate an anisotropic radiation field. The modelling of the dust emission suggests that the radiation field intensity is at least three times the value of the normal interstellar radiation field.
Conclusions: The inter-comparison between the extinction, light scattering, and dust emission provides very stringent constraints on the cloud structure, the illuminating radiation field, and the grain properties.
The interstellar clouds are hierarchical structures where gravitationally bound prestellar cores are found on the smallest scales. The study of the dense clouds is largely motivated by this connection to star formation. The filamentary nature of interstellar clouds and the possible connection between filaments and star formation has been known for a long time (e.g. Barnard 1919; Fessenkov 1952; Elmegreen & Elmegreen 1979; Schneider & Elmegreen 1979; Bally et al. 1987). The recent ground-based studies and Herschel surveys have shown cloud filaments in exquisite detail (Miville-Deschênes et al. 2010; André et al. 2010; Arzoumanian et al. 2011; Hill et al. 2011; Juvela et al. 2010) and have demonstrated that pre-stellar cores and protostars are preferentially located along these structures (André et al. 2010; Men’shchikov et al. 2010; Juvela et al. 2012). Filaments are a natural outcome of interstellar turbulence (e.g. Padoan & Nordlund 2011), with further contributions from immediate triggering by supernova explosions and the radiation and stellar winds from massive stars. However, some filaments may also be formed directly by gravitational processes, as have been modelled in numerous cosmological simulations, and in other theoretical studies (e.g. Burkert & Hartmann 2004). The filaments should fragment as dictated by the local Jeans condition, a process addressed by many theoretical studies (Inutsuka & Miyama 1997; Myers 2009). In addition to the pure compressive instability, also the so-called sausage instability may sometimes play an additional role (see Fischera & Martin 2012; McLeman et al. 2012). Although the morphology of all the star-forming clouds is not predominantly filamentary (Juvela et al. 2012), a global picture of star formation is emerging, where the turbulence creates filaments, the filaments become gravitationally unstable and subsequently fragment forming the cores that may still fed by material flowing in along the filaments. Low mass stars can be found in individual filaments while high mass stars are preferentially born at the intersections of several filaments. This scenario is supported by the observations made within the Herschel Gould Belt survey (e.g. André et al. 2010; Könyves et al. 2010; Men’shchikov et al. 2010) and other Herschel programs (Molinari et al. 2010; Schneider et al. 2010; Nguyen Luong et al. 2011; Hill et al. 2011; Schneider et al. 2012; Juvela et al. 2012) as well as by numerical simulations (Padoan & Nordlund 2011; Vázquez-Semadeni et al. 2011; Klessen 2011; Bonnell et al. 2011).
To determine the initial conditions for the star formation to occur, we need to measure the physical properties of the clouds, the filaments, and the cores. The mass distributions of the clouds can be measured with a number of methods. This is fortunate because each method suffers from different sources of uncertainty.
The dense clouds have been traditionally mapped in molecular lines. Line observations are invaluable because of the information they give on the physical state, kinematics, and chemistry of the clouds (Bergin & Tafalla 2007). On the other hand, they do not always provide reliable or consistent estimates of the mass distribution. The line intensity depends on the local physical conditions, mainly density and temperature, and the chemical abundances. If the lines are not optically thin, the radiative transfer effects introduce additional uncertainty.
For the above reasons – and because of the development in detector technology and the appearance of new ground-based, balloon-borne, and space-borne facilities – the thermal dust emission has become important as a tracer of the densest clouds. Observations of the dust at submm wavelengths are of particular interest because they are sensitive to the emission of the cold dust that may have been missed in earlier far-infrared studies. With the knowledge of dust temperature and dust opacity, the column density can be estimated. However, the colour temperature obtained from observations is known to be a biased estimator of the average dust temperature (Shetty et al. 2009b; Juvela & Ysard 2012b) and additionally the value of the dust opacity is uncertain. Furthermore, there are clear indications that the dust opacity is not constant, the variations being a likely consequence of grain coagulation and aggregation processes (Cambrésy et al. 2001; del Burgo et al. 2003; Kramer et al. 2003; Lehtinen et al. 2007). Similarly, the spectral index appears to vary from region to region and this increases the uncertainty of the colour temperature and column density estimates. The observations point to a negative correlation between the colour temperature and the observed spectral index (Dupac et al. 2003; Désert et al. 2008; Anderson et al. 2010; Paradis et al. 2010; Veneziani et al. 2010; Planck Collaboration et al. 2011b; Arab et al. 2012). Such variations can affect the accuracy of the column density estimates derived from the dust emission. However, the evaluation of errors is complicated by the observational noise (also through its effect on the relation, see Shetty et al. (2009a); Juvela & Ysard (2012a) and the unknown line-of-sight temperature variations (Shetty et al. 2009b; Juvela & Ysard 2012a). One of the main advantages of submm observations of dust emission is the large range of column densities probed. With the Herschel satellite SPIRE instrument (Griffin et al. 2010) one can probe regions from to with spatial resolution of . However, the accuracy will depend on the line-of-sight temperature variations and other opacity effects (Malinen et al. 2011).
Complementary information can be obtained from near-infrared (NIR) observations, from reddening of the light of background stars, and from the measurement of the scattered light. The observed near-infrared colour excesses can be converted to column density or to provide estimates of the visual extinction using techniques such as the NICER algorithm (Lombardi & Alves 2001). The calculations are based on the assumption that the average intrinsic colour of the stars is constant and, therefore, reliable estimates are obtained only as an average over a number of stars. Therefore, the extinction must be estimated for larger pixels, as an average over 10 or more background stars. This reduces the resolution of the extinction maps and can cause bias in the presence of strong gradients. With dedicated observations it is possible to achieve a spatial resolution of some tens of arc seconds, which is comparable to the resolution of the submm surveys. The range of values probed extends from to , depending on the stellar density and the depth of the observations.
If the NIR observations are made using an ON-OFF mode the diffuse signal is preserved and it is possible to measure the light scattered by the same dust particles that are responsible for the NIR extinction. If the cloud is optically thin at the observed wavelengths and if the cloud does not include local radiation sources, the signal will be proportional to the column density (Juvela et al. 2006). The main advantage of the surface brightness measurement is the high, potentially spatial resolution. Even in normal NIR observations at mid-Galactic latitudes, the scattered light will provide better resolution than similar observations of the background stars. The method remains useful in the range –. At higher column densities the signal saturates and the surface brightness will depend on the structure of the cloud. With multiwavelength measurements (e.g., J, H, and K bands) it becomes possible to use the intensity ratios to measure and to partially correct for the saturation of the signal (Juvela et al. 2006, 2008, hereafter Paper I and Paper II, respectively). The observed surface brightness depends on the grain properties and upon the intensity of the radiation field at the location of the cloud. Therefore, combining surface brightness data with other column density tracers one can also provide constraints to these parameters.
In this paper we will examine a filament in the northern end of the Corona Australis molecular cloud. The filament contains a dense clump with a central previously estimated to be at least . The source has been studied with the help of NIR observations using background stars and the scattered light (Juvela et al. 2006, 2008, 2009). The region has now been mapped with the Herschel satellite and these data provide the opportunity to measure the structure of the central part of the filament where the column density was too high to enable reliable estimates either with the NICER or with the surface brightness (scattered light) methods. This provides an excellent opportunity to compare the performance of the three column density tracers and to use the comparison to look for indications for changes in the dust properties.
The structure of the paper is the following: In Sect. 2 we describe the observations and the main properties of the data. In Sect. 3 we present the direct results derived from each data set and make the first comparison between the column density estimates. We continue by constructing radiative transfer models for the NIR data and the sub-millimetre data separately. These calculations are described and compared in Sect. 4. Our final conclusions are given in Sect. 5.
2 The observations
2.1 Near-infrared data
The NIR observations were made with the SOFI instrument on the NTT telescope, Chile. All observations were carried out in ON-OFF mode. The northern part, an area of was observed in 2006. The integration times on the ON-fields were 1.0 hour in the J band, 1.5 hours in the H band, and 4.5 hours in the Ks band. The map was extended in June 2007 with two fields that cover the southern side of the filament. The integration times were 1.5 hours for the Ks band and 0.5 hours in the J and H bands. On the basis of the 2MASS data (Skrutskie et al. 2006), the OFF field was estimated to have an extinction below . For the details of the observations, see Juvela et al. (2008) and Juvela et al. (2009).
The observations were calibrated using the 2MASS stellar photometry and that calibration was carried over to the surface brightness. For the final surface brightness map the southern ON fields were mosaiced together with the northern part, this requiring a small additive adjustment. The observations do not provide an absolute measure of the surface brightness. Therefore, we present the results relative to the signal within a diameter aperture centred at coordinates =19043.5, =-364936 (see Fig. 1d). Based on the 2MASS data, the visual extinction in this reference area is (see Juvela et al. 2009).
In Juvela et al. (2008) the NICER method (Lombardi & Alves 2001) was used to convert the measured stellar colour excesses to estimates of extinction. This extinction map is reproduced in Fig. 1a with values relative to the above described reference area. The frames b and c show for the same area the 250 m surface brightness and the column density derived from the Herschel data.
2.2 Herschel observations
The Corona Australis cloud has been mapped as part of the Herschel Gould Belt Survey (André et al. 2010). The entire cloud is covered by parallel mode maps that also cover the northern filament. The observation were made with the PACS instrument111PACS has been developed by a consortium of institutes led by MPE (Germany) and including UVIE (Austria); KU Leuven, CSL, IMEC (Belgium); CEA, LAM (France); MPIA (Germany); INAF-IFSI/OAA/OAP/OAT, LENS, SISSA (Italy); IAC (Spain). This development has been supported by the funding agencies BMVIT (Austria), ESA-PRODEX (Belgium), CEA/CNES (France), DLR (Germany), ASI/INAF (Italy), and CICYT/MCYT (Spain). (Poglitsch et al. 2010) at wavelengths 70 m and 160 m and with the SPIRE222SPIRE has been developed by a consortium of institutes led by Cardiff University (UK) and including Univ. Lethbridge (Canada); NAOC (China); CEA, LAM (France); IFSI, Univ. Padua (Italy); IAC (Spain); Stockholm Observatory (Sweden); Imperial College London, RAL, UCL-MSSL, UKATC, Univ. Sussex (UK); and Caltech, JPL, NHSC, Univ. Colorado (USA). This development has been supported by national funding agencies: CSA (Canada); NAOC (China); CEA, CNES, CNRS (France); ASI (Italy); MCINN (Spain); SNSB (Sweden); STFC (UK); and NASA (USA) instrument (Griffin et al. 2010) at wavelengths 250 m, 350 m, and 500 m. The observation consist of two orthogonal scans with observation identifiers 1342206677 and 1342206678. We use data provided by the Herschel Gould Belt Consortium. We use the SPIRE maps produced with the naive mapping routine in HIPE software (Ott 2010). The 160 m PACS map is made using Scanamorphos v. 16 (Roussel 2011, submitted), starting from data processed to level 1 in HIPE. The 70 m data were not used, since at that wavelength the emission is not produced by the large dust grains that are responsible for the longer wavelength emission. For the Herschel data the zero point of the intensity scale was estimated by comparison with data from the Planck satellite survey (see Bernard et al. 2010). These zero points are used in the derivation of the column density estimates. However, in the plots (e.g., Fig.1b) and in comparison with the other tracers, the zero level is set using the reference area described in Sect. 2.1.
2.3 Laboca 870 m data
The Corona Australis filament was mapped in 2008 with the APEX/LABOCA instrument. This instrument operates at a wavelength of 870 m, and has the beam size . Observations consist of long scans that were made perpendicular to the filament, and at up to 30 degree angles relative to the normal of the filament. The observations were reduced using the Boa program, version 1.11 and the calibration checked with planet observations and with the help of the source S CrA that is located immediate south-east of the area covered by the NIR data. The final noise level in the central part of the map is better than 20 mJy/beam. For details of the Laboca observations, see Juvela et al. (2009).
2.4 Mid-infrared data
The area was mapped in Spitzer guaranteed time programs with both the IRAC instrument (3.6, 4.5, 5.8, and 8.0 m) and the MIPS instrument (24, 70, and 160 m). The MIPS data were already compared to the LABOCA and NIR data in Juvela et al. (2009). The images of the four IRAC bands are shown in Fig.2. The 5.8 and 8.0 m surface brightness is affected by the nearby stars and does not trace the column density. The 3.6 and 4.5 m data will be used in Sect. 3.3 check for signs of enhanced mid-infrared scattering.
We will use also data from the WISE333The Wide-field Infrared Survey Explorer is a joint project of the University of California, Los Angeles, and the Jet Propulsion Laboratory/California Institute of Technology, funded by the National Aeronautics and Space Administration. survey that have recently become public. The WISE data includes the wavelengths of 3.4, 4.6, 12, and 22 m (Wright et al. 2010).
3 The results
3.1 The column densities
3.1.1 Column density derived from the Herschel data
The Herschel measurements of submm dust emission were converted to column density using the usual procedure that assumes a constant dust temperature along the line of sight and a constant dust opacity everywhere within the mapped area. For each map pixel, the observed intensities at 160 m, 250 m, 350 m, and 500 m were fitted with a modified black body curve to determine the colour temperature . With knowledge of the dust opacity (and the spectral index ) the fit determines the column density. The column density of the molecular hydrogen can be written formally as:
where is the assumed temperature, here equal to the colour temperature , is the average particle mass per hydrogen molecule, and the dust opacity is given relative to the gas mass.
In the analysis, the surface brightness maps are convolved to a common resolution of resolution before the determination of and (H). The spectra are modelled assuming a constant value of the spectral index, and a dust opacity 0.1 cm/g (/1000 GHz) that is applicable to high density environments (Hildebrand 1983; Beckwith et al. 1990). The same opacity law has been adopted in several papers on Herschel results (André et al. 2010; Könyves et al. 2010; Arzoumanian et al. 2011; Planck Collaboration et al. 2011a; Juvela et al. 2012). However, for any particular sources, the values of these parameters are not known to high precision and both the and the values introduce an uncertainty in the column density that may be several tens of per cent.
The calculated column density map is shown in Fig. 1c and the colour temperature map is shown separately in Fig. 3. The minimum temperature at the location of the densest clump is K. This value is obtained using the total intensity (i.e., the original zero point of the Herschel maps) without subtracting the signal in the local reference area (see Sect. 2.2).
3.1.2 The maps from NIR data
The extinction derived from the colour excesses of the background stars was calculated using the NICER method (Lombardi & Alves 2001). The resulting maps of visual extinction, have a resolution of 30. A similar map with a 20 resolution has already been presented by Juvela et al. (2008). Here we settle for a more reliable but a slightly lower initial resolution of 30 because the data will later be compared to lower resolution column density maps derived from the Herschel data. The extinction map covers the whole area included in the NIR observations (see Fig.1a). The maximum extinction values are 33 but because no background stars are visible at the centre of the dense clumps within the filament, this is only a lower limit. Furthermore, because of the strong gradients the values near the filament centre may be biased downwards (see Juvela et al. 2008).
In Juvela et al. (2008) and Juvela et al. (2009) the surface brightness, assumed to consist only of scattered light, was converted to column density assuming analytical formulas for the relation between the surface brightness and the extinction. This resulted in an extinction map that was consistent with the values up to but gave higher values in the range of 10–15. It was not clear which of the two extinction estimates was more robust although simulations indicated that most of the difference could be explained by bias in . When the extinction approaches 20 the surface brightness is strongly saturated and the simple analysis method adopted in Juvela et al. (2008, 2009) leads to diverging values. The saturation of the surface brightness and the dip at the centre of the dense clump are evident even in the band (see Fig.1d) where, compared to the band the optical depth is lower by close to a factor of ten. This alone confirms that the central must be higher than 20. As in Juvela et al. (2008), the spatial resolution of this initial map is 10.
3.1.3 The comparison of the column density maps
In Fig. 4 we compare the column density estimates from the submm emission, the reddening of the background stars (the NICER method by Lombardi & Alves (2001)), and from the NIR surface brightness data (Juvela et al. 2009).
For easier comparison, the Herschel estimates of (H) were also converted to units of visual extinction, . The conversion was done using the canonical relation (H)= (Bohlin et al. 1978). That relation was determined for more diffuse lines of sight and may not be representative of the present field, but does at least provide a convenient point of reference. The map of (a scaled version of Fig. 1c) is shown in Fig. 4a. The map of was convolved to the resolution of 40. We have masked the areas near the borders of the observed area where the result of the convolution is not well defined (where more than 10% of the convolving beam falls outside observed region). The map is shown in Fig. 4b. In the case of the NIR surface brightness, values remain undefined in the central part of the filament where the simple analytical solution diverges. In the area where the estimates exceed 25, before making the convolution, the data were replaced with values read from the NICER map. After the convolution, those areas are masked from the subsequent analysis (Fig. 4c). In the plots Fig. 4a–c the map zero points were set using the common reference region (see 2.1).
The three maps are morphologically very similar at low and intermediate extinctions. In the lowest regions the absolute value of the visual extinction is of order of 1–3. At those diffuse lines-of-sight the correspondence is best between the two NIR-derived maps because one is close to the sensitivity of the sub-millimetre data (or because the sub-millimetre maps may be affected by minor baseline uncertainties).
Figure 4d confirms the good correspondence between the tracers. Below the correlation coefficients coefficient is 0.96 for the – relation with a least squares line =-0.11+1.76. For the relation – the correlation coefficient is 0.86 but the presence of two different relations is evident. In the northern half of the map (blue points in Fig. 4d) the least squares fit gives a relation =0.15+1.05 with a correlation coefficient of 0.95. The slopes are consistent with the results presented by Juvela et al. (2008) where the was below the values up to while the situation was reversed at higher . On the southern side the relation is =0.81+2.03 with an even higher correlation coefficient, =0.98. The main difference compared to the northern side is the much steeper slope of the relation. The northern and the southern fields were calibrated separately, but a factor of two difference in the calibration is very improbable. Furthermore, the surface brightness values decrease to zero both at the southern and northern edges of the maps and the surface brightness values were matched at the centre when the fields were mosaiced together. Therefore, there is no room for an over 60% multiplicative error in the relative calibration. The difference is an indication of a change in the local NIR radiation field. Assuming that the grain properties are the same, and ignoring the fact that the scattering is mostly in the forward direction, this would suggest that the radiation field intensity is at least 50% higher on the southern side of the filament. This is qualitatively consistent with the conclusions of Juvela et al. (2009) that were based on the morphology of the NIR and the Spitzer dust emission maps.
The average ratios / and / are or higher. The interpretation of is affected by the uncertainty of the intensity of the interstellar radiation field (ISRF) but that does not influence . It would thus appear that ISRF is not very different from Mathis et al. (1983) values north of the filament but are clearly elevated on the southern side. The values appear to be too low by almost a factor of two. Part of the difference could be caused the assumed scaling between N(H) and , but it is clear that especially the value of the dust opacity has a large uncertainty. The values were derived with an opacity of 0.1 cm/g (/1000 GHz) that is higher than the value found in diffuse regions. Boulanger et al. (1996) derived for high latitude clouds a value of =0.005 cm g which is 40% of the value given by the previous formula. Thus, with the diffuse medium the values would rise above the values. This suggests that the actual value of should be somewhere in between. However, it is also known that the colour temperature used in Eq. 1 overestimates the mass averaged dust temperature and may lead to the underestimation of the column densities (see Shetty et al. 2009b; Malinen et al. 2011; Ysard et al. 2012; Juvela & Ysard 2012b). If this effect can be quantified, it becomes possible to increase the lower limit the dust opacity. This is possible only with modelling and we will return to this question in Sect. 4.
3.2 The NIR and sub-millimetre profiles of the filament
We measure the profile of the filament using the three tracers available to us. The sub-millimetre emission is the only one that, at least in principle, is capable of probing the column density distribution across the whole filament. Good estimates are missing for a small central part where no background stars are visible (roughly the white area in Fig. 4a). The NIR surface brightness data are reliable only in regions of less than . Because we want to compare all the tracers, we concentrate on the narrow region that is marked with dashed lines in Fig. 1d. We calculate the one-dimensional profile along the main axis of that region, averaging the data in the perpendicular direction.
Figure 5 shows the resulting column density profiles that, like in Fig. 4, are converted to units of visual extinction. The first frames show the data at the original resolution (cf Fig. 4). There is an offset between the maps. In the northern end the decreases below 0.5 while the other maps are, based on the large scale 2MASS extinction maps, were set at . This is not likely to be caused by an uncertainty in the zero point of the Herschel surface brightness data. The missing corresponds to a 250 m surface brightness of MJy sr, this assuming a temperature of K and a dust opacity of m) cm g (see Sect. 3.1.3). Such a surface brightness offset would be comparable to the minimum signal found in the northern part of examined area. Figure 5 again shows the discrepancy in the estimates. On the northern side these follow closely the shape of the curve but they are significantly higher in the southern end. In the plot we have not included the values for the central filament, , because of their large uncertainty near the region where the surface brightness saturates.
The frame of Fig. 5 shows the same data convolved to a common resolution of 40. The value at the offset of has been subtracted, and the values have been scaled by a factor of 1.3. With this scaling, the column density estimates are very consistent on the northern slope. For and the scaled version of the match is good up to . However, the peak of is almost one minute of arc south of the maximum. The profile is strongly skewed towards the north and rough agreement is found again in the south end of the map. If we trust the NICER map, both and show deviations on the south side of the filament and qualitatively both features could be affected by the asymmetry of the radiation field. Figure 5a shows as dotted lines the minimum and maximum along the stripe for which the other curves show the average value. These indicate that there is a strong gradient perpendicular to the stripe and parallel to the filament, especially close to the high column density clump. This could cause significant bias in (see Juvela et al. 2008) but cannot explain all the differences between and , neither in the levels nor in the profiles.
3.3 Mid-infrared profiles
Although we will not attempt to model the mid-infrared data of the Corona Australis cloud (Sect. 4) it is also useful to examine the Spitzer IRAC and the WISE satellite observations of the filament. Figure 6 compares the 3.6 m and 4.5 m Spitzer IRAC profiles with the scattered light and the extinction derived from the Herschel data.
The mid-infrared profiles show some resemblance to the NIR bands but follow more closely the column density profile of the filament as derived from the sub-millimetre observations. In particular, the 4.5 m emission peak coincides with the location of the column density maximum derived from the Herschel data. Compared to the estimated column density, both the 3.6 m and 4.5 m data show significantly higher intensity levels on the south side. We compared the Spitzer and WISE data at the wavelengths of 3.4 m and 4.6 m. The WISE measurements are in complete agreement with the Spitzer data at the corresponding wavelengths except for the 3.4 m intensity that is slightly lower in WISE at the location column density peak (0.18 MJy/sr vs. 0.21 MJy/sr) so that the intensity remains almost flat between offsets 3.0-5.5. Further outside the filament, the Spitzer ratio 3.6 m/4.5 m and the WISE ratio 3.4 m/4.6 m show some increase. This could be caused by the effect the optically thick filament has on the local radiation field but the quantification of these effects would require separate modelling.
The mid-infrared data have revealed in many dense cloud cores the presence of the “coreshine” phenomenon where the emission at 3.5 m becomes brighter relative to the emission at 4.5 m (Steinacker et al. 2010; Pagani et al. 2010). The changes are interpreted as a sign of an increase of dust grain sizes that leads to enhanced light scattering at wavelengths beyond 3 m. The first detection of coreshine was made in the cloud LDN 183 using Spitzer data (Steinacker et al. 2010). Further detections have subsequently been made with WISE observations (Juvela et al. 2012) but LDN 183 remains the best example of the phenomenon. Figure 6 shows significant increase of the 3.6 m/4.5 m ratio but only south of the column density peak. This is thus probably caused by a change of the radiation field rather than by grain growth. Regarding the overall asymmetry with respect to the filament centre, the 3.5 m profile is rather similar to that of the band.
|Radiation field||Dust model||Figures||Notes|
|isotropic,||MWD||Figs. 8, 13|
|isotropic,||OH||Figs. 8, 9, 13, 16|
|anisotropic,||OH||Figs. 10, 11, 14|
Ratio 3:1 between isotropic background and additional radiation from the south-eastern direction.
Densities twice the values obtained from the modelling of sub-millimetre emission.
4 Radiative transfer models
4.1 Model of the sub-millimetre emission
We construct a three-dimensional model to explain the Herschel observations of 160 m, 250 m, 350 m, 500 m. In the present paper we examine only the main effects, including the possibly anisotropic radiation field, using two dust models. The first model corresponds to =5.5 and is described in Draine (2003)444The data files describing the dust properties are available at http://www.astro.princeton.edu/ draine/dust/. In the following, this dust model is called MWD. As a point of comparison, we use the Ossenkopf & Henning (1994) dust model (in the following, the model OH) for coagulated grains with thin ice mantles that have accreted in years at a density of cm. At wavelengths below 1 m we adopt the short wavelength extension discussed in Stamatellos & Whitworth (2003). The dust models differ with respect to the sub-millimetre spectral index that, measured between 250 m and 500 m is for the MWD dust model and 1.76 for the OH model. The optical depth ratios for extinction are 1760 and 2080 for MWD and OH, respectively. Detailed studies of the possible variations of the dust properties are deferred to a later paper. Here the two models are used only to test the sensitivity of the results to the actual dust properties. The dust emission is calculated with our radiative transfer program (Juvela & Padoan 2003; Juvela 2005). The program uses Monte Carlo simulation to determine the radiation field intensity at each position within the model cloud. This information is used to determine the distribution of dust temperatures. The line-of-sight integration of the radiative transfer equation then results in synthetic surface brightness maps that are calculated for each of the observed wavelengths.
To avoid edge effects near the boundaries of the area covered by the NIR data, we include in the model a larger area than that shown in Fig. 1. The modelled area is , centred at coordinates =19054.3, =-365522. The model will be adjusted to reproduce the observed 350 m data. This is accomplished by having as free parameters the column densities corresponding to each map pixel. The optimisation will determine the column density for each line-of-sight and thus the mass distribution in the plane of the sky. In principle, the solution should be unique for any combination of the external radiation field and dust properties. However, the line-of-sight density distribution must be fixed for this modelling. According to Fig. 5 in the plane of the sky the FWHM of the column density distribution is perpendicular to the filament. For the assumed distance of 130 pc (Marraco & Rydgren 1981) this corresponds to a filament width of pc. We assume that the filament is approximately cylindrical. Thus, in the model the density distribution along the line-of-sight is set to be Gaussian with FWHM equal to 0.09 pc. The sensitivity of the results on this assumption is discussed further in Sect. 4. The model is modelled by a Cartesian grid of 60 cells with the cell size corresponding to 10. The surface brightness maps are calculated towards one principal axis, the simulated map consisting of pixels. When the model is compared with observations, the model results are convolved to the resolution of the corresponding observation. In particular, at 350 m the observations have a resolution of while the model discretisation corresponds to an angular resolution of 10.
We start by assuming that the filament is illuminated by an isotropic ISRF with the spectrum given by Mathis et al. (1983). The solution is obtained iteratively, calculating the model prediction of the 350 m surface brightness and scaling independently the column density corresponding to each of the map pixels. The densities of cells corresponding to a single map pixel are scaled with the same number, thus preserving the original shape of the line-of-sight density profile. The first realisation is that the CrA cannot be modelled using the standard value of the ISRF. The column densities increase without limit but the models never reach the surface brightness observed at the centre of the filament. There are two solutions to this problem. Either the sub-millimetre opacity must be increased by at least factor of 2 (or more considering the associated decrease in dust temperatures). The other alternative is to increase the radiation field.
Figure 7 shows the fit residuals in the case of the OH dust model and the scaling of the radiation field by factors =2 and 3. With the surface brightness never reaches the observed values, the error at 350 m remaining above 50%. With the 350 m can be fitted but only when the peak values are above 100 over a significant area of the map. Apart from being incompatible with the data, the model fails to produce the correct shape of the emission spectrum. For example, the 160 m surface brightness is too low typically by more than 20%.
More satisfactory solutions are found by scaling the ISRF values by a factor of 4. Figure 8 shows the observed surface brightness maps and the fit residuals for both the MWD and OH dust models. Because the model is fitted only based on the 350 m intensity, and the product of the radiation field and column density is capable of producing high enough surface brightness, the rms error at this wavelength is now below 1%.
Concentrating on the central area, the fit is better with the OH model. Apart from the noise affecting the low column density part, the relative errors are 10% at both 250 m and 500 m. There is no strong overall bias even at 160 m. This indicates that the average shape of the SED is well reproduced. There are, however, some systematic effects. The observed level of the 500 m emission is well reproduced in the filament but is overestimated is the lower density regions. The error is small, 10%, and could be associated with the uncertainty of the intensity zero point of the observations. A stronger effect is seen at the location of the eastern clump, near the boundary of the area shown. There the observed 500 m surface brightness is up to 15% lower than the model prediction. The percentage error is not large but it is very clear in the map and significant considering that the nearby wavelength of 350 m is perfectly fitted. In the same area the model fails to produce sufficiently high surface brightness at short wavelengths, the error increasing to 30% at 160 m. One possible interpretation is that the actual radiation field is higher in this area. By underestimating the heating, the model would require a column density that is too high to explain the 350 m map (hence the high 500 m signal) but does not include enough warm dust to reproduce the shorter wavelengths. The 160 m (and 250 m) residuals show a general gradient towards south-east that can indicate an increase of the radiation field intensity. This would be consistent with the higher level of NIR scattered light seen in that direction (see Sect. 3.2). An alternative explanation for the behaviour of the eastern clump involves the dust opacity spectral index. A higher value might be needed to correct the ratio between the 350 m and 500 m data that are relatively insensitive to variations in the dust temperature.
The error maps of the MWD model look qualitatively similar to the OH results but the value (considering maps other than the 350 m) is worse by more than a factor of two. Outside the filament the 160 m and the 250 m signals are mainly overestimated and thus the emission appears too warm. However, the 250–500 m spectral index of MWD is 2.1, i.e., higher than the value =1.76 of the OH model.Therefore, compared to the OH dust model, the best fit would also be expected to correspond to a somewhat lower value of the radiation field. In the eastern clump the 500 m emission is again overestimated and the m emission underestimated but, possibly because of the different values, the errors are smaller. There is still a gradient towards south-east consistent with an increase in of the dust temperature. However, as already noted, the sign of the 160 m errors has changed compared to the OH model, and the MWD model would produce insufficient intensity at the southern side of the filament.
We show in Fig. 9 the map and the filament profile for the OH model where the Mathis et al. (1983) radiation field was scaled by a factor of four. Because of the isotropy of the radiation field, the shape of the recovered profile is similar to that of . However, the predicted values are not just higher than the previous estimates, which is possible for a number of reasons, but also higher than the measurements of . This is an indication that either the assumption of the NIR extinction curve (used in the conversion of NIR colour excesses to ) is incorrect or, more likely, the ratio used in the modelling is incorrect. An increase in the sub-millimetre opacity would decrease the model . If this modification were restricted to the centre of the filament, the radiation field would not have to be changed and the correct shape of the SED would still be recovered at low and intermediate column densities. In the same figure we also show again the minimum and maximum profiles along the selected stripe (dotted lines). The shape and magnitude of in again seen to be in agreement with the maximum profile rather than the average profile.
The previous study of the NIR observations and Spitzer FIR data suggested an anisotropy with a stronger radiation field in the south (Juvela et al. 2009). Figure 5 lead to a similar conclusion. Because the Galactic plane is located north of Corona Australis cloud, there must be a more local cause for the asymmetry. The star S CrA is located on the south side of the filament, immediately east of the area covered by our NIR data. The Herbig Ae/Be star R Corona Australis is farther in the east (Gray et al. 2003). The stars appear to affect the dust emission at least in the eastern part of our maps (as seen in Fig. 7). The precise spectral type of S CrA is not known (the main component S CrA A was classified as G5Ve by Carmona et al. (2007)). Based on the observed K band magnitude of S CrA, , the NIR intensity produced by S CrA is comparable to the general ISRF up to distances corresponding to several arc minutes. R CrA may have similar contribution to the NIR radiation field although, because of large intervening extinction, its effect on dust heating is probably smaller. The other bright stars in the region, HD 176269 and HD 176270, are near the south end of the area covered by our NIR observations (see Fig. 2). The distance estimates of these B9V stars are uncertain but they may both be associated with the cloud (see Peterson et al. 2011).
We examine the situation by illuminating the model cloud with an isotropic ISRF component and adding a second radiation source in the southern direction, at a position angle of 20 degrees from south to east. The source is assumed to be at a large distance. Both radiation field components are assumed to have the same spectral shape as the normal ISRF and their combined energy input is kept the same as in Fig. 9, four times that of the normal ISRF.
When the southern source stands for 25% of the total radiation field energy, the north-south gradients in the residuals of the 160 m and 500 m fits disappear (Fig. 10). Thus the north-south radiation field anisotropy is now correct as far as can be concluded from the FIR and sub-mm emission. The result assumes that the additional radiation is coming directly from the SE direction. If the direction for the incoming radiation is not perpendicular to the line-of-sight (e.g., more on the front side of the cloud), the actual anisotropy of the radiation field may be larger. Because the dust is now warmer on the southern side of the filament the maximum of the column density should move towards north thus increasing the discrepancy with that peaks south of . However, as shown in Fig. 11, the effect on the shape of the model column density profile is not significant. The gradient in the east-west direction still remains unexplained by the model. More specifically, the clump near the east boundary of the NIR map shows, compared to our model, an excess at 250 m and low surface brightness at 500 m. As already noted, this can be a sign of a further temperature gradient (natural considering the vicinity of the R CrA region) possibly combined with an increase in the dust spectral index.
If the ISRF level is further increased by 20%, the column densities decrease but the values are higher by more than a factor of two, mainly because the 160 m intensities are overestimated throughout the field. This indicates that, apart from the inherent uncertainties involved in the modelling, the total intensity of the radiation field is well constrained. The models cannot produce sufficient surface brightness with lower level of heating and with increased ISRF the spectrum becomes incorrect. However, if the dust opacity in the cloud centre increases by a factor of a few, the lower limit of the ISRF could be relaxed correspondingly.
4.2 Modelling of the NIR scattered light
We calculate the NIR surface brightness maps for the previous models where the radiation field corresponded to four times the normal ISRF (Mathis et al. 1983) and the field is either isotropic or anisotropic. In the models of the sub-millimetre dust emission, the density distributions are similar for the MWD and the OH dust, the column densities being only % higher for the MWD. Therefore, in the following we consider only the density distributions derived from the fits of the Herschel data with the OH dust model. With this fixed density distribution, the scattered light can be calculated using both the MWD and the OH dusts. The OH model does not specify the scattering function (i.e., the directional distribution of scattered photons) and we use the same scattering function as in the case of MWD.
The differences between the optical and near-infrared properties of the two dust models are important as can be seen in Fig. 12. For MWD the scattering is still dominant in the band while for the OH model the efficiencies and are almost equal. In band the albedos are 0.68 and 0.58 for MWD and OH, respectively. In case of multiple scattering, the effect on the surface brightness can be expected to be noticeable.
The surface brightness images are shown in Fig. 13. For the MWD dust, the predicted surface brightness levels are of the right order of magnitude. The band levels are very close to the right value, the band values are 30% higher than the observed values and the band values are lower by a similar amount. This could point to a NIR spectrum of the ISRF that is redder than in the Mathis et al. (1983) model. Juvela et al. (2008) obtained similar results noting that in particular the and band values were twice as high as expected in the case of the normal ISRF. In the present study, the modelling also takes into account the effect that the optically thick filament has on the radiation field. This shadowing is probably the main reason why the modelling prefers higher values for the intensity of the external radiation field.
The NIR calculations fail to reproduce the surface brightness dip that is evident in the observations. If the modelling of the dust emission underestimates the column density at the centre of the filament, this could also partly explain the difference in the NIR colours. Because the NIR scattering takes place is mostly in the forward direction, an increase in the column density would make the spectrum of the scattered light redder. Higher column densities would appear to be incompatible with the data but, because no stars are seen through the centre of the filament (see Fig. 7 in Juvela et al. 2008), this cannot be excluded on the basis of alone. The present model is based on the modelling of dust emission and already has a peak extinction that is twice the maximum of the map. If the radiation field intensity is decreased by a factor of two, the column density could increase almost without limit and there would be no constraints on the maximum column density. This is excluded only by the fact that a lower intensity of the radiation field would alter the SED of the dust emission in a way that is incompatible with the observations.
The OH dust model (Fig. 13, bottom row) reproduces the morphology of the NIR data better than the MWD dust model. The dip in the band is now noticeable at the location of the dense clump. The depth of the depression is about half of the peak values around the central clump. This is still far from the almost complete absence of scattered light in the observations but is clearly a step in the right direction. A small depression is seen even in the band where the MWD dust produced a clear peak. The maximum band intensity is again slightly higher than the observed values. On the other hand, the band signal is very low, only one quarter of the observed.
We show in Fig. 14 similar calculations for the model of Fig. 10. The scattered light is calculated consistently assuming that, in addition to an isotropic radiation field component, 25% of the total radiation comes from the southern direction. The resulting surface brightness asymmetry between the south and north side of the filament is already slightly too large. Although the ratio of 3:1 between the isotropic and the anisotropic radiation field components was appropriate for the distribution of dust emission, the NIR surface brightness at least precludes larger degree of anisotropy, assuming that the spectrum of the anisotropic source is the same as for the isotropic component.
The major shortcoming of the models is the absence of a sufficiently strong surface brightness dip at the location of the main clump. This could be caused either by the column densities being underestimated (unlikely given the extinction data) or by some modification of the dust properties. To check the first alternative and to check that the previous differences between the OH and MWD models were not just the effect of different opacities (instead of the albedo), we recalculated the MWD model after scaling its densities by a factor of two. This results only a very minor improvement (see Fig. 15) although the maximum of the map is already well over . The effect on the morphology of the NIR maps is smaller than the difference between the two dust models.
The models predict higher 870 m intensities than those observed with LABOCA (see Fig. 16). However, a large fraction of the difference can be explained by the spatial filtering that removes large scale structure from the LABOCA map. In the NIR observations the band surface brightness goes almost to zero at the centre of the main clump. This is difficult to explain because some scattered light is always observed from the outer cloud layers. Thus, the explanation may require a very special geometry (e.g., a smaller line-of-sight extent of the cloud at the location of the clump) and/or further modifications to the dust models as suggested by the difference between the MWD and OH model results. Although the NIR signal saturates by , the data still set very strong constraints both on the cloud structure and the dust properties.
4.3 Sensitivity to the cloud shape
In the modelling the column density structure is determined directly by the observations. For fixed dust properties and a fixed external radiation field the result should be unique. However, there are still factors connected with the density distribution that can affect the results. The actual three-dimensional shape of the cloud is unknown. Although the morphology suggests the presence of a single cylindrical filament, the cloud could still be flattened or elongated along the line-of-sight. Some effects can also result from small scale inhomogeneities that allow short wavelength radiation to penetrate deeper into the cloud. We tested the sensitivity of the results to these factors using the model presented in Fig. 10.
The clumpiness was implemented by multiplying the density of each cell in the model by a random number uniformly distributed between zero and one. This results in a fair amount of inhomogeneity considering that the cloud radius is only 30 cells. All densities were then rescaled so that the model again reproduced the observed 350 m surface brightness. Compared to the original model with a smooth density distribution, this resulted in less than 10% changes in the column density along individual lines-of-sight. More importantly, this level of inhomogeneity did not introduce significant systematic change in the column density nor in the intensity of NIR scattered light produced. Of course, the effects would be stronger if the cloud included large continuous cavities extending to the centre of the cloud.
The overall cloud shape was varied by changing the FWHM of the line-of-sight density distribution by 30%. The column densities were again adjusted according to the 350 m data. When the cloud was 30% more extended along the line-of-sight, the resulting column density towards the model centre was lower by 20%. When the cloud was flattened by 30%, the percentual increase of the model column density was closer to 40%. These numbers do not yet take into account the fact that a for the best fit to the observed SED also the intensity of the external radiation field might have to be readjusted. This could reduce the net change of the column density. However, with the radiation field of Fig. 10 the 160 m surface brightness was already correct to within 10%.
We conclude that the uncertainty of the three-dimensional shape of the cloud translates to an uncertainty of some tens of per cent in the column density. These uncertainties affect the accuracy to which the dust opacity can be determined. However, the modelling is better constrained when independent extinction measurements are available.
We have examined the structure of the northern filament of the Corona Australis cloud using the combination of Herschel sub-millimetre data and near-infrared observations. The study has lead to the following conclusions:
At the three column density estimators, the thermal dust emission, the NIR reddening of the background stars, and the surface brightness caused by the NIR scattered light, are all linearly related to each other. The differences in the actual numbers can be explained by the uncertainty of the dust opacity and of the radiation field intensity.
At low column densities ( or below) the scattered light provides the best column density map in terms of resolution. Its predictions are very tightly correlated with the estimates obtained from the dust emission. However, the slope of the relation depends on the intensity of the radiation field.
Based on the dust emission, the north-south column density profile was observed to be skewed with a sharp drop in the column density on the northern side. The profile obtained from the reddening of the background stars is more symmetric but is strongly affected by the lack of stars visible through the central filament.
The modelling of the dust emission suggests that the radiation field intensity at the location of the cloud filament is about four times the value of the normal ISRF. The value can be lowered only by assuming that the dust sub-millimetre opacity has increased at the centre of the filament. However, the models of the NIR scattered light are also consistent with an elevated radiation field intensity.
According to the models, the radiation field is anisotropic with an approximate ratio of 3:1 between the isotropic component and additional radiation coming from the southern direction.
The models are unable to reproduce the deep dip in the NIR surface brightness at the centre of the filament. One possible explanation is a change in the dust properties (as indicated by the differences between the dust models examined) that has lowered the NIR albedo.
Acknowledgements.MJ acknowledges the support of the Academy of Finland Grants No. 127015 and 250741.
- Anderson et al. (2010) Anderson, L. D., Zavagno, A., Rodón, J. A., et al. 2010, A&A, 518, L99+
- André et al. (2010) André, P., Men’shchikov, A., Bontemps, S., et al. 2010, A&A, 518, L102
- Arab et al. (2012) Arab, H., Abergel, A., Habart, E., et al. 2012, A&A, 541, A19
- Arzoumanian et al. (2011) Arzoumanian, D., André, P., Didelon, P., et al. 2011, A&A, 529, L6+
- Bally et al. (1987) Bally, J., Lanber, W. D., Stark, A. A., & Wilson, R. W. 1987, ApJ, 312, L45
- Barnard (1919) Barnard, E. E. 1919, ApJ, 49, 1
- Beckwith et al. (1990) Beckwith, S. V. W., Sargent, A. I., Chini, R. S., & Guesten, R. 1990, AJ, 99, 924
- Bergin & Tafalla (2007) Bergin, E. A. & Tafalla, M. 2007, ARA&A, 45, 339
- Bernard et al. (2010) Bernard, J.-P., Paradis, D., Marshall, D. J., et al. 2010, A&A, 518, L88
- Bohlin et al. (1978) Bohlin, R. C., Savage, B. D., & Drake, J. F. 1978, ApJ, 224, 132
- Bonnell et al. (2011) Bonnell, I. A., Smith, R. J., Clark, P. C., & Bate, M. R. 2011, MNRAS, 410, 2339
- Boulanger et al. (1996) Boulanger, F., Abergel, A., Bernard, J., et al. 1996, A&A, 312, 256
- Burkert & Hartmann (2004) Burkert, A. & Hartmann, L. 2004, ApJ, 616, 288
- Cambrésy et al. (2001) Cambrésy, L., Boulanger, F., Lagache, G., & Stepnik, B. 2001, A&A, 375, 999
- Carmona et al. (2007) Carmona, A., van den Ancker, M. E., & Henning, T. 2007, A&A, 464, 687
- del Burgo et al. (2003) del Burgo, C., Laureijs, R. J., Ábrahám, P., & Kiss, C. 2003, MNRAS, 346, 403
- Désert et al. (2008) Désert, F., Macías-Pérez, J. F., Mayet, F., et al. 2008, A&A, 481, 411
- Draine (2003) Draine, B. T. 2003, ApJ, 598, 1017
- Dupac et al. (2003) Dupac, X., Bernard, J., Boudet, N., et al. 2003, A&A, 404, L11
- Elmegreen & Elmegreen (1979) Elmegreen, D. M. & Elmegreen, B. G. 1979, AJ, 84, 615
- Fessenkov (1952) Fessenkov, V. G. 1952, Trans. IAU, 8, 707
- Fischera & Martin (2012) Fischera, J. & Martin, P. G. 2012, ArXiv e-prints
- Gray et al. (2003) Gray, R. O., Corbally, C. J., Garrison, R. F., McFadden, M. T., & Robinson, P. E. 2003, AJ, 126, 2048
- Griffin et al. (2010) Griffin, M. J., Abergel, A., Abreu, A., et al. 2010, A&A, 518, L3
- Hildebrand (1983) Hildebrand, R. H. 1983, QJRAS, 24, 267
- Hill et al. (2011) Hill, T., Motte, F., Didelon, P., et al. 2011, A&A, 533, A94
- Inutsuka & Miyama (1997) Inutsuka, S.-I. & Miyama, S. M. 1997, ApJ, 480, 681
- Juvela (2005) Juvela, M. 2005, A&A, 440, 531
- Juvela & Padoan (2003) Juvela, M. & Padoan, P. 2003, A&A, 397, 201
- Juvela et al. (2006) Juvela, M., Pelkonen, V.-M., Padoan, P., & Mattila, K. 2006, A&A, 457, 877
- Juvela et al. (2008) Juvela, M., Pelkonen, V.-M., Padoan, P., & Mattila, K. 2008, A&A, 480, 445
- Juvela et al. (2009) Juvela, M., Pelkonen, V.-M., & Porceddu, S. 2009, A&A, 505, 663
- Juvela et al. (2010) Juvela, M., Ristorcelli, I., Montier, L. A., et al. 2010, A&A, 518, L93+
- Juvela et al. (2012) Juvela, M., Ristorcelli, I., Pagani, L., et al. 2012, A&A, 541
- Juvela & Ysard (2012a) Juvela, M. & Ysard, N. 2012a, A&A, 541, A33
- Juvela & Ysard (2012b) Juvela, M. & Ysard, N. 2012b, A&A, 539, A71
- Klessen (2011) Klessen, R. S. 2011, in EAS Publications Series, Vol. 51, EAS Publications Series, ed. C. Charbonnel & T. Montmerle, 133–167
- Könyves et al. (2010) Könyves, V., André, P., Men’shchikov, A., et al. 2010, A&A, 518, L106
- Kramer et al. (2003) Kramer, C., Richer, J., Mookerjea, B., Alves, J., & Lada, C. 2003, A&A, 399, 1073
- Lehtinen et al. (2007) Lehtinen, K., Juvela, M., Mattila, K., Lemke, D., & Russeil, D. 2007, A&A, 466, 969
- Lombardi & Alves (2001) Lombardi, M. & Alves, J. 2001, A&A, 377, 1023
- Malinen et al. (2011) Malinen, J., Juvela, M., Collins, D. C., Lunttila, T., & Padoan, P. 2011, A&A, 530, A101+
- Marraco & Rydgren (1981) Marraco, H. G. & Rydgren, A. E. 1981, AJ, 86, 62
- Mathis et al. (1983) Mathis, J. S., Mezger, P. G., & Panagia, N. 1983, A&A, 128, 212
- McLeman et al. (2012) McLeman, J. A., Wang, C. H.-T., & Bingham, R. 2012, ArXiv e-prints
- Men’shchikov et al. (2010) Men’shchikov, A., André, P., Didelon, P., et al. 2010, A&A, 518, L103
- Miville-Deschênes et al. (2010) Miville-Deschênes, M.-A., Martin, P. G., Abergel, A., et al. 2010, A&A, 518, L104
- Molinari et al. (2010) Molinari, S., Swinyard, B., Bally, J., et al. 2010, A&A, 518, L100
- Myers (2009) Myers, P. C. 2009, ApJ, 700, 1609
- Nguyen Luong et al. (2011) Nguyen Luong, Q., Motte, F., Hennemann, M., et al. 2011, A&A, 535, A76
- Ossenkopf & Henning (1994) Ossenkopf, V. & Henning, T. 1994, A&A, 291, 943
- Ott (2010) Ott, S. 2010, in Astronomical Society of the Pacific Conference Series, Vol. 434, Astronomical Data Analysis Software and Systems XIX, ed. Y. Mizumoto, K.-I. Morita, & M. Ohishi, 139
- Padoan & Nordlund (2011) Padoan, P. & Nordlund, Å. 2011, ApJ, 741, L22
- Pagani et al. (2010) Pagani, L., Steinacker, J., Bacmann, A., Stutz, A., & Henning, T. 2010, Science, 329, 1622
- Paradis et al. (2010) Paradis, D., Veneziani, M., Noriega-Crespo, A., et al. 2010, A&A, 520, L8
- Peterson et al. (2011) Peterson, D. E., Caratti o Garatti, A., Bourke, T. L., et al. 2011, ApJS, 194, 43
- Planck Collaboration et al. (2011a) Planck Collaboration, Ade, P. A. R., Aghanim, N., et al. 2011a, A&A, 536, A22
- Planck Collaboration et al. (2011b) Planck Collaboration, Ade, P. A. R., Aghanim, N., et al. 2011b, A&A, 536, A23
- Poglitsch et al. (2010) Poglitsch, A., Waelkens, C., Geis, N., et al. 2010, A&A, 518, L2
- Schneider et al. (2010) Schneider, N., Csengeri, T., Bontemps, S., et al. 2010, A&A, 520, A49
- Schneider et al. (2012) Schneider, N., Csengeri, T., Hennemann, M., et al. 2012, A&A, 540, L11
- Schneider & Elmegreen (1979) Schneider, S. & Elmegreen, B. G. 1979, ApJS, 41, 87
- Shetty et al. (2009a) Shetty, R., Kauffmann, J., Schnee, S., & Goodman, A. A. 2009a, ApJ, 696, 676
- Shetty et al. (2009b) Shetty, R., Kauffmann, J., Schnee, S., Goodman, A. A., & Ercolano, B. 2009b, ApJ, 696, 2234
- Skrutskie et al. (2006) Skrutskie, M., Cutri, R., Stiening, R., et al. 2006, AJ, 131, 1163
- Stamatellos & Whitworth (2003) Stamatellos, D. & Whitworth, A. P. 2003, A&A, 407, 941
- Steinacker et al. (2010) Steinacker, J., Pagani, L., Bacmann, A., & Guieu, S. 2010, A&A, 511, A9+
- Vázquez-Semadeni et al. (2011) Vázquez-Semadeni, E., Banerjee, R., Gómez, G. C., et al. 2011, MNRAS, 414, 2511
- Veneziani et al. (2010) Veneziani, M., Ade, P. A. R., Bock, J. J., et al. 2010, ApJ, 713, 959
- Wright et al. (2010) Wright, E. L., Eisenhardt, P. R. M., Mainzer, A. K., et al. 2010, AJ, 140, 1868
- Ysard et al. (2012) Ysard, N., Juvela, M., Demyk, K., et al. 2012, A&A, 542, A21