Thistle

Satellite orbit propagation and data generation built on SGP4 and Skyfield. Handles automatic TLE switching for long-duration propagation, and generates orbit data, ground site ranges, and satellite events from a single propagation pass.

Installation

pip install thistle

Requires Python >= 3.9.

Optional extras:

pip install 'thistle[cli]'   # the thistle command-line tool (Python >= 3.10)
pip install 'thistle[plot]'  # CLI plus plot/groundtrack commands
pip install 'thistle[db]'    # CLI plus thistle-db database fallback (Python >= 3.11)

With thistle[db] installed and a thistle-db database configured (THISTLE_DB_* environment variables or a config.toml named by THISTLE_DB_CONFIG), CLI commands resolve NORAD catalog IDs from the database whenever no matching file is found via THISTLE_TLE_DIR. See thistle config for the effective settings.

Quick start

import numpy as np
from thistle import Propagator, read_tle, generate

tles = read_tle("satellite_data.tle")
prop = Propagator(tles, method="midpoint")

times = np.arange(
    np.datetime64("2024-01-15T00:00"),
    np.datetime64("2024-01-16T00:00"),
    np.timedelta64(60, "s"),
)

data = generate(times, prop, ["eci", "lla"])
# data["eci_x"], data["lat"], data["lon"], ...

Propagator

A single TLE degrades in accuracy as you propagate further from its epoch. For tracking a satellite over days or weeks you need multiple TLEs and a strategy for switching between them.

from thistle import Propagator, read_tle

tles = read_tle("satellite_data.tle")
prop = Propagator(tles, method="midpoint")

Switching strategies

You can also pass a strategy instance directly:

from thistle import Propagator, EpochSwitchStrategy, read_tle

tles = read_tle("satellite_data.tle")
prop = Propagator(tles, method=EpochSwitchStrategy([]))

This is useful for subclassing SwitchingStrategy to implement custom switching logic.

Lookup methods

time = np.datetime64("2024-01-15T12:00:00")

satellite = prop.find_satellite(time)  # Skyfield EarthSatellite
satrec = prop.find_satrec(time)        # sgp4 Satrec
tle = prop.find_tle(time)              # (line1, line2) strings

Direct propagation

The Propagator can be used as a drop-in replacement for a Skyfield EarthSatellite:

from skyfield.api import load
ts = load.timescale()

geo = prop.at(ts.utc(2024, 1, 15, 12, 0, 0))

Data generation

generate() propagates once and extracts multiple data groups in a single pass.

from thistle import generate

data = generate(times, prop, ["eci", "lla", "keplerian", "sunlight"])

Groups

All values are returned as NumPy arrays. Angles and dimensionless quantities are float32; positions, velocities, and range data are float64.

Ground site range

Pass sites to generate() to compute slant range and range rate without a second propagation:

data = generate(times, prop, ["lla"], sites={"ksc": (28.57, -80.65)})
# data["range_ksc"]       -> slant range in meters
# data["range_rate_ksc"]  -> range rate in m/s

Sites can also be a list (keys become range_0, range_1, ...):

data = generate(times, prop, ["lla"], sites=[(28.57, -80.65), (34.05, -118.24)])

Each site tuple is (lat, lon) or (lat, lon, alt_m).

The standalone generate_range() function is also available:

from thistle import generate_range

rng = generate_range(times, prop, sites={"ksc": (28.57, -80.65)})

Doppler shift

from thistle import doppler_shift

doppler_hz = doppler_shift(data["range_rate_ksc"], freq=437e6)

Events

Find satellite events within a time window. All event functions accept either an EarthSatellite or a Propagator and return lists of dicts.

Passes

from thistle import find_passes

passes = find_passes(start, stop, prop, lat=28.57, lon=-80.65, min_elevation=10.0)
for p in passes:
    print(p["start"], p["stop"], p["peak_elevation"])

Returns dicts with keys: start, stop, peak_time, peak_elevation.

Node crossings

from thistle import find_node_crossings

crossings = find_node_crossings(start, stop, prop)
for c in crossings:
    print(c["start"], c["longitude"], c["ascending"])

Returns dicts with keys: start, stop, longitude, ascending.

Sunlit and eclipse periods

from thistle import find_sunlit_periods, find_eclipse_periods

sunlit = find_sunlit_periods(start, stop, prop)
eclipse = find_eclipse_periods(start, stop, prop)

Ascending and descending periods

from thistle import find_ascending_periods, find_descending_periods

ascending = find_ascending_periods(start, stop, prop)
descending = find_descending_periods(start, stop, prop)

All period functions return dicts with keys: start, stop.

Visibility circle

Compute the ground footprint where a satellite at a given altitude is visible above a minimum elevation angle:

from thistle import visibility_circle

lats, lons = visibility_circle(28.57, -80.65, alt=0.0, sat_alt=408_000, min_el=10.0)

Accuracy

TLE propagation

SGP4 propagation from TLEs is inherently limited to roughly 1-10 km near epoch, degrading to 10-100 km over days. For sub-kilometer accuracy, use precision ephemerides (SP3) instead of TLEs.

Extrapolation warnings

When requested times are more than 7 days from the nearest loaded TLE epoch -- before the first TLE, after the last, or inside a coverage gap -- Propagator emits a thistle.TLEExtrapolationWarning, since SGP4 results that far from epoch are questionable. Tune or disable it per propagator, or with the standard warnings machinery:

import warnings
from datetime import timedelta
from thistle import Propagator, TLEExtrapolationWarning

prop = Propagator(tles, warn_threshold=timedelta(days=14))  # widen
prop = Propagator(tles, warn_threshold=None)                # disable

warnings.filterwarnings("ignore", category=TLEExtrapolationWarning)  # global
warnings.simplefilter("error", TLEExtrapolationWarning)  # make it fatal

The CLI commands that propagate (propagate, find-tle, groundtrack) print the same warning once per satellite on stderr and accept --warn-days N (0 disables).

Time scale

Thistle treats UTC as UT1, introducing ~0.2-0.9 seconds of time offset (~1-7 m position error for LEO). This is negligible compared to TLE propagation errors and avoids the ~12x overhead of proper UTC/leap-second handling.

Coordinate system

ECI outputs use the ICRF (International Celestial Reference Frame), equivalent to J2000 for most applications.

License

MIT. See LICENSE for details.