Decoding MVT Tiles to GeoJSON in Python
To turn a .pbf vector tile back into GeoJSON in Python, decode the protobuf with the mapbox_vector_tile package and then apply an affine transform that maps the tile’s local 0..extent grid onto the EPSG:4326 bounding box implied by its (z, x, y) address. The decode is one function call; the coordinate rescale is the part that trips people up, because the tile bytes carry no geographic coordinates of their own.
When to Use This
Reach for a Python decoder when you need to know what actually shipped in a tile rather than what your source data or Tippecanoe logs claim. Two situations dominate:
- Debugging production tiles. A layer is missing, an attribute is absent, or a feature is clipped. Decoding the exact
.pbfthe CDN served answers “is the data in the tile?” definitively, isolating generation bugs from styling bugs. - CI assertions on tile contents. After a nightly Tippecanoe build, a test decodes a handful of representative tiles and asserts that expected layers exist, feature counts are within tolerance, and required attribute keys survived filtering — catching regressions before they reach clients.
For a broad structural walkthrough of the fields you are decoding, read the MVT encoding internals guide first; this page is the runnable Python counterpart.
Specification Detail
mapbox_vector_tile.decode(raw_bytes) returns a dict keyed by layer name. Each layer is a dict with extent, version, and a features list; each feature has geometry (already GeoJSON-shaped, but in tile-local coordinates) and properties (the keys[]/values[] tables already joined for you).
| Element | Shape after decode() |
Notes |
|---|---|---|
| top level | { "roads": {...}, "water": {...} } |
one entry per Layer.name |
layer["extent"] |
int, usually 4096 |
the divisor for rescaling — never hard-code it |
layer["features"] |
list[dict] |
each has geometry, properties, id, type |
feature["geometry"] |
GeoJSON geometry | coordinates are 0..extent, y-axis up by default |
feature["properties"] |
dict |
decoded typed attributes |
The affine transform from tile-local coordinates to lon/lat has two stages. First, the tile’s Web Mercator bounds come from its address: at zoom z there are 2**z tiles per axis, so tile (x, y) spans a known slice of the mercator square. Second, within the tile, a local coordinate cx in 0..extent is a simple linear interpolation across that slice. Because MVT’s local grid has its origin at the tile’s top-left with y increasing downward, but mapbox_vector_tile yields geometry with y increasing upward, you flip y during rescale. The longitude is linear in mercator x; the latitude requires the inverse Web Mercator (atan(sinh(...))) because latitude is non-linear in the projection.
One consequence of the linear-per-tile interpolation is worth internalizing: the latitude edges (north, south) are computed once per tile with the non-linear inverse, and everything inside the tile is then interpolated linearly between them. That approximation is exact enough at the scale of a single tile — a tile spans a tiny fraction of a degree at anything past zoom 8 — that the error is well below MVT’s own quantization grid. Do not be tempted to run the inverse-mercator per vertex; it is both slower and no more accurate than interpolating within the pre-computed tile bounds.
The properties dict also preserves attribute types faithfully. MVT’s values[] table stores typed scalars — string, int, uint, sint, float, double, bool — and mapbox_vector_tile restores them to the matching Python types. So an assertion like assert feat["properties"]["population"] > 0 works directly against an int, and a value that unexpectedly arrives as a str is a strong signal that mixed source types forced a string coercion during generation.
Production Command
The block below fetches one tile over HTTP, handles gzip, decodes it, rescales every coordinate to EPSG:4326, and writes one GeoJSON FeatureCollection per layer. It depends only on httpx and mapbox-vector-tile.
import gzip
import json
import math
from pathlib import Path
import httpx
import mapbox_vector_tile
def tile_lonlat_bounds(z: int, x: int, y: int) -> tuple[float, float, float, float]:
"""Return (west, south, east, north) in EPSG:4326 for an XYZ tile."""
n = 2.0 ** z
west = x / n * 360.0 - 180.0
east = (x + 1) / n * 360.0 - 180.0
north = math.degrees(math.atan(math.sinh(math.pi * (1 - 2 * y / n))))
south = math.degrees(math.atan(math.sinh(math.pi * (1 - 2 * (y + 1) / n))))
return west, south, east, north
def rescale(coords, extent, bounds):
"""Recursively rescale tile-local (0..extent) coords to lon/lat."""
west, south, east, north = bounds
if isinstance(coords[0], (list, tuple)):
return [rescale(c, extent, bounds) for c in coords]
cx, cy = coords
lon = west + (cx / extent) * (east - west)
lat = south + (cy / extent) * (north - south) # y is already up after decode
return [lon, lat]
def fetch_tile(url: str) -> bytes:
resp = httpx.get(url, headers={"Accept-Encoding": "gzip"})
resp.raise_for_status()
body = resp.content
# httpx transparently gunzips Content-Encoding; guard against double-gzip
# (tiles stored pre-gzipped and served without the header) via magic bytes.
if body[:2] == b"\x1f\x8b":
body = gzip.decompress(body)
return body
def tile_to_geojson(url: str, z: int, x: int, y: int) -> dict[str, dict]:
raw = fetch_tile(url)
decoded = mapbox_vector_tile.decode(raw)
bounds = tile_lonlat_bounds(z, x, y)
collections: dict[str, dict] = {}
for layer_name, layer in decoded.items():
extent = layer["extent"]
features = []
for feat in layer["features"]:
geom = feat["geometry"]
geom = {
**geom,
"coordinates": rescale(geom["coordinates"], extent, bounds),
}
features.append({
"type": "Feature",
"geometry": geom,
"properties": feat.get("properties", {}),
"id": feat.get("id"),
})
collections[layer_name] = {
"type": "FeatureCollection",
"features": features,
}
return collections
if __name__ == "__main__":
Z, X, Y = 10, 512, 384
url = f"https://tiles.example.com/v2/{Z}/{X}/{Y}.pbf"
for name, fc in tile_to_geojson(url, Z, X, Y).items():
out = Path(f"{name}_{Z}_{X}_{Y}.geojson")
out.write_text(json.dumps(fc))
print(f"{name}: {len(fc['features'])} features -> {out}")
Run it against any XYZ endpoint. For a (z, x, y) you extracted from a local archive, serve the archive first (pmtiles serve roads.pmtiles --port 8080) and point the URL at http://localhost:8080.
Interaction Effects
Gzip and Content-Encoding. mapbox_vector_tile.decode() expects raw protobuf. If bytes still carry the gzip magic header 1f 8b, decode raises immediately. httpx and requests decompress transparently when the server sets Content-Encoding: gzip, but tiles are frequently stored pre-gzipped on disk and served with Content-Type: application/x-protobuf and no encoding header — so the client never decompresses. The magic-byte guard in fetch_tile() handles both cases. This is the same double-compression hazard that surfaces on the delivery side; see the parent MVT encoding internals troubleshooting notes.
Layer name equals style source-layer. The keys of the decode() result are the exact Layer.name values a MapLibre style must reference as source-layer. A quick print(list(decoded.keys())) is the fastest way to confirm the contract that keeps styles and tiles in sync.
Decoding straight from an archive. When the tiles live in a PMTiles or MBTiles container rather than behind an HTTP endpoint, skip the network entirely. For MBTiles, remember the tiles table is TMS-indexed, so flip y before decoding:
import sqlite3, gzip, mapbox_vector_tile
con = sqlite3.connect("roads.mbtiles")
z, x, y = 10, 512, 384
tms_y = (2 ** z - 1) - y # MBTiles stores rows in TMS order
row = con.execute(
"SELECT tile_data FROM tiles WHERE zoom_level=? AND tile_column=? AND tile_row=?",
(z, x, tms_y),
).fetchone()
blob = row[0]
if blob[:2] == b"\x1f\x8b":
blob = gzip.decompress(blob)
decoded = mapbox_vector_tile.decode(blob)
Feeding these bytes into the same rescale() and tile_lonlat_bounds() helpers — using the XYZ y, not tms_y — produces identical GeoJSON without a running tile server.
Cross-checking with tippecanoe-decode. For a second opinion, decode the same tile with the Tippecanoe tool and diff feature counts:
tippecanoe-decode tile.pbf 10 512 384 | python3 -c "import sys,json; d=json.load(sys.stdin); print(sum(len(l['features']) for l in d['layers']))"
If the Python and tippecanoe-decode counts disagree, you are almost certainly decoding stale or double-compressed bytes.
Performance Impact
Decoding is cheap per tile — a typical 50–150 KB tile decodes in single-digit milliseconds — but it is pure-Python object construction, so cost scales with feature count, not byte count. A dense z14 building tile with 20,000 features can take 50–100 ms and allocate several megabytes of dicts. The rescale recursion adds a modest constant factor per vertex.
The practical guidance is to sample, not sweep. A CI job that decodes every tile in a continental archive will run for hours and dwarf the Tippecanoe build itself. Instead, decode a fixed set of representative tiles — one dense urban tile, one sparse rural tile, one coastline tile per zoom band — which catches the regressions that matter (missing layers, dropped attributes, clipped coastlines) in seconds. When you do need bulk decoding, stream tiles out of the archive with pmtiles or a direct SQLite read rather than fetching each over HTTP, which removes network latency from the loop entirely.
Common Mistakes
Forgetting to gunzip. The error is blunt:
google.protobuf.message.DecodeError: Error parsing message
or a mapbox_vector_tile failure on the first field. It means gzip bytes reached the decoder. Check the first two bytes; if they are 1f 8b, gzip.decompress() before decoding. This is the single most common failure and the reason the guard above exists.
Not rescaling — coordinates stuck in 0..4096. If your output GeoJSON has coordinates like [2048, 3711] instead of [-73.98, 40.71], you skipped the affine transform. The geometry is correct in tile space; it just has not been mapped to lon/lat. Apply rescale() with the tile’s bounds, and read extent from the layer rather than assuming 4096 — a tile generated with a custom --extent will otherwise be scaled wrong.
Y-flip and TMS confusion. Two independent y-axis conventions collide here. First, if the archive is TMS-indexed (y origin at the south, common in older MBTiles), convert to XYZ with xyz_y = 2**z - 1 - tms_y before computing bounds, or every tile lands in the wrong hemisphere. Second, MVT’s internal grid has y increasing downward while mapbox_vector_tile returns y increasing upward — the code above already accounts for the latter, so do not flip a second time. The full coordinate derivation lives in converting lat/lon to slippy map tile numbers, the inverse of the bounds math used here.
Parent: MVT Encoding Internals: The Protobuf Tile Structure
Related
- Converting Lat/Lon to Slippy Map Tile Numbers — the forward
(lon, lat) → (z, x, y)math that inverts the tile-bounds function used here to rescale coordinates. - How to Inspect PMTiles Metadata with CLI Tools — locating and extracting the exact
.pbftile to feed this decoder from a PMTiles archive.