InsightSoftwareConsortium / ITKElastix

An ITK Python interface to elastix, a toolbox for rigid and nonrigid registration of images
Apache License 2.0
201 stars 22 forks source link

External initial transform causes a failure in ConvertToItkTransform #245

Open dzenanz opened 1 year ago

dzenanz commented 1 year ago
elastix_object = itk.ElastixRegistrationMethod.New(fixed_image, moving_image)
elastix_object.SetExternalInitialTransform(initial_transform)
...
comb_transform = elastix_object.GetCombinationTransform()
itk_transform = elastix_object.ConvertToItkTransform(comb_transform)  # error

produces:

RuntimeError: D:\a\im\_skbuild\win-amd64-3.9\cmake-build\_deps\elx-src\Core\Main\itkElastixRegistrationMethod.hxx:1019:
ITK ERROR: Failed to convert transform object SimilarityTransformElastix (000001A9A7B0E0B0)
  RTTI typeinfo:   class elastix::SimilarityTransformElastix<class elastix::ElastixTemplate<class itk::Image<float,3>,class itk::Image<float,3> > >
  Reference Count: 12
  Modified Time: 50477
  Debug: Off
  Object Name: 
  Observers: 
    none

The code does not crash if elastix_object.SetInitialTransform(initial_transform) or elastix_object.SetInitialTransformParameterFileName(initial_transform_txt_path) is invoked instead.

dzenanz commented 1 year ago

If this is expected behavior, perhaps documentation could be improved.

What is the difference between SetExternalInitialTransform and SetInitialTransform?

N-Dekker commented 1 year ago

Thanks for reporting this issue, @dzenanz Can you please post the complete code example? Or at least tell me how initial_transform was initialized?

I just tried passing an itk::Similarity2DTransform<double> as ExternalInitialTransform to a 2D registration, and it seems to work fine, without any RuntimeError or ITK ERROR.

dzenanz commented 1 year ago

read_slicer_fiducials is from here: https://github.com/KitwareMedical/HASI/blob/d58acf4cda836bdedd1b7531dddc8ca7b20b396a/src/hasi/mouse_femur_tibia_ct_morphometry.py#L39-L67 Here is how I construct initial transform:

def register_landmarks(fixed_landmarks, moving_landmarks):
    transform_type = itk.Transform[itk.D, 3, 3]
    landmark_transformer = itk.LandmarkBasedTransformInitializer[transform_type].New()
    rigid_transform = itk.Similarity3DTransform[itk.D].New()
    landmark_transformer.SetFixedLandmarks(fixed_landmarks)
    landmark_transformer.SetMovingLandmarks(moving_landmarks)
    landmark_transformer.SetTransform(rigid_transform)
    landmark_transformer.InitializeTransform()
    return rigid_transform

def main():
    print("Fiducial registration")
    fixed_fiducials = read_slicer_fiducials(fixed_pts_fpath)
    moving_fiducials = read_slicer_fiducials(moving_pts_fpath)
    fiducial_transform = register_landmarks(fixed_fiducials, moving_fiducials)
    # write to file to work around this bug
    itk_fiducial_filename = str(output_dir / "fiducial_transform.tfm")
    elx_fiducial_filename = str(output_dir / "fiducial_transform.txt")
    itk.transformwrite([fiducial_transform], itk_fiducial_filename)
    out_elastix_transform = open(elx_fiducial_filename, "w")
    out_elastix_transform.writelines(['(Transform "File")\n',
                                      '(TransformFileName "./fiducial_transform.tfm")\n'])
    out_elastix_transform.close()
    # invoke ITKElastix
N-Dekker commented 1 year ago

perhaps documentation could be improved. What is the difference between SetExternalInitialTransform and SetInitialTransform?

Thanks @dzenanz I made that an issue:

N-Dekker commented 1 year ago

@dzenanz Sorry I still haven't been able to reproduce the error you reported, "Failed to convert transform object SimilarityTransformElastix". Would it be possible for you to post a minimal example to reproduce the issue, that I could simply copy-and-paste into a Jupyter Notebook?

In the mean time, Marius (@mstaring) already pinpointed some limitations of the use of an external initial transform to me. It appears that some very specific metrics components may not support an external initial transform. Specifically "TransformBendingEnergyPenalty" (and maybe some other metrics as well). To be investigated. However, when the specified metric does not allow using an external initial transform, a different error message will appear, saying "Not implemented for AdvancedTransformAdapter" To be continued...

dzenanz commented 1 year ago

Below is the code which triggers the error. To be used from a jupyter notebook or a script in examples directory. I based it on https://github.com/InsightSoftwareConsortium/ITKElastix/blob/4fbdc9c9ca32aa740d8842c0892b817cf8ef0f41/examples/ITK_Example04_InitialTransformAndMultiThreading.ipynb.

---------------------------------------------------------------------------
RuntimeError                              Traceback (most recent call last)
Cell In[1], line 27
     25 result_image = elastix_object.GetOutput()
     26 comb_transform = elastix_object.GetCombinationTransform()
---> 27 itk_transform = elastix_object.ConvertToItkTransform(comb_transform)
     28 elastix_transform_parameters = elastix_object.GetTransformParameterObject()

RuntimeError: D:\a\im\_skbuild\win-amd64-3.9\cmake-build\_deps\elx-src\Core\Main\itkElastixRegistrationMethod.hxx:1019:
ITK ERROR: Failed to convert transform object SimilarityTransformElastix (000001D6E5473250)
  RTTI typeinfo:   class elastix::SimilarityTransformElastix<class elastix::ElastixTemplate<class itk::Image<float,2>,class itk::Image<float,2> > >
  Reference Count: 9
  Modified Time: 36004
  Debug: Off
  Object Name: 
  Observers: 
    none
import itk

# Import Images
fixed_image = itk.imread('data/CT_2D_head_fixed.mha', itk.F)
moving_image = itk.imread('data/CT_2D_head_moving.mha', itk.F)

# Import Default Parameter Map
parameter_object = itk.ParameterObject.New()
parameter_map_rigid = parameter_object.GetDefaultParameterMap('rigid')
parameter_map_rigid['Transform'] = ['SimilarityTransform']
parameter_object.AddParameterMap(parameter_map_rigid)

# Load Elastix Image Filter Object with initial transfrom and number of threads
elastix_object = itk.ElastixRegistrationMethod.New(fixed_image,moving_image)
elastix_object.SetParameterObject(parameter_object)
initial_transform = itk.Similarity2DTransform[itk.D].New()
elastix_object.SetExternalInitialTransform(initial_transform)
# Set additional options
elastix_object.SetLogToConsole(False)

# Update filter object (required)
elastix_object.UpdateLargestPossibleRegion()

# Results of Registration
result_image = elastix_object.GetOutput()
comb_transform = elastix_object.GetCombinationTransform()
itk_transform = elastix_object.ConvertToItkTransform(comb_transform)
elastix_transform_parameters = elastix_object.GetTransformParameterObject()
N-Dekker commented 1 year ago

Thansk @dzenanz Looks like it can be further reduced to:

import itk
fixed_image = itk.imread('data/CT_2D_head_fixed.mha', itk.F)
moving_image = itk.imread('data/CT_2D_head_moving.mha', itk.F)
parameter_object = itk.ParameterObject.New()
parameter_map_rigid = parameter_object.GetDefaultParameterMap('rigid')
parameter_object.AddParameterMap(parameter_map_rigid)
elastix_object = itk.ElastixRegistrationMethod.New(fixed_image,moving_image)
elastix_object.SetParameterObject(parameter_object)
initial_transform = itk.Similarity2DTransform[itk.D].New()
elastix_object.SetExternalInitialTransform(initial_transform)
elastix_object.Update()
comb_transform = elastix_object.GetCombinationTransform()
itk_transform = elastix_object.ConvertToItkTransform(comb_transform)

Right? Which produces:

RuntimeError: D:\a\im_skbuild\win-amd64-3.8\cmake-build_deps\elx-src\Core\Main\itkElastixRegistrationMethod.hxx:1020: ITK ERROR: Failed to convert transform object EulerTransformElastix (00000202831CCF10) RTTI typeinfo: class elastix::EulerTransformElastix<class elastix::ElastixTemplate<class itk::Image<float,2>,class itk::Image<float,2> > >

N-Dekker commented 1 year ago

I just tried to do the very same in C++:

using ImageType = itk::Image<float>;
const auto fixed_image = itk::ReadImage<ImageType>("data/CT_2D_head_fixed.mha");
const auto moving_image = itk::ReadImage<ImageType>("data/CT_2D_head_moving.mha");
const auto parameter_object = elx::ParameterObject::New();
const auto parameter_map_rigid = elx::ParameterObject::GetDefaultParameterMap("rigid");
parameter_object->AddParameterMap(parameter_map_rigid);
const auto elastix_object = itk::ElastixRegistrationMethod<ImageType, ImageType>::New();
elastix_object->SetFixedImage(fixed_image);
elastix_object->SetMovingImage(moving_image);
elastix_object->SetParameterObject(parameter_object);
const auto initial_transform = itk::Similarity2DTransform<double>::New();
elastix_object->SetExternalInitialTransform(initial_transform);
elastix_object->Update();
const auto comb_transform = elastix_object->GetCombinationTransform();
const auto itk_transform = elastix_object->ConvertToItkTransform(*comb_transform);

But in C++, it seems to work fine: no error, no exception! Do you have a clue?

dzenanz commented 1 year ago

Nothing better than "SetExternalInitialTransform might need special care in wrapping", but that is unlikely 😞

N-Dekker commented 1 year ago

I think this issue should have been fixed by:

Could it be that that PR is not yet included with your ITKElastix version?

dzenanz commented 1 year ago

Without that PR, there is an error "there is no SetExternalInitialTransform symbol". To create minimal repro, I did pip install --upgrade itk-elastix and got itk-elastix-0.18.0, the latest one.

N-Dekker commented 1 year ago

PR https://github.com/SuperElastix/elastix/pull/949 was merged to SuperElastix/elastix/main on August 22. While itk-elastix-0.18.0 was released before, on August 17: https://github.com/InsightSoftwareConsortium/ITKElastix/releases/tag/v0.18.0

I think the elastix PR still has to be included with ITKElastix. The next release...!

dzenanz commented 1 year ago

Hopefully you can update elastix via a PR? Bumping the version number in setup.py at the same time would be great. After that is merged I could tag it if you want, as that is easy.

dzenanz commented 11 months ago

Using itk-elastix 0.19.0, fixes the above problem. However, adding the following line to it:

parameter_object.WriteParameterFile(elastix_transform_parameters, "rigid_transform.txt")

causes this error:

---------------------------------------------------------------------------
RuntimeError                              Traceback (most recent call last)
Cell In[2], line 30
     27 itk_transform = elastix_object.ConvertToItkTransform(comb_transform)
     28 elastix_transform_parameters = elastix_object.GetTransformParameterObject()
---> 30 parameter_object.WriteParameterFile(elastix_transform_parameters, "rigid_transform.txt")

RuntimeError: D:\a\im\_skbuild\win-amd64-3.9\cmake-build\_deps\elx-src\Core\Main\elxParameterObject.cxx:314:
ITK ERROR: ParameterObject(000001FF4E3C2DE0): Error writing to disk: The number of parameter maps (2) does not match the number of provided filenames (1). Please call WriteParameterFiles instead, and provide a vector of filenames.

Is different code needed now? Or do we need a fix in elastix and to update version within itk-elastix?

N-Dekker commented 11 months ago

If elastix_transform_parameters contains multiple maps, my initial guess is:

Does that help already? Otherwise I'd have to try it myself.

dzenanz commented 11 months ago

I tried the plural version, and it has its own problem: https://github.com/InsightSoftwareConsortium/ITKElastix/issues/246#issuecomment-1753154752.