Transit light curve routines

Contents:

The Mandel & Agol (2002) transit light curve equations.

FUNCTIONS:

occultuniform() – uniform-disk transit light curve

occultquad() – quadratic limb-darkening

occultnonlin() – full (4-parameter) nonlinear limb-darkening

occultnonlin_small() – small-planet approximation with full

nonlinear limb-darkening.

t2z() – convert dates to transiting z-parameter for circular

orbits.

REQUIREMENTS:

numpy

scipy.special

NOTES:

Certain values of p (<0.09, >0.5) cause some routines to hang; your mileage may vary. If you find out why, please let me know!

Cursory testing suggests that the Python routines contained within

are slower than the corresponding IDL code by a factor of 5-10.

For occultquad() I relied heavily on the IDL code of E. Agol and J. Eastman.

Function appellf1() comes from the mpmath compilation, and is adopted (with modification) for use herein in compliance with its BSD license (see function documentation for more details).

REFERENCE:

The main reference is that seminal work by Mandel and Agol (2002).

LICENSE:

Created by Ian Crossfield at UCLA. The code contained herein may be reused, adapted, or modified so long as proper attribution is made to the original authors.

REVISIONS:
2011-04-22 11:08 IJMC: Finished, renamed occultation functions.

Cleaned up documentation. Published to website.

2011-04-25 17:32 IJMC: Fixed bug in ellpic_bulirsch().

2012-03-09 08:38 IJMC: Several major bugs fixed, courtesy of
  1. Aigrain at Oxford U.
2012-03-20 14:12 IJMC: Fixed modeleclipse_simple based on new

format of :func:`occultuniform. `

transit.analyzetransit_general(params, NL, NP, time, data, weights=None, dopb=False, domcmc=False, gaussprior=None, ngaussprior=None, nsigma=4, maxiter=10, parinfo=None, nthread=1, nstep=2000, nwalker_factor=20, xtol=9.9999999999999998e-13)[source]

Fit transit to data, and estimate uncertainties on the fit.

Inputs :
params : sequence

A guess at the best-fit model transit parameters, to be passed to modeltransit_general().

NL : int

number of limb-darkening parameters (cf. modeltransit_general())

NP : int

number of normalizing polynomial coefficients (cf. modeltransit_general())

time : sequence

time values (e.g., BJD_TDB - BJD_0)

data : sequence

photometric values (i.e., the transit light curve) to be fit to.

weights : sequence

weights to the photometric values. If None, weights will be set equal to the inverse square of the residuals to the best-fit model. In either case, extreme outliers will be de-weighted in the fitting process. This will not change the values of the input ‘weights’.

nsigma : scalar

Residuals beyond this value of sigma-clipped standard deviations will be de-weighted.

dopb : bool

If True, run prayer-bead (residual permutation) error analysis.

domcmc : bool

If True, run Markov Chain Monte Carlo error analysis (requires EmCee)

Outputs :

(some object with useful fields)

transit.appelf1_ac(a, b1, b2, c, z1, z2, **kwargs)[source]

Analytic continuations of the Appell hypergeometric function of 2 variables.

Reference :Olsson 1964, Colavecchia et al. 2001
transit.appellf1(a, b1, b2, c, z1, z2, **kwargs)[source]

Give the Appell hypergeometric function of two variables.

Inputs :

six parameters, all scalars.

Options :

eps – scalar, machine tolerance precision. Defaults to 1e-12.

Notes :

Adapted from the mpmath module, but using the scipy (instead of mpmath) Gauss hypergeometric function speeds things up.

License :

MPMATH Copyright (c) 2005-2010 Fredrik Johansson and mpmath contributors. All rights reserved.

Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:

  1. Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer.
  2. Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in the documentation and/or other materials provided with the distribution.
  3. Neither the name of mpmath nor the names of its contributors may be used to endorse or promote products derived from this software without specific prior written permission.

THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS “AS IS” AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.

transit.ellke(k)[source]

Compute Hasting’s polynomial approximation for the complete elliptic integral of the first (ek) and second (kk) kind.

Inputs :k – scalar or Numpy array
Outputs :ek, kk
Notes :Adapted from the IDL function of the same name by J. Eastman (OSU).
transit.ellke2(k, tol=2.2204460492503131e-14, maxiter=100)[source]

Compute complete elliptic integrals of the first kind (K) and second kind (E) using the series expansions.

transit.ellpic_bulirsch(n, k, tol=2.2204460492503131e-14, maxiter=10000.0)[source]

Compute the complete elliptical integral of the third kind using the algorithm of Bulirsch (1965).

Inputs :

n – scalar or Numpy array

k– scalar or Numpy array

Notes :

Adapted from the IDL function of the same name by J. Eastman (OSU).

transit.fiteclipse(data, sv, ords, tlc, edata=None, index=None, dotransit=True, dopb=True)[source]

data: time series to fit using least-squares.

sv: state vectors (e.g., various instrumental parameters)

ords: orders to raise each sv vector to: e.g., [1, [1,2], 3]

tlc: eclipse light curve

edata: error on the data (for chisq ONLY! No weighted fits.)

index: array index to apply to data, sv, and tlc

dopb: do prayer-bead uncertainty analysis

dotransit: include tlc in the fitting; otherwise, leave it out.

transit.integral_smallplanet_nonlinear(z, p, cn, lower, upper)[source]

Return the integral in I*(z) in Eqn. 8 of Mandel & Agol (2002). – Int[I(r) 2r dr]_{z-p}^{1}, where:

Inputs :
z = scalar or array. Distance between center of star &

planet, normalized by the stellar radius.

p = scalar. Planet/star radius ratio.

cn = 4-sequence. Nonlinear limb-darkening coefficients,

e.g. from Claret 2000.

lower, upper – floats. Limits of integration in units of mu

Returns :

value of the integral at specified z.

transit.mcmc_eclipse(flux, sigma, t, func, params, tparams, stepsize, numit, nstep=1, posdef=None, holdfixed=None)[source]

MCMC for 3-parameter eclipse function with KNOWN orbit

Inputs :
flux : 1D array

Contains dependent data

sigma : 1D array

Contains standard deviation (uncertainties) of flux data

t : 1D array

Contains independent data: timing info

func : function

Function to model eclipse (e.g., transit.occultuniform())

params : parameters to be fit:
EITHER:

[T_center, depth, Fstar]

OR:

[c0, ..., c13, T_center, depth, Fstar]

params : 4 KNOWN, CONSTANT orbital parameters

[b, Rstar/a, Rp/Rstar, period]

stepsize : 1D array

Array of 1-sigma change in parameter per iteration

numit : int

Number of iterations to perform

nstep : int

Saves every “nth” step of the chain

posdef : None, ‘all’, or sequences of indices.

Which elements should be restricted to positive definite? If indices, it should be of the form (e.g.): [0, 1, 4]

holdfixed : None, or sequences of indices.

Which elements should be held fixed in the analysis? If indices, it should be of the form (e.g.): [0, 1, 4]

Returns :
allparams : 2D array

Contains all parameters at each step

bestp : 1D array

Contains best paramters as determined by lowest Chi^2

numaccept: int

Number of accepted steps

chisq: 1D array

Chi-squared value at each step

References :

Numerical Recipes, 3rd Edition (Section 15.8); Wikipedia

transit.mcmc_transit_single(flux, sigma, t, per, func, params, stepsize, numit, nstep=1, posdef=None, holdfixed=None)[source]

MCMC for 5-parameter eclipse function of transit with KNOWN period

Inputs :
flux : 1D array

Contains dependent data

sigma : 1D array

Contains standard deviation (uncertainties) of flux data

t : 1D array

Contains independent data: timing info

per : scalar

Known orbital period (same units as t)

func : function

Function to model transit (e.g., transit.occultuniform)

params : 5+N parameters to be fit
[T_center, b, Rstar/a, Rp/Rstar, Fstar] + (limb-darkening parameters?)

#[Fstar, t_center, b, v (in Rstar/day), p (Rp/Rs)]

stepsize : 1D or 2D array

if 1D: array of 1-sigma change in parameter per iteration if 2D: array of covariances for new parameters

numit : int

Number of iterations to perform

nstep : int

Saves every “nth” step of the chain

posdef : None, ‘all’, or sequences of indices.

Which elements should be restricted to positive definite? If indices, it should be of the form (e.g.): [0, 1, 4]

holdfixed : None, or sequences of indices.

Which elements should be held fixed in the analysis? If indices, it should be of the form (e.g.): [0, 1, 4]

Returns :
allparams : 2D array

Contains all parameters at each step

bestp : 1D array

Contains best paramters as determined by lowest Chi^2

numaccept: int

Number of accepted steps

chisq: 1D array

Chi-squared value at each step

References :

Numerical Recipes, 3rd Edition (Section 15.8); Wikipedia

transit.mcmc_transit_single14(flux, sigma, t, per, func, params, stepsize, numit, nstep=1, posdef=None, holdfixed=None)[source]

MCMC for 5-parameter eclipse function of transit with KNOWN period

Inputs :
flux : 1D array

Contains dependent data

sigma : 1D array

Contains standard deviation (uncertainties) of flux data

t : 1D array

Contains independent data: timing info

per : scalar

Known orbital period (same units as t)

func : function

Function to model transit (e.g., transit.occultuniform)

params : 14+5+N parameters to be fit

[c0,...,c13] + [T_center, b, Rstar/a, Rp/Rstar, Fstar] + (limb-darkening parameters?)

stepsize : 1D or 2D array

If 1D: 1-sigma change in parameter per iteration If 2D: covariance matrix for parameter changes.

numit : int

Number of iterations to perform

nstep : int

Saves every “nth” step of the chain

posdef : None, ‘all’, or sequences of indices.

Which elements should be restricted to positive definite? If indices, it should be of the form (e.g.): [0, 1, 4]

holdfixed : None, or sequences of indices.

Which elements should be held fixed in the analysis? If indices, it should be of the form (e.g.): [0, 1, 4]

Returns :
allparams : 2D array

Contains all parameters at each step

bestp : 1D array

Contains best paramters as determined by lowest Chi^2

numaccept: int

Number of accepted steps

chisq: 1D array

Chi-squared value at each step

References :

Numerical Recipes, 3rd Edition (Section 15.8); Wikipedia

transit.modeleclipse(params, func, per, t)[source]

Model an eclipse light curve of arbitrary type to a flux time series, assuming zero eccentricity and a fixed, KNOWN period.

Inputs :
params – (6-or-7)-sequence with the following:

the time of conjunction for each individual eclipse (Tc),

the impact parameter (b = a cos i/Rstar)

the stellar radius in units of orbital distance (Rstar/a),

planet-to-star radius ratio (Rp/Rstar),

eclipse depth (dimensionless),

stellar flux (F0),

orbital period (OPTIONAL!)

func – function to fit to data; presumably transit.occultuniform()

per – float.

Orbital period, OR

None, if period is included in params

t – numpy array.

Time of observations (same units as Tc and per)

transit.modeleclipse_simple(params, tparams, func, t)[source]

Model an eclipse light curve of arbitrary type to a flux time series, assuming zero eccentricity and a fixed, KNOWN orbit.

Inputs :
params – (3)-sequence with eclipse parameters to FIT:

the time of conjunction for each individual eclipse (Tc),

eclipse depth (dimensionless),

stellar flux (F0),

tparams – (4)-sequence of transit parameters to HOLD FIXED:

the impact parameter (b = a cos i/Rstar)

the stellar radius in units of orbital distance (Rstar/a),

planet-to-star radius ratio (Rp/Rstar),

orbital period (same units as Tc and t)

func – function to fit to data; presumably transit.occultuniform()

t – numpy array. Time of observations.

transit.modeleclipse_simple14(params, tparams, func, t)[source]

Model an eclipse light curve of arbitrary type to a flux time series, assuming zero eccentricity and a fixed, KNOWN orbit.

Inputs :
params – (14+3)-sequence with eclipse parameters to FIT:

the multiplicative sensitivity effects (c0, ..., c13), which affect each bit of data as (1. + c_j) * ... HOWEVER, to keep these from becoming degenerate with the overall stellar flux level, only 13 of these are free parameters: the first (c0) will always be set such that the product PROD_j(1 + c_j) = 1.

the time of conjunction for each individual eclipse (Tc),

eclipse depth (dimensionless),

stellar flux (F0),

tparams – (4)-sequence of transit parameters to HOLD FIXED:

the impact parameter (b = a cos i/Rstar)

the stellar radius in units of orbital distance (Rstar/a),

planet-to-star radius ratio (Rp/Rstar),

orbital period (same units as Tc and t)

func – function to fit to data; presumably transit.occultuniform()

t – numpy array. Time of observations.

Must either be of size (14xN), or if a 1D vector then t.reshape(14, N) must correctly reformat the data into data streams at 14 separate positions.

transit.modellightcurve(params, t, tfunc=<function occultuniform at 0x65b5030>, nlimb=0, nchan=0)[source]

Model a full planetary light curve: transit, eclipse, and (sinusoidal) phase variation. Accept independent eclipse and transit times-of-center, but otherwise assume a circular orbit (and thus symmetric transits and eclipses).

Inputs :

params – (M+10+N)-sequence with the following:

OPTIONALLY:

sensitivity variations for each of M channels (e.g., SST/MIPS). This assumes the times in ‘t’ are in the order (T_{1,0}, ... T_{1,M-1}, ... T_{2,0}, ...). The parameters affect the data multiplicatively as (1 + c_i), with the constraint that Prod_i(1+c_i) = 1.

the time of conjunction for each individual transit (T_t),

the time of conjunction for each individual eclipse (T_e),

the orbital period (P),

the impact parameter (b = a cos i/Rstar)

the stellar radius in units of orbital distance (Rstar/a),

planet-to-star radius ratio (Rp/Rstar),

stellar flux (F0),

maximum (hot-side) planet flux (Fbright),

minimum (cold-side) planet flux (Fdark),

phase curve offset (phi_0; 0 implies maximum flux near eclipse)

OPTIONALLY:

limb-darkening parameters (depending on tfunc):

EITHER:

gamma1, gamma2 – quadratic limb-darkening coefficients

OR:

c1, c2, c3, c4 – nonlinear limb-darkening coefficients

t – numpy array. Time of observations; same units as orbital

period and ephemerides. If nchan>0, t should be of shape (nchan, L), or a .ravel()ed version of that.

Options :
tfunc : model transit function

One of occultuniform(), occultnonlin_small(), occultquad(), or occultnonlin()

nlimb : int

number of limb-darkening parameters; these should be the last values of params.

nchan : int

number of photometric channel sensitivity perturbations; these should be the first ‘nchan’ values of params.

Example :

TBW

transit.modeltransit(params, func, per, t)[source]

Model a transit light curve of arbitrary type to a flux time series, assuming zero eccentricity and a fixed, KNOWN period.

Inputs :
params – (5+N)-sequence with the following:

the time of conjunction for each individual transit (Tc),

the impact parameter (b = a cos i/Rstar)

the stellar radius in units of orbital distance (Rstar/a),

planet-to-star radius ratio (Rp/Rstar),

stellar flux (F0),

the limb-darkening parameters u1 and u2:

EITHER:

gamma1, gamma2 – quadratic limb-darkening coefficients

OR:

c1, c2, c3, c4 – nonlinear limb-darkening coefficients

OR:

Nothing at all (i.e., only 5 parameters).

func – function to fit to data, e.g. transit.occultquad

per – float. Orbital period, in days.

t – numpy array. Time of observations.

transit.modeltransit14(params, func, per, t)[source]

Model a transit light curve of arbitrary type to a flux time series, assuming zero eccentricity and a fixed, KNOWN period, and assuming MIPS-type data with 14 separate sensitivity dependencies.

Inputs :

params – (14+5+N)-sequence with the following:

the multiplicative sensitivity effects (c0, ..., c13), which affect each bit of data as (1. + c_j) * ... HOWEVER, to keep these from becoming degenerate with the overall stellar flux level, only 13 of these are free parameters: the first (c0) will always be set such that the product PROD_j(1 + c_j) = 1.

the time of conjunction for each individual transit (Tc),

the impact parameter (b = a cos i/Rstar)

the stellar radius in units of orbital distance (Rstar/a),

planet-to-star radius ratio (Rp/Rstar),

stellar flux (F0),

the limb-darkening parameters u1 and u2:

EITHER:

gamma1, gamma2 – quadratic limb-darkening coefficients

OR:

c1, c2, c3, c4 – nonlinear limb-darkening coefficients

OR:

Nothing at all (i.e., only 5 parameters).

func – function to fit to data, e.g. transit.occultquad

per – float. Orbital period, in days.

t – numpy array. Time of observations.

Must either be of size (14xN), or if a 1D vector then t.reshape(14, N) must correctly reformat the data into data streams at 14 separate positions.

See ALSO:

modeltransit()

transit.modeltransit_general(params, t, NL, NP=1, smallplanet=True)[source]

Model a transit light curve of arbitrary type to a flux time series, assuming zero eccentricity.

Inputs :
params – (6 + NL + NP)-sequence with the following:

Tc, the time of conjunction for each individual transit,

P, the orbital period (in units of “t”)

i, the orbital inclination (in degrees; 90 is edge-on)

R*/a, the stellar radius in units of orbital distance,

Rp/R*, planet-to-star radius ratio,

the NP polynomial coefficients to normalize the data.

EITHER:

F0 – stellar flux _ONLY_ (set NP=1)

OR:

[p_1, p_2, ..., p_(NP)] – coefficients for polyval, to be used as: numpy.polyval([p_1, ...], t)

the NL limb-darkening parameters (cf. Claret+2011):

EITHER:

u – linear limb-darkening (set NL=1)

OR:

a, b – quadratic limb-darkening (set NL=2)

OR:

c, d – root-square limb-darkening (set NL= -2)

OR:

a1, a2, a3, a4 – nonlinear limb-darkening (set NL=4)

OR:

Nothing at all – uniform limb-darkening (set NL=0)

t – numpy array. Time of observations.

smallplanet : bool

This only matters if NL>2. If “smallplanet” is True, use occultnonlin_small(). Otherwise, use occultnonlin()

Notes :

If quadratic or linear limb-darkening (L.D.) is used, the sum of the L.D. coefficients cannot exceed 1. If they do, this routine normalizes the coefficients [g1,g2] to: g_i = g_i / (g1 + g2).

If “Rp/R*”, or “R*/a” are < 0, they will be set to zero.

If “P” < 0.01, it will be set to 0.01.

If “inc” > 90, it will be set to 90.

transit.occultnonlin(z, p0, cn)[source]

Nonlinear limb-darkening light curve; cf. Section 3 of Mandel & Agol (2002).

Inputs :

z – sequence of positional offset values

p0 – planet/star radius ratio

cn – four-sequence. nonlinear limb darkening coefficients

Example :
# Reproduce Figure 2 of Mandel & Agol (2002):
from pylab import *
import transit
z = linspace(0, 1.2, 50)
cns = vstack((zeros(4), eye(4)))
figure()
for coef in cns:
    f = transit.occultnonlin(z, 0.1, coef)
    plot(z, f)
See ALSO:

t2z(), occultnonlin_small(), occultuniform(), occultquad()

Notes :

Scipy is much faster than mpmath for computing the Beta and Gauss hypergeometric functions. However, Scipy does not have the Appell hypergeometric function – the current version is not vectorized.

transit.occultnonlin_small(z, p, cn)[source]

Nonlinear limb-darkening light curve in the small-planet approximation (section 5 of Mandel & Agol 2002).

Inputs :

z – sequence of positional offset values

p – planet/star radius ratio

cn – four-sequence nonlinear limb darkening coefficients. If

a shorter sequence is entered, the later values will be set to zero.

Note :

I had to divide the effect at the near-edge of the light curve by pi for consistency; this factor was not in Mandel & Agol, so I may have coded something incorrectly (or there was a typo).

Example :
# Reproduce Figure 2 of Mandel & Agol (2002):
from pylab import *
import transit
z = linspace(0, 1.2, 100)
cns = vstack((zeros(4), eye(4)))
figure()
for coef in cns:
    f = transit.occultnonlin_small(z, 0.1, coef)
    plot(z, f, '--')
See ALSO:

t2z()

transit.occultquad(z, p0, gamma, retall=False, verbose=False)[source]

Quadratic limb-darkening light curve; cf. Section 4 of Mandel & Agol (2002).

Inputs :

z – sequence of positional offset values

p0 – planet/star radius ratio

gamma – two-sequence.

quadratic limb darkening coefficients. (c1=c3=0; c2 = gamma[0] + 2*gamma[1], c4 = -gamma[1]). If only a single gamma is used, then you’re assuming linear limb-darkening.

Options :
retall – bool.

If True, in addition to the light curve return the uniform-disk light curve, lambda^d, and eta^d parameters. Using these quantities allows for quicker model generation with new limb-darkening coefficients – the speed boost is roughly a factor of 50. See the second example below.

Example :
# Reproduce Figure 2 of Mandel & Agol (2002):
from pylab import *
import transit
z = linspace(0, 1.2, 100)
gammavals = [[0., 0.], [1., 0.], [2., -1.]]
figure()
for gammas in gammavals:
    f = transit.occultquad(z, 0.1, gammas)
    plot(z, f)
# Calculate the same geometric transit with two different
#    sets of limb darkening coefficients:
from pylab import *
import transit
p, b = 0.1, 0.5
x = (arange(300.)/299. - 0.5)*2.
z = sqrt(x**2 + b**2)
gammas = [.25, .75]
F1, Funi, lambdad, etad = transit.occultquad(z, p, gammas, retall=True)

gammas = [.35, .55]
F2 = 1. - ((1. - gammas[0] - 2.*gammas[1])*(1. - F1) + 
   (gammas[0] + 2.*gammas[1])*(lambdad + 2./3.*(p > z)) + gammas[1]*etad) / 
   (1. - gammas[0]/3. - gammas[1]/6.)
figure()
plot(x, F1, x, F2)
legend(['F1', 'F2'])
See ALSO:

t2z(), occultnonlin_small(), occultuniform()

Notes :

In writing this I relied heavily on the occultquad IDL routine by E. Agol and J. Eastman, especially for efficient computation of elliptical integrals and for identification of several apparent typographic errors in the 2002 paper (see comments in the source code).

From some cursory testing, this routine appears about 9 times slower than the IDL version. The difference drops only slightly when using precomputed quantities (i.e., retall=True). A large portion of time is taken up in ellpic_bulirsch() and ellke(), but at least as much is taken up by this function itself. More optimization (or a C wrapper) is desired!

transit.occultuniform(z, p, complement=False)[source]

Uniform-disk transit light curve (i.e., no limb darkening).

Inputs :
z – scalar or sequence; positional offset values of planet in

units of the stellar radius.

p – scalar; planet/star radius ratio.

complement : bool

If True, return (1 - occultuniform(z, p))

See ALSO:

t2z(), occultquad(), occultnonlin_small()

transit.smallplanet_nonlinear(*arg, **kw)[source]

Placeholder for backwards compatibility with my old code. The function is now called occultnonlin_small().

transit.t2z(tt, per, inc, hjd, ars, ecc=0, longperi=0)[source]

Convert HJD (time) to transit crossing parameter z.

Inputs :

tt – scalar. transit ephemeris

per – scalar. planetary orbital period

inc – scalar. orbital inclination (in degrees)

hjd – scalar or array of times, typically heliocentric or

barycentric julian date.

ars – scalar. ratio a/Rs, orbital semimajor axis over stellar radius

ecc – scalar. orbital eccentricity.

longperi=0 scalar. longitude of periapse (in radians)

Algorithm :

At zero eccentricity, z relates to physical quantities by:

z = (a/Rs) * sqrt(sin[w*(t-t0)]**2+[cos(i)*cos(w*[t-t0])]**2)

transit.uniform(*arg, **kw)[source]

Placeholder for my old code; the new function is called occultuniform().

transit.z2dt_circular(per, inc, ars, z)[source]

Convert transit crossing parameter z to a time offset for circular orbits.

Inputs :

per – scalar. planetary orbital period

inc – scalar. orbital inclination (in degrees)

ars – scalar. ratio a/Rs, orbital semimajor axis over stellar radius

z – scalar or array; transit crossing parameter z.

Returns :
|dt| – magnitude of the time offset from transit center at

which specified z occurs.

Table Of Contents

Previous topic

Plotting and analysis tools

Next topic

Transit/eclipse prediction charts

This Page