© ESO, 2012

1. Introduction

Blue compact dwarfs (BCDs) are low-mass galaxies that are experiencing a starburst. They are usually characterized by small physical sizes (~2–3 kpc), low metallicities (0.2 ≲ Z/Z_⊙ ≲ 0.02), and relatively large amounts of gas (M_HI/L_B ≳ 1). The question has arisen as to whether they are young galaxies undergoing their first burst of star formation (Searle & Sargent 1972), but several studies based on surface brightness and colour profiles (e.g. Gil de Paz & Madore 2005), as well as colour–magnitude diagrams of resolved stellar populations (e.g. Tosi 2009) have demonstrated that BCDs also contain old stars, with ages  >2–3 Gyr. The star-formation histories of the nearby BCDs, as derived using colour–magnitude diagrams (e.g. Tosi 2009), show that the starburst is a short-lived phenomenon, typically sustained for  ~100 Myr (e.g. McQuinn et al. 2010). Thus, BCDs are transition-type dwarfs, but it is unclear whether there are evolutionary connections with dwarf irregulars (dIrrs), spheroidals (dSphs) and/or ellipticals (dEs) (e.g. Papaderos et al. 1996; van Zee et al. 2001). In addition, the mechanisms that trigger, sustain, and quench the starburst activity are not understood.

Various studies of the H i  distribution and kinematics (e.g. van Zee et al. 1998a, 2001) have highlighted two striking properties of BCDs: i) they have strong concentrations of H i  within the starburst region near the galaxy centre; ii) they usually have steep central velocity gradients. Both properties are not observed in more quiescent dIrrs (e.g. Swaters et al. 2002). This suggests that there is a close connection between the starburst, the compact distribution of baryons (gas and stars), and the H i  kinematics. The nature of the steep velocity gradients is unclear and two main interpretations are possible: i) rapid rotation (e.g. van Zee et al. 2001), ii) gaseous inflows/outflows (e.g. Kobulnicky & Skillman 2008, Cannon et al. 2004). Rapid rotation in the inner regions could imply that there is a strong concentration of mass that may be either luminous or dark. Gaseous inflows/outflows might be linked to fuelling processes and/or feedback mechanisms.

On larger scales, BCDs usually have extended and diffuse H i  structures, which may form reservoirs for fueling the starburst. In general, two different kinds of structures are observed: i) extended H i  disks in regular rotation, e.g. NGC 2915 (Elson et al. 2010) and NGC 2366 (Oh et al. 2008); ii) complex filamentary structures, e.g. II Zw 40 (van Zee et al. 1998a) and NGC 5253 (Kobulnicky & Skillman 2008). The study of these extended H i  structures can provide key information about the triggering mechanism (external or internal processes), the properties of the progenitor galaxies (gas-rich dIrrs or gas-poor dEs/dSphs), and the possible presence of massive gas inflows/outflows.

We present an H i  study of I Zw 18, the BCD prototype (e.g. Zwicky 1966; Searle & Sargent 1972), which is one of the most metal-poor galaxies known (12+log(O/H)  ~  7.2, Izotov & Thuan 1999) and has long been considered a candidate young galaxy in the Local Universe, formed within the past 0.5 Gyr (e.g. Papaderos et al. 2002; Izotov & Thuan 2004). However, Aloisi et al. (2007), using HST observations, detected stars older than 1–2 Gyr and ruled out the possibility that I Zw 18 is a truly primordial galaxy (as suggested by e.g. Aloisi et al. 1999, Östlin & Mouhcine 2005). The intense star-formation activity started only  ~20 Myr ago (e.g. Aloisi et al. 1999). Two key questions remain open: What triggered the starburst? Why is the metallicity so low?

Previous H i  studies (Viallefond et al. 1987; van Zee et al. 1998b) showed that I Zw 18 is characterized by a strong central concentration of H i  and a steep velocity gradient, as is typical of BCDs. In addition, the optical galaxy is surrounded by an extended H i  envelope that was described by van Zee et al. (1998b) as “a fragmenting H i  cloud in the early stages of galaxy evolution”.

We analysed archival H i  data to: i) investigate the relation between the gas distribution and the starburst activity; ii) clarify the nature of the velocity gradient; and iii) study the structure, kinematics, and origin of the extended gas.

2. Data reduction and analysis

We analysed public H i  data taken from the VLA archive. The observations were carried out between 1993 and 2004, using the VLA in all four configurations (see Table 1). Data from the B, C, and D configurations were presented by van Zee et al. (1998b). In this new analysis, we also included data taken in 2004 with the high-resolution A-array configuration. The correlator was used in 2AD mode, with a total bandwidth of 0.8 MHz (~165 km s^-1). An on-line Hanning taper was applied to the data, producing 127 spectral line channels with a width of 6.3 kHz (~1.3 km s^-1).

The raw UV data were interactively flagged, calibrated, and combined using the AIPS package and following standard VLA procedures. The UV data were mapped using a robust weighting technique (Briggs 1995) and various Gaussian baseline tapers to attenuate the longest baselines. We built three datacubes with different spatial resolutions and pixel sizes by using different combinations of the robust parameter and the taper FWHM (see Table 2). After various trials, we chose the combinations that minimize sidelobes and wings in the beam profiles.

After the Fourier transform, the datacubes were analysed using the Groningen Imaging Processing SYstem (GIPSY) (van der Hulst et al. 1992). Continuum maps were constructed by averaging line-free channels. Owing to the narrow bandwidth of the observations, a few line-free channels were available and the resulting continuum-subtracted datacubes displayed correlated noise in the spectral direction. Thus, we constructed a continuum map by using a mask, defining the area of H i  emission in every channel and averaging, for each spatial pixel, all the channels without any H i  signal. The masks were constructed by smoothing the datacubes both in velocity (by a factor 4) and spatially (by a factor  ~3, i.e. at 45′′, 10′′, and 5′′ for the low, intermediate, and high resolution data, respectively) and clipping at 2.5σ_s (where σ_s is the noise in the smoothed cubes). The masks were inspected channel by channel and any remaining noise peaks were blotted out.

The use of a mask for the continuum subtraction may have the disadvantage that the noise is no longer uniform across the channel maps, as a different number of channels is used at every pixel to build the continuum map. Thus, we built signal-to-noise ratio (S/N) maps for every channel (similarly to Verheijen & Sancisi 2001) and calculated a pseudo-1σ contour by averaging the values for the pixels with 0.75 < S/N < 1.25. The resulting pseudo-1σ level is close to that obtained by calculating the noise in a box without signal, suggesting that the noise is still almost uniform.

The channel maps were cleaned (Högbom 1974) down to 0.5σ, using the masks to define the search areas for the clean-components, which were then restored with a Gaussian beam of the same FWHM as the antenna pattern. In order to improve the S/N, the cubes were smoothed in velocity to a resolution of 5.2 km s^-1 (10.4 km s^-1 for the high-resolution data) and spatially to 20′′, 5′′, and 2′′ for the low, intermediate, and high resolution data, respectively. Table 2 summarizes the properties of the cubes.

Total H i  maps were constructed by summing the signal inside the clean-masks. A pseudo-3σ contour was calculated following Verheijen & Sancisi (2001). Velocity fields were derived by fitting a Gaussian function to the H i  line profiles. Fitted Gaussians with a peak intensity less than 2.5σ and a FWHM smaller than 5.2 km s^-1 were discarded; the remaining noise in the velocity fields (i.e. signal outside the pseudo-3σ contour of the total H i  maps) was blotted out. The H i  line profiles are quite broad and asymmetric, thus the velocity fields must be considered as only a rough indication of the global kinematics. Our kinematical analysis is based on position–velocity diagrams (Sect. 3) and on three-dimensional (3D) models of the observations (Sect. 4.1).

3. H I distribution and kinematics

We describe the overall H i  structure of I Zw 18. We adopt the standard nomenclature introduced by Davidson et al. (1989). The main body is designated as I Zw 18 A (Fig. 2, top-right) and is characterized by two starburst regions: one to the north-west (NW) and one to the south-east (SE). The light concentrations denoted by Davidson et al. (1989) as B, D, and E are background galaxies. The stellar complex to the NW is named I Zw 18 C or C-component. We assume a distance of 18.2 Mpc, as derived from the tip of the red giant branch (Aloisi et al. 2007) and confirmed by observations of Cepheids (Fiorentino et al. 2010).

We use data at three different resolutions (see Table 2) to probe different spatial scales and H i  column densities.

3.1. The low-resolution view

The low-resolution data (with FWHM = 20′′    ~ 1.8 kpc and 3σ column density sensitivity N_H   I(3σ) ~ 4 × 10¹⁸ atoms cm^-2 per channel) illustrate the large-scale overall structure of I Zw 18.

Figure 1 (top) shows the channel maps, superimposed with two isophotes of a B-band image (from Gil de Paz et al. 2003). The central H i  emission presents a velocity gradient at a position angle PA  ≃  140°–150°. At velocities from  ~750 to  ~720 km s^-1, extended emission also appears to the south.

Figure 2 (top-left) shows the integrated H i  map, superimposed on a V-band HST image (from Aloisi et al. 2007). I Zw 18 A is associated with a strong concentration of gas, while diffuse emission extends beyond the optical galaxy, covering an angular size of  ~3$\stackrel{\mathrm{\prime }}{.}$5 (~18.5 kpc). The H i  gas to the south of I Zw 18 A displays a tail-like morphology that extends over  ~2$\stackrel{\mathrm{\prime }}{.}$5 (~13.5 kpc).

Figure 2 (middle-left) shows the velocity field. The main body (I Zw 18 A) is associated with the central velocity gradient. The southern “tail” does not seem to be kinematically connected to the SE region of I Zw 18 A, as the gas velocity changes abruptly from  ~790 km s^-1 to  ~720 km s^-1. Moreover, at the junction between I Zw 18 A and the “tail”, the H i  line profiles are double peaked, suggesting that there are two distinct components, which are possibly well-separated in space but projected on the same location on the sky.

Figure 2 (bottom-left) shows a position–velocity (PV) diagram, obtained from the 20′′ datacube following the “tail” (the dashed line overlaid on the velocity field). The central velocity gradient (associated with I Zw 18 A) is very steep and there is a spatial broadening towards the NW direction between  ~700 and  ~780 km s^-1. The gas to the south forms a coherent kinematical structure at velocities between  ~710 and  ~760 km s^-1. Strikingly, the broadened part of the PV-diagram and the “tail” are almost at the same velocities, suggesting that they may be physically connected. This hypothesis is investigated further in Sect. 5.

3.2. The intermediate-resolution view

The intermediate-resolution data (FWHM = 5′′ ~ 440 pc and N_H   I(3σ) ~ 4 × 10¹⁹ atoms cm^-2 per channel) illustrate the H i   emission associated with I Zw 18 A and I Zw 18 C, and their possible connections.

Figure 1 (bottom) shows the channel maps, superimposed with two isophotes of a B-band image. The H i  emission to the NW is spatially resolved but still visible, whereas the southern “tail” is completely resolved out, indicating the diffuse nature of this gas. Between  ~770 and  ~740 km s^-1, there are H i  clumps near the C-component. Their association with I Zw 18 C is very likely, because the H i  clumps are at the same velocities as the Hα emission (V_{sys,   Hα} [C ]  = 751 ± 5 km s^-1, Dufour et al. 1996).

Figure 2 (top-right) shows the integrated H i  map, overlaid on a V-band HST image. The main body is characterized by two H i  peaks, roughly corresponding to the NW and SE starburst regions. With respect to I Zw 18 A, the H i  is more extended to the north-west, in the direction of I Zw 18 C.

Figure 2 (middle-right) shows the velocity field, while Fig. 2 (bottom-right) shows a PV-diagram obtained by following the dashed line superimposed on the velocity field. The steep velocity gradient is aligned approximately along the two H i  peaks. The H i  emission to the NW shows a shallow velocity gradient from I Zw 18 A to I Zw 18 C and seems to connect the two stellar bodies. This connection is also visible in Hα at velocities similar to those observed in H i  (Dufour & Hester 1990; Dufour et al. 1996). Along I Zw 18 A, the Hα velocity gradient, instead, shows a “wiggly” behaviour that is not observed in H i. This may be caused by an Hα superbubble (see Martin 1996 and Sect. 5.2). In addition, there is H i  emission to the west of the main body, that shows a velocity gradient and seems to have an Hα counterpart (see Sect. 5.2).

3.3. The high-resolution view

The high-resolution data (FWHM = 2′′ ~ 180 pc and N_H   I(3σ) ~ 2 × 10²⁰ atoms cm^-2 per channel) show in detail the H i  emission associated with the NW and SE starburst regions.

Figure 3 shows the total H i  map in grayscale (left) and in contours overlaid on a Hα image (middle). The two H i  peaks are spatially resolved. The H i  clump to the SE coincides with a complex of H ii  regions and has a strong peak, where the H i  column density reaches  ~100 M_⊙ pc^-2 (~1.2 × 10²² atoms cm^-2). In-between the two clumps, an H i  hole is associated with a strong H ii  region, suggesting that the neutral gas has been consumed, ionized, and/or blown out by young stars. The H i  clump to the NW coincides with a Hα shell, which surrounds the bulk of the young stars (see Fig. 1 of Cannon et al. 2002). This shell is probably connected to the high-velocity Hα emission detected by Dufour et al. (1996) and Martin (1996) at  ± 200 km s^-1 with respect to the systemic velocity.

The velocity field at 2′′ resolution is shown in Fig. 3 (right). This velocity field is very uncertain because of the clumpy H i  distribution, the asymmetric line profiles, and the low S/N of the data at this high angular resolution. However, it shows a clear velocity gradient from the SE to the NW region, as already observed at lower spatial resolution. The H i  depression is in the approaching NW side of the galaxy.

4. Dynamics of I Zw 18 A

In Sect. 3, we described the overall structure of I Zw 18. Two important findings need to be explained: i) the steep velocity gradient associated with I Zw 18 A; ii) the extended H i  emission to the south and to the north-west of the galaxy. In this section, we focus on the dynamics of I Zw 18 A, while in Sect. 5 we study the large-scale gas emission.

4.1. Kinematical models

There is some controversy in the literature about the interpretation of the velocity gradient of I Zw 18 A. Viallefond et al. (1987), Petrosian et al. (1997), and van Zee et al. (1998b) analysed velocity fields and interpreted the gradient as rotation, whereas Skillman & Kennicutt (1993) and Dufour et al. (1996) obtained long-slit spectroscopy and argued that the gradient may result from the merger of two (or more) gaseous clouds. The velocity gradient is along the optical major axis of the galaxy (Fig. 4) and the velocity field displays a pattern that is typical of rotation: this strongly suggests that there is a rotating disk. Here we present 3D kinematical models that demonstrate that the H i  disk is differentially rotating and has a global inflow/outflow motion.

The disk is modelled by a set of gas rings with fixed values of centre, systemic velocity, position angle, inclination, surface density, thickness, velocity dispersion, and rotation velocity. The centre, the systemic velocity, and the position angle were estimated by eye using both optical and H i  data (see Table 3). The centre is between the NW and the SE starburst regions (see Fig. 4, right). For the radial distribution, we used the H i  surface density profile derived from the total H i  map at 2′′ resolution by azimuthally averaging over ellipses (Fig. 6, top). For the vertical distribution, we assumed an exponential law exp( − z/z₀). We built a set of models with different values for the inclination i, the scale height z₀, and the velocity dispersion σ_H   I, assuming that each of these parameters is constant with radius. The inclination and the scale height are constrained by the observed H i  map; their values are slightly degenerate but do not strongly affect the final result; we assumed that i = 70° and z₀ = 100 pc (see Fig. 4, right). The mean velocity dispersion is constrained by the shape of different PV-diagrams and values larger than  ~10 km s^-1 are ruled out; we assumed σ_H   I = 7.5 km s^-1. A thickness of 100 pc and a H i  velocity dispersion of 7.5 km s^-1 are values typical of a H i  disk.

The actual H i  distribution of I Zw 18 A is clearly not axisymmetric (see Fig. 3). Thus, once we had fixed the structural and geometrical parameters of the disk, we built models with axisymmetric kinematics but a clumpy H i  distribution, i.e. the surface density varies with position as in the observed H i  map. The procedure is as follows: we built a disk model with a uniform density distribution (fixing z₀, σ_H   I and the rotation curve), projected it on the sky, and then renormalized the H i  line profiles at every spatial position to reproduce the H i  density distribution observed at 2′′ resolution. For the rotation curve, we tried two extreme cases: solid body (slowly rising rotation curve) and differential (steeply rising and flat rotation curve).

Figure 5 (top) shows PV-diagrams obtained from both different models and observations at a resolution of 2′′ and 5′′. The slices are taken along the major and minor axes, as indicated in Fig. 4 (left). The velocity gradient along the major axis is grossly reproduced by all the models. This demonstrates that: i) a rotating disk is a good representation of the data; ii) the asymmetry between the NW and the SE region is mostly caused by the clumpy H i  distribution. Moreover, it is possible to discriminate between solid body and differential rotation. The observed PV-diagram along the major axis indicates that there is H i  emission close to the galaxy centre (R < 5′′) at high rotational velocities (~810 and  ~730 km s^-1). The solid body model does not reproduce this emission, as the gas in the inner radii is mostly concentrated near the systemic velocity. The differentially rotating model, instead, correctly reproduces the high-velocity gas. This is clearly illustrated by the channel maps at 2′′ resolution at receding velocities (Fig. 5, bottom): the solid body model (dashed line) is not extended enough towards the galaxy centre, whereas the differentially rotating model (solid line) provides a good match of the observations. The approaching NW side of the galaxy is not reproduced as well as the receding SE side, perhaps because the NW starburst region is more active than the SE one (Contreras Ramos et al. 2011). Nonetheless, on this side, a differentially rotating disk is also preferable to a solid body one.

A simple rotating disk (independently of the assumed rotation curve) cannot reproduce the observed PV-diagrams along the minor axis (Fig. 5, top), because of the kinematic asymmetry and the presence of H i  emission at velocities forbidden by circular motions. This effect is also visible in the velocity field at 5′′ resolution (Fig. 4, left): the kinematic minor axis, defined by the contours close to the systemic velocity, is not orthogonal to the kinematic major axis. This is usually attributed to radial motions (e.g. Fraternali et al. 2002). Alternatively, the non-orthogonality between the minor and major axes may be due to a bar-like or oval distortion (e.g. Bosma 1978), although the optical images of I Zw 18 provide no evidence of these strong distortions. We improved the differentially rotating model by adding a global radial motion of 15 km s^-1. The resulting model reproduces the H i  emission at forbidden velocities. The model, however, cannot reproduce all the details present in the observed PV-diagram. The non-circular motions are not uniform across the disk and small variations (of the order of 3–4 km s^-1) could account for the observed discrepancies. It is impossible to discriminate between inflow and outflow, as it is unclear which side of the disk is the nearest to the observer. Vertical motions with roughly the same speed as the radial ones may also be present. Non-circular motions in excess of 20 km s^-1 are ruled out, confirming that the disk kinematics is dominated by rotation.

The rotation curve used to build our best model was not derived from a standard tilted-ring fit to the velocity field (Begeman 1987), but by building 3D models by trial-and-error. The uncertainties in the rotation velocities are difficult to quantify, but we conservatively estimated them by assigning an error equal to ΔV/2.35 = 4.4 km s^-1, where ΔV is the velocity resolution of the 2′′ datacube. The first point of the rotation curve is the most uncertain, as the velocity dispersion in the inner ring may be larger than the mean value of 7.5 km s^-1. For example, if σ_H   I were 10 km s^-1 larger than the mean value, the rotation velocity would decrease by  ~5 km s^-1. Following Meurer et al. (1996), we calculated the asymmetric drift correction, which is smaller than the errors.

4.2. Mass models

In Sect. 4.1, we showed that I Zw 18 A has a rotating H i  disk. The rotation curve is uncertain, but has an inner steep rise and an outer flat part, which indicate that there is a strong central concentration of mass that may be either luminous or dark. Using this rotation curve, we built mass models to estimate the relative contributions of luminous and dark matter to the gravitational potential, following Begeman (1987).

The contribution of the gaseous disk was computed using the surface density profile derived from the total H i  map at 2′′ resolution (Fig. 6, top), multiplied by a factor of 1.4 to take into account the presence of helium. The possible gravitational effect of the H i  components outside the disk (i.e. the “tail” to the south and the extensions to the north-west and to the west) was not taken into account. Molecular gas was not explicitly considered in the mass model because its amount is very uncertain (Leroy et al. 2007). However, if molecules are distributed in a similar way to the stars, their contribution is reflected as an increase in the stellar mass-to-light ratio (M_∗/L). In agreement with the models in Sect. 4.1, we assumed an exponential vertical distribution with z₀ = 100 pc.

The contribution of the stars was computed using the R-band surface brightness profile from Papaderos et al. (2002) (Fig. 6, middle), which was derived from an HST image after subtracting the nebular emission (dominated by the Hα line). The colour profiles of Papaderos et al. (2002) (their Fig. 11) show that, after subtracting the nebular emission, the colour of I Zw 18 A is almost constant with radius. Thus, it makes sense to use a constant value of M_∗/L_R. We assumed a stellar disk with a vertical density distribution given by ρ(z) = sech²(z/z₀) (van der Kruit & Searle 1981), with z₀ = 100 pc.

For the dark matter distribution, we assumed a pseudo-isothermal halo described by the equation ${\mathit{\rho }}_{\mathrm{ISO}}\mathrm{\left(}\mathit{r}\mathrm{\right)}\mathrm{=}\frac{{\mathit{\rho }}_{\mathrm{0}}}{\mathrm{1}\mathrm{+}\mathrm{(}\mathit{r}\mathit{/}{\mathit{r}}_{\mathrm{c}}{\mathrm{)}}^{\mathrm{2}}}\mathit{,}$(1)where the central density ρ₀ and the core radius r_c are both free parameters of the mass models.

Figure 6 (bottom) shows the “maximum disk” decomposition of the rotation curve. The resulting stellar mass-to-light ratio is M _∗ /L_R ~ 1.5. In the maximum disk hypothesis, the baryons dominate the gravitational potential in the inner regions of the galaxy, while the dark matter halo dominates in the outer parts. The parameters of the halo are uncertain because the sampling of the rotation curve is quite course. The halo shown in Fig. 6 (bottom) has ρ₀ = 833 × 10^-3 M_⊙ pc^-3 and r_c = 0.13 kpc.

A M_∗/L_R ~ 1.5 implies a stellar mass of  ~9 × 10⁷ M_⊙. According to Aloisi et al. (1999), the starburst started  ~20 Myr ago with a star-formation rate (SFR) of 6 × 10^-2 M_⊙ yr^-1, giving a mass in young stars of  ~10⁶ M_⊙¹. Thus, the newly formed stars and the concentration of H i  cannot explain the steep rise in the rotation curve, implying that the mass concentration consists of either old stars, molecules, or dark matter. Old stars were detected by Aloisi et al. (2007), and their total mass can be constrained by deriving the galaxy star-formation history from color–magnitude diagrams. The maximum-disk value requires a mean SFR of  ~7 × 10^-3 M_⊙ yr^-1 over the past 13 Gyr, that cannot be ruled out. For the molecules, the upper limit to the H₂ mass within  ~400 pc is  ~7 × 10⁵ M_⊙ (Leroy et al. 2007)², using a Galactic CO-to-H₂ conversion factor (X_CO). However, Leroy et al. (2007) argued that, in I Zw 18, X_CO may be 10^-2 times the Galactic value. The same result is found by extrapolating the relation between X_CO and metallicity by Boselli et al. (2002) down to the metallicity of I Zw 18. Thus, the H₂ mass within  ~400 pc may be dynamically important and be as high as  ~7 × 10⁷ M_⊙.

We also used MOdified Newtonian Dynamics (MOND) (Milgrom 1983, Sanders & McGaugh 2002) to fit the rotation curve. We assumed that a₀ = 1.21 × 10^-8 cm s^-2 (Begeman et al. 1991) and the distance D = 18.2 Mpc (Aloisi et al. 2007), thus the only free parameter is M_∗/L_R. Using MOND, acceptable fits are possible for M_∗/L_R = 1.5 using the “standard” interpolation function (Milgrom 1983) and M_∗/L_R = 1 using the “simple” one (Famaey & Binney 2005).

Following McGaugh (2011), we check the position of I Zw 18 A on the baryonic Tully-Fisher relation. The galaxy follows the correlation within the observed scatter.

5. The extended emission

In this section we study the extended H i  emission. This may provide some clues to the mechanism that triggered the starburst. In addition, we compare the large-scale H i  and Hα emission to investigate the possible presence of outflows.

5.1. The C-component and the H I tail

In Sect. 3.1, we reported two puzzling results (see Fig. 2, left): i) an H i  “tail” at line-of-sight velocities of 710–760 km s^-1 extending to the south of I Zw 18 A and kinematically disconnected from the south-east side of the central rotating disk; ii) a broadening of the H i  emission in the PV-diagram at almost the same velocities (700–780 km s^-1) to the north-west (in the direction of I Zw 18 C). This is clearly shown in Fig. 7, where the PV-diagrams at resolutions of 20′′ and 5′′ (with different column density sensitivities) are plotted on the same scale.

To study these components and their possible connection in more detail, we subtracted the compact H i  disk of I Zw 18 A from the surrounding extended H i  emission. We used the high-resolution data to define the emission from the disk and subtracted this emission from the low-resolution datacube³. Subsequently, we smoothed the residual datacube to both 30′′ and 10.4 km s^-1 and used it to obtain: i) a total H i  map by summing the channels in the velocity range  ~700–780 km s^-1; ii) a velocity field by estimating an intensity-weighted mean velocity; iii) a PV-diagram by following the tail (dashed line in Fig. 8).

Figure 8 shows the results of the subtraction. Interestingly, extended H i  emission is centered on I Zw 18 C (left panel) and forms a coherent kinematical structure with velocities ranging between  ~700 and 800 km s^-1 (right panel). The physical association of the C-component with this surrounding H i  emission is likely because the Hα systemic velocity of I Zw 18 C is  ~751  ±  5 km s^-1 (Dufour et al. 1996) (see star in Fig. 8, right).

The southern H i  tail seems to be connected in both space and velocity to the H i  structure around I Zw 18 C, as shown by the velocity field (Fig. 8, middle). The connection may be either behind or in front of I Zw 18 A. Possible interpretations of the extended H i  emission are discussed in Sect. 6.1.

5.2. The connection between the H I and Hα emission

In Sect. 3.3, we pointed out the relative H i  and Hα distributions in the inner regions of I Zw 18 A. Here we compare the distribution and kinematics of neutral and ionized gas on larger scales.

Figure 9 shows an Hα image (from Gil de Paz et al. 2003) superimposed with the H i  emission (black and white contours). The Hα emission extends well beyond the stellar body and is almost perpendicular to the H i  disk, suggesting that it traces an outflow. This interpretation agrees with that of Martin (1996), who used Hα long-slit spectroscopy and detected a bipolar superbubble expanding with velocities of  ± 60 km s^-1 out to  ~2 kpc from I Zw 18 A. The H i  to the north-east of the main body (see black contour) seems to border the Hα emission and is almost at the same velocity, suggesting that part of the diffuse H i  may be associated with the outflow.

The Hα emission also presents a prominent arc to the west of I Zw 18 A. This feature was first identified by Dufour & Hester (1990) and interpreted as a radiation-bound ionization front driven into the ISM. Petrosian et al. (1997), instead, argued that the Hα arc also contains stellar emission and suggested that it is a structure with stars that are able to ionize the gas “in situ”. The H i  to the west of I Zw 18 A is associated with the Hα arc. To show this, we summed the channel maps at $\mathrm{2}\stackrel{\mathrm{″}}{.}\mathrm{8}\mathrm{\times }\mathrm{3}\stackrel{\mathrm{″}}{.}\mathrm{3}$ resolution in the velocity range 800–770 km s^-1 (see Fig. 9, white contours). Moreover, the Hα velocity field of Petrosian et al. (1997) shows a gradient along the arc similar to the one observed in Fig. 2 (middle-right), confirming the physical association between H i  and Hα. The hypothesis of a radiation-bound ionization front is difficult to reconcile with the presence along the Hα arc of a high-density (~10²¹ atoms cm^-2) neutral gas.

6. Discussion

6.1. Observational evidence and interpretation

This H i  study of I Zw 18 has shown that:

• I Zw 18 A has a compact rotating disk with very high H i  densities. The rotation curve is flat with a steep inner rise, indicating that there is a strong concentration of mass. A global inflow/outflow motion is also present.

• I Zw 18 C is located in the direction of the major axis of I Zw 18 A and is almost at the same velocities as its approaching side. Gas emission with a smooth velocity gradient connects the two stellar bodies. I Zw 18 C appears to be at the centre of a diffuse H i  structure.

• An H i  tail extends to the south of I Zw 18 A out to  ~13.5 kpc. The tail has a coherent kinematical structure and seems to be connected with the H i  emission to the north-west.

Studies of the resolved stellar populations (Aloisi et al. 2007; Contreras Ramos et al. 2011) have shown that: i) the two starburst regions in I Zw 18 A (NW and SE) are embedded in a common envelope of old stars with ages  >1 Gyr; ii) I Zw 18 A and I Zw 18 C are two completely separate stellar bodies and there are no stars between them; iii) I Zw 18 C also contains both old (>1 Gyr) and young (~10 Myr) stars, but its current star formation rate (SFR) is lower than that of I Zw 18 A.

For the interpretation, we consider first the hypothesis of an interaction/merger of two (or more) gas-rich dwarfs. It is well-known that interactions/mergers can produce tidal tails (e.g. Toomre & Toomre 1972). In addition, numerical simulations (e.g. Hibbard & Mihos 1995) suggest that mergers can lead to gas inflows, produce strong gas concentrations, and trigger intense star-formation. Thus, an interaction/merger may provide an explanation of: i) the concentration of H i, ii) the ongoing starburst, and iii) the southern H i  tail. The C-component may be either a “relic” of the interaction or a dwarf galaxy that is interacting/merging with I Zw 18 A. The two objects are at a projected distance of  ~2.2 kpc, the difference between their systemic velocities is  ~12 km s^-1, and are connected by H i  emission with a smooth velocity gradient. The ratio of the R-band luminosities of I Zw 18 A to I Zw 18 C is  ~14, thus this would be classified as a minor merger. The merger hypothesis may also explain the extremely low metallicity, as we now discuss.

Bekki (2008) argued that BCDs with low nebular metallicity are the results of mergers between gas-rich dwarfs with extended H i  disks. According to his simulations, the central starburst is fuelled with metal-poor gas transferred from the outer regions of the extended disks, where the star formation and the chemical enrichment were inefficient owing to the low H i  densities. Our results are consistent with this picture, as the C-component is surrounded by an extended H i  structure, that does not have a stellar counterpart and has probably not been efficiently enriched by SN explosions. Therefore, this H i  structure may provide “fresh” unprocessed gas into the starburst regions of I Zw 18 A. A similar mechanism of metal dilution was proposed by Ekta & Chengalur (2010) to explain why I Zw 18 and the other extremely metal-deficient BCDs are outliers of the mass-metallicity relation.

The extended H i  emission, including the tail, has been considered above as supporting evidence of a merger, but might this instead be the result of a blowout from the starburst? The Hα observations of I Zw 18, indeed, suggest that there is an outflow (Sect. 5.2). In addition, numerical simulations predict that starbursting dwarfs undergo massive outflows because they have a shallow gravitational potential (e.g. Mac Low & Ferrara 1999). The rate of the outflowing gas dM_out/dt can be roughly estimated as $\frac{\mathrm{d}{\mathit{M}}_{\mathrm{out}}}{\mathrm{d}\mathit{t}}\mathrm{=}\frac{\mathrm{2}\mathrm{\times }\mathit{\epsilon }\mathrm{\times }\mathit{SNR}\hspace{0.17em}\mathrm{\times }\overline{{\mathit{E}}_{\mathrm{SN}}}}{{\mathit{V}}_{\mathrm{esc}}^{\mathrm{2}}}\mathit{,}$(2)where SNR is the rate of supernovae (SN), $\overline{{\mathit{E}}_{\mathrm{SN}}}$ is the mean energy of a SN, V_esc is the escape velocity, and ε is the efficiency of the SN feedback. Thus, M_out = dM_out/dt × ΔT, where ΔT is the duration of the starburst. Assuming that $\overline{{\mathit{E}}_{\mathrm{SN}}}\mathrm{=}\mathrm{1.2}\mathrm{\times }{\mathrm{10}}^{\mathrm{51}}$ erg, ${\mathit{V}}_{\mathrm{esc}}\mathrm{=}\sqrt{\mathrm{2}}\mathrm{\times }{\mathit{V}}_{\mathrm{rot}}$, ε = 0.15, ΔT = 20 Myr (Aloisi et al. 1999), and SNR = 0.01 × SFR with SFR = 0.06 M_⊙ yr^-1 (Aloisi et al. 1999), we find that M_out ~ 5−8 × 10⁷ M_⊙. The total mass of the extended gas is  ~1.6 × 10⁸ M_⊙ (corrected for the presence of He). Thus, a massive outflow may explain all the diffuse gas if only slightly higher values of SFR, ΔT, and ε are assumed. Since the extended H i  emission is entirely at approaching velocities, any outflow should be highly asymmetric and confined.

Finally, we would also like to mention the hypothesis of a “fragmenting H i  cloud in the early stages of galaxy evolution”, which was suggested by van Zee et al. (1998b). This picture can explain the extremely low metallicity of I Zw 18, but is in contrast to the results of Aloisi et al. (2007), who concluded that I Zw 18 has old stars and is not a young galaxy in formation. Alternatively, I Zw 18 A and I Zw 18 C may be old stellar systems that are accreting cold gas from the inter-galactic medium and are now forming new stars. This may be in line with some simulations of dwarf galaxy formation (e.g. Kereš et al. 2005; Dekel & Birnboim 2006).

6.2. Comparison with other dwarf galaxies

The evolution of BCDs is still not understood. In particular, it is unclear what objects can be identified as their progenitors and descendants (Papaderos et al. 1996; van Zee et al. 2001). It is useful, therefore, to compare their properties with those of other types of dwarf galaxies.

In Fig. 10, we compare I Zw 18 A with a typical dwarf irregular (UGC 7232) taken from the sample of Swaters et al. (2009). The H i  observations of these two objects have almost the same linear resolution (~200 pc), making it possible to compare H i  surface densities and velocity gradients. The two galaxies have approximately the same H i  size and the same rotation velocity at the last measured point, thus they have roughly the same dynamical mass. However, their structural properties are very different:

• 1.

  The azimuthally averaged H i  surface densities of I Zw 18 A are a factor  ~2 higher than those of UGC 7232 (Fig 10, top). In addition, the H i  distribution of I Zw 18 A is clumpy with H i  column densities as high as  ~50–100 M_⊙ pc^-2 (Fig. 3).

• 2.

  The stellar component of I Zw 18 A is much more compact than that of UGC 7232 (Fig. 10, middle).

• 3.

  The rotation curves have a completely different shape (Fig. 10, bottom). UGC 7232 has a slowly rising rotation curve indicating a smooth mass distribution, whereas I Zw 18 A has a flat rotation curve with an inner steep rise, indicating that there is a strong concentration of mass.

In Sect. 4.2, we showed that the central mass concentration cannot be explained by the newly formed stars and/or by the concentration of HI, but can be identified with the old stars and/or dark matter. As to the molecules, their amount is highly uncertain and it is unclear whether they are dynamically important or not. The now uncovered concentration of mass is unique among dwarf irregular galaxies and must be tightly linked to the starburst. Furthermore, this result sheds new light on the question of the evolution of BCDs and their descendants. It is clear that, unless a significant re-distribution of mass takes place, a steeply rising rotation curve, as found for I Zw 18, would be the distinctive signature that would make the descendants recognizable.

For other BCDs, it is known that their underlying stellar component, which consists of old stars, is generally more compact than common dEs and dIrrs (e.g. Papaderos et al. 1996; Gil de Paz & Madore 2005). In particular, the old stellar component of the majority of BCDs has a typical central surface brightness μ₀ ~ 21 mag arcsec^-2 in the B-band (e.g. Gil de Paz & Madore 2005), that is similar to those of high surface brightness (HSB) disk galaxies (e.g. van der Kruit & Freeman 2011). If the distribution of mass is strongly coupled to the distribution of light (Sancisi 2004), we expect BCDs to have a similar dynamical behaviour of HSB spiral galaxies, i.e. steeply rising rotation curves that can be described under the maximum disk hypothesis. This seems to be the case for I Zw 18 A. In all these respects, I Zw 18 A resembles a “miniature” HSB disk galaxy. There are already indications that BCDs may have “steeper rotation curves than similar luminosity, low surface brightness dwarf galaxies” (van Zee et al. 2001), but a detailed dynamical study is needed to derive reliable rotation curves and determine the relative contributions of gas, stars, and dark matter to the gravitational potential.

7. Conclusions

We have analysed H i  observations of the blue compact dwarf galaxy I Zw 18. Our main results can be summarized as follows:

• The H i  gas associated with the starburst region (I Zw 18 A) is in a compact rotating disk. The H i  column densities are very high, up to  ~50–100   M_⊙ pc^-2 (~0.6–1.2    ×    10²² atoms cm^-2).

• The disk has a flat rotation curve with an inner steep rise. This indicates that there is a strong concentration of mass that may be either luminous or dark. Baryons may dominate the gravitational potential in the inner regions.

• The disk has a radial inflow/outflow motion of  ~15 km s^-1.

• The stellar concentration to the north-west (I Zw 18 C) is surrounded by extended H i  emission, which is smoothly connected with I Zw 18 A.

• An H i  tail extends to the south of I Zw 18 A out to  ~13.5 kpc. It has a coherent kinematical structure and seems to be connected to the H i  emission to the north-west.

I Zw 18 A appears structurally different from a typical dIrr in terms of H i  distribution, stellar distribution, and dynamics. In particular, it has a strong central concentration of mass. It may be considered as a “miniature” HSB disk galaxy. The H i  concentration and the dynamical properties must be tightly linked to the starburst. They are also crucial to address the question of the progenitors/descendants of BCDs.

The mechanism that triggered the starburst is most likely to have been an interaction/merger between gas-rich dwarf galaxies.

---

Acknowledgments

We thank F. Annibali, G. Fiorentino, and M. Tosi for helpful discussions about the stellar populations of BCDs and for providing the HST image of I Zw 18. We are grateful to J. van Gorkom for stimulating discussions. We also thank J. Cannon for kindly making his HST images available to us.

References

All Tables

All Figures