How to Convert Latitude and Longitude to UTM in Python

Problem statement

Universal Transverse Mercator coordinates are what field surveys, military grids, many national datasets and every "metres, not degrees" calculation use: an easting and a northing in metres within one of 60 zones, each 6° of longitude wide. Converting from latitude and longitude is one call with pyproj, but you have to choose the zone, handle the southern hemisphere's false northing, and decide what to do with data that spans two zones, where coordinates from different zones cannot be mixed.

This page converts the USGS Geographic Names Information System's features for Hawaii to UTM with pyproj, picks zones automatically, checks a known point both ways, shows what happens when points outside a zone are forced into it, and handles southern-hemisphere points with real airports from OurAirports.

Quick answer

Find the zone from the longitude, build the EPSG code, and transform with always_xy=True:

from pyproj import Transformer

def utm_epsg(lat, lon):
    zone = int((lon + 180) // 6) + 1
    return (32600 if lat >= 0 else 32700) + zone

lat, lon = 19.8206, -155.4681                          # Mauna Kea summit area
epsg = utm_epsg(lat, lon)
tr = Transformer.from_crs(4326, epsg, always_xy=True)
e, n = tr.transform(lon, lat)
print(f"EPSG:{epsg} (zone {epsg % 100}{'N' if epsg < 32700 else 'S'}): easting {e:,.1f} m, northing {n:,.1f} m")

EPSG:32605 (zone 5N): easting 241,470.6 m, northing 2,193,518.0 m

Mauna Kea's summit area, at 155.47°W, is in zone 5N, EPSG:32605: 241,470.6 m east and 2,193,518.0 m north. The easting is less than 500,000 m because the point is west of zone 5's central meridian at 153°W.

Zones across a dataset, round trips, forced and wrong zones

1. Zones across a dataset

Hawaii straddles the boundary between zones 4 and 5 at 156°W. Assign each feature its own zone and count:

G["zone"] = ((G.lon + 180) // 6).astype(int) + 1
print(G.zone.value_counts().sort_index().to_dict())
print(f"zone 4 spans {-180 + 3 * 6}° to {-180 + 4 * 6}°; zone 5 spans {-180 + 4 * 6}° to {-180 + 5 * 6}°")
gs = gpd.GeoSeries(gpd.points_from_xy(G.lon, G.lat), crs=4326)
print(f"estimate_utm_crs() for the whole set: {gs.estimate_utm_crs().name}")

{4: 4746, 5: 2557}
zone 4 spans -162° to -156°; zone 5 spans -156° to -150°
estimate_utm_crs() for the whole set: WGS 84 / UTM zone 4N

4,746 features fall in zone 4, which covers Kauaʻi, Oʻahu and Maui, and 2,557 in zone 5, most of the Island of Hawaiʻi. GeoPandas' estimate_utm_crs() picks a single zone for the whole set, zone 4N, the zone of the data's centre.

2. Convert and convert back

Transform every feature in its own zone and back again, and check the round trip:

out = []
for z, part in G.groupby("zone"):
    t = Transformer.from_crs(4326, 32600 + z, always_xy=True)
    e, n = t.transform(part.lon.values, part.lat.values)
    lo, la = Transformer.from_crs(32600 + z, 4326, always_xy=True).transform(e, n)
    out.append(pd.DataFrame({"zone": z, "e": e, "n": n, "err_m": np.hypot((lo - part.lon.values) * 111_000 * np.cos(np.radians(la)), (la - part.lat.values) * 111_000)}, index=part.index))
U = pd.concat(out)
print(U.groupby("zone").agg(easting_min=("e", "min"), easting_max=("e", "max"), northing_min=("n", "min"), northing_max=("n", "max")).round(0).to_string())
print(f"worst round-trip error: {U.err_m.max() * 1000:.4f} mm")

      easting_min  easting_max  northing_min  northing_max
zone                                                      
4        340426.0     814588.0     2174141.0     2458824.0
5        185587.0     310354.0     2093123.0     2301497.0
worst round-trip error: 0.0000 mm

In their own zones, every easting lies between 185,587 and 814,588 m, within the range a zone is designed for, and the round trip back to latitude and longitude is exact to well under a millimetre. Note that the two zones' eastings overlap numerically, 340,426 to 814,588 in zone 4 and 185,587 to 310,354 in zone 5, yet the same easting in the two zones describes two different places: the zone number is part of the coordinate.

3. Forcing one zone on everything

Using a single zone for a dataset that spans two keeps all coordinates comparable but distorts the ones far from the zone's central meridian. Measure the scale distortion of zone 4 at the eastern tip of Hawaii:

z4 = pyproj.CRS.from_epsg(32604)
far = G.loc[G.lon.idxmax()]
f = pyproj.Proj(z4).get_factors(far.lon, far.lat)
f_own = pyproj.Proj(pyproj.CRS.from_epsg(utm_epsg(far.lat, far.lon))).get_factors(far.lon, far.lat)
print(f"{far.feature_name} at {far.lon:.3f}°: scale in zone 4 {f.meridional_scale:.6f}, in its own zone {f_own.meridional_scale:.6f}")
print(f"distortion: {(f.meridional_scale - 1) * 1e6:,.0f} ppm in zone 4 against {(f_own.meridional_scale - 1) * 1e6:,.0f} ppm in zone {utm_epsg(far.lat, far.lon) % 100}")

Kings Pillars at -154.807°: scale in zone 4 1.001996, in its own zone 1.000045
distortion: 1,996 ppm in zone 4 against 45 ppm in zone 5

At Kings Pillars, the easternmost feature, zone 4 stretches distances by 1,996 parts per million, 2 m per kilometre, because the point is 4.2° beyond zone 4's central meridian, outside the zone. In its own zone 5 the distortion is 45 ppm.

4. Reading coordinates back in the wrong zone

An easting and a northing do not say which zone they belong to, and every zone uses the same easting range, so coordinates read back in the wrong zone look perfectly valid. Convert Mauna Kea in zone 5N and read the numbers back in its own zone, the zone to the west, and the southern zone of the same number:

z5 = Transformer.from_crs(4326, 32605, always_xy=True)
e, n = z5.transform(mk.lon, mk.lat)                                   # Mauna Kea in zone 5N
for epsg in (32605, 32604, 32705):
    lon, lat = Transformer.from_crs(epsg, 4326, always_xy=True).transform(e, n)
    off = geod.inv(mk.lon, mk.lat, lon, lat)[2] / 1000
    print(f"read back as EPSG:{epsg}: {lat:.4f}, {lon:.4f}; {off:,.1f} km from the true point")

read back as EPSG:32605: 19.8206, -155.4681; 0.0 km from the true point
read back as EPSG:32604: 19.8206, -161.4681; 628.6 km from the true point
read back as EPSG:32705: -70.2333, -159.8622; 9,993.5 km from the true point

Read in zone 4N, the same numbers land at exactly the same latitude 6° farther west, 628.6 km away in the ocean. Read in zone 5S, they land at 70°S in Antarctica, 9,993.5 km away.

z5f = G[G.zone == 5]
e, n = z5.transform(z5f.lon.values, z5f.lat.values)
lon4, lat4 = Transformer.from_crs(32604, 4326, always_xy=True).transform(e, n)
d = geod.inv(z5f.lon.values, z5f.lat.values, lon4, lat4)[2] / 1000
print(f"zone 5 features read as zone 4: {len(z5f):,}; shift {d.min():,.0f} to {d.max():,.0f} km; longitude shift {np.median(z5f.lon.values - lon4):.2f}°")

zone 5 features read as zone 4: 2,557; shift 625 to 632 km; longitude shift 6.00°

Every one of the 2,557 zone 5 features moves by exactly 6.00° of longitude when read in zone 4, 625 to 632 km depending on latitude: a neighbouring zone is one zone width away.

The southern hemisphere's false northing makes the other mistake larger. Convert every southern airport in its proper southern zone and read it back as northern:

south = AP[AP.latitude_deg < 0].head(200)
rows = []
for _, r in south.iterrows():
    z = int((r.longitude_deg + 180) // 6) + 1
    e, n = Transformer.from_crs(4326, 32700 + z, always_xy=True).transform(r.longitude_deg, r.latitude_deg)
    lon, lat = Transformer.from_crs(32600 + z, 4326, always_xy=True).transform(e, n)
    rows.append((r.name, r.latitude_deg, lat))
S = pd.DataFrame(rows, columns=["iata", "true_lat", "read_lat"])
print(f"southern airports: {len(S)}; read as northern: latitudes {S.read_lat.min():.1f} to {S.read_lat.max():.1f}; non-finite: {(~np.isfinite(S.read_lat)).sum()}")
print(S.head(3).round(3).to_string(index=False))

southern airports: 182; read as northern: latitudes 37.3 to 89.3; non-finite: 0
iata  true_lat  read_lat
 HIR    -9.428    80.624
 WSI   -33.883    56.333
 LAE    -6.568    83.510

All 182 southern airports land between 37°N and 89°N. Honiara (HIR), at 9°S, appears at 81°N in the Arctic, and Western Sydney International (WSI), at 34°S, at 56°N: the 10,000,000 m false northing read as a real distance.

A wrong zone can be caught from the numbers. Eastings far outside about 166,000 to 834,000 m mean a point is not in the zone it claims, and northings near 10,000,000 m in a northern zone mean a southern point lost its hemisphere. Forcing all of Hawaii into zone 4:

e4, n4 = Transformer.from_crs(4326, 32604, always_xy=True).transform(G.lon.values, G.lat.values)
out = (e4 < 166_000) | (e4 > 834_000)
print(f"Hawaii in zone 4: eastings {e4.min():,.0f} to {e4.max():,.0f} m; outside 166,000-834,000: {out.sum():,} features")
print(f"of those, features that belong to zone 5: {G.zone.values[out].tolist().count(5):,} of {out.sum():,}")

Hawaii in zone 4: eastings 340,426 to 940,212 m; outside 166,000-834,000: 1,875 features
of those, features that belong to zone 5: 1,875 of 1,875

Forced into zone 4, Hawaiian eastings run up to 940,212 m, and 1,875 features are beyond 834,000 m. Every one of them belongs to zone 5, so here the range check identifies the misplaced points exactly.

The cure is to keep the zone — better, the EPSG code — with every coordinate, in a column or in the file's CRS metadata, so each row can be transformed back with the zone it was made in.

UTM in the southern hemisphere

The data used on this page

The USGS Geographic Names Information System domestic names for Hawaii, with decimal coordinates on NAD83, which for this purpose matches WGS84 to within about a metre, and the large airports of OurAirports, a public-domain database, for the southern hemisphere.

import io
import zipfile
from pathlib import Path

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

Path("data").mkdir(exist_ok=True)
def cached(path, url):
    p = Path(path)
    if not p.exists():
        r = requests.get(url, timeout=300)
        r.raise_for_status()
        p.write_bytes(r.content)
    return p
with zipfile.ZipFile(cached("data/DomesticNames_HI_Text.zip", "https://prd-tnm.s3.amazonaws.com/StagedProducts/GeographicNames/DomesticNames/DomesticNames_HI_Text.zip")) as z:
    G = pd.read_csv(io.BytesIO(z.read("Text/DomesticNames_HI.txt")), sep="|", encoding="utf-8-sig")
G = G[(G.state_name == "Hawaii") & G.prim_long_dec.between(-161, -154)].rename(columns={"prim_lat_dec": "lat", "prim_long_dec": "lon"}).reset_index(drop=True)
AP = pd.read_csv(cached("data/airports.csv", "https://davidmegginson.github.io/ourairports-data/airports.csv"), keep_default_na=False, na_values=[""])
AP = AP[(AP.type == "large_airport") & AP.iata_code.notna()].set_index("iata_code")
mk = G[G.feature_name == "Mauna Kea"].iloc[0]
geod = pyproj.Geod(ellps="WGS84")
print(f"{len(G):,} features in the main Hawaiian Islands; {len(AP):,} large airports")

7,303 features in the main Hawaiian Islands; 1,172 large airports

7,303 GNIS features lie in the main Hawaiian Islands, between 161°W and 154°W; the setup also loads 1,172 large airports for the southern hemisphere example.

Example 1: the southern hemisphere

South of the equator, UTM adds a false northing of 10,000,000 m so northings stay positive. Convert Sydney's airport with the north and the south zone of the same number:

syd = AP.loc["SYD"]
for epsg in (32656, 32756):
    e, n = Transformer.from_crs(4326, epsg, always_xy=True).transform(syd.longitude_deg, syd.latitude_deg)
    print(f"SYD in EPSG:{epsg}: easting {e:,.0f} m, northing {n:,.0f} m")
print(f"utm_epsg() picks EPSG:{utm_epsg(syd.latitude_deg, syd.longitude_deg)}")

SYD in EPSG:32656: easting 331,533 m, northing -3,757,676 m
SYD in EPSG:32756: easting 331,533 m, northing 6,242,324 m
utm_epsg() picks EPSG:32756

Sydney's airport has the same easting in both, 331,533 m, but its northing is -3,757,676 m in the northern zone and 6,242,324 m in the southern one, which adds the 10,000,000 m false northing. The function picks the southern zone, EPSG:32756.

Sixty zones, each its own projection

UTM divides the world between 80°S and 84°N into 60 zones of 6° of longitude, each with its own transverse Mercator projection centred on the zone's central meridian. A scale factor of 0.9996 on the central meridian spreads the distortion: about −400 parts per million in the middle of a zone, rising to about +1,000 ppm at its edges, 3° away, on the equator, and less at higher latitudes; beyond the edge it grows quickly. Eastings are offset by 500,000 m so they are positive across the zone; in the southern hemisphere northings are offset by 10,000,000 m so they are positive down to the South Pole.

Two things follow. Coordinates from different zones are in different coordinate systems: subtracting a zone 4 easting from a zone 5 easting is meaningless, even though both are in metres. And a dataset that spans a zone boundary needs a decision: transform everything into one zone and accept growing distortion away from it, use a projection designed for the whole area, such as a state plane or a custom transverse Mercator, or keep zones separate and only compare within them.

Norway, polar regions and axis order

• Norway and Svalbard: zones 32 and 31–37 have irregular boundaries there; estimate_utm_crs handles them, the formula above does not.
• Polar regions: UTM stops at 80°S and 84°N; the polar stereographic UPS system covers the poles.
• Datum: EPSG 326xx and 327xx are WGS84; NAD83 and ETRS89 have their own UTM codes (269xx, 258xx).
• always_xy: without it, pyproj follows EPSG axis order, latitude first for EPSG:4326, and swaps your inputs.
• Latitude band letters in MGRS-style references (C to X) are not hemisphere flags: band "S" is in the northern hemisphere.

Internal links

• How to convert degrees, minutes and seconds to decimal degrees in Python — the same features in DMS
• Latitude and longitude explained — the input coordinates

FAQ

How do I convert latitude and longitude to UTM in Python?

Find the zone as floor((lon + 180) / 6) + 1, use EPSG 326 + zone in the north or 327 + zone in the south, and transform with pyproj.Transformer.from_crs(4326, epsg, always_xy=True).

How do I find the UTM zone for a coordinate?

Zone = floor((longitude + 180) / 6) + 1, except in parts of Norway and Svalbard; GeoPandas' estimate_utm_crs() returns it as a CRS.

Why are my UTM northings around 10 million?

The point is in the southern hemisphere, where northings include a false northing of 10,000,000 m.

Can I mix coordinates from two UTM zones?

No: each zone is its own coordinate system; transform to one zone or to a common projection first.

Why are my UTM coordinates hundreds of kilometres from where they should be?

They were converted in one zone and read back in another: a neighbouring zone moves a point one zone width, 625 to 632 km at Hawaii's latitude, and a southern point read with a northern zone moves 10,000 km. Keep the zone or EPSG code with the coordinates.