A MULTIWAVELENGTH STUDY OF YOUNG MASSIVE STAR-FORMING REGIONS. III. MID-INFRARED EMISSION

We detected MIR emission toward the observed regions in all the observed filters. Figures 1 and 2 show gray-scale images of the smoothed emission in the N-band (8.7 μm or 11.7 μm) toward 12 regions. We do not show images of the emission toward IRAS 14095−6102 and 16524−4300, which have a low signal-to-noise ratio (S/N). In seven regions, we detected a single MIR source. Most of these have a cometary-like morphology, exhibiting a bright compact component at the head and a weak trailing extended emission. In three regions (IRAS 15502−5302, 17016−4124, and 17008−4040) we detected the presence of two distinct MIR objects, and in two regions (IRAS 16128−5109 and 17271−3439) the MIR emission has a complex morphology. The morphologies of the sources detected in both the N and Q bands are very similar. The only exception is in the IRAS 14593−5852 region, where we detected a source at 17.7 μm (14593−5852 II) which is not seen in the N filter (see Figure 3).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Among the compact components, no one is strictly pointlike: whereas the average FWHM size of the standard stars is 09 in the 8.7 μm filter, 12 at 11.7 μm, and 13 at 17.7 μm, the compact components have FWHM sizes ≳14 and ≳15 at 11.7 μm and 17.7 μm, respectively. However, we do not discard the presence of point sources in a few cases, due to the uncertainties in the FWHM determination, which was made using a radial profile fitting.

According to Saviane & Doublier (2005), the pointing accuracy of the ESO 3.6 m telescope is ≃5'' which prevents a direct comparison of the MIR and radio continuum spatial distributions. We used the Spitzer–IRAC data to determine accurate astrometry of our images. The bright compact part of the smoothed TIMMI2 image was fitted to the corresponding peak in the longest wavelength IRAC image available (usually the 8.0 μm image). In most cases, the morphology of the emission exhibited in the IRAC and our images are very similar, allowing us to match the positions without any confusion. The offsets between the original TIMMI2 coordinates and the new adjusted positions are typically ≃11''. Since the GLIMPSE point-source accuracy is ≃03 (Meade et al. 2005), we estimate that our final TIMMI2 coordinates have an error ≲1'', produced by centering uncertainties (e.g., corresponding regions not compact enough or peaks not well defined).

Figures 1 and 2 also show contour maps of the radio continuum emission (Paper I, 4.8 or 8.6 GHz) overlaid on the N-band images. These figures show that in most sources the spatial distribution of the radio and MIR emissions are highly correlated, both in extended (IRAS 12383−6128, 13291−6249, 14593−5852, 15502−5302, 16128−5109, and 17009−4042) and compact structures (the same sources, together with IRAS 15520−5234, 16458−4512, and 17016−4124 I). This suggests an intimate association between the ionized gas and the warm dust within these regions (see further discussion in Section 4.2). Although in five (of 10) compact sources the corresponding radio and MIR peaks are slightly shifted, the displacement never reaches more than ≃1'' and therefore does not affect the global morphological match. Since these shifts are small, they probably can be explained uniquely by positional uncertainties of the data. The exceptional cases in which an MIR source is not associated with radio continuum emission will be discussed in Section 4.4.

Table 2 lists the source properties computed directly from the TIMMI2 images (positions, angular sizes, and flux densities). For each source, the given position corresponds to the centroid of the compact component; the brightest compact object is taken as a reference when the region is complex. The angular size was calculated as $\theta = \sqrt{\mathrm{FWHM}^2 - \mathrm{FWHM_{std}}^2}$, where FWHM is the full width at half-maximum of the compact component and FWHM_std is the averaged value of the observed standard stars, both calculated by fitting the radial profile on the N-band image (11.7 μm when available, 8.7 μm otherwise). We computed two flux densities: a "point-source" flux density F^point, calculated on compact sources using a 20 radius aperture photometry (diameter of about twice the smoothed resolution of 18), and a "total" flux density F^tot, calculated on a region (circular in almost all cases) enclosing the whole emission from each source. Photometric error was computed adding by quadrature the statistical noise of the image and the flux calibration error, which is mainly determined by the atmospheric fluctuations. We estimated this contribution as the standard deviation of the calibration factor along the three nights of observation. The resulting total photometric error is ≃10% for the 8.7 and 11.7 μm filters and 22%–34% for the 17.7 μm filter, which is affected by the presence of many atmospheric features within the filter wavelength range.

Within the 14 regions observed with TIMMI2, we can identify at least 15 simple bright compact sources with N-band deconvolved FWHM sizes in the range 07–33 (θ in Table 2). Using the kinematical distances, the physical diameters of these objects are in the range 0.008–0.18 pc, with an average value of 0.03 pc. The physical sizes of the compact MIR components associated with radio continuum sources are in good agreement with each other, when measured at a comparable angular resolution (at 8.6 GHz). This supports the idea that the radio and MIR emission trace similar structures.

Spitzer and Q-band images of all sources are shown in a dedicated Web page.⁶ In this Web page we further discuss the characteristics of the MIR emission toward each of the individual regions of the survey and summarize what is already known about them in the literature. In general, the IRAC images show similar structures to those seen in the TIMMI2 images, although the former are much more sensitive to diffuse extended emission, at 5.8 and 8.0 μm. To illustrate this, we show in Figures 3 and 4 two examples of three-color images made combining the individual 3.6, 4.5, and 8.0 μm IRAC bands, together with the corresponding N- and Q-band images from TIMMI2. These two regions are those containing the most promising high-mass protostellar object (HMPO) candidates (see Section 4.4). The diffuse emission seen in the 5.8 and 8.0 μm IRAC filters, not detected with TIMMI2, is not associated with radio continuum emission from ionized gas and probably contains a considerable contribution from polycyclic aromatic hydrocarbon (PAH) emission bands at 6.2, 7.7, and 8.6 μm (see Section 4.3). In addition, due to its higher sensitivity and shorter wavelengths available, the IRAC images reveal that in most cases the TIMMI2 source is the brightest one of a cluster of MIR objects.

Zoom In Zoom Out Reset image size

Some of the three-color images show the presence of IR features that appear greenish (the so-called green fuzzies), most of which are extended and not associated with radio continuum emission. This type of emission is found toward IRAS 15394−5358, 15520−5234, 15596−5301, 16272−4837, 16547−4247, 17008−4040 (Figure 4), and 17016−4124. Beuther et al. (2005) and Rathborne et al. (2005) first suggested that the bright and extended IRAC 4.5 μm emission might trace shocked gas by outflowing material ramming in the ambient interstellar medium. Indeed, numerical simulations predict that the contribution of H₂ shock-excited line emission in the 4.5 μm band is an order of magnitude brighter than in the other IRAC bands (Smith & Rosen 2005). Observational support for this interpretation has been provided by recent Spitzer surveys (e.g., Noriega-Crespo et al. 2004; Smith et al. 2006), which show that molecular outflows are particularly strong in the 4.5 μm band. We also found that most of the sources which exhibit green fuzzies (all but IRAS 15596−5301 and 17008−4040) are embedded in large-scale 8.0 μm extinction features known as infrared dark clouds (IRDCs), which host the earliest stages of high-mass star formation (see Rathborne et al. 2005, and references therein).

All the sources of our multiwavelength survey covered by MIPSGAL present very bright emission at 24 μm, quite saturated at the core of the PSF (with the exception of IRAS 16272−4837); in particular, the TIMMI2 objects lie just within these saturated cores. In some regions, interesting large-scale structures could be seen at 24 μm, hinting the presence of warm dust probably excited by our target sources beyond the FOV studied with TIMMI2. We tried to overcome the saturation and rescue the lost MIPS fluxes using the publicly available imageworks code⁷ written by T. Jarrett, which recovers the saturated pixels by applying a PSF model that is fit to the nonsaturated wings of the source. However, since most of our sources are embedded, the PSF wings are often contaminated by diffuse emission and the recovered flux was overestimated by the code. We therefore only consider valid those fluxes which are consistent with our spectral energy distribution (SED) models, as described in Section 4.1; they are listed in the last column of Table 3. For the remaining sources, we give lower limit fluxes, obtained by integrating the emission of the PSF wings within an aperture of 35'' and applying the aperture correction of Engelbracht et al. (2007).

Table 3 gives also the results of the photometry obtained from the IRAC images. In order to compare properly with TIMMI2 data, we integrated the emission across the same regions as for the TIMMI2 photometry, even though the emission in the IRAC images is usually more extended. This resulted, however, in a good correlation between Midcourse Space Experiment (MSX) 8.3 μm and IRAC 8.0 μm fluxes (see Section 4.1), indicating that the regions selected for photometry were appropriate. For the sources not observed with TIMMI2, we chose similar apertures (≃6''–9'' radius) enclosing the brightest source(s). According to Meade et al. (2005), magnitude uncertainties are estimated to be less than 0.2 for the GLIMPSE Catalog, corresponding to a flux error of ≲18%. Considering that our sources are not pointlike, we estimated our photometric uncertainty as 20%. For the saturated sources, we proceeded as for MIPS photometry, applying the imageworks code and considering the recovered flux only if it is consistent with the SED model; otherwise, we give a lower limit flux integrated within the TIMMI2-defined region.

Figure 5 presents the 8–13 μm spectra of the eight observed sources, extracted within small angular apertures of ≃4''–8'' around the central positions of the spectral images. The north–south oriented slit was centered in the brightest part of each source. The spectra integrated in the whole emission regions, not shown here, exhibit practically the same shape and features, but logically with a different flux scale. We truncated the spectra at 8.2 μm due to calibration uncertainties at the lower wavelength edge. Error bars are plotted as filled regions around the main curve; they only represent the statistical noise derived from the original spectral image and do not include other causes of error.

Zoom In Zoom Out Reset image size

All the spectra exhibit a broad—and in most cases prominent—silicate band seen in absorption, at ≃9.7 μm. Seven sources (all but IRAS 17008−4040 I, Figure 5(f)) present [Ne II] line emission at ≃12.8 μm, which is spatially extended. This is consistent with the association of these sources with radio continuum emission of ionized gas. Column 2 of Table 4 gives the [Ne II] line flux measured on the integrated spectrum of each source; the statistical uncertainties derived directly from the spectrum are ≲2%. Emission in other two atomic fine-structure lines, [Ar III] at 9.0 μm and [S IV] at 10.5 μm, was observed toward IRAS 16128−5109 (Figure 5(d), but this is more evident in the integrated spectrum). The [Ar III] is also present in IRAS 13291−6249 (Figure 5(b)). The PAH features at 11.3 μm and 8.6 μm were detected in only one region (IRAS 12383−6128, Figure 5(a)). The faint emission feature at ≃9.7 μm in IRAS 13291−6249, 15520−5234, and 17016−4124 I (Figures 5(a), (c), and (h), respectively) could correspond to the H₂ 0–0 S(3) line, but it is located in a wavelength range highly contaminated by telluric ozone features, which make noisier this part of the spectrum (note that this is in general not properly represented by the statistical noise shown in the spectra). We therefore label this line with a question mark in Figure 5. If present, the H₂ 0–0 S(3) line could be tracing a photodissociation region or the existence of shock events, since it can be either radiatively or collisionally excited (see, e.g., Morris et al. 2004).

We fitted the continuum part of the TIMMI2 spectrum (i.e., removing first the line features) using a simple model (Pascucci et al. 2004) which assumes that the MIR continuum flux arises from a small region of hot dust, whose emission can be described by a gray body with temperature T_d, column density N(H₂)_hot, and solid angle Ω_d, surrounded by a larger cloud of cold dust with column density N(H₂)_cold. Under these assumptions, the flux density as a function of frequency is given by

$$\begin{equation} F_\nu = \Omega _{\mathrm{d}} B_\nu (T_{\mathrm{d}})(1 - e^{-\sigma _\nu N(\mathrm{H}_2)_{\mathrm{hot}}}) e^{-\sigma _{\nu } N(\mathrm{H}_2)_{\mathrm{cold}}/2}, \end{equation} \tag{ 1 }$$

where

$$\begin{eqnarray} \sigma _{\nu } &=& R_{\mathrm{dg}}\mu m_{\mathrm{H}}\, \kappa _{\nu }. \end{eqnarray} \tag{ 2 }$$

R_dg is the dust-to-gas mass ratio, μ is the mean molecular weight per hydrogen molecule, m_H is the hydrogen mass, and κ_ν is the absorption opacity per mass of dust, whose tabulated values were taken from the Weingartner & Draine (2001) dust model, with R_V = 5.5. Chapman et al. (2009) found that this model is more consistent with the MIR extinction law of high-density regions than the Weingartner & Draine (2001) model with R_V = 3.1, because it includes larger dust grains and therefore accounts for possible grain growth. The fitting procedure uses Equation (1) with Ω_d, T_d, and N(H₂)_cold as free parameters. We constrain the parameter N(H₂)_hot by imposing the condition that

$$\begin{equation} N(\mathrm{H}_{2})_{\mathrm{hot}} = \sqrt{\frac{\Omega _{\mathrm{d}}}{\Omega _{\mathrm{c}}}}N(\mathrm{H}_{2})_{\mathrm{cold}}, \end{equation}\vspace*{-6pt} \tag{ 3 }$$

where Ω_c is the solid angle subtended by the cold component. This condition, which arises by assuming constant density and spherical symmetry for both the hot and cold components, is imposed in order to overcome the degeneracy problem of the parameters N(H₂)_hot and Ω_d, as was done for the SEDs fitting (see Section 4.1). We use as angular size of the cold cloud that determined from the 1.2 mm emission (Paper II).

Columns 3 and 4 of Table 4 list the parameters T_d and N(H₂)_cold obtained from the best fit to the spectrum extracted from the central part of each source. The best fit for each source is shown in Figure 5 as a dashed line. We assumed R_dg = 1/125 (Draine 2003) and μ = 2.8 (mean molecular weight, adopting 10% abundance of He with respect to H). The uncertainties derived from the covariance matrix computed in the least-square minimization are <1%, but we suspect that they are significantly underestimated, since the fitting only considers the statistical noise in the spectra, which is also low.

To assess the robustness of the derived values of the parameters we investigated the effects produced by changes in the assumptions. First, we used the Weingartner and Draine model for the opacities with R_V = 3.1. Since for the TIMMI2 spectral range, the R_V = 5.5 and R_V = 3.1 models do not present significant differences, we found that the fitting procedure is not sensitive to such change, and the adjusted parameters only differ by ≃5% in T_d and ≃3% in N(H₂)_cold.

The average optical depths between 8 and 13 μm inferred from the hot-component column densities N(H₂)_hot derived using the fiducial model, are in the range 0.01 ≲ 〈τ_ν〉_8-13 μm≲ 0.1, suggesting that the MIR emission from hot dust is indeed optically thin in these regions. Similar values for the optical depths at 11.7 μm were derived by De Buizer et al. (2005) from the MIR flux densities of their sample of massive star-forming regions, which are mostly ≲0.3. Hence, we also investigated the changes produced by assuming that the MIR emission is optically thin, i.e., using

$$\begin{equation} F_\nu = A \sigma _\nu B_\nu (T_{\mathrm{d}}) e^{-\sigma _{\nu } N(\mathrm{H}_2)_{\mathrm{cold}}/2}, \end{equation} \tag{ 4 }$$

where A = N(H₂)_hotΩ_d. Note that here the parameters N(H₂)_hot and Ω_d cannot be fitted simultaneously, but this allows us to avoid the constant-density condition (Equation (3)). Under this assumption we derived almost the same column densities N(H₂)_cold and temperatures than those given by the fiducial model (differences <1%).

Since in all eight sources the MIR emission is located at the peak of the millimeter core, the fitted value of the cold dust cloud column density, N(H₂)_cold, can be compared with the column density, N(H₂)_1.2 mm, derived from the peak value of the 1.2 mm flux density, S^peak_1.2 mm (Paper II). Assuming that the dust emission is optically thin at 1.2 mm, then

$$\begin{equation} N(\mathrm{H}_{2})_{\mathrm{1.2\, mm}} = \frac{S^{\mathrm{peak}}_{\mathrm{1.2\, mm}}}{\Omega _{\mathrm{B}}B_\mathrm{1.2\,mm}(T_{\mathrm{c}})\sigma _{\mathrm{1.2\, mm}}}, \end{equation} \tag{ 5 }$$

where Ω_B is the beam (24'') and T_c is the cold dust temperature derived in Section 4.1. We use κ_1.2 mm = 1 cm² g⁻¹ from Ossenkopf & Henning (1994), since the effect of coagulation and presence of dirty ice mantles on dust grains in dense cold cores is likely to be important at this wavelength, where the Weingartner & Draine (2001) R_V = 5.5 model might not represent it properly (in particular, it ignores ice mantles). The column densities calculated using this relation are given in Column 5 of Table 4. The column densities derived from both methods are in good agreement, within a factor of ≃2, particularly considering the uncertainties inherent to both methods (e.g., in the opacities and dust-to-gas ratio), as well as the differences in angular resolution between the TIMMI2 and 1.2 mm observations.

We constructed SEDs for all 18 sources of our survey, collecting the flux densities at 1.2 mm (SIMBA data, Paper II), at 12, 25, 60, and 100 μm (IRAS data), and at 8.3, 12.1, 14.7, and 21.3 μm, computed from the images of the MSX Survey of the Galactic Plane (Price et al. 2001). We also included the IRAC (nonsaturated) and TIMMI2 F^tot fluxes listed in Tables 2 and 3, summing the flux of the two components for double sources. In most cases, a simple model consisting of emission from three uniform gray bodies with different temperatures and sizes was able to fit the SED. Since the regions exhibit generally complex MIR structures and probably present temperature gradients, this is indeed a coarse simplification. However, this model allows us to determine the average dust temperatures representative of each wavelength range.

The MSX fluxes were computed just over the corresponding compact sources in the MSX images, in order to reduce the effect of the different resolution (≃18'') with respect to the IRAC data (≃18). Despite that difference, the SEDs show (Figure 6) that the MSX fluxes are consistent with the IRAC and TIMMI2 fluxes (see, e.g., the good correlation between MSX 8.3 μm fluxes and IRAC 8.0 μm fluxes), indicating that the major contribution to the "compact" MSX flux arises from the brightest source(s) within its beam. The IRAS data, however, have a much coarser resolution, so they include emission from a more extended region. Consequently, we often ignored the IRAS 12 μm flux for the fitting, and at 1.2 mm we considered the whole flux to make it representative of the same region that IRAS accounts. Conservative flux uncertainties of 20% were taken for SIMBA, IRAS, MSX, and IRAC data. For TIMMI2 data, we used the computed errors.

Zoom In Zoom Out Reset image size

The model assumes the presence of three spherically symmetric components: a cold (c) extended component, a warm (w) inner component, and a hot (h) compact component. We consider them as radiating gray bodies and include the absorption of the radiation by the enclosing components. The total flux density is then given by,

$$\begin{equation} F_{\nu } = F_{\nu }^{\mathrm{h}} + F_{\nu }^{\mathrm{w}} + F_{\nu }^{\mathrm{c}}, \end{equation} \tag{ 6 }$$

where

$$\begin{eqnarray} F_{\nu }^{\mathrm{h}} &=& \Omega _{\mathrm{h}}B_{\nu }(T_{\mathrm{h}})(1 - e^{-N_{\mathrm{h}} \sigma _{\nu }})e^{-(N_{\mathrm{c}} + N_{\mathrm{w}})\sigma _{\nu }/2}, \end{eqnarray} \tag{ 7 }$$

$$\begin{eqnarray} F_{\nu }^{\mathrm{w}} &=& \Omega _{\mathrm{w}}B_{\nu }(T_{\mathrm{w}})(1 - e^{-N_{\mathrm{w}} \sigma _{\nu }})e^{-N_{\mathrm{c}}\sigma _{\nu }/2}, \end{eqnarray} \tag{ 8 }$$

$$\begin{eqnarray} F_{\nu }^{\mathrm{c}} &=& \Omega _{\mathrm{c}}B_{\nu }(T_{\mathrm{c}})(1 - e^{-(\nu /\nu _0)^{\beta }}), \end{eqnarray} \tag{ 9 }$$

where the superscripts c, w, and h denote quantities for the cold, warm, and hot components, respectively; Ω is the solid angle; T the dust temperature, and N the H₂ column density. The cross section per hydrogen molecule, σ_ν, is given by Equation (2), where R_dg = 1/125 and μ = 2.8, and we adopted the Weingartner & Draine (2001) R_V = 5.5 model for dust opacities (see Section 3.3). For the fit we considered as free parameters ν₀, β, Ω_h, Ω_w, T_h, and T_w. The parameters of the cold component, T_c, Ω_c, and N_c were taken as fixed, and equal to the values listed in Paper II. Further, we adopted as upper limits for Ω_h and Ω_w, the values derived, respectively, from the N-band observations (Table 2) and MSX E-band observations. Because the column densities for the hot and warm components are small, the expression $(1 - e^{-N \sigma _{\nu }})$ acts as a linear term in N, producing degeneracy between the parameters N and Ω. To overcome this problem, we constrained the N_h and N_w parameters assuming constant density for all the temperature components, which implies $N_{\mathrm{h,w}} = \sqrt{\Omega _{\mathrm{h,w}}/\Omega _{\mathrm{c}}}\,N_{\mathrm{c}}$ for the hot and warm regions. Given the simplicity of the SED model and the differences in the data used, a more realistic density distribution is unnecessary.

Table 5 lists the parameters of the hot and warm components obtained from the best fit to the SED using the model described above. The FWHM angular sizes θ_h and θ_w were computed from the fitted solid angles assuming Gaussian emission distributions (Ω = 1.133 × θ²). Covariance-matrix errors are about 5%, 15%, 12%, and 30% for T_h, T_w, θ_h, and θ_w, respectively, but they could be underestimated. Most of the hot dust temperatures fall in the range ≃250–450 K, whereas warm dust temperatures fall in the range ≃60–150 K.

Figure 6 shows the SEDs and their corresponding fits for all sources. In most cases the fit is excellent. For the objects with the largest cold cloud column densities (N_c > 3 × 10²³ cm⁻²), which gives rise to a very deep silicate absorption feature at ≃10 μm (e.g., IRAS 17271−3439, Figure 6(r)), the inclusion of absorption in the model in the MIR range became difficult. Because of this, particularly unstable fits were obtained for IRAS 17009−4042 and 17016−4124 (Figures 6(o) and (p), respectively). Hence, for these objects we set an upper limit of 400 K for T_h in order to obtain reasonable output parameters. The SEDs of the sources IRAS 16272−4837 and 16547−4247 (Figures 6(j) and (m) exhibit a peculiar tendency in the IRAC fluxes, resulting in high hot-dust temperatures (591 K and 504 K, respectively). In these sources, the IRAC 4.5 μm flux density is similar (only slightly lower) to the 5.8 μm flux, which is higher than the 8.0 μm flux. This trend is consistent with the high absorption and the green fuzzies found in these regions.

For each source, the best-fit SED model was used as a validity check for the saturated IRAC/MIPS fluxes, which were recovered by applying the imageworks code (see Section 3.2), but were not used by the SED fitting. We list in Table 3, and show in Figure 6 as filled triangles, only those saturated fluxes that satisfy the condition

$$\begin{equation} \left(\frac{F_{\lambda _0}^{\mathrm{sat}} - F_{\lambda _0}^{\mathrm{model}}}{\Delta F_{\lambda _0}}\right)^2 \le \frac{1}{N} \sum _{i=1}^{N} \left(\frac{F_{\lambda _{i}} - F_{\lambda _{i}}^{\mathrm{model}}}{\Delta F_{\lambda _{i}}}\right)^2, \end{equation} \tag{ 10 }$$

where $F_{\lambda _0}^{\mathrm{sat}}$ is the recovered flux at the saturated wavelength λ₀ = 4.5, 5.8, 8.0, or 24 μm, F^model_λ is the flux density given by the SED model, N is the number of data points being fitted, $\lbrace F_{\lambda _{i}}\rbrace$ are the observed fluxes, and ΔF_λ represents the corresponding error (it was taken as 30% for the saturated flux). In other words, we request that $F_{\lambda _0}^{\mathrm{sat}}$ did not deviate from the SED model more than the χ²-per-data point value of the fit. Using this criterion, it was found that the imageworks tool had about 50% of success for embedded sources (e.g., seven out of 16 sources with saturated MIPS data).

The high correlation between the spatial distribution of the MIR and radio emissions may suggest that the MIR emission is due to the ionized gas. However, extrapolating the observed radio continuum flux of the optically thin H II regions (Paper I) to the MIR, assuming a ν^−0.1 frequency dependence, we find that the free–free emission at MIR wavelengths is negligible in comparison to the observed MIR fluxes (less than ≃6% at 8.0, 11.7, and 17.7 μm). We conclude that the MIR emission is tracing warm dust, located either within the H II region or in a thin shell around it. Kraemer et al. (2003) already found well-correlated MIR and radio distributions, and discarded the possibility that the emitting dust is distributed in a thin shell, since in that case the MIR emission would show evidence of limb-brightened morphology, which was not observed. None of the regions of our sample exhibit limb-brightened MIR emission.

We conclude that the MIR emission of the sources associated with radio-continuum emission arises from warm dust mixed with the ionized gas inside the H II region. Although dust grains and gas particles are likely to be dynamically coupled in H II regions, the heating of the dust grains is dominated by the UV radiation from the exciting source(s) and the ionized nebula, and not by the collisional exchange of energy with the gas (see, e.g., Spitzer 1978). This makes possible the coexistence of a plasma with a kinetic temperature of ≃7000 K and solid particles with a temperature of ≃250–450 K (Section 4.1). These dust temperatures are consistent with the models by Natta & Panagia (1976) for dusty H II regions.

The N-band TIMMI2 images (and probably the Q-band images too) are sensitive mainly to the emission of warm dust intimately associated with the ionized gas, or close to the exciting sources when there is no radio emission. They do not trace the presence of PAH, as indicated by the 8–13 μm spectroscopy toward seven of the eight sources of our sample (all but IRAS 12383−6128), where the emission bands in this wavelength range (the more prominent at 8.6 and 11.3 μm) are not seen.

Within the H II region or very close to a massive protostar, hard UV radiation destroys these molecules. For example, Povich et al. (2007) calculated the photodestruction rate of PAHs by extreme-UV (EUV) photons and found that this destruction mechanism accounts for the absence of PAHs in the M17 H II region. Because this region is one of the most energetic in the Galaxy, we expect that for the UC sources of our sample, the PAH destruction edge is located much closer to the exciting star(s). IRAS 12383−6128, which is the only one that presents PAH emission bands at 8.6 and 11.3 μm, has probably a softer UV spectrum than the others.

It is possible that a part of the outer diffuse emission observed in the TIMMI2 images is due to PAHs, but since the slit was centered in the bright compact MIR sources, the PAH radiation was not strong enough to be detected. This explanation is supported by other MIR-spectroscopic observations found in the literature, with larger apertures, toward some of our sources, in which PAH emission is detected. In particular, the IRAS LSR spectrum (aperture of 6' × 5') of IRAS 13291−6249 (Jourdain de Muizon et al. 1990) and the ISO spectrum (aperture of 14'' × 20'' in the wavelength range of TIMMI2) of IRAS 16128−5109 (Peeters et al. 2002) reveal the presence of PAH emission bands within the wavelength range of the TIMMI2 spectra (at 8.6 and 11.3 μm). For IRAS 16128−5109, the PAH bands at 3.3, 6.2, 7.7, and 12.7 μm are also present (Peeters et al. 2002). PAH emission was also detected toward other sources of our sample which were not observed spectroscopically by us: IRAS 14593−5852 (IRAS LSR, Zavagno et al. 1992), IRAS 15502−5302 (ISO, Peeters et al. 2002), IRAS 15596−5301, and IRAS 17271−3439 (IRAS LSR; Jourdain de Muizon et al. 1990).

The latter result strongly suggests that the extended diffuse emission seen in most regions by the Spitzer–IRAC images beyond the radio continuum contours include an important contribution of PAH emission bands. The [3.6], [5.8], and [8.0] IRAC filters cover the PAH features at 3.3, 6.2, and 7.7 μm, respectively. The [8.0] filter also covers the PAH band at 8.6 μm. This extended diffuse emission is, however, more clearly seen in the [5.8] and [8.0] IRAC images than in the [3.6] images. A possible explanation (Povich et al. 2007) is that the 6.2 μm feature might have a greater fractional contribution to the [5.8] filter than the 3.3 μm band to the [3.6] filter (and the same for the [8.0] and [3.6] filters), or simply that the 3.3 μm band is intrinsically fainter. As for the TIMMI2 images, the Spitzer MIR emission at the location of the UC H II regions is more likely to correspond to warm dust continuum radiation with little or any contribution of PAH features.

The high angular resolution MIR observations presented here revealed the presence of bright compact MIR sources in the vicinity of, but not coincident with, UC H II regions, which are likely to signpost HMPOs. We identify a total of five new HMPO candidates within the regions observed with TIMMI2: IRAS 14593−5852 II (only detected at 17.7 μm, see Figure 3), 17008−4040 I and II (Figure 4), 17016−4124 II (Figure 2(d)), and the compact component of IRAS 17158−3901 (Figure 2(e)).

In order to estimate the spectral types of these HMPO candidates, we derived MIR luminosities from the measured N- and Q-band fluxes, following a method based on the work by De Buizer et al. (2000). This technique provides good lower limits to the true bolometric luminosities of the stellar sources (see De Buizer et al. 2002). We note that estimates made using the IRAS fluxes are inappropriate since they include the contribution from the exciting star(s) of the nearby UC H II regions, as well as from other IR sources within the large IRAS beam. First, a dust color temperature is computed assuming that the MIR emission from warm dust is optically thin, as we suggested in Section 3.3. Our SED fits are also consistent with optically thin emission from the hot-dust component (using the derived N_h, the optical depths at 10 μm are in the range ≃0.01–0.3). We adopt a standard expression for the flux density (analog to Equation (4)),

$$\begin{equation} F_{\nu } = A \sigma _{\nu } B_{\nu }(T_{\mathrm{color}}) e^{-N_{\mathrm{c}}\sigma _{\nu }/2}. \end{equation} \tag{ 11 }$$

Here, the factor A = N_hΩ is constant, σ_ν is given by Equation (2), and N_c by Equation (5), but truncated to the maximum column density derived from TIMMI2 spectra if it is higher than that value (17.1 × 10²² cm⁻², see Table 4), in order to moderate the effect of dust absorption (see below). The color temperature T_color is derived from the ratio of the integrated fluxes at the N- and Q-band (lower limit N-band flux for IRAS 14593−5852 II):

$$\begin{equation} \frac{e^{hc/(k\lambda _2 T_{\mathrm{color}})} - 1}{e^{hc/(k\lambda _1 T_{\mathrm{color}})} - 1} = \frac{F_{\nu _1}}{F_{\nu _2}} \left(\frac{\lambda _1}{\lambda _2}\right)^3 \frac{\kappa _{\nu _2}}{\kappa _{\nu _1}} e^{-N_{\mathrm{c}}(\sigma _{\nu _2} - \sigma _{\nu _1})/2}. \end{equation} \tag{ 12 }$$

The MIR luminosity is then obtained integrating an expression of the form F_ν = Aσ_νB_ν(T_color) from 1 to 1000 μm  using the kinematical distance and the factor A to scale the flux density to the observed ones. Columns 2 and 3 of Table 6 give the computed color temperatures T_color and MIR luminosities L_MIR, respectively. We find that the derived spectral types, using the calibration of Panagia (1973) and the MIR luminosities, are earlier than B3 for four of our HMPO candidates (Column 4 of Table 6).

The spectral types derived by us are earlier than the typical spectral types obtained by De Buizer et al. (2000) using a similar method. A substantial difference between our technique and theirs is that we have considered the absorption of a cold cloud to model the observed fluxes, and corrected by extinction when computing the MIR luminosities. We found, however, that the employed method is quite sensitive to assumed column density N_c of the absorbing cloud, and therefore, in order to obtain strict lower limits for the spectral types, we recomputed the MIR luminosities under the assumption that there is no extinction (N_c = 0). This results now in two sources earlier than B3, namely, B0 for IRAS 14593−5852 II and B2 for 17008−4040 I, which thus might be considered as genuine HMPOs.

Figures 3 and 4 present the whole set of MIR images for the regions containing the two bona fide O-type high-mass protostar candidates: IRAS 14593−5852 II and 17008−4040 I. The Spitzer-IRAC images toward IRAS 14593−5852 (Figure 3(a)) show only diffuse emission toward the HMPO, bright in red (8.0 μm) and thus probably corresponding to the environment's PAH emission, which dominates the large-scale structure. IRAS 14593−5852 II is then likely to be a deeply embedded object. IRAS 17008−4040 I is bright in all available MIR images, even saturated in the IRAC filters (Figure 4(a)). The H II region is displaced by ≃30'' to the west of the TIMMI2 sources and associated with extended diffuse emission at 8.0 μm. In addition to its high MIR luminosity and the lack of a detectable H II region (as indicated by the lack of both radio continuum emission and of the [Ne II] line in its N-band spectrum), there is additional evidence supporting the nature of IRAS 17008−4040 I as a high-mass protostar: (i) the Spitzer–IRAC image shows the presence of extended filamentary green fuzzies (see Section 3.2) to its northeast, suggesting outflowing activity; (ii) it is located at the center of a massive and dense molecular core (peak of 1.2 mm emission, Paper II), as predicted by all the massive star formation theories; and (iii) it is associated with different types of maser emission (Figure 4(b)), all of which are known to be associated with high-mass star-forming regions: OH masers at 1665/1667 MHz (e.g., Forster & Caswell 1989; Caswell 1998) and 6.7 GHz methanol masers (e.g., Walsh et al. 1998) just coincident with the position of IRAS 17008−4040 I, and water vapor maser emission at 22 GHz located ≃4'' to its southeast (Forster & Caswell 1989). In particular, 6.7 GHz methanol masers are only detected toward high-mass star-forming regions (Minier et al. 2003).

Finally, although the MIR luminosities given in Table 6 should only be considered as rough estimates, it is interesting to know where the HMPO candidates are located in the L_bol–M_env (bolometric luminosity–envelope mass) diagram presented in the work of Molinari et al. (2008, their Figure 9). They successfully used it to follow the pre-main-sequence evolution of massive young stellar objects, by classifying a numerous sample of objects according to SED modeling and then comparing with the predictions of the turbulent core model (McKee & Tan 2003). We can do such exercise only for IRAS 17008−4040 I, since it is the only HMPO candidate that is far enough from nearby UC H II regions and thus the mass of the associated 1.2 mm core can be considered as its envelope mass (neglecting also the secondary object IRAS 17008−4040 II, which is much less luminous). Then, assuming L_bol = 40900 L_☉ (MIR luminosity) and M_env = 1200 M_☉ (Paper II), IRAS 17008−4040 I lies in a zone dominated by MM-P objects, below that occupied by IR-P objects in the notation of Molinari et al. (2008). They suggested that the IR-P stage may correspond to the arrival of the protostar on the zero-age main sequence (ZAMS), and the MM-P objects represent an earlier stage. The location of IRAS 17008−4040 I in the L_bol–M_env plot would then categorize it as a pre-ZAMS object.