keplertools package

Submodules

keplertools.CyKeplerSTM module

keplertools.CyKeplerSTM.CyKeplerSTM(ndarray x0, DTYPE_t dt, ndarray mus, DTYPE_t epsmult)

Kepler State Transition Matrix

Cythonized version of the Kepler State Transition Matrix Calculations. This provides a single method to iteratively call the backend KeplerSTM_C function on inputs.

Parameters:
  • x0 (ndarray) – 6n x 1 vector of stacked positions and velocities for n planets

  • dt (float) – Time step

  • mus (ndarray) – n x 1 vector of standard gravitational parameters mu = G(m+m_s) where m is the planet mass, m_s is the star mass and G is the gravitational constant

  • epsmult (float) – default multiplier on floating point precision, used as convergence metric. Higher values mean faster convergence, but sacrifice precision.

Returns:

Propagated orbital values (equivalent dimension to x0)

Return type:

x1 (ndarray)

Notes

All units must be complementary (i.e., if position is AU and velocity is AU/day, mu must be in AU^3/day^2).

keplertools.CyRV module

keplertools.CyRV.CyRV_from_E(ndarray rv, ndarray E, ndarray sinE, ndarray cosE, DTYPE_t e, DTYPE_t w, DTYPE_t K)

Finds radial velocity for a single object at the desired epochs, uses method from orvara that relies on eccentric anomaly and sine and cosine of the eccentric anomaly.

Parameters:
  • rv (ndarray) – Preexisting radial velocities, can also be zeros (rad)

  • E (ndarray) – Eccentric anomaly (rad)

  • sinE (ndarray) – Sine of eccentric anomaly (rad)

  • cosE (ndarray) – Cosine of eccentric anomaly (rad)

  • e (float) – Eccentricity

  • I (float) – Inclination (rad)

  • w (float) – Argument of periapsis (rad)

  • K (float) – RV semi-amplitude (m/s)

Returns:

Radial velocities

Return type:

rv (ndarray)

keplertools.CyRV.CyRV_from_time(ndarray rv, ndarray t, ndarray tp, ndarray per, ndarray e, ndarray w, ndarray K)

Finds radial velocity for a single object at the desired epochs

Parameters:
  • rv (ndarray) – Preexisting radial velocities, can also be zeros (rad)

  • t (ndarray) – Times of to calculate RV at (jd)

  • tp (float) – Time of periastron

  • per (float) – Period

  • e (float) – Eccentricity

  • w (float) – Argument of periapsis (rad)

  • K (float) – RV semi-amplitude (m/s)

Returns:

Radial velocities

Return type:

rv (ndarray)

keplertools.Cyeccanom module

keplertools.Cyeccanom.Cyeccanom(ndarray M, ndarray e, DTYPE_t epsmult, ITYPE_t maxIter)

Finds eccentric anomaly from mean anomaly and eccentricity

This method uses Newton-Raphson iteration to find the eccentric anomaly from mean anomaly and eccentricity, assuming a closed (0<e<1) orbit.

Parameters:
  • M (ndarray) – mean anomaly (rad)

  • e (ndarray) – eccentricity

  • epsmult (float) – Precision of convergence (multiplied by precision of floating data type).

  • maxiter (int) – Maximum numbr of iterations.

Returns:

eccentric anomaly (rad)

Return type:

E (float or ndarray)

Notes

M and e must be the same size.

keplertools.Cyeccanom.Cyeccanom_orvara(ndarray M, DTYPE_t e)

Finds eccentric anomaly E, sinE, and cosE from mean anomaly and eccentricity

This uses the method described the orvara paper which uses a 5th order polynomial to approximate E and does a single Newton-Raphson iteration to refine it.

Parameters:
Returns:

E (numpy.ndarray):

eccentric anomaly (rad)

sinE (numpy.ndarray):

Sine of eccentric anomaly

cosE (numpy.ndarray):

Cosine of eccentric anomaly

Return type:

tuple

Note

Currently only works for a single orbit since it relies on creating a lookup table for a single orbit. A loop can be added.

keplertools.angutils module

keplertools.fun module

keplertools.fun.RV_from_time_py(t: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes], tp: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes], per: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes], e: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes], w: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes], K: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes]) _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes][source]

Calculate radial velocities using the standard method.

Parameters:
Returns:

Array of radial velocities at each epoch (n,).

Return type:

rv (numpy.ndarray)

keplertools.fun.c2c3(psi: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes]) Tuple[ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]]][source]

Calculate the c2, c3 coefficients for the universal variable

Parameters:

psi (iterable or float) – psi = chi^2/a for universal variable chi and semi-major axis a

Returns:

c2 (numpy.ndarray):

c2 coefficients (same size as input)

c3 (numpy.ndarray):

c3 coefficients (same size as input)

Return type:

tuple

keplertools.fun.calcAB(a: ndarray[tuple[int, ...], dtype[float64]], e: ndarray[tuple[int, ...], dtype[float64]], O: ndarray[tuple[int, ...], dtype[float64]], I: ndarray[tuple[int, ...], dtype[float64]], w: ndarray[tuple[int, ...], dtype[float64]]) Tuple[ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]]][source]

Calculate inertial frame components of perifocal frame unit vectors scaled by orbit semi-major and semi-minor axes.

Note that these quantities are closely related to the Thiele-Innes constants

Parameters:
  • a (ndarray) – Semi-major axes

  • e (ndarray) – eccentricities

  • O (ndarray) – longitudes of ascending nodes (rad)

  • I (ndarray) – inclinations (rad)

  • w (ndarray) – arguments of pericenter (rad)

Returns:

A (ndarray):

Components of eccentricity vector scaled by a

B (ndarray):

Components of q vector (orthogonal to e and h) scaled by b (=asqrt{1-e^2})

Return type:

tuple

Notes

All inputs must be of same size. Outputs are 3xn for n input points. See Vinti (1998) for details on element/coord sys defintions.

keplertools.fun.calc_RV_from_M(M: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes], e: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes], w: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes], K: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes])[source]

Calculate the combined radial velocity of a system of n objects at m epochs.

Parameters:
  • M (numpy.ndarray) – Mean anomalies of the objects at desired epochs (n x m) (rad)

  • e (numpy.ndarray) – Eccentricities of the objects (n x 1)

  • w (numpy.ndarray) – Argument of periapsis of the objects (n x 1) (rad)

  • K (numpy.ndarray) – Semi-amplitudes of the objects (n x 1) (m/s)

Returns:

System radial velocities at desired epochs

Return type:

numpy.ndarray

keplertools.fun.calc_RV_from_time(t: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes], tp: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes], per: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes], e: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes], w: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes], K: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes], noc: bool = True) _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes][source]

Calculate the combined radial velocity of a system of n objects at m epochs.

Parameters:
  • t (numpy.ndarray) – Epochs in jd (m x 1)

  • tp (numpy.ndarray) – Time of periastrons of the objects (n x 1)

  • per (numpy.ndarray) – Period

  • e (numpy.ndarray) – Eccentricities of the objects (n x 1)

  • w (numpy.ndarray) – Argument of periapsis of the objects (n x 1) (rad)

  • K (numpy.ndarray) – Semi-amplitudes of the objects (n x 1) (m/s)

  • noc (bool) – Use the Cython implementation if True, otherwise use the pure Python implementation.

Returns:

System radial velocities at desired epochs

Return type:

numpy.ndarray

keplertools.fun.eccanom(M: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes], e: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes], epsmult: float = 4.01, maxIter: int = 100, returnIter: bool = False, noc: bool = False, verb: bool = False) Tuple[ndarray[tuple[int, ...], dtype[float64]], int] | ndarray[tuple[int, ...], dtype[float64]][source]

Finds eccentric anomaly from mean anomaly and eccentricity

This method uses Newton-Raphson iteration to find the eccentric anomaly from mean anomaly and eccentricity, assuming a closed (0<e<1) orbit.

Parameters:
  • M (float or ndarray) – mean anomaly (rad)

  • e (float or ndarray) – eccentricity (eccentricity may be a scalar if M is given as an array, but otherwise must match the size of M.)

  • epsmult (float) – Precision of convergence (multiplied by precision of floating data type). Optional, defaults to 4.01.

  • maxIter (int) – Maximum numbr of iterations. Optional, defaults to 100.

  • returnIter (bool) – Return number of iterations (defaults false, only available in python version, ignored if using C version)

  • noc (bool) – Don’t use C version even if it can be loaded.

  • verb (bool) – Print exactly which version (C or Python is being used)

Returns:

E (float or ndarray):

eccentric anomaly (rad)

numIter (int):

Number of iterations (returned only if returnIter=True)

Return type:

tuple

Notes

If either M or e are scalar, and the other input is an array, the scalar input will be expanded to the same size array as the other input. So, a scalar M and array e will result in the calculation of the eccentric anomaly for one mean anomaly at a variety of eccentricities, and a scalar e and array M input will result in the calculation of eccentric anomalies for one eccentricity at a variety of mean anomalies. If both inputs are arrays then they are matched element by element.

keplertools.fun.eccanom_orvara(M: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes], e: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes]) Tuple[ndarray[tuple[int, ...], dtype[float64]], int][source]

Finds eccentric anomaly E, sinE, cosE from mean anomaly and eccentricity

This uses the method described in the orvara paper which uses a 5th order polynomial to approximate E and does a single Newton-Raphson iteration to refine it.

Parameters:
  • M (float or ndarray) – mean anomaly (rad)

  • e (float or ndarray) – eccentricity (eccentricity may be a scalar if M is given as an array, but otherwise must match the size of M.)

Returns:

E (float or ndarray):

eccentric anomaly (rad)

sinE (float or ndarray):

Sine of eccentric anomaly (rad)

cosE (float or ndarray):

Cosine of eccentric anomaly (rad)

Return type:

tuple

Notes

If either M or e is scalar, and the other input is an array, the scalar input will be expanded to the same size array as the other input. If both inputs are arrays then they are matched element by element.

keplertools.fun.forcendarray(x: float | ndarray[tuple[int, ...], dtype[float64]]) ndarray[tuple[int, ...], dtype[float64]][source]

Convert any numerical value into 1-D ndarray

Parameters:

x (float or numpy.ndarray) – Input

Returns:

Same size as input but in ndarray form

Return type:

numpy.ndarray

keplertools.fun.invKepler(M: float | ndarray[tuple[int, ...], dtype[float64]], e: float | ndarray[tuple[int, ...], dtype[float64]], tol: float | None = None, E0: float | ndarray[tuple[int, ...], dtype[float64]] | None = None, maxIter: int = 100, return_nu: bool = False, convergence_error: bool = True) Tuple[ndarray[tuple[int, ...], dtype[float64]], ...][source]

Finds eccentric/hyperbolic/parabolic anomaly from mean anomaly and eccentricity

This method uses Newton-Raphson iteration to find the eccentric anomaly from mean anomaly and eccentricity, assuming a closed (0<e<1) orbit.

Parameters:
  • M (float or ndarray) – mean anomaly (rad)

  • e (float or ndarray) – eccentricity (eccentricity may be a scalar if M is given as an array, but otherwise must match the size of M.)

  • tol (float) – Convergence of tolerance. Defaults to eps(2*pi)

  • E0 (float or ndarray) – Initial guess for iteration. Defaults to Taylor-expansion based value for closed orbits and Vallado-derived heuristic for open orbits. If set, must match size of M.

  • maxIter (int) – Maximum numbr of iterations. Optional, defaults to 100.

  • return_nu (bool) – Return true anomaly (defaults false)

  • convergence_error (bool) – Raise error on convergence failure. Defaults True. If false, throws a warning.

Returns:

E (ndarray):

eccentric/parabolic/hyperbolic anomaly (rad)

numIter (ndarray):

Number of iterations

nu (ndarray):

True anomaly (returned only if return_nu=True)

Return type:

tuple

Notes

If either M or e are scalar, and the other input is an array, the scalar input will be expanded to the same size array as the other input. So, a scalar M and array e will result in the calculation of the eccentric anomaly for one mean anomaly at a variety of eccentricities, and a scalar e and array M input will result in the calculation of eccentric anomalies for one eccentricity at a variety of mean anomalies. If both inputs are arrays then they are matched element by element.

keplertools.fun.kepler2orbstate(a: float | ndarray[tuple[int, ...], dtype[float64]], e: float | ndarray[tuple[int, ...], dtype[float64]], O: float | ndarray[tuple[int, ...], dtype[float64]], I: float | ndarray[tuple[int, ...], dtype[float64]], w: float | ndarray[tuple[int, ...], dtype[float64]], mu: float | ndarray[tuple[int, ...], dtype[float64]], nu: float | ndarray[tuple[int, ...], dtype[float64]]) Tuple[ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]]][source]

Calculate orbital state vectors from Keplerian elements

Parameters:
Returns:

r (numpy.ndarray):

Components of orbital radius (n x 3)

v (numpy.ndarray):

Components of orbital velocity (n x 3)

Return type:

tuple

Notes

r.flatten() and v.flatten() will automatically stack elements in the proper order in a 1D array

keplertools.fun.orbElem2vec(E: ndarray[tuple[int, ...], dtype[float64]], mus: float | ndarray[tuple[int, ...], dtype[float64]], orbElem: Tuple[ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]]] | None = None, AB: Tuple[ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]]] | None = None, returnAB: bool = False) Tuple[ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]], Tuple[ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]]]] | Tuple[ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]]][source]

Convert Keplerian orbital elements to position and velocity vectors

Parameters:
  • E (ndarray) – nx1 array of eccentric anomalies (rad)

  • mus (ndarray or float) – nx1 array of gravitational parameters (G*m_i) where G is the gravitational constant and m_i is the mass of the ith body. if all vectors represent the same body, mus may be a scalar.

  • orbElem (tuple) – (a,e,O,I,w) Exact inputs to calcAB. Either this or AB input must be set

  • AB (tuple) – (A,B) Exact outpus from calcAB

  • returnAB (bool) – Default False. If True, returns (A,B) as thrid output.

Returns:

rs (ndarray):

3 x n stacked position vectors

vs (ndarray):

3 x n stacked velocity vectors

AB (tuple):

(A,B)

Return type:

tuple

Notes

All units are complementary, i.e., if mus are in AU^3/day^2 then positions will be in AU, and velocities will be AU/day.

Possible combinations or inputs are:

  1. E scalar, mu scalar - single body, single position. A, B should be 3x1 (or orbElem should be all scalars).

  2. E vector, mu scalar - single body, many orbital positions. A, B should be 3x1 (or orbElem should be all scalars).

  3. E vector, mu vector - multiple bodies at varying orbital positions. A, B should be 3xn where E.size==n (or all orbElem should be size n) and mus.size must equal E.size.

keplertools.fun.orbstate2kepler(r: ndarray[tuple[int, ...], dtype[float64]], v: ndarray[tuple[int, ...], dtype[float64]], mu: float | ndarray[tuple[int, ...], dtype[float64]]) Tuple[ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]]][source]

Calculate Keplerian elements given orbital state vectors

Parameters:
  • r (numpy.ndarray) – Components of orbital radius. 3n elements in 1D as [r1(1);r1(2);r1(3);r2(1);r2(2)r2(3);…;rn(1);rn(2);rn(3)] or in 2D as nx3 or 3xn

  • v (numpy.ndarray) – Components of orbital velocity. Same stacking as r

  • mu (float or numpy.ndarray) – Gravitational parameters. If float, assuming all state vectors belong to the same system.

Returns:

a (ndarray):

Semi-major axis (or semi-parameter where e = 1)

e (ndarray):

eccentricity

O (ndarray):

longitude of ascending node (rad)

I (ndarray):

inclination (rad)

w (ndarray):

arguments of periapsis (rad)

tp (ndarray):

time of periapsis passage

Return type:

tuple

keplertools.fun.orvara_vector_helper(M_val, e_val)[source]

Wraps the Cyeccanom_orvara function to handle single M value as array.

keplertools.fun.true2eccanom(nu: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes], e: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes]) ndarray[tuple[int, ...], dtype[float64]][source]

Finds eccentric anomaly from true anomaly and eccentricity

Parameters:
  • nu (float or ndarray) – true anomaly (rad)

  • e (float or ndarray) – eccentricity

Returns:

true anomaly (rad)

Return type:

ndarray

Notes

If either nu or e are scalar, and the other input is an array, the scalar input will be expanded to the same size array as the other input.

keplertools.fun.trueanom(E: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes], e: _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes]) ndarray[tuple[int, ...], dtype[float64]][source]

Finds true anomaly from eccentric anomaly and eccentricity

The implemented method corresponds to Eq. 6.28 in Green assuming a closed (0<e<1) orbit.

Parameters:
  • E (float or ndarray) – eccentric anomaly (rad)

  • e (float or ndarray) – eccentricity

Returns:

true anomaly (rad)

Return type:

ndarray

Notes

If either E or e are scalar, and the other input is an array, the scalar input will be expanded to the same size array as the other input.

keplertools.fun.unitvector(vec: ndarray[tuple[int, ...], dtype[float64]], mag: ndarray[tuple[int, ...], dtype[float64]]) ndarray[tuple[int, ...], dtype[float64]][source]

Return the unit vectors of an array of vectors

Parameters:
Returns:

Unit vectors in the same layout as input

Return type:

numpy.ndarray

keplertools.fun.universalfg(r0: ndarray[tuple[int, ...], dtype[float64]], v0: ndarray[tuple[int, ...], dtype[float64]], mu: float | ndarray[tuple[int, ...], dtype[float64]], dt: float | ndarray[tuple[int, ...], dtype[float64]], maxIter: int = 100, return_counter: bool = False, convergence_error: bool = True) Tuple[ndarray[tuple[int, ...], dtype[float64]], ...][source]

Propagate orbital state vectors by delta t via universal variable-based f and g

Parameters:
  • r0 (numpy.ndarray) – Components of orbital radius. 3n elements in 1D as [r1(1);r1(2);r1(3);r2(1);r2(2)r2(3);…;rn(1);rn(2);rn(3)] or in 2D as nx3 or 3xn

  • v0 (numpy.ndarray) – Components of orbital velocity. Same stacking as r

  • mu (float or numpy.ndarray) – Gravitational parameters. If float, assuming all state vectors belong to the same system.

  • dt (float or numpy.ndarray) – Propagation time. If float, assuming all states are propagated for the same time

  • maxIter (int) – Maximum numbr of iterations. Optional, defaults to 100.

  • return_counter (bool) – If True, returns the number of iterations for each input state. Defaults False.

  • convergence_error (bool) – Raise error on convergence failure if True. Defaults True.

Returns:

r (numpy.ndarray):

Components of orbital radius (n x 3)

v (numpy.ndarray):

Components of orbital velocity (n x 3)

counter(numpy.ndarray):

Number of required iterations (size n). Only returned if return_counter is True

Return type:

tuple

Note

r.flatten() and v.flatten() will automatically stack elements in the proper order in a 1D array

keplertools.fun.validateOrbitalStateInputs(r: ndarray[tuple[int, ...], dtype[float64]], v: ndarray[tuple[int, ...], dtype[float64]], mu: float | ndarray[tuple[int, ...], dtype[float64]]) Tuple[ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]]][source]

Validate and standardize dimensionality of orbital state vector inputs

Parameters:
  • r (numpy.ndarray) – Components of orbital radius. 3n elements in 1D as [r1(1);r1(2);r1(3);r2(1);r2(2)r2(3);…;rn(1);rn(2);rn(3)] or in 2D as nx3 or 3xn

  • v (numpy.ndarray) – Components of orbital velocity. Same stacking as r

  • mu (float or numpy.ndarray) – Gravitational parameters. If float, assuming all state vectors belong to the same system.

Returns:

r (numpy.ndarray):

Components of orbital radius. (n x 3)

v (numpy.ndarray):

Components of orbital velocity. (n x 3)

mu (numpy.ndarray):

Gravitational parameters. (size 1 or n)

Return type:

tuple

keplertools.fun.vec2orbElem(rs: ndarray[tuple[int, ...], dtype[float64]], vs: ndarray[tuple[int, ...], dtype[float64]], mus: float | ndarray[tuple[int, ...], dtype[float64]]) Tuple[ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]]][source]

Convert position and velocity vectors to Keplerian orbital elements

Implements the (corrected) algorithm from Vinti

Parameters:
  • rs (ndarray) – 3n x 1 stacked initial position vectors: [r1(1);r1(2);r1(3);r2(1);r2(2)r2(3);…;rn(1);rn(2);rn(3)] or 3 x n or n x 3 matrix of position vectors.

  • vs (ndarray) – 3n x 1 stacked initial velocity vectors or 3 x n or n x3 matrix

  • mus (ndarray or float) – nx1 array of gravitational parameters (G*m_i) where G is the gravitational constant and m_i is the mass of the ith body. if all vectors represent the same body, mus may be a scalar.

Returns:

a (ndarray):

Semi-major axes

e (ndarray):

eccentricities

E (ndarray):

eccentric anomalies

O (ndarray):

longitudes of ascending nodes (rad)

I (ndarray):

inclinations (rad)

w (ndarray):

arguments of pericenter (rad)

P (ndarray):

orbital periods

tau (ndarray):

time of periapsis crossing

Return type:

tuple

Notes

All units must be complementary, i.e., if positions are in AU, and time is in days, vs must be in AU/day, mus must be in AU^3/day^2

keplertools.fun.vec2orbElem2(rs: ndarray[tuple[int, ...], dtype[float64]], vs: ndarray[tuple[int, ...], dtype[float64]], mus: float | ndarray[tuple[int, ...], dtype[float64]]) Tuple[ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]]][source]

Convert position and velocity vectors to Keplerian orbital elements

Implements the algorithm from Vallado

Parameters:
  • rs (ndarray) – 3n x 1 stacked initial position vectors: [r1(1);r1(2);r1(3);r2(1);r2(2)r2(3);…;rn(1);rn(2);rn(3)] or 3 x n or n x 3 matrix of position vectprs.

  • vs (ndarray) – 3n x 1 stacked initial velocity vectors or 3 x n or n x3 matrix

  • mus (ndarray or float) – nx1 array of gravitational parameters (G*m_i) where G is the gravitational constant and m_i is the mass of the ith body. if all vectors represent the same body, mus may be a scalar.

Returns:

a (ndarray):

Semi-major axes

e (ndarray):

eccentricities

E (ndarray):

eccentric anomalies

O (ndarray):

longitudes of ascending nodes (rad)

I (ndarray):

inclinations (rad)

w (ndarray):

arguments of pericenter (rad)

P (ndarray):

orbital periods

tau (ndarray):

time of periapsis crossing

Return type:

tuple

Notes

All units must be complementary, i.e., if positions are in AU, and time is in days, vs must be in AU/day, mus must be in AU^3/day^2

keplertools.keplerSTM module

class keplertools.keplerSTM.planSys(x0, mu, epsmult=4.0, prefVallado=False, noc=False)[source]

Bases: object

Kepler State Transition Matrix

Class container for defining a planetary system (or group of planets in multiple systems) via their gravitational parameters and state vectors. Contains methods for propagating state vectors forward in time via the Kepler state transition matrix.

Parameters:
  • x0 (numpy.ndarray) – 6n vector of stacked positions and velocities for n planets

  • mu (numpy.ndarray) – n vector of standard gravitational parameters mu = G(m+m_s) where m is the planet mass, m_s is the star mass and G is the gravitational constant

  • epsmult (float) – default multiplier on floating point precision, used as convergence metric. Higher values mean faster convergence, but sacrifice precision.

  • prefVallado (bool) – If True, always try the Vallado algorithm first, otherwise try Shepherd first. Defaults False.

  • noc (bool) – Do not attempt to use cythonized code. Defaults False.

Note

All units must be complementary (i.e., if position is AU and velocity is AU/day, mu must be in AU^3/day^2.

Note

Two algorithms are implemented, both using Batting/Goodyear universal variables. The first is from Shepperd (1984), using continued fraction to solve the Kepler equation. The second is from Vallado (2004), using Newton iteration to solve the time equation. One algorithm is used preferentially, and the other is called only in the case of convergence failure on the first. All convergence is calculated to machine precision of the data type and variable size, scaled by a user-selected multiple.

calcSTM(dt)[source]

Compute STM for input time

Parameters:

dt (float) – Time step

Returns:

6x6 STM

Return type:

ndarray(float)

calcSTM_vallado(dt)[source]

Compute STM for input time

Parameters:

dt (float) – Time step

Returns:

6x6 STM

Return type:

ndarray(float)

contFrac(x, a=5.0, b=0.0, c=2.5)[source]

Compute continued fraction

Parameters:
Returns:

converged iterant

Return type:

ndarray(float)

psi2c2c3(psi0)[source]

Compute c_2 and c_3 values given psi

Parameters:

psi0 (float) – Input psi value

Returns:

c2 (float):

c2 value

c3 (float):

c3 value

Return type:

tuple

takeStep(dt)[source]

Propagate state by input time

Parameters:

dt (float) – Time step

updateState(x0)[source]

Update internal state variable and associated constants

Parameters:

x0 (ndarray(float)) – 6n vector of stacked positions and velocities for n planets

Module contents