ANIMA  4.0
animaApproximateMCMSmoothingCostFunction.h
Go to the documentation of this file.
1 #pragma once
2 
3 #include <itkSingleValuedCostFunction.h>
4 #include <AnimaMCMSimilarityExport.h>
6 #include <animaMCMConstants.h>
7 
8 namespace anima
9 {
10 
11 class ANIMAMCMSIMILARITY_EXPORT ApproximateMCMSmoothingCostFunction : public itk::SingleValuedCostFunction
12 {
13 public:
16  typedef itk::SingleValuedCostFunction Superclass;
17  typedef itk::SmartPointer<Self> Pointer;
18  typedef itk::SmartPointer<const Self> ConstPointer;
19 
20  itkNewMacro(Self)
21 
22 
23  itkTypeMacro(ApproximateMCMSmoothingCostFunction, itk::SingleValuedCostFunction)
24 
25  typedef Superclass::MeasureType MeasureType;
26  typedef Superclass::DerivativeType DerivativeType;
27  typedef Superclass::ParametersType ParametersType;
28 
30  typedef MCModelType::Pointer MCMPointer;
31  typedef MCModelType::Vector3DType GradientType;
34 
35  virtual MeasureType GetValue(const ParametersType &parameters) const ITK_OVERRIDE;
36  virtual void GetDerivative(const ParametersType &parameters, DerivativeType &derivative) const ITK_OVERRIDE;
37 
38  void SetReferenceModels(const std::vector <MCMPointer> &refModels, const std::vector <GradientType> &gradients,
39  const double &smallDelta, const double &bigDelta, const std::vector <double> &gradientStrengths);
40  void SetMovingModels(const std::vector<MCMPointer> &movingModels, const std::vector <GradientType> &gradients,
41  const double &smallDelta, const double &bigDelta, const std::vector <double> &gradientStrengths);
42 
43  void SetGradientStrengths(const std::vector <double> &val);
44  void SetGradientDirections(const std::vector <GradientType> &val);
45 
46  void SetBValueWeightIndexes(const std::vector <unsigned int> &val);
47  void SetSphereWeights(const std::vector <double> &val);
48 
49  itkSetMacro(ParameterScale, double)
50  void SetSmallDelta(double val);
51  void SetBigDelta(double val);
52 
53  unsigned int GetNumberOfParameters() const ITK_OVERRIDE {return 1;}
54 
55 protected:
57  {
58  m_UpdatedData = false;
59  m_ConstantTerm = 0;
60  m_ParameterScale = 1.0e-3;
61  m_SmallDelta = anima::DiffusionSmallDelta;
62  m_BigDelta = anima::DiffusionBigDelta;
63  }
64 
66 
67 private:
68  ApproximateMCMSmoothingCostFunction(const Self&); //purposely not implemented
69  void operator=(const Self&); //purposely not implemented
70 
71  std::vector < std::vector <double> > m_ReferenceModelSignalValues;
72  std::vector < std::vector <double> > m_MovingModelSignalValues;
73 
74  std::vector <unsigned int> m_BValueWeightIndexes;
75  std::vector <double> m_GradientStrengths, m_SphereWeights;
76  std::vector <GradientType> m_GradientDirections;
77 
78  mutable bool m_UpdatedData;
79  mutable double m_ConstantTerm;
80  double m_ParameterScale;
81  double m_SmallDelta, m_BigDelta;
82 };
83 
84 } // end namespace anima
STL namespace.
const double DiffusionBigDelta
Default big delta value (classical values)
const double DiffusionSmallDelta
Default small delta value (classical values)
MultiCompartmentModel: holds several diffusion compartments, ordered by type It also handles weights ...