4
$\begingroup$

When the ISS orbit the Earth, its path traced on a 2D Mercator Earth Map is similar to a sinusoidal wave. I am wondering is there any way to determine the trig function of the wave, knowing the orbital height (which decides its speed) and direction in relation to Earth's rotational axis. This will possibility allow you to determine the position of a satellite at an specific time on a 2D map similarly to here: https://observablehq.com/@jake-low/satellite-ground-track-visualizer

Edit: Julius gave an awesome answer using planar geometric laws. How does it affect spherical triangles?

$\endgroup$
1
  • $\begingroup$ It's certainly possible, VSOP87 is a popular algorithm for computing positions of the planets, which kinda does it that way. But it's not very compact, Kepler's equation is pretty compact, but not that accurate. Especially for satellites near the Earth which experience gravitational effects from the moon, sun, and Earth non-symetric figure, and atmospheric drag. There really is no substitute for numerical integration for accuracy. But a Fourier Series could give you at least something that looks like a propper ground path over a short period. $\endgroup$ Commented Jun 29, 2022 at 0:09

2 Answers 2

6
$\begingroup$

If your goal is to plot on a 2D map the position of the satellite, what you can do is draw its pairs longitude-latitude $(Lo,La)$ as a function of time. To find these values you first calculate the pairs $(\alpha,\delta)$ named right ascension and declination of the satellite. You'll need the satellite orbital elements and eventually its position in an Earth centered reference system.

You can start by looking at the following diagram: I made it on paper cause I learned this in university but I can't find anything similar on astrodynamics textbooks. We're on a geocentric-equatorial reference frame with $I,J,K$ axis.

The satellite is at time $t$ on point $L$, you can see a portion of its orbit in red starting from its ascending node on point $N$, on the equatorial plane. From point $N$ to point $L$ the satellite swept an angle $u(t)$ (argument of latitude at time $t$). If you look at the spherical triangle $NML_F$ you can write the sine law: $$\frac{\sin\delta}{\sin i}=\frac{\sin u(t)}{\sin \pi/2}$$

Where $i$ is the orbit's inclination. The declination is equal to the latitude $La=\delta$, so you find the first coordinate to plot on the 2D map by solving the non linear equation $$sin\delta =\sin i\sin u(t)$$

The longitude is found in a similar way: from the spherical triangle $NLH_F$ you write the following sine law $$\frac{\sin\alpha'}{\sin u(t)}=\frac{\sin(\pi/2-1)}{\sin(\pi/2-\delta)}$$

where $\alpha'$ is the angle between the ascending node and the projection of the current satellite position on the equatorial plane (point $M$). So you find this angle by solving the following equation: $$sin\alpha'=\sin u(t)\frac{\cos i}{\cos\delta}$$ Then to get $\alpha$ you'll sum to this result right longitude of ascending node angle $$\alpha=\Omega+\alpha'$$

Problem: the $Lo$ is the angle between the line throught point $M$ and the Greenwhich meridian. Thus you need to subtract to $\alpha$ the angle $\alpha_G$ between the Greenwich meridian and the $I$ axis. enter image description here This angle is time dependent because of Earth's rotation. While Earth rotates, the $IJK$ reference frame doesn't. Assuming you know how to find this angle (you can use a GMST - Greenwhich Mean Sidereal Time calculator online like this one http://neoprogrammics.com/sidereal_time_calculator/) then the longitude will finally be determinated: $$Lo = \alpha-\alpha_G$$

If you know $\alpha_{G0}$ at a specific time $t_0$ you'll still be able to compute the longitude as $$Lo = \alpha-[\alpha_{G0}+\omega_E(t-t0)]$$

where $\omega_E$ is Earth's rotational speed. Another maybe easier way to compute the longitude requires the knowledge of the satellite position in the $IJK$ frame, so the vector $r(t)$ in the picture: $$\alpha = \arctan(r_J/r_I)$$

Use the arctan2 function to get solutions between $0$ and $2\pi$. I don't know your confidence with different reference systems in astrodynamics and the orbital elements calculations, but some useful references migth be:

  • Bate, Mueller, White: Fundamentals of Astrodynamics
  • Curtis: Orbital Mechanics for Engineering Students
$\endgroup$
4
  • 1
    $\begingroup$ Nice work! The hand-drawn graphics are lovely. The most detailed textbook treatment of coordinate systems, including pictures like these and many others, is probably Vallado, Fundamentals of Astrodynamics and Applications. $\endgroup$
    – Ryan C
    Commented Jul 5, 2022 at 14:22
  • 1
    $\begingroup$ @Ryan Vallado's looks like a pretty complete book, I checked it and it even has a "Converting IJK to latitude and longitude" chapter that might be useful to Dian $\endgroup$
    – JuliusC
    Commented Jul 5, 2022 at 14:33
  • $\begingroup$ Nice work, but I don't think it's what the OP is looking for. The way I read the question, given a set of lat,lon coodinates, fit a trig fuction to them to predict future coordinates. $\endgroup$ Commented Jul 6, 2022 at 15:49
  • $\begingroup$ @JuliusC Thanks! But I am not sure if sin(𝜋/2−i) can be used as the interior angle of a spherical triangle is less than 180 degrees. In addition, the approximation of the sin rule drastically deviates from the true calculations as the sides become comparable to the radius of the celestial sphere, such as at around the longitude of 180 degrees from starting point. $\endgroup$
    – Dian Sheng
    Commented Nov 9, 2022 at 22:18
3
$\begingroup$

I'll repost my answer to Analytical expression for the ground track of the International Space Station to supplement these answers, and also mention this answer to Why does the ISS track appear to be sinusoidal?

(I) use a parametric equation.

If the Earth were not rotating, then we would have something like

\begin{align} x & = \cos \omega (t-t_0)\\ y & = \sin \omega (t-t_0) \ \cos i\\ z & = \sin \omega (t-t_0) \ \sin i\\ \end{align}

where the radius of the orbit is 1, $\omega$ is $2 \pi/T$ and $T$ is the orbital period, and $i$ is the inclination of the orbit.

Then we would have

\begin{align} lon & = \arctan2(y, x) + const\\ lat & = \arcsin(z)\\ \end{align}

If the Earth is rotating then

$$lon = \arctan2(y, x) - \omega_E (t-t_0) + const$$

where $\omega_E$ is $2 \pi/T_D$ and $T_D$ is a sidereal day (23h, 56m, 4s roughly).

Solving this for longitude as a function of latitude looks like some serious work and I am not sure there is an analytical solution.

Instead you can try the parametric equation approach where you first make a hidden table of times, and then solve for $lon(t)$ and $lat(t)$ and plot $lat$ vs $lon$

Here is a plot, I haven't adjusted $t_0$ or $const$ and just used rough values for $\omega$, $\omega_E$ and $i$ but it should be enough to get you stared.

$t_0$ and $const$ represent the known starting conditions of the spacecraft that you are trying to plot; $t_0$ is the time at which it crosses the equator going north, and $const$ is the longitude on the Earth below the spacecraft at that time.

Here is some further reading:

ISS ground track simulation

Python script:

import numpy as np
import matplotlib.pyplot as plt

twopi = 2*np.pi
to_degs, to_rads = 180/np.pi, np.pi/180.

omega = twopi/(92*60)
omega_E = twopi/(23*3600 + 56*60 + 4)

time = 60 * np.arange(101.) # 100 minutes

t0 = 1000. # arbitrary, you can fit this later
inc = 51.
const = 1.0  # arbitrary, you can fit this later

x = np.cos(omega * (time-t0))
y = np.sin(omega * (time-t0)) * np.cos(to_rads*inc)
z = np.sin(omega * (time-t0)) * np.sin(to_rads*inc)

lon = np.arctan2(y, x) - omega_E * (time-t0) + const
lat = np.arcsin(z)

if True:
    plt.figure()
    plt.plot(to_degs*lon, to_degs*lat, '.k')
    plt.xlim(-180, 180)
    plt.ylim(-60, 60)
    #plt.gca().set_aspect('equal')
    plt.show()
$\endgroup$

Not the answer you're looking for? Browse other questions tagged or ask your own question.