MeshGlyph — Unstructured Mesh Visualization¶

This notebook demonstrates cleopatra.mesh_glyph.MeshGlyph, a class for visualizing UGRID-style unstructured mesh data using matplotlib.

Features covered:

1. Building a triangular mesh and plotting face-centered data
2. Node-centered (interpolated contour) plotting
3. Wireframe / outline rendering
4. Mixed-element meshes (triangles + quads)
5. Overlaying data on a wireframe
6. Generating a realistic river-channel mesh

In [ ]:

import matplotlib.pyplot as plt
import matplotlib.tri as mtri
import numpy as np

from cleopatra.mesh_glyph import MeshGlyph

import matplotlib.pyplot as plt import matplotlib.tri as mtri import numpy as np from cleopatra.mesh_glyph import MeshGlyph

1. Simple Triangular Mesh — Face-Centered Data¶

Create a small triangular mesh from random points using Delaunay triangulation, then color each face by its centroid elevation.

In [ ]:

rng = np.random.default_rng(42)
n_points = 80
node_x = rng.uniform(0, 10, n_points)
node_y = rng.uniform(0, 8, n_points)

# matplotlib's Triangulation computes Delaunay internally
tri = mtri.Triangulation(node_x, node_y)
face_nodes = tri.triangles  # (n_faces, 3)

# Compute a synthetic "water depth" at each face centroid
cx = node_x[face_nodes].mean(axis=1)
cy = node_y[face_nodes].mean(axis=1)
face_depth = np.sin(cx * 0.6) * np.cos(cy * 0.5) + 1.5

mg = MeshGlyph(node_x, node_y, face_nodes)
print(f"Mesh: {mg.n_nodes} nodes, {mg.n_faces} faces")

rng = np.random.default_rng(42) n_points = 80 node_x = rng.uniform(0, 10, n_points) node_y = rng.uniform(0, 8, n_points) # matplotlib's Triangulation computes Delaunay internally tri = mtri.Triangulation(node_x, node_y) face_nodes = tri.triangles # (n_faces, 3) # Compute a synthetic "water depth" at each face centroid cx = node_x[face_nodes].mean(axis=1) cy = node_y[face_nodes].mean(axis=1) face_depth = np.sin(cx * 0.6) * np.cos(cy * 0.5) + 1.5 mg = MeshGlyph(node_x, node_y, face_nodes) print(f"Mesh: {mg.n_nodes} nodes, {mg.n_faces} faces")

In [ ]:

fig, ax = mg.plot(
    face_depth,
    location="face",
    cmap="RdYlBu_r",
    title="Face-Centered Water Depth",
    edgecolor="grey",
    figsize=(10, 7),
)
ax.set_xlabel("x [m]")
ax.set_ylabel("y [m]")
plt.tight_layout()
plt.show()

fig, ax = mg.plot( face_depth, location="face", cmap="RdYlBu_r", title="Face-Centered Water Depth", edgecolor="grey", figsize=(10, 7), ) ax.set_xlabel("x [m]") ax.set_ylabel("y [m]") plt.tight_layout() plt.show()

2. Node-Centered Data — Smooth Contours¶

When data is defined at nodes rather than faces, MeshGlyph.plot uses tricontourf for smooth interpolated contours.

In [ ]:

# Synthetic elevation field at each node
node_elevation = np.sin(node_x * 0.5) * np.cos(node_y * 0.4) * 3.0 + 5.0

fig, ax = mg.plot(
    node_elevation,
    location="node",
    cmap="terrain",
    levels=15,
    title="Node-Centered Elevation (Contour Fill)",
    figsize=(10, 7),
)
ax.set_xlabel("x [m]")
ax.set_ylabel("y [m]")
plt.tight_layout()
plt.show()

# Synthetic elevation field at each node node_elevation = np.sin(node_x * 0.5) * np.cos(node_y * 0.4) * 3.0 + 5.0 fig, ax = mg.plot( node_elevation, location="node", cmap="terrain", levels=15, title="Node-Centered Elevation (Contour Fill)", figsize=(10, 7), ) ax.set_xlabel("x [m]") ax.set_ylabel("y [m]") plt.tight_layout() plt.show()

3. Wireframe Outline¶

plot_outline renders the mesh edges as a lightweight wireframe, useful for inspecting mesh quality or overlaying on top of data plots.

In [ ]:

fig, ax = mg.plot_outline(color="steelblue", linewidth=0.5, figsize=(10, 7))
ax.plot(node_x, node_y, "k.", markersize=3)
ax.set_title("Mesh Wireframe with Node Locations")
ax.set_xlabel("x [m]")
ax.set_ylabel("y [m]")
plt.tight_layout()
plt.show()

fig, ax = mg.plot_outline(color="steelblue", linewidth=0.5, figsize=(10, 7)) ax.plot(node_x, node_y, "k.", markersize=3) ax.set_title("Mesh Wireframe with Node Locations") ax.set_xlabel("x [m]") ax.set_ylabel("y [m]") plt.tight_layout() plt.show()

4. Overlay — Data + Wireframe on the Same Axes¶

Pass an existing ax to combine a data plot with a wireframe overlay.

In [ ]:

fig, ax = mg.plot(
    face_depth,
    location="face",
    cmap="Blues",
    colorbar=True,
    title="Water Depth with Wireframe Overlay",
    figsize=(10, 7),
)
# Overlay the wireframe on the same axes
mg.plot_outline(ax=ax, color="black", linewidth=0.3)
ax.set_xlabel("x [m]")
ax.set_ylabel("y [m]")
plt.tight_layout()
plt.show()

fig, ax = mg.plot( face_depth, location="face", cmap="Blues", colorbar=True, title="Water Depth with Wireframe Overlay", figsize=(10, 7), ) # Overlay the wireframe on the same axes mg.plot_outline(ax=ax, color="black", linewidth=0.3) ax.set_xlabel("x [m]") ax.set_ylabel("y [m]") plt.tight_layout() plt.show()

5. Mixed-Element Mesh (Triangles + Quads)¶

Many real-world UGRID meshes combine quadrilateral and triangular elements. MeshGlyph handles these via fan triangulation — quads are decomposed into 2 triangles internally. Use fill_value=-1 to pad rows for faces with fewer nodes.

In [ ]:

# Build a structured quad grid (3x2 cells) with triangular transitions
#
#   6---7---8---9
#   |   |   | / |
#   3---4---5   |
#   |   |   | \ |
#   0---1---2--10
#
mixed_x = np.array([0, 1, 2, 0, 1, 2, 0, 1, 2, 3, 3], dtype=float)
mixed_y = np.array([0, 0, 0, 1, 1, 1, 2, 2, 2, 2, 0], dtype=float)

# 4 quads + 2 triangles — padded with -1
face_conn = np.array(
    [
        [0, 1, 4, 3],   # quad
        [1, 2, 5, 4],   # quad
        [3, 4, 7, 6],   # quad
        [4, 5, 8, 7],   # quad
        [2, 10, 5, -1],  # triangle (bottom-right)
        [5, 10, 9, -1],  # triangle (right side)
    ]
)
face_values = np.array([1.0, 2.0, 3.0, 4.0, 2.5, 3.5])

mg_mixed = MeshGlyph(mixed_x, mixed_y, face_conn, fill_value=-1)
print(f"Mixed mesh: {mg_mixed.n_nodes} nodes, {mg_mixed.n_faces} faces")
print(f"Nodes per face: {mg_mixed.nodes_per_face}")

# Build a structured quad grid (3x2 cells) with triangular transitions # # 6---7---8---9 # | | | / | # 3---4---5 | # | | | \ | # 0---1---2--10 # mixed_x = np.array([0, 1, 2, 0, 1, 2, 0, 1, 2, 3, 3], dtype=float) mixed_y = np.array([0, 0, 0, 1, 1, 1, 2, 2, 2, 2, 0], dtype=float) # 4 quads + 2 triangles — padded with -1 face_conn = np.array( [ [0, 1, 4, 3], # quad [1, 2, 5, 4], # quad [3, 4, 7, 6], # quad [4, 5, 8, 7], # quad [2, 10, 5, -1], # triangle (bottom-right) [5, 10, 9, -1], # triangle (right side) ] ) face_values = np.array([1.0, 2.0, 3.0, 4.0, 2.5, 3.5]) mg_mixed = MeshGlyph(mixed_x, mixed_y, face_conn, fill_value=-1) print(f"Mixed mesh: {mg_mixed.n_nodes} nodes, {mg_mixed.n_faces} faces") print(f"Nodes per face: {mg_mixed.nodes_per_face}")

In [ ]:

fig, axes = plt.subplots(1, 2, figsize=(14, 5))

# Left: face data
mg_mixed.plot(face_values, location="face", ax=axes[0], cmap="coolwarm",
              edgecolor="black", title="Mixed Mesh — Face Data")
axes[0].set_xlabel("x")
axes[0].set_ylabel("y")

# Right: wireframe only
mg_mixed.plot_outline(ax=axes[1], color="darkgreen", linewidth=1.5)
for i, (x, y) in enumerate(zip(mixed_x, mixed_y)):
    axes[1].annotate(str(i), (x, y), fontsize=9, ha="center", va="bottom",
                     fontweight="bold", color="darkgreen")
axes[1].set_title("Mixed Mesh — Wireframe with Node IDs")
axes[1].set_xlabel("x")
axes[1].set_ylabel("y")

plt.tight_layout()
plt.show()

fig, axes = plt.subplots(1, 2, figsize=(14, 5)) # Left: face data mg_mixed.plot(face_values, location="face", ax=axes[0], cmap="coolwarm", edgecolor="black", title="Mixed Mesh — Face Data") axes[0].set_xlabel("x") axes[0].set_ylabel("y") # Right: wireframe only mg_mixed.plot_outline(ax=axes[1], color="darkgreen", linewidth=1.5) for i, (x, y) in enumerate(zip(mixed_x, mixed_y)): axes[1].annotate(str(i), (x, y), fontsize=9, ha="center", va="bottom", fontweight="bold", color="darkgreen") axes[1].set_title("Mixed Mesh — Wireframe with Node IDs") axes[1].set_xlabel("x") axes[1].set_ylabel("y") plt.tight_layout() plt.show()

6. Realistic Example — River Channel Bathymetry¶

Generate a curved river-channel mesh using Delaunay triangulation on a set of points that follow a sinusoidal centreline. The face data represents simulated bed elevation.

In [ ]:

rng = np.random.default_rng(123)

# River centreline: sinusoidal path
s = np.linspace(0, 4 * np.pi, 200)
cx_line = s * 50 / (4 * np.pi)  # 0..50 m along x
cy_line = 4 * np.sin(s)         # sinusoidal meander

# Scatter points around the centreline with a channel half-width of ~3 m
n_river = 600
idx = rng.integers(0, len(s), n_river)
offsets = rng.normal(0, 1.2, n_river)
rx = cx_line[idx] + rng.normal(0, 0.5, n_river)
ry = cy_line[idx] + offsets

# Triangulate using matplotlib (Delaunay)
tri_river = mtri.Triangulation(rx, ry)
faces_river = tri_river.triangles

# Compute face-centroid bed elevation: deeper in the centre, higher at banks
fcx = rx[faces_river].mean(axis=1)
fcy = ry[faces_river].mean(axis=1)
dist_from_centre = np.abs(fcy - 4 * np.sin(fcx * 4 * np.pi / 50))
bed_elevation = -2.0 + 1.5 * (dist_from_centre / dist_from_centre.max()) + rng.normal(0, 0.1, len(fcx))

mg_river = MeshGlyph(rx, ry, faces_river)
print(f"River mesh: {mg_river.n_nodes} nodes, {mg_river.n_faces} faces")

rng = np.random.default_rng(123) # River centreline: sinusoidal path s = np.linspace(0, 4 * np.pi, 200) cx_line = s * 50 / (4 * np.pi) # 0..50 m along x cy_line = 4 * np.sin(s) # sinusoidal meander # Scatter points around the centreline with a channel half-width of ~3 m n_river = 600 idx = rng.integers(0, len(s), n_river) offsets = rng.normal(0, 1.2, n_river) rx = cx_line[idx] + rng.normal(0, 0.5, n_river) ry = cy_line[idx] + offsets # Triangulate using matplotlib (Delaunay) tri_river = mtri.Triangulation(rx, ry) faces_river = tri_river.triangles # Compute face-centroid bed elevation: deeper in the centre, higher at banks fcx = rx[faces_river].mean(axis=1) fcy = ry[faces_river].mean(axis=1) dist_from_centre = np.abs(fcy - 4 * np.sin(fcx * 4 * np.pi / 50)) bed_elevation = -2.0 + 1.5 * (dist_from_centre / dist_from_centre.max()) + rng.normal(0, 0.1, len(fcx)) mg_river = MeshGlyph(rx, ry, faces_river) print(f"River mesh: {mg_river.n_nodes} nodes, {mg_river.n_faces} faces")

In [ ]:

fig, axes = plt.subplots(2, 1, figsize=(14, 10))

# Top: bed elevation with face colouring
mg_river.plot(
    bed_elevation,
    location="face",
    ax=axes[0],
    cmap="YlOrBr_r",
    vmin=-2.5,
    vmax=0.0,
    title="River Channel — Bed Elevation (face-centered)",
)
axes[0].set_xlabel("Along-channel distance [m]")
axes[0].set_ylabel("Cross-channel [m]")

# Bottom: node-centered velocity magnitude (synthetic)
node_vel = np.exp(-0.3 * np.abs(ry - 4 * np.sin(rx * 4 * np.pi / 50)))
mg_river.plot(
    node_vel,
    location="node",
    ax=axes[1],
    cmap="plasma",
    levels=20,
    title="River Channel — Velocity Magnitude (node-centered contour)",
)
mg_river.plot_outline(ax=axes[1], color="white", linewidth=0.15)
axes[1].set_xlabel("Along-channel distance [m]")
axes[1].set_ylabel("Cross-channel [m]")

plt.tight_layout()
plt.show()

fig, axes = plt.subplots(2, 1, figsize=(14, 10)) # Top: bed elevation with face colouring mg_river.plot( bed_elevation, location="face", ax=axes[0], cmap="YlOrBr_r", vmin=-2.5, vmax=0.0, title="River Channel — Bed Elevation (face-centered)", ) axes[0].set_xlabel("Along-channel distance [m]") axes[0].set_ylabel("Cross-channel [m]") # Bottom: node-centered velocity magnitude (synthetic) node_vel = np.exp(-0.3 * np.abs(ry - 4 * np.sin(rx * 4 * np.pi / 50))) mg_river.plot( node_vel, location="node", ax=axes[1], cmap="plasma", levels=20, title="River Channel — Velocity Magnitude (node-centered contour)", ) mg_river.plot_outline(ax=axes[1], color="white", linewidth=0.15) axes[1].set_xlabel("Along-channel distance [m]") axes[1].set_ylabel("Cross-channel [m]") plt.tight_layout() plt.show()

7. Explicit Edge Connectivity — Wireframe with Pre-built Edges¶

When edge-node connectivity is available (common in UGRID NetCDF files), pass it to MeshGlyph for faster vectorized wireframe rendering.

In [ ]:

# Build a simple structured quad grid with explicit edges
nx, ny = 6, 4
gx, gy = np.meshgrid(np.arange(nx, dtype=float), np.arange(ny, dtype=float))
gx = gx.ravel()
gy = gy.ravel()

# Build quad faces: each cell (i,j) -> [i*nx+j, i*nx+j+1, (i+1)*nx+j+1, (i+1)*nx+j]
quads = []
for i in range(ny - 1):
    for j in range(nx - 1):
        n0 = i * nx + j
        quads.append([n0, n0 + 1, n0 + nx + 1, n0 + nx])
quad_faces = np.array(quads)

# Build explicit edges (horizontal + vertical)
edges = []
for i in range(ny):
    for j in range(nx - 1):
        edges.append([i * nx + j, i * nx + j + 1])
for i in range(ny - 1):
    for j in range(nx):
        edges.append([i * nx + j, (i + 1) * nx + j])
edge_arr = np.array(edges)

mg_quad = MeshGlyph(gx, gy, quad_faces, edge_node_connectivity=edge_arr)
print(f"Quad grid: {mg_quad.n_nodes} nodes, {mg_quad.n_faces} faces, {mg_quad.n_edges} edges")

# Build a simple structured quad grid with explicit edges nx, ny = 6, 4 gx, gy = np.meshgrid(np.arange(nx, dtype=float), np.arange(ny, dtype=float)) gx = gx.ravel() gy = gy.ravel() # Build quad faces: each cell (i,j) -> [i*nx+j, i*nx+j+1, (i+1)*nx+j+1, (i+1)*nx+j] quads = [] for i in range(ny - 1): for j in range(nx - 1): n0 = i * nx + j quads.append([n0, n0 + 1, n0 + nx + 1, n0 + nx]) quad_faces = np.array(quads) # Build explicit edges (horizontal + vertical) edges = [] for i in range(ny): for j in range(nx - 1): edges.append([i * nx + j, i * nx + j + 1]) for i in range(ny - 1): for j in range(nx): edges.append([i * nx + j, (i + 1) * nx + j]) edge_arr = np.array(edges) mg_quad = MeshGlyph(gx, gy, quad_faces, edge_node_connectivity=edge_arr) print(f"Quad grid: {mg_quad.n_nodes} nodes, {mg_quad.n_faces} faces, {mg_quad.n_edges} edges")

In [ ]:

quad_data = np.sin(gx[quad_faces].mean(axis=1)) * np.cos(gy[quad_faces].mean(axis=1))

fig, axes = plt.subplots(1, 2, figsize=(14, 5))

mg_quad.plot(quad_data, location="face", ax=axes[0], cmap="PiYG",
             edgecolor="black", title="Quad Grid — Face Data")
axes[0].set_xlabel("x")
axes[0].set_ylabel("y")

mg_quad.plot_outline(ax=axes[1], color="navy", linewidth=1.0)
axes[1].set_title("Quad Grid — Wireframe (explicit edges)")
axes[1].set_xlabel("x")
axes[1].set_ylabel("y")

plt.tight_layout()
plt.show()

quad_data = np.sin(gx[quad_faces].mean(axis=1)) * np.cos(gy[quad_faces].mean(axis=1)) fig, axes = plt.subplots(1, 2, figsize=(14, 5)) mg_quad.plot(quad_data, location="face", ax=axes[0], cmap="PiYG", edgecolor="black", title="Quad Grid — Face Data") axes[0].set_xlabel("x") axes[0].set_ylabel("y") mg_quad.plot_outline(ax=axes[1], color="navy", linewidth=1.0) axes[1].set_title("Quad Grid — Wireframe (explicit edges)") axes[1].set_xlabel("x") axes[1].set_ylabel("y") plt.tight_layout() plt.show()

8. Side-by-Side Comparison — Face vs Node¶

Compare the two rendering modes on the same mesh and data to see the difference between flat face colouring and smooth node interpolation.

In [ ]:

fig, axes = plt.subplots(1, 2, figsize=(14, 6))

# Face-centered (flat colours per triangle)
mg.plot(face_depth, location="face", ax=axes[0], cmap="RdYlBu_r",
        edgecolor="lightgrey", title="Face-Centered (tripcolor)")
axes[0].set_xlabel("x [m]")
axes[0].set_ylabel("y [m]")

# Node-centered (smooth contour interpolation)
mg.plot(node_elevation, location="node", ax=axes[1], cmap="RdYlBu_r",
        levels=15, title="Node-Centered (tricontourf)")
axes[1].set_xlabel("x [m]")
axes[1].set_ylabel("y [m]")

plt.tight_layout()
plt.show()

fig, axes = plt.subplots(1, 2, figsize=(14, 6)) # Face-centered (flat colours per triangle) mg.plot(face_depth, location="face", ax=axes[0], cmap="RdYlBu_r", edgecolor="lightgrey", title="Face-Centered (tripcolor)") axes[0].set_xlabel("x [m]") axes[0].set_ylabel("y [m]") # Node-centered (smooth contour interpolation) mg.plot(node_elevation, location="node", ax=axes[1], cmap="RdYlBu_r", levels=15, title="Node-Centered (tricontourf)") axes[1].set_xlabel("x [m]") axes[1].set_ylabel("y [m]") plt.tight_layout() plt.show()