Anyway, in these first months of the project my focus has been on just trying to learn the basics of orbital mechanics and it is going rather well. For some time now I have been working on simulating what a ground station would measure when a satellite passes above it, and I wanted to talk a bit about it.
In order to simulate the passing of a satellite we need to go through a few steps. Firstly, we need to figure out the initial state of our satellite and then propagate this forward in time. Secondly, we need to find out when the satellite passes our station. And lastly, we have to determine what the pass will look like from the station.
Propagating the orbit
The propagation of an orbit is most easily done by using an existing orbit propagator, employing models such as the widely used sgp4 . However, in order to better understand the concepts behind this process I began writing my own propagator instead. This does have some additional disadvantages in terms of accuracy. The sgp4 model is a perturbations model and takes into account effects such as solar radiation, atmospheric drag, Earth's non-uniform gravitaional field etc. Adding these effects to my own propagator would be a daunting task, so instead I was limited to an unperturbed propagator, which assumes that the only force acting on the satellite is that of Earth's gravity acting uniformly from the centre of the Earth. This will obviously add some amount of error in the results, but we will be disregarding these errors for now.
To start doing anything at all, we need some initial information about the satellite. Luckily, NORAD tracks all (or at least most) satellites in orbit, and for each of these satellites (except the secret ones) they regularly release something called a two-line element set (TLE), which looks like this:
1 32789U 08021G 15341.86593944 +.00006209 +00000-0 +48321-3 0 9990
2 32789 097.6272 039.4006 0012136 181.1494 178.9704 15.01549889413023
The TLE contains among other things the satellite's catalogue number, and the date and time of the TLE, i.e. the epoch. But more interesting it also contains six independent orbital parameters called Keplerian elements. Two of these elements (semi-major axis, which has to be derived, and eccentricity) describe the size and shape of the elliptical orbit, one (mean anomaly) defines the satellites location in this orbit, and the remaining three (right-ascension of ascending node, inclination, and argument of periapsis) describe how the orbit is oriented in space. With some basic theory of astrodynamics and the two-body problem, these six parameters can then be used to derive a state vector, which describes the location and velocity of the satellite in Cartesian coordinates.
$$\mathbf{S}=[x,y,z,V_x,V_y,V_z]$$
Now that we have a vector describing all the initial conditions of the satellite in a useful way, we can start propagating the orbit from the epoch and forward in time. For this I decided to use the Runge-Kutta iterative method known as RK4.
\begin{align*}\mathbf{S}_{n+1} = \mathbf{S}_n + \frac{h}{6}(k_1 + 2k_2 + 2k_3 + k4)\end{align*}
\begin{align*}k_1 =& f(\mathbf{S}_n)\\ k_2 =& f(\mathbf{S}_n+\frac{h}{2}k_1) \\ k_3 =& f(\mathbf{S}_n+\frac{h}{2}k_2) \\ k_4 =& f(\mathbf{S}_n+k_3)\end{align*}
\begin{align*}k_1 =& f(\mathbf{S}_n)\\ k_2 =& f(\mathbf{S}_n+\frac{h}{2}k_1) \\ k_3 =& f(\mathbf{S}_n+\frac{h}{2}k_2) \\ k_4 =& f(\mathbf{S}_n+k_3)\end{align*}
where $h$ is the time step, and the function,
\begin{align*} f(\mathbf{S})=\mathbf{S}'=[V_x,V_y,V_z,a_x,a_y,a_z]\end{align*}
is the change in the state vector, including the acceleration due to Earth's gravity.
\begin{align*} f(\mathbf{S})=\mathbf{S}'=[V_x,V_y,V_z,a_x,a_y,a_z]\end{align*}
is the change in the state vector, including the acceleration due to Earth's gravity.
Propagating the orbit forward in time with a time step of 1-20 seconds, we can create a large dataset with the state vectors at each time step. This dataset then contains all the necessary information about the orbit
We also need to implement a ground station from which we can observe the satellite. The ground station is defined by its latitude, longitude, and height above sea level, which we transform into a Cartesian state vector similar to the satellite's. We now encounter a problem regarding our reference system. The satellite is located in an Earth-centred inertial (ECI) reference system which does not rotate with the earth. And the station is located in an Earth-centred, Earth-fixed (ECEF) reference system which does rotate with the Earth. By transforming the state vectors of the satellite from ECI to ECEF, we ensure all of our data is in a single reference system.
We also need to implement a ground station from which we can observe the satellite. The ground station is defined by its latitude, longitude, and height above sea level, which we transform into a Cartesian state vector similar to the satellite's. We now encounter a problem regarding our reference system. The satellite is located in an Earth-centred inertial (ECI) reference system which does not rotate with the earth. And the station is located in an Earth-centred, Earth-fixed (ECEF) reference system which does rotate with the Earth. By transforming the state vectors of the satellite from ECI to ECEF, we ensure all of our data is in a single reference system.
Finding all the passes
With our data set of state vectors we can start doing some interesting things. For example, we can easily plot the orbit, since our state vectors already contains the Cartesian coordinates of the satellite at every individual time step.
However, we are more interested in simulating passes of the satellite over our ground station. For this we have to compare the respective state vectors of the station with those of the satellite. With some basic vector operations we establish the times when the satellite is visible to the station, i.e. above the horizon. This can be seen in the following graph where a value of 1 means the satellite is visible to the station.
In this specific case there seems to be a certain periodic nature with series of three and four passes and with an overall period of approximately one day. We can also find the elevation of the satellite with respect to the station for each state vector, and from this the maximum elevation of each pass. The following output is then produced.
However, we are more interested in simulating passes of the satellite over our ground station. For this we have to compare the respective state vectors of the station with those of the satellite. With some basic vector operations we establish the times when the satellite is visible to the station, i.e. above the horizon. This can be seen in the following graph where a value of 1 means the satellite is visible to the station.
In this specific case there seems to be a certain periodic nature with series of three and four passes and with an overall period of approximately one day. We can also find the elevation of the satellite with respect to the station for each state vector, and from this the maximum elevation of each pass. The following output is then produced.
| Start | Max elevation | End
Date | Time | Time El | Time
----------------------------------------------------
2015-06-13 | 15:49 | 15:54 12 | 15:59
2015-06-13 | 17:24 | 17:30 86 | 17:36
2015-06-13 | 18:59 | 19:04 13 | 19:09
2015-06-14 | 02:42 | 02:47 17 | 02:53
2015-06-14 | 04:16 | 04:22 63 | 04:28
2015-06-14 | 15:49 | 15:54 13 | 15:59
2015-06-14 | 17:24 | 17:30 78 | 17:36
2015-06-14 | 19:00 | 19:04 12 | 19:09
2015-06-15 | 02:42 | 02:47 18 | 02:53
2015-06-15 | 04:16 | 04:22 57 | 04:28
Date | Time | Time El | Time
----------------------------------------------------
2015-06-13 | 15:49 | 15:54 12 | 15:59
2015-06-13 | 17:24 | 17:30 86 | 17:36
2015-06-13 | 18:59 | 19:04 13 | 19:09
2015-06-14 | 02:42 | 02:47 17 | 02:53
2015-06-14 | 04:16 | 04:22 63 | 04:28
2015-06-14 | 15:49 | 15:54 13 | 15:59
2015-06-14 | 17:24 | 17:30 78 | 17:36
2015-06-14 | 19:00 | 19:04 12 | 19:09
2015-06-15 | 02:42 | 02:47 18 | 02:53
2015-06-15 | 04:16 | 04:22 57 | 04:28
This summary presents the start time and end time of each pass, as well as the maximum elevation and time hereof. All passes with a maximum elevation of less than 10 degrees have been excluded, which is why a few passes from the graph are not present in the summary. The data set of state vectors is split into each individual pass for further use.
What does a pass look like?
Now that we finally have state vector data for each of the passes as well as the ground station, we can simulate how the pass would look when viewed from the station. More specifically we are interested in elevation and range-rate, with range-rate being the velocity component of the satellite in the direction of the ground station. The range rate is the only parameter that can be measured using the Doppler-shift tracking that Doptrack employs. This is done by measuring the Doppler-shifted frequency of the satellites radio signal.
The following figures show the range-rate and elevation of the first two passes from the above summary.
These two passes have widely different elevations of 12 and 86 degrees and we can see some of the effects of this on the range-rate. A pass with higher maximum elevation has a greater change in range rate at the middle of the pass. The maximum and minimum range-rate, as well as the duration of the pass, are also greater.
With the range-rate data we calculate what our ground station would measure in terms of frequency and the associated frequency shift. The following relation describes the measured frequency.
$$ f = \left( 1 - \frac{\dot{\mathbf{\rho}}}{c} \right) f_0 $$
where $\dot{\mathbf{\rho}}$ is the range-rate, $c$ is the speed of light, and $f_0$ is the frequency transmitted by the satellite. The following graphs show the frequency shift of the Delfi-C3, with a transmitted frequency of 145.87 MHz, and a satellite with a transmitted frequency half of that.
The graph of the frequency shift also shows what we should expect when we actually measure the frequency of the satellite at our ground station. The following figures show the plot of such measurements and they clearly contain curves that resemble our simulations. The plots also contain side lobes of the frequency curve, but we are only interested in the main lobe at the centre. With our new knowledge we can quickly conclude that the first plot shows the data of a pass with a higher elevation than the second plot due to the shape of the frequency curve.
This then concludes the initial simulation of satellite passes. The coming weeks I will continue my work by looking into the atmospheric effects on the frequency in hope of even more accurate results.
