a Code for the Combination of Indirect and Direct Constraints on High Energy Physics Models Logo
MonteCarlo.cpp
Go to the documentation of this file.
1 /*
2  * Copyright (C) 2012 HEPfit Collaboration
3  *
4  *
5  * For the licensing terms see doc/COPYING.
6  */
7 
8 #include "MonteCarlo.h"
9 #include "ModelFactory.h"
10 #include "ThObsFactory.h"
11 #include <BAT/BCParameter.h>
12 #include <TSystem.h>
13 #include <BAT/BCAux.h>
14 #include <BAT/BCLog.h>
15 #ifdef _MPI
16 #include <mpi.h>
17 #endif
18 #include <fstream>
19 #include <sstream>
20 #include <stdexcept>
21 #include <ctime>
22 #include <iostream>
23 #include <iomanip>
24 
25 
26 MonteCarlo::MonteCarlo(
27  ModelFactory& ModelF, ThObsFactory& ThObsF,
28  const std::string& ModelConf_i,
29  const std::string& MonteCarloConf_i,
30  const std::string& OutFile_i,
31  const std::string& JobTag_i)
32 : myInputParser(ModelF, ThObsF), MCEngine(ModPars, Obs, Obs2D, CGO, CGP) {
33  ModelConf = ModelConf_i;
34  MCMCConf = MonteCarloConf_i;
35  JobTag = JobTag_i;
36  if (OutFile_i.compare("") == 0) OutFile = "MCout" + JobTag + ".root";
37  else OutFile = OutFile_i + JobTag + ".root";
38  ObsDirName = "Observables" + JobTag;
39  FindModeWithMinuit = false;
40  RunMinuitOnly = false;
41  CalculateNormalization = "false";
42  NIterationNormalizationMC = 0;
43  PrintAllMarginalized = false;
44  PrintCorrelationMatrix = false;
45  PrintKnowledgeUpdatePlots = false;
46  PrintParameterPlot = false;
47  PrintTrianglePlot = false;
48  WritePreRunData = false;
49  checkrun = false;
50  writechains = false;
51 }
52 
53 MonteCarlo::~MonteCarlo()
54 {
55  Obs.clear();
56  boost::ptr_vector<Observable>().swap(Obs);
57 }
58 
59 void MonteCarlo::TestRun(int rank) {
60  if (checkrun == true) {
61  if (rank == 0) throw std::runtime_error("ERROR: MonteCarlo::TestRun() cannot be called after calling MonteCarlo::Run().\n");
62  else sleep (2);
63  } else {
64  checkrun = true;
65  }
66 
67  if (rank == 0) {
68  ModelName = myInputParser.ReadParameters(ModelConf, rank, ModPars, Obs, Obs2D, CGO, CGP);
69  std::map<std::string, double> DP;
70  if (Obs.size() == 0 && CGO.size() == 0) throw std::runtime_error("\nMonteCarlo::TestRun(): No observables or correlated Gaussian observables defined in SomeModel.conf file\n");
71 
72  for (std::vector<ModelParameter>::iterator it = ModPars.begin(); it < ModPars.end(); it++) {
73  if (it->IsCorrelated()) {
74  for (unsigned int i = 0; i < CGP.size(); i++) {
75  if (CGP[i].getName().compare(it->getCgp_name()) == 0) {
76  std::string index = it->getname().substr(CGP[i].getName().size());
77  long int lindex = strtol(index.c_str(), NULL, 10);
78  if (lindex > 0)
79  DP[CGP[i].getPar(lindex - 1).getname()] = CGP[i].getPar(lindex - 1).getave();
80  else {
81  std::stringstream out;
82  out << it->getname();
83  throw std::runtime_error("MonteCarlo::Run(): " + out.str() + "seems to be part of a CorrelatedGaussianParameters object, but I couldn't find the corresponding object");
84  }
85  }
86  }
87  } else
88  DP[it->getname()] = it->getave();
89  }
90 
91  if (!myInputParser.getModel()->Init(DP)) {
92  if (myInputParser.getModel()->getmissingModelParameters().size() > 0) {
93  std::cout << "\nPlease set the following parameters in the model configuration files:\n" << std::endl;
94  for (std::vector<std::string>::iterator it = myInputParser.getModel()->getmissingModelParameters().begin(); it != myInputParser.getModel()->getmissingModelParameters().end(); it++) {
95  std::cout << "ModelParameter\t" << *it << std::endl;
96  }
97  std::cout << std::endl;
98  myInputParser.getModel()->getmissingModelParameters().clear();
99  }
100  if (rank == 0) throw std::runtime_error("ERROR: Parameter(s) missing in model initialization. \n");
101  else sleep (2);
102  }
103  std::vector<std::string> unknownParameters = myInputParser.getModel()->getUnknownParameters();
104  if (unknownParameters.size() > 0 && rank == 0) {
105  std::cout << "\n" << std::endl;
106  for (std::vector<std::string>::iterator it = unknownParameters.begin(); it != unknownParameters.end(); it++)
107  std::cout << "WARNING: unknown parameter " << *it << " not added to MCMC Test Run" << std::endl;
108  }
109  if (Obs.size() > 0) std::cout << "\nObservables: \n" << std::endl;
110  for (boost::ptr_vector<Observable>::iterator it = Obs.begin(); it < Obs.end(); it++) {
111  double th = it->computeTheoryValue();
112  std::cout << it->getName() << " = " << th << std::endl;
113  }
114 
115  if (CGO.size() > 0) std::cout << "\nCorrelated Gaussian Observables: \n" << std::endl;
116  for (std::vector<CorrelatedGaussianObservables>::iterator it1 = CGO.begin(); it1 < CGO.end(); it1++) {
117  std::cout << it1->getName() << ":" << std::endl;
118  std::vector<Observable> ObsInCGO = it1->getObs();
119  for (std::vector<Observable>::iterator it2 = ObsInCGO.begin(); it2 < ObsInCGO.end(); it2++) {
120  double th = it2->computeTheoryValue();
121  std::cout << it2->getName() << " = " << th << std::endl;
122  }
123  std::cout << std::endl;
124  }
125  }
126 }
127 
128 void MonteCarlo::Run(const int rank) {
129  if (checkrun == true) {
130  if (rank == 0) throw std::runtime_error("ERROR: MonteCarlo::Run() cannot be called after calling MonteCarlo::TestRun().\n");
131  else sleep(2);
132  } else {
133  checkrun = true;
134  }
135 
136  try {
137 
138  /* set model parameters */
139  ModelName = myInputParser.ReadParameters(ModelConf, rank, ModPars, Obs, Obs2D, CGO, CGP);
140  int buffsize = 0;
141  std::map<std::string, double> DP;
142  for (std::vector<ModelParameter>::iterator it = ModPars.begin(); it < ModPars.end(); it++) {
143  if (it->geterrg() > 0. || it->geterrf() > 0.)
144  buffsize++;
145  if (it->IsCorrelated()) {
146  for (unsigned int i = 0; i < CGP.size(); i++) {
147  if (CGP[i].getName().compare(it->getCgp_name()) == 0) {
148  std::string index = it->getname().substr(CGP[i].getName().size());
149  long int lindex = strtol(index.c_str(), NULL, 10);
150  if (lindex > 0)
151  DP[CGP[i].getPar(lindex - 1).getname()] = CGP[i].getPar(lindex - 1).getave();
152  else {
153  std::stringstream out;
154  out << it->getname();
155  throw std::runtime_error("MonteCarlo::Run(): " + out.str() + "seems to be part of a CorrelatedGaussianParameters object, but I couldn't find the corresponding object");
156  }
157  }
158  }
159  } else
160  DP[it->getname()] = it->getave();
161  }
162  if (buffsize == 0){
163  if (rank == 0) throw std::runtime_error("No parameters being varied. Aborting MCMC run.\n");
164  else sleep(2);
165  }
166  buffsize++;
167  if (!myInputParser.getModel()->Init(DP)) {
168  if (myInputParser.getModel()->getmissingModelParameters().size() > 0) {
169  if (rank == 0) std::cout << "\nPlease set the following parameters in the model configuration files:\n" << std::endl;
170  for (std::vector<std::string>::iterator it = myInputParser.getModel()->getmissingModelParameters().begin(); it != myInputParser.getModel()->getmissingModelParameters().end(); it++) {
171  if (rank == 0) std::cout << "ModelParameter\t" << *it << std::endl;
172  }
173  std::cout << std::endl;
174  myInputParser.getModel()->getmissingModelParameters().clear();
175  }
176 
177  if (rank == 0) throw std::runtime_error("\nERROR: " + ModelName + " cannot be initialized.\n");
178  else sleep(2);
179  }
180  int discardbuffer = 0;
181  std::vector<std::string> unknownParameters = myInputParser.getModel()->getUnknownParameters();
182  for (std::vector<ModelParameter>::iterator it = ModPars.begin(); it < ModPars.end(); it++) {
183  if (std::find(unknownParameters.begin(), unknownParameters.end(), it->getname()) != unknownParameters.end()) {
184  if (it->geterrg() > 0. || it->geterrf() > 0.) {
185  discardbuffer++;
186  if (rank == 0) std::cout << "\nWARNING: " << it->getname() << " has a flat or Gaussian error in the configuration file even though it is an unknown parameter!!!\n";
187  }
188  }
189  }
190  buffsize = buffsize - discardbuffer;
191  if (buffsize == 1){
192  if (rank == 0) throw std::runtime_error("No relevant parameters being varied. Aborting MCMC run.\n");
193  else sleep(2);
194  }
195 
196  /* create a directory for outputs */
197  if (rank == 0) {
198  FileStat_t info;
199  if (gSystem->GetPathInfo(ObsDirName.c_str(), info) != 0) {
200  if (gSystem->MakeDirectory(ObsDirName.c_str()) == 0) std::cout << ObsDirName << " directory has been created." << std::endl;
201  else throw std::runtime_error("ERROR: " + ObsDirName + " directory cannot be created.\n");
202  }
203  }
204 
205  MCEngine.SetName(ModelName.c_str());
206  MCEngine.Initialize(myInputParser.getModel());
207  MCEngine.SetInitialPositionAttemptLimit(10000);
208 
209 #ifdef _MPI
210  double *recvbuff = new double[buffsize];
211 #endif
212 
213  if (rank != 0) {
214 #ifdef _MPI
215  double **sendbuff = new double*[1];
216  sendbuff[0] = new double[1];
217 
218  std::vector<double> pars;
219  int obsbuffsize = Obs.size() + 2 * Obs2D.size();
220  for (std::vector<CorrelatedGaussianObservables>::iterator it1 = CGO.begin(); it1 < CGO.end(); it1++)
221  obsbuffsize += it1->getObs().size();
222 
223  while (true) {
224  MPI_Scatter(sendbuff[0], buffsize, MPI_DOUBLE, recvbuff, buffsize, MPI_DOUBLE, 0, MPI_COMM_WORLD);
225 
226  if (recvbuff[0] == 0.) { // do nothing and return ll
227  double ll = log(0.);
228  MPI_Gather(&ll, 1, MPI_DOUBLE, sendbuff[0], 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);
229  } else if (recvbuff[0] == 1.) { //compute log likelihood
230  pars.assign(recvbuff + 1, recvbuff + buffsize);
231  double ll = MCEngine.LogEval(pars);
232  MPI_Gather(&ll, 1, MPI_DOUBLE, sendbuff[0], 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);
233  } else if (recvbuff[0] == 2.) { // compute observables
234  double sbuff[obsbuffsize];
235  std::map<std::string, double> DPars;
236  pars.assign(recvbuff + 1, recvbuff + buffsize);
237  MCEngine.setDParsFromParameters(pars,DPars);
238  myInputParser.getModel()->Update(DPars);
239 
240  int k = 0;
241  for (boost::ptr_vector<Observable>::iterator it = Obs.begin(); it < Obs.end(); it++) {
242  sbuff[k++] = it->computeTheoryValue();
243  }
244  for (std::vector<Observable2D>::iterator it = Obs2D.begin(); it < Obs2D.end(); it++) {
245  sbuff[k++] = it->computeTheoryValue();
246  sbuff[k++] = it->computeTheoryValue2();
247  }
248 
249  for (std::vector<CorrelatedGaussianObservables>::iterator it1 = CGO.begin(); it1 < CGO.end(); it1++) {
250  std::vector<Observable> pino(it1->getObs());
251  for (std::vector<Observable>::iterator it = pino.begin(); it != pino.end(); ++it)
252  sbuff[k++] = it->computeTheoryValue();
253  }
254  MPI_Gather(sbuff, obsbuffsize, MPI_DOUBLE, sendbuff[0], obsbuffsize, MPI_DOUBLE, 0, MPI_COMM_WORLD);
255  } else if (recvbuff[0] == 3.) { // do not compute observables, but gather the buffer
256  double sbuff[obsbuffsize];
257  MPI_Gather(sbuff, obsbuffsize, MPI_DOUBLE, sendbuff[0], obsbuffsize, MPI_DOUBLE, 0, MPI_COMM_WORLD);
258  } else if (recvbuff[0] == -1.)
259  break;
260  else {
261  std::cout << "recvbuff = " << recvbuff[0] << " rank " << rank << std::endl;
262  throw std::runtime_error("MonteCarlo::Run(): error in MPI message!\n");
263  }
264  }
265  delete sendbuff[0];
266  delete [] sendbuff;
267 #endif
268  } else {
269  ParseMCMCConfig(MCMCConf);
270  if (!RunMinuitOnly) std::cout << std::endl << "\nRunning in MonteCarlo mode...\n" << std::endl;
271  else std::cout << std::endl << "\nRunning Minuit Minimizer..." << std::endl;
272  std::cout << std::endl;
273  if (ModPars.size() > 0) std::cout << ModPars.size() - unknownParameters.size() << " parameters defined." << std::endl;
274  if (Obs.size() > 0) std::cout << Obs.size() << " observables defined." << std::endl;
275  if (Obs2D.size() > 0) std::cout << Obs2D.size() << " 2D observables defined." << std::endl;
276  if (CGO.size() > 0) std::cout << CGO.size() << " correlated gaussian observables defined:" << std::endl;
277  for (std::vector<CorrelatedGaussianObservables>::iterator it1 = CGO.begin();
278  it1 != CGO.end(); ++it1)
279  std::cout << " " << it1->getName() << " containing "
280  << it1->getObs().size() << " observables." << std::endl;
281  std::cout << std::endl;
282 
283  if (RunMinuitOnly) {
284  BCLog::OpenLog(("MinuitMinimizationResults" + JobTag + ".txt").c_str(), BCLog::debug, BCLog::debug);
285  RedirectHandle_t * rHandle = NULL;
286  gSystem->RedirectOutput(("MinuitMinimizationResults" + JobTag + ".txt").c_str(), "a", rHandle);
287  std::time_t ti = std::time(NULL);
288 #if __GNUC__ >= 5 || defined __clang__
289  std::cout << "\nMinuit Minimization started at " << std::put_time(std::localtime(&ti), "%c %Z") << " (" << ti << "s since the Epoch)" << "\n" << std::endl;
290 #else
291  char mitime[128];
292  strftime(mitime,sizeof(mitime),"%c %Z", std::localtime(&ti));
293  std::cout << "\nMinuit Minimization started at " << mitime << " (" << ti << "s since the Epoch)" << "\n" << std::endl;
294 #endif
295  std::cout << " " <<std::endl;
296  std::cout << "\t\t*** Minimizer Options ***" << std::endl;
297  MCEngine.GetMinuit().Options().Print();
298  std::cout << " " <<std::endl;
299  MCEngine.FindMode();
300  std::time_t tf = std::time(NULL);
301 #if __GNUC__ >= 5 || defined __clang__
302  std::cout << "\nMinuit Minimization ended at " << std::put_time(std::localtime(&tf), "%c %Z") << " (" << tf << "s since the Epoch)" << "\n" << std::endl;
303 #else
304  char mftime[128];
305  strftime(mitime,sizeof(mftime),"%c %Z", std::localtime(&tf));
306  std::cout << "\nMinuit Minimization ended at " << mftime << " (" << tf << "s since the Epoch)" << "\n" << std::endl;
307 #endif
308  gSystem->RedirectOutput(0);
309  BCLog::OutSummary(("Minuit results printed in MinuitMinimizationResults" + JobTag + ".txt").c_str());
310  BCLog::SetPrefix(false);
311  BCLog::OutSummary("");
312  BCLog::SetPrefix(true);
313  BCLog::CloseLog();
314  MCEngine.GetMinuit().PrintResults();
315  std::cout << std::endl;
316  return;
317  }
318 
319  MCEngine.CreateHistogramMaps();
320 
321  /* Open root file for storing data. */
322  if (writechains) {
323  MCEngine.WriteMarkovChain(OutFile, "RECREATE", true, false); /*Run: true, PreRun: false*/
324  MCEngine.AddChains();
325  } else {
326  MCEngine.WriteMarkovChain(OutFile, "RECREATE", false, false);
327  MCEngine.InitializeMarkovChainTree();
328  MCEngine.WriteMarkovChainRun(false);
329  MCEngine.WriteMarkovChainPreRun(false);
330  if (MCEngine.getchainedObsSize() > 0 || MCEngine.getWriteLogLikelihoodChain() || MCEngine.getWriteParametersChain()) MCEngine.AddChains();
331  }
332 
333  std::time_t ti = std::time(NULL);
334 #if __GNUC__ >= 5 || defined __clang__
335  std::cout << "\nMCMC run started at " << std::put_time(std::localtime(&ti), "%c %Z") << " (" << ti << "s since the Epoch)" << "\n" << std::endl;
336 #else
337  char mitime[128];
338  strftime(mitime,sizeof(mitime),"%c %Z", std::localtime(&ti));
339  std::cout << "\nMCMC run started at " << mitime << " (" << ti << "s since the Epoch)" << "\n" << std::endl;
340 #endif
341  // open log file
342  BCLog::OpenLog(("log" + JobTag + ".txt").c_str(), BCLog::debug, BCLog::debug);
343  // run the MCMC and marginalize w.r.t. to all parameters
344  MCEngine.BCIntegrate::SetNbins(NBINSMODELPARS);
345  MCEngine.SetMarginalizationMethod(BCIntegrate::kMargMetropolis);
346  MCEngine.MarginalizeAll();
347  std::time_t tf = std::time(NULL);
348 #if __GNUC__ >= 5 || defined __clang__
349  std::cout << "\nMCMC run ended at " << std::put_time(std::localtime(&tf), "%c %Z") << " (" << tf << "s since the Epoch)" << std::endl;
350 #else
351  char mftime[128];
352  strftime(mftime,sizeof(mftime),"%c %Z", std::localtime(&tf));
353  std::cout << "\nMCMC run ended at " << mftime << " (" << tf << "s since the Epoch)" << "\n" << std::endl;
354 #endif
355  int hour = (tf-ti)/3600;
356  int min = ((tf-ti)%3600) / 60;
357  int sec = ((tf-ti)%3600) % 60;
358  std::cout << "Time taken by the MCMC run: ";
359  std::cout << std::setfill('0') << std::setw(2) << hour << ":";
360  std::cout << std::setfill('0') << std::setw(2) << min << ":" ;
361  std::cout << std::setfill('0') << std::setw(2) << sec << "\n" << std::endl;
362 
363  // find mode using the best fit parameters as start values
364  if (FindModeWithMinuit) MCEngine.FindMode(MCEngine.GetBestFitParameters());
365 
366  // draw all marginalized distributions into a pdf file
367  if (PrintAllMarginalized) MCEngine.PrintAllMarginalized(("MonteCarlo_plots" + JobTag + ".pdf").c_str());
368 
369  // print histograms
370  MCEngine.PrintHistogram(OutFile, ObsDirName);
371 
372  // draw the correlation matrix into a pdf file and make the tex files
373  if (PrintCorrelationMatrix) {
374  MCEngine.PrintCorrelationMatrix(("CorrelationMatrix" + JobTag + ".pdf").c_str());
375  MCEngine.PrintCorrelationMatrix(("CorrelationMatrix" + JobTag + ".tex").c_str());
376  MCEngine.PrintCorrelationMatrixToLaTeX(("ParameterCorrelations" + JobTag + ".tex").c_str());
377  }
378 
379  // print comparisons of the prior knowledge to the posterior knowledge
380  // for all parameters into a pdf file
381  if (PrintKnowledgeUpdatePlots) MCEngine.PrintKnowledgeUpdatePlots(("ParamUpdate" + JobTag + ".pdf").c_str());
382 
383  // draw an overview plot of the parameters into a pdf file
384  if (PrintParameterPlot) MCEngine.PrintParameterPlot(("ParamSummary" + JobTag + ".pdf").c_str());
385 
386  MCEngine.WriteMarginalizedDistributions(OutFile, "UPDATE");
387 
388  // print logs for the histograms of the observables into a text file
389  std::ofstream outHistoLog;
390  outHistoLog.open((ObsDirName + "/HistoLog" + JobTag + ".txt").c_str(), std::ios::out);
391  outHistoLog << MCEngine.getHistoLog();
392  outHistoLog.close();
393 
394  // print statistics for the theory values of the observables into a text file
395  std::ofstream outStatLog;
396  outStatLog.open((ObsDirName + "/Statistics" + JobTag + ".txt").c_str(), std::ios::out);
397  outStatLog << MCEngine.computeStatistics();
398  outStatLog.close();
399 
400  // print global mode and scale factors for the 1st chain into a text file
401  if (WritePreRunData) {
402  std::ofstream outPreRun;
403  outPreRun.open(("PreRun" + JobTag + ".txt").c_str(), std::ios::out);
404  outPreRun << MCEngine.writePreRunData();
405  outPreRun.close();
406  }
407 
408  /* Number of events */
409  std::stringstream ss;
410  ss << "";
411  BCLog::SetPrefix(false);
412  BCLog::OutSummary(ss.str().c_str());
413  BCLog::SetPrefix(true);
414  ss << "Number of used events: " << MCEngine.getNumOfUsedEvents();
415  BCLog::OutSummary(ss.str().c_str());
416  ss.str("");
417  ss << "Number of discarded events: " << MCEngine.getNumOfDiscardedEvents() << "\n";
418  BCLog::OutSummary(ss.str().c_str());
419  // close log file
420  BCLog::CloseLog();
421 
422  // print results of the analysis into a text file
423  BCLog::SetPrefix(false);
424  BCLog::OpenLog(("MonteCarlo_results" + JobTag + ".txt").c_str(), BCLog::summary, BCLog::nothing);
425  MCEngine.PrintSummary();
426  BCLog::CloseLog();
427 
428  // Calculate and print normalization
429  if (CalculateNormalization.compare("false") != 0) {
430  std::ofstream outStatLog;
431  outStatLog.open((ObsDirName + "/Statistics" + JobTag + ".txt").c_str(), std::ios::app);
432 
433  if (CalculateNormalization.compare("LME") == 0) {
434  // Make sure mode is found by MINUIT for normalization computation.
435  if (!FindModeWithMinuit) MCEngine.FindMode(MCEngine.GetBestFitParameters());
436  outStatLog << "\nNormalization for " << ModelName.c_str() << ": " << MCEngine.computeNormalizationLME() << "\n" << std::endl;
437  }
438  else if (CalculateNormalization.compare("MC") == 0) {
439  std::vector<double> normalizationMC = MCEngine.computeNormalizationMC(NIterationNormalizationMC);
440  outStatLog << "\nNormalization for " << ModelName.c_str() << ": " << normalizationMC[0] << " +- " << normalizationMC[1] << "\n" << std::endl;
441  }
442  else
443  throw std::runtime_error(("\n ERROR: Normalization method" + CalculateNormalization + " not implemented.\n").c_str());
444 
445  outStatLog.close();
446  }
447 
448  // Print the triangle plot (the size of the file is quite large for large number of parameters).
449  // this is done last for now since BAT sometimes crashes while doing this. Did not look for the fix.
450  if (PrintTrianglePlot) {
451  BCLog::OpenLog(("log_CorrelationPlots" + JobTag + ".txt").c_str(), BCLog::debug, BCLog::debug);
452  MCEngine.PrintCorrelationPlot("CorrelationPlots" + JobTag + ".pdf");
453  BCLog::CloseLog();
454  }
455 
456 
457 #ifdef _MPI
458  double ** sendbuff = new double *[MCEngine.procnum];
459  sendbuff[0] = new double[MCEngine.procnum * buffsize];
460  for (int il = 1; il < MCEngine.procnum; il++) {
461  sendbuff[il] = sendbuff[il - 1] + buffsize;
462  sendbuff[il][0] = -1.; //Exit command
463  }
464  MPI_Scatter(sendbuff[0], buffsize, MPI_DOUBLE, recvbuff, buffsize, MPI_DOUBLE, 0, MPI_COMM_WORLD);
465  delete sendbuff[0];
466  delete [] sendbuff;
467 #endif
468  }
469 #ifdef _MPI
470  delete [] recvbuff;
471 #endif
472  } catch (std::string message) {
473  std::cerr << message << std::endl;
474  exit(EXIT_FAILURE);
475  }
476 }
477 
478 void MonteCarlo::ParseMCMCConfig(std::string file)
479 {
480  std::ifstream ifile(file.c_str());
481  if (!ifile.is_open())
482  throw std::runtime_error("\nERROR: " + MCMCConf + " does not exist. Make sure to specify a valid Monte Carlo configuration file.\n");
483  std::string line;
484  bool IsEOF = false;
485  do {
486  IsEOF = getline(ifile, line).eof();
487  if (*line.rbegin() == '\r') line.erase(line.length() - 1); // for CR+LF
488  if (line.empty() || line.at(0) == '#')
489  continue;
490  boost::char_separator<char> sep(" \t");
491  boost::tokenizer<boost::char_separator<char> > tok(line, sep);
492  boost::tokenizer<boost::char_separator<char> >::iterator beg = tok.begin();
493  if (beg->compare("NChains") == 0) {
494  ++beg;
495  if (isdigit(beg->at(0)) && atoi((*beg).c_str()) > 0) MCEngine.setNChains(atoi((*beg).c_str()));
496  else if (isdigit(beg->at(0)) && atoi((*beg).c_str()) == 0)
497  if (MCEngine.getMPIWorldSize() > 1) MCEngine.setNChains(MCEngine.getMPIWorldSize());
498  else throw std::runtime_error("\nERROR: NChains in the MonteCarlo configuration file: " + MCMCConf + " can only be an integer > 0.\n");
499  else
500  throw std::runtime_error("\nERROR: NChains in the MonteCarlo configuration file: " + MCMCConf + " can only be an integer >= 0.\n");
501  } else if (beg->compare("PrerunMaxIter") == 0) {
502  ++beg;
503  if (isdigit(beg->at(0)) && atoi((*beg).c_str()) > 0) MCEngine.SetNIterationsPreRunMax(atoi((*beg).c_str()));
504  else
505  throw std::runtime_error("\nERROR: PrerunMaxIter in the MonteCarlo configuration file: " + MCMCConf + " can only be an integer > 0.\n");
506  } else if (beg->compare("NIterationsUpdateMax") == 0) {
507  ++beg;
508  if (isdigit(beg->at(0)) && atoi((*beg).c_str()) > 0) MCEngine.SetNIterationsPreRunCheck(atoi((*beg).c_str()));
509  else
510  throw std::runtime_error("\nERROR: NIterationsUpdateMax in the MonteCarlo configuration file: " + MCMCConf + " can only be an integer > 0.\n");
511  } else if (beg->compare("Seed") == 0) {
512  ++beg;
513  if (!isdigit(beg->at(0))) throw std::runtime_error("\nERROR: Seed in the MonteCarlo configuration file: " + MCMCConf + " can only be a number.\n");
514  int seed = atoi((*beg).c_str());
515  if (seed != 0)
516  MCEngine.SetRandomSeed(seed);
517  } else if (beg->compare("Iterations") == 0) {
518  ++beg;
519  if (isdigit(beg->at(0))) MCEngine.SetNIterationsRun(atoi((*beg).c_str()));
520  else
521  throw std::runtime_error("\nERROR: Iterations in the MonteCarlo configuration file: " + MCMCConf + " can only be an integer.\n");
522  } else if (beg->compare("MinimumEfficiency") == 0) {
523  ++beg;
524  double efficiency = atof((*beg).c_str());
525  if (efficiency > 0. && efficiency <= 1.) MCEngine.SetMinimumEfficiency(efficiency);
526  else
527  throw std::runtime_error("\nERROR: MinimumEfficiency in the MonteCarlo configuration file: " + MCMCConf + " can only be a real number greater than 0.0 and less than or equal to 1.0.\n");
528  } else if (beg->compare("MaximumEfficiency") == 0) {
529  ++beg;
530  double efficiency = atof((*beg).c_str());
531  if (efficiency > 0. && efficiency <= 1.) MCEngine.SetMaximumEfficiency(efficiency);
532  else
533  throw std::runtime_error("\nERROR: MaximumEfficiency in the MonteCarlo configuration file: " + MCMCConf + " can only be a real number greater than 0.0 and less than or equal to 1.0.\n");
534  } else if (beg->compare("RValueForConvergence") == 0) {
535  ++beg;
536  double rvalue = atof((*beg).c_str());
537  if (rvalue > 1. ) MCEngine.SetRValueParametersCriterion(rvalue);
538  else
539  throw std::runtime_error("\nERROR: RValueForConvergence in the MonteCarlo configuration file: " + MCMCConf + " can only be a real number greater than 1.0.\n");
540  } else if (beg->compare("WriteChain") == 0) {
541  ++beg;
542  if (beg->compare("true") == 0 || beg->compare("false") == 0) writechains = (beg->compare("true") == 0);
543  else
544  throw std::runtime_error("\nERROR: WriteChain in the MonteCarlo configuration file: " + MCMCConf + " can only be 'true' or 'false'.\n");
545  } else if (beg->compare("FindModeWithMinuit") == 0) {
546  ++beg;
547  if (beg->compare("true") == 0 || beg->compare("false") == 0) FindModeWithMinuit = (beg->compare("true") == 0);
548  else
549  throw std::runtime_error("\nERROR: FindModeWithMinuit in the MonteCarlo configuration file: " + MCMCConf + " can only be 'true' or 'false'.\n");
550  } else if (beg->compare("RunMinuitOnly") == 0) {
551  ++beg;
552  if (beg->compare("true") == 0 || beg->compare("false") == 0) RunMinuitOnly = (beg->compare("true") == 0);
553  else
554  throw std::runtime_error("\nERROR: RunMinuitOnly in the MonteCarlo configuration file: " + MCMCConf + " can only be 'true' or 'false'.\n");
555  } else if (beg->compare("CalculateNormalization") == 0) {
556  ++beg;
557  if (beg->compare("true") == 0 || beg->compare("false") == 0) CalculateNormalization = *beg;
558  else
559  throw std::runtime_error("\nERROR: CalculateNormalization in the MonteCarlo configuration file: " + MCMCConf + " can only be 'true' or 'false'.\n");
560  } else if (beg->compare("NIterationNormalizationMC") == 0) {
561  ++beg;
562  if (isdigit(beg->at(0))) NIterationNormalizationMC = atoi((*beg).c_str());
563  else
564  throw std::runtime_error("\nERROR: NIterationNormalizationMC in the MonteCarlo configuration file: " + MCMCConf + " can only be an integer.\n");
565  } else if (beg->compare("PrintAllMarginalized") == 0) {
566  ++beg;
567  if (beg->compare("true") == 0 || beg->compare("false") == 0) PrintAllMarginalized = (beg->compare("true") == 0);
568  else
569  throw std::runtime_error("\nERROR: PrintAllMarginalized in the MonteCarlo configuration file: " + MCMCConf + " can only be 'true' or 'false'.\n");
570  } else if (beg->compare("PrintCorrelationMatrix") == 0) {
571  ++beg;
572  if (beg->compare("true") == 0 || beg->compare("false") == 0) PrintCorrelationMatrix = (beg->compare("true") == 0);
573  else
574  throw std::runtime_error("\nERROR: PrintCorrelationMatrix in the MonteCarlo configuration file: " + MCMCConf + " can only be 'true' or 'false'.\n");
575  } else if (beg->compare("PrintKnowledgeUpdatePlots") == 0) {
576  ++beg;
577  if (beg->compare("true") == 0 || beg->compare("false") == 0) PrintKnowledgeUpdatePlots = (beg->compare("true") == 0);
578  else
579  throw std::runtime_error("\nERROR: PrintKnowledgeUpdatePlots in the MonteCarlo configuration file: " + MCMCConf + " can only be 'true' or 'false'.\n");
580  } else if (beg->compare("PrintParameterPlot") == 0) {
581  ++beg;
582  if (beg->compare("true") == 0 || beg->compare("false") == 0) PrintParameterPlot = (beg->compare("true") == 0);
583  else
584  throw std::runtime_error("\nERROR: PrintParameterPlot in the MonteCarlo configuration file: " + MCMCConf + " can only be 'true' or 'false'.\n");
585  } else if (beg->compare("PrintTrianglePlot") == 0) {
586  ++beg;
587  if (beg->compare("true") == 0 || beg->compare("false") == 0) PrintTrianglePlot = (beg->compare("true") == 0);
588  else
589  throw std::runtime_error("\nERROR: PrintTrianglePlot in the MonteCarlo configuration file: " + MCMCConf + " can only be 'true' or 'false'.\n");
590  } else if (beg->compare("WritePreRunData") == 0) {
591  ++beg;
592  if (beg->compare("true") == 0 || beg->compare("false") == 0) WritePreRunData = (beg->compare("true") == 0);
593  else
594  throw std::runtime_error("\nERROR: WritePreRunData in the MonteCarlo configuration file: " + MCMCConf + " can only be 'true' or 'false'.\n");
595  } else if (beg->compare("ReadPreRunData") == 0) {
596  ++beg;
597  ReadPreRunData(*beg);
598  } else if (beg->compare("MultivariateProposal") == 0) {
599  ++beg;
600  if (beg->compare("true") == 0 || beg->compare("false") == 0) MCEngine.SetProposeMultivariate((beg->compare("true") == 0));
601  else
602  throw std::runtime_error("\nERROR: MultivariateProposal in the MonteCarlo configuration file: " + MCMCConf + " can only be 'true' or 'false'.\n");
603  } else if (beg->compare("Histogram1DSmooth") == 0) {
604  ++beg;
605  if (beg->compare("true") == 0) {
606  MCEngine.setSmooth(1);
607  } else if (beg->compare("false") == 0) {
608  MCEngine.setSmooth(0); /* Default */
609  } else if (isdigit(beg->at(0))) {
610  if (atoi((*beg).c_str()) >= 0 && atoi((*beg).c_str()) <= 5) MCEngine.setSmooth(atoi((*beg).c_str()));
611  } else
612  throw std::runtime_error("\nERROR: Histogram1DSmooth in the MonteCarlo configuration file: " + MCMCConf + " can only be 'true', 'false' or an integer from 0 to 5.\n");
613  } else if (beg->compare("Histogram2DType") == 0) {
614  ++beg;
615  if (!isdigit(beg->at(0)))
616  throw std::runtime_error("\nERROR: Histogram2DType in the MonteCarlo configuration file : " + MCMCConf + "can only be an integer amongst 1001 -> Lego, 101 -> Filled, 1 -> Contour.\n");
617  int type = atoi((*beg).c_str());
618  if (type == 1 || type == 101 || type == 1001) {
619  MCEngine.setHistogram2DType(type);
620  } else
621  throw std::runtime_error("\nERROR: Histogram2DType in the MonteCarlo configuration file : " + MCMCConf + "can only be an integer amongst 1001 -> Lego, 101 -> Filled, 1 -> Contour.\n");
622  } else if (beg->compare("MCMCInitialPosition") == 0) {
623  ++beg;
624  if (beg->compare("Center") == 0) {
625  MCEngine.SetInitialPositionScheme(BCEngineMCMC::kInitCenter); /* Default */
626  } else if (beg->compare("RandomUniform") == 0) {
627  MCEngine.SetInitialPositionScheme(BCEngineMCMC::kInitRandomUniform);
628  } else if (beg->compare("RandomPrior") == 0) {
629  MCEngine.SetInitialPositionScheme(BCEngineMCMC::kInitRandomPrior);
630  } else
631  throw std::runtime_error("\nERROR: MCMCInitialPosition in MonteCarlo config file: " + MCMCConf + " can only be 'Center', 'RandomUniform' or 'RandomPrior'.\n");
632  } else if (beg->compare("PrintLogo") == 0) {
633  ++beg;
634  if (beg->compare("true") == 0 || beg->compare("false") == 0) MCEngine.setPrintLogo((beg->compare("true") == 0));
635  else
636  throw std::runtime_error("\nERROR: PrintLogo in the MonteCarlo configuration file: " + MCMCConf + " can only be 'true' or 'false'.\n");
637  } else if (beg->compare("NoHistogramLegend") == 0) {
638  ++beg;
639  if (beg->compare("true") == 0 || beg->compare("false") == 0) MCEngine.setNoLegend((beg->compare("true") == 0));
640  else
641  throw std::runtime_error("\nERROR: PrintLogo in the MonteCarlo configuration file: " + MCMCConf + " can only be 'true' or 'false'.\n");
642  } else if (beg->compare("PrintLoglikelihoodPlots") == 0) {
643  ++beg;
644  if (beg->compare("true") == 0 || beg->compare("false") == 0) MCEngine.setPrintLoglikelihoodPlots((beg->compare("true") == 0));
645  else
646  throw std::runtime_error("\nERROR: PrintLoglikelihoodPlots in the MonteCarlo configuration file: " + MCMCConf + " can only be 'true' or 'false'.\n");
647  if (beg->compare("true") == 0) gSystem->MakeDirectory((ObsDirName + "/LogLikelihoodPlots").c_str());
648  } else if (beg->compare("WriteLogLikelihoodChain") == 0) {
649  ++beg;
650  if (beg->compare("true") == 0 || beg->compare("false") == 0) MCEngine.setWriteLogLikelihoodChain((beg->compare("true") == 0));
651  else
652  throw std::runtime_error("\nERROR: WriteLogLikelihoodChain in the MonteCarlo configuration file: " + MCMCConf + " can only be 'true' or 'false'.\n");
653  } else if (beg->compare("WriteParametersChains") == 0) {
654  ++beg;
655  if (beg->compare("true") == 0 || beg->compare("false") == 0) MCEngine.setWriteParametersChain((beg->compare("true") == 0));
656  else
657  throw std::runtime_error("\nERROR: WriteParametersChain in the MonteCarlo configuration file: " + MCMCConf + " can only be 'true' or 'false'.\n");
658  } else if (beg->compare("Histogram2DAlpha") == 0) {
659  ++beg;
660  double alpha = atof((*beg).c_str());
661  if (alpha > 0. && alpha <= 1.) MCEngine.setAlpha2D(alpha);
662  else
663  throw std::runtime_error("\nERROR: Histogram2DAlpha in the MonteCarlo configuration file: " + MCMCConf + " can only be a real number greater than 0.0 and less than or equal to 1.0.\n");
664  } else if (beg->compare("NBinsHistogram1D") == 0) {
665  ++beg;
666  unsigned int nBins1D = atoi((*beg).c_str());
667  if (nBins1D >= 0) MCEngine.setNBins1D(nBins1D);
668  else throw std::runtime_error("\nERROR: NBinsHistogram1D in the MonteCarlo configuration file: " + MCMCConf + " can only be a integer greater than 0 or 0 to set to default value (100).\n");
669  } else if (beg->compare("NBinsHistogram2D") == 0) {
670  ++beg;
671  unsigned int nBins2D = atoi((*beg).c_str());
672  if (nBins2D >= 0) MCEngine.setNBins2D(nBins2D);
673  else throw std::runtime_error("\nERROR: NBinsHistogram2D in the MonteCarlo configuration file: " + MCMCConf + " can only be a integer greater than 0 or 0 to set to default value (100).\n");
674  } else if (beg->compare("InitialPositionAttemptLimit") == 0) {
675  ++beg;
676  unsigned int max_tries = atoi((*beg).c_str());
677  if (max_tries > 0) MCEngine.SetInitialPositionAttemptLimit(max_tries);
678  else if (max_tries == 0) MCEngine.SetInitialPositionAttemptLimit(10000);
679  else throw std::runtime_error("\nERROR: InitialPositionAttemptLimit in the MonteCarlo configuration file: " + MCMCConf + " can only be a integer greater than 0 or 0 to set to default value (100).\n");
680  } else if (beg->compare("SignificantDigits") == 0) {
681  ++beg;
682  unsigned int significants = atoi((*beg).c_str());
683  if (significants >= 0) MCEngine.setSignificants(significants);
684  else throw std::runtime_error("\nERROR: SignificantDigits in the MonteCarlo configuration file: " + MCMCConf + " can only be a integer greater than 0 or 0 to set to default values.\n");
685  } else if (beg->compare("HistogramBufferSize") == 0) {
686  ++beg;
687  unsigned int histogramBufferSize = atoi((*beg).c_str());
688  if (histogramBufferSize >= 0) MCEngine.setHistogramBufferSize(histogramBufferSize);
689  else throw std::runtime_error("\nERROR: HistogramBufferSize in the MonteCarlo configuration file: " + MCMCConf + " can only be a integer greater than 0 or 0 to set to default values.\n");
690  } else
691  throw std::runtime_error("\nERROR: Wrong keyword in MonteCarlo config file: " + MCMCConf + "\n Make sure to specify a valid Monte Carlo configuration file.\n");
692  } while (!IsEOF);
693  if (FindModeWithMinuit || RunMinuitOnly)
694  MCEngine.SetOptimizationMethod(BCIntegrate::kOptMinuit);
695  if (MCEngine.GetMaximumEfficiency() <= MCEngine.GetMinimumEfficiency())
696  throw std::runtime_error("\nERROR: MaximumEfficiency (default 0.5) must be greater than MaximumEfficiency (default 0.15) in the MonteCarlo configuration file: " + MCMCConf + ".\n");
697  if (CalculateNormalization.compare("MC") == 0 && NIterationNormalizationMC <= 0)
698  throw std::runtime_error(("\nMonteCarlo ERROR: CalculateNormalization cannot be set to MC without setting NIterationNormalizationMC > 0 in " + MCMCConf + " .\n").c_str());
699 }
700 
701 void MonteCarlo::ReadPreRunData(std::string file)
702 {
703  std::ifstream ifile(file.c_str());
704  if (!ifile.is_open())
705  throw std::runtime_error("\nMonteCarlo::ReadPreRunData ERROR: " + file + " does not exist.\n");
706  std::string line;
707  bool IsEOF = false;
708  std::vector<double> mode;
709  std::vector<double> scale;
710  do {
711  IsEOF = getline(ifile, line).eof();
712  if (line.empty())
713  continue;
714  boost::char_separator<char> sep(" ");
715  boost::tokenizer<boost::char_separator<char> > tok(line, sep);
716  boost::tokenizer<boost::char_separator<char> >::iterator beg = tok.begin();
717  ++beg;
718  mode.push_back(atof((*beg).c_str()));
719  ++beg;
720  scale.push_back(atof((*beg).c_str()));
721  } while (!IsEOF);
722  if (mode.size() != MCEngine.GetNParameters())
723  throw std::runtime_error("\nMonteCarlo::ReadPreRunData ERROR: wrong data size.\n");
724  std::vector<double> mode_all;
725  std::vector<double> scale_all;
726  for (unsigned int i = 0; i < MCEngine.GetNChains(); i++){
727  mode_all.insert(mode_all.end(), mode.begin(), mode.end());
728  scale_all.insert(scale_all.end(), scale.begin(), scale.end());
729  }
730  MCEngine.SetInitialPositions(mode_all);
731  MCEngine.SetInitialScaleFactors(scale_all);
732 }
733 
734 void MonteCarlo::addCustomObservableType(const std::string name, boost::function<Observable*() > funct) {
735  myInputParser.addCustomObservableType(name, funct);
736 }
MonteCarlo::TestRun
void TestRun(int rank)
This member is used for test runa to generate a single event.
Definition: MonteCarlo.cpp:59
ModelFactory
A class for.
Definition: ModelFactory.h:25
gslpp::sec
complex sec(const complex &z)
Definition: gslpp_complex.cpp:447
MonteCarlo::checkrun
bool checkrun
A check to make sure TestRun()and Run() are not called consecutively.
Definition: MonteCarlo.h:178
MonteCarlo::ObsDirName
std::string ObsDirName
String for the output directory name.
Definition: MonteCarlo.h:167
MonteCarlo::addCustomObservableType
void addCustomObservableType(const std::string name, boost::function< Observable *() > funct)
Definition: MonteCarlo.cpp:734
MonteCarloEngine::setAlpha2D
void setAlpha2D(double alpha)
A set method to toggle the printing of legends in 1D and 2D histograms.
Definition: MonteCarloEngine.h:359
MonteCarlo::CalculateNormalization
std::string CalculateNormalization
< Flag for calculating the evidence.
Definition: MonteCarlo.h:170
MonteCarloEngine::PrintHistogram
void PrintHistogram(std::string &OutFile, Observable &it, const std::string OutputDir)
Overloaded from PrintHistogram(BCModelOutput&, const std::string) to print histogram for observables.
Definition: MonteCarloEngine.cpp:786
MonteCarlo::RunMinuitOnly
bool RunMinuitOnly
Flag for running Minuit only.
Definition: MonteCarlo.h:169
MonteCarloEngine::CreateHistogramMaps
void CreateHistogramMaps()
Creation of the histogram maps for Observables, Observable2D and Correlated Gaussian Observable.
Definition: MonteCarloEngine.cpp:121
MonteCarlo::myInputParser
InputParser myInputParser
An object of the InputParser class.
Definition: MonteCarlo.h:156
Model::getmissingModelParameters
std::vector< std::string > getmissingModelParameters()
Definition: Model.h:237
MonteCarloEngine::setHistogramBufferSize
void setHistogramBufferSize(unsigned int histogramBufferSize_i)
A set method to set the size of the buffer used by the histograms.
Definition: MonteCarloEngine.h:413
MonteCarlo::JobTag
std::string JobTag
String for the optional JobTag argument to be passes to the executable.
Definition: MonteCarlo.h:166
MonteCarloEngine::PrintCorrelationMatrixToLaTeX
void PrintCorrelationMatrixToLaTeX(const std::string filename)
This member generates the correlation matrix using BCH2D from the BAT libraries.
Definition: MonteCarloEngine.cpp:973
MonteCarlo::PrintAllMarginalized
bool PrintAllMarginalized
Flag for printing all Marginalized distributions to be passed on to the BAT routines.
Definition: MonteCarlo.h:172
MonteCarloEngine::getWriteLogLikelihoodChain
bool getWriteLogLikelihoodChain()
A get method to get the value of the bool the writing of Loglikelihood in the ROOT tree.
Definition: MonteCarloEngine.h:333
gslpp::log
complex log(const complex &z)
Definition: gslpp_complex.cpp:342
MonteCarlo::NIterationNormalizationMC
int NIterationNormalizationMC
< Number of iterations for MC integral done to compute normalization of a model
Definition: MonteCarlo.h:171
ModelFactory.h
MonteCarloEngine::AddChains
void AddChains()
A method to add the observable values and weights to the chain information.
Definition: MonteCarloEngine.cpp:884
MonteCarloEngine::setPrintLogo
void setPrintLogo(bool print)
A set method to toggle the printing of logo on the histogram pdf.
Definition: MonteCarloEngine.h:298
MonteCarloEngine::setNBins1D
void setNBins1D(unsigned int nbins)
A set method to set the number of bins for a 1D histograms.
Definition: MonteCarloEngine.h:386
MonteCarloEngine::getMPIWorldSize
int getMPIWorldSize()
A get method to get the number of MPI world size.
Definition: MonteCarloEngine.h:422
ThObsFactory.h
MonteCarlo::Obs2D
std::vector< Observable2D > Obs2D
Vector for the Observables2D defined in SomeModel.conf.
Definition: MonteCarlo.h:160
InputParser::ReadParameters
std::string ReadParameters(const std::string filename_i, const int rank, std::vector< ModelParameter > &ModelPars, boost::ptr_vector< Observable > &Observables, std::vector< Observable2D > &Observables2D, std::vector< CorrelatedGaussianObservables > &CGO, std::vector< CorrelatedGaussianParameters > &CGP)
The member that parses the Observable2D directives from SomeModel.conf file.
Definition: InputParser.cpp:35
ThObsFactory
A class for.
Definition: ThObsFactory.h:26
StandardModel::Init
virtual bool Init(const std::map< std::string, double > &DPars)
A method to initialize the model parameters.
Definition: StandardModel.cpp:159
MonteCarloEngine::setNBins2D
void setNBins2D(unsigned int nbins)
A set method to set the number of bins for a 2D histograms.
Definition: MonteCarloEngine.h:395
MonteCarloEngine::getchainedObsSize
int getchainedObsSize() const
A get method to access the number of Observable chains requested.
Definition: MonteCarloEngine.h:195
MonteCarlo::CGP
std::vector< CorrelatedGaussianParameters > CGP
Vector for the Correlated Gaussian Parameters defined in SomeModel.conf.
Definition: MonteCarlo.h:162
MonteCarlo::PrintTrianglePlot
bool PrintTrianglePlot
Flag for printing the triangle plot.
Definition: MonteCarlo.h:176
MonteCarloEngine::setDParsFromParameters
void setDParsFromParameters(const std::vector< double > &parameters, std::map< std::string, double > &DPars_i)
A method to rotate the diagonalized parameters to the original basis for correlated parameters.
Definition: MonteCarloEngine.cpp:261
MonteCarlo::OutFile
std::string OutFile
String for the name of the output root file without the .root extension.
Definition: MonteCarlo.h:165
MonteCarloEngine::getNumOfUsedEvents
int getNumOfUsedEvents() const
A get method to access the number of events used in the MCMC run.
Definition: MonteCarloEngine.h:233
MonteCarlo::writechains
bool writechains
Flag for writing the chains of paramters and observables during the MCMC run.
Definition: MonteCarlo.h:179
MonteCarlo::ParseMCMCConfig
void ParseMCMCConfig(std::string file)
Definition: MonteCarlo.cpp:478
MonteCarlo::PrintCorrelationMatrix
bool PrintCorrelationMatrix
Flag for printing the correlation matrix.
Definition: MonteCarlo.h:173
MonteCarloEngine::computeNormalizationLME
double computeNormalizationLME()
A method to calculate the normalization based on the Laplace-Metropolis Estimator.
Definition: MonteCarloEngine.cpp:1335
MonteCarlo::Obs
boost::ptr_vector< Observable > Obs
Vector for the observables defined in SomeModel.conf.
Definition: MonteCarlo.h:159
MonteCarloEngine::computeStatistics
std::string computeStatistics()
A get method to compute the mean and rms of the computed observables.
Definition: MonteCarloEngine.cpp:1012
MonteCarlo::Run
void Run(const int rank)
This member is responsible for setting the Monte Carlo run parameters and conducting the Monte Carlo ...
Definition: MonteCarlo.cpp:128
InputParser::getModel
StandardModel * getModel() const
A get method to access the pointer to the object of the StandardModel class.
Definition: InputParser.h:109
MonteCarloEngine::setHistogram2DType
void setHistogram2DType(int type)
A set method to set the band fill type for 2D histograms.
Definition: MonteCarloEngine.h:377
MonteCarlo::ModPars
std::vector< ModelParameter > ModPars
Vector for the model parameters defined in SomeModel.conf.
Definition: MonteCarlo.h:158
MonteCarlo::WritePreRunData
bool WritePreRunData
Flag for printing the overview parameter plots.
Definition: MonteCarlo.h:177
MonteCarloEngine::writePreRunData
std::string writePreRunData()
A method to write in a text file the best fit parameters and the prerun scale factors.
Definition: MonteCarloEngine.cpp:1306
MonteCarloEngine::setNChains
void setNChains(unsigned int i)
A set method to fix the number of chains.
Definition: MonteCarloEngine.cpp:178
Observable
A class for observables.
Definition: Observable.h:28
MonteCarlo::~MonteCarlo
virtual ~MonteCarlo()
The default destructor.
Definition: MonteCarlo.cpp:53
MonteCarlo::ReadPreRunData
void ReadPreRunData(std::string file)
Definition: MonteCarlo.cpp:701
MonteCarlo::CGO
std::vector< CorrelatedGaussianObservables > CGO
Vector for the Correlated Gaussian Observables defined in SomeModel.conf.
Definition: MonteCarlo.h:161
MonteCarlo::MCEngine
MonteCarloEngine MCEngine
An object of the MonteCarloEngine class.
Definition: MonteCarlo.h:157
MonteCarlo::PrintKnowledgeUpdatePlots
bool PrintKnowledgeUpdatePlots
Flag for printing plots to compare prior vs. posterior knowledge of parameters.
Definition: MonteCarlo.h:174
MonteCarloEngine::getHistoLog
std::string getHistoLog() const
A get method to access the stream that stores the log messages coming from histogram printing and che...
Definition: MonteCarloEngine.h:186
MonteCarlo::PrintParameterPlot
bool PrintParameterPlot
Flag for printing the overview parameter plots.
Definition: MonteCarlo.h:175
InputParser::addCustomObservableType
void addCustomObservableType(const std::string name, boost::function< Observable *() > funct)
Definition: InputParser.cpp:250
MonteCarloEngine::setNoLegend
void setNoLegend(bool legend)
A set method to toggle the printing of legends in 1D and 2D histograms.
Definition: MonteCarloEngine.h:307
MonteCarlo::FindModeWithMinuit
bool FindModeWithMinuit
Flag for using Minuit libraries.
Definition: MonteCarlo.h:168
MonteCarloEngine::computeNormalizationMC
std::vector< double > computeNormalizationMC(int NIterationNormalizationMC)
A method to calculate the normalization based on the Monte Carlo Simulation.
Definition: MonteCarloEngine.cpp:1319
StandardModel::Update
virtual bool Update(const std::map< std::string, double > &DPars)
The update method for StandardModel.
Definition: StandardModel.cpp:183
MonteCarloEngine::setSmooth
void setSmooth(int int_N)
A set method to set the number of smoothing passes for ROOT in 1D histograms.
Definition: MonteCarloEngine.h:368
MonteCarlo::MonteCarlo
MonteCarlo(ModelFactory &ModelF, ThObsFactory &ThObsF, const std::string &ModelConf_i, const std::string &MonteCarloConf_i, const std::string &OutFile_i, const std::string &JobTag_i)
Constructor.
Definition: MonteCarlo.cpp:26
MonteCarlo::ModelConf
std::string ModelConf
String for the name of the SomeModel.conf file.
Definition: MonteCarlo.h:163
MonteCarloEngine::setSignificants
void setSignificants(unsigned int significants_i)
A set method to set the number of significant digits in the output.
Definition: MonteCarloEngine.h:404
MonteCarloEngine::getWriteParametersChain
bool getWriteParametersChain()
A get method to get the value of the bool the writing of parameters in the ROOT tree.
Definition: MonteCarloEngine.h:350
MonteCarloEngine::setWriteParametersChain
void setWriteParametersChain(bool LL)
A set method to toggle the writing of parameters in the ROOT tree.
Definition: MonteCarloEngine.h:342
MonteCarlo::MCMCConf
std::string MCMCConf
String for the name of the MonteCarlo.conf file.
Definition: MonteCarlo.h:164
MonteCarlo.h
MonteCarloEngine::Initialize
void Initialize(StandardModel *Mod_i)
Initialization of the Monte Carlo Engine.
Definition: MonteCarloEngine.cpp:76
MonteCarloEngine::setWriteLogLikelihoodChain
void setWriteLogLikelihoodChain(bool LL)
A set method to toggle the writing of Loglikelihood in the ROOT tree.
Definition: MonteCarloEngine.h:325
MonteCarlo::ModelName
std::string ModelName
The name of the model.
Definition: MonteCarlo.h:155
QCD::getUnknownParameters
std::vector< std::string > getUnknownParameters()
A method to get the vector of the parameters that have been specified in the configuration file but n...
Definition: QCD.h:467
MonteCarloEngine::getNumOfDiscardedEvents
int getNumOfDiscardedEvents() const
A get method to access the number of events discarded due to failure to update model....
Definition: MonteCarloEngine.h:224
MonteCarloEngine::setPrintLoglikelihoodPlots
void setPrintLoglikelihoodPlots(bool plot)
A set method to toggle the printing of Parameters vs. Loglikelihood.
Definition: MonteCarloEngine.h:316