Below is the complete description of the kima and pykima APIs.

kima classes and modules

Class RVData

Variables

: bool indicator_correlations
: int number_indicators
: vector< string > indicator_names
: string datafile
: vector< string > datafiles
: string dataunits
: int dataskip
: bool datamulti
: int number_instruments
: double M0_epoch

docs for M0_epoch

Methods

int N()

Get the total number of RV points.

int Ninstruments()

Get the number of instruments.

double get_RV_max()

Get the maximum RV.

double get_RV_min()

Get the mininum RV.

double get_RV_span()

Get the RV span.

double get_RV_std()

Get the standard deviation of the RVs.

double get_RV_var()

Get the variance of the RVs.

const vector< vector< double > > & get_actind()

Get the array of activity indictators.

double get_adjusted_RV_span()

Get the RV span, adjusted for multiple instruments.

double get_adjusted_RV_std()

Get the RV standard deviation, adjusted for multiple instruments.

double get_adjusted_RV_var()

Get the RV variance, adjusted for multiple instruments.

get_instance()
double get_max_RV_span()

Get the maximum RV span.

const vector< int > & get_obsi()

Get the array of instrument identifiers.

const vector< double > & get_sig()

Get the array of errors.

const vector< double > & get_sig2()
const vector< double > & get_t()

Get the array of times.

double get_t_max()

Get the maximum (ending) time.

double get_t_middle()

Get the middle time.

double get_t_min()

Get the mininum (starting) time.

double get_t_span()
double get_timespan()

Get the timespan.

int get_trend_magnitude(int degree)

Order of magnitude of trend coefficient (of degree) given the data.

const vector< double > & get_y()

Get the array of RVs.

const vector< double > & get_y2()
double get_y2_max()

Get the maximum y2.

double get_y2_min()

Get the mininum y2.

double get_y2_span()

Get the y2 span.

void load(const string filename, const string units, int skip=2, const vector< string > &indicators=vector< string >())

Load RV data from a file.

void load_multi(const string filename, const string units, int skip=2)

Load RV data from a multi-instrument file.

double topslope()

Get the maximum slope allowed by the data.

Class RVmodel

Priors

: ContinuousDistribution Cprior

Prior for the systemic velocity.

: ContinuousDistribution Jprior

Prior for the extra white noise (jitter).

: ContinuousDistribution slope_prior

Prior for the slope.

: ContinuousDistribution quadr_prior

Prior for the quadratic coefficient of the trend.

: ContinuousDistribution cubic_prior

Prior for the cubic coefficient of the trend.

: ContinuousDistribution offsets_prior

(Common) prior for the between-instruments offsets.

: vector < ContinuousDistribution > individual_offset_prior
: ContinuousDistribution betaprior

no doc.

: ContinuousDistribution sigmaMA_prior

no doc.

: ContinuousDistribution tauMA_prior

no doc.

: ContinuousDistribution log_eta1_prior

Prior for the logarithm of $\eta_1$, the GP "amplitude".

: ContinuousDistribution eta2_prior

Prior for $\eta_2$, the GP correlation timescale.

: ContinuousDistribution eta3_prior

Prior for $\eta_3$, the GP period.

: ContinuousDistribution log_eta4_prior

Prior for the logarithm of $\eta_4$, the recurrence timescale.

: ContinuousDistribution alpha_prior

Prior for the Rational Quadratic shape parameter.

: ContinuousDistribution eta5_prior

Prior for the "amplitude" of the cosine term of the QPC kernel.

: vector < ContinuousDistribution > KO_Pprior

Prior for the KO orbital period(s)

: vector < ContinuousDistribution > KO_Kprior

Prior for the KO semi-amplitude(s)

: vector < ContinuousDistribution > KO_eprior

Prior for the KO eccentricity(ies)

: vector < ContinuousDistribution > KO_phiprior

Prior for the KO mean anomaly(ies)

: vector < ContinuousDistribution > KO_wprior

Prior for the KO argument(s) of pericenter.

: ContinuousDistribution nu_prior

Prior for the degrees of freedom $\nu$ of the Student t likelihood.

Methods

void add_known_object()
void calculate_C()

Fill the GP covariance matrix.

void calculate_mu()

Calculate the full RV model.

std::string description()
void from_prior(DNest4::RNG &rng)

Generate a point from the prior.

Data & get_data()
void initialise()
int is_stable()
double log_likelihood()
class T make_prior(Args &&... args)

Assign a prior distribution.

double perturb(DNest4::RNG &rng)

Do Metropolis-Hastings proposals.

void print(std::ostream &out)
void remove_known_object()
void save_setup()
void setPriors()
Class RVFWHMmodel

Priors

: ContinuousDistribution Vprior

Prior for the systemic velocity.

: ContinuousDistribution C2prior
: ContinuousDistribution Jprior

Prior for the extra white noise (jitter) in the RVs.

: ContinuousDistribution J2prior

Prior for the extra white noise (jitter) in the FWHM.

: ContinuousDistribution slope_prior

Prior for the slope.

: ContinuousDistribution quadr_prior

Prior for the quadratic coefficient of the trend.

: ContinuousDistribution cubic_prior

Prior for the cubic coefficient of the trend.

: ContinuousDistribution offsets_prior

(Common) prior for the between-instruments RV offsets

: ContinuousDistribution offsets2_prior

(Common) prior for the between-instruments FWHM offsets

: ContinuousDistribution betaprior

no doc.

: ContinuousDistribution sigmaMA_prior

no doc.

: ContinuousDistribution tauMA_prior

no doc.

: ContinuousDistribution eta1_1_prior

Prior for $\eta_1$, the GP "amplitude".

: ContinuousDistribution eta2_1_prior

Prior for $\eta_2$, the GP correlation timescale.

: ContinuousDistribution eta3_1_prior

Prior for $\eta_3$, the GP period.

: ContinuousDistribution eta4_1_prior

Prior for $\eta_4$, the recurrence timescale.

: ContinuousDistribution alpha_1_prior

Prior for the Rational Quadratic shape parameter.

: ContinuousDistribution eta5_1_prior

Prior for the "amplitude" of the cosine term of the QPC kernel.

: ContinuousDistribution eta1_2_prior

Same as $\eta_1$ but for the FWHM.

: ContinuousDistribution eta2_2_prior

Same as $\eta_2$ but for the FWHM.

: ContinuousDistribution eta3_2_prior

Same as $\eta_3$ but for the FWHM.

: ContinuousDistribution eta4_2_prior

Same as $\eta_4$ but for the FWHM.

: ContinuousDistribution eta5_2_prior

Same as $\eta_5$ but for the FWHM.

: vector < ContinuousDistribution > KO_Pprior

Prior for the KO orbital period(s)

: vector < ContinuousDistribution > KO_Kprior

Prior for the KO semi-amplitude(s)

: vector < ContinuousDistribution > KO_eprior

Prior for the KO eccentricity(ies)

: vector < ContinuousDistribution > KO_phiprior

Prior for the KO mean anomaly(ies)

: vector < ContinuousDistribution > KO_wprior

Prior for the KO argument(s) of pericenter.

: ContinuousDistribution nu_prior

Prior for the degrees of freedom $\nu$ of the Student t likelihood.

Variables

: bool share_eta2

Whether $\eta_2$ is shared between RVs and FWHM.

: bool share_eta3

Whether $\eta_3$ is shared between RVs and FWHM.

: bool share_eta4

Whether $\eta_4$ is shared between RVs and FWHM.

: bool share_eta5

Whether $\eta_5$ is shared between RVs and FWHM.

Methods

void add_known_object()
void calculate_C_1()

Fill the GP covariance matrix.

void calculate_C_2()

Fill the GP covariance matrix.

void calculate_mu()

Calculate the full RV model.

void calculate_mu_2()

Calculate the FWHM model.

std::string description()
double ecc_anomaly(double time, double prd, double ecc, double peri_pass)
double eps3(double e, double M, double x)
void from_prior(DNest4::RNG &rng)

Generate a point from the prior.

get_data()
void initialise()
int is_stable()
double keplerstart3(double e, double M)
double log_likelihood()
class T make_prior(Args &&... args)

Assign a prior distribution.

double perturb(DNest4::RNG &rng)

Do Metropolis-Hastings proposals.

void print(std::ostream &out)
void remove_known_object()
void save_setup()
void setPriors()
double true_anomaly(double time, double prd, double ecc, double peri_pass)


pykima modules

Module pykima.results

This module defines the KimaResults class to hold results from a run.

Classes

class KimaResults (options, data_file=None, save_plots=False, return_figs=True, verbose=False, hyperpriors=None, trend=None, GPmodel=None, posterior_samples_file='posterior_sample.txt')

A class to hold, analyse, and display the results from kima

Static methods

def load(filename=None, diagnostic=False)

Load a KimaResults object from the current directory, a pickle file, or a zip file.

Instance variables

var data_properties
var eta1
var eta2
var eta3
var eta4
var instruments
var parameter_priors

A list of priors which can be indexed using self.indices

var ratios
var star

Methods

def apply_cuts_period(self, samples, pmin=None, pmax=None, return_mask=False)

apply cuts in orbital period

def burst_model(self, sample, t, v=None)

For models with multiple instruments, this function "bursts" the computed RV into n_instruments individual arrays. This is mostly useful for plotting the RV model together with the original data.

Arguments

sample : ndarray
One posterior sample, with shape (npar,)
t : ndarray
Times at which to evaluate the model
v : ndarray
Pre-computed RVs, if None we call self.eval_model
def corner_planet_parameters(self, fig=None, pmin=None, pmax=None)

Corner plot of the posterior samples for the planet parameters

def eval_model(self, sample, t=None, include_planets=True, include_known_object=True, include_trend=True, single_planet=None, except_planet=None)

Evaluate the deterministic part of the model at one posterior sample. If t is None, use the observed times. Instrument offsets are only added if t is None, but the systemic velocity is always added. To evaluate at all posterior samples, consider using np.apply_along_axis(self.eval_model, 1, self.posterior_sample)

Note: this function does not evaluate the GP component of the model.

Arguments

sample : ndarray
One posterior sample, with shape (npar,)
t : ndarray
Times at which to evaluate the model, or None to use observed times
include_planets : bool, default True
Whether to include the contribution from the planets
include_known_object : bool, default True
Whether to include the contribution from the known object planets
include_trend : bool, default True
Whether to include the contribution from the trend
single_planet : int
Index of a single planet to include in the model, starting at 1. Use positive values (1, 2, …) for the Np planets and negative values (-1, -2, …) for the known object planets.
except_planet : int or list
Index (or list of indices) of a single planet to exclude from the model, starting at 1. Use positive values (1, 2, …) for the Np planets and negative values (-1, -2, …) for the known object planets.
def from_prior(self, n=1)

Generate n samples from the priors for all parameters.

def full_model(self, sample, t=None, **kwargs)

Evaluate the full model at one posterior sample, including the GP. If t is None, use the observed times. Instrument offsets are only added if t is None, but the systemic velocity is always added. To evaluate at all posterior samples, consider using np.apply_along_axis(self.full_model, 1, self.posterior_sample)

Arguments

sample : ndarray
One posterior sample, with shape (npar,)
t : ndarray (optional)
Times at which to evaluate the model, or None to use observed times
**kwargs
Keyword arguments passed directly to eval_model
def get_marginals(self)

Get the marginal posteriors from the posterior_sample matrix. They go into self.T, self.A, self.E, etc

def get_medians(self)

return the median values of all the parameters

def get_sorted_planet_samples(self, full=True)
def hist_MA(self)

Plot the histogram of the posterior for the MA parameters

def hist_correlations(self)

Plot the histogram of the posterior for the activity correlations

def hist_jitter(self, show_prior=False, **kwargs)

Plot the histogram of the posterior for the additional white noise

def hist_nu(self, show_prior=False, **kwargs)

Plot the histogram of the posterior for the Student-t degrees of freedom

def hist_trend(self, per_year=True, show_prior=False, ax=None)

Plot the histogram of the posterior for the coefficients of the trend

def hist_vsys(self, show_offsets=True, specific=None, show_prior=False, **kwargs)

Plot the histogram of the posterior for the systemic velocity and for the between-instrument offsets (if show_offsets is True and the model has multiple instruments). If specific is not None, it should be a tuple with the name of the datafiles for two instruments (matching self.data_file). In that case, this function works out the RV offset between the specific[0] and specific[1] instruments.

def log_posterior(self, sample, separate=False)
def log_prior(self, sample)
def make_plot1(self, **kwargs)

Plot the histogram of the posterior for Np

def make_plot2(self, nbins=100, bins=None, plims=None, logx=True, density=False, kde=False, kde_bw=None, show_peaks=False, show_prior=False, show_year=True, show_timespan=True, separate_colors=False, **kwargs)

Plot the histogram (or the kde) of the posterior for the orbital period(s). Optionally provide the number of histogram bins, the bins themselves, the limits in orbital period, or the kde bandwidth. If both kde and show_peaks are true, the routine attempts to plot the most prominent peaks in the posterior.

def make_plot3(self, mask=None, include_known_object=False, points=True, gridsize=50, **kwargs)

Plot the 2d histograms of the posteriors for semi-amplitude and orbital period and eccentricity and orbital period. If points is True, plot each posterior sample, else plot hexbins

def make_plot4(self, *args, **kwargs)

Plot histograms for the GP hyperparameters. If Np is not None, highlight the samples with Np Keplerians.

def make_plot5(self, include_jitters=False, show=True, ranges=None)

Corner plot for the GP hyperparameters

def make_plots(self, options, save_plots=False)
def map_sample(self, Np=None, mask=None, printit=True, cache=True)
def maximum_likelihood_sample(self, from_posterior=False, Np=None, printit=True, mask=None)

Get the maximum likelihood sample. By default, this is the highest likelihood sample found by DNest4. If from_posterior is True, this returns instead the highest likelihood sample from those that represent the posterior. The latter may change, due to random choices, between different calls to "showresults". If Np is given, select only samples with that number of planets.

def median_sample(self, Np=None, printit=True)

Get the median posterior sample. If Np is given, select only from posterior samples with that number of planets.

def phase_plot(self, sample, highlight=None, **kwargs)

Plot the phase curves given the solution in sample

def planet_model(self, sample, t=None, include_known_object=True, single_planet=None, except_planet=None)

Evaluate the planet part of the model at one posterior sample. If t is None, use the observed times. To evaluate at all posterior samples, consider using np.apply_along_axis(self.planet_model, 1, self.posterior_sample)

Note: this function does not evaluate the GP component of the model nor the systemic velocity and instrument offsets.

Arguments

sample : ndarray
One posterior sample, with shape (npar,)
t : ndarray
Times at which to evaluate the model, or None to use observed times
include_known_object : bool, default True
Whether to include the contribution from the known object planets
single_planet : int
Index of a single planet to include in the model, starting at 1. Use positive values (1, 2, …) for the Np planets and negative values (-1, -2, …) for the known object planets.
except_planet : int or list
Index (or list of indices) of a single planet to exclude from the model, starting at 1. Use positive values (1, 2, …) for the Np planets and negative values (-1, -2, …) for the known object planets.
def plot1(self, **kwargs)

Plot the histogram of the posterior for Np

def plot2(self, nbins=100, bins=None, plims=None, logx=True, density=False, kde=False, kde_bw=None, show_peaks=False, show_prior=False, show_year=True, show_timespan=True, separate_colors=False, **kwargs)

Plot the histogram (or the kde) of the posterior for the orbital period(s). Optionally provide the number of histogram bins, the bins themselves, the limits in orbital period, or the kde bandwidth. If both kde and show_peaks are true, the routine attempts to plot the most prominent peaks in the posterior.

def plot3(self, mask=None, include_known_object=False, points=True, gridsize=50, **kwargs)

Plot the 2d histograms of the posteriors for semi-amplitude and orbital period and eccentricity and orbital period. If points is True, plot each posterior sample, else plot hexbins

def plot4(self, *args, **kwargs)

Plot histograms for the GP hyperparameters. If Np is not None, highlight the samples with Np Keplerians.

def plot5(self, include_jitters=False, show=True, ranges=None)

Corner plot for the GP hyperparameters

def plot6(self, **kwargs)

Display the RV data together with curves from the posterior predictive. A total of ncurves random samples are chosen, and the Keplerian curves are calculated covering 100 + over% of the data timespan. If the model has a GP component, the prediction is calculated using the GP hyperparameters for each of the random samples.

def plot_posterior_np(self, **kwargs)

Plot the histogram of the posterior for Np

def plot_random_planets(self, **kwargs)

Display the RV data together with curves from the posterior predictive. A total of ncurves random samples are chosen, and the Keplerian curves are calculated covering 100 + over% of the data timespan. If the model has a GP component, the prediction is calculated using the GP hyperparameters for each of the random samples.

def print_sample(self, p, star_mass=1.0, show_a=False, show_m=False, mass_units='mjup', squeeze=False)
def residual_rms(self, sample, weighted=True, printit=True)
def residuals(self, sample, full=False)
def save_pickle(self, filename, verbose=True)

Pickle this KimaResults object into a file.

def save_zip(self, filename, verbose=True)

Save this KimaResults object and the text files into a zip.

def show_kima_setup(self)
def simulate_from_sample(self, sample, times, add_noise=True, errors=True, append_to_file=False)
def stochastic_model(self, sample, t=None, return_std=False, derivative=False, **kwargs)

Evaluate the stochastic part of the model (GP) at one posterior sample. If t is None, use the observed times. Instrument offsets are only added if t is None, but the systemic velocity is always added. To evaluate at all posterior samples, consider using np.apply_along_axis(self.stochastic_model, 1, self.posterior_sample)

Arguments

sample : ndarray
One posterior sample, with shape (npar,)
t : ndarray (optional)
Times at which to evaluate the model, or None to use observed times
return_std : bool (optional, default False)
Whether to return the st.d. of the predictive
derivative : bool (optional, default False)
Return the first time derivative of the GP prediction instead
Module pykima.analysis

Functions

def column_dynamic_ranges(results)

Return the range of each column in the posterior file

def columns_with_dynamic_range(results)

Return the columns in the posterior file which vary

def detection_limits(results, star_mass=1.0, Np=None, bins=200, plot=True, semi_amplitude=False, semi_major_axis=False, logX=True, sorted_samples=True, return_mask=False, remove_nan=True, add_jitter=False, show_prior=False, show_eccentricity=False, smooth=False, smooth_degree=3)

Calculate detection limits using samples with more than Np planets. By default, this function uses the value of Np which passes the posterior probability threshold.

Arguments

star_mass : float or tuple
Stellar mass and optionally its uncertainty [in solar masses].
Np : int
Consider only posterior samples with more than Np planets.
bins : int
Number of bins at which to calculate the detection limits. The period ranges from the minimum to the maximum orbital period in the posterior.
plot : bool (default True)
Whether to plot the detection limits
semi_amplitude : bool (default False)
Show the detection limits for semi-amplitude, instead of planet mass
semi_major_axis : bool (default False)
Show semi-major axis in the x axis, instead of orbital period
logX : bool (default True)
Should X-axis bins be logarithmic ?
sorted_samples : bool
undoc
return_mask : bool
undoc
remove_nan : bool
remove bins with no counts (no posterior samples)
add_jitter …
show_prior …
show_prior : bool (default False)
If true, color the posterior samples according to orbital eccentricity
smooth : bool (default False)
Smooth the binned maximum with a polynomial
smooth_degree : int (default 3)
Degree of the polynomial used for smoothing

Returns

P, K, E, M : ndarray
Orbital periods, semi-amplitudes, eccentricities, and planet masses used in the calculation of the detection limits. These correspond to all the posterior samples with more than Np planets
s : DLresult, namedtuple
Detection limits result, with attributes max and bins. The max array is in units of Earth masses, bins is in days.
def find_outliers(results, sample, threshold=10, verbose=False)

Estimate which observations are outliers, for a model with a Student t likelihood. This function first calculates the residuals given the parameters in sample. Then it computes the relative probability of each residual point given a Student-t (Td) and a Gaussian (Gd) likelihoods. If the probability Td is larger than Gd (by a factor of threshold), the point is flagged as an outlier. The function returns an "outlier mask".

def full_model_table(res, sample, instruments=None, star_mass=1.0)
def most_probable_np(results)

Return the value of Np with the highest posterior probability.

Arguments

results : KimaResults
A results instance
def parameter_clusters(results, n_clusters=None, method='KMeans', include_ecc=True, scale=False, plot=True, downsample=1, **kwargs)
def passes_threshold_np(results, threshold=150)

Return the value of Np supported by the data, considering a posterior ratio threshold (default 150).

Arguments

results : KimaResults
A results instance
threshold : float
Posterior ratio threshold, default 150.
def planet_parameters(results, star_mass=1.0, sample=None, printit=True)
def sort_planet_samples(res, planet_samples, byP=True, byK=False)
Module pykima.utils

Small (but sometimes important) utility functions

Functions

def apply_argsort(arr1, arr2, axis=-1)

Apply arr1.argsort() on arr2, along axis.

def chdir(dir)

A simple context manager to switch directories temporarily

def clipped_mean(arr, min, max)

Mean of arr between min and max

def clipped_std(arr, min, max)

Standard deviation of arr between min and max

def date_to_jd(year, month, day)

Convert a date (year, month, day) to Julian Day.

Algorithm from 'Practical Astronomy with your Calculator or Spreadsheet', 4th ed., Duffet-Smith and Zwart, 2011.

Parameters

year : int
Year as integer. Years preceding 1 A.D. should be 0 or negative. The year before 1 A.D. is 0, 10 B.C. is year -9.
month : int
Month as integer, Jan = 1, Feb. = 2, etc.
day : float
Day, may contain fractional part.

Returns

jd : float
Julian Day

Examples

Convert 6 a.m., February 17, 1985 to Julian Day

>>> date_to_jd(1985, 2, 17.25)
2446113.75
def datetime_to_jd(date)

Convert a datetime.datetime object to Julian Day.

Parameters

date : datetime.datetime instance
 

Returns

jd : float
Julian day.

Examples

>>> d = datetime.datetime(1985, 2, 17, 6)
>>> d
datetime.datetime(1985, 2, 17, 6, 0)
>>> jdutil.datetime_to_jd(d)
2446113.75
def days_to_hmsm(days)

Convert fractional days to hours, minutes, seconds, and microseconds. Precision beyond microseconds is rounded to the nearest microsecond.

Parameters

days : float
A fractional number of days. Must be less than 1.

Returns

hour : int
Hour number.
min : int
Minute number.
sec : int
Second number.
micro : int
Microsecond number.

Raises

ValueError
If days is >= 1.

Examples

>>> days_to_hmsm(0.1)
(2, 24, 0, 0)
def find_prior_limits(prior)

Find lower and upper limits of a prior from the kima_model_setup.txt file.

def find_prior_parameters(prior)

Find parameters of a prior from the kima_model_setup.txt file.

def get_instrument_name(data_file)

Find instrument name

def get_planet_mass(P, K, e, star_mass=1.0, full_output=False, verbose=False)

Calculate the planet (minimum) mass Msini given orbital period P, semi-amplitude K, eccentricity e, and stellar mass. If star_mass is a tuple with (estimate, uncertainty), this (Gaussian) uncertainty will be taken into account in the calculation.

Units

P [days] K [m/s] e [] star_mass [Msun]

Returns

if P is float:
if star_mass is float:
Msini [Mjup], Msini [Mearth]
if star_mass is tuple:
(Msini, error_Msini) [Mjup], (Msini, error_Msini) [Mearth]
if P is array:
if full_output
mean Msini [Mjup], std Msini [Mjup], Msini [Mjup] (array) else: mean Msini [Mjup], std Msini [Mjup], mean Msini [Mearth], std Msini [Mearth]
def get_planet_mass_and_semimajor_axis(P, K, e, star_mass=1.0, full_output=False, verbose=False)

Calculate the planet (minimum) mass Msini and the semi-major axis given orbital period P, semi-amplitude K, eccentricity e, and stellar mass. If star_mass is a tuple with (estimate, uncertainty), this (Gaussian) uncertainty will be taken into account in the calculation.

Units

P [days] K [m/s] e [] star_mass [Msun]

Returns

(M, A) where M is the output of get_planet_mass A is the output of get_planet_semimajor_axis

def get_planet_mass_latex(P, K, e, star_mass=1.0, earth=False, **kargs)
def get_planet_semimajor_axis(P, K, star_mass=1.0, full_output=False, verbose=False)

Calculate the semi-major axis of the planet's orbit given orbital period P, semi-amplitude K, and stellar mass.

Units

P [days] K [m/s] star_mass [Msun]

Returns

if P is float
a [AU]
if P is array:
if full_output
mean a [AU], std a [AU], a [AU] (array) else: mean a [AU], std a [AU]
def get_planet_semimajor_axis_latex(P, K, star_mass=1.0, earth=False, **kargs)
def get_planet_teq(Tstar=5777, Rstar=1.0, a=1.0, A=0, f=1)

Tstar = stellar effective temperature in kelvins (K) Rstar = stellar radius in solar radii a = distance from the star in AU (change to solRad) A = bond albedo f = redistribution factor Returns equilibrium temperature in Kelvins

def get_prior(prior)

Return a scipt.stats-like prior from a kima_model_setup.txt prior

def get_star_name(data_file)

Find star name (usually works for approx. standard filenames)

def get_transit_probability(Rstar=1.0, a=1.0)

Transit probability, simple. Eq 6.1 in Exoplanet Handbook Rstar: stellar radius [Rsun] a: semi-major axis [au]

def hmsm_to_days(hour=0, min=0, sec=0, micro=0)

Convert hours, minutes, seconds, and microseconds to fractional days.

Parameters

hour : int, optional
Hour number. Defaults to 0.
min : int, optional
Minute number. Defaults to 0.
sec : int, optional
Second number. Defaults to 0.
micro : int, optional
Microsecond number. Defaults to 0.

Returns

days : float
Fractional days.

Examples

>>> hmsm_to_days(hour=6)
0.25
def hyperprior_samples(size)
def jd_to_date(jd)

Convert Julian Day to date.

Algorithm from 'Practical Astronomy with your Calculator or Spreadsheet', 4th ed., Duffet-Smith and Zwart, 2011.

Parameters

jd : float
Julian Day

Returns

year : int
Year as integer. Years preceding 1 A.D. should be 0 or negative. The year before 1 A.D. is 0, 10 B.C. is year -9.
month : int
Month as integer, Jan = 1, Feb. = 2, etc.
day : float
Day, may contain fractional part.

Examples

Convert Julian Day 2446113.75 to year, month, and day.

>>> jd_to_date(2446113.75)
(1985, 2, 17.25)
def jd_to_mjd(jd)

Convert Julian Day to Modified Julian Day

Parameters

jd : float
Julian Day

Returns

mjd : float
Modified Julian Day
def lighten_color(color, amount=0.5)

Lightens the given color by multiplying (1-luminosity) by the given amount. Input can be a matplotlib color string, hex string, or RGB tuple.

Examples:

>>> lighten_color('g', 0.3)
>>> lighten_color('#F034A3', 0.6)
>>> lighten_color((.3, .55, .1), 0.5)
def mjd_to_jd(mjd)

Convert Modified Julian Day to Julian Day.

Parameters

mjd : float
Modified Julian Day

Returns

jd : float
Julian Day
def percentile68_ranges(a, min=None, max=None)
def percentile68_ranges_latex(a, min=None, max=None)
def percentile_ranges(a, percentile=68, min=None, max=None)
def percentile_ranges_latex(a, percentile, min=None, max=None)
def priors_latex(results, title='Priors', label='tab:1')
def read_big_file(filename)
def read_datafile(datafile, skip)

Read data from datafile for multiple instruments. Can be str, in which case the 4th column is assumed to contain an integer identifier of the instrument. Or list, in which case each element will be one different filename containing three columns each.

def read_datafile_rvfwhm(datafile, skip)

Read data from datafile for multiple instruments and RV-FWHM data. Can be str, in which case the 4th column is assumed to contain an integer identifier of the instrument. Or list, in which case each element will be one different filename containing three columns each.

def read_model_setup(filename='kima_model_setup.txt')
def rms(array)

Root mean square of array

def show_tips()

Show a few tips on how to use kima (but only sometimes)

def wrms(array, weights)

Weighted root mean square of array, given weights

Classes

class Fixed (value: float)

Not a real distribution, just a parameter fixed at a given value

Methods

def logpdf(self, x: float) ‑> float

Logarithm of the probability density function

def pdf(self, x: float) ‑> float

Probability density function: 1.0 if x==value, otherwise 0.0

def rvs(self) ‑> float

Random sample (always returns value)

class ZeroDist

A dummy probability distribution which always returns 0.0

Methods

def logpdf(self, x: float) ‑> float
def pdf(self, x: float) ‑> float
def rvs(self) ‑> float