11institutetext:  IAASARS, National Observatory of Athens, 15326 Penteli,
Greece 22institutetext: National and Kapodistrian University of Athens,
Department of Physics, Panepistimiopolis, Zografos, GR15784, Greece
33institutetext: Institute of Astronomy, KU Leuven, Celestijnlaan 200D, 3001
Leuven, Belgium 44institutetext: Key Laboratory of Space Astronomy and
Technology, National Astronomical Observatories, Chinese Academy of Sciences,
Beijing 100101, People’s Republic of China 55institutetext: Institute of
Astrophysics FORTH, GR-71110, Heraklion, Greece 66institutetext: Carnegie
Observatories, Las Campanas Observatory, Colina El Pino, Casilla 601, La
Serena, Chile 77institutetext: University of Liège, Allée du 6 Août 19c (B5C),
B-4000 Sart Tilman, Liège, Belgium 88institutetext: Department of Astronomy,
University of Geneva, Chemin Pegasi 51, CH-1290 Versoix, Switzerland

# Properties of luminous red supergiant stars
in the Magellanic Clouds

S. de Wit 1122 A.Z. Bonanos 11 F. Tramper 1133 M. Yang 1144 G. Maravelias 1155
K. Boutsia 66 N. Britavskiy 77 E. Zapartas 1188

###### Abstract

Context. There is evidence that some red supergiants (RSGs) experience short
lived phases of extreme mass loss, producing copious amounts of dust. These
episodic outburst phases help to strip the hydrogen envelope of evolved
massive stars, drastically affecting their evolution. However, to date, the
observational data of episodic mass loss is limited.

Aims. This paper aims to derive surface properties of a spectroscopic sample
of fourteen dusty sources in the Magellanic Clouds using the Baade telescope.
These properties may be used for future spectral energy distribution fitting
studies to measure the mass loss rates from present circumstellar dust
expelled from the star through outbursts.

Methods. We apply marcs models to obtain the effective temperature ($T_{\rm
eff}$) and extinction ($A_{V}$) from the optical TiO bands. We use a
$\chi^{2}$ routine to determine the best fit model to the obtained spectra. We
compute the $T_{\rm eff}$ using empirical photometric relations and compare
this to our modelled $T_{\rm eff}$.

Results. We have identified a new yellow supergiant and spectroscopically
confirmed eight new RSGs and one bright giant in the Magellanic Clouds.
Additionally, we observed a supergiant B[e] star and found that the spectral
type has changed compared to previous classifications, confirming that the
spectral type is variable over decades. For the RSGs, we obtained the surface
and global properties, as well as the extinction $A_{V}$.

Conclusions. Our method has picked up eight new, luminous RSGs. Despite
selecting dusty RSGs, we find values for $A_{V}$ that are not as high as
expected given the circumstellar extinction of these evolved stars. The most
remarkable object from the sample, LMC3, is an extremely massive and luminous
evolved massive star and may be grouped amongst the largest and most luminous
RSGs known in the Large Magellanic Cloud (log(L∗/L${}_{\odot})\sim$5.5 and
$R=1400\,\ \textrm{R}_{\odot}$).

###### Key Words.:

stars: massive - stars: supergiants - stars: fundamental parameters - stars:
atmospheres - stars: late-type - Magellanic Clouds

##  1 Introduction

The red supergiant (RSG) phase is the last evolutionary phase for the majority
of massive stars (Levesque 2017). RSGs are subject to intense mass loss
through stellar winds and outbursts (i.e. short lived phases of extreme mass
loss, see Smith 2014, for a review). These outbursts remain poorly
constrained, despite increasing evidence supporting their existence (Prieto et
al. 2008; Mauerhan et al. 2013; Bruch et al. 2021; Humphreys & Jones 2022).
The ability to potentially strip massive stars of their hydrogen envelope may
drastically affect the course of evolution of massive stars in the latest
stages before supernova collapse. At present, stellar evolutionary models for
single stars (e.g. Eggenberger et al. 2008, 2021; Brott et al. 2011; Ekström
et al. 2012; Köhler et al. 2015) or for binary stars (e.g. BPASS; Stanway &
Eldridge 2018) adopt a continuous mass-loss rate, strongly influencing the
type of supernova event. One approach to understand episodic mass loss in
evolved massive stars, is to study the stellar properties of these stars
empirically and test whether the derived properties can be recovered by
existing evolutionary predictions. Several recent works have studied the mass
loss properties (e.g. Beasor et al. 2020, 2021) and stellar properties (e.g.
Massey et al. 2021; González-Torà et al. 2021) of stars in the RSG phase.
Interestingly, Beasor et al. (2020, 2021) have found the mass-loss rates in
the RSG phase to be lower than that of classical recipes (i.e. de Jager et al.
1988; van Loon et al. 2005), predicting that most single stars fail to remove
their hydrogen rich envelope and therefore should explode as a type IIP
supernova. This is inconsistent with the lack of observed high mass type IIP
progenitor stars, known as the ”red supergiant problem” (Smartt et al. 2009;
Smartt 2015; Davies & Beasor 2020). The solution to this contradiction may
either reside in poorly understood short-lived phases ($\sim\\!\\!10^{4}$ yrs;
Beasor et al. 2020) of extreme mass loss, in addition to outbursts in the last
years before the supernova explosion, resulting in a type IIn supernova (Smith
et al. 2009; Smith 2014) or in stripping of the envelope due to binary
interactions (Podsiadlowski et al. 1992; Eldridge et al. 2013; Zapartas et al.
2017).
The mass-loss rates of RSGs can be measured using SED fitting techniques (e.g.
DUSTY; Ivezic & Elitzur 1997). RSGs can reveal an excess emission due to
circumstellar dust as a result of slow and thick stellar winds in the RSG
phase. When an object is surrounded by circumstellar dust, one can derive the
contribution from the central component and the dust component by fitting the
shape of the full SED (Goldman et al. 2017, Yang et al. 2022, in prep.).
Reliable mid-infrared (IR) photometry to model the dust emission bumps (i.e.
10$\mu$m $\&$ 18$\mu$m silicate emission and 11.3$\mu$m PAH emission; Jones et
al. 2017) is vital, but such data are scarce due to instrument limitations,
especially for more distant galaxies. Constraining the properties of the
central source through spectroscopy eliminates potential degeneracies and
uncertainties in the SED fitting procedure, greatly reducing the number of
viable models for fitting the SED. To assess the mass loss, we therefore first
need to derive accurate stellar parameters, such as the effective temperature
of the central source. We approach this by targeting RSGs in the Large
Magellanic Cloud (LMC) and Small Magellanic Cloud (SMC), using selection
criteria from the near-IR and mid-IR (e.g. Bonanos et al. 2009, 2010;
Britavskiy et al. 2014), to look for potential dust-obscured RSGs. In recent
decades, various methods have been applied to derive the properties of RSGs.
Using the marcs model atmospheres (Gustafsson et al. 2008), Davies et al.
(2013) and González-Torà et al. (2021) have fitted the shape of the SED, while
Davies et al. (2015) and Tabernero et al. (2018) have used these models to fit
spectral lines in the $J$ and $I-$band, respectively. In our work however, as
we had access to optical spectra, we used a similar approach to Levesque et
al. (2005, 2006), where we measure the effective temperature from the depths
of the optical TiO bands.
In Section 2 we discuss the selection of targets and describe the
spectroscopic observations, data reduction and spectral classification of the
spectra. We present the results of spectral modelling in Section 3 and compare
them to empirical photometric relations and stellar evolutionary models. In
Section 4, we compare the methodology and results to other, similar studies,
and present the conclusions in Section 5.

##  2 Target selection and observations

###  2.1 Target selection

We based our target selection on foreground-cleaned, multi-wavelength
photometric catalogues for the SMC (Yang et al. 2019) and LMC (Yang et al.
2021). These extensive catalogues comprise many photometric data sets,
including Gaia (Gaia Collaboration et al. 2016, 2018), Spitzer (Werner et al.
2004) and 2MASS (Cutri et al. 2003; Skrutskie et al. 2006) photometry for
numerous evolved massive stars in the Magellanic Clouds. To properly select
dusty supergiant candidates, we used a set of criteria to distinguish them
from the asymptotic giant branch (AGB) population. Following a similar
approach to Britavskiy et al. (2014), we used the following set of criteria:

i) M[3.6] ¡ $-$8 mag (see Britavskiy et al. 2014);
ii) J$-$[$3.6$] ¿ 1 mag (see Bonanos et al. 2009, 2010);
iii) detection in [24]

The above criteria allow us to select the most luminous (i.e. massive and
evolved) and red sources in the near-IR and mid-IR bands. RSGs are bright
sources in the near-IR due to their dramatically increased radii and low
effective temperatures, while disks around sgB[e] stars are regions in which
circumstellar dust may form, making them appear bright and red in the near and
mid-IR. Other potential contaminants may be Extreme-AGBs. These overlap with
the region of the mid-IR CMD in which yellow supergiant stars (YSGs) and
sgB[e] stars are found.
A detection in [$24$] indicates a cooler dusty circumstellar environment,
revealing a potential preceding phase of enhanced mass loss. In Table 1 we
list the coordinates and photometry for all selected targets. The first four
columns indicate the name and coordinates of the selected targets. The
following nine columns show the magnitudes used to construct the CMDs and to
derive effective temperatures. Gaia magnitudes were taken from Gaia DR2,
$V_{\textrm{mean}}$ from ASAS-SN, J $\&$ Ks from 2MASS, [3.6] $\&$ [4.5] from
IRAC and [24] from MIPS. When available, we added time series data from the
All-Sky Automated Survey for Supernovae (ASAS-SN; Shappee et al. 2014;
Jayasinghe et al. 2018), to highlight the type and amplitude of the
variability for our sources. The type (i.e. semi-regular; SR, Mira variable; M
or slow-irregular; L) and amplitude of the variability are shown in the last
two columns of Table 1. We obtained the amplitude of the variability from the
difference between the two most extreme points in the ASAS-SN light curves.
The variability information improved our understanding of the sources as well
as the interpretation of the uncertainties on the derived parameters.
Fig. 1 shows colour-magnitude diagrams (CMDs) in the near-IR to mid-IR for the
SMC and LMC. We plot spectroscopically verified RSG populations from Levesque
et al. (2006) and Davies et al. (2013) for comparison. These two works have
used similar methods to derive the properties of RSGs and therefore provided
the best possible comparison to our results. We then illustrate the criteria
used to select evolved massive stars with the blue shaded box. The coloured
markers in Fig. 1 make up the final sample after applying these criteria. From
the top right CMD in Fig. 1 it is apparent that we have selected redder
sources compared to other studies, implying that these sources may be dust-
obscured RSGs with high extinction. Furthermore, we verified that most of our
sources were clustered at the tip of the RSG branch of the optical CMD. In
this CMD, AGB stars are expected to be found extending towards the right,
while other contaminants are expected to be at the fainter end of the diagram
(see Fig. 2 bottom panel). Next, we searched the literature for spectral
classifications for sources fitting the selection criteria and discarded those
with accurate existing spectral classifications from the target list. The
final sample of observed stars consisted of ten previously understudied
objects, from which eight objects are located in the LMC and two in the SMC.
Most of these objects had been assigned a spectral type in the past, but no
luminosity class, implying they have not been previously (spectroscopically)
verified as RSGs. One object was added to the sample (SMC1) due to its rare
class (sgB[e]; e.g. Kraus 2019), to investigate evidence of variations in its
spectral type, completing the sample to eleven stars. We also retained three
sources that partially satisfied our selection criteria, to investigate
potential contaminants; i.e. altogether we observed fourteen sources.

Figure 1: Color magnitude diagrams for the SMC (left) and LMC (right). The
background points are sources from the catalogues used for target selection
(Yang et al. 2019, 2021). The colour bar indicates the number density of the
stellar population. We show RSGs classified by Levesque et al. (2006) in
black. Inverted grey triangles indicate contaminating objects that were
observed (i.e. SMC4, SMC5 and LMC9). Our final sample is indicated with the
coloured symbols indicated by the legend. Top: M[3.6] vs. $J-[3.6]$ CMD
showing the criteria used to select supergiant stars. In the panel on the
right, seemingly only 7 inverted red triangles are plotted, which is due to
the nearly identical position of LMC1 and LMC5 in this panel. Bottom: M[3.6]
vs. ${\rm[3.6]-[4.5]}$ CMD to visually inspect the locations of the objects
using IRAC bands. The numbers correspond to the object IDs of the verified
supergiants.

Figure 2: Similar to Fig. 1, but for [24] vs. $\textrm{G}-[24]$ (top) and G
vs. $\textrm{G}_{\textrm{BP}}-\textrm{G}_{\textrm{RP}}$ (bottom). The top CMD
was used for the third selection criterion, while the bottom CMD was used to
check whether the selected sources cluster in regions where they are expected.
The RSG branch extends upwards in the middle, with the AGB branch extending to
the right. Table 1: Photometric data of the sample.

ID | SIMBAD Name | R.A. | DEC | $G$ | $G_{\textrm{BP}}$ | $G_{\textrm{RP}}$ | $V_{\textrm{mean}}$ | $J$ | $K_{S}$ | $[3.6]$ | $[4.5]$ | $[24]$ | Var. | $V_{\textrm{amp}}$
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---
|  | (J2000) | (J2000) | (mag) | (mag) | (mag) | (mag) | (mag) | (mag) | (mag) | (mag) | (mag) | Type | (mag)
SMC1 | LHA 115-S 23 | 00 55 53.81 | -72 08 59.51 | 13.32 | 13.39 | 13.16 | 13.36 | 13.01 | 12.18 | 10.63 | 9.87 | 7.34 | - | -
SMC2 |  $[$MA93$]$ 1810 | 01 15 21.22 | -73 30 15.10 | 12.82 | 12.87 | 12.71 | 12.80 | 12.67 | 11.69 | 9.94 | 9.45 | 6.92 | - | -
SMC3 | PMMR 64 | 00 53 53.08 | -71 42 44.37 | 12.10 | 13.29 | 11.06 | 12.89 | 9.67 | 8.65 | 8.46 | 8.33 | 6.06 | L | 1.50
SMC4 | IRAS 00350-7436 | 00 36 59.56 | -74 19 50.26 | 14.17 | 14.92 | 13.25 | 14.56 | 11.32 | 9.13 | 7.67 | 7.12 | 3.58 | - | 0.19
SMC5 | SV* HV 859 | 01 10 26.92 | -72 35 48.57 | 14.18 | 16.52 | 12.66 | 15.26 | 10.28 | 9.08 | 8.72 | 8.58 | 7.86 | M | 3.80
LMC1 | WOH S 57 | 04 53 14.80 | -69 12 18.00 | 12.34 | 13.98 | 11.11 | 13.59 | 9.14 | 7.75 | 7.77 | 7.62 | 4.99 | SR | 1.14
LMC2 | WOH S 374 | 05 31 47.42 | -66 03 40.58 | 13.01 | 14.78 | 11.78 |  14.401 | 9.80 | 8.67 | 7.84 | 7.39 | 3.27 |  Var:1 | -
LMC3 |  $[$W60$]$ B90 | 05 24 19.31 | -69 38 49.37 | 12.29 | 14.47 | 10.95 | 15.11 | 8.36 | 6.81 | 6.29 | 6.29 | 3.35 | L | 0.38
LMC4 | UCAC2 2674864 | 05 44 13.77 | -66 16 44.62 | 12.86 | 14.99 | 11.54 | 14.28 | 8.96 | 7.49 | 7.01 | 7.46 | 4.71 |  SR2 |  0.313
LMC5 | SP77 28-2 | 04 57 26.36 | -66 23 25.78 | 11.91 | 12.99 | 10.61 | 13.03 | 9.09 | 7.88 | 7.72 | 7.73 | 4.92 | SR | 1.54
LMC6 | SV* HV 12185 | 05 09 43.58 | -65 21 59.18 | 11.74 | 13.30 | 10.59 | 12.97 | 8.83 | 7.69 | 7.41 | 7.36 | 4.12 | SR | 1.49
LMC7 | SV* HV 12793 | 05 23 43.61 | -65 41 59.84 | 12.23 | 13.98 | 11.01 | 13.57 | 9.01 | 7.74 | 7.39 | 7.53 | 3.94 | SR | 1.11
LMC8 | W61 19-24 | 05 32 20.19 | -67 32 42.36 | 12.00 | 13.45 | 10.85 | 12.30 | 9.22 | 8.01 | 7.60 | 7.72 | 4.57 | SR | 0.71
LMC9 | MSX LMC 806 | 05 32 03.42 | -67 42 25.70 | 13.47 | 13.91 | 12.75 | 13.69 | 11.82 | 10.04 | 8.15 | 7.54 | 0.37 | - | -

  * •

Notes: The variability type indicates the following classes: semi regular
(SR), mira variable (M) and slow irregular (L) variables.
1Magnitude and variability designation as listed by Samus’ et al. (2017).
2Designation from Gaia DR2.
3Semi-amplitude from Groenewegen & Sloan (2018).

###  2.2 Spectroscopic observations

We obtained optical spectra for the fourteen selected targets with the 6.5m
Baade telescope at Las Campanas Observatory, Chile. The targets were observed
with the Magellan Echellette (MagE) spectrograph (Marshall et al. 2008) using
the 1” science slit. MagE is a medium resolution (resolving power $\rm{R}\sim
4100$ for the 1” slit) single-object spectrograph and provides spectral
coverage in the entire optical domain ($\lambda$ $3200-10000\,\AA$). Multiple
exposures were taken for each object. The observations were completed in four
separate observing runs between December 2018 and March 2020. A log of the
observations is presented in the first three columns of Table 2.

Table 2: Observational properties.

ID | Obs. Date | Exp. Time | S/N | RV | New Spectral | Previous Spectral | References
---|---|---|---|---|---|---|---
| (UT) | (s) |  |  (km s-1) | Type | Type |
SMC1 | 19-Dec-2018 | 3x200 | 80 |  155 $\pm$ 4 | A0 I[e] | B8 I[e], A1 I[e] |  Zickgraf et al. (1992), Kraus et al. (2008)
SMC2 | 19-Dec-2018 | 3x180 | 85 |  189 $\pm$ 10 | F8 I | … | …
SMC3 | 13-Nov-2019 | 3x100 | 95 |  178 $\pm$ 7 | M2 I | K/M | Sanduleak (1989)
SMC4 | 12-Sep-2019 | 3x600 | 65 | … | C0-2,2e | C3,2e | Groenewegen & Blommaert (1998)
SMC5 | 13-Nov-2019 | 3x300 | 20 | … | M7e III | K5/M7e | Wood et al. (1983)
LMC1 | 12-Nov-2019 | 3x200 | 60 |  272 $\pm$ 8 | M2 I | M0: | Westerlund et al. (1981)
LMC2 | 13-Nov-2019 | 3x150 | 45 |  306 $\pm$ 11 | M4-5 II-III | M6 | Reid et al. (1990)
LMC3 | 08-Mar-2020 | 3x180,1x240 | 55 |  260 $\pm$ 6 | M3 I | M2 | Chen et al. (2016)
LMC4 | 08-Mar-2020 | 3x180 | 50 |  296 $\pm$ 6 | M0 I | M0.5 | Chen et al. (2016)
LMC5 | 08-Mar-2020 | 3x100 | 55 |  290 $\pm$ 5 | M1 I | M0 | Westerlund et al. (1981)
LMC6 | 08-Mar-2020 | 3x100 | 75 |  318 $\pm$ 5 | M1 I | M1 I: | Samus’ et al. (2017)
LMC7 | 08-Mar-2020 | 3x180 | 70 |  316 $\pm$ 9 | M2 I | M3/M4 | Reid et al. (1990)
LMC8 | 08-Mar-2020 | 3x120,1x100 | 65 |  304 $\pm$ 5 | M3 I | M2 | Egan et al. (2001)
LMC9 | 08-Mar-2020 | 3x180 | 70 | … |  H ii \+ LBVc | YSO: | Seale et al. (2009)

  * •

Notes: The signal-to-noise ratio (S/N) of the combined spectrum is measured in
the red ($\sim$8500$\AA$). No previous luminosity classes were indicated for
most stars. The previous classification of LMC6 included a luminosity class,
but this was an uncertain classification.

###  2.3 Data reduction

The spectroscopic data were reduced using the MagE Spectral Extractor (MASE;
Bochanski et al. 2009). The MASE pipeline offers the full reduction process
from raw data to a 1D spectrum, including the bias subtraction, flat fielding,
sky subtraction, cosmic ray removal, flux and wavelength calibration and
finally the extraction of the echelle orders to 1D-spectra for all exposures.
After carefully transforming the wavelength units from vacuum to air
wavelengths, which is required for spectral modelling, we merged the
individual exposures using the median flux in every wavelength bin. The 1D
echelle orders were then connected using linear weighting at overlapping
wavelengths to assemble the full 1D spectrum.
Artificial bumps were present in the merged spectrum, caused by the small
offset between orders. We have corrected for this by scaling down the
individual orders so that their relative fluxes matched at overlapping
wavelength sections. Once the spectrum was cleaned for these bumps, we scaled
the flux back up in units of absolute flux. The spectra were not corrected for
tellurics. The S/N presented in the fourth column of Table 2 has been derived
by calculating the median S/N in small chunks of pseudo-continuum in the
I-band and serves as a general indicator of the spectral quality (S/N $\geq
50$ in most cases).

###  2.4 RV measurements

Before proceeding with the spectral modelling, the spectra were corrected for
the radial velocity (RV) shift. Specifically, the Ca ii triplet
($\lambda$8498, $\lambda$8542, $\lambda$8662), the Na i doublet
($\lambda$5890, $\lambda$5896) and two resolved single metal lines close to
the Ca ii triplet (Ti i $\lambda$8426 and Fe i $\lambda$8612) were used to
estimate the RV111See table C.1 in Dorda et al. (2016) for a complete atlas of
spectral lines for cool supergiants in the $i$ and $z$-bands.. We measured the
central wavelength for the selected metal lines by fitting a Gaussian to the
line profile. Then, we calculated the offset from the rest-wavelength for
these lines to determine a velocity shift for each line profile. Precisely
determined values for the rest-wavelength of spectral lines were taken from
the NIST atomic database (Kramida et al. 2020). The final RV was determined as
the weighted mean of the shifts of individual lines. We then applied this RV
shift to the spectrum. Towards mid-late M-types, prominent metal lines
disappear due to the increasing strength of TiO $\lambda$8432-8442-8452. In
our sample this only applied to LMC2, for which we have derived the velocity
shift solely from the Ca ii triplet. The mean radial velocity for each target
is presented in the sixth column of Table 2. We note that the radial
velocities we have found are consistent with the systemic velocity of the LMC
and SMC (260 $\rm km\,s^{-1}$ and 150 $\rm km\,s^{-1}$, respectively,
McConnachie 2012; van der Marel & Kallivayalil 2014) confirming their
membership to these systems. When a RV value from Gaia DR2 was available, we
compared our results and confirmed them to be consistent. Additionally, the
improved Gaia EDR3 parallax (Gaia Collaboration et al. 2021) was inspected. We
found the parallaxes of all objects, considering their errors and the
systematic correction from Lindegren et al. (2021) to be consistent with zero,
which is in turn consistent with the extragalactic origin of the sources.

###  2.5 Spectral classification

####  2.5.1 Classification of RSGs

The abundant presence of molecular bands in the optical are a clear indicator
of late-type stars. Depending on the metal content of the star, the relative
strengths of the molecular bands are typically used to determine the spectral
sub-type of M-stars. Despite being derived from strong and abundant molecular
bands, the optical classification of M0-M3 stars posed difficulties, as the
relative strengths of the TiO bands change only slightly in this range, while
metal lines at $\sim$8400-8800$\AA$ remain unchanged due to their
insensitivity to temperature (Solf 1978). For later types (M4-M6), Solf (1978)
indicated the blending of the Ti ii $\lambda$8435 doublet with the triple TiO
band-head at $\lambda$8432-8442-8452 as the primary classifier.
Given these spectral classification criteria, we then proceeded with the
classification of our spectra by visual inspection. Nine stars were identified
as M-type stars due to the presence of strong molecular TiO bands at
$\lambda$6150 and $\lambda$7050\. Fig. 3 presents the optical spectra of these
stars. The most important spectral features are indicated and a spectral type
was added. For seven of these stars, we identified similar TiO band strengths
amongst them. For this group, we sorted the spectra from weakest TiO
$\lambda$7050 band (M1) to strongest (M3). Intermediate objects were
classified as M2. The remaining two stars (LMC2 and LMC4) were noticeably
different from the rest of the sample. Fig. 4 highlights spectral lines in the
Ca ii triplet region, for example Ti ii $\lambda$8435\. For LMC2, the TiO
$\lambda$8432 band is present and blends with Ti ii $\lambda$8435, but does
not dominate the metal lines in this region. When the spectrum is depleted of
metal line features due to the increasing TiO $\lambda$8432 band strength, the
star is considered as spectral type M6 (Negueruela et al. 2012). Therefore, we
classified LMC2 as M4-5. For LMC4, the TiO bands at $\lambda$6150 and
$\lambda$7050 were present, but were noticeably weaker compared to the other
stars, indicating an early M type. We classified this object as M0.
For the luminosity class we used the strength of the Ca ii triplet
($\lambda$8498, $\lambda$8542 $\&$ $\lambda$8662) as the main criterion. The
Ca ii triplet is sensitive to changes in surface gravity due to pressure
broadening, but insensitive to $T_{\rm eff}$ and $Z$ changes (Massey 1998).
Therefore, a strong Ca ii triplet is expected for supergiant stars. A strong
Ca ii triplet was indeed detected in all stars, with the exception of LMC2
(see Fig. 4). Apart from the Ca ii triplet, the ratio between Fe i
$\lambda$8514 and Ti i $\lambda$8518 was indicative (except LMC2, where
neither line clearly dominates the other) of luminosity class I (Britavskiy et
al. 2014). LMC2 was classified as M4-5 II-III, as it does not satisfy the
criteria to be classified as a luminosity class I star. To our knowledge,
criteria for luminosity class II are not well established. Therefore, we
cannot reject this possibility and present a range of luminosity class II-III
for this object.
The variability types presented in Table 1 support the spectral
classification, as stars of type semi-regular (SR) and slow-irregular (L) are
often attributed to the RSG or red-giant class. Furthermore, most of our
classified RSGs show minimum to modest (up to 1.5 magnitudes) variability in
the optical, which is expected for RSGs. The final classification of our
objects, along with a comparison to older classifications and references for
these, is presented in the final columns of Table 2.

Figure 3: Spectra from M-type sources in our sample (an arbitrary offset has
been applied for illustration purposes). The Ca ii triplet and key TiO bands
are indicated.

Figure 4: A zoom-in on the Ca ii triplet and prominent metal lines (see Dorda
et al. 2016) used for the classification of the objects and RV analysis.
Features marked with a * are blends and only the most significant line is
indicated. For LMC2, Ti i $\lambda$8435 is relatively strong and the ratio
between Fe i $\lambda$8514 and Ti i $\lambda$8518 line strengths is close to
unity.

####  2.5.2 Classification of other supergiants

During the observing campaign, two additional evolved massive stars were
observed. We present spectra for these in Fig. 6, 10, 11 and 12. The nature of
their spectra was such that they could be classified, but needed higher
spectral resolution for adequate modelling (i.e. the models constrain
properties from narrow metal line ratios), and this was not pursued. We
briefly describe their spectral types here.
The spectrum of SMC1 displayed strong double peaked emission for the Balmer
series and several forbidden [O i] and [Fe ii] emission lines (see Fig. 10).
These features originate from a circumstellar disk or ring structure
(Maravelias et al. 2018). Despite being previously classified as a B8 I[e]
(Zickgraf et al. 1992, spectrum taken in December 1989) and as A1 I[e] (Kraus
et al. 2008, spectrum taken in October 2000), we determined a spectral type of
A0 I[e] (spectrum taken in December 2018) from the Mg ii $\lambda$4481 to He i
$\lambda$4471 line ratio (see Fig. 11). If Mg ii $\lambda 4481\geq$ He i
$\lambda$4471, the spectral type is later than B9, while when the ratio of Ca
K/(H$\epsilon+$Ca H) $\leq 0.33$, the spectral type is A0 (Evans et al. 2004),
which is the case for our spectrum. We also note that the He i lines
reappeared in our spectrum, after being absent in the spectrum analysed by
Kraus et al. (2008), making the spectrum appear slightly hotter than it was in
Kraus et al. (2008). Assuming the spectra have been classified correctly in
past studies, the spectral type of SMC1 apparently varied on the scale of
decades. For the radial velocity measurement, several metal lines in the
optical wavelength range were used, such as Mg ii $\lambda$4481, He i
$\lambda$4471, 5875, Si ii $\lambda$6347 and Fe ii $\lambda$4174, 4179.
Modelling of this source was not pursued due to the complexity of modelling
the disk emission features and their contaminating effect on key Balmer and He
i lines.
SMC2 is a yellow supergiant star that displayed a strong H$\alpha$ emission
component, indicating a mass-loss component surrounding the star. The spectrum
is characterized by a series of Ti and Fe lines in the
$4000-5000\,\textrm{\AA}$ region and strong Ca ii H and K absorption (see Fig.
12). We classified this star as F8 I, due to the absence of the G-band (¡ G0)
and the strong Ti-Fe forest lines relative to the Balmer lines at
$Z=0.2\,Z_{\odot}$ (¿ F5). The sharpness of the Balmer lines and the strengths
of the metal lines suggests that this star is of luminosity class I. For the
radial velocity, several resolved metal lines were used, such as Mg ii
$\lambda$4481, Si ii $\lambda$6347 and Na i $\lambda$5890, 5896. For yellow
supergiants, spectral models are scarce. Models for a $T_{\rm eff}$ estimate
based on Fe line ratios have recently been published (Kourniotis et al. 2022),
however they require us to resolve narrow Fe lines. Our current spectrum is
not able to resolve these lines and therefore a higher resolution spectrum is
needed to be able to derive a $T_{\rm eff}$ empirically.

Figure 5: Top: Spectrum of SMC1, an A0 I[e] star. The spectrum is
characterized by several strong emission lines, which are indicated. Bottom:
Spectrum of SMC2, a F8I star. Strong H$\alpha$ emission was detected
indicating expanding circumstellar material.

####  2.5.3 Classification of other sources

Three more objects fitting the selection criteria were observed (see Fig. 1
and 2, inverted grey triangles), but were not classified as evolved massive
stars. We found SMC4 and SMC5 to have characteristic features of giant stars,
while LMC9 shows the characteristics of an ionized nebula.
The spectrum of SMC4 contains carbon absorption bands at $\lambda$4737 and
$\lambda$5165 and deep P-Cygni profiles indicative of strong mass loss
(Whitelock et al. 1989), which could explain its extreme magnitude at [24].
SMC4 is a spectroscopic binary, with pairs of shifted spectral lines. We
classified the primary star as a carbon star with moderately strong Swan C2
bands at $\lambda$4737 and $\lambda$5165, but with the absence of the C2
$\lambda$5585 band. Due to the strength of the Balmer lines, the absence of Ca
i $\lambda$4426 and the presence of the G-band, albeit weak compared to
H$\gamma$, we conclude this star is of G-type, and its carbon star equivalent
should therefore be in the range C0-2. We were not able to obtain a more
secure classification due to the contaminating features of the companion star.
SMC5 is characterized by extreme TiO absorption bands. The relative strengths
of the TiO bands are indicative of a late M-star, beyond M6. The strength of
VO $\lambda$7900 suggests a spectral type of M7. The strength of the Ca ii
triplet indicates the star is of luminosity class III.
The spectrum of LMC9 was characterized by strong, narrow emission lines, most
notably the [O iii] feature at $\lambda$5007 and a saturated H$\alpha$
emission component. Most of the emission lines were identified as
recombination lines from low ionisation states, typically present in a Hii
region. A stellar continuum has been observed underneath the sharp emission
lines, revealing a bright blue object. From the broadened Balmer line profiles
and Fe ii emission lines present throughout the spectrum we conclude that a
Luminous Blue Variable candidate (LBVc) may be present inside this Hii region,
although further investigation is needed to verify this object as such.

Figure 6: Top: Spectrum of SMC4, a carbon star. The spectrum is characterized
by two molecular carbon bands and several strong, deep P-Cygni profiles.
Middle: Spectrum of SMC5, a late M giant star. Strong TiO bands are revealed
in the red part of the spectrum, indicating a low effective temperature. From
the Ca iitriplet, only $\lambda$8662 is present, albeit almost completely
blended in with overarching TiO band. Bottom: Spectrum of LMC9, an ionized
nebula. The overall spectrum is dominated by strong emission lines coming from
a nebula. Strong [O iii] was detected and the H$\alpha$ emission was
saturated. The emission line profile of H$\alpha$ shows significant broadening
and other circumstellar lines (Fe ii) were detected in the left wing of
H$\alpha$, revealing a LBVc.

##  3 Spectral modelling and resulting parameters

###  3.1 Fitting the marcs models

We used a grid of alpha-poor, spherical marcs model atmospheres (Gustafsson et
al. 2008) to fit the M-type stars in our sample. The computed models have a
mass of $15\,\textrm{M}_{\odot}$. Given the expected mass range for stars that
evolve towards the RSG phase (8M⊙ $\leq$ M∗ $\leq$ 25M⊙) and the availability
of models with limited discrete stellar masses
($0.5,1,2,5\,\&\,15\,\textrm{M}_{\odot}$), the 15$\,\textrm{M}_{\odot}$ model
was the most suitable for our studies. Furthermore, we assume that a single
mass to represent the entire range of masses for RSGs is justified given that
the geometrical thickness, and thus the atmospheric structure is largely
unaffected in this mass range (Davies et al. 2010). The microturbulent
velocity of the models was fixed to $\xi=5\,\textrm{km s}^{-1}$. González-Torà
et al. (2021) have shown that changes in $\xi$ had little effect on the final
result. We used marcs models with log$(Z/Z_{\odot})=-0.35\,\textrm{dex}$ for
the LMC and log$(Z/Z_{\odot})=-0.55\,\textrm{dex}$ for the SMC, which
adequately represent the average metallicity derived empirically from a set of
RSGs in the Magellanic Clouds (Davies et al. 2015). The surface gravities vary
between -0.5 and +$0.5\,\textrm{dex}$ in steps of $0.1\,\textrm{dex}$. Only
models with a surface gravity of +0.5, +0.0 and -$0.5\,\textrm{dex}$ were
available from the marcs platform, hence we interpolated the flux of adjacent
models linearly to make the grid more dense. The effective temperature ranges
from $3300-4500\,\textrm{K}$, with $10\,\textrm{K}$ steps for the range of
$3300-4000\,\textrm{K}$ and $25\,\textrm{K}$ steps for the range of
$4000-4500\,\textrm{K}$ using similar interpolation strategies.
We first degraded the resolution of the models from $R=20000$ to $R=4000$, so
that the model spectrum matches the spectral resolution of MagE. We then
applied Fitzpatrick (1999) reddening laws with varying $A_{V}$ to the models
to derive the best fit extinction factor, assuming a typical total-to-
selective extinction of $R_{V}=2.74$ and $R_{V}=3.41$ for the SMC and LMC,
respectively (Gordon et al. 2003). We note that the values for $R_{V}$ may be
higher, considering the grain size distribution in the circumstellar
environment of RSGs (Massey et al. 2005). A recent study by González-Torà et
al. (2021) however, have indicated that changes in $R_{V}$ and the type of
extinction law used, do not largely contribute to temperature changes of the
best fit model. To obtain the best fit model to our spectra, we computed the
reduced chi-squared ($\chi^{2}_{\rm red}$):

$\displaystyle\chi^{2}_{\textrm{red}}=\frac{1}{\nu}\sum^{\nu}_{i}\left(\frac{F_{i}-F_{i,\textrm{mod}}}{\sigma_{i}}\right)^{2},$
(1)

where $\nu$ is the degrees of freedom, $F_{i}$ and $F_{i,\textrm{mod}}$ are
the fluxes of the spectrum and model, respectively, and $\sigma_{i}$ is the
uncertainty on the flux in bin $i$. The degrees of freedom are set by the
amount of bins minus the amount of free parameters ($T_{\rm eff}$ and
$A_{V}$). The accepted best-fit model was the model on the 2D grid with lowest
$\chi^{2}_{\rm red}$ (i.e. $\chi^{2}_{\rm min}$). As in González-Torà et al.
(2021), we chose to smooth the spectra to determine the $\chi^{2}_{\rm red}$,
given that there may be uncertainties in the molecular transitions in the
marcs models. Bin by bin fitting of the single wavelength bins may therefore
yield unrealistic $\chi^{2}_{i}$, which could heavily impact the
$\chi^{2}_{\textrm{red}}$. For both the model and the observed spectrum, we
grouped smaller wavelength bins into larger bins of $50\,\textrm{\AA}$ and set
the mean flux of the new bin as the corresponding data point. We then
calculated the $\chi^{2}_{i}$ of the larger bin using the mean flux and mean
error (i.e. the standard deviation of the flux measurements). We fitted
selected wavelength bins between $5400\,\textrm{\AA}$ and
$8800\,\textrm{\AA}$, which include the most prominent and temperature
sensitive TiO bands as well as the bluest and reddest $100\,\textrm{\AA}$ to
describe the slope of the spectrum to estimate $A_{V}$. Upon inspecting the
$\chi^{2}_{\rm red}$ values after the initial run, we chose to discard the TiO
band at $6150\,\textrm{\AA}$ from the calculation, as the molecular
transitions in the models do not accurately match the features of real RSG
spectra (see Fig. 7) and therefore skew the results. As our spectra were not
corrected for telluric contamination, we avoided bins including telluric
bands. Isolated telluric lines may still be present (Catanzaro 1997), but
given the use of the aforementioned smoothing technique, this does not
significantly affect the $\chi^{2}_{\rm red}$ calculation.
We approached the modelling as follows: we first fit the Ca ii triplet to
derive the surface gravity (log $g$). The Ca ii triplet is the primary
diagnostic sensitive to the pressure broadening in the optical, so fitting
this narrow spectral domain allowed us to get a precise measurement of the
surface gravity. We proceeded by fixing the derived value of log $g$ for the
remainder of the fitting routine. We then fitted a 2D grid of models with
variable $A_{V}$ and $T_{\rm eff}$ to the spectrum and located the point on
the grid with $\chi^{2}_{\rm min}$. An example fit is presented in Fig. 7, in
which we show the observed spectrum in black and the model in red (see Fig. 14
in the Appendix for the remainder of objects). From here on, we combine the
extracted best fit model parameters with photometry to derive the luminosity
and subsequent radius of the star. For this, we used the $K-$band magnitude
($m_{K}$), an appropriate bolometric correction for the $K-$band based on the
spectral type (Davies et al. 2018) and the distance modulus ($\mu_{\rm LMC}$
$=18.477\pm 0.030$ and $\mu_{\rm SMC}$ $=18.95\pm 0.07$ mag; Pietrzyński et
al. 2019; Graczyk et al. 2014) to derive $M_{\rm bol}$ ($M_{\rm bol}=$
$m_{K}-\mu-A_{K}+$BCK). Using the magnitude-luminosity relation, we then
derived L∗ ($L=L_{0}\times 10^{(M_{\rm bol})/2.5}$, with $L_{0}$ being the
zero point luminosity of the sun $\sim 3\times 10^{28}\,\ \textrm{W}$). The
radius (R∗) was derived from the Stefan-Boltzmann equation ($L=4\pi\sigma
R^{2}T^{4}_{\rm eff}$). The parameters corresponding to the best fit model of
each star (i.e. $\chi^{2}_{\rm red,min}$, $T_{\rm eff,TiO}$, $A_{V}$, $M_{\rm
bol}$, log(L∗/L⊙) and R∗) are presented in Table 3.
To propagate the uncertainties, similar to Davies et al. (2010), we
constructed a 2D $\chi^{2}$ map for all fitted stars (see Fig. 8). The
$\chi^{2}$ map shows the derived $\chi^{2}$ for every point on the 2D grid.
Using a $\chi^{2}_{\textrm{min}}+\Delta\chi^{2}$ (where $\Delta\chi^{2}=2.3$,
6.2 and 11.8 for 1, 2 and 3$\sigma$), we were able to derive uncertainties
from the contours around the best fit solution and take the extreme points as
the inferred uncertainties. Following standard error propagation rules, we
then estimated uncertainties for $M_{\textrm{bol}}$, log(L∗/L⊙) and R∗ and
listed the uncertainties in Table 3.

Figure 7: Best fit marcs model (red dotted line) to the spectrum of LMC1
(black solid line). Regions included in the $\chi^{2}_{\textrm{red,min}}$
calculation are indicated with shaded areas. The best fit properties are
indicated in the top left corner. Tellurics are indicated with a T, the Ca ii
triplet with green shades, the regions used for $A_{V}$ in yellow shades and
the TiO bands are indicated with blue shades. The TiO band at $\lambda$6150
(hatched) was excluded from the fit.

Figure 8: A map of $A_{V}$ and $T_{\rm eff}$ for LMC1, with darker shades
indicating parts of the grid with lower $\chi^{2}$ values. The $\chi^{2}$
values of the color bar are with respect to the minimum $\chi^{2}$. Contours
at 1$\sigma$, 2$\sigma$ and 3$\sigma$ are indicated to estimate the
uncertainties. Table 3: Parameters derived for stars of spectral class M.

ID | $\chi^{2}_{\rm red}$ | $T_{\rm eff,TiO}$ | $A_{V}$ | $\rm log\,\it g$ | $M_{\textrm{bol}}$ |  log(L∗/L⊙) |  R∗
---|---|---|---|---|---|---|---
|  | (K) | (mag) | (dex) | (mag) | (dex) |  (R⊙)
SMC3 | 1.0 |  3570${}^{+21}_{-32}$ |  0.15${}^{+0.04}_{-0.08}$ |  +0.3 $\pm$ 0.1 |  -7.50${}^{+0.20}_{-0.21}$ |  4.90${}^{+0.08}_{-0.08}$ |  730${}^{+75}_{-65}$
LMC1 | 1.0 |  3600${}^{+68}_{-37}$ |  0.40${}^{+0.13}_{-0.12}$ |  +0.4 $\pm$ 0.1 |  -7.92${}^{+0.15}_{-0.15}$ |  5.06${}^{+0.06}_{-0.06}$ |  875${}^{+70}_{-60}$
LMC2 | 5.2 |  3410${}^{+32}_{-52}$ |  0.30${}^{+0.20}_{-0.19}$ | ¿ +1.0 |  -6.90${}^{+0.23}_{-0.22}$ |  4.65${}^{+0.09}_{-0.09}$ |  610${}^{+75}_{-60}$
LMC3 | 2.9 |  3570${}^{+59}_{-46}$ |  1.00${}^{+0.14}_{-0.14}$ |  +0.5 $\pm$ 0.1 |  -8.89${}^{+0.18}_{-0.18}$ |  5.45${}^{+0.07}_{-0.07}$ |  1390${}^{+130}_{-110}$
LMC4 | 1.1 |  3710${}^{+216}_{-87}$ |  0.90${}^{+0.18}_{-0.18}$ |  +0.1 $\pm$ 0.1 |  -8.32${}^{+0.21}_{-0.21}$ |  5.22${}^{+0.08}_{-0.08}$ |  990${}^{+115}_{-100}$
LMC5 | 1.6 |  3620${}^{+53}_{-55}$ |  0.35${}^{+0.10}_{-0.13}$ |  +0.3 $\pm$ 0.1 |  -7.82${}^{+0.13}_{-0.15}$ |  5.03${}^{+0.05}_{-0.06}$ |  825 ${}^{+60}_{-60}$
LMC6 | 1.9 |  3640${}^{+59}_{-51}$ |  0.35${}^{+0.08}_{-0.14}$ | ¿ +0.5 |  -8.01${}^{+0.11}_{-0.16}$ |  5.10${}^{+0.05}_{-0.06}$ |  890 ${}^{+55}_{-65}$
LMC7 | 3.5 |  3600${}^{+25}_{-49}$ |  0.50${}^{+0.05}_{-0.14}$ | ¿ +0.5 |  -7.94${}^{+0.10}_{-0.16}$ |  5.07${}^{+0.04}_{-0.06}$ |  880 ${}^{+45}_{-65}$
LMC8 | 2.5 |  3540${}^{+17}_{-48}$ |  0.30${}^{+0.05}_{-0.18}$ |  +0.5 $\pm$ 0.1 |  -7.61${}^{+0.12}_{-0.21}$ |  4.94${}^{+0.05}_{-0.08}$ |  780 ${}^{+50}_{-70}$

  * •

Notes: The errors presented for $T_{\rm eff,TiO}$ and $A_{V}$ correspond to
the edges of the 1-$\sigma$ contour. For LMC6 and LMC7, we give a lower limit
on the surface gravity. By eye we noticed that +0.5 came close to fitting
these stars, but the real surface gravity is likely in the range from
[+0.5,+1.0], which could not be fitted with the models we have available. For
LMC2, it is evident that the surface gravity is higher than +1.0, and
therefore too compact to be a RSG.

###  3.2 Comparison to $T_{\rm eff}$($J-K$)

The effective temperature one derives depends on the region in the RSG
atmosphere probed, resulting in systematic offsets between different methods.
Numerous approaches have been developed in the last two decades to derive
$T_{\rm eff}$ for RSGs. Each of these approaches uses a different set of
spectral or photometric features, so that depending on the data available, for
both optical and near-IR wavelengths, methods have been established to derive
$T_{\rm eff}$ for RSGs. One approach is to convert photometric data into an
effective temperature through empirical relations. As near-IR wavelengths are
not very susceptible to extinction and variability (van Loon et al. 1999;
Whitelock et al. 2003; Mauron & Josselin 2011) we opted to use the $T_{\rm
eff}$($J-K$) relations (Tabernero et al. 2018; Britavskiy et al. 2019a, b).
These relations probe a temperature closer to the $\tau_{\lambda}=2/3$
continuum region, allowing us to discuss potential weaknesses of the TiO
method as the molecular absorption happens in layers where $\tau_{\lambda}\leq
2/3$. Britavskiy et al. (2019a, b) fitted a sample of several hundreds of RSGs
from Tabernero et al. (2018). In this study, Tabernero et al. (2018) derived
effective temperatures from spectral lines in the $8400-8800\,\textrm{\AA}$
region by averaging the $T_{\rm eff}$ derived from plane-parallel Kurucz
models (Mészáros et al. 2012) and spherical marcs models. Britavskiy et al.
(2019a, b) then obtained a linear relation between the $T_{\rm eff}$ and the
($J-K$) colour for both Magellanic Clouds. However, $T_{\rm eff}$($J-K$) does
not serve as a direct comparison to the $T_{\rm eff,TiO}$ due to the nature of
the derived spectroscopic $T_{\rm eff}$. A direct comparison could only be
made if a relation between the $J-K$ colours and the $T_{\rm eff,TiO}$ were to
be established in the future, or if one would construct a relation between the
$J-K$ colours and the temperature of the marcs models, which will be attempted
in a future work.
We used the relations presented in Britavskiy et al. (2019a, b), which are
applicable for the LMC and SMC, respectively;

$\displaystyle\rm T_{eff}(\it J-K)\rm=-791\cdot(\it J-K\rm_{S})_{0}+4741$ (2)

$\displaystyle\rm T_{eff}(\it J-K)\rm=-1432\cdot(\it J-K\rm_{S})_{0}+5549$ (3)

First we correct for the reddening by calculating $(J-K_{S})_{0}$ from the
relation $(J-K_{S})_{0}=(J-K_{S})-E(J-K)$, where $E(J-K)=0.535E(B-V)$
(Schlegel et al. 1998) and $E(B-V)=\frac{A_{V}}{R_{V}}$, with $R_{V}=2.74$ and
$R_{V}=3.41$ for the SMC and LMC, respectively, so that ultimately the
reddening correction depends on the derived $A_{V}$ from Section 3.1. Although
this affects the integrity of the comparison, as we use our derived $A_{V}$ to
deredden the independent $J-K$ colour, using an independent $A_{V}$ (i.e. from
an extinction map) would not affect the $T_{\rm eff}$($J-K$) much ($\pm
50\,\textrm{K}$ per 0.5 mag, which is arguably small compared to the adopted
error of $\pm 140\,\textrm{K}$). The resulting $T_{\rm eff}$($J-K$) using Eq.
2 $\&$ 3 are presented in the last column of Table 4. In most cases, $T_{\rm
eff}$($J-K$) is higher than $T_{\rm eff,TiO}$, as expected. It is noted that
Eq. 2 $\&$ 3 are valid between $0.8\leq J-K_{S}\leq 1.4\,\ \textrm{mag}$, yet
two of our objects had colours slightly beyond the upper limit, so we
extrapolated the relations accordingly. Britavskiy et al. (2019a, b) indicated
a typical error of $\pm 140\,\textrm{K}$ as the statistical uncertainty on
their linear fit, which we adopted.

Table 4: Spectroscopic and photometric $T_{\rm eff}$ for the M supergiants. ID | $T_{\rm eff,TiO}$ |  $T_{\rm eff}$($J-K$)
---|---|---
| (K) | (K)
SMC3 |  3570${}^{+21}_{-32}$ |  4020 $\pm$ 140
LMC1 |  3600${}^{+68}_{-37}$ |  3685 $\pm$ 140
LMC3 |  3570${}^{+59}_{-46}$ |  3630 $\pm$ 140
LMC4 |  3710${}^{+216}_{-87}$ |  3685 $\pm$ 140
LMC5 |  3620${}^{+53}_{-55}$ |  3830 $\pm$ 140
LMC6 |  3640${}^{+59}_{-51}$ |  3880 $\pm$ 140
LMC7 |  3600${}^{+25}_{-49}$ |  3795 $\pm$ 140
LMC8 |  3540${}^{+17}_{-48}$ |  3805 $\pm$ 140

##  4 Discussion

###  4.1 Hertzsprung-Russell Diagram analysis

We have used evolutionary models to verify the evolutionary status through the
derived $T_{\rm eff}$ and $L$ of our stars. We use MIST models (Dotter 2016;
Choi et al. 2016) for rotating stars ($v/v_{crit}=0.4$) in the range of 8 to
$25\,\textrm{M}_{\odot}$. We proceeded with the MIST models over the Geneva
models, as they better reproduced the RSG population from Yang et al. (2021)
at higher masses. The goal of this section is to provide a crude comparison of
our sample to the tracks, not to derive properties. The tracks terminate at
carbon core depletion, thus approaching the imminent collapse as a supernova,
effectively mapping the entire evolution of a massive star from the main-
sequence to the final evolutionary stages.

We overlay evolutionary tracks on our data in Fig. 9. We also added sources
from Levesque et al. (2006) and Davies et al. (2013), who presented similarly
derived $T_{\rm eff,TiO}$ results, allowing us to make a direct comparison
between the two samples. We note that a considerable fraction of the RSGs lie
in the forbidden zone of the evolutionary tracks right of the Hayashi limit.
Recall that the $T_{\rm eff}$ of these points is determined from the TiO
bands, while the $T_{\rm eff}$ of the models probe the (hotter) continuum
temperature. Indeed, when we construct a HRD using the $T_{\rm eff}$($J-K$)
approach (see Fig. 9, right panel), we see that the discrepancy becomes
minimal and only a few outliers remain. Binary interactions cannot be used to
explain the position in the HRD for these outliers, as the location of the
RSGs in the HRD is largely unaffected by the amount of envelope mass of the
RSGs, except for very low or very high envelope masses (Farrell et al. 2020;
Beasor et al. 2020), it is difficult to infer these scenarios. For moderate
envelope masses, an increase in envelope mass would make an object appear
slightly hotter.
Most of our objects are in agreement with predictions from single star
evolutionary tracks. We identified three potential outliers: LMC2, LMC3 and
SMC3.
i) LMC2 (M4–5 II–III) is located at the lower luminosity end of the population
compared to the RSGs studied in Levesque et al. (2006); Davies et al. (2013).
This is in agreement with the luminosity classification we have assigned.
Using the derived $T_{\rm eff,TiO}$, it appears that LMC2 is too cool to be
explained by any of the evolutionary tracks (see Fig. 9, left panel). This
source is well within the range of tracks if we consider $T_{\rm eff}$($J-K$)
however (see right panel of Fig. 9). Despite the spectral classification, and
considering that the $T_{\rm eff}$($J-K$) is closer to the real temperature of
the photosphere, we find that the properties of this star are consistent with
RSGs, which is in conflict with the spectral classification. It is possible
however that there is a circumstellar contribution to the $K-$band magnitude,
which results in a higher luminosity.
ii) LMC3 (M3 I) is brighter than other M-type supergiants. Since the star
appears to be saturated in archival optical HST images, we could not find
direct evidence of a contaminated $K-$band magnitude by an unresolved binary,
which could potentially affect the derived luminosity. If this star is truly
single, it may be grouped amongst the largest supergiants ($\rm R=1390\rm\,\
R_{\odot}$) known in the LMC, albeit still trailing WOH G64 in size ($\rm
R=1540\rm\,\ R_{\odot}$, log(L∗/L${}_{\odot})\sim$5.45; Levesque et al. 2009).
Even though LMC3 is slightly brighter than WOH G64, we note that the
luminosity of WOH G64 has been carefully corrected for contaminating
contributions in the $K-$band magnitude, while LMC3 is not. Further
investigation is needed to carefully constrain the stellar radius and
luminosity of LMC3 with more certainty. Considering the RSG sample of Davies
et al. (2018), this object is at the high luminosity end of what has been
observed and what is predicted by the Geneva models at Z = 0.006 from
Eggenberger et al. (2021), indicating that this star is a high mass RSG at the
latest stages of evolution.
iii) SMC3 (M2 I) is one of the latest type RSGs in the SMC and lies in the
forbidden zone as can be seen in the left panel of Fig. 15. SMC3 does not
cluster with the sources from similar studies and appears to be cooler and
less luminous. However, if we employ a different methodology to derive the
effective temperature (i.e. $T_{\rm eff}$($J-K$), see right panel of Fig. 15),
this source is not an outlier. We also note that if one considers the
variability of the stars, LMC2 and SMC3 may not strictly be outliers.

Figure 9: Top: HRD indicating the locations of our LMC targets with inverted
red triangles. The $T_{\rm eff}$ for all data points was derived through the
TiO method. Smaller light grey squares and stars are objects from Levesque et
al. (2006) and Davies et al. (2013), respectively. For the two outliers we
have extended the uncertainty assuming a shift of 0.3 mag in the $K-$band
(dotted vertical error bar) instead of only the propagated uncertainty, to
visualize the effect of intrinsic variability. The colour map represents the
central ${{}^{12}{\rm C}}$ mass fraction, while the nodes on the track again
indicate a step of $10^{4}$ years. Bottom: Same as top, but for $T_{\rm eff}$
derived using the $T_{\rm eff}$($J-K$) method. We show the general RSG
population in the background (Yang et al. 2021, black points).

###  4.2 Metallicity dependence

As stated in Section 3.1, we are guided by empirically derived metallicities
from Davies et al. (2015) (log$(Z/Z_{\odot})=-0.35\,\textrm{dex}$ for the LMC
and log$(Z/Z_{\odot})=-0.55\,\textrm{dex}$ for the SMC RSG populations,
respectively) and we proceeded with the spectral modelling using the most
appropriate models available for these values. Given that i) the models and
empirical values are slightly different, ii) we have fitted metal-dependent
molecular bands and iii) the metallicity may not strictly be uniform within
the Magellanic Clouds, we tested if the metallicity could possibly affect the
results. We argue that the use of one single metallicity per galaxy is
justified: i) The general metallicity gradient ([Fe/H] abundance) in the SMC
is flat and shows a small linear decrease outwards for the LMC (Cioni 2009).
ii) Davies et al. (2013) discussed that small variations in metallicity do not
affect the results significantly. Indeed, changing the value of
log$(Z/Z_{\odot})$ significantly (by $\pm$ $0.25\,\textrm{dex}$), did not
affect the effective temperature on a $50\,\textrm{K}$ scale in their study,
and is therefore not expected to affect our results much.

###  4.3 $T_{\rm eff,TiO}$ \- $T_{\rm eff}$($J-K$) discrepancy

The $T_{\rm eff,TiO}$ is generally systematically lower in comparison to
continuum temperature methods (i.e. $T_{\rm eff}$($J-K$) and the $T_{\rm eff}$
derived through SED fitting), due to the region at which these TiO molecules
are formed. Similar to $T_{\rm eff}$, the extinction $A_{V}$ is remarkably
low. By selecting dusty, red sources in the CMDs presented in Section 2, we
expected to find high extinction objects. Using this fitting method however,
we find extinction factors $A_{V}<1$, despite the fact that the extinction
should be the sum of the contributions from foreground extinction, extinction
from the star forming region and circumstellar extinction.
The SED method has been indicated by many authors to be a more reliable method
to extract $T_{\rm eff}$ compared to the TiO bands (Davies et al. 2013;
González-Torà et al. 2021; Davies & Plez 2021). This is due to the fact that
the TiO molecular absorption happens at a region in which $\tau_{\lambda}\neq
2/3$, not near the photosphere. Here, we stress that we have derived the
$T_{\rm eff,TiO}$, not the true $T_{\rm eff}$ of the RSG, which is technically
defined to be at $\tau_{\lambda}=2/3$. In Davies et al. (2013) it has been
discussed extensively that the marcs models do not correctly represent the
radial temperature structure of a real RSG, as the line formation zone of the
TiO bands is pushed to higher altitudes and lower temperatures. To fit the
optical molecular bands properly, detailed 3D models that correctly take into
account convection are needed. Davies et al. (2013) conclude that the
temperature of the layer at which the continuum forms is consistent between 1
and 3D models, hence they argue that, for 1D models, using the continuum is
more reliable. Finally, they suggest that the TiO method may be further
hampered by the fact that TiO molecules may form in the wind region and the
TiO absorption strengths therefore also depend on the wind properties.
Some studies however, provide opposing results. The recent study by Massey et
al. (2021) suggested that Geneva single star (Eggenberger et al. 2008; Ekström
et al. 2012) and BPASS binary (Stanway & Eldridge 2018) models fit observed
RSG/WR ratios the best if the adopted effective temperature is the one derived
through spectroscopy, not photometry (i.e. $T_{\rm eff}$($J-K$)). They argue
that using the spectroscopic derived $T_{\rm eff}$ of stars yields a better
match of their position in the HRD with respect to evolutionary tracks.
However, Massey et al. (2021) also mention that it is uncertain whether the
$T_{\rm eff}$ adopted in the models can be interpreted as equal to what is
derived from a real RSG spectrum.

##  5 Summary and conclusions

We set out to find dusty, luminous evolved massive stars using near and mid-IR
selection criteria. We have spectroscopically confirmed eight new dusty,
luminous RSGs, a new bright giant (LMC2) and a new YSG (SMC2, with strong
H$\alpha$ emission) and three contaminants (i.e. two giants and a LBVc + Hii
region) through the analysis of newly obtained MagE spectra. In addition, we
found a different spectral type for the sgB[e] (SMC1, A0 I[e]), compared to
the previous classifications reported in the literature (B8 I[e] and A1 I[e]),
suggesting this star has changed spectral type more than once over a thirty
year period, although we remain speculative over the physical mechanism
modifying the photosphere of the central star. The RSGs were analysed and
modelled using the marcs models to derive key properties. We approached the
modelling by fitting the optical TiO bands to a precision of $50\,\textrm{K}$,
and have used evolutionary tracks to determine the current evolutionary stage
of these objects. We have derived the physical properties of these red
supergiants. Most sources are very luminous (log(L∗/L${}_{\odot})\geq
5\textrm{\,\ dex}$) and have large radii (R${}_{*}\sim 700-1400$ R⊙). The
brightest RSG (LMC3) appears disconnected from the bulk of the RSG population
in the HRD. This object is located at the tail of the empirical and
theoretical RSG luminosity distribution and is arguably one of the brightest
and largest stars known in the LMC, with log(L∗/L${}_{\odot})\sim
5.45\textrm{\,\ dex}$ and R${}_{*}\sim 1400$ R⊙, resembling WOH G64 in size
and luminosity (Levesque et al. 2009).

We have compared our results and methodology with respect to other studies,
computed the $T_{\rm eff}$($J-K$) of our RSGs and discussed the implications.
The majority of our RSGs have a $T_{\rm eff,TiO}$ that is arguably different
from what is expected through $T_{\rm eff}$($J-K$), although hotter
temperatures are expected for $T_{\rm eff}$($J-K$) as this temperature probes
a deeper layer, closer to $\tau_{\lambda}\sim 2/3$.

Despite the aim to select dusty supergiants, using the methodology as
presented in this paper has resulted in extinction coefficients ($A_{V}<1$
mag) lower than expected. Davies et al. (2013) demonstrated that the derived
extinction from the optical was in almost all cases much lower than that
derived using the near-IR. We share their conclusion that using the optical
TiO bands to derive the extinction is not ideal. We suggest the use of near-IR
spectra or SED fitting for a more reliable result.

We conclude that the selection criteria were successful in selecting cool,
luminous RSGs. However, in the near future we aim to apply the machine
learning algorithm by Maravelias et al. (2022) to improve the selection of
candidate evolved massive stars. The next step is to employ DUSTY models
(Ivezic & Elitzur 1997) to derive mass loss rates from the SED (e.g. Yang et
al. 2022, in prep.). The derivation of the central star properties constrains
the parameter space used with DUSTY. This significantly decreases the amount
of viable models, which reduces degeneracies and computation time. Finally, we
emphasize that LMC3 is a remarkable source and we aim to present a more
detailed analysis of this star in a future study.

## Acknowledgements

We thank the referee, Ben Davies, for the careful reading of the manuscript
and the insightful comments and suggestions that significantly improved the
manuscript. SdW, AZB, FT, GM, MY, EZ acknowledge funding support from the
European Research Council (ERC) under the European Union’s Horizon 2020
research and innovation program (Grant agreement No. 772086). NB acknowledges
support from the postdoctoral program (IPD-STEMA) of Liege University. EZ also
acknowledges support by the Swiss National Science Foundation Professorship
grant (project number PP00P2 176868; PI Tassos Fragos). This paper includes
data gathered with the 6.5 meter Magellan Telescopes located at Las Campanas
Observatory, Chile. This research has made use of NASA’s Astrophysics Data
System. This research has made use of the SIMBAD database, operated at CDS,
Strasbourg, France. This research has made use of the VizieR catalogue access
tool, CDS, Strasbourg, France.

## References

  * Beasor et al. (2021) Beasor, E. R., Davies, B., & Smith, N. 2021, ApJ, 922, 55
  * Beasor et al. (2020) Beasor, E. R., Davies, B., Smith, N., et al. 2020, MNRAS, 492, 5994
  * Bochanski et al. (2009) Bochanski, J. J., Hennawi, J. F., Simcoe, R. A., et al. 2009, PASP, 121, 1409
  * Bonanos et al. (2010) Bonanos, A. Z., Lennon, D. J., Köhlinger, F., et al. 2010, AJ, 140, 416
  * Bonanos et al. (2009) Bonanos, A. Z., Massa, D. L., Sewilo, M., et al. 2009, AJ, 138, 1003
  * Britavskiy et al. (2019a) Britavskiy, N., Lennon, D. J., Patrick, L. R., et al. 2019a, A&A, 624, A128
  * Britavskiy et al. (2019b) Britavskiy, N. E., Bonanos, A. Z., Herrero, A., et al. 2019b, A&A, 631, A95
  * Britavskiy et al. (2014) Britavskiy, N. E., Bonanos, A. Z., Mehner, A., et al. 2014, A&A, 562, A75
  * Brott et al. (2011) Brott, I., de Mink, S. E., Cantiello, M., et al. 2011, A&A, 530, A115
  * Bruch et al. (2021) Bruch, R. J., Gal-Yam, A., Schulze, S., et al. 2021, ApJ, 912, 46
  * Catanzaro (1997) Catanzaro, G. 1997, Ap&SS, 257, 161
  * Chen et al. (2016) Chen, R., Luo, A., Liu, J., & Jiang, B. 2016, AJ, 151, 146
  * Choi et al. (2016) Choi, J., Dotter, A., Conroy, C., et al. 2016, ApJ, 823, 102
  * Cioni (2009) Cioni, M. R. L. 2009, A&A, 506, 1137
  * Cutri et al. (2003) Cutri, R. M., Skrutskie, M. F., van Dyk, S., et al. 2003, VizieR Online Data Catalog, II/246
  * Davies & Beasor (2020) Davies, B. & Beasor, E. R. 2020, MNRAS, 493, 468
  * Davies et al. (2018) Davies, B., Crowther, P. A., & Beasor, E. R. 2018, MNRAS, 478, 3138
  * Davies et al. (2010) Davies, B., Kudritzki, R.-P., & Figer, D. F. 2010, MNRAS, 407, 1203
  * Davies et al. (2015) Davies, B., Kudritzki, R.-P., Gazak, Z., et al. 2015, ApJ, 806, 21
  * Davies et al. (2013) Davies, B., Kudritzki, R.-P., Plez, B., et al. 2013, ApJ, 767, 3
  * Davies & Plez (2021) Davies, B. & Plez, B. 2021, MNRAS, 508, 5757
  * de Jager et al. (1988) de Jager, C., Nieuwenhuijzen, H., & van der Hucht, K. A. 1988, A&AS, 72, 259
  * Dorda et al. (2016) Dorda, R., González-Fernández, C., & Negueruela, I. 2016, A&A, 595, A105
  * Dotter (2016) Dotter, A. 2016, ApJS, 222, 8
  * Egan et al. (2001) Egan, M. P., Van Dyk, S. D., & Price, S. D. 2001, AJ, 122, 1844
  * Eggenberger et al. (2021) Eggenberger, P., Ekström, S., Georgy, C., et al. 2021, A&A, 652, A137
  * Eggenberger et al. (2008) Eggenberger, P., Meynet, G., Maeder, A., et al. 2008, Ap&SS, 316, 43
  * Ekström et al. (2012) Ekström, S., Georgy, C., Eggenberger, P., et al. 2012, A&A, 537, A146
  * Eldridge et al. (2013) Eldridge, J. J., Fraser, M., Smartt, S. J., Maund, J. R., & Crockett, R. M. 2013, MNRAS, 436, 774
  * Evans et al. (2004) Evans, C. J., Howarth, I. D., Irwin, M. J., Burnley, A. W., & Harries, T. J. 2004, MNRAS, 353, 601
  * Farrell et al. (2020) Farrell, E. J., Groh, J. H., Meynet, G., & Eldridge, J. J. 2020, MNRAS, 494, L53
  * Fitzpatrick (1999) Fitzpatrick, E. L. 1999, PASP, 111, 63
  * Gaia Collaboration et al. (2018) Gaia Collaboration, Brown, A. G. A., Vallenari, A., et al. 2018, A&A, 616, A1
  * Gaia Collaboration et al. (2021) Gaia Collaboration, Brown, A. G. A., Vallenari, A., et al. 2021, A&A, 649, A1
  * Gaia Collaboration et al. (2016) Gaia Collaboration, Prusti, T., de Bruijne, J. H. J., et al. 2016, A&A, 595, A1
  * Goldman et al. (2017) Goldman, S. R., van Loon, J. T., Zijlstra, A. A., et al. 2017, MNRAS, 465, 403
  * González-Torà et al. (2021) González-Torà, G., Davies, B., Kudritzki, R.-P., & Plez, B. 2021, MNRAS, 505, 4422
  * Gordon et al. (2003) Gordon, K. D., Clayton, G. C., Misselt, K. A., Landolt, A. U., & Wolff, M. J. 2003, ApJ, 594, 279
  * Graczyk et al. (2014) Graczyk, D., Pietrzyński, G., Thompson, I. B., et al. 2014, ApJ, 780, 59
  * Groenewegen & Blommaert (1998) Groenewegen, M. A. T. & Blommaert, J. A. D. L. 1998, A&A, 332, 25
  * Groenewegen & Sloan (2018) Groenewegen, M. A. T. & Sloan, G. C. 2018, A&A, 609, A114
  * Gustafsson et al. (2008) Gustafsson, B., Edvardsson, B., Eriksson, K., et al. 2008, A&A, 486, 951
  * Humphreys & Jones (2022) Humphreys, R. M. & Jones, T. J. 2022, AJ, 163, 103
  * Ivezic & Elitzur (1997) Ivezic, Z. & Elitzur, M. 1997, MNRAS, 287, 799
  * Jayasinghe et al. (2018) Jayasinghe, T., Kochanek, C. S., Stanek, K. Z., et al. 2018, MNRAS, 477, 3145
  * Jones et al. (2017) Jones, O. C., Woods, P. M., Kemper, F., et al. 2017, MNRAS, 470, 3250
  * Köhler et al. (2015) Köhler, K., Langer, N., de Koter, A., et al. 2015, A&A, 573, A71
  * Kourniotis et al. (2022) Kourniotis, M., Kraus, M., Maryeva, O., Borges Fernandes, M., & Maravelias, G. 2022, MNRAS, 511, 4360
  * Kramida et al. (2020) Kramida, A., Yu. Ralchenko, Reader, J., & and NIST ASD Team. 2020, NIST Atomic Spectra Database (ver. 5.8), [Online]. Available: https://physics.nist.gov/asd [2017, April 9]. National Institute of Standards and Technology, Gaithersburg, MD.
  * Kraus (2019) Kraus, M. 2019, Galaxies, 7, 83
  * Kraus et al. (2008) Kraus, M., Borges Fernandes, M., Kubát, J., & de Araújo, F. X. 2008, A&A, 487, 697
  * Levesque (2017) Levesque, E. M. 2017, Astrophysics of Red Supergiants
  * Levesque et al. (2005) Levesque, E. M., Massey, P., Olsen, K. A. G., et al. 2005, ApJ, 628, 973
  * Levesque et al. (2006) Levesque, E. M., Massey, P., Olsen, K. A. G., et al. 2006, ApJ, 645, 1102
  * Levesque et al. (2009) Levesque, E. M., Massey, P., Plez, B., & Olsen, K. A. G. 2009, AJ, 137, 4744
  * Lindegren et al. (2021) Lindegren, L., Bastian, U., Biermann, M., et al. 2021, A&A, 649, A4
  * Maravelias et al. (2022) Maravelias, G., Bonanos, A. Z., Tramper, F., et al. 2022, arXiv e-prints, arXiv:2203.08125
  * Maravelias et al. (2018) Maravelias, G., Kraus, M., Cidale, L. S., et al. 2018, MNRAS, 480, 320
  * Marshall et al. (2008) Marshall, J. L., Burles, S., Thompson, I. B., et al. 2008, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 7014, Ground-based and Airborne Instrumentation for Astronomy II, 701454
  * Massey (1998) Massey, P. 1998, ApJ, 501, 153
  * Massey et al. (2021) Massey, P., Neugent, K. F., Dorn-Wallenstein, T. Z., et al. 2021, ApJ, 922, 177
  * Massey et al. (2005) Massey, P., Plez, B., Levesque, E. M., et al. 2005, ApJ, 634, 1286
  * Mauerhan et al. (2013) Mauerhan, J. C., Smith, N., Filippenko, A. V., et al. 2013, MNRAS, 430, 1801
  * Mauron & Josselin (2011) Mauron, N. & Josselin, E. 2011, A&A, 526, A156
  * McConnachie (2012) McConnachie, A. W. 2012, AJ, 144, 4
  * Mészáros et al. (2012) Mészáros, S., Allende Prieto, C., Edvardsson, B., et al. 2012, AJ, 144, 120
  * Negueruela et al. (2012) Negueruela, I., Marco, A., González-Fernández, C., et al. 2012, A&A, 547, A15
  * Pietrzyński et al. (2019) Pietrzyński, G., Graczyk, D., Gallenne, A., et al. 2019, Nature, 567, 200
  * Podsiadlowski et al. (1992) Podsiadlowski, P., Joss, P. C., & Hsu, J. J. L. 1992, ApJ, 391, 246
  * Prieto et al. (2008) Prieto, J. L., Kistler, M. D., Thompson, T. A., et al. 2008, ApJ, 681, L9
  * Reid et al. (1990) Reid, N., Tinney, C., & Mould, J. 1990, ApJ, 348, 98
  * Samus’ et al. (2017) Samus’, N. N., Kazarovets, E. V., Durlevich, O. V., Kireeva, N. N., & Pastukhova, E. N. 2017, Astronomy Reports, 61, 80
  * Sanduleak (1989) Sanduleak, N. 1989, AJ, 98, 825
  * Schlegel et al. (1998) Schlegel, D. J., Finkbeiner, D. P., & Davis, M. 1998, ApJ, 500, 525
  * Seale et al. (2009) Seale, J. P., Looney, L. W., Chu, Y.-H., et al. 2009, ApJ, 699, 150
  * Shappee et al. (2014) Shappee, B. J., Prieto, J. L., Grupe, D., et al. 2014, ApJ, 788, 48
  * Skrutskie et al. (2006) Skrutskie, M. F., Cutri, R. M., Stiening, R., et al. 2006, AJ, 131, 1163
  * Smartt (2015) Smartt, S. J. 2015, PASA, 32, e016
  * Smartt et al. (2009) Smartt, S. J., Eldridge, J. J., Crockett, R. M., & Maund, J. R. 2009, MNRAS, 395, 1409
  * Smith (2014) Smith, N. 2014, ARA&A, 52, 487
  * Smith et al. (2009) Smith, N., Hinkle, K. H., & Ryde, N. 2009, AJ, 137, 3558
  * Solf (1978) Solf, J. 1978, A&AS, 34, 409
  * Stanway & Eldridge (2018) Stanway, E. R. & Eldridge, J. J. 2018, MNRAS, 479, 75
  * Tabernero et al. (2018) Tabernero, H. M., Dorda, R., Negueruela, I., & González-Fernández, C. 2018, MNRAS, 476, 3106
  * van der Marel & Kallivayalil (2014) van der Marel, R. P. & Kallivayalil, N. 2014, ApJ, 781, 121
  * van Loon et al. (2005) van Loon, J. T., Cioni, M. R. L., Zijlstra, A. A., & Loup, C. 2005, A&A, 438, 273
  * van Loon et al. (1999) van Loon, J. T., Groenewegen, M. A. T., de Koter, A., et al. 1999, A&A, 351, 559
  * Werner et al. (2004) Werner, M. W., Roellig, T. L., Low, F. J., et al. 2004, ApJS, 154, 1
  * Westerlund et al. (1981) Westerlund, B. E., Olander, N., & Hedin, B. 1981, A&AS, 43, 267
  * Whitelock et al. (1989) Whitelock, P. A., Feast, M. W., Menzies, J. W., & Catchpole, R. M. 1989, MNRAS, 238, 769
  * Whitelock et al. (2003) Whitelock, P. A., Feast, M. W., van Loon, J. T., & Zijlstra, A. A. 2003, MNRAS, 342, 86
  * Wood et al. (1983) Wood, P. R., Bessell, M. S., & Fox, M. W. 1983, ApJ, 272, 99
  * Yang et al. (2021) Yang, M., Bonanos, A. Z., Jiang, B., et al. 2021, A&A, 646, A141
  * Yang et al. (2019) Yang, M., Bonanos, A. Z., Jiang, B.-W., et al. 2019, A&A, 629, A91
  * Yang et al. (2020) Yang, M., Bonanos, A. Z., Jiang, B.-W., et al. 2020, A&A, 639, A116
  * Zapartas et al. (2017) Zapartas, E., de Mink, S. E., Van Dyk, S. D., et al. 2017, ApJ, 842, 125
  * Zickgraf et al. (1992) Zickgraf, F. J., Stahl, O., & Wolf, B. 1992, A&A, 260, 205

##  Appendix A Supplementary Figures

Figure 10: Prominent disk emission lines of SMC1 near the double peaked
H$\alpha$ emission, which is shown in the inset at the top right.

Figure 11: Zoom in on the metal lines used for the A0 I[e] spectral
classification of SMC1.

Figure 12: Zoom in on specific lines selected for the classification for SMC2
(spectral type F8 I). Balmer lines are moderately strong, metal lines are
abundant, and the G-band is weak.

Figure 14: Same as Fig. 7, but for SMC3 and LMC2-8.

Figure 15: Left: same as top panel in Fig. 9, but for the SMC. Right: same as
bottom panel in Fig. 9, but we show the general SMC RSG population in the
background from Yang et al. (2020).