How to obtain the complete center lines using VMTK?

Thank you @lassoan! I was not aware of these indications for quantitative analysis of the vessel’s characteristics.

I have been experimenting with several configurations for the resampling of the binary labelmap (I can’t change the resolution of the source image since I am segmenting it using a trained neural net) and I have been able to find configurations that, along with a smoothing filter, work reasonably well (although visually I can’t reach the same quality of the original segmentation).

However, I still run into the same problem of same problem of straight “centerlines” traveling outside the model in some cases (I attach a picture below).

I have tried to snap the endpoints to the closest Voronoi Diagram point as you mention in this issue, but it does not seem to change anything. I attach the code I used for the modification of the endpoints:

# Instance Extract Centerline Widget
extractCenterlineWidget = slicer.modules.extractcenterline.widgetRepresentation().self()

# Set up parameter node
parameterNode = slicer.mrmlScene.GetSingletonNode("ExtractCenterline", "vtkMRMLScriptedModuleNode")
extractCenterlineWidget.setParameterNode(parameterNode)

# Set up widget
extractCenterlineWidget.setup()

# Update from GUI to get segmentationNode as inputSurfaceNode
extractCenterlineWidget.updateParameterNodeFromGUI()
# Set network node reference to new empty node
extractCenterlineWidget._parameterNode.SetNodeReferenceID("InputSurface", segmentationNode.GetID())

# Autodetect endpoints
extractCenterlineWidget.onAutoDetectEndPoints()

extractCenterlineWidget.updateGUIFromParameterNode()

# Create new Surface model node for the centerline model
voronoiDiagramNode = slicer.mrmlScene.AddNewNodeByClass("vtkMRMLModelNode")
# Set centerline node reference to new empty node
extractCenterlineWidget._parameterNode.SetNodeReferenceID("VoronoiDiagram", voronoiDiagramNode.GetID())
extractCenterlineWidget.onApplyButton()

endpointsNode = slicer.util.getNode(extractCenterlineWidget._parameterNode.GetNodeReferenceID("EndPoints"))

voronoiDiagramPointsArray = np.ndarray([voronoiDiagramNode.GetPolyData().GetNumberOfPoints() // 10, 3])
for idx in range(voronoiDiagramNode.GetPolyData().GetNumberOfPoints() // 10):
    voronoiDiagramPointsArray[idx] = voronoiDiagramNode.GetPolyData().GetPoints().GetPoint(idx * 10)

endpointsVoronoiNormArray = np.ndarray([endpointsNode.GetNumberOfMarkups(), voronoiDiagramNode.GetPolyData().GetNumberOfPoints() // 10])

for idx in range(endpointsNode.GetNumberOfMarkups()):
    for idx2 in range(voronoiDiagramNode.GetPolyData().GetNumberOfPoints() // 10):
        endpointsVoronoiNormArray[idx][idx2] = np.linalg.norm(endpointsNode.GetCurvePoints().GetPoint(idx) - voronoiDiagramPointsArray[idx2])
    
endpointsArgVoronoiNormArray = np.argmin(endpointsVoronoiNormArray, axis=1)

for idx in range(endpointsNode.GetNumberOfMarkups()):
    endpointsNode.GetCurvePoints().SetPoint(idx, voronoiDiagramNode.GetPolyData().GetPoints().GetPoint(endpointsArgVoronoiNormArray[idx] * 10))

extractCenterlineWidget.updateGUIFromParameterNode()

# Set network node reference to new empty node
extractCenterlineWidget._parameterNode.SetNodeReferenceID("InputSurface", segmentationNode.GetID())
# extractCenterlineWidget._parameterNode.SetNodeReferenceID("InputSurface", segmentationNode.GetID())
# Create new Surface model node for the centerline model
centerlineModelNode = slicer.mrmlScene.AddNewNodeByClass("vtkMRMLModelNode")
# Set centerline node reference to new empty node
extractCenterlineWidget._parameterNode.SetNodeReferenceID("CenterlineModel", centerlineModelNode.GetID())
extractCenterlineWidget.onApplyButton()

Note that I only check one in every 10 points from the Voronoi diagram to speed up the computation. Do you see why the code would not work? I can provide more context if necessary. This basically assumes that you have a volume which you have segmented (and stored the segmentation in the segmentationNode).

Thank you again for your help.