Transit duration and contact point times

Transit and in- and egress duration for circular orbit

Functions to estimate the duration of the transit and the in- and egress for circular orbits based on input in solar system units (AU, solar and Jovian radius)

or parameters expressed in terms of the stellar radius

Example: Estimate duration of Earth’s transit and ingress (solar system units)

from __future__ import print_function, division
from PyAstronomy import pyasl
from PyAstronomy import constants as pc

# Earth radius expressed in Jovian radii
reJ = pc.REarth/pc.RJ
print("Earth radius in Jovian units: ", reJ)

# Estimate the duration of Earth's transit
td = pyasl.transitDuration(1.0, reJ, 1.0, 90.0, 365.0)
print("The transit of Earth lasts about: %5.3f days" % td)

# Estimate the duration of in- and egress
ti = pyasl.ingressDuration(1.0, reJ, 1.0, 90.0, 365.0)
print("The in- and egress last: %6.4f days or %4.2f hours" % (ti, ti*24))

Example: Estimate transit and ingress duration of HD 189733 b (stellar units)

from __future__ import print_function, division
from PyAstronomy import pyasl

# Semi-major axis in units of stellar radius
sma = 8.8
# Radius ratio (Rp/Rs)
rprs = 0.16

# Estimate the duration of Earth's transit
td = pyasl.transitDuration_Rs(sma, rprs, 85.7, 2.21858)
print("The transit of HD 189733 b lasts about: %5.2f hours" % (td*24.))

# Estimate the duration of in- and egress
ti = pyasl.ingressDuration_Rs(sma, rprs, 85.7, 2.21858)
print("The in- and egress of HD 189733 b lasts: %5.2f hours" % (ti*24.))

Function documentation

PyAstronomy.pyasl.transitDuration(sma, rp, rs, inc, period, exact=False)

Calculate the transit duration.

The routine calculates the transit duration assuming a spherical star and planet and a circular planetary orbit.

It evaluates the expression

\[T_D = \frac{P}{\pi} \arcsin\left(\frac{\sqrt{(R_s+R_p)^2 - b^2}}{a} \right)\]

where P is the orbital period, b the impact parameter (\(b=a \cos(i)\)), and a the semi-major axis.

If the exact flag is set True, the slightly more accurate expression

\[T_D = \frac{P}{\pi} \arcsin\left(\frac{\sqrt{(R_s+R_p)^2 - b^2}}{a \sin(i)} \right)\]

is evaluated, where the additional sine term accounts for orbit curvature.

Note

The units of the transit duration are the same as the units of the input orbital period.

Parameters
smafloat

The semi-major axis in AU.

rpfloat

The planetary radius in Jovian radii.

rsfloat

The stellar radius in solar radii.

incfloat

The orbital inclination in degrees.

periodfloat

The orbital period.

exactboolean, optional

If True, a slightly more accurate expression is evaluated. The default is False.

Returns
Transit durationfloat

The duration of the transit (same units as period).

PyAstronomy.pyasl.transitDuration_Rs(sma, rprs, inc, period, exact=False)

Calculate transit duration

Invokes transitDuration() after parameter transformation.

Parameters
smafloat

Semi-major axis [stellar radii]

rprsfloat

Planet-to-star radius ratio (Rp/Rs)

incfloat

Orbital inclination [deg]

periodfloat

Orbital period

exactboolean, optional

If True, a slightly more accurate expression is evaluated. The default is False.

Returns
Transit durationfloat

The duration of the transit (same units as period).

PyAstronomy.pyasl.ingressDuration(sma, rp, rs, inc, period)

Calculate the ingress (and egress) duration of the transit.

The routine calculates the transit in- and egress duration assuming a spherical star and planet and a circular planetary orbit.

It evaluates the expression

\[T_I = \frac{P}{2 \pi} \left( \arcsin\left(\frac{\sqrt{(R_s+R_p)^2 - b^2}}{a \sin(i)} \right) - \arcsin\left(\frac{\sqrt{(R_s-R_p)^2 - b^2}}{a \sin(i)} \right) \right)\]

is evaluated, where b is the impact parameter (\(b=a \cos(i)\)).

Note

The units of the ingress duration are the same as the units of the input orbital period.

Parameters
smafloat

The semi-major axis in AU.

rpfloat

The planetary radius in Jovian radii.

rsfloat

The stellar radius in solar radii.

incfloat

The orbital inclination in degrees.

periodfloat

The orbital period.

Returns
Ingress durationfloat

The duration of in- and egress (same units as period).

PyAstronomy.pyasl.ingressDuration_Rs(sma, rprs, inc, period)

Calculate transit duration

Invokes transitDuration() after parameter transformation.

Parameters
smafloat

Semi-major axis [stellar radii]

rprsfloat

Planet-to-star radius ratio (Rp/Rs)

incfloat

Orbital inclination [deg]

periodfloat

Orbital period

Returns
Transit durationfloat

The duration of the transit (same units as period).

Estimate times of contact (T1-T4) for eccentric orbit

Numerically estimate contact points for primary and secondary eclipse.

Example

from __future__ import print_function
from PyAstronomy import pyasl

# SMA in stellar radii
sma = 5.67
# Rp/Rs
rprs = 0.15
# Orbital inclination
inc = 89.2
# Orbital period (time units are arbitrary but must be consistent)
p = 2.0
# Eccentricity
e = 0.63
# Argument of periastron (planetary orbit)
w = 155.
# Time of periastron passage
tau = 2412345.346

# Contact times for primary transit
pts = pyasl.transit_T1_T4_ell(sma, rprs, inc, p, tau, e, w, transit="p")
# Contact times for secondary transit
sts = pyasl.transit_T1_T4_ell(sma, rprs, inc, p, tau, e, w, transit="s")

print("Transit times at arbitrary epoch (N*period may added)")
print("Primary transit T1-4: ", pts)
print("Secondary trabnsit T1-4: ", sts)
print()
print("Duration of primary and secondary transit: %5.3f, %5.3f " % (pts[3]-pts[0], sts[3]-sts[0]))

Function documentation

PyAstronomy.pyasl.transit_T1_T4_ell(sma, rprs, inc, period, tau, e, w, transit='p')

Calculate first to fourth contact times for elliptical orbit.

The contact times are numerically estimated using the algorithms provided by KeplerEllipse.

Parameters
smafloat

Semi-major axis [stellar radii]

rprsfloat

Planet-to-star radius ratio (Rp/Rs)

incfloat

Orbital inclination [deg]

periodfloat

Orbital period

taufloat

Time of periastron passage

efloat

Eccentricity

wfloat

Argument of periastron

transitstring, {p, s}, optional

Determines whether contact point times for primary (p, default) or secondary (s) transit are to be calculated.

Returns
T1-T4tuple of floats

Contact times (arbitrary epoch). Note that contact time may not exist in the case of grazing transits or no transit at all. None is returned for undefined contact times.