Accepted for publication in The Astrophysical Journal Thermal Emission and Tidal Heating of the Heavy and Eccentric Planet XO-3b Pavel Machalek1,2, Tom Greene1, Peter R. McCullough3, Adam Burrows4, 0 Christopher J. Burke5, Joseph L. Hora5, Christopher M. Johns-Krull6, Drake L. Deming7 1 0 [email protected] 2 n a ABSTRACT J 3 1 We determined the flux ratios of the heavy and eccentric planet XO-3b to its ] parent star in the four IRAC bands of the Spitzer Space Telescope: 0.101% ± P E 0.004% at 3.6 micron; 0.143% ± 0.006% at 4.5 micron; 0.134% ± 0.049% at . 5.8 micron and 0.150% ± 0.036% at 8.0 micron. The flux ratios are within [-2.2, h p 0.3, -0.8, -1.7]-σ of the model of XO-3b with a thermally inverted stratosphere - o in the 3.6 micron, 4.5 micron, 5.8 micron and 8.0 micron channels, respectively. r t XO-3b has a high illumination from its parent star (F ∼(1.9 - 4.2) × 109 ergs s p a cm−2 s−1) and is thus expected to have a thermal inversion, which we indeed [ observe. When combined with existing data for other planets, the correlation 1 v between the presence of an atmospheric temperature inversion and the substellar 9 flux is insufficient to explain why some high insolation planets like TrES-3 do not 1 3 have stratospheric inversions and some low insolation planets like XO-1b do have 2 . inversions. Secondary factors such as sulfur chemistry, atmospheric metallicity, 1 0 amounts of macroscopic mixing in the stratosphere or even dynamical weather 0 effects likely play a role. Using the secondary eclipse timing centroids we deter- 1 : mined the orbital eccentricity of XO-3b as e = 0.277±0.009. The model radius- v i age trajectories for XO-3b imply that at least some amount of tidal-heating is X r a 1NASA Ames Research Center, MS 245-6, Moffett Field, CA 94035 2Bay Area Environmental Research Institute, 560 Third St West, Sonoma, CA 95476 3Space Telescope Science Institute, 3700 San Martin Dr., Baltimore MD 21218 4Department of Astrophysical Sciences, Princeton University, Princeton, NJ 08544 5Harvard-Smithsonian Center for Astrophysics, 60 Garden St., Cambridge, MA 02138 6DeptartmentofPhysicsandAstronomy,RiceUniversity,6100MainStreet,MS-108,Houston,TX77005 7Planetary Systems Laboratory, NASA/GSFC, Code 693.0, Greenbelt, MD 20771 – 2 – required to inflate the radius of XO-3b, and the tidal heating parameter of the planet is constrained to Q (cid:46) 106. p Subject headings: stars:individual(XO-3) — binaries:eclipsing — infrared:stars — planetary systems 1. Introduction The study of hot Jupiter atmospheres is maturing. In particular, low resolution spectra and broadband spectral energy distributions have been assembled from high precision pho- tometry of Hot-Jupiter’s day and night sides using the Spitzer Space Telescope’s InfraRed Array Camera (IRAC) (Knutson et al. 2008, 2009b; Tinetti et al. 2007; Charbonneau et al. 2008; Machalek et al. 2008, 2009; O’Donovan et al. 2009; D´esert et al. 2009; Todorov et al. 2010; Fressin et al. 2009; Christiansen et al. 2009), Infrared Spectrograph (IRS) (Grillmair et al. 2008) and Multi Band Imaging Spectrometer (MIPS) (Knutson et al. 2009a) as well the Hubble Space Telescope (Swain et al. 2008a, 2009a,b). Upper atmospheres of hot Jupiters are currently thought to be split into two classes depending on the stellar insolation at their substellar points: planets with substellar flux higher than F (cid:38) 109 erg cm−2 s−1 should posses temperature inversions in their stratosphere p as the intense stellar radiation is absorbed by upper atmospheric gaseous absorbing species (Hubeny et al. 2003; Burrows et al. 2008; Fortney et al. 2006, 2008; Spiegel et al. 2009). Planets with insolation fluxes F ∼ 0.5-1.0 × 109 ergs cm−2 s−1 like XO-2b, HAT-P-1, p OGLE-TR-113, and WASP-2 are in a transition zone between atmospheres with or without a stratosphere. Secondary effects like sulfur chemistry and atmospheric metallicity (Zahnle et al. 2009), amounts of macroscopic mixing in the stratosphere (Spiegel et al. 2009) or even dynamical weather effects (Showman et al. 2009; Rauscher & Menou 2009) could determine the stratospheric temperature profiles of these transition planets. XO-3b is a hot Jupiter with a high mass M = 11.79 ± 0.59 M (Winn et al. 2008; p Jup Johns-Krull et al. 2008), which is close to the deuterium burning limit and has one of the highest observed surface gravities, g = 209 m s−2 amongst the known transiting planets. Its 3.1915239 day long orbit around the parent star XO-3 (spectral type F5V, d = 260 ± 23 pc, Johns-Krull et al. (2008) ) has significant eccentricity e = 0.287± 0.005 (H´ebrard et al. 2008)), which causes stellar irradiance to vary three-fold over the entire orbit and causes the secondary eclipse to shift in time from half-phase. Furthermore Liu et al. (2008) estimated the amount of tidal energy dissipation rate – 3 – 3000 2500 4.5-µm 5.8-µm 8.0-µm 3.6-µm e (K) 2000 ur at er p em 1500 T 1000 P=0.3, κ=0.0 cm2/g Pn=0.2, κe=0.2 cm2/g n e 500 10-6 10-4 10-2 100 102 104 106 Pressure (bar) Fig. 1.— (top) Temperature / Pressure profiles for the atmosphere of XO-3b following the methodology of Burrows et al. (2007b, 2008); Spiegel et al. (2009) for heat redistribution parameter P =0.3 with no upper atmospheric optical absorber (dot-dashed line) and a cor- n responding model with a uniform upper atmospheric absorber (solid line) (absorption coeffi- cient κ =0.2 cm2/g and heat re-distribution parameter P =0.2) with depths corresponding e n to emission in the IRAC channels denoted. Both Temperature/Pressure profiles are calcu- lated for XO-3b orbital distance a = 0.0454 AU (Winn et al. 2008) and stellar insolation F ∼2.01× 109 ergs cm−2 s−1. See §3.2 for more details. p contributing to the inflated radius of XO-3b (R = 1.217 ± 0.073 R ; Winn et al. 2008) p Jup assuming the age of XO-3b t= 2.82 +0.58 Gyr (Winn et al. 2008). Liu et al. (2008) concluded −0.82 that the radius-age relationship for XO-3b is consistent to within 1.0-σ with no internal heat source (i.e. no tidal heating) or tidal heating dissipation with dimensionless tidal heating parameter Q (cid:38) 106 as defined by Goldreich & Soter (1966). By determining the exact p timing of the secondary eclipse in our 4 infrared light curves obtained with Spitzer Space Telescope IRAC we will refine the orbital eccentricity of XO-3b and constrain the amount of tidal heating (if any) responsible for inflating the planetary radius. In addittion to its high mass and significant orbital eccentricity, XO-3b was also the first planet with detected and confirmed non-zero sky projection angle λ = 37.3 ± 3.7 deg between the orbital axis and stellar rotation axis obtained from the Rositter-McLaughlin effect (H´ebrard et al. 2008; Winn et al. 2009), currently thought to be a result of planet-planet scattering (Nagasawa et al. 2008; Juri´c & Tremaine 2008). The substellar point flux at XO-3b is F ∼(1.9 - 4.2) × 109 ergs cm−2 s−1. The exact p value depends on the adopted stellar and planetary mass and radius, which are still un- certain Liu et al. (2008), as well as the changing distance distance from the star due to an – 4 – eccentric orbit. However this range of substellar point flux is clearly consistent with a promi- nent thermal inversion in the stratosphere. Figure 1 shows Temperature- Pressure models of Burrows et al. (2007b, 2008); Spiegel et al. (2009) and the predicted thermal inversion inthestratosphereandthenegativetemperaturegradientintheupperatmosphereofXO-3b. By obtaining the light curve of XO-3b in the 4 IRAC (4-8 microns) channels on the Spitzer Space Telescope and determining the depth and timing of the secondary eclipse in multiple wavelengths, we will be able to constrain the upper atmospheric temperature structure of XO-3b, refine the orbital eccentricity of the planet from the secondary eclipse timing centroids and hence its tidal heating rate, which could be responsible for inflating the radius of XO-3b. The Cold Spitzer IRAC observations in this work will provide a firm observational and theoretical foothold on the properties of the XO-3b atmosphere during the secondary eclipse and serve as comparison for future full orbit observation of XO-3b with Warm Spitzer similar to previous extended duration phase curves of hot Jupiters (Knutson et al. 2007, 2009a,c; Laughlin et al. 2009). Since there is a strong water band near IRAC 5.8 micron, coverage in all four IRAC bands will test for transitions between water in emission and in absorption, which can not be observed with Warm Spitzer. Furthermore as Fig. 1 illustrates, thereisasteeptemperaturegradientbetweendepthscorrespondingtoemissionin IRAC 5.8/8.0 micron channels, which can be uniquely studied with Cold Spitzer or otherwise with JWST in the future. The 5.8 and 8.0 micron channel planet/star flux ratios will further be correlated with the 3.6/ 4.5 micron flux ratio to test the two signatures of stratospheres. 2. Observations & Data Analysis The InfraRed Array Camera (IRAC; Fazio et al. 2004) has a field of view of 5.2(cid:48) × 5.2(cid:48) in each of its four bands. Two adjacent fields are imaged in pairs (3.6 and 5.8 microns; 4.5 and 8.0 microns). The detector arrays each measure 256 × 256 pixels, with a pixel size of approximately 1.22(cid:48)(cid:48) × 1.22(cid:48)(cid:48). We closely repeat the data analysis of Machalek et al. (2008, 2009) with modifications and improvements mentioned in the text. WehaveobservedXO-3 systeminall4channelsintwoseparateAstronomicalObserving Requests (AORs) in two different sessions: the 3.6 and 5.8 micron channels for 6.9 hours (with 2.9 hour long secondary eclipse) on UT 2009 March 17 (AOR 31618560) and the 4.5 and 8.0 micron channels for 6.9 hours on UT 2009 April 21 (AOR 31618816) with a 30- minute preflash on a bright uniform part of NGC1569. We used the full array 2s+2s/12s frame time in the stellar mode in which the 3.6 micron and 4.5 micron bands are exposed – 5 – for two consecutive 2s exposures while the 5.8 micron and 8.0 micron bands are integrating for 12s to prevent detector saturation. The 4.5 micron and 8.0 micron time series has been preflashed with a bright uniform extended target to prevent the initial “ramp-up” effect (Charbonneau et al. 2005; Deming et al. 2005; Knutson et al. 2008; Machalek et al. 2008, 2009), consequently no data points were removed from the beginning of the time series. The 3.6 micron and 5.8 micron time series however were obtained with no pre-flashing and hence exhibit an initial charge build up which is consistently removed during our detector effect removal. 2.1. InSb Detectors We have repeated our methodology from Machalek et al. (2009) by performing aperture photometry on the 3.6 micron and 4.5 micron time series with radii between 2.5 and 6.0 pixels. In order to test whether our secondary eclipse depths and centroid timings depend on apertureradius,wehaverepeatedtheentiredatareductionforapertureradiibetween2.5and 6.0 pixels in 0.5 pixel increments and obtained consistent results for different apertures. We have improved our photometry pipeline by obtaining the stellar centroids from flux-weighted position of a 5 × 5 pixel square centered on the peak stellar pixel (method suggested by Sean Carrey, private communication). Since our starting point was the BCD images produced by the pipeline version 18.7, cosmic rays were already rejected. The heliocentric modified Julian date at Spitzer spacecraft position recorded in the header keyword “HMJD OBS” did not necessitate our previous calculations of spacecraft positions (Machalek et al. 2008, 2009). We have chosen the aperture radii based on the RMS of residuals after detector effects and the secondary eclipse were removed. We used an aperture of radius 3.0 pixels for the 3.6 micron time series of XO-3, which had an RMS 0.0034 for out of transit points after decorelation. This is essentially Poisson noise limited, being only 1.01 higher than the predicted noise based on source brightness, detector read noise and gain. Similarly the 4.5 micron time series of XO-3 was obtained from 3.0 pixel radius aperture photometry which had the lowest RMS of 0.0049 which is 1.08 times higher than the predicted noise. The appropriate aperture corrections were applied to the photometry as specified by the Spitzer Data Handbook. As is evident from Fig. 3, the 3.6 micron time series exhibits a prominent flux variation with magnitude of ∼ 0.8 %, which is a well studied instrumental effect (Charbonneau et al. 2005; Morales-Caldero´n et al. 2006; Machalek et al. 2008; Knutson et al. 2009b; Machalek et al. 2009; D´esert et al. 2009) due to sub-pixel sensitivity variations caused by spacecraft – 6 – position drift of 0.1 - 0.3 arcsec over a period of ∼ 3000 seconds, which makes the star move onthepixel. The4.5microntimeseries, however, hasnegligiblefluxvariations, probablydue to a chance positioning on a pixel phase with a flat response curve (pixel reference: 126.46; 128.78). This pixel could be useful in planning for extended duration observations with Warm Spitzer. D´esert et al. (2009) has noted a similar pixel with a flat response function at pixel coordinates [147.20; 198.25]. Ourremovalofthesystematiceffectsandeclipsecurvefittingcloselyfollowsthemethod- ology of Machalek et al. (2009). The subpixel intensity variations in the 3.6 micron and 4.5 micron time series are detrended as a linear function of subpixel positions of the stellar centroid x, y, x2, y2, a linear function of time t, plus a constant for each of the two InSb channels: I = 1.0+ b x+ b y + b x2 + b y2 + b t, (1) 3.6micron 1 2 3 4 5 I = 1.0+ b x+ b y, (2) 4.5micron 1 2 We tried adding higher order terms of x and y, a cross terms of x×y, and a linear term linearinttothe4.5micron timeseriesdecorelation. However, addingtermsdidnotdecrease theχ2 orchangethesecondaryeclipsedepthorcentroidtiminginthe4.5micron timeseries, so we chose only two degrees of freedom (Eq. 2) for the 4.5 micron time series decorelation. Furthermore as can be seen from Fig. 2, the binned residuals in the decorrelated and fitted 4.5 micron light curve of XO-3 scale as N−1/2, where N is the number of points per bin. Since the binning of the residuals scales as N−1/2 we can conclude that negligible systamtic errors remain in the decorelated light curve. We fit the secondary eclipse with the formalism of Mandel & Agol (2002) with no stellar limb darkening and adopt the stellar and planetary parameters of Winn et al. (2008): R = 1.38+0.08 R , M = 11.79+0.59 M , R = 1.22+0.07 R 1, i = 84.20+0.54 degrees, (cid:63) −0.08 (cid:12) p −0.59 Jup p −0.07 Jup −0.54 and a = 0.0454 ± 0.0008 AU with ephemeris: T (E) = 2,454,449.86816(HJD)+E(3.1915239 days). (3) c We fit the 5 baseline parameters of Eq. 1 and the 2 baseline fitting parameters Eq. 2 concurrently with the secondary eclipse depth ∆F and the phase of the eclipse centroid Φ 11 R = 71,492 km. Jup – 7 – 0.0060 0.0040 s al u d si e of r s m r 0.0020 0.0000 1 10 100 1000 Fig. 2.—Rootmeansquareofbinnedlightcurvepointsafterdetectorremovalandsecondary eclipse fitting for the 4.5 micron photometry of XO-3 vs. the number of points per bin. The sub-pixelphasesensitivityvariationswereremovedusingEq. 2. Thesolidlineisproportional to the number of residual points per bin N−1/2. for a total of 7 and 4 fitting parameters, respectively. This was done to properly account for the way in which systematic effects removal affects the secondary eclipse fitting. The best parameter solutions were obtained by using a Monte Carlo Markov Chain (MCMC) with 105 iterations(Gregory2005;Markwardt2009)withratioofjumpsbetween20-40%. Thebestfit parameter values were obtained by discarding the first 20% of the iterations to prevent initial conditions from influencing the results and adopting the median of the distribution of each parameter as the best fit value. These values are reported in Table 1 with errors obtained from symmetric 66.8% contours around the median of the posterior probability distribution of the MCMC runs. The decorelated best-fit light curves are depicted in Fig. 4 binned in 3.5 minute intervals. Note, however, that all our analysis is performed on the unbinned data. We find that the XO-3 3.6 micron time series shows a linear flux increase with a slope of b = 0.015% ± 0.002% per hour which is consistently removed from our photometry, 4 but inconsistent with the slope of XO-2 at 3.6 micron of b = -0.011% ± 0.005% per hour 4 (Machalek et al. 2009). This flux decrease has been attributed by Machalek et al. (2009) and Knutson et al. (2009b) to an instrumental effect on the In:Sb detectors. When we added a linear time term b t to the decorelation of the 4.5 micron time series of XO-3 in Eq. 2, 3 its value was consistent with zero. Thus we omitted a linear time term b t from the final 3 analysis. – 8 – 2.2. Si:As Detectors The 5.8 micron and 8.0 micron time series is recorded with Si:As detectors, which have a different set of systematic effects from the 3.6 micron and 4.5 micron InSb detectors. We have performed aperture photometry on the 5.8 micron and 8.0 micron images with aperture radii ranging from 3.0 to 6.0 pixels, choosing the aperture radius with the lowest RMS of the residuals after systematic effects and the secondary eclipse were removed. This resulted in an aperture of radius 3.5 pixels for the 5.8 micron time series with a detrended RMS of 0.0055 (42% higher than Poisson noise) and an aperture radius of 4.5 pixels for the 8.0 micron time series with a detrended RMS of 0.0049 (60% higher than Poisson noise). Nopointswereremovedfromthebeginningofeitherthe5.8micronor8.0micron timeseries. A well studied instrumental effect of the Si:As arrays is the gain variations of individual pixels over time, which result in flux decrease/increase in the light curve (e.g. Deming et al. 2005; Knutson et al. 2007, 2008; Machalek et al. 2008; D´esert et al. 2009), quite un- like the pixel position dependent flux effect in the InSb 3.6 micron and 4.5 micron arrays. Machalek et al. (2009) and Laughlin et al. (2009) have reported that the gain variations in the 5.8 micron and 8.0 micron channels and resultant flux trends in the light curves differ for the two components of a binary star, which have the same brightness and similar colors, suggesting that relative placement of the stellar centroid with respect to the edges of the pixels determines the Si:As detector pixel response. The gain variations can be clearly seen in the 5.8 micron and 8.0 micron light series in Fig. 3: a nonlinear decrease in brightness in the 5.8 micron light series and a nonlinear flux increase in the 8.0 micron time series. To remove the nonlinear flux variation inherent to the Si:As detector, we fit the sec- ondary eclipse depth ∆F along with the eclipse centroid phase Φ concurrently with the 3 “ramp” decorelation coefficients as follows: I = a +a × ln(∆t+0.05)+a × ln(∆t+0.05)2 (4) model 1 2 3 where I is the normalized model flux and ∆t is the time in days since the beginning model of the integration (constant of +0.05 inserted to avoid singularity at ∆t=0.). We fit the 5 parameters (2 for the eclipse and 3 for the “ramp” in Eq. 4) for the 5.8 micron and 8.0 micron time series concurrently using 105 MCMC runs with errors adopted as the 66.8 % contours around the median of the posterior distribution of the MCMC runs for each parameter. To ensure that our results are not dependent on the aperture radius we have repeated the MCMC runs for all aperture radii between 3.0 and 6.0 piels in 0.5 pixel in- – 9 – crements and found the timing centroids to be consistent. The secondary eclipse depths were, however, found to vary by about 1-σ for photometry with aperture radii between 3.0 and 6.0 pixels. Hence to be conservative, as stated above, we have adopted the secondary eclipse depths from the aperture photometry with the lowest RMS of residuals after eclipse removal. These aperture radii were 3.5 pixels for the 5.8 micron time series and 4.5 pixels for the 8.0 micron time series. We adopted uncertainties as the upper and lower envelope of the eclipse depths with their uncertainties for photometry with aperture radii between 3.0 and 6.0 pixels. Note, however, that these large, conservative uncertainties of the 5.8 micron and 8.0 micron eclipse depths ( ∆F = 0.134% ± 0.049% and ∆F = 0.150% ± 5.8 micron 8.0 micron 0.036%) still allow us to distinguish between the two models for the upper atmospheric tem- perature structure of XO-3b (see Fig. 5). The final results are reported in Table 1, and the decorelated time series was binned in 3.5-minute bins for viewing clarity in Fig. 4. 3. Discussion The IR light curves presented in this work allow for the determination of the exact timing of the secondary eclipse centroid such that the orbital eccentricity can be refined more accurately than from the radial velocity (RV) curve. The temperature structure of the upper atmosphere of XO-3b can also be determined from the light curves by comparing the secondary eclipse depths (i.e. the planet/star contrastratios) to atmospheric models. A refined eccentricity determines the rate of tidal heating of the planet and can help explain the inflated radius of XO-3b. 3.1. Tidal heating rate and the radius of the planet A subset of transiting extrasolar giant planets (EGPs) have radii larger than standard models can accommodate (Guillot et al. 1996; Bodenheimer et al. 2001, 2003; Chabrier et al. 2004; Ibgui & Burrows 2009). Numerous explanations have been suggested as sources of the inflated radii of EGPs (see Fortney & Nettelmann (2009) for a review). Working in opposition to any inflation mechanism, heavy-element inner cores lead to smaller planetary radii compared to pure H/He objects (Burrows et al. 2007a; Fortney et al. 2007; Baraffe et al. 2008). Tidal inflation has been a popular explanation, as the dissipation of orbital energy into the inner regions of a planet can lead to inflated radii (Bodenheimer et al. 2001, 2003; Liu et al. 2008). Radius-age trajectories for extra-solar giant planets (EGPs) are presented by Liu et al. (2008) to explain the inflated radius of several planets, including XO-3b, which is larger than theoretical predictions. – 10 – Liu et al. (2008) have investigated the radius of XO-3b as a function of planetary age t, planetary radius R , planetary mass M , orbital eccentricity e, planetary metallicity [Fe/H], p p and tidal heating parameter Q . They conclude that for the parameters adopted from the p photometric follow-up of XO-3b by Winn et al. (2008), which are used in our study (see §2.1), the radius R = 1.22+0.07 R adopted by Winn et al. (2008) is consistent within p −0.07 Jup 1-sigma to either no internal heat source or tidal energy dissipation with tidal heating pa- rameter Q (cid:38) 106.0. This is the heating parameter for the adopted eccentricity based upon p RV measurements from Winn et al. (2008), e = 0.260± 0.017. We have also adopted these parameters of Winn et al. (2008) up to this point in this study (see §2.1). We refine the eccentricity e of XO-3b using the weighted average of the secondary eclipse timing centroids from Table 1 using the displacement from half orbital phase as a measurement of eccentricity (e.g. Kopal (1959) Eq. 9.23): 2πΦ (cid:39) π +2e×cos(ω)(1+csc2(i))+... (5) where e is the eccentricity, Φ is the orbital phase of the time centroid of the secondary eclipse, ω is the longitude of periastron, and i is the planetary orbit inclination. Using argument of pericenter ω=345.8 deg ± 7.3 deg and inclination i = 84.20 deg ± 0.54 deg from Winn et al. (2008), we derive a refined eccentricity of the XO-3b system from our secondary eclipse timings: e = 0.277±0.009 (6) with uncertainties formally propagated through Eq. 5. Taken individually the March 2009 secondary eclipse phase centroids from Table 1 imply an eccentricity of 0.278±0.010 and the April 2009 secondary eclipse phase centroids imply eccentricity of 0.276 ± 0.009, which are consistent with each other. The intriguing possibility of eccentricity changing on a timescale of months will be further studied during the Warm Spitzer mission phase obser- vationsofXO-3binspring2010,whenboththetransitandsecondaryeclipsewillbeobserved. The refined value of XO-3b eccentricity e = 0.277± 0.009 is 1.0-σ higher than the Winn et al. (2008); Johns-Krull et al. (2008) eccentricity e = 0.260 ± 0.017. It is also 2.0-σ lower than the radial velocity derived eccentricity e = 0.287± 0.005, H´ebrard et al. (2008). The tidal heating of XO-3b, which is a strong function of eccentricity, can inflate the radius of