Main MRPT website > C++ reference
MRPT logo
CMonteCarlo.h
Go to the documentation of this file.
1 /* +---------------------------------------------------------------------------+
2  | The Mobile Robot Programming Toolkit (MRPT) C++ library |
3  | |
4  | http://www.mrpt.org/ |
5  | |
6  | Copyright (C) 2005-2012 University of Malaga |
7  | |
8  | This software was written by the Machine Perception and Intelligent |
9  | Robotics Lab, University of Malaga (Spain). |
10  | Contact: Jose-Luis Blanco <jlblanco@ctima.uma.es> |
11  | |
12  | This file is part of the MRPT project. |
13  | |
14  | MRPT is free software: you can redistribute it and/or modify |
15  | it under the terms of the GNU General Public License as published by |
16  | the Free Software Foundation, either version 3 of the License, or |
17  | (at your option) any later version. |
18  | |
19  | MRPT is distributed in the hope that it will be useful, |
20  | but WITHOUT ANY WARRANTY; without even the implied warranty of |
21  | MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the |
22  | GNU General Public License for more details. |
23  | |
24  | You should have received a copy of the GNU General Public License |
25  | along with MRPT. If not, see <http://www.gnu.org/licenses/>. |
26  | |
27  +---------------------------------------------------------------------------+ */
28 #ifndef _MRPT_MONTE_CARLO_H_
29 #define _MRPT_MONTE_CARLO_H_
30 
31 #include <map>
32 #include <cmath>
33 #include <vector>
34 #include <numeric>
35 #include <algorithm>
36 #include <stdexcept>
37 #include <functional>
38 //#include <mrpt/gui/CDisplayWindowPlots.h>
39 #include <mrpt/random.h>
40 #include <mrpt/utils/CTicTac.h>
41 
42 namespace mrpt { namespace math {
43 
44  /** Montecarlo simulation for experiments in 1D.
45  Template arguments are:
46  - T: base type, i.e., if an experiment needs to generate random points, then T may be a TPoint3D, and so on.
47  - NUM: the numeric type used to represent the error. Usually, double.
48  - OTHER: an intermediate type, used especially when testing inverse functions. Leave as int or double if you don't use it.
49 
50  HOW TO USE THIS CLASS:
51  - Create an instance of the class.
52  - Refill the "valueGenerator" member with an appropriate function.
53  - If your experiment calculates the error directly from the base value, then refill the "errorFun1" member.
54  - Otherwise, if your experiment involves the calculation of some value whom with the experimental function is compared, refill "intermediateFun" and "errorFun2".
55  - Refill only on of the alternatives.
56  * \ingroup mrpt_base_grp
57  */
58  template<typename T,typename NUM,typename OTHER> class CMonteCarlo {
59  private:
60  mrpt::random::CRandomGenerator gen;
61  class CStatisticalAnalyzer {
62  private:
63  Eigen::Matrix<NUM,Eigen::Dynamic,1> data;
64  public:
65  template<typename VEC> inline CStatisticalAnalyzer(const VEC &v1):data(v1.begin(),v1.end()) {}
66  template<typename VEC> inline void setData(const VEC &v1) {
67  data.assign(v1.begin(),v1.end());
68  }
69  template<typename VEC> inline void getData(VEC &v1) const {
70  v1.assign(data.begin(),data.end());
71  }
72  template<typename VEC1,typename VEC2> inline void getDistribution(VEC1 &vx,VEC2 &vy,const NUM width=1.0) const {
73  mrpt::vector_double vvx,vvy;
74  getDistribution(vvx,vvy,width);
75  vx.assign(vvx.begin(),vvx.end());
76  vy.assign(vvy.begin(),vvy.end());
77  }
78  // Function overload, not specialization (GCC complains otherwise):
79  inline void getDistribution(mrpt::vector_double &vx,mrpt::vector_double &vy,const NUM width=1.0) const {
80  CHistogram hist(CHistogram::createWithFixedWidth(0,*max_element(data.begin(),data.end()),width));
81  hist.add(data);
82  hist.getHistogram(vx,vy);
83  }
84 
85  };
86  public:
87  //TODO: use templates for function types.
88  //Random generator.
89  T (*valueGenerator)(mrpt::random::CRandomGenerator &);
90  //Computes automatically the error (without an intermediate type)
91  NUM (*errorFun1)(const T &);
92 
93  OTHER (*intermediateFun)(const T &);
94  NUM (*errorFun2)(const T &,const OTHER &);
95  inline CMonteCarlo():gen(),valueGenerator(NULL),errorFun1(NULL),intermediateFun(NULL),errorFun2(NULL) {}
96  NUM doExperiment(size_t N,double &time,bool showInWindow=false) {
97  if (!valueGenerator) throw std::logic_error("Value generator function is not set.");
98  std::vector<T> baseData(N);
99  std::vector<NUM> errorData(N);
100  mrpt::utils::CTicTac meter;
101  for (size_t i=0;i<N;++i) baseData[i]=valueGenerator(gen);
102  if (errorFun1) {
103  meter.Tic();
104  std::transform(baseData.begin(),baseData.end(),errorData.begin(),errorFun1);
105  time=meter.Tac();
106  } else {
107  if (!intermediateFun||!errorFun2) throw std::logic_error("Experiment-related functions are not set.");
108  std::vector<OTHER> intermediate(N);
109  transform(baseData.begin(),baseData.end(),intermediate.begin(),intermediateFun);
110  meter.Tic();
111  for (size_t i=0;i<N;++i) errorData[i]=errorFun2(baseData[i],intermediate[i]);
112  time=meter.Tac();
113  }
114  NUM res=accumulate(errorData.begin(),errorData.end(),NUM(0))/errorData.size();
115  //if (showInWindow) {
116  // CStatisticalAnalyzer st(errorData);
117  // mrpt::gui::CDisplayWindowPlots wnd("Error results from Monte Carlo simulation");
118  // std::vector<NUM> errorX,errorY;
119  // st.getDistribution(errorX,errorY,0.1);
120  // wnd.plot(errorX,errorY,"b-","Plot1");
121  // NUM maxVal=*std::max_element(errorY.begin(),errorY.end());
122  // std::vector<NUM> dx(2,res),dy(mrpt::utils::make_vector<2,NUM>(0,maxVal));
123  // wnd.plot(dx,dy,"r-","Plot2");
124  // while (wnd.isOpen());
125  //}
126  return res;
127  }
128  };
129 
130 }} //End of namespaces
131 #endif


Page generated by Doxygen 1.8.3 for MRPT 0.9.6 SVN: at Fri Feb 15 22:05:02 EST 2013