System

class orbitize.system.System(num_secondary_bodies, data_table, stellar_mass, plx, mass_err=0, plx_err=0, restrict_angle_ranges=None, tau_ref_epoch=58849, fit_secondary_mass=False, results=None)[source]

A class to store information about a system (data & priors) and calculate model predictions given a set of orbital parameters.

Parameters
  • num_secondary_bodies (int) – number of secondary bodies in the system. Should be at least 1.

  • data_table (astropy.table.Table) – output from orbitize.read_input.read_file()

  • stellar_mass (float) – mean mass of the primary, in M_sol. See fit_secondary_mass docstring below.

  • plx (float) – mean parallax of the system, in mas

  • mass_err (float, optional) – uncertainty on stellar_mass, in M_sol

  • plx_err (float, optional) – uncertainty on plx, in mas

  • restrict_angle_ranges (bool, optional) – if True, restrict the ranges of the position angle of nodes to [0,180) to get rid of symmetric double-peaks for imaging-only datasets.

  • tau_ref_epoch (float, optional) – reference epoch for defining tau (MJD). Default is 58849 (Jan 1, 2020).

  • fit_secondary_mass (bool, optional) – if True, include the dynamical mass of the orbiting body as a fitted parameter. If this is set to False, stellar_mass is taken to be the total mass of the system. (default: False)

  • results (list of orbitize.results.Results) – results from an orbit-fit will be appended to this list as a Results class

Users should initialize an instance of this class, then overwrite priors they wish to customize.

Priors are initialized as a list of orbitize.priors.Prior objects, in the following order:

semimajor axis 1, eccentricity 1, inclination 1,
argument of periastron 1, position angle of nodes 1,
epoch of periastron passage 1,
[semimajor axis 2, eccentricity 2, etc.],
[parallax, [mass1, mass2, ..], total_mass/m0]

where 1 corresponds to the first orbiting object, 2 corresponds to the second, etc. Mass1, mass2, … correspond to masses of secondary bodies. If fit_secondary_mass is set to True, the last element of this list is initialized to the mass of the primary. If not, it is initialized to the total system mass.

Written: Sarah Blunt, Henry Ngo, Jason Wang, 2018

add_results(results)[source]

Adds an orbitize.results.Results object to the list in system.results

Parameters

results (orbitize.results.Results object) – add this object to list

clear_results()[source]

Removes all stored results

compute_model(params_arr)[source]

Compute model predictions for an array of fitting parameters.

Parameters

params_arr (np.array of float) – RxM array of fitting parameters, where R is the number of parameters being fit, and M is the number of orbits we need model predictions for. Must be in the same order documented in System() above. If M=1, this can be a 1d array.

Returns

Nobsx2xM array model predictions. If M=1, this is a 2d array, otherwise it is a 3d array.

Return type

np.array of float

convert_data_table_radec2seppa(body_num=1)[source]

Converts rows of self.data_table given in radec to seppa. Note that self.input_table remains unchanged.

Parameters

body_num (int) – which object to convert (1 = first planet)

orbitize.system.radec2seppa(ra, dec, mod180=False)[source]

Convenience function for converting from right ascension/declination to separation/ position angle.

Parameters
  • ra (np.array of float) – array of RA values, in mas

  • dec (np.array of float) – array of Dec values, in mas

  • mod180 (Bool) – if True, output PA values will be given in range [180, 540) (useful for plotting short arcs with PAs that cross 360 during observations) (default: False)

Returns

(separation [mas], position angle [deg])

Return type

tulple of float