TLE Fitting¶

Two-Line Element sets (TLEs) are the standard format for distributing satellite orbital elements, but they degrade in accuracy over time. When higher-fidelity state vectors are available (e.g., from GPS or precision orbit determination), it is useful to fit a new TLE to those states.

satkit.TLE.fit_from_states performs this fit using Levenberg-Marquardt non-linear least-squares optimization. It tunes the TLE orbital parameters to minimize the difference between the input state positions and the SGP4-predicted positions. Input states are assumed to be in the GCRF frame and are internally rotated to the TEME frame used by SGP4.

In [1]:

# Imports
import satkit as sk
import numpy as np
import math as m

# Create a high-precision state
# Altitude for circular orbit
altitude = 450e3

# Radius & velocity
r0 = altitude + sk.consts.earth_radius
v0 = m.sqrt(sk.consts.mu_earth / r0)

# Inclination
inclination = 15 * m.pi / 180.0

# Create the state (3D position in meters, 3D velocity in meters / second)
state0 = np.array([r0, 0, 0, 0, v0 * m.cos(inclination), v0 * m.sin(inclination)])
# Make up an epoch
time0 = sk.time(2024, 3, 15, 13, 0, 0)

# Propagate the state forward by a day with high-precision propagator
res = sk.propagate(state0, time0, time0 + sk.duration(days=1.0))

# Get interpolated states every 10 minutes
times = [time0 + sk.duration(minutes=i) for i in range(0, 1440, 10)]
states = [res.interp(t) for t in times]

# Fit the TLE
(tle, fitresults) = sk.TLE.fit_from_states(states, times, time0 + sk.duration(days=0.5))  # type: ignore

# Print the result
print(tle)
print(fitresults["success"])

# Imports import satkit as sk import numpy as np import math as m # Create a high-precision state # Altitude for circular orbit altitude = 450e3 # Radius & velocity r0 = altitude + sk.consts.earth_radius v0 = m.sqrt(sk.consts.mu_earth / r0) # Inclination inclination = 15 * m.pi / 180.0 # Create the state (3D position in meters, 3D velocity in meters / second) state0 = np.array([r0, 0, 0, 0, v0 * m.cos(inclination), v0 * m.sin(inclination)]) # Make up an epoch time0 = sk.time(2024, 3, 15, 13, 0, 0) # Propagate the state forward by a day with high-precision propagator res = sk.propagate(state0, time0, time0 + sk.duration(days=1.0)) # Get interpolated states every 10 minutes times = [time0 + sk.duration(minutes=i) for i in range(0, 1440, 10)] states = [res.interp(t) for t in times] # Fit the TLE (tle, fitresults) = sk.TLE.fit_from_states(states, times, time0 + sk.duration(days=0.5)) # type: ignore # Print the result print(tle) print(fitresults["success"])

            TLE: none
                         NORAD ID: 00000,
                      Launch Year: 2000,
                            Epoch: 2024-03-16T01:00:00.000000Z,
                  Mean Motion Dot: 0 revs / day^2,
              Mean Motion Dot Dot: 0 revs / day^3,
                             Drag: -0.00015328359493676168,
                      Inclination: 14.99137382650794 deg,
                             RAAN: 355.99766801858095 deg,
                            eccen: 0.001361340632230509,
                   Arg of Perigee: 199.3159526612913 deg,
                     Mean Anomaly: 62.01104282934201 deg,
                      Mean Motion: 15.408691622420529 revs / day
                            Rev #: 0
        

---------------------------------------------------------------------------
KeyError                                  Traceback (most recent call last)
Cell In[1], line 34
     30 (tle, fitresults) = sk.TLE.fit_from_states(states, times, time0 + sk.duration(days=0.5))  # type: ignore
     31 
     32 # Print the result
     33 print(tle)
---> 34 print(fitresults["success"])

KeyError: 'success'

Generate Test Data¶

To demonstrate the fitting, we create a synthetic truth trajectory by propagating a circular orbit at 450 km altitude with satkit's high-precision propagator. We then sample position and velocity states every 10 minutes over one day, and fit a TLE to those samples.

In [2]:

# Compute position errors (differences between TLE & state)

# Get the positions from sgp4
(pteme, vteme) = sk.sgp4(tle, times)
# Rotate positions from TEME to GCRF frame
pgcrf = [sk.frametransform.qteme2gcrf(t) * p for t, p in zip(times, pteme)]
# Take difference between state vector and SGP4 positions, and compute norm
pdiff = [p - s[0:3] for p, s in zip(pgcrf, states)]
pdiff = np.array([np.linalg.norm(p) for p in pdiff])


# Plot position errors
import matplotlib.pyplot as plt
import scienceplots  # noqa: F401
plt.style.use(["science", "no-latex", "../satkit.mplstyle"])
%config InlineBackend.figure_formats = ['svg']

fig, ax = plt.subplots(figsize=(10, 5))
ax.plot([t.as_datetime() for t in times], pdiff, color="black", linewidth=2)
ax.set_xlabel("Time")
ax.set_ylabel("Position Error (m)")
ax.set_title("TLE Fitting Position Errors")
fig.autofmt_xdate()
plt.tight_layout()
plt.show()

# Compute position errors (differences between TLE & state) # Get the positions from sgp4 (pteme, vteme) = sk.sgp4(tle, times) # Rotate positions from TEME to GCRF frame pgcrf = [sk.frametransform.qteme2gcrf(t) * p for t, p in zip(times, pteme)] # Take difference between state vector and SGP4 positions, and compute norm pdiff = [p - s[0:3] for p, s in zip(pgcrf, states)] pdiff = np.array([np.linalg.norm(p) for p in pdiff]) # Plot position errors import matplotlib.pyplot as plt import scienceplots # noqa: F401 plt.style.use(["science", "no-latex", "../satkit.mplstyle"]) %config InlineBackend.figure_formats = ['svg'] fig, ax = plt.subplots(figsize=(10, 5)) ax.plot([t.as_datetime() for t in times], pdiff, color="black", linewidth=2) ax.set_xlabel("Time") ax.set_ylabel("Position Error (m)") ax.set_title("TLE Fitting Position Errors") fig.autofmt_xdate() plt.tight_layout() plt.show()

Evaluate Fit Quality¶

Compare the fitted TLE against the original states by propagating the TLE with SGP4, rotating from TEME to GCRF, and computing position differences. Since TLEs are a simplified analytical model, some residual error is expected even with a perfect fit.