🔧 Fix Linear referencing Data cleaning Linear Referencing · step 20 Notebook By · Updated

Event lengths do not add up to the route length

On this page

Problem statement

The route is 69.322 miles. The line events on it declare 69.322 miles and only 66.958 of those have a matching route segment underneath — or they declare 74 miles, which is worse. Either way something is wrong. Under-coverage means stretches of road with no value; over-coverage means events describing the same stretch twice, and whichever one the code happens to find first wins.

Neither shows up on a map. A drawn event table with a gap looks continuous, and a drawn event table with an overlap looks normal. The only way to see it is arithmetic: span, covered, overlap, and the residual between them.

Quick answer

Three numbers, printed every time an event table is built:

import numpy as np, pandas as pd
from lrs import query, routes

shn = routes("County='SLO' AND Route=101 AND AlignCode='Right'").sort_values("bPM")
spans = [((s.PMPrefix or ""), s.bPM, s.ePM) for s in shn.itertuples()]

aadt = pd.DataFrame([f["attributes"] for f in
                     query("Traffic_AADT", "CNTY='SLO' AND RTE='101'")])
aadt["PM"] = pd.to_numeric(aadt.PM, errors="coerce")
aadt["AHEAD_AADT"] = pd.to_numeric(aadt.AHEAD_AADT, errors="coerce")
aadt["prefix"] = aadt.PM_PFX.fillna("")
st = (aadt.drop_duplicates(subset=["prefix", "PM", "DESCRIPTION"])
      .sort_values("PM").dropna(subset=["PM"]).reset_index(drop=True))

events = [(a.prefix, a.PM, b.PM) for a, b in zip(st.itertuples(), st.shift(-1).itertuples())
          if not pd.isna(b.PM) and b.PM > a.PM]

def covered_measure(prefix, lo, hi):
    """How much of [lo, hi) has a route segment of this prefix under it."""
    return sum(max(0.0, min(hi, b) - max(lo, a)) for p, a, b in spans if p == prefix)

span = st.PM.max() - st.PM.min()
declared = sum(hi - lo for _, lo, hi in events)
covered = sum(covered_measure(*e) for e in events)
print(f"events span      {span:7.3f} miles of postmile")
print(f"event lengths    {declared:7.3f} miles declared")
print(f"route underneath {covered:7.3f} miles carried by a matching segment")
print(f"unaccounted      {declared - covered:7.3f} miles")
events span       69.322 miles of postmile
event lengths     69.322 miles declared
route underneath  66.958 miles carried by a matching segment
unaccounted        2.364 miles

Two and a third miles of declared event length has no route segment under it. That is not a rounding error and it is not visible on any map of this data.

A route with event bars above it showing a gap, an overlap and the residual arithmetic. span, declared, overlap, gapoverlapgapthe routespanfirst start to last enddeclaredsum of event lengthsoverlapcounted twicegapscovered by nothingdeclared − overlap + gaps should equal the span, exactly spatialworkflow.io
Four quantities, and the interesting one is the difference between two of them.

Span, declared length and covered length, then sweep

1. Compute the span, the declared length and the covered length

They are three different things, and any two of them agreeing proves nothing about the third.

2. Find the overlaps by sweeping in order

Sort by start measure, track the furthest end reached, and accumulate whenever an event starts before it.

3. Find the gaps the same way

Accumulate whenever an event starts after the furthest end reached. Gaps inside the span are the ones that matter.

4. Attribute the residual

Prefix mismatches, numbering gaps, events off the end of the route, duplicated source rows. Each has a different fix.

5. Deduplicate the source before building events

A table with every row twice produces zero-length events and doubled coverage, and neither is obvious downstream.

6. Publish the reconciliation with the data

Four numbers in the metadata cost nothing and tell the next reader exactly how complete the layer is.

A single ordered sweep that finds both gaps and overlaps. one pass in measure order1sort the events by start measurenot by row order2track the furthest end reachednot the previous event’s end3starts before it -> overlapcatches containment too4starts after it -> gapinside the span, it matters5reconcile the four numbersthe residual should be zeroComparing each event only with its predecessor misses an event contained in an earlier one spatialworkflow.io
One pass in measure order finds both faults; the trick is tracking the furthest end, not the previous one.

Where the missing miles are

The data used on this page

Caltrans's state highway network and its event layers — traffic counts (Traffic_AADT), bridges (State_Highway_Bridges) and functional classification — from the Caltrans ArcGIS feature services. lrs.py packages the download code from the reading guide as a module: query() pages through any layer, routes() returns centrelines with their postmiles, projected by the server to California Albers (EPSG:3310), and point_at_measure() interpolates a position along one line.

%%writefile lrs.py
"""Caltrans linear referencing: the reading guide's fetch and to_frame, as a module."""
import hashlib
import urllib.parse

import geopandas as gpd
import numpy as np
import requests
from shapely.geometry import LineString

SERVICE = "https://caltrans-gis.dot.ca.gov/arcgis/rest/services/CHhighway/{layer}/FeatureServer/0/query"
# Every page this guide requests, as Caltrans served it on 9 October 2026 (the layers are revised as
# counts are updated): FROZEN holds one file per request, named by a SHA-1 of the live URL.
FROZEN = "https://data.spatialworkflow.io/caltrans-2026-10-09/{key}.json"


def query(layer, where, out_sr=4326, want_m=False):
    """Every feature of a Caltrans layer that matches where, following the service's pages."""
    features, offset = [], 0
    while True:
        params = {"where": where, "outFields": "*", "f": "json", "outSR": out_sr,
                  "resultOffset": offset, "resultRecordCount": 2000}
        if want_m:
            params["returnM"] = "true"                   # without this there is no measure
        url = SERVICE.format(layer=layer) + "?" + urllib.parse.urlencode(params)
        key = hashlib.sha1(url.encode()).hexdigest()[:16]   # the live request: requests.get(url)
        page = requests.get(FROZEN.format(key=key), timeout=180).json()
        features += page["features"]
        if not page.get("exceededTransferLimit"):         # the service says when there is more
            return features
        offset += len(page["features"])


def to_frame(features, crs):
    """One row per path, the measures kept as an array per row."""
    rows, geoms, measures = [], [], []
    for f in features:
        for path in f["geometry"]["paths"]:
            a = np.asarray(path, dtype="float64")
            rows.append(f["attributes"])
            geoms.append(LineString(a[:, :2]))
            measures.append(a[:, 2])
    frame = gpd.GeoDataFrame(rows, geometry=geoms, crs=crs)
    frame["m"] = measures
    return frame


def routes(where, out_sr=3310):
    """State highway centrelines with postmiles, projected by the server (California Albers)."""
    frame = to_frame(query("SHN_Lines", where, out_sr=out_sr, want_m=True), out_sr)
    frame["pm_route"] = [f"{r.County}-{r.Route}{r.PMPrefix or ''}{r.PMSuffix or ''}-{r.AlignCode}"
                         for r in frame.itertuples()]         # the full route key
    return frame


def point_at_measure(line, m, target):
    """The position at a measure along one line, or None if the line does not carry it."""
    xy, m = np.asarray(line.coords), np.asarray(m)
    if not m.min() - 1e-9 <= target <= m.max() + 1e-9:
        return None
    distance = np.concatenate([[0.0], np.cumsum(np.hypot(*np.diff(xy, axis=0).T))])
    order = np.argsort(m)                                # np.interp needs increasing x
    return line.interpolate(float(np.interp(target, m[order], distance[order])))
Writing lrs.py

Example 1 — where the missing miles are

missing = []
for prefix, lo, hi in events:
    gap = (hi - lo) - covered_measure(prefix, lo, hi)
    if gap > 0.001:
        missing.append({"prefix": prefix or "(none)", "from_pm": lo, "to_pm": hi,
                        "missing_mi": round(gap, 3)})
missing = pd.DataFrame(missing)
print(f"{len(missing)} events are not fully carried by the route")
print(missing.to_string(index=False))
print(f"\ntotal {missing.missing_mi.sum():.3f} miles")
7 events are not fully carried by the route
prefix  from_pm  to_pm  missing_mi
(none)    7.851 12.521       0.090
(none)   17.756 19.812       0.274
     R   24.296 25.911       0.094
(none)   30.360 37.863       1.640
(none)   48.331 49.319       0.004
(none)   49.319 50.644       0.006
(none)   67.228 67.712       0.256

total 2.364 miles

The largest is an unprefixed event from postmile 30.360 to 37.863 running across a stretch numbered R between 32.976 and 34.567. The event is real, the road is real, and the two are keyed differently — a conflation problem wearing an arithmetic disguise.

Example 2 — overlaps from a source not deduplicated

raw = aadt.dropna(subset=["PM"]).sort_values("PM").reset_index(drop=True)
raw_events = [(a.prefix, a.PM, b.PM)
              for a, b in zip(raw.itertuples(), raw.shift(-1).itertuples())
              if not pd.isna(b.PM) and b.PM > a.PM]

def sweep(rows):
    ordered = sorted(rows, key=lambda r: r[1])
    overlap = gap = 0.0
    reach = None
    for _, lo, hi in ordered:
        if reach is None:
            reach = hi
            continue
        if lo < reach - 1e-9:
            overlap += min(reach, hi) - lo
        elif lo > reach + 1e-9:
            gap += lo - reach
        reach = max(reach, hi)
    return overlap, gap

dropped = sum(1 for a, b in zip(raw.itertuples(), raw.shift(-1).itertuples())
              if not pd.isna(b.PM) and b.PM <= a.PM)
print(f"raw table: {len(raw)} rows -> {len(raw_events)} events, "
      f"{dropped} pairs silently discarded as zero-length or backwards")
print(f"deduplicated: {len(st)} rows -> {len(events)} events")
for name, rows in [("raw", raw_events), ("deduplicated", events)]:
    overlap, gap = sweep(rows)
    total = sum(hi - lo for _, lo, hi in rows)
    print(f"{name:14} declared {total:7.3f} mi, overlap {overlap:6.3f}, gap {gap:6.3f}")
raw table: 96 rows -> 47 events, 48 pairs silently discarded as zero-length or backwards
deduplicated: 48 rows -> 47 events
raw            declared  69.322 mi, overlap  0.000, gap  0.000
deduplicated   declared  69.322 mi, overlap  0.000, gap  0.000

The two runs agree, and that is the trap: the raw table has each station twice, so 48 of its 95 consecutive pairs are zero-length and are silently discarded by the b.PM > a.PM guard. The totals come out identical to the deduplicated run, so nothing in the output suggests the source had a problem. A table that needs deduplicating should be deduplicated explicitly, where it can be counted, rather than filtered by accident inside a comprehension.

Example 3 — a reusable audit

def audit(rows, name):
    ordered = sorted(rows, key=lambda r: r[1])
    overlap, gap = sweep(ordered)
    lo, hi = ordered[0][1], max(r[2] for r in ordered)
    declared = sum(b - a for _, a, b in ordered)
    return {"table": name, "events": len(ordered), "span": round(hi - lo, 3),
            "declared": round(declared, 3), "overlap": round(overlap, 3),
            "gaps": round(gap, 3), "residual": round(declared - overlap - (hi - lo) + gap, 3)}

fclass = []
for f in query("CRS_Functional_Classification", "RouteID='SHS_101._P'", want_m=True):
    for path in f["geometry"]["paths"]:
        fclass.append(("", min(path[0][2], path[-1][2]), max(path[0][2], path[-1][2])))

report = pd.DataFrame([audit(events, "AADT volumes (postmile)"),
                       audit(fclass, "functional class (odometer)")])
print(report.to_string(index=False))
                      table  events    span  declared  overlap  gaps  residual
    AADT volumes (postmile)      47  69.322    69.322      0.0 0.000       0.0
functional class (odometer)    2285 806.626   806.624      0.0 0.002       0.0

The functional-class table covers its 806-mile span with 2,285 pieces, two thousandths of a mile of gaps and no overlap at all: a table built by a system that enforces the invariant. The volume table is built from point stations by a convention, and its arithmetic shows it.

Track the furthest end, and publish the residual

Why the sweep must track the furthest end

Comparing each event with its predecessor misses an event that is entirely contained in an earlier long one — it starts after the previous event's start and ends before its end, so a pairwise check sees no overlap. Tracking the maximum end reached so far catches containment as well as partial overlap, in one pass.

Why gaps and overlaps need different responses

A gap is missing information: the right answer is usually to leave it null and say so. An overlap is contradictory information: two rows claim the same stretch, and something has to decide which wins — most recent, most specific, highest priority. That decision is domain knowledge, and it should be written down rather than left to row order.

Why the residual is the number to publish

Span, declared length, overlap and gaps reconcile exactly when the table is well formed. When they do not, the size of the discrepancy is the honest summary of the layer's quality — better than a percentage complete, because it is in the same units as everything else in the analysis.

Reconciliation for two event tables on the same network. the same five columns for every tabletableeventsspan / declaredoverlap / gapsAADT volumes (postmile)4769.322 / 69.3220.000 / 0.000functional class (odometer)2,285806.626 / 806.6240.000 / 0.002volumes against the route4769.322 / 66.958 carried2.364 unaccountedboth tables reconcile internally; only the cross-check against the route does notA table can be perfectly consistent and still not fit the network it names spatialworkflow.io
The same five columns for every event table; the last one should be zero.

Zero-length events and direction-specific tables

  • Zero-length events are legal in some registers (a point event stored as a line) and should be counted separately.
  • Floating point: compare with a tolerance the register's precision justifies, not with ==.
  • Direction-specific events double the declared length legitimately; audit each direction separately.
  • Temporal events overlap in measure and not in time; filter to one date before auditing.
  • Ramps and connectors belong to their own routes; events on them do not reconcile against the mainline.
  • Publish the numbers, not a boolean: "2.364 miles unaccounted" is actionable, "failed validation" is not.

FAQ

Why do my event lengths not add up to the route length?

Gaps, overlaps, prefix mismatches or a source table that was never deduplicated.

How do I find overlaps?

Sort by start measure and sweep, tracking the furthest end reached rather than the previous event's end.

Is a gap always an error?

No — it may be a genuine stretch with no value. Leave it null and report the mileage.

What should I do about overlapping events?

Decide which wins by a documented rule, not by row order.

Should the declared lengths equal the span?

Only when the table is gapless and non-overlapping. Publish all four numbers and let a reader see.

How much residual is acceptable?

That is a judgement, and it can only be made if the number is published.