#Coordinate Precision Reduction for Smaller GeoParquet

Snap every ordinate to a grid matched to the data’s real positional accuracy, and a geometry column typically loses 40 to 60 per cent of its compressed size — with a bounded, documented, and uniform positional error. The saving is not from storing fewer bytes per coordinate; a float64 is always eight bytes. It is from making the low-order mantissa bits identical across the column, so the compressor collapses what was previously incompressible noise. This page is the operational procedure behind geometry encoding and coordinate precision.

#Quick Reference

#Why the Grid, Not the Type

A compressor is a predictor. It shrinks data by noticing that a byte is likely given the bytes around it, and it cannot shrink data it cannot predict. In a column of float64 coordinates spanning one country, the sign bit, the exponent, and the top of the mantissa are nearly constant — highly predictable, heavily compressed. The bottom three bytes of every ordinate are effectively uniform random, and uniform random data compresses to itself.

Rounding to a grid sets those bits to a fixed pattern. The column is the same width, every value is still a float64, and the compressed column is dramatically smaller because the unpredictable part has been removed. This is why the saving only appears after compression, and why teams who measure the raw column before and after conclude that nothing happened.

#The Procedure

# Requires: shapely>=2.0, geopandas>=1.0, numpy>=1.26  (Python 3.10+)
from __future__ import annotations

from dataclasses import dataclass

import geopandas as gpd
import numpy as np
from shapely import get_coordinates, set_precision
from shapely.errors import GEOSException

METRES_PER_DEGREE = 111_320.0


@dataclass(frozen=True)
class ReductionReport:
    grid_crs_units: float
    accuracy_m: float
    features: int
    collapsed: int
    invalid_added: int
    max_shift_m: float
    p99_shift_m: float

    def summary(self) -> str:
        return (
            f"grid {self.grid_crs_units:g} · max shift {self.max_shift_m:.3f} m · "
            f"p99 {self.p99_shift_m:.3f} m · {self.collapsed} collapsed · "
            f"{self.invalid_added} newly invalid"
        )


def reduce_precision(
    frame: gpd.GeoDataFrame,
    accuracy_m: float,
    *,
    max_collapse_fraction: float = 0.0,
) -> tuple[gpd.GeoDataFrame, ReductionReport]:
    """Snap to an absolute grid and prove the result is still the same data."""
    if frame.crs is None:
        raise ValueError("a CRS is required — grids are expressed in CRS units")
    if frame.empty:
        raise ValueError("nothing to reduce")

    geographic = bool(frame.crs.is_geographic)
    grid = accuracy_m / METRES_PER_DEGREE if geographic else accuracy_m
    scale = METRES_PER_DEGREE if geographic else 1.0

    before = frame.geometry.to_numpy()
    invalid_before = int((~frame.geometry.is_valid).sum())

    try:
        # An absolute grid anchored at the CRS origin — identical for every
        # geometry, so two features that shared a vertex still share it.
        after = set_precision(before, grid_size=grid)
    except GEOSException as exc:
        raise RuntimeError(f"precision reduction failed at grid {grid}: {exc}") from exc

    out = frame.copy()
    out.geometry = after

    collapsed = int(sum(1 for g in after if g is None or g.is_empty))
    if collapsed > max_collapse_fraction * len(frame):
        raise RuntimeError(
            f"{collapsed} feature(s) collapsed at grid {grid} — smaller than one "
            f"grid cell. Use a finer grid or filter them deliberately."
        )

    invalid_added = int((~out.geometry.is_valid).sum()) - invalid_before
    if invalid_added > 0:
        raise RuntimeError(f"grid {grid} introduced {invalid_added} invalid geometries")

    shifts = []
    for old, new in zip(before, after):
        if old is None or new is None or new.is_empty:
            continue
        a, b = get_coordinates(old), get_coordinates(new)
        shifts.append(grid if a.shape != b.shape else float(np.abs(a - b).max()))
    arr = (np.asarray(shifts) if shifts else np.zeros(1)) * scale

    return out, ReductionReport(
        grid_crs_units=grid,
        accuracy_m=accuracy_m,
        features=len(frame),
        collapsed=collapsed,
        invalid_added=max(invalid_added, 0),
        max_shift_m=float(arr.max()),
        p99_shift_m=float(np.percentile(arr, 99)),
    )

#Validation

Two things must be true after reduction: the geometry is still valid, and the positional error is inside the budget you declared. Both are cheap to assert and neither is optional.

# Requires: geopandas>=1.0, shapely>=2.0 — prove the reduction was within budget
import geopandas as gpd
from shapely import get_coordinates
import numpy as np

original = gpd.read_parquet("parcels.full.parquet")
reduced = gpd.read_parquet("parcels.10cm.parquet")

assert len(original) == len(reduced), "feature count changed — something collapsed"
assert reduced.geometry.is_valid.all(), "invalid geometry after reduction"

shift = max(
    float(np.abs(get_coordinates(a) - get_coordinates(b)).max())
    for a, b in zip(original.geometry, reduced.geometry)
    if a is not None and b is not None
    and get_coordinates(a).shape == get_coordinates(b).shape
)
print(f"max positional shift: {shift * 111_320:.3f} m")   # expect ≤ the declared grid

Expected results at a 10 cm grid on a metre-based CRS: maximum shift at or just under 0.10 m, p99 shift around 0.06 m, zero collapsed features on any layer whose smallest feature exceeds a metre, and a compressed geometry column 55–60% smaller.

#Edge Cases and Caveats

Features smaller than one grid cell. A 60 cm inspection chamber on a 1 m grid becomes a point or disappears. The implementation refuses to proceed by default rather than silently dropping features; if some collapse is genuinely acceptable, raise the threshold explicitly so the decision is recorded in code review rather than in a log line nobody read.

Snapping between projections. Reduce precision in the CRS you intend to store, as the last step before writing. A grid applied before a re-projection becomes an irregular, latitude-dependent spacing afterwards, which defeats both the size argument and the accuracy guarantee.

Layers in a family reduced at different grids. Parcels at 1 cm and the roads that bound them at 1 m will not share edges, and every overlay produces slivers. Record the grid in the dataset metadata alongside the CRS, the way metadata preservation records the projection, so a later reader can check compatibility before joining.

#Frequently Asked Questions

How much smaller does a geometry column actually get?

On a typical national polygon layer, 40 to 45 per cent at a one-centimetre grid and 55 to 60 per cent at ten centimetres, measured after ZSTD level 3. Point layers gain less because they carry fewer ordinates per feature relative to their attributes, and dense coastline or contour data gains more. The saving comes entirely from making the low mantissa bits constant so the entropy coder can collapse them.

What grid should I use if I do not know the data’s accuracy?

Find out before you snap. That said, there is a safe diagnostic: compute the distribution of distances between consecutive vertices. If the smallest gaps cluster around a value, that value is close to the digitising resolution and a grid at or below it is safe. If they extend smoothly down to sub-millimetre, the data has been through a re-projection and the tail is arithmetic noise rather than measurement.

Can precision reduction be undone?

No. The discarded bits are gone, and no later processing recovers them. This is why the recommended pattern is to keep one archival copy of survey-grade source data at full precision and apply reduction only to derived delivery copies. Storage for a single master is cheap relative to the cost of discovering, months later, that a boundary dataset was quantized past what a legal process requires.

#Related

• Geometry Encoding and Coordinate Precision — parent guide: the accuracy budget and the three encoding levers
• WKB vs Native Arrow Geometry Encoding — the container choice that decides whether delta encoding can act
• Simplifying Geometries Before Compression — removing vertices instead of digits
• Optimal ZSTD Levels for Vector vs Raster Data — the stage that runs after this one

← Back to Geometry Encoding and Coordinate Precision