MonteCarlo Class Reference

A class for Monte Carlo. More...

#include <MonteCarlo.h>

Collaboration diagram for MonteCarlo:
[legend]

Detailed Description

A class for Monte Carlo.

Author
HEPfit Collaboration

This class is responsible for the MCMC runs using the MCMCEngine class together with input parameters from the InputParser class. The Monte Carlo mode is the default run mode if the –noMC flag is not specified. The settings for the Monte Carlo run can be fixed in in the format specified in the MonteCarlo.conf file, the name of this configuration file should be specified as the second argument to the executable.

Definition at line 35 of file MonteCarlo.h.

Public Member Functions

void addCustomObservableType (const std::string name, boost::function< Observable *() > funct)
 
std::map< std::string, BCH1D * > getHistograms1D () const
 
std::map< std::string, BCH2D * > getHistograms2D () const
 
 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. More...
 
void Run (const int rank)
 This member is responsible for setting the Monte Carlo run parameters and conducting the Monte Carlo run including initiating all output generation. More...
 
void TestRun (int rank)
 The default destructor. More...
 

Private Member Functions

void ReadPreRunData (std::string file)
 

Private Attributes

bool CalculateNormalization
 < Flag for calculating the evidence. More...
 
std::vector< CorrelatedGaussianObservablesCGO
 Vector for the Correlated Gaussian Observables defined in SomeModel.conf. More...
 
std::vector< CorrelatedGaussianParametersCGP
 Vector for the Correlated Gaussian Parameters defined in SomeModel.conf. More...
 
bool checkrun
 A check to make sure TestRun()and Run() are not called consecutively. More...
 
bool FindModeWithMinuit
 Flag for using Minuit libraries. More...
 
std::string JobTag
 String for the optional JobTag argument to be passes to the executable. More...
 
MonteCarloEngine MCEngine
 An object of the MonteCarloEngine class. More...
 
std::string MCMCConf
 String for the name of the MonteCarlo.conf file. More...
 
std::string ModelConf
 String for the name of the SomeModel.conf file. More...
 
std::string ModelName
 The name of the model. More...
 
std::vector< ModelParameterModPars
 Vector for the model parameters defined in SomeModel.conf. More...
 
InputParser myInputParser
 An object of the InputParser class. More...
 
double normalization
 A variable to store the evidence of a model. More...
 
boost::ptr_vector< ObservableObs
 Vector for the observables defined in SomeModel.conf. More...
 
std::vector< Observable2DObs2D
 Vector for the Observables2D defined in SomeModel.conf. More...
 
std::string ObsDirName
 String for the output directory name. More...
 
std::string OutFile
 String for the name of the output root file without the .root extension. More...
 
bool PrintAllMarginalized
 Flag for printing all Marginalized distributions to be passed on to the BAT routines. More...
 
bool PrintCorrelationMatrix
 Flag for printing the correlation matrix. More...
 
bool PrintKnowledgeUpdatePlots
 Flag for printing plots to compare prior vs. posterior knowledge of parameters. More...
 
bool PrintParameterPlot
 Flag for printing the overview parameter plots. More...
 
bool WritePreRunData
 Flag for printing the overview parameter plots. More...
 

Constructor & Destructor Documentation

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.

The default constructor sets the names of the configuration files and the names of the output root file and specific output directory if a job ID is specified as an argument to the executable. The constructor also sets the default values of some flags related to the Monte Carlo run. The boolean flags (FindModeWithMinuit, PrintAllMarginalized, PrintCorrelationMatrix, PrintKnowledgeUpdatePlots, PrintParameterPlot) attain their run time values after the MonteCarlo.conf file is parsed by MonteCarlo::Run().

Parameters
[in]ModelF
[in]ThObsF
[in]ModelConf_ithe name of the input configuration file for the model name, the model parameters and observables to be calculated (Observables, Observables2D, Model Parameters vs. Observables and Correlated Gaussian Observables)
[in]MonteCarloConf_ithe name of the Monte Carlo configuration file that specifies the parameters of the Monte Carlo run like no. of chains, no. of pre run iterations etc.
[in]OutFile_ithe name of the root output file to be given without the .root extention
[in]JobTag_ioptional job tag that might be specified

Definition at line 26 of file MonteCarlo.cpp.

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  CalculateNormalization = false;
41  PrintAllMarginalized = false;
42  PrintCorrelationMatrix = false;
44  PrintParameterPlot = false;
45  WritePreRunData = false;
46  checkrun = false;
47  normalization = 0.;
48 }
bool checkrun
A check to make sure TestRun()and Run() are not called consecutively.
Definition: MonteCarlo.h:175
double normalization
A variable to store the evidence of a model.
Definition: MonteCarlo.h:176
std::string ObsDirName
String for the output directory name.
Definition: MonteCarlo.h:167
bool CalculateNormalization
< Flag for calculating the evidence.
Definition: MonteCarlo.h:169
bool PrintParameterPlot
Flag for printing the overview parameter plots.
Definition: MonteCarlo.h:173
MonteCarloEngine MCEngine
An object of the MonteCarloEngine class.
Definition: MonteCarlo.h:157
std::string OutFile
String for the name of the output root file without the .root extension.
Definition: MonteCarlo.h:165
bool PrintAllMarginalized
Flag for printing all Marginalized distributions to be passed on to the BAT routines.
Definition: MonteCarlo.h:170
bool WritePreRunData
Flag for printing the overview parameter plots.
Definition: MonteCarlo.h:174
std::vector< CorrelatedGaussianParameters > CGP
Vector for the Correlated Gaussian Parameters defined in SomeModel.conf.
Definition: MonteCarlo.h:162
std::vector< ModelParameter > ModPars
Vector for the model parameters defined in SomeModel.conf.
Definition: MonteCarlo.h:158
bool PrintKnowledgeUpdatePlots
Flag for printing plots to compare prior vs. posterior knowledge of parameters.
Definition: MonteCarlo.h:172
bool PrintCorrelationMatrix
Flag for printing the correlation matrix.
Definition: MonteCarlo.h:171
std::vector< Observable2D > Obs2D
Vector for the Observables2D defined in SomeModel.conf.
Definition: MonteCarlo.h:160
boost::ptr_vector< Observable > Obs
Vector for the observables defined in SomeModel.conf.
Definition: MonteCarlo.h:159
std::vector< CorrelatedGaussianObservables > CGO
Vector for the Correlated Gaussian Observables defined in SomeModel.conf.
Definition: MonteCarlo.h:161
InputParser myInputParser
An object of the InputParser class.
Definition: MonteCarlo.h:156
bool FindModeWithMinuit
Flag for using Minuit libraries.
Definition: MonteCarlo.h:168
std::string MCMCConf
String for the name of the MonteCarlo.conf file.
Definition: MonteCarlo.h:164
std::string JobTag
String for the optional JobTag argument to be passes to the executable.
Definition: MonteCarlo.h:166
std::string ModelConf
String for the name of the SomeModel.conf file.
Definition: MonteCarlo.h:163

Member Function Documentation

void MonteCarlo::addCustomObservableType ( const std::string  name,
boost::function< Observable *() >  funct 
)

Definition at line 489 of file MonteCarlo.cpp.

489  {
491 }
void addCustomObservableType(const std::string name, boost::function< Observable *() > funct)
InputParser myInputParser
An object of the InputParser class.
Definition: MonteCarlo.h:156
std::map<std::string, BCH1D * > MonteCarlo::getHistograms1D ( ) const
inline

Definition at line 142 of file MonteCarlo.h.

143  {
144  return MCEngine.getHistograms1D();
145  }
MonteCarloEngine MCEngine
An object of the MonteCarloEngine class.
Definition: MonteCarlo.h:157
std::map< std::string, BCH1D * > getHistograms1D() const
std::map<std::string, BCH2D * > MonteCarlo::getHistograms2D ( ) const
inline

Definition at line 147 of file MonteCarlo.h.

148  {
149  return MCEngine.getHistograms2D();
150  }
MonteCarloEngine MCEngine
An object of the MonteCarloEngine class.
Definition: MonteCarlo.h:157
std::map< std::string, BCH2D * > getHistograms2D() const
void MonteCarlo::ReadPreRunData ( std::string  file)
private

Definition at line 456 of file MonteCarlo.cpp.

457 {
458  std::ifstream ifile(file.c_str());
459  if (!ifile.is_open())
460  throw std::runtime_error("\nMonteCarlo::ReadPreRunData ERROR: " + file + " does not exist.\n");
461  std::string line;
462  bool IsEOF = false;
463  std::vector<double> mode;
464  std::vector<double> scale;
465  do {
466  IsEOF = getline(ifile, line).eof();
467  if (line.empty())
468  continue;
469  boost::char_separator<char> sep(" ");
470  boost::tokenizer<boost::char_separator<char> > tok(line, sep);
471  boost::tokenizer<boost::char_separator<char> >::iterator beg = tok.begin();
472  ++beg;
473  mode.push_back(atof((*beg).c_str()));
474  ++beg;
475  scale.push_back(atof((*beg).c_str()));
476  } while (!IsEOF);
477  if (mode.size() != MCEngine.GetNParameters())
478  throw std::runtime_error("\nMonteCarlo::ReadPreRunData ERROR: wrong data size.\n");
479  std::vector<double> mode_all;
480  std::vector<double> scale_all;
481  for (unsigned int i = 0; i < MCEngine.MCMCGetNChains(); i++){
482  mode_all.insert(mode_all.end(), mode.begin(), mode.end());
483  scale_all.insert(scale_all.end(), scale.begin(), scale.end());
484  }
485  MCEngine.MCMCSetInitialPositions(mode_all);
486  MCEngine.MCMCSetTrialFunctionScaleFactor(scale_all);
487 }
MonteCarloEngine MCEngine
An object of the MonteCarloEngine class.
Definition: MonteCarlo.h:157
void MonteCarlo::Run ( const int  rank)

This member is responsible for setting the Monte Carlo run parameters and conducting the Monte Carlo run including initiating all output generation.

The algorithm implemented by this member is as follows:

  • Initiate InputParser::ReadParameters() which read the SomeModel.conf file for setting the model parameters and the observables to be generated. This call also passes on the name of the model to the private member ModelName.
  • Map the model parameter mean values with the map DP and calculate buffsize which can be used for implementing MPI runs. The variable buffsize is incremented only for those model parameters that are varied in the Monte Carlo run.
  • Initiate the model using the map DP checking for consistency between required model parameters and the ones supplied in the SomeModel.conf file
  • Set the model name and initialize the MCMC engine with the model at hand.
  • For MPI runs, send the evaluates of the loglikelihood for the prerun to the slave processes (rank != 0)
  • The master in an MPI run or the process in a serial run then parses the MonterCarlo.conf file to read the parameters for the Monte Carlo run.
  • The root output file is handed over to the object out of type BCModelOutput and some BAT options are set for the log and output histograms.
  • The MCMC is run and all the out put is generated according to the options set in the MonteCarlo.conf file.
  • For a MPI run, the final call to the MPI class prepares the processes for finalizing.

The details for BCSummaryTool, BCLog, BCAux and BCModelOutput can be found in the BAT website. These are used mainly to generate logs and output.

The details of the object MCEngine of type MonteCarloEngine which overloads the BCEngineMCMC class can be found in our documentation of the former class.

Parameters
[in]rank= MPI::COMM_WORLD.Get_rank(), specifies the rank of the process. This carries a non zero value only when the executable is compiled with the parallelized version of BAT and run as parallel processes with MPI.

Definition at line 108 of file MonteCarlo.cpp.

108  {
109  if (checkrun == true) {
110  if (rank == 0) throw std::runtime_error("ERROR: MonteCarlo::Run() cannot be called after calling MonteCarlo::TestRun().\n");
111  else sleep(2);
112  } else {
113  checkrun = true;
114  }
115 
116  try {
117 
118  /* set model parameters */
120  int buffsize = 0;
121  std::map<std::string, double> DP;
122  for (std::vector<ModelParameter>::iterator it = ModPars.begin(); it < ModPars.end(); it++) {
123  if (it->geterrg() > 0. || it->geterrf() > 0.)
124  buffsize++;
125  if (it->IsCorrelated()) {
126  for (unsigned int i = 0; i < CGP.size(); i++) {
127  if (CGP[i].getName().compare(it->getCgp_name()) == 0) {
128  std::string index = it->getname().substr(CGP[i].getName().size());
129  long int lindex = strtol(index.c_str(), NULL, 10);
130  if (lindex > 0)
131  DP[CGP[i].getPar(lindex - 1).getname()] = CGP[i].getPar(lindex - 1).getave();
132  else {
133  std::stringstream out;
134  out << it->getname();
135  throw std::runtime_error("MonteCarlo::Run(): " + out.str() + "seems to be part of a CorrelatedGaussianParameters object, but I couldn't find the corresponding object");
136  }
137  }
138  }
139  } else
140  DP[it->getname()] = it->getave();
141  }
142  if (buffsize == 0){
143  if (rank == 0) throw std::runtime_error("No parameters being varied. Aborting MCMC run.\n");
144  else sleep(2);
145  }
146  buffsize++;
147  if (!myInputParser.getModel()->Init(DP)){
148  if (rank == 0) throw std::runtime_error("\nERROR: Model cannot be initialization.\n");
149  else sleep(2);
150  }
151 
152  if (rank == 0) std::cout << std::endl << "Running in MonteCarlo mode...\n" << std::endl;
153 
154  /* create a directory for outputs */
155  if (rank == 0) {
156  FileStat_t info;
157  if (gSystem->GetPathInfo(ObsDirName.c_str(), info) != 0) {
158  if (gSystem->MakeDirectory(ObsDirName.c_str()) == 0)
159  std::cout << ObsDirName << " directory has been created." << std::endl;
160  else
161  throw std::runtime_error("ERROR: " + ObsDirName + " directory cannot be created.\n");
162  }
163  }
164 
165  MCEngine.SetName(ModelName.c_str());
167 
168 #ifdef _MPI
169  double *recvbuff = new double[buffsize];
170 #endif
171 
172  if (rank != 0) {
173 #ifdef _MPI
174  double **sendbuff = new double*[1];
175  sendbuff[0] = new double[1];
176 
177  std::vector<double> pars;
178  int obsbuffsize = Obs.size() + 2 * Obs2D.size();
179  for (std::vector<CorrelatedGaussianObservables>::iterator it1 = CGO.begin(); it1 < CGO.end(); it1++)
180  obsbuffsize += it1->getObs().size();
181 
182  while (true) {
183  MPI::COMM_WORLD.Scatter(sendbuff[0], buffsize, MPI::DOUBLE,
184  recvbuff, buffsize, MPI::DOUBLE, 0);
185 
186  if (recvbuff[0] == 0.) { // do nothing and return ll
187  double ll = log(0.);
188  MPI::COMM_WORLD.Gather(&ll, 1, MPI::DOUBLE, sendbuff[0], 1, MPI::DOUBLE, 0);
189  } else if (recvbuff[0] == 1.) { //compute log likelihood
190  pars.assign(recvbuff + 1, recvbuff + buffsize);
191  double ll = MCEngine.LogEval(pars);
192  MPI::COMM_WORLD.Gather(&ll, 1, MPI::DOUBLE, sendbuff[0], 1, MPI::DOUBLE, 0);
193  } else if (recvbuff[0] == 2.) { // compute observables
194  double sbuff[obsbuffsize];
195  std::map<std::string, double> DPars;
196  pars.assign(recvbuff + 1, recvbuff + buffsize);
197  MCEngine.setDParsFromParameters(pars,DPars);
198  myInputParser.getModel()->Update(DPars);
199 
200  int k = 0;
201  for (boost::ptr_vector<Observable>::iterator it = Obs.begin(); it < Obs.end(); it++) {
202  sbuff[k++] = it->computeTheoryValue();
203  }
204  for (std::vector<Observable2D>::iterator it = Obs2D.begin(); it < Obs2D.end(); it++) {
205  sbuff[k++] = it->computeTheoryValue();
206  sbuff[k++] = it->computeTheoryValue2();
207  }
208 
209  for (std::vector<CorrelatedGaussianObservables>::iterator it1 = CGO.begin(); it1 < CGO.end(); it1++) {
210  std::vector<Observable> pino(it1->getObs());
211  for (std::vector<Observable>::iterator it = pino.begin(); it != pino.end(); ++it)
212  sbuff[k++] = it->computeTheoryValue();
213  }
214  MPI::COMM_WORLD.Gather(sbuff, obsbuffsize, MPI::DOUBLE, sendbuff[0], obsbuffsize, MPI::DOUBLE, 0);
215  } else if (recvbuff[0] == 3.) { // do not compute observables, but gather the buffer
216  double sbuff[obsbuffsize];
217  MPI::COMM_WORLD.Gather(sbuff, obsbuffsize, MPI::DOUBLE, sendbuff[0], obsbuffsize, MPI::DOUBLE, 0);
218  } else if (recvbuff[0] == -1.)
219  break;
220  else {
221  std::cout << "recvbuff = " << recvbuff[0] << " rank " << rank << std::endl;
222  throw "MonteCarlo::Run(): error in MPI message!\n";
223  }
224  }
225  delete sendbuff[0];
226  delete [] sendbuff;
227 #endif
228  } else {
229 
230  bool writechains = false;
231  std::cout << std::endl;
232  std::cout << std::endl;
233  if (ModPars.size() > 0) std::cout << ModPars.size() << " parameters defined." << std::endl;
234  if (Obs.size() > 0) std::cout << Obs.size() << " observables defined." << std::endl;
235  if (Obs2D.size() > 0) std::cout << Obs2D.size() << " 2D observables defined." << std::endl;
236  if (CGO.size() > 0) std::cout << CGO.size() << " correlated gaussian observables defined";
237  if (CGO.size() == 0)
238  std::cout << "." << std::endl;
239  else
240  std::cout << ":" << std::endl;
241  for (std::vector<CorrelatedGaussianObservables>::iterator it1 = CGO.begin();
242  it1 != CGO.end(); ++it1)
243  std::cout << " " << it1->getName() << " containing "
244  << it1->getObs().size() << " observables." << std::endl;
245  //MonteCarlo configuration parser
246  std::ifstream ifile(MCMCConf.c_str());
247  if (!ifile.is_open())
248  throw std::runtime_error("\nERROR: " + MCMCConf + " does not exist. Make sure to specify a valid Monte Carlo configuration file.\n");
249  std::string line;
250  bool IsEOF = false;
251  do {
252  IsEOF = getline(ifile, line).eof();
253  if (*line.rbegin() == '\r') line.erase(line.length() - 1); // for CR+LF
254  if (line.empty() || line.at(0) == '#')
255  continue;
256  boost::char_separator<char> sep(" \t");
257  boost::tokenizer<boost::char_separator<char> > tok(line, sep);
258  boost::tokenizer<boost::char_separator<char> >::iterator beg = tok.begin();
259  if (beg->compare("NChains") == 0) {
260  ++beg;
261  MCEngine.setNChains(atoi((*beg).c_str()));
262  } else if (beg->compare("PrerunMaxIter") == 0) {
263  ++beg;
264  MCEngine.MCMCSetNIterationsMax(atoi((*beg).c_str()));
265  } else if (beg->compare("NIterationsUpdateMax") == 0) {
266  ++beg;
267  MCEngine.MCMCSetNIterationsUpdateMax(atoi((*beg).c_str()));
268  } else if (beg->compare("Seed") == 0) {
269  ++beg;
270  int seed = atoi((*beg).c_str());
271  if (seed != 0)
272  MCEngine.MCMCSetRandomSeed(seed);
273  } else if (beg->compare("Iterations") == 0) {
274  ++beg;
275  MCEngine.MCMCSetNIterationsRun(atoi((*beg).c_str()));
276  } else if (beg->compare("MinimumEfficiency") == 0) {
277  ++beg;
278  MCEngine.MCMCSetMinimumEfficiency(atof((*beg).c_str()));
279  } else if (beg->compare("WriteChain") == 0) {
280  ++beg;
281  if (beg->compare("true") == 0)
282  writechains = true;
283  } else if (beg->compare("FindModeWithMinuit") == 0) {
284  ++beg;
285  if (beg->compare("true") == 0) {
286  FindModeWithMinuit = true;
287  }
288  } else if (beg->compare("CalculateNormalization") == 0) {
289  ++beg;
290  if (beg->compare("true") == 0) {
291  CalculateNormalization = true;
292  }
293  } else if (beg->compare("PrintAllMarginalized") == 0) {
294  ++beg;
295  if (beg->compare("true") == 0) {
296  PrintAllMarginalized = true;
297  }
298  } else if (beg->compare("PrintCorrelationMatrix") == 0) {
299  ++beg;
300  if (beg->compare("true") == 0) {
301  PrintCorrelationMatrix = true;
302  }
303  } else if (beg->compare("PrintKnowledgeUpdatePlots") == 0) {
304  ++beg;
305  if (beg->compare("true") == 0) {
307  }
308  } else if (beg->compare("PrintParameterPlot") == 0) {
309  ++beg;
310  if (beg->compare("true") == 0) {
311  PrintParameterPlot = true;
312  }
313  } else if (beg->compare("WritePreRunData") == 0) {
314  ++beg;
315  if (beg->compare("true") == 0) {
316  WritePreRunData = true;
317  }
318  } else if (beg->compare("ReadPreRunData") == 0) {
319  ++beg;
320  ReadPreRunData(*beg);
321  } else if (beg->compare("OrderParameters") == 0) {
322  ++beg;
323  if (beg->compare("false") == 0) {
324  MCEngine.MCMCSetFlagOrderParameters(false);
325  }
326  } else if (beg->compare("StatOverFlow") == 0) {
327  ++beg;
328  if (beg->compare("true") == 0) MCEngine.setHistogramOverFlow(true);
329  if (beg->compare("false") == 0) MCEngine.setHistogramOverFlow(false);
330  else throw std::runtime_error("\nERROR: Wrong value in MonteCarlo config file for keyword StatOverFlow. It can be either true or false.\n");
331  } else
332  throw std::runtime_error("\nERROR: Wrong keyword in MonteCarlo config file: " + *beg + "\n Make sure to specify a valid Monte Carlo configuration file.\n");
333  } while (!IsEOF);
334 
335  BCModelOutput out(&MCEngine, OutFile.c_str());
336  if (writechains) {
337  out.WriteMarkovChain(true);
339  }
340 
341  // set nicer style for drawing than the ROOT default
342  BCAux::SetStyle();
343 
344  // open log file
345  BCLog::OpenLog(("log" + JobTag + ".txt").c_str());
346  BCLog::SetLogLevel(BCLog::debug);
347 
348  // run the MCMC and marginalize w.r.t. to all parameters
349  MCEngine.BCIntegrate::SetNbins(NBINSMODELPARS);
350  MCEngine.SetMarginalizationMethod(BCIntegrate::kMargMetropolis);
351  std::time_t ti = std::time(NULL);
352  char mbstr[100];
353  if (std::strftime(mbstr, sizeof(mbstr), "%H:%M:%S %d/%m/%y ", std::localtime(&ti)))
354  std::cout << "MCMC Run started at " << mbstr << std::endl;
355  MCEngine.MarginalizeAll();
356  std::time_t tf = std::time(NULL);
357  if (std::strftime(mbstr, sizeof(mbstr), "%H:%M:%S %d/%m/%y ", std::localtime(&tf)))
358  std::cout << "MCMC Run ended at " << mbstr << std::endl;
359 
360  // find mode using the best fit parameters as start values
361  if (FindModeWithMinuit)
362  MCEngine.FindMode(MCEngine.GetBestFitParameters());
363 
365 
366  // draw all marginalized distributions into a pdf file
368  MCEngine.PrintAllMarginalized(("MonteCarlo_plots" + JobTag + ".pdf").c_str());
369 
370  // print results of the analysis into a text file
371  MCEngine.PrintResults(("MonteCarlo_results" + JobTag + ".txt").c_str());
372 
373  // print histograms
375 
376  BCSummaryTool myBCSummaryTool(&MCEngine);
377 
378  // draw the correlation matrix into a pdf file
380  myBCSummaryTool.PrintCorrelationMatrix(("ParamCorrelations" + JobTag + ".pdf").c_str());
381 
382  // print the correlation matrix into a tex file
384  MCEngine.PrintCorrelationMatrix(("ParamCorrelations" + JobTag + ".tex").c_str());
385 
386  // print comparisons of the prior knowledge to the posterior knowledge
387  // for all parameters into a pdf file
389  myBCSummaryTool.PrintKnowledgeUpdatePlots(("ParamUpdate" + JobTag + ".pdf").c_str());
390 
391  // draw an overview plot of the parameters into a pdf file
392  if (PrintParameterPlot)
393  myBCSummaryTool.PrintParameterPlot(("ParamSummary" + JobTag + ".pdf").c_str());
394 
395  // print a LaTeX table of the parameters into a tex file
396  //myBCSummaryTool.PrintParameterLatex(("ParamSummary" + JobTag + ".tex").c_str());
397 
398  out.WriteMarginalizedDistributions();
399  out.FillAnalysisTree();
400  out.Close();
401 
402  // print logs for the histograms of the observables into a text file
403  std::ofstream outHistoLog;
404  outHistoLog.open((ObsDirName + "/HistoLog" + JobTag + ".txt").c_str(), std::ios::out);
405  outHistoLog << MCEngine.getHistoLog();
406  outHistoLog.close();
407 
408  // print statistics for the theory values of the observables into a text file
409  std::ofstream outStatLog;
410  outStatLog.open((ObsDirName + "/Statistics" + JobTag + ".txt").c_str(), std::ios::out);
411  if (CalculateNormalization) outStatLog << "Normalization for " << ModelName.c_str() << ": " << normalization << "\n" << std::endl;
412  outStatLog << MCEngine.computeStatistics();
413  outStatLog.close();
414 
415  // print global mode and scale factors for the 1st chain into a text file
416  if (WritePreRunData) {
417  std::ofstream outPreRun;
418  outPreRun.open(("PreRun" + JobTag + ".txt").c_str(), std::ios::out);
419  outPreRun << MCEngine.writePreRunData();
420  outPreRun.close();
421  }
422 
423  /* Number of events */
424  std::stringstream ss;
425  ss << "Number of used events: " << MCEngine.getNumOfUsedEvents();
426  BCLog::OutSummary(ss.str().c_str());
427  ss.str("");
428  ss << "Number of discarded events: " << MCEngine.getNumOfDiscardedEvents();
429  BCLog::OutSummary(ss.str().c_str());
430 
431  // close log file
432  BCLog::CloseLog();
433 
434 #ifdef _MPI
435  double ** sendbuff = new double *[MCEngine.procnum];
436  sendbuff[0] = new double[MCEngine.procnum * buffsize];
437  for (int il = 1; il < MCEngine.procnum; il++) {
438  sendbuff[il] = sendbuff[il - 1] + buffsize;
439  sendbuff[il][0] = -1.; //Exit command
440  }
441  MPI::COMM_WORLD.Scatter(sendbuff[0], buffsize, MPI::DOUBLE,
442  recvbuff, buffsize, MPI::DOUBLE, 0);
443  delete sendbuff[0];
444  delete [] sendbuff;
445 #endif
446  }
447 #ifdef _MPI
448  delete [] recvbuff;
449 #endif
450  } catch (std::string message) {
451  std::cerr << message << std::endl;
452  exit(EXIT_FAILURE);
453  }
454 }
bool checkrun
A check to make sure TestRun()and Run() are not called consecutively.
Definition: MonteCarlo.h:175
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:28
double normalization
A variable to store the evidence of a model.
Definition: MonteCarlo.h:176
virtual bool Update(const std::map< std::string, double > &DPars)
The update method for StandardModel.
std::string ObsDirName
String for the output directory name.
Definition: MonteCarlo.h:167
std::string computeStatistics()
A get method to compute the mean and rms of the computed observables.
std::string getHistoLog() const
A get method to access the stream that stores the log messages coming from histogram printing and che...
int getNumOfDiscardedEvents() const
A get method to access the number of events discarded due to failure to update model. These events are not used for the MCMC run.
bool CalculateNormalization
< Flag for calculating the evidence.
Definition: MonteCarlo.h:169
void ReadPreRunData(std::string file)
Definition: MonteCarlo.cpp:456
bool PrintParameterPlot
Flag for printing the overview parameter plots.
Definition: MonteCarlo.h:173
MonteCarloEngine MCEngine
An object of the MonteCarloEngine class.
Definition: MonteCarlo.h:157
void PrintHistogram(BCModelOutput &out, Observable &it, const std::string OutputDir)
Overloaded from PrintHistogram(BCModelOutput&, const std::string) to print histogram for observables...
std::string OutFile
String for the name of the output root file without the .root extension.
Definition: MonteCarlo.h:165
bool PrintAllMarginalized
Flag for printing all Marginalized distributions to be passed on to the BAT routines.
Definition: MonteCarlo.h:170
bool WritePreRunData
Flag for printing the overview parameter plots.
Definition: MonteCarlo.h:174
void AddChains()
A method to add the observable values and weights to the chain information.
std::vector< CorrelatedGaussianParameters > CGP
Vector for the Correlated Gaussian Parameters defined in SomeModel.conf.
Definition: MonteCarlo.h:162
virtual bool Init(const std::map< std::string, double > &DPars)
A method to initialize the model parameters.
std::vector< ModelParameter > ModPars
Vector for the model parameters defined in SomeModel.conf.
Definition: MonteCarlo.h:158
void setDParsFromParameters(const std::vector< double > &parameters, std::map< std::string, double > &DPars_i)
bool PrintKnowledgeUpdatePlots
Flag for printing plots to compare prior vs. posterior knowledge of parameters.
Definition: MonteCarlo.h:172
#define NBINSMODELPARS
bool PrintCorrelationMatrix
Flag for printing the correlation matrix.
Definition: MonteCarlo.h:171
std::string ModelName
The name of the model.
Definition: MonteCarlo.h:155
void setNChains(unsigned int i)
A set method to fix the number of chains.
std::string writePreRunData()
A method to write in a text file the best fit parameters and the prerun scale factors.
std::vector< Observable2D > Obs2D
Vector for the Observables2D defined in SomeModel.conf.
Definition: MonteCarlo.h:160
void setHistogramOverFlow(bool overflow)
boost::ptr_vector< Observable > Obs
Vector for the observables defined in SomeModel.conf.
Definition: MonteCarlo.h:159
void PrintCorrelationMatrix(const std::string filename)
This member generates the correlation matrix using BCH2D from the BAT libraries.
complex log(const complex &z)
std::vector< CorrelatedGaussianObservables > CGO
Vector for the Correlated Gaussian Observables defined in SomeModel.conf.
Definition: MonteCarlo.h:161
InputParser myInputParser
An object of the InputParser class.
Definition: MonteCarlo.h:156
bool FindModeWithMinuit
Flag for using Minuit libraries.
Definition: MonteCarlo.h:168
void Initialize(Model *Mod_i)
Initialization of the Monte Carlo Engine.
StandardModel * getModel() const
A get method to access the pointer to the object of the StandardModel class.
Definition: InputParser.h:110
std::string MCMCConf
String for the name of the MonteCarlo.conf file.
Definition: MonteCarlo.h:164
int getNumOfUsedEvents() const
A get method to access the number of events used in the MCMC run.
std::string JobTag
String for the optional JobTag argument to be passes to the executable.
Definition: MonteCarlo.h:166
std::string ModelConf
String for the name of the SomeModel.conf file.
Definition: MonteCarlo.h:163
void MonteCarlo::TestRun ( int  rank)

The default destructor.

This member is used for test runa to generate a single event

The values of the parameters that are used to generate the single event is the average value passed from the SomeModel.conf file.

The algorithm implemented by this member is as follows:

  • Initiate InputParser::ReadParameters() which read the SomeModel.conf file for setting the model parameters and the observables to be generated. This call also passes on the name of the model to the private member ModelName.
  • Map the model parameter mean values with the map DP and calculate buffsize which can be used for implementing MPI runs. The variable buffsize is incremented only for those model parameters that are varied in the Monte Carlo run.
  • Initiate the model using the map DP checking for consistency between required model parameters and the ones supplied in the SomeModel.conf file
  • Generate a single set of observables and correlated Gaussian observables as listed in the SomeModel.conf file. The output is passes to the standard output.
Parameters
[in]rank= MPI::COMM_WORLD.Get_rank(), specifies the rank of the process. This carries a non zero value only when the executable is compiled with the parallelized version of BAT and run as parallel processes with MPI.

Definition at line 52 of file MonteCarlo.cpp.

52  {
53  if (checkrun == true) {
54  if (rank == 0) throw std::runtime_error("ERROR: MonteCarlo::TestRun() cannot be called after calling MonteCarlo::Run().\n");
55  else sleep (2);
56  } else {
57  checkrun = true;
58  }
59 
60  if (rank == 0) {
62  std::map<std::string, double> DP;
63  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");
64 
65  for (std::vector<ModelParameter>::iterator it = ModPars.begin(); it < ModPars.end(); it++) {
66  if (it->IsCorrelated()) {
67  for (unsigned int i = 0; i < CGP.size(); i++) {
68  if (CGP[i].getName().compare(it->getCgp_name()) == 0) {
69  std::string index = it->getname().substr(CGP[i].getName().size());
70  long int lindex = strtol(index.c_str(), NULL, 10);
71  if (lindex > 0)
72  DP[CGP[i].getPar(lindex - 1).getname()] = CGP[i].getPar(lindex - 1).getave();
73  else {
74  std::stringstream out;
75  out << it->getname();
76  throw std::runtime_error("MonteCarlo::Run(): " + out.str() + "seems to be part of a CorrelatedGaussianParameters object, but I couldn't find the corresponding object");
77  }
78  }
79  }
80  } else
81  DP[it->getname()] = it->getave();
82  }
83 
84  if (!myInputParser.getModel()->Init(DP)) {
85  if (rank == 0) throw std::runtime_error("ERROR: Parameter(s) missing in model initialization. \n");
86  else sleep (2);
87  }
88 
89  if (Obs.size() > 0) std::cout << "\nOservables: \n" << std::endl;
90  for (boost::ptr_vector<Observable>::iterator it = Obs.begin(); it < Obs.end(); it++) {
91  double th = it->computeTheoryValue();
92  std::cout << it->getName() << " = " << th << std::endl;
93  }
94 
95  if (CGO.size() > 0) std::cout << "\nCorrelated Gaussian Oservables: \n" << std::endl;
96  for (std::vector<CorrelatedGaussianObservables>::iterator it1 = CGO.begin(); it1 < CGO.end(); it1++) {
97  std::cout << it1->getName() << ":" << std::endl;
98  std::vector<Observable> ObsInCGO = it1->getObs();
99  for (std::vector<Observable>::iterator it2 = ObsInCGO.begin(); it2 < ObsInCGO.end(); it2++) {
100  double th = it2->computeTheoryValue();
101  std::cout << it2->getName() << " = " << th << std::endl;
102  }
103  std::cout << std::endl;
104  }
105  }
106 }
bool checkrun
A check to make sure TestRun()and Run() are not called consecutively.
Definition: MonteCarlo.h:175
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:28
std::vector< CorrelatedGaussianParameters > CGP
Vector for the Correlated Gaussian Parameters defined in SomeModel.conf.
Definition: MonteCarlo.h:162
virtual bool Init(const std::map< std::string, double > &DPars)
A method to initialize the model parameters.
std::vector< ModelParameter > ModPars
Vector for the model parameters defined in SomeModel.conf.
Definition: MonteCarlo.h:158
std::string ModelName
The name of the model.
Definition: MonteCarlo.h:155
std::vector< Observable2D > Obs2D
Vector for the Observables2D defined in SomeModel.conf.
Definition: MonteCarlo.h:160
boost::ptr_vector< Observable > Obs
Vector for the observables defined in SomeModel.conf.
Definition: MonteCarlo.h:159
std::vector< CorrelatedGaussianObservables > CGO
Vector for the Correlated Gaussian Observables defined in SomeModel.conf.
Definition: MonteCarlo.h:161
InputParser myInputParser
An object of the InputParser class.
Definition: MonteCarlo.h:156
StandardModel * getModel() const
A get method to access the pointer to the object of the StandardModel class.
Definition: InputParser.h:110
std::string ModelConf
String for the name of the SomeModel.conf file.
Definition: MonteCarlo.h:163

Member Data Documentation

bool MonteCarlo::CalculateNormalization
private

< Flag for calculating the evidence.

Definition at line 169 of file MonteCarlo.h.

std::vector<CorrelatedGaussianObservables> MonteCarlo::CGO
private

Vector for the Correlated Gaussian Observables defined in SomeModel.conf.

Definition at line 161 of file MonteCarlo.h.

std::vector<CorrelatedGaussianParameters> MonteCarlo::CGP
private

Vector for the Correlated Gaussian Parameters defined in SomeModel.conf.

Definition at line 162 of file MonteCarlo.h.

bool MonteCarlo::checkrun
private

A check to make sure TestRun()and Run() are not called consecutively.

Definition at line 175 of file MonteCarlo.h.

bool MonteCarlo::FindModeWithMinuit
private

Flag for using Minuit libraries.

Definition at line 168 of file MonteCarlo.h.

std::string MonteCarlo::JobTag
private

String for the optional JobTag argument to be passes to the executable.

Definition at line 166 of file MonteCarlo.h.

MonteCarloEngine MonteCarlo::MCEngine
private

An object of the MonteCarloEngine class.

Definition at line 157 of file MonteCarlo.h.

std::string MonteCarlo::MCMCConf
private

String for the name of the MonteCarlo.conf file.

Definition at line 164 of file MonteCarlo.h.

std::string MonteCarlo::ModelConf
private

String for the name of the SomeModel.conf file.

Definition at line 163 of file MonteCarlo.h.

std::string MonteCarlo::ModelName
private

The name of the model.

Definition at line 155 of file MonteCarlo.h.

std::vector<ModelParameter> MonteCarlo::ModPars
private

Vector for the model parameters defined in SomeModel.conf.

Definition at line 158 of file MonteCarlo.h.

InputParser MonteCarlo::myInputParser
private

An object of the InputParser class.

Definition at line 156 of file MonteCarlo.h.

double MonteCarlo::normalization
private

A variable to store the evidence of a model.

Definition at line 176 of file MonteCarlo.h.

boost::ptr_vector<Observable> MonteCarlo::Obs
private

Vector for the observables defined in SomeModel.conf.

Definition at line 159 of file MonteCarlo.h.

std::vector<Observable2D> MonteCarlo::Obs2D
private

Vector for the Observables2D defined in SomeModel.conf.

Definition at line 160 of file MonteCarlo.h.

std::string MonteCarlo::ObsDirName
private

String for the output directory name.

Definition at line 167 of file MonteCarlo.h.

std::string MonteCarlo::OutFile
private

String for the name of the output root file without the .root extension.

Definition at line 165 of file MonteCarlo.h.

bool MonteCarlo::PrintAllMarginalized
private

Flag for printing all Marginalized distributions to be passed on to the BAT routines.

Definition at line 170 of file MonteCarlo.h.

bool MonteCarlo::PrintCorrelationMatrix
private

Flag for printing the correlation matrix.

Definition at line 171 of file MonteCarlo.h.

bool MonteCarlo::PrintKnowledgeUpdatePlots
private

Flag for printing plots to compare prior vs. posterior knowledge of parameters.

Definition at line 172 of file MonteCarlo.h.

bool MonteCarlo::PrintParameterPlot
private

Flag for printing the overview parameter plots.

Definition at line 173 of file MonteCarlo.h.

bool MonteCarlo::WritePreRunData
private

Flag for printing the overview parameter plots.

Definition at line 174 of file MonteCarlo.h.


The documentation for this class was generated from the following files: