20 min read

Data Conversion & Migration Pipelines for Cloud-Native Geospatial Storage

The geospatial data landscape is undergoing a fundamental architectural shift. Legacy vector formats — Shapefiles, GeoJSON, and KML — were engineered for desktop GIS workflows and single-threaded processing. At enterprise scale they introduce severe bottlenecks: excessive I/O overhead, poor compression ratios, lack of concurrent read support, and no efficient spatial filtering path. Cloud-native storage formats — GeoParquet, FlatGeobuf, Cloud Optimized GeoTIFF (COG), and Zarr — solve these constraints through columnar encoding, predicate pushdown, and HTTP range request optimisation.

Transitioning terabytes of spatial assets to these formats requires more than a single ogr2ogr invocation. Production-grade migration pipelines must enforce strict schema validation, preserve critical spatial metadata, handle coordinate reference system (CRS) normalisation, and guarantee idempotency across distributed compute environments.


Geospatial Data Migration Pipeline A three-stage medallion pipeline: source data (Shapefile, GeoJSON, PostGIS) enters an ingestion and validation stage; valid records proceed to transformation and encoding while invalid records route to a dead-letter queue; transformed data is published to object storage as GeoParquet, FlatGeobuf, or COG, and a manifest records row counts, bounding boxes, and checksums. Sources Shapefile GeoJSON PostGIS Ingestion & Validation CRS check Geometry validity Attribute completeness valid invalid Dead-letter queue Transformation & Encoding Schema mapping Hilbert sort WKT2 CRS normalise Publication & Verify GeoParquet FlatGeobuf / COG Checksum + row parity Manifest log

Four Core Problems Every Migration Must Solve

A migration pipeline is not a single operation — it is the composition of four distinct engineering problems. Each one has dedicated patterns and trade-offs.

1. Geometry Encoding & Format Selection

The choice of output format is a workload-aware decision, not a one-size-fits-all default. GeoParquet stores geometries as Well-Known Binary (WKB) columns inside a columnar Parquet file, making it ideal for analytical queries that scan millions of features and apply bounding-box or attribute filters. FlatGeobuf, compared head-to-head with GeoParquet, uses a packed Hilbert R-tree index that makes single-feature fetches and tile-serving requests very fast, because a client can calculate the byte offset of any spatial region without reading the entire file.

Cloud Optimized GeoTIFF embeds internal overviews and a tiling structure in the TIFF image file directory (IFD), allowing web mapping clients to request only the zoom-appropriate tile via a single HTTP range request. Zarr is preferred for multidimensional raster stacks (climate model outputs, time-series satellite imagery) where chunk-aligned access patterns dominate.

2. CRS Normalisation & Metadata Integrity

Coordinate reference system mismatches are the most dangerous silent failure in geospatial migrations. A dataset that survives ingestion, transformation, and publication with the wrong CRS produces geometrically corrupt outputs that pass row-count and checksum validation but are spatially wrong.

Schema mapping for legacy to modern formats must standardise all CRS representations to WKT2:2019 or explicit EPSG identifiers as the first transformation step, before any geometry processing. Legacy datasets regularly mix PROJ4 strings, WKT1, and EPSG codes — sometimes in the same layer, when data has been concatenated from multiple sources.

3. Spatial Indexing & Partitioning

Columnar formats unlock efficient analytical queries only when features are physically ordered by their spatial location. A naively ordered GeoParquet file — rows in their original ingestion sequence — forces a query engine to scan every row group to answer a bounding-box filter. Apply a Hilbert curve sort before serialisation: features that are geographically proximate are stored in contiguous row groups, so a spatial filter touches a small fraction of the file. See spatial partitioning with quadtree indexes for a comparison of Hilbert, Z-order, and quadtree strategies.

Partitioning complements indexing at the file level. For continental-scale datasets, partition by administrative grid or UTM zone so that a query targeting a single region reads only the relevant Parquet files. For time-series raster data, partition by year/month. Target file sizes of 256 MB–1 GB after compression to balance metadata overhead against parallelism.

4. Compression & Layout

ZSTD compression levels for geospatial data have a material effect on both storage cost and query latency. Level 3–6 is the operational sweet spot for GeoParquet vector data: compression ratios of 4–8× over raw GeoJSON with minimal CPU overhead. Raster data encoded in COG typically benefits from level 1–3, because raster tiles are already partially decorrelated by the overview pyramid and high compression levels yield diminishing returns at significant CPU cost.

Row group sizing strategies interact directly with both compression and query performance. Row groups that are too small (under 64k rows) produce excessive footer metadata; row groups that are too large (over 1 M rows) reduce filter selectivity because a predicate must match the full row group before it can be skipped.


Format & Technology Landscape

The table below maps the major output formats to their primary workload class, query pattern, and key trade-off.

Format Primary Use Case Query Pattern Spatial Index Key Trade-off
GeoParquet Analytical / OLAP Column scan + predicate pushdown Hilbert sort (embedded) Requires columnar-aware engine (DuckDB, Athena, BigQuery)
FlatGeobuf Tile serving / web mapping Single-feature or bbox fetch Packed Hilbert R-tree Not suited for columnar analytical scans
Cloud Optimized GeoTIFF Raster visualisation / tiling Zoom-level tile request IFD overview pyramid Image only; no vector attribute storage
Zarr Multidimensional raster stacks Chunk-aligned n-D slice None (chunk grid) Chunk alignment must match access pattern
GeoJSON API payloads / interchange Full-file read None Serialisation overhead scales with feature count
Shapefile Desktop GIS interchange Full-file read None 2 GB limit, 10-char fields, multi-file fragility

Pipeline Architecture: Medallion Design

A production migration pipeline follows a three-stage medallion pattern — bronze (raw ingest), silver (validated and normalised), gold (partitioned, indexed, and compressed for serving) — decoupling each stage’s compute and storage.

Bronze (Ingestion & Validation)
Source data is pulled from legacy data lakes, PostGIS instances, or network shares into an object storage landing zone without transformation. Geometry topology (self-intersections, ring orientation), CRS consistency, and attribute completeness are validated here. Records that fail validation are routed to a dead-letter queue rather than halting the batch.

Silver (Transformation & Encoding)
Data is CRS-normalised, schema-mapped to explicit Arrow types, geometry-flattened to homogeneous tables, and spatially sorted. This is the compute-intensive stage. Orchestration frameworks — Airflow, Prefect, or Dagster — manage task dependencies and resource quotas. Apache Arrow and DuckDB enable in-memory vectorised operations that bypass Python’s GIL bottleneck during batch transformation.

Gold (Publication & Verification)
Converted files are written to a staging prefix first, then checksums, row counts, and bounding box extents are compared against the source manifest. Only after parity is confirmed does an atomic rename promote the files to the final destination prefix. A manifest record is written capturing source_path, output_path, conversion_timestamp, row_count, bbox, and checksum.


Performance & Trade-off Analysis

Compression Ratio vs. Query Latency

The relationship between compression level and query latency is not monotonic. Higher compression reduces I/O time but increases decompression CPU cost. The optimal operating point depends on whether your workload is I/O-bound (slow object storage, high egress cost) or CPU-bound (in-memory analytical cluster).

Format Compression vs. GeoJSON Typical Latency (DuckDB bbox filter) Best For
GeoParquet + ZSTD-3 6–8× smaller 80–200 ms (10 M features) Analytical batch queries
GeoParquet + Snappy 4–5× smaller 60–150 ms Latency-sensitive dashboards
FlatGeobuf 3–4× smaller 5–20 ms (single bbox) Tile server, map API
COG (raster, ZSTD-1) 2–4× vs. uncompressed TIFF <50 ms (tile request) Web map raster serving

Spatial Indexing Strategy Comparison

Strategy Primary Use Case Build Cost Range Query Speed Supported Formats
Hilbert curve sort Analytical bbox scans Low (sort + write) Excellent GeoParquet
Z-order (Morton) curve Lakehouse partition pruning Low (bit-interleave) Good GeoParquet, Delta Lake
Packed Hilbert R-tree Single-feature / tile fetch Moderate (tree build) Excellent for point fetches FlatGeobuf
Quadtree partition Tile-aligned web mapping Moderate Good for tile-aligned access Any partitioned format

Format Selection Decision Path

The diagram below reduces the format decision to four questions answered in priority order: raster vs. vector, then the dominant access pattern. Resolve these before tuning compression or indexing — the format constrains every parameter downstream.

Output Format Selection Decision Path Starting from the source data type, the flow asks whether the data is raster or vector. Raster branches on dimensionality: single-band or RGB imagery becomes Cloud Optimized GeoTIFF, while multidimensional stacks become Zarr. Vector branches on access pattern: per-feature or tile fetches become FlatGeobuf, while columnar analytical scans become GeoParquet sorted on a Hilbert curve. Source dataset Raster or vector? raster vector Multi-band / time stack? Per-feature / tile fetch? no Cloud Optimized GeoTIFF yes Zarr yes FlatGeobuf no GeoParquet (Hilbert sort)

Production Implementation

The following Python snippet demonstrates the core silver-layer transformation: read a legacy Shapefile, sort features along a Hilbert curve, apply ZSTD compression, and write GeoParquet to object storage with an idempotent staging pattern.

python
# Requires: geopandas>=0.14, pyarrow>=14, lonboard>=0.4, s3fs>=2024.2
from __future__ import annotations

import hashlib
import json
import logging
import tempfile
from pathlib import Path

import geopandas as gpd
import pyarrow as pa
import pyarrow.parquet as pq
import s3fs

logger = logging.getLogger(__name__)


def spatial_sort(gdf: gpd.GeoDataFrame) -> gpd.GeoDataFrame:
    """Sort GeoDataFrame rows along a space-filling curve for spatial locality.

    Uses a Z-order (Morton) curve here for clarity; swap ``interleave`` for a
    true Hilbert encoder (e.g. ``hilbertcurve``) when contiguity matters more.
    """
    from shapely.geometry import box
    import numpy as np

    bounds = gdf.total_bounds  # [minx, miny, maxx, maxy]
    cx = (gdf.geometry.bounds["minx"] + gdf.geometry.bounds["maxx"]) / 2
    cy = (gdf.geometry.bounds["miny"] + gdf.geometry.bounds["maxy"]) / 2

    # Normalise to [0, 1] range
    x_range = bounds[2] - bounds[0] or 1.0
    y_range = bounds[3] - bounds[1] or 1.0
    nx = ((cx - bounds[0]) / x_range * 65535).astype(np.uint32)
    ny = ((cy - bounds[1]) / y_range * 65535).astype(np.uint32)

    # Interleave bits (simple Z-order as Hilbert approximation)
    def interleave(x: np.ndarray, y: np.ndarray) -> np.ndarray:
        result = np.zeros(len(x), dtype=np.uint64)
        for i in range(16):
            result |= (x >> i & 1).astype(np.uint64) << (2 * i)
            result |= (y >> i & 1).astype(np.uint64) << (2 * i + 1)
        return result

    order = interleave(nx, ny).argsort()
    return gdf.iloc[order].reset_index(drop=True)


def convert_to_geoparquet(
    source_path: str,
    dest_s3_prefix: str,
    layer: str | None = None,
    zstd_level: int = 3,
    row_group_size: int = 131_072,
) -> dict[str, object]:
    """
    Convert a legacy vector file to GeoParquet and write to S3 with idempotent staging.

    Returns a manifest dict with output_path, row_count, bbox, and sha256.
    """
    logger.info("Reading source: %s", source_path)
    gdf = gpd.read_file(source_path, layer=layer)

    if gdf.crs is None:
        raise ValueError(f"Source has no CRS: {source_path}")

    # Normalise to EPSG:4326 (WGS84) for interoperability; adjust if needed
    if gdf.crs.to_epsg() != 4326:
        logger.info("Reprojecting from %s to EPSG:4326", gdf.crs)
        gdf = gdf.to_crs(epsg=4326)

    gdf = spatial_sort(gdf)

    table = pa.Table.from_pandas(gdf, preserve_index=False)
    compression = pq.ParquetCompression("ZSTD", compression_level=zstd_level)

    with tempfile.NamedTemporaryFile(suffix=".parquet", delete=False) as tmp:
        tmp_path = Path(tmp.name)

    pq.write_table(
        table,
        str(tmp_path),
        compression=compression,
        row_group_size=row_group_size,
    )

    sha256 = hashlib.sha256(tmp_path.read_bytes()).hexdigest()
    bbox = list(gdf.total_bounds)  # [minx, miny, maxx, maxy]
    row_count = len(gdf)

    # Atomic staging write to S3
    fs = s3fs.S3FileSystem()
    stem = Path(source_path).stem
    staging_key = f"{dest_s3_prefix.rstrip('/')}/_staging/{stem}.parquet"
    final_key = f"{dest_s3_prefix.rstrip('/')}/{stem}.parquet"

    logger.info("Uploading to staging: s3://%s", staging_key)
    fs.put(str(tmp_path), staging_key)
    fs.mv(staging_key, final_key)  # atomic rename within same bucket
    tmp_path.unlink(missing_ok=True)

    manifest = {
        "source_path": source_path,
        "output_path": f"s3://{final_key}",
        "row_count": row_count,
        "bbox": bbox,
        "sha256": sha256,
        "zstd_level": zstd_level,
        "row_group_size": row_group_size,
    }
    logger.info("Manifest: %s", json.dumps(manifest))
    return manifest

Schema Governance & Compliance

Type Normalisation Rules

Legacy GIS datasets rarely conform to strict typing. Shapefiles silently truncate strings; GeoJSON’s verbose text structure has no explicit type declarations; mixed geometry types (Point, MultiPolygon, GeometryCollection) frequently coexist in a single layer.

Modern columnar formats require deterministic schemas. The key mappings to enforce:

  • stringlarge_string (Arrow) to support attribute values over 2 GB in aggregate
  • int32int64 for feature ID columns with large datasets to prevent overflow
  • mixed geometry → homogeneous tables via a geometry_type discriminator column or separate layers
  • All nullable fields must be declared explicitly; implicit nullability causes serialisation failures in strict Parquet writers

Schema mapping for legacy to modern formats provides a systematic methodology for defining transformation rules, validating type compatibility, and generating schema manifests that downstream consumers can rely on. Validate all output schemas against the Apache Parquet specification to ensure compatibility across query engines like Athena, BigQuery, and DuckDB.

CRS Management

Standardise to WKT2:2019 using PROJ’s CRSType.GEOGRAPHIC_2D or emit explicit EPSG identifiers in the GeoParquet geometry_columns metadata block. For raster COGs, populate GeoKeys GTModelTypeGeoKey, GTRasterTypeGeoKey, and ProjectedCRSGeoKey fully — partial GeoKey blocks are the most common cause of silent reprojection errors in web mapping clients.

Schema Evolution

When migrating datasets that will be updated incrementally, use column-level backwards-compatible evolution: add new columns at the end, never remove or rename existing ones. Parquet’s schema evolution rules mean that readers with an older schema simply ignore unknown columns, while writers with an older schema produce null for new required columns — test both directions in your CI validation step.


Strategic Selection Framework

Workload Dataset Size Query Pattern Recommended Format Secondary Option
Analytical OLAP Any Attribute + bbox filter GeoParquet + ZSTD-3 Parquet + Snappy
Web tile serving Any Single-bbox or feature fetch FlatGeobuf COG (raster)
Raster visualisation Large (>1 TB) Zoom-level tile COG + ZSTD-1 Zarr
Multidimensional raster Very large n-D chunk slice Zarr COG stack
Desktop GIS interchange Small (<50 MB) Full-layer read GeoPackage Shapefile
API payloads Small (<10 MB) Full deserialise GeoJSON FlatGeobuf
Streaming ingestion Continuous Append-only write GeoParquet (partitioned) Parquet + Delta Lake

Failure Handling & Operational Safeguards

Distributed batch conversions encounter transient failures: network timeouts, OOM errors, corrupted source files, and storage quota exhaustion. Resilient pipelines implement layered failure handling.

Exponential backoff and retries. Configure orchestration tasks (Airflow operators, Prefect tasks) to retry transient errors with randomised jitter to avoid thundering-herd bursts against object storage rate limits.

Dead-letter queues. Route irrecoverable records or malformed files to isolated storage for manual inspection without halting the batch. Track the DLQ path and error reason in the manifest log.

Partial success handling. When a single partition fails, isolate it and continue processing the remaining chunks. Mark the job partial_success and trigger a downstream alert so operators can re-run only the failed shard.

Fallback routing for failed migration jobs covers how to automatically switch to secondary compute pools or legacy processing paths when primary infrastructure degrades — critical for migration windows with SLA commitments.

Immutable source snapshots. Maintain write-protected snapshots of source data during migration windows. Version output directories with date-stamped prefixes so any batch can be re-run against the original inputs.


Production Deployment & Observability

Once pipelines are validated in staging, deployment moves from ad-hoc scripts to managed infrastructure-as-code. Containerise conversion workers using lightweight base images (python:3.12-slim with GDAL/OGR bindings compiled against PROJ 9+). Deploy via Kubernetes Jobs or serverless batch services (AWS Batch, GCP Cloud Run Jobs) with auto-scaling policies tied to queue depth.

Instrument pipelines with three observability layers:

  • Metrics: rows processed per second, conversion latency by format, error rate, and storage egress volume. Expose as Prometheus counters from within the conversion worker.
  • Distributed tracing: propagate a run_id and partition_id through ingestion, transformation, and publication so any failed record can be traced to its source file and offset.
  • Data quality dashboards: track schema compliance rates, geometry validity percentages, CRS consistency, and spatial coverage gaps over time. Alert when validity drops below a threshold.

Cost optimisation should run parallel to performance tuning. Leverage spot or preemptible instances for stateless conversion workers (they can be retried safely). Implement lifecycle policies to transition older output partitions to infrequent-access or cold storage tiers. Monitor compute-to-storage ratios monthly to prevent over-provisioning — conversion workloads are typically CPU-intensive for a short window, then storage-dominated long-term.


Frequently Asked Questions

When should I keep Shapefiles instead of converting?
Keep Shapefiles only when interoperating with desktop GIS tools that ingest nothing else, or when datasets are under 50 MB and shared ad hoc. At any meaningful analytical scale, the 2 GB file limit, 10-character field name restriction, and fragmented multi-file layout create enough operational pain to justify conversion. See the full breakdown in why Shapefiles fail at scale.

GeoParquet or FlatGeobuf for tile serving?
FlatGeobuf is the better fit for tile serving because its spatial index enables efficient single-feature or bounding-box fetches via HTTP range requests without deserialising the whole file. GeoParquet is optimised for columnar analytical queries across many features. The detailed GeoParquet vs FlatGeobuf performance comparison covers measured latency at different feature counts.

How large should partitions be in a migrated GeoParquet dataset?
Target individual Parquet files between 256 MB and 1 GB after compression. Files smaller than 128 MB accumulate metadata overhead that degrades query planner performance; files larger than 2 GB limit the parallelism available to distributed query engines. Run a compaction pass to merge small residual files after initial migration.

How do I preserve CRS information across format conversions?
Standardise to WKT2:2019 before writing. GeoParquet encodes CRS in its geometry_columns metadata block; COG stores it in TIFF GeoKeys. The danger zone is format-hopping through intermediate representations — PROJ strings, WKT1, and EPSG codes can describe the same datum differently and cause silent coordinate shifts. The metadata preservation guide covers the full metadata survival checklist.

What causes idempotency failures in batch migration jobs?
The most common cause is writing directly to the final output path rather than a staging prefix, then crashing mid-write. Partial files are indistinguishable from complete ones on the next run. Fix this by writing to a temp prefix, validating checksums and row counts, then performing an atomic copy or rename to the destination.

Can DuckDB query GeoParquet directly without a full download?
Yes. DuckDB’s httpfs extension fetches only the footer and the row groups relevant to a query predicate via HTTP range requests. A bounding-box filter that matches 2% of rows causes DuckDB to download roughly 2% of the file, not the whole thing — the core advantage of Hilbert-sorted columnar layout.


← Back to Home

Continue exploring