Transit Timing Variations for Inclined and Retrograde Exoplanetary Systems
aa r X i v : . [ a s t r o - ph . E P ] M a r Draft version October 31, 2018
Preprint typeset using L A TEX style emulateapj v. 08/22/09
TRANSIT TIMING VARIATIONS FOR INCLINED AND RETROGRADE EXOPLANETARY SYSTEMS
Matthew J. Payne, Eric B. Ford and Dimitri Veras
Department of Astronomy, University of Florida, 211 Bryant Space Science Center, P.O. Box 112055, Gainesville, FL 32611-2055
Draft version October 31, 2018
ABSTRACTWe perform numerical calculations of the expected transit timing variations (TTVs) induced on aHot-Jupiter by an Earth-mass perturber. Motivated by the recent discoveries of retrograde transitingplanets, we concentrate on an investigation of the effect of varying relative planetary inclinations,up to and including completely retrograde systems. We find that planets in low order (E.g. 2:1)mean-motion resonances (MMRs) retain approximately constant TTV amplitudes for 0 < ◦ i < ◦ ,only reducing in amplitude for i > ◦ . Systems in higher order MMRs (E.g. 5:1) increase in TTVamplitude as inclinations increase towards 45 ◦ , becoming approximately constant for 45 < i < i > ◦ . Planets away from resonance slowly decrease in TTV amplitudeas inclinations increase from 0 to 180, where-as planets adjacent to resonances can exhibit a hugerange of variability in TTV amplitude as a function of both eccentricity and inclination. For highlyretrograde systems (135 ◦ < i ≤ ◦ ), TTV signals will be undetectable across almost the entiretyof parameter space, with the exceptions occurring when the perturber has high eccentricity or is veryclose to a MMR. This high inclination decrease in TTV amplitude (on and away from resonance) isimportant for the analysis of the known retrograde and multi-planet transiting systems, as inclinationeffects need to be considered if TTVs are to be used to exclude the presence of any putative planetarycompanions: absence of evidence is not evidence of absence. Subject headings: celestial mechanics — planetary systems — methods: numerical INTRODUCTION
Of the known extra-Solar planets, more than 60 transit the host star. Of these systems, at least four show evidencefor an external companion (GJ 436, Maness et al. 2007; HAT-P-13, Bakos et al. 2009; HAT-P-7, Winn et al. 2009 andCoRoT-7 Queloz et al. 2009).If the transiting planet were the only planet in the system, then the period between each successive transit wouldbe constant (neglecting complicating effects such as general relativity, stellar oblateness and tides). The presence ofadditional planets in the system (which themselves may or may not transit the star) can cause perturbations to theorbit of the transiting planet, leading to detectable transit timing variations (TTVs) of the known transiting planet(Miralda-Escud´e 2002; Holman and Murray 2005; Agol et al. 2005). Other studies have extended the TTV method todemonstrate the feasibility of detecting planetary moons (Kipping 2009 a , b ) and Trojan companions (Ford and Holman2007) in extra-solar planetary systems.However, practical use of TTVs as a detection tool requires the solution of the difficult “inverse problem”, i.e. givena particular TTV profile, can one reconstruct (or at least restrict) the mass and orbit of the unseen perturber? Thisproblem is non-trivial, as numerous different perturber mass-orbit configurations can lead to degenerate TTV solutions(e.g. Nesvorn´y and Morbidelli 2008).Several investigations have considered (and ruled out) planets more massive than the Earth in certain orbits / regionsof parameter space close to many of the known transiting planets (e.g. Steffen and Agol 2005; Agol and Steffen2007; Alonso et al. 2008; Bean and Seifahrt 2008), while Nesvorn´y and Morbidelli (2008) and Nesvorn´y (2009) havedeveloped an approximate analytic method to try and tackle the inverse problem in a more general manner.However, most previous investigations have focused on analyzing the effect of prograde, coplanar companions. Evenwhere explicit investigations of inclination effects have been conducted (e.g. Nesvorn´y 2009), the investigations havebeen restricted to prograde cases, examining only relative inclinations significantly lower than 90 ◦ .There are several reasons to consider TTVs of inclined systems. Observationally, it is possible to measure the sky-projected angle between the spin vector of a star, and the orbital angular momentum vector of a planet transiting thatstar via the Rossiter-McLaughlin effect. Measurements over the past year have revealed that a number of the knowntransiting planets are strongly inclined, or even retrograde (e.g. HAT-P-7b, Narita et al. 2009 and Winn et al. 2009,as well as WASP-17b, Anderson et al. 2009).Although the mechanism(s) driving the creation of such retrograde planets remains unclear, models used to ex-plain high eccentricity exo-planets via planet-planet scattering (Rasio and Ford 1996; Weidenschilling and Marzari1996; Lin and Ida 1997; Levison et al. 1998; Papaloizou et al. 2001; Moorhead and Adams 2005; Chatterjee et al. 2008;Juri´c and Tremaine 2008; Ford and Rasio 2008) and/or Kozai oscillations (E.g. Kozai 1962; Takeda and Rasio 2005;Nagasawa et al. 2008) can naturally excite large orbital inclinations (although the quantitative details of the distri- Electronic address: [email protected]
Payne, Ford & Verasbutions can differ greatly). Other mechanisms such as inclination-pumping during migration (e.g. Lee and Thommes2009) may also contribute.In addition, previous dynamical studies have suggested that systems of prograde and retrograde planets are more stable than standard prograde-prograde cases (Gayon and Bois 2008; Gayon et al. 2009; Smith and Lissauer 2009),essentially because the planets spend less time “close together” and thus perturbations are smaller.Given the observational evidence for highly inclined and retrograde systems, and the suggestion of enhanced stabilityin such systems, we investigate and quantify the hypothesis that retrograde systems will have a significantly reducedTTV profile compared to a standard prograde case. METHODOLOGY
Integration and TTVs
We consider a fiducial system consisting of two planets: a transiting Hot-Jupiter (with mass M = M Jupiter ,eccentricity e = 0 and semi-major axis a = 0 . AU ) and a small non-transiting, Earth-mass planet which perturbsthe transit times of the larger planet. We evolve each system forward in time for ∼
10 years, equating to ∼ , Mercury (Chambers 1999). We use a barycentric coordinate system and limit the time steps to no more than 0.04 times theorbital period. After each time step, we test whether the projected star-planet separation (∆) passed through a localminimum and the planet in question is closer to the observer than the star. Each time these conditions are met, wefind the nearby time that minimizes ∆ via Newton-Raphson iteration and increment an index i . If the minimum ∆ isless than the stellar radius, then we record the mid-time of the transit, t i . Observations of transit times are perturbedby the light travel time, δ t ltt ( i ) ≃ − ( r p ( t i ) · r los ) /c , where r p ( t i ) is the barycentric vector of the planet at time t i , r los is a unit vector pointing to the observer, and c is the speed of light. The observable transit time variations arecalculated as δ t ( i ) = t i + δ t llt ( i ) − i × ˆ P − ˆ t , where the constant ˆ P and ˆ t are determined by linear least squaresminimization of P i ( δ t ( i )) . We neglect any motion of the stellar center between the time of light emission and thetime of transit. TTV Maps In § i , e ) and ( a a , e ) parameter-space. These plots are producedby fixing the mass and the initial orbit of the transiting hot-Jupiter ( e = 0 & a = 0 . a , eccentricity, e , & inclination, i of the perturbing planet- we conduct 5 simulations randomising the other orbital elements (argument of pericentre, ω , longitude of ascendingnode, Ω, & initial mean anomaly M ).We then calculate the RMS of δt ( i ) for each a, e, i configuration. The required contour plots are then plotted usingan approximately logarithmic color scale (shown in Figs. 2 and 3).We restrict our simulations to the region defined by the prograde coplanar three-body Hill stability criterion ofGladman (1993). This ensures that all of the i = 0 ◦ simulations are stable, and ensures for the i > ◦ systems thatwe have covered all of the guaranteed stability zone and more (see Veras and Armitage (2004) for a discussion of thestability limit as a function of inclination). As an additional check on systems with extremely large TTVs close to theHill-stability criterion boundary, we conducted 1,000 longer-term Mercury simulations. In most cases the Earth-massplanet experiences significant changes in semi-major axis within 1,000 years, due to the “kicks” at pericentre passage.Thus, although these planets do not collide, the less massive planet can experience some orbit variation. RESULTS
Example TTVs For a Hot-Jupiter System
In Fig 1a we show a circular hot-Jupiter at 0.05 au being perturbed by an Earth-mass companion with e = 0 . ◦ ≤ i ≤ ◦ )have TTV profiles that differ significantly, their overall amplitudes are similar. In contrast, the amplitude of the i = 180 ◦ signal is very low ( ∼ i = 180 ◦ signal is extremely low, thereare now significant differences in the behaviour at smaller angles: the i = 0 ◦ and i = 180 ◦ now have similarly lowamplitudes, with an approximate maximum occurring at i = 90 ◦ .Away from resonance (Figs 1d, 1e & 1f), we again find that the i = 180 ◦ signal is generally significantly lowerthan that observed for cases with i < ◦ , but that the detailed dependence on inclination is far from obvious, withsome cases (Figs 1d & 1f) having the highest amplitude at i = 0 ◦ , while the low eccentricity plot of Fig. 1e has aparticularly low amplitude signal at i = 0 ◦ .We note that there are a number of frequencies evident in Fig 1, with periods from a few months to >>
10 years.In general, the short-term oscillations tend to be driven by variations in e , leading to variations in the distance ofclosest approach between the two planets. The longer term quadratic trends are primarily due to outward drifts in a /a as a result of kicks at close approach.The large number of different frequencies evident in the plots also serves to highlight the importance of consideringthe sampling period over which observations and simulations are to be conducted and compared.TVs: Inclined & Retrograde Systems 3 -200-150-100-50 0 50 100 150 0 1 2 3 4 5 6 7 8 9 10a) 2:1 Resonance, e=0.025, i=0 -150-100-50 0 50 100 150 0 2 4 6 8 10b) 3:1 Resonance, e=0.12 -150-100-50 0 50 100 0 2 4 6 8 10c) 5:1 Resonance, e=0.05i=0i=45i=90i=180-10-5 0 5 10 15 20 0 1 2 3 4 5 6 7 8 9 10d) Adjacent to 2:1 Resonance, e=0.05 -0.15-0.1-0.05 0 0.05 0.1 0.15 0 2 4 6 8 10e) Non-Resonant, e=0.15 -20-15-10-5 0 5 10 15 20 0 1 2 3 4 5 6 7 8 9 10f) Non-Resonant, e=0.55 R M S TT V ( s e c o nd s ) Time (years)
Fig. 1.—
The effect of varying inclination at fixed semi-major axes and eccentricities. Prograde orbits with: i = 0 ◦ have solid red lines; i = 45 ◦ have dotted green lines; 90 ◦ have dashed blue lines; while retrograde i = 180 ◦ orbits have dot-dashed black lines. In (a) weshow the TTVs for a system on/near the 2:1 MMR ( a /a = 1 . P /P = 2 . e = 0 . a /a = 2 . P /P = 3 . e = 0 . a /a = 2 . P /P = 5 . e = 0 . a /a = 1 . P /P = 2 . e = 0 . a /a = 3 . P /P = 6 . e = 0 . a /a = 3 . P /P = 6 . e = 0 . i = 180 ◦ case tends to be significantly lower in amplitude than for other inclinations. Furthermore, TTV profiles may change in amplitude and frequency by several orders of magnitude simply due toa variation in the initial values of the mean anomalies (Veras and Ford 2009). This means that while approximateexpressions for TTV amplitude dependencies on mass, semi-major axis separation, etc, can serve as a useful guide,individual system variations can be far removed from the mean. We do not intend this statement as a criticism ofprevious work, but rather as a cautionary note against over-reliance on approximations when analyzing individualsystems which can exhibit great sensitivity to initial conditions. We do not seek in this paper to address in detail allsuch issues, but instead defer a more detailed consideration to a companion paper
Veras et al. 2010, in prep. , where amore thorough examination of the coplanar prograde case will be presented, allowing a much more detailed expositionof the huge number of variables which can affect the TTV signal.
Eccentricity-Inclination Contour Plots
To better understand the TTV amplitudes as a function of inclination, for each of the semi-major axis plots in Fig 1we perform a search of the ( e , i ) parameter space, producing contour plots in Fig 2 showing the RMS TTV amplitudevariation at that semi-major axis separation as the eccentricity and inclination of the outer planet is varied.In the 2:1 case of Fig 2a the resulting contour plot is rather uniform, with the only area of significant reductionbeing in the e < . i ∼ ◦ region. An approximately similar result can be seen for the 3:1 case.For the 5:1 case (Fig 2c), the situation changes: the region of low amplitude close to i = 180 ◦ significantly expands,while another low amplitude region emerges at e < . i < ◦ .For Fig 2d the perturber is located very close to (but not on) the 2:1 MMR. The plotted behavior is particularlyrich, showing huge variations in the predicted TTV amplitude as functions of both eccentricity and inclination.Finally, the non-resonant case in Fig 2e with a /a = 3 .
6, shows a fairly simple decline in TTV amplitude as afunction of inclination.
TTV Maps for Varying a a We now look at a different projection of the data, looking in the ( a a , e ) plane at the TTV structure as a functionof inclination.In Fig 3a we reproduce the coplanar prograde results of Agol et al. (2005), showing the pronounced “flames of reso-nance” (Veras and Ford 2009). The majority of the non-resonant parameter-space has TTV signals of comparativelylow amplitude ( <
30 seconds). It is only when eccentricities become very high, or the planets are in / near MMRs,that the TTV signal amplitudes become large.As the relative inclination of the planetary orbits is increased, for low inclinations changes in TTV amplitude are Payne, Ford & Veras
RMS TTV (seconds)0-11-33-1010-3030-100100-300300-10001000+No Simulation E cc e n t r i c i t y o f t h e O u t e r B o d y ( e ) Relative Inclination of the Outer Body ( i ). Fig. 2.—
Plots at fixed ( a /a ) showing the median RMS TTV amplitude as a function of e & i . The contour key is displayed in thebottom right-hand corner. The semi-major axes selected for plots (a) - (e) are the same as those of Fig 1a - Fig 1e. We can see significantdifferences in behaviour as a function of inclination depending on (i) position on / off resonance, (ii) the type of any resonance, and (iii)the eccentricity of the the perturber. relatively subtle; for the i = 30 ◦ & i = 45 ◦ cases (Fig 3c & 3d) the RMS TTV signal is reduced in some high-eccentricity locations, corresponding to the reduction previously seen in the close-to-resonance map of Fig 2d. Inaddition, we can also see that spikes around higher resonances start to drop lower down, i.e. the low eccentricityamplitude increases as the inclination decreases, an effect noted in § a a , e ) param-eter space (see the i = 90 ◦ plot in Fig 3d). The reduction in high amplitude signals at high eccentricities is the mostmarked decrease, but at any given eccentricity one tends to find a reduction in amplitude (shifting of the contours)compared to the i = 0 ◦ plot.When the systems are pushed into retrograde orbits, the difference becomes more pronounced (E.g. i = 135 ◦ in Fig3e and i = 180 in Fig 3f). For these inclinations we find a reduction in amplitude across large swathes of the plot,with significant amplitude remaining only near MMRs. In addition, the TTV amplitude and the width of regions withlarge TTVs (near MMRs) is markedly reduced. External 2:1 MMR
Since highly inclined systems only have large RMS TTVs in / near MMRs (Fig 3f), we examine in more detailone such resonant region. The external 2:1 MMR is known to be of importance for exo-planet systems from bothTVs: Inclined & Retrograde Systems 5 E cc e n t r i c i t y o f t h e O u t e r B o d y ( e ) Semi-major Axis Ratio ( a /a ). Fig. 3.—
Median RMS TTV amplitudes in the ( a /a , e ) plane for various i . The contour key is shown at the top-left in plot (a). Ineach ( a /a , e ) plot, we see a strong increase in amplitude (i) at high eccentricity, and (ii) close to MMRs. As the inclination is increasedthere is an obvious general decrease in TTV amplitude, particularly pronounced for the retrograde orbits. theoretical (e.g. S´andor et al. 2007) and observational (e.g. Tinney et al. 2006) studies.On and adjacent to the 2:1 MMR region, the TTV signal is extremely sensitive to changes in semi-major axis and/ or eccentricity. In the co-planar prograde case (Fig 4a), a fractional increase of just a few percent in either a or e causes the TTV amplitude to jump by more than an order of magnitude.As the inclination increases to 30 ◦ and 45 ◦ , we find that the central region of the MMR maintains an approximatelyconstant TTV amplitude. The regions just outside the resonance feel the strongest effects. E.g. regions near a /a =1 .
54 show TTV amplitudes reduced by an order of magnitude for many eccentricities (compare Fig 4a and Fig 4c).This decreases the width of the region of parameter space around the 2:1 MMR resonance with TTV signals > a /a ∼ . ◦ and 135 ◦ , TTV signals >
30 seconds are restricted purely to the resonanceregions, with regions away from resonance having signals <
30 seconds almost everywhere, even at high eccentricity.As the inclination reaches 180 ◦ , and the planets are now again coplanar but with completely counter-rotating orbits,we find TTV amplitude over 100 seconds only in very narrow regions around the first- and second-order MMRs.Finally, we note that as we move from Fig 4e to Fig 4f, the amplitude in regions slightly offset from resonance appearsto increase slightly at high eccentricities. However this only occurs at low order resonances (small semi-major axisseparations - see Fig 3e and Fig 3f) as well as Fig 2. More generally, away from such regions (see Fig 3f), retrogradeorbits tend to produce the smallest TTVs. Payne, Ford & Veras E cc e n t r i c i t y o f t h e O u t e r B o d y ( e ) Semi-major Axis Ratio ( a /a ). Fig. 4.—
Close-up of the 2:1 MMR region, demonstrating extreme sensitivity to variations in a /a , e & i . The contour key is thesame as that used in Figs. 2 & 3. Internal Perturbers
Finally, for completeness we consider an Earth-mass planet on an orbit interior to the Jupiter-mass planet, keepingthe more massive planet at a = 0 . AU and e = 0, and then monitor the TTVs of the massive planet due to theinterior perturber.From Fig 5 it seems that the background (non-resonant) regions now retain an approximately constant TTV signalirrespective of inclination, while there is a tendency for the TTV signal to become stronger around the MMR regionsas the inclination is increased towards 90 ◦ , but then drop away to almost nothing as the inclination further increasestowards i = 180 ◦ . CONCLUSIONS
We have investigated the TTV signals for systems of highly inclined and retrograde planets. We find that: • In the vicinity of exterior MMRs the inclination dependence is complex: low order resonances maintain a highTTV amplitude for all regions i < ◦ , declining in amplitude only for low eccentricity cases close to i = 180 ◦ ,where-as higher order resonances display an increase in TTV amplitude as inclinations rise from 0 to 45 degrees.Moreover, the regions immediately adjacent to MMRs show extreme sensitivity to changes in perturber a , e and i . • Exterior perturbers away from resonances tend to show a slow decrease in TTV amplitude with increasingTVs: Inclined & Retrograde Systems 7 E cc e n t r i c i t y o f t h e O u t e r B o d y ( e ) Semi-major Axis Ratio ( a /a ). Fig. 5.—
Internal perturber parameter space. Top; a general view. Bottom; a close-up of the 1:2 MMR. The contour key is the sameas in Fig 4. We find that the behaviour in the internal region differs from the external region: the background stays constant and there isnow a strengthening of the resonant regions as the inclination increases towards 90 ◦ . inclination, although regions adjacent to resonances can show remarkably complex behaviour. • Perturbing planets on interior orbits display a slightly different behaviour: Away from resonance the amplituderemains approximately constant with inclination, but around MMRs the perturbations become stronger as theinclination increases towards 90 ◦ before decreasing again beyond 90 ◦ .We note that: • Absence of evidence is not evidence of absence: Planets in retrograde orbits should be expected to producemarkedly reduced TTV signals as compared to the standard prograde case. For an Earth-mass perturber in ananti-aligned orbit, almost the entirety of the sample parameter space would result in a very small or undetectableTTV signal (unless the planet happened to be fortuitously located precisely on a strong MMR). • Retrograde orbits may be a natural way to explain transiting systems in which little or no TTV signal is observed,but in which the radial velocity observations point towards the existence of an additional planetary companions(E.g. GJ 436, HAT-P-13, etc). • In addition to the TTV considerations in this work, it is important to acknowledge that inclined orbits inmulti-planet systems will precess, leading to variations in transit duration. More work is required to try andunderstand whether a combination of TTV signals with transit duration variation signals could remove some ofthe degeneracies inherent in this problem. • Many additional dependencies (E.g. sampling period, perturber mass, relative mean anomaly, etc, etc) can signif-icantly alter the expected TTV signal for a particular system. We defer the provision of a detailed investigationof such matters to a companion paper,
Veras et al. 2010, in prep. • Future work could also investigate in more detail the regions above the Hill stability curve, seeking to identifythe (likely large) TTVs for any stable systems which exist in that region. Payne, Ford & VerasThis material is based upon work supported by the National Science Foundation under Grant No. 0707203 andis also based on work supported by the National Aeronautics and Space Administration under grant NNX08AR04Gissued through the Kepler Participating Scientist Program. We are very grateful to Eric Agol and our anonymousreferee for providing detailed comments which helped to significantly clarify sections of this letter.
REFERENCESAgol E and Steffen J H 2007
MNRAS , 941–948.Agol E, Steffen J, Sari R and Clarkson W 2005
MNRAS , 567–579.Alonso R, Barbieri M, Rabus M, Deeg H J, Belmonte J A andAlmenara J M 2008
A&A , L5–L8.Anderson D R, Hellier C, Gillon M, Triaud A H M J, Smalley B,Hebb L, Collier Cameron A, Maxted P F L, Queloz D, WestR G, Bentley S J, Enoch B, Horne K, Lister T A, Mayor M,Parley N R, Pepe F, Pollacco D, S´egransan D, Udry S andWilson D M 2009
ArXiv e-prints .Bakos G A, Howard A W, Noyes R W, Hartman J, Torres G,Kovacs G, Fischer D A, Latham D W, Johnson J A, MarcyG W, Sasselov D D, Stefanik R P, Sipocz B, Kovacs G,Esquerdo G A, Pal A, Lazar J and Papp I 2009
ArXiv e-prints .Bean J L and Seifahrt A 2008
A&A , L25–L28.Chatterjee S, Ford E B, Matsumura S and Rasio F A 2008
ApJ , 580–602.Ford E B and Holman M J 2007
ApJ , L51–L54.Ford E B and Rasio F A 2008
ApJ , 621–636.Gayon J and Bois E 2008
A&A , 665–672.Gayon J, Bois E and Scholl H 2009
Celestial Mechanics andDynamical Astronomy , 267–279.Gladman B 1993
Icarus , 247–+.Holman M J and Murray N W 2005
Science , 1288–1291.Juri´c M and Tremaine S 2008
ApJ , 603–620.Kipping D M 2009 a MNRAS , 181–189.Kipping D M 2009 b MNRAS , 1797–1804.Kozai Y 1962 AJ , 591–+.Lee M H and Thommes E W 2009 ApJ , 1662–1672.Levison H F, Lissauer J J and Duncan M J 1998 AJ , 1998–2014. Lin D N C and Ida S 1997 ApJ , 781–+.Maness H L, Marcy G W, Ford E B, Hauschildt P H, Shreve A T,Basri G B, Butler R P and Vogt S S 2007
PASP , 90–101.Miralda-Escud´e J 2002
ApJ , 1019–1023.Moorhead A V and Adams F C 2005
Icarus , 517–539.Nagasawa M, Ida S and Bessho T 2008
ApJ , 498–508.Narita N, Sato B, Hirano T and Tamura M 2009
PASJ , L35–L40.Nesvorn´y D 2009 ApJ , 1116–1122.Nesvorn´y D and Morbidelli A 2008
ApJ , 636–646.Papaloizou J C B, Nelson R P and Masset F 2001
A&A , 263–275.Queloz D, Bouchy F, Moutou C, Hatzes A, Hegrard G, Alonso R,Auvergne M, Baglin A, Barbieri M and Barge P 2009
A+Apreprints .Rasio F A and Ford E B 1996
Science , 954–956.S´andor Z, Kley W and Klagyivik P 2007
A&A , 981–992.Smith A W and Lissauer J J 2009
Icarus , 381–394.Steffen J H and Agol E 2005
MNRAS , L96–L100.Takeda G and Rasio F A 2005
ApJ , 1001–1010.Tinney C G, Butler R P, Marcy G W, Jones H R A, Laughlin G,Carter B D, Bailey J A and O’Toole S 2006
ApJ , 594–599.Veras D and Armitage P J 2004
MNRAS , 613–624.Veras D and Ford E B 2009 in ‘IAU Symposium’ Vol. 253 of IAUSymposium pp. 486–489.Weidenschilling S J and Marzari F 1996
Nature , 619–621.Winn J N, Johnson J A, Albrecht S, Howard A W, Marcy G W,Crossfield I J and Holman M J 2009