Notebook

Creating a Printable Model from a 3D Medical Image¶

A Tutorial on dicom2stl.py¶

heads

In [ ]:

import SimpleITK as sitk
%matplotlib notebook

Digital Imaging and Communications in Medicine (DICOM)¶

DICOM is the standard for the communication and management of medical imaging information and related data.

DICOM is most commonly used for storing and transmitting medical images enabling the integration of medical imaging devices such as scanners, servers, workstations, printers, network hardware, and picture archiving and communication systems (PACS) from multiple manufacturers

https://en.wikipedia.org/wiki/DICOM

Imaging Modalities¶

CT (computed tomography)
MRI (magnetic resonance imaging)
ultrasound
X-ray
fluoroscopy
angiography
mammography
breast tomosynthesis
PET (positron emission tomography)
SPECT (single photon emission computed tomography)
Endoscopy
microscopy and whole slide imaging
OCT (optical coherence tomography).

In [ ]:

ct_image = sitk.ReadImage('Data/ct_example.nii.gz')
mri_image = sitk.ReadImage('Data/mri_t1_example.nii.gz')

import gui
gui.MultiImageDisplay(image_list=[ct_image, mri_image], title_list=['CT Head', 'MRI T1 Head'])

CT Houndsfield Units¶

Hounsfield units (HU) are a dimensionless unit universally used in computed tomography (CT) scanning to express CT numbers in a standardized and convenient form. Hounsfield units are obtained from a linear transformation of the measured attenuation coefficients 1

Water is 0 HU
Air is -1000 HU
Very dense bone is 2000 HU
Metal is 3000 HU

Houndsfield Wikipedia page

Image Segmentation¶

The process of partitioning an image into multiple segments.

Typically used to locate objects and boundaries in images.

We use thresholding (selecting a range of image intesities), but SimpleITK has a variety of algorithms

SimpleITK Notebooks

In [ ]:

from myshow import myshow, myshow3d

ct_bone = ct_image>200

# To visualize the labels image in RGB with needs a image with 0-255 range
ct255_image = sitk.Cast(sitk.IntensityWindowing(ct_bone,0,500.0,0.,255.), 
                        sitk.sitkUInt8)

ct255_bone = sitk.Cast(ct_bone, sitk.sitkUInt8)

myshow(sitk.LabelOverlay(ct255_image, ct255_bone), "Basic Thresholding")

Iso-surface extraction¶

Extract a polygonal surface from a 3D image. The most well known algorithm is Marching Cubes (Lorenson & Cline, SIGGRAPH 1987). The 2D version is Marching Squares, shown below

Marching Squares

Marching Cubes¶

And here is the lookup table for Marching Cubes

Marching Cubes

dicom2stl.py processing pipeline¶

SimpleITK image processing pipeline

Shrink the volume to 256^3
Apply anisotripic smoothing
Threshold
- Preset tissue types: skin, bone, fat, soft tissue
- User specified iso-value
Median filter
Pad the volume with black

VTK mesh pipeline

Run Marching Cubes to extract surface
Apply CleanMesh filter to merge vertices
Apply SmoothMesh filter
Run polygon reduction
Write STL

In [ ]:

import itkwidgets
head = sitk.ReadImage("Data/ct_head.nii.gz")
itkwidgets.view(head)

In [ ]:

import sys, os

# download dicom2stl if it's not here already
if not os.path.isdir('dicom2stl'):
    !{'git clone https://github.com/dave3d/dicom2stl.git'}
    
!{sys.executable} dicom2stl/dicom2stl.py -h

In [ ]:

!{sys.executable} dicom2stl/dicom2stl.py -i 400 -o bone.stl Data/ct_head.nii.gz

In [ ]:

from dicom2stl.utils import vtkutils
mesh = vtkutils.readMesh('bone.stl')
itkwidgets.view(head, geometries=[mesh])

In [ ]: