PHOTOMETRIC FOLLOW-UP OBSERVATIONS OF THE TRANSITING NEPTUNE-MASS PLANET GJ 436b

, , , , , , and

Published 2009 March 25 © 2009. The American Astronomical Society. All rights reserved.
, , Citation Avi Shporer et al 2009 ApJ 694 1559 DOI 10.1088/0004-637X/694/2/1559

0004-637X/694/2/1559

ABSTRACT

This paper presents multiband photometric follow-up observations of the Neptune-mass transiting planet GJ 436b, consisting of five new ground-based transit light curves obtained in 2007 May. Together with one already published light curve, we have at hand a total of six light curves, spanning 29 days. The analysis of the data yields an orbital period P = 2.64386 ± 0.00003 days, midtransit time Tc [HJD] = 2454235.8355 ± 0.0001, planet mass Mp = 23.1 ± 0.9 M = 0.073 ± 0.003 MJup, planet radius Rp = 4.2 ± 0.2 R = 0.37 ± 0.01 RJup, and stellar radius Rs = 0.45 ± 0.02 R. Our typical precision for the midtransit timing for each transit is about 30 s. We searched the data for a possible signature of a second planet in the system through transit timing variations (TTV) and variation of the impact parameter. The analysis could not rule out a small, of the order of a minute, TTV and a long-term modulation of the impact parameter, of the order of +0.2 yr−1.

Export citation and abstract BibTeX RIS

1. INTRODUCTION

Of the almost 300 extrasolar planets discovered to date6, the ∼35 transiting planets are the only ones allowing a measurement of their mass and radius and the study of their atmospheres (e.g., Charbonneau et al. 2007; Fortney 2008; Guillot 2008; Pont et al. 2008). Of those, GJ 436b serves as a unique opportunity to study a planet with mass and radius as small as Neptune. Although several planets have already been detected with minimum masses similar to the mass of GJ 436b, for example GJ 581b, c, and d, (Bonfils et al. 2005; Udry et al. 2007), HD 4308b (Udry et al. 2006), and Cnc 55e (Fischer et al. 2008), GJ 436b is the only transiting Neptune-mass planet discovered so far. Moreover, it is the only known transiting planet orbiting an M-type star, presenting an opportunity to help constrain stellar models in this mass range.

GJ 436b was initially discovered by Butler et al. (2004) through a radial velocity (RV) modulation of its host star. Maness et al. (2007) have refined its orbital elements, specifically the eccentricity, e = 0.160 ± 0.019, and identified a linear velocity trend of 1.36 ± 0.4 m s−1 yr−1. Those authors suggested that the velocity trend might result from the presence of a long-period second planet in the system, a planet that could induce a periodic modulation of the orbital eccentricity of GJ 436b. This type of effect was suggested already for triple stellar systems (Mazeh & Shaham 1979, see also Mazeh 1990), and for a planet in a binary system for 16 Cygni B (Mazeh et al. 1997; Holman et al. 1997). Maness et al. (2007) pointed out that this interpretation of the eccentricity looks especially attractive if the circularization timescale of GJ 436b is shorter than the age of the system, because then we need to explain why the orbit of the planet has not been circularized.

Butler et al. (2004) obtained photometric observations of GJ 436 at the expected time of possible transits and concluded that a transit was unlikely. Nevertheless, the transiting nature of GJ 436b was recently discovered by Gillon et al. (2007a) who measured a planetary radius of ≃4 R and mass of 22.6 M. Soon after the discovery, a transit and secondary eclipse events were observed by Spitzer Space Telescope at 8 μm (Gillon et al. 2007b; Deming et al. 2007; Demory et al. 2007), further constraining the system parameters and reducing the uncertainties by a factor of about 10 for the planet-to-star radius ratio and midtransit timing. Those observations resulted in a somewhat increased planetary radius of 4.2 ± 0.2 R. Torres (2007) was able to refine further the planet parameters by deriving more accurate constraints on the host star and obtained a slightly larger planetary mass, 23.17 ± 0.79 M, and a more precise planetary radius of 4.22+0.09−0.10R.

Recently, Ribas et al. (2008) restudied the system and suggested that the observed radial velocities of the system are consistent with an additional small, super-Earth planet in a close orbit around GJ 436 in a 1:2 mean-motion resonance with the known planet. Similar to Maness et al. (2007), they suggested that this planet could pump eccentricity into the orbit. Ribas et al. (2008) further suggested that such an additional planet can induce a precession of the orbital plane of GJ 436b, if the orbital planes of motion of the two planets are inclined relative to each other. Such a precession should change the inclination relative to our line of sight, and therefore might explain why Butler et al. (2004) failed to observe a transit in 2004.

Immediately after the detection of the transits of GJ 436b by Gillon et al. (2007a) we launched a ground-based (GB) observational campaign to obtain high-quality transit light curves in different filters. Our motivation was to obtain the first group of light curves which will be used by future studies to look for a possible variation of the impact parameter and to search for transit timing variations (TTV; Agol et al. 2005; Holman & Murray 2005). Our observations are described in Section 2. In our analysis, described in Section 3, we simultaneously analyze six complete GB transit light curves: five light curves obtained in our campaign and the light curve from Gillon et al. (2007a). In Section 4, we discuss our results.

2. OBSERVATIONS

In addition to the transit of UT 2007 May 2, which was reported by Gillon et al. (2007a) and was the first complete transit observation of GJ 436b, we obtained five complete transit light curves of four different transit events. The following paragraphs briefly describe these observations.

We observed the transits of UT 2007 May 4 and May 25 with the Wise Observatory 1 m telescope, located in Israel. We used a Princeton Instruments (PI) VersArray camera and a 1340 × 1300 pixel back-illuminated CCD, with a pixel scale of 0farcs58 per pixel, giving a field of view of 13farcm0 × 12farcm6. On May 4 we used a Johnson V filter and on May 25 a Bessel R filter. The flexible scheduling of the Wise 1 m telescope was an important factor in obtaining the transit light curve of May 4, less than three days after the first observation of a complete transit.

The transit of UT 2007 May 4 was simultaneously observed with the Wise Observatory 0.46 m telescope operated remotely. The camera was a 2148 × 1472 pixel Santa Barbara Instrument Group (SBIG) ST-10 XME CCD detector with a field of view of 40farcm5 × 27farcm3 and a scale of 1farcs1 per pixel. This camera has no filters. For a detailed description of this telescope and instrument, see Brosch et al. (2008).

On UT 2007 May 23 and May 31 we used the 1.2 m telescope at the Fred L. Whipple Observatory (FLWO) on Mount Hopkins, Arizona. The camera was the KeplerCam, which is a 40962 Fairchild 486 CCD with a square field of view 23farcm1 on a side (Szentgyorgyi et al. 2005). On both nights, the filter was Sloan z and a 2 × 2 binning was used, giving a pixel scale of 0farcs68 per binned pixel.

Basic data reduction procedures, including bias and flat-field correction, were applied to all images using standard IRAF7 routines. We then performed aperture photometry of GJ 436 and several comparison stars of similar brightness showing no significant variability. The light curve of GJ 436 was calibrated by dividing it by the summed flux of the comparison stars. Next, we fitted a polynomial of degree one or two to the out-of-transit (OOT) points versus time, and divided all points by this polynomial. The amplitude of these polynomial corrections in mmag per transit duration was in the range of 0.1–1.0 mmag. As a final step, we divided the entire light curve by the OOT median intensity and normalized the measurement uncertainties so the median uncertainty OOT is equal to the OOT RMS. Table 1 lists all photometric measurements of the five light curves obtained in this campaign.

Table 1. Photometry of GJ 436

Observatory + Telescope Filter HJD Relative Flux Relative Flux Error
Wise 0.46 m clear 2454225.229864 1.0015 0.0020
Wise 1 m V 2454225.218116 0.9997 0.0010
FLWO 1.2 m z 2454243.658856 1.0002 0.0020
Wise 1 m R 2454246.309414 1.0024 0.0017
FLWO 1.2 m z 2454251.641074 0.9987 0.0024

Note. This table includes the five light curves obtained in this work.

Only a portion of this table is shown here to demonstrate its form and content. A machine-readable version of the full table is available.

Download table as:  DataTypeset image

Table 2 lists for each light curve the UT start date and time, observatory and telescope, the filter used, limb-darkening coefficients used in its analysis, average exposure time, average cadence (in minute−1), duration of the entire observation, start and end airmass, the β factor (see below), RMS residuals, and the Photometric Noise Rate (PNR). The PNR is a quantity which takes into account the RMS and the cadence, defined as ${\rm RMS}\times \sqrt{{\rm cadence}}$, in units of mmag per minute (see also Burke et al. 2008, Section 2.2).

Table 2. List of Light Curves Analyzed in this Work

Start Date UT Start Time (hh:mm) UT Observatory + Telescope Filter u1 u2 Exp. Time (s) Cadence (minute−1) Duration (hr) Start AM End AM β RMS (mmag) PNRa (mmag minute−1)
2007 May  2 00:08 Euler 1.2 m V 0.340 0.444 80 0.44 3.86 2.10 2.24 1.0 1.2 1.8
2007 May  4 17:39 Wise 0.46 m clear 0.343 0.398 20 1.23 3.19 1.03 1.14 1.6 2.0 1.8
2007 May  4 17:22 Wise 1 m V 0.340 0.444 60 0.65 2.90 1.05 1.07 1.5 0.9 1.1
2007 May 23 03:56 FLWO 1.2 m z 0.088 0.522 10 2.44 3.91 1.02 1.99 1.3 2.0 1.3
2007 May 25 19:33 Wise 1 m R 0.343 0.398 40 0.88 3.60 1.15 3.79 1.3 2.4 2.6
2007 May 31 03:30 FLWO 1.2 m z 0.088 0.522 10 2.38 3.20 1.02 1.64 2.1 2.1 1.4

Note. aPhotometric noise ${\rm rate} = {\rm RMS} \times \sqrt{ {\rm cadence}}$.

Download table as:  ASCIITypeset image

3. DATA ANALYSIS AND RESULTS

We decided not to include the publicly available Spitzer 8 μm light curve in our simultaneous fitting since the shape of that light curve is dependent on how the Spitzer "ramp" was modeled (Gillon et al. 2007b; Deming et al. 2007). As the Spitzer measurements are much more precise than the GB measurements, this ramp could have systematically influenced our results. For comparison, the Spitzer light curve PNR is a factor of ∼3 smaller than that of our ground-based light curves.

Accounting for correlated noise (Pont et al. 2006) was done similarly to the "time-averaging" method of Winn et al. (2008). After a preliminary analysis, we binned the residual light curves using bin sizes close to the duration of ingress and egress. The presence of correlated noise in the data was quantified as the ratio between the binned residual light-curves standard deviation and the expected standard deviation in the absence of correlated noise. For each light curve, we estimated β as the largest ratio among the bin sizes we used. We then multiplied the relative flux errors by β. Values of β are in the range of 1.0–2.1 and are listed in Table 2.

3.1. Simultaneous Fitting of All Parameters

We fitted the system parameters to all six GB light curves simultaneously. Our transit light-curve model consisted of nine parameters: the orbital period, P; a particular midtransit time, Tc; planet-to-star radius ratio, Rp/Rs; semimajor axis scaled by the stellar radius, a/Rs; impact parameter, b (see Equation (3) of Winn et al. 2007 for the formula of the impact parameter in an eccentric orbit); two limb-darkening coefficients u1 and u2, for a quadratic limb-darkening law; orbital eccentricity, e; and longitude of periastron, ω. The latter four parameters were held fixed in the fitting process, hence our model had five free parameters.

For transiting planets on eccentric orbits, the midpoint, in time, between transit start and end is different than the time of sky-projected star–planet closest approach. For GJ 436b, this difference is comparable to the typical uncertainty on the midtransit times obtained here. In our model, the midtransit time is defined as the time of closest approach, which is also the time of minimum light.

For e and ω, we used the values given by Maness et al. (2007). Although other authors (e.g., Demory et al. 2007; Deming et al. 2007) report more precise values for e, their uncertainty on ω is either much higher than that of Maness et al. (2007) or not reported. As the errors in two orbital parameters are correlated, we chose to adopt the reference giving the smallest uncertainty on both. To determine u1 and u2, we adopted a model of Teff = 3500 K (Maness et al. 2007), log g = 4.843 (Torres 2007), and [Fe/H] = −0.03 (Bonfils et al. 2005) and interpolated on the Claret (2000, 2004) grids. For the light curve obtained with the Wise Observatory 0.46 m telescope, which has no filter, we used limb-darkening coefficients corresponding to the R filter, as its CCD response resembles a "wide-R" filter (Brosch et al. 2008). The coefficients used for each light curve are listed in Table 2.

We used the formulas of Mandel & Agol (2002) to determine the relative flux as a function of the planet–star sky-projected separation and the Monte Carlo Markov Chain (MCMC) algorithm (e.g., Tegmark et al. 2004; Ford 2005; Gregory 2005) to determine the parameters that best fit the data, along with their uncertainties. MCMC algorithms have already been used extensively in the literature for fitting transit light curves (e.g., Holman et al. 2006; Collier Cameron et al. 2007; Gillon et al. 2008; Burke et al. 2007). Assuming a uniform prior for the five fitted parameters, the algorithm can be viewed as a random walk in the five-dimensional parameter space of the model. Starting from an initial guess, at each step the algorithm examines a new point in the parameter space which is reached by adding a Gaussian random value to each of the parameters. The algorithm decides whether to accept the new point in the parameter space and move there, or repeat the current point in the chain depending on the resulting posterior probability. If the new point has a higher relative posterior probability, the algorithm will move there, and if not, it will move to the new point with a relative probability of exp(−Δχ2/2), where Δχ2 is the χ2 difference between the two points in parameter space. The distribution of parameter values in all points is used to derive the best-fit value and its uncertainties. Here, we took the distribution median to be the best-fit value and the values at the 84.13 and 15.87 percentiles to be the +1σ and −1σ confidence uncertainties, respectively.

The widths of the Gaussian distributions used to find a new point are a crucial part of the MCMC algorithm and they affect the fraction of accepted points, f. In order to have f close to 25% (Gregory 2005; Holman et al. 2006), we divided each long chain, of 500, 000 steps, into smaller chains (minichains), of 1000 steps, and re-evaluated the distribution widths after the execution of each minichain. This re-evaluation was done by scaling them according to the relative difference between f of the previous minichain, and the target value of 25%. The same scaling was applied to all parameters.

To estimate the relative size of the distribution widths for each one of the model parameters, we initially performed a sequence of 10 minichains, where only that parameter was allowed to vary and the target f was 50% (Gregory 2005). Thus, a single MCMC run was comprised of several sequences of 10 minichains and a long chain of 500 minichains. For determining the final result, we considered only the long chain while ignoring the first 20% of its steps.

We performed 10 MCMC runs as described above while changing each time the initial guess and initial distribution widths. Results of those runs were highly consistent with each other and our final result is based on all runs.

Before describing our results, we briefly describe how we tested our analysis:

  • 1.  
    To test the adaptive form of our MCMC algorithm, we reran it while fixing the distribution widths during the long chain part. The differences in the results for the fitted parameters did not exceed 0.15σ, and the differences in the estimated uncertainties were no more than 10%.
  • 2.  
    We repeated our analysis while skipping the polynomial corrections mentioned in Section 2. The results for the fitted parameters were less than 0.2σ from the original ones, and the estimated uncertainties were similar at the 20% level.
  • 3.  
    We applied our analysis on the Spitzer 8 μm light curve alone, and derived results similar to Gillon et al. (2007b); Deming et al. (2007); Southworth (2008), at the ∼1σ level.
  • 4.  
    To further test the validity of our results, we reran MCMC while adopting a Gaussian prior probability for e and ω. The central values and standard deviations were those given by Maness et al. (2007). The results were not modified significantly.

Figure 1 presents the six GB light curves, with our fitted model overplotted in solid lines. The three light curves obtained in the beginning of 2007 May are of better quality, especially the Euler 1.2 m V and Wise 1 m V light curves, with RMS residuals of 1.2 and 0.9 mmag, respectively. The increased scatter in the three light curves obtained at the end of 2007 May is the result of the high airmass of these observations, done close to the end of the observational season.

Figure 1.

Figure 1. Six light curves analyzed in this work. The overplotted solid line is our best-estimate model. For each light curve, residuals from the model are plotted, centered on a relative flux of 0.98 for clarity. Within each panel the UT date, observatory, telescope, and filter used are given.

Standard image High-resolution image

Our estimates for the fitted parameters, while fitting all six light curves simultaneously, are listed in Table 3. For comparison, this table includes results from previous studies whenever these parameters were fitted directly. In case they were not, but can be calculated from other fitted parameters, they are given without an error. Gillon et al. (2007a) based their parameter determinations on an analysis of the Euler 1.2 m V light curve from UT 2007 May 2, which is included in this work. The other three studies referred to in Table 3 (Gillon et al. 2007b; Deming et al. 2007; Southworth 2008) are based on the Spitzer 8 μm light curve. Results obtained by Torres et al. (2008) are not presented since for the light-curve parameters of GJ 436b, these authors gave the weighted average of previous works.

Table 3. Light-Curve Fitted Parameters

Reference P (day) Tc-2454200 (HJD) Rp/Rs a/Rs b
This Work 2.64386 35.8355 0.085 13.6 0.85
  ±0.00003 ±0.0001 ±0.001 ±0.5 ±0.01
Gillon et al. (2007a)a     0.082    
      ±0.005    
Gillon et al. (2007b)a   80.78148 0.0830   0.849
    +0.00015     +0.010
    −0.00008     −0.013
Deming et al. (2007)a   80.78149 0.0839 13.2 0.85
    ±0.00016 ±0.0005 ±0.6 +0.03
          −0.02
Southworth (2008)a   80.78174 0.08284 13.68  
    ±0.00011 ±0.00090 ±0.51  

Note. aPeriod was fixed at the Maness et al. (2007) value of P = 2.64385 ±  0.00009 days.

Download table as:  ASCIITypeset image

Table 4 lists our results for the physical parameters of the system, along with the corresponding values obtained previously. Assuming the mass of the star to be Ms = 0.452 ± 0.013 M (Torres 2007), and adopting the orbital parameters of Maness et al. (2007; specifically K, e, and ω), we derived the planet mass and semimajor axis. Using Newton's version for Kepler's third law, we derive the star and planet radius.

Table 4. Results for the Physical Parameters

Reference a (AU) Mp (M) Rs (R) Rp (R)
This Work 0.02872 +0.00030 23.1 ±0.9 0.45 ±0.02 4.2 ±0.2
    −0.00025            
Gillon et al. (2007a)     22.6 ±1.9 0.44 ±0.04 3.95 +0.41
                −0.28
Gillon et al. (2007b)         0.463 +0.022 4.19 +0.21
            −0.017   −0.16
Deming et al. (2007) 0.0291  ±0.0004  22.24 ±0.95 0.47 ±0.02 4.33 ±0.18
Torres (2007) 0.02872 ±0.00027 23.17 ±0.79 0.464 +0.009 4.22 +0.09
            −0.011   −0.10
Southworth (2008)a         0.452 ±0.017 4.08 ±0.16

Note. aAssuming Ms = 0.452+0.014−0.012M (Torres 2007), orbital parameters of Maness et al. (2007), and the period derived here.

Download table as:  ASCIITypeset image

3.2. Fitting Transit Timing and Impact Parameter for Each Light Curve

As mentioned in Section 1, Ribas et al. (2008) suggested the presence of a second planet orbiting GJ 436, in an outer orbit close to a 2:1 mean-motion resonance with the transiting planet. The gravitational interaction between the two planets may cause a detected TTV signal (Holman & Murray 2005; Agol et al. 2005). Ribas et al. (2008) suggest that the second planet is responsible for changing the inclination angle i of the transiting planet orbit. The presence of such a planet may be detected through measuring transit timing and the impact parameter of many transit light curves. As the light curves obtained here are the first group of such light curves, we performed two additional MCMC analyses described below.

First, we fitted Rp/Rs, a/Rs, and b to all light curves and an independent Tc to each transit light curve to estimate the midtransit time of each event. In this analysis, the period does not determine the time intervals between transit events. It is used only for determining the true anomaly of each measurement within each light curve, so it is only weakly constrained. Hence, in this analysis P was fixed at the value derived in our original analysis. The resulting six Tc values were later used to derive the transit ephemeris by a linear fit. The new period is ≈0.35σ from the period fitted in our original MCMC run. Repeating the process described in this paragraph with the new period results in a variation of only 0.02σ, so no further iterations were done.

Table 5 contains the midtransit timings fitted independently to each light curve. Note that light curves with the smallest RMS residuals do not necessarily have the smallest midtransit time uncertainty. For example, the FLWO 1.2 m light curve, at E = 3, has relatively high RMS residuals although a small midtransit time uncertainty. This is simply the effect of the short KeplerCam cadence, as the increased number of points act to better constrain the midtransit time. This effect is quantified by the PNR (see Table 2).

Table 5. Midtransit Time and Impact Parameter Fitted Independently to Each Light Curvea

    ΔTcb   Δbc
E Tc-2454200 (HJD) (minute) (σ) b   (σ)
−5 22.61612 ± 0.00037 −0.15 −0.28 0.838 ± 0.017 −0.013 −0.61
−4 25.26002 ± 0.00038 −0.10 −0.18 0.828 ± 0.020 −0.023 −0.96
−4 25.26052 ± 0.00030 0.62 1.45 0.843 ± 0.016 −0.008 −0.37
3 43.76657 ± 0.00026 −0.71 −1.88 0.846 ± 0.018 −0.005 −0.23
4 46.40982 ± 0.00040 −1.58 −2.74 0.843 ± 0.014 −0.008 −0.41
6 51.69956 ± 0.00030 1.34 3.12 0.854 ± 0.016 0.003 0.13

Notes. aThe table gives the result of two separate MCMC analyses. In each analysis, only a single parameter (Tc or b) was fitted separately to every light curve and all other parameters were fitted to all light curves simultaneously, except for P which was fixed when fitting an individual Tc to each light curve. bDifference between the measured Tc and the expected midtransit time according to the fitted ephemeris (see Table 4). The difference is given in minutes and also in units of the uncertainty on each midtransit time. cDifference between the measured impact parameter while fitting each light curve separately, and the result of the simultaneous fit (see Table 4). The difference is also given in units of the uncertainty on each impact parameter.

Download table as:  ASCIITypeset image

By fitting a linear function,

Equation (1)

to the transit epoch, E, we derived P = 2.64386 ± 0.00003 days and Tc(0) = 2454235.8355 ± 0.0001 [HJD]. The transit event with E = 0 was chosen to be right in the middle of our observed events, although we did not observe that particular event. We take this ephemeris as our final result, listed in Table 3, and note that it is consistent with P and Tc derived from our original MCMC fitting.

The residuals from the linear fit are presented in Figure 2, known as the OC diagram. The figure shows some indication for a variability of the period, reflected by the high χ2 value of the fit, of 23 for 4 degrees of freedom. However, the small number of points does not allow to assign any significance to the detection of this variability, especially because some systematic effects could shift some of the points.

Figure 2.

Figure 2. Observed minus calculated (OC) transit timing of the six light curves included in this study. The graph shows the residuals from a linear fit to the transit timing as a function of the transit epoch. Table 5 lists the actual transit timings.

Standard image High-resolution image

Second, we performed an MCMC analysis where we looked for a transit impact parameter variation. The four parameters P, Tc, Rp/Rs, and a/Rs were fitted to the entire data while the impact parameter, b, was fitted independently to each light curve. The derived impact parameters are listed in Table 5 and plotted in Figure 3. As a reference, we overplotted in Figure 3 the value of our original run, in a solid line, and the upper and lower 1σ confidence limits in dashed lines. No significant variation can be seen in the data, although a long-term trend of the order of 0.2 yr−1 cannot be ruled out.

Figure 3.

Figure 3. Impact parameter b fitted independently to each light curve, as a function of the transit epoch. The solid line marks the best estimate derived in our original MCMC run, while fitting the impact parameter to all light curves simultaneously. The dashed lines mark the upper and lower 1σ confidence limits. Table 5 lists the actual values of the impact parameters presented here.

Standard image High-resolution image

The two separate analyses described in this section were repeated without applying the polynomial corrections mentioned in Section 2. The results were consistent with the original ones. In addition, for both the OC transit timing residuals and the impact parameters derived here, we looked for a possible correlation with several external parameters, including midtransit airmass, filter central wavelength, and the light-curves polynomial correction amplitude. All correlations were below the 1.5σ confidence level.

4. DISCUSSION AND SUMMARY

This paper presents an analysis of six complete transit light curves of GJ 436b, spanning 29 days, or 11 orbital periods.

Our period is consistent with the period given by Maness et al. (2007), of 2.64385 ± 0.00009 days, and is more precise. Our value for Rp/Rs (see Table 3) is larger by 1σ–2σ, than the radii ratio derived for the Spitzer 8 μm light curve by Gillon et al. (2007b), Deming et al. (2007), and Southworth (2008). If this difference is real, it may be induced by stellar spots, or wavelength-dependent opacities in the planetary atmosphere (Barman 2007; Pont et al. 2008, 2009). Interestingly, a similar effect, i.e., an increased radius ratio in the optical relative to the IR, was already observed for HD 189733b (Pont et al. 2007).

Our derived physical parameters, mainly the planetary radius and mass, are consistent with those derived previously (see Table 4).

We examined the individual transit timing of each light curve (TTV, see Figure 2) and derived the impact parameter of each light curve independently (see Figure 3). The transit timing residuals show a possible hint for variability, supported by the high χ2 of the fit. However, many more light curves are needed to establish the variability and to explore its nature. The impact parameter plot shows a possible trend, although a linear fit to the six data points gives a slope consistent with zero.

Comparing our derived transit ephemeris to transit times published recently (Alonso et al. 2008; Bean et al. 2008) shows that the predicted transit times are somewhat earlier than the measured times. The difference is at the 1σ–2σ level, so it is not highly significant. However, if real, it supports the claim made by Bean et al. (2008) for a long-term drift in the transit times. Together with the long-term slope in GJ 436 RVs (Maness et al. 2007), it suggests the existence of another object in the system with a long period. Future photometric and spectroscopic data will be able to study this possibility better.

As mentioned in the introduction, the impact parameter of GJ 436b can be changed by a precession of the line of nodes, induced by a second planet whose orbital plane is inclined relative to the orbit of GJ 436b (e.g., Miralda-Escudé 2002). Such an effect was already suggested by Soderhjelm (1975) and Mazeh & Shaham (1976) for stellar triple systems. The precession period is

Equation (2)

(Miralda-Escudé 2002), where Mp,2 and P2 are the mass and period of the second planet, respectively. For the planet suggested by Ribas et al. (2008), this precession is of the order of 10 yr.

For eccentric orbits, the impact parameter can also be modified by the precession of the periastron because the impact parameter depends on the orientation of the line of apsides (see Equation (3) of Winn et al. 2007). The apsidal motion can be driven by the stellar quadrupole moment, or by a second planet. The precession driven by the stellar quadrupole moment could be quite slow, and therefore probably does not contribute to the variation of the impact parameter. The rate of the apsidal precession induced by another planet is of the same order of magnitude as the rate of the nodal precession, estimated in Equation (2), and therefore can be of the order of years, depending on the parameters of the unseen additional planet in the system.

To estimate the possible implication of the apsidal precession on the impact parameter, we plot in Figure 4 the impact parameter of GJ 436b as a function of the argument of the periastron ω, at fixed orbital eccentricity and inclination. We can see that b is modulated between ∼1.0 and 0.7 over the apsidal precession period. The upper dashed line of the figure indicates that at some phase of the precession period the transit of GJ 436b will even lose its flat-bottom shape. The figure suggests that if the apsidal motion is of the order of 10 years, we should be able to observe soon a change of the impact parameter. In case this modulation is due to apsidal motion, we will then be able to confirm the change of the line of apsides direction by timing the secondary eclipse or by precise RV measurements.

Figure 4.

Figure 4. Impact parameter, b, vs. ω, for the eccentricity and inclination of GJ 436b (the solid line). The filled circle (green) marks the current position, error bars are too small to be shown. The upper and lower dashed lines (red) mark the impact parameter when the stellar radius is multiplied by a factor of 1 − Rp/Rs and 1 + Rp/Rs, respectively. When multiplying the stellar radius by 1 − Rp/Rs (the upper dashed line), the impact parameter will be smaller than one only for nongrazing transits. When multiplying by 1 + Rp/Rs (the lower dashed line), the impact parameter will equal one when the minimal star–planet sky-projected distance equals the radii sum, i.e., for b > 1 there will be no transit at all.

Standard image High-resolution image

We are deeply thankful to Michael Gillon for his photometric processing of the data obtained on 2007 May 4 and May 25. We thank the anonymous referee for his thorough reading of the paper and his useful comments which helped improve the paper. T. Mazeh and A. Shporer thank Elia Leibowitz and Liliana Formiggini for allowing the use of their telescope time at the Wise Observatory 1 m telescope on 2007 May 4. This work was partly supported by Grant No. 2006234 from the United States–Israel Binational Science Foundation (F), by the Kepler Mission under NASA Cooperative Agreement NCC 2-1390 with the Smithsonian Astrophysical Observatory, and by NASA Origins grant NNG06GH69G (to M. J. H.).

Footnotes

  • For an updated list, see: http://exoplanet.eu/.

  • The Image Reduction and Analysis Facility (IRAF) is distributed by the National Optical Astronomy Observatories, which is operated by the Association of Universities for Research in Astronomy, Inc., under cooperative agreement with the National Science Foundation.

Please wait… references are loading.
10.1088/0004-637X/694/2/1559