# PRÁCTICA GPS (2017)

Pràctica InglésUniversidad | Universidad Politécnica de Cataluña (UPC) |

Grado | Ingeniería de Aeronavegación - 3º curso |

Asignatura | Navegación, cartografía y cosmografía |

Año del apunte | 2017 |

Páginas | 11 |

Fecha de subida | 25/06/2017 |

Descargas | 3 |

Subido por | areig |

### Vista previa del texto

GROUND TRACK OF
GPS SATELLITES
PRACTICE 1
Air Navigation, Cartography and Cosmology
Alba Martín
Anna Reig
6GX31
17.04.2017
Table of contents
1.

Introduction ......................................................................................................................3
1.1
Aim of the practice .................................................................................................3
1.2
Almanac explanation .............................................................................................3
2.

ECEF coordinates .............................................................................................................5
3.

LLA coordinates ...............................................................................................................7
4.

Complete satellite ground track .................................................................................8
5.

All satellites ground track in the next hour ............................................................9
6.

Conclusions ....................................................................................................................11
1. Introduction
1.1 Aim of the practice
The main objective of this project is to plot the satellite ground track for any given
satellite or constellation of it. To do so, real-time ephemerides are available on
Internet for known Keplerian elements, what makes it easy to work with GPS
constellation.

1.2 Almanac explanation
First it is crucial to understand ephemerides and Keplerian elements. Satellite
ephemerides are on GPS almanac, which includes information of the entire GPS
satellite constellation for every orbit. The files used can be found on
http://celestrak.com/GPS/almanac/Yuma/2017/.

Every almanac contains the following data and Keplerian elements about each
satellite that must be defined and correctly understand:
1.

ID: PRN (Pseudo Random Noise) ID of the SVN (Space Vehicle NAVSTAR)
2.

Health: Indicates the operational state of the SV (Space Vehicle). The value
“000” means that the SV is usable.

3.

Eccentricity: It shows the amount of the orbit deviation from a circular orbit. It
is also the distance between the foci divided by the length of the semi-major
axis (the NAVSTAR orbits are almost circular).

4.

Time of Applicability (ToA): The time instant (seconds within the GPS week)
for which the almanac has been computed (ephemeris epoch).

5.

Orbital Inclination: The angle of the SV orbit plane with respect to the equator
(GPS is at approx. 55 degrees) in radians. The SV ground track will not rise
above approximately 55 degrees of latitude.

6.

Rate of Right Ascension: Rate of change of the angle of right ascension in
rad/s.

7.

SQRT(A) Square Root of Semi-Major Axis: The semi-major axis is the distance
from the centre of the ellipse to either the point of apogee or the point of
perigee.

8.

Right Ascension at GPS week epoch: Geographic longitude (in radians, with
respect to Greenwich meridian) of the ascending node of the orbit at the GPS
week epoch (Saturday to Sunday midnight).

9.

Argument of Perigee at ToA: An angular measurement (in radians) along the
orbital path measured from the ascending node to the point of perigee. It is
measured in the direction of the SV's motion.

10.

Mean Anomaly at ToA: Angle (in radians) travelled past the perigee at ToA.

When the SV has passed perigee and heading towards apogee, the mean
anomaly is positive. After the point of apogee, the mean anomaly value will be
negative to the point of perigee.

11.

Af(0): SV clock bias in seconds
12.

Af(1): SV clock Drift in seconds per seconds
13.

Week: GPS week number (0000-1023). Counter increased (module 1024)
every week since 0h UTC of Jan 6th, 1980 (the GPS epoch = 2444244.5 JD). First
roll-out happened at 23:59:47 (UTC) on Aug 21, 1999 (roll epoch =
2451412.499850 JD).

Once defined all the parameters, the almanac.m function, already created, will
download current GPS almanac and create a matrix within the ephemerides.

2. ECEF coordinates
ECEF (Earth-Centred, Earth-Fixed) coordinates
is a geographic and Cartesian coordinate
system that represents positions as X, Y and Z
coordinates. Its axes are aligned with the
International Reference Pole and International
Reference Meridian that are fixed with respect
to the earth. In other words, it has its origin at
the Earth’s centre of mass with its x-axis
extend through the intersection of the prime
meridian (0º longitude) and the equator (0º
Figure 1 ECEF Coordinates system
latitude). The z-axis extends through the true
north pole (i.e., coincident with the Earth spin axis). The y-axis completes the righthanded coordinate system, passing through the equator and 90º longitude.

Consequently, the created function ECEF.m in Matlab will be able to compute ECEF
coordinates of any satellite at a given time, and will have the possibility of plotting
ground track in three dimensions in meters for a given period of time.

As exposed before, the orbits are characterized by the six Keplerian parameters.

Using these parameters, the ECEF coordinates will be computed for any satellite at
current time without taking into account the slow changes in the elements. To obtain
the mentioned coordinates, the following steps must be done:
1.

Computation of the time between the ToA of the ephemerides and the current
time:
𝑡𝑘 = 𝑡 − 𝑡0
2.

Getting the semi-major axis:
𝑎 = (√𝑎)
3.

2
Getting the mean motion (corresponding to the angular speed):
𝑛=√
𝐺·𝑀
𝑎3
Where G in the gravitational constant and M the earth mass.

4.

Computation of the mean anomaly:
𝑀𝑘 = 𝑀0 + 𝑛 · 𝑡𝑘
5.

Iteratively for Ek , find the eccentric anomaly:
𝑀𝑘 = 𝐸𝑘 − 𝑒 · sin(𝐸𝑘 )
𝐸𝑘 (𝑛) = 𝑀𝑘 + 𝑒 · sin(𝐸𝑘 (𝑛 − 1))
6.

Computation of the true anomaly without ambiguity (converting sinus and
cosine into complex numbers and using then the function angle):
sin(𝑣𝑘 ) = √1 − 𝑒 2 ·
cos(𝑣𝑘 ) =
7.

sin(𝐸𝑘 )
,
1 − 𝑒 · cos(𝐸𝑘 )
cos(𝐸𝑘 ) − 𝑒
1 − 𝑒 · cos(𝐸𝑘 )
Computation of the argument of latitude:
𝑢𝑘 = 𝑣𝑘 + 𝜔
8.

Find the orbit radius, that corresponds to the distance to Earth center):
𝑟𝑘 = 𝑎 · (1 − 𝑒 · cos(𝐸𝑘 ))
9.

Current longitude of the ascending node:
Ω𝑘 = Ω0 + Ω̇0 · 𝑡𝑘 − Ω̇𝑒 · 𝑡
Where Ω̇𝑒 = 7.2921151467 · 10−5 𝑟𝑎𝑑/𝑠
10.

Computation of the X coordinate within the orbital plane:
𝑋𝑝 = 𝑟𝑘 · cos(𝑢𝑘 )
11.

Computation of the Y coordinate within the orbital plane:
𝑌𝑝 = 𝑟𝑘 · sin(𝑢𝑘 )
Once done step eleven, it is time to compute the X, Y and Z coordinates for ECEF.

12.

ECEF X-coordinate:
𝑋 = 𝑋𝑝 · cos(Ω𝑘 ) − 𝑌𝑝 · cos(𝑖0 ) · sin(Ω𝑘 )
13.

ECEF Y-coordinate:
𝑌 = 𝑋𝑝 · sin(Ω𝑘 ) + 𝑌𝑝 · cos(𝑖0 ) · cos(Ω𝑘 )
14.

ECEF Z-coordinate:
𝑍 = 𝑌𝑝 · sin(𝑖0)
This type of system is useful to know
where we are in the Earth or the
coordinates of a point fixed on the
surface of the earth but if one wants to
plot the satellite ground track it is
necessary to convert them into LLA
coordinates as the ECEF coordinates
moves with the Earth unlike the system
of LLA coordinates.

Figure 2 ECEF Coordinates representation
3. LLA coordinates
LLA (Latitude, Longitude and Altitude)
coordinates are chosen such that one of
the
numbers
represents
a
vertical
position, and two of the numbers
represent horizontal position. Enables
every location on Earth to be specified
by a set of numbers, letters or symbols.

Consequently, the
created
function
LLA.m in Matlab will be able to compute
Figure 3 LLA Coordinates system
LLA coordinates of any satellite at a
given time, and will have the possibility of plotting ground track using degrees for
the latitude and the longitude and meters for the height for a given period of time.

As it is a two dimensions’ representation, the altitude will not be considered in the
plots.

So that, once having the conversion to ECEF Coordinates, they must be converted to
LLA following the algorithm shown below so to have a plot in degrees representing
the ground track of a certain satellite correctly.

1.

𝑟 = √𝑋 2 + 𝑌 2
2.

𝐸2 = 𝑎2 − 𝑏 2
3.

𝐹 = 54 · 𝑏 2 · 𝑍 2
4.

𝐺 = 𝑟 2 + (1 − 𝑒 2 ) · 𝑍 2 − 𝑒 2 · 𝐸 2
5.

𝑠 = √1 + 𝑐 + √𝑐 2 + 2 · 𝑐
6.

𝑃=
7.

𝑄 = √1 + 2 · 𝑒 4 · 𝑃
8.

𝑟0 = −
9.

𝑈 = √𝑍 2 + (𝑟 − 𝑟0 · 𝑒 2 )2
3
10.

11.

12.

13.

𝐹
2
1
𝑠
3·(𝑠+ +1) ·𝐺 2
𝑃·𝑒 2 ·𝑟
1+𝑄
1
1
𝑃·(1−𝑒 2 )·𝑍 2
2
𝑄
𝑄(1+𝑄)
+ √ · 𝑎2 · (1 + ) −
𝑏 2 ·𝑍
𝑌
𝜆 = 𝑎𝑟𝑐𝑡𝑔 ( )
𝑋
Where 𝑒 and 𝑎 are given by WGS-84
and have the value of:
𝑎·𝑉
ℎ = 𝑈 · (1 −
𝜙 = 𝑎𝑟𝑐𝑡𝑔 (
2
14.

𝑉 = √𝑍 2 · (1 − 𝑒 2 ) + (𝑟 − 𝑟0 · 𝑒 2 )2
𝑍0 =
1
− · 𝑃 · 𝑟2
𝑏2
𝑎·𝑉
𝑎 = 6378137𝑚
)
𝑍+𝑒 ′2 ·𝑍0
𝑟
𝑒 2 = 0.00669437999014
)
𝑐
𝑐
𝑎
𝑏
And also: 𝑎2 = 𝑏 2 + 𝑐 2 , 𝑒 = → 𝑒 ′ =
4. Complete satellite ground track
The complete satellite ground track has been obtained by creating a loop in the
Matlab program called Main_One_Satellite.m. This loop iterates for each 300
seconds until it reaches 24 hours and, at each iteration, it computes the ECEF
coordinates for the chosen satellite and then with these obtained data, the LLA
coordinates. As the period of the satellite orbit is ½ of a sideral day (12h), it takes
two consecutive orbits (periods) for the satellite to travel a complete ground track.

After selecting the satellite 12 of the 31 available and by means of LLA coordinates
the latitude and longitude can be plotted with in a world map in the background,
has shown below.

Figure 4 Ground track of the satellite number 12
5. All satellites ground track in the next hour
To know the current subsatellite point of all the GPS satellites and their
corresponding ground track for the next hour, it is needed to develop a new Matlab
program called Main_All_Satellits.m. This program has a loop that computes the
ECEF coordinates and then the LLA coordinates for each time interval until it reaches
1 hour for a certain satellite as the Main_One_Satellite.m but also has another loop
which englobes the previous one and computes the next 60 minutes for the 31
available satellites.

Figure 5 GPS satellite ground tracks in the next hour
Figure 6 Online satellites ground track for past hour
The resulting plot must be validated to assure their reliability. In the Web-page
http://www.nstb.tc.faa.gov/rt_waassatellitestatus.htm it is possible to download the
updated map on real-time but taking into account that this plot (Figure 6)
corresponds to the past hour, while the one extracted from the developed program
Main_All_Satellits.m is computing the ground track for the next hour. Looking to
both figures 6 and 7, it can be appreciated the match between the profile followed
by each satellite in the past hour and the profile that each satellite will perform in
the next hour.

6. Conclusions
Computing the GPS satellites position is a hard task that requires the use of many
equations involving the Keplerian elements. By means of the provided information
of the ephemerides, it is needed to change the coordinates system to ECEF
coordinates and then to LLA coordinates to obtain a better and more understanding
representation of the satellites positions and ground tracks.

Once developed the functions changing the coordinates system and the almanac, it
is just needed some looping iterations to represent the ground tracks of the desired
satellites.

If one analyses Figure 4, the satellite is covering different areas of the earth, both
the northern hemisphere and the south hemisphere. Even though, this is not enough
to cover the whole surface, that’s why it is needed at least four satellites to give the
position of a device in a precise way. As it has been seen, all over the world, there
are 31 satellites that are always covering the GPS service.

In figure 5, it is possible to verify how, thanks to these 31 satellites which compose
the GPS system, the covered Earth area is much greater. Although it is only possible
to see their positions in the next hour, it can be seen how the orbits represented on
the map of the earth would weave a kind of net, covering almost the whole terrestrial
surface. This fact explains the efficiency of the GPS system as well as the fact that
sometimes, when someone tries to get GPS positioning, the device takes a while as
it tries to get information from 4 of the 31 satellites.

...