diff --git a/clinica/utils/exceptions.py b/clinica/utils/exceptions.py index 2dfb13601..bcea3420f 100644 --- a/clinica/utils/exceptions.py +++ b/clinica/utils/exceptions.py @@ -34,6 +34,10 @@ def __init__(self, dataset_folder_path): ) +class ClinicaImageDimensionError(ClinicaException): + """Base class for errors linked to image dimensions.""" + + class ClinicaParserError(ClinicaException): """Base class for parser errors.""" diff --git a/clinica/utils/image.py b/clinica/utils/image.py index 8b04fe97a..c9c41e44c 100644 --- a/clinica/utils/image.py +++ b/clinica/utils/image.py @@ -300,10 +300,26 @@ def __repr__(self): MNI_CROP_BBOX = Bbox3D.from_coordinates(12, 181, 13, 221, 0, 179) +def _is_bbox_within_array(array: np.ndarray, bbox: Bbox3D) -> bool: + if not (0 <= bbox.x_slice.start and bbox.x_slice.end <= array.shape[0]): + return False + if not (0 <= bbox.y_slice.start and bbox.y_slice.end <= array.shape[1]): + return False + if not (0 <= bbox.z_slice.start and bbox.z_slice.end <= array.shape[2]): + return False + return True + + def _crop_array(array: np.ndarray, bbox: Bbox3D) -> np.ndarray: # TODO: When Python 3.10 is dropped, replace with 'return array[*bbox.get_slices()]' - x, y, z = bbox.get_slices() - return array[x, y, z] + from clinica.utils.exceptions import ClinicaImageDimensionError + + if _is_bbox_within_array(array, bbox): + x, y, z = bbox.get_slices() + return array[x, y, z] + raise ClinicaImageDimensionError( + f"Cannot use the bounding box {bbox} to crop the provided array of shape {array.shape}." + ) def _get_file_locally_or_download( @@ -387,7 +403,7 @@ def _get_mni_template_flair() -> Path: ) -def crop_nifti(input_image: Path, output_dir: Optional[Path] = None) -> Path: +def crop_nifti(input_image_path: Path, output_dir: Optional[Path] = None) -> Path: """Crop input image. The function expects a 3D anatomical image and will crop it @@ -397,9 +413,13 @@ def crop_nifti(input_image: Path, output_dir: Optional[Path] = None) -> Path: into the cropped template (located in clinica/resources/masks/ref_cropped_template.nii.gz). + If the bounding box falls outside of the image limits, the + function will perform resampling of the input image onto the + reference cropped template. + Parameters ---------- - input_image : Path + input_image_path : Path The path to the input image to be cropped. output_dir : Path, optional @@ -413,25 +433,42 @@ def crop_nifti(input_image: Path, output_dir: Optional[Path] = None) -> Path: Raises ------ - ValueError: - If the input image is not 3D. + ClinicaImageDimensionError: + If the input image is not 3D or if the output image has unexpected dimension. """ - from nilearn.image import new_img_like + from nilearn.image import new_img_like, resample_to_img + from clinica.utils.exceptions import ClinicaImageDimensionError from clinica.utils.filemanip import get_filename_no_ext + from clinica.utils.stream import log_and_warn - filename_no_ext = get_filename_no_ext(input_image) - input_image = nib.load(input_image) + filename_no_ext = get_filename_no_ext(input_image_path) + input_image = nib.load(input_image_path) + reference_image = nib.load(get_mni_cropped_template()) if len(input_image.shape) != 3: - raise ValueError( + raise ClinicaImageDimensionError( "The function crop_nifti is implemented for anatomical 3D images. " f"You provided an image of shape {input_image.shape}." ) output_dir = output_dir or Path.cwd() - crop_img = new_img_like( - nib.load(get_mni_cropped_template()), - _crop_array(input_image.get_fdata(), MNI_CROP_BBOX), - ) + try: + cropped_array = _crop_array(input_image.get_fdata(), MNI_CROP_BBOX) + crop_img = new_img_like(reference_image, cropped_array) + except ClinicaImageDimensionError: + log_and_warn( + ( + f"The image {input_image_path} has dimensions {input_image.get_fdata().shape} and cannot be " + f"cropped using the bounding box {MNI_CROP_BBOX}. The `crop_nifti` function will try to resample the " + f"input image to the reference template {get_mni_cropped_template()} instead of cropping." + ), + UserWarning, + ) + crop_img = resample_to_img(input_image, reference_image, force_resample=True) + if crop_img.shape != reference_image.shape: + raise ClinicaImageDimensionError( + f"The cropped image has shape {crop_img.shape} different from the expected shape " + f"{reference_image.shape} of the reference template {get_mni_cropped_template()}." + ) output_img = output_dir / f"{filename_no_ext}_cropped.nii.gz" crop_img.to_filename(output_img) diff --git a/test/unittests/utils/test_image.py b/test/unittests/utils/test_image.py index 6877049db..230868bb0 100644 --- a/test/unittests/utils/test_image.py +++ b/test/unittests/utils/test_image.py @@ -5,6 +5,7 @@ import pytest from numpy.testing import assert_array_equal +from clinica.utils.exceptions import ClinicaImageDimensionError from clinica.utils.testing_utils import ( assert_nifti_equal, build_test_image_cubic_object, @@ -192,6 +193,24 @@ def test_slice(): assert s.get_slice() == slice(3, 16) +@pytest.mark.parametrize( + "array,coords,expected", + [ + (np.ones((10, 10, 10)), (3, 6, 0, 10, 6, 6), True), + (np.ones((10, 10, 10)), (3, 11, 0, 10, 6, 6), False), + (np.ones((10, 10, 10)), (3, 6, -1, 10, 6, 6), False), + (np.ones((10, 10, 10)), (1, 1, 1, 1, 1, 1), True), + (np.ones((10, 10, 10)), (11, 11, 11, 11, 11, 11), False), + ], +) +def test_is_bbox_within_array( + array: np.ndarray, coords: tuple[int, int, int, int, int, int], expected: bool +): + from clinica.utils.image import Bbox3D, _is_bbox_within_array # noqa + + assert _is_bbox_within_array(array, Bbox3D.from_coordinates(*coords)) is expected + + def test_mni_cropped_bbox(): from clinica.utils.image import MNI_CROP_BBOX # noqa @@ -253,7 +272,7 @@ def test_get_mni_template(tmp_path, mocker): assert img.shape == expected_shape -def test_crop_nifti_error(tmp_path): +def test_crop_nifti_input_image_not_3d_error(tmp_path): from clinica.utils.image import crop_nifti nib.Nifti1Image(np.random.random((10, 10, 10, 10)), np.eye(4)).to_filename( @@ -261,7 +280,7 @@ def test_crop_nifti_error(tmp_path): ) with pytest.raises( - ValueError, + ClinicaImageDimensionError, match=re.escape( "The function crop_nifti is implemented for anatomical 3D images. " "You provided an image of shape (10, 10, 10, 10)." @@ -270,6 +289,26 @@ def test_crop_nifti_error(tmp_path): crop_nifti(tmp_path / "test.nii.gz") +def test_crop_nifti_with_resampling(tmp_path): + from clinica.utils.image import ( + crop_nifti, + get_mni_cropped_template, + get_mni_template, + ) + + nib.load(get_mni_template("flair")).to_filename(tmp_path / "mni.nii.gz") + with pytest.warns( + UserWarning, + match=re.escape( + f"The image {tmp_path / 'mni.nii.gz'} has dimensions (182, 218, 182) and cannot " + "be cropped using the bounding box ( ( 12, 181 ), ( 13, 221 ), ( 0, 179 ) ). " + "The `crop_nifti` function will try to resample the input image to the reference template" + ), + ): + cropped = crop_nifti(tmp_path / "mni.nii.gz", output_dir=tmp_path) + assert nib.load(cropped).shape == nib.load(get_mni_cropped_template()).shape + + def test_crop_nifti(tmp_path): from clinica.utils.image import ( crop_nifti,