A class for handling numerical operations for models. More...
#include <BCIntegrate.h>
Inherits BCEngineMCMC.
Inherited by BCModel.
Detailed Description
A class for handling numerical operations for models.
- Version:
- 1.0
- Date:
- 08.2008 This is a base class for a model class. It contains numerical methods to carry out the integration, marginalization, peak finding etc.
Definition at line 45 of file BCIntegrate.h.
Member Enumeration Documentation
An enumerator for Cuba integration method
Definition at line 71 of file BCIntegrate.h.
{ kCubaVegas, kCubaSuave, kCubaDivonne, kCubaCuhre };
An enumerator for the integration algorithm
Definition at line 55 of file BCIntegrate.h.
{ kIntMonteCarlo, kIntImportance, kIntMetropolis, kIntCuba, NIntMethods };
An enumerator for the marginalization algorithm
Definition at line 59 of file BCIntegrate.h.
{ kMargMonteCarlo, kMargMetropolis, NMargMethods };
An enumerator for the mode finding algorithm
Definition at line 63 of file BCIntegrate.h.
{ kOptSA, kOptMetropolis, kOptMinuit, NOptMethods };
An enumerator for the Simulated Annealing schedule
Definition at line 67 of file BCIntegrate.h.
{ kSACauchy, kSABoltzmann, kSACustom, NSAMethods };
Constructor & Destructor Documentation
| BCIntegrate::BCIntegrate | ( | ) |
The default constructor
Definition at line 31 of file BCIntegrate.cxx.
: BCEngineMCMC() , fNvar(0) , fNbins(100) , fNSamplesPer2DBin(100) , fMarkovChainStepSize(0.1) , fMarkovChainAutoN(true) , fDataPointLowerBoundaries(0) , fDataPointUpperBoundaries(0) , fFillErrorBand(false) , fFitFunctionIndexX(-1) , fFitFunctionIndexY(-1) , fErrorBandContinuous(true) , fErrorBandNbinsX(100) , fErrorBandNbinsY(500) , fMinuit(0) , fFlagIgnorePrevOptimization(false) , fFlagWriteMarkovChain(false) , fMarkovChainTree(0) , fSAT0(100) , fSATmin(0.1) , fTreeSA(0) , fFlagWriteSAToFile(false) , fNiterPerDimension(100) #ifdef HAVE_CUBA_H , fIntegrationMethod(BCIntegrate::kIntCuba) #else , fIntegrationMethod(BCIntegrate::kIntMonteCarlo) #endif , fMarginalizationMethod(BCIntegrate::kMargMetropolis) , fOptimizationMethod(BCIntegrate::kOptMinuit) , fOptimizationMethodMode(BCIntegrate::kOptMinuit) , fSASchedule(BCIntegrate::kSACauchy) , fNIterationsMax(1000000) , fNIterations(0) , fRelativePrecision(1e-3) , fAbsolutePrecision(1e-12) , fCubaIntegrationMethod(BCIntegrate::kCubaVegas) , fCubaMinEval(0) , fCubaMaxEval(2000000) , fCubaVerbosity(0) , fCubaVegasNStart(25000) , fCubaVegasNIncrease(25000) , fCubaSuaveNNew(10000) , fCubaSuaveFlatness(50) , fError(-999.) { fMinuitArglist[0] = 20000; fMinuitArglist[1] = 0.01; }
| BCIntegrate::BCIntegrate | ( | BCParameterSet * | par | ) |
A constructor
Definition at line 84 of file BCIntegrate.cxx.
: BCEngineMCMC() , fNvar(0) , fNbins(100) , fNSamplesPer2DBin(100) , fMarkovChainStepSize(0.1) , fMarkovChainAutoN(true) , fDataPointLowerBoundaries(0) , fDataPointUpperBoundaries(0) , fFillErrorBand(false) , fFitFunctionIndexX(-1) , fFitFunctionIndexY(-1) , fErrorBandContinuous(true) , fErrorBandNbinsX(100) , fErrorBandNbinsY(500) , fMinuit(0) , fFlagIgnorePrevOptimization(false) , fFlagWriteMarkovChain(false) , fMarkovChainTree(0) , fSAT0(100) , fSATmin(0.1) , fTreeSA(0) , fFlagWriteSAToFile(false) , fNiterPerDimension(100) #ifdef HAVE_CUBA_H , fIntegrationMethod(BCIntegrate::kIntCuba) #else , fIntegrationMethod(BCIntegrate::kIntMonteCarlo) #endif , fMarginalizationMethod(BCIntegrate::kMargMetropolis) , fOptimizationMethod(BCIntegrate::kOptMinuit) , fOptimizationMethodMode(BCIntegrate::kOptMinuit) , fSASchedule(BCIntegrate::kSACauchy) , fNIterationsMax(1000000) , fNIterations(0) , fRelativePrecision(1e-3) , fAbsolutePrecision(1e-12) , fCubaIntegrationMethod(BCIntegrate::kCubaVegas) , fCubaMinEval(0) , fCubaMaxEval(2000000) , fCubaVerbosity(0) , fCubaVegasNStart(25000) , fCubaVegasNIncrease(25000) , fCubaSuaveNNew(10000) , fCubaSuaveFlatness(50) , fError(-999.) { SetParameters(par); fMinuitArglist[0] = 20000; fMinuitArglist[1] = 0.01; }
| BCIntegrate::BCIntegrate | ( | const BCIntegrate & | bcintegrate | ) |
The copy constructor
Definition at line 139 of file BCIntegrate.cxx.
: BCEngineMCMC(bcintegrate) { fNvar = bcintegrate.fNvar; fNbins = bcintegrate.fNbins; fNSamplesPer2DBin = bcintegrate.fNSamplesPer2DBin; fMarkovChainStepSize = bcintegrate.fMarkovChainStepSize; fMarkovChainNIterations = bcintegrate.fMarkovChainNIterations; fMarkovChainAutoN = bcintegrate.fMarkovChainAutoN; if (bcintegrate.fDataPointLowerBoundaries) fDataPointLowerBoundaries = new BCDataPoint(*bcintegrate.fDataPointLowerBoundaries); else fDataPointLowerBoundaries = 0; if (bcintegrate.fDataPointUpperBoundaries) fDataPointUpperBoundaries = new BCDataPoint(*bcintegrate.fDataPointUpperBoundaries); else fDataPointUpperBoundaries = 0; fDataFixedValues = bcintegrate.fDataFixedValues; fBestFitParameters = bcintegrate.fBestFitParameters; fBestFitParameterErrors = bcintegrate.fBestFitParameterErrors; fBestFitParametersMarginalized = bcintegrate.fBestFitParametersMarginalized; for (int i = 0; i < int(bcintegrate.fHProb1D.size()); ++i) { if (bcintegrate.fHProb1D.at(i)) fHProb1D.push_back(new TH1D(*(bcintegrate.fHProb1D.at(i)))); else fHProb1D.push_back(0); } for (int i = 0; i < int(bcintegrate.fHProb2D.size()); ++i) { if (bcintegrate.fHProb2D.at(i)) fHProb2D.push_back(new TH2D(*(fHProb2D.at(i)))); else fHProb2D.push_back(0); } fFillErrorBand = bcintegrate.fFillErrorBand; fFitFunctionIndexX = bcintegrate.fFitFunctionIndexX; fFitFunctionIndexY = bcintegrate.fFitFunctionIndexY; fErrorBandX = bcintegrate.fErrorBandX; if (bcintegrate.fErrorBandXY) fErrorBandXY = new TH2D(*(bcintegrate.fErrorBandXY)); else fErrorBandXY = 0; fErrorBandNbinsX = bcintegrate.fErrorBandNbinsX; fErrorBandNbinsY = bcintegrate.fErrorBandNbinsY; fMinuit = new TMinuit(); // debugKK // *fMinuit = *(bcintegrate.fMinuit); fMinuitArglist[0] = bcintegrate.fMinuitArglist[0]; fMinuitArglist[1] = bcintegrate.fMinuitArglist[1]; fMinuitErrorFlag = bcintegrate.fMinuitErrorFlag; fFlagIgnorePrevOptimization = bcintegrate.fFlagIgnorePrevOptimization; fFlagWriteMarkovChain = bcintegrate.fFlagWriteMarkovChain; fMarkovChainTree = bcintegrate.fMarkovChainTree; fMCMCIteration = bcintegrate.fMCMCIteration; fSAT0 = bcintegrate.fSAT0; fSATmin = bcintegrate.fSATmin; // debugKK fTreeSA = 0; fFlagWriteSAToFile = bcintegrate.fFlagWriteSAToFile; fSANIterations = bcintegrate.fSANIterations; fSATemperature = bcintegrate.fSATemperature; fSALogProb = bcintegrate.fSALogProb; fSAx = bcintegrate.fSAx; if (bcintegrate.fx) fx = new BCParameterSet(*(bcintegrate.fx)); else fx = 0; fMin = new double[fNvar]; fMax = new double[fNvar]; fVarlist = new int[fNvar]; fMin = bcintegrate.fMin; fMax = bcintegrate.fMax; fVarlist = bcintegrate.fVarlist; fNiterPerDimension = bcintegrate.fNiterPerDimension; fIntegrationMethod = bcintegrate.fIntegrationMethod; fMarginalizationMethod = bcintegrate.fMarginalizationMethod; fOptimizationMethod = bcintegrate.fOptimizationMethod; fOptimizationMethodMode = bcintegrate.fOptimizationMethodMode; fSASchedule = bcintegrate.fSASchedule; fNIterationsMax = bcintegrate.fNIterationsMax; fNIterations = bcintegrate.fNIterations; fRelativePrecision = bcintegrate.fRelativePrecision; fAbsolutePrecision = bcintegrate.fAbsolutePrecision; fCubaIntegrationMethod = bcintegrate.fCubaIntegrationMethod; fCubaMinEval = bcintegrate.fCubaMinEval; fCubaMaxEval = bcintegrate.fCubaMaxEval; fCubaVerbosity = bcintegrate.fCubaVerbosity; fCubaVegasNStart = bcintegrate.fCubaVegasNStart; fCubaVegasNIncrease = bcintegrate.fCubaVegasNIncrease; fCubaSuaveNNew = bcintegrate.fCubaSuaveNNew; fCubaSuaveFlatness = bcintegrate.fCubaSuaveFlatness; fError = bcintegrate.fError; fNmetro = bcintegrate.fNmetro; fNacceptedMCMC = bcintegrate.fNacceptedMCMC; fXmetro0 = bcintegrate.fXmetro0; fXmetro1 = bcintegrate.fXmetro1; fMarkovChainValue = bcintegrate.fMarkovChainValue; }
| BCIntegrate::~BCIntegrate | ( | ) | [virtual] |
The default destructor
Definition at line 338 of file BCIntegrate.cxx.
Member Function Documentation
| int BCIntegrate::CubaIntegrand | ( | const int * | ndim, | |
| const double | xx[], | |||
| const int * | ncomp, | |||
| double | ff[], | |||
| void * | userdata | |||
| ) | [static] |
Integrand for the Cuba library.
- Parameters:
-
ndim The number of dimensions to integrate over xx The point in parameter space to integrate over (scaled to 0 - 1 per dimension) ncomp The number of components of the integrand (usually 1) ff The function value
- Returns:
- The integral
Definition at line 2025 of file BCIntegrate.cxx.
{
#ifdef HAVE_CUBA_H
// scale variables
double jacobian = 1.0;
std::vector<double> scaled_parameters;
for (int i = 0; i < *ndim; i++) {
double range = global_this->fx->at(i)->GetUpperLimit() - global_this->fx->at(i)->GetLowerLimit();
// multiply range to jacobian
jacobian *= range;
// get the scaled parameter value
scaled_parameters.push_back(global_this->fx->at(i)->GetLowerLimit() + xx[i] * range);
}
// call function to integrate
ff[0] = global_this->Eval(scaled_parameters);
// multiply jacobian
ff[0] *= jacobian;
// multiply fudge factor
ff[0] *= 1e99;
// remove parameter vector
scaled_parameters.clear();
#else
BCLog::OutError("!!! This version of BAT is compiled without Cuba.");
BCLog::OutError(" Use other integration methods or install Cuba and recompile BAT.");
#endif
return 0;
}
| double BCIntegrate::CubaIntegrate | ( | BCIntegrate::BCCubaMethod | method, | |
| std::vector< double > | parameters_double, | |||
| std::vector< double > | parameters_int | |||
| ) |
Calculate integral using the Cuba library. For details see documentation.
- Parameters:
-
method A short cut for the method parameters_double A vector of parameters (double) parameters_int A vector of parameters (int)
- Returns:
- The integral
Definition at line 2106 of file BCIntegrate.cxx.
{
#ifdef HAVE_CUBA_H
const int NDIM = int(fx ->size());
const int NCOMP = 1;
const int USERDATA = 0;
const int SEED = 0;
const int NBATCH = 1000;
const int GRIDNO = -1;
const char*STATEFILE = "";
const double EPSREL = parameters_double[0];
const double EPSABS = parameters_double[1];
const int VERBOSE = int(parameters_int[0]);
const int MINEVAL = int(parameters_int[1]);
const int MAXEVAL = int(parameters_int[2]);
int neval;
int fail;
int nregions;
double integral[NCOMP];
double error[NCOMP];
double prob[NCOMP];
global_this = this;
integrand_t an_integrand = &BCIntegrate::CubaIntegrand;
if (method == 0) {
// set VEGAS specific parameters
const int NSTART = int(parameters_int[3]);
const int NINCREASE = int(parameters_int[4]);
// call VEGAS integration method
Vegas(NDIM, NCOMP, an_integrand, USERDATA,
EPSREL, EPSABS, VERBOSE, SEED,
MINEVAL, MAXEVAL, NSTART, NINCREASE, NBATCH,
GRIDNO, STATEFILE,
&neval, &fail, integral, error, prob);
// interface for Cuba version 1.5
/*
Vegas(NDIM, NCOMP, an_integrand,
EPSREL, EPSABS, VERBOSE, MINEVAL, MAXEVAL,
NSTART, NINCREASE,
&neval, &fail, integral, error, prob);
*/
}
else if (method == 1) {
// set SUAVE specific parameters
// const int LAST = int(parameters_int[3]);
const int NNEW = int(parameters_int[3]);
const double FLATNESS = parameters_double[2];
// call SUAVE integration method
Suave(NDIM, NCOMP, an_integrand,
USERDATA, EPSREL, EPSABS, VERBOSE, SEED,
MINEVAL, MAXEVAL,
NNEW, FLATNESS,
&nregions, &neval, &fail, integral, error, prob);
// interface for Cuba version 1.5
/*
Suave(NDIM, NCOMP, an_integrand,
EPSREL, EPSABS, VERBOSE | LAST, MINEVAL, MAXEVAL,
NNEW, FLATNESS,
&nregions, &neval, &fail, integral, error, prob);
*/
}
else if (method == 2) {
// set DIVONNE specific parameters
const int KEY1 = int(parameters_int[3]);
const int KEY2 = int(parameters_int[4]);
const int KEY3 = int(parameters_int[5]);
const int MAXPASS = int(parameters_int[6]);
const int BORDER = int(parameters_int[7]);
const int MAXCHISQ = int(parameters_int[8]);
const int MINDEVIATION = int(parameters_int[9]);
const int NGIVEN = int(parameters_int[10]);
const int LDXGIVEN = int(parameters_int[11]);
const int NEXTRA = int(parameters_int[12]);
const int FLAGS = 0;
// call DIVONNE integration method
Divonne(NDIM, NCOMP, an_integrand, USERDATA,
EPSREL, EPSABS, FLAGS, SEED, MINEVAL, MAXEVAL,
KEY1, KEY2, KEY3, MAXPASS, BORDER, MAXCHISQ, MINDEVIATION,
NGIVEN, LDXGIVEN, NULL, NEXTRA, NULL,
&nregions, &neval, &fail, integral, error, prob);
// interface for Cuba version 1.5
/*
Divonne(NDIM, NCOMP, an_integrand,
EPSREL, EPSABS, VERBOSE, MINEVAL, MAXEVAL,
KEY1, KEY2, KEY3, MAXPASS, BORDER, MAXCHISQ, MINDEVIATION,
NGIVEN, LDXGIVEN, NULL, NEXTRA, NULL,
&nregions, &neval, &fail, integral, error, prob);
*/
}
else if (method == 3) {
// set CUHRE specific parameters
// const int LAST = int(parameters_int[3]);
const int KEY = int(parameters_int[4]);
const int FLAGS = 0;
// call CUHRE integration method
Cuhre(NDIM, NCOMP, an_integrand, USERDATA,
EPSREL, EPSABS, FLAGS, MINEVAL, MAXEVAL, KEY,
&nregions, &neval, &fail, integral, error, prob);
// interface for Cuba version 1.5
/*
Cuhre(NDIM, NCOMP, an_integrand,
EPSREL, EPSABS, VERBOSE | LAST, MINEVAL, MAXEVAL, KEY,
&nregions, &neval, &fail, integral, error, prob);
*/
}
else {
std::cout << " Integration method not available. " << std::endl;
integral[0] = -1e99;
}
if (fail != 0) {
std::cout << " Warning, integral did not converge with the given set of parameters. "<< std::endl;
std::cout << " neval = " << neval << std::endl;
std::cout << " fail = " << fail << std::endl;
std::cout << " integral = " << integral[0] << std::endl;
std::cout << " error = " << error[0] << std::endl;
std::cout << " prob = " << prob[0] << std::endl;
}
return integral[0] / 1e99;
#else
BCLog::OutError("!!! This version of BAT is compiled without Cuba.");
BCLog::OutError(" Use other integration methods or install Cuba and recompile BAT.");
return 0.;
#endif
}
| double BCIntegrate::CubaIntegrate | ( | ) |
Definition at line 2064 of file BCIntegrate.cxx.
{
#ifdef HAVE_CUBA_H
std::vector<double> parameters_double;
std::vector<double> parameters_int;
parameters_double.push_back(fRelativePrecision);
parameters_double.push_back(fAbsolutePrecision);
parameters_int.push_back(fCubaVerbosity);
parameters_int.push_back(fCubaMinEval);
parameters_int.push_back(fCubaMaxEval);
switch (fCubaIntegrationMethod) {
case BCIntegrate::kCubaSuave:
parameters_int.push_back(fCubaSuaveNNew);
parameters_double.push_back(fCubaSuaveFlatness);
break;
case BCIntegrate::kCubaDivonne:
break;
case BCIntegrate::kCubaCuhre:
break;
default: // if unknown method run Vegas. Shouldn't ever happen anyway
case BCIntegrate::kCubaVegas:
parameters_int.push_back(fCubaVegasNStart);
parameters_int.push_back(fCubaVegasNIncrease);
}
// print to log
BCLog::OutDebug( Form(" --> Running Cuba/%s integation over %i dimensions.", DumpCubaIntegrationMethod().c_str(), fNvar));
BCLog::OutDebug( Form(" --> Maximum number of iterations: %i", fCubaMaxEval));
BCLog::OutDebug( Form(" --> Aimed relative precision: %e", fRelativePrecision));
return CubaIntegrate(fCubaIntegrationMethod, parameters_double, parameters_int);
#else
BCLog::OutError("!!! This version of BAT is compiled without Cuba.");
BCLog::OutError(" Use other integration methods or install Cuba and recompile BAT.");
return 0.;
#endif
}
| void BCIntegrate::DeleteVarList | ( | ) |
| std::string BCIntegrate::DumpCubaIntegrationMethod | ( | BCIntegrate::BCCubaMethod | type | ) |
Definition at line 2372 of file BCIntegrate.cxx.
{
switch(type) {
case BCIntegrate::kCubaVegas:
return "Vegas";
case BCIntegrate::kCubaSuave:
return "Suave";
case BCIntegrate::kCubaDivonne:
return "Divonne";
case BCIntegrate::kCubaCuhre:
return "Cuhre";
default:
return "Undefined";
}
}
| std::string BCIntegrate::DumpCubaIntegrationMethod | ( | ) | [inline] |
Definition at line 849 of file BCIntegrate.h.
{ return DumpCubaIntegrationMethod(fCubaIntegrationMethod); }
| std::string BCIntegrate::DumpIntegrationMethod | ( | BCIntegrate::BCIntegrationMethod | type | ) |
Definition at line 2327 of file BCIntegrate.cxx.
{
switch(type) {
case BCIntegrate::kIntMonteCarlo:
return "Sampled Mean Monte Carlo";
case BCIntegrate::kIntImportance:
return "Importance Sampling";
case BCIntegrate::kIntMetropolis:
return "Metropolis";
case BCIntegrate::kIntCuba:
return "Cuba";
default:
return "Undefined";
}
}
| std::string BCIntegrate::DumpIntegrationMethod | ( | ) | [inline] |
Definition at line 835 of file BCIntegrate.h.
{ return DumpIntegrationMethod(fIntegrationMethod); }
| std::string BCIntegrate::DumpMarginalizationMethod | ( | ) | [inline] |
Definition at line 839 of file BCIntegrate.h.
{ return DumpMarginalizationMethod(fMarginalizationMethod); }
| std::string BCIntegrate::DumpMarginalizationMethod | ( | BCIntegrate::BCMarginalizationMethod | type | ) |
Definition at line 2344 of file BCIntegrate.cxx.
{
switch(type) {
case BCIntegrate::kMargMonteCarlo:
return "Monte Carlo Integration";
case BCIntegrate::kMargMetropolis:
return "Metropolis MCMC";
default:
return "Undefined";
}
}
| std::string BCIntegrate::DumpOptimizationMethod | ( | BCIntegrate::BCOptimizationMethod | type | ) |
Definition at line 2357 of file BCIntegrate.cxx.
{
switch(type) {
case BCIntegrate::kOptSA:
return "Simulated Annealing";
case BCIntegrate::kOptMetropolis:
return "Metropolis MCMC";
case BCIntegrate::kOptMinuit:
return "Minuit";
default:
return "Undefined";
}
}
| std::string BCIntegrate::DumpOptimizationMethod | ( | ) | [inline] |
Definition at line 843 of file BCIntegrate.h.
{ return DumpOptimizationMethod(fOptimizationMethod); }
| std::string BCIntegrate::DumpUsedOptimizationMethod | ( | ) | [inline] |
Definition at line 845 of file BCIntegrate.h.
{ return DumpOptimizationMethod(fOptimizationMethodMode); }
| double BCIntegrate::Eval | ( | std::vector< double > | x | ) | [virtual] |
Evaluate the unnormalized probability at a point in parameter space. Method needs to be overloaded by the user.
- Parameters:
-
x The point in parameter space
- Returns:
- The unnormalized probability
Reimplemented in BCModel.
Definition at line 505 of file BCIntegrate.cxx.
{
BCLog::OutWarning( "BCIntegrate::Eval. No function. Function needs to be overloaded.");
return 0;
}
| double BCIntegrate::EvalPrint | ( | std::vector< double > | x | ) |
Evaluate the un-normalized probability at a point in parameter space and prints the result to the log.
- Parameters:
-
x The point in parameter space
- Returns:
- The un-normalized probability
- See also:
- Eval(std::vector <double> x)
Definition at line 532 of file BCIntegrate.cxx.
{
double val=Eval(x);
BCLog::OutDebug(Form("BCIntegrate::EvalPrint. Value: %g.", val));
return val;
}
| double BCIntegrate::EvalSampling | ( | std::vector< double > | x | ) | [virtual] |
Evaluate the sampling function at a point in parameter space. Method needs to be overloaded by the user.
- Parameters:
-
x The point in parameter space
- Returns:
- The value of the sampling function
Reimplemented in BCModel.
Definition at line 519 of file BCIntegrate.cxx.
{
BCLog::OutWarning( "BCIntegrate::EvalSampling. No function. Function needs to be overloaded.");
return 0;
}
| void BCIntegrate::FCNLikelihood | ( | int & | npar, | |
| double * | grad, | |||
| double & | fval, | |||
| double * | par, | |||
| int | flag | |||
| ) | [static] |
Definition at line 1942 of file BCIntegrate.cxx.
{
// copy parameters
std::vector <double> parameters;
int n = global_this->GetNvar();
for (int i = 0; i < n; i++)
parameters.push_back(par[i]);
fval = - global_this->LogEval(parameters);
// delete parameters
parameters.clear();
}
| void BCIntegrate::FindModeMCMC | ( | ) |
Does the mode finding using Markov Chain Monte Carlo (prerun only!)
Definition at line 1961 of file BCIntegrate.cxx.
{
// call PreRun
// if (flag_run == 0)
if (!fMCMCFlagPreRun)
MCMCMetropolisPreRun();
// find global maximum
// double probmax = (MCMCGetMaximumLogProb()).at(0);
double probmax = 0;
if ( int(fBestFitParameters.size()) == fNvar) {
probmax = Eval( fBestFitParameters );
// fBestFitParameters = MCMCGetMaximumPoint(0);
}
// loop over all chains and find the maximum point
for (int i = 0; i < fMCMCNChains; ++i) {
double prob = exp( (MCMCGetMaximumLogProb()).at(i));
// copy the point into the vector
if ( (prob >= probmax && !fFlagIgnorePrevOptimization) || fFlagIgnorePrevOptimization) {
probmax = prob;
fBestFitParameters.clear();
fBestFitParameterErrors.clear();
fBestFitParameters = MCMCGetMaximumPoint(i);
for (int j = 0; j < fNvar; ++j)
fBestFitParameterErrors.push_back(-1.); // error undefined
fOptimizationMethodMode = BCIntegrate::kOptMetropolis;
}
}
}
| void BCIntegrate::FindModeMinuit | ( | std::vector< double > | start = std::vector<double>(0), |
|
| int | printlevel = 1 | |||
| ) | [virtual] |
Does the mode finding Does the mode finding using Minuit. If starting point is not specified, finding will start from the center of the parameter space.
- Parameters:
-
start point in parameter space from which the mode finding is started. printlevel The print level.
Reimplemented in BCModel.
Definition at line 1457 of file BCIntegrate.cxx.
{
bool have_start = true;
// check start values
if (int(start.size()) != fNvar)
have_start = false;
// set global this
global_this = this;
// define minuit
if (fMinuit)
delete fMinuit;
fMinuit = new TMinuit(fNvar);
// set print level
fMinuit->SetPrintLevel(printlevel);
// set function
fMinuit->SetFCN(&BCIntegrate::FCNLikelihood);
// set UP for likelihood
fMinuit->SetErrorDef(0.5);
// set parameters
int flag;
for (int i = 0; i < fNvar; i++) {
double starting_point = (fMin[i]+fMax[i])/2.;
if(have_start)
starting_point = start[i];
fMinuit->mnparm(i,
fx->at(i)->GetName().data(),
starting_point,
(fMax[i]-fMin[i])/100.,
fMin[i],
fMax[i],
flag);
}
// do mcmc minimization
// fMinuit->mnseek();
// do minimization
fMinuit->mnexcm("MIGRAD", fMinuitArglist, 2, flag);
// improve search for local minimum
// fMinuit->mnimpr();
// copy flag
fMinuitErrorFlag = flag;
// check if mode found by minuit is better than previous estimation
double probmax = 0;
bool valid = false;
if ( int(fBestFitParameters.size()) == fNvar) {
valid = true;
for (int i = 0; i < fNvar; ++i)
if (fBestFitParameters.at(i) < fMin[i] || fBestFitParameters.at(i) > fMax[i])
valid= false;
if (valid)
probmax = Eval( fBestFitParameters );
}
std::vector<double> tempvec;
for (int i = 0; i < fNvar; i++) {
double par;
double parerr;
fMinuit->GetParameter(i, par, parerr);
tempvec.push_back(par);
}
double probmaxminuit = Eval( tempvec );
// set best fit parameters
if ( (probmaxminuit > probmax && !fFlagIgnorePrevOptimization) || !valid || fFlagIgnorePrevOptimization) {
fBestFitParameters.clear();
fBestFitParameterErrors.clear();
for (int i = 0; i < fNvar; i++) {
double par;
double parerr;
fMinuit->GetParameter(i, par, parerr);
fBestFitParameters.push_back(par);
fBestFitParameterErrors.push_back(parerr);
}
// set optimization method used to find the mode
fOptimizationMethodMode = BCIntegrate::kOptMinuit;
}
// delete minuit
delete fMinuit;
fMinuit = 0;
return;
}
| void BCIntegrate::FindModeSA | ( | std::vector< double > | start = std::vector<double>(0) |
) |
Does the mode finding using Simulated Annealing. If starting point is not specified, finding will start from the center of the parameter space.
- Parameters:
-
start point in parameter space from thich the mode finding is started.
Definition at line 1573 of file BCIntegrate.cxx.
{
// note: if f(x) is the function to be minimized, then
// f(x) := - LogEval(parameters)
bool have_start = true;
std::vector<double> x, y, best_fit; // vectors for current state, new proposed state and best fit up to now
double fval_x, fval_y, fval_best_fit; // function values at points x, y and best_fit (we save them rather than to re-calculate them every time)
int t = 1; // time iterator
// check start values
if (int(start.size()) != fNvar)
have_start = false;
// if no starting point is given, set to center of parameter space
if ( !have_start ) {
start.clear();
for (int i = 0; i < fNvar; i++)
start.push_back((fMin[i]+fMax[i])/2.);
}
// set current state and best fit to starting point
x.clear();
best_fit.clear();
for (int i = 0; i < fNvar; i++) {
x.push_back(start[i]);
best_fit.push_back(start[i]);
}
// calculate function value at starting point
fval_x = fval_best_fit = LogEval(x);
// run while still "hot enough"
while ( SATemperature(t) > fSATmin ) {
// generate new state
y = GetProposalPointSA(x, t);
// check if the proposed point is inside the phase space
// if not, reject it
bool is_in_ranges = true;
for (int i = 0; i < fNvar; i++)
if (y[i] > fMax[i] || y[i] < fMin[i])
is_in_ranges = false;
if ( !is_in_ranges )
; // do nothing...
else {
// calculate function value at new point
fval_y = LogEval(y);
// is it better than the last one?
// if so, update state and chef if it is the new best fit...
if (fval_y >= fval_x) {
x.clear();
for (int i = 0; i < fNvar; i++)
x.push_back(y[i]);
fval_x = fval_y;
if (fval_y > fval_best_fit) {
best_fit.clear();
for (int i = 0; i < fNvar; i++)
best_fit.push_back(y[i]);
fval_best_fit = fval_y;
}
}
// ...else, only accept new state w/ certain probability
else {
if (fRandom->Rndm() <= exp( (fval_y - fval_x) / SATemperature(t) )) {
x.clear();
for (int i = 0; i < fNvar; i++)
x.push_back(y[i]);
fval_x = fval_y;
}
}
}
// update tree variables
fSANIterations = t;
fSATemperature = SATemperature(t);
fSALogProb = fval_x;
fSAx.clear();
for (int i = 0; i < fNvar; ++i)
fSAx.push_back(x[i]);
// fill tree
if (fFlagWriteSAToFile)
fTreeSA->Fill();
// increate t
t++;
}
// check if mode found by minuit is better than previous estimation
double probmax = 0;
bool valid = false;
if ( int(fBestFitParameters.size()) == fNvar) {
valid = true;
for (int i = 0; i < fNvar; ++i)
if (fBestFitParameters.at(i) < fMin[i] || fBestFitParameters.at(i) > fMax[i])
valid= false;
if (valid)
probmax = Eval( fBestFitParameters );
}
double probmaxsa = Eval( best_fit);
if ( (probmaxsa > probmax && !fFlagIgnorePrevOptimization) || !valid || fFlagIgnorePrevOptimization) {
// set best fit parameters
fBestFitParameters.clear();
fBestFitParameterErrors.clear();
for (int i = 0; i < fNvar; i++) {
fBestFitParameters.push_back(best_fit[i]);
fBestFitParameterErrors.push_back(-1.); // error undefined
}
// set optimization moethod used to find the mode
fOptimizationMethodMode = BCIntegrate::kOptSA;
}
return;
}
| virtual double BCIntegrate::FitFunction | ( | std::vector< double > | , | |
| std::vector< double > | ||||
| ) | [inline, virtual] |
Defines a fit function.
- Parameters:
-
parameters A set of parameter values x A vector of x-values
- Returns:
- The value of the fit function at the x-values given a set of parameters
Reimplemented in BCEfficiencyFitter, BCGraphFitter, and BCHistogramFitter.
Definition at line 597 of file BCIntegrate.h.
{ return 0.; }
| double BCIntegrate::GetAbsolutePrecision | ( | ) | [inline] |
- Returns:
- Requested absolute precision of the numerical integation
Definition at line 196 of file BCIntegrate.h.
{ return fAbsolutePrecision; }
| BCCubaMethod BCIntegrate::GetCubaIntegrationMethod | ( | ) | [inline] |
- Returns:
- Cuba Integration method
Definition at line 201 of file BCIntegrate.h.
{ return fCubaIntegrationMethod; }
| int BCIntegrate::GetCubaMaxEval | ( | ) | [inline] |
- Returns:
- Maximum number of evaluations in Cuba integration
Definition at line 211 of file BCIntegrate.h.
{ return fCubaMaxEval; }
| int BCIntegrate::GetCubaMinEval | ( | ) | [inline] |
- Returns:
- Minimum number of evaluations in Cuba integration
Definition at line 206 of file BCIntegrate.h.
{ return fCubaMinEval; }
| double BCIntegrate::GetCubaSuaveFlatness | ( | ) | [inline] |
- Returns:
- Flatness for Cuba Suave
Definition at line 236 of file BCIntegrate.h.
{ return fCubaSuaveFlatness; }
| int BCIntegrate::GetCubaSuaveNNew | ( | ) | [inline] |
- Returns:
- Number of new integrand evaluations in each subdivision for Cuba Suave
Definition at line 231 of file BCIntegrate.h.
{ return fCubaSuaveNNew; }
| int BCIntegrate::GetCubaVegasNIncrease | ( | ) | [inline] |
- Returns:
- Increase in number of evaluations per iteration for Cuba Vegas
Definition at line 226 of file BCIntegrate.h.
{ return fCubaVegasNIncrease; }
| int BCIntegrate::GetCubaVegasNStart | ( | ) | [inline] |
- Returns:
- Initial number of evaluations per iteration for Cuba Vegas
Definition at line 221 of file BCIntegrate.h.
{ return fCubaVegasNStart; }
| int BCIntegrate::GetCubaVerbositylevel | ( | ) | [inline] |
- Returns:
- Verbosity level of Cuba integration
Definition at line 216 of file BCIntegrate.h.
{ return fCubaVerbosity; }
| double BCIntegrate::GetError | ( | ) | [inline] |
- Returns:
- The uncertainty in the most recent Monte Carlo integration
Definition at line 241 of file BCIntegrate.h.
{ return fError; }
| TH1D * BCIntegrate::GetH1D | ( | int | parIndex | ) |
- Parameters:
-
parIndex1 Index of parameter
- Returns:
- Pointer to 1D histogram (TH1D) of marginalized distribution wrt. parameter with given index.
Definition at line 1191 of file BCIntegrate.cxx.
{
if(fHProb1D.size()==0) {
BCLog::OutWarning("BCModel::GetH1D. MarginalizeAll() has to be run prior to this to fill the distributions.");
return 0;
}
if(parIndex<0 || parIndex>fNvar) {
BCLog::OutWarning(Form("BCIntegrate::GetH1D. Parameter index %d is invalid.",parIndex));
return 0;
}
return fHProb1D[parIndex];
}
| TH2D * BCIntegrate::GetH2D | ( | int | parIndex1, | |
| int | parIndex2 | |||
| ) |
- Parameters:
-
parIndex1 Index of first parameter parIndex2 Index of second parameter, with parIndex2>parIndex1
- Returns:
- Pointer to 2D histogram (TH2D) of marginalized distribution wrt. parameters with given indeces.
Definition at line 1243 of file BCIntegrate.cxx.
{
if(fHProb2D.size()==0) {
BCLog::OutWarning("BCModel::GetH2D. MarginalizeAll() has to be run prior to this to fill the distributions.");
return 0;
}
int hindex = GetH2DIndex(parIndex1, parIndex2);
if(hindex==-1)
return 0;
if(hindex>(int)fHProb2D.size()-1) {
BCLog::OutWarning("BCIntegrate::GetH2D. Got invalid index from GetH2DIndex(). Something went wrong.");
return 0;
}
return fHProb2D[hindex];
}
| int BCIntegrate::GetH2DIndex | ( | int | parIndex1, | |
| int | parIndex2 | |||
| ) |
- Parameters:
-
parIndex1 Index of first parameter parIndex2 Index of second parameter, with parIndex2>parIndex1
- Returns:
- Index of the distribution in the vector of 2D distributions, which corresponds to the combination of parameters with given indeces
Definition at line 1207 of file BCIntegrate.cxx.
{
if(parIndex1>fNvar-1 || parIndex1<0) {
BCLog::OutWarning(Form("BCIntegrate::GetH2DIndex. Parameter index %d is invalid", parIndex1));
return -1;
}
if(parIndex2>fNvar-1 || parIndex2<0) {
BCLog::OutWarning(Form("BCIntegrate::GetH2DIndex. Parameter index %d is invalid", parIndex2));
return -1;
}
if(parIndex1==parIndex2) {
BCLog::OutWarning(Form("BCIntegrate::GetH2DIndex. Parameters have equal indeces: %d , %d", parIndex1, parIndex2));
return -1;
}
if(parIndex1>parIndex2) {
BCLog::OutWarning("BCIntegrate::GetH2DIndex. First parameters must be smaller than second (sorry).");
return -1;
}
int index=0;
for(int i=0;i<fNvar-1;i++)
for(int j=i+1;j<fNvar;j++) {
if(i==parIndex1 && j==parIndex2)
return index;
index++;
}
BCLog::OutWarning("BCIntegrate::GetH2DIndex. Invalid index combination.");
return -1;
}
| BCIntegrate::BCIntegrationMethod BCIntegrate::GetIntegrationMethod | ( | ) | [inline] |
- Returns:
- The integration method
Definition at line 108 of file BCIntegrate.h.
{ return fIntegrationMethod; }
| BCIntegrate::BCMarginalizationMethod BCIntegrate::GetMarginalizationMethod | ( | ) | [inline] |
- Returns:
- The marginalization method
Definition at line 113 of file BCIntegrate.h.
{ return fMarginalizationMethod; }
| std::vector<double>* BCIntegrate::GetMarkovChainPoint | ( | ) | [inline] |
Returns the actual point in the markov chain
Definition at line 263 of file BCIntegrate.h.
{ return &fXmetro1; }
| TTree* BCIntegrate::GetMarkovChainTree | ( | ) | [inline] |
- Returns:
- The ROOT tree containing the Markov chain
Definition at line 258 of file BCIntegrate.h.
{ return fMarkovChainTree; }
| double* BCIntegrate::GetMarkovChainValue | ( | ) | [inline] |
Returns the value of the loglikelihood at the point fXmetro1
Definition at line 273 of file BCIntegrate.h.
{ return &fMarkovChainValue; }
| int* BCIntegrate::GetMCMCIteration | ( | ) | [inline] |
Returns the iteration of the MCMC
Definition at line 268 of file BCIntegrate.h.
{ return &fMCMCIteration; }
| TMinuit * BCIntegrate::GetMinuit | ( | ) |
- Returns:
- Minuit used for mode finding
Definition at line 1447 of file BCIntegrate.cxx.
| int BCIntegrate::GetMinuitErrorFlag | ( | ) | [inline] |
Definition at line 253 of file BCIntegrate.h.
{ return fMinuitErrorFlag; }
| int BCIntegrate::GetNbins | ( | ) | [inline] |
- Returns:
- number of bins per dimension for the marginalized distributions
Definition at line 246 of file BCIntegrate.h.
{ return fNbins; }
| int BCIntegrate::GetNIterations | ( | ) | [inline] |
- Returns:
- The number of iterations for the most recent Monte Carlo integration
Definition at line 186 of file BCIntegrate.h.
{ return fNIterations; }
| int BCIntegrate::GetNIterationsMax | ( | ) | [inline] |
- Returns:
- The number of maximum iterations for Monte Carlo integration
Definition at line 181 of file BCIntegrate.h.
{ return fNIterationsMax; }
| int BCIntegrate::GetNiterationsPerDimension | ( | ) | [inline] |
- Returns:
- The number of iterations per dimension for the Monte Carlo integration
Definition at line 166 of file BCIntegrate.h.
{ return fNiterPerDimension; }
| int BCIntegrate::GetNSamplesPer2DBin | ( | ) | [inline] |
- Returns:
- Number of samples per 2D bin per variable in the Metropolis marginalization.
Definition at line 171 of file BCIntegrate.h.
{ return fNSamplesPer2DBin; }
| int BCIntegrate::GetNvar | ( | ) | [inline] |
- Returns:
- The number of variables to integrate over
Definition at line 176 of file BCIntegrate.h.
{ return fNvar; }
| BCIntegrate::BCOptimizationMethod BCIntegrate::GetOptimizationMethod | ( | ) | [inline] |
- Returns:
- The current optimization method
Definition at line 118 of file BCIntegrate.h.
{ return fOptimizationMethod; }
| BCIntegrate::BCOptimizationMethod BCIntegrate::GetOptimizationMethodMode | ( | ) | [inline] |
- Returns:
- The optimization method used to find the mode
Definition at line 123 of file BCIntegrate.h.
{ return fOptimizationMethodMode; }
| std::vector< double > BCIntegrate::GetProposalPointSA | ( | std::vector< double > | x, | |
| int | t | |||
| ) |
Generates a new state in a neighbourhood around x that is to be accepted or rejected by the Simulated Annealing algorithm. Delegates the generation to the appropriate method according to fSASchedule.
- Parameters:
-
x last state. t time iterator to determine current temperature.
Definition at line 1738 of file BCIntegrate.cxx.
{
// do we have Cauchy (default), Boltzmann or custom annealing schedule?
if (fSASchedule == BCIntegrate::kSABoltzmann)
return GetProposalPointSABoltzmann(x, t);
else if (fSASchedule == BCIntegrate::kSACauchy)
return GetProposalPointSACauchy(x, t);
else
return GetProposalPointSACustom(x, t);
}
| std::vector< double > BCIntegrate::GetProposalPointSABoltzmann | ( | std::vector< double > | x, | |
| int | t | |||
| ) |
Generates a new state in a neighbourhood around x that is to be accepted or rejected by the Simulated Annealing algorithm. This method is used for Boltzmann annealing schedule.
- Parameters:
-
x last state. t time iterator to determine current temperature.
Definition at line 1751 of file BCIntegrate.cxx.
{
std::vector<double> y;
y.clear();
double new_val, norm;
for (int i = 0; i < fNvar; i++) {
norm = (fMax[i] - fMin[i]) * SATemperature(t) / 2.;
new_val = x[i] + norm * fRandom->Gaus();
y.push_back(new_val);
}
return y;
}
| std::vector< double > BCIntegrate::GetProposalPointSACauchy | ( | std::vector< double > | x, | |
| int | t | |||
| ) |
Generates a new state in a neighbourhood around x that is to be accepted or rejected by the Simulated Annealing algorithm. This method is used for Cauchy annealing schedule.
- Parameters:
-
x last state. t time iterator to determine current temperature.
Definition at line 1768 of file BCIntegrate.cxx.
{
std::vector<double> y;
y.clear();
if (fNvar == 1) {
double cauchy, new_val, norm;
norm = (fMax[0] - fMin[0]) * SATemperature(t) / 2.;
cauchy = tan(3.14159 * (fRandom->Rndm() - 0.5));
new_val = x[0] + norm * cauchy;
y.push_back(new_val);
}
else {
// use sampling to get radial n-dim Cauchy distribution
// first generate a random point uniformly distributed on a
// fNvar-dimensional hypersphere
y = SAHelperGetRandomPointOnHypersphere();
// scale the vector by a random factor determined by the radial
// part of the fNvar-dimensional Cauchy distribution
double radial = SATemperature(t) * SAHelperGetRadialCauchy();
// scale y by radial part and the size of dimension i in phase space
// afterwards, move by x
for (int i = 0; i < fNvar; i++)
y[i] = (fMax[i] - fMin[i]) * y[i] * radial / 2. + x[i];
}
return y;
}
| std::vector< double > BCIntegrate::GetProposalPointSACustom | ( | std::vector< double > | x, | |
| int | t | |||
| ) | [virtual] |
Generates a new state in a neighbourhood around x that is to be accepted or rejected by the Simulated Annealing algorithm. This is a virtual method to be overridden by a user-defined custom proposal function.
- Parameters:
-
x last state. t time iterator to determine current temperature.
Definition at line 1803 of file BCIntegrate.cxx.
{
BCLog::OutError("BCIntegrate::GetProposalPointSACustom : No custom proposal function defined");
return std::vector<double>(fNvar);
}
| double BCIntegrate::GetRandomPoint | ( | std::vector< double > & | x | ) |
Fills a vector of (flat) random numbers in the limits of the parameters and returns the probability at that point
- Parameters:
-
x A vector of doubles
- Returns:
- The (unnormalized) probability at the random point
Definition at line 1263 of file BCIntegrate.cxx.
| double BCIntegrate::GetRandomPointImportance | ( | std::vector< double > & | x | ) |
Fills a vector of random numbers in the limits of the parameters sampled by the sampling function and returns the probability at that point
- Parameters:
-
x A vector of doubles
- Returns:
- The (unnormalized) probability at the random point
Definition at line 1274 of file BCIntegrate.cxx.
{
double p = 1.1;
double g = 1.0;
while (p>g) {
GetRandomVector(x);
for(int i=0;i<fNvar;i++)
x[i] = fMin[i] + x[i] * (fMax[i]-fMin[i]);
p = fRandom->Rndm();
g = EvalSampling(x);
}
return Eval(x);
}
| void BCIntegrate::GetRandomPointMetro | ( | std::vector< double > & | x | ) |
Fills a vector of random numbers in the limits of the parameters sampled by the probality function and returns the probability at that point (Metropolis)
- Parameters:
-
x A vector of doubles
Definition at line 1317 of file BCIntegrate.cxx.
{
// get new point
GetRandomVectorMetro(fXmetro1);
// scale the point to the allowed region and stepsize
int in=1;
for(int i=0;i<fNvar;i++) {
fXmetro1[i] = fXmetro0[i] + fXmetro1[i] * (fMax[i]-fMin[i]);
// check whether the generated point is inside the allowed region
if( fXmetro1[i]<fMin[i] || fXmetro1[i]>fMax[i] )
in=0;
}
// calculate the log probabilities and compare old and new point
double p0 = fMarkovChainValue; // old point
double p1 = 0; // new point
int accept=0;
if(in) {
p1 = LogEval(fXmetro1);
if(p1>=p0)
accept=1;
else {
double r=log(fRandom->Rndm());
if(r<p1-p0)
accept=1;
}
}
// fill the return point after the decision
if(accept) {
// increase counter
fNacceptedMCMC++;
for(int i=0;i<fNvar;i++) {
fXmetro0[i]=fXmetro1[i];
x[i]=fXmetro1[i];
fMarkovChainValue = p1;
}
}
else
for(int i=0;i<fNvar;i++) {
x[i]=fXmetro0[i];
fXmetro1[i] = x[i];
fMarkovChainValue = p0;
}
fNmetro++;
}
| void BCIntegrate::GetRandomPointSamplingMetro | ( | std::vector< double > & | x | ) |
Fills a vector of random numbers in the limits of the parameters sampled by the sampling function and returns the probability at that point (Metropolis)
- Parameters:
-
x A vector of doubles
Definition at line 1372 of file BCIntegrate.cxx.
{
// get new point
GetRandomVectorMetro(fXmetro1);
// scale the point to the allowed region and stepsize
int in=1;
for(int i=0;i<fNvar;i++) {
fXmetro1[i] = fXmetro0[i] + fXmetro1[i] * (fMax[i]-fMin[i]);
// check whether the generated point is inside the allowed region
if( fXmetro1[i]<fMin[i] || fXmetro1[i]>fMax[i] )
in=0;
}
// calculate the log probabilities and compare old and new point
double p0 = LogEvalSampling(fXmetro0); // old point
double p1=0; // new point
if(in)
p1= LogEvalSampling(fXmetro1);
// compare log probabilities
int accept=0;
if(in) {
if(p1>=p0)
accept=1;
else {
double r=log(fRandom->Rndm());
if(r<p1-p0)
accept=1;
}
}
// fill the return point after the decision
if(accept)
for(int i=0;i<fNvar;i++) {
fXmetro0[i]=fXmetro1[i];
x[i]=fXmetro1[i];
}
else
for(int i=0;i<fNvar;i++)
x[i]=fXmetro0[i];
fNmetro++;
}
| void BCIntegrate::GetRandomVector | ( | std::vector< double > & | x | ) |
| void BCIntegrate::GetRandomVectorMetro | ( | std::vector< double > & | x | ) | [virtual] |
Definition at line 1433 of file BCIntegrate.cxx.
{
double * randx = new double[fNvar];
fRandom->RndmArray(fNvar, randx);
for(int i=0;i<fNvar;i++)
x[i] = 2.0 * (randx[i] - 0.5) * fMarkovChainStepSize;
delete[] randx;
randx = 0;
}
| double BCIntegrate::GetRelativePrecision | ( | ) | [inline] |
- Returns:
- Requested relative precision of the numerical integation
Definition at line 191 of file BCIntegrate.h.
{ return fRelativePrecision; }
| BCIntegrate::BCSASchedule BCIntegrate::GetSASchedule | ( | ) | [inline] |
- Returns:
- The Simulated Annealing schedule
Definition at line 128 of file BCIntegrate.h.
{ return fSASchedule; }
| double BCIntegrate::GetSAT0 | ( | ) | [inline] |
Returns the Simulated Annealing starting temperature.
Definition at line 278 of file BCIntegrate.h.
{ return fSAT0; }
| double BCIntegrate::GetSATmin | ( | ) | [inline] |
Returns the Simulated Annealing threshhold temperature.
Definition at line 283 of file BCIntegrate.h.
{ return fSATmin; }
| TTree* BCIntegrate::GetSATree | ( | ) | [inline] |
Getter for the tree containing the Simulated Annealing chain.
Definition at line 529 of file BCIntegrate.h.
{ return fTreeSA; }
| void BCIntegrate::InitializeSATree | ( | ) |
Initialization of the tree for the Simulated Annealing
Definition at line 1557 of file BCIntegrate.cxx.
{
if (fTreeSA)
delete fTreeSA;
fTreeSA = new TTree("SATree", "SATree");
fTreeSA->Branch("Iteration", &fSANIterations, "iteration/I");
fTreeSA->Branch("NParameters", &fNvar, "parameters/I");
fTreeSA->Branch("Temperature", &fSATemperature, "temperature/D");
fTreeSA->Branch("LogProbability", &fSALogProb, "log(probability)/D");
for (int i = 0; i < fNvar; ++i)
fTreeSA->Branch(TString::Format("Parameter%i", i), &fSAx[i], TString::Format("parameter %i/D", i));
}
| void BCIntegrate::InitMetro | ( | ) |
Initializes the Metropolis algorithm (for details see manual)
Definition at line 1294 of file BCIntegrate.cxx.
{
fNmetro=0;
if (fXmetro0.size() <= 0) {
// start in the center of the phase space
for(int i=0;i<fNvar;i++)
fXmetro0.push_back((fMin[i]+fMax[i])/2.0);
}
fMarkovChainValue = LogEval(fXmetro0);
// run metropolis for a few times and dump the result... (to forget the initial position)
std::vector <double> x;
x.assign(fNvar,0.);
for(int i=0;i<1000;i++)
GetRandomPointMetro(x);
fNmetro = 0;
}
| double BCIntegrate::IntegralImportance | ( | std::vector< double > | x | ) |
Perfoms the importance sampling Monte Carlo integration. For details see documentation.
- Parameters:
-
x An initial point in parameter space
- Returns:
- The integral
Definition at line 767 of file BCIntegrate.cxx.
{
// print debug information
BCLog::OutDebug(Form("BCIntegrate::IntegralImportance. Integate over %i dimensions.", fNvar));
// get total number of iterations
double Niter = pow(fNiterPerDimension, fNvar);
// print if more than 100,000 iterations
if(Niter>1e5)
BCLog::OutDebug(Form("BCIntegrate::IntegralImportance. Total number of iterations: %d.", (int)Niter));
// reset sum
double sumI = 0;
std::vector <double> randx;
randx.assign(fNvar,0.);
// prepare maximum value
double pmax = 0.0;
// loop over iterations
for(int i = 0; i <= Niter; i++) {
// get random point from sampling distribution
GetRandomPointImportance(randx);
// calculate probability at random point
double val_f = Eval(randx);
// calculate sampling distributions at that point
double val_g = EvalSampling(randx);
// add ratio to sum
if (val_g > 0)
sumI += val_f / val_g;
// search for maximum probability
if (val_f > pmax) {
// set new maximum value
pmax = val_f;
// delete old best fit parameter values
fBestFitParameters.clear();
// write best fit parameters
for(int i = 0; i < fNvar; i++)
fBestFitParameters.push_back(randx.at(i));
}
// write intermediate results
if((i+1)%100000 == 0)
BCLog::OutDebug(Form("BCIntegrate::IntegralImportance : Iteration %d, integral: %g.", i+1, sumI/double(i+1)));
}
// calculate integral
double result=sumI/Niter;
// print debug information
BCLog::OutDebug(Form("BCIntegrate::IntegralImportance : Integral %g in %i iterations.", result, (int)Niter));
return result;
}
| double BCIntegrate::IntegralMC | ( | std::vector< double > | x, | |
| int * | varlist | |||
| ) |
Perfoms the Monte Carlo integration. For details see documentation.
- Parameters:
-
x An initial point in parameter space varlist A list of variables
- Returns:
- The integral
Definition at line 601 of file BCIntegrate.cxx.
{
SetVarList(varlist);
return IntegralMC(x);
}
| double BCIntegrate::IntegralMC | ( | std::vector< double > | x | ) |
Definition at line 608 of file BCIntegrate.cxx.
{
// count the variables to integrate over
int NvarNow=0;
for(int i = 0; i < fNvar; i++)
if(fVarlist[i])
NvarNow++;
// print to log
BCLog::OutDebug(Form(" --> Running MC integation over %i dimensions.", NvarNow));
BCLog::OutDebug(Form(" --> Maximum number of iterations: %i", GetNIterationsMax()));
BCLog::OutDebug(Form(" --> Aimed relative precision: %e", GetRelativePrecision()));
// calculate D (the integrated volume)
double D = 1.0;
for(int i = 0; i < fNvar; i++)
if (fVarlist[i])
D = D * (fMax[i] - fMin[i]);
// reset variables
double pmax = 0.0;
double sumW = 0.0;
double sumW2 = 0.0;
double integral = 0.0;
double variance = 0.0;
double relprecision = 1.0;
std::vector <double> randx;
randx.assign(fNvar, 0.0);
// reset number of iterations
fNIterations = 0;
// iterate while precision is not reached and the number of iterations is lower than maximum number of iterations
while ((fRelativePrecision < relprecision && fNIterationsMax > fNIterations) || fNIterations < 10) {
// increase number of iterations
fNIterations++;
// get random numbers
GetRandomVector(randx);
// scale random numbers
for(int i = 0; i < fNvar; i++) {
if(fVarlist[i])
randx[i]=fMin[i]+randx[i]*(fMax[i]-fMin[i]);
else
randx[i]=x[i];
}
// evaluate function at sampled point
double value = Eval(randx);
// add value to sum and sum of squares
sumW += value;
sumW2 += value * value;
// search for maximum probability
if (value > pmax) {
// set new maximum value
pmax = value;
// delete old best fit parameter values
fBestFitParameters.clear();
// write best fit parameters
for(int i = 0; i < fNvar; i++)
fBestFitParameters.push_back(randx.at(i));
}
// calculate integral and variance
integral = D * sumW / fNIterations;
if (fNIterations%10000 == 0) {
variance = (1.0 / double(fNIterations)) * (D * D * sumW2 / double(fNIterations) - integral * integral);
double error = sqrt(variance);
relprecision = error / integral;
BCLog::OutDebug(Form("BCIntegrate::IntegralMC. Iteration %i, integral: %e +- %e.", fNIterations, integral, error));
}
}
fError = variance / fNIterations;
// print to log
BCLog::OutDebug(Form(" --> Result of integration: %e +- %e in %i iterations.", integral, sqrt(variance), fNIterations));
BCLog::OutDebug(Form(" --> Obtained relative precision: %e. ", sqrt(variance) / integral));
return integral;
}
| double BCIntegrate::IntegralMetro | ( | std::vector< double > | x | ) |
Perfoms the Metropolis Monte Carlo integration. For details see documentation.
- Parameters:
-
x An initial point in parameter space
- Returns:
- The integral
Definition at line 701 of file BCIntegrate.cxx.
{
// print debug information
BCLog::OutDebug(Form("BCIntegrate::IntegralMetro. Integate over %i dimensions.", fNvar));
// get total number of iterations
double Niter = pow(fNiterPerDimension, fNvar);
// print if more than 100,000 iterations
if(Niter>1e5)
BCLog::OutDebug(Form(" --> Total number of iterations in Metropolis: %d.", (int)Niter));
// reset sum
double sumI = 0;
// prepare Metropolis algorithm
std::vector <double> randx;
randx.assign(fNvar,0.);
InitMetro();
// prepare maximum value
double pmax = 0.0;
// loop over iterations
for(int i = 0; i <= Niter; i++) {
// get random point from sampling distribution
GetRandomPointSamplingMetro(randx);
// calculate probability at random point
double val_f = Eval(randx);
// calculate sampling distributions at that point
double val_g = EvalSampling(randx);
// add ratio to sum
if (val_g > 0)
sumI += val_f / val_g;
// search for maximum probability
if (val_f > pmax) {
// set new maximum value
pmax = val_f;
// delete old best fit parameter values
fBestFitParameters.clear();
// write best fit parameters
for(int i = 0; i < fNvar; i++)
fBestFitParameters.push_back(randx.at(i));
}
// write intermediate results
if((i+1)%100000 == 0)
BCLog::OutDebug(Form("BCIntegrate::IntegralMetro. Iteration %d, integral: %g.", i+1, sumI/double(i+1)));
}
// calculate integral
double result=sumI/Niter;
// print debug information
BCLog::OutDebug(Form(" --> Integral is %g in %g iterations.", result, Niter));
return result;
}
| double BCIntegrate::Integrate | ( | ) |
Does the integration over the un-normalized probability.
- Returns:
- The normalization value
Definition at line 564 of file BCIntegrate.cxx.
{
std::vector <double> parameter;
parameter.assign(fNvar, 0.0);
BCLog::OutSummary(
Form("Running numerical integration using %s (%s)",
DumpIntegrationMethod().c_str(),
DumpCubaIntegrationMethod().c_str()));
switch(fIntegrationMethod)
{
case BCIntegrate::kIntMonteCarlo:
return IntegralMC(parameter);
case BCIntegrate::kIntMetropolis:
return IntegralMetro(parameter);
case BCIntegrate::kIntImportance:
return IntegralImportance(parameter);
case BCIntegrate::kIntCuba:
#ifdef HAVE_CUBA_H
return CubaIntegrate();
#else
BCLog::OutError("!!! This version of BAT is compiled without Cuba.");
BCLog::OutError(" Use other integration methods or install Cuba and recompile BAT.");
#endif
default:
BCLog::OutError(
Form("BCIntegrate::Integrate : Invalid integration method: %d", fIntegrationMethod));
}
return 0;
}
| int BCIntegrate::IntegrateResetResults | ( | ) |
Reset all information on the best fit parameters.
- Returns:
- An error code
Definition at line 553 of file BCIntegrate.cxx.
{
fBestFitParameters.clear();
fBestFitParameterErrors.clear();
fBestFitParametersMarginalized.clear();
// no errors
return 1;
}
| double BCIntegrate::LogEval | ( | std::vector< double > | x | ) | [virtual] |
Evaluate the natural logarithm of the Eval function. For better numerical stability, this method should also be overloaded by the user.
- Parameters:
-
x The point in parameter space
- Returns:
- log(Eval(x))
Reimplemented from BCEngineMCMC.
Reimplemented in BCModel.
Definition at line 512 of file BCIntegrate.cxx.
{
// this method should better also be overloaded
return log(Eval(x));
}
| double BCIntegrate::LogEvalSampling | ( | std::vector< double > | x | ) |
Evaluate the natural logarithm of the EvalSampling function. Method needs to be overloaded by the user.
- Parameters:
-
x The point in parameter space
- Returns:
- log(EvalSampling(x))
Definition at line 526 of file BCIntegrate.cxx.
{
return log(EvalSampling(x));
}
| TH1D * BCIntegrate::Marginalize | ( | BCParameter * | parameter | ) |
Performs the marginalization with respect to one parameter.
- Parameters:
-
parameter The parameter w.r.t. which the marginalization is performed
- Returns:
- A histogram which contains the marginalized probability distribution (normalized to 1)
Definition at line 831 of file BCIntegrate.cxx.
{
BCLog::OutDebug(Form(" --> Marginalizing model wrt. parameter %s using %s.", parameter->GetName().data(), DumpMarginalizationMethod().c_str()));
switch(fMarginalizationMethod)
{
case BCIntegrate::kMargMonteCarlo:
return MarginalizeByIntegrate(parameter);
case BCIntegrate::kMargMetropolis:
return MarginalizeByMetro(parameter);
default:
BCLog::OutError(
Form("BCIntegrate::Marginalize. Invalid marginalization method: %d. Return 0.", fMarginalizationMethod));
}
return 0;
}
| TH2D * BCIntegrate::Marginalize | ( | BCParameter * | parameter1, | |
| BCParameter * | parameter2 | |||
| ) |
Performs the marginalization with respect to two parameters.
- Parameters:
-
parameter1 The first parameter w.r.t. which the marginalization is performed parameter2 The second parameter w.r.t. which the marginalization is performed
- Returns:
- A histogram which contains the marginalized probability distribution (normalized to 1)
Definition at line 853 of file BCIntegrate.cxx.
{
switch(fMarginalizationMethod)
{
case BCIntegrate::kMargMonteCarlo:
return MarginalizeByIntegrate(parameter1,parameter2);
case BCIntegrate::kMargMetropolis:
return MarginalizeByMetro(parameter1,parameter2);
default:
BCLog::OutError(
Form("BCIntegrate::Marginalize. Invalid marginalization method: %d. Return 0.", fMarginalizationMethod));
}
return 0;
}
| int BCIntegrate::MarginalizeAllByMetro | ( | const char * | name = "" |
) |
Performs the marginalization with respect to every single parameter as well as with respect all combinations to two parameters using the Metropolis algorithm.
- Parameters:
-
name Basename for the histograms (e.g. model name)
- Returns:
- Total number of marginalized distributions
Definition at line 1016 of file BCIntegrate.cxx.
{
int niter=fNbins*fNbins*fNSamplesPer2DBin*fNvar;
BCLog::OutDetail(Form(" --> Number of samples in Metropolis marginalization: %d.", niter));
// define 1D histograms
for(int i=0;i<fNvar;i++)
{
double hmin1 = fx->at(i)->GetLowerLimit();
double hmax1 = fx->at(i)->GetUpperLimit();
TH1D * h1 = new TH1D(
TString::Format("h%s_%s", name, fx->at(i)->GetName().data()),"",
fNbins, hmin1, hmax1);
fHProb1D.push_back(h1);
}
// define 2D histograms
for(int i=0;i<fNvar-1;i++)
for(int j=i+1;j<fNvar;j++) {
double hmin1 = fx->at(i)->GetLowerLimit();
double hmax1 = fx->at(i)->GetUpperLimit();
double hmin2 = fx->at(j)->GetLowerLimit();
double hmax2 = fx->at(j)->GetUpperLimit();
TH2D * h2 = new TH2D(
TString::Format("h%s_%s_%s", name, fx->at(i)->GetName().data(), fx->at(j)->GetName().data()),"",
fNbins, hmin1, hmax1,
fNbins, hmin2, hmax2);
fHProb2D.push_back(h2);
}
// get number of 2d distributions
int nh2d = fHProb2D.size();
BCLog::OutDetail(Form(" --> Marginalizing %d 1D distributions and %d 2D distributions.", fNvar, nh2d));
// prepare function fitting
double dx = 0.;
double dy = 0.;
if (fFitFunctionIndexX >= 0) {
dx = (fDataPointUpperBoundaries->GetValue(fFitFunctionIndexX) -
fDataPointLowerBoundaries->GetValue(fFitFunctionIndexX))
/ double(fErrorBandNbinsX);
dx = (fDataPointUpperBoundaries->GetValue(fFitFunctionIndexY) -
fDataPointLowerBoundaries->GetValue(fFitFunctionIndexY))
/ double(fErrorBandNbinsY);
fErrorBandXY = new TH2D(
TString::Format("errorbandxy_%d",BCLog::GetHIndex()), "",
fErrorBandNbinsX,
fDataPointLowerBoundaries->GetValue(fFitFunctionIndexX) - 0.5 * dx,
fDataPointUpperBoundaries->GetValue(fFitFunctionIndexX) + 0.5 * dx,
fErrorBandNbinsY,
fDataPointLowerBoundaries->GetValue(fFitFunctionIndexY) - 0.5 * dy,
fDataPointUpperBoundaries->GetValue(fFitFunctionIndexY) + 0.5 * dy);
fErrorBandXY->SetStats(kFALSE);
for (int ix = 1; ix <= fErrorBandNbinsX; ++ix)
for (int iy = 1; iy <= fErrorBandNbinsX; ++iy)
fErrorBandXY->SetBinContent(ix, iy, 0.0);
}
// prepare Metro
std::vector <double> randx;
randx.assign(fNvar, 0.0);
InitMetro();
// reset counter
fNacceptedMCMC = 0;
// run Metro and fill histograms
for(int i=0;i<=niter;i++) {
GetRandomPointMetro(randx);
// save this point to the markov chain in the ROOT file
if (fFlagWriteMarkovChain) {
fMCMCIteration = i;
fMarkovChainTree->Fill();
}
for(int j=0;j<fNvar;j++)
fHProb1D[j]->Fill( randx[j] );
int ih=0;
for(int j=0;j<fNvar-1;j++)
for(int k=j+1;k<fNvar;k++) {
fHProb2D[ih]->Fill(randx[j],randx[k]);
ih++;
}
if((i+1)%100000==0)
BCLog::OutDebug(Form("BCIntegrate::MarginalizeAllByMetro. %d samples done.",i+1));
// function fitting
if (fFitFunctionIndexX >= 0) {
// loop over all possible x values ...
if (fErrorBandContinuous) {
double x = 0;
for (int ix = 0; ix < fErrorBandNbinsX; ix++) {
// calculate x
x = fErrorBandXY->GetXaxis()->GetBinCenter(ix + 1);
// calculate y
std::vector <double> xvec;
xvec.push_back(x);
double y = FitFunction(xvec, randx);
xvec.clear();
// fill histogram
fErrorBandXY->Fill(x, y);
}
}
// ... or evaluate at the data point x-values
else {
int ndatapoints = int(fErrorBandX.size());
double x = 0;
for (int ix = 0; ix < ndatapoints; ++ix) {
// calculate x
x = fErrorBandX.at(ix);
// calculate y
std::vector <double> xvec;
xvec.push_back(x);
double y = FitFunction(xvec, randx);
xvec.clear();
// fill histogram
fErrorBandXY->Fill(x, y);
}
}
}
}
// normalize histograms
for(int i=0;i<fNvar;i++)
fHProb1D[i]->Scale( 1./fHProb1D[i]->Integral("width") );
for (int i=0;i<nh2d;i++)
fHProb2D[i]->Scale( 1./fHProb2D[i]->Integral("width") );
if (fFitFunctionIndexX >= 0)
fErrorBandXY->Scale(1.0/fErrorBandXY->Integral() * fErrorBandXY->GetNbinsX());
if (fFitFunctionIndexX >= 0) {
for (int ix = 1; ix <= fErrorBandNbinsX; ix++) {
double sum = 0;
for (int iy = 1; iy <= fErrorBandNbinsY; iy++)
sum += fErrorBandXY->GetBinContent(ix, iy);
for (int iy = 1; iy <= fErrorBandNbinsY; iy++) {
double newvalue = fErrorBandXY->GetBinContent(ix, iy) / sum;
fErrorBandXY->SetBinContent(ix, iy, newvalue);
}
}
}
BCLog::OutDebug(Form("BCIntegrate::MarginalizeAllByMetro done with %i trials and %i accepted trials. Efficiency is %f",fNmetro, fNacceptedMCMC, double(fNacceptedMCMC)/double(fNmetro)));
return fNvar+nh2d;
}
| TH1D * BCIntegrate::MarginalizeByIntegrate | ( | BCParameter * | parameter | ) |
Performs the marginalization with respect to one parameter using the simple Monte Carlo technique.
- Parameters:
-
parameter The parameter w.r.t. which the marginalization is performed
- Returns:
- A histogram which contains the marginalized probability distribution (normalized to 1)
Definition at line 872 of file BCIntegrate.cxx.
{
// set parameter to marginalize
ResetVarlist(1);
int index = parameter->GetIndex();
UnsetVar(index);
// define histogram
double hmin = parameter->GetLowerLimit();
double hmax = parameter->GetUpperLimit();
double hdx = (hmax - hmin) / double(fNbins);
TH1D * hist = new TH1D("hist","", fNbins, hmin, hmax);
// integrate
std::vector <double> randx;
randx.assign(fNvar, 0.);
for(int i=0;i<=fNbins;i++) {
randx[index] = hmin + (double)i * hdx;
double val = IntegralMC(randx);
hist->Fill(randx[index], val);
}
// normalize
hist->Scale( 1./hist->Integral("width") );
return hist;
}
| TH2D * BCIntegrate::MarginalizeByIntegrate | ( | BCParameter * | parameter1, | |
| BCParameter * | parameter2 | |||
| ) |
Performs the marginalization with respect to two parameters using the simple Monte Carlo technique.
- Parameters:
-
parameter1 The first parameter w.r.t. which the marginalization is performed parameter2 The second parameter w.r.t. which the marginalization is performed
- Returns:
- A histogram which contains the marginalized probability distribution (normalized to 1)
Definition at line 903 of file BCIntegrate.cxx.
{
// set parameter to marginalize
ResetVarlist(1);
int index1 = parameter1->GetIndex();
UnsetVar(index1);
int index2 = parameter2->GetIndex();
UnsetVar(index2);
// define histogram
double hmin1 = parameter1->GetLowerLimit();
double hmax1 = parameter1->GetUpperLimit();
double hdx1 = (hmax1 - hmin1) / double(fNbins);
double hmin2 = parameter2->GetLowerLimit();
double hmax2 = parameter2->GetUpperLimit();
double hdx2 = (hmax2 - hmin2) / double(fNbins);
TH2D * hist = new TH2D(Form("hist_%s_%s", parameter1->GetName().data(), parameter2->GetName().data()),"",
fNbins, hmin1, hmax1,
fNbins, hmin2, hmax2);
// integrate
std::vector <double> randx;
randx.assign(fNvar, 0.0);
for(int i=0;i<=fNbins;i++) {
randx[index1] = hmin1 + (double)i * hdx1;
for(int j=0;j<=fNbins;j++) {
randx[index2] = hmin2 + (double)j * hdx2;
double val = IntegralMC(randx);
hist->Fill(randx[index1],randx[index2], val);
}
}
// normalize
hist->Scale(1.0/hist->Integral("width"));
return hist;
}
| TH1D * BCIntegrate::MarginalizeByMetro | ( | BCParameter * | parameter | ) |
Performs the marginalization with respect to one parameter using the Metropolis algorithm.
- Parameters:
-
parameter The parameter w.r.t. which the marginalization is performed
- Returns:
- A histogram which contains the marginalized probability distribution (normalized to 1)
Definition at line 946 of file BCIntegrate.cxx.
{
int niter = fMarkovChainNIterations;
if (fMarkovChainAutoN == true)
niter = fNbins*fNbins*fNSamplesPer2DBin*fNvar;
BCLog::OutDetail(Form(" --> Number of samples in Metropolis marginalization: %d.", niter));
// set parameter to marginalize
int index = parameter->GetIndex();
// define histogram
double hmin = parameter->GetLowerLimit();
double hmax = parameter->GetUpperLimit();
TH1D * hist = new TH1D("hist","", fNbins, hmin, hmax);
// prepare Metro
std::vector <double> randx;
randx.assign(fNvar, 0.0);
InitMetro();
for(int i=0;i<=niter;i++) {
GetRandomPointMetro(randx);
hist->Fill(randx[index]);
}
// normalize
hist->Scale(1.0/hist->Integral("width"));
return hist;
}
| TH2D * BCIntegrate::MarginalizeByMetro | ( | BCParameter * | parameter1, | |
| BCParameter * | parameter2 | |||
| ) |
Performs the marginalization with respect to two parameters using the Metropolis algorithm.
- Parameters:
-
parameter1 The first parameter w.r.t. which the marginalization is performed parameter2 The second parameter w.r.t. which the marginalization is performed
- Returns:
- A histogram which contains the marginalized probability distribution (normalized to 1)
Definition at line 980 of file BCIntegrate.cxx.
{
int niter=fNbins*fNbins*fNSamplesPer2DBin*fNvar;
// set parameter to marginalize
int index1 = parameter1->GetIndex();
int index2 = parameter2->GetIndex();
// define histogram
double hmin1 = parameter1->GetLowerLimit();
double hmax1 = parameter1->GetUpperLimit();
double hmin2 = parameter2->GetLowerLimit();
double hmax2 = parameter2->GetUpperLimit();
TH2D * hist = new TH2D(Form("hist_%s_%s", parameter1->GetName().data(), parameter2->GetName().data()),"",
fNbins, hmin1, hmax1,
fNbins, hmin2, hmax2);
// prepare Metro
std::vector <double> randx;
randx.assign(fNvar, 0.0);
InitMetro();
for(int i=0;i<=niter;i++) {
GetRandomPointMetro(randx);
hist->Fill(randx[index1],randx[index2]);
}
// normalize
hist->Scale(1.0/hist->Integral("width"));
return hist;
}
| void BCIntegrate::MCMCFillErrorBand | ( | ) | [private] |
Fill error band histogram for curreent iteration. This method is called from MCMCIterationInterface()
Definition at line 2260 of file BCIntegrate.cxx.
{
if (!fFillErrorBand)
return;
// function fitting
if (fFitFunctionIndexX < 0)
return;
// loop over all possible x values ...
if (fErrorBandContinuous) {
double x = 0;
for (int ix = 0; ix < fErrorBandNbinsX; ix++) {
// calculate x
x = fErrorBandXY->GetXaxis()->GetBinCenter(ix + 1);
// calculate y
std::vector <double> xvec;
xvec.push_back(x);
// loop over all chains
for (int ichain = 0; ichain < MCMCGetNChains(); ++ichain) {
// calculate y
double y = FitFunction(xvec, MCMCGetx(ichain));
// fill histogram
fErrorBandXY->Fill(x, y);
}
xvec.clear();
}
}
// ... or evaluate at the data point x-values
else {
int ndatapoints = int(fErrorBandX.size());
double x = 0;
for (int ix = 0; ix < ndatapoints; ++ix) {
// calculate x
x = fErrorBandX.at(ix);
// calculate y
std::vector <double> xvec;
xvec.push_back(x);
// loop over all chains
for (int ichain = 0; ichain < MCMCGetNChains(); ++ichain) {
// calculate y
double y = FitFunction(xvec, MCMCGetx(ichain));
// fill histogram
fErrorBandXY->Fill(x, y);
}
xvec.clear();
}
}
}
| void BCIntegrate::MCMCIterationInterface | ( | ) | [private, virtual] |
Method executed for every iteration of the MCMC, overloaded from BCEngineMCMC.
Reimplemented from BCEngineMCMC.
Reimplemented in BCRooInterface.
Definition at line 2247 of file BCIntegrate.cxx.
{
// what's within this method will be executed
// for every iteration of the MCMC
// fill error band
MCMCFillErrorBand();
// do user defined stuff
MCMCUserIterationInterface();
}
| virtual void BCIntegrate::MCMCUserIterationInterface | ( | ) | [inline, virtual] |
Method executed for every iteration of the MCMC. User's code should be provided via overloading in the derived class
Reimplemented in BCGoFTest, and BCTemplateFitter.
Definition at line 826 of file BCIntegrate.h.
{}
| BCIntegrate & BCIntegrate::operator= | ( | const BCIntegrate & | bcintegrate | ) |
Defaut assignment operator
Definition at line 237 of file BCIntegrate.cxx.
{
BCEngineMCMC::operator=(bcintegrate);
fNvar = bcintegrate.fNvar;
fNbins = bcintegrate.fNbins;
fNSamplesPer2DBin = bcintegrate.fNSamplesPer2DBin;
fMarkovChainStepSize = bcintegrate.fMarkovChainStepSize;
fMarkovChainNIterations = bcintegrate.fMarkovChainNIterations;
fMarkovChainAutoN = bcintegrate.fMarkovChainAutoN;
if (bcintegrate.fDataPointLowerBoundaries)
fDataPointLowerBoundaries = new BCDataPoint(*bcintegrate.fDataPointLowerBoundaries);
else
fDataPointLowerBoundaries = 0;
if (bcintegrate.fDataPointUpperBoundaries)
fDataPointUpperBoundaries = new BCDataPoint(*bcintegrate.fDataPointUpperBoundaries);
else
fDataPointUpperBoundaries = 0;
fDataFixedValues = bcintegrate.fDataFixedValues;
fBestFitParameters = bcintegrate.fBestFitParameters;
fBestFitParameterErrors = bcintegrate.fBestFitParameterErrors;
fBestFitParametersMarginalized = bcintegrate.fBestFitParametersMarginalized;
for (int i = 0; i < int(bcintegrate.fHProb1D.size()); ++i) {
if (bcintegrate.fHProb1D.at(i))
fHProb1D.push_back(new TH1D(*(bcintegrate.fHProb1D.at(i))));
else
fHProb1D.push_back(0);
}
for (int i = 0; i < int(bcintegrate.fHProb2D.size()); ++i) {
if (bcintegrate.fHProb2D.at(i))
fHProb2D.push_back(new TH2D(*(fHProb2D.at(i))));
else
fHProb2D.push_back(0);
}
fFillErrorBand = bcintegrate.fFillErrorBand;
fFitFunctionIndexX = bcintegrate.fFitFunctionIndexX;
fFitFunctionIndexY = bcintegrate.fFitFunctionIndexY;
fErrorBandX = bcintegrate.fErrorBandX;
if (bcintegrate.fErrorBandXY)
fErrorBandXY = new TH2D(*(bcintegrate.fErrorBandXY));
else
fErrorBandXY = 0;
fErrorBandNbinsX = bcintegrate.fErrorBandNbinsX;
fErrorBandNbinsY = bcintegrate.fErrorBandNbinsY;
fMinuit = new TMinuit();
// debugKK
// *fMinuit = *(bcintegrate.fMinuit);
fMinuitArglist[0] = bcintegrate.fMinuitArglist[0];
fMinuitArglist[1] = bcintegrate.fMinuitArglist[1];
fMinuitErrorFlag = bcintegrate.fMinuitErrorFlag;
fFlagIgnorePrevOptimization = bcintegrate.fFlagIgnorePrevOptimization;
fFlagWriteMarkovChain = bcintegrate.fFlagWriteMarkovChain;
fMarkovChainTree = bcintegrate.fMarkovChainTree;
fMCMCIteration = bcintegrate.fMCMCIteration;
fSAT0 = bcintegrate.fSAT0;
fSATmin = bcintegrate.fSATmin;
// debugKK
fTreeSA = 0;
fFlagWriteSAToFile = bcintegrate.fFlagWriteSAToFile;
fSANIterations = bcintegrate.fSANIterations;
fSATemperature = bcintegrate.fSATemperature;
fSALogProb = bcintegrate.fSALogProb;
fSAx = bcintegrate.fSAx;
if (bcintegrate.fx)
fx = new BCParameterSet(*(bcintegrate.fx));
else
fx = 0;
fMin = new double[fNvar];
fMax = new double[fNvar];
fVarlist = new int[fNvar];
fNiterPerDimension = bcintegrate.fNiterPerDimension;
fIntegrationMethod = bcintegrate.fIntegrationMethod;
fMarginalizationMethod = bcintegrate.fMarginalizationMethod;
fOptimizationMethod = bcintegrate.fOptimizationMethod;
fOptimizationMethodMode = bcintegrate.fOptimizationMethodMode;
fSASchedule = bcintegrate.fSASchedule;
fNIterationsMax = bcintegrate.fNIterationsMax;
fNIterations = bcintegrate.fNIterations;
fRelativePrecision = bcintegrate.fRelativePrecision;
fAbsolutePrecision = bcintegrate.fAbsolutePrecision;
fCubaIntegrationMethod = bcintegrate.fCubaIntegrationMethod;
fCubaMinEval = bcintegrate.fCubaMinEval;
fCubaMaxEval = bcintegrate.fCubaMaxEval;
fCubaVerbosity = bcintegrate.fCubaVerbosity;
fCubaVegasNStart = bcintegrate.fCubaVegasNStart;
fCubaVegasNIncrease = bcintegrate.fCubaVegasNIncrease;
fCubaSuaveNNew = bcintegrate.fCubaSuaveNNew;
fCubaSuaveFlatness = bcintegrate.fCubaSuaveFlatness;
fError = bcintegrate.fError;
fNmetro = bcintegrate.fNmetro;
fNacceptedMCMC = bcintegrate.fNacceptedMCMC;
fXmetro0 = bcintegrate.fXmetro0;
fXmetro1 = bcintegrate.fXmetro1;
fMarkovChainValue = bcintegrate.fMarkovChainValue;
// return this
return *this;
}
| void BCIntegrate::ResetVarlist | ( | int | v | ) |
Sets all values of the variable list to a particular value The value
Definition at line 498 of file BCIntegrate.cxx.
| double BCIntegrate::SAHelperGetRadialCauchy | ( | ) |
Generates the radial part of a n-dimensional Cauchy distribution. Helper function for Cauchy annealing.
Definition at line 1856 of file BCIntegrate.cxx.
{
// theta is sampled from a rather complicated distribution,
// so first we create a lookup table with 10000 random numbers
// once and then, each time we need a new random number,
// we just look it up in the table.
double theta;
// static vectors for theta-sampling-map
static double map_u[10001];
static double map_theta[10001];
static bool initialized = false;
static int map_dimension = 0;
// is the lookup-table already initialized? if not, do it!
if (!initialized || map_dimension != fNvar) {
double init_theta;
double beta = SAHelperSinusToNIntegral(fNvar - 1, 1.57079632679);
for (int i = 0; i <= 10000; i++) {
init_theta = 3.14159265 * (double)i / 5000.;
map_theta[i] = init_theta;
map_u[i] = SAHelperSinusToNIntegral(fNvar - 1, init_theta) / beta;
}
map_dimension = fNvar;
initialized = true;
} // initializing is done.
// generate uniform random number for sampling
double u = fRandom->Rndm();
// Find the two elements just greater than and less than u
// using a binary search (O(log(N))).
int lo = 0;
int up = 10000;
int mid;
while (up != lo) {
mid = ((up - lo + 1) / 2) + lo;
if (u >= map_u[mid])
lo = mid;
else
up = mid - 1;
}
up++;
// perform linear interpolation:
theta = map_theta[lo] + (u - map_u[lo]) / (map_u[up] - map_u[lo]) * (map_theta[up] - map_theta[lo]);
return tan(theta);
}
| std::vector< double > BCIntegrate::SAHelperGetRandomPointOnHypersphere | ( | ) |
Generates a uniform distributed random point on the surface of a fNvar-dimensional Hypersphere. Used as a helper to generate proposal points for Cauchy annealing.
Definition at line 1811 of file BCIntegrate.cxx.
{
std::vector<double> rand_point(fNvar);
// This method can only be called with fNvar >= 2 since the 1-dim case
// is already hard wired into the Cauchy annealing proposal function.
// To speed things up, hard-code fast method for 2 and dimensions.
// The algorithm for 2D can be found at
// http://mathworld.wolfram.com/CirclePointPicking.html
// For 3D just using ROOT's algorithm.
if (fNvar == 2) {
double x1, x2, s;
do {
x1 = fRandom->Rndm() * 2. - 1.;
x2 = fRandom->Rndm() * 2. - 1.;
s = x1*x1 + x2*x2;
}
while (s >= 1);
rand_point[0] = (x1*x1 - x2*x2) / s;
rand_point[1] = (2.*x1*x2) / s;
}
else if (fNvar == 3) {
fRandom->Sphere(rand_point[0], rand_point[1], rand_point[2], 1.0);
}
else {
double s = 0.,
gauss_num;
for (int i = 0; i < fNvar; i++) {
gauss_num = fRandom->Gaus();
rand_point[i] = gauss_num;
s += gauss_num * gauss_num;
}
s = sqrt(s);
for (int i = 0; i < fNvar; i++)
rand_point[i] = rand_point[i] / s;
}
return rand_point;
}
| double BCIntegrate::SAHelperSinusToNIntegral | ( | int | dim, | |
| double | theta | |||
| ) |
Returns the Integral of sin^dim from 0 to theta. Helper function needed for generating Cauchy distributions.
Definition at line 1913 of file BCIntegrate.cxx.
{
if (dim < 1)
return theta;
else if (dim == 1)
return (1. - cos(theta));
else if (dim == 2)
return 0.5 * (theta - sin(theta) * cos(theta));
else if (dim == 3)
return (2. - sin(theta) * sin(theta) * cos(theta) - 2. * cos(theta)) / 3.;
else
return - pow(sin(theta), (double)(dim - 1)) * cos(theta) / (double)dim
+ (double)(dim - 1) / (double)dim
* SAHelperSinusToNIntegral(dim - 2, theta);
}
| void BCIntegrate::SAInitialize | ( | ) |
Definition at line 2320 of file BCIntegrate.cxx.
| double BCIntegrate::SATemperature | ( | double | t | ) |
Temperature annealing schedule for use with Simulated Annealing. Delegates to the appropriate method according to fSASchedule.
- Parameters:
-
t iterator for lowering the temperature over time.
Definition at line 1703 of file BCIntegrate.cxx.
{
// do we have Cauchy (default), Boltzmann or custom annealing schedule?
if (fSASchedule == BCIntegrate::kSABoltzmann)
return SATemperatureBoltzmann(t);
else if (fSASchedule == BCIntegrate::kSACauchy)
return SATemperatureCauchy(t);
else
return SATemperatureCustom(t);
}
| double BCIntegrate::SATemperatureBoltzmann | ( | double | t | ) |
Temperature annealing schedule for use with Simulated Annealing. This method is used for Boltzmann annealing schedule.
- Parameters:
-
t iterator for lowering the temperature over time.
Definition at line 1716 of file BCIntegrate.cxx.
{
return fSAT0 / log((double)(t + 1));
}
| double BCIntegrate::SATemperatureCauchy | ( | double | t | ) |
Temperature annealing schedule for use with Simulated Annealing. This method is used for Cauchy annealing schedule.
- Parameters:
-
t iterator for lowering the temperature over time.
Definition at line 1723 of file BCIntegrate.cxx.
{
return fSAT0 / (double)t;
}
| double BCIntegrate::SATemperatureCustom | ( | double | t | ) | [virtual] |
Temperature annealing schedule for use with Simulated Annealing. This is a virtual method to be overridden by a user-defined custom temperature schedule.
- Parameters:
-
t iterator for lowering the temperature over time.
Definition at line 1730 of file BCIntegrate.cxx.
{
BCLog::OutError("BCIntegrate::SATemperatureCustom : No custom temperature schedule defined");
return 0.;
}
| void BCIntegrate::SetAbsolutePrecision | ( | double | absprecision | ) | [inline] |
Set absolute precision of the numerical integation
Definition at line 362 of file BCIntegrate.h.
{ fAbsolutePrecision = absprecision; }
| void BCIntegrate::SetCubaIntegrationMethod | ( | BCIntegrate::BCCubaMethod | type | ) |
Set Cuba integration method
Definition at line 1999 of file BCIntegrate.cxx.
{
#ifdef HAVE_CUBA_H
switch(type) {
case BCIntegrate::kCubaVegas:
case BCIntegrate::kCubaSuave:
fCubaIntegrationMethod = type;
return;
case BCIntegrate::kCubaDivonne:
case BCIntegrate::kCubaCuhre:
BCLog::OutError(TString::Format(
"BAT does not yet support global setting of Cuba integration method to %s. "
"To use this method use explicit call to CubaIntegrate() with arguments.",
DumpCubaIntegrationMethod(type).c_str()));
return;
default:
BCLog::OutError(TString::Format("Integration method of type %d is not defined for Cuba",type));
return;
}
#else
BCLog::OutWarning("!!! This version of BAT is compiled without Cuba.");
BCLog::OutWarning(" Setting Cuba integration method will have no effect.");
#endif
}
| void BCIntegrate::SetCubaMaxEval | ( | int | n | ) | [inline] |
Set maximum number of evaluations in Cuba integration
Definition at line 376 of file BCIntegrate.h.
{ fCubaMaxEval = n; }
| void BCIntegrate::SetCubaMinEval | ( | int | n | ) | [inline] |
Set minimum number of evaluations in Cuba integration
Definition at line 371 of file BCIntegrate.h.
{ fCubaMinEval = n; }
| void BCIntegrate::SetCubaSuaveFlatness | ( | double | p | ) | [inline] |
Set flatness for Cuba Suave
Definition at line 401 of file BCIntegrate.h.
{ fCubaSuaveFlatness = p; }
| void BCIntegrate::SetCubaSuaveNNew | ( | int | n | ) | [inline] |
Set number of new integrand evaluations in each subdivision for Cuba Suave
Definition at line 396 of file BCIntegrate.h.
{ fCubaSuaveNNew = n; }
| void BCIntegrate::SetCubaVegasNIncrease | ( | int | n | ) | [inline] |
Set increase in number of evaluations per iteration for Cuba Vegas
Definition at line 391 of file BCIntegrate.h.
{ fCubaVegasNIncrease = n; }
| void BCIntegrate::SetCubaVegasNStart | ( | int | n | ) | [inline] |
Set initial number of evaluations per iteration for Cuba Vegas
Definition at line 386 of file BCIntegrate.h.
{ fCubaVegasNStart = n; }
| void BCIntegrate::SetCubaVerbosityLevel | ( | int | n | ) | [inline] |
Set verbosity level of Cuba integration
Definition at line 381 of file BCIntegrate.h.
{ fCubaVerbosity = n; }
| void BCIntegrate::SetDataPointLowerBoundaries | ( | BCDataPoint * | datasetlowerboundaries | ) | [inline] |
Sets the data point containing the lower boundaries of possible data values
Definition at line 446 of file BCIntegrate.h.
{ fDataPointLowerBoundaries = datasetlowerboundaries; }
| void BCIntegrate::SetDataPointLowerBoundary | ( | int | index, | |
| double | lowerboundary | |||
| ) | [inline] |
Sets the lower boundary of possible data values for a particular variable
Definition at line 458 of file BCIntegrate.h.
{ fDataPointLowerBoundaries -> SetValue(index, lowerboundary); }
| void BCIntegrate::SetDataPointUpperBoundaries | ( | BCDataPoint * | datasetupperboundaries | ) | [inline] |
Sets the data point containing the upper boundaries of possible data values
Definition at line 452 of file BCIntegrate.h.
{ fDataPointUpperBoundaries = datasetupperboundaries; }
| void BCIntegrate::SetDataPointUpperBoundary | ( | int | index, | |
| double | upperboundary | |||
| ) | [inline] |
Sets the upper boundary of possible data values for a particular variable
Definition at line 464 of file BCIntegrate.h.
{ fDataPointUpperBoundaries -> SetValue(index, upperboundary); }
| void BCIntegrate::SetErrorBandHisto | ( | TH2D * | h | ) | [inline] |
| void BCIntegrate::SetFillErrorBand | ( | bool | flag = true |
) | [inline] |
Turn on or off the filling of the error band during the MCMC run.
- Parameters:
-
flag set to true for turning on the filling
Definition at line 418 of file BCIntegrate.h.
{ fFillErrorBand=flag; }
| void BCIntegrate::SetFitFunctionIndexX | ( | int | index | ) | [inline] |
Sets index of the x values in function fits.
- Parameters:
-
index Index of the x values
Definition at line 430 of file BCIntegrate.h.
{ fFitFunctionIndexX = index; }
| void BCIntegrate::SetFitFunctionIndexY | ( | int | index | ) | [inline] |
Sets index of the y values in function fits.
- Parameters:
-
index Index of the y values
Definition at line 436 of file BCIntegrate.h.
{ fFitFunctionIndexY = index; }
| void BCIntegrate::SetFitFunctionIndices | ( | int | indexx, | |
| int | indexy | |||
| ) | [inline] |
Definition at line 439 of file BCIntegrate.h.
{ this -> SetFitFunctionIndexX(indexx);
this -> SetFitFunctionIndexY(indexy); }
| void BCIntegrate::SetFlagIgnorePrevOptimization | ( | bool | flag | ) | [inline] |
Definition at line 295 of file BCIntegrate.h.
{ fFlagIgnorePrevOptimization = flag; }
| void BCIntegrate::SetFlagWriteSAToFile | ( | bool | flag | ) | [inline] |
Definition at line 519 of file BCIntegrate.h.
{ fFlagWriteSAToFile = flag; }
| void BCIntegrate::SetIntegrationMethod | ( | BCIntegrate::BCIntegrationMethod | method | ) |
- Parameters:
-
method The integration method
Definition at line 541 of file BCIntegrate.cxx.
{
#ifdef HAVE_CUBA_H
fIntegrationMethod = method;
#else
BCLog::OutWarning("!!! This version of BAT is compiled without Cuba.");
BCLog::OutWarning(" Monte Carlo Sampled Mean integration method will be used.");
BCLog::OutWarning(" To be able to use Cuba integration, install Cuba and recompile BAT.");
#endif
}
| void BCIntegrate::SetMarginalizationMethod | ( | BCIntegrate::BCMarginalizationMethod | method | ) | [inline] |
- Parameters:
-
method The marginalization method
Definition at line 318 of file BCIntegrate.h.
{ fMarginalizationMethod = method; }
| void BCIntegrate::SetMarkovChainAutoN | ( | bool | flag | ) | [inline] |
Sets a flag for automatically calculating the number of iterations
Definition at line 497 of file BCIntegrate.h.
{ fMarkovChainAutoN = flag; }
| void BCIntegrate::SetMarkovChainInitialPosition | ( | std::vector< double > | position | ) |
Sets the initial position for the Markov chain
Definition at line 401 of file BCIntegrate.cxx.
| void BCIntegrate::SetMarkovChainNIterations | ( | int | niterations | ) | [inline] |
Sets the number of iterations in the markov chain
Definition at line 491 of file BCIntegrate.h.
{ fMarkovChainNIterations = niterations;
fMarkovChainAutoN = false; }
| void BCIntegrate::SetMarkovChainStepSize | ( | double | stepsize | ) | [inline] |
Sets the step size for Markov chains
Definition at line 486 of file BCIntegrate.h.
{ fMarkovChainStepSize = stepsize; }
| void BCIntegrate::SetMarkovChainTree | ( | TTree * | tree | ) | [inline] |
Sets the ROOT tree containing the Markov chain
Definition at line 477 of file BCIntegrate.h.
{ fMarkovChainTree = tree; }
| void BCIntegrate::SetMinuitArlist | ( | double * | arglist | ) | [inline] |
Definition at line 291 of file BCIntegrate.h.
{ fMinuitArglist[0] = arglist[0];
fMinuitArglist[1] = arglist[1]; }
| void BCIntegrate::SetMode | ( | std::vector< double > | mode | ) |
Sets mode
Definition at line 1931 of file BCIntegrate.cxx.
{
if((int)mode.size() == fNvar) {
fBestFitParameters.clear();
for (int i = 0; i < fNvar; i++)
fBestFitParameters.push_back(mode[i]);
}
}
| void BCIntegrate::SetNbins | ( | int | nbins, | |
| int | index = -1 | |||
| ) |
- Parameters:
-
n Number of bins per dimension for the marginalized distributions. Default is 100. Minimum number allowed is 2. Set the number of bins for the marginalized distribution of a parameter. nbins Number of bins (default = 100) index Index of the parameter.
Definition at line 430 of file BCIntegrate.cxx.
{
if (fNvar == 0)
return;
// check if index is in range
if (index >= fNvar) {
BCLog::OutWarning("BCIntegrate::SetNbins : Index out of range.");
return;
}
// set for all parameters at once
else if (index < 0) {
for (int i = 0; i < fNvar; ++i)
SetNbins(nbins, i);
return;
}
// sanity check for nbins
if (nbins <= 0)
nbins = 10;
fMCMCH1NBins[index] = nbins;
return;
// if(n<2) {
// BCLog::OutWarning("BCIntegrate::SetNbins. Number of bins less than 2 makes no sense. Setting to 2.");
// n=2;
// }
// MCMCSetH1NBins(n, -1);
// fNbins=n;
// fMCMCH1NBins = n;
// fMCMCH2NBinsX = n;
// fMCMCH2NBinsY = n;
}
| void BCIntegrate::SetNIterationsMax | ( | int | niterations | ) | [inline] |
- Parameters:
-
niterations The maximum number of iterations for Monte Carlo integration
Definition at line 351 of file BCIntegrate.h.
{ fNIterationsMax = niterations; }
| void BCIntegrate::SetNiterationsPerDimension | ( | int | niterations | ) | [inline] |
- Parameters:
-
niterations Number of iterations per dimension for Monte Carlo integration.
Definition at line 340 of file BCIntegrate.h.
{ fNiterPerDimension = niterations; }
| void BCIntegrate::SetNSamplesPer2DBin | ( | int | n | ) | [inline] |
- Parameters:
-
n Number of samples per 2D bin per variable in the Metropolis marginalization. Default is 100.
Definition at line 346 of file BCIntegrate.h.
{ fNSamplesPer2DBin = n; }
| void BCIntegrate::SetOptimizationMethod | ( | BCIntegrate::BCOptimizationMethod | method | ) | [inline] |
- Parameters:
-
method The mode finding method
Definition at line 323 of file BCIntegrate.h.
{ fOptimizationMethod = method; }
| void BCIntegrate::SetOptimizationMethodMode | ( | BCIntegrate::BCOptimizationMethod | method | ) | [inline] |
- Parameters:
-
method The mode finding method that was used to find the current best estimate of the parameters at the mode
Definition at line 329 of file BCIntegrate.h.
{ fOptimizationMethodMode = method; }
| void BCIntegrate::SetParameters | ( | BCParameterSet * | par | ) |
- Parameters:
-
par The parameter set which gets translated into array needed for the Monte Carlo integration
Definition at line 361 of file BCIntegrate.cxx.
{
DeleteVarList();
fx = par;
fNvar = fx->size();
fMin = new double[fNvar];
fMax = new double[fNvar];
fVarlist = new int[fNvar];
ResetVarlist(1);
for(int i=0;i<fNvar;i++) {
fMin[i]=(fx->at(i))->GetLowerLimit();
fMax[i]=(fx->at(i))->GetUpperLimit();
}
fXmetro0.clear();
for(int i=0;i<fNvar;i++)
fXmetro0.push_back((fMin[i]+fMax[i])/2.0);
fXmetro1.clear();
fXmetro1.assign(fNvar,0.);
fMCMCBoundaryMin.clear();
fMCMCBoundaryMax.clear();
for(int i=0;i<fNvar;i++) {
fMCMCBoundaryMin.push_back(fMin[i]);
fMCMCBoundaryMax.push_back(fMax[i]);
fMCMCFlagsFillHistograms.push_back(true);
}
for (int i = int(fMCMCH1NBins.size()); i<fNvar; ++i)
fMCMCH1NBins.push_back(100);
fMCMCNParameters = fNvar;
}
| void BCIntegrate::SetRelativePrecision | ( | double | relprecision | ) | [inline] |
- Parameters:
-
relprecision The relative precision envisioned for Monte Carlo integration
Definition at line 357 of file BCIntegrate.h.
{ fRelativePrecision = relprecision; }
| void BCIntegrate::SetSASchedule | ( | BCIntegrate::BCSASchedule | schedule | ) | [inline] |
- Parameters:
-
method The Simulated Annealing schedule
Definition at line 335 of file BCIntegrate.h.
{ fSASchedule = schedule; }
| void BCIntegrate::SetSAT0 | ( | double | T0 | ) | [inline] |
- Parameters:
-
T0 new value for Simulated Annealing starting temperature.
Definition at line 511 of file BCIntegrate.h.
{ fSAT0 = T0; }
| void BCIntegrate::SetSATmin | ( | double | Tmin | ) | [inline] |
- Parameters:
-
Tmin new value for Simulated Annealing threshold temperature.
Definition at line 516 of file BCIntegrate.h.
{ fSATmin = Tmin; }
| void BCIntegrate::SetSATree | ( | TTree * | tree | ) | [inline] |
Sets the tree containing the Simulated Annealing chain.
Definition at line 524 of file BCIntegrate.h.
{ fTreeSA = tree; }
| void BCIntegrate::SetVar | ( | int | index | ) | [inline] |
- Parameters:
-
index The index of the variable to be set
Definition at line 309 of file BCIntegrate.h.
{fVarlist[index]=1;}
| void BCIntegrate::SetVarList | ( | int * | varlist | ) |
- Parameters:
-
varlist A list of parameters
Definition at line 491 of file BCIntegrate.cxx.
| void BCIntegrate::UnsetFillErrorBand | ( | ) | [inline] |
Turn off filling of the error band during the MCMC run. This method is equivalent to SetFillErrorBand(false)
Definition at line 424 of file BCIntegrate.h.
{ fFillErrorBand=false; }
| void BCIntegrate::UnsetVar | ( | int | index | ) | [inline] |
Set value of a particular integration variable to 0.
- Parameters:
-
index The index of the variable
Definition at line 553 of file BCIntegrate.h.
{ fVarlist[index] = 0; }
| void BCIntegrate::WriteMarkovChain | ( | bool | flag | ) | [inline] |
Flag for writing Markov chain to ROOT file (true) or not (false)
Definition at line 470 of file BCIntegrate.h.
{ fFlagWriteMarkovChain = flag;
fMCMCFlagWriteChainToFile = flag;
fMCMCFlagWritePreRunToFile = flag; }
Member Data Documentation
double BCIntegrate::fAbsolutePrecision [private] |
Requested relative precision of the integation
Definition at line 1031 of file BCIntegrate.h.
std::vector<double> BCIntegrate::fBestFitParameterErrors [protected] |
Definition at line 892 of file BCIntegrate.h.
std::vector<double> BCIntegrate::fBestFitParameters [protected] |
A vector of best fit parameters estimated from the global probability and the estimate on their uncertainties
Definition at line 891 of file BCIntegrate.h.
std::vector<double> BCIntegrate::fBestFitParametersMarginalized [protected] |
A vector of best fit parameters estimated from the marginalized probability
Definition at line 896 of file BCIntegrate.h.
Cuba integration method
Definition at line 1034 of file BCIntegrate.h.
int BCIntegrate::fCubaMaxEval [private] |
Maximum number of evaluations in Cuba integration
Definition at line 1040 of file BCIntegrate.h.
int BCIntegrate::fCubaMinEval [private] |
Minimum number of evaluations in Cuba integration
Definition at line 1037 of file BCIntegrate.h.
double BCIntegrate::fCubaSuaveFlatness [private] |
Flatness for Cuba Suave
Definition at line 1055 of file BCIntegrate.h.
int BCIntegrate::fCubaSuaveNNew [private] |
Number of new integrand evaluations in each subdivision for Cuba Suave
Definition at line 1052 of file BCIntegrate.h.
int BCIntegrate::fCubaVegasNIncrease [private] |
Increase in number of evaluations per iteration for Cuba Vegas
Definition at line 1049 of file BCIntegrate.h.
int BCIntegrate::fCubaVegasNStart [private] |
Initial number of evaluations per iteration for Cuba Vegas
Definition at line 1046 of file BCIntegrate.h.
int BCIntegrate::fCubaVerbosity [private] |
Verbosity level of Cuba integration
Definition at line 1043 of file BCIntegrate.h.
std::vector<bool> BCIntegrate::fDataFixedValues [protected] |
Definition at line 886 of file BCIntegrate.h.
BCDataPoint* BCIntegrate::fDataPointLowerBoundaries [protected] |
data point containing the lower boundaries of possible data values
Definition at line 880 of file BCIntegrate.h.
BCDataPoint* BCIntegrate::fDataPointUpperBoundaries [protected] |
data point containing the upper boundaries of possible data values
Definition at line 884 of file BCIntegrate.h.
double BCIntegrate::fError [private] |
The uncertainty in the most recent Monte Carlo integration
Definition at line 1059 of file BCIntegrate.h.
bool BCIntegrate::fErrorBandContinuous [protected] |
A flag for single point evaluation of the error "band"
Definition at line 917 of file BCIntegrate.h.
int BCIntegrate::fErrorBandNbinsX [protected] |
Number of X bins of the error band histogram
Definition at line 926 of file BCIntegrate.h.
int BCIntegrate::fErrorBandNbinsY [protected] |
Nnumber of Y bins of the error band histogram
Definition at line 930 of file BCIntegrate.h.
std::vector<double> BCIntegrate::fErrorBandX [protected] |
Definition at line 918 of file BCIntegrate.h.
TH2D* BCIntegrate::fErrorBandXY [protected] |
The error band histogram
Definition at line 922 of file BCIntegrate.h.
bool BCIntegrate::fFillErrorBand [protected] |
Flag whether or not to fill the error band
Definition at line 908 of file BCIntegrate.h.
int BCIntegrate::fFitFunctionIndexX [protected] |
The indices for function fits
Definition at line 912 of file BCIntegrate.h.
int BCIntegrate::fFitFunctionIndexY [protected] |
Definition at line 913 of file BCIntegrate.h.
double BCIntegrate::fFlagIgnorePrevOptimization [protected] |
Flag for ignoring older results of minimization
Definition at line 941 of file BCIntegrate.h.
bool BCIntegrate::fFlagWriteMarkovChain [protected] |
Flag for writing Markov chain to file
Definition at line 945 of file BCIntegrate.h.
bool BCIntegrate::fFlagWriteSAToFile [protected] |
Flag deciding whether to write SA to file or not.
Definition at line 969 of file BCIntegrate.h.
std::vector<TH1D *> BCIntegrate::fHProb1D [protected] |
Vector of TH1D histograms for marginalized probability distributions
Definition at line 900 of file BCIntegrate.h.
std::vector<TH2D *> BCIntegrate::fHProb2D [protected] |
Vector of TH2D histograms for marginalized probability distributions
Definition at line 904 of file BCIntegrate.h.
Current integration method
Definition at line 1000 of file BCIntegrate.h.
Current marginalization method
Definition at line 1004 of file BCIntegrate.h.
bool BCIntegrate::fMarkovChainAutoN [protected] |
Definition at line 876 of file BCIntegrate.h.
int BCIntegrate::fMarkovChainNIterations [protected] |
Definition at line 874 of file BCIntegrate.h.
double BCIntegrate::fMarkovChainStepSize [protected] |
Step size in the Markov chain relative to min and max
Definition at line 872 of file BCIntegrate.h.
TTree* BCIntegrate::fMarkovChainTree [protected] |
ROOT tree containing the Markov chain
Definition at line 949 of file BCIntegrate.h.
double BCIntegrate::fMarkovChainValue [private] |
A double containing the log likelihood value at the point fXmetro1
Definition at line 1076 of file BCIntegrate.h.
double* BCIntegrate::fMax [private] |
Array containing the upper boundaries of the variables to integrate over.
Definition at line 988 of file BCIntegrate.h.
int BCIntegrate::fMCMCIteration [protected] |
Iteration of the MCMC
Definition at line 953 of file BCIntegrate.h.
double* BCIntegrate::fMin [private] |
Array containing the lower boundaries of the variables to integrate over.
Definition at line 984 of file BCIntegrate.h.
TMinuit* BCIntegrate::fMinuit [protected] |
Minuit
Definition at line 934 of file BCIntegrate.h.
double BCIntegrate::fMinuitArglist[2] [protected] |
Definition at line 936 of file BCIntegrate.h.
int BCIntegrate::fMinuitErrorFlag [protected] |
Definition at line 937 of file BCIntegrate.h.
int BCIntegrate::fNacceptedMCMC [private] |
Definition at line 1064 of file BCIntegrate.h.
int BCIntegrate::fNbins [protected] |
Number of bins per dimension for the marginalized distributions
Definition at line 863 of file BCIntegrate.h.
int BCIntegrate::fNIterations [private] |
Number of iterations in the most recent Monte Carlo integation
Definition at line 1025 of file BCIntegrate.h.
int BCIntegrate::fNIterationsMax [private] |
Maximum number of iterations
Definition at line 1021 of file BCIntegrate.h.
int BCIntegrate::fNiterPerDimension [private] |
Number of iteration per dimension for Monte Carlo integration.
Definition at line 996 of file BCIntegrate.h.
int BCIntegrate::fNmetro [private] |
The number of iterations in the Metropolis integration
Definition at line 1063 of file BCIntegrate.h.
int BCIntegrate::fNSamplesPer2DBin [protected] |
Number of samples per 2D bin per variable in the Metropolis marginalization.
Definition at line 868 of file BCIntegrate.h.
int BCIntegrate::fNvar [protected] |
Number of variables to integrate over.
Definition at line 859 of file BCIntegrate.h.
Current mode finding method
Definition at line 1008 of file BCIntegrate.h.
Method with which the global mode was found (can differ from fOptimization method in case more than one algorithm is used).
Definition at line 1013 of file BCIntegrate.h.
double BCIntegrate::fRelativePrecision [private] |
Requested relative precision of the integation
Definition at line 1028 of file BCIntegrate.h.
double BCIntegrate::fSALogProb [protected] |
Definition at line 973 of file BCIntegrate.h.
int BCIntegrate::fSANIterations [protected] |
Definition at line 971 of file BCIntegrate.h.
Current Simulated Annealing schedule
Definition at line 1017 of file BCIntegrate.h.
double BCIntegrate::fSAT0 [protected] |
Starting temperature for Simulated Annealing
Definition at line 957 of file BCIntegrate.h.
double BCIntegrate::fSATemperature [protected] |
Definition at line 972 of file BCIntegrate.h.
double BCIntegrate::fSATmin [protected] |
Minimal/Threshold temperature for Simulated Annealing
Definition at line 961 of file BCIntegrate.h.
std::vector<double> BCIntegrate::fSAx [protected] |
Definition at line 974 of file BCIntegrate.h.
TTree* BCIntegrate::fTreeSA [protected] |
Tree for the Simulated Annealing
Definition at line 965 of file BCIntegrate.h.
int* BCIntegrate::fVarlist [private] |
List of variables containing a flag whether to integrate over them or not.
Definition at line 992 of file BCIntegrate.h.
BCParameterSet* BCIntegrate::fx [private] |
Set of parameters for the integration.
Definition at line 980 of file BCIntegrate.h.
std::vector<double> BCIntegrate::fXmetro0 [private] |
A vector of points in parameter space used for the Metropolis algorithm
Definition at line 1068 of file BCIntegrate.h.
std::vector<double> BCIntegrate::fXmetro1 [private] |
A vector of points in parameter space used for the Metropolis algorithm
Definition at line 1072 of file BCIntegrate.h.
The documentation for this class was generated from the following files:
