Upload centerline_extraction.py

centerline_extraction.py

@@ -142,7 +142,7 @@ def compute_centerlines_network(
     # Extract main aorta→iliacs centerline via Voronoi method
     try:
         print("Extracting main centerline...")
-        centerline_polydata = extract_centerline_from_surface(surface_polydata_ras, aorta_iliac_endpoints_ras, curve_sampling_distance, smooth_iterations=0)
+        centerline_polydata = extract_centerline_from_surface(surface_polydata_ras, aorta_iliac_endpoints_ras, curve_sampling_distance, smooth_iterations=20)
         print(f"Main centerline: {centerline_polydata.GetNumberOfCells()} cells, {centerline_polydata.GetNumberOfPoints()} points")
     except Exception as e:
         print(f"Main centerline extraction failed: {type(e).__name__}: {e}")
@@ -157,7 +157,7 @@ def compute_centerlines_network(
         except Exception as e:
             print(f"Branch extraction failed: {type(e).__name__}: {e}")
 
-    # Merge main + branch centerlines, and finally smooth all the centerlines
+    # Merge main + branch centerlines
     if branch_centerlines:
         try:
             merged = merge_main_and_branch_centerlines(centerline_polydata, branch_centerlines)
@@ -726,7 +726,7 @@ def smooth_centerline_polydata(centerline_polydata, iterations, relaxation_facto
 # Process all branches, generating their centerlines, computing their new bifurcation points on actual voronoi centerline (not skeleton)
 def process_all_branches(branch_information, main_centerline_polydata, surface_polydata, curve_sampling_distance):
     if not branch_information:
-        return []
+        return [], []
 
     # Convert branch info from LPS to RAS
     branch_info_ras = []
@@ -768,7 +768,8 @@ def extract_branch_centerline(branch_info, main_centerline_polydata, surface_pol
 
     # Extract the branch centerline
     try:
-        branch_centerline = extract_centerline_from_surface(clipped_surface, branch_endpoints, curve_sampling_distance, smooth_iterations=0)
+        branch_centerline = extract_centerline_from_surface(clipped_surface, branch_endpoints, curve_sampling_distance, smooth_iterations=10)
+        branch_centerline = snap_branch_to_contact(branch_centerline, contact_point) 
         return branch_centerline, contact_point
     except Exception as e:
         print(f"\nWarning: branch centerline extraction failed: {e}")
@@ -812,11 +813,25 @@ def clip_surface_around_points(surface_polydata, point1, point2, radius_multipli
     cleaner.Update()
     return cleaner.GetOutput()
 
+# Forces the branch centerline endpoint nearest to contact_point to land at EXACTLY that coordinate, guaranteeing vtkCleanPolyData can merge it with the main centerline
+def snap_branch_to_contact(branch_cl, contact_point):
+    pts = branch_cl.GetPoints()
+    n = pts.GetNumberOfPoints()
+    if n < 2:
+        return branch_cl
+    p_first = pts.GetPoint(0)
+    p_last  = pts.GetPoint(n - 1)
+    d_first = vtk.vtkMath.Distance2BetweenPoints(p_first, contact_point)
+    d_last  = vtk.vtkMath.Distance2BetweenPoints(p_last,  contact_point)
+    snap_id = 0 if d_first < d_last else n - 1
+    pts.SetPoint(snap_id, contact_point[0], contact_point[1], contact_point[2])
+    pts.Modified()
+    return branch_cl
+
 # Merge the main aorta→iliacs centerline with all branch centerlines into one polydata.
 def merge_main_and_branch_centerlines(main_centerline_polydata, branch_centerlines):
     if not branch_centerlines:
-        # Still smooth the main-only result
-        return smooth_centerline_polydata(main_centerline_polydata, iterations=20, relaxation_factor=0.1)
+        return main_centerline_polydata
 
     append_filter = vtk.vtkAppendPolyData()
     append_filter.AddInputData(main_centerline_polydata)
@@ -828,14 +843,11 @@ def merge_main_and_branch_centerlines(main_centerline_polydata, branch_centerlin
     cleaner = vtk.vtkCleanPolyData()
     cleaner.SetInputData(append_filter.GetOutput())
     cleaner.SetToleranceIsAbsolute(True)
-    cleaner.SetAbsoluteTolerance(0.5)  # 0.5mm tolerance to merge junction points from separately-extracted branch centerlines
+    cleaner.SetAbsoluteTolerance(0.1)  # 0.1mm tolerance to merge junction points from separately-extracted branch centerlines
     cleaner.Update()
 
-    # Smooth ONCE on the fully-connected network — endpoints and junction nodes are fixed by BoundarySmoothingOff
     merged = cleaner.GetOutput()
-    smoothed = smooth_centerline_polydata(merged, iterations=20, relaxation_factor=0.1)
-    print(f"Merged+smoothed: {smoothed.GetNumberOfCells()} cells, {smoothed.GetNumberOfPoints()} points")
-    return smoothed
+    return merged
 
 # Main Pipeline for extracting centerline data, labeling bifurcations and arteries, computing ostium and reference point positions
 def extract_centerline_data(centerlines, branch_junctions=None):