© ESO, 2012

1. Introduction

The formation and early evolution of massive stars are not well understood. Direct observations of these processes are hampered because massive stars form on a much shorter time-scale than lower mass stars. In their earliest infancy they are hidden behind huge amounts of optical extinction. The nearest massive young stars are also more distant than their lower-mass counterparts. The much higher luminosity and intense UV radiation field of young massive stars strongly affects their environment at very early times. The geometry and composition of the gas and dust around the star may differ from that seen in lower mass stars, as may the time-scale on which this circumstellar material dissipates. Massive stars are also known to have prominent stellar winds (e.g. Vink et al. 2001), which complicates the picture even more. The lack of a long pre-main-sequence phase, during which the circumstellar material slowly dissipates and a planetary system may form, marks an important qualitative difference between massive and lower mass young stars. This being said, the upper mass for stars to form planetary systems is not known.

Herbig Be (HBe) stars are a subclass of Herbig Ae/Be stars, intermediate-mass young stars first described as a group by Herbig (1960). They represent the more massive stars in the group, with zero-age main-sequence (ZAMS) masses up to roughly 20 M_⊙. Many observational studies have described the properties of Herbig Ae/Be stars (see e.g. Natta et al. 2000). A picture has emerged in which the lower mass Herbig Ae (HAe) pre-main-sequence stars have disks that can be considered as scaled-up versions of disks found around the solar mass T Tau pre-main-sequence stars. These disks are formed during the collapse and proto-stellar phase as a result of angular momentum conservation, and are the site of planet formation. The dust in these gas-rich disks has undergone substantial settling and grain growth, which can be inferred from the large mass in cold, millimeter sized grains detected in many systems (Natta et al. 2007). The warm disk surface layers in many HAe stars also show evidence for grain growth and crystallization in the mid-infrared (IR) (e.g. van Boekel et al. 2003; van Boekel et al. 2005; Schegerer et al. 2006; Kessler-Silacci et al. 2006).

There are observational indications that this qualitative picture no longer applies for higher mass stars. Searches for a reservoir of large, cold grains in more massive HBe stars have shown that the mass in cold grains is substantially lower than that in lower mass stars (Fuente et al. 2003). In a more recent study, Alonso-Albi et al. (2009) show that the disk mass in some HBe stars is only on the order of 10 per cent of the total envelope mass; they propose that the strong UV radiation field of the central star causes rapid photo-evaporation of the outer disk, leaving a relatively small inner disk. Observations of the inner disk structure also suggest differences between high- and low-mass stars. Monnier et al. (2005) showed on the basis of near-IR interferometry that the inner radius of disks in lower mass stars is set by the dust sublimation radius and that the inner cavity is optically thin. For the high-luminosity B0-B3 stars, however, smaller inner dust radii are found, possibly implying optically thick gas in the inner regions. A detailed study of the geometry of the disk surrounding the B1.5 star MWC 297 showed that the disk is very compact and may extend to (close to) the stellar surface, with little or no evidence for an inner cavity (Acke et al. 2008).

Theoretical studies of disk evolution in massive young stars predict that the outer disk evaporates on short time-scales while the inner disk can survive longer (e.g. Hollenbach et al. 1994; Gorti & Hollenbach 2009). Disk evaporation time-scales are found not to vary with stellar mass for masses lower than  ~3 M_⊙ and disk lifetimes are on the order of 10⁶ yr. However, for the more massive stars the intense UV radiation field evaporates the disks on time-scales on the order of 10⁵ yr. Translating this mass into a spectral type at the ZAMS, substantial changes in disk lifetime and disk geometry are expected in the mid-B spectral range.

What then are the observed properties of circumstellar matter in HBe stars? Is there evidence for grain growth and processing in HBe stars, similar to that seen in HAe stars? Is the spectral energy distribution of HBe stars at IR wavelengths different from that of lower mass stars, and how does this translate into (inner) disk geometry? To answer these questions and put HBe stars into perspective, observations at mid-IR wavelengths are needed at high angular resolution. So far, many studies were based on IR data obtained with large beams, which caused contamination with ambient molecular cloud emission, which is prominent for massive stars because of their strong ultraviolet (UV) radiation field. In the present study, we employ the subarcsecond spatial resolution in imaging and spectroscopy that VISIR at the VLT offers to study the immediate environment of HBe stars on spatial scales comparable to the disks they may have.

This paper is organized as follows: in Sect. 2 we discuss the selection of our sample of HBe stars and their basic properties. Sections 3 and 4 describe the observations and data reduction procedur es, and we analyze the images and spectra in Sects. 5 and 6. In Sect. 7 we discuss the results and compare the HBe stars to their lower mass HAe counterparts. Section 8 contains the main conclusions.

2. Sample

2.1. Selection

We have taken targets from the Herbig Ae/Be catalogs of Thé et al. (1994) and Vieira et al. (2003) and supplemented this list with Herbig Ae/Be stars discovered by Perez et al. (1989), de Winter et al. (1997), and van den Ancker et al. (1998). We selected the targets based on several criteria. They have a spectral type B or O, with an uncertainty in spectral type of less than three subclasses. They are known to have an IR excess from near-IR photometry. They are associated with a star-formation region. They are bright enough to observe with VISIR in the mid-IR according to their IRAS 12 μm flux or according to the flux expected from an extrapolation of the photospheric flux at optical wavelengths. They are observable with the VLT in Chile and consequently mostly located in the southern hemisphere. While our approach is optimized to select young massive stars, it is possible that our sample is contaminated with stars in other evolutionary stages. We address this problem in Sect. 7.4. Furthermore, recent improvements of the spectral types have shifted two of the targets to the HAe class (LkHα 208 and LkHα 339). These were kept in the sample for the sake of comparison and completeness.

2.2. SEDs

To establish the pre-main-sequence nature of our targets, we pinpointed their position in the Hertzsprung-Russell (HR) diagram. To do this we gathered all available photometry and constructed spectral energy distributions (SEDs) for all targets that were successfully observed with the VISIR spectrograph (see Fig. 1). We gathered spectral types, effective temperatures, and distances from the literature (see Table 1). If the effective temperature was not stated by the authors, we used de Jager & Nieuwenhuijzen (1987) to find a value corresponding to the spectral type. With the effective temperature and distance we were able to assign a theoretical stellar SED using the models from Kurucz (1991). By comparing this with the observed photometry, we obtained the optical extinction, A_V, and the luminosity of the star L _∗ . In this comparison, we used the Steenman & The (1991) interstellar extinction law characterized by a total to selective extinction R_V = 3.1. For three targets no distance estimates are specified in the literature. For these we used the spectroscopic parallax to determine the stellar luminosity. For Hen 3-847 we used log    g = 3.84 (see van den Ancker et al. 1998 and references therein), for PDS 477 and PDS 543 we derived log g = 4.0.

In Fig. 1 we plot the constructed SEDs. Square symbols denote the de-reddened photometry. Curves are the Kurucz models. The plots show that the stellar photosphere of a number of sources is not properly constrained, which increases the uncertainty on the luminosity somewhat, especially when the uncertainty on the spectral type is large as well (e.g. PDS 204). For the other targets the goodness of fit of the Kurucz model to the de-reddened photometry give confidence in the reliability of the derived luminosities.

2.3. Stellar parameters

In Table 1 we summarize the stellar parameters for our sample targets. We only derived optical extinctions and luminosities for the targets that were successfully observed with the VISIR spectrograph. Masses for the stars below 7 M_⊙ were found by interpolating in the HR diagram (see Fig. 2) between PMS evolution tracks of Siess et al. (2000). For the stars above 7 M_⊙ we interpolated between MS evolution tracks of Meynet et al. (1994). We adopted 20% errors on the masses. The high-mass Herbig stars (M ≳ 8 M_⊙) are known to become visible only when they reach the ZAMS (Palla & Stahler 1993). Therefore we would expect a lot of young stars close to the ZAMS. However, the optical selection of the targets has a bias toward low A_V and consequently it favors stars that have dispersed their natal clouds. These likely are young main-sequence stars.

To calculate the uncertainties on the luminosities and radii (see Table 1), we did not take the uncertainties on the effective temperatures into account, but accounted for the uncertainties in the distances only. Most of the distances are based on associations with star-formation regions and clusters and their uncertainties assume the associations are correct. When no error on the distance was stated in the literature, we assumed an error of 20%. Another source of uncertainty that was not taken into account is in the choice of R_V. Some studies (e.g. Hernández et al. 2004) argue that R_V = 5.0 is more appropriate for star-forming regions. Such a higher R_V value leads to considerably higher luminosities. In Table 2 we provide values available in the literature and compare these to values that were derived with R_V = 5.0 and R_V = 3.1. The variation in the luminosities in Table 2 illustrates the remaining uncertainties.

3. Observations

3.1. Imaging

From October 2006 to May 2007 N-band imaging was obtained with the VLT Imager and Spectrometer for mid Infrared (VISIR, see Lagage et al. 2004). The targets were observed for typically  ~3 min. in the small-field mode (pixel field of view = 0.075′′) and in each of the three filters: SiC, PAH2, and PAH2_2. These focus on the 10 μm silicate feature, the 11.2 μm PAH feature, and the continuum redward of the PAH feature respectively (see Table 3). For the science targets we employed parallel chopping and nodding to remove the atmospheric background. A maximum chopper throw of 30″ was chosen to avoid any confusion with surrounding nebulosity. Originally we chose chopping rotations, which were offset from the standard north-south orientation for the same reason of avoiding confusion, because this was expected from optical and IRAS imaging of the targets. In the next observing period we did use the standard north-south orientation, because almost all targets turned out to be nearly point-like in the mid-IR. Observations were carried out in ESO’s service mode at airmass  <1.6 and optical seeing  <1.2″, corresponding to a seeing in the N-band of  ≲0.5″. Photometric standard stars were observed in all filters less than three hours before or after the science targets. Table 4 contains a detailed log of the observations. It lists the target names, the imaging filters, the image orientation, the chopper throw, the date and time, the airmass and optical seeing, the integration time and the target category (calibration or science). Four sources were not detected with our VISIR imaging. They have been indicated with a † symbol in Table 1.

3.2. Spectroscopy

VISIR long-slit low-resolution (R ≈ 300) spectra were also obtained in ESO’s service mode from January 2007 to February 2008. The complete N-band was mapped with the sequence of four separate observations in the 8.5, 9.8, 11.4, and 12.2 μm settings. Each observation had a typical integration time of  ~5 min, a slit-width of 0.4″, a north-south orientation and a chopping throw of 8″. Based on our imaging results, some observations were performed using different choices of slit orientation and chop throw. In addition to the four sources that were not detected in imaging at all (indicated with † in Table 1), six sources turned out to be too faint for spectroscopic follow up (#) and one was not covered owing to lack of observation time (Hen 3-1444). Several sources displayed apparent extended emission; for these the slit was rotated to select the major axis of the target. For AS 310 and PDS 204 wider chopper throws were chosen because of their larger mid-IR sizes. Bright spectroscopic standard stars were observed directly before and after the science observation for the purpose of telluric correction, but also to determine the spectroscopic point spread function (PSF). Table 5 contains a detailed log of the observations. The table lists the target names, the wavelength settings, the orientation of the slit, the chopper throw, the date and time, the airmass and optical seeing, the integration time and the target category (calibration or science).

4. Data reduction

4.1. Image reduction

The raw VISIR images and spectra were reduced with a dedicated pipeline that we developed (see Pantin 2010). The basics are explained in our study of HAe stars (see Verhoeff 2009). The resulting images are displayed in Fig. 3 for the PAH2 filter. For some targets with detached emission patches we display both the PAH2 filter (Fig. 4) and the PAH2_2 filter (Fig. 5). The images are displayed on a logarithmic scale and the color contrast ranges from 1σ of the background (white) to 0.9 of the peak flux (black). When the orientation deviates from the standard north up and east left, a compass in the bottom right of the images is supplied. For the observations that were followed-up with spectroscopy, the slit orientation is indicated.

4.2. Telluric correction

The VISIR spectra of our targets were calibrated immediately before and after the science observation, which in most cases allowed us to do the telluric correction with an airmass interpolation of the extinction. The principle is explained in detail in van Boekel et al. (2004) and it was applied for VISIR spectra of HAe stars in Verhoeff (2009). In some cases there was only one suitable calibration measurement. For those cases we followed the formalism described in Verhoeff (2009). The flux calibration was achieved with the photometry from the SiC filter (see Sect. 5.1). The resulting spectra of our HBe stars are plotted in Fig. 6. For comparison we overplotted ISO spectra from the ISO archive and Spitzer spectra taken from Acke et al. (2010) and Van den Ancker (unpublished; Spitzer PID 50180).

4.3. FWHM extraction

Because VISIR employs a long slit, spatial information of the targets is maintained along the slit. To characterize this information, we fitted Gaussians to the spatial direction of the 2D spectra to extract the spatial FWHM in all wavelength bins. The resulting plots can be found in Fig. 7. The combined measurement error and fitting error was estimated by evaluating the standard deviation of the residuals after subtraction of a spline fit. The point spread function (PSF) of the combined system of VISIR and the VLT has been proven to be quite unstable (see e.g. Pantin 2010). To determine the PSF during science measurements we constructed a “statistical calibrator” from a sample of  ~ 130 calibration measurements. A detailed description of our analysis can be found in Verhoeff (2009), below we provide a brief summary.

We identified three components that contribute to the PSF: the diffraction limit, the atmospheric seeing, and a stochastic telescopic broadening. To establish a reliable PSF and its uncertainty we determined the dependence on the atmospheric optical seeing and the magnitude of the random distribution. This was achieved by plotting all calibration measurements in a diagram with their IR FWHM size against their optical seeing, after a correction for their airmass. The fit of these data defines our statistical calibrator and the scatter around it gives the uncertainty. In practice a PSF and its uncertainty can now be found for a random measurement by supplying an airmass-corrected optical seeing value, reading the corresponding IR FWHM size from the fit, and transforming that value back to the airmass of the measurement. In our analysis we will also use the calibrator nearest in time as a reference.

5. Analysis of the images

5.1. Photometry

Most of our sources are surrounded by nebulosity, which contributes significantly to the IRAS flux. Imaging with VISIR allows one to zoom in to an area of  ~1″ around the star. The flux levels that are obtained in this way give new information about the star and disk.

From inspection of the morphology of the stars in the obtained imaging we conclude that a standard circular photometry aperture is suitable for all sources. More specific, we measured the source counts within a circle of radius r = 13 pixels (~1.0″) and the background within an annulus within radii r = 18 to r = 23 pixels (~1.4″ to  ~1.7″). For a typical PSF the source aperture encompasses the first and second airy rings, the third airy ring is ignored and the background aperture runs from the fourth to the fifth airy ring.

The apertures were applied to obtain the background-subtracted flux levels in counts per second of all science and calibration measurements. To obtain the tellurically corrected flux levels of our science targets we used (1)where I is the intrinsic flux, S the measured background-subtracted flux, m the airmass of the measurement and R and A the responsivity of the system and the atmospheric extinction per unit airmass. To find the responsivity of the system, we looked at the ensemble of all calibrators of all nights and plotted the natural logarithm of their observed counts over their model fluxes against the corresponding airmass values. In this way a least-squares linear fit to the data points allows us to obtain the responsivity: (2)The intersection with the y-axis represents the responsivity and the slope is determined by the extinction per unit airmass. The latter is not a constant, it depends on the humidity and temperature of the atmosphere. Nevertheless, when we fit the ensemble of calibrators, the uncertainty of the responsivity becomes small compared to the other errors.

Now that we know the value of the responsivity, we can calculate a value of the extinction per unit airmass for every single calibration measurement through an inversion of Eq. (1). To obtain the extinction per unit airmass at the time of the science measurement, we used the calibrator nearest in time. To estimate the accompanying error, we used the variation in extinction as derived from all calibrators during the same night.

This procedure was followed for all three filters (see Table 3) and all science targets; the results are summarized in Table 6. The listed errors are obtained by differential propagation of the measurement errors and the errors in the responsivity and extinction. For the undetected sources Hilt 641, AS 267, NGC 6530-45, and LkHa 112, we assert an upper limit in the SiC filter of 0.1 Jy.

5.2. Spatial emission profile

We analyzed the spatial emission profile of the images by extracting horizontal cuts through the source maxima. In the vertical direction the profiles are corrupted by detector ghosts, caused by the fairly high brightness of some sources (see Pantin 2010). Given the high amplitude of the chopping throw (30″) and the observational mode used (parallel nodding mode), the limited accuracy of superimposition of the two central beams results in a potential degradation of the spatial resolution. The chopped-only images were therefore analyzed separately. The two profiles corresponding to the two nod positions are shown in Fig. 8. These profiles were averaged over three pixel rows centered on the peak intensity of the sources. Overplotted in green are the calibration measurements that were observed within two hours before or after the science target and had minimum airmass differences. The average of these two measurements is a reliable representation of the PSF. The difference between the profiles of the different nodding positions provides an indication of the reliability of the measurements. For the sources with low signal-to-noise ratio (S/N), such as PDS 241, we were forced to use a total stack of the data, which thus results in a single profile curve. Hen 3-1418 and MCW 1264 were too faint for a meaningful analysis.

To obtain a size measurement of our targets, we fitted the extracted emission profiles with a Gaussian. Where we had multiple profiles, we averaged the resulting FWHM values. The nearest calibration image was then taken as the PSF and was quadratically subtracted from the science measurement. The error in the PSF (σ_PSF) was estimated by taking the standard deviation of all calibration measurements in the same filter of the program. The measurement errors were estimated to be 0.3*σ_PSF. For high seeing conditions we increased this error continuously with 1.0 σ_PSF/arcsec over an optical seeing of 0.7″. The final errors were obtained by numerical propagation of the measurement and PSF errors. This rather conservative error budgeting was chosen because the imaging was meant to provide a quick grasp on morphology and photometry and the PSF was therefore relatively underdetermined. The analysis was performed for every target and every filter. The results are displayed in Table 7.

The last two columns of Table 7 give a description of the spatial emission profile at low intensities. After inspection of the images we selected targets that contain detached emission components, indicated with the label “cloud”. We display these separately in Fig. 4. We assume that these components are caused by illuminated natal cloud remnants. These sources also show associated optical nebulosity in the Digitized Sky Survey 0.65 μm imaging (see “Neb” column in Table 6), which is consistent with this hypothesis. In the continuum filter (PAH2_2) the features are nearly absent (see Fig. 5), indicating that the nature of the emission is most likely fluorescence from PAH molecules. After comparing the emission profiles in Fig. 8 with their calibrators, we were able to discriminate between targets that are “point”-like and those that have spatially extended emission. We assigned the label “extended” when the profiles of both nod positions were at any intensity level more extended than the most extended profile of the calibrator. Note, however, that these qualifications are slightly speculative because of the underdetermined PSF. We stressed this by adding a question mark where there were seeing variations between science and calibration measurements. Finally, we indicate the presence of the diagonally striped pattern in the background of the PAH images (Fig. 3) with a tick-mark. These features appear when observing a fairly extended (15″) background, because such a background causes a telescope residual, which is rotated in two consecutive nod positions (see Pantin 2010). This means that around these targets we can infer the presence of a homogeneous plateau emission on scales larger than 15″.

6. Analysis of the spectra

6.1. Extinction correction

In order to determine the composition of the circumstellar dust responsible for the 10 μm feature, we first need to correct the spectra for interstellar extinction. To convert the optical extinction A_V to an extinction in the N-band, we need to determine the shape of the extinction curve, which is quite different for the diffuse interstellar medium (ISM) and for dense molecular clouds. Chiar et al. (2007) set the maximum A_V for reliable use of the diffuse ISM extinction curve at 12 mag. None of the A_V of the sample stars reach that maximum. We therefore adopted the diffuse ISM extinction law given by Whittet (2003): A_V = 6.2·E(J − K_s), where E(J − K_s) is the near-IR color excess. A linear fit given by Chiar et al. (2007) of the relation between the optical depth of the 9.7 μm silicate absorption feature, τ_9.7, and E(J − K_s) provides the conversion to the N-band: τ_9.7 = 0.35·E(J − K_s) = 0.056·A_V. We multiplied the spectra with the function e^τ_λ, with τ_λ given by the ISM silicate feature scaled to the calculated 9.7 μm optical depth. The opacity of the ISM dust was taken from Min et al. (2007).

6.2. Compositional fits

In order to allow a direct comparison, the extinction corrected spectra were fitted with the same routine and the same spectral templates as used by van Boekel et al. (2005) for HAe stars. The reference for amorphous silica was incorrectly stated in that paper (priv. comm.) and should be Henning & Mutschke (1997). The routine fits the most commonly found dust species in circumstellar material that show spectral structure in the 10 μm region. The only difference in our version of the routine is that it fits the three different PAH features at 8.6, 11.2, and 12.7 μm independently. The extinction corrected spectra and their best fits are displayed in Fig. 9. The resulting mass fractions are given in Table 8. The small and large grains have volume equivalent radii of 0.1 and 1.5 μm respectively. The targets LkHα 215 and Hen 3-1110 were omitted from the analysis because of their low S/N.

There are a few remarkable characteristics appearing in the compositional analysis. The 10 μm features contain little PAH emission, there is a fair amount of crystalline silicates and there is a notable presence of amorphous silica at 9 μm. In Sect. 7 we will put these characteristics into context when we compare the HBe sample to the HAe sample of van Boekel et al. (2005). R Mon will be excluded from the discussion of the silicates because of the large error bars, probably because of a substantial underestimation of the IR extinction. To check the influence of our extinction correction on the other sources, we also obtained compositional fits to the uncorrected spectra. We found that an underestimate of the extinction correction could lead to a slight overestimate of the crystallinity. PDS 27 turns out to be unusual, because its 10 μm feature does not show spectral structure. This most likely indicates the absence of grains smaller than a few micron.

We note strong [Ne II] emission at 12.81 μm in the Spitzer and VISIR spectra of MWC 878 (Fig. 6). This line in the VISIR spectrum shows that it must originate in the immediate vicinity of the star. [Ne II] 12.81 μm emission is also detected in a significant fraction (~20%) of T Tauri stars, where it has been interpreted as caused by X-ray excitation on the surface of a flared circumstellar disk (e.g. Lahuis et al. 2007; Pascucci et al. 2007). We note that MWC 878 is one of the earliest-type stars in our sample (O9/B0), suggesting that X-ray excitation may also be important for the line emission in these stars.

6.3. Spatial extent

The observed sizes as a function of wavelength as displayed in Fig. 7 are quite similar to the results found with the same method for HAe stars (see Verhoeff 2009). The overall slopes are quite flat and no signature of the silicate feature can be identified. For some targets the size measured in the continuum is larger than the statistical PSF. An increase in size can be found in some targets at the wavelengths of the PAH emission features (PDS 204, MWC 147, Hen 3-331, Hen 3-847). Notable is the upturn to the left of 8.3 μm, which is associated with the 7.7–7.9 μm PAH complex. There are also slight increases in extent at 8.6 and 11.2 μm, which are associated with the C-H in plane and C-H out of plane bending modes of the PAH molecule, respectively.

To obtain the intrinsic sizes of the sources, we quadratically subtracted of the statistical PSF for every setting and in every wavelength bin. The error for the absolute scale of every setting was then derived by numerical propagation of the measurement errors and the error on the statistical PSF. Subsequently, we defined the continuum size of a particular target as the mean of the median values of all settings. In this way we used all settings as semi-independent measurements of the N-band size and filtered out the size increases caused by PAH emission.

For the sources that contained a significant increase in FWHM at the wavelength of the 8.6 and 11.2 μm PAH features, we retrieved the intrinsic scale of the PAH emission by decomposing the spatial emission profile at the central wavelength of the features. We first determined the underlying continuum profile from the 2D spectrum. This was achieved by extracting the spatial emission profiles just outside the wavelength range of the PAH features. We averaged these profiles and then scaled the flux level to the continuum flux level as determined with the VISIR flux spectra. We then subtracted this continuum profile from the spatial profile at the central wavelength of the PAH feature. We fitted a Gaussian to the remaining “pure” PAH profiles to obtain the FWHM size. Finally, we quadratically subtracted the statistical PSF. In Table 9 we report the results of all analyses. For sources with a significant PAH flux feature (Table 8), but no significant change in FWHM, we state the continuum size for that feature. For a more detailed description of the applied methods we refer the reader to Verhoeff (2009).

6.4. Differential extractions

For a few objects it was possible to extract more than one spectrum. This provides information about the spatial distribution of different dust species and PAH molecules. For PDS 204, we extracted a five pixel wide trace through the central component and a trace through the attached emission patch 1.7′′ south-southeast of the target. Figure 10 shows them overplotted. The emission patch  ~1700 AU away from the star generates the same amount of flux in the continuum-subtracted 11.2 μm PAH feature.

The large spatial scale (diameter  ~0.05 pc) of the continuum emission of AS 310 is very peculiar (see Fig. 5). The nature of this emission is not clear. In Fig. 11 we display five-pixel-wide traces through the primary source, the background 1′′ north of the primary, and the secondary component 4.6′′ north of AS 310. The primary has a distinct silicate feature, whose shape suggests significant grain processing (crystallization and/or grain growth), which is commonly observed in disks. However, the large continuum size in the imaging and in the spectrum (Table 9) seems to suggest we are looking at a huge cloud of dust. The secondary component does not seem to display a silicate feature. The deficit of PAH emission in either of the components or the background is also puzzling (compare Fig. 10). Note that in the optical, this source also has quite prominent nebulosity. How the distribution of material can generate the observed mid-IR features is not clear.

7. Discussion

7.1. IR excess

The HAe stars are relatively isolated, but for the HBe stars the IRAS, MSX, and 2MASS photometry contains contributions from the surrounding molecular cloud ( ≳ 1″). Therefore a classification of the HBe objects based on wide-field photometry is not always meaningful. The VISIR observations help to isolate the mid-IR emission from the immediate vicinity of the star and provides more clarity on the nature of the objects. Of our original sample (see Table 1) four sources were not even detected in the VISIR imaging (indicated with a †). These sources (Hilt 641, AS 267, NGC 6530-45, and LkHα 112) do not have a compact IR excess. The emission detected by IRAS is either caused by an unrelated nearby source, or by diffuse interstellar matter heated by the star. Another six sources (#) had a 12 μm excess as seen by VISIR that was too faint for spectroscopic follow-up (PDS 241, MCW 1264, Hen 3-1418, PDS 152, MWC 593, and LkHα 115). These sources are unlikely to possess a significant amount of circumstellar dust. This is also valid for LkHα 215 and Hen 3-1110, albeit to a lesser extent. In Sect. 7.4 we will discuss the nature of the sources in more detail. For now we focus on the sources that do have a compact IR excess and for which we were therefore able to obtain VISIR spectra (see Fig. 6).

The IR excess emission of all of these sources starts in the near-IR and extends to the far-IR. Figure 1 displays the wide variety of shapes that is observed. There are declining SEDs (e.g. Hen 3-331), flat SEDs (e.g. PDS 27), and rising SEDs (e.g. LkHα 339). Some SEDs show a broad bump that peaks around 10 μm (e.g. Hen 3-847). We must keep in mind, though, that most of these targets are still shrouded in their natal clouds, which significantly contributes to the mid- and far-IR fluxes. With the higher spatial resolution of VISIR we are in most cases capable of isolating the emission of the circumstellar disk. The contrast is illustrated in Fig. 1. The red triangles represent our VISIR photometric results at 12 μm (see Sect. 5.1). Sources with an obvious substantial environmental contribution include LkHα 339, LkHα 215, LkHα 241, and AS 310. Even with the VISIR photometry it is difficult to classify the disk in terms of flaring or flat disk geometries (Meeus et al. 2001). This requires additional high spatial resolution photometry at longer wavelengths. In the last column of Table 1 we indicated our classification into group I, II, G, and U. The “G” stands for gas-disk (see Sect. 7.4) and the “U” for uncertain. This latter classification is needed because it is not possible to discriminate between a flaring disk and a flat disk plus cloud contribution. Note that as a consequence, we cannot identify group I sources in our sample.

The collection of shapes of the IR excess of the HBe stars resembles those that we see for the lower mass HAe stars (see Meeus et al. 2001). A noticeable difference, however, is that the average IR excess is weaker for the HBes, and that the sample contains a higher diversity. The median near-IR excess in the HAe sample of Verhoeff (2009) is a factor 8 higher than the median value of our HBe sample. The spread is also notably different; it is twice as extended in the HBe sample. The far-IR excess in the IRAS 60 and 100 μm bands, on the other hand is higher on average for the HBe stars, and their spread is much wider. This can be explained with cloud contributions, as mentioned before. For the comparison of the IR excess at 12 μm we are able to consider the emission from the disk only. In Fig. 12 we plotted the 12 μm luminosity, L₁₂, against the stellar luminosity for the HAe sample of van Boekel et al. (2005) and for our HBe sample. The L₁₂ was calculated from the distance, d, and the IRAS 12 μm flux for HAe stars and the VISIR SiC flux for the HBe stars (F₁₂) as follows: L₁₂ = 4π·d²·F₁₂·Δν, where we set Δν equal to the bandwidth in Hz of the SiC filter (see Table 3) for both samples. The plot shows that the average L₁₂ increases with stellar luminosity, but the typical reprocessing power, L₁₂/L_∗, decreases with stellar luminosity. The median reprocessing at 12 μm is a factor 5 smaller for the HBe stars than for the HAe stars. This is also true for the total IR reprocessing power. From the literature photometry we obtain a median L_IR/L_∗ of 45% for the HAes and only 15% for the HBes. A correction for cloud contribution would lower the latter even more. This result was already observed for one particular HBe star: Acke et al. (2008) found that the total IR excess of MWC 297 is only L_IR/L_∗ ~  5%.

We interpret the differences in the IR SEDs of the HBe and HAe stars as geometrical differences of their disks. The lower total IR reprocessing power (see Fig. 1) is directly related to a narrower covering angle of the disk. The lower near-IR reprocessing corresponds to a lesser contribution from the hot inner part of the disk. This suggests that the inner disks of many HBe stars are flatter. This has also been demonstrated with IR interferometers for MWC 147 (Kraus et al. 2008) and MWC 297 (Acke et al. 2008). Corrected for cloud contribution, we would expect to find a lower far-IR reprocessing efficiency as well, corresponding to a lesser contribution from the outer disk. This would imply that the outer disk is also flatter and/or that the magnitude of the outer radius has been reduced. The latter is consistent with predictions from photoevaporation theory (Gorti & Hollenbach 2009), which suggest that the higher FUV flux from the central star decreases the dissipation time-scale of the outer disk considerably. The disks around HBe stars are apparently flatter and are probably truncated on the outside.

7.2. PAHs

The comparison of VISIR and Spitzer/ISO spectra illustrates nicely how observations taken with wide apertures can report a much higher PAH flux (see e.g. MWC 147 in Fig. 6). This raises the question whether the PAH emission comes from a circumstellar disk, from a surrounding nebula, or both. Combining the results of the analysis of the spectra and images, we can locate the origin of the PAH emission. We indicate the results in the last two columns of Table 9. Where the data were not conclusive, we state added a question mark.

From the contrast with wide aperture spectra in Fig. 6 we know that PDS 204, MWC 137, MWC 147, Hen 3-847, V921 Sco, and MWC 297 display PAH from surrounding clouds, while LkHα 339, Hen 3-331, and MWC 878 do not. A comparison of the VISIR and IRAS photometry in Table 6 excludes MWC 789, LkHα 208, PDS 27, and MWC 878 from having any cloud emission in the mid-IR. A comparison of the VISIR imaging in the PAH2 filter and in the PAH2_2 filter confirms and adds a couple of sources to the list of cloud PAH emitters. For R Mon, PDS 27, V921 Sco, PDS 477, and LkHα 260 the imaging shows a weak cloud PAH component. PDS 204, MWC 147, and AS 310 display strong cloud PAH components (compare Figs. 4 and 5).

The VISIR spectra were extracted with a weight-map that emphasizes the spatially centered emission. They are therefore a relatively good tracer for the emission from disk-scales. We recall, though, that at 2 kpc the 0.4″ slit still encompasses 800 AU. By taking spectra at different locations we compared the PAH contribution from the central component and the surrounding cloud (see Figs. 10–11). The PAH components from the centers of these sources are not stronger than those of their surrounding media. Furthermore, the deconvolved sizes of the PAH emission have quite high values for these sources. This could mean that the PAH emission at the center of these sources is foreground emission from an envelope and not from a circumstellar disk. Of course, observations of even higher spatial resolution are required to clarify this matter. For now we will assume all the PAH emission in the VISIR spectra has a disk origin. The individual PAH bands were not all significantly detected. However, the total flux in the 8.6, 11.2, and 12.7 μm PAH bands together turned out to be significant for all sources in Table 8, with the exception of LkHα 260.

In Fig. 13 we plot the logarithm of F_PAH/F_∗, the total flux in the 8.6, 11.2, and 12.7 μm PAH bands divided by the stellar flux as a function of the stellar luminosity. The plot includes the HAe sample of van Boekel et al. (2005; open blue squares). It shows a trend of a decreasing PAH reprocessing power with increasing luminosity of the star. However, when we plot the F_PAH relative to F₁₂, we see that the ratio of the PAH and IR reprocessing powers is roughly constant (see Fig. 14). The lower PAH contribution from the disks of HBe stars could entirely be ascribed to smaller outer radii and/or more narrow covering angles of the HBe disks that we argued for in Sect. 7.1. The outer disk is known to deliver the dominant contribution to the PAH emission of the integrated flux spectrum. Observations of HAe stars also show that flatter disks have weaker PAH reprocessing powers (Meeus et al. 2001; Acke et al. 2004; Habart et al. 2004; Acke et al. 2010). We conclude that there is no clear evidence for a decreasing PAH abundance with increasing stellar luminosity.

7.3. Silicates

We have compared the silicate composition of our HBe sources as displayed in Table 8 (R Mon excluded) with the HAe sample of van Boekel et al. (2005). From Figs. 6 and 9 it is already clear that the overall silicate composition of the HBe stars is very similar to that of the HAe stars. Here we address the subtle differences that we did find.

Figure 15 shows that for the HAe stars the crystallinity is typically low (~5%). However, those sources for which the large grain component is particularly abundant (>80%) can have much higher crystallinity fractions. The HBe stars do not show this behavior. The average crystallinity fraction is slightly lower in the large-grain domain than in the small-grain domain. On the whole, the average crystallinity fraction is similar (8% versus 12% for the HAes). Within the margin of the errors this conclusion is independent of the extinction correction.

We find that for a typical HBe star the mass fraction of enstatite, amorphous silica, and forsterite is similar to that of a typical HAe star. At low forsterite abundances (<4%) we do see a difference, however, (see Fig. 16). The amorphous silica abundance is seven times higher in that range for the HBe stars. More generally it seems that the positive correlation between silica and forsterite for HAe stars, found by van Boekel et al. (2005), is no longer present for the HBe stars. The presence of this correlation was predicted by experiments that show that amorphous silica can be formed as a by-product in the annealing process of enstatite to forsterite (see e.g. Rietmeijer et al. 1986; Fabian et al. 2000). In this scenario, enstatite (MgSiO₃), silica (SiO₂), and forsterite (Mg₂SiO₄) are connected by the chemical equilibrium reaction (3)Sargent et al. (2009) suggested that in protoplanetary disks this reaction takes place via incongruent melting. In this context it is very interesting that we do not find such a correlation between silicate and forsterite for the HBe stars (see Fig. 16). Maybe the silica in HBe stars is predominantly formed as a primary condensation product of MgO-SiO₂ vapor in non-equilibrium conditions (Rietmeijer et al. 1999). This suggests that the high silica fractions in the HBe disks with low forsterite abundances reflect a history of dust vaporization followed by condensation.

7.4. Relation to B[e] and classical Be stars

Several classes of hot B-type stars exist that show permitted and/or forbidden emission lines in their optical spectra, and an IR excess caused by circumstellar gas and/or dust. Apart from the HBe stars discussed in this paper, the B[e] supergiants and the classical Be stars can show observational characteristics that are similar to those of HBe stars. We refer to Lamers et al. (1998) for a more complete discussion of this topic. B[e] supergiants are evolved, probably rapidly rotating high-mass stars with highly stratified outflows, a fast ionized wind at the poles and a slow, dusty outflow at the equator (Zickgraf et al. 1985). Their IR SEDs show a wide range of dust temperatures, similar to what is expected from disks surrounding young pre-main-sequence stars. The luminosities of B[e] supergiants in the LMC have been well established, and are typically in the range 10⁵ to 10⁶ L_⊙, with a few cases at lower luminosities (Zickgraf et al. 1986). The luminosities estimated for our sample of HBe stars are significantly lower, and so we conclude that it is unlikely that our sample is strongly contaminated by B[e] supergiants.

Classical Be stars are rapidly rotating main-sequence stars or giants with circumstellar gas located in a dense equatorial disk (Porter & Rivinius 2003). No forbidden emission lines are seen in the optical, and the IR excess is caused by free-free and bound-free emission from the gas in the disk. The mid-IR spectra of these stars typically show a power-law continuum with hydrogen recombination line emission superposed (Hony et al. 2000). There is not much evidence for dust in disks surrounding classical Be stars. The density in the disk is probably insufficient to form molecules and dust. The IR excess of classical Be stars at 12 μm is modest, typically up to 2.5 mag (Waters 1986), and much weaker than in disks surrounding HBe stars. In our sample we found some stars with a strong 12 μm excess based on IRAS photometry, but a much lower flux based on high spatial resolution VISIR mid-IR imaging, corresponding to a much weaker IR excess (see Table 6). These are LkHα 215, PDS 241, Hen 3-1110, MCW 1264, Hen 3-1418, PDS 152, MWC 593, and LkHα 115. Considering the SEDs, these stars might be classical Be stars rather than HBe stars (see e.g. LkHα 215 and PDS 241 in Fig. 1). Unfortunately, these sources are too weak in the mid-IR to conclusively establish the nature of their IR excess with VISIR spectroscopy. The shape of the feeble spectra of LkHα 215 and Hen 3-1110 are compatible, however, with free-free and bound-free emission from gas rather than the more typical dust-dominated spectra seen in the other studied targets (hence the “G” label in Table 1). It would be interesting to confirm that these are indeed classical Be stars in star-forming regions, because it would show that B type stars can rotate close to their break-up speed at a very young age (see also Acke et al. 2008).

A final remark concerns massive interacting binaries. These systems may also show a dusty disk around one of the stars in the system, or around the binary. However, these evolutionary phases are rare and we do not expect that our sample is significantly contaminated with such systems.

7.5. Emerging picture

The mid-IR imaging of our HBe stars revealed compact circumstellar material, most likely distributed in disks. The SEDs of the HBe stars display a higher diversity of IR excess shapes than the HAe stars. Also, the total IR excess is less prominent on average than that of the HAe stars. We interpret these differences as a variation in cloud contributions and geometrical differences of the disks. The disks around HBe stars are flatter and likely truncated on the outside. The typical disk size was constrained with the VISIR spectra to a FWHM upper limit of  ~500 AU. For comparison, Verhoeff (2009) found a typical FWHM of  ~50 AU for the HAes using the same method. These observations question the notion that the disks of HBe are scaled up versions of HAe stars. The observed 10 μm features show that the mineralogy of these disks is quite similar to that of their lower mass Ae counterparts. Just as for the HAes, we detect amorphous silicates, crystalline silicates, and emission from PAHs. Noteworthy differences are the IR extinction (e.g. MWC 297), the absence of a relation between the abundance of amorphous silica and forsterite, the lower contribution of PAH emission (e.g. MWC 147), and an absence of PAH-dominated sources.

The flux-levels, the spatial compactness and the underlying dust composition all point to the fact that, as for the HAe stars, the IR emission of Herbig Be stars is typically characterized by a circumstellar disk that efficiently reprocesses a substantial portion of the stellar flux. The picture described above of the HBe disks as being truncated is consistent with modeling results of Gorti & Hollenbach (2009), who find that FUV photoevaporation predominantly removes gas from the outer disk and that more luminous stars lose their disks rapidly. This could also explain the absence of a silica-forsterite relation, because this would require a fair deal of dust evaporation followed by condensation. Furthermore, it can explain the lower PAH flux contribution we found, because this emission is known to come from radii larger than the thermal silicate emission (see e.g. Habart et al. 2004; van Boekel et al. 2004; Lagage et al. 2006; Visser et al. 2007; Verhoeff 2009).

8. Conclusions

We have composed SEDs from literature photometry and have derived stellar parameters for a sample of HBe stars. We performed N-band imaging and long-slit spectroscopy with VISIR at the VLT. We analyzed these data to determine the composition and geometry of the circumstellar dust. The main conclusions we can draw from the results of our study can be summarized as follows:

• The HBe stars are surrounded by nebulae that manifestthemselves in increased optical extinction, a difference betweenVISIR and IRAS photometry, strong far-IR excess, high PAHflux of wide field spectroscopy, and large spatial extent of the PAHemission in VISIR spectroscopy.

• The IR emission of the HBe stars can be characterized by a circumstellar disk that efficiently reprocesses a substantial portion of the stellar flux.

• The disks around HBe stars appear to be comparable to those around the lower mass HAe counterparts in terms of flux levels, mineralogy, and spatial compactness.

• Compared to HAe stars, HBes typically have a weaker overall IR excess, and the sample contains a higher variety of IR excess shapes.

• The 12 μm reprocessing power, L₁₂/L _∗ , of Herbig Ae/Be stars goes down with stellar luminosity.

• The 10 μm spectra show that HBes have similar silicate abundances to those of the HAes.

• The disks around HBe stars generally have a lower PAH reprocessing power, F_PAH/F _∗ , than HAe disks.

• The abundance of PAH molecules does not have to be different for disks around HBe and HAe stars.

• The 10 μm spectra of the HBes do not show a relation between the mass fraction in crystalline grains and the mass fraction in large grains, as is the case for HAe stars.

• The 10 μm spectra of the HBes do not show a relation between the abundance of amorphous silica and forsterite, suggesting a history of dust vaporization followed by dust condensation.

• In the spectroscopic mode we spatially resolved 5 out of 17 HBe stars in either the continuum or in the 8.6 or 11.2 μm PAH-bands.

• We find upper limits to the mid-IR FWHM size of typically  ~ 500 AU.

• On the basis of the low mid-IR fluxes in our high spatial resolution photometry of LkHα 215, PDS 241, Hen 3-1110, MCW 1264, Hen 3-1418, PDS 152, MWC 593, and LkHα 115, it is possible that these stars are classical Be stars rather than HBe stars.

• The shape of the IR SED, the lower 12 μm reprocessing power, the weaker PAH emission, and the lack of a silica-forsterite relation can all be explained by assuming that HBes have flatter disks than the HAes and that they are truncated on the outside by photoevaporation.

Online material

Acknowledgments

This research was sponsored by NWO under grant number 614.000.411 and was partly supported by the Agence Nationale de la Recherche (ANR) of France through contract ANR-07-BLAN-0221. M. Min acknowledges financial support from the Netherlands Organization for Scientific Research (NWO) through a Veni grant.

References

All Tables

All Figures