I try to plot satellite ground track. I found lots of formulas but hardly find algorithm which describe how to uses them one-by-one and step-by-step. As I understand the algorithm will be:

1) Get orbital elements:

a -  semimajor axis
i - inclination
Ω - right ascension of the ascending node
ω  - argument of perigee
v - true anomaly
Note: I need to plot ground track knowing only this elements and I don’t get them from TLE.

2) Get time variables:

tCur  – current time
tBegin – begin time of propagation period
tEnd – end time of propagation period
deltaT – time step

3) Calculate supporting variables like:

Ra = Ha + Re – radius of semimajor axis
Rp = Hp + Re – radius of semiminor axis
e = (Ra-Rp)/(Ra+Rp) – eccentricity
T = 2PI*sqrt(a^3/MU) – period
n = 2PI/T – mean motion
Ha – satellite apogee height 
Hp – satellite perigee height
MU – Earth gravity constant = 398602

4) Calculate deltas of Ω and ω:

dΩ = -(1.5*sqrt(MU)*J2*(R3*R3)/((1-e^2)* a^3.5))*cos(i)
dω = dΩ *(2.5*(sin(i)*sin(i)) - 2)/cos(i)
J2 – constant harmonic coefficient  = 0.0010826257

5) Before loop calculate init values of:

E = 2*atan(tan(v/2)*sqrt((1-e)/(1+e))) – eccentricity anomaly
M = E - e*sin(E) – mean anomaly
    I find another formulas for E, but I don’t know which is better for accuracy.

6) Begin loop and continue while tCur < tEnd

7) Calculate:

E = KeplerEquatation(M,e,1E-09) – calculate eccentricity anomaly by solving Kepler equitation
M = E-e*sin(E) – change step in mean anomaly, also I find formula:
M = n*(tCur-t0) – Where: t0 – time when satellite cross perigee point
v = atan2(sqrt(1-e^2)*sin(E), cos(E)-e) – calculate true anomaly
u = v + ω – latitude argument? Show angle position respective to center of Earth, but I don’t actually know what is it.
r = a*(1-e*cos(E))

8) Calculate state vector movement coordinates XYZeci (ECI – Earth Centered Inertial) – position of satellite in 3D ?:

X = (cos(Ω)*cos(u)-sin(Ω)*sin(u)*cos(i)) * r
Y = (sin(Ω)*cos(u)+cos(Ω)*sin(u)*cos(i)) * r
Z = (sin(u)*sin(i)) * r

9) Rotate XYZ in ECI to ECEF(Earth-Centered Eartg-Fixed) with rotation

XYZecef = Rot(GST)* XYZeci
        cos(GST) sin(GST)  0
Rot =  -sin(GST) cos(GST)  0
       0         0         1
GST – greenwich sidereal time (in radians)

10) Translate ECEF to LLA(longitude-latitude-altitude) to find coordinates of satellite for plot ground track on 2D map.

11) Change time step:

tCur = tCur + deltaT

12) Repeat from 6.

Is this algorithm correct for calculating satellite ground track and orbit position?

  • 1
    $\begingroup$ You are asking for someone to go line-by-line through all of this looking for errors, and that's not a good question to ask in Stack Exchange. You should try it yourself! If you don't know how, then THAT would be an excellent Stack Exchange question because it requires a thoughtful answer. Of course, you could get lucky and someone will do exactly what your asking for, so I won't up/down vote or flag. Let's see what happens! $\endgroup$ – uhoh Feb 6 '18 at 16:42

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.