Atmosphere of Betelgeuse before and during the Great Dimming event revealed by tomography

Despite being the best studied red supergiant star in our Galaxy, the physics behind the photometric variability and mass loss of Betelgeuse is poorly understood. Moreover, recently the star has experienced an unusual fading with its visual magnitude reaching a historical minimum. We investigate the nature of this event with the help of a recently developed tomographic technique. Tomography allows us to probe different depths in the stellar atmosphere and to recover the corresponding disk-averaged velocity field. We apply the tomographic method to a few-year time-series high-resolution spectroscopic observations of Betelgeuse in order to relate its atmospheric dynamics to the photometric variability. Our results show that a sudden increase of the molecular opacity in the cooler upper atmosphere of Betelgeuse is likely the reason of the observed unusual decrease of the star's brightness.


Introduction
Red supergiant (RSG) stars represent the late stage in the evolution of stars with initial masses larger than 8 M 1 before they end their lives in spectacular supernova explosions. Despite the fact that RSGs were extensively studied during the last decades, important properties such as photometric variability and mass loss are still poorly constrained. Understanding these properties is crucial for a broad range of astrophysical questions including the chemical enrichment of the Galaxy, supernova progenitors, and the extragalactic distance scale (Levesque 2017).
The reduced HERMES spectra are only available in electronic form at the CDS via anonymous ftp to cdsarc.u-strasbg.fr (130.79.128.5) or via http://cdsweb.u-strasbg.fr/cgi-bin/ qcat?J/A+A/ 1 The 8 M threshold corresponds to the minimum mass required for carbon ignition in the stellar core. This value is not well established and depends on the treatment of convection in stellar evolution models (Höfner & Olofsson 2018). Betelgeuse is one of the closest (222 +48 −34 pc, Harper et al. 2017a), brightest (0.0 -1.6 mag in V-band 2 ), and best studied RSGs in our Galaxy. Its atmospheric parameters are listed in Table B.1. Interferometric observations of Betelgeuse by Haubois et al. (2009), Chiavassa et al. (2010), Ohnaka et al. (2011), and Montargès et al. (2016) have shown evidence for the presence of large convective cells in its atmosphere, which is in accordance with predictions by Schwarzschild (1975). The photometric periodicity of Betelgeuse is characterized by two main periods: a short one of about 400 days and a long one of about 2000 days (Kiss et al. 2006). Mechanisms such as atmospheric pulsations in the fundamental or low-overtone modes and oscillations excited by convection cells were invoked to explain the short photometric period (Kiss et al. 2006). The long photometric variations were proposed to be due to binarity or magnetic activity (Wood et al. 2004), or turnover of giant convection cells (Stothers 2010). 2 According to the American Association of Variable Star Observers (AAVSO) database (Kafka 2018); https://aavso.org/ Article number, page 1 of 12 arXiv:2104.08105v1 [astro-ph.SR] 16 Apr 2021 A&A proofs: manuscript no. 39801corr Josselin & Plez (2007) and Gray (2008) analyzed highresolution (R ∼ 40 000 and 100 000, respectively) spectroscopic observations of Betelgeuse and detected time-variable Doppler shifts and line strengths. In particular, Gray (2008) reported a hysteresis-like behavior of the temperature (represented by a line-depth ratio) as a function of line-core velocity. It was shown that hysteresis loops turn counterclockwise, that is to say an increase in temperature is followed by a rise of the matter, its subsequent cooling, and falling down. This behavior was explained by the presence of large convective cells in the stellar atmosphere, which, in turn, dominate the 400-day variations of Betelgeuse. However, according to Kravchenko et al. (2019, Paper II), stationary convection cannot account for the phase shifts between temperature and velocity variations observed by Gray (2008). Paper II shows that hysteresis loops are linked to the propagating shock waves, which originate from acoustic waves generated by a disturbance in the convective flow in the stellar interior.
Between October 2019 and February 2020, Betelgeuse has experienced an unusual decrease in its visual brightness by about 1 mag (relative to the most recent maximum in September 2019; , see also the ESO press release 3 ), which is referred to as the Great Dimming event. Different hypotheses have been proposed to interpret the nature of this unprecedented fading. They include the formation of a dust cloud due to a recent mass-loss episode (Levesque & Massey 2020;Dupree et al. 2020;Cotton et al. 2020, Montargès et al., submitted), thermal changes in the photosphere accompanying a temperature drop of ∼200 K or the presence of a dark spot covering ∼50% of the stellar surface (Dharmawardena et al. 2020), or a critical transition in Betelgeuse seen as a complex dynamical system undergoing a regime shift (George et al. 2020).
The present letter aims to relate the atmospheric dynamics of Betelgeuse to its photometric variability, with the main focus on the dimming event. For this purpose, the tomographic method developed by Alvarez et al. (2001) and Kravchenko et al. (2018, hereafter Paper I) was applied to time series of high-resolution spectra of Betelgeuse.
The letter is structured as follows. Section 2 describes the spectral time series of Betelgeuse. The application to Betelgeuse of the tomographic method, which is summarized in Appendix A, is presented in Sect. 3. The resulting interpretation of the dimming event and our conclusions are discussed in Sect. 4.

Observations
Betelgeuse was observed with the high-resolution fiber-fed cross-dispersed echelle spectrograph HERMES (Raskin et al. 2011) mounted on the 1.2m Mercator telescope at the Roque de Los Muchachos Observatory, La Palma (Spain). The spectral resolution of HERMES is R = 86000, and the wavelength coverage is from 3800 to 9000 Å. In total, 37 high-resolution, high signal-to-noise ratio (S/N ∼ 200-400 at 5000 Å) spectra were obtained between November 2015 and September 2020, that is, covering a time span of approximately 5 years. The epochs of the HERMES observations may be located on the Betelgeuse AAVSO light curve, which can be seen in the top panel of Fig. 1. The HERMES spectra were reduced with an automated pipeline, merging the different orders and correcting for the blaze function of the echelle grating as well as for the Earth motion around the solar-system barycenter. The long-term stability of the HER-MES velocity frame is on the order of 80 m s −1 , which was estimated from a long-term monitoring of RV standard stars (Jorissen et al. 2016

Results
In this section, we apply the tomographic method of Paper I (summarized in Appendix A) to the high-resolution HERMES spectra of Betelgeuse to investigate its atmospheric motions and to search for a possible link with photometric variability, with the main focus on the recent dimming event. The results are compared to those obtained by Paper II, where the tomographic method of Paper I was applied for the first time to time-series observations of another RSG star (µ Cep).

Cross-correlation functions (CCFs)
First, we cross-correlated the HERMES spectra of Betelgeuse with the set of five tomographic masks (see Table A.1) constructed in Paper II. Among these, the mask denoted C1 probes the innermost photospheric layer, whereas mask C5 probes the outermost layer. The resulting CCFs (displayed in Figs. B.3 -B.8) show asymmetries in all masks. Moreover, at a few epochs (marked with a star symbol in the top panel of Fig. 1), the CCF profiles in masks C1 and C2 are characterized by two components shifted by 10-15 km s −1 with respect to each other (Fig. B.1). To be considered meaningful, the second component must have a depth larger than 3σ, where σ is the standard deviation in the flat part of the CCFs measuring the correlation noise.
For instance, at JD 2458172 (i.e., February 22, 2018), the CCF becomes asymmetric with a faint extra component appearing on the blue wing of the main peak in masks C1 and C2. After an episode where the CCF in mask C1 shows no correlation peak (JD 2458372.7), the CCF reappears around JD 2458492 (i.e., January 8, 2019), but this time with a secondary peak on the red wing of the main peak. The additional CCF component only appears in the innermost masks and remains weak, never equalling the contrast of the main CCF peak. According to Paper II, the intensity of a CCF component is related to the size of the corresponding emitting surface on the star. For example, if a large fraction of the stellar surface is covered by rising (falling) material, then the CCF is characterized by a strong blue-shifted (red-shifted) peak. A similar scenario is observed in long-period variable (LPV) stars of Mira-type, where line doubling happens around maximum light (Alvarez et al. 2000). It is associated with the appearance and passage of a spherically-symmetric shock wave across the atmosphere (the Schwarzschild scenario 4 , Schwarzschild 1954; Alvarez et al. 2000). In our observations of Betelgeuse, line-doubling events in February 2018 and January 2019 (also marked as star symbols in the second panel of Fig. 1) occur on the rising part of the light curve (top panel of Fig. 1). Paper II shows that, despite the absence of a Schwarzschild scenario in the temporal and spatial evolution of CCFs as in Mira stars, this behavior nevertheless illustrates the emergence of shock waves in the atmosphere of RSGs, but with lower amplitudes (compared to Miras) and more asymmetric shapes. The absence of similar line-doubling events at earlier epochs of our observations (i.e., during light cycles 1 and 2 in Fig. 1) means that those shock waves -if present -had amplitudes that were too low to be distinguishable on CCFs (see Figs. B.3 and B.4).

Radial velocities
As a next step, we derived radial velocities (RVs) from the CCFs using the two methods described in Paper II. First, we fit the CCFs with a single-or multiple-component Gaussian profile using the detection of extrema (DOE) tool (Merle et al. 2017). In this way, the RV information is extracted for each CCF compo-nent separately. The resulting RVs for each tomographic mask are shown in the second panel of Fig. 1. There, RVs from the single-peaked CCFs and the strongest CCF components of double-peaked CCFs are connected with lines. The unconnected symbols correspond to the RVs of the secondary components of double-peaked CCFs. CCFs in mask C1 are often either too shallow or too asymmetric to warrant a reliable Gaussian fit (see Figs. B.5 and B.6). Therefore, the analysis below is, in several instances, restricted to masks C2-C5. Second, RVs were computed as the first (raw) moment M 1 of CCFs (considering the CCF as a probability distribution, see Paper II). The resulting M 1 velocities for masks C2-C5 are shown in the third panel of Fig. 1. Both RV and M 1 velocities are shown in Fig. 1 on a relative scale, with the zero-point set at v CoM = 20.7 km s −1 corresponding to the center-of-mass (CoM) velocity of Betelgeuse (±0.3 km s −1 ), which is represented by the horizontal line in Fig. 1. This CoM velocity is the average between the values of Harper et al. (2017b) and Kervella et al. (2018a).
As seen in the second and third panels in Fig. 1, the maximum peak-to-peak amplitude of the RV and M 1 velocity variations is slightly larger in the innermost masks. A comparison with the light curve reveals that for the cycles preceding the dimming event, it amounts to 3 − 5 km s −1 . However, during dimming, the velocity amplitude of Betelgeuse reached ∼ 6 km s −1 in mask C5 and even ∼ 10 km s −1 in masks C2-C3. To interpret this behavior (which cannot be ascribed to instrumental errors, which are considerably smaller, as discussed in Sect. 2), we defined three atmospheric velocity gradients between mask C5 and masks C2, C3, and C4 as of atmospheric regions over which the gradient is taken. The derived velocity gradients are displayed in the fourth panel in Fig. 1 and reveal the following. First, the velocity gradients are generally nonzero, thus revealing the existence of compression or expansion within the photosphere, and they are stronger for the inner atmospheric layers. This explains the differences in the peak-to-peak amplitudes of the RV and M 1 velocity variations between the masks. Second, the velocity gradients very closely follow the light curve variations, with expansion corresponding to dimming. In particular, during the dimming event, the velocity gradients are the strongest in all masks, which in turn explains the increase in the RV and M 1 amplitudes. Finally, during the dimming event, the velocity gradients reach values that are at least a factor of 2 larger (negative) compared to any other epoch of our observations in all masks, indicating rapid expansion of a portion of the atmosphere of Betelgeuse. Finally, the fact that the velocity gradients recovered their normal values in September 2020 (JD 2459100) clearly indicates that the dynamical consequences of the dimming event were then definitely over. Figure 1 also reveals that, as in µ Cep, the RV and M 1 velocities of Betelgeuse vary with the same periodicity as the light curve, but with a phase shift such that the light-curve maxima occur at later epochs than the velocity maxima. The phase shift amounts to about 0.2 for mask C5 and 0.3 for mask C2, and thus decreases toward the outer atmospheric layers. In particular, the difference between masks C2 and C5 is about 40 days.

The FWHM of CCFs
The fifth panel of Fig. 1 displays the full width at half maximum (FWHM) of CCFs in different masks, which reflects the amplitude of the velocity dispersion inside a given atmospheric layer. The FWHM values increase by a factor of more than two when going from the inner to the outer atmospheric layers. The average FWHM in mask C5 is ∼ 39 km s −1 while it is ∼ 15 km s −1 in mask C1, whereas their uncertainty is only 0.1 -0.2 km s −1 .
The FWHMs of CCFs do not show a clear correlation with the light curve variations. Nevertheless, there appears to be a phase lag between the FWHM variations in the different masks. This behavior is especially apparent between JD 2458492 (January 8, 2019) and JD 2458603 (April 29, 2019) in Fig. 1, thus somewhat before the dimming event, where a strong increase in the FWHM first occurs in mask C2 and then gradually appears in the subsequent masks. The phase lag between the masks is the same as observed in the M 1 velocity variations (i.e., 40 days between masks C2 and C5).
The strong increase in the FWHM visible at the beginning of Cycle 4 results from the apparition of line doubling (January 2019) in mask C1 (see the second panel of Fig. 1) at earlier epochs. The line doubling in mask C1 indicates the emergence of a shock wave at the base of the atmosphere. As the shock propagates outward, it gives rise to broader CCF profiles in the outer atmospheric layers. Thus, the phase lag of the strong increase in the FWHM during Cycle 4 along the different masks reflects the shock propagation. This interpretation is supported by the fact that the line doubling in mask C1 appears when the stellar brightness increases. In addition, the comparison with the temporal evolution of velocity M 1 reveals that, for a given mask, the maximum of the FWHM occurs at the exact time when the M 1 velocity turns negative (indicating rising material), as can be seen by comparing panels 3 and 5 of Fig. 1 around January-March 2019. Finally, the identical phase lag between the M 1 velocity and the FWHM variations points at the same mechanism being responsible for their behavior. As pointed out in Paper II, physical processes linked to shocks are responsible for the M 1 velocity variations.
Interestingly, the above-described FWHM behavior occurring about 6 months before the Great Dimming episode was pre-

TiO-band temperatures
For each HERMES observation, we derived a temperature 5 for Betelgeuse using the method described in Van Eck et al. (2017) and Paper II. In short, the technique is based on the computation of the strengths of the TiO bands having their bandheads at 5847, 6158, 6187, 7054, and 7126 Å. First, the relation between the TiO-band strength and the effective temperature of a model was derived for a grid of 1D MARCS model atmospheres (Gustafsson et al. 2008). Then, the strengths of the TiO bands were computed for the HERMES spectra and the temperatures were deduced using the 1D-model calibration. The absolute error of our temperature determination is about 50 K, corresponding to the half-step of the effective temperature grid for MARCS models. However, the relative precision is better as may be judged from 5 According to Davies et al. (2013), temperatures derived using the 1D hydrostatic model spectral fit of the TiO bands are not identical to the effective temperature of a star, which is best derived from a fit of the whole spectral energy distribution (SED). Indeed, the temperature in the outer layers (where the TiO lines form) of more realistic 3D models is lower than in their 1D counterparts, leading to increased TiO absorption. Thus, the term "temperature" used in this letter refers only to the "TiOband" temperature. the smooth variations observed in the top panel of Fig. 1, with steps much smaller than 50 K.
The temperatures derived in the manner described above are shown in the top panel of Fig. 1 (see also Table B.3). A comparison with the light curve reveals that the variations of the TiOband temperatures are generally in phase with the photometric variations. Thus, the phase shifts observed between the photometric and velocity (and the FWHM of CCF) variations apply to the temperature variations as well. The only discrepancy between the temperature and photometric variations is observed at epochs JD 2 458 900.5 and JD 2 458 907.5, when the star recovers from the dimming event, but the temperature continues to decrease.
We compared our temperatures with those obtained by previous studies. Levesque et al. (2005) and Levesque & Massey (2020), using moderate-resolution observations (4 − 8 Å, corresponding to R = 625 − 1250 at 5000 Å), proceeded in a similar way, but by using TiO bands at 4761, 4954, 5167, 5448, 5847, and 6158 Å bands, thus the reddest ones were lacking. These authors found that the temperature of Betelgeuse amounts to 3650 ± 25 K and 3600 ± 25 K when V ∼ 0.5 and V ∼ 1.6, respectively. Guinan & Wasatonic (2020) measured the Betelgeuse surface temperature from the strength of the TiO 7190 Å bandhead and obtained 3650 K and 3565 K for V ∼ 0.3 and V ∼ 1.6, respectively. Harper et al. (2020) derived TiO-indices from the Wing narrow-band photometric filters centered at the 7190 Å TiO band and two continuum regions. Their measured TiO-temperatures amount to 3645 ± 25 K and 3520 ± 25 K in September 2019 (V ∼ 0.6) and February 2020 (V ∼ 1.6), respectively. According to the bottom panel of Fig. 1, for Cycle 4, at V ∼ 0.3, 0.5, and 1.6, we find temperatures of 3645 K, 3620 K, and 3542 K, respectively. These are in good agreement (within the errorbars) with those derived by , Harper et al. (2020), and Levesque & Massey (2020). The small discrepancy with the latter authors can probably be attributed to the different set of TiO bands used for the analysis, since HER-MES spectra give access to the stronger red bands at 6187, 7054, and 7126 Å. These bands are actually the most sensitive to temperature, as shown in Fig. B.2. The comparison of our temperatures with the photometric variations shown in the top panel of Fig. 1 reveals that the peakto-peak amplitude of the temperature variations is the largest for the light cycles 1 and 4. In particular, during the rising part of Cycle 1, the temperature increases by about 81 K. Then during Cycle 4, corresponding to the dimming event, the temperature decreases by about 105 K, which is comparable to the 125 K temperature drop found by Harper et al. 2020. Such a change in TiO temperature, if accompanied by a similar change in T eff , is not enough to explain the photometric variation of Betelgeuse. Indeed, the synthetic photometry of MARCS models reveals that a temperature drop of ∆T = 100 K corresponds to a decrease of ∆V = 0.5 mag at most, whereas the observed variation in V was more than twice as large (1.3 mag). in Fig. 2 resemble those obtained earlier by Gray (2008) and those for µ Cep from Paper II. According to Table B.2, the characteristic timescale of hysteresis loops closely matches the short photometric period of Betelgeuse (400 days, Kiss et al. 2006). A similar match between the timescales of the hysteresis loops and the light variations was reported for µ Cep in Paper II and for Betelgeuse by Gray (2008). Finally, the hysteresis loops from Fig. 2, from Gray (2008) and from Paper II show qualitative similarity in terms of their velocity and temperature ranges (see Table B.2 in this paper and Table 4 in Paper II). This further confirms that the mechanism linked to shock-wave propagation is responsible for the velocity, temperature, and photometric variations in Betelgeuse.

Discussion and conclusions
According to the results presented in Sect. 3, the main dynamical specificities prior to the dimming event are as follows: (i) the line-doubling episodes of February 2018 and January 2019 (second panel of Fig. 1); (ii) a concomittant increase in all CCF FWHMs (fifth panel of Fig. 1); and (iii) the relative velocity M 1 −v CoM turning negative during the duration of Cycle 4, which is the signature of an outflow (second and third panels of Fig. 1, as confirmed by Dupree et al. 2020, who also report UV line and continuum emissions in September -November 2019 at-tributed to this strong outflow). Assuming that the RV stayed at −7 km s −1 for 240 d in 2019 (second panel in Fig. 1 and Fig. 2 of Dupree et al. 2020), the distance travelled by the outflowing matter amounts to 210 R , which represents a small fraction of the photospheric radius of ∼ 1000 R (Ohnaka et al. 2011, combined with Gaia DR2), which is certainly too small for the outflowing matter to reach the dust-formation region at > ∼ 2 stellar radii (counted from the stellar center). Between January and March 2019, the M 1 − v CoM velocity turned negative for each mask at the exact time the FWHM maximum reached the corresponding layer. Moreover, the upward propagation of the broad CCFs during Cycle 4 (fifth panel of Fig. 1) reveals the propagation of the shock wave through the atmosphere. Thus, it is clear that a strong shock has appeared at the bottom of the photosphere in January 2019 and it has moved upward, which first led to a strong compression (v C 5 − v C i > 0) and then, as the dimming started, to a strong expansion (v C 5 − v C i < 0) of the atmosphere. Observational evidence is strong to link this compression (expansion) in the photosphere with the increase (decrease) in the TiO temperature. In particular, during the strong compression episode (seen in January -May 2019 in panel 4 of Fig. 1), the temperature shows a clear maximum as expected when a perfect gas is being compressed (e.g., Fig. 6 of Liljegren et al. 2016). Conversely, during the expansion episode, the temperature decreases as expected. The emergence of a similar shock in February 2018, with even larger upward velocities than in January 2019, is likely; although, it is unfortunately not documented by observations as well.
Thus the extraordinary strong upflow occurring during the duration of Cycle 4 results from the conjunction of several effects: (i) the succession of two shocks along our line-of-sight, the second one amplifying the effect of the first one; (ii) a specific convective pattern with mostly rising matter on the visible hemisphere of Betelgeuse, with a large cell 6 resulting from random fluctuations in the stellar interior; and (iii) the previous phenomena amplified by the coincidence with outward motion present at this phase of the ∼ 400 d pulsation cycle. We interpret the Great Dimming event as a consequence of this strong outflow which has caused a cooling of the outer layers, resulting in a sudden increase in molecular opacity. This phenomenon has been described in the literature as "molecular plumes" rising from the photosphere of supergiants (Kervella et al. 2018b, Freytag et al., in prep.), or "molecular reservoirs" (Harper et al. 2020), a nonspherical version of the "molsphere" observed around AGB stars (Freytag et al. 2017).
This scenario is supported by the VLT/SPHERE measurements (see footnote 3), which revealed an obscuration in the southern half of Betelgeuse disk in the optical in January 2020. Whether or not this obscuration was caused by dust is a crucial point. Contrary to expectations, when new dust is produced, no significant flux increase has been recorded in the IR (Harper et al. 2020, and references therein). Moreover this obscuration event has been recorded (20% dimming) in the submillimeter range (Dharmawardena et al. 2020); however, new dust is not expected to significantly change the spectrum at wavelengths larger than 100 µm. Therefore, similarly to Harper et al. (2020), we conclude that dust does not appear to be a crucial ingredient to account for the Great Dimming.

Appendix A: The tomographic method
The tomographic method allows one to probe different geometrical depths in the stellar atmosphere and to recover the surface-averaged line-of-sight velocity at those various geometrical depths (Alvarez et al. 2001, and Paper I). The technique is based on sorting different spectral lines according to their formation depth, which is provided by the maximum of the line-depth contribution function (CF; see Fig. 1 of Paper I). The formation depth of spectral lines is expressed in an optical-depth scale computed at a reference wavelength (e.g., λ = 5000 Å) built from continuum opacities only according to Eq. 3 of Paper I. This reference optical depth serves as a proxy for geometrical depth. A set of spectral templates (so-called masks) was constructed so as to identify spectral lines forming inside a given range of reference optical depths (or equivalently, geometrical depths) in the atmosphere, as listed in Table A.1. Cross-correlating these masks with a stellar spectrum provides information on the velocity field, the average shape, and strength of lines for a series of atmospheric slices. The capability of the tomographic method to correctly recover the line-of-sight velocity at different depths in a stellar atmosphere was validated in Paper I using state-of-the-art 3D RHD simulations of convection in an RSG atmosphere.
The tomographic method was first applied to the spectroscopic observations of Mira stars in order to interpret their atmospheric dynamics (Alvarez et al. 2001). It was shown that, around maximum light, the line profiles of Mira stars appear double-peaked and follow the Schwarzschild scenario, reflecting the upward propagation of a spherically-symmetric shock front across the photosphere. The tomographic method was then applied by Josselin & Plez (2007) to time-series spectroscopic observations of RSG stars, which revealed complex time-variable upward and downward motions in their atmospheres. Paper II applied an improved version of the tomographic method to the RSG star µ Cep. This study revealed that the spatial and temporal evolution of µ Cep CCFs does not follow the Schwarzschild scenario as in Mira variables due to smaller-amplitude velocity fields, the absence of spherically-symmetric shock waves, and more extended atmospheres.
Finally, Kravchenko et al. (2020) applied the tomographic technique to the high-resolution spectro-interferometric VLTI/AMBER observations of the Mira-type star S Ori and validated the capability of tomography to probe different geometrical depths in the stellar atmosphere. Moreover, this study derived, for the first time, a quantitative relation between opticaland geometrical-depth scales for S Ori.
The various applications of the tomographic method reported above use different sets of masks, constructed from 1D static model atmospheres matching the atmospheric parameters of the stars under consideration. The set of masks used here (Table A.1) is the same as that built in Paper II for RSGs, based on a 1D MARCS model with T eff = 3400 K and log g = −0.4, which falls in the parameter range of RSGs (Levesque et al. 2005). The model was extended to log τ 0 = −4.6 with respect to the standard version of MARCS models (limited to log τ 0 = −3.5) in order to properly handle the lines forming in the outermost layers (mask C5). For the sake of testing the sensitivity of the masks to the adopted model, another set of masks was built from a model with T eff = 3600 K and log g = 0. The vast majority of lines (> 88% per mask) are assigned to the same mask in both models. Therefore, the small mismatch between the atmospheric parameters of Betelgeuse (Table B.1) and those of the model that were used to build the set of masks from Paper II (Table A.1) does not jeopardize the conclusions reached in the present work.  −4.6 < log τ 0 ≤ −4.0 378 * τ 0 is the reference optical depth computed with the continuum opacity at λ = 5000 Å and it serves as a proxy for geometrical depth.    The figure clearly shows that the TiO 7125 Å band is the most sensitive to T eff , whereas the TiO 5847 Å band is the least sensitive.