DETECTION OF PLANETARY EMISSION FROM THE EXOPLANET TrES-2 USING SPITZER/IRAC

There has been a recent dramatic increase in the number of extrasolar planets within 300 pc whose structures and atmospheric compositions can be probed using the Spitzer Space Telescope (Werner et al. 2004). These are the nearby transiting exoplanets. By measuring with Spitzer the decrease in light as one of these planets passes behind its star in an event known as a secondary eclipse, we can estimate the flux emitted by the planet. Detection of this emission from planetary atmospheres is made possible by taking advantage of the enhanced contrast between stars and their planets in the infrared wavelengths observable with Spitzer. Combining several of these flux measurements allows us to characterize the shape of the planet's emission spectrum, which tells us about the properties of its dayside atmosphere. (For a discussion of extrasolar planetary atmospheres, see Charbonneau et al. 2007 and Marley et al. 2007.)

The planets HD 209458b and HD 189733b have been the optimal choices for Spitzer studies of extrasolar planetary atmospheres because of their early discovery and the relative brightness of their stellar hosts. We have detected infrared emission from these exoplanets, both photometrically (Deming et al. 2005, 2006; Knutson et al. 2007, 2008) and with low-resolution spectroscopy (Richardson et al. 2007; Grillmair et al. 2007, 2008; Swain et al. 2008a, 2008b, 2009a, 2009b; Charbonneau et al. 2008). Another early success was the report by Charbonneau et al. (2005) of the measurement of infrared light from TrES-1. More recently, Harrington et al. (2007), Deming et al. (2007a), Demory et al. (2007), and Machalek et al. (2008) announced the results from Spitzer observations of the transiting exoplanets HD 149026b, GJ 436b, and XO-1b.

There has been a flurry of activity to reconcile atmospheric models with this limited number of infrared measurements. While several attempts have been made to explain the infrared observations (see, e.g., Barman et al. 2005; Burrows et al. 2005, 2007a; Fortney et al. 2005, 2006; Seager et al. 2005), the models are not entirely in agreement and no single model can explain every observation. Recently, Burrows et al. (2007b), Burrows et al. (2008), and Fortney et al. (2008) supplied a possible piece of the puzzle by proposing that the very highly irradiated hot Jupiters such as HD 209458b and TrES-2 (see Figure 1 of Fortney et al. 2008) will exhibit water emission rather than the expected water absorption at the IRAC wavelengths, a result of a temperature inversion in their atmospheres. The large stellar insolation that TrES-2 experiences may permit the presence of TiO and/or VO molecules in the hot planetary atmosphere that would condense in a cooler atmosphere. These opaque molecules would then cause the inversion. However, while this hypothesis can explain the emission features seen in the spectrum of HD 209458b (Burrows et al. 2007b; Knutson et al. 2008), it is incomplete. The relatively lower insolation experienced by the planet XO-1b is similar to that of TrES-1, and hence both atmospheres would be predicted to have water absorption features, as is indeed supported by the infrared observations of TrES-1 by Charbonneau et al. (2005). Nevertheless, Machalek et al. (2008) have shown that XO-1b displays the contrary, with evidence of water emission in its atmosphere. More recently, Fressin et al. (2009) did not find any evidence for the expected thermal inversion in the atmosphere of the highly irradiated planet TrES-3. Zahnle et al. (2009) have shown that sulfur photochemistry may also lead to the formation of inversions. In their models, this photochemistry is more or less temperature independent between 1200 K and 2000 K, and is also relatively insensitive to atmospheric metallicity.

Although over 400 extrasolar planets are known, it is only for these nearby transiting exoplanets that we can measure the planetary radii and true planetary masses precisely enough to provide useful constraints for theoretical models. There have been problems reconciling the observed planetary masses and radii with models (see Laughlin et al. 2005, Charbonneau et al. 2007, Liu et al. 2008, and Ibgui & Burrows 2009, for a review), namely that there are planets such as TrES-2 whose radii are larger than predicted by the standard models for an irradiated gas giant planet. Several explanations have been proposed for bloated planets like TrES-2, mainly related to some additional source of energy to combat planetary contraction. One such possible energy source is the tidal damping of a non-zero eccentricity (Bodenheimer et al. 2001, 2003): the planetary orbit of a hot Jupiter is expected to be circular, unless it is gravitationally affected by an unseen planetary companion (see, e.g., Rasio & Ford 1996). We can constrain the likelihood of a tidal damping energy source by measuring the timing of the secondary eclipse of a transiting hot Jupiter like TrES-2 and comparing these timings to predictions based on transit timings and the hypothesis of a circular planetary orbit (see, e.g., Deming et al. 2005).

The transiting hot Jupiter TrES-2 (O'Donovan et al. 2006) is one of the known hot Jupiters with a radius larger than expected from current models. The atmosphere of TrES-2 experiences similar levels of irradiation from its host star as the atmosphere of the likewise bloated planet HD 209458b, and hence we expect to find evidence of a thermal inversion for TrES-2, according to the predictions of Burrows et al. (2007b, 2008) and Fortney et al. (2008). Here, we present the first Spitzer observations of TrES-2 (Section 2). From our analysis (Section 3), we have detected thermal emission from the transiting planet. We found no evidence for timing offsets of the secondary eclipses and deduced the possible presence of a thermal inversion in the atmosphere of TrES-2 (Section 4).

We monitored TrES-2 using Spitzer during the time of two secondary eclipses, employing a different pair of the four bandpasses available on the Infrared Array Camera (IRAC; Fazio et al. 2004) during each of the eclipses. We took care to position TrES-2 (2MASS J19071403+4918590: J = 10.232 mag, J − K_s = 0.386 mag) away from array regions impaired by bad pixels or scattered light. We also kept the corresponding IRAC stray light avoidance zones free of stars that are bright in the infrared. We observed a 52 × 52 field of view (FOV) containing TrES-2 during two eclipses, using the stellar mode of the IRAC instrument. On UT 2006 November 30 (starting at HJD 2, 454, 069.956), we obtained 1073 images of this FOV at 4.5 μm and 8.0 μm with an effective integration time of 10.4 s for a total observing time of 3.9 hr. Our 3.6 μm and 5.8 μm observations of TrES-2 were taken on UT 2007 August 16. The observations began at HJD 2, 454, 324.436 and lasted 4.0 hr, during which we acquired 2130 and 1065 images in the respective channels (of effective exposure time 1.2 s and 10.4 s, respectively).

As part of the Spitzer Science Center pipeline for IRAC data (ver. S15.0.5 for the 2006 observations, ver. S16.1.0 for the 2007 observations), the images were corrected for dark current, flat-field variations, and certain detector nonlinearities. Each header of these Basic Calibrated Data (BCD) images contains the time and date of observation and the effective integration time. We used these to compute the Julian date corresponding to mid-exposure of each observation. In order to convert these dates to Heliocentric Julian Dates, we calculated the corresponding light travel time between the Sun and Spitzer (ranging from 0.5 to 2.5 minutes), using NASA JPL's HORIZONS⁹ service. We then computed the orbital phases using an updated orbital period and transit epoch (see Table 1) derived from observations (Holman et al. 2007) made as part of the Transit Light Curve project.

Using the nearest integer pixel (x₀, y₀) as an initial estimate of the position of TrES-2 on the array, we computed the flux-weighted centroid (x, y) of TrES-2 in each BCD image. The intra-pixel position of TrES-2 is then (x' = x − x₀, y' = y − y₀). We measured the flux from our target using circular apertures ranging from 2 to 10 pixels, and subtracted the background signal assessed in a sky annulus with inner and outer radii of 20 and 30 pixels, respectively. We normalized the fluxes for a given channel and aperture size by dividing each time series by the median of its values outside the times of eclipse. We examined the variation of the residual scatter of the out-of-eclipse data and found that an aperture radius of 3 pixels produced the smallest residual scatter for all four channels. Figure 1 shows the four light curves we obtained using this photometric aperture.

Zoom In Zoom Out Reset image size

For each of the observed secondary eclipses, we created a model of the two IRAC light curves obtained during that eclipse. We first accounted for various detector effects known to be present in IRAC data. There is a known correlation between the IRAC 3.6 μm or 4.5 μm flux from a source and the intra-pixel position (x', y') on the detector (Reach et al. 2005; Charbonneau et al. 2005, 2008; Knutson et al. 2008): the sensitivity of an individual pixel varies depending on the location of the stellar point-spread function, with higher fluxes measured near the center of the pixel and lower fluxes near the edges. Figure 1 shows this effect in our TrES-2 data. Data from these two channels also demonstrate a linear trend with time (as previously observed by Knutson et al. 2009 in IRAC observations of TrES-4), with a positive trend at 3.6 μm and a negative slope at 4.5 μm. We removed these trends to obtain f₀, the actual stellar flux:

$$\begin{equation} f_{0} = f / ([c_{1} + c_{2}x^{\prime } + c_{3}y^{\prime } + c_{4}(x^{\prime })^2 + c_{5}(y^{\prime })^2] \times [1 + Cdt]), \end{equation} \tag{ 1 }$$

where f is the measured stellar flux, (x', y') is the intra-pixel position, dt is the amount of time from the first observation, and (c_1-5, C) are free parameters in our model. The 5.8 μm and 8.0 μm data (see Figure 1) display the "ramp" associated with these IRAC channels first noticed by Charbonneau et al. (2005) and expanded upon by Harrington et al. (2007; supplementary information). Both data sets showed an overall increase in flux with dt. For these two data sets, we removed this detector effect by including the following correction in the model:

$$\begin{equation} f_{0} = f / [d_{1} + d_{2}(\ln {dt^{\prime }}) + d_{3}(\ln {dt^{\prime }})^2], \end{equation} \tag{ 2 }$$

where d₁–d₃ are free parameters and dt' = dt + 0.02 days. Here our substitution of dt' for dt is simply to avoid division by infinity. We then modeled simultaneously the two corrected time series using the eclipse light curve equation for a uniform source from Mandel & Agol (2002). We obtained the required system parameters (see Table 1) from Holman et al. (2007): the planetary orbital period, impact parameter, orbital inclination, and the radius ratio between the planet and the star. Based on this ephemeris, we calculated the predicted eclipse epoch (T_s; see Table 1). The three free parameters for the eclipse model were the timing offset (Δt), the eclipse depth at the shorter wavelength (Δf_l), and the depth at the longer wavelength (Δf_h). The model of the two light curves therefore has 12 free parameters: [c₁, c₂, c₃, c₄, c₅, C, d₁, d₂, d₃, Δt, Δf_l, and Δf_h].

For a particular instance of this model, we compute the χ² measure of the goodness of fit to the two relevant data sets as follows. For each light curve, we compute an initial uncertainty for the normalized flux values as their standard deviation. We exclude 5σ outliers in flux from further consideration. We then exclude large outliers in intra-pixel position (max [|x' − x'_m|, |y' − y'_m|]>0.15, where [x'_m, y'_m] are the median intra-pixel values) on the detector. We compute the χ² of the two light curves separately, and rescale each χ² to reflect the reduction of the number of data points by the above outlier exclusion. We then sum the resulting values to compute the χ² of the overall model.

To find an initial estimate of the best-fit parameters for each model of two light curves, we first used the AMOEBA algorithm (Press et al. 1992) to minimize the χ² of the fit. Using this initial estimate as a starting point, we applied the Markov Chain Monte Carlo method (see, e.g., Ford 2005; Winn et al. 2007), computing the χ² at each of the 10⁶ steps of the chain. We then calculated the median of the 10⁶ χ² values, and excluded from further analysis all the steps prior to the occurrence of the first value lower than this median. For the ith free parameter of our model, we derived the best-fit value p_i as the median of the remaining distribution of values for that parameter. We computed the (possibly unequal) lower (σ_−,i) and upper (σ_+,i) errors in this value such that the ranges [p_i − σ_−,i, p_i] and [p_i, p_i + σ_+,i] each contains 68%/2 of the values less than or greater than, respectively, the best-fit value. We computed the flux residuals after dividing out the best-fit model, and then derived an updated uncertainty as the new standard deviation of the residuals from the model. We again computed and rescaled the χ² of the two light curves separately, and then summed the resulting values to compute the χ² of the best-fit.

We tabulate for both pairs of channels the best-fit values for the free parameters of the eclipse models and the reduced χ² in Table 2. We have overplotted in Figure 1 the above corrective functions using the best-fit parameters we derived from the Markov chains. In Figure 2, we plot the corrected fluxes from TrES-2 at the four wavelengths and overplot the best-fit eclipse models. In both figures, the error bar shown for each binned data point is the standard deviation of the flux values in that bin, divided by the square root of the number of points in the bin.

Zoom In Zoom Out Reset image size

From the best-fit values (see Table 2) for the timing offset of the two observed secondary eclipses, we see that their weighted average (Δt = 1.2 ± 2.3 minutes) is consistent with no offset from the predicted epochs for the eclipses. Note that we have not accounted for the light travel delay time (see Loeb 2005) of 37 s across the TrES-2 system because of the relatively large size of the errors in these timing offsets compared to this delay time. An upper limit for the orbital eccentricity of a transiting planet can be computed from the timing offset Δt, using e cos ω ≃ π Δt/2 P, where ω is the unknown longitude of periastron and P is the known orbital period (Charbonneau et al. 2005, Equation (4)). The 3σ upper limit for TrES-2 is therefore 0.0036, consistent with a negligible orbital eccentricity, unless w ≃ 90°. Tidal damping of orbital eccentricity (Bodenheimer et al. 2001, 2003) is therefore unlikely to be a sizable contribution to the internal energy of this bloated exoplanet.

We now turn to a discussion of the TrES-2 planet–star contrasts. We first emphasize how well the data can be fit by a blackbody spectrum with no molecular band features. The black filled circles with error bars in Figure 3 show the Spitzer IRAC data points from this work. The blue dashed line is a 1500 K blackbody flux divided by the Kurucz stellar model with stellar parameters (T_eff = 5750 K, Z = 0.0, log g = 4.5) closest to those (T_eff = 5850 ± 50 K, Z = −0.15 ± 0.10, log g = 4.4 ± 0.1) derived by Torres et al. (2008). It is clear that the blackbody fits all the data points reasonably well, given the large error bars. However, a blackbody spectrum is only a nominal guideline, since the actual planetary spectrum is influenced by the myriad contributions due to molecular band features, collision induced opacities, and temperature gradients.

Zoom In Zoom Out Reset image size

Model atmospheric spectra for TrES-2 are also shown in Figure 3. The red and green spectra show models with and without a thermal inversion on the planet dayside, respectively. The red and green circles (enclosed in black circles) show the corresponding model points obtained by integrating the spectra over the Spitzer IRAC bandpasses. The corresponding model thermal profiles are shown in Figure 4. The spectra were generated using the hot Jupiter atmosphere model developed in Madhusudhan & Seager (2009). We consider a cloudless atmosphere, and the molecular species are assumed to be well mixed. The stellar spectrum was represented by the appropriate Kurucz model. The model spectrum without a thermal inversion has T_eff = 1634 K, and the model spectrum with a thermal inversion has T_eff = 1459 K. Both the models allow for extremely efficient day-night energy redistribution. At face value, we find that both the models fit the data almost equally well.

Zoom In Zoom Out Reset image size

However, the two models require different molecular compositions, which helps us determine the more probable model. The model without a thermal inversion has uniform molecular mixing ratios of 10⁻⁴ for H₂O, 10⁻⁶ for CO, and 10⁻⁶ for CH₄. While the mixing ratio of H₂O is plausible, the low mixing ratio of CO is surprising. In a hot atmosphere, with T_eff = 1500 K, CO is expected to be highly abundant. On the other hand, the model with a thermal inversion has uniform molecular mixing ratios of 10⁻⁴ for H₂O, 10⁻⁴ for CO, 5 × 10⁻⁵ for CH₄, and 2 × 10⁻⁶ for CO₂. These compositions show a relatively high abundance of CO, as expected in a hot Jupiter atmosphere like that of TrES-2. Based on these compositions, the presence of a thermal inversion seems like a more favorable scenario. Additionally, the requirement of 2 × 10⁻⁶ of CO₂ for the inversion model suggests an enhanced metallicity of [M/H] ∼ 0.7 (Zahnle et al. 2009). However, a thorough exploration of the parameter space is needed to place constraints on the thermal inversion in conjunction with the molecular compositions (Madhusudhan & Seager 2010).

Thus, we see some evidence that TrES-2 is part of the group of planets that display excess emission at wavelengths longer than 4 μm. The presence of a thermal inversion in the atmosphere of TrES-2 would support the theory of Fortney et al. (2008) and Burrows et al. (2008) that planets with substellar fluxes greater than or equal to approximately 10⁹ erg s^-1 cm⁻² should show evidence for a temperature inversion in their atmospheres.

Our theoretical models for the planetary flux from TrES-2 are currently constrained only at infrared wavelengths. Such an infrared signal appears to be relatively insensitive to cloud cover (Burrows et al. 2005). However, the amount of starlight reflected off the planetary atmosphere (at optical wavelengths) is highly dependent on the presence and size of upper-atmospheric condensates such as MgSiO₃, Fe, and Al₂O₃ (Marley et al. 1999; Green et al. 2003). Hence, a tighter constraint of atmospheric models for TrES-2 would be derived from the combination of measurements of the optical planetary flux from observations of secondary eclipses with the results from our infrared observations. Two likely sources for such observations of the optical flux from TrES-2 are the Kepler mission (Borucki et al. 1997) and the EPOXI mission (Deming et al. 2007b). Kepler has already demonstrated the ability to monitor the variation in thermal emission and reflected light from HAT-P-7 (Borucki et al. 2009). We eagerly await the comparison of our infrared observations with future optical observations from these two missions.

Although our observations of the atmospheric emission from TrES-2 at 3.6 μm, 4.5 μm, 5.8 μm, and 8.0 μm provide limited spectral coverage, we have been able to deduce the probable presence of gaseous molecules with high opacities in the atmosphere that result in emission in the 4.5 μm band. The highly irradiated gas giant TrES-2 thus may provide additional evidence for the correlation between the occurrence of thermal inversion in the atmosphere of a planet and the level of stellar insolation experienced by the planet.

Our observations of two secondary eclipses by this exoplanet occurred at the time predicted using the time of transit and the assumption of a circular orbit, within the errors. From this we conclude that the source of additional energy required to inflate the planetary radius to its bloated size is unlikely to be tidal heating caused by the circularization of an eccentric orbit.

This work is based on observations made with the Spitzer Space Telescope, which is operated by the Jet Propulsion Laboratory, California Institute of Technology under a contract with the National Aeronautics and Space Administration (NASA). This research was supported in part by NASA under grant NNG05GJ29G (issued through the Origins of Solar Systems Program) and also by an appointment to the NASA Postdoctoral Program at the Goddard Space Flight Center (administered by Oak Ridge Associated Universities through a contract with NASA).

Facilities: Spitzer (IRAC) - Spitzer Space Telescope satellite