GeoAlchemy2 vs Raw WKB: When to Use Each

Knowing when to lean on the Type Coercion and Serialization machinery in GeoAlchemy2 versus when to bypass it and handle raw Well-Known Binary directly is the difference between an AIS backend that is both ergonomic and fast and one that is only one of the two. GeoAlchemy2 maps every geometry column to a WKBElement and gives you object identity, relationships, and safe conversion; raw WKB through psycopg gives you a firehose with no per-object overhead. This page draws the line with concrete ais_positions and port_zones examples in SRID 4326, so each path is used where it actually pays off.

Why the naive approach fails

Two opposite mistakes appear in real codebases. The first forces the ORM through a bulk export, paying WKBElement construction for millions of rows nobody will ever treat as objects:

python
# SLOW at scale: builds a WKBElement per row for a pure data dump
def export_positions_orm(session, day):
    rows = session.scalars(
        select(AisPosition).where(AisPosition.recorded_at >= day)
    ).all()                                  # millions of ORM objects
    return [
        {"mmsi": r.mmsi, "wkt": to_shape(r.geom).wkt}  # decode each one
        for r in rows
    ]

The second overreacts and hand-parses WKB with struct on an ordinary read path, reinventing what Shapely already does correctly — and getting the byte order or the EWKB SRID prefix subtly wrong:

python
# FRAGILE: hand-rolled WKB parsing that ignores the EWKB SRID flag
import struct

def parse_point_wrong(wkb: bytes):
    byte_order = wkb[0]
    # Assumes plain WKB; PostGIS emits EWKB with an SRID block that this
    # offset math skips, so x/y land on the wrong bytes for SRID geometries.
    x, y = struct.unpack_from("<dd", wkb, 5)
    return x, y

The right rule is simpler than either instinct: use GeoAlchemy2 wherever the code is genuinely object-shaped, and drop to raw WKB only on wide bulk paths where per-object cost dominates.

Production-ready implementation

The single copy-paste block below implements both paths against the same schema and shows exactly where each belongs. The ORM path serves a single vessel’s recent track as GeoJSON; the raw path bulk-loads and bulk-exports millions of points with psycopg and COPY.

python
import psycopg
from datetime import datetime, timezone
from sqlalchemy import BigInteger, DateTime, Integer, func, select
from sqlalchemy.orm import DeclarativeBase, Mapped, mapped_column, Session
from sqlalchemy.dialects.postgresql import JSONB
from geoalchemy2 import Geometry
from geoalchemy2.shape import to_shape, from_shape
from shapely import wkb as shp_wkb
from shapely.geometry import Point


class Base(DeclarativeBase):
    pass


class AisPosition(Base):
    __tablename__ = "ais_positions"
    id:   Mapped[int] = mapped_column(Integer, primary_key=True)
    mmsi: Mapped[int] = mapped_column(BigInteger, index=True)
    recorded_at: Mapped[datetime] = mapped_column(DateTime(timezone=True))
    geom: Mapped[object] = mapped_column(
        Geometry(geometry_type="POINT", srid=4326, spatial_index=True)
    )


# ── PATH A — GeoAlchemy2 for ORM-shaped work ─────────────────────────────
# Use when you want object identity, a small result set, or geometric logic.

def recent_track_geojson(session: Session, mmsi: int, limit: int = 500):
    """Serialize in the DB: ST_AsGeoJSON ships text, no per-row decode."""
    stmt = (
        select(
            AisPosition.recorded_at,
            func.ST_AsGeoJSON(AisPosition.geom).cast(JSONB).label("geometry"),
        )
        .where(AisPosition.mmsi == mmsi)
        .order_by(AisPosition.recorded_at.desc())
        .limit(limit)
    )
    return [
        {"t": r.recorded_at.isoformat(), "geometry": r.geometry}
        for r in session.execute(stmt)
    ]


def snap_to_point(session: Session, mmsi: int):
    """Use to_shape when you need a Shapely object for client-side geometry."""
    pos = session.scalars(
        select(AisPosition)
        .where(AisPosition.mmsi == mmsi)
        .order_by(AisPosition.recorded_at.desc())
        .limit(1)
    ).first()
    pt = to_shape(pos.geom)              # WKBElement -> shapely Point (SRID 4326)
    return pt.x, pt.y


def insert_one(session: Session, mmsi: int, lon: float, lat: float):
    """ORM insert: from_shape attaches the SRID to the WKBElement."""
    session.add(AisPosition(
        mmsi=mmsi,
        recorded_at=datetime.now(timezone.utc),
        geom=from_shape(Point(lon, lat), srid=4326),
    ))
    session.commit()


# ── PATH B — Raw WKB with psycopg for bulk firehose paths ────────────────
# Use when moving millions of geometries; skip WKBElement entirely.

DSN = "postgresql://ais:secret@db.internal:5432/ais_tracking"

def bulk_load_positions(rows: list[tuple]) -> None:
    """rows: (mmsi, epoch, lon, lat). COPY raw WKB straight into the table."""
    with psycopg.connect(DSN) as conn:
        with conn.cursor() as cur:
            with cur.copy(
                "COPY ais_positions (mmsi, recorded_at, geom) "
                "FROM STDIN (FORMAT BINARY)"
            ) as copy:
                copy.set_types(["int8", "timestamptz", "geometry"])
                for mmsi, epoch, lon, lat in rows:
                    ts = datetime.fromtimestamp(epoch, tz=timezone.utc)
                    geom = shp_wkb.dumps(Point(lon, lat), srid=4326, hex=False)
                    copy.write_row((mmsi, ts, geom))
        conn.commit()


def bulk_export_coords(mmsi: int):
    """Stream ST_AsBinary and decode with shapely — no ORM, no WKBElement."""
    with psycopg.connect(DSN) as conn:
        with conn.cursor(name="export", binary=True) as cur:  # server-side
            cur.itersize = 5000
            cur.execute(
                "SELECT ST_AsBinary(geom) FROM ais_positions WHERE mmsi = %s",
                (mmsi,),
            )
            for (raw,) in cur:
                pt = shp_wkb.loads(bytes(raw))   # correct EWKB/WKB decode
                yield pt.x, pt.y

The two paths are complementary, not competing. Path A gives the API layer clean GeoJSON and Shapely objects for the handful of vessels a request touches. Path B moves the daily firehose of AIS points without ever constructing a WKBElement. Notice that even Path B leans on Shapely’s wkb.loads rather than manual struct parsing — Shapely already handles PostGIS EWKB, including the SRID prefix, correctly.

Tradeoffs at a glance

Concern GeoAlchemy2 (WKBElement) Raw WKB / WKT with psycopg
Object identity & relationships Yes No
Per-row Python overhead One WKBElement per row None until you decode
Serialization to_shape or ST_AsGeoJSON ST_AsBinary / ST_AsGeoJSON you handle
Bulk insert throughput Moderate (ORM unit of work) High (COPY, executemany)
SRID safety Enforced by the Geometry type You attach SRID in SQL or EWKB
Best for API reads, writes, geometric logic Exports, ingestion, ETL firehoses

The decision collapses to one question: does this code path treat geometries as domain objects, or as bytes flowing past? Object-shaped work belongs in GeoAlchemy2; byte-shaped work belongs in raw WKB.

Configuration and tuning knobs

  • ST_AsGeoJSON precision — pass a coordinate-precision argument, e.g. ST_AsGeoJSON(geom, 6), to cap AIS coordinates at ~0.1 m resolution and shrink payloads; the default 9 digits is wasteful for lon/lat.
  • itersize on the raw cursor — the server-side cursor above fetches 5000 rows per round-trip. Raise it for throughput, lower it to cap memory; it is the raw-path analogue of yield_per.
  • COPY (FORMAT BINARY) vs text — binary COPY avoids hex encoding of WKB and is the fastest bulk-load path; set_types([..., "geometry"]) lets psycopg send the EWKB bytes directly.
  • from_shape(..., srid=4326) — always pass the SRID on the ORM write path; omitting it yields SRID 0 and a mismatch error against a 4326 column.
  • JSONB cast — casting ST_AsGeoJSON to JSONB in the query lets the driver hand you a ready dict, skipping a json.loads in Python for feature-collection endpoints.

Verification steps

Prove that both paths encode and decode the same coordinate identically before trusting either in production:

python
from shapely import wkb as shp_wkb
from shapely.geometry import Point
from geoalchemy2.shape import from_shape, to_shape

lon, lat = 4.4792, 51.9080     # Port of Rotterdam, SRID 4326

# Raw WKB round-trip
raw = shp_wkb.dumps(Point(lon, lat), srid=4326, hex=False)
assert shp_wkb.loads(raw).x == lon

# GeoAlchemy2 round-trip
element = from_shape(Point(lon, lat), srid=4326)
assert to_shape(element).y == lat
print("both paths agree:", (lon, lat))

Then confirm the database serialization matches the SRID you expect:

sql
-- Both serializers must report SRID 4326 and identical coordinates
SELECT ST_SRID(geom)            AS srid,
       ST_AsGeoJSON(geom, 6)    AS geojson,
       encode(ST_AsBinary(geom), 'hex') AS wkb_hex
FROM   ais_positions
WHERE  mmsi = 244660000
ORDER  BY recorded_at DESC
LIMIT  1;
-- Expect srid = 4326 and the geojson coordinates matching the decoded WKB.

If the GeoJSON and the decoded WKB disagree, an axis-order or SRID problem is present — recheck that the column is declared geometry(Point,4326) and that inserts attach the SRID, rather than assuming the serializer is at fault.

Gotchas checklist

  • Forcing the ORM through a bulk export materialises a WKBElement per row for data nobody treats as an object. Switch that path to raw WKB with COPY or a server-side cursor.
  • Hand-parsing WKB with struct ignores the EWKB SRID prefix PostGIS emits and lands coordinates on the wrong bytes. Use shapely.wkb.loads, which handles EWKB correctly.
  • Serializing large collections with to_shape per row is slower than pushing ST_AsGeoJSON into the query, where PostgreSQL builds the JSON once.
  • Dropping the SRID on a raw write stores SRID 0 and breaks later ST_DWithin or ST_Transform. Embed it via wkb.dumps(..., srid=4326) or ST_GeomFromWKB(%s, 4326).
  • Over-precise GeoJSON ships nine decimal places for AIS lon/lat that is accurate to a few metres. Pass an explicit precision to ST_AsGeoJSON to cut payload size without losing meaningful accuracy.