11institutetext: Instituto de Astrofísica de Andalucía (IAA-CSIC), Glorieta de
la Astronomía s/n, 18008 Granada, Spain 22institutetext: Facultad de
Ingeniería y Ciencias, Universidad Adolfo Ibáñez, Av. Diagonal las Torres
2640, Peñalolén, Santiago, Chile 33institutetext: Millennium Institute for
Astrophysics, Chile; 33email<EMAIL_ADDRESS>44institutetext: Leiden
Observatory, Leiden University, Postbus 9513, 2300 RA, Leiden, The Netherlands
55institutetext: Hamburger Sternwarte, Universität Hamburg, Gojenbergsweg 112,
21029 Hamburg, Germany 66institutetext: Instituto de Astrofísica de Canarias
(IAC), 38200 La Laguna, Tenerife, Spain 77institutetext: Deptartamento de
Astrofísica, Universidad de La Laguna (ULL), 38206 La Laguna, Tenerife, Spain
88institutetext: Universitäts-Sternwarte, Ludwig-Maximilians-Universität
München, Scheinerstrasse 1, 81679 München, Germany 99institutetext:
Exzellenzcluster Origins, Boltzmannstraße 2, 85748 Garching, Germany
1010institutetext: Landessternwarte, Zentrum für Astronomie der Universität
Heidelberg, Königstuhl 12, 69117 Heidelberg, Germany 1111institutetext: Centro
de Astrobiología (CSIC-INTA), ESAC, Camino bajo del castillo s/n, 28692
Villanueva de la Cañada, Madrid, Spain 1212institutetext: Institut für
Astrophysik, Georg-August-Universität, FriedrichHund-Platz 1, 37077 Göttingen,
Germany 1313institutetext: Institut de Ciències de l’Espai (ICE, CSIC), Campus
UAB, C/Can Magrans s/n, 08193 Bellaterra, Spain 1414institutetext: Institut
d’Estudis Espacials de Catalunya (IEEC), 08034 Barcelona, Spain

# Moderately misaligned orbit of the warm sub-Saturn HD 332231 b

E. Sedaghati 112233 A. Sánchez-López 44 S. Czesla 55 M. López-Puertas 11 P. J.
Amado 11 E. Palle 6677
K. Molaverdikhani 88991010 J. A. Caballero 1111 L. Nortmann 1212 A.
Quirrenbach 1010 I. Ribas 13131414 A. Reiners 1212

(Received xxx; accepted xxx)

Measurements of exoplanetary orbital obliquity angles, for different classes
of planets, are an essential tool in testing various planet formation
theories. Measurements for those transiting planets on relatively large
orbital periods ($P$ $>$ 10 d) present a rather difficult observational
challenge. Here we present the obliquity measurement for the warm sub-Saturn
planet HD 332231 b, discovered through TESS photometry of sectors 14 and 15.
Through a joint analysis of previously obtained spectroscopic data and our
newly obtained CARMENES transit observations, we estimate the spin-orbit
misalignment angle, , to be $-42.0^{+11.3}_{-10.6}$ $\deg$, which challenges
Laplacian ideals of planet formation. Through the addition of these new data
points obtained with CARMENES, we also derive marginal improvements on other
orbital and bulk parameters for the planet, as compared to previously
published values. We show the robustness of the obliquity measurement through
model comparison with an aligned orbit. Finally, we demonstrate the inability
of the obtained data to probe any possible extended atmosphere of the planet,
due to a lack of precision, and place the atmosphere in the context of a
parameter detection space.

###### Key Words.:

planetary systems – planets and satellites: individual: HD 332231 b – planets
and satellites: atmospheres – methods: observational – techniques:
spectroscopic – techniques: radial velocities

##  1 Introduction

Relative to the invariable plane of Laplace (de Laplace, 1796), one that is
orthogonal to the angular momentum vector and passes through the barycentre,
both terrestrial planets and gas giants in our Solar system exist within well-
aligned orbits. This obliquity angle is within $\sim$ 2 $\deg$ for all
planets, with the exception of Mercury (${}_{\Mercury}\approx 6.3$ $\deg$;
Souami & Souchay, 2012). This picture, however, is drastically different for
the minor bodies in the Solar system. This coplanarity observed in our own
backyard paved the way for understanding the formation and evolution of
planetary systems, tied to the acceptance of the heliocentric model (Kuiper,
1951; Gingerich, 1973).

Discovery of exoplanetary systems has presented a rather more complex picture
of planetary architectures. Transiting exoplanets, those that, from our
vantage point, cross the visible disk of their host stars, permit the
measurement of, among other crucial characteristics, the spin-orbit
misalignment between the orbital plane and the stellar equatorial plane,
projected onto the sky (Triaud, 2018). This is what is referred to as the sky-
projected obliquity angle ( hereafter). Its measurement is performed through
the observations of the Rossiter-McLaughlin effect (Rossiter, 1924;
McLaughlin, 1924), which entails the measurements of the stellar radial
velocity (RV) during the planetary transit (that is when exoplanet transits
are concerned). The anomaly in the measured RV values arises from the
deformation of absorption lines from which they are determined, caused by the
transiting planet occulting either the blue- or red-shifted portion of the
spinning stellar disk. The measurement of this effect is possible thanks to
precision, high dispersion spectrographs at large telescopes that allow for
obtaining high resolution and large Signal to Noise (S/N) spectra at
relatively short temporal sampling. The measurement of has been performed for
a large number of transiting exoplanets, full details of which can be found in
TEPCAT111www.astro.keele.ac.uk/jkt/tepcat/obliquity.html (Southworth, 2011).
These have revealed a surprising diversity in the orbital alignments (e.g.
Queloz et al., 2000; Winn et al., 2005, 2006, 2009; Triaud et al., 2009;
Gandolfi et al., 2010; Mancini et al., 2018; Yu et al., 2018; Lendl et al.,
2020; Sedaghati et al., 2021), contrasting the Laplacian ideals of planets
forming inside a flat disk, coplanar with the stellar equator and staying
there (de Laplace, 1796). A surprising picture that has emerged is that a
significant fraction of those close-in hot-Jupiter regime planets are on
misaligned orbits, as is evident in panel (b) of Fig. 1 (Albrecht et al.,
2012; Dawson, 2014). Furthermore, the spectral type of the host also appears
to play a role, whereby giant planets around hot stars seem to exist on more
oblique orbits (panel (a) of Fig. 1), perhaps pointing to a different, more
chaotic formation history, as compared to their cooler counterparts (Winn et
al., 2010). Hébrard et al. (2011) also point out a lack of planets with masses
¿ 3M${}_{\textrm{jup}}$ on retrograde orbits, the distribution for which is
shown in panel (c) of Fig. 1. Tidal interactions with the host star are also
expected to realign orbits of close-in, massive planets (Zahn, 1977), but lack
of age estimations for the current sample prevents this picture becoming clear
in the parameter space, as is evident in panel (d) of Fig. 1.

Figure 1: Orbital obliquity distributions of those transiting exoplanets with
values determined from the Rossiter-McLaughlin anomaly. (a) shows the
distribution of with the temperature of the host star, where it is quite
evident that giant planets around hot stars are almost entirely on misaligned
orbits. (b) represents the distribution of with the scaled semi-major axis,
whereby giant planets whose orbital semi-major axes are $\lesssim
10\,R_{\star}$ tend to be on highly misaligned orbits. (c) shows again the
same distribution but for the mass of the planet, where the more massive
planets present a larger sample of misaligned orbits. (d) represents the
distribution of as a function of the estimated stellar age, for those targets
where an estimation is available (from NASA exoplanet archive September 2021).
In all panels HD 332231 b is represented as the purple point.

The underlying causes of the aforementioned misalignment for the hot-Jupiter
sample are a subject of debate. One hypothesis suggests that during the high
eccentricity migration, through which hot Jupiters are formed (Petrovich,
2015; Dawson & Johnson, 2018), dynamical interactions contribute to increasing
the obliquity of planetary orbit, pushing it away from the initial equatorial
plane. Such dynamical interactions are attributed to phenomena such as Kozai-
Lidov cycles (Katz et al., 2011; Storch & Lai, 2014) or secular planet-planet
interactions (Naoz et al., 2011). There exists, however, another school of
thought that does not invoke any mechanisms related to planetary migration.
For instance, torques from wide-orbiting companions (Batygin, 2012; Huber et
al., 2013), or a chaotic star formation environment (Bate et al., 2010), among
others, have been suggested as responsible for the observed misalignment of
the hot Jupiter class planets.

Giant planets on wider orbits, those typically with orbital periods ranging
from 10 to 100 days, are classified as “warm”. Warm Jupiters and Saturns are
typically subject to much weaker tidal interactions as compared to their hot
counterparts. Dong et al. (2014) suggest that planets in this regime form
through high eccentricity migration, overcoming the precession caused by
general relativistic effects. Measurement of obliquity angles for a large
sample of giant planets in the warm regime, can somewhat help differentiate
between the two competing schools of thought presented above. Namely, if the
misalignments are due to high eccentricity migration alone, then the majority
of warm giant planets are expected to be on aligned orbits. This would favour
arguments invoking migration mechanisms. However, if the same pattern of
misalignment is observed for both classes, then explanations not involving
planetary migration are favoured. Once a large sample of values for the warm
class is obtained, the interpretation of its distribution might be simpler, as
the orbital architectures are expected to be more pristine. This is due to the
fact that for the close-in hot Jupiters, this information is most likely lost
because of significant dynamical interactions with the host star.

From an observational point of view, measurement of for warm giant planets
presents a more difficult challenge, as compared to the hot class. This is due
to the fact that they are less probable to transit, offer relatively fewer
observable transits for a given period of time, and their long transit
durations make the observations of a complete transit from the ground rather
cumbersome. It must be noted that Rossiter-McLaughlin observations, for now,
are only feasible from ground-based observatories as they require very high
resolution, precision spectroscopy.

In this study we present the obliquity measurement of the warm sub-Saturn
planet, HD 332231 b (Dalba et al., 2020), discovered from NASA’s Transiting
Exoplanet Survey Satellite (TESS; Ricker et al., 2015) photometry of sectors
14 and 15. It has a mass and radius of 0.244 M${}_{\textrm{Jup}}$ and 0.867
R${}_{\textrm{Jup}}$, respectively, and orbits a main-sequence F8, 8.56
m${}_{\textrm{V}}$ star on an 18.71 day circular orbit.

In what follows we present the observations and data preparation, including
telluric correction of the spectra, in section 2, calculation of the RV
values, modelling the Rossiter-McLaughlin effect to determine the obliquity of
the orbit, as well as transmission spectroscopy to search for atmospheric
signals in section 3, discuss the possible implications in section 4, and
finally conclude the study in section 5.

##  2 Observations and data preparation

We observed a single transit of HD 332231 b on the night of 18-10-2020
spanning the entire night, and consequently out of transit observations were
obtained during the two nights on either side. The observations were performed
with the high dispersion échelle spectrograph CARMENES (Quirrenbach et al.,
2014) installed at the 3.5m telescope at the Calar Alto Observatory, in
southern Spain. The instrument consists of two spectrographs, covering the
visible (VIS) and the near-infrared (NIR) channels separately. The VIS channel
covers the range of 0.52 to 0.96 µm, with $\mathcal{R}$ $\sim$ 94 600,
encompassing 55 spectral orders. The NIR channel covers 0.96 to 1.71 µm within
28 orders at $\mathcal{R}$ $\sim$ 80 400. The instrument offers two fibres for
light injection, where fibre A is placed on the astronomical source and fibre
B is typically used for calibration purposes. For our observations we placed
fibre B on sky to monitor and possibly correct for atmospheric emission lines,
as well as remove any potential lunar light contamination. The observations
were performed at 400s exposure time in the VIS channel and 406s in the NIR,
to account for different detector readout times in the two channels. 43
exposures were taken on the night of 18-10, of which only the last 4 were out
of transit, as well as 15 and 14 exposures on the nights before and after,
respectively. The on-sky configuration of the observations are shown in the
top panel of Fig. 2.

Figure 2: Geometry and conditions for CARMENES observations of HD 332231 b
transit. Top panel presents the altitude of the target at the observatory,
whereby the red lines present the duration of the observations, with the cyan
line representing the full exoplanetary transit. Middle panel shows the S/N of
the spectra calculated by the pipeline in the VIS and NIR channels, at the
given specific wavelengths. Bottom panel shows the retrieved water vapour
column density by molecfit in mm, with the 1 confidence shaded in light blue.
In all panels the dark and light shaded regions are the astronomical night and
twilight, respectively, whereby the time axis has been broken for better
presentation and midnight is in UT. Please note that no atmospheric seeing
values were recorded at the observatory for the observations and the cause of
the discrepancy in the S/N values of the first night observations in the VIS
and NIR channels is not entirely clear.

###  2.1 Data reduction

The recorded spectra were subsequently reduced using the caracal reduction
pipeline (CARMENES Reduction And Calibration; Caballero et al., 2016), which
considers all standard astronomical data reduction steps, including a flat-
relative optimal extraction algorithm (Zechmeister et al., 2014). The
wavelength calibration is performed using Th-Ne and U-Ne lamps for the VIS and
NIR channels, respectively, both in combination with a Fabry-Pérot etalon. The
wavelengths are given in the observatory frame and in vacuum. The median S/N
values for the in-transit spectra (the night of 18-10) are 98 ($\sim$ 0.74 µm)
and 90 ($\sim$ 1.23 µm), for the VIS and NIR channels, respectively. The
variations of the S/N values for both channels are presented in the middle
panel of Fig. 2.

###  2.2 Telluric correction

Stellar spectra obtained with ground-based facilities are always inherently
affected by the presence of the atmosphere. This telluric impact manifests
itself as both absorption and emission features, which significantly depend on
the location of the observatory and local meteorological conditions.
Particularly, the atmospheric humidity and temperature profiles along the line
of sight of the observations dictate the specific shapes and depths of those
telluric lines, which are highly time-dependent. For our observations, the
variations of humidity along the line of sight, obtained from the modelling
routine presented below, are shown in the bottom panel of Fig. 2.

In this study, we present results from the initial set of spectra, as well as
those that have been corrected for the telluric absorption features. The
emission features, most prominent in the NIR channel, are simply masked out.
For the correction of the telluric absorption features we use ESO’s molecfit
routines (Kausch et al., 2015; Smette et al., 2015), version 1.5.9, which
synthetically model the lines using a line-by-line radiative transfer model.
In order to model line shapes and depths, the code requires the atmospheric
pressure profile, the initial guess for which it calculates from the humidity
and temperature profiles that it obtains from the
GDAS222https://www.ncei.noaa.gov/products/weather-climate-models/global-data-
assimilation (Global Data Assimilation System) database. The optimal
atmospheric profiles are then obtained through minimization routines, given
the recorded spectra. As molecfit requires one-dimensional spectra, we stitch
all the spectral orders from a single exposure together, creating one spectrum
for each pair of VIS and NIR observation.

For the fitting algorithm of molecfit, we include the main atmospheric optical
and NIR absorbers of O2, H2O, CO2 and CH4, allow the continuum to be modelled
with a local polynomial of the 4th order and relative convergence criterion of
$10^{-10}$. The line shapes are modelled using a Voigt profile, which is
assumed as variable to account for the instrumental profile variation with
wavelength. Example regions for all the fitted species are shown in Fig. 3.

Figure 3: Process of telluric absorption correction of the CARMENES spectra
using molecfit. The panels show zooms into regions where the main absorbers
are, those being O2, H2O, CO2 and CH4, from top to bottom, respectively. The
modelled telluric transmission spectrum is presented as a solid black line,
shifted vertically for clarity, where the zero flux level is indicated with
the dashed line. In the bottom panel some sky emission lines are clearly
visible, which are masked at the analysis stage.

In order to remove the telluric absorption signature from the original two-
dimensional spectra (échelle order-by-order data matrices), the corrected
spectrum obtained from molecfit is resampled on to the wavelength grids of the
individual échelle orders and stored as fits files for further analysis.

In addition to correcting the telluric absorption with molecfit, we also used
the SYSREM algorithm, which is equivalent to a Principal Component Analysis
(PCA) algorithm for unequal uncertainties (Tamuz et al., 2005), to remove the
pseudo-stationary telluric and stellar absorption lines from the spectral
matrices. This approach has been shown to be effective, in a number of high
resolution exoplanetary atmospheric studies, in removing the unwanted star and
Earth atmospheric signals from the recorded spectra (e.g. Birkby et al., 2013;
Sánchez-López et al., 2020). The conclusions reached from this approach are
similar to those obtained from those molecfit-corrected spectra. However, we
refrain from presenting further this analysis, since the relatively low
$K_{\textrm{p}}$ for the planet means that SYSREM quite likely removes a
relatively large portion of any potential atmospheric absorption.

##  3 Analysis

###  3.1 RV calculation

We obtain the RV values of the star for the full set of observations using the
serval pipeline (Zechmeister et al., 2018), applying it to both sets of raw
and telluric corrected spectra, from the VIS and NIR channels. This pipeline
calculates the RVs through a template matching approach (Butler et al., 1996;
Zechmeister et al., 2018), which it constructs by coadding all available
spectra for the target. It provides values which are corrected for
instrumental drift, as well as secular acceleration. Additionally, the
pipeline provides a whole host of spectral diagnostics, some of which are
plotted in Fig. 4. Further to the pipeline calculated values, we also measured
the RVs by fitting a Gaussian peak to the cross-correlation function (CCF;
Pepe et al., 2002) of the observed spectra with a phoenix stellar model
(Husser et al., 2013) matching best the parameters of HD 332231, and observe
no statistically significant differences between the two sets of results (all
pairs of values within 1). Subsequently, throughout this analysis we will
present RV results calculated by the pipeline through the template-matching
approach. These calculated RV values are presented in Fig. 5, whereby a
comparison between the values estimated from the raw and telluric corrected
spectra is made. Considering only the out of transit spectra, the RV jitter
for the target reduces by $\sim$ 0.44 m/s; namely reduced from 8.11 m/s for
the RVs measured from raw spectra to 7.67 m/s for those estimated from the
telluric corrected spectra. The improvements made to stellar RV measurements
through telluric correction of high resolution spectra, have been demonstrated
in a number of studies previously (Artigau et al., 2014; Cunha et al., 2014;
Figueira et al., 2016; Kimeswenger et al., 2021). The calculated RV values
from both sets of spectra are given in Table 2. Furthermore, we ignore the RV
values determined from the NIR channel, which present significantly larger
uncertainties and jitter as compared to the VIS channel data (c.f. Fig. 5).
This is rather expected due to the lower S/N spectra in this channel owing to
the spectral type of the star, as well as the presence of significant telluric
absorption residuals and sky emission lines.

Figure 4: Stellar activity diagnostics calculated by serval. From top to
bottom: Ca II IRT index calculated from the middle line of the triplet at
$\pm$15 km/s; the chromatic RV index which represents the dependence of RV on
wavelength, and is clearly affected by the transiting planet; H activity index
calculated for bin a of $\pm$40 km/s centred on the line core; and the
differential Line Width, again impacted by line deformation caused by the
transiting planet. Overall the star presents no significant sign of
variability due to activity. Figure 5: Radial velocity of HD 332231 measured
during the exoplanet transit, as well as pre- and post-transit epochs. The two
VIS sets of measurements are made for the raw pipeline-corrected spectra
(blue) and telluric-corrected spectra (black). The raw spectra have been
shifted along the time-axis for clarity. Additionally, the RVs measured from
the NIR channel data have also been plotted in the background as red
triangles, which have also been shifted along the x-axis in the opposite
direction.

###  3.2 Obliquity measurement

In order to possibly refine the orbital parameters obtained by Dalba et al.
(2020), as well as model the RM effect, we define a composite Keplerian model
for a circular orbit of a transiting planet, i.e. one that also includes the
RM effect. We assume a circular orbit, as Dalba et al. (2020) obtain a value
for eccentricity consistent with $\sim$ 0 ($0.032^{+0.030}_{-0.022}$). This is
performed using the modelSuite sub-routine from the PyAstronomy python package
(Czesla et al., 2019), which employs the formulation of Ohta et al. (2005) for
the analytical description of the RM effect. In order to estimate best fitting
values for the free parameters, we ran MCMC simulations with 106 steps,
burning in the first 105, through which we sample from the posterior
probability distributions. The best fit model, fitted to the three data sets
used (HIRES/Keck, APF/Lick (Dalba et al., 2020) and CARMENES (this work)), is
shown in panel (a) of Fig. 6, with the phase-folded RV curve, better
representing the RM fit to the in-transit data, shown in panel (b) of the same
figure. The best fit parameters and their uncertainties, the prior
assumptions, as well as a comparison to previously published values are given
in Table 3. It must be noted that the last 4 observations performed during the
night of 18-10, which are out-of-transit (OOT), present large and significant
residuals for the fitted RM model. The nature of this offset is not entirely
clear; it is most likely due to the fact that they were taken for the star at
very high airmass ($1.75\rightarrow 1.98$). In order to investigate whether
these data points bias any of the results for the fitted parameters, we
performed the analysis twice; once with all the data included and a second
time with these data points masked. The results from fitting data with masked
OOT points are given in Figs. 6 and 7, and their equivalents with all the data
points included in Figs. 11 and 12. A comparison between the posterior
distributions, shows that the inclusion of these 4 data points only marginally
impacts the final outcome, and the results are not sensitive to them in any
statistically significant way.

Figure 6: Joint RM+RV analysis of all available RV values for HD 332231,
excluding the last 4 OOT points of CARMENES 18-10 observations, shown as red
points. (a) shows the joint analysis of previously obtained RV data from
spectrographs at Keck (HIRES) and Lick (APF) observatories (Dalba et al.,
2020), together with the data obtained with CARMENES. (b) represents the same
data phase-folded for our derived system parameters from the fit above, where
the inset shows a zoom into the transit region and the fitted RM function to
the CARMENES data. The thin gray lines in the inset are 300 random
realizations drawn from 1 regions of the posteriors. In all panels, the red
line represents the best fit model.

For the fitting of the model, we fixed the scaled semi-major axis
($a/R_{\star}$) and the inclination of the planetary orbit ($i$) to those
values found by Dalba et al. (2020) from fitting the TESS transit light curves
with an analytical model. This method provides much more stringent solutions
as compared to what could be obtained from fitting the RM effect. Furthermore,
the inclination of the stellar rotation axis ($I_{\star}$) is also fixed to 90
$\deg$, since leaving it as a free parameter did not result in convergence of
the posterior. This fact, of course, means that our fitted value for the
projected rotational velocity of the star ($\nu\sin I_{\star}$) is only an
upper limit. All other parameters are taken as free, with the assumed prior
distributions, as well as the best fit values derived from the analysis of
posterior probabilities, given in Table 3. These posterior probability
distributions are plotted in Fig. 7.

Figure 7: Posterior distributions for the fitted parameters of the RV+RM
model, to the data where the last 4 data points of CARMENES 18-10 observations
have been masked. The median values are quoted as solutions, with the 16th and
84th quartiles given as uncertainties, at the top of each column, which are
represented as vertical gray dashed lines graphically. All posteriors are
modelled well with a Gaussian pdf, with the exception of the mid-transit time,
where a double-degenerate solution is derived, as well as and $\nu\sin
I_{\star}$, where a slightly skewed normal distribution is used. In the
parameter-parameter plots, 1-, 2- and 3 levels are presented as solid blue
lines. Table 1: Orbital and physical parameters of HD 332231 b system from
simultaneous RV and RM modelling of all data available.

Parameter | Prior | This study | Dalba et al. (2020)
---|---|---|---
Period [d] |  $\mathcal{U}$(P${}_{\mathrm{D}}-\frac{1}{24}$, P${}_{\mathrm{D}}+\frac{1}{24}$) | $18.70173\pm 0.00028$ |  $18.71204\pm 0.00043$ (PD)
T0 [BJD] |  $\mathcal{U}$(T${}_{\mathrm{D}}-\frac{0.5}{24}$, T${}_{\mathrm{D}}+\frac{0.5}{24}$) | $2459141.34721^{+0.00497}_{-0.00764}$ | 2459141.35287 (TD)
$R_{\textrm{p}}$/$R_{\star}$ |  $\mathcal{N}$(0.06976, 0.012) | $0.06893^{+0.01027}_{-0.00998}$ | $0.06976^{+0.00041}_{-0.00039}$
[$\deg$] |  $\mathcal{U}$($-90$, 90) | $-42.0^{+11.3}_{-10.6}$ | …
$\nu\sin I_{\star}$ [km/s] |  $\mathcal{U}$(0, 40) |  $16.3^{+6.9}_{-4.4}$aa$a$This is an upper limit for the sky-projected stellar rotation velocity and cannot directly be compared to the values from Dalba et al. (2020), as they are derived from spectral analysis. | ($5.3-7.0$) $\pm 1.0$
$K_{\star}$ [m/s] |  $\mathcal{U}$(10, 30) | $17.8\pm 0.5$ | $17.3\pm 1.2$
RV0 [m/s] |  $\mathcal{U}$($-10$, 10) | $-0.91\pm 0.31$ | …
$u_{1}$ |  $\mathcal{U}$(0, 1) | $0.70^{+0.18}_{-0.24}$ | …bb$b$In fitting the transit light curve, Dalba et al. (2020) implement the quadratic limb darkening law, whereas the analytical RM model of Ohta et al. (2005) employs the linear formulation.
$a/R_{\star}$ | … | $<24.21>$ | $24.21^{+0.62}_{-0.78}$
$i$ [$\deg$] | … | $<89.68>$ | $89.68^{+0.22}_{-0.28}$
$I_{\star}$ [$\deg$] | … | $<90>$ | …
$M_{p}$ [$M_{\mathrm{Jup}}$] | … | $0.251\pm 0.017$ | $0.244\pm 0.021$
$\rho_{p}$ [g cm-3] | … | $0.478^{+0.047}_{-0.044}$ | $0.464^{+0.054}_{-0.052}$

333For the period and mid-transit time, values from Dalba et al. (2020) are
used to construct prior distributions. $<>$ symbols indicate fixed parameters.

From this analysis we obtain a value of $-42.0^{+11.3}_{-10.6}$ $\deg$ for the
orbital obliquity angle () of HD 332231 b, with $\nu\sin
I_{\star}=16.3^{+6.9}_{-4.4}$ km/s. The value of has been plotted in the four
parameter spaces presented earlier in Fig. 1. The best fit RM model to the in-
transit data has been highlighted in the inset of the bottom panel in Fig. 6,
where 300 random realizations from the inner 1- posterior distributions have
also been drawn, representing the uncertainty in the final fitted model.

###  3.3 Transmission spectroscopy

We search for possible atmospheric signatures from individually resolved lines
emanating from the exoplanet. This is performed through two distinct, albeit
fundamentally equivalent, approaches, namely; 1) searching for absorption from
strong singular transition lines from species such as Na, K, H or He, or 2)
systematically adding absorption signals from a multitude of weaker lines
through the cross-correlation technique.

The first of the two aforementioned methods has previously been employed in a
number of studies to detect a host of atomic species in the upper atmospheres
of a number of transiting exoplanets. Those detections include neutral Na from
the strong D doublet lines in KELT-9b (Hoeijmakers et al., 2019), KELT-20b
(Casasayas-Barris et al., 2018, 2019; Hoeijmakers et al., 2020a), WASP-12b
(Jensen et al., 2018), WASP-49b (Wyttenbach et al., 2017), WASP-69b
(Casasayas-Barris et al., 2017), WASP-76b (Seidel et al., 2019; Tabernero et
al., 2021), WASP-121b (Cabot et al., 2020; Hoeijmakers et al., 2020b; Borsa et
al., 2021) and WASP-166b (Seidel et al., 2020); the redder line of the K
doublet at 7701 Å in WASP-76b (Tabernero et al., 2021) and WASP-121b (Borsa et
al., 2021); multiple Balmer lines in KELTP-9b (Wyttenbach et al., 2020),
KELT-20b (Casasayas-Barris et al., 2018, 2019), WASP-33b (Borsa et al., 2021;
Yan et al., 2021) and H in WASP-12b (Jensen et al., 2018); and IRT lines of
neutral He in GJ 3470b (Palle et al., 2020), HAT-P11b (Allart et al., 2018),
HD 189733b (Salz et al., 2018), HD 209458b (Alonso-Floriano et al., 2019),
WASP-69b (Nortmann et al., 2018) and WASP-107b (Allart et al., 2019; Kirk et
al., 2020; Spake et al., 2021). We perform this analysis of HD 332231b for all
four of these neutral species, the results for which are shown in Fig. 8. It
is quite evident that the precision at which the data probes the atmosphere is
far below what is needed for any possible detection in this atmosphere, due to
its relatively low equilibrium temperature. Furthermore, a comparison to
atmospheric models suggests that placing any upper limits on the temperature
or abundances from these data is not particularly informative. These models
are calculated with the petitRADTRANS package (Mollière et al., 2019), with
equilibrium chemistry considered through FastChem (Stock et al., 2018),
examples of which are shown in two panels of Fig. 8.

Figure 8: Narrowband transmission spectroscopy of strong transition lines. The
upper panels present the composite out-of-transit stellar spectrum. It must be
noted that, with the exception of He I, all the absorption lines present in
the stellar spectrum correspond to the species being probed, labelled above
each column. The bottom panels present the transmission spectra. The thin gray
lines are the individual in-transit residuals shifted to the planetary rest
frame. The black lines are the weighted mean of all those shifted residuals.
The rest-frame wavelength of each individual transition probed are indicated
with vertical blue dashed lines. For Na and K, two atmospheric models of solar
metallicity and C/O, with equilibrium temperatures of 876 and 1100 K have been
over-plotted.

In addition to this narrow-band approach, we search for the possible presence
of H2O in the atmosphere, through the addition of absorption signals from a
multitude of individually resolved lines in various strong bands of water in
the visible and the near-IR wavelengths. This is performed through the
calculation of the cross-correlation function (CCF) of the telluric-corrected
spectra in the stellar rest-frame with atmospheric templates. The weighted CCF
is defined as:

$C(v,t)=\frac{\sum_{i}^{N}x_{i}(t)\leavevmode\nobreak\
T_{i}(v)}{\sum_{i}^{N}T_{i}(v)}$ (1)

where $T(\nu)$ is the atmospheric template Doppler shifted to velocity $\nu$,
and $x(t)$ is the residual planet spectrum observed at time $t$. It must be
noted that in order to perform the operation in the numerator, both $T$ and
$x$ must necessarily be sampled along the same wavelength grid, which is
typically achieved by calculating the template at very high resolution, and
re-sampling it onto the observational wavelength grid. The cross-correlation
$C$ is subsequently a two-dimensional matrix, with rows representing epochs of
observations and the columns spanning the velocity space being probed. During
this analysis, for a slightly more comprehensive assessment of data quality
and detection limits, we calculate a grid of models varying in equilibrium
temperature and metallicity, and inject them into the residual planet spectra
through multiplication. The abundance profiles for this range of models,
calculated with the equilibrium chemistry code FastChem, are presented in Fig.
9.

Figure 9: H2O abundance profiles (in mass fraction, ) for a range of
atmospheric equilibrium temperatures, $825\leq T_{eq}\leq 1100$ K, and
metallicities, $-2\leq\textrm{[Fe/H]}\leq 2$, calculated for HD 332231 b ,
with C/O fixed to the solar value of 0.55. These profiles, determined with
FastChem, are subsequently used to calculate the models injected into the
data, retrieval of which yields the detection space presented in Fig. 10.

The cross-correlation maps, with the injected signals included, are then
systematically summed for a range of planetary radial velocities to obtain the
so-called ”velocity-velocity” maps, from which detection significances are
determined. We perform this analysis for H2O only, taking advantage of the
deep absorption bands in the near-IR wavelengths, as shown in the top panel of
Fig. 10. The FastChem equilibrium chemistry model for the calculated planetary
$T_{eq}$, solar C/O and metallicity, for a wide range of species is given in
Fig. 13.

As mentioned above, atmospheric models, calculated using petitRADTRANS and
FastChem, spanning in temperature and metallicity are injected into the data
and their corresponding peaks in the velocity maps are compared to the
variance of the entire map to determine the significance at which each point
in the parameter space is detected. A smoothed map of this detection space is
shown in the main panel of Fig. 10, where it is evident, given the quality of
the data, that for the expected equilibrium temperature of the planet, any
presence of H2O in the atmosphere would not be detectable. In fact much higher
temperatures and abundances are necessary for any potential signal to be
accessible. One important caveat to note with this analysis is that we have
assumed a C/O ratio equal to that of the solar value, as well as a cloudless,
hydrostatic and symmetric atmosphere in chemical equilibrium. All such
assumptions of course have impacts, to varying degrees, upon the detectability
space.

Figure 10: Assessing the quality of CARMENES observations for the potential
detection of H2O in the atmosphere of HD 332231 b. Top panel shows an example
of a transmission model calculated with petitRADTRANS (Teq = 876 K and [Fe/H]
= 0), normalised to a baseline estimated from the continuum species only, and
convolved with the instrumental profile. Bottom panel shows the detection
significances for the injected models, smoothed with a Gaussian filter of = 1.
The expected location of HD 332231 b’s atmosphere for its estimated
equilibrium temperature, and assumed solar matallicity, has been indicated as
the red data point.

##  4 Discussion

### Robustness of measurement

From the analysis of the RV values measured for HD 332231, we have presented a
moderately misaligned orbit for the warm sub-Saturn planetary companion of
this star. This obliquity angle is measured through the fitting of an RV+RM
model to all the available data points, whereby the in-transit data are
exclusively obtained through our observations with CARMENES. This angle is
measured at $=-42.0^{+11.3}_{-10.6}$ $\deg$, which in the context of the
relatively small ”warm” planet sample (presented in panel (b) of Fig. 1),
presents a rather rare moderately misaligned orbit.

In order to assess the significance of this measurement, we perform a
statistical comparison between the currently fitted model and one in which is
fixed to 0 $\deg$. We fit all the available data again with an identical
procedure to the one described in section 3.2, whereby the only difference is
that the planetary orbit is assumed to be completely aligned with the stellar
rotation axis ( = 0). We then compare these two models to the data obtained
with CARMENES only (excluding the masked points), to determine the
significance of fit improvement by keeping free. Making the approximation that
the residuals of both models are normally distributed, we calculate the
likelihood ratio, LR, as:

$LR=2\ln{\left(\mathcal{L}(\theta_{\lambda=0})-\mathcal{L}(\theta_{\lambda=free})\right)}=14.38$

where $\mathcal{L}(\theta)$ are the log likelihood values of the models with
fixed and taken as a free parameter. This ratio corresponds to $p=0.00015$,
with the critical value at $P=0.005$ for the degree of freedom 1 being 7.879.
Therefore, we reject the null hypothesis (H0: no statistically significant
difference between the variances of the residuals) at high significance and
conclude that the model with free describes best the data. The F-ratio, which
in this case is a monotone transformation of the likelihood ratio also results
in an identical conclusion.

### Tidal effects

Albrecht et al. (2012) estimate the realignment time scale for planets with
host stars with convective envelopes ($\tau_{CE}$) as:

$\frac{1}{\tau_{CE}}=\frac{1}{10^{10}\textrm{yr}}\left(\frac{M_{p}}{M_{\star}}\right)^{2}\left(\frac{40}{a/R_{\star}}\right)^{6}$
(2)

which is estimated through a calibration of relations derived by Zahn (1977),
with binary star observations. We choose to use this relation, as opposed to
its equivalent for those stars with radiative envelopes, because HD 332231 is
a main sequence star with a measured effective temperature of
$6089^{+97}_{-96}$ K, i.e. spectral type F8 with a radiative core and a
convective envelope. This time scale for HD 332231 b is therefore estimated as
$\sim$ $10^{16}$ yr, far exceeding the age of the universe. Therefore, tidal
interactions with the host star certainly have had no significant impact upon
the obliquity of the planetary orbit.

##  5 Conclusions

In this study we presented high resolution spectroscopic transit observations
of the warm sub-Saturn HD 332231 b, using the CARMENES spectrograph. The
correction of the telluric absorption features using molecfit improves the
precision of the RV values obtained from those recorded spectra. Measurements
of various activity diagnostics indicated that the observations are not likely
to be affected by intrinsic variability of the star, which is further
corroborated by the stability of the long term photometry of the star obtained
by TESS (Dalba et al., 2020).

We measured a projected spin-orbit obliquity angle () of
$-42.0^{+11.3}_{-10.6}$ $\deg$ for HD 332231 b, which presents a moderately
misaligned planet, with a relatively large orbital period. In the context of
planet formation theory presented in section 1, such finding would suggest
that those formation theories not invoking planetary migration are favoured.
However, for such conclusions to be made concrete, certainly a much larger
sample of obliquity measurements for this class of warm transiting planets is
required. This is, of course, made rather difficult by the rarity and duration
of the transits of those planets. We also discussed the robustness of the
measured obliquity value through a statistical model comparison and showed
that the inferred measurement significantly improves the residuals of the
fitted model, as compared to an aligned orbit. Trivially, it was shown that
tidal effects of the host star have not played a role in altering the
alignment of the planetary orbital orientation.

As an aside, we also analysed the spectra for atmospheric signatures through a
couple of distinct methods; namely narrow-band transmission spectroscopy of
strong, singular transition lines and cross-correlation analysis with models
containing a myriad of atmospheric lines of H2O. We showed that given the
precision of the obtained spectra, it is not feasible to detect any
atmospheric characteristics of this planet, especially considering its
relatively low estimated equilibrium temperature of 876K. Through an injection
and retrieval routine, we presented a slice through the detection parameter
space of the atmosphere for a range of metallicities and equilibrium
temperatures. In order for such data to be amenable to atmospheric studies,
higher precision instrumentation coupled to much larger aperture telescopes
are required, a problem which spectrographs such as HIRES at the ELT (Marconi
et al., 2021) will address.

###### Acknowledgements.

CARMENES is an instrument at the Centro Astronómico Hispano-Alemán (CAHA) at
Calar Alto (Almería, Spain), operated jointly by the Junta de Andalucía and
the Instituto de Astrofísica de Andalucía (CSIC). The authors wish to express
their sincere thanks to all members of the Calar Alto staff for their expert
support of the instrument and telescope operation. We acknowledge financial
support from the Agencia Estatal de Investigación of the Ministerio de
Ciencia, Innovación y Universidades through projects Ref.
PID2019-110689RB-I00/AEI/10.13039/501100011033 and the Centre of Excellence
“Severo Ochoa” award to the Instituto de Astrofísica de Andalucía
(SEV-2017-0709). We also acknowledge the use of the ExoAtmospheres database
during the preparation of this work and thank the team at the IAC for their
excellent work preparing such a useful database. ES acknowledges support from
ANID - Millennium Science Initiative - ICN12_009. KM acknowledges support from
the Excellence Cluster ORIGINS, which is funded by the Deutsche
Forschungsgemeinschaft (DFG, German Research Foundation) under Germany’s
Excellence Strategy - EXC-2094 - 390783311. Beside those explicitly mentioned
in the manuscript, throughout this work the following python packages have
been used: astropy (Astropy Collaboration et al., 2018), corner (Foreman-
Mackey, 2016), matplotlib (Hunter, 2007), NumPy (Harris et al., 2020), pandas
(pandas development team, 2020), PyMC (Salvatier et al., 2016) and SciPy
(Virtanen et al., 2020).

## References

  * Albrecht et al. (2012) Albrecht, S., Winn, J. N., Johnson, J. A., et al. 2012, ApJ, 757, 18
  * Allart et al. (2019) Allart, R., Bourrier, V., Lovis, C., et al. 2019, A&A, 623, A58
  * Allart et al. (2018) Allart, R., Bourrier, V., Lovis, C., et al. 2018, Science, 362, 1384
  * Alonso-Floriano et al. (2019) Alonso-Floriano, F. J., Snellen, I. A. G., Czesla, S., et al. 2019, A&A, 629, A110
  * Artigau et al. (2014) Artigau, É., Astudillo-Defru, N., Delfosse, X., et al. 2014, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 9149, Observatory Operations: Strategies, Processes, and Systems V, ed. A. B. Peck, C. R. Benn, & R. L. Seaman, 914905
  * Astropy Collaboration et al. (2018) Astropy Collaboration, Price-Whelan, A. M., Sipőcz, B. M., et al. 2018, AJ, 156, 123
  * Bate et al. (2010) Bate, M. R., Lodato, G., & Pringle, J. E. 2010, MNRAS, 401, 1505
  * Batygin (2012) Batygin, K. 2012, Nature, 491, 418
  * Birkby et al. (2013) Birkby, J. L., de Kok, R. J., Brogi, M., et al. 2013, MNRAS, 436, L35
  * Borsa et al. (2021) Borsa, F., Allart, R., Casasayas-Barris, N., et al. 2021, A&A, 645, A24
  * Butler et al. (1996) Butler, R. P., Marcy, G. W., Williams, E., et al. 1996, PASP, 108, 500
  * Caballero et al. (2016) Caballero, J. A., Guàrdia, J., López del Fresno, M., et al. 2016, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 9910, Observatory Operations: Strategies, Processes, and Systems VI, ed. A. B. Peck, R. L. Seaman, & C. R. Benn, 99100E
  * Cabot et al. (2020) Cabot, S. H. C., Madhusudhan, N., Welbanks, L., Piette, A., & Gandhi, S. 2020, MNRAS, 494, 363
  * Casasayas-Barris et al. (2017) Casasayas-Barris, N., Palle, E., Nowak, G., et al. 2017, A&A, 608, A135
  * Casasayas-Barris et al. (2018) Casasayas-Barris, N., Pallé, E., Yan, F., et al. 2018, A&A, 616, A151
  * Casasayas-Barris et al. (2019) Casasayas-Barris, N., Pallé, E., Yan, F., et al. 2019, A&A, 628, A9
  * Cunha et al. (2014) Cunha, D., Santos, N. C., Figueira, P., et al. 2014, A&A, 568, A35
  * Czesla et al. (2019) Czesla, S., Schröter, S., Schneider, C. P., et al. 2019, PyA: Python astronomy-related packages
  * Dalba et al. (2020) Dalba, P. A., Gupta, A. F., Rodriguez, J. E., et al. 2020, AJ, 159, 241
  * Dawson (2014) Dawson, R. I. 2014, ApJ, 790, L31
  * Dawson & Johnson (2018) Dawson, R. I. & Johnson, J. A. 2018, ARA&A, 56, 175
  * de Laplace (1796) de Laplace, P. S. 1796, Exposition du système du monde
  * Dong et al. (2014) Dong, S., Katz, B., & Socrates, A. 2014, ApJ, 781, L5
  * Figueira et al. (2016) Figueira, P., Adibekyan, V. Z., Oshagh, M., et al. 2016, A&A, 586, A101
  * Foreman-Mackey (2016) Foreman-Mackey, D. 2016, The Journal of Open Source Software, 1, 24
  * Gandolfi et al. (2010) Gandolfi, D., Hébrard, G., Alonso, R., et al. 2010, A&A, 524, A55
  * Gingerich (1973) Gingerich, O. 1973, Proceedings of the American Philosophical Society, 117, 513
  * Harris et al. (2020) Harris, C. R., Millman, K. J., van der Walt, S. J., et al. 2020, Nature, 585, 357
  * Hébrard et al. (2011) Hébrard, G., Evans, T. M., Alonso, R., et al. 2011, A&A, 533, A130
  * Hoeijmakers et al. (2020a) Hoeijmakers, H. J., Cabot, S. H. C., Zhao, L., et al. 2020a, A&A, 641, A120
  * Hoeijmakers et al. (2019) Hoeijmakers, H. J., Ehrenreich, D., Kitzmann, D., et al. 2019, A&A, 627, A165
  * Hoeijmakers et al. (2020b) Hoeijmakers, H. J., Seidel, J. V., Pino, L., et al. 2020b, A&A, 641, A123
  * Huber et al. (2013) Huber, D., Carter, J. A., Barbieri, M., et al. 2013, Science, 342, 331
  * Hunter (2007) Hunter, J. D. 2007, Computing in Science & Engineering, 9, 90
  * Husser et al. (2013) Husser, T. O., Wende-von Berg, S., Dreizler, S., et al. 2013, A&A, 553, A6
  * Jensen et al. (2018) Jensen, A. G., Cauley, P. W., Redfield, S., Cochran, W. D., & Endl, M. 2018, AJ, 156, 154
  * Katz et al. (2011) Katz, B., Dong, S., & Malhotra, R. 2011, Phys. Rev. Lett., 107, 181101
  * Kausch et al. (2015) Kausch, W., Noll, S., Smette, A., et al. 2015, A&A, 576, A78
  * Kimeswenger et al. (2021) Kimeswenger, S., Rainer, M., Przybilla, N., & Kausch, W. 2021, AJ, 161, 66
  * Kirk et al. (2020) Kirk, J., Alam, M. K., López-Morales, M., & Zeng, L. 2020, AJ, 159, 115
  * Kuiper (1951) Kuiper, G. P. 1951, Proceedings of the National Academy of Sciences of the United States of America, 37, 1
  * Lendl et al. (2020) Lendl, M., Csizmadia, S., Deline, A., et al. 2020, A&A, 643, A94
  * Mancini et al. (2018) Mancini, L., Esposito, M., Covino, E., et al. 2018, A&A, 613, A41
  * Marconi et al. (2021) Marconi, A., Abreu, M., Adibekyan, V., et al. 2021, The Messenger, 182, 27
  * McLaughlin (1924) McLaughlin, D. B. 1924, ApJ, 60, 22
  * Mollière et al. (2019) Mollière, P., Wardenier, J. P., van Boekel, R., et al. 2019, A&A, 627, A67
  * Naoz et al. (2011) Naoz, S., Farr, W. M., Lithwick, Y., Rasio, F. A., & Teyssandier, J. 2011, Nature, 473, 187
  * Nortmann et al. (2018) Nortmann, L., Pallé, E., Salz, M., et al. 2018, Science, 362, 1388
  * Ohta et al. (2005) Ohta, Y., Taruya, A., & Suto, Y. 2005, ApJ, 622, 1118
  * Palle et al. (2020) Palle, E., Nortmann, L., Casasayas-Barris, N., et al. 2020, A&A, 638, A61
  * pandas development team (2020) pandas development team, T. 2020, pandas-dev/pandas: Pandas
  * Pepe et al. (2002) Pepe, F., Mayor, M., Galland, F., et al. 2002, A&A, 388, 632
  * Petrovich (2015) Petrovich, C. 2015, ApJ, 805, 75
  * Queloz et al. (2000) Queloz, D., Eggenberger, A., Mayor, M., et al. 2000, A&A, 359, L13
  * Quirrenbach et al. (2014) Quirrenbach, A., Amado, P. J., Caballero, J. A., et al. 2014, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 9147, Ground-based and Airborne Instrumentation for Astronomy V, ed. S. K. Ramsay, I. S. McLean, & H. Takami, 91471F
  * Ricker et al. (2015) Ricker, G. R., Winn, J. N., Vanderspek, R., et al. 2015, Journal of Astronomical Telescopes, Instruments, and Systems, 1, 014003
  * Rossiter (1924) Rossiter, R. A. 1924, ApJ, 60, 15
  * Salvatier et al. (2016) Salvatier, J., Wiecki, T. V., & Fonnesbeck, C. 2016, PeerJ Computer Science, 2, e55
  * Salz et al. (2018) Salz, M., Czesla, S., Schneider, P. C., et al. 2018, A&A, 620, A97
  * Sánchez-López et al. (2020) Sánchez-López, A., López-Puertas, M., Snellen, I. A. G., et al. 2020, A&A, 643, A24
  * Sedaghati et al. (2021) Sedaghati, E., MacDonald, R. J., Casasayas-Barris, N., et al. 2021, MNRAS, 505, 435
  * Seidel et al. (2020) Seidel, J. V., Ehrenreich, D., Bourrier, V., et al. 2020, A&A, 641, L7
  * Seidel et al. (2019) Seidel, J. V., Ehrenreich, D., Wyttenbach, A., et al. 2019, A&A, 623, A166
  * Smette et al. (2015) Smette, A., Sana, H., Noll, S., et al. 2015, A&A, 576, A77
  * Souami & Souchay (2012) Souami, D. & Souchay, J. 2012, A&A, 543, A133
  * Southworth (2011) Southworth, J. 2011, MNRAS, 417, 2166
  * Spake et al. (2021) Spake, J. J., Oklopčić, A., & Hillenbrand, L. A. 2021, arXiv e-prints, arXiv:2107.08999
  * Stock et al. (2018) Stock, J. W., Kitzmann, D., Patzer, A. B. C., & Sedlmayr, E. 2018, MNRAS, 479, 865
  * Storch & Lai (2014) Storch, N. I. & Lai, D. 2014, MNRAS, 438, 1526
  * Tabernero et al. (2021) Tabernero, H. M., Zapatero Osorio, M. R., Allart, R., et al. 2021, A&A, 646, A158
  * Tamuz et al. (2005) Tamuz, O., Mazeh, T., & Zucker, S. 2005, MNRAS, 356, 1466
  * Triaud (2018) Triaud, A. H. M. J. 2018, The Rossiter-McLaughlin Effect in Exoplanet Research, ed. H. J. Deeg & J. A. Belmonte, 2
  * Triaud et al. (2009) Triaud, A. H. M. J., Queloz, D., Bouchy, F., et al. 2009, A&A, 506, 377
  * Virtanen et al. (2020) Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, Nature Methods, 17, 261
  * Winn et al. (2010) Winn, J. N., Fabrycky, D., Albrecht, S., & Johnson, J. A. 2010, ApJ, 718, L145
  * Winn et al. (2009) Winn, J. N., Johnson, J. A., Albrecht, S., et al. 2009, ApJ, 703, L99
  * Winn et al. (2006) Winn, J. N., Johnson, J. A., Marcy, G. W., et al. 2006, ApJ, 653, L69
  * Winn et al. (2005) Winn, J. N., Noyes, R. W., Holman, M. J., et al. 2005, ApJ, 631, 1215
  * Wyttenbach et al. (2017) Wyttenbach, A., Lovis, C., Ehrenreich, D., et al. 2017, A&A, 602, A36
  * Wyttenbach et al. (2020) Wyttenbach, A., Mollière, P., Ehrenreich, D., et al. 2020, A&A, 638, A87
  * Yan et al. (2021) Yan, F., Wyttenbach, A., Casasayas-Barris, N., et al. 2021, A&A, 645, A22
  * Yu et al. (2018) Yu, L., Zhou, G., Rodriguez, J. E., et al. 2018, AJ, 156, 250
  * Zahn (1977) Zahn, J. P. 1977, A&A, 500, 121
  * Zechmeister et al. (2014) Zechmeister, M., Anglada-Escudé, G., & Reiners, A. 2014, A&A, 561, A59
  * Zechmeister et al. (2018) Zechmeister, M., Reiners, A., Amado, P. J., et al. 2018, A&A, 609, A12

##  Appendix A Calculated RVs

RV values calculated by the serval pipeline from both raw and telluric
corrected spectra.

Table 2: Caracal calculated RVs. Time [BJD] | RV
---|---
| Raw spectra | Telluric-corrected
| [m/s] | [m/s]
2459140.42205 |  13.57 $\pm$ 7.51 |  6.61 $\pm$ 8.40
2459140.42713 |  1.74 $\pm$ 7.17 |  -4.04 $\pm$ 7.63
2459140.43232 |  13.80 $\pm$ 6.18 |  14.42 $\pm$ 7.15
2459140.43736 |  -1.89 $\pm$ 7.19 |  -2.90 $\pm$ 7.71
2459140.44257 |  6.84 $\pm$ 7.36 |  -1.43 $\pm$ 6.93
2459140.45008 |  0.41 $\pm$ 5.76 |  -3.16 $\pm$ 7.45
2459140.45508 |  4.62 $\pm$ 6.62 |  1.90 $\pm$ 7.40
2459140.46005 |  5.41 $\pm$ 8.15 |  -5.29 $\pm$ 7.08
2459140.46531 |  5.83 $\pm$ 7.48 |  5.88 $\pm$ 9.07
2459140.47045 |  12.57 $\pm$ 6.42 |  4.61 $\pm$ 6.91
2459140.47712 |  2.29 $\pm$ 8.34 |  -2.46 $\pm$ 8.81
2459140.48229 |  2.61 $\pm$ 7.70 |  -4.38 $\pm$ 8.78
2459140.48746 |  8.01 $\pm$ 8.23 |  0.56 $\pm$ 8.60
2459140.49260 |  -6.95 $\pm$ 6.87 | -13.61 $\pm$ 7.66
2459140.49770 |  8.59 $\pm$ 7.15 |  4.00 $\pm$ 7.36
2459141.26203 |  58.04 $\pm$ 8.72 |  62.54 $\pm$ 8.00
2459141.26715 |  39.72 $\pm$ 8.88 |  46.71 $\pm$ 8.71
2459141.27239 |  37.27 $\pm$ 8.62 |  41.86 $\pm$ 8.10
2459141.28010 |  40.74 $\pm$ 8.69 |  42.53 $\pm$ 8.83
2459141.28502 |  44.57 $\pm$ 8.92 |  48.66 $\pm$ 8.88
2459141.29021 |  33.54 $\pm$ 8.74 |  36.68 $\pm$ 8.72
2459141.29538 |  33.16 $\pm$ 8.62 |  39.14 $\pm$ 8.12
2459141.30233 |  22.48 $\pm$ 8.22 |  31.19 $\pm$ 7.93
2459141.30736 |  17.60 $\pm$ 7.72 |  25.73 $\pm$ 8.50
2459141.31257 |  17.35 $\pm$ 8.54 |  22.31 $\pm$ 7.87
2459141.31760 |  23.95 $\pm$ 9.73 |  29.82 $\pm$ 8.17
2459141.32402 |  12.00 $\pm$ 9.16 |  18.55 $\pm$ 8.53
2459141.32908 |  16.88 $\pm$ 8.37 |  21.29 $\pm$ 7.66
2459141.33426 |  10.93 $\pm$ 7.09 |  19.36 $\pm$ 6.48
2459141.33930 |  -3.67 $\pm$ 9.50 |  4.87 $\pm$ 9.23
2459141.34620 |  3.01 $\pm$ 8.30 |  12.06 $\pm$ 8.25
2459141.35135 |  -5.85 $\pm$ 8.45 |  2.23 $\pm$ 8.47
2459141.35644 | -18.26 $\pm$ 8.94 |  -5.37 $\pm$ 7.96
2459141.36158 | -12.92 $\pm$ 8.68 |  -7.08 $\pm$ 8.23
2459141.36812 | -14.23 $\pm$ 8.86 |  -3.88 $\pm$ 8.40
2459141.37326 | -15.08 $\pm$ 8.14 |  -6.91 $\pm$ 8.37
2459141.37834 | -21.00 $\pm$ 7.86 | -11.15 $\pm$ 7.60
2459141.38351 | -20.50 $\pm$ 8.74 | -11.91 $\pm$ 8.63
2459141.39060 | -31.54 $\pm$ 7.81 | -22.62 $\pm$ 8.45
2459141.39564 | -39.09 $\pm$ 7.41 | -27.28 $\pm$ 7.36
2459141.40081 | -40.27 $\pm$ 8.92 | -30.50 $\pm$ 9.05
Table A.1: continued. Time [BJD] | RV
---|---
| Raw spectra | Telluric-corrected
| [m/s] | [m/s]
2459141.40591 | -25.35 $\pm$ 9.21 | -17.93 $\pm$ 9.42
2459141.41225 | -35.80 $\pm$ 6.61 | -26.76 $\pm$ 8.08
2459141.41732 | -31.45 $\pm$ 7.90 | -25.13 $\pm$ 9.21
2459141.42243 | -36.51 $\pm$ 8.27 | -34.07 $\pm$ 8.79
2459141.42755 | -40.56 $\pm$ 8.60 | -33.83 $\pm$ 9.21
2459141.43373 | -40.67 $\pm$ 9.00 | -33.90 $\pm$ 8.68
2459141.43889 | -36.43 $\pm$ 10.79 | -29.48 $\pm$ 11.23
2459141.44400 | -38.88 $\pm$ 10.55 | -33.39 $\pm$ 11.35
2459141.44920 | -40.02 $\pm$ 10.58 | -35.86 $\pm$ 10.39
2459141.45685 | -23.75 $\pm$ 12.46 | -15.97 $\pm$ 11.48
2459141.46185 | -29.52 $\pm$ 11.09 | -21.37 $\pm$ 10.61
2459141.46801 | -33.96 $\pm$ 11.97 | -30.09 $\pm$ 11.43
2459141.47317 | -26.48 $\pm$ 10.66 | -22.84 $\pm$ 11.26
2459141.47954 | -29.19 $\pm$ 11.41 | -19.92 $\pm$ 13.34
2459141.48453 | -27.43 $\pm$ 12.16 | -25.61 $\pm$ 10.99
2459141.49077 | -33.93 $\pm$ 11.18 | -32.91 $\pm$ 10.22
2459141.49586 | -27.93 $\pm$ 13.34 | -36.14 $\pm$ 12.77
2459142.25492 |  19.34 $\pm$ 10.82 |  2.72 $\pm$ 9.18
2459142.26009 |  18.24 $\pm$ 10.84 |  2.39 $\pm$ 8.65
2459142.26518 |  7.27 $\pm$ 10.22 |  -7.42 $\pm$ 8.81
2459142.27040 |  22.76 $\pm$ 9.49 |  3.06 $\pm$ 7.90
2459142.27841 |  9.91 $\pm$ 9.84 | -14.59 $\pm$ 10.22
2459142.28337 |  3.55 $\pm$ 11.71 | -11.81 $\pm$ 10.83
2459142.28829 |  12.84 $\pm$ 8.58 |  -1.98 $\pm$ 9.59
2459142.29390 |  11.85 $\pm$ 11.16 |  -6.01 $\pm$ 9.07
2459142.29980 |  7.30 $\pm$ 10.20 |  -5.24 $\pm$ 9.22
2459142.30499 |  -3.88 $\pm$ 10.28 | -17.11 $\pm$ 9.04
2459142.31003 |  -5.55 $\pm$ 9.38 | -14.60 $\pm$ 9.45
2459142.31533 |  5.89 $\pm$ 11.91 | -13.03 $\pm$ 10.55
2459142.32160 |  -7.84 $\pm$ 11.65 | -13.08 $\pm$ 10.13
2459142.32665 | -11.09 $\pm$ 11.95 | -13.09 $\pm$ 12.30

##  Appendix B Model fit to all data

Here we present the fitted RV+RM model, as well as the posterior
distributions, where all the data points are considered. Namely, the four OOT
data points have not been masked.

Figure 11: Same as Fig. 6 but with all the data points considered. Figure 12:
Same as Fig. 7 but for the fitting routine where all the data points are
considered.

##  Appendix C FastChem abundances

In this appendix we present the equilibrium chemistry model calculated using
the FastChem code, for an atmosphere with solar C/O and metallicity, and
$T_{eq}$ equal to that estimated for HD 332231 b, 876 K.

Figure 13: Equilibrium chemistry model for HD 332231 b calculated with
FastChem, for an atmosphere with $T_{eq}$ of 876 K, and solar C/O and
metallicity. For clarity the most abundant species are included, with the
molecular species plotted as solid lines and those ionic/atomic species
plotted with dashed lines.