Inverse parametermap calculation using SimpleElastix (SimpleITK extension)

inverse-parammap-simpleelastix.py

# This gist is a portion of a code that doesn't work. Right now, since I
# cant use SimpleElastix easily, I switched to FSL for registration,
# which has a CLI.

import SimpleITK as sitk
import os

# Inverse parametermap calculation. It gives a SIGSEV error. Check updates in
# SimpleElastix and see if it's solved.
# https://github.com/SuperElastix/SimpleElastix/issues/167
# Although I try do it manually, the SIGSEV error appears.

invreg_filter = sitk.ElastixImageFilter()

# For the inverse transform calculation we need to set the original fixed image
# as both fixed and moving image
# (section 6.1.6 of the elastix manual)
original_fixed_image = sitk.ReadImage(fixed_image_path)
invreg_filter.SetFixedImage(original_fixed_image)
invreg_filter.SetMovingImage(original_fixed_image)

# Take the parameter map of the previous registration filter_name
inv_par_map = sitk.ReadParameterFile(PARAMETER_MAP_PATH)

# This metric_str makes elastix find the inverse.
inv_par_map['Metric'] = ['DisplacementMagnitudePenalty']

invreg_filter.SetParameterMap(inv_par_map)
invreg_filter.SetInitialTransformParameterFileName(
    reg_param_file)  # Provide the previous registration as initialization

invreg_filter.Execute()  # Get the inverse transform

inv_par_map = invreg_filter.GetTransformParameterMap(0)
# inverse_reg_matrix = list(inv_par_map['TransformParameters'])  # This is the
# inverse matrix, but with a wrong displacement

# Registration matrix
# reg_matrix = transform_parameter_map['TransformParameters']

# Change the displecement to minus the registration displacement
# inverse_reg_matrix[9:] = list(-1 * np.array(reg_matrix[9:], dtype=np.float))

# Write the inverse matrix to the parameter map.
# inv_par_map['TransformParameters'] = [str(f) for f in list(inverse_reg_matrix)]

# Manually set the basic image parameters
inv_par_map['Size'] = [str(f) for f in list(original_size)]
inv_par_map['Spacing'] = [str(f) for f in list(original_spacing)]
inv_par_map['Origin'] = [str(f) for f in list(original_origin)]
inv_par_map['Direction'] = [str(f) for f in list(original_direction)]
inv_par_map['Direction'] = [str(f) for f in list(original_direction)]
inv_par_map['DefaultPixelValue'] = [str(-1000)]

inv_par_map["InitialTransformParametersFileName"] = [
    "NoInitialTransform"]  # As mentioned in the elastix manual

# Save the inverse transform
reg_param_file = os.path.join(self.save_path,
                              os.path.split(self.f_name)[1] + "_invreg.txt")
reg_filter.WriteParameterFile(inv_par_map, reg_param_file)