text:
# The Dependence of the Galaxy Luminosity Function on Large-Scale Environment

## 1 Introduction

In a hierarchical cosmogony like the cold dark matter (CDM) model, galaxies are assumed to form in dark matter haloes (e.g. White & Rees 1978). A generic prediction of such cosmogony is that the properties of the galaxy population must be closely related to that of the halo population. If the cosmological density field is Gaussian and if the power of density perturbations extends to large scales, as is the case in the current CDM cosmogony, halo properties are expected to be correlated to some degree with the large-scale structure, and so the properties of the galaxy population are also expected to change with large-scale environment. With high-resolution numerical simulations much has been learned about the halo population. It turns out that a strong correlation with the large-scale environment exists only for halo mass; all other important halo properties are at most weakly correlated with the large-scale environment at any given time (e.g. Lemson & Kauffmann, 1999). This result, combined with our current model of galaxy formation, implies that the environmental dependence of the galaxy population is mainly due to the change of the halo mass function with large-scale environment. We can therefore hope to understand the environmental dependence of the galaxy population by understanding how galaxies occupy dark haloes of different masses. This is indeed the approach taken by the so-called halo occupation models, in which the number of galaxies per halo is assumed to depend on halo mass only (e.g. Jing, Mo & Börner 1998; Peacock & Smith 2000; Seljak 2000; Scoccimarro et al. 2001; Benson 2001; Bullock, Wechsler & Sommerville 2002; Berlind & Weinberg 2002; Kang et al. 2002; Zheng et al. 2002; Yang, Mo & van den Bosch 2003a; van den Bosch, Yang & Mo 2003; Magliocchetti & Porciani 2003; Scranton 2002; Zehavi et al. 2003; Yan, Madgwick & White 2003; Kravtsov et al. 2003). Note that this assumption is non-trivial, as it implies that galaxy formation is largely a local process in individual dark matter haloes. This assumption was also made in many of the semi-analytic models based on Monte-Carlo halo merging trees (e.g. Kauffmann, White & Guiderdoni 1993; Cole et al. 1994; Sommerville et al. 1999), because in such models the merging histories are statistically the same for all haloes of the same mass at a given time Berlind et al. (2003) found that this assumption is consistent with their numerical simulations.

To test the validity of this assumption, we use the halo occupation model recently developed by Yang et al. (2003a), which is based on the conditional luminosity function (hereafter CLF) of galaxies in dark haloes of given mass, to predict how the galaxy luminosity function (hereafter LF) changes with large-scale environment. In this model, the change of the LF with large-scale environment is entirely due to the change of the halo mass function with large-scale environment. Since the conditional mass function of dark matter haloes can be accurately obtained from cosmological $`N`$-body simulation, the CLF model can be used to make accurate predictions for how the galaxy LF changes with large-scale density field. We also show that, when averaged over a large volume, the predicted galaxy luminosity density is closely correlated with the underlying mass density. Therefore, the luminosity density can be used to rank the mass density, and the model predictions presented here can be tested using large galaxy redshift surveys, such as the two-degree Field Galaxy Redshift Survey (2dFGRS; Colless et al. 2001) and the Sloan Digital Sky Survey (SDSS; York et al. 2000).

The paper is arranged as follows. In Section 2 we review the CLF model and describe how it can be used to calculate the LF of galaxies in different environments. Section 3 presents the model predictions for how the galaxy LF changes with large-scale density field. In Section 4 we use our model to study the correlation between the galaxy luminosity density and the underlying mass density, and discuss how the predictions of our CLF formalism can be tested with observations. We summarize our results in Section 5.

## 2 From the conditional luminosity function to galaxy luminosity function

In Yang, Mo & van den Bosch (2003a, hereafter Paper I), we developed a formalism, based on the conditional luminosity function $`\mathrm{\Phi }(L|M)`$, to link the distribution of galaxies to that of dark matter haloes. We introduced a parameterized form for $`\mathrm{\Phi }(L|M)`$ which we constrained using the LF and the correlation lengths as function of luminosity. In van den Bosch, Yang & Mo (2003a, hereafter Paper II), we extended our model by constructing separate CLFs for the early- and late-type galaxies. The CLF formalism developed in these two papers has subsequently been used to investigate galaxy clustering as function of luminosity and type (Yang et al. 2003b), to put constraints on cosmological parameters (van den Bosch, Mo & Yang 2003b), to constrain redshift-evolution in the halo model (Yan, Madgwick & White 2003), and to characterize the population of satellite galaxies (van den Bosch et al. 2004). In this paper, we use the CLF formalism to predict how the luminosity function of galaxies depends on large-scale environment. For completeness, we briefly summarize in the following the main ingredients of the CLF formalism, and we refer the reader to papers I and II for more details.

The conditional luminosity function $`\mathrm{\Phi }(L|M)\mathrm{d}L`$ gives the average number of galaxies with luminosities in the range $`L\pm \mathrm{d}L/2`$ that reside in haloes of mass $`M`$. It is parameterized by a Schechter function:

$$\mathrm{\Phi }(L|M)\mathrm{d}L=\frac{\stackrel{~}{\mathrm{\Phi }}^{}}{\stackrel{~}{L}^{}}\left(\frac{L}{\stackrel{~}{L}^{}}\right)^{\stackrel{~}{\alpha }}\mathrm{exp}(L/\stackrel{~}{L}^{})\mathrm{d}L,$$
(1)
where $`\stackrel{~}{L}^{}=\stackrel{~}{L}^{}(M)`$, $`\stackrel{~}{\alpha }=\stackrel{~}{\alpha }(M)`$ and $`\stackrel{~}{\mathrm{\Phi }}^{}=\stackrel{~}{\mathrm{\Phi }}^{}(M)`$ are all functions of halo mass $`M`$. Following Papers I and II, we write the average total mass-to-light ratio of a halo with mass $`M`$ as

$$\frac{M}{L}(M)=\frac{1}{2}\left(\frac{M}{L}\right)_0\left[\left(\frac{M}{M_1}\right)^{\gamma _1}+\left(\frac{M}{M_1}\right)^{\gamma _2}\right],$$
(2)
which has four free parameters: a characteristic mass $`M_1`$, for which the mass-to-light ratio is equal to $`(M/L)_0`$, and two slopes, $`\gamma _1`$ and $`\gamma _2`$, that specify the behavior of $`M/L`$ at the low and high mass ends, respectively. A similar parameterization is adopted for the characteristic luminosity $`\stackrel{~}{L}^{}(M)`$:

$$\frac{M}{\stackrel{~}{L}^{}(M)}=\frac{1}{2}\left(\frac{M}{L}\right)_0f(\stackrel{~}{\alpha })\left[\left(\frac{M}{M_1}\right)^{\gamma _1}+\left(\frac{M}{M_2}\right)^{\gamma _3}\right],$$
(3)
with

$$f(\stackrel{~}{\alpha })=\frac{\mathrm{\Gamma }(\stackrel{~}{\alpha }+2)}{\mathrm{\Gamma }(\stackrel{~}{\alpha }+1,1)}.$$
(4)
Here $`\mathrm{\Gamma }(x)`$ is the Gamma function and $`\mathrm{\Gamma }(a,x)`$ the incomplete Gamma function. This parameterization has two additional free parameters: a characteristic mass $`M_2`$ and a power-law slope $`\gamma _3`$. For $`\stackrel{~}{\alpha }(M)`$ we adopt a simple linear function of $`\mathrm{log}(M)`$,

$$\stackrel{~}{\alpha }(M)=\alpha _{15}+\eta \mathrm{log}(M_{15}),$$
(5)
with $`M_{15}`$ the halo mass in units of $`10^{15}h^1\mathrm{M}_{}`$, $`\alpha _{15}=\stackrel{~}{\alpha }(M_{15}=1)`$, and $`\eta `$ describes the change of the faint-end slope $`\stackrel{~}{\alpha }`$ with halo mass. Note that once $`\stackrel{~}{\alpha }(M)`$ and $`\stackrel{~}{L}^{}(M)`$ are given, the normalization of the conditional LF, $`\stackrel{~}{\mathrm{\Phi }}^{}(M)`$, is obtained through equations (1) and (2), using the fact that the total (average) luminosity in a halo of mass $`M`$ is

$$L(M)=_0^{\mathrm{}}\mathrm{\Phi }(L|M)LdL=\stackrel{~}{\mathrm{\Phi }}^{}\stackrel{~}{L}^{}\mathrm{\Gamma }(\stackrel{~}{\alpha }+2).$$
(6)
Finally, we introduce the mass scale $`M_{\mathrm{min}}`$ below which we set the CLF to zero; i.e., we assume that no stars form inside haloes with $`M<M_{\mathrm{min}}`$. Motivated by reionization considerations (see Paper I for details) we adopt $`M_{\mathrm{min}}=10^9h^1\mathrm{M}_{}`$ throughout.

In order to split the galaxy population in early and late types, we follow Paper II and introduce the function $`f_{\mathrm{late}}(L,M)`$, which specifies the fraction of galaxies with luminosity $`L`$ in haloes of mass $`M`$ that are late-type. The CLFs of late- and early-type galaxies are then given by

$$\mathrm{\Phi }_{\mathrm{late}}(L|M)\mathrm{d}L=f_{\mathrm{late}}(L,M)\mathrm{\Phi }(L|M)\mathrm{d}L$$
(7)
and

$$\mathrm{\Phi }_{\mathrm{early}}(L|M)\mathrm{d}L=\left[1f_{\mathrm{late}}(L,M)\right]\mathrm{\Phi }(L|M)\mathrm{d}L.$$
(8)
As with the CLF for the entire population of galaxies, $`\mathrm{\Phi }_{\mathrm{late}}(L|M)`$ and $`\mathrm{\Phi }_{\mathrm{early}}(L|M)`$ are constrained by 2dFGRS measurements of the LFs and the correlation lengths as function of luminosity. We assume that $`f_{\mathrm{late}}(L,M)`$ has a quasi-separable form

$$f_{\mathrm{late}}(L,M)=g(L)h(M)q(L,M).$$
(9)
Here

$$q(L,M)=\{\begin{array}{ccc}1\hfill & \text{if }g(L)h(M)1\hfill & \\ \frac{1}{g(L)h(M)}\hfill & \text{if }g(L)h(M)>1\hfill & \end{array}$$
(10)
is to ensure that $`f_{\mathrm{late}}(L,M)1`$. We adopt

$$g(L)=\frac{\widehat{\mathrm{\Phi }}_{\mathrm{late}}(L)}{\widehat{\mathrm{\Phi }}(L)}\frac{_0^{\mathrm{}}\mathrm{\Phi }(L|M)n(M)dM}{_0^{\mathrm{}}\mathrm{\Phi }(L|M)h(M)n(M)dM}$$
(11)
where $`n(M)`$ is the halo mass function (Sheth & Tormen 1999; Sheth, Mo & Tormen 2001), $`\widehat{\mathrm{\Phi }}_{\mathrm{late}}(L)`$ and $`\widehat{\mathrm{\Phi }}(L)`$ correspond to the observed LFs of the late-type and entire galaxy samples, respectively, and

$$h(M)=\mathrm{max}(0,\mathrm{min}[1,\left(\frac{\mathrm{log}(M/M_a)}{\mathrm{log}(M_b/M_a)}\right)])$$
(12)
with $`M_a`$ and $`M_b`$ two additional free parameters, defined as the masses at which $`h(M)`$ takes on the values $`0`$ and $`1`$, respectively. As shown in Paper II, this parameterization allows the population of galaxies to be split in early- and late-types such that their respective LFs and clustering properties are well fitted.

In Papers I and II we presented a number of different CLFs for different cosmologies and different assumptions regarding the free parameters. In what follows we focus on the flat $`\mathrm{\Lambda }`$CDM cosmology with $`\mathrm{\Omega }_m=0.3`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$ and $`h=H_0/(100\mathrm{km}\mathrm{s}^1\mathrm{Mpc}^1)=0.7`$ and with initial density fluctuations described by a scale-invariant power spectrum with normalization $`\sigma _8=0.9`$. These cosmological parameters are in good agreement with a wide range of observations, including the recent WMAP results (Spergel et al. 2003), and in what follows we refer to it as the “concordance” cosmology. Finally, we adopt the CLF with the following parameters: $`M_1=10^{10.94}h^1\mathrm{M}_{}`$, $`M_2=10^{12.04}h^1\mathrm{M}_{}`$, $`M_a=10^{17.26}h^1\mathrm{M}_{}`$, $`M_b=10^{10.86}h^1\mathrm{M}_{}`$, $`(M/L)_0=124h(\mathrm{M}/\mathrm{L})_{}`$, $`\gamma _1=2.02`$, $`\gamma _2=0.30`$, $`\gamma _3=0.72`$, $`\eta =0.22`$ and $`\alpha _{15}=1.10`$. As shown in paper II, this model (referred to as model D) yields excellent fits to the observed LFs and the observed correlation lengths as function of both luminosity and type. We emphasize, however, that our results are not sensitive to uncertainties in the CLF; had we chosen models A, B or C in paper II, instead of model D, the results presented below would have been virtually identical.

Since $`\mathrm{\Phi }(L|M)`$ is the average number of galaxies per unit luminosity in a halo with mass $`M`$, one can obtain the galaxy LF, $`\mathrm{\Phi }(L)`$, from the halo mass function, $`n(M)`$, according to

$$\mathrm{\Phi }(L)=\mathrm{\Phi }(L|M)n(M)dM.$$
(13)
Fig.1 shows how haloes of different masses contribute to the total luminosity function. Note that the shape of the CLF changes significantly with halo mass, and that the faint end of the LF is dominated by galaxies hosted by low-mass haloes.

If we make the assumption that the CLF is statistically independent of the large scale environment, the galaxy LF in a region of overdensity $`\delta `$ follows from

$$\mathrm{\Phi }(L|\delta )=\mathrm{\Phi }(L|M)n(M|\delta )dM.$$
(14)
Here $`n(M|\delta )`$ is the conditional mass function of dark matter haloes, which gives the number density of haloes as a function of halo mass in an environment with average mass overdensity $`\delta [\overline{}]/\overline{}`$ (with $``$ the total mass in volume $`V`$, and $`\overline{}`$ the mean mass in all volumes $`V`$). Thus, in the CLF formalism, the $`\delta `$-dependence of the galaxy LF enters only though the conditional mass function $`n(M|\delta )`$.

Fig. 2 shows the conditional mass functions for haloes in several representative environments. These functions are derived from a high-resolution $`N`$-body simulation which follows the motions of $`512^3`$ particles with a P3M code in a $`100h^1\mathrm{Mpc}`$ box (see Jing & Suto 2002 for details), assuming the $`\mathrm{\Lambda }`$CDM ‘concordance’ cosmology specified above. Dark matter haloes are identified using the standard friends-of-friends algorithm with a linking length of $`0.2`$ times the mean inter-particle separation. As one can see, the shape of the conditional mass function is quite independent of $`\delta `$ at the low mass end, while the shape at the high-mass end depends significantly on $`\delta `$. The break in the conditional mass function at the high-mass end occurs at a lower mass for lower $`\delta `$. This reflects the fact that the formation of massive haloes is suppressed in low density regions due to bias (e.g. Mo & White 1996; Gottloeber et al. 2003). For a given $`\delta `$, the break occurs at smaller mass if the volume used to define $`\delta `$ is smaller.

Ideally, if we have an accurate model for the conditional mass function, we can combine it with the conditional luminosity function to construct an analytical model for the $`\delta `$-dependence of the galaxy LF. Unfortunately, the model based on peak-background splitting (e.g. Cole & Kaiser 1989; Mo & White 1996) is not sufficiently accurate for our purpose, even if ellipsoidal collapse is taken into account (Sheth, Mo & Tormen 2001; Sheth & Tormen 2002). The reason for this is that, when the total mass contained in a volume becomes comparable to the mass of individual halos in consideration, as is the case for massive halos in low-density volumes, the peak-background splitting becomes inaccurate (e.g. Sheth & Tormen 2002). We therefore decided to use the simulated samples constructed in Yang et al. (2003b) for our purpose. These simulated samples were obtained by populating dark matter haloes in $`N`$-body simulations with galaxies according to the CLF model described above. Here we use the sample constructed from one of the $`100h^1\mathrm{Mpc}`$ simulations, which is complete down to a $`b_j`$-band luminosity of $`L_{b_j}10^8h^2\mathrm{L}_{}`$ (see Yang et al. 2003b for details).

## 3 Galaxy luminosity function in different environments

Fig. 3 shows the LFs of galaxies in regions of different mass overdensities $`\delta `$ (defined over spheres with radius $`R=8h^1\mathrm{Mpc}`$). As one can see, the shape of the LF at the faint end is quite independent of $`\delta `$. This is due to the fact that most faint galaxies reside in low mass haloes, for which the shape of the conditional mass function $`n(M|\delta )`$ is only weakly dependent on $`\delta `$ (cf. Fig. 2). The most pronounced difference between LFs in different mass density environments is the break at the bright end, which occurs at fainter luminosities in lower-density regions. This, in turn, is simply a reflection of the fact that most bright galaxies reside in massive haloes.

Our CLF model also predicts a pronounced difference for the environment dependence of the LFs of early- and late-type galaxies. As shown in Fig. 4, the shape of the LF of early-type galaxies depends rather strongly on $`\delta `$. For the late-type galaxies, however, no such pronounced shape-dependence is predicted. Once again, this behaviour is easy to understand from the conditional halo mass function and the CLF: many of the fainter early-types actually reside in clusters (see Paper II), whose abundance depends strongly on $`\delta `$. The majority of the faint late-type galaxies, on the other hand, reside in relatively low mass haloes, for which the shape of the conditional halo mass function is roughly independent of environment (cf. Fig.2).

To further quantify how the LF depends on the mean density of the environment, we fit each of the luminosity functions by a Schechter form. For most of our discussion, we only fit over the luminosity range $`L10^9h^2\mathrm{L}_{}`$ ($`M_{b_j}5\mathrm{log}h17.2`$). In real redshift surveys, such as the 2dFGRS, galaxies fainter than this are observed only within a relatively small local volume, making it difficult to study their large-scale environment dependencies. As shown in Figs. 3 and 4, the Schechter function provides in general a reasonable fit to the LF over a large range in luminosity. This is not a trivial result, as the shape of the CLF varies significantly with halo mass and the conditional mass function at the high-mass end changes significantly with $`\delta `$.

The left panel of Fig. 5 shows the characteristic luminosity $`L^{}`$ and the faint end slope $`\alpha `$ as functions of $`\delta `$. In all cases, the characteristic luminosity increases with $`\delta `$, reflecting the fact that massive haloes, which host bright galaxies, are preferentially located in high-density regions. The increase is only moderate: for the total sample it is about one magnitude from the lowest to the highest densities probed. The increase is more significant for early-type galaxies than for late-type galaxies. For late-type galaxies, the increase becomes insignificant for $`\delta 1`$.

For the total population the faint-end slope $`\alpha `$ is roughly a constant, with a value about $`1.1`$. The faint-end slope is also quite independent of $`\delta `$ for late-type galaxies, with $`\alpha 1.25`$. In contrast, the faint-end slope for early-type galaxies changes from $`\alpha 0.3`$ for low-density regions to $`\alpha 0.8`$ for $`\delta 10`$. As mentioned above, this is mainly due to the fact that a large number of relatively faint early types reside in massive haloes whose mass function depends strongly on $`\delta `$.

If we fit the luminosity function over the entire range with $`L10^{7.5}h^2\mathrm{L}_{}`$ ($`M_{b_j}5\mathrm{log}h13.5`$), the faint-end slope becomes steeper for low-density regions (cf. open and solid circles in left-panel of Fig. 5. Note that in this case, the Schechter function is no longer a good fit (see Fig. 3) for low-density regions; the faint-end slope is much steeper than that of the extrapolation from the luminosity function at the brighter end. This departure starts at $`L10^9h^2\mathrm{L}_{}`$ and is more significant for early-type galaxies.

We have also made similar analyses for spherical volumes with a radius of $`4h^1\mathrm{Mpc}`$, the results of which are shown in the right-hand panel of Fig. 5. A noticeable difference with using volumes of $`8h^1\mathrm{Mpc}`$ radius is that the characteristic luminosity increases faster with $`\delta `$. The reason is that the conditional mass function at the high-mass end reveals a stronger dependence on $`\delta `$ for smaller radius $`R`$ (see Fig. 2).

## 4 Predictions for observational tests

So far we have investigated how the galaxy LF depends on the mass density averaged over spherical volumes of $`8h^1\mathrm{Mpc}`$ radius. Unfortunately, it is not easy observationally to obtain accurate measures for the mass density field in the Universe. Therefore, in order to facilitate an observational test of the predictions presented above, we need to find a quantity that is (i) easy to derive from observations, and (ii) that can be used to rank the mass density.

When averaged within a large volume, the number density, $`\delta _\mathrm{g}`$ of galaxies (with luminosities above some value) and the luminosity density, $`\delta _L`$, are both expected to be correlated with the underlying mass density. Such relations can be derived directly from the simulated catalogs constructed with the CLF (Yang et al. 2003b). Fig. 6 shows the results. As one can see, there is a tight correlation between the luminosity overdensity $`\delta _L`$ and mass overdensity $`\delta `$. Here $`\delta _L`$ is calculated using all galaxies with $`LL_{\mathrm{min}}=10^7h^2\mathrm{L}_{}`$, but it is not sensitive to this lower limit as long as $`L_{\mathrm{min}}10^9h^2\mathrm{L}_{}`$. Thus the luminosity density can be used to rank the mass density. The use of galaxy number density is trickier, because the scatter is quite large for low-density volumes (for example, $`\delta _\mathrm{g}`$, here defined by all galaxies with $`L10^9h^2\mathrm{L}_{}`$, becomes zero at $`\delta <0.8`$ in the right panel of Fig. 6) and because it depends on the luminosity range used in defining the galaxy overdensity $`\delta _\mathrm{g}`$ (since the number density of galaxies is dominated by the faint ones).

In order to show that $`\delta _L`$ can be used to rank $`\delta `$, we first fit the mean $`\rho _L`$-$`\rho `$ relation by the following function

$$y=Ae^x_{}x^\beta e^{x_{}/x},$$
(15)
where $`x=1+\delta `$, $`y=1+\delta _L`$, and $`A`$, $`x_{}`$ and $`\beta `$ are fitting parameters. For high-density cells ($`xx_{}`$), $`yx^\beta `$, while for $`xx_{}`$, $`y`$ decreases exponentially with decreasing $`x`$. For $`R=8h^1\mathrm{Mpc}`$, we obtain $`A=1.0203`$, $`x_{}=0.3053`$ and $`\beta =0.8035`$. We then use this mean relation to convert the distribution function of the mass density into a distribution function for the luminosity density and compare it with the luminosity-density distribution function derived directly from the simulation. The result is shown in Fig. 7. Clearly, the luminosity-density distribution function is recovered quite accurately with the mean relation (15), suggesting that the luminosity density can be used to represent the mass density.

Fig. 8 shows how $`L^{}`$ and $`\alpha `$ change with the luminosity density. The results are similar to those shown in Fig. 5, except that the range of $`\delta _L`$ is stretched at the low-density end relative to $`\delta `$, because of the non-linear relationship between $`\delta _L`$ and $`\delta `$. As comparison, we also show $`L^{}`$ and $`\alpha `$ as functions of $`\delta _g`$, where $`\delta _g`$ is based on the number of galaxies with $`L10^9h^2\mathrm{L}_{}`$ (Fig. 9). With such a luminosity limit, the results based on $`\delta _g`$ are very close to those based on $`\delta _L`$.

## 5 Discussion

We have shown that the halo occupation model, which assumes that the luminosity distribution in a dark halo only depends on its mass, makes specific predictions about how the galaxy LF changes with large-scale environment. The model predicts that the LFs for relatively bright galaxies can be fit reasonably well by a Schechter function, independent of environment. In all cases the characteristic luminosity $`L^{}`$ increases with local mass density. For late-type galaxies, the faint-end slope is quite independent of large-scale environment while for early-types the value of $`\alpha `$ changes from $`0`$ in low density regions to $`1`$ in high-density regions. The predicted LF for low-density regions shows significant departure from the Schechter form at $`L10^9h^2\mathrm{L}_{}`$, and this departure is more significant for early-type galaxies. All these predictions have simple explanations based on the change of halo mass function with large-scale environment, and the change of halo occupation with halo mass.

We have also shown that the luminosity density is tightly correlated with the mass density on large scales. Thus, our predictions can be checked using large galaxy redshift surveys. At the moment, the study of the dependence of galaxy LF on local environment has mainly focused on comparing the LF of cluster galaxies with that of field galaxies (see e.g. De Propris et al. 2003 and references therein). More related to the model predictions presented here are the analyses made by Hütsi et al. (2002) and Bromley et al. (1998). Hütsi et al. estimated galaxy LFs in spherical volumes with a radius $`10h^1\mathrm{Mpc}`$ (in redshift space) within which the overdensities of galaxy number are respectively $`\delta _g<0`$, $`0<\delta _g<1`$, and $`\delta _g>1`$, and found that the faint-end slope $`\alpha `$ is about $`1.1`$ for all the three cases, but that the characteristic luminosity brightens by about 0.4 magnitude from the lowest-density sample to the highest-density sample. These results are consistent with our model prediction. Bromley et al. found that the faint-end slope of the LF of early-type galaxies (defined according to spectral type) depends strongly on local density, with $`\alpha `$ increasing from $`0.4`$ in high density regions to $`0.2`$ in low density regions. This trend is also consistent with our model prediction. Unfortunately, current observational results are not yet sufficiently accurate to give a stringent constraint on the model. The situation will soon change. With the use of large redshift surveys of galaxies, such as the 2dFGRS and SDSS 111During the final stages of this project, Hoyle et al (2003) published a paper based on SDSS data that adresses the environment-dependence of the LF. Their results are qualitatively in good agreement with our predictions, though we delay a more detailed comparison to a future paper., one can estimate to high accuracy the galaxy luminosity functions in regions of different luminosity densities or galaxy number densities (Croton et al. in preperation), facilitating an accurate comparison between model and observations

If our model prediction agrees in detail with observational data, it will add strong support to the assumptions made in the theory: the ‘concordance’ cosmology (which determines the properties of the halo population) and that galaxy properties only depend on halo mass. If we assume that the halo population in the real universe is not very different from that predicted by the ‘concordance’ cosmology, any significant difference between our model prediction and observation would have important implications for the formation of galaxies in the cosmic density field. For example, the prediction that the faint-end slope of the LF should not depend significantly on large-scale environment is due to the fact that most of the faint galaxies in our CLF model are in low-mass haloes. A significant discrepancy between model prediction and observation would therefore require a CLF model in which many faint galaxies are in massive haloes. This change of CLF is, however, strongly constrained by the observed weak clustering of faint galaxies. Thus, any significant discrepancy would strongly indicate the assumption that the galaxy properties only depend on halo mass is incorrect, invalidating the whole halo occupation approach. The other prediction of the model, that the characteristic luminosity increases mildly with local density, is due to the fact that in our CLF model bright galaxies are preferentially found in massive halos. Here any discrepancy with observational data would therefore require a change in the fraction of bright galaxies in massive haloes. Such change will lead to a change in the correlation length as a function of luminosity, and is subject to sringent observational constraints. As mentioned above, the other CLF models we have obtained give very similar results, and so any significant discrepancy would mean that either galaxy properties depend significantly on other properties of halos in addition to the mass, or our parameterization of the CLF is not sufficiently general. As discussed in Paper I, our parameterization is motivated by physical considerations, and is quite general as long as the CLF is assumed to be a smooth function of halo mass. Thus, should a sigificant change in the CLF be required by the observation, it would have a significant impact on our understanding of galaxy formation in dark matter haloes.

As far as an accurate comparison with observation is concerned, there are several effects to be taken into account. The first is redshift distortion. Since bright galaxies in massive clusters are expected to have large velocity dispersion, some of these galaxies may appear in ‘low-density’ regions in redshift space, weakening the dependence on $`\delta `$. This effect was seen in our experiment where $`\delta `$ was defined in redshift space. The second is concerned with the definition of overdensity. In a magnitude-limited sample, volumes at large distances from us contain only bright galaxies, and so the estimate of the overdensity in such a cell has to rely on a small number of bright galaxies. Since the shape of the luminosity function changes with $`\delta `$, the correction of the incompleteness due to the magnitude limit in individual volumes must take this effect into account. Finally, since the Schechter function is not a perfect model, the fitting parameters depend on the luminosity range used for the fit. All these effects can be taken into account properly with the use of mock catalogs by applying the same analyses to both the simulated and observed samples. We plan to present such a detailed comparison with observational data in a forthcoming paper.

It is also interesting to compare the predictions presented here, based on the CLF, with those based on semi-analytical models for galaxy formation (e.g. Kauffmann et al. 1999; Somerville & Primack 1999; Benson et al. 2002; Mathis et al. 2002). Mathis & White (2002) have shown that the characteristic luminosity of the galaxy luminosity function in their semi-analytical model increases progressively from low-density to high-density environments while the faint-end slope becomes slightly shallower. These results are in qualitative agreement with our predictions presented here. As already emphasized in paper II, this suggests that the halo occupation statistics as described by our CLF nicely fit within the standard framework of galaxy formation.

## Acknowledgement

We are grateful to Darren Croton, Glennys Farrar, and Guinevere Kauffmann for useful discussions. Numerical simulations used in this paper were carried out at the Astronomical Data Analysis Center (ADAC) of the National Astronomical Observatory, Japan. YPJ is supported in part by NKBRSF (G19990754) and by NSFC.


format: text

text:
# The Ofpe/WN9 stars in M 33 1footnote 11footnote 1Based on observations with the Hubble Space Telescope which is operated by AURA, Inc. under NASA contract NAS5-26555

## 1 Introduction

There is a general consensus about the Wolf-Rayet stars being formed from O stars through high mass-loss phases, but many specific questions remain still open. The Ofpe/WN9 stars (“slash stars”) class, introduced by Walborn (1977, 1982a), shows properties intermediate between those of Of and WN stars: the emission is stronger than in the most extreme Of stars, but HeII photospheric absorption indicates a less dense wind than seen in any Wolf-Rayet stars (see Bohannan & Walborn 1989). The implication is that Ofpe/WN9 stars are transition objects between the Of stars and the Wolf-Rayet in the “Conti scenario” (Conti 1976; Conti & Bohannan 1989). Based on the reduced surface H fraction observed in these stars, Pasquali et al (1997) suggest instead the sequence O $``$ Of $``$ H-rich WNL $``$ Ofpe/WN9, for initial masses less than 100 M. One of the prototypes of the Ofpe/WN9 stars, R127 (HDE 269858) underwent an “S Dor like” eruption in 1982 (Walborn 1982b; Stahl et al. 1983) revealing Luminous Blue Variable (LBV) characteristics. Other evidence suggests that at least some of the Ofpe/WN9 stars have progenitors that are considerably lower in mass than the WNLs; see Massey, DeGioia-Eastwood, & Waterhouse (2001). LBV’s irregular episodes of extreme mass-loss may play an important role in the formation of W-R stars by substantial shedding of the progenitor’s outer layers (e.g. Maeder & Conti 1994).

Ten “slash stars” are known in the LMC (Bohannan & Walborn 1989). Seven LMC Ofpe/WN9 stars were investigated in the optical (Nota et al. 1996) and in the UV with Faint Object Spectrograph (FOS) spectra (Pasquali et al. 1997), and with IR spectroscopy by Morris et al (1996). The spectral morphology confirms the transition nature of these objects and their resemblance to LBVs. Crowther et al. (1977) also analysed optical spectra of LMC Ofpe/WN9 stars and re-classified them as WNL. Krabbe et al. (1991) observed a number of HeI emission line stars in the Galactic Centre, whose HeI line ratios and line-to-continuum ratio are consistent with Ofpe/WN9 stars. An Ofpe/WN9 star in NGC 300 was reported by Bresolin et al. (2002).

With follow-up, ground-based (KPNO, MMT) spectroscopy of about 400 UV-brightest stellar sources found in the Ultraviolet Imaging Telescope (UIT) far- and near-UV images of M 33, we discovered six Ofpe/WN9 stars in this galaxy (Massey et al. 1996, hereafter MBHS). The “slash stars” are UV bright. In our study of the UIT UV-brightest sources in M 33, these stars are relatively well separated from other stars in the $`UB`$ vs. $`FUVNUV`$ color-color diagram. They lie near the LBV candidates. Their colors are consistent with the fact that the Balmer jump is in emission in Ofpe/WN9 stars, which requires both high UV flux and an extended atmosphere, consistent with the strong P Cygni profiles (e.g. HeI $`\lambda `$ 4471) characteristic of the Ofpe/WN9 class.

We present Space Telescope Imaging Spectrograph  (STIS) ultraviolet spectra of the six known Ofpe/WN9 stars in M 33. In section 2 details are given about the observations and the reduction, in section 3 the STIS UV spectra are analysed with non-LTE, line blanketed, hydrodynamical, spherical model atmospheres to derive physical stellar parameters. The results are interpreted in terms of stellar evolution by comparison with evolutionary model predictions (Section 4).

## 2 The programs stars. Data and reduction

### 2.1 The program stars

In Bianchi’s HST program GO8207, we obtained UV spectra of the six Ofpe/WN9 stars known in M33. The spectra are shown in Figure 1. The objects were classified as Ofpe/WN9 from ground-based (KPNO, MMT) classification spectra of the UV-brightest stellar sources measured in the UIT far-UV and near-UV images (MBHS). The optical spectra, taken several years earlier (1993 - 1995), are shown in Figure 2 for the purpose of discussion in the following sections. A seventh object (M33-UIT339) was included in our program because of its photometric properties (H$`\alpha `$ emission-line source ). An additional object, M33-UIT340, was included serendipitously in the STIS long slit during the M33-UIT339 observations by chosing an appropriate orientation. Its spectrum has a UV flux level comparable to the primary target, but much stronger UV lines, in spite of the similar optical colors. Another serendipitous star was observed due to failure of the acquisition on a repeated observation of M33-UIT349.

Data from our previous ground-based photometry are compiled in Table 1. In the table we give the star name, coordinates and V-mag, from the UIT sources list of MBHS. Also, for the position of every star we computed the de-projected galactocentric distance, using PA = 22 for the position angle of the semi-major axis of M 33, an inclination of the disk i=54 and a systemic velocity for M33 of Vsys=-180 km s-1 (from Warner et al. 1973). Warner et al. (1973) also provide a model for the rotation of the M33 disk from HI measurements. We computed the expected disk velocity, projected along the line of sight, at the position of each star. Because the targets are likely to lie in the M 33 disk, we used this velocity to shift the spectra to a rest-wavelength frame at the star. The position of the stellar lines indicates that this velocity correction is appropriate in all cases, although the difference from just using the M33 systemic velocity is not really appreciable at the resolution of our spectra. Both the galactocentric distance and the recessional velocity are also given in Table 1. Finally, the STIS datasets are listed.

### 2.2 STIS spectroscopy and data reduction

The observations were completed over a long time span, because the HST entered a safe mode event just after our program (GO8207) began, and most of the observations were delayed because we missed the window for the desired orientation constraints. To assure homogeneous data quality and consistent flux calibration, we used the most recent on-the-fly calibration for all of the data after the program was completed.

For each star the same observing pattern was repeated: after target acquisition, a 24 min. long exposure was obtained with the STIS G230L grating (range 1570-3180 Å, scale 1.58Å/ pixel), followed by an 8 min. exposure in the same orbit with the grating G140L (range 1150-1730Å, scale 0.60Å/pixel). The whole subsequent orbit was devoted to a long (57 min) G140L exposure. The total exposure time was therefore 3300sec (G140L grating) and 1440sec (G230L grating) for every target. The STIS 0.2″ x 52 ″ slit was used. For target M33-UIT349, which is fainter than the others in the sample, we repeated this pattern twice, aiming at doubling the total time and achieve a S/N comparable to the other sample objects. The two sets of observations were taken about a week apart (Dec. 19 and Dec. 27, 2000). The spectra are significantly different on the two dates, both in flux level and wind lines, as can be seen in Fig 1. After examining the STIS CCD50 acquisition images, we discovered that the intended target (M33-UIT349) was actually observed on Dec. 27, while a nearby star was centered in the aperture on Dec. 19, during the fine centering stage, in spite the right target was exactly centered (on the same pixel) after coarse aquisition on both dates, as shown in Figure 3. We measured photometry on both dates, on the STIS CCD acquisition images: the two stars are equally bright in the CCD50 filter. Additionally, the target (M33-UIT349) appears slightly elongated in the acquisition images (Figure 3), but it is not clearly resolved into multiple sources. A cut of the long slit 2D spectra along the spatial direction, by summing several columns across the wavelength dispersion (Figure 4), shows a hint of a faint extension in M33-UIT349 spectrum (Dec. 27, datasets O5CO07010,20,30), although the two peaks, or asymmetric intensity profile, may just be an effect of the low countrate. The STIS MAMA spatial scale is 0.025”/pxl, and the width of the entire spectrum (see Figure 4) is about 10 pixels or .25”, with the two brightest components less than .1” apart. Therefore it appears that, although the object is not resolved into multiple stars in the STIS CCD, which has a factor of two lower resolution than the MAMA, there may be some \[marginal\] extra light into the STIS slit from other fainter sources. The spectra are centered in the y-direction (perpendicular to the dispersion) in the nominal position. Just based on the elongation seen on the finding chart, and the two apparent peaks in the spectrum intensity profile, we can speculate: could the object be a binary (multiple?) system, or is there a chance superposition of a faint source? A separation of $``$0.1” at the distance of M33 (840,000 pc, Magrini et al. 2000) would correspond to $``$ 0.4pc, extremely large for a bound binary. An additional source seem to be excluded by the quantitative modeling of the spectra.

We will call the serendipitous star which was in the slit on Dec. 19, M33-UIT349-B hereafter. In the acquisition image, it is just 1.33′′ east of the target, thus the coordinates given in Table 1 are relative to the ground-based coordinates of the source (+0.10 seconds in R.A.). Most likely, the ground-based coordinates of M33-UIT349 were an intermediate position between the two sources, whose separation of 1.3′′ is not resolved in the ground-based image, and is comparable to the coordinate accuracy. The two stars can uniquely be identified from their relative positions and our finding chart (Figure 3).

For every target star we examined the two G140L exposures taken in subsequent orbits, and found that they perfectly match in absolute flux and the strong line features reproduce very well - in spite of the very different S/N due to the different exposure time. Also, we verified that the G230L and G140L fluxes match in the wavelength region of overlap. With the exception of M33-UIT339 (see later), we then combined the two G140L exposures, taking into account the different exposure times, to obtain the maximum possible total S/N over the range. The STIS flux calibration is expected to be accurate for sources in the broad 2′′ wide slit to 2% (Bohlin 2000; Bohlin, Dickinson, & Calzetti 2001). For the 0.2′′ slit used in this program, photometric precision of one wavelength with respect to a distant wavelength is 4.5%rms (Bohlin & Hartig 1998). There are no known non-linearities in the response of the STIS MAMA detectors. However, additional uncertainties arise from the Poisson counting statistics and from possible errors in background subtraction for the faintest sources.

An additional target, M33-UIT339, was included in our program GO8207 as it appeared to be an H$`\alpha `$ emission-line source, however it was not subsequently classified as an Ofpe/WN9 star. With an appropriate orientation, the STIS long slit centered on M33-UIT339 included also M33-UIT340, separated by about 4$`\mathrm{}`$ . The pipeline automatic extractions for the two G140L spectra, taken in subsequent orbits, are centered on the two differemt stars respectively, therefore the extracted spectra from the archive should not be used. The automatic pipeline extraction has centered the extraction slit on the target in one case (dataset O5CO03020), and on the serendipitous star in the other case (dataset O5CO03030). Therefore we manually re-extracted the spectra of both stars in all three exposures (two G140L and one G230L) and coadded them separately. While there is no acquisition image for the second orbit (only reacquisition of previous guide stars was performed), we did obtaine a parallel WFPC2 image with the same filter (F170W) in each orbit, as well as with different filters. Our F170W WFPC2 images with 0.1$`\mathrm{}`$ resolution show that there is no shift between the two pointings down to the pixel level. The continuum level is the same over most of the G140L range, but the wind features are completely different for the two stars (Figure 1).

The spectra of the serendipitous stars are shown in Figure 1 for completeness, and for clarification of all aquisition and multiple-source issues to potential archive users. They will be analyzed separately.

Because the targets were selected from ground-based photometry and spectroscopy, we examined all the STIS acquisition images to check for multiplicity that might be revealed at the HST STIS imaging resolution but not in the ground-based data. All targets appear point-like.

## 3 Analysis of the UV spectra

### 3.1 UV line morphology

The short wavelength ($`<`$ 1750Å) portion of the spectra for all the program stars and the serendipitous stars is shown in Fig. 1. This wavelength region is richest in lines from the stellar wind (and the interstellar medium). The main wind lines are shown in detail in Figure 5, where their extent in velocity (largely varying across the sample) can be appreciated. The spectra are arranged from top to bottom of the figures in a progression of decreasing strength of the UV lines and decreasing wind velocity (width of the P Cygni absorption). Previous ground-based spectroscopy from MBHS also showed a range of line strengths within the sample. However, the progression is different from what the recent UV spectra indicates. The optical spectra from MBHS are shown in Figure 2. For example, M33-UIT008 was a weak-liner in the optical spectra, but displays some of the strongest UV features in the STIS spectra. Viceversa, M33-UIT349 has fairly strong optical lines (Figure 2), but not in the UV spectra. As another example, M33-UIT045 has weak lines both in the old optical spectra and in the UV spectra.

The UV line strength progression indicates a range of $`T_{\mathrm{eff}}`$(luminosities) and mass-loss rates among the sample. It is difficult to relate the progression to an evolutionary sequence, because the relative line strengths are different between the optical spectra ($``$ 1993 -1995) and the recent (2000-2001) UV spectra, suggesting variability. One stellar parameter which can plausibly vary on very short time-scales is the mass-loss rate, if these objects were single stars. Or are they evolved binaries? The quantitative spectral modeling described below is aimed at claryfying the evolutionary status, by measuring the stellar parameters.

### 3.2 NTLE modeling with WM-basic

We performed detailed quantitative modeling of the spectra with the WM-basic code (Pauldrach et al. 2001). The code solves the wind hydrodynamics equations for specified parameters and includes non-LTE metal line blanketing in spherical geometry. It also allows us to include the effects of shocks (from clumps) on the wind ionization in the radiative transfer calculations, which are extremely important for a correct modeling of hot massive stars (e.g. Bianchi & Garcia 2002, 2003). Below we discuss the determination of the physical parameters for each individual star and how the uncertanties were estimated. More discussion on the uncertainties can be found in Section 4. The parameters from the best fit models are listed in table 2. The process is quite complex, given the many physical parameters which affect the emerging spectrum. The analysis of the first object is thus described in detail as a sample case. The effects of the physical parameters on the emergent spectrum are similar for the other objects. Initial guesses for the terminal velocities were taken by measuring the maximum extent of the P Cygni profile absorptions, which largely varies across the sample as can be appreciated (on an expanded scale, and independent of any modeling) in Figure 5.

#### 3.2.1 M33-UIT236

The UV spectrum of this object indicates that it might be the hottest of the sample, and that it has a high mass-loss rate, because it has the strongest UV wind lines. We ran a grid of WM-basic models, with initial parameters covering a range of values for $`T_{\mathrm{eff}}`$, luminosity and mass-loss rate inferred from a wider grid computed for other purposes. We computed synthetic model spectra for both solar and LMC-type metallicities, and then we also varied individual abundances of the CNO elements. In Table 2 we give the resulting stellar parameters and their uncertainty range, and here we discuss how they were constrained.

Before we discuss the strong P Cygni profiles that primarily constrain the main stellar parameters, we note that all models of the initial grid computed for LMC-type metallicity (z=0.008) are “flatter” than the observed spectra in the some spectral regions which contain no single strong lines, but a multitude of fainter lines, mainly from Fe III and Fe IV transitions (e.g. from 1420 to 1500Å). In these regions, the emerging flux in the low (LMC-type) metallicity models is close to a pure continuum, while the observed spectrum definitely shows shallow but significant depressions. These can be reproduced better in models with solar metallicity (z=0.02). A precise determination of z is not possible. Therefore we adopted solar metallicity in the subsequent computation of the finer grid which includes the best fit adopted solutions, discussed hereafter. The result is not surprising given the metallicity gradient in M33 and that M33-UIT236 is relatively close to the galaxy center. The model spectra are sensitive to this effect only for large wind velocities.

The STIS spectrum of M33-UIT-236 (Figure 1) shows strong P Cygni profiles of both Si IV and N V, as well as a strong C IV $`\lambda \lambda `$ 1548,1551 doublet and C III multiplet around 1174Å. To reproduce the observed strength of both N V and Si IV consistently, we need to take into account the X-rays produced by shocks in the wind as a relevant ionization process. Because there is no prior information, to our knowledge, about soft X-ray emission from these objects, we tried a range of values that would be appropriate for early type luminous stars as observed e.g. in our Galaxy (e.g. Bianchi 1982a,b, Bianchi & Garcia 2002a,b). We found that for $`T_{\mathrm{eff}}`$= 32000 K or lower, $`L_x`$/$`L_{}`$$``$ -7 is appropriate to reproduce the strength and shape of the N V doublet and the other “hot” transitions relative to the lines from lower ionization potential stages such as Si IV . A lower X-ray luminosity of $`L_x`$/$`L_{}`$=-7.5 is found to be adequate at a higher $`T_{\mathrm{eff}}`$(34,000K). The Si IV $`\lambda \lambda `$1394,1403 doublet, as well as N IV $`\lambda `$ 1718, is not very sensitive to the effect of shocks, but is instead very sensitive to $`T_{\mathrm{eff}}`$(especially in the range 34000 - 30000 K) and to the mass-loss rate. HeII $`\lambda `$ 1640 is very sensitive to the mass-loss rate, in this temperature regime, thus cannot be used to constrain the helium abundance. The relative strength of the C IV to C III lines is again best reproduced by models between $`T_{\mathrm{eff}}`$= 32000 K and $`T_{\mathrm{eff}}`$= 34000 K. However, both of these C lines are always too strong for all parameter combinations which can reproduce the overall spectral features, in the models with solar abundances. We have therefore progressively decreased the C abundance (and changed correspondingly the abundance of N and O as expected in different evolutionary scenarios) and finally achieved a consistent fit of all the strongest wind lines with carbon 20 times less abundant than the solar value and nitrogen 10 times overabundant. These values correspond to the abundances predicted by evolutionary models for e.g. a star with initial mass of $``$30-50$`M_{}`$ reaching the WNL stage. The evolutionary models also predict oxygen to be underabundant by a factor of ten with respect to solar, and He/H up to 0.5 for evolved massive objects reaching such C and N abundances (e.g. Maeder 1987, Maeder & Meynet 1994, 2000, 2003). We thus adopted also \[O\]=0.1 the solar value, for consistency with our findings of the C and N abundances, although there are no strong oxigen lines in the range covered by our STIS spectra for stars in this temperature range to directly measure the oxigen abundance.

We varied He/H from 0.1 to 0.5 (by number) and in both cases we can achieve an acceptable fit of the overall spectrum, by slightly tuning the other parameters. Best fits are obtained with $`T_{\mathrm{eff}}`$= 32000 K if He/H = 0.1 and with $`T_{\mathrm{eff}}`$= 34000 K for He/H = 0.5. Changing the helium abundance does not only affect the strength of the He II $`\lambda `$ 1640, but also influences the relative strengths of other lines as well. Therefore we consider the actual helium abundance a source of uncertainty and adopt the above range for $`T_{\mathrm{eff}}`$.

Once the strength of all P Cygni lines in the spectrum was matched by the models and the acceptable ranges constrained, we fine-tuned the trough of the absorptions, especially the Si IV doublet which is partly resolved. We varied both the terminal velocity and the rotational velocity. While the terminal velocity significantly affects the separation of the two absorption components of the Si IV doublet, and its uncertainty range is easily determined (see Fig. 5 and Table 2), rotational velocities of V$`sini`$ = 60 km s-1and twice as high V$`sini`$ = 120 km s-1 were tested in the model calculations: the change only slightly affects the shape of the absorption troughs. This analysis is not as sensitive to this parameter (V$`sini`$) as the optical photospheric lines, since lines at UV wavelengths are mainly produced in the wind thus broad and asymmetric, and the STIS spectra resolution is too low for this purpose. The fit of the C III and Si IV lines, favours the high rotation. Again, a quantitative determination of this parameter is not possible; we can only infer from the comparison an indication of possible high rotation.

All the best fit models shown are computed with log g =3.2. The low gravity value was initially driven by the intention to compute models with reasonable mass/radius/luminosity values for this type of stars. We then tried higher gravities and obtained worse line fits. Note the good match of both the absolute flux values (a further confirmation that our choice of the Radius in the best fit models is adequate) and of the slope of the spectrum, independently confirming the temperature derived above by fitting the line strengths, and the extinction value (which is also consistent with the previous work of MBHS). Because the spectral slope depends on $`T_{\mathrm{eff}}`$(and reddening) and the absolute flux depends on R2 and $`T_{\mathrm{eff}}`$4, the overall match to both the line profiles and the absolute flux level confirms the consistency of the parameters derived from our analysis.

The derived C/N ratio (by mass) of 0.02 is significantly lower than the cosmic value (4.8) and in the range of the WNL type stars (see e.g. Massey 2003). The C/N value is consistent, according to the stellar evolution calculations of Maeder (1987), Maeder & Meynet (1994), Meynet & Maeder (2000, 2003) with the surface abundances of a massive star during the stage from B supergiant to WNL. The current mass, $`T_{\mathrm{eff}}`$ and luminosity derived from our spectral modeling for M33-UIT236 agree with the evolutionary models for an initial mass of about 50 $`M_{}`$.

#### 3.2.2 M33-UIT008

If we superimpose the spectrum of this star to the previous one (M33-UIT236), the two spectra, normalized to the respective continua, are basically identical except for much weaker nitrogen lines in M33-UIT008. In particular, the C IV and C III lines are identical, and so is the Si IV doublet, except for a clearer separation of the two absorption components, indicating slightly lower terminal velocity or mass-loss in M33-UIT008. He II $`\lambda `$ 1640 is also identical. However, the N IV and N V lines are both significantly weaker than in the previous star. We attempted to explain this difference by either a lower nitrogen abundance or by a lower $`T_{\mathrm{eff}}`$ or $`\dot{M}`$. In our grid of models for solar CNO abundances, the N IV $`\lambda `$ 1720 line becomes very weak (comparable to the one seen in the M33-UIT008 spectrum) when $`T_{\mathrm{eff}}`$ is cooler than 28000 K, for mass-loss rates of 3-5 10-6 M yr-1, but the strength of this line also decreases significantly with decreasing mass-loss rate. However, these models with lower $`T_{\mathrm{eff}}`$ produce too much Si IV and too little C III, respect to the observed spectrum, if we keep the other parameters similar to M33-UIT236, including the CNO “WNL-type” abundances. In short, a good fit of the other lines cannot be achieved at such low $`T_{\mathrm{eff}}`$.

The strength of the N V doublet significantly decreases in models with lower $`L_x`$/$`L_{}`$, however N IV is practically not affected by this change and the Si IV doublet is significantly enhanced. Therefore, varying this parameter does not produce a good fit either.

Because both N IV and N V transitions are weaker in UIT008 than in UIT236, but the spectra are otherwise identical, we computed models with the same parameters as the best fit of M33-UIT236, but lowered the nitrogen abundance from “WNL” values down to the solar value. This produces slightly lower (not lower enough) nitrogen lines but enhances the Si IV and C IV emissions enormously. Probably because we keep other parameters (such as mass-loss rate) the same, the opacities are redistributed. We also tried to decrease the mass-loss rate, and alter other parameters including different combinations of abundances, and the amount of shocks.

The best fit, in terms of the relative line strength and match to the flux below 1700Å, is achieved with the parameters given in Table 2. The results indicate that the main difference with respect to M33-UIT236 is a lower mass-loss rate.

We caution that Bianchi & Garcia (2002), Garcia & Bianchi (2003) encountered some difficulties in fitting weak N V lines in pop. I Milky Way stars with the WM-basic code, and better matches could be achieved in some cases using Hillier’s CMFGEN code, although the resulting stellar parameters were not changed (e.g. Bianchi et al. 2003a, b). Therefore, the derived abundances, in particular of N in this case, must be taken as indicative, as stressed in the previous section.

#### 3.2.3 M33-UIT104

The line spectrum suggests lower wind velocity and a lower stellar temperature than the previous stars examined, the “cooler” transitions (C III and N IV) being enhanced and the higher ions (N V and C IV) weaker. A similar effect on the line strength, to some extent, may be produced by a lower amount of X-rays from shocks. All these effects were quantified by the modeling (Table 2). A CIII line, visible as a small but detectable absorption in the emission of N V, is the most unambiguous indication of $`T_{\mathrm{eff}}`$ $``$32000 K according to our grid of models (see Bianchi & Garcia 2002), at least for solar abundances. Therefore, this line provides a strong constraint on $`T_{\mathrm{eff}}`$.

In this $`T_{\mathrm{eff}}`$ regime, the C III 1175Å line becomes much more sensitive to the C abundance than to the mass-loss rate (while C IV is saturated so rather insensitive to all parameters). This is very fortunate because it allowed us to constrain one parameter (the C abundance) relatively independently from the other ones in this case. Along the stellar evolution, the C-depletion develops concurrently with the N-enrichment (according to Maeder 1987), therefore we tried different abundance ratios for C and N, corresponding to a typical predicted evolutionary sequence. The best fit is shown in figure 6, and was achieved with \[C\]=0.1 x solar and \[N\]=2 x solar. In mass fraction, this means C/N=0.9.

Again the match to both the line strength and to the absolute flux level, and the continuum shape, is very good shortwards of 1700Å.

#### 3.2.4 M33-UIT003, M33-UIT349 and M33-UIT045

Continuing towards progressively weaker P Cygni profiles and lower terminal velocities, M33-UIT003, M33-UIT045 and M33-UIT349 (Dec.27) show a smooth decrease in the lines of Si IV and C IV, compared to the previous targets. Instead, the N V $`\lambda `$1240 doublet abruptly desappears in the spectra of M33-UIT045 and of M33-UIT349. In modeling the spectrum of M33-UIT045, it has been impossible to obtain a synthetic spectrum which matches the line features, and the absolute flux, with nitrogen abundance solar or higher. The best fit spectrum (shown in Figure 6) was obtained by lowering the nitrogen abundance by a factor of ten with respect to solar. A nitrogen underabundance in such evolved massive star is hard to explain and extremely suspicious. Again, we caution that model calculations may be less reliable for such extreme parameters than for the previous cases.

Similarly to M33-UIT045, it has been extremely difficult to produce a synthetic spectrum without N V to match the observed spectrum of M33-UIT349. The fit was achieved with an extremely low mass loss rate, about two order of magnitude lower than the other stars in the sample. The result is consistent with the fact that the typical wind lines in the UV spectrum are almost absent. However, the low mass-loss rate is surprising since this star has photospheric parameters not too different from the rest of the sample. We must keep in mind that additional faint sources may contribute -although marginally - to the spectrum, but they are not resolved nor confirmed. It is also puzzling that the wind lines in the M33-UIT349 spectrum seems wider rather than narrower, compared to the spectra with stronger lines. Typically, a progression from stronger to weaker P Cygni profiles (i.e. decreasing mass-loss) also corresponds to decreasing wind velocities, as it evident e.g. in the sequence (from the top down) of spectra in Figures 1 and 5. The large velocity may suggest a physical binary, i.e. two stars with high ($``$ 1000km s-1) velocity relative to each other, both having weak lines. Or, more simply, the fact that the lines are so weak (and the limited S/N) prevents to resolve blends and makes the analysis inconclusive in this case.

Finally, for two “weak-lined” spectra, modeling was achieved with a carbon underabundance (respect to solar) of a factor of two, and not twenty like the hotter stars with stronger winds, and lower mass values. The derived parameters $`T_{\mathrm{eff}}`$, $`L_{\mathrm{bol}}`$, Mass, are again consistent with the evolutionary models of Maeder, since lower mass stars produce less enhancement of N and less depletion of C, than more massive stars.

Looking at the absolute flux level and slope at longer wavelengths (Figure 7), we note a slight mismatch between model and observed spectrum for M33-UIT008, M33-UIT045 and M33-UIT003. While the best fit model reproduces well the continuum and all the lines in the shorter wavelengths range (where the flux from the hot star dominates), the slope of the entire observed spectrum cannot be matched perfectly by any single star model for any combination of reddening (using known reddening types). The mismatch between model and observed spectra, however, is of the order of 0.1 - 0.2 10-14ergs cm-2 s-1, or $``$ up to 5% beyond $``$ 1500Å, thus comparable to the calibration uncertainty. The flux level and the line spectrum are matched well in the short wavelength range, pointing to a consistent and plausible model fit for the hot star. In the finding chart produced from the STIS CCD acquisition image and the spectral intensity cuts shown in Figure 4, no extra faint components can be clearly resolved. Therefore, the mismatch could indicate either an unresolved chance superposition of a fainter cooler star (statistically unlikely), or a faint binary component with separation less than $``$ 0.2pc (unresolved), or an UV extinction slightly different than the types known so far in the MW or LMC (see next section). The level of mismatch is however of the order of the calibration accuracy (see Section 2.2.), therefore it could just be an instrumental effect.

## 4 Discussion and Conclusions

We analyzed STIS Ultraviolet spectra of six Ofpe/WN9 stars in M33. Previous optical spectra of our sample (MBHS) showed a progression from “weak-lined” slash stars, similar to BE470 in the LMC (with strong P Cygni HeI and Balmer emission but very weak NIII $`\lambda `$ 4634,42 and HeII $`\lambda `$ 4686) to “strong lined” slash stars, similar to BE381 in the LMC (with strong NIII and HeII). The UV line spectra show a progression (Figure 1) of line intensities which does not correlate with the earlier optical spectra (Figure 2), suggesting possible variability on short time-scale (years) of the mass-loss rate.

One target, M33-UIT349, was observed twice, 8 days apart, but successfully only on Dec. 27. Its spectrum seems to include a very faint extension (Figure 4), separated spatially by $`<`$0.1” (thus if physically associated, by $`<`$ 0.4pc), but the source is not resolved into separate components at the STIS CCD resolution (Figure 3). In any case, the apparent secondary component has negligible flux with respect to the main source (figure 4), as confirmed by the good fit to the spectral distribution. The spectrum shows extremely weak wind lines.

For three stars (half of the sample) we find C/N=0.02 by mass, consistent with predictions from stellar evolution models for an evolved massive star (initial mass $``$40-50 $`M_{}`$) in the phase transitioning from Blue supergiant to WNL (i.e. towards the W-R stage) at constant luminosity in the H-R diagram. The other physical parameters ($`T_{\mathrm{eff}}`$, mass, luminosity) are also entirely consistent with the scenario of massive stars approaching the WNL stage (Maeder 1987, Maeder & Meynet 1994, Meynet & Maeder 2000, 2003). It must be kept in mind that the available spectroscopic range provides several wind lines but no photospheric lines, therefore derivation of abundances is not very accurate as wind line strengths depend also on the mass-loss rate and other factors. However, the consistency of the entire modeling (lines, absolute flux, velocity) supports the results as quite reliable, and all the parameters indicate a consistent evolutionary picture. For the other three stars, C was found to be (by number), only about half the solar value. For M33-UIT104, the analysis indicates N to be overabundant by a factor of two, i.e. C/N = 0.9 by mass. For the other two stars with extremely weak lines, the nitrogen abundances and the overall modeling is way more uncertain, for the reasons explained in the previous section. Either these stars are in an earlier evolutionary stage, evolving off the ZAMS, still towards lower $`T_{\mathrm{eff}}`$’s, or their progenitors has lower initial masses (by comparison with Maeder’s evolutionary calculations). Variability of the mass loss rate is also possibile, as mentioned, complicating the interpretation.

Studies of coeval clusters in the Milky Way (Massey, DeGioia-Eastwood, & Waterhouse 2001) and Magellanic Clouds (Massey, Waterhouse, & DeGioia-Eastwood 2000) suggest that LBVs are descendent from the most massive stars ($`>90_{}`$). The two Ofpe/WN9 stars in their sample, however, come from stars of much lower mass, $`>`$ 25-35$`_{}`$, calling into question the evolutionary link between LBVs and Ofpe/WN9s. The range of progenitor masses inferred by comparing our derived current masses, abundances, $`T_{\mathrm{eff}}`$ and luminosities with evolutionary calculations ($``$ 30-50$`M_{}`$), is consistent with their findings. On the other hand, a link between LBVs and Ofpe/WN9 has been observationally established in the cases e.g. of AG Car, R127, S61, S119 (Stahl 1987, Stahl et al. 1983, Pasquali et al. 1999, Nota et al. 1995) although these studies do not provide enough statistics to quantify the correlation in terms of evolution. Langer et al. (1998) and Lamers et al. (2001) discuss the effect of stellar rotation, on the evolution of massive stars, and show that rotation ”extends” the occurrance of the LBV phase to lower progenitor masses.

One important parameter, the mass-loss rate, varies from 8 10-6 M yr-1to 2 10-6 M yr-1among the first four stars in table 2, which appear to have similar masses (within the uncertainty of our analysis) but a range of temperatures. The uncertanties on the mass-loss rates are less than 30% according to our model calculations. Although the difference in $`T_{\mathrm{eff}}`$ within our sample covers a fairly narrow range, and the sample is numerically limited, there is a trend of mass-loss increasing as the star evolves towards the W-R stage and becomes hotter. Mass-loss rate increases with luminosity (e.g. Bianchi & Garcia 2002, Nugis & Lamers 2000, Vink et al. 2001), but the small range of luminosity variation among the sample does not explain the large mass-loss rate spread (see below).

The values of $`\dot{M}`$ that we derive for the “strong-lines” stars (the first four in Table 2) are lower than the high mass-loss rates of W-R stars, and higher than mass-loss rates expected for O stars with similar parameters. Maeder’s evolutionary models assume an average $`\dot{M}`$= 4 10-5 M yr-1for WNL stars, and the Nieuwenhuijzen & de Jager (1990) parametrization for the earlier stages (this parametrization is based on a compilation of observed mass loss rates). Our mass-loss rates are also lower than Nieuwenhuijzen & de Jager (1990) parametrization which, for the physical parameters given in table 2, would predict for the first four stars of the sample $`\dot{M}`$ between 2.1 and 1.5 10-5 M yr-1. Recent works revised the mass-loss rates for Wolf-Rayet stars, taking into account “clumping” in the wind. These mass-loss rates are generally lower than earlier determinations, of the order of $``$ 1 10-5 M yr-1for WN-type stars of comparable luminosities to our targets (Nugis & Lamers 2000; Nugis, Crowther & Willis 1998; Graefner, Koesterke & Hamann 2002), thus higher than what we measure for the Ofpe/WN9 stars. To compare with O-type stars, we use the recent “recipe” (based on theoretical monte-carlo simulations) from Vink et al. (2000, 2001), for pop.I stars, and find that the predicted mass-loss rates are of the order of 2 10-6 M yr-1for the first five stars (2.4 to 1.7 10-6 M yr-1from M33-UIT236 to M33-UIT045), for solar metallicity, and less for metallicity lower than solar (e.g. for z=0.1 $`\times `$ solar, $`\dot{M}`$ would be lower by a factor of seven). In summary: the values of mass-loss rates derived for our sample are lower than typical W-R mass-loss rates, but in comparison to predictions from radiation-pressure wind theory for pop.I stars, mass loss is enhanced for the hottest stars (M33-UIT236, M33-UIT008, M33-104), is comparable to the predictions for M33-UIT003, and lower than predicted for M33-UIT045 and M33-UIT349. However, as we stressed before, $`\dot{M}`$ derivations are less reliable for the two latter stars, as the lines are extremely weak. Another result emerging by these comparisons is that mass-loss rate varies within the sample more than we would expect from the variation of the other physical parameters ($`T_{\mathrm{eff}}`$, $`L_{\mathrm{bol}}`$, mass). For this fact, and because the UV line strengths do not correlate with optical line strengths observed several years earlier, we suggest that mass-loss rate may vary in time as well. LBV stars also display variable mass-loss rates, which Vink & de Koter (2002) explain in terms of changes in the line driving efficiency. The large spread in mass-loss rates among our sample may also relate to evolutionary effects. The physical conditions of the radiation-pressure driven wind change during the evolutionary phases from O-type to LBV to W-R, and these classes of objects display very different mass loss rates and wind velocities, although the stellar luminosity remains basically constant during the transition. According to Lamers & Nugis (2002), the wind changes are mainly due to change in stellar parameters (radius, gravity) and to a lesser extent to change in surface chemical abundances. The ensemble of our findings consistently points out that the first four objects in Table 2 are transitioning towards the W-R (WNL) stage through an enhanced mass-loss phase.

All spectra are very well fitted in both continuum flux level and line strength at the short wavelengths (G140L range), while a mismatch between observed flux and model is seen longwards for most stars. All stars appear to be point-like sources at the STIS spatial resolution. The amount of the mismatch (up to 5% at the longer wavelengths) is comparable to the absolute flux confidence level (section 2.2), and to the uncertainty in the extinction amount and extinction law. A variation in E(B-V) of 0.01 would have an effect comparable to the mismatch seen beyond 1700Å. However, it would create a larger mismatch at the shorter wavelength, which are more sensitive to the reddening, thus was excluded. The amount of extinction ($`E(BV)`$) adopted as our best fit to the spectra is given in table 2, and the adopted extinction law is a combination of Galactic extinction (foreground) and LMC-type extinction. This was found to be generally the preferrable solution, after examining different combinations with all known extinction laws (for other galaxies), concurrently with varying the stellar model parameters. The UV extinction law in M 33 is currently being investigated by us with another STIS program, and the results may refine this issue. We attempted to add different black-body components to the stellar-model fits. While we may visually improve the fit to the observed spectra at the long wavelengths, removing the mismatch, the results would not be significant given that the flux level of a postulated additional component would be up to a few % that of the luminous hot star. Thefore, we cannot conclude at this point, whether a faint excess flux, peaking at about 1800Å or longwards, is present in addition to the brighter hot luminous star spectrum, or if the extinction curve slightly differs from the MW and LMC ones.

Finally, we point out again two remaining possible sources of uncertainty in the present analysis. The first may come from some inconsistency in the WM-basic calculations, for the cases of extreme parameters such as those of M33-UIT045 and M33-UIT349, where the flux indicates a hot, luminous object but the lack of wind lines (or just of NV in the case of M33-UIT045) indicates an unusually low mass-loss rate. These caveats are discussed at length by Bianchi & Garcia (2002,2003b), Garcia & Bianchi (2003). Second, the amount of shocks in the winds and the effect of the related X-rays on the ionization - which bares on the $`T_{\mathrm{eff}}`$ results, cannot be precisely estimated when only the UV range is available. The FUSE range, at shorter wavelengths, contains transitions much more sensitive to this parameter than those in our STIS spectra, primarily the O VI doublet (Bianchi & Garcia 2002, Bianchi et al. 2003). We plan to extend the analysis of the sample to the FUSE range as a next step, although the fluxes of these distant objects are at the limit of detection for FUSE.

It is interesting to compare our findings with similar objects in other galaxies, to verify metallicity effects on the evolution. To our knowledge, no similar modeling of UV spectra has been performed for Ofpe/WN9 stars in the MW or other galaxies. However, previous works, mostly based on optical spectroscopy, provide very useful comparisons. Similarly to our STIS spectra, UV FOS spectra of all the LMC “slash stars” show P Cygni profiles of CIII $`\lambda `$ 1176, SiIV$`\lambda `$$`\lambda `$ 1394,1403 and NIII $`\lambda `$$`\lambda `$ 1748, 1752. Other strong lines, CIV$`\lambda `$$`\lambda `$ 1548,1550, NV $`\lambda `$$`\lambda `$1238,1243, HeII $`\lambda `$ 1640, NIV$`\lambda `$1719 and AlIII $`\lambda `$$`\lambda `$ 1855,1863, vary from pure absorption to strong P Cygni profiles in the sample (Pasquali et al. 1997). these authors derive terminal velocities of the order of 400 km s-1 , much lower than the wind velocities measured in our sample, and mass-loss rates (from H$`\alpha `$ equivalent widths) of $``$ 2 - 5 10-5 M yr-1, much higher than what we found for the M33 counterparts. The discrepancy may be due to the H$`\alpha `$ equivalent width method, where the correction for the underlying photospheric absorption is a very large source of uncertainty (Bianchi & Scuderi 1999). A non-LTE analysis of the optical spectrum of an Ofpe/WN9 star in NGC 300 was performed by Bresolin et al. (2001) using the Hiller & Miller (1998) code. The $`T_{\mathrm{eff}}`$ is lower, but the luminosiy and the mass loss rate higher, than the range of values found in our sample. Crowther et al. (1977) performed a non-LTE analysis of optical spectra of LMC’s Ofpe/WN9 stars, which they however re-classified as WNL (WN9-11). They find a range of $`T_{\mathrm{eff}}`$ which overlaps with the lower end of our range, most of their sample (WN9-11) stars having lower $`T_{\mathrm{eff}}`$’s than ours, and higher $`\dot{M}`$ values. We caution again that mass loss rates derived from UV and from optical spectra give often discrepant results (e.g. Crowther et al. 2002, Bianchi & Garcia 2002), thus the comparison has to be taken with caution. A more conclusive comparison would require similar analysis of similar data-sets.

Support for this work (proposal GO 8207) was provided by NASA through grant GO8207-0197A from the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS5-26555. We are very grateful to Miriam Garcia for computing several WM-basic models of the grid used in this work, and to Jorick Vink, Lars Koesterke and Wolf-Rainer Hamann for useful discussions.