SuperElastix / SimpleElastix

Multi-lingual medical image registration library
http://simpleelastix.github.io
Apache License 2.0
509 stars 149 forks source link

Set interpolator for TransformixImageFilter #409

Closed ctorti closed 3 years ago

ctorti commented 3 years ago

Hi all,

I registered two 3D DICOM images using sitk.ElastixImageFilter() and used the image filter to transform a binary 3D image using sitk.TransformixImageFilter():

TxImFilt = sitk.TransformixImageFilter()

TxParamMap = RegImFilt.GetParameterMap() 

TxImFilt.SetTransformParameterMap(RegImFilt.GetTransformParameterMap())

TxImFilt.SetMovingImage(Im)

TxImFilt.Execute()

TxIm = TxImFilt.GetResultImage()

where RegImFilt is the image filter used to apply the registration, and Im is the binary image I wish to transform.

However TxIm has many (non-binary) values. To put some numbers to this, the input image Im has two unique values (0 and 1) and TxIm has 337, 653 unique values ranging from -0.19 to 1.08 (mean = 0.001). It is my expectation that a pixel is either transformed somewhere or it is not.

Question 1: Should multiple values be expected to arise from the transformation of a binary image?

I wondered if I could set the interpolator that is used so that new "labels" are not introduced (akin to setting the interpolator for ResampleImageFilter when resampling a binary image:

Resampler = sitk.ResampleImageFilter()
Resampler.SetInterpolator(sitk.sitkNearestNeighbor)

So I tried adding the line TxImFilt.SetInterpolator(sitk.sitkNearestNeighbor) to the transformation code above but it didn't accept the command.

Question 2: Is it possible to force TransformixImageFilter to provide a binary output (e.g. by setting an interpolator)?

Question 3: If the answers to 1 and 2 are "yes" and "no", I suppose the only ways around this are to either binarise TxIm (but I'm not sure I can trust TxIm), or rather than transforming the binary image, I could convert it to a list of indices, create a text file with

index {number of indices} index0 index1 index2 ...

use TxImFilt.SetFixedPointSetFileName('filename.txt') to transform the indices, then re-create the 3D labelmap from the list of transformed indices.

It's not elegant but it's the only work around I can think of. Does anyone have any other suggestions?

NHPatterson commented 3 years ago
  1. No. There are some other issues on the GH you can check out where this is addressed.

  2. Yes, see below.

  3. This may be problematic as the points must be specified in fixed image space, which is probably not your goal, you would need to invert the transform to go from moving to fixed. There are a lot of issues on this too.

You can access the transform parameter map directly in python as you would a dict and change the resampling interpolator to what is appropriate. All settings for transformation are set in this parameter map rather than the TransformixImageFilter object.

# n.b. if you used multiple models, you may need to loop over them to set the interpolator for each
# or at the very least for the last one
transform_map = RegImFilt.GetTransformParameterMap()
transform_map["ResampleInterpolator"] = ["FinalNearestNeighborInterpolator"]
TxImFilt.SetTransformParameterMap(transform_map)
# continue as previously

Let us know if that works

ctorti commented 3 years ago

Thanks @NHPatterson. It worked with only a slight modification (#273) to your suggested code:

transform_map[0]["ResampleInterpolator"] = ["FinalNearestNeighborInterpolator"]

Thanks to you I've no need to perform a binary thresholding of the transformed labelmap! 👍

Closing this now.