# Planetary phase curve routines¶

Contents:

Functions for fitting phase curves.

Use of ‘errfunc’ is encouraged.

REQUIREMENTS: Transit light curve routines
phasecurves.devfunc(*arg, **kw)[source]

Generic function to give the weighted residuals on a function or functions:

INPUTS: OPTIONAL INPUTS: (fitparams, function, arg1, arg2, ... , depvar, weights) OR: (fitparams, function, arg1, arg2, ... , depvar, weights, kw) OR: (allparams, (args1, args2, ..), npars=(npar1, npar2, ...)) where allparams is an array concatenation of each functions input parameters. If the last argument is of type dict, it is assumed to be a set of keyword arguments: this will be added to resfunc’s direct keyword arguments, and will then be passed to the fitting function **kw. This is necessary for use with various fitting and sampling routines (e.g., kapteyn.kmpfit and emcee.sampler) which do not allow keyword arguments to be explicitly passed. So, we cheat! Note that any keyword arguments passed in this way will overwrite keywords of the same names passed in the standard, Pythonic, way. jointpars – list of 2-tuples. For use with multi-function calling (w/npars keyword). Setting jointpars=[(0,10), (0,20)] will always set params=params and params=params. gaussprior – list of 2-tuples (or None values), same length as “fitparams.” The i^th tuple (x_i, s_i) imposes a Gaussian prior on the i^th parameter p_i by adding ((p_i - x_i)/s_i)^2 to the total chi-squared. Here in devfunc(), we _scale_ the error-weighted deviates such that the resulting chi-squared will increase by the desired amount. uniformprior – list of 2-tuples (or ‘None’s), same length as “fitparams.” The i^th tuple (lo_i, hi_i) imposes a uniform prior on the i^th parameter p_i by requiring that it lie within the specified “high” and “low” limits. We do this (imprecisely) by multiplying the resulting deviates by 1e9 for each parameter outside its limits. ngaussprior – list of 3-tuples of Numpy arrays. Each tuple (j_ind, mu, cov) imposes a multinormal Gaussian prior on the parameters indexed by ‘j_ind’, with mean values specified by ‘mu’ and covariance matrix ‘cov.’ This is the N-dimensional generalization of the ‘gaussprior’ option described above. Here in devfunc(), we _scale_ the error-weighted deviates such that the resulting chi-squared will increase by the desired amount. For example, if parameters 0 and 3 are to be jointly constrained (w/unity means), set: jparams = np.array([0, 3]) mu = np.array([1, 1]) cov = np.array([[1, .9], [9., 1]]) ngaussprior=[[jparams, mu, cov]] # Double brackets are key!
EXAMPLE:
```from numpy import *
import phasecurves
def sinfunc(period, x): return sin(2*pi*x/period)
snr = 10
x = arange(30.)
y = sinfunc(9.5, x) + randn(len(x))/snr
guess = 8.
period = optimize.fmin(phasecurves.errfunc,guess,args=(sinfunc,x, y, ones(x.shape)*snr**2))
```
phasecurves.domodfit(profile, fdatc, wmat, xmat)[source]

Helper function for fitsin

Generates takes a matrix of variables and adds (1) a phase profile to the first row and (2) a flux conservation constraint to the last column.

phasecurves.doubleexp(params, t)[source]

Model Agol et al. 2010’s double-exponential IRAC ramp function.

params - 4- or 5-sequence, defining the function as shown below.

t - sequence. Input time (presumed: since start of observations)

Functional form:
if len(params)==4:
F’/F = 1 - p0*exp(-t/p2) - p1*exp(-t/p3)
elif len(params)>4:
F’/F = p4 * (1 - p0*exp(-t/p2) - p1*exp(-t/p3))
phasecurves.doubleexp2(params, t)[source]

Model a different double-exponential ramp function.

INPUTS: Functional_form: params : 3- or 4-sequence Parameters that define the function as shown below. t : sequence. Input time (presumed: since start of observations) if len(params)==3: F’/F = (1 - p0*exp(-t/p1)) * exp(-t/p2) elif len(params)>3: F’/F = (1 - p0*exp(-t/p1)) * exp(-t/p2) * p3
phasecurves.doubleexp214(params, t)[source]

Model a different double-exponential ramp function w/14 positional offsets.

INPUTS: Functional_form: params : sequence of length [(3 or 4) + 14] Parameters that define the function as shown below. t : sequence. Input time (presumed: since start of observations). If not of shape (14 x N), will be reshaped to that. if len(params)==3: F = (1 - p0*exp(-t/p1)) * exp(-t/p2) elif len(params)>3: F = (1 - p0*exp(-t/p1)) * exp(-t/p2) * p3 return (F * (1 + p[3/4::]).reshape(14,1) )
phasecurves.eclipse14_single(param, tparam, t, xyord=None, x=None, y=None)[source]

compute 3-parameter eclipse function of a single event, and account for 14 different possible flux offsets and X/Y positional motions.

param : parameters to be fit
[Fstar, Fplanet, t_center, c0, c1, ... c13]
tparam : parameters to be held constant (from transit)
[b, v (in Rstar/day), p (Rp/Rs)]

Input data ‘t’,’x’,’y’ must therefore be of size (14xN); if not, it will be reshaped into that.

If xyord==None or 0:
f[i, j] = (1. + p[3+i]) * [eclipse light curve]

“xyord” determines the linear order of the polynomial in x and y. If xyord==1, then: f[i,j] = (...) * (1.+ p[3+i] +p[17+i]*x + p[31+i]*y)

[note that the the individual offsets will be subjected to the constraint: param = 1./(1.+param[4:17]).prod() - 1. ]

phasecurves.eclipse_single(param, tparam, t)[source]

compute 3-parameter eclipse function of a single event

param : 3 parameters to be fit
[Fstar, Fplanet, t_center]
tparam : 3 parameters to be held constant (from transit)
[b, v (in Rstar/day), p (Rp/Rs)]
Input data ‘t’ must be of size (14xN); if not, it will be reshaped
into that.
phasecurves.errfunc(*arg, **kw)[source]
Generic function to give the chi-squared error on a generic
function or functions:
INPUTS: OPTIONAL INPUTS: (fitparams, function, arg1, arg2, ... , depvar, weights) OR: (fitparams, function, arg1, arg2, ... , depvar, weights, kw) OR: (allparams, (args1, args2, ..), npars=(npar1, npar2, ...)) where allparams is an array concatenation of each functions input parameters. If the last argument is of type dict, it is assumed to be a set of keyword arguments: this will be added to errfunc2’s direct keyword arguments, and will then be passed to the fitting function **kw. This is necessary for use with various fitting and sampling routines (e.g., kapteyn.kmpfit and emcee.sampler) which do not allow keyword arguments to be explicitly passed. So, we cheat! Note that any keyword arguments passed in this way will overwrite keywords of the same names passed in the standard, Pythonic, way. jointpars – list of 2-tuples. For use with multi-function calling (w/npars keyword). Setting jointpars=[(0,10), (0,20)] will always set params=params and params=params. gaussprior – list of 2-tuples (or None values), same length as “fitparams.” The i^th tuple (x_i, s_i) imposes a Gaussian prior on the i^th parameter p_i by adding ((p_i - x_i)/s_i)^2 to the total chi-squared. Here in devfunc(), we _scale_ the error-weighted deviates such that the resulting chi-squared will increase by the desired amount. uniformprior – list of 2-tuples (or ‘None’s), same length as “fitparams.” The i^th tuple (lo_i, hi_i) imposes a uniform prior on the i^th parameter p_i by requiring that it lie within the specified “high” and “low” limits. We do this (imprecisely) by multiplying the resulting deviates by 1e9 for each parameter outside its limits. ngaussprior – list of 3-tuples of Numpy arrays. Each tuple (j_ind, mu, cov) imposes a multinormal Gaussian prior on the parameters indexed by ‘j_ind’, with mean values specified by ‘mu’ and covariance matrix ‘cov.’ This is the N-dimensional generalization of the ‘gaussprior’ option described above. Here in devfunc(), we _scale_ the error-weighted deviates such that the resulting chi-squared will increase by the desired amount. For example, if parameters 0 and 3 are to be jointly constrained (w/unity means), set: jparams = np.array([0, 3]) mu = np.array([1, 1]) cov = np.array([[1, .9], [9., 1]]) ngaussprior=[[jparams, mu, cov]] # Double brackets are key! scaleErrors – bool If True, fit for the measurement uncertainties. In this case, the first element of ‘fitparams’ (“s”) is used to rescale the measurement uncertainties. Thus weights –> weights/s^2, and chi^2 –> 2 N log(s) + chi^2/s^2 (for N data points).
EXAMPLE:
```from numpy import *
import phasecurves
def sinfunc(period, x): return sin(2*pi*x/period)
snr = 10
x = arange(30.)
y = sinfunc(9.5, x) + randn(len(x))/snr
guess = 8.
period = optimize.fmin(phasecurves.errfunc,guess,args=(sinfunc,x, y, ones(x.shape)*snr**2))
```
phasecurves.errfunc14xymult_cfix(*arg, **kw)[source]

Generic function to give the chi-squared error on a generic function:

INPUTS:
(fitparams, arg1, arg2, ... indepvar, depvar, weights)
phasecurves.fitcurve(pdat0, fdat0, efdat0, i2, xmat, phi=30, bsind=None, func=<ufunc 'cos'>, retall=False, args=None, nsigma=inf)[source]

Decorrelate data with 14 XYs, 14 offsets, and a curve model of arbitrary amplitude.

INPUTS: phase data (14xN) flux data (14xN) error on flux data (14xN) boolean time index array (14) matrix generated for testing

phi – phi values to test bsind – (N). bootstrap indices. Indices, 0-N inclusive and with

repetitions allowed.
func – phase function, returning amplitude for an input phase
phasecurves.fithelper(param, func, pdat0, fdat0, efdat0, fdatc, wmatc, xmat, i2, i3, bsind=None, nsigma=5, retall=False)[source]

Helper funtion for fitcurve – XXX

Param – either [phi_0] or [phi_0, inclination]

phasecurves.fitsin(pdat0, fdat0, efdat0, i2, xmat, phi=30, bsind=None)[source]

Decorrelate data with 14 XYs, 14 offsets, and a sinusoid model.

INPUTS: phase data (14xN) flux data (14xN) error on flux data (14xN) boolean time index array (14) matrix generated for testing

phi – either phi values to test, or number of evenly space phi
values to test.
bsind – (N). bootstrap indices. Indices, 0-N inclusive and with
repetitions allowed.
phasecurves.lam4fit(param, x)[source]

Compute labertian phase function with a fixed period=1, and x in units of orbital cycles.

param = [DC pedestal, AC amplitude, inclination (radians)]

phasecurves.lam4fit2(param, x)[source]

Compute labertian phase function with a fixed period=1, and x in units of orbital cycles – this time with a variable phase offset.

phasecurves.lam4fit_noinc(param, inc, x)[source]

Compute labertian phase function with a fixed period=1, and x in units of orbital cycles – a variable phase offset but FIXED inclination (in radians)

param = [DC pedestal, AC amplitude, phase offset (radians)]

phasecurves.lambertian(ophase, inc=1.5707963267948966)[source]

Return a lambertian phase function with peak-to-valley amplitude unity.

INPUTS:
ophase (seq) an orbital phase (in radians). Secondary eclipse
(or ‘opposition’ for non-transiting planets) occurs at pi; Primary transit (or ‘conjuction’) occurs at 0 or 2*pi
inc (float) system inclination angle (also in radians). Edge-on
is pi/2, face-on is 0.
phasecurves.lambertian_amplitude(inc, n=5000)[source]

Return amplitude of a nominally unity-amplitude lambertian with a given inclination angle, using function ‘lambertian’. inc is in radians.

phasecurves.lambertian_mean(inc, n=5000)[source]

Return mean of a nominally unity-amplitude lambertian with a given inclination angle, using function ‘lambertian’. inc is in radians.

phasecurves.lnprobfunc(*arg, **kw)[source]

Return natural logarithm of posterior probability: i.e., -chisq/2.

Inputs are the same as for errfunc().

phasecurves.makexmat(xpos, ypos, constraint=True)[source]

Generate matrix for least-squares MIPS photometric detrending.

Generate the LSq dependent-variable matrix. The rows are: 0–constant pedestal, (ii*3+1)–dither position-dependent pedestal offset, (ii*3+2)–dither position-dependent x-position correlation, (ii*3+3)–dither position-dependent y-position correlation.

EXAMPLE:
xmat = makexmat() cleanedphot = dot(coef, xmat)

Take a 2D planet map, and convert it into an eclipse light curve.

INPUTS:
map : 2D Numpy array

[M x M] square map of planet surface brightness distribution. Note that map.sum() corresponds to the fractional eclipse depth.

z : 1D Numpy array

length-N planet crossing parameter z (i.e., distance between planet and stellar geocenters in units of stellar radii).

k : scalar

[nslice x M x M] boolean maps of planet occultation at various stages, from prepEclipseMap()

[nslice x N] boolean maps of which slices correspond to which z-indices, from prepEclipseMap()

RETURNS: lightcurve : 1D Numpy array length-N light curve, normalized to unity in eclipse.
phasecurves.mcmc14xy(z, t, x, y, sigma, params, stepsize, numit, nstep=1, sumtol=1e+20, xyord=1, prodtol=1e-06, fixc0=False)[source]
Applies Markov Chain Monte Carlo model fitting using the Metropolis-Hastings algorithm.
INPUTS: z : 1D array Contains dependent data t,x,y : 1D array Contains independent data: phase, x- and y-positions sigma : 1D array Standard deviation of dependent (y) data params : 1D array Initial guesses for parameters 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 sumtol : float Tolerance error, used as: abs(params[3:17].sum()/sumtol)**2 xyord : int “xyord” determines the linear order of the polynomial in x and y. 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 Numerical Recipes, 3rd Edition (Section 15.8); Wikipedia 2008-05-02 Kevin Stevenson, UCF kevin218@knights.ucf.edu Started converting MCMC from IDL to Python while making upgrades 2008-06-21 Kevin Stevenson Finished updating, current stable version 2008-11-17 Kevin Stevenson Updated docstring Simplified routine for AST3110 project and AST5765 class demo 2010-06-12 14:21 IJC: Modified for phasecurve test 2010-06-30 11:43 IJC: Added nstep option 2010-07-01 16:47 IJC: Added sum constraint for parameters 3-17. 2010-07-14 11:26 IJC: Added product constraints for parameters 3-17
phasecurves.mcmc_eclipse14_single(z, t, x, y, sigma, params, tparam, stepsize, numit, nstep=1, xyord=None)[source]

MCMC for 17-parameter eclipse function of a single event

INPUTS: z : 1D array Contains dependent data t,x,y : 1D array Contains independent data: phase, x- and y-positions sigma : 1D array Contains standard deviation of dependent (z) data params : 17 parameters to be fit [Fstar, Fplanet, t_center, c0, ... , c13] tparam : 3 parameters to be held constant (from transit) [b, v (in Rstar/day), p (Rp/Rs)] 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 xyord : int Highest order in x/y motions 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 Numerical Recipes, 3rd Edition (Section 15.8); Wikipedia
phasecurves.mcmc_eclipse_single(z, t, sigma, params, tparam, stepsize, numit, nstep=1)[source]

MCMC for 3-parameter eclipse function of a single event

INPUTS: z : 1D array Contains dependent data t : 1D array Contains independent data: phase, x- and y-positions sigma : 1D array Contains standard deviation of dependent (z) data params : 3 parameters to be fit [Fstar, Fplanet, t_center] tparam : 3 parameters to be held constant (from transit) [b, v (in Rstar/day), p (Rp/Rs)] 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 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 Numerical Recipes, 3rd Edition (Section 15.8); Wikipedia
phasecurves.model_fixed_param(varparam, fixedparam, fixedindex, func, *arg, **kw)[source]

Allow modeling with some parameters held constant.

INPUTS: varparam : sequence Primary parameters (which can be varied). fixedparam : sequence Secondary parameters (which should be held fixed.) fixedindex : sequence of ints Indices of parameters which should be held fixed, when passed to func : function Modeling function. Arguments Thus if param = [10, 20, 50, 70] and holdfixed = , one would set varparam = [10, 50, 70] and fixedparam = . *arg : tuple Arguments to be passed to func **kw : dict Keyword (optional) arguments to be passed to func func(param, *arg, **kw)
phasecurves.nsin(param, x)[source]

Compute phase function with unknown period, assuming n sinusoids:

p(x) = param - param[i+1] * cos(2*pi*x/param[i+2] + param[i+3]) - param[i+3+1] * cos(2*pi*x/param[i+3+2] + param[i+3+3]) - ...

phasecurves.phase2sin(param, x, absamp=True)[source]

compute double-sinusoid with a fixed period=1 offset (in radians):

p(x) = param - param * cos(2*pi*x + param) - param * cos(4*pi*x + param)

ampsamp: if True, take the absolute value of param & param

phasecurves.phaselamb14(param, x)[source]

compute phase function with a fixed period=1, assuming a sinusoid, and account for 14 different possible flux offsets.

Input data ‘x’ must therefore be of size (14xN); if not, it will be reshaped into that.

p[i,j] = (1. + param[3+j]) * (param + param*lambertian(2*pi*x[i,j]+param))

[note that the first c-parameter (individual offset) will be
constrained such that: prod(1. + param[3::]) = 1.]
phasecurves.phaselinsin(param, x)[source]

compute phase function with a linear drift and period=1:

p(x) = param+param*x - param*cos(2*pi*x + param)

phasecurves.phaselinsin2(param, phaseoffset, x)[source]

compute phase function with a linear drift and fixed period=1 and fixed phase offset (in radians):

p(x) = param+param*x - param*cos(2*pi*x + phaseoffset)

phasecurves.phasepoly14(param, x)[source]

compute phase function, assuming a polynomial, and account for 14 different possible flux offsets.

Input data ‘x’ must be in units of orbital phase, and must be of size (14xN); if not, it will be reshaped into that.

For an order-N polynomial:
p[i,j] = (1. + param[N+j]) * (numpy.polyval(param[0:N], x))
[note that the first c-parameter (individual offset) will be
constrained such that: prod(1. + param[N::]) = 1.]
phasecurves.phasesin(param, x)[source]

compute phase function with a fixed period=1, assuming a sinusoid:

p(x) = param - param * cos(2*pi*x + param)

phasecurves.phasesin14(param, x)[source]

compute phase function with a fixed period=1, assuming a sinusoid, and account for 14 different possible flux offsets.

Input data ‘x’ must therefore be of size (14xN); if not, it will be reshaped into that.

p[i,j] = (1. + param[3+j]) * (param - param*cos(2*pi*x[i,j]+param))

[note that the first c-parameter (individual offset) will be
constrained such that: prod(1. + param[3::]) = 1.]
phasecurves.phasesin14xymult(param, xyord, crossord, t, x, y)[source]

compute phase function with a fixed period=1, assuming a sinusoid, and account for 14 different possible flux offsets and X/Y positional motions.

Input data ‘t’,’x’,’y’ must therefore be of size (14xN); if not, it will be reshaped into that.

“xyord” determines the linear order of the polynomial in x and y. If xyord==1, then: f[i,j] = p - p*cos(2*pi*x[i,j]+p) + p[3+i] +p[17+i]*x + p[31+i]*y

If xyord==2, then: f[i,j] = p + p[3+i] +p[17+i]*x + p[31+i]*y +p[45+i]*x**2 + p[59+i]*y**2

• p*cos(2*pi*x[i,j]+p)
If crossord==1, then the cross-terms (x*y) will be included using
14 coefficients. crossord>1 still just uses the single cross-terms; I haven’t generalized it yet to higher orders.

[note that the the individual offsets will be subjected to the constraint: param[3::] -= param[3::].mean() ]

phasecurves.phasesin14xymult_cfix(param, xyord, crossord, t, x, y)[source]

compute phase function with a fixed period=1, assuming a sinusoid, and account for 14 different possible flux offsets and X/Y positional motions.

Input data ‘t’,’x’,’y’ must therefore be of size (14xN); if not, it will be reshaped into that.

“xyord” determines the linear order of the polynomial in x and y. If xyord==1, then: f[i,j] = (p - p*cos(2*pi*x[i,j]+p)) * (1.+ p[3+i] +p[17+i]*x + p[31+i]*y)

If xyord==2, then: f[i,j] = (p - p*cos(2*pi*x[i,j]+p)) * (1.+ p[3+i] +p[17+i]*x + p[31+i]*y)

f[i,j] = (p - p*cos(2*pi*x[i,j]+p)) * (1.+ p[3+i] +p[17+i]*x + p[31+i]*y +p[45+i]*x**2 + p[59+i]*y**2

If crossord==1, then the cross-terms (x*y) will be included using
14 coefficients. crossord>1 still just uses the single cross-terms; I haven’t generalized it yet to higher orders.

[note that the the individual offsets will be subjected to the constraint: param = 1./(1.+param[4:17]).prod() - 1. ]

phasecurves.phasesin2(param, phaseoffset, x)[source]

compute phase function with a fixed period=1 and phase offset (in radians):

p(x) = param - param * cos(2*pi*x + phaseoffset)

phasecurves.phasesinsin14(param, x)[source]

compute phase function with a fixed period=1, assuming a sinusoid and first harmonic, and account for 14 different possible flux offsets.

Input data ‘x’ must therefore be of size (14xN); if not, it will be reshaped into that.

p[i,j] = (1. + param[5+j]) * [ param - param*cos(2*pi*x[i,j]+param) + param*cos(4*pi*x[i,j]+param) ]

[note that the first c-parameter (individual offset) will be
constrained such that: prod(1. + param[5::]) = 1.]
NOTE: The magnitude of the amplitudes will always be taken; they cannot be negative.
phasecurves.phasesinsin14_2(param, x)[source]

compute phase function with a fixed period=1, assuming a sinusoid and first harmonic, and account for 14 different possible flux offsets.

Input data ‘x’ must therefore be of size (14xN); if not, it will be reshaped into that.

p[i,j] = (1. + param[5+j]) * [ param - param*cos(2*pi*x[i,j]) - param*sin(2*pi*x[i,j]) -
param*cos(4*pi*x[i,j]) - param*sin(4*pi*x[i,j]) ]
[note that the first c-parameter (individual offset) will be
constrained such that: prod(1. + param[5::]) = 1.]
phasecurves.prepEclipseMap(nslice, npix0, k, b, z=None)[source]

Set up variables necessary for (DISCRETE!) eclipse-mapping.

INPUTS: nslice : int Number of slices in model npix : int Number of pixels across the digitized planet face. k : scalar Planet/star radius ratio. b : scalar Transit impact parameter [always < (1+k) ] z : 1D Numpy array Transit crossing parameter (i.e., separation between geometric planet and stellar centers) at epochs of interest. slicemasks : Numpy array [N x M x M] boolean array, one for each map cell. timemasks : Numpy array [N x T] boolean array. Equals true for the data indices which contribute to the corresponding slice mask. cumulative_ingress_masks, cumulative_egress_masks : 3D Numpy arrays [nslice x M x M] boolean maps of planet occultation at various stages, from prepEclipseMap() ingress_zmasks, egress_zmasks : 2D Numpy arrays [nslice x N] boolean maps of which slices correspond to which z-indices, from prepEclipseMap()
phasecurves.prolatesize(f, phase)[source]

The projected size of a prolate ellipsoidal planet, viewed edge-on, at a given orbital phase.

INPUTS:
f : scalar

Ratio of large and small axes of the ellipsoid

phase : scalar or NumPy array

OUTPUTS: a_scale Value by which to scale up the area of a circle (f=1). Vickers 1996: http://dx.doi.org/10.1016/0032-5910(95)03049-2
phasecurves.ramp10(params, phase, args={'guess': [1, 0.2], 'n': 2})[source]

Model Ramp Eq. 10 from Stevenson et al. (2011).

params: 2-sequence
parameters that define the function, as shown below.
phase: NumPy array.
Orbital phase (or more generally, ‘time’)
Functional form:
ramp = r * (1. + r * (phase - 0.5))
phasecurves.ramp11(params, phase, args={'guess': [1, 0.14, -1.9], 'n': 3})[source]

Model Ramp Eq. 11 from Stevenson et al. (2011).

params: 3-sequence
parameters that define the function, as shown below.
phase: NumPy array.
Orbital phase (or more generally, ‘time’)
Functional form:
ramp = r * (1. + r * (phase - 0.5) + r * (phase - 0.5)**2)
phasecurves.ramp2n(params, phase, args={'guess': [1, 26.6, 7.8], 'n': 3})[source]

Model Ramp Eq. 2 (negative) from Stevenson et al. (2011).

params: 3-sequence
parameters that define the function, as shown below.
phase: NumPy array.
Orbital phase (or more generally, ‘time’)
Functional form:
ramp = r * (1. - np.exp(-r*phase + r))
phasecurves.ramp2p(params, phase, args={'guess': [1, -0.16, 4.2], 'n': 3})[source]

Model Ramp Eq. 2 (positive) from Stevenson et al. (2011).

params: 3-sequence
parameters that define the function, as shown below.
phase: NumPy array.
Orbital phase (or more generally, ‘time’)
Functional form:
ramp = r * (1. + np.exp(-r*phase + r))
phasecurves.ramp3n(params, phase, args={'guess': [1, 141, 57.7, 0.123], 'n': 4})[source]

Model Ramp Eq. 3 (negative) from Stevenson et al. (2011).

params: 4-sequence
parameters that define the function, as shown below.
phase: NumPy array.
Orbital phase (or more generally, ‘time’)
Functional form:
ramp = r * (1. - np.exp(-r*t + r) + r * (t - 0.5))
phasecurves.ramp3p(params, phase, args={'guess': [1, -0.16, 4.2, 0.1], 'n': 4})[source]

Model Ramp Eq. 3 (positive) from Stevenson et al. (2011).

params: 4-sequence
parameters that define the function, as shown below.
phase: NumPy array.
Orbital phase (or more generally, ‘time’)
Functional form:
ramp = r * (1. + np.exp(-r*t + r) + r * (t - 0.5))
phasecurves.ramp4n(params, phase, args={'guess': [1, -0.00037, -0.94, 0.087, -1.08], 'n': 5})[source]

Model Ramp Eq. 4 (negative) from Stevenson et al. (2011).

params: 5-sequence
parameters that define the function, as shown below.
phase: NumPy array.
Orbital phase (or more generally, ‘time’)
Functional form:
ramp = r * (1. - np.exp(-r*phase + r) + r * (phase - 0.5) + r * (phase - 0.5)**2)
phasecurves.ramp4p(params, phase, args={'guess': [1, -0.068, 2.33, 0.933, -20.5], 'n': 5})[source]

Model Ramp Eq. 4 (positive) from Stevenson et al. (2011).

params: 5-sequence
parameters that define the function, as shown below.
phase: NumPy array.
Orbital phase (or more generally, ‘time’)
Functional form:
ramp = r * (1. + np.exp(-r*phase + r) + r * (phase - 0.5) + r * (phase - 0.5)**2)
phasecurves.ramp5n(params, phase, args={'guess': [1.0, 20, 83, 8.1, -0.1], 'n': 5})[source]

Model Ramp Eq. 5 (negative) from Stevenson et al. (2011).

params: 5-sequence
parameters that define the function, as shown below.
phase: NumPy array.
Orbital phase (or more generally, ‘time’)
Functional form:
ramp = r * (1. - np.exp(-r*phase + r) - np.exp(-r*phase + r))
phasecurves.ramp5p(params, phase, args={'guess': [1, -0.32, 2, -0.08, 2], 'n': 5})[source]

Model Ramp Eq. 5 (positive) from Stevenson et al. (2011).

params: 5-sequence
parameters that define the function, as shown below.
phase: NumPy array.
Orbital phase (or more generally, ‘time’)
Functional form:
ramp = r * (1. + np.exp(-r*phase + r) + np.exp(-r*phase + r))
phasecurves.ramp6(params, phase, args={'guess': [1, 0.053, 0.004, 0.4], 'n': 4})[source]

Model Ramp Eq. 6 from Stevenson et al. (2011).

params: 4-sequence
parameters that define the function, as shown below.
phase: NumPy array.
Orbital phase (or more generally, ‘time’)
Functional form:
ramp = r * (1. + r * (phase - 0.5) + r * np.log(phase - r))
phasecurves.ramp7(params, phase, args={'guess': [1, 0.034, 0.35, 0.005, 0.35], 'n': 5})[source]

Model Ramp Eq. 7 from Stevenson et al. (2011).

params: 5-sequence
parameters that define the function, as shown below.
phase: NumPy array.
Orbital phase (or more generally, ‘time’)
Functional form:
ramp = r * (1. + r * (phase - 0.5) + r * (phase - 0.5)**2 + r * np.log(phase - r))
phasecurves.ramp8(params, phase, args={'guess': [1, 0.0096, 0.35, 0.00053], 'n': 4})[source]

Model Ramp Eq. 8 from Stevenson et al. (2011).

params: 4-sequence
parameters that define the function, as shown below.
phase: NumPy array.
Orbital phase (or more generally, ‘time’)
Functional form:
ramp = r * (1. + r * np.log(phase - r) + r * np.log(phase - r)**2)
phasecurves.ramp9(params, phase, args={'guess': [1, 0.003, 0.6, 0.009, 0.35, 0.0004], 'n': 6})[source]

Model Ramp Eq. 9 from Stevenson et al. (2011).

params: 6-sequence
parameters that define the function, as shown below.
phase: NumPy array.
Orbital phase (or more generally, ‘time’)
Functional form:
ramp = r * (1. + r * (phase - 0.5) + r * (phase - 0.5)**2 + r * np.log(phase - r) + r * np.log(phase - r)**2)
phasecurves.resfunc(*arg, **kw)[source]

Generic function to give the error-weighted deviates on a function or functions:

INPUTS: OPTIONAL INPUTS: (fitparams, function, arg1, arg2, ... depvar, errs) OR: (allparams, (args1, args2, ..), npars=(npar1, npar2, ...)) where allparams is an array concatenation of each functions input parameters. jointpars – list of 2-tuples. For use with multi-function calling (w/npars keyword). Setting jointpars=[(0,10), (0,20)] will always set params=params and params=params.

EXAMPLE:

resfunc()
phasecurves.rotmod(param, airmass, rotang, phase=None)[source]

Model the Bean & Seifert rotation angle effect.

INPUTS: param : 3- or 4-sequence airmass : NumPy array rotang : NumPy array Instrument rotator angle, in radians. phase : NumPy array, or None Orbital phase (or time, etc.) param + param * airmass * np.cos(param + rotang) + param * (phase - phase.mean())
phasecurves.sin2_errs(params, eparams, nphi=10000.0, ntrials=10000.0)[source]

Estimate the uncertainties from a double-sinusoid fit.

FUNCTION: p - p*cos(phi + p) + p*cos(2*phi + p) params : 5-sequence parameters for the function, as defined immediately above. eparams : 5-sequence 1-sigma uncertainties on parameters ntrials : float number of Monte Carlo trials to run nphi : float number of points in phase curve (0-1, inclusive) (visibilities, peak_offset (rad), trough_offset (rad), true_vals) phasesinsin14()
phasecurves.singleexp(params, t)[source]

Model a simple, single-exponential function.

params: 2- or 3-sequence, defining the function as shown below.

t: sequence. Input time (presumed: since start of observations)

Functional form:
if len(params)==2:
F’/F = 1 - p0*exp(-t/p1)
elif len(params)>2:
F’/F = p2 * (1 - p0*exp(-t/p1))
phasecurves.singleexp14(params, t)[source]

Model a simple, single-exponential function.

params: 16- or 17-sequence, defining the function as shown below.

t: sequence. Input time (presumed: since start of observations)

Functional form:
if len(params)==2:
F’/F = 1 - p0*exp(-t/p1)
elif len(params)>2:

F’/F = p2 * (1 - p0*exp(-t/p1))

... with fourteen additional sensitivity parameters.

phasecurves.slicemodel(param, xi)[source]

Compute a slice model via Cowan & Agol (2008).

xi is from 0 to 2*pi

phasecurves.subfit_kw(params, input, i0, i1)[source]

Parse keyword inputs (for errfunc(), etc.) for multiple concatenated inputs (e.g., with input key ‘npars’ set).

INPUTS: params : 1D NumPy array of parameters input : dict of keywords i0 : int, first index of current parameters i1 : int, last index of current parameters (e.g.: params[i0:i1])
phasecurves.transit_single(param, t)[source]

compute 6+L-parameter eclipse function of a single event

velocity, impact parameter, stellar flux, planet/star radius ratio, time of center transit, period, limb darkening

param : parameters to be fit:

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

Up to two additional parameter can be concatenated onto the end, respectively representing linear and quadratic limb-darkening.

Input data ‘t’ must be of size (14xN); if not, it will be reshaped
into that.
phasecurves.tsin(param, x)[source]

Compute phase function with unknown period, assuming a sinusoid:

p(x) = param - param * cos(2*pi*x/param + param)

phasecurves.unwrap_joint_params(params, jfw_indices)[source]

Unwrap parameters that are jointly constrained.

INPUTS: params – 1D NumPy array The P non-redundant parameters to the input function – i.e., any parameters which are to be jointly fit (i.e., both held to the same, floating, value) are included only once. jfw_indices : sequence of scalars and sequences A length-P sequence of scalars and sequences. Each element jfw_indices[i] indicates the indices in the unwrapped set of parameters that will be assigned the value of params[i]. The final value of jfw_indices should be an integer equal to the length of the final set of unwrapped parameters. ```import tools import numpy as np npts = 100 snr = 20. params = [1, 0.5, 1, 1] x = np.linspace(0, 1, npts) y = np.polyval(params, x) + np.random.randn(npts)/snr jointpars = [(0, 2), (0, 3)] joint_guess = np.array([1, 0.5]) jfw_indices = [[0, 2, 3], , 4] full_params = tools.unwrap_joint_params(joint_guess, jfw_indices) ``` wrap_joint_params()
INPUTS:
visitcoef : 1D NumPy array

offsets from unity for each of N HST visits

masks : 2D NumPy array, N x M

Set of boolean masks for each of N visits (not orbit!), assuming M total observations.

NOTES: Note that visitcoef will be set so that the quantity (1. + visitcoef).prod() always equals unity.
phasecurves.wrap_joint_params(params, jointpars)[source]

Wrap parameters that are jointly constrained.

INPUTS: params – 1D NumPy array All parameters, some of which may be jointly constrained. jointpars – list of 2-tuples. For use with multi-function calling (w/npars keyword). Setting jointpars=[(0,10), (0,20)] will always set params=params and params=params. ```import tools import numpy as np npts = 100 snr = 20. params = [1, 0.5, 1, 1] x = np.linspace(0, 1, npts) y = np.polyval(params, x) + np.random.randn(npts)/snr jointpars = [(0, 2), (0, 3)] all_params, joint_indices = tools.wrap_joint_params(full_params, jointpars) wrapped_params = tools.unwrap_joint_params(all_params, joint_indices) ``` unwrap_joint_params()

#### Previous topic

Gaussian Processes

Planetary maps