Convert .nii.gz warp to .h5 (with proper header)

Hi!

I have been modifying an anatomical processing workflow that largely uses ANTs. I have been trying to incorporate FreeSurfer’s synthmorph to calculate a registration, which outputs in RAS convention. I have been able to convert this deformation to a LPS (like ITK/ANTs expects) NIFTI file fine with FreeSurfer’s mri_warp_convert.

However, for the final step of the pipeline when a set of registrations are combined and applied, I would need a .h5 file with the proper formatting/header information. FreeSurfer doesn’t know what kind of formatting is expected, so simply specifying the output of .h5 in mri_warp_convert doesn’t work, as I am met with the following error when trying to apply it with antsApplyRegistration:

Description: ITK ERROR: TransformFileReaderTemplate(0x55565df6abb0): Could not create Transform IO object for reading file /home/smeisler/projects/fmriprep_ants/single_subject_bids/derivatives/synthmorph_test/sub-NDARINV1Y40DZT8/ses-baselineYear1Arm1/anat/sub-NDARINV1Y40DZT8_ses-baselineYear1Arm1_rec-normalized_from-T1w_to-MNI152NLin2009cAsym_mode-image_xfm.h5
	  Tried to create one of the following:
	    HDF5TransformIOTemplate
	    HDF5TransformIOTemplate
	    MINCTransformIOTemplate
	    MINCTransformIOTemplate
	    MatlabTransformIOTemplate
	    MatlabTransformIOTemplate
	    TxtTransformIOTemplate
	    TxtTransformIOTemplate
	  You probably failed to set a file suffix, or
	    set the suffix to an unsupported type.

I am hoping Convert3D will help, but I am not sure exactly how I would implement this. Any thoughts would be appreciated.

Thank you,
Steven

As an update, I have tried the following to manually convert,

   def convert_nifti_to_ants_h5_deformation(self, nifti_file, h5_file):
        # Load the NIfTI file
        img = nib.load(nifti_file)
        data = img.get_fdata()

        # Ensure the data is in float32 format for compatibility
        if data.dtype != np.float32:
            data = data.astype(np.float32)

        # Create the HDF5 file with ANTs specific structure
        with h5py.File(h5_file, 'w') as f:
            # Create the Transform group
            xform_group = f.create_group('/TransformGroup/0')
            
            # Write out the transformation parameters
            xform_group.create_dataset('TransformParameters', data=data.flatten(order='C'))
            
            # Write additional attributes and datasets expected by ANTs
            xform_group.attrs['TransformType'] = 'DisplacementFieldTransform'
            xform_group.attrs['TransformDisplacementField'] = 'DisplacementField'
            xform_group.attrs['TransformFixedParameters'] = img.affine.flatten().tolist()
            
            f.attrs['ITKTransformType'] = 'DisplacementFieldTransform_double_3_3'
            f.attrs['ITKVersion'] = '4.11.0'
            f.attrs['HDFVersion'] = '1.8.15'

and now ANTs returns the following errors when I go to apply the registration:

Stderr:
	HDF5-DIAG: Error detected in HDF5 (1.14.3) thread 1:
	  #000: H5D.c line 403 in H5Dopen2(): unable to synchronously open dataset
	    major: Dataset
	    minor: Can't open object
	  #001: H5D.c line 364 in H5D__open_api_common(): unable to open dataset
	    major: Dataset
	    minor: Can't open object
	  #002: H5VLcallback.c line 1980 in H5VL_dataset_open(): dataset open failed
	    major: Virtual Object Layer
	    minor: Can't open object
	  #003: H5VLcallback.c line 1947 in H5VL__dataset_open(): dataset open failed
	    major: Virtual Object Layer
	    minor: Can't open object
	  #004: H5VLnative_dataset.c line 331 in H5VL__native_dataset_open(): unable to open dataset
	    major: Dataset
	    minor: Can't open object
	  #005: H5Dint.c line 1418 in H5D__open_name(): not found
	    major: Dataset
	    minor: Object not found
	  #006: H5Gloc.c line 421 in H5G_loc_find(): can't find object
	    major: Symbol table
	    minor: Object not found
	  #007: H5Gtraverse.c line 816 in H5G_traverse(): internal path traversal failed
	    major: Symbol table
	    minor: Object not found
	  #008: H5Gtraverse.c line 596 in H5G__traverse_real(): traversal operator failed
	    major: Symbol table
	    minor: Callback failed
	  #009: H5Gloc.c line 381 in H5G__loc_find_cb(): object 'TransformType' doesn't exist
	    major: Symbol table
	    minor: Object not found
	Transform reader for /home/smeisler/projects/fmriprep_ants/single_subject_bids/derivatives/synthmorph_test/sub-NDARINV1Y40DZT8/ses-baselineYear1Arm1/anat/sub-NDARINV1Y40DZT8_ses-baselineYear1Arm1_rec-normalized_from-T1w_to-MNI152NLin2009cAsym_mode-image_xfm.h5 caught an ITK exception:

	itk::ExceptionObject (0x560251e0ed50)
	Location: "unknown" 
	File: /home/conda/feedstock_root/build_artifacts/libitk_1717078286681/work/Modules/IO/TransformHDF5/src/itkHDF5TransformIO.cxx
	Line: 365
	Description: ITK ERROR: HDF5TransformIOTemplate(0x560251d81890): H5Dopen2 failed

Is this manual setting a good direction to go down or are there more automated tools for this?

Thanks,
Steven

Have you tried putting your warp field into a DisplacementFieldTransform? And then save that via itk.transformwrite()?

2 Likes

Hi @dzenanz, thanks for responding!

Just to make sure, do you mean something like:

import SimpleITK as sitk
import nibabel as nib

# Load nifti warp
img = nib.load(nifti_file)
data = img.get_fdata()

# Define and save ITK warp
displacement_image = sitk.GetImageFromArray(data, isVector=True)
tx = sitk.DisplacementFieldTransform(displacement_image)
sitk.WriteTransform(tx, "out_warp.h5");

Best,
Steven

That is close, but you’ll need to use SimpleITK/ITK to read the NIFTI deformation field to preserve all the metadata.

You can see an example here using the sitk.DisplacementFieldTransform with a displacement field image:
https://simpleitk.readthedocs.io/en/master/link_DemonsRegistration1_docs.html

1 Like

Great, I’ll give this a shot and report back, thanks @dzenanz and @blowekamp !

Hi @blowekamp and @dzenanz,

I tried the following:

import SimpleITK as sitk
in_nii="path_to_LPS_warp.nii.gz"
out_h5="path_to_LPS_warp.h5"
displacement_image = sitk.ReadImage(in_nii, sitk.sitkFloat64, imageIO="NiftiImageIO")
tx = sitk.DisplacementFieldTransform(displacement_image)
sitk.WriteTransform(tx,out_h5)

and am met with the following error:

File "/opt/conda/envs/fmriprep/lib/python3.11/site-packages/niworkflows/interfaces/itk.py", line 346, in _run_interface
  tx = sitk.DisplacementFieldTransform(displacement_image)
       ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
File "/opt/conda/envs/fmriprep/lib/python3.11/site-packages/SimpleITK/SimpleITK.py", line 5608, in __init__
  _SimpleITK.DisplacementFieldTransform_swiginit(self, _SimpleITK.new_DisplacementFieldTransform(*args))
                                                       ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
RuntimeError: Exception thrown in SimpleITK new_DisplacementFieldTransform: /tmp/SimpleITK/Code/Common/src/sitkDisplacementFieldTransform.cxx:66:
sitk::ERROR: Expected input displacement field image must be of pixel type: vector of 64-bit float

I thought that hard-coding that precision with sitk.sitkFloat64 would get past this error. Do you have any suggestions about how I may proceed?

Thanks,
Steven

The above line converts the vector pixels to scalars, the correct type would be sitk.sitkVectorFloat64. However the above specificity is not needed. Why not just do:

displacement_image = sitk.ReadImage(in_nii)

P.S. If the error continues do a print(displacement_image.ToString()) to get the details about the image that was read.

1 Like

Thanks @blowekamp,

I got the same error when just doing displacement_image = sitk.ReadImage(in_nii), but seem to have gotten past it when hard coding the correct sitk.sitkVectorFloat64 precision!

Waiting until my pipeline finishes to see if registration looks like it should, but thank you so much for helping me get to this point!

Best,
Steven

1 Like

The displacement field may have been save as 32-bit float ( the print statement should show), white the ITK Transform requires 64-bits.

1 Like