Skip to content

Instantly share code, notes, and snippets.

@pangyuteng
Last active April 27, 2024 16:49
Show Gist options
  • Save pangyuteng/0493816e0398cdb59fde48a946745b31 to your computer and use it in GitHub Desktop.
Save pangyuteng/0493816e0398cdb59fde48a946745b31 to your computer and use it in GitHub Desktop.
conversion from stl to nifti and vice versa
  • for conversion of nifti to vtk stl, see nifti2stl.py below

  • for conversion from stl to nifgi, see stl2nifti.py below


wget https://discourse.itk.org/uploads/default/original/1X/ced4efe8c1eb78342afe474f53ad4409ad5dfc0c.obj
mv ced4efe8c1eb78342afe474f53ad4409ad5dfc0c.obj bunny.obj

wget https://upload.wikimedia.org/wikipedia/commons/4/43/Stanford_Bunny.stl


docker build -t stl2nifti .
docker run -it -u $(id -u):$(id -g) -v $PWD:$PWD -w $PWD stl2nifti bash
python stl2nifti.py bunny.obj bunny.nii.gz

python stl2nifti.py Stanford_Bunny.stl Stanford_Bunny.nii.gz


https://discourse.itk.org/t/trianglemeshtobinaryimagefilter-in-python/1604/7

https://examples.itk.org/src/core/mesh/converttrianglemeshtobinaryimage/documentation

https://github.com/InsightSoftwareConsortium/ITKIOMeshSTL
https://stackoverflow.com/questions/42281881/get-stl-format-3d-mesh-from-binary-mask-segmentation

https://discourse.slicer.org/t/itk-image-to-mesh-to-vtk-polydata/12717

https://slicer.readthedocs.io/en/latest/developer_guide/script_repository.html#rasterize-a-model-and-save-it-to-a-series-of-image-files

https://github.com/InsightSoftwareConsortium/ITK/issues/2416

https://stackoverflow.com/questions/31430171/convert-stl-into-numpy-array

https://github.com/InsightSoftwareConsortium/ITK/issues/2884

https://stackoverflow.com/questions/65112987/convert-stl-object-to-vtk-geometry-in-python/68584066#68584066

DICOM RT-STRUCT

https://zhangresearch.org/post/parse-dicom-rtstruct-into-binary-masks
https://github.com/qurit/rt-utils


FROM python:3.8-bullseye
COPY requirements.txt /tmp/requirements.txt
RUN pip install -r /tmp/requirements.txt
'''
http://vtk.1045678.n5.nabble.com/How-to-write-vtkDelaunay3D-into-vtk-file-td5741818.html
https://gist.github.com/pangyuteng/facd430d0d9761fc67fff4ff2e5fffc3
'''
import sys
import itk
import numpy as np
import vtk
nifti_file = sys.argv[1]
output_stl_file = sys.argv[2]
reader = vtk.vtkNIFTIImageReader()
reader.SetFileName(nifti_file)
reader.Update()
threshold = vtk.vtkImageThreshold ()
threshold.SetInputConnection(reader.GetOutputPort())
threshold.ThresholdByLower(2000) #th
threshold.ReplaceInOn()
threshold.SetInValue(0) # set all values below th to 0
threshold.ReplaceOutOn()
threshold.SetOutValue(1) # set all values above th to 1
threshold.Update()
dmc = vtk.vtkDiscreteMarchingCubes()
dmc.SetInputConnection(threshold.GetOutputPort())
dmc.GenerateValues(1, 1, 1)
dmc.Update()
mapper = vtk.vtkPolyDataMapper()
mapper.SetInputConnection(dmc.GetOutputPort())
mapper.Update()
smooth = vtk.vtkSmoothPolyDataFilter()
smooth.SetInputConnection(dmc.GetOutputPort())
smooth.SetNumberOfIterations(100)
smooth.SetRelaxationFactor(1)
smooth.Update()
if output_stl_file.endswith('.vtk'):
writer = vtk.vtkPolyDataWriter()
writer.SetInputData(smooth.GetOutput())
writer.SetFileName('bunny.vtk')
writer.Write()
if output_stl_file.endswith('.ply'):
writer = vtk.vtkPLYWriter()
writer.SetInputConnection(dmc.GetOutputPort())
writer.SetFileName('bunny.ply')
writer.Write()
if output_stl_file.endswith('.stl'):
writer = vtk.vtkSTLWriter()
writer.SetInputConnection(dmc.GetOutputPort())
writer.SetFileTypeToBinary()
writer.SetFileName("bunny.stl")
writer.Write()
itk==5.3.0
itk-iomeshstl==2.3.0
import sys
import itk
import numpy as np
mesh_file = sys.argv[1] # input stl or obj filepath
nifti_file = sys.argv[2] # output nifti filepath
Dimension = 3
PixelType = itk.SS
ImageType = itk.Image[PixelType, Dimension]
meshType = itk.Mesh[itk.SS, 3]
meshReader = itk.MeshFileReader[meshType].New()
if mesh_file.endswith('.obj'):
meshIO = itk.OBJMeshIO.New()
elif mesh_file.endswith('.stl'):
meshIO = itk.STLMeshIO.New()
else:
raise NotImplementedError()
meshReader.SetMeshIO(meshIO)
meshReader.SetFileName(mesh_file)
meshReader.Update()
print(meshReader.GetOutput().GetBoundingBox())
origin = meshReader.GetOutput().GetBoundingBox().GetMinimum()
mesh_size = np.array(meshReader.GetOutput().GetBoundingBox().GetMaximum()) -\
np.array(meshReader.GetOutput().GetBoundingBox().GetMinimum())
image_size = [512,512,512]
sp = np.max(mesh_size/np.array(image_size))
spacing = [sp]*3
print('mesh_origin',origin)
print('mesh_size',mesh_size)
print('image_size',image_size)
print('spacing',spacing)
region = itk.ImageRegion[Dimension]()
region.SetSize(image_size)
region.SetIndex([0,0,0])
image = ImageType.New()
image.SetRegions(region)
image.Allocate()
image.SetOrigin(origin)
image.SetSpacing(spacing)
myfilter = itk.TriangleMeshToBinaryImageFilter.New()
myfilter.SetInput(meshReader.GetOutput())
myfilter.SetInfoImage(image)
myfilter.SetInsideValue(1)
myfilter.Update()
print(myfilter.GetSize())
print(myfilter.GetSpacing())
print(myfilter.GetOrigin())
print(myfilter.GetDirection())
itk.imwrite(myfilter.GetOutput(),nifti_file)
@pangyuteng
Copy link
Author

pangyuteng commented Jul 14, 2023

rendering of input bunny.obj (right), and output bunny.nii.gz (left).

image

@pangyuteng
Copy link
Author

@pangyuteng
Copy link
Author

smoothing_iterations = 15
pass_band = 0.001
feature_angle = 120.0
smoother = vtk.vtkWindowedSincPolyDataFilter()
smoother.SetInputConnection(dmc.GetOutputPort())
smoother.SetNumberOfIterations(smoothing_iterations)
smoother.BoundarySmoothingOff()
smoother.FeatureEdgeSmoothingOff()
smoother.SetFeatureAngle(feature_angle)
smoother.SetPassBand(pass_band)
smoother.NonManifoldSmoothingOn()
smoother.NormalizeCoordinatesOn()
smoother.Update()

@pangyuteng
Copy link
Author

# https://examples.vtk.org/site/Python/Meshes/Decimation
reduction = 0.9
decimate = vtkDecimatePro()
decimate.SetInputData(dmc.GetOutputPort())
decimate.SetTargetReduction(reduction)
decimate.PreserveTopologyOn()
decimate.Update()
decimated = vtkPolyData()
decimated.ShallowCopy(decimate.GetOutput())

@pangyuteng
Copy link
Author

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment