Follow me! https://github.com/Juanlu001/
poliastro is a collection of Python functions useful in Astrodynamics and Orbital Mechanics, focusing on interplanetary applications.
import matplotlib.pyplot as plt
plt.ion()
from plotly.offline import init_notebook_mode
init_notebook_mode(connected=False)
from astropy import units as u
from astropy.time import Time
from poliastro.bodies import Earth
from poliastro.twobody import Orbit
r = [-6045, -3490, 2500] * u.km
v = [-3.457, 6.618, 2.533] * u.km / u.s
ss = Orbit.from_vectors(Earth, r, v, Time.now())
ss
7283 x 10293 km x 153.2 deg orbit around Earth (♁)
from poliastro.plotting import plot, plot3d
plot(ss, label="Sample orbit");
plot3d(ss, label="Sample orbit");
from astropy.coordinates import solar_system_ephemeris
solar_system_ephemeris.set("jpl")
<ScienceState solar_system_ephemeris: 'jpl'>
Orbit.from_body_ephem(Earth)
1 x 1 AU x 23.4 deg orbit around Sun (☉)
from poliastro.neos import neows
florence = neows.orbit_from_name("Florence")
florence
1 x 3 AU x 22.1 deg orbit around Sun (☉)
from poliastro.neos import dastcom5
halley_1835 = dastcom5.orbit_from_name('1P')[30]
halley_1835
1 x 35 AU x 162.3 deg orbit around Sun (☉)
frame = OrbitPlotter()
# Florence close approach
frame.plot(Orbit.from_body_ephem(Earth, epoch=Time("2017-09-01 12:05:50", scale="tdb")), label=Earth)
for planet in Mercury, Venus, Mars:
frame.plot(Orbit.from_body_ephem(planet), label=planet)
frame.plot(halley_1835, label='Halley', color='#666666')
frame.plot(florence, label='Florence', color='#000000')
plt.title("Inner solar system + Florence + Halley")
plt.xlim(-.3e9, .3e9)
plt.ylim(-.3e9, .3e9)
(-300000000.0, 300000000.0)
halleys = dastcom5.orbit_from_name('1P')
frame = OrbitPlotter(num_points=200)
frame.plot(halleys[0], label='Halley')
frame.plot(halleys[5], label='Halley')
frame.plot(halleys[10], label='Halley')
frame.plot(halleys[20], label='Halley')
frame.plot(halleys[30], label='Halley')
[<matplotlib.lines.Line2D at 0x7fab65df2470>, <matplotlib.lines.Line2D at 0x7fab65e5add8>]
from poliastro.twobody.propagation import propagate, cowell
from poliastro.util import norm, time_range
from numba import njit
@njit
def accel(t0, state, k):
"""Constant acceleration aligned with the velocity. """
v_vec = state[3:]
norm_v = (v_vec * v_vec).sum() ** .5
return 1e-5 * v_vec / norm_v
def custom_propagator(orbit, tof, rtol, accel=accel):
# Workaround for https://github.com/poliastro/poliastro/issues/328
if tof == 0:
return orbit.r.to(u.km).value, orbit.v.to(u.km / u.s).value
else:
# Use our custom perturbation acceleration
return cowell(orbit, tof, rtol, ad=accel)
from poliastro.plotting import OrbitPlotter3D
frame = OrbitPlotter3D()
_, rr_earth = ss.sample(
time_range(ss.epoch, end=ss.epoch + 3 * u.day, periods=300),
method=custom_propagator
)
frame.set_attractor(Earth)
frame.plot_trajectory(rr_earth, label=Earth)
frame.show()
The wrooooooooooooooong way:
def sinc(x):
return np.sin(x) / x
import pytest
@pytest.mark.parametrize("x", [0, 1, 10])
def test_sinc(x):
assert sinc(x) == np.sin(x) / x
0.1 + 0.2 == 0.3
False
0.2 + 0.3 == 0.5
True
What's the good way?
...If you're really interested, go read my Final Masters Project: https://github.com/juanlu001/pfc-uc3m
So... let's make our code Fortran-esque!
High level API:
Dangerous™ algorithms:
astropy.units
performance¶astropy.units
does a lot of introspection, which makes it slowI believe the choice of license is an important one, and I advocate a BSD-style license. In my experience, the most important commodity an open source project needs to succeed is users.
-- John Hunter † http://nipy.org/nipy/faq/johns_bsd_pitch.html
Go find your users!