HAT-P-67b: An Extremely Low Density Saturn Transiting an F-Subgiant Confirmed via Doppler Tomography

TL;DR

HAT-P-67b is an extremely low-density Saturn-sized exoplanet orbiting a rapidly rotating F-subgiant, confirmed through Doppler tomography, with its properties challenging traditional mass measurement techniques due to stellar rotation.

Contribution

This study reports the discovery and confirmation of HAT-P-67b using Doppler tomography, demonstrating a method to validate planets around rapidly rotating stars where radial velocities are limited.

Findings

HAT-P-67b has a radius of approximately 2.09 Jupiter radii.

The planet's mass is constrained to be less than 0.59 Jupiter masses.

The planet's orbit is aligned within 12 degrees of the host star's spin.

Abstract

We report the discovery of HAT-P-67b, a hot-Saturn transiting a rapidly rotating F-subgiant. HAT-P-67b has a radius of Rp = 2.085 -0.071/+0.096 RJ,, orbiting a M* = 1.642 -0.072/+0.155 Msun, R* = 2.546 -0.084/+0.099 Rsun host star in a ~4.81-day period orbit. We place an upper limit on the mass of the planet via radial velocity measurements to be Mp < 0.59 MJ, and lower limit of > 0.056 MJ by limitations on Roche lobe overflow. Despite being a subgiant, the host star still exhibits relatively rapid rotation, with a projected rotational velocity of v sin I* = 35.8 +/- 1.1 km/s, making it difficult to precisely determine the mass of the planet using radial velocities. We validated HAT-P-67b via two Doppler tomographic detections of the planetary transit, which eliminated potential eclipsing binary blend scenarios. The Doppler tomographic observations also confirmed that HAT-P-67b has an orbit that is aligned to within 12 degrees, in projection, with the spin of its host star. HAT-P-67b receives strong UV irradiation, and is amongst the one of the lowest density planets known, making it a good candidate for future UV transit observations to search for an extended hydrogen exosphere.

Tables

Table 1. Table 1 Differential photometry of HAT-P-67

BJD	Mag (Raw)   aa This table is available in a machine-readable form in the online journal. A portion is shown here for guidance regarding its form and content. Raw, EPD, and TFA magnitudes are presented for HATNet light curves. The detrending and potential blending may cause the HATNet transit to be shallower than the true transit in the EPD and TFA light curves. This is accounted for in the global modelling by the inclusion of a third light factor. Follow-up light curves have been treated with EPD simultaneous to the transit fitting. Pre-EPD magnitudes are presented for the follow-up light curves.	Mag (EPD)	Mag (TFA)	$\sigma$ Mag	Instrument	Filter
2454521.99042	9.39859	9.72414	9.72271	0.00303	HATNet	R
2454521.99445	9.39572	9.72048	9.73067	0.00325	HATNet	R
2454521.99854	9.39262	9.72005	9.71695	0.00276	HATNet	R
2454522.00264	9.38941	9.7199	9.72327	0.00344	HATNet	R
2454522.00673	9.41475	9.74011	9.73327	0.00328	HATNet	R

Table 2. Table 2 Summary of photometric observations

Facility	Date(s)	Number of Images   aa Outlying exposures have been discarded.	Cadence (s)   bb Median time difference between points in the light curve. Uniform sampling was not possible due to visibility, weather, pauses.	Filter
HATNet	2005 Jan – 2005 Jul	4050	328	$I$
HATNet	2008 Feb – 2008 Aug	4518	246	Cousins $R$
FLWO 1.2 m KeplerCam	2011 Apr 15	730	24	Sloan-$i$
FLWO 1.2 m KeplerCam	2011 May 19	509	44	Sloan-$i$
FLWO 1.2 m KeplerCam	2011 Jun 07	801	29	Sloan-$i$
FLWO 1.2 m KeplerCam	2012 May 28	730	34	Sloan-$i$
FLWO 1.2 m KeplerCam	2013 Apr 25	960	24	Sloan-$i$
FLWO 1.2 m KeplerCam	2013 May 24	361	24	Sloan-$z$

Table 3. Table 3 KECK-HIRES relative radial velocities and bisector span measurements of HAT-P-67

BJD	RV   aa Internal errors excluding the component of astrophysical/instrumental jitter considered in Section 3. Bisector spans (BS) are given where available.	$\sigma$ RV	BS	$\sigma$ BS
(UTC)	$\left(m{s}^{-1}\right)$	$\left(m{s}^{-1}\right)$	$\left(m{s}^{-1}\right)$	$\left(m{s}^{-1}\right)$
2455696.8366	-105	28	24	11
2455696.88382	-151	31	67	17
2455697.833	31	25	-58	9
2455698.92918	93	23	-53	8
2455699.83162	63	25	46	10
2455700.88206	1	22	30	13
2455704.84352	18	23	-2	13
2455705.86007	100	22	-34	8
2455706.83882	8	22	14	11
2455707.85238	35	22	-17	13
2455853.70871	-6	28	…	…
2455945.15236	-89	24	…	…
2455997.02884	-17	38	…	…
2455017.0082	58	22	…	…
2455042.88956	-65	26	…	…
2455043.9989	-25	29	…	…
2455044.94822	-63	26	…	…
2455048.86097	-112	30	…	…
2455107.71733	241	40	-3	28

Table 4. Table 4 Summary of spectroscopic observations

Telescope/Instrument	Date Range	Number of Observations	Resolution	Observing Mode
NOT 2.5 m/FIES	2009 Aug 4 – 2009 Oct 10	5	67000	RECON RV
FLWO 1.5 m/TRES	2011 Apr 17 – 2011 Apr 20	2	44000	RECON RV
KECK 10 m/HIRES	2009 Jul 04 – 2012 Mar 10	19	55000	RVaa High resolution spectra to obtain stellar atmospheric parameters and high precision radial velocities
FLWO 1.5 m/TRES	2016 Apr 17	14	44000	Transitbb High resolution in-transit spectra to detect the Doppler tomographic signal of the planet
FLWO 1.5 m/TRES	2016 May 16	16	44000	Transitbb High resolution in-transit spectra to detect the Doppler tomographic signal of the planet

Table 5. Table 5 Stellar parameters for HAT-P-67

Parameter	Circular Fit Geneva	Eccentric Fit Geneva	Circular Fit Dartmouth	Eccentric Fit Dartmouth
Catalogue Information
Tycho-2	3084-533-1
GSC	03084-00533
2MASS	J17062656+4446371
GAIA	1358614978835493120
GAIA RA (J2015)	17:06:26.574
GAIA DEC (J2015)	+44:46:36.794
GAIA ${\mu }_{\alpha }$ $\left(\mathrm{mas}{\mathrm{yr}}^{-1}\right)$	$9.32\pm 0.88$
GAIA ${\mu }_{\delta }$ $\left(\mathrm{mas}{\mathrm{yr}}^{-1}\right)$	$18.5\pm 1.2$
GAIA Parallaxaa A correction of $-0.325\pm 0.062$ mas has been applied to the GAIA DR1 parallax, as per Stassun & Torres (2016). $\left(\mathrm{mas}\right)$	$2.92\pm 0.23$
Spectroscopic properties bb Derived from the global modelling described in Section 3.4, co-constrained by spectroscopic stellar parameters and the GAIA parallax. cc These stellar parameters are heavily constrained by Gaussain priors about their derived values from the Keck-HIRES iodine-free spectrum using the Stellar Parameter Classification (SPC) pipeline (Buchhave et al., 2012).
$T_{\mathrm{eff}\star }$ (K)	${6406}_{-61}^{+65}$	${6408}_{-65}^{+63}$	${6406}_{-63}^{+58}$	${6414}_{-59}^{+69}$
$\left[\mathrm{Fe}/H\right]$	$-0.08\pm 0.05$	$-0.08\pm 0.05$	$-{0.07}_{-0.05}^{+0.04}$	$-0.08\pm 0.05$
$v\sin I_{\star }$ ($\mathrm{km}{s}^{-1}$)	$35.8\pm 1.1$	$35.8\pm 1.1$	${33.2}_{-1.2}^{+1.5}$	${33.9}_{-1.3}^{+1.2}$
Photometric properties
GALEX FUV (AB mag)	$19.759\pm 0.137$
GALEX NUV (AB mag)	$14.251\pm 0.007$
GAIA $g$ (mag)	9.94
APASS $B$ (mag)	$10.682\pm 0.010$
APASS $g^{\prime }$ (mag)	$10.351$
APASS $V$ (mag)	$10.069\pm 0.016$
APASS $r^{\prime }$ (mag)	$10.010$
TASS $I$ (mag)	$9.518\pm 0.048$
2MASS $J$ (mag)	$9.145\pm 0.021$
2MASS $H$ (mag)	$8.961\pm 0.019$
2MASS $K_s$ (mag)	$8.900\pm 0.019$
Derived propertiesbb Derived from the global modelling described in Section 3.4, co-constrained by spectroscopic stellar parameters and the GAIA parallax.
$M_{\star }$ ($M_{☉}$)	${1.642}_{-0.072}^{+0.155}$	${1.73}_{-0.13}^{+0.21}$	$1.43\pm 0.05$	${1.38}_{-0.05}^{+0.05}$
$R_{\star }$ ($R_{☉}$)	${2.546}_{-0.084}^{+0.099}$	${2.71}_{-0.39}^{+0.48}$	${2.389}_{-0.038}^{+0.040}$	${2.13}_{-0.14}^{+0.17}$
$\log g_{\star }$ (cgs)	${3.854}_{-0.023}^{+0.014}$	${3.800}_{-0.080}^{+0.106}$	${3.837}_{-0.011}^{+0.009}$	${3.932}_{-0.060}^{+0.035}$
$L_{\star }$ ($L_{☉}$)	${8.68}_{-0.86}^{+1.50}$	${8.3}_{-1.9}^{+4.0}$	${8.62}_{-0.50}^{+0.57}$	${6.8}_{-0.9}^{+1.2}$
$M_V$ (mag)	${2.50}_{-0.23}^{+0.13}$	${2.57}_{-0.37}^{+0.29}$	${2.403}_{-0.063}^{+0.083}$	${2.67}_{-0.17}^{+0.15}$
$M_K$ (mag,ESO)	${1.26}_{-0.34}^{+0.15}$	${1.36}_{-0.39}^{+0.25}$	${1.304}_{-0.045}^{+0.046}$	${1.56}_{-0.18}^{+0.14}$
$A_V$ (mag)
Age (Gyr)	${1.24}_{-0.22}^{+0.27}$	${1.00}_{-0.41}^{+0.21}$	${2.83}_{-0.19}^{+0.22}$	${3.04}_{-0.27}^{+0.31}$
Distance (pc)	${320}_{-14}^{+48}$	${322}_{-19}^{+35}$	${335}_{-7}^{+7}$	${297}_{-18}^{+26}$

Table 6. Table 6 Orbital and planetary parameters

Parameter	Circular Fit Geneva	Eccentric Fit Geneva	Circular Fit Dartmouth	Eccentric Fit Dartmouth
Light curve parameters
$P$ (days)	${4.8101025}_{-0.00000033}^{+0.00000043}$	${4.8101038}_{-0.00000037}^{+0.00000054}$	${4.8101017}_{-0.00000030}^{+0.00000034}$	${4.8101082}_{-0.00000051}^{+0.00000052}$
$T_c$ ($\mathrm{BJD}$) aa $T_c$: Reference epoch of mid transit that minimizes the correlation with the orbital period. BJD is calculated from UTC. $T_{14}$: total transit duration, time between first to last contact; $T_{12}=T_{34}$: ingress/egress time, time between first and second, or third and fourth contact.	${2455961.38467}_{-0.00064}^{+0.00076}$	${2455961.38472}_{-0.00082}^{+0.00090}$	${2455961.38465}_{-0.00065}^{+0.00074}$	${2455961.3852}_{-0.0010}^{+0.0008}$
$T_{14}$ (days) aa $T_c$: Reference epoch of mid transit that minimizes the correlation with the orbital period. BJD is calculated from UTC. $T_{14}$: total transit duration, time between first to last contact; $T_{12}=T_{34}$: ingress/egress time, time between first and second, or third and fourth contact.	$0.2912\pm 0.0019$	${0.308}_{-0.031}^{+0.029}$	$0.2910\pm 0.0015$	${0.257}_{-0.010}^{+0.019}$
$T_{12}=T_{34}$ (days) aa $T_c$: Reference epoch of mid transit that minimizes the correlation with the orbital period. BJD is calculated from UTC. $T_{14}$: total transit duration, time between first to last contact; $T_{12}=T_{34}$: ingress/egress time, time between first and second, or third and fourth contact.	$0.0229\pm 0.0010$	$0.0246\pm 0.0027$	${0.02330}_{-0.00055}^{+0.00030}$	${0.0213}_{-0.0014}^{+0.0015}$
$a/R_{\star }$	${5.691}_{-0.124}^{+0.057}$	${5.34}_{-0.46}^{+0.61}$	${5.659}_{-0.061}^{+0.066}$	${6.34}_{-0.42}^{+0.30}$
$R_p/R_{\star }$	$0.0834\pm 0.0017$	${0.084}_{-0.0020}^{+0.0019}$	${0.0821}_{-0.0009}^{+0.0013}$	${0.0846}_{-0.0018}^{+0.0016}$
$b\equiv a\cos i/R_{\star }$	${0.12}_{-0.08}^{+0.12}$	${0.12}_{-0.08}^{+0.12}$	${0.214}_{-0.045}^{+0.023}$	${0.20}_{-0.12}^{+0.11}$
$i$ (deg)	${88.8}_{-1.3}^{+1.1}$	$88.9\pm 1.6$	${88.37}_{-0.57}^{+0.61}$	${88.2}_{-1.1}^{+1.3}$
$\left|\lambda \right|$ (deg)	${2.9}_{-4.9}^{+6.4}$	${2.5}_{-4.6}^{+5.8}$	$-{1.6}_{-4.6}^{+3.9}$	${2.3}_{-6.4}^{+6.6}$
Limb-darkening coefficients bb Values for a quadratic law given separately for each of the filters with which photometric observations were obtained. These values were adopted from the tabulations by Claret & Bloemen (2011) according to the spectroscopic (SPC) parameters listed in Table 5. The limb darkening coefficients are held fixed during the global modelling.
$a_r$ (linear term)	0.2497
$b_r$ (quadratic term)	0.3765
$a_I$	0.1701
$b_I$	0.3744
$a_i$	0.1897
$b_i$	0.3747
$a_z$	0.1397
$b_z$	0.3661
RV parameters
$K$ ($m{s}^{-1}$)
$e\cos \omega$		$-{0.21}_{-0.14}^{+0.15}$		$-{0.03}_{-0.22}^{+0.20}$
$e\sin \omega$		${0.027}_{-0.11}^{+0.10}$		$-{0.150}_{-0.055}^{+0.075}$
$e$		$0.24\pm 0.12$		${0.22}_{-0.08}^{+0.12}$
$\omega$		${172}_{-43}^{+31}$		${105}_{-66}^{+46}$
RV jitter ($m{s}^{-1}$)cc This jitter was added linearly to the RV uncertainties for each instrument such that ${\chi }^2/\mathrm{dof}=1$ for the observations from that instrument.	59	58	59	59
Systemic RV ($\mathrm{km}{s}^{-1}$)dd The systemic RV for the system as measured relative to the telluric lines	$-1.4\pm 0.5$
Planetary parameters
$M_p$ ($M_J$)ee The mass measurement is quoted as the median of the posterior, with the uncertainties defined as the 68 percentile region.	${0.34}_{-0.19}^{+0.25}$	${0.49}_{-0.22}^{+0.27}$	${0.33}_{-0.17}^{+0.22}$	${0.29}_{-0.19}^{+0.24}$
$R_p$ ($R_J$)	${2.085}_{-0.071}^{+0.096}$	${2.25}_{-0.23}^{+0.20}$	${1.975}_{-0.038}^{+0.045}$	${1.78}_{-0.10}^{+0.14}$
${\rho }_p$ ($g{\mathrm{cm}}^{-3}$)	${0.052}_{-0.028}^{+0.039}$	${0.058}_{-0.025}^{+0.039}$	${0.058}_{-0.030}^{+0.039}$	${0.065}_{-0.044}^{+0.062}$
$\log g_p$ (cgs)	${2.32}_{-0.34}^{+0.24}$	${2.41}_{-0.25}^{+0.20}$	${3.837}_{-0.011}^{+0.009}$	${2.36}_{-0.47}^{+0.28}$
$a$ (AU)	${0.06505}_{-0.00079}^{+0.00273}$	${0.0663}_{-0.0014}^{+0.0016}$	${0.062844}_{-0.00049}^{+0.00053}$	${0.061994}_{-0.00072}^{+0.00068}$
$T_{\mathrm{eq}}$ (K)	$1903\pm 25$	${1963}_{-99}^{+85}$	${1903}_{-21}^{+19}$	${1803}_{-43}^{+62}$
$\Theta$ff The Safronov number is given by $\Theta =\frac{1}{2}{\left(V_{\mathrm{esc}}/V_{\mathrm{orb}}\right)}^2=\left(a/R_p\right)\left(M_p/M_{\star }\right)$ (see Hansen & Barman, 2007).	${0.0138}_{-0.0075}^{+0.0099}$	${0.0178}_{-0.0077}^{+0.0098}$	${0.015}_{-0.015}^{+0.025}$	${0.015}_{-0.010}^{+0.013}$
$⟨F⟩$ (${10}^9$$\mathrm{erg}{s}^{-1}{\mathrm{cm}}^{-2}$) gg Incoming flux per unit surface area, averaged over the orbit.	${2.74}_{-0.17}^{+0.19}$	${2.57}_{-0.46}^{+0.45}$	${2.98}_{-0.13}^{+0.14}$	${2.41}_{-0.25}^{+0.33}$

Full text

1. Introduction

Finding well-characterized planets in a variety of environments is key to understanding the processes that govern planet formation and evolution. Planets orbiting high mass stars are likely born in high mass protoplanetary disks (e.g. Muzerolle et al., 2003; Natta et al., 2006), environments that may yield higher planet occurrence rates (e.g. Johnson et al., 2010; Bowler et al., 2010) and higher-mass planets (e.g. Lovis & Mayor, 2007; Jones et al., 2014) than around solar type stars. Planets around early type stars also receive higher incident flux over their lifetimes, which in turn make them anchor-points in the planet mass-radius-equilibrium temperature relationships (e.g. Béky et al., 2011; Enoch et al., 2012).

However, only 1% of known transiting planets orbit stars more massive than $1.5\,M_{\sun}$. Early type stars have larger radii, resulting in shallower transit depths for any planets; they are also more likely to have rotationally blended spectral lines due to the lack of magnetic braking over the main-sequence lifetime, making traditional radial-velocity confirmation techniques more difficult. One successful strategy is to conduct radial-velocity surveys of ‘retired A-stars’ – stars that have evolved off the main-sequence and spun down enough to exhibit sharp spectroscopic lines that enable precise radial-velocity measurements. These surveys have been extremely successful, yielding 122 planetary systems to date111Choosing host stars with $\log{g}<4.0$, from NASA Exoplanet Archive, July 2016 (e.g. Johnson et al., 2007; Wittenmyer et al., 2011; Jones et al., 2014). Recently, transit surveys have also been successful in discovering planets around high mass stars. These include planets around subgiants and giants whose shallow transits were identified by Kepler (e.g. Kepler-56b,c Huber et al. 2013, Kepler-96b Lillo-Box et al. 2014, Kepler-432b Quinn et al. 2015; Ciceri et al. 2015, KOI-206b, KOI-680b Almenara et al. 2015, and K2-39b Van Eylen et al. 2016), and hot-Jupiters around main-sequence A-stars confirmed via Doppler tomography (WASP-33b Collier Cameron et al. 2010b, Kepler-13b Szabó et al. 2011; Shporer et al. 2011; Johnson et al. 2014, HAT-P-57b Hartman et al. 2015, and KELT-17b Zhou et al. 2016b).

In this paper, we present the discovery of HAT-P-67b, a Saturn-mass planet found to transit an F-subgiant by the HATNet survey (Bakos et al., 2004). Despite the evolved status of HAT-P-67, the host star still exhibits a rapid rotation rate of $v\sin{I_{\star}}=35.8\pm 1.1\,\rm km\,s^{-1}$, making precise radial velocities difficult to obtain. Eventual confirmation was achieved via a detection of the Doppler tomographic shadow of the planet during transit. When a planet transits a rapidly rotating star, it successively blocks parts of the rotating stellar disk, causing an asymmetry in the observed spectral line profiles. At low rotational velocities, the asymmetry can be measured by the Holt-Rossiter-McLaughlin effect (Holt, 1893; Rossiter, 1924; McLaughlin, 1924). At higher rotational velocities, the shadow of the planet can be resolved in the broadened stellar spectroscopic lines (e.g. Collier Cameron et al., 2010a, b). A detection of the Doppler tomographic signal, at a depth and width that are in agreement with the photometric transit, eliminates eclipsing binary blend scenarios that may mimic transiting planet signals. Further radial velocity measurements can then provide an upper-limit mass constraint of the orbiting companion. If the mass can be constrained to less than that of brown dwarfs, the transiting object is confirmed to be a planet.

2. Observations

2.1. Photometry

The transits of HAT-P-67b were first detected with the HATNet survey (Bakos et al., 2004). HATNet employs a network of small, wide field telescopes located at the Fred Lawrence Whipple Observatory (FLWO) in Arizona, and the Mauna Kea Observatory (MKO) in Hawaii, USA, to photometrically monitor selected $8\times 8^{\circ}$ fields of the sky. A total of 4050 $I$ band observations were taken by HAT-5 and HAT-8 over 2005 Jan – July, and an additional 4518 observations were obtained in the Cousins $R$ band using HAT-5, HAT-7, and HAT-8 telescopes between 2008 Feb – Aug. The data reduction follows Bakos et al. (2010). Light curves were produced via aperture photometry, and detrended with External Parameter Decorrelation (EPD, Bakos et al., 2007) and Trend Filtering Algorithm (TFA, Kovács et al., 2005). The Box-Fitting Least Squares (BLS, Kovács et al., 2002) analysis revealed the periodic transits of the planet candidate. The discovery light curve of HAT-P-67b is shown in Figure 1, and the photometry presented in Table 1.

To better characterize the planetary properties, follow-up photometry of the transits were obtained using KeplerCam on the FWLO 1.2 m telescope. KeplerCam is a $4\mathrm{K}\times 4\mathrm{K}$ CCD camera with a pixel scale of $0\farcs 672\,\mathrm{pixel}^{-1}$ at $2\times 2$ pixel binning. The photometry was reduced as per Bakos et al. (2010). A full transit was observed in the Sloan-$i$ band on 2012 May 28, and five partial transits observed on 2011 Apr 15, 2011 May 19, 2011 Jun 07, 2013 Apr 25 in the Sloan-$i$, and 2013 May 24 in the Sloan-$z$ band. The light curves and best fit models are shown in Figure 2, and the data presented in Table 1.

2.2. Spectroscopy

Spectroscopic observations of HAT-P-67 were carried out using the FIber-fed Echelle Spectrograph (FIES), the Tillinghast Reflector Echelle Spectrograph (TRES), and the High Resolution Echelle Spectrometer (HIRES). The observations are summarized in Table 4 and described below.

Initial spectroscopic characterization of HAT-P-67b was obtained with the FIES instrument (Telting et al., 2014) on the 2.5 m Nordic Optical Telescope (NOT). FIES is a fiber fed high resolution echelle spectrograph with a resolution of $\lambda/\Delta\lambda\equiv R=67000$ and spectral coverage of $3700-7300\,\AA$. Four FIES radial velocities were obtained over the 2009 Aug – Oct period. The observations were obtained and reduced as per the procedure from Buchhave et al. (2010). No radial velocity variation was detected with the FIES observations, with a scatter of $200\,\mathrm{m\,s}^{-1}$ over the four observations.

Additional observations were obtained with the TRES instrument (Fűrész, 2008) on the FLWO 1.5 m telescope. TRES is a fiber fed echelle with a spectral resolution of $R=44000$, over the spectral region of $3850-9100\,\AA$. Radial velocities and spectral classifications are measured from each spectrum as per Buchhave et al. (2012). Each TRES observation consists of three exposures combined together for cosmic-ray removal, and wavelength-calibrated by Th-Ar lamp exposures that bracket each set of three exposures. Two TRES observations at phase quadrature were taken in 2011 Apr 17 and 2011 Apr 20, with signal-to-noise at the Mg b lines of $\sim 100$ per resolution element. The velocity difference between the two observations was $80\,\mathrm{m\,s}^{-1}$, with a per-point uncertainty of $100\,\mathrm{m\,s}^{-1}$. As such, the the FIES and TRES observations showed that any companion orbiting HAT-P-67 must be sub-brown dwarf in mass.

In addition, we observed two partial spectroscopic transits of HAT-P-67b, on 2016 Apr 17 and 2016 May 16, with TRES to detect the Doppler tomographic shadow of the planet. These observations were performed as per the strategy described in Zhou et al. (2016a). A set of time series spectra, at 900 s cadence, where collected on both nights. The Doppler tomographic analysis for these two transit sets are described in Section 3.1.

To constrain the mass of the companion, we obtained spectroscopic observations from HIRES on the 10 m KECK telescope (Vogt et al., 1994) at MKO over the 2009 Jul – 2012 Mar period. A total of 19 observations were obtained through the $I_{2}$ cell to provide precise radial velocities. An additional $I_{2}$-free observation was obtained to provide a template for the radial-velocity measurements. The instrument was set up to use the C2 decker, which provides a $14\arcsec\times 0\farcs 861$ slit, yielding a spectral resolution of $R=48000$. The radial velocities were measured as per Butler et al. (1996), and the bisector spans calculated as per Torres et al. (2007). The high signal-to-noise HIRES observations provide the best constraints on the radial velocities of HAT-P-67, and were used in the global analysis in Section 3.4. The radial velocities from HIRES are plotted in Figure 3 and presented in Table 3.

3. Analysis

3.1. Doppler tomographic detection of the planetary transit

The significant rotational broadening of HAT-P-67 allows us to detect the spectroscopic transit of the planet via Doppler tomography (Collier Cameron et al., 2010a, b). Two sets of transit spectroscopy were obtained for HAT-P-67b with TRES. The TRES spectra were processed as per the procedure laid out in Zhou et al. (2016a): the broadening profiles were derived via a least-squares deconvolution (LSD) of the observed spectra against a non-rotating stellar template (as per Donati et al., 1997). Synthetic template spectra were generated using the SPECTRUM (Gray & Corbally, 1994) spectral synthesis program, using the ATLAS9 model atmospheres (Castelli & Kurucz, 2004). The synthetic templates were generated at the same $T_{\rm eff}$, $\log{g}$ and $\rm[Fe/H]$ as HAT-P-67, with no line broadening imposed. A broadening profile was derived for each spectrum and subtracted from the average out-of-transit profile, revealing the planetary transit signal (Figure 4). We model the rotational profiles and the planetary signal as part of our global analysis, described in Section 3.4.

3.2. Stellar parameters

Stellar atmospheric parameters of HAT-P-67 were derived from the 32 TRES spectra using the Stellar Parameter Classification pipeline (SPC, Buchhave et al., 2012). We first run SPC to retrieve an initial estimate of the stellar atmospheric parameters. These are then incorporated in a first run of the global modeling and isochrone retrieval analysis described later in Section 3.4. We then re-run SPC with the stellar surface gravity $\log{g}$ fixed to that measured from the transit duration in the global analysis (Section 3.4, with $\log{g}=3.854_{-0.023}^{+0.014}$) to provide updated, and better constrained $T_{\rm eff}$ and $\rm[Fe/H]$ values. We find that HAT-P-67 is consistent with an F-subgiant, of effective temperature $T_{\mathrm{eff}}=6406\pm 62\,\mathrm{K}$, metallicity $\mathrm{[m/H]}=-0.08\pm 0.05$, and projected rotational velocity $v\sin I_{\star}=36.5\pm 0.3\,\mathrm{km\,s}^{-1}$. Similarly, running SPC on the four FIES spectra yield $T_{\mathrm{eff}}=6380\pm 50\,\mathrm{K}$, $\log g=3.91\pm 0.10$, $\mathrm{[m/H]}=-0.05\pm 0.08$, $v\sin I_{\star}=38\,\mathrm{km\,s}^{-1}$, consistent with the interpretation that HAT-P-67 is an F subgiant. Since an accurate $v\sin{I_{\star}}$ measurement is vital to correctly modeling the Doppler tomographic signal, we also use the set of time-series TRES spectra to measure the $v\sin{I_{\star}}$ of HAT-P-67. Following Zhou et al. (2016a), the broadening kernel for each spectrum is modeled by a rotational kernel, with width of $v\sin{I_{\star}}$, and a Gaussian kernel to account for macroturbulence and instrumental broadening, finding $v\sin{I_{\star}}=30.9\pm 2.0\,\rm km\,s^{-1}$, and macroturbulence of $9.22\pm 0.5\,\mathrm{km\,s}^{-1}$. The uncertainties are estimated from the standard deviation scatter between exposures. The difference between the $v\sin{I_{\star}}$ measured via SPC and that from the rotational profile can be partially attributed to the inclusion of macroturbulence.

3.3. GAIA parallax

HAT-P-67 is included in the Tycho-GAIA-Astrometric-Catalogue in the first data release (DR1) of GAIA (Lindegren et al., 2016), which measured a parallax of $2.60\pm 0.23$ mas. Several literature investigations have pointed out a systematic under-estimation in the DR1 parallaxes, as per separate studies via eclipsing binaries (Stassun & Torres, 2016), close-by Cepheids (Casertano et al., 2016), asteroseismic distances (Silva Aguirre et al., 2016), and comparison with existing parallaxes of solar neighborhood stars (Jao et al., 2016). We adopt the correction offered in Stassun & Torres (2016) of $-0.325\pm 0.062$ mas to the DR1 parallax of HAT-P-67, arriving at an adopted parallax value of $2.92\pm 0.23$ mas, and corresponding astrometric distance measurement of $342\pm 27$ pc. This parallax measurement is used to co-constrain the stellar parameters during the global modeling in Section 3.4.

3.4. Global fitting and derived planet parameters

We perform a global analysis of the HATNet discovery light curves, follow-up transit light curves, KECK-HIRES $I_{2}$ radial velocities, and the TRES Doppler tomographic signal, co-constrained by stellar isochrones and the GAIA distance measurement. The transits are modelled according to Mandel & Agol (2002), with the transit shape defined by the transit centroid time $T_{0}$, star-planet distance $a/R_{\star}$, planet-star radius ratio $R_{p}/R_{\star}$, and transit inclination $i$. Individual quadratic limb darkening parameters are assigned to each light curve (interpolated from Claret & Bloemen, 2011), and fixed throughout the fitting. Separate dilution factors are allowed for the HATNet $I$ and $R_{C}$ band light curves to account for any distortions to the light curve shape from the TFA detrending process. The follow-up light curves are simultaneously detrended against instrumental parameters describing the $X$, $Y$ pixel centroids of the target star, background flux, and target airmass. The radial velocities are described by an arbitrary offset $\gamma$ and orbital semi-amplitude $K$. The orbital eccentricity parameters $e\cos\omega$ and $e\sin\omega$ are also included when eccentricity is allowed to vary. The Doppler tomographic signal is modeled as per Zhou et al. (2016b), via a 2D integration of the stellar surface covered by the planet. The free parameters describing the Doppler tomography effect include the projected spin-orbit angle $\lambda$, and the projected rotational broadening velocity $v\sin I_{\star}$. Note that we do not account for the broadening of the planetary shadow due to the motion of the planet during an exposure; the blurring of the planetary shadow during an exposure ($2\,\rm km\,s^{-1}$) is smaller than the width of the shadow ($7.2\,\rm km\,s^{-1}$), but is not an insignificant effect. We also allow the effective temperature $T_{\mathrm{eff}}$, metallicity [M/H], and the apparent $K$-band magnitude to be iterated, though heavily constrained about their spectroscopic and photometric values. At each step, we derive a stellar density $\rho_{\star}$ from the transit duration as per Seager & Mallén-Ornelas (2003); Sozzetti et al. (2007), and query the stellar isochrones to derive a distance modulus. Isochrone interpolation is performed at each step using the gradient boosting regression algorithm implemented in scikit-learn. This distance modulus is compared to the actual distance as measured from the GAIA parallax, with the difference applied as a penalty on the likelihood function.

The rapid rotation rate of HAT-P-67 can introduce a bias in the isochrone-derived parameters for the system. For stars with radiative envelopes, the convective core overshoot and mixing length parameters are different to that of non-rotating stellar models, with the overall effect of lengthening the main-sequence lifetime (e.g. Meynet & Maeder, 2000). We adopt the Geneva 2D stellar evolution models (Ekström et al., 2012), which account for the effects of rotation, for our analysis. For the isochrone fitting, we introduce the added dimension of equatorial velocity $v_{eq}$ into our interpolation. The $v_{eq}$ distribution is calculated from the measured $v\sin I_{\star}$ value, scaled by a uniform distribution of orientations sampled in $\cos I_{\star}$.

To compare our Geneva isochrone results to fittings with more traditional 1D isochrones, we also present the results from analyses using the Dartmouth isochrones (Dotter et al., 2008).

The parameter space is explored with a Markov chain Monte Carlo (MCMC) analysis, using the affine-invariant ensemble sampler emcee (Foreman-Mackey et al., 2013). The observations are fitted for twice, with the per-point uncertainties for each dataset inflated such that the reduced $\chi^{2}$ is at unity for the second run. A $\cos i$ prior is imposed on the transit inclination, while a Gaussian prior is imposed on $T_{\mathrm{eff}}=6406\pm 64\,\mathrm{K}$, $\mathrm{[M/H]}=0.08\pm 0.05$, and $v\sin{I_{\star}}=30.9\pm 2.0\mathrm{km\,s}^{-1}$ based on the spectroscopic values outlined in Section 3.2. We note that the derived posterior $v\sin{I_{\star}}$ ($35.8\pm 1.1$ $\rm km\,s^{-1}$) is offset with the prior by $\sim 2\sigma$. Resetting the prior to $35.8\pm 1.1$ $\rm km\,s^{-1}$did not change the system parameters significantly, with a derived $\lambda$ $1\sigma$ upper limit of $<11^{\circ}$ (compared to $<14^{\circ}$ from our adopted results). The $K$-band magnitude is also constrained by a Gaussian prior about its 2MASS value (Skrutskie et al., 2006). The GAIA parallax is also heavily constrained by a Gaussian prior about our adopted value of $2.92\pm 0.23$ mas as described in Section 3.2. A $\beta$ distribution prior is imposed on the eccentricity, following the prescription for short period planets set out in Kipping (2013). Uniform priors are imposed on all other parameters.

Due to the large radius of HAT-P-67b, potential solutions in the MCMC chain lead to the planet overflowing its Roche lobe (e.g. Lecavelier des Etangs et al., 2004). We can use this to place a lower limit on the mass of the planet by assuming no Roche lobe overflow. For each link of the MCMC chain, we calculate the corresponding Roche lobe radius using equation A5 of Hartman et al. (2011). Links with $R_{p}/a$ overflowing the Roche lobe are eliminated. For the circular orbit fit, the Roche lobe provides a weak lower limit on the mass of the planet of $0.056\,M_{\mathrm{J}}$. The posterior distribution for planet mass is plotted in Figure 5. The final mass measurement we report is the 68% confidence interval for the Roche lobe constrained posterior distribution.

We present four sets of solutions in Tables 5 and 6 for the circular and eccentric orbit scenarios from the Geneva and Dartmouth isochrone fits. The circular orbit Geneva isochrone fit solution is preferred, favored over the eccentric solution with a Bayesian Information Criterion $\Delta\mathrm{BIC}$ of 212. That is, the increased degrees of freedom in an eccentric orbit fit do not justify the improvements in the goodness of fit over that of a circular orbit model.

The evolutionary stage of HAT-P-67 is shown in Figure 6 on the Hertzsprung-Russell diagram, along with evolutionary tracks of various stellar masses and rotation rates marked for context. The derived stellar and planetary parameters are presented in Tables 5 and 6, respectively.

3.5. Eccentricity constraint

We can constrain the eccentricity of the system via the photometric light curves despite a lack of detection of the radial velocity orbit, since the GAIA parallax provides a good constraint on the stellar radius and transit duration (Kipping, 2008; Dawson & Johnson, 2012). The eccentricity posterior, as constrained primarily from this ‘photo-eccentric’ effect, is shown in Figure 7. The eccentricity $2\sigma$ upper limit is 0.43, with a posterior median and 64% confidence region of $ecc=0.24\pm 0.12$.

The parallax we choose to adopt has an effect on our best-fit solutions. If we choose to adopt the GAIA parallax of $2.60\pm 0.23$ mas ($385\pm 34$pc) from Lindegren et al. (2016) without the systematic correction offered by Stassun & Torres (2016), we would have a modest $1.3\,\sigma$ tension between the best fit isochrone distance and the parallax distance. Adopting a distance of $385\pm 34$pc yields an eccentric orbit of $e=0.356_{-0.077}^{+0.072}$. The tidal circularization time scale for the system is $<500$ Myr (Dobbs-Dixon et al., 2004), so the likelihood of such an eccentric orbit is low for the system.

3.6. Transit timing variations and additional companions

To check for potential transit timing variations that may be indicative of additional orbiting companions, we re-fit the follow-up transit observations, allowing for individual transit centroids for each epoch. The timing residuals are shown in Figure 8. The transit geometry parameters $a/R_{\star}$, $R_{p}/R_{\star}$, and inclination, are heavily constrained by Gaussian priors about their best fit values from the global analysis (adopted as the circular orbit fit in Table 6). We find no convincing evidence for transit timing variations, but also note that the $\sim 7$ hr transit duration makes it difficult for us to capture full transits via ground-based follow-up, and partial transits provide poorer transit timing measurements. In addition, we find no evidence for long term radial velocity trend, with the quadratic and linear fits to the radial velocity data consistent with flat slopes.

3.7. Imaging Constraints on Resolved Neighbors

In order to detect possible neighboring stars which may be diluting the photometric transits, we obtained optical and near infrared imaging of HAT-P-67 using the Clio2 near-IR imager (Freed et al., 2004) on the MMT 6.5 m telescope on Mt. Hopkins, in AZ, together with the Differential Speckle Survey Instrument (DSSI; Howell et al., 2011; Horch et al., 2012, 2011) and the WIYN High-Resolution Infrared Camera (WHIRC), both on the WIYN 3.5 m telescope222The WIYN Observatory is a joint facility of the University of Wisconsin-Madison, Indiana University, the National Optical Astronomy Observatory and the University of Missouri. at Kitt Peak National Observatory in Arizona.

The Clio2 images were obtained on the night of UT 2011 June 22. Observations in $H$-band and $L^{\prime}$-band were made using the adaptive optics (AO) system. A possible neighbor was detected $4\farcs 9$ to the southeast of HAT-P-67 with a relative magnitude difference of $\Delta H=7.4\pm 0.5$ mag, but no closer objects are seen. The neighbor was blended with HAT-P-67 in the HATNet survey observations, but was fully resolved by all subsequent follow-up observations. Figure 9 shows the $H$-band magnitude contrast curve for HAT-P-67 based on these observations. This curve was calculated using the method and software described by Espinoza et al. (2016). The band shown in this image represents the variation in the contrast limit depending on the position angle of the putative neighbor. We can rule out other neighbors with a magnitude difference of $\Delta H<2$ mag, down to a separation of $0\farcs 3$, and $\Delta H<6$ mag, down to a separation of $0\farcs 8$. The $L^{\prime}$ observations suffered from high thermal background, and the $4\farcs 9$ neighbor was not detected. Meaningful constraints could not be placed on closer neighbors in $L^{\prime}$ based on these observations.

$J$-band snapshot images of HAT-P-67 were obtained with WHIRC on the night of 2016 April 24, with a seeing of $\sim 0\farcs 9$. The images were collected at four nod positions, and were calibrated, background-subtracted, registered and median-combined using the same tools that we used for reducing the KeplerCam images. The $4\farcs 9$ neighbor was not detected, and we concluded that it must have $\Delta J>7$ mag. The closest neighbor detected in these observations was at a separation of $9\farcs 3$ to the northwest, and has a relative magnitude difference of $\Delta J=4.96\pm 0.01$ mag compared to HAT-P-67. Figure 9 shows the $J$-band magnitude contrast curve computed in a similar manner to the $H$-band contrast curve.

The DSSI observations were gathered between the nights of UT 26 September 2015 and UT 3 October 2015. A dichroic beamsplitter was used to obtain simultaneous imaging through 692 nm and 880 nm filters. Each observation consisted of a sequence of 1000 40 ms exposures read-out on $128\times 128$ pixel ($2\farcs 8\times 2\farcs 8$) subframes, which were reduced to reconstructed images following Horch et al. (2011). These images were searched for companions, with none detected. Based on this, the $5\sigma$ lower limits on the differential magnitude between a putative companion and the primary star were determined as a function of angular separation as described in Horch et al. (2011). Based on these observations we exclude neighbors with $\Delta m<2.56$ at 692 nm, or $\Delta m<2.80$ at 880 nm, down to a limiting separation of $0\farcs 2$ (Figure 10).

3.8. Blend analysis

Blend scenarios are eliminated by the detection of the planetary Doppler tomographic transit signal. In the cases where an eclipsing binary blended with a foreground star is the cause of the transit signal, the Doppler tomographic shadow will be significantly diluted with respect to the photometric transit signal.

The flux under the shadow of the planet, as a fraction of the total flux under the rotational broadening kernel, describes the area blocked by the planet. This directly corresponds to a ‘transit light curve’ over the broadband of the TRES spectrum (following Zhou et al., 2016a). We plot this Doppler tomographic light curve in Figure 4 (bottom). We also plot the model transit light curve as per the global best fit solution. The spectroscopic transit depth is consistent with that of the photometric transit depth, confirming the lack of any significant dilution by background stars. The elimination of blend scenarios and the mass upper limit determined from HIRES radial velocities validates HAT-P-67b as a planet. We can also place strict upper limits on any third light contamination from background stars by modelling the line broadening profiles from the LSD analysis. A high signal-to-noise broadening profile was derived by averaging the 32 TRES spectra obtained for HAT-P-67. By modeling this profile as two stars, we place an upper limit on the flux ratio of any potential companion to be $<0.004$, or within 6 magnitudes of the primary star, with the caveat that any potential blended companion exhibits no radial velocity variation.

4. Discussion

We presented the discovery of HAT-P-67b, a hot-Saturn transiting an F-subgiant. HAT-P-67b has a radius of $2.085_{-0.071}^{+0.096}\,R_{\rm J}$, and a mass constrained by radial velocity measurements to be $M_{p}<0.59\,M_{J}$ at $1\,\sigma$. Confirmation of the planetary nature of HAT-P-67b involved numerous high precision follow-up transit light curves, radial-velocity constraints on its mass, and two Doppler tomographic transits that eliminated potential blended eclipsing binary scenarios.

The mass, radius, and densities of HAT-P-67b are plotted in Figure 11, along with selected parts of the gas giant population. HAT-P-67b is one of the largest, and one of the lowest density planets known $(\rho_{p}=0.052_{-0.028}^{+0.039}\,\rm g\,cm^{-3})$. A number of other inflated gas giants have been discovered around subgiants (KOI-680b Almenara et al. 2015, EPIC 206247743b Van Eylen et al. 2016, KELT-8b Fulton et al. 2015, KELT-11b Pepper et al. 2016, HAT-P-65b and HAT-P-66b Hartman et al. 2016) and giants (e.g. EPIC 211351816b Grunblatt et al. 2016). One hypothesis is that these gas giants are re-inflated by the evolved host star (Lopez & Fortney, 2016). In this scenario, as the host star evolves off the main sequence, ‘warm Jupiters’ are subjected to higher incident flux and stronger tidal heating (assuming a non-zero initial eccentricity). The heating reaches deep enough into the planetary interior to inflate the planet radius. Hartman et al. (2016)) found empirical evidence that the level of planet inflation is correlated with the fractional age of the host star, further supporting the idea of re-inflation. Figure 12 shows the evolution in the incident flux received by HAT-P-67b over its lifetime. Currently HAT-P-67b receives $\sim 2\times$ the incident-flux of a zero-age-main-sequence HAT-P-67, potentially inducing an inflation of the planetary radius.

Figure 12 also plots the incident flux received by the hot-Jupiter distribution against their planet masses. There is a paucity of low mass planets that receive high incident flux – a sharp envelope that likely resulted from the evaporation of Saturn and Neptune mass planets in close-in orbits (e.g. Lecavelier Des Etangs, 2007; Ehrenreich & Désert, 2011; Owen & Wu, 2013). HAT-P-67b lies on the edge of the envelope – unlike planets of similar masses that receive high incident irradiation, HAT-P-67b did not ‘boil-off’, but survived to the present day. The high incident flux may also have halted contraction early on, leading to its current radius. Since there is a lack of inflated Saturn-mass planets in high incident flux environments, HAT-P-67b is an important point in the mass-radius-flux relationship.

The low density and high irradiation of HAT-P-67b also results in a bloated atmosphere, with a large scale height of $\sim 500$ km (assuming an $H_{2}$ atmosphere), making the planet a good candidate for transmission spectroscopy follow-up studies.

In addition, X-ray and EUV-driven hydrodynamic escape play an especially important role in low density, low mass planets (e.g. Lecavelier Des Etangs, 2007; Murray-Clay et al., 2009; Ehrenreich & Désert, 2011; Owen & Jackson, 2012; Owen & Wu, 2013). For hot-Jupiters, X-ray and EUV photoionizes the upper atmosphere, causing it to heat up and expand, resulting in escaping flows. Atmospheric escape has been observed for HD 209458 b (Vidal-Madjar et al., 2003, 2004) and HD 189733 b (Lecavelier Des Etangs et al., 2010), where the Lyman-$\alpha$ radii of the planets are $\sim 10$ times larger than their optical radii, extending beyond the Roche sphere. Since the mass loss rate is largely proportional to the incident UV and X-ray flux received by the planets (e.g. Murray-Clay et al., 2009), we checked for existing X-ray and UV measurements of HAT-P-67. While no EUV or X-ray flux measurements exist, HAT-P-67 is identified as a source by GALEX, with flux measurements in the FUV (1344–1786Å) and NUV (1771–2831Å) bands. In Figure 13, we compile all transiting planet systems with GALEX FUV and NUV measurements (Bianchi et al., 2011), as well as GAIA parallaxes and updated system parameters from Stassun et al. (2016). To examine the potential mass loss rate of HAT-P-67b in the context of existing systems, we plot the UV fluxes received by each planet (normalized to that received by HD 209458 b) against their escape velocities. HAT-P-67b receives 24 times the FUV, and 10 times the NUV flux of HD 209458 b, and has one of the lowest escape velocities of known transiting planets ($25\,\mathrm{km\,s}^{-1}$, compared to $43\,\mathrm{km\,s}^{-1}$ for HD 209458 b). As such, it should be an excellent target for Lyman-$\alpha$ transit observations to measure its extended hydrogen exosphere.

HAT-P-67 has the highest $v\sin{I_{\star}}$ amongst all the evolved planet hosts ($35.8\pm 1.1\,\rm km\,s^{-1}$), allowing us to spectroscopically measure its projected spin alignment angle. The only previous Holt-Rossiter-McLaughlin measurement of a planet around a subgiant was WASP-71b, (Smith et al., 2013), which was found to be in a well-aligned orbit. Two additional planet-hosting evolved stars have had their line-of-sight stellar inclination measured via asteroseismology. Quinn et al. (2015) found the super-Jupiter Kepler-432b to be spin-orbit aligned in a eccentric, 53-day period orbit around a red-giant. Huber et al. (2013) found two co-planar planets residing in inclined orbits around their red-giant host Kepler-56.

As a host star expands over its post main-sequence evolution, star-planet tidal interaction should increase in strength, modifying the orbit of the planet. radial-velocity searches have found that eccentric planets are rarer around evolved hosts than around dwarfs (e.g. Jones et al., 2014). However, tidal interactions are likely weak during the lifetime of HAT-P-67b: the characteristic timescale for orbital decay is on the order of $\sim 10^{11}$ years (adopting eq.2 of Hansen, 2012, assuming an effective stellar dissipation coefficient of $\sigma_{\star}=10^{-8}$). Similarly, the tidal synchronization timescale – the time taken to synchronize the stellar spin vector with the orbit normal vector, leading to spin-orbit alignment and spin-orbit synchronization, is $\sim 10^{13}$ years (eq. 3 of Hansen, 2012), so the stellar spin is unlikely to have been modified over the lifetime of the system due to planet-star tidal interactions. The same cannot be said of WASP-71b, which orbits at a shorter period of 2.9 days, and a characteristic synchronization timescale of $\sim 10^{10}$ years, short enough that tides likely played a role in modifying the spin of the host star.

HAT-P-67b is part of a growing list of planets confirmed via Doppler tomography, including WASP-33b (Collier Cameron et al., 2010b; Johnson et al., 2015), Kepler-13b (Johnson et al., 2014), KELT-7b (Bieryla et al., 2015; Zhou et al., 2016a) HAT-P-57b (Hartman et al., 2015), KOI-12b (Bourrier et al., 2015), KELT-17b (Zhou et al., 2016b), and XO-6b (Crouzet et al., 2016). The routine success of the Doppler tomography technique yields exciting prospects to fill the paucity of transiting planets around early type stars.

We thank the referee for their insightful contributions to this manuscript. HATNet operations have been funded by NASA grants NNG04GN74G and NNX13AJ15G. Follow-up of HATNet targets has been partially supported through NSF grant AST-1108686. G.Á.B., Z.C., and K.P. acknowledge partial support from NASA grant NNX09AB29G. J.H. acknowledges support from NASA grant NNX14AE87G. K.P. acknowledges support from NASA grant NNX13AQ62G. We acknowledge partial support also from the Kepler Mission under NASA Cooperative Agreement NCC2-1390 (D.W.L., PI). Data presented in this paper are based on observations obtained at the HAT station at the Submillimeter Array of SAO, and the HAT station at the Fred Lawrence Whipple Observatory of SAO. This research has made use of the NASA Exoplanet Archive, which is operated by the California Institute of Technology, under contract with the National Aeronautics and Space Administration under the Exoplanet Exploration Program. Data presented herein were obtained at the WIYN Observatory from telescope time allocated to NN-EXPLORE through the scientific partnership of the National Aeronautics and Space Administration, the National Science Foundation, and the National Optical Astronomy Observatory. This work was supported by a NASA WIYN PI Data Award, administered by the NASA Exoplanet Science Institute. We acknowledge of J.A. Johnson in supporting the Keck HIRES observations. The authors wish to recognize and acknowledge the very significant cultural role and reverence that the summit of Mauna Kea has always had within the indigenous Hawaiian community. We are most fortunate to have the opportunity to conduct observations from this mountain.