PythonでVTKファイルの3次元データの断面（スライス）を取得する方法（スカラー値）

pip install vtk

plane = vtk.vtkPlane()
plane.SetOrigin(0, 0, 25)  # 通る点（平面の位置）
plane.SetNormal(0, 1, 1)   # 法線ベクトル（平面の傾き）

cutter = vtk.vtkCutter() # vtkCutterクラスのインスタンスを作成する
cutter.SetCutFunction(plane) # 定義した平面をセットする
cutter.SetInputData(dataset) # カットする元データ（3次元データ）をセットする
cutter.Update()
sliced_output = cutter.GetOutput()

if dataset.IsA("vtkDataSet"):
    print("このデータはvtkCutterでスライス可能です！")
else:
    print("このデータはvtkCutterには使えません。")

import vtk
import numpy as np

# 構造格子のサイズ
nx, ny, nz = 50, 50, 50

# スカラー値（例：中心からの距離の2乗）
scalars = np.zeros((nx, ny, nz))
for k in range(nz):
    for j in range(ny):
        for i in range(nx):
            x = i - nx / 2
            y = j - ny / 2
            z = k - nz / 2
            scalars[i, j, k] = x**2 + y**2 + z**2

# vtkImageDataに変換
image_data = vtk.vtkImageData()
image_data.SetDimensions(nx, ny, nz)
image_data.SetSpacing(1.0, 1.0, 1.0)

# VTK用スカラー配列に変換
flat_scalars = scalars.ravel(order="F")  # VTKはFortran順
vtk_array = vtk.vtkFloatArray()
vtk_array.SetName("Distance2")
vtk_array.SetNumberOfValues(len(flat_scalars))
for i, val in enumerate(flat_scalars):
    vtk_array.SetValue(i, val)

# スカラーをセット
image_data.GetPointData().SetScalars(vtk_array)

# ファイルに保存
writer = vtk.vtkStructuredPointsWriter()
writer.SetFileName("sample.vtk")
writer.SetInputData(image_data)
writer.Write()

import vtk

# VTKファイルのパスを指定
vtk_file_path = "sample.vtk"

# 可視化に使うスカラー値のフィールド名（VTK内に定義された名前）
field_name = "Distance2"

# スカラー値の定義場所：点（point）またはセル（cell）に対する値かを指定
# サンプルファイルでは、点にスカラー値が定義されている
field_position = "point"
# field_position = "cell"

# VTKファイルを読み込むためのリーダーを作成し、ファイルを読み込む
reader = vtk.vtkDataSetReader()
reader.SetFileName(vtk_file_path)
reader.Update()
dataset = reader.GetOutput() # 読み込まれたデータセット（メッシュ＋属性）を取得

# 指定したフィールド名に対応するスカラー配列を取得
if field_position == "cell":
    field = dataset.GetCellData().GetArray(field_name)
elif field_position == "point":
    field = dataset.GetPointData().GetArray(field_name)
else:
    assert 0 # エラーチェック（"point" か "cell" 以外は不正）

# スカラー値の最小値・最大値（カラーマップの範囲）を取得
field_range = field.GetRange()

# スライスする平面を定義：ここではZ=25の位置で、(0,1,1)方向に垂直な面
plane = vtk.vtkPlane()
plane.SetOrigin(0, 0, 25)  # 平面が通る座標
plane.SetNormal(0, 1, 1)  # 平面の傾き（法線ベクトルで指定）

# vtkCutter を使って、定義した平面でデータをスライス
cutter = vtk.vtkCutter()
cutter.SetCutFunction(plane) # 定義した平面をセットする
cutter.SetInputData(dataset) # カットする元データ（3次元データ）をセットする
cutter.Update()
sliced_output = cutter.GetOutput()

# スカラー値に応じて色をつけるマッパー（描画設定）
mapper = vtk.vtkDataSetMapper()
mapper.SetInputData(sliced_output)
mapper.ScalarVisibilityOn()  # スカラー値による色分けを有効に

# スカラー値の位置（点 or セル）に応じてマッピング設定を切り替え
if field_position == "cell":
    mapper.SetScalarModeToUseCellData()
elif field_position == "point":
    mapper.SetScalarModeToUsePointData()
else:
    assert 0
# カラーマップの範囲を指定（最小値～最大値）
mapper.SetScalarRange(field_range)
mapper.SetColorModeToMapScalars() # スカラー値を色にマッピングするように設定

# Actorを設定
actor = vtk.vtkActor()
actor.SetMapper(mapper) # マッパーを関連づける
actor.RotateX(0.0) # 必要に応じて方向調整も可能
actor.RotateY(0.0)
actor.RotateZ(0.0)

# レンダラー（描画担当）を作成し、背景色と表示内容を設定
renderer = vtk.vtkRenderer()
renderer.SetBackground(1, 1, 1) # 背景は白
renderer.AddActor(actor)

# レンダリングウィンドウを作成し、サイズを指定
window = vtk.vtkRenderWindow()
window.AddRenderer(renderer)
window.SetSize(600, 500)

# インタラクション（マウス操作など）を可能にするインタラクターを作成
interactor = vtk.vtkRenderWindowInteractor()
interactor.SetRenderWindow(window)

window.Render() # 初期描画を実行
interactor.Start() # インタラクション（回転・ズームなど）を開始