Kalman filter vs Savitzky-Golay for real-time vs batch GPS smoothing

The choice between a Kalman filter and a Savitzky-Golay smoother is really a choice about latency, not accuracy. A Kalman filter is causal — it produces a smoothed estimate from each observation the moment it arrives, using only the past — which makes it the correct tool for streaming pipelines. A Savitzky-Golay filter is centred — it fits a local polynomial using points on both sides of the current sample — which makes it cleaner and cheaper for batch post-processing, but unusable on a live stream because the future points it needs do not exist yet. This guide compares the two across the dimensions that matter in production and gives a complete, runnable implementation of each, both operating on projected metric coordinates.

Why the distinction matters

Both filters remove the same high-frequency positional noise that GNSS receivers inject under multipath and poor satellite geometry, the sources catalogued in the parent GPS precision and error handling reference. Where they differ is when they can act. A live fleet-tracking dashboard must show a smoothed position within a second of the ping arriving, so it cannot wait for future samples; a nightly analytics job has the entire trajectory in memory and can exploit both directions of context.

Picking the wrong tool is a silent quality failure. Run Savitzky-Golay in a stream and you either introduce half-a-window of latency or smooth against a truncated, asymmetric window that distorts the newest points. Run a causal Kalman filter offline and you accept a peak-lagging estimate when a symmetric smoother would have preserved the turn geometry that drift-correction pipelines work hard to keep. The related problem of filling gaps rather than smoothing noise is treated separately in interpolating missing GPS points with Kalman filters.

Comparison

Dimension Kalman filter (causal) Savitzky-Golay (centred)
Latency Zero — emits an estimate per observation Half a window; unusable live
Direction of context Past only Past and future (symmetric)
Best fit Streaming, real-time dashboards Batch, offline analytics
Peak / turn preservation Lags peaks unless Q is raised Preserves peak shape well
Sensor fusion Native — extend the state and matrices Not supported
Tuning surface Process noise Q, measurement noise R Window length, polynomial order
Handles irregular cadence Yes — dt enters the model per step Assumes uniform cadence
Compute cost Per-point matrix ops; cheap incrementally One vectorized C pass; very cheap in bulk
Failure mode Diverges if Q/R mis-scaled Over-smooths peaks if window too wide
Smoothing Filter Decision Path Starting from projected metric coordinates, a latency question branches to a causal Kalman filter for streaming pipelines or a centred Savitzky-Golay smoother for batch pipelines. Projected metric coordinates (UTM) Need smoothing as points arrive? Yes -> Kalman filter causal, streaming, fusible No -> Savitzky-Golay centred, batch, peak-safe

Core smoothing pipeline

Both paths share the first two steps and diverge only at the filter.

  1. Project to a metric CRS. Transform WGS84 to a UTM zone derived from the data centroid so the filter works in metres. Both a Kalman process model and a Savitzky-Golay window are defined in physical distance; angular degrees break both.
  2. Regularize where needed. Savitzky-Golay assumes a uniform cadence, so resample first. The Kalman filter accepts the raw irregular cadence because the time step dt enters its prediction directly.
  3. Apply the matching filter. Causal Kalman for streaming; centred Savitzky-Golay for batch.
  4. Validate peak preservation. Confirm speed peaks and turn angles survived.

Both implementations

Both functions take projected metre coordinates and are fully vectorized where the algorithm allows. The Kalman filter is inherently sequential (each step depends on the last) but uses small fixed-size matrix operations per point; the Savitzky-Golay path is a single SciPy call over the whole array.

PYTHON 
from __future__ import annotations

import numpy as np
import pandas as pd
from scipy.signal import savgol_filter


def smooth_savgol(
    x_m: np.ndarray,
    y_m: np.ndarray,
    window_length: int = 11,
    polyorder: int = 3,
) -> tuple[np.ndarray, np.ndarray]:
    """
    Batch (centred) Savitzky-Golay smoothing of projected coordinates.

    x_m, y_m are easting/northing in a metric CRS (metres) — NEVER raw
    WGS84 degrees, or the window no longer maps to a physical distance.

    Requires a uniform time cadence (resample upstream). Uses both past
    and future samples, so it is offline-only.

    Raises
    ------
    ValueError
        If arrays mismatch, are empty, or are shorter than window_length.
    """
    x_m = np.asarray(x_m, dtype=float)
    y_m = np.asarray(y_m, dtype=float)
    if x_m.shape != y_m.shape:
        raise ValueError("x_m and y_m must share the same shape")
    if x_m.size == 0:
        raise ValueError("Empty coordinate arrays")
    if window_length % 2 == 0:
        window_length += 1  # SG requires an odd window
    if x_m.size < window_length:
        raise ValueError(
            f"Need >= {window_length} points, got {x_m.size}"
        )
    xs = savgol_filter(x_m, window_length, polyorder, mode="interp")
    ys = savgol_filter(y_m, window_length, polyorder, mode="interp")
    return xs, ys


def smooth_kalman_cv(
    x_m: np.ndarray,
    y_m: np.ndarray,
    dt_s: np.ndarray,
    meas_std_m: float = 5.0,
    accel_std_ms2: float = 1.5,
) -> tuple[np.ndarray, np.ndarray]:
    """
    Causal constant-velocity Kalman filter on projected coordinates.

    State is [x, y, vx, vy] in metres and metres/second — a metric CRS is
    mandatory, so project WGS84 to UTM before calling this.

    Parameters
    ----------
    x_m, y_m : easting/northing in metres.
    dt_s     : per-step time delta in seconds (same length as x_m); the
               first element is unused. Handles irregular cadence natively.
    meas_std_m   : receiver horizontal accuracy (1 sigma, metres) -> R.
    accel_std_ms2: expected acceleration disturbance -> process noise Q.

    Returns
    -------
    (x_filt, y_filt) causal position estimates, past-only context.

    Raises
    ------
    ValueError
        If arrays mismatch or are empty.
    """
    x_m = np.asarray(x_m, dtype=float)
    y_m = np.asarray(y_m, dtype=float)
    dt_s = np.asarray(dt_s, dtype=float)
    n = x_m.size
    if not (n == y_m.size == dt_s.size):
        raise ValueError("x_m, y_m, dt_s must share length")
    if n == 0:
        raise ValueError("Empty coordinate arrays")
    if n == 1:
        return x_m.copy(), y_m.copy()  # nothing to filter

    r = meas_std_m ** 2
    R = np.eye(2) * r
    H = np.array([[1.0, 0, 0, 0], [0, 1.0, 0, 0]])

    # Initialise state at the first measurement, zero velocity.
    state = np.array([x_m[0], y_m[0], 0.0, 0.0])
    P = np.eye(4) * (r * 10.0)

    xf = np.empty(n)
    yf = np.empty(n)
    xf[0], yf[0] = x_m[0], y_m[0]

    q = accel_std_ms2 ** 2
    for k in range(1, n):
        dt = max(dt_s[k], 1e-3)  # guard against zero/negative deltas
        F = np.array([
            [1, 0, dt, 0],
            [0, 1, 0, dt],
            [0, 0, 1, 0],
            [0, 0, 0, 1],
        ])
        # Discrete white-noise acceleration process covariance Q.
        dt2, dt3, dt4 = dt * dt, dt ** 3, dt ** 4
        Q = q * np.array([
            [dt4 / 4, 0, dt3 / 2, 0],
            [0, dt4 / 4, 0, dt3 / 2],
            [dt3 / 2, 0, dt2, 0],
            [0, dt3 / 2, 0, dt2],
        ])
        # Predict
        state = F @ state
        P = F @ P @ F.T + Q
        # Update with the measurement z = (x, y)
        z = np.array([x_m[k], y_m[k]])
        y_res = z - H @ state
        S = H @ P @ H.T + R
        K = P @ H.T @ np.linalg.inv(S)
        state = state + K @ y_res
        P = (np.eye(4) - K @ H) @ P
        xf[k], yf[k] = state[0], state[1]
    return xf, yf

Validation block

After either filter, confirm it removed jitter without erasing genuine manoeuvres. Compare the total path length and the retained turn energy against the raw projected track.

PYTHON 
def validate_smoothing(
    x_raw: np.ndarray, y_raw: np.ndarray,
    x_smooth: np.ndarray, y_smooth: np.ndarray,
) -> dict[str, float]:
    """
    Compare raw vs smoothed projected tracks (all inputs in metres).
    Raises AssertionError on over-smoothing; returns summary metrics.
    """
    def path_len(x, y):
        return float(np.sum(np.hypot(np.diff(x), np.diff(y))))

    raw_len = path_len(x_raw, y_raw)
    sm_len = path_len(x_smooth, y_smooth)
    reduction = (raw_len - sm_len) / raw_len * 100 if raw_len else 0.0

    # A healthy noise-removal pass shortens the jittery path by 5-30%.
    # Beyond ~50% the filter is eating real geometry.
    assert reduction < 50, (
        f"Over-smoothing: path length cut by {reduction:.1f}%"
    )
    return {
        "raw_path_m": round(raw_len, 1),
        "smoothed_path_m": round(sm_len, 1),
        "length_reduction_pct": round(reduction, 1),
    }

When to choose which

  • Streaming / live dashboards: Kalman. It is causal, delivers zero-latency estimates, and handles irregular ping intervals through dt.
  • Sensor fusion (IMU, odometer, wheel speed): Kalman. Extend the state vector and measurement matrix; Savitzky-Golay cannot ingest a second sensor.
  • Batch analytics on complete tracks: Savitzky-Golay. Symmetric context preserves peak and turn geometry at the lowest cost.
  • Both: run Kalman online for the live view, then re-smooth in batch for the analytics store once the trajectory is complete.

Common mistakes and gotchas

  • Smoothing in WGS84 degrees. Both filters assume metres. A Kalman Q calibrated in m/s² is meaningless applied to degrees, and a Savitzky-Golay window of “11 samples” no longer maps to a fixed ground distance. Project to a UTM zone first, exactly as CRS transformation best practices require.
  • Using an even Savitzky-Golay window. SciPy raises on even windows. The implementation above auto-increments, but a direct call will fail.
  • Mis-scaling Kalman Q and R. If R (measurement noise) is far too small the filter trusts every noisy fix and barely smooths; if Q (process noise) is far too small it lags every turn. Derive R from the receiver accuracy and tune Q against a ground-truth track per transport mode.
  • Feeding a centred smoother into a stream. Applying Savitzky-Golay to only the trailing points available live produces an asymmetric window that distorts the newest, most important sample. Use Kalman online.
  • Skipping gap handling. Neither filter should smooth across a long signal outage. Split the trajectory at gaps and smooth each segment independently, then rejoin.

Back to GPS Precision & Error Handling

Related in this section