API reference

This page holds poliastro’s API documentation, which might be helpful for final users or developers to create their own poliastro-based utilities. Among the different sub-packages and modules, we might differentiate two big categories: core utilities and high-level ones.

High level API

The high level API of poliastro allows you to do most common tasks (propagate an osculating orbit, sampling an ephemerides, compute maneuvers) in a straightforward way. All the methods expect Astropy units.

The most important high level objects and methods are poliastro.twobody.Orbit, poliastro.ephem.Ephem, and poliastro.maneuver.Maneuver. Here is a summarized reference of commonly used methods:

class poliastro.twobody.Orbit(state, epoch)

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.

propagate(value, method=FarnocchiaPropagator())

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.

  • method (function, optional) – Method used for propagation, default to farnocchia.

Returns

New orbit after propagation.

Return type

Orbit

to_ephem(strategy=TrueAnomalyBounds())

Samples Orbit to return an ephemerides.

New in version 0.17.0.

classmethod from_vectors(attractor, r, v, epoch=J2000, plane=Planes.EARTH_EQUATOR)

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_classical(attractor, a, ecc, inc, raan, argp, nu, epoch=J2000, plane=Planes.EARTH_EQUATOR)

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_sbdb(name, **kwargs)

Return osculating Orbit by using SBDB from Astroquery.

Parameters
  • name (str) – Name of the body to make the request.

  • **kwargs – Extra kwargs for astroquery.

Returns

ss – Orbit corresponding to body_name

Return type

poliastro.twobody.orbit.Orbit

Examples

>>> from poliastro.twobody.orbit import Orbit
>>> apophis_orbit = Orbit.from_sbdb('apophis')  
class poliastro.ephem.Ephem(coordinates, epochs, plane)

Time history of position and velocity of some object at particular epochs.

Instead of creating Ephem objects directly, use the available classmethods.

Parameters
classmethod from_body(body, epochs, *, attractor=None, plane=Planes.EARTH_EQUATOR)

Return Ephem for a SolarSystemPlanet at certain epochs.

Parameters
  • body (SolarSystemPlanet) – Body.

  • epochs (Time) – Epochs to sample the body positions.

  • attractor (SolarSystemPlanet, optional) – Body to use as central location, if not given the Solar System Barycenter will be used.

  • plane (Planes, optional) – Fundamental plane of the frame, default to Earth Equator.

classmethod from_horizons(name, epochs, *, attractor=None, plane=Planes.EARTH_EQUATOR, id_type=None)

Return Ephem for an object using JPLHorizons module of Astroquery.

Parameters
  • name (str) – Name of the body to query for.

  • epochs (Time) – Epochs to sample the body positions.

  • attractor (SolarSystemPlanet, optional) – Body to use as central location, if not given the Solar System Barycenter will be used.

  • plane (Planes, optional) – Fundamental plane of the frame, default to Earth Equator.

  • id_type (NoneType or str, optional) – Use “smallbody” for Asteroids and Comets and None (default) to first search for Planets and Satellites.

classmethod from_orbit(orbit, epochs, plane=Planes.EARTH_EQUATOR)

Return Ephem from an Orbit at certain epochs.

Parameters
  • orbit (Orbit) – Orbit.

  • epochs (Time) – Epochs to sample the orbit positions.

  • plane (Planes, optional) – Fundamental plane of the frame, default to Earth Equator.

sample(epochs=None, *, interpolator=SplineInterpolator())

Returns coordinates at specified epochs.

Parameters
  • epochs (Time, optional) – Epochs to sample the ephemerides, if not given the original one from the object will be used.

  • interpolator (BaseInterpolator, optional) – Interpolation method to use for epochs outside of the original ones, default to splines.

Returns

Sampled coordinates with velocities.

Return type

CartesianRepresentation

rv(epochs=None, **kwargs)

Position and velocity vectors at given epochs.

Parameters
  • epochs (Time, optional) – Epochs to sample the ephemerides, default to now.

  • **kwargs – Extra kwargs for interpolation method.

class poliastro.maneuver.Maneuver(*impulses)

Class to represent a Maneuver.

Each Maneuver consists on a list of impulses \(\Delta v_i\) (changes in velocity) each one applied at a certain instant \(t_i\). You can access them directly indexing the Maneuver object itself.

>>> man = Maneuver((0 * u.s, [1, 0, 0] * u.km / u.s),
... (10 * u.s, [1, 0, 0] * u.km / u.s))
>>> man[0]
(<Quantity 0. s>, <Quantity [1., 0., 0.] km / s>)
>>> man.impulses[1]
(<Quantity 10. s>, <Quantity [1., 0., 0.] km / s>)
classmethod impulse(dv)

Single impulse at current time.

Parameters

dv (numpy.ndarray) – Velocity components of the impulse.

classmethod hohmann(orbit_i, r_f)

Compute a Hohmann transfer between two circular orbits.

Parameters
  • orbit_i (poliastro.twobody.orbit.Orbit) – Initial orbit

  • r_f (astropy.unit.Quantity) – Final orbital radius

classmethod bielliptic(orbit_i, r_b, r_f)

Compute a bielliptic transfer between two circular orbits.

Parameters
  • orbit_i (poliastro.twobody.orbit.Orbit) – Initial orbit

  • r_b (astropy.unit.Quantity) – Altitude of the intermediate orbit

  • r_f (astropy.unit.Quantity) – Final orbital radius

classmethod lambert(orbit_i, orbit_f, method=lambert_izzo, **kwargs)

Computes Lambert maneuver between two different points.

Parameters
  • orbit_i (Orbit) – Initial orbit

  • orbit_f (Orbit) – Final orbit

  • method (function) – Method for solving Lambert’s problem

  • **kwargs – Extra kwargs for Lambert method.

You can read the complete reference of the high level API here:

Core API

The core API is a low level layer that contains simple functions. They are accelerated using Numba, a Just-in-Time compiler for Python, to achieve good performance. However, they take raw NumPy arrays and Python scalars, so they will not protect you from dimensional errors.