Spatial Data Archival · Cold Storage

Partitioning GeoParquet by H3 Spatial Index

Partitioning a GeoParquet archive by H3 cell turns a directory of opaque columnar files into a spatially prunable dataset, where a query engine reads only the partitions whose hexagons intersect a bounding box instead of scanning every object in cold storage. This guide is for data engineers and GIS archivists who need a deterministic, cold-storage-aware partition scheme: choosing an H3 resolution that balances predicate-pushdown selectivity against a cloud provider’s minimum-object-size billing, deriving stable cell keys from geometry centroids, and writing a Hive-style h3_cell= layout that DuckDB, GDAL, and pyarrow can prune without opening files. Default writers ignore this entirely — they emit a flat file set with no partition column, so every spatial filter degrades into a full-archive scan. It operates under the Format Conversion & Pipeline Automation framework and extends the batch GeoParquet Migration Workflows reference with a partition-key strategy tuned for archival economics.

Why H3 for Archive Partitioning

H3 is a hierarchical hexagonal grid: every cell has exactly one parent at the next-coarser resolution, so a partition scheme built on it is reprojectable up the hierarchy without recomputing keys from geometry. Hexagons also have uniform neighbor distance — unlike a quadkey grid, adjacent cells are equidistant, which keeps a radius query from fetching a lopsided set of partitions. The single decision that governs everything downstream is resolution: it fixes cell area, and therefore the feature count and byte size of each partition file.

H3 resolution versus partition file size trade-off A four-row matrix mapping H3 resolution to average cell area, pruning granularity, and typical partition file size. Coarse resolutions yield oversized files with weak pruning; fine resolutions prune well but produce tiny objects that fall below cold-storage minimum billing. Resolution 4 is highlighted as the balanced default. Resolution Avg cell area Pruning granularity Typical partition file res 3 ~12,400 km² coarse · whole regions multi-GB · weak pushdown res 4 ~1,770 km² regional · balanced 200 MB–1 GB · default res 5 ~252 km² fine · metro-scale 30–200 MB res 6 ~36 km² very fine <30 MB · file-count bloat

The resolution you pick is a storage-economics decision as much as a query one, so make it against the same cold-tier minimums that drive Hot/Warm/Cold Tier Design for Geospatial Data — pick the coarsest resolution whose median file still clears the tier’s minimum billable object size.

Sizing the Resolution Against Cold-Storage Minimums

Cold classes bill a minimum object size and a minimum retention duration regardless of the bytes you actually store. AWS Glacier Deep Archive rounds every object up to a 128 KB billable floor and adds a per-object overhead, so a partition scheme that shatters a dataset into millions of sub-megabyte files pays overhead on every one and cripples restore throughput. The countervailing pressure is pushdown: coarser cells mean fewer, larger files, so a tight bounding-box query still drags in gigabytes of irrelevant features.

Compute the resolution empirically rather than guessing. Sample the feature density of the source, project a per-cell byte estimate at each candidate resolution, and pick the coarsest one whose median file lands in your target range (a common archival target is 200 MB–1 GB per object, large enough to amortize per-object overhead, small enough to parallelize restore).

import h3
import geopandas
import numpy as np

# Sample 1% of a large source to estimate per-cell feature density.
sample = geopandas.read_parquet("s3://spatial-archive/staging/parcels_sample.parquet")
centroids = sample.geometry.to_crs(4326).centroid
avg_bytes_per_feature = 320  # measured: WKB geometry + attribute columns, post-ZSTD

for res in (3, 4, 5, 6):
    cells = [h3.latlng_to_cell(pt.y, pt.x, res) for pt in centroids]
    counts = np.array(list({c: cells.count(c) for c in set(cells)}.values()))
    # scale the 1% sample up to the full population
    est_features = counts * 100
    est_mb = est_features * avg_bytes_per_feature / 1_048_576
    print(f"res={res}: partitions={len(counts):>6d}  "
          f"median_file={np.median(est_mb):8.1f} MB  p95={np.percentile(est_mb,95):8.1f} MB")

Read the output as a curve: resolution is right when the median file clears your cold-tier floor and the p95 file stays under the restore-parallelism ceiling you can tolerate. If the p95 is an order of magnitude above the median, the data is spatially skewed (dense urban cells, empty ocean cells) and a single global resolution will not fit — see the mitigation below.

Resolution, Cell Area, and What Fits in One Object

H3 resolutions differ by a factor of seven in area at each step, so the choice between two adjacent resolutions changes partition size by nearly an order of magnitude. Reading the ladder against a known feature density makes the decision arithmetic rather than intuition.

H3 resolution against cell area and partition size Five H3 resolutions with edge length, cell area and the resulting partition size for a parcel dataset at 340 features per square kilometre. Resolution 6 lands at about 320 megabytes and is marked as the target. resedge lengthcell areafeatures/cellpartition size 422.6 km1,770 km² 602,00015.4 GB 58.5 km253 km² 86,0002.2 GB 63.2 km36.1 km² 12,300320 MB — target 71.2 km5.16 km² 1,75046 MB 80.46 km0.74 km² 2506.6 MB — too many objects Worked at 340 features/km² and 26 KB per feature; substitute your own density and the target resolution falls out directly.

Density is rarely uniform, so treat the table as a starting point and then check the tail: compute the actual per-cell feature count across the dataset, and if the busiest cell exceeds the object-size ceiling, subdivide only those cells to the next resolution rather than moving the whole archive down a step.

Deriving Stable Cell Keys and Writing the Layout

Derive the H3 index from the geometry centroid in EPSG:4326, because H3 is defined on the WGS84 sphere and latlng_to_cell expects geographic coordinates. Deriving from a projected CRS silently places features in the wrong cell. Pin the CRS explicitly — the same discipline enforced by CRS Synchronization in Pipelines — and store the resolution in the archive manifest so future writers reproduce identical keys.

import json
import geopandas
import pyarrow as pa
import pyarrow.parquet as pq

H3_RES = 4  # pinned in the archive manifest; never change without a full repartition

def write_h3_partitioned(src: str, dst_root: str, crs_epsg: int = 4326):
    gdf = geopandas.read_parquet(src)
    # Centroid in EPSG:4326 is the ONLY correct input to H3.
    centroids = gdf.geometry.to_crs(4326).centroid
    gdf["h3_cell"] = [h3.latlng_to_cell(p.y, p.x, H3_RES) for p in centroids]

    geo_meta = {
        "version": "1.1.0",
        "primary_column": "geometry",
        "columns": {"geometry": {"encoding": "WKB",
                                  "geometry_types": sorted(gdf.geom_type.unique().tolist()),
                                  "crs": f"EPSG:{crs_epsg}"}},
    }
    gdf["geometry"] = gdf.geometry.to_wkb()

    for cell, part in gdf.groupby("h3_cell"):
        table = pa.Table.from_pandas(part, preserve_index=False)
        table = table.replace_schema_metadata({b"geo": json.dumps(geo_meta).encode()})
        # Hive-style directory: the engine reads h3_cell from the path, not the file.
        out = f"{dst_root}/h3_cell={cell}/part-0000.parquet"
        pq.write_table(table, out, compression="zstd", compression_level=3,
                       row_group_size=100_000)

The written tree is Hive-partitioned: s3://spatial-archive/parcels/h3_cell=8428309ffffffff/part-0000.parquet. The h3_cell value lives in the path, not inside the file, so a partition-pruning engine skips whole directories before it opens a single Parquet footer. Note the h3_cell column is dropped from the file body by convention — it is redundant with the path — which also shrinks each object.

For skewed data where one resolution cannot serve both dense and sparse regions, partition at a coarse base resolution and let row-group sizing handle intra-partition selectivity instead of pushing to a finer H3 level; the Row Group Sizing Strategies reference covers picking a row-group size that keeps in-file pruning effective on multi-gigabyte partitions.

Verifying Predicate Pushdown on the Partition Column

Writing the layout is only half the job; confirm that engines actually prune it. DuckDB with hive_partitioning=1 exposes h3_cell as a virtual column and skips non-matching directories.

duckdb -c "
  SET enable_object_cache=true;
  EXPLAIN ANALYZE
  SELECT count(*) FROM read_parquet(
      's3://spatial-archive/parcels/**/*.parquet', hive_partitioning=1)
  WHERE h3_cell IN ('8428309ffffffff','842830bffffffff');"

The EXPLAIN ANALYZE output must show that only the matching partitions were read. Look for the Files Read counter — it should equal the number of cells in your IN list, not the archive total:

┌─────────────────────────────────────┐
│           PARQUET_SCAN              │
│   Files Read: 2 / 18443             │  ← pruned to the two matching cells
│   Rows Scanned: 214,880            │
│   Total Time: 0.31s                │
└─────────────────────────────────────┘

If Files Read equals the archive total, pruning failed — the filter did not reference the partition column, or hive_partitioning was off. To resolve a bounding box to the set of cells to query, use h3.polygon_to_cells on the query envelope and pass the result as the IN list; this is how you translate a map viewport into a partition predicate.

Cell Keys That Stay Stable Across Reloads

A partition key must be a pure function of the geometry, or a reload will file the same feature somewhere new and leave the archive holding two copies under different keys. Three ordinary pipeline choices break that purity.

Unstable versus stable derivations of an H3 partition key Three unstable key derivations — centroid after reprojection, first vertex, and bounding-box centre — each with the change that moves the key, contrasted with a stable derivation from a point-on-surface in a fixed reference system. centroid computed after reprojection moves when the PROJ pipeline or grid inventory changes first vertex of the ring moves when a tool normalises ring order or winding centre of the bounding box moves when any single outlying vertex is edited point-on-surface in a pinned CRS stable — deterministic, inside the polygon, independent of vertex order The key must be a pure function of the geometry, or reloads split a feature across partitions

Pin the reference system the key is computed in, record it next to the resolution in the manifest, and compute the key once at ingest rather than recomputing it in each downstream job. Where a key does change legitimately — because a feature’s geometry was corrected — treat it as a delete in the old partition and an insert in the new one, so the archive never holds the same identifier twice.

Troubleshooting Partition Pruning

Symptom Cause Fix
Query scans every file despite a spatial filter Filter uses ST_Intersects on geometry, not the h3_cell path column Resolve the query envelope to cells with h3.polygon_to_cells and filter WHERE h3_cell IN (...)
Millions of sub-megabyte partition files Resolution too fine for the data’s spatial density Drop one H3 resolution and re-derive keys via h3.cell_to_parent(cell, res-1) — no geometry recompute needed
One partition is 40× larger than the median Spatial skew (dense metro cell) at a coarse resolution Keep the coarse base and rely on row-group pruning inside the file rather than splitting the whole archive finer
Features land in the wrong cell Centroid computed in a projected CRS instead of EPSG:4326 Always to_crs(4326) before latlng_to_cell; H3 is defined on WGS84

Operational Execution Checklist

Frequently Asked Questions

Should the H3 cell be stored as a column as well as a directory level?

Yes. The directory level prunes at planning time, before any file is opened; the column lets an engine filter within a file and lets a reader recover the key without parsing the path. The duplication costs eight bytes per row as an integer and is recovered many times over by dictionary encoding, since the value is constant within each partition.

What happens to features whose geometry spans several cells?

Assign each to the cell containing its representative point and rely on row-group bounding-box statistics for extent queries that clip a neighbour. Duplicating a feature into every intersecting cell makes counts unreliable and forces every reader to deduplicate; it is only worth it for datasets where features are large relative to the cell and exact per-cell completeness is a stated requirement.

Does H3 partitioning survive a resolution change?

Partially, and in one direction. A finer partition can be rolled up because every cell has a well-defined parent, so a resolution-7 archive can be presented as resolution-6 without recomputing keys. Going finer requires recomputing from geometry, because a parent cell does not determine which child each feature belongs to. Choose the finest resolution you are willing to store, and aggregate upward when queries want coarser.

Part of the Spatial Data Archival knowledge base.

Related pages

Converting Legacy Shapefiles to GeoParquet at Scale Incremental GeoParquet Updates Without a Full Rewrite Writing GeoParquet Bounding-Box Metadata for Predicate Pushdown