🔧 Fix DuckDB & analytic SQL Files & formats DuckDB for Spatial Data · step 20 of 20 Updated

Fixing a Geometry Column Lost on Export from DuckDB

Problem statement

The query is correct, the geometry is right in DuckDB, and the exported file is not spatial:

copy places to 'places.csv' (format csv, header true);
NAME,geom
Colonia del Sacramento,POINT (-57.836116004496425 -34.469787716602944)
Trinidad,POINT (-56.9009966 -33.5439989)

The geometry survived — as text. It has no CRS, no binary encoding, and nothing downstream will recognise the file as spatial. The same happens with JSON:

{"NAME":"Colonia del Sacramento","geom":"POINT (-57.836116004496425 -34.4697877)"}

And there is a second, sharper failure: the obvious way to write a spatial format does not work at all.

copy places to 'places.parquet' with (format gdal, driver 'Parquet');
-- Binder Error: Could not find GDAL driver: Parquet

Three exports, three different outcomes, and only one of them produces a file another GIS tool will open correctly.

Quick answer

Use the native Parquet writer for Parquet, and the GDAL driver for everything else:

-- GeoParquet: native writer, keeps the CRS in the file's `geo` metadata
copy places to 'places.parquet' (format parquet, compression zstd);

-- GeoPackage, GeoJSON, shapefile: the GDAL driver, with the SRS stated
copy places to 'places.gpkg'    with (format gdal, driver 'GPKG');
copy places to 'places.geojson' with (format gdal, driver 'GeoJSON',
                                      srs 'EPSG:4326');

Verified: the native Parquet export round-tripped through GeoPandas with crs intact, and the GPKG export came back as 3 features in EPSG:4326.

For CSV, convert the geometry deliberately so the loss is a decision rather than an accident:

copy (select name, st_x(geom) as lon, st_y(geom) as lat from places)
to 'places.csv' (format csv, header true);
Triage table of five DuckDB export forms and what each does to geometry.
The CSV export succeeds, which is exactly why it is the dangerous one.

Step-by-step solution

1. Understand what each export actually writes

Export Geometry becomes CRS kept Readable as spatial
(format csv) WKT text no only after parsing
(format json) WKT text no no
(format parquet) WKB + geo metadata yes yes, as GeoParquet
with (format gdal, driver 'GPKG') native geometry yes yes
with (format gdal, driver 'GeoJSON', srs '…') GeoJSON geometry yes, if srs given yes
with (format gdal, driver 'Parquet') — — fails: no such driver

The last row is the trap: the GDAL-driver syntax is the general answer for spatial formats, and Parquet is the one format where it does not apply because the bundled GDAL build has no Parquet driver.

2. Use the native Parquet writer, and check the metadata 

copy places to 'places.parquet' (format parquet, compression zstd);

DuckDB writes the geometry as WKB and records the CRS in the file's geo key-value metadata, which is the GeoParquet convention:

import pyarrow.parquet as pq
import json

meta = pq.read_metadata("places.parquet").metadata
print([k.decode() for k in meta])                  # ['geo']
geo = json.loads(meta[b"geo"])
print(geo["columns"]["geom"]["crs"]["id"])          # {'authority': 'EPSG', 'code': 4326}

Measured, GeoPandas read that file back with the CRS present as full PROJJSON.

3. State the SRS when the driver needs it

Some GDAL drivers take the CRS from the geometry and some need telling. Passing srs explicitly costs nothing and removes the doubt:

copy places to 'out.geojson' with (format gdal, driver 'GeoJSON', srs 'EPSG:4326');
copy places to 'out.shp'     with (format gdal, driver 'ESRI Shapefile',
                                   srs 'EPSG:27700');

This matters more in DuckDB than elsewhere because the CRS is a column-type annotation that ST_Transform discards — so by export time the engine frequently does not know it.

4. Check the driver exists before relying on it 

select short_name, long_name from st_drivers() where can_create order by 1;
ESRI Shapefile, MapInfo File, S57, DGN, Memory, CSV, GML, GPX, KML,
GeoJSON, GeoJSONSeq, OGR_GMT, …

The bundled GDAL is a smaller build than a system GDAL. If the driver you want is missing, write Parquet or GeoPackage from DuckDB and convert with another tool.

5. Make the CSV loss explicit

There are legitimate reasons to write CSV — a spreadsheet, a system that only accepts it. Do it deliberately:

-- coordinates as numbers: readable everywhere, no CRS
copy (select name, round(st_x(geom), 6) as lon, round(st_y(geom), 6) as lat,
             'EPSG:4326' as crs
      from places)
to 'places.csv' (format csv, header true);

-- or WKT, explicitly, with the CRS in a column
copy (select name, st_astext(geom) as wkt, 'EPSG:4326' as crs from places)
to 'places.csv' (format csv, header true);

Six decimal places is about 11 cm at the equator, which is beyond the accuracy of nearly all source data.

6. Verify the export by reading it back 

import geopandas as gpd

gdf = gpd.read_parquet("places.parquet")
assert gdf.crs is not None, "the CRS did not survive the export"
assert gdf.geometry.notna().all(), "geometry is missing"
print(f"{len(gdf):,} features, CRS {gdf.crs.to_string()}")

Reading the file in the tool that will consume it is the only check that matters, and it takes three lines.

Table of what DuckDB writes into a GeoParquet file and where.
Pass srs explicitly to GDAL drivers; by export time the type annotation may be gone.

Code examples

Example 1 — an export function that cannot silently lose geometry 

import duckdb


SPATIAL_FORMATS = {
    ".parquet": ("native", "format parquet, compression zstd"),
    ".gpkg":    ("gdal", "GPKG"),
    ".geojson": ("gdal", "GeoJSON"),
    ".json":    ("gdal", "GeoJSON"),
    ".shp":     ("gdal", "ESRI Shapefile"),
    ".fgb":     ("gdal", "FlatGeobuf"),
}


def export_spatial(con, relation, path, crs="EPSG:4326", overwrite=True):
    """Write a spatial format, or refuse and say why."""
    import os

    ext = os.path.splitext(path)[1].lower()
    if ext in (".csv", ".tsv"):
        raise ValueError(
            f"{ext} cannot hold geometry. Either write coordinates explicitly:\n"
            f"  copy (select *, st_x(geom) lon, st_y(geom) lat from {relation}) "
            f"to '{path}' (format csv)\n"
            f"or choose a spatial format: {', '.join(sorted(SPATIAL_FORMATS))}")

    if ext not in SPATIAL_FORMATS:
        raise ValueError(f"no known writer for {ext!r}")

    kind, spec = SPATIAL_FORMATS[ext]
    if overwrite and os.path.exists(path):
        os.remove(path)

    if kind == "native":
        con.execute(f"copy {relation} to '{path}' ({spec})")
    else:
        available = {row[0] for row in con.execute(
            "select short_name from st_drivers() where can_create").fetchall()}
        if spec not in available:
            raise ValueError(
                f"the bundled GDAL has no '{spec}' driver. "
                f"Write .parquet or .gpkg and convert with another tool.")
        con.execute(f"copy {relation} to '{path}' "
                    f"with (format gdal, driver '{spec}', srs '{crs}')")

    print(f"wrote {path} ({os.path.getsize(path) / 1e6:,.2f} MB) as {kind}/{spec}")
    return path

Example 2 — verifying the round trip 

def verify_export(path, expected_rows=None, expected_crs=None):
    """Read the file back the way its consumer will."""
    import geopandas as gpd
    import os

    reader = gpd.read_parquet if path.endswith(".parquet") else gpd.read_file
    gdf = reader(path)

    problems = []
    if gdf.geometry.isna().any():
        problems.append(f"{gdf.geometry.isna().sum()} null geometries")
    if gdf.crs is None:
        problems.append("no CRS — the export lost it")
    elif expected_crs and gdf.crs.to_string() != expected_crs:
        problems.append(f"CRS is {gdf.crs.to_string()}, expected {expected_crs}")
    if expected_rows and len(gdf) != expected_rows:
        problems.append(f"{len(gdf):,} rows, expected {expected_rows:,}")

    for problem in problems:
        print(f"  ! {problem}")
    if not problems:
        print(f"  {os.path.basename(path)}: {len(gdf):,} features, "
              f"CRS {gdf.crs.to_string()}, geometry intact")
    return not problems

Example 3 — turning a WKT CSV back into geometry 

def recover_from_wkt_csv(con, csv_path, crs="EPSG:4326", wkt_column="geom"):
    """Rescue a CSV that was exported with geometry as text."""
    con.execute("load spatial")
    return con.execute(f"""
        select * exclude {wkt_column},
               st_geomfromtext({wkt_column}) as geom
        from read_csv('{csv_path}')
    """).df(), crs

The CRS has to be supplied from outside, because the CSV never carried it. That is the whole cost of the accidental export: a value that has to be remembered rather than read.

Explanation

Why CSV writes WKT rather than failing

CSV has one type: text. Faced with a geometry column, DuckDB does the most useful thing available and writes its textual representation, which is at least reversible.

The alternative — refusing — would break the legitimate case of exporting WKT deliberately. So the behaviour is defensible and it is silent, which is why an export that "worked" produces a file no GIS tool recognises.

Why the Parquet GDAL driver does not exist here

GDAL has a Parquet driver, but it depends on Arrow and is not compiled into the static GDAL bundled with the DuckDB spatial extension. Hence Binder Error: Could not find GDAL driver: Parquet.

This is not a limitation in practice, because DuckDB's native Parquet writer is better for the job: it is faster, it writes proper GeoParquet metadata, and it is the same writer used for every other Parquet export.

The trap is only that the GDAL syntax is the general answer for spatial formats and Parquet is the exception.

Why the CRS is so easily lost

DuckDB carries the CRS on the geometry column's type, and the annotation does not survive most operations — ST_Transform returns a bare GEOMETRY, and so does any constructed geometry.

By the time a query result is exported, the engine frequently has no CRS to write. Passing srs explicitly to the GDAL driver, or knowing that the native Parquet writer preserves what the type still holds, is the difference between a file that documents its coordinate system and one that does not.

Why verifying the round trip is worth three lines

An export is only correct with respect to the tool that will read it. GeoPandas, QGIS, PostGIS and a web map each have different tolerances for a missing CRS, an unusual geometry type or a mixed-geometry column.

Reading the file back with the consumer's library, checking the row count, the CRS and the geometry validity, catches all of that at the point of writing rather than in somebody else's session a week later.

Checklist of four export verification steps and the anti-pattern of trusting the write.
Catching it here beats catching it in somebody else’s session next week.

Edge cases or notes

  • format parquet is native; format gdal, driver 'Parquet' does not exist in the bundled build.
  • COPY ... TO 'x.csv' writes WKT — reversible, but with no CRS.
  • Pass srs to GDAL drivers explicitly. The type annotation may already be gone.
  • Shapefile truncates field names to 10 characters and has a 2 GB limit.
  • GeoJSON should be EPSG:4326 by specification; writing another CRS is legal and confusing.
  • Mixed geometry types are fine in DuckDB and GeoPackage, and rejected by shapefile.
  • st_drivers() lists what this build can write — check before relying on a format.
  • Round coordinates when writing text. Six decimal places is about 11 cm.

FAQ

Why is my exported CSV's geometry just text?

Because CSV has no geometry type. DuckDB writes WKT, which is reversible but carries no CRS and is not recognised as spatial by any GIS tool.

How do I write GeoParquet from DuckDB?

copy tbl to 'x.parquet' (format parquet) — the native writer. It stores WKB plus the geo metadata key, and GeoPandas reads the CRS back.

Why does format gdal, driver 'Parquet' fail?

The bundled GDAL build has no Parquet driver. Use the native Parquet writer instead, which is faster and writes proper GeoParquet metadata anyway.

How do I write a GeoPackage?

copy tbl to 'x.gpkg' with (format gdal, driver 'GPKG'). Verified: it round-trips through GeoPandas with the geometry and EPSG:4326 intact.

Why did my export lose the CRS?

DuckDB tracks CRS on the column type and drops it in most operations, including ST_Transform. Pass srs explicitly to the GDAL driver.

How do I check the export worked?

Read it back with the library that will consume it and assert the row count, the CRS and that no geometry is null. Three lines, at the point of writing.