SuperElastix / elastix

Official elastix repository
http://elastix.dev
Apache License 2.0
483 stars 118 forks source link

Inverting transformationparameters with `DisplacementMagnitudePenalty` does not seem to work #500

Open jakob1379 opened 3 years ago

jakob1379 commented 3 years ago

I have made successfully made a transform from CBCT to CT with Elastix and now I want to invert this. I have read the manual and followed suggestions in issue #149 which made sure I did not reach inf. I am using the exact same parameter files but have changed the metrid to be DisplacementMagnitudePenalty instead of the original mutual information.

The situation is depicted in this image image The yellow body is the gray body registered to the planning-ct I have delineations of tumours and organs which I need the inverse transform to transfer them to CBCT-space. The brown patches are the original delineations from in ct space transformed with what should have been the inverse transform, but it seems like the affine part is missing as the location is still in the image and not on the CBCT. The call to find the inverse transform is the following:

elastix \
    -f cbct1.nii.gz \
    -m cbct1.nii.gz \
    -p "Par0003.affine.inverted.txt" \
    -p "Par0003.bs-R6-ug.inverted.txt" \
    -t0 cbct1_271138_to_planning_ct/TransformParameters.1.txt \
    -threads 4 \
    -out planning_ct_to_cbct1_271138

Where Par0003.affine.inverted.txt is:

//ImageTypes
(FixedInternalImagePixelType "float")
(FixedImageDimension 3)
(MovingInternalImagePixelType "float")
(MovingImageDimension 3)
(UseBinaryFormatForTransformationParameters "false")

//Components
(Registration "MultiResolutionRegistration")
(FixedImagePyramid "FixedSmoothingImagePyramid")
(MovingImagePyramid "MovingSmoothingImagePyramid")
(Interpolator "BSplineInterpolator")
(Metric "DisplacementMagnitudePenalty")
(Optimizer "AdaptiveStochasticGradientDescent")
(ResampleInterpolator "FinalBSplineInterpolator")
(Transform "AffineTransform")
(UseDirectionCosines "false")
// ********** Pyramid

// Total number of resolutions
(NumberOfResolutions 6)

// ********** Transform

(AutomaticTransformInitialization "true")
(AutomaticTransformInitializationMethod "CenterOfGravity")
(AutomaticScalesEstimation "true")
(HowToCombineTransforms "Compose")

// ********** Optimizer

// Maximum number of iterations in each resolution level:
(MaximumNumberOfIterations 1000)

//SP: Param_a in each resolution level. a_k = a/(A+k+1)^alpha
(SP_a 500.0)

//SP: Param_alpha in each resolution level. a_k = a/(A+k+1)^alpha
(SP_alpha 0.602)

//SP: Param_A in each resolution level. a_k = a/(A+k+1)^alpha
(SP_A 50.0)

// ********** Metric

//Number of grey level bins in each resolution level:
(NumberOfHistogramBins 32)
(FixedLimitRangeRatio 0.0)
(MovingLimitRangeRatio 0.0)
(FixedKernelBSplineOrder 1)
(MovingKernelBSplineOrder 3)

// ********** Several

(WriteTransformParametersEachIteration "false")
(WriteTransformParametersEachResolution "false")
(WriteResultImage "true")
(ShowExactMetricValue "false")
(ErodeFixedMask "false")
(ErodeMovingMask "false")
(UseDifferentiableOverlap "false")

// ********** ImageSampler

//Number of spatial samples used to compute the mutual information in each resolution level:
(ImageSampler "Random")
(NumberOfSpatialSamples 2048)
(NewSamplesEveryIteration "true")

// ********** Interpolator and Resampler

//Order of B-Spline interpolation used in each resolution level:
(BSplineInterpolationOrder 1)

//Order of B-Spline interpolation used for applying the final deformation:
(FinalBSplineInterpolationOrder 3)

//Default pixel value for pixels that come from outside the picture:
(DefaultPixelValue -1000)

// Resampler specific
(Resampler "DefaultResampler")
(ResultImageFormat "nii")
(ResultImagePixelType "float")
(CompressResultImage "true")

and Par0003.bs-R6-ug.inverted.txt is:

// ********** Image Types

(FixedInternalImagePixelType "float")
(FixedImageDimension 3)
(MovingInternalImagePixelType "float")
(MovingImageDimension 3)

// ********** Components

(Registration "MultiResolutionRegistration")
(FixedImagePyramid "FixedSmoothingImagePyramid")
(MovingImagePyramid "MovingSmoothingImagePyramid")
(Interpolator "BSplineInterpolator")
(Metric "DisplacementMagnitudePenalty")
(Optimizer "AdaptiveStochasticGradientDescent")
(ResampleInterpolator "FinalBSplineInterpolator")
(Resampler "DefaultResampler")
(Transform "BSplineTransform")

// ********** Pyramid

// Total number of resolutions
(NumberOfResolutions 6)
// default schedule: isotropic upsampling with factor 2

// ********** Transform

(FinalGridSpacingInPhysicalUnits 12.0 12.0 12.0)
(GridSpacingSchedule 16.0 16.0 8.0 4.0 2.0 1.0)
(HowToCombineTransforms "Compose")

// ********** Optimizer

// Maximum number of iterations in each resolution level:
(MaximumNumberOfIterations 1000)

//SP: Param_a in each resolution level. a_k = a/(A+k+1)^alpha
(SP_a 1000.0)

//SP: Param_alpha in each resolution level. a_k = a/(A+k+1)^alpha
(SP_alpha 0.602)

//SP: Param_A in each resolution level. a_k = a/(A+k+1)^alpha
(SP_A 50.0)

// ********** Metric

//Number of grey level bins in each resolution level:
(NumberOfHistogramBins 32)
(FixedLimitRangeRatio 0.0)
(MovingLimitRangeRatio 0.0)
(FixedKernelBSplineOrder 1)
(MovingKernelBSplineOrder 3)

// ********** Several

(WriteTransformParametersEachIteration "false")
(WriteTransformParametersEachResolution "true")
(WriteResultImageAfterEachResolution "false")
(WriteResultImage "true")
(ShowExactMetricValue "false")
(ErodeMask "false")

// ********** ImageSampler

//Number of spatial samples used to compute the mutual information in each resolution level:
(ImageSampler "Random")
(NumberOfSpatialSamples 2000)
(NewSamplesEveryIteration "true")
(UseRandomSampleRegion "false")

// ********** Interpolator and Resampler

//Order of B-Spline interpolation used in each resolution level:
(BSplineInterpolationOrder 1)

//Order of B-Spline interpolation used for applying the final deformation:
(FinalBSplineInterpolationOrder 3)

//Default pixel value for pixels that come from outside the picture:
(DefaultPixelValue 0)

// Resampler specific
(ResultImageFormat "nii")
(ResultImagePixelType "float")
(CompressResultImage "true")

and the t0 is the transform used to transform the gray image to the yellow image which is aligning with the planning ct almost perfectly.

I am running elastix version: 5.000 on manjaro

N-Dekker commented 1 year ago

@ntatsisk Please have a look