DuckDB Spatial Analytical SQL for GIS

Choosing a Projected CRS for Distance Work

Every distance, area and buffer you compute inherits the distortion of the frame it is computed in, and Web Mercator — the frame most data arrives in — is the worst common choice for all three. This walkthrough, part of the CRS mapping and transformations reference, sets out how to pick a frame from the extent of the data rather than from habit, what each family of projections preserves, and how to check the error you have accepted.

Root-Cause Analysis: why the wrong frame is chosen so consistently

  • Web Mercator is the default of the web. Tiles, basemaps and most JavaScript mapping libraries use EPSG:3857, so data acquires it in transit and keeps it. It is a display projection: it preserves angles and destroys area, and its scale error grows with latitude.
  • The data arrived in degrees. EPSG:4326 is the interchange default, and planar distance in degrees is meaningless — a “distance” of 0.01 is about 1,100 m in longitude at the equator and about 700 m at 50° north, in the same dataset.
  • The frame was never stated. DuckDB stores no SRID, so nothing forces the question. A pipeline can compute areas for years without anyone naming the frame those areas are in.
  • One frame was chosen for a continent. A single projected frame is accurate over a limited extent. Applying a national grid to a continental dataset produces errors that grow with distance from its central meridian.
  • The error was never measured. Distortion is smooth and plausible: a distance that is four per cent wrong looks exactly like a distance that is right. Nothing surfaces it except a deliberate check.

The distinguishing question is what the numbers are for. A tolerance of a few per cent is irrelevant for a heat map and unacceptable for a land-area calculation, and the same frame can be right for one and wrong for the other.

What each projection family preserves Conformal preserves angles, equal-area preserves area, equidistant preserves distance from one origin, and a national grid is conformal with negligible error over its own extent. FAMILY PRESERVES USE IT FOR conformal (Mercator, UTM) local angles and shape maps; UTM for measurement in-zone equal-area (Albers, LAEA) area, exactly land-use and density statistics equidistant (azimuthal) distance from one origin catchments around a single point a national grid both, negligibly, in-extent the default where one exists

Only the last row is right by default, and only inside the extent it was designed for.

Deterministic Configuration 

INSTALL spatial; LOAD spatial;
SET memory_limit = '6GB';
SET threads = 8;

-- Reproject once into a materialised column rather than inside predicates:
-- a transform in a join is evaluated per candidate pair, and it also hides
-- the geometry column from the optimizer so no index can serve it.

Optimized Execution Pattern

The pattern is to reproject once at ingest into a frame chosen from the extent, and to keep the original only if it has to be re-derived.

-- ANTI-PATTERN: measuring in Web Mercator because that is what arrived.
-- The number is plausible and, at 55° north, about 3% too large.
SELECT ST_Area(geom) AS area_m2 FROM parcels;   -- geom is EPSG:3857
-- PATTERN: project once, into a frame appropriate to the extent, at ingest.
CREATE OR REPLACE TABLE parcels AS
SELECT parcel_id, zone_id,
       ST_Transform(geom, 'EPSG:4326', 'EPSG:27700') AS geom   -- British National Grid
FROM staging;

CREATE INDEX idx_parcels_geom ON parcels USING RTREE (geom);

Choosing the frame at ingest rather than per query also removes a whole class of inconsistency: two queries that reproject differently produce two different answers to the same question, and nothing in either says so.

Choosing the frame from the extent A national grid inside one country, a UTM zone inside one zone, an equal-area projection for a continent, and spherical maths for global distance. EXTENT OF THE DATA CHOOSE WHY one country with a national grid that grid designed for this extent inside one UTM zone that zone sub-metre over ~6° of longitude a country or continent, area matters an equal-area projection area is what the answer is about global, distance matters no planar frame compute on the sphere instead

The extent picks the row. The tolerance decides whether the row is good enough.

What Web Mercator actually costs

The scale factor in Web Mercator is the secant of the latitude, so linear distances are overstated by that factor and areas by its square. At the equator the error is nil. At 40° north a distance is about 30% too long and an area about 70% too large. At 60° north an area is four times too large. These are not subtle errors, and they are entirely invisible in the output.

The reason they survive is that most uses of Web Mercator are relative rather than absolute: a heat map, a proportional symbol, a ranking. Those are unaffected as long as everything is in the same frame, which is why the habit persists. The moment a number is reported in metres or square metres, the frame stops being a display choice and becomes a measurement claim.

Diagnostic Queries & Plan Validation

The check that settles the question is a known distance. Pick two points whose separation is known independently, compute it in the candidate frame, and compare.

-- A known baseline, measured in the candidate frame. The relative error is
-- the number to compare against the tolerance you wrote down earlier.
WITH baseline AS (
    SELECT ST_Transform(ST_Point(-0.0014, 51.4778), 'EPSG:4326', 'EPSG:27700') AS a,
           ST_Transform(ST_Point(-0.1246, 51.5007), 'EPSG:4326', 'EPSG:27700') AS b,
           8_600.0 AS known_metres          -- from an independent source
)
SELECT ST_Distance(a, b)                              AS computed_metres,
       abs(ST_Distance(a, b) - known_metres) / known_metres AS relative_error
FROM baseline;

A relative error inside your stated tolerance means the frame is adequate for this dataset. One outside it means the extent is larger than the frame was designed for, and the answer is a different frame rather than a correction factor.

Web Mercator area error by latitude Nil at the equator, about 13% at 20°, about 70% at 40°, and about 300% at 60°. LATITUDE AREA ERROR A 1 KM² PARCEL REPORTS AS 0° none 1.00 km² 20° ~13% too large 1.13 km² 40° ~70% too large 1.70 km² 60° ~300% too large 4.00 km²

Smooth, monotonic and invisible, which is exactly why it survives.

When degrees are actually fine

Not every question needs a projected frame. A bounding-box filter is a comparison of coordinates against coordinates and is exact in any frame, so a spatial index and an envelope test work perfectly well in degrees. A containment test is topological rather than metric and is likewise unaffected. What needs a projected frame is anything that produces a number in ground units — distance, length, area, buffer radius, or a density expressed per square kilometre.

That distinction is worth applying deliberately, because reprojecting a large layer is not free and doing it for a query that did not need it is pure cost. The workable rule is to keep the working frame projected when measurement is the point, and to accept degrees when the pipeline only ever filters and joins.

Geometry Validation & Fallback Routing

Where no single frame covers the extent, the fallback is to partition the computation by zone rather than to accept the error.

-- Per-zone reprojection: each row is measured in the frame appropriate to
-- where it is, and the results are comparable because all of them are metres.
SELECT parcel_id,
       ST_Area(ST_Transform(geom, 'EPSG:4326', 'EPSG:' || (32600 + utm_zone))) AS area_m2
FROM (
    SELECT parcel_id, geom,
           floor((ST_X(ST_Centroid(geom)) + 180) / 6)::INT + 1 AS utm_zone
    FROM parcels_wgs84
);

Frequently Asked Questions

Is Web Mercator ever the right frame for measurement?

Only at the equator, and only by coincidence. Its scale factor is the secant of the latitude, so distances are overstated by that factor and areas by its square — about 70% too large at 40° north and four times too large at 60°. It is a display projection, and using it for a reported number turns a rendering choice into a measurement claim.

Can I just measure in degrees and convert?

No, because a degree is not a fixed distance. A degree of latitude is roughly constant, but a degree of longitude shrinks toward the poles, so any conversion factor is only correct at one latitude. Within a single dataset spanning a few hundred kilometres north to south, the factor changes measurably.

Which frame should I use if the data spans several UTM zones?

If area is the answer, an equal-area projection centred on the region. If distance is the answer and the extent is large, compute per zone and combine, or use spherical distance directly. A single UTM zone applied outside its own extent degrades quickly and the degradation is not uniform.

How do I record which frame a table is in?

As metadata your pipeline owns — a companion table, a naming convention, or a configuration entry. DuckDB stores no SRID, so nothing enforces it, and a frame that lives only in someone’s memory is a frame that will be wrong within a year of them leaving.

Does the choice affect index behaviour?

Not directly — an R-tree over envelopes works in any frame. What it affects is the meaning of a radius: ST_DWithin(a, b, 500) means 500 of whatever unit the coordinates are in, so the same query changes meaning entirely between a geographic and a projected frame while continuing to return rows.

Is reprojecting the whole layer worth it?

Whenever measurement is the point, yes, and it is a one-off cost recovered on every later query — a transform inside a join predicate is evaluated per candidate pair rather than per row, which is routinely an order of magnitude more work. Where the pipeline only filters and joins, degrees are fine and the reprojection is pure cost.

Up: CRS Mapping & Transformations in DuckDB Spatial