Basis
- class orbitize.basis.Basis(stellar_or_system_mass, mass_err, plx, plx_err, num_secondary_bodies, fit_secondary_mass, angle_upperlim=6.283185307179586, hipparcos_IAD=None, rv=False, rv_instruments=None)[source]
Abstract base class for different basis sets. All new basis objects should inherit from this class. This class is meant to control how priors are assigned to various basis sets and how conversions are made from the basis sets to the standard keplarian set.
Author: Tirth, 2021
- set_default_mass_priors(priors_arr, labels_arr)[source]
Adds the necessary priors for the stellar and/or companion masses.
- Parameters
priors_arr (list of orbitize.priors.Prior objects) – holds the prior objects for each parameter to be fitted (updated here)
labels_arr (list of strings) – holds the name of all the parameters to be fitted (updated here)
- set_hip_iad_priors(priors_arr, labels_arr)[source]
Adds the necessary priors relevant to the hipparcos data to ‘priors_arr’ and updates ‘labels_arr’ with the priors’ corresponding labels.
- Parameters
priors_arr (list of orbitize.priors.Prior objects) – holds the prior objects for each parameter to be fitted (updated here)
labels_arr (list of strings) – holds the name of all the parameters to be fitted (updated here)
- set_rv_priors(priors_arr, labels_arr)[source]
Adds the necessary priors if radial velocity data is supplied to ‘priors_arr’ and updates ‘labels_arr’ with the priors’ corresponding labels. This function assumes that ‘rv’ data has been supplied and a secondary mass is being fitted for.
- Parameters
priors_arr (list of orbitize.priors.Prior objects) – holds the prior objects for each parameter to be fitted (updated here)
labels_arr (list of strings) – holds the name of all the parameters to be fitted (updated here)
- class orbitize.basis.Period(stellar_or_system_mass, mass_err, plx, plx_err, num_secondary_bodies, fit_secondary_mass, angle_upperlim=6.283185307179586, hipparcos_IAD=None, rv=False, rv_instruments=None)[source]
Modification of the standard basis, swapping our sma for period: (per, ecc, inc, aop, pan, tau).
- Parameters
stellar_or_system_mass (float) – mass of the primary star (if fitting for dynamical masses of both components) or total system mass (if fitting using relative astrometry only) [M_sol]
mass_err (float) – uncertainty on ‘stellar_or_system_mass’, in M_sol
plx (float) – mean parallax of the system, in mas
plx_err (float) – uncertainty on ‘plx’, in mas
num_secondary_bodies (int) – number of secondary bodies in the system, should be at least 1
fit_secondary_mass (bool) – if True, include the dynamical mass of orbitting body as fitted parameter, if False, ‘stellar_or_system_mass’ is taken to be total mass
angle_upperlim (float) – either pi or 2pi, to restrict the prior range for ‘pan’ parameter (default: 2pi)
hipparcos_IAD (orbitize.HipparcosLogProb object) – if not ‘None’, then add relevant priors to this data (default: None)
rv (bool) – if True, then there is radial velocity data and assign radial velocity priors, if False, then there is no radial velocity data and radial velocity priors are not assigned (default: False)
rv_instruments (np.array) – array of unique rv instruments from the originally supplied data (default: None)
- construct_priors()[source]
Generates the parameter label array and initializes the corresponding priors for each parameter that’s to be sampled. For the standard basis, the parameters common to each companion are: per, ecc, inc, aop, pan, tau. Parallax, hipparcos (optional), rv (optional), and mass priors are added at the end.
- Returns
list: list of strings (labels) that indicate the names of each parameter to sample
list: list of orbitize.priors.Prior objects that indicate the prior distribution of each label
- Return type
tuple
- to_period_basis(param_arr)[source]
Convert parameter array from standard basis to period basis by swapping out the semi-major axis parameter to period for each companion. This function is used primarily for testing purposes.
- Parameters
param_arr (np.array of float) – RxM array of fitting parameters in the standard basis, where R is the number of parameters being fit, and M is the number of orbits. If M=1 (for MCMC), this can be a 1D array.
- Returns
- modifies ‘param_arr’ to contain the period for each companion
in each orbit rather than semi-major axis. Shape of ‘param_arr’ remains the same.
- Return type
np.array of float
- to_standard_basis(param_arr)[source]
Convert parameter array from period basis to standard basis by swapping out the period parameter to semi-major axis for each companion.
- Parameters
param_arr (np.array of float) – RxM array of fitting parameters in the period basis, where R is the number of parameters being fit, and M is the number of orbits. If M=1 (for MCMC), this can be a 1D array.
- Returns
- modifies ‘param_arr’ to contain the semi-major axis for each companion
in each orbit rather than period. Shape of ‘param_arr’ remains the same.
- Return type
np.array of float
- class orbitize.basis.SemiAmp(stellar_or_system_mass, mass_err, plx, plx_err, num_secondary_bodies, fit_secondary_mass, angle_upperlim=6.283185307179586, hipparcos_IAD=None, rv=False, rv_instruments=None)[source]
Modification of the standard basis, swapping our sma for period and additionally sampling in the stellar radial velocity semi-amplitude: (per, ecc, inc, aop, pan, tau, K).
Note
Ideally, ‘fit_secondary_mass’ is true and rv data is supplied.
- Parameters
stellar_or_system_mass (float) – mass of the primary star (if fitting for dynamical masses of both components) or total system mass (if fitting using relative astrometry only) [M_sol]
mass_err (float) – uncertainty on ‘stellar_or_system_mass’, in M_sol
plx (float) – mean parallax of the system, in mas
plx_err (float) – uncertainty on ‘plx’, in mas
num_secondary_bodies (int) – number of secondary bodies in the system, should be at least 1
fit_secondary_mass (bool) – if True, include the dynamical mass of orbitting body as fitted parameter, if False, ‘stellar_or_system_mass’ is taken to be total mass
angle_upperlim (float) – either pi or 2pi, to restrict the prior range for ‘pan’ parameter (default: 2*pi)
hipparcos_IAD (orbitize.HipparcosLogProb object) – if not ‘None’, then add relevant priors to this data (default: None)
rv (bool) – if True, then there is radial velocity data and assign radial velocity priors, if False, then there is no radial velocity data and radial velocity priors are not assigned (default: False)
rv_instruments (np.array) – array of unique rv instruments from the originally supplied data (default: None)
- compute_companion_mass(period, ecc, inc, semi_amp, m0)[source]
Computes a single companion’s mass given period, eccentricity, inclination, stellar rv semi-amplitude, and stellar mass. Uses scipy.fsolve to compute the masses numerically.
- Parameters
period (np.array of float) – the period values for each orbit for a single companion (can be float)
ecc (np.array of float) – the eccentricity values for each orbit for a single companion (can be float)
inc (np.array of float) – the inclination values for each orbit for a single companion (can be float)
semi_amp (np.array of float) – the stellar rv-semi amp values for each orbit (can be float)
m0 (np.array of float) – the stellar mass for each orbit (can be float)
- Returns
the companion mass values for each orbit (can also just be a single float)
- Return type
np.array of float
- compute_companion_sma(period, m0, m_n)[source]
Computes a single companion’s semi-major axis using Kepler’s Third Law for each orbit.
- Parameters
period (np.array of float) – the period values for each orbit for a single companion (can be float)
m0 (np.array of float) – the stellar mass for each orbit (can be float)
m_n (np.array of float) – the companion mass for each orbit (can be float)
- Returns
the semi-major axis values for each orbit
- Return type
np.array of float
- construct_priors()[source]
Generates the parameter label array and initializes the corresponding priors for each parameter that’s to be sampled. For the semi-amp basis, the parameters common to each companion are: per, ecc, inc, aop, pan, tau, K (stellar rv semi-amplitude). Parallax, hipparcos (optional), rv (optional), and mass priors are added at the end.
The mass parameter will always be m0.
- Returns
list: list of strings (labels) that indicate the names of each parameter to sample
list: list of orbitize.priors.Prior objects that indicate the prior distribution of each label
- Return type
tuple
- func(x, lhs, m0)[source]
Define function for scipy.fsolve to use when computing companion mass.
- Parameters
x (float) – the companion mass to be calculated (Msol)
lhs (float) – the left hand side of the rv semi-amplitude equation (Msol^(1/3))
m0 (float) – the stellar mass (Msol)
- Returns
- the difference between the rhs and lhs of the rv semi-amplitude equation, ‘x’ is a
good companion mass when this difference is very close to zero
- Return type
float
- to_semi_amp_basis(param_arr)[source]
Convert parameter array from standard basis to semi-amp basis by swapping out the semi-major axis parameter to period for each companion and computing the stellar rv semi-amplitudes for each companion.
- Parameters
param_arr (np.array of float) – RxM array of fitting parameters in the period basis, where R is the number of parameters being fit, and M is the number of orbits. If M=1 (for MCMC), this can be a 1D array.
- Returns
- modifies ‘param_arr’ to contain the semi-major axis for each companion
in each orbit rather than period, appends stellar rv semi-amplitude parameters, and removes companion masses
- Return type
np.array of float
- to_standard_basis(param_arr)[source]
Convert parameter array from semi-amp basis to standard basis by swapping out the period parameter to semi-major axis for each companion and computing the masses of each companion.
- Parameters
param_arr (np.array of float) – RxM array of fitting parameters in the period basis, where R is the number of parameters being fit, and M is the number of orbits. If M=1 (for MCMC), this can be a 1D array.
- Returns
- modifies ‘param_arr’ to contain the semi-major axis for each companion
in each orbit rather than period, removes stellar rv semi-amplitude parameters for each companion, and appends the companion masses to ‘param_arr’
- Return type
np.array of float
- class orbitize.basis.Standard(stellar_or_system_mass, mass_err, plx, plx_err, num_secondary_bodies, fit_secondary_mass, angle_upperlim=6.283185307179586, hipparcos_IAD=None, rv=False, rv_instruments=None)[source]
Standard basis set based upon the 6 standard Keplarian elements: (sma, ecc, inc, aop, pan, tau).
- Parameters
stellar_or_system_mass (float) – mass of the primary star (if fitting for dynamical masses of both components) or total system mass (if fitting using relative astrometry only) [M_sol]
mass_err (float) – uncertainty on ‘stellar_or_system_mass’, in M_sol
plx (float) – mean parallax of the system, in mas
plx_err (float) – uncertainty on ‘plx’, in mas
num_secondary_bodies (int) – number of secondary bodies in the system, should be at least 1
fit_secondary_mass (bool) – if True, include the dynamical mass of orbitting body as fitted parameter, if False, ‘stellar_or_system_mass’ is taken to be total mass
angle_upperlim (float) – either pi or 2pi, to restrict the prior range for ‘pan’ parameter (default: 2pi)
hipparcos_IAD (orbitize.HipparcosLogProb object) – if not ‘None’, then add relevant priors to this data (default: None)
rv (bool) – if True, then there is radial velocity data and assign radial velocity priors, if False, then there is no radial velocity data and radial velocity priors are not assigned (default: False)
rv_instruments (np.array) – array of unique rv instruments from the originally supplied data (default: None)
- construct_priors()[source]
Generates the parameter label array and initializes the corresponding priors for each parameter that’s to be sampled. For the standard basis, the parameters common to each companion are: sma, ecc, inc, aop, pan, tau. Parallax, hipparcos (optional), rv (optional), and mass priors are added at the end.
- Returns
list: list of strings (labels) that indicate the names of each parameter to sample
list: list of orbitize.priors.Prior objects that indicate the prior distribution of each label
- Return type
tuple
- to_standard_basis(param_arr)[source]
For standard basis, no conversion needs to be made.
- Parameters
param_arr (np.array of float) – RxM array of fitting parameters in the standard basis, where R is the number of parameters being fit, and M is the number of orbits. If M=1 (for MCMC), this can be a 1d array.
- Returns
param_arr
without any modification- Return type
np.array of float
- class orbitize.basis.XYZ(stellar_or_system_mass, mass_err, plx, plx_err, num_secondary_bodies, fit_secondary_mass, data_table, best_epoch_idx, epochs, angle_upperlim=6.283185307179586, hipparcos_IAD=None, rv=False, rv_instruments=None)[source]
Defines an orbit using the companion’s position and velocity components in XYZ space (x, y, z, xdot, ydot, zdot). The conversion algorithms used for this basis are defined in the following paper: http://www.dept.aoe.vt.edu/~lutze/AOE4134/9OrbitInSpace.pdf
Note
Does not have support with sep,pa data yet.
Note
Does not work for all multi-body data.
- Parameters
stellar_or_system_mass (float) – mass of the primary star (if fitting for dynamical masses of both components) or total system mass (if fitting using relative astrometry only) [M_sol]
mass_err (float) – uncertainty on ‘stellar_or_system_mass’, in M_sol
plx (float) – mean parallax of the system, in mas
plx_err (float) – uncertainty on ‘plx’, in mas
num_secondary_bodies (int) – number of secondary bodies in the system, should be at least 1
fit_secondary_mass (bool) – if True, include the dynamical mass of orbitting body as fitted parameter, if False, ‘stellar_or_system_mass’ is taken to be total mass
input_table (astropy.table.Table) – output from ‘orbitize.read_input.read_file()’
best_epoch_idx (list) – indices of the epochs corresponding to the smallest uncertainties
epochs (list) – all of the astrometric epochs from ‘input_table’
angle_upperlim (float) – either pi or 2pi, to restrict the prior range for ‘pan’ parameter (default: 2*pi)
hipparcos_IAD (orbitize.HipparcosLogProb object) – if not ‘None’, then add relevant priors to this data (default: None)
rv (bool) – if True, then there is radial velocity data and assign radial velocity priors, if False, then there is no radial velocity data and radial velocity priors are not assigned (default: False)
rv_instruments (np.array) – array of unique rv instruments from the originally supplied data (default: None)
Author: Rodrigo
- construct_priors()[source]
Generates the parameter label array and initializes the corresponding priors for each parameter that’s to be sampled. For the xyz basis, the parameters common to each companion are: x, y, z, xdot, ydot, zdot. Parallax, hipparcos (optional), rv (optional), and mass priors are added at the end.
The xyz basis describes the position and velocity vectors with reference to the local coordinate system (the origin of the system is star).
- Returns
list: list of strings (labels) that indicate the names of each parameter to sample
list: list of orbitize.priors.Prior objects that indicate the prior distribution of each label
- Return type
tuple
- standard_to_xyz(epoch, elems, tau_ref_epoch=58849, tolerance=1e-09, max_iter=100)[source]
Converts array of orbital elements from the regular base of Keplerian orbits to positions and velocities in xyz Uses code from orbitize.kepler
- Parameters
epoch (float) – Date in MJD of observation to calculate time of periastron passage (tau).
elems (np.array of floats) – Orbital elements (sma, ecc, inc, aop, pan, tau, plx, mtot). If more than 1 set of parameters is passed, the first dimension must be the number of orbital parameter sets, and the second the orbital elements.
- Returns
Orbital elements in xyz (x-coordinate [au], y-coordinate [au], z-coordinate [au], velocity in x [km/s], velocity in y [km/s], velocity in z [km/s], parallax [mas], total mass of the two-body orbit
(M_* + M_planet) [Solar masses])
- Return type
np.array
- to_standard_basis(param_arr)[source]
Makes a call to ‘xyz_to_standard’ to convert each companion’s xyz parameters to the standard parameters an returns the updated array for conversion.
- Parameters
param_arr (np.array of float) – RxM array of fitting parameters in the period basis, where R is the number of parameters being fit, and M is the number of orbits. If M=1 (for MCMC), this can be a 1D array.
- Returns
Orbital elements in the standard basis for all companions.
- Return type
np.array
- to_xyz_basis(param_arr)[source]
Makes a call to ‘standard_to_xyz’ to convert each companion’s standard keplerian parameters to the xyz parameters an returns the updated array for conversion.
- Parameters
param_arr (np.array of float) – RxM array of fitting parameters in the period basis, where R is the number of parameters being fit, and M is the number of orbits. If M=1 (for MCMC), this can be a 1D array.
- Returns
Orbital elements in the xyz for all companions.
- Return type
np.array
- verify_params()[source]
For now, additionally throws exceptions if data is supplied in sep/pa or if the best epoch for each body is one of the last two (which would inevtably mess up how the priors are imposed).
- xyz_to_standard(epoch, elems, tau_ref_epoch=58849)[source]
Converts array of orbital elements in terms of position and velocity in xyz to the standard basis.
- Parameters
epoch (float) – Date in MJD of observation to calculate time of periastron passage (tau).
elems (np.array of floats) – Orbital elements in xyz basis (x-coordinate [au], y-coordinate [au], z-coordinate [au], velocity in x [km/s], velocity in y [km/s], velocity in z [km/s], parallax [mas], total mass of the two-body orbit (M_* + M_planet) [Solar masses]). If more than 1 set of parameters is passed, the first dimension must be the number of orbital parameter sets, and the second the orbital elements.
- Returns
- Orbital elements in the standard basis
(sma, ecc, inc, aop, pan, tau, plx, mtot)
- Return type
np.array
- orbitize.basis.switch_tau_epoch(old_tau, old_epoch, new_epoch, period)[source]
Convert tau to another tau that uses a different referench epoch
- Parameters
old_tau (float or np.array) – old tau to convert
old_epoch (float or np.array) – old reference epoch (days, typically MJD)
new_epoch (float or np.array) – new reference epoch (days, same system as old_epoch)
period (float or np.array) – orbital period (years)
- Returns
new taus
- Return type
float or np.array
- orbitize.basis.tau_to_manom(date, sma, mtot, tau, tau_ref_epoch)[source]
Gets the mean anomlay. Wrapper for kepler.tau_to_manom, kept here for backwards compatibility.
- Parameters
date (float or np.array) – MJD
sma (float) – semi major axis (AU)
mtot (float) – total mass (M_sun)
tau (float) – epoch of periastron, in units of the orbital period
tau_ref_epoch (float) – reference epoch for tau
- Returns
mean anomaly on that date [0, 2pi)
- Return type
float or np.array
- orbitize.basis.tau_to_tp(tau, ref_epoch, period, after_date=None)[source]
Convert tau (epoch of periastron in fractional orbital period after ref epoch) to t_p (date in days, usually MJD, but works with whatever system ref_epoch is given in)
- Parameters
tau (float or np.array) – value of tau to convert
ref_epoch (float or np.array) – date (in days, typically MJD) that tau is defined relative to
period (float or np.array) – period (in years) that tau is noralized with
after_date (float) – tp will be the first periastron after this date. If None, use ref_epoch.
- Returns
corresponding t_p of the taus
- Return type
float or np.array
- orbitize.basis.tp_to_tau(tp, ref_epoch, period)[source]
Convert t_p to tau
- Parameters
tp (float or np.array) – value to t_p to convert (days, typically MJD)
ref_epoch (float or np.array) – reference epoch (in days) that tau is defined from. Same system as tp (e.g., MJD)
period (float or np.array) – period (in years) that tau is defined by
- Returns
corresponding taus
- Return type
float or np.array