Regridding API

These are the offline conservative-regridding construction, persistence, and application helpers. See Regridding for the physical contract.

• build_regridder returns an IdentityRegrid for equivalent meshes, a NestedCubedSphereRegrid (Restriction or Prolongation) for nested cubed spheres, and a ConservativeRegridding.jl Regridder otherwise.

• apply_regridder!(dst, r, src) maps an intensive field; apply_regridder!(dst, r, src, kind) takes the quantity kind (IntensiveCellField(), ExtensiveCellField() for per-cell totals), regrid_conservation_report compares the global totals, and source_coverage gives the covered area of each source cell.

• coarsen_sum! and coarsen_area_weighted! are the cubed-sphere block kernels of the transport-binary coarseners.

• default_regrid_cache_dir is the shared weight cache.

AtmosTransport.Regridding Module

Regridding (src)

Offline conservative regridding between mesh types, built on ConservativeRegridding.jl.

Build a sparse weights matrix once per (source_mesh, target_mesh) pair, optionally cache it to disk, and reuse it. Equivalent meshes get an IdentityRegrid and nested cubed spheres (one definition, integer Nc ratio) an exact NestedCubedSphereRegrid, neither of which needs weights. Fields are regridded by quantity kind (apply_regridder!(dst, r, src, kind)): intensive fields directly, per-cell totals through their density. The module also owns the cubed-sphere block kernels (coarsen_sum!, coarsen_area_weighted!) that the transport-binary coarseners use. Callers are Preprocessing (met data to target grids), Models.InitialConditionIO (initial conditions and surface-flux inventories at run setup) and Visualization (cubed-sphere to lon-lat rasters). The operators and the stepping loop never call into this module — TransportBinaryReader consumes binaries that are already on the target grid.

Workflow

using AtmosTransport: LatLonMesh, CubedSphereMesh
using AtmosTransport.Quantities: ExtensiveCellField
using AtmosTransport.Regridding

src = LatLonMesh(Nx=1440, Ny=721)
dst = CubedSphereMesh(Nc=90)

# Build weights (expensive — minutes for C90, cached to disk)
r = build_regridder(src, dst; cache_dir="/tmp/atmos_regrid_cache")

# Apply to a 2D field (intensive: a value per area)
src_field = zeros(src.Nx * src.Ny)
dst_field = zeros(6 * dst.Nc * dst.Nc)
apply_regridder!(dst_field, r, src_field)

# Per-cell totals (kg per cell) are regridded through their density
apply_regridder!(dst_field, r, src_field, ExtensiveCellField())
regrid_conservation_report(r, src_field, dst_field, ExtensiveCellField())

# Persist / export
save_regridder("weights.jld2", r)
save_esmf_weights("weights_esmf.nc", r)

Supported mesh types

Source/Destination	Tree strategy	Spatial acceleration
LatLonMesh	CellBasedGrid from (Nx+1)×(Ny+1) face corners → TopDownQuadtreeCursor	O(log(Nx·Ny))
CubedSphereMesh	6 per-panel convention-aware grids from CS corners → CubedSphereToplevelTree	O(log Nc²) per panel
ReducedGaussianMesh	Sectorized per-ring CellBasedGrid trees → MultiTreeWrapper	O(nrings · log nlon)

All three produce SphericalCap extents at every tree level, which is required by CR.jl's spherical dual-DFS intersection search.

Known limitations

• CubedSphereMesh uses analytical coordinates from src/Grids/cs_mesh_coordinates.jl. GEOSNativePanelConvention includes the GEOS-FP/GEOS-IT panel order, native orientation, and global -10° longitude offset used by GEOS grid files. Left-handed GEOS panels are wound correctly for tree traversal while preserving file-order indices.

• ReducedGaussianMesh clamps polar face latitudes by 0.001° to avoid degenerate polygons at the poles. Also, its cell polygons use great-circle edges to approximate curved latitude boundaries. This approximation is visible on unusually coarse synthetic rings and converges as ring resolution increases. Where adjacent rings have different longitude counts the two great-circle polylines of their shared latitude leave a sliver that no cell covers (6.3e-4 of the sphere on O24); source_coverage gives the covered part of each source cell.

Architecture reference

See docs/src/preprocessing/regridding.md for the public API, the weights cache, and the mass-consistency correction applied after regridding, and src/Regridding/README.md for the file map and invariants.

source

AtmosTransport.Regridding.AbstractNestingDirection Type

AbstractNestingDirection

Direction of a NestedCubedSphereRegrid: Restriction (fine → coarse) or Prolongation (coarse → fine).

source

AtmosTransport.Regridding.IdentityRegrid Type

IdentityRegrid{M}

Passthrough sentinel returned by build_regridder when source and target meshes are structurally equivalent (per meshes_equivalent). Carries the shared mesh and its cell areas [m²] as src_areas and dst_areas, in the flat cell order of apply_regridder!, so callers that convert between densities and cell totals use it as they use a conservative Regridder. It has no intersection matrix: it cannot be saved as weights, and normalize does not scale its areas.

source

AtmosTransport.Regridding.NestedCubedSphereRegrid Type

NestedCubedSphereRegrid{D<:AbstractNestingDirection}(ratio, src_areas, dst_areas)
NestedCubedSphereRegrid(src_mesh::CubedSphereMesh, dst_mesh::CubedSphereMesh)

Exact conservative remap between nested cubed spheres (meshes_nested), returned by build_regridder for such a pair. ratio is the integer ratio R of the edge sizes; src_areas and dst_areas are the cell areas [m²] of the two meshes in the flat cell order of apply_regridder!, (p − 1)Nc² + i + (j − 1)Nc. The mesh constructor evaluates them in Float64 from the mesh corners (Girard's theorem), whatever the mesh precision, so Float32 and Float64 meshes give the same regridder.

apply_regridder!(dst, r, src) remaps an intensive field (per area or per point) as a conservative regridder with intersection matrix W does, $d_t=\sum _s{W}_{ts}{s}_s/A_t^{\mathrm{d}\mathrm{s}\mathrm{t}}$. For nested meshes W needs no polygon clipping: fine cell f lies in coarse cell c(f).

• Restriction (fine → coarse), $W_{cf}=A_f$:

  $$d_c=\frac{1}{A_c}\sum _{f\in c}{A}_f{s}_f$$

• Prolongation (coarse → fine), $W_{fc}=A_c{A}_f/\sum _{f^{\prime }\in c}{A}_{f^{\prime }}$:

  $$d_f=s_c\frac{A_c}{\sum _{f^{\prime }\in c}{A}_{f^{\prime }}}$$

Both satisfy $\sum _t{W}_{ts}=A_s^{\mathrm{s}\mathrm{r}\mathrm{c}}$, so the extensive path of apply_regridder! (ExtensiveCellField: divide by src_areas, remap, multiply by dst_areas) conserves the total exactly up to rounding: a block sum for restriction, a split by area fraction $m_f=m_c{A}_f/\sum A_{f^{\prime }}$ for prolongation. The coarse mesh area and the block sum of fine areas differ only by round-off (6.3e-13 relative at C180/C90 on the GMAO grid), so a uniform intensive field stays uniform to that level, and the result equals a CR.jl Regridder of the same pair to about 1e-12. Sums run in Float64 in a fixed order (ii fastest within a block); the result is rounded once to the element type of dst.

Nothing is cached on disk: the regridder is built from the meshes in O(Nc²). It has no intersection matrix, so it cannot be saved as weights, and normalize does not scale its areas.

source

AtmosTransport.Regridding.Prolongation Type

Prolongation <: AbstractNestingDirection

Coarse → fine remap between nested cubed spheres: each coarse cell is split among its R × R fine cells by area fraction.

source

AtmosTransport.Regridding.Restriction Type

Restriction <: AbstractNestingDirection

Fine → coarse remap between nested cubed spheres: area-weighted block sum.

source

AtmosTransport.Regridding.apply_regridder! Method

apply_regridder!(dst, r, src, kind::QuantityKind; scratch = nothing) -> dst

Regrid src onto dst with regridder r according to the kind of field (see Quantities). Axis 1 of dst and src is the flat horizontal index; trailing axes (levels, times, …) are regridded column by column.

• IntensiveCellField(): a value per area or per point (ps, T, mixing ratios, flux densities). Same as apply_regridder!(dst, r, src): $d_t=\sum _s{W}_{ts}{s}_s/A_t^{\mathrm{d}\mathrm{s}\mathrm{t}}$.

• ExtensiveCellField(): a per-cell total (m in kg, a flux in kg s⁻¹ per cell). Regridded as a density and multiplied back,

  $$d_t=A_t^{\mathrm{d}\mathrm{s}\mathrm{t}}\mathcal{R}{\left(\frac{s}{A^{\mathrm{s}\mathrm{r}\mathrm{c}}}\right)}_t=\sum _s\frac{W_{ts}}{A_s^{\mathrm{s}\mathrm{r}\mathrm{c}}}{s}_s,$$

  with A^src = r.src_areas and A^dst = r.dst_areas. Since $\sum _t{W}_{ts}=A_s^{\mathrm{s}\mathrm{r}\mathrm{c}}$ for a target that covers the source, $\sum _t{d}_t=\sum _s{s}_s$ whatever polygon areas the regridder uses. Regridding a per-cell total as intensive instead distorts it in proportion to the source cell area (the ~12× polar mass deficit of early C180 ERA5 binaries).

• HorizontalVectorField(), HorizontalFluxField(): not scalar regrids; throw an ArgumentError. Winds are regridded component by component as intensive and rotated (rotate_winds_to_panel_local!); face fluxes are rebuilt from regridded winds (reconstruct_cs_fluxes!).

src is not modified. The density $s/A^{\mathrm{s}\mathrm{r}\mathrm{c}}$ of one horizontal column at a time is formed in the element type of src and stored in scratch; the product with $A^{\mathrm{d}\mathrm{s}\mathrm{t}}$ is stored in the element type of dst. Both steps are broadcasts, not scalar loops.

scratch is a vector of length length(r.src_areas) with element type eltype(src) that does not share memory with src or dst; its contents are overwritten. With nothing (the default) an ExtensiveCellField regrid allocates one such vector per call, similar(src, n_src). Callers that regrid repeatedly pass a reusable buffer, as the cubed-sphere preprocessing workspace does (density_scratch), so its steady state allocates nothing. The other kinds accept and ignore scratch, so one call site serves every kind.

source

AtmosTransport.Regridding.apply_regridder! Method

apply_regridder!(dst, regridder, src) -> dst

Apply a regridder to field data. Thin wrapper over ConservativeRegridding.regrid! that handles n-dimensional arrays whose first dimension is the flattened horizontal index with any number of trailing dimensions (levels, time windows, tracers, etc).

dst and src may be:

• plain Vectors of length n_dst / n_src, or

• multi-dimensional arrays with a leading horizontal axis of length n_dst / n_src and matching trailing shape.

For the multi-dim case the horizontal regridding is applied column-by-column using CR.jl's temporary work arrays — no extra allocations per column.

The field is treated as intensive (a value per area or per point). Per-cell totals go through apply_regridder!(dst, r, src, ExtensiveCellField()) (quantity_regrid.jl). IdentityRegrid and NestedCubedSphereRegrid have their own three-argument methods with the same contract.

source

AtmosTransport.Regridding.build_regridder Method

build_regridder(src_mesh, dst_mesh; normalize=false, cache_dir=nothing, kwargs...)
    -> IdentityRegrid | NestedCubedSphereRegrid | Regridder

Construct the regridder that maps source-mesh fields to destination-mesh fields, choosing in this order:

1. IdentityRegrid when meshes_equivalent(src, dst);

2. NestedCubedSphereRegrid when meshes_nested(src, dst) (cubed spheres of one definition whose Nc differ by an integer factor): exact block remapping, no polygon intersection;

3. otherwise a ConservativeRegridding.Regridder built by true spherical polygon intersection.

The first two take no weights from or to disk and ignore normalize and kwargs. If cache_dir is supplied, a conservative Regridder is cached on disk using a SHA-1 fingerprint of the (src_mesh, dst_mesh) pair (see default_regrid_cache_dir for the shared default). Subsequent calls with identical meshes load the cached weights instead of recomputing them. Custom construction kwargs cannot be combined with cache_dir, because arbitrary operators do not have a stable, portable cache fingerprint.

normalize=false matches the convention used by xESMF: the sparse intersections matrix contains raw area values (m²), and the per-cell normalization by dst_areas happens inside regrid!. Pass normalize=true to scale intersections, src_areas, and dst_areas by max(intersections) — useful for display but not for weight-level comparisons.

Any extra kwargs are forwarded to ConservativeRegridding.Regridder.

source

AtmosTransport.Regridding.coarsen_area_weighted! Method

coarsen_area_weighted!(dst, src, area, ratio::Int) -> dst

Area-weighted block mean of the fine panel array src onto the coarse panel array dst, for an intensive quantity (K, m s⁻¹, kg kg⁻¹, kg m⁻² s⁻¹, …):

$$d_{i,j,k}={\mathrm{fl}}_{FT}\left(\frac{\sum _{f\in B(i,j)}{A}_f{s}_{f,k}}{\sum _{f\in B(i,j)}{A}_f}\right),B(i,j)=\{((i-1)R+ii,(j-1)R+jj):1\le ii,jj\le R\},$$

with R = ratio and A = area, the (R·Nc, R·Nc) fine cell areas of one panel (every panel has the same areas). Numerator and denominator accumulate in Float64 (ii fastest) and the quotient is rounded once to FT = eltype(dst), so a Float32 result is the correctly rounded block mean up to the Float64 accumulation error, and a uniform field stays uniform. area may be Float64 while the data are Float32. dst/src are both 2-D (Nc, Nc) or both 3-D (Nc, Nc, Nz); the 3-D result equals the 2-D kernel applied level by level, bit for bit.

The denominator is the block sum of fine areas, not the coarse mesh area: nested meshes agree to round-off only (relative 6.3e-13 at C180 → C90 on the GMAO grid), and the block sum keeps constants exact. NestedCubedSphereRegrid divides by the coarse mesh area instead, so that its extensive path conserves exactly.

source

AtmosTransport.Regridding.coarsen_sum! Method

coarsen_sum!(dst, src, ratio::Int) -> dst

Block-sum the fine panel array src onto the coarse panel array dst, for an extensive quantity (kg per cell, kg per substep, …):

$$d_{i,j,k}=\sum _{jj=1}^R\sum _{ii=1}^R{s}_{(i-1)R+ii,(j-1)R+jj,k}$$

with R = ratio. dst and src are both 2-D (Nc, Nc) or both 3-D (Nc, Nc, Nz) with size(src) = (R·Nc, R·Nc[, Nz]), of the same element type FT. The sum accumulates in FT, ii fastest, which is the order of the two earlier kernel copies (Preprocessing's nested binary coarsener and GEOS block coarsener), so their outputs are unchanged bit for bit.

See also coarsen_area_weighted! for intensive fields.

source

AtmosTransport.Regridding.cubed_sphere_face_corners Method

cubed_sphere_face_corners(mesh::CubedSphereMesh)
    -> NTuple{6, Matrix{UnitSphericalPoint{Float64}}}

Return a 6-tuple of (Nc+1) × (Nc+1) corner-point matrices, one per panel, with corners on the unit sphere. Panel ordering follows mesh.convention; for each user-visible panel p, corner (i, j) matches panel_cell_corner_lonlat(mesh, p, Float64). Corners are evaluated in Float64 for every mesh precision, so a Float32 run regrids with the same weights as a Float64 run (the cache key does not include the precision).

Used by Trees.treeify(::Spherical, ::CubedSphereMesh) to assemble a CubedSphereToplevelTree. Cached externally (the regridder cache includes the mesh hash) so this is only called once per (Nc, convention).

source

AtmosTransport.Regridding.default_regrid_cache_dir Method

default_regrid_cache_dir() -> String

Directory of the shared conservative-regridding weight cache: the value of ATMOSTR_REGRID_CACHE_DIR if set and not empty (the golden-output harness points it into each case, so weights are recomputed), else ~/.cache/AtmosTransport/cr_regridding, with ~ expanded. Read at each call, so the precompiled package carries no home directory. Runtime surface fluxes pass it as cache_dir to build_regridder.

source

AtmosTransport.Regridding.load_regridder Method

load_regridder(path) -> Regridder

Inverse of save_regridder. Rebuilds the temporary work arrays from the area vectors.

source

AtmosTransport.Regridding.meshes_equivalent Method

meshes_equivalent(src, dst) -> Bool

True when the two meshes describe the same horizontal grid up to floating-point equality of the geometry-defining fields. Used by build_regridder to skip building a conservative regridder when the source and target are the same grid.

source

AtmosTransport.Regridding.meshes_nested Method

meshes_nested(a, b) -> Bool

True when a and b are cubed-sphere meshes with the same definition (tag, coordinate law, centre law, panel convention, longitude offset) and radius, different edge sizes Nc with an integer ratio, and corners that coincide: corner (i, j) of each panel of the coarser mesh is corner ((i − 1)R + 1, (j − 1)R + 1) of the finer one, to 1e-12 on the unit sphere (evaluated, not assumed). Every fine cell then lies in exactly one coarse cell. Meshes with equal Nc are equivalent (meshes_equivalent), not nested; meshes of other types are never nested.

source

AtmosTransport.Regridding.regrid_conservation_report Method

regrid_conservation_report(r, src, dst, kind::QuantityKind) -> (; src_total, dst_total, rel_error)

Global totals of a field before (src) and after (dst) regridding with r, and the relative error $|T_{\mathrm{d}\mathrm{s}\mathrm{t}}-T_{\mathrm{s}\mathrm{r}\mathrm{c}}|/|T_{\mathrm{s}\mathrm{r}\mathrm{c}}|$ (0 when both totals are 0, Inf when only src_total is). Axis 1 is the flat horizontal index, of length length(r.src_areas) for src and length(r.dst_areas) for dst (else a DimensionMismatch, for both kinds); totals include every trailing column.

• ExtensiveCellField(): $T=\sum _c{x}_c$.

• IntensiveCellField(): $T=\sum _c{x}_c{A}_c$, with r.src_areas for src and r.dst_areas for dst.

Sums are taken in Float64 in a fixed order. Vector and flux kinds throw an ArgumentError.

source

AtmosTransport.Regridding.save_esmf_weights Method

save_esmf_weights(path, regridder;
                  src_shape=nothing, dst_shape=nothing,
                  src_grid_name="source", dst_grid_name="destination") -> path

Export the regridder weights to an ESMF offline-weights NetCDF file.

Format (ESMF convention):

variable	shape	meaning
S	(n_s,)	weight value S[k] = intersections[row[k], col[k]] / dst_areas[row[k]]
row	(n_s,)	destination cell index (1-indexed)
col	(n_s,)	source cell index (1-indexed)
frac_a	(n_a,)	fraction of source cell covered by the destination grid
frac_b	(n_b,)	fraction of destination cell covered by the source grid
area_a	(n_a,)	source cell areas (mesh units, typically m²)
area_b	(n_b,)	destination cell areas

This format is what xESMF / ESMF_RegridWeightGen / GCHP's tile-file loader all consume, which lets the CR.jl output be diffed against xESMF fixtures at weight level and handed off to ESMF-based tooling.

Optional src_shape / dst_shape tuples (e.g. (Nx, Ny), (Nc, Nc, 6)) are stored as global attributes for provenance.

source

AtmosTransport.Regridding.save_regridder Method

save_regridder(path, regridder) -> path

Write the regridder's intersection matrix and area vectors to path as a JLD2 file. Temporary work arrays are not serialized (they're rebuilt on load). The file is published with AtomicFiles.atomic_publish: a reader sees the previous file or the complete new one, and the file gets mode 0o666 & ~umask like any file the user creates, so a shared cache stays readable by the group.

source

AtmosTransport.Regridding.source_coverage Method

source_coverage(r) -> Vector{Float64}

The area of each source cell that the destination cells of r cover, $C_s=\sum _t{W}_{ts}$, in the flat source order of r.src_areas. When the destination tiles the source, $C_s=A_s^{\mathrm{s}\mathrm{r}\mathrm{c}}$ to the round-off of the polygon intersections, and an ExtensiveCellField() regrid keeps $\sum _s{s}_s$. IdentityRegrid and NestedCubedSphereRegrid cover every source cell exactly and return r.src_areas.

A ReducedGaussianMesh destination does not tile the sphere where adjacent rings have different longitude counts: its cell edges are great circles, so the latitude two such rings share is two different polylines, and the sliver between them belongs to no destination cell. There $C_s<A_s^{\mathrm{s}\mathrm{r}\mathrm{c}}$ (by 17 % for 5° source cells on the boundary of a 12-cell and a 24-cell ring; the octahedral grids change the count at every ring). A per-cell total divided by $C_s$ instead of $A_s^{\mathrm{s}\mathrm{r}\mathrm{c}}$ before an intensive regrid, and multiplied by r.dst_areas after it, keeps its total on such destinations: $\sum _t\sum _s{W}_{ts}{s}_s/C_s=\sum _s{s}_s$.

source