© ESO, 2015

1. Introduction

CW Leo (also known as IRC +10216) is the nearest carbon-rich asymptotic giant branch (AGB) star at a distance of about 120 to 150 pc (Men’shchikov et al. 2001; Groenewegen et al. 2012, and reference therein). During the AGB phase, the evolution is determined by mass loss, with the mass-loss rate being significantly higher than the rate at which mass is comsumed by nuclear burning (Lagadec et al. 2008). The generally accepted idea is that mass loss proceeds in a two-step process: (1) pulsations lift photospheric layers to temperatures that are cool and dense enough for molecules to nucleate and form dust species; (2) radiation pressure on the dust grains pushes the material outwards, and a wind is generated. However, even for the nearest carbon-rich AGB star it is still not well established whether the dominant behaviour of the mass loss occurs by homogeneous, isotropic processes, or if small-scale irregularities or instabilities dictate the overall appearance.

Specifically, for CW Leo it has been shown that at milli-arcsecond scale dust clumps that formed in the dust formation region are moving as a result of the expanding wind of the star, with clump velocity motions varying between 7.9 and 17.5 km s^-1 (Tuthill et al. 2000; Weigelt et al. 2002; Menut et al. 2007). Formation of new dust clumps at irregular time intervals might explain the fading illumination of some dust features during some periods. At arcsecond-scale, a bipolar morphology is seen in optical data with a bipolar axis at a position angle estimated at about 8° (Mauron & Huggins 2000) to 22° (Leão et al. 2006; with respect to north). This feature probably indicates that scattering from stellar photons is more efficient in the polar direction. Evidence was provided from optical and near-infrared data for an equatorial density enhancement in the form of a geometrically thick disk (Dyck et al. 1987) or dusty torus (Skinner et al. 1998; Murakawa et al. 2005; Jeffers et al. 2014) seen almost edge-on. Beyond 1′′, multiple almost concentric shells (or arcs) are detected on top of the smooth extended envelope (Mauron & Huggins 1999, 2000; Leão et al. 2006; Dinh-V-Trung & Lim 2008; Fong et al. 2003; Decin et al. 2011). The nested shells are composed of thinner elongated arcs with a typical width of ~1.6′′ and separated by ~3′′–20′′, corresponding to ~150–1000 yr. These partial shells have a typical angular extent of ~30° to 90° and are detected up to a distance of 320′′ (Decin et al. 2011). Around 450′′, Herschel data have revealed a bow shock (Ladjal et al. 2010).

In this paper, we present the first ALMA band 9 observations of the inner envelope of CW Leo at a spatial resolution of ~0.2′′ (~10 R_⋆, with 1 R_⋆~ 5 × 10¹³ cm; De Beck et al. 2012). In Sect. 2 we describe the ALMA observations and data reduction including imaging and image fidelity. Section 3 reports on the astrometry and continuum properties. The ALMA spectra are described in Sect. 4, and the channel maps and PV-diagrams are shown in Sect. 5. A qualitative interpretation of the ALMA ¹³CO J = 6–5 channel maps is given in Sect. 6. Preliminary morpho-kinematical modelling of the ALMA data is presented in Sect. 7. The results are discussed in Sect. 8, and the conclusions are given in Sect. 9.

2. ALMA observations and data reduction

2.1. ALMA observations

IRC +10216 was observed by ALMA on 1 December 2012 for proposal code 2011.0.00277.S (Cycle 0). The total time on IRC +10216 was 17.5 min spread over 35 min. The precipitable water vapour was 0.55 mm. The 18 antennas provided good data, with minimum and maximum baselines of 25–340 m.

Four 1.875-GHz frequency bands (spectral windows, spw) were used for science, centred at 643.6645, 646.4827, 658.0587, and 661.1205 GHz (ALMA Band 9). Each spectral window was divided into 3840 channels of 488 kHz. After Hanning smoothing, the effective velocity resolution is 0.455 km s^-1. The usual ALMA instrumental corrections were applied, including those derived from water vapour radiometry and system temperature measurements. The Cycle 0 flux scale was derived with respect to Titan using a resolved, frequency-dependent model derived from its ephemeris¹, as appropriate for the baseline lengths whose flux density was greater than the noise (maximum ~6.535 Jy). The flux density of the QSO 3C 279 was derived with respect to Titan and found to be (9.39, 9.28, 8.82, 8.70) Jy in spw 0–3. QSO 3C 279 is point-like on all baselines and was used in turn to scale the antenna gains for all other sources. The flux density scale is formally accurate to 20% in band 9; the accuracy of our data might be even better because the flux densities derived for other calibration sources were consistent with values measured during other observations. The compact, bright QSO 3C 279 and J0854+201 were used as bandpass calibration and phase reference sources, respectively; both were also used to calibrate the pointing.

2.2. Data reduction

Normal data reduction procedures were performed in CASA (McMullin et al. 2007)². Inspection of the phase reference source data showed that the phase of the first few scans could be connected smoothly, without wrapping, but there were apparent discontinuities in the second half of the data. We applied the phase-reference solutions to the target data and inspected the first half of the corrected data to identify the line-free channels. We used them to make a continuum image, which provided a model for self-calibration of all the target continuum data. The residual scatter in the final self-calibration solutions was about 5° in phase and 10% in amplitude. Using the expressions from Perley (1999), this corresponds to a dynamic range limitation of ~150–300, depending on whether the errors are correlated between scans. We applied the solutions to all the channels of the target data. Finally, we made a preliminary cube to refine the selection of continuum channels, subtracted them, and made the final line and continuum images.

We used the software CASA for the Fourier transformation and to clean each channel to produce image cubes for all spectral windows. During this process, we corrected for the motion of Earth relative to the local standard of rest. For the identified lines, we extracted the channels covering that line and set the velocity to zero at its rest frequency. The noise in the image cubes is 0.2–0.3 Jy beam^-1 per 488 kHz channel³. The noise is highest near the 658 GHz atmospheric water line and is also higher in bright channels because of the dynamic range limitations and because it is slightly increased by continuum subtraction and residual calibration errors.

The short duration of the observations provided sparse visibility plane coverage, most noticeable at the highest spectral resolution. The full-resolution images were made with partial uniform weighting to give a synthesised beam 0.̋42  ×  0.̋24 for a position angle of 117° (measured anticlockwise from north). This is similar to the natural beam position angle, but improves the resolution without increasing the noise appreciably⁴. A beam of 0.̋4 × 0.̋4 was used for the three-channel averaged images, which have a resolution of ~0.7 km s^-1 and a noise value of 0.13–0.18 Jy beam^-1 (see Sect. 5.1).

The shortest baselines in our observations impose a limit of ~3′′ on the size of any structure that can be imaged accurately. Emission on scales 3′′–6′′ is imperfectly sampled, and emission that is smoothly extended over >6′′ is undetectable, even if it is much brighter than our limit (~1 Jy beam^-1 at full resolution) for compact emission.

These limitations (dynamic range and sparse visibility plane coverage) produce two kinds of artefacts. The dynamic range limitations are apparent as residual sidelobes of the synthesised beam: the intensity is proportional to the channel peak and is placed symmetrically w.r.t. the channel peak at a constant spacing (scaled with frequency). We took care to check that any apparently symmetric structures were not in fact co-incident with sidelobes. The effects of poorly sampled emission are more complicated because they are the convolution of the actual distribution with the sidelobes. If the true distribution would for example be a disk >6′′ across in a given channel, the inner ~3′′ would be faithfully imaged, but at larger separations negative and positive ring-like artefacts (known as the “cereal-bowl” effect) might be visible. Nevertheless, smaller and larger rings (or patches of emission) seen in the channel maps with a size in the 3′′–6′′ range result from actual structures with a smaller or larger size, respectively, which allows one to at least partly reconstruct the actual morphological structures.

2.3. Position measurements

We measured the positions and extents of emission by fitting 2D Gaussian components to the continuum peak and to the zeroth moment of each line (see Sect. 3). The noise-based uncertainty (for a sparsely filled snapshot) is given by the beam size divided by the signal-to-noise ratio (S/N; before self-calibration), and gives the uncertainty in comparing positions within the same data set. However, if the actual flux distribution is non-Gaussian, the uncertainties are larger.

3. Imaging results

In this section, we discuss the astrometry and stellar and dust properties of CW Leo as deduced from the ALMA images in band 9.

3.1. Astrometry

The position of the 650 GHz continuum peak is 09:47:57.4553 +13:16:43.749 (J2000), epoch 2012.92. Four factors affect the astrometry. Uncertainties in the antenna positions and the phase reference source coordinates are smaller than a few mas. The noise-based uncertainty (for a sparsely filled snapshot) is given by the beam size divided by the signal-to-noise ratio (before self-calibration), ~ for the continuum peak. The uncertainty in measuring the extent of emission is about twice the position uncertainty. The phase reference source is 14° from IRC +10216, and the timescales of phase changes due to the atmosphere lead to the dominant uncertainty of . Menten et al. (2012) used the VLA in 2006.16 to measure the position and proper motion of IRC +10216, which would give 09:47:57.4417 +13:16:43.896 in 2012.92 (uncertainty ). The discrepancy from the position measured by ALMA is (, ), which is within the uncertainties.

The position fitted to the zeroth moment of each line apart from ¹³CO is within 30 mas (the combined position uncertainties) of the continuum peak; the mean position of all these lines is within 6 mas. The ¹³CO peak is offset by (−67, −37) mas, but this is probably due to the relatively poor fit of a Gaussian component to the complex and extended CO emission (see Sect. 5.1).

3.2. Stellar and dust properties

The total flux density within the 3σ contour is 5.66 Jy, see Fig. 1. The continuum peak is 4.362 Jy beam^-1 with a noise uncertainty of 0.023 Jy beam^-1 and a calibration uncertainty <0.87 Jy beam^-1. The size (deconvolved from the beam) is 120 ± 4 mas by 108 ± 5 mas at a position angle (PA) of 128 ± 20°. The ALMA total flux density at 650 GHz is slightly higher, although within the uncertainty limits, than the continuum emission flux density derived by Young et al. (2004) using the SMA at 680 GHz, who found a compact unresolved component in their ~2′′ full width at half maximum (FWHM) beam with a flux density of 3.9 ± 1.2 Jy.

The stellar photospheric emission contributes part of the total flux density at 650 GHz. Following Jenness et al. (2002), the ALMA data were taken near maximum light (at phase ϕ = 0.98), which would imply a stellar luminosity of ~15 600 L_⊙ according to De Beck et al. (2012). For an effective temperature of 2330 K (De Beck et al. 2012), this yields a stellar diameter of 48 mas. Groenewegen et al. (1997) reported a stellar diameter of 70.2 mas at 243 GHz for an effective temperature of T_eff of 2000 K. Menten et al. (2012) found a diameter of 83 mas at 43 GHz. Thus, the star is probably unresolved at our spatial resolution.

We used two approaches to estimate the stellar flux contribution at 650 GHz. First, using the IRAM Plateau de Bure Interferometer (PdBI), Lucas & Guélin (1999) found a “point source” at 89 GHz (3.5 mm) and 242 GHz (1.3 mm) with flux densities of 65 ± 7 mJy and 487 ± 70 mJy, respectively. They identified the point source as being the photospheric emission. Using a spectral index of 1.96 (±0.04) as derived by Menten et al. (2006), this yields a stellar contribution of 3.3 Jy at 650 GHz. The uncertainty on this value is estimated to be about 20% and arises from the variability of the pulsating star and the uncertainties on the spectral index. Second, the stellar contribution at 650 GHz can also be calculated from the effective temperature and stellar diameter, under the assumption of a black-body spectrum. Using the values obtained by De Beck et al. (2012) at maximum luminosity (see previous paragraph), the stellar flux at 650 GHz is 1.3 Jy, while the values derived by Groenewegen et al. (1997) yield a stellar flux of 2.3 Jy. The difference in stellar flux estimates mainly arises from the uncertainties in the angular diameter; the uncertainty on the effective temperature only gives an uncertainty in the flux of 0.2 Jy. These values for the stellar contribution are a factor ~1.4–2.5 lower than the value obtained by extrapolating the “point-source” values of Lucas & Guélin (1999). This might mean that part of the “point-source” flux densities at 89 GHz (3.5 mm) and 243 GHz (1.3 mm) might not come from the central star, resulting in too high an estimate for the stellar contribution at 650 GHz when using the data of Lucas & Guélin (1999). Following Reid & Menten (1997), Groenewegen (1997), Groenewegen et al. (1997), we suggest that free-free emission at wavelengths beyond 800 μm might (at least partly) contribute to the high flux density values of Lucas & Guélin (1999).

The estimated stellar flux density at 650 GHz of ~2.3 Jy (adopting the parameters of Groenewegen et al. 1997) is significantly lower than the flux density of 5.66 Jy obtained from the ALMA data. Since free-free emission is negligible at 650 GHz (Groenewegen 1997), we conclude that the extended emission seen in the ALMA data is due to dust emission alone.

We estimated the 650 GHz emission from dust by subtracting a 2.3 Jy point source, convolved with the natural synthesised beam, from the continuum image, leaving a residual 3.36 Jy (see Fig. 2). The main sources of uncertainty are the 20% uncertainty on the ALMA data and the estimated stellar flux density. The residual emission is extended and clearly shows an asymmetric distribution, with a central peak flux density of ~2 Jy. The lowest three-sigma contour in Fig. 2 corresponds to a position uncertainty ~1/3 of the beam size, so although the apparent southern bifurcation may be unreliable, the overall extension is real. The inner emission appears to be elongated in a direction similar to the ~128° position angle of the equatorial density enhancement detected in scattered light and polarimetry data (Skinner et al. 1998; Murakawa et al. 2005; Jeffers et al. 2014). However, this is close to the direction of elongation of the natural synthesised beam and may be coincidence. The faintest emission is, nonetheless, elongated in an almost orthogonal direction, at a position angle of ~20° and 200°, corresponding to the biconical cone detected by Skinner et al. (1998), Mauron & Huggins (2000; see Sect. 6).

The dust nucleation zone starts around 5 R_⋆ (or 100–150 mas, Decin et al. 2010b). This is similar to the resolution of these images, but we cannot distinguish faint details this close to the strong stellar peak. The outer radius of the dust detectable by ALMA is set by the radius at which it becomes too diffuse or too cool to be detectable. Subtracting any realistic estimate for the stellar contribution (1.3–3.3 Jy) leaves an approximately flat or centrally peaked flux distribution within the inner ~160 mas radius, which means that the dust could be close to optically thick. If the stellar contribution were ≳4 Jy, the dust emission would be characterized by a central hollow.

Using the continuum image after subtracting a 2.3 Jy star, assuming that the dust is optically thin and for an opacity κ_λ at 650 GHz of 90 cm²/g (Mennella et al. 1998) and a dust temperature distributed as derived by Decin et al. (2010a), we obtain a dust mass within 0.8′′ of 3–6 × 10^-7M_⊙, corresponding to a dust mass-loss rate of 1.1–2.2 × 10^-8M_⊙/yr. This estimate for the dust mass-loss rate is a factor ~4 lower than derived by Decin et al. (2010a) from modelling the spectral energy distribution. The large error margin is due to the flux scale uncertainty, the unknown stellar contribution, and the unknown composition of the dust and geometry in the inner regions. Assuming a constant dust velocity equal to the terminal velocity of 14.5 km s^-1, this yields a dust optical depth of ~0.1. Since the dust velocity might be overestimated in the wind acceleration region, the derived value for the dust optical depth might be an underestimate by a factor of a few. Future multi-wavelength high-resolution ALMA observations will provide a spectral index map and thus a better measurement of the optical depth.

4. Spectral line results

In this section, we present the ALMA spectra of CW Leo in band 9 around 650 GHz (Sect. 4.1). For each detected emission line, we determine the line strength, line width, and spatial extension. The line strengths are compared with the line strengths as deduced from Herschel/HIFI observations in the same frequency window (Sect. 4.2). The line widths and spatial extension are used to determine the velocity structure of CW Leo in the inner wind region, where the wind is accelerated from the sound velocity to the terminal velocity (Sect. 4.3).

4.1. ALMA spectra of CW Leo in band 9

Figure 3 shows the spectra of the four spectral bands extracted using a circular aperture of 1.2′′ centred on the dust continuum peak. In the ALMA data, 25 molecular line transitions are detected, of which 7 are unidentified (see Table 1). The identified emission lines arise from ¹³CO, C¹⁸O, SiO, ²⁹SiO, SiS, ²⁹SiS, ³⁰SiS, Si³³S, Si³⁴S, ²⁹Si³⁴S, SiC₂, ¹³CS, and AlF.

Following De Beck et al. (2010), we applied a minimization routine to numerically fit the profiles of the molecular emission lines and extract some of the physical wind properties. We used the soft-parabola function (Olofsson et al. 1993), which provides the flux density F as a function of velocity v: (1)where F_peak is the peak flux in the line centre, v_LSR the systemic velocity with respect to the local standard of rest, v_line the half width of the profile at continuum level, and γ is a shape parameter. The results are listed in Table 1. The soft-parabola fit to the ALMA data yields a mean difference between the rest frequencies and the ALMA frequencies of −57.4 MHz or −26.41 km s^-1, in accordance with the local standard of rest velocity derived by De Beck et al. (2012) from the high-spectral resolution Herschel/HIFI data. The value of the half line width, v_line, measures the gas velocity at the line formation region and will equal the terminal wind velocity (v_∞) for lines being formed in the region where the wind has achieved its full (terminal) velocity. The half-line width is determined by the expansion velocity, v_exp(r), and turbulent velocity, v_turb(r), with r being the radial distance. The turbulent velocity in stellar winds is not well established. Keady et al. (1988) proposed a turbulent velocity in the inner wind region of CW Leo of ~5 km s^-1 and a terminal turbulent velocity of ~1.0 km s^-1; De Beck et al. (2012) derived a value of 1.5 km s^-1 from a fit to high-resolution Herschel/HIFI CO lines. The derived values for the half-line width range between ~3.5 km s^-1 and ~15 km s^-1. The velocity structure as derived from the ALMA data is discussed in detail in Sect. 4.3.

To determine the spatial extension, we computed the first moment for the inner 10 km s^-1 of each line and measured the spatial FWHM by fitting a 2D Gaussian component and deconvolving the restoring beam (see Table 1). We only performed the spatial fitting to the inner 10 km s^-1 to avoid line blending and since for an expanding shell the emission at line centre would be most spatially extended, whilst the extreme velocity emission comes from closer to the line of sight to the star. All lines have a high enough signal-to-noise ratio to be spatially resolved, except for the faintest unidentified line (“U” line) at 645.507 GHz. They are all compact enough for a good fit except for the ¹³CO line, for which its apparent spatial FWHM, s_FWHM, represents only the inner peak. There is likely also to be extended C¹⁸O emission below our detection threshold. Table 1 also lists the (not deconvolved) largest angular size (LAS), that is, the maximum width of the area enclosed by the central contour at 2 Jy km s^-1 per beam (approximately 3σ_rms). For ¹³CO J = 6–5 there is also non-contiguous emission extended throughout the field of view (see Sect. 5.1). As illustrated in Fig. 4, the extension scales inversely with the upper state energy, E_upper. The anomalies in this plot are due to blending and the fact that the ALMA data resolve out a lot of the ¹³CO flux (see Sect. 4.2).

4.2. Comparison between ALMA and Herschel

The frequency range observed with ALMA was also observed by the Herschel/HIFI instrument (de Graauw et al. 2010). A survey in all HIFI bands was carried out on May 11–15 2010. The data reduction and first results are described in Cernicharo et al. (2010). The HIFI beam size around 640 GHz is 33.2′′ and the main-beam antenna efficiency is 0.75 (Roelfsema et al. 2012). The data were taken in double beam-switching mode at a frequency resolution of 1.1 MHz. To increase the S/N of the weaker lines, some spectral regions were rebinned. The rms noise after averaging all scans is given in Table 2. To convert from antenna temperatures (in K) to fluxes in Jy, we used the point-source sensitivity of 466 Jy/K (Roelfsema et al. 2012). The source size of all transitions listed in Table 2 can be approximated with a point source within the Herschel beam, with the exception of the ¹³CO J = 6–5 line. Based on observations with the JCMT and CSO (Crosas & Menten 1997) and our modelling of the ¹³CO emission presented in Decin et al. (2010b) and De Beck et al. (2012), we derive a source size of ~10′′ for the ¹³CO J = 6–5 line, yielding a 3% correction w.r.t. point-source sensitivity.

Table 2 lists the integrated line intensities measured with Herschel and ALMA. Note that the ALMA integrated line intensities are higher than the values listed in Table 1 because the spatial fitting used for the values in Table 1 was restricted to the inner 10 km s^-1 (see previous section). As can be seen from Table 2, only a few lines are not plagued by blends with other molecular line transitions. The half-line width at continuum level, v_line, of all common unblended lines agrees for both instruments. The integrated flux density of the unblended lines measured from the Herschel data is typically 1.2–1.7 times higher than that measured using ALMA (note that the overall ALMA flux scale is uncertain by up to 20%). The ratio is lower for the Si³³S J = 36–35 transition, but this line has a low signal-to-noise ratio in the Herschel data, and the ratio is much higher for ¹³CO J = 6–5 line. The 20–50% additional flux detected by Herschel compared with ALMA can be explained by two factors.

Firstly, there are observational differences. Herschel is fundamentally more sensitive to large-scale emission. Its 5σ sensitivity for these observations, per 0.7 km s^-1, was ~20 Jy beam^-1, compared with <0.2 Jy beam^-1 for ALMA. However, the Herschel beam is ~10⁴ times the area of the ALMA synthesised beam. A 20 Jy point source would be very bright in the ALMA data, but would be spread out over an area of more than 33.̋2 in diameter, the surface brightness per ALMA beam would be below the detection limit. In addition, even brighter emission that is smooth on scales >3′′ would be poorly imaged by ALMA, and not detected at all, however bright, if it were >6′′, as explained in Sect. 2.2. Our results are consistent with this because the fraction of the Herschel flux detected by ALMA is higher for higher excitation temperature lines and lower for more spatially extended species, and we detect a similar fraction of the flux with ALMA for similar transitions. Specifically, for the transitions listed in Table 2, the estimated emission region for the SiS isotopologues (J = 37–36 or J = 36–35) in the ground-vibrational state is <6.̋5 (based on modelling by Decin et al. 2010b). For the ¹³CS J = 14–13 line the main emission region is within ~7′′ (based on the ¹²CS J = 7–6 and ¹³CS J = 7–6 line observed by Williams & White 1992; Patel et al. 2011), and the ¹²CS J = 14–13 SMA data of Young et al. (2004) even indicate an emission region within ~4′′. As mentioned above, the main emission region for the ¹³CO J = 6–5 transition is ~10′′ and the ALMA observations have resolved out at least 95% of the flux. No flux will be resolved out for the high-excitation SiS v = 2J = 36–35 line with an upper state energy of 2724 K; the difference between ALMA and Herschel integrated line flux values is within the calibration uncertainties of both instruments.

Secondly, the variability of thermal line emission due to the variability of the central star is different. Following Jenness et al. (2002), the Herschel and ALMA data were taken at phase 0.46 and 0.98, respectively. The changes in the radiation field during the pulsation period can influence the molecular excitation. The infrared radiation field is expected to mainly affect molecules with strong vibrational transitions, with stronger modulations for lines excited close to the star. Using Herschel data, Teyssier et al. (2013) have recently shown that the ¹²CS J = 14–13 has a maximum amplification factor in the line intensities of ~1.3, while for SiS J = 37–36 and J = 36–35 this is ~1.2. The low-excitation CO lines seem quite insensitive to a change in the infrared radiation field.

Predicting whether a line will have a higher or lower amplification factor depending on the variability phase at the time of the observations is complex, however. The ALMA data are taken near light maximum and the Herschel data near light minimum, which indicates that higher excitation levels might be more populated at the time of the ALMA observations. However, the exact effect is dependent on the transition-rate probabilities with the connecting levels.

4.3. Gas velocity structure in the inner wind region

Table 1 and Fig. 5 show that the wind velocity of the lines falls into two groups. All but two of the identified lines trace a velocity higher than 13.5 km s^-1 at a LAS radius (representing the not-deconvolved largest detectable scale) beyond 600 mas radius (or 30 R_⋆) and a spatial FWHM (illustrating the main line formation region) beyond 200 mas (10 R_⋆), indicating that the wind has reached its terminal velocity (~14.5 km s^-1) within ~10 R_⋆ from the central star, in agreement with the results of Keady et al. (1988). The other group has line width velocities <8 km s^-1 within 800 mas LAS radius (40 R_⋆). This contains the two identified lines in v = 2 excited vibrational states and all unblended “U” lines, suggesting that these also come from high-excitation states. The narrowest line is the “U” line at 646.31 GHz with a velocity of only 3.6 km s^-1, slightly higher than the sound velocity of ~3 km s^-1 at the base of the dust nucleation region.

Plotting the velocities as a function of the spatial FWHM (see Fig. 5) shows that the measured velocities rapidly increase from the sound velocity to ~8 km s^-1 between ~100 to 200 mas (5–10 R_⋆). Assuming that the turbulent velocity does not vary substantially in this small region, this indicates that carbonaceous material condenses around 5 R_⋆, resulting in a rapid acceleration of the gas until a velocity of ~8 km s^-1. In the four ALMA spectral windows, no lines are detected with a velocity between 8 and 13.5 km s^-1, but the data suggest a rapid increase of the gas velocity from ~8 km s^-1 to ~14.5 km s^-1 around ~10–12 R_⋆. This acceleration profile is not captured by the solution of the wind momentum equation derived by Decin et al. (2010a), who assumed instantaneous condensation of amorphous carbon, iron, silicon carbide, and magnesium sulfide dust species at a radius of 5.6 R_⋆. Keady et al. (1988) suggested a two-step acceleration process from an analysis of high-resolution near-infrared ro-vibrational CO spectra. They suggested an increase around 3 and 11 R_⋆, with the first increase up to a velocity of 11 km s^-1 being due to the condensation of carbon and a second increase that might be linked to dust condensation involving Mg and/or S. This idea was used and adapted by Fonfría et al. (2008), who modelled a large set of C₂H₂ lines. They suggested a strong velocity gradient from 5 to 11 km s^-1 around 4.7–5.7 R_⋆ and a second velocity increase to 14.5 km s^-1 around 20.75–21.75 R_⋆ (see Fig. 5). Based on the ALMA band 9 data, the wind acceleration profile can be further refined. The data can be fitted assuming a linear velocity increase from around 2–3 km s^-1 at ~5 R_⋆ reaching terminal velocity around 11 R_⋆ (see the full grey line in the right panel of Fig. 5). Assuming a two-step acceleration process, we obtain (2)that is, a first jump around 5.6 R_⋆ and a second one around 10 R_⋆ (see the dotted grey line in the right panel of Fig. 5). To improve upon the solution of the momentum equation as derived by Decin et al. (2010a), different possibilities exist. (1) Either the dust absorption efficiencies need to increase in the near-infrared, where the bulk of the stellar photons are emitted. (2) Or larger grains (≳0.2 μm) should be formed more efficiently close to the star, thus enhancing the wind acceleration via scattering processes (Höfner 2008). (3) A third possibility might exist by adapting the velocity at the start of the dust condensation region, which is now assumed to be the sound velocity (Decin et al. 2006). However, in that case, the dust-to-gas ratio would need to be decreased to avoid overshooting the observed terminal velocity.

5. Channel maps and PV diagrams

In Sects. 5.1 and 5.2, we show and discuss the channel maps and position-velocity (PV) diagrams of one representative high-excitation line with an upper state energy (E_upper) around 2700 K, of one medium-excitation line (E_upper ~ 450 K), and of one low-excitation line (E_upper ~ 100 K). The vibrationally excited line traces the wind acceleration region, molecules in the inner envelope with excitation energies of about a few hundred K are slightly resolved with the ALMA beam of ~0.2′′, while the ¹³CO J = 6–5 line shows a complex emission pattern.

5.1. Channel maps and morphology

As expected, the channel map of a very high-excitation line of a molecule in the inner wind region just above the photosphere, such as SiS v = 2J = 36–35 (Fig. 6), shows that the emission of this line comes from the very inner regions of the envelope, close to the stellar surface, where the wind material has a velocity of about 4 km s^-1 and the gas kinetic temperature is typically around 900–1300 K. The morphology in the region <1′′ (or diameter of 50 R_⋆) is still more or less circular (for a circular beam of 0.̋4).

The channel maps of molecular transitions formed in the inner wind region and with upper state energies of about a few hundred K, such as the SiC₂ 26(4, 22)–25(4, 21) line displayed in Fig. 7, show that the emission of this line is extended and displays a non-homogeneous distribution of the molecular density. The typical spatial FWHM of these lines is around a few hundred milli-arcsec (see Table 1), and these lines are formed in a region of the wind where the velocity has already reached the terminal velocity.

The ¹³CO J = 6–5 channel map (Fig. 8) displays a complex structure. At different velocities, parts of a spiral-like arm (arcs) can be distinguished: around −17.5 km s^-1 an arc extends to the north at ~1.5′′ offset from the central source, around −24.6 km s^-1 at 0.7′′ offset to the south-east, around −27.5 km s^-1 at 1′′ to the south-west, and around −37.4 km s^-1 at 2′′ to the south-west. Local density enhancements, the change in (projected) velocity for different arcs, gaps in the UV-coverage (see Sect. 2.2) and the fact that we resolve out flux for scales larger than ~3′′ can explain why the spiral arms are not complete. For example, at −37.4 km s^-1, the extended emission has a distinctive curve that is unrelated to the beam sidelobe pattern. The negative artefact is due to the “bowl” effect mentioned in Sect. 2.2, resulting from emission on scales just larger than are sampled by the interferometer. Although the detailed flux distribution at radii ≳ is unreliable, it does represent the general shape of the brighter emission at each velocity.

The local density enhancements or spiral arcs are seen best in the maps of ¹³CO because this is the lowest-excitation and most abundant transition in our ALMA observations, mitigating the over-resolution of its extended emission. Some of the less-abundant, intermediate-excitation lines have a similar LAS for integrated emission in the ALMA data, but are too faint to resolve extended details, spectrally or spatially. Figure 9 compares the integrated ¹³CO J = 6–5 emission (zeroth moment) and dust emission. Bright arcs at 1′′ to the north and 1′′ and 1.6′′ to the south-west can clearly be discerned in the ¹³CO emission, but the dust emission at 650 GHz is below our detection threshold at radii beyond 1′′. This is the first time that we can confirm that the extended continuum and molecular emission are both centred on the continuum peak position.

5.2. Position-velocity diagrams and morphology

Using the casa task impv, position-velocity diagrams were calculated for all transitions. In the standard set-up, 201 slits with a width of 1 pixel (=0.060′′) and a length of 101 pixels (=6.06′′) were taken. Each slit was then centred perpendicular to the direction slice (i.e. in this case 3.03′′ above and 3.03′′ below the central pixel), and the average for each slit was taken. A PV along the right ascension axis hence has PA = 90°, and along the declination axis PA = 0°. The strength of PV plots is that they correlate structure at different spatial offsets, which is otherwise a challenging task from channel maps alone.

As discussed in Sect. 2.2, the PV diagrams are not free from artefacts. (1) For lines with an integrated line peak flux higher than 0.4 Jy/beam, sidelobe effects are clearly visible at −4′′ and +4′′ offset (see right panel in Fig. 10). (2) The instrument resolves out some flux for most of the lines, which is especially true for the emission of the CO isotopologues (see Sect. 4.2). In some PV diagrams, emission extends beyond the −10 to −40 km s^-1 range (covering the radially outflowing wind), which is due to blending with molecular line transitions at slightly different frequencies (see right panel in Fig. 10).

For all transitions, the channel maps and PV diagrams indicate that (some part of) the emission originates from the region just above the stellar photosphere, that is, within ~0.2′′ or 10 R_⋆ (see Appendix A). The ALMA data prove that SiO is formed in the inner wind region, indicating an active photochemical region where shock-driven chemistry (Cherchneff 2011) and/or photodissociation of molecules in a clumpy medium (Agúndez et al. 2010) dictate the chemical balance.

Figure 11 shows the PV diagrams of the ¹³CO J = 6–5 channel maps for a direction slice in right ascension or in declination. These PV diagrams display a more complex morpho-kinematical structure than the PV diagrams in Fig. 10, or any other line. This is expected because the ¹³CO J = 6–5 line is the strongest in the frequency range we covered, as discussed in Sect 5.1. It is clear that the curved distribution of the ¹³CO PV maps represents the genuine morpho-kinematical shape, but that the actual distribution of the emission on scales larger than ~3′′ is artificially fragmented as a result of the gaps in UV spacing (see Sect. 2.2). By changing the position angle of the PV diagram, changes in the morpho-kinematical structure can be diagnosed. The most revealing ¹³CO J = 6–5 PV-diagram is obtained for position angles around 10–25° to the north-east (see Fig. 11). To show and enhance the effect of the small-scale structure, we also calculated these PV diagrams for a slit of 51 pixels (instead of 101 pixels in the standard setup); see bottom panels in Fig. 12. In all panels of Fig. 12, one can clearly distinguish correlated curved structures ranging from −40 to −12 km s^-1 and the distribution of the brightest emission around 0′′ offset displays an “S”-shape.

6. Qualitative interpretation of the ¹³CO J = 6–5 ALMA data

The ¹³CO J = 6–5 PV diagrams in Fig. 12 show a very characteristic shape with correlated structure at an interval of ~1′′ (around the systemic velocity) for position angles of about 10–20°. These structures are almost absent from the PV diagrams for PA ~ 112.5°. In almost all PV diagrams a typical “S”-shape is seen around zero offset. Based on the simulations of Kim et al. (2013), we postulate that this type of morphology implies a spiral structure induced by a binary companion. In contrast to R Scl, for which ALMA data have recently shown a spiral seen almost face-on (Maercker et al. 2012), ALMA data of CW Leo are reminiscent of a spiral seen almost edge-on, with the orbital axis at an angle of ~22.5° north-east (see Sect. 8.1). Our ALMA images capture fragments of the spiral arms, seen almost edge-on, with a typical width of ~300 mas. We cannot be sure whether the fragmentation is an instrumental artefact, but for CIT 6 the 0.̋7-resolution data analysed by Kim et al. (2013) show that the structure in the spiral arms is rough. High spatial resolution observations (Tuthill et al. 2000; Weigelt et al. 2002; Menut et al. 2007) and hydrodynamical models (Woitke 2006) have shown that gas and dust clumps are created in the inner envelope by the wind formation mechanism. These density inhomogeneities will remain in the inner spiral-arm structures until they are eventually dispersed.

No direct evidence has yet been found for a binary companion. However, the current ALMA data in combination with a plethora of other observational diagnostics (schematized in Fig. 13) suggest that CW Leo is part of a binary system. (1) Guélin et al. (1993) suggested that CW Leo is a binary since the clumpy shells observed in MgNC, C₄H, and other carbon-chain molecules and radicals, appear off-centre from the continuum source by 2–3′′ east. The cause of such a drift could be an acceleration of the star due to a companion. (2) As mentioned in Sect. 1, the central regions around CW Leo show an axisymmetric peanut- or bipolar-like structure with the major axis lying at a position angle of ~8–22° (see Fig. 13; Le Bertre et al. 1989; Skinner et al. 1998; Mauron & Huggins 2000; Kastner & Weintraub 1994). This is indirect evidence of a possible equatorial density enhancement of dust and gas, with the equatorial plane being located perpendicular to the major axis of the nebula, enhancing scattered stellar light emission through the biconical openings. A direct detection of an equatorial dust lane seen almost edge-on was recently presented by Jeffers et al. (2014), corroborating the results of Murakawa et al. (2005). The inferred radius of the dust lane is ~0.5–1′′. The elongated east-west structure seen in our ALMA continuum images at a PA of ~128° with a length of 1.8′′ across (Figs. 2 and 9) might be the signature of this equatorial density enhancement⁵, while the extended emission at a PA of ~20° and ~200° might reflect the biconical cones. This type of morphology can be explained in terms of a binary system in which angular momentum provides a natural way to define the equatorial plane and polar axis. Interaction between the secondary and expanding red giant primary star will create a bipolar nebula, and the dissipative tidal interaction might eventually lead to a merger of the two stars in a common-envelope system (Morris 1987). (3) Optical and infrared data proved that there are multiple almost concentric shells (or arcs) out to 320′′ from the central star (Mauron & Huggins 1999; Leão et al. 2006; Fong et al. 2003; Decin et al. 2011). The shells are incomplete and cover ~30°–90° in azimuth. These shells might represent the limb-brightened spiral arms. The shell-intershell density contrast deduced from the observations is about a factor 3–10, in agreement with recent hydrodynamical simulations for binary-induced spiral patterns as calculated by Kim & Taam (2012c). (4) Previous observations by Guélin et al. (1993), Lucas & Guélin (1999), and Dinh-V-Trung & Lim (2008) of molecular shells show a lack of molecular emission at position angles of ~0°–30°. This position angle aligns with our inferred model rotational axis of the spiral structure (see Sect. 8.1). As shown by Mastrodemos & Morris (1999), a binary-induced spiral structure has a lower density along the orbital axis, which might result in a lack of molecular emission along the inferred rotational axis.

Mauron & Huggins (2000) argued against the binary scenario based on the irregular arc spacing seen in the optical and infrared images. However, the shells are projections in the plane of the sky of a complex 3D structure. A binary companion in a smooth envelope structure would indeed create a spiral structure with regular arm-spacing. However, the wind and envelope creation does not occur in a homogeneous way, adding complexity to the morpho-kinematical structure. Magnetic cool spots during the active phase (Soker & Clayton 1999) and/or large-scale photospheric convective cells (Freytag & Höfner 2008) can locally reduce the gas temperature, which might locally enhance the dust formation and thus lead to gas and dust clumps. The magnetic cool spots cover only a fraction of the stellar surface, and this inhomogeneity might explain the clumpiness seen on smaller spatial scales and in the arcs⁶. Stellar spots have not yet been observed on AGB stars, although recent observations of linear polarization in CO J = 3–2, SiS J = 19–18 and CS J = 7–6 suggest a complex magnetic field configuration in the wind of CW Leo with a strength of about 50–300 mG at ~3′′ offset of the central source (Girart et al. 2012), and recent interferometric data show strong evidence for inhomogeneities in the molecular photospheric layers of AGB stars (Wittkowski et al. 2011). In addition, modulations in the density structure might be a natural outcome from the dust formation process during which a complex non-linear interplay between gas-grain drift, grain nucleation, radiation pressure, and envelope hydrodynamics create density irregularities. Moreover, flow instabilities (such as Rayleigh-Taylor instabilities) have time to fragment the outward-moving density structures and can produce numerous small-scale cloud-like sub-structures (Woitke 2006).

7. Modelling a binary-induced spiral structure

A full reconstruction of the intensities, the 3D shell pattern and its kinematics using the ALMA ¹³CO channel maps and PV diagrams requires detailed hydrodynamical simulations and radiative transfer calculations. This is beyond the scope of the current paper, but will be presented in a forthcoming paper (Homan et al., in prep.). In this section, we wish to demonstrate how different density fluctuations affect a PV diagram and how the morphology in the current ALMA PV data implies the influence of a binary companion in the shaping of the wind structure around CW Leo. The simulations were obtained using the Shape software developed by Steffen & López (2006) and Steffen et al. (2011).

In this section, we first briefly describe the effect of a binary companion on the wind envelope structure (Sect. 7.1). Then, we present Shape simulations, where we gradually increased the complexity of the envelope structure (Sect. 7.2). In Sect. 8.1, we demonstrate the similarities between the ALMA ¹³CO PV diagrams and PV diagrams based on spiral structures.

7.1. Binary-induced spiral structures

As shown in the seminal work of Mastrodemos & Morris (1999), a binary companion might produce a spiral pattern in the circumstellar wind material. Kim & Taam (2012b) and Kim & Taam (2012c) have studied the separate effects of the orbital motions of the individual stars and have shown that two types of spiral patterns are created. Firstly, a more flattened spiral pattern confined within a very limited height from the orbital plane is created by the companion’s motion (Kim & Taam 2012b). Secondly, a spiral-shell-shaped pattern is created by the orbital motion of the mass-losing star around the centre of gravity (Kim & Taam 2012c). These two spiral patterns are different in several ways. (1) The companion’s wake is attached to the companion, while the pattern due to the motion of the mass-losing star has a stand-off radius defined by the orbital and wind velocity (see Eq. (7) in Kim & Taam 2012c). (2) Secondly, the propagation speeds of the patterns determining the shape are different. However, for AGB stars, this difference is very small since the propagation speed of all patterns is close to the wind speed, which dominates the orbital and sound speed. (3) Thirdly, while the arc pattern due to the reflex motion of the mass-losing star nearly reaches the orbital axis and introduces an oblate-shaped flattening of the circumstellar envelope density, the direct effect of the companion results in a spiral structure confined toward the orbital plane. When combining both types of spiral structures, the hydrodynamical models show clumpy structures within the vertical extension limit of the companion’s wake as a result of shocks.

Kim et al. (2013) presented PV diagrams at different inclinations for simulations including both types of spirals for the carbon star CIT 6. Note that a spherical central region of 2′′ was carved out in their simulations to mimic the central hole seen in the observed HC₃N data. The PV diagrams in their simulations are dominated by the spiral induced by the reflex motion. In Sect. 7.2, we show that for spiral structures quite confined towards the orbital plane and for position angles significantly different from the orbital axis or plane, an S-type feature around offset zero appears in a PV diagram around offset zero when seen almost edge-on; this feature is also seen in the ALMA data (see Fig. 11).

7.2. Morpho-kinematical simulations

To understand the complex morpho-kinematical structure seen in the ¹³CO PV diagrams and link the ALMA data to a 3D shell pattern, we used the Shape modelling tool (Steffen & López 2006; Steffen et al. 2011). This is a flexible interactive 3D morpho-kinematical modelling application for astrophysics that is publicly available. By interactively defining 3D structures, it calculates intensity maps, PV diagrams, channel maps, and spectra. While eulerian 3D grid-based hydrodynamic simulations are possible, in this work we used a purely mathematical description of the object structure (see below). A 3D mesh was constructed that served as a container of the emissivity and velocity field. If the radiation transfer is computed, the information in the mesh is transferred to a regular Cartesian 3D grid, which is used to compute the radiation transfer. The advantage of Shape is that it is computationally very fast, which facilitates a first broad screening of the huge 3D kinematical and morphological parameter space. The deduced model parameters can then be used as input for a more detailed (hydrodynamical) simulation. Recently, a non-local-thermodynamic equilibrium (non-LTE) radiative transfer solver (Shapemol) has been added to compute molecular excitation levels (Santander-García et al. 2012, 2015). This solver is based on the well-known LVG (Large Velocity Gradient, Castor 1970) approximation, which significantly simplifies the radiative transfer problem. The current version of Shapemol still has limitations when calculating accurately the molecular level populations in an AGB wind: it only includes collisional rate constants for temperatures up to 1000 K and H₂ densities between 1 × 10¹³ to 1 × 10⁸ m^-3, rotational levels in the vibrational excited states are not included, and the tabulated densities are sometimes too low for the high densities encountered in the inner wind⁷. Nonetheless, it can still be used for a first morpho-kinematical interpretation of the structure seen in the ALMA PV data of ¹³CO.

Based on the results of Decin et al. (2010b) and De Beck et al. (2012), we chose as input for our basic models a distance of 150 pc, a constant gas mass-loss rate of 1.5 × 10^-5M_⊙/yr, a fractional abundance [CO/H₂] of 6 × 10^-4, an isotopologue ratio ¹²CO/¹³CO of 30, a stellar radius of 4.55 × 10¹³ cm (or 0.02′′ at 150 pc), and a wind velocity of 14.5 km s^-1. The gas kinetic temperature is assumed to follow a power law (3)with T_eff the effective temperature of 2330 K and ζ ≈ 0.5. The molecular line data (energy levels, Einstein A coefficients and collisional rates) were taken from the LAMDA database (Schöier et al. 2005).

We produced models of the inner 6′′ (i.e., total width of 12′′), at ~47 mas resolution for spectral channels separated by 0.4 km s^-1. We simulated the effect of the ALMA setup using the CASA task simobserve, with antenna positions and atmospheric conditions matching those used for our observations. To obtain the same total duration (and range of uv spacings sampled), the simulations had a slightly longer on-target duration. The simulated model images and PV plots therefore provide sampling, including missing spacings, similar to the original observations. The noise in the simulated images is slightly lower, not only because of the longer time on-source, but because although atmospheric and instrumental noise affecting the target are included, the limitations due to other imperfections in phase calibration are not. The plots of simulated data do have similar dynamic range limitations and artefacts due to missing flux (where the input model has a large-scale structure) as the observed images. For example, the absolute value of the most negative emission in the simulated data in the right panels of Fig. 14 is ~30% of the peak. The look-up tables in some plots were curtailed to highlight the most relevant features for comparison. We used the same clean parameters as for our observations to produce image cubes and position-velocity plots.

7.2.1. PV diagram for a radially outflowing wind

We first simulated two simple examples for a radially outflowing wind at a constant expansion velocity of 14.5 km s^-1. In the first case, the mass-loss rate was assumed to be constant (or density ρ ∝ r^-2; see top left panel in Fig. 14); while in the second example (bottom left and right panels in Fig. 14) we demonstrate the effect of density-enhanced shells (Cordiner & Millar 2009; Decin et al. 2011) with a shell-intershell density contrast of a factor 10. The wind region close to the stellar surface is responsible for the bright bar at zero offset in the PV diagrams, while the outer density-enhanced shells create well-distinguished correlated curved structures at larger offset in the PV diagram. Including the ALMA instrumental effects results in a (virtual) fragmentation of the shells (see right panels in Fig. 14). Owing to the short hour-angle coverage of the simulated observations in Fig. 14 (right panels), the simulated PV diagrams at different position angles are not identical⁸.

The simulated PV diagrams for the density-enhanced shells taking the ALMA instrumental effects into account show some resemblance with the ALMA ¹³CO J = 6–5 PV diagram for a position angle around 22° north-east (Fig. 12). However, a spherical symmetric geometry and isotropic density wind structure as applied in the current setup cannot explain the different morphologies seen in the PV diagram with position angles 90° apart from each other and the asymmetry when reflecting around zero-offset (the “S”-shape as seen, e.g., for a PA = 0° and PA = 90° in Fig. 12).

7.2.2. PV diagram for binary-induced spiral wind structure

In this section, we simulate the effect of a spiral structure on the PV diagram. Both the spiral structure due to the companion’s wake and that due to the reflex motion of the mass-losing AGB star are simulated using a simple analytical description. To do so, we followed the results described by Kim & Taam (2012b) and Kim & Taam (2012c) by defining an Archimedean spiral r = a_pθ, with r the radial coordinate, θ the longitudinal coordinate, and a_p the factor of proportionality between the angle θ and the radial coordinate r.

The spiral structure is generated by a thin shell mesh with its position given as the position vector r in spherical coordinates (r,θ,φ) as a mapping of the unit sphere (1, θ, φ), with φ the latitude (φ = 90° being the equator). The coordinates of the mesh vertices are given by(4)with r in units of meter, θ between 0 and 360°, and φ between 0 and 180°. The spacing between two successive turnings in the plane of the sky, a, is proportional to a_p, with a defined in units of arcsec for a distance of 150 pc.

The overall radial density was taken to be ∝r^−p (with p ~ 2, Kim & Taam 2012c), which was then modified by the longitudinal and latitudinal dependencies. To simulate both types of spiral arms (confined toward the orbital plane and extending almost completely toward the orbital axis), we define the latitudinal dependence of the density as(5)in which case low values of w (in units of degrees) generate a spiral arm confined toward the orbital plane. Hence the density in the spiral arms is defined as(6)with(7)and ρ₀ the number density at the stellar radius (r₀) as derived from the equation of mass conservation for a mass-loss rate of 1.5 × 10^-5 M_⊙/yr. The thickness β of the spiral arms is given as a function of distance r from the centre as β(r) = br + c, where b and c are parameters that determine the growth rate and the initial thickness of the spiral shell.

The different free parameters (a,p,b,c, inclination i, and w) were changed to study their effect with the aim of understanding the morpho-kinematical structure visible in the ¹³CO J = 6–5 PV diagrams. The inclination angle i is defined as the angle between the spiral axis and the line of sight, or, equivalently, the angle between the equator and the plane of the sky. Homan et al. (in prep.) will present a much more detailed study. Although limited by the short duration of the observations, we estimate a value for the spiral-arm spacing a in the range of 1–2′′ using Fig. 8.

As setup for our basic spiral model, we chose as free parameters i = 80° (Jeffers et al. 2014), a = 1.45′′, p = 2, b = 0, c = 8 × 10¹² cm (or 16 R_⋆), and w = 50° (see Table 3). Each of the parameters was then consecutively changed to demonstrate its effect (see Fig. 15). The PV diagram for the basic model of a one-armed spiral at a position angle along the orbital plane (second panel in Fig. 15) shows very regular structures due to the formation of higher density arcs in a plane perpendicular to the orbital plane (see also Fig. 11 in Mastrodemos & Morris 1999), although its signature is slightly more asymmetric than the PV diagram for the radial outflowing wind with density-enhanced shells (first panel in Fig. 15). At position angles not aligned with the orbital axis or orbital plane, the difference with the radial outflowing wind is much more pronounced, as can be seen in the first panel of Fig. 16. At a PA along the orbital plane, the morphological difference between the models for an inclination i = 80° and i = 20° is mainly seen around zero offset, for which a higher intensity contrast is reached for i = 20°. For position angles not aligned with the orbital axis or orbital plane, a more distinct morphological difference is seen for different inclination angles (see Fig. 16): at inclination angle ≳60° an “S-type” structure around zero-offset is seen in the PV-diagram. This S-type signature is also seen in the ¹³CO J = 6–5 ALMA images of CW Leo (see Fig. 11); this feature cannot be reproduced by a model assuming a radial outflowing wind with density-enhanced shells (see Fig. 14). The difference in PV diagrams for different position angles at different inclinations is used below to constrain the inclination angle of the spiral arms in CW Leo (Sect. 8.1).

As demonstrated by Kim & Taam (2012c), the mean density in the binary-induced stellar wind is proportional to r^−p, with p ~ 2. Figure 15 shows the effect for increasing the value of p from 2 to 2.4 (although the mass is not conserved for this latter simulation). For p ~ 2.4, the density is more concentrated in the inner arm spirals, clearly visible in the PV diagram.

When the parameter b is larger than zero, the thickness of the spiral arms increases with radial distance from the star, simulating the effect of dispersion. This leads to a lower mean density in the outer spiral arms. In contrast, lowering the initial thickness of the spiral shells from 8 × 10¹² cm to 3 × 10¹² cm results in a higher mean density in the complete spiral structure because of the conservation of mass. Decreasing the parameter a leads to a decrease in the distance between the successive arm turnings of the spiral.

A value for the parameter w around 50 results in the mean density being rather uniformly distributed in latitude, while the much lower value of w = 5 results in the latitudinal dependence of the density being more concentrated toward the orbital plane, hence mimicking the effect of gravitational focusing. While at a PA along the orbital axis the difference is not so large, the latitudinal dependence is more pronounced at PA = 145° (see Fig. 16).

As discussed by Mastrodemos & Morris (1999), the effect of a binary companion is to focus some fraction of the material into a spiral structure, the remainder still forming a quasi-spherical envelope. This smooth wind structure is well-visible in the ALMA PV diagrams as the dominant maximum around zero spatial offset covering the full velocity range (see Fig. 17). As will be shown by Homan et al. (in prep.), the flux contrast between the correlated arcs and bright inner bar can be used to deduce the density contrast between the smooth wind and spiral structure.

8. Discussion

8.1. Morphological comparison between ALMA ¹³CO PV diagrams and the theoretical simulations

With an on-source integration time of only 17 min, an in-depth understanding of the ¹³CO J = 6–5 channel maps (Fig. 8) is challenging. However, it is clear that the different PV diagrams (Figs. 11 and 12) display correlated structures that cannot be interpreted with a simple spherical symmetric wind model. While the wind model with the density-enhanced shell structures shows some resemblance to the ALMA PV data, the asymmetry when reflecting around zero-offset seen in the ALMA PV diagrams (Fig. 12) cannot be explained by this simple model, neither can it explain the difference in the PV diagrams with 90° difference in PA. A possible scenario to explain the structures seen in the ALMA PV diagrams is a binary-induced spiral structure. This interpretation is corroborated by other data and analyses already presented in the literature (see Sect. 6). Using the simple mathematical prescription for a binary-induced spiral structure described in Sect. 7.2, we derive from the ALMA PV diagrams at PA ~ 10–20° that the spiral parameter a is ~1.1′′. Different combinations for p, b, and c give a good fit to the data, but the current limitations in the ALMA sensitivity of the current data and the limited UV-coverage inhibit a strict constraint. We therefore opted to use p = 2 (mass conservation) and b = 0 (the thickness of the spiral arms remains constant). A good fit to correlated structures in the ALMA PV data is obtained for the spiral parameters a = 1.12″, p = 2, b = 0, c = 3 × 10¹² cm, and w = 30°, and i = 60−80° (see Fig. 18). The best constrained parameter is the spiral arm distance a, with an uncertainty of only ~10%. The general morphology is reproduced quite nicely, but the fit to the relative intensities can be improved. The latter is expected taking into account the simple approach for modelling the spiral structure. An in-depth study describing the results of a large parameter grid, including a variation in radial temperature structure of the spiral, density contrast, etc. will be presented in Homan et al. (in prep.).

It is clear that the Shape simulations shown in Fig. 18 do not capture the clumpiness seen in the ALMA data. Part of the fragmentation is artificial and comes from the poor UV coverage of the data (see Fig. 14). However, part of the clumpiness seems also inherent to the wind structure of CW Leo and indicates a rough gas (and dust) density distribution within the spiral arms. Magnetic cool spots, photospheric convective motions, flow instabilities, etc. may result in localized gas and dust clumps (see the discussion in Sect. 6). The fact that we see the arc structure extending toward the orbital axis (see illustration in Fig. 13) implies that the reflex motion of the primary AGB star (CW Leo itself) acts as a cause for the arc-like signature detected in the ALMA PV, Herschel, and Hubble Space Telescope images.

8.2. Orbital period and binary separation

The model-fitting of the ALMA ¹³CO J = 6–5 PV data indicates a spiral-arm spacing of a ~ 1.12′′, which agrees with the arc spacing seen in the HST wide-V band image of the inner 10′′ region of CW Leo (see Fig. 3 in Mauron & Huggins 2000). The arm-spacing is determined by the product of the binary orbital period T_p and the pattern propagation speed in the orbital plane, which following Kim & Taam (2012a) is given by(8)with ⟨ V_w ⟩ the wind velocity, V_p the orbital velocity, and r_p the orbital radius of the primary, and the orbital period T_p given by the last term (2πr_p/V_p). The first term in the right-hand part of Eq. (8), , is the pattern propagation speed throughout the orbital plane, which is close to the wind speed if V_w dominates over the orbital and sound speeds (which is the case for almost all AGB binary simulations, see Kim & Taam 2012c, and Kim priv. comm.). At a distance of 1′′, the wind has already reached its terminal velocity of 14.5 km s^-1 (see Sect. 4.3). The derived arm-spacing hence results in an orbital period of ~55 yr.

The derived orbital period is lower than the orbital periods derived for some other carbon-rich AGB stars, of which a binary companion is thought to be the cause of the detected spiral arm structure. Maercker et al. (2012) derived for R Scl a period of 350 yr on the basis of CO J = 3–2 ALMA data, for CIT 6 Dinh-V.-Trung & Lim (2009) found a period of ~600 yr, and for AFGL 3068 an orbital period of 830 yr was obtained by Mauron & Huggins (2006). The longer orbital periods for these AGB stars indicate that the companion must be on a wide orbit around the AGB star; typically around 50–70 au. According to Kepler’s third law and assuming a primary initial mass of 4 M_⊙ for CW Leo (Guélin et al. 1995) and the mass of the secondary being lower than that of the primary, the mean binary separation, d, for CW Leo is ~25 ± 2 au (or ~8.2 R_⋆ = 0.17′′), that is, a binary system that can reside within the inferred dust lane of ~0.5–1′′ radius (see Fig. 13). If the primary has already lost some 2 M_⊙ through its stellar wind (Decin et al. 2011), the binary separation would reduce to ~19 au.

8.3. Companion mass

According to the method outlined by Kim & Taam (2012a), one needs the (projected) separation between the primary star and its companion to derive the binary mass ratio. However, even after an extensive search of the literature, no direct identification of the potential binary companion was possible. We are hence lacking diagnostics to constrain the mass of the companion. However, because the ALMA PV data around PA ~ 10–20° might be explained by a spiral shock caused by the reflex motion of the mass-losing primary star implies that the mass of the secondary cannot be very low, meaning that the secondary cannot be a Jupiter-like planet or a brown dwarf. Although the simulations of Mastrodemos & Morris (1999) were not fine-tuned toward the specific situation of CW Leo (i.e., the radius of CW Leo is larger than the input radius for the model simulations, and the estimated mass of CW Leo is around 4 M_⊙, while the input model mass for the primary, M_p, is 1.5 M_⊙), we can use their results as a guideline to estimate the mass-ratio, M_s/M_p, with M_s the mass of the secondary. Their model sequence M10 → M17 → M18 shows that for a decrease of the mass of the secondary (1 M_⊙ → 0.5 M_⊙→  0.25 M_⊙) the morphology changes from bipolar → elliptical → quasi-spherical, meaning that the latter case resembles the signature in our ALMA data and multiple shells seen in the optical and infrared images (Mauron & Huggins 2000; Decin et al. 2011). This would imply a mass ratio of about 1/6 or a mass for the secondary around 0.6 M_⊙. Using Eq. (2) in Eggleton (1983), this implies that the effective radius of the Roche lobe, r_L, is 0.54 times the orbital separation, or r_L is ~13.5 au (4.5 R_⋆).

An alternative method was offered by Huggins et al. (2009), who used the models of Mastrodemos & Morris (1999) to develop a prescription for the observed envelope shapes in terms of the binary parameters. The companion modifies the mass loss by gravitationally focusing the wind towards the orbital plane and thereby determines the shape of the envelope at large distances from the star. They defined the envelope shape parameter, K_n, as the density contrast between the pole and the equator(9)with n_eq and n_po being the density at the equator and pole. Using the results of Mastrodemos & Morris (1999), Huggins et al. (2009) developed an empirical relation between the envelope shape parameter K_n and the binary parameters, being the masses of the primary and secondary (M_p and M_s), the binary separation (d) and the velocity of the wind at the orbit of the secondary (V_s). From the numerical simulations, they derived following equation (Eq. (4) in their paper): (10)independent of the mass of the primary M_p, and with M_s in units of solar mass, V_s in km s^-1 and d in astronomical units (au). Although the models were not fine-tuned toward the specific situation of CW Leo (i.e., the radius of CW Leo is larger than the input radius for the model simulations), we can use this relation and the value of w, characterizing the latitudinal dependence (see Sect. 7.2.2), to derive a first crude estimate of the companion mass.

In Sect. 8.1 we derive a best-fit value for w of ~30°, resulting in a value of K_n of about 10 (using Eq. (6)). Using Eq. (10), we derive that(11)For a binary separation d of 20 au, and the velocity of the wind at the orbit of the secondary, V_s, being ~8 km s^-1 (see Fig. 5) the derived mass for the secondary, M_s, is ~1.1 M_⊙. However, the uncertainty on the estimated value for M_s is significant because of the square dependence on the highly uncertain value of V_s and the uncertainty of ~15° in the value of w.

The derived rough estimate for the companion mass places the companion in the category of a white dwarf or an unevolved low-mass main sequence star. However, a binary with a white dwarf companion could result in a so-called dusty symbiotic system. This type of system harbours a very hot and active accretion disk, with clear signatures in the UV and many forbidden atomic lines of [O III], [Ne III], and [Fe IV] in the optical spectrum, for example. Since no signs of activity are seen in the optical/UV spectra of CW Leo, we tentatively postulate that the companion of CW Leo is an unevolved M-type dwarf.

8.4. Other scenarios

One might wonder if localized mass-loss ejections radially driven outward, without the action of a binary companion, might result in the type of PV diagrams seen in the ALMA data. Dinh-V-Trung & Lim (2008) suggested that the shells seen in the VLA data of HC₃N and HC₅N at a distance of ~15′′ from the central star cover some 10% of the stellar surface at the time of ejection. This scenario can be simulated using the simple parametric description for density-enhanced shells (Sect. 7.2.1), but constraining randomly the longitudinal and latitudinal coordinate for each of the five shells so that each segment covers only ~10% of the surface. This suggestion provides much (parametric) freedom to model the wind of CW Leo. An example of such a PV diagram is shown in the left panel of Fig. 19. While this example shows some resemblance with the ALMA PV diagrams, this scenario would need much fine-tuning to obtain a model showing almost symmetrical signatures as seen in the ALMA PV diagram for PA ~ 10–20° (Fig. 12), but for which almost no structure is left in a PV diagram taken at a PA 90° apart.

A way to constrain the huge parameter space entailed by the model for the localized mass-loss ejections is offered by the model presented by Le Bertre et al. (1989), who suggested that ejected material also moves as a result of the stellar rotation, creating latitudinal bands. This type of model was recently used by Jeffers et al. (2014) to explain the ExPo images of CW Leo. The Shape simulations for this scenario are also based on the model for the density-enhanced shells (Sect. 7.2.1), but the latitudinal coordinates of each shell only cover ~5°. These simulations do not rule out the ejection of clumps with initial extent ~10% of the stellar surface, nor the possible effects of stellar rotation, but the addition of the effects of the companion in creating a spiral provides a more robust explanation for the asymmetry in the PV diagrams when reflecting around zero-offset, for the differences seen in the PV diagrams 90° apart and for the other diagnostics discussed in Sect. 6 and shown in Fig. 13.

8.5. Further evolution

There are numerous similarities between CW Leo and the young post-AGB star CRL 2688 (also called the Egg Nebula; Sahai et al. 1998). Both objects are carbon-rich, suggesting an initial mass higher than 3 M_⊙, Although CRL 2688 is more advanced in its evolution, both objects have an optically thick circumstellar core with biconical cavities, and shell-like density enhancements are observed in scattered light as incomplete arcs. In CRL 2688 high-speed collimated jets are present with velocities of about 300 km s^-1, which seem to be only powered during the last 200 years. According to Soker & Livio (1994), different mechanisms can explain the collimated outflows in CRL 2688, each of them requiring the central star to be a binary, with the collimated outflow resulting from an accretion disk. The discovery of a spiral-arm structure in the ALMA data of CW Leo, which might be induced by a binary companion, supports the suggestion that CW Leo will evolve into an object similar to CRL 2688.

9. Conclusions

We have presented the first ALMA band 9 data of the inner wind region of CW Leo at a spatial resolution of 0.̋42 × 0.̋24 and with a sensitivity of 0.2–0.3 Jy beam^-1 per 488 kHz channel. We detected 25 emission lines, all of them centred on the dust continuum peak, which has a total (star+dust) flux density within the 3σ contour of 5.66 Jy. The images prove that the vibrational lines are excited just above the stellar photosphere and that SiO is a parent molecule, formed close to the stellar surface, probably by shock-induced non-equilibrium chemistry. The velocity traced by the line width of the emission lines implies a steep increase of the wind velocity starting at around 5 R_⋆, almost reaching the terminal velocity at ~11 R_⋆.

Both the dust emission and the emission of the brightest lines show a clear asymmetric distribution. We note, however, that the incomplete UV-coverage partly limits the quality of the data. The PV maps at different position angles of the ¹³CO J = 6–5 line display correlated arc-like structures, which can be explained by a spiral arm. This spiral arm can be caused by a binary companion. Using the Shape modeling tool, we modelled the ¹³CO J = 6–5 PV diagrams. We deduce that the orbital axis lies at a position angle of ~10–20° north-east and that the spiral arm spacing is ~1.12′′. At a distance of 150 pc, this leads to an orbital period of 55 yr and a binary separation of ~20–25 au (or ~6–8 R_⋆). We tentatively suggest that the companion is an unevolved low-mass main-sequence star.

The scenario of a binary system can explain (1) the spiral arm structure seen in the ALMA PV data; (2) the signature of a bipolar structure seen at arcsecond scales; and (3) the multiple non-concentric shells detected in the outer wind, which might represent the limb-brightened edges of the spiral arms seen almost edge-on. Granted ALMA cycle 1, time will increase the S/N and UV-coverage, hence facilitating the in-depth interpretation of the ALMA data.

---

Acknowledgments

This paper makes use of the following ALMA data: ADS/JAO.ALMA#2011.0.00277.S. ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan), together with NRC (Canada) and NSC and ASIAA (Taiwan), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO and NAOJ. W.S. acknowledges support by grant UNAM-PAPIIT 101014. The authors thank Nicholas Koning and Miguel Santander-García for technical support with Shape.

References

Online material

Appendix A: Position-velocity diagrams

All Tables

All Figures