I haven't had much time form Orbiter of late, but I've been musing on how to make use of PyKEP / PyGMO output within Orbiter. PyKEP / PyGMO is a free mission planning tool from the European Space Agency. However, it knows nothing about Orbiter and produces optimal trajectory plans in ways that are, most decidedly, Orbiter-unfriendly. To make practical use of PyKEP / PyGMO within Orbiter, then, one needs a way of converting PyKEP optimal trajectory plans into something that can be implemented in Orbiter using widely used MFDS - such as TransX and IMFD.
The purpose of this note, then, is to capture a few thoughts on these transformations. This is as much as an 'aide de memoire' for me as anything else since on more than one occasion, I have had recourse to rifling through my own posts to remind myself of useful equations and conversions.
Why convert from PyKEP-speak to Orbiter-speak?
PyKEP works by optimising a series of 'lambert arcs' that connect a series of planetary encounters. For example, for a simple Hohmann transfer from Earth to Jupiter, say, the two planetary encounters are: a departure from Earth; and an arrival at Jupiter. When PyKEP has completed its optimisation task, usually by minimising total [MATH]\Delta V[/MATH] requirements, it reports back the optimal times of the planetary encounters - e.g., the departure date from Earth; and the arrival date at Jupiter. However, PyKEP doesn't provide any useful information about how to set up the Earth 'escape plan'. Now, one can use TransX / IMFD and, by a process of trial-and-error, one can back-solve the encounter dates to work out what the 'escape plan' has to be to arrive at Jupiter on the planned arrival date. For more complicated trajectory plans involving a long sequence of slingshots and DSMs, using this trial-and-error technique to re-construct PyKEP's optimal trajectory solution in Orbiter can be challenging, if not impossible. Surely, there must be a better way than trial-and-error.
What output does PyKEP / PyGMO produce?
Let's start with the information that PyKEP does produce. First and foremost, it yields the times of the encounters with the planetary bodies. But PyKEP also has an ephemeris engine so it knows the position and velocity (the 'state vectors') of the planetary bodies. Now, PyKEP uses JPL's Low Precision ephemeris whereas Orbiter uses the VSOP87 ephemeris solution. Although these ephemerides are not exactly the same, for practical purposes the two ephemerides can be treated as equivalent.
Given that PyKEP knows 'where' the planetary bodies are, there is a function in PyKEP that allows one to re-construct the lambert arc connecting each encounter in the sequence. The lambert-arc is generally just a ballistic elliptical trajectory that connects the two planetary encounters and so it corresponds to the trajectory that a spacecraft needs to take in order to go from one encounter to the next in the specified amount of time. Part of the reconstruction of the Lambert arc is information about the velocity that the spacecraft must have when it leaves one planetary encounter; and the velocity that the spacecraft will have when it arrives at the next planetary encounter.
So, PyKEP provides the following information:
* the position of the spacecraft (and planetary body) at each encounter;
* the velocity of the planetary body at each encounter;
* the velocity of the spacecraft as it arrives at a planetary encounter; and
* the velocity of the spacecraft as it leaves a planetary encounter.
From this information we need to construct information that can be used in Orbiter's mission tools - e.g. TransX and IMFD. (In this note, I'm going to focus on TransX, since that tool is better suited to setting up slingshots than IMFD).
As an example from a PyKEP optimisation of a VEGA-style trajectory (Venus-Earth Gravity Assist), after extracting the relevant information from the lambert arcs, PyKEP tells me that at the Venus encounter (MJD 58002.0765):
* the position of Venus and the spacecraft (in a heliocentric co-ordinate system) is:
[MATH]r_v=(5070159035.023689, 107556835822.169, 1182577911.1819444)[/MATH]
* the velocity of Venus is:
[MATH] v_v=(-35101.3332441319, 1466.4390915935328, 2045.7339135134034)[/MATH]
* the spacecraft's approach velocity vector (heliocentric coordinates again) is:
[MATH] v_{in}=(-36817.31294777, -5469.38165681, -195.38174684)[/MATH]
* the spacecraft's departure velocity vector (heliocentric coordinates) is:
[MATH] v_{out} = (-41823.78381252, -1170.11804522, 63.30610391)[/MATH]
(In case you are wondering, the positions are the x-y-z coordinates in metres; and the velocities are the x-y-z coordinates in metres per second.)
From this information, we can accurately and reliably construct the 'escape plan' (if we are on Venus and leaving it); and the slingshot (if we are arriving from somewhere else and want to head of to another planetary destination). How do we do that? Well, let's begin by talking about TransX's coordinate system.
TransX's coordinate system - escape plans
Rather than use the x-y-z coordinate system of the underlying ephemeris, TransX works in another coordinate system consisting of 'prograde', 'outward' and 'plane'. These three quantities can be thought of as a new x-y-z coordinate system - in that the coordinates are all at right-angles to each other and form a dextral coordinate system in the same way that the underlying heliocentric x-y-z coordinate system does. To set up an escape plan, we need first to construct the unit vectors corresponding to TransX's 'prograde', 'outward' and 'plane' triplet of directions. Let's call these unit vectors [MATH]\hat{f}[/MATH], [MATH]\hat{o}[/MATH] and [MATH]\hat{p}[/MATH]. Here 'f' (or 'forward') stands for prograde to avoid confusion with using 'p' for 'plane'. These vectors are constructed as follows:
[MATH] \hat{f} = \frac{v_v}{\left\| v_v\right\| }[/MATH]
[MATH] \hat{p} = \frac{v_v \times r_v}{\left\| v_v \times r_v \right\|}[/MATH]
[MATH] \hat{o} = \hat{p} \times \hat{f}[/MATH]
So, using the values given above:
[MATH] \hat{f} = (-0.99743888, 0.04167031, 0.05813154) [/MATH]
[MATH] \hat{o} = ( 0.04239308, 0.99903761, 0.01125561) [/MATH]
[MATH] \hat{p} = (-0.05760657, 0.01369116, -0.99824548) [/MATH]
Having calculated these unit vectors, we can now calculate the information needed to form the escape plan. First, we need to work out the departure speed of the spacecraft relative to Venus:
[MATH] v_{out}' = v_{out} - v_v [/MATH]
In the above example, we calculate the spacecraft's departure velocity - relative to Venus - to be:
[MATH] v_{out}' = (-6722.45056839 -2636.55713681 -1982.4278096) [/MATH]
Now, we can always write $v_{out}'$ as:
[MATH] v_{out}' = v_f\,\hat{f} + v_o\,\hat{o} + v_p\,\hat{p} [/MATH]
where the values [MATH]v_f[/MATH], [MATH]v_o[/MATH] and [MATH]v_p[/MATH] are just the components of our Venus escape plan. We calculate these quantities by taking the 'dot product' of [MATH]v_{out}'[/MATH] with each of the TransX unit vectors as follows:
[MATH] v_f = v_{out}'.\hat{f} [/MATH]
[MATH] v_o = v_{out}'.\hat{o} [/MATH]
[MATH] v_p = v_{out}'.\hat{p} [/MATH]
If we apply this to our example, we calculate (in m/s) the following TransX escape plan:
[MATH] v_f = 6480.126 [/MATH]
[MATH] v_o = -2941.319 [/MATH]
[MATH] v_p = 2330.109 [/MATH]
And, to be sure, if starting from Venus one sets up this escape plan for MJD 58002.0765 we should find that TransX will predict a close encounter with Earth at the next scheduled encounter date of MJD 58572.5770.
The calculations here are quite general and can be applied to any PyKEP trajectory solution.
TransX's coordinate system - slingshots
Escape plans are generally used at the start of a sequence of encounters in TransX. In TransX, other planetary encounters (prior to arrival at the destination) are called 'slingshots' - aka 'gravity assists' and 'flybys'. Rather than using the cartesian coordinate system consisting of 'prograde', 'outward' and 'plane', TransX uses a spherical coordinate system for slingshots. Again, three components are needed: [MATH]v_\infty[/MATH] ('hyperbolic escape velocity'), [MATH]\theta[/MATH] ('outward angle') and [MATH]\phi[/MATH] ('inclination'). The 'hyperbolic escape velocity' is the speed of the spacecraft (relative to the planet) once it has escaped the planet's gravity well; the 'outward angle' is the angle between the prograde direction and the direction of departure from the planet (in the orbital plane of planet); and 'inclination' is the angle that the direction of departure makes with the planet's orbital plane. This sounds more complicated than it is in reality.
How do we calculate these quantities? Using the same notation as above:
[MATH] v_\infty = \left\| v_{out}' \right\|[/MATH]
[MATH]\theta = \tan ^{-1}\left(v_x,v_y\right)[/MATH]
[MATH]\phi = \tan^{-1}\left(\sqrt{v_x^2+v_y^2},v_z\right)[/MATH]
(Here, the arctan functions are of the 'atan2' variety - requiring both an 'x' and a 'y' input.) Normally in TransX, the slingshot 'inherits' the departure hyperbolic excess velocity from the arrival hyperbolic excess velocity. However, and in particular for powered flybys, one can specify a hyperbolic excess velocity that is faster or slower than the arrival hyperbolic excess velocity.
In our example, if we work though the numbers, we find that:
[MATH] v_\infty = 7488.177 \quad \text{m/s}[/MATH]
[MATH]\theta = -24.413 \quad \text{deg}[/MATH]
[MATH]\phi = 18.130 \quad \text{deg}[/MATH]
And, upon approach to Venus, if one inserts these numbers into TransX's slingshot plan (again for Venus arrival at MJD 58002.0765. The departure trajectory should again coincide with Earth arrival around MJD 58572.5770. However, this isn't as easy to set up *de novo* than the escape plan.
Setting up slingshots in this way is especially important if there is a DSM mid-way between one encounter and another. Without an ability to construct a planet departure trajectory in this fashion is becomes almost impossible (using TransX and IMFD) to ensure that the spacecraft is heading towards the DSM point.
PyKEP's full optimal VEGA trajectory for Earth departure MJD 57634.5001
For those interested, I've translated PyKEP's optimised VEGA-style trajectory plan for Earth departure MJD 57634.5001. This is based on a list of Jupiter trajectories provided in an earlier post by 'dgatsoulis'.
The plan is as follows:
This plan is perfectly flyable in Orbiter now. Total delta-V cost of getting to Jupiter (and being captured by it) is 4,507.7 m/s (starting from a 300 km x 300 km orbit around Earth). It involves a ballistic encounter with Venus and a powered flyby of Earth. In comparison, an EGA-style Jupiter approach will typically cost around 5,300 m/s; and a Hohmann-stye transfer to Jupiter close to 7,000 m/s. Journey time is a little over five years.
The purpose of this note, then, is to capture a few thoughts on these transformations. This is as much as an 'aide de memoire' for me as anything else since on more than one occasion, I have had recourse to rifling through my own posts to remind myself of useful equations and conversions.
Why convert from PyKEP-speak to Orbiter-speak?
PyKEP works by optimising a series of 'lambert arcs' that connect a series of planetary encounters. For example, for a simple Hohmann transfer from Earth to Jupiter, say, the two planetary encounters are: a departure from Earth; and an arrival at Jupiter. When PyKEP has completed its optimisation task, usually by minimising total [MATH]\Delta V[/MATH] requirements, it reports back the optimal times of the planetary encounters - e.g., the departure date from Earth; and the arrival date at Jupiter. However, PyKEP doesn't provide any useful information about how to set up the Earth 'escape plan'. Now, one can use TransX / IMFD and, by a process of trial-and-error, one can back-solve the encounter dates to work out what the 'escape plan' has to be to arrive at Jupiter on the planned arrival date. For more complicated trajectory plans involving a long sequence of slingshots and DSMs, using this trial-and-error technique to re-construct PyKEP's optimal trajectory solution in Orbiter can be challenging, if not impossible. Surely, there must be a better way than trial-and-error.
What output does PyKEP / PyGMO produce?
Let's start with the information that PyKEP does produce. First and foremost, it yields the times of the encounters with the planetary bodies. But PyKEP also has an ephemeris engine so it knows the position and velocity (the 'state vectors') of the planetary bodies. Now, PyKEP uses JPL's Low Precision ephemeris whereas Orbiter uses the VSOP87 ephemeris solution. Although these ephemerides are not exactly the same, for practical purposes the two ephemerides can be treated as equivalent.
Given that PyKEP knows 'where' the planetary bodies are, there is a function in PyKEP that allows one to re-construct the lambert arc connecting each encounter in the sequence. The lambert-arc is generally just a ballistic elliptical trajectory that connects the two planetary encounters and so it corresponds to the trajectory that a spacecraft needs to take in order to go from one encounter to the next in the specified amount of time. Part of the reconstruction of the Lambert arc is information about the velocity that the spacecraft must have when it leaves one planetary encounter; and the velocity that the spacecraft will have when it arrives at the next planetary encounter.
So, PyKEP provides the following information:
* the position of the spacecraft (and planetary body) at each encounter;
* the velocity of the planetary body at each encounter;
* the velocity of the spacecraft as it arrives at a planetary encounter; and
* the velocity of the spacecraft as it leaves a planetary encounter.
From this information we need to construct information that can be used in Orbiter's mission tools - e.g. TransX and IMFD. (In this note, I'm going to focus on TransX, since that tool is better suited to setting up slingshots than IMFD).
As an example from a PyKEP optimisation of a VEGA-style trajectory (Venus-Earth Gravity Assist), after extracting the relevant information from the lambert arcs, PyKEP tells me that at the Venus encounter (MJD 58002.0765):
* the position of Venus and the spacecraft (in a heliocentric co-ordinate system) is:
[MATH]r_v=(5070159035.023689, 107556835822.169, 1182577911.1819444)[/MATH]
* the velocity of Venus is:
[MATH] v_v=(-35101.3332441319, 1466.4390915935328, 2045.7339135134034)[/MATH]
* the spacecraft's approach velocity vector (heliocentric coordinates again) is:
[MATH] v_{in}=(-36817.31294777, -5469.38165681, -195.38174684)[/MATH]
* the spacecraft's departure velocity vector (heliocentric coordinates) is:
[MATH] v_{out} = (-41823.78381252, -1170.11804522, 63.30610391)[/MATH]
(In case you are wondering, the positions are the x-y-z coordinates in metres; and the velocities are the x-y-z coordinates in metres per second.)
From this information, we can accurately and reliably construct the 'escape plan' (if we are on Venus and leaving it); and the slingshot (if we are arriving from somewhere else and want to head of to another planetary destination). How do we do that? Well, let's begin by talking about TransX's coordinate system.
TransX's coordinate system - escape plans
Rather than use the x-y-z coordinate system of the underlying ephemeris, TransX works in another coordinate system consisting of 'prograde', 'outward' and 'plane'. These three quantities can be thought of as a new x-y-z coordinate system - in that the coordinates are all at right-angles to each other and form a dextral coordinate system in the same way that the underlying heliocentric x-y-z coordinate system does. To set up an escape plan, we need first to construct the unit vectors corresponding to TransX's 'prograde', 'outward' and 'plane' triplet of directions. Let's call these unit vectors [MATH]\hat{f}[/MATH], [MATH]\hat{o}[/MATH] and [MATH]\hat{p}[/MATH]. Here 'f' (or 'forward') stands for prograde to avoid confusion with using 'p' for 'plane'. These vectors are constructed as follows:
[MATH] \hat{f} = \frac{v_v}{\left\| v_v\right\| }[/MATH]
[MATH] \hat{p} = \frac{v_v \times r_v}{\left\| v_v \times r_v \right\|}[/MATH]
[MATH] \hat{o} = \hat{p} \times \hat{f}[/MATH]
So, using the values given above:
[MATH] \hat{f} = (-0.99743888, 0.04167031, 0.05813154) [/MATH]
[MATH] \hat{o} = ( 0.04239308, 0.99903761, 0.01125561) [/MATH]
[MATH] \hat{p} = (-0.05760657, 0.01369116, -0.99824548) [/MATH]
Having calculated these unit vectors, we can now calculate the information needed to form the escape plan. First, we need to work out the departure speed of the spacecraft relative to Venus:
[MATH] v_{out}' = v_{out} - v_v [/MATH]
In the above example, we calculate the spacecraft's departure velocity - relative to Venus - to be:
[MATH] v_{out}' = (-6722.45056839 -2636.55713681 -1982.4278096) [/MATH]
Now, we can always write $v_{out}'$ as:
[MATH] v_{out}' = v_f\,\hat{f} + v_o\,\hat{o} + v_p\,\hat{p} [/MATH]
where the values [MATH]v_f[/MATH], [MATH]v_o[/MATH] and [MATH]v_p[/MATH] are just the components of our Venus escape plan. We calculate these quantities by taking the 'dot product' of [MATH]v_{out}'[/MATH] with each of the TransX unit vectors as follows:
[MATH] v_f = v_{out}'.\hat{f} [/MATH]
[MATH] v_o = v_{out}'.\hat{o} [/MATH]
[MATH] v_p = v_{out}'.\hat{p} [/MATH]
If we apply this to our example, we calculate (in m/s) the following TransX escape plan:
[MATH] v_f = 6480.126 [/MATH]
[MATH] v_o = -2941.319 [/MATH]
[MATH] v_p = 2330.109 [/MATH]
And, to be sure, if starting from Venus one sets up this escape plan for MJD 58002.0765 we should find that TransX will predict a close encounter with Earth at the next scheduled encounter date of MJD 58572.5770.
The calculations here are quite general and can be applied to any PyKEP trajectory solution.
TransX's coordinate system - slingshots
Escape plans are generally used at the start of a sequence of encounters in TransX. In TransX, other planetary encounters (prior to arrival at the destination) are called 'slingshots' - aka 'gravity assists' and 'flybys'. Rather than using the cartesian coordinate system consisting of 'prograde', 'outward' and 'plane', TransX uses a spherical coordinate system for slingshots. Again, three components are needed: [MATH]v_\infty[/MATH] ('hyperbolic escape velocity'), [MATH]\theta[/MATH] ('outward angle') and [MATH]\phi[/MATH] ('inclination'). The 'hyperbolic escape velocity' is the speed of the spacecraft (relative to the planet) once it has escaped the planet's gravity well; the 'outward angle' is the angle between the prograde direction and the direction of departure from the planet (in the orbital plane of planet); and 'inclination' is the angle that the direction of departure makes with the planet's orbital plane. This sounds more complicated than it is in reality.
How do we calculate these quantities? Using the same notation as above:
[MATH] v_\infty = \left\| v_{out}' \right\|[/MATH]
[MATH]\theta = \tan ^{-1}\left(v_x,v_y\right)[/MATH]
[MATH]\phi = \tan^{-1}\left(\sqrt{v_x^2+v_y^2},v_z\right)[/MATH]
(Here, the arctan functions are of the 'atan2' variety - requiring both an 'x' and a 'y' input.) Normally in TransX, the slingshot 'inherits' the departure hyperbolic excess velocity from the arrival hyperbolic excess velocity. However, and in particular for powered flybys, one can specify a hyperbolic excess velocity that is faster or slower than the arrival hyperbolic excess velocity.
In our example, if we work though the numbers, we find that:
[MATH] v_\infty = 7488.177 \quad \text{m/s}[/MATH]
[MATH]\theta = -24.413 \quad \text{deg}[/MATH]
[MATH]\phi = 18.130 \quad \text{deg}[/MATH]
And, upon approach to Venus, if one inserts these numbers into TransX's slingshot plan (again for Venus arrival at MJD 58002.0765. The departure trajectory should again coincide with Earth arrival around MJD 58572.5770. However, this isn't as easy to set up *de novo* than the escape plan.
Setting up slingshots in this way is especially important if there is a DSM mid-way between one encounter and another. Without an ability to construct a planet departure trajectory in this fashion is becomes almost impossible (using TransX and IMFD) to ensure that the spacecraft is heading towards the DSM point.
PyKEP's full optimal VEGA trajectory for Earth departure MJD 57634.5001
For those interested, I've translated PyKEP's optimised VEGA-style trajectory plan for Earth departure MJD 57634.5001. This is based on a list of Jupiter trajectories provided in an earlier post by 'dgatsoulis'.
The plan is as follows:
Code:
Date of Earth departue: 2016-Sep-03 12:00:08.640002
Date of Venus encounter: 2017-Sep-06 01:50:05.290876
Date of Earth encounter: 2019-Mar-30 13:50:54.766147
Date of Jupiter arrival: 2021-Sep-24 07:39:52.717700
Transfer time from Earth to Venus: 367.58 days
Transfer time from Venus to Earth: 570.5 days
Transfer time from Earth to Jupiter: 908.74 days
Total mission duration: 1846.82 days
TransX escape plan - Earth escape
--------------------------------------
MJD: 57634.5001
Prograde: -3461.680 m/s
Outward: -291.288 m/s
Plane: -298.955 m/s
Hyp. excess velocity: 3486.754 m/s
Earth escape burn: 3743.016 m/s
Venus encounter
--------------------------------------
MJD: 58002.0765
Approach velocity: 7488.177 m/s
Departure velocity: 7488.177 m/s
Outward angle: -24.413 deg
Inclination: 18.130 deg
Turning angle: 52.332 deg
Periapsis altitude: 1292.475 km
dV needed: 0.000 m/s
Earth encounter
--------------------------------------
MJD: 58572.5770
Approach velocity: 12080.335 m/s
Departure velocity: 12425.696 m/s
Outward angle: 54.099 deg
Inclination: 4.666 deg
dV needed for flyby: 442.676 m/s
Jupiter arrival
--------------------------------------
MJD: 59481.3194
Hyp. excess velocity: 6193.941 m/s
Capture burn 322.023 m/s - (C3 = 0)
Total fuel cost: 4507.714 m/sThis plan is perfectly flyable in Orbiter now. Total delta-V cost of getting to Jupiter (and being captured by it) is 4,507.7 m/s (starting from a 300 km x 300 km orbit around Earth). It involves a ballistic encounter with Venus and a powered flyby of Earth. In comparison, an EGA-style Jupiter approach will typically cost around 5,300 m/s; and a Hohmann-stye transfer to Jupiter close to 7,000 m/s. Journey time is a little over five years.
Last edited:
