Im trying to get the historical Lon/Lat/Alt of a starlink satellite. I'm using TLE's from space-track and Skyfield to calculate latitude, longitude and altitude of satellites for a given date range. The methodology I'm using:
- Get TLE for date range 06/01 - 06/07 (say you get 30 data points). EDIT: got 19 historical TLE's.
- Use Skyfield to generate the LLA (lat,lon,alt) for those 30 data points + timestamps between those 30 TLE's.
- Plot the altitude as predicted using Skyfield (that uses SGP4 under the hood?) over time.
Here's what I got for STARLINK-55775 (y-axis in kilometres).
I'm not familiar with astrodynamics, just trying to learn something new.
--EDITS--
Code:
from sgp4.api import Satrec
from sgp4.api import SGP4_ERRORS
from skyfield.api import EarthSatellite
import requests
from skyfield.api import load, wgs84
import pandas as pd
import datetime
from dateutil import parser, tz
import pytz
import plotly.express as px
def create_session():
user_name = 'username'
password = 'password'
base_url = 'https://www.space-track.org/'
request_login = "/ajaxauth/login"
request_cmd_action = "/basicspacedata/query"
request_find_starlinks = "/class/tle_latest/NORAD_CAT_ID/53871/ORDINAL/1/OBJECT_NAME/STARLINK~~/format/json/orderby/NORAD_CAT_ID%20asc"
session = requests.Session()
# need to log in first. note that we get a 200 to say the web site got the data, not that we are logged in
resp = session.post(base_url + request_login, data={'identity': user_name, 'password': password})
if resp.status_code != 200:
raise MyError(resp, "POST fail on login")
# this query picks up all Starlink satellites from the catalog. Note - a 401 failure shows you have bad credentials
resp = session.get(base_url + request_cmd_action + request_find_starlinks)
if resp.status_code != 200:
print(resp)
raise MyError(resp, "GET fail on request for Starlink satellites")
else:
print("########### SUCCESSFULLY SET UP SESSION ###############")
return session
def get_tle(start_dt, end_dt, session, id=53871, ):
url = f'https://www.space-track.org/basicspacedata/query/class/gp_history/NORAD_CAT_ID/{id}/orderby/EPOCH ASC/EPOCH/{start_dt}--{end_dt}/format/json'
result = session.get(url) #logs into space-track to get data.
if result.status_code == 200:
docs = result.json()
return_doc = []
for doc in docs:
_t = {}
_t['timestamp'] = doc['EPOCH']
_t['TLE_LINE1'] = doc['TLE_LINE1']
_t['TLE_LINE2'] = doc['TLE_LINE2']
return_doc.append(_t)
return return_doc
else:
print(result.text)
def get_satellite_lla(tle_1, tle_2, timestamp, **kwargs):
ts = load.timescale()
if type(timestamp) == str:
t = ts.utc(parser.parse(timestamp))
else:
t = ts.utc(timestamp)
satellite = EarthSatellite(tle_1, tle_2, 'StarLink')
geocentric = satellite.at(t)
lat, lon = wgs84.latlon_of(geocentric)
alt = wgs84.height_of(geocentric).km
session.close()
return pd.Series([lat.degrees, lon.degrees, alt], index=['lat', 'lon', 'alt'])
session = create_session()
id = 56529
tle_data = get_tle(id=id,start_dt='2023-05-31', end_dt='2023-06-07', session=session)
df = pd.DataFrame.from_dict(tle_data)
df.drop_duplicates(subset='timestamp', inplace=True)
resampled = df.copy()
resampled.timestamp = resampled.timestamp.apply(lambda x: parser.parse(x + ' UTC'))
resampled = resampled.set_index('timestamp').resample('5T').ffill().reset_index().dropna() # forward filling the TLE data for newly generated timestamps.
resampled[['lat','lon', 'alt']] = resampled.apply(lambda x: get_satellite_lla(tle_1=x.TLE_LINE1, tle_2=x.TLE_LINE2, timestamp=x.timestamp), axis=1)
resampled.to_csv(f"Starlink_{id}.csv", index=False)
fix = px.line(resampled, x='timestamp', y='alt')
fix.show()
session.close()
Chose a different satellite, which isn't being raised (but hasn't been put into station yet) and changed the upsampling from 30minutes to 5 minutes.
The satellite still seems to swing up and down by close to 10km!
