Atmospheric characterization of the ultra-hot Jupiter MASCARA-2b/KELT-20b

TL;DR

This study uses high-resolution spectroscopy to analyze the atmosphere of the ultra-hot Jupiter MASCARA-2b/KELT-20b, detecting multiple atomic species and confirming ionized gas presence consistent with theoretical models.

Contribution

First detection of CaII triplet lines in transmission spectrum of an exoplanet, providing new insights into ultra-hot Jupiter atmospheric composition.

Findings

Detected CaII, FeII, NaI, Hα, Hβ in the atmosphere

Measured consistent absorption depths across instruments

Confirmed presence of ionized gases as predicted by models

Abstract

Ultra hot Jupiters orbit very close to their host star and consequently receive strong irradiation that makes their atmospheric chemistry different from the common gas giants. Here, we study the atmosphere of one of these particular hot planets, MASCARA-2b/KELT-20b, using four transit observations with high resolution spectroscopy facilities. Three of these observations were performed with HARPS-N and one with CARMENES. We simultaneously observed one of the transits with MuSCAT2 to monitor possible spots in the stellar surface. At high resolution, the transmission residuals show the effects of Rossiter-McLaughlin and center-to-limb variations from the stellar lines profiles, which we correct to finally extract the transmission spectra of the planet. We clearly observe absorption of CaII, FeII, NaI, H$\alpha$ and H$\beta$ in the atmosphere of MASCARA-2b, and indications of H$\gamma$ and MgI at low SNR. The results obtained with CARMENES and HARPS-N are consistent, measuring an H$\alpha$ absorption depth of $0.68\pm0.05\%$ and $0.59\pm0.07\%$, and NaI absorption of $0.11\pm0.04\%$ and $0.09\pm0.05\%$ ($0.75\r{A}$ passband), in both instruments respectively. For H$\beta$ only HARPS-N covers this wavelength range, measuring $0.28\pm0.06\%$ of absorption. The CaII is only covered with CARMENES, measuring an average absorption of $0.32\pm0.05\%$. Three additional FeII lines are observed in the transmission spectrum by HARPS-N, with mean absorption depth of $0.08\pm0.04\%$. The results presented here are consistent with theoretical models of ultra hot Jupiters atmospheres, suggesting the emergence of an ionised gas on the day-side of such planets. Ca and Fe, with other elements, are expected to be singly ionised at these temperatures and be more numerous than its neutral state. The CaII triplet lines are detected here for the first time in transmission in an exoplanet atmosphere.

Tables

Table 1. Table 1: Physical and orbital parameters of MASCARA-2b.

Description	Symbol	Value
Stellar parameters
Identifiers	-	KELT-20, HD 185603
V-band magnitude	$m_V$	$7.6$
Effective temperature	$T_{\mathrm{eff}}$	${8980}_{-130}^{+90}K$
Projected rotation speed	$v\sin i_{\star }$	$114\pm 3\mathrm{km}{s}^{-1}$
Surface gravity	$\log g$	$4.31\pm 0.02\mathrm{cgs}$
Metallicity	[Fe/H]	$-0.02\pm 0.07$
Stellar mass	$M_{\star }$	${1.89}_{-0.05}^{+0.06}{M}_{\odot }$
Stellar radius	$R_{\star }$	$1.60\pm 0.06R_{\odot }$
Planet parameters
Planet mass(11footnotemark: 1)	$M_p$
Planet radius	$R_p$	$1.83\pm 0.07R_{\mathrm{Jup}}$
Equilibrium temperature	$T_{eq}$	$2260\pm 50K$
Surface gravity(11footnotemark: 1)	$\log g_P$
System parameters
Right ascension	-	19h38m38.73s
Declination	-	+31o13’09.2”
Epoch	$T_c$	${2457909.5906}_{-0.0002}^{+0.0003}\mathrm{BJD}$
Period(11footnotemark: 1)	$P$	$3.4741070\pm 0.0000019\mathrm{days}$
Transit duration(11footnotemark: 1)	$T_{14}$	${0.14882}_{-0.00090}^{+0.00092}\mathrm{days}$
Ingress/Egress duration(11footnotemark: 1)	$\tau$	${0.01985}_{-0.00079}^{+0.00082}\mathrm{days}$
Semimajor axis(11footnotemark: 1)	$a$	${0.0542}_{-0.0021}^{+0.0014}\mathrm{AU}$
Inclination(11footnotemark: 1)	$i$	${86.15}_{-0.27}^{+0.28}\mathrm{deg}$
Eccentricity	e	$0$ (fixed)
Systemic velocity	$\gamma$	$-21.07\pm 0.03\mathrm{km}{s}^{-1}$
Projected obliquity	$\lambda$	$0.6\pm 4\mathrm{deg}$

Table 2. Table 2: Observing log of the MASCARA-2b transit observations.

Night	Tel.	Instrument	Date of	Start	End	Filter	$T_{\exp }$	$N_{\mathrm{obs}}$	Airmass	Aperture
			observation	UT	UT		[s]			[px]
1	TNG	HARPS-N	2017-08-16	21:21	03:56	$-$	200	90	1.089$\to$1.001$\to$2.089	Spec.
2	TNG	GIARPS	2018-07-12	21:27	05:15	$-$	200	116	1.604$\to$1.001$\to$1.527	Spec.
3	TNG	GIARPS	2018-07-19	21:25	04:23	$-$	300	78	1.442$\to$1.000$\to$1.389	Spec.
4	CA 3.5 m	CARMENES	2017-08-23	21:09	02:36	$-$	192	74	1.012$\to$1.006$\to$1.903	Spec.
2	TCS	MuSCAT2	2018-07-12	23:49	05:41	$g$	0.7	3365	1.066$\to$1.001$\to$1.767	14,23,29
2	TCS	MuSCAT2	2018-07-12	23:49	05:41	$r$	0.9	3262	1.066$\to$1.001$\to$1.767	14,14,23
2	TCS	MuSCAT2	2018-07-12	23:49	05:41	$i$	0.9	3205	1.066$\to$1.001$\to$1.767	15,23,29
2	TCS	MuSCAT2	2018-07-12	23:49	05:41	$z_s$	0.9	3117	1.066$\to$1.001$\to$1.767	14,14,23
5	TCS	MuSCAT2	2018-09-30	20:16	00:30	$g$	2	2321	1.002$\to$1.001$\to$1.801	29,32,37
5	TCS	MuSCAT2	2018-09-30	20:16	00:30	$r$	2	2322	1.002$\to$1.001$\to$1.801	27,32,37
5	TCS	MuSCAT2	2018-09-30	20:16	00:30	$i$	4	2322	1.002$\to$1.001$\to$1.801	28,30,35
5	TCS	MuSCAT2	2018-09-30	20:16	00:30	$z_s$	10	1320	1.002$\to$1.001$\to$1.801	29,36,40

Table 3. Table 3: Best fit parameters from MuSCAT2 transit light curve.

Parameter	Symbol	Units	Value
Transit epoch	$T_0$	BJD(44footnotemark: 4)	2457503.120120 $\pm$ 0.00018
Period	$P$	days	3.47410196 $\pm$ 0.00000106
Inclination	$i$	deg	85.61 $\pm$ 0.19
Radius ratio	$R_p/R_{\star }$	$-$	0.1176 $\pm$ 0.0014
	${\left(R_p/R_{\star }\right)}^2$	$-$	0.01382 $\pm$ 0.00033
Semimajor axis	$a/R_{\star }$	$-$	7.254 $\pm$ 0.096
Mid-transit 1(11footnotemark: 1)	$T_{mid}^1$	BJD	2458312.58566 $\pm$ 0.00022
Mid-transit 2(33footnotemark: 3)	$T_{mid}^2$	BJD	2458392.48932 $\pm$ 0.00042

Table 4. Table 4: Best fit parameters and 1- σ 𝜎 \sigma error bars from the MCMC analysis of the H α 𝛼 \alpha line for each individual night and combined HARPS-N observations.

	$h$	FWHM	Kp	vwind	Rλ (11footnotemark: 1)	R${}_{\lambda }{}^{}\left(h\right)$ (44footnotemark: 4)
	[$\%$]	[km/s]	[km/s]	[km/s]	[Rp]	[Rp]
HARPS-N 1	$-{0.75}_{-0.11}^{+0.10}$	${34.9}_{-6.4}^{+7.7}$	${140.3}_{-42.1}^{+39.0}$	${1.3}_{-2.9}^{+3.0}$	${1.21}_{-0.05}^{+0.04}$	${1.24}_{-0.06}^{+0.06}$
HARPS-N 2	$-{0.52}_{-0.11}^{+0.10}$	${15.6}_{-3.9}^{+3.7}$	${101.5}_{-33.5}^{+33.7}$	$-{3.9}_{-2.4}^{+2.1}$	${1.08}_{-0.03}^{+0.03}$	${1.17}_{-0.07}^{+0.06}$
HARPS-N 3	$-{0.69}_{-0.11}^{+0.10}$	${16.4}_{-2.4}^{+2.6}$	${196.0}_{-32.0}^{+30.1}$	$-{2.8}_{-2.2}^{+2.2}$	${1.15}_{-0.04}^{+0.04}$	${1.22}_{-0.06}^{+0.06}$
CARMENES 1	$-{0.85}_{-0.03}^{+0.03}$	${22.6}_{-0.9}^{+0.9}$	${166.2}_{-7.4}^{+7.3}$	$-{4.5}_{-0.5}^{+0.5}$	${1.21}_{-0.01}^{+0.01}$	${1.27}_{-0.02}^{+0.02}$
HARPS-N$x3$	$-{0.68}_{-0.06}^{+0.06}$	${19.0}_{-1.6}^{+1.7}$	${165.6}_{-16.6}^{+16.7}$	$-{3.0}_{-1.2}^{+1.2}$	${1.13}_{-0.02}^{+0.02}$	${1.22}_{-0.04}^{+0.04}$

Table 5. Table 5: Summary of the results obtained in the analysis presented here, for the different species and instruments. Hx3 corresponds to the analysis of three HARPS-N nights combined, and Cx1 to the analysis of the unique transit observed with CARMENES.

		$h$ [$\%$]	$R_{\lambda }\left(h\right)$ [Rp]	$v_{wind}$ [km/s]	FWHM [km/s]	$K_p$ [km/s]	AD${}_{}{}^{1.5\AA }{}_{\mathrm{TS}}{}^{}$ [$\%$] (11footnotemark: 1)	AD${}_{}{}^{1.5\AA }{}_{\mathrm{TLC}}{}^{}$ [$\%$] (33footnotemark: 3)	AD${}_{}{}^{1.5\AA }{}_{\mathrm{EMC}}{}^{}$ [$\%$] (55footnotemark: 5)
H$\alpha$	Hx3	$-0.68\pm 0.06$	$1.22\pm 0.04$	$-3.0\pm 1.2$	$19.0\pm 1.6$	$165.6\pm 16.7$	$0.33\pm 0.05$	$0.30\pm 0.05$	$0.36\pm 0.14$
	Cx1	$-0.85\pm 0.03$	$1.27\pm 0.02$	$-4.5\pm 0.5$	$22.6\pm 0.9$	$166.2\pm 7.4$	$0.44\pm 0.04$	$0.43\pm 0.05$	$0.40\pm 0.14$
H$\beta$	Hx3	$-0.45\pm 0.05$	$1.15\pm 0.03$	$-1.2\pm 1.4$	$19.4\pm 2.5$	$136.2\pm 18.6$	$0.13\pm 0.04$	$0.17\pm 0.04$	$0.14\pm 0.08$
H$\gamma$	Hx3	$-0.38\pm 0.08$	$1.13\pm 0.05$	$-2.3\pm 2.7$	$16.6\pm 4.2$	$135.0\pm 34.8$	$0.09\pm 0.05$	$0.07\pm 0.09$	$0.08\pm 0.04$
CaII	Cx1	$-0.56\pm 0.05$	$1.18\pm 0.03$	$-1.9\pm 0.6$	$9.2\pm 1.0$	$157.7\pm 8.2$	$0.16\pm 0.04$	$0.14\pm 0.04$	$0.11\pm 0.09$
NaI	Hx3	$-0.34\pm 0.05$	$1.11\pm 0.03$	$-3.1\pm 0.9$	$9.2\pm 2.0$	$182.5\pm 14.3$	$0.05\pm 0.03$	$0.08\pm 0.03$	$0.07\pm 0.05$
	Cx1	$-0.29\pm 0.04$	$1.10\pm 0.03$	$-3.2\pm 0.7$	$8.0\pm 1.2$	$176.6\pm 11.7$	$0.04\pm 0.03$	$0.03\pm 0.04$	$0.08\pm 0.08$
FeII	Hx3	$-0.33\pm 0.05$	$1.11\pm 0.03$	$-2.8\pm 0.8$	$7.2\pm 1.2$	$174.4\pm 14.0$	$0.04\pm 0.03$	$0.01\pm 0.03$	$0.04\pm 0.04$

Table 6. Table 6: Best fit parameters and 1- σ 𝜎 \sigma error bars from the MCMC analysis.

		h	FWHM	Kp	vwind	Rλ (11footnotemark: 1)	R${}_{\lambda }{}^{}\left(h\right)$ (44footnotemark: 4)
		[$\%$]	[km/s]	[km/s]	[km/s]	[Rp]	[Rp]
H$\beta$	HARPS-N 1	$-{0.20}_{-0.19}^{+0.16}$	${14.7}_{-9.4}^{+16.0}$	${110.6}_{-71.7}^{+85.5}$	$-{5.1}_{-6.6}^{+8.0}$	${1.12}_{-0.03}^{+0.04}$	${1.07}_{-0.13}^{+0.11}$
	HARPS-N 2	$-{0.57}_{-0.08}^{+0.07}$	${17.9}_{-2.6}^{+2.8}$	${131.8}_{-18.4}^{+19.3}$	$-{1.5}_{-1.7}^{+1.6}$	${0.89}_{-0.03}^{+0.03}$	${1.19}_{-0.05}^{+0.04}$
	HARPS-N 3	$-{0.38}_{-0.09}^{+0.08}$	${24.2}_{-6.0}^{+7.0}$	${166.3}_{-43.7}^{+55.6}$	$+{4.2}_{-3.5}^{+3.6}$	${0.86}_{-0.03}^{+0.03}$	${1.13}_{-0.06}^{+0.05}$
	HARPS-Nx3	$-{0.45}_{-0.05}^{+0.05}$	${19.4}_{-2.4}^{+2.7}$	${136.2}_{-18.5}^{+18.7}$	$-{1.2}_{-1.5}^{+1.4}$	${0.91}_{-0.02}^{+0.02}$	${1.15}_{-0.03}^{+0.03}$
H$\gamma$	HARPS-N 1	$-{0.29}_{-0.26}^{+0.13}$	${29.5}_{-15.8}^{+56.3}$	${189.7}_{-105.8}^{+43.7}$	$-{1.4}_{-6.3}^{+7.2}$	${1.01}_{-0.05}^{+0.05}$	${1.10}_{-0.17}^{+0.08}$
	HARPS-N 2	$-{0.43}_{-0.12}^{+0.10}$	${16.9}_{-4.6}^{+5.6}$	${101.5}_{-105.8}^{+43.7}$	$-{4.0}_{-3.0}^{+2.8}$	${0.86}_{-0.03}^{+0.03}$	${1.15}_{-0.08}^{+0.07}$
	HARPS-N 3	$-{0.28}_{-0.23}^{+0.17}$	${18.5}_{-11.2}^{+18.6}$	${104.6}_{-65.8}^{+81.0}$	$-{2.8}_{-7.5}^{+8.5}$	${0.78}_{-0.05}^{+0.05}$	${1.13}_{-0.06}^{+0.05}$
	HARPS-Nx3	$-{0.38}_{-0.09}^{+0.08}$	${16.6}_{-4.0}^{+4.5}$	${135.0}_{-35.6}^{+34.1}$	$-{2.3}_{-3.0}^{+2.5}$	${0.85}_{-0.03}^{+0.02}$	${1.13}_{-0.06}^{+0.05}$
NaI D2	HARPS-N 1	$-{0.37}_{-0.16}^{+0.14}$	${5.93}_{-2.2}^{+4.0}$	${151.4}_{-36.4}^{+35.4}$	$-{0.9}_{-2.6}^{+2.5}$	${0.92}_{-0.05}^{+0.05}$	${1.13}_{-0.10}^{+0.09}$
	HARPS-N 2	$-{0.12}_{-0.02}^{+0.02}$	${96.4}_{-6.0}^{+2.7}$	${137.1}_{-91.1}^{+78.3}$	${1.4}_{-9.6}^{+9.0}$	${0.97}_{-0.04}^{+0.04}$	${1.04}_{-0.01}^{+0.01}$
	HARPS-N 3	$-{0.28}_{-0.10}^{+0.08}$	${14.9}_{-6.2}^{+8.5}$	${103.3}_{-26.4}^{+26.8}$	$-{1.3}_{-2.4}^{+2.2}$	${0.94}_{-0.04}^{+0.04}$	${1.10}_{-0.07}^{+0.05}$
	HARPS-Nx3	$-{0.33}_{-0.06}^{+0.06}$	${8.9}_{-1.9}^{+2.7}$	${171.1}_{-12.8}^{+15.2}$	$-{2.4}_{-0.9}^{+0.9}$	${0.94}_{-0.03}^{+0.03}$	${1.11}_{-0.04}^{+0.04}$
	CARMENES 1	$-{0.29}_{-0.04}^{+0.04}$	${8.3}_{-1.2}^{+1.3}$	${176.5}_{-12.3}^{+12.8}$	$-{2.8}_{-0.8}^{+0.8}$	${0.91}_{-0.02}^{+0.02}$	${1.10}_{-0.03}^{+0.03}$
NaI D1	HARPS-N 1	$-{0.29}_{-0.11}^{+0.11}$	${10.2}_{-3.2}^{+5.1}$	${177.8}_{-41.4}^{+39.9}$	$-{2.2}_{-2.4}^{+2.4}$	${0.97}_{-0.05}^{+0.06}$	${1.10}_{-0.07}^{+0.07}$
	HARPS-N 2	$-{0.12}_{-0.02}^{+0.01}$	${95.8}_{-7.5}^{+3.2}$	${143.3}_{-93.1}^{+75.2}$	$-{7.4}_{-5.7}^{+8.6}$	${0.90}_{-0.03}^{+0.04}$	${1.04}_{-0.01}^{+0.01}$
	HARPS-N 3	$-{0.09}_{-0.02}^{+0.02}$	${81.7}_{-21.8}^{+12.9}$	${196.5}_{-79.3}^{+40.9}$	$-{6.4}_{-6.1}^{+9.6}$	${0.97}_{-0.04}^{+0.04}$	${1.03}_{-0.02}^{+0.01}$
	HARPS-Nx3	$-{0.35}_{-0.05}^{+0.05}$	${9.6}_{-1.5}^{+1.9}$	${194.0}_{-14.4}^{+15.1}$	$-{3.9}_{-0.9}^{+0.9}$	${0.91}_{-0.03}^{+0.03}$	${1.12}_{-0.03}^{+0.03}$
	CARMENES 1	$-{0.29}_{-0.04}^{+0.04}$	${7.8}_{-1.0}^{+1.2}$	${176.8}_{-10.6}^{+10.8}$	$-{3.6}_{-0.7}^{+0.7}$	${0.90}_{-0.02}^{+0.02}$	${1.10}_{-0.03}^{+0.03}$
CaII $\lambda 8498\AA$	CARMENES 1	$-{0.52}_{-0.05}^{+0.05}$	${8.9}_{-0.9}^{+1.4}$	${173.2}_{-7.5}^{+7.8}$	$-{3.2}_{-0.6}^{+0.6}$	${0.87}_{-0.03}^{+0.03}$	${1.17}_{-0.03}^{+0.03}$
CaII $\lambda 8542\AA$	CARMENES 1	$-{0.60}_{-0.04}^{+0.04}$	${11.7}_{-0.9}^{+1.0}$	${139.2}_{-9.2}^{+8.9}$	$-{1.3}_{-0.6}^{+0.6}$	${1.09}_{-0.02}^{+0.02}$	${1.20}_{-0.03}^{+0.03}$
CaII $\lambda 8662\AA$	CARMENES 1	$-{0.55}_{-0.06}^{+0.06}$	${6.9}_{-0.9}^{+1.1}$	${160.7}_{-7.8}^{+8.0}$	$-{1.1}_{-0.5}^{+0.5}$	${0.95}_{-0.03}^{+0.02}$	${1.18}_{-0.04}^{+0.04}$
FeII $\lambda 5018\AA$	HARPS-Nx3	$-{0.26}_{-0.05}^{+0.04}$	${7.9}_{-1.3}^{+1.3}$	${181.9}_{-18.7}^{+18.3}$	$-{2.7}_{-1.1}^{+1.1}$	${1.19}_{-0.03}^{+0.03}$	${1.09}_{-0.03}^{+0.03}$
FeII $\lambda 5169\AA$	HARPS-Nx3	$-{0.47}_{-0.05}^{+0.05}$	${8.1}_{-1.0}^{+1.1}$	${170.6}_{-7.8}^{+8.3}$	$-{3.6}_{-0.5}^{+0.5}$	${1.18}_{-0.02}^{+0.02}$	${1.16}_{-0.03}^{+0.03}$
FeII $\lambda 5316\AA$	HARPS-Nx3	$-{0.26}_{-0.06}^{+0.05}$	${5.6}_{-1.1}^{+1.4}$	${170.8}_{-16.6}^{+14.4}$	$-{2.1}_{-0.8}^{+0.8}$	${1.18}_{-0.03}^{+0.03}$	${1.09}_{-0.04}^{+0.03}$
	CARMENES 1	$-{0.55}_{-0.09}^{+0.08}$	${10.3}_{-2.0}^{+2.3}$	${139.6}_{-22.6}^{+32.3}$	$+{0.2}_{-1.3}^{+1.7}$	${1.24}_{-0.06}^{+0.05}$	${1.18}_{-0.05}^{+0.05}$
MgI $\lambda 5173$	HARPS-Nx3	$-{0.23}_{-0.04}^{+0.03}$	${12.8}_{-2.0}^{+2.4}$	${159.6}_{-23.8}^{+30.9}$	$-{0.9}_{-1.9}^{+1.7}$	${1.23}_{-0.02}^{+0.02}$	${1.08}_{-0.02}^{+0.02}$

Table 7. Table 7: Absorption depth (in % percent \% ) measured for 0.75 ​ Å 0.75 Å 0.75\leavevmode\nobreak\ \mathrm{\AA} and 1.5 ​ Å 1.5 Å 1.5\leavevmode\nobreak\ \mathrm{\AA} bandwidths calculated in the final transmission spectra (TS) and transmission light curves (TLC).

		TS		TLC
		$0.75\AA$	$1.5\AA$	$0.75\AA$	$1.5\AA$
H$\alpha$	HARPS-N 1	$0.69\pm 0.10$	$0.44\pm 0.07$	$0.69\pm 0.11$	$0.45\pm 0.10$
	HARPS-N 2	$0.41\pm 0.08$	$0.20\pm 0.06$	$0.37\pm 0.10$	$0.16\pm 0.08$
	HARPS-N 3	$0.61\pm 0.10$	$0.24\pm 0.07$	$0.45\pm 0.12$	$0.10\pm 0.10$
	HARPS-N$x3$	$0.59\pm 0.07$	$0.33\pm 0.05$	$0.58\pm 0.06$	$0.30\pm 0.05$
	CARMENES 1	$0.68\pm 0.05$	$0.44\pm 0.04$	$0.68\pm 0.06$	$0.43\pm 0.05$
H$\beta$	HARPS-N 1	$0.30\pm 0.09$	$0.17\pm 0.06$	$0.35\pm 0.10$	$0.22\pm 0.08$
	HARPS-N 2	$0.35\pm 0.07$	$0.14\pm 0.05$	$0.37\pm 0.08$	$0.17\pm 0.07$
	HARPS-N 3	$0.23\pm 0.08$	$0.10\pm 0.06$	$0.28\pm 0.10$	$0.15\pm 0.08$
	HARPS-N$x3$	$0.28\pm 0.06$	$0.13\pm 0.04$	$0.31\pm 0.05$	$0.17\pm 0.04$
H$\gamma$	HARPS-N 1	$0.19\pm 0.10$	$0.09\pm 0.07$	$0.18\pm 0.11$	$0.09\pm 0.10$
	HARPS-N 2	$0.28\pm 0.08$	$0.14\pm 0.06$	$0.29\pm 0.09$	$0.15\pm 0.08$
	HARPS-N 3	$0.13\pm 0.10$	$0.04\pm 0.07$	$0.08\pm 0.13$	$0.00\pm 0.11$
	HARPS-N$x3$	$0.21\pm 0.07$	$0.09\pm 0.05$	$0.17\pm 0.06$	$0.07\pm 0.05$
NaI D2	HARPS-N 1	$0.07\pm 0.07$	$0.02\pm 0.05$	$0.06\pm 0.08$	$0.01\pm 0.07$
	HARPS-N 2	$0.12\pm 0.06$	$0.05\pm 0.04$	$0.15\pm 0.07$	$0.09\pm 0.06$
	HARPS-N 3	$0.10\pm 0.06$	$0.06\pm 0.05$	$0.10\pm 0.08$	$0.06\pm 0.07$
	HARPS-N$x3$	$0.09\pm 0.05$	$0.05\pm 0.03$	$0.15\pm 0.04$	$0.07\pm 0.04$
	CARMENES 1	$0.10\pm 0.05$	$0.03\pm 0.03$	$0.08\pm 0.05$	$0.01\pm 0.05$
NaI D1	HARPS-N 1	$0.16\pm 0.07$	$0.10\pm 0.05$	$0.14\pm 0.04$	$0.08\pm 0.07$
	HARPS-N 2	$0.09\pm 0.06$	$0.04\pm 0.04$	$0.05\pm 0.07$	$0.01\pm 0.06$
	HARPS-N 3	$0.11\pm 0.06$	$0.05\pm 0.05$	$0.10\pm 0.08$	$0.04\pm 0.07$
	HARPS-N$x3$	$0.10\pm 0.05$	$0.06\pm 0.03$	$0.16\pm 0.04$	$0.09\pm 0.04$
	CARMENES 1	$0.11\pm 0.04$	$0.06\pm 0.03$	$0.09\pm 0.05$	$0.04\pm 0.04$
Na D21 (11footnotemark: 1)	HARPS-N 1	$0.11\pm 0.07$	$0.06\pm 0.05$	$0.10\pm 0.07$	$0.05\pm 0.06$
	HARPS-N 2	$0.10\pm 0.06$	$0.04\pm 0.04$	$0.10\pm 0.06$	$0.05\pm 0.05$
	HARPS-N 3	$0.10\pm 0.06$	$0.05\pm 0.05$	$0.10\pm 0.07$	$0.05\pm 0.06$
	HARPS-N$x3$	$0.09\pm 0.05$	$0.05\pm 0.03$	$0.15\pm 0.04$	$0.08\pm 0.03$
	CARMENES 1	$0.11\pm 0.04$	$0.04\pm 0.03$	$0.09\pm 0.04$	$0.03\pm 0.04$
CaII $\lambda 8498\AA$	CARMENES 1	$0.28\pm 0.05$	$0.15\pm 0.04$	$0.28\pm 0.06$	$0.14\pm 0.05$
CaII $\lambda 8542\AA$	CARMENES 1	$0.41\pm 0.05$	$0.18\pm 0.04$	$0.40\pm 0.06$	$0.17\pm 0.05$
CaII $\lambda 8662\AA$	CARMENES 1	$0.27\pm 0.06$	$0.14\pm 0.04$	$0.26\pm 0.07$	$0.15\pm 0.04$
CaII Combined(33footnotemark: 3)	CARMENES 1	$0.32\pm 0.05$	$0.16\pm 0.04$	$0.31\pm 0.05$	$0.14\pm 0.04$
FeII $\lambda 5018\AA$	HARPS-Nx3	$0.09\pm 0.04$	$0.04\pm 0.03$	$0.03\pm 0.04$	$0.02\pm 0.03$
FeII $\lambda 5169\AA$	HARPS-Nx3	$0.09\pm 0.04$	$0.04\pm 0.03$	$0.12\pm 0.04$	$0.05\pm 0.03$
FeII $\lambda 5316\AA$	HARPS-Nx3	$0.07\pm 0.04$	$0.04\pm 0.03$	$0.04\pm 0.03$	$0.00\pm 0.03$
FeII Combined(33footnotemark: 3)	HARPS-Nx3	$0.08\pm 0.04$	$0.04\pm 0.03$	$0.06\pm 0.03$	$0.01\pm 0.03$
MgI $\lambda 5173\AA$ HARPS-Nx3	HARPS-Nx3	$0.07\pm 0.04$	$0.02\pm 0.03$	$0.09\pm 0.04$	$0.04\pm 0.03$

Full text

1 Introduction

The study of Hot Jupiters atmospheres has revealed a broad composition and structural diversity in gas-rich planets (Jensen et al. 2012, Sing et al. 2012, Crossfield 2015, Deming & Seager 2017, Madhusudhan et al. 2016). When these planets are under extreme conditions, in particular, located very close to their host stars and consequently being strongly irradiated, their temperatures increase and lead to different atmospheric chemistry between the day and night-sides of atmosphere (Arcangeli et al. 2018, Bell & Cowan 2018, Helling & Rimmer 2019) in comparison with cooler planets. This creates a new type of exoplanets called ultra hot Jupiters (hereafter UHJ).

We define UHJs as planets with day-side temperatures higher than $2200\leavevmode\nobreak\ \mathrm{K}$ (Parmentier et al., 2018). One of the most significant differences between hot and ultra hot Jupiters is perhaps the presence of water vapour in their atmospheres, which in hot Jupiters has become a common species but is missing in the day-side of ultra hot Jupiters, even though the ingredients to form water vapor are present. The high irradiation received by the day-side of such UHJs causes the atmospheric gas temperature to increase to further more than $3000$K, as for for example demonstrated for WASP-18b (Helling et al., 2019). The dominating gas species in the upper, day-side atmosphere is therefore atomic hydrogen (HI) instead of molecular hydrogen (H2). Water (H2O) is by $\approx$ 5 orders of magnitude less abundant than on the night-side. Such high day-side gas temperatures further result in an increased number of H- compared to cooler giant gas planets like HD 189733b. In the case of WASP-18b, in the upper low-pressure atmosphere, the day-side has more than $10$ orders of magnitude more H- and H+ than the night-side. Therefore, H- opacity at the day-side plays an important role (Lothringer et al., 2018).

Several ultra hot Jupiters have been studied: WASP-33b (Haynes et al., 2015), WASP-121b (Evans et al., 2017), WASP-103b (Kreidberg et al., 2018), WASP-18b (Sheppard et al. 2017, Arcangeli et al. 2018, Helling et al. 2019) , HAT-P-7b (Armstrong et al., 2016), most of them presenting thermal inversions. Focusing on transmission spectroscopy at high resolution, a study of the atmosphere of KELT-9b, the hottest planet known to date ($T_{eq}=4050\leavevmode\nobreak\ \mathrm{K}$ and day-side temperature of $4600\leavevmode\nobreak\ \mathrm{K}$), by Yan & Henning (2018) shows a detection of an extended atmosphere of H$\alpha$ produced by irradiation. On the other hand, Hoeijmakers et al. (2018) studied the atmosphere of this same planet, detecting Fe, Fe+ and Ti+ and concluding that at this high temperatures metals are preponderantly in their ionic forms. Cauley et al. (2019) further reported the detection of H$\beta$ and the optical MgI triplet, and found that planetary rotational broadening is needed to reproduce the Balmer line transmission profile shapes.

Here, we present the study of the atmosphere of MASCARA-2b (Talens et al., 2017), also known as KELT-20b (Lund et al., 2017), using four transit observations, three of them observed with the HARPS-North spectrograph and one observed with the CARMENES spectrograph. MASCARA-2b is an ultra hot Jupiter with an equilibrium temperature of $2260\leavevmode\nobreak\ \mathrm{K}$, located in a $3.5\leavevmode\nobreak\ \mathrm{day}$ orbit around an A-type star ($m_{V}=7.6$) with an effective temperature ($T_{\mathrm{eff}}$) of $8980\leavevmode\nobreak\ \mathrm{K}$ and a rotation ($v\sin i$) of $116\leavevmode\nobreak\ \mathrm{km/s}$. MASCARA-2b orbits very close to its host star ($a=0.057\leavevmode\nobreak\ \mathrm{AU}$) and is strongly irradiated. Its mass, because radial-velocity measurements of rapid-rotators are very challenging, has not been measured, and only an upper limit is available: $M_{p}<3.510\leavevmode\nobreak\ M_{J}$. The details about the system are summarized in Table 1. In Figure 1 our target can be contextualized among the other known ultra hot Jupiters. MASCARA-2b orbits around the second hottest host star of the UHJ population, although its equilibrium temperature is typical for UHJs. Given that the effective temperature of MASCARA-2 is larger than WASP-18, the effect of thermal ionisation of the day-side can be expected to be stronger for MASCARA-2b/KELT-20b due to a stronger irradiation by the host star.

Due to its brightness and the fast rotation of the host star, MASCARA-2b is an ideal target for transmission spectroscopy. Its atmosphere was previously studied in Casasayas-Barris et al. (2018), where only the first HARPS-North transit on 16 August 2017 was analyzed, finding evidence of H$\alpha$ and NaI absorption features and noisy features at H$\beta$ position. Here, combining four transits, we achieve higher signal-to-noise ratio (SNR) being able to clearly distinguish the center-to-limb variation (CLV) and the Rossiter-McLaughlin effect (RME), whose treatment will become very important for atmospheric studies of this kind of planets when observing with ESPRESSO-like facilities.

This paper is organized in seven sections. In Section 1 we present a general introduction of transmission spectroscopy and ultra hot Jupiters. In Section 2 we detail the observations. Methods to extract the color-photometry light curves and high resolution transmission spectra are explained in Section 3. In Section 4 we explain how the CLV and RME are modeled. In Section 5 we present the results obtained in the analysis of high resolution transmission spectra. The analysis of the systematic effects are presented in Section 6 and the conclusions in Section 7.

2 Observations

We observed a total of four transits of MASCARA-2b using the HARPS-N and CARMENES high resolution spectrographs. One of these transits was simultaneously observed with MuSCAT2 (Multicolor Simultaneous Camera for studying Atmospheres of Transiting exoplanets 2), a 4-color simultaneous imager, in order to monitor possible stellar activity. One additional epoch was also observed with MuSCAT2 to reproduce the results of the first observation. Details about the observations can be found in Table 2 and Figure 3.

2.1 HARPS-N observations

We observed three transits of MASCARA-2b using the HARPS-North spectrograph mounted on the $3.58$-m TNG telescope, located at the Observatorio del Roque de los Muchachos (ORM, La Palma). The observations were on 16 August 2017, 12 and 19 July 2018. For the 2 last observations we used the GIARPS mode, which permits simultaneous high resolution spectroscopy observations in the near-infrared and in the optical wavelength regions, using GIANO-B and HARPS-North, respectively.

We exposed continuously before, during and after the transit, using the fiber A on the target and fiber B on the sky in order to monitor possible sky emissions. For the first and second night we used an exposure time of $200$ s while in the third night we increased it to $300$ s, to obtain a higher signal to noise ratio (SNR).

2.2 CARMENES observations

One more transit was observed with CARMENES spectrograph (Calar Alto high-Resolution search for M dwarfs with Exoearths with Near-infrared and optical Echelle Spectrographs) located at the Calar Alto Observatory, on 23 August 2017. The strategy of these observations was the same as with HARPS-North. For these observations an exposure time of $192$ s was used, obtaining a SNR of around $100$ near the H$\alpha$ line.

2.3 Photometric observations with MuSCAT2

The transit of MASCARA-2b on 12 July 2018 was observed simultaneously using MuSCAT2 (Narita et al., 2018) located at the $1.52\leavevmode\nobreak\ \rm{m}$ Carlos Sánchez Telescope (TCS) at the Teide Observatory. Another transit on 30 September 2018 was observed using MuSCAT2 only. Both observations were heavily defocussed to avoid saturation. While the first observation covered the whole transit event, the second observation had to stop at $\sim$43 minutes after the mid-transit point due to low elevation of the target. Both nights were mostly clear. Although defocussed, the profile of point spread function (PSF) was not stable for both observations. The observing log is presented in Table 2.

3 Methods

In this section we present the data processing applied to extract the four color light curves observed with MuSCAT2 (photometry) and the transmission spectra using the high resolution observations (spectroscopy).

3.1 Photometry

Two independent approaches were used for the MuSCAT2 photometry, light curve detrending, and transit modelling to search for possible spot-crossing events.

The first approach uses the dedicated MuSCAT2 photometry and transit modelling pipeline: the photometry follows the standard aperture photometry practices, the transit modelling is carried out with PyTransit (Parviainen, 2015) and the priors for the quadratic limb darkening coefficients are calculated using LDTk (Parviainen & Aigrain, 2015). The final transit light curve modelling was carried out jointly for the two nights and all passbands, using a Gaussian process (GP) based systematics model calculated with george (Ambikasaran et al., 2014) with the sky level, airmass, seeing, and CCD centroid shifts as covariates.

The second approach followed closely the approach described in Chen et al. (2014). The photometry was carried out with the DAOPHOT package in the NASA IDL Astronomy User’s Library, the transit modelling was done with batman (Kreidberg, 2015), Celerite (Foreman-Mackey et al., 2017) was used to implement a Gaussian process-based systematics model, and the method by Espinoza & Jordán (2015) was used to obtain priors on the quadratic limb darkening coefficients.

In both cases, emcee (Foreman-Mackey et al., 2013a) was used to estimate the model parameter posteriors. The parameter posterior estimates are presented in Table 3, and the MuSCAT2 light curves are shown in Figure 2. The transit light curves do not feature evident spot-crossing events, nor do they show significant colour-dependent depth variations.

3.2 Spectroscopy

The observations were reduced with the HARPS-North Data reduction Software (DRS), version $3.7$. The DRS extracts the spectra order-by-order, which are then flat-fielded using the daily calibration set. A blaze correction and the wavelength calibration are applied to each spectral order and, finally, all the spectral orders from each two-dimensional echelle spectrum are combined and resampled into a one-dimensional spectrum. The resulting one-dimensional spectra cover the optical range between $3800\leavevmode\nobreak\ \mathrm{\AA}$ and $6900\leavevmode\nobreak\ \mathrm{\AA}$ and has a wavelength step of $0.01\leavevmode\nobreak\ \mathrm{\AA}$. The spectra are referred to the Baricentric rest frame and the wavelengths are given in the air.

All nights present stable seeing at ${\sim 1\leavevmode\nobreak\ \mathrm{arcsec}}$ with the Moon at more than $90\leavevmode\nobreak\ \mathrm{deg}$ from our target. HARPS-N Night 1 and CARMENES observations are observed to be very stable in terms of SNR during the night. For Night 2 and 3, however, we noticed variable SNR and variations of the continuum similar to a bad blaze function correction. This was caused by a problem with the Atmospheric Dispersion Corrector (ADC) of the telescope, that causes PSF distortions and variations of the continuum. For these both nights, re-reducing the data with the DRS did not produced any improvement. The way of correcting this variation is discussed in detail in Section 5.1. On the other hand, during the Night 2 we observe a drop in the SNR of some spectra before the transit, possibly caused by the presence of a cloud. We discard the spectra from 22:58 UT until 23:24 UT ($8$ spectra), where the SNR is lower than $35$ (around $6600\leavevmode\nobreak\ \mathrm{\AA}$). Additionally, for consistency with the SNR established for the cloud region, for this night we also removed the spectra taken at 01:28UT which also present SNR smaller than $35$ (i.e. a total of $9$ spectra are discarded).

CARMENES observations are processed with the CARMENES pipeline CARACAL (CARMENES Reduction And Calibration; Caballero et al. (2016)) which considers bias, flat-relative optimal extraction (Zechmeister et al., 2014), cosmic ray correction, and a wavelength calibration described in Bauer et al. (2015). The reduced spectra are referred to the terrestrial rest frame and the wavelengths are given in the vacuum. At the beginning of the night some technical problems lead to some very short time exposures with very low SNR which are directly discarded. The spectra, SNR and airmass evolution of each night can be observed in Figure 3.

3.2.1 Telluric correction

Taking advantage that the observations were carried out using fiber B monitoring the sky, we checked for possible atmospheric emissions features during the night. We did not observe any telluric emission feature that needs to be corrected, and since subtracting the sky to the target spectra introduces additional noise, we did not apply this correction.

Following Allart et al. (2017), we correct the telluric absorption contamination using version 1.2.0 of Molecfit (Smette et al. 2015 and Kausch et al. 2015). Molecfit is an ESO tool built to correct Earth atmospheric features in ground-based spectra, which uses a line-by-line radiative transfer model (LBLRTM). HARPS-North spectra are given in the Solar System baricentric rest frame, while Molecfit models are computed in the terrestrial rest frame. For this reason, before fitting the telluric contamination we shift the spectra to the terrestrial rest frame using the Baricentric Earth Radial Velocity (BERV) information. We shift the spectra using the IRAF (Tody, 1986) routine dopcor. The parameters used to compute the high-resolution telluric spectrum are detailed in Allart et al. (2017), and the wavelength regions with strong lines used here are very similar to those used in this same paper, but considering only those regions out of stellar features.

In the case of CARMENES data, the spectra after CARACAL reduction are already given in the terrestrial rest frame. The Molecfit correction applied is slightly different, and follows the method discussed in Nortmann et al. (2018) and Salz et al. (2018). The telluric uncorrected and corrected spectra in the atomic sodium (NaI) region for Night 4 obtained with CARMENES are shown in Figure 4.

3.2.2 General transmission spectrum extraction

Once the telluric correction has been applied, we cut the spectra in the wavelength region we want to study and normalize them by the continuum. This is done recursively over the whole spectral range covered by both spectrographs to screen for all possible planetary signals. Then, as the spectra after the Molecfit correction is referred to the terrestrial rest frame, we shift them to the stellar rest frame considering the system velocity ($\gamma$) and the BERV information. The stellar reflex motion is not considered because the mass of the planet ($M_{p}$) has not yet been determined, and the radial velocity measurements in rapid-rotating stars are difficult to determine. If we assume $K_{s}=322.51\leavevmode\nobreak\ \mathrm{m/s}$ (i.e. the mass upper limit $M_{p}=3.510\leavevmode\nobreak\ M_{J}$; Lund et al. 2017), the maximum stellar radial-velocity change is $\sim\pm 0.3\leavevmode\nobreak\ \mathrm{km/s}$. Since one pixel in HARPS is around $0.8\leavevmode\nobreak\ \mathrm{km/s}$ and the stellar lines present a large broadening caused by the fast stellar rotation, considering or not considering this reflex motion contribution to the radial velocity correction does not have any significant impact to the results, as presented in Casasayas-Barris et al. (2018). In addition, the spectra could be also affected by stellar activity and instrumental radial-velocity drifts. However, an A-type star like MASCARA-2 is not expected to be active, and the instrumental effects would similarly affect to all spectra, resulting in global wavelength offsets from the rest frame position (which after the daily calibration are expected to be of the order of $10\leavevmode\nobreak\ \mathrm{m/s}$ for the instruments used here). For CARMENES observations, when shifting the spectra to the stellar rest frame, we rebin the data to the same wavelength grid using the IDL routine rebinw from the PINTofALE package, which ensures flux conservation (Drake & Kashyap, 2010).

For the next steps, we follow Wyttenbach et al. (2015): we compute the out-of-transit master spectrum (combination of all the stellar spectra taken when the planet is not transiting) and we then divide each spectrum by this master spectrum to remove the stellar lines. The resulting residuals contain the variations of the stellar lines profile produced by the presence of the planet and additional systematic effects. These variations include the absorption produced by the exoplanet atmosphere, the center-to-limb variation (CLV) and the Rossiter-McLaughlin efect (RME). At this step we find the best fit model of the absorption, CLV and RME, and the resulting best fit parameters are used to compute the best fit CLV and RME model, which is then removed from the residuals (see Section 4 for a detailed explanation of the models). Once these effects have been corrected, we shift the remaining residuals to the planet rest frame and we then combine the in-transit residuals to find the transmission spectrum for each night.

HARPS-N Night 2 and Night 3 are affected by an ADC problem, presenting weird continuum variations similar to a bad blaze function correction. In order to solve this problem, we tried to re-reduce the raw data with the DRS pipeline, without any improvement. This effect introduces differences in the line wings changing in time, and which are more significant in very broad spectral lines, specially in the H$\alpha$ line. This differences become very important after dividing the spectra by the out-of-transit master spectrum, as the resulting spectra present broad variations that partially hide the narrow residuals we are looking for behind this broad effect. For this reason, after computing the residuals spectra we applied a very broad median filter using medfilt from scipy library with a kernel size of $1505$ to protect the planetary signals from being diluted and remove only the low frequency variations produced by the variation of the continuum, and the final transmission spectrum was normalised with UnivariateSplines from scipy in order to remove remaining undulations. This additional correction is only applied to HARPS-N Night 2 and Night 3.

Finally, we note that the parameters used during the whole process are those presented in Table 1 except for the epoch ($T_{0}$), orbital period ($P$), and inclination ($i$), which are from the recalculated parameters of Table 3.

4 Modeling and fitting the CLV and RME

The Rossiter-McLaughlin effect (RME) strongly affects the stellar lines profile of a fast-rotating star such as MASCARA-2b, which has a $v\sin{i}=116\leavevmode\nobreak\ \rm{km/s}$. As the transit progresses, the light from different regions of the stellar disk is blocked by the planet, producing a lack of flux in different parts of the stellar lines which depends on the geometry of the system. On the other hand, the center-to-limb variation (CLV) also affects the stellar lines profile (Czesla et al., 2015), which become broader and/or deeper than the integrated flux spectrum depending on the stellar properties of the region obscured by the planet.

The stellar spectrum is modeled using VALD3 (Ryabchikova et al., 2015) line lists and Kurucz ATLAS9 models, computed with Spectroscopy Made Easy tool (SME; Valenti & Piskunov 1996). As presented in Yan et al. (2017), using SME we compute the stellar spectrum at $21$ different limb-darkening angles ($\mu$) without considering the rotation broadening. For the stellar model computation we use local thermodynamic equilibrium (LTE) and solar abundance for the hydrogen lines. For the NaI, the solar abundance model strongly differs from the data. For this reason, we fit the NaI abundance using the SME tool. In this process, the spectra used to fit the abundance is the combination of all the out-of-transit data of Night 4 with CARMENES, which have higher SNR. After checking the fitted [Na/H] results for a range of [Fe/H] values from literature ([Fe/H]$=-0.02\pm 0.07$ in Talens et al. 2017) we noticed that the NaI abundance change significantly when varying the [Fe/H] value (see below). The most reasonable NaI abundance is obtained in the $1-\sigma$ upper limit value i.e. assuming [Fe/H]$=+0.05$ finding [Na/H]$=0.98$. As we do not consider non-LTE effects that could be important in this stellar type (Kubát & Korcáková, 2004). These [Fe/H], [Na/H] values together with the parameters given in Table 1 are those used to compute the stellar models. We note here that for our computation of these models, we need a good description of the line profile and the actual abundance value is not the focus of our analysis. The maximum relative difference obtained between models when assuming [Fe/H]$=-0.29$ (Lund et al., 2017) (for which we obtain the [Na/H]$=1.39$) or [Fe/H]$=+0.05$ ([Na/H]$=0.98$) is ${\sim 0.08\%}$.

We divide, then, the stellar disk in cells of size 0.01\leavevmode\nobreak\ R_{\star}\leavevmode\nobreak\x$\leavevmode\nobreak\ 0.01\leavevmode\nobreak\ R_{\star}$. Each of these cells contain its own properties, such as the rotation velocity due to $v\sin i$, which depends on its distance to the stellar rotation axis (i.e. the spin-orbit angle, $\lambda$), its $\mu$ angle value, and could be obscured by the planet at a given time. At the time of each exposure we can compute the position of the planet in the stellar disk using EXOFAST (Eastman et al., 2013) and compute the integrated stellar disk flux considering the obscured cells, their proper radial velocity shift caused by the rotation and the correspondent $\mu$ spectrum (which are linearly interpolated using the $21$ $\mu$ angle reference spectra). The only difference between the HARPS-N and CARMENES modeled spectra is the spectral resolution, $\mathcal{R}=115\leavevmode\nobreak\ 000$ and $\mathcal{R}=94\leavevmode\nobreak\ 600$, respectively.

The modeled stellar spectra at different orbital phases depend on the geometry of the system but also the planetary radius, $R_{p}$. For a particular wavelength, and depending on the amount of material, the atmosphere becomes opaque at a given pressure level or radius of the planet, $R_{\lambda}$. In order to consider this radius change when fitting the data, the modeled stellar spectra is computed for a grid of $R_{\lambda}$, increasing from $0.7\leavevmode\nobreak\ R_{p}$ ($1\leavevmode\nobreak\ R_{p}=1.83\leavevmode\nobreak\ R_{J}$; Talens et al. (2017)) to $2.5\leavevmode\nobreak\ R_{p}$ in steps of $0.1\leavevmode\nobreak\ R_{p}$. For a given $R_{\lambda}$, we then calculate the stellar spectrum by linearly interpolating these grid of models. Finally, the modeled residuals are computed dividing all the modeled stellar spectra by the out-of-transit model.

Finally, and similarly to Yan & Henning (2018), in the modeled in-transit residuals we add a Gaussian profile to describe the absorption of the planet. This profile, depends on the center of the Gaussian ($v_{0}$), the full width at half maximum ($FWHM$) and constrast ($h$). The center position is defined as $v_{0}=v_{p}+v_{wind}$ where $v_{p}=K_{p}\sin(2\pi\varphi)$ is the planet radial velocity in each in-transit exposure time ($\varphi$ is the orbital phase) and $v_{wind}$ is the radial velocity of the atmospheric wind.

With this, we perform a Markov chain Monte Carlo (MCMC) analysis to fit the observed residuals using emcee (Foreman-Mackey et al., 2013b) code. The free parameters are $R_{\lambda}$, $h$, $v_{wind}$, $K_{p}$ and $FWHM$. We use $14$ walkers and $1$x$10^{5}$ steps. The best fit $R_{\lambda}$ value is used to compute the CLV and RME model, which is then subtracted to the data, remaining only the true planetary absorption. Before combining the data residuals, this absorption is shifted to the planet rest frame using the best fit $K_{p}$ value. The MCMC analysis is based on binned spectra ($0.05\leavevmode\nobreak\ \mathrm{\AA}$ and $0.002$ orbital phase bins) and only applied to the fully in-transit data (i.e. excluding the ingress and egress). The planetary absorption spectrum in the ingress and in the egress could present different geometries (Yan & Henning, 2018), as it was recently shown in Salz et al. (2018).

During the modelling, we observe that in some cases the best-fit value of $R_{\lambda}$ is smaller than $1\leavevmode\nobreak\ R_{p}$. This is possibly due to the intrinsic error of the stellar models, which are computed considering LTE and solar abundance in case of the hydrogen lines (the NaI abundance is fitted). The non-LTE effects in these stellar-types can become strong (Kubát & Korcáková, 2004) and, as calculated previously, the abundance of the lines also influences the final model. For two different values of [Fe/H] inside the error bars presented in Lund et al. (2017) and Talens et al. (2017), we can obtain the same intensity of the RME for two different $R_{\lambda}$. Additionally, if we compare two models computed with the same abundance, one calculated for $0.85\leavevmode\nobreak\ R_{p}$ and one for $1.0\leavevmode\nobreak\ R_{p}$, in the mid-transit time their relative difference in the RME centre is ${\sim 0.2\%}$. Thus, there are different assumptions that can produce variations of the final $R_{p}$ values. However, our main goal here with the model fitting is to remove the strong RME and CLV effects, which are easily observed in the residuals. For this reason, considering the $R_{p}$ change permits us to find the best model that reproduce the effects, even if the best-fit $R_{\lambda}$ value is not exactly the expected. This is one of the reasons we are not linking the increase of radius ($R_{\lambda}$) with the absorption that we measure ($h$) in the models, and we are leaving both free, independently. We do not expect this $R_{\lambda}$ values to have any effect to our results and we use the models as a tool to discard the RME and CLV effects from the absorption signatures. When comparing models, we observe that in the mid-transit time the RME is more constraining for the $R_{\lambda}$ parameter. In the ingress and egress regions, where the CLV become important, both CLV and RME are constraining $R_{\lambda}$. However, we stress that only the the fully in-transit data is used in the fitting procedure i.e. where the RME is the more constraining for the $R_{\lambda}$ parameter.

5 Results

We present here the results obtained in the different wavelength regions. In all cases, we have analysed each individual night separately and, additionally, the three HARPS-N data together, sorting the spectra in time from the mid-transit center. We note that in this latter case all the spectra were telluric corrected first.

Here, we show the 2D maps of the final residuals obtained after applying the differential spectroscopy, the final transmission spectrum and the transmission light curves around several spectral lines detectable in MASCARA-2b’s transmission spectrum. The 2D maps of individual nights are presented in the Appendix A, their best-fit values in Table 6 and MCMC correlation diagrams in Appendix D.1. All light curves are also presented in Appendix C. The absorption depth values measured in the transmission spectra and light curves can be found in Table 7.

The transmission light curves are obtained calculating the equivalent width of the residuals after dividing each spectrum by the master out-of-transit spectrum, and then moved to the planet rest frame. The integration to obtain the equivalent width is computed using the scipy library (Jones et al., 2001) from Python2.7.

Hereafter, in Tables and Figures, we refer to the different nights as HARPS-N 1 (Night 1), HARPS-N 2 (Night 2), HARPS-N 3 (Night 3) and CARMENES 1 (Night 4) in order to not create confusion between instruments. The combined results of all three HARPS-N nights are presented as HARPS-N$x3$.

5.1 H$\alpha$

The residuals after computing the ratio between each spectra and the master out-of-transit spectrum, the best-fit models and the transmission spectra around H$\alpha$ ($6562.81\leavevmode\nobreak\ \mathrm{\AA}$) of each individual night can be observed in Figure 5. The best fit values and the 1-$\sigma$ error bars obtained with the MCMC analysis are presented in Table 4.

In Figure 5 we can see that the three nights observed with the HARPS-N are in general more noisy than the one with CARMENES, with residuals showing larger standard deviation. No obvious reason for this fact can be derived from checking the weather or seeing logs, only the ADC problem in HARPS-N Night 2 and 3 are noticeable. During the transit, the RME is observed as a red trail (positive relative flux) moving in the $\pm 116\leavevmode\nobreak\ \mathrm{km/s}$ ($v\sin i$) range of velocities and the absorption as a blue trail (negative relative flux) around $\pm 24\leavevmode\nobreak\ \mathrm{km/s}$, approximately.

Due to the noise, specially for Night 2 and Night 3 for which it increases during the transit, we are not able to totally fit the planetary absorption for the HARPS-N nights, which becomes underestimated by the model, and can not be recovered in the region were it intersects with the RME around mid-transit. However, for Night 4 where the SNR is higher, the absorption during the whole transit is recovered. In particular, this underestimation of the absorption can be observed in the residuals of the Gaussian profiles derived from the best fit model, which seem to present smaller contrast ($h$) than the absorption observed in the transmission spectrum of HARPS-N (see fourth and fifth row of Figure 5). When we fit all HARPS-N data together (see Figure 6), the central absorption still remains within the noise level and the absorption is slightly underestimated. This is possibly caused by the decrease in SNR during the transit. We also analyzed the HARPS-N data considering only Night 1 and Night 3, with results closer to those obtained for CARMENES in Night 4 but the mid-transit region still remains within the noise level. The final transmission spectra of CARMENES 1 and the analysis combining all three HARPS-N nights, together with the Gaussian profiles derived from the best-fit parameters, can be observed in Figure 7.

The observed transmission light curves before and after correcting the CLV and RME effects can be observed in Figure 8. The models are also obtained with the best fit parameters from the general fitting shown in Figure 5. In this light curves the RME effect and the time regions where the noise is higher can be clearly observed. The error bars are obtained by the error propagation from the photon noise. In order to include the ingress and egress shape of the transmission model we use the PyLDTk (Parviainen & Aigrain, 2015) and the PyTransit (Parviainen, 2015) Python Packages. The PyLDTk estimates the limb-darkening coefficients for a given filter and stellar properties using the library of PHOENIX stellar atmospheres (Husser et al., 2013). Then, PyTransit allows us to obtain the transit shape with the same absorption ($h$) measured during the fitting procedure and the limb-darkening coefficients assumed from PyLDTk.

The absorption that we measure from the best fit model in the line core is $0.85\pm 0.03\leavevmode\nobreak\ \%$ for CARMENES data and $0.68\pm 0.06\leavevmode\nobreak\ \%$ for HARPS-N data combined, meaning an effective radius ($R_{\lambda}(h)$) of $\sim 1.20\leavevmode\nobreak\ R_{p}$ as presented in Casasayas-Barris et al. (2018) where only HARPS-N Night 1 was available. We measure a mean standard deviation of $1.3$ and $2.3\leavevmode\nobreak\ \%$ in the residuals and a SNR of ${\sim 350}$ and ${\sim 300}$ around the H$\alpha$ position ($\pm 180\leavevmode\nobreak\ \mathrm{km/s}$) for CARMENES and HARPS-Nx3 analysis, respectively. We note that the results from both instruments are not totally consistent, as explained above. HARPS-N Night 1 and CARMENES observations are the most stable nights, with consistent values, while HARPS-N Night 2 and Night 3 present SNR variations during the transit that are then propagated to the final result when combining all three HARPS-N nights. These two nights, specially Night 3, present strong SNR variations during the observation due to the ADC problem. Although this effect is mainly removed during the analysis, in this process we apply additional steps that could introduce small over-corrections and cause variation of the final absorption level in these two nights. On the other hand, for these nights, the RME is not totally removed (see transmission light curves in Figure 8), and when combining the mid-transit signals we are averaging regions with RME residuals that decrease the absorption depth.

The differences observed in the effective radius values obtained when fitting the CLV and RME models ($R_{\lambda}$) to the 2D residual maps and those derived from the absorption in the core of the lines ($R_{\lambda}\leavevmode\nobreak\ (h)$) are possibly caused by intrinsic errors of the models. In particular, in case of hydrogen atoms, we are considering Solar abundance and LTE, but for MASCARA2 stellar-type we expect the non-LTE effects to become significant (Kubát & Korcáková, 2004). There is also the low and variable SNR difficulty added when fitting the data. We can observe that for CARMENES data, which have higher SNR, the RME measurement is better.

We measure, in all cases, a shift of the absorption line with respect the theoretical position of H$\alpha$, which could be indicative of planetary winds. However, as discussed in Section 3.2, we do not use the radial-velocity information from the instrument pipeline, and consequently we are not considering possible instrumental offsets in radial-velocity. For this reason, the winds here measured can be influenced by this effect. In all cases, before combining the residuals to obtain the transmission spectrum, we shift these residuals to the planet rest frame using the best fit $K_{p}$ parameter of each night. This value has large error bars for the noisiest nights, thus, we checked the differences when combining the spectra with different $K_{p}$ values within these $1-\sigma$ error bars, but only differences inside the noise level are noticeable. From Lund et al. (2017) we estimate that $K_{p}$ should be about $170\leavevmode\nobreak\ \pm 7\leavevmode\nobreak\ \mathrm{km/s}$. We expect all nights and lines present consistent values (within the error bars) close to this estimation, however, caused by irregularities in the maps due to small residuals and the low S/N that makes difficult for the model to distinguish between CLV, RME and absorption regions, the measured $K_{p}$ values can be far from the expectation for a single night but become consistent after the combination of more nights. Small deviations of this measurements are introduced by the not corrected reflex motion from the host star.

The absorption depths measured for $1.5\leavevmode\nobreak\ \mathrm{\AA}$ passband in the final transmission spectra and light curves are presented in Table 5. We can observe that the results of CARMENES data and the combined HARPS-N nights are consistent, while individual HARPS-N nights give more dispersed values. We also note the slightly different results obtained in the analysis presented here and in Casasayas-Barris et al. (2018) for HARPS-N Night 1. The difference comes because here we are only considering the spectra fully in-transit (the ingress and egress files are not considered when computing the in-transit sample). Additionally, we are using different methods of telluric correction, $K_{p}$ values used when shifting the spectra to the planet rest-frame, and slightly improved corrections of the CLV and RME effects.

5.2 H$\beta$ and H$\gamma$

The method applied to study the H$\beta$ ($4861.28\leavevmode\nobreak\ \mathrm{\AA}$) and H$\gamma$ ($4340.46\leavevmode\nobreak\ \mathrm{\AA}$) lines is exactly the same as for H$\alpha$. In Figure 9 we present the results for the joint analysis of all nights observed with HARPS-N. The best-fit values are shown in Table 5. Note that CARMENES data do not cover this wavelength range, and thus only HARPS-N data is presented. Same as the H$\alpha$ line, the models were computed for a grid of $[0.7-2.5]\leavevmode\nobreak\ R_{p}$ with $0.1\leavevmode\nobreak\ R_{p}$ step.

In Casasayas-Barris et al. (2018) we observed noisy features in the transmission spectrum around the H$\beta$ line analysing only Night 1. Here, with three transits and after correcting the strong CLV and RME, we observe absorption centred at the H$\beta$ laboratory value. We note that some CLV and RME residuals remain after removing the best fit model (see the residuals in Figure 9). The absorption is observable in the analysis of the individual nights, except for Night 1 where the MCMC fitting is not able to distinguish between the noise and the absorption i.e. the parameters related to the absorption profile ($h$, $K_{p}$ and $FWHM$) are not well determined. However, when combining all three transits the absorption is well recovered. The absorption measured in this case corresponds to $-0.45\pm 0.05\leavevmode\nobreak\ \%$ in the line core, which is consistent with the value measured in the transmission light curves during the transit (see Table 5). It is noticeable that the $R_{\lambda}$ values that we obtain in the fitting procedure are lower than $1\leavevmode\nobreak\ R_{p}$ in some cases.

On the other hand, we find a H$\gamma$ absorption of $-0.4\pm 0.1\leavevmode\nobreak\ \%$ in the line core when combining all nights. However, the in-transit trail in the 2D maps is not clear, possibly because the signal is less intense and remains partially hidden in the noise. This absorption is only observed when combining all data and in Night 2 when analysing each individual night. Additionally, the final transmission spectrum presents some residuals near the line.

In Table 5 we give the absorption depth values measured in the spectra and the transmission light curves of H$\beta$ and H$\gamma$ for $1.5\leavevmode\nobreak\ \mathrm{\AA}$ passband. The results of each individual night and light curves are presented in Appendix A and  C.

5.3 The Ca$\mathrm{II}$ triplet

The CaII triplet lines are located at $8498.02\leavevmode\nobreak\ \mathrm{\AA}$, $8542.09\leavevmode\nobreak\ \mathrm{\AA}$ and $8662.14\leavevmode\nobreak\ \mathrm{\AA}$. Only the transit observed with CARMENES (Night 4) is used to analyze these lines, as HARPS-N data do not cover this wavelength range. The telluric contamination observed in this region is low. The second and third lines (ordered in wavelength) are in the limit of two orders i.e. both lines are covered by two different CARMENES orders. The analysis shown here is performed using the orders where these lines are located in a most central region of the order. Although the detection of the lines, at lower SNR, is also possible in the other orders.

In Figure 10 we show the 2D maps of these lines, where the planetary absorption trail is clear, specially for the second line ($\lambda 8542\leavevmode\nobreak\ \mathrm{\AA}$) which is the strongest of the triplet. We measure absorption of $0.52\pm 0.05\leavevmode\nobreak\ \%$, $0.60\pm 0.04\leavevmode\nobreak\ \%$ and $0.55\pm 0.06\leavevmode\nobreak\ \%$ in the lines core, respectively, which correspond to an effective radius of ${\sim 1.2\leavevmode\nobreak\ R_{P}}$. The transmission spectrum around $\lambda 8662\leavevmode\nobreak\ \mathrm{\AA}$ show an smaller second absorption peak in the right side of the main absorption feature. Even the telluric contamination is low in this region, it can be caused by a small telluric residual. The transmission light curves of each individual CaII line and the transmission light curve after combining of all three lines (see Appendix C) present consistent absorption depth values in comparison with those measured in the transmission spectra (see Table 5). The transmission light curves and MCMC correlation diagrams can be observed in Appendix C and D.

As absorption residuals are observed at the CaII triplet lines, we also checked the CaII H & K lines, which are only covered with HARPS-N observations. Unfortunately, these lines are located at the beginning of the first echelle order of the CCD and the spectra are too noisy in that region to be of any value. These lines are chromospherically sensitive lines, which could present variations caused by the contrast between different regions of the stellar disk (Cauley et al., 2018). However, as expected for this stellar-type, we do not notice variations or emissions in the line cores that could be produced by stellar activity.

5.4 Na $\mathrm{I}$ doublet

Due to the rotational broadening, the NaI stellar lines expand over a large wavelength range, and a large number of telluric lines fall into the NaI lines, as we showed in Figure 4. Additionally, interstellar medium (ISM) NaI absorption at the stellar lines core is observed for all nights, affecting the region where we would expect the planetary signal. We used the LISM Kinematic Calculator (Redfield & Linsky, 2008) in order to determine the origin of the interstellar NaI lines, which gives us information about the interstellar medium clouds that our line-of-sight traverses while observing a target. Although we expect this ISM NaI absorption at the same spectral position and to be stable during the night, any small variation of this absorption during the transit could mimic a planetary atmospheric absorption. As first approximation, assuming that the ISM remains totally stable, and not correcting by the stellar reflex motion around the Barycenter of the system, the ISM effects would be removed when dividing each spectrum by the master-out spectrum.

Because of the large telluric contamination of the NaI stellar lines, we include an additional correction. After computing the ratio between each spectrum and the out-of-transit master spectrum, similarly to Snellen et al. (2010), we remove the remaining telluric residuals by fitting the time variation in each pixel by, in this case, linear splines using the UnivariateSpline class from scipy to refuse outliers. The resulting fit is compared to the data, and those points differing more than $5\leavevmode\nobreak\ \%$ from the fit values are set to $1.0$ (we note that the residuals spectra are normalized to the unit). This would not affect the possible planetary signals because they are expected to be smaller than $5\leavevmode\nobreak\ \%$, as observed in Casasayas-Barris et al. (2018).

As explained in Section 4, the stellar models computed with SME in the NaI lines region assuming solar abundance were observed to be very different compared with the observational data in the same region. For this reason, we fit the abundance, and the final models were computed using the resulting NaI abundance value ([Na/H]$=0.98$). We then fit the residuals after dividing each spectrum by the master out spectrum as we did for the H$\alpha$ in section 5.1. The MCMC fitting procedure is applied to each NaI D2 and D1 lines separately (for which we take reference positions at $5889.95\leavevmode\nobreak\ \mathrm{\AA}$ and $5895.92\leavevmode\nobreak\ \mathrm{\AA}$, respectively) and the grid of models is computed for $[0.7-2.5]\leavevmode\nobreak\ R_{p}$ with steps of $0.1\leavevmode\nobreak\ R_{p}$. The results for CARMENES and combined HARPS-N data analysis are shown in Figure 10, and the best-fit parameters are given in Table 5.

In the 2D maps, the in-transit trail seems to be observable for both D2 and D1 lines and both CARMENES and HARPS-N analysis (see Figure 10), being clearer in HARPS-N. The MCMC analysis is able to distinguish and characterise the absorption, finding a consistent value of $K_{p}$ with the one theoretically predicted. When we analyse the nights separately, the D2 line absorption can be reproduced every night except for HARPS-N Night 2, and the D1 line can be reproduced in Night 1, CARMENES data, and when all HARPS-N data are combined. The measured absorption in the line core ($h$) is consistent for both analysis, with $\sim-0.30\leavevmode\nobreak\ \%$, which corresponds to an effective radius of $\sim 1.1\leavevmode\nobreak\ R_{p}$ of absorption in both lines of the NaI doublet. These results are also consistent with the results obtained for Night 1 in Casasayas-Barris et al. (2018) applying a different analysis.

The light curves and transmission spectra present consistent absorption depth (see values in Table 5). See 2D maps of each night and light curves in Appendix A and  C. The transmission light curves are computed for each line of the doublet separately and are then combined using the weighted mean.

5.5 FeII

Three other transmission signals are observed when scanning the residuals in the wavelength range covered by both spectrographs. The position of these lines are $5018.38\pm 0.01\leavevmode\nobreak\ \mathrm{\AA}$, $5168.96\pm 0.01$, and $5316.57\pm 0.01\leavevmode\nobreak\ \mathrm{\AA}$, which we attribute to FeII $\lambda 5018.43$, $\lambda 5169.03$ and $\lambda 5316.61$, respectively. However, we note that due to the large rotation velocity of the star, it is difficult to distinguish a single line within the blended stellar lines. The FeII $\lambda 5316\leavevmode\nobreak\ \mathrm{\AA}$ line is also covered by CARMENES at very low S/N (the mean S/N of this order is around $10$). The transmission spectrum obtained shows absorption in the same place as HARPS-N ($5316.61\pm 0.02\leavevmode\nobreak\ \mathrm{\AA}$), but we also observe additional features which are not observed in HARPS-N results, and whose origin is not clear, possibly caused by the low S/N (see this results in Appendix A).

In the same way as for the other lines, we remove the RM and CLV, and find the best-fit parameters (see Figure 11 and Table 5). A faint trail during the transit is observed for each line 2D map, but stronger for FeII $\lambda 5169$ line, which present an absorption of $0.47\pm 0.05$ in the core. We note that near this line, we observe a second and fainter residual at $5172.66\pm 0.02$. The unique specie that overlaps this position is MgI at $5172.68$, but the significance of this absorption signal is very low ($1.7\sigma$) with an absorption depth of $0.07\pm 0.04\%$ in the transmission spectrum for a $0.75\mathrm{\AA}$ passband, and ${\sim-0.2\%}$ of absorption in the line core (see the results for this MgI residual in Figure 11).

We compute the transmission light curves of each FeII and MgI lines, and the FeII combined transmission light curve (see Appendix C). We then measure the absorption depth in the transmission spectra and light curves, with values presented in Table 5. The results from the transmission spectra and light curves are consistent.

6 Systematic effects

As presented in Redfield et al. (2008), we can measure the stability of the individual observations by quantifying the systematic effects with the Empirical Monte Carlo (EMC) or bootstrapping diagnostic. The EMC analysis helps us to quantify the systematic errors, since the absorption depth measurements error bars come from the the propagation of the random photon noise from Poisson statistics through the process. This analysis involves the random selection of a number of individual exposures to build the ”in-transit” sample and another random selection for the ”out-of-transit” sample. With these two samples we then extract the transmission spectrum and we measure the absorption depth of this result. This process is applied $10\leavevmode\nobreak\ 000$ times with different random samples in order to have statistical significance. With this, we are able to check how likely it is that the signal we are observing has planetary origin, or is caused by a random combination of the data.

We assume three different scenarios, as presented in Redfield et al. (2008) and other recent studies such as Wyttenbach et al. (2015) and Wyttenbach et al. (2017). The first scenario is called ”out-out” where we randomly split the out-of-transit observations in two samples with the same number of individual observations; one of these two samples will be used as the ”in-transit” sample and the other as the ”out-of-transit” sample. The second scenario is called ”in-in” and constructed in the same way, but now splitting the in-transit observations. The final scenario is the ”in-out” and represent the atmospheric absorption case. Here, the ”in-transit” and ”out-of-transit” samples are built with the in- and out-of-transit observations, respectively. For the ”in-in” and ”out-out” scenarios, we select the same number of exposures included in the ”in-transit” and in the ”out-of-transit” samples. This means that in the ”in-in” scenario this value is half the total number of real in-transit observations and in the ”out-out” is half the real out-of-transit observations. For the ”in-out” scenario the number of exposures included in each sample changes randomly, never being less than half the total number of real in- or out-of-transit exposures and always accomplishing the number ratio of real in-transit to the out-of-transit observations. We note that each sample including in-transit data is corrected by the CLV and RME effects, computed considering the best fit parameters of each line presented in Section 7.

This analysis is applied to each line here analysed. The distributions for a $1.5\leavevmode\nobreak\ \mathrm{\AA}$ passband can be observed in Figure 32. In all cases the ”out-out” and ”in-in” scenarios show distributions centered to $0\leavevmode\nobreak\ \%$ absorption while the ”in-out” scenario is centered to a different position. The clearest case is the H$\alpha$ line, where the absorption scenario is centered to $0.36\pm 0.14\leavevmode\nobreak\ \%$ and $0.40\pm 0.14\leavevmode\nobreak\ \%$ (the error corresponds to the standard deviation of the distribution) for HARPS-N and CARMENES analysis, respectively, consistent with the absorption depth measured from the final transmission spectrum (see Table 5). In case of the CaII we observe the ”in-out” samples centered on $0.10\pm 0.08$, $0.12\pm 0.11$ and $0.11\pm 0.09$ absorption, for the $\lambda 8498\leavevmode\nobreak\ \mathrm{\AA}$, $\lambda 8542\leavevmode\nobreak\ \mathrm{\AA}$ and $\lambda 8662\leavevmode\nobreak\ \mathrm{\AA}$ lines, respectively.

For the NaI doublet case, the ”in-out” scenario distributions clearly present different characteristics in comparison with the other scenarios for HARPS-N analysis. Both D2 and D1 distributions are centered to $0.07\pm 0.05\leavevmode\nobreak\ \%$. In case of CARMENES, the three different scenarios give closer results; in this case, the distributions are centered on $0.05\pm 0.09\leavevmode\nobreak\ \%$ for the D2 line and $0.12\pm 0.08\%$ for the D1 line. The absorption scenario of H$\beta$ presents a distribution centered to $0.14\pm 0.08\leavevmode\nobreak\ \%$. Finally, for H$\gamma$ the scenarios for all three distribution overlapped, being the ”in-out” distribution centered on $0.04\pm 0.09\leavevmode\nobreak\ \%$ absorption.

The FeII lines obtained with HARPS-N data show ”in-out” distributions centered to $0.03\pm 0.03$, $0.04\pm 0.04$ and $0.04\pm 0.03$, for the $\lambda 5018$, $\lambda 5169$ and $\lambda 5316$, respectively, close to the ”in-in” and ”out-out” scenarios. Similarly to H$\gamma$, the ”in-out” scenario distributions obtained for FeII $\lambda 5316$ using CARMENES data and MgI $\lambda 5173$ using HARPS-N data, show mainly indistinguishable distributions with respect ”in-in” and ”out-out” scenarios, centered to $0.09\pm 0.15$ and $0.03\pm 0.05$, respectively.

7 Discussion and conclusions

Casasayas-Barris et al. (2018) presented the analysis of one HARPS-N transit observation of MASCARA-2b, Night 1 in this paper, and reported the detection of H$\alpha$ and NaI in its atmosphere, and a hint of H$\beta$ detection. Here we have repeated the analysis after combining two more transit observations of MASCARA-2b with HARPS-N and one more with the CARMENES spectrograph. After analysing the 2D residual transit maps, the transmission spectra, the transmission light curves and EMC results, we conclude that we have detected a hint to the presence of a possible ionosphere around MASCARA2-b, including the detection of CaII, FeII, NaI, H$\alpha$, H$\beta$, and H$\gamma$ transitions. We note that, with respect to the previous study, we are using here different telluric correction and only the spectra taken in the fully in-transit time are combined when computing the transmission spectra (the ingress and egress data are excluded from this sample). We summarise the results obtained in Table 5 and Figure 12.

Focusing on the Balmer lines, the H$\alpha$ absorption by the planetary atmosphere can be observed by eye in the residuals of each single night, appearing during the transit and ranging different radial-velocities than the RME, which is also clearly observed. The planet orbital radial-velocity semi-amplitude ($K_{p}$) values measured are consistent (considering the error bars) with those from literature. For H$\beta$, on the other hand, the absorption by the planetary atmosphere is not clear if we do not combine the three transits observed with HARPS-N. But when doing so, the absorption in the transmission spectrum and the light curves is statistically significant ($5\sigma$ and $6\sigma$, respectively, for a $0.75\leavevmode\nobreak\ \mathrm{\AA}$ passband). Finally, the H$\gamma$ analysis of each individual night show mainly flat transmission spectra and light curves, but a faint H$\gamma$ signal moving at the planetary speed for MASCARA-2b results when combining all the HARPS-N data ($3\sigma$ for a $0.75\leavevmode\nobreak\ \mathrm{\AA}$ passband).

Recently, Yan & Henning (2018) detected H$\alpha$ in the atmosphere of KELT-9b, the ultra hot Jupiter with highest day-side temperature ($4600\leavevmode\nobreak\ \mathrm{K}$). They measure $1.15\pm 0.05\%$ of absorption which corresponds to an hydrogen atmosphere around $1.64$ times the radius of the planet, close to the Roche Lobe ($1.91^{+0.22}_{-0.26}R_{p}$), and an estimated mass loss rate of ${\sim 10}^{12}\mathrm{g\leavevmode\nobreak\ s^{-1}}$. In extreme cases, photoevaporation causes a transonic planetary wind and the mass loss rate could affect the planetary evolution (Salz et al., 2015). Cauley et al. (2019) reproduced the H$\alpha$ detection in KELT-9b and also reported significant absorption in H$\beta$. Here, for MASCARA-2b, the measurement of the H$\alpha$ absorption is ${\sim 0.7\%}$, which corresponds to an atmospheric extension of $1.2\leavevmode\nobreak\ R_{p}$. Considering that the lower-limit of the scale height ($H$) is $281\leavevmode\nobreak\ \mathrm{km}$ (calculated with $\mu=2.3$ and the gravity upper-limit from Lund et al. (2017)), the atmosphere is extended by $545\leavevmode\nobreak\ H$. The atmosphere is expanded due to the gas heating caused by the irradiation received from the host star. For H$\beta$, on the other hand, the $0.45\%$ of absorption corresponds to $1.15\leavevmode\nobreak\ R_{p}$ ($525\leavevmode\nobreak\ H$). Considering the mass upper-limit from Lund et al. (2017), we estimate that the Roche Lobe radius (Eggleton, 1983) of MASCARA-2b is $<3.6\leavevmode\nobreak\ R_{p}$. If we now assume a mass of 2\leavevmode\nobreak\MJ this radius decreases to $3\leavevmode\nobreak\ R_{p}$. The effective radius measured for the different atoms/ions do not reach the Roche Lobe in any of those cases. The measured FWHM of H$\alpha$ line is ${\sim 20}$ km/s which is significantly lower than that of KELT-9b (${\sim 51}$ km/s). The large FWHM of H$\alpha$ in KELT-9b is a result of high density of excited hydrogen atoms that causes the absorption to be optically thick at the line centre ($\tau\approx 57$). According to our estimation using the analytic equation (3) in Huang et al. (2017), the optical depth at the H$\alpha$ line centre of MASCARA-2b is small ($\tau\approx 2$), indicating that the density of excited hydrogen is much lower than in KELT-9b. Such a lower density also explains why the effective radius of H$\alpha$ line does not reach a high altitude like in KELT-9b. Thus, we conclude the Balmer lines absorption in MASCARA-2b are from relatively low altitude atmosphere which is far away from the Roche Lobe and therefore can not be used to estimate the mass loss rate.

For CaII, the individual lines of the triplet are clearly detected in the CARMENES night that covers this wavelength range. The trail of the strongest line ($8\sigma$ of absorption significance for $0.75\leavevmode\nobreak\ \mathrm{\AA}$) at $\lambda 8542\leavevmode\nobreak\ \mathrm{\AA}$ can be observed in the 2D maps tracing the radial-velocity change of the planet during the transit, while the other two lines of the triplet are also significant but fainter (around $5\sigma$). We note that this results from the analysis of only one transit. Using the combination of three HARPS-N nights, we also observe the presence of FeII ions in MASCARA-2b atmosphere at the $2\sigma$ level measured in the absorption depth for $0.75\leavevmode\nobreak\ \mathrm{\AA}$ passband, and an averaged significance around $7\sigma$ in the lines core. Finally, the traces of NaI planetary absorption are observable in each individual night, although the true absorption values retrieved here could be affected by the strong telluric and interstellar contamination in this region of the spectrum. The significance of these lines is around $7\sigma$ in the lines core and $2-3\sigma$ when averaging the absorption in a $0.75\leavevmode\nobreak\ \mathrm{\AA}$ passband, for both instruments.

Absorption excess in neutral Ca lines was firstly reported by Astudillo-Defru & Rojo (2013) in HD 209458b atmosphere, and the Ca+ triplet lines variation during a transit was studied for the first time by Khalafinejad et al. (2018) as indicator of stellar activity. Hoeijmakers et al. (2018) detected Fe, Fe+ and Ti+ atoms in the KELT-9b atmosphere, using the cross-correlation technique. For MASCARA-2b, we are able to resolved individual lines of Ca+ and Fe+ in transmission. The individual Fe+ lines and neutral Mg at $\lambda 5173$ observed here were also observed in the atmosphere of KELT-9b by Cauley et al. (2019). Fe+, together with Ti+, Ca+, is expected to be more abundant than its neutral form (Helling et al., 2019) due to thermal ionisation. In the high altitude atmosphere of MASCARA-2b, the temperature could be very high causing the Fe+ more abundant than neutral Fe. Our detection of these species, in a planet with equilibrium temperature almost $2000$ K colder that KELT-9b, could indicate that either a large difference between actual and equilibrium temperatures or that there are strong mechanism of transport to move atmospheric species from day-side to the planet’s terminator. Our detection of neutral Na falls in line also with the recent detection of this species in the atmosphere of WASP-76b, a hot planet with equilibrium temperature similar to MASCARA-2b (Seidel et al., 2019).

For all species, the resulting absorption signatures are blueshifted with respect their reference values, with velocities larger than $1\leavevmode\nobreak\ \mathrm{km/s}$, while the instrumental drifts for the spectrographs used here are expected to be around $10\leavevmode\nobreak\ \mathrm{m/s}$. These blueshifts can be associated to atmospheric winds. For H$\alpha$ and NaI, the blueshift measured in CARMENES and HARPS-N data are consistent. For the H$\beta$ and H$\gamma$ lines, the significance is low ($<0.9\sigma$) to determine if the measured blueshift is real. On the other hand, the FWHMs measured for the different lines are consistent for those measured with HARPS-N and CARMENES. The Balmer lines are measured to be broader than CaII, NaI and FeII lines. The line spread function (LSF) of HARPS-N is $2.7\leavevmode\nobreak\ \mathrm{km/s}$ while for CARMENES is $3.5\leavevmode\nobreak\ \mathrm{km/s}$. In all cases, the FWHM measurements are larger than the resolution element of those instruments, and are probably associated to kinematic and/or thermal broadening. In Figure 13 we relate the velocity shifts ($v_{wind}$) and FWHM measurements with respect the absorption ($R_{\lambda}(h)$). In the $v_{wind}$ versus absorption relation no obvious correlation is observed, with a mean $v_{wind}$ value around $-2.4\pm 1.0\leavevmode\nobreak\ \mathrm{km/s}$. For FWHM versus $R_{\lambda}(h)$, however, each element relation could be explained by a linear regression. Koll & Komacek (2018) studied the magnitude of the wind speeds expected in hot Jupiters, as a result of the incident flux that they receive and tidally locked spin states. Considering this study, the average wind measurement of $-2.4\pm 1.0\leavevmode\nobreak\ \mathrm{km/s}$ is in agreement with the estimated values, for a planet with the equilibrium temperature of MASCARA-2b, assuming that the dissipation is caused by shear instabilities (equation 15).

Ultra hot Jupiters lie in the temperature transition region between gas-giants and stars, some of them having very similar temperatures as the coldest stars. As such they are useful case tests to study the atmospheric chemistry, understand planetary mass loss and provide constrains for planetary theories of formation and evolution. As presented in several studies (Parmentier et al. 2018, Arcangeli et al. 2018 and Bell & Cowan 2018), we expect these extremely hot planets to present close to stellar-like atmospheric temperatures. The stellar radiation heats the atmospheric gas to temperatures considerably higher than $3000$K causing the day-side atmosphere to be composed of atoms rather than molecules in addition to an increased fraction of ions (see discussion in Helling & Rimmer (2019). Such a high irradiation leads to large day-night temperature differences and therefore we expect these planets to present hot days and cold nights, in particular for those planets with short orbital periods. For WASP-18b, for example, a planet with an effective temperature very close to MASCARA-2b, detailed theoretical studies show that we expect two very different sides, a cloudy night-side depleted of elements (which are then left to form atoms, ions, molecules and will make up the atmosphere of the planet) and a cloud-free day-side that forms a thermal ionosphere (Helling et al., 2019). As the day-side is $\approx 2500$K hotter than the night-side on WASP-18b, hydrogen is present in its atomic form (HI) such that n(HI)¿n(H2)$\gg$n(H2O) in the low-pressure, upper atmospheric regions. The same pattern holds for most atoms, including Na, K, Ca, Ti, Al, Fe, Mg, and Si. Only at higher pressures, molecule like H2, SiO, AlH are the dominating species of their element. Na, K, Ca, Ti, Al, Fe, Mg, and Si are singly ionised, with Na+, Ca+, K+, Al+, and Ti+ being more abundant than their neutral atomic/ionic form. On the night-side, the thermal ionisation is smaller and the elements remain in their neutral state of ionization or are bound in molecules. The detection of HI, NaI, Ca+ (CaII) and Fe+ (FeII) in MASCARA-2b using transmission spectroscopy data from the planet’s terminator is consistent with this theoretical works. We note, however, that transmission spectroscopy does not probe the day-side directly but rather indirectly through the observation of the terminator regions. The terminator regions are transition regions between the day- and the night-side where the temperature differences will not be as large as between the day- and night-side. The effect on the chemistry will therefore also be less pronounced. The observation of HI, NaI, Ca+ and Fe+ is therefore indicative for rather warm terminator regions. The blueshifted HI, NaI, Ca+ and Fe+ absorption suggestes that strong winds emerge from the day-side transporting hot and considerably ionised material into the respective terminator region, hence, providing us with information about the chemistry on the day-side, too. We would therefore expect the evening terminator to be geometrically more extended than the cooler morning terminator of MASCARA-2b. We note, however, that with the observations presented here we can not conclude if the heat is coming from the day-side or we are observing at very high altitudes where the atmosphere is hotter.

Thus, the theoretical picture of UHJs possessing extended ionospheres and two differentiated day-side/night-side chemistry’s is consistent with our MASCARA-2b results. However, more observational evidence is desirable to underpin this hypothesis such as CCF studies. For this planet the detection of several species at the same time has become possible with 3-4 meter class telescopes. UHJs present an exciting opportunity to study their detailed composition, with state-of-the-art spectrographs such as ESPRESSO, which will be able to characterize the H$\gamma$ features observed in MASCARA-2b and possibly other fainter elements predicted by the theoretical models. In the near future high-resolution spectrographs, such as HIRES for the ELT, will be able to disentangle the complex picture of UHJ atmospheres through the simultaneous detection of multiple species in their atmospheres.

Acknowledgements.

The authors would like to thank Dr Sven Buder for his useful comments and discussion about the Na abundances. Based on observations made with the Italian Telescopio Nazionale Galileo (TNG) operated on the island of La Palma by the Fundación Galileo Galilei of the INAF (Istituto Nazionale di Astrofisica) at the Spanish Observatorio del Roque de los Muchachos of the Instituto de Astrofisica de Canarias. CARMENES is an instrument for the Centro Astronómico Hispano-Alemán de Calar Alto (CAHA, Almería, Spain) funded by the German Max-Planck-Gesellschaft (MPG), the Spanish Consejo Superior de Investigaciones Científicas (CSIC), the European Union through FEDER/ERF FICTS-2011-02 funds, and the members of the CARMENES Consortium. This work is partly financed by the Spanish Ministry of Economics and Competitiveness through project ESP2016-80435-C2-2-R. G.C. acknowledges the support by the National Natural Science Foundation of China (Grant No. 11503088, 11573073, 11573075) and by the project “Technology of Space Telescope Detecting Exoplanet and Life” from National Defense Science and Engineering Bureau civil spaceflight advanced research project (D030201). F.Y. acknowledges the support of the DFG priority program SPP 1992 ”Exploring the Diversity of Extrasolar Planets (RE 1664/16-1). A.F. acknowledges the support from JSPS KAKENHI Grant Number JP17H04574. This work is partly supported by JSPS KAKENHI Grant Numbers JP18H01265 and 18H05439, and JST PRESTO Grant Number JPMJPR1775. This article is based on observations made with the MuSCAT2 instrument, developed by ABC, at Telescopio Carlos Sánchez operated on the island of Tenerife by the IAC in the Spanish Observatorio del Teide. This work made use of PyAstronomy and of the VALD database, operated at Uppsala University, the Institute of Astronomy RAS in Moscow, and the University of Vienna.

Appendix A Individual transmission spectra

We present here the transmission spectra obtained for each individual night.

A.1 H$\beta$

A.2 H$\gamma$

A.3 NaI D2

A.4 NaI D1

A.5 FeII

Appendix B Best fit parameters and absorption depth values

We present here the best-fit parameters obtained with the MCMC fitting procedure and the absorption depth values measured in the transmission spectra and light curves.

B.1 MCMC best-fit parameters

B.2 Absorption depth values

Appendix C Individual transmission light curves

We present here the transmission light curves obtained for each specie analysed here.

C.1 H$\beta$

C.2 H$\gamma$

C.3 CaII triplet

C.4 NaI doublet

C.5 FeII

C.6 MgI

Appendix D Correlation diagrams of the best-fit models

We present here the correlation diagrams of the best-fit CLV, RME and absorption model for the different nights and wavelength regions.

D.1 H$\alpha$

D.2 H$\beta$ and H$\gamma$

D.3 CaII

D.4 NaI

D.5 FeII

D.6 MgI

Appendix E Empirical Monte-Carlo distributions