Medical Imaging Interaction Toolkit  2018.4.99-12ad79a3
Medical Imaging Interaction Toolkit
mitkMRNormTwoRegionBasedFilter.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 
15 #include "mitkImageToItk.h"
16 #include "mitkImageAccessByItk.h"
17 
18 #include "itkImageRegionIterator.h"
19 // MITK
20 #include <mitkITKImageImport.h>
21 #include <mitkImageCast.h>
22 #include <mitkImageAccessByItk.h>
23 // ITK
24 #include <itkLabelStatisticsImageFilter.h>
25 #include <itkMinimumMaximumImageCalculator.h>
26 
28  m_Area1(MRNormTwoRegionsBasedFilter::MEDIAN),
29  m_Area2(MRNormTwoRegionsBasedFilter::MEDIAN)
30 {
31  this->SetNumberOfIndexedInputs(3);
32  this->SetNumberOfRequiredInputs(1);
33 }
34 
36 {
37 }
38 
40 {
41  // Process object is not const-correct so the const_cast is required here
42  auto* nonconstMask = const_cast< mitk::Image * >( mask );
43  this->SetNthInput(1, nonconstMask );
44 }
45 
47 {
48  // Process object is not const-correct so the const_cast is required here
49  auto* nonconstMask = const_cast< mitk::Image * >( mask );
50  this->SetNthInput(2, nonconstMask );
51 }
52 
54 {
55  return this->GetInput(1);
56 }
57 
59 {
60  return this->GetInput(2);
61 }
62 
64 {
65  Superclass::GenerateInputRequestedRegion();
66 
67  mitk::Image* input = this->GetInput();
68 
70 }
71 
73 {
74  mitk::Image::ConstPointer input = this->GetInput();
75  mitk::Image::Pointer output = this->GetOutput();
76 
77  itkDebugMacro(<<"GenerateOutputInformation()");
78 
79  output->Initialize(input->GetPixelType(), *input->GetTimeGeometry());
80  output->SetPropertyList(input->GetPropertyList()->Clone());
81 }
82 
83 template < typename TPixel, unsigned int VImageDimension >
84 void mitk::MRNormTwoRegionsBasedFilter::InternalComputeMask(itk::Image<TPixel, VImageDimension>* itkImage)
85 {
86  // Define all necessary Types
87  typedef itk::Image<TPixel, VImageDimension> ImageType;
88  typedef itk::Image<int, VImageDimension> MaskType;
89  typedef itk::LabelStatisticsImageFilter<ImageType, MaskType> FilterType;
90  typedef itk::MinimumMaximumImageCalculator<ImageType> MinMaxComputerType;
91 
92  typename MaskType::Pointer itkMask0 = MaskType::New();
93  typename MaskType::Pointer itkMask1 = MaskType::New();
94  mitk::CastToItkImage(this->GetMask1(), itkMask0);
95  mitk::CastToItkImage(this->GetMask2(), itkMask1);
96 
97  typename ImageType::Pointer outImage = ImageType::New();
98  mitk::CastToItkImage(this->GetOutput(0), outImage);
99 
100  typename MinMaxComputerType::Pointer minMaxComputer = MinMaxComputerType::New();
101  minMaxComputer->SetImage(itkImage);
102  minMaxComputer->Compute();
103 
104  typename FilterType::Pointer labelStatisticsImageFilter = FilterType::New();
105  labelStatisticsImageFilter->SetUseHistograms(true);
106  labelStatisticsImageFilter->SetHistogramParameters(256, minMaxComputer->GetMinimum(),minMaxComputer->GetMaximum());
107  labelStatisticsImageFilter->SetInput( itkImage );
108 
109  labelStatisticsImageFilter->SetLabelInput(itkMask0);
110  labelStatisticsImageFilter->Update();
111  double median0 = labelStatisticsImageFilter->GetMedian(1);
112  double mean0 = labelStatisticsImageFilter->GetMean(1);
113  double modulo0=0;
114  {
115  auto histo = labelStatisticsImageFilter->GetHistogram(1);
116  double maxFrequency=0;
117  for (auto hIter=histo->Begin();hIter!=histo->End();++hIter)
118  {
119  if (maxFrequency < hIter.GetFrequency())
120  {
121  maxFrequency = hIter.GetFrequency();
122  modulo0 = (histo->GetBinMin(0,hIter.GetInstanceIdentifier()) + histo->GetBinMax(0,hIter.GetInstanceIdentifier())) / 2.0;
123  }
124  }
125  }
126  labelStatisticsImageFilter->SetLabelInput(itkMask1);
127  labelStatisticsImageFilter->Update();
128  double median1 = labelStatisticsImageFilter->GetMedian(1);
129  double mean1 = labelStatisticsImageFilter->GetMean(1);
130  double modulo1 = 0;
131  {
132  auto histo = labelStatisticsImageFilter->GetHistogram(1);
133  double maxFrequency=0;
134  for (auto hIter=histo->Begin();hIter!=histo->End();++hIter)
135  {
136  if (maxFrequency < hIter.GetFrequency())
137  {
138  maxFrequency = hIter.GetFrequency();
139  modulo1 = (histo->GetBinMin(0,hIter.GetInstanceIdentifier()) + histo->GetBinMax(0,hIter.GetInstanceIdentifier())) / 2.0;
140  }
141  }
142  }
143 
144  double value0=0;
145  double value1=0;
146  switch (m_Area1)
147  {
149  value0=mean0; break;
151  value0=median0; break;
153  value0=modulo0; break;
154  }
155  switch (m_Area2)
156  {
158  value1=mean1; break;
160  value1=median1; break;
162  value1=modulo1; break;
163  }
164 
165  double offset = std::min(value0, value1);
166  double scaling = std::max(value0, value1) - offset;
167  if (scaling < 0.0001)
168  return;
169 
170  itk::ImageRegionIterator<ImageType> inIter(itkImage, itkImage->GetLargestPossibleRegion());
171  itk::ImageRegionIterator<ImageType> outIter(outImage, outImage->GetLargestPossibleRegion());
172  while (! inIter.IsAtEnd())
173  {
174  TPixel value = inIter.Value();
175  outIter.Set((value - offset) / scaling);
176  ++inIter;
177  ++outIter;
178  }
179 }
180 
182 {
184 }
void SetRequestedRegionToLargestPossibleRegion() override
void InternalComputeMask(itk::Image< TPixel, VImageDimension > *itkImage)
itk::Image< unsigned char, 3 > ImageType
static Vector3D offset
Image class for storing images.
Definition: mitkImage.h:72
#define AccessByItk(mitkImage, itkImageTypeFunction)
Access a MITK image by an ITK image.
static T max(T x, T y)
Definition: svm.cpp:56
static T min(T x, T y)
Definition: svm.cpp:53
InputImageType * GetInput(void)
void MITKCORE_EXPORT CastToItkImage(const mitk::Image *mitkImage, itk::SmartPointer< ItkOutputImageType > &itkOutputImage)
Cast an mitk::Image to an itk::Image with a specific type.
mitk::Image::Pointer mask
OutputType * GetOutput()
Get the output data of this image source object.