Move medial tire drawing script into lib

papers/medial_tire_decompositions_of_plane_triangulations/experiments/draw_random_medial_tire_decompositions.py

@@ -1,622 +1,15 @@
-"""Draw medial tire decompositions of random 5-connected triangulations.
+"""Compatibility wrapper for the medial tire decomposition drawing script."""
-
-The source graphs come from ``plantri -c5`` in graph6 format.  For each sampled
-30-vertex triangulation, this script chooses a random source vertex, builds the
-BFS depth-component tire tree, and draws both the tire tree and the medial
-tread model for each depth component.
-"""
 
 from __future__ import annotations
 
-import argparse
-import math
 import os
-import random
-import subprocess
 import sys
-from collections import defaultdict
-from dataclasses import dataclass
-from pathlib import Path
 
-PAPER_DIR = Path(__file__).resolve().parents[1]
+PAPER_DIR = os.path.dirname(os.path.dirname(os.path.abspath(__file__)))
-REPO_ROOT = PAPER_DIR.parents[1]
+if PAPER_DIR not in sys.path:
-os.environ.setdefault(
+    sys.path.insert(0, PAPER_DIR)
-    "MPLCONFIGDIR", str(PAPER_DIR / "experiments" / ".matplotlib-cache")
-)
-os.environ.setdefault("XDG_CACHE_HOME", str(PAPER_DIR / "experiments" / ".cache"))
 
-import matplotlib
+from lib.draw_random_medial_tire_decompositions import main
-
-matplotlib.use("Agg")
-import matplotlib.pyplot as plt
-from matplotlib.lines import Line2D
-import networkx as nx
-
-if str(PAPER_DIR) not in sys.path:
-    sys.path.insert(0, str(PAPER_DIR))
-
-from lib.medial_tire_decomposition import (
-    annular_cycle_components,
-    ekey,
-    medial_tire_facemodel,
-)
-
-
-@dataclass(frozen=True)
-class TreadNode:
-    idx: int
-    depth: int
-    face_indices: tuple[int, ...]
-    annular: frozenset
-    up: frozenset
-    down: frozenset
-    bites: frozenset
-    medial: nx.Graph
-    annular_cycles: tuple[tuple, ...]
-
-
-@dataclass(frozen=True)
-class Augmentation:
-    graph: nx.Graph
-    added_vertices: tuple[int, ...]
-    filled_faces: tuple[tuple[int, tuple[int, int, int], int], ...]
-
-
-def sample_plantri_graphs(n: int, count: int, seed: int, scan_limit: int) -> list[nx.Graph]:
-    cmd = [str(REPO_ROOT / "plantri" / "plantri"), "-g", "-c5", str(n)]
-    rng = random.Random(seed)
-    sample: list[tuple[int, nx.Graph]] = []
-    seen = 0
-    with subprocess.Popen(cmd, stdout=subprocess.PIPE, stderr=subprocess.PIPE) as proc:
-        assert proc.stdout is not None
-        for raw in proc.stdout:
-            line = raw.strip()
-            if not line or line.startswith(b">>"):
-                continue
-            graph = nx.from_graph6_bytes(line)
-            if nx.node_connectivity(graph) < 5:
-                continue
-            seen += 1
-            if len(sample) < count:
-                sample.append((seen, graph))
-            else:
-                j = rng.randrange(seen)
-                if j < count:
-                    sample[j] = (seen, graph)
-            if seen >= scan_limit:
-                proc.terminate()
-                break
-        proc.wait(timeout=10)
-    if len(sample) < count:
-        raise RuntimeError(f"only found {len(sample)} graphs after scanning {seen}")
-    return [graph for _ordinal, graph in sample]
-
-
-def triangular_faces(g: nx.Graph):
-    ok, emb = nx.check_planarity(g)
-    if not ok:
-        raise ValueError("not planar")
-    seen = set()
-    faces = []
-    for u, v in list(emb.edges()):
-        if (u, v) in seen:
-            continue
-        face = tuple(emb.traverse_face(u, v, mark_half_edges=seen))
-        faces.append(face)
-    return faces
-
-
-def edge_face_data(faces):
-    face_edges = []
-    edge_faces: dict[tuple, list[int]] = defaultdict(list)
-    for i, face in enumerate(faces):
-        edges = {
-            ekey(face[0], face[1]),
-            ekey(face[1], face[2]),
-            ekey(face[2], face[0]),
-        }
-        face_edges.append(edges)
-        for edge in edges:
-            edge_faces[edge].append(i)
-    return face_edges, edge_faces
-
-
-def augment_same_level_faces(g: nx.Graph, source: int) -> Augmentation:
-    """Stack a new vertex into every facial triangle with one BFS level.
-
-    If a triangular face has all three vertices at level d, the new vertex is
-    adjacent to those three vertices and therefore has level d + 1.  This turns
-    the same-level region into three adjacent-level tread faces before the tire
-    decomposition is extracted.
-    """
-    levels = nx.single_source_shortest_path_length(g, source)
-    faces = triangular_faces(g)
-    augmented = g.copy()
-    next_vertex = max(augmented.nodes()) + 1
-    added = []
-    filled = []
-
-    for face in faces:
-        face_levels = {levels[v] for v in face}
-        if len(face_levels) != 1:
-            continue
-        new_vertex = next_vertex
-        next_vertex += 1
-        augmented.add_node(new_vertex)
-        augmented.add_edges_from((new_vertex, v) for v in face)
-        added.append(new_vertex)
-        filled.append((new_vertex, tuple(face), next(iter(face_levels))))
-
-    return Augmentation(
-        graph=augmented,
-        added_vertices=tuple(added),
-        filled_faces=tuple(filled),
-    )
-
-
-def depth_components(faces, face_edges, edge_faces, levels):
-    depths = [min(levels[v] for v in face) for face in faces]
-    dual_adj: dict[int, set[int]] = defaultdict(set)
-    for incident in edge_faces.values():
-        for a in range(len(incident)):
-            for b in range(a + 1, len(incident)):
-                dual_adj[incident[a]].add(incident[b])
-                dual_adj[incident[b]].add(incident[a])
-
-    comps = []
-    seen = [False] * len(faces)
-    for start in range(len(faces)):
-        if seen[start]:
-            continue
-        depth = depths[start]
-        stack = [start]
-        comp = []
-        seen[start] = True
-        while stack:
-            face = stack.pop()
-            comp.append(face)
-            for other in dual_adj[face]:
-                if not seen[other] and depths[other] == depth:
-                    seen[other] = True
-                    stack.append(other)
-        comps.append((depth, tuple(sorted(comp))))
-    return comps, depths, dual_adj
-
-
-def tread_from_component(faces, levels, face_indices):
-    tread_faces = [faces[i] for i in face_indices]
-    if not tread_faces:
-        return None
-    depth = min(min(levels[v] for v in face) for face in tread_faces)
-    annular, up, down = set(), set(), set()
-    face_of_down = defaultdict(int)
-    for face in tread_faces:
-        for x, y in ((face[0], face[1]), (face[1], face[2]), (face[2], face[0])):
-            e = ekey(x, y)
-            lx, ly = levels[x], levels[y]
-            if {lx, ly} == {depth, depth + 1}:
-                annular.add(e)
-            elif lx == ly == depth:
-                up.add(e)
-            elif lx == ly == depth + 1:
-                down.add(e)
-                face_of_down[e] += 1
-    if len(annular) < 3:
-        return None
-    return {
-        "tread_faces": tread_faces,
-        "annular": annular,
-        "up": up,
-        "down": down,
-        "bites": {e for e in down if face_of_down[e] == 2},
-    }
-
-
-def build_tire_tree(g: nx.Graph, source: int, augment: bool = True):
-    augmentation = augment_same_level_faces(g, source) if augment else Augmentation(g, (), ())
-    work_graph = augmentation.graph
-    faces = triangular_faces(work_graph)
-    face_edges, edge_faces = edge_face_data(faces)
-    levels = nx.single_source_shortest_path_length(work_graph, source)
-    comps, depths, dual_adj = depth_components(faces, face_edges, edge_faces, levels)
-    comp_of_face = {}
-    for comp_idx, (_depth, face_indices) in enumerate(comps):
-        for face in face_indices:
-            comp_of_face[face] = comp_idx
-
-    nodes: list[TreadNode] = []
-    comp_to_node = {}
-    for comp_idx, (depth, face_indices) in enumerate(comps):
-        tread = tread_from_component(faces, levels, face_indices)
-        if tread is None or len(tread["up"]) < 3:
-            continue
-        mt = medial_tire_facemodel(tread["tread_faces"])
-        annular_cycles = tuple(annular_cycle_components(mt, tread["annular"]))
-        if not annular_cycles:
-            continue
-        node = TreadNode(
-            idx=len(nodes),
-            depth=depth,
-            face_indices=face_indices,
-            annular=frozenset(tread["annular"]),
-            up=frozenset(tread["up"]),
-            down=frozenset(tread["down"]),
-            bites=frozenset(tread["bites"]),
-            medial=mt,
-            annular_cycles=annular_cycles,
-        )
-        comp_to_node[comp_idx] = node.idx
-        nodes.append(node)
-
-    tree_edges = set()
-    for comp_idx, (depth, face_indices) in enumerate(comps):
-        if comp_idx not in comp_to_node:
-            continue
-        child = comp_to_node[comp_idx]
-        parent_candidates = set()
-        for face in face_indices:
-            for other in dual_adj[face]:
-                other_comp = comp_of_face[other]
-                if depths[other] == depth - 1 and other_comp in comp_to_node:
-                    parent_candidates.add(comp_to_node[other_comp])
-        for parent in parent_candidates:
-            tree_edges.add((parent, child))
-    return augmentation, faces, levels, nodes, sorted(tree_edges)
-
-
-def graph_layout(g: nx.Graph):
-    try:
-        return nx.planar_layout(g)
-    except nx.NetworkXException:
-        return nx.spring_layout(g, seed=0)
-
-
-def draw_base_graph(ax, g, levels, source, added_vertices=()):
-    pos = graph_layout(g)
-    max_level = max(levels.values())
-    cmap = plt.get_cmap("viridis", max_level + 1)
-    node_colors = [cmap(levels[v]) for v in g.nodes()]
-    nx.draw_networkx_edges(g, pos, ax=ax, edge_color="#cbd5e1", width=0.8)
-    added_set = set(added_vertices)
-    nx.draw_networkx_nodes(
-        g,
-        pos,
-        ax=ax,
-        node_color=node_colors,
-        node_size=[
-            150 if v == source else 96 if v in added_set else 72
-            for v in g.nodes()
-        ],
-        edgecolors=[
-            "#dc2626" if v == source else "#7c3aed" if v in added_set else "#111827"
-            for v in g.nodes()
-        ],
-        linewidths=[
-            1.8 if v == source else 1.2 if v in added_set else 0.45
-            for v in g.nodes()
-        ],
-    )
-    labels = {v: str(v) for v in g.nodes()}
-    nx.draw_networkx_labels(g, pos, labels=labels, ax=ax, font_size=5)
-    ax.set_title(
-        f"Augmented G, source {source}; vertex levels 0..{max_level}",
-        fontsize=10,
-    )
-    ax.set_aspect("equal")
-    ax.axis("off")
-
-
-def tree_positions(nodes: list[TreadNode], tree_edges):
-    children: dict[int, list[int]] = defaultdict(list)
-    has_parent = set()
-    for parent, child in tree_edges:
-        children[parent].append(child)
-        has_parent.add(child)
-    roots = [node.idx for node in nodes if node.idx not in has_parent]
-    for child_list in children.values():
-        child_list.sort(key=lambda idx: (nodes[idx].depth, idx))
-
-    x_counter = 0
-    pos = {}
-
-    def place(idx, depth):
-        nonlocal x_counter
-        if not children[idx]:
-            pos[idx] = (x_counter, -depth)
-            x_counter += 1
-            return pos[idx][0]
-        xs = [place(child, depth + 1) for child in children[idx]]
-        x = sum(xs) / len(xs)
-        pos[idx] = (x, -depth)
-        return x
-
-    for root in sorted(roots, key=lambda idx: (nodes[idx].depth, idx)):
-        place(root, 0)
-        x_counter += 1
-    return pos
-
-
-def draw_tire_tree(ax, nodes: list[TreadNode], tree_edges):
-    pos = tree_positions(nodes, tree_edges)
-    for parent, child in tree_edges:
-        x0, y0 = pos[parent]
-        x1, y1 = pos[child]
-        ax.plot([x0, x1], [y0, y1], color="#374151", lw=1.0, zorder=1)
-    for node in nodes:
-        x, y = pos[node.idx]
-        ax.text(
-            x,
-            y,
-            f"T{node.idx}\nd={node.depth}\n{len(node.annular_cycles)} cycle(s)",
-            ha="center",
-            va="center",
-            fontsize=8,
-            bbox={
-                "boxstyle": "round,pad=0.32",
-                "facecolor": "#fef3c7",
-                "edgecolor": "#111827",
-                "linewidth": 0.9,
-            },
-            zorder=3,
-        )
-    ax.set_title("Depth-component tire tree", fontsize=10)
-    if pos:
-        xs = [p[0] for p in pos.values()]
-        ys = [p[1] for p in pos.values()]
-        ax.set_xlim(min(xs) - 1.0, max(xs) + 1.0)
-        ax.set_ylim(min(ys) - 0.7, max(ys) + 0.7)
-    ax.axis("off")
-
-
-def vertex_xy(k: int, n: int, radius: float) -> tuple[float, float]:
-    angle = math.pi / 2 - 2 * math.pi * k / n
-    return radius * math.cos(angle), radius * math.sin(angle)
-
-
-def edge_midpoint_angle(i: int, n: int) -> float:
-    return math.pi / 2 - 2 * math.pi * (i + 0.5) / n
-
-
-def draw_tread_model(ax, node: TreadNode):
-    cycle_count = len(node.annular_cycles)
-    offsets = [3.25 * (i - (cycle_count - 1) / 2) for i in range(cycle_count)]
-    apex_positions: dict[tuple, tuple[float, float]] = {}
-    apex_corners: dict[tuple, list[tuple[float, float]]] = defaultdict(list)
-    ann_positions: dict[tuple, tuple[float, float]] = {}
-
-    for cycle_idx, order in enumerate(node.annular_cycles):
-        n = len(order)
-        dx = offsets[cycle_idx]
-        ann = {
-            vertex: (dx + x, y)
-            for vertex, (x, y) in zip(order, [vertex_xy(k, n, 1.0) for k in range(n)])
-        }
-        ann_positions.update(ann)
-
-        cyc_x = [ann[v][0] for v in order] + [ann[order[0]][0]]
-        cyc_y = [ann[v][1] for v in order] + [ann[order[0]][1]]
-        ax.plot(cyc_x, cyc_y, color="black", lw=1.3, zorder=2)
-
-        for i, a in enumerate(order):
-            b = order[(i + 1) % n]
-            apexes = [
-                w for w in set(node.medial.neighbors(a)) & set(node.medial.neighbors(b))
-                if w not in node.annular
-            ]
-            for apex in apexes:
-                apex_corners[apex].extend([ann[a], ann[b]])
-                if apex in apex_positions:
-                    continue
-                angle = edge_midpoint_angle(i, n)
-                if apex in node.up:
-                    radius = 1.42
-                else:
-                    radius = 0.58
-                apex_positions[apex] = (
-                    dx + radius * math.cos(angle),
-                    radius * math.sin(angle),
-                )
-
-    for apex, corners in apex_corners.items():
-        if apex in node.bites and corners:
-            cx = sum(p[0] for p in corners) / len(corners)
-            cy = sum(p[1] for p in corners) / len(corners)
-            center_x = sum(offsets) / len(offsets) if offsets else 0.0
-            apex_positions[apex] = (
-                center_x + 0.82 * (cx - center_x),
-                0.82 * cy,
-            )
-        pos = apex_positions[apex]
-        for corner in corners:
-            ax.plot([pos[0], corner[0]], [pos[1], corner[1]], color="#9ca3af", lw=0.5)
-
-    for apex, pos in apex_positions.items():
-        if apex in node.up:
-            color, size, edgecolor = "#2563eb", 13, "none"
-        elif apex in node.bites:
-            color, size, edgecolor = "#7f1d1d", 24, "black"
-        else:
-            color, size, edgecolor = "#dc2626", 13, "none"
-        ax.scatter(
-            [pos[0]],
-            [pos[1]],
-            s=size,
-            color=color,
-            edgecolors=edgecolor,
-            linewidths=0.4,
-            zorder=3,
-        )
-
-    if ann_positions:
-        ax.scatter(
-            [p[0] for p in ann_positions.values()],
-            [p[1] for p in ann_positions.values()],
-            s=9,
-            color="black",
-            zorder=4,
-        )
-
-    singleton_down = set(node.down) - set(node.bites)
-    ax.set_title(
-        f"T{node.idx} d={node.depth}: {len(node.annular_cycles)} annular cycle(s)\n"
-        f"ann={len(node.annular)} up={len(node.up)} down={len(singleton_down)} "
-        f"bite={len(node.bites)}",
-        fontsize=6.4,
-        pad=1.5,
-    )
-    pad = 1.7
-    ax.set_xlim(min(offsets, default=0.0) - pad, max(offsets, default=0.0) + pad)
-    ax.set_ylim(-1.65, 1.65)
-    ax.set_aspect("equal")
-    ax.axis("off")
-
-
-def draw_medial_tire_grid(fig, outer_spec, nodes):
-    if not nodes:
-        ax = fig.add_subplot(outer_spec)
-        ax.text(0.5, 0.5, "No medial treads extracted", ha="center")
-        ax.axis("off")
-        return
-    cols = min(3, max(1, math.ceil(math.sqrt(len(nodes)))))
-    rows = math.ceil(len(nodes) / cols)
-    sub = outer_spec.subgridspec(rows, cols, wspace=0.08, hspace=0.35)
-    for i in range(rows * cols):
-        ax = fig.add_subplot(sub[i // cols, i % cols])
-        if i < len(nodes):
-            draw_tread_model(ax, nodes[i])
-        else:
-            ax.axis("off")
-
-
-def write_index(
-    path: Path,
-    graph_idx: int,
-    source: int,
-    original_graph: nx.Graph,
-    augmentation: Augmentation,
-    nodes,
-    tree_edges,
-):
-    g = augmentation.graph
-    lines = [
-        f"# Random medial tire decomposition {graph_idx}",
-        "",
-        f"- original vertices: {original_graph.number_of_nodes()}",
-        f"- original edges: {original_graph.number_of_edges()}",
-        f"- original node connectivity: {nx.node_connectivity(original_graph)}",
-        f"- augmented vertices: {g.number_of_nodes()}",
-        f"- augmented edges: {g.number_of_edges()}",
-        f"- same-level faces filled: {len(augmentation.added_vertices)}",
-        f"- source vertex: {source}",
-        f"- tire-tree nodes: {len(nodes)}",
-        f"- tire-tree edges: {len(tree_edges)}",
-        "",
-        "| node | depth | faces | annular cycles | annular | up | singleton down | bite apexes |",
-        "|--:|--:|--:|--:|--:|--:|--:|--:|",
-    ]
-    for node in nodes:
-        singleton_down = set(node.down) - set(node.bites)
-        lines.append(
-            f"| T{node.idx} | {node.depth} | {len(node.face_indices)} | "
-            f"{len(node.annular_cycles)} | {len(node.annular)} | {len(node.up)} | "
-            f"{len(singleton_down)} | {len(node.bites)} |"
-        )
-    path.write_text("\n".join(lines) + "\n")
-
-
-def draw_case(out_dir: Path, graph_idx: int, g: nx.Graph, source: int, augment: bool = True):
-    augmentation, _faces, levels, nodes, tree_edges = build_tire_tree(g, source, augment=augment)
-    work_graph = augmentation.graph
-    fig = plt.figure(figsize=(17, 10))
-    spec = fig.add_gridspec(2, 2, width_ratios=[1.15, 1.0], height_ratios=[1.0, 1.25])
-    ax_graph = fig.add_subplot(spec[0, 0])
-    ax_tree = fig.add_subplot(spec[1, 0])
-    draw_base_graph(ax_graph, work_graph, levels, source, augmentation.added_vertices)
-    draw_tire_tree(ax_tree, nodes, tree_edges)
-    draw_medial_tire_grid(fig, spec[:, 1], nodes)
-    fig.suptitle(
-        f"Random 5-connected maximal planar graph {graph_idx}: "
-        f"n={g.number_of_nodes()} (+{len(augmentation.added_vertices)}), "
-        f"source={source}",
-        fontsize=13,
-    )
-    legend = [
-        Line2D([0], [0], marker="o", color="w", label="source",
-               markerfacecolor="#fde68a", markeredgecolor="#dc2626", markersize=8),
-        Line2D([0], [0], marker="o", color="w", label="inserted vertex",
-               markerfacecolor="#fde68a", markeredgecolor="#7c3aed", markersize=8),
-        Line2D([0], [0], color="black", lw=1.3, label="annular cycle A(T)"),
-        Line2D([0], [0], marker="o", color="w", label="up tooth",
-               markerfacecolor="#2563eb", markersize=6),
-        Line2D([0], [0], marker="o", color="w", label="down tooth",
-               markerfacecolor="#dc2626", markersize=6),
-        Line2D([0], [0], marker="o", color="w", label="bite apex",
-               markerfacecolor="#7f1d1d", markeredgecolor="black", markersize=6),
-    ]
-    fig.legend(handles=legend, loc="lower center", ncol=5, fontsize=9)
-    fig.subplots_adjust(left=0.03, right=0.99, top=0.92, bottom=0.08, wspace=0.08, hspace=0.16)
-
-    png = out_dir / f"random_c5_n30_medial_tire_decomposition_{graph_idx}.png"
-    pdf = out_dir / f"random_c5_n30_medial_tire_decomposition_{graph_idx}.pdf"
-    fig.savefig(png, dpi=180)
-    fig.savefig(pdf)
-    plt.close(fig)
-    write_index(
-        out_dir / f"random_c5_n30_medial_tire_decomposition_{graph_idx}.md",
-        graph_idx,
-        source,
-        g,
-        augmentation,
-        nodes,
-        tree_edges,
-    )
-    return png, pdf, len(nodes), sum(len(node.annular_cycles) for node in nodes)
-
-
-def run(args: argparse.Namespace):
-    out_dir = Path(args.out_dir)
-    out_dir.mkdir(parents=True, exist_ok=True)
-
-    if args.graph6:
-        graphs = [nx.from_graph6_bytes(args.graph6.encode())]
-        if args.source is None:
-            raise ValueError("--source is required with --graph6")
-        sources = [args.source]
-    else:
-        graphs = sample_plantri_graphs(args.n, args.count, args.seed, args.scan_limit)
-        rng = random.Random(args.seed + 101)
-        sources = [rng.choice(list(graph.nodes())) for graph in graphs]
-
-    for i, (graph, source) in enumerate(zip(graphs, sources), start=1):
-        png, pdf, node_count, annular_cycle_count = draw_case(
-            out_dir, i, graph, source, augment=not args.no_augment_same_level_faces
-        )
-        print(
-            f"case {i}: source={source}, connectivity={nx.node_connectivity(graph)}, "
-            f"tire nodes={node_count}, annular cycles={annular_cycle_count}"
-        )
-        print(f"  wrote {png}")
-        print(f"  wrote {pdf}")
-
-
-def main():
-    parser = argparse.ArgumentParser(description=__doc__)
-    parser.add_argument("--n", type=int, default=30)
-    parser.add_argument("--count", type=int, default=2)
-    parser.add_argument("--seed", type=int, default=20260615)
-    parser.add_argument("--scan-limit", type=int, default=500)
-    parser.add_argument("--graph6", help="draw this graph6 graph instead of sampling")
-    parser.add_argument("--source", type=int, help="source vertex for --graph6")
-    parser.add_argument(
-        "--no-augment-same-level-faces",
-        action="store_true",
-        help="skip the same-level-face vertex insertion step",
-    )
-    parser.add_argument(
-        "--out-dir",
-        default=str(PAPER_DIR / "experiments" / "random_medial_tire_decompositions"),
-    )
-    run(parser.parse_args())
 
 
 if __name__ == "__main__":

papers/medial_tire_decompositions_of_plane_triangulations/lib/draw_random_medial_tire_decompositions.py

@@ -0,0 +1,706 @@
+"""Draw medial tire decompositions of random 5-connected triangulations.
+
+The source graphs come from ``plantri -c5`` in graph6 format.  For each sampled
+30-vertex triangulation, this script chooses a random source vertex, builds the
+BFS depth-component tire tree, and draws both the tire tree and the medial
+tread model for each depth component.
+"""
+
+from __future__ import annotations
+
+import argparse
+import math
+import os
+import random
+import subprocess
+import sys
+from collections import defaultdict
+from dataclasses import dataclass
+from pathlib import Path
+
+PAPER_DIR = Path(__file__).resolve().parents[1]
+REPO_ROOT = PAPER_DIR.parents[1]
+os.environ.setdefault(
+    "MPLCONFIGDIR", str(PAPER_DIR / "experiments" / ".matplotlib-cache")
+)
+os.environ.setdefault("XDG_CACHE_HOME", str(PAPER_DIR / "experiments" / ".cache"))
+
+import matplotlib
+
+matplotlib.use("Agg")
+import matplotlib.pyplot as plt
+from matplotlib.lines import Line2D
+import networkx as nx
+
+if str(PAPER_DIR) not in sys.path:
+    sys.path.insert(0, str(PAPER_DIR))
+
+from lib.medial_tire_decomposition import (
+    annular_cycle_components,
+    ekey,
+    medial_tire_facemodel,
+)
+
+
+@dataclass(frozen=True)
+class TreadNode:
+    idx: int
+    depth: int
+    face_indices: tuple[int, ...]
+    annular: frozenset
+    up: frozenset
+    down: frozenset
+    bites: frozenset
+    medial: nx.Graph
+    annular_cycles: tuple[tuple, ...]
+
+
+@dataclass(frozen=True)
+class Augmentation:
+    graph: nx.Graph
+    added_vertices: tuple[int, ...]
+    filled_faces: tuple[tuple[int, tuple[int, int, int], int], ...]
+
+
+def sample_plantri_graphs(n: int, count: int, seed: int, scan_limit: int) -> list[nx.Graph]:
+    cmd = [str(REPO_ROOT / "plantri" / "plantri"), "-g", "-c5", str(n)]
+    rng = random.Random(seed)
+    sample: list[tuple[int, nx.Graph]] = []
+    seen = 0
+    with subprocess.Popen(cmd, stdout=subprocess.PIPE, stderr=subprocess.PIPE) as proc:
+        assert proc.stdout is not None
+        for raw in proc.stdout:
+            line = raw.strip()
+            if not line or line.startswith(b">>"):
+                continue
+            graph = nx.from_graph6_bytes(line)
+            if nx.node_connectivity(graph) < 5:
+                continue
+            seen += 1
+            if len(sample) < count:
+                sample.append((seen, graph))
+            else:
+                j = rng.randrange(seen)
+                if j < count:
+                    sample[j] = (seen, graph)
+            if seen >= scan_limit:
+                proc.terminate()
+                break
+        proc.wait(timeout=10)
+    if len(sample) < count:
+        raise RuntimeError(f"only found {len(sample)} graphs after scanning {seen}")
+    return [graph for _ordinal, graph in sample]
+
+
+def triangular_faces(g: nx.Graph):
+    ok, emb = nx.check_planarity(g)
+    if not ok:
+        raise ValueError("not planar")
+    seen = set()
+    faces = []
+    for u, v in list(emb.edges()):
+        if (u, v) in seen:
+            continue
+        face = tuple(emb.traverse_face(u, v, mark_half_edges=seen))
+        faces.append(face)
+    return faces
+
+
+def edge_face_data(faces):
+    face_edges = []
+    edge_faces: dict[tuple, list[int]] = defaultdict(list)
+    for i, face in enumerate(faces):
+        edges = {
+            ekey(face[0], face[1]),
+            ekey(face[1], face[2]),
+            ekey(face[2], face[0]),
+        }
+        face_edges.append(edges)
+        for edge in edges:
+            edge_faces[edge].append(i)
+    return face_edges, edge_faces
+
+
+def augment_same_level_faces(g: nx.Graph, source: int) -> Augmentation:
+    """Stack a new vertex into every facial triangle with one BFS level.
+
+    If a triangular face has all three vertices at level d, the new vertex is
+    adjacent to those three vertices and therefore has level d + 1.  This turns
+    the same-level region into three adjacent-level tread faces before the tire
+    decomposition is extracted.
+    """
+    levels = nx.single_source_shortest_path_length(g, source)
+    faces = triangular_faces(g)
+    augmented = g.copy()
+    next_vertex = max(augmented.nodes()) + 1
+    added = []
+    filled = []
+
+    for face in faces:
+        face_levels = {levels[v] for v in face}
+        if len(face_levels) != 1:
+            continue
+        new_vertex = next_vertex
+        next_vertex += 1
+        augmented.add_node(new_vertex)
+        augmented.add_edges_from((new_vertex, v) for v in face)
+        added.append(new_vertex)
+        filled.append((new_vertex, tuple(face), next(iter(face_levels))))
+
+    return Augmentation(
+        graph=augmented,
+        added_vertices=tuple(added),
+        filled_faces=tuple(filled),
+    )
+
+
+def depth_components(faces, face_edges, edge_faces, levels):
+    depths = [min(levels[v] for v in face) for face in faces]
+    dual_adj: dict[int, set[int]] = defaultdict(set)
+    for incident in edge_faces.values():
+        for a in range(len(incident)):
+            for b in range(a + 1, len(incident)):
+                dual_adj[incident[a]].add(incident[b])
+                dual_adj[incident[b]].add(incident[a])
+
+    comps = []
+    seen = [False] * len(faces)
+    for start in range(len(faces)):
+        if seen[start]:
+            continue
+        depth = depths[start]
+        stack = [start]
+        comp = []
+        seen[start] = True
+        while stack:
+            face = stack.pop()
+            comp.append(face)
+            for other in dual_adj[face]:
+                if not seen[other] and depths[other] == depth:
+                    seen[other] = True
+                    stack.append(other)
+        comps.append((depth, tuple(sorted(comp))))
+    return comps, depths, dual_adj
+
+
+def tread_from_component(faces, levels, face_indices):
+    tread_faces = [faces[i] for i in face_indices]
+    if not tread_faces:
+        return None
+    depth = min(min(levels[v] for v in face) for face in tread_faces)
+    annular, up, down = set(), set(), set()
+    face_of_down = defaultdict(int)
+    for face in tread_faces:
+        for x, y in ((face[0], face[1]), (face[1], face[2]), (face[2], face[0])):
+            e = ekey(x, y)
+            lx, ly = levels[x], levels[y]
+            if {lx, ly} == {depth, depth + 1}:
+                annular.add(e)
+            elif lx == ly == depth:
+                up.add(e)
+            elif lx == ly == depth + 1:
+                down.add(e)
+                face_of_down[e] += 1
+    if len(annular) < 3:
+        return None
+    return {
+        "tread_faces": tread_faces,
+        "annular": annular,
+        "up": up,
+        "down": down,
+        "bites": {e for e in down if face_of_down[e] == 2},
+    }
+
+
+def build_tire_tree(g: nx.Graph, source: int, augment: bool = True):
+    augmentation = augment_same_level_faces(g, source) if augment else Augmentation(g, (), ())
+    work_graph = augmentation.graph
+    faces = triangular_faces(work_graph)
+    face_edges, edge_faces = edge_face_data(faces)
+    levels = nx.single_source_shortest_path_length(work_graph, source)
+    comps, depths, dual_adj = depth_components(faces, face_edges, edge_faces, levels)
+    comp_of_face = {}
+    for comp_idx, (_depth, face_indices) in enumerate(comps):
+        for face in face_indices:
+            comp_of_face[face] = comp_idx
+
+    nodes: list[TreadNode] = []
+    comp_to_node = {}
+    for comp_idx, (depth, face_indices) in enumerate(comps):
+        tread = tread_from_component(faces, levels, face_indices)
+        if tread is None or len(tread["up"]) < 3:
+            continue
+        mt = medial_tire_facemodel(tread["tread_faces"])
+        annular_cycles = tuple(annular_cycle_components(mt, tread["annular"]))
+        if not annular_cycles:
+            continue
+        node = TreadNode(
+            idx=len(nodes),
+            depth=depth,
+            face_indices=face_indices,
+            annular=frozenset(tread["annular"]),
+            up=frozenset(tread["up"]),
+            down=frozenset(tread["down"]),
+            bites=frozenset(tread["bites"]),
+            medial=mt,
+            annular_cycles=annular_cycles,
+        )
+        comp_to_node[comp_idx] = node.idx
+        nodes.append(node)
+
+    tree_edges = set()
+    for comp_idx, (depth, face_indices) in enumerate(comps):
+        if comp_idx not in comp_to_node:
+            continue
+        child = comp_to_node[comp_idx]
+        parent_candidates = set()
+        for face in face_indices:
+            for other in dual_adj[face]:
+                other_comp = comp_of_face[other]
+                if depths[other] == depth - 1 and other_comp in comp_to_node:
+                    parent_candidates.add(comp_to_node[other_comp])
+        for parent in parent_candidates:
+            tree_edges.add((parent, child))
+    return augmentation, faces, levels, nodes, sorted(tree_edges)
+
+
+def graph_layout(g: nx.Graph):
+    try:
+        return nx.planar_layout(g)
+    except nx.NetworkXException:
+        return nx.spring_layout(g, seed=0)
+
+
+def draw_base_graph(ax, g, levels, source, added_vertices=()):
+    pos = graph_layout(g)
+    max_level = max(levels.values())
+    cmap = plt.get_cmap("viridis", max_level + 1)
+    node_colors = [cmap(levels[v]) for v in g.nodes()]
+    nx.draw_networkx_edges(g, pos, ax=ax, edge_color="#cbd5e1", width=0.8)
+    added_set = set(added_vertices)
+    nx.draw_networkx_nodes(
+        g,
+        pos,
+        ax=ax,
+        node_color=node_colors,
+        node_size=[
+            150 if v == source else 96 if v in added_set else 72
+            for v in g.nodes()
+        ],
+        edgecolors=[
+            "#dc2626" if v == source else "#7c3aed" if v in added_set else "#111827"
+            for v in g.nodes()
+        ],
+        linewidths=[
+            1.8 if v == source else 1.2 if v in added_set else 0.45
+            for v in g.nodes()
+        ],
+    )
+    labels = {v: str(v) for v in g.nodes()}
+    nx.draw_networkx_labels(g, pos, labels=labels, ax=ax, font_size=5)
+    ax.set_title(
+        f"Augmented G, source {source}; vertex levels 0..{max_level}",
+        fontsize=10,
+    )
+    ax.set_aspect("equal")
+    ax.axis("off")
+
+
+def tree_positions(nodes: list[TreadNode], tree_edges):
+    children: dict[int, list[int]] = defaultdict(list)
+    has_parent = set()
+    for parent, child in tree_edges:
+        children[parent].append(child)
+        has_parent.add(child)
+    roots = [node.idx for node in nodes if node.idx not in has_parent]
+    for child_list in children.values():
+        child_list.sort(key=lambda idx: (nodes[idx].depth, idx))
+
+    x_counter = 0
+    pos = {}
+
+    def place(idx, depth):
+        nonlocal x_counter
+        if not children[idx]:
+            pos[idx] = (x_counter, -depth)
+            x_counter += 1
+            return pos[idx][0]
+        xs = [place(child, depth + 1) for child in children[idx]]
+        x = sum(xs) / len(xs)
+        pos[idx] = (x, -depth)
+        return x
+
+    for root in sorted(roots, key=lambda idx: (nodes[idx].depth, idx)):
+        place(root, 0)
+        x_counter += 1
+    return pos
+
+
+def draw_tire_tree(ax, nodes: list[TreadNode], tree_edges):
+    pos = tree_positions(nodes, tree_edges)
+    for parent, child in tree_edges:
+        x0, y0 = pos[parent]
+        x1, y1 = pos[child]
+        ax.plot([x0, x1], [y0, y1], color="#374151", lw=1.0, zorder=1)
+    for node in nodes:
+        x, y = pos[node.idx]
+        ax.text(
+            x,
+            y,
+            f"T{node.idx}\nd={node.depth}\n{len(node.annular_cycles)} cycle(s)",
+            ha="center",
+            va="center",
+            fontsize=8,
+            bbox={
+                "boxstyle": "round,pad=0.32",
+                "facecolor": "#fef3c7",
+                "edgecolor": "#111827",
+                "linewidth": 0.9,
+            },
+            zorder=3,
+        )
+    ax.set_title("Depth-component tire tree", fontsize=10)
+    if pos:
+        xs = [p[0] for p in pos.values()]
+        ys = [p[1] for p in pos.values()]
+        ax.set_xlim(min(xs) - 1.0, max(xs) + 1.0)
+        ax.set_ylim(min(ys) - 0.7, max(ys) + 0.7)
+    ax.axis("off")
+
+
+def vertex_xy(k: int, n: int, radius: float) -> tuple[float, float]:
+    angle = math.pi / 2 - 2 * math.pi * k / n
+    return radius * math.cos(angle), radius * math.sin(angle)
+
+
+def edge_midpoint_angle(i: int, n: int) -> float:
+    return math.pi / 2 - 2 * math.pi * (i + 0.5) / n
+
+
+def annular_cycle_edges(node: TreadNode) -> set[tuple]:
+    edges = set()
+    for order in node.annular_cycles:
+        for i, a in enumerate(order):
+            b = order[(i + 1) % len(order)]
+            edges.add(tuple(sorted((a, b))))
+    return edges
+
+
+def draw_compound_tread_model(ax, node: TreadNode):
+    """Draw a compound tread using a planar layout of its actual medial graph."""
+    try:
+        pos = nx.planar_layout(node.medial)
+    except nx.NetworkXException:
+        pos = nx.spring_layout(node.medial, seed=node.idx)
+
+    cycle_edges = annular_cycle_edges(node)
+    non_cycle_edges = [
+        edge for edge in node.medial.edges()
+        if tuple(sorted(edge)) not in cycle_edges
+    ]
+
+    nx.draw_networkx_edges(
+        node.medial,
+        pos,
+        edgelist=non_cycle_edges,
+        ax=ax,
+        edge_color="#cbd5e1",
+        width=0.7,
+    )
+    nx.draw_networkx_edges(
+        node.medial,
+        pos,
+        edgelist=list(cycle_edges),
+        ax=ax,
+        edge_color="black",
+        width=1.4,
+    )
+
+    annular = set(node.annular)
+    singleton_down = set(node.down) - set(node.bites)
+    categories = [
+        (annular, "black", 13, "none"),
+        (set(node.up) - annular, "#2563eb", 18, "none"),
+        (singleton_down - annular, "#dc2626", 18, "none"),
+        (set(node.bites) - annular, "#7f1d1d", 28, "black"),
+    ]
+    for vertices, color, size, edgecolor in categories:
+        drawn = [v for v in vertices if v in pos]
+        if not drawn:
+            continue
+        ax.scatter(
+            [pos[v][0] for v in drawn],
+            [pos[v][1] for v in drawn],
+            s=size,
+            color=color,
+            edgecolors=edgecolor,
+            linewidths=0.4,
+            zorder=3,
+        )
+
+    xs = [p[0] for p in pos.values()]
+    ys = [p[1] for p in pos.values()]
+    xpad = max(0.05, (max(xs) - min(xs)) * 0.12)
+    ypad = max(0.05, (max(ys) - min(ys)) * 0.12)
+    ax.set_xlim(min(xs) - xpad, max(xs) + xpad)
+    ax.set_ylim(min(ys) - ypad, max(ys) + ypad)
+    ax.set_aspect("equal")
+    ax.axis("off")
+
+
+def draw_tread_model(ax, node: TreadNode):
+    if len(node.annular_cycles) > 1:
+        draw_compound_tread_model(ax, node)
+        singleton_down = set(node.down) - set(node.bites)
+        ax.set_title(
+            f"T{node.idx} d={node.depth}: {len(node.annular_cycles)} annular cycle(s)\n"
+            f"ann={len(node.annular)} up={len(node.up)} down={len(singleton_down)} "
+            f"bite={len(node.bites)}",
+            fontsize=6.4,
+            pad=1.5,
+        )
+        return
+
+    cycle_count = len(node.annular_cycles)
+    offsets = [3.25 * (i - (cycle_count - 1) / 2) for i in range(cycle_count)]
+    apex_positions: dict[tuple, tuple[float, float]] = {}
+    apex_corners: dict[tuple, list[tuple[float, float]]] = defaultdict(list)
+    ann_positions: dict[tuple, tuple[float, float]] = {}
+
+    for cycle_idx, order in enumerate(node.annular_cycles):
+        n = len(order)
+        dx = offsets[cycle_idx]
+        ann = {
+            vertex: (dx + x, y)
+            for vertex, (x, y) in zip(order, [vertex_xy(k, n, 1.0) for k in range(n)])
+        }
+        ann_positions.update(ann)
+
+        cyc_x = [ann[v][0] for v in order] + [ann[order[0]][0]]
+        cyc_y = [ann[v][1] for v in order] + [ann[order[0]][1]]
+        ax.plot(cyc_x, cyc_y, color="black", lw=1.3, zorder=2)
+
+        for i, a in enumerate(order):
+            b = order[(i + 1) % n]
+            apexes = [
+                w for w in set(node.medial.neighbors(a)) & set(node.medial.neighbors(b))
+                if w not in node.annular
+            ]
+            for apex in apexes:
+                apex_corners[apex].extend([ann[a], ann[b]])
+                if apex in apex_positions:
+                    continue
+                angle = edge_midpoint_angle(i, n)
+                if apex in node.up:
+                    radius = 1.42
+                else:
+                    radius = 0.58
+                apex_positions[apex] = (
+                    dx + radius * math.cos(angle),
+                    radius * math.sin(angle),
+                )
+
+    for apex, corners in apex_corners.items():
+        if apex in node.bites and corners:
+            cx = sum(p[0] for p in corners) / len(corners)
+            cy = sum(p[1] for p in corners) / len(corners)
+            center_x = sum(offsets) / len(offsets) if offsets else 0.0
+            apex_positions[apex] = (
+                center_x + 0.82 * (cx - center_x),
+                0.82 * cy,
+            )
+        pos = apex_positions[apex]
+        for corner in corners:
+            ax.plot([pos[0], corner[0]], [pos[1], corner[1]], color="#9ca3af", lw=0.5)
+
+    for apex, pos in apex_positions.items():
+        if apex in node.up:
+            color, size, edgecolor = "#2563eb", 13, "none"
+        elif apex in node.bites:
+            color, size, edgecolor = "#7f1d1d", 24, "black"
+        else:
+            color, size, edgecolor = "#dc2626", 13, "none"
+        ax.scatter(
+            [pos[0]],
+            [pos[1]],
+            s=size,
+            color=color,
+            edgecolors=edgecolor,
+            linewidths=0.4,
+            zorder=3,
+        )
+
+    if ann_positions:
+        ax.scatter(
+            [p[0] for p in ann_positions.values()],
+            [p[1] for p in ann_positions.values()],
+            s=9,
+            color="black",
+            zorder=4,
+        )
+
+    singleton_down = set(node.down) - set(node.bites)
+    ax.set_title(
+        f"T{node.idx} d={node.depth}: {len(node.annular_cycles)} annular cycle(s)\n"
+        f"ann={len(node.annular)} up={len(node.up)} down={len(singleton_down)} "
+        f"bite={len(node.bites)}",
+        fontsize=6.4,
+        pad=1.5,
+    )
+    pad = 1.7
+    ax.set_xlim(min(offsets, default=0.0) - pad, max(offsets, default=0.0) + pad)
+    ax.set_ylim(-1.65, 1.65)
+    ax.set_aspect("equal")
+    ax.axis("off")
+
+
+def draw_medial_tire_grid(fig, outer_spec, nodes):
+    if not nodes:
+        ax = fig.add_subplot(outer_spec)
+        ax.text(0.5, 0.5, "No medial treads extracted", ha="center")
+        ax.axis("off")
+        return
+    cols = min(3, max(1, math.ceil(math.sqrt(len(nodes)))))
+    rows = math.ceil(len(nodes) / cols)
+    sub = outer_spec.subgridspec(rows, cols, wspace=0.08, hspace=0.35)
+    for i in range(rows * cols):
+        ax = fig.add_subplot(sub[i // cols, i % cols])
+        if i < len(nodes):
+            draw_tread_model(ax, nodes[i])
+        else:
+            ax.axis("off")
+
+
+def write_index(
+    path: Path,
+    graph_idx: int,
+    source: int,
+    original_graph: nx.Graph,
+    augmentation: Augmentation,
+    nodes,
+    tree_edges,
+):
+    g = augmentation.graph
+    lines = [
+        f"# Random medial tire decomposition {graph_idx}",
+        "",
+        f"- original vertices: {original_graph.number_of_nodes()}",
+        f"- original edges: {original_graph.number_of_edges()}",
+        f"- original node connectivity: {nx.node_connectivity(original_graph)}",
+        f"- augmented vertices: {g.number_of_nodes()}",
+        f"- augmented edges: {g.number_of_edges()}",
+        f"- same-level faces filled: {len(augmentation.added_vertices)}",
+        f"- source vertex: {source}",
+        f"- tire-tree nodes: {len(nodes)}",
+        f"- tire-tree edges: {len(tree_edges)}",
+        "",
+        "| node | depth | faces | annular cycles | annular | up | singleton down | bite apexes |",
+        "|--:|--:|--:|--:|--:|--:|--:|--:|",
+    ]
+    for node in nodes:
+        singleton_down = set(node.down) - set(node.bites)
+        lines.append(
+            f"| T{node.idx} | {node.depth} | {len(node.face_indices)} | "
+            f"{len(node.annular_cycles)} | {len(node.annular)} | {len(node.up)} | "
+            f"{len(singleton_down)} | {len(node.bites)} |"
+        )
+    path.write_text("\n".join(lines) + "\n")
+
+
+def draw_case(out_dir: Path, graph_idx: int, g: nx.Graph, source: int, augment: bool = True):
+    augmentation, _faces, levels, nodes, tree_edges = build_tire_tree(g, source, augment=augment)
+    work_graph = augmentation.graph
+    fig = plt.figure(figsize=(17, 10))
+    spec = fig.add_gridspec(2, 2, width_ratios=[1.15, 1.0], height_ratios=[1.0, 1.25])
+    ax_graph = fig.add_subplot(spec[0, 0])
+    ax_tree = fig.add_subplot(spec[1, 0])
+    draw_base_graph(ax_graph, work_graph, levels, source, augmentation.added_vertices)
+    draw_tire_tree(ax_tree, nodes, tree_edges)
+    draw_medial_tire_grid(fig, spec[:, 1], nodes)
+    fig.suptitle(
+        f"Random 5-connected maximal planar graph {graph_idx}: "
+        f"n={g.number_of_nodes()} (+{len(augmentation.added_vertices)}), "
+        f"source={source}",
+        fontsize=13,
+    )
+    legend = [
+        Line2D([0], [0], marker="o", color="w", label="source",
+               markerfacecolor="#fde68a", markeredgecolor="#dc2626", markersize=8),
+        Line2D([0], [0], marker="o", color="w", label="inserted vertex",
+               markerfacecolor="#fde68a", markeredgecolor="#7c3aed", markersize=8),
+        Line2D([0], [0], color="black", lw=1.3, label="annular cycle A(T)"),
+        Line2D([0], [0], marker="o", color="w", label="up tooth",
+               markerfacecolor="#2563eb", markersize=6),
+        Line2D([0], [0], marker="o", color="w", label="down tooth",
+               markerfacecolor="#dc2626", markersize=6),
+        Line2D([0], [0], marker="o", color="w", label="bite apex",
+               markerfacecolor="#7f1d1d", markeredgecolor="black", markersize=6),
+    ]
+    fig.legend(handles=legend, loc="lower center", ncol=5, fontsize=9)
+    fig.subplots_adjust(left=0.03, right=0.99, top=0.92, bottom=0.08, wspace=0.08, hspace=0.16)
+
+    png = out_dir / f"random_c5_n30_medial_tire_decomposition_{graph_idx}.png"
+    pdf = out_dir / f"random_c5_n30_medial_tire_decomposition_{graph_idx}.pdf"
+    fig.savefig(png, dpi=180)
+    fig.savefig(pdf)
+    plt.close(fig)
+    write_index(
+        out_dir / f"random_c5_n30_medial_tire_decomposition_{graph_idx}.md",
+        graph_idx,
+        source,
+        g,
+        augmentation,
+        nodes,
+        tree_edges,
+    )
+    return png, pdf, len(nodes), sum(len(node.annular_cycles) for node in nodes)
+
+
+def run(args: argparse.Namespace):
+    out_dir = Path(args.out_dir)
+    out_dir.mkdir(parents=True, exist_ok=True)
+
+    if args.graph6:
+        graphs = [nx.from_graph6_bytes(args.graph6.encode())]
+        if args.source is None:
+            raise ValueError("--source is required with --graph6")
+        sources = [args.source]
+    else:
+        graphs = sample_plantri_graphs(args.n, args.count, args.seed, args.scan_limit)
+        rng = random.Random(args.seed + 101)
+        sources = [rng.choice(list(graph.nodes())) for graph in graphs]
+
+    for i, (graph, source) in enumerate(zip(graphs, sources), start=1):
+        png, pdf, node_count, annular_cycle_count = draw_case(
+            out_dir, i, graph, source, augment=not args.no_augment_same_level_faces
+        )
+        print(
+            f"case {i}: source={source}, connectivity={nx.node_connectivity(graph)}, "
+            f"tire nodes={node_count}, annular cycles={annular_cycle_count}"
+        )
+        print(f"  wrote {png}")
+        print(f"  wrote {pdf}")
+
+
+def main():
+    parser = argparse.ArgumentParser(description=__doc__)
+    parser.add_argument("--n", type=int, default=30)
+    parser.add_argument("--count", type=int, default=2)
+    parser.add_argument("--seed", type=int, default=20260615)
+    parser.add_argument("--scan-limit", type=int, default=500)
+    parser.add_argument("--graph6", help="draw this graph6 graph instead of sampling")
+    parser.add_argument("--source", type=int, help="source vertex for --graph6")
+    parser.add_argument(
+        "--no-augment-same-level-faces",
+        action="store_true",
+        help="skip the same-level-face vertex insertion step",
+    )
+    parser.add_argument(
+        "--out-dir",
+        default=str(PAPER_DIR / "experiments" / "random_medial_tire_decompositions"),
+    )
+    run(parser.parse_args())
+
+
+if __name__ == "__main__":
+    main()