ANIMA  4.0
animaClassificationStrategy.h
Go to the documentation of this file.
1 #pragma once
2 
5 #include <itkGaussianMembershipFunction.h>
6 #include <map>
7 
8 namespace anima
9 {
10 
14 template <typename TInputImage, typename TMaskImage>
15 class ClassificationStrategy: public itk::ProcessObject
16 {
17 public:
18 
21  typedef itk::ProcessObject Superclass;
22  typedef itk::SmartPointer <Self> Pointer;
23  typedef itk::SmartPointer <const Self> ConstPointer;
24 
26  itkNewMacro(Self);
27 
29  itkTypeMacro(ClassificationStrategy, itk::ProcessObject);
30 
32  typedef TInputImage InputImageType;
33 
35  typedef TMaskImage MaskImageType;
36 
39 
40  typedef itk::VariableLengthVector<double> MeasurementVectorType;
41  typedef itk::Statistics::GaussianMembershipFunction< MeasurementVectorType > GaussianFunctionType;
42 
46  void Update() ITK_OVERRIDE;
47 
52  void SetEstimator( GaussianREMEstimatorPointerType theValue ){m_Estimator=theValue;}
53 
54 
59 
60 
64 
69  bool GetSolutionMap(std::map< double,std::vector<GaussianFunctionType::Pointer> > &solution, std::map< double, std::vector<double> > &solutionAlpha, int step=-1);
70 
71  void SetStrategy( std::vector< unsigned int >& ems, std::vector<unsigned int> &iters );
72 
73  itkSetMacro(EM_Mode, bool);
74  itkGetMacro(EM_Mode, bool);
75 
76  void SetNumberOfIterations(std::vector<unsigned int> NumberOfIterations){m_NumberOfIterations=NumberOfIterations;}
77  std::vector<unsigned int> GetNumberOfIterations() const {return m_NumberOfIterations;}
78 
79  void SetNumberOfEstimators(std::vector<unsigned int> NumberOfEstimators){m_NumberOfEstimators=NumberOfEstimators;}
80  std::vector<unsigned int> GetNumberOfEstimators() const {return m_NumberOfEstimators;}
81 
82 
83 protected:
84 
85  ClassificationStrategy(const Self&); //purposely not implemented
86  void operator=(const Self&); //purposely not implemented
87 
89  {
90  m_NumberOfIterations.resize(1,100);
91  m_NumberOfEstimators.resize(1,1);
92  m_EM_Mode=true;
93  }
94 
96 
97 
101  bool sameModel( std::vector<GaussianFunctionType::Pointer> &mod1,std::vector<GaussianFunctionType::Pointer> &mod2);
102 
105  GaussianREMEstimatorPointerType m_Estimator;
106 
110 
114  std::vector< std::map<double, std::vector<GaussianFunctionType::Pointer> > > m_ListGaussianModels;
115  std::vector< std::map<double, std::vector<double> > > m_ListAlphas;
116 
119  std::vector<unsigned int> m_NumberOfIterations;
120 
123  std::vector<unsigned int> m_NumberOfEstimators;
124 
127  bool m_EM_Mode;
128 
129 };
130 
131 }
132 
133 
bool sameModel(std::vector< GaussianFunctionType::Pointer > &mod1, std::vector< GaussianFunctionType::Pointer > &mod2)
compare two models
GaussianREMEstimatorType::Pointer GaussianREMEstimatorPointerType
std::vector< unsigned int > GetNumberOfIterations() const
std::vector< std::map< double, std::vector< GaussianFunctionType::Pointer > > > m_ListGaussianModels
vector with all the solutions found in each step For EM&#39;s the double value stores the -likelihood...
void SetInitializer(RandomInitializer::Pointer theValue)
set the m_RandomInitializer method
std::vector< unsigned int > GetNumberOfEstimators() const
Class performing a robust expectation-maximation (REM) algorithm. Allow finding the 3-classes NABT mo...
itk::Statistics::GaussianMembershipFunction< MeasurementVectorType > GaussianFunctionType
std::vector< std::map< double, std::vector< double > > > m_ListAlphas
void operator=(const Self &)
bool m_EM_Mode
modification of processing if an em is being used
itk::SmartPointer< const Self > ConstPointer
void SetNumberOfIterations(std::vector< unsigned int > NumberOfIterations)
void SetStrategy(std::vector< unsigned int > &ems, std::vector< unsigned int > &iters)
void SetEstimator(GaussianREMEstimatorPointerType theValue)
set the estimation algorithm to use The algorithm must be already configured with the joint histogram...
bool GetSolutionMap(std::map< double, std::vector< GaussianFunctionType::Pointer > > &solution, std::map< double, std::vector< double > > &solutionAlpha, int step=-1)
return the solutions map for each step
void SetNumberOfEstimators(std::vector< unsigned int > NumberOfEstimators)
Class initializing ramdomly a gaussian model.
itk::VariableLengthVector< double > MeasurementVectorType
itk::SmartPointer< Self > Pointer
RandomInitializer::Pointer m_RandomInitializer
Initialization object.
std::vector< unsigned int > m_NumberOfEstimators
vector with the number of random estimations for each step
GaussianREMEstimator< InputImageType, MaskImageType > GaussianREMEstimatorType
GaussianREMEstimatorPointerType m_Estimator
estimation algorithm object
std::vector< unsigned int > m_NumberOfIterations
vector with the iterations for each step
RandomInitializer * getInitializer()
return m_RandomInitializer object
itk::SmartPointer< Self > Pointer
void Update() ITK_OVERRIDE
executes the strategy