How to Map-Match a GPS Track to a Street Network

Problem statement

Map matching answers "which streets did this track use", and the obvious implementation — snap every point to its nearest edge — produces a sequence of streets nobody could have driven.

Measured on 46 real GPS tracks from Brighton, every fix snapped to its nearest edge of OpenStreetMap's drive network for the same area (Example 3):

snap distance p50   8.3 m, p90  80.5 m, max 1276.2 m
11766 edge changes, 459 to an edge sharing no node, 401 of them more than 50 m of road away

Across 55,316 fixes the "street" changed 11,766 times, and 401 of those changes jumped to a road more than 50 m of network away. Matched as a sequence instead, the five longest tracks went from 3,053 street changes to 492 (Example 4): most of what the nearest-edge snap reports as movement between streets is flicker, because each point was matched on its own.

Quick answer

Match the sequence, not the points. Score candidate edges by both distance and connectivity:

import numpy as np
import networkx as nx
import osmnx as ox


def match_sequence(G, points, radius=50, max_candidates=5):
    """Viterbi-style match: emission by distance, transition by route length."""
    candidates = [nearby_edges(G, p, radius, max_candidates) for p in points]

    scores = {c: emission(points[0], c) for c in candidates[0]}
    back = [{}]
    for i in range(1, len(points)):
        new_scores, pointers = {}, {}
        for c in candidates[i]:
            best, best_prev = -np.inf, None
            for prev, prev_score in scores.items():
                s = prev_score + transition(G, prev, c, points[i - 1], points[i]) \
                    + emission(points[i], c)
                if s > best:
                    best, best_prev = s, prev
            new_scores[c], pointers[c] = best, best_prev
        scores, back = new_scores, back + [pointers]

    path = [max(scores, key=scores.get)]
    for pointers in reversed(back[1:]):
        path.append(pointers[path[-1]])
    return list(reversed(path))

The transition term is what removes the teleports: a candidate edge is penalised when the route from the previous edge to it is far longer than the distance the vehicle could have travelled.

46 tracks, 55,316 points; drive network 8,720 edges, walk 35,294

The right network, candidates, two scores, Viterbi

1. Use the right network for the mode

The p90 snap distance above is 81 m, which is far more than GPS error. The reason is that many of these tracks are on foot — the median speed was 4.4 km/h — and they were being matched against a drive network that has no footpaths. On the walk network the same five tracks snapped at a median of 3.3 m and a p90 of 8.0 m, against 6.4 m and 20.4 m on the drive network (Example 4).

G = ox.graph_from_bbox(bbox, network_type="walk")   # not "drive"

Choosing the wrong network type produces large, systematic snap distances that look like GPS error and are not. Check the distribution: a p50 of 8 m is plausible for GPS; a p90 of 81 m and a maximum over a kilometre are not.

2. Project everything before measuring

Snapping compares distances. In degrees those distances are wrong, and wrong differently in each direction.

G = ox.project_graph(G, to_crs="EPSG:32630")

3. Generate candidates, not a single match

For each fix, take every edge within a radius — 30–50 m for good urban data, more where accuracy is poor. Keeping five candidates per point is usually plenty.

A single nearest edge throws away the information that resolves ambiguity: on a dual carriageway the two carriageways are metres apart, and only the sequence says which one you were on.

4. Score each candidate two ways

• Emission: how well does this edge explain this observation? Usually a Gaussian on the perpendicular distance, with a standard deviation set to the GPS error.
• Transition: how plausible is moving from the previous candidate to this one? Compare the shortest-path distance along the network with the straight-line distance between the two fixes. A ratio near 1 is plausible; a ratio of 10 means the route required a detour the timing does not allow.

5. Solve the sequence with Viterbi

The two scores define a hidden Markov model where the hidden state is "which edge am I on". The Viterbi algorithm finds the most likely sequence in time linear in the number of points and quadratic in the number of candidates.

6. Check the output for the failure you were fixing

near = lambda a: nx.single_source_dijkstra_path_length(G, a[1], cutoff=50, weight="length")
jumps = sum(1 for a, b in zip(path, path[1:])
            if a != b and b[0] not in near(a))
print(f"{jumps} steps to an edge more than 50 m of road away")

A good match has few, and those only across genuine data gaps. Test road distance, not shared nodes: at a fix a second, a car crosses city-centre edges shorter than one step.

Candidates, emission and transition on Brighton traces

The data used on this page

Public GPS traces from the OpenStreetMap trackpoints API for a box over central Brighton and Hove, −0.20 to −0.10° E, 50.81 to 50.86° N: the first 12 pages of 5,000 points. The setup loads them the way How to load GPS tracks in Python does — one track per <trkseg>, exact duplicates dropped, elements with more than 10% repeated timestamps set aside, tracks under ten points dropped — and projects them to UTM zone 30N, adding x, y and the interval dt in seconds. It also downloads OpenStreetMap's drive and walk networks for the same box with OSMnx, projected the same way, as G_drive and G_walk. The API serves the newest uploads first, so the setup reads the pages as saved on 9 October 2026. Nothing on this page prints a location or a username.

import glob
import io
import os
import xml.etree.ElementTree as ET
import zipfile

import geopandas as gpd
import numpy as np
import pandas as pd
import requests

# OpenStreetMap as Overpass returned it on 9 October 2026, saved as osmnx's own cache so the counts
# below do not drift as OSM is edited; delete the cache folder to query live instead
if not os.path.exists("cache"):
    reply = requests.get("https://data.spatialworkflow.io/osm-cache-map-match-gps-track-python-2026-10-09.zip", timeout=300)
    reply.raise_for_status()
    with zipfile.ZipFile(io.BytesIO(reply.content)) as archive:
        archive.extractall("cache")

# The OSM trackpoints API for this box, pages 0-11, as saved on 9 October 2026. The API serves the
# newest uploads first, so a live download changes whenever someone uploads a trace here. Live:
#   requests.get("https://api.openstreetmap.org/api/0.6/trackpoints",
#                params={"bbox": "-0.20,50.81,-0.10,50.86", "page": page})
if not os.path.exists("traces/page_11.gpx"):
    reply = requests.get("https://data.spatialworkflow.io/osm-gps-traces-brighton-2026-10-09-v2.zip",
                         timeout=120)
    reply.raise_for_status()
    with zipfile.ZipFile(io.BytesIO(reply.content)) as archive:
        archive.extractall("traces")

NS = {"g": "http://www.topografix.com/GPX/1/0"}
rows = []
for n, path in enumerate(sorted(glob.glob("traces/*.gpx"))):
    for t, trk in enumerate(ET.parse(path).getroot().findall("g:trk", NS)):
        for s, seg in enumerate(trk.findall("g:trkseg", NS)):
            rows += [(f"p{n}:{t}-{s}", float(p.get("lat")), float(p.get("lon")),
                      p.findtext("g:time", namespaces=NS)) for p in seg.findall("g:trkpt", NS)]
raw = pd.DataFrame(rows, columns=["track_id", "lat", "lon", "time"]).dropna()
raw["time"] = pd.to_datetime(raw["time"], utc=True, format="ISO8601")
raw = raw.drop_duplicates(["track_id", "time", "lat", "lon"]).sort_values(["track_id", "time"], kind="stable")
zero = raw.groupby("track_id")["time"].apply(lambda s: (s.diff().dt.total_seconds() == 0).mean())
size = raw.groupby("track_id").size()
keep = zero.index[(zero <= 0.1) & (size >= 10)]

df = gpd.GeoDataFrame(raw[raw.track_id.isin(keep)], crs=4326,
                      geometry=gpd.points_from_xy(raw.lon, raw.lat)[raw.track_id.isin(keep).to_numpy()])
df = df.to_crs(32630).reset_index(drop=True)
df["x"], df["y"] = df.geometry.x, df.geometry.y
df["dt"] = df.groupby("track_id")["time"].diff().dt.total_seconds()
import osmnx as ox

ox.settings.overpass_url = "https://maps.mail.ru/osm/tools/overpass/api"
ox.settings.requests_timeout = 600
ox.settings.max_query_area_size = 50e9
graphs = {}
for kind in ("drive", "walk"):
    if not os.path.exists(f"brighton_{kind}.graphml"):
        g = ox.graph_from_bbox((-0.20, 50.81, -0.10, 50.86), network_type=kind)
        ox.save_graphml(ox.project_graph(g, to_crs="EPSG:32630"), f"brighton_{kind}.graphml")
    graphs[kind] = ox.load_graphml(f"brighton_{kind}.graphml")
G_drive, G_walk = graphs["drive"], graphs["walk"]
print(f"{df.track_id.nunique()} tracks, {len(df):,} points; "
      f"drive network {G_drive.number_of_edges():,} edges, walk {G_walk.number_of_edges():,}")

46 tracks, 55,316 points; drive network 8,720 edges, walk 35,294

Example 1 — candidates and emission scores

import numpy as np
import osmnx as ox
from shapely.geometry import Point, LineString


def edge_geometry(G, u, v, k):
    edge = G.edges[u, v, k]
    if "geometry" in edge:
        return edge["geometry"]
    return LineString([(G.nodes[u]["x"], G.nodes[u]["y"]),
                       (G.nodes[v]["x"], G.nodes[v]["y"])])


def candidates(G, edges_gdf, sindex, x, y, radius=50.0, limit=5):
    """Every edge within radius of the fix, with its perpendicular distance."""
    point = Point(x, y)
    idx = list(sindex.intersection(point.buffer(radius).bounds))
    found = []
    for i in idx:
        row = edges_gdf.iloc[i]
        distance = point.distance(row.geometry)
        if distance <= radius:
            found.append((row.name, float(distance)))
    found.sort(key=lambda t: t[1])
    return found[:limit]


def emission(distance, sigma=10.0):
    """Log-likelihood that this edge produced this observation."""
    return -0.5 * (distance / sigma) ** 2


edges_drive = ox.graph_to_gdfs(G_drive, nodes=False)
sindex_drive = edges_drive.sindex
for i, (x0, y0) in enumerate(zip(df.x, df.y)):
    found = candidates(G_drive, edges_drive, sindex_drive, x0, y0)
    if found:                                    # a fix far from any drivable road has none
        break
print(f"{len(edges_drive):,} drive edges; fix {i} is the first with candidates: {len(found)} within 50 m")
print([((int(u), int(v)), round(distance, 1), round(emission(distance), 2)) for (u, v, _), distance in found])

8,720 drive edges; fix 5 is the first with candidates: 3 within 50 m
[((25345924, 14594929), 49.3, -12.16), ((8105829302, 14594929), 49.3, -12.16), ((14594929, 8105829302), 49.3, -12.16)]

All three candidates score the same: distance alone cannot choose between the two directions of one street, which is what the transition score is for.

Setting sigma to the actual GPS error matters. Too small and the match becomes a nearest-edge snap with extra steps; too large and distance stops discriminating at all.

Example 2 — the transition score

import networkx as nx


def transition(G, prev_edge, next_edge, prev_point, next_point,
               beta=10.0, cache=None):
    """Penalise routes much longer than the straight-line step."""
    straight = ((next_point[0] - prev_point[0]) ** 2 +
                (next_point[1] - prev_point[1]) ** 2) ** 0.5
    if next_edge == prev_edge:
        return 0.0                         # still on the same street: no detour at all

    key = (prev_edge[1], next_edge[0])
    if cache is not None and key in cache:
        route = cache[key]
    else:
        try:
            route = nx.shortest_path_length(G, key[0], key[1], weight="length")
        except nx.NetworkXNoPath:
            route = float("inf")
        if cache is not None:
            cache[key] = route

    if not np.isfinite(route):
        return -1e6                        # unreachable: effectively forbidden
    return -abs(route - straight) / beta

Caching the shortest path between node pairs is what makes this tractable. With five candidates per fix, consecutive steps ask for the same node pairs repeatedly.

The -1e6 for unreachable pairs means an impossible transition can never win, whatever its emission score. The next_edge == prev_edge line matters as much: without it, staying on the same street is costed as a route from the edge's end back to its start, and the matcher hops off every street it is on. The first version of this page had exactly that bug, and it roughly tripled the number of street changes instead of cutting them.

Example 3 — a diagnostic before and after

import networkx as nx
import numpy as np


def match_quality(G, points, matched, max_gap=50.0):
    """Snap distance and connectivity, before trusting the result."""
    distances = []
    for (x, y), edge in zip(points, matched):
        geom = edge_geometry(G, *edge)
        distances.append(Point(x, y).distance(geom))
    distances = np.array(distances)

    steps = [(a, b) for a, b in zip(matched, matched[1:]) if a != b]
    unjoined = [(a, b) for a, b in steps if not set(a[:2]) & set(b[:2])]
    # sharing no node is not yet a jump: a car covers 15 m a second, and city-centre edges are
    # often shorter, so check how much road actually separates the two edges
    jumps = sum(1 for a, b in unjoined
                if b[0] not in nx.single_source_dijkstra_path_length(G, a[1], cutoff=max_gap, weight="length"))

    print(f"  snap distance p50 {np.median(distances):5.1f} m, "
          f"p90 {np.percentile(distances, 90):5.1f} m, "
          f"max {distances.max():5.1f} m")
    print(f"  {len(steps)} edge changes, {len(unjoined)} to an edge sharing no node, "
          f"{jumps} of them more than {max_gap:.0f} m of road away")
    if np.percentile(distances, 90) > 40:
        print("  ! large snap distances — check the network type matches the mode")
    return {"p50": float(np.median(distances)),
            "p90": float(np.percentile(distances, 90)),
            "jumps": jumps, "unjoined": len(unjoined), "changes": len(steps)}


# before: every fix snapped on its own to the nearest drive edge
pts = list(zip(df.x, df.y))
nearest = ox.distance.nearest_edges(G_drive, df.x.to_numpy(), df.y.to_numpy())
before = match_quality(G_drive, pts, nearest)

  snap distance p50   8.3 m, p90  80.5 m, max 1276.2 m
  11766 edge changes, 459 to an edge sharing no node, 401 of them more than 50 m of road away
  ! large snap distances — check the network type matches the mode

That is the before. The warning is the useful part: it points at the network type rather than the algorithm, which Example 4 confirms.

Example 4 — the sequence match, before and after

import time


def match_track(G, edges_gdf, sindex, points, radius=50.0, sigma=10.0):
    """Viterbi over candidate edges; points with no candidate are skipped."""
    cache, kept, cands = {}, [], []
    for p in points:
        c = candidates(G, edges_gdf, sindex, *p, radius=radius)
        if c:
            kept.append(p)
            cands.append(c)
    scores = {e: emission(d, sigma) for e, d in cands[0]}
    back = []
    for i in range(1, len(kept)):
        new, ptr = {}, {}
        for e, d in cands[i]:
            best, arg = max((sc + transition(G, pe, e, kept[i - 1], kept[i], cache=cache), pe)
                            for pe, sc in scores.items())
            new[e], ptr[e] = best + emission(d, sigma), arg
        scores = new
        back.append(ptr)
    path = [max(scores, key=scores.get)]
    for ptr in reversed(back):
        path.append(ptr[path[-1]])
    return kept, path[::-1]


# the five longest tracks, every fix: at 1 s consecutive fixes cannot skip a street
longest = df.groupby("track_id").size().nlargest(5).index
sample = df[df.track_id.isin(longest)]
for name, G in (("drive", G_drive), ("walk", G_walk)):
    edges = ox.graph_to_gdfs(G, nodes=False)
    started = time.perf_counter()
    kept_all, path_all, near_all = [], [], []
    for _, sub in sample.groupby("track_id"):
        kept, path = match_track(G, edges, edges.sindex, list(zip(sub.x, sub.y)))
        kept_all += kept
        path_all += path
        near_all += list(ox.distance.nearest_edges(G, [p[0] for p in kept], [p[1] for p in kept]))
    print(f"{name} network, {len(kept_all):,} of {len(sample):,} fixes had a candidate within 50 m "
          f"({time.perf_counter() - started:.0f} s)")
    print("  nearest edge:")
    match_quality(G, kept_all, near_all)
    print("  sequence match:")
    match_quality(G, kept_all, path_all)

drive network, 15,872 of 19,371 fixes had a candidate within 50 m (8 s)
  nearest edge:
  snap distance p50   6.4 m, p90  20.4 m, max  50.0 m
  3053 edge changes, 57 to an edge sharing no node, 50 of them more than 50 m of road away
  sequence match:
  snap distance p50   6.6 m, p90  20.7 m, max  50.0 m
  492 edge changes, 100 to an edge sharing no node, 28 of them more than 50 m of road away
walk network, 19,209 of 19,371 fixes had a candidate within 50 m (24 s)
  nearest edge:
  snap distance p50   3.3 m, p90   8.0 m, max  49.7 m
  6028 edge changes, 325 to an edge sharing no node, 243 of them more than 50 m of road away
  sequence match:
  snap distance p50   3.5 m, p90   8.3 m, max  49.7 m
  1073 edge changes, 174 to an edge sharing no node, 51 of them more than 50 m of road away

On both networks the sequence match cuts street changes by more than four-fifths, and the jumps that matter — to an edge more than 50 m of road away — from 50 to 28 on the drive network and from 243 to 51 on the walk network. Counting only edges that share no node would say the opposite on the drive network, 57 against 100: between two fixes a second apart the matched route often crosses a short third edge, which is travel, not a jump. It does not remove them all: fixes with no edge within 50 m are skipped, and the step across such a hole can land on a street the route only reaches further on. Snap distances barely move — the matcher trades a slightly worse fit at single points for a sequence that could have been travelled.

Why nearest-edge snapping fails

Why nearest-edge snapping fails

Each point is matched independently, so nothing prevents consecutive points landing on unrelated roads. The measurement is unambiguous: across 55,316 fixes the nearest edge changed 11,766 times, and 401 of those changes went to an edge more than 50 m of road from the previous one.

In dense networks the nearest edge flips between parallel streets, between a road and a service road behind it, and between the two carriageways of a dual carriageway. Each flip is a plausible local decision and the sequence is nonsense.

Why the transition term is the whole algorithm

The insight behind every modern matcher is that movement is constrained. Between two fixes three seconds apart, a vehicle can travel maybe fifty metres along the network — so any candidate requiring a five-hundred-metre detour is wrong regardless of how close it is.

Comparing network distance with straight-line distance captures this in one number. It also handles the U-turn case correctly: doubling back is possible but requires travelling much further than the straight line, so it is penalised rather than forbidden.

Why the network type matters more than the algorithm

The measured p90 snap distance of 81 m on the drive network is not GPS error — on the walk network the five longest tracks snapped at a p90 of 8.0 m. It is the distance from a footpath to the nearest road.

Matching pedestrians to a drive network gives a confident, smooth, entirely fictional route along the roads beside the paths they used. No matcher can detect this, because the input contains nothing that says the pedestrian was not in the road.

Choose network_type from the mode. Where the mode is unknown, match against all and infer the mode from which classes of edge the match used.

Why matching cannot always be right

Map matching is inference, not measurement. Where the network is dense relative to the position error and the sampling is coarse, several routes explain the observations equally well.

Good matchers return a confidence — the margin between the best path and the next best. A low margin means the answer is a guess, and downstream analysis should be able to see that rather than treating every matched route as fact.

Network type, emission sigma and one-way streets

• Choose the network type from the mode. A p90 snap distance above about 40 m usually means the wrong network.
• Project before snapping. Degrees are not metres.
• Split at data gaps. Matching across a fifty-minute gap invents a route.
• Cache shortest paths between node pairs, or the transition term dominates the runtime.
• Keep three to five candidates per fix. More rarely changes the answer and is quadratically more expensive.
• Set the emission sigma to the real GPS error, not to a round number.
• One-way streets matter. Match on the directed graph, or the route will include illegal moves.
• Report a confidence, not just a route.

Internal links

• GPS error explained: why your track wanders — why the emission sigma is what it is
• How to snap points to a street network in Python — the single-point case
• Street networks as graphs explained — the graph being matched against
• How to build a street network graph in Python with OSMnx — building the network
• How to clean a GPS track: outliers, duplicates and jumps — cleaning before matching
• How to calculate the shortest path along a street network — the routine the transition term uses
• Stops, trips and segmentation explained — segment before matching
• How to aggregate movement into flows between zones — what matched routes feed into

FAQ

How do I map-match a GPS track in Python?

Generate several candidate edges per fix, score each by perpendicular distance and by whether the network route from the previous candidate matches the distance travelled, and solve the sequence with Viterbi.

Why does nearest-edge snapping give a broken route?

Because each point is matched independently. Measured on 46 real tracks, the nearest edge changed 11,766 times across 55,316 fixes, and 401 of those changes jumped to a road more than 50 m of network away.

What snap distance should I expect?

A few metres at the median for good GPS on the right network — 3.3 m here on the walk network. A p90 of 81 m, as measured on the drive network, means the network type does not match the travel mode.

Which network type should I use?

The one matching the mode: walk for pedestrians, drive for vehicles, bike for cycling. Matching pedestrians to a drive network produces a confident, fictional route.

How do I handle data gaps?

Split the track and match each leg separately. Matching across a gap forces the algorithm to invent a route.

How do I know if the match is right?

Check that consecutive matched edges share a node, and report the margin between the best path and the second best. A small margin means the answer is a guess.

Do I need a library, or is this enough?

The sketch here is a working Viterbi matcher. Established libraries add one-way handling, better candidate pruning and calibrated parameters — worth using in production, but the mechanism is exactly this.