Open Access
Issue
A&A
Volume 638, June 2020
Article Number A87
Number of page(s) 20
Section Planets and planetary systems
DOI https://doi.org/10.1051/0004-6361/201937316
Published online 18 June 2020

© A. Wyttenbach et al. 2020

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1 Introduction

Observations of exoplanet atmospheres allow us to characterize several of these remote worlds. Indeed, the measurement of planetary spectra are rich in information about these planets’ physical and chemical conditions. Thus, the study of exoplanetary atmospheres are motivated by the possibility that they can provide crucial hints about their origin and evolution.

When observing the primary transit of an exoplanet as a function of the wavelength, we measure a transmission spectrum. The transmission spectrum of an exoplanet bears the signature of its atmospheric constituents. The possibility of detecting atomic or molecular species through transmission spectroscopy was theorized early on (Seager & Sasselov 2000; Brown 2001). Since giant close-in planets (Mayor & Queloz 1995), with hot and extended atmospheres, are easy targets to probe in transmission, early observations of hot Jupiter atmospheres were achieved with sodium detections (Charbonneau et al. 2002; Redfield et al. 2008; Snellen et al. 2008). Recent results suggest that high-resolution (typically λ∕Δλ ~ 105) transmission spectra of hot Jupiters probe their thermospheres, a key region of the upper atmosphere, in order to understand irradiated planets and allow us to achieve precise temperature and wind measurements because the lines (e.g., the Na I D doublet, the K I D doublet, the H I Balmer lines Hα, β, γ, the He I (23S) triplet) are resolved and studies benefit from large signals (e.g., Wyttenbach et al. 2015, 2017; Casasayas-Barris et al. 2017, 2018, 2019; Khalafinejad et al. 2017, 2018; Seidel et al. 2019; Yan & Henning 2018; Cauley et al. 2019; Allart et al. 2018, 2019; Nortmann et al. 2018; Salz et al. 2018; Alonso-Floriano et al. 2019; Keles et al. 2019; Chen et al. 2020). Another landmark achievement obtained in the optical domain is the detection, with cross-correlation, of iron and other metals in the atmosphere of several ultra-hot Jupiters (Hoeijmakers et al. 2018, 2019; Borsa et al. 2019; Bourrier et al. 2020a; Cabot et al. 2020; Gibson et al. 2020; Nugroho et al. 2020; Stangret et al. 2020; Ehrenreich et al. 2020).

Ultra-hot Jupiters are the hottest and the most irradiated hot Jupiters (with a planetary equilibrium temperatureabove 2000–2500 K, Fig. 1). Well-known members of that category include KELT-9 b (Gaudi et al. 2017, the hottest exoplanet discovered so far around a main-sequence star, with Teq ≳ 4000 K, see Table 1), the newly discovered MASCARA-4 b/bRing-1 b (Dorval et al. 2020), as well as targets often observed due to their bright host stars, namely, WASP-189 b, WASP-33 b, MASCARA-1 b, WASP-12 b, WASP-103 b, WASP-18 b, WASP-121 b, KELT-20 b/MASCARA-2 b, HAT-P-7 b, WASP-76 b, etc. Hot exoplanet atmospheres are expected to be near chemical equilibrium from T ≳ 2500 K (Miguel & Kaltenegger 2014; Lothringer et al. 2018; Kitzmann et al. 2018). Under such assumptions, planet atmospheres are expected to be cloud-free on the day-side, to have strong inverted temperature-pressure profiles, and to be nearly barren of molecules, apart from H2, CO, and possible traces of H2O, TiO, and OH. Because of thermal dissociation and ionization, mostly neutral and ionized atomic hydrogen and helium, along with such metals as iron, sodium, magnesium, titanium, etc., are expected to be found (Kitzmann et al. 2018; Lothringer et al. 2018; Lothringer & Barman 2019; Mansfield et al. 2020). In this type of atmosphere, the important source of H bound-free and free–free opacities acts as a dominant continuum opacity (Arcangeli et al. 2018; Parmentier et al. 2018; Bourrier et al. 2020b). These predictions can change due to heat transport and circulation, depending on whether we observe the day or the night side, or the terminator (Parmentier et al. 2018; Bell & Cowan 2018). For example, in emission spectroscopy, these features of ultra-hot Jupiters are mixed up with the H2 O and TiO emission features (that may be absent due to dissociation), which make emission spectra appear as a continuum (see Parmentier et al. 2018; Arcangeli et al. 2018; Mansfield et al. 2018; Kreidberg et al. 2018). In transmission spectroscopy, the probed pressures are lower, thus increasing the presence of atomic species because there is a higher rate of thermal dissociation. These higher rates are produced even without taking disequilibrium processes into account. In particular, the effect of the H opacity increases, playing a more crucial role in the modeled transmission spectra (Kitzmann et al. 2018).

Recent observations of KELT-9 b in transmission spectroscopy bring in key aspects for our understanding of ultra-hot Jupiters. The detections of Fe I, Fe II, Ti II (but not of Ti I; Hoeijmakers et al. 2018) are confirming the predicted presence of neutral and ionized species in ultra-hot Jupiters. In Cauley et al. (2019); Hoeijmakers et al. (2019); Yan et al. (2019); Turner et al. (2020), more neutral and ionized atomic species are reported, such as Na I, Ca I, Ca II, Mg I, Cr II, Sc II, and Y II. In KELT-9 b, ionized species are more abundant than neutral species, which is expected for such high temperatures. However, the absorptions measured for ionized species are not reproduced by the models (Hoeijmakers et al. 2019; Yan et al. 2019). Those measurements probe altitudes between 1.01 and 1.1 planetary radii (1–15 scale-height or 10−3 –10−6 bar), thus, in the lower thermosphere. A complete exploration of species absorbing at optical wavelengths in ultra-hot Jupiter is very constraining in terms of aeronomic conditions. Indeed, the additional detections of Hα by Yan & Henning (2018); Turner et al. (2020), and of Hα, Hβ by Cauley et al. (2019), also bring in important clues about the KELT-9 system because measuring absorption from specific strong opacity lines, such as the sodium Na I D doublet or the hydrogen Balmer lines, allow us to probe higher layers (≳ 1.2 planetary radii) in the thermosphere (Christie et al. 2013; Heng et al. 2015; Huang et al. 2017). So far Hα has only been observed in six exoplanets, four of which are classed as ultra-hot Jupiters (Jensen et al. 2012, 2018; Cauley et al. 2015, 2016, 2017a,b,c, 2019; Barnes et al. 2016; Casasayas-Barris et al. 2018, 2019; Cabot et al. 2020; Chen et al. 2020). The detection of an optically thick Hα absorption in the atmosphere of KELT-9 b, which takes place at up to ~1.6 planetary radii (at ~10 000 K in the thermosphere), provides hints that the hydrogen is filling up the planetary Roche lobe (1.95 Rp) and is escaping from the planet at a rate of ~ 1012 g s−1. This confirms that such planets undergo a process of evaporation (Sing et al. 2019). However, it is important to get a more precise idea of the mass loss rate of such objects to understand the role of the stellar irradiation and its processes of interaction with the planet atmosphere (Fossati et al. 2018; Allan & Vidotto 2019; García Muñoz & Schneider 2019).

This paper reports on the results of the Sensing Planetary Atmospheres with Differential Echelle Spectroscopy (SPADES) program (Bourrier et al. 2017b, 2018; Hoeijmakers et al. 2018, 2019). This program, performed with the HARPS-N spectrograph, is a twin of the Hot Exoplanetary Atmospheres Resolved with Transit Spectroscopy (HEARTS) program based on the use of the HARPS spectrograph (Wyttenbach et al. 2017; Seidel et al. 2019; Bourrier et al. 2020a). Throughout this work, we introduce the python CHESS code (the CHaracterization of Exoplanetary and Stellar Spectra code) and its sub-modules, which allow us to measure, separate, and model the stellar and exoplanet spectra.

The HARPS-N transit observations1 are described in Sect. 2. In Sect. 3, we propose a comprehensive analysis of spectroscopic transit data and present an improvedmethod to extract transmission spectra. Furthermore, we show a stellar pulsation analysis and a Rossiter-McLaughlin analysis, and we present the Balmer lines signatures we found in the transmission spectrum of KELT-9 b. In Sect. 4, we present a new model and retrieve fundamental atmospheric properties from the high-resolution data observations. In Sect. 5, we discuss the constraints we find on the KELT-9b atmospheric temperature, mass loss rate, and local thermodynamic properties. We present our conclusions in Sect. 6.

thumbnail Fig. 1

Distribution of planetary masses as a function of the total power received by the planet from the incident bolometric stellar radiation, expressed in units of power received on Earth at the top of the atmosphere (P = 174 PW = 1.74 × 1024 erg s−1). The markersize scales with the planet density, and the color with the planetary equilibrium temperature. Physical properties of exoplanets were extracted from the Extrasolar Planets Encyclopedia (exoplanet.eu) in December 2019. The lower left corner represents the evaporation desert. KELT-9 (HD 195689), encircled in black, stand out of all other planetary population in term of its received power and equilibrium temperature.

Table 1

Adopted values for the orbital and physical parameters of the KELT-9 system.

Table 2

Log of HARPS-N observations.

2 The HARPS-N dataset

We observed two transits of the ultra-hot Jupiter KELT-9 b around its bright A0V host star with the HARPS-N spectrograph (Cosentino et al. 2012)mounted on the 3.58 m Telescopio Nazionale Galileo (TNG) telescope in the Roque de los Muchachos observatory (La Palma, Spain)2.

We dedicated a full night of continuous observations to each transit of KELT-9 b on 31 July 2017 and 20 July 2018 (see Table 2). One fiber recorded the target flux (fiber A) and the other monitored the sky (fiber B). In total, we registered 95 spectra, each with an exposure time of 600 s. The time series cover each transit in full and cover baselines of ~2–3 h before and ~1.5 h after each transit, used to build accurate out-of-transit master spectra. The typical signal-to-noise ratio (S/N) calculated in the continuum around 400 and 600 nm is between 27–131 and 70–190, respectively (no spectra were discarded in the analysis; however, the spectra collected on the second night are affected by a failure in the ADC, making their blue part much noisier than in the first night). Finally, the CHESS.KING module (KINematic and Geometric properties of the system) handles all the system parameters (see Table 1) needed during the data analysis.

We identified 44 spectra (22 in each night) fully in transit (i.e., entirely exposed between the first and fourth contacts) and we considered the 51 remaining spectra out of transit. The HARPS-N raw frames are automatically reduced with the HARPS-N Data Reduction Software (DRS version 3.5). The DRS pipeline consists of an order-by-order optimal extraction of the two-dimensional echelle spectra. Then, for each exposure, the daily calibration set is used to flat-field, blaze correct, and wavelength-calibrate all the 69 spectral orders. Finally, the orders are merged into a one-dimensional spectrum and resampled (with flux conservation) between 387 and 690 nm with a 0.001 nm wavelength step (the wavelengths are given in the air and the reference frame is the solar system barycenter one, see Fig. 2). Our main analysis is conducted on the final HARPS-N merged spectra that have a spectral resolution of λ∕Δλ = 115 000 or 2.7 km s−1. In the present study, we took advantage of the cross-correlation functions (CCFs) that are computed with the DRS. The CCFs are constructed by correlating each spectral order with a weighted stellar mask (a line list containing the main lines of the spectral type of interest, Pepe et al. 2002). Then, the data are further prepared with the CHESS.ROOK module (spectRoscopic data tOOl Kit). For the KELT-9 observations, we use a customized A0 mask (as in Anderson et al. 2018) containing about 600 lines of typically Fe i, Fe ii, and other metal species. The A0 mask is sensitive to both the stellar and planetary spectra. Since the stellar spectrum is rotationally broadened, we computed the CCFs over a 700 km s−1 window around the systemic velocity and with an oversampling step of 0.082 km s−1 (a tenth of the pixel size). Then, we conducted the CCF analysis on the sum of the CCF of all the 69 orders3 (Fig. 2). The use of a custom A0 mask boost the stellar CCFs contrast by a factor of ~10 in comparison to the standard DRS reduction, allowing us to have a precision similar to the independent reduction presented in Borsa et al. (2019).

The correction of the telluric absorption lines in each spectrum is the purpose of the CHESS.KNIGHT module (Knowledge of the NIGHtly Telluric lines). In this study, we used the molecfit software, version 1.5.7, provided by ESO (Smette et al. 2015; Kausch et al. 2015). Following the method presented in Allart et al. (2017) and applied in Hoeijmakers et al. (2018, 2019); Seidel et al. (2019), synthetic telluric spectra were computed with the real-time atmosphericprofiles measured above the TNG/HARPS-N location and then fitted to each observed spectrum prior to dividing them. The correction of the changing telluric lines is measured to be at the noise level, especially for the (4ν + δ) water band in the Hα region. We note that the Hβ to Hζ regions are not affected by telluric lines.

thumbnail Fig. 2

Left panel: normalized out-of-transit master spectrum of KELT-9 in the 387–687 nm optical region. The stellar A0 spectral type makes the spectrum dominated by strong Balmer lines, though we noticethe presence of numerous shallow metal lines, as shown in the inset. The Balmer lines (Hα, β, γ, δ, ɛ, ζ), that are covered by the HARPS-N spectrograph wavelength range, are identified by the vertical dashed lines. Right panel: normalized out-of-transit master cross-correlation function (CCF) obtained with a dedicated A0 mask. The averaged stellar line has a measured projected velocity broadening of v sin i* = 112.6 km s−1. The vertical dashed line represent the systemic velocity (−17.74 km s−1).

3 A comprehensive analysis of spectroscopic transit data

Spectroscopic transit observations (at high-resolution) have been widely used to characterize exoplanetary systems. On the one hand, they allow us to measure exoplanetary spin-orbit alignment, via the Rossiter-McLaughlin (RM) effect, with classical velocimetricmeasurements or with the tomographic and reloaded methods (Queloz et al. 2000; Collier Cameron et al. 2010; Cegla et al. 2016). On the other hand, the same data also allow us to measure exoplanetary transmission spectra (Wyttenbach et al. 2015). These two signals can be entangled between each other or with any non-homogeneous stellar surface effects, suchas the center-to-limb variation (CLV) of spectral lines, or stellar activity (Czesla et al. 2015; Barnes et al. 2016; Cauley et al. 2017a). The method presented here aims to separate and to accurately measure the stellar and planetary spectra (see Fig. 3). Even though we are only taking some major effects into account, our method can be, in principle, adapted to correct for any stellar feature. An upcoming paper by Bourrier et al. (in prep.) will establish a rigorous formalism and present an accompanying pipeline (but see also Bourrier et al. 2020a; Ehrenreich et al. 2020, who already used the same method). Our method is implemented as in on our previous work (Wyttenbach et al. 2015, 2017), using the CHESS modules BISHOP (Basic Insights on the Stellar pHotosphere Occulted by the Planet) and QUEEN (QUantitative Extraction of Exoplanet atmospheres in transmissioN), and applied to the KELT-9 HARPS-N data.

3.1 On the extraction of a transmission spectrum

Because of the fluctuations of the Earth’s atmospheric transmission, ground-based transit spectroscopy observations at high-resolution have to be self-normalized. The stellar spectra f(λ, t) are divided by their respective values in a reference wavelength band ⟨λref⟩ chosen in the continuum near the region of interest to get the normalized spectra: . The master spectra in- and out-of-transit Fin(λ) = ∑t∈inf(λ, t) and Fout (λ) = ∑t∈outf(λ, t) are self-normalized in the same manner, giving and . We now base our approach from the schematic shown in Fig. 3 and will proceed step by step.

We first treat an ideal case where there is no exoplanetary atmosphere, and we only consider the stellar rotation and the limb-darkening (and ignore other stellar effects). We are interested in extracting the Rossiter-McLaughlin effect. A possibility is to use the reloaded method (Cegla et al. 2016; Bourrier et al. 2017a, 2018): (1)

where 1 − δ(t) represents the value of the photometric light curve (white light curve), and δ(t) is the transit depth at time t. We compute the light curve with a Mandel & Agol (2002) model using the batman code (Kreidberg 2015) and the system parameters from Table 1. While the observed light-curve contains errors, these latter have a negligible effect on the normalization of the spectroscopic observations (Cegla et al. 2016; Bourrier et al. 2017a). In this ideal case, fres (λ, tin) represents the stellar flux from the planet-occulted photosphere regions, that is to say the local stellar spectrum slocal (λ, t) emitted from the region behind the planet at time t. The measurement and modeling of fres(λ, tin) over time during the transit allow us to determine the sky-projected spin-orbit angle (Cegla et al. 2016). The Rossiter-McLaughlin measurements are mostly done with CCFs in the velocity domain (because of the boost in signal). However, we underline that this measurement is also possible in a single spectral line.

We continue by treating a second ideal case where an exoplanet with an atmosphere transits a star deprived of its rotation (and of any other stellar features, except the limb-darkening). Since the star is not rotating, we know that the shape of the local stellar spectra slocal (λ, t) is not time-dependent and is equal to the out-of-transit spectrum . In this case, the only difference between the shape of the in- and out-of-transit spectra is the transmission through the planet atmosphere,implying that the stellar flux is multiplied by the planet apparent size. Then, the formula proposed by Brown (2001) and implemented by considering the change in radial velocity of the planet in Wyttenbach et al. (2015), is applied in a straightforward manner: (2)

The −′(λ) denotes the normalized spectrum ratio (averaged over Nin in-transit spectra), or the averaged transmission spectrum, and is equal to . The |p expresses the spectrum Doppler shift according to the planet radial velocity and places the observer in the planet rest frame. As for the Eq. (1), we also use the light curve normalization 1 − δ(t) to put the in-transit spectra to their correct flux levels compared to the out-of-transit spectra. This correction is also necessary for the planetary atmosphere since the absorption depends on the relative stellar flux that is passing through the atmosphere (Pino et al. 2018). This means that atmospheric absorptions were slightly overestimated in previous works. The calculation of Eq. (2) remains precise enough for slowly rotating stars (Wyttenbach et al. 2017).

We then present a comprehensive case that considers the presence of the planet atmosphere and an inhomogeneous stellar surface. The local stellar spectra slocal(λ, t) vary as a function of the stellar rotation and other effects (center-limb variation, convective blueshift, pulsations, gravity darkening, active regions, etc.). In transit configuration, we note that the light source transmitted through the planetary atmosphere is emitted from the region behind the atmosphere limb annulus4. We assume that the light source slocal,atm(λ, t) coming from this annulus (with a surface ratio Σatm which is several times ) is equal in shape to the local stellar spectrum slocal(λ, t) (that comesfrom the surface behind the planet), but differs in intensity from it by a factor Σatmδ, meaning that slocal,atm∕Σatm = slocalδ. The light source slocal,atm is absorbed by the exoplanet atmosphere and constitute slocal,atmstr∕Σatm, str (λ, t) being the transmission spectrum wesearch to measure5. Owing to the absorption by the exoplanet atmosphere, the residuals fres (λ, t), computed with Eq. (1), are not equal to the local stellar spectrum slocal (λ, t), but contain the sum of slocal(λ, t) and of str(λ, t) ⋅ slocal,atm∕Σatm. Hence, the spectrum ratio, (3)

is equal to the surface of the atmosphere limb (alone) times its wavelength dependent absorption, meaning that is the transmission spectrum. Taking into account that the shape of the out-of-transit flux is relative to the continuum reference band ⟨λref⟩, we get the final transmission spectrum formula relative to ⟨λref⟩: (4)

Placing the observer in the planet rest frame is the last performed operation before the sum. We note that by setting slocal (λ, t) = δ(t) ∀λ, we can simplify Eq. (4) to get Eq. (2).

We first acknowledge that the presence of the planet atmosphere, as anticipated by Snellen (2004), magnifies the RM effect since the amplitude of the RM effect is proportional to the transit depth (see also Dreizler et al. 2009). Additionally, at high-resolution, the planetary orbital velocity plays a major role. First, the change in shape of a stellar line due to the planet atmosphere happens when the planetary radial velocity is smaller than the full line broadening (which is dominated by the stellar v sin i*). Second, the line shape anomaly due to the RM effect is also altered by the planetary atmosphere when the two signals overlap in terms of velocity. Hence, the presence of a planet atmosphere biases the RM effect, in a way that is not a simple amplification, but rather by creating new anomalies that depend on the planetary radial velocity and of the spin-orbit angle. Despite Snellen (2004) and Czesla et al. (2012) analyzed specific stellar lines to find RM anomalies, these latter were not modeled with the effect of the planet orbital velocity. This effect has been precisely taken into account only recently by Borsa et al. (2019), who called it the “atmospheric Rossiter-McLaughlin” effect. This effect is particularly important for ultra-hot Jupiters who share many spectral lines with their host stars (see more details in Bourrier et al. 2020a; Ehrenreich et al. 2020).

Similarly to the RM effect being modified by the planet atmosphere, the planetary atmosphere absorption depends on the stellar flux that is passing through it. For example, when the local stellar and the atmospheric spectra overlap, the atmospheric absorption is lowered because less flux is coming from the star. This effect is automatically taken into account and corrected for in Eq. (4) thanks to the division by slocal (λ, t).

This method relies on the knowledge of the local stellar profile slocal (λ, t). One approach is to model slocal(λ, t) by using stellar spectrum model at different limb-darkening angles (e.g., Yan et al. 2017; Yan & Henning 2018; Casasayas-Barris et al. 2019). Our approach is to empirically fit the observed local stellar profile (see Sect. 3.2.1). We emphasize the caveat, particularly during the overlap, of the uncertainty whether the local stellar profile remains constant.

thumbnail Fig. 3

Schematic of a spectroscopic transit observation of an exoplanetary system. The difference between an out-of-transit spectrum (left) and an in-transit spectrum (middle) gives us two essential quantities: the exoplanetary atmosphere transmission spectrum and the local stellar spectrum of the surface that is hidden behind the exoplanet (right). The temporal evolution of the stellar local spectrum will mostly depend on the projected stellar rotation (represented by the black arrow going through the star from pole to pole), the planetary orbit parameters (represented with the green arrows) and particularly the projected spin-orbit alignment angle. The temporal evolution of the exoplanet atmosphere spectrum will mostly depend on the atmosphere size, the planetary orbital projected velocity, and the local stellar spectrum that is transmitted through. Hence, it is necessary to separate and to measure accurately each component to retrieve the unknown parameters. The present figure is partially inspired by figures in Dravins et al. (2015, 2017); Dumusque et al. (2014); Cegla et al. (2016); Bourrier et al. (2017a).

3.2 Analysis of the local stellar photospheric CCFs

The analysis of the stellar CCFs is conducted on the CCFs constructed with the aforementioned A0 mask (similarly to Anderson et al. 2018; Borsa et al. 2019). The CCF residuals are computed with Eq. (1) in the velocity domain (see Fig. 4). They clearly show the local stellar CCFs (the CCF of the spectrum of the stellar surface that is hidden behind the exoplanet, see Sect. 3.2.2), the CCFs of some stellar pulsations (see Sect. 3.2.1) and the planetary atmosphere CCFs trace. The planet atmospheric trace is an independent confirmation of the metal (mainly Fe I, Fe II) detection of Hoeijmakers et al. (2018, 2019) and Borsa et al. (2019), and is not further discussed in this study. The main goal of this section is to focus on the pulsation signal, and then on the local stellar CCFs to have an updated value of the spin-orbit angle λ and to use this value as a prior in our Balmer lines analysis.

thumbnail Fig. 4

Map of the residual CCF flux during the transit of KELT-9 b. The observations of two nights are binned by 2 in phase. The 1st and 4th contacts are traced with the horizontal dashed white lines; the 0 km s−1 and ± v sin i* are traced with vertical dashed white lines. Top panel: local stellar spectrum tracing the spin-orbit misalignment through the Rossiter-McLaughlin effect. The solid black line shows our best fit. The exoplanet transmission CCF, that has a slanted shape following the orbital velocity (dashed black line), traces the metal detections by Hoeijmakers et al. (2018, 2019). Bottom panel: residual CCF flux after the removal of the local stellar and of the planetary atmosphere spectra. The residuals are dominated by the stellar pulsations, whose strengths are 10 × lower than the local stellar signal (see the color scale). The pulsations are due to pressure-mode in a δ Scuti-type star.

3.2.1 Spectroscopic detection of stellar pulsation

Prior to discussing the reloaded RM study, we first mention the spectroscopic detection of stellar pulsation. As all the different signals can be disentangled, we study the CCFs cleaned from the local stellar and the planetary atmosphere CCFs, but still keep the pulsation signals (Fig. 4). A rotationally broadened model (including convolution with the instrumental profile and with a local stellar CCFs width) are fitted to each CCF (Gray 2005; Anderson et al. 2018; Borsa et al. 2019; Dorval et al. 2020). We measure an averaged v sin i* value of 112.60 ± 0.04 km s−1. We note that our CCFs (see Fig. 2) and v sin i* values are Hβoxsimilar to those published in the results of Borsa et al. (2019). The instrumental profile (2.7 km s−1) with the local stellar CCFs (see Sect. 3.2.2) have an averaged broadening of ξ = 9.17 ± 0.05 km s−1, making the total broadening ~113 km s−1. The stellar local profiles are fully described by Gaussians. Indeed, the local stellar profile have a small broadening, despite the high disk integrated rotational broadening (ξv sin i*). This is because of the high observing cadence, the small impact parameter, and the inclined orbit (see Sect. 3.2.2). We measure a sinusoidal oscillation in the time series of the stellar v sin i* around the averaged value. The oscillation signal is not in phase with the planetary transit. However, this signal has a period P = 7.54 ± 0.12 h, in phase within the two nights separated by one year (Fig. 5). This period corresponds to the stellar pulsation signal seen in the KELT-9 TESS light curve of Ppuls = 7.5868 ± 0.0009 h (Wong et al. 2019), and is compatible with pressure modes (p-mode) in a δ Scuti-type star. This is therefore an independent spectroscopic detection of the KELT-9 δ Scuti-type pulsation. The first and second CCF moments are clearly described by a single sine function with a period Ppuls (the CCF radial velocities are the first CCF moment, and the v sin i* are proportional to the second CCF moment), hinting to a sectoral mode oscillation (the azimuthal order of the mode m is equal to ±l, the degree of the mode, see e.g., Balona 1986; Aerts et al. 1992, 2010). A complete analysis of the pulsation mode is out of the scope of thispaper, but a follow-up study could bring complementary information to the spin-orbit analysis, in particular, for the stellar inclination.

thumbnail Fig. 5

Stellar pulsation signal seen from variation in the stellar v sin i*. The stellar CCF are cleaned from the local stellar and planetary atmosphere spectra, and fit with a stellar rotation profile. The pulsation period is Ppuls = 7.54 ± 0.12 h, and is consistent between the two nights (in dark blue and green, respectively) that are separated by one year. The averaged v sin i* is reported for each night with the horizontal dashed lines, and has an averaged value of v sin i* = 112.60 ± 0.04 km s−1.

3.2.2 “Reloaded” Rossiter-McLaughlin analysis

For the reloaded RM analysis, we treated the CCF residuals in a similar way to that of recent studies by Bourrier et al. (2020a); Ehrenreich et al. (2020). We first removed the signals from the stellar pulsations in the out-of-transit phases by fitting and removing as many Gaussians as there are of individual pulsation signals, allowing us to recompute a new out-of-transit master, from which we recompute again the CCF residuals. At this stage, we remove the stellar pulsation in-transit and the planet atmosphere signal with Gaussian fitting. By doing so, we are left with the local stellar CCF alone (slocal; see e.g., Cegla et al. (2016); Wyttenbach et al. (2017); Bourrier et al. (2020a) for more descriptions of the local stellar CCF), from which we measure each local velocity (Fig. 6). We discard the velocities coming from the overlapping region with the atmospheric signal, because they are biased (Borsa et al. 2019). Local stellar CCFs are considered to be contaminated by the planet atmosphere if the local CCF full width at half maximum (FWHM) overlaps the planetary atmosphere CCF FWHM (10 data points between phase −0.025 and 0.000). We noticed a posteriori that discarding or keeping the overlapping region changes our results by about 1σ. We also discard the outliers points at the egress, for which the CCFs (and velocities) are not detected following the S/N criteria of Cegla et al. (2016) and Bourrier et al. (2018).

We model the residual CCF velocities with the formalism of Cegla et al. (2016), that computes the brightness-weighted average rotational velocity behind the planet during each exposure. The model parameters are the spin-orbit angle, and the stellar velocity field. The latter could be described by a rigid rotation (with v sin i* as a parameter) or a solar-like differential rotation law, which could be expected for A-type star (Reiners & Royer 2004; Balona & Abedigamba 2016). In this case, we have the three parameters veq, i*, and α which describe a solar-like differential rotation law (Cegla et al. 2016). We use emcee (Foreman-Mackey et al. 2013) to explore the parameter space and to find the best models and their parameter confidence intervals (CI). For each of the two nights and the combination of the two nights, we use 200 walkers for 1500 steps and a burn-in size of 500 steps. We use the Bayesian information criterion (BIC) to compare between models. For the solid rotation case (BIC = 151.5), we found v sin i* = 122.5 ± 0.6 km s−1, and λ = −85.65° ± 0.09° (see Fig. A.1). The v sin i* differs between the two nights and is also bigger than our spectral fit, and the one of Borsa et al. (2019). Despite this, the spin-orbit angle is consistent with Gaudi et al. (2017) and Borsa et al. (2019). The differential rotation case (BIC = 147.4) has a Δ BIC = 4.1 in its favor,which is interpreted as “positive evidence” but not a confirmed detection (σ ~ 2.2). In this scenario, we found , km s−1, , and λ = −85.01° ± 0.23° (best solution, see Fig. A.2). We note that λ is still compatible (at 3σ) with the rigid rotation case, and that the corresponding v sin i* = 116.9 ± 1.8 km s−1 is more consistent with our spectral fit. Another differential solution, with a ΔBIC = 3.8, was found with α = 0.09 ± 0.03, veq = 155 ± 18 km s−1, , and λ = −85.10° ± 0.25°, with a verysimilar local stellar trace to the previous solution, but with a higher km s−1 (see Fig. A.3). Our two differential rotation solutions give a true spin-orbit angle of ψ = 89.6° ± 0.7 and , respectively.The best fit, shown in Figs. 4 and 6, is later used for our Balmer lines extraction.

The differential rotation solutions give equatorial velocities that are typical of A-type stars (Zorec & Royer 2012); furthermore an A-type star is likely to show differential rotation (Zorec et al. 2017). The different solutions have a solar-type differential rotation (α = 0.08), and an “anti-solar”-type (α = −0.1), respectively, meaning that in the latter case, that the high latitudes are rotating slightly faster than the equator. The equatorial velocities found imply that the stellar rotation periods are Prot = 0.7 ± 0.1 day (~17.1 h) and Prot = 0.56 ± 0.09 day (~13.6 h), respectively. We note that these periods do not correspond to the pulsation period of ~7.6 h. This means that the star is turning at ~33–56% of its theoretical critical Keplerian break-up velocity (vcrit ~ 416 km s−1). With such an equatorial velocity, the angular rotation parameter ω = Ωrot∕Ωcrit is about 0.45–0.75 (Maeder 2009). This implies that gravity darkening effects become significant and measurable, for example by TESS or high quality photometry (as already done by Wong et al. 2019; Barnes et al. 2011; Ahlers et al. 2015, 2020). For example, with ω = 0.45–0.75, the star becomes oblate and the equatorial radius becomes up to ~4–12% bigger than the polar radius (Georgy et al. 2014). In a gravity darkened star, the polar regions are hotter and more luminous than the equatorial regions. Because of the stellar inclination measured in the differential rotation scenario, we are mainly measuring light from the polar region of the star. Hence, at ω = 0.45–0.75, the stellar luminosity is overestimated up to ~5–15%, and the effective temperature up to ~1–3%, depending on the stellar inclination toward the observer (Georgy et al. 2014). This may be the source of the differences in retrieved parameters between Gaudi et al. (2017) and Borsa et al. (2019). Thus, a follow-up study with high-resolution spectroscopy and precise photometry, which take care of gravity darkening effects, must be undertaken, and will allow us to refine our solution and to choose between different degenerate geometries (e.g., Zhou et al. 2019; Dorval et al. 2020; Ahlers et al. 2020). In particular, the gravity darkening provides a third possibility to measure the stellar inclination alongside to the Rossiter-McLaughlin effect and the stellar pulsation mode analysis. Finally, because of the polar orbit, the small impact parameter and the stellar inclination, the planet may be transiting near the pole during the transit (Wong et al. 2019). Such a transit configuration is favorable to measure star-planet interaction because of anisotropic stellar flux and winds (“gravity-darkened seasons”, Ahlers et al. 2020).

thumbnail Fig. 6

Velocities of the local stellar surface hidden behind the planet during its transit. The observations of the two nights are shown in dark blue and green, respectively. The first to fourth contacts are traced with the vertical dashed lines. The overlap region between the local stellar and atmospheric signals (between phases −0.025 and 0) is shown in red and is not used for the fit. The spin-orbit misalignment is modeled with a solid body rotation (in orange), and with a differential rotation (in light blue). The differential rotation model is slightly preferred.

3.3 The transmission spectrum of KELT-9 b in the Balmer lines regions

To extract the planetary atmosphere information from the spectra, we follow the procedure that leads to Eq. (4). We assume the system parameters from Table 1. In particular, for each Balmer line we use the same white light curve from Gaudi et al. (2017) to normalize the spectra because it is the only one available6. To correct the transmission spectra from the local stellar photospheric signal, we must know slocal(λ, t) for every Balmer line. We assume that the local stellar signal in the Balmer lines follows the CCF one in terms of velocity. For each Balmer line, and for each phase, we fit a Gaussian to the local stellar spectrum, with the velocity centers following our best reloaded RM model from Sect. 3.2. The use of the differential rotation model instead of the rigid rotation model does not change our final line extraction. Indeed, the two models are similar and the velocity errors from the CCFs are much smaller than the one present in a single Balmer line. Thus, we can fix the velocities to the best reloaded RM model in our fit. The Gaussian offset values follow the light curve δ(t). We therefore only fit for the contrast and the FWHM values. In a second step, we averaged all the contrast and FWHM values, and fitted again a Gaussian (for each line, for each phase) with starting positions coming from the averaged values, and by not allowing us the fit to exceed ±1σ value of the averaged contrast and FWHM. Whenever the local stellar spectrum signal is too weak for the Gaussian fit to converge (or to be significant), we assume the averaged contrast and FWHM values. We note that the stellar pulsations, despite being expected in the Balmer lines, are too weak to be detectable given the S/N of a single line. Thus, the pulsation contaminations in the Balmer lines are at the noise level. Also, we removed low-frequency trends in the residual continuum by fitting polynomials (Seidel et al. 2019). The knowledge of slocal(λ, t), and the Kp from Hoeijmakers et al. (2019) allows us to compute the weighted averaged transmission spectra with Eq. (4)7. We discard the spectra where the atmospheric signal overlaps with the local stellar signal (−0.027 ≤ Phase ≤ 0.0015), because we cannot fit slocal(λ, t) with a Gaussian for these phases. Figure B.1 shows a representative example of the effect of the local stellar photospheric signal correction (RM effect correction). We measure the RM effect to have the same relative amplitude for the different Balmer lines and to be corrected down to the noise level. We compute the errors on the absorption lines by considering the systematic noise, empirically measured by the standard deviation into the continuum. We also consider the stellar line contrasts that decrease the S/N in the absorption region (see Table 2). This is particularly relevant for the blue part of the detector. Each line center takes the systemic velocity into account (e.g., 6562.408 Å for Hα); the reference passbands for the normalization of every spectrum are centered 6 Å on each side of the line center, and are 6 Å wide (e.g., [6553.408:6559.408]∪[6565.408:6571.408] for Hα). For the different Balmer lines, our main results are:

Hα

The Hα line is detected with a high confidence of 24.3σ, which makes the planetary trail visible for each phase (Fig. 7). Figure B.1 shows that the magnitude of the RM effect on the Hα line is about 0.3% in the blue part of the line (around 30% of the line core contrast), and about ≲ 0.05–0.1% in the other parts of the line. The RM effect correction is made down to the noise level. We note that our results for Hα are in conflict with previous literature: by ~3σ with Yan & Henning (who used the CARMENES spectrograph 2018), by ~4.5σ with Cauley et al. (2019, PEPSI) and by ~1σ with Turner et al. (2020, CARMENES). The difference could come from our new extraction method, the different resolutions, or different normalizations. However, we note that Casasayas-Barris et al. (2019) reported Hα detections in KELT-20 b/MASCARA-2 b, that are different for the HARPS-N and CARMENES instruments, despite using the same data reduction. It is thus likely that there are instrumental systematics that are not yet taken into account (detector non-linearity, blaze and bias removal, spectrum interpolations, moonlight pollution, etc.). A follow-up study exploring such effects must be undertaken.

thumbnail Fig. 7

Left: detection of the planetary Hα line (significance of 24.3σ). Top left panel: map of the residual spectra ±200 km s−1 around the Hα line in the stellar rest frame (Eq. (1)). The observations of two nights are binned by 2 in phase. The planetary atmosphere Hα absorption follows the orbital velocity (solid black line). The local Hα stellar line follows the reloaded RM model (dashed black line). The first and fourth contacts are traced with the horizontal dashed white lines, the 0 km s−1 and ± v sin i* are traced with vertical dashed white lines. Bottom left panel: weighted averaged Hα transmission spectrum (corrected from the local stellar photospheric signal, Eq. (4)) in the planet rest frame (light gray), also binned by 15× in black circles (the continuum is set to the white light curve transit depth δ = 0.00677). The verticalblack dashed line shows the line center taking the systemic velocity into account. We show the Gaussian fit (blue) with contrast of 0.91 ± 0.04% and FWHM of 44.3 ± 1.8 km s−1. Right: same as for Hα, but for the Hβ line (detection significance of 16.1σ, see Table 3).

Hβ

The Hβ line and its trail are also clearly detected (16.1σ, Fig. 7). Our Hβ detection also disagrees with the result of Cauley et al. (2019, by ~2σ) and is subject to the same caution. Another explanation for the Hα and Hβ differences between studies could be intrinsic variation in the atmosphere. We also note that the line absorption is stronger in the second part of the transit than in the first part, hence the average measured line could arise from different weights given to each spectrum (due to the observation conditions).

Hγ

For Hγ the atmospheric trail is barely visible, since the total spectrum 7.5σ detection is distributed to ~1σ in each spectrum (Fig. 8). Hβ and Hγ have only been reported in the atmosphere of three exoplanets: HD 189733 b, KELT-9 b, and KELT-20 b/MASCARA-2 b (Cauley et al. 2016, 2019; Casasayas-Barris et al. 2019).

Hδ

The Hδ detection in the KELT-9 b atmosphere (4.6σ) is the first of its kind to our knowledge (Fig. 8). Despite its significance, the line is probably slightly biased because the local stellar signal becomes less visible and is more difficult to correct for. This may be due to the lack of signal in the bluer part of the detector.

Hɛ

The Hɛ absorption (2.9σ) is just below the significance level of 3σ. Therefore it should be treated as a hint of detection (Fig. 9). Additional data are needed to better constrain the Hδ and Hɛ line shapes.

Hζ

Hζ is not detected in our data set (Fig. 9). The most probable reason is due to the lack of flux in the very blue part of the HARPS-N wavelength range. Our upper limit value is hinting that the Hζ absorption is smaller than all the other Balmer lines, as expected.

Table 3 presents a summary of our spectral analysis of the Balmer series. From the contrasts of the KELT-9 b Balmer Hα, Hβ, and Hγ lines, we know that neutral hydrogen is present at altitudes between ~1.2 RP and ~1.5 RP (see Figs. 7 and 8). Since this is close to the planetary Roche radius (1.95 ± 0.2 RP), the strong Balmer lines absorptions are a hint of atmospheric escape (Yan & Henning 2018). We note that no lineis significantly blueshifted or redshifted in the line of sight. With an average shift of 0.2 ± 0.6 km s−1, we do not infer any strong winds in the upper atmosphere. This is compatible with the high circulation drag expected in very hot atmospheres(Koll & Komacek 2018, but also see Wong et al. (2019) who found a relatively high recirculation efficiency). Finally, we note that the Hα, Hβ, and Hγ lines have a similar FWHM, with an averaged value of 44.9 ± 1.4 km s−1. For the rest of our analysis, we investigate the strong and consistent Hα, Hβ, and Hγ signals that contain most of the information.

4 Interpretation of the Balmer lines of KELT-9 b

The measured FWHMs of the KELT-9 b Balmer Hα, Hβ, and Hγ lines cannot be explained only by thermal broadening. Indeed, for hydrogen, such a FWHM would require a thermal broadening with a temperature around 40 000 K, far above any theoretical prediction for a thermosphere temperature (e.g., Salz et al. 2016). Since the temperature cannot be measured directly, a rough estimate of the atmospheric properties of KELT-9 b from the Balmer series can be made following Huang et al. (2017, Sect. 2) as done in previous studies (Yan & Henning 2018; Turner et al. 2020). The estimation of Huang et al. (2017) assumes an isothermal hydrostatic atmosphere with a mean molecular weight μ = 1.3 and a temperature T = 5000 K. For KELT-9 b, a better assumption is an atmosphere dominated by ionized hydrogen (μ = 0.66) and a thermosphere temperature of T = 10 000 K (Fossati et al. 2018; García Muñoz & Schneider 2019). With these numbers, we find that to explain the FWHM8 of our Hα detection, one needs an optically thick layer (τ ~ 40) of excited hydrogen H I(2) with a number density of n ~ 8.2 × 103 cm−3 along a typical line path of L13.8 × 109 cm, where Nsc ~ 6 is the number of absorbing scale-heights (Madhusudhan et al. 2014; Madhusudhan & Redfield 2015; Huang et al. 2017). This means that at an altitude of NscH ~ 0.3 RP (H ~ 6000 km), which is about half of the line core contrast, the Hα line core have τ ≫ 1. Since the line is optically thick in most of the layers around that altitude, the line FWHM is larger than the thermal width (Huang et al. 2017). Also, this is consistent with our observations since the Hα core is observed at ~1.5 RP, where τ becomes close to unity (Lecavelier Des Etangs et al. 2008; de Wit & Seager 2013; Heng et al. 2015; Heng & Kitzmann 2017).

Another estimate of the atmospheric conditions can be done with the Lecavelier Des Etangs et al. (2008) formalism, which again assumes an isothermal hydrostatic atmosphere. Indeed, one can measure the atmospheric scale-height H thanks to the linear correlation between the line contrasts (or the altitudes z) and the hydrogen Balmer line oscillator strengths (see the log10gf values in Sect. 4.1). Knowing that Δz = HΔln(gf) = Hln(10)Δlog10(gf), the slope of the z–versus– ln(gf) linear fit gives us H. We computed that H = 9600 ± 1400 km around the average altitude of absorption (1.3 RP). When the gravity is evaluated at 1.3 RP (~ 0.59 gsurf), we have that Tμ ≃ 13 500 ± 2000 [K/u]. Assuming a μ between 0.66 and 1.26 (atmospheric composition dominated by a mixture of ionized and neutral hydrogen and helium), the upper atmosphere temperature is constrained between 7500 and 19 500 K. The lower half temperature range is more probable because the thermal ionization of hydrogen decreases μ9.

Finally, our line contrast measurements and the aforementioned assumptions lead to the same conclusions as in Yan & Henning (2018) of a Jeans mass loss rate of ~1012 g s−1, although hydrodynamic escape is more likely the mass loss mechanism in KELT-9 b (Fossati et al. 2018).

While this kind of estimation is useful for getting a qualitative idea of the planet aeronomic conditions, the data quality calls for a quantitative estimate. A few forward self-consistent models of the hydrogen Hα line exist (Christie et al. 2013; Huang et al. 2017; Allan & Vidotto 2019; García Muñoz & Schneider 2019). They explore different physical questions and have different atmospheric structures, but none of them are linked to any retrieval algorithms. Retrieval approaches (on similar absorption lines) have been made by Cauley et al. (2019); Fisher & Heng (2019); Seidel et al. (2020), but without the mass loss rate being a parameter. Thus, we decided to follow our own retrieval approach with a model adapted to our data set.

thumbnail Fig. 8

Same as for Fig. 7, but for the Hγ line (left) and the Hδ line (right). The Hγ and Hδ lines are detected with a significance of 7.5σ and 4.6σ, respectively (see Table 3).

thumbnail Fig. 9

Same as for Fig. 7, but for the Hɛ line (left), and the Hζ line (right). The Hɛ and Hζ lines are not statistically detected (see Table 3).

Table 3

Summary of the hydrogen Balmer line detections.

4.1 The PAWN model

We model the planet atmosphere in transmission with the customized model PAWN (PArker Winds and Saha-BoltzmanN atmospheric model) which is included in our CHESS framework as a tool for modeling exoplanetary thermospheric lines in transmission spectroscopy. Our purpose is to constrain fundamental parameters of the thermosphere regions such as the temperature, the mass loss rate and the hydrogen number density (local thermodynamical condition), if possible. Another goal is to have a sufficiently fast model to feed a retrieval algorithm. Hence, our modeling holds several assumptions.

We model a 1D planetary atmosphere, with a special interest in the upper atmosphere part (thermosphere). A logarithmic grid with nr radius r layers is chosen between 1 and rup RP. For the hydrogen Balmer lines, rup = 2 ensures the convergence of the models. We model the atmospheric structure ρ(r) with a hydrostatic or a hydrodynamic structure. For the hydrostatic solution, we have the base density ρ0 (r = 1) as a free parameter, and consider that the pressure scale-height H = kBTμg can change with the altitude according to T(r),  μ(r),  g(r) (see e.g., theπη code in Ehrenreich et al. 2006; Pino et al. 2018; Seidel et al. 2020). In the hydrostatic case, we can estimate a lower limit to the mass loss rate, by computing the Jeans escape rate at the planetary Roche radius, rRoche (Eggleton 1983; Ridden-Harper et al. 2016), or at the exobase radius, rexobase, if below (Lecavelier des Etangs et al. 2004; Tian et al. 2005; Salz et al. 2016; Heng 2017): (5)

where FJeans is the Jeans escape flux (Shizgal & Arkos 1996; Hunten et al. 1989; Chamberlain & Hunten 1987). For the hydrodynamical solution, we use the isothermal Parker wind solution following the formalism of Heng (2017)10, Oklopčić & Hirata (2018), and Lampón et al. (2020). Our free parameter for the atmospheric structure becomes the atmosphericmass loss rate: (6)

where ρ(r) and v(r) solve the Parker equations (Parker 1958; Lamers & Cassinelli 1999) from the sonic point (rs,   vs,  ρs,  μs) down to the lowest radius in the grid. We note that is constant through the atmosphere. Since the Parker wind solution is isothermal, we keep an isothermal profile for both the hydrostatic and hydrodynamic structure with the temperature T as a free parameter. We note that the assumption of an isothermal profile for the upper thermosphere (typically between 1.1 and 2–3 RP for hot Jupiters) is a good approximation of aeronomy models (e.g., Lammer et al. 2003; Lecavelier des Etangs et al. 2004; Yelle 2004; García Muñoz 2007; Murray-Clay et al. 2009; Koskinen et al. 2013; Salz et al. 2016). Indeed, above the thermobase where the temperature changes dramatically, the temperature encounters a maximum and is varying much slower through several orders of magnitude in pressure than in the lower regions.

For the atmospheric constituent structure, we use a pre-computed chemical equilibrium table obtained with the chemical equilibrium code of Mollière et al. (2017). The chemical equilibrium (of a variety of neutral and singly ionized atomic, molecular, and condensed species) is computed, with a Gibbs energy minimization, into a grid of pressure (10−15–102 bar) and temperature (1000–20 000 K). The mean molecular weight μ(ρ(r), T(r)), the mass fraction mι (r), and volume mixing ratio nι (r) are given by the chemical grid, where ι is a chosen element and ionization of interest (H I, He I, Na I, K I, Ca II, etc.). For example, in the hydrogen case, we know the ratio of ionized to neutral hydrogen nH IInH I from reading and interpolating the chemicalgrid with the correct T and ρ from our atmosphericstructure. When the Parker wind solution is in use, we link it with the chemical grid at the sonic point down to the lowest radius. We note that for hydrogen, the Gibbs energy minimization is equivalent to computing the Saha law with an electron number density ne coming from the ionization of the hydrogen and from all the other atomic species in the grid. By default, the chemical grid is always used, but it is also possible to assume a constant (free parameter) μ through the atmosphere. In that case the nH IInH I can be constant or can follow the Saha equation with ne accounting only for the electron coming from the hydrogen ionization. In any case, we are ignoring photo-ionization effects that can be important in the upper atmosphere.

Once we know the vertical distribution of our element of interest, we still need to know what is the number density of the electronic level ι(n) from which the transitions we observed are arising. In the case of H I, we are observing the Balmer series in absorption, which arise from photons being absorbed by the n = 2 electronic state, that is to say the H I(2) state. If the atmosphere is in local thermodynamic equilibrium (LTE), the distribution of the electronic states follow the Boltzmann equation (for hydrogen): (7)

where E1 = 13.606 eV, and υ(T) is the hydrogen partition function. The Boltzmann equilibrium is assumed by default in our code. However, because we are probing high in the thermosphere, the collisions between particles become less numerous and the electronic state distribution can depart from the HβoxBoltzmann distribution, the atmosphere is in non-LTE or NLTE (Barman et al. 2002; Fortney et al. 2003; Christie et al. 2013; Huang et al. 2017; Oklopčić & Hirata 2018; Fisher & Heng 2019; Allan & Vidotto 2019; García Muñoz & Schneider 2019; Lampón et al. 2020). For the NLTE case, we define the departure coefficients (e.g., Mihalas 1970; Salem & Brocklehurst 1979; Barman et al. 2002; Draine 2011; Hubeny & Mihalas 2014) as: (8)

Then, we can let the different βH I(n) be free parameters (or equivalently all the nH I(n) are free parameters).Another possibility is to force all the nH I(n) to depart by the same amount from LTE, meaning that only one constant free parameter β is used for all nH I(n). We note that a constant nH I(n)nH I through the thermosphere is an assumption justified by Huang et al. (2017) and García Muñoz & Schneider (2019).

When an element and its transitions are chosen (in our case the Hydrogen Balmer series), the opacity κ [cm2 g−1] are computed in each layer of the atmosphere according to the line list of Kurucz (1979, 1992), to the temperature, and to the number density nH I(n). We do assume that the line intensity arises from photon absorption and from stimulated emission (Sharp & Burrows 2007). For the hydrogen Balmer transition from electronic level i to j, we have: (9)

where e is the elementary charge, me is the electron mass, and c the speed of light in vacuum. The fij, gi, gj are the oscillator strength, and the statistical weights of the electronic levels due to their degeneracies, respectively. For the hydrogen Balmer series, we have log10g2f2j = 0.71, −0.02, and −0.447 for j = 3, 4, and 5, respectively. A Voigt profile ΦV, centered on the line transition center νij (or λij), is assumed for the line shape. The Lorentzian wing component considers the natural broadening (for the hydrogen Balmer series, the Einstein A-coefficients are log10(Aj2 [s−1]) = 8.76, 8.78, 8.79 for j = 3, 4, and 5, respectively). There is no pressure broadening. The Gaussian core component takes the thermal broadening into account. We can also add velocity components due to rotation, circulation, winds in the atmosphere (see below). For this exploratory study, we decided not to consider the opacity arising from H (the main source of background continuum opacity, Kitzmann et al. 2018), as a rough estimate shows that its opacity (~102 cm2 g−1) stays smaller than the one from each Balmer line (103 –108 cm2 g−1 in the Gaussian core). Hence, the modeled absorptions are generally robust for the line cores that we observe (Hoeijmakers et al. 2019). The spectral resolution of our model is chosen to be 0.001 Å (λ∕Δλ ~ 106), and the wavelength coverage is 20 Å around each line.

Next, we perform the radiative transfer in a grazing geometry following the method of Ehrenreich et al. (2006); Mollière et al. (2019). For each impact parameter, the optical depth , where is the column mass of the element ι in [g cm−2], is computed through the line of sight (LOS). The transmission function is computedaccording to a pure absorption in each LOS. The transmission spectrum is computed assuming a spherically symmetric atmosphere, weighted along the impact parameter surface annuli, and normalized to the stellar surface.

Within the previous two steps, our line profile can be broadened according to different atmospheric patterns (see Brogi et al. 2016; Seidel et al. 2020). First, we assume a solid body planetary rotation (perpendicular to the orbital plane). Since the planet rotation is axisymmetric, this breaks the spherical symmetry assumed before. In that case, we divide each hemisphere of the planet in a grid of mθ circular sector. The sector grid is evenly separated in cosθ, where θ is the angle going from the equator to the pole back to the equator (anti-clockwise). If the rotation is included, our code becomes essentially a 2D atmospheric code. The planet rotation period PP is the parameter that defines the solid rotation. PP can be a free parameter, or more simply be fixed to the planetary orbital period by assuming that the planet is tidally locked (which is our default choice). It is interesting to note that in this configuration, each LOS has a fixed velocity along the cord, which is given by: (10)

In the case of KELT-9 b, the planetary rotation have a broadening effect of ~ 6.1 km s−1 at 1 RP for each hemisphere (~ 12.2 km s−1 in total).

Second, when the atmosphere is hydrodynamical (Parker winds), a spherically symmetric expansion of the atmosphere occurs with a velocity v(r). This expansion changes the projected velocity of the different components of each LOS, thus creating another broadening source. This Parker wind broadening can be considered in our model, but not simultaneously to the rotation (because the symmetry is not the same, the structure would be 3- and the computational cost would become too large). The LOS velocity is given by: (11)

In practice, for KELT-9b, for most of the mass loss rates (except for unphysical large rates), the upward velocity v(r) is smaller than 0.01–1 km s−1 in the probed layers (i.e., the sonic point rs, where the velocity becomes several km s−1 is far above the Roche radius). This leads to a negligible (much smaller than the instrumental resolution) parker wind broadening (though see Seidel et al. 2020, for other possible scenarios). Thus, for the rest of our study, we decided to ignore this effect in favor of the planetary rotation. We note that these two aspects have been previously explored in Seidel et al. (2020). They were able to keepa 1D atmospheric structure at the price of several approximations in the atmospheric velocity patterns (pseudo-2D/3D).

In the next step, our line profile is modified according to the observational set-up. First, the orbital blurring due to the planet changing velocity during an exposure time is computed with a box-shape convolution (of a width of ~6.8 km s−1 for the KELT-9 b orbit and for the 600 s exposure time). The spectrum is then convolved with the instrument profile, which is assumed to a be a Gaussian (with a width of ~2.7 km s−1 for HARPS-N).Then, the spectrum is binned down to match the data sampling. Finally, the differential normalization (due to the loss of the absolute depth in ground-based observations, see Wyttenbach et al. 2015, 2017; Pino et al. 2018) is performed according to the aforementioned reference passbands of each line. The models, as for the data, are normalized to the white light curve transit depth of KELT-9b (δ = 0.00677). The models are now in a form which is comparable with the data. Our PAWN model is finally linked with the Markov chain Monte Carlo (MCMC) posterior sampling code emcee (Foreman-Mackey et al. 2013). Typically, for one line, one model without rotation is computed in ≲0.1 s (nr = 100), while with rotation, one model takes ~0.5 s (nr = 50, mθ = 14). We note that a model with nr = 50 and mθ = 14 is ≲1% equal to a model with nr = 250 and mθ = 28.

In the next section, we apply our PAWN model to the Balmer lines detection in the KELT-9 b atmosphere.

4.2 Interpretation of the Balmer series in the KELT-9 b atmosphere with the PAWN model

In this section, we investigate the Hα, Hβ, and Hγ signals based on their line shape. All lines are fitted together in a range of 20 Å (6000 data points, notably to cover the normalization passbands). Our retrieval approach uses our PAWN model with the emcee sampler to explore the parameter space and to find the best models and their parameter confidence intervals (CI). For each different set of models, we use ten walkers for each dimension (parameter) during 2500 steps and with a burn-in size of 500 steps. We use the bayesian information criterion (BIC) to compare between models.

For every model, each line center position is a free parameter (three parameters for Hα, Hβ, and Hγ). The line center positions are sensitive to circulation winds in the upper atmosphere (Snellen et al. 2010; Brogi et al. 2016; Seidel et al. 2020). Uniform priors between ±20 km s−1 are used to ensure that our prior cover the full posteriors. We noticed a posteriori that these parameters are not correlated with any of the other parameters. Also, we found that our MCMC results are confirming the line shifts and errors quoted in Table 3, thus we do not show the results for these parameters. The second set of parameters defines the atmospheric structure with the temperature T, and the base density ρ0(r = 1) or the mass loss rate depending on if the structure is hydrostatic or hydrodynamic, respectively (for these latter, we use their log10 values). We use a uniform prior on T between 1000 and 20 000 K. We use a log-uniform prior between −18 and −4 for the log10(ρ0 [g cm−3]). We use a log-uniform prior between 8 and 17 for the . In our base models, we use the Boltzmann equilibrium and there are no other parameters. We first explore this set-up, compare the hydrostatic and hydrodynamic models, before exploring some NLTE models.

4.2.1 LTE models with a hydrostatic or hydrodynamic structure

We first asses how our simplest models are able to fit the observed lines. These two set of models are hydrostatic and hydrodynamic, respectively. They have two free parameters (apart of the line centers), and assume a Boltzmann equilibrium. For the hydrostatic model, we find a solution with a thermosphere temperature of K, and a base density of log10(ρ0 [g cm−3]) = −10.74 ± 0.05 (BIC = 5016.1, Fig. 10). The corresponding Jeans mass loss rate at the Roche radius is computed to be Jeans ~ 6.6 × 1011 g s−1. Interestingly, our hydrodynamical model, which has a BIC = 5017.2, is found to have a similar thermosphere temperature of K, while the mass loss rate is tending towards to the higher value of log10(PW [g s−1]) = 12.8 ± 0.3 (i.e. PW ≃ 5.9 × 1012 g s−1, Fig. 11). Both models are comparable in term of statistics and none can be preferred; this is because the density profiles are not different enough at the altitudes probed. The mass loss rate is in both cases around 1012 g s−1, confirming the result of Yan & Henning (2018). Also, we can safely conclude that the thermosphere temperature, if in LTE, must be near 13 000 K (in line with our first estimate using the Lecavelier Des Etangs et al. 2008 equation). Despite the high temperature, the line FWHMs (44.9 ± 1.4 km s−1) are not dominated by thermal broadening. Indeed, for T = 13 200 K, the hydrogen thermal broadening is ~ 25 km s−1, while the other source of broadening are ~16, ~ 7, and ~ 2.7 km s−1 for the rotational (at r = 1.3 RP), orbital blurring, and instrumental broadening, respectively. This means that the first source of broadening (~ 33 km s−1) comes from optically thick absorptions through most of the lower layers (Huang et al. 2017).

We notice that the line contrasts and FWHMs are both enhanced when the temperature or the mass loss rate are increasing. For the temperature, this is caused by the atmosphere being more extended, and by thermal broadening. For the mass loss rate, this is caused by larger densities present at higher altitudes that increase the LOS opticaldepth, thus boosting the broadening and the contrast. Therefore, the shape of the correlation in the hydrodynamical case (Fig. 11) is explained by the following mechanism: when the temperature increases, the number density of the hydrogen second electronic level builds up. However, as for the stellar case, the total number density of neutral hydrogen decreases because of ionization. This last effect eventually dominates at high temperature. So, to keep a constant absorption corresponding to our observations, a higher density is needed. This is achieved by increasing the mass loss rate (Eq. (6)). Thus, we get the positive correlation between T and PW, which shapeis given by the Boltzmann equation.

4.2.2 NLTE models with a hydrodynamic structure

We explore the NLTE case for the hydrodynamical atmosphere case. We first run MCMC chains by letting all the electronic level number density be free parameters. The prior on each was set to be a log-uniform distribution between −40 and 0. With this test, we found that the model hardly converges, and that the parameters are not Hβoxwell-constrained (see Fig. C.1). The explanation lies in two different families of solutions. One family is distributed around the Boltzmann equilibrium, while the other is away from it, and has a lower temperature and a higher mass loss rate. We therefore further explore these two families.

For the first family, we limit our prior ranges on the between −9, −12, −16, −20 and 0, for n = 2, 3, 4, 5, respectively. The results are shown in Fig. C.2, and the model results in a BIC = 5031.4. The increased number of parameters clearly exclude this model (ΔBIC = −14.2 < −5), despite a better best fit. We still notice that the retrieved temperature lies between 12 600 and 16 600 K, while the mass loss rate lies between 1012 and 1016 g s−1, encompassing the LTE result.

For the second family of solutions in NLTE, we notice that it is mainly dominated by the departure of the second electronic level. We explore this family, by reducing the number of parameters for the number densities. Instead of having all the be free parameters, we decide to use a unique departure coefficient β for all the electronic levels, that is βH I(n) = β for all H I(n). For this last MCMC run, a log-uniform prior for log10(β) was set between −10 and 5. The result is shown in Fig. 12. We confirm that this set-up explores a parameter space close to the second aforementioned family, where the temperature is lower and the mass loss rate is higher. In this case, the best fit is again better, and have a BIC = 5014.8. The ΔBIC = 2.4 from the hydrodynamical-Boltzmann model, meaning a small positive evidence against the LTE model (and only a smaller ΔBIC = 1.3 from the hydrostatic Boltzmann model). In any case, the ΔBIC is not sufficient to conclude in favor of any NLTE models (see the comparison in Fig 13). Nevertheless it is interesting to note a few particularities in this scenario. First, the β parameter is strongly correlated with the temperature and mass loss rate, giving a conservative idea of the uncertainties on these two parameters. Second, the mass loss rate is g s−1 (between 8 × 1012 and 3 × 1015 g s−1). Such a massloss rate could have dramatic consequences (see discussion in Sect. 5). Third, the departure from LTE is (in the 10−8.7–10−2.2 range at 3σ), making the LTE excluded in this framework. Interestingly, our best value is far (but still compatible at ~2.5σ) from the theoretical work by García Muñoz & Schneider (2019) who predicted a departure from Boltzmann to be at the level of ~ 10−3. Our finding could hint to another possible departure from LTE than only the one from the Boltzmann equilibrium (see discussion in Sect. 5).

thumbnail Fig. 10

Correlation diagrams in the case of a hydrostatic atmosphere in the Boltzmann equilibrium.

thumbnail Fig. 11

Correlation diagrams in the case of a hydrodynamic atmosphere in the Boltzmann equilibrium.

thumbnail Fig. 12

Correlation diagrams in the case of a hydrodynamic atmosphere with a departure β from LTE.

5 Discussion

With our PAWN framework, we are able to retrieve a good fit to the data (Fig. 13), and show that the thermosphere temperature (in LTE) is near 13 000 K and the mass loss rate is around 1012 g s−1. If the hydrodynamic expansion of the thermosphere is a valid assumption, then the thermosphere temperature is constrained to be K, and the mass loss rate to be g s−1. This mass loss rate is higher than the energy-limited escape rate of 1010 –1011 g s−1 found by Fossati et al. (2018, see Fig. 14). It is important to note that this study only considers the X-FUV and NUV irradiation from the star, which may be small in the case of a A0-type star without chromosphere. However, García Muñoz & Schneider (2019) presented a new escape mechanism that applies to ultra-hot Jupiters orbiting early-type stars. This mechanism shows that it is important to also consider the NUV and optical Balmer continuum and Balmer series irradiation from the star. In this scenario, the thermosphere is first heated by Lyman photons (Lyman Hβoxcontinuum and especially Lyα), which increases the number density of excited hydrogen in the high atmosphere (Lyman irradiation and de-excitation put the population levels out of Boltzmann equilibrium), thus increasing the Balmer opacity. At a certain point, the excited hydrogen atoms are able to absorb significantly much more flux from the stellar Balmer photons, which constitutes an enormous source of energy in early-type stars. This increases even more the thermosphere temperature, which reaches an equilibrium temperature of ~12 000–15 000 K (3000 K more than with only the X-FUV-NUV flux). In this case the mass loss rate increases to 1012 –1013 g s−1, the atmosphere undergoes a “Balmer-driven” escape, named so because the Balmer energy deposition dominates over all other. This escape rate is well in line with our findings (Fig. 14). Knowing that the age of KELT-9 is ~300–500 Myr (Gaudi et al. 2017; Borsa et al. 2019), a mass loss of this rate over the planetary life time represents a total mass loss of ~0.03–0.05 MJup (~1.5–2% of the planet mass).

In the thermosphere, the excited hydrogen atoms are expected to be in departure from LTE. The main reason being the lower frequency of collisions at high altitudes, making the hydrogen energy levels depart from the Boltzmann equilibrium, typically with a β factor of 10−3 (García Muñoz & Schneider 2019). Despite the data being too low in S/N to show clear evidences for a departure from the LTE, our NLTE model hints at a lower temperature of T = 9600 ± 1200 K, and a higher mass loss rate of . If true, the planet may have already lost 1–2 MJup, and would undergo a catastrophic evaporation that would make the planet barely survive the next ~500–700 My, when its host star will enter the red-giant phase. Such a mass loss rate is not yet physically motivated (though see García Muñoz & Schneider 2019). One possibility would be to explore the real effect of the stellar gravity and magnetic field. Indeed, such effects could push down the sonic point to much lower altitude, forcing the density and velocity to increase at high altitude, hence increasing the mass loss rate (Murray-Clay et al. 2009). Another possibility would be to investigate the effect of the high rotation of KELT-9, which is also relatively young. Indeed, both effects potentially increase the mass loss rate (Allan & Vidotto 2019).

In our NLTE model, we found that the departure coefficient β is on the order of 10−6, while theory predicts 10−3 (only for Boltzmann departure). Thus, our much smaller β coefficient could hint to some other possible LTE departure. Indeed, because β is eventually normalized to the total amount of neutral hydrogen, our finding may signify that there is less neutral hydrogen than in the chemical (Saha) equilibrium by a factor of 103. This impliesthat the LTE departure may also be caused by photo-ionization (Saha equilibrium departure, e.g., Fortney et al. 2003; Sing et al. 2008; Vidal-Madjar et al. 2011; Lavvas et al. 2014). In this case, the ratio of neutral to ionized hydrogen would be 10−3 lower. The other atoms would probably also undergo photo-ionization, making the ionized species even more dominant in the high atmosphere,increasing the importance of ionospheric phenomenons (Helling et al. 2019).

Moreover, our NLTE scenario may be able to explain the strong absorption by ionized atomic species, especially the Ca II and Fe ii, by Hoeijmakers et al. (2018, 2019); Cauley et al. (2019); Yan et al. (2019); Turner et al. (2020). It is interesting that Yan et al. (2019) found a different line ratio between the calcium doublet and triplet for the two ultra-hot Jupiter, WASP-33 b and KELT-9 b. Hence, this type of data is also sensing the local thermodynamical conditions, in complement to the Balmer lines. We can predict absorptions of the calcium doublet Ca II H&K, and infrared triplet (IRT) lines, by taking the two main atmospheric structures of our best fit models of the Balmer lines. To model the calcium lines, we assume the calcium abundance to be solar (χCa = 6.34, Asplund et al. 2009), and rup = 10 (for good convergence). In our code, the calcium population is dominated by ionized Ca II at the temperatures and pressures of our models (see also Kitzmann et al. 2018; Lothringer et al. 2018; Yan et al. 2019, although for very high temperatures, the Ca III and further ionization levels are becoming important). In this case, our model is compatible with the detections of Yan et al. (2019); Turner et al. (2020) only if we assume a departure from LTE on the order of ~10−2, or ~10−4.5 for our two main structures, respectively (Fig. 15). Furthermore, the Ca II opacities are much larger than the Balmer opacities (by ~100 times), but the measured absorptions have about the same contrasts, which hints at the presence of some Ca II NLTE. Because the calcium H&K lines are arising from the ground state of the singly ionized calcium, Boltzmann equilibrium departure is not likely to explain this, because the lack of collision is mostly populating the ground state level (the smaller calcium IRT absorption may be explained by non-Boltzmann distribution since this triplet arises from the first excited level). The β departure in the Ca II case may also be due to a non-solar metallicity, but it is unlikely that the planet atmosphere metallicity is far below the solar one (KELT-9 has a solar metallicity, and planets are expected to generally have a higher metallicity than their host star). Thus, the departure from LTE may be better explained by a second ionization of calcium (through thermal and photo-ionization), which would lower the singularly ionized calcium population. The departure level for Ca II is in line with the extra departure found for hydrogen (we note that the ionization energy of the H I and Ca II have similar value of 13.6 eV and 11.9 eV, respectively). More investigations have to be done to explore if the hypothesis of chemical equilibrium in the thermosphere of (ultra) hot Jupiter is still justified.

One caveat of our framework is that our set of assumptions may be too strong when interpreting, together or not, the hydrogen and calcium absorptions. For example, the assumption of an isothermal atmosphere is probably too strong because temperature gradients and temperature inversion are predicted (Kitzmann et al. 2018; Lothringer et al. 2018; Mansfield et al. 2020). Actually, such a temperature inversion has been detected for the first time by Pino et al. (2020) in the KELT-9 b atmosphere by observing neutral iron lines in emission on the day side. Their detection shows that the temperature inversion happens from the 1 bar region up to the 10−3–10−5 bar region. However, Hoeijmakers et al. (2019) infers that the maximum probable pressure in transmission is in the mbar regime, due to the H opacity in the optical. The consequence for our model is that defining the lower boundary for the radius integration at 1 RP (the measured white light radius) is justified because this radius is set by the, effectively gray, H opacity. At no wavelength should the high-resolution observations probe smaller radii. However, because we assume that the atmosphere is isothermal, at high temperatures this can lead to pressures at 1 RP that are much smaller than mbar. Forcing a fit to continue to mbar pressures increases the line core to continuum amplitude unrealistically. Such a continuation to larger pressures at fixed temperature is also not correct because the atmosphere must transition from the high thermosphere temperature to the planetary equilibrium temperature. Accounting for the decrease in atmosphericscale-height, the line core to continuum amplitude would then stay compatible with our current model. Another effect, not considered in our model, is a potential difference in the vertical distribution of the individual atmospheric gases that are determined by their respective masses, whose effect is to give birth to a heterosphere on top of the homosphere (Shizgal & Arkos 1996). For example, the presence of a heterosphere may also explain the weak calcium absorption, since these atoms would stay in the lower thermosphere. However, this effect of atomic diffusion may be counterbalanced by the hydrodynamic outflow that is able to lift heavy atomic species, that has been observed in exospheres (Vidal-Madjar et al. 2004). A last caveat is the modeling of the exosphere itself. Because the Roche radius of KELT-9 b is low (1.95 ± 0.2 RP), it is expected that hydrogen and metals escape easily. The kinetic broadening in the exosphere would change the shape of the lines, making them asymmetric (Ehrenreich et al. 2015; Bourrier et al. 2015, 2016; Allart et al. 2019). This effect has not yet been seen for the KELT-9 b Balmer lines. Thus, additional data gathering, and modeling will be necessary to get a coherent image of the KELT-9 b atmosphere consistent with all the observations.

thumbnail Fig. 13

Hα, Hβ, and Hγ transmission spectra (light gray, and binned by 10 × in black circles) in the planet rest frame and in velocity space (the continuums are set to the white light transit depth δ = 0.00677). The best PAWN models, in the case of a hydrodynamically expanding thermosphere, and in the Boltzmann equilibrium, are shown in blue. The retrieved temperature is K, and the mass loss rate is g s−1. The best PAWN models, in the case of a hydrodynamically expanding thermosphere, with a LTE departure, are shown in dashed red. The retrieved temperature is K, the mass lossrate is g s−1, and the common LTE departure coefficient is .

thumbnail Fig. 14

Summary of the mass loss rates and their errors found in the LTE case compared to the expected mass loss rates from the theoretical studies of Fossati et al. (2018) and García Muñoz & Schneider (2019). Our finding shows that the KELT-9 b thermosphere probably undergoes a “Balmer-driven” expansion, that leads to a high escape rate. In this case, the total mass loss through the planetary life time (500 My) represents sim 0.05 MJup, which is sim 2% of the planet mass.

thumbnail Fig. 15

Example of PAWN models of the Ca II H line detection by Yan et al. (2019). Our model is compatible with their detection only when assuming some LTE departure. For the Ca II, the departure is likely explained by further thermal and/or photo-ionization not yet taken in our code into account. Additional data gathering, and modeling will be necessary to confirm our interpretation.

6 Conclusion

We analyzed two optical high-resolution transit observations of KELT-9 b gathered with HARPS-N, with our CHESS framework. We found a spin-orbit misalignment of λ = −85.01° ± 0.23°, and positiveevidence in favor of a differential rotation model. We also confirmed the presence of stellar pulsations with a period of Ppuls = 7.54 ± 0.12 h. We use the reloadedRM analysis to correct the Balmer lines from the local stellar photospheric spectrum effect in transmission spectroscopy. We detected theHα, Hβ, Hγ, and Hδ absorptions coming from the planet’s upper atmosphere, with high confidence (between 24.3 and 4.6σ) during both nights. The Hα, Hβ, Hγ lines probe altitudes up to ~1.2 to ~1.5 planetary radii; they do not show any significant velocity shifts, and share a large FWHM of 44.9 ± 1.4 km s−1 that is dominated by optically thick absorptions. We presented a new model, CHESS.PAWN, to interpret the absorption in the Balmer series in exoplanet atmospheres that is linked with a retrieval algorithm. We constrained the thermosphere temperature to be K, and the mass loss rate to be g s−1, under the assumptions of a hydrodynamic expansion and of LTE. This escape rate is fully explained by the “Balmer-driven” escape mechanism proposed by García Muñoz & Schneider (2019). NLTE effects are not yet significantly detected with this data set.

Our work shows that with the infrared metastable Helium triplet (λ10830 Å, Seager & Sasselov 2000; Oklopčić & Hirata 2018; Spake et al. 2018; Allart et al. 2018; Nortmann et al. 2018), the Balmer series is an alternative new tool for probing exoplanet evaporation, alongside with classical UV observations (Vidal-Madjar et al. 2003). Furthermore, these two probes seems to be complementary since the helium lines are mainly sensitive to G and K-type exoplanetary systems (Oklopčić 2019), while the Balmer lines are more sensitive to early-type (A-type to G-type) systems (Jensen et al. 2012, 2018; Yan & Henning 2018; Casasayas-Barris et al. 2018, 2019; Cabot et al. 2020; Chen et al. 2020). Finally, such observations will benefit from the ESPRESSO spectrograph (Pepe et al. 2013, 2014). On the one hand, ESPRESSO will increase the S/N of the absorption lines (Ehrenreich et al. 2020; Chen et al. 2020). On the other hand, in combination with near-infrared spectrographs (CARMENES, SPIRou, NIRPS), it will be possible to detect more absorption lines arising from different electronic levels of a given atomic species. For example, the combination of the Balmer and Paschen (or Brackett) series for the hydrogen, of the infrared triplet, and of the D3 line for the helium, or of the H&K doublet and of the infrared triplet for the calcium, will allow us to infer the local thermodynamic state. This will eventually allow us to put constraints on several chemical and physical aeronomic processes.

Acknowledgements

A.W. acknowledges the financial support of the SNSF by grant number P2GEP2_178191 and P400P2_186765. This work has been carried out within the frame of the National Centre for Competence in Research “PlanetS” supported by the Swiss National Science Foundation (SNSF). This project has received funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation program (grant agreement number 724427; project FourAces). I.S. and P.M acknowledge funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation program (grant agreement No. 694513 and 832428; ExoplanetBio). L.P. acknowledges that the research leading to these results has received funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation program (grant agreement No. 679633; Exo-Atmos). We thank the anonymous referee for the reading and the pertinent comments and questions. A.W. thanks J. Alonso-Floriano for discussions about transmission spectroscopy, C. Aerts for discussions about stellar pulsations, and A. Vidotto, A. Garcìa-Muñoz, M. Bergemann, and Y. Miguel for discussions about NLTE in hot Jupiter atmospheres.

Appendix A Rossiter-McLaughlin MCMC results

Here, we show the correlation diagrams for the reloaded Rossiter-McLaughlin study obtained with the emcee algorithm.

thumbnail Fig. A.1

Correlation diagrams in the case of a stellar solid body rotation.

thumbnail Fig. A.2

Correlation diagrams for the case of a stellar differential rotation, with an “anti-solar” law.

thumbnail Fig. A.3

Correlation diagrams for the case of a stellar differential rotation, with a solar-like law.

Appendix B: Rossiter-McLaughlin correction for transmission spectroscopy

Here, we show the effect of the Rossiter-McLaughlin correction on the Hα line presented in Sect. 3.3.

thumbnail Fig. B.1

Hα transmission spectrum in the planet rest frame (gray line and black circles), uncorrected (top) and corrected (bottom) from the Rossiter-McLaughlin effect. The vertical black dashed line shows the line center taking the systemic velocity into account. We show the HβoxGaussian fit (blue) to the corrected line in both sub-panels for the purposes of comparison.

Appendix C Transmission spectroscopy MCMC results

Here, we show complementary correlation diagrams for the Balmer lines NLTE study presented in Sect. 4.2.2.

thumbnail Fig. C.1

Correlation diagrams in the case of a complete NLTE atmosphere (each neutral hydrogen electronic level is a separated free parameter). Two families of solution are visible, one being compatible and not so far from the Boltzmann equilibrium (shown in red dashed lines), the other being significantly away from the Boltzmann equilibrium. These two families are studied in the main text.

thumbnail Fig. C.2

Correlation diagrams in the case of a complete NLTE atmosphere (each neutral hydrogen electronic level is a separated free parameter). Exploration of the high temperature, near Boltzmann equilibrium (red dashed lines) family solution (see main text).

References

  1. Aerts, C., de Pauw, M., & Waelkens, C. 1992, A&A, 266, 294 [NASA ADS] [Google Scholar]
  2. Aerts, C., Christensen-Dalsgaard, J., & Kurtz, D. W. 2010, Asteroseismology (Berlin: Springer) [Google Scholar]
  3. Ahlers, J. P., Barnes, J. W., & Barnes, R. 2015, ApJ, 814, 67 [NASA ADS] [CrossRef] [Google Scholar]
  4. Ahlers, J. P., Kruse, E., Colón, K. D., et al. 2020, ApJ, 888, 63 [Google Scholar]
  5. Allan, A., & Vidotto, A. A. 2019, MNRAS, 490, 3760 [Google Scholar]
  6. Allart, R., Lovis, C., Pino, L., et al. 2017, A&A, 606, A144 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  7. Allart, R., Bourrier, V., Lovis, C., et al. 2018, Science, 362, 1384 [NASA ADS] [CrossRef] [Google Scholar]
  8. Allart, R., Bourrier, V., Lovis, C., et al. 2019, A&A, 623, A58 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  9. Alonso-Floriano, F. J., Snellen, I. A. G., Czesla, S., et al. 2019, A&A, 629, A110 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  10. Anderson, D. R., Temple, L. Y., Nielsen, L. D., et al. 2018, MNRAS, submitted [arXiv:1809.04897] [Google Scholar]
  11. Arcangeli, J., Désert, J.-M., Line, M. R., et al. 2018, ApJ, 855, L30 [NASA ADS] [CrossRef] [Google Scholar]
  12. Asplund, M., Grevesse, N., Sauval, A. J., & Scott, P. 2009, ARA&A, 47, 481 [NASA ADS] [CrossRef] [Google Scholar]
  13. Balona, L. A. 1986, MNRAS, 219, 111 [Google Scholar]
  14. Balona, L. A., & Abedigamba, O. P. 2016, MNRAS, 461, 497 [NASA ADS] [CrossRef] [Google Scholar]
  15. Barman, T. S., Hauschildt, P. H., Schweitzer, A., et al. 2002, ApJ, 569, L51 [NASA ADS] [CrossRef] [Google Scholar]
  16. Barnes, J. R., Haswell, C. A., Staab, D., & Anglada-Escudé, G. 2016, MNRAS, 462, 1012 [Google Scholar]
  17. Barnes, J. W., Linscott, E., & Shporer, A. 2011, ApJS, 197, 10 [NASA ADS] [CrossRef] [Google Scholar]
  18. Bell, T. J., & Cowan, N. B. 2018, ApJ, 857, L20 [NASA ADS] [CrossRef] [Google Scholar]
  19. Borsa, F., Rainer, M., Bonomo, A. S., et al. 2019, A&A, 631, A34 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  20. Bourrier, V., Ehrenreich, D., & Lecavelier des Etangs, A. 2015, A&A, 582, A65 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  21. Bourrier, V., Lecavelier des Etangs, A., Ehrenreich, D., Tanaka, Y. A., & Vidotto, A. A. 2016, A&A, 591, A121 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  22. Bourrier, V., Cegla, H. M., Lovis, C., & Wyttenbach, A. 2017a, A&A, 599, A33 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  23. Bourrier, V., Ehrenreich, D., Allart, R., et al. 2017b, A&A, 602, A106 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  24. Bourrier, V., Lovis, C., Beust, H., et al. 2018, Nature, 553, 477 [NASA ADS] [CrossRef] [Google Scholar]
  25. Bourrier, V., Ehrenreich, D., Lendl, M., et al. 2020a, A&A, 635, A205 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  26. Bourrier, V., Kitzmann, D., Kuntzer, T., et al. 2020b, A&A, 637, A36 [CrossRef] [EDP Sciences] [Google Scholar]
  27. Brogi, M., de Kok, R. J., Albrecht, S., et al. 2016, ApJ, 817, 106 [NASA ADS] [CrossRef] [Google Scholar]
  28. Brown, T. M. 2001, ApJ, 553, 1006 [NASA ADS] [CrossRef] [Google Scholar]
  29. Cabot, S. H. C., Madhusudhan, N., Welbanks, L., Piette, A., & Gandhi, S. 2020, MNRAS, 494, 363 [NASA ADS] [CrossRef] [Google Scholar]
  30. Casasayas-Barris, N., Palle, E., Nowak, G., et al. 2017, A&A, 608, A135 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  31. Casasayas-Barris, N., Pallé, E., Yan, F., et al. 2018, A&A, 616, A151 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  32. Casasayas-Barris, N., Pallé, E., Yan, F., et al. 2019, A&A, 628, A9 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  33. Cauley, P. W., Redfield, S., Jensen, A. G., et al. 2015, ApJ, 810, 13 [NASA ADS] [CrossRef] [Google Scholar]
  34. Cauley, P. W., Redfield, S., Jensen, A. G., & Barman, T. 2016, AJ, 152, 20 [NASA ADS] [CrossRef] [Google Scholar]
  35. Cauley, P. W., Redfield, S., & Jensen, A. G. 2017a, AJ, 153, 217 [NASA ADS] [CrossRef] [Google Scholar]
  36. Cauley, P. W., Redfield, S., & Jensen, A. G. 2017b, AJ, 153, 81 [Google Scholar]
  37. Cauley, P. W., Redfield, S., & Jensen, A. G. 2017c, AJ, 153, 185 [NASA ADS] [CrossRef] [Google Scholar]
  38. Cauley, P. W., Shkolnik, E. L., Ilyin, I., et al. 2019, AJ, 157, 69 [NASA ADS] [CrossRef] [Google Scholar]
  39. Cegla, H. M., Lovis, C., Bourrier, V., et al. 2016, A&A, 588, A127 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  40. Chamberlain, J. W., & Hunten, D. M. 1987, Theory of Planetary Atmospheres. An Introduction to Their Physics and Chemistry (Cambridge: Academic Press), 36 [Google Scholar]
  41. Charbonneau, D., Brown, T. M., Noyes, R. W., & Gilliland, R. L. 2002, ApJ, 568, 377 [NASA ADS] [CrossRef] [Google Scholar]
  42. Chen, G., Casasayas-Barris, N., Pallé, E., et al. 2020, A&A, 635, A171 [CrossRef] [Google Scholar]
  43. Christie, D., Arras, P., & Li, Z.-Y. 2013, ApJ, 772, 144 [NASA ADS] [CrossRef] [Google Scholar]
  44. Collier Cameron, A., Bruce, V. A., Miller, G. R. M., Triaud, A. H. M. J., & Queloz, D. 2010, MNRAS, 403, 151 [NASA ADS] [CrossRef] [Google Scholar]
  45. Cosentino, R., Lovis, C., Pepe, F., et al. 2012, Proc. SPIE, 8446, 84461V [Google Scholar]
  46. Czesla, S., Schröter, S., Wolter, U., et al. 2012, A&A, 539, A150 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  47. Czesla, S., Klocová, T., Khalafinejad, S., Wolter, U., & Schmitt, J. H. M. M. 2015, A&A, 582, A51 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  48. de Wit, J., & Seager, S. 2013, Science, 342, 1473 [NASA ADS] [CrossRef] [Google Scholar]
  49. Dorval, P., Talens, G. J. J., Otten, G. P. P. L., et al. 2020, A&A, 635, A60 [Google Scholar]
  50. Draine, B. T. 2011, Physics of the Interstellar and Intergalactic Medium (Princeton: Princeton University Press) [Google Scholar]
  51. Dravins, D., Ludwig, H.-G., Dahlen, E., & Pazira, H. 2015, in Cambridge Workshop on Cool Stars, Stellar Systems, and the Sun, 18th Cambridge Workshop on Cool Stars, Stellar Systems, and the Sun, 18, 853 [NASA ADS] [Google Scholar]
  52. Dravins, D., Ludwig, H.-G., Dahlén, E., & Pazira, H. 2017, A&A, 605, A90 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  53. Dreizler, S., Reiners, A., Homeier, D., & Noll, M. 2009, A&A, 499, 615 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  54. Dumusque, X., Boisse, I., & Santos, N. C. 2014, ApJ, 796, 132 [NASA ADS] [CrossRef] [Google Scholar]
  55. Eggleton, P. P. 1983, ApJ, 268, 368 [NASA ADS] [CrossRef] [Google Scholar]
  56. Ehrenreich, D., Bourrier, V., Wheatley, P. J., et al. 2015, Nature, 522, 459 [Google Scholar]
  57. Ehrenreich, D., Lovis, C., Allart, R., et al. 2020, Nature, 580, 597 [CrossRef] [Google Scholar]
  58. Ehrenreich, D., Tinetti, G., Lecavelier Des Etangs, A., Vidal-Madjar, A., & Selsis, F. 2006, A&A, 448, 379 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  59. Fisher, C., & Heng, K. 2019, ApJ, 881, 25 [NASA ADS] [CrossRef] [Google Scholar]
  60. Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306 [CrossRef] [Google Scholar]
  61. Fortney, J. J., Sudarsky, D., Hubeny, I., et al. 2003, ApJ, 589, 615 [Google Scholar]
  62. Fossati, L., Koskinen, T., Lothringer, J. D., et al. 2018, ApJ, 868, L30 [NASA ADS] [CrossRef] [Google Scholar]
  63. García Muñoz, A. 2007, Planet. Space Sci., 55, 1426 [NASA ADS] [CrossRef] [Google Scholar]
  64. García Muñoz, A., & Schneider, P. C. 2019, ApJ, 884, L43 [CrossRef] [Google Scholar]
  65. Gaudi, B. S., Stassun, K. G., Collins, K. A., et al. 2017, Nature, 546, 514 [NASA ADS] [Google Scholar]
  66. Georgy, C., Granada, A., Ekström, S., et al. 2014, A&A, 566, A21 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  67. Gibson, N. P., Merritt, S., Nugroho, S. K., et al. 2020, MNRAS, 493, 2215 [NASA ADS] [CrossRef] [Google Scholar]
  68. Gray, D. F. 2005, The Observation and Analysis of Stellar Photospheres (Cambridge: Cambridge University Press) [CrossRef] [Google Scholar]
  69. Helling, C., Gourbin, P., Woitke, P., & Parmentier, V. 2019, A&A, 626, A133 [CrossRef] [Google Scholar]
  70. Heng, K. 2017, Exoplanetary Atmospheres: Theoretical Concepts and Foundations (Princeton: Princeton University Press) [Google Scholar]
  71. Heng, K., & Kitzmann, D. 2017, MNRAS, 470, 2972 [NASA ADS] [CrossRef] [Google Scholar]
  72. Heng, K., Wyttenbach, A., Lavie, B., et al. 2015, ApJ, 803, L9 [Google Scholar]
  73. Hoeijmakers, H. J., Ehrenreich, D., Heng, K., et al. 2018, Nature, 560, 453 [NASA ADS] [CrossRef] [Google Scholar]
  74. Hoeijmakers, H. J., Ehrenreich, D., Kitzmann, D., et al. 2019, A&A, 627, A165 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  75. Huang, C., Arras, P., Christie, D., & Li, Z.-Y. 2017, ApJ, 851, 150 [NASA ADS] [CrossRef] [Google Scholar]
  76. Hubeny, I., & Mihalas, D. 2014, Theory of Stellar Atmospheres (Princeton: Princeton University Press) [Google Scholar]
  77. Hunten, D. M., Donahue, T. M., Walker, J. C. G., & Kasting, J. F. 1989, Escape of Atmospheres and Loss of Water., eds. S. K. Atreya, J. B. Pollack, & M. S. Matthews (Tucson: University of Arizona Press), 386 [Google Scholar]
  78. Jensen, A. G., Redfield, S., Endl, M., et al. 2012, ApJ, 751, 86 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  79. Jensen, A. G., Cauley, P. W., Redfield, S., Cochran, W. D., & Endl, M. 2018, AJ, 156, 154 [NASA ADS] [CrossRef] [Google Scholar]
  80. Kausch, W., Noll, S., Smette, A., et al. 2015, A&A, 576, A78 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  81. Keles, E., Mallonn, M., von Essen, C., et al. 2019, MNRAS, 489, L37 [NASA ADS] [CrossRef] [Google Scholar]
  82. Khalafinejad, S., Salz, M., Cubillos, P. E., et al. 2018, A&A, 618, A98 [Google Scholar]
  83. Khalafinejad, S., von Essen, C., Hoeijmakers, H. J., et al. 2017, A&A, 598, A131 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  84. Kitzmann, D., Heng, K., Rimmer, P. B., et al. 2018, ApJ, 863, 183 [NASA ADS] [CrossRef] [Google Scholar]
  85. Koll, D. D. B., & Komacek, T. D. 2018, ApJ, 853, 133 [NASA ADS] [CrossRef] [Google Scholar]
  86. Koskinen, T. T., Harris, M. J., Yelle, R. V., & Lavvas, P. 2013, Icarus, 226, 1678 [Google Scholar]
  87. Kreidberg, L. 2015, PASP, 127, 1161 [NASA ADS] [CrossRef] [Google Scholar]
  88. Kreidberg, L., Line, M. R., Parmentier, V., et al. 2018, AJ, 156, 17 [NASA ADS] [CrossRef] [Google Scholar]
  89. Kurucz, R. L. 1979, ApJS, 40, 1 [NASA ADS] [CrossRef] [Google Scholar]
  90. Kurucz, R. L. 1992, IAU Symp., 149, 225 [Google Scholar]
  91. Lamers, H. J. G. L. M., & Cassinelli, J. P. 1999, Introduction to Stellar Winds (Cambridge: Cambridge University Press) [CrossRef] [Google Scholar]
  92. Lammer, H., Selsis, F., Ribas, I., et al. 2003, ApJ, 598, L121 [NASA ADS] [CrossRef] [Google Scholar]
  93. Lampón, M., López-Puertas, M., Lara, L. M., et al. 2020, A&A, 636, A13 [CrossRef] [EDP Sciences] [Google Scholar]
  94. Lavvas, P., Koskinen, T., & Yelle, R. V. 2014, ApJ, 796, 15 [NASA ADS] [CrossRef] [Google Scholar]
  95. Lecavelier des Etangs, A., Vidal-Madjar, A., McConnell, J. C., & Hébrard, G. 2004, A&A, 418, L1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  96. Lecavelier des Etangs, A., Pont, F., Vidal-Madjar, A., & Sing, D. 2008, A&A, 481, L83 [NASA ADS] [CrossRef] [Google Scholar]
  97. Lothringer, J. D., & Barman, T. 2019, ApJ, 876, 69 [NASA ADS] [CrossRef] [Google Scholar]
  98. Lothringer, J. D., Barman, T., & Koskinen, T. 2018, ApJ, 866, 27 [NASA ADS] [CrossRef] [Google Scholar]
  99. Madhusudhan, N., Knutson, H., Fortney, J. J., & Barman, T. 2014, Protostars and Planets VI (Tucson: University of Arizona Press), 739 [Google Scholar]
  100. Madhusudhan, N., & Redfield, S. 2015, Int. J. Astrobiol., 14, 177 [Google Scholar]
  101. Maeder, A. 2009, Physics, Formation and Evolution of Rotating Stars (Berlin: Springer) [Google Scholar]
  102. Mandel, K., & Agol, E. 2002, ApJ, 580, L171 [NASA ADS] [CrossRef] [Google Scholar]
  103. Mansfield, M., Bean, J. L., Line, M. R., et al. 2018, AJ, 156, 10 [NASA ADS] [CrossRef] [Google Scholar]
  104. Mansfield, M., Bean, J. L., Stevenson, K. B., et al. 2020, ApJ, 888, L15 [NASA ADS] [CrossRef] [Google Scholar]
  105. Mayor, M., & Queloz, D. 1995, Nature, 378, 355 [Google Scholar]
  106. Miguel, Y., & Kaltenegger, L. 2014, ApJ, 780, 166 [NASA ADS] [CrossRef] [Google Scholar]
  107. Mihalas, D. 1970, Stellar Atmospheres (Wales: W.H. Freeman & Co Ltd;) [Google Scholar]
  108. Mollière, P., van Boekel, R., Bouwman, J., et al. 2017, A&A, 600, A10 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  109. Mollière, P., Wardenier, J. P., van Boekel, R., et al. 2019, A&A, 627, A67 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  110. Murray-Clay, R. A., Chiang, E. I., & Murray, N. 2009, ApJ, 693, 23 [Google Scholar]
  111. Nortmann, L., Pallé, E., Salz, M., et al. 2018, Science, 362, 1388 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  112. Nugroho, S. K., Gibson, N. P., de Mooij, E. J. W., et al. 2020, MNRAS, submitted [arXiv:2003.04856] [Google Scholar]
  113. Oklopčić, A., & Hirata, C. M. 2018, ApJ, 855, L11 [NASA ADS] [CrossRef] [Google Scholar]
  114. Oklopčić, A. 2019, ApJ, 881, 133 [NASA ADS] [CrossRef] [Google Scholar]
  115. Parker, E. N. 1958, ApJ, 128, 664 [Google Scholar]
  116. Parmentier, V., Line, M. R., Bean, J. L., et al. 2018, A&A, 617, A110 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  117. Pepe, F., Mayor, M., Galland, F., et al. 2002, A&A, 388, 632 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  118. Pepe, F., Cristiani, S., Rebolo, R., et al. 2013, The Messenger, 153, 6 [NASA ADS] [Google Scholar]
  119. Pepe, F., Ehrenreich, D., & Meyer, M. R. 2014, Nature, 513, 358 [NASA ADS] [CrossRef] [Google Scholar]
  120. Pino, L., Ehrenreich, D., Wyttenbach, A., et al. 2018, A&A, 612, A53 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  121. Pino, L., Désert, J. M., Brogi, M., et al. 2020, ApJ, 894, L27 [NASA ADS] [CrossRef] [Google Scholar]
  122. Queloz, D., Eggenberger, A., Mayor, M., et al. 2000, A&A, 359, L13 [NASA ADS] [Google Scholar]
  123. Redfield, S., Endl, M., Cochran, W. D., & Koesterke, L. 2008, ApJ, 673, L87 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
  124. Reiners, A., & Royer, F. 2004, A&A, 415, 325 [Google Scholar]
  125. Ridden-Harper, A. R., Snellen, I. A. G., Keller, C. U., et al. 2016, A&A, 593, A129 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  126. Salem, M., & Brocklehurst, M. 1979, ApJS, 39, 633 [Google Scholar]
  127. Salz, M., Czesla, S., Schneider, P. C., & Schmitt, J. H. M. M. 2016, A&A, 586, A75 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  128. Salz, M., Czesla, S., Schneider, P. C., et al. 2018, A&A, 620, A97 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  129. Seager, S., & Sasselov, D. D. 2000, ApJ, 537, 916 [NASA ADS] [CrossRef] [Google Scholar]
  130. Seidel, J. V., Ehrenreich, D., Wyttenbach, A., et al. 2019, A&A, 623, A166 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  131. Seidel, J. V., Ehrenreich, D., Pino, L., et al. 2020, A&A, 633, A86 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  132. Sharp, C. M., & Burrows, A. 2007, ApJS, 168, 140 [NASA ADS] [CrossRef] [Google Scholar]
  133. Shizgal, B. D., & Arkos, G. G. 1996, Rev. Geophys., 34, 483 [Google Scholar]
  134. Sing, D. K., Vidal-Madjar, A., Lecavelier des Etangs, A., et al. 2008, ApJ, 686, 667 [NASA ADS] [CrossRef] [Google Scholar]
  135. Sing, D. K., Lavvas, P., Ballester, G. E., et al. 2019, AJ, 158, 91 [NASA ADS] [CrossRef] [Google Scholar]
  136. Smette, A., Sana, H., Noll, S., et al. 2015, A&A, 576, A77 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  137. Snellen, I. A. G. 2004, MNRAS, 353, L1 [NASA ADS] [CrossRef] [Google Scholar]
  138. Snellen, I. A. G., Albrecht, S., de Mooij, E. J. W., & Le Poole, R. S. 2008, A&A, 487, 357 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  139. Snellen, I. A. G., de Kok, R. J., de Mooij, E. J. W., & Albrecht, S. 2010, Nature, 465, 1049 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  140. Spake, J. J., Sing, D. K., Evans, T. M., et al. 2018, Nature, 557, 68 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  141. Stangret, M., Casasayas-Barris, N., Pallé, E., et al. 2020, A&A, 638, A26 [CrossRef] [EDP Sciences] [Google Scholar]
  142. Tian, F., Toon, O. B., Pavlov, A. A., & De Sterck, H. 2005, ApJ, 621, 1049 [NASA ADS] [CrossRef] [Google Scholar]
  143. Turner, J. D., de Mooij, E. J. W., Jayawardhana, R., et al. 2020, ApJ, 888, L13 [NASA ADS] [CrossRef] [Google Scholar]
  144. Vidal-Madjar, A., Lecavelier des Etangs, A., Désert, J.-M., et al. 2003, Nature, 422, 143 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  145. Vidal-Madjar, A., Désert, J.-M., Lecavelier des Etangs, A., et al. 2004, ApJ, 604, L69 [NASA ADS] [CrossRef] [Google Scholar]
  146. Vidal-Madjar, A., Sing, D. K., Lecavelier Des Etangs, A., et al. 2011, A&A, 527, A110 [Google Scholar]
  147. Wong, I., Shporer, A., Morris, B. M., et al. 2019, AJ, submitted [arXiv:1910.01607] [Google Scholar]
  148. Wyttenbach, A., Ehrenreich, D., Lovis, C., Udry, S., & Pepe, F. 2015, A&A, 577, A62 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  149. Wyttenbach, A., Lovis, C., Ehrenreich, D., et al. 2017, A&A, 602, A36 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  150. Yan, F., & Henning, T. 2018, Nat. Astron., 2, 714 [NASA ADS] [CrossRef] [Google Scholar]
  151. Yan, F., Pallé, E., Fosbury, R. A. E., Petr-Gotzens, M. G., & Henning, T. 2017, A&A, 603, A73 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  152. Yan, F., Casasayas-Barris, N., Molaverdikhani, K., et al. 2019, A&A, 632, A69 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  153. Yelle, R. V. 2004, Icarus, 170, 167 [NASA ADS] [CrossRef] [Google Scholar]
  154. Zhou, G., Bakos, G. Á., Bayliss, D., et al. 2019, AJ, 157, 31 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
  155. Zorec, J., & Royer, F. 2012, A&A, 537, A120 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  156. Zorec, J., Rieutord, M., Espinosa Lara, F., et al. 2017, A&A, 606, A32 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

1

This data set is the same as in Hoeijmakers et al. (2018, 2019). Those previous studies focused on metal detections (species between atomic numbers 3 and 78), while our present study focuses only on the hydrogen lines, making it a complementary analysis. We note that no Helium (λ10 830 Å) has been found in Nortmann et al. (2018).

2

Our observational strategy of recording several transits – to assess the reproducibility of the detected spectroscopic signatures – was made possible thanks to a Director’s Discretionary Time proposal asked for in summer 2017 (A35DDT4) and a standard OPTICON proposal in the general frame of our SPADES survey (2018A/038, PI: Ehrenreich).

3

We discarded (by visual inspection) the order number 1, 2, 50–58, 61, 64–66, 68 and 69 because they contain no or very few stellar lines, or because they are strongly affected by interstellar medium (ISM) or telluric line absorptions.

4

We ignore any refraction in the exoplanetary atmosphere.

5

We define the planetary absorption A to be A = str∕Σatm.

6

This has little impact since the expected relative change of the line contrasts would be less than 0.5% (from the prediction of Fig. 14 of Yan et al. 2019), that is smaller than our individual error bars.

7

In practice, Eq. (4) is used with a weighted mean that take the changing data quality into account. The weights w at time t are the inverse of the squared flux uncertainties, .

8

To compare our FWHM, we removed the additional width contribution by other effects (see Sects. 4.1 and 4.2.1).

9

A composition dominated by molecular hydrogen and helium (μ = 2.34) is probably excluded by an unphysical T ≳ 30 000 K.

10

Equation (13.14) of Heng (2017, Chap. 13) should be corrected as follows: where the Mach number , and x = rrs, with cs, rs being the sound speed, and the sonic point radius, respectively.

All Tables

Table 1

Adopted values for the orbital and physical parameters of the KELT-9 system.

Table 2

Log of HARPS-N observations.

Table 3

Summary of the hydrogen Balmer line detections.

All Figures

thumbnail Fig. 1

Distribution of planetary masses as a function of the total power received by the planet from the incident bolometric stellar radiation, expressed in units of power received on Earth at the top of the atmosphere (P = 174 PW = 1.74 × 1024 erg s−1). The markersize scales with the planet density, and the color with the planetary equilibrium temperature. Physical properties of exoplanets were extracted from the Extrasolar Planets Encyclopedia (exoplanet.eu) in December 2019. The lower left corner represents the evaporation desert. KELT-9 (HD 195689), encircled in black, stand out of all other planetary population in term of its received power and equilibrium temperature.

In the text
thumbnail Fig. 2

Left panel: normalized out-of-transit master spectrum of KELT-9 in the 387–687 nm optical region. The stellar A0 spectral type makes the spectrum dominated by strong Balmer lines, though we noticethe presence of numerous shallow metal lines, as shown in the inset. The Balmer lines (Hα, β, γ, δ, ɛ, ζ), that are covered by the HARPS-N spectrograph wavelength range, are identified by the vertical dashed lines. Right panel: normalized out-of-transit master cross-correlation function (CCF) obtained with a dedicated A0 mask. The averaged stellar line has a measured projected velocity broadening of v sin i* = 112.6 km s−1. The vertical dashed line represent the systemic velocity (−17.74 km s−1).

In the text
thumbnail Fig. 3

Schematic of a spectroscopic transit observation of an exoplanetary system. The difference between an out-of-transit spectrum (left) and an in-transit spectrum (middle) gives us two essential quantities: the exoplanetary atmosphere transmission spectrum and the local stellar spectrum of the surface that is hidden behind the exoplanet (right). The temporal evolution of the stellar local spectrum will mostly depend on the projected stellar rotation (represented by the black arrow going through the star from pole to pole), the planetary orbit parameters (represented with the green arrows) and particularly the projected spin-orbit alignment angle. The temporal evolution of the exoplanet atmosphere spectrum will mostly depend on the atmosphere size, the planetary orbital projected velocity, and the local stellar spectrum that is transmitted through. Hence, it is necessary to separate and to measure accurately each component to retrieve the unknown parameters. The present figure is partially inspired by figures in Dravins et al. (2015, 2017); Dumusque et al. (2014); Cegla et al. (2016); Bourrier et al. (2017a).

In the text
thumbnail Fig. 4

Map of the residual CCF flux during the transit of KELT-9 b. The observations of two nights are binned by 2 in phase. The 1st and 4th contacts are traced with the horizontal dashed white lines; the 0 km s−1 and ± v sin i* are traced with vertical dashed white lines. Top panel: local stellar spectrum tracing the spin-orbit misalignment through the Rossiter-McLaughlin effect. The solid black line shows our best fit. The exoplanet transmission CCF, that has a slanted shape following the orbital velocity (dashed black line), traces the metal detections by Hoeijmakers et al. (2018, 2019). Bottom panel: residual CCF flux after the removal of the local stellar and of the planetary atmosphere spectra. The residuals are dominated by the stellar pulsations, whose strengths are 10 × lower than the local stellar signal (see the color scale). The pulsations are due to pressure-mode in a δ Scuti-type star.

In the text
thumbnail Fig. 5

Stellar pulsation signal seen from variation in the stellar v sin i*. The stellar CCF are cleaned from the local stellar and planetary atmosphere spectra, and fit with a stellar rotation profile. The pulsation period is Ppuls = 7.54 ± 0.12 h, and is consistent between the two nights (in dark blue and green, respectively) that are separated by one year. The averaged v sin i* is reported for each night with the horizontal dashed lines, and has an averaged value of v sin i* = 112.60 ± 0.04 km s−1.

In the text
thumbnail Fig. 6

Velocities of the local stellar surface hidden behind the planet during its transit. The observations of the two nights are shown in dark blue and green, respectively. The first to fourth contacts are traced with the vertical dashed lines. The overlap region between the local stellar and atmospheric signals (between phases −0.025 and 0) is shown in red and is not used for the fit. The spin-orbit misalignment is modeled with a solid body rotation (in orange), and with a differential rotation (in light blue). The differential rotation model is slightly preferred.

In the text
thumbnail Fig. 7

Left: detection of the planetary Hα line (significance of 24.3σ). Top left panel: map of the residual spectra ±200 km s−1 around the Hα line in the stellar rest frame (Eq. (1)). The observations of two nights are binned by 2 in phase. The planetary atmosphere Hα absorption follows the orbital velocity (solid black line). The local Hα stellar line follows the reloaded RM model (dashed black line). The first and fourth contacts are traced with the horizontal dashed white lines, the 0 km s−1 and ± v sin i* are traced with vertical dashed white lines. Bottom left panel: weighted averaged Hα transmission spectrum (corrected from the local stellar photospheric signal, Eq. (4)) in the planet rest frame (light gray), also binned by 15× in black circles (the continuum is set to the white light curve transit depth δ = 0.00677). The verticalblack dashed line shows the line center taking the systemic velocity into account. We show the Gaussian fit (blue) with contrast of 0.91 ± 0.04% and FWHM of 44.3 ± 1.8 km s−1. Right: same as for Hα, but for the Hβ line (detection significance of 16.1σ, see Table 3).

In the text
thumbnail Fig. 8

Same as for Fig. 7, but for the Hγ line (left) and the Hδ line (right). The Hγ and Hδ lines are detected with a significance of 7.5σ and 4.6σ, respectively (see Table 3).

In the text
thumbnail Fig. 9

Same as for Fig. 7, but for the Hɛ line (left), and the Hζ line (right). The Hɛ and Hζ lines are not statistically detected (see Table 3).

In the text
thumbnail Fig. 10

Correlation diagrams in the case of a hydrostatic atmosphere in the Boltzmann equilibrium.

In the text
thumbnail Fig. 11

Correlation diagrams in the case of a hydrodynamic atmosphere in the Boltzmann equilibrium.

In the text
thumbnail Fig. 12

Correlation diagrams in the case of a hydrodynamic atmosphere with a departure β from LTE.

In the text
thumbnail Fig. 13

Hα, Hβ, and Hγ transmission spectra (light gray, and binned by 10 × in black circles) in the planet rest frame and in velocity space (the continuums are set to the white light transit depth δ = 0.00677). The best PAWN models, in the case of a hydrodynamically expanding thermosphere, and in the Boltzmann equilibrium, are shown in blue. The retrieved temperature is K, and the mass loss rate is g s−1. The best PAWN models, in the case of a hydrodynamically expanding thermosphere, with a LTE departure, are shown in dashed red. The retrieved temperature is K, the mass lossrate is g s−1, and the common LTE departure coefficient is .

In the text
thumbnail Fig. 14

Summary of the mass loss rates and their errors found in the LTE case compared to the expected mass loss rates from the theoretical studies of Fossati et al. (2018) and García Muñoz & Schneider (2019). Our finding shows that the KELT-9 b thermosphere probably undergoes a “Balmer-driven” expansion, that leads to a high escape rate. In this case, the total mass loss through the planetary life time (500 My) represents sim 0.05 MJup, which is sim 2% of the planet mass.

In the text
thumbnail Fig. 15

Example of PAWN models of the Ca II H line detection by Yan et al. (2019). Our model is compatible with their detection only when assuming some LTE departure. For the Ca II, the departure is likely explained by further thermal and/or photo-ionization not yet taken in our code into account. Additional data gathering, and modeling will be necessary to confirm our interpretation.

In the text
thumbnail Fig. A.1

Correlation diagrams in the case of a stellar solid body rotation.

In the text
thumbnail Fig. A.2

Correlation diagrams for the case of a stellar differential rotation, with an “anti-solar” law.

In the text
thumbnail Fig. A.3

Correlation diagrams for the case of a stellar differential rotation, with a solar-like law.

In the text
thumbnail Fig. B.1

Hα transmission spectrum in the planet rest frame (gray line and black circles), uncorrected (top) and corrected (bottom) from the Rossiter-McLaughlin effect. The vertical black dashed line shows the line center taking the systemic velocity into account. We show the HβoxGaussian fit (blue) to the corrected line in both sub-panels for the purposes of comparison.

In the text
thumbnail Fig. C.1

Correlation diagrams in the case of a complete NLTE atmosphere (each neutral hydrogen electronic level is a separated free parameter). Two families of solution are visible, one being compatible and not so far from the Boltzmann equilibrium (shown in red dashed lines), the other being significantly away from the Boltzmann equilibrium. These two families are studied in the main text.

In the text
thumbnail Fig. C.2

Correlation diagrams in the case of a complete NLTE atmosphere (each neutral hydrogen electronic level is a separated free parameter). Exploration of the high temperature, near Boltzmann equilibrium (red dashed lines) family solution (see main text).

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.