Catalog

Catalog construction has two separable concerns:

1. Fetch — query a STAC endpoint (CMR-STAC) for what / when / where → a Catalog (a stac-geoparquet table of granule metadata, reusable across many grids).
2. Shard map — take a Catalog plus an output grid → a ShardMap: the work-distribution manifest mapping shard keys to granules.

The CLI chains them, building the output grid from the same pipeline config the aggregator uses, so a shard map can never be built against a different grid than the run (enforced at run time via grid.signature()).

Building a shard map (CLI)

# HEALPix grid from atl06.yaml, an ICESat-2 cycle, Antarctic polygon:
python -m zagg.catalog --config atl06.yaml --short-name ATL06 --cycle 22 \
    --polygon antarctica.geojson

# Rectilinear (UTM) grid from a config, explicit dates, over a bbox:
python -m zagg.catalog --config serc_atl03.yaml --short-name ATL03 \
    --start-date 2025-01-01 --end-date 2025-12-31 \
    --bbox=-76.62107,38.84504,-76.50583,38.93512

# Persist the fetched Catalog too (reusable for other grids):
python -m zagg.catalog --config atl06.yaml --short-name ATL06 --cycle 22 \
    --polygon antarctica.geojson --catalog-out cycle22.parquet

--polygon drives both the CMR query bbox and the coverage mask; --bbox gives the query box directly (coverage falls back to that rectangle). The geometry backend (--backend) defaults to auto: exact-S2 spherely if the spherely fork (with SpatialIndex) is installed separately, else mortie (HEALPix) / shapely (rectilinear).

Endpoint selection (S3 vs HTTPS) is not made here — each granule record keeps both hrefs, and the aggregator picks one at run time via data_source.driver.

The footprint_cells convention column

The per-granule footprint→cells cover is identical for every shard map built against the same catalog, so a catalog can carry it (issue #396):

# index while fetching, so the saved catalog ships pre-indexed
python -m zagg.catalog --config atl06.yaml --short-name ATL06 --cycle 22 \
    --polygon antarctica.geojson --catalog-out cycle22.parquet \
    --index-footprints 9

cat = Catalog.from_geoparquet("cycle22.parquet").index_footprints(9)
cat.to_geoparquet("cycle22_indexed.parquet")

footprint_cells is a ragged large_list of morton MOC words — mortie's morton_index extension type, the same typed-morton convention as morton arrow — one entry per table row. This is a zagg convention column, not a stac-geoparquet upstream field: an indexed catalog is still an ordinary stac-geoparquet file, and any other reader ignores the extra column. It is not part of the store specification, which governs the Zarr store the aggregator writes, not the catalog it reads.

It cannot go stale. The column rides in the same file as the geometry it was covered from, so there is no second artifact to keep in sync and no version skew to detect: subsetting the catalog (filter_bbox) carries the column with the rows, and rewriting the geometry means rewriting the file.

ShardMap.build against an indexed catalog does no geometry at all — each granule's stored MOC is intersected (moc_and) with the AOI's own shard MOC — and records metadata["footprint_cells"] = True when it took that path, False when the catalog is indexed but the build took the geometry path anyway. The catalog's own footprint_cells_order rides into the manifest either way and only says the column exists, so read footprint_cells for the verdict; a manifest with neither key came from a catalog that was never indexed. The fast path engages only where the stored cover is exactly what the build asked for: the mortie backend, footprint="swath", and no caller-pinned mortie_order. Anything else takes the geometry path unchanged, so an exact-S2 spherely run is never silently swapped for a MOC one.

The assignment it produces is the geometry path's, with one intended exception: a MultiPolygon footprint. index_footprints covers the union of the rings in each blob, while granule_records reads only the largest part's exterior ring — so for a multi-part granule the column is a superset and the index can place it in shards the geometry path misses. No CMR ATL03/06 footprint is multi-part today; antimeridian-split STAC footprints are the natural producer.

What indexing buys, after #445

An unindexed mortie swath build on a HEALPix grid now covers the geometry column itself — the same from_wkb cover, at the grid's parent_order, thrown away when the build ends — and runs the same intersection (issue #445). The column no longer buys a different code path; it buys skipping the cover, which is most of a first build. Index when many builds share one catalog; skip it for a one-shot build and pay the cover once, in memory. The MultiPolygon superset above and the null-geometry refusal are properties of that shared cover, so the unindexed path has them too. footprint="beams", the spherely backend, rectilinear grids and paired builds still cover from decoded records, and the footprint_cells metadata verdict still means "the stored column answered this build" — an unindexed build that covered its own is not that.

Choosing the order

Index at the grid's parent_order. A MOC answers any shard order coarser than or equal to its own; a build against a grid whose parent_order is finer than the column is refused outright, because answering it would refine every cell onto all its descendants and put ~every granule in ~every shard (issue #92). The other end is mortie's order-18 coverage cap, which index_footprints refuses above — the same bound ShardMap.build clamps its own derived MOC order to.

Going finer than you need is expensive in both directions — words per granule roughly double per order:

order	words/granule (88S)	column, 35,639 granules	index pass	build (35,639 granules)
9	288	17 MB parquet	2.6 s	1.6 s (vs 4.8 s from geometry)
13	7,207	560 MB parquet	53.3 s	22.9 s (vs 59.3 s from geometry)

At the shard order the index pays for itself on the first build. At the chunk order it does not: the pass costs about as much as the build it replaces, and the column is two orders of magnitude larger for resolution the shard cells cannot see.

Fetch

zagg.catalog.sources.Query dataclass

A spatiotemporal metadata query: what, when, where.

Parameters:

• short_name (str) –

  Product short name (e.g. "ATL03").

• version (str) –

  Product version (e.g. "007").

• start_date (str) –

  Inclusive date bounds, YYYY-MM-DD.

• end_date (str) –

  Inclusive date bounds, YYYY-MM-DD.

• region (tuple or str) –

  Either a (lon_min, lat_min, lon_max, lat_max) bbox or a path to a GeoJSON file (its bounding box is used for the STAC query).

• provider (str, default: 'NSIDC_CPRD' ) –

  CMR provider / STAC sub-catalog. Default "NSIDC_CPRD".

collection property

collection: str

CMR-STAC collection id, {short_name}_{version}.

zagg.catalog.sources.CMRSource

Fetch granule metadata from NASA's CMR-STAC endpoint.

Parameters:

• provider (str, default: None ) –

  Overrides the query provider for the STAC sub-catalog URL.

• timeout (int, default: 60 ) –

  Per-request timeout in seconds.

fetch

fetch(
    query: Query, *, preserve_thumbnails: bool = False, limit: int = 2000
) -> "Catalog"

Run query against CMR-STAC and return a Catalog.

Parameters:

• query (Query) –

  What/when/where to fetch.

• preserve_thumbnails (bool, default: False ) –

  Keep thumbnail_*/browse assets (default drops them).

• limit (int, default: 2000 ) –

  Page size hint; CMR clamps it and paging follows rel=next.

Returns:

• Catalog –

zagg.catalog.sources.Catalog dataclass

Fetched granule metadata: a stac-geoparquet table + provenance.

Reusable across many ShardMap builds. Endpoint-neutral -- each granule carries both its S3 and HTTPS .h5 hrefs.

Parameters:

• table (Table) –

  stac-geoparquet table (one row per granule).

• metadata (dict, default: dict() ) –

  Query provenance (product, version, bbox, dates, ...).

cover_footprints

cover_footprints(order: int) -> tuple

Cover every footprint at order and return the MOCs -- persisting nothing.

The shared core of :meth:index_footprints and of ShardMap.build's unindexed mortie path (issue #445): one mortie.arrow.from_wkb call over the geometry column, screened by :meth:granule_row_mask, in the same row-aligned ragged layout :meth:footprint_cells returns. The only difference between the two callers is what they do with it -- index_footprints writes it into a column, build intersects it and drops it -- so they cannot drift on the cover itself. "Indexed" is therefore a statement about persistence, not about how a build assigns.

Parameters:

• order (int) –

  HEALPix order to cover at. Unvalidated here: callers own the bound (index_footprints refuses above MORTIE_MOC_ORDER_CAP, build covers at the grid's parent_order).

Returns:

• tuple –

  (values, offsets, rows). values is every covered row's morton words concatenated (uint64); offsets are arrow list offsets over the whole table (int64, length table.num_rows + 1), so row i's MOC is values[offsets[i]:offsets[i + 1]] and a screened-out row carries a zero-length run; rows is the covered table rows (int64, ascending) -- np.flatnonzero of the screen, which is exactly the order :meth:granule_records emits.

Notes

Row-aligned rather than compact so both callers index it by table row: the column index_footprints writes is one entry per row by contract, and build's intersection is positional against the same table. The empty runs cost 8 bytes each.

Memory is :meth:index_footprints' documented posture, because it is this method: from_wkb bounds its own peak at roughly the result size, and the screen's WKB copy plus live shapely objects is the term it does not bound (~1 GB over the parquet read on the 555,867-row ATL03 clone -- a peak, not a leak, issue #429).

filter_bbox

filter_bbox(boxes) -> 'Catalog'

Subset to granules whose bbox overlaps any of boxes (superset cut).

A columnar prefilter over the stac-geoparquet bbox column — no geometry runs here. The exact footprint-vs-shard intersection happens in ShardMap.build (mortie / spherely backends); this cut exists so a large catalog (e.g. a full-mission clone) hands the exact backend thousands of candidates instead of the whole archive.

Parameters:

• boxes (tuple or list of tuple) –

  One (lon_min, lat_min, lon_max, lat_max) box, or a list of them — one per scattered AOI part, so a multipart AOI is cut per part rather than by its (possibly continental) union box. Pair with grid.coverage_bbox for shard-complete cuts.

Returns:

• Catalog –

  New catalog with the subset table; metadata carried verbatim.

footprint_cells

footprint_cells()

The stored footprint index as (values, offsets, order), or None.

None when the catalog was never indexed (:meth:index_footprints), which is what keeps the build fast path opt-in: an ordinary catalog simply has no column and takes the geometry path.

Returns:

• tuple or None –

  (values, offsets, order) where values is the concatenated uint64 morton words of every table row (not every granule_records entry -- build aligns the two by row position, via np.flatnonzero(granule_row_mask()), since the screened rows carry a zero-length run rather than being absent) and offsets are arrow list offsets into it, so row i's MOC is values[offsets[i]:offsets[i + 1]].

from_geoparquet classmethod

from_geoparquet(path: str) -> 'Catalog'

Load a catalog from a stac-geoparquet file (CMR or user-supplied).

granule_records

granule_records() -> list[dict]

Decode the table into per-granule dicts for ShardMap building.

Returns:

• list of dict –

  Each: {"id", "s3", "https", "lats", "lons"} where lats/ lons are the footprint exterior-ring coordinate arrays (WGS84) and s3/https are the canonical data-asset hrefs (either may be None). Records with no canonical data asset (raster sources,

  218) additionally carry assets ({key: href} for every

  non-canonical asset), datetime (ISO acquisition time), and time_key (the acquisition-group property named by the catalog's time_key metadata, when present); any record with a data/data_s3 asset -- every CMR record, including preserve_thumbnails=True -- keeps its exact pre-#218 shape, except that any record whose catalog carries STAC start_datetime/end_datetime also gains time_start/ time_end (ISO acquisition range, issue #246) for the per-window dispatch subsetting.

granule_row_mask

granule_row_mask() -> ndarray

Boolean mask over table rows of the rows :meth:granule_records emits.

:meth:granule_records skips rows whose geometry is empty or not polygonal, so record i is table row np.flatnonzero(mask)[i] and mask.sum() is the record count -- both without decoding a single record. Same predicate as the row-wise loop, applied with one batched shapely.from_wkb instead of one call per row, and refusing on the one input where "same predicate" would otherwise be false (see Raises).

Two callers need the alignment without the records: :meth:index_footprints, which gives screened rows an empty MOC so the column stays one entry per table row, and ShardMap.build's stored-index fast path, which intersects on the column before materializing anything (issues #396, #439).

Returns:

• ndarray –

  bool, length table.num_rows.

Raises:

• ValueError –

  When any row's geometry is null. shapely.from_wkb maps a null to None, whose get_type_id is -1 -- so screening it out here would be silent, while :meth:granule_records' row-wise geom.is_empty raises AttributeError on the same row. A granule that vanishes from an indexed build but crashes an unindexed one is worse than either, so both paths refuse: this one loudly, naming the count and the first offending row.

Notes

Time is small (0.02 s over the 35,639-granule 88S catalog, ~0.3 s over the 555,867-row ATL03 clone) but memory is not: the to_numpy WKB copy plus the live shapely objects put the clone's RSS ~1 GB over the parquet read. It is a peak, not a leak -- the objects die with the call.

index_footprints

index_footprints(order: int) -> 'Catalog'

Precompute the footprint_cells morton MOC column (issue #396).

Covers every row's geometry WKB once, at order, and returns a new catalog carrying the result as a ragged column plus the order in its metadata. The per-granule footprint cover is identical for every ShardMap.build against this catalog, so paying it once here turns each later build into set intersection with no geometry work at all -- see ShardMap.build's fast path.

Parameters:

• order (int) –

  HEALPix order to cover at. Choose the shard order (the grid's parent_order), not the finer chunk order: coverage words per granule roughly double per order, so a full ATL03 clone indexes to ~270-420 MB of parquet at order 9 but ~9-13 GB at order 13, and the extra resolution is invisible to order-9 shard cells. The column serves every grid whose parent_order is at most order; a finer grid is refused by build rather than answered coarsely. Bounded above by MORTIE_MOC_ORDER_CAP (mortie's order-18 coverage cap), the same bound ShardMap.build clamps its own MOC order to, so the two paths agree on what "too fine" means. build clamps because its order is derived; this one is the operator's own number, so it raises rather than silently indexing at another order than the one recorded in footprint_cells_order.

Raises:

• ValueError –

  When order is above MORTIE_MOC_ORDER_CAP.

Returns:

• Catalog –

  New catalog, same rows and metadata plus footprint_cells and footprint_cells_order. Re-indexing at another order replaces the column rather than appending a second one.

Notes

mortie.arrow.from_wkb (mortie >= 0.9.5, espg/mortie#157/#163) takes the geometry column across the Python/Rust boundary once, with the GIL released and chunking that bounds peak at roughly the result size. It covers the union of the rings inside each blob, where :meth:granule_records reads the largest part's exterior ring only, so the two agree exactly on single-part footprints (every CMR ATL03/06 granule) and the column is a superset for a MultiPolygon -- it keeps the smaller parts granule_records drops.

Rows :meth:granule_records would skip -- empty or non-polygonal geometry -- are screened out with the same shapely predicate it uses and get an empty MOC, so the column stays one entry per table row and a catalog carrying a stray Point indexes rather than raising (mortie's coverage refuses a point outright, naming the blob). The screen is :meth:granule_row_mask -- one vectorised shapely.from_wkb, shared with ShardMap.build's fast path; on the 35,639-granule 88S catalog it is 0.02 s against 2.65 s for the order-9 cover, so it costs under 1% of a pass that runs once per catalog. A null geometry is the one row the screen refuses instead of skipping, so this method raises where it would once have indexed -- see :meth:granule_row_mask's Raises.

The screen is cheap in time but it is this pass's peak in memory, and it is the term from_wkb's chunking does not bound: on the 555,867-row ATL03 clone RSS goes 835 MB after the parquet read -> 1,169 MB after to_numpy (a full WKB copy) -> 1,794 MB with the shapely objects live. They die with granule_row_mask's frame, so that stays a peak rather than stacking with the cover, but a whole-clone index wants headroom for it. :meth:cover_footprints' keep.all() short-circuit avoids a second ~334 MB WKB copy in the case that actually occurs (nothing screened out -- every catalog in the tree). Reading the geometry-type word straight out of the WKB, or chunking the from_wkb call, would drop the screen's peak entirely; not done here because it trades the shared shapely predicate for a hand-rolled one.

to_geoparquet

to_geoparquet(path: str) -> None

Write the catalog to a stac-geoparquet file.

stac_geoparquet rewrites schema metadata with only the GeoParquet geo key, so we reopen and merge zagg provenance back in (keeping geo intact) before the final write.

Shard map

zagg.catalog.shardmap.ShardMap dataclass

Work-distribution manifest: shard key -> granules, tied to one grid.

Parameters:

• grid_signature (dict) –

  grid.spatial_signature() at build time -- the spatial layout only (#89). The runner checks it against the run grid's spatial signature so a map can't be paired with a mismatched spatial grid, while staying reusable across configs that differ only in aggregation fields. (Kept as grid_signature for back-compat; old maps carry the full signature and still validate via a spatial-subset projection.)

• shard_keys (list of int) –

  Sorted shard keys with at least one granule.

• granules (list of list of dict) –

  Parallel to shard_keys. Each granule is {"id", "s3", "https"} (option C -- self-contained, endpoint-neutral).

• metadata (dict, default: dict() ) –

  Provenance copied from the Catalog plus backend/timing info.

aoi_mask class-attribute instance-attribute

aoi_mask: List[List[int]] | None = None

Optional strict-AOI per-shard mask payload (issue #101), parallel to shard_keys. None when output.aoi_mask is off (the default) — the manifest then carries no extra key and is byte-identical to a pre-feature map. Each entry is a JSON int list the grid expands to a per-cell bool over the shard's children(): a compact MOC (HEALPix) or the True-cell indices into children order (rectilinear).

build classmethod

build(
    catalog,
    grid,
    *,
    region=None,
    aoi=None,
    backend: str = "auto",
    mortie_order: int | None = None,
    footprint: str = "swath",
    sibling_catalog=None,
    sibling_asset: str = "l2a",
) -> "ShardMap"

Build a ShardMap from a Catalog and an output grid.

Parameters:

• catalog (Catalog) –

  Fetched granule metadata (provides granule_records()). A catalog indexed by Catalog.index_footprints additionally carries every granule's morton MOC, which takes the geometry-free fast path described in the Notes (issue #396).

• grid (OutputGrid) –

  Output grid (provides coverage, shard_footprint, spatial_signature).

• region (list of (lats, lons), default: None ) –

  Coverage mask in WGS84. Defaults to the catalog bbox rectangle.

• aoi (AOIGeometry | bytes | str | list of (lats, lons), default: None ) –

  Strict-AOI polygon for the optional output.aoi_mask (issue #101), supplied as an :class:~zagg.grids.aoi.AOIGeometry, WKB bytes, WKT str, or (lats, lons) ring parts. None (default) reuses region (or the bbox rectangle), so a ring run is unchanged. Only consulted when output.aoi_mask is on — a flag-off run never builds it and stays byte-identical.

• backend (('auto', 'spherely', 'mortie'), default: "auto" ) –

  Geometry backend. "auto" -> spherely when importable, else mortie for HEALPix grids (non-HEALPix grids require spherely and raise an ImportError with an install pointer when it is absent).

• mortie_order (int, default: None ) –

  MOC order for the mortie backend. None (default) pins it to the grid's inner-chunk order grid.chunk_order (the chunk_inner order, defaulting to parent_order when unset), clamped to mortie's order-18 coverage cap -- the dispatch chunk's own resolution, enough to keep moc_to_order from upsampling a footprint onto neighbor shards (#92) at near-minimal compute. Raises if the resolved order is coarser than parent_order.

• footprint (('swath', 'beams'), default: "swath" ) –

  Granule footprint used for intersection. "swath" (default) uses the raw CMR polygon. "beams" decomposes ICESat-2 ATL03/06 swaths into per-beam-pair corridors so granules stop being assigned to shards their beams never cross (issue #65); non-beam products fall back to the swath ring.

  .. deprecated:: The "beams" corridor mechanism is a stopgap (see beams.py); remove it once native per-beam CMR geometry, the memory-handling robustness in #66, or data virtualization (#97) lands.

• sibling_catalog (Catalog, default: None ) –

  Paired-asset sibling product (issue #425), e.g. the GEDI02_A catalog beside a GEDI01_B primary. Granules are joined at build time on :func:sibling_join_key (the shared orbital id core, pinned within a product generation); each paired entry carries the sibling's hrefs under assets[sibling_asset], spatial intersection runs once on the primary footprints (the products share them). Pairless granules are excluded and reported: metadata["pairless"] lists {id, missing} for every primary without a sibling, sibling without a primary, and sibling shadowed by an earlier record on the same join key ("duplicate-key"), and a build-time warning fires when the list is non-empty. The sibling catalog is taken as authoritative, so it must be queried over the same AOI and time window as the primary — a narrower one reports genuinely-paired acquisitions as missing and thins the product; past _PAIRLESS_ALERT_FRACTION of the primary catalog the warning escalates to name that cause. A paired build always takes the eager/geometry path: pairing filters the record list, which the footprint_cells fast path cannot follow (its alignment is positional over the raw table -- issue #439), so the stored-index plan is skipped even on an indexed catalog.

• sibling_asset (str, default: 'l2a' ) –

  Asset key the sibling's hrefs are stored under (default "l2a"); matches the data_source.assets name the reader joins on.

Returns:

• ShardMap –

Notes

Indexed catalogs skip the geometry entirely (issue #396). When the catalog carries the footprint_cells column, the backend resolves to mortie, footprint="swath" and mortie_order is left at its default, the intersection becomes moc_and of each stored granule MOC with the AOI's own shard MOC -- no WKB parse, no coverage walk. The result is the mortie backend's, unchanged on the single-part footprints every CMR ATL03/06 granule has; metadata["footprint_cells"] records that the index answered the build. The one intended divergence is a MultiPolygon footprint, where the column is a superset: it covers every part, while :meth:~zagg.catalog.Catalog.granule_records reads only the largest part's exterior ring, so the index can place such a granule in shards the geometry path misses. A column coarser than the grid's parent_order raises rather than answering (see :func:_footprint_cells_plan).

That path also materializes no granule records until after the intersection (issue #439). Decoding every row's WKB and asset map ran ~25 s of a 29.5 s build over the 555,867-granule ATL03 clone against a California AOI, to feed an intersection that assigned 2,357 of them; the stored column is already table-row-ordered, so the build intersects on it and then decodes only the rows it kept. The manifest is unchanged -- total_granules still counts the records considered (the whole post-screen catalog), which :meth:~zagg.catalog.Catalog.granule_row_mask supplies without decoding any.

An unindexed mortie build covers first too (issue #445). "Indexed" now means only that the cover was persisted: a mortie swath build on a HEALPix grid with no usable column covers the WKB column itself (:meth:~zagg.catalog.Catalog.cover_footprints) and runs the same intersection, instead of decoding every record to cover from its rings. On the clone/California case that is 86.9 s -> ~33 s. Records still follow the intersection, and total_granules is still the screen popcount. The cover runs at the grid's parent_order -- the order the map's shard membership is stated in, and the order beyond which the real-catalog sweep in :func:_resolve_mortie_order measures granules per shard flat -- so an unpinned unindexed build records mortie_order = parent_order where it used to record the chunk order. Same assignment at the orders production runs (byte-identical over the 555,867-granule clone, and 0/200 rows differing in a randomized sweep at parent/chunk 9/13, 11/13, 8/12 and 9/11), and superset-safe in general: mortie's coverage is conservative per order, so the coarser cover can admit a boundary shard that the finer cover refined down misses -- up to 0.71% extra cells at 3/7 -- and never drops one. See :func:_live_cells_plan. A caller-pinned mortie_order is still honored literally. Unchanged paths: footprint="beams" (the cover is the swath, not the corridors), the spherely backend, rectilinear grids, and paired builds (issue #425).

Two deltas ride on it, both inherited from the cover this path now shares with the index. MultiPolygon footprints get the same superset as above -- the cover takes every part where :meth:~zagg.catalog.Catalog.granule_records reads only the largest part's exterior ring -- so a multi-part granule may assign to more shards than it did pre-#445; every CMR ATL03/06 granule is single-part, so no production build moves. And a null geometry is refused by name (:meth:~zagg.catalog.Catalog.granule_row_mask) where the record loop raised AttributeError on the same row: both refuse, one legibly. Memory is the index's documented posture, screen peak included (~1 GB over the parquet read at clone scale, issue #429).

A pinned mortie_order covers live whether or not the catalog is indexed. A persisted column cannot restate an arbitrary order, so the stored plan declines a pin (:func:_footprint_cells_plan) -- and since issue #445 what catches it is the ephemeral cover, not the geometry path. So an indexed catalog built with an explicit mortie_order is a third population whose assignment can move: it inherits the same two deltas above (measured on a two-part MultiPolygon: 8 shards on the records path, 13 on the cover). The pin is honored either way, and metadata["footprint_cells"] still records False there -- it means "the stored column did not answer this build", which stays the thing the operator can act on. Flagged for review on PR #447.

from_json classmethod

from_json(path: str) -> 'ShardMap'

Load a manifest from JSON.

from_parquet classmethod

from_parquet(path: str) -> 'ShardMap'

Load a manifest from the parquet form written by :meth:to_parquet.

Importing :mod:mortie.arrow first registers the morton_index extension type, so the shard_keys column rehydrates typed; the words are pulled over the C Data Interface (import_c_array) regardless.

reproject

reproject(target_grid, catalog=None) -> 'ShardMap'

Derive a ShardMap at target_grid's parent_order (issue #294).

HEALPix nesting means a shard map at one order is derivable from another without touching the source catalog again in the coarsen direction, and with only a scoped (per-shard) re-intersection in the refine direction -- either is far cheaper than a full :meth:build over the whole catalog.

Parameters:

• target_grid (HealpixGrid) –

  Grid to reproject onto. Must share child_order (the leaf resolution) with the source grid -- reprojecting across different DGGS resolutions isn't meaningful, only the shard (dispatch) order changes.

• catalog (Catalog, default: None ) –

  Required only when refining (target_grid.parent_order > this map's parent_order): the shard map itself stores only {"id", "s3", "https"} per granule, not footprint geometry, so recovering which finer cell each granule falls in needs the granule records (catalog.granule_records()) back.

Returns:

• ShardMap –

  New map at target_grid.parent_order. metadata records the source order and reproject: {method: "coarsen"|"refine"|"noop"}.

Notes

Coarsen (target_order < source_order): pure regroup, exact, no geometry. Each source shard's key coarsens via mortie.clip2order(target_order, shard_key); shards sharing a coarse parent are grouped and their granule lists unioned, deduplicated by granule id (a granule spanning multiple children counts once in the parent). Exact because footprint-to-cell assignment nests: a granule intersects a coarse cell iff it intersects one of its finer children.

Refine (target_order > source_order): cannot be a pure regroup -- the coarse map never recorded which child cell a granule fell in. Instead, for each source shard, its own granules are re-intersected at target_order (the same polygons_to_morton_mocs machinery :meth:build uses), restricted to that shard's own descendant cells (generate_morton_children(shard_key, target_order)). In the interior this reproduces the direct :meth:build at the finer order; at a region/AOI boundary it may over-include child shards a region-restricted direct build would drop (the #101 whole-shard overhang class -- reproject applies no region clip), but it never drops a real intersection. Costs only this shard's granules, not the whole catalog. Refine always re-intersects via the mortie MOC path regardless of the source's backend (a spherely-built source is not reproduced by it), so the derived map records backend="mortie".

to_json

to_json(path: str) -> None

Write the manifest as JSON.

to_parquet

to_parquet(path: str) -> None

Write the manifest as parquet with a TYPED morton shard_keys column.

The Arrow-native sibling of :meth:to_json (issue #135): shard_keys carries mortie's morton_index pyarrow extension type (registered by mortie on import), so any Arrow-aware consumer sees morton words, not anonymous ints. granules (and aoi_mask when present) ride as per-shard JSON strings — the same self-contained payloads the JSON form stores — and metadata/grid_signature live in the schema metadata, mirroring the Catalog geoparquet convention (sources.py).

Requires pyarrow (the off-Lambda catalog extra); the worker path never calls this — the runner dispatches from the JSON manifest.

Convenience

zagg.catalog.make_shardmap

make_shardmap(
    query,
    grid,
    *,
    region=None,
    aoi=None,
    backend="auto",
    catalog_out=None,
    footprint="swath",
)

Fetch a Catalog and build a ShardMap in one call (concerns 1+2 chained).

Parameters:

• query (Query) –

  What/when/where to fetch.

• grid (OutputGrid) –

  Output grid (typically from_config(config)).

• region (list of (lats, lons), default: None ) –

  Coverage mask. Defaults to the query bbox rectangle.

• aoi (AOIGeometry | bytes | str | list of (lats, lons), default: None ) –

  Strict-AOI polygon for the optional output.aoi_mask (issue #101) — WKB bytes, WKT str, an :class:~zagg.grids.aoi.AOIGeometry, or ring parts. None reuses region. Forwarded to :meth:ShardMap.build.

• backend (str, default: 'auto' ) –

  Geometry backend for the shard map.

• catalog_out (str, default: None ) –

  If given, persist the fetched Catalog to this geoparquet path.

• footprint (('swath', 'beams'), default: "swath" ) –

  Granule footprint for intersection; "beams" tightens ICESat-2 ATL03/06 assignment to per-beam-pair corridors (issue #65).

  .. deprecated:: The "beams" corridor mechanism is a stopgap. Remove it once a better fix lands -- native per-beam CMR geometry, the memory-handling robustness in #66, or data virtualization tracked in #97.

Returns:

• ShardMap –

Temporal / spatial helpers

zagg.catalog.cycle_to_dates

cycle_to_dates(cycle: int) -> tuple[datetime, datetime]

Convert an ICESat-2 repeat cycle number to a (start, end) date range.

Parameters:

• cycle (int) –

  ICESat-2 cycle number (1-based).

Returns:

• tuple of (datetime, datetime) –

zagg.catalog.load_polygon

load_polygon(geojson_path: str) -> list[tuple]

Load polygon(s) from a GeoJSON file as (lats, lons) parts.

Supports Feature, FeatureCollection, Polygon, and MultiPolygon geometries.

Parameters:

• geojson_path (str) –

  Path to a GeoJSON file.

Returns:

• list of (lats, lons) –

  One coordinate-array pair per polygon ring (WGS84).

zagg.catalog.polygon_to_bbox

polygon_to_bbox(parts: list[tuple]) -> tuple[float, float, float, float]

Compute a (lon_min, lat_min, lon_max, lat_max) bbox from polygon parts.

Parameters:

• parts (list of (lats, lons)) –

Returns:

• tuple of (lon_min, lat_min, lon_max, lat_max) –

zagg.catalog.load_antarctic_basins

load_antarctic_basins(filepath=None) -> list[tuple]

Load Antarctic drainage basin polygons as (lats, lons) parts.

Parameters:

• filepath (str, default: None ) –

  Path to the basin polygon file. Defaults to the file shipped with mortie.

Returns:

• list of (lats, lons) –

  One pair per basin.

Closest-observation pairing (issue #509)

One raster store (e.g. Sentinel-2 L2A) can serve several point-cloud reference stores (ATL03 + GEDI): for every reference epoch a shard's stores actually observed, ingest the single nearest acquisition from the raster catalog. The pairing is a property of the ingest query, not the store schema — the raster store stays a plain raster store, which granules were ingested is the pairing, and coincidence at read time is toc intersection.

Epochs are store-derived, never catalog-derived: each reference store's coverage.toc sibling (spec §10.5) records per-shard word-set covers of the data that actually landed, quantized at temporal order 18 (2^45 ns ≈ 9.77 h buckets). The builder expands each cover word into its constituent buckets and takes one epoch per bucket midpoint — good to ±4.9 h against Sentinel-2's ~4.3-day revisit. Granule catalogs would inherit the CMR-hull over-assignment (~70 assigned granules vs 49 contributing pass-days on a measured Californian shard); covers reflect contribution, not assignment.

from zagg.catalog.closest_obs import closest_obs_shardmap
from zagg.catalog.sources import Catalog
from zagg.grids import HealpixGrid
import numpy as np

grid = HealpixGrid(9, 13)  # parent_order must equal the covers' shard order
s2 = Catalog.from_geoparquet("catalog_s2_ca.parquet")

# Size the run first — the dry-run builds nothing and prices the fan-out:
est = closest_obs_shardmap(
    s2,
    ["s3://bucket/atl03_store", "s3://bucket/gedi_store"],
    grid=grid,
    aoi="california.geojson",
    max_time_offset=np.timedelta64(3, "D"),
    max_granules_per_shard=200,  # the same gate the build below applies
    estimate=True,
)
# violations is [] unless the gate is passed here too -- estimate returns
# before the build's raise, so this is the safe way to size it.
est["histogram"], est["max_cost_usd"], est["violations"]

# Then build the map; dispatch consumes it like any other ShardMap:
sm = closest_obs_shardmap(
    s2,
    ["s3://bucket/atl03_store", "s3://bucket/gedi_store"],
    grid=grid,
    aoi="california.geojson",
    max_time_offset=np.timedelta64(3, "D"),
    max_granules_per_shard=200,
)
sm.to_json("s2_closest_obs.json")

Everything refuses or records loudly, never silently: a reference store with no readable coverage.toc raises (sweep the store first); an epoch whose nearest acquisition lies beyond max_time_offset selects nothing and is recorded per-epoch in metadata["closest_obs"]["dropped"] with its near-miss offset; a shard past max_granules_per_shard raises naming the worst shards (estimate=True reports the violations instead, so the gate can be sized first); a cover block coarsened below the §10.5 pin is warned about and reported in coarsened_orders — and under a max_time_offset, epochs whose coarse-bucket half-span exceeds the stated offset cannot be paired to that precision, so they drop into the ledger as their own category (epochs_dropped_low_resolution, rows naming the block's effective order; espg tolerance ruling 2026-08-24). Selected granule entries carry paired_epochs / epoch_offsets_ns provenance so the paired product is reconstructable from the manifest alone. Epochs are bucket midpoints — size max_time_offset with ReferenceEpochs.tolerance()'s half-bucket slack in mind.

zagg.catalog.closest_obs.reference_epochs

reference_epochs(
    reference_stores, *, aoi=None, **store_kwargs
) -> ReferenceEpochs

Per-shard epochs from the reference stores' coverage.toc covers.

For each store root: fetch the §10.5 sibling (:func:zagg.coverage_toc.read_cover), strict-decode its per-shard word sets (:func:zagg.coverage_toc.cover_words), and expand each word into its constituent buckets' midpoints (:func:_word_midpoints). Per shard the result is the union across stores, deduplicated (espg ruling: one raster store serves both sensors, epochs are the union across the reference stores).

The union is canonical at the bucket level, not the word level. Words are post-toc_normalize runs whose extent depends on that store's other data, so two stores sharing one pass routinely emit overlapping-but-unequal range words — a raw np.unique over words would keep both and represent the shared pass twice, at two displaced midpoints. Expanding to buckets first removes that degree of freedom: the bucket grid is fixed by the order alone, so a shared pass contributes the same bucket midpoint from every store and dedupes exactly.

§10.5 lets a block coarsen below the object's pinned temporal order to fit the cover cap, and its landed order is recorded in the block itself. :func:zagg.coverage_toc.cover_words returns only the words, so this reads each block's temporal_order straight from the grammar, expands that block's words at its order, logs a warning whenever a block sits below the pin (the read half of §10.5's "widening only, loudly recorded"), and reports the coarsest contributing order per shard on :attr:ReferenceEpochs.orders.

A store that carries no readable cover refuses loudly — this builder is cover-driven by design (store-derived epochs, never the granule catalogs), so "no cover yet" means "sweep the store first", not "fall back silently". Likewise a shard-order mismatch between stores: D1 ids at two orders are not comparable.

Parameters:

• reference_stores (str or sequence of str) –

  Store roots whose covers drive the epochs (e.g. ATL03 + GEDI).

• aoi (optional, default: None ) –

  Restrict the shard set: a mortie.Moc, a GeoJSON path, or [(lats, lons), ...] ring parts (see :func:_aoi_shard_set).

• **store_kwargs –

  Forwarded to the object-store open (region, credentials, ...).

Returns:

• ReferenceEpochs –

zagg.catalog.closest_obs.ReferenceEpochs dataclass

Per-shard reference epochs decoded from one or more store covers.

Attributes:

• order (int) –

  The covers' shard (dispatch) order — every contributing store must agree on it, and the raster grid the epochs pair against must match.

• epochs (dict of int -> np.ndarray) –

  Shard key (packed morton word, the canonical in-memory form D1 renders as a decimal string externally) -> sorted unique UTC datetime64[ns] epoch midpoints, the union across the contributing stores' covers. Shards whose cover block decodes to an empty word set are omitted (they claim nothing).

• stores (list of str) –

  The store roots that contributed, in the order given (provenance).

• orders (dict of int -> int) –

  Shard key -> the coarsest temporal order that contributed to that shard's epochs. Normally the pinned :data:~zagg.coverage_toc.TEMPORAL_COVER_ORDER, but §10.5 lets a block coarsen below the pin to fit the cover cap, and a coarser order means a wider midpoint bound (2**(62 - order) ns). Phase 2's max_time_offset reasoning reads this rather than assuming the pin.

epoch_orders class-attribute instance-attribute

epoch_orders: dict[int, ndarray] = field(default_factory=dict)

Shard key -> int64 array row-aligned with epochs: each epoch's own effective temporal order (espg tolerance ruling, 2026-08-24). :attr:orders is the per-shard COARSEST — right for a headline warning, too blunt for a precision gate: one shard can mix a pinned store's epochs with a coarsened store's, and only the coarse ones fail a stated max_time_offset. A midpoint claimed at two orders keeps the finest.

total property

total: int

Total epoch count across shards (post-union, post-dedupe).

tolerance

tolerance(shard: int) -> timedelta64

Half a bucket at shard's effective order — its epoch bound.

zagg.catalog.closest_obs.nearest_acquisitions

nearest_acquisitions(epochs, times, *, max_time_offset=None)

Nearest acquisition per epoch — the vectorized closest-1 core.

The phase-2 selection of the closest-observation join (issue #509): for each reference epoch, the single nearest acquisition — closest-1, never a two-sided bracket (espg ruling); several epochs bracketing one acquisition each select it, and the builder dedupes granules downstream.

Parameters:

• epochs (array-like of datetime64[ns]) –

  One shard's reference epochs (:func:reference_epochs). Any order. NaT refuses — see the note below.

• times (array-like of datetime64[ns]) –

  The shard's acquisition times, in catalog record order. Any order — the returned selection indexes THIS array's positions. NaT refuses — see the note below.

• max_time_offset (timedelta64 or int, default: None ) –

  An epoch whose nearest acquisition lies further than this selects nothing. Exactly-at selects; one ns past does not. A bare int means nanoseconds. None (default) always selects the nearest, however far. Negative refuses, and so does a duration that will not convert to nanoseconds (np.timedelta64(1000, "Y") overflows int64 ns). Callers gating against the epochs' own precision should widen by :meth:ReferenceEpochs.tolerance — a cover epoch is a bucket midpoint, good to half a bucket, not exact.

Returns:

• selection ( np.ndarray of int64 ) –

  Per epoch, the index into times of the selected acquisition, or -1 where the epoch selects nothing (no acquisitions at all, or nearest beyond max_time_offset).

• offsets ( np.ndarray of timedelta64[ns] ) –

  Per epoch, the SIGNED offset times[nearest] - epoch of the nearest acquisition — positive when the acquisition follows the epoch. Reported for every epoch, dropped ones included (the loud record a drop rides — the builder's report needs the near-miss distance, not just the fact of the drop). NaT when times is empty, and for the offset that will not fit timedelta64[ns] — selection is exact over the whole datetime64[ns] span (~584 years), but a difference past ~292 years is unrepresentable, so it saturates to NaT rather than reporting a wrapped duration. Such an epoch selects nothing under any cap (no finite tolerance reaches it) and still selects its nearest with max_time_offset=None.

Raises:

• ValueError –

  If either input carries NaT. A missing instant has no nearest anything: NaT sorts last but casts to iinfo(int64).min, so it would leave ts unsorted (silently mis-pairing its neighbours) and wrap the offset subtraction. :func:reference_epochs never emits NaT and the catalog's time parser refuses a missing acquisition time, so NaT here is caller debris — refused loudly, this module's posture.

Notes

A tie — an epoch exactly equidistant between two acquisitions — selects the EARLIER acquisition, deterministically. Equal acquisition times (one Sentinel-2 datatake stamps many granules with the same instant) are broken by catalog record order: the FIRST record of the equal-time run, whichever flank the epoch approaches it from.

zagg.catalog.closest_obs.closest_obs_shardmap

closest_obs_shardmap(
    s2_catalog,
    reference_stores,
    *,
    grid,
    aoi=None,
    max_time_offset=None,
    max_granules_per_shard=None,
    estimate=False,
    backend="auto",
    bytes_per_granule=None,
    **store_kwargs,
)

Closest-observation ingest map: covers -> epochs -> nearest granules.

The issue-#509 builder. Epochs come from the reference stores' coverage.toc covers (:func:reference_epochs); the raster catalog is spatially assigned to shards through the existing machinery (:meth:~zagg.catalog.shardmap.ShardMap.build — the stored-index / batch-cover fast paths included); per shard, each epoch then selects its single nearest acquisition (:func:nearest_acquisitions) and the selected granules, deduplicated, become the shard's ingest list. The result is a standard :class:~zagg.catalog.shardmap.ShardMap — JSON round-trip, total_pairs bookkeeping — so dispatch consumes it unchanged; the pairing is a property of this query, never of the raster store's schema (espg ruling).

Parameters:

• s2_catalog (Catalog or str) –

  The raster acquisition catalog (stac-geoparquet), or a path to one (Catalog.from_geoparquet). Every record must carry an acquisition time (datetime, or STAC start_datetime).

• reference_stores (str or sequence of str) –

  Store roots whose covers drive the epochs (e.g. ATL03 + GEDI).

• grid (HealpixGrid) –

  The raster store's output grid. Its parent_order must equal the covers' shard order — the map's shard keys and the covers' D1 ids live on the same grid, and the emitted grid_signature is what the ingest run validates against.

• aoi (optional, default: None ) –

  Restrict the shard set (mortie.Moc, GeoJSON path, or ring parts); intersected with the store coverage (:func:reference_epochs). Ring parts / a GeoJSON path additionally ride as the spatial build's region=, scoping the raster intersection; a mortie.Moc has no parts form, so it restricts the epoch side only and the spatial build still covers the whole catalog bbox.

• max_time_offset (timedelta64 or int(ns), default: None ) –

  An epoch whose nearest acquisition lies beyond this selects nothing — recorded per epoch in metadata["closest_obs"]["dropped"] and warned about, never silent. It is also a precision bar on the epochs themselves (espg tolerance ruling, 2026-08-24): an epoch from a cover block coarsened far enough that its bucket half-span exceeds this offset cannot be paired to the stated precision, and is dropped into the ledger as its own category (rows carrying temporal_order/cover_half_span_ns; counted in epochs_dropped_low_resolution). Half-span exactly at the offset stays pairable — the same side the selection gate's exactly-at rule pins. None always selects the nearest, at whatever resolution the covers offer (a single warning names the effective resolution when any block sits below the §10.5 pin).

• max_granules_per_shard (int, default: None ) –

  Cost gate: a shard exceeding this REFUSES loudly (ValueError naming the worst shards) — never truncates. estimate=True reports violations instead of raising, so the gate can be sized first.

• estimate (bool, default: False ) –

  Dry-run: return the per-shard histogram + cost estimate dict (see Notes) WITHOUT building the map. The espg-operated ingest runs are cost-gated through this.

• backend (('auto', 'spherely', 'mortie'), default: "auto" ) –

  Forwarded to :meth:ShardMap.build.

• bytes_per_granule (int, default: None ) –

  Per-granule ingest volume for the estimate byte figure; without it est_bytes is None (the catalog does not carry sizes).

• **store_kwargs –

  Forwarded to the reference stores' object-store opens.

Returns:

• ShardMap or dict –

  The ingest map — or, with estimate=True, a dict: {"shards", "granules", "pairs", "epochs_total", "epochs_paired", "epochs_dropped", "epochs_dropped_low_resolution", "per_shard" (decimal -> granule count), "histogram" (granule count -> shard count), "est_bytes", "max_cost_usd", "violations"}. max_cost_usd is the :func:zagg.dispatch.max_cost_usd ceiling at the production worker size (one invoke per shard, 900 s timeout) — a bound, not a forecast.

Notes

A selected entry whose record carries no datetime (STAC's null-datetime + start_datetime form) is emitted with datetime backfilled from time_start — the same instant the pairing used, and the key raster dispatch requires. Selected granule entries also gain two provenance keys so the eventual paired product is reconstructable from the manifest alone: paired_epochs (ISO instants of every epoch that selected the granule) and epoch_offsets_ns (row-aligned SIGNED acquisition - epoch ns). :meth:ShardMap.reproject carries both through every arm (issue #517), so a derived map keeps its per-entry provenance; because these keys are accumulated per shard, a coarsen unions the rows of every child a granule reached rather than keeping one child's. What reproject still cannot fix is metadata["closest_obs"], which rides through describing the SOURCE map — after a coarsen its shard ids name a shard set that no longer exists. Rebuild at the target grid rather than reprojecting whenever that ledger has to be true of the map carrying it. metadata["closest_obs"] records the query: the reference stores, the epoch totals, every dropped epoch with its near-miss offset (None where there was no acquisition to measure against, so epochs_total == epochs_paired + epochs_dropped holds), shards whose epochs found no acquisition at all, and any cover blocks coarsened below the §10.5 pin. Epochs are bucket midpoints, good to half a bucket (:meth:ReferenceEpochs.tolerance) — size max_time_offset with that slack in mind.

The strict-AOI per-shard mask (output.aoi_mask, issue #101) is the SPATIAL map's payload and is not carried onto the emitted map, so the aoi_mask metadata claim is dropped with it: an output.aoi_mask ingest run must compute its mask at run time rather than read one off this manifest.