Logging Feature Counts and EPSG Codes per Task

Two fields on every stage’s log record — the feature count and the CRS — are enough to reconstruct what a vector pipeline did to a dataset and where it went wrong. The count exposes silent drops: a join that halved the features, a filter that emptied a layer, a repair that discarded slivers. The CRS exposes the transform that happened where none was intended, or the one that did not happen where it was. Both are one cheap read per stage, and together they turn “the output looks wrong” into a table showing exactly which step changed what.

When to Use This Pattern

  • A vector chain has several stages, so a discrepancy at the end could have originated at any of them.
  • Features have gone missing before and nobody could say which step lost them.
  • Multiple CRSs are in play, which they are whenever a source and a target differ.
  • Structured logging is in place, since these are fields rather than sentences.

Complete Working Example

A tiny helper reads both facts from a layer’s header, and every stage logs them on entry and exit.

from __future__ import annotations

from dataclasses import dataclass
from pathlib import Path

import fiona
import structlog
from pyproj import CRS

log = structlog.get_logger()


@dataclass(frozen=True)
class LayerFacts:
    """The two numbers that describe a vector layer's state."""

    features: int
    epsg: int | None
    crs_text: str
    geometry_type: str

    @classmethod
    def read(cls, path: Path, layer: str = "features") -> "LayerFacts":
        # Header only: len(src) uses the driver's feature count where available,
        # which is O(1) for GeoPackage and FlatGeobuf.
        with fiona.open(path, layer=layer) as src:
            crs = CRS.from_user_input(src.crs) if src.crs else None
            return cls(
                features=len(src),
                epsg=crs.to_epsg() if crs else None,
                crs_text=crs.to_string() if crs else "undeclared",
                geometry_type=src.schema["geometry"],
            )


def log_stage(stage: str, before: LayerFacts, after: LayerFacts) -> None:
    """One record per stage, with the deltas computed rather than implied."""
    delta = after.features - before.features
    log.info(
        "stage_completed",
        stage=stage,
        features_in=before.features,
        features_out=after.features,
        features_delta=delta,
        # A ratio is what an alert or a dashboard can threshold on; a raw
        # delta means nothing without knowing the size of the input.
        features_ratio=round(after.features / before.features, 4) if before.features else None,
        epsg_in=before.epsg,
        epsg_out=after.epsg,
        crs_changed=before.epsg != after.epsg,
        geometry_type_changed=before.geometry_type != after.geometry_type,
    )

Wiring it into a chain is one line per stage, and the assertion is what turns a log into a guard:

def run_chain(source: Path, scratch: Path, target_crs: str = "EPSG:25832") -> Path:
    raw = extract(source, scratch / "01_raw.gpkg")
    facts_raw = LayerFacts.read(raw)
    log.info("chain_started", features=facts_raw.features, epsg=facts_raw.epsg)

    clean = repair(raw, scratch / "02_valid.gpkg")
    facts_clean = LayerFacts.read(clean)
    log_stage("repair", facts_raw, facts_clean)

    projected = reproject(clean, target_crs, scratch / "03_projected.gpkg")
    facts_proj = LayerFacts.read(projected)
    log_stage("reproject", facts_clean, facts_proj)
    # A reprojection must not change the feature count. Ever.
    if facts_proj.features != facts_clean.features:
        raise ValueError(
            f"reprojection changed the feature count: "
            f"{facts_clean.features}{facts_proj.features}"
        )

    joined = join_attributes(projected, scratch / "04_joined.gpkg")
    facts_join = LayerFacts.read(joined)
    log_stage("join", facts_proj, facts_join)
    return joined
Where the features wentExtract produces four hundred and twelve thousand features. Repair drops nine. Reprojection changes nothing. The join drops a hundred and eighty-four thousand, which is where the problem is.features after each stageextract412 004repair411 995reproject411 995join228 331 — ratio 0.554An inner join on a key with mismatched types. The output loaded successfully and looked plausible.
Without a per-stage count the only visible fact is that the final table has 228 331 rows, which is a number nobody has a reason to distrust.

The invariants differ per stage, and writing them down is most of the value here. A reprojection must never change the feature count — if it does, something reprojected into an invalid area and the driver dropped rows. A repair may reduce the count and must never increase it. A join may reduce it substantially and legitimately, which is why it gets a threshold rather than an assertion. And an extract sets the baseline, so it has no invariant at all, only a plausibility check that it is not zero. Four stages, four different rules, and each one catches a failure the others cannot.

One invariant per stageExtract sets the baseline and is only checked for being non-zero. Repair may reduce the count. Reprojection must not change it at all. The join may reduce it and is monitored against its own history.stageinvariantenforcementextractcount > 0plausibilityrepairout ≤ in, drop < 1%thresholdreprojectout == in, alwaysraisejoinratio near its own normmonitor
Only one stage gets a hard assertion, and it is the one where any change at all is definitionally a bug. Applying the same strictness everywhere would fail on every legitimate join.

Parameter & Option Reference

Field Source Spatial notes
features_in / features_out len(src) O(1) for GeoPackage and FlatGeobuf; a full scan for GeoJSON, which is a reason to avoid it as an intermediate.
features_ratio computed What a threshold can be set on. A raw delta is uninterpretable without the input size.
epsg_in / epsg_out CRS.to_epsg() None for a custom or undeclared CRS, which is itself the signal.
crs_text CRS.to_string() Keeps the full description when to_epsg() returns nothing.
crs_changed boolean Makes “which stage transformed?” a filter rather than a comparison across records.
geometry_type_changed boolean Repairs change polygons into multipolygons and collections; the load will care.
assertion per stage The stages that must not change the count should say so in code, not only in a log.

Verification & Testing

The tests worth writing are about the invariants each stage promises, which the fields make checkable.

def test_reprojection_preserves_the_feature_count(sample_layer, tmp_path) -> None:
    before = LayerFacts.read(sample_layer)
    out = reproject(sample_layer, "EPSG:3857", tmp_path / "p.gpkg")
    after = LayerFacts.read(out)
    assert after.features == before.features
    assert after.epsg == 3857 and before.epsg != 3857


def test_repair_may_drop_but_must_not_add(invalid_layer, tmp_path) -> None:
    before = LayerFacts.read(invalid_layer)
    after = LayerFacts.read(repair(invalid_layer, tmp_path / "r.gpkg"))
    assert after.features <= before.features
    assert (before.features - after.features) / before.features < 0.01


def test_join_ratio_is_logged(sample_layer, tmp_path, log_capture) -> None:
    run_chain(sample_layer, tmp_path)
    join_record = next(r for r in log_capture if r.get("stage") == "join")
    assert "features_ratio" in join_record
    assert 0.0 <= join_record["features_ratio"] <= 1.0

Once the records exist, the query that catches a silent drop is a single aggregation, and it is worth putting on the dashboard rather than running only after someone complains:

-- Stages whose feature ratio moved from its own norm, over the last week.
SELECT stage,
       round(avg(features_ratio)::numeric, 4)      AS avg_ratio,
       round(stddev(features_ratio)::numeric, 4)   AS sd,
       min(features_ratio)                          AS worst
FROM   stage_log
WHERE  logged_at > now() - interval '7 days'
GROUP  BY stage
HAVING stddev(features_ratio) > 0.02
ORDER  BY sd DESC;
The ratio is stable until it is notThe join stage's feature ratio sits at just under one for eleven days, then drops sharply to zero point five five when an upstream key format changed.1.00.5join stage, features_ratio, two weeksthe source switched to zero-padded key stringsEvery run after the change succeeded, published, and contained 45% fewer features than it should have.
A ratio with a stable history is one of the cheapest anomaly detectors available, and it requires no model — only that the number was recorded every time.

Two weeks of history is enough to make this detector useful, which is worth knowing because it means the value arrives quickly. There is no model to train and no baseline to configure: the ratio either has a narrow historical distribution, in which case a departure from it is a signal, or it does not, in which case the stage is genuinely variable and the field is still useful for post-hoc investigation. Either way the cost is one number per stage per run, and the first time it catches a silent join failure it repays every log line it has ever produced.

Common Pitfalls

  • Logging a count without the input count. “228 331 features” is not actionable. The ratio is what makes a drop visible without knowing what the input was supposed to be.
  • Assuming len() is cheap. It is O(1) for GeoPackage and FlatGeobuf and a full scan for GeoJSON and CSV. If a stage’s intermediate is one of the latter, the count costs a pass — which is another reason those formats are poor intermediates.
  • Recording only the EPSG code. A custom or compound CRS returns None from to_epsg(), so a pipeline logging only the code loses the information entirely. Keep the full string alongside it.
  • Not asserting the invariants. A reprojection that changes the feature count is always a bug, and a log record that shows it is only useful if someone reads it. Where the invariant is absolute, raise.
  • Ignoring the geometry type. A repair turning polygons into a geometry collection is invisible in the count and fatal at the load. One boolean field makes it visible at the stage that caused it.
  • Logging at the wrong granularity. These are per-stage facts, not per-feature. A record per feature is four hundred thousand lines and answers nothing the aggregate does not.

Frequently Asked Questions

Should the counts be metrics instead?

Both, for different questions. A metric answers “is the join ratio trending down across all runs”, which is a monitoring question. The log record answers “what happened to this dataset in this run”, which is a debugging question and needs the per-run identity that metric labels cannot carry. The two are cheap enough that choosing between them is a false economy.

What about counting geometries rather than features?

They differ when a feature’s geometry is null, which is more common than it should be in vector deliveries. If your loader treats a null geometry as a failure, count both and log the difference — a dataset where the two diverge is telling you something about the source that no other field will.

How do I set a threshold on the ratio?

From its own history, per stage, as in the query above. A join whose ratio has been 0.998 for a month should alert at 0.95; one that legitimately varies between 0.6 and 0.9 needs a wider band or a different signal. Stages that must not change the count at all get an assertion instead of a threshold, which is stricter and simpler.

Does this replace the contract in the ETL chain?

No — it complements it. The contract is enforced at the handoff and fails the run; these fields are recorded whether or not anything failed, which is what makes the historical comparison possible. In practice the contract carries the same values, and logging them is one line at the point where they are already in hand.

What should be logged for raster stages?

The same shape with different fields: pixel dimensions, band count, nodata value and CRS instead of feature count. The invariants differ — a warp legitimately changes dimensions — but the principle is identical: record what the stage received and what it produced, and the deltas become queryable.

Structured Logging for Geospatial Flows