Recover a receiver's velocity and clock drift from one epoch of Doppler or pseudorange-rate measurements against a precise (SP3) ephemeris source.

This is the velocity counterpart to Orbis.GNSS.Positioning: where single-point positioning inverts pseudoranges for position and clock bias, this module inverts pseudorange rates (or equivalently Doppler shifts) for the receiver velocity and clock drift, given a known receiver position. It reuses the forward model in Orbis.GNSS.Observables for the satellite geometry and velocity, and the explicit 4x4 inverse in Orbis.GNSS.Geometry.inv4/1 for the normal-equation solution.

Observation model (standard Doppler / range-rate least squares)

For satellite i with line-of-sight unit vector e_i (receiver -> satellite, ECEF), satellite velocity v_sat_i, receiver velocity v_rx, the measured pseudorange rate is

rho_dot_i = e_i . (v_sat_i - v_rx) + c * (rx_clock_drift - sat_clock_drift_i)

Collecting the four unknowns into x = [v_rx_x, v_rx_y, v_rx_z, c * rx_clock_drift] (the clock term carried in length-rate units, m/s), the model rearranges to a linear system

e_i . v_sat_i - rho_dot_i - c * sat_clock_drift_i = e_i . v_rx - c * rx_clock_drift

which is one row H_i . x = y_i with

H_i = [-e_ix, -e_iy, -e_iz, 1]
y_i = rho_dot_i - (e_i . v_sat_i) + c * sat_clock_drift_i

Stacking the visible satellites and solving the (unweighted) normal equations gives x = (H^T H)^-1 H^T y, from which the receiver velocity is {x0, x1, x2} and the clock drift is x3 / c (seconds per second).

How the geometry terms are obtained

Every quantity comes from Orbis.GNSS.Observables.predict/5 evaluated at the known receiver position, which is treated as static for the forward model:

• los_unit is e_i, the receiver -> satellite ECEF unit vector;
• range_rate_m_s for a static receiver is exactly e_i . v_sat_i (the los . (v_sat - v_rx) projection with v_rx = 0), already expressed in the consistent receive-epoch ECEF frame with light-time and Sagnac corrections applied. The satellite velocity is finite-differenced inside predict/5, so it is never accessed directly here.

Satellite clock drift

The per-satellite term c * sat_clock_drift_i is a known correction folded into y_i. By default it is taken as zero: the SP3 forward model exposes the satellite clock offset but not its time derivative, and a constant common satellite-clock-drift bias is absorbed by the estimated receiver clock drift, exactly as a common clock offset is absorbed by the receiver clock bias in single-point positioning. Callers who have per-satellite drift estimates may supply them through opts[:sat_clock_drift] (a %{sat => drift_s_s} map or a one-argument function of the satellite id); the same value is then subtracted from y_i.

Doppler observations

A Doppler shift relates to the pseudorange rate by doppler_hz = -rho_dot * f / c (the convention in Orbis.GNSS.Observables), so with opts[:observable] = :doppler each measurement is converted via rho_dot = -doppler_hz * c / f before the solve, with the carrier f taken from opts[:carrier_hz] (default L1, 1575.42 MHz).

Result map

%{
  velocity_m_s:    {vx, vy, vz},   # receiver velocity, ECEF m/s
  speed_m_s:       float(),        # norm of velocity_m_s
  clock_drift_s_s: float(),        # receiver clock drift, s/s
  residuals_m_s:   %{sat => r},    # post-fit range-rate residuals, m/s
  used_sats:       [sat],          # satellites contributing, build order
  n_satellites:    non_neg_integer()
}