Interpreting Reprojection Error Histograms

Every reconstruction report quotes a mean reprojection error, and it is close to useless on its own. A block with a 0.4 px mean can be excellent, can be a well-converged fit to partly wrong data, or can be a solver that discarded a third of its observations to reach that number. The distribution’s shape distinguishes them, and it costs one plot.

This page reads that shape: what a healthy distribution looks like, the three abnormal shapes and what each means, and why the spatial arrangement of the errors matters at least as much as the histogram.

What the distribution should look like

Reprojection error is the distance in pixels between where an observation was measured and where the solved geometry says it should be. On a sound block those residuals are dominated by measurement noise in the keypoint localisation, which is approximately Gaussian in each image axis — so their magnitude follows a Rayleigh-like shape: zero at zero, a single peak below one pixel, and a tail that decays quickly.

Three deviations from that shape are diagnostic.

A heavy tail — a visible population beyond three or four pixels — means outliers survived into the solution. Under a squared loss those observations dominate the cost, which is the divergence mechanism described in fixing OpenSfM bundle adjustment divergence.

A second mode — a distinct bump at two or three pixels — means two populations of observations are being fitted at once, typically two camera models or two flights whose intrinsics differ.

A distribution that is too tight — almost everything under 0.1 px — means the model has more freedom than the data constrains. That happens when self-calibration is estimating distortion parameters a small block cannot determine, and the price is paid at the block edges where nothing was measured.

Four reprojection-error distributions and what each indicates Four residual distributions on a shared axis of reprojection error in pixels. The healthy shape rises from zero to a single peak below one pixel and decays quickly to nothing by three pixels. The heavy-tailed shape has the same peak plus a long population extending past six pixels, indicating surviving outliers. The bimodal shape has a second distinct bump near two and a half pixels, indicating two populations such as two camera models. The over-fitted shape is a narrow spike below a fifth of a pixel with nothing beyond, indicating more model freedom than the data constrains. 0 1 2 3 5 7 px reprojection error healthy heavy tail — outliers survived bimodal — two populations too tight — over-fitted All four can report a mean between 0.3 and 0.6 pixels. The mean is the one statistic that does not distinguish them.

Figure 1 — Four distributions that a single mean cannot tell apart. Reading the shape takes a glance and rules out three quite different faults.

Minimal reproducible solution

Compute the residuals per observation, then report the statistics that describe shape rather than centre.

import numpy as np


def residual_stats(residuals_px: np.ndarray) -> dict:
    """Shape-describing statistics for a reprojection residual set."""
    r = np.asarray(residuals_px, dtype=float)
    r = r[np.isfinite(r)]
    if r.size < 100:
        raise ValueError("too few observations to describe a distribution")

    p50, p95, p99 = np.percentile(r, [50, 95, 99])
    return {
        "n": int(r.size),
        "median_px": float(p50),
        "p95_px": float(p95),
        "p99_px": float(p99),
        "max_px": float(r.max()),
        # Tail heaviness: on a Rayleigh-like shape this ratio sits near 2.
        "tail_ratio": float(p99 / max(p50, 1e-9)),
        # Fraction beyond a threshold a sound block should almost never reach.
        "beyond_3px": float((r > 3.0).mean()),
    }


def classify_shape(stats: dict) -> str:
    if stats["median_px"] < 0.08:
        return "over-fitted: residuals are below plausible keypoint localisation noise"
    if stats["beyond_3px"] > 0.01:
        return f"heavy tail: {stats['beyond_3px']:.1%} of observations beyond 3 px"
    if stats["tail_ratio"] > 5.0:
        return f"tail ratio {stats['tail_ratio']:.1f} — outliers or a second population"
    return "healthy"

The over-fitting bound deserves justification. Keypoint localisation on real imagery is good to roughly a tenth of a pixel at best; a median residual below that means the model is reproducing measurement noise, which it can only do by having spare parameters. On a small block that is almost always self-calibration estimating distortion terms the geometry does not determine.

tail_ratio is the cheapest single number for the shape. On a Rayleigh distribution the ratio of the 99th percentile to the median is about 2.2; anything above 5 means the tail is not the tail of a single population.

Edge-case matrix

Distribution Median Beyond 3 px Meaning
Healthy 0.3–0.8 px < 0.1% Nothing to do
Heavy tail 0.3–0.8 px 1–5% Outliers survived; use a robust loss
Bimodal 0.5–1.2 px < 1% Two camera models or two flights
Over-fitted < 0.1 px ~0% Too many free parameters
Uniformly high 1.5–3 px 5–20% Wrong intrinsics, or blurred imagery
Tight with a few extremes 0.4 px < 0.1%, max > 20 px A handful of gross mismatches
Bimodal at 0 and 2 px — — Some observations were fixed, not solved
Heavy tail on one flight only varies localised Flight-specific, not block-specific

The last two rows are worth separating out. Observations pinned at exactly zero residual are not being adjusted at all — usually control points weighted so heavily they are effectively fixed, which makes their residual meaningless as a quality measure. And a tail confined to one flight is a data problem in that flight rather than a solver problem in the block.

Verification snippet

The histogram describes the population; the map describes where it lives. A block can have a textbook distribution and concentrate every large residual in one corner, which is a geometry problem the histogram cannot see.

import numpy as np


def assert_residuals_spatially_uniform(residuals_px: np.ndarray,
                                       positions_xy: np.ndarray,
                                       cells: int = 4, ratio: float = 2.0) -> None:
    """No region of the block may carry a much larger median residual.

    Bins observations on a coarse grid over the block and compares each cell's
    median against the block median, which catches a corner the histogram hides.
    """
    x, y = positions_xy[:, 0], positions_xy[:, 1]
    xi = np.clip(((x - x.min()) / (np.ptp(x) + 1e-9) * cells).astype(int), 0, cells - 1)
    yi = np.clip(((y - y.min()) / (np.ptp(y) + 1e-9) * cells).astype(int), 0, cells - 1)

    overall = float(np.median(residuals_px))
    worst_cell, worst_val = None, 0.0
    for i in range(cells):
        for j in range(cells):
            sel = (xi == i) & (yi == j)
            if sel.sum() < 30:
                continue
            m = float(np.median(residuals_px[sel]))
            if m > worst_val:
                worst_cell, worst_val = (i, j), m

    assert worst_val <= overall * ratio, (
        f"cell {worst_cell} has a median residual of {worst_val:.2f} px against "
        f"{overall:.2f} px overall — the error is concentrated, not distributed")

A concentrated residual field has a small number of causes and they are all worth knowing about: a region flown at a different altitude, an area of repetitive texture producing systematically wrong matches, or the edge of the block where camera poses are weakly constrained. All three are visible on the map and invisible in the histogram.

Two blocks with identical histograms and different residual maps Two plan views of a survey block, shaded by median reprojection residual per cell. The first is uniform: every cell sits near the block median, which is what a sound block looks like. The second has the same overall distribution but concentrates its larger residuals in one corner, which indicates a local cause such as a region flown at a different altitude or an area of repetitive texture. A note states that both blocks report the same histogram and the same mean, so only the spatial view distinguishes them. uniform — sound every cell near the block median concentrated — a local cause two cells carry three times the median Identical histograms, identical means, different blocks. A concentrated field points at a region: a different altitude, repetitive texture, or a weakly constrained edge.

Figure 2 — Why the map matters as much as the histogram. The distribution answers “how large are the errors”; only the map answers “where”, and the second question is the one with an actionable answer.

What to record per run

Three numbers and one plot are enough to make this a routine check rather than an investigation, and all four come free from data the solver already produced.

Record the median, the 99th percentile and the fraction beyond three pixels in the run manifest alongside the engine version and the worker count. Those three describe the shape well enough that a comparison between two runs is meaningful, which a mean does not support: two blocks with the same mean can differ by an order of magnitude in their tail.

Write the per-cell residual map as a small raster next to the outputs. It costs a few kilobytes and it is the artefact that answers “was this region always like that?” when a client queries a particular corner of a deliverable months later.

Accumulated across a fleet, those numbers become a baseline. The useful property is not any single run’s value but the distribution of values across runs with the same camera and the same settings — a survey whose tail fraction is five times the fleet norm has something specific wrong with it, and that is visible immediately rather than after somebody thinks to look.

When to escalate

  • The distribution is healthy and the checkpoint accuracy is poor. Reprojection error measures internal consistency, not accuracy. A block can be perfectly self-consistent and systematically displaced; that is what independent checkpoints exist to catch, as described in control points versus checkpoints.
  • Residuals tighten every time self-calibration is given more freedom. That is the over-fitting signature. Fix the intrinsics from a calibration and solve for pose only; if the residuals rise substantially, the extra parameters were absorbing something real and the calibration is wrong.
  • A second mode appears only after adding a flight. The two flights have different effective intrinsics — a different camera, a changed zoom, or a firmware update that altered the reported focal length. Reconstruct them with separate camera models rather than forcing one.

← Optimizing Bundle Adjustment with Python

Residual against free parameters, inside and outside the block A chart of error against the number of free camera parameters allowed to self-calibration. The reprojection residual measured at the observations falls steadily as more parameters are freed, approaching zero. The checkpoint error measured at withheld ground points falls at first, reaches a minimum, and then rises as the extra parameters begin fitting noise. The gap between the two curves beyond the minimum is labelled as over-fitting, and a note states that the residual alone would recommend the rightmost setting, which is the worst one for accuracy. f only f, cx, cy + k1 k2 + k3 p1 p2 all free camera parameters → reprojection residual checkpoint minimum checkpoint error rises Optimising on residual alone selects the rightmost setting, which is the worst one for delivered accuracy.

Figure 3 — The reason a smaller residual is not automatically better. The two curves diverge at the point where extra parameters stop describing the camera and start describing the noise.