Orbit

This is probably one of the main important modules of poliastro since it enables the user to create the poliastro.twobody.orbit.Orbit objects.

User can create or initialize an orbit in a very different ways but the most common ones are the following:

digraph {
   "poliastro.twobody.orbit.Orbit" -> "from_vectors", "from_classical", "from_equinoctial", "from_body_ephem";
}

Howeverm it is also possible to create ‘special’ orbits such us poliastro.twobody.orbit.Orbit.circular or even poliastro.twobody.orbit.Orbit.parabolic.

Furthermore, orbits can be sampled by poliastro.twobody.orbit.Orbit.sample or propagated in time by poliastro.twobody.orbit.Orbit.propagate. But even it is possible to apply a maneuver to an orbit with the poliastro.twobody.orbit.Orbit.apply_maneuver.

exception poliastro.twobody.orbit.TimeScaleWarning
exception poliastro.twobody.orbit.OrbitSamplingWarning
class poliastro.twobody.orbit.Orbit(state, epoch, plane)

Position and velocity of a body with respect to an attractor at a given time (epoch).

Regardless of how the Orbit is created, the implicit reference system is an inertial one. For the specific case of the Solar System, this can be assumed to be the International Celestial Reference System or ICRS.

attractor

Main attractor.

epoch

Epoch of the orbit.

plane

Fundamental plane of the frame.

frame

Reference frame of the orbit.

New in version 0.11.0.

r

Position vector.

v

Velocity vector.

a

Semilatus rectum.

p

Semimajor axis.

r_p

Radius of pericenter.

r_a

Radius of apocenter.

ecc

Eccentricity.

inc

Inclination.

raan

Right ascension of the ascending node.

argp

Argument of the perigee.

nu

True anomaly.

f

Second modified equinoctial element.

g

Third modified equinoctial element.

h

Fourth modified equinoctial element.

k

Fifth modified equinoctial element.

L

True longitude.

period

Period of the orbit.

n

Mean motion.

energy

Specific energy.

e_vec

Eccentricity vector.

h_vec

Specific angular momentum vector.

arglat

Argument of latitude.

classmethod from_vectors(attractor, r, v, epoch=<Time object: scale='tt' format='jyear_str' value=J2000.000>, plane=<Planes.EARTH_EQUATOR: 'Earth mean Equator and Equinox of epoch (J2000.0)'>)

Return Orbit from position and velocity vectors.

Parameters:
  • attractor (Body) – Main attractor.
  • r (Quantity) – Position vector wrt attractor center.
  • v (Quantity) – Velocity vector.
  • epoch (Time, optional) – Epoch, default to J2000.
  • plane (Planes) – Fundamental plane of the frame.
classmethod from_coords(attractor, coord, plane=<Planes.EARTH_EQUATOR: 'Earth mean Equator and Equinox of epoch (J2000.0)'>)

Creates an Orbit from an attractor and astropy SkyCoord or BaseCoordinateFrame instance.

This method accepts position and velocity in any reference frame unlike Orbit.from_vector which can accept inputs in only inertial reference frame centred at attractor. Also note that the frame information is lost after creation of the orbit and only the inertial reference frame at body centre will be used for all purposes.

Parameters:
  • attractor (Body) – Main attractor
  • coord (astropy.coordinates.SkyCoord or BaseCoordinateFrame) – Position and velocity vectors in any reference frame. Note that coord must have a representation and its differential with respect to time.
  • plane (Planes, optional) – Final orbit plane, default to Earth Equator.
classmethod from_classical(attractor, a, ecc, inc, raan, argp, nu, epoch=<Time object: scale='tt' format='jyear_str' value=J2000.000>, plane=<Planes.EARTH_EQUATOR: 'Earth mean Equator and Equinox of epoch (J2000.0)'>)

Return Orbit from classical orbital elements.

Parameters:
  • attractor (Body) – Main attractor.
  • a (Quantity) – Semi-major axis.
  • ecc (Quantity) – Eccentricity.
  • inc (Quantity) – Inclination
  • raan (Quantity) – Right ascension of the ascending node.
  • argp (Quantity) – Argument of the pericenter.
  • nu (Quantity) – True anomaly.
  • epoch (Time, optional) – Epoch, default to J2000.
  • plane (Planes) – Fundamental plane of the frame.
classmethod from_equinoctial(attractor, p, f, g, h, k, L, epoch=<Time object: scale='tt' format='jyear_str' value=J2000.000>, plane=<Planes.EARTH_EQUATOR: 'Earth mean Equator and Equinox of epoch (J2000.0)'>)

Return Orbit from modified equinoctial elements.

Parameters:
  • attractor (Body) – Main attractor.
  • p (Quantity) – Semilatus rectum.
  • f (Quantity) – Second modified equinoctial element.
  • g (Quantity) – Third modified equinoctial element.
  • h (Quantity) – Fourth modified equinoctial element.
  • k (Quantity) – Fifth modified equinoctial element.
  • L (Quantity) – True longitude.
  • epoch (Time, optional) – Epoch, default to J2000.
  • plane (Planes) – Fundamental plane of the frame.
classmethod from_body_ephem(body, epoch=None)

Return osculating Orbit of a body at a given time.

classmethod from_horizons(name, epoch=None, plane=<Planes.EARTH_EQUATOR: 'Earth mean Equator and Equinox of epoch (J2000.0)'>, id_type='smallbody')

Return osculating Orbit of a body using JPLHorizons module of Astroquery.

Parameters:
  • name (string) – Name of the body to query for.
  • epoch (Time, optional) – Epoch, default to None.
  • plane (Planes) – Fundamental plane of the frame.
  • id_type (string, optional) – Use “smallbody” for Asteroids and Comets, and “majorbody” for Planets and Satellites.
classmethod circular(attractor, alt, inc=<Quantity 0. deg>, raan=<Quantity 0. deg>, arglat=<Quantity 0. deg>, epoch=<Time object: scale='tt' format='jyear_str' value=J2000.000>, plane=<Planes.EARTH_EQUATOR: 'Earth mean Equator and Equinox of epoch (J2000.0)'>)

Return circular Orbit.

Parameters:
  • attractor (Body) – Main attractor.
  • alt (Quantity) – Altitude over surface.
  • inc (Quantity, optional) – Inclination, default to 0 deg (equatorial orbit).
  • raan (Quantity, optional) – Right ascension of the ascending node, default to 0 deg.
  • arglat (Quantity, optional) – Argument of latitude, default to 0 deg.
  • epoch (Time, optional) – Epoch, default to J2000.
  • plane (Planes) – Fundamental plane of the frame.
classmethod geostationary(attractor, angular_velocity=None, period=None, hill_radius=None)

Return the geostationary orbit for the given attractor and its rotational speed.

Parameters:
  • attractor (Body) – Main attractor.
  • angular_velocity (Quantity) – Rotational angular velocity of the attractor.
  • period (Quantity) – Attractor’s rotational period, ignored if angular_velocity is passed.
  • hill_radius (Quantity) – Radius of Hill sphere of the attractor (optional). Hill sphere radius(in contrast with Laplace’s SOI) is used here to validate the stability of the geostationary orbit, that is to make sure that the orbital radius required for the geostationary orbit is not outside of the gravitational sphere of influence of the attractor. Hill SOI of parent(if exists) of the attractor is ignored if hill_radius is not provided.
classmethod parabolic(attractor, p, inc, raan, argp, nu, epoch=<Time object: scale='tt' format='jyear_str' value=J2000.000>, plane=<Planes.EARTH_EQUATOR: 'Earth mean Equator and Equinox of epoch (J2000.0)'>)

Return parabolic Orbit.

Parameters:
  • attractor (Body) – Main attractor.
  • p (Quantity) – Semilatus rectum or parameter.
  • inc (Quantity, optional) – Inclination.
  • raan (Quantity) – Right ascension of the ascending node.
  • argp (Quantity) – Argument of the pericenter.
  • nu (Quantity) – True anomaly.
  • epoch (Time, optional) – Epoch, default to J2000.
  • plane (Planes) – Fundamental plane of the frame.
represent_as(representation)

Converts the orbit to a specific representation.

New in version 0.11.0.

Parameters:representation (BaseRepresentation) – Representation object to use. It must be a class, not an instance.

Examples

>>> from poliastro.examples import iss
>>> from astropy.coordinates import CartesianRepresentation, SphericalRepresentation
>>> iss.represent_as(CartesianRepresentation)
<CartesianRepresentation (x, y, z) in km
    (859.07256, -4137.20368, 5295.56871)
 (has differentials w.r.t.: 's')>
>>> iss.represent_as(CartesianRepresentation).xyz
<Quantity [  859.07256, -4137.20368,  5295.56871] km>
>>> iss.represent_as(CartesianRepresentation).differentials['s']
<CartesianDifferential (d_x, d_y, d_z) in km / s
    (7.37289205, 2.08223573, 0.43999979)>
>>> iss.represent_as(CartesianRepresentation).differentials['s'].d_xyz
<Quantity [7.37289205, 2.08223573, 0.43999979] km / s>
>>> iss.represent_as(SphericalRepresentation)
<SphericalRepresentation (lon, lat, distance) in (rad, rad, km)
    (4.91712525, 0.89732339, 6774.76995296)
 (has differentials w.r.t.: 's')>
to_icrs()

Creates a new Orbit object with its coordinates transformed to ICRS.

Notice that, strictly speaking, the center of ICRS is the Solar System Barycenter and not the Sun, and therefore these orbits cannot be propagated in the context of the two body problem. Therefore, this function exists merely for practical purposes.

New in version 0.11.0.

rv()

Position and velocity vectors.

classical()

Classical orbital elements.

pqw()

Perifocal frame (PQW) vectors.

propagate(value, method=<function mean_motion>, rtol=1e-10, **kwargs)

Propagates an orbit a specified time.

If value is true anomaly, propagate orbit to this anomaly and return the result. Otherwise, if time is provided, propagate this Orbit some time and return the result.

Parameters:
  • value (Quantity, Time, TimeDelta) – Scalar time to propagate.
  • rtol (float, optional) – Relative tolerance for the propagation algorithm, default to 1e-10.
  • method (function, optional) – Method used for propagation
  • **kwargs – parameters used in perturbation models
Returns:

New orbit after propagation.

Return type:

Orbit

propagate_to_anomaly(value)

Propagates an orbit to a specific true anomaly.

Parameters:value (Quantity) –
Returns:Resulting orbit after propagation.
Return type:Orbit
sample(values=100, *, min_anomaly=None, max_anomaly=None, method=<function mean_motion>)

Samples an orbit to some specified time values.

New in version 0.8.0.

Parameters:
  • values (int) – Number of interval points (default to 100).
  • max_anomaly (min_anomaly,) – Anomaly limits to sample the orbit. For elliptic orbits the default will be \(E \in \left[0, 2 \pi \right]\), and for hyperbolic orbits it will be \(\nu \in \left[-\nu_c, \nu_c \right]\), where \(\nu_c\) is either the current true anomaly or a value that corresponds to \(r = 3p\).
  • method (function, optional) – Method used for propagation
Returns:

positions – Array of x, y, z positions, with proper times as the frame attributes if supported.

Return type:

BaseCoordinateFrame

Notes

When specifying a number of points, the initial and final position is present twice inside the result (first and last row). This is more useful for plotting.

Examples

>>> from astropy import units as u
>>> from poliastro.examples import iss
>>> iss.sample()  # doctest: +ELLIPSIS
<GCRS Coordinate ...>
>>> iss.sample(10)  # doctest: +ELLIPSIS
<GCRS Coordinate ...>
apply_maneuver(maneuver, intermediate=False)

Returns resulting Orbit after applying maneuver to self.

Optionally return intermediate states (default to False).

Parameters:
  • maneuver (Maneuver) – Maneuver to apply.
  • intermediate (bool, optional) – Return intermediate states, default to False.
plot(label=None, use_3d=False)

Plots the orbit as an interactive widget.

Parameters:
  • label (str, optional) – Label for the orbit, defaults to empty.
  • use_3d (bool, optional) – Produce a 3D plot, default to False.