Resampling PET, CT and binary mask images to resolution [1,1,1].

Beginner Questions
, ,
1

Hello SimpleITK community,

I am new to this community and python image processing using SimpleITK in general, so i hope that i have asked the question in the correct way, and i apologise if my question is trivial.

As part of my masters project, i have been provided with a dataset of PET/CT images, with corresponding contours delineating prostate lesions. I have successfully imported the PET and CT images into python, and used the following code to resample the PET images to the same spacing and size as the CT images.

image_dirs = [data_path_pet, data_path_ct, data_path_mask]

def read_first_series_in_dir(dir):
    series_reader = sitk.ImageSeriesReader()
    series_IDs = sitk.ImageSeriesReader.GetGDCMSeriesIDs(dir)
    series_file_names = sitk.ImageSeriesReader.GetGDCMSeriesFileNames(dir, series_IDs[0])
    series_reader.SetFileNames(series_file_names)    
    return series_reader.Execute()
    
images = [read_first_series_in_dir(image_dir) for image_dir in image_dirs]

pet_image = images[0]
ct_image = images[1]
mask_image = images[2]

pet_image_resampled = sitk.Resample(pet_image, ct_image)

Following this, i am left with all three images (PET,CT and mask) with spacing and size as follows:

image size: (512, 512, 221)
image spacing: (1.5234375, 1.5234375, 5.0)

This is great, but the ultimate goal is to feed all these images into a Machine Learning model, and so i would like to resample all PET/CT/Mask combinations to an isotropic spacing of 1x1x1. I have attempted to do this with the following code, but i continue getting a type error with respect to the ResampleImageFilter. The code i am using goes as follows:

images_resampled = [pet_image_resampled, ct_image, mask_image]

dimension = images_resampled[1].GetDimension()

reference_physical_size = np.zeros(dimension)

for img in images_resampled:
    reference_physical_size[:] = [(sz-1)*spc if sz*spc>mx  else mx for sz,spc,mx in zip(img.GetSize(), img.GetSpacing(), reference_physical_size)]
    
reference_origin = np.zeros(dimension)
reference_direction = np.identity(dimension).flatten()
reference_size_x = 512
reference_spacing = [reference_physical_size[0]/(reference_size_x-1)]*dimension
reference_size = [int(phys_sz/(spc) + 1) for phys_sz,spc in zip(reference_physical_size, reference_spacing)]

reference_image = sitk.Image(reference_size, images_resampled[0].GetPixelIDValue())
reference_image.SetOrigin(reference_origin)
reference_image.SetSpacing(reference_spacing)
reference_image.SetDirection(reference_direction)

reference_center = np.array(reference_image.TransformContinuousIndexToPhysicalPoint(np.array(reference_image.GetSize())/2.0))

final_data = []

for img in images_resampled:
    # Transform which maps from the reference_image to the current img with the translation mapping the image
    # origins to each other.
    transform = sitk.AffineTransform(dimension)
    transform.SetMatrix(img.GetDirection())
    transform.SetTranslation(np.array(img.GetOrigin()) - reference_origin)
    # Modify the transformation to align the centers of the original and reference image instead of their origins.
    centering_transform = sitk.TranslationTransform(dimension)
    img_center = np.array(img.TransformContinuousIndexToPhysicalPoint(np.array(img.GetSize())/2.0))
    centering_transform.SetOffset(np.array(transform.GetInverse().TransformPoint(img_center) - reference_center))
    centered_transform = sitk.Transform(transform)
    centered_transform.AddTransform(centering_transform)
    
    # Set all the output image parameters
    resampler = sitk.ResampleImageFilter()
    resampler.SetReferenceImage(reference_image)
    resampler.SetDefaultPixelValue(img.GetPixelIDValue)
    resampler.SetTransform(centered_transform)
    resampler.SetSize(reference_image.GetSize())
    resampler.SetOutputSpacing(reference_image.GetSpacing())
    resampler.SetOutputOrigin(reference_image.GetOrigin())
    resampler.SetOutputDirection(reference_image.GetDirection())
    resampler.SetInterpolator(sitk.sitkLinear)
    resampled_img = resampler.Execute(img)
    
    final_data.append(resampled_img)

Can anyone offer any suggestions as to why what i am doing is not working, and/or if there is a different/easier method of resampling these images to [1,1,1]?

2

Hello @JakeKendrick,
Welcome to SimpleITK!

It is a bit hard to identify the issue without having access to your data, and without additional information (e.g. In what line are you getting the type error? For which image.)
Several observations:

  1. You have images of three types CT (pixel type is sitk.sitkInt16), PET (pixel type is likely sitk.sitkFloat32), mask (pixel type is likely sitk.sitkUInt8). There may be a mismatch between the ResampleImageFilter internal information and the image passed to it.
  2. After setting the reference image for the resampler:
resampler.SetReferenceImage(reference_image)

You don’t need to set:

resampler.SetOutputSpacing(reference_image.GetSpacing())
resampler.SetOutputOrigin(reference_image.GetOrigin())
resampler.SetOutputDirection(reference_image.GetDirection())

Possibly the default pixel value is also taken from the reference image, but I don’t remember. Also, no need to recreate the resampler in the loop, create it once and modify in the loop.

1 Like
3

Hi Zivy,

Thank you for the welcome, glad to be aboard :slight_smile:

The comments are much appreciated, i have removed the superfluous code from my loop. I have managed to solve the problem regarding the type error (I was missing a single pair of brackets, which i have now added!). The error was on the following line, and the error has now been corrected.

resampler.SetDefaultPixelValue(img.GetPixelIDValue())

So the code is running nicely now. Thank you for the help!

1 Like