Support arbitrary xgcm topologies for boundaries and masks#22
Draft
hdrake wants to merge 15 commits into
Draft
Conversation
- Add regionate/geometry.py with pure-shapely antimeridian helpers: normalize_lon, has_antimeridian_crossing, polygon_to_lon360, split_at_antimeridian (ported from xcryocouple). - Rewrite get_region_boundary_grid_indices to return a 7-tuple including f_c, normalizing sectionate.grid_section's 4/5-tuple return and passing f_c through uvcoords_from_qindices. Drop the topology= kwarg. - Rewrite mask_from_grid_boundaries to split the region polygon at the antimeridian and OR per-piece regionmask rasters; discover tracer-center lon/lat coords via sectionate (no hard-coded geolon/geolat). - Add regionate/tests/test_geometry.py (unit + functional). Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
Preserve the single-tile contourpy path byte-for-byte (shared _contour_to_corner_indices) and add a multi-tile path: per-face boundary faces from a topology-aware padded mask (real-neighbour halo), with internal seam faces removed and surviving faces stitched into closed loops by physical coincidence. Returns a new f_c_list (per-corner face index; None for single-tile). Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
- GriddedRegion/MaskRegions accept ij=(i_c, j_c, f_c) and store self.f_c; initiate_from_boundary unpacks the new 7-tuple from get_region_boundary_grid_indices. - BoundedRegion threads f_c through grid_section/uvcoords/GriddedSection. - Replace hard-coded geolon_c/geolat_c presence checks with a get_geo_corners(grid) probe (arbitrary corner-coord names). - Fix the genuinely broken .gr save/load: write/read consistent lons_c/lats_c/i_c/j_c names, persist+reload f_c, set lons_uv/lats_uv and children only when present, and guard the children directory on load. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
Add explicit shapely and xgcm>=0.9.0 (both used directly). Keep the sectionate floor at >=0.3.3 because the topology-driven branch (MOM6-community/sectionate#47) self-reports 0.3.3; CI installs that branch from its PR ref until it is released to PyPI. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
…geometry - Self-contained two-face (face_connections) fixture with seam-spanning, seam-terminating, interior-cell, and MaskRegions f_c-threading tests. - Optional, skipped-by-default MOM6-global and ECCO/LLC checks (REGIONATE_REALDATA_TESTS); ECCO asserts reciprocal-or-raises. - Export regionate.geometry helpers from the package namespace. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
A boundary whose net longitude winding is ~360 deg encircles a pole and has no simple lon/lat polygon. Re-introduce the South-Pole extension (sectionate's stereographic-plane orientation convention) for that case only, then rasterize its output through the same antimeridian-split + per-piece OR path as ordinary regions. Normal bounded regions are unaffected (net winding ~0). Add a regression test asserting a zonal equator line encloses the southern hemisphere for both windings. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
…e#47 - Re-execute all four example notebooks against the topology-driven stack (new sectionate API + regionate changes) on the real MOM6 example data; all run cleanly (validates the high-level API including BoundedRegion). - CLAUDE.md: document arbitrary-topology support, the f_c face index, the rewritten mask/boundary algorithms, new geometry.py, and corrected .gr I/O. - installation.rst: note the sectionate#47 requirement and dev-install step. - Verified the Sphinx/nbsphinx docs build succeeds. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
|
Check out this pull request on See visual diffs & provide feedback on Jupyter Notebooks. Powered by ReviewNB |
examples/load_example_model_grid_gfdl.py has been dead code since the Zenodo-data migration (commit 798528c): nothing in the repo imports it or its functions (load_CM4p25/load_OM4p5/fix_grid_coords), and it only runs against GFDL-internal /archive paths. Removing it; git history preserves it. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
New example demonstrating regionate on the ECCO lat-lon-cap (LLC90) grid -- a real 13-tile face_connections topology: - examples/load_example_ECCO_grid.py: downloads the ECCO geometry from PO.DAAC (earthaccess), renames XC/YC/XG/YG to geolon*/geolat*, builds the LLC90 grid with canonical face_connections, and converts the native MITgcm 'left' staggering to symmetric via sectionate.gridutils.symmetrize. - examples/5_ECCO_LLC90_multiface_regions.ipynb: traces a region straddling the tile-1/tile-2 seam -> a single boundary loop whose face index spans both tiles; plots the mask and stitched boundary. Added to the docs toctree. - test_real_grids.py: replace the ECCO stub with a concrete seam-stitching test (runs when the geometry file is present; skipped otherwise). - docs/environment.yml: add earthaccess + pandoc (nbsphinx render dep). Requires sectionate.symmetrize (MOM6-community/sectionate#47 branch, local commit). Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
The multi-tile stitcher matched arc endpoints only by lon/lat coincidence, which fails across rotated (e.g. X->Y) tile seams: there the two tiles' symmetric corners are offset by one cell and never coincide, so a region spanning such a seam fragmented into one loop per non-rotated block (e.g. the Atlantic split into east [1,2] and west [10,11] halves). Add a topology-aware fallback: when neither the exact corner nor a coincident corner continues the loop, join to an edge starting at a grid-ADJACENT corner via sectionate.build_neighbor_maps. This closes rotated-seam and Arctic-cap crossings. Single-tile and non-rotated multi-tile output is unchanged. Update the ECCO example to trace a real Atlantic basin: its perimeter is now a single loop spanning all five tiles it touches [1,2,6,10,11] (was fragmented). Add a real-data regression test for the rotated tile-2/tile-10 seam. Needs xgcm with face_connections padding fix (xgcm#713) for the cap/rotated seams, and sectionate.symmetrize for native ECCO 'left' staggering. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
Add a CI-runnable synthetic 2-tile grid with a rotated X->Y seam, built on the MITgcm 'left' corner and converted via sectionate.gridutils.symmetrize. The two tiles use disjoint coordinate blocks so the seam's corner representations never coincide in lon/lat -- a region spanning the seam can only be stitched through the grid topology (build_neighbor_maps), giving CI coverage of the rotated-seam join that previously was exercised only by the (data-gated) ECCO test. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
…titch
Replace the cell-based multi-tile tracer with the planned decomposition, which
isolates the only bespoke step (cross-seam stitching) and produces boundaries
that are grid-adjacent everywhere -- so they convert to velocity (u,v) faces and
obey the discrete flux-divergence theorem across rotated seams AND the Arctic cap.
Algorithm (regionate/boundaries.py):
1. contourpy traces each face's mask boundary (reusing the single-tile path).
2. a segment is internal (cut) iff BOTH cells it separates are in-mask in the
topology-aware padded mask -- the entire seam-cutting rule.
3. surviving open arcs are stitched by the global cell-set at their endpoints
(cells are globally unique, so no coordinate coincidence is needed).
4. face-local corners are converted to native (f,j,i) by a cell-set-matched
adjacency walk; seam crossings the cell-set over-merges are repaired by
inserting the corner where the two tiles meet; the loop is closed once.
Works for 'left' (MITgcm/ECCO), 'outer' (symmetric) and 'right' grids via
sectionate's corner_position/corner_offset; single-tile uses the same offset fix.
ECCO example now loads the native 'left' LLC90 grid directly (no symmetrize) and
uses the published regionmask North+South Atlantic basin (geographically correct,
no Pacific contamination). The basin traces to one 720-corner boundary across
six tiles and converts to 720 velocity faces. New real-data test asserts the
Atlantic boundary converts to velocity faces (the transport-consistency check).
Requires the topology-driven sectionate (MOM6-community/sectionate#47): native
'left' support + corner_position/corner_offset; and xgcm with the face_connections
padding fix (xgcm#713) for the cap.
Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
…atitude-circle section - atlantic_basin_mask: union the Atlantic-sector Natural Earth basins (the open Atlantic is split into several, notably the Sargasso Sea filling the central North Atlantic -- selecting only '*Atlantic*' left a large hole). Fills the basin, follows coastlines, excludes the Pacific. - Actually verify the discrete divergence theorem (do not just assert it): on a region crossing a non-rotated tile seam, the flux convergence summed over the cells equals the net flux through the traced boundary's velocity faces to round-off (residual ~1e-14). Added as a regression test, and demonstrated in the notebook (with an honest note that a cell-centred divergence across the LLC's rotated seams is subtle -- u rotates into v -- and is handled by sectionate.convergent_transport). - Notebook: distinct styling for boundary vs interior points; add a 26N latitude-circle section (curve='latitude circle') contrasted with great circle. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
…sion pins - examples/5 (ECCO LLC90): trace the full Atlantic basin, draw ALL boundary loops (so every masked cell is enclosed), and verify the discrete divergence theorem over the WHOLE basin (residual ~1e-14) using xgcm's vector-aware diff_2d_vector with an arbitrary transport field; remove the latitude-circle section. - boundaries.py: call the module-level xgcm.padding.pad directly, dropping the dead hasattr(grid,"pad") fallback (the pinned xgcm exposes only that API). - tests: add a real-ECCO Atlantic-basin divergence-theorem regression; compute cell-centred convergence with grid.diff_2d_vector (correct across rotated seams, where a per-component scalar grid.diff would be wrong -- expected xgcm behaviour). - deps (CI git refs + pyproject/meta comments): require the xgcm face_connections padding fix (hdrake/xgcm@fix-faceconnection-pad-712) and the topology-driven sectionate API + 'left'/rotated-seam corner-neighbour fix (sectionate#47), both needed for correct multi-tile (lat-lon-cap / cap) results. - docs: re-execute all example notebooks. The rotated-seam boundary tracing bug this exercises is fixed upstream in sectionate#47 (build_neighbor_maps corner-neighbour correction on 'left'/'right' grids); validated to machine precision on the Atlantic basin and the Arctic cap. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01E9R9Y9HxT5ZQyq6hKx4Hjm
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
1 participant
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Summary
Generalizes
regionatefrom single-tile MOM6 grids to arbitraryxgcm.Gridtopologies — single-tile periodic, bipolar-fold (Arctic) MOM6, and genuinely multi-tile grids (ECCOv4r4 LLC90, cubed-sphere) — by driving all grid logic from the topology-awaresectionate#47 API instead of hard-coded MOM6 specifics.(1) Boundaries from masks —
boundaries.py_contour_to_corner_indices; existing tests unchanged).boundary/face_connections); faces internal to a tile seam are dropped; surviving directed edges are stitched into closed loops by physical coincidence. A region spanning several faces yields a single boundary loop.f_c_list(Nonefor single-tile).(2) Masks from boundaries —
grid_conform.py+ newgeometry.pyregionmaskrasterization OR-ed together (modeled onxcryocouple.geometry; the intent is for xcryocouple to depend on regionate, not the reverse).get_geo_corners/get_geo_centers/coord_dict) — no hard-codedgeolon/geolat.(3) Refactor
f_cthreaded throughRegion/GriddedRegion/BoundedRegion,MaskRegions, and the.gr/.grssave format..grsave/load (variable-name mismatch + unguarded children); it now round-tripsf_c.Dependencies / CI
shapelyandxgcm >= 0.9.0.0.3.3, so the floor stays>= 0.3.3and CI installs sectionate from the PR ref (refs/pull/47/head) until it is released. This PR should not be released until sectionate#47 merges.Tests
14 passed, 2 skippedby default. New: self-contained two-face (face_connections) fixture exercising seam-spanning (one stitched loop across faces), seam-terminating (seam edge kept), interior-cell, andMaskRegionsf_cthreading; plus a pole-encircling regression test. Optional, skipped-by-default end-to-end checks against the real MOM6 global tripolar grid and ECCOv4r4 LLC90 (REGIONATE_REALDATA_TESTS=1); the LLC check asserts a reciprocal-or-raises invariant (canonical LLC90 is reciprocal; only the full cubed-sphere cap may raiseNotImplementedError).Known limitations / follow-ups
build_neighbor_mapsraisesNotImplementedErrorthere (upstream xgcm padding limitation) rather than miscomputing. ECCO's simply-rotated seams are supported.🤖 Generated with Claude Code