Simple ITK Problem when saving deformation field : the deformation field is empty

Question

I'm looking for help with simpleITK, I'm working with NII files in order to calculate the deformation field. I use Python and simpleITK (simpleitk 1.2.0rc2.dev1167+gd4cf2). I can easily register my images or volumes but I cannot get the deformation field!

When I print it in the console simpleITK tells me that the vector is an image as predicted but when I try to save it with "WriteImage" which is a sitk member function, the written image is full of null values

the code used to do that:

def RegistrationGel(pathImage1,pathImage2):

elastixImageFilter = sitk.ElastixImageFilter()
elastixImageFilter.SetFixedImage(sitk.ReadImage(pathImage1))
elastixImageFilter.SetMovingImage(sitk.ReadImage(pathImage2))

parameterMapVector = sitk.VectorOfParameterMap()
parameterMapVector.append(sitk.GetDefaultParameterMap("affine"))
#parameterMapVector.append(sitk.GetDefaultParameterMap("bspline"))
elastixImageFilter.SetParameterMap(parameterMapVector)
elastixImageFilter.Execute()
sitk.WriteImage(elastixImageFilter.GetResultImage(),"C:/Users/Armand/Desktop/RecalageGel/Output/resultImageGel.nii")

transformixImageFilter = sitk.TransformixImageFilter()
transformixImageFilter.SetMovingImage(sitk.ReadImage(pathImage2))
transformixImageFilter.ComputeDeformationFieldOn()
#transformixImageFilter.SetOutputDirectory("C:/Users/Armand/Desktop/RecalageGel/Output/DeformationFIELD.nii")
transformixImageFilter.SetTransformParameterMap(elastixImageFilter.GetTransformParameterMap())
transformixImageFilter.Execute()
sitk.WriteImage(transformixImageFilter.GetDeformationField(),"C:/Users/Armand/Desktop/RecalageGel/Output/OutputDef.nii")

I also tried to set the output directory of transformixImageFilter but without success too. I really do not understand what's happening here... Does anyone have a solution or the same issue? I have read some related topics but without any success

Answer 1

I had the same problem but strangely now it seemed to worked.. here is my code:

moving = sitk.ReadImage("path_to_my_moving_image")
fixed = sitk.ReadImage("path_to_my_fixed_image")

elastixImageFilter = sitk.ElastixImageFilter()
parameterMapVector = sitk.VectorOfParameterMap()
    
elastixImageFilter.SetFixedImage(fixed)
elastixImageFilter.SetMovingImage(moving)
parameterMapVector = sitk.VectorOfParameterMap()
parameterMapVector.append(sitk.GetDefaultParameterMap("rigid"))
parameterMapVector.append(sitk.GetDefaultParameterMap("affine"))
parameterMapVector.append(sitk.GetDefaultParameterMap("bspline"))
elastixImageFilter.SetParameterMap(parameterMapVector)
elastixImageFilter.Execute()

   
transformixImageFilter = sitk.TransformixImageFilter()
transformixImageFilter.SetMovingImage(moving)  
    transformixImageFilter.SetTransformParameterMap(elastixImageFilter.GetTransformParameterMap())
transformixImageFilter.ComputeDeformationFieldOn()
transformixImageFilter.LogToConsoleOn()
transformixImageFilter.SetOutputDirectory("path_to_my_outputfolder/")
transformixImageFilter.Execute()
sitk.WriteImage(transformixImageFilter.GetDeformationField(), "path_to_my_outputfolder/transformix_rigid_affine_bspline.nii.gz")

the output in the terminal for the transformix part looked like this:

ELASTIX version: 5.000
Command line options from ElastixBase:
-out      path_to_my_outputfolder/
-threads  unspecified, so all available threads are used
-def      all
-jac      unspecified, so no det(dT/dx) computed
-jacmat   unspecified, so no dT/dx computed

Reading input image ...
  Reading input image took 0.000005 s
Calling all ReadFromFile()'s ...
WARNING: The parameter "UseBinaryFormatForTransformationParameters", requested at entry number 0, does not exist at all.
  The default value "false" is used instead.
WARNING: The parameter "UseBinaryFormatForTransformationParameters", requested at entry number 0, does not exist at all.
  The default value "false" is used instead.
WARNING: The parameter "UseBinaryFormatForTransformationParameters", requested at entry number 0, does not exist at all.
  The default value "false" is used instead.
  Calling all ReadFromFile()'s took 0.041007 s
Transforming points ...
  The transform is evaluated on all points. The result is a deformation field.
  Transforming points done, it took 8.95s
Compute determinant of spatial Jacobian ...
  The command-line option "-jac" is not used, so no det(dT/dx) computed.
  Computing determinant of spatial Jacobian done, it took 0.00s
Compute spatial Jacobian (full matrix) ...
  The command-line option "-jacmat" is not used, so no dT/dx computed.
  Computing spatial Jacobian done, it took 0.00s
Resampling image and writing to disk ...
  Resampling took 25.14s