{'fraction_non_alphanumeric': 0.03599351975769529, 'fraction_numerical': 0.03157357188138339, 'mean_word_length': 4.462581762216237, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
The Rossiter-McLaughlin effect Revolutions: An ultra-short period planet and a warm mini-Neptune on perpendicular orbits

Comparisons of the alignment of exoplanets with a common host star can be used to distinguish among concurrent evolution scenarios. However, multi-planet systems usually host mini-Neptunes and super-Earths, whose size make orbital architecture measurements challenging. We introduce the Rossiter-McLaughlin effect Revolutions technique, which can access spin-orbit angles of small planets by exploiting the full information contained in spectral transit time series. We validated the technique on published HARPS-N data of the mini-Neptune HD3167c, refining its high sky-projected spin-orbit angle (-108.9+5.4-5.5 deg), and we applied it to new ESPRESSO observations of the super-Earth HD3167b, revealing an aligned orbit (-6.6+6.6-7.9 deg). Surprisingly different variations in the contrast of the stellar lines occulted by the planets can be reconciled with a latitudinal dependence of the stellar line shape. In this scenario, a joint fit to both datasets constrains the inclination of the star (111.6+3.1-3.3 deg) and the 3D spin-orbit angles of HD3167b (29.5+7.2-9.4 deg) and HD3167c (107.7+5.1-4.9 deg). The projected spin-orbit angles do not depend on the model for the line contrast variations, and so, with a mutual inclination of 102.3+7.4-8.0 deg, we conclude that the two planets are on perpendicular orbits. This could be explained by HD3167b being strongly coupled to the star and retaining its primordial alignment, whereas HD3167c would have been brought to a nearly polar orbit via secular gravitational interactions with an outer companion. Follow-up observations and dynamical evolution simulations are required to search for this companion and explore this scenario. HD3167b is the smallest exoplanet with a confirmed spectroscopic Rossiter-McLaughlin signal. Our new technique opens the way to determining the orbital architectures of the super-Earth and Earth-sized planet populations.


Introduction
Measurements of star-planet alignments are essential to improving our understanding of exoplanet formation and evolution, especially for planets that migrated close to their star. Disk-driven migration is expected to conserve the initial alignment between the angular momentums of the protoplanetary disk and of the planetary orbits (e.g., Winn & Fabrycky 2015). A variety of formation and migration pathways can, however, lead to misalignments between the spins of a star and of its planets' orbits, namely: a primordial tilt of the star or the protoplanetary disk (e.g., Lai et al. 2011, Bate et al. 2010, Fielding et al. 2015; a tidal torque on the protoplanetary disk from a neighboring star (e.g., Lai 2014, Batygin 2012Zanazzi & Lai 2018); a tidal torque on the inner planetary system from an outer companion (e.g., Huber et al. 2013); and scattering or secular interactions between the inner planets and an outer planetary or stellar companion (e.g., Wu & Murray 2003, Chatterjee et al. 2008, Fabrycky & Tremaine 2007, Teyssandier et al. 2013. In this context, determining the orbital architecture of multi-planet systems is of great relevance as the mutual inclinations between the planets, and with the star, can help distinguish among various formation and migration scenarios. Early measurements revealed coplanar orbits that are well aligned with the stellar equator in multi-planet systems (e.g. Figueira et al. 2012;Hirano et al. 2012;Sanchis-Ojeda et al. 2012;Albrecht et al. 2013, Chaplin et al. 2013, Van Eylen et al. 2014, standing in contrast to the broad distribution of misalignments observed for hot Jupiters (e.g., Naoz et al. 2012, Albrecht et al. 2012, Davies et al. 2014) and suggesting that their orbital planes still trace a primordial alignment with the protoplanetary disk. Interestingly, spin-orbit measurements in young systems (ă100 Myr) have shown prograde and aligned orbits (Palle et al. 2020;Zhou et al. 2020;Mann et al. 2020). Other studies have unveiled substantial misalignments in multi-planet systems (Huber et al. 2013, Walkowicz & Basri 2013Hirano et al. 2014, Hjorth et al. 2021, Zhang et al. 2021), hinting at more complex formation processes (e.g., Spalding & Batygin 2014;Spalding & Millholland 2020) or dynamical histories (e.g. Gratia & Fabrycky 2017).
Analyses of an exoplanet transit at high temporal cadence and with high spectral resolution can be used to resolve the path of the planet across the star and thus measure its obliquity, defined as the angle between the normal to the planetary orbit and the rotation axis of the star (usually only known through its projection λ in the plane of the sky). The occultation of the rotating stellar photosphere by the planet creates a signature in the disk-integrated lines of the star, whose spectral position is linked with the spatial position of the occulted stellar regions (the so-called Rossiter-McLaughlin, or RM, effect;Rossiter 1924;McLaughlin 1924). This distortion of the stellar lines results in an anomalous deviation of their RV centroid from the Keplerian curve, whose analysis has traditionally been used to study the RM effect and derive  (see review by Triaud 2018). Alternative techniques consist of directly fitting the disk-integrated lines with a spectral model accounting for the signature of the planet (Doppler tomography; e.g., Collier Bourrier et al. 2015;Temple et al. 2019) or retrieving the stellar line occulted by the planet and measuring its local RV centroid directly (Reloaded RM effect; e.g., Cegla et al. 2016b;Bourrier et al. 2020, Kunovac Hodžić et al. 2021. By using the planet as a probe of the stellar surface, the latter technique provides direct access to its properties and reduces biases in the derived obliquity (Cegla et al. 2016b, Bourrier et al. 2017).
In the current paper, we investigate the orbital architecture of the HD 3167 system, host to three known planets (Sect. 2). We used high-resolution optical spectroscopy to analyze the Rossiter-McLaughlin effect of HD 3167b and c (Sect. 3), using a new technique that we describe and apply in Sect. 4. We discuss the orbital architecture and possible history of the HD 3167 system in light of our results in Sect. 5 and we present our conclusions in Sect. 6.

HD 3167 system
HD 3167 is a nearby (47 pc) and bright (V = 9) K0-type star hosting three known planets (Vanderburg et al. 2016, Gandolfi et al. 2017, Christiansen et al. 2017). The innermost planet, HD 3167b, is a super-Earth (R p = 1.70`0 .18 0.15 R C , M p = 5.02˘0.38 M C ) on an ultra-short period of 0.96 d and the outermost one, HD 3167c, is a mini-Neptune (R p = 3.01`0 .42 0.28 R C , M p = 9.80`1 .30 1.24 M C ) on a 29.84 d orbit. While HD 3167b and HD 3167c transit, this is not the case for HD 3167d (M p sin i "6.90˘0.71 M C ), which orbits in between them on a 8.51 d orbit. HD 3167c is on a highly misaligned orbit around the star and dynamical analysis has shown that its orbit is likely nearly coplanar with that of HD 3167d, albeit with a mutual inclination, likely ranging between 2.3 and 21˝ (Dalal et al. 2019).
Properties of the HD 3167 system useful to our analysis were set to the values derived by Christiansen et al. 2017 (hereafter, C17). Precise values for the ephemeris, however, are critical to the analysis of RM signals (e.g., Casasayas-Barris et al. 2021); C17's ephemeris yield 1σ uncertainties of 19 min for HD 3167b and 16 min for HD 3167c at the time of our observations. As is common for close-in planets, these uncertainties are dominated by the error propagated from the period, which increases linearly with the number of orbital revolutions. To address this issue, we exploited mid-transit times derived from recent transit observations of HD 3167b with the CHEOPS satellite (T 0 = 2459141.90513`0 .00048 0.00050 ; Bourrier et al. in prep.) and of HD 3167c with the Hubble Space Telescope (Guilluy et al. 2021). The large number of planetary transits elapsed between the mid-transit time T C17 0 derived by C17 and those recent T new 0 measurements allow for a substantial Article number, page 3 of 22 A&A proofs: manuscript no. article_RM_HD3167 refinement of the ephemeris as: where P indicates the orbital period, σ the uncertainty on the considered parameter (C17 and new parameters are assumed to be uncorrelated), and n " n new`nC17 (with n C17 and n new the number of orbital revolutions elapsed between a given RM observation and the C17 and new mid-transit times, respectively). We derive P b = 0.9596550˘6.0ˆ10´7 d and P c = 29.846449˘1.9ˆ10´5 d and we decrease to below 1 min the uncertainties on the mid-transit time at the epochs of the HARPS-N (T RM 0 = 2457663.59662˘4.4ˆ10´4 BJD) and ESPRESSO (T RM 0 = 2458534.44470˘6.7ˆ10´4 BJD) observations.

Transit observations
We exploited two datasets on the HD 3167 system, one obtained with ESPRESSO during the transit of HD 3167 b and the other obtained with HARPS-N during the transit of HD 3167 c. Observing conditions for the two transits are given in Table 1. Observations outside of the transit are essential for analyses of the RM effect: to measure the stable stellar flux, to determine the precise RV zero point at the epoch of the observations, and to search for variations due to the Earth's atmosphere and the instrument, which can then be corrected over the entire data set. About half of the exposures were obtained during transit for both HD 3167 b and HD 3167 c, with out-of-transit exposures obtained both before and after the transits (e.g., Fig. 1).
ESPRESSO ) is a fiber-fed, ultrastabilized high-resolution echelle spectrograph installed at the Very Large Telescope (VLT) at ESO's Paranal site. The light is dispersed on 85 spectral orders (each order being covered by two independent slices on the detectors) from 380 to 788 nm. HD 3167 was observed in the frame of Guaranteed Time Observations (Prog: 1104.C-0350(D)) using one of the four 8 m Unit Telescopes (UT3) of the VLT. A constant exposure time was maintained throughout the night, yielding 39 exposures with an average signal-to-noise ratio (S/N) per pixel of 50 at 550 nm (Table 1). We used fiber A to observe the target and fiber B to monitor the sky. The sky signature was not corrected for, as no detrimental contamination was identified. Spectra were extracted from the detector frames, corrected for, and calibrated by version 2.2.8 of the Data Reduction Software (DRS) pipeline. They were then passed through weighted cross-correlation (Baranne et al. 1996;Pepe et al. 2002) with a custom-built numerical mask (Sect. 3.2), yielding cross-correlation functions (CCFs) with a step corresponding to the ESPRESSO pixel width (500 m s´1). Observations were obtained using the HR21 mode (R"140 000), which yields an instrumental resolution of "2.1 km s´1.
HARPS-N (Cosentino et al. 2012) is a fiber-fed, crossdispersed, environmentally stabilized echelle spectrograph installed at the 3.58 m Telescopio Nazionale Galileo (TNG, La Palma, Spain). The light is dispersed on 69 spectral orders from 383 to 690 nm. HD 3167 was observed in the frame of a Director's Discretionary Time program (PI: Hebrard), which was used to study the RM effect of HD 3167 c (Dalal et al. 2019). Here we use these archival data to validate our new technique and refine the spinorbit angle of HD 3167 c. Following Dalal et al. (2019), we discarded the second night of the program obtained in poor weather conditions and with little out-of-transit baseline. Our analysis was performed on the 35 first exposures in the first night, obtained with an exposure time of 900 s and an average S/N of 87 at 550 nm ( Table 1). The last exposure was discarded because of its low exposure time (427 s), high airmass (4.2), and low resulting S/N (33). The remaining four last exposures still have airmass larger than 2, but their S/N remains high, they only contribute to the reference spectrum for the star (Sect. 4.2), and they do not show any spurious features (Fig. 9). Fiber A was set on the star and fiber B on a thorium-argon lamp to ensure a higher RV precision. The data were extracted and CCFs calculated using version 2.2.8 of the ESPRESSO pipeline adapted to HARPS-N , using a custom-built mask (Sect. 3.2) and a step corresponding to the pixel width (820 m s´1). HARPS-N provides a resolving power R"115000, corresponding to an instrumental resolution of "2.6 km s´1.
The variability of the extinction caused by Earth's atmosphere modifies the spectral flux balance of the stellar spectra in a different way throughout each night. In turn, this so-called color effect impacts the shape of the CCFs and can bias the measurements derived from transit time series (e.g., Bourrier & Hébrard 2014, Bourrier et al. 2018b, Wehbe et al. 2020. To avoid these biases, the flux balance of HD 3167 spectra was reset to a stellar spectrum template representative of its spectral type by the ESPRESSO and HARPS-N DRS.

Improved CCF masks
CCFs calculated from the observed spectra correspond to the stellar lines averaged over the disk of the star (hereafter CCF DI , for disk-integrated CCFs). Each CCF DI is automatically fitted by the DRS with a Gaussian model, yielding time series of their contrast, full-width at half-maximum (FWHM), and radial velocity centroid (RV). After subtracting a Keplerian model calculated with the properties of HD 3167 b, c, and d from C17 from the RVs, we assessed the quality and stability of the measurements via these three properties (Fig. 1). In-transit CCF DI were excluded from this analysis since they show anomalous deviations due to the planets occulting local regions of the stellar disk. For both HARPS-N and ESPRESSO, the dispersion of the out-of-transit CCF DI properties was found to be larger than their averaged photon-noise uncertainty, suggesting the presence of systematic variations (Table 2). Order-byorder and line-by-line analyses (Dumusque 2018) revealed such increased dispersions in the bluest orders and in specific orders over the green band. Short wavelengths exhibit many blended stellar lines, which are more sensitive to activity (Anglada-Escudé & Butler 2012). At longer wavelengths, several orders contain stellar lines with deep narrow cores and wide pressure-broadened wings, such as the sodium doublet, whose profiles may not be well captured by a Gaussian model, thus decreasing the stability of the CCF DI they contribute to. We hypothesized that the current DRS masks were amplifying the contribution of these variable lines to the dispersion of the CCF DI properties by overweighting them. Indeed, the weight of a mask line is defined as the square of its depth, which corresponds to the line RV precision assuming that all lines have the same FWHM (Lovis & Fischer 2010). In the case of an inaccurate determination of the continuum level for blended lines or lines with broad wings, the weight in the mask will not correctly reflect the RV precision of the lines.
We therefore developed new masks with weights more representative of the photonic error on the line positions, and with a line selection optimized for the present observations (as in Bourrier et al. 2020). Out-of-transit spectra were first normalized with the RASSINE code (Cretignier et al. 2020b) and aligned to a common rest frame using their individual RV measurements. Master spectra were then obtained by taking the median of the spectra time series, for the HARPS-N and ESPRESSO data independently. Because upper envelop algorithms such as the "Jarvis March" used in RASSINE can yield inaccurate continuum levels, the master continuum was corrected for using a stellar template from the POLLUX database (Palacios et al. 2020), selected for a star with similar effective temperature and gravity as HD 3167. Such an absolute correction of the continuum level is described in the appendix of Cretignier et al. 2021 and can be understood as a matching of the master spectra and stellar template upper envelopes. More accurate line depths can be calculated with this corrected continuum, and we rejected shallow lines with contrast smaller than 10%.
Line positions, identified by local minima surrounded by two local maxima, were fitted with a second order polynomial function of the core of spectral lines (Cretignier et al. 2020a). Weights on individual lines were defined using the formula from Bouchy et al. (2001), integrating the line profiles between˘5 km s´1 around the fitted line center. We used a standard telluric template (Smette et al. 2015, Cretignier et al. 2020a) to identify and reject stellar lines whose relative depth with contaminant telluric lines is larger than 2%. Optimized line selection indeed requires a trade-off (sometimes referred to as a noise-bias trade-off) between rejecting lines affected by systematics and keeping as many lines as possible to increase the S/N. Telluric contamination was computed for the barycentric Earth radial velocity (BERV) span of the observations, 3 km s´1 to account for the width of the tellurics lines. This selection may thus slightly depend on the observing epoch, but this approach allows excluding the smallest possible range of telluric lines contaminating the data on the observed nights, rather than a wider range accounting for all possible BERVs, as in generic masks.
A second cutoff was applied to reject lines unsuitable for the determination of precise RV measurements. We first measured the line-by-line RVs in each exposure (Dumusque 2018), and then calculated the weighted rms for each line over the out-of-transit time series. Lines presenting a rms larger than 50 m s´1 (for HARPS-N) and 80 m s´1 (for ESPRESSO) were rejected, this cutoff again being determined as the best noise-bias trade-off. The larger cutoff for the ESPRESSO data is linked to the larger dispersion of the line-by-line RVs, itself related to the lower S/N of this dataset ( Table 2). The underlying assumption behind the cutoff is that out-of-transit residual RVs are expected to be constant on average, the observed dispersion arising mainly from lines affected by stellar activity. This second line selection can thus also differ between instruments on different telescopes, as their sensitivity to bias depend on the photon noise level. We note that using the rms as a metric does not allow us to know whether rejected lines are affected by systematics or are intrinsically less precise, but the gain brought by rejecting lines affected by systematics is higher than the loss of lines with a high level of white noise. Figure 2 presents a comparison of the distributions of weights between the generic and custom ESPRESSO masks. The generic mask contains a large number of shallow lines that even when combined together bring more noise than information to the CCF DI ; thus they are excluded from the custom mask. We can also see how the generic mask overweighs a substantial fraction of stellar lines, particularly toward shorter wavelengths. Besides an improved line weighing, the RV dispersion cutoff we applied naturally excludes both the low-precision and strongly active lines at the extreme blue of the spectrum. Together, the tellurics and RV dispersion cutoffs decrease the total weight of the line selection by an average of 20%, which is equivalent to an increase in the RV uncertainties by 12%. However, the resulting reduction in bias considerably improves the dispersion of the CCF DI constrast, FWHM, and RVs (Fig. 1). The custom masks yield shallower but substantially narrower CCF DI (Table 2) because they do not include many blended lines in the blue part of the spectrum, which artificially broaden the FWHM of the CCF. The custom masks further yield much more stable FWHM measurements for both the HARPS-N and ESPRESSO data (possibly because it is a known activity tracer and the new weights reduce the contribution of active lines) and a remarkable stabilisation of the contrast in the ESPRESSO data. The dispersion of all properties is closer to their photonic error, especially for the contrast and RV residuals, with a substantial reduction of the dispersion on the RV residuals for both HARPS-N (by 34%) and ESPRESSO (by 25%).

Motivation
In recent years, a new method has proven more efficient at analyzing the Rossiter-McLaughlin effect of exoplanets at higher precision than the classical RM technique (e.g., Ohta et al. 2005, Giménez 2006, Hirano et al. 2011, Boué et al. 2013. The reloaded RM (RRM) technique (Cegla et al. 2016a, see also Bourrier et al. 2017Bourrier et al. , 2018bBourrier et al. , 2020Ehrenreich et al. 2020;Allart et al. 2020;Kunovac Hodžić et al. 2021) isolates the local CCFs (heareafter CCF loc ) from the regions of the stellar photosphere that are occulted by a planet during its transit. The RVs of the planet-occulted regions can then be derived as the centroids of the average stellar line contained in each CCF loc . Fitting those RVs with a realistic model (Cegla et al. 2016a;Bourrier et al. 2020) yields constraints on the stellar surface velocity field and since they directly trace the path of the planet across the star, on the orbital architecture of the system as well.
In its current implementation, however, this approach relies on the possibility to detect and fit the average stellar line in an individual CCF loc with a simple χ 2 minimization. Often, this is not possible in low-flux CCF loc at the limbs of the star, which are darker at optical wavelengths and only partially occulted by the planet (e.g., Bourrier et al. 2018b). More generally, a RRM analysis cannot be carried   (cm s´1) 109 72 111 83 Notes: σ indicates standard deviations of the contrast, FWHM, and Keplerian RV residuals with respect to their mean out-of-transit value (ăxą). The σ rel dispersions have been normalized by this mean to allow for a direct comparison between masks. e x indicate the mean out-of-transit error on x.
out if the stellar line is too shallow, the host star too faint, or the transit depth too small. Up to now, this caveat has stood in the way of determining the orbital architecture of the smallest planets transiting all but the brightest stars (Kunovac Hodžić et al. 2021).
To address these limitations, we present a new approach to the analysis of the RM effect, which we call the RM effect Revolutions (RMR) method. The first step (Sect. 4.2), up to the extraction of the planet-occulted CCFs, is the same as in the traditionnal RRM approach. The second step (Sect. 4.3) improves on the analysis of the average stellar line in individual CCF loc by using MCMC exploration. The third step (Sect. 4.4) exploits the full temporal and spectral information contained in the CCF loc time series to derive the stellar surface and orbital architecture properties with a higher precision.

Step 1: Extraction of planet-occulted CCFs
The CCF DI reduced by the DRS are first aligned by correcting for the Keplerian motion of the star induced by all planets in the system. The relative flux level of the CCF DI is lost due to extinction by Earth atmosphere (Sect. 3.1). The continuum of CCF DI is thus scaled to reflect the planetary disk absorption during the studied transit, using light curves computed with the batman package (Kreidberg 2015). The Keplerian and transit models are calculated using properties from C17 and our refined ephemeris. The CCF DI outside of the transit are co-added to build a master-out CCF DI representative of the unocculted star. By construction, the masks built for HD 3167 are defined in the star rest frame (Sect. 3.2), but we corrected all CCF DI for a small residual offset of 37 m s´1 measured via a fit to the masters-out of each dataset. Residual CCF loc are obtained by subtracting the scaled CCF DI from their corresponding master-out. We note that the new versions of the HARPS-N and ESPRESSO DRS propagate errors from the raw spectra to the CCF DI , allowing us to further propagate the errors on the CCF loc and removing the need to use continuum Values are colored as a function of their normalized flux and plotted as a function of RV in the stellar rest frame (in abscissa) and orbital phase (in ordinate). HD 3167c is large enough that the stellar line can be distinguished in individual exposures, while the small transit depth of HD 3167b makes the line undetectable in most exposures even with the VLT/ESPRESSO data. The green solid lines show the best-fit models for the stellar surface RVs derived with the joint RMR fit to both data sets. dispersion as estimate for flux uncertainties (e.g., Cegla et al. 2016a, Bourrier et al. 2017. Furthermore, rather than working with CCF loc as in previous RRM studies, we reset them to a common flux level by dividing their continuum by the flux scaling applied to the CCF DI . These "intrinsic" CCF intr , shown in Fig. 3, are thus independent of planetary occultation and of the effect of limb-darkening on the continuum flux. They only trace variations in the local stellar line profiles, allowing for a more direct comparison.

Step 2: Analysis of individual exposures
In the traditional RRM approach, each CCF loc is fitted individually with a stellar line model via Levenberg-Marquardt least-squares minimization. A single criterion (e.g., Bourrier et al. 2018b) is used to determine whether the local stellar line is considered detected and whether the best-fit properties can subsequently be used for analysis. This approach has several limitations; for example, the detectability of the stellar line and its exploitation depends on an arbitrary threshold calculated using the fit properties. If the data have systematic errors, the χ 2 minimization makes it difficult to check that the model was fitted to the stellar line rather than to a spurious feature. If the properties of the stellar line model do not have Gaussian uncertainty distributions, the χ 2 minimization can then bias their best-fit estimates and errors.
In the RMR approach, we improve on the analysis of individual planet-occulted CCFs by using a Bayesian approach. A model of the stellar line is fitted to the CCF intr in each exposure, sampling the posterior distributions of the model parameters using emcee MCMC (Foreman-Mackey et al. 2013) and allowing for a variety of priors. In the present study we use a simple Gaussian model with jump parameters set to the contrast, FWHM, and RV centroid, and we chose uninformative priors whose main role is to prevent the MCMC walkers from diverging. A uniform prior distribution between [-5 , 5] km s´1 is set on the RV centroid, considering that stellar surface RVs cannot be larger than the maximum projected stellar rotational velocity (about 2.5 km s´1 for HD 3167, as derived from RM analysis or stellar spectroscopy; Dalal et al. 2019). The flux in the core of stellar absorption lines is bound by zero and the local continuum (corresponding to a contrast range of [0 , 1]), but we set a uniform prior distribution on the line contrast over a large enough range ([-2 , 2]) that we can assess the impact of noise and potentially noisy features on the line detection. A uniform prior distribution between [0 , 20] km s´1 is set on the FWHM, about three times larger than the average FWHM of the local stellar line (see below). Walkers are initialized randomly over the prior ranges to ensure a thorough exploration of the parameter space. In the present analysis, we empirically chose (checking for convergence of the chains) to use 100 walkers running for 2000 steps, with a burn-in phase of 500 steps. We tested the model against the data over [-50 , 50] km s´1 in the star rest frame. The median of the posterior probability distributions (PDF) are considered as the derived values for the model parameters. To better account for multimodal and asymmetric PDFs we define 1σ uncertainty ranges using highest density intervals (HDIs), which contains 68.3% of the PDF mass such that no point outside the interval has a higher density than any point within it.

Single-exposure analysis for HD 3167c
Results from the MCMC exploration of the HD 3167c transit are displayed in Fig. 4. Most exposures show well-defined PDFs, narrow enough that they constrain the model parameters within the prior ranges. In such cases, the local stellar line can be detected in individual CCF intr via χ 2 minimization, as in the RRM study performed by Dalal et al. (2019). Nonetheless, several exposures show PDFs with asymmetrical profiles, highlighting the interest of a Bayesian approach and the combination of median + HDI to interpret the fit results and prevent biases. Furthermore, this approach reveals that the last three exposures have a contrast and FWHM PDFs that are consistent with null values, which shows that the flux noise in these CCF intr is too large for the stellar line to be detected. The S/N of the spectra decreased over the transit as airmass increased, dropping below the average S/N in the night for the last three in-transit exposures. Hereafter, we thus exclude these exposures from all analysis, as they bring more noise than constraints to the interpretation of the transit time series.
Local RVs from individual exposures have small enough uncertainties that they clearly show the signature of the RM effect from HD 3167c, which allowed Dalal et al. (2019) to fit them using the traditional RRM approach (Fig. 5). The local RV series is negative and slightly decreases with phase, showing that HD 3167c crosses the stellar hemisphere rotating toward us, and moves farther away from the projected stellar spin axis over the course of its transit. This configuration yields the redshifted anomaly visible in the RV centroids of the CCF DI (Fig. 1). An analysis of the derived properties in Fig. 5 suggests possible center-to-limb variations in local contrast (see e.g., Cavallini et al. 1985, Cavallini et al. 1985Löhner-Böttcher et al. 2019 for the Sun), while uncertainties on the local FWHM are too large to draw any conclusions on center-to-limb variations.

Single-exposure analysis for HD 3167b
Despite the use of VLT/ESPRESSO to observe the transit of HD 3167b, the planet is so small that the classical RRM approach fails at detecting and constraining the local stellar line in individual exposures. This can be understood when looking at the results of the MCMC exploration (Fig. 6). Most exposures show ill-defined PDFs, or broad with respect to the prior ranges, especially for the local RVs. Even the few exposures near the center of the transit, which have the largest S/N on the night and show better-defined PDFs, still poorly constrain the local RVs (Fig. 7). Nonetheless, information on the local stellar line can still be present even when it is hidden by noise. Binning the consecutive exposures could unveil the local line by increasing the S/N, but at the cost of blurring the signal. The motivation behind the RMR approach, as detailed in the next section, is to circumvent our inability to measure the stellar line properties in individual exposures by fitting all the exposures together, without modifications.
The MCMC fit to individual exposures, along with the analysis of the associated PDFs, remains a useful step in assessing the quality of each CCF intr and determining which model best fits the local line profile variations. Here, the last exposure does an especially poor job in constraining the model parameters (Fig. 6), likely because it is only partly in-transit, and thus it was excluded from further analysis.

Step 3: Global transit analysis
The main novelty of the RMR approach lies in fitting all CCF intr together with a joint model, rather than fitting them individually. The stellar line profile can be described by a variety of analytical models (e.g., Gaussian; double Gaussian, Bourrier et al. 2018bBourrier et al. , 2020. The properties describing the line profile can be kept constant in all exposures or be set by parametric models as a function of a given coordinate parameter. The RV centroids of the theoretical lines are set by the surface RV model described in Cegla et al. (2016a) and Bourrier et al. (2017). It accounts for the projected stellar rotational velocity field (plus convec- tive blueshift and differential rotation when relevant) and its blur over the region occulted by the planet for a given exposure. The time series of theoretical stellar lines is convolved with the instrumental response and then fitted to all CCF intr together, using emcee MCMC to sample the PDFs of the model parameters. Derived values and uncertainties for the model parameters are again defined using medians and HDIs. We note that modeling the stellar line before convolution allows for a more direct comparison between different instruments.
The main interest of the RMR technique is that it boosts the S/N of the occulted stellar line, at first order by the number of in-transit exposures, thus increasing the possibility of constraining the planet path even when the line is not detectable in individual exposures. The RMR approach further exploits the full extent of information contained in the data, as it fits the stellar line profiles rather than just their centroids (as in the traditional RRM approach). Below, we first show how we validated the RMR approach by comparing its results on HD 3167 c with those from the traditional RRM approach, and we then demonstrate its importance for small planets via its application to HD 3167 b.
We used a Gaussian profile to model the local stellar lines occulted by both planets. Given that the analysis of individual exposures in Step 2 (Sect. 4.3) hinted at possible changes in the line shapes along the transit chords, we explored polynomial variations of the contrast or FWHM with µ = cos θ (where θ is the center-to-limb angle) and used the Bayesian information criterion (BIC, Schwarz 1978, Liddle 2007 to determine the best polynomial order (see Kass & Raftery 1995 for the details of the BIC comparison). The data are not precise enough to search for second-order RV Article number, page 9 of 22 A&A proofs: manuscript no. article_RM_HD3167 effects caused by the stellar surface motion, and we thus assumed solid-body rotation for the star. Jump parameters for the MCMC are the polynomial coefficients for the contrast and FWHM, the sky-projected obliquity λ and stellar rotational velocity, v sin i˚. Uniform priors are set on the model parameters, which are uninformative for the contrast and FWHM coefficients (i.e., over a much larger range than their range of variations), between [0 -10] km s´1 for v sin i˚, and over its definition range ([-180 , 180]˝) for λ. The model was fitted to the CCF intr over the range of [-50 , 50] km s´1 in the star rest frame.

Global analysis for HD 3167c
The PDFs for the best-fit parameters are shown in Fig. 8. They are well-defined with symmetrical profiles, so that HDIs are here equivalent to quantile-based credible intervals. Stellar lines along the transit chord of HD 3167c are best reproduced with a linear variation of the contrast with µ and a constant FWHM. This model improves upon the BIC by ten points compared to a constant contrast, and by six points compared to a linear variation of the FWHM with µ. Contrast and FWHM of the best-fit line model convolved with the HARPS-N response are shown in Fig. 5. They are in agreement with the properties derived from individual CCF intr , showing that the RMR approach captures the line variations along the transit chord well. This is also visible in the residual map between the CCF intr and their best-fit models (Fig. 9). The best-fit FWHM is 5.46˘0.45 km s´1, corresponding to 6.05˘0.41 km s´1 after accounting for instrumental convolution. The model contrast increases from the center of the transit chord toward the limbs. Given the highly misaligned orbit of HD 3167c this could imply a deepening of the lines with increasing latitude, textiti.e. from stellar equator to stellar poles (e.g., due to the stellar magnetic field, see a similar case in Bourrier et al. 2017). The exact dependence of this contrast variation with stellar latitude is, however, linked to the unkown inclination of the stellar spin axis with respect to the observer.
We find that λ and v sin i˚are not correlated with the contrast and FWHM of the fitted lines, so that their bestfit values do not depend on the line model. We derive λ = -102.2`7 .5 8.7˝a nd v sin i˚= 2.37˘0.42 km s´1. These results are consistent within 1σ with the values derived by Dalal et al. (2019), using three different techniques, including the traditional RRM approach, thus validating the RMR approach. In RM analysis λ and v sin i˚can be degenerate with the scaled semi-major axis a{R˚and orbital inclination i p . We allowed these parameters to vary (see Bourrier et al. 2020), using values from C17 as priors. The resulting fit, however, increases the BIC by 15 points and yields PDFs for a{R˚and i p that match their priors, showing that these parameters cannot be better constrained by the present data set.

Global analysis for HD 3167b: Detection
We exploited results from the HD 3167c fit to refine the prior on v sin i˚, using a normal distribution with center 2.37 km s´1 and width 0.42 km s´1. Residuals between the CCF intr and the best-fit line models are shown in Fig. 11. The corresponding PDFs for the model parameters are shown in Fig. 10. Apart from λ, these are well defined with mostly symmetrical profiles. Stellar lines along the transit chord are best reproduced with a linear variation in contrast with µ and a constant FWHM (Fig. 7), improving the BIC by six points compared to a constant contrast or to a linear variation in FWHM with µ. The best-fit FWHM is 5.00`0 .78 0.91 km s´1 (corresponding to 5.38`0 .72 0.85 km s´1 for the line observed with ESPRESSO), in good agreement with the value derived from the HARPS-N transit even though the ESPRESSO and HARPS-N CCF masks are not the same. In contrast to HD 3167c, the derived line contrast decreases from the center of HD 3167b transit chord toward the limbs. This could indicate that the contrast variations observed for the two planets do not depend primarily on µ, but on stellar latitude instead.
Indeed, one of the main results from our study is that HD 3167b and HD 3167c transit along nearly perpendicular tracks across the stellar surface, and thus occult very different stellar latitudes. The RMR approach allows us to detect the RM signal of HD 3167 b, as shown by the PDF of the model obliquity, with λ = 5`2 8 38˝ ( Fig. 10). This is the first clear detection of a spectroscopic RM signal for a planet smaller than 2 R C (the second smallest being Pi  Men c, 2.1 R C , Kunovac Hodžić et al. 2021), showing the potential of the RMR approach to unveil the orbital architectures of the smallest exoplanets. As with HD 3167c, we note that λ is not correlated with the line shape properties (Fig. 10) and, accordingly, its value does not depend on the contrast and FWHM models.
The ESPRESSO data do not constrain much the stellar rotational velocity compared to the prior from the HARPS-N fit, with v sin i˚= 2.16`0 .41 0.39 km s´1. They also do not constrain the scaled semi-major axis a{R˚and orbital inclination i p of HD 3167b better than the priors from C17, as letting these parameters free to vary increased the BIC by 16 points with no change in their distributions.

Global analysis for HD 3167b: Validation
We carried out several tests to assess the validity of the HD 3167b RM detection.
First, we tested the null hypothesis, namely, the absence of a detectable stellar line from the regions occulted by HD 3167b. The corresponding model, a constant flux in all CCF intr , yields a BIC larger than the best-fit model by 80 points, unambiguously displaying the presence of a signal in the data.
Then, to ensure that the detected signal is not a spurious feature arising from noise, we performed the same fits as in Sect. 4.4.2, but in the blue and in the red continuum of the CCF intr . The nominal fit was carried out over the range of [-50 ; 50] km s´1 in the star rest frame. For the "blue" test, we artificially set the star rest frame at -160 km s´1, that is, we shifted all CCF intr by 160 km s´1 to perform the fit over the original range [-210 ; -110] km s´1. For the "red" test, all CCF intr were shifted by -160 km s´1 , so that the fit was performed over [110 ; 210] km s´1. We show in Appendix A, the correlation diagrams for these two fits. They converge toward model stellar lines with contrast and A&A proofs: manuscript no. article_RM_HD3167 We thus emphasize that the solid black lines models were not fitted to the properties loosely constrained by each individual exposure, but were derived from the global fit to the CCFintr time series.
FWHM consistent with null values, and a roughly uniform PDF for v sin i˚and λ. This is the expected result when fitting a region of constant flux with photon noise and no spurious features, giving us confidence that there are no systematics in the CCF intr that could mimic the detected signal.
Finally, to further ensure that the detected signal is not a spurious feature we performed a fit on the residual CCF loc time series obtained after the transit, which covers a phase range of same duration as the fitted transit window. As can be seen in Fig. A.3, the fit again converges toward a model stellar line with null contrast and FWHM, as expected from a constant flux with white noise.

Multi-planet fit
While we cannot exclude this possibility, we found no evidence for an instrumental origin for the contrast variations measured for HD 3167b and HD 3167c. The transit light curves used to scale the CCF DI (Sect. 4.2) have an impact on the derived local contrasts, but varying the transit depth within its uncertainties does not suppress the contrast variations. The different trends in local line contrast for the two planets can, however, be reconciled within the frame of an intriguing scenario. Due to the different obliquities of their orbits, HD 3167b and HD 3167c effectively occult regions of very different latitudes across the stellar surface. While no meridional variations in line contrast are observed over the solar surface, HD 3167 is a K-type star with a super-Earth orbiting within its corona, at only four times the stellar radius. The measured line contrast variations can be explained by a common trend if we assume that the stellar surface properties change with latitude, either due to the star itself or because of interactions with its close-in companion (e.g., Bourrier et al. 2018a). In the framework of this scenario, we applied a joint RMR fit to the two transit datasets, using a common intrinsic line model. The line positions are set by the theoretical transit chords of each planet, controlled by their respective sky-projected obliquities, and by the theoretical surface velocity field of the star, now controlled independently by the stellar equatorial velocity, v eq , and the stellar inclination, i ‹ . We fit a common FWHM for the intrinsic Gaussian line in the ESPRESSO and HARPS-N data, since the independent fits to each dataset yielded widths consistent within their uncertainties (Sect. 4.4). Inspired by differential rotation models, we assume that the intrinsic line contrast varies as a polynomial function of sinpθ lat q 2 , where θ lat is the latitude relative to the stellar equator (we also assume that variations of the contrast with µ are negligible compared to its latitudinal variations). Since this model depends on absolute latitude, the stellar inclination is fitted within [0 , 180]˝(where i ‹ is the angle between the line-of-sight and stellar spin axis, see Cegla et al. 2016a for details on the coordinate system). We set uniform priors between [0 -10] km s´1 for v eq , and between [-1,1] for cospi ‹ q (used as jump parameter for the stellar inclination). The PDFs for the model parameters, sampled as in Sect. 4.4,are shown in Fig. 12. We find that the different line contrast variations measured for both planets can be explained by a common quadratic function of sinpθ lat q 2 (Fig. 13), providing a better BIC than a linear function by eight points, and than a third-order polynomial by five points. In this scenario, the stellar spin has to be slightly inclined (i ‹ = 111.6`3 .1 3.3˝) , so that HD 3167b transits closer to the stellar pole (Fig. 15). Breaking the degeneracy on v sin i˚allows us to derive the stellar equatorial rotational velocity v eq = 2.65`0 .47 0.42 km s´1. Accounting for the uncertainties on R ‹ = 0.850˘0.020 R d (derived from a spectroscopic analysis of the ESPRESSO data, following the approach from da Silva et al. 2006), this result further yields the equatorial rotation period of the star as P eq = 16.2`2 .8 2.7 d. As with the fits to each dataset, the obliquities are not correlated with the line contrast, giving us confidence in the derived λ b = -6.6`6 .6 7.9˝a nd λ c = -108.9`5 .4 5.5˝. These values are consistent with those derived from each dataset, but their precision was improved by fitting the two datasets together. We note that local line contrasts from some of the regions occulted by HD 3167b are outside [0 ; 1], although the largest deviation remains within 2σ from this range (Fig. 13). This could be due to uncertainties in the model light curve used to scale in-transit CCF DI , as a transit depth larger than its actual value would increase the local line contrasts artificially. We checked, however, that varying the transit depth of both planets within their 1σ uncertainties has no impact on the derived obliquities. We show in Fig. 14 the average of the CCF intr measured along each transit chord, after they were aligned using the bestfit surface RV models. Gaussian fits to these average lines yield contrasts of about 52% for HD 3167c and 86% for HD 3167b, illustrating the fact that the two planets occult stellar regions with different properties.

Orbital architecture of the HD 3167 system in context
The main architectural properties we derived on the HD 3167 system are reported in Table 3. The observables in our study, namely, the RV position and contrast of the local stellar lines occulted by the planets, leave a degeneracy between two possible orbital configurations shown in Fig. 15: A (i b , λ b , i c , λ c , i ‹ ) and B (π-i b , -λ b , π-i c , -λ c , π-i ‹ ). These two configurations, however, yield the same values for the 3D spin-orbit angles of the planets ψ = arccos(sin i ‹ cos λ sin i p + cos i ‹ cos i p ), and for the mutual inclination between their orbital planes i mut = arccos(cos i b cos i c + cos Ω sin i b sin i c ). Here, Ω is the difference between the longitudes of the planets' ascending nodes, corresponding to λ c -λ b (A) or -λ c +λ b (B).
We combined the probability distributions of i p , λ, and i˚to calculate the PDFs of ψ and i mut (Fig. 16). We derive ψ b = 29.5`7 8.0˝. Thus, HD 3167b is close to orbiting within the stellar equatorial plane, while HD 3167c is on a near-polar orbit around the star, and the orbital planes of HD 3167b and HD 3167c are nearly perpendicular.
Dynamical analysis by C17 and Dalal et al. (2019) showed that HD 3167c and d are not coplanar but have low mutual inclinations, likely between 2.3 and 21˝. The high spin-orbit angle of HD 3167c thus implies that both planets are on near-polar orbits around the star. Dalal et al. (2019) further predicted that the stellar spin could not be aligned with the angular momentum of these two planets, which is confirmed by our measurement of ψ for HD 3167c. These authors concluded that if they formed on aligned orbits, HD 3167c and HD 3167d are unlikely to have reached such misaligned orbits by interacting with the star. They proposed instead that the orbital planes of these planets were misaligned via gravitational interactions with an outer inclined companion (Boué & Fabrycky 2014). This companion should have been massive enough and close enough to HD 3167c and HD 3167d that their dynamical evolution was not dominated by gravitational interactions with the star, but not too massive and close that it entered a Kozai resonance with HD 3167c and destabilized the outer planets. In the intermediate regime, the dynamical evolution of HD 3167c and HD 3167d would have been controlled by gravitational interactions with the companion, which gently excited their orbits to high inclinations with respect to the star over secular timescales.
At such short orbital distance (only four times the stellar radius), the present-day orbit of HD 3167b is governed by interactions with the star and the dynamical stability of the system does not depend on this planet's inclination (Dalal et al. 2019). We propose that HD 3167b formed with a low enough mass and far enough away from planets c and d that their dynamical evolution was always decoupled, and that HD 3167b was not influenced by the outer companion. While HD 3167b is close enough to the star that it might have reached its present orbit due to tidal interactions, the three planets likely initially migrated via disk-driven interactions, remaining close to the stellar equatorial plane. As the star spinned down HD 3167c and HD 3167d became more sensitive to the gravitational influence of the outer companion (Dalal et al. 2019) and were progressively misaligned, while HD 3167b remained strongly coupled to the star and kept its primordial alignment. We note that we consider HD 3167b as aligned following the limit of Ψ "30˝commonly used in the exoplanet literature (e.g., Triaud et al. 2010;Albrecht et al. 2012). Albrecht et al. (2021) recently confirmed that exoplanet systems do not span the full range of obliquities, but show a preference for nearly-perpendicular ("misaligned"; Ψ = 80-125˝) or low-obliquity ("aligned"; Ψ ă " 35˝) orbits. deg It is interesting to compare HD 3167 with other multiplanet systems with known architectures. HD 3167 provides an interesting counterpoint to the multi-planet systems Kepler-56 and K2-290 A. Kepler-56 hosts two inner giant planets on nearly coplanar orbits, but that are misaligned with their star (Huber et al. 2013). Dynamical simulations by Gratia & Fabrycky (2017) show that late scattering with the third more massive outer planet in the system (Otor et al. 2016) could have misaligned the inner system with the star. K2-290 A also hosts a mini-Neptune and a Jupitersize planet on coplanar orbits (Hjorth et al. 2021). Their Rossiter-McLaughlin signatures, analyzed with the classical RV technique, are however consistent with mutually aligned planetary orbits both retrograde with respect to the star. Combined with the presence of a stellar companion on a wide orbit, this finding led Hjorth et al. (2021) to conclude that the most likely scenario for the K2-290 A system is a primordial misalignment of the protoplanetary disk with the primary star, caused by the gravitational torque from the stellar companion. The dynamical evolution of the HD 3167's outer planetary system could be similar to what was proposed for Kepler-129 by Zhang et al. (2021), also in the framework developed by Boué & Fabrycky (2014). The orbits of the two sub-Neptunes around Kepler-129 could have become misaligned as a result of precessing around  the orbital normal of a long-period giant planet discovered in the system. The HD 3167 system -with outer planets that were likely tilted by secular interactions with an outer companion but an inner planet still nearly aligned with the star -offers yet another example of a complex dynamical history. These various results illustrate how the analysis of RM signatures in transiting multi-planet systems allows us to distinguish between various primordial or secular misalignment scenarios.

Conclusion
Rossiter-McLaughlin studies of transiting exoplanets yield precious information about their orbital architecture. Recent improvements in analysis techniques have unveiled the signal of smaller exoplanets (e.g., Kunovac Hodžić et al. 2021) around slower rotators (e.g., Bourrier et al. 2018b). Yet even with next-generation spectrographs such as ESPRESSO, the RM signal of super-Earths and smaller planets remain difficult to reach. Building on the latest technique (Cegla et al. 2016a), we developed the Rossiter-McLaughling effect Revolutions (RMR) approach to exploit even further the information contained in transit time series. This technique relies on fitting together the spectral profiles from all planet-occulted regions with a joint model, with the possibility to combine transit data sets from multiple planets obtained with different instruments. We first validate the RMR approach on HARPS-N observations of the mini-Neptune HD 3167c, confirming its known RM signal and then apply it to ESPRESSO observations of the super-Earth HD 3167b, detecting its RM signal for the first time. We emphasize the critical importance of precise ephemeris for RM studies of close-in planets, which requires re-observing their transit regularly (e.g., Casasayas-Barris et al. 2021). HD 3167 b is the smallest exoplanet (1.7 R Earth ), with a clear detection of a spectroscopic RM signature (see the cases of 55 Cnc e, Bourrier & Hébrard 2014, López-Morales et al. 2014; and of the TRAPPIST-1 system, Hirano et al. 2020). The RMR approach thus opens a new window into the orbital architecture of Earth-size planets, up to now mainly probed via photometric techniques (e.g., Kamiaka et al. 2019).
Curiously, the contrast of the stellar lines occulted by HD 3167b and HD 3167c display different variations from the center to the limb of the star. While we cannot exclude an instrumental origin or a bias due to the photometric Fig. 11. Residual maps during the HD 3167b visit (same details as in Fig. 9).
scaling of the data, it is remarkable that the contrast variations can be explained by a common dependence of the stellar line shape with latitude. In that scenario, a joint RMR fit to the ESPRESSO and HARPS-N datasets allow us to constrain the stellar inclination, thus breaking the sky-projection degeneracy and measuring the 3D spin-orbit angles of HD 3167b (29.5`7 .2 9.4˝) and HD 3167c (107.7`5 .1 4.9 ). HD 3167b is close to orbiting within the stellar equatorial plane, while HD 3167c orbits above the poles of the star, and the orbits of the two planets are nearly perpendicular (mutual inclination 102.3`7 .4 8.0˝) . This unusual orbital architecture likely traces different dynamical pathways for the planets. Building on the results from Dalal et al. (2019), we propose that the evolution of HD 3167b has always been dominated by its interactions with the star, so that its orbit retained its primordial alignment, while the evolution of HD 3167c and d was influenced by secular interactions with an outer companion, so that their orbits were gently excited to their present high inclination. Follow-up observations of the system via velocimetry, transits, and direct imaging are needed to refine the properties of the planets, in particular their eccentricity, and to confirm the presence of the companion hinted in RV observations of the system (Dalal et al. 2019). Exoplanets on perpendicular orbits, resulting from gravitational interactions with outer companions, might actually be a more common occurrence than previously considered ). Fig. 13. Local contrast of the stellar lines occulted by HD 3167b (green) and HD 3167c (orange), derived from the fits to individual exposures (Sect. 4.3,Figs. 5,and 7). Solid black lines show the best-fit model for the contrast variations as a function of the squared latitude coordinate sinpθ lat q 2 (ranging from 0 at the stellar equator to 1 at the poles), after convolution by the HARPS-N and ESPRESSO LSFs. We note that the models were derived from the RMR fit to the combined datasets.   15. Projection of HD 3167 in the plane of sky for the best-fit orbital architecture (the top and bottom panels show configurations A and B, respectively; see text). The stellar spin axis is displayed as a black arrow extending from the north pole. The stellar equator is represented as a black line. The stellar disk is colored as a function of its surface RV field. Normals to the orbital planes of HD 3167b and HD 3167c are shown as green and brown arrows, respectively. Thick solid curves with corresponding colors represent the best-fit orbital trajectories. The thin lines surrounding them show orbits obtained for orbital inclination, semi-major axis and sky-projected obliquity values drawn randomly within 1σ from their probability distributions. The star, the planets HD 3167b and HD 3167c (black disks), and their orbits are shown to scale. HD 3167d, which orbits in between HD 3167b and HD 3167c, is not shown because of its unknown inclination.