Medical Imaging Interaction Toolkit  2018.4.99-12ad79a3
Medical Imaging Interaction Toolkit
mitkGIFCooccurenceMatrix2.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 
20 // ITK
22 #include <itkShapedNeighborhoodIterator.h>
23 #include <itkImageRegionConstIterator.h>
24 
25 // STL
26 #include <sstream>
27 #include <cmath>
28 
29 namespace mitk
30 {
31  struct CoocurenceMatrixHolder
32  {
33  public:
34  CoocurenceMatrixHolder(double min, double max, int number);
35 
36  int IntensityToIndex(double intensity);
37  double IndexToMinIntensity(int index);
38  double IndexToMeanIntensity(int index);
39  double IndexToMaxIntensity(int index);
40 
41  double m_MinimumRange;
42  double m_MaximumRange;
43  double m_Stepsize;
44  int m_NumberOfBins;
45  Eigen::MatrixXd m_Matrix;
46 
47  };
48 
49  struct CoocurenceMatrixFeatures
50  {
51  CoocurenceMatrixFeatures() :
52  JointMaximum(0),
53  JointAverage(0),
54  JointVariance(0),
55  JointEntropy(0),
56  RowMaximum(0),
57  RowAverage(0),
58  RowVariance(0),
59  RowEntropy(0),
60  FirstRowColumnEntropy(0),
61  SecondRowColumnEntropy(0),
62  DifferenceAverage(0),
63  DifferenceVariance(0),
64  DifferenceEntropy(0),
65  SumAverage(0),
66  SumVariance(0),
67  SumEntropy(0),
68  AngularSecondMoment(0),
69  Contrast(0),
70  Dissimilarity(0),
71  InverseDifference(0),
72  InverseDifferenceNormalised(0),
73  InverseDifferenceMoment(0),
74  InverseDifferenceMomentNormalised(0),
75  InverseVariance(0),
76  Correlation(0),
77  Autocorrelation(0),
78  ClusterTendency(0),
79  ClusterShade(0),
80  ClusterProminence(0),
81  FirstMeasureOfInformationCorrelation(0),
82  SecondMeasureOfInformationCorrelation(0)
83  {
84  }
85 
86  public:
87  double JointMaximum;
88  double JointAverage;
89  double JointVariance;
90  double JointEntropy;
91  double RowMaximum;
92  double RowAverage;
93  double RowVariance;
94  double RowEntropy;
95  double FirstRowColumnEntropy;
96  double SecondRowColumnEntropy;
97  double DifferenceAverage;
98  double DifferenceVariance;
99  double DifferenceEntropy;
100  double SumAverage;
101  double SumVariance;
102  double SumEntropy;
103  double AngularSecondMoment;
104  double Contrast;
105  double Dissimilarity;
106  double InverseDifference;
107  double InverseDifferenceNormalised;
108  double InverseDifferenceMoment;
109  double InverseDifferenceMomentNormalised;
110  double InverseVariance;
111  double Correlation;
112  double Autocorrelation;
113  double ClusterTendency;
114  double ClusterShade;
115  double ClusterProminence;
116  double FirstMeasureOfInformationCorrelation;
117  double SecondMeasureOfInformationCorrelation;
118  };
119 
120 }
121 
122 static
123 void MatrixFeaturesTo(mitk::CoocurenceMatrixFeatures features,
124  std::string prefix,
125  mitk::GIFCooccurenceMatrix2::FeatureListType &featureList);
126 static
127 void CalculateMeanAndStdDevFeatures(std::vector<mitk::CoocurenceMatrixFeatures> featureList,
128  mitk::CoocurenceMatrixFeatures &meanFeature,
129  mitk::CoocurenceMatrixFeatures &stdFeature);
130 static
131 void NormalizeMatrixFeature(mitk::CoocurenceMatrixFeatures &features,
132  std::size_t number);
133 
134 
135 
136 
137 mitk::CoocurenceMatrixHolder::CoocurenceMatrixHolder(double min, double max, int number) :
138 m_MinimumRange(min),
139 m_MaximumRange(max),
140 m_NumberOfBins(number)
141 {
142  m_Matrix.resize(number, number);
143  m_Matrix.fill(0);
144  m_Stepsize = (max - min) / (number);
145 }
146 
147 int mitk::CoocurenceMatrixHolder::IntensityToIndex(double intensity)
148 {
149  int index = std::floor((intensity - m_MinimumRange) / m_Stepsize);
150  return std::max(0, std::min(index, m_NumberOfBins - 1));
151 }
152 
153 double mitk::CoocurenceMatrixHolder::IndexToMinIntensity(int index)
154 {
155  return m_MinimumRange + index * m_Stepsize;
156 }
157 double mitk::CoocurenceMatrixHolder::IndexToMeanIntensity(int index)
158 {
159  return m_MinimumRange + (index+0.5) * m_Stepsize;
160 }
161 double mitk::CoocurenceMatrixHolder::IndexToMaxIntensity(int index)
162 {
163  return m_MinimumRange + (index + 1) * m_Stepsize;
164 }
165 
166 template<typename TPixel, unsigned int VImageDimension>
167 void
168 CalculateCoOcMatrix(itk::Image<TPixel, VImageDimension>* itkImage,
169  itk::Image<unsigned short, VImageDimension>* mask,
170  itk::Offset<VImageDimension> offset,
171  int range,
172  mitk::CoocurenceMatrixHolder &holder)
173 {
174  typedef itk::Image<TPixel, VImageDimension> ImageType;
175  typedef itk::Image<unsigned short, VImageDimension> MaskImageType;
176  typedef itk::ShapedNeighborhoodIterator<ImageType> ShapeIterType;
177  typedef itk::ShapedNeighborhoodIterator<MaskImageType> ShapeMaskIterType;
178  typedef itk::ImageRegionConstIterator<ImageType> ConstIterType;
179  typedef itk::ImageRegionConstIterator<MaskImageType> ConstMaskIterType;
180 
181 
182  itk::Size<VImageDimension> radius;
183  radius.Fill(range+1);
184  ShapeIterType imageOffsetIter(radius, itkImage, itkImage->GetLargestPossibleRegion());
185  ShapeMaskIterType maskOffsetIter(radius, mask, mask->GetLargestPossibleRegion());
186  imageOffsetIter.ActivateOffset(offset);
187  maskOffsetIter.ActivateOffset(offset);
188  ConstIterType imageIter(itkImage, itkImage->GetLargestPossibleRegion());
189  ConstMaskIterType maskIter(mask, mask->GetLargestPossibleRegion());
190  // iterator.GetIndex() + ci.GetNeighborhoodOffset()
191  auto region = mask->GetLargestPossibleRegion();
192 
193 
194  while (!maskIter.IsAtEnd())
195  {
196  auto ciMask = maskOffsetIter.Begin();
197  auto ciValue = imageOffsetIter.Begin();
198  if (maskIter.Value() > 0 &&
199  ciMask.Get() > 0 &&
200  imageIter.Get() == imageIter.Get() &&
201  ciValue.Get() == ciValue.Get() &&
202  region.IsInside(maskOffsetIter.GetIndex() + ciMask.GetNeighborhoodOffset()))
203  {
204  int i = holder.IntensityToIndex(imageIter.Get());
205  int j = holder.IntensityToIndex(ciValue.Get());
206  holder.m_Matrix(i, j) += 1;
207  holder.m_Matrix(j, i) += 1;
208  }
209  ++imageOffsetIter;
210  ++maskOffsetIter;
211  ++imageIter;
212  ++maskIter;
213  }
214 }
215 
217  mitk::CoocurenceMatrixHolder &holder,
218  mitk::CoocurenceMatrixFeatures & results
219  )
220 {
221  auto pijMatrix = holder.m_Matrix;
222  auto piMatrix = holder.m_Matrix;
223  auto pjMatrix = holder.m_Matrix;
224  double Ng = holder.m_NumberOfBins;
225  int NgSize = holder.m_NumberOfBins;
226  pijMatrix /= pijMatrix.sum();
227  piMatrix.rowwise().normalize();
228  pjMatrix.colwise().normalize();
229 
230  for (int i = 0; i < holder.m_NumberOfBins; ++i)
231  for (int j = 0; j < holder.m_NumberOfBins; ++j)
232  {
233  if (pijMatrix(i, j) != pijMatrix(i, j))
234  pijMatrix(i, j) = 0;
235  if (piMatrix(i, j) != piMatrix(i, j))
236  piMatrix(i, j) = 0;
237  if (pjMatrix(i, j) != pjMatrix(i, j))
238  pjMatrix(i, j) = 0;
239  }
240 
241  Eigen::VectorXd piVector = pijMatrix.colwise().sum();
242  Eigen::VectorXd pjVector = pijMatrix.rowwise().sum();
243  double sigmai = 0;;
244  for (int i = 0; i < holder.m_NumberOfBins; ++i)
245  {
246  double iInt = i + 1;// holder.IndexToMeanIntensity(i);
247  results.RowAverage += iInt * piVector(i);
248  if (piVector(i) > 0)
249  {
250  results.RowEntropy -= piVector(i) * std::log(piVector(i)) / std::log(2);
251  }
252  }
253  for (int i = 0; i < holder.m_NumberOfBins; ++i)
254  {
255  double iInt = i + 1; // holder.IndexToMeanIntensity(i);
256  results.RowVariance += (iInt - results.RowAverage)*(iInt - results.RowAverage) * piVector(i);
257  }
258  results.RowMaximum = piVector.maxCoeff();
259  sigmai = std::sqrt(results.RowVariance);
260 
261  Eigen::VectorXd pimj(NgSize);
262  pimj.fill(0);
263  Eigen::VectorXd pipj(2*NgSize);
264  pipj.fill(0);
265 
266 
267  results.JointMaximum += pijMatrix.maxCoeff();
268 
269  for (int i = 0; i < holder.m_NumberOfBins; ++i)
270  {
271  for (int j = 0; j < holder.m_NumberOfBins; ++j)
272  {
273  //double iInt = holder.IndexToMeanIntensity(i);
274  //double jInt = holder.IndexToMeanIntensity(j);
275  double iInt = i + 1;// holder.IndexToMeanIntensity(i);
276  double jInt = j + 1;// holder.IndexToMeanIntensity(j);
277  double pij = pijMatrix(i, j);
278 
279  int deltaK = (i - j)>0?(i-j) : (j-i);
280  pimj(deltaK) += pij;
281  pipj(i + j) += pij;
282 
283  results.JointAverage += iInt * pij;
284  if (pij > 0)
285  {
286  results.JointEntropy -= pij * std::log(pij) / std::log(2);
287  results.FirstRowColumnEntropy -= pij * std::log(piVector(i)*pjVector(j)) / std::log(2);
288  }
289  if (piVector(i) > 0 && pjVector(j) > 0 )
290  {
291  results.SecondRowColumnEntropy -= piVector(i)*pjVector(j) * std::log(piVector(i)*pjVector(j)) / std::log(2);
292  }
293  results.AngularSecondMoment += pij*pij;
294  results.Contrast += (iInt - jInt)* (iInt - jInt) * pij;
295  results.Dissimilarity += std::abs<double>(iInt - jInt) * pij;
296  results.InverseDifference += pij / (1 + (std::abs<double>(iInt - jInt)));
297  results.InverseDifferenceNormalised += pij / (1 + (std::abs<double>(iInt - jInt) / Ng));
298  results.InverseDifferenceMoment += pij / (1 + (iInt - jInt)*(iInt - jInt));
299  results.InverseDifferenceMomentNormalised += pij / (1 + (iInt - jInt)*(iInt - jInt)/Ng/Ng);
300  results.Autocorrelation += iInt*jInt * pij;
301  double cluster = (iInt + jInt - 2 * results.RowAverage);
302  results.ClusterTendency += cluster*cluster * pij;
303  results.ClusterShade += cluster*cluster*cluster * pij;
304  results.ClusterProminence += cluster*cluster*cluster*cluster * pij;
305  if (iInt != jInt)
306  {
307  results.InverseVariance += pij / (iInt - jInt) / (iInt - jInt);
308  }
309  }
310  }
311  results.Correlation = 1 / sigmai / sigmai * (-results.RowAverage*results.RowAverage+ results.Autocorrelation);
312  results.FirstMeasureOfInformationCorrelation = (results.JointEntropy - results.FirstRowColumnEntropy) / results.RowEntropy;
313  if (results.JointEntropy < results.SecondRowColumnEntropy)
314  {
315  results.SecondMeasureOfInformationCorrelation = sqrt(1 - exp(-2 * (results.SecondRowColumnEntropy - results.JointEntropy)));
316  }
317  else
318  {
319  results.SecondMeasureOfInformationCorrelation = 0;
320  }
321 
322  for (int i = 0; i < holder.m_NumberOfBins; ++i)
323  {
324  for (int j = 0; j < holder.m_NumberOfBins; ++j)
325  {
326  //double iInt = holder.IndexToMeanIntensity(i);
327  //double jInt = holder.IndexToMeanIntensity(j);
328  double iInt = i + 1;
329  double pij = pijMatrix(i, j);
330 
331  results.JointVariance += (iInt - results.JointAverage)* (iInt - results.JointAverage)*pij;
332  }
333  }
334 
335  for (int k = 0; k < NgSize; ++k)
336  {
337  results.DifferenceAverage += k* pimj(k);
338  if (pimj(k) > 0)
339  {
340  results.DifferenceEntropy -= pimj(k) * log(pimj(k)) / std::log(2);
341  }
342  }
343  for (int k = 0; k < NgSize; ++k)
344  {
345  results.DifferenceVariance += (results.DifferenceAverage-k)* (results.DifferenceAverage-k)*pimj(k);
346  }
347 
348 
349  for (int k = 0; k <2* NgSize ; ++k)
350  {
351  results.SumAverage += (2+k)* pipj(k);
352  if (pipj(k) > 0)
353  {
354  results.SumEntropy -= pipj(k) * log(pipj(k)) / std::log(2);
355  }
356  }
357  for (int k = 0; k < 2*NgSize; ++k)
358  {
359  results.SumVariance += (2+k - results.SumAverage)* (2+k - results.SumAverage)*pipj(k);
360  }
361 
362  //MITK_INFO << std::endl << holder.m_Matrix;
363  //MITK_INFO << std::endl << pijMatrix;
364  //MITK_INFO << std::endl << piMatrix;
365  //MITK_INFO << std::endl << pjMatrix;
366 
367  //for (int i = 0; i < holder.m_NumberOfBins; ++i)
368  //{
369  // MITK_INFO << "Bin " << i << " Min: " << holder.IndexToMinIntensity(i) << " Max: " << holder.IndexToMaxIntensity(i);
370  //}
371  //MITK_INFO << pimj;
372  //MITK_INFO << pipj;
373 
374 }
375 
376 template<typename TPixel, unsigned int VImageDimension>
377 void
378 CalculateCoocurenceFeatures(itk::Image<TPixel, VImageDimension>* itkImage, mitk::Image::Pointer mask, mitk::GIFCooccurenceMatrix2::FeatureListType & featureList, mitk::GIFCooccurenceMatrix2::GIFCooccurenceMatrix2Configuration config)
379 {
380  typedef itk::Image<unsigned short, VImageDimension> MaskType;
381  typedef itk::Neighborhood<TPixel, VImageDimension > NeighborhoodType;
382  typedef itk::Offset<VImageDimension> OffsetType;
383 
385  double rangeMin = config.MinimumIntensity;
386  double rangeMax = config.MaximumIntensity;
387  int numberOfBins = config.Bins;
388 
389  typename MaskType::Pointer maskImage = MaskType::New();
390  mitk::CastToItkImage(mask, maskImage);
391 
392  //Find possible directions
393  std::vector < itk::Offset<VImageDimension> > offsetVector;
394  NeighborhoodType hood;
395  hood.SetRadius(1);
396  unsigned int centerIndex = hood.GetCenterNeighborhoodIndex();
397  OffsetType offset;
398  for (unsigned int d = 0; d < centerIndex; d++)
399  {
400  offset = hood.GetOffset(d);
401  bool useOffset = true;
402  for (unsigned int i = 0; i < VImageDimension; ++i)
403  {
404  offset[i] *= config.range;
405  if (config.direction == i + 2 && offset[i] != 0)
406  {
407  useOffset = false;
408  }
409  }
410  if (useOffset)
411  {
412  offsetVector.push_back(offset);
413  }
414  }
415  if (config.direction == 1)
416  {
417  offsetVector.clear();
418  offset[0] = 0;
419  offset[1] = 0;
420  offset[2] = 1;
421  }
422 
423  std::vector<mitk::CoocurenceMatrixFeatures> resultVector;
424  mitk::CoocurenceMatrixHolder holderOverall(rangeMin, rangeMax, numberOfBins);
425  mitk::CoocurenceMatrixFeatures overallFeature;
426  for (std::size_t i = 0; i < offsetVector.size(); ++i)
427  {
428  if (config.direction > 1)
429  {
430  if (offsetVector[i][config.direction - 2] != 0)
431  {
432  continue;
433  }
434  }
435 
436 
437  offset = offsetVector[i];
438  mitk::CoocurenceMatrixHolder holder(rangeMin, rangeMax, numberOfBins);
439  mitk::CoocurenceMatrixFeatures coocResults;
440  CalculateCoOcMatrix<TPixel, VImageDimension>(itkImage, maskImage, offset, config.range, holder);
441  holderOverall.m_Matrix += holder.m_Matrix;
442  CalculateFeatures(holder, coocResults);
443  resultVector.push_back(coocResults);
444  }
445  CalculateFeatures(holderOverall, overallFeature);
446  //NormalizeMatrixFeature(overallFeature, offsetVector.size());
447 
448 
449  mitk::CoocurenceMatrixFeatures featureMean;
450  mitk::CoocurenceMatrixFeatures featureStd;
451  CalculateMeanAndStdDevFeatures(resultVector, featureMean, featureStd);
452 
453  std::ostringstream ss;
454  ss << config.range;
455  std::string strRange = ss.str();
456 
457  MatrixFeaturesTo(overallFeature, config.prefix + "Overall ", featureList);
458  MatrixFeaturesTo(featureMean, config.prefix + "Mean ", featureList);
459  MatrixFeaturesTo(featureStd, config.prefix + "Std.Dev. ", featureList);
460 }
461 
462 
463 static
464 void MatrixFeaturesTo(mitk::CoocurenceMatrixFeatures features,
465  std::string prefix,
466  mitk::GIFCooccurenceMatrix2::FeatureListType &featureList)
467 {
468  featureList.push_back(std::make_pair(prefix + "Joint Maximum", features.JointMaximum));
469  featureList.push_back(std::make_pair(prefix + "Joint Average", features.JointAverage));
470  featureList.push_back(std::make_pair(prefix + "Joint Variance", features.JointVariance));
471  featureList.push_back(std::make_pair(prefix + "Joint Entropy", features.JointEntropy));
472  featureList.push_back(std::make_pair(prefix + "Difference Average", features.DifferenceAverage));
473  featureList.push_back(std::make_pair(prefix + "Difference Variance", features.DifferenceVariance));
474  featureList.push_back(std::make_pair(prefix + "Difference Entropy", features.DifferenceEntropy));
475  featureList.push_back(std::make_pair(prefix + "Sum Average", features.SumAverage));
476  featureList.push_back(std::make_pair(prefix + "Sum Variance", features.SumVariance));
477  featureList.push_back(std::make_pair(prefix + "Sum Entropy", features.SumEntropy));
478  featureList.push_back(std::make_pair(prefix + "Angular Second Moment", features.AngularSecondMoment));
479  featureList.push_back(std::make_pair(prefix + "Contrast", features.Contrast));
480  featureList.push_back(std::make_pair(prefix + "Dissimilarity", features.Dissimilarity));
481  featureList.push_back(std::make_pair(prefix + "Inverse Difference", features.InverseDifference));
482  featureList.push_back(std::make_pair(prefix + "Inverse Difference Normalized", features.InverseDifferenceNormalised));
483  featureList.push_back(std::make_pair(prefix + "Inverse Difference Moment", features.InverseDifferenceMoment));
484  featureList.push_back(std::make_pair(prefix + "Inverse Difference Moment Normalized", features.InverseDifferenceMomentNormalised));
485  featureList.push_back(std::make_pair(prefix + "Inverse Variance", features.InverseVariance));
486  featureList.push_back(std::make_pair(prefix + "Correlation", features.Correlation));
487  featureList.push_back(std::make_pair(prefix + "Autocorrelation", features.Autocorrelation));
488  featureList.push_back(std::make_pair(prefix + "Cluster Tendency", features.ClusterTendency));
489  featureList.push_back(std::make_pair(prefix + "Cluster Shade", features.ClusterShade));
490  featureList.push_back(std::make_pair(prefix + "Cluster Prominence", features.ClusterProminence));
491  featureList.push_back(std::make_pair(prefix + "First Measure of Information Correlation", features.FirstMeasureOfInformationCorrelation));
492  featureList.push_back(std::make_pair(prefix + "Second Measure of Information Correlation", features.SecondMeasureOfInformationCorrelation));
493  featureList.push_back(std::make_pair(prefix + "Row Maximum", features.RowMaximum));
494  featureList.push_back(std::make_pair(prefix + "Row Average", features.RowAverage));
495  featureList.push_back(std::make_pair(prefix + "Row Variance", features.RowVariance));
496  featureList.push_back(std::make_pair(prefix + "Row Entropy", features.RowEntropy));
497  featureList.push_back(std::make_pair(prefix + "First Row-Column Entropy", features.FirstRowColumnEntropy));
498  featureList.push_back(std::make_pair(prefix + "Second Row-Column Entropy", features.SecondRowColumnEntropy));
499 }
500 
501 static
502 void CalculateMeanAndStdDevFeatures(std::vector<mitk::CoocurenceMatrixFeatures> featureList,
503  mitk::CoocurenceMatrixFeatures &meanFeature,
504  mitk::CoocurenceMatrixFeatures &stdFeature)
505 {
506 #define ADDFEATURE(a) \
507  if ( ! (featureList[i].a == featureList[i].a)) featureList[i].a = 0; \
508  meanFeature.a += featureList[i].a;stdFeature.a += featureList[i].a*featureList[i].a
509 #define CALCVARIANCE(a) stdFeature.a =sqrt(stdFeature.a - meanFeature.a*meanFeature.a)
510 
511  for (std::size_t i = 0; i < featureList.size(); ++i)
512  {
513  ADDFEATURE(JointMaximum);
514  ADDFEATURE(JointAverage);
515  ADDFEATURE(JointVariance);
516  ADDFEATURE(JointEntropy);
517  ADDFEATURE(RowMaximum);
518  ADDFEATURE(RowAverage);
519  ADDFEATURE(RowVariance);
520  ADDFEATURE(RowEntropy);
521  ADDFEATURE(FirstRowColumnEntropy);
522  ADDFEATURE(SecondRowColumnEntropy);
523  ADDFEATURE(DifferenceAverage);
524  ADDFEATURE(DifferenceVariance);
525  ADDFEATURE(DifferenceEntropy);
526  ADDFEATURE(SumAverage);
527  ADDFEATURE(SumVariance);
528  ADDFEATURE(SumEntropy);
529  ADDFEATURE(AngularSecondMoment);
530  ADDFEATURE(Contrast);
531  ADDFEATURE(Dissimilarity);
532  ADDFEATURE(InverseDifference);
533  ADDFEATURE(InverseDifferenceNormalised);
534  ADDFEATURE(InverseDifferenceMoment);
535  ADDFEATURE(InverseDifferenceMomentNormalised);
536  ADDFEATURE(InverseVariance);
537  ADDFEATURE(Correlation);
538  ADDFEATURE(Autocorrelation);
539  ADDFEATURE(ClusterShade);
540  ADDFEATURE(ClusterTendency);
541  ADDFEATURE(ClusterProminence);
542  ADDFEATURE(FirstMeasureOfInformationCorrelation);
543  ADDFEATURE(SecondMeasureOfInformationCorrelation);
544  }
545  NormalizeMatrixFeature(meanFeature, featureList.size());
546  NormalizeMatrixFeature(stdFeature, featureList.size());
547 
548  CALCVARIANCE(JointMaximum);
549  CALCVARIANCE(JointAverage);
550  CALCVARIANCE(JointVariance);
551  CALCVARIANCE(JointEntropy);
552  CALCVARIANCE(RowMaximum);
553  CALCVARIANCE(RowAverage);
554  CALCVARIANCE(RowVariance);
555  CALCVARIANCE(RowEntropy);
556  CALCVARIANCE(FirstRowColumnEntropy);
557  CALCVARIANCE(SecondRowColumnEntropy);
558  CALCVARIANCE(DifferenceAverage);
559  CALCVARIANCE(DifferenceVariance);
560  CALCVARIANCE(DifferenceEntropy);
561  CALCVARIANCE(SumAverage);
562  CALCVARIANCE(SumVariance);
563  CALCVARIANCE(SumEntropy);
564  CALCVARIANCE(AngularSecondMoment);
565  CALCVARIANCE(Contrast);
566  CALCVARIANCE(Dissimilarity);
567  CALCVARIANCE(InverseDifference);
568  CALCVARIANCE(InverseDifferenceNormalised);
569  CALCVARIANCE(InverseDifferenceMoment);
570  CALCVARIANCE(InverseDifferenceMomentNormalised);
571  CALCVARIANCE(InverseVariance);
572  CALCVARIANCE(Correlation);
573  CALCVARIANCE(Autocorrelation);
574  CALCVARIANCE(ClusterShade);
575  CALCVARIANCE(ClusterTendency);
576  CALCVARIANCE(ClusterProminence);
577  CALCVARIANCE(FirstMeasureOfInformationCorrelation);
578  CALCVARIANCE(SecondMeasureOfInformationCorrelation);
579 
580 #undef ADDFEATURE
581 #undef CALCVARIANCE
582 }
583 
584 static
585 void NormalizeMatrixFeature(mitk::CoocurenceMatrixFeatures &features,
586  std::size_t number)
587 {
588  features.JointMaximum = features.JointMaximum / number;
589  features.JointAverage = features.JointAverage / number;
590  features.JointVariance = features.JointVariance / number;
591  features.JointEntropy = features.JointEntropy / number;
592  features.RowMaximum = features.RowMaximum / number;
593  features.RowAverage = features.RowAverage / number;
594  features.RowVariance = features.RowVariance / number;
595  features.RowEntropy = features.RowEntropy / number;
596  features.FirstRowColumnEntropy = features.FirstRowColumnEntropy / number;
597  features.SecondRowColumnEntropy = features.SecondRowColumnEntropy / number;
598  features.DifferenceAverage = features.DifferenceAverage / number;
599  features.DifferenceVariance = features.DifferenceVariance / number;
600  features.DifferenceEntropy = features.DifferenceEntropy / number;
601  features.SumAverage = features.SumAverage / number;
602  features.SumVariance = features.SumVariance / number;
603  features.SumEntropy = features.SumEntropy / number;
604  features.AngularSecondMoment = features.AngularSecondMoment / number;
605  features.Contrast = features.Contrast / number;
606  features.Dissimilarity = features.Dissimilarity / number;
607  features.InverseDifference = features.InverseDifference / number;
608  features.InverseDifferenceNormalised = features.InverseDifferenceNormalised / number;
609  features.InverseDifferenceMoment = features.InverseDifferenceMoment / number;
610  features.InverseDifferenceMomentNormalised = features.InverseDifferenceMomentNormalised / number;
611  features.InverseVariance = features.InverseVariance / number;
612  features.Correlation = features.Correlation / number;
613  features.Autocorrelation = features.Autocorrelation / number;
614  features.ClusterShade = features.ClusterShade / number;
615  features.ClusterTendency = features.ClusterTendency / number;
616  features.ClusterProminence = features.ClusterProminence / number;
617  features.FirstMeasureOfInformationCorrelation = features.FirstMeasureOfInformationCorrelation / number;
618  features.SecondMeasureOfInformationCorrelation = features.SecondMeasureOfInformationCorrelation / number;
619 }
620 
622 m_Range(1.0)
623 {
624  SetShortName("cooc2");
625  SetLongName("cooccurence2");
626  SetFeatureClassName("Co-occurenced Based Features");
627 }
628 
629 mitk::GIFCooccurenceMatrix2::FeatureListType mitk::GIFCooccurenceMatrix2::CalculateFeatures(const Image::Pointer & image, const Image::Pointer &mask)
630 {
631  InitializeQuantifier(image, mask);
632 
633  FeatureListType featureList;
634 
636  config.direction = GetDirection();
637  config.range = m_Range;
638 
639  config.MinimumIntensity = GetQuantifier()->GetMinimum();
640  config.MaximumIntensity = GetQuantifier()->GetMaximum();
641  config.Bins = GetQuantifier()->GetBins();
642  config.prefix = FeatureDescriptionPrefix();
643 
644  AccessByItk_3(image, CalculateCoocurenceFeatures, mask, featureList,config);
645 
646  return featureList;
647 }
648 
650 {
651  FeatureNameListType featureList;
652  return featureList;
653 }
654 
655 
656 
657 
659 {
660  std::string name = GetOptionPrefix();
661 
662  parser.addArgument(GetLongName(), name, mitkCommandLineParser::Bool, "Use Co-occurence matrix", "calculates Co-occurence based features (new implementation)", us::Any());
663  parser.addArgument(name+"::range", name+"::range", mitkCommandLineParser::String, "Cooc 2 Range", "Define the range that is used (Semicolon-separated)", us::Any());
664  AddQuantifierArguments(parser);
665 }
666 
667 void
668 mitk::GIFCooccurenceMatrix2::CalculateFeaturesUsingParameters(const Image::Pointer & feature, const Image::Pointer &, const Image::Pointer &maskNoNAN, FeatureListType &featureList)
669 {
670  auto parsedArgs = GetParameter();
671  std::string name = GetOptionPrefix();
672 
673  if (parsedArgs.count(GetLongName()))
674  {
675  InitializeQuantifierFromParameters(feature, maskNoNAN);
676  std::vector<double> ranges;
677 
678  if (parsedArgs.count(name + "::range"))
679  {
680  ranges = SplitDouble(parsedArgs[name + "::range"].ToString(), ';');
681  }
682  else
683  {
684  ranges.push_back(1);
685  }
686 
687  for (std::size_t i = 0; i < ranges.size(); ++i)
688  {
689  MITK_INFO << "Start calculating coocurence with range " << ranges[i] << "....";
690  this->SetRange(ranges[i]);
691  auto localResults = this->CalculateFeatures(feature, maskNoNAN);
692  featureList.insert(featureList.end(), localResults.begin(), localResults.end());
693  MITK_INFO << "Finished calculating coocurence with range " << ranges[i] << "....";
694  }
695  }
696 }
697 
698 
700 {
701  std::ostringstream ss;
702  ss << m_Range;
703  std::string strRange = ss.str();
704  return QuantifierParameterString() + "_Range-" + ss.str();
705 }
#define AccessByItk_3(mitkImage, itkImageTypeFunction, arg1, arg2, arg3)
float k(1.0)
void CalculateCoOcMatrix(itk::Image< TPixel, VImageDimension > *itkImage, itk::Image< unsigned short, VImageDimension > *mask, itk::Offset< VImageDimension > offset, int range, mitk::CoocurenceMatrixHolder &holder)
#define MITK_INFO
Definition: mitkLogMacros.h:18
itk::Image< unsigned char, 3 > ImageType
virtual void SetLongName(std::string _arg)
DataCollection - Class to facilitate loading/accessing structured data.
static void MatrixFeaturesTo(mitk::CoocurenceMatrixFeatures features, std::string prefix, mitk::GIFCooccurenceMatrix2::FeatureListType &featureList)
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)
virtual void SetShortName(std::string _arg)
static void CalculateMeanAndStdDevFeatures(std::vector< mitk::CoocurenceMatrixFeatures > featureList, mitk::CoocurenceMatrixFeatures &meanFeature, mitk::CoocurenceMatrixFeatures &stdFeature)
static Vector3D offset
std::string FeatureDescriptionPrefix()
Returns a string that encodes the feature class name.
#define ADDFEATURE(a)
virtual std::string GetLongName() const
void CalculateCoocurenceFeatures(itk::Image< TPixel, VImageDimension > *itkImage, mitk::Image::Pointer mask, mitk::GIFCooccurenceMatrix2::FeatureListType &featureList, mitk::GIFCooccurenceMatrix2::GIFCooccurenceMatrix2Configuration config)
Definition: usAny.h:163
void AddArguments(mitkCommandLineParser &parser) override
itk::Image< unsigned char, 3 > MaskImageType
Definition: CLBrainMask.cpp:32
FeatureListType CalculateFeatures(const Image::Pointer &image, const Image::Pointer &feature) override
Calculates the Cooccurence-Matrix based features for this class.
virtual int GetDirection() const
static T max(T x, T y)
Definition: svm.cpp:56
mitk::Image::Pointer image
void AddQuantifierArguments(mitkCommandLineParser &parser)
static T min(T x, T y)
Definition: svm.cpp:53
mitkClassMacro(AbstractGlobalImageFeature, BaseData) typedef std typedef std::vector< std::string > FeatureNameListType
std::vector< double > SplitDouble(std::string str, char delimiter)
virtual void SetFeatureClassName(std::string _arg)
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
FeatureListType CalculateFeatures(const Image::Pointer &image, const Image::Pointer &feature) override
Calculates the Cooccurence-Matrix based features for this class.
virtual IntensityQuantifier::Pointer GetQuantifier()
void InitializeQuantifierFromParameters(const Image::Pointer &feature, const Image::Pointer &mask, unsigned int defaultBins=256)
FeatureNameListType GetFeatureNames() override
Returns a list of the names of all features that are calculated from this class.
void CalculateFeaturesUsingParameters(const Image::Pointer &feature, const Image::Pointer &mask, const Image::Pointer &maskNoNAN, FeatureListType &featureList) override
Calculates the feature of this abstact interface. Does not necessarily considers the parameter settin...
std::string GetCurrentFeatureEncoding() override
Adds an additional Separator to the name of the feature, which encodes the used parameters.
void InitializeQuantifier(const Image::Pointer &feature, const Image::Pointer &mask, unsigned int defaultBins=256)
#define CALCVARIANCE(a)
virtual ParameterTypes GetParameter() const
static void NormalizeMatrixFeature(mitk::CoocurenceMatrixFeatures &features, std::vcl_size_t number)
virtual void SetRange(double _arg)