Modeling GPS Multipath Error in Urban Canyons

Synthetic urban trajectories perturbed with independent white Gaussian jitter look nothing like real downtown GPS, because the dominant error source between tall buildings is multipath — a slowly wandering, position-dependent bias that white noise cannot reproduce.

Part of Noise Injection & Stochastic Drift: that stage degrades clean generated paths into sensor-realistic ones, and this page addresses the specific failure that appears the moment a trajectory enters a dense urban core — a noise model calibrated on open-sky data undershoots the magnitude, the correlation, and the spatial structure of the error a receiver actually produces in a canyon.

Root Cause: Reflected Signals Produce Correlated, Biased Error

Open-sky GNSS error is close to zero-mean and nearly white: each fix is a fresh, roughly independent draw around the true position, so adding independent Gaussian noise per sample is a fair approximation. In an urban canyon the geometry breaks both assumptions. Tall buildings block the direct line of sight to low-elevation satellites and reflect their signals off glass and stone, so the receiver locks onto a path that travelled farther than the straight-line distance. That extra path length enters the range measurement as a positive bias, and because the reflecting surfaces and blocked satellites change slowly as the vehicle moves down a street, the resulting position error is temporally correlated and position-dependent — it drifts smoothly over seconds and tens of metres rather than resampling every fix.

White Gaussian noise $\varepsilon_t \sim \mathcal{N}(0, \sigma^2)$ with $\mathbb{E}[\varepsilon_t \varepsilon_{t-1}] = 0$ therefore gets three things wrong at once: it is zero-mean (multipath is biased), it is uncorrelated (multipath persists across fixes), and it is spatially uniform (multipath is severe on a narrow street and mild in an open plaza). A first-order Gauss-Markov process captures the first two:

$$b_t = \phi\, b_{t-1} + \eta_t, \qquad \eta_t \sim \mathcal{N}(0, \sigma_\eta^2), \qquad 0 < \phi < 1,$$

where the bias state $b_t$ decays toward zero with correlation coefficient $\phi$ (larger $\phi$ means longer-lived drift) and is driven by process noise $\eta_t$. The stationary variance is $\sigma_b^2 = \sigma_\eta^2 / (1 - \phi^2)$, and the correlation time is $\tau = -\Delta t / \ln\phi$. Making the driving variance and correlation depend on where the vehicle is — high on a narrow canyon segment, low in the open — adds the spatial structure. This is a strict generalization of the open-sky model that adding realistic GPS noise to synthetic vehicle trajectories establishes: set $\phi \to 0$ and the spatial term to a constant and the Gauss-Markov process collapses back to white noise.

Minimal Reproducer: Show White Noise Fails the Autocorrelation Test

Confirm the gap by measuring the autocorrelation of a real urban error series and comparing it to what white Gaussian noise produces. Pin the toolchain so the RNG stream is reproducible.

Copy
# requirements.txt — pin majors, let patch releases float
numpy==1.26.*
scipy==1.11.*
pyproj==3.*

pythonCopy
import numpy as np

def lag1_autocorr(err: np.ndarray) -> float:
    """Lag-1 autocorrelation of a scalar error series (e.g. cross-track offset)."""
    e = err - err.mean()
    denom = np.sum(e * e)
    return float(np.sum(e[1:] * e[:-1]) / denom) if denom > 0 else 0.0

rng = np.random.default_rng(7)
white = rng.normal(0.0, 8.0, size=600)          # open-sky-style jitter, sigma 8 m
# Real downtown cross-track error is biased and sticky:
real_urban_lag1 = 0.93                            # measured from a canyon drive log

print(round(lag1_autocorr(white), 3))            # ~ 0.0  -> uncorrelated
print(real_urban_lag1)                            # ~ 0.93 -> strongly correlated
print(round(white.mean(), 2))                     # ~ 0.0  -> zero-mean, but urban error is biased

White noise reports a lag-1 autocorrelation near zero and a mean near zero; a canyon drive log reports autocorrelation above 0.9 and a non-zero mean. Any model whose synthetic error autocorrelation does not approach the measured value will not fool a detector that inspects error dynamics.

Fix: Position-Dependent Gauss-Markov Bias with Occasional Jumps

Model the urban error as the sum of a slowly drifting Gauss-Markov bias whose intensity depends on local canyon severity, plus a small white measurement term, plus rare discrete jumps when the satellite set changes. Work in a metric CRS so the noise is in metres, never in EPSG:4326 degrees.

Step 1 — A canyon-severity field

Drive the noise from a per-position severity in $[0, 1]$: high on narrow streets flanked by tall buildings, low in open space. In production this comes from a building-height raster or sky-view factor; here it is any callable over the path.

pythonCopy
import numpy as np

def severity_along(path_xy: np.ndarray, sky_view: np.ndarray) -> np.ndarray:
    """Map a per-sample sky-view factor (1=open, 0=fully enclosed) to canyon severity."""
    sky_view = np.clip(sky_view, 0.0, 1.0)
    return 1.0 - sky_view                              # 0 open plaza ... 1 deep canyon

Step 2 — Correlated, position-dependent bias process

Integrate a first-order Gauss-Markov bias whose correlation $\phi$ and driving standard deviation scale with severity, so the error is stickier and larger in the canyon and relaxes back to open-sky behaviour in the open.

pythonCopy
import numpy as np

def multipath_bias(severity: np.ndarray, dt: float, *, seed: int,
                   phi_open: float = 0.6, phi_canyon: float = 0.97,
                   sigma_open: float = 1.5, sigma_canyon: float = 12.0
                   ) -> np.ndarray:
    """Return an (N, 2) correlated bias in metres, driven by local canyon severity."""
    rng = np.random.default_rng(seed)
    n = len(severity)
    b = np.zeros((n, 2))
    for t in range(1, n):
        s = severity[t]
        phi = phi_open + (phi_canyon - phi_open) * s        # longer memory in canyon
        sigma_eta = sigma_open + (sigma_canyon - sigma_open) * s
        eta = rng.normal(0.0, sigma_eta, size=2)
        b[t] = phi * b[t - 1] + eta                         # Gauss-Markov update
    # Bias points preferentially away from the reflecting facade: add a mild along-severity
    # offset so deep-canyon segments carry a persistent positive component, not just zero-mean drift.
    b += 0.4 * severity[:, None] * sigma_canyon * np.array([1.0, 0.0])
    return b

Step 3 — Satellite-set jumps and the final observation

Add rare discrete jumps where the visible constellation changes (turning a corner, passing a gap), then combine bias, jumps, and a small white term into the observed position. Keep the timestamps untouched — resampling and clock alignment belong to the temporal synchronization for moving objects stage, not here.

pythonCopy
import numpy as np

def apply_urban_gps(clean_xy: np.ndarray, severity: np.ndarray, dt: float, *,
                    seed: int, jump_rate: float = 0.01, sigma_white: float = 1.0
                    ) -> np.ndarray:
    """clean_xy: (N,2) metric truth from the physics stage -> observed urban positions."""
    rng = np.random.default_rng(seed)
    bias = multipath_bias(severity, dt, seed=seed)
    jumps = np.zeros_like(clean_xy)
    active = np.zeros(2)
    for t in range(len(clean_xy)):
        if rng.random() < jump_rate * (0.5 + severity[t]):   # more re-locks in canyon
            active = rng.normal(0.0, 6.0 + 10.0 * severity[t], size=2)
        jumps[t] = active                                    # step persists until next re-lock
    white = rng.normal(0.0, sigma_white, size=clean_xy.shape)
    return clean_xy + bias + jumps + white

Because the bias is driven off the clean geometry from the physics-based path generation stage, run this perturbation after the kinematic envelope is enforced — the correlated drift will inflate finite-difference derivatives, so a path that sat on its jerk limit before noise will breach it after.

Verification Step: Gate Autocorrelation and Spatial Structure

Assert that synthetic urban error reproduces both the temporal correlation and the position dependence measured from real drives, within tolerance.

pythonCopy
import numpy as np

def test_urban_error_is_correlated_and_position_dependent(
        clean_xy, observed_xy, severity, *, target_lag1=0.90, tol=0.08):
    err = observed_xy - clean_xy
    # 1. Temporal correlation must approach the measured urban value, not ~0.
    e = err[:, 0] - err[:, 0].mean()
    lag1 = float(np.sum(e[1:] * e[:-1]) / np.sum(e * e))
    assert lag1 >= target_lag1 - tol, f"error not sticky enough: lag1={lag1:.2f}"

    # 2. Error magnitude must be larger in canyon segments than open ones.
    mag = np.linalg.norm(err, axis=1)
    hi = mag[severity > 0.7].mean()
    lo = mag[severity < 0.3].mean()
    assert hi > 1.8 * lo, f"no position dependence: canyon {hi:.1f} vs open {lo:.1f}"

    # 3. Sanity: not degenerate / not exploding.
    assert np.isfinite(mag).all() and mag.max() < 150.0, "implausible error magnitude"

The autocorrelation assertion is the load-bearing one: white Gaussian noise cannot pass it at any variance, so a green gate proves the model captures error dynamics and not merely error magnitude. The position-dependence assertion prevents a well-tuned but spatially uniform model from sneaking through.

Edge Cases & Gotchas

Degrees are not metres — project first. Applying a metre-valued bias to EPSG:4326 coordinates adds roughly 8 metres of longitude error at the equator but wildly more near the poles, because a degree of longitude shrinks with latitude. Always transform the clean path into a UTM zone or EPSG:6933 before adding noise, then transform back; a bias generated in degrees produces error that is anisotropic and latitude-dependent for entirely the wrong reason.

Correlation time versus sample rate. The correlation coefficient $\phi$ is tied to the sampling interval through $\tau = -\Delta t / \ln \phi$. A $\phi = 0.97$ calibrated at 1 Hz implies a very different physical memory at 10 Hz, so if you resample the trajectory you must recompute $\phi$ from the target correlation time, not carry the raw coefficient across. Carrying it unchanged makes the drift look implausibly persistent or implausibly brief after a rate change.

Antimeridian and UTM-zone seams. A vehicle crossing a UTM zone boundary or the antimeridian while a bias state is active will have that bias applied in two different projected frames, injecting a spurious jump at the seam. Keep the whole canyon segment in a single local projected CRS whose central meridian sits inside the region, and only stitch frames at points where the bias state has relaxed near zero.

Frequently Asked Questions

Related

• Noise Injection & Stochastic Drift — the perturbation stage this page extends, from open-sky jitter through correlated urban error.
• Adding Realistic GPS Noise to Synthetic Vehicle Trajectories — the open-sky baseline this multipath model generalizes by adding correlation and position dependence.
• Physics-Based Path Generation — the clean-geometry stage whose output this noise model degrades, and whose kinematic limits must be enforced first.
• Temporal Synchronization for Moving Objects — resample and align the noisy fixes to a common clock without touching the spatial error model.