© J. M. Winters et al. 2022

1 Introduction

Mass-loss in red giants is due to a combination of stellar pulsations and radiation pressure on dust forming in dense shocked regions in the outer stellar atmosphere (e.g., Höfner & Olofsson 2018). Even if the basic principles are understood, a fully consistent picture – including the role of convection, the time-dependent chemistry, and a consistent description of dust formation – still needs to be developed. In particular, the contribution of transparent grains to the acceleration of matter close to the stellar photosphere (Norris et al. 2012) still needs to be assessed.

The mechanisms shaping circumstellar environments around asymptotic giant branch (AGB) stars are vividly debated. Among them, magnetic fields (Matt et al. 2000; Duthu et al. 2017), binarity (Theuns & Jorissen 1993; Mastrodemos & Morris 1999; Decin et al. 2020), stellar rotation (Dorfi & Höfner 1996), and common-envelope evolution (Olofsson et al. 2015; Glanz & Perets 2018) have been considered.

A major difficulty is to explain the observed velocity field in axi-symmetrical sources, with larger velocities at high latitudes than at low latitudes (Hoai et al. 2014; Nhung et al. 2015b). Also, recent observations of rotating structures and streams bring additional conundrums (Tuan-Anh et al. 2019; Hoai et al. 2019).

We have concentrated our efforts on two relatively close (d ~ 150 pc) sources that show composite profiles in CO rotational lines (Winters et al. 2003): EP Aqr (Winters et al. 2007, hereafter referred to as W2007), and RS Cnc (Libert et al. 2010). Data obtained at IRAM show that these two sources have an axi-symmetrical structure with a low-velocity (~ 2 km s⁻¹) wind close to the equatorial plane, and faster (~ 8 km s⁻¹) outflows around the polar axes (Hoai et al. 2014; Nhung et al. 2015b). For EP Aqr, W2007 find a radial dependence of the density showing intermediate maxima. Additional data obtained with ALMA (Nhung et al. 2019b; Homan et al. 2018b) reveal a spiral structure explaining the earlier W2007 results.

RS Cnc is one of the best examples of the interaction between the stellar wind from an AGB star and the surrounding interstellar medium (Hoai et al. 2014). Its high declination makes RS Cnc an ideal target for the Northern Extended Millimeter Array (NOEMA). Previous studies based on IRAM data show that it is a twin of EP Aqr, but observed at a different angle, with a polar axis inclined at about 30° with respectto the line of sight (Libert et al. 2010; Hoai et al. 2014; Nhung et al. 2015b). This is favorable for studying polar and equatorial structures simultaneously, whereas the different viewing angle between EP Aqr and RS Cnc can be exploited to discriminate between different models in explaining the observed composite CO line profiles (Le Bertre et al. 2016). In contrast to EP Aqr, technetium is detected in the atmosphere of RS Cnc (Lebzelter & Hron 1999), proving that it is evolving along the thermal pulsing asymptotic giant branch (TP-AGB) in the Hertzsprung-Russell (HR) diagram. From a chemical point of view, RS Cnc is in a slightly more advanced evolutionary stage on the AGB, as indicated by its spectral classification as an MS star (see below) and by a higher photospheric ratio of ¹²C/¹³C (~ 35; Smith & Lambert (1986), but see Sect. 4.1 for an improved evaluation based on CO rotational lines from the circumstellar environment).

RS Cnc is a semi-regular variable star with periods of ≈ 122 d and ≈ 248 days (Adelman & Dennis 2005), located at a distance of ≈150 pc (Gaia Collaboration 2021; Bailer-Jones et al. 2021). It is listed as S-star CSS 589 in Stephenson (1984) based on its spectral classification M6S given in Keenan (1954). With its weak ZrO bands, its chemical type is intermediate between M and S (Keenan 1954). The stellar temperature is estimated to T_* ≈ 3200 K and its luminosity is L_* ≈ 4950 L_⊙ (Dumm & Schild 1998). From CO rotational line observations, two circumstellar wind components were identified: an equatorial structure expanding at about 2 km s⁻¹ and a bipolar outflow reaching a terminal velocity of v_exp ≈ 8 km s⁻¹ (Libert et al. 2010; Hoai et al. 2014), carrying mass-loss rates of 4 × 10⁻⁸ M_⊙ yr⁻¹ and 8 × 10⁻⁸ M_⊙ yr⁻¹, respectively (see Sect. 4.1 for an improved value of the mass-loss rate derived here). Lines of ¹²CO, ¹³CO, SiO, and HI were detected from previous observations at millimeter (mm) and radio wavelengths (Nyman et al. 1992; Danilovich et al. 2015; de Vicente et al. 2016; Gérard & Le Bertre 2003; Matthews & Reid 2007).

NOEMA was recently equipped with the wide band correlator PolyFiX, covering a total bandwidth of 15.6 GHz and therefore offering the potential to observe several lines from different species simultaneously. In this paper we present new data obtained with NOEMA in D- and A-configuration, complemented by short spacing observations obtained at the IRAM 30m telescope. Observational details are summarized in Sect. 2 and our results are presented in Sect. 3. Section 4 contains a discussion of the morphological structures and compares them to similar structures found in EP Aqr. Our conclusions are summarized in Sect. 5.

2 Observations

New observations of RS Cnc have been obtained in CO(2–1) with NOEMA/WideX in the (extended) nine-antenna A-configuration in December 2016 (Nhung et al. 2018) and with NOEMA/PolyFiX in the (compact) nine-antenna D-configuration during the science verification phase of PolyFiX in December 2017 and in the ten-antenna A-configuration in February 2020. The WideX correlator covered an instantaneous bandwidth of 3.8 GHz in two orthogonal polarizations with a channel spacing of 2 MHz. Additionally, up to eight high-spectral resolution units could be placed on spectral lines, providing channel spacings down to 39 kHz. WideX was decommissioned in September 2017 and replaced in December 2017 by the new correlator PolyFiX. This new correlator simultaneously covers 7.8 GHz in two sidebands and for both polarizations, and provides a channel spacing of 2 MHz throughout the 15.6 GHz total bandwidth. In addition, up to 128 high-spectral-resolution “chunks” can be placed in the 15.6 GHz-wide frequency range covered by PolyFiX for both polarizations, each providing a fixed channel spacing of 62.5 kHz over their 64 MHz bandwidth.

RS Cnc was observed with two individual frequency setups covering a total frequency range of ≈ 32 GHz in the 1.3 mm atmospheric window (see Fig. 1). We used the two quasars J0923+282 and 0923+392 as phase and amplitude calibrators; these were observed every ≈20 min. Pointing and focus of the telescopes was checked about every hour, and corrected when necessary. The bandpass was calibrated on the strong quasars 3C84 and 3C273, and the absolute flux scale was fixed on MWC349 and LkHa101, respectively. The accuracy of the absolute flux calibration at 1.3 mm is estimated to be better than 20%.

In order to add the short spacing information filtered out by the interferometer, in May and July 2020 we observed at the IRAM 30m telescope maps of 1′ by 1′ using the On-The-Fly (OTF) mode. This turned out to be necessary for the ¹²CO(2–1) and ¹³CO(2–1) lines but was not needed for the SiO lines, whose emitting region was found to be smaller than ~ 3′′. In the case of the ¹²CO(2–1) and ¹³CO(2–1) lines, the interferometer filters out large-scale structures that account for about two-thirds and three-quarters, respectively, of the total line flux, information that is recovered by adding the short spacing data from the OTF map. A comparison of the respective line profiles is shown in Fig. A.1.

The data were calibrated and imaged within the GILDAS¹ suite of software packages using CLIC for the NOEMA data calibration and the uv table creation, CLASS for calibrating the OTF maps, and the MAPPING package for merging and subsequent uv fitting, imaging, and self-calibration of the combined data sets. Continuum data were extracted for each sideband of the two frequency setups individually by filtering out spectral lines, and then averaging over 400 MHz bins to properly rescale the uv coordinates to the mean frequency of each bin. Phase self-calibration was performed on the corresponding continuum data. The gain table containing the self-calibration solutions was then applied to the spectral line uv tables using the SELFCAL procedures provided in MAPPING.

The resulting data sets were imaged applying either natural weighting, or, on the high-signal-to-noise (S/N) cubes, by applying robust weighting with a threshold of 0.1 to increase the spatial resolution by typically a factor 2. The resulting dirty maps were then CLEANed using the Hogbom algorithm (Högbom 1974).

The beam characteristics and sensitivities of the individual combined data sets from A- and D-configuration (and including the pseudo-visibilities from the OTF maps, where appropriate) are listed in Table 1 for all detected lines.

3 Results

The PolyFiX data, covering the frequency ranges 213–221 GHz (setup1, LSB), 228-236 GHz (setup1, USB), 243–251 GHz (setup2, LSB), and 258–266 GHz (setup2, USB) with two setups (see Fig. 1), showed different lines of CO and SiO, and, for the first time, many lines of species like SO, SO₂, HCN, and PN and some of their isotopologs. Furthermore, the data confirmed the H₂O line at 232.687 GHz already detected serendipitously with WideX in 2016, with a second H₂O line at 263.451 GHz seen for the first time in RS Cnc.

All linescovered by the same setup (1 or 2, see Fig. 1) share the same phase-, amplitude-, and flux calibration. All 32 detected lines are listed in Table 1.

3.1 Continuum

Figure 2 shows the continuum map from A-configuration only, using robust weighting to increase the spatial resolution to 0.39′′ × 0.22′′ at PA 28°. After self-calibration, S/N = 492 is obtained. The continuum source is unresolved, a point source fit results in a flux at ≈ 247 GHz of 23.65 ± 4.7 mJy (where the quoted error accounts for the accuracy of the absolute flux calibration of 20%) and a source position at RA = 09:10:38.780 and Dec = 30:57:46.62 in February 2020. All line data cubes discussed in the remainder of this paper are re-centered on this continuum position.

The source position is offset from the J2000 coordinates by −0.26′′ in RA and by −0.68′′ in Dec, consistent with the proper motion of RS Cnc (−10.72 mas yr⁻¹ in RA and −33.82 mas yr⁻¹ in Dec, Gaia Collaboration 2021; Bailer-Jones et al. 2021). From the PolyFiX data, spanning a total frequency range of about 53 GHz, we determine a spectral index of 1.99 ± 0.09 for RS Cnc in the 1 mm range, which is fully consistent with a black body spectrum of the continuum (see also Libert et al. 2010).

Fig. 1 Overview of the frequency ranges observed with PolyFiX using two spectral setups (setup1: red and setup2: blue, respectively). Lower diagrams: zoom onto the individual spectra covering 7.8 GHz each. Upper row: setup1, lower row: setup2. Thecentral 20 MHz at the border between inner and outer baseband are blanked out, i.e., set to zero, as this region is contaminated by the LO2 separation of the 8 GHz-wide IF in the IF processor (“LO2 zone”).

3.2 Detected molecules and lines

Within the total frequency coverage of about 32 GHz, we detect 32 lines of 13 molecules and isotopologs, including several transitions from vibrationally excited states. All these lines are listed in Table 1 and are presented in the following sections. The peak flux and FWHP of the line-emitting regions, as listed in Table 1, are determined by circular Gaussian fits in the uv-plane to the central channel (if the source is (partially) spatially resolved) or by point-source fits to the central channel (if the source is unresolved). All line profiles shown in the following sections in Fig. 3 and Figs. 5 through 14 are integrated over square apertures whose sizes are given in each figure caption. Two-component profiles are seen in CO and ¹³CO only, and not in any other of the lines detected here.

We looked for but did not detect the vibrationally excited ¹²CO(v = 1, 2–1) line, nor do we detect C¹⁸O(2–1), resulting in 3σ upper limits for the line peaks of 6 mJy beam⁻¹ and 3 mJy beam⁻¹, respectively (the ¹²CO(v = 1, 2–1) line was not covered in our A-configuration data).

3.2.1 CO

The profiles of ¹²CO(2–1), ¹³CO(2–1) (see Fig. 3), and ¹²CO(1–0) (see Libert et al. 2010) show a very distinct shape composed of a broad component that extends out to v_lsr,* ± 8 km s⁻¹ and a narrow component indicating velocities of ±2 km s⁻¹ with respect to v_lsr,* = 7 km s⁻¹. Velocity-integrated intensity maps of CO are shown in Fig. 18, indicating a clear kinematic structure in the north–south direction. InFig. 4, we present a schematic representation of the geometrical structure of RS Cnc as implied by the data; see Sect. 4.1. The CO emitting region is spatially extended, consisting of a dense equatorial structure that corresponds to the low-velocity expansion and an inclined, bipolar structure corresponding to an outflow at a projected velocity of 8 km s⁻¹. These structures were discussed in Hoai et al. (2014) based on Plateau de Bure data obtained on ¹²CO(2–1) and ¹²CO(1–0) that had a spatial resolution ofabout 1′′. The model built by these latter authors was later refined by Nhung et al. (2018) based on ¹²CO(2–1) data obtained with the WideX correlator in NOEMA’s A-configuration, providing a spatial resolution of 0.44′′ × 0.28′′. Nhung et al. (2018) find a position angle of the projected bipolar outflow axis of ω =7° (measured counter-clockwise from north) and an inclination angle of the outflow axis with respect to the line of sight of i = 30°. The CO distribution is further investigated in Sect. 4.1 below.

Such a structure had already been found in the S-type star π¹ Gru (Sahai 1992), which was later confirmed by higher spatial resolution observations using ALMA (Doan et al. 2017). This object has a G0V companion (Feast 1953) and possibly a second, much closer companion (Homan et al. 2020). In Hoai et al. (2014), we reported for RS Cnc the possible presence of a companion seen in the ¹²CO(1–0) channel maps at velocities around 6.6 km s⁻¹ and located about 1′′ west-northwest of the continuum source. The new data allow for a more detailed study of this feature, which is presented in Sect. 4.1.

Fig. 2 Continuum map around 247 GHz from A-configuration. Contours are plotted in 100σ steps, where 1σ is 47.6 μJy beam−1. The synthesized beam is indicated in the lower left corner.

Fig. 3 CO line profiles, showing a two-component structure. Left: 12CO(2–1). Right: 13CO(2–1). A-configuration and D-configuration are merged, OTF data are added, and the spectral resolution is 0.5 km s−1. The CO emission is integrated over the central 22′′ × 22′′, i.e., over thefull field of view of the NOEMA antennas at 230 GHz.

Fig. 4 Sketch of the geometrical structure of the wind components as inferred from the current data (see Sect. 4.1). The sketch is not to scale: there is a smooth transition between the equatorial enhancement and the polar outflows.

Fig. 5 Profiles of SiO ground-state and first vibrationally excited state lines. Left: SiO(6–5): upper: v = 1, lower: v = 0. A-configuration and D-configuration merged. Right: SiO(5–4): upper: v = 1, D-configuration (black) and A-configuration (red), lower: v = 0, A-configuration and D-configuration merged. The spectral resolution is 1 km s−1 for (v = 1) and 0.5 km s−1 for the (v = 0) lines, respectively. The emission is integrated over the central 5′′ × 5′′ aperture.

3.2.2 SiO

We detect a suite of ²⁸Si¹⁶O (henceforth SiO) transitions, including the vibrational ground-state lines of SiO(5–4) and SiO(6–5), the first and second vibrationally excited state of SiO(5–4), and the first vibrationally excited state of SiO(6–5). All SiO profiles are shown in Figs. 5 and 6. The spatial region emitting the vibrational ground-state lines extends out to about 2′′ from the continuum peak (see Table 1, Fig. 18, and Sect. 4.2). Interestingly, we detect a strong maser component on the SiO(v = 1, J = 5–4) line at v_lsr ≈ 14 km s⁻¹ in the data obtained in December 2017, which had completely disappeared when we re-observed RS Cnc in February 2020 (see Fig. 5). Such behavior is well known for pulsating AGB stars, and lends support to the idea that the SiO masers are excited by infrared pumping as opposed to collisional pumping (see, e.g., Pardo et al. 2004).

The SiO(v = 2, J = 5–4) line is detected above the 3σ level of 3 mJy beam⁻¹ over a broad range of Doppler velocities from at least −5 to 18 km s⁻¹ (Fig. 6). Given its high excitation energy (~3500 K), we expect this line to trace exclusively the innermost region around RS Cnc, as was the case in o Cet, where SiO(v = 2) absorption and emission was spatially resolved by ALMA (Wong et al. 2016). Its broad line width suggests that it may trace the same high-velocity wings seen in other detected SiO lines (Sect. 3.3). However, our detection is too weak to allow for a detailed study of the morpho-kinematics of the emission.

At the upper edge of the LSB of setup 2 at 250.744 GHz, we serendipitously detect a strong line that we identify as ground-state Si¹⁷O(6–5) at 250.7446954 GHz (Müller et al. 2013) from the Cologne Database of Molecular Spectroscopy (CDMS², Müller et al. 2005); the profile is shown in Fig. 7. This line and other transitions of Si¹⁷O have already been detected in a number of well-studied objects, such as the S-type star W Aql (De Beck & Olofsson 2020), the M-type star R Dor (De Beck & Olofsson 2018), and the evolved, high-mass-loss-rate oxygen-rich star IK Tau (Velilla Prieto et al. 2017). No other Si¹⁷O transitions are covered in our setups, but there is a highly excited H₂O line at 250.7517934 GHz (v₂ = 2, J(K_a, K_c) = 9(2,8)–8(3,5); E_u ∕k = 6141 K) listed in the JPL catalog³ and predicted by Yu et al. (2012) from the Bending-Rotation approach analysis. If the detected line was H₂O emission, it would be redshifted from the systemic velocity by about 9 km s⁻¹. As indicated by the modeling of Gray et al. (2016), the 250.752 GHz line may exhibit strong maser action in regions of hot gas (T_kin = 1500 K) with cool dust (T_d ≤ 1000 K). While we cannot unequivocally exclude some contamination from a potential new, redshifted H₂O maser, we consider Si¹⁷O a more likely identification of the 250.744 GHz emission. From the respective integrated line intensities of Si¹⁶O(6–5) and Si¹⁷O(6–5), which are ~ 163 Jy km s⁻¹ and ~ 3 Jy km s⁻¹, and taking the difference of the Einstein coefficients of the transitions into account, we estimate the isotopolog ratio ¹⁶O/¹⁷O ~ 50, assuming equal excitation conditions for both transitions and optically thin emission of both lines. This value is much lower than the solar isotopic ratio of ~2700 (Lodders et al. 2009) due to dredge-up events (Karakas & Lattanzio 2014; Hinkle et al. 2016) and is broadly consistent with those obtained in the M-type star R Dor and the S-type star W Aql (61–74; De Beck & Olofsson 2018, 2020). The initial mass of RS Cnc is about 1.5 M_⊙ (Libert et al. 2010)⁴, which is in the same range as R Dor (1.4 M_⊙; De Beck & Olofsson 2018) and W Aql (1.6 M_⊙; De Nutte et al. 2017) that gives a ¹⁶O/¹⁷O ratio of <1000 (Hinkle et al. 2016). However, we note that the oxygen isotopic ratio (¹⁶O/¹⁷O) derived from the line intensity ratio is likely underestimated if the Si¹⁶O line is not optically thin, as has been shown in De Beck & Olofsson (2018), who obtained a value of ~ 400 in R Dor with radiative transfer modeling. Indeed, we demonstrate in Sect. 4.2 that the Si¹⁶O emission in RS Cnc is optically thick, especially within a projected radius of ~ 1′′. A photospheric ¹⁶O/¹⁷O ratio of 710 in RS Cnc (=HR 3639) was estimated by Smith & Lambert (1990) from the spectra of near-infrared overtone band transitions of C¹⁶O and C¹⁷O, which is probably a more realistic ratio. We do not cover C¹⁷O(2–1) in our setups and therefore cannot give an independent estimate of the ¹⁶O/¹⁷O ratio. As Si¹⁸O(6–5) and C¹⁸O(2–1) are either not covered or not detected, there is not enough information from our data to obtain a meaningful constraint on the initial stellar mass from oxygen isotopic ratios (e.g. from the ¹⁷O/¹⁸O ratio; De Nutte et al. 2017).

We detect a line at 247.482 GHz at low S/N that might be identified as ²⁹Si¹⁷O(v = 0,J = 6–5) at 247.4815250 GHz based on the line list by Müller et al. (2013) and used in the CDMS (see Fig. 7). However, in contrast to Si¹⁷O(6–5), ²⁹Si¹⁷O(6–5) has never been detected; only higher-J lines of ²⁹Si¹⁷O have been tentatively detected in R Dor (J = 7–6 and J = 8–7, De Beck & Olofsson 2018). More specifically, the 247.482 GHz line is seen with an integrated line intensity of ~ 0.08 Jy km s⁻¹ in our D-configuration data only, observed in December 2017, but it does not show up in the A-configuration data, taken in February 2020. This may largely be due to the much reduced brightness sensitivity in the A-configuration, which is a factor of approximately 15 smaller because of the smaller synthesized beam area rather than some variable maser action in this line. Based on the D-configuration data, the source position of the 247.482 GHz emission appears slightly offset toward the northwest direction from the Si¹⁷O(6–5) emission. Further data on ²⁹Si¹⁷O, possibly covering the J = 6–5, J = 7–6, and J = 8–7 transitions, would be needed to draw any firm conclusion.

Fig. 6 Profile around the SiO(v = 2, J = 5–4) line frequency. A-configuration and D-configuration merged. The spectral resolution is 1 km s−1 and the emission is integrated over the central 1′′ × 1′′ aperture.

Fig. 7 Line profiles of SiO isotopologs. Upper left: profile of the 247.482 GHz line, possibly 29Si17O(6–5); D-configuration, only. Lower left: Si17O(6–5): D-configuration, only (line was not covered in A-configuration). Right: 29SiO(5–4); A-configuration and D-configuration merged. The spectral resolution in all cases is 3 km s−1 and the emission is integrated over the central 5′′ × 5′′ aperture.

Fig. 8 HCN line profiles. Left: HCN(3–2); A-configuration and D-configuration merged with a spectral resolution of 0.5 km s−1. Right: H13CN(3–2); A-configuration and D-configuration merged with a spectral resolution of 3 km s−1. The emission of both lines is integrated over the central 2′′ × 2′′.

3.2.3 HCN

We clearly detect the HCN(3–2) and H¹³CN(3–2) lines; the profiles are displayed in Fig. 8, and velocity-integrated intensity maps of both species are shown in Fig. B.1. Both lines are slightly spatially resolved and a circular Gaussian fit to HCN(3–2) gives a peak flux of 1.12 Jy and a FWHP size of 0.76′′ on the merged data. To our knowledge, this is the first detection of HCN and H¹³CN in RS Cnc (see Sect. 4.4). From the first-moment map (shown in Fig. 17, left), a clear velocity pattern is evident that indicates possible rotation in the HCN-emitting region (see Sect. 3.4). Also, the velocity-integrated intensity maps presented in Fig. B.1 show a clear kinematic structure in the east–west direction.

Formation of the HCN molecule in oxygen-rich environments is further discussed in Sect. 4.4. A modeling using the 1D local thermodynamic equilibrium (LTE) radiative transfer code XCLASS (Möller et al. 2017, see Appendix D) gives a column density for HCN in RS Cnc of N_HCN = 1.6 × 10¹⁵ cm⁻², corresponding to an abundance of X(HCN/H₂) = 6.6 × 10⁻⁷. This value is well within the range found for other M- and S-type stars as modeled by Schöier et al. (2013), who find X(HCN/H₂) equal to a few times 10⁻⁷ (for more details see Sect. 4.4 and Appendix D).

3.2.4 H₂O

The WideX spectrum obtained in A-configuration in December 2016 serendipitously revealed a line at 232.687 GHz that we ascribe to the J(K_a,K_c) = 5(5,0)–6(4,3) transition ofo-H₂O in the v₂ = 1 vibrational state. The H₂O source is weak and seems still unresolved within the synthesized beam of 0.5′′ × 0.34′′ obtained in the A-configuration in February 2020, consistent with its high upper-state energy of 3462 K. The line profile is shown in Fig. 9. With the follow-up observations employing PolyFiX in D-configuration and A-configuration we also covered and detected the 263.451 GHz o-H₂O v₂ = 1, J(K_a,K_c) = 7(7,0)–8(6,3) line (Fig. 9, right; E_u∕k = 4475 K). Both lines are resampled to a resolution of 3 km s⁻¹, data are merged from A-configuration and D-configuration, and the emission is integrated over an aperture of 1′′ × 1′′. Intensity maps of both lines are shown in Fig. B.2, testifying to the compactness of the H₂O-emitting region.

These are the first detections of millimeter vibrationally excited H₂O emission in RS Cnc. We note that the 22 GHz H₂O maser in the ground state was tentatively detected by Szymczak &Engels (1995) in one of the two epochs they covered, but the 22 GHz line is not detected in other observations (Dickinson et al. 1973; Lewis 1997; Han et al. 1995; Yoon et al. 2014). RS Cnc also shows clear photospheric H₂O absorption at 2.7 μm (Merrill & Stein 1976; Noguchi & Kobayashi 1993), and at 1.3 μm (7500 cm⁻¹; Joyce et al. 1998), although the H₂O band near 900 nm is not detected (Spinrad et al. 1966).

Both the 232 and 263 GHz water lines have upper levels belonging to the so-called transposed backbone in the v₂ = 1 vibrationally excited state of H₂O, that is K_a = J and K_c = 0 or 1 (see Fig. 1 of Alcolea & Menten 1993). The 232 GHz line was first detected in evolved stars together with the 96 GHz line from another transposed backbone upper level by Menten & Melnick (1989) toward the red supergiant VY CMa and the AGB star W Hya. The latter is an M-type star with a similar mass-loss rate to RS Cnc. The authors find that the 232 GHz line emission in both stars may be of (quasi-)thermal nature while the 96 GHz line clearly showed maser action. The (unpublished) detection of the 263 GHz line was mentioned in Alcolea & Menten (1993), who also described a mechanism that may lead to a systematic overpopulation of the transposed backbone upper levels in the v₂ = 1 state of H₂O in the innerregion of circumstellar envelopes. If the vibrational decay routes (to the ground state) of the transposed backbone upper levels become more optically thick than the lower levels in the v₂ = 1 state, then differential radiative trapping may cause population inversion of these lines. Additional vibrationally excited H₂O emission lines from transposed backbone upper levels were predicted and later detected in VY CMa by Menten et al. (2006) and Kamiński et al. (2013). We observed the 232 GHz line in RS Cnc at three epochs (December 2016, December 2017, and February 2020) and the 263 GHz line at the latter two epochs, and the emission appears to be stable in time for both lines. The profiles appear to be very similar, both are broad, even broader than the (ground-state) lines of other species reported here, and there is no sign for any narrow component in either of the two profiles at any of the epochs. As the lines should arise from a region very close to the star – compatible with their broad widths; see Sect. 3.3 – one might expect to see time variations due to the varying density and radiation field caused by the stellar pulsation, in particular if the emission were caused by maser action, as seen on the SiO(v = 1, 5–4) line observed in December 2017 (see Fig. 5). Also, the modeling of Gray et al. (2016) shows only very little inversion of the involved level populations for the 263 GHz H₂O transition. We therefore think that both lines could be thermally excited. A definite assessment of the nature of the vibrationally excited H₂O emission would however require some detailed modeling of the emission, together with high-sensitivity monitoring of the line profiles with high spectral resolution, possibly including other H₂O lines from transposed backbone upper levels and/or known maser lines for comparison, which is beyond the scope of the present paper.

Fig. 9 H2O line profiles. Left: H2O line at 232.687 GHz. Right: H2O line at 263.451 GHz. Data are merged from A-configuration and D-configuration, the spectral resolution is 3 km s−1, and the emission of both lines is integrated over the central 1′′ × 1′′ aperture.

Fig. 10 Profiles of the four detected SO lines, with A-configuration and D-configuration merged. The spectral resolution is 3 km s−1 and the emission is integrated over the central 2′′ × 2′′ aperture.

Fig. 11 Profiles of the two 34SO lines detected here, with A-configuration and D-configuration merged. The spectral resolution is 3 km s−1 and the emission is integrated over the central 2′′ × 2′′ aperture.

3.2.5 SO

Four lines of SO are detected (see Fig. 10) along with two lines of the isotopolog ³⁴SO (Fig. 11). These represent the first detections of SO and ³⁴SO in RS Cnc. SO has beenobserved in several M-type stars, including R Dor and W Hya, (Danilovich et al. 2016)), but remains undetected in S-type stars (e.g., W Aql, Decin et al. 2008; De Beck & Olofsson 2020). All SO lines detected here are slightly spatially resolved with a FWHP around 0.8′′ and therefore seem to be emitted from the same region as HCN. Velocity-integrated intensity maps of SO are shown in Fig. B.3. The SO lines show the same velocity pattern (indicating rotation) as HCN, although the velocity resolution of the SO lines is only 3 km s⁻¹; see Fig. B.3 and the first-moment map in the right panel of Fig. 17.

Using the integrated line strengths of SO(6(5)–5(4)) and ³⁴SO(6(5)–5(4)) found here (~ 4.69 Jy km s⁻¹ and ~ 0.20 Jy km s⁻¹, respectively) and taking the difference of the Einstein coefficients of the transitions into account, we estimate the isotopolog ratio ³²SO/³⁴SO ~ 23, assuming equal excitation conditions for both transitions and optically thin emission of both lines. This value is in good agreement with the values of 21.6 ± 8.5 and 18.5 ± 5.8 derived from the radiative transfer models for M-type stars by Danilovich et al. (2016, 2020), respectively. We note that, for the S-type star W Aql, an Si³²S/Si³⁴S isotopolog ratio of 10.6 ± 2.6 was derived by De Beck & Olofsson (2020). As ³²S is mainly produced by oxygen burning in massive stars and, to a lesser extent, in type Ia supernovae, and as ³⁴S is formed by subsequent neutron capture (e.g., Nomoto et al. 1984; Wilson & Matteucci 1992; Timmes et al. 1995; Hughes et al. 2008), the ³²S/³⁴S isotopic ratio remains virtually unaltered during AGB evolution (see, e.g. tables in the FRUITY⁵ database, Cristallo et al. 2011) and therefore should reflect the chemical initial conditions of the natal cloud from which the star has formed. The spread in the isotopic ratio seen among the different AGB stars mentioned above would then rather be indicative of the Galactic environment in which the star has formed (see, e.g., Chin et al. 1996; Humire et al. 2020) instead of reflecting any evolutionary effect. For the low-mass-loss-rate M-type stars R Dor and W Hya, Danilovich et al. (2016) reproduce their observed line profiles best with centrally peaked SO (and SO₂) distributions, consistent with the maps presented in Fig. B.3.

3.2.6 SO₂

In SO₂, 11 lines are detected; their parameters are summarized in Table 2, and all profiles are shown in Fig. C.1. These are the first detections of SO₂ in RS Cnc. A previous survey with the IRAM 30 m telescope by Omont et al. (1993) did not detect SO₂ in RS Cnc with an rms noise of 0.052 K (or ~0.25 Jy at 160.8 GHz). As an example, we show the SO₂ (14(0,14)–13(1,13)) line at 244.3 GHz, only in Fig. 12. A first-moment map of the SO₂ (14(0,14)–13(1,13)) line is shown in Fig. 17 in the middle left panel. Although the source remains barely resolved (source size ~ 0.43′′) by the beam(0.69′′ × 0.49′′), there is a signature of a rotating structure in SO₂, as was also seen in EP Aqr (Homan et al. 2018b; Tuan-Anh et al. 2019). Integrated intensity maps of three SO₂ lines (SO₂(9(3, 7)–9(2, 8)), which has the lowest upperlevel energy of the SO₂ lines detected here (E_u = 64 K); SO₂ (14(0,14)–13(1,13)), the strongest line, and SO₂(34(4,30)–34(3,31)), which has the highest upper level energy of the detected lines, E_u = 595 K) are shown in Fig. B.4. All lines show kinematic structure in the E–W direction, approximately orthogonal to the outflow structure seen in CO and SiO, cf. Fig. 18.

We derive the rotational temperature and column density of the SO₂ -emitting region with a population diagram analysis (Sect. 3.6) and by an XCLASS modeling (Appendix D). Both methods give a similar rotational temperature of ~ 320−350 K and a column density of ~ 3.5 × 10¹⁵ cm⁻².

3.2.7 PN

We detect a line at 234.936 GHz that we ascribe to the PN molecule, which would be the first detection of PN in RS Cnc. PN has been detected in several M-type stars (e.g., De Beck et al. 2013; Ziurys et al. 2018), and in the C-rich envelopes of IRC +10216 and CRL 2688 (Guélin et al. 2000; Cernicharo et al. 2000; Milam et al. 2008). The presence of PN in an MS-type star therefore does not seem to come as a surprise. However, RS Cnc appears to be the source with lowest mass-loss rate in which this molecule has been reported so far. The PN line profile is shown in Fig. 13. The line is spatially resolved at 0.8′′, which places it in about the same region as HCN and SO. The first-moment map of this line also shows signatures of rotation but due to the weakness of the line, the evidence is low. An integrated intensity map of PN is presented in Fig. B.5, showing that the line-emitting region is slightly spatially resolved. The 3σ feature seen about 1.5′′ south of the phase center should not be considered as a detection but rather as a noise peak, as long as this structure is not confirmed by higher sensitivity observations.

Fig. 12 Profile of SO2 (14(0,14)–13(1,13)) with A-configuration and D-configuration merged, a spectral resolution of 3 km s−1, and emissionintegrated over the central 2′′ × 2′′ aperture.

Fig. 13 Profile of PN(N = 5–4, J = 6–5) with A-configuration and D-configuration merged, a spectral resolution of 3 km s−1, and emissionintegrated over the central 2′′ × 2′′ aperture.

Fig. 14 Linewings in SiO(5–4) and SiO(6–5) compared to CO(2–1). The emission is integrated over the central 5′′ × 5′′ aperture.

Fig. 15 High velocities close to the line of sight as seen in SiO. PV maps are shown in the Vz vs. R plane for SiO(5–4) (left) and SiO(6–5) (right). The horizontal black line indicates the wind terminal velocity as traced in CO and the white scale bar indicates the spatial resolution. $R=\sqrt{{(\delta \text{Dec})}^2+{(\delta \text{RA})}^2}$, |Vz | = |vlsr − vlsr,*|.

3.3 High-velocity wings in SiO, and in other molecules

In SiO, five lines in three different vibrational states (v = 0,1,2) are detected (see Figs. 5 and 6). The vibrational ground-state lines clearly indicate the presence of material at velocities much higher than the wind terminal velocity of ~ 8 km s⁻¹ as traced by CO lines at this stellar latitude (see Sect. 4.1). This is illustrated in Fig. 14, and in Fig. 15 where we define v_z = v_lsr − v_lsr,*, the Doppler velocity relative to the star. The high-velocity region is centered on the line of sight and is confined to the inner ≈ 0.3′′; see Fig. 15. A similar feature was seen in high-spatial-resolution observations of other oxygen-rich, low-mass-loss-rate AGB stars, such as W Hya (Vlemmings et al. 2017), EP Aqr (Tuan-Anh et al. 2019), o Cet (Hoai et al. 2020), R Dor (Decin et al. 2018; Nhung et al. 2019a, 2021), and in 15 out of 17 sources observed in the ALMA Large Program ATOMIUM (Decin et al. 2020; Gottlieb et al. 2022), calling for a common mechanism causing high-velocity wings in this type of object. In the case of EP Aqr, where the bipolar outflow axis almost coincides with the line of sight (with an inclination angle of i ≈ 10°), the high-velocity wings were interpreted in terms of narrow polar jets. For R Dor and o Cet, which do not show obvious signs of axial symmetry in their winds, such an interpretation could not be retained and it was argued instead that the high-velocity wings were caused by (a mixture of) turbulence, thermal broadening, and some effect of shocks, acting at distances below some 10 to 15 AU from the central star. The presence of broad wings in the SiO lines emitted from RS Cnc, whose symmetry axis is inclined by ≈ 30° with respect to the line of sight (see Sect. 4.1), lends support to the latter type of interpretation and casts serious doubts on the polar jet interpretation proposed earlier for EP Aqr, which shows a morpho-kinematics similar to that of RS Cnc (Nhung et al. 2015b). Indeed, if the broad line widths are present regardless of the orientation of a possible symmetry axis with the line of sight, they must be caused by a mechanism of nondirectional (accounting for the resolving beam) nature. A possible candidate, whose action is limited to the close vicinity of the star, is pulsation-driven shocks that dissipate their energy relatively close to the star and imply positive and negative velocities in the shocked region that can be much higher than the terminal outflow velocity of the wind. Such structures could be explained by the B-type models discussed in Winters et al. (2000b) as presented in Winters et al. (2002); see their Fig. 3. Recent 3D model calculations that self-consistently describe convection and fundamental-mode radial pulsations in the stellar mantle would provide the physical mechanism that leads to the development of such shocks close to the star surface (e.g., Freytag et al. 2017) and could therefore replace the simplified inner boundary condition (the so-called “piston approximation”) that was used in the earlier 1D models mentioned above.

In the data presented here, wings at high Doppler velocity are seen in nearly all lines detected with sufficient sensitivity to probe the profile over at least v_lsr,* ± 10 km s⁻¹. This is illustrated in Fig. 16, where v_z profiles are integrated over a circle of radius 0.2′′ centered on the star. Gaussian profiles centered at the origin are shown as visual references (not fits), showing how absorption produces asymmetric profiles. A major difference is seen between vibrational ground-state lines, which have a Gaussian FWHM of ~ 10 km s⁻¹, and vibrationally excited-state lines, which have a Gaussian FWHM of ~ 14 km s⁻¹. Such a difference is not surprising, assuming that the high-velocity wings are formed in the inner layer of the circumstellar envelope (CSE), which is preferentially probed by the (v = 1) lines. In this context, we note that Rizzo et al. (2021) recently reported the detection of a narrow SiO(v = 1, 1–0) maser line in RS Cnc at a velocity of +14 km s⁻¹ with respectto the star’s lsr velocity. The effect of shocks on line profiles was first observed in the near-infrared range on CO ro-vibrational lines, probing the stellar photosphere and the innermost circumstellar region within ~ 10R_* (e.g., χ Cyg, an S-type star, Hinkle et al. 1982). Very-high-angular-resolution observations obtained over the past decade using VLT, VLTI, and ALMA show that the effect of shocks from pulsations and convection cell ejections is confined within some 10 AUfrom the star (see, e.g., Khouri et al. 2018; Höfner & Olofsson 2018; Ohnaka et al. 2019, and references therein). Rotation, when observed, is instead found to extend beyond this distance, typically up to ~ 20 AU (e.g., Vlemmingset al. 2018; Homan et al. 2018a; Nhung et al. 2021). The angular resolution of the present data is insufficient to detect such differences directly; however, the effect of rotation and shocks on lines of sight contained within a beam centered on the star depends on the region probed by each specific line: lines that probe the inner layers exclusively, such as the (v = 1) lines, are mostly affected by shocks, and somewhat by rotation; CO lines, for which the probed region extends very far out, see little effects of rotation and even less effects of shocks because the emission from the inner envelope provides too small a fraction of the total emission. Between these two extremes, the relative importance of the contributions of shocks and rotation depends on the radial extent of the region probed by the line. Such an interpretation is consistent with the data displayed in Fig. 16.

Fig. 16 Line profiles of different molecules on a logarithmic intensity scale. Gaussian profiles are shown for comparison, FWHM = 10 km s−1 for the ground-state lines of all molecules, and FWHM = 14 km s−1 for the (v = 1) lines of SiO. All observed profiles are integrated over R < 0.2′′.

3.4 Rotation

In Fig. 17, we present first-moment maps of HCN(3–2) (left), SO₂ (14(0,14)–13(1,13)) (middle left), SiO(v = 1, 6–5) (middle right), and SO(7(6)–6(5)) (right). At projected distances from the star not exceeding 0.5′′, all four tracers display approximate anti-symmetry with respect to a line at PA ~ 10°. This is suggestive of the presence of rotation in the inner CSE layer around an axis that projects on this line in the plane of the sky. Such a morpho-kinematic structure has also been observed in other stars, notably R Dor (Vlemmings et al. 2018; Homan et al. 2018a; Nhung et al. 2021). The angular resolution of the present data does not allow for a detailed exploration of this region, which prevents us from commenting on its possible cause. Nevertheless, the anti-symmetry axis of the velocity pattern projected on the plane of the sky at a PA that approximately coincides with the projected symmetry axis of the polar outflows (see Sect. 4.1) is remarkable and suggests that rotation is taking place about this same polar axis in the inner CSE layer.

The line-of-sight velocities of these structures are small, on the order of the velocities derived from CO for the equatorial region, and we interpret them here as possible signs of rotation (rather than indicating another bipolar outflow oriented perpendicular to the larger scale outflow traced in CO and SiO (v = 0) lines). We note that out of these four lines, the HCN(3–2) line is detected with the highest S/N (S/N = 233 in the line peak, cf. Table 1).

The mean Doppler velocity ⟨v_z⟩, averaged over the inner 0.5′′, of the HCN line can be fit in position angle ω, measured counter-clockwise from north, by $$\text{nn left}<v_z\text{ right}{>}_{\text{HCN}}=-0.19\text{km}{\text{s}}^{-1}+1.0\text{km}{\text{s}}^{-1}\sin (\omega -{19}^{°}),$$(1)

whereas the SiO(v = 1, 6–5) velocity is well fit by $$\text{nn left}<v_z\text{ right}{>}_{\text{SiO}(v=1,6-5)}=-0.37\text{km}{\text{s}}^{-1}+0.46\text{km}{\text{s}}^{-1}\sin (\omega -{26}^{°}).$$(2)

The small offsets of ~ − 0.3 km s⁻¹ on average are within the uncertainty attached to the measurement of the star’s LSR velocity. The coefficients of the sine terms measure the projected rotation velocity, namely the rotation velocity divided by the sine of the angle made by the rotation axis with the line of sight. Assuming that the rotation axis is the axi-symmetry axis of the CSE, this angle is i ~ 30° (see Sect. 4.1), meaning rotation velocities of ~2 and ~ 1 km s⁻¹ for HCN and SiO respectively. Observations of higher angular resolution are needed to confirm the presence of rotation within a projected distance of 0.5′′ from the star and we prefer to summarize the results presented in this section in the form of an upper limit to the mean rotation velocity of a few km s⁻¹.

Fig. 17 First-moment maps of different lines, indicating a possibly rotating structure (see Sect. 3.4). Left: HCN(3–2), middle left: SO2 (14(0,14)–13(1,13)), middle right: SiO(v = 1,6–5), right: SO(7(6)–6(5)). The black ellipses indicate the synthesized beam.

3.5 Global outflow structure traced by CO and SiO

The detailed structure of the morpho-kinematics of the CSE has been studied using observations of the ¹²CO(1–0) and ¹²CO(2–1) molecular line emission. The analyses of Hoai et al. (2014) and Nhung et al. (2015b) confirmed the interpretation of the two-component nature of the Doppler velocity spectrum originally given by Libert et al. (2010). The CSE is axi-symmetric about an axis making an angle of i ~ 30° with the line of sight and projecting on the plane of the sky at a position angle ω ~7° east of north (see also the sketch in Fig. 4). The expansion velocity reaches ~ 8 to 9 km s⁻¹ along the axis – we refer to this part of the CSE as bipolar outflow – and ~ 3 to 4 km s⁻¹ in the plane perpendicular to the axis – we refer to this part of the CSE as equatorial enhancement. The transition from the equator to the poles of the CSE is smooth. Section 4.1 below, using observations of the ¹²CO(2–1) and ¹³CO(2–1) molecular lines, confirms and significantly refines this picture. The right panels of Fig. 20 show projections of the CSE on the plane containing the axis and perpendicular to the plane of the sky, which give a good qualitative idea of the global structure.

Velocity-integrated channel maps of the CO(2–1) and SiO(6–5) observations analyzed in the present article are displayed in Fig. 18. They clearly show the bipolar outflows, inclined toward the observer in the north and receding in the south. We note that the red wings are brighter than the blue wings as a result of absorption (see Sects. 4.1 and 4.2) The SiO-emitting region is seen to be significantly more compact than the CO-emitting region; this is in conformity with observations of many other oxygen-rich AGB stars and is generally interpreted as the result of SiO molecules condensing on dust grains and being ultimately dissociated by the interstellar radiation at some 200 AU from the star, well before CO molecules are dissociated (see e.g., Schöier et al. 2004).

3.6 Temperature and SO₂ abundance

In this section, we use the 11 detected SO₂ lines to derive an approximate temperature and column density of the SO₂ -emitting region by means of a population diagram. Following Goldsmith & Langer (1999), in the optically thin case, the column density of the upper level population N_u of a transition u -> l can be expressed as $$N_u=\frac{8\pi \text{k}{\nu }^{\text{2}}}{\text{h}{c}^3{A}_{ul}}{\displaystyle \int T_b}\text{d}v.$$(3)

N_u is the column density of the upper level population of the transition, k and h are the Boltzmann and Planck constant, respectively, ν is the line frequency, c the speed of light, A_ul is the Einstein coefficient for spontaneous emission of the transition, and ∫ T_b dv is the velocity-integrated main-beam brightness temperature. The latter is converted to the surface brightness distribution of the source S_ν per beam, measured by the interferometer, by means of $$T_b=\frac{{\lambda }^2}{2k{\Omega }_b}{S}_{\nu },$$(4)

where $\lambda =\frac{c}{\nu }$ is the observing wavelength, and ${\Omega }_b=\frac{\pi {\theta }_a{\theta }_b}{4\ln 2}$ with θ_a and θ_b being the major and minor axis of the synthesized beam.

We determine ∫ S(v)dv =: W_I from a circular Gaussian fit to the velocity-integrated emission in the uv-plane, where the integration is taken from (v_lsr,* − 4.5) km s⁻¹ to (v_lsr,* + 4.5) km s⁻¹, that is over the three central channels of the SO₂ lines.

For the population diagram, we then get $$\ln \left(\frac{N_u}{g_u}\right)=\ln \left(\frac{{\zeta }_u{W}_I}{g_u}\right)=\ln \left(\frac{N_{{\text{SO}}_2}}{Q_{{\text{SO}}_2,\text{rot}}}\right)-\frac{E_u}{kT},$$(5)

where we define ζ_u as $${\zeta }_u=\frac{4.784\times {10}^{-7}}{\text{h}c{A}_{ul}}\frac{1}{{\theta }_a{\theta }_b}.$$(6)

In Eq. (5), W_I is expressed in Jy × km s⁻¹, θ_a and θ_b are given in arcsec, E_u is the upper level energy, g_u the statistical weight of the upper level, and T is the excitation temperature. All relevant parameters of the SO₂ transitions used here are listed in Table 2.

In fitting a straight line to the population diagram (shown in Fig. 19) to determine a rotational temperature according to Eq. (5), we assume that the SO₂ level populations are dominated by collisions, i.e., that LTE holds for the rotational excitation of SO₂, and that the lines are optically thin. The assumption of LTE populations may be questionable for SO₂ (see Danilovich et al. 2016), and essentially could result in underestimation of the kinetic gas temperature in the SO₂ -emitting region (cf. the discussion in Goldsmith & Langer 1999, their Sect. 5). The effect on the derived column density is more difficult to assess without a detailed non-LTE modeling. However, we note that our result is consistent with the SO₂ abundance derived for similar objects by means of a comprehensive non-LTE description (see below). On the other hand, the assumption that the lines are optically thin is justified by the results of our XCLASS modeling; see below and Appendix D.

With the temperature of T ≈ 320 K resulting from a linear fit to the population diagram shown in Fig. 19, we get the partition function $Q_{{\text{SO}}_2,\text{rot}}$ (interpolated from values given in the CDMS), from which an effective SO₂ column density of $N_{{\text{SO}}_2}=3.5\times {10}^{15}$ cm⁻² is determined by means of Eq. (5). The SO₂ source is compact (see Table 1 and Fig. B.4, FWHP ≈ 0.5′′, corresponding to 75 AU, or ≈70 R_*), which is consistent with the estimated temperature in this inner (possibly rotating) region. Assuming a mass-loss rate of 1 × 10⁻⁷M_⊙ yr⁻¹ in the equatorial region (cf. Sect. 4.1), an outflow velocity of at maximum 8 km s⁻¹ as indicated bythe line widths (but excluding the high-velocity wings discussed in Sect. 3.3), an inner radius of ≈ 10¹⁴ cm, and an outer radius of that region of ≈10¹⁵ cm, corresponding to 0.5′′, we estimate upper limits of the SO₂ abundance of X(SO₂∕<H>) = 7.3 × 10⁻⁷, or, if hydrogen were completely bound in H₂, X(SO₂/H₂) = 1.5 × 10⁻⁶. For comparison, Danilovich et al. (2016) find an SO₂ abundance of ~ 5 × 10⁻⁶ for the low-mass-loss-rate (~ 1–2 × 10⁻⁷ M_⊙ yr⁻¹) oxygen-rich stars R Dor and W Hya, about a factor of approximately 3 higher than the value found here for the MS-type star RS Cnc. We note that Decin et al. (2008) and De Beck & Olofsson (2020) do not detect SO₂ in the S-type star W Aql. Our result therefore appears to be consistent with an intermediate chemical state of RS Cnc.

In Appendix D, an XCLASS modeling is presented that results in the same SO₂ column density of 3.5 × 10¹⁵ cm⁻² as derived from the population diagram, and in a slightly higher rotational temperature of 350 K, still within the uncertainties of ~ ± 50 K derived here (see Fig. 19). The XCLASS modeling of the 11 SO₂ lines results in an average optical depth of 0.1, confirming that the lines are optically thin.

Fig. 18 Velocity-integrated intensity maps of the 12CO(2–1) line (upper row) and the vibrational ground-state line of SiO(6–5) (lower row), covering three velocity intervals. Left: blue line wing [vlsr,* − 10,vlsr,* − 2] km s−1, middle: line center [vlsr,* − 2,vlsr,* + 2] km s−1, right: red line wing [vlsr,* + 2,vlsr,* + 10] km s−1. North is up and east is to the left. We note the different color scales. Contours are plotted every 5σ for CO and every 20σ for SiO, where (from left to right) 1σ = 14.6, 22.0, 19.2 mJy beam−1 km s−1 for CO(2–1) and 1σ = 11.1, 16.7, 16.4 mJy beam−1 km s−1 for SiO(6–5). The black ellipse in the lower left corner indicates the synthesized beam.

Fig. 19 Population diagram for SO2. The three data points in brackets correspond to the three lowest frequency SO2 lines observed with setup 1, which, due to a different uv coverage, resulted in a comparably larger beam than the setup 2 observations (see Table 2 and Eq. (6)).

Fig. 20 CO data and 3D reconstruction. Left panels: flux density averaged over 0.4′′ < R <2.5′′ in the vz vs. ω plane for 12CO(2–1) (left) and 13CO(2–1) (center-left). Right panels: de-projected flux density projected on the plane perpendicular to the sky and containing the polar outflow axis, shown as black lines, for 12CO(2–1) (center-right) and 13CO(2–1) (right).

4 Discussion

4.1 ¹²CO(2–1) and ¹³CO(2–1): morpho-kinematics of the circumstellar envelope

Rotational CO lines in the vibrational ground state are known to probe the circumstellar envelope of AGB stars up to distances on the 1000 AU scale, where the CO molecules are dissociated by the interstellar UV radiation (Mamon et al. 1988). In this section, we use observations of the ¹²CO(2–1) and ¹³CO(2–1) lines to probe the morpho-kinematics of the wind at distances from the star in excess of ~ 50 AU, in a region where the wind is expected to evolve smoothly toward the constant expansion regime.

Both ¹²CO(2–1) and ¹³CO(2–1) data have a Doppler velocity (v_z) spectrum that clearly shows a two-component structure, where the blueshifted wing of the ¹²CO line is partly absorbed, very similar to the corresponding situation in EP Aqr (Tuan-Anh et al. 2019). CO rotational lines from higher J levels (5–4 and 9–8) observed with Herschel (Danilovich et al. 2015) show triangular profiles consistent with the double-component wind structure that we are proposing.

In the spatial distribution of the CO emission, a clear separation is observed between the equatorial region and the polar outflows. This is illustrated in Fig. 20, where we show the maps of the flux density in the Doppler velocity (v_z) versus position angle (ω) plane. The flux density is integrated over an interval of projected distance (R) from the star between 0.4′′ and 2.5′′. The equatorial region is seen as an intense oscillation at low values of |v_z | while the polar outflows are seen as emission at larger |v_z| values, in phase opposition. From these maps, we estimate, in agreement with earlier findings, that the outflow axis projects on the plane of the sky at a positionangle of ω = 7° ± 5°, at which the equatorial and bipolar outflow oscillations are maximal and minimal, respectively. Figure 20 also shows de-projections of the data cubes projected on the plane perpendicular to the plane of the sky that contains thepolar outflow axis. These are drawn assuming a dependence of the wind velocity on stellar latitude α of the form v_term = 4 + 5sin⁴α km s⁻¹, increasing from 4 to 9 km s⁻¹ from the equator to the poles. The velocity field is assumed to be radial and to have reached its terminal value in the R interval considered here.

In order to interpret these observations, in particular possible differences in the ¹²CO/¹³CO ratio between the equatorial region and polar outflows, we need to take into account the effect of absorption. We choose to consider only regions outside 0.4′′ ~ 60 AU from the star, which allows for major simplifications: we can ignore the very complex kinematics that governs the inner CSE layers, including shocks and possible rotation. As we cannot expect the data to strongly constrain the temperature distribution, we take it of the form T =T₀ exp(−r∕r₀) and fit the observed data cubes by means of a χ² minimization scheme for different values of the (T₀, r₀) pair consistent with earlier estimates (Nhung et al. 2015a). We find that the morpho-kinematics is best described by modeling the polar outflows and equatorial region separately, with a separation at stellar latitude α₀ ~ 30°. We parametrize the radial dependence of the wind velocity as v = v_termr∕(r + r_1∕2): it increases from 0 at r = 0 to v_term at r= ∞, reaching $\frac{1}{2}{v}_{\text{term}}$ at r= r_1∕2. In order to account for different opening angles of the polar outflows and different flaring angles of the equatorial enhancement, we allow, within each region, for latitudinal dependencies of the velocity described as Gaussian functions in sin α in the equatorial region and in cosα in the polar outflows with adjustable Gaussian widths σ.

We do the same for the number densities, which vary radially as ρ₀ (v_term∕v)r⁻² with different ρ₀ values for the equatorial region and the polar outflows, and of course for ¹²CO and ¹³CO. The radiative transfer equation is integratedin a sphere of 10′′ radius and the calculated flux densities are then smeared with the synthesized beam of the observations. Reasonably good fits of the data cubes are obtained given the crudeness of the model (see Fig. 21).

We find that the best fits to the morpho-kinematics of the two lines are relatively insensitive to the exact form assumed for the temperature structure. For r₀ = 2′′, we find that values of T₀ between 70 K and 170 K in the equatorial region and between 100 K and 200 K in the polar outflows are acceptable. Changing the value of r₀ modifies the values of T₀ accordingly but does not improve the quality of the fit. Parameters of the best-fit model include an angle of the polar outflow axis with the line of sight i ≈ 30° as expected, a separation in latitude between polar and equatorial regions, α₀ = 29 ± 4°, terminal velocities of 3–4 km s⁻¹ in the equatorial region, and 8–9 km s⁻¹ in the polar outflows, outflow opening angle and equatorial flaring angle of ~ 70° and ~ 30°, respectively (FWHM = 2.35$\times \arccos \left(\sigma (\rho )\right)$ and $2.35\times \arcsin \left(\sigma (\rho )\right)$, respectively); these are listed in Table 3. The wind is still being accelerated (the escape velocity of a 1.5M_⊙ star at a distance of 150 AU (≈ 1′′ for RS Cnc) is 4.2 km s⁻¹) – having reached half terminal velocity at the inner edge of the observed radial range – earlier in the poles than in the equatorial region. The number densities of CO molecules correspond to mass-loss rates of 1.0 × 10⁻⁷ M_⊙ yr⁻¹ in the equatorial region and 2.0 × 10⁻⁷ M_⊙ yr⁻¹ in the outflows⁶ for a CO/H₂ ratio of 2 × 10⁻⁴. The ¹²CO/¹³CO ratio is measured to be ~20 on average, but larger in the polar outflows (24 ± 2) than in the equatorial region (19 ± 3). This is a barely significant difference, but a similar asymmetry seems to be present in EP Aqr (Tuan-Anh et al. 2019): such a result is unexpected and needs to be confirmed by higher sensitivity observations before being accepted. Indeed, if it were confirmed, this might suggest that the polar outflows are fed in part from material freshly produced in the 3-α process and mixed into the atmosphere by the third dredge-up following a He shell flash. This process would not only increase the ¹²C abundance in the atmosphere, but also the ¹²C/¹³C isotope ratio, as indicated for example by Smith & Lambert (1990); see in particular their Fig. 9.

The presence of an equatorial density enhancement with a rather small flaring angle, ~ 30° FWHM, suggests that it may rather be a disk, which might be expected to be rotating and to have an inner rim. However, the size of the beam is too large to study this reliably. From a close inspection of the ⟨ v_z ⟩ distribution near the star, using the method described in Sect. 3.4 for HCN and SiO, we infer a rotation velocity at r ~0.5′′ of v_rot = |v_z| ∕sini ~ 2.5 km s⁻¹. However, this is a very crude estimate given the size of the beam and the lack of precise knowledge of the morpho-kinematics in the innermost radial range.

The blob of enhanced ¹²CO(2–1) line emission that had been identified in Hoai et al. (2014) as possibly suggesting the presence of a companion is seen on the channel maps in the v_lsr = 6 to 7 km s⁻¹ range as an elongation in the west–northwest to east–southeast direction (Fig. 22). Projections of the data cube on different planes (v_z vs. ω, v_z vs. R, and ω vs. R) in its neighborhood show that it can be described as a pair of elongations at position angles of ~ 120° and ~ 270°, the latter being significantly more intense than the former and covering a broad range of R between 1 and 2′′. While these features provide no justification for a possible identification of a companion, they cannot be used either as arguments against the presence of an unobserved companion.

Fig. 21 Best-fit results for CO. Left panels: same as in Fig. 20, but for the model best fit. Right panels: modeled vz spectra (red) compared with observations (black) for 12CO(2–1) (center-right) and 13CO(2–1) (right).

4.2 SiO(5–4) and SiO(6–5): evidence for strong absorption

In the present section, we compare the SiO(5–4) and SiO(6–5) line emission with the ¹²CO results described in the previous section. In contrast to CO, the SiO emission does not resolve the equatorial region from the polar outflows, as illustrated in the leftmost panels of Fig. 23. Part of the reason for this is the much smaller radial range being probed, as illustrated in the right panel of Fig. 23. As mentioned earlier in Sect. 3.5, the radial extent of the SiO emission is often significantly smaller than that of the CO emission, which is usually interpreted as evidence for the progressive condensation of SiO molecules on dust grains (e.g., Schöier et al. 2004). This would cause the progressive decline of the SiO/CO ratio observed in the right panel of Fig. 23 up to R ~ 1.5′′, followed by a more abrupt cut-off around R ~2′′ caused by the dissociation of the SiO molecules by the interstellar UV radiation.

Another striking difference is the presence of high-velocity wings in the SiO data close to the star, that are mostly absent in the CO data (see Figs. 24 and 16). This may be understood as the SiO data probing the close neighborhood of the star much more efficiently, where pulsation shocks and possibly convection cells are known to play an important role, as demonstrated by a host of observations at shorter wavelengths, from far-infrared to near-UV, and in agreement with theoretical modeling (see, e.g., Hinkle et al. 1982; Winters et al. 2000a; Richter et al. 2003; Nowotny et al. 2005; Freytag et al. 2017; Montez et al. 2017; Höfner & Olofsson 2018, and references therein). Also worth noting in Fig. 24 is the significant absorption in front of the stellar disk seen in the SiO data at ~ − 6 km s⁻¹, reminiscent of similar observations in stars such as R Dor (Nhung et al. 2021), W Hya (Takigawa et al. 2017), and o Ceti (Wong et al. 2016). Of these three stars, two are found to be surrounded by an optically thick SiO layer that extends well beyond the stellar disk, while in the third, o Ceti, SiO emission decreases abruptly beyond some 10 to 50 AU from the center of the star. In order to understand the mechanism of condensation of SiO molecules on dust grains, it is therefore important to study the SiO emission region in oxygen-rich AGB stars. Accordingly, we devote the remainder of the present section to this task. In a first part, we make some qualitative comments that help with unraveling the general picture; in a second part, we account for radiative transfer using the model of the morpho-kinematics developed for CO emission in the preceding section. The morpho-kinematics of the SiO(5–4) and SiO(6–5) lines are so similar that we cannot expect their comparison to be very sensitive to temperature. At temperature T, in the optically thin LTE approximation, the ratio R_T of the SiO(5–4)/SiO(6–5) line emission depends only on the ratio of the Einstein coefficients, 5.2 × 10⁻⁴∕9.1 × 10⁻⁴ = 0.57, and the ratio of the level populations,(2J + 1)(exp(−E_u∕T), where the upper level energies E_u are 31.3 K and 43.8 K, respectively(line parameters from the CDMS, based on Müller et al. 2013).

Hence, we have $$R_T=0.48\exp (-31.3/T)/\exp (-43.8/T),$$(7)

from which we obtain $$T[\text{K}]=12.5/(\ln (R_T)+\mathrm{0.74.}$$(8)

Figure 25 (left) displays the observed dependence of R_T on R. If the SiO layer is optically thin, the value of R_T provides a direct measure of the temperature; if it is optically thick, the emission probes only the outer part of the SiO layer and measures a temperature that is representative of larger values of r. When R increases, the observed value of R_T remains nearly constant at ~0.65 up to R ~1′′ and then increases to a value of ~0.8 at 1.5′′, corresponding to temperatures of ~40 K and ~ 24 K, respectively. The value of 40 K, associated with R <1′′, is close to the lower value of the temperature obtained from the CO model for r =1′′ (T_CO[K] = 70exp(−r∕2′′) = 42 K for r =1′′) in Sect. 4.1. This suggests that for low values of the projected distance from the star, R, the distance in space from the star, r, effectively probed by the SiO emission stays approximately constant and of the order of 1′′. Indeed, the strong self-absorption causes the region probed by the SiO emission to be confined beyond r ~1′′ (corresponding to τ_SiO ~ 1), but not to probe the volume beneath this surface. For larger values of the projected distance R >1′′, we expect the distance in 3D space, r, of the effective SiO-emitting surface to increase progressively as a result of the elongation of the emission volume along the line of poles, a pure geometrical effect (Fig. 20).

Qualitatively, we describe this trend in Fig. 25 (left) by assuming that r stays at the 1′′ level for R <1′′ and increases approximately from 1′′ to 2.5′′ when the projected R increases from 1′′ to 2′′. Such a trend, when translated in terms of temperatures using the CO model relation T_CO[K] = 70exp (− r∕2′′) gives a fair description of the observed dependence of R_T on R given the important approximations that have been made (Fig. 25 left).

The ratio ε∕τ of the emissivity to optical depth is a measure of the lower limit of the emission of a self-absorbing layer. In the LTE approximation and to first order in ΔE∕T, where ΔE is the energy of the transition, (ε∕τ)∕(Tν²) is a constant, where ν is the frequency of the emission. Therefore, at the same temperature, (ε∕τ)[SiO(5–4)] = 0.88(ε∕τ)[¹²CO(2–1)]. The value of(ε∕τ) is T∕26.3 Jy arcsec⁻² for CO and T∕29.5 Jy arcsec⁻² for SiO. Taking as reference T = 100 K for SiO and 50 K for CO, the corresponding values of (ε∕τ) are 1.9 Jy arcsec⁻² for CO and 3.4 Jy arcsec⁻² for SiO. The observed values (Fig. 25 center-right) are ~ 2.5 and 6 Jy arcsec⁻², respectively at R ~ 0. This is less than a factor 2 above the reference values, showing that the optical thickness is close to that of the self-absorption regime. A similar result has been observed for other AGB stars with mass-loss rates on the order of 10⁻⁷ M_⊙ yr⁻¹ (R Dor, Nhung et al. 2021, W Hya, Takigawa et al. 2017).

To gain further insight into this issue, having obtained qualitative evidence for strong absorption, we need to properly account for radiative transfer. To do so, we use the model developed in Sect. 4.1 to describe the morpho-kinematics of the CO component of the CSE: in Fig. 26 we compare Doppler velocity spectra of SiO and CO in octants of position angle (the first one defined as 7° < ω < 52°, then counter-clockwise for increasing numbers, as indicated on Fig. 23) in a ring 1′′ < R < 1.5′′. This range of R is far enough from the star to be relatively independent of modeling the star neighborhood, which contributes to the high-velocity wings of SiO spectra but our current data do not offer significant constraints on this region. We only show SiO(5–4) but the results are essentially the same for SiO(6–5). The CO data are seen to trace the polar outflows from the equatorial region, the SiO data donot. To adapt the CO model to a description of the SiO data, we keep all parameters as listed in Table 3 with the only exception being the radial profile: we divide the CO density parameter ρ₀ by three, and we further multiply it by exp(−r∕1.6), and set it to zero beyond r = 1.6′′, in qualitative accordance with the distributions displayed in the right panel of Fig. 23. These modifications result in a ratio of SiO/CO = 0.18 at r = 1′′, in agreement with the estimate of Van de Sande et al. (2018a), who obtain a SiO/CO ratio of 0.17–0.19 for this type of star.

Apart from this major modification of the radial distribution of the density, we use the same morpho-kinematics and same temperaturedistribution as for the best fit to the CO data. In order to adapt to SiO, we change the parameters defining the transition, in particular with Einstein coefficients three orders of magnitude larger than for the CO transitions. The result, displayed in Fig. 26, describes the data surprisingly well given the crudeness of the exercise; in addition to offering an understanding of the large difference between CO and SiO emission, it provides confirmation of the validity of the description of the morpho-kinematics given by the model. What happens is that the SiO emission in the modeled annular ring between 1′′ and 1.5′′ probes only the outer layer of the (confined) SiO volume because of the very large absorption, while CO probes the lower density gas at much larger distances from the star where the equatorial region and the polar outflows are well separated. We refrain from attempting to adjust parameters to obtain a better fit: modeling the inner region could not be done reliably with the limited angular resolution of the present data.

Another illustration of the stronger absorption of the SiO line compared to the CO line is displayed in the right panel of Fig. 25. It compares the (projected) radial distribution of the ratio ²⁸SiO(5–4)/²⁹SiO(5–4) with that of the ratio ¹²CO(2–1)/¹³CO(2–1). The ²⁹SiO data are obtained with a similar beam size to the ²⁸SiO and CO data, but with a spectral resolution of only 3 km s⁻¹, preventing us from a detailed study of their Doppler velocity spectra. We expect the ²⁸SiO/²⁹SiO isotopic ratio to be ~12 or lower (De Beck & Olofsson 2018, 2020), but this line ratio is reached only at large enough values of R, where the optical depth has decreased sufficiently for absorption to be less important. In the limit of complete self-absorption (i.e., emission at the ε∕τ level), the line ratio is unity. In contrast, the CO ratio is much less affected by absorption, staying at ~ 80% of the value evaluated in Sect. 4.1 (~20) over the whole radial range.

Fig. 22 Enhanced 12CO(2–1) emission at vlsr ~ 6.75 km s−1. Upper panels: channel maps of the 12CO(2–1) line emissionbetween 6.25 and 7.25 km s−1. Lower panels: projection of the flux density multiplied by R on the vlsr vs. R plane for 250° < ω <330° (left), on thevlsr vs. ω plane for 0.5 < R < 2 arcsec (middle), and on the ω vs. R plane for 5.75 <vlsr < 7.25 km s−1 (right).

Fig. 23 vz vs. ω maps for 0.5′′ < R <1.5′′. Left: SiO(5–4), middle: SiO(6–5), octant intervals in ω as discussed in the text are indicated. Right: intensity ratio SiO(5–4)/12CO(2–1) (black) and SiO(6–5)/12CO(2–1) (red) as a function of R for |vz | < 8 km s−1.

Fig. 24 Doppler velocity spectra integrated over 0.2 < R < 2.5 arcsec (black) and within R < 0.2 arcsec (red) for SiO(5–4), SiO(6–5), and 12CO(2–1) from left to right.

Fig. 25 Analysis of SiO and comparison with CO. Left: ratio RT of SiO(5–4) to SiO(6–5) line emitted fluxes as a function of R for |vz | < 8 km s−1. Corresponding temperature values using the relation T[K] = 12.5∕(ln(RT) + 0.74 are shown in red on the scale of the ordinate. The red dashed line indicates constant RT of 0.64 (T ~ 42 K) for R ≤1′′ and a linear increase of RT to 0.89 (T ~ 20 K) from 1′′ to 2′′. Center-left: Doppler velocity spectrum of the SiO(5–4) data for R <0.1′′. Center-right: R distribution of the SiO(5–4) (red) and 12CO(2–1) (black) emission averaged over position angle ω, measured in Jy arcsec−2, and averaged over − 2 < vz < 6 km s−1 as shown by the arrow in the center-left panel. Right: R dependence of the 12CO(2–1)/13CO(2–1) ratio (black) and of the 28SiO(5–4)/29SiO(5–4) ratio (red).

Fig. 26 12CO(2–1) and 28SiO(5–4) Doppler velocity spectra in the ring 1′′ < R <1.5′′. The data are shown in black and the best-fit model results are shown in red (see text). Octant intervals in position angle ω are numbered as 1 = [7°, 52°], 2 = [52°, 97°], …8 = [322°, 7°], these are indicated on the left and central panels of Fig. 23.

4.3 MHD wind models

An important feature of the morpho-kinematics of RS Cnc is the latitudinal dependence of the velocity, which is larger towards the pole than along the equatorial plane. On the other hand, the modeling of Hoai et al. (2014) leads to an almost spheroidal distribution of matter (their Fig. 6), which implies an isotropic radiation field within the circumstellar shell (assuming position coupling between dust and gas). In these conditions, a driver other than radiation is needed to accelerate matter preferentially along the polar axis.

Axi-symmetric models have been developed on the assumption that a companion disturbs the gravitational field (Theuns & Jorissen 1993; Mastrodemos & Morris 1999; Decin et al. 2020). These models have been successful in explaining several observed features, such as the presence of spiral structures, although the preferential acceleration along the polar axis seems more difficult to reproduce.

In RS Cnc, we do not find evidence for the presence of a companion. While this negative result does not necessarily mean that there is no companion, we are considering here another mechanism that could play a role: the presence of a magnetic field. In this section, we argue that magnetic fields are one possible candidate to explain a significant pole-to-equator velocity contrast of > 1.

Indeed, magnetic fields are an appealing way to explain some of the axi-symmetric shapes of planetary nebulae, such as polar jets or equatorial winds. Early simulations of winds with dipole fields and rotation showed that excess magnetic pressure could be able to repel the wind towards both the equator and/or the poles (Washimi & Shibata 1993). Similar simulations also showed the formation of equatorial disks with enhanced outflow velocities in the equatorial region (Matt et al. 2000). In a more toroidal configuration, Matt & Balick (2004) showed that both jets and equatorial disks could be produced depending on the magnitude of the rotation. Finally, García-Segura et al. (2005) showed that magnetically driven winds yield strongly anisotropic outflows with highly collimated polar jets.

Observational polarimetric studies (Greaves 2002) revealed ordered magnetic fields in planetary nebulae, with variousdegrees of toroidal configurations depending on the target. Using a handful of SiO masers (Herpin et al. 2006)or CN Zeeman measurements, Duthu et al. (2017) later attempted to further characterize the magnitude and the radial dependence of the magnetic field in the winds of AGB stars. These latter authors concluded that the field is on the order of a few Gauss near the stellar surface, and consistent with a 1∕r dependence on the distance r from the star (see their Fig. 6).

In this section, we report results of simplistic calculations which integrate magnetized fluid parcel trajectories from the surface of an AGB star up to 10 stellar radii. We integrate the acceleration due to gravitational and Lorentz forces over time (pressure gradients become quickly negligible after the sonic point is crossed): $$\ddot{r}=(\Gamma -1)\frac{G{M}_{*}}{r^3}r+\frac{1}{\rho }J\times B,$$(9)

where r is the position of the fluid parcel, Γ is the ratio between the radiative and the gravitational force, G is the universal gravity constant, M_* is the stellar mass, and $J=\frac{1}{4\pi }\nabla \times B$ is the current vector. The mass density ρ and the magnetic field B are prescribed, while we solve for the position and velocity of the fluid parcels. We parametrize our equations with $(\Gamma -1)G{M}_{*}=v_{\infty }{R}_{*}^2$ and $\rho =\frac{\dot{M}}{4\pi r^2{\text{v}}_{\infty }}$. We assume Ṁ = 1.24 × 10⁻⁷ M_⊙ yr⁻¹, and R_* = 1.6 × 10¹³ cm as reasonable values for RS Cnc (Hoai et al. 2014). Our choice of parametrization allows us to easily explore nondimensional values of the parameters independently of the absolute observational constraints. After we obtain a suitable contrast between polar and equatorial velocities, we retrieve the physical value for the velocity scale – here v_∞ = 5.6 km s⁻¹ – in order to obtain a given polar outflow velocity of 8 km s⁻¹ at r = 10 R_*. The initial velocity vector is set with a small uniform radial velocity and we probe starting trajectories at the surface with a uniformly distributed initial latitude.

This simple setup allows us to quickly investigate various magnetic field configurations. Figure 27 displaysthe trajectories in the meridional plane obtained for a toroidal magnetic field with a 1∕r^1.1 decline from B_* = 0.5 G at the stellar surface: $B_{\varphi }=B_{*}\cos (\theta )/{(r/R_{*})}^{1.1}$ where θ is the latitude (1∕r^1.1 gives a more pronounced velocity contrast between pole and equator than 1∕r). The Lorentz force J × B in this case is directed toward the symmetry axis and acts as a focusing agent. All trajectories eventually end up on the axis where they would presumably launch a jet: this focuses mass loss towards the poles. Figure 28 shows the resulting “terminal” velocity at 10 stellar radii, where we have separated the latitudinal and the radial components. The flow velocities at this radius are dominated by their radial component, but with a clearly slower wind at the equator compared to the pole, as indicated by observations of, for instance, RS Cnc and EP Aqr. We note that rotation tends to produce the opposite effect (faster radial flow at equator compared to poles, due to centrifugal acceleration). We find similar behavior for toroidal magnetic field configurations closer to what Matt & Balick (2004) found ($B_{\varphi }=3B_{*}\cos (\theta )\sin {(\theta )}^2/{(r/R_{*})}^2$ with a 1∕r² dependence and a concentration of B at intermediate latitudes). We investigated additional dipolar fields which are able to generate some amount of rotation in the wind. Thanks to the versatility of the present setting, we can quickly explore various configurations, linear combinations between them, several magnetic field decay exponents, but have not yet investigated them systematically. The purpose of our investigation here is simply to show that a magnetic field is a valid candidate to produce pole-to-equator velocity ratios of significantly greater than 1. Finally, we note that our termination radius of r = 10 R_* is arbitrary, and a pertinent match to the observations could be considered at various distances depending on where the given tracer is expected to be concentrated. These crude models are still far from quantitatively matching the observational constraints from RS Cnc or EP Aqr, which require a broader polar outflow and a thinner and denser equatorial disk. This could be adjusted by providing a sharper toroidal magnetic barrier to funnel the wind at the appropriate places. However, such fine-tuning would overcome the limits of this crude exercise which still lacks self-consistency as the density profile remains radial and unaffected by the magnetic constraints (themselves blind to the wind), and shocks generated by the crossings of trajectories at the polar axis are not accounted for. We plan to investigate further with such simplified models in future work as they might provide a useful means to constrain the magnetic field configuration from the morpho- kinetics.

Fig. 27 Magnetized fluid parcel trajectories from a simplified model of RS Cnc (see text). Red trajectories meet on the polar axis, where they will likely generate a jet.

Fig. 28 Radial velocity at 10 stellar radii depending on the latitude for the same simplified model of magnetized wind as shown in Fig. 27.

4.4 HCN in M- to S-type stars

Cherchneff (2006) recognized that the formation of CN/HCN depends on the high activation barrier of the H + C₂ → CH + C reaction, followed by rapid CN formation via N + CH →CN+H. Their abundance therefore depends on thermal excursions in shocks, or inhomogeneities of temperature. In addition, both the total rate of formation of the pair CN/HCN and the respective share between CN and HCN depend on the H/H₂ ratio which itself depends on out-of-equilibrium chemistry because of the slow conversion between H and H₂. Cherchneff (2006) was thus able to show that the shocks produced by the pulsations close to the stellar photosphere could produce highly increased yields of the HCN molecule in M or S-type stars, despite their high O/C ratio. We note here that magnetic fields produce shocks on the symmetry axis which might also help to boost HCN production in the polar jet.

HCN has long been detected and surveyed in M-type and S-type stars (e.g., Deguchi & Goldsmith 1985; Lindqvist et al. 1988; Bujarrabal et al. 1994; Olofsson et al. 1998; Schöier et al. 2013). RS Cnc on the other hand has never been detected in ground-state or vibrationally excited HCN despite various observational efforts with different telescopes (Lucas et al. 1988; Sopka et al. 1989; Nercessian et al. 1989; Lindqvist et al. 1992; Bujarrabal et al. 1994; Bieging & Latter 1994; Olofsson et al. 1998). Adopting a mass-loss rate of 3 × 10⁻⁷ M_⊙ yr⁻¹, Bujarrabal et al. (1994) estimated an upper limit to the HCN abundance in RS Cnc of 4.5 × 10⁻⁷. This upper limit becomes 1.35 × 10⁻⁶ if we adopt Ṁ = 1 × 10⁻⁷ M_⊙ yr⁻¹. Schöier et al. (2013) presented a comprehensive analysis of the HCN abundance in a sample of 59 AGB stars, including 25 carbon-rich, 19 S-type, and 25 M-type stars, by means of a non-LTE radiative transfer modeling. For M-type and S-type stars, these authors derived a median HCN/H₂ abundance of order 1 × 10⁻⁷ and 7 × 10⁻⁷, respectively, with a large spread between ~5 × 10⁻⁸ and ~5 × 10⁻⁶.

By an XCLASS modeling of the HCN(3–2) line detected here (see Appendix D), and assuming a rotational temperature of 350 K, we derive an HCN column density of ~ 1.6 × 10¹⁵ cm⁻²; see Fig. 29. With the same assumptions as those made in Sect. 3.6, this translates to HCN abundances of X(HCN/⟨H ⟩) = 3.3 × 10⁻⁷, or, if hydrogen were completely bound in H₂, X(HCN/H₂) = 6.6 × 10⁻⁷. This abundance perfectly fits in the range found by Schöier et al. (2013) for M- to S-type stars and is also consistent with the upper limit derived by Bujarrabal et al. (1994) corrected for the mass-loss rate.

Clumpy and porous winds help UV photons to penetrate closer to the star, thus photodissociating CO and N₂ to release moreof the C/NO and the N/CS pairs of reactants, which both produce CN: this process was shown to considerably enhance the HCN abundances close to the star (Van de Sande et al. 2018b, 2020). Rather than uniformly distributed random clumpiness and porosity, one can also imagine ordered density distributions in the wind, which could let UV photons penetrate through low-density channels. These structures could sometimes be hard to witness due to line of sight confusion. In fact, our simplified magnetized wind models (see Fig. 27 with scarcity of trajectories around the equator) or more sophisticated magnetized wind models (e.g., Matt et al. 2000; Washimi & Shibata 1993) allow for lower column-density channels at certain angles. These may result in increased HCN abundance, but chemical post-processing in magnetized models will be necessary to assess whether this is a viable interpretation of the observations.

Fig. 29 Map of the HCN column density as derived by an XCLASS modeling (see Appendix D). The black ellipse in the lower left corner indicates the synthesized beam.

5 Conclusions

Using NOEMA equipped with PolyFiX, we obtained high-spatial resolution (~ 0.3′′) images of RS Cnc in several lines of different molecules. We detect, and in most cases are able to map, 32 lines of 13 molecules and isotopologs (CO, ¹³CO, SiO, ²⁹SiO, SO, ³⁴SO, SO₂, H₂O, HCN, H¹³CN, PN), including several transitions from vibrationally excited states, and a tentative identification of Si¹⁷O and possibly ²⁹Si¹⁷O. HCN as well as millimeter vibrationally excited H₂O, SO, SO₂, and PN and their isotopologs are first detections in RS Cnc. From their first-moment maps, some of the lines, SiO(v = 1,6–5), HCN, SO, SO₂, show signs of rotation in the close vicinity of the star. A population diagram analysis for the 11 observed SO₂ lines provides a rotational temperature of about 320 K in the region that shows signs of rotation. Temperatures of this order are also found from an XCLASS modeling of the SO₂ lines. For SO₂ and HCN, we find column densities from the XCLASS modeling of $N_{{\text{SO}}_2}~3.5\times {10}^{15}$ cm⁻² and N_HCN ~ 1.6 × 10¹⁵ cm⁻², which translate to abundance ratios of X(SO₂∕H₂) = 1.5 × 10⁻⁶ and X(HCN/H₂) = 6.6 × 10⁻⁷, respectively, well within the range expected for an MS-type star.

We find broad wings in the spectral line profiles of vibrational ground-state transitions of SiO and SO and in first vibrationally excited transitions of SiO, which indicate radial velocities of about twice the terminal outflow velocity as probed by CO. As high velocities very close to the star are also seen in similar objects, such as EP Aqr, o Cet, and R Dor, the presence of these broad line wings calls for a mechanism common to the class of pulsating AGB stars. We interpret these high-velocity line wings as the imprints of pulsation shocks acting in the very inner region around these stars.

The spatially resolved images allow us to trace the morpho-kinematics of the wind around RS Cnc at different scales. In the inner part (<0.5′′, or 75 AU), we find a rotating structure well traced by the less abundant molecules (HCN, SO, SO₂), and by SiO in (v = 1) lines. Outside 75 AU, we find an expanding axi-symmetric outflow, with velocities ~ 4 km s⁻¹ in the equatorial plane, and ~ 9 km s⁻¹ along the polar axis. This polar axis coincides with the axis of the internal rotating structure. A model that fits the data cubes obtained on the ¹²CO(2–1), ¹³CO(2–1) and SiO(v = 0, 5–4 and 6–5) lines gives a mass-loss rate of 1 × 10⁻⁷ M_⊙ yr⁻¹ for the equatorial region (latitude <30°) and of 2 × 10⁻⁷ M_⊙ yr⁻¹ for the polar outflows (latitude >30°). The ¹²CO/¹³CO ratio is measured to be ~20 on average, 24 ± 2 in the polar outflows and 19 ± 3 in the equatorial region.

Although we cannot exclude the possibility that an unseen stellar or substellar companion shapes the circumstellar environmentof RS Cnc, we also consider the possibility of a magnetic field playing this role. In particular, a toroidal magnetic field configuration would provide a mechanism able to produce the significant velocity contrast between high polar-outflow velocities and low expansion velocities in the equatorial region that is observed in RS Cnc and other similar stars.

Acknowledgements

We thank the staff at the NOEMA and Pico Veleta observatories for their support of these observations. The authorsare grateful to the anonymous referee for a very detailed and valuable report that helped improving the presentation of the material. This work is based on observations carried out under project numbers W16BE, D17AE, W19AX with the IRAM NOEMA interferometer and under project ID 136-19 with the IRAM 30m telescope. IRAM is supported by INSU/CNRS (France), MPG (Germany) and IGN (Spain). The Ha Noi team acknowledges financial support from the World Laboratory, the Odon Vallet Foundation and VNSC. This research is funded in part by the Vietnam National Foundation for Science and Technology Development (NAFOSTED) under grant number 103.99-2019.368. This work was supported by the Programme National “Physique et Chimie du Milieu Interstellaire” (PCMI) of CNRS/INSU with INC/INP co-funded by CEA and CNES. This work has made use of data from the European Space Agency (ESA) mission Gaia (https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing and Analysis Consortium (DPAC, https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has been provided by national institutions, in particular the institutions participating in the Gaia Multilateral Agreement.

Appendix A Resolved out flux

The effect on the flux of filtering out large-scale structure with the interferometer is shown in Fig. A.1 for the CO and ¹³CO lines. These two are the only lines discussed in this paper that are affected by the short-spacing problem.

Fig. A.1 Effect of missing short spacing information on the CO lines. Left: 12CO(2–1). Right: 13CO(2-1). A-configuration and D-configuration are merged, the spectral resolution is 0.5 km s−1 and the CO emission is integrated over the central 22′′ × 22′′, i.e., over the full field of view of the NOEMA antennas at 230 GHz. Black profiles: OTF data (i.e., the short-spacing information that is filtered out by the interferometer) are added. Red profiles: A+D configuration interferometer data, only.

Appendix B Intensity maps

Here we present velocity-integrated intensity maps in three velocity ranges for HCN and H¹³CN (Fig. B.1), the four detected SO lines (Fig. B.3), and three out of the 11 SO₂ lines (Fig. B.4). These are the SO₂ line with the lowest upper level energy, the strongest SO₂ line detectedhere, and the SO₂ line with the highest upper level energy, respectively. All these nine lines display kinematic structure ineast-west direction. Also shown are zeroth moment maps for the (unresolved) H₂O lines (Fig. B.2) and for the (weak) PN line (Fig. B.5).

Fig. B.1 Velocity-integrated intensity maps of the HCN(3-2) (upper row) and H13CN(3-2) lines (lower row), covering three velocity intervals. Left: Blue line wing [vlsr,* − 10,vlsr,* − 2] km s−1, Middle: line center [vlsr,* − 2,vlsr,* + 2] km s−1, Right: Red line wing [vlsr,* + 2,vlsr,* + 10] km s−1. North is up and east is to the left. We note the different color scales. Contours are plotted every 10σ for HCN andevery 3σ for H13CN, where (from left to right) 1σ = 5.7, 7.2, 7.2 mJy/beam km s−1 for HCN(3-2) and 1σ = 5.2, 5.0, 5.6 mJy/beam km s−1 for H13CN(3-2). The black ellipse in the lower left corner indicates the synthesized beam. We note that the HCN maps were produced using robust weighting, whereas for the H13CN maps, we applied natural weighting.

Fig. B.2 Zeroth moment maps of the H2O 232 GHz (left) and H2O 263 GHz (right) lines. North is up and east is to the left. We note the different color scales. Contours are plotted every 5σ where 1σ = 13.8 mJy/beam⋅ km s−1 for the 232 GHz line and 1σ = 13.92 mJy/beam⋅ km s−1 for the 263 GHz line. The black ellipse in the lower left corner indicates the synthesized beam.

Fig. B.3 Velocity-integrated intensity maps of the SO(6(6)-5(5)) (upper row), SO(5(5)-4(4)) (second row), SO(6(5)-5(4)) (third row) and SO(7(6)-6(5)) (lower row) lines covering three velocity intervals. Left: Blue line wing [vlsr,* − 10,vlsr,* − 2] km s−1, Middle: Linecenter [vlsr,* − 2,vlsr,* + 2] km s−1, Right: Red line wing [vlsr,* + 2,vlsr,* + 10] km s−1. North is up and east is to the left. We note the different color scales. Contours are plotted every 10σ, where (from left to right) 1σ = 10.0, 10.4, 10.3 mJy/beam⋅ km s−1 for SO(6(6)-5(5)), 1σ = 7.3, 7.6, 7.4 mJy/beam⋅ km s−1 for SO(5(5)-4(4)), 1σ = 7.6, 8.1, 7.4 mJy/beam⋅ km s−1 for SO(6(5)-5(4)), and 1σ = 8.7, 9.6, 8.7 mJy/beam⋅ km s−1 for SO(7(6)-6(5)). The black ellipse in the lower left corner indicates the synthesized beam.

Fig. B.4 Velocity-integrated intensity maps of SO2 (9(3, 7)-9(2, 8)) (lowest Eu, upper row), SO2 (14(0,14)-13(1,13)) (strongest SO2 line, second row), and SO2 (34(4,30)-34(3,31)) (highest Eu, lower row), covering three velocity intervals. Left: Blue line wing [vlsr,* − 10,vlsr,* − 2] km s−1, Middle: line center [vlsr,* − 2,vlsr,* + 2] km s−1, Right: Red line wing [vlsr,* + 2,vlsr,* + 10] km s−1. North is up and east is to the left. We note the different color scales. Contours are plotted every 3σ, where (from left to right) 1σ = 6.5, 6.7, 6.6 mJy/beam⋅ km s−1 for SO2(9(3, 7)- 9(2, 8)) 1σ = 5.9, 6.2, 6.3 mJy/beam⋅ km s−1 for SO2(14(0,14)-13(1,13)), 1σ = 7.8, 7.8, 7.7 mJy/beam⋅ km s−1 for SO2 (34(4,30)-34(3,31)). The black ellipse in the lower left corner indicates the synthesized beam.

Fig. B.5 Zeroth moment map of the PN(N=5-4,J=6-5) line. Contours are plotted every 3σ, where 1σ = 9.6 mJy/beam⋅ km s−1. North is up and east is to the left. The black ellipse in the lower left corner indicates the synthesized beam.

Appendix C SO₂ line profiles

In this section, we present the line profiles of all the 11 SO₂ lines detected with our setups (see Fig. C.1).

Fig. C.1 Profiles of all the 11 detected SO2 lines. A-configuration and D-configuration data are merged, the spectral resolution is 3 km s−1, and the emission is integrated over the central 2′′ × 2′′ square aperture. The lines are ordered by decreasing frequency from top to bottom and left to right.

Appendix D XCLASS modeling of HCN and SO₂

The HCN(3-2) line and the 11 detected SO₂ emission lines were modeled using the eXtended CASA Line Analysis Software Suite (XCLASS⁷, Möller et al. 2017). XCLASS models and fits molecular lines by solving the 1D radiative transfer equation with the assumption of LTE and of an isothermal source. Here, 1D means that the radiative transfer equation is integrated along the line of sight. Spectral lines are fitted with Gaussian profiles, and optical depth effects and source size are considered in the calculations. Molecular properties (e.g., Einstein coefficients, partition functions, etc.) are taken from an embedded SQLite database containing entries fromthe Cologne Database for Molecular Spectroscopy (CDMS, Müller et al. 2001, 2005) and from the Jet Propulsion Laboratory database (JPL, Pickett et al. (1998)) using the Virtual Atomic and Molecular Data Center (VAMDC, Endres et al. 2016). The fit parameter set for each line component consists of the source size θ_source, the rotation temperatureT_rot, the total column density N_tot, the line width Δv, and the velocity offset v_off (given here in the LSR system).

The XCLASS package offers various algorithms to find the best-fit parameters by minimizing the χ² value, and here we utilized the Levenberg-Marquardt (LM) method. To obtain maps of the physical parameters, we use the myXCLASSMapFIt function to fit HCN(3-2) and the 11 detected SO₂ emission lines (see Fig. C.1) pixel by pixel.

For SO₂, we modeled and fitted 11 lines simultaneously with a threshold of 18σ and a single Gaussian component. All fit parameters are regarded as free parameters in the fitting process. The XCLASS models for SO₂ result in a temperature of T_rot ~ 350 K (Fig. D.1), somewhat higher than the result from our population diagram analysis, but within the error bars (see Sect. 3.6). It also results in an average line optical depth of 0.1, confirming that the assumption of the lines being optically thin is justified when constructing the population diagram. Assuming LTE, the derived rotation temperature equals the kinetic gas temperature T_kin. For the SO₂ column density, the XCLASS modeling results in a value of $N_{{\text{SO}}_2}~3.5\times {10}^{15}$ cm⁻², almost identical to the result from the population diagram analysis presented in Sect. 3.6, resulting in an abundance X(SO₂∕H₂) = 1.5 × 10⁻⁶. The respective results are shown in Figs. D.1 and D.2.

For HCN(3-2), we applied the same threshold of 18σ as for SO₂ and fitted Gaussian components to the line profiles on each pixel, taking into account the hyperfine structure of the line. The velocity map and line widths displayed on Fig. D.3 are due to the intrinsic (thermal and rotational) broadening only, whereas the profile shown in Fig. D.4 represents the sum of the hyperfine components of HCN(3-2). As only one HCN rotational line is available, the HCN rotational temperature cannot be determined. We therefore fixed the rotation temperature at a value of 350 K, the same temperature as we find from the SO₂ modeling based on the similar emission region of SO₂ and HCN, see Figs. 17, B.1, and B.4. The HCN results are displayed in Fig. D.3. In particular, we derive an HCN column density of ~ 1.6 × 10¹⁵ cm⁻², which translates to an HCN abundance of X(HCN/H₂) = 6.6 × 10⁻⁷.

Fig. D.1 Maps of total column density, rotation temperature, velocity offset (in the LSR system), and line width for SO2 are derived with a threshold of 18σ and fitting one Gaussian component to the line profiles.

Fig. D.2 Profiles of the 11 SO2 emission lines (in black) extracted on the central pixel and XCLASS modeled lines (in red) at the position of the continuum at RA = 09:10:38.77 and Dec = +30:57:46.68 (see Sect. 3.1) as marked on the upper diagram.

Fig. D.3 Maps of total column density, velocity offset (= vlsr), and line width of HCN(3-2) (see text for details). The black ellipse in the lower left corner of each map indicates the synthesized beam.

Fig. D.4 Profile of the HCN(3-2) emission line (in black) extracted on the central pixel and XCLASS modeled line (in red) at the position of the continuum at RA = 09:10:38.77 and Dec = +30:57:46.68 (see Sect. 3.1) as marked on the upper diagram of Fig. D.3.

References

All Tables

All Figures

Fig. 1 Overview of the frequency ranges observed with PolyFiX using two spectral setups (setup1: red and setup2: blue, respectively). Lower diagrams: zoom onto the individual spectra covering 7.8 GHz each. Upper row: setup1, lower row: setup2. Thecentral 20 MHz at the border between inner and outer baseband are blanked out, i.e., set to zero, as this region is contaminated by the LO2 separation of the 8 GHz-wide IF in the IF processor (“LO2 zone”).

In the text

	Fig. 2 Continuum map around 247 GHz from A-configuration. Contours are plotted in 100σ steps, where 1σ is 47.6 μJy beam−1. The synthesized beam is indicated in the lower left corner.
In the text

	Fig. 3 CO line profiles, showing a two-component structure. Left: 12CO(2–1). Right: 13CO(2–1). A-configuration and D-configuration are merged, OTF data are added, and the spectral resolution is 0.5 km s−1. The CO emission is integrated over the central 22′′ × 22′′, i.e., over thefull field of view of the NOEMA antennas at 230 GHz.
In the text

	Fig. 4 Sketch of the geometrical structure of the wind components as inferred from the current data (see Sect. 4.1). The sketch is not to scale: there is a smooth transition between the equatorial enhancement and the polar outflows.
In the text

Fig. 5 Profiles of SiO ground-state and first vibrationally excited state lines. Left: SiO(6–5): upper: v = 1, lower: v = 0. A-configuration and D-configuration merged. Right: SiO(5–4): upper: v = 1, D-configuration (black) and A-configuration (red), lower: v = 0, A-configuration and D-configuration merged. The spectral resolution is 1 km s−1 for (v = 1) and 0.5 km s−1 for the (v = 0) lines, respectively. The emission is integrated over the central 5′′ × 5′′ aperture.

In the text

	Fig. 6 Profile around the SiO(v = 2, J = 5–4) line frequency. A-configuration and D-configuration merged. The spectral resolution is 1 km s−1 and the emission is integrated over the central 1′′ × 1′′ aperture.
In the text

Fig. 7 Line profiles of SiO isotopologs. Upper left: profile of the 247.482 GHz line, possibly 29Si17O(6–5); D-configuration, only. Lower left: Si17O(6–5): D-configuration, only (line was not covered in A-configuration). Right: 29SiO(5–4); A-configuration and D-configuration merged. The spectral resolution in all cases is 3 km s−1 and the emission is integrated over the central 5′′ × 5′′ aperture.

In the text

	Fig. 8 HCN line profiles. Left: HCN(3–2); A-configuration and D-configuration merged with a spectral resolution of 0.5 km s−1. Right: H13CN(3–2); A-configuration and D-configuration merged with a spectral resolution of 3 km s−1. The emission of both lines is integrated over the central 2′′ × 2′′.
In the text

	Fig. 9 H2O line profiles. Left: H2O line at 232.687 GHz. Right: H2O line at 263.451 GHz. Data are merged from A-configuration and D-configuration, the spectral resolution is 3 km s−1, and the emission of both lines is integrated over the central 1′′ × 1′′ aperture.
In the text

	Fig. 10 Profiles of the four detected SO lines, with A-configuration and D-configuration merged. The spectral resolution is 3 km s−1 and the emission is integrated over the central 2′′ × 2′′ aperture.
In the text

	Fig. 11 Profiles of the two 34SO lines detected here, with A-configuration and D-configuration merged. The spectral resolution is 3 km s−1 and the emission is integrated over the central 2′′ × 2′′ aperture.
In the text

	Fig. 12 Profile of SO2 (14(0,14)–13(1,13)) with A-configuration and D-configuration merged, a spectral resolution of 3 km s−1, and emissionintegrated over the central 2′′ × 2′′ aperture.
In the text

	Fig. 13 Profile of PN(N = 5–4, J = 6–5) with A-configuration and D-configuration merged, a spectral resolution of 3 km s−1, and emissionintegrated over the central 2′′ × 2′′ aperture.
In the text

	Fig. 14 Linewings in SiO(5–4) and SiO(6–5) compared to CO(2–1). The emission is integrated over the central 5′′ × 5′′ aperture.
In the text

	Fig. 15 High velocities close to the line of sight as seen in SiO. PV maps are shown in the Vz vs. R plane for SiO(5–4) (left) and SiO(6–5) (right). The horizontal black line indicates the wind terminal velocity as traced in CO and the white scale bar indicates the spatial resolution. $R=\sqrt{{(\delta \text{Dec})}^2+{(\delta \text{RA})}^2}$, |Vz | = |vlsr − vlsr,*|.
In the text

	Fig. 16 Line profiles of different molecules on a logarithmic intensity scale. Gaussian profiles are shown for comparison, FWHM = 10 km s−1 for the ground-state lines of all molecules, and FWHM = 14 km s−1 for the (v = 1) lines of SiO. All observed profiles are integrated over R < 0.2′′.
In the text

	Fig. 17 First-moment maps of different lines, indicating a possibly rotating structure (see Sect. 3.4). Left: HCN(3–2), middle left: SO2 (14(0,14)–13(1,13)), middle right: SiO(v = 1,6–5), right: SO(7(6)–6(5)). The black ellipses indicate the synthesized beam.
In the text

Fig. 18 Velocity-integrated intensity maps of the 12CO(2–1) line (upper row) and the vibrational ground-state line of SiO(6–5) (lower row), covering three velocity intervals. Left: blue line wing [vlsr,* − 10,vlsr,* − 2] km s−1, middle: line center [vlsr,* − 2,vlsr,* + 2] km s−1, right: red line wing [vlsr,* + 2,vlsr,* + 10] km s−1. North is up and east is to the left. We note the different color scales. Contours are plotted every 5σ for CO and every 20σ for SiO, where (from left to right) 1σ = 14.6, 22.0, 19.2 mJy beam−1 km s−1 for CO(2–1) and 1σ = 11.1, 16.7, 16.4 mJy beam−1 km s−1 for SiO(6–5). The black ellipse in the lower left corner indicates the synthesized beam.

In the text

	Fig. 19 Population diagram for SO2. The three data points in brackets correspond to the three lowest frequency SO2 lines observed with setup 1, which, due to a different uv coverage, resulted in a comparably larger beam than the setup 2 observations (see Table 2 and Eq. (6)).
In the text

	Fig. 20 CO data and 3D reconstruction. Left panels: flux density averaged over 0.4′′ < R <2.5′′ in the vz vs. ω plane for 12CO(2–1) (left) and 13CO(2–1) (center-left). Right panels: de-projected flux density projected on the plane perpendicular to the sky and containing the polar outflow axis, shown as black lines, for 12CO(2–1) (center-right) and 13CO(2–1) (right).
In the text

	Fig. 21 Best-fit results for CO. Left panels: same as in Fig. 20, but for the model best fit. Right panels: modeled vz spectra (red) compared with observations (black) for 12CO(2–1) (center-right) and 13CO(2–1) (right).
In the text

	Fig. 22 Enhanced 12CO(2–1) emission at vlsr ~ 6.75 km s−1. Upper panels: channel maps of the 12CO(2–1) line emissionbetween 6.25 and 7.25 km s−1. Lower panels: projection of the flux density multiplied by R on the vlsr vs. R plane for 250° < ω <330° (left), on thevlsr vs. ω plane for 0.5 < R < 2 arcsec (middle), and on the ω vs. R plane for 5.75 <vlsr < 7.25 km s−1 (right).
In the text

	Fig. 23 vz vs. ω maps for 0.5′′ < R <1.5′′. Left: SiO(5–4), middle: SiO(6–5), octant intervals in ω as discussed in the text are indicated. Right: intensity ratio SiO(5–4)/12CO(2–1) (black) and SiO(6–5)/12CO(2–1) (red) as a function of R for |vz | < 8 km s−1.
In the text

	Fig. 24 Doppler velocity spectra integrated over 0.2 < R < 2.5 arcsec (black) and within R < 0.2 arcsec (red) for SiO(5–4), SiO(6–5), and 12CO(2–1) from left to right.
In the text

Fig. 25 Analysis of SiO and comparison with CO. Left: ratio RT of SiO(5–4) to SiO(6–5) line emitted fluxes as a function of R for |vz | < 8 km s−1. Corresponding temperature values using the relation T[K] = 12.5∕(ln(RT) + 0.74 are shown in red on the scale of the ordinate. The red dashed line indicates constant RT of 0.64 (T ~ 42 K) for R ≤1′′ and a linear increase of RT to 0.89 (T ~ 20 K) from 1′′ to 2′′. Center-left: Doppler velocity spectrum of the SiO(5–4) data for R <0.1′′. Center-right: R distribution of the SiO(5–4) (red) and 12CO(2–1) (black) emission averaged over position angle ω, measured in Jy arcsec−2, and averaged over − 2 < vz < 6 km s−1 as shown by the arrow in the center-left panel. Right: R dependence of the 12CO(2–1)/13CO(2–1) ratio (black) and of the 28SiO(5–4)/29SiO(5–4) ratio (red).

In the text

	Fig. 26 12CO(2–1) and 28SiO(5–4) Doppler velocity spectra in the ring 1′′ < R <1.5′′. The data are shown in black and the best-fit model results are shown in red (see text). Octant intervals in position angle ω are numbered as 1 = [7°, 52°], 2 = [52°, 97°], …8 = [322°, 7°], these are indicated on the left and central panels of Fig. 23.
In the text

	Fig. 27 Magnetized fluid parcel trajectories from a simplified model of RS Cnc (see text). Red trajectories meet on the polar axis, where they will likely generate a jet.
In the text

	Fig. 28 Radial velocity at 10 stellar radii depending on the latitude for the same simplified model of magnetized wind as shown in Fig. 27.
In the text

	Fig. 29 Map of the HCN column density as derived by an XCLASS modeling (see Appendix D). The black ellipse in the lower left corner indicates the synthesized beam.
In the text

Fig. A.1 Effect of missing short spacing information on the CO lines. Left: 12CO(2–1). Right: 13CO(2-1). A-configuration and D-configuration are merged, the spectral resolution is 0.5 km s−1 and the CO emission is integrated over the central 22′′ × 22′′, i.e., over the full field of view of the NOEMA antennas at 230 GHz. Black profiles: OTF data (i.e., the short-spacing information that is filtered out by the interferometer) are added. Red profiles: A+D configuration interferometer data, only.

In the text

Fig. B.1 Velocity-integrated intensity maps of the HCN(3-2) (upper row) and H13CN(3-2) lines (lower row), covering three velocity intervals. Left: Blue line wing [vlsr,* − 10,vlsr,* − 2] km s−1, Middle: line center [vlsr,* − 2,vlsr,* + 2] km s−1, Right: Red line wing [vlsr,* + 2,vlsr,* + 10] km s−1. North is up and east is to the left. We note the different color scales. Contours are plotted every 10σ for HCN andevery 3σ for H13CN, where (from left to right) 1σ = 5.7, 7.2, 7.2 mJy/beam km s−1 for HCN(3-2) and 1σ = 5.2, 5.0, 5.6 mJy/beam km s−1 for H13CN(3-2). The black ellipse in the lower left corner indicates the synthesized beam. We note that the HCN maps were produced using robust weighting, whereas for the H13CN maps, we applied natural weighting.

In the text

	Fig. B.2 Zeroth moment maps of the H2O 232 GHz (left) and H2O 263 GHz (right) lines. North is up and east is to the left. We note the different color scales. Contours are plotted every 5σ where 1σ = 13.8 mJy/beam⋅ km s−1 for the 232 GHz line and 1σ = 13.92 mJy/beam⋅ km s−1 for the 263 GHz line. The black ellipse in the lower left corner indicates the synthesized beam.
In the text

Fig. B.3 Velocity-integrated intensity maps of the SO(6(6)-5(5)) (upper row), SO(5(5)-4(4)) (second row), SO(6(5)-5(4)) (third row) and SO(7(6)-6(5)) (lower row) lines covering three velocity intervals. Left: Blue line wing [vlsr,* − 10,vlsr,* − 2] km s−1, Middle: Linecenter [vlsr,* − 2,vlsr,* + 2] km s−1, Right: Red line wing [vlsr,* + 2,vlsr,* + 10] km s−1. North is up and east is to the left. We note the different color scales. Contours are plotted every 10σ, where (from left to right) 1σ = 10.0, 10.4, 10.3 mJy/beam⋅ km s−1 for SO(6(6)-5(5)), 1σ = 7.3, 7.6, 7.4 mJy/beam⋅ km s−1 for SO(5(5)-4(4)), 1σ = 7.6, 8.1, 7.4 mJy/beam⋅ km s−1 for SO(6(5)-5(4)), and 1σ = 8.7, 9.6, 8.7 mJy/beam⋅ km s−1 for SO(7(6)-6(5)). The black ellipse in the lower left corner indicates the synthesized beam.

In the text

Fig. B.4 Velocity-integrated intensity maps of SO2 (9(3, 7)-9(2, 8)) (lowest Eu, upper row), SO2 (14(0,14)-13(1,13)) (strongest SO2 line, second row), and SO2 (34(4,30)-34(3,31)) (highest Eu, lower row), covering three velocity intervals. Left: Blue line wing [vlsr,* − 10,vlsr,* − 2] km s−1, Middle: line center [vlsr,* − 2,vlsr,* + 2] km s−1, Right: Red line wing [vlsr,* + 2,vlsr,* + 10] km s−1. North is up and east is to the left. We note the different color scales. Contours are plotted every 3σ, where (from left to right) 1σ = 6.5, 6.7, 6.6 mJy/beam⋅ km s−1 for SO2(9(3, 7)- 9(2, 8)) 1σ = 5.9, 6.2, 6.3 mJy/beam⋅ km s−1 for SO2(14(0,14)-13(1,13)), 1σ = 7.8, 7.8, 7.7 mJy/beam⋅ km s−1 for SO2 (34(4,30)-34(3,31)). The black ellipse in the lower left corner indicates the synthesized beam.

In the text

	Fig. B.5 Zeroth moment map of the PN(N=5-4,J=6-5) line. Contours are plotted every 3σ, where 1σ = 9.6 mJy/beam⋅ km s−1. North is up and east is to the left. The black ellipse in the lower left corner indicates the synthesized beam.
In the text

	Fig. C.1 Profiles of all the 11 detected SO2 lines. A-configuration and D-configuration data are merged, the spectral resolution is 3 km s−1, and the emission is integrated over the central 2′′ × 2′′ square aperture. The lines are ordered by decreasing frequency from top to bottom and left to right.
In the text

	Fig. D.1 Maps of total column density, rotation temperature, velocity offset (in the LSR system), and line width for SO2 are derived with a threshold of 18σ and fitting one Gaussian component to the line profiles.
In the text

	Fig. D.2 Profiles of the 11 SO2 emission lines (in black) extracted on the central pixel and XCLASS modeled lines (in red) at the position of the continuum at RA = 09:10:38.77 and Dec = +30:57:46.68 (see Sect. 3.1) as marked on the upper diagram.
In the text

	Fig. D.3 Maps of total column density, velocity offset (= vlsr), and line width of HCN(3-2) (see text for details). The black ellipse in the lower left corner of each map indicates the synthesized beam.
In the text

	Fig. D.4 Profile of the HCN(3-2) emission line (in black) extracted on the central pixel and XCLASS modeled line (in red) at the position of the continuum at RA = 09:10:38.77 and Dec = +30:57:46.68 (see Sect. 3.1) as marked on the upper diagram of Fig. D.3.
In the text