Medical Imaging Interaction Toolkit  2018.4.99-12ad79a3
Medical Imaging Interaction Toolkit
mitkGIFIntensityVolumeHistogramFeatures.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 // MITK
16 #include <mitkITKImageImport.h>
17 #include <mitkImageCast.h>
18 #include <mitkImageAccessByItk.h>
19 #include <mitkPixelTypeMultiplex.h>
20 
21 // ITK
22 #include <itkImageRegionIteratorWithIndex.h>
23 #include <itkNeighborhoodIterator.h>
24 // STL
25 #include <limits>
26 
27 struct GIFIntensityVolumeHistogramFeaturesParameters
28 {
30  std::string prefix;
31 };
32 
33 
34 template<typename TPixel, unsigned int VImageDimension>
35 static void
36 CalculateIntensityPeak(itk::Image<TPixel, VImageDimension>* itkImage, mitk::Image::Pointer mask, GIFIntensityVolumeHistogramFeaturesParameters params, mitk::GIFIntensityVolumeHistogramFeatures::FeatureListType & featureList)
37 {
38  typedef itk::Image<TPixel, VImageDimension> ImageType;
39  typedef itk::Image<unsigned short, VImageDimension> MaskType;
40 
41  typename MaskType::Pointer itkMask = MaskType::New();
42  mitk::CastToItkImage(mask, itkMask);
43 
44  mitk::IntensityQuantifier::Pointer quantifier = params.quantifier;
45  std::string prefix = params.prefix;
46 
47  itk::ImageRegionConstIterator<ImageType> iter(itkImage, itkImage->GetLargestPossibleRegion());
48  itk::ImageRegionConstIterator<MaskType> iterMask(itkMask, itkMask->GetLargestPossibleRegion());
49 
50  MITK_INFO << "Quantification: " << quantifier->GetMinimum() << " to " << quantifier->GetMaximum() << " with " << quantifier->GetBins()<< " bins";
51 
52  iter.GoToBegin();
53  iterMask.GoToBegin();
54  std::vector<double> hist;
55  hist.resize(quantifier->GetBins() , 0);
56 
57  int count = 0;
58  while (!iter.IsAtEnd())
59  {
60  if (iterMask.Get() > 0)
61  {
62  double value = iter.Get();
63  //std::size_t index = std::floor((value - minimum) / (maximum - minimum) * (bins-1));
64  std::size_t index = quantifier->IntensityToIndex(value);
65  ++count;
66  hist[index] += 1.0;// / count;
67  }
68  ++iterMask;
69  ++iter;
70  }
71 
72  bool notFoundIntenstiy010 = true;
73  bool notFoundIntenstiy090 = true;
74 
75  double intensity010 = -1;
76  double intensity090 = -1;
77  double fraction = 0;
78  double auc = 0;
79  bool firstRound = true;
80  for (int i = quantifier->GetBins()-1; i >= 0; --i)
81  {
82  hist[i] /= count;
83  hist[i] += fraction;
84  fraction = hist[i];
85  if (!firstRound)
86  {
87  auc += 0.5 * (hist[i] + hist[i+1]) / (quantifier->GetBins()-1);
88  }
89  firstRound = false;
90 
91  if (notFoundIntenstiy010 && fraction > 0.1)
92  {
93  intensity010 = quantifier->IndexToMeanIntensity(i + 1);
94  notFoundIntenstiy010 = false;
95  }
96  if (notFoundIntenstiy090 && fraction > 0.9)
97  {
98  intensity090 = quantifier->IndexToMeanIntensity(i + 1);
99  notFoundIntenstiy090 = false;
100  }
101  }
102 
103  unsigned int index010 = std::ceil(quantifier->GetBins() * 0.1);
104  unsigned int index090 = std::floor(quantifier->GetBins() * 0.9);
105 
106  featureList.push_back(std::make_pair(prefix + "Volume fraction at 0.10 intensity", hist[index010]));
107  featureList.push_back(std::make_pair(prefix + "Volume fraction at 0.90 intensity", hist[index090]));
108  featureList.push_back(std::make_pair(prefix + "Intensity at 0.10 volume", intensity010));
109  featureList.push_back(std::make_pair(prefix + "Intensity at 0.90 volume", intensity090));
110  featureList.push_back(std::make_pair(prefix + "Difference volume fraction at 0.10 and 0.90 intensity", std::abs<double>(hist[index010] - hist[index090])));
111  featureList.push_back(std::make_pair(prefix + "Difference intensity at 0.10 and 0.90 volume", std::abs<double>(intensity090 - intensity010)));
112  featureList.push_back(std::make_pair(prefix + "Area under IVH curve", auc));
113  //featureList.push_back(std::make_pair("Local Intensity Global Intensity Peak", globalPeakValue));
114 }
115 
116 
117 mitk::GIFIntensityVolumeHistogramFeatures::GIFIntensityVolumeHistogramFeatures()
118 {
119  SetLongName("intensity-volume-histogram");
120  SetShortName("ivoh");
121  SetFeatureClassName("Intensity Volume Histogram");
122 }
123 
124 mitk::GIFIntensityVolumeHistogramFeatures::FeatureListType mitk::GIFIntensityVolumeHistogramFeatures::CalculateFeatures(const Image::Pointer & image, const Image::Pointer &mask)
125 {
126  InitializeQuantifier(image, mask, 1000);
127  FeatureListType featureList;
128  GIFIntensityVolumeHistogramFeaturesParameters params;
129  params.quantifier = GetQuantifier();
130  params.prefix = FeatureDescriptionPrefix();
131  AccessByItk_3(image, CalculateIntensityPeak, mask, params, featureList);
132  return featureList;
133 }
134 
135 mitk::GIFIntensityVolumeHistogramFeatures::FeatureNameListType mitk::GIFIntensityVolumeHistogramFeatures::GetFeatureNames()
136 {
137  FeatureNameListType featureList;
138  return featureList;
139 }
140 
141 
142 void mitk::GIFIntensityVolumeHistogramFeatures::AddArguments(mitkCommandLineParser &parser)
143 {
144  std::string name = GetOptionPrefix();
145 
146  parser.addArgument(GetLongName(), name, mitkCommandLineParser::Bool, "Use Local Intensity", "calculates local intensity based features", us::Any());
147  AddQuantifierArguments(parser);
148 }
149 
150 void
151 mitk::GIFIntensityVolumeHistogramFeatures::CalculateFeaturesUsingParameters(const Image::Pointer & feature, const Image::Pointer &mask, const Image::Pointer &, FeatureListType &featureList)
152 {
153  InitializeQuantifierFromParameters(feature, mask, 1000);
154 
155  auto parsedArgs = GetParameter();
156  if (parsedArgs.count(GetLongName()))
157  {
158  MITK_INFO << "Start calculating local intensity features ....";
159  auto localResults = this->CalculateFeatures(feature, mask);
160  featureList.insert(featureList.end(), localResults.begin(), localResults.end());
161  MITK_INFO << "Finished calculating local intensity features....";
162  }
163 }
164 
165 std::string mitk::GIFIntensityVolumeHistogramFeatures::GetCurrentFeatureEncoding()
166 {
167  return QuantifierParameterString();
168 }
169 
#define AccessByItk_3(mitkImage, itkImageTypeFunction, arg1, arg2, arg3)
#define MITK_INFO
Definition: mitkLogMacros.h:18
itk::Image< unsigned char, 3 > ImageType
void addArgument(const std::string &longarg, const std::string &shortarg, Type type, const std::string &argLabel, const std::string &argHelp=std::string(), const us::Any &defaultValue=us::Any(), bool optional=true, bool ignoreRest=false, bool deprecated=false, mitkCommandLineParser::Channel channel=mitkCommandLineParser::Channel::None)
Definition: usAny.h:163
mitk::Image::Pointer image
void MITKCORE_EXPORT CastToItkImage(const mitk::Image *mitkImage, itk::SmartPointer< ItkOutputImageType > &itkOutputImage)
Cast an mitk::Image to an itk::Image with a specific type.
static void CalculateIntensityPeak(itk::Image< TPixel, VImageDimension > *itkImage, mitk::Image::Pointer mask, GIFIntensityVolumeHistogramFeaturesParameters params, mitk::GIFIntensityVolumeHistogramFeatures::FeatureListType &featureList)
mitk::Image::Pointer mask
void CalculateFeatures(mitk::CoocurenceMatrixHolder &holder, mitk::CoocurenceMatrixFeatures &results)