Medical Imaging Interaction Toolkit  2018.4.99-a3d2e8fb
Medical Imaging Interaction Toolkit
mitkShapeBasedInterpolationAlgorithm.cpp
Go to the documentation of this file.
1 /*============================================================================
2 
3 The Medical Imaging Interaction Toolkit (MITK)
4 
5 Copyright (c) German Cancer Research Center (DKFZ)
6 All rights reserved.
7 
8 Use of this source code is governed by a 3-clause BSD license that can be
9 found in the LICENSE file.
10 
11 ============================================================================*/
12 
14 #include "mitkImageAccessByItk.h"
15 #include "mitkImageCast.h"
16 #include <mitkITKImageImport.h>
17 
18 #include <itkFastChamferDistanceImageFilter.h>
19 #include <itkInvertIntensityImageFilter.h>
20 #include <itkIsoContourDistanceImageFilter.h>
21 #include <itkSubtractImageFilter.h>
22 
24  Image::ConstPointer lowerSlice,
25  unsigned int lowerSliceIndex,
26  Image::ConstPointer upperSlice,
27  unsigned int upperSliceIndex,
28  unsigned int requestedIndex,
29  unsigned int /*sliceDimension*/, // commented variables are not used
30  Image::Pointer resultImage,
31  unsigned int /*timeStep*/,
32  Image::ConstPointer /*referenceImage*/)
33 {
34  mitk::Image::Pointer lowerDistanceImage = mitk::Image::New();
35  AccessFixedDimensionByItk_1(lowerSlice, ComputeDistanceMap, 2, lowerDistanceImage);
36 
37  mitk::Image::Pointer upperDistanceImage = mitk::Image::New();
38  AccessFixedDimensionByItk_1(upperSlice, ComputeDistanceMap, 2, upperDistanceImage);
39 
40  // calculate where the current slice is in comparison to the lower and upper neighboring slices
41  float ratio = (float)(requestedIndex - lowerSliceIndex) / (float)(upperSliceIndex - lowerSliceIndex);
43  resultImage, InterpolateIntermediateSlice, 2, upperDistanceImage, lowerDistanceImage, ratio);
44 
45  return resultImage;
46 }
47 
48 template <typename TPixel, unsigned int VImageDimension>
49 void mitk::ShapeBasedInterpolationAlgorithm::ComputeDistanceMap(const itk::Image<TPixel, VImageDimension> *binaryImage,
50  mitk::Image::Pointer &result)
51 {
52  typedef itk::Image<TPixel, VImageDimension> DistanceFilterInputImageType;
53 
54  typedef itk::FastChamferDistanceImageFilter<DistanceFilterImageType, DistanceFilterImageType> DistanceFilterType;
55  typedef itk::IsoContourDistanceImageFilter<DistanceFilterInputImageType, DistanceFilterImageType> IsoContourType;
56  typedef itk::InvertIntensityImageFilter<DistanceFilterInputImageType> InvertIntensityImageFilterType;
57  typedef itk::SubtractImageFilter<DistanceFilterImageType, DistanceFilterImageType> SubtractImageFilterType;
58 
59  typename DistanceFilterType::Pointer distanceFilter = DistanceFilterType::New();
60  typename DistanceFilterType::Pointer distanceFilterInverted = DistanceFilterType::New();
61  typename IsoContourType::Pointer isoContourFilter = IsoContourType::New();
62  typename IsoContourType::Pointer isoContourFilterInverted = IsoContourType::New();
63  typename InvertIntensityImageFilterType::Pointer invertFilter = InvertIntensityImageFilterType::New();
64  typename SubtractImageFilterType::Pointer subtractImageFilter = SubtractImageFilterType::New();
65 
66  // arbitrary maximum distance
67  int maximumDistance = 100;
68 
69  // this assumes the image contains only 1 and 0
70  invertFilter->SetInput(binaryImage);
71  invertFilter->SetMaximum(1);
72 
73  // do the processing on the image and the inverted image to get inside and outside distance
74  isoContourFilter->SetInput(binaryImage);
75  isoContourFilter->SetFarValue(maximumDistance + 1);
76  isoContourFilter->SetLevelSetValue(0);
77 
78  isoContourFilterInverted->SetInput(invertFilter->GetOutput());
79  isoContourFilterInverted->SetFarValue(maximumDistance + 1);
80  isoContourFilterInverted->SetLevelSetValue(0);
81 
82  distanceFilter->SetInput(isoContourFilter->GetOutput());
83  distanceFilter->SetMaximumDistance(maximumDistance);
84 
85  distanceFilterInverted->SetInput(isoContourFilterInverted->GetOutput());
86  distanceFilterInverted->SetMaximumDistance(maximumDistance);
87 
88  // inside distance should be negative, outside distance positive
89  subtractImageFilter->SetInput2(distanceFilter->GetOutput());
90  subtractImageFilter->SetInput1(distanceFilterInverted->GetOutput());
91  subtractImageFilter->Update();
92 
93  result = mitk::GrabItkImageMemory(subtractImageFilter->GetOutput());
94 }
95 
96 template <typename TPixel, unsigned int VImageDimension>
97 void mitk::ShapeBasedInterpolationAlgorithm::InterpolateIntermediateSlice(itk::Image<TPixel, VImageDimension> *result,
98  const mitk::Image::Pointer &lower,
99  const mitk::Image::Pointer &upper,
100  float ratio)
101 {
102  typename DistanceFilterImageType::Pointer lowerITK = DistanceFilterImageType::New();
103  typename DistanceFilterImageType::Pointer upperITK = DistanceFilterImageType::New();
104 
105  CastToItkImage(lower, lowerITK);
106  CastToItkImage(upper, upperITK);
107 
108  itk::ImageRegionConstIteratorWithIndex<DistanceFilterImageType> lowerIter(lowerITK,
109  lowerITK->GetLargestPossibleRegion());
110 
111  lowerIter.GoToBegin();
112 
113  if (!lowerITK->GetLargestPossibleRegion().IsInside(upperITK->GetLargestPossibleRegion()) ||
114  !lowerITK->GetLargestPossibleRegion().IsInside(result->GetLargestPossibleRegion()))
115  {
116  // TODO Exception etc.
117  MITK_ERROR << "The regions of the slices for the 2D interpolation are not equally sized!";
118  return;
119  }
120 
121  float weight[2] = {1.0f - ratio, ratio};
122 
123  while (!lowerIter.IsAtEnd())
124  {
125  typename DistanceFilterImageType::PixelType lowerPixelVal = lowerIter.Get();
126  typename DistanceFilterImageType::PixelType upperPixelVal = upperITK->GetPixel(lowerIter.GetIndex());
127 
128  typename DistanceFilterImageType::PixelType intermediatePixelVal =
129  (weight[0] * upperPixelVal + weight[1] * lowerPixelVal > 0 ? 0 : 1);
130 
131  result->SetPixel(lowerIter.GetIndex(), static_cast<TPixel>(intermediatePixelVal));
132 
133  ++lowerIter;
134  }
135 }
#define MITK_ERROR
Definition: mitkLogMacros.h:20
#define AccessFixedDimensionByItk_3(mitkImage, itkImageTypeFunction, dimension, arg1, arg2, arg3)
Image::Pointer GrabItkImageMemory(itk::SmartPointer< ItkOutputImageType > &itkimage, mitk::Image *mitkImage=nullptr, const BaseGeometry *geometry=nullptr, bool update=true)
Grabs the memory of an itk::Image (with a specific type) and puts it into an mitk::Image.The memory is managed by the mitk::Image after calling this function. The itk::Image remains valid until the mitk::Image decides to free the memory.
Image::Pointer Interpolate(Image::ConstPointer lowerSlice, unsigned int lowerSliceIndex, Image::ConstPointer upperSlice, unsigned int upperSliceIndex, unsigned int requestedIndex, unsigned int sliceDimension, Image::Pointer resultImage, unsigned int timeStep, Image::ConstPointer referenceImage) override
static Pointer New()
#define AccessFixedDimensionByItk_1(mitkImage, itkImageTypeFunction, dimension, arg1)
void MITKCORE_EXPORT CastToItkImage(const mitk::Image *mitkImage, itk::SmartPointer< ItkOutputImageType > &itkOutputImage)
Cast an mitk::Image to an itk::Image with a specific type.