{ "cells": [ { "cell_type": "markdown", "id": "5f835a6a", "metadata": {}, "source": [ "# **multimoda-rs**: CCTA Module Tutorial\n", "\n", "This notebook demonstrates the full CCTA–intravascular fusion pipeline:\n", "\n", "1. Read in and label a CCTA geometry from an STL file and centerline CSV files\n", "2. Prepare centerlines, detect branches, and discretize the vessel tree\n", "3. Load an intravascular geometry and fine-align it to the CCTA point cloud\n", "4. Label the anomalous (intramural) region within the CCTA mesh\n", "5. Compute radial scaling factors for the proximal, distal, and aortic regions\n", "6. Morph the CCTA mesh to match the intravascular geometry\n", "7. Remove the intramural region and stitch the CCTA to the intravascular geometry\n", "8. Remesh and smooth the stitched geometry\n", "9. Visualise labeled regions and export section STL files\n", "10. Re-label the final stitched geometry\n", "\n", "> **Note**: Example data is loaded automatically from `examples/data/` when run from the repository root." ] }, { "cell_type": "markdown", "id": "8e4dfe8f", "metadata": {}, "source": [ "Load required packages and set the working directory. Performance is typically faster in the console than in a notebook.\n" ] }, { "cell_type": "code", "execution_count": 1, "id": "3e81c481", "metadata": { "execution": { "iopub.execute_input": "2026-05-04T15:10:57.242463Z", "iopub.status.busy": "2026-05-04T15:10:57.242255Z", "iopub.status.idle": "2026-05-04T15:11:15.856292Z", "shell.execute_reply": "2026-05-04T15:11:15.855584Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Working directory: /mnt/c/WorkingData/Documents/2_Coding/Rust/multimodars/examples/data\n" ] } ], "source": [ "import os\n", "from pathlib import Path\n", "import numpy as np\n", "import multimodars as mm\n", "\n", "cwd = Path.cwd()\n", "for candidate in [cwd, cwd.parent, cwd.parent.parent]:\n", " if (candidate / \"examples\" / \"data\").exists():\n", " os.chdir(candidate / \"examples\" / \"data\")\n", " break\n", " elif (candidate / \"data\").exists():\n", " os.chdir(candidate / \"data\")\n", " break\n", "print(f\"Working directory: {os.getcwd()}\")\n" ] }, { "cell_type": "code", "execution_count": 2, "id": "3e8428c5", "metadata": { "execution": { "iopub.execute_input": "2026-05-04T15:11:15.858537Z", "iopub.status.busy": "2026-05-04T15:11:15.858168Z", "iopub.status.idle": "2026-05-04T15:11:23.984837Z", "shell.execute_reply": "2026-05-04T15:11:23.983925Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Note: you may need to restart the kernel to use updated packages.\n" ] } ], "source": [ "# Install if needed\n", "%pip --disable-pip-version-check install trimesh plotly scipy pymeshlab | grep -v 'already satisfied'\n", "\n", "import trimesh\n", "import plotly.graph_objects as go\n", "import plotly.io as pio\n", "pio.renderers.default = \"notebook_connected\"\n", "\n", "_DARK = dict(\n", " scene=dict(aspectmode='data', bgcolor='black'),\n", " margin=dict(r=0, l=0, b=0, t=30),\n", " legend=dict(yanchor=\"top\", y=0.99, xanchor=\"left\", x=0.01),\n", " paper_bgcolor=\"black\",\n", ")\n", "\n", "\n", "def trimesh_to_mesh3d(mesh, color, name, opacity=0.6):\n", " verts, faces = mesh.vertices, mesh.faces\n", " return go.Mesh3d(\n", " x=verts[:,0], y=verts[:,1], z=verts[:,2],\n", " i=faces[:,0], j=faces[:,1], k=faces[:,2],\n", " color=color, opacity=opacity, name=name, flatshading=True,\n", " )\n", "\n", "\n", "def _to_arr(pts):\n", " return np.array(pts, dtype=np.float64) if pts else np.empty((0, 3))\n", "\n", "\n", "def plot_labeled_geometry(results):\n", " \"\"\"Visualise initial geometry labeling (aorta / RCA / LCA / removed).\"\"\"\n", " mesh = results[\"mesh\"]\n", " fig = go.Figure()\n", " fig.add_trace(go.Mesh3d(\n", " x=mesh.vertices[:,0], y=mesh.vertices[:,1], z=mesh.vertices[:,2],\n", " i=mesh.faces[:,0], j=mesh.faces[:,1], k=mesh.faces[:,2],\n", " opacity=0.1, color='lightgray', name='Mesh', hoverinfo='skip',\n", " ))\n", " for key, color, label in [\n", " (\"aorta_points\", \"yellow\", \"Aorta\"),\n", " (\"rca_points\", \"blue\", \"RCA\"),\n", " (\"lca_points\", \"green\", \"LCA\"),\n", " (\"rca_removed_points\", \"red\", \"RCA removed (intramural)\"),\n", " (\"lca_removed_points\", \"red\", \"LCA removed (intramural)\"),\n", " ]:\n", " arr = _to_arr(results.get(key, []))\n", " if len(arr):\n", " fig.add_trace(go.Scatter3d(\n", " x=arr[:,0], y=arr[:,1], z=arr[:,2],\n", " mode='markers', marker=dict(size=2, color=color), name=label,\n", " ))\n", " fig.update_layout(**_DARK)\n", " fig.show()\n", "\n", "\n", "def plot_anomalous_regions(results):\n", " \"\"\"Visualise proximal / anomalous / distal sub-regions.\"\"\"\n", " mesh = results[\"mesh\"]\n", " fig = go.Figure()\n", " fig.add_trace(go.Mesh3d(\n", " x=mesh.vertices[:,0], y=mesh.vertices[:,1], z=mesh.vertices[:,2],\n", " i=mesh.faces[:,0], j=mesh.faces[:,1], k=mesh.faces[:,2],\n", " opacity=0.1, color='lightgray', name='Mesh', hoverinfo='skip',\n", " ))\n", " for key, color, label in [\n", " (\"proximal_points\", \"cyan\", \"Proximal\"),\n", " (\"anomalous_points\", \"orange\", \"Anomalous\"),\n", " (\"distal_points\", \"magenta\", \"Distal\"),\n", " ]:\n", " arr = _to_arr(results.get(key, []))\n", " if len(arr):\n", " fig.add_trace(go.Scatter3d(\n", " x=arr[:,0], y=arr[:,1], z=arr[:,2],\n", " mode='markers', marker=dict(size=2, color=color), name=label,\n", " ))\n", " fig.update_layout(**_DARK)\n", " fig.show()\n", "\n", "\n", "def plot_results_key(results):\n", " \"\"\"Plotly replacement for mm.plot_results_key — all region labels on one figure.\"\"\"\n", " mesh = results[\"mesh\"]\n", " fig = go.Figure()\n", " fig.add_trace(go.Mesh3d(\n", " x=mesh.vertices[:,0], y=mesh.vertices[:,1], z=mesh.vertices[:,2],\n", " i=mesh.faces[:,0], j=mesh.faces[:,1], k=mesh.faces[:,2],\n", " opacity=0.15, color='lightgray', name='Mesh', hoverinfo='skip',\n", " ))\n", " for key, color, label in [\n", " (\"aorta_points\", \"yellow\", \"Aorta\"),\n", " (\"rca_points\", \"blue\", \"RCA\"),\n", " (\"lca_points\", \"green\", \"LCA\"),\n", " (\"rca_removed_points\", \"red\", \"Intramural wall\"),\n", " (\"proximal_points\", \"cyan\", \"Proximal\"),\n", " (\"distal_points\", \"magenta\", \"Distal\"),\n", " (\"anomalous_points\", \"orange\", \"Anomalous\"),\n", " ]:\n", " arr = _to_arr(results.get(key, []))\n", " if len(arr):\n", " fig.add_trace(go.Scatter3d(\n", " x=arr[:,0], y=arr[:,1], z=arr[:,2],\n", " mode='markers', marker=dict(size=2, color=color), name=label,\n", " ))\n", " fig.update_layout(**_DARK)\n", " fig.show()\n", "\n", "\n", "def plot_centerline_branches(rca_cl, lca_cl, results_dict=None):\n", " \"\"\"Visualise centerline branch assignments for RCA and LCA.\n", "\n", " Each branch gets its own toggleable legend entry. When *results_dict*\n", " is provided, labelled surface-mesh points are overlaid so you can verify\n", " that branch labels transferred correctly to the geometry.\n", "\n", " Colour coding: blues = RCA (dark = main, lighter = side branches);\n", " reds/oranges = LCA.\n", " \"\"\"\n", " _RCA = [\"#1f77b4\", \"#17becf\", \"#9467bd\", \"#2ca02c\", \"#7f7f7f\"]\n", " _LCA = [\"#d62728\", \"#ff7f0e\", \"#e377c2\", \"#bcbd22\", \"#8c564b\"]\n", "\n", " fig = go.Figure()\n", "\n", " def add_cl(cl, colors, vessel):\n", " from collections import defaultdict\n", " by_branch = defaultdict(list)\n", " for p in cl.points:\n", " by_branch[p.branch_id].append(p.contour_point)\n", " for bid in sorted(by_branch):\n", " pts = by_branch[bid]\n", " label = f\"{vessel} main\" if bid == 0 else f\"{vessel} branch {bid}\"\n", " fig.add_trace(go.Scatter3d(\n", " x=[p.x for p in pts], y=[p.y for p in pts], z=[p.z for p in pts],\n", " mode='markers', marker=dict(color=colors[bid % len(colors)], size=3),\n", " name=label,\n", " ))\n", "\n", " def add_surf(key, color, label):\n", " arr = _to_arr(results_dict.get(key, []))\n", " if len(arr):\n", " fig.add_trace(go.Scatter3d(\n", " x=arr[:,0], y=arr[:,1], z=arr[:,2],\n", " mode='markers', marker=dict(color=color, size=1.5, opacity=0.3),\n", " name=label,\n", " ))\n", "\n", " add_cl(rca_cl, _RCA, \"RCA\")\n", " add_cl(lca_cl, _LCA, \"LCA\")\n", "\n", " if results_dict is not None:\n", " add_surf(\"rca_points_main\", _RCA[0], \"RCA main pts\")\n", " i = 1\n", " while f\"rca_points_side_{i}\" in results_dict:\n", " add_surf(f\"rca_points_side_{i}\", _RCA[i % len(_RCA)], f\"RCA side {i} pts\")\n", " i += 1\n", " add_surf(\"lca_points_main\", _LCA[0], \"LCA main pts\")\n", " i = 1\n", " while f\"lca_points_side_{i}\" in results_dict:\n", " add_surf(f\"lca_points_side_{i}\", _LCA[i % len(_LCA)], f\"LCA side {i} pts\")\n", " i += 1\n", "\n", " fig.update_layout(title=\"Centerline Branch Assignment\", **_DARK)\n", " fig.show()\n", "\n", "\n", "def plot_vessel_tree(tree, pts_per_contour=24):\n", " \"\"\"Visualise a discretized vessel tree.\n", "\n", " Renders contour rings (subsampled to *pts_per_contour* points),\n", " yellow centroid dots tracing each vessel axis, and orientation\n", " reference triplets at the ostium and bifurcation sites.\n", "\n", " Colour coding: silver = aorta; steel-blue = RCA main; shades of blue =\n", " RCA side branches; coral = LCA main; shades of orange = LCA side branches.\n", " Red × = main reference; orange ▲ = CCW reference; magenta ▼ = CW reference.\n", " \"\"\"\n", " _RCA = [\"#1f77b4\", \"#4fa3e0\", \"#7ec8e3\", \"#a8d8ea\"]\n", " _LCA = [\"#d62728\", \"#e07f4f\", \"#e3a87e\", \"#eac0a8\"]\n", "\n", " fig = go.Figure()\n", "\n", " def ring_xyz(contour):\n", " pts = contour.points\n", " m = len(pts)\n", " if m == 0:\n", " return [], [], []\n", " step = max(1, m // pts_per_contour)\n", " s = pts[::step]\n", " return ([p.x for p in s] + [s[0].x],\n", " [p.y for p in s] + [s[0].y],\n", " [p.z for p in s] + [s[0].z])\n", "\n", " def add_contours(contours, color, name, opacity=0.7):\n", " rx, ry, rz, cx, cy, cz = [], [], [], [], [], []\n", " for c in contours:\n", " xs, ys, zs = ring_xyz(c)\n", " if not xs:\n", " continue\n", " rx += xs + [None]; ry += ys + [None]; rz += zs + [None]\n", " if c.centroid:\n", " cx.append(c.centroid[0]); cy.append(c.centroid[1]); cz.append(c.centroid[2])\n", " if rx:\n", " fig.add_trace(go.Scatter3d(x=rx, y=ry, z=rz, mode='lines',\n", " name=name, line=dict(color=color, width=1), opacity=opacity))\n", " if cx:\n", " fig.add_trace(go.Scatter3d(x=cx, y=cy, z=cz, mode='markers',\n", " marker=dict(color='yellow', size=2), showlegend=False,\n", " name=f\"{name} centroids\"))\n", "\n", " def add_refs(refs, label):\n", " mx, my, mz, ccx, ccy, ccz, clx, cly, clz = [], [], [], [], [], [], [], [], []\n", " for main_ref, cc_ref, clock_ref in refs:\n", " mx.append(main_ref[0]); my.append(main_ref[1]); mz.append(main_ref[2])\n", " ccx.append(cc_ref[0]); ccy.append(cc_ref[1]); ccz.append(cc_ref[2])\n", " clx.append(clock_ref[0]); cly.append(clock_ref[1]); clz.append(clock_ref[2])\n", " if mx:\n", " fig.add_trace(go.Scatter3d(x=mx, y=my, z=mz, mode='markers',\n", " name=f\"{label} main ref\", marker=dict(color='red', symbol='x', size=5)))\n", " fig.add_trace(go.Scatter3d(x=ccx, y=ccy, z=ccz, mode='markers',\n", " name=f\"{label} CCW ref\", marker=dict(color='orange', symbol='diamond', size=5)))\n", " fig.add_trace(go.Scatter3d(x=clx, y=cly, z=clz, mode='markers',\n", " name=f\"{label} CW ref\", marker=dict(color='magenta', symbol='square', size=5)))\n", "\n", " add_contours(tree.discretized_aorta, \"silver\", \"Aorta\", opacity=0.3)\n", " add_contours(tree.discretized_rca_main, \"steelblue\", \"RCA main\")\n", " for i, branch in enumerate(tree.rca_branches):\n", " add_contours(branch, _RCA[i % len(_RCA)], f\"RCA branch {i+1}\")\n", " add_contours(tree.discretized_lca_main, \"coral\", \"LCA main\")\n", " for i, branch in enumerate(tree.lca_branches):\n", " add_contours(branch, _LCA[i % len(_LCA)], f\"LCA branch {i+1}\")\n", " add_refs(tree.rca_references, \"RCA\")\n", " add_refs(tree.lca_references, \"LCA\")\n", "\n", " fig.update_layout(title=\"Discretized Vessel Tree\", **_DARK)\n", " fig.show()\n", "\n", "def plot_centerline_edges(cl, cos_threshold=0.0, title=\"Centerline Edges\"):\n", " \"\"\"Show centerline branches with sharp-angle positions as red x markers.\n", "\n", " Each branch gets a distinct colour. Large red x markers indicate positions\n", " from find_sharp_angles — use them to decide where to call\n", " split_branch or merge_branches.\n", " \"\"\"\n", " _PALETTE = [\n", " \"#1f77b4\", \"#ff7f0e\", \"#2ca02c\", \"#d62728\", \"#9467bd\",\n", " \"#8c564b\", \"#e377c2\", \"#7f7f7f\", \"#bcbd22\", \"#17becf\",\n", " ]\n", " from collections import defaultdict\n", " by_branch = defaultdict(list)\n", " for p in cl.points:\n", " by_branch[p.branch_id].append(p.contour_point)\n", "\n", " fig = go.Figure()\n", " for bid in sorted(by_branch):\n", " pts = by_branch[bid]\n", " color = _PALETTE[bid % len(_PALETTE)]\n", " label = \"main\" if bid == 0 else f\"branch {bid}\"\n", " fig.add_trace(go.Scatter3d(\n", " x=[p.x for p in pts], y=[p.y for p in pts], z=[p.z for p in pts],\n", " mode='lines+markers', name=label,\n", " line=dict(color=color, width=2),\n", " marker=dict(color=color, size=3),\n", " ))\n", " sharp_pos = cl.find_sharp_angles(bid, cos_threshold)\n", " if sharp_pos:\n", " sharp_pts = [pts[i] for i in sharp_pos if i < len(pts)]\n", " fig.add_trace(go.Scatter3d(\n", " x=[p.x for p in sharp_pts], y=[p.y for p in sharp_pts], z=[p.z for p in sharp_pts],\n", " mode='markers', name=f\"{label} sharp angles\",\n", " marker=dict(color='red', size=8, symbol='x'),\n", " ))\n", " fig.update_layout(title=title, **_DARK)\n", " fig.show()" ] }, { "cell_type": "markdown", "id": "68985b25", "metadata": {}, "source": [ "## 1. Read in and prepare CCTA geometries and centerlines\n", "\n", "The entry point for the CCTA module is `mm.label_geometry`, which reads a triangulated surface mesh\n", "(STL) together with centerline CSV files for the aorta, RCA, and LCA. It returns a labeled `results`\n", "dictionary and three `PyCenterline` objects used in all subsequent steps.\n", "\n", "The centerline CSV files must contain three columns (no header): `x`, `y`, `z` in mm.\n", "\n", "**Key parameters:**\n", "- `bounding_sphere_radius_mm`: radius of the rolling sphere for the initial vessel labeling pass.\n", " Larger values capture more distant vertices; smaller values are more conservative.\n", "- `n_points_intramural`: number of centerline points defining the end of the intramural segment.\n", "- `anomalous_rca` / `anomalous_lca`: when `True`, re-assigns incorrectly labeled intramural points.\n", " Set to `False` for normal coronary anatomy.\n", "- `control_plot`: open an interactive scene to inspect labeling. Useful when tuning the parameters above.\n", "\n", "The returned `results` dictionary contains keys `\"mesh\"`, `\"aorta_points\"`, `\"rca_points\"`,\n", "`\"lca_points\"`, `\"rca_removed_points\"`, and `\"lca_removed_points\"`.\n" ] }, { "cell_type": "code", "execution_count": 3, "id": "b33a0ed3", "metadata": { "execution": { "iopub.execute_input": "2026-05-04T15:11:23.987453Z", "iopub.status.busy": "2026-05-04T15:11:23.987223Z", "iopub.status.idle": "2026-05-04T15:11:35.123454Z", "shell.execute_reply": "2026-05-04T15:11:35.122613Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Loaded mesh: 25171 vertices, 50338 faces\n", "Loaded aorta centerline: 75 points\n", "Loaded LCA centerline: 1092 points\n", "Loaded RCA centerline: 788 points\n", "\n", "RCA points found: 4820\n", "LCA points found: 5899\n", "Applying occlusion removal for anomalous RCA...\n", "Total faces to exclude: 232\n", "Excluded 232 faces, removed 110 points (filtered from 4820 to 4710 points)\n", "RCA: relabeled 110 points in intramual course\n", "\n", "Removing LCA and RCA island points...\n", "length before: 5899\n", "length after: 5879\n", "\n", "Applying final reclassification based on adjacency map...\n", "aorta_points:14472\n", "rca_points:4710\n", "lca_points:5879\n", "rca_removed_points:110\n", "lca_removed_points:0\n" ] }, { "data": { "text/html": [ " \n", " \n", " " ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/html": [ "
\n", "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "rca_cl_raw = np.genfromtxt(\"./centerline_rca_short.csv\", delimiter=',')\n", "lca_cl_raw = np.genfromtxt(\"./centerline_lca.csv\", delimiter=',')\n", "aorta_cl_raw = np.genfromtxt(\"./centerline_aorta.csv\", delimiter=',')\n", "rca_cl = mm.numpy_to_centerline(rca_cl_raw)\n", "lca_cl = mm.numpy_to_centerline(lca_cl_raw)\n", "aorta_cl = mm.numpy_to_centerline(aorta_cl_raw)\n", "\n", "results, (rca_cl, lca_cl, ao_cl) = mm.label_geometry(\n", " path_ccta_geometry=\"./NARCO_119_noside.stl\",\n", " path_centerline_aorta=\"./centerline_aorta.csv\",\n", " path_centerline_rca=\"./centerline_rca_short.csv\",\n", " path_centerline_lca=\"./centerline_lca.csv\",\n", " bounding_sphere_radius_mm=3.0,\n", " n_points_intramural=100,\n", " anomalous_rca=True,\n", " anomalous_lca=False,\n", " control_plot=False,\n", ")\n", "\n", "plot_labeled_geometry(results)\n" ] }, { "cell_type": "markdown", "id": "61175f0c", "metadata": {}, "source": [ "## 2. Discretize CCTA and map reference points\n", "\n", "`prepare_centerlines` detects branches on both coronary centerlines, validates topology, and labels\n", "the surface-mesh points by branch — all in one call. The returned `results` dictionary gains keys\n", "`rca_points_main`, `rca_points_side_1`, …, `lca_points_main`, `lca_points_side_1`, ….\n", "\n", "After verifying the branch assignment in the interactive plot, `discretize_vessel_tree` slices each\n", "vessel into evenly-spaced cross-sectional contours and computes orientation reference triplets at\n", "the ostium and every side-branch bifurcation.\n", "\n", "**Key parameters:**\n", "- `branch_sigma`: Gaussian smoothing radius (mm) for branch detection — increase if the algorithm\n", " over-segments noisy centerlines.\n", "- `step_size`: arc-length spacing between consecutive cross-sections (mm).\n", "- `n_points`: number of evenly-spaced points per output contour ring.\n", "- `b_spline`: replace each discretized contour with a closed periodic B-spline fit — useful for\n", " smoothing jagged contours before computing reference points. Control strength with\n", " `bspline_smoothing` (0 = exact; ≈ `n_points` = gentle; ≈ `5 × n_points` = strong)." ] }, { "cell_type": "code", "execution_count": 4, "id": "eee97270", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "Branch labeling for 'rca_points' (branch_ids=[0]):\n", " rca_points_main: 3749\n", " rca_points_side: 961\n", " rca_points_side_1: 515\n", " rca_points_side_2: 423\n", " rca_points_side_3: 23\n", "\n", "Branch labeling for 'lca_points' (branch_ids=[0]):\n", " lca_points_main: 5228\n", " lca_points_side: 651\n", " lca_points_side_1: 248\n", " lca_points_side_2: 277\n", " lca_points_side_3: 114\n", " lca_points_side_4: 12\n" ] }, { "data": { "text/html": [ "
\n", "
" ] }, "metadata": {}, "output_type": "display_data" }, { "name": "stdout", "output_type": "stream", "text": [ "\n", "Branch labeling for 'lca_points' (branch_ids=[0]):\n", " lca_points_main: 2623\n", " lca_points_side: 3256\n", " lca_points_side_1: 496\n", " lca_points_side_2: 277\n", " lca_points_side_3: 114\n", " lca_points_side_4: 2617\n" ] }, { "data": { "text/html": [ "
\n", "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "rca_cl, lca_cl, results = mm.prepare_centerlines(\n", " rca_cl, lca_cl, results,\n", " branch_sigma=2.0,\n", ")\n", "\n", "# Inspect each centerline for sharp angles — red × marks a position that may\n", "# need a split_branch / merge_branches fix. Tune cos_threshold toward negative\n", "# values (e.g. -0.3) to suppress false positives on gently curved segments.\n", "plot_centerline_edges(lca_cl, cos_threshold=0.0, title=\"LCA edges\")\n", "\n", "# Manual corrections (uncomment and adapt as needed):\n", "list_edges = lca_cl.find_sharp_angles(branch_id=0, cos_threshold=0.0)\n", "lca_cl = lca_cl.split_branch(0, list_edges[4])\n", "lca_cl = lca_cl.merge_branches(0, 4)\n", "lca_cl = lca_cl.check_centerline()\n", "results = mm.label_branches(lca_cl, results, results_key=\"lca_points\")\n", "\n", "# Verify final branch assignment after any corrections.\n", "plot_centerline_branches(rca_cl, lca_cl, results)\n" ] }, { "cell_type": "code", "execution_count": 5, "id": "edd07aa7", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "tree = mm.discretize_vessel_tree(\n", " ao_cl, rca_cl, lca_cl,\n", " results,\n", " step_size=1.0,\n", " n_points=100,\n", " b_spline=True, # set True + tune bspline_smoothing to smooth noisy contours\n", " bspline_smoothing=5.0,\n", ")\n", "\n", "plot_vessel_tree(tree)" ] }, { "cell_type": "markdown", "id": "28ae6ef4", "metadata": {}, "source": [ "## 3. Load and align intravascular geometry\n", "\n", "Load the intravascular segmentation with `mm.from_file_singlepair` (see the intravascular tutorial\n", "for the full range of loading options and parameter tuning). Then align it to the CCTA centerline\n", "and point cloud with `mm.align_combined`, which first performs a coarse three-point alignment using\n", "anatomical landmarks and then refines the rotation by minimising Hausdorff distances.\n", "\n", "**Landmark points** (in mm, CCTA coordinate system):\n", "- Aortic reference — centre of vessel on the aortic side\n", "- Superior reference — proximal end of the intramural segment\n", "- Inferior reference — distal end of the intramural segment\n", "\n", "**Key parameters for `align_combined`:**\n", "- `angle_range_deg`: angular search window (±degrees) for Hausdorff refinement.\n", "- `write` / `watertight` / `output_dir`: export OBJ meshes when `write=True`." ] }, { "cell_type": "code", "execution_count": 6, "id": "fdb35527", "metadata": { "execution": { "iopub.execute_input": "2026-05-04T15:11:35.139859Z", "iopub.status.busy": "2026-05-04T15:11:35.139627Z", "iopub.status.idle": "2026-05-04T15:11:45.784349Z", "shell.execute_reply": "2026-05-04T15:11:45.783619Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "✅ Successfully built geometry from path\n", "-----------------------------------------\n", "✅ Lumen\n", "❌ Eem\n", "❌ Calcification\n", "❌ Sidebranch\n", "✅ Catheter\n", "-----------------------------------------\n", "Label: aligned_dia\n", "Diastole phase: Yes\n", "\n", "\n", "✅ Successfully built geometry from path\n", "-----------------------------------------\n", "✅ Lumen\n", "❌ Eem\n", "❌ Calcification\n", "❌ Sidebranch\n", "✅ Catheter\n", "-----------------------------------------\n", "Label: aligned_sys\n", "Diastole phase: No\n", "\n", "\n", "+--------------------------------------------------------------------+\n", "| ✅ Finished aligning 'aligned_sys' (anomalous: true) |\n", "+---------+------------+---------------+-------+-------+-------------+\n", "| Contour | Matched To | Rotation (°) | Tx | Ty | Centroid |\n", "+---------+------------+---------------+-------+-------+-------------+\n", "| 1 | 0 | -2.50 | -0.19 | -0.01 | (3.63,5.12) |\n", "| 2 | 1 | 12.00 | -0.54 | 0.15 | (3.63,5.12) |\n", "| 3 | 2 | -10.00 | -0.44 | 0.53 | (3.63,5.12) |\n", "| 4 | 3 | -4.00 | -0.72 | 0.94 | (3.63,5.12) |\n", "| 5 | 4 | 1.00 | -0.75 | 0.26 | (3.63,5.12) |\n", "| 6 | 5 | 5.00 | -0.79 | 0.65 | (3.63,5.12) |\n", "| 7 | 6 | -6.50 | -1.21 | 0.69 | (3.63,5.12) |\n", "| 8 | 7 | 15.00 | -1.08 | 0.81 | (3.63,5.12) |\n", "| 9 | 8 | -35.00 | -1.03 | 0.24 | (3.63,5.12) |\n", "| 10 | 9 | 25.50 | -1.46 | 0.28 | (3.63,5.12) |\n", "| 11 | 10 | 40.00 | -2.40 | 0.73 | (3.63,5.12) |\n", "| 12 | 11 | -21.00 | -1.95 | 0.68 | (3.63,5.12) |\n", "| 13 | 12 | 5.00 | -2.06 | 0.66 | (3.63,5.12) |\n", "| 14 | 13 | -7.00 | -2.14 | 0.52 | (3.63,5.12) |\n", "| 15 | 14 | 23.00 | -2.27 | 0.46 | (3.63,5.12) |\n", "| 16 | 15 | 7.50 | -2.28 | 0.53 | (3.63,5.12) |\n", "+---------+------------+---------------+-------+-------+-------------+\n", "\n", "+--------------------------------------------------------------------+\n", "| ✅ Finished aligning 'aligned_dia' (anomalous: true) |\n", "+---------+------------+---------------+-------+-------+-------------+\n", "| Contour | Matched To | Rotation (°) | Tx | Ty | Centroid |\n", "+---------+------------+---------------+-------+-------+-------------+\n", "| 1 | 0 | 7.50 | 0.12 | -0.36 | (3.72,5.25) |\n", "| 2 | 1 | 1.00 | 0.04 | -0.16 | (3.72,5.25) |\n", "| 3 | 2 | 6.00 | -0.28 | 0.03 | (3.72,5.25) |\n", "| 4 | 3 | -13.00 | -0.16 | 0.45 | (3.72,5.25) |\n", "| 5 | 4 | -24.00 | -0.32 | 0.89 | (3.72,5.25) |\n", "| 6 | 5 | -5.00 | -0.48 | 1.10 | (3.72,5.25) |\n", "| 7 | 6 | -29.50 | -0.94 | 0.82 | (3.72,5.25) |\n", "| 8 | 7 | -18.00 | -1.29 | 0.87 | (3.72,5.25) |\n", "| 9 | 8 | 1.50 | -1.36 | 0.86 | (3.72,5.25) |\n", "| 10 | 9 | -1.00 | -1.68 | 1.27 | (3.72,5.25) |\n", "| 11 | 10 | 47.00 | -1.54 | 1.87 | (3.72,5.25) |\n", "| 12 | 11 | 3.50 | -1.45 | 1.95 | (3.72,5.25) |\n", "| 13 | 12 | 1.00 | -1.51 | 2.08 | (3.72,5.25) |\n", "| 14 | 13 | 30.50 | -0.96 | 2.12 | (3.72,5.25) |\n", "| 15 | 14 | -7.00 | -1.11 | 1.93 | (3.72,5.25) |\n", "| 16 | 15 | 19.50 | -0.68 | 2.47 | (3.72,5.25) |\n", "| 17 | 16 | 3.50 | -0.50 | 2.35 | (3.72,5.25) |\n", "| 18 | 17 | -1.00 | -0.58 | 2.24 | (3.72,5.25) |\n", "| 19 | 18 | -11.00 | -0.95 | 2.37 | (3.72,5.25) |\n", "+---------+------------+---------------+-------+-------+-------------+\n", "\n", "✅ Aligned geometry 'aligned_sys' to 'aligned_dia'\n", "-----------------------------------------\n", "Applied initial translation: (0.09, 0.13, -3.94) mm\n", "Found best rotation of 3.50° with parameters: \n", "range: 90.00° \n", "step size: 0.5°\n", "Applied final translation: ( 0, 0.00, 0.00) mm\n", "-----------------------------------------\n", "\n", "Saving files for 'aligned_dia - aligned_sys' to 'output/rest'\n", "LUMEN .obj files: 2/2 written successfully\n", "CATHETER .obj files: 2/2 written successfully\n", "WALL .obj files: 2/2 written successfully\n", "\n", "Step 1: Finding initial rotation via three-point method\n", "---------------------Centerline alignment: Finding optimal rotation---------------------\n", "✅ Best angle found: 219.00°\n", "Step 2: Refining with Hausdorff distance\n", "---------------------Refining alignment with Hausdorff---------------------\n", "Initial rotation: 0.00°, Initial CL index: 22\n", "Refined rotation: 5.00°, Refined CL index: 21, Hausdorff: 5.12\n", "---------------------Applying final transformation---------------------\n", "Total rotation (initial + delta): 224.00°\n", "Moving ostium by 1 centerline points\n", "\n", "Saving files for 'None' to 'test'\n", "LUMEN .obj files: 2/2 written successfully\n", "CATHETER .obj files: 2/2 written successfully\n", "WALL .obj files: 2/2 written successfully\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ "sidebranch file not found, skipping: \"ivus_rest/branch_diastolic_contours.csv\"\n", "eem file not found, skipping: \"ivus_rest/eem_diastolic_contours.csv\"\n", "process_directory: unknown mapping name 'catheter', skipping\n", "calcification file not found, skipping: \"ivus_rest/calcium_diastolic_contours.csv\"\n", "process_directory: unknown mapping name 'catheter', skipping\n", "eem file not found, skipping: \"ivus_rest/eem_systolic_contours.csv\"\n", "sidebranch file not found, skipping: \"ivus_rest/branch_systolic_contours.csv\"\n", "calcification file not found, skipping: \"ivus_rest/calcium_systolic_contours.csv\"\n", "resample_centerline_by_contours: centroid_count=14, centroid_mean_spacing=Some(1.2955758552631584), centerline_length=185.4959091541807, spacing=1.295576\n", "resample_centerline_by_contours: produced 144 points\n" ] } ], "source": [ "rest, (dia_logs, sys_logs) = mm.from_file_singlepair(\n", " input_path=\"ivus_rest\",\n", " labels=[\"aligned_dia\", \"aligned_sys\"],\n", " output_path=\"output/rest\",\n", ")\n", "\n", "ref_points = tree.rca_references[0]\n", "\n", "rca_cl_main = rca_cl.get_branch(0) # alignment needs single-branch CL\n", "aligned, resampled_cl = mm.align_combined(\n", " rca_cl_main,\n", " rest,\n", " ref_points[0], # aortic reference point\n", " ref_points[1], # superior reference point\n", " ref_points[2], # inferior reference point\n", " results['rca_points'], # CCTA point cloud for Hausdorff refinement\n", " angle_range_deg=30.0,\n", " write=True,\n", " watertight=False,\n", " output_dir=\"test\",\n", " align_wall_anomalous=True,\n", ")\n" ] }, { "cell_type": "markdown", "id": "edd186c3", "metadata": {}, "source": [ "## 4. Label the anomalous region\n", "\n", "`mm.label_anomalous_region` subdivides the RCA points into three sub-regions — proximal, anomalous\n", "(intramural), and distal — based on spatial overlap between the aligned intravascular frames and the\n", "CCTA mesh. The `results` dictionary is extended with `\"proximal_points\"`, `\"anomalous_points\"`,\n", "and `\"distal_points\"`.\n", "\n", "Set `debug_plot=True` to open an interactive scene that shows the sub-region boundaries — useful\n", "when a boundary appears misplaced." ] }, { "cell_type": "code", "execution_count": 7, "id": "88731123", "metadata": { "execution": { "iopub.execute_input": "2026-05-04T15:11:45.786463Z", "iopub.status.busy": "2026-05-04T15:11:45.786247Z", "iopub.status.idle": "2026-05-04T15:11:46.267849Z", "shell.execute_reply": "2026-05-04T15:11:46.267175Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "Applying anomalous labeling based on aligned intravascular frames...\n", "proximal_points: 36\n", "distal_points: 4273\n", "anomalous_points: 401\n" ] }, { "data": { "text/html": [ "
\n", "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "results = mm.label_anomalous_region(\n", " centerline=rca_cl,\n", " frames=aligned.geom_a.frames,\n", " results=results,\n", " results_key='rca_points',\n", " debug_plot=False,\n", ")\n", "\n", "plot_anomalous_regions(results)\n" ] }, { "cell_type": "markdown", "id": "d45f88a6", "metadata": {}, "source": [ "## 5. Compute scaling factors\n", "\n", "Before morphing, the optimal radial scaling factor for each region is computed independently.\n", "Each function searches for the scale that minimises distance between the CCTA mesh and the\n", "corresponding portion of the aligned intravascular geometry. All return values are signed floats\n", "in mm (+ expands, − contracts).\n", "\n", "| Function | Region |\n", "|---|---|\n", "| `find_distal_and_proximal_scaling` | proximal and distal segments |\n", "| `find_aorta_scaling` | aortic vertices, minimising distance to outer intravascular wall |\n", "| `find_aortic_wall_scaling` | aortic *wall* vertices, targeting the intramural→free-segment transition |\n", "\n", "> **Note:** `find_aortic_wall_scaling` raises `ValueError` if no frame with elliptic ratio < 1.3\n", "> is found (nearly circular geometry). Omit this step in that case." ] }, { "cell_type": "code", "execution_count": 8, "id": "041ba58d", "metadata": { "execution": { "iopub.execute_input": "2026-05-04T15:11:46.276799Z", "iopub.status.busy": "2026-05-04T15:11:46.276584Z", "iopub.status.idle": "2026-05-04T15:11:46.509792Z", "shell.execute_reply": "2026-05-04T15:11:46.509126Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "Finding best proximal/distal radial scaling factors...\n", "Proximal scaling: -0.4 mm\n", "Distal scaling: 0.4 mm\n", "\n", "Finding best aortic radial scaling factor...\n", "Aortic scaling: -0.4 mm\n", "\n", "Finding best aortic wall radial scaling factor...\n", "elliptic ratio <1.3 for frame index 7\n", "Aortic wall scaling: 1.16 mm\n", "Proximal scaling: -0.400 mm\n", "Distal scaling: 0.400 mm\n", "Aortic scaling: -0.400 mm\n", "Aortic wall scaling: 1.158 mm\n" ] } ], "source": [ "prox_scaling, distal_scaling = mm.find_distal_and_proximal_scaling(\n", " frames=aligned.geom_a.frames,\n", " centerline=rca_cl,\n", " results=results,\n", ")\n", "\n", "aortic_scaling = mm.find_aorta_scaling(\n", " frames=aligned.geom_a.frames,\n", " cl_aorta=ao_cl,\n", " results=results,\n", ")\n", "\n", "aortic_wall_scaling = mm.find_aortic_wall_scaling(\n", " frames=aligned.geom_a.frames,\n", " cl_aorta=ao_cl,\n", " results=results,\n", ")\n", "\n", "print(f\"Proximal scaling: {prox_scaling:.3f} mm\")\n", "print(f\"Distal scaling: {distal_scaling:.3f} mm\")\n", "print(f\"Aortic scaling: {aortic_scaling:.3f} mm\")\n", "print(f\"Aortic wall scaling: {aortic_wall_scaling:.3f} mm\")\n" ] }, { "cell_type": "markdown", "id": "10647148", "metadata": {}, "source": [ "## 6. Morph CCTA to intravascular geometry\n", "\n", "`mm.scale_region_centerline_morphing` moves mesh vertices radially along the local centerline normal\n", "so that the region diameter changes by `diameter_adjustment_mm`. After each call,\n", "`mm.sync_results_to_mesh` updates all coordinate lists in `results` to reflect the new vertex\n", "positions. **Always sync before the next morphing step.**\n", "\n", "The three regions are morphed in order: distal → aortic → proximal." ] }, { "cell_type": "code", "execution_count": 9, "id": "39693f96", "metadata": { "execution": { "iopub.execute_input": "2026-05-04T15:11:46.511863Z", "iopub.status.busy": "2026-05-04T15:11:46.511671Z", "iopub.status.idle": "2026-05-04T15:11:47.408124Z", "shell.execute_reply": "2026-05-04T15:11:47.407427Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "Scaling 4273 vertices around Centerline(len=786, spacing=0.57 mm)\n", "Diameter adjustment: 0.4 mm\n", "\n", "Scaling 14582 vertices around Centerline(len=75, spacing=0.75 mm)\n", "Diameter adjustment: -0.4 mm\n", "\n", "Scaling 36 vertices around Centerline(len=786, spacing=0.57 mm)\n", "Diameter adjustment: -0.4 mm\n" ] }, { "data": { "text/html": [ "
\n", "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# 1. Scale the distal segment along the RCA centerline\n", "scaled_distal = mm.scale_region_centerline_morphing(\n", " mesh=results['mesh'],\n", " region_points=results['distal_points'],\n", " centerline=rca_cl,\n", " diameter_adjustment_mm=distal_scaling,\n", ")\n", "results = mm.sync_results_to_mesh(results, results['mesh'], scaled_distal)\n", "\n", "# 2. Scale the aortic region (aorta + intramural wall) along the aortic centerline\n", "scaled_distal_aortic = mm.scale_region_centerline_morphing(\n", " mesh=results['mesh'],\n", " region_points=results['aorta_points'] + results['rca_removed_points'],\n", " centerline=aorta_cl,\n", " diameter_adjustment_mm=aortic_scaling,\n", ")\n", "results = mm.sync_results_to_mesh(results, results['mesh'], scaled_distal_aortic)\n", "\n", "# 3. Scale the proximal segment along the RCA centerline\n", "scaled_proximal = mm.scale_region_centerline_morphing(\n", " mesh=results['mesh'],\n", " region_points=results['proximal_points'],\n", " centerline=rca_cl,\n", " diameter_adjustment_mm=prox_scaling,\n", ")\n", "results = mm.sync_results_to_mesh(results, results['mesh'], scaled_proximal)\n", "\n", "# Overlay scaled CCTA with the aligned intravascular geometry\n", "anomaly_mesh = trimesh.load(\"test/lumen_000_None.obj\")\n", "anomaly_wall_mesh = trimesh.load(\"test/wall_000_None.obj\")\n", "fig = go.Figure(data=[\n", " trimesh_to_mesh3d(results['mesh'], \"darkred\", \"CCTA scaled\", opacity=0.5),\n", " trimesh_to_mesh3d(anomaly_mesh, \"royalblue\", \"Intravascular lumen\", opacity=0.7),\n", " trimesh_to_mesh3d(anomaly_wall_mesh, \"lightblue\", \"Intravascular wall\", opacity=0.3),\n", "])\n", "fig.update_layout(scene=dict(aspectmode=\"data\"), margin=dict(l=0, r=0, t=30, b=0))\n", "fig.show()\n" ] }, { "cell_type": "markdown", "id": "98b16fc6", "metadata": {}, "source": [ "## 7. Remove intramural region and stitch geometries\n", "\n", "`mm.remove_labeled_points_from_mesh` deletes the anomalous and proximal vertices, opening a boundary\n", "ring at the proximal end of the intravascular segment. The `results` dictionary gains a\n", "`\"boundary_points\"` key.\n", "\n", "`mm.stitch_ccta_to_intravascular` then triangulates a patch connecting this boundary ring to the\n", "proximal contour of the aligned intravascular geometry.\n", "\n", "**`prox_start_mode`:**\n", "- `\"nearest_iv\"` (default) — rotate to the boundary vertex closest to intravascular point 0.\n", "- `\"highest_z\"` — rotate to the boundary vertex with the largest z-coordinate; prefer for straight\n", " intramural segments where the pullback axis is nearly aligned with the image z-axis." ] }, { "cell_type": "markdown", "id": "4eae8a22", "metadata": {}, "source": [ "### Anomalous ostium clamping (`clamp_overshoot`)\n", "\n", "For anomalous coronary arteries the intramural course runs roughly **perpendicular** to the aortic\n", "wall. This means the plane of the intravascular (IV) ring at the ostium is nearly at 90° to the\n", "plane of the aortic boundary ring. Without correction, some aortic boundary vertices end up on the\n", "wrong side of the ostium plane (i.e. inside the coronary lumen), creating a jagged seam.\n", "\n", "When the angle between the two planes exceeds 45°, the following three-step correction is applied\n", "automatically:\n", "\n", "1. **Clamp** — any boundary vertex that has crossed to the intravascular side of the ostium plane\n", " is projected back onto that plane.\n", "2. **Minimum gap** (`clamp_overshoot`, default 0.5 mm) — every boundary vertex (including freshly\n", " clamped ones) is pushed at least `clamp_overshoot` mm away from the ostium plane on the aortic\n", " side. This prevents the stitching patch from meeting the aortic wall at a razor-sharp angle at\n", " the ostium edge.\n", "3. **Layer propagation** — the two rings of aortic mesh vertices adjacent to the boundary ring are\n", " shifted radially *outward* from the coronary axis (within the ostium plane), by 0.1 mm and\n", " 0.2 mm respectively. This prevents those second-layer vertices from sitting closer to the\n", " coronary axis than the boundary ring itself, which would otherwise create a visible ridge just\n", " outside the seam.\n", "\n", "`clamp_overshoot` can be tuned per case. Set it to `0` to disable steps 2 and 3 (pure clamping\n", "only)." ] }, { "cell_type": "code", "execution_count": 10, "id": "7cf465ff", "metadata": { "execution": { "iopub.execute_input": "2026-05-04T15:11:47.420847Z", "iopub.status.busy": "2026-05-04T15:11:47.420626Z", "iopub.status.idle": "2026-05-04T15:11:52.411733Z", "shell.execute_reply": "2026-05-04T15:11:52.411021Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Applying removal of '['anomalous_points', 'proximal_points']'\n", "Removed 437\n", "Created 56 boundary points\n", "Stitching: 100/100 triangles created (n_boundary=37, n_iv=100, step=2, remainder=26)\n", "Stitching: 100/100 triangles created (n_boundary=19, n_iv=100, step=5, remainder=5)\n", "Raw stitched mesh exported → prefixed_mesh.stl\n" ] } ], "source": [ "updated_results = mm.remove_labeled_points_from_mesh(\n", " results,\n", " [\"anomalous_points\", \"proximal_points\"],\n", ")\n", "\n", "stitched = mm.stitch_ccta_to_intravascular(\n", " aligned.geom_a,\n", " updated_results['mesh'],\n", " updated_results,\n", " prox_start_mode=\"highest_z\",\n", " clamp_overshoot=0.5,\n", ")\n", "stitched['mesh'].export(\"prefixed_mesh.stl\")\n", "print(\"Raw stitched mesh exported → prefixed_mesh.stl\")\n" ] }, { "cell_type": "code", "execution_count": 11, "id": "281f85b9", "metadata": { "execution": { "iopub.execute_input": "2026-05-04T15:11:52.413737Z", "iopub.status.busy": "2026-05-04T15:11:52.413536Z", "iopub.status.idle": "2026-05-04T15:11:52.617452Z", "shell.execute_reply": "2026-05-04T15:11:52.616852Z" } }, "outputs": [ { "data": { "text/html": [ "
\n", "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# Visualise the stitching seam:\n", "# boundary ring (red→blue by index) + IV frame-0 lumen (red→blue by index)\n", "boundary_pts = np.array(stitched['prox_boundary_points'], dtype=np.float64)\n", "n_bnd = len(boundary_pts)\n", "\n", "iv_viz = aligned.geom_a.downsample(100).sort_frame_points()\n", "iv_pts = iv_viz.frames[0].lumen.points\n", "iv_xyz = np.array([[p.x, p.y, p.z] for p in iv_pts])\n", "n_iv = len(iv_xyz)\n", "\n", "fig = go.Figure()\n", "fig.add_trace(trimesh_to_mesh3d(stitched['mesh'], \"lightgray\", \"Stitched mesh\", opacity=0.2))\n", "fig.add_trace(go.Scatter3d(\n", " x=boundary_pts[:,0], y=boundary_pts[:,1], z=boundary_pts[:,2],\n", " mode='markers',\n", " marker=dict(size=4, color=list(range(n_bnd)), colorscale='RdBu', showscale=False),\n", " name='CCTA boundary ring',\n", "))\n", "fig.add_trace(go.Scatter3d(\n", " x=iv_xyz[:,0], y=iv_xyz[:,1], z=iv_xyz[:,2],\n", " mode='markers',\n", " marker=dict(size=5, color=list(range(n_iv)), colorscale='RdBu', showscale=False),\n", " name='IV frame-0 lumen',\n", "))\n", "fig.update_layout(scene=dict(aspectmode='data'), margin=dict(l=0, r=0, t=30, b=0))\n", "fig.show()\n" ] }, { "cell_type": "markdown", "id": "c607c287", "metadata": {}, "source": [ "## 8. Remesh and smooth\n", "\n", "`mm.fix_and_remesh_stitched_mesh` applies a three-step repair pipeline:\n", "non-manifold repair → hole filling → isotropic remesh (via `pymeshlab`).\n", "Taubin smoothing then reduces surface noise while preserving overall shape.\n", "\n", "> **Requires:** `pip install 'multimodars[meshlab]'`\n", "\n", "**Key parameters:**\n", "- `target_edge_length_mm`: desired edge length after remeshing (`None` uses the 25th-percentile edge\n", " length of the input mesh).\n", "- `verbose`: print per-step vertex/face counts and watertightness." ] }, { "cell_type": "code", "execution_count": 12, "id": "c2dc01b1", "metadata": { "execution": { "iopub.execute_input": "2026-05-04T15:11:52.623933Z", "iopub.status.busy": "2026-05-04T15:11:52.623734Z", "iopub.status.idle": "2026-05-04T15:12:12.428681Z", "shell.execute_reply": "2026-05-04T15:12:12.428049Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[input ] verts= 26,134 faces= 52,264 watertight=False\n", " non-manifold edges/vertices repaired\n", " holes closed\n", "[after hole fill ] verts= 26,134 faces= 52,264 watertight=True\n", " target edge=0.5000 mm (0.2586% of bbox diag=193.34 mm)\n", "[after remesh ] verts= 58,004 faces=116,004 watertight=True\n", "Watertight? True\n", "Remeshed and smoothed geometry exported → fixed_mesh.stl\n" ] } ], "source": [ "remeshed = stitched.copy()\n", "remeshed['mesh'] = mm.fix_and_remesh_stitched_mesh(\n", " stitched['mesh'],\n", " target_edge_length_mm=0.5,\n", " verbose=True,\n", ")\n", "print(f\"Watertight? {remeshed['mesh'].is_watertight}\")\n", "\n", "trimesh.smoothing.filter_taubin(remeshed['mesh'], lamb=0.6)\n", "remeshed['mesh'].export(\"fixed_mesh.stl\")\n", "print(\"Remeshed and smoothed geometry exported → fixed_mesh.stl\")\n" ] }, { "cell_type": "markdown", "id": "67bef5a2", "metadata": {}, "source": [ "## 9. Visualise and export sections\n", "\n", "Inspect all labeled regions of the stitched result, then export individual anatomical sections as\n", "STL files with `mm.export_section_stl`.\n", "\n", "| Colour | Region |\n", "|--------|--------|\n", "| Yellow | Aorta |\n", "| Blue | RCA |\n", "| Green | LCA |\n", "| Red | Intramural wall (removed) |\n", "| Cyan | Proximal |\n", "| Magenta | Distal |\n", "| Orange | Anomalous |" ] }, { "cell_type": "code", "execution_count": 13, "id": "6c6881dc", "metadata": { "execution": { "iopub.execute_input": "2026-05-04T15:12:12.430819Z", "iopub.status.busy": "2026-05-04T15:12:12.430626Z", "iopub.status.idle": "2026-05-04T15:12:13.062948Z", "shell.execute_reply": "2026-05-04T15:12:13.062257Z" } }, "outputs": [ { "data": { "text/html": [ "
\n", "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "mesh = stitched['mesh']\n", "verts, faces = mesh.vertices, mesh.faces\n", "\n", "intensity = mesh.vertex_normals[:, 2]\n", "\n", "fig = go.Figure(data=[go.Mesh3d(\n", " x=verts[:,0], y=verts[:,1], z=verts[:,2],\n", " i=faces[:,0], j=faces[:,1], k=faces[:,2],\n", " intensity=intensity,\n", " colorscale=\"RdBu\",\n", " opacity=1.0,\n", " flatshading=False,\n", " name=\"Stitched geometry\",\n", ")])\n", "fig.update_layout(scene=dict(aspectmode=\"data\"), margin=dict(l=0, r=0, t=30, b=0))\n", "fig.show()" ] }, { "cell_type": "code", "execution_count": 14, "id": "08feda4b", "metadata": { "execution": { "iopub.execute_input": "2026-05-04T15:12:13.072955Z", "iopub.status.busy": "2026-05-04T15:12:13.072738Z", "iopub.status.idle": "2026-05-04T15:12:13.678383Z", "shell.execute_reply": "2026-05-04T15:12:13.677336Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Exported: all.stl aorta.stl rca.stl lca.stl\n" ] } ], "source": [ "mm.export_section_stl(stitched, \"all\")\n", "mm.export_section_stl(stitched, \"aorta\")\n", "mm.export_section_stl(stitched, \"rca\")\n", "mm.export_section_stl(stitched, \"lca\")\n", "print(\"Exported: all.stl aorta.stl rca.stl lca.stl\")\n" ] }, { "cell_type": "markdown", "id": "b8e37352", "metadata": {}, "source": [ "## 10. Re-label the final geometry\n", "\n", "After remeshing and smoothing, vertex coordinates have changed and the stored point lists are no\n", "longer valid. Re-running `mm.label_geometry` on the exported fixed mesh produces a fresh, consistent\n", "labeling for downstream biomechanical simulation or further analysis." ] }, { "cell_type": "code", "execution_count": 15, "id": "b2a0d2af", "metadata": { "execution": { "iopub.execute_input": "2026-05-04T15:12:13.680963Z", "iopub.status.busy": "2026-05-04T15:12:13.680768Z", "iopub.status.idle": "2026-05-04T15:12:58.067406Z", "shell.execute_reply": "2026-05-04T15:12:58.066666Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Loaded mesh: 58004 vertices, 116004 faces\n", "Loaded aorta centerline: 75 points\n", "Loaded LCA centerline: 1092 points\n", "Loaded RCA centerline: 788 points\n", "\n", "RCA points found: 13674\n", "LCA points found: 14130\n", "Applying occlusion removal for anomalous RCA...\n", "Total faces to exclude: 322\n", "Excluded 322 faces, removed 191 points (filtered from 13674 to 13483 points)\n", "RCA: relabeled 191 points in intramual course\n", "\n", "Removing LCA and RCA island points...\n", "length before: 14130\n", "length after: 14087\n", "\n", "Applying final reclassification based on adjacency map...\n", "aorta_points:30243\n", "rca_points:13483\n", "lca_points:14087\n", "rca_removed_points:191\n", "lca_removed_points:0\n" ] }, { "data": { "text/html": [ "
\n", "
" ] }, "metadata": {}, "output_type": "display_data" }, { "name": "stdout", "output_type": "stream", "text": [ "Final labeled sections exported.\n" ] } ], "source": [ "results_final, (rca_cl_f, lca_cl_f, ao_cl_f) = mm.label_geometry(\n", " path_ccta_geometry=\"fixed_mesh.stl\",\n", " path_centerline_aorta=\"../data/centerline_aorta.csv\",\n", " path_centerline_rca=\"../data/centerline_rca_short.csv\",\n", " path_centerline_lca=\"../data/centerline_lca.csv\",\n", " bounding_sphere_radius_mm=3.0,\n", " n_points_intramural=100,\n", " anomalous_rca=True,\n", " anomalous_lca=False,\n", " control_plot=False,\n", ")\n", "\n", "plot_labeled_geometry(results_final)\n", "\n", "mm.export_section_stl(results_final, \"all\")\n", "mm.export_section_stl(results_final, \"aorta\")\n", "mm.export_section_stl(results_final, \"lca\")\n", "mm.export_section_stl(results_final, \"rca\")\n", "print(\"Final labeled sections exported.\")\n" ] } ], "metadata": { "kernelspec": { "display_name": ".venv", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.10.12" } }, "nbformat": 4, "nbformat_minor": 5 }