Trajectory propagation example

This example shows how to propagated and plot a trajectory using the GODOT cosmos module.

First we load some python packages for plotting and linear algebra.

In [1]:
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

Then we load godot specific libraries

In [2]:
from godot.core import tempo
from godot import cosmos

# optionally avoid verbose logging messages
import godot.core.util as util
util.suppressLogger()

Load the universe configuration and create the universe object

In [3]:
uniConfig = cosmos.util.load_yaml('universe.yml')
uni = cosmos.Universe(uniConfig)

Load the trajectory configuration and create the trajectory object using the universe object

In [4]:
traConfig = cosmos.util.load_yaml('trajectory.yml')
tra = cosmos.Trajectory(uni, traConfig)

The trajectory can be evaluated: in this phase the timeline elements are processed and propagation arcs are computed

In [5]:
tra.compute(partials = False)

After evaluation the trajectory, we can request the solutions of specific timeline elements

In [6]:
sol = tra.getTimelineSolution()
for entry in sol:
    for item in entry:
        print('%24s%15.6f' % (item.name, item.epoch.mjd()))

                  launch       0.000000
              man1_start       0.514660
                man1_end       0.622055
             match1_left       0.900795
            match1_right       1.207151
                man2_end       1.455459
              man2_start       1.501376
                 arrival       2.000000

For each propagation arc, we generate some plotting data for the position vector, propagated mass and accumulated delta-V

In [7]:
data = list()
for entry in sol:
    for left, right in zip(entry, entry[1:]):
        grid = tempo.EpochRange(left.epoch, right.epoch).createGrid(1e-3 * 86400)
        t = [e.mjd() for e in grid]
        x = np.asarray([uni.frames.vector3('Earth', 'GeoSat_center', 'ICRF', e) for e in grid])
        m = [uni.evaluables.get('GeoSat_mass').eval(e) for e in grid]
        v = [uni.evaluables.get('GeoSat_dv').eval(e) for e in grid]
        data.append((t, x, m, v))

Let's plot the X-Y components of the evolution of the position vector

In [8]:
plt.grid()
plt.xlabel('X pojection [km]')
plt.ylabel('Y pojection [km]')
for i, entry in enumerate(data):
    c = [u'b', u'g', u'r', u'c', u'm', u'y', u'k'][i]
    plt.plot([entry[1][0, 0]], [entry[1][0, 1]], 'o', color=c)
    plt.plot(entry[1][:, 0], entry[1][:, 1], '-', color=c)
    plt.plot([entry[1][-1, 0]], [entry[1][-1, 1]], 'o', color=c)

Let's plot the 3D Earth-centered trajectory

In [9]:
fig = plt.figure()
ax = fig.add_subplot(1, 1, 1, projection='3d')
ax.grid()
for i, entry in enumerate(data):
    c = [u'b', u'g', u'r', u'c', u'm', u'y', u'k'][i]
    ax.plot([entry[1][0, 0]], [entry[1][0, 1]], [entry[1][0, 2]], 'o', color=c)
    ax.plot(entry[1][:, 0], entry[1][:, 1], entry[1][:, 2], '-', color=c)
    ax.plot([entry[1][-1, 0]], [entry[1][-1, 1]], [entry[1][-1, 2]], 'o', color=c)

Let's plot the mass evolution

In [10]:
plt.grid()
plt.xlabel('Time [mjd]')
plt.ylabel('Mass [kg]')
for i, entry in enumerate(data):
    c = [u'b', u'g', u'r', u'c', u'm', u'y', u'k'][i]
    plt.plot(entry[0][0], entry[2][0], 'o', color=c)
    plt.plot(entry[0], entry[2], '-', color=c)
    plt.plot(entry[0][-1], entry[2][-1], 'o', color=c)

Let's plot the accumulated delta-V evolution

In [11]:
plt.grid()
plt.xlabel('Time [mjd]')
plt.ylabel('Delta-V [km/s]')
for i, entry in enumerate(data):
    c = [u'b', u'g', u'r', u'c', u'm', u'y', u'k'][i]
    plt.plot(entry[0][0], entry[3][0], 'o', color=c)
    plt.plot(entry[0], entry[3], '-', color=c)
    plt.plot(entry[0][-1], entry[3][-1], 'o', color=c)