Satellite clock drift uses t instead of (t − t_0c)

[giove-a]

Summary

calc_satellite_clock_drift evaluates the clock-drift polynomial at t instead of (t − t_0c), so it does not match the derivative of the clock-correction polynomial used in correct_clock.

Location

src/sat_time.jl (the calc_satellite_clock_drift method):

function calc_satellite_clock_drift(decoder::GNSSDecoder.GNSSDecoderState, t)
    decoder.data.a_f1 +
    decoder.data.a_f2 * t * 2
end

Why this is wrong

The satellite clock correction in correct_clock is the standard polynomial about the clock reference epoch t_0c:

Δt = a_f0 + a_f1·(t − t_0c) + a_f2·(t − t_0c)² + Δt_rel

Its time derivative (the clock drift) is therefore:

d(Δt)/dt = a_f1 + 2·a_f2·(t − t_0c) + d(Δt_rel)/dt

The current code computes a_f1 + 2·a_f2·t — i.e. it uses t rather than (t − t_0c). The a_f2 term is offset by a spurious constant 2·a_f2·t_0c.

Impact

a_f2 is a real broadcast parameter — GPS L1 C/A carries it in subframe 1 (8-bit, two's complement, LSB 2⁻⁵⁵ s/s²; GNSSDecoder.GPSL1Data.a_f2), and Galileo E1B carries it too.

• In practice, operators almost always broadcast a_f2 = 0 for GPS, so the spurious term is exactly zero and there is no observable effect today.
• Worst case, though, |a_f2| can reach ~3.5e-15 s/s² (≈128·LSB) and t_0c runs up to ~604800 s (end of week), giving a spurious drift of 2·a_f2·t_0c ≈ 4.3e-9 s/s, i.e. c·4.3e-9 ≈ 1.3 m/s of range-rate. Since each satellite has its own a_f2/t_0c, this is not common-mode, so it leaks into the estimated velocity, not just the clock drift, via calc_user_velocity_and_clock_drift.

So the bug is masked in typical data but can be ~1 m/s in the worst case — worth fixing regardless, since it's an unambiguous inconsistency with the model correct_clock differentiates.

Suggested fix

function calc_satellite_clock_drift(decoder::GNSSDecoder.GNSSDecoderState, t)
    decoder.data.a_f1 + 2 * decoder.data.a_f2 * (t - decoder.data.t_0c)
end

Related (secondary)

calc_satellite_clock_drift also omits the derivative of the relativistic periodic term Δt_rel = F·e·√A·sin(E) that correct_clock does include. Its rate, F·e·√A·cos(E)·Ė, is roughly ~1 mm/s and is currently neglected in the drift. Worth considering for completeness, though it is a separate (modeling-completeness) question from the t_0c fix above.