© ESO, 2012

1. Introduction

Wolf-Rayet (WR) stars are the evolved descendants of massive OB stars (M > 25   M_⊙), and they are believed to constitute the previous phase in their evolution before supernova explosion. Many of them have been shown to be surrounded by nebular emission in a structure of shells (ring nebula). Johnson & Hogg (1965) were the first to suggest that these nebulae were due to material ejected from the star sweeping up the surrounding interstellar medium (ISM). Pikel’Ner & Shcheglov (1969) and Weaver et al. (1977) both made models deriving analytical solutions about dynamic evolution of shock bubbles created by interaction between the ISM and the stellar wind. One of the most comprehensive models describing the basic optical morphology of bubbles in WR stars was produced by Garcia-Segura et al. (1996). Their hydrodynamical simulations take the evolution and mass-loss history into account over the lifetime of the central star from the main sequence (MS) to the WR phase. In summary, stellar winds from the MS star sweep up the ISM, forming a hot bubble surrounded by a thin shell; later, dense material from the red supergiant (RSG) winds partially fills up the bubble creating another internal shell; finally, in the WR stage, the star develops another mass-loss episode with a fast wind that reaches the RSG shell and can eventually break up the MS bubble, creating a wind-blown bubble (WBB).

NGC 6888 (S105) is a proto-typical WBB associated with the WN6 star WR 136 (HD 192163) (Johnson & Hogg 1965). It was first reported by Sharpless (1959), and its relative proximity and large optical angular size has made it one of the best-studied examples of this class of objects, which are observed in a wide wavelength range. The basic parameters of star and nebula are summarized in Table 1.

The ellipsoidal shell appears to be geometrically prolate and to have a highly filamentary structure (Lozinskaya 1970). The expansion velocity of the shell varies in the range 55−110 km s^-1 (Lozinskaya 1970; Marston & Meaburn 1988; Esteban & Vilchez 1992, among others). Spectroscopy of the nebular shell shows that the ionized gas is enriched with nitrogen and helium (Wendker et al. 1975) and is slightly underabundant in oxygen (Kwitter 1981; Esteban & Vilchez 1992). This could be explained as a consequence of the mixing of interstellar and stellar wind material (Kwitter 1981) and the transport from the core of the star to the shell during the RSG phase (Esteban & Vilchez 1992). Marston (1995) presented IRAS images obtaining two concentric shells associated with WR 136. With the inferred masses, he proposed that material of the inner shell comes from the mass lost by the star in the RSG phase and the outer shell formed from the early O star phase. Gruendl et al. (2000) observed Hα and [OIII]λ5007 Å emission from NGC 6888, finding an offset between both images and suggest that these offsets are due to the shock from the WR bubble expanding into the circumstellar medium. This contrast between the Hα and [OIII]λ5007 Å morphology is shown better in the images from the WFPC2 on the HST (Moore et al. 2000).

Bubbles blown by massive stars are filled by the fast stellar wind that interacts strongly with the ISM, and some fraction of the kinetic energy is converted into thermal energy heating the gas at X-ray emitting temperatures (Weaver et al. 1977). NGC 6888 was the first WBB that has been detected in X-rays, it was made with Einstein by Bochkarev (1988). The Wrigge et al. (1994) X-ray observation did not detect a filled bubble emitting as models predicted, but a filamentary X-ray structure was instead found following the Hα features of the nebula. In addition, spectra in the shell were fit to a double-temperature plasma model with a high-temperature component at T ~ 8.5 × 10⁶ K and a low-temperature component at T ~ 1.3 × 10⁶ K without detecting significant temperature variation between different areas of the nebula (Wrigge et al. 2005). In the last X-ray study of this nebula, Zhekov & Park (2011) have found a bi-modal distribution of the X-ray emitting plasma, and made models to explain how the soft-X rays come from evaporation of clumps. They also discuss whether the heating conduction is efficient in the hot bubble.

NGC 6888 represents one of the most suitable candidates for searching for the observational footprint of the interactions between ISM and stellar winds. In addition, it can give us clues about the mass-loss history of massive stars, as well as about the enrichment of the ISM (Esteban et al. 1992; Monreal-Ibero et al. 2010). Despite all the studies already done on NGC 6888, the spatial extent of the chemical enrichment and the homogeneity of the physical properties and structure are still unknown, since no 2D measures of their spectroscopic properties (and imaging, simultaneously) have been done so far. To study the morphology and chemodynamics in 2D, we decided to include NGC 6888 as a target in our programme of integral field spectroscopy (IFS) to perform the spatial analysis of spectral maps and compare them with 1D studies of selected integrated spectra.

The paper is organized as follows. The observations and data reduction are described in Sect. 2. Description, analysis, and results of the physical conditions of the 2D study are presented in Sect. 3, and Sect. 4 is devoted to study selected 1D spectra. The discussion and the summary of the main conclusions are given in Sects. 5 and 6, respectively.

2. Observations and data reduction

2.1. Observations

Observations were carried out in 2005 July at the 3.5-m telescope of the Centro Astronómico Hispano Alemán (CAHA) at the observatory of Calar Alto (Almera, Spain) with the Potsdam Multi-Aperture Spectrograph instrument (PMAS) (Roth et al. 2005) in PPAK mode (Kelz et al. 2006). The PPAK configuration consists of 331 science fibres with a diameter of 2.7 arcsec covering a hexagonal field of view (FoV) of 74  ×  65 arcsec². The surrounding sky is sampled by 36 additional fibres distributed in six bundles located following a circle at about 90′′ from the centre. Furthermore, there are 15 fibres for calibration purposes (see Fig. 5 in Kelz et al. 2006).

We chose the zones of NGC 6888 to be mapped based mainly on the narrow band images observed by our group in the Isaac Newton Telescope (INT) with the Wide Field Camera (WFC). Figure 1 shows an RGB composited image of the nebula; Hα and [OIII] emissions represented by red and green colours respectively, the [OII] emission, in blue, is fainter but noticeable as the pink shade in the mix with Hα. Here, we present IFS observations of NGC 6888 in four zones of the nebula (see Fig. 1): the first zone (1) in the northeast area, was chosen because it can be a good candidate for presenting signs of shocked gas, and it is known that this area emits in X-rays. Two other pointings, (2) and (3), were chosen in the outer zones of NGC 6888 where Hα emission is very weak but an [OIII] skin is observed in Fig. 1. Finally, another zone (4) was selected to study the centre of the nebula where we can observe a striking dark region in all the wavelengths or a knot that could be associated to our object. We have called all these regions: X-ray zone (1), edge (2), mini-bubble(3), and bullet (4) zones, as marked in Fig. 1.

We used three different gratings to obtain information for the most important emission lines in the optical spectral range. The R1200 grating used was centred in two wavelengths covering the spectral range from 6220 to 6870 Å and from 4320 to 5060 Å with a dispersion of 0.3 Å/pix, giving a spectral resolution of  ~2.4 Å FWHM at 6300 Å. A V1200 grating with a dispersion of 0.35 Å/pix was used that covered the spectral range from 3660 to 4430 Å with a spectral resolution of FWHM ~ 2.7 Å at 4047 Å. This setup overlaps 100 Å with the wavelength range of the R1200 blue grating. Finally, we used a low-resolution configuration with the V300 grating, covering the spectral range from 3660 to 7040 Å with a dispersion of 1.67 Å/pix, giving a spectral resolution of  ~8.7 Å FWHM at 5577 Å. Table 2 shows the observational log for NGC 6888 in the four regions and three exposures were taken in each pointing. Moreover, in the X-ray zone with R1200 gratings, we adopted a dithering scheme with three positions having offsets measured from the centre of Δα =  −1.35″, Δδ =  −0.78″ for the second pointing and Δα =  −1.35″, Δδ =  + 0.78″ for the third. This pattern allowed us to fill up the gap between fibres.

The weather was photometric throughout the observations with the typical seeing subarcsecond. Two spectrophotometric standard stars, Hz 44 and BD+28 4211, were observed for all the setups. The necessary bias frames, continuum, and arcs were acquired.

2.2. Data reduction

The data were reduced using the R3D software (Sánchez 2006) in combination with the Euro3D packages (Sánchez 2004) and IRAF¹. The reduction consisted of the standard steps for fibre-fed IFS observations.

A master bias frame for each night of observations was created by averaging the corresponding bias frames observed during the night and subtracted from all the images. The different exposures taken at the same position on the sky were combined to reject the cosmic rays using the package IMCOMBINE of IRAF.

Using continuum-illuminated exposures taken before the science exposures, we determined the location of spectra on the detector for each pixel along the dispersion axis. Then, each spectrum was extracted from the science and standard star frames, co-adding the flux within an aperture on the location of the spectral peak in the raw data using the tracing information, and storing in a 2D image called row-stacked-spectrum (RSS) (Sánchez 2004). Furthermore, we checked that the contamination from flux coming from adjacent fibres using this aperture was negligible (Bershady et al. 2004; Sánchez 2006) . For a given aperture and FWHM ~ (0.5 ×    aperture) we found a level of cross-talk that is always  <10%. This seems to be an acceptable compromise between maximizing the recovered flux and minimizing the cross-talk.

Wavelength calibration was performed in a two-step procedure using arc calibration-lamp exposures (a neon lamp for the R1200 grating in the range from 6220 to 6870 Å and He for the rest of the pointings). In the first step, distortion correction was carried out to avoid the curvature. Then, the wavelength coordinate system was determined by identifying the lines. Finally, the distortion and dispersion solutions were applied over the science data. The accuracy achieved was better than  ~0.1 Å (rms) for the arc exposures in every spectral range. We performed a further wavelength “re-centring” using sky emission lines present in the observed wavelength ranges.

Corrections to minimize the differences between fibre-to-fibre transmission throughput were also applied, creating a fibre-flat with an exposure of a continuum. All the science frames were divided by the obtained fibre-flat.

Observations of spectrophotometric standard stars were used to perform the flux calibration. Only Hz 44 was observed all the nights, so we used it for our primary calibration. With the E3D software a spectrum for the calibration star was created each night by co-adding the central fibres of the standard star frames. Then, with IRAF tasks, a sensitivity function was determined. Subsequently, we checked that the flux-calibration for the two standard stars was self-consistent for all the nights observed and with the library values. Finally, flux calibration was applied to all the science frames.

To achieve an accurate sky subtraction is an important issue in data reduction. Although PPAK has additional fibres for sampling the sky, our object is very extended, and all the sky fibres are located within an area containing signals from the nebula. Extra frames to sample the background near NGC 6888 were not available. For each night, we adopted different sky-subtraction schemes. For the night when we used the V1200 grating, an external sky-frame was observed. Since the sky and object exposures were taken within a few minutes of each other in the same night and under similar weather conditions, we decided to estimate a sky spectrum using this extra frame. To do so, the spectra of all individual fibres of the sky were combined with a mean in a single spectrum and a 2D spectrum was created with the value in each fibre of the combined spectrum. This sky was subsequently subtracted from every science spectrum. This method was used in observations with the high-resolution gratings (R1200 and V1200), and no satisfactory sky frame was available for the night of V300 observations.

Considering the short wavelength range and the airmass of the observations, and using Filippenko (1982), we estimated that the offsets due to the differential atmospheric refraction (DAR) were always smaller than one third of the diameter fibre. Correction for DAR was not necessary in our data.

3. Two-dimensional study

Once data had been reduced, we first carried out a basic 2D study to describe the observed morphology of the four pointings in several emission lines. Later, in zones with a dithering scheme, we performed a more detailed analysis.

3.1. Morphology in the emission lines

Using E3D software, we generated interpolated images of the four observed regions in three wavelength ranges including the lines: [OII]λλ3726, 3728 Å, [OIII]λ5007 Å, and [NII]λλ6548, 6584 Å+Hα. Figure 2 displays these images. It can be seen two different patterns in all the regions, [OII] and Hα+[NII] present similar behaviour, while [OIII] shows a different structure. This can be understood as the signature of the ionization structure of the nebula.

In the X-ray emitting region, most of the emission is located in the northwest and seems spread towards the northeast and centre. The high ionization emission (i.e. [OIII]) shows two main peaks, one in the northwest and the other in the southwest.

The second zone (Edge) was chosen to be observed because is placed in a big bubble of [OIII] that appears surrounding all the nebula. It can be observed in Fig. 2, that in the [OIII]λ5007 Å emission-line, a discontinuity separates the nebula into two different regions, while it appears uniform in the low ionization emission-line images.

The little bubble in the northern region of the nebula (Mini-bubble) is cleary seen in the [OIII] images, but appears very diffuse in Hα in Fig. 1. In the PPAK observations of this region different features can be distinguished, a bright zone to the south in all the wavelengths and a semi-circle structure only seen in the [OIII] image.

Finally, the last row of Fig. 2 corresponds to a darker region placed near the centre of NGC 6888 that could be a sort of “bullet”  crossing the nebula. This region appears quite similar in the three wavelength maps. The dark zone in the north is very faint in all the emission lines, though again low ionization lines show slight differences from [OIII].

3.2. Data-cube analysis of the X-ray zone

The bi-dimensional structure of the three dithering observations was studied in detail by means of the creation of data cubes. We performed this kind of analysis only for the X-ray zone as observed with R1200 gratings in two wavelength ranges because for the other pointings dithering was not available. In the next Sects. 3.2.1–3.2.3, all the results refer to the X-ray zone.

As explained above, the 331 science fibres lead the light from different positions of the target to the pseudo-slit. When dithering exposures (to fill holes) were performed, it is important to combine them when knowing the overlapping areas and registering the spectra to their positions in the sky. With this aim, we used the R3D software to construct a grid with a pixel scale about one third of the original fibres (~1 arcsec), obtaining easier-to-use data cubes with two spatial and one spectral dimension. Therefore, two cubes of the X-ray region in the wavelength ranges 6220–6870 Å and 4320−5060 Å were built.

One of our objectives was to obtain maps from the datacubes of the main emission lines, which could then be analysed to describe the 2D ionization structure and physical properties of the X-ray emitting region of NGC 6888. This procedure was performed using our own routine, which fits each line with a Gaussian function. We introduce four parameters as a reference to the fitting (the central wavelength, the FWHM, the line wavelength range, and the corresponding continuum range) and maps are returned with useful information in each pixel for each selected line (i.e. flux, FWHM, central wavelength, error of flux, continuum, etc.).

As a first step, we fit all the lines by a single Gaussian for the two data cubes obtained. For the grating placed in the red range (from 6200 to 6870 Å), seven lines were fit: [SIII]λ 6312 Å, Hα, [NII]λλ6548, 6584 Å, HeIλ6678 Å, and [SII]λλ6717, 6731 Å. In the other spectral range (from 4320 to 5060 Å), another seven lines were fit: Hγ, Hβ, [OIII]λλ4363, 4959, 5007 Å, and HeIλλ4471,5015 Å. In addition, to prevent contamination by low signal-to-noise (S/N) data, we masked all pixels with an S/N lower than 5 out. Under these conditions, maps were used in the analysis that we present in the next sections.

Two main problems appeared in the maps creation. First, with those fainter lines placed in zones with high noise, the fits of the automatic routine were not good in all the spaxels, and some adjacent signal was added to the fit, giving a wrong measurement of the line flux. It happened with three lines: [SIII]λ6312 Å, HeIλ4471 Å, and HeIλ5015 Å. With the obtained maps for these lines, masks were created in order to not take bad pixels into account in the analysis. The second problem that we found was fitting the [OIII]λ4363 Å. This line is located very close to the strong HgIλ4358 Å night sky line, in particular, both lines are overlapped in our object, and we cannot resolved them, so we could not obtain workable maps of [OIII]λ4363 Å for this nebula. Without the measure of this line, it was not possible to obtain the electron temperature map. Furthermore, no other auroral line, like [NII]λ5755 Å, was measured in the cubes with a high S/N (but see Sect. 4.1).

3.2.1. Extinction and electron density maps

The reddening coefficient c(Hβ) was derived from the Hα/Hβ line ratio. We used an intrinsic Balmer emission line ratio of Hα/Hβ = 2.81 obtained from the public software of Storey & Hummer (1995) assuming Case B recombination with an electron temperature of T_e ~ 10⁴ K and an electron density of n_e ~ 10² cm^-3.

Although other Balmer lines were detected in our wavelength spectral range, Hα and Hβ were the lines with the highest S/N. We performed the estimation of c(Hβ) including Hγ. The average value obtained was similar but the errors increased, and the S/N in the Hγ flux was lower than 5 in 10% of the spaxels, so finally we only used Hα and Hβ. However, we checked that in the spaxels with good S/N in Hγ both derivations were consistent. The derived c(Hβ) map is shown in Fig. 3, where a non-uniform, rather filamentary structure can be clearly distinguished, with values ranging from c(Hβ) = 0.3 to 1.00, with a mean value of 0.57. Adopting E(B − V) = 0.692 × c(Hβ) and using the Cardelli et al. (1989) extinction law with R_v = 3.1, E(B − V) was computed with values in the range 0.21 − 0.88. In general terms, the structure presented in this map is in good agreement with values reported in the literature. Using the flux density at 21 cm and the flux in Hα, Wendker et al. (1975) estimated an absorption over the whole nebula of A_V ~ 2.05 mag (E(B − V) = A_V/3.1 = 0.66). Esteban & Vilchez (1992) obtained values of the c(Hβ) reddening coefficient of 0.28 for the blueshifted component and 0.50 for the redshifted one, correspond to E(B − V) of 0.19 and 0.34, respectively. Some of these values are slightly lower than ours, but their slit positions do no spatially coincide with our observed zones.

Finally, all the emission-line maps were reddening-corrected using the derived c(Hβ) map. Observing their appearance we checked the morphology explained in Sect. 3.1 but now with more lines: Hα, Hβ, [SII], [NII], and He I with similar patterns as examples of low ionization and the doublet of [OIII] as high ionization.

The electron density (n_e) map of the X-ray region was produced from the [SII]λλ6717/6731 ratio using the package TEMDEN of IRAF based in a five-level statistical equilibrium model (De Robertis et al. 1987; Shaw & Dufour 1995). The spatial distribution of the electron density can be seen in Fig. 4. The map shows a rather filamentary aspect that is more or less uniform with a mean value of 178 cm^-3. The main significant feature in this spatial distribution is the peak density (~400 cm^-3) coincident with the maximum flux emission of all the species detected. Our highest n_e values agree with the density obtained by Kwitter (1981) and Esteban & Vilchez (1992).

3.2.2. Diagnostic diagrams

The excitation mechanisms and the ionization structure have been studied using diagnostic diagrams based on the ratio of strong emission lines: the classic BPT analysis (Baldwin et al. 1981). Using the data from the 2D study of the X-ray zone, we generated maps for the necessary line ratios (with S/N > 5) and plotted all the spaxels in the three diagrams displayed in Fig. 5.

The first diagram (Fig. 5a) shows [NII]λλ6584/Hα vs. [SII]λλ6731/Hα (NvsS) and presents a bivaluate behaviour. We performed a reference fit to the points, paying attention to which spaxels are inside/outside the nominal limits of  ± 3σ from this fit. It can be seen that the relation between these two line ratios shows a significant scatter, especially for the points below the  − 3σ line. We notice that spaxels above the +3σ line show a linear-like correlation as could be expected from the ionization structure given the ionization degree of these ions. To understand these behaviours, we located the points of the diagram in the map differentiating between points inside, above and below the 3σ limits of the global fit, to identify their spatial locations (see Fig. 6). We found that points with a similar tendency in the BPT appear grouped in the PPAK map. The second diagram (Fig. 5b) presents [OIII]λλ5007/Hβ vs. [SII]λλ6731/Hα (OvsS), it shows a larger scatter with two tendencies. To understand how the zones defined above in Fig. 6 behave in this diagram, we identify pixels belonging to these zones with the same colours as in Fig. 5a. Points with large scatter in NvsS (blue crosses) show scatter also in this diagnostic diagram, and points with a correlation in NvsS (red crosses) appear here anti-correlated. In the last diagram (Fig. 5c) we can see the relation between [OIII]λλ5007/Hβ and [NII]λλ6584/Hα (OvsN), where we find a strong anti-correlation. Once again, we study the behaviour of the three regions defined in Fig. 6 in this diagram with the same colour code. In this case, the points appear anti-correlated for all the zones.

Taking the analysis and results explained above into account, we have defined three spatial zones in this pointing: Zone A to the northwestern part of this pointing, corresponding to all the pixels above +3σ in NvsS. We notice that in all the diagrams these spaxels show a normal and expected tendency. In addition, in this area the maximum emission in low excitation species is found. Zone B is located to the southwest and corresponds to points below  − 3σ in NvsS. In this case, the relations between the line ratios in the two first diagrams were rather anomalous and had a large dispersion. In the map this zone corresponds to the peak of the [OIII]λ5007 Å emission. Finally, we defined Zone C by the points inside  ± 3σ limits of the fit in the diagram NvsS, which are located over the rest of the field and present a behaviour that is a mixture between the two other zones.

3.2.3. Kinematics

The radial velocity of the gas was derived using the central wavelength of the Gaussian fit in every pixel. All the measured radial velocities were translated into the local standard of rest (LSR) and corrected for the Earth’s motions using IRAF. In addition, a zero-velocity map was created using a “sky” line, and this map gives us realistic information about the error in the wavelength calibration and the accuracy in the velocity that can be achieved. We concluded that differences in radial velocity below 8 km s^-1 cannot be covered with our data. Figure 7 shows the velocity field derived from the Hα emission line. The map has a mean value  ~84 km s^-1, but it is not spatially homogeneous, because the southwest zone moves faster, reaching 110 km s^-1, while in the northwest zone the velocity is lower than 60 km s^-1.

Taking advantage of the large number of pixels in our maps, we have produced statistical frequency distributions of the radial velocity. When we represent velocities from all the pixels, the histograms present a Gaussian distribution as expected. But, if pixels from Zones A and B are represented separately, a bimodal distribution is found in some relevant emission-lines like [NII]λ6584 Å. Figure 8 shows two histograms generated with a binning of 3 km s^-1 for two emission lines differentiating the two defined zones. The [NII]λ6584 Å line presents two patterns: while Zone A shows a distribution with a peak redshifted by 10 km s^-1, Zone B shows a bimodal distribution indicating a shell-like behaviour with pixels showing gas that is moving away and other pixels coming towards us. On the other hand, in the second histogram, the [OIII]λ5007 Å line presents a single peak distribution in both zones and the velocity difference between the two peaks is  ~15 km s^-1.

With this kinematic statistical study we support the idea that NGC 6888 is a non-uniform nebula with shells and filaments moving with different velocities. Errors in the wavelength calibration and the instrument resolution of the observations prevent us from carrying out a more exhaustive analysis of the kinematics of the nebula in this work.

4. One-dimensional study

A different analysis of the IFS data can be carried out by combining the observed spectra in several fibres to produce integrated spectra in different regions of NGC 6888. With this method we can use data without dithering exposures to obtain information over a broader spectral range and derive some physical properties that could not be obtained using a 2D approach. To decide which fibres are to be combined, we examined the appearance of each zone in emission lines at several wavelengths using E3D as shown in Fig. 9 (see also Fig. 2):

• X-ray zone: two integrated spectra were created in this pointing. The first spectrum (X1) corresponds to spaxels where the [NII]λ5755 Å line is visible. We have chosen this criterion because we were looking for emission lines useful for deriving the electron temperature. In addition, spaxels where [NII]λ5755 Å is bright enough coincide with the brightest low-excitation lines. To define the second region (X2), we selected spaxels placed in the area where we had a peak in the high ionization line emission ([OIII]). Moreover, X1 and X2 are in Zones A and B, respectively.

• Edge zone: we were interested in checking the different conditions and parameters on both sides of the nebular discontinuity observed in [OIII]λ5007 Å, so the first integrated spectrum was performed in the brighter zone to the left (E1) and the second one in the fainter part (E2).

• Mini-bubble: we created integrated spectra over three of the features seen. Firstly, we chose fibres emitting in [NII]λ5755 Å (MB1). Secondly, a bright zone in [OIII]λ5007 Å but without [NII]λ5755 Å emission (MB2). Finally, we defined the zone with low level of emission at all the wavelengths (MB3).

• Bullet zone: to check the nature of the bullet candidate, we created one integrated spectrum over the dark zone (B1). In addition, as a representative example of the emission of the inner area of the nebula, we extracted a spectrum (B2) in this brighter region but outside of the bullet.

Using these definitions, several sets of fibres were selected and combined for each grating in all the zones to obtain their integrated spectra, which were analysed to derive their spectroscopic properties in the whole spectral range of the observations. For observations where a dithering scheme was adopted we combined the fibres of all the exposures. At the end, we obtained 19 1D spectra considering the different gratings. Figure 10 shows an example of integrated spectra in zone X1.

4.1. Integrated properties

The emission-line fluxes were measured using the SPLOT task in IRAF. Lines show a degree of asymmetry in their profiles, so we did not perform Gaussian fits The fluxes were derived measuring areas below the lines, however, partially blended lines, such as Hα and [NII]λλ6584, 6584 Å in the V300 grating, were deblended with two or three Gaussian profiles to measure the line fluxes separately. The statistical errors associated with the observed fluxes were calculated using the expression (1)where σ₁ represents the error in the observed line flux, N is the number of pixels used to measure the line, EW the line equivalent width, σ_c the standard deviation of the continuum in a box near the line, and Δ represents the dispersion in Å/pix (Pérez-Montero & Díaz 2003).

Table 3 lists the emission lines measured in every zone labelled with its standard identification. The estimated fluxes and errors were normalized to F(Hβ) = 100. Line intensities were reddening-corrected using the same procedure as described in Sect. 3.2.1. In this case, other HI Balmer lines with S/N > 5 were measured in the integrated spectra and therefore used to derive the reddening coefficient. In all the spectra, Hα, Hβ, Hγ, and Hδ were measured with enough S/N and used to estimate c(Hβ), with the exception of the pointing B1 in which Hγ and Hδ were very faint. The reddening coefficients derived in zones X1 and X2 agree with the values obtained in the 2D-study in the X-ray pointing using only Hα and Hβ. The values obtained for c(Hβ) are presented in the last row of Table 3. Errors in the emission-line intensities were derived by propagating the observational errors in the fluxes and the reddening constant uncertainties.

In regions X1, X2, E1, and E2, some lines were measured more than once since they appear in exposures with different gratings. In these cases, we calculated the corresponding intensity as the mean weighted by the error in each zone.

To obtain the physical conditions of the gas, we performed an iterative process until achieving agreement between electron temperature and electron density. Table 4 shows the values for T_e and n_e derived for each zone.

With the [SII]λλ6717/6731 line ratio, electron densities (n_e) were obtained using the package TEMDEN of IRAF. We derived differences in the electron density of up to 350 cm^-3 between the studied zones. The results from the X-ray zones are coherent with the values obtained in the bi-dimensional analysis where densities were higher in the northwest of our field of view that corresponds to the integrated spectrum X1. On the other hand, values obtained on the edge of the nebula (E1 and E2) and in the less bright zones of the mini-bubble (MB2 and MB3) are lower than the densities derived in other regions (always n_e < 100 cm^-3).

Electron temperature, T_e, can be derived using the line ratios (2)and (3)The [OIII]λ4363 Å line was blended with the “sky” line HgIλ4358 Å, and therefore we cannot compute T_e([OIII]) directly. In contrast, the [NII]λ5755 Å auroral line was measured with enough confidence in zones X1, MB1, and B2, and direct estimates of T_e([NII]) were obtained from the R_N2 ratio. In the rest of zones (X2, E1, E2, MB2, MB3, and B1), for which this line was not bright enough to be measured with confidence, we resorted to the empirical flux-flux relations from Pilyugin (2005) that relate the faint auroral lines with other strong nebular lines. To find out what was the best relation to our particular case, we used several fits proposed by Pilyugin over the zone X1, for which R_N2 is well measured, and selected which relation represented our results better. We chose the relation Q_NII = f(R₂,P) from Pilyugin (2007) to estimate the electron temperature T_e([NII]) and checked that this relation gives consistent results with our measurements of the line in the X1 zone. Later, from T_e([NII]), we derived T_e([OIII]) and T_e([OII]) in each zone by using the relations based on models proposed by Pérez-Montero & Contini (2009) and Pérez-Montero & Díaz (2003), respectively. Derived temperatures range from  ~7700 K to  ~10 200 K over the nebula, which is consistent with the previous studies of Kwitter (1981) and Esteban & Vilchez (1992).

4.2. Chemical abundances

Considering the electron density and temperatures estimated for each ionization zone, we obtained ionic abundances from the forbidden-to-hydrogen emission-line ratios using the functional forms given by Hägele et al. (2008), which are based on the package IONIC of IRAF. To determine the singly ionized helium abundance, we used the equations proposed by Olive & Skillman (2004). We defined two temperatures: T_e([NII]) as temperature representative for the low ionization ions, S⁺, N⁺, and O⁺, and T_e([OIII]) for the high ionization species, O²⁺ and Ne²⁺, and for deriving the abundance of He⁺. All the ionic abundances are presented in Table 4 with their corresponding errors.

Furthermore, in Table 4 we show the total abundances derived as follows. The O/H has been obtained by directly adding the two ionic abundances (O⁺ and O²⁺). The same procedure was used for He/H but only with He⁺ because He²⁺ is not seen in our observations, the value chosen of HeI was a mean weighted by the errors of the three HeI observed lines (HeIλ4471 Å, HeIλ5876 Å, and HeIλ6678 Å). The N/O can be approximated to N⁺/O⁺. Finally, to derive the total Ne/H abundance, ionization correction factors (ICF) are required to take the unseen ionization stages into account, we used the ICF expression for Ne from Pérez-Montero et al. (2007). Taking Izotov et al. (2007) as reference, we estimated ICF(He)  ~ 1 for all the zones, so correction for the unseen neutral helium was not necessary.

Mainly, three results in the total abundances are interesting. First, we see that nearly all the zones present an oxygen abundance slightly below the solar value, taking the solar abundance from Asplund et al. (2009) as reference. This effect was reported above for the areas observed by Kwitter (1981) and Esteban & Vilchez (1992). On the other hand, these authors derived an overabundance of nitrogen. For some of our zones (E1, E2, MB2, MB3, and B1), we derived a nitrogen abundance similar to the solar value, however in other zones the N/H appears enhanced by a factor of 6, or even 8 in the X1 zone. Helium presents an enrichment in most of the integrated zones too. (MB3, E2, and B1 do not show the HeI line.) One of the most striking results is the important differences in the N/O between the different zones, mainly due to the large N/H variation and not to differences in O/H. This could be understood by assuming that we are seeing structures of the nebula formed in different stages of the evolution of the WR star, each one with its own ionization degree, physical properties, and chemical patterns. This idea is analysed in Sect. 5. Furthermore, we calculated the expected values in the ISM at the galactocentric radius of NGC 6888, considering the galactic chemical abundance radial gradients. Taking Balser et al. (2011) as reference, we derived 12 + log (O/H) = 8.6 ± 0.05 and log (N/O) =  −1.26    ±    0.26 from Carigi et al. (2005). Our inferred abundances for E2, MB3, and B1 are consistent with these expected values (within the errors).

4.3. Photoionization models

To investigate and consolidate both ionization parameter (U) and chemical abundances in the different pointings of the nebula, we performed a grid of photoionization models made with the code Cloudy v.8.0 (Ferland et al. 1998). These models are based on the grid described by Pérez-Montero & Contini (2009), using one single stellar atmosphere of WR, in particular, a WN star model (Smith et al. 2002). Our grids assume an effective temperature of the star of 50 000 K (Rosa & Mathis 1990; Esteban et al. 1993), with a metallicity of Z = 0.008. Taking the solar abundances reported by Asplund et al. (2009) (12 + log (O/H) = 8.69, Z = 0.013) as a reference, the assumed metallicity implies an oxygen abundance 12 + log (O/H)) = 8.47, which is the closest value to the O/H abundances derived in this work for NGC 6888. These models vary the ionization parameter and the nitrogen-to-oxygen ratio as free parameters, which are very relevant for our purpose. We considered values log (U) =  −2.5, −3.0, and  −3.5, and log (N/O) from  −1.5 to 0.0 (in steps of 0.25).

We compared selected emission-line ratios of the integrated spectra with models in two diagnostic diagrams: [OIII]λλ5007/Hβ vs. [NII]λλ6584/Hα and [OIII]λλ5007/Hβ vs. [NII]λ6584/[OII]λλ3726,3729 (Dopita et al. 2000). As can be seen in Fig. 11, for each zone, the values of log (N/O) obtained with the comparison always agree with abundances calculated with our measures. For each integrated spectra, the ionization parameter is similar in the two diagrams: zones X2, E1, and MB2 have the greater U; zones X1 and B2 have intermediate values; and the paremeters are the lowest in zones E2, MB1, MB3, and B1².

5. Discussion

To develop a scheme of the internal structure of NGC 6888, we put together all the information derived from the 2D maps and the nine integrated spectra. Furthermore, we used the narrow-band images from the INT-WFC shown in Sect. 1 to have a global overview of the nebula and to locate regions of different observed conditions. The proposed scheme consists of a shell structure (close to “onion-like”) but allowing irregular and/or broken shells. We found, at least, two shells with different physical, kinematical, and chemical properties: an outer shell (OS) and an inner shell (IS). In addition, a more external and fainter structure is seen surrounding the whole nebula like a thin skin; we call it surrounding skin (SS). A basic sketch of this scheme is displayed in Fig. 12.

Our interpretation of the structure of NGC 6888 is supported by comparison to models of WR ring evolution. In particular, there is very good agreement with the simulation by Garcia-Segura et al. (1996, hereinafter GS96) for the interaction of a 35 M_⊙ star with the circumstellar medium over the entire stellar lifetime in the case of the slow RSG wind (~15 km s^-1). Their simulations start at the MS stage, when the stellar winds sweep up the ISM, forming a hot bubble surrounded by a thin shell. Later, in the red supergiant (RSG) phase, the star is characterized by a slow and very dense wind that shocks with the MS bubble, creating a shell (RSG shell). Finally, a faster wind from the WR phase sweeps up the RSG wind material into a new shell; however, the propagation of this WR shell is faster, and it soon arrives at the RSG shell. After the collision, clumps and dense material of the WR shell are propelled by their high inertia, crossing the RSG shell as bullets forming blow-outs. The final stage in the GS96 simulations (Fig. 7 in that paper) shows that the bubble formed in the main sequence is, eventually, broken out as consequence of the RSG-shell shocked by the WR shell and its expansion.

5.1. The inner shell (IS)

In the proposed scheme, the central part of NGC 6888 is an ellipsoidal shell surrounding the WR star with a high surface brightness, as can be seen in the INT image (Fig. 1). This shape agrees with the scheme proposed previously by Lozinskaya et al. (1988). Based on interferometric observations, they derived a spatial model in which NGC 6888 is a prolate ellipsoid of revolution (a/b    ~ 0.5) whose major axis is inclined 20°–40° to the projection plane.

In our data, the information about this shell was obtained from the integrated spectra X1, X2, MB1, and B2. From the analysis of the properties of these four zones, we propose that this shell is not a complete shell, but a broken structure, with holes, filaments and bubbles. One evidence pointing to this are the variations of the electron density always  >100 cm^-3 ranging from 106 to 360 cm^-3 along the defined IS. Also, differences in the electron temperature and in the chemical composition show inhomogeneities between different zones of this shell.

One of the more important results for the IS is the strong overabundance in N/H by a factor between 5 and 9 with respect to the solar abundance and by a factor up to  ~11 with respect to other regions of the nebula. The enriched N/H found varies from 12 + log (N/H) = 8.29 to 8.72 in the integrated spectra included in this shell, supporting the idea of a patchy and incomplete structure. On the other hand, the O/H abundance is found to be slightly deficient, suggesting that some oxygen could have been processed into N mainly via the ON cycle. In the case of helium, we obtained a value of up to He/H = 0.18 in this shell, indicating a strong overabundance of this element when we compare it with the solar value He/H = 0.09 (Asplund et al. 2009). These facts give us clues to the stage in the star evolution in which this shell originated.

We have additional information about the IS thanks to the observations with IFS because the X-ray pointing (with dithering and full datacubes) is over this region. In Sect. 3.2.2 we defined two zones, Zones A and B, associated with the emission of the low and high ionization ions, respectively, and showing very different trends in the diagnostic diagrams. Later, we found two behaviours in the radial velocities associated to these zones. In particular, the statistical frequency distributions of velocity (Fig. 8) are consistent with a broken shell scheme: the [NII]λ6584 Å velocity histogram in Zone A suggests only one redshifted shell going away from us, while in Zone B the bimodal velocity distribution found can be understood if the shell is moving away in some places and in others is coming towards us. Rarely are the two behaviours observed simultaneously, because we do not see a central and dominant velocity component in the histograms. This kinematical study supports the idea of an intrinsically inhomogeneous shell in expansion, more obvious in the nitrogen than in the oxygen emission. Figure 13 shows a sketch of how the scheme of broken-shell works and what we would see for each zone.

Although in our work we cannot consider the 2D study of the whole nebula, including all the spatial information, we propose that the most internal shell of NGC 6888 should present a structure with inhomogeneities in the physical properties and abundances. Differences in the chemical composition could be understood by attending to the phase in the evolution of the central star when the shell was formed; in particular and in comparison with GS96 simulations, we argue that the IS material corresponds to the red supergiant (RSG) and WR shells shocked expanding in a low-density medium. This would explain its fragmented and filamentary structure and the high emission of Hα. Also, the unusual abundances found indicate that this shell is mostly composed of stellar material ejected from the central star when it became a WR star and later probably mixed with material from the RSG stage. This idea is also supported by the recent simulation performed by Toalá & Arthur (2011). They find that the excess of N/O in NGC 6888 can be explained and reproduced by a stellar evolution model with rotation (Meynet & Maeder 2003) for an initial stellar mass of 37 M_⊙.

5.1.1. Can we detect the footprint of shocks?

The previous observations of X ray emission detected in this region and the theoretical simulations (GS96) lead us to consider that on the north edge of the major axis of this IS, a collision between the WR and RSG shells happened during the star’s lifetime. To check that there are spectral signatures of shocks, we resorted to the models presented by Allen et al. (2008). Their code, MAPPINGS III, gives an extensive library of radiative shock models covering a wide range of abundances and velocities. First, we had the maps of the emission lines ratios created in our 2D study; then, we generated a grid of shock models with different parameters using libraries from Allen et al. (2008); finally, we overplotted all the spaxels from the maps on the theoretical grids in the diagnostic diagrams to find out the consistency with models.

In a first attempt, we plotted all the pixels of our maps in two diagnostic diagrams ([OIII]λλ5007/Hβ vs. [NII]λλ6548,6584/Hα and [OIII]λλ5007/Hβ vs. [SII]λλ6717,6731/Hα), but we failed to find a model that represented all the observational data in both diagrams well. Later, we studied only the region with a striking behaviour in Figs. 5 and 6 (called Zone B). When we plotted only the spaxels that belong to this region, we found that they can be fit with one of the grids of models. Figure 14 represents the two diagnostic diagrams showing the pixels from the Zone B, together with the grid of shock models.

We found that our data of Zone B can be represented well by a model in which shock and precursor are involved assuming n = 1000 cm^-3 and with twice solar abundances. In both diagrams we find that the models that fit better have shock velocities from 250 to 400 km s^-1. The agreement of our data with the predictions of shock models with twice solar abundances can be understood by attending to the strong overabundance of nitrogen. Nevertheless, we should bear in mind that in our 1D study the line ratios of this region (integrated spectra X2) are also consistent with the tailored photoionization model presented in Sect. 4.3. An alternative scenario can be considered that is similar to the ionization of “pillars of creations”  structures (Ercolano et al. 2012). In this scenario, altered line ratios could be produced as a consequence of the partial observation of ionized to neutral surfaces in the nebula.

As a conclusion, considering the shock models, the X-ray emission and the stellar evolution simulations, we cannot rule out the presence of a shock, so if this event has happened, the analysis of our data indicates that the most likely place to find it is Zone B of the X-ray pointing. It would be interesting to check this with new observations of the other side of the major axis (southwest of NGC 6888).

5.2. The outer shell (OS)

Studying the integrated spectra E1 and MB2, we found several differences in the physical properties in comparison with other zones that lead us to think that they belong to a different shell. In our nebular scheme, we propose an intermediate shell placed between the SS and the IS. The morphological study made in the bi-dimensional analysis (Sect. 3.1) supports the existence of different shells or spatial inhomogeneities in the structure. In particular, the images in the second row of Fig. 2, corresponding to the Edge pointing, show a clear discontinuity in the spatial distribution of the emission of the [OIII]λ5007 Å line that is not visible in Hα nor [OII]λ3727 Å. The integrated spectrum E1 was created from the emission to the right of this discontinuity. Something similar happens in the mini-bubble pointing (third row in Fig. 2) where a well-shaped semicircular structure (where the MB2 spectrum was extracted) is seen only in the [OIII] emission line.

Taking the previous observations of NGC 6888 into account (Gruendl et al. 2000; Moore et al. 2000), it is very likely that the OS has an apparent spherical shape, at least in the [OIII]λ5007 Å emission line. The INT-WFC image (Fig. 1) agrees with this morphology: it can be seen that the oxygen distribution is not as elliptical as the other elements, showing some extended emission beyond the end of the minor axis of the IS.

This bubble is surrounding the inner shell and looks dominant in double ionized oxygen. Intensities measured of the emission lines are weaker here than in the inner shell; nevertheless, the [OIII]λ5007 Å line has the maximum flux in comparison with the rest of the integrated spectra. Here we estimate the highest electron temperatures of the whole sampled nebula (>10 000 K). The expansion of the shell could explain the electron density found here, always lower than 100 cm^-3. The chemical analysis shows that the He/H is enriched with respect to the solar values, while N/H and O/H appear not to be enhanced, or they present a slight deficiency. It is interesting to notice that the O/H abundance is similar to the one of the IS, but N/H is lower by a factor  ~8. It was not possible to estimate the neon abundance here because these lines were not bright enough in the integrated spectra.

We inferred the gas distribution in this shell studying data from the maps in the X-ray zone, because, although in our scheme the X-ray pointing does not target the OS, it gives us other interesting information. Remembering the statistical frequency velocity distributions in Fig. 8, it can be seen that the overall range in radial velocity covered by [OIII]λ5007 Å and [NII]λ6584 Å emissions is similar (~80 km s^-1), but the [OIII]λ5007 Å velocity distribution shows a single peak for the two defined zones (A and B), whereas the [NII]λ6584 Å emission presents a bi-modal distribution in Zone B and a redshifted peak in Zone A. This behaviour can indicate that the O²⁺ appears to be distributed, forming a rather filled expanding bubble.

Analysing the kinematics, abundances, and physical properties, we propose an origin for this shell that could explain why it is only observed in the direction of the minor axis of the ellipse defined above. We argue that ejections of material in the RSG stage were non-spherical (Humphreys 2010), forming an ellipsoidal RSG shell denser towards its major axis. When the WR wind collides with the RSG shell, along the minor axis, the collision does not find strong resistance, fragmenting it and breaking up the MS-surrounded bubble (as predicted in GS96 models). The shock then propagates in a tenuous medium where the high temperature favours the [OIII]λ5007 Å emission seen and where a lower cooling rate could prevent the material from the WR stage to be observed in the optical. This idea is basically consistent with previous findings of Gruendl et al. (2000).

5.3. Surrounding skin (SS)

Finally, we propose the existence of an external and faint shell around the whole nebula, the surrounding skin. This shell is barely noticeable after close inspection of the INT-WFC images, but analysis of the integrated spectra E2 and MB3 shows that a faint emission exists beyond the outer shell. Although the zones where both spectra were extracted do not present any apparent emission in Fig. 9, non negligible emission-line fluxes were measured on the spectra.

Properties obtained from these two zones reveal that they belong to a different region than the rest of integrated spectra. Both appear to have similar electron temperatures and densities, and we find here the lowest N/H abundance (lower than solar). Neither helium nor neon lines were detected in these integrated spectra. The N/O and O/H in these regions are in the range of the expected values of the ISM at R_G = 10 kp (NGC 6888 galactocentric radius) considering the chemical abundance gradients. Also, I(Hβ) has its lowest intensity in E2 (except for the bullet, as explained below). We checked that the emission in this shell is intense enough not to be confused with general diffuse ISM emission, but it is not strong enough to belong to the OS. The derived oxygen and nitrogen abundances appear consistent with the ISM expected values at this galactocentric distance.

Following GS96, we think that this skin should represent the early interaction between the bubble created by winds from the MS star with the circumstellar medium.

5.4. The bullet

As we explained in Sect. 1, observations in the pointing called Bullet were made to shed light on the nature of a dark object near the centre of the nebula. Measures in this zone (integrated B1) indicate that it is the region where the Hβ flux (and other emission lines) is lowest, but still not negligible. In this region, abundances of N and O are slightly lower than for the solar values of Asplund et al. (2009), but in the same range as the outer shell. The physical properties inferred (T_e and n_e) are not remarkable in relation to other regions of the nebula.

Attending to the measured intensities, we propose that this object is not behind the nebula because in that case all the nebular emission would be measured in front of it, so it is probably located between the main body of NGC 6888 and us, but not out of the nebula, because in that case all the emission lines would be blocked. With this analysis and considering the physical properties of B1, we suggest that the bullet is placed in the surrounding skin defined above, and it should only be blocking light coming from the internal shells (see Fig. 15).

The observed bubble can also be explained using the GS96 scheme. Maybe such a clump formed in the early collision between MS winds and ISM, or perhaps as an accumulation of molecules and dust of the very cool RSG wind, later impelled by the WR shock. In any case, new observations in other wavelengths, such as radio and infrared, would be necessary to understand the nature and composition of this object.

6. Summary and conclusions

In this article, we have presented an exhaustive study of the wind-blown bubble NGC 6888 with integral field spectroscopic observations carried out with PPAK. The most relevant results presented in this paper can be summarized as follows.

• 1.

  We obtained interpolated images of four regions of the nebula in several wavelengths. We observe two different patterns for the low and high ionization ions. Some features are only observed in the [OIII]λ5007 Å emission line, indicating clear differences in the excitation conditions.

• 2.

  In the X-ray emitting zone of the nebula, a dithering scheme was adopted in the observations. This allowed us to generate cubes to carry out a thorough bi-dimensional analysis of this region:

  • (a)

    We obtained a map for the reddening coefficient, c(Hβ), finding a clear non-uniform structure, rather filamentary with a mean value of  ~0.57.

  • (b)

    We derived the electron density map, n_e, showing a distribution with a peak density (~400 cm^-3) spatially coincident with the zone with the highest emission-line fluxes in the nebula.

  • (c)

    We generated radial velocity maps. The velocity field derived from the Hα emission line is not spatially homogeneous and is moving faster in its southwest part.

  • (d)

    We studied the ionization conditions plotting pixels from intensity line ratios maps in three diagnostic diagrams. From these, we generated a map defining two main spatially well delimited zones: Zone A showing correlations between lines as could be expected from ionization structure, given the different ionization degree of the ions, and where the maximum emission in low excitations species was found; and Zone B where the relations between the line ratios show significant scatter in diagnostic diagrams that is located in the area corresponding to the peak of the [OIII]λ5007 Å emission.

  • (e)

    We compared the statistical frequency distributions of the radial velocity for different lines. The [OIII]λ5007 Å histogram shows a single and centred peak in both zones (A and B). The [NII]λ6584 Å presents one redshifted peak in Zone A and a bimodal distribution in Zone B. This result could be explained as the two components of an expanding and broken shell.

  • (f)

    We evaluated the possible presence of shocks to explain the observed behaviours in Zone B, while resorting to the models presented by Allen et al. (2008). We found that our data of Zone B can be represented well by models in which shock and a precursor are involved, assuming n = 1000 cm^-3, twice solar abundances, and shock velocities from 250 to 400 km s^-1.

• 3.

  We generated nine integrated spectra extracted from the four pointings of the nebula. We measured line intensities to obtain physical parameters and chemical abundances:

  • (a)

    We derived electron temperatures, T_e, using the R_N2 = [NII] λλ6548,6584/[NII] λ5755 ratio, obtaining values ranging from  ~ 7700 K to  ~10 200 K. From the [SII]λλ6717/6731 line ratios, electron densities n_e, were calculated by deriving differences of up to 350 cm^-3 between the studied zones.

  • (b)

    We inferred that nearly all the zones present a slightly underabundance of oxygen with respect to the solar abundances. We derived a nitrogen abundance similar to the solar value in five regions; however, in other zones N/H appears enhanced up to a factor of 6 or even 8. Helium abundance presents an enrichment in most of the integrated zones, too.

  • (c)

    We compared selected emission-line ratios of the integrated spectra with predictions from photoionization models on diagnostic diagrams. The values obtained with these models for log (N/O) agree with the abundance derived.

• 4.

  We provided a scenario for the evolution of the central WN6 star to explain the features observed in NGC 6888. The proposed scheme consists of a shell structure, but allows irregular and/or broken shells with different physical, kinematical, and chemical properties:

  • (a)

    An inner, elliptical and broken shell (IS), with a strong overabundance in N/H and slightly deficient in O/H. We argue that the IS material corresponds to the RSG and WR shells that were shocked and expand in a low-density medium.

  • (b)

    An outer and spherical shell (OS) that is weaker in all the observed emission lines except in [OIII]λ5007 Å. The helium appears enriched here, while N/H and O/H do not appear enhanced or do not present any slight deficiency. We argue that this shell is the early MS bubble broken up as a consequence of the collision between RSG and WR shells.

  • (c)

    We propose the existence of an external and faint shell around the whole nebula, the surrounding skin (SS). Measured fluxes in this regions are the lowest but not negligible. We think that this skin should represent the early interaction between the bubble created by winds from the MS star with the interstellar medium. The abundances inferred for this component appear to be typical of the local ISM.

  • (d)

    We checked that the dark region (bullet) observed in the centre of NGC 6888 appears to be associated with the nebula and is probably located in its most external shell. Observations in other wavelengths, such as radio and infrared, would be necessary to fully understand the nature of this object, but we speculate that it could be a dense clump impelled by the WR winds crossing the shells like a bullet.

---

Acknowledgments

This work is supported by the Spanish Ministry of Science and Innovation (MICINN) under the grant BES-2008-008120. This work has been partially funded by the projects: AYA2007-67965-C01 and AYA2010-21887-C04-01 of the Spanish PNAYA; CSD2006-00070 “1st Science with GTC”  from the CONSOLIDER 2010 programme of the Spanish MICINN and TIC114 of the Junta de Andalucía. S.F.S. thanks the Viabilidad, Diseno, Acceso y Mejora funding programme, ICTS-2009-10, and the Plan Nacional de Investigacion y Desarrollo funding programme, AYA2010-22111-C03-03, of the Spanish Ministerio de Ciencia e Innovación, for the support given to this project. We thank Mercedes Mollá for her useful help with the models of the WR stellar atmospheres.

References

Online material

All Tables

All Figures