API Reference

poliastro.twobody package

class poliastro.twobody.State

Class to represent the position of a body wrt to an attractor.

static from_vectors(attractor, r, v, epoch=<Time object: scale='utc' format='jyear_str' value=J2000.000>)

Return State object 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.
static from_classical(attractor, a, ecc, inc, raan, argp, nu, epoch=<Time object: scale='utc' format='jyear_str' value=J2000.000>)

Return State object from classical orbital elements.

Parameters:
  • attractor (Body) – Main attractor.
  • p (Quantity) – Semilatus rectum.
  • 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.
static from_equinoctial(attractor, p, f, g, h, k, L, epoch=<Time object: scale='utc' format='jyear_str' value=J2000.000>)

Return State object 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.
classmethod circular(attractor, alt, inc=<Quantity 0.0 deg>, raan=<Quantity 0.0 deg>, arglat=<Quantity 0.0 deg>, epoch=<Time object: scale='utc' format='jyear_str' value=J2000.000>)

Return State corresponding to a 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.
classmethod parabolic(attractor, p, inc, raan, argp, nu, epoch=<Time object: scale='utc' format='jyear_str' value=J2000.000>)

Return State corresponding to 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.
r

Position vector.

v

Velocity vector.

a

Semimajor axis.

p

Semilatus rectum.

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.

rv()

Position and velocity vectors.

coe()

Classical orbital elements.

pqw()

Perifocal frame (PQW) vectors.

propagate(time_of_flight, rtol=1e-10)

Propagate this State some time and return the result.

apply_maneuver(maneuver, intermediate=False)

Returns resulting State 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.
to_vectors()

Converts to position and velocity vector representation.

to_classical()

Converts to classical orbital elements representation.

to_equinoctial()

Converts to modified equinoctial elements representation.

poliastro.twobody.propagation module

Propagation algorithms.

poliastro.twobody.propagation.func_twobody(t0, u, k, ad)

Differential equation for the initial value two body problem.

This function follows Cowell’s formulation.

Parameters:
  • t0 (float) – Time.
  • u (ndarray) – Six component state vector [x, y, z, vx, vy, vz] (km, km/s).
  • k (float) – Standard gravitational parameter.
  • ad (function(t0, u, k)) – Non Keplerian acceleration (km/s2).
poliastro.twobody.propagation.cowell(k, r0, v0, tof, ad=None, nsteps=1000, rtol=1e-10)

Propagates orbit using Cowell’s formulation.

Parameters:
  • k (float) – Gravitational constant of main attractor (km^3 / s^2).
  • r0 (array) – Initial position (km).
  • v0 (array) – Initial velocity (km).
  • ad (function(t0, u, k), optional) – Non Keplerian acceleration (km/s2), default to None.
  • tof (float) – Time of flight (s).
  • nsteps (int, optional) – Maximum number of internal steps, default to 1000.
  • rtol (float, optional) – Maximum relative error permitted, default to 1e-10.
Raises:

RuntimeError – If the algorithm didn’t converge.

Notes

This method uses a Dormand & Prince method of order 8(5,3) available in the scipy.integrate.ode module.

poliastro.twobody.propagation.kepler(k, r0, v0, tof, numiter=35, rtol=1e-10)

Propagates Keplerian orbit.

Parameters:
  • k (float) – Gravitational constant of main attractor (km^3 / s^2).
  • r0 (array) – Initial position (km).
  • v0 (array) – Initial velocity (km).
  • tof (float) – Time of flight (s).
  • numiter (int, optional) – Maximum number of iterations, default to 35.
  • rtol (float, optional) – Maximum relative error permitted, default to 1e-10.
Raises:

RuntimeError – If the algorithm didn’t converge.

Notes

This algorithm is based on Vallado implementation, and does basic Newton iteration on the Kepler equation written using universal variables. Battin claims his algorithm uses the same amount of memory but is between 40 % and 85 % faster.

poliastro.twobody.angles module

Angles and anomalies.

poliastro.twobody.angles.nu_to_E(nu, ecc)

Eccentric anomaly from true anomaly.

New in version 0.4.0.

Parameters:
  • nu (float) – True anomaly (rad).
  • ecc (float) – Eccentricity.
Returns:

E – Eccentric anomaly.

Return type:

float

poliastro.twobody.angles.E_to_nu(E, ecc)

True anomaly from eccentric anomaly.

New in version 0.4.0.

Parameters:
  • E (float) – Eccentric anomaly (rad).
  • ecc (float) – Eccentricity.
Returns:

nu – True anomaly (rad).

Return type:

float

poliastro.twobody.angles.M_to_E(M, ecc)

Eccentric anomaly from true anomaly.

New in version 0.4.0.

Parameters:
  • M (float) – Mean anomaly (rad).
  • ecc (float) – Eccentricity.
Returns:

E – Eccentric anomaly.

Return type:

float

poliastro.twobody.angles.E_to_M(E, ecc)

Mean anomaly from eccentric anomaly.

New in version 0.4.0.

Parameters:
  • E (float) – Eccentric anomaly (rad).
  • ecc (float) – Eccentricity.
Returns:

M – Mean anomaly (rad).

Return type:

float

poliastro.twobody.angles.M_to_nu(M, ecc)

True anomaly from mean anomaly.

New in version 0.4.0.

Parameters:
  • M (float) – Mean anomaly (rad).
  • ecc (float) – Eccentricity.
Returns:

nu – True anomaly (rad).

Return type:

float

Examples

>>> nu = M2nu(np.radians(30.0), 0.06)
>>> np.rad2deg(nu)
33.673284930211658
poliastro.twobody.angles.nu_to_M(nu, ecc)

Mean anomaly from true anomaly.

New in version 0.4.0.

Parameters:
  • nu (float) – True anomaly (rad).
  • ecc (float) – Eccentricity.
Returns:

M – Mean anomaly (rad).

Return type:

float

poliastro.twobody.angles.fp_angle(nu, ecc)

Flight path angle.

New in version 0.4.0.

Parameters:
  • nu (float) – True anomaly (rad).
  • ecc (float) – Eccentricity.

Notes

Algorithm taken from Vallado 2007, pp. 113.

poliastro.maneuver module

Orbital maneuvers.

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>)
__init__(*impulses)

Constructor.

Parameters:impulses (list) – List of pairs (delta_time, delta_velocity)

Notes

TODO: Fix docstring, *args convention

classmethod impulse(dv)

Single impulse at current time.

classmethod hohmann(ss_i, r_f)

Compute a Hohmann transfer between two circular orbits.

classmethod bielliptic(ss_i, r_b, r_f)

Compute a bielliptic transfer between two circular orbits.

get_total_time()

Returns total time of the maneuver.

get_total_cost()

Returns total cost of the maneuver.

poliastro.bodies module

Bodies of the Solar System.

Contains some predefined bodies of the Solar System:

  • Sun (☉)
  • Earth (♁)

and a way to define new bodies (Body class).

TODO: Add more parameters (e.g. J2)

class poliastro.bodies.Body(k, name=None, symbol=None, R=<Quantity 0.0 km>)

Class to represent a body of the Solar System.

__init__(k, name=None, symbol=None, R=<Quantity 0.0 km>)

Constructor.

Parameters:
  • k (Quantity) – Standard gravitational parameter
  • name (str, optional) – Name of the body, default to None.
  • symbol (str, optional) – Symbol for the body, default to None.
  • R (Quantity, optional) – Radius of the body, default to 0 km.

poliastro.plotting module

Plotting utilities.

poliastro.plotting.plot(state, label=None)

Plots a State.

For more advanced tuning, use the OrbitPlotter class.

class poliastro.plotting.OrbitPlotter(ax=None, num_points=100)

OrbitPlotter class.

This class holds the perifocal plane of the first State plotted in it using plot(), so all following plots will be projected on that plane. Alternatively, you can call set_frame() to set the frame before plotting.

__init__(ax=None, num_points=100)

Constructor.

Parameters:
  • ax (Axes) – Axes in which to plot. If not given, new ones will be created.
  • num_points (int, optional) – Number of points to use in plots, default to 100.
set_frame(p_vec, q_vec, w_vec)

Sets perifocal frame if not existing.

Raises:
  • ValueError – If the vectors are not a set of mutually orthogonal unit vectors.
  • NotImplementedError – If the frame was already set up.
plot(state, osculating=True, label=None)

Plots state and osculating orbit in their plane.

poliastro.iod module

poliastro.ephem module

poliastro.stumpff module

Stumpff functions.

poliastro.stumpff.c2

Second Stumpff function.

For positive arguments:

\[c_2(\psi) = \frac{1 - \cos{\sqrt{\psi}}}{\psi}\]
poliastro.stumpff.c3

Third Stumpff function.

For positive arguments:

\[c_3(\psi) = \frac{\sqrt{\psi} - \sin{\sqrt{\psi}}}{\sqrt{\psi^3}}\]

poliastro.util module

Function helpers.

poliastro.util.rotate(vector, angle, axis='z', unit=None)

Rotates a vector around axis a right-handed positive angle.

This is just a convenience function around astropy.coordinates.angles.rotation_matrix.

Parameters:
  • vector (array_like) – Dimension 3 vector.
  • angle (convertible to Angle) – Angle of rotation.
  • axis (str or 3-sequence) – Either ‘x’,’y’, ‘z’, or a (x,y,z) specifying an axis to rotate about. If ‘x’,’y’, or ‘z’, the rotation sense is counterclockwise looking down the + axis (e.g. positive rotations obey left-hand-rule).
  • unit (UnitBase, optional) – If angle does not have associated units, they are in this unit. If neither are provided, it is assumed to be degrees.

Notes

This is just a convenience function around astropy.coordinates.angles.rotation_matrix. This performs a so-called active or alibi transformation: rotates the vector while the coordinate system remains unchanged. To do the opposite operation (passive or alias transformation) call the function as rotate(vec, ax, -angle, unit) or use the convenience function transform(), see [1].

References

[1]http://en.wikipedia.org/wiki/Rotation_matrix#Ambiguities
poliastro.util.transform(vector, angle, axis='z', unit=None)

Rotates a coordinate system around axis a positive right-handed angle.

Notes

This is a convenience function, equivalent to rotate(vec, -angle, axis, unit). Refer to the documentation of rotate() for further information.

poliastro.util.check_units(quantities, units)

Check if list of quantities is consistent with given list of units.

poliastro.util.norm(vec)

Norm of a Quantity vector that respects units.