Towards consistent mapping of distant worlds: secondaryeclipse scanning of the exoplanet HD 189733b^{⋆}
^{1} Department of Earth, Atmospheric and Planetary Sciences, MIT, 77 Massachusetts Avenue, Cambridge, MA 02139, USA
email: jdewit@mit.edu
^{2} Faculté des Sciences Appliquées, Université de Liège, Grande Traverse 12, 4000 Liège, Belgium
^{3} Institut d’Astrophysique et de Géophysique, Université de Liège, Allée du 6 Août 17, 4000 Liège, Belgium
^{4} Department of Physics and Kavli Institute for Astrophysics and Space Research, MIT, 77 Massachusetts Avenue, Cambridge, MA 02138, USA
Received: 16 February 2012
Accepted: 31 October 2012
Context. Mapping distant worlds is the next frontier for exoplanet infrared (IR) photometry studies. Ultimately, constraining spatial and temporal properties of an exoplanet atmosphere (e.g., its temperature) will provide further insight into its physics. For tidallylocked hot Jupiters that transit and are eclipsed by their host star, the first steps are now possible.
Aims. Our aim is to constrain an exoplanet’s (1) shape, (2) brightness distribution (BD) and (3) system parameters from its phase curve and eclipse measurements. In particular, we rely on the secondaryeclipse scanning which is obtained while an exoplanet is gradually masked by its host star.
Methods. We use archived Spitzer/IRAC 8μm data of HD 189733 (six transits, eight secondary eclipses, and a phase curve) in a global Markov chain Monte Carlo (MCMC) procedure for mitigating systematics. We also include HD 189733’s outoftransit radial velocity (RV) measurements to assess their incidence on the inferences obtained solely from the photometry.
Results. We find a 6σ deviation from the expected occultation of a uniformlybright disk. This deviation emerges mainly from a largescale hot spot in HD 189733b’s atmosphere, not from HD 189733b’s shape. We indicate that the correlation of the exoplanet orbital eccentricity, e, and BD (“uniform time offset”) does also depend on the stellar density, ρ_{⋆}, and the exoplanet impact parameter, b (“ebρ_{⋆}BD correlation”). For HD 189733b, we find that relaxing the eccentricity constraint and using more complex BDs lead to lower stellar/planetary densities and a more localized and latitudinallyshifted hot spot. We, therefore, show that the light curve of an exoplanet does not constrain uniquely its brightness peak localization. Finally, we obtain an improved constraint on the upper limit of HD 189733b’s orbital eccentricity, e ≤ 0.011 (95% confidence), when including HD 189733’s RV measurements.
Conclusions. Reanalysis of archived HD 189733’s data constrains HD 189733b’s shape and BD at 8 μm. Our study provides new insights into the analysis of exoplanet light curves and a proper framework for future eclipsescanning observations. In particular, observations of the same exoplanet at different wavelengths could improve the constraints on HD 189733’s system parameters while ultimately yielding a largescale timedependent 3D map of HD 189733b’s atmosphere. Finally, we discuss the perspective of extending our method to observations in the visible (e.g., Kepler data), in particular to better understand exoplanet albedos.
Key words: Eclipses / planets and satellites: individual: HD 189733b / techniques: photometric
Movies are available in electronic form at http://www.aanda.org
© ESO, 2012
1. Introduction
Among the hundreds of exoplanets detected to date^{1}, the transiting ones can be most extensively characterized with current technology. Beyond derivation of the exoplanet orbital inclination and density (e.g., Winn 2010), transiting exoplanets are key objects because their atmospheres are observationally accessible through transittransmission and occultationemission spectrophotometry (see e.g., Seager & Deming 2010, and references therein). In addition, the phasedependence of exoplanet IR flux provides observational constraints on their atmospheric temperature distribution (see the component labeled “phase” in Fig. 1) and, consequently, further insight into their atmospheric structures. For instance, from phasecurve IR measurements of the exoplanet HD 189733b, Knutson et al. (2007) derived the first longitudinal^{2} brightness distribution (BD) of an exoplanet; it indicates an offset hot spot in agreement with the predictions for hot Jupiters with shifted atmospheric temperature distributions because of equatorial jets (e.g., Showman & Guillot 2002).
Fig. 1 Schematic description of the different scanning processes observable for an occulted exoplanet. The green dotted lines indicate the scanning processes during the exoplanet occultation ingress/egress. The red dashed line indicates the scanning process that results from the exoplanet rotation and produces the exoplanet phase curve – this scanning appears longitudinal for an observer as long as the exoplanet spin is close to the projection plane, e.g., for a transiting and synchronized exoplanet. The component labeled “combined” shows the specific grid generated by these three scanning processes. 

Open with DEXTER 
Since then, significant theoretical developments have been achieved in modeling exoplanet atmospheres by combining hydrodynamic flow with thermal forcing (e.g., DobbsDixon et al. 2010; Rauscher & Menou 2010, 2012; Showman et al. 2009) and/or with ohmic dissipation (e.g., Batygin et al. 2011; Heng 2012; Menou 2012). However, other possible forcing factors are yet unmodeled. For example, the magnetic starplanet interactions can also have a significant role in this matter – although, to our knowledge no studies discuss in detail their impact on the exoplanet atmospheric structure so far. Nevertheless, magnetic interactions have to date been only observed at the stellar surface, in the form of chromospheric hot spots rotating synchronously with the companions (e.g., Lanza 2009; Shkolnik et al. 2005).
The observation of specific spatial features within an exoplanet atmosphere, such as hot spots or cold vortices, is essential for constraining its structure and for gaining further insight into its physics. Eclipses have proved to be powerful tools for “spatially resolving” distant objects, including binary stars (e.g., Warner et al. 1971) and accretion disks (e.g., Horne 1985). Previous theoretical studies introduced the potential of eclipse scanning^{3} (see Fig. 1) for exoplanets in order to disentangle atmospheric circulation regimes (e.g., Rauscher et al. 2007; Williams et al. 2006).
Ideally, the light curve of a transiting and occulted exoplanet with a nonzero impact parameter can enable a 2D surface brightness map of its day side. In fact, as represented in Fig. 1, such an exoplanet is scanned through several processes along its orbit. First, the exoplanet is gradually masked/unmasked by its host star during occultation ingress/egress. Secondly, the exoplanet rotation provides its phasedependent hemisphereintegrated flux (i.e., its phase curve), which constrains its BD in longitudinal slices – as long as the exoplanet spin is close to the projection plane, e.g., for a transiting and synchronized exoplanet. In particular, the phase curve is modulated by the orbital period as long as the exoplanet is tidally locked. The three scanning processes (ingress, egress and phase curve) provide thus complementary pieces of information that could ultimately constrain the BD over a specific “grid” (e.g., see the component labeled “combined” in Fig. 1). In this way, only a “uniform time offset”^{4} has been detected so far; Agol et al. (2010) showed that, assuming HD 189733b’s orbit to be circularized, its occultation was offset by 38 ± 11 s. This time offset is in agreement with the expected effect in occultation of the offset hot spot indicated by HD 189733b’s longitudinal map from Knutson et al. (2007).
The subset of exoplanet thermal IR observations that aims to characterize a planetary BD is growing. Currently, the Spitzer Space Telescope (Werner et al. 2004) has observed thermal phase curves for a dozen different exoplanets as well as the IR occultations of over thirty exoplanets. Among these, HD 189733b (Bouchy et al. 2005) is arguably the most favorable transiting exoplanet for detailed observational atmospheric studies; in particular, because its Kdwarf host is the closest star to Earth with a transiting hot Jupiter. This means the star is bright and the eclipses are relatively deep yielding favorable signaltonoise ratio (SNR). As such, HD 189733b represents a “Rosetta Stone” for the field of exoplanetology with one of the highest SNR secondary eclipses (Charbonneau et al. 2008; Deming et al. 2006), phasecurve observations (Knutson et al. 2007, 2009, 2012) and, consequently, numerous atmospheric characterizations (e.g., Deroo et al. 2010; Désert et al. 2009; Gibson et al. 2011; Grillmair et al. 2007; Huitson et al. 2012; Madhusudhan & Seager 2009; Pont et al. 2007; Redfield et al. 2008; Sing et al. 2011; Swain et al. 2008; Tinetti et al. 2007). Although HD 189733b’s atmospheric models are in qualitative agreement with observations, important discrepancies remain between simulated and observed light curves as well as between emission spectra (see e.g., Showman et al. 2009, Figs. 8 and 10). In addition, discrepancies exist between several published inferences – in particular molecular detections – which emphasize the impact of data reduction and analysis procedures (e.g., Désert et al. 2009; Gibson et al. 2011). Hence, we undertake a global analysis of all HD 189733’s public photometry obtained with the Spitzer Space Telescope for assessing the validity of published inferences (de Wit & Gillon, in prep.).
AOR’s description.
In this paper, we present the first secondaryeclipse scanning of an exoplanet showing a deviation from the occultation of a uniformlybright disk at the 6σ level, which we obtain from the archived Spitzer/IRAC 8μm data of the star HD 189733. In addition, we propose a new methodology of analysis for secondary eclipse scanning to disentangle the possible contributing factors of such deviations. As a result, we perform a new step toward mapping distant worlds by constraining consistently (i.e., simultaneously) HD 189733b’s shape, BD at 8 μm and system parameters.
At the time of submission, we learned about a similar study by Majeau et al. (2012), hereafter M12, focusing on the derivation of HD 189733b’s 2D eclipse map, using the same data but different frameworks for data reduction and analysis. Our study differs from M12 in three main ways. First, we find a deviation from the occultation of a uniformlybright disk at the 6σ level in contrast to the ≲ 3.5σ level deviation in the phasefolded light curves from Agol et al. (2010, see their Fig.12) used in M12. Secondly, this deviation has multiple possible contributing factors (i.e., not only a nonuniform BD but also the exoplanet shape or biased orbital parameters). Our study provides a framework for constraining consistently these contributing factors. Thirdly, and related to the second point, we do not constrain a priori the system parameters to the bestfit of a conventional analysis, nor the orbital eccentricity to zero; instead, we estimate the system parameters simultaneously with the BD. We compare the procedures and the results in Sect. 6.6.
We begin with a summary description of the Spitzer 8 μm data. In Sect. 2 we present our data reduction and conventional data analysis, i.e., for which we model HD 189733b as a uniformly bright disk. We present and assess the robustness of the structure detected in HD 189733b’s occultation in Sect. 3. In Sect. 4, we describe our new methodology for disentangling the possible contributing factors of such a structure. We present our results in Sect. 5. We, then, discuss in Sect. 6 the robustness of our results and the perspectives of our method. Finally, we conclude in Sect. 7.
2. Data reduction and conventional analysis
We present in this section the conventional data reduction and global data analysis performed for determining HD 189733’s system parameters (i.e., the orbital and physical parameters of the star and its companion) based on the Infrared Array Camera (IRAC: Fazio et al. 2004) 8μm eclipse photometry. We emphasize that we simultaneously analyze the whole data set in what we term a “global analysis”, instead of combining each separatelyanalyzed eclipse events. The global analysis approach helps to mitigate the effects of noise by extracting simultaneously common information between multiple light curves. Therefore, a global procedure also enables detection of previously unnoticed effects in the data. We begin with a summary description of the data sets; then, we introduce the analysis method and the physical models used for the parameter determination.
2.1. Data description and reduction
The eight secondary eclipses and the six transits of HD 189733b^{5} used in this study are described (by Astronomical Observation Requests, hereafter AOR) in Table 1. The data were obtained from November 2005 to June 2008 with IRAC at 8 μm. These are calibrated by the Spitzer pipeline version S18.18.0. The new S18.18.0 version enables improvements in the quality of data reduction over the original published data sets that used older Spitzer pipeline versions^{6}. Each AOR is composed of data sets; each of which corresponds to 64 individual subarray images of 32 × 32 pixels.
The data reduction consists in converting each AOR into a light curve; for that purpose, we follow a procedure (see Gillon et al. 2010b, hereafter G10, and references therein) that is performed individually for each AOR (de Wit & Gillon, in prep.). For each AOR, we first convert fluxes from Spitzer units of specific intensity (MJy/sr) to photon counts; then, we perform aperture photometry on each image with the IRAF^{7}/DAOPHOT software (Stetson 1987). We estimate the PSF center by fitting a Gaussian profile to each image. We estimate the best aperture radius (see Table 1) based on the instrument pointspread function (PSF^{8}), HD189733b’s, HD189733’s, HD 189733B’s and the skybackground flux contributions. For each image, we correct the sky background by subtracting from the measured flux its mean contribution in an annulus extending from 10 to 16 pixels from the PSF center. Then, we discard:

the first 30 min of each AOR for allowing the detectorsensitivitystabilization;

the few significant outliers to the bulk of the xy distribution of the PSF centers using a 10σ median clipping;

and, for each subset of 64 subarray images, the few measurements with discrepant flux values, background and PSF center positions using a 10σ median clipping.
Finally, for each AOR, the resulting light curve is the time series of the flux values averaged across each subset of 64 subarray images; while the photometric errors are the standard deviations on the averaged flux from each subset.
2.2. Photometry data analysis
We describe now the procedure used for constraining HD 189733’s system parameters from the light curves obtained after data reduction. We present the analysis method and the models – eclipse and systematic models – used to fit the light curves. In addition, we describe the procedure used for taking into account the correlated noise for each light curve (for details on the effect of correlated noise see, e.g., Pont et al. 2006).
2.2.1. Analysis method
We use an adaptive Markov Chain Monte Carlo (MCMC; see e.g. Ford 2006; Gregory 2005) algorithm. MCMC is a Bayesian inference method based on stochastic simulations that sample the posterior probability distribution (PPD) of adjusted parameters for a given fitting model. We use here the implementation presented in detail in Gillon et al. (2009, 2010a,b). More specifically, this implementation uses Keplerian orbits and models the eclipse photometry using the model of Mandel & Agol (2002). In addition, the simulated eclipse photometry is multiplied by a baseline, different for each timeseries, to take into account the systematics (see below).
We used in our first analysis the following jump parameters^{9}; the planet/star area ratio (R_{p}/R_{⋆})^{2}, the orbital period P, the transit duration (from the first to last contact) W, the time of minimum light T_{0}, , and the impact parameter b = a/R_{⋆}cosi (where a is the exoplanet semimajor axis). We assumed a uniform prior distribution for all these jump parameters and draw at each step a random stellar mass, M_{⋆}, based on the Gaussian prior: M_{⋆} = 0.84 ± 0.06 M_{⊙} (Southworth 2010). We estimate these parameters similarly to Agol et al. (2010), setting e = 0 based on the small inferred value of { ecosω,esinω } and theoretical predictions advocating for HD 189733b’s orbital circularization (e.g., Fabrycky 2010). We discuss the incidence of this assumption (hereafter COA for circularizedorbit assumption) in Sect. 3 and tackle it in detail in Sect. 5.
Fig. 2 IRAC 8μm HD 189733b’s transit and occultation ingress/egress photometry binned per 1 min and corrected for the systematics (black dots) with their 1σ error bars (red triangles) and the bestfitting eclipse model superimposed (in red). The green dots present the residuals from Agol et al. (2010, their Fig. 12), obtained using an average of the bestfit models from 7 individual eclipse analyses. Left: phasefolded transits show no significant deviation to the transit of a disk during ingress/egress. Right: phasefolded occultation ingress/egress deviate from the eclipse of a uniformlybright disk, highlighting the secondaryeclipse scanning of HD 189733b’s dayside. Top: phasefolded and corrected eclipse photometry. Bottom: phasefolded residuals. 

Open with DEXTER 
2.2.2. Eclipse models & limbdarkening
We model the transit assuming a quadratic limbdarkening law for the star, and the secondary eclipse assuming the exoplanet to be a uniformlybright disk. We draw the limbdarkening coefficients from the theoretical tables of Claret & Bloemen (2011), u_{1} = 0.0473 ± 0.0032 and u_{2} = 0.0991 ± 0.0036, based on the spectroscopic parameters of HD 189733 (T_{eff} = 5050 ± 50 K, log g = 4.61 ± 0.03 and , see Southworth 2010, and references therein). We add this a priori knowledge as a Bayesian penalty to our merit function, using as additional jump parameters the combinations c_{1} = 2u_{1} + u_{2} and c_{2} = u_{1} − 2u_{2}, as described in G10. The incidence of this coefficient choice does not affect our results (see Sect. 3).
2.2.3. Systematic correction models
IRAC instrumental systematic variations of the observed flux, such as the pixelphase or the detector ramp, are welldocumented (e.g., Désert et al. 2009, and references therein). At 8 μm, Si:Asbased detector showed a uniform intrapixel sensitivity (i.e., negligible pixelphase effect) but temporal evolution of their pixels gain (i.e., detector ramp). Although Agol et al. (2010) advocate using a doubleexponential for modeling the detector ramp, we find that the most adequate baseline model for the present AORs is the quadratic function of log (dt) introduced by Charbonneau et al. (2008). We discuss this matter in Sect. 3.
We also take into account possible lowfrequency noise sources (e.g., instrumental and/or stellar) with a secondorder timedependent polynomial.
The use of linear baseline models enables to determine their coefficients by linear leastsquares minimization at each step of the MCMC. For this purpose, we employed the singular value decomposition (SVD) method (Press et al. 1992).
2.2.4. Correlated noise
We take into account the correlated noise following a procedure similar to Winn et al. (2008) for obtaining reliable error bars on our parameters. For each light curve, we estimate the standard deviation of the bestfitting solution residuals for time bins ranging from 3.5 to 30 min in order to assess their deviation to the behavior of white noise with binning. For that purpose, the following factor β_{red} is determined for each time bin: (1)where N is the mean number of points in each bin, M is the number of bins, and σ_{1} and σ_{N} are respectively the standard deviation of the unbinned and binned residuals. The largest value obtained with the different time bins is used to multiply the error bars of the measurements.
3. Secondary eclipse: anomalous ingress/egress
3.1. Significance
We have detected an anomalous structure in the HD 189733b occultation ingress/egress residuals (see Fig. 2, bottomright panel); which shows that the observations deviate from the occultation of uniformlybright disk with a 6.2σ significance^{10}.
The IRAC 8μm photometry of the eclipse ingress/egress corrected for the systematics and binned per 1 min are shown in Fig. 2, with the bestfitting eclipse model superimposed (solid red line). The error bars (red triangles) are rescaled by β_{red} (~1.2) to take into account the correlated noise effects on our detection. In addition, we take advantage of the MCMC framework to account for the uncertainty induced by the systematic correction on the phasefolded light curves. For that purpose, we use the posterior distribution of the accepted baselines to estimate their median instead of using the bestfit model – which has no particular significance. Furthermore, we propagate the uncertainty of the systematic correction using the standard deviation of each bin from this posterior distribution; and increase their error bars accordingly (up to 20%).
We include in Fig. 2 the residuals from Agol et al. (2010, their Fig. 12). These are obtained using an average of the bestfit models from 7 individual eclipse analyses, not from a global analysis. For further comparison, the structure detected in HD 189733b’s occultation leads to a uniform time offset of 37 ± 6 s (~6σ), in agreement with Agol et al. (2010) estimate of 38 ± 11 s (~3.5σ), light travel time deduced.
3.2. Robustness
We test the robustness of our results against various effects including the baseline models, the limb darkening, HD 189733b’s circularizedorbit assumption (COA) and the AORs. Although below we mainly discuss the robustness of the structure in occultation, the robustness of system parameter estimates is assessed simultaneously but tackled in detail in Sect. 5.1.
As mentioned in Sect. 2, Agol et al. (2010) advocate using a doubleexponential for modeling the detector ramp; but we use the quadratic function of log (dt) introduced by Charbonneau et al. (2008). The reason is that by taking advantage of our Bayesian framework we show that the quadratic function of log (dt) is the most adequate for correcting the present AORs. In particular, we use two different information criteria (the BIC and the AIC see, e.g., Gelman et al. 2004) that prevent from overfitting based on the likelihood function and on a penalty term related to the number of parameters in the fitting model. (Note that the penalty term is larger in the BIC than in the AIC.) We obtain both a higher BIC (Δ BIC ~ 90) and a higher AIC (Δ AIC ~ 1.3) with the doubleexponential; this means the additional parameters do not improve the fit enough, according to both criteria. The most adequate ramp model is thus the less complex quadratic function of log (dt). Nevertheless, we assess the robustness of our results to different baseline models including the doubleexponential ramp, phasepixel corrections and sinusoidal terms. These different MCMC simulations do not significantly affect the anomalous shape found in occultation ingress/egress to within 0.5σ. The reason is that the time scale of the baseline models are much larger than the structure detected in occultation ingress/egress; which, therefore, has no incidence on the baseline models, and vice versa.
We assess the incidence of the priors on the limbdarkening coefficients; for that purpose, we perform MCMC simulations with no priors on u_{1} and u_{2} (see Sect. 2). Again, we observe no significant incidence to within 0.5σ. However, note that we outline in the next subsection the necessity of precise and independent constraints on the limbdarkening coefficient, even if these might appear to be wellconstrained by highSNR transit photometry.
We assess the incidence of assuming HD 189733b’s orbital circularization. For that purpose, we perform MCMC simulations with a free eccentricity. We observe no significant incidence to within 0.5σ for the jump parameters, except for the transit duration, the impact parameter, and (correlation detailed in Sect. 5.1). In addition, we observe a net drop of the anomalousshape significance which is explained by a compensation of the “uniform time offset”, enabled when relaxing the constraint on e (see next subsection and Sect. 5.1).
Finally, we also validate the independence of our results to inclusion or not of AOR subsets by analyzing different subsets of seven out of the eight eclipses. We observe relative significance decrease of ~; which are consistent with a significance drop due to a reduction of the sample. In particular, it shows that the anomalous shape in occultation ingress/egress is not due to one specific AOR.
Fig. 3 Schematic description of the anomalous occultation ingress/egress induced by the shape or the brightness distribution of an exoplanet. The red curve indicates the occultation photometry for a nonuniformlybright disk (hot spot in red). The yellow curve indicates the occultation photometry for an oblate exoplanet (yellow ellipse). Both synthetic scenarios show specific deviations from the occultation photometry of uniformlybright disk (black curve) in the occultation ingress/egress. 

Open with DEXTER 
3.3. Possible contributing factors
The anomalous shape detected in occultation ingress/egress is in agreement with the expected signature of the offset thermal pattern indicated by HD 189733b’s phase curve from Knutson et al. (2007). We discuss here the possible contributing factors – or origins – of this anomalous shape to propose further a consistent methodology for analysing secondaryeclipse scanning (see Sect. 4).
An anomalous occultation ingress/egress possibly emerges from: (1) a noncircular projection of the exoplanet at conjunctions; (2) a biased eccentricity estimate (i.e., uniform time offset); and/or (3) a nonuniform brightness distribution (BD). We present in Fig. 3 a schematic description of the anomalous occultation ingress/egress induced by a noncircular projection of the exoplanet at conjunctions (yellow) and a nonuniform BD (red). Both synthetic scenarios show specific deviations from the occultation photometry of uniformlybright disk (black curve) in the occultation ingress/egress.
(1) We reject the noncircularprojection concept because the transit residuals show no anomalous structure (Fig. 2, bottomleft panel). This is in agreement with current constraints on the HD 189733b oblateness (projected oblateness below 0.056, 95% confidence, Carter & Winn 2010) and winddriven shape (expected to introduce a light curve deviation below 10 ppm, see Barnes et al. 2009). In particular, our transit residuals constrain HD 189733b’s oblateness (95% confidence) below 0.0267 in case of a projected obliquity of 45° and below 0.147 in case of a projected obliquity of 0° – based in Fig. 1 of Carter & Winn (2010). No significant deviation in transit ingress/egress means that HD 189733b’s shapeinduced effects in occultation ingress/egress are negligible, because it is expected to be about one order of magnitude lower than in transit (i.e., effect ratio ∝ I_{p}/I_{⋆}, with I_{p} and I_{⋆} respectively the exoplanet and the star mean intensities at 8 μm). We therefore assume the exoplanet to be spherical, for a further analysis of the present data.
We, therefore, emphasize that for constraining an exoplanet shape, one needs independent a priori knowledge on the host star limbdarkening (e.g., Claret & Bloemen 2011) to avoid overfitting its possible signature in transit ingress/egress.
We highlight that the phase curve also constrains a companion shape, similarly to the ellipsoidal light variations caused by a tidallydistorted star (see Kopal 1959; Russell & Merrill 1952). The contributions of the shape and the BD of an exoplanet to its phase curve may be disentangle because of their different periods (see e.g., Faigler & Mazeh 2011). Indeed, for a synchronized exoplanet, the shapeinduced modulation has a period of P/2 – twice the same projection an orbital period, while the brightnessinduced modulation as a period of P. HD 189733b’s phasecurves show mainly a Pmodulation; what is in agreement with our previous constraint on HD 189733b’s shape.
(2) and (3) Williams et al. (2006) introduced the concept of uniform time offset to emphasize that a small variation of eccentricity can partially mimic the anomalous occultation ingress/egress of a nonuniformlybright exoplanet. Therefore, an anomalous occultation ingress/egress can be partially compensated leading to biased estimates of and . In fact, we show in Sect. 5 that not only e enables such partial compensations.
We present in Sect. 4 the methodology proposed to disentangle the factors possibly responsible for partial compensations.
4. Analysis of HD 189733b’s scans
Our second analysis aims to investigate the contributing factors of the anomalous occultation ingress/egress we detect at the 6σ level (see Sect. 3). We require additional information to constrain the role of the remaining possible contributing factors. We, therefore, take advantage of the phase curve – in addition to the secondaryeclipse scanning – to constrain the exoplanet BD simultaneously with the system parameters.
We describe below the framework used for analysing HD 189733b’s secondaryeclipse scanning. We begin with an introduction of the data. Then we describe the analysis method and the models used for constraining HD 189733b’s BD at 8 μm; we recall that HD 189733b’s shape has been constrained in Sect. 3.
Fit properties for different fitting models of HD 189733’s photometry in the Spitzer/IRAC 8μm channel.
4.1. HD 189733b’s scans
In this second analysis, we use the corrected and phasefolded light curves from Sect. 3 (Fig. 2) and the phase curve from (Knutson et al. 2007) corrected for stellar variability (see Agol et al. 2010, Fig. 11). As discussed in Sect. 3, the transit and the phase curve enable to constrain respectively HD 189733b’s shape at conjunctions and its BD. Note that we discard the first third of the phase curve since it is strongly affected by systematics (i.e., detector stability).
To constrain the system parameters, we choose to use HD 189733b’s phasefolded transits instead of using the system parameter estimates from Sect. 2 (see Table 2, Col. 2) for two reasons. First, these estimates are under the form of 1D Gaussian distribution while light curves provide a complex posterior probability distribution (PPD) over the whole parameter space. Secondly, these estimates are affected by HD 189733b’s COA (see implications in Sect. 5). For those reasons, it is relevant to simultaneously analyse HD 189733b’s scans and estimate the system parameters, to be consistent with our methodology of a “global” analysis and to avoid propagation of bias through inadequate priors, i.e., inadequate assumptions.
4.2. Analysis method
The procedure we develop for analysing HD 189733b’s scans uses a second, new, MCMC implementation. This implementation differs from the one introduced in Sect. 2 in that it models a nonuniformlybright exoplanet. This implementation uses Keplerian orbits and models the transit photometry using the model from Mandel & Agol (2002). For that reason, it uses the jump parameters of the conventionalanalysis method together with new jump parameters for the exoplanet brightness models (see below). However, note that we choose to use the stellar density, ρ_{⋆}, instead of W as jump parameter; because it relates directly to the orbital parameters (Seager & MallénOrnelas 2003). Furthermore, because we use the phasefolded light curves, the main constraint on the orbital period is missing. We, therefore, use a uniform prior on P, centered on the conventionalanalysis estimate. In particular, we use a large prior (with an arbitrary symmetric extension of 10 σ_{P}, where σ_{P} is the estimated uncertainty on P, see Table 2) to prevent our results to be affected by the assumption underlying this estimate. Note however that this has no incidence on our further results – the reason is that the orbital period is highlyconstrained by the transit epochs and, therefore, is not affected by effects on the occultation.
4.3. Nonuniformlybright exoplanet: light curve models
We model the phase curve and the secondary eclipse by performing a numerical integration of the observed exoplanet flux. Model approximations include: ignoring the time variability of the target atmosphere (in line with atmospheric models, e.g., Cooper & Showman 2005) – especially because we use here phasefolded light curves; ignoring the planet limb darkening – expected to be weak for hot Jupiters, in particular at 8 μm – and assuming HD 189733b’s rotational period to be synchronized with its orbital period – synchronization occurs over ~10^{6} yr for hot Jupiters, e.g., Winn (2010).
From the exoplanet BD and the system orbits, we model the flux temporal evolution of the exoplanet by sampling its surface with a grid of 2N points in longitude (φ) and N + 1 points in latitude (θ). We fix N to 100 to mitigate numerical effects up to 10^{3} the secondaryeclipse depth, i.e. below 4 ppm; in comparison, the data photometric precision is ~130 ppm. (We validate the independence of our inferences to higher grid resolutions.)
Our brightness models (described in the next subsection) are composed of several modes/degrees, e.g., spherical harmonics. At each step of the MCMC, we (1) simulate the system orbits, (2) estimate the light curve corresponding to each of the brightness model modes and (3) determine the mode amplitude based on their simulated light curve by leastsquares minimization using the SVD method. Practically, to employ consistently the SVD method within a stochastic framework, we use the covariance matrix to perturb the modeamplitude estimates.
We model here the transit as in our first analysis (see Sect. 2) because we show that the projection of HD 189733b at conjunctions does not deviate significantly from a disk (see Sect. 3).
4.4. Nonuniform brightness models
To model the broad patterns of HD 189733b’s BD, we use two groups of models: the spherical harmonics and toy models. For the spherical harmonics, we use two additional jump parameters per degree for direction (i.e., longitude and latitude). The generalized formulation of the BD, Γ_{SH,d}(φ,θ), can thus be written as; (2)where is the real spherical harmonic of degree l and order 0, directed to Δφ_{l} in longitude and Δθ_{l} in latitude. I_{l} is the lthmode amplitude that is estimated at each step of the MCMC using the perturbed SVD method. For instance, the additional parameters for a dipolar fitting model compared to the conventionalanalysis method are: Δφ_{1} and Δθ_{1} as jump parameters and I_{1} as a linear coefficient; I_{0} is the amplitude of the uniformlybright mode, included in the conventionalanalysis method.
On the other hand, we use toy models that enable modeling a thermal pattern – a hot or cold spot – of various shape. Their expression can be written simply as; where φ_{°} and θ_{°} are respectively the longitude and latitude relative to the position of the model extremum, i.e. φ_{°} = f(φ,Δφ,α,β,γ) and θ_{°} = f(θ,Δθ,α,β,γ). Δφ and Δθ are respectively the longitudinal and the latitudinal shift of the model peak from the substellar point. α,β and γ parametrize the shape of the hot/cold spot. These toy models add to the conventional analysis method five jump parameters (Δφ,Δθ,α,β,γ) and a linear coefficient (I_{1}).
Finally, we emphasize that the use of several brightness models is critical to assess the modeldependence of our results. We discuss this important matter in Sect. 6.
Fig. 4 Incidence of assuming an exoplanet to be uniformly bright and its orbit to be circularized. a) Marginal posterior probability distribution (PPD, 68% and 95%confidence intervals) of and that shows an unusual correlation. b) Marginal PPD (95%confidence intervals) of ρ_{⋆} and b that highlights the increase of adequate solutions enabled by the additional dimensions of the parameter space probed when relaxing the circularized orbit assumption. c) Deviations in occultation ingress/egress from the medianfit model for the individual perturbations of b, , and ρ_{⋆} by their estimated uncertainty (see Table 2, Col. 2). It outlines that the system parameters enable compensation of an anomalous occultation that emerges from, e.g., a nonuniformlybright exoplanet (see Fig. 3) leading to biased estimates of the system parameters. 

Open with DEXTER 
5. Results
For HD 189733’s system, we find that relaxing the eccentricity constraint and using more complex brightness distributions (BDs) lead to lower stellar/planetary densities and a more localized and latitudinallyshifted hot spot. In particular, we find that the more complex HD 189733b’s brightness model, the larger the eccentricity, the lower the densities, the larger the impact parameter and the more localized and latitudinallyshifted the hot spot estimated. This “ebρ_{⋆}BD correlation” is of primary importance for data of sufficient quality. We present in this section our results for increasing model complexity to gain insight into the incidence of the model underlying assumptions – e.g., the COA and the uniformlybright exoplanet assumption (hereafter UBEA). In particular, we relate their incidence to the ebρ_{⋆}BD correlation.
We gather the system parameter estimates for different fitting models in Table 2; it shows the median values and the 68% probability interval for our jump parameters. We compute our estimates based on the posterior probability distribution (PPD) of global MCMC simulations, i.e., not as a weighted mean of individual transit or eclipse analyses. Our conventional analysis estimates are in good agreement with previous studies (e.g., Agol et al. 2010; Triaud et al. 2009; Winn et al. 2007). We discuss further the system parameter estimates obtained from our second analysis.
Fig. 5 Phasefolded IRAC 8μm HD 189733b’s photometry binned per 1 min and corrected for the systematics (black dots) with their 1σ error bars (red triangles) and the bestfitting eclipse models superimposed. Left: phasefolded transits. Right: phase curve and phasefolded occultations that show the benefit of using nonuniform brightness model. 

Open with DEXTER 
Fig. 6 Incidence of the brightness model complexity on the systemparameter posterior probability distribution (PPD), using unipolar models. a), b) Marginal PPDs (68% and 95%confidence intervals) of and for two unipolar brightness models, respectively, with a fixed and large structure (Γ_{SH,1}) and with freeconfinement structure (Γ_{2}). A comparison with Fig. 4a shows that is constrained closer to zero. The reason is that nonuniform brightness models enable the exploration of additional dimensions of the parameter space and, therefore, provide additional adequate combinations of the contributing factors to compensate an anomalous occultation. In particular, the uniform time offset is now mainly compensated by a nonuniform BD, rather than by as in conventional analysis. It shows also the evolution of the marginal PPD toward larger when using more complex brightness models – which enable more localized brightness structure. c), d) Marginal PPDs of ρ_{⋆} and for, respectively, Γ_{SH,1} and Γ_{2}. These show the impact of the brightness distribution on the retrieved system parameters (i.e., ebρ_{⋆}BD correlation); in particular, it outlines the possible overestimation of ρ_{⋆} by 5% (i.e., at 6σ of the conventional estimate) when using extended brightness models, e.g., a uniform brightness distribution. 

Open with DEXTER 
Fig. 7 Estimate of HD 189733b’s global brightness distribution in the IRAC 8μm channel using the Γ_{SH,1} brightness model. Left: relative brightness distribution at HD 189733b’s dayside. Right: dayside standard deviation. Because of its fixed and large structure, the Γ_{SH,1} brightness model is wellconstrained in amplitude (by the occultation depth) and in longitudinal localization (by the phase curve). However, it is less constrained in latitude (by the secondaryeclipse scanning) than more confined model (schematic description in Fig. 3), e.g., Γ_{2} (see Fig. 8). These modelinduced constraints are observable on the dayside standard deviation; which is significantly lower than for more complex brightness models (see Figs. 8, 11a and b) and is related to its gradient, because the BDs accepted along the MCMC simulations differ from each other mainly in latitudinal orientation (see Fig. A.1a). 

Open with DEXTER 
Fig. 8 Estimate of HD 189733b’s global brightness distribution in the IRAC 8μm channel using the Γ_{2} brightness model. Left: relative brightness distribution at HD 189733b’s dayside. Right: dayside standard deviation. Because of its increased complexity, the Γ_{2} brightness model enable more localized structure that are less constrained in amplitude (by the occultation depth) and in longitudinal localization (by the phase curve) than largeandfixedstructure model, e.g., the Γ_{SH,1} model (see Fig. 7). However, it is wellconstrained in latitude by the secondaryeclipse scanning that is sensitive to confined brightness structure. These modelinduced constraints are observable on the dayside standard deviation that shows a maximum at the brightness peak localization and extended wings towards west and east along the planetary equator, because the BDs accepted along the MCMC simulations mainly affect the former by their amplitude change and the latter by their structure change (see Fig. A.1b). 

Open with DEXTER 
Fig. 9 Dependence and significance of the brightness peak localization. a) Marginal PPD (68% and 95%confidence intervals) of the brightness peak longitude, Δφ, and ecosω for the Γ_{SH,1} brightness model – the simplest nonuniform brightness model used in this study. This shows the correlation between the brightness model and the orbital eccentricity. This correlation emerges from enabling compensation of the anomalous occultation with a larger set of contributing factors, i.e., by including nonuniform brightness models. In particular, the uniform time offset is now mainly compensated by an nonuniformlybright model, rather than by ecosω as in conventional analysis (see Fig. 3, Eq. (6) and Fig. 4c). b) Marginal PPDs (68%confidence intervals) of the brightness peak localization for the Γ_{SH,1} and Γ_{2} brightness models. It shows that the brightness peak localization is modeldependent. For example, the longitudinal Γ_{SH,1} peak localization is constrained by the phase curve because of its large and constant extension; while the free extension of the Γ_{2} model relaxes this longitudinal constraint (see Sect. 5). Therefore, the light curve of an exoplanet does not constrain uniquely its brightness peak localization; furthermore, the brightness peak localization is not an adequate parameter to characterize complex exoplanet brightness distributions. 

Open with DEXTER 
Fig. 10 Incidence of the brightness model extension on the system parameter posterior probability distribution (PPD), using multipolar models. a), b) Marginal PPDs (68% and 95%confidence intervals) of and for the quadrupolar and octupolar brightness models, respectively. c), d) Marginal PPDs of ρ_{⋆} and for the quadripolar and octopolar brightness models, respectively. These confirm the trend observed in Sect. 5.1 (see Fig. 6) towards larger and lower stellar density when increasing the complexity of the brightness distribution, i.e., when enabling more localized structures. 

Open with DEXTER 
Fig. 11 Estimate of HD 189733b’s global brightness distribution in the IRAC 8μm channel using multipolar brightness models. Left: relative brightness distribution at HD 189733b’s dayside. Right: dayside standard deviation. a) Estimate using the Γ_{SH,2} brightness model. b) Estimate using the Γ_{SH,3} brightness model. These confirm the trend toward a more localized and latitudinallyshifted hot spot when increasing the brightnessmodel complexity. Therefore, together with Fig. 10 these outline that the more complex HD 189733b’s brightness model, the larger the eccentricity, the lower the densities, the larger the impact parameter and the more localized and latitudinallyshifted the hot spot estimated. 

Open with DEXTER 
5.1. Assuming HD 189733b to be uniformly bright
We show here the similar effects of the COA and the UBEA. Both assumptions prevent from exploring dimensions of the parameter space and, therefore, lead to more localized, and possibly biased, PPD.
We first present an unusual correlation between and (see Fig. 4a). It emerges from the partial compensation of the anomalous occultation ingress/egress enabled primarily by adequate combinations of and , for conventional analyses. In particular, these parameters enable both to shift the occultation by (6)(i.e., uniform time offset) and to change its duration (7)what is sufficient to partially compensate, e.g., the effect of a hot spot – compare the red and the black curves in Fig. 3. However, ρ_{⋆} and b are also affected by the occultation shape and timing, as they constrain respectively a/R_{⋆} and i. We emphasize this point in Fig. 4c. Figure 4c shows the deviations in occultation ingress/egress induced by individual perturbations of b, , and ρ_{⋆} by their estimated uncertainty (see Table 2, Col. 3). This indicates that adequate combinations of these parameters mimic partially an anomalous occultation ingress/egress (see Fig. 2, bottomright panel, and Fig. 4c). Therefore, the COA and the UBEA also affect the marginal PPD of {ρ_{⋆},b} by inhibiting the exploration of dimensions of the parameter space. We present in Fig. 4b the extension of the {ρ_{⋆},b} marginal PPD that results from relaxing the COA. Similarly, we expect that the relaxation of the UBEA would significantly affect the systemparameter PPD – compare Figs. 3 and 4c. This motivates our global approach that constrains simultaneously the possible contributing factors of an anomalous occultation ingress/egress and, therefore, prevents the incidence of supplementary assumptions – e.g., the COA and the UBEA.
5.2. Assuming HD 189733b to be nonuniformly bright
We show in this subsection the incidence of relaxing the UBEA on the fit improvement and on the systemparameter PPD. We present in Fig. 5 the phasefolded IRAC 8μm photometry of HD 189733b, corrected for the systematics with the bestfitting eclipse models for a uniformly (blue) and a nonuniformly (green) bright exoplanet superimposed. The bestfitting nonuniformlybright eclipse model is shown for the Γ_{SH,1} model – chosen arbitrarily, as the nonuniformlybright models provide similar fits (see Table 2). In particular, these models are significantly more adequate according to both the BIC and the AIC, see Table 2 (odds ratio: ~10^{36}).
First, we introduce the results obtained using unipolar (i.e., with one spot on the planetary dayside) BDs to gain insight into the incidence of relaxing the UBEA. Then, we introduce the results obtained using multipolar BDs to assess the validity of trends observed in the unipolarmodel results.
5.2.1. Unipolar brightness distribution
We first present the incidence of relaxing the UBEA on the system parameters. For that purpose, we show in Figs. 6a and c respectively the marginal PPDs of and for the Γ_{SH,1} model, and the ones for the Γ_{2} model in Figs. 6b and d. We superimpose in Figs. 6c and d the 95%confidence interval obtained for the uniformlybright model to extend our previous observations regarding the incidence of relaxing supplementary assumptions, e.g., the COA (see Sect. 5.1; Fig. 4b). We observe the increases of and b and the decrease of ρ_{⋆} while is constrained closer to zero (see Table 2). The reason is that the compensation of HD 189773b’s anomalous occultation is now also enabled by the nonuniform brightness models; which provide a better compensation than e solely – with (i.e., uniform time offset) for conventional analysis. Therefore, numerous combinations of e and BDbased compensations are adequate. In other words, e and the BD are correlated as highlighted by the PPD in Fig. 9a. Finally, we note a progressive evolution of the systemparameter PPD with the brightnessmodel complexity (from uniform to Γ_{2}). We assess further the validity of these observations, using spherical harmonics of higher degree.
We now turn to HD 189733b’s BD. We show the dayside estimates for the Γ_{SH,1} and Γ_{2} models with their corresponding uncertainties in Figs. 7 and 8 respectively; we focus on HD 189733b’s dayside as it is effectively constrained by the combination of the phase curve and the secondaryeclipse scanning. In particular, note that Figs. 7 and 8 present HD 189733b’s brightness relative to HD 189733’s hemisphereaveraged brightness in the IRAC 8μm channel. In addition, the figures are timeaveraged; our estimates aim to approach the global pattern of HD 189733b’s BD based on eight snapshots taken from November 2005 to June 2008. Finally, these estimates correspond to the median and standard deviation of the map trials accepted along the MCMC simulations, similarly to our approach for the corrected and phasefolded light curves (see Sect. 3).
Both models retrieve a spatial feature in HD 189733b’s BD; which corresponds to a hot spot. The Γ_{SH,1} model retrieves a hot spot shifted to the east of the substellar point, see Fig. 7. The Γ_{2} model retrieves a hot spot shifted to the east of the substellar point but also away from the equator, see Fig. 8. However, we cannot discuss the direction of this latitudinal shift due to a NorthSouth ambiguity (Agol, priv. comm.).
The BD estimates shown in Figs. 7 and 8 are significantly different both in pattern and in intensity. These differences are due to the estimate modeldependence; which motivates the use of different fitting models to enable a thorough discussion. For example, brightness models with nonconstant structure (“complex”, i.e., in opposition to a dipole) are less constrained by a phase curve that is only dependent on the hemisphereintegrated brightness. To emphasize these modelinduced constraints, we present in Fig. A.1 animations showing compilations of dayside BDs accepted along the MCMC simulations for the Γ_{SH,1} and Γ_{2} models. These compilations show that (1) the amplitude and (2) the longitudinal localization for the Γ_{SH,1}BD are more constrained than for the Γ_{2} model (by the occultation depth and by the phase curve, respectively) because of its fixed and large structure. However, (3) the Γ_{SH,1} model is less constrained in latitude (by the secondaryeclipse scanning) than the more complex Γ_{2} model which enables more confined structures that induce larger deviations in occultation ingress/egress. For that reason, the brightness peak localization for the Γ_{2} model is wellconstrained in latitude (see Fig. 9b), while for the Γ_{SH,1} model it is wellconstrained in longitude.
This shows that the light curve of an exoplanet does not constrain uniquely its brightness peak localization without a priori assumption (e.g., assuming a dipolar BD). Therefore, we will further refer to our brightnessdistribution estimates instead of the brightness peak localization; which is not representative of complex BDs, in addition to being modeldependent. Nevertheless, we propose in Sect. 6.4 another unidimensional parameter to replace the brightness peak localization.
Finally, note that these modelinduced constraints are also observable on the dayside standard deviation; which is significantly lower for Γ_{SH,1} model than for more complex models. In particular, the standard deviation distribution for the Γ_{SH,1} model is related to its gradient – with a larger variation from the brightness peak localization along the latitude axis than along the longitude axis, because the BDs accepted along the MCMC simulations differ from each other mainly in (latitudinal) orientation, see Fig. A.1a. This is in contrast with the standard deviation distribution for the Γ_{2} model that shows a maximum at the brightness peak localization and extended wings towards west and east along the equator; because the BDs accepted along the MCMC simulations mainly affect the former by their amplitude change and the latter by their structure change (see Fig. A.1b).
5.2.2. Multipolar brightness distribution
We observe an evolution of our inferences (systemparameter PPD and BD) when increasing the complexity of our fitting model. To assess the validity of this observation, we present here the results obtained when using spherical harmonics up to the degrees 2 (quadrupole) and 3 (octupole).
We present in Figs. 10a and c respectively the marginal PPDs of and for the Γ_{SH,2} model, and in Figs. 10b and d for the Γ_{SH,3} model. These PPDs appear as intermediate steps between the results obtained with the Γ_{SH,1} and Γ_{2} models (see Fig. 6). This, therefore, confirms the strong incidence of HD 189733b’s brightness model on the retrieved system parameters, because of the ebρ_{⋆}BD correlation.
We show the dayside brightness estimates for the Γ_{SH,2} and Γ_{SH,3} models with their uncertainty in Figs. 11a and b. As observed for the systemparameter PPDs, the BDs also appear as intermediate steps. The major evolutions are the shrinking of the structure retrieved and its shift away from the equator.
These progressive evolutions of both the systemparameter PPD and the retrieved brightness structure show that, for HD 189733’s system, relaxing the eccentricity constraint and using more complex BDs lead to lower stellar/planetary densities and more localized and latitudinallyshifted hot spot. In particular, we find that the more complex HD 189733b’s brightness model, the larger the eccentricity, the lower the densities, the larger the impact parameter and the more localized and latitudinallyshifted the hot spot estimated. We discuss the significance of these suggestions in the next section.
6. Discussion
6.1. The most adequate model
We present in Sect. 5 the results obtained using several fittingmodels for assessing the modeldependence of our inferences. As a result, this enables us to reveal the significant correlation between the system parameters and the BD of an exoplanet (i.e., “ebρ_{⋆}BD correlation”). In particular, we show a progressive evolution of the systemparameter PPD and the BD estimate when increasing the complexity of the fitting model. For that reason, we discuss below the relevance of our modelcomplexity increase, which may ultimately lead to overfitting the data.
We take advantage of our Bayesian framework using the BIC and the AIC. Both information criteria are in favor of models that relax the assumptions of circularized orbit and uniformly bright HD 189733b. In particular, the Γ_{SH,1} and Γ_{SH,3} models are favored; the BIC insignificantly favors Γ_{SH,3} (odds ratio ~ 1.02) while the AIC significantly favors Γ_{SH,3} (odds ratio ~ 3.5). Because not decisive, the information criteria only suggest that the Γ_{SH,3} model provides the most adequate constraints on HD 189733’s system. (Note that theoretical studies favor the AIC, e.g., Burnham & Anderson 2002; Yang 2005).
In particular, it suggests that HD 189733b’s hot spot is shifted both east of the substellar point and away from the equator and HD 189733b’s density has been overestimated by 3.6%. Furthermore, it suggests that HD 189733b’s orbit is possibly not fully circularized (), although its eccentricity is consistent with zero. This emphasizes that the assumption of circularized orbit has to be continuously assessed, because of data of constantly increasing quality; including for old hot Jupiters that may show a hint of eccentricity (e.g., CoRoT16b, see Ollivier et al. 2012). Finally, for dataquality reason, the interpretation of HD 189733b’s BD has to focus on global trends: the presence of an asymmetrical hot spot.
HD 189733b’s dayside presents a shifted hot spot. The eastward shift is in agreement with the literature: (1) with previous derivations, from HD 189733b’s phase curve (Knutson et al. 2007) and an eclipse timing constraint (Agol et al. 2010); and (2) atmospheric models suggesting a superrotating equatorial jet (e.g., Showman et al. 2009). In opposition, the suggested shift away from the equator is new. The smallscale origin of this latitudinal asymmetry remains unconstrained because we use largescale brightness models to be consistent with the data quality. For that reason, additional observations would be required to improve our understanding of HD 189733b’s atmosphere (see Sect. 6.5); in particular, its yet unmodeled interaction with HD 189733 (see e.g., Lecavelier des Etangs et al. 2012) that could induce unexpected thermal patterns, e.g., asymmetric patterns in its BD. For example, magnetic starplanet interactions may lead to energy dissipation due to the stellar field penetration into the exoplanet envelope (e.g., Laine et al. 2008) and to extensive energy injections into the auroral zones of the exoplanet from magnetic reconnections (e.g., Ip et al. 2004) – similarly to the JupiterIo flux tube (e.g., Bigg 1964). However, such magnetic reconnections have so far been only observed at the stellar surface, in the form of chromospheric hot spots rotating synchronously with the companions (e.g., Lanza 2009; Shkolnik et al. 2005).
6.2. Adequacy of conventional analyses
We discuss here the possible limitation of conventional analyses for interpreting light curves of “sufficient” data quality. We highlight the ebρ_{⋆}BD correlation and demonstrate in this context the significant impact on the system parameter estimates of assuming HD 189733b to be uniformly bright – conceptually similar to neglecting the limb darkening of the host star for transitphotometry analysis. Nevertheless, the significance of this impact is related to the significance of HD 189733b’s secondaryeclipse scanning; which is enabled by the highSNR Spitzer/IRAC 8μm photometry. This detection, therefore, motivates using more complex fitting models and outlines the limitation of conventional analyses because of the ebρ_{⋆}BD correlation; which would have remained hidden if the photometric precision was less. In this context, a “sufficient” photometric precision requires resolving the occultation ingress/egress; therefore, it has to be about one order of magnitude less than the occultation depth, for a time bin about one order of magnitude less than the occultation ingress/egress duration. For less photometric precisions, modelling the eclipse using a uniformlybright disk, therefore, is adequate, so are the conventional analyses.
Finally, we briefly outline that the conventional assumption of a uniformlybright exoplanet could also affect the inferred planetary interior models; because ρ_{⋆} is possibly affected and, therefore, ρ_{p} is too (e.g., suggested 3.6%overestimation for HD 189733b’s).
6.3. Relevance of reconstruction methods
We discuss here the relevance of mapping methods based on direct reconstruction from an exoplanet “slicing” (e.g., the slice method introduced in M12). We emphasize that a 2D map is not directly accessible over a specific “grid” (see Fig. 1). The main reason is that such method requires the use of fixed system parameters, which are derived from a conventional analysis. Therefore, it leads to biased and overlyprecise estimates by neglecting the correlation between the system parameters and the exoplanet BD – and shape, see Sect. 5. Note that in addition, it also neglects the exoplanet rotation during its occultation (~12° for HD 189733b) and does not take advantage of its phase curve (see M12). Finally, the degeneracy induced by the limited data (N slices by scanning processes while roughly N^{2} cells over the map, see M12) implies the use of a priori constraint on the BD.
6.4. Reducing brightness distributions to 1Dparameters
We have shown in Sect. 5.2.1 that the light curve of an exoplanet does not constrain uniquely its brightness peak localization. We discuss here the reasons why we strongly advocate discussing the BD estimates and, if necessary, using with care the dayside barycenter as a representative parameter.
Discussing the significance of a hot spot with complex structure from a 2D map is difficult; but may be simplified using a unidimensional parameter. We show in Sect. 5.2.1 that the brightness peak localization poorly represents complex BDs, in addition to being modeldependent. Therefore, we investigate using the dayside barycenter to replace the brightness peak localization. The reason is that the dayside brightness barycenter weights the BD according to the geometrical configuration at superior conjunction (i.e., it contains partial 2D information).
We show in Fig. 12 the marginal PPDs (68%confidence intervals) of the brightness peak localization for the Γ_{SH,1} and Γ_{2} brightness models. A comparison with the marginal PPDs of the brightness peak localization (Fig. 9b) shows the reduced modeldependence of the dayside barycenter. In particular, it shows a lessextended PPD for the Γ_{2}model barycenter; because this extension for the brightnesspeaklocalization PPD emerges from the model wings – weighted by the dayside barycenter. In addition, it shows the shift and slight shrinking of the Γ_{SH,1} PPDs that reflect the barycenter weighting according to the geometrical configuration at superior conjunction; map cells closer to the substellar point have more weight. This emphasizes the primary drawback of the dayside barycenter that is to attenuate the offset of BDs. This recalls that unidimensional parameters cannot stand adequately for complex BDs and, therefore, have to be used complementary to BD estimates.
6.5. Perspectives
We show that conventional analysis – which assume uniformlybright exoplanets – may lead to biased estimates of the system parameters. We introduce below two of the possible applications of our method in this context: (1) IR multiwavelength observations of HD 189733b for improving the constraints on the system parameters and, ultimately, for yielding a timedependent 3D map and (2) observations, in the visible, of targets with apparently high albedos.
(1) The atmospheric layers scanned during an occultation ingress/egress are wavelengthdependent. Therefore, additional high quality occultations at different wavelengths (i.e., optical depths) would contain the same information about the system parameters while different information about the BD. These would enable us to improve our constraints on the system parameters while gaining insights into the atmospheric physics of HD 189733b from its thermal patterns at different depths.
In particular, the origin of HD 189733’s hot spot could be revealed from observations at shorter wavelengths. These observations would refer to a higher temperature zone whose localization and extent could disentangle between its possible origins, e.g., radiative or magnetic.
Fig. 12 Reducing brightness distributions to unidimensional parameters. Marginal PPDs (68%confidence intervals) of the dayside barycenter brightness peak localization for the Γ_{SH,1} and Γ_{2} brightness models. A comparison with the marginal PPDs of the brightness peak localization (Fig. 9b) shows the reduced modeldependence of the dayside barycenter. In particular, it shows a lessextended PPD for the Γ_{2}model barycenter; because this extension for the brightnesspeaklocalization PPD emerges from the model wings, weighted by the dayside barycenter. In addition, it shows the shift and slight shrinking of the Γ_{SH,1} PPD that reflect the barycenter weighting based on to the geometrical configuration at superior conjunction; map cells closer to the substellar point have more weight. 

Open with DEXTER 
In addition, the time variability of these spatial features could also be targeted. Indeed, as discussed in Rauscher et al. (2007), JWST/NIRspec performance should allow high significance detection of anomalous occultation ingress/egress based on a unique eclipse (e.g., for the grating centered at 4 μm). The time variability would then be assessed from the BDs derived for different occultations. In other words, the future exoplanetaryatmosphere investigations of spacedbased facilities like JWST (Clampin 2010) and EChO (Tinetti et al. 2012) could ultimately lead to timedependent 3D maps of distant worlds.
(2) Extension of our method to observations in the visible would also be valuable, both for ground and spacebased observations. For example, the Kepler mission provides additional targets that would be relevant for similar analyses, thanks to its outstanding photometry quality. Kepler provides broadband visible photometry and is, therefore, capable of observing the reflected light from exoplanet atmospheres. Recent studies found high geometric albedos for two hot Jupiters (Berdyugina et al. 2011; Demory et al. 2011; Kipping & Bakos 2011). In addition, Demory et al. (2011) suggest hazes/cloud and Rayleigh scattering as plausible explanations for such high albedos. Recently, Madhusudhan & Burrows (2012) provide a theoretical framework for interpreting geometrical albedos from phase curves, which is is indicative of the scattering and absorptive properties of the atmosphere. We emphasize that the albedo origins could make the exoplanet nonuniformly bright and modify its occultation ingress/egress, similarly to the effect of an hot spot. Therefore, our study could constrain independently the geometrical albedo by scanning an exoplanet dayside, what would be a complementary asset to Madhusudhan & Burrows (2012)’s framework. Furthermore, our study will be necessary to constrain consistently the system parameters from highSNR light curves, see Sect. 6.2.
6.6. Comparison with Majeau et al. (2012)
We obtain qualitatively similar results as M12: an offset hot spot. Nevertheless, significant differences exist between both studies. (1) The starting point of our second analysis is our detection of HD 189733b’s eclipse scanning at the 6σ level (in contrast to their ~3.5σ). (2) Eclipse scanning has multiple possible contributing factors (i.e., not only a nonuniform BD). Our study provides a framework for consistently constraining their contribution, e.g., HD 189733b’s shape. (3) In addition, and related to the second point, we do not constrain a priori the system parameters to the bestfit of a conventional analysis nor the orbital eccentricity to zero; instead, we estimate the system parameters simultaneously with the BD. (4) We also investigate the modeldependence of our inferences, in contrast to M12 who focused on a dipolar brightness model (see related discussion in Sect. 5.2.1). In particular, we show their impact on the systemparameter PPD (see Sect. 5). Finally, M12 proposed an additional method (the slice method) to generate an exoplanet map from its occultation ingress/egress only; we discuss the relevance of similar methods in Sect. 6.3.
As a consequence, M12 estimate the brightness peak localization with narrow error bars (21.8 ± 1.5° east and 3.1 ± 9.4° away from the equator) for a circularized HD 189733b’s orbit; while we show that the brightness peak localization, as well as the systemparameter PPD, are modeldependent because of the ebρ_{⋆}BD correlation. In particular, we show that the light curve of an exoplanet does not constrain uniquely its brightness peak localization. Nevertheless, for a direct comparison to M12’s estimate of the brightness peak localization for a dipolar brightness model, we estimate it to 11.5 ± 4.3° east and 3.1 ± 11.4° away from the equator, see Sect. 6, Fig. 9b.
Fig. 13 Incidence of HD 189733’s variability on individual transitdepth estimates in the Spitzer/IRAC 8μm channel. Transitdepth estimated for six transits in individual (markers) and global (solid line, the dashed lines refer to the 1σ error bar) analyses. Our estimates (green) show no significant transitdepth variation with the relative eclipse phase (blue numbers) in opposition to Agol et al. (2010) (red). Similarly, we observe no significant variation of the transit parameters from one individual analysis to another. In addition, we observe no pattern specific to the occultation of a star spot (i.e., similar to, e.g., “Features A and B” in Fig. 1 of Pont et al. 2007). Therefore, we consider that our timeaveraged inferences (see Sect. 5) are not biased by HD 189733’s activity. 

Open with DEXTER 
6.7. Complementary analysis
We emphasize in our study the necessity of global approaches for fitting consistently the available data, in particular by highlighting the ebρ_{⋆}BD correlation. Regarding this correlation, we discuss in this subsection the effect of the stellar activity (through the planetary transit analysis) and of radial velocity (RV) measurements (i.e., for constraining ).
Fig. 14 Incidence of HD 189733’s RV measurements on inferences obtained from the Spitzer/IRAC 8μm photometry. a) Overall Keplerian fit (green) of the outoftransit RV data (black dots) with their 1σ error bars (red triangles) from Winn et al. (2006) and Boisse et al. (2009). b) Marginal PPD (68% and 95%confidence intervals) of and obtained solely from the RV data. This shows that the RV data do not constrain further HD 189733b’s eccentricity than the Spitzer/IRAC 8μm photometry for low complexity brightness models (see Figs. 4a, 6a and 10a). However, for more complex models that favor a localized hot spot and larger , HD 189733’s RV measurements may affect our inferences by rejecting the solutions involving ≳ 0.15 – especially in the context of the ebρ_{⋆}BD correlation (see Sect. 5.2.1). 

Open with DEXTER 
6.7.1. Effect of HD 189733’s activity
HD 189733 presents highactivity levels that may affect the transit parameters – incl., (R_{p}/R_{⋆})^{2}, due to occulted/unocculted star spots (Pont et al. 2007; Sing et al. 2011). Therefore, treating coherently spots of active stars is necessary in the present context – global approach and ebρ_{⋆}BD correlation. However, while important in the optical, the stellar activity may be negligible at 8 μm, regarding the data quality. We assess this statement performing individual analysis of the 6 transits used in this study (see Table 1). We present in Fig. 13 our individual transitdepth estimates. These estimates show no significant temporal variation (in opposition to Agol et al. 2010, who attributed these to stellar activity). Similarly, we observe no significant variation of the transit parameters from one individual analysis to another. In addition, we observe no pattern specific to the occultation of a star spot (i.e., similar to, e.g., “Features A and B” in the Fig. 1 of Pont et al. 2007). Therefore, we consider that our timeaveraged inferences (see Sect. 5) are not biased by HD 189733’s activity.
6.7.2. Complementary input of HD 189733’s RV measurements
RV measurements may help constraining and, therefore, providing a complementary insight into the ebρ_{⋆}BD correlation. Therefore, we analysed HD 189733’s outoftransit RV data (Boisse et al. 2009; Winn et al. 2006) while using priors on independent system parameters; (R_{p}/R_{⋆})^{2}, P, W, and b. This enable us to assess the constraint derived solely from the RV data on the eccentricity. We present in Fig. 14 our overall Keplerian fit and the marginal PPD of the parameters and . This shows that the RV data does not constrain HD 189733b’s eccentricity further than the Spitzer/IRAC 8μm photometry for low complexity brightness models (see Figs. 4a, 6a and 10a). However HD 189733’s RV measurements may affect our inferences for more complex models that favor a localized hot spot and larger eccentricity (see Figs. 10b and 6b), by rejecting the solutions involving ≳ 0.15. Therefore, we implement the analysis of RV data in our global method.
Fig. 15 Incidence of the RV measurements on the systemparameter PPD estimated using complex brightness models that suggest a nonzero eccentricity from the Spitzer/IRAC 8μm photometry (see Figs. 10b and 6b). a), b) Marginal PPDs (68% and 95%confidence intervals) of and for the Γ_{SH,3} and the Γ_{2} brightness models, respectively. Conceptually these correspond to the marginalized product of the PPDs estimated using solely the photometry (see Figs. 10b and 6b) and using solely the RV data (Fig. 14a). In particular, these highlight that the solutions involving ≳ 0.15 are strongly rejected by HD 189733’s RV data, leading to a redistribution of the probability density. c), d) Marginal PPDs of ρ_{⋆} and for the Γ_{SH,3} and Γ_{2} brightness models, respectively. These show the impact of RV data in the context of the ebρ_{⋆}BD correlation (see Sect. 5.2.1); by rejecting the solutions involving ≳ 0.15, HD 189733’s RV data favors solutions consistent with the inferences obtained with less complex brightness models (see Fig. 6). This emphasizes the necessity of global approaches to consistently probe the highlycorrelated parameter space of exoplanetary data. 

Open with DEXTER 
We observe no change of the system parameters PPDs for simple brightness models (see Table 3), as expected from the comparison of Fig. 14b and Figs. 6a and 10a. However, we observe the incidence of rejecting the solutions involving – in the context of the ebρ_{⋆}BD correlation − for complex brightness models (Γ_{SH,3} and Γ_{2}). In particular, the evolution of the parameter estimates – toward lower densities, larger impact parameter and a more localized and latitudinallyshifted hot spot – providing larger eccentricity are attenuated. We present the incidence on the systemparameter PPD in Fig. 15. The PPDs (Figs. 15a and b) conceptually correspond to the marginalized product of the PPDs estimated using solely the photometry (see Figs. 10b and 6b) and using solely the RV data (Fig. 14a). These highlight the redistribution of the probability density that results from the rejection of solutions involving ≳ 0.15 by HD 189733’s RV data. As a consequence of this probability redistribution in , the set of acceptable combinations to compensate HD 189733b’s anomalous occultation is affected too. In other words, the b, ρ_{⋆} and BD estimates are affected. On the one hand, the system parameter estimates are consistent with those obtained for less complex brightness models (see Table 3 and compare Figs. 15c and d to 6c and 10c). On the other hand, HD 189733b’s dayside brightness estimates present less confined patterns (compare Figs. 11b and 8 with 16a and b, respectively). This emphasizes the necessity of global approaches to consistently probe the highlycorrelated parameter space of exoplanetary data.
Finally, we present our improved constraint on the upper limit of HD 189733b’s orbital eccentricity, e ≤ 0.011 (95% confidence), based on our global analysis of the Spitzer/IRAC 8μm photometry and the RV measurements of HD 189733.
Parameter estimates obtained from the photometry and the RV measurements.
Fig. 16 Incidence of HD 189733’s RV measurements on HD 189733b’s dayside brightness distribution, estimated using complex brightness models that suggest nonzero eccentricity from the Spitzer/IRAC 8μm photometry (see Figs. 10b and 6b). Left: relative brightness distribution at HD 189733b’s dayside. Right: dayside standard deviation. a) Estimate using the Γ_{SH,3} brightness model. b) Estimate using the Γ_{2} brightness model. The incidence of HD 189733’s RV data is an extension of the brightness patterns (compare with Figs. 11b and 8, respectively) associated with a decrease of the brightness peak – conservation of the hemisphere integrated flux. These show that localized brightness pattern that are favored by the photometry are rejected by the RV data in the context of the ebρ_{⋆}BD correlation; because these are associated with . 

Open with DEXTER 
7. Summary
We have performed two analyses of the HD 189733b public eclipse data obtained at 8 μm with Spitzer/IRAC. The first, conventional analysis provides HD 189733’s system parameters and the corrected and phasefolded light curves; which deviate from the occultation of a uniformlybright disk at the 6σ level. The second analysis investigates the possible contributing factors of this deviation following a new, general methodology.
We demonstrate that this deviation emerges mainly from a largescale brightness structure in HD 18733b’s atmosphere; we reject HD 189733b’s shape as a possible contributing factor based on the transit residuals. It indicates a hot spot within the atmospheric layers probed at 8 μm as well as a hemispheric asymmetry. In addition, we show a correlation between the system parameters and HD 189733b’s BD – what we term the “ebρ_{⋆}BD correlation”. This enables us to outline that for data of sufficient quality the conventional assumption of a uniformlybright exoplanet leads to an underestimated uncertainty on (and possibly biased estimates of) system parameters, in particular on ρ_{⋆} and ρ_{p}. Notably, we find that the more complex HD 189733b’s brightness model, the larger the eccentricity, the lower the densities, the larger the impact parameter and the more localized and latitudinallyshifted the hot spot estimated. In this context, we show that the light curve of an exoplanet does not constrain uniquely its brightness peak localization. Therefore, we propose to use the dayside barycenter – when a unidimensional parameter is necessary – because it contains partial twodimensional information from weighting the BD according to the geometrical configuration at superior conjunction.
Because of the highlighted ebρ_{⋆}BD correlation, we assess the incidence of HD 189733’s RV measurements on our inferences obtained solely from Spitzer/IRAC 8μm photometry. We observe redistributions of the probability density for the PPDs obtained with complex brightness models, which result from the rejection by HD 189733’s RV data of solutions involving ≳ 0.15, favored by the photometry in case of a localized hot spot – i.e., for complex brightness models. This emphasizes the necessity of global approaches to consistently probe the highlycorrelated parameter space of exoplanetary data.
As a final result, we show that, for present data quality, the system parameter estimates for complex brightness model are similar to those obtained for a dipolar brightness models. In particular, we present our constraint on the upper limit of HD 189733b’s orbital eccentricity, e ≤ 0.011 (95% confidence).
We discuss the perspectives of applying our methods to observations in different spectral bands than in the IRAC 8μm channel. In particular, observations of HD 189733b at other wavelengths could improve the constraints on the system parameters while ultimately yielding a largescale 3D map – because different optical depths would be probed under the same orbital configuration. In addition, applications to observations in the visible (ground and spacebased, e.g., Kepler data) could also be valuable for targets nonuniformly bright because of albedo (e.g., cloud/hazes or Rayleigh scattering). We emphasize that the albedo origins could make the exoplanet nonuniformly bright and modify its occultation ingress/egress, similarly to the effect of an hot spot. Therefore, our study could constrain independently the geometrical albedo by scanning an exoplanet dayside, what would be a complementary asset to Madhusudhan & Burrows (2012)’s framework. We discuss the relevance of mapping methods based on direct reconstruction from an exoplanet “slicing”. We emphasize the necessity of independent constraints on the host star to consistently constrain the exoplanet shape, in particular on its limbdarkening coefficient – although we emphasize that the phase curve also constrains the exoplanet shape. Finally, we suggest that the future exoplanetaryatmosphere investigations of spacedbased facilities like JWST and EChO could ultimately lead to timedependent 3D maps of distant worlds.
For an uptodate list see the Extrasolar Planets Encyclopaedia: http://exoplanet.eu
The uniform time offset measures the time lag between the observed secondary eclipse and that predicted by a planet with spatially uniform emission (defined by Williams et al. 2006).
Data available in the form of Basic Calibrated Data (BCD) on the Infrared Science Archive: http://sha.ipac.caltech.edu/applications/Spitzer/SHA//
Acknowledgments
This study is based on observations made with the Spitzer Space Telescope, which is operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract to NASA. We thank the people who made and keep the Spitzer Space Telescope such a wonderful scientific adventure. We thank H. Knutson for sharing HD 189733b’s phase curve data. We are grateful to A. Zsom, V. Stamenkovic, A. Triaud and E. Agol for helpful discussions and to the anonymous referee for her/his comments that strengthen this manuscript. J.d.W. acknowledges the hospitality of the Institut d’Astrophysique et de Géophysique, Université de Liège where portions of this study were completed and thanks A. Lanotte and S. Sohy for their help. J.d.W. acknowledges support from the Grayce B. Kerr Foundation and the WBI (WallonieBruxelles International) in the form of fellowships. Finally, J.d.W. also acknowledges support provided by the Odissea Prize initiated by the Belgian Senate.
References
 Agol, E., Cowan, N. B., Knutson, H. A., et al. 2010, ApJ, 721, 1861 [NASA ADS] [CrossRef] [Google Scholar]
 Barnes, J. W., Cooper, C. S., Showman, A. P., & Hubbard, W. B. 2009, ApJ, 706, 877 [NASA ADS] [CrossRef] [Google Scholar]
 Batygin, K., Stevenson, D. J., & Bodenheimer, P. H. 2011, ApJ, 738, 1 [NASA ADS] [CrossRef] [Google Scholar]
 Berdyugina, S. V., Berdyugin, A. V., Fluri, D. M., & Piirola, V. 2011, ApJ, 728, L6 [NASA ADS] [CrossRef] [Google Scholar]
 Bigg, E. K. 1964, Nature, 203, 1008 [NASA ADS] [CrossRef] [Google Scholar]
 Boisse, I., Moutou, C., VidalMadjar, A., et al. 2009, A&A, 495, 959 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Bouchy, F., Udry, S., Mayor, M., et al. 2005, A&A, 444, L15 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Burnham, K., & Anderson, D. 2002, Model Selection and Multimodel Inference: A Practical Informationtheoretic Approach (Springer) [Google Scholar]
 Carter, J. A., & Winn, J. N. 2010, ApJ, 709, 1219 [NASA ADS] [CrossRef] [Google Scholar]
 Charbonneau, D., Knutson, H. A., Barman, T., et al. 2008, ApJ, 686, 1341 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
 Clampin, M. 2010, in Proceedings of the conference In the Spirit of Lyot 2010 [Google Scholar]
 Claret, A., & Bloemen, S. 2011, A&A, 529, A75 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Cooper, C. S., & Showman, A. P. 2005, ApJ, 629, L45 [NASA ADS] [CrossRef] [Google Scholar]
 Deming, D., Harrington, J., Seager, S., & Richardson, L. J. 2006, ApJ, 644, 560 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
 Demory, B.O., Seager, S., Madhusudhan, N., et al. 2011, ApJ, 735, L12 [NASA ADS] [CrossRef] [Google Scholar]
 Deroo, P., Swain, M. R., & Vasisht, G. 2010, MNRAS, submitted [arXiv:1011.0476] [Google Scholar]
 Désert, J.M., Lecavelier des Etangs, A., Hébrard, G., et al. 2009, ApJ, 699, 478 [NASA ADS] [CrossRef] [Google Scholar]
 DobbsDixon, I., Cumming, A., & Lin, D. N. C. 2010, ApJ, 710, 1395 [NASA ADS] [CrossRef] [Google Scholar]
 Fabrycky, D. C. 2010, NonKeplerian Dynamics of Exoplanets, ed. S. Seager, 217 [Google Scholar]
 Faigler, S., & Mazeh, T. 2011, MNRAS, 415, 3921 [NASA ADS] [CrossRef] [Google Scholar]
 Fazio, G. G., Hora, J. L., Allen, L. E., et al. 2004, ApJS, 154, 10 [NASA ADS] [CrossRef] [Google Scholar]
 Ford, E. B. 2006, ApJ, 642, 505 [NASA ADS] [CrossRef] [Google Scholar]
 Gelman, A., Carlin, J., Stern, H., & Rubin, D. 2004, Bayesian Data Analysis (Chapman & Hall/CRC) [Google Scholar]
 Gibson, N. P., Pont, F., & Aigrain, S. 2011, MNRAS, 411, 2199 [NASA ADS] [CrossRef] [Google Scholar]
 Gillon, M., Demory, B.O., Triaud, A. H. M. J., et al. 2009, A&A, 506, 359 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Gillon, M., Hatzes, A., Csizmadia, S., et al. 2010a, A&A, 520, A97 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Gillon, M., Lanotte, A. A., Barman, T., et al. 2010b, A&A, 511, A3 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Gregory, P. C. 2005, Bayesian Logical Data Analysis for the Physical Sciences: A Comparative Approach with Mathematica Support (Cambridge University Press) [Google Scholar]
 Grillmair, C. J., Charbonneau, D., Burrows, A., et al. 2007, ApJ, 658, L115 [NASA ADS] [CrossRef] [Google Scholar]
 Heng, K. 2012, ApJ, 748, L17 [NASA ADS] [CrossRef] [Google Scholar]
 Horne, K. 1985, MNRAS, 213, 129 [NASA ADS] [Google Scholar]
 Huitson, C. M., Sing, D. K., VidalMadjar, A., et al. 2012, MNRAS, 422, 2477 [NASA ADS] [CrossRef] [Google Scholar]
 Ip, W.H., Kopp, A., & Hu, J.H. 2004, ApJ, 602, L53 [NASA ADS] [CrossRef] [Google Scholar]
 Kipping, D., & Bakos, G. 2011, ApJ, 730, 50 [NASA ADS] [CrossRef] [Google Scholar]
 Knutson, H. A., Charbonneau, D., Allen, L. E., et al. 2007, Nature, 447, 183 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
 Knutson, H. A., Charbonneau, D., Cowan, N. B., et al. 2009, ApJ, 690, 822 [NASA ADS] [CrossRef] [Google Scholar]
 Knutson, H. A., Lewis, N., Fortney, J. J., et al. 2012, ApJ, 754, 22 [NASA ADS] [CrossRef] [Google Scholar]
 Kopal, Z. 1959, Close binary systems [Google Scholar]
 Laine, R. O., Lin, D. N. C., & Dong, S. 2008, ApJ, 685, 521 [NASA ADS] [CrossRef] [Google Scholar]
 Lanza, A. F. 2009, A&A, 505, 339 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Lecavelier des Etangs, A., Bourrier, V., Wheatley, P. J., et al. 2012, A&A, 543, L4 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Madhusudhan, N., & Burrows, A. 2012, ApJ, 747, 25 [NASA ADS] [CrossRef] [Google Scholar]
 Madhusudhan, N., & Seager, S. 2009, ApJ, 707, 24 [NASA ADS] [CrossRef] [Google Scholar]
 Majeau, C., Agol, E., & Cowan, N. B. 2012, ApJ, 747, L20 [NASA ADS] [CrossRef] [Google Scholar]
 Mandel, K., & Agol, E. 2002, ApJ, 580, L171 [NASA ADS] [CrossRef] [Google Scholar]
 Menou, K. 2012, ApJ, 745, 138 [NASA ADS] [CrossRef] [Google Scholar]
 Ollivier, M., Gillon, M., Santerne, A., et al. 2012, A&A, 541, A149 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Pont, F., Zucker, S., & Queloz, D. 2006, MNRAS, 373, 231 [NASA ADS] [CrossRef] [Google Scholar]
 Pont, F., Gilliland, R. L., Moutou, C., et al. 2007, A&A, 476, 1347 [NASA ADS] [CrossRef] [EDP Sciences] [MathSciNet] [Google Scholar]
 Press, W. H., Flannery, B. P., Teukolsky, S. A., & Vetterling, W. T. 1992, Numerical Recipes in FORTRAN: The Art of Scientific Computing, 2nd edn. (Cambridge University Press), 51 [Google Scholar]
 Rauscher, E., & Menou, K. 2010, ApJ, 714, 1334 [NASA ADS] [CrossRef] [Google Scholar]
 Rauscher, E., & Menou, K. 2012, ApJ, 750, 96 [NASA ADS] [CrossRef] [Google Scholar]
 Rauscher, E., Menou, K., Seager, S., et al. 2007, ApJ, 664, 1199 [NASA ADS] [CrossRef] [Google Scholar]
 Redfield, S., Endl, M., Cochran, W. D., & Koesterke, L. 2008, ApJ, 673, L87 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
 Russell, H. N., & Merrill, J. E. 1952, The determination of the elements of eclipsing binaries [Google Scholar]
 Seager, S., & Deming, D. 2010, ARA&A, 48, 631 [NASA ADS] [CrossRef] [Google Scholar]
 Seager, S., & MallénOrnelas, G. 2003, ApJ, 585, 1038 [NASA ADS] [CrossRef] [Google Scholar]
 Shkolnik, E., Walker, G. A. H., Bohlender, D. A., Gu, P.G., & Kürster, M. 2005, ApJ, 622, 1075 [NASA ADS] [CrossRef] [Google Scholar]
 Showman, A. P., & Guillot, T. 2002, A&A, 385, 166 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Showman, A. P., Fortney, J. J., Lian, Y., et al. 2009, ApJ, 699, 564 [NASA ADS] [CrossRef] [Google Scholar]
 Sing, D. K., Pont, F., Aigrain, S., et al. 2011, MNRAS, 416, 1443 [NASA ADS] [CrossRef] [Google Scholar]
 Southworth, J. 2010, MNRAS, 408, 1689 [NASA ADS] [CrossRef] [Google Scholar]
 Stetson, P. B. 1987, PASP, 99, 191 [NASA ADS] [CrossRef] [Google Scholar]
 Swain, M. R., Vasisht, G., & Tinetti, G. 2008, Nature, 452, 329 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
 Tinetti, G., VidalMadjar, A., Liang, M.C., et al. 2007, Nature, 448, 169 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
 Tinetti, G., Beaulieu, J. P., Henning, T., et al. 2012, Exper. Astron., 34, 311 [NASA ADS] [CrossRef] [Google Scholar]
 Triaud, A. H. M. J., Queloz, D., Bouchy, F., et al. 2009, A&A, 506, 377 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Warner, B., Robinson, E. L., & Nather, R. E. 1971, MNRAS, 154, 455 [NASA ADS] [Google Scholar]
 Werner, M. W., Roellig, T. L., Low, F. J., et al. 2004, ApJS, 154, 1 [NASA ADS] [CrossRef] [Google Scholar]
 Williams, P. K. G., Charbonneau, D., Cooper, C. S., Showman, A. P., & Fortney, J. J. 2006, ApJ, 649, 1020 [NASA ADS] [CrossRef] [Google Scholar]
 Winn, J. N. 2010, Exoplanet Transits and Occultations, ed. S. Seager, 55 [Google Scholar]
 Winn, J. N., Johnson, J. A., Marcy, G. W., et al. 2006, ApJ, 653, L69 [NASA ADS] [CrossRef] [Google Scholar]
 Winn, J. N., Holman, M. J., Henry, G. W., et al. 2007, AJ, 133, 1828 [NASA ADS] [CrossRef] [Google Scholar]
 Winn, J. N., Holman, M. J., Torres, G., et al. 2008, ApJ, 683, 1076 [NASA ADS] [CrossRef] [Google Scholar]
 Yang, Y. 2005, Biometrika, 92, 937 [CrossRef] [Google Scholar]
Appendix A: Movies
Fig. A.1. (Movies) Insight into the brightness distribution estimates. Animations showing compilations of HD 189733b’s dayside brightness distributions accepted along the MCMC simulations for the Γ_{SH,1} and Γ_{2} models, a) and b) respectively. These animations show that (1) the amplitude and (2) the longitudinal localization for the Γ_{SH,1} brightness model is more constrained than for the Γ_{2} model (by the occultation depth and by the phase curve, respectively) because of its fixed and large structure. However, (3) the Γ_{SH,1} model is less constrained in latitude (by the secondaryeclipse scanning) than the more complex Γ_{2} model which enable more confined structures that induce larger deviation in occultation ingress/egress (schematic description in Fig. 3).
Online material
Movie 1a (Access here)
Movie 1b (Access here)
All Tables
Fit properties for different fitting models of HD 189733’s photometry in the Spitzer/IRAC 8μm channel.
All Figures
Fig. 1 Schematic description of the different scanning processes observable for an occulted exoplanet. The green dotted lines indicate the scanning processes during the exoplanet occultation ingress/egress. The red dashed line indicates the scanning process that results from the exoplanet rotation and produces the exoplanet phase curve – this scanning appears longitudinal for an observer as long as the exoplanet spin is close to the projection plane, e.g., for a transiting and synchronized exoplanet. The component labeled “combined” shows the specific grid generated by these three scanning processes. 

Open with DEXTER  
In the text 
Fig. 2 IRAC 8μm HD 189733b’s transit and occultation ingress/egress photometry binned per 1 min and corrected for the systematics (black dots) with their 1σ error bars (red triangles) and the bestfitting eclipse model superimposed (in red). The green dots present the residuals from Agol et al. (2010, their Fig. 12), obtained using an average of the bestfit models from 7 individual eclipse analyses. Left: phasefolded transits show no significant deviation to the transit of a disk during ingress/egress. Right: phasefolded occultation ingress/egress deviate from the eclipse of a uniformlybright disk, highlighting the secondaryeclipse scanning of HD 189733b’s dayside. Top: phasefolded and corrected eclipse photometry. Bottom: phasefolded residuals. 

Open with DEXTER  
In the text 
Fig. 3 Schematic description of the anomalous occultation ingress/egress induced by the shape or the brightness distribution of an exoplanet. The red curve indicates the occultation photometry for a nonuniformlybright disk (hot spot in red). The yellow curve indicates the occultation photometry for an oblate exoplanet (yellow ellipse). Both synthetic scenarios show specific deviations from the occultation photometry of uniformlybright disk (black curve) in the occultation ingress/egress. 

Open with DEXTER  
In the text 
Fig. 4 Incidence of assuming an exoplanet to be uniformly bright and its orbit to be circularized. a) Marginal posterior probability distribution (PPD, 68% and 95%confidence intervals) of and that shows an unusual correlation. b) Marginal PPD (95%confidence intervals) of ρ_{⋆} and b that highlights the increase of adequate solutions enabled by the additional dimensions of the parameter space probed when relaxing the circularized orbit assumption. c) Deviations in occultation ingress/egress from the medianfit model for the individual perturbations of b, , and ρ_{⋆} by their estimated uncertainty (see Table 2, Col. 2). It outlines that the system parameters enable compensation of an anomalous occultation that emerges from, e.g., a nonuniformlybright exoplanet (see Fig. 3) leading to biased estimates of the system parameters. 

Open with DEXTER  
In the text 
Fig. 5 Phasefolded IRAC 8μm HD 189733b’s photometry binned per 1 min and corrected for the systematics (black dots) with their 1σ error bars (red triangles) and the bestfitting eclipse models superimposed. Left: phasefolded transits. Right: phase curve and phasefolded occultations that show the benefit of using nonuniform brightness model. 

Open with DEXTER  
In the text 
Fig. 6 Incidence of the brightness model complexity on the systemparameter posterior probability distribution (PPD), using unipolar models. a), b) Marginal PPDs (68% and 95%confidence intervals) of and for two unipolar brightness models, respectively, with a fixed and large structure (Γ_{SH,1}) and with freeconfinement structure (Γ_{2}). A comparison with Fig. 4a shows that is constrained closer to zero. The reason is that nonuniform brightness models enable the exploration of additional dimensions of the parameter space and, therefore, provide additional adequate combinations of the contributing factors to compensate an anomalous occultation. In particular, the uniform time offset is now mainly compensated by a nonuniform BD, rather than by as in conventional analysis. It shows also the evolution of the marginal PPD toward larger when using more complex brightness models – which enable more localized brightness structure. c), d) Marginal PPDs of ρ_{⋆} and for, respectively, Γ_{SH,1} and Γ_{2}. These show the impact of the brightness distribution on the retrieved system parameters (i.e., ebρ_{⋆}BD correlation); in particular, it outlines the possible overestimation of ρ_{⋆} by 5% (i.e., at 6σ of the conventional estimate) when using extended brightness models, e.g., a uniform brightness distribution. 

Open with DEXTER  
In the text 
Fig. 7 Estimate of HD 189733b’s global brightness distribution in the IRAC 8μm channel using the Γ_{SH,1} brightness model. Left: relative brightness distribution at HD 189733b’s dayside. Right: dayside standard deviation. Because of its fixed and large structure, the Γ_{SH,1} brightness model is wellconstrained in amplitude (by the occultation depth) and in longitudinal localization (by the phase curve). However, it is less constrained in latitude (by the secondaryeclipse scanning) than more confined model (schematic description in Fig. 3), e.g., Γ_{2} (see Fig. 8). These modelinduced constraints are observable on the dayside standard deviation; which is significantly lower than for more complex brightness models (see Figs. 8, 11a and b) and is related to its gradient, because the BDs accepted along the MCMC simulations differ from each other mainly in latitudinal orientation (see Fig. A.1a). 

Open with DEXTER  
In the text 
Fig. 8 Estimate of HD 189733b’s global brightness distribution in the IRAC 8μm channel using the Γ_{2} brightness model. Left: relative brightness distribution at HD 189733b’s dayside. Right: dayside standard deviation. Because of its increased complexity, the Γ_{2} brightness model enable more localized structure that are less constrained in amplitude (by the occultation depth) and in longitudinal localization (by the phase curve) than largeandfixedstructure model, e.g., the Γ_{SH,1} model (see Fig. 7). However, it is wellconstrained in latitude by the secondaryeclipse scanning that is sensitive to confined brightness structure. These modelinduced constraints are observable on the dayside standard deviation that shows a maximum at the brightness peak localization and extended wings towards west and east along the planetary equator, because the BDs accepted along the MCMC simulations mainly affect the former by their amplitude change and the latter by their structure change (see Fig. A.1b). 

Open with DEXTER  
In the text 
Fig. 9 Dependence and significance of the brightness peak localization. a) Marginal PPD (68% and 95%confidence intervals) of the brightness peak longitude, Δφ, and ecosω for the Γ_{SH,1} brightness model – the simplest nonuniform brightness model used in this study. This shows the correlation between the brightness model and the orbital eccentricity. This correlation emerges from enabling compensation of the anomalous occultation with a larger set of contributing factors, i.e., by including nonuniform brightness models. In particular, the uniform time offset is now mainly compensated by an nonuniformlybright model, rather than by ecosω as in conventional analysis (see Fig. 3, Eq. (6) and Fig. 4c). b) Marginal PPDs (68%confidence intervals) of the brightness peak localization for the Γ_{SH,1} and Γ_{2} brightness models. It shows that the brightness peak localization is modeldependent. For example, the longitudinal Γ_{SH,1} peak localization is constrained by the phase curve because of its large and constant extension; while the free extension of the Γ_{2} model relaxes this longitudinal constraint (see Sect. 5). Therefore, the light curve of an exoplanet does not constrain uniquely its brightness peak localization; furthermore, the brightness peak localization is not an adequate parameter to characterize complex exoplanet brightness distributions. 

Open with DEXTER  
In the text 
Fig. 10 Incidence of the brightness model extension on the system parameter posterior probability distribution (PPD), using multipolar models. a), b) Marginal PPDs (68% and 95%confidence intervals) of and for the quadrupolar and octupolar brightness models, respectively. c), d) Marginal PPDs of ρ_{⋆} and for the quadripolar and octopolar brightness models, respectively. These confirm the trend observed in Sect. 5.1 (see Fig. 6) towards larger and lower stellar density when increasing the complexity of the brightness distribution, i.e., when enabling more localized structures. 

Open with DEXTER  
In the text 
Fig. 11 Estimate of HD 189733b’s global brightness distribution in the IRAC 8μm channel using multipolar brightness models. Left: relative brightness distribution at HD 189733b’s dayside. Right: dayside standard deviation. a) Estimate using the Γ_{SH,2} brightness model. b) Estimate using the Γ_{SH,3} brightness model. These confirm the trend toward a more localized and latitudinallyshifted hot spot when increasing the brightnessmodel complexity. Therefore, together with Fig. 10 these outline that the more complex HD 189733b’s brightness model, the larger the eccentricity, the lower the densities, the larger the impact parameter and the more localized and latitudinallyshifted the hot spot estimated. 

Open with DEXTER  
In the text 
Fig. 12 Reducing brightness distributions to unidimensional parameters. Marginal PPDs (68%confidence intervals) of the dayside barycenter brightness peak localization for the Γ_{SH,1} and Γ_{2} brightness models. A comparison with the marginal PPDs of the brightness peak localization (Fig. 9b) shows the reduced modeldependence of the dayside barycenter. In particular, it shows a lessextended PPD for the Γ_{2}model barycenter; because this extension for the brightnesspeaklocalization PPD emerges from the model wings, weighted by the dayside barycenter. In addition, it shows the shift and slight shrinking of the Γ_{SH,1} PPD that reflect the barycenter weighting based on to the geometrical configuration at superior conjunction; map cells closer to the substellar point have more weight. 

Open with DEXTER  
In the text 
Fig. 13 Incidence of HD 189733’s variability on individual transitdepth estimates in the Spitzer/IRAC 8μm channel. Transitdepth estimated for six transits in individual (markers) and global (solid line, the dashed lines refer to the 1σ error bar) analyses. Our estimates (green) show no significant transitdepth variation with the relative eclipse phase (blue numbers) in opposition to Agol et al. (2010) (red). Similarly, we observe no significant variation of the transit parameters from one individual analysis to another. In addition, we observe no pattern specific to the occultation of a star spot (i.e., similar to, e.g., “Features A and B” in Fig. 1 of Pont et al. 2007). Therefore, we consider that our timeaveraged inferences (see Sect. 5) are not biased by HD 189733’s activity. 

Open with DEXTER  
In the text 
Fig. 14 Incidence of HD 189733’s RV measurements on inferences obtained from the Spitzer/IRAC 8μm photometry. a) Overall Keplerian fit (green) of the outoftransit RV data (black dots) with their 1σ error bars (red triangles) from Winn et al. (2006) and Boisse et al. (2009). b) Marginal PPD (68% and 95%confidence intervals) of and obtained solely from the RV data. This shows that the RV data do not constrain further HD 189733b’s eccentricity than the Spitzer/IRAC 8μm photometry for low complexity brightness models (see Figs. 4a, 6a and 10a). However, for more complex models that favor a localized hot spot and larger , HD 189733’s RV measurements may affect our inferences by rejecting the solutions involving ≳ 0.15 – especially in the context of the ebρ_{⋆}BD correlation (see Sect. 5.2.1). 

Open with DEXTER  
In the text 
Fig. 15 Incidence of the RV measurements on the systemparameter PPD estimated using complex brightness models that suggest a nonzero eccentricity from the Spitzer/IRAC 8μm photometry (see Figs. 10b and 6b). a), b) Marginal PPDs (68% and 95%confidence intervals) of and for the Γ_{SH,3} and the Γ_{2} brightness models, respectively. Conceptually these correspond to the marginalized product of the PPDs estimated using solely the photometry (see Figs. 10b and 6b) and using solely the RV data (Fig. 14a). In particular, these highlight that the solutions involving ≳ 0.15 are strongly rejected by HD 189733’s RV data, leading to a redistribution of the probability density. c), d) Marginal PPDs of ρ_{⋆} and for the Γ_{SH,3} and Γ_{2} brightness models, respectively. These show the impact of RV data in the context of the ebρ_{⋆}BD correlation (see Sect. 5.2.1); by rejecting the solutions involving ≳ 0.15, HD 189733’s RV data favors solutions consistent with the inferences obtained with less complex brightness models (see Fig. 6). This emphasizes the necessity of global approaches to consistently probe the highlycorrelated parameter space of exoplanetary data. 

Open with DEXTER  
In the text 
Fig. 16 Incidence of HD 189733’s RV measurements on HD 189733b’s dayside brightness distribution, estimated using complex brightness models that suggest nonzero eccentricity from the Spitzer/IRAC 8μm photometry (see Figs. 10b and 6b). Left: relative brightness distribution at HD 189733b’s dayside. Right: dayside standard deviation. a) Estimate using the Γ_{SH,3} brightness model. b) Estimate using the Γ_{2} brightness model. The incidence of HD 189733’s RV data is an extension of the brightness patterns (compare with Figs. 11b and 8, respectively) associated with a decrease of the brightness peak – conservation of the hemisphere integrated flux. These show that localized brightness pattern that are favored by the photometry are rejected by the RV data in the context of the ebρ_{⋆}BD correlation; because these are associated with . 

Open with DEXTER  
In the text 