How to crop a NIFTI image while keeping its origin?

I have a 512x512x460 NIFTI image. I am trying to crop the image to 512x216x460 while retaining its origin.

Here is my code trying to crop the image

import nibabel as nib
import numpy as np

def cropNifti(img_path, out_path, crop_width=256, crop_height=512):
        
    # Load the NIfTI file
    img = nib.load(img_path)

    # Get the image data
    data = img.get_fdata()

    # Calculate the dimensions for the central crop
    original_shape = data.shape
    central_x = original_shape[0] // 2
    central_y = original_shape[1] // 2
    start_x = max(0, central_x - crop_width // 2)
    end_x = min(original_shape[0], central_x + crop_width // 2)
    start_y = max(0, central_y - crop_height // 2)
    end_y = min(original_shape[1], central_y + crop_height // 2)

    # Calculate the centering offsets
    offset_x = (crop_width - (end_x - start_x)) // 2
    offset_y = (crop_height - (end_y - start_y)) // 2

    # Crop each slice to the central region and center it
    cropped_data = np.zeros((crop_height, crop_width, original_shape[2]), dtype=data.dtype)
    for i in range(original_shape[2]):
        cropped_data[offset_y:(end_y - start_y + offset_y), offset_x:(end_x - start_x + offset_x), i] = data[start_y:end_y, start_x:end_x, i]

    # Save the cropped and centered image
    cropped_img = nib.Nifti1Image(cropped_data, img.affine, img.header)
    nib.save(cropped_img, out_path)    

img_path = r"some\path"
out_path = r"some\other\path"

cropNifti(img_path, out_path)

When I view both the original image vs the cropped image, I expected the cropped image to overlap on the center of the original image. However, the result is that the cropped image is on the side of the original image (the left edge of the cropped image is aligned to the center of the original image). I have tried modifying the code a bit and added img.affine and img.header in the parameters in the hope that it would keep its origin information but that doesn't seem to have worked.

I have also tried cropping the image using SimpleITK. Although the result of the crop is centered, it seem to have rescaled the image or voxel sizes so that I can no longer run the result on one of my deep learning pipeline.

I am relatively new to image processing, thank you for any comments and help.

Solution

The image affine transforms the IJK indices into RAS coordinates via the formula RAS = img.affine @ IJK. When you remove rows, columns or slices, then what used to be IJK of [0, 0, 0, 1] becomes [di, dj, dk, 1]. To preserve the RAS coordinates of that first voxel (and thus all voxels), we need some transform RAS = img.affine @ xfm @ IJK' = img.affine @ IJK + img.affine @ [di, dj, dk, 1].

This transform is just:

[[1, 0, 0, di],
 [0, 1, 0, dj],
 [0, 0, 1, dk],
 [0, 0, 0,  1]]

So you can update your affine with:

cropped_affine = nib.affines.from_matvec(np.eye(3), [start_x, start_y, 0])

This affine update is handled by nibabel images through the img.slicer property:

import nibabel as nib

img = nib.load(img_path)
new_img = img.slicer[i_start:i_end, j_start:j_end, k_start:k_end]

So in your case:

import nibabel as nib

def cropNifti(img_path, out_path, crop_width=256, crop_height=512):
    img = nib.load(img_path)

    i_start = max(0, (crop_height - img.shape[0]) // 2)
    j_start = max(0, (crop_width - img.shape[1]) // 2)

    cropped_img = img.slicer[i_start:i_start + crop_height, j_start:j_start + crop_width]
    nib.save(cropped_img, out_path)