8 #include <itkSingleValuedCostFunction.h> 11 #include <itkCostFunction.h> 12 #include <itkNonLinearOptimizer.h> 20 template <
class InputPixelType,
class OutputPixelType>
92 itkSetMacro(SmallDelta,
double)
93 itkSetMacro(BigDelta,
double)
95 itkSetMacro(B0Threshold,
double)
99 itkSetMacro(AbsoluteCostChange,
double)
105 itkSetMacro(ModelWithFreeWaterComponent,
bool)
106 itkSetMacro(ModelWithStationaryWaterComponent,
bool)
107 itkSetMacro(ModelWithRestrictedWaterComponent,
bool)
108 itkSetMacro(ModelWithStaniszComponent,
bool)
114 itkSetMacro(NumberOfCoils,
unsigned int)
115 itkGetMacro(NumberOfCoils,
unsigned int)
116 itkSetMacro(NumberOfCompartments,
unsigned int)
117 itkSetMacro(FindOptimalNumberOfCompartments,
bool)
119 itkSetMacro(UseConstrainedDiffusivity,
bool)
120 itkSetMacro(UseConstrainedFreeWaterDiffusivity,
bool)
121 itkSetMacro(UseConstrainedIRWDiffusivity,
bool)
122 itkSetMacro(UseConstrainedStaniszDiffusivity,
bool)
123 itkSetMacro(UseConstrainedStaniszRadius,
bool)
125 itkSetMacro(UseConstrainedOrientationConcentration,
bool)
126 itkSetMacro(UseConstrainedExtraAxonalFraction,
bool)
127 itkSetMacro(UseCommonConcentrations,
bool)
128 itkSetMacro(UseCommonExtraAxonalFractions,
bool)
130 itkSetMacro(UseCommonDiffusivities,
bool)
135 itkGetMacro(SmallDelta,
double)
136 itkGetMacro(BigDelta,
double)
152 itkSetMacro(AxialDiffusivityValue,
double)
153 itkSetMacro(StaniszDiffusivityValue,
double)
154 itkSetMacro(IRWDiffusivityValue,
double)
155 itkSetMacro(RadialDiffusivity1Value,
double)
156 itkSetMacro(RadialDiffusivity2Value,
double)
158 itkSetMacro(XTolerance,
double)
159 itkSetMacro(FTolerance,
double)
160 itkSetMacro(MaxEval,
unsigned int)
167 m_SigmaSquareVolume = 0;
170 m_GradientStrengths.clear();
171 m_GradientDirections.clear();
173 m_NumberOfDictionaryEntries = 500;
174 m_Optimizer =
"bobyqa";
175 m_AbsoluteCostChange = 0.01;
179 m_ModelWithFreeWaterComponent =
true;
180 m_ModelWithStationaryWaterComponent =
true;
181 m_ModelWithRestrictedWaterComponent =
true;
182 m_ModelWithStaniszComponent =
true;
188 m_NumberOfCompartments = 2;
189 m_FindOptimalNumberOfCompartments =
true;
191 m_UseConstrainedDiffusivity =
false;
192 m_UseConstrainedFreeWaterDiffusivity =
true;
193 m_UseConstrainedIRWDiffusivity =
true;
194 m_UseConstrainedStaniszDiffusivity =
true;
195 m_UseConstrainedStaniszRadius =
true;
196 m_UseCommonDiffusivities =
false;
198 m_UseConstrainedOrientationConcentration =
false;
199 m_UseConstrainedExtraAxonalFraction =
false;
200 m_UseCommonConcentrations =
false;
201 m_UseCommonExtraAxonalFractions =
false;
203 m_AxialDiffusivityValue = 1.71e-3;
204 m_StaniszDiffusivityValue = 1.71e-3;
205 m_IRWDiffusivityValue = 7.5e-4;
206 m_RadialDiffusivity1Value = 1.9e-4;
207 m_RadialDiffusivity2Value = 1.5e-4;
209 m_NumberOfImages = 0;
210 m_ExternalMoseVolume =
false;
222 for (
unsigned int i = 0;i < m_MCMCreators.size();++i)
223 delete m_MCMCreators[i];
224 m_MCMCreators.clear();
241 double &aiccValue,
double &b0Value,
double &sigmaSqValue);
245 std::vector <double> &observedSignals, itk::ThreadIdType threadId,
246 double &aiccValue,
double &b0Value,
double &sigmaSqValue);
250 std::vector <double> &observedSignals, itk::ThreadIdType threadId,
251 double &aiccValue,
double &b0Value,
double &sigmaSqValue);
255 itk::ThreadIdType threadId,
double &aiccValue,
double &b0Value,
double &sigmaSqValue);
259 itk::Array<double> &upperBounds);
263 std::vector<double> &observedSignals, itk::ThreadIdType threadId);
295 void InitializeDictionary();
297 double m_SmallDelta, m_BigDelta;
298 std::vector <double> m_GradientStrengths;
299 std::vector< GradientType > m_GradientDirections;
306 std::vector <MCMCreatorType *> m_MCMCreators;
308 std::string m_Optimizer;
311 vnl_matrix <double> m_SparseSticksDictionary;
312 unsigned int m_NumberOfDictionaryEntries;
313 std::vector < std::vector <double> > m_DictionaryDirections;
315 double m_B0Threshold;
316 unsigned int m_NumberOfImages;
318 double m_AbsoluteCostChange;
321 bool m_ModelWithFreeWaterComponent, m_ModelWithStationaryWaterComponent, m_ModelWithRestrictedWaterComponent, m_ModelWithStaniszComponent;
326 unsigned int m_NumberOfCoils;
327 unsigned int m_NumberOfCompartments;
328 bool m_FindOptimalNumberOfCompartments;
330 bool m_UseConstrainedDiffusivity;
331 bool m_UseConstrainedFreeWaterDiffusivity;
332 bool m_UseConstrainedIRWDiffusivity;
333 bool m_UseConstrainedStaniszRadius;
334 bool m_UseConstrainedStaniszDiffusivity;
335 bool m_UseCommonDiffusivities;
337 bool m_UseConstrainedOrientationConcentration;
338 bool m_UseConstrainedExtraAxonalFraction;
339 bool m_UseCommonConcentrations;
340 bool m_UseCommonExtraAxonalFractions;
342 double m_AxialDiffusivityValue;
343 double m_IRWDiffusivityValue;
344 double m_StaniszDiffusivityValue;
345 double m_RadialDiffusivity1Value;
346 double m_RadialDiffusivity2Value;
348 bool m_ExternalMoseVolume;
350 unsigned int m_MaxEval;
355 std::vector < std::vector <double> > m_ValuesCoarseGrid;
void TensorCoarseGridInitialization(MCMPointer &mcmUpdateValue, CostFunctionBasePointer &cost, MCMType::ListType &workVec, ParametersType &p)
Coarse grid initialization of tensor model.
OptimizerType::ParametersType ParametersType
anima::MultiCompartmentModelCreator::CompartmentType CompartmentType
std::vector< GradientType > & GetGradientDirections()
std::vector< double > & GetGradientStrengths()
itk::SmartPointer< Self > Pointer
anima::BaseCompartment BaseCompartmentType
virtual void InitializeModelFromSimplifiedOne(MCMPointer &simplifiedModel, MCMPointer &complexModel)
Performs initialization from simplified model with the same number of compartments.
anima::MultiCompartmentModelCreator MCMCreatorType
MoseImageType * GetMoseVolume()
void EstimateFreeWaterModel(MCMPointer &mcmValue, std::vector< double > &observedSignals, itk::ThreadIdType threadId, double &aiccValue, double &b0Value, double &sigmaSqValue)
Specific method for N=0 compartments estimation (only free water)
double PerformSingleOptimization(ParametersType &p, CostFunctionBasePointer &cost, itk::Array< double > &lowerBounds, itk::Array< double > &upperBounds)
Performs an optimization of the supplied cost function and parameters using the specified optimizer(s...
itk::Image< unsigned char, 3 > MoseImageType
itk::NonLinearOptimizer OptimizerType
void SetGradientStrengths(std::vector< double > &mb)
virtual void SparseInitializeSticks(MCMPointer &complexModel, bool authorizeNegativeB0Value, std::vector< double > &observedSignals, itk::ThreadIdType threadId)
Performs initialization from single DTI.
BaseCompartment::ListType ListType
MoseImageType::Pointer MoseImagePointer
SignalNoiseType
Denotes noise type on input signal. Based on this, a different cost function should be created...
MCMEstimatorImageFilter< InputPixelType, OutputPixelType > Self
void ComputeExtraAxonalAndKappaCoarseGrids()
Computes extra axonal and kappa coarse grids (used for NODDI and DDI initialization) ...
Superclass::InputImageRegionType InputImageRegionType
void SetMoseVolume(MoseImageType *img)
const double DiffusionBigDelta
Default big delta value (classical values)
OutputScalarImageType * GetB0Volume()
void AddGradientDirection(unsigned int i, GradientType &grad)
const double DiffusionSmallDelta
Default small delta value (classical values)
CostFunctionBaseType::Pointer CostFunctionBasePointer
itk::Image< InputPixelType, 3 > TInputImage
double ComputeAICcValue(MCMPointer &mcmValue, double costValue)
Compute AICc value from a cost function value and model.
void GetProfiledInformation(CostFunctionBasePointer &cost, MCMPointer &mcm, double &b0Value, double &sigmaSqValue)
Utility function to get profiled data from the cost function into an MCM model.
OutputImageType::Pointer OutputImagePointer
itk::Image< OutputPixelType, 3 > OutputScalarImageType
DiffusionModelCompartmentType
MCMType::ListType MCMVectorType
void ModelEstimation(MCMPointer &mcmValue, bool authorizedNegativeB0Value, std::vector< double > &observedSignals, itk::ThreadIdType threadId, double &aiccValue, double &b0Value, double &sigmaSqValue)
Doing estimation, calling initialization procedure until ball and zeppelin, returns AICc value...
itk::CostFunction CostFunctionBaseType
void OptimizeNonIsotropicCompartments(MCMPointer &mcmValue, unsigned int currentNumberOfCompartments, std::vector< double > &observedSignals, itk::ThreadIdType threadId, double &aiccValue, double &b0Value, double &sigmaSqValue)
Doing estimation of non isotropic compartments (for a given number of anisotropic compartments) ...
void ComputeTensorRadialDiffsAndAzimuthCoarseGrids()
Computes extra axonal and kappa coarse grids (used for tensor final initaialization) ...
VectorImageType::Pointer VectorImagePointer
void GenerateOutputInformation() ITK_OVERRIDE
void WriteMCMOutput(std::string fileName)
OptimizerType::Pointer OptimizerPointer
anima::MCMImage< OutputPixelType, 3 > ModelImageType
itk::VectorImage< OutputPixelType, 3 > VectorImageType
itk::SmartPointer< const Self > ConstPointer
Superclass::OutputImageRegionType OutputImageRegionType
void ExtraAxonalAndKappaCoarseGridInitialization(MCMPointer &mcmUpdateValue, CostFunctionBasePointer &cost, MCMType::ListType &workVec, ParametersType &p)
Coarse grid initialization of NODDI and DDI models.
MultiCompartmentModel: holds several diffusion compartments, ordered by type It also handles weights ...
OutputScalarImageType * GetAICcVolume()
vnl_vector_fixed< double, 3 > GradientType
TInputImage InputImageType
anima::MaskedImageToImageFilter< TInputImage, TOutputImage > Superclass
virtual ~MCMEstimatorImageFilter()
OptimizerPointer CreateOptimizer(CostFunctionBasePointer &cost, itk::Array< double > &lowerBounds, itk::Array< double > &upperBounds)
Create an optimizer following the optimizer type and estimation mode.
MaximumLikelihoodEstimationMode
OutputImageType::PixelType VariableLengthVectorType
void InitialOrientationsEstimation(MCMPointer &mcmValue, bool authorizedNegativeB0Value, unsigned int currentNumberOfCompartments, std::vector< double > &observedSignals, itk::ThreadIdType threadId, double &aiccValue, double &b0Value, double &sigmaSqValue)
Doing estimation only of multiple orientations.
Superclass::MaskImageType MaskImageType
virtual MCMCreatorType * GetNewMCMCreatorInstance()
MCMType::Pointer MCMPointer
void SetOptimizer(std::string &opt)
anima::MCMImage< OutputPixelType, 3 > TOutputImage
MCMCreatorType::MCMType MCMType
itk::Image< InputPixelType, 4 > Image4DType
double GetCostValue(CostFunctionBasePointer &cost, ParametersType &p)
Utility function to get a value from a cost function.
virtual void BeforeThreadedGenerateData() ITK_OVERRIDE
void DynamicThreadedGenerateData(const OutputImageRegionType &outputRegionForThread) ITK_OVERRIDE
void CheckComputationMask() ITK_OVERRIDE
InputImageType::Pointer InputImagePointer
OutputScalarImageType::Pointer OutputScalarImagePointer
MCMEstimatorImageFilter()
MCMCreatorType::MCMPointer MCMPointer
OutputScalarImageType * GetSigmaSquareVolume()
MCMCreatorType * GetMCMCreator(unsigned int i)
std::string GetOptimizer()
Really this class is some simplified factory that creates the MCM that it knows.
virtual CostFunctionBasePointer CreateCostFunction(std::vector< double > &observedSignals, MCMPointer &mcmModel)
Create a cost function following the noise type and estimation mode.