Elements module

digraph {
   "poliastro.core.elements" -> "rv_pqw", "coe2rv", "coe2mee", "rv2coe";
}

The poliastro.core.elements module contains a set of functions related to orbital elements conversion.

This module contains a set of functions that can be used to convert between different elements that define the orbit of a body.

poliastro.core.elements.rv_pqw

Returns r and v vectors in perifocal frame.

\[ \begin{align}\begin{aligned}\begin{split}\vec{r} = \frac{h^2}{\mu}\frac{1}{1 + e\cos(\theta)}\begin{bmatrix} \cos(\theta)\\ \sin(\theta)\\ 0 \end{bmatrix} \\\\\\\end{split}\\\begin{split}\vec{v} = \frac{h^2}{\mu}\begin{bmatrix} -\sin(\theta)\\ e+\cos(\theta)\\ 0 \end{bmatrix}\end{split}\end{aligned}\end{align} \]
Parameters:
  • k (float) – Standard gravitational parameter (km^3 / s^2).
  • p (float) – Semi-latus rectum or parameter (km).
  • ecc (float) – Eccentricity.
  • nu (float) – True anomaly (rad).
Returns:

  • r (ndarray) – Position. Dimension 3 vector
  • v (ndarray) – Velocity. Dimension 3 vector

Examples

>>> from poliastro.core.elements import rv_pqw
>>> from poliastro.constants import GM_earth
>>> k = GM_earth #Earth gravitational parameter
>>> ecc = 0.3 #Eccentricity
>>> h = 60000e6 #Angular momentum of the orbit [m^2]/[s]
>>> nu = np.deg2rad(120) #True Anomaly [rad]
>>> p = h**2 / k #Parameter of the orbit
>>> r, v = rv_pqw(k, p, ecc, nu)
>>> #Printing the results
r = [-5312706.25105345  9201877.15251336    0] [m]
v = [-5753.30180931 -1328.66813933  0] [m]/[s]

Note

These formulas can be checked at Curtis 3rd. Edition, page 110. Also the example proposed is 2.11 of Curtis 3rd Edition book.

poliastro.core.elements.coe2rv

Converts from classical orbital elements to vectors.

Parameters:
  • k (float) – Standard gravitational parameter (km^3 / s^2).
  • p (float) – Semi-latus rectum or parameter (km).
  • ecc (float) – Eccentricity.
  • inc (float) – Inclination (rad).
  • omega (float) – Longitude of ascending node (rad).
  • argp (float) – Argument of perigee (rad).
  • nu (float) – True anomaly (rad).
poliastro.core.elements.coe2mee

Converts from classical orbital elements to modified equinoctial orbital elements.

The definition of the modified equinoctial orbital elements is taken from [Walker, 1985].

The modified equinoctial orbital elements are a set of orbital elements that are useful for trajectory analysis and optimization. They are valid for circular, elliptic, and hyperbolic orbits. These direct modified equinoctial equations exhibit no singularity for zero eccentricity and orbital inclinations equal to 0 and 90 degrees. However, two of the components are singular for an orbital inclination of 180 degrees.

\[\begin{split}\begin{align} p &= a(1-e^2) \\ f &= e\cos(\omega + \Omega) \\ g &= e\sin(\omega + \Omega) \\ h &= \tan(\frac{i}{2})\cos(\Omega) \\ k &= \tan(\frac{i}{2})\sin(\Omega) \\ L &= \Omega + \omega + \theta \\ \end{align}\end{split}\]
Parameters:
  • k (float) – Standard gravitational parameter (km^3 / s^2).
  • p (float) – Semi-latus rectum or parameter (km).
  • ecc (float) – Eccentricity.
  • inc (float) – Inclination (rad).
  • omega (float) – Longitude of ascending node (rad).
  • argp (float) – Argument of perigee (rad).
  • nu (float) – True anomaly (rad).

Note

The conversion equations are taken directly from the original paper.

poliastro.core.elements.rv2coe

Converts from vectors to classical orbital elements.

  1. First the angular momentum is computed:
    \[\vec{h} = \vec{r} \times \vec{v}\]
  2. With it the eccentricity can be solved:
    \[\begin{split}\begin{align} \vec{e} &= \frac{1}{\mu}\left [ \left ( v^{2} - \frac{\mu}{r}\right ) \vec{r} - (\vec{r} \cdot \vec{v})\vec{v} \right ] \\ e &= \sqrt{\vec{e}\cdot\vec{e}} \\ \end{align}\end{split}\]
  3. The node vector line is solved:
    \[\begin{split}\begin{align} \vec{N} &= \vec{k} \times \vec{h} \\ N &= \sqrt{\vec{N}\cdot\vec{N}} \end{align}\end{split}\]
  4. The rigth ascension node is computed:
    \[\begin{split}\Omega = \left\{ \begin{array}{lcc} cos^{-1}{\left ( \frac{N_{x}}{N} \right )} & if & N_{y} \geq 0 \\ \\ 360^{o} -cos^{-1}{\left ( \frac{N_{x}}{N} \right )} & if & N_{y} < 0 \\ \end{array} \right.\end{split}\]
  5. The argument of perigee:
    \[\begin{split}\omega = \left\{ \begin{array}{lcc} cos^{-1}{\left ( \frac{\vec{N}\vec{e}}{Ne} \right )} & if & e_{z} \geq 0 \\ \\ 360^{o} -cos^{-1}{\left ( \frac{\vec{N}\vec{e}}{Ne} \right )} & if & e_{z} < 0 \\ \end{array} \right.\end{split}\]
  6. And finally the true anomaly:
    \[\begin{split}\nu = \left\{ \begin{array}{lcc} cos^{-1}{\left ( \frac{\vec{e}\vec{r}}{er} \right )} & if & v_{r} \geq 0 \\ \\ 360^{o} -cos^{-1}{\left ( \frac{\vec{e}\vec{r}}{er} \right )} & if & v_{r} < 0 \\ \end{array} \right.\end{split}\]
Parameters:
  • k (float) – Standard gravitational parameter (km^3 / s^2)
  • r (array) – Position vector (km)
  • v (array) – Velocity vector (km / s)
  • tol (float, optional) – Tolerance for eccentricity and inclination checks, default to 1e-8
Returns:

  • p (float) – Semi-latus rectum of parameter (km)
  • ecc (float) – Eccentricity
  • inc (float) – Inclination (rad)
  • raan (float) – Right ascension of the ascending nod (rad)
  • argp (float) – Argument of Perigee (rad)
  • nu (float) – True Anomaly (rad)

Examples

>>> from poliastro.core.elements import rv2coe
>>> from poliastro.constants import GM_earth
>>> from astropy import units as u
>>> import numpy as np
>>> k = GM_earth.to(u.km**3 / u.s**2) #Earth gravitational parameter
>>> r = [-6045, -3490, 2500] * u.km
>>> v = [-3.457, 6.618, 2.533] * u.km / u.s
>>> p, ecc, inc, raan, argp, nu = rv2coe(k, r, v)
>>> print("p:", p, "[km]")
p: 8530.47436396927 [km]
>>> print("ecc:", ecc)
ecc: 0.17121118195416898
>>> print("inc:", np.rad2deg(inc), "[deg]")
inc: 153.2492285182475 [deg]
>>> print("raan:", np.rad2deg(raan), "[deg]")
raan: 255.27928533439618 [deg]
>>> print("argp:", np.rad2deg(argp), "[deg]")
argp: 20.068139973005362 [deg]
>>> print("nu:", np.rad2deg(nu), "[deg]")
nu: 28.445804984192108 [deg]

Note

This example is a real exercise from Orbital Mechanics for Engineering students by Howard D.Curtis. This exercise is 4.3 of 3rd. Edition, page 200.

poliastro.core.elements.mee2coe

Converts from modified equinoctial orbital elements to classical orbital elements.

The definition of the modified equinoctial orbital elements is taken from [Walker, 1985].

Note

The conversion is always safe because arctan2 works also for 0, 0 arguments.