VTK 使用 vtkDiscreteMarchingCubes 將二元遮罩影像轉為 Mesh 網格教學與範例

介紹如何在 Python 中使用 VTK 的 vtkDiscreteMarchingCubes 將二元遮罩影像轉為 Mesh 網格，並進行平滑化處理。

vtkDiscreteMarchingCubes 建立網格

VTK 的 vtkDiscreteMarchingCubes 可以用來將離散型的標示影像（label image）或遮罩影像（mask image）轉為 Mesh 網格：

import vtk
import itkwidgets

# 讀取二元遮罩影像
reader = vtk.vtkMetaImageReader()
reader.SetFileName("structure_477.mha")

# 將二元遮罩影像轉為 Mesh 網格
discrete = vtk.vtkDiscreteMarchingCubes()
discrete.SetInputConnection(reader.GetOutputPort())
discrete.GenerateValues(1, 1, 1)

vtkWindowedSincPolyDataFilter 網格平滑處理

由於使用二元遮罩影像透過 vtkDiscreteMarchingCubes 所產生的網格會有鋸齒狀，因此通常會搭配 vtkWindowedSincPolyDataFilter 將網格做平滑處理：

# 將 Mesh 網格做平滑處理
sinc = vtk.vtkWindowedSincPolyDataFilter()
sinc.SetInputConnection(discrete.GetOutputPort())

# 平滑處理疊代次數
sinc.SetNumberOfIterations(50)

# 設定 Pass Band（值越小越平滑）
sinc.SetPassBand(0.003)

# 其他參數設定
#sinc.BoundarySmoothingOn()
#sinc.FeatureEdgeSmoothingOn()
#sinc.SetFeatureAngle(120.0)
#sinc.NonManifoldSmoothingOn()
#sinc.NormalizeCoordinatesOn()

# 降低 Mesh 網格中的三角形數量
decimate = vtk.vtkQuadricDecimation()
decimate.SetInputConnection(sinc.GetOutputPort())
decimate.SetVolumePreservation(True)

# 降低比例
decimate.SetTargetReduction(0.95)

decimate.Update()

最後以 itkwidgets 顯示網格：

# 顯示結果
itkwidgets.view(geometries = decimate.GetOutput())

vtkSmoothPolyDataFilter 網格平滑處理

除了 vtkWindowedSincPolyDataFilter 之外，亦可使用 vtkSmoothPolyDataFilter 來進行網格的平滑處理，其使用拉普拉斯平滑（Laplacian smoothing）演算法來調整每個頂點（vertex）的座標位置，在每一次疊代中，對於每一個頂點 v，先找出所有連接 v 的頂點，將這些頂點座標做平均之後，根據 relaxation factor 所指定的比例來修正 v 的座標。

# 將 Mesh 網格做平滑處理
smoother = vtk.vtkSmoothPolyDataFilter()
smoother.SetInputConnection(discrete.GetOutputPort())

# 平滑處理疊代次數
smoother.SetNumberOfIterations(30)

# 設定 Relaxation Factor
smoother.SetRelaxationFactor(0.6)

# 其他參數設定
#smoother.BoundarySmoothingOn()
#smoother.FeatureEdgeSmoothingOn()
#smoother.SetFeatureAngle(120.0)

# 降低 Mesh 網格中的三角形數量
decimate = vtk.vtkQuadricDecimation()
decimate.SetInputConnection(smoother.GetOutputPort())
decimate.SetVolumePreservation(True)
decimate.SetTargetReduction(0.95)
decimate.Update()

# 顯示結果
itkwidgets.view(geometries = decimate.GetOutput())

在使用 vtkSmoothPolyDataFilter 做平滑處理時，通常以較小的 relaxation factor 搭配較大的疊代次數會比較穩定。

參考資料：VTKExamples：SmoothDiscreteMarchingCubes