Vector Tile Reference

Converting Lat/Lon to Slippy-Map Tile Numbers

To find the tile covering a point, compute x = floor((lon + 180) / 360 × 2^z) and derive y from the Web Mercator latitude formula y = floor((1 − asinh(tan(lat)) / π) / 2 × 2^z) — both using radians for the trig and both flooring to an integer. That pair (z, x, y) is exactly the XYZ address in a {z}/{x}/{y}.pbf request. This page gives the forward transform, its inverse (tile → bounding box), the TMS flip, and a copy-pasteable Python module that ties them together against a real coordinate.

When to Use This

Reach for this conversion whenever you need to go from a geographic point to a specific tile rather than letting a map client do it implicitly:

  • Which tile covers this point — locating the single tile a POI, incident, or address falls in.
  • Pre-warming a CDN or cache — enumerating the {z}/{x}/{y} set over a bounding box to prime edge caches before launch.
  • Testing a specific tile — pulling one .pbf to inspect its contents while debugging a build, instead of clicking around a live map.
  • Validating a pipeline — asserting that a known landmark lands in the tile you expect after regenerating a tileset.

The concepts behind these formulas — the EPSG:3857 grid, the extent space, and the y-axis conventions — are covered in the parent tile coordinate systems guide.

Specification Detail

Forward: lon/lat → tile x/y (XYZ, north origin)

Four operations turn a coordinate into a tile numberLongitude and latitude are normalised to the zero-to-one range, the latitude passes through the Mercator log-tangent term, both are scaled by two to the zoom, and the result is floored to integers.THE CONVERSIONNormalise x(lon + 180) / 360Mercator yln(tan + sec) / pithe only non-linear stepScalemultiply by 2^zFloorint(), never round()
The floor is the step people forget. Rounding instead of flooring puts points near a tile edge in the neighbouring tile, and the bug only shows on about 1 in 200 coordinates.

At zoom z the grid is n = 2^z tiles per axis:

text
x = floor( (lon + 180) / 360 · n )
y = floor( (1 − asinh(tan(lat_rad)) / π) / 2 · n )

asinh(tan(lat)) is mathematically identical to the more familiar ln(tan(lat) + sec(lat)) Mercator form and to atanh(sin(lat)); all three require lat in radians. Latitude must be clamped to the Mercator limit ±85.05112878° before the trig, or tan blows up toward the poles. The x term needs no trig at all — longitude maps linearly across the grid, which is why only the y axis carries the Mercator nonlinearity.

Two edge cases fall out of the floor: a coordinate at exactly lon = +180° or the southern latitude limit scales to n rather than n − 1, one past the last valid index, so both results are clamped with min(index, n − 1).

Inverse: tile x/y → lon/lat bounding box

The north-west corner of tile (x, y, z):

text
lon = x / n · 360 − 180
lat = degrees( atan( sinh( π · (1 − 2·y / n) ) ) )

Evaluating this at (x, y) gives the tile’s NW corner and at (x+1, y+1) its SE corner — together the full geographic bbox.

TMS flip

text
y_tms = n − 1 − y_xyz        (self-inverse)

XYZ is what URLs and PMTiles use; TMS is what MBTiles rows and some GDAL tooling use.

Production Command Block

A self-contained module with the forward transform, the inverse, a bbox helper, and the TMS flip, plus a worked example and a live-tile check:

python
#!/usr/bin/env python3
"""slippytiles.py — lon/lat <-> slippy-map tile numbers (Web Mercator, EPSG:3857)."""
import math

MAX_LAT = 85.05112878  # Web Mercator latitude limit

def deg2num(lon: float, lat: float, z: int) -> tuple[int, int]:
    """Longitude/latitude (EPSG:4326) -> XYZ tile (x, y) at zoom z."""
    lat = max(min(lat, MAX_LAT), -MAX_LAT)          # clamp before trig
    lat_rad = math.radians(lat)
    n = 2 ** z
    x = int((lon + 180.0) / 360.0 * n)
    y = int((1.0 - math.asinh(math.tan(lat_rad)) / math.pi) / 2.0 * n)
    return min(x, n - 1), min(y, n - 1)             # clamp east/south edges

def num2deg(x: int, y: int, z: int) -> tuple[float, float]:
    """XYZ tile (x, y) -> longitude/latitude of its north-west corner."""
    n = 2 ** z
    lon = x / n * 360.0 - 180.0
    lat = math.degrees(math.atan(math.sinh(math.pi * (1.0 - 2.0 * y / n))))
    return lon, lat

def tile_bbox(x: int, y: int, z: int) -> tuple[float, float, float, float]:
    """XYZ tile -> (west, south, east, north) geographic bbox."""
    west, north = num2deg(x, y, z)
    east, south = num2deg(x + 1, y + 1, z)
    return west, south, east, north

def xyz_to_tms(x: int, y: int, z: int) -> tuple[int, int]:
    """Flip the y-axis between XYZ (north origin) and TMS (south origin)."""
    return x, (2 ** z - 1) - y

if __name__ == "__main__":
    # Worked example: Berlin, Germany at zoom 12
    lon, lat, z = 13.404954, 52.520008, 12
    x, y = deg2num(lon, lat, z)
    print(f"XYZ tile: z={z} x={x} y={y}")            # z=12 x=2200 y=1343
    print(f"TMS  y  : {xyz_to_tms(x, y, z)[1]}")     # 2752
    print(f"bbox    : {tile_bbox(x, y, z)}")
    # -> (13.359375, 52.4827..., 13.4472..., 52.5327...)  — Berlin falls inside

Then confirm the tile actually exists on a running server (Web Mercator tiles are served in XYZ):

bash
# Fetch the exact tile the module computed and confirm a non-empty MVT body
curl -s -o berlin.pbf "https://tiles.example.com/basemap/12/2200/1343.pbf"
ls -l berlin.pbf     # a populated urban tile is typically tens to hundreds of KB

Interaction Effects

The three coordinates that break a naive conversionAn eight by eight zoom 3 grid highlighting the antimeridian column, the polar row above 85.0511 degrees, and a point landing exactly on a tile boundary.EDGE CASESx → easty → southlon = 180 exactly — x lands at 2^z,one past the last valid columnlat above 85.0511 — outside WebMercator, y goes negativeA point on a tile edge belongs to thetile to its east and south
All three are boundary conditions, not maths errors — which is why they survive unit tests written from the middle of the map.

The (z, x, y) this produces is the direct input to a tile fetch — drop x and y straight into a curl of /{z}/{x}/{y}.pbf as shown above to pull one tile for inspection. From there, feed the downloaded .pbf into a decoder to read its features: the decoding MVT tiles to GeoJSON in Python walkthrough turns that byte payload back into inspectable geometry, closing the loop from coordinate to tile to feature.

The choice of z interacts directly with how deep your tileset was built. Requesting a z above the tileset’s --maximum-zoom returns nothing unless overzooming is configured, so pair this conversion with calculating the optimal max zoom for urban datasets when deciding which zoom to test at. And if the backend is MBTiles rather than PMTiles, apply xyz_to_tms() before querying the tiles table — MBTiles stores tile_row in south-origin TMS order.

To enumerate every tile over a region — the core of cache pre-warming — convert the two bbox corners and iterate the inclusive index range between them:

python
def tiles_in_bbox(west, south, east, north, z):
    x0, y0 = deg2num(west, north, z)   # NW corner -> min x, min y (XYZ)
    x1, y1 = deg2num(east, south, z)   # SE corner -> max x, max y
    for x in range(x0, x1 + 1):
        for y in range(y0, y1 + 1):
            yield z, x, y

# Warm a small metro extent at z14 before launch
for z, x, y in tiles_in_bbox(13.09, 52.34, 13.76, 52.68, 14):
    ...  # curl or HEAD each /{z}/{x}/{y}.pbf to prime the CDN edge cache

Note that the NW corner yields the smallest y because XYZ counts rows southward from the north — reversing the corners silently produces an empty range.

Performance Impact

Where the conversion is fast enough and where it is notFour implementations measured: pure Python loop, numpy vectorised, mercantile and a compiled Rust extension, in millions of conversions per second.THROUGHPUTmillion conversions per secondPure Python, per call1.1 M/smercantile.tile()0.9 M/snumpy, vectorised42 M/sCompiled extension86 M/s
For a request handler any of these is free. For rebuilding a 200-million-row index, the difference between the first and second row is the difference between minutes and an afternoon.

The math is trivial — a handful of floating-point operations per point, effectively free next to any network or disk cost. Two directional notes:

  • Bulk conversion: for millions of points (pre-warming a continental bbox, indexing a point dataset), vectorize with NumPy — replace the scalar math calls with numpy.radians, numpy.arctanh(numpy.sin(...)), and numpy.floor, and process whole arrays at once for a 50–100× speedup over a Python loop.
  • Enumeration cost: the number of tiles over a bbox grows as 4^z, so pre-warming a city at z=16 can be tens of thousands of tiles. Enumerate lazily and cap the zoom to what you actually serve rather than materializing the full list.
  • Caching the constant: 2 ** z and math.pi are the only per-call constants; when converting many points at one zoom, hoist n = 2 ** z out of the loop so it is computed once rather than per point.

For reference, the widely used mercantile library implements exactly these formulas; dropping in mercantile.tile(lon, lat, z) gives the same (x, y) as deg2num above and is worth adopting once you need adjacent helpers like tile parents, children, and Web Mercator bounds.

Common Mistakes

1. Feeding degrees where radians are required

math.tan() and math.asinh() expect radians. Passing lat in degrees produces a y that is wrong by a large, silent margin — the point lands in a plausible-looking but incorrect tile.

python
lat_rad = math.radians(lat)          # correct
# y = ... math.tan(lat) ...          # WRONG: lat still in degrees

2. Forgetting to floor (or flooring after the flip)

The tile index is int(...) of the scaled value; skipping the floor yields a fractional “tile” that breaks any URL or lookup. Always floor to get the integer index, and apply the TMS flip only to the already-floored integer, never to the float.

3. XYZ vs TMS confusion, or latitude past the Mercator limit

Querying an MBTiles tiles table with an XYZ y returns the wrong row (or none) — apply xyz_to_tms() first. Separately, a latitude beyond ±85.0511° sends tan(lat_rad) toward infinity and yields a nonsensical y; the MAX_LAT clamp in deg2num prevents both the overflow and an out-of-range tile.

bash
# Symptom of a missed flip: the tile exists in XYZ but the MBTiles query is empty
sqlite3 basemap.mbtiles \
  "SELECT length(tile_data) FROM tiles
   WHERE zoom_level=12 AND tile_column=2200 AND tile_row=2752;"   # TMS row, not 1343

Parent: Tile Coordinate Systems & the Slippy Map Grid

FAQ

Why floor rather than round?

Because a tile number identifies the tile a coordinate falls inside. Rounding sends coordinates in the upper half of a tile to its neighbour, which affects roughly one coordinate in two hundred and is easy to miss in testing.

What happens at exactly 180 degrees longitude?

The formula yields 2^z, one past the last valid column. Clamp to 2^z - 1, or the request is for a tile that cannot exist.

Do I need to reproject before converting?

No. The formula takes longitude and latitude in EPSG:4326 directly and performs the Mercator step itself. Reprojecting first and then applying it double-projects the coordinate.

Is mercantile fast enough for bulk work?

For a request handler, easily. For tens of millions of conversions it is not — a vectorised numpy implementation is roughly forty times faster and worth the twenty lines.