Read CT volume compatible to nibabel orientation using SimpleITK

Beginner Questions
,
1

I want to read the CT volume from simpleITK with compatible to nibabel.

When I read the CT volume using nibabel, the orientation (i.e. LAS), image origin, spacing and start, and ending physical points are output as seen below.

img = nibabel.load('./MCR.nii.gz')
print('\nImage shape:', img.shape)

# check orientation of the 3D CT
x, y, z = nibabel.aff2axcodes(img.affine)
print('\nImage orientation:', x, y, z)

voxel_spacing = img.affine[:3, :3]
image_origin = img.affine[:3, 3]

print('\nVoxel spacing:\n', voxel_spacing)

print('\nImage origin:\n', image_origin)

starting_physical_point = voxel_spacing.dot([0, 0, 0]) + image_origin
print('\nstarting_physical_point:', starting_physical_point)

ending_physical_point = voxel_spacing.dot([511, 511, 104]) + image_origin
print('\nending_physical_point:', ending_physical_point)

The output is as below:

Image shape: (512, 512, 105)

Image orientation: L A S

Voxel spacing:
[[-0.9765625 0. 0. ]
[-0. 0.9765625 0. ]
[ 0. -0. 2. ]]

Image origin:
[ 250. -280.9234314 -68. ]

starting_physical_point: [ 250. -280.9234314 -68. ]

ending_physical_point: [-249.0234375 218.1000061 140. ]

However, if I use simpleITK, the output is different as shown below.

read_sitk_mcr_img = sitk.ReadImage('./MCR.nii.gz')
print('\nImage shape:', read_sitk_mcr_img.GetSize())
print('\nVoxel spacing:', read_sitk_mcr_img.GetSpacing())
print('\nImage origin:', read_sitk_mcr_img.GetOrigin())
print('\nDirection:', read_sitk_mcr_img.GetDirection())
print('\nstarting_physical_point:', read_sitk_mcr_img.TransformIndexToPhysicalPoint([0,0,0]))
print('\nending_physical_point:', read_sitk_mcr_img.TransformIndexToPhysicalPoint([sz-1 for sz in read_sitk_mcr_img.GetSize()]))

The output is as below:

Image shape: (512, 512, 105)

Voxel spacing: (0.9765625, 0.9765625, 2.0)

Image origin: (-250.0, 280.9234313964844, -68.0)

Direction: (1.0, 0.0, 0.0, 0.0, -1.0, 0.0, 0.0, 0.0, 1.0)

starting_physical_point: (-250.0, 280.9234313964844, -68.0)

ending_physical_point: (249.0234375, -218.10000610351562, 140.0)

How can I use SimpleITK to convert or read a CT volume in nifti format that is compatible with the nibabel output ? Any help/code would be appreciated.

2

You will probably find this Jupyter notebook helpful:

This was written for ITKPython, not SimpleITK, but it might be directly usable for SimpleITK images too. TorchIO implementation was written for SimpleITK images, so you could use that directly.

1 Like
3

Thank you and it is working. I just modified the code you provided below to generate the updated simpleITK image.

def get_nibabel_compatible_img_from_itk(itk_image) -> np.ndarray:
    spacing = np.array(itk_image.GetSpacing())
    direction_lps = np.array(itk_image.GetDirection())
    origin_lps = np.array(itk_image.GetOrigin())
    rotation_lps = direction_lps.reshape(3, 3)
    
    FLIPXY_33 = np.diag([-1, -1, 1])
    rotation_ras = np.dot(FLIPXY_33, rotation_lps)
    rotation_ras_zoom = rotation_ras * spacing
    translation_ras = np.dot(FLIPXY_33, origin_lps)
    affine = np.eye(4)
    affine[:3, :3] = rotation_ras_zoom
    affine[:3, 3] = translation_ras
    
    voxel_spacing = np.sqrt(np.sum(rotation_ras_zoom * rotation_ras_zoom, axis=0))
    
    image_out = itk_image
    #setup other image characteristics
    image_out.SetOrigin(tuple(translation_ras))
    image_out.SetSpacing(tuple(voxel_spacing))
    #set direction
    image_out.SetDirection(tuple(rotation_ras.ravel()))
    
    return image_out

It now perfectly matches the nibabel output. The result is as follows.

Image shape: (512, 512, 105)

Voxel spacing: (0.9765625, 0.9765625, 2.0)

Image origin: (250.0, -280.9234313964844, -68.0)

Direction: (-1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0)

starting_physical_point: (250.0, -280.9234313964844, -68.0)

ending_physical_point: (-249.0234375, 218.10000610351562, 140.0)

2 Likes
4

@isuruwi I am using this code, it’s changing the results as you mentioned, but when I plot the figure it’s still the same as before applying the function.
can you please help?