SuperElastix / SimpleElastix

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

Transformed labelmap image is cropped #419

Open ctorti opened 3 years ago

ctorti commented 3 years ago

Hello,

I'm getting what appears to be a cropped binary image following use of Transformix. I would appreciate some insight as to whether the effect is some artefact due to improper use, an error in my interpretation, or other.

I begin with a DICOM SEG that relates to a "Source" 3D DICOM image, SrcIm. The objective is to resample the Source 3D segment to a "Target" 3D DICOM image's (TrgIm) grid. The SEG is parsed to obtain a Numpy pixel array, SrcPixArr and based on the ReferencedSOPInstanceUIDs, a list of indices, SrcF2Sinds (frame-to-slice indices), that link each frame in SrcPixArr to the slice number within SrcIm.

SrcPixArr has shape (F, R, C), where F = the number of frames that make up the segment, R/C = the number of rows/cols of SrcIm.

Next I convert SrcPixArr to what I call a "labelmap array", which is simply a zero-padded representation of SrcPixArr. SrcLabmapArr has shape (S, R, C), where S = the number of slices in SrcIm. The indices SrcF2Sinds determine where the F 2D masks are mapped to within the S frames in SrcLabmapArr.

Now that I have a pixel array that matches the dimensions of SrcIm, I convert SrcLabmapArr to a SimpleITK image. I call this a "labelmap image" - not to be confused with ITK's Labelmap class. SrcLabmapIm not only has the same size, but also origin, spacing and direction as SrcIm.

I register SrcIm (= moving image) toTrgIm (= fixed image), then use Transformix to transform SrcLabmapIm to get TxSrcLabmapIm.

I can provide the code used to perform the above, but in the interest of keeping this readable I have only copied the bits that I suspect are relevant to ascertain what's the cause of the cropping.

Following is how I perform the registration:

RegImFilt = sitk.ElastixImageFilter()
RegImFilt.SetFixedImage(TrgIm)
RegImFilt.SetMovingImage(SrcIm)

ParamMap = sitk.GetDefaultParameterMap('affine')

ParamMap['AutomaticTransformInitialization'] = ['true']
ParamMap['AutomaticTransformInitializationMethod'] = ['GeometricalCenter']
ParamMap['AutomaticScalesEstimation'] = ['true']
ParamMap['BSplineInterpolatorOrder'] = ['1']
ParamMap['FixedInternalImagePixelType'] = ['float']
ParamMap['MovingInternalImagePixelType'] = ['float']
ParamMap['HowToCombineTransforms'] = ['Compose']
ParamMap['MaximumNumberOfIterations'] = ['500']
ParamMap['NumberOfHistogramBins'] = ['32']
ParamMap['UseDirectionCosines'] = ['true']
ParamMap['WriteIterationInfo'] = ['true']

RegImFilt.SetParameterMap(ParamMap)

RegImFilt.Execute()

And the transformation of SrcLabmapIm:

TxImFilt = sitk.TransformixImageFilter()
TxImFilt.SetMovingImage(SrcLabmapIm)

TxParamMap = RegImFilt.GetParameterMap()

TxMap = RegImFilt.GetTransformParameterMap()

""" See [1]: """
TxMap[0]["ResampleInterpolator"] = ["FinalNearestNeighborInterpolator"]

""" The image is binary so set the FinalBSplineInterpolationOrder to 0 to avoid "overshoot"-related 
"garbage" as described in the elastix manual (see 4.3 The transform parameter file): """
TxMap[0]["FinalBSplineInterpolationOrder"] = ["0"]

TxImFilt.SetTransformParameterMap(TxMap)

TxImFilt.Execute()

TxSrcLabmapIm = TxImFilt.GetResultImage()

[1] As per @NHPatterson 's advice #409

I've uploaded a jpg of a plot of the segmentations present for Source (spanning slices 108 to 122), and the transformed segmentations for Target (spanning slices 19 to 25).

Test_SR10_Ventricles_3_6_rerun

If it's of any use, here I have put some NIFTI files of SrcIm, TrgIm, SrcLabmapIm and TxSrcLabmapIm (ignore the "New" part of the filename), as well as the source SEG, and folders containing the Source and Target DICOMs. Thanks!

print(SrcIm):

Image (00000246F9D33D60) RTTI typeinfo: class itk::Image<float,3> Reference Count: 2 Modified Time: 965 Debug: Off Object Name: Observers: none Source: (none) Source output name: (none) Release Data: Off Data Released: False Global Release Data: Off PipelineMTime: 943 UpdateMTime: 964 RealTimeStamp: 0 seconds LargestPossibleRegion: Dimension: 3 Index: [0, 0, 0] Size: [512, 512, 192] BufferedRegion: Dimension: 3 Index: [0, 0, 0] Size: [512, 512, 192] RequestedRegion: Dimension: 3 Index: [0, 0, 0] Size: [512, 512, 192] Spacing: [0.5, 0.5, 1] Origin: [-129.808, -109.976, -117.14] Direction: 1 0 0 0 1 0 0 0 1

IndexToPointMatrix: 0.5 0 0 0 0.5 0 0 0 1

PointToIndexMatrix: 2 0 0 0 2 0 0 0 1

Inverse Direction: 1 0 0 0 1 0 0 0 1

PixelContainer: ImportImageContainer (00000246FA83D070) RTTI typeinfo: class itk::ImportImageContainer<unsigned __int64,float> Reference Count: 1 Modified Time: 961 Debug: Off Object Name: Observers: none Pointer: 000002468000A040 Container manages memory: true Size: 50331648 Capacity: 50331648

print(TrgIm):

Image (00000246FB0768B0) RTTI typeinfo: class itk::Image<float,3> Reference Count: 2 Modified Time: 1156 Debug: Off Object Name: Observers: none Source: (none) Source output name: (none) Release Data: Off Data Released: False Global Release Data: Off PipelineMTime: 1133 UpdateMTime: 1155 RealTimeStamp: 0 seconds LargestPossibleRegion: Dimension: 3 Index: [0, 0, 0] Size: [270, 320, 40] BufferedRegion: Dimension: 3 Index: [0, 0, 0] Size: [270, 320, 40] RequestedRegion: Dimension: 3 Index: [0, 0, 0] Size: [270, 320, 40] Spacing: [0.71875, 0.71875, 4] Origin: [-63.5501, -156.401, -57.1759] Direction: 0.972062 -0.19324 -0.133243 0.217328 0.955414 0.199881 0.0886773 -0.223254 0.970718

IndexToPointMatrix: 0.69867 -0.138891 -0.532972 0.156205 0.686704 0.799524 0.0637368 -0.160464 3.88287

PointToIndexMatrix: 1.35243 0.30237 0.123377 -0.268855 1.32927 -0.310615 -0.0333107 0.0499702 0.242679

Inverse Direction: 0.972062 0.217328 0.0886773 -0.19324 0.955414 -0.223254 -0.133243 0.199881 0.970718

PixelContainer: ImportImageContainer (00000246FAE4AFC0) RTTI typeinfo: class itk::ImportImageContainer<unsigned __int64,float> Reference Count: 1 Modified Time: 1152 Debug: Off Object Name: Observers: none Pointer: 000002468C011040 Container manages memory: true Size: 3456000 Capacity: 3456000

print(SrcLabmapIm):

Image (00000246FB076E30) RTTI typeinfo: class itk::Image<unsigned int,3> Reference Count: 1 Modified Time: 1351 Debug: Off Object Name: Observers: none Source: (none) Source output name: (none) Release Data: Off Data Released: False Global Release Data: Off PipelineMTime: 0 UpdateMTime: 0 RealTimeStamp: 0 seconds LargestPossibleRegion: Dimension: 3 Index: [0, 0, 0] Size: [512, 512, 192] BufferedRegion: Dimension: 3 Index: [0, 0, 0] Size: [512, 512, 192] RequestedRegion: Dimension: 3 Index: [0, 0, 0] Size: [512, 512, 192] Spacing: [0.5, 0.5, 1] Origin: [-129.808, -109.976, -117.14] Direction: 1 0 0 0 1 0 0 0 1

IndexToPointMatrix: 0.5 0 0 0 0.5 0 0 0 1

PointToIndexMatrix: 2 0 0 0 2 0 0 0 1

Inverse Direction: 1 0 0 0 1 0 0 0 1

PixelContainer: ImportImageContainer (00000246FAE48D30) RTTI typeinfo: class itk::ImportImageContainer<unsigned __int64,unsigned int> Reference Count: 1 Modified Time: 1348 Debug: Off Object Name: Observers: none Pointer: 00000246BABB8040 Container manages memory: true Size: 50331648 Capacity: 50331648

print(TxSrcLabmapIm):

Image (00000246FB0789B0) RTTI typeinfo: class itk::Image<float,3> Reference Count: 1 Modified Time: 64477 Debug: Off Object Name: Observers: none Source: (none) Source output name: (none) Release Data: Off Data Released: False Global Release Data: Off PipelineMTime: 0 UpdateMTime: 0 RealTimeStamp: 0 seconds LargestPossibleRegion: Dimension: 3 Index: [0, 0, 0] Size: [270, 320, 40] BufferedRegion: Dimension: 3 Index: [0, 0, 0] Size: [270, 320, 40] RequestedRegion: Dimension: 3 Index: [0, 0, 0] Size: [270, 320, 40] Spacing: [0.71875, 0.71875, 4] Origin: [-63.5501, -156.401, -57.1759] Direction: 0.972062 -0.19324 -0.133243 0.217328 0.955414 0.199881 0.0886773 -0.223254 0.970718

IndexToPointMatrix: 0.69867 -0.138891 -0.532972 0.156204 0.686704 0.799524 0.0637368 -0.160464 3.88287

PointToIndexMatrix: 1.35243 0.30237 0.123377 -0.268855 1.32927 -0.310615 -0.0333107 0.0499702 0.24268

Inverse Direction: 0.972062 0.217329 0.0886773 -0.193239 0.955414 -0.223254 -0.133243 0.199881 0.970719

PixelContainer: ImportImageContainer (00000246F9CA4DD0) RTTI typeinfo: class itk::ImportImageContainer<unsigned __int64,float> Reference Count: 1 Modified Time: 64478 Debug: Off Object Name: Observers: none Pointer: 0000024694D6D040 Container manages memory: true Size: 3456000 Capacity: 3456000