Representing Orbits with Classical (Keplerian) Elements

Introduction to Keplerian Orbital Elements

Perhaps the most well-known representation of the orbit is the Classical or Keplerian Elements. These elements are based on 6 parameters (conic geometry and its orientation in space) and are equivalent to the 6 element cartesian coordinates (position and velocity).

A long and detailed definition is well beyond the scope of this documentation. More information can be found in the TLE page in Wikipedia.

However, it must be noted that, these elements have singularities and ill-defined parameters for some special but common cases (such as true anomaly and argument of periapsis when the orbit is equatorial), and the implementation generally handles them gracefully. That said, parabolic orbits (eccentricity = 1.0) cannot be handled.

Mean vs Osculating Elements

A not very obvious fact to those learning about the Keplerian Orbital Elements is that it is just a format to represent the orbit; without further definition, it is not very useful. This is similar to cartesian coordinates needing a coordinate system definition to actually represent a position in space.

The actual implementation of the Keplerian Elements rests on the force model under which they are generated. For example, the TLEs are Mean Keplerian Orbital Elements generated with a specific theory that includes linearised J2, J3 and J4 zonal harmonics and the elements are defined in True Equator, Mean Equinox coordinate frame. They would be incompatible with other Mean Elements generated with other force models and theories. Osculating Keplerian Orbital Elements are generated with a two-body force model.

Therefore, over a trajectory generated with a two-body force model, the Osculating Keplerian Orbital Elements (apart from true anomaly) should stay constant, limited by the accuracy of the trajectory generation algorithm. On the other hand, osculating orbital elements should be used with care, particularly in contexts where instantaneous orbital elements are computed on a trajectory generated by a non-two-body model (e.g. including geopotentials). The orbital elements will not stay constant along the trajectory, simply because the force model over successive points are not strictly two-body.

Initialising and Using the Keplerian Orbital Elements

There is an Abstract Base Class AbstractKeplerianOrbitElements, from which the concrete implementation OsculatingKeplerianOrbElems derive. It represents the Osculating Orbital Elements for coordinates in the local inertial frame (e.g. GCRS for the Earth).

The usual way to initialise the class is via initial orbital elements (note the initialisation with units):


from astropy.time import Time
from astropy import units as u
from satmad.core.celestial_bodies_lib import EARTH
from satmad.propagation.classical_orb_elems import OsculatingKeplerianOrbElems

time = Time("2020-01-11T11:00:00.000", scale="utc")
central_body = EARTH

sm_axis = 7055.95320378 * u.km
ecc = 0.0020835 * u.dimensionless_unscaled
incl = 1.71234602 * u.rad
raan = 4.42402394 * u.rad
arg_perigee = 5.23982923 * u.rad
true_an = 1.5 * u.rad

orb_elems = OsculatingKeplerianOrbElems(
    time, sm_axis, ecc, incl, raan, arg_perigee, true_an, central_body
)

The second way to initialise the orbital elements is via cartesian coordinates, through the OsculatingKeplerianOrbElems.from_cartesian() method. The initial cartesian coordinates can be in any frame, they are automatically converted to the inertial frame of the central body. Once the Osculating Keplerian Elements are initialised, it is possible to query the parameters, and some derived parameters like period or radius of periapsis and apoapsis.


from astropy.time import Time
from astropy import units as u
from astropy.coordinates import GCRS, CartesianDifferential, CartesianRepresentation
from satmad.core.celestial_bodies_lib import EARTH
from satmad.coordinates.frames import init_pvt
from satmad.propagation.classical_orb_elems import OsculatingKeplerianOrbElems

time = Time("2020-01-11T11:00:00.000", scale="utc")
central_body = EARTH

r = CartesianRepresentation([7.213392947764267e+03, 8.523654531348812e+01, 2.783146976770290e-16], unit=u.km)
v = CartesianDifferential([5.902225938368851e-02, 7.421779936019859e+00, 1.595360086373873e-18], unit=u.km / u.s)
rv_init = init_pvt(GCRS, time, r.with_differentials(v))

orb_elems = OsculatingKeplerianOrbElems.from_cartesian(rv_init, central_body)

# Query parameters
print(f"Semimajor axis: {orb_elems.sm_axis}")
print(f"Period: {orb_elems.period}")
print(f"Radius of Periapsis: {orb_elems.periapsis}")
print(f"Period: {orb_elems.true_anomaly.to(u.deg)}")

Another useful output of the Osculating Keplerian Elements to convert them to the cartesian position and velocity in the inertial coordinate frame belonging to the orbital elements. Just appending the following line to the first example will yield the cartesian coordinates in a SkyCoord object.

>>> rv = orb_elems.to_cartesian()
>>> print(rv)

The conversions to and from the cartesian coordinates are based on GMAT [OM3]. Note that the classical orbital elements have a number of singularities for a lot of common orbits (e.g. circular and/or equatorial). GMAT Mathematical Specifications handles these cases gracefully, but care must be taken when interpreting the results. In some cases, return-trip testing may not be successful. However, the code is tested against GMAT and is working as expected.

Reference/API

Classical Orbital Elements definitions.

class satmad.propagation.classical_orb_elems.AbstractKeplerianOrbitElements(epoch, sm_axis, eccentricity, inclination, raan, arg_periapsis, true_anomaly, central_body=<satmad.core.celestial_body.CelestialBody object>)

Classical (Keplerian) Orbital Elements that should be overridden with the specific type (Mean, Osculating etc.).

Parameters
  • central_body (CelestialBody) – celestial body (e.g. Earth) around which the orbit is defined

  • epoch (Time) – Epoch Time corresponding to the orbital elements

  • sm_axis (Quantity) – semimajor axis of the orbit [km]

  • inclination (Quantity) – inclination of the orbit [rad]

  • raan (Quantity) – right ascension of ascending node (RAAN) of the orbit [rad]

  • eccentricity (Quantity or float) – eccentricity of the orbit [dimensionless]

  • arg_periapsis (Quantity) – argument of periapsis [rad]

  • true_anomaly (Quantity) – true anomaly of the orbit [rad]

Raises

ValueError – Parabolic orbits or singularity

property apoapsis

Computes the apoapsis distance [km].

property arg_periapsis

Gets and sets the argument of periapsis [rad].

Argument of periapsis should be in range 0 <= argp < 2*PI. Any input outside this range will be forced into this range.

property eccentricity

Gets and sets the eccentricity of the orbit.

Input can be a Quantity (dimensionless) or float. Getter value is a float.

Eccentricity should be in range 0 <= e < 1.0. Raises a ValueError otherwise.

property epoch

Returns the epoch time associated with the orbital parameters. Setter replicates the time with copy=False.

abstract classmethod from_cartesian(init_coords, central_body=<satmad.core.celestial_body.CelestialBody object>)

Generates Keplerian Elements from cartesian coordinates in the inertial frame of the celestial body.

If the initial coordinate is in a different frame than the inertial frame of the celestial body, it is automatically converted to the proper frame for conversion.

Parameters
  • init_coords (SkyCoord) – Initial coordinates (the first value is used)

  • central_body (CelestialBody) – celestial body (e.g. Earth) around which the orbit is defined

Returns

Osculating classical (or Keplerian) orbital elements

Return type

OsculatingKeplerianOrbElems

Raises

ValueError – Parabolic orbits or singularity

property inclination

Gets and sets the inclination of the orbit [rad].

Inclination should be in range 0 <= om < PI. Raises a ValueError otherwise.

property mean_motion

Computes the mean motion of the orbit [rad/sec].

Returns

mean motion in u.rad / u.s

Return type

Quantity

property periapsis

Computes the periapsis distance [km].

property period

Computes the satellite period [sec].

property raan

Gets and sets the right ascension of ascending node (RAAN) of the orbit [rad].

RAAN should be in range 0 <= om < 2*PI. Any input outside this range will be forced into this range.

property sm_axis

Gets the semimajor axis [km].

abstract to_cartesian()

Converts the orbital elements to the cartesian coordinates in the local inertial frame of the central body.

Returns

cartesian coordinates in the local inertial frame of the central body

Return type

SkyCoord

Raises

ValueError – Parabolic orbits or singularity

property true_anomaly

Gets and sets the true anomaly of the orbit [rad].

True Anomaly should be in range 0 <= true_anomaly < 2*PI. Any input outside this range will be forced into this range.

class satmad.propagation.classical_orb_elems.OsculatingKeplerianOrbElems(epoch, sm_axis, eccentricity, inclination, raan, arg_periapsis, true_anomaly, central_body=<satmad.core.celestial_body.CelestialBody object>)

Osculating Classical (Keplerian) Orbital Elements in the local inertial frame.

By definition, this uses a two-body potential for computations. Therefore, this is not a “mean elements” model such as a TLE. Over a trajectory generated with a two-body force model, the orbital elements (apart from true anomaly) should stay constant, limited by the accuracy of the trajectory generation algorithm.

Osculating orbital elements should be used with care, particularly in contexts where instantaneous orbital elements are computed on a trajectory generated by a non-two-body model (e.g. including geopotentials). The orbital elements will not stay constant along the trajectory, simply because the force model over successive points are not strictly two-body.

Parameters
  • central_body (CelestialBody) – celestial body (e.g. Earth) around which the orbit is defined

  • epoch (Time) – Epoch Time corresponding to the orbital elements

  • sm_axis (Quantity) – semimajor axis of the orbit [km]

  • inclination (Quantity) – inclination of the orbit [rad]

  • raan (Quantity) – right ascension of ascending node (RAAN) of the orbit [rad]

  • eccentricity (Quantity or float) – eccentricity of the orbit [dimensionless]

  • arg_periapsis (Quantity) – argument of periapsis [rad]

  • true_anomaly (Quantity) – true anomaly of the orbit [rad]

Raises

ValueError – Parabolic orbits or singularity

classmethod from_cartesian(init_coords, central_body=<satmad.core.celestial_body.CelestialBody object>)

Generates Osculating Keplerian Elements from cartesian coordinates in the inertial frame of the celestial body.

If the initial coordinate is in a different frame than the inertial frame of the celestial body, it is automatically converted to the proper frame for conversion.

Parameters
  • init_coords (SkyCoord) – Initial coordinates (the first value is used)

  • central_body (CelestialBody) – celestial body (e.g. Earth) around which the orbit is defined

Returns

Osculating classical (or Keplerian) orbital elements

Return type

OsculatingKeplerianOrbElems

Raises

ValueError – Parabolic orbits or singularity

to_cartesian()

Converts the orbital elements to the cartesian coordinates in the local inertial frame of the central body.

Returns

cartesian coordinates in the local inertial frame of the central body

Return type

SkyCoord

Raises

ValueError – Parabolic orbits or singularity