How to Encode Coordinates as Geohashes in Python
On this page
Problem statement
A geohash turns a latitude and longitude into a short string, such as gcpsv for the area around Heathrow, by repeatedly halving the world and recording which half the point is in. Points that share a prefix are near each other, so geohashes make coordinates easy to index in a database, group in a groupby, use as cache keys or partition a large table by area, with no GIS library at all. Using them well means knowing how big a cell each length gives, how to find a point's neighbours, and where the trick of "same prefix means nearby" fails.
This page implements geohash encoding and decoding in plain Python, measures the cell size at each precision, uses geohashes to count and group 5,281 real airports from OurAirports, finds a cell's neighbours, and shows the edge effect that makes prefix searches miss close points.
Quick answer
Interleave longitude and latitude bisections into bits and map every five bits to a base-32 character:
B32 = "0123456789bcdefghjkmnpqrstuvwxyz"
def encode(lat, lon, precision=7):
lat_r, lon_r, bits, out, even = [-90.0, 90.0], [-180.0, 180.0], 0, "", True
for i in range(precision * 5):
rng, val = (lon_r, lon) if even else (lat_r, lat)
mid = (rng[0] + rng[1]) / 2
bits = bits * 2 + (val >= mid)
rng[0 if val >= mid else 1] = mid
even = not even
if i % 5 == 4:
out, bits = out + B32[bits], 0
return out
lhr = AP.loc["LHR"]
for p in (3, 5, 7, 9):
print(f"precision {p}: {encode(lhr.latitude_deg, lhr.longitude_deg, p)}")precision 3: gcp
precision 5: gcpsv
precision 7: gcpsv9c
precision 9: gcpsv9cng
Heathrow is gcp at precision 3, a cell covering much of south-east England, and gcpsv9cng at precision 9, a cell a few metres across. Every longer geohash starts with every shorter one: that is the property that makes prefix grouping work.
Decode to a box, cell sizes, group by cell
1. Decode back to a box
Decoding reverses the bisections and gives the cell's bounds; the centre is the usual decoded point, and half the cell is the error:
def decode_box(gh):
lat_r, lon_r, even = [-90.0, 90.0], [-180.0, 180.0], True
for ch in gh:
v = B32.index(ch)
for k in range(4, -1, -1):
rng = lon_r if even else lat_r
mid = (rng[0] + rng[1]) / 2
rng[0 if (v >> k) & 1 else 1] = mid
even = not even
return lat_r, lon_r
(la, lo) = decode_box(encode(lhr.latitude_deg, lhr.longitude_deg, 7))
print(f"cell: lat {la[0]:.5f} to {la[1]:.5f}, lon {lo[0]:.5f} to {lo[1]:.5f}")
print(f"centre {(la[0] + la[1]) / 2:.5f}, {(lo[0] + lo[1]) / 2:.5f}; airport {lhr.latitude_deg:.5f}, {lhr.longitude_deg:.5f}")cell: lat 51.46957 to 51.47095, lon -0.46005 to -0.45868
centre 51.47026, -0.45937; airport 51.47075, -0.45991
At precision 7 Heathrow's cell spans 0.00138° of latitude and 0.00137° of longitude. The decoded centre is 55 m south and 37 m east of the airport's reference point: the coordinate is known only to within the cell.
2. How big is a cell?
Each character adds five bits, alternately splitting longitude and latitude. Measure the cell size in kilometres at the equator and at 51.5°N:
for p in range(1, 10):
la, lo = decode_box(encode(0.0001, 0.0001, p))
h = (la[1] - la[0]) * 111.32
for lat0 in (0, 51.5):
w = (lo[1] - lo[0]) * 111.32 * np.cos(np.radians(lat0))
if lat0 == 0: line = f"precision {p}: {h:9.3f} km tall, {w:9.3f} km wide at the equator"
else: line += f", {w:9.3f} km wide at 51.5°N"
print(line)precision 1: 5009.400 km tall, 5009.400 km wide at the equator, 3118.425 km wide at 51.5°N
precision 2: 626.175 km tall, 1252.350 km wide at the equator, 779.606 km wide at 51.5°N
precision 3: 156.544 km tall, 156.544 km wide at the equator, 97.451 km wide at 51.5°N
precision 4: 19.568 km tall, 39.136 km wide at the equator, 24.363 km wide at 51.5°N
precision 5: 4.892 km tall, 4.892 km wide at the equator, 3.045 km wide at 51.5°N
precision 6: 0.611 km tall, 1.223 km wide at the equator, 0.761 km wide at 51.5°N
precision 7: 0.153 km tall, 0.153 km wide at the equator, 0.095 km wide at 51.5°N
precision 8: 0.019 km tall, 0.038 km wide at the equator, 0.024 km wide at 51.5°N
precision 9: 0.005 km tall, 0.005 km wide at the equator, 0.003 km wide at 51.5°N
Each extra character divides the cell by 32: alternately into 4 × 8 and 8 × 4 pieces. Odd precisions give cells as tall as they are wide at the equator, 4.9 km at precision 5 and 153 m at precision 7; even precisions give cells twice as wide as they are tall. At 51.5°N every cell is 38% narrower than at the equator, because the cells are fixed in degrees.
3. Group points by cell
A geohash column is an ordinary string, so grouping by area is a groupby on a prefix:
A["gh"] = [encode(a, b, 6) for a, b in zip(A.latitude_deg, A.longitude_deg)]
by3 = A.gh.str[:3].value_counts()
print(f"airports: {len(A):,}; distinct cells at precision 3: {len(by3):,}; at 6: {A.gh.nunique():,}")
print("busiest precision-3 cells:", ", ".join(f"{k} ({v})" for k, v in by3.head(5).items()))airports: 5,281; distinct cells at precision 3: 3,181; at 6: 5,279
busiest precision-3 cells: gcp (12), c28 (11), 9q5 (11), u15 (10), xn7 (10)
At precision 3, about 156 km cells, the 5,281 airports fall into 3,181 cells; at precision 6, about 0.6 by 1.2 km, almost every airport has its own cell. The busiest precision-3 cell is gcp itself, London and south-east England, with 12 large and medium airports; 9q5, around Los Angeles, has 11.
Neighbours and the edge problem
The data used on this page
Airport coordinates from OurAirports, a public-domain database: the large and medium airports. The geohash functions are written out in full so they run without a geohash library; pygeohash and python-geohash implement the same algorithm.
from pathlib import Path
import numpy as np
import pandas as pd
import requests
Path("data").mkdir(exist_ok=True)
p = Path("data/airports.csv")
if not p.exists():
r = requests.get("https://davidmegginson.github.io/ourairports-data/airports.csv", timeout=300)
r.raise_for_status()
p.write_bytes(r.content)
A = pd.read_csv(p, keep_default_na=False, na_values=[""])
A = A[A.type.isin(["large_airport", "medium_airport"])].reset_index(drop=True)
AP = A.dropna(subset=["iata_code"]).set_index("iata_code")
print(f"{len(A):,} large and medium airports, {len(AP):,} with an IATA code")5,281 large and medium airports, 4,568 with an IATA code
The data has 5,281 large and medium airports, 4,568 of them with an IATA code.
Example 1: neighbours and the edge problem
Two points can be a few kilometres apart yet share no prefix if a cell boundary runs between them. Find each airport's nearest other airport and check how often the pair shares a precision-4 cell:
from scipy.spatial import cKDTree
xyz = np.c_[np.cos(np.radians(A.latitude_deg)) * np.cos(np.radians(A.longitude_deg)),
np.cos(np.radians(A.latitude_deg)) * np.sin(np.radians(A.longitude_deg)), np.sin(np.radians(A.latitude_deg))]
d, j = cKDTree(xyz).query(xyz, k=2)
km = 2 * np.arcsin(d[:, 1] / 2) * 6371.0088
close = km < 20
same4 = A.gh.str[:4].values == A.gh.str[:4].values[j[:, 1]]
print(f"airports whose nearest neighbour is under 20 km away: {close.sum()}")
print(f"of those, pairs in the same precision-4 cell (about 39 x 20 km): {same4[close].mean():.0%}; in different cells: {(~same4[close]).sum()}")airports whose nearest neighbour is under 20 km away: 722
of those, pairs in the same precision-4 cell (about 39 x 20 km): 40%; in different cells: 436
722 airports have another airport within 20 km. Only 40% of those pairs share a precision-4 cell, although the cells are 20 by up to 39 km: in 436 cases a cell boundary runs between two airports that are close. A prefix search for nearby points has to look at the neighbouring cells too.
A geohash is a Z-order curve
A geohash is a Z-order curve: by interleaving the bits of longitude and latitude, it maps two dimensions onto one line so that most points close in space are close in the sorted order. Each character carries five bits, so precision 5 has 25 bits, 13 for longitude and 12 for latitude, and cells are about 4.9 km by 4.9 km at the equator. Odd precisions give roughly square cells, even precisions cells twice as wide as they are tall, and all cells narrow towards the poles because they are fixed in degrees.
The weakness is the boundaries. Neighbouring cells can have completely different strings when they sit either side of a major split, for example either side of the equator or the prime meridian, so a search that only looks for one prefix misses points just across the line. The standard fix is to search the cell and its eight neighbours. For analysis rather than indexing, hexagonal systems such as H3 have cells of nearly equal area and uniform neighbours, which geohashes do not.
Precision, privacy and unequal areas
- Precision and privacy: a 7-character geohash pins a point to about 150 m; truncate to a shorter prefix before sharing locations.
- Letters a, i, l, o: the base-32 alphabet leaves them out to avoid confusion with digits.
- Database indexes: a B-tree index on a geohash column supports prefix queries with
LIKE 'gcps%'. - Not equal area: a precision-5 cell covers about 24 km² at the equator but about 15 km² at 51.5°N.
Internal links
- How to find points within a radius in Python — exact radius search
- How to find the nearest coordinates with a BallTree in Python — nearest neighbours without cell edges
- How many decimal places do coordinates need? — precision in another form
FAQ
What is a geohash?
A string encoding of a latitude and longitude made by repeatedly halving the longitude and latitude ranges and recording each half as a bit; every five bits become one base-32 character.
How accurate is a geohash?
It depends on its length: 5 characters is a cell about 4.9 km across, 7 characters about 150 m, 9 characters about 5 m.
Do points with the same geohash prefix lie close together?
Yes, but the reverse is not guaranteed: close points either side of a cell boundary can have different prefixes.
How do I decode a geohash in Python?
Reverse the bisections, as in step 1, or use the pygeohash package's decode.