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:
- An error code
Definition at line 2029 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.");
return 1;
#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 2112 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 {
BCLog::OutError(" Integration method not available. Set integral to -1e99.");
integral[0] = -1e99;
}
if (fail != 0) {
BCLog::OutWarning("Warning, integral did not converge with the given set of parameters. ");
BCLog::OutWarning(Form(" nevel : %d", neval));
BCLog::OutWarning(Form(" fail : %d", fail));
BCLog::OutWarning(Form(" integral[0] : %f", integral[0]));
BCLog::OutWarning(Form(" error[0] : %f", error[0]));
BCLog::OutWarning(Form(" prob[0] : %f", prob[0]));
}
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 | ( | ) |
Calculate integral using the Cuba library. For details see documentation.
- Returns:
- The integral
Definition at line 2070 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 | ) |
Return string with the name for a given Cuba integration type.
- Parameters:
-
type code for the Cuba integration type
- Returns:
- string containing the name of the Cuba integration type
Definition at line 2378 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] |
Return string with the name for the currently set Cuba integration type.
- Returns:
- string containing the name of the Cuba integration type
Definition at line 887 of file BCIntegrate.h.
{ return DumpCubaIntegrationMethod(fCubaIntegrationMethod); }
| std::string BCIntegrate::DumpIntegrationMethod | ( | BCIntegrate::BCIntegrationMethod | type | ) |
Return string with the name for a given integration type.
- Parameters:
-
type code for the integration type
- Returns:
- string containing the name of the integration type
Definition at line 2333 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] |
Return string with the name for the currently set integration type.
- Returns:
- string containing the name of the integration type
Definition at line 845 of file BCIntegrate.h.
{ return DumpIntegrationMethod(fIntegrationMethod); }
| std::string BCIntegrate::DumpMarginalizationMethod | ( | ) | [inline] |
Return string with the name for the currently set marginalization type.
- Returns:
- string containing the name of the marginalization type
Definition at line 857 of file BCIntegrate.h.
{ return DumpMarginalizationMethod(fMarginalizationMethod); }
| std::string BCIntegrate::DumpMarginalizationMethod | ( | BCIntegrate::BCMarginalizationMethod | type | ) |
Return string with the name for a given marginalization type.
- Parameters:
-
type code for the marginalization type
- Returns:
- string containing the name of the marginalization type
Definition at line 2350 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 | ) |
Return string with the name for a given optimization type.
- Parameters:
-
type code for the optimization type
- Returns:
- string containing the name of the optimization type
Definition at line 2363 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] |
Return string with the name for the currently set optimization type.
- Returns:
- string containing the name of the optimization type
Definition at line 869 of file BCIntegrate.h.
{ return DumpOptimizationMethod(fOptimizationMethod); }
| std::string BCIntegrate::DumpUsedOptimizationMethod | ( | ) | [inline] |
Return string with the name for the optimization type used to find the current mode.
- Returns:
- string containing the name of the optimization type
Definition at line 875 of file BCIntegrate.h.
{ return DumpOptimizationMethod(fOptimizationMethodMode); }
| double BCIntegrate::Eval | ( | const 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::EvalSampling | ( | const 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 1946 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 1965 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 1473 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
// 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 1588 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 | ( | const std::vector< double > & | , | |
| const 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 593 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 1207 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 1259 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 1223 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 265 of file BCIntegrate.h.
{ return &fXmetro1; }
| TTree* BCIntegrate::GetMarkovChainTree | ( | ) | [inline] |
- Returns:
- The ROOT tree containing the Markov chain
Definition at line 260 of file BCIntegrate.h.
{ return fMarkovChainTree; }
| double* BCIntegrate::GetMarkovChainValue | ( | ) | [inline] |
Returns the value of the loglikelihood at the point fXmetro1
Definition at line 275 of file BCIntegrate.h.
{ return &fMarkovChainValue; }
| int* BCIntegrate::GetMCMCIteration | ( | ) | [inline] |
Returns the iteration of the MCMC
Definition at line 270 of file BCIntegrate.h.
{ return &fMCMCIteration; }
| TMinuit * BCIntegrate::GetMinuit | ( | ) |
- Returns:
- Minuit used for mode finding
Definition at line 1463 of file BCIntegrate.cxx.
| int BCIntegrate::GetMinuitErrorFlag | ( | ) | [inline] |
- Returns:
- Error flag from Minuit run
Definition at line 255 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 | ( | const 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 1748 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 | ( | const 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 1760 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 | ( | const 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 1776 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 | ( | const 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 1810 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 1279 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 1290 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 1333 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 1388 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 1449 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 280 of file BCIntegrate.h.
{ return fSAT0; }
| double BCIntegrate::GetSATmin | ( | ) | [inline] |
Returns the Simulated Annealing threshhold temperature.
Definition at line 285 of file BCIntegrate.h.
{ return fSATmin; }
| TTree* BCIntegrate::GetSATree | ( | ) | [inline] |
Getter for the tree containing the Simulated Annealing chain.
Definition at line 533 of file BCIntegrate.h.
{ return fTreeSA; }
| void BCIntegrate::InitializeSATree | ( | ) |
Initialization of the tree for the Simulated Annealing
Definition at line 1572 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 1310 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 | ( | const 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 779 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 | ( | const std::vector< double > & | x, | |
| int * | varlist | |||
| ) |
Perfoms a 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 594 of file BCIntegrate.cxx.
{
SetVarList(varlist);
return IntegralMC(x);
}
| double BCIntegrate::IntegralMC | ( | const std::vector< double > & | x | ) |
Perfoms a Monte Carlo integration. For details see documentation.
- Parameters:
-
x An initial point in parameter space
- Returns:
- The integral
Definition at line 601 of file BCIntegrate.cxx.
{
// count the variables to integrate over
// and calculate D (integrated volume)
int NvarNow=0;
double D=1.;
for(int i = 0; i < fNvar; i++)
if(fVarlist[i]) {
NvarNow++;
D *= fMax[i] - fMin[i];
}
// print to log
BCLog::LogLevel level=BCLog::summary;
if(fNvar!=NvarNow) {
level=BCLog::detail;
BCLog::OutDetail(Form("Running MC integation over %i dimensions out of %i.", NvarNow, fNvar));
BCLog::OutDetail(" --> Fixed parameters:");
for(int i = 0; i < fNvar; i++)
if(!fVarlist[i])
BCLog::OutDetail(Form(" %3i : %g", i, x[i]));
}
else
BCLog::OutSummary(Form("Running MC integation over %i dimensions.", NvarNow));
BCLog::Out(level, Form(" --> Maximum number of iterations: %i", GetNIterationsMax()));
BCLog::Out(level, Form(" --> Aimed relative precision: %e", GetRelativePrecision()));
// reset variables
double pmax = 0.;
double sumW = 0.;
double sumW2 = 0.;
double integral = 0.;
double variance = 0.;
double relprecision = 1.;
std::vector<double> randx;
randx.assign(fNvar, 0.);
// how often to print out the info line to screen
int nwrite = fNIterationsMax/10;
if(nwrite < 10000)
nwrite=1000;
else if(nwrite < 100000)
nwrite=10000;
else if(nwrite < 1000000)
nwrite=100000;
else
nwrite=1000000;
// 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));
}
if (fNIterations%1000 == 0 || fNIterations%nwrite == 0) {
// calculate integral and variance
integral = D * sumW / fNIterations;
variance = (1.0 / double(fNIterations)) * (D * D * sumW2 / double(fNIterations) - integral * integral);
double error = sqrt(variance);
relprecision = error / integral;
if (fNIterations%nwrite == 0)
BCLog::OutDetail(Form("BCIntegrate::IntegralMC. Iteration %i, integral: %e +- %e.", fNIterations, integral, error));
}
}
integral = D * sumW / fNIterations;
fError = variance / fNIterations;
// print to log
BCLog::OutSummary(Form(" --> Result of integration: %e +- %e in %i iterations.", integral, sqrt(variance), fNIterations));
BCLog::OutSummary(Form(" --> Obtained relative precision: %e. ", sqrt(variance) / integral));
return integral;
}
| double BCIntegrate::IntegralMetro | ( | const 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 713 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 555 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.");
break;
#endif
default:
BCLog::OutError(
Form("BCIntegrate::Integrate : Invalid integration method: %d", fIntegrationMethod));
break;
}
return 0;
}
| int BCIntegrate::IntegrateResetResults | ( | ) |
Reset all information on the best fit parameters.
- Returns:
- An error code
Definition at line 544 of file BCIntegrate.cxx.
{
fBestFitParameters.clear();
fBestFitParameterErrors.clear();
fBestFitParametersMarginalized.clear();
// no errors
return 1;
}
| double BCIntegrate::LogEval | ( | const 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 | ( | const 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 843 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));
break;
}
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 866 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));
break;
}
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 1032 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 886 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);
// remember i = 0 => underflow bin
hist->SetBinContent(i+1, 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 918 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);
// remember i = 0 => underflow bin
hist->SetBinContent(i+1, j+1, 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 962 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 996 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 2266 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 2253 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, BCTemplateFitter, and BCMTF.
Definition at line 828 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 1861 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 1817 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 1917 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 | ( | ) |
Initializes the Simulated Annealing algorithm (for details see manual)
Definition at line 2326 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 1717 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 1729 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 1735 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 1741 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 368 of file BCIntegrate.h.
{ fAbsolutePrecision = absprecision; }
| void BCIntegrate::SetCubaIntegrationMethod | ( | BCIntegrate::BCCubaMethod | type | ) |
Set Cuba integration method
Definition at line 2003 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 382 of file BCIntegrate.h.
{ fCubaMaxEval = n; }
| void BCIntegrate::SetCubaMinEval | ( | int | n | ) | [inline] |
Set minimum number of evaluations in Cuba integration
Definition at line 377 of file BCIntegrate.h.
{ fCubaMinEval = n; }
| void BCIntegrate::SetCubaSuaveFlatness | ( | double | p | ) | [inline] |
Set flatness for Cuba Suave
Definition at line 407 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 402 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 397 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 392 of file BCIntegrate.h.
{ fCubaVegasNStart = n; }
| void BCIntegrate::SetCubaVerbosityLevel | ( | int | n | ) | [inline] |
Set verbosity level of Cuba integration
Definition at line 387 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 451 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 463 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 457 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 469 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 419 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 431 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 437 of file BCIntegrate.h.
{ fFitFunctionIndexY = index; }
| void BCIntegrate::SetFitFunctionIndices | ( | int | indexx, | |
| int | indexy | |||
| ) | [inline] |
Sets indices of the x and y values in function fits.
- Parameters:
-
indexx Index of the x values indexy Index of the y values
Definition at line 444 of file BCIntegrate.h.
{ SetFitFunctionIndexX(indexx);
SetFitFunctionIndexY(indexy); }
| void BCIntegrate::SetFlagIgnorePrevOptimization | ( | bool | flag | ) | [inline] |
Flag whether or not to ignore result of previous mode finding
Definition at line 301 of file BCIntegrate.h.
{ fFlagIgnorePrevOptimization = flag; }
| void BCIntegrate::SetFlagWriteSAToFile | ( | bool | flag | ) | [inline] |
Definition at line 523 of file BCIntegrate.h.
{ fFlagWriteSAToFile = flag; }
| void BCIntegrate::SetIntegrationMethod | ( | BCIntegrate::BCIntegrationMethod | method | ) |
- Parameters:
-
method The integration method
Definition at line 532 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 324 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 501 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 495 of file BCIntegrate.h.
{ fMarkovChainNIterations = niterations;
fMarkovChainAutoN = false; }
| void BCIntegrate::SetMarkovChainStepSize | ( | double | stepsize | ) | [inline] |
Sets the step size for Markov chains
Definition at line 490 of file BCIntegrate.h.
{ fMarkovChainStepSize = stepsize; }
| void BCIntegrate::SetMarkovChainTree | ( | TTree * | tree | ) | [inline] |
Sets the ROOT tree containing the Markov chain
Definition at line 481 of file BCIntegrate.h.
{ fMarkovChainTree = tree; }
| void BCIntegrate::SetMinuitArlist | ( | double * | arglist | ) | [inline] |
pointer to list of doubles to be passed as arguments to Minuit
Definition at line 295 of file BCIntegrate.h.
{ fMinuitArglist[0] = arglist[0];
fMinuitArglist[1] = arglist[1]; }
| void BCIntegrate::SetMode | ( | std::vector< double > | mode | ) |
Sets mode
Definition at line 1935 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 | |||
| ) |
Set the number of bins for the marginalized distribution of a parameter.
- Parameters:
-
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 357 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 346 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 352 of file BCIntegrate.h.
{ fNSamplesPer2DBin = n; }
| void BCIntegrate::SetOptimizationMethod | ( | BCIntegrate::BCOptimizationMethod | method | ) | [inline] |
- Parameters:
-
method The mode finding method
Definition at line 329 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 335 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 363 of file BCIntegrate.h.
{ fRelativePrecision = relprecision; }
| void BCIntegrate::SetSASchedule | ( | BCIntegrate::BCSASchedule | schedule | ) | [inline] |
- Parameters:
-
method The Simulated Annealing schedule
Definition at line 341 of file BCIntegrate.h.
{ fSASchedule = schedule; }
| void BCIntegrate::SetSAT0 | ( | double | T0 | ) | [inline] |
- Parameters:
-
T0 new value for Simulated Annealing starting temperature.
Definition at line 515 of file BCIntegrate.h.
{ fSAT0 = T0; }
| void BCIntegrate::SetSATmin | ( | double | Tmin | ) | [inline] |
- Parameters:
-
Tmin new value for Simulated Annealing threshold temperature.
Definition at line 520 of file BCIntegrate.h.
{ fSATmin = Tmin; }
| void BCIntegrate::SetSATree | ( | TTree * | tree | ) | [inline] |
Sets the tree containing the Simulated Annealing chain.
Definition at line 528 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 315 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 425 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 557 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 474 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 1069 of file BCIntegrate.h.
std::vector<double> BCIntegrate::fBestFitParameterErrors [protected] |
Definition at line 930 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 929 of file BCIntegrate.h.
std::vector<double> BCIntegrate::fBestFitParametersMarginalized [protected] |
A vector of best fit parameters estimated from the marginalized probability
Definition at line 934 of file BCIntegrate.h.
Cuba integration method
Definition at line 1072 of file BCIntegrate.h.
int BCIntegrate::fCubaMaxEval [private] |
Maximum number of evaluations in Cuba integration
Definition at line 1078 of file BCIntegrate.h.
int BCIntegrate::fCubaMinEval [private] |
Minimum number of evaluations in Cuba integration
Definition at line 1075 of file BCIntegrate.h.
double BCIntegrate::fCubaSuaveFlatness [private] |
Flatness for Cuba Suave
Definition at line 1093 of file BCIntegrate.h.
int BCIntegrate::fCubaSuaveNNew [private] |
Number of new integrand evaluations in each subdivision for Cuba Suave
Definition at line 1090 of file BCIntegrate.h.
int BCIntegrate::fCubaVegasNIncrease [private] |
Increase in number of evaluations per iteration for Cuba Vegas
Definition at line 1087 of file BCIntegrate.h.
int BCIntegrate::fCubaVegasNStart [private] |
Initial number of evaluations per iteration for Cuba Vegas
Definition at line 1084 of file BCIntegrate.h.
int BCIntegrate::fCubaVerbosity [private] |
Verbosity level of Cuba integration
Definition at line 1081 of file BCIntegrate.h.
std::vector<bool> BCIntegrate::fDataFixedValues [protected] |
Definition at line 924 of file BCIntegrate.h.
BCDataPoint* BCIntegrate::fDataPointLowerBoundaries [protected] |
data point containing the lower boundaries of possible data values
Definition at line 918 of file BCIntegrate.h.
BCDataPoint* BCIntegrate::fDataPointUpperBoundaries [protected] |
data point containing the upper boundaries of possible data values
Definition at line 922 of file BCIntegrate.h.
double BCIntegrate::fError [private] |
The uncertainty in the most recent Monte Carlo integration
Definition at line 1097 of file BCIntegrate.h.
bool BCIntegrate::fErrorBandContinuous [protected] |
A flag for single point evaluation of the error "band"
Definition at line 955 of file BCIntegrate.h.
int BCIntegrate::fErrorBandNbinsX [protected] |
Number of X bins of the error band histogram
Definition at line 964 of file BCIntegrate.h.
int BCIntegrate::fErrorBandNbinsY [protected] |
Nnumber of Y bins of the error band histogram
Definition at line 968 of file BCIntegrate.h.
std::vector<double> BCIntegrate::fErrorBandX [protected] |
Definition at line 956 of file BCIntegrate.h.
TH2D* BCIntegrate::fErrorBandXY [protected] |
The error band histogram
Definition at line 960 of file BCIntegrate.h.
bool BCIntegrate::fFillErrorBand [protected] |
Flag whether or not to fill the error band
Definition at line 946 of file BCIntegrate.h.
int BCIntegrate::fFitFunctionIndexX [protected] |
The indices for function fits
Definition at line 950 of file BCIntegrate.h.
int BCIntegrate::fFitFunctionIndexY [protected] |
Definition at line 951 of file BCIntegrate.h.
double BCIntegrate::fFlagIgnorePrevOptimization [protected] |
Flag for ignoring older results of minimization
Definition at line 979 of file BCIntegrate.h.
bool BCIntegrate::fFlagWriteMarkovChain [protected] |
Flag for writing Markov chain to file
Definition at line 983 of file BCIntegrate.h.
bool BCIntegrate::fFlagWriteSAToFile [protected] |
Flag deciding whether to write SA to file or not.
Definition at line 1007 of file BCIntegrate.h.
std::vector<TH1D *> BCIntegrate::fHProb1D [protected] |
Vector of TH1D histograms for marginalized probability distributions
Definition at line 938 of file BCIntegrate.h.
std::vector<TH2D *> BCIntegrate::fHProb2D [protected] |
Vector of TH2D histograms for marginalized probability distributions
Definition at line 942 of file BCIntegrate.h.
Current integration method
Definition at line 1038 of file BCIntegrate.h.
Current marginalization method
Definition at line 1042 of file BCIntegrate.h.
bool BCIntegrate::fMarkovChainAutoN [protected] |
Definition at line 914 of file BCIntegrate.h.
int BCIntegrate::fMarkovChainNIterations [protected] |
Definition at line 912 of file BCIntegrate.h.
double BCIntegrate::fMarkovChainStepSize [protected] |
Step size in the Markov chain relative to min and max
Definition at line 910 of file BCIntegrate.h.
TTree* BCIntegrate::fMarkovChainTree [protected] |
ROOT tree containing the Markov chain
Definition at line 987 of file BCIntegrate.h.
double BCIntegrate::fMarkovChainValue [private] |
A double containing the log likelihood value at the point fXmetro1
Definition at line 1114 of file BCIntegrate.h.
double* BCIntegrate::fMax [private] |
Array containing the upper boundaries of the variables to integrate over.
Definition at line 1026 of file BCIntegrate.h.
int BCIntegrate::fMCMCIteration [protected] |
Iteration of the MCMC
Definition at line 991 of file BCIntegrate.h.
double* BCIntegrate::fMin [private] |
Array containing the lower boundaries of the variables to integrate over.
Definition at line 1022 of file BCIntegrate.h.
TMinuit* BCIntegrate::fMinuit [protected] |
Minuit
Definition at line 972 of file BCIntegrate.h.
double BCIntegrate::fMinuitArglist[2] [protected] |
Definition at line 974 of file BCIntegrate.h.
int BCIntegrate::fMinuitErrorFlag [protected] |
Definition at line 975 of file BCIntegrate.h.
int BCIntegrate::fNacceptedMCMC [private] |
Definition at line 1102 of file BCIntegrate.h.
int BCIntegrate::fNbins [protected] |
Number of bins per dimension for the marginalized distributions
Definition at line 901 of file BCIntegrate.h.
int BCIntegrate::fNIterations [private] |
Number of iterations in the most recent Monte Carlo integation
Definition at line 1063 of file BCIntegrate.h.
int BCIntegrate::fNIterationsMax [private] |
Maximum number of iterations
Definition at line 1059 of file BCIntegrate.h.
int BCIntegrate::fNiterPerDimension [private] |
Number of iteration per dimension for Monte Carlo integration.
Definition at line 1034 of file BCIntegrate.h.
int BCIntegrate::fNmetro [private] |
The number of iterations in the Metropolis integration
Definition at line 1101 of file BCIntegrate.h.
int BCIntegrate::fNSamplesPer2DBin [protected] |
Number of samples per 2D bin per variable in the Metropolis marginalization.
Definition at line 906 of file BCIntegrate.h.
int BCIntegrate::fNvar [protected] |
Number of variables to integrate over.
Definition at line 897 of file BCIntegrate.h.
Current mode finding method
Definition at line 1046 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 1051 of file BCIntegrate.h.
double BCIntegrate::fRelativePrecision [private] |
Requested relative precision of the integation
Definition at line 1066 of file BCIntegrate.h.
double BCIntegrate::fSALogProb [protected] |
Definition at line 1011 of file BCIntegrate.h.
int BCIntegrate::fSANIterations [protected] |
Definition at line 1009 of file BCIntegrate.h.
Current Simulated Annealing schedule
Definition at line 1055 of file BCIntegrate.h.
double BCIntegrate::fSAT0 [protected] |
Starting temperature for Simulated Annealing
Definition at line 995 of file BCIntegrate.h.
double BCIntegrate::fSATemperature [protected] |
Definition at line 1010 of file BCIntegrate.h.
double BCIntegrate::fSATmin [protected] |
Minimal/Threshold temperature for Simulated Annealing
Definition at line 999 of file BCIntegrate.h.
std::vector<double> BCIntegrate::fSAx [protected] |
Definition at line 1012 of file BCIntegrate.h.
TTree* BCIntegrate::fTreeSA [protected] |
Tree for the Simulated Annealing
Definition at line 1003 of file BCIntegrate.h.
int* BCIntegrate::fVarlist [private] |
List of variables containing a flag whether to integrate over them or not.
Definition at line 1030 of file BCIntegrate.h.
BCParameterSet* BCIntegrate::fx [private] |
Set of parameters for the integration.
Definition at line 1018 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 1106 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 1110 of file BCIntegrate.h.
The documentation for this class was generated from the following files:
