API Reference

poliastro.twobody module

class poliastro.twobody.State

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

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

Return State object from orbital elements.

Parameters:
  • attractor (Body) – Main attractor.
  • a (Quantity) – Semimajor 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.
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) – Semi-latus 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.
elements

Classical orbital elements.

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.

period

Period of the orbit.

n

Mean motion.

h

Specific angular momentum.

e_vec

Eccentricity vector.

h_vec

Specific angular momentum vector.

rv()

Position and velocity vectors.

pqw()

Perifocal frame (PQW) vectors.

propagate(time_of_flight)

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.

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)

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)

Plots state and osculating orbit in their plane.

poliastro.iod module

Initial orbit determination.

poliastro.iod.lambert(k, r0, r, tof, short=True, numiter=35, rtol=1e-08)

Solves the Lambert problem.

New in version 0.3.0.

Parameters:
  • k (Quantity) – Gravitational constant of main attractor (km^3 / s^2).
  • r0 (Quantity) – Initial position (km).
  • r (Quantity) – Final position (km).
  • tof (Quantity) – Time of flight (s).
  • short (boolean, optional) – Find out the short path, default to True. If False, find long path.
  • numiter (int, optional) – Maximum number of iterations, default to 35.
  • rtol (float, optional) – Relative tolerance of the algorithm, default to 1e-8.
Raises:

RuntimeError – If it was not possible to compute the orbit.

Note

This uses the universal variable approach found in Battin, Mueller & White with the bisection iteration suggested by Vallado. Multiple revolutions not supported.

poliastro.ephem module

Planetary ephemerides.

poliastro.ephem.select_kernel()

Selects appropriate kernel.

Returns DE421 if available in data directory, else the first kernel found, else None.

New in version 0.3.0.

poliastro.ephem.download_kernel(name)

Downloads SPICE kernel by name.

The function will try the top SPK path first, and then the old versions path in case the .bsp file is not found.

New in version 0.3.0.

poliastro.ephem.planet_ephem(body, epoch, kernel=None)

Position and velocity vectors of a given planet at a certain time.

The vectors are computed with respect to the Solar System barycenter.

New in version 0.3.0.

Parameters:
  • body (int) – Planetary body.
  • epoch (astropy.time.Time) – Computation time. Can be scalar or vector.
  • kernel (jplephem.spk.SPK, optional) – jplephem SPK kernel to make the computation, if not given a default one will be used.
Returns:

r, v – Position and velocity vectors.

Return type:

Quantity

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.

poliastro.util.dot

Returns the dot product of two vectors.