Adding Realistic GPS Noise to Synthetic Vehicle Trajectories

When a synthetic vehicle path is perturbed with naive isotropic Gaussian noise, it fails downstream validation because real GNSS error is temporally correlated, anisotropic, and satellite-geometry dependent — this page gives the per-receiver error budgets and verification benchmarks that make injected noise indistinguishable from a real receiver’s output.

Part of Noise Injection & Stochastic Drift: the parent area defines the ordered transform — project, perturb, re-validate, inverse-project — and this page fills in the numbers that transform needs: which σ and θ to use for which receiver class, and how to prove the result is statistically faithful before it reaches a training set.

Root Cause: Why White Noise Is Detectable

GNSS positional error is a composite signal driven by atmospheric propagation delay, receiver clock bias, satellite ephemeris inaccuracy, and local multipath reflection. Treating these as independent white noise destroys the autocorrelation structure that downstream tracking filters — Kalman, particle, or spline-based — depend on. A model trained against white noise learns nothing about the failure modes that dominate in production: urban-canyon navigation, tunnel exits, and cold-start drift.

The error must be decomposed into three orthogonal components, each a distinct stochastic process with its own time constant:

  1. Slow-varying bias — ionospheric and tropospheric delay plus receiver clock drift. Low-frequency wander over minutes to hours, modeled as a random walk (integrated Wiener process).
  2. Correlated spatial jitter — multipath reflection and signal attenuation in urban canyons, strongly dependent on vehicle heading and surrounding geometry. Best represented by an Ornstein-Uhlenbeck process with environment-dependent reversion speed θ.
  3. High-frequency measurement noise — thermal noise in the receiver front-end and quantization error. Zero-mean Gaussian with variance scaled by satellite geometry (HDOP/VDOP).

Critically, noise must be applied in a local tangent plane (ENU or NED), never directly in latitude/longitude. Applying metric noise to geographic degrees introduces latitude-dependent anisotropy — the same silent-corruption class as coordinate reference system drift, where numbers stay finite and the map still renders so the bug survives review. The injection layer should also operate on the kinematic state vector — position, velocity, heading — so it never violates the constraints enforced by physics-based path generation upstream.

Three-component GNSS noise decomposition summed onto a clean baseline path Top row: a clean, kinematically valid baseline path. Second row adds slow bias drift modeled as a random walk scaled by sigma_b in metres per minute, wandering over minutes to hours. Third row adds correlated multipath jitter modeled as an Ornstein-Uhlenbeck process with reversion rate theta and a 10 to 60 second time constant. Fourth row adds a high-frequency, zero-mean Gaussian measurement floor with standard deviation sigma_base times HDOP, applied per fix. The bottom row sums all four into the final perturbed trace, shown hugging a dashed copy of the clean baseline, which clears the lag-1 autocorrelation and one-over-f power-spectral-density verification gates. Clean baseline path kinematically valid · physics-checked in a local metric frame (ENU/UTM) + Bias drift random walk · σ_b (m/min) τ ≈ minutes → hours · iono/tropo/clock + Multipath jitter Ornstein-Uhlenbeck · reversion θ (s⁻¹) τ ≈ 10–60 s · heading-aligned, correlated + Measurement floor Gaussian · σ = σ_base · HDOP white, per-fix · thermal + quantization = Perturbed trace clears lag-1 autocorr (r₁>0.2) and 1/f PSD shape gates - - - clean baseline (reference)

Per-Receiver Error Budgets

The parameters below anchor each process to a hardware class. Use these as defaults, then calibrate against a real error trace from the target receiver.

Receiver class Bias σ_b (m/min) OU reversion θ (s⁻¹) Measurement σ_base (m) Multipath jump prob.
Survey / RTK (L1+L2) 0.05 2.0 0.3 <0.1%
Consumer phone (open sky) 0.2 1.0 1.5 0.5%
Consumer phone (urban canyon) 0.4 0.2 3.0 1.0%
Low-cost / IMU-fused 0.5 0.1 3.0 1.5%

The measurement-noise variance is then scaled by simulated satellite geometry, σpos2=σbase2⋅HDOP2\sigma_{pos}^2 = \sigma_{base}^2 \cdot \text{HDOP}^2, so the same receiver shows tighter fixes under good geometry and wider scatter when satellites cluster.

Prerequisite Check

Confirm the projection stack is present and that pyproj carries its own PROJ data — a mismatched transform pipeline is a silent reproducibility hazard.

python
import numpy as np
import pyproj

print("numpy", np.__version__)        # 1.26.x or newer
print("pyproj", pyproj.__version__)   # pin to 3.x
print("PROJ", pyproj.proj_version_str)  # 9.x bundled with pyproj 3.x

# A metre must be a metre before any noise is applied:
assert pyproj.CRS.from_epsg(4326).is_geographic   # degrees — never perturb here

Fix: Three-Component Noise Injection

The injector consumes a clean path in EPSG:4326, a receiver profile, and a deterministic seed, then returns perturbed geodetic coordinates. Each stochastic process is seeded from one base seed so the whole fleet reconstructs byte-for-byte.

python
import numpy as np
from pyproj import Transformer

def inject_gnss_noise(lon, lat, t, profile, seed):
    """Perturb a clean path with correlated GNSS error. Arrays are 1-D, time in seconds."""
    rng = np.random.default_rng(seed)
    dt = np.diff(t, prepend=t[0] - 1.0)              # per-fix interval

    # --- Step 1: project EPSG:4326 -> local UTM so noise is metric ---
    utm = pyproj.database.query_utm_crs_info(
        "WGS 84",
        pyproj.aoi.AreaOfInterest(lon.min(), lat.min(), lon.max(), lat.max()),
    )[0]
    fwd = Transformer.from_crs("EPSG:4326", utm.code, always_xy=True)
    inv = Transformer.from_crs(utm.code, "EPSG:4326", always_xy=True)
    e, n = fwd.transform(lon, lat)

    # --- Step 2: slow bias drift (random walk, m/min -> m/s) ---
    sb = profile["sigma_b"] / 60.0
    bias_e = np.cumsum(sb * rng.standard_normal(len(t)) * np.sqrt(dt))
    bias_n = np.cumsum(sb * rng.standard_normal(len(t)) * np.sqrt(dt))

    # --- Step 3: OU multipath jitter, exact discretization ---
    theta, sigma_j = profile["theta"], profile["sigma_base"]
    J_e = np.zeros(len(t)); J_n = np.zeros(len(t))
    for k in range(1, len(t)):
        a = np.exp(-theta * dt[k])
        sd = sigma_j * np.sqrt((1 - a**2) / (2 * theta))
        J_e[k] = a * J_e[k-1] + sd * rng.standard_normal()
        J_n[k] = a * J_n[k-1] + sd * rng.standard_normal()

    # --- Step 4: HDOP-scaled high-frequency measurement noise ---
    hdop = profile.get("hdop", 1.0)
    meas = profile["sigma_base"] * hdop
    noise_e = meas * rng.standard_normal(len(t))
    noise_n = meas * rng.standard_normal(len(t))

    # --- Step 5: rare heavy-tailed multipath jumps ---
    jump = rng.random(len(t)) < profile["jump_prob"]
    jump_mag = rng.standard_cauchy(len(t)) * meas  # heavy tail
    noise_e += jump * jump_mag
    noise_n += jump * rng.standard_cauchy(len(t)) * meas * jump

    e_p = e + bias_e + J_e + noise_e
    n_p = n + bias_n + J_n + noise_n
    lon_p, lat_p = inv.transform(e_p, n_p)
    return lon_p, lat_p

python
URBAN = {"sigma_b": 0.4, "theta": 0.2, "sigma_base": 3.0,
         "jump_prob": 0.01, "hdop": 1.8}
lon_p, lat_p = inject_gnss_noise(lon, lat, t, URBAN, seed=20260402)

Verification Step

Injected noise must clear three CI gates before it enters any training or routing pipeline. Each gate catches a specific way naive noise betrays itself.

python
def verify_noise(clean_e, pert_e, dt=1.0):
    err = pert_e - clean_e

    # 1. Autocorrelation decay — white noise drops to ~0 at lag 1.
    r1 = np.corrcoef(err[:-1], err[1:])[0, 1]
    assert r1 > 0.2, f"no temporal correlation (r1={r1:.3f}) — white noise leaked in"

    # 2. Spectral shape — real GNSS error is 1/f at low frequency, not flat.
    psd = np.abs(np.fft.rfft(err - err.mean())) ** 2
    lo, hi = psd[1:5].mean(), psd[len(psd)//2:].mean()
    assert lo > 2 * hi, "flat PSD — error has no low-frequency drift"

    # 3. Kinematic continuity — derived accel must stay vehicle-plausible.
    v = np.gradient(pert_e, dt)
    a = np.gradient(v, dt)
    assert np.abs(a).max() < 4.0, "acceleration spike — perturbed beyond vehicle dynamics"
    return {"r1": r1, "psd_ratio": lo / hi, "amax": np.abs(a).max()}

Wire verify_noise into the same gate that runs on every regression build. The autocorrelation assertion is the one that matters most: it is the single check that makes it impossible for a refactor to silently regress the OU process back to independent Gaussian draws. Calibrate the real GNSS time constant against the target hardware — position error typically decays exponentially with a 10–60 s constant.

Edge Cases & Gotchas

  • Antimeridian and UTM-zone straddling. A trajectory that crosses ±180° longitude, or one that spans two UTM zones, will pick a single zone via query_utm_crs_info and accumulate distortion at the far edge. For long-haul routes, perturb in a local ENU tangent plane anchored at the path centroid instead of a fixed UTM zone, or segment the path per zone and stitch the perturbed pieces.
  • Null Island and unset coordinates. A (0.0, 0.0) fix from a failed lock sits off the African coast and will be perturbed as if real, injecting plausible-looking garbage. Drop or interpolate sentinel (0, 0) and NaN fixes before injection — never let the noise layer dignify a dropout with realistic drift.
  • Float precision at low latitude vs. high latitude. Inverse-projecting metre-scale perturbations back to degrees loses precision differently near the poles, where a degree of longitude collapses. Keep all arithmetic in float64 and inverse-project only once, at export — repeated round-trips through EPSG:4326 accumulate sub-metre bias that the PSD gate will eventually flag.

Frequently Asked Questions

Why not just add Gaussian noise to every fix?

Because independent Gaussian draws have zero temporal autocorrelation, so any Kalman filter or moving average downstream erases them instantly and the model sees effectively clean paths. Real GNSS error drifts and reverts; reserve plain Gaussian noise for the thin HDOP-scaled measurement floor layered under the OU drift, never as the whole model.

How do I choose the OU reversion rate θ for a given vehicle?

Match θ to how long the receiver’s bias persists. A survey-grade unit corrects quickly (large θ, error reverts within a second or two); a phone in an urban canyon drifts for tens of seconds (small θ). Fit the OU autocorrelation decay to a real error trace from the target hardware class, then confirm the synthetic lag-1 autocorrelation lands in the same band.

How do I scale measurement noise with satellite geometry?

Multiply the base variance by HDOP², so σpos2=σbase2⋅HDOP2\sigma_{pos}^2 = \sigma_{base}^2 \cdot \text{HDOP}^2. Drive HDOP from a simulated sky-occlusion model — high values when buildings block satellites, low values in open sky — so the same receiver tightens and widens realistically along the route.

Does adding GPS noise make a trajectory dataset private?

No. Perturbation reduces re-identification risk but does not bound it. Pair it with formal differential privacy mechanisms and a spatial k-anonymity gate, and confirm the perturbation magnitude clears the minimum safe distance for every distinctive trace before export.