Go to the documentation of this file.00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028 #ifndef _MRPT_MONTE_CARLO_H_
00029 #define _MRPT_MONTE_CARLO_H_
00030
00031 #include <map>
00032 #include <cmath>
00033 #include <vector>
00034 #include <numeric>
00035 #include <algorithm>
00036 #include <stdexcept>
00037 #include <functional>
00038
00039 #include <mrpt/random.h>
00040 #include <mrpt/utils/CTicTac.h>
00041
00042 namespace mrpt { namespace math {
00043
00044
00045
00046
00047
00048
00049
00050
00051
00052
00053
00054
00055
00056
00057 template<typename T,typename NUM,typename OTHER> class CMonteCarlo {
00058 private:
00059 mrpt::random::CRandomGenerator gen;
00060 class CStatisticalAnalyzer {
00061 private:
00062 Eigen::Matrix<NUM,Eigen::Dynamic,1> data;
00063 public:
00064 template<typename VEC> inline CStatisticalAnalyzer(const VEC &v1):data(v1.begin(),v1.end()) {}
00065 template<typename VEC> inline void setData(const VEC &v1) {
00066 data.assign(v1.begin(),v1.end());
00067 }
00068 template<typename VEC> inline void getData(VEC &v1) const {
00069 v1.assign(data.begin(),data.end());
00070 }
00071 template<typename VEC1,typename VEC2> inline void getDistribution(VEC1 &vx,VEC2 &vy,const NUM width=1.0) const {
00072 mrpt::vector_double vvx,vvy;
00073 getDistribution(vvx,vvy,width);
00074 vx.assign(vvx.begin(),vvx.end());
00075 vy.assign(vvy.begin(),vvy.end());
00076 }
00077
00078 inline void getDistribution(mrpt::vector_double &vx,mrpt::vector_double &vy,const NUM width=1.0) const {
00079 CHistogram hist(CHistogram::createWithFixedWidth(0,*max_element(data.begin(),data.end()),width));
00080 hist.add(data);
00081 hist.getHistogram(vx,vy);
00082 }
00083
00084 };
00085 public:
00086
00087
00088 T (*valueGenerator)(mrpt::random::CRandomGenerator &);
00089
00090 NUM (*errorFun1)(const T &);
00091
00092 OTHER (*intermediateFun)(const T &);
00093 NUM (*errorFun2)(const T &,const OTHER &);
00094 inline CMonteCarlo():gen(),valueGenerator(NULL),errorFun1(NULL),intermediateFun(NULL),errorFun2(NULL) {}
00095 NUM doExperiment(size_t N,double &time,bool showInWindow=false) {
00096 if (!valueGenerator) throw std::logic_error("Value generator function is not set.");
00097 std::vector<T> baseData(N);
00098 std::vector<NUM> errorData(N);
00099 mrpt::utils::CTicTac meter;
00100 for (size_t i=0;i<N;++i) baseData[i]=valueGenerator(gen);
00101 if (errorFun1) {
00102 meter.Tic();
00103 std::transform(baseData.begin(),baseData.end(),errorData.begin(),errorFun1);
00104 time=meter.Tac();
00105 } else {
00106 if (!intermediateFun||!errorFun2) throw std::logic_error("Experiment-related functions are not set.");
00107 std::vector<OTHER> intermediate(N);
00108 transform(baseData.begin(),baseData.end(),intermediate.begin(),intermediateFun);
00109 meter.Tic();
00110 for (size_t i=0;i<N;++i) errorData[i]=errorFun2(baseData[i],intermediate[i]);
00111 time=meter.Tac();
00112 }
00113 NUM res=accumulate(errorData.begin(),errorData.end(),NUM(0))/errorData.size();
00114
00115
00116
00117
00118
00119
00120
00121
00122
00123
00124
00125 return res;
00126 }
00127 };
00128
00129 }}
00130 #endif