Go to the documentation of this file.00001 #ifndef __BCENGINEMCMC__H
00002 #define __BCENGINEMCMC__H
00003
00023
00024
00025 #include <vector>
00026
00027
00028 class TH1D;
00029 class TH2D;
00030 class TTree;
00031 class TRandom3;
00032
00033
00034
00035 class BCEngineMCMC
00036 {
00037
00038 public:
00039
00044 enum Precision{ kLow, kMedium, kHigh, kVeryHigh };
00045
00052 BCEngineMCMC();
00053
00057 BCEngineMCMC(int n);
00058
00061 BCEngineMCMC(const BCEngineMCMC & enginemcmc);
00062
00065 virtual ~BCEngineMCMC();
00066
00073 BCEngineMCMC & operator = (const BCEngineMCMC & engineMCMC);
00074
00081 int MCMCGetNParameters()
00082 { return fMCMCNParameters; }
00083
00086 int MCMCGetNChains()
00087 { return fMCMCNChains; }
00088
00091 int MCMCGetNLag()
00092 { return fMCMCNLag; }
00093
00096 std::vector <int> MCMCGetNIterations()
00097 { return fMCMCNIterations; }
00098
00101 int MCMCGetCurrentIteration()
00102 { return fMCMCCurrentIteration; }
00103
00106 int MCMCGetCurrentChain()
00107 { return fMCMCCurrentChain; }
00108
00112 int MCMCGetNIterationsConvergenceGlobal()
00113 { return fMCMCNIterationsConvergenceGlobal; }
00114
00117 bool MCMCGetFlagConvergenceGlobal()
00118 { return fMCMCFlagConvergenceGlobal; }
00119
00122 int MCMCGetNIterationsMax()
00123 { return fMCMCNIterationsMax; }
00124
00127 int MCMCGetNIterationsRun()
00128 { return fMCMCNIterationsRun; }
00129
00132 int MCMCGetNIterationsPreRunMin()
00133 { return fMCMCNIterationsPreRunMin; }
00134
00137 int MCMCGetNIterationsUpdate()
00138 { return fMCMCNIterationsUpdate; }
00139
00142 int MCMCGetNIterationsUpdateMax()
00143 { return fMCMCNIterationsUpdateMax; }
00144
00147 std::vector <int> MCMCGetNTrialsTrue()
00148 { return fMCMCNTrialsTrue; }
00149
00152 std::vector <int> MCMCGetNTrialsFalse()
00153 { return fMCMCNTrialsFalse; }
00154
00158 std::vector <double> MCMCGetprobMean()
00159 { return fMCMCprobMean; }
00160
00164 std::vector <double> MCMCGetVariance()
00165 { return fMCMCprobVar; }
00166
00169 std::vector <double> MCMCGetTrialFunctionScaleFactor()
00170 { return fMCMCTrialFunctionScaleFactor; }
00171
00175 std::vector <double> MCMCGetTrialFunctionScaleFactor(int ichain);
00176
00181 double MCMCGetTrialFunctionScaleFactor(int ichain, int ipar);
00182
00185 std::vector <double> MCMCGetx()
00186 { return fMCMCx; }
00187
00191 std::vector <double> MCMCGetx(int ichain);
00192
00197 double MCMCGetx(int ichain, int ipar);
00198
00201 std::vector <double> MCMCGetLogProbx()
00202 { return fMCMCprob; }
00203
00207 double MCMCGetLogProbx(int ichain);
00208
00211 int MCMCGetPhase()
00212 { return fMCMCPhase; }
00213
00216 int MCMCGetCycle()
00217 { return fMCMCCycle; }
00218
00221 std::vector <double> MCMCGetMaximumPoints()
00222 { return fMCMCxMax; }
00223
00227 std::vector <double> MCMCGetMaximumPoint(int i);
00228
00231 std::vector <double> MCMCGetMaximumLogProb()
00232 { return fMCMCprobMax; }
00233
00236 int MCMCGetFlagInitialPosition()
00237 { return fMCMCFlagInitialPosition; }
00238
00241 double MCMCGetRValueCriterion()
00242 { return fMCMCRValueCriterion; }
00243
00246 double MCMCGetRValueParametersCriterion()
00247 { return fMCMCRValueParametersCriterion; }
00248
00251 double MCMCGetRValue()
00252 { return fMCMCRValue; }
00253
00257 double MCMCGetRValueParameters(int i)
00258 { return fMCMCRValueParameters.at(i); }
00259
00261 bool MCMCGetRValueStrict()
00262 { return fMCMCRValueUseStrict; }
00263
00266 bool MCMCGetFlagRun()
00267 { return fMCMCFlagRun; }
00268
00273 TTree * MCMCGetMarkovChainTree(int i)
00274 { return fMCMCTrees.at(i); }
00275
00280 TH1D * MCMCGetH1Marginalized(int i);
00281
00287 TH2D * MCMCGetH2Marginalized(int i, int j);
00288
00293 TRandom3 * MCMCGetTRandom3()
00294 { return fRandom; }
00295
00303 void MCMCSetTrialFunctionScaleFactor(std::vector <double> scale)
00304 { fMCMCTrialFunctionScaleFactorStart = scale; }
00305
00308 void MCMCSetNChains(int n);
00309
00312 void MCMCSetNLag(int n)
00313 { fMCMCNLag = n; }
00314
00317 void MCMCSetNIterationsMax(int n)
00318 { fMCMCNIterationsMax = n; }
00319
00322 void MCMCSetNIterationsRun(int n)
00323 { fMCMCNIterationsRun = n; }
00324
00327 void MCMCSetNIterationsPreRunMin(int n)
00328 { fMCMCNIterationsPreRunMin = n; }
00329
00334 void MCMCSetNIterationsUpdate(int n)
00335 { fMCMCNIterationsUpdate = n; }
00336
00342 void MCMCSetNIterationsUpdateMax(int n)
00343 { fMCMCNIterationsUpdateMax = n; }
00344
00347 void MCMCSetMinimumEfficiency(double efficiency)
00348 { fMCMCEfficiencyMin = efficiency; }
00349
00352 void MCMCSetMaximumEfficiency(double efficiency)
00353 { fMCMCEfficiencyMax = efficiency; }
00354
00357 void MCMCSetWriteChainToFile(bool flag)
00358 { fMCMCFlagWriteChainToFile = flag; }
00359
00362 void MCMCSetWritePreRunToFile(bool flag)
00363 { fMCMCFlagWritePreRunToFile = flag; }
00364
00368 void MCMCSetInitialPositions(std::vector<double> x0s);
00369
00373 void MCMCSetInitialPositions(std::vector< std::vector<double> > x0s);
00374
00377 void MCMCSetFlagInitialPosition(int flag)
00378 { fMCMCFlagInitialPosition = flag; }
00379
00383 void MCMCSetFlagOrderParameters(bool flag)
00384 { fMCMCFlagOrderParameters = flag; }
00385
00387 void MCMCSetFlagFillHistograms(bool flag);
00388
00390 void MCMCSetFlagFillHistograms(int index, bool flag);
00391
00393 void MCMCSetFlagPreRun(bool flag)
00394 { fMCMCFlagPreRun = flag; }
00395
00398 void MCMCSetRValueCriterion(double r)
00399 { fMCMCRValueCriterion = r; }
00400
00403 void MCMCSetRValueParametersCriterion(double r)
00404 { fMCMCRValueParametersCriterion = r; }
00405
00407 void MCMCSetRValueStrict(bool strict=true)
00408 { fMCMCRValueUseStrict = strict; }
00409
00412 void MCMCSetMarkovChainTrees(std::vector <TTree *> trees);
00413
00416 void MCMCInitializeMarkovChainTrees();
00417
00422 int SetMarginalized(int index, TH1D * h);
00423
00429 int SetMarginalized(int index1, int index2, TH2D * h);
00430
00433 void MCMCSetValuesDefault();
00434
00437 void MCMCSetValuesQuick();
00438
00441 void MCMCSetValuesDetail();
00442
00445 void MCMCSetPrecision(BCEngineMCMC::Precision precision);
00446
00447
00457 int MCMCAddParameter(double min, double max);
00458
00465 virtual void MCMCTrialFunction(int ichain, std::vector <double> &x);
00466
00474 virtual double MCMCTrialFunctionSingle(int ichain, int ipar);
00475
00481 bool MCMCGetProposalPointMetropolis(int chain, std::vector <double> &x);
00482
00488 bool MCMCGetProposalPointMetropolis(int chain, int parameter, std::vector <double> &x);
00489
00493 bool MCMCGetNewPointMetropolis(int chain = 0);
00494 bool MCMCGetNewPointMetropolis(int chain, int parameter);
00495
00498 void MCMCInChainCheckMaximum();
00499
00502 void MCMCInChainUpdateStatistics();
00503
00506 void MCMCInChainFillHistograms();
00507
00510 void MCMCInChainTestConvergenceAllChains();
00511
00532 double MCMCCalculateRValue(const std::vector<double> & chainMeans,
00533 const std::vector<double> & chainVariances, double mean, double B,
00534 double W, int nPoints, int iParam = -1);
00535
00538 void MCMCInChainWriteChains();
00539
00543 virtual double LogEval(std::vector <double> parameters);
00544
00547 int MCMCMetropolis();
00548
00551 int MCMCMetropolisPreRun();
00552
00555 void MCMCResetRunStatistics();
00556
00559 void MCMCInitializeMarkovChains();
00560
00564 int MCMCInitialize();
00565
00569 int MCMCResetResults();
00570
00577 virtual void MCMCIterationInterface()
00578 {}
00579
00588 virtual void MCMCCurrentPointInterface(std::vector <double> & point, int ichain, bool accepted)
00589 {
00590
00591
00592 (void)point;
00593 (void)ichain;
00594 (void)accepted;
00595 }
00596
00599 private:
00600
00603 typedef bool (BCEngineMCMC::*MCMCPointerToGetProposalPoint) (int chain, std::vector <double> xnew, std::vector <double> xold) const;
00604
00607 MCMCPointerToGetProposalPoint fMCMCPointerToGetProposalPoint;
00608
00609 protected:
00610
00613 int fMCMCNParameters;
00614
00617 std::vector <double> fMCMCBoundaryMin;
00618 std::vector <double> fMCMCBoundaryMax;
00619
00622 std::vector <bool> fMCMCFlagsFillHistograms;
00623
00626 int fMCMCNChains;
00627
00630 int fMCMCNLag;
00631
00635 std::vector<int> fMCMCNIterations;
00636
00640 int fMCMCCurrentIteration;
00641
00645 int fMCMCCurrentChain;
00646
00649 int fMCMCNIterationsUpdate;
00650
00653 int fMCMCNIterationsUpdateMax;
00654
00658 int fMCMCNIterationsConvergenceGlobal;
00659
00662 bool fMCMCFlagConvergenceGlobal;
00663
00666 int fMCMCNIterationsMax;
00667
00670 int fMCMCNIterationsRun;
00671
00674 int fMCMCNIterationsPreRunMin;
00675
00679 std::vector<int> fMCMCNTrialsTrue;
00680
00684 std::vector<int> fMCMCNTrialsFalse;
00685
00688 bool fMCMCFlagWriteChainToFile;
00689
00692 bool fMCMCFlagWritePreRunToFile;
00693
00697 std::vector <double> fMCMCTrialFunctionScaleFactor;
00698
00699
00702 std::vector <double> fMCMCTrialFunctionScaleFactorStart;
00703
00706 bool fMCMCFlagPreRun;
00707
00710 bool fMCMCFlagRun;
00711
00717 std::vector <double> fMCMCInitialPosition;
00718
00721 double fMCMCEfficiencyMin;
00722
00725 double fMCMCEfficiencyMax;
00726
00731 int fMCMCFlagInitialPosition;
00732
00736 bool fMCMCFlagOrderParameters;
00737
00740 bool fMCMCFlagFillHistograms;
00741
00746 int fMCMCPhase;
00747
00751 int fMCMCCycle;
00752
00758 std::vector <double> fMCMCx;
00759
00764 std::vector <double> fMCMCxMax;
00765
00769 std::vector <double> fMCMCxMean;
00770
00774 std::vector <double> fMCMCxVar;
00775
00778 std::vector <double> fMCMCxLocal;
00779
00783 std::vector<double> fMCMCprob;
00784
00788 std::vector <double> fMCMCprobMax;
00789
00793 std::vector <double> fMCMCprobMean;
00794
00798 std::vector <double> fMCMCprobVar;
00799
00806 bool fMCMCRValueUseStrict;
00807
00810 double fMCMCRValueCriterion;
00811
00814 double fMCMCRValueParametersCriterion;
00815
00818 double fMCMCRValue;
00819
00821 std::vector <double> fMCMCRValueParameters;
00822
00825 TRandom3 * fRandom;
00826
00829 std::vector<int> fMCMCH1NBins;
00830
00833 std::vector <TH1D *> fMCMCH1Marginalized;
00834 std::vector <TH2D *> fMCMCH2Marginalized;
00835
00839 std::vector<TTree *> fMCMCTrees;
00840 };
00841
00842
00843
00844 #endif