Resampling an Image onto Another's Physical Space

The purpose of this Notebook is to demonstrate how the physical space described by the meta-data is used when resampling onto a reference image.


In [ ]:
from __future__ import print_function
import matplotlib.pyplot as plt
%matplotlib inline

import SimpleITK as sitk
print(sitk.Version())
from myshow import myshow
# Download data to work on
from downloaddata import fetch_data as fdata

OUTPUT_DIR = "Output"

Load the RGB cryosectioning of the Visible Human Male dataset. The data is about 1GB so this may take several seconds, or a bit longer if this is the first time the data is downloaded from the midas repository.


In [ ]:
fixed = sitk.ReadImage(fdata("vm_head_rgb.mha"))

In [ ]:
moving = sitk.ReadImage(fdata("vm_head_mri.mha"))

In [ ]:
print(fixed.GetSize())
print(fixed.GetOrigin())
print(fixed.GetSpacing())
print(fixed.GetDirection())

In [ ]:
print(moving.GetSize())
print(moving.GetOrigin())
print(moving.GetSpacing())
print(moving.GetDirection())

In [ ]:
import sys
resample = sitk.ResampleImageFilter()
resample.SetReferenceImage(fixed)
resample.SetInterpolator(sitk.sitkBSpline)
resample.AddCommand(sitk.sitkProgressEvent, lambda: print("\rProgress: {0:03.1f}%...".format(100*resample.GetProgress()),end=''))
resample.AddCommand(sitk.sitkProgressEvent, lambda: sys.stdout.flush())
out = resample.Execute(moving)

Because we are resampling the moving image using the physical location of the fixed image without any transformation (identity), most of the resulting volume is empty. The image content appears in slice 57 and below.


In [ ]:
myshow(out)

In [ ]:
#combine the two images using a checkerboard pattern:
  #because the moving image is single channel with a high dynamic range we rescale it to [0,255] and repeat 
  #the channel 3 times
vis = sitk.CheckerBoard(fixed,sitk.Compose([sitk.Cast(sitk.RescaleIntensity(out),sitk.sitkUInt8)]*3), checkerPattern=[15,10,1])

In [ ]:
myshow(vis)

Write the image to the Output directory: (1) original as a single image volume and (2) as a series of smaller JPEG images which can be constructed into an animated GIF.


In [ ]:
import os

sitk.WriteImage(vis, os.path.join(OUTPUT_DIR, "example_resample_vis.mha"))

temp = sitk.Shrink(vis,[3,3,2])
sitk.WriteImage(temp, [os.path.join(OUTPUT_DIR,"r{0:03d}.jpg".format(i)) for i in range(temp.GetSize()[2])])

In [ ]: