How satpass works

This document covers the chain from a two-line element set to "the ISS is at magnitude −3.4, 52° up in the north-west". It goes through the propagator, the reference frames, illumination, brightness, pass search and accuracy, with worked numbers from the bundled data.

TLE text ──parse──▶ mean elements ──SGP4/SDP4──▶ r,v (TEME) ──GMST──▶ r,v (ECEF)
                                                     │                    │
                                     Sun (MOD≈TEME) ─┤                    ├─ observer (WGS-84) ─▶ SEZ ─▶ az/el/range/rate
                                                     ▼                    ▼
                                         shadow cone, phase angle    sub-point lat/lon/alt
                                                     │
                                         visual magnitude, visibility

1. TLE decoding (tle.py)

A TLE has two 69-column lines and may carry a name line before them. Column 69 holds a modulo-10 checksum: you add up every digit in columns 1–68 and count each - as 1.

1 25544U 98067A   26288.82192346  .00005580  00000+0  10980-3 0  9999
2 25544  51.6317  63.1136 0006755 265.8514  94.1703 15.48879905590452

• Epoch 26288.82192346: years 57–99 belong to the 1900s and 00–56 to the 2000s. Day 288.82192346 of 2026 is 2026-10-15 19:43:34.187 UTC. The Julian date comes from "January 0" of that year (JD 2461040.5) plus the day number, which gives 2461329.32192346.
• B* 10980-3 uses an implied decimal point: 0.10980 × 10⁻³ = 1.098e-4 per Earth radius.
• Eccentricity 0006755 also has an implied leading decimal point: 0.0006755.
• Mean motion is 15.48879905 rev/day. SGP4 wants radians per minute, so it becomes 15.48879905 × 2π / 1440 = 0.067582 rad/min.

2. SGP4/SDP4 (sgp4.py, deepspace.py, lunisolar.py)

The implementation follows Vallado, Crawford, Hujsak & Kelso, Revisiting Spacetrack Report #3 (AIAA 2006-6753), in the "improved" operation mode with WGS-72 constants (μ = 398600.8 km³/s², Rₑ = 6378.135 km, J2/J3/J4). TLEs are fitted with WGS-72, so propagating them with WGS-84 constants would add error rather than remove it.

Initialisation

1. Un-Kozai the mean motion. TLE mean motion follows Kozai's convention, but SGP4 works with Brouwer's mean motion. The code takes one δ-correction to get a₁, refines it, and then computes n₀'' = n₀ / (1 + δ₀).
2. Atmosphere parameters. The density function uses s = 78 km and q₀ = 120 km. When perigee is below 156 km, s is lowered: it becomes perigee − 78 km, with a floor of 20 km for perigees below 98 km. This is how the reference code handles re-entering objects.
3. Secular rates. It computes the J2/J4 rates of mean anomaly, argument of perigee and node (mdot, argpdot, nodedot) and the drag coefficients C1…C5 and D2…D4. When perigee is below 220 km, the "simple" drag model is used (isimp = 1) and the higher-order terms are dropped.
4. Choose the theory. A period of 225 min or more switches to SDP4 (method = 'd'): * dscom computes the solar and lunar coefficients from the Sun's and Moon's mean elements at the epoch. * dsinit classifies resonance:
   • irez = 1 covers 24-hour geosynchronous orbits (n between 0.0034906585 and 0.0052359877 rad/min).
   • irez = 2 covers 12-hour Molniya-type orbits (n between 8.26e-3 and 9.24e-3 rad/min with e ≥ 0.5). This case uses the G-function polynomials in e, split at e = 0.65 and e = 0.715.

Propagation (per time step t minutes from epoch)

• Secular gravity and drag come first. Drag makes a shrink by (1 − C1 t − D2 t² …)² and adds the t²…t⁵ terms to the mean longitude.
• Deep space (dspace) adds the secular lunar-solar rates. For resonant orbits it also integrates the resonance terms (xli, xni) with a fixed 720-minute Euler-Maclaurin step. That integrator keeps its state between calls, as the reference code does, and restarts from the epoch whenever the direction of time changes. A test checks that propagating forwards and backwards gives the same answers.
• Long-period periodics (dpper) add the lunar-solar perturbations. Below 0.2 rad of inclination the Lyddane form is used to avoid dividing by sin i. When the perturbed inclination goes negative, the code flips it and shifts the node and argument of perigee by π.
• Kepler's equation is solved in its equinoctial form, U = E + aᵧₙ cos E − aₓₙ sin E, with Newton iterations. Each step is capped at 0.95 rad and the loop stops at 10 iterations or |Δ| < 1e-12.
• Short-period J2 terms are added last. The result is rotated into TEME with the orientation vectors U and V.

The code raises SGP4Error with the reference error codes:

Verification

tests/data embeds 22 satellites from Vallado's SGP4-VER.TLE together with the matching lines of tcppver.out. Together they cover:

• near-Earth objects with normal drag
• objects with perigee below 220 km (simple drag) and below 98 km (the s* floor)
• SDP4 cases without resonance, with 12-hour resonance at e < 0.65, between 0.65 and 0.715, and above 0.715, and with 24-hour GEO resonance
• the Lyddane low-inclination cases
• a decayed object that has to fail with error 6

Every position matches to better than 1 m and every velocity to better than 1 mm/s; the actual agreement is about 1e-9 km. Modern TLEs from the bundled files are also checked against the sgp4 PyPI package, which is used only as a test oracle, at random times from −3 to +10 days.

Worked example, satellite 00005 at t = 0:

>>> parse_tle(l1, l2).propagator.propagate(0.0)
((7022.46529266, -1400.08296755, 0.03995155), (1.893841015, 6.405893759, 4.534807250))

These are exactly the first line of tcppver.out.

3. Frames (frames.py)

TEME → ECEF. SGP4 outputs positions in the True Equator, Mean Equinox frame. Rotating about z by the Greenwich Mean Sidereal Time (IAU-82 polynomial, as in Vallado's gstime) gives the Pseudo-Earth-Fixed frame:

θ = 67310.54841 s + (876600 h + 8640184.812866 s)·T + 0.093104·T² − 6.2e-6·T³    (T = Julian centuries UT1)
r_ecef = Rz(θ) r_teme
v_ecef = Rz(θ) v_teme − ω⊕ × r_ecef            (ω⊕ = 7.292115146706979e-5 rad/s)

Two simplifications are made on purpose:

• Polar motion is ignored. It is about 0.3″, or roughly 10 m at the Earth's surface.
• UT1 is taken as UTC. They differ by less than 0.9 s, which is at most about 0.4 km along the equator.

Both are smaller than the error in any TLE more than a few hours old, and they let the tool run without an Earth-orientation table. Vallado's TEME example (r = 5094.18, 6127.64, 6380.34 km on 2004-04-06 07:51:28.386 UTC) reproduces the published ITRF vector to within 30 m when the published UT1 is passed in.

Geodetic → ECEF. This uses the WGS-84 ellipsoid (a = 6378.137 km, f = 1/298.257223563) and N = a/√(1 − e² sin²φ). The inverse used for the sub-satellite point iterates on latitude until the change is below 1e-12 rad. Round trips are exact to about 1 µm, up to GEO altitudes and near the poles.

Topocentric SEZ. The range vector ρ = r_sat − r_obs is rotated into South-East-Zenith coordinates:

S =  sinφ cosλ ρx + sinφ sinλ ρy − cosφ ρz
E = −sinλ ρx + cosλ ρy
Z =  cosφ cosλ ρx + cosφ sinλ ρy + sinφ ρz
el = asin(Z/|ρ|),  az = atan2(E, −S)  (north = 0°, east = 90°),  range-rate = ρ·v_ecef/|ρ|

4. Sun, shadow and twilight (sun.py)

Sun position. This is Vallado's low-precision algorithm, good to about 0.01°. It uses the mean longitude, mean anomaly, ecliptic longitude with two equation-of-centre terms, a distance series in AU, and the obliquity of the ecliptic. The output is a mean-of-date equatorial vector, which is treated as TEME because the two frames differ only by nutation (under 20″). Example 5-1 (2006-04-02 0h) gives (0.9771945, 0.1924424, 0.0834308) AU, matching the book to 2e-5 AU.

Earth shadow, conical model. Seen from the satellite, the Sun is a disc of angular radius a = asin(R☉/|r☉ − r|) and the Earth is a disc of angular radius b = asin(Rₑ/|r|). Their centres are separated by c. The lit fraction is one minus the overlapping area of the two discs divided by the Sun's disc:

This handles both the umbra cone and the penumbra cone without working out where the cone apex is. For the ISS at 420 km, b = 69.8° and a = 0.27°, so the satellite crosses the penumbra in about 8–10 s; the tests check that the fraction falls monotonically across it. The classic cylindrical test is also provided: the satellite is shadowed if it is behind the Earth and less than one Earth radius from the Sun–Earth axis. A test checks that the two models agree everywhere except in the penumbra band.

A satellite counts as sunlit when at least half of the solar disc is visible. The observer twilight check uses the Sun's geometric altitude: "dark enough" means below −6° (the end of civil twilight). The now command also names nautical and astronomical twilight.

5. Visual magnitude (magnitude.py)

The satellite is modelled as a diffuse (Lambertian) sphere:

m = m_std + 5·log10(range / 1000 km) − 2.5·log10( F(φ) · lit_fraction ),
F(φ) = sin φ + (π − φ)·cos φ        (normalised so that F(90°) = 1)

Here φ is the phase angle, the angle at the satellite between the Sun and the observer. m_std is the standard magnitude at 1000 km range and 90° phase (the McCants / Heavens-Above convention). Built-in values:

The 4.0 default is typical of the large rocket bodies in CelesTrak's "visual" group.

Worked example: the ISS at 600 km range and 60° phase.

5·log10(0.6)            = −1.109
F(60°) = 0.866 + 2.094·0.5 = 1.913  →  −2.5·log10(1.913) = −0.704
m = −1.8 − 1.109 − 0.704 = −3.6

The model ignores specular flares (such as Iridium-style glints or solar-panel glints) and the satellite's orientation. Expect ±1 mag scatter against what you actually see.

6. Pass search (passes.py, rootfind.py)

1. Coarse scan. Elevation is sampled every 60 s for LEO objects, scaled down to 30 s for very short periods, and every 600 s for deep-space objects. A 60 s step cannot miss a pass above 10°: such an ISS pass lasts more than 4 minutes. Only grazing passes of under a minute near 0° could slip through.
2. AOS/LOS refinement. Each change of sign in elevation brackets a horizon crossing. Brent's method (Brent-Dekker: inverse quadratic interpolation, secant, bisection) refines it to 0.1 s in a handful of evaluations.
3. Culmination. The highest coarse sample and its two neighbours bracket the maximum. Brent's parabolic/golden-section maximiser refines it to 0.05 s. The tests compare against a 1-second brute-force scan driven by the oracle sgp4 package and agree to within 1.5 s on AOS/LOS and 0.01° on peak elevation.
4. Visibility windows. Inside each pass, the code samples every 5–30 s. A point counts as visible when the satellite is sunlit and the Sun is below −6°. The transitions are refined by bisection to about 0.01 s; bisection is used here because the predicate is boolean. Each window records its brightest magnitude and its highest elevation.

If the Sun is above −6° at AOS, at culmination and at LOS, the window search is skipped for that pass. This makes flyover-list fast enough to scan all 155 objects in the "visual" group for a whole night in under a second. 5. Clipping. Passes still in progress at the start or end of the search window are flagged with < and >. Geostationary satellites that never set appear as one clipped pass.

flyover-list --tonight first finds the next dark interval: it scans the Sun's altitude in 10-minute steps and refines the −6° crossings with Brent's method. It then lists every visible pass in that interval that is brighter than --max-mag.

7. Accuracy and TLE age

For an observer, along-track error shows up as a timing error: the ISS moves at 7.66 km/s, so 10 km is about 1.3 s. But reboosts and changes in solar activity alter drag in ways a TLE cannot predict. In practice, use elements less than about three days old for minute-level planning. That is why satpass warns at ±7 days and loudly at ±30 days.

The bundled sample TLEs have epochs of 2026-10-13…16. Predictions for late October 2026 are fine; for any later date, run satpass tle update.