ITK  4.9.0
Insight Segmentation and Registration Toolkit
itkPowellOptimizer.h
Go to the documentation of this file.
1 /*=========================================================================
2  *
3  * Copyright Insight Software Consortium
4  *
5  * Licensed under the Apache License, Version 2.0 (the "License");
6  * you may not use this file except in compliance with the License.
7  * You may obtain a copy of the License at
8  *
9  * http://www.apache.org/licenses/LICENSE-2.0.txt
10  *
11  * Unless required by applicable law or agreed to in writing, software
12  * distributed under the License is distributed on an "AS IS" BASIS,
13  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14  * See the License for the specific language governing permissions and
15  * limitations under the License.
16  *
17  *=========================================================================*/
18 #ifndef itkPowellOptimizer_h
19 #define itkPowellOptimizer_h
20 
21 #include "itkVector.h"
22 #include "itkMatrix.h"
24 #include "ITKOptimizersExport.h"
25 
26 namespace itk
27 {
62 class ITKOptimizers_EXPORT PowellOptimizer:
64 {
65 public:
71 
74 
76  itkNewMacro(Self);
77 
80 
84 
86  itkSetMacro(Maximize, bool);
87  itkBooleanMacro(Maximize);
88  itkGetConstReferenceMacro(Maximize, bool);
90 
92  itkSetMacro(MaximumIteration, unsigned int);
93  itkGetConstReferenceMacro(MaximumIteration, unsigned int);
95 
97  itkSetMacro(MaximumLineIteration, unsigned int);
98  itkGetConstMacro(MaximumLineIteration, unsigned int);
100 
103  itkSetMacro(StepLength, double);
104  itkGetConstReferenceMacro(StepLength, double);
106 
109  itkSetMacro(StepTolerance, double);
110  itkGetConstReferenceMacro(StepTolerance, double);
112 
116  itkSetMacro(ValueTolerance, double);
117  itkGetConstReferenceMacro(ValueTolerance, double);
119 
121  itkGetConstReferenceMacro(CurrentCost, MeasureType);
122  MeasureType GetValue() const { return this->GetCurrentCost(); }
124 
126  itkGetConstReferenceMacro(CurrentIteration, unsigned int);
127 
129  itkGetConstReferenceMacro(CurrentLineIteration, unsigned int);
130 
132  virtual void StartOptimization() ITK_OVERRIDE;
133 
137  void StopOptimization()
138  { m_Stop = true; }
139 
140  itkGetConstReferenceMacro(CatchGetValueException, bool);
141  itkSetMacro(CatchGetValueException, bool);
142 
143  itkGetConstReferenceMacro(MetricWorstPossibleValue, double);
144  itkSetMacro(MetricWorstPossibleValue, double);
145 
146  virtual const std::string GetStopConditionDescription() const ITK_OVERRIDE;
147 
148 protected:
149  PowellOptimizer();
150  PowellOptimizer(const PowellOptimizer &);
151  virtual ~PowellOptimizer();
152  virtual void PrintSelf(std::ostream & os, Indent indent) const ITK_OVERRIDE;
153 
154  itkSetMacro(CurrentCost, double);
155 
158  void SetLine(const ParametersType & origin,
159  const vnl_vector< double > & direction);
160 
164  double GetLineValue(double x) const;
165 
166  double GetLineValue(double x, ParametersType & tempCoord) const;
167 
170  void SetCurrentLinePoint(double x, double fx);
171 
174  void Swap(double *a, double *b) const;
175 
178  void Shift(double *a, double *b, double *c, double d) const;
179 
189  virtual void LineBracket(double *ax, double *bx, double *cx,
190  double *fa, double *fb, double *fc);
191 
192  virtual void LineBracket(double *ax, double *bx, double *cx,
193  double *fa, double *fb, double *fc,
194  ParametersType & tempCoord);
195 
201  virtual void BracketedLineOptimize(double ax, double bx, double cx,
202  double fa, double fb, double fc,
203  double *extX, double *extVal);
204 
205  virtual void BracketedLineOptimize(double ax, double bx, double cx,
206  double fa, double fb, double fc,
207  double *extX, double *extVal,
208  ParametersType & tempCoord);
209 
210  itkGetMacro(SpaceDimension, unsigned int);
211  void SetSpaceDimension(unsigned int dim)
212  {
213  this->m_SpaceDimension = dim;
214  this->m_LineDirection.set_size(dim);
215  this->m_LineOrigin.set_size(dim);
216  this->m_CurrentPosition.set_size(dim);
217  this->Modified();
218  }
219 
220  itkSetMacro(CurrentIteration, unsigned int);
221 
222  itkGetMacro(Stop, bool);
223  itkSetMacro(Stop, bool);
224 
225 private:
226  unsigned int m_SpaceDimension;
227 
229  unsigned int m_CurrentIteration;
231 
233  unsigned int m_MaximumIteration;
235 
238 
241 
243  double m_StepLength;
245 
247  vnl_vector< double > m_LineDirection;
248 
250 
252  MeasureType m_CurrentCost;
253 
258  bool m_Stop;
259 
260  std::ostringstream m_StopConditionDescription;
261 }; // end of class
262 } // end of namespace itk
263 
264 #endif
unsigned int m_MaximumIteration
SmartPointer< const Self > ConstPointer
This class is a base for the CostFunctions returning a single value.
Light weight base class for most itk classes.
MeasureType GetValue() const
This class is a base for the Optimization methods that optimize a single valued function.
unsigned int m_MaximumLineIteration
CostFunctionType::Pointer CostFunctionPointer
SingleValuedCostFunction CostFunctionType
std::ostringstream m_StopConditionDescription
unsigned int m_CurrentLineIteration
vnl_vector< double > m_LineDirection
Implements Powell optimization using Brent line search.
SmartPointer< Self > Pointer
SingleValuedNonLinearOptimizer::ParametersType ParametersType
ParametersType m_LineOrigin
SingleValuedNonLinearOptimizer Superclass
Control indentation during Print() invocation.
Definition: itkIndent.h:49
unsigned int m_CurrentIteration