The origin of mid-infrared emission in massive young stellar objects: multi-baseline VLTI observations of W33A

W. J. de Wit; M. G. Hoare; R. D. Oudmaijer; S. L. Lumsden

doi:10.1051/0004-6361/200913209

Free Access

Issue		A&A Volume 515, June 2010


Article Number		A45
Number of page(s)		13
Section		Interstellar and circumstellar matter
DOI		https://doi.org/10.1051/0004-6361/200913209
Published online		08 June 2010

A&A 515, A45 (2010)

The origin of mid-infrared emission in massive young stellar objects: multi-baseline VLTI observations of W33A

W. J. de Wit - M. G. Hoare - R. D. Oudmaijer - S. L. Lumsden

School of Physics & Astronomy, University of Leeds, Woodhouse Lane, Leeds LS2 9JT, UK

Received 30 August 2009 / Accepted 9 December 2009

Abstract
Aims. In this paper we aim to determine the structure on 100 AU scales of the massive young stellar object W33A, using interferometric observations in the mid-infrared. This emission could be caused by a variety of elements, for example, the inner protostellar envelope, outflow cavity walls, or a dusty or gaseous accretion disk.
Methods. We used the Unit Telescopes of the VLT Interferometer in conjunction with the MIDI instrument to obtain spectrally dispersed visibilities in the N-band on 4 baselines with an angular resolution between 25 and 60 milli-arcsec (equivalent to 95 and 228 AU at 3.8 kpc). The visibility spectra and spectral energy distribution were compared to 2D-axi-symmetric dust radiative transfer models with a geometry that includes a rotationally flattened envelope and outflow cavities. We assumed an O 7.5 ZAMS star as the central source, consistent with the observed bolometric luminosity. The observations were compared to models with and without (dusty and gaseous) accretion disks.
Results. The visibilities are between 5% and 15%, and the non-spherically symmetric emitting structure has a typical size of 100 AU. A satisfactory model is constructed to reproduce the visibility spectra for each (u,v) point. It fits the N-band flux spectrum, the mid-infrared slope, the far-infrared peak, and the (sub)mm regime of the SED. It produces a 350 $\mu$ m morphology consistent with the observations.
Conclusions. The mid-infrared emission of W33A on 100 AU scales is dominated by the irradiated walls of the cavity sculpted by the outflow. The protostellar envelope has an equivalent mass infall rate of $7.5 \times 10^{-4}~\mbox{$M_{\odot}$ ~yr$^{-1}$ }$ , and an outflow opening angle of $2\theta=20\hbox{$^\circ$ }$ . The visibilities rule out the presence of any dust disk with total (gas and dust) mass more than $0.01~M_{\odot}$ . Within the model, this implies a disk $\dot{M}_{\rm acc}$ of less than $1.1 \times 10^{-7} (\alpha/0.01)~\mbox{$M_{\odot}$ ~yr$^{-1}$ }$ , where $\alpha$ is the viscosity of the Shakura-Sunyaev prescription. However, optically thick accretion disks, which are inside the dust sublimation radius, are allowed to accrete at rates equalling the envelope's mass infall rate (up to $10^{-3}~\mbox{$M_{\odot}$ ~yr$^{-1}$ }$ ) without substantially affecting the visibilities due to the extinction by the extremely massive envelope of W33A.

Key words: stars: formation - stars: early-type - ISM: jets and outflows - accretion, accretion disks - techniques: interferometric

1 Introduction

Circumstellar accretion disks are an essential element in the formation of stars. The relatively long time scales involved in low-mass star formation (SF) and the proximity of the regions where they form allow the disks to be imaged at millimetre (Dutrey et al. 1996), near-infrared (near-IR; Padgett et al. 1999) or even optical wavelengths (Grady et al. 1999) and studied in great detail. The short time scales on which high-mass stars are formed and the distance of massive SF regions amongst others render the study of the accretion process particularly difficult. Whether or not a high-mass star grows in mass through the accretion of matter delivered to the stellar surface by means of a circumstellar disk is therefore still left unanswered. Recent 3-D radiation hydrodynamic simulations demonstrate that accretion through a disk can continue in the presence of strong radiation pressure (Krumholz et al. 2009). The pressure is greatly reduced by radiation escaping through the outflow cavities (Krumholz et al. 2005). Observationally a growing number of direct and indirect pieces of evidence point to large-scale (10 000 AU) rotating toroids (e.g. Beltrán et al. 2004; Beltrán et al. 2005; Furuya et al. 2008), but also to intriguingly disk-like structures on scales of several 100 AU in a small number of comparatively nearby examples (e.g. Hoare 2006; Jiménez-Serra et al. 2007, for a review see Cesaroni et al. 2007).

The reality of a disk accretion scenario in massive star formation can be assessed by observations of massive stars during the main mass-gathering phase at significant angular resolution. Prime examples of such young massive stars constitutes the class of object referred to as massive YSO (MYSO), but also as high-mass protostellar object, and BN-object. The class is characterised by luminous (> $10^{4}~\mbox{$L_{\odot}$ }$ ) infrared objects with spectral energy distributions that peak around 100 $\mu$ m. The luminosity indicates a single ZAMS star that would be able to ionise its surroundings, although too little or no radio emission is observed. This absence has been used as an argument to claim ongoing mass infall onto the central object. As a result the star is extended hence relatively cool (Hosokawa & Omukai 2009) and therefore unable to ionise its surroundings (Hoare & Franco 2007). Alternatively, high-mass infall rates (Walmsley 1995) or gravitational entrapment (Keto 2002) could quench the development of an H II that dominates the radio continuum emission. These alternative scenarios invoke, however, very high emission-measure (EM) H II regions to make them optically thick in the radio, e.g. $\rm\tau_{1~cm}= 3 \times 10^{-10} \times EM~(cm^{-6}~pc)$ (Osterbrock 1989). However, it is unlikely in this picture that the region would still be optically thick in the near-IR H I recombination lines Br $\alpha$ and Br $\gamma$ , as these have line-centre optical depths that are 10^5-6 times less than $\rm\tau_{1~cm}$ (Hummer & Storey 1987). This would imply that these objects should have strong near-IR H I lines, and yet they only show weak, broad lines indicating a stellar wind origin (Bunn et al. 1995). Nonetheless, the MYSO class has the potential to test the disk-accretion scenario. The actual accretion region is probably found on AU to tens of AU scales near the central stellar object. This region is, however, buried by hundreds of $A_{\rm v}$ , the dust that makes up the extended ( $\sim$ 0.1 pc) protostellar envelope (van der Tak et al. 2000). Only long-wavelength radiation is capable of carrying away information on the physical processes in play, limiting the attainable spatial resolution with single-dish apertures. Interferometers can overcome the problem of limited resolution and the near and mid-infrared wavelength regions currently deliver the highest angular resolution by means of e.g. the Very Large Telescope Interferometer (VLTI).

The inner parts of the dusty MYSO protostellar envelope is heated to a few 100 K and emits strongly at mid-IR wavelengths. The 10 $\mu$ m emission is usually unresolved at the 0.3 $^{\prime \prime }$ resolution of single 8 m dish telescopes, but it is partially resolved at 24.5 $\mu$ m (de Wit et al. 2009). Here, we present new observations made of the MYSO W33A with the MID-infrared Interferometric (MIDI) instrument at the VLTI. MIDI is a two-beam combiner that operates in the N-band (8-13 $\mu$ m) and that delivers spectrally dispersed interferometric observables (see Leinert et al. 2003). W33A is one of the most intensely studied MYSO thanks to its IR spectrum rich in ice features (Gibb et al. 2000). Its physical nature has received much less attention. The object has an IRAS luminosity of $10^{5}~\mbox{$L_{\odot}$ }$ (Faúndez et al. 2004) for a kinematical distance of 3.8 kpc (Bronfman et al. 1996). The object shows weak, compact, and optically thick radio continuum emission (Rengarajan & Ho 1996; van der Tak & Menten 2005) and broad ( $\sim$ 100 km s^-1), single-peaked H I recombination emission consistent with an ionised stellar wind origin (Bunn et al. 1995). Modelling of the protostellar envelope with spherical dust radiative transfer (RT) model at constant density by Gürtler et al. (1991) required the presence of an inner dust-free cavity of 135 AU (35 milli-arcsec at 3.8 kpc), a size that would be well resolved by MIDI. In de Wit et al. (2007, hereafter Paper I) we presented a single MIDI baseline for W33A in order to probe this hundred AU scale region. Via 1-D modelling we found an inner dust-free radius at 7 milli-arcsec (26.5 AU), corresponding to the dust sublimation radius. We concluded that the dispersed MIDI visibilities and SED can be reproduced by a spherical dust model provided that the density law is relatively flat. In this paper we present three additional VLTI/MIDI observations of W33A with baseline position angles close to perpendicular and with different projected lengths. This allows us to constrain the distribution of the emitting material further using 2-D axi-symmetric RT models.

We describe and discuss the new and archived MIDI observations in Sects. 2 and 3. The basic features of the RT model, our modelling approach, and the final preferred model are explained in Sect. 4. The implications for W33A and their effect for MYSOs in general are presented in Sect. 5. We conclude our work in Sect. 6

2 Observations and data reduction

Table 1: List of MIDI observations of W33A.

2.1 MIDI observations

W33A was observed with the VLTI and the MIDI instrument in service mode on four different occasions using the UTs as aperture stations (see Table 1). We aimed at covering different position angles on the shortest possible UT baselines. The baselines are presented in Fig. 1 projected onto the near-IR reflection nebula of W33A as imaged by UKIDSS (Casali et al. 2007; Warren et al. 2007; Lucas et al. 2008). The shortest VLTI baselines are attained with the ATs, but W33A cannot be observed in this mode because of its relatively low brightness. The obtained angular resolution for the data set presented here stretches from 25 to 60 milli-arcsec. Detailed descriptions of the MIDI observation procedure are given in Przygodda et al. (2003), Chesneau et al. (2005) and Leinert et al. (2004).

$\begin{figure} \par\includegraphics[height=9cm,width=8cm]{13209fg1.eps} \vspace*{-3mm} \end{figure}$	Figure 1: UKIDSS K-band observations of W33A. Projected VLTI/MIDI baselines are indicated by white bars.
Open with DEXTER

The W33A observations were executed in the so-called High-Sens MIDI mode, which uses all the incoming light for beam combination and fringe tracking. The incoming beams were combined producing two complementary interferometric channels that have by definition a phase difference of $\pi$ radians. The uncorrelated flux spectrum is measured immediately after the interferometric observation in order to obtain final visibilities. The high-sens mode is advantageous when observing faint targets or low visibilities, however the accuracy of the final visibility spectrum is limited by the sky brightness variation between the interferometric and photometric measurement. This can amount to 10-15% uncertainty, depending on the atmospheric stability. During the nights, the unresolved interferometric calibrator star W33A was observed in order to determine the instrumental visibility and correct for it. We also inspected all other observed calibrators from different programmes for each night (a total of 18 objects), which allowed us to estimate the stability and accuracy of the instrumental visibility. We found it to be uncertain by not more than 5%, much smaller than the uncertainty in the determined photometric measurement (20%, see next section). The calibrator star also provides an approximate flux calibration (Cohen et al. 1999) in addition to the visibility calibration. A prism with a spectral resolution of 30 was employed to disperse the incoming beams.

The correlated flux spectra were estimated by the MIA+EWS software package (version 1.6; see Jaffe 2004; Koehler 2005). Interferograms were coherently added to maximise the signal-to-noise. First, fringe spectra had to be corrected, because they are affected by the instantaneous atmospheric and instrumental piston. The fringe scan was used to estimate the group delay due to the atmosphere. Removing the atmospheric and the (known) instrumental group delays constitutes a correction to the dispersed fringe signal, and straightens the dispersed fringe spectra; i.e., the phase is independent of wavelength. Next, the phase offset due to varying water refraction index between the time of recording of the fringe spectra had to be accounted for. In principle, all spectra can then be added to a final fringe spectrum and the amplitude of the power at each frequency estimated, which will give the correlated flux.

The total flux spectra were taken in a chopping mode and each observation resulted in two spectra as MIDI splits the beams. The spectra were extracted from the raw data in the EWS mode. An estimate of the sky background was made and subtracted. The area corresponding to the mask with which the correlated spectra were extracted is used here. The counts were then simply summed in the y-direction. The final spectrum is the sum of the two geometric means of the four individual spectra: $\sqrt{A1 \times B1}+\sqrt{A2 \times B2}$ . This quantity is also what was obtained for the correlated flux after beam combination and thus ensures consistency in deriving visibilities.

2.2 JCMT/HARP-B observations

W33A was observed on the night of 5 June 2007 in ¹²CO, ¹³CO, and C¹⁸O using the HARP-B instrument on the James Clerk Maxwell Telescope. The ¹²CO data have a resolution of 488 kHz, and the other two lines were rebinned to the same resolution, giving a velocity resolution of approximately 0.4 km s^-1. The maps had a linear baseline removed. The core emission was presumed to be traced by the C¹⁸O emission (which has line profile very close to a Gaussian). We used this as a model to trace the emission in the ¹²CO and ¹³CO maps that arises from the actual outflow. Full details of this work, which is part of a larger project studying outflows from the Red MSX Source Survey (see e.g. Urquhart et al. 2008), will be published separately in due course.

3 Description of observational data of W33A

3.1 N-band from MIDI

The measured MIDI observables (flux, correlated flux) and derived visibilities of W33A are presented in Fig. 3. Each column represents a VLTI configuration, and we have included the observations already presented in Paper I (configuration D). The top panels show the flux-calibrated N-band spectra. They show a variation of around 20% as expected for data taken in the observational high-sens mode. The errorbars represent two uncertainties that are added in quadrature: (1) the rms estimates of the flux; and (2) the systematic error determined from a total of five W33A MIDI spectra. (More than the standard one flux spectrum per interferometric measurement were taken). The distinct feature of W33A in the N-band is the extremely deep silicate absorption, the actual depth is not detected at the central wavelength where it drops to less than 1% times the level of the pseudo-continuum.

The second row of panels shows the correlated flux spectra. The detector noise level was estimated by processing that part of the detector directly adjacent to the W33A fringes, indicative of the background noise. The noise produces a signal a few times 0.1 Jy around 10 $\mu$ m, despite the absence of incoming flux at the base of the silicate absorption. Similar noise levels have been found previously (e.g. Jaffe 2004; Matsuura et al. 2006). The errorbars of the W33A correlated flux spectra reflect the uncertainties due to spurious correlated noise level. We discard any signal below this noise level in further analysis.

The third row of figures presents the resulting visibilities for the wavelengths not affected by low flux levels. The errorbars are the uncertainties from correlated and total flux added in quadrature. This implies that the systematic uncertainty inherent to the high-sens mode are taken into account in the visibility uncertainties. We find that W33A has visibility values between 5% and 15% on the used baselines. The spectral shape of the visibilities reveal a clear trend toward decreasing visibilities on the blue wing of the silicate absorption feature. The lower visibilities imply a larger emitting region, hence that the increasing optical depth due to the silicate absorption comes from structures on larger size scales than the pseudo-continuum. The visibilities corresponding to the two similar configurations (B and C) are practically the same, as are their correlated flux levels. It gives an independent indication of the stability of the set-up and of the relatively good weather. Configuration A shows a correlated flux level twice higher than the ones of configuration B and C, despite the relatively little difference in total flux. Given that the weather during the observations was stable, we believe that the difference in visibilities between the various configurations is real. If we translate the visibilities at 8 $\mu$ m to Gaussian emission regions, then the FWHM decreases from 115 AU to about 95 AU rotating in PA from 15 $^\circ$ (configuration A) to $\sim$ $120\hbox{$^\circ$ }$ (B and C) and changing the baseline length accordingly.

$\begin{figure} \par\includegraphics[height=9cm,width=7cm,angle=90]{13209fg2.eps} \end{figure}$	Figure 2: $\rm ^{12}CO~(3-2)$ map taken with JCMT/HARP illustrating the molecular outflow of W33A. The PA and orientation of the blue and red-shifted lobes of the outflow are consistent with the K-band reflection nebula of Fig. 1. Coordinates are in J2000, and they correspond to the UKIDSS near-IR source and accord with the Capps et al. (1978) main IR source.
Open with DEXTER

3.2 The molecular outflow and near-IR reflection nebula

The star formation activity in W33A produces an archetypical bipolar molecular outflow (Fig. 2). The southeastern lobe constitutes the approaching, blue-shifted flow. Detailed millimetre line maps that would provide constraints on the outflow opening angle and the inclination of the system are currently unavailable. However, outflow activity produces reflection nebulae particularly well traced by near-IR scattered light (see e.g. Tamura et al. 1991). The near-IR image of W33A already presented in Fig. 1 is consistent with the CO map in position angle and orientation. We used this image therefore in order to constrain these system parameters, which will further help in modelling the MIDI visibilities and SED later on. The full extent of the $\sim$ $50\hbox{$^{\prime\prime}$ }$ long K-band nebula is shown in Fig. 1. The nebula has a position angle of $\sim$ $145\hbox{$^\circ$ }$ , but patchy extinction distorts the morphology on scales of 1 $^{\prime \prime }$ . A fainter extension can also be seen in the northwestern direction. This lends support to $145\hbox{$^\circ$ }$ position angle to be the actual orientation of the system. The high-angular resolution IFU images presented in Davies et al. (2010) confirms the reflected nature of the K-band emission at the base of the outflow and position angle. The outer regions of the nebula may have contributions from shocked $\rm H_{2}$ emission as traced by extended 4.5 $\mu$ m emission in Spitzer IRAC images (Cyganowski et al. 2008).

In Davies et al. (2010), spectro-astrometry of the Br $\gamma$ transition is also used to show that the near-IR source is at the base of the bipolar outflow on scales of 200 $\mu$ -arcseconds. The near-IR source most likely corresponds to the central object in W33A, as was assumed during the execution of the MIDI observations. This source shows a K-band profile in the UKIDSS image that is broader than the unresolved stars in the nearby field. In the H-band, the scattering nebula is much weaker (see Fig. 4). Given that photons are more efficiently scattered at shorter wavelengths the relatively dim nature of the H-band nebula argues for strong extinction in the line-of-sight towards W33A. Moreover, the central source is absent in the H-band.

$\begin{figure} \par\includegraphics[height=18cm,width=14.5cm,angle=90]{13209fg3.ps} \end{figure}$	Figure 3: MIDI observations of W33A on four different occasions (see Table 1). Top: flux spectra. Middle: correlated flux spectra (points with errorbars), and the detector noise level (dashed curve). Bottom: derived visibilities.
Open with DEXTER

Most of the emission seen at K-band is dust scattered. Contribution to the emission by molecular hydrogen gas is small as proven by intermediate-resolution spectra, obtained as part of the rms survey. In the H-band a foreground star is visible projected onto the nebula, located 2.5 $^{\prime \prime }$ south and 1.0 $^{\prime \prime }$ east of the main K-band source. We re-measured the flux densities coming from the full extended nebula as seen on the UKIDSS images, after first removing all the interloping stars from the images and replacing the pixels values using the adjacent regions. We find that in the K-band the nebula is 0.39 Jy ( $K \sim 8.0$ ) and in the H-band 0.025 Jy ( $H \sim 11.5$ ).

4 Modelling

In this section we present the simultaneous modelling of the visibility and spectral energy distribution of W33A. Historically, the SEDs of MYSOs have been interpreted as hot stars embedded in a spherically symmetric dusty envelope (e.g. Guertler & Henning 1986). In Paper I we adopted the approach of modelling the observables with a 1-D spherical symmetric radiative transfer code. We concluded that models with shallow density laws are preferred, provided that a relatively low effective temperature for the central star is adopted. Here we construct a fiducial model of a hot star surrounded by a protostellar envelope with outflow cavities. Our initial model aims to be as simple as possible, and we begin by investigating the imprint on N-band visibilities by a protostellar envelope that includes outflow cavities. The walls of the cavities are suggested to be the source of mid-IR emission seen on 1 $^{\prime \prime }$ scales in the nearly edge-on MYSO W33A (De Buizer 2006) and in W33A (de Wit et al. 2009).

4.1 Description of the radiative transfer code and input

We adopt the 2D-axi-symmetric dust RT model by Whitney et al. (see for details Whitney et al. 2003a,b). The model calculates radiation transfer (absorption, re-emission, and scattering) through a dusty structure. The inner rim of the structure is the dust sublimation radius that follows from a fit to models with different stellar temperatures (see Whitney et al. 2004: $R_{\rm subl}/R_{*}=(T_{\rm subl}/T_{*})^{-2.1}$ ); space is empty within the inner rim except for the star itself. The structure can include three distinct geometrical elements: an accretion disk, a protostellar envelope, and low-density polar outflow cavities. The protostellar envelope is described by the analytical TSC solution (Ulrich 1976; Terebey et al. 1984) of a simultaneously rotating and collapsing spherical structure. This solution is governed by the infall parameter ( $\dot{M}_{\rm infall}$ ) and the centrifugal radius $R_{\rm c}$ , i.e. the radius where rotational motion dominates infall. The envelope mass infall rate is a parameterization of the density distribution and does not necessarily represent an actual determination of infalling material. The presence of an accretion disk is optional in the code. When present, the disk structure follows the one for a flared $\alpha$ -type disk prescription, and the released accretion luminosity within the dust inner rim is added to the central luminosity source. The outflow cavity geometry can have two possible shapes, either a streamline or polynomial geometry. The streamline is conical on large scales and the polynomial has a $180\hbox{$^\circ$ }$ opening angle at the stellar surface. Each of these geometrical elements are potential mid-IR emitters and could contribute to the visibilities. For most model runs, the $R_{\rm subl}$ is comparatively small (tens of AU) and unresolved by the employed VLTI baselines in the N-band.

$\begin{figure} \par\includegraphics[height=12cm,width=7cm]{13209fg4.eps} \end{figure}$

Figure 4:

Left: UKIDSS H-band ( upper panel) and K-band ( middle panel) observations of the base of the W33A outflow scattering nebula. Images are logarithmically scaled, in units of MJy/sr and rotated by $35\hbox {$^\circ $ }$ so that the main lobe is pointing downwards. H-band contours are at 7.5%, 15%, 40% and 70% of peak flux; the K-band contours are at 1%, 2.5%, 10%, and 40% of peak flux. Secondary stellar sources have been removed from the images. The lower panel shows the H-K colour images (in magnitudes). Right: model images produced with the RT scattering code, see Sect. 4.1.2. The model images are convolved with 2-D Gaussian function with a FWHM of 1.6 $^{\prime \prime }$ . All images cover an area of $25\hbox {$^{\prime \prime }$ }$ on the side.

Open with DEXTER

A model run is necessarily characterised by a large number of free parameters. A fraction of these can be constrained by the known properties of W33A. The spatial information at milli-arcsecond angular resolution provided by MIDI allows the RT model to be seriously tested on its correctness in describing an MYSO.

4.1.1 Fixed input parameters

The central star parameters are chosen close to an O7.5 ZAMS according to the scale of Martins et al. (2005). The bolometric luminosity is $10^{5}~L_{\odot}$ (Faundez et al. 2004), and it is located at the kinematical distance of 3.8 kpc (Bronfman et al. 1996). The luminosity is based on IRAS flux density, and could be overestimated. We take the opportunity here to correct our 70 $\mu$ m Spitzer MIPS photometry presented in Paper I, which was not corrected for detector nonlinearity. Taking this correction into account we find a value of $1.9 \times 10^{3}$ Jy, quite similar to the IRAS measurement at 60 $\mu$ m. This indicates that the luminosity based on IRAS is probably reasonable. We think it is unlikely that a significant fraction of the IR luminosity is because of the secondary millimetre source noted in van der Tak et al. (2000) at 6 $^{\prime \prime }$ from the main source. The W33A MIPS 24 $\mu$ m and 70 $\mu$ m data show no sign of a secondary object, any more than a new 24.5 $\mu$ m image taken with Subaru/COMICS at an angular resolution of 0.6 $^{\prime \prime }$ (in prep.). This suggests that the IR luminosity is dominated by a single object.

A number of the important parameters that determine the properties of the TSC envelope are fixed for the following reasons. We have chosen a single type of dust that consists of an interstellar mixture of ``warm silicates'' (Ossenkopf et al. 1992) and of amorphous carbon with a size distribution according to Mathis et al. (MRN; 1977). The type of silicates provides a better fit to the observed shape of the absoprtion feature than MRN DL grain-mixture (Draine & Lee 1984), in particular to the blue wing (see Paper I, and see also Capps et al. 1978; Guertler & Henning 1986). The centrifugal radius is set at 33 AU which is comparable to the dust sublimation radius. The TSC envelope makes a transition from a pure infall, -1.5 radial density law to a shallower rotating, -0.5 law at the centrifugal radius. A TSC envelope with a centrifugal radius smaller than the sublimation radius is dominated by the infalling part of the envelope's structure solution; i.e., the envelope is quite spherically symmetric apart from the outflow cavities. The shape of the outflow cavities is polynomial (see Whitney et al. 2003a), which in the presence of a stellar wind could be a more appropriate geometry than a conical streamline. The density in the outflow cavity is chosen to be typical of MYSO outflows, i.e. between $5 \times 10^{3}$ and a few times 10⁴ cm^-3 (see Beuther et al. 2002). The cavity material consists of gas and dust. We initially base the outer edge of the envelope on 350 $\mu$ m observations by van der Tak et al. (2000). The radial intensity profile shows a downturn in submm emission at a radial distance of $\sim$ $50\hbox{$^{\prime\prime}$ }$ , or $2 \times 10^{5}$ AU. We require a posteriori that the model densities at the outer edge are similar to that of the large-scale molecular cloud material in which W33A is embedded, viz. 10⁴ cm^-3.

4.1.2 Cavity opening angle

One can exploit the observed shape of the near-IR scattering nebula to constrain the opening angle of the outflow. To this end we use the less computationally intensive scattering-only version of the dust RT code (Whitney & Hartmann 1992; Whitney & Hartmann 1993). The scattering code provides images at ten inclination angles $\mu$ simultaneously, which one can readily compare to the observations. We model the near-IR nebula under the assumption that the emission is free of direct thermal emission. We put particular emphasis on matching the outflow opening angle with the observed one and constraining the inclination angle. Model input parameters conform to the fixed parameters discussed in the previous subsection, except for the dust properties. The detailed dust properties are not important for the goal to constrain the geometrical parameters of the W33A system from the near-IR nebula. For completeness, the scattering dust used here is coated with a thin layer of water ice. The dust sizes follow a distribution that fits the extinction curve from (Cardelli et al. 1989) for a prescription of $R_{\rm V} = 4.0$ . This distribution is not a single power law, and it is presented in (Kim et al. 1994). For more details on the dust properties see Fig. 6 in Whitney et al. (2003a).

$\begin{figure} \par\mbox{\includegraphics[height=9cm,width=7cm,angle=90]{13209fg... ...} \includegraphics[height=9cm,width=7cm,angle=90]{13209fg6.eps} }\end{figure}$

Figure 5:

Left: SED of the preferred model overplot W33A flux densities. Small filled circles trace the ISO SWS spectrum, small open circles are the MIDI spectrum, the filled square is the Spitzer MIPS point, open square is the IRAS measurement at 60 $\mu$ m. In the (sub)mm range, the integrated fluxes at 350 $\mu$ m, 850 $\mu$ m, and 1.2 mm are from van der Tak et al. (2000), Di Francesco et al. (2008), and Faúndez et al. (2004), respectively. Right: zoom in onto the Nband spectrum. Shown are the MIDI flux spectrum with errorbars and the model fluxes appropriate for the MIDI slit width of 0.6 $^{\prime \prime }$ .

Open with DEXTER

$\begin{figure} \par\includegraphics[height=18cm,width=15cm,angle=90]{13209fg7.eps} \end{figure}$	Figure 6: Model (full line and filled circles) fits to the observed visibilities. Model parameters are listed in Table 2. Corresponding images are presented in Fig. 8. Note the discrepancy with the spherically symmetric model (dotted lines) used to fit configuration D in Paper I.
Open with DEXTER

Figure 4 compares the H and K-band images of a model run with the W33A scattering nebula. The images of the observed nebula have been rotated by $35\hbox {$^\circ $ }$ . The observed and model images are in the same units (MJy sr^-1). Model images are convolved with a Gaussian point spread function that simulates the one derived from the UKIDSS images, and Gaussian noise derived from the observed images is added to the model ones. The outflow angle measured from the images directly would be $45\hbox{$^\circ$ }$ ; however, the presented model has an envelope geometry with a full cavity opening angle of less than half that, viz. $2\theta=20\hbox{$^\circ$ }$ . In our fitting procedure, attention is devoted to the absence of an H-band counterpart of the central source. Whether or not a central source is detectable depends on the system inclination. The model that is presented in Fig. 4 has an inclination of $50\hbox{$^\circ$ }$ . An inclination less than $40\hbox{$^\circ$ }$ or more than $70\hbox{$^\circ$ }$ are excluded. In the former case, the central star is visible in the H-band, whereas the redshifted outflow lobe becomes prominent in the latter. These features are not present in the UKIDSS near-IR images.

Although the shape on the sky can be reproduced reasonably well, not a single model run reproduces the observed H-K colour of $\sim$ 3.5. Photometry done on the model images shows that they produce colours that are much bluer instead. The model images of Fig. 4 match the colour of the observed images, but additional foreground extinction of $A_{\rm V} \approx 25$ is required. The model H and K-band fluxes need to be scaled up accordingly. The central source visible in the K-band only constitutes a minor contribution to the total observed flux density at this wavelength. We note that MYSO models invariably have problems matching the short wavelength fluxes of MYSO scattering nebulae (e.g. Alvarez et al. 2004). The deviant model colours might be rectified by including a clumpy structure within the outflow cavity (Indebetouw et al. 2006) and/or very small dust particles that are transiently heated to high temperatures. We find that introducing an outflow cavity with a smooth density distribution into an infalling envelope cannot be the only physical element determining the observed near-IR emission of MYSOs. Regardless of the solution, a deeper investigation of this problem is beyond the purpose of this paper. Nonetheless, a comparison of near-IR images with scattering RT models has led to relatively stringent constraints on the outflow opening angle and comparatively less stringent limits on the inclination.

4.1.3 Methodology

The prime model parameters that remain to be determined are the density as parametrized by the envelope mass infall rate, the presence of an accretion disk (and its properties), and the outer radius. They determine the peak of the SED and, depending on inclination angle, the depth of the silicate feature, and the dominant emission regions on scales of 100 AU. The RT code is computationally demanding which excludes a grid search for the optimum value for the decisive model parameters. Building on a simple fiducial model, we try to find a model that fits the SED and MIDI visibilities simultaneously. One run delivers the SEDs at 10 inclinations simultaneously, evenly spaced in $\cos(\mu)$ . Model images are produced for a single inclination at a time but for any number of desired wavelengths/filters. A total of $4 \times 10^{7}$ ``photons'' are used to produce images, for SED-only-runs the code uses four times less. We divided the N-band up in 6 square filters of 1 $\mu$ m width, running from 7.5 $\mu$ m to 12.5 $\mu$ m. The resulting narrow-band model images are multiplied with a Gaussian model for the VLT-UT Airy disk. This represents the aperture field of view with an FWHM appropriate to the particular wavelength. The resulting image is then Fourier-transformed and the visibilities extracted for the projected interferometer baseline and position angle. The model images have a pixel size of $0.01\hbox{$^{\prime\prime}$ }$ , which is less than half the angular resolution (0.024 $^{\prime \prime }$ ) at the longest baseline for the shortest wavelength. The SED is matched from $\sim$ 8 $\mu$ m to 1 mm. The input data are the MIDI flux spectrum, the slope of the SED around 40 $\mu$ m determined from the ISO-SWS spectrum (see Paper I) and the MIPS 70 $\mu$ m flux density. At the longer wavelengths, the integrated fluxes from the following observations are used: 350 $\mu$ m from with the CSO/SHARC instrument (van der Tak et al. 2000), 850 $\mu$ m with JCMT/SCUBA (Di Francesco et al. 2008), and 1.2 mmm with SEST/SIMBA (Faúndez et al. 2004). The far-IR and (sub)mm observations are compared to model SEDs representing all the emitted energy at each wavelength, whilst the MIDI flux spectrum is compared to the model SED with an aperture of 0.6 $^{\prime \prime }$ corresponding to the MIDI slit width. Matching of the SED and visibilities is performed by eye.

Table 2: Parameters for the preferred model discussed in Sect. 4.2.

$\begin{figure} \par\mbox{\includegraphics[height=9cm,width=8cm,angle=90]{13209fg... ...} \includegraphics[height=9cm,width=8cm,angle=90]{13209fg9.ps} } \end{figure}$	Figure 7: Logarithmically scaled slices through the spatial distribution of the density ( left) and temperature ( right) of the preferred model perpendicular to the outflow axis. The cavity has a very small, but non-zero density consisting of gas and dust (see text).
Open with DEXTER

$\begin{figure} \par\includegraphics[width=18cm,angle=0]{13209fg10.eps}\par\includegraphics[width=18cm,angle=0]{13209fg11.eps} \vspace*{4mm} \end{figure}$

Figure 8:

Logarithmically scaled images ( top row) and corresponding Fourier transforms ( bottom row) of the preferred model (Table 2) at wavelengths 7.5 $\mu$ m to 12.5 $\mu$ m in steps of 1 $\mu$ m. Image field of view is $0.6\hbox {$^{\prime \prime }$ }$ (2280 AU) on a side, corresponding to the MIDI slit width. Pixel scale is 0.01 $^{\prime \prime }$ . The model images have been multiplied with a wavelength dependent 2-D Gaussian function, simulating the VLT-UT airy disk. Warm cavity walls dominate the emission at 10 $\mu$ m.

Open with DEXTER

4.2 Preferred model fit

With the parameter values described in the previous section, a model is obtained that produces a reasonable match to both the SED and the visibilities. We ran a total of some 30 models in order to narrow down the parameter values by varying the mass infall rate, inclination angle, and the outer radius. A satisfactory model is presented in Figs. 5 and 6. The parameter values of this particular model can be found in Table 2. The figures show that the model fits both the SED and visibilities reasonably. A discrepancy between model and observed SED is found in the near-IR, as expected (see Sect. 4.1.2). The model reproduces the spectral trend of smaller visibilities in the wings of the silicate feature, especially clear on the blue wing quite well. Figure 6 also displays the visibilities from the spherically symmetric solution that matches the SED and visibilities of configuration D only (from Paper I). Each panel shows the spherically symmetric model visibilities calculated for the appropriate baseline length. It is clear that they provide a considerably worse fit to the visibilities obtained on the additional configurations A, B, and C, illustrating the need for multiple baselines.

The actual model images of the preferred model are shown in Fig. 8, along with a slice through the spatial temperature and density distribution in Fig. 7. They make clear that the dominant emission regions in the N-band at a few 100 AU scale are the envelope regions close to the surface of the outflow cavity. These envelope regions are warmed up because of irradiation by the star. We refer to these regions as ``cavity walls'' for brevity. The density distribution in Fig. 7 illustrates that the density of the paraboloidal cavity is a few orders of magnitude less than that of the envelope. The cavity itself is hotter than the bulk of the envelope but contributes very little to the N-band emission due to the low density. Figure 8 shows that the emission is relatively extended on the probed scales and produces consequently very low visibilities. An increase in optical depth, by either probing a higher inclination angle or increasing the mass infall rate, extinguishes the warmest, inner part of the irradiated walls. As a result visibilities generally decrease with an optical depth increase. However the extremely deep silicate absorption feature of W33A requires a large optical depth. This is the trade off for obtaining a decent fit. We find that for an inclination angle larger than $35\hbox {$^\circ $ }$ , the mass infall rate needs to be lower than $1 \times 10^{-3}~\mbox{$M_{\odot}$ ~yr$^{-1}$ }$ in order to see the warm cavity wall regions and avoid spatially resolving the N-band emission out. The increase in optical depth because of the silicate absorption has a similar effect. Figure 8 makes clear that the warm inner regions of the cavity walls become progressively fainter with wavelength causing an increase in the typical emitting size, hence lower visibilities.

The same effect is found when one increases the centrifugal radius. The preferred model with an increased $R_{\rm c}$ to a value of 300 AU shows that the density in the equatorial region increases, but the increase is confined to within a fraction of a degree of the equator. The density actually decreases away from the equator allowing a larger area of the cavity walls to be heated up. Visibilities drop therefore with increasing $R_{\rm c}$ . Further fine-tuning of the model parameters would require additional observations at preferentially shorter baselines.

In summary, MIDI visibilities and SED of W33A can be matched with a simple geometry that consists of a protostellar envelope with a -1.5 radial density distribution that has polar outflow cavities with a $2\theta=20\hbox{$^\circ$ }$ opening angle. The modelling of the MIDI data is not very sensitive to the exact density power law in the envelope within the framework of the applied RT model. The density power law in this case is effectively governed by the centrifugal radius. Spherical models as applied to W33A in de Wit et al. (2007) are inadequate for modelling structure on 100 AU scales, as Fig. 5 illustrates. On the other hand, resolved observations of MYSOs at $\sim$ 25 $\mu$ m at average inclination angles (de Wit et al. 2009) and at $350~\mu$ m (van der Tak et al. 2000) provide insight into the density law starting at radii of 1000 AU and indicate a radial density law of -1.0. Matching images of W33A with model images at these respective wavelengths will provide much stronger constraints on the centrifugal radius. For now we show the 350 $\mu$ m morphology of our model in Fig. 9. We convolved the model image with a 10 $^{\prime \prime }$ beam appopriate for comparison with the CSO observations presented in van der Tak et al. (2000). These observations show a structure that is more flattened than the model image; however, the FWHM of the model image of 16 $^{\prime \prime }$ is bracketed by the observed FWHM range of 14.4 $^{\prime \prime }$ -19.4 $^{\prime \prime }$ , showing a good correspondance. The density at the outer radius of the model complies with expected densities of molecular clouds. Observations in the submm give a lower limit to the size, while in reality the outer parts of the envelope may become even steeper than -1.5(e.g. Mueller et al. 2002). The outer radius is therefore uncertain by a factor of 2. The best-matching RT model presented does not need contributions to the emission from other components within the MYSO environment; instead, N-band emission on a 100 AU scale is completely dominated by the warm cavity walls.

5 Discussion

In several studies the accretion based formation scenario for massive stars has been forwarded partially based on observations with 1 $^{\prime \prime }$ angular resolution and/or SED studies of MYSOs (e.g. De Buizer et al. 2005). If MYSOs are to be considered young massive stars still in their accretion phase, and this accretion occurs in a similar fashion to low-mass stars, then circumstellar accretion disks ought to be present. The MYSO accretion disk is expected to be luminous (approaching that of the star itself), and is estimated to have a size a few times 100 AU. Our analysis raises the question to what degree an accretion disk can contribute to the mid-IR emission.

5.1 The contribution of an accretion disk

The dust RT code allows the inclusion of a dust disk truncated at the 25 AU dust sublimation radius. As a first step we consider the effect on the model fit by adding such a dusty disk to the preferred model presented previously (Table 2). We emphasize that the disk is directly irradiated by the star, as the space between disk and star is empty in the model. The mass of the dusty disk is initially chosen to be 1% of the stellar mass, i.e. 0.25 $M_{\odot}$ . The radial density law of the disk agrees with that expected of an $\alpha$ -type disk (i.e. a power law exponent of 1.875, Frank et al. 1985). Because of the large inner truncation radius, the actual accretion luminosity generated by the dust disk is marginal and about 1% of the stellar luminosity (see Eq. (6) in Whitney et al. 2003a). The formal accretion rate (see Eq. (5) in Whitney et al. 2003a) for $\alpha=0.01$ is then $2.8 \times 10^{-6} (\alpha/0.01)~\mbox{$M_{\odot}$ ~yr$^{-1}$ }$ . The accretion luminosity inside the dust sublimation radius, where no actual disk structure is present in the RT code, is emitted with the stellar spectrum (see Whitney et al. 2003a). The disk has a constant opening angle; i.e., the scaleheight increases linearly with radius. The model disk's outer radius is arbitrary because it does not emit in the N-band and is better constrained by millimetre data; for now, it is set to 500 AU.

The resulting model visibilities of a protostellar envelope plus dust disk are shown in Fig. 10. It is clear that such a model violates the observed MIDI visibilities shortward of the silicate absorption feature. At these wavelengths the dust disk contributes $\sim$ 50% to the total flux density, and one can deduce that the disk has a visibility of approximately 70%. For this model to be consistent with the observations, we computed the fractional contribution to the total 8 $\mu$ m flux by the dust disk as probably at most 1%, if the disk emission at this wavelength is fully unresolved. If it is marginally resolved (as it is), the 1% upper limit becomes less. By running models with increasingly lower disk masses, we find that this flux fraction is reached for a total mass of the dusty disk of 0.01 $M_{\odot}$ , or equivalently $1.1 \times 10^{-7} (\alpha/0.01)~\mbox{$M_{\odot}$ ~yr$^{-1}$ }$ . The SED of W33A is still fitted well by the dust disk+envelope model. We thus find that the disk mass is very low compared to the envelope mass and, for reasonable values of $\alpha$ , that the dust disk accretion rate is more than three orders of magnitude lower than the envelope mass infall rate. These results on the disk masses are in accordance with other studies of MYSOs (e.g. Murakawa et al. 2008) and also with some evolved early-type Herbig Be stars (Alonso-Albi et al. 2009). Were the disk outer radius much smaller than the assumed 500 AU and similar to the small outer radii found for the Herbig Be stars, then the disk mass should be scaled down accordingly: a 50 AU outer radius results in a disk mass of $5 \times 10^{-4}~\mbox{$M_{\odot}$ }$ . A small disk mass may imply that such a dust disk is not present, or that it is simply extremely faint in the mid-IR because of the absence of gas and dust opacities in the model disk, which would lead to strong shielding from the stellar radiation field.

$\begin{figure} \par\includegraphics[height=8.5cm,width=8cm,angle=90]{13209fg12.eps} \end{figure}$	Figure 9: Morphology at 350 $\mu$ m predicted by the preferred model. Contours are at 15%, 30%, 50%, 70% and 90% of the peak value. Model image has been convolved with a 10 $^{\prime \prime }$ beam. See van der Tak et al. (2000) for the observed 350 $\mu$ m image.
Open with DEXTER

$\begin{figure} \par\includegraphics[height=8.5cm,width=7cm,angle=90]{13209fg13.eps} \end{figure}$	Figure 10: Model predictions for baseline configuration A. Represented are the preferred model (full line), with the addition of a dust disk ( $M=0.25~\mbox{$M_{\odot}$ }$ , dashed line), and an optically thick disk interior to the dust sublimation radius (dotted line).
Open with DEXTER

Most of the accretion luminosity, however, would be released within the 25 AU dust sublimation radius by the hot gaseous extension of the dust disk. Self-shielding means that such a disk probably contains dust at smaller radii than the 25 AU sublimation radius determined from direct irradiation. Inclusion of a self-luminous accretion disk within the dust sublimation radius is required to fit the SED and the near-IR interferometric data of the Herbig B6e star W33A (Kraus et al. 2008). Accretion luminosity generated by an accretion rate of $7\times 10^{-6}~\mbox{$M_{\odot}$ ~yr$^{-1}$ }$ for this particular object is, however, not essential for modelling the star's mid-IR interferometry, which is already reproduced reasonably well by a passive disk. We explored to what extent the mid-IR interferometry of W33A is affected by including an optically thick accretion disk. We approximated the reslulting fluxes and visibilities with an optically thick, geometrically thin $\alpha$ disk model (see Malbet & Bertout 1995; Malbet et al. 2007) and added these to the preferred envelope model in Table 2. The disk is heated only by stationary accretion and for an accretion rate equalling the mass infall rate ( $7.5 \times 10^{-4}~\mbox{$M_{\odot}$ ~yr$^{-1}$ }$ ) the generated accretion luminosity is $\sim$ $7 \times 10^{4}~\mbox{$L_{\odot}$ }$ for our adopted central star. The spatial properties of the optically thick disk emission is such that 99% of the N-band is emitted within 50 AU. Most of the 8 $\mu$ m flux is emitted at 20 AU. The disk thus remains largely unresolved on the employed VLTI baselines with MIDI. Disks accreting at a slower rate than the one assumed are cooler and the mid-IR emitting region more compact. We chose the inner radius of the disk to be at 1 R_*. The disk temperature falls below 1500 K at 2.5 AU, beyond which the disk opacity is given by both gas and dust.

As in the previous paragraph, if the preferred model visibilities are to be affected by the optically thick disk, the ratio of disk emission to the emerging envelope emission should be at least 1% at 8 $\mu$ m. We calculated the resulting disk spectrum for the said disk accretion rate and applied an extinction of $A_{\rm V} = 230$ due to the overlying envelope, as found in the line-of-sight to W33A from the preferred model. We show the resulting visibilities in Fig. 10. The figure shows that an $\alpha$ disk model embedded in a protostellar envelope is compatible with the observed MIDI visibilities. The SED for this particular model does not differ from the SED with only an envelope contributing to the emission, provided that the total luminosity of each model is chosen the same. The $\alpha$ disk contributes less than 1% to the total 8 $\mu$ m flux, despite its high luminosity. The envelope of W33A is simply too opaque for any disk emission to be directly observed. MIDI observations of MYSOs with less massive envelopes or at smaller inclinations, however, may well disclose the presence of an accretion disk in such systems. The expected visibilities would be much higher than the ones observed for W33A (depending on distance and accretion rate), considering the compactness and brightness of the 8 $\mu$ m emission in accretion disks.

The popular YSO SED model grid developed by Robitaille et al. (2006) allows an estimate of the extent to which MYSO SEDs are dominated by dust disk emission. This YSO SED grid makes use of the same 2-D axi-symmetric RT calculations as we use in this paper (Whitney et al. 2003a), and the disks in the SED model grid are therefore dust disks. We extract all (135) SEDs corresponding to ZAMS hot stars with a luminosity of $10^{3} < L/L_{\odot}<10^{6}$ and an envelope mass infall rate of at least $10^{-4}~M_{\odot}~{\rm yr}^{-1}$ . Inspection of the fractional disk flux contribution at 8 $\mu$ m (just short of the extended wings of the silicate absorption feature) shows that all except three of the extracted models have a fractional contribution over 15%, and half of the models have a contribution over 70%. The geometry of the accretion disk determines the total visibility at these wavelengths. The SED model grid thus contains a preponderance of MYSO models with relatively high visibilities, similar to the envelope plus dust disk model presented in Fig. 10. At present, the MIDI observations of MYSOs demonstrate the opposite trend, however. Visibilities presented for W33A are in line with the low N-band visibilities observed for other MYSOs (see Linz et al. 2008; Vehoff et al. 2008; de Wit et al. in prep.). The emission of these objects could equally well be explained partially or fully in the context of cavity wall emission as proposed here for W33A.

Table 3: Models from the YSO SED model grid (Robitaille et al. 2006), obtained using the web-based fit procedure.

5.2 SED model grid fits

To what extent are models obtained from SED fits consistent with the spatial information on milli arcsecond scales? To address this question, we applied the web-based SED fitter procedure developed by Robitaille et al. (2007) to the W33A SED data shown in Fig. 12. The procedure returns a ranking on the basis of a $\chi^{2}$ minimization method, and we discuss the grid models following this ranking.

$\begin{figure} \par\includegraphics[height=8cm,width=8cm]{13209fg14.ps} \end{figure}$	Figure 11: Logarithmically scaled K-band image for grid model 3004217 to be compared with the observations in Fig. 4. The contours are at 1%, 40%, 50%, and 65% of the peak flux. The image is convolved with a 2-D Gaussian function with an FWHM of 1.6 $^{\prime \prime }$ .
Open with DEXTER

The procedure delivers a model consisting of a central hot star and a dust disk as the best fit (model grid number 3002520, see Table 3). We see from Fig. 12 that the grid model fit to the SED is acceptable, as expected. This particular grid model has a system inclination of $50\hbox{$^\circ$ }$ , consistent with the one of the reflection nebula. The best-fit model consists of a protostellar envelope with a relatively large cavity opening angle compared to the one derived from the scattering nebula. Moreover, a fit of the model to the observed SED requires a distance of 1 kpc, which is a factor $\sim$ 4 too short. The discrepancy between the modelled and observed SED is mainly at the depth of the silicate feature. This depth is better accounted for by dust partially consisting of ``warm silicates'' (see Ossenkopf et al. 1992) as discussed previously, while the appropriate set of optical constants is not used in the SED model grid. When it is broken down in its various contributing components, we find that the grid model SED is completely dominated by disk emission at wavelengths shortward of the silicate absorption feature. The emission accounts for the relatively high flux levels observed shortward of the silicate feature. This is a common property of MYSO SEDs. Spherical and indeed 2D-axisymmetric RT codes for the protostellar envelope are unable to reproduce these high flux levels (see e.g. de Wit et al. 2009). Contrary to the good SED fit, Fig. 13 shows that the corresponding visibilities are far too high. The visibilities are obtained by producing images with the dust RT code using the same parameter file as produced the SED. The grid model visibilities are inconsistent simply because of the dominance of the (compact) disk emission. Clearly, a good model fit to the SED does not guarantee a proper description of the object in consideration, because the best-fit grid model violates the spatial information delivered by MIDI.

If we use the known properties of W33A and limit the model distance to be between 3.0 and 4.6 kpc with a maximum foreground extinction of A_v=50, then we find that the returned fitting models have a strong preference for small inclinations. In the best 10 models, 7 have an inclination of 18.19 degrees (the smallest present within the grid), also the best-fit grid model (number 3004217). In addition, the visibilities fit the MIDI observations reasonably well (Fig. 13). The model is characterized by a central cool star with a large stellar radius. Such a description of the MYSO M8E-IR is preferred in Linz et al. (2009). In the case of W33A, the small inclination angles are inconsistent with the allowed range we derived from the scattering nebula. The predicted K-band morphology of this model is shown in Fig. 11, and it is strongly dominated by the central source in contradiction with the observations of W33A (Fig. 4). The ranked fourth model however has a more reasonable inclination angle of $50\hbox{$^\circ$ }$ (number 3009737). The latter model consists of a hot star and an accretion disk and encounters the previously discussed problems in fitting the visibilities (Fig. 13). If we applied all the W33A fiducial parameters ( $40\hbox{$^\circ$ }<i<70\hbox{$^\circ$ }$ , cavity opening angle $10<2\theta<30$ ) in selecting a model from the grid for a hot star ( $T_{\rm eff}>20~000~K$ ) without a dust disk, then the model complying with these constraints is ranked 245th (model number 3003047).

$\begin{figure} \par\includegraphics[height=8.2cm,width=8.35cm,angle=90]{13209fg15.ps} \end{figure}$	Figure 12: Predicted SEDs of the discussed models obtained from the SED web fit procedure listed in Table 3.
Open with DEXTER

$\begin{figure} \par\includegraphics[height=8.2cm,width=8.35cm,angle=90]{13209fg16.ps} \end{figure}$	Figure 13: Predicted visibilities along the major (lower curve for each model) and minor axis (top curve for each model) at 8 $\mu$ m as a function of baseline for the discussed models obtained from the SED web fit procedure (Table 3). MIDI visibilities for the four baselines are also shown.
Open with DEXTER

It is clear that the SED model grid can deliver models that fit the SED, but all the SED fitting models presented in this section violate the MIDI visibilities. Regarding W33A, the SED grid does not provide a fit to both SED and visibilities simply because the appropriate models are not present in the grid (see also the discussion in Robitaille 2008; Linz et al. 2009).

6 Conclusions

Observations with MIDI at the VLTI of the MYSO W33A indicate that N-band emission on scales of 100 AU stems from the warm parts of the envelope close to the outflow cavity, i.e. the ``cavity walls''. When described by a TSC type geometry, the overlying envelope is characterised by an equivalent mass infall rate of order $10^{-3}~\mbox{$M_{\odot}$ ~yr$^{-1}$ }$ and is quite massive. We find that the presence of a dust disk located beyond the dust sublimation radius violates the observed visibilities. To make the model consistent with the observations, the dust disk must be marginal in mass. In contrast, the addition of an optically thick accretion disk inside the dust sublimation radius has little effect on the visibilities, even when the disk dominates the luminosity. The disk emission is highly extincted by the massive overlying protostellar envelope. MYSOs viewed at smaller inclination and/or with less massive protostellar envelopes have the potential of revealing such accretion disks.

Acknowledgements

We would like to thank O. Chesneau for his prompt advice, and T. Fujiyoshi for the W33A COMICS data. It is a pleasure to thank B. Whitney for discussions of the subject, and to an anonymous referee for the clear remarks that improved this manuscript.

References

Alonso-Albi, T., Fuente, A., Bachiller, R., et al. 2009, A&A, 497, 117 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Alvarez, C., Hoare, M., & Lucas, P. 2004, A&A, 419, 203 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Beltrán, M. T., Cesaroni, R., Neri, R., et al. 2004, ApJ, 601, L187 [NASA ADS] [CrossRef] [Google Scholar]
Beltrán, M. T., Cesaroni, R., Neri, R., et al. 2005, A&A, 435, 901 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Beuther, H., Schilke, P., Sridharan, T. K., et al. 2002, A&A, 383, 892 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Bronfman, L., Nyman, L.-A., & May, J. 1996, A&AS, 115, 81 [NASA ADS] [CrossRef] [EDP Sciences] [MathSciNet] [PubMed] [Google Scholar]
Bunn, J. C., Hoare, M. G., & Drew, J. E. 1995, MNRAS, 272, 346 [NASA ADS] [Google Scholar]
Capps, R. W., Gillett, F. C., & Knacke, R. F. 1978, ApJ, 226, 863 [NASA ADS] [CrossRef] [Google Scholar]
Cardelli, J. A., Clayton, G. C., & Mathis, J. S. 1989, ApJ, 345, 245 [NASA ADS] [CrossRef] [Google Scholar]
Casali, M., Adamson, A., Alves de Oliveira, C., et al. 2007, A&A, 467, 777 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Cesaroni, R., Galli, D., Lodato, G., Walmsley, C. M., & Zhang, Q. 2007, in Protostars and Planets V, ed. B. Reipurth, D. Jewitt, & K. Keil, 197 [Google Scholar]
Chesneau, O., Verhoelst, T., Lopez, B., et al. 2005, A&A, 435, 563 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Cohen, M., Walker, R. G., Carter, B., et al. 1999, AJ, 117, 1864 [NASA ADS] [CrossRef] [Google Scholar]
Cyganowski, C. J., Whitney, B. A., Holden, E., et al. 2008, AJ, 136, 2391 [NASA ADS] [CrossRef] [Google Scholar]
Davies, B., Lumsden, S. L., Hoare, M. G., Oudmaijer, R. D., & de Wit, W. J. 2010, MNRAS, 402, 1504 [NASA ADS] [CrossRef] [Google Scholar]
De Buizer, J. M. 2006, ApJ, 642, L57 [NASA ADS] [CrossRef] [Google Scholar]
De Buizer, J. M., Osorio, M., & Calvet, N. 2005, ApJ, 635, 452 [NASA ADS] [CrossRef] [Google Scholar]
de Wit, W. J., Hoare, M. G., Fujiyoshi, T., et al. 2009, A&A, 494, 157 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
de Wit, W. J., Hoare, M. G., Oudmaijer, R. D., et al. 2007, ApJ, 671, L169 [NASA ADS] [CrossRef] [Google Scholar]
Di Francesco, J., Johnstone, D., Kirk, H., MacKenzie, T., & Ledwosinska, E. 2008, ApJS, 175, 277 [NASA ADS] [CrossRef] [Google Scholar]
Draine, B. T., & Lee, H. M. 1984, ApJ, 285, 89 [NASA ADS] [CrossRef] [Google Scholar]
Dutrey, A., Guilloteau, S., Duvert, G., et al. 1996, A&A, 309, 493 [NASA ADS] [Google Scholar]
Faúndez, S., Bronfman, L., Garay, G., et al. 2004, A&A, 426, 97 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Frank, J., King, A. R., & Raine, D. J. 1985, Accretion power in astrophysics, ed. J. Frank, A. R. King, & D. J. Raine [Google Scholar]
Furuya, R. S., Cesaroni, R., Takahashi, S., et al. 2008, ApJ, 673, 363 [NASA ADS] [CrossRef] [Google Scholar]
Gibb, E. L., Whittet, D. C. B., Schutte, W. A., et al. 2000, ApJ, 536, 347 [NASA ADS] [CrossRef] [Google Scholar]
Grady, C. A., Woodgate, B., Bruhweiler, F. C., et al. 1999, ApJ, 523, L151 [NASA ADS] [CrossRef] [Google Scholar]
Guertler, J., & Henning, T. 1986, Ap&SS, 128, 163 [NASA ADS] [CrossRef] [Google Scholar]
Hoare, M. G. 2006, ApJ, 649, 856 [NASA ADS] [CrossRef] [Google Scholar]
Hoare, M. G., & Franco, J. 2007, Massive Star Formation, ed. T. W. Hartquist, J. M. Pittard, & S. A. E. G. Falle, 61 [Google Scholar]
Hosokawa, T., & Omukai, K. 2009, ApJ, 691, 823 [NASA ADS] [CrossRef] [Google Scholar]
Hummer, D. G., & Storey, P. J. 1987, MNRAS, 224, 801 [NASA ADS] [CrossRef] [Google Scholar]
Indebetouw, R., Whitney, B. A., Johnson, K. E., et al. 2006, ApJ, 636, 362 [NASA ADS] [CrossRef] [Google Scholar]
Jaffe, W. J. 2004, in SPIE Conf. Ser., ed. W. A. Traub, 5491, 715 [Google Scholar]
Jiménez-Serra, I., Martín-Pintado, J., Rodríguez-Franco, A., et al. 2007, ApJ, 661, L187 [NASA ADS] [CrossRef] [Google Scholar]
Keto, E. 2002, ApJ, 580, 980 [NASA ADS] [CrossRef] [Google Scholar]
Kim, S., Martin, P. G., & Hendry, P. D. 1994, ApJ, 422, 164 [NASA ADS] [CrossRef] [Google Scholar]
Koehler, R. 2005, Astron. Nachr., 326, 563 [NASA ADS] [Google Scholar]
Kraus, S., Preibisch, T., & Ohnaka, K. 2008, ApJ, 676, 490 [NASA ADS] [CrossRef] [Google Scholar]
Krumholz, M. R., Klein, R. I., McKee, C. F., Offner, S. S. R., & Cunningham, A. J. 2009, Science, 323, 754 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
Krumholz, M. R., McKee, C. F., & Klein, R. I. 2005, ApJ, 618, L33 [NASA ADS] [CrossRef] [Google Scholar]
Leinert, C., Graser, U., Richichi, A., et al. 2003, The Messenger, 112, 13 [NASA ADS] [Google Scholar]
Leinert, C., van Boekel, R., Waters, L. B. F. M., et al. 2004, A&A, 423, 537 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Linz, H., Henning, T., Stecklum, B., et al. 2008, in Massive Star Formation: Observations Confront Theory, ed. H. Beuther, H. Linz, & T. Henning, ASP Conf. Ser., 387, 132 [Google Scholar]
Linz, H., Henning, T., Feldt, M., et al. 2009, A&A, 505, 655 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Lucas, P. W., Hoare, M. G., Longmore, A., et al. 2008, MNRAS, 391, 136 [NASA ADS] [CrossRef] [Google Scholar]
Malbet, F., Benisty, M., de Wit, W.-J., et al. 2007, A&A, 464, 43 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Malbet, F., & Bertout, C. 1995, A&AS, 113, 369 [NASA ADS] [CrossRef] [EDP Sciences] [MathSciNet] [PubMed] [Google Scholar]
Martins, F., Schaerer, D., & Hillier, D. J. 2005, A&A, 436, 1049 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Mathis, J. S., Rumpl, W., & Nordsieck, K. H. 1977, ApJ, 217, 425 [NASA ADS] [CrossRef] [Google Scholar]
Matsuura, M., Chesneau, O., Zijlstra, A. A., et al. 2006, ApJ, 646, L123 [NASA ADS] [CrossRef] [Google Scholar]
Mueller, K. E., Shirley, Y. L., Evans, II, N. J., et al. 2002, ApJS, 143, 469 [NASA ADS] [CrossRef] [Google Scholar]
Murakawa, K., Preibisch, T., Kraus, S., et al. 2008, A&A, 490, 673 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Ossenkopf, V., Henning, T., & Mathis, J. S. 1992, A&A, 261, 567 [NASA ADS] [Google Scholar]
Osterbrock, D. E. 1989, Astrophysics of gaseous nebulae and active galactic nuclei, ed. D. E. Osterbrock [Google Scholar]
Padgett, D. L., Brandner, W., Stapelfeldt, K. R., et al. 1999, AJ, 117, 1490 [NASA ADS] [CrossRef] [Google Scholar]
Przygodda, F., Chesneau, O., Graser, U., Leinert, C., & Morel, S. 2003, Ap&SS, 286, 85 [NASA ADS] [CrossRef] [Google Scholar]
Rengarajan, T. N., & Ho, P. T. P. 1996, ApJ, 465, 363 [NASA ADS] [CrossRef] [Google Scholar]
Robitaille, T. P. 2008, in Massive Star Formation: Observations Confront Theory, ed. H. Beuther, H. Linz, & T. Henning, ASP Conf. Ser., 387, 290 [Google Scholar]
Robitaille, T. P., Whitney, B. A., Indebetouw, R., Wood, K., & Denzmore, P. 2006, ApJS, 167, 256 [NASA ADS] [CrossRef] [Google Scholar]
Robitaille, T. P., Whitney, B. A., Indebetouw, R., et al. 2007, ApJS, 169, 328 [NASA ADS] [CrossRef] [Google Scholar]
Tamura, M., Gatley, I., Joyce, R. R., et al. 1991, ApJ, 378, 611 [NASA ADS] [CrossRef] [Google Scholar]
Terebey, S., Shu, F. H., & Cassen, P. 1984, ApJ, 286, 529 [NASA ADS] [CrossRef] [Google Scholar]
Ulrich, R. K. 1976, ApJ, 210, 377 [NASA ADS] [CrossRef] [Google Scholar]
Urquhart, J. S., Hoare, M. G., Lumsden, S. L., Oudmaijer, R. D., & Moore, T. J. T. 2008, in Massive Star Formation: Observations Confront Theory, ed. H. Beuther, H. Linz, & T. Henning, ASP Conf. Ser., 387, 381 [Google Scholar]
van der Tak, F. F. S., van Dishoeck, E. F., Evans, II, N. J., et al. 2000, ApJ, 537, 283 [NASA ADS] [CrossRef] [Google Scholar]
van der Tak, F. F. S., & Menten, K. M. 2005, A&A, 437, 947 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Vehoff, S., Nürnberger, D. E. A., Hummel, C. A., et al. 2008, in Massive Star Formation: Observations Confront Theory, ed. H. Beuther, H. Linz, & T. Henning, ASP Conf. Ser., 387, 444 [Google Scholar]
Walmsley, M. 1995, in Rev. Mex. Astron. Astrof., Conf. Ser. 1, Circumstellar Disks, Outflows and Star Formation, Cozumel, Mexico, 137 [Google Scholar]
Warren, S. J., Hambly, N. C., Dye, S., et al. 2007, MNRAS, 375, 213 [NASA ADS] [CrossRef] [Google Scholar]
Whitney, B. A., & Hartmann, L. 1992, ApJ, 395, 529 [NASA ADS] [CrossRef] [Google Scholar]
Whitney, B. A., & Hartmann, L. 1993, ApJ, 402, 605 [NASA ADS] [CrossRef] [Google Scholar]
Whitney, B. A., Wood, K., Bjorkman, J. E., et al. 2003a, ApJ, 598, 1079 [NASA ADS] [CrossRef] [Google Scholar]
Whitney, B. A., Wood, K., Bjorkman, J. E., et al. 2003b, ApJ, 591, 1049 [NASA ADS] [CrossRef] [Google Scholar]
Whitney, B. A., Indebetouw, R., Bjorkman, J. E., et al. 2004, ApJ, 617, 1177 [NASA ADS] [CrossRef] [Google Scholar]

Footnotes

... W33A: Based on observations with the VLTI, proposal 381.C-0602.
... survey: http://www.leeds.ac.uk/RMS
...Robitaille et al. (2006): http://caravan.astro.wisc.edu/protostars/index.php

All Tables

Table 1: List of MIDI observations of W33A.

Table 2: Parameters for the preferred model discussed in Sect. 4.2.

Table 3: Models from the YSO SED model grid (Robitaille et al. 2006), obtained using the web-based fit procedure.

All Figures

$\begin{figure} \par\includegraphics[height=9cm,width=8cm]{13209fg1.eps} \vspace*{-3mm} \end{figure}$	Figure 1: UKIDSS K-band observations of W33A. Projected VLTI/MIDI baselines are indicated by white bars.
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[height=9cm,width=7cm,angle=90]{13209fg2.eps} \end{figure}$	Figure 2: $\rm ^{12}CO~(3-2)$ map taken with JCMT/HARP illustrating the molecular outflow of W33A. The PA and orientation of the blue and red-shifted lobes of the outflow are consistent with the K-band reflection nebula of Fig. 1. Coordinates are in J2000, and they correspond to the UKIDSS near-IR source and accord with the Capps et al. (1978) main IR source.
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[height=18cm,width=14.5cm,angle=90]{13209fg3.ps} \end{figure}$	Figure 3: MIDI observations of W33A on four different occasions (see Table 1). Top: flux spectra. Middle: correlated flux spectra (points with errorbars), and the detector noise level (dashed curve). Bottom: derived visibilities.
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[height=12cm,width=7cm]{13209fg4.eps} \end{figure}$	Figure 4: Left: UKIDSS H-band ( upper panel) and K-band ( middle panel) observations of the base of the W33A outflow scattering nebula. Images are logarithmically scaled, in units of MJy/sr and rotated by $35\hbox {$^\circ $ }$ so that the main lobe is pointing downwards. H-band contours are at 7.5%, 15%, 40% and 70% of peak flux; the K-band contours are at 1%, 2.5%, 10%, and 40% of peak flux. Secondary stellar sources have been removed from the images. The lower panel shows the H-K colour images (in magnitudes). Right: model images produced with the RT scattering code, see Sect. 4.1.2. The model images are convolved with 2-D Gaussian function with a FWHM of 1.6 $^{\prime \prime }$ . All images cover an area of $25\hbox {$^{\prime \prime }$ }$ on the side.
Open with DEXTER
In the text

$\begin{figure} \par\mbox{\includegraphics[height=9cm,width=7cm,angle=90]{13209fg... ...} \includegraphics[height=9cm,width=7cm,angle=90]{13209fg6.eps} }\end{figure}$	Figure 5: Left: SED of the preferred model overplot W33A flux densities. Small filled circles trace the ISO SWS spectrum, small open circles are the MIDI spectrum, the filled square is the Spitzer MIPS point, open square is the IRAS measurement at 60 $\mu$ m. In the (sub)mm range, the integrated fluxes at 350 $\mu$ m, 850 $\mu$ m, and 1.2 mm are from van der Tak et al. (2000), Di Francesco et al. (2008), and Faúndez et al. (2004), respectively. Right: zoom in onto the Nband spectrum. Shown are the MIDI flux spectrum with errorbars and the model fluxes appropriate for the MIDI slit width of 0.6 $^{\prime \prime }$ .
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[height=18cm,width=15cm,angle=90]{13209fg7.eps} \end{figure}$	Figure 6: Model (full line and filled circles) fits to the observed visibilities. Model parameters are listed in Table 2. Corresponding images are presented in Fig. 8. Note the discrepancy with the spherically symmetric model (dotted lines) used to fit configuration D in Paper I.
Open with DEXTER
In the text

$\begin{figure} \par\mbox{\includegraphics[height=9cm,width=8cm,angle=90]{13209fg... ...} \includegraphics[height=9cm,width=8cm,angle=90]{13209fg9.ps} } \end{figure}$	Figure 7: Logarithmically scaled slices through the spatial distribution of the density ( left) and temperature ( right) of the preferred model perpendicular to the outflow axis. The cavity has a very small, but non-zero density consisting of gas and dust (see text).
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[width=18cm,angle=0]{13209fg10.eps}\par\includegraphics[width=18cm,angle=0]{13209fg11.eps} \vspace*{4mm} \end{figure}$	Figure 8: Logarithmically scaled images ( top row) and corresponding Fourier transforms ( bottom row) of the preferred model (Table 2) at wavelengths 7.5 $\mu$ m to 12.5 $\mu$ m in steps of 1 $\mu$ m. Image field of view is $0.6\hbox {$^{\prime \prime }$ }$ (2280 AU) on a side, corresponding to the MIDI slit width. Pixel scale is 0.01 $^{\prime \prime }$ . The model images have been multiplied with a wavelength dependent 2-D Gaussian function, simulating the VLT-UT airy disk. Warm cavity walls dominate the emission at 10 $\mu$ m.
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[height=8.5cm,width=8cm,angle=90]{13209fg12.eps} \end{figure}$	Figure 9: Morphology at 350 $\mu$ m predicted by the preferred model. Contours are at 15%, 30%, 50%, 70% and 90% of the peak value. Model image has been convolved with a 10 $^{\prime \prime }$ beam. See van der Tak et al. (2000) for the observed 350 $\mu$ m image.
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[height=8.5cm,width=7cm,angle=90]{13209fg13.eps} \end{figure}$	Figure 10: Model predictions for baseline configuration A. Represented are the preferred model (full line), with the addition of a dust disk ( $M=0.25~\mbox{$M_{\odot}$ }$ , dashed line), and an optically thick disk interior to the dust sublimation radius (dotted line).
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[height=8cm,width=8cm]{13209fg14.ps} \end{figure}$	Figure 11: Logarithmically scaled K-band image for grid model 3004217 to be compared with the observations in Fig. 4. The contours are at 1%, 40%, 50%, and 65% of the peak flux. The image is convolved with a 2-D Gaussian function with an FWHM of 1.6 $^{\prime \prime }$ .
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[height=8.2cm,width=8.35cm,angle=90]{13209fg15.ps} \end{figure}$	Figure 12: Predicted SEDs of the discussed models obtained from the SED web fit procedure listed in Table 3.
Open with DEXTER
In the text

$\begin{figure} \par\includegraphics[height=8.2cm,width=8.35cm,angle=90]{13209fg16.ps} \end{figure}$	Figure 13: Predicted visibilities along the major (lower curve for each model) and minor axis (top curve for each model) at 8 $\mu$ m as a function of baseline for the discussed models obtained from the SED web fit procedure (Table 3). MIDI visibilities for the four baselines are also shown.
Open with DEXTER
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.