[Neuroimaging] Slicing images on different planes

[Christopher Markiewicz]

Hi Carl,

It sounds like what you want to do is to interpolate the image in a space that's more conveniently rotated so that image dimensions correspond to meaningful directions in the world and then visualize slices. The way to go about that is going to depend on how well the orientation is described in the image affine matrix.

If the affine correctly describes rotations, then you can simply resample the image with nibabel.processing.resample_to_output() (https://nipy.org/nibabel/reference/nibabel.processing.html#resample-to-output). That will produce an image in-memory that you can work with as you show in your code, or save and open in another viewer.

If your images don't have reliable rotation information stored (e.g., img.affine[:3, :3] contains 6 zeroes), then you're going to need to find a rotation so that you can do something like the above. I have to imagine that tools exist that allow you to place axes in a 3D image and determine the rotation matrix, but I don't know one to suggest.

Best,
Chris

________________________________________
From: Neuroimaging <neuroimaging-bounces+markiewicz=stanford.edu at python.org> on behalf of Buist, Carl R. <cbuist at ou.edu>
Sent: Friday, March 19, 2021 3:30 PM
To: neuroimaging at python.org
Subject: [Neuroimaging] Slicing images on different planes

Hello all,

I would like to preference this question by saying that I am new to the world of python in general and apologize if this is a silly question and if my code is suboptimal. I am also not using medical data, but microCT scanned images of fossils, but I think that the principles for slicing the images are pretty similar.

I have recently written a very basic code to slice a NIfTI (.nii) file along the major axes or some variation of those axes. Now I would like to slice the file along a plane that is not along one of these major axes, for example maybe a plane that is rotated 45 degrees off the x axis in any direction or something like that. My initial thought was to define a plane and use that in place of one of the planes I am currently displaying but am running into trouble.

I tried to look through the archive to see if this was a problem for others but didn't see anything. I would really appreciate any help or if someone could point me in the correct direction. Below is the simple code that I made and attached is an image of the output.

Thank you so much for your help,
Carl

def display_views(file_path):
    image_axis = 2
    medical_image = nib.load(file_path)
    image = medical_image.get_fdata()
    image = np.squeeze(image)
    print(image.shape)

    # Infer dimensions on each Cartesian axis
    x_dim = len(image[:,0,0])
    y_dim = len(image[0,:,0])
    z_dim = len(image[0,0,:])

    # Midpoint layers, as integers
    x_midpoint = int(x_dim/2)
    y_midpoint = int(y_dim/2)
    z_midpoint = int(z_dim/2)

    sagital_image = image[x_midpoint, :, :]
    coronal_image = image[:, y_midpoint, :]
    axial_image = image[:, :, z_midpoint]
    axial_imagep10 = image[:, :, int(z_midpoint + 100)]

    plt.figure(figsize=(20, 10))
    plt.style.use('grayscale')

    plt.subplot(141)
    plt.imshow(np.rot90(sagital_image))
    plt.title('Sagital Plane')
    plt.axis('off')

    plt.subplot(142)
    plt.imshow(np.rot90(coronal_image))
    plt.title('Coronal Plane')
    plt.axis('off')

    plt.subplot(143)
    plt.imshow(np.rot90(axial_image))
    plt.title('Axial Plane')
    plt.axis('off')

    plt.subplot(144)
    plt.imshow(np.rot90(axial_imagep10))
    plt.title('Axial Plane moved')
    plt.axis('off')

    plt.show()

#fpath = sys.argv[1]
fpath = '/content/drive/MyDrive/3D_Forams/Foram_NH1108sit_gdr_7_1__IR_rec480-1180_testcrop.nii'

display_views(fpath)