Utilizing the Official SGP4 DLLs in Python
Interacting with TLEs the Right Way
SAA DLL
AFSPC (Air Force SPace Command) distributes the Standard Astrodynamic Library Dynamic Link Library (SAA DLL) in binary form. This library implements the official SGP4 propagation algorithms used by AFSPC to interpret Two Line Element sets (TLEs).
Of course, there are other SGP4 propagation tools. The mathematics behind this system are outlined in Spacetrack Report No 3, and can be implemented by any competent programmer that understands the math. The most famous Python3 model is maintained by Brandon Rhodes and is known simply as sgp4.
However, if you want to exactly match the system maintained by AFSPC and NORAD (and you do want to match them exactly when working with real satellites) your best bet is to utilize the official SAA DLL distributed at Space Track..
Strangely, there is no publically maintained Python3 driver for this purpose. This tutorial explains the basics of how to interact with these DLLs using Python3. This tutorial aims for a "quick and dirty" approach that explains the fundamentals, but I've started an open source project on github.
Note: As usual, this will be a Linux based tutorial. Windows users will need to adjust their approach slightly.
Obtaining SAA DLL
First, the SAA DLLs are regulated under a US Federal Law known as "Export Administration Regulations," which is administrated by the Bureau of Industry and Security. The documentation and executable forms of SAA DLL Version 7 (but not the Version 7 betas or previous versions) are approved for public and unlimited distribution. Even so, the federal government is fickle, so while I will guide you to where to get the SAA DLL, I will not distribute it myself.
If you followed along with my previous blog entry, you already have the DLLs. If not, go to space-track.org, make an account, then go to the sgp4 download page and download the appropriate package for your operating system. For me, this is SGP4_small_V7.9_LINUX64.tar.gz. You should also download the windows package, SGP4_small_V7.9_WIN64.zip, because it provides additional Python2 driver examples.
Reference Documentation
Extract the archives. There are several important reference documents inside.
./README FIRST - STANDARDIZED ASTRODYNAMIC ALGORITHMS (V7.8.1).pdf describes the content and layout of the package. It also describes the validation protocol, which you should complete to make sure your computer is capable of using the DLLs included. More instructions are provided in my prior blog post on this topic.
./Driver_Examples/* contains several examples of how to interface with the library. The Windows package has many more examples than the Linux package, and I felt the Windows Python2 package was the most helpful of the bunch.
./Documents/* contains the API reference documents for the DLLs. As a Linux user, you will need to install a program called xchm to open these .chm files. It should be available via your software respository.
The Theory of DLLs
A DLL is a "Dynamic Link Library." This is a special kind of executable object that other programs can interact with. The "Dynamic" part means the dependent executables only interact with the DLL at run time. This makes it a lot easier to simply swap out the DLLs as needed.
The term "DLL" is generally associated with Windows. In Linux, we use .so, or "Shared Object" files. I have tested both DLLs and SOs and they both seem to work interchangeably. However, it is good practice to point Unix software to .so files and Windows software to .dll files because AFSPC may choose to implement slightly different behavior on a per-platform basis.
DLLs solve three problems for AFSPC.
- AFSPC wants to keep its source code a secret. By distributing compiled DLLs, AFSPC only needs to secure an export license for the compiled DLLs and associated documentation.
- AFSPC wants to avoid duplication of effort. By distributing these DLLs, all of their astrodynamic applications can be built on the same code base. This means that customers, contractors, and government agencies are all on the same page software-wise. If a bug is found, it can be repaired for all users and all applications simultaneously.
- AFSPC needs to distribute bug-fixes and keep a wide variety of code up to date. By using DLLs, the end-users don't need to recompile code to gain the benefits of the bugfixes and updates. They need only swap out the DLL files. In fact, there are installer programs that can manage this for the end-user.
Project Setup
This tutorial aims for the simplest, most easily understood usage of the SAA DLLs. Make a directory called pysaademo/ and structure it like this:
pysaademo/
├── libdll/
├── pysgp4demo.py
├── run.log
└── set_env
Copy the entire contents of the SAA DLL directory located at SGP4_small_V7.9_LINUX64/Lib/ into pysaademo/libdll/.
Environment Setup
By default, Linux expects to find DLLs in specific locations, like /usr/local/lib. When developing, it is necessary to set an environment variable called LD_LIBRARY_PATH to add a new search location for these DLLs. Unfortunately, it is impossible to do this from within the python script, so you will need to add that environment variable manually for every session.
Edit set_env and add the following single line:
export LD_LIBRARY_PATH='./libdll'
From now on, you can set the environment variable by simply running source set_env from the terminal. Note that this is a relative path, so it will only work when your current working directory contains the libdll/ directory.
Initialization Sequence for the DLLs
To interact with a DLL from the Python3 interpreter, you'll need to use the ctypes library. Each time you want to interface with an API, you'll need to use ctypes to tell the Python3 interpreter exactly how to communicate through the API.
You will also need to fulfill all the requirements demanded by the SAA DLLs.
Sgp4PropWrapper.chm
If you read the "SGP4 Architecture" page in Sgp4PropWrapper.chm, you will find the above chart. What it means is that you need to intialize the DLLs in a specific sequence for them to work. There are 6 DLLs with specific interactions that must be understood.
DLLMain must be initialized first. It will provide a "handle" that is needed for the other DLLs to intercommunicate.
Next you'll need TimeFunc, then you'll need TLE. TLE depends on TimeFunc to understand Epoch. TLE will help you manage SatKey objects, which are codes used by the DLLs to uniquely identify TLE sets.
After that, you should initialize EnvConst and AstroFunc. EnvConst contains general numerical constants and accessor methods for the other DLLs. AstroFunc depends on EnvConst to get a variety of constants like the gravitational parameter, \(\mu\). However, TimeFunc also has some unmentioned dependencies on EnvConst, so make sure to initialize EnvConst before doing any serious computations.
Starting the pysaademo.py script
The full pysaademo.py script is available on github.
Start the pysaademo.py script as follows:
#! /usr/bin/env python3
import os #used to set environment variable
import ctypes as c #module that interacts with DLLs
import pdb
# Python can't set this environment variable for ctypes,
# so tell the user if they need to set it before starting
# the script
if 'LD_LIBRARY_PATH' not in os.environ:
print("Set LD_LIBRARY_PATH before running")
exit()
Importing os and checking that LD_LIBRARY_PATH is a good way to catch basic user error while developing.
pdb is the 'python debugger'. If you get stuck, put pdb.set_trace() on the line before whichever line throws an error to trigger the interactive debugger. You can read more about pdb here.
ctypes is the python C interface library. It is fully documented here and is the primary tool for accessing the SAA DLLs.
Init MainDLL
The next section of pysaademo.py is all about initializing the DLLs.
# Init Main
maindll = c.CDLL('libdllmain.so')
maindll.DllMainInit.restype = c.c_int64
maindll_handle = maindll.DllMainInit()
This block initializes maindll.
maindll = c.CDLL('libdllmain.so') tells the Python ctypes library to go retrieve libdllmain.so. If this fails, it is likely because LD_LIBRARY_PATH isn't set correctly.
We want to run the API DllMainInit() method, but before doing that we need to tell Python exactly what kind of data DLLMainInit() accepts as arguments and what kind of data it returns. By default, ctypes will convert integers, bytes, and strings, but don't rely on this implicit behavior. You will get unpredictable results.
We define the return type by setting the .restype property to be a valid ctypes type object. In the case of DllMainInit, the returned object is a 64 bit integer which is a "handle," or a special kind of memory address used as a common meeting place for the other DLLs to find each other.
maindll_handle = maindll.DllMainInit()
If you reference DllMainWrapper.chm, the documentation states that the main init function returns an item of type "longlong". However, if you reference the documentation for the init functions of the other DLLs, they accept an item of type "int64". For 64 bit systems, these data types should be equivalent. On long time scales, it is a good idea to use the less ambiguous "int64" data type. Who knows, maybe you'll prevent a bug on a 128 bit computer in the future!
Start a log file
Main DLL has a method called OpenLogFile() which will tell all of the DLLs to write error messages to a common location. This project would have been impossible if I hadn't used this feature.
The documentation states that the returned value from this function is an integer value of 0 on success, and various small numbers on failure.
We also need to set the argument for this. C and Fortran have no native concept of "strings". They have a data type for "character" called "char", and to build a string the programmer creates an array that ends with a null character. This is the so-called "null-terminated character array."
To tell this function to save the data to "./run.log", we'll need to make an array that looks like this:
['.','/','r','u','n','.','l','o','g','\0']
To add to our difficulties, C and Fortran are incapable of directly passing arrays. To receive the above array, we need to pass a pointer to that array. A pointer is a memory address. Directly manipulating pointers is very powerful and very error-prone.
Fortunately, the ctypes library will do this automatically if we feed a bytes formatted string into the ctypes.c_char_p() constructor. This is handled implicitly by ctypes because we specified the .argtypes property for this method.
maindll.OpenLogFile.restype = c.c_int
maindll.OpenLogFile.argtypes = [c.c_char_p]
retcode = maindll.OpenLogFile(b"./run.log")
In this section, I set the return type by defining the OpenLogFile.restype to be equal to c.c_int.
I then set the argument types. This will always be represented with a list of type constructors. In our case, there was only one argument, a pointer to a null-terminated character array, so our argument type is c.c_char_p.
Third, I ran the OpenLogFile function with an appropriately formatted null-terminated string. My python string, b"./run.log", is automatically fed into the c.c_char_p() constructor because of how I specified this method's argtypes.
If all goes well, error messages will now be written to the ./run.log file.
Init the TimeFunc, TLE, EnvConst, and Astro DLLs
After initing the main dll, these are fairly straight forward. In each case, we specify the name of the "shared object" (so) or "dynamic link library" (DLL) file. Then we specify the return type and the argument type. In each case, the documentation specified an int64 input and an integer return. Finally, we initialize the DLLs using their various init functions.
# Init TimeFunc
timedll = c.CDLL('libtimefunc.so')
timedll.TimeFuncInit.restype = c.c_int
timedll.TimeFuncInit.argtypes = [c.c_int64]
retcode = timedll.TimeFuncInit(maindll_handle)
# Init TLE
tledll = c.CDLL('libtle.so')
tledll.TleInit.restype = c.c_int
tledll.TleInit.argtypes = [c.c_int64]
retcode = tledll.TleInit(maindll_handle)
# Init EnvConst
envdll = c.CDLL('libenvconst.so')
envdll.EnvInit.restype = c.c_int
envdll.EnvInit.argtypes = [c.c_int64]
retcode = envdll.EnvInit(maindll_handle)
# Init Astro
astrodll = c.CDLL('libastrofunc.so')
astrodll.AstroFuncInit.restype = c.c_int
astrodll.AstroFuncInit.argtypes = [c.c_int64]
retcode = astrodll.AstroFuncInit(maindll_handle)
Keep in mind that this isn't a very robust way to initialize DLLs. It would be wise to test that retcode is zero, then pull the last known error using the maindll function GetLastErrMsg() if retcode isn't 0. Since this is a quick and dirty tutorial, I won't dwell on that.
Init SGP4
The SGP4 propagation DLL is a bit trickier to initialize. It requires that a license file be correctly detected prior to initialization. This behavior is documented in the Sgp4SetLicFilePath() and Sgp4Init() methods in Sgp4PropWrapper.chm.
# Init SGP4
sgp4dll = c.CDLL('libsgp4prop.so')
# Gotta get that license file
sgp4dll.Sgp4SetLicFilePath.argtypes = [c.c_char_p]
sgp4dll.Sgp4SetLicFilePath(c.c_char_p(b"./libdll/"))
# Init the DLL
sgp4dll.Sgp4Init.restype = c.c_int
sgp4dll.Sgp4Init.argtypes = [c.c_int64]
retcode = sgp4dll.Sgp4Init(maindll_handle)
Pass By Reference
C and Fortran are quite old languages. Most modern languages exist in a function based paradigm. The idea is as follows: Given a function \(f(x) = y\), your programming functions should accept an input \(x\) and return an output \(y\).
In C and Fortran, this behavior is the exception, not the norm. Instead, it was quite common to provide the address of where the function should write its results. This is more of a subroutine based paradigm.
I've written the following short C program to highlight these differences.
#include<stdio.h>
// headers
int add_fun(int num1, int num2);
void add_ref(int * num1, int * num2, int * ans);
int main() {
int num1 = 1;
int num2 = 2;
int ans = 0;
//perform functional addition
ans = add_fun(num1, num2);
printf("add_fun: %i\n", ans); //prints '3'
//perform reference addition
ans = 0;
add_ref(&num1, &num2, &ans);
printf("add_ref: %i\n", ans); //prints '3'
//safe exit
return (0);
}
// functional programming
int add_fun(int num1, int num2)
{
int num3;
num3 = num1 + num2;
return (num3);
}
// pass by reference
void add_ref(int * num1, int * num2, int * ans)
{
*ans = *num1 + *num2;
}
The add_fun function behaves like a modern programming function. It has two inputs, it computes and returns an output.
The add_ref structure is very different. The memory addresses of the three numbers are passed to the function, then they are directly modified using the de-reference operator, *.
The next several SAA DLL methods we are using will employ this memory based subroutine approach. The return values will be integers that equal 0 if there are no errors.
Retrieving the license path
I'm going to demonstrate how to use ctypes to pass by reference from Python3.
The SGP4 DLL documentation states that Sgp4GetLicFilePath() is used to retrieve the current path to the SGP4 license file. This is a void function, which means nothing is returned. There is only one argument, and it is specified to be a 512 byte wide character array pointer. This is a bit different than a null-terminated string, because the length of 512 is guaranteed.
Sgp4GetLicFilePath() will modify the 512 byte character array that is passed to it.The modified value will be the current search location for the license file.
sgp4dll.Sgp4GetLicFilePath.argtypes = [c.c_char_p]
byte512 = c.c_char*512
path = byte512()
sgp4dll.Sgp4GetLicFilePath(path)
print('-'*20)
print('SAA Licence Path:')
print(path.value)
print('-'*20)
The first line sets the argument type, which is a pointer to a character array.
The second line initiates a 512 character long ctypes array. The ctypes documentation recommends this method of creating array data types. I named my 512 byte wide array type object byte512.
The third line initializes the variable path that we will be passing by reference. Because this is an array, we don't have to do anything special to get the memory address: An array is its memory address. In fact, if you add one to the value, you get the address of the second entry, add two and you get the address of the third entry, and so on.
The fourth line calls the Sgp4GetLicFilePath() method with the path argument. Remember that we specified that path would be a c.c_char_p type object.
The next several lines print the result. Because path is a ctypes object, it has a .value property which is a python compatible string. Printing it will reveal that it is a 512 character long string starting with ./libdll/. This matches the spec.
--------------------
SAA Licence Path:
b'./libdll/ '
--------------------
Reading a TLE using TleAddSatFrLines()
Now that we know how to pass by reference, we should have all the knowledge we need to use this API.
The next step is to read a TLE into the TLE DLL. The objective is to register this TLE in the system, and retrieve the SatKey object which is used to uniquely identify this TLE element.
tledll.TleAddSatFrLines.restype = c.c_int64
tledll.TleAddSatFrLines.argtypes = [c.c_char_p, c.c_char_p]
line1 = c.c_char_p(b"1 90001U SGP4-VAL 93 51.47568104 .00000184 0 0 00000-4 814")
line2 = c.c_char_p(b"2 90001 0.0221 182.4922 0000720 45.6036 131.8822 1.00271328 1199")
SatKey = tledll.TleAddSatFrLines(line1, line2)
The TLE library has the TleAddSatFrLines() method, which allows me to directly input two strings of data, one for each line of the TLE. After doing so, it will return a SatKey int64 object, which is similar to the MainDLL handle in its purpose. This SatKey number will be used by the SAA DLLs to uniquely identify this instance of a TLE.
Initializing a Satellite Object in the SGP4PROP DLL
The DLLs use an object oriented methodology that requires the satellite be initialized from the TLEs we just read. This should look familiar by now.
# init an Sgp4InitSat object
sgp4dll.Sgp4InitSat.restype = c.c_int
sgp4dll.Sgp4InitSat.argtypes = [c.c_int64]
retcode = sgp4dll.Sgp4InitSat(SatKey)
Propagating The TLE
I've chosen to propagate the TLEs by using the number of minutes past Epoch. This is captured in the Sgp4PropMse() method, which uses a particularly complicated form of the "pass by reference" model explained before. You can think of this as the "boss fight" of this tutorial.
First, go read the documentation on this method carefully.
The returned value is a simple integer, 0 on success.
The input values are, in order:
- SatKey
- Minutes past Epoch (
mse) - A reference to a double which will contain our ds50UTC value, which is decimal days since 1950.
- An array composed of three doubles representing the three components of position in the ECI reference frame
- An array composed of three doubles representing the three components of velocity in the ECI reference frame
- An array composed of three doubles representing the geodedic latitude, longitude, and height
Items 1 and 2 are the only true inputs. The remaining items are being passed by reference to be modified.
vector = c.c_double * 3
r_ECI = vector()
v_ECI = vector()
llh = vector()
ds50UTC = c.c_double()
sgp4dll.Sgp4PropMse.restype = c.c_int
sgp4dll.Sgp4PropMse.argtypes = [c.c_int64, c.c_double, c.POINTER(c.c_double), vector, vector, vector]
The first line defines the "vector" data type to be a double array with three elements.
The second, third, and fourth line define a position, velocity, and LLH vector which can be passed by reference.
The fifth line defines a ds50UTC double which can be passed by reference.
The seventh and eighth lines define the return value and inputs to Sgp4PropMse(). There is something new here: The ctypes library doesn't have a "pointer to a double" data type, but it does have a generic pointer constructor called POINTER(). In order to pass a double object by reference, I used the POINTER() method to construct a ctypes pointer object.
After we compute r_ECI, v_ECI, and llh, we're going to want to print the resulting vectors. This function will print those vectors to standard out. This is plain 'ol python, so no tricks here.
def printvector(name, vector):
print("{:5} = < {:13.7f}, {:14.7f}, {:14.7f} >".format(name, vector[0], vector[1], vector[2]))
Now, let us propagate 0 minutes. This effectively converts our TLE into position and velocity, but does no propagation. This is a useful way to do conversions.
# Do one run at 0 minutes past Epoch
mse = c.c_double(0)
retcode = sgp4dll.Sgp4PropMse(SatKey, mse, c.byref(ds50UTC), r_ECI, v_ECI, llh)
print('Sgp4PropMse Return Code: {:d}'.format(retcode))
print("ds50UTC: {0:.7f}".format(ds50UTC.value))
printvector('r_ECI', r_ECI)
printvector('v_ECI', v_ECI)
printvector('llh', llh)
print('-'*20)
Line 2 sets our minutes past Epoch to zero. Line 3 performs the propagation and passes by reference. The output is printed to the shell.
I want to take a moment to talk about the byref function. The C language has a concept called a "reference," which is like a "pointer" but it has more rules attached to it. Some people insist that a reference isn't a pointer, and some insist that it is. In practice, references can often be used like pointers, as is the case here. If you decide to learn some C, it is worthwhile to really dig down into the differences.
The next propagation we want to do is 2700 minutes. That is identical except that mse will of course be set to 2700.
# Do one run at 2700 minutes past Epoch
mse = c.c_double(2700)
retcode = sgp4dll.Sgp4PropMse(SatKey, mse, c.byref(ds50UTC), r_ECI, v_ECI, llh)
print('Sgp4PropMse Return Code: {:d}'.format(retcode))
print("ds50UTC: {0:.7f}".format(ds50UTC.value))
printvector('r_ECI', r_ECI)
printvector('v_ECI', v_ECI)
printvector('llh', llh)
print('-'*20)
Results will be printed to the shell.
Validation
Now that we've done all this work, we want to know our results are correct. To prepare for this, I used the same TLE and time steps used in the validation set included with the SAA DLL package.
If you did everything right, the following should be printed to the terminal:
--------------------
SAA Licence Path:
b'./libdll/ '
--------------------
Sgp4PropMse Return Code: 0
ds50UTC: 15757.4756810
r_ECI = < 42166.3724464, 12.9531593, -9.8621994 >
v_ECI = < -0.0007811, 3.0745638, -0.0017028 >
llh = < -0.0134144, 38.3680760, 35788.2405904 >
--------------------
Sgp4PropMse Return Code: 0
ds50UTC: 15759.3506810
r_ECI = < 30767.5128077, -28829.9475673, 11.4616939 >
v_ECI = < 2.1025258, 2.2435062, -0.0017973 >
llh = < 0.0155908, 38.3644662, 35785.9000608 >
--------------------
These values can be validated by locating them in the baseline data set located in the SAA package at these locations:
SGP4_small_V7.9_WIN64/Verify/Baseline/sgp4_val_C_LatLonHeight.txt
TSINCE (MIN) LAT(DEG) LON (DEG) HT (KM) X (KM) Y (KM) Z (KM)
1 90001U SGP4-VAL 93051.47568104 +.00000184 00000 0 00000-4 0 0814
2 90001 0.0221 182.4922 0000720 45.6036 131.8822 1.0027132801199
0.0000000 -0.0134144 38.3680760 35788.2405904 42166.3724464 12.9531591 -9.8621994
2700.0000000 0.0155908 38.3644662 35785.9000608 30767.5128078 -28829.9475673 11.4616939
SGP4_small_V7.9_WIN64/Verify/Baseline/sgp4_val_C_OscState.txt
TSINCE (MIN) X (KM) Y (KM) Z (KM) XDOT (KM/S) YDOT(KM/S) ZDOT (KM/SEC)
1 90001U SGP4-VAL 93051.47568104 +.00000184 00000 0 00000-4 0 0814
2 90001 0.0221 182.4922 0000720 45.6036 131.8822 1.0027132801199
0.0000000 42166.3724464 12.9531591 -9.8621994 -0.0007811 3.0745638 -0.0017028
2700.0000000 30767.5128078 -28829.9475673 11.4616939 2.1025258 2.2435062 -0.0017973
As you can see, our values match the values provided in the baseline set exactly.
Conclusion
By using the SAA DLLs, we can match AFSPC's satellite propagation solution exactly. If there is a major bugfix or update, we can update our software quickly and easily by copying over the relevant DLLs.
My next step in this process is to build a robust Python3 module and license it for open source use. Given Python3's popularity, I'm frankly surprised I haven't found one already exists.