Orbit Routines

The psoap.orbit module is the backbone for the main PSOAP tasks of computing radial velocities to shift the wavelength vectors. For a description of the various radial velocity orbital models, see Models. There is also an additional module included within the PSOAP packaged, psoap.orbit_astrometry, which is not used for the main Gaussian process routines, but includes functions for jointly modeling radial velocity measurements and relative astrometric measurements.

psoap.orbit

Orbital conventions derived according to

\[\begin{split}v_{r,1} = K_1 (\cos (\omega_1 + f) + e \cos \omega_1) \\ v_{r,2} = K_2 (\cos (\omega_2 + f) + e \cos \omega_2)\end{split}\]

where a positive \(v_r\) indicates a redshifted source and \(\omega_2 = \omega_1 + \pi\). In an RV-only paper, typically authors report \(\omega = \omega_1\).

In general, quantities like vA, vB, and vC refer to radial velocities relative to the heliocentric frame of the earth, i.e., quantities to directly compare to actual radial velocity measurements.

(Private) quantities like v1, v2, v1_in, and v1_out have more subtle meanings—these are generally partial velocities used in the process of constructing vA, etc...

class psoap.orbit.SB1(K, e, omega, P, T0, gamma, obs_dates=None, **kwargs)[source]

Orbital model for a single-line spectroscopic binary.

Parameters:
  • K (float) – velocity semi-amplitude
  • e (float) – eccentricity (must be between [0.0, 1.0))
  • omega (float) – argument of periastron for the primary star (degrees)
  • P (float) – period (days)
  • T0 (float) – epoch of periastron (JD)
  • gamma (float) – systemic velocity (km/s)
  • obs_dates (1D np.array) – dates of observation (JD)
get_velocities(dates=None)[source]

Calculate the velocities of the primary (vA) for all dates provided.

Parameters:dates (optional) – if provided, calculate velocities at this new vector of dates, rather than the one provided when the object was initialized.
Returns:A (1, nvels) shape array of the primary velocities
Return type:np.array
class psoap.orbit.SB2(q, K, e, omega, P, T0, gamma, obs_dates=None, **kwargs)[source]

Orbital model for a double-line spectroscopic binary.

Parameters:
  • q (float) – mass ratio M_B / M_A
  • K (float) – velocity semi-amplitude
  • e (float) – eccentricity (must be between [0.0, 1.0))
  • omega (float) – argument of periastron for primary star (degrees)
  • P (float) – period (days)
  • T0 (float) – epoch of periastron (JD)
  • gamma (float) – systemic velocity (km/s)
  • obs_dates (1D np.array) – dates of observation (JD)
get_velocities(dates=None)[source]

Calculate the velocities of the primary (vA) and secondary (vB) for all dates provided.

Parameters:dates (optional) – if provided, calculate velocities at this new vector of dates, rather than the one provided when the object was initialized.
Returns:A (2, nvels) shape array of the primary and secondary velocities
Return type:np.array
class psoap.orbit.ST1(K_in, e_in, omega_in, P_in, T0_in, K_out, e_out, omega_out, P_out, T0_out, gamma, obs_dates=None, **kwargs)[source]

Orbital model for a single-line hierarchical spectroscopic triple.

Parameters:
  • K_in (float) – velocity semi-amplitude for inner orbit
  • e_in (float) – eccentricity for inner orbit (must be between [0.0, 1.0))
  • omega_in (float) – argument of periastron for inner orbit (primary star) (degrees)
  • P_in (float) – inner period (days)
  • T0_in (float) – epoch of periastron for inner orbit (JD)
  • K_out (float) – velocity semi-amplitude for outer orbit
  • e_out (float) – eccentricity for outer orbit (must be between [0.0, 1.0))
  • omega_out (float) – argument of periastron for outer orbit (primary stars) (degrees)
  • P_out (float) – outer period (days)
  • T0_out (float) – epoch of periastron for outer orbit (JD)
  • gamma (float) – systemic velocity (km/s)
  • obs_dates (1D np.array) – dates of observation (JD)
get_component_velocities(dates=None)[source]

Routine to grab v1_inner and v1_outer.

Parameters:dates (optional) – if provided, calculate velocities at this new vector of dates, rather than the one provided when the object was initialized.
Returns:A (1, nvels) shape array of the primary velocities for the outer motion
Return type:np.array
get_velocities(dates=None)[source]

Calculate the velocities of the primary (vA) for all dates provided.

Parameters:dates (optional) – if provided, calculate velocities at this new vector of dates, rather than the one provided when the object was initialized.
Returns:A (1, nvels) shape array of the primary velocities
Return type:np.array
class psoap.orbit.ST2(q_in, K_in, e_in, omega_in, P_in, T0_in, K_out, e_out, omega_out, P_out, T0_out, gamma, obs_dates=None, **kwargs)[source]

Orbital model for a double-line (inner orbit) hierarchical spectroscopic triple.

Parameters:
  • q_in (float) – mass ratio M_B / M_A
  • K_in (float) – velocity semi-amplitude for inner orbit
  • e_in (float) – eccentricity for inner orbit (must be between [0.0, 1.0))
  • omega_in (float) – argument of periastron for inner orbit (degrees)
  • P_in (float) – inner period (days)
  • T0_in (float) – epoch of periastron for inner orbit (JD)
  • K_out (float) – velocity semi-amplitude for outer orbit
  • e_out (float) – eccentricity for outer orbit (must be between [0.0, 1.0))
  • omega_out (float) – argument of periastron for outer orbit (degrees)
  • P_out (float) – outer period (days)
  • T0_out (float) – epoch of periastron for outer orbit (JD)
  • gamma (float) – systemic velocity (km/s)
  • obs_dates (1D np.array) – dates of observation (JD)
get_component_velocities(dates=None)[source]

Routine to grab v1_inner and v1_outer.

Parameters:dates (optional) – if provided, calculate velocities at this new vector of dates, rather than the one provided when the object was initialized.
Returns:A (1, nvels) shape array of the primary velocities for the outer motion
Return type:np.array
get_velocities(dates=None)[source]

Calculate the velocities of the primary (vA) and secondary (vB) for all dates provided.

Parameters:dates (optional) – if provided, calculate velocities at this new vector of dates, rather than the one provided when the object was initialized.
Returns:A (2, nvels) shape array of the primary and secondary velocities
Return type:np.array
class psoap.orbit.ST3(q_in, K_in, e_in, omega_in, P_in, T0_in, q_out, K_out, e_out, omega_out, P_out, T0_out, gamma, obs_dates=None, **kwargs)[source]

Orbital model for a double-line (inner orbit) hierarchical spectroscopic triple.

Parameters:
  • q_in (float) – mass ratio M_B / M_A
  • K_in (float) – velocity semi-amplitude for inner orbit
  • e_in (float) – eccentricity for inner orbit (must be between [0.0, 1.0))
  • omega_in (float) – argument of periastron for inner orbit (degrees)
  • P_in (float) – inner period (days)
  • T0_in (float) – epoch of periastron for inner orbit (JD)
  • q_out (float) – mass ratio M_C / (M_A + M_B)
  • K_out (float) – velocity semi-amplitude for outer orbit
  • e_out (float) – eccentricity for outer orbit (must be between [0.0, 1.0))
  • omega_out (float) – argument of periastron for outer orbit (degrees)
  • P_out (float) – outer period (days)
  • T0_out (float) – epoch of periastron for outer orbit (JD)
  • gamma (float) – systemic velocity (km/s)
  • obs_dates (1D np.array) – dates of observation (JD)
get_component_velocities(dates=None)

Routine to grab v1_inner and v1_outer.

Parameters:dates (optional) – if provided, calculate velocities at this new vector of dates, rather than the one provided when the object was initialized.
Returns:A (1, nvels) shape array of the primary velocities for the outer motion
Return type:np.array
get_velocities(dates=None)[source]

Calculate the velocities of the primary (vA), secondary (vB), and tertiary (vC) for all dates provided.

Parameters:dates (optional) – if provided, calculate velocities at this new vector of dates, rather than the one provided when the object was initialized.
Returns:A (3, nvels) shape array of the primary, secondary, and tertiary velocities
Return type:np.array

psoap.orbit_astrometry

This is an alternate parameterization suited for joint astrometric and radial velocity orbital modeling.

Orbital conventions derived according to

\[\begin{split}\Delta \delta = X = \rho \cos \theta \\ \Delta \alpha \cos \delta = Y = \rho \sin \theta\end{split}\]

The positive X-axis points North, and the positive Y-axis points East. Rotation angle \(\theta\) is measured in degrees East of North, i.e.

\[\theta = \arctan(Y / X).\]

In our implementation, we use np.atan2(Y/X) to resolve the quadrant ambiguity.

The ascending node is defined as the point where the secondary crosses the plane of the sky receeding from the observer. It is measured in degrees East of North.

Inclination is defined such that \(0 < i < \pi/2\) yield prograde orbits (counter-clockwise, such that \(\theta\) increases), and \(\pi/2 < i < \pi\) yield retrograde orbits (clockwise, such that \(\theta\) decreases).

These lead to the following equations of motion

\[\begin{split}X = r (\cos \Omega \cos(\omega + f) - \sin(\Omega) \sin(\omega + f) \cos(i)) \\ Y = r (\sin \Omega \cos(\omega + f) + \cos(\Omega) \sin(\omega + f) \cos(i)) \\ Z = - r \sin(\omega + f) \sin(i)\end{split}\]

and

\[\begin{split}v_{r,1} = K_1 (\cos (\omega_1 + f) + e \cos \omega_1) \\ v_{r,2} = K_2 (\cos (\omega_2 + f) + e \cos \omega_2).\end{split}\]
class psoap.orbit_astrometry.Binary(a, e, i, omega, Omega, T0, M_tot, M_2, gamma, obs_dates=None, **kwargs)[source]

Binary orbital model that can deliver absolute astrometric position, relative astrometric position (B relative to A), and radial velocities for A and B.

Parameters:
  • a (float) – semi-major axis [AU]
  • e (float) – eccentricity (must be between [0.0, 1.0))
  • i (float) – inclination [deg]
  • omega (float) – argument of periastron of the primary, i.e. \(\omega_1\) [degrees]
  • Omega (float) – position angle of the ascending node (going into sky) [deg] east of north
  • T0 (float) – epoch of periastron passage [JD]
  • M_tot (float) – sum of the masses \(M_\mathrm{A} + M_\mathrm{B}\) [\(M_\odot\)]
  • M_2 (float) – mass of B [\(M_\odot\)]
  • gamma (float) – systemic velocity (km/s)
  • obs_dates (1D np.array) – dates of observation (JD)
get_full_orbit(dates=None)[source]

Deliver the full set of astrometric and radial velocity quantities, namely the radial velocities vA, vB, the position of A and B relative to the center of mass in the plane of the sky (XY_A and XY_B, respectively), the position of B relative to the position of A in the plane of the sky (XY_AB), the position of A and B in the plane of the orbit (xy_A and xy_B, respectively), and the position of B relative to the position of A in the plane of the orbit (xy_AB), for all dates provided. All positions are given in units of AU, and so must be converted to arcseconds after assuming a distance to the system.

Parameters:dates (optional) – if provided, calculate quantities at this new vector of dates, rather than the one provided when the object was initialized.
Returns:A dictionary with items of "vAs", "vBs", "XYZ_As", "XYZ_Bs", "XYZ_ABs", "xy_As", "xy_Bs", "xy_ABs"
Return type:dict
get_orbit(dates=None)[source]

Deliver only the main quantities useful for performing a joint astrometric + RV fit to real data, namely the radial velocities vA, vB, the relative offsets \(\rho\), and relative position angles \(\theta\), for all dates provided. Relative offsets are provided in AU, and so must be converted to arcseconds after assuming a distance to the system. Relative position angles are given in degrees east of north.

Parameters:dates (optional) – if provided, calculate quantities at this new vector of dates, rather than the one provided when the object was initialized.
Returns:A dictionary with items "vAs", "vBs", "rhos", "thetas".
Return type:dict
class psoap.orbit_astrometry.Triple(a_in, e_in, i_in, omega_in, Omega_in, T0_in, a_out, e_out, i_out, omega_out, Omega_out, T0_out, M_1, M_2, M_3, gamma, obs_dates=None, **kwargs)[source]

Triple orbital model that can deliver absolute astrometric position, relative astrometric position (B relative to A, and C relative to A), and radial velocities for A, B, and C.

Parameters:
  • a_in (float) – semi-major axis for inner orbit [AU]
  • e_in (float) – eccentricity for inner orbit (must be between [0.0, 1.0))
  • i_in (float) – inclination for inner orbit [deg]
  • omega_in (float) – argument of periastron for inner orbit, for the primary (\(\omega_1\)) [degrees]
  • Omega_in (float) – position angle of the ascending node [deg] east of north for inner orbit
  • T0_in (float) – epoch of periastron passage for inner orbit [JD]
  • a_out (float) – semi-major axis for outer orbit [AU]
  • e_out (float) – eccentricity for outer orbit (must be between [0.0, 1.0))
  • i_out (float) – inclination for outer orbit [deg]
  • omega_out (float) – argument of periastron for outer orbit, for the “primary” (\(\omega_{12}\)) [degrees]
  • Omega_out (float) – position angle of the ascending node [deg] east of north for outer orbit
  • T0_out (float) – epoch of periastron passage for outer orbit [JD]
  • M_1 (float) – mass of A [\(M_\odot\)]
  • M_2 (float) – mass of B [\(M_\odot\)]
  • M_3 (float) – mass of C [\(M_\odot\)]
  • gamma (float) – systemic velocity (km/s)
  • obs_dates (1D np.array) – dates of observation (JD)
get_full_orbit(dates=None)[source]

Deliver the full set of astrometric and radial velocity quantities, namely the radial velocities \(v_{r,A}\), \(v_{r,B}\), \(v_{r,C}\), the position of A, B, and C relative to the center of mass in the plane of the sky (XYZ_A, XYZ_B, and XYZ_C, respectively), the absolute position of the center of mass of (AB), (XYZ_AB), the position of A relative to the center of mass of AB (XYZ_A_loc), the position of B relative to the center of mass of (AB) (XYZ_B_loc), the absolute position of C in the plane of the orbit (xy_C), the absolute positon of the center of mass of AB in the plane of the orbit (xy_AB), the position of A in the plane of the orbit, relative to the center of mass of AB (xy_A_locs), and the position of B in the plane of the orbit, relative to the center of mass of AB (xy_B_locs), for all dates provided. All positions are given in units of AU, and so must be converted to arcseconds after assuming a distance to the system.

Parameters:dates (optional) – if provided, calculate quantities at this new vector of dates, rather than the one provided when the object was initialized.
Returns:A dictionary with keys vAs, vBs, vCs, XYZ_As, XYZ_Bs, XYZ_Cs, XYZ_ABs, XYZ_A_locs, XYZ_B_locs, xy_Cs, xy_ABs, xy_A_locs, xy_B_locs
Return type:dict
get_orbit(dates=None)[source]

Deliver only the main quantities useful for performing a joint astrometric + RV fit to real data, namely the radial velocities \(v_{r,A}\), \(v_{r,B}\), \(v_{r,C}\), the relative offsets of B to A \(\rho_\mathrm{AB}\), and relative position angles \(\theta_\mathrm{AB}\), and the relative offsets of C to A \(\rho_\mathrm{AC}\) and \(\theta_\mathrm{AC}\) for all dates provided. Relative offsets are provided in AU, and so must be converted to arcseconds after assuming a distance to the system. Relative position angles are given in degrees east of north.

Parameters:dates (optional) – if provided, calculate quantities at this new vector of dates, rather than the one provided when the object was initialized.
Returns:A dictionary with keys "vAs", "vBs", "vCs", "rho_ABs", "theta_ABs", "rho_ACs", "theta_ACs"
Return type:dict