Spatial Dashboard Streamlit · Panel · Production
Data Reference

Reprojecting Between EPSG:4326 and EPSG:3857 Without Distortion

Reproject with to_crs(3857) only for display, clip geometries to ±85.06° latitude first, and never measure distance or area in Web Mercator — switch to an equal-area or UTM CRS for any calculation.

Why this matters

Almost every web tile map — Folium, Leaflet, Deck.gl, pydeck — renders in Web Mercator (EPSG:3857), while most data is stored and exchanged in WGS 84 (EPSG:4326). So the single most common reprojection in a spatial dashboard is the hop between these two, and it is also the one most often done wrong. Two mistakes dominate: measuring area or distance in 3857, where the numbers are inflated by a latitude-dependent factor, and reprojecting geometries that extend past the projection’s ±85.06° limit, where coordinates blow up to infinity. Doing this hop correctly is the display half of the CRS & coordinate systems reference; the measurement half is covered by looking up EPSG codes for your region.

Prerequisites

Step-by-step solution

Step 1 — Understand why measurement in 3857 is wrong

Web Mercator is conformal: it preserves the shape of small features by stretching the map vertically as you move away from the equator. The vertical scale factor is 1 / cos(latitude), which means a shape at 60° latitude covers roughly four times the area on the map that the same shape covers at the equator. The diagram makes the distortion concrete.

Web Mercator distortion and the valid latitude bandThree equal-area cells on the globe are drawn at the equator, at 45 degrees, and at 70 degrees latitude. In Web Mercator the equator cell keeps its size, the 45 degree cell is enlarged, and the 70 degree cell is enlarged much more, showing that area inflates with latitude. A dashed line marks the 85.06 degree clipping limit beyond which the projection is undefined.EQUAL-AREA TRUTHWEB MERCATOR (3857)0° equator45° lat70° lat×1.0 area≈×2 area≈×8 areavalid band ends at ±85.06° — beyond this, 3857 is undefined

The takeaway: 3857 is a display CRS. Any number you need — a length, an area, a buffer radius — is computed in a different CRS.

Step 2 — Reproject the GeoDataFrame with to_crs

The transform itself is one call. Start from a frame you know is in 4326.

python
import geopandas as gpd
from shapely.geometry import Point

# Plausible points over Berlin and Hamburg (lon, lat)
gdf = gpd.GeoDataFrame(
    {"city": ["Berlin", "Hamburg"]},
    geometry=[Point(13.40, 52.52), Point(9.99, 53.55)],
    crs="EPSG:4326",
)

gdf_web = gdf.to_crs("EPSG:3857")
print(gdf_web.geometry.iloc[0])   # POINT (1491772.6 6894111.9) — metres

Berlin’s (13.40, 52.52) becomes roughly (1,491,773, 6,894,112) in Web Mercator metres. Those large easting/northing values are correct — they are metres from the projection origin, not a bug.

The scale factor, and why area is worse than lengthWeb Mercator's scale factor is one over the cosine of latitude, plotted from the equator to eighty degrees. It is exactly one at the equator, 1.15 at thirty degrees, 1.41 at forty-five, 2.0 at sixty and 2.92 at seventy — so a kilometre of ground measures nearly three kilometres in projected units at seventy degrees north. Areas scale with the square of that factor, so the same latitudes give area errors of one, thirty-three percent, double, quadruple and over eight times. Two markers show what this means concretely: a length measured near Milan is out by forty-one percent, and an area measured near Oslo is out by a factor of four. The projection is not defective — it is conformal, which is exactly what a tile pyramid needs, and preserving angle at the cost of area is the trade it was designed to make. It is simply not a system to measure in.SCALE FACTOR = 1 / cos(latitude)1×2×3×4×45° — a length is 41% too long60° — an area is 4× too large0°45°60°80°Length errors follow the curve; area errors follow its square, which is why a map that looks only slightly stretchedreports totals that are wrong by multiples. The projection is not defective — it is conformal, which is what a tile pyramid needs. ### Step 3 — Handle the ±85.06° latitude clipping

Web Mercator’s y coordinate tends to infinity as latitude approaches ±90°, so the projection is cut off at ±85.06°. A dataset with an Arctic vertex — a global coastline, an ice-extent polygon — will produce infinite coordinates after to_crs. Clip to the valid band before reprojecting.

python
from shapely.geometry import box

MERC_LAT = 85.06
valid_band = box(-180, -MERC_LAT, 180, MERC_LAT)

gdf_clipped = gpd.clip(gdf, valid_band)     # drop/trim beyond ±85.06°
gdf_web = gdf_clipped.to_crs("EPSG:3857")

# guard against any surviving infinities
assert gdf_web.geometry.is_valid.all()
assert gdf_web.total_bounds[3] < 2.0e7      # northing stays finite
Three systems, three jobs, one requestA single request is traced through the coordinate systems it touches. Storage and interchange happen in the geographic system, because it is what every source and every GeoJSON consumer agrees on. Rendering happens in Web Mercator, because that is what the tile pyramid is cut in and what the map component expects. Measurement — buffers, areas, distances, nearest-neighbour radii — happens in a third system chosen for the extent, either the local UTM zone or an equal-area projection, and its results are scalars rather than geometry, so they come back as numbers and never need reprojecting. The mistake the diagram is drawn to prevent is collapsing three jobs into two: measuring in the render system because the frame is already there. The reprojection for measurement costs milliseconds; the wrong total costs a decision.THREE JOBS, THREE SYSTEMSEPSG:4326 — store and exchangewhat every source and every GeoJSONconsumer already agrees onEPSG:3857 — renderwhat the tile pyramid is cut in andwhat the map component expectsUTM or equal-area — measurebuffers, areas, distances, radii —chosen for the extent, not for the mapresults are scalarsnumbers, not geometryThe collapse to avoidMeasuring in the render system,because the frame is already sittingthere in 3857 and the call returns aplausible-looking number.The reprojection costs milliseconds.The wrong total costs a decision.The three systems are not interchangeable and the conversions between them are cheap, so the only reason to reuseone for another's job is that nothing stopped you. ### Step 4 — Measure in an equal-area or UTM CRS

When a filter or tooltip needs an actual area or distance, reproject a copy into the right measurement CRS — never read it from the 3857 frame.

python
# area: use an equal-area CRS
gdf["area_km2"] = gdf.to_crs("EPSG:6933").area / 1_000_000

# distance between the two cities: use a local UTM zone (metres)
metric = gdf.to_crs("EPSG:25832")
d_km = metric.geometry.iloc[0].distance(metric.geometry.iloc[1]) / 1000
print(f"{d_km:.0f} km")   # ~255 km Berlin–Hamburg

Verification

Confirm the round-trip is lossless within the valid band — the fastest proof the transform is correct.

python
import geopandas as gpd
from shapely.geometry import Point

gdf = gpd.GeoDataFrame(
    geometry=[Point(13.40, 52.52)], crs="EPSG:4326"   # Berlin
)
roundtrip = gdf.to_crs("EPSG:3857").to_crs("EPSG:4326")

orig = gdf.geometry.iloc[0]
back = roundtrip.geometry.iloc[0]
assert abs(orig.x - back.x) < 1e-6, "longitude drifted on round-trip"
assert abs(orig.y - back.y) < 1e-6, "latitude drifted on round-trip"
print("round-trip clean")

Edge cases and gotchas

  • Poles clipping: any geometry crossing ±85.06° must be clipped first (Step 3), or a single vertex sends the whole layer’s bounds to infinity and the map fails to fit.
  • Axis order: if points arrive as (lat, lon) instead of (lon, lat), they will land in the wrong hemisphere. Construct Point(lon, lat) and, when using pyproj.Transformer directly, pass always_xy=True.
  • Area in 3857: calling .area on the reprojected frame returns an inflated value that grows with latitude. Always reproject to EPSG:6933, EPSG:3035, or EPSG:5070 for area (Step 4).

Round-tripping is not free, but it is usually fine

Converting from the geographic system to Web Mercator and back does not return the exact bits you started with. The forward and inverse transforms are both well conditioned away from the poles, but they run in floating point, so the returned coordinates differ from the originals in the last few decimal places — on the order of nanometres of ground distance in the middle latitudes.

That is far below anything a dashboard renders or measures, so a single round trip is harmless. What is worth avoiding is making round trips part of a loop: reprojecting on every rerun, or storing the reprojected result and reprojecting it again on the next interaction, accumulates the error and, much more importantly, accumulates the cost, since a transform over a large frame is not cheap. Reproject once at the boundary, keep both representations if both are genuinely needed, and treat the geographic version as the source of truth that the others are derived from rather than as one of several equal copies that drift apart.

FAQ

Why does a polygon look bigger near the poles in Web Mercator?

Web Mercator is a conformal projection that preserves local angles by stretching the map vertically as latitude increases. The scale factor is 1 / cos(latitude), so at 60° a shape is exaggerated about twofold in each direction and roughly fourfold in area. This is a display property, not a data error, but it means you must never compute area directly in EPSG:3857.

What happens to points beyond 85° latitude when reprojecting to 3857?

Web Mercator is only defined up to about ±85.06°, where the projected y coordinate would otherwise go to infinity. Points beyond that latitude produce infinite or extreme values after to_crs. Clip or filter geometries to the valid band before reprojecting so a single polar vertex does not break the whole layer.

Is reprojecting 4326 to 3857 and back lossless?

Within the valid latitude band the round trip is accurate to floating-point tolerance, so an assert with a small tolerance passes. Loss only appears if geometries were clipped at the ±85.06° limit or if coordinates were rounded to reduce precision between the two transforms.

Back to CRS & Coordinate Systems Reference

Related