Skip to content

Instantly share code, notes, and snippets.

@mdsumner

mdsumner/myriahedral_connections.md

Created March 10, 2026 21:53
Show Gist options

  • Select an option

    Learn more about clone URLs
  • Save mdsumner/64d52743af1ffd5aae2affa07bdc4ad8 to your computer and use it in GitHub Desktop.

Select an option

Learn more about clone URLs
Save mdsumner/64d52743af1ffd5aae2affa07bdc4ad8 to your computer and use it in GitHub Desktop.
Download ZIP
Raw
myriahedral_connections.md

Myriahedral Projections: Implementations

Myriahedral projections (van Wijk, 2008) are not projections in the PROJ sense — there's no forward/inverse formula. They're a mesh-algorithm-projection pipeline:

  1. Build a fine triangular mesh on the sphere (icosahedral subdivision, graticule grid, or geography-adaptive)
  2. Assign edge weights (land/ocean crossings, graticule alignment, etc.)
  3. Compute a minimum spanning tree of the dual graph to decide where to cut
  4. Unfold the tree of faces flat using per-face gnomonic projections

The result has negligible area and angle distortion at the cost of many interrupts.

Implementations

temporaer/myriaworld — C++ with Google S2. The most complete implementation of the full pipeline: recursive icosahedral subdivision, edge weighting from Natural Earth shapefiles (via OGR), MST, and unfolding to SVG/PNG. Self-described as "still a mess"; depends on a modified S2 build. The author tried Boost.Geometry and GDAL before settling on S2 as the only library with adequate spherical mesh support.

jameslaneconkling/myriahedron — Pure JavaScript. Generates GeoJSON myriahedra by recursive icosahedral subdivision at arbitrary depth. More of a geodesic binning tool than a full unfolding projection. Outputs to GeoJSON/TopoJSON.

philogb WebGL demo — Interactive WebGL (luma.gl) demo by Nicolas Garcia Belmonte celebrating the 10th anniversary. Covers graticule-based meshes, Platonic/Archimedean solid subdivisions, and geography-aligned meshes with a live unfold slider. Source is embedded in the author's GitHub Pages repo, not published standalone.

d3-geo-projection #67 — Open feature request since 2016, no implementation.

Algorithm notes

The key insight is that the "projection" is really just per-face gnomonic applied after a graph-theoretic decision about mesh topology. The hard parts are spherical mesh construction, the spatial join of mesh edges against geography to assign weights, and the MST-based unfolding. Adjacent tools and libraries: RTriangle and decido for triangulation, rgl for mesh visualization, the s2 R package for spherical geometry, and PROJ for the per-face projection step.

Connections to DGGS and unstructured grids

The myriahedral approach shares deep structural DNA with Discrete Global Grid Systems (DGGS). Both start from the same premise: project the sphere onto a polyhedron (icosahedron, cube, rhombic dodecahedron), then recursively subdivide. The difference is that DGGS keeps the faces attached and uses them as an indexing/binning system, while myriahedral projections cut the dual graph open and flatten it. They're two answers to the same question about how to tessellate the sphere with minimal distortion.

Specific connections in the xarray/Pangeo ecosystem:

  • xdggs — Xarray extension providing a unified API across DGGS backends (H3, HEALPix, ISEA etc.), with cell-ID-based indexing. The mesh subdivision and cell hierarchy machinery is exactly what myriahedral projections need for step 1 (mesh construction), just without the cutting/unfolding.
  • UXarray — Xarray extension for unstructured grids (UGRID, MPAS, ICON, SCRIP). Stores full mesh topology: face-node connectivity, edge-face connectivity, dual meshes. MPAS uses icosahedral Voronoi meshes; ICON uses icosahedral triangular grids. The connectivity tables are exactly the graph structure needed for the MST step. UXarray's edge_face_connectivity is the dual graph adjacency list.
  • HEALPix — Equal-area pixelisation based on a rhombic dodecahedron (H=4 variant). The Calabretta & Roukema (2007) paper on HEALPix projections explicitly discusses unfolding the facets and notes the connection to Tegmark's icosahedral projection. ICON model output on HEALPix grids is accessed via hierarchical Zarr through intake catalogues. HEALPix is also a PROJ-supported projection (+proj=healpix).
  • GEOS Cubed-Sphere — Gnomonic projection onto 6 cube faces, each subdivided into a regular grid. UXarray supports this format. This is essentially a low-face-count myriahedral projection (6 faces, gnomonic per face) frozen into a standard grid.

The broader pattern is that the climate/weather modelling community has converged on the same polyhedra and subdivision schemes that van Wijk used, but for computational grids rather than cartographic display. The mesh topology data structures in UXarray (UGRID conventions) would serve as a natural interchange format for myriahedral projection algorithms.

References

@mdsumner
Copy link
Author

mdsumner commented Mar 10, 2026 •
edited
Loading

Myriahedral Projections in R: Feasibility Notes

Core insight

The algorithm decomposes into steps that already have good R implementations. No C++, JS, or Rust porting is required for a working package. The existing codebases (temporaer/myriaworld in C++/S2, jameslaneconkling/myriahedron in JS) are useful as algorithmic references but not as code to wrap or port.

Algorithm steps and R coverage

Step 1: Mesh construction on the sphere

Three mesh types from the van Wijk paper:

Icosahedral subdivision. The icosa package (Adam Kocsis) is the strongest existing R foundation here. icosa::trigrid(c(4,4)) gives you the triangular mesh directly, with @vertices (XYZ coordinates), @faces (which vertices form which faces), and @edges (vertex pairs) slots. Crucially, icosa::gridgraph() already builds an igraph face-adjacency graph — this is exactly the dual graph needed for step 3. The package uses Rcpp internally, handles both triangular (trigrid) and penta-hexagonal (hexagrid) grids, supports sf output via newsf(), and can do locate() for point-in-cell lookup. It also already imports igraph. The tessellation is controlled by a vector (e.g. c(2,3,4) for successive subdivision levels), giving finer control than simple recursive bisection. This package was optimised for coarser resolutions and gradual cell-size changes, making it a natural fit for myriahedral-scale meshes.

Alternatively, rgl::icosahedron3d() + rgl::subdivision3d() works, or ~50 lines of pure R index arithmetic (the jameslaneconkling JS code is a clean reference for the subdivision indexing).

Graticule-based meshes. Construct a lon/lat grid at some resolution (e.g. 5°), triangulate each quad into two triangles. This is trivial — no external dependency needed beyond base R. Vertex coordinates are just the graticule intersections on the sphere.

ISEA hexagonal grids (via dggridR). dggridR wraps Kevin Sahr's DGGRID C++ library via Rcpp, providing ISEA3H (aperture 3 hexagonal), ISEA4T (aperture 4 triangular), ISEA4D (diamond), and Fuller projection variants at resolutions 0–30. It gives you cell polygons (dgcellstogeo()), cell-to-SEQNUM lookup (dgGEO_to_SEQNUM()), and the ISEA/Fuller inverse icosahedral projections. The package is heavier (bundled C++ from DGGRID, proj4lib, shapelib, clipper) but gives access to the full zoo of DGGS grid types and the proper equal-area icosahedral projections that underpin them. Currently maintained by Sebastian Krantz on a fork from Richard Barnes' original. For myriahedral use, the key question is whether you can extract the face/edge topology from dggridR's cell polygons — the package is oriented around cell IDs and polygon generation rather than mesh connectivity, so there's an impedance mismatch. You'd likely need to reconstruct adjacency from the cell polygons, which is doable but less direct than icosa's built-in graph support.

Geography-adaptive meshes. Constrained triangulation incorporating coastline vertices. RTriangle or decido for the triangulation. This is the most complex mesh type and could be deferred to a later version.

Step 2: Edge weighting

For each edge in the mesh, assign a weight reflecting whether it should be cut (high weight = prefer to keep, low weight = prefer to cut). Weight sources:

  • Land/ocean crossing: does the edge cross a coastline? Use s2::s2_intersects() to test mesh edges against coastline polygons. Edges crossing land get high weight (don't cut there). This is the only step that needs the s2 package, and only for geography-aware meshes.
  • Graticule alignment: edges aligned with parallels or meridians get lower weight (prefer to cut along graticule lines). Pure arithmetic on vertex lon/lat.
  • Gaussian smoothing: myriaworld applies smoothing to edge weights to get cleaner cut patterns. Worth understanding from that codebase — it prevents noisy/jagged cuts.

Step 3: MST on the dual graph

The dual graph has one node per mesh face and an edge between faces that share a mesh edge. If using icosa, this is already done: icosa::gridgraph() returns an igraph object of the face adjacency. Then igraph::mst() on the weighted dual graph. Edges in the MST are folds (faces stay connected); edges not in the MST are cuts (faces separate). For non-icosa meshes, construct the dual from the face-node table by finding shared edges (index manipulation). igraph is mature and fast for this.

Step 4: Unfolding (the novel part)

Walk the spanning tree from a root face, laying each face flat in 2D:

  1. Pick a root face, project it to 2D via gnomonic projection centred on its centroid (the gnomonic formula is ~10 lines).
  2. For each child face in the tree, find the shared fold edge with its parent.
  3. Place the child face adjacent to its parent along that shared edge, computing the rotation needed to unfold it flat.
  4. Recurse through the tree.

The per-face gnomonic projection is trivial (available via reproj::reproj() or direct formula). The layout is a tree-walk with 2D rotations — maybe ~100 lines of careful coordinate geometry. The myriaworld C++ source is the best reference for this step, particularly for handling the rotation/placement logic.

Self-intersection caveat: for regular meshes (icosahedral, graticule) the unfolding is well-behaved. For geography-adaptive meshes with irregular topology, the layout can self-intersect. The graticule weight smoothing in step 2 mitigates this. A robust solution would need overlap detection, but that's an enhancement, not a blocker.

Proposed API

# Build a spherical mesh
mesh <- myriahedron_mesh(type = "icosahedral", depth = 5)
mesh <- myriahedron_mesh(type = "graticule", resolution = 5) # 5-degree

# Assign edge weights
mesh <- myriahedron_weight(mesh,
  coastline = NULL,        # sf polygon, optional
  graticule_weight = 3.0,  # weight for graticule-aligned edges
  sigma = 0.4              # Gaussian smoothing
)

# Compute cuts via MST
mesh <- myriahedron_cut(mesh)

# Unfold to 2D
unfolded <- myriahedron_unfold(mesh, root_face = 1L)
# Returns: 2D vertex coordinates + face table + face-ID linking to sphere

The output of myriahedron_unfold() is a 2D mesh: vertex xy coordinates, face indices (triangles), and a face-ID column linking each 2D face back to its spherical original. Directly plottable, and usable as input for warping/rasterization.

Dependencies

Package Role Required?
icosa Icosahedral mesh construction, face/edge/vertex tables, dual graph via gridgraph(), sf output, point-in-cell lookup Recommended (strongest foundation for step 1 + 3)
igraph MST on dual graph Yes (already an icosa dependency)
dggridR ISEA/Fuller DGGS grids, cell polygon generation, equal-area icosahedral projections Optional (alternative mesh source, heavier dependency)
s2 Coastline spatial join for edge weighting Only for geography-aware meshes
RTriangle Constrained triangulation Only for geography-adaptive meshes
decido Ear-clipping triangulation Alternative to RTriangle
rgl 3D visualisation, alternative mesh construction Optional
PROJ (or sf, terra project) Per-face gnomonic projection Optional (can use direct formula)

A minimal first version using icosa for mesh construction + graph gets igraph for free (it's already imported by icosa). That means the only hard dependency is icosa itself.

icosa vs dggridR for myriahedral use

Both packages build icosahedral grids, but they're oriented very differently:

icosa exposes the mesh topology directly. The trigrid S4 class gives you @vertices (XYZ), @faces (vertex indices), @edges (vertex pairs), @faceCenters (XYZ), and gridgraph() builds the igraph face-adjacency graph. Tessellation depth is controlled by a vector (e.g. c(2,3,4)) giving gradual resolution steps. It supports both triangular and penta-hexagonal grids, sf output, grid rotation, and locate() for point lookup. It uses Rcpp but the compiled code is small. The package was designed for coarse-to-medium resolution analysis (paleobiology, macroecology), which aligns well with myriahedral projections where you want maybe hundreds to low thousands of faces, not millions. This is the natural starting point for a myriahedral package — you get mesh + graph + sf in one dependency.

dggridR wraps Kevin Sahr's DGGRID C++ engine. It gives you the full ISEA and Fuller projection families (ISEA3H, ISEA4T, ISEA4D, ISEA43H, FULLER3H, FULLER4T, etc.) at up to resolution 30 (billions of cells). The API is cell-ID-oriented: you go from lon/lat → SEQNUM → cell polygon. It does not directly expose mesh topology (face-node connectivity, edge-face adjacency). To use it for myriahedral projections, you'd need to reconstruct adjacency from the cell polygons, which is feasible but adds a layer of work. The package also bundles its own proj4lib and shapelib, making it heavier. The ISEA equal-area icosahedral projection that DGGRID implements is interesting in its own right — it's the projection applied per icosahedral face, which is closely related to the per-face gnomonic used in myriahedral unfolding. Now maintained by Sebastian Krantz (v3.1.0, imports sf, s2, collapse).

Recommendation: Start with icosa::trigrid() for mesh construction. If you later want to connect to DGGS standards (ISEA3H cell IDs, hierarchical resolution levels, interop with xdggs/UXarray), dggridR provides the bridge, but it's not needed for the core myriahedral algorithm. The deeper opportunity with dggridR is that its bundled ISEA projection code could potentially be extracted to provide the per-face projection step — but that's a future investigation, not a blocker.

What to read from existing code

temporaer/myriaworld (C++/S2, 51 commits): Read for the edge weighting heuristics (graticular weight function, land/ocean classification, Gaussian smoothing) and the unfold tree-walk geometry. The S2 and caching infrastructure is not portable but the algorithmic logic is clear. The --depth, --sigma, --wlon parameters map directly to the proposed API.

jameslaneconkling/myriahedron (JS, 57 commits): Read for the icosahedral subdivision indexing only. Clean implementation of recursive subdivision → GeoJSON. Does not do cutting or unfolding.

philogb WebGL demo: Read the Medium article for intuition about the different mesh types and their visual properties. Source not separately published.

Van Wijk (2008) paper: The essential reference. Describes the weight function for graticule meshes, the MST algorithm on the dual graph, and the three mesh types. Available at https://www.win.tue.nl/~vanwijk/myriahedral/.

Connections to DGGS / xarray ecosystem

DGGS and myriahedral projections share the same substrate: project the sphere onto a polyhedron, recursively subdivide. DGGS keeps faces attached for indexing; myriahedral projections cut and flatten. Relevant projects:

  • xdggs — Xarray extension for DGGS (H3, HEALPix, ISEA). The mesh hierarchy is the same as step 1.
  • UXarray — Xarray for unstructured grids (UGRID, MPAS, ICON). Stores face-node and edge-face connectivity tables — literally the dual graph adjacency needed for step 3.
  • HEALPix — Equal-area pixelisation on a rhombic dodecahedron. Has a proper PROJ projection (+proj=healpix), bridging the DGGS/cartographic worlds. ICON climate model outputs on HEALPix grids via hierarchical Zarr.
  • GEOS Cubed-Sphere — Gnomonic projection onto 6 cube faces, each with a regular grid. A 6-face myriahedral projection frozen into a standard format.

Estimation

A working package with icosahedral and graticule mesh types, MST cutting, and 2D unfolding is probably 500–800 lines of R. Geography-adaptive meshes with coastline-aware weighting add maybe 200 more. The main risk is getting the unfold tree-walk geometry right (rotation placement of child faces), which is fiddly but mathematically straightforward.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment