Skip to content

Instantly share code, notes, and snippets.

@fepegar
Last active April 14, 2024 07:31
Show Gist options
  • Save fepegar/a8814ff9695c5acd8dda5cf414ad64ee to your computer and use it in GitHub Desktop.
Save fepegar/a8814ff9695c5acd8dda5cf414ad64ee to your computer and use it in GitHub Desktop.
Python function to create maximum intensity projections using SimpleITK

Requirements:

  • Python
  • SimpleITK
  • NumPy

Install

conda install numpy -y
conda install -c simpleitk simpleitk -y
#!/usr/bin/env python
# -*- coding: utf-8 -*-
"""
Created on Fri Nov 10 16:51:47 2017
@author: fernando
"""
import os
import numpy as np
import SimpleITK as sitk
def make_mips(image_path, output_dir):
image = sitk.ReadImage(image_path)
image_size = image.GetSize()
basename = os.path.basename(image_path)
if not os.path.isdir(output_dir):
os.makedirs(output_dir)
for dim in range(3):
projection = sitk.MaximumProjection(image, dim)
if image_size[dim] % 2: # odd number
voxel = [0, 0, 0]
voxel[dim] = (image_size[dim] - 1) / 2
origin = image.TransformIndexToPhysicalPoint(voxel)
else: # even
voxel1 = np.array([0, 0, 0], int)
voxel2 = np.array([0, 0, 0], int)
voxel1[dim] = image_size[dim] / 2 - 1
voxel2[dim] = image_size[dim] / 2
point1 = np.array(image.TransformIndexToPhysicalPoint(voxel1.tolist()))
point2 = np.array(image.TransformIndexToPhysicalPoint(voxel2.tolist()))
origin = np.mean(np.vstack((point1, point2)), 0)
projection.SetOrigin(origin)
projection.SetDirection(image.GetDirection())
proj_basename = basename.replace('.nii.gz', '_mip_{}.nii.gz'.format(dim))
sitk.WriteImage(projection, os.path.join(output_dir, proj_basename))
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment