go home Home | Main Page | Modules | Namespace List | Class Hierarchy | Alphabetical List | Data Structures | File List | Namespace Members | Data Fields | Globals | Related Pages
itkGenericConjugateGradientOptimizer.h
Go to the documentation of this file.
00001 /*======================================================================
00002 
00003   This file is part of the elastix software.
00004 
00005   Copyright (c) University Medical Center Utrecht. All rights reserved.
00006   See src/CopyrightElastix.txt or http://elastix.isi.uu.nl/legal.php for
00007   details.
00008 
00009      This software is distributed WITHOUT ANY WARRANTY; without even
00010      the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
00011      PURPOSE. See the above copyright notices for more information.
00012 
00013 ======================================================================*/
00014 
00015 
00016 #ifndef __itkGenericConjugateGradientOptimizer_h
00017 #define __itkGenericConjugateGradientOptimizer_h
00018 
00019 #include "itkScaledSingleValuedNonLinearOptimizer.h"
00020 #include "itkLineSearchOptimizer.h"
00021 #include <vector>
00022 #include <map>
00023 
00024 namespace itk
00025 {
00037   class GenericConjugateGradientOptimizer :
00038     public ScaledSingleValuedNonLinearOptimizer
00039   {
00040   public:
00041 
00042     typedef GenericConjugateGradientOptimizer     Self;
00043     typedef ScaledSingleValuedNonLinearOptimizer  Superclass;
00044     typedef SmartPointer<Self>                    Pointer;
00045     typedef SmartPointer<const Self>              ConstPointer;
00046 
00047     itkNewMacro(Self);
00048     itkTypeMacro(GenericConjugateGradientOptimizer,
00049       ScaledSingleValuedNonLinearOptimizer);
00050 
00051     typedef Superclass::ParametersType            ParametersType;
00052     typedef Superclass::DerivativeType            DerivativeType;
00053     typedef Superclass::CostFunctionType          CostFunctionType;
00054     typedef Superclass::ScaledCostFunctionType    ScaledCostFunctionType;
00055     typedef Superclass::MeasureType               MeasureType;
00056     typedef Superclass::ScalesType                ScalesType;
00057 
00058     typedef LineSearchOptimizer                   LineSearchOptimizerType;
00059     typedef LineSearchOptimizerType::Pointer      LineSearchOptimizerPointer;
00060 
00063     typedef double (Self::*                       ComputeBetaFunctionType
00064       )( const DerivativeType & ,
00065          const DerivativeType & ,
00066          const ParametersType & );
00067     typedef std::string                           BetaDefinitionType;
00068     typedef std::map< BetaDefinitionType,
00069       ComputeBetaFunctionType >                   BetaDefinitionMapType;
00070 
00071     typedef enum {
00072       MetricError,
00073       LineSearchError,
00074       MaximumNumberOfIterations,
00075       GradientMagnitudeTolerance,
00076       ValueTolerance,
00077       InfiniteBeta,
00078       Unknown }                                   StopConditionType;
00079 
00080     virtual void StartOptimization(void);
00081     virtual void ResumeOptimization(void);
00082     virtual void StopOptimization(void);
00083 
00085     itkGetConstMacro(CurrentIteration, unsigned long);
00086     itkGetConstMacro(CurrentValue, MeasureType);
00087     itkGetConstReferenceMacro(CurrentGradient, DerivativeType);
00088     itkGetConstMacro(InLineSearch, bool);
00089     itkGetConstReferenceMacro(StopCondition, StopConditionType);
00090     itkGetConstMacro(CurrentStepLength, double);
00091 
00093     itkSetObjectMacro(LineSearchOptimizer, LineSearchOptimizerType);
00094     itkGetObjectMacro(LineSearchOptimizer, LineSearchOptimizerType);
00095 
00097     itkGetConstMacro(MaximumNumberOfIterations, unsigned long);
00098     itkSetClampMacro(MaximumNumberOfIterations, unsigned long,
00099       1, NumericTraits<unsigned long>::max());
00100 
00107     itkGetConstMacro(GradientMagnitudeTolerance, double);
00108     itkSetMacro(GradientMagnitudeTolerance, double)
00109 
00110     
00117     itkGetConstMacro(ValueTolerance, double);
00118     itkSetMacro(ValueTolerance, double);
00119 
00123     virtual void SetMaxNrOfItWithoutImprovement(unsigned long arg);
00124     itkGetConstMacro(MaxNrOfItWithoutImprovement, unsigned long);
00125 
00127     virtual void SetBetaDefinition(const BetaDefinitionType & arg);
00128     itkGetConstReferenceMacro(BetaDefinition, BetaDefinitionType);
00129 
00130   protected:
00131     GenericConjugateGradientOptimizer();
00132     virtual ~GenericConjugateGradientOptimizer(){};
00133 
00134     void PrintSelf( std::ostream & os, Indent indent ) const;
00135 
00136     DerivativeType                m_CurrentGradient;
00137     MeasureType                   m_CurrentValue;
00138     unsigned long                 m_CurrentIteration;
00139     StopConditionType             m_StopCondition;
00140     bool                          m_Stop;
00141     double                        m_CurrentStepLength;
00142 
00145     bool                          m_UseDefaultMaxNrOfItWithoutImprovement;
00146 
00148     bool                          m_InLineSearch;
00149     itkSetMacro(InLineSearch, bool);
00150 
00155     bool                          m_PreviousGradientAndSearchDirValid;
00156 
00158     BetaDefinitionType            m_BetaDefinition;
00159 
00162     BetaDefinitionMapType         m_BetaDefinitionMap;
00163 
00169     virtual void AddBetaDefinition(
00170       const BetaDefinitionType & name,
00171       ComputeBetaFunctionType function);
00172 
00183     virtual void ComputeSearchDirection(
00184       const DerivativeType & previousGradient,
00185       const DerivativeType & gradient,
00186       ParametersType & searchDir);
00187 
00192     virtual void LineSearch(
00193       const ParametersType searchDir,
00194       double & step,
00195       ParametersType & x,
00196       MeasureType & f,
00197       DerivativeType & g );
00198 
00203     virtual bool TestConvergence(bool firstLineSearchDone);
00204 
00206     virtual double ComputeBeta(
00207       const DerivativeType & previousGradient,
00208       const DerivativeType & gradient,
00209       const ParametersType & previousSearchDir);
00210 
00214     double ComputeBetaSD(
00215       const DerivativeType & previousGradient,
00216       const DerivativeType & gradient,
00217       const ParametersType & previousSearchDir);
00219     double ComputeBetaFR(
00220       const DerivativeType & previousGradient,
00221       const DerivativeType & gradient,
00222       const ParametersType & previousSearchDir);
00224     double ComputeBetaPR(
00225       const DerivativeType & previousGradient,
00226       const DerivativeType & gradient,
00227       const ParametersType & previousSearchDir);
00229     double ComputeBetaDY(
00230       const DerivativeType & previousGradient,
00231       const DerivativeType & gradient,
00232       const ParametersType & previousSearchDir);
00234     double ComputeBetaHS(
00235       const DerivativeType & previousGradient,
00236       const DerivativeType & gradient,
00237       const ParametersType & previousSearchDir);
00239     double ComputeBetaDYHS(
00240       const DerivativeType & previousGradient,
00241       const DerivativeType & gradient,
00242       const ParametersType & previousSearchDir);
00243 
00244 
00245   private:
00246     GenericConjugateGradientOptimizer(const Self&); //purposely not implemented
00247     void operator=(const Self&); //purposely not implemented
00248 
00249     unsigned long                 m_MaximumNumberOfIterations;
00250     double                        m_ValueTolerance;
00251     double                        m_GradientMagnitudeTolerance;
00252     unsigned long                 m_MaxNrOfItWithoutImprovement;
00253 
00254     LineSearchOptimizerPointer    m_LineSearchOptimizer;
00255 
00256   }; // end class GenericConjugateGradientOptimizer
00257 
00258 
00259 } // end namespace itk
00260 
00261 
00262 #endif //#ifndef __itkGenericConjugateGradientOptimizer_h
00263 


Generated on 24-05-2012 for elastix by doxygen 1.7.6.1 elastix logo