Coordinate error when outputting VTK model

Development
,
1

I am writing a code to read a nifiti file with python code, convert it to vtk format, and save it again. When I display the output result in a slicer, the location is different. Below is the output result (original in red) and the code. I would appreciate your help.

image
image643×543 48.2 KB

import nibabel as nib
import numpy as np
import vtk
import vmtk.vtkvmtkMiscPython as vtkvmtkMisc
from vtk.util import numpy_support
import matplotlib.pyplot as plt

input_nifti = 'artery1_obj_1.nii.gz'
msk_name='output_test.vtk'

reader = vtk.vtkNIFTIImageReader()
reader.SetFileName(input_nifti)
reader.Update()

cast = vtk.vtkImageCast()
cast.SetInputData(reader.GetOutput())
cast.SetOutputScalarTypeToFloat()
cast.Update()

mc = vtk.vtkMarchingCubes()
mc.SetInputData(cast.GetOutput())
mc.ComputeNormalsOn()
mc.SetValue(0,0.5)
mc.Update()

transform = vtk.vtkTransform()
if reader.GetQFormMatrix():
    transform.Concatenate(reader.GetQFormMatrix())
elif reader.GetSFormMatrix():
    transform.Concatenate(reader.GetSFormMatrix())


tpd = vtk.vtkTransformPolyDataFilter()
tpd.SetInputData(mc.GetOutput())
tpd.SetTransform(transform)
tpd.Update()

writer = vtk.vtkPolyDataWriter()
writer.SetInputData(tpd.GetOutput())
writer.SetFileName(msk_name)
writer.Update()
2

You have saved your mesh to file in RAS coordinate system. This is highly unusual, as DICOM standard uses LPS and all files storing DICOM-derived data is expected to use LPS coordinate system (unless another coordinate system is explicitly specified in the file).

The solution is to add an RAS->LPS transformation (diag([-1,-1,1,1]) matrix) in your polydata transform.

Note that it is not appropriate to use vtkMarchingCubes on a label image. You either need to apply low-pass filter (e.g., Gaussian image smoothing) on the image before passing it to contour filter, or you need to apply low-pass filter (e.g., windowed sinc polydata smoothing) on the contour filter output. Also, using the classic marching cubes algorithm is not the best solution anymore. Flying edges or surface nets provide results at higher quality and/or faster.