OTB  9.0.0
Orfeo Toolbox
otbMRFOptimizerMetropolis.h
Go to the documentation of this file.
1 /*
2  * Copyright (C) 2005-2022 Centre National d'Etudes Spatiales (CNES)
3  *
4  * This file is part of Orfeo Toolbox
5  *
6  * https://www.orfeo-toolbox.org/
7  *
8  * Licensed under the Apache License, Version 2.0 (the "License");
9  * you may not use this file except in compliance with the License.
10  * You may obtain a copy of the License at
11  *
12  * http://www.apache.org/licenses/LICENSE-2.0
13  *
14  * Unless required by applicable law or agreed to in writing, software
15  * distributed under the License is distributed on an "AS IS" BASIS,
16  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
17  * See the License for the specific language governing permissions and
18  * limitations under the License.
19  */
20 
21 #ifndef otbMRFOptimizerMetropolis_h
22 #define otbMRFOptimizerMetropolis_h
23 
24 #include "otbMRFOptimizer.h"
25 #include "otbMath.h"
26 #include "itkNumericTraits.h"
27 #include "itkMersenneTwisterRandomVariateGenerator.h"
28 
29 namespace otb
30 {
51 class ITK_EXPORT MRFOptimizerMetropolis : public MRFOptimizer
52 {
53 public:
56  typedef itk::SmartPointer<Self> Pointer;
57  typedef itk::SmartPointer<const Self> ConstPointer;
58  typedef Superclass::ParametersType ParametersType;
59 
60  typedef itk::Statistics::MersenneTwisterRandomVariateGenerator RandomGeneratorType;
61 
62  itkNewMacro(Self);
63 
64  itkTypeMacro(MRFOptimizerMetropolis, MRFOptimizer);
65 
67  void SetSingleParameter(double parameterVal)
68  {
69  this->m_Parameters.SetSize(1);
70  this->m_Parameters.Fill(parameterVal);
71  this->Modified();
72  }
74 
75  inline bool Compute(double deltaEnergy) override
76  {
77  if (deltaEnergy < 0)
78  {
79  return true;
80  }
81  if (deltaEnergy == 0)
82  {
83  return false;
84  }
85  else
86  {
87  double proba = std::exp(-(deltaEnergy) / this->m_Parameters[0]);
88  if ((m_Generator->GetIntegerVariate() % 10000) < proba * 10000)
89  {
90  return true;
91  }
92  }
93  return false;
94  }
95 
97  void InitializeSeed(int seed)
98  {
99  m_Generator->SetSeed(seed);
100  }
102  {
103  m_Generator->SetSeed();
104  }
106 
107 protected:
109  {
110  this->m_NumberOfParameters = 1;
111  this->m_Parameters.SetSize(1);
112  this->m_Parameters[0] = 1.0;
113  m_Generator = RandomGeneratorType::GetInstance();
114  m_Generator->SetSeed();
115  }
117  {
118  }
119  RandomGeneratorType::Pointer m_Generator;
120 };
121 }
122 
123 #endif
otb::MRFOptimizerMetropolis::MRFOptimizerMetropolis
MRFOptimizerMetropolis()
Definition: otbMRFOptimizerMetropolis.h:108
otb::MRFOptimizerMetropolis::m_Generator
RandomGeneratorType::Pointer m_Generator
Definition: otbMRFOptimizerMetropolis.h:119
otb::MRFOptimizerMetropolis::Self
MRFOptimizerMetropolis Self
Definition: otbMRFOptimizerMetropolis.h:54
otb::MRFOptimizerMetropolis::RandomGeneratorType
itk::Statistics::MersenneTwisterRandomVariateGenerator RandomGeneratorType
Definition: otbMRFOptimizerMetropolis.h:60
otb::MRFOptimizerMetropolis::Pointer
itk::SmartPointer< Self > Pointer
Definition: otbMRFOptimizerMetropolis.h:56
otb::MRFOptimizerMetropolis::ParametersType
Superclass::ParametersType ParametersType
Definition: otbMRFOptimizerMetropolis.h:58
otb::MRFOptimizerMetropolis::Superclass
MRFOptimizer Superclass
Definition: otbMRFOptimizerMetropolis.h:55
otb::MRFOptimizerMetropolis
This is the optimizer class implementing the Metropolis algorithm.
Definition: otbMRFOptimizerMetropolis.h:51
otbMath.h
otb
The "otb" namespace contains all Orfeo Toolbox (OTB) classes.
Definition: otbJoinContainer.h:32
otb::MRFOptimizerMetropolis::~MRFOptimizerMetropolis
~MRFOptimizerMetropolis() override
Definition: otbMRFOptimizerMetropolis.h:116
otb::MRFOptimizerMetropolis::InitializeSeed
void InitializeSeed(int seed)
Definition: otbMRFOptimizerMetropolis.h:97
otb::MRFOptimizerMetropolis::SetSingleParameter
void SetSingleParameter(double parameterVal)
Definition: otbMRFOptimizerMetropolis.h:67
otb::MRFOptimizerMetropolis::ConstPointer
itk::SmartPointer< const Self > ConstPointer
Definition: otbMRFOptimizerMetropolis.h:57
otb::MRFOptimizer
This is the base class for optimizer used in the MRF framework.
Definition: otbMRFOptimizer.h:43
otb::MRFOptimizerMetropolis::InitializeSeed
void InitializeSeed()
Definition: otbMRFOptimizerMetropolis.h:101
otb::MRFOptimizerMetropolis::Compute
bool Compute(double deltaEnergy) override
Definition: otbMRFOptimizerMetropolis.h:75
otbMRFOptimizer.h