It might be better to read the data with vtkNIFTIImageReader.
The vtkNIFTIImageReader converts the image volume to RAS when it reads it, because some VTK filters, mappers, and writers don’t work properly with LAS. The vtkNIFTIImageWriter can convert the volume back to LAS when it saves it.
The workflow would be something like the following. First, read the image:
reader = vtk.vtkNIFTIImageReader()
reader.SetFileName(ref_image_volume_path)
reader.Update()
imageMatrix = reader.GetSFormMatrix()
if imageMatrix is None:
imageMatrix = reader.GetQFormMatrix()
The imageMatrix
contains the offset (converted from LAS to RAS) as well as the orientation.
Next, read the meshes and transform them into the same coords as the vtkImageData:
refMeshSurfaceSourcePoints = getSurfaceMeshPoints(ref_volumetric_mesh_path)
deformedMeshSurfaceTargetPoints = getSurfaceMeshPoints(predicted_deformed_mesh_path)
pointsToImage = vtk.vtkTransform()
pointsToImage.Concatenate(imageMatrix)
pointsToImage.Inverse()
newSourcePoints = vtk.vtkPoints()
pointsToImage.TransformPoints(refMeshSurfaceSourcePoints, newSourcePoints)
newTargetPoints = vtk.vtkPoints()
pointsToImage.TransformPoints(deformedMeshSurfaceTargetPoints, newTargetPoints)
Now the newSourcePoints
and newTargetPoints
are in the same coordinate space as the vtkImageData, so they can be used to build the reslice transform.
Finally, you can write the image. Use SetQFac(-1)
to write as LAS instead of RAS:
writer = vtk.vtkNIFTIImageWriter()
writer.SetFileName(deformed_image_volume_path)
writer.SetQFac(reader.GetQFac())
writer.SetQFormMatrix(imageMatrix)
writer.SetSFormMatrix(imageMatrix)
writer.SetInputData(reslice.GetOutput())
writer.Write()
As a sanity check, the Bounds of all the data sets should be similar:
print("input bounds", reader.GetOutput().GetBounds())
print("output bounds", reslice.GetOutput().GetBounds())
print("source bounds", newSourcePoints.GetBounds())
print("target bounds", newTargetPoints.GetBounds())
All of the above code is completely untested, but it gives an overview of how I might transform a nifti file with VTK. This kind of code is very tricky, and I usually use other packages to do things like this, for example the ANTs package has WarpImageMultiTransform for warping NIFTI files.