DuckDB Spatial Analytical SQL for GIS

Sorting Writes with Hilbert Curves

Row-group statistics can only skip a block whose value range excludes your predicate, and a two-dimensional predicate has no value range unless the write imposed one — this walkthrough, part of the partitioning and file layout reference, covers using a space-filling curve as the sort key so that a bounding-box query lands on a contiguous run of row groups instead of scattering across all of them.

Root-Cause Analysis: why spatial filters do not prune by default

An attribute filter prunes because the writer sorted by that attribute, so each row group covers a narrow slice of its range. A spatial filter has no equivalent by default, and the reasons are worth separating because each has a different remedy.

  • Geometry has no total order. There is no meaningful “less than” between two polygons, so ORDER BY geom either errors or sorts by serialised bytes, which correlates with nothing geographic. The statistics that result are technically present and completely useless.
  • Sorting by a single coordinate helps in one dimension only. ORDER BY ST_X(ST_Centroid(geom)) produces row groups that are narrow in longitude and full-width in latitude, so a compact bounding-box query still touches every group that overlaps its longitude band — which on a national dataset is most of them.
  • Statistics on a blob column are not usable. Even where a geometry column carries min/max statistics, they are over the serialised bytes and mean nothing spatially. The pruning has to happen on plain numeric columns.
  • Arrival order correlates with the wrong thing. Data usually arrives ordered by ingestion date or by source file, both of which mix geography thoroughly. Every row group ends up spanning the whole extent.

A space-filling curve addresses the first two directly: it maps two dimensions into one in a way that keeps nearby points nearby in the ordering, so a compact region in space becomes a small number of short runs in the sort.

Four write orders, and how many row groups a bbox query must read From 200 row groups: arrival order reads all 200, longitude order about 40, Z-order about 14, Hilbert order about 9. WRITE ORDER GROUPS READ (of 200) WHY arrival order 200 — all of them every group spans the extent by longitude alone ~40 narrow in x, full height in y Z-order (Morton) key ~14 locality in both, with jumps Hilbert key ~9 the curve never makes a long jump

Same query, same data, same file size. Only the order in which the rows were written differs.

Deterministic Configuration

Sorting on write is an external sort over the whole dataset, which is the most memory-hungry step in the entire pipeline.

INSTALL spatial; LOAD spatial;

-- The sort dominates. It will spill on anything worth sorting, so the
-- question is not whether to give it a spill target but how fast that
-- target is.
SET memory_limit = '12GB';
SET threads = 6;
SET temp_directory = '/var/tmp/duckdb_sort';   -- fast local disk, never a network mount
SET max_temp_directory_size = '120GB';
SET preserve_insertion_order = false;          -- lets the writer stream the sorted output

Optimized Execution Pattern

The pattern has two halves that only work together: derive a curve key and sort by it, and materialise the envelope as plain numeric columns so statistics have something to prune on.

-- ANTI-PATTERN: sorted spatially, but there is no numeric column for the
-- statistics to describe, so the ordering buys nothing at read time.
COPY (SELECT * FROM parcels ORDER BY ST_Hilbert(geom, ST_Extent(geom) OVER ()))
  TO 'parcels.parquet' (FORMAT PARQUET);
-- PATTERN: sort by the curve AND write the envelope as four DOUBLEs, so the
-- row-group statistics describe something a bbox predicate can compare against.
CREATE OR REPLACE TABLE parcels_keyed AS
WITH extent AS (SELECT ST_Extent(geom) AS e FROM parcels)
SELECT
    p.*,
    ST_XMin(p.geom) AS bbox_xmin, ST_YMin(p.geom) AS bbox_ymin,
    ST_XMax(p.geom) AS bbox_xmax, ST_YMax(p.geom) AS bbox_ymax,
    ST_Hilbert(p.geom, (SELECT e FROM extent)) AS hkey
FROM parcels p;

COPY (SELECT * EXCLUDE (hkey) FROM parcels_keyed ORDER BY hkey)
  TO 'parcels.parquet' (FORMAT PARQUET, ROW_GROUP_SIZE 30000);

-- The read then prunes on plain numbers, before any WKB is decoded.
SELECT count(*) FROM read_parquet('parcels.parquet')
WHERE bbox_xmax >= 2.20 AND bbox_xmin <= 2.45
  AND bbox_ymax >= 48.80 AND bbox_ymin <= 48.95;

The EXCLUDE (hkey) matters more than it looks. The key exists to impose an order, not to be queried, and carrying it into the file adds a column to every row for no benefit — while also tempting a later reader to filter on it, which would couple every query to the extent the key was computed against.

Diagnostic Queries & Plan Validation

The measurement that matters is how many row groups a representative bounding-box query actually reads, compared with how many exist.

-- Groups in the file, versus groups the query touched. If they are equal,
-- the ordering did not help — check that the bbox columns were written and
-- that the predicate is on them rather than on the geometry.
SELECT count(DISTINCT row_group_id) AS groups_total FROM parquet_metadata('parcels.parquet');

EXPLAIN ANALYZE
SELECT count(*) FROM read_parquet('parcels.parquet')
WHERE bbox_xmax >= 2.20 AND bbox_xmin <= 2.45
  AND bbox_ymax >= 48.80 AND bbox_ymin <= 48.95;

Diagnostic — how tight the ordering actually is: the min/max of the bbox columns per row group tells you directly whether the curve did its job. A well-ordered file has groups whose extents are small and roughly square; a badly ordered one has groups whose extents each cover most of the dataset.

-- Per-group extent as a fraction of the whole. Values near 1 mean the group
-- spans the dataset and will never be pruned; values near 0 mean tight
-- spatial clustering and effective pruning.
SELECT row_group_id,
       (max_value::DOUBLE - min_value::DOUBLE) AS x_span
FROM parquet_metadata('parcels.parquet')
WHERE path_in_schema = 'bbox_xmin'
ORDER BY x_span DESC
LIMIT 10;
Z-order and Hilbert compared as write-order keys Z-order is cheaper to compute and jumps at power-of-two boundaries; Hilbert is extent-dependent, never jumps, and prunes better. PROPERTY Z-ORDER (MORTON) HILBERT cost to compute a bit interleave — trivial a real transformation long jumps in the order at power-of-two boundaries never depends on the extent no, if quantised consistently yes — stated at build time pruning achieved good better — usually by 30–50% Either is a large improvement on arrival order. The choice between them is a refinement, not the decision. When a spatial sort repays its cost Write-once-query-many with spatial filters repays it fully; write-once-query-once does not; frequent rewrites should compact periodically instead; attribute-only workloads should sort by the attribute. WORKLOAD DOES THE SORT REPAY? WHAT TO DO written once, queried often, spatially fully sort on write — the ideal case written once, queried once no skip it entirely rewritten every load partly move the sort to periodic compaction queried by attribute only no sort by that attribute instead

The ratio of queries to writes decides it, and whether those queries are spatial at all.

Geometry Validation & Fallback Routing

Two properties have to hold before a curve key is meaningful, and both fail quietly. The geometry must be in a projected frame, because a curve computed over degrees inherits the same latitude distortion that makes degree-based grids unequal. And the extent used to build the key must cover the whole dataset — a key built from a subset’s extent maps out-of-range geometries to the boundary, which piles them all into one place in the ordering.

-- Guard both before computing the key. A layer still in degrees, or an
-- extent that does not contain every row, produces a key that orders badly
-- without ever failing.
SELECT
    count(*) FILTER (WHERE abs(ST_XMin(geom)) <= 180) AS looks_like_degrees,
    count(*) FILTER (WHERE geom IS NULL OR ST_IsEmpty(geom)) AS unkeyable_rows
FROM parcels;

If either count is non-zero, fix it before sorting rather than after: rows with no geometry cannot be given a key at all and will sort to one end, and a layer in degrees will produce an ordering whose row groups are tall and thin in ground terms even though they look square in coordinate terms.

Frequently Asked Questions

Why do I need the bbox columns if I have already sorted spatially?

Because statistics are computed per column, and the geometry column’s statistics are over serialised bytes, which mean nothing spatially. The sort puts nearby features in nearby rows; the bbox columns are what lets the reader see that, because they are plain doubles whose min and max per row group describe a real region.

Is Hilbert better than Z-order?

Usually by a moderate margin — typically thirty to fifty per cent fewer row groups read — because a Hilbert curve never makes a long jump, so a compact query region stays contiguous in the ordering. Z-order is cheaper to compute and jumps at power-of-two boundaries, splitting some regions across distant runs. Either is a large improvement on arrival order, which is the decision that actually matters.

How much does the sort cost?

One full external sort over the dataset, which will spill and will take minutes on anything large. It is paid once at write time and recovered on every subsequent query, so it is right for a table read many times and wrong for a staging table read once and dropped.

Does the key have to be stored in the file?

No, and it usually should not be. Its job is to impose an order during the write; once the rows are written in that order the key has done its work. Storing it adds a column to every row and invites a later query to filter on it, which would couple that query to the extent the key was built against.

Up: Partitioning and File Layout