Python VTKを使ってSTLをVoxelメッシュに変換する

はじめに

STLのCADデータからを符号付距離関数(Signed Distance Function: SDF)を使ってVoxelメッシュ作成します。使ったライブラリはPythonのVTKです。なお、コードは全面的に以下の記事を参考にしました。

• https://qiita.com/torisan_piyopiyo/items/588e607c578b4b00903a

環境

Windows 10 home

Anaconda(Python 3.7.6)

VTK 8.2

Paraview 5.4.1(STLとvoxelメッシュの表示に利用）

CADデータ（STLデータ）

スタンフォードドラゴン。

ソースコード

import vtk

filename = "./dragon.stl"
#stlファイルを読む
reader = vtk.vtkSTLReader()
reader.SetFileName(filename)
reader.Update()

#Poly dataを取り出す
poly_data = reader.GetOutput()

#メッシュサイズを決める
cell_dims = [50, 50, 50]#x, y, z

#メッシュグリッドを作成する
#x_min:0 x_max:1, y_min:2,y_max:3,z_min:4,z_max:5
bounds = poly_data.GetBounds()
mesh_pitch = [(bounds[1] - bounds[0])/cell_dims[0],
              (bounds[3] - bounds[2])/cell_dims[1],
              (bounds[5] - bounds[4])/cell_dims[2]]
mins = [bounds[0], bounds[2], bounds[4]]

points = vtk.vtkPoints()
for ix in range(cell_dims[0]+1):
    for iy in range(cell_dims[1]+1):
        for iz in range(cell_dims[2]+1):
            x = ix * mesh_pitch[0] + mins[0]
            y = iy * mesh_pitch[1] + mins[1]
            z = iz * mesh_pitch[2] + mins[2]
            points.InsertNextPoint(x,y,z)
structured_base_mesh = vtk.vtkStructuredGrid()
structured_base_mesh.SetExtent(0, cell_dims[0], 0, cell_dims[1], 0, cell_dims[2])
structured_base_mesh.SetPoints(points)

#構造格子データを非構造格子データに変換する
append = vtk.vtkAppendFilter()
append.AddInputData(structured_base_mesh)
append.Update()
base_mesh = append.GetOutput()

#Voxelの中心を座標を求める
cell_centers = vtk.vtkCellCenters()
cell_centers.SetInputData(base_mesh)
cell_centers.Update()

#STLからVoxelメッシュを作成する 
center_points = cell_centers.GetOutput().GetPoints()
cell_list = vtk.vtkIdList()
sdf = vtk.vtkImplicitPolyDataDistance()
sdf.SetInput(poly_data)
distance_sdf = vtk.vtkDoubleArray()
distance_sdf.SetName("sdf")
for idx in range(center_points.GetNumberOfPoints()):
    current_center = center_points.GetPoint(idx)
    distance = sdf.FunctionValue(current_center)
    distance_sdf.InsertNextValue(distance)
    if distance <= 0:
        cell_list.InsertNextId(idx)

base_mesh.GetCellData().SetScalars(distance_sdf)#SDFの値を追加する

extract_cells = vtk.vtkExtractCells()
extract_cells.SetInputData(base_mesh)
extract_cells.SetCellList(cell_list)
extract_cells.Update()

#VoxelデータをVTKファイルで出力する
writer = vtk.vtkXMLDataSetWriter()
writer.SetFileName("closed_voxel.vtu")
writer.SetInputData(extract_cells.GetOutput())
writer.Update()

#SDFの値を入力した格子データをVTKファイルで出力する
#writer = vtk.vtkXMLDataSetWriter()
#writer.SetFileName("all_voxel.vtu")
#writer.SetInputData(base_mesh)
#writer.Update()