A&A 415, 577-594 (2004)
DOI: 10.1051/0004-6361:20034622

Mid-IR observations of Galactic H II regions: Constraining ionizing spectra of massive stars and the nature of the observed excitation sequences

C. Morisset ^1,2 - D. Schaerer ^3,4 - J.-C. Bouret ² - F. Martins ^4,3

1 - Instituto de Astronomía, Universidad Nacional Autónoma de México, Apdo. postal 70-264, Ciudad Universitaria, México DF 04510, México
2 - Laboratoire d'Astrophysique de Marseille, CNRS, BP 8, 13376 Marseille Cedex 12, France
3 - Observatoire de Genève, 51 Ch. des Maillettes, 1290 Sauverny, Switzerland
4 - Laboratoire d'Astrophysique, UMR 5572, Observatoire Midi-Pyrénées, 14 Av. E. Belin, 31400 Toulouse, France

Received 3 June 2003 / Accepted 1 October 2003

Abstract
Extensive photoionization model grids for single star H II regions using state-of-the-art stellar atmosphere models have been computed to test their predicted ionizing spectra against recent ISO mid-IR observations of Galactic H II regions. Particular care has been paid to examining in detail the dependences of the nebular properties on the numerous nebular parameters which are generally unconstrained. Provided the ionization parameter U is fairly constant on average and the atomic data is correct these comparisons show the following:

• Both recent non-LTE codes including line blanketing and stellar winds (WM-Basic and CMFGEN) show a reasonable agreement with the observations, although non-negligible differences between their predicted ionizing spectra are found. Recurrently none of the models can be preferred over the other.
• The softening of the ionizing spectra with increasing metallicity predicted by the WM-Basic models is found to be too strong.
• We confirm earlier indications that the CoStar atmospheres, including an approximate treatment of line blanketing, overpredict somewhat the ionizing flux at high energies.
• Both LTE and non-LTE plane parallel hydrostatic atmosphere codes predict ionizing spectra that are too soft, especially over the energy range between 27.6, 35.0, and 41.1 eV and above. The inclusion of wind effects is crucial for accurate predictions of ionizing fluxes.

These conclusions are found to be robust to effects such as changes of U, stellar metallicity changes, and the inclusion of dust. Uncertainties due to atomic data (especially for Ar) are discussed. We also discuss the difficulties in estimating absolute stellar temperatures from mid-IR line ratios. Finally we have examined which parameters are chiefly responsible for the observed mid-IR excitation sequences. The galactic gradient of metallicity changing the shape of the stellar emission is found to be one of the drivers for the excitation sequence of Galactic H II regions, the actual contribution of this effect being finally atmosphere model dependent. The observed excitation scatter can be explained by effects due to statistical sampling of the IMF leading to a ${T}_{{\rm eff}}$ dispersion plus additional dispersion of U.

Key words: ISM: abundances - ISM: dust, extinction - ISM: HII regions - ISM: lines and bands - atomic data - stars: atmospheres

   
1 Introduction

Despite their paucity, hot massive stars are prominent contributors to the chemical and dynamical evolution of their host galaxies. Because of their intense nucleosynthesis, they process large amounts of material, on very short time scales. Furthermore, in addition to type II supernovae, of which they are progenitors, massive stars drive the dynamics and energetics of the ISM through their supersonic massive winds, thus affecting the subsequent star formation process in their surrounding environment. Their strong UV radiative fluxes ionize the ISM and create H II regions. The ionization structure of the latter is therefore, for the most part, controlled by the EUV radiation field of their massive stars content. In order to determine the properties of H II regions, it is therefore essential to understand the physical properties of massive stars and most importantly, to constrain their FUV and EUV (H-ionizing continuum) flux distribution. Yet, this part of the stellar spectrum is generally unaccessible to direct observations and it is crucial to find indirect tests to constrain it. In this context, nebular observations of H II regions combined with extensive grids of photoionization models including state-of-the-art model atmospheres offer the best opportunity to achieve this goal.

A large number of galactic H II regions have been observed with the ISO satellite (see e.g. Martín-Hernández 2002, and references therein). These spectra provide a wealth of spectral information, through fine-structure lines of ions whose ionization/excitation threshold are located below 912 Å. The shape of the SED in the EUV, and more specifically the number of ionizing photons in this region, is directly probed by ratios of successive ionization states such as [Ar III] 8.98 $\mu$m/[Ar II] 6.98 $\mu$m, [N III] 57.3 $\mu$m/[N II] 121.8 $\mu$m, [S IV] 10.5 $\mu$m/[S III] 18.7 $\mu$m, and [Ne III] 15.5 $\mu$m/[Ne II] 12.8 $\mu$m. Building line ratios diagrams for these species that are very sensitive to different parts of the flux distribution below the Lyman threshold allow one to derive informations on the actual spectral energy distribution at wavelengths usually unaccessible to direct observations. This not only provides valuable informations on the physical properties of the H II regions but on their stellar content as well. As a matter of fact, it is nowadays often used to estimate the spectral type of the ionizing source of single star H II regions, and offers a useful counterparts to more classical techniques of typing, based on optical or near-infrared absorption features (Mathys 1988; Kaper et al. 2002; Hanson et al. 1996; Watson & Hanson 1997).

On the other hand, modeling tools to analyze the photosphere and winds of hot, massive stars with a high level of accuracy and reliability have become available in recent years. In particular, major progress has been achieved modeling the stellar photosphere and stellar wind in a unified approach incorporating also a treatment of non-LTE line blanketing for the major opacity sources (Lanz & Hubeny 2003a,b; Hubeny & Lanz 1995; Hillier & Miller 1998; Pauldrach et al. 2001).

The impact of the first generation of atmosphere models including stellar winds and non-LTE line blanketing on nebular diagnostics was studied by Stasinska & Schaerer (1997) using the CoStar atmosphere models of Schaerer & de Koter (1997). This study showed already several improvements with respect to the widely used LTE models of Kurucz (1991). More recently Martín-Hernández et al. (2004); Martín-Hernández (2002) have investigated the metallicity dependence of the spectral energy distribution of O stars and the ionization structure of H II regions, using the CMFGEN code by Hillier & Miller (1998). They also compared the EUV fluxes from CMFGEN to those of the CoStar (Schaerer & de Koter 1997) and WM-Basic (Pauldrach et al. 2001) codes. They concluded that different treatment of line-blanketing between CoStar on the one hand and WM-Basic and CMFGEN on the other hand results in significant differences in the predicted EUV SEDs and ionizing fluxes.

In this context, it is of special interest to investigate how the different models available nowadays compare to each other in predicting nebular lines ratios. Similarly, it is of importance to test the role that a handful of various nebular parameters might have on the line ratios diagrams provided by ISO observations. The parameters influencing the ionization structure of a photoionized region are: 1) the geometry, the density distribution, the metallicity of the gas, and the possible absorption of the ionizing radiation by dust, 2) any physical quantity affecting the shape of the ionizing flux like, for example, the effective temperature of the ionizing star, its metallicity, the presence of a wind and the characteristics of the latter, 3) the hypothesis used to model the atmosphere like the number of elements taken into account, the treatment of the line-blanketing, etc. in summary: the physical ingredients and the related assumptions used to model the emitting atmosphere.

The present paper describes photoionization models performed with various atmosphere models, separating the effects of all these parameters. The paper is structured as follows: The various adopted atmosphere models are briefly described and compared in Sect. 2. The ionizing spectra from these models are then used as input to our photoionization code for the calculation of extended grids of nebular models (Sect. 3). The compilation of ISO observations of H II regions is described in Sect. 4. In Sect. 5 we compare our photoionization models to the observations and discuss the effect of changing parameters one by one on the excitation diagnostics. In Sect. 6 we test the validity of the different excitation diagnostics and softness radiation parameters for the determination of  ${T}_{{\rm eff}}$. The discussion takes place in Sect. 7. Our main conclusions are summarized in Sect. 8.

2 O star atmosphere models

 Among the key ingredients for the description of O star model atmospheres are the treatment of non-LTE effects, the inclusion of stellar winds, and a treatment of line blanketing (see e.g. Kudritzki et al. 1988; Gabler et al. 1989; Abbott & Hummer 1985). In recent years considerable improvements have been made in the modeling of these processes and model grids computed with various sophisticated atmosphere codes have become available (see e.g. the recent conference on "Stellar atmosphere modeling'', Hubeny et al. 2003). For our photoionization models, we adopt the ionizing spectra predicted from five different codes (Kurucz, TLUSTY, CoStar, WM-Basic, CMFGEN) briefly described hereafter. With the exception of the TLUSTY and Kurucz models, which assume a plane parallel geometry and thus no wind, all models describe the photosphere and winds in spherical geometry, in a unified manner.

Except mentioned otherwise, we have used atmosphere models computed for solar abundances: He, C, N, O, Ne, Si, S, Ar and Fe being 0.1, 4.7 $~\times 10^{-4}$, 9.8 $~\times 10^{-5}$, 8.3 $~\times 10^{-4}$, 5.4 $~\times 10^{-5}$, 4 $~\times 10^{-5}$, 1.6 $~\times 10^{-5}$, 6.8 $~\times 10^{-6}$ and 4 $~\times 10^{-5}$ resp.

   
2.1 Kurucz models

We use the well-known plane parallel LTE line blanketed models of Kurucz (1994,1991). Computations were done for models with ${T}_{{\rm eff}}$ (and $\log~(g)$) between 26 and 50 kK (3.0 and 5.0). For stability reasons, the available high ${T}_{{\rm eff}}$ models are restricted to cases of high gravity. The employed Kurucz models are therefore representative of dwarfs rather than supergiants mostly considered for the other model atmospheres (cf. below).

   
2.2 TLUSTY

A grid of plane-parallel non-LTE line blanketed models based on the TLUSTY code of Hubeny & Lanz (1995) has recently been calculated using a super-level approach and an Opacity distribution function or a modified opacity sampling (Lanz & Hubeny 2003a,b). About 100 000 individual atomic levels have been included, for more than 40 ions of H, He, C, N, O, Ne, Si, P, S, Fe and Ni, using a superlevel approach. For all the models, a standard microturbulent velocity $V_{\rm turb}=10$ km s^-1 has been used. The parameter space of the grid covers $27~500~{\rm K} \leq T_{\rm eff} \leq 55~000$ K with 2500 K steps and  $3.0 \leq \log~(g) \leq 4.75$ with 0.25 dex steps. Up to 10 different metallicities, from 2 times solar to metal free chemical composition, have been considered by Lanz & Hubeny (2003a,b). We extracted from this database models with  ${T}_{{\rm eff}}$ ($\log~(g)$) ranging from 30 to 50 kK (3.0 to 4.0), with solar metallicity.

   
2.3 CoStar

The CoStar models of Schaerer & de Koter (1997) include stellar winds, treat H-He in full non-LTE, and include line blanketing effects with an opacity sampling method based on Monte-Carlo simulations (Schmutz 1991). The impact of these effects on the ionizing fluxes and nebular diagnostics of H II regions has been discussed in detail by Stasinska & Schaerer (1997).

For our computations we use CoStar models with the lowest value for $\log~(g)$, i.e. models D5, D4, D3, E3, F3, F2 and F1 from the CoStar grid of Schaerer & de Koter (1997). The ${T}_{{\rm eff}}$ (and $\log~(g)$) range from 27 kK (2.9) to 53 kK (4.1).

   
2.4 WM-Basic

The WM-Basic models of Pauldrach et al. (2001) treat a large number of ions in non-LTE and include their line blocking effect by means of an opacity sampling technique. The atmospheric structure is computed from the hydrodynamic equations including radiative acceleration from numerous metal-lines and continua. We used the grid available on the web and described in Pauldrach et al. (2001) for Supergiant models with  ${T}_{{\rm eff}}$ (and $\log~(g)$) ranging from 30 kK (3.0) to 50 kK (3.9).

   
2.5 CMFGEN

We have constructed spherically symmetric wind models, using the non-LTE, comoving frame code CMFGEN (Hillier & Miller 1998). This code solves the radiative transfer equation, together with the statistical equilibrium equations, and line blanketing is self-consistently taken into account, using a super-level formulation. The chemical elements included in our model calculations are H, He, C, N, O, S, Si and Fe. For the 28 ions explicitly treated, a total of 2292 levels distributed in 819 superlevels are included, representing 22 762 bound-bound transitions. Atomic data for Fe IV and Fe V have been slightly improved, compared to those first introduced in CMFGEN (Hillier & Miller 1998) and made consistent with those used in TLUSTY. As shown in Bouret et al. (2003), this was required to get a very good agreement in the determination of iron abundances, when fitting lines of these two successive ionization stages in individual O stars in NGC 346, the largest H II region in the SMC. The temperature structure is calculated from the assumption of radiative equilibrium. The atmospheric structure consists of the wind, parameterized by the classical $\beta$-law, which is connected to hydrostatic layers obtained from the ISA-Wind code of de Koter et al. (1996), such that at the connecting point both the velocity and velocity gradient match. We assume a constant Doppler profile of 20 km s^-1 for all lines. As shown by Martins et al. (2002) for dwarfs and by additional test calculations this assumption leads to negligible changes of emergent spectrum. The stellar parameters, including the abundances, used to compute the CMFGEN grid of supergiants are identical to those used by WM-Basic and described in Pauldrach et al. (2001).

   
2.6 Atmosphere models for Dwarf stars

The main results of the present paper are obtained for Supergiant stars. Nevertheless, we also computed grids of photoionization models using Dwarf stellar atmosphere models as ionizing spectrum, to check the effect of $\log~(g)$ on the excitation of the nebula. In this purpose, WM-Basic D models from Pauldrach et al. (2001) have been used. We have also computed CMFGEN models using the same set of parameters than those used for the WM-Basic D models.

The models using Dwarf stellar atmosphere are discussed in Sect. 5.3.

   
2.7 Rebinning of the ionizing spectra

For subsequent use in our photoionization code NEBU (described in Sect. 3) the different atmosphere models have to be rebinned to the energy grid used in this code. The SEDs are first converted to the units used in NEBU (number of photons/eV/s/cm²). The SED is then interpolated on the NEBU grid, such as to preserve the integrated number of photons in each energy interval in NEBU. For most of the energy intervals, the number of points in the original stellar atmosphere domain is some tens to some hundreds, giving a good accuracy for the rebinning. Note that despite the much lower number of points used to describe the ionizing spectrum in NEBU, the results are reliable, as the most important quantities are the number of photons able to ionize the different ions. The grid points actually fully takes into account the discontinuities at the ionisation thersholds of the differents ions.

   
2.8 Comparing the EUV spectra

	Figure 1: Position in a $\log~(g)$ versus ${T}_{{\rm eff}}$ diagram of the Supergiant models and the Kurucz dwarf models used in this paper. See Fig. 2 for the line symbols. CoStar models are labeled according to the original naming convention. CMFGEN models have the same parameters than  WM-Basic models and are not drawn here.
Open with DEXTER

Figure 1 present all the Supergiant models used in this paper in a $\log~(g)$ versus ${T}_{{\rm eff}}$ diagram. The values for $\log~(g)$ at a given  ${T}_{{\rm eff}}$ are very close together, with the exception of the Kurucz models, which have a systematic higher value for $\log~(g)$, up to be even higher than the value adopted for Dwarf models (see also Fig. 1 in Schaerer & de Koter 1997).

Figure 2 illustrates the differences in the SED obtained from different atmosphere models after the rebinning procedure described above, for the same ${T}_{{\rm eff}}$, here 35 and 40 kK, with the exception of CoStar model for which no value at 35 kK is available in the Supergiant subset of models used here, model D4 at 32.2 kK is plotted. While the five models agree quite well in the domain of energies lower than 20 eV (and very well in the optical and IR range, not shown here), their differences can be as big as 4 orders of magnitudes just before 4 Rydberg. In this paper, we will use IR lines to trace the SED between 27, 35 and 41 eV (see next section), where the models differences already reach 1 to 2 orders of magnitude.

Of more interest for the analysis of the behavior of the fine-structure lines is the distribution of the ionizing photons at each energy. This is quantified by Q_E, which is the number of photons with energy greater than E, shown in right panels of Fig. 2. More precisely, the relevant quantity determining the nebular structure and properties would be the photon output weighted by the photoionization cross section. In the range of energy traced by the excitation diagnostics, 27-41 eV, the behavior of the four models is very different. We will discuss this further in Sects. 5.1 and 6.1.

Figure 2: Comparison between the 6 stellar atmosphere models: CoStar (solid), WM-Basic (dotted), CMFGEN (dashed), TLUSTY (dash dot), Kurucz (dash dot dot) and the Blackbody (long dashes, left panels only), for the same ${T}_{{\rm eff}}$ of 35 kK (upper plots), except for CoStar model (see text), and 40 kK (lower plots). The left panels show the Spectral Energy Distribution and the right panels show, for any energy E (eV), the number QE of photons with energy greater than E, relative to the corresponding number for the Blackbody emission, all the spectra having the same value for Q13.6. Vertical lines are plotted at 13.6 eV (solid) and 27.6, 35.0 and 41.1 eV (dotted), corresponding to the ionization potentials of the ions considered in this paper (Ar+, S++, and Ne+ resp.).

Open with DEXTER

   
3 Grid of photoionization models

Extensive grids of photoionization models were computed with the NEBU code (Morisset & Péquignot 1996; Morisset et al. 2002; Péquignot et al. 2001) in order to evaluate in detail the dependence of the mid-IR atomic fine-structure line emission of Galactic H II regions on the atmosphere models, and the stellar and nebular properties. Our main aims are a) to derive constraints on the stellar ionizing spectra and b) to examine the origin of the observed excitation gradients in (compact) Galactic H II regions.

In principle nebular emission line properties depend on fairly a large number of parameters, namely:

• the shape of the stellar ionizing spectrum, determined (mostly) by the stellar temperature ${T}_{{\rm eff}}$, gravity, and metallicity $Z_\star$;
• the ionization parameter $U(r)= Q_{13.6}/N_{\rm e} 4\pi r^2c$. As the geometry of the H II regions modeled here is an empty cavity surrounded by an ionized shell, we prefer to use hereafter the mean ionization parameter $\bar{U}=U(\bar{r})$, computed following Evans & Dopita (1985) at a distance from the ionizing star $\bar{r} = r_{\rm empty} + \Delta R$/2, where $r_{\rm empty}$ is the size of the empty cavity and $\Delta R$ the thickness of the H II shell. $\bar{U}$ is essentially given by the ratio of the ionizing photon density over the nebular particle density, i.e. properties of the ionizing source and the nebular geometry;
• the nebular abundances/metallicity  $Z_{\rm gas}$;
• atomic parameters driving the ionization equilibrium and line emissivities, e.g. photoionization cross section, recombination coefficients (radiative and dielectronic), collisional excitation cross sections, etc.;
• other secondary parameters like the presence of dust.

There is no doubt the existence of a systematic metallicity variation among the observed sources considered below. On the other hand, as for most of these objects (compact/ultra-compact H II regions) the properties of the ionizing source(s) and their geometry are not known, it is imperative to assess the impact of all parameters on the observables and to establish that the conclusions drawn from comparisons with observations are robust in this respect. The dependence of the emission line properties on the above parameters is examined in Sect. 5 with the help of photoionization models computed for a wide range of model parameters.

The bulk of "standard'' models were computed for the following cases/assumptions. The ionizing spectra from the five atmosphere models described in Sect. 2 and plotted in Fig. 2 plus blackbody spectra are adopted. Stellar ${T}_{{\rm eff}}$ ranging from 30 to 50 kK were used. This range in ${T}_{{\rm eff}}$ is likely to describe the physical conditions of the sample of H II regions (Morisset 2003). For most cases we assume a solar composition for the nebular and stellar abundances. Metallicity variations are considered in Sect. 5.5. For each of these stellar atmosphere, series of photoionization models were computed for the following nebular parameters. We set the electron density to 10³ cm^-3, one order of magnitude below the lowest critical density of the lines under consideration (cf. Martín-Hernández et al. 2002a). An empty cavity of radius 3 $~\times 10^{17}$ cm is assumed. The luminosity of the ionizing star is adjusted to lead to a constant number of Lyman continuum photons ( ${Q}_{13.6}=1.5 \times 10^{49}$ s^-1) corresponding to an ionization parameter $\log~(U) = -1.5$. Additional models quantifying the effect of variations of $\bar{U}$ are presented in Sects. 5.2, 5.4, and 6.2.

The effect of dust can be included in the photoionization computation, with two different optical properties corresponding to graphite or astronomical silicate (see Sect. 5.6).

The observables predicted from these extensive model grids will be compared to observations in Sect. 5.1.

   
4 The ISO H II regions catalogs

Infra-red spectra between 2.3 and 196 $\mu$m were taken from a sample of 43 compact H II regions using the two spectrometers (SWS and LWS) on board ISO (Peeters et al. 2002). Details about the data reduction and a first analysis of the ionic lines in terms of abundances can be found in Martín-Hernández et al. (2002a). Error bars on the lines intensities are within 10% to 20%. Note that a detailed study of one source has been achieved in Morisset et al. (2002).

Giveon et al. (2002b) published a catalog of 112 H II regions observed by ISO SWS spectrometer. Some of the sources are common with the Martín-Hernández et al. (2002a) catalog. The two catalogs are very coherent in terms of line intensities, as concluded by Giveon et al. (2002a), and are therefore included in our analysis. The effect of local and interstellar attenuation, even if lower in the IR range used in this work than for the optical domain, can be important and need to be corrected for. We correct the observed line intensities from the reddening using the extinction law described in Table 2 of Giveon et al. (2002a).

In the SWS and LWS spectral domain, 4 fine-structure line ratios are sensitive to the ionizing flux distribution: [Ar III] 8.98 $\mu$m/[Ar II] 6.98 $\mu$m, [N III] 57.3 $\mu$m/ [N II] 121.8 $\mu$m, [S IV] 10.5 $\mu$m/[S III] 18.7 $\mu$m, and [Ne III] 15.5 $\mu$m/[Ne II] 12.8 $\mu$m, hereafter [Ar III/II]/, [N III/II]/, [S IV/III]/, and [Ne III/II]/ respectively.

The excitation ratio implying nitrogen lines will not be used in the next discussion, since: 1) the two nitrogen lines are observed by LWS spectrometer, with a larger aperture size than the SWS: some nitrogen emission can arise from regions not seen in the other lines; 2) Giveon et al. (2002b) have observations only with SWS and then the number of observational constraints strongly decrease when using only Martín-Hernández et al. (2002a) data; 3) the critical densities of the nitrogen lines are very low compared to the one of the other lines (see Martín-Hernández et al. 2002a) and will not be emitted by medium density gas which can still emit the other lines, and 4) the ionization potential of N⁺⁺ is very close to the one of Ar⁺⁺ (29.7 and 27.6 eV resp.), so the main conclusions regarding the 30 eV energy domain will be obtained from argon lines.

Depending on the element, the number of sources for which we have finite value for the corrected excitation ratios is ranging from 45 to 51. Error bars on the line intensities are approximately 10 to 20%.

    
5 Excitation diagnostics

Figure 3: Dereddened observed values for the excitation sensitive line ratio [Ar III/II] versus [S IV/III] (the corresponding ionization potentials are also given). Source with a galactocentric distance lower than 7 kpc are symbolized with a +, otherwise with an X. Results from the photoionization model grid are line plotted using the same codes as in Fig. 2. The plot have been done such as the lowest ionization potential (indicated in braces) is always on the y-axis. Models obtained with 35 and 40 kK stars are shown using filled diamonds and empty squares respectively (except for CoStar model at 32.2 kK, empty diamond, see text). The y=x line is also drawn.

Open with DEXTER

	Figure 4: Same as Fig. 3 for the excitation sensitive line ratio [Ar III/II] versus [Ne III/II].
Open with DEXTER

	Figure 5: Same as Fig. 3 for the excitation sensitive line ratio [S IV/III] versus [Ne III/II].
Open with DEXTER

Figures 3 to 5 show the main results of the photoionization models using different atmosphere models for the ionizing star (lines), and the deredenned observed values, for the two merged catalogs (Giveon et al. 2002b; Martín-Hernández et al. 2002a). As we consider 3 diagnostic ratios, 3 plots can be drawn. The models obtained with ${T}_{{\rm eff}}$ = 35 and 40 kK are symbolized by filled diamonds and open squares respectively. The open diamond indicates a CoStar model at 32.2 kK as no model at 35 kK is available.

In principle the position of a model in such diagrams depends on all the following parameters: the hardness of the stellar SED (parametrized here for each set of model atmospheres by ${T}_{{\rm eff}}$ for a fixed stellar metallicity) and the main nebular parameters, i.e. the ionization parameter $\bar{U}$ and the nebular composition. Let us consider first the case of constant (solar) metallicity and constant $\bar{U}$ (but see Sects. 5.2, 5.4 and 5.5). In this case the location of a point on such diagnostic diagrams depends only on a) the global excitation of the gas and b) the "slope'' of the ionizing photon distribution between the two corresponding ionization potentials. For constant $\bar{U}$ and a given set of atmosphere models the excitation (a) is determined by  ${T}_{{\rm eff}}$. In other words, when ${T}_{{\rm eff}}$ increases, the number of ionizing photons at all the energies traced by the observed ions increases, and the points in the excitation diagrams will essentially move along the diagonal (y=x) direction. Note, however, that different atmosphere models with the same  ${T}_{{\rm eff}}$ predict fairly different absolute positions in these plots. This simply reflects the differences in the predicted number of ionizing photons above the relevant energy (cf. Fig. 2).

For a given line ratio the other line ratios depend to the first order on the "slope'' of the ionizing spectra (b). More precisely, the relevant quantity is the slope of the cumulative number of ionizing photon flux Q_E between the corresponding ionization potentials (see right panels in Fig. 2). For example, TLUSTY and Kurucz models show in general the softest spectra (i.e. steepest slopes) between 27.6 and 41.1 eV. For a given [Ne III/II] these models therefore show the highest [Ar III/II] values.

For the assumptions made here (constant $\bar{U}$ and metallicity) each location of the model results in the three excitation diagrams can be approximately understood in terms of the ionizing photon distributions Q_E. The correspondence is not always exact, as some competitive processes take place in the use of the ionizing photons, but the overall trends can be simply understood from the shapes of the spectra. Other additional assumptions (e.g. on the luminosity class of the exciting sources, the presence of dust, and uncertainties of the atomic data) also affect the predicted excitation diagrams. These effects are discussed below.

The observed excitations, correlated between the three excitation ratios [Ar III/II], [S IV/III], and [Ne III/II], can be decomposed into two components: an excitation sequence showing a global increase of the excitation ratios over $\sim$2 orders of magnitude, following to first order a trend parallel to y=x in the excitation diagrams, and a superposed excitation scatter of typically $\sim$0.5-1 order of magnitude around the mean excitation (cf. Figs. 3 to 5, but see also Figs. 19 and 20, where excitations versus metallicity are plotted).

   
5.1 First comparison with the observations

The zero-th order trend of the observations plotted in Figs. 3 to 5 is reproduced by the models: the excitation of the ionized gas, traced by the Xⁱ⁺¹/Xⁱ ratios, are well correlated. A ${T}_{{\rm eff}}$ sequence from $\sim$30 to 45 kK succeeds in reproducing the entire range of gas excitations. Note however, that, as discussed below (Sect. 7.2), this does not imply that the ionizing stars of our objects indeed cover this range of  ${T}_{{\rm eff}}$.

Fairly large differences are found in the predicted excitation diagnostic diagrams (Figs. 3 to 5) when using different atmosphere models. As expected from the intrinsic SEDs, the largest differences are found in Fig. 4, which traces the largest energy domain ($\sim$28 to 41 eV) corresponding to the [Ar III/II] and [Ne III/II] ratios.

When taken literally (i.e. assuming a fixed constant value of $\bar{U}$ for all atmosphere model sets and a fixed solar metallicity) Figs. 3 to 5 indicate the following concerning the shape of the ionizing spectra.

1.

Both recent codes, WM-Basic and CMFGEN, show a reasonable agreement with the observations. Given their different behavior in the three excitation diagnostics depending on luminosity class (cf. Sect. 5.3), and other remaining uncertainties discussed subsequently, it seems quite clear that none of the models can be preferred on this basis.

Despite these similarities we note, however, that an important offset is found in the excitation ratios predicted by these codes for a given absolute value of ${T}_{{\rm eff}}$ (cf. Sect. 6.1).

2.

Both plane parallel hydrostatic codes (Kurucz, TLUSTY) predict spectra which are too soft, especially over the energy range between 27.6 and 41.1 eV and above. For the Kurucz LTE models this problem is just a manifestation of the well-known "Ne  III'' problem (Simpson et al. 1995; Stasinska & Schaerer 1997; Rubin et al. 1988). This problem is persistently found with the non-LTE TLUSTY models. Although not completely clear, this "softness'' is likely due the neglect of wind effects which are known to alter the ionizing spectrum (cf. Stasinska & Schaerer 1997; Gabler et al. 1989) albeit in fairly complex way involving line blanketing from large numbers of metal-lines.

3.

As already found in other investigations (cf. Schaerer 2000; Oey et al. 2000, and references therein) we see that the CoStar models (showing the hardest spectra among the ones considered here) overpredict somewhat the ionizing flux at high energies. This likely overestimate of the ionizing flux at high energy (cf. below) must be due to the approximate and incomplete treatment of line blanketing (see also Martín-Hernández et al. 2004; Crowther et al. 1999).

4.

Interestingly blackbodies reproduce best the observed excitation diagrams, which indicates a posteriori that the ionizing spectra should have relative ionizing photon flux productions (Q_E at energies 27.6, 35.0 and 41.1 eV) close to that of blackbody spectra. This will be discussed in more detail in Sect. 7.

We note that the results concerning CMFGEN, WM-Basic, CoStar, and Kurucz models presented here confirm and support those presented by Stasinska et al. (2002) in terms of position of the models between each other in the [Ar III/II] versus [Ne III/II] excitation diagram and distance to the H II regions observations (cf. below).

   
5.2 Main dependences: Stellar temperature, ionization parameter, metallicity

For clarity it is useful to discuss first the dependences of the excitation diagnostics on the main parameters, i.e. the stellar temperature, ionization parameter, and metallicity. This is illustrated here somewhat schematically for the case of the [Ne III/II] versus [Ar III/II] diagnostic. Qualitatively the same results are found for the other excitation diagrams.

An increase of the stellar temperature ${T}_{{\rm eff}}$ or ionization parameter $\bar{U}$ or a decrease of the metallicity all lead overall to a higher excitation of the nebula, which e.g. manifests itself by larger [Ne III/II] and [Ar III/II] line ratios. However, although both line ratios change in similar ways, their effect is distinguishable to some extent. This is illustrated in Fig. 6, which shows for blackbody (and WM-Basic, see Sect. 5.5) spectra the implied shift in the [Ne III/II] versus [Ar III/II] excitation diagnostics due to a change of  ${T}_{{\rm eff}}$, $\bar{U}$ and  $Z_{\rm gas}$(consistent changes of both the stellar and nebular metallicity are discussed in Sect. 5.5). From this figure we see that ${T}_{{\rm eff}}$ and $\bar{U}$ variations are not completely degenerate (i.e. "parallel'').

Also, an increase of the nebular metallicity leads to a (very) small decrease of the excitation diagnostics, quite parallel to the variations induced by changing $\bar{U}$. The effect of changing coherently the metallicity of both the H II region and the star (which is principally acting on excitation diagnostics) are considered in Sect. 5.5.

The effect of continuum absorption by dust inside the H II region on the excitation diagnostics is also parallel to changes of  ${T}_{{\rm eff}}$, and is discussed in Sect. 5.6. Although quantitatively these variations depend e.g. on the adopted SED (and of the point in the ${T}_{{\rm eff}}$-$\bar{U}$-Z space chosen to compute the partial derivatives traced by the arrows in Fig. 6), these qualitative distinctions remain valid for the entire parameter space considered in the present paper and will be useful for the discussions below.

Note that earlier investigations have considered that changes of $\bar{U}$ are completely degenerate with ${T}_{{\rm eff}}$ (Giveon et al. 2002b; Martín-Hernández et al. 2002b), or have simply assumed an arbitrarily fixed, constant value of $\bar{U}$ (Martín-Hernández et al. 2004).

Figure 6: Increase of the excitation diagnostics [Ar III/II] versus [Ne III/II], when $\log(\bar{U})$ is increased from -2.6 to -0.8 (solid arrows), when the metallicity of the gas is increased from half solar to twice solar (dashed arrows), and when  ${T}_{{\rm eff}}$ is increasing from 35 to 40 kK (dotted arrows). Upper set of arrows are for BlackBody models, and lower set for WM-Basic models. In the latter case the metallicity is changes coherently for both the gas and the star. Codes are the same than in Figs. 2, 3.

Open with DEXTER

As apparent from Fig. 6, [Ar III/II] varies little with $\bar{U}$ in comparison with other mid-IR excitation ratios. This behavior is easily understood, as explained in the following brief digression.

 In ${T}_{{\rm eff}}$ and $\bar{U}$ domain used in this work, the [Ar III/II] ratio is controlled by the helium equilibrium: the IPs of Ar⁺ and He⁰ are closed together (27.7 and 24.6 eV resp.); in this energy domain the photon dominant predator is He⁰. The Ar⁺⁺ region (Ar⁺) is then cospatial with the He⁺ (He⁰) region (Some Ar⁺ can also be present in the He⁺ region, depending on the ionization parameter). The H II regions modeled here are all radiation bounded, the size (and the emission) of the He⁺ and Ar⁺⁺ region is mainly proportional to Q_24.6, while the size of the He⁰ and Ar⁺ region is controlled by the size of the H II region removing the He⁺ region. [Ar III/II] is then mainly controlled by  Q_24.6/Q_13.6. The previous argumentation is valid only if the recombination of Ar⁺⁺ remains quite small, which is not the case when strong dielectronic recombination occurs. In this case, the Ar⁺ region penetrates inside the He⁺ region, and the [Ar III/II] is decreased (see Sect. 5.7 for the effects of dielectronic recombinations of Ar⁺⁺.) Nevertheless, using atmosphere models instead of BlackBody leads to a more important increase of [Ar III/II] while increasing $\bar{U}$, as seen in Fig. 6.

The extreme correlation between He⁺/H and Ar⁺⁺/Ar⁺ can be verified in Fig. 7, where He I 5876 Å/H$_\beta$ versus [Ar III/II] is plotted, for all the atmosphere models. While the He I 5876 Å/H$_\beta$ ratio saturate at a value between 0.1 and 0.2, the [Ar III/II] excitation diagnostic still evolve with  ${T}_{{\rm eff}}$. This will be discussed further in Sect. 6.1.

	Figure 7: Correlation between He I 5876 Å/H$_\beta$ and [Ar III/II]. The line styles are the same as in Fig. 2. Arrows have the same meaning as in Fig. 6, for WM-Basic models only.
Open with DEXTER

   
5.3 Comparing results of Supergiant and Dwarf stars

The results shown in previous sections are obtained for Supergiant stars. Models were also performed using Dwarfs stars (see description in Sect. 2.6). The decrease of $\log~(g)$ changes the shape of the ionizing radiation, as shown in Fig. 8, where Q_E is shown for Supergiants and Dwarfs atmosphere models obtained using CMFGEN and WM-Basic, all at 35 kK. The main effect of increasing $\log~(g)$, observed on the two models, is to decrease strongly Q_E at 41 eV, and to increase the slope of Q_E between 27 and 41 eV. This overall hardening of the ionizing flux for stars with lower gravity is due mostly to an increased ionization in the continuum forming layers, the latter effect resulting from the increased wind density (mass-loss rate).

Figure 9 shows the excitation diagram [Ne III/II] versus [Ar III/II] performed using the models described above, comparing the Supergiant (light curves) and Dwarf (bold curves) results for both CMFGEN (dashed) and WM-Basic (dotted).

As expected from the increased hardness of the ionizing fluxes for supergiants, the use of dwarf atmospheres leads in general to an excitation decrease which is more important at the highest energies (i.e. [Ne III/II] decreases more rapidly than [Ar III/II]). Overall the differences between supergiant and dwarf spectra do not importantly affect our conclusions.

	Figure 8: Comparison between energy distribution (same as right panels of Fig. 2) between Supergiants (light curves) and Dwarfs (bold curves) of WM-Basic (dotted) and CMFGEN (dashed) models, all at 35 kK.
Open with DEXTER

   
5.4 Effect of changing the mean ionization parameter $\bar{U}$

Given that the ionizing sources and the nebular geometry of the observed objects are essentially unknown and $\bar{U}$ therefore as well, it is important to examine how robust the above results are with respect to changes of $\bar{U}$. E.g. is it possible to reconcile the discrepant predictions using the Kurucz and TLUSTY atmosphere models (i.e. to increase the predicted [Ne III/II] ratio of a given [Ar III/II]) by invoking a larger ionization parameter toward the low excitation end of the observed sequence? As shown in Sect. 5.2 and by detailed model calculations, variations of $\bar{U}$ imply changes nearly parallel to the "standard'' sequences for the Kurucz and TLUSTY atmosphere models, similar to the case of WM-Basic models shown in Fig. 6.

	Figure 9: Comparison between excitation diagnostic [Ar III/II] obtained with Supergiants (light curves) and Dwarfs stars (bold curves), for WM-Basic (dotted) and CMFGEN (dashed) atmosphere models. Models at 35 and 40 kK are shown by filled diamonds and empty squares respectively.
Open with DEXTER

We can therefore quite safely conclude that variations of $\bar{U}$ cannot reconcile the Kurucz and TLUSTY atmosphere models with the observation. However, for the other model atmospheres showing more "curved'' predictions in the excitation diagrams, variations of $\bar{U}$ might be invoked to improve/alter the predicted sequences for constant $\bar{U}$.

   
5.5 Effect of adopting a metallicity gradient, for both the stellar atmosphere and the ionized gas

It is well known that the metallicity decreases in the Galaxy when the distance to the center increases (e.g. Martín-Hernández et al. 2002a; Giveon et al. 2002a, and references therein). The metallicity varies approximately by a factor 4 from the center out to 15 kpc, where the most external regions used in this work are located (cf. Fig. 19).

Figure 6 shows the effect of changing the metallicity Z coherently in the ionizing star's SED computation and of the ionized gas in the photoionization computation, from half solar to twice solar. A metallicity increase in the atmospheres leads to a stronger blanketing with the effect of softening the EUV spectra of early type stars (e.g. Pauldrach et al. 2001), leading thus to a lower nebular excitation. As seen by comparing the Z-arrows in Fig. 6 and by additional test calculations, the increase of the nebular abundance plays a minor role in the resulting excitation shift.

In fact the WM-Basic models used here show that the ionizing spectra soften too strongly with increasing metallicity, leading to a stronger reduction of [Ne III/II] compared to [Ar III/II]. This results in a progressive shift away from the observed sequence toward higher metallicity. This discrepant trend, also found by Giveon et al. (2002b), was actually used by these authors to argue that the observed excitation sequence was mostly driven by ${T}_{{\rm eff}}$ variations. However, as abundances of these sources are known to vary by approximately the same factor as the Z variations considered here for the WM-Basic models, metallicity cannot be neglected. Therefore the discrepancy between the observations and the expected changes of [Ne III/II] and [Ar III/II] show that the predicted softening of the WM-Basic ionizing spectra with metallicity at high energies ($\ga$41 eV) is probably incorrect. Alternate solutions to this puzzle include postulating a increase of $\bar{U}$ toward higher Z (cf. above), or processes currently not accounted for in the WM-Basic models altering the high energy part of the SED (cf. Sect. 5.8), or changes in atomic physics parameters (cf. Sect. 5.7).

   
5.6 Effect of the dust

The effect of the presence of dust inside an H II region is firstly to decrease the global amount of ionizing photons from the point of view of the ionized gas, the effect being then to reduce the ionization parameter. On the other hand, the efficiency of dust in absorbing of the ionizing photons inside the H II region decreases with the energy of the photons after about 18 eV (e.g. Aannestad 1989; Mathis 1985), the global effect being to increase the excitation of the gas when increasing the amount of dust, for a given ionization parameter. As already pointed out by Morisset et al. (2002) for the case of ISO observations of G29.96-0.02, if dust is present in H II regions, quite the same excitation of the gas will be recovered using an higher ionization parameter and a lower  ${T}_{{\rm eff}}$. As illustration, inclusion of Graphite and Astronomical Silicate dusts, in proportion of 2.5 $~\times 10^{-3}$relative to hydrogen, for each type of dust, leads to an increase of all the excitation diagnostic ratios by a factor close to 2. The excitation increase due to dust is found to be "parallel'' to a ${T}_{{\rm eff}}$ increase to reasonable accuracy, when keeping $\bar{U}$ constant.

   
5.7 Effect of the atomic data

Could uncertainties in the atomic data affect the results? Indeed, from Figs. 3 to 5, we could suspect the argon ionization equilibrium to be wrong, favoring the emission of [Ar III]. This could e.g. be due to an overestimation of the ionizing flux at 27.6 eV with respect to higher energies, or to a systematic error in the observed intensities of one of the two lines involved in the [Ar III/II] ratio ([Ar III] 8.98 $\mu$m is affected by silicate band). However, we cannot exclude also the effect of atomic data in the photoionization computation. Collisional rates are generally believed to be accurate within 20%, while our knowledge of recombination coefficients are less probant. Dielectronic recombination coefficients for the elements of the third row of the periodic Table are poorly known, and even the new computations done today are usually only for the first and second rows, corresponding to highly charged elements of the third rows, which is not the case for Ar⁺⁺ (see e.g. Savin & Laming 2002). Dielectronic recombination coefficients have been computed for less charged third row elements by e.g. Mazzotta et al. (1998), but only for high electron temperature (coronal gas), which is not the case for H II regions. Very recently, new dielectronic recombination rates have been computed for Ar⁺⁺ (Zatsarinny et al. 2003, private communication), but the results these authors obtain have still to be checked (dielectronic recombination rates reach values as 10³ times the classical recombination rates for electron temperature closed to 10⁴ K!).

To simulate the effect of dielectronic recombination and charge transfer reactions we have multiplied the classical recombination coefficient for Ar⁺⁺ by factors up to 20. Figure 10 shows the effects of multiplying this coefficient arbitrarily by 10, on the excitation diagnostics [Ar III/II] versus [Ne III/II]. The figure shows results using WM-Basic and CMFGEN models, but the same effect can be observed with any of the atmosphere models used in this paper.

The increase of the recombination coefficient improves overall the agreement with the observed sequence, although new discrepancies appear now at the high excitation end. However, no dependence of the dielectronic recombination coefficient on the electron temperature $T_{\rm e}$ has been taken into account here. As $T_{\rm e}$ is known to vary along the excitation sequence this could in fact "twist'' the global shape of the predicted excitation sequence. Currently both the exact "direction'' and importance of this effect remain, however, unknown. Note, that the uncertainties due to atomic data of Ar were already pointed out by Stasinska et al. (2002). We join these authors in encouraging atomic physicists to improve our knowledge of such data.

	Figure 10: Variation of the excitation diagnostics [Ne III/II] versus [Ar III/II] for the same atmosphere models (here WM-Basic and  CMFGEN, dotted and dashed thin lines respectively), multiplying the effective recombination coefficient for only Ar++ by 10 (bold lines).
Open with DEXTER

   
5.8 Effects not included in atmosphere models

What could be the limitations of the most sophisticated atmosphere models currently available, capable of altering the excitation diagnostics discussed here? Although an exhaustive discussion is obviously not possible, one can suspect one major process, namely the presence of X-rays, to alter in a non-negligible way the ionizing spectra of O stars. This has been shown clearly by Macfarlane et al. (1994), and has been discussed later e.g. by Schaerer & de Koter (1997). The relative importance of X-rays compared to normal photospheric emission is expected to increase for stars with weaker winds and toward later spectral types. In late O types their contribution can be non-negligible down to energies $\ga$30 eV, see e.g. Macfarlane et al. (1994) and a model at  $T_{\rm eff}=30$ kK by Pauldrach et al. (2001), with obvious consequences on nebular diagnostics. Regrettably few models treating the X-ray emission in O stars exist, their impact on the overall emergent spectrum including the EUV has hardly been studied with complete non-LTE codes including winds and blanketing, and their dependence with stellar parameters (wind density, stellar temperature, even metallicity?!) remains basically unknown.

For now, we can only qualitatively expect the inclusion of X-rays to harden the ionizing spectra, probably down to the energy range probed by (some) mid-IR diagnostics. While this could in principle improve some difficulties observed by the CMFGEN and WM-Basic models (e.g. increasing [Ne III/II]) their precise effect remains open.

   
6 ${T}_{{\rm eff}}$ diagnostics and "second-order'' diagnostic diagrams

   
6.1 Using excitation diagnostic to determine ${T}_{{\rm eff}}$

	Figure 11: Variation of the excitation diagnostic [Ar III/II] according to the ${T}_{{\rm eff}}$ of the ionizing star. Codes are the same as in Figs. 2, 3. The set of arrows have the same meaning as in Fig. 6, for WM-Basic models only.
Open with DEXTER

Since the excitation of the gas increases with ${T}_{{\rm eff}}$, it is tempting to infer stellar temperatures from excitation diagnostic ratios. However, such an approach is intrinsically highly uncertain, as the nebular excitation is also strongly dependent on other parameters (see Sect. 3), such as the ionization parameter $\bar{U}$, which remain in most cases poorly known, cf. Mathis (1982) for optical lines and Schaerer & Stasinska (1999) for mid-IR ratios. These cautionary remarks should be kept in mind when e.g. using single line ratios or even several line ratios (e.g. Takahashi et al. 2000; Okamoto et al. 2001,2003), but see also Morisset (2003), to estimate stellar properties of individual objects from nebular observations. Tailored photoionization models including numerous constraints can lead to substantially different results and should clearly be the preferred method (see e.g. Morisset et al. 2002).

For illustration we show in Fig. 11 the dependence of the [Ar III/II] excitation ratio on  ${T}_{{\rm eff}}$ for a fixed $\bar{U}$ and metallicity. Other mid-IR excitation diagnostics show similar behaviors as can already be seen from various figures above, except that their dependence upon $\bar{U}$ are higher than for [Ar III/II], as already discussed in Sect. 5.2. The discussion of  ${T}_{{\rm eff}}$ and $\bar{U}$ determinations using mid-IR excitation diagnostics and the H II regions metallicities is developed in Morisset (2003).

	Figure 12: Variation of the optical excitation diagnostic He I 5876 Å/H$_\beta$ with ${T}_{{\rm eff}}$ for supergiant stars. Codes are the same as in Figs. 2, 3. The set of arrows have the same meaning as in Fig. 6, for WM-Basic models only.
Open with DEXTER

The most important conclusion from Fig. 11 is the important difference of the excitation of the gas ionized by WM-Basic (dotted line) and CMFGEN (dashed line) stars, for the same  ${T}_{{\rm eff}}$ and $\bar{U}$, even if the two types of models are showing the same behavior in Figs. 3 to 5. In other words, taking for example a value of 10. For [Ar III/II], we can determine a  ${T}_{{\rm eff}}$ of 39 kK using CMFGEN and a value of 45 kK using WM-Basic. This behavior can easily be understood when comparing the Q_E distribution, as shown in Fig. 2 and discussed in Sect. 2.8.

Two of the classical ways to constrain  ${T}_{{\rm eff}}$ from optical observations are to use the He I 5876 Å/H$_\beta$ ratio (e.g. Kennicutt et al. 2000) or the [O III] 5007 Å/[O II] 3727,29 Å ratio (Dors & Copetti 2003). The predictions for He I 5876 Å/H$_\beta$ are shown in Fig. 12. As for [Ar III/II], this line ratio is fairly independent of $\bar{U}$. This diagnostic line ratio saturates above $T_{\rm eff} \ga 40$ kK, when helium is completely ionized to He⁺. Note that even among WM-Basic and CMFGEN, which treat very similar physics, some differences in this  ${T}_{{\rm eff}}$ indicator remain. Furthermore note that the predicted  Q_24.6/Q_13.6 and hence He I 5876 Å/H$_\beta$ vary non-negligibly between dwarfs and supergiants (see e.g. Pauldrach et al. 2001; Smith et al. 2002). The predictions for [O III] 5007 Å/[O II] 3727,29 Å are presented in Fig. 13. As for the He I 5876 Å/H$_\beta$ ratio, the results differ strongly from one atmosphere model to another one. The Oxygen excitation diagnostic is also strongly sensitive to $\bar{U}$ and Z, the effect of Z being slightly more important than found by Dors & Copetti (2003), but these authors used WM-Basic dwarfs atmosphere models.

	Figure 13: Variation of the optical excitation diagnostic [O III] 5007 Å/ [O II] 3727,29 Å with ${T}_{{\rm eff}}$ for supergiant stars. Codes are the same as in Figs. 2, 3. The set of arrows have the same meaning as in Fig. 6, for WM-Basic models only.
Open with DEXTER

We can conclude from the above examples that any determination of  ${T}_{{\rm eff}}$ based on diagnostic ratios can be reliable only if the metallicity of the star is coherently taken into account, and if $\bar{U}$ is also determined at the same time, as shown in Morisset (2003).

   
6.2 Radiation softness parameters $\eta$

Following Vilchez & Pagel (1988), radiation softness parameters can be defined combining the excitation diagnostic ratios, namely: $\eta_{\rm Ar-Ne} = \frac{[Ar {\sc iii/ii}]}{[Ne {\sc iii/ii}]}$and so on for the other diagnostic ratios.

To first order (but see the discussion on [Ar III/II] in Sect. 5.2) an excitation ratio Xⁱ⁺¹/Xⁱ depends on the ionization parameter $\bar{U}$ and the hardness of the ionizing radiation and is given by

$$\frac{X^{i+1}}{X^i} \propto \bar{U} \frac{\int_{\nu(X^i)}^\in... ... {\rm d}\nu} = \bar{U} \frac{{Q}_{E(X^{i})}}{{Q}_{13.6}}\cdot$$ (1)

Therefore one has:

$$\eta_{X-Y} \propto \frac{{Q}_{E(X^{i})}}{{Q}_{E(Y^{i})}},$$ (2)

where E(Xⁱ) is the ionization potential of the ion Xⁱ. $\eta$ is thus in principle independent of the ionization parameter $\bar{U}$, and a measure of the "slope'' of the ionizing spectrum between the ionization energies  E(Xⁱ) and  E(Yⁱ) respectively. Therefore such quantities - used so far for optical lines only - have often been thought to be good estimators of stellar effective temperatures, see e.g. Garnett (1989), Kennicutt et al. (2000) but also Skillman (1989), Bresolin et al. (1999), Oey et al. (2002). For this reason we here explore whether the observed mid-IR fine structure lines offer, both from the empirical and theoretical standpoint, such a diagnostic power.

	Figure 14: Variation of the radiation softness parameter $\eta _{\rm Ar-Ne}$versus [Ar III/II]. The codes are the same than in Figs. 2, 3. Note that the model predictions for both axis depend also on the assumed ionization parameter $\bar{U}$. Arrows as in Fig. 6 for WM-Basic models.
Open with DEXTER

Figure 14 shows the variation of  $\eta_{\rm Ne-Ar}$versus [Ar III/II]. No correlation between these two observables is found. While a galactic gradient is found for [Ar III/II] (Giveon et al. 2002b; Martín-Hernández et al. 2002a), $\eta_{\rm Ne-Ar}$ and none of the other mid-IR $\eta$'s which can be constructed show a galactic gradient. This can also be seen from the correlations between the various excitation ratios plotted by Martín-Hernández et al. (2002a), their Fig. 10. Already this finding indicates empirically that mid-IR softness parameters do not carry important information on the stellar ionizing sources.

How do the models compare with the observed softness parameters? Given the considerable spread between photoionization models using different stellar atmospheres and the various discrepancies already found earlier, it is not surprising that overall a large spread is also found here (Fig. 14). Compared to the observations the blackbody SED seems again to fit best. The WM-Basic, CMFGEN and CoStar models are marginally compatible with the observations each one on one side, while the TLUSTY and Kurucz results are definitively far away from the observed values. Note, however, that both theoretical quantities plotted here (including $\eta$) depend also on the ionization parameter.

	Figure 15: Variation of the radiation softness parameter  $\eta _{\rm Ar-Ne}$according to the ${T}_{{\rm eff}}$ of the ionizing star. Codes are the same as in Figs. 2, 3. Arrows have the same meaning as in Fig. 6, for WM-Basic models.
Open with DEXTER

	Figure 16: Variation of the radiation softness parameter  $\eta _{\rm S-O}$ (defined by Vilchez & Pagel 1988, see note 8), according to the ${T}_{{\rm eff}}$ of the ionizing star. Codes are the same as in Figs. 2, 3. Arrows have the same meaning as in Fig. 6, for WM-Basic models, for WM-Basic models only.
Open with DEXTER

In Fig. 15 we illustrate theoretical predictions of mid-IR softness parameters as a function of  ${T}_{{\rm eff}}$ for the case of  $\eta_{\rm Ne-Ar}$. Note that most of the model atmospheres predict that  $\eta_{\rm Ne-Ar}$ becomes insensitive to ${T}_{{\rm eff}}$ above a certain value, here of the order of  $T_{\rm eff}\ga 35~000$ K (the exact value also depending on $\bar{U}$). Similar dependencies on the adopted model atmosphere and "saturation effects'' have also been found for the traditional optical softness parameter (cf. Kennicutt et al. 2000; Oey et al. 2002).

Furthermore model calculations also show that all the mid-IR $\eta$ depend quite strongly on the ionization parameter $\bar{U}$ and on metallicity, as shown here for the case of the WM-Basic models (see Fig. 17). From the theoretical point of view, and without tailored photoionization modeling including constraints on $\bar{U}$ and Z, the use of mid-IR $\eta$'s appears therefore highly compromised. Again, similar difficulties have also been found for the optical softness parameter, cf. Skillman (1989), Bresolin et al. (1999), Oey et al. (2002), but Kennicutt et al. (2000) and hereafter.

	Figure 17: Variation of the radiation softness parameter $\eta _{\rm S-Ne}$versus $\eta _{\rm Ar-Ne}$, when changing $\bar{U}$, ${T}_{{\rm eff}}$ or Z. Arrows codes are the same as in Fig. 6.
Open with DEXTER

Last but not least, as each $\eta$ implies 4 line intensities, the softness parameters are more sensitive to any observational uncertainty (as the attenuation correction or detector calibration), possible collisional effects, uncertainties in the atomic data etc. In view of all these considerations and the model results presented above, applications of softness parameters in the mid-IR appear therefore to be of very limited use.

For completeness we show in Fig. 16 the behavior of the traditional optical softness parameter  $\eta_{\rm O-S}$ for all model atmospheres. Again a considerable spread between different atmosphere models and a dependence on $\bar{U}$ and Z is found.

7 Discussion 

From the comparisons of the extensive model calculations with the observed excitation diagnostics we can draw some rather general conclusions about the shape of the ionizing spectra of early type stars and the on the nature of the ionizing sources of the Galactic H II regions.

   
7.1 Implications on stellar SEDs

Overall the "best fit'' SED to the observed excitation diagrams (Figs. 3 to 4) was found with blackbody spectra. Does this imply that the ionizing fluxes of hot stars are best described by the Planck function? The overall answer is no, but. The observations probe (to 1st order) the relative number of ionizing photons with energies above the relevant ionization potentials, i.e. what was called the "slope'' in Sect. 5. Therefore if we considered that all observed 3 line ratios are correctly reproduced by a blackbody this would imply that the 27.6-35.0 and 27.6-41.1 eV slopes (and hence also 35.0-41.1 eV) be equal to that of the Planck function (bb) of the same  ${T}_{{\rm eff}}$, i.e. ${Q}_{27.6}/{Q}_{35.0}= {Q}_{27.6}^{\rm bb}/{Q}_{35.0}^{\rm bb}$and  ${Q}_{27.6}/{Q}_{41.1}= {Q}_{27.6}^{\rm bb}/{Q}_{41.1}^{\rm bb}$. This would therefore represent three "integral'' constraints on the stellar SED. This is rather strong, but still leaves room for the detailed shape of the SED between these energies. Actually the agreement between blackbody spectra and [Ar III/II] vs. [Ne III/II] is better than diagrams involving [S IV/III]. This indicates that the constraint on the Q_27.6/Q_41.1 is better than at intermediate energies.

However, it is important to remember that these constraints on the underlying SED can be deduced only if the observed excitation diagram (Figs. 3 to 5) is essentially driven by a temperature sequence (i.e. $\bar{U}={\rm const}.$). Note also the effects of uncertainties on the atomic data, as discussed in Sect. 5.7.

   
7.2 On the origin of the observed excitation sequences

The origin of the observed mid-IR excitation sequences and their correlation with galactocentric distance has been discussed recently by Giveon et al. (2002b); Martín-Hernández et al. (2002b). As both studies present somewhat limited arguments a more general discussion is appropriate here. The basic question is whether variations of the stellar effective temperature or metallicity variations are responsible for the observed decrease of excitation toward the Galactic Center?

Giveon et al. (2002b) have noted from photoionization models using WM-Basic spectra that [Ne III/II] is predicted to decrease more rapidly than [Ar III/II] with increasing metallicity - a finding also confirmed here (cf. Fig. 6). However, as the observed sequence does NOT follow this trend, they conclude that the decrease of excitation must be due to a reduction of  ${T}_{{\rm eff}}$ as opposed to a softening of the stellar SED with increasing metallicity. Obviously as such this conclusion cannot be upheld as the same mid-IR observations, allowing fairly accurate abundance determinations, clearly establish the existence of a metallicity gradient (Giveon et al. 2002a). In fact, the apparent contradiction between the observed trend with metallicity and the one predicted with WM-Basic model atmospheres indicates quite likely that the stellar SEDs soften too quickly with increasing metallicity and/or that the ionization parameter in regions at small galactocentric distance must be larger than assumed (cf. Sect. 5.5).

In contrast to the above study, Martín-Hernández et al. (2002b) show a loose correlation between excitation and metal abundances (e.g. between [Ne III/II] and Ne/H), stress the importance of metallicity effects on the SED (see also Martín-Hernández et al. 2004), and conclude that at least partly the observed decrease of [Ne III/II] must be due to a softening of the stellar SEDs with increasing metallicity.

7.2.1 Stellar evolution effects

Figure 18: Average stellar effective temperature and dispersion as a function of time for metallicities Z=0.008 ( $1/2.5~Z_\odot$, solid lines) and 0.04 ($2~Z_\odot$, dashed) predicted for an ensemble of single star H II regions for a Salpeter IMF with $M_{\rm up}=120~M_\odot$ and a minimum Lyman continuum flux $\log~({Q}_{\rm lim}) \ge 48.$. Note the rather small differences in average ${T}_{{\rm eff}}$ and the large intrinsic dispersion for each given metallicity.

Open with DEXTER

From what we know, three effects are related to metallicity and must all be taken into account. First, higher metallicity is known in stellar evolution to lead to a cooler zero age main sequence and to an overall shift to cooler temperatures (e.g. Schaller et al. 1992). Second, blanketing effects in the atmospheres become stronger with increasing metallicity and lead to softer ionizing spectra (e.g. Sect. 5.5; Martín-Hernández et al. 2004; Pauldrach et al. 2001). Finally, an increased nebular abundance leads also to a somewhat lower excitation of the gas in the H II regions (cf. Sect. 5.5). The remaining questions are then a) which of these effects dominate and b) whether taken together they can indeed quantitatively reproduce the entire range of observed excitation variations in the Galactic H II regions.

Although the content of stellar ionizing sources of the objects considered is not known (but see Morisset 2003) we can estimate the ${T}_{{\rm eff}}$ variations expected from stellar evolution, e.g. by assuming a single ionizing source. We then perform Monte Carlo simulations of single star H II regions of different metallicities assuming that the ionizing stars are drawn from a Salpeter IMF with a given upper mass cut-off  $M_{\rm up}$. In order to compute the mean properties of these stars, such as their average ${T}_{{\rm eff}}$, the predicted variation with metallicity and their dispersion, the equivalent of a lower mass limit must also be specified. This is done by imposing a lower limit on the total Lyman continuum photon flux  ${Q}_{\rm lim}$. Only stars with $Q >{Q}_{\rm lim}$ at a given age are retained for this computation. In practice we use the Meynet et al. (1994) stellar tracks, we consider metallicities between $\sim$1/2 and 2 times solar, as indicated by the observed range of Ne/H or Ar/H abundances and we adopt $\log~({Q}_{\rm lim})=48.$, corresponding a typical lower limit for the H II regions of Martín-Hernández et al. (2002a). Very massive stars entering the Wolf-Rayet phase already on the main sequence are also excluded.

The resulting average and spread of ${T}_{{\rm eff}}$ as a function of age is shown in Fig. 18 for Z=0.04 and 0.008 respectively. This figure shows the following: first, the reduction of the average stellar temperature due to a metallicity increase by a factor 4 is rather modest, of the order of $\sim$3-4 kK. Second, for a given metallicity the predicted dispersion in  ${T}_{{\rm eff}}$ is larger; typically of the order of $\sim$3-9 kK for reasonable ages. Although the absolute values of these ${T}_{{\rm eff}}$ depend for obvious reasons on the exact choice of  ${Q}_{\rm lim}$, the ${T}_{{\rm eff}}$ differences and dispersion depend little on this value. From all the models considered above, the decrease of the mean  ${T}_{{\rm eff}}$ due to stellar evolution effects appears to be too small to explain the full range of excitations.

On the other hand the fairly large ${T}_{{\rm eff}}$ dispersion at a given Zwill induce quite important variations in the excitation. This effect probably dominates the observed excitation scatter (defined at the end of Sect. 5), as also suggested by the observations of a large spread of the excitation diagnostics Figs. 19 and 20 for a given metallicity measured here by Ne/Ne$_\odot$. Indeed the observed spread of [Ar III/II] and [Ne III/II] at a given Z is somewhat larger than or similar to the decrease of the mean excitation with increasing Z. Taken together these findings imply that while undeniably metallicity effects on stellar evolution and nebular abundances must be present, statistical fluctuations of the effective stellar temperature due to the IMF are likely the dominant source of scatter for the observed mid-IR excitation sequence of Galactic H II regions, while the excitation sequences must be predominantly driven by other effects which we will discuss now.

	Figure 19: [Ar III/II] versus metallicity, measured here by the Ne abundance (see text). The line is a linear fit to the observations in log-log space. Arrows as in Fig. 6, for WM-Basic models.
Open with DEXTER

	Figure 20: [Ne III/II] versus metallicity, as in Fig. 19. A very similar plot is obtained for [S IV/III].
Open with DEXTER

7.2.2 Stellar atmospheres and nebular effects

As apparent from our modeling (see Figs. 19 and 20) the effects of metallicity on the shape of the stellar ionizing spectra strongly alter the predicted excitation. The magnitude of the predicted effect is found to be comparable to the observed variation. Both these findings and the above results concerning stellar evolution effects indicate the Z dependence of the ionizing spectra is the main driver for the correlation of the excitation with galactocentric distance.

As this result depends on a specific set of model atmospheres (WM-Basic) a few words of caution are, however, necessary here. First, we note that the effects of  ${T}_{{\rm eff}}$ and Z are not exactly the same for [Ar III/II], [Ne III/II], and [S IV/III]. The predicted excitation variation (shown here for a change of  ${T}_{{\rm eff}}$ from 35 to 40 kK) tends to be somewhat larger (smaller) for [Ne III/II] and [S IV/III] ([Ar III/II]) than the observed variations. Second, it must be remembered that the WM-Basic atmosphere models employed here could predict too strong a softening with increasing Zas suspected from Fig. 20 and already discussed in Sect. 5.5. Despite these imperfections there is little doubt that the above result remains valid.

Finally we may also comment on excitation changes related to the ionization parameter. Again, as for  ${T}_{{\rm eff}}$, the above Monte Carlo simulations show small differences between the average ionizing photon flux Q_13.6 with metallicity, but a considerable dispersion for each Z. Combined with the observational fact of fairly similar nebular densities in our objects this could be a justification for a constant ionization parameter, at least on average. Random variations of $\bar{U}$ are, however, expected to contribute to the excitation scatter at a given metallicity.

In conclusion we see little doubt that the observed excitation sequence of Galactic H II regions is shaped by the joint effects of metallicity on stellar evolution, atmospheric line blanketing, and cooling of the ISM. From our investigations it seems, however, that metallicity effects on ionizing stellar flux is the dominant effect causing the excitation gradients while statistical fluctuations of  ${T}_{{\rm eff}}$ and $\bar{U}$ are likely the dominant source of scatter in the observed excitations. A more detailed study of possible systematic  ${T}_{{\rm eff}}$ and ionization parameter gradients is presented in Morisset (2003): no clear gradient of  ${T}_{{\rm eff}}$ nor $\bar{U}$ versus the galactocentric distance are found by this author, but some trends of increase (decrease) of  ${T}_{{\rm eff}}$ ($\bar{U}$) with the metallicity are observed. Importance of taking into account the effect of the metallicity on the stellar spectral shape is addressed.

8 Summary and conclusion 

We have presented results from extensive photoionization model grids for single star H II regions using a variety of recent state-of-the-art stellar atmosphere models such as CMFGEN, WM-Basic, TLUSTY, CoStar, and Kurucz models. Even among the two recent non-LTE line blanketed codes including stellar winds (WM-Basic and CMFGEN) the predicted ionizing spectra differ by amounts leading to observable differences in nebular spectra.

The main aim of this investigation was to compare these model predictions to recent catalogs of ISO mid-IR observations of Galactic H II regions, which present rich spectra probing the ionizing spectrum between $\sim$24 to 41 eV thanks to the measurements of [Ar III/II], [Ne III/II], and [S IV/III] line ratios. Particular care has been paid to examining in detail the dependences of the nebular properties on the numerous nebular parameters (ionization parameter $\bar{U}$, abundances, dust etc.) which are generally unconstrained for the objects considered here.

Most excitation diagnostics are found to be fairly degenerate, but not completely so, with respect to increases of  ${T}_{{\rm eff}}$, $\bar{U}$, a change from dwarf to supergiant spectra, a decrease of the nebular metallicity (Sects. 5.2 and 5.3), and the presence of dust in the H II region (Sect. 5.6). Each of these parameters increases the overall excitation of the gas, and in absence of constraints on them, a derivation of such a parameter, e.g. an estimate of the stellar  ${T}_{{\rm eff}}$ of the ionizing source, is intrinsically uncertain. In consequence, while for sets of objects with similar gas properties statistical inferences are probably meaningful, such estimates for individual objects must be taken with care.

Provided the ionization parameter is fairly constant on average and the atomic data is correct (but cf. below) the comparisons between the photoionization model predictions and the observations allow us to conclude the following concerning the different stellar atmosphere models (Sect. 5.1):

• Both recent non-LTE codes including line blanketing and stellar winds (WM-Basic and CMFGEN) show a reasonable agreement with the observations. Given their different behavior in the three excitation diagnostics, depending on luminosity class, and other remaining uncertainties, it appears that none of the models can be preferred on this basis.

• The plane parallel hydrostatic codes (Kurucz, TLUSTY) predict spectra which are too soft, especially over the energy range between 27.6, 35.0, and 41.1 eV and above. Although a good agreement is found for UV to optical spectra predicted by the hydrostatic TLUSTY code and the photosphere-wind code CMFGEN (Bouret et al. 2003; Hillier 2003) important differences are found in the EUV range probed by the present observations and photoionization models. Apparently the full non-LTE treatment of numerous elements accounted for by TLUSTY is insufficient to accurately predict the ionizing spectra at these energies, and the inclusion of stellar winds is imperative.

• We confirm the finding of earlier investigations (e.g. Oey et al. 2000) showing that the CoStar models overpredict somewhat the ionizing flux at high energies.

• Interestingly blackbodies reproduce best the observed excitation diagrams, which indicates that the ionizing spectra of our observed object should have relative ionizing photon flux productions Q_E at energies 27.6, 35.0 and 41.1 eV close to that of blackbody spectra. Although this integral constraint on the SED remains approximate, it should still be useful to guide future improvements in atmosphere modeling.

• Finally, the softening of the ionizing spectra with increasing metallicity predicted by the WM-Basic models is found to be too strong. As already apparent from observed correlations between excitation diagnostics probing various energies, the observed softening of the radiation field (in part due to metallicity) affects fairly equally the range between $\sim$27 and 41 eV (Martín-Hernández et al. 2002a) in contrast to the atmosphere model predictions, which soften most at the highest energies.

These conclusions are found to be fairly robust to effects such as changes of $\bar{U}$, nebular and stellar metallicity changes, and the inclusion of dust. We suggest that the main uncertainty which could alter the above conclusions is the poorly known atomic data for Ar⁺⁺(especially dielectronic recombination coefficients) as also pointed out by Stasinska et al. (2002). Reliable computations for such data are strongly needed. From the perspective of atmosphere codes probably the most important step toward improving the reliability of ionizing fluxes resides in a quantitative exploration of the influence of X-rays on the emergent spectra at lower energy.

The potential of mid-IR line ratios or "softness parameters'', defined in analogy to the well known $\eta$ parameter for optical emission lines, has been explored (Sect. 6.1). The following main results have been obtained:

• Given the non-negligible differences between the various atmosphere models it is not surprising that individual line ratios (e.g. [Ar III/II], [Ne III/II]) show quite different dependences on  ${T}_{{\rm eff}}$. We find that [Ar III/II] depends little on the ionization parameter as the ionization of Ar⁺ is closely coupled to that of He. This suggests that [Ar III] 8.98 $\mu$m/[Ar II] 6.98 $\mu$m should in principle be a fairly robust temperature indicator, provided the atmosphere models are sufficiently accurate up to $\sim$24-27 eV and the atomic data is reliable (cf. above).

  In comparison to He I/H indicators (e.g. He I 6678 Å/H$\beta$ or He I 2.06 $\mu$m/Br$\gamma$ in Lumsden et al. 2003; Kennicutt et al. 2000, and references therein), [Ar III/II] does not show a saturation effect but remain sensitive to  ${T}_{{\rm eff}}$ up to the highest temperature examined here ( $T_{\rm eff}\sim 50$ kK). Due to the large uncertainties of dielectronic recombination coefficient of Ar⁺⁺, Morisset (2003) prefer to use [S IV/III] and [Ne III/II] to determine ${T}_{{\rm eff}}$ and $\bar{U}$ simultaneously, for the H II regions used in this work.

• Both empirically and theoretically the mid-IR softness parameters which can be constructed from [Ar III/II], [S IV/III], and [Ne III/II] are found to provide little if any information on stellar temperatures if not used to determine $\bar{U}$ at the same time (Morisset 2003). Observationally little / no correlation is found between the $\eta$'s and excitation ratios. Furthermore, our photoionization models show a considerable dependence of $\eta$ on the ionization parameter. We therefore conclude that mid-IR $\eta$'s appear to be of limited diagnostic power.

Finally we have examined which parameter(s) is (are) chiefly responsible for the observed mid-IR excitation sequences. Combining the results from our extensive photoionization model grids with Monte Carlo simulations of ensembles of single star H II regions of different metallicity and age we conclude the following (Sect. 7). While metallicity effects on stellar evolution, atmospheres and the nebulae all have an undeniable influence, they are probably of minor importance compared to the fairly large dispersion of  ${T}_{{\rm eff}}$ expected at each metallicity from a simple statistical sampling of the IMF. The ${T}_{{\rm eff}}$ scatter plus additional scatter in the ionization parameter are probably the dominant driver for the observed mid-IR excitation scatter of Galactic H II regions, while the effect of metallicity on the shape of the ionizing spectra is partially responsible of the global excitation sequence, the proportion of this effect being strongly dependent of the reliability of the atmosphere models (Morisset 2003).

Acknowledgements

We wish to thank various persons who have contributed directly and indirectly to the shape of this paper, over its rather long gestation time. Among those are in particular Grazyna Stasinska and Daniel Péquignot who have been consulted on photoionization codes and atomic physics, and Leticia Martín-Hernandez. Thierry Lanz kindly provided model atmosphere results prior to publication. We are also grateful to John Hillier for assistance in adapting and modifying his atmosphere code to our purpose. We thank Ryszard Szczerba for helping introducing dust in NEBU code. DS thanks the CNRS, the French Programme National de Galaxies, and the Swiss National Fund for Research for support. DS and FM thank the CALMIP and IDRIS centers for generous allocation of computing time. CM thanks the IA/UANM in Mexico for offering him the opportunity to finish this work and investigate many others in the future. This work was partly supported by the CONACyT (México) grant 40096-F.

References

 
• Aannestad, P. A. 1989, ApJ, 338, 162 In the text NASA ADS  
• Abbott, D. C., & Hummer, D. G. 1985, ApJ, 294, 286 In the text NASA ADS  
• Bouret, J., Lanz, T., Hillier, D. J., et al. 2003, ApJ, 595, 1182 In the text NASA ADS  
• Bresolin, F., Kennicutt, R. C., & Garnett, D. R. 1999, ApJ, 510, 104 In the text NASA ADS  
• Crowther, P. A., Pasquali, A., de Marco, O., et al. 1999, A&A, 350, 1007 In the text NASA ADS  
• de Koter, A., Lamers, H. J. G. L. M., & Schmutz, W. 1996, A&A, 306, 501 In the text NASA ADS  
• Dors, O. L., & Copetti, M. V. F. 2003, A&A, 404, 969 In the text NASA ADS  
• Evans, I. N., & Dopita, M. A. 1985, ApJS, 58, 125 In the text NASA ADS  
• Gabler, R., Gabler, A., Kudritzki, R. P., Puls, J., & Pauldrach, A. 1989, A&A, 226, 162 In the text NASA ADS  
• Garnett, D. R. 1989, ApJ, 345, 282 In the text NASA ADS  
• Giveon, U., Morisset, C., & Sternberg, A. 2002a, A&A, 392, 501 In the text NASA ADS  
• Giveon, U., Sternberg, A., Lutz, D., Feuchtgruber, H., & Pauldrach, A. W. A. 2002b, ApJ, 566, 880 In the text NASA ADS  
• Hanson, M. M., Conti, P. S., & Rieke, M. J. 1996, ApJS, 107, 281 In the text NASA ADS  
• Hillier, D. J. 2003, in Stellar Atmosphere Modeling, ed. I. Hubeny, D. Mihalas, & K. Werner, ASP Conf. Ser., 288 In the text  
• Hillier, D. J., & Miller, D. L. 1998, ApJ, 496, 407 In the text NASA ADS  
• Hubeny, I., & Lanz, T. 1995, ApJ, 439, 875 In the text NASA ADS  
• Hubeny, I., Mihalas, D., & Werner, K. 2003, Stellar Atmosphere Modeling, ASP Conf. Ser., 288 In the text  
• Kaper, L., Bik, A., Hanson, M. M., & Comerón, F. 2002, in Hot Star Workshop III, The Earliest Stages of Massive Star Birth, ed. P. A. Crowther, ISBN: 1-58381-107-9 (San Francisco: Astronomical Society of the Pacific), ASP Conf. Proc., 267, 95 In the text  
• Kennicutt, R. C., Bresolin, F., French, H., & Martin, P. 2000, ApJ, 537, 589 In the text NASA ADS  
• Kudritzki, R. P., Yorke, H. W., & Frisch, H. 1988, Research supported by the Swiss Academy of Sciences and DFG. Sauverny, Switzerland, Observatoire de Genève, 458 p., 89, 50499 In the text  
• Kurucz, R. 1991, in Stellar Atmospheres: Beyond Classical Models, ed. L. Crivellari, I. Hubeny, & D. Hummer, NATO ASI Series C, 341, 441 In the text  
• Kurucz, R. 1994, Solar abundance model atmospheres for 0, 1, 2, 4, 8 km/s. Kurucz CD-ROM No. 19. Cambridge, Mass.: Smithsonian Astrophysical Observatory, 1994, 19 In the text  
• Lanz, T., & Hubeny, I. 2003a, ApJS, 146, 417 In the text NASA ADS  
• Lanz, T., & Hubeny, I. 2003b, ApJS, 147, 225, erratum In the text NASA ADS  
• Lumsden, S. L., Puxley, P. J., Hoare, M. G., Moore, T. J. T., & Ridge, N. A. 2003, MNRAS, 340, 799 In the text NASA ADS  
• Macfarlane, J. J., Cohen, D. H., & Wang, P. 1994, ApJ, 437, 351 In the text NASA ADS  
• Martín-Hernández, N., Peeters, E., Damour, F., et al. 2002a, A&A, 381, 606 In the text NASA ADS  
• Martín-Hernández, N. L. 2002, Ph.D. Thesis, Rijksuniversiteit Groningen In the text  
• Martín-Hernández, N. L., Mokiem, M., de Koter, A., & Tielens, A. G. G. M. 2004, A&A, submited In the text  
• Martín-Hernández, N. L., Vermeij, R., Tielens, A. G. G. M., van der Hulst, J. M., & Peeters, E. 2002b, A&A, 389, 286 In the text NASA ADS  
• Martins, F., Schaerer, D., & Hillier, D. J. 2002, A&A, 382, 999 In the text NASA ADS  
• Mathis, J. S. 1982, ApJ, 261, 195 In the text NASA ADS  
• Mathis, J. S. 1985, ApJ, 291, 247 In the text NASA ADS  
• Mathys, G. 1988, A&AS, 76, 427 In the text NASA ADS  
• Mazzotta, P., Mazzitelli, G., Colafrancesco, S., & Vittorio, N. 1998, A&AS, 133, 403 In the text NASA ADS  
• Meynet, G., Maeder, A., Schaller, G., Schaerer, D., & Charbonnel, C. 1994, A&AS, 103, 97 In the text NASA ADS  
• Morisset, C. 2003, ApJ, in press [astro-ph/0310275] In the text  
• Morisset, C., & Péquignot, D. 1996, A&A, 312, 135 In the text NASA ADS  
• Morisset, C., Schaerer, D., Martín-Hernández, N. L., et al. 2002, A&A, 386, 558 In the text NASA ADS  
• Oey, M. S., Dopita, M. A., Shields, J. C., & Smith, R. C. 2000, ApJS, 128, 511 In the text NASA ADS  
• Oey, M. S., Shields, J. C., Dopita, M. A., & Smith, R. C. 2002, in Rev. Mex. Astron. Astrofis. Conf. Ser., 77 In the text  
• Okamoto, Y. K., Kataza, H., Yamashita, T., Miyata, T., & Onaka, T. 2001, ApJ, 553, 254 In the text NASA ADS  
• Okamoto, Y. K., Kataza, H., Yamashita, T., et al. 2003, ApJ, 584, 368 In the text NASA ADS  
• Péquignot, D., Ferland, G., Netzer, H., et al. 2001, in Spectroscopic Challenges of Photoionized Plasmas, ASP Conf. Ser., 247, 533 In the text  
• Pauldrach, A. W. A., Hoffmann, T. L., & Lennon, M. 2001, A&A, 375, 161 In the text NASA ADS  
• Peeters, E., Martín-Hernández, N. E., Damour, F., et al. 2002, A&A, 381, 571 In the text NASA ADS  
• Rubin, R. H., Simpson, J. P., Erickson, E. F., & Haas, M. R. 1988, ApJ, 327, 377 In the text NASA ADS  
• Savin, D. W., & Laming, J. M. 2002, ApJ, 566, 1166 In the text NASA ADS  
• Schaerer, D. 2000, in Stars, Gas and Dust in Galaxies: Exploring the Links, ASP Conf. Ser., 221, 99 In the text  
• Schaerer, D., & de Koter, A. 1997, A&A, 322, 598 In the text NASA ADS  
• Schaerer, D., & Stasinska, G. 1999, in ESA SP-427: The Universe as Seen by ISO, 751 In the text  
• Schaller, G., Schaerer, D., Meynet, G., & Maeder, A. 1992, A&AS, 96, 269 In the text NASA ADS  
• Schmutz, W. 1991, in NATO ASIC Proc. 341: Stellar Atmospheres - Beyond Classical Models, 191 In the text  
• Simpson, J. P., Colgan, S. W. J., Rubin, R. H., Erickson, E. F., & Haas, M. R. 1995, ApJ, 444, 721 In the text NASA ADS  
• Skillman, E. D. 1989, ApJ, 347, 883 In the text NASA ADS  
• Smith, L. J., Norris, R. P. F., & Crowther, P. A. 2002, MNRAS, 337, 1309 In the text NASA ADS  
• Stasinska, G., & Schaerer, D. 1997, A&A, 322, 615 In the text NASA ADS  
• Stasinska, G., Schaerer, D., & Zeipen, C. 2002, in A Massive Star Odyssey, from Main Sequence to Supernovae, ed. K. van der Hucht, A. Herrero, & C. Esteban, IAU Symp., 212 In the text  
• Takahashi, H., Matsuhara, H., Watarai, H., & Matsumoto, T. 2000, ApJ, 541, 779 In the text NASA ADS  
• Vilchez, J. M., & Pagel, B. E. J. 1988, MNRAS, 231, 257 In the text NASA ADS  
• Watson, A. M., & Hanson, M. M. 1997, ApJ, 490, L165 In the text NASA ADS