BCModel Class Reference

The base class for all user-defined models. More...

#include <BCModel.h>

Inherits BCIntegrate.

Inherited by BCEfficiencyFitter, BCGoFTest, BCGraphFitter, and BCHistogramFitter.

Collaboration diagram for BCModel:

Collaboration graph
[legend]

List of all members.

Public Member Functions

Member functions (miscellaneous methods)


int AddParameter (BCParameter *parameter)
int AddParameter (const char *name, double lowerlimit, double upperlimit)
double APrioriProbability (std::vector< double > parameters)
BCH1DCalculatePValue (std::vector< double > par, bool flag_histogram=false)
int CheckParameters (std::vector< double > parameters)
double ConditionalProbabilityEntry (BCDataPoint *datapoint, std::vector< double > parameters)
virtual void CorrelateDataPointValues (std::vector< double > &x)
double Eval (std::vector< double > parameters)
double EvalSampling (std::vector< double > parameters)
void FindMode (std::vector< double > start=std::vector< double >(0))
void FixDataAxis (unsigned int index, bool fixed)
BCH2DGetMarginalized (const char *name1, const char *name2)
BCH2DGetMarginalized (BCParameter *parameter1, BCParameter *parameter2)
BCH1DGetMarginalized (const char *name)
BCH1DGetMarginalized (BCParameter *parameter)
double GetPValue ()
double GetPvalueFromChi2 (std::vector< double > par, int sigma_index)
double HessianMatrixElement (BCParameter *parameter1, BCParameter *parameter2, std::vector< double > point)
double Likelihood (std::vector< double > parameter)
virtual double LogAPrioriProbability (std::vector< double > parameters)
virtual double LogConditionalProbabilityEntry (BCDataPoint *datapoint, std::vector< double > parameters)
double LogEval (std::vector< double > parameters)
virtual double LogLikelihood (std::vector< double > parameter)
double LogProbability (std::vector< double > parameter)
double LogProbabilityNN (std::vector< double > parameter)
int MarginalizeAll ()
double Normalize ()
int PrintAllMarginalized (const char *file, unsigned int hdiv=1, unsigned int ndiv=1)
int PrintAllMarginalized1D (const char *filebase)
int PrintAllMarginalized2D (const char *filebase)
void PrintHessianMatrix (std::vector< double > parameters)
void PrintResults (const char *file)
void PrintShortFitSummary ()
void PrintSummary ()
double Probability (std::vector< double > parameter)
double ProbabilityNN (std::vector< double > parameter)
int ReadErrorBandFromFile (const char *file)
int ReadMarginalizedFromFile (const char *file)
int ReadMode (const char *file)
virtual double SamplingFunction (std::vector< double > parameters)
void SetGoFNChains (int n)
void SetGoFNIterationsMax (int n)
void SetGoFNIterationsRun (int n)
void WriteMode (const char *file)
Constructors and destructors


 BCModel (const char *name)
 BCModel ()
virtual ~BCModel ()
Member functions (get)


double GetBestFitParameter (unsigned int index)
double GetBestFitParameterMarginalized (unsigned int index)
std::vector< double > GetBestFitParameters ()
std::vector< double > GetBestFitParametersMarginalized ()
BCDataPointGetDataPoint (unsigned int index)
BCDataPointGetDataPointLowerBoundaries ()
double GetDataPointLowerBoundary (unsigned int index)
BCDataPointGetDataPointUpperBoundaries ()
double GetDataPointUpperBoundary (unsigned int index)
BCDataSetGetDataSet ()
std::vector< double > GetErrorBand (double level)
TGraph * GetErrorBandGraph (double level1, double level2)
TH2D * GetErrorBandXY ()
TH2D * GetErrorBandXY_yellow (double level=.68, int nsmooth=0)
TGraph * GetFitFunctionGraph (std::vector< double > parameters, double xmin, double xmax, int n=1000)
TGraph * GetFitFunctionGraph ()
TGraph * GetFitFunctionGraph (std::vector< double > parameters)
bool GetFixedDataAxis (unsigned int index)
bool GetFlagBoundaries ()
int GetIndex ()
double GetModelAPosterioriProbability ()
double GetModelAPrioriProbability ()
std::string GetName ()
int GetNDataPoints ()
unsigned int GetNDataPointsMaximum ()
unsigned int GetNDataPointsMinimum ()
double GetNormalization ()
unsigned int GetNParameters ()
BCParameterGetParameter (const char *name)
BCParameterGetParameter (int index)
BCParameterSetGetParameterSet ()
Member functions (set)


void SetDataBoundaries (unsigned int index, double lowerboundary, double upperboundary, bool fixed=false)
void SetDataSet (BCDataSet *dataset)
void SetErrorBandContinuous (bool flag)
void SetIndex (int index)
void SetModelAPosterioriProbability (double probability)
void SetModelAPrioriProbability (double probability)
void SetName (const char *name)
void SetNbins (const char *parname, int nbins)
void SetNDataPointsMaximum (unsigned int maximum)
void SetNDataPointsMinimum (unsigned int minimum)
void SetNormalization (double norm)
void SetParameterSet (BCParameterSet *parset)
void SetSingleDataPoint (BCDataSet *dataset, unsigned int index)
void SetSingleDataPoint (BCDataPoint *datapoint)

Protected Attributes

BCDataSetfDataSet
int fGoFNChains
int fGoFNIterationsMax
int fGoFNIterationsRun
int fIndex
bool flag_ConditionalProbabilityEntry
double fModelAPosteriori
double fModelAPriori
std::string fName
unsigned int fNDataPointsMaximum
unsigned int fNDataPointsMinimum
BCParameterSetfParameterSet
double fPValue

Private Member Functions

int CompareStrings (const char *string1, const char *string2)
BCDataPointVectorToDataPoint (std::vector< double > data)

Private Attributes

double fNormalization


Detailed Description

The base class for all user-defined models.

Author:
Daniel Kollar

Kevin Kröninger

Version:
1.0
Date:
08.2008 This class represents a model. It contains a container of parameters, their prior distributions and the conditional probabilities given those parameters. The methods which implement the prior and conditional probabilities have to be overloaded by the user in the user defined model class which will inherit from this class.

Definition at line 50 of file BCModel.h.


Constructor & Destructor Documentation

BCModel::BCModel (  ) 

The default constructor.

Definition at line 61 of file BCModel.cxx.

00061                  : BCIntegrate()
00062 {
00063    fNormalization = -1.0;
00064    fDataSet = 0;
00065    fParameterSet = new BCParameterSet();
00066 
00067    fIndex = -1;
00068    fPValue = -1;
00069 
00070    fName = "model";
00071    fDataPointUpperBoundaries = 0;
00072    fDataPointLowerBoundaries = 0;
00073 
00074    flag_ConditionalProbabilityEntry = true;
00075 
00076    fGoFNChains = 5;
00077    fGoFNIterationsMax = 100000;
00078    fGoFNIterationsRun = 2000;
00079 }

BCModel::BCModel ( const char *  name  ) 

A constructor.

Parameters:
name The name of the model

Definition at line 37 of file BCModel.cxx.

00037                                   : BCIntegrate()
00038 {
00039    fNormalization = -1.0;
00040    fDataSet = 0;
00041    fParameterSet = new BCParameterSet;
00042 
00043    fIndex = -1;
00044    fPValue = -1;
00045 
00046    fName = (char *) name;
00047    flag_ConditionalProbabilityEntry = true;
00048 
00049    fDataPointUpperBoundaries = 0;
00050    fDataPointLowerBoundaries = 0;
00051 
00052    fErrorBandXY = 0;
00053 
00054    fGoFNChains = 5;
00055    fGoFNIterationsMax = 100000;
00056    fGoFNIterationsRun = 2000;
00057 }

BCModel::~BCModel (  )  [virtual]

The default destructor.

Definition at line 83 of file BCModel.cxx.

00084 {
00085    delete fParameterSet;
00086 
00087    if (fDataPointLowerBoundaries)
00088       delete fDataPointLowerBoundaries;
00089 
00090    if (fDataPointUpperBoundaries)
00091       delete fDataPointUpperBoundaries;
00092 }


Member Function Documentation

int BCModel::AddParameter ( BCParameter parameter  ) 

Adds a parameter to the model.

Parameters:
parameter A model parameter
See also:
AddParameter(const char * name, double lowerlimit, double upperlimit);

Definition at line 436 of file BCModel.cxx.

00437 {
00438    // check if parameter set exists
00439    if (!fParameterSet)
00440    {
00441       BCLog::OutError("BCModel::AddParameter : Parameter set does not exist");
00442       return ERROR_PARAMETERSETDOESNOTEXIST;
00443    }
00444 
00445    // check if parameter with same name exists
00446    int flag_exists = 0;
00447    for (unsigned int i = 0; i < this -> GetNParameters(); i++)
00448       if (this -> CompareStrings(parameter -> GetName().data(), this -> GetParameter(i) -> GetName().data()) == 0)
00449          flag_exists = -1;
00450 
00451    if (flag_exists < 0)
00452    {
00453       BCLog::OutError(
00454             Form("BCModel::AddParameter : Parameter with name %s exists already. ", parameter -> GetName().data()));
00455       return ERROR_PARAMETEREXISTSALREADY;
00456    }
00457 
00458    // define index of new parameter
00459    unsigned int index = fParameterSet -> size();
00460    parameter -> SetIndex(index);
00461 
00462    // add parameter to parameter container
00463    fParameterSet -> push_back(parameter);
00464 
00465    // add parameters to integation methods
00466    this -> SetParameters(fParameterSet);
00467 
00468    return 0;
00469 }

int BCModel::AddParameter ( const char *  name,
double  lowerlimit,
double  upperlimit 
)

Adds a parameter to the parameter set

Parameters:
name The name of the parameter
lowerlimit The lower limit of the parameter values
upperlimit The upper limit of the parameter values
See also:
AddParameter(BCParameter* parameter);

Definition at line 422 of file BCModel.cxx.

00423 {
00424    // create new parameter
00425    BCParameter * parameter = new BCParameter(name, lowerlimit, upperlimit);
00426 
00427    int flag_ok = this -> AddParameter(parameter);
00428    if (flag_ok)
00429       delete parameter;
00430 
00431    return flag_ok;
00432 }

double BCModel::APrioriProbability ( std::vector< double >  parameters  )  [inline]

Returns the prior probability.

Parameters:
parameters A set of parameter values
Returns:
The prior probability p(parameters)
See also:
GetPrior(std::vector <double> parameters)

Definition at line 329 of file BCModel.h.

00330          { return exp( this->LogAPrioriProbability(parameters) ); };

BCH1D * BCModel::CalculatePValue ( std::vector< double >  par,
bool  flag_histogram = false 
)

Definition at line 1172 of file BCModel.cxx.

01173 {
01174    BCH1D * hist = 0;
01175 
01176    // print log
01177    BCLog::OutSummary("Do goodness-of-fit-test");
01178 
01179    // create model test
01180    BCGoFTest * goftest = new BCGoFTest("modeltest");
01181 
01182    // set this model as the model to be tested
01183    goftest -> SetTestModel(this);
01184 
01185    // set the point in parameter space which is tested an initialize
01186    // the model testing
01187    if (!goftest -> SetTestPoint(par))
01188       return 0;
01189 
01190    // set parameters of the MCMC for the GoFTest
01191    goftest -> MCMCSetNChains(fGoFNChains);
01192    goftest -> MCMCSetNIterationsMax(fGoFNIterationsMax);
01193    goftest -> MCMCSetNIterationsRun(fGoFNIterationsRun);
01194 
01195    // get p-value
01196    fPValue = goftest -> GetCalculatedPValue(flag_histogram);
01197 
01198    // get histogram
01199    if (flag_histogram)
01200    {
01201       hist = new BCH1D();
01202       hist -> SetHistogram(goftest -> GetHistogramLogProb());
01203    }
01204 
01205    // delete model test
01206    delete goftest;
01207 
01208    // return histogram
01209    return hist;
01210 }

int BCModel::CheckParameters ( std::vector< double >  parameters  ) 

Checks if a set of parameters values is within the given range.

Parameters:
parameters A set of parameter values
Returns:
Error code (0: OK, -1 length of parameters not correct, -2 values not within range)

Definition at line 565 of file BCModel.cxx.

00566 {
00567    // check if vectors are of equal size
00568    if (!parameters.size() == fParameterSet -> size())
00569       return  ERROR_INVALIDNUMBEROFPARAMETERS;
00570 
00571    // check if parameters are within limits
00572    for (unsigned int i = 0; i < fParameterSet -> size(); i++)
00573    {
00574       BCParameter * modelparameter = fParameterSet -> at(i);
00575 
00576       if (modelparameter -> GetLowerLimit() > parameters.at(i) ||
00577             modelparameter -> GetUpperLimit() < parameters.at(i))
00578       {
00579          BCLog::OutError(
00580                 Form("BCModel::CheckParameters : Parameter %s not within limits.", fParameterSet -> at(i) -> GetName().data()));
00581          return ERROR_PARAMETERNOTWITHINRANGE;
00582       }
00583    }
00584 
00585    return 0;
00586 }

int BCModel::CompareStrings ( const char *  string1,
const char *  string2 
) [private]

Compares to strings

Definition at line 1577 of file BCModel.cxx.

01578 {
01579    int flag_same = 0;
01580 
01581    if (strlen(string1) != strlen(string2))
01582       return -1;
01583 
01584    for (int i = 0; i < int(strlen(string1)); i++)
01585       if (string1[i] != string2[i])
01586          flag_same = -1;
01587 
01588    return flag_same;
01589 }

double BCModel::ConditionalProbabilityEntry ( BCDataPoint datapoint,
std::vector< double >  parameters 
) [inline]

Returns a conditional probability. Method needs to be overloaded by the user.

Parameters:
datapoint A data point
parameters A set of parameter values
Returns:
The conditional probability p(datapoint|parameters)
See also:
GetConditionalEntry(BCDataPoint* datapoint, std::vector <double> parameters)

Definition at line 389 of file BCModel.h.

00390          { return exp( this->LogConditionalProbabilityEntry(datapoint, parameters) ); };

void BCModel::CorrelateDataPointValues ( std::vector< double > &  x  )  [virtual]

Constrains a data point

Parameters:
x A vector of double

Definition at line 1214 of file BCModel.cxx.

01215 {
01216    // ...
01217 }

double BCModel::Eval ( std::vector< double >  parameters  )  [inline, virtual]

Overloaded function to evaluate integral.

Reimplemented from BCIntegrate.

Definition at line 411 of file BCModel.h.

00412          { return exp( this->LogEval(parameters) ); };

double BCModel::EvalSampling ( std::vector< double >  parameters  )  [virtual]

Overloaded function to evaluate integral.

Reimplemented from BCIntegrate.

Definition at line 524 of file BCModel.cxx.

00525 {
00526    return this -> SamplingFunction(parameters);
00527 }

void BCModel::FindMode ( std::vector< double >  start = std::vector<double>(0)  ) 

Do the mode finding using a method set via SetOptimizationMethod. Default is Minuit. The mode can be extracted using the GetBestFitParameters() method.

A starting point for the mode finding can be specified for Minuit. If not specified, Minuit default will be used (center of the parameter space).

If running mode finding after the MCMC it is a good idea to specify the mode obtained from MCMC as a starting point for the Minuit minimization. MCMC will find the location of the global mode and Minuit will converge to the mode precisely. The commands are:

			model -> MarginalizeAll();
			model -> FindMode( model -> GetBestFitParameters() );
			
startinf point of Minuit minimization

Definition at line 590 of file BCModel.cxx.

00591 {
00592 // this implementation is CLEARLY not good we have to work on this.
00593 
00594    BCLog::OutSummary(Form("Model \'%s\': Finding mode", this -> GetName().data()));
00595 
00596    // synchronize parameters in BCIntegrate
00597    this -> SetParameters(fParameterSet);
00598 
00599    switch(this -> GetOptimizationMethod())
00600    {
00601       case BCIntegrate::kOptSimulatedAnnealing:
00602          BCLog::OutError("BCModel::FindMode : Simulated annaeling not yet implemented");
00603          return;
00604 
00605       case BCIntegrate::kOptMinuit:
00606          this -> FindModeMinuit(start);
00607          return;
00608 
00609       case BCIntegrate::kOptMetropolis:
00610          this -> MarginalizeAll();
00611          return;
00612       }
00613 
00614    BCLog::OutError(
00615       Form("BCModel::FindMode : Invalid mode finding method: %d",
00616          this->GetOptimizationMethod()));
00617 
00618    return;
00619 }

void BCModel::FixDataAxis ( unsigned int  index,
bool  fixed 
)

Definition at line 1268 of file BCModel.cxx.

01269 {
01270    // check if index is within range
01271    if (index < 0 || index > fDataSet -> GetDataPoint(0) -> GetNValues())
01272    {
01273       BCLog::OutWarning("BCModel::FixDataAxis. Index out of range.");
01274       return;
01275    }
01276 
01277    if (fDataFixedValues.size() == 0)
01278       fDataFixedValues.assign(fDataSet -> GetDataPoint(0) -> GetNValues(), false);
01279 
01280    fDataFixedValues[index] = fixed;
01281 }

double BCModel::GetBestFitParameter ( unsigned int  index  )  [inline]

Returns the value of a particular parameter (defined by index) at the global mode of the posterior pdf.

Parameters:
index The index of the parameter.
Returns:
The best fit parameter.

Definition at line 177 of file BCModel.h.

00178          { return fBestFitParameters[index]; };

double BCModel::GetBestFitParameterMarginalized ( unsigned int  index  )  [inline]

Returns the value of a particular parameter (defined by index) at the modes of the marginalized posterior pdfs.

Parameters:
index The index of the parameter.
Returns:
The best fit parameter

Definition at line 192 of file BCModel.h.

00193          { return fBestFitParametersMarginalized[index]; };

std::vector<double> BCModel::GetBestFitParameters (  )  [inline]

Returns the set of values of the parameters at the global mode of the posterior pdf.

Returns:
The best fit parameters

Definition at line 184 of file BCModel.h.

00185          { return fBestFitParameters; };

std::vector<double> BCModel::GetBestFitParametersMarginalized (  )  [inline]

Returns the set of values of the parameters at the modes of the marginalized posterior pdfs.

Returns:
The best fit parameters.

Definition at line 199 of file BCModel.h.

00200          { return fBestFitParametersMarginalized; };

BCDataPoint * BCModel::GetDataPoint ( unsigned int  index  ) 

Parameters:
index The index of the data point.
Returns:
The data point in the current data set at index

Definition at line 112 of file BCModel.cxx.

00113 {
00114    if (fDataSet)
00115       return fDataSet -> GetDataPoint(index);
00116 
00117    BCLog::OutWarning("BCModel::GetDataPoint. No data set defined.");
00118    return 0;
00119 }

BCDataPoint* BCModel::GetDataPointLowerBoundaries (  )  [inline]

Returns:
The lower boundaries of possible data values.

Definition at line 108 of file BCModel.h.

00109          { return fDataPointLowerBoundaries; };

double BCModel::GetDataPointLowerBoundary ( unsigned int  index  )  [inline]

Parameters:
index The index of the variable.
Returns:
The lower boundary of possible data values for a particular variable.

Definition at line 119 of file BCModel.h.

00120          { return fDataPointLowerBoundaries -> GetValue(index); };

BCDataPoint* BCModel::GetDataPointUpperBoundaries (  )  [inline]

Returns:
The upper boundaries of possible data values.

Definition at line 113 of file BCModel.h.

00114          { return fDataPointUpperBoundaries; };

double BCModel::GetDataPointUpperBoundary ( unsigned int  index  )  [inline]

Parameters:
index The index of the variable.
Returns:
The upper boundary of possible data values for a particular variable.

Definition at line 125 of file BCModel.h.

00126          { return fDataPointUpperBoundaries -> GetValue(index); };

BCDataSet* BCModel::GetDataSet (  )  [inline]

Returns:
The data set.

Definition at line 103 of file BCModel.h.

00104          { return fDataSet; };

std::vector< double > BCModel::GetErrorBand ( double  level  ) 

Returns a vector of y-values at a certain probability level.

Parameters:
level The level of probability
Returns:
The vector of y-values

Definition at line 180 of file BCModel.cxx.

00181 {
00182    std::vector <double> errorband;
00183 
00184    if (!fErrorBandXY)
00185       return errorband;
00186 
00187    int nx = fErrorBandXY -> GetNbinsX();
00188    errorband.assign(nx, 0.0);
00189 
00190    // loop over x and y bins
00191    for (int ix = 1; ix <= nx; ix++)
00192    {
00193       TH1D * temphist = fErrorBandXY -> ProjectionY("temphist", ix, ix);
00194 
00195       int nprobSum = 1;
00196       double q[1];
00197       double probSum[1];
00198       probSum[0] = level;
00199 
00200       temphist -> GetQuantiles(nprobSum, q, probSum);
00201 
00202       errorband[ix-1] = q[0];
00203    }
00204 
00205    return errorband;
00206 }

TGraph * BCModel::GetErrorBandGraph ( double  level1,
double  level2 
)

Definition at line 210 of file BCModel.cxx.

00211 {
00212    if (!fErrorBandXY)
00213       return 0;
00214 
00215    // define new graph
00216    int nx = fErrorBandXY -> GetNbinsX();
00217 
00218    TGraph * graph = new TGraph(2 * nx);
00219    graph -> SetFillStyle(1001);
00220    graph -> SetFillColor(kYellow);
00221 
00222    // get error bands
00223    std::vector <double> ymin = this -> GetErrorBand(level1);
00224    std::vector <double> ymax = this -> GetErrorBand(level2);
00225 
00226    for (int i = 0; i < nx; i++)
00227    {
00228       graph -> SetPoint(i,      fErrorBandXY -> GetXaxis() -> GetBinCenter(i + 1), ymin.at(i));
00229       graph -> SetPoint(nx + i, fErrorBandXY -> GetXaxis() -> GetBinCenter(nx - i), ymax.at(nx - i - 1));
00230    }
00231 
00232    return graph;
00233 }

TH2D* BCModel::GetErrorBandXY (  )  [inline]

Returns:
The 2-d histogram of the error band.

Definition at line 204 of file BCModel.h.

00205          { return fErrorBandXY; };

TH2D * BCModel::GetErrorBandXY_yellow ( double  level = .68,
int  nsmooth = 0 
)

Definition at line 237 of file BCModel.cxx.

00238 {
00239    if (!fErrorBandXY)
00240       return 0;
00241 
00242    int nx = fErrorBandXY -> GetNbinsX();
00243    int ny = fErrorBandXY -> GetNbinsY();
00244 
00245    // copy existing histogram
00246    TH2D * hist_tempxy = (TH2D*) fErrorBandXY -> Clone(TString::Format("%s_sub_%f.2",fErrorBandXY->GetName(),level));
00247    hist_tempxy -> Reset();
00248    hist_tempxy -> SetFillColor(kYellow);
00249 
00250    // loop over x bins
00251    for (int ix = 1; ix < nx; ix++)
00252    {
00253       BCH1D * hist_temp = new BCH1D();
00254 
00255       TH1D * hproj = fErrorBandXY -> ProjectionY("temphist", ix, ix);
00256       if(nsmooth>0)
00257          hproj->Smooth(nsmooth);
00258 
00259       hist_temp -> SetHistogram(hproj);
00260 
00261       TH1D * hist_temp_yellow = hist_temp -> GetSmallestIntervalHistogram(level);
00262 
00263       for (int iy = 1; iy <= ny; ++iy)
00264          hist_tempxy -> SetBinContent(ix, iy, hist_temp_yellow -> GetBinContent(iy));
00265 
00266       delete hist_temp_yellow;
00267       delete hist_temp;
00268    }
00269 
00270    return hist_tempxy;
00271 }

TGraph * BCModel::GetFitFunctionGraph ( std::vector< double >  parameters,
double  xmin,
double  xmax,
int  n = 1000 
)

Definition at line 302 of file BCModel.cxx.

00303 {
00304    // define new graph
00305    TGraph * graph = new TGraph(n+1);
00306 
00307    double dx = (xmax-xmin)/(double)n;
00308 
00309    // loop over x values
00310    for (int i = 0; i <= n; i++)
00311    {
00312       double x = (double)i*dx+xmin;
00313       std::vector <double> xvec;
00314       xvec.push_back(x);
00315       double y = this -> FitFunction(xvec, parameters);
00316 
00317       xvec.clear();
00318 
00319       graph -> SetPoint(i, x, y);
00320    }
00321 
00322    return graph;
00323 }

TGraph* BCModel::GetFitFunctionGraph (  )  [inline]

Definition at line 219 of file BCModel.h.

00220          { return this -> GetFitFunctionGraph(this -> GetBestFitParameters()); };

TGraph * BCModel::GetFitFunctionGraph ( std::vector< double >  parameters  ) 

Definition at line 275 of file BCModel.cxx.

00276 {
00277    if (!fErrorBandXY)
00278       return 0;
00279 
00280    // define new graph
00281    int nx = fErrorBandXY -> GetNbinsX();
00282    TGraph * graph = new TGraph(nx);
00283 
00284    // loop over x values
00285    for (int i = 0; i < nx; i++)
00286    {
00287       double x = fErrorBandXY -> GetXaxis() -> GetBinCenter(i + 1);
00288 
00289       std::vector <double> xvec;
00290       xvec.push_back(x);
00291       double y = this -> FitFunction(xvec, parameters);
00292       xvec.clear();
00293 
00294       graph -> SetPoint(i, x, y);
00295    }
00296 
00297    return graph;
00298 }

bool BCModel::GetFixedDataAxis ( unsigned int  index  ) 

Definition at line 1285 of file BCModel.cxx.

01286 {
01287    // check if index is within range
01288    if (index < 0 || index > fDataSet -> GetDataPoint(0) -> GetNValues())
01289    {
01290       BCLog::OutWarning("BCModel::GetFixedDataAxis. Index out of range.");
01291       return false;
01292    }
01293 
01294    return fDataFixedValues.at(index);
01295 }

bool BCModel::GetFlagBoundaries (  ) 

Definition at line 327 of file BCModel.cxx.

00328 {
00329    if (!fDataPointLowerBoundaries)
00330       return false;
00331 
00332    if (!fDataPointUpperBoundaries)
00333       return false;
00334 
00335    if (fDataPointLowerBoundaries -> GetNValues() != fDataSet -> GetDataPoint(0) -> GetNValues())
00336       return false;
00337 
00338    if (fDataPointUpperBoundaries -> GetNValues() != fDataSet -> GetDataPoint(0) -> GetNValues())
00339       return false;
00340 
00341    return true;
00342 }

int BCModel::GetIndex (  )  [inline]

Returns:
The index of the model.

Definition at line 83 of file BCModel.h.

00084          { return fIndex; };

BCH2D* BCModel::GetMarginalized ( const char *  name1,
const char *  name2 
) [inline]

Definition at line 493 of file BCModel.h.

00494          { return this -> GetMarginalized(this -> GetParameter(name1), this -> GetParameter(name2)); };

BCH2D * BCModel::GetMarginalized ( BCParameter parameter1,
BCParameter parameter2 
)

If MarginalizeAll method was used, the individual marginalized distributions with respect to otwo parameters can be retrieved using this method.

Parameters:
parameter1 First parameter
parameter2 Second parameter
Returns:
2D marginalized probability

Definition at line 1104 of file BCModel.cxx.

01105 {
01106 // if(fMCMCH2Marginalized.size()==0)
01107 // {
01108 //    BCLog::Out(BCLog::warning, BCLog::warning,
01109 //           "BCModel::GetMarginalized. MarginalizeAll() has to be run prior to this.");
01110 //    return 0;
01111 // }
01112 
01113    int index1 = parameter1 -> GetIndex();
01114    int index2 = parameter2 -> GetIndex();
01115 
01116    if (index1 == index2)
01117    {
01118       BCLog::OutError("BCModel::GetMarginalized : Provided parameters are identical. Distribution not available.");
01119       return 0;
01120    }
01121 
01122    if (index1 > index2)
01123    {
01124       int itmp = index1;
01125       index1 = index2;
01126       index2 = itmp;
01127    }
01128 
01129    // get histogram
01130    TH2D * hist = this -> MCMCGetH2Marginalized(index1, index2);
01131 
01132    if(hist==0)
01133       return 0;
01134 
01135    BCH2D * hprob = new BCH2D();
01136 
01137    // set axis labels
01138    hist -> SetName(Form("hist_%s_%s_%s", this -> GetName().data(), parameter1 -> GetName().data(), parameter2 -> GetName().data()));
01139    hist -> SetXTitle(Form("%s", parameter1 -> GetName().data()));
01140    hist -> SetYTitle(Form("%s", parameter2 -> GetName().data()));
01141    hist -> SetStats(kFALSE);
01142 
01143    double gmode[] = { fBestFitParameters.at(index1), fBestFitParameters.at(index2) };
01144    hprob->SetGlobalMode(gmode);
01145 
01146    // set histogram
01147    hprob -> SetHistogram(hist);
01148 
01149    return hprob;
01150 }

BCH1D* BCModel::GetMarginalized ( const char *  name  )  [inline]

Definition at line 482 of file BCModel.h.

00483          { return this -> GetMarginalized(this -> GetParameter(name)); };

BCH1D * BCModel::GetMarginalized ( BCParameter parameter  ) 

If MarginalizeAll method was used, the individual marginalized distributions with respect to one parameter can be retrieved using this method.

Parameters:
parameter Model parameter
Returns:
1D marginalized probability

Definition at line 746 of file BCModel.cxx.

00747 {
00748 // if(fMCMCH1Marginalized.size()==0)
00749 // {
00750 //    BCLog::Out(BCLog::warning, BCLog::warning,
00751 //          "BCModel::GetMarginalized. MarginalizeAll() has to be run prior to this.");
00752 //    return 0;
00753 // }
00754 
00755    int index = parameter -> GetIndex();
00756 
00757    // get histogram
00758    TH1D * hist = this -> MCMCGetH1Marginalized(index);
00759    if(!hist)
00760       return 0;
00761 
00762    BCH1D * hprob = new BCH1D();
00763 
00764    // set axis labels
00765    hist -> SetName(Form("hist_%s_%s", this -> GetName().data(), parameter -> GetName().data()));
00766    hist -> SetXTitle(parameter -> GetName().data());
00767    hist -> SetYTitle(Form("p(%s|data)", parameter -> GetName().data()));
00768    hist -> SetStats(kFALSE);
00769 
00770    // set histogram
00771    hprob -> SetHistogram(hist);
00772 
00773    // set best fit parameter
00774    double bestfit = hprob -> GetMode();
00775 
00776    if (fBestFitParametersMarginalized.size() == 0)
00777       for (unsigned int i = 0; i < this -> GetNParameters(); i++)
00778          fBestFitParametersMarginalized.push_back(0.0);
00779 
00780    fBestFitParametersMarginalized[index] = bestfit;
00781 
00782    hprob->SetGlobalMode(fBestFitParameters.at(index));
00783 
00784    return hprob;
00785 }

double BCModel::GetModelAPosterioriProbability (  )  [inline]

Returns:
The a posteriori probability.

Definition at line 93 of file BCModel.h.

00094          { return fModelAPosteriori; };

double BCModel::GetModelAPrioriProbability (  )  [inline]

Returns:
The a priori probability.

Definition at line 88 of file BCModel.h.

00089          { return fModelAPriori; };

std::string BCModel::GetName (  )  [inline]

Returns:
The name of the model.

Definition at line 78 of file BCModel.h.

00079          { return fName; };

int BCModel::GetNDataPoints (  ) 

Returns:
The number of data points in the current data set.

Definition at line 96 of file BCModel.cxx.

00097 {
00098    int npoints = 0;
00099    if (fDataSet)
00100       npoints = fDataSet -> GetNDataPoints();
00101    else
00102    {
00103       BCLog::OutWarning("BCModel::GetNDataPoints(). No data set defined.");
00104       return ERROR_NOEVENTS;
00105    }
00106 
00107    return npoints;
00108 }

unsigned int BCModel::GetNDataPointsMaximum (  )  [inline]

Returns:
The maximum number of data points.

Definition at line 149 of file BCModel.h.

00150          { return fNDataPointsMaximum; };

unsigned int BCModel::GetNDataPointsMinimum (  )  [inline]

Returns:
The minimum number of data points.

Definition at line 144 of file BCModel.h.

00145          { return fNDataPointsMinimum; };

double BCModel::GetNormalization (  )  [inline]

Returns:
The normalization factor of the probability

Definition at line 98 of file BCModel.h.

00099          { return fNormalization; };

unsigned int BCModel::GetNParameters (  )  [inline]

Returns:
The number of parameters of the model.

Definition at line 154 of file BCModel.h.

00155          { return fParameterSet ? fParameterSet -> size() : 0; };

BCParameter * BCModel::GetParameter ( const char *  name  ) 

Parameters:
name The name of the parameter in the parameter set.
Returns:
The parameter.

Definition at line 140 of file BCModel.cxx.

00141 {
00142    if (!fParameterSet)
00143       return 0;
00144 
00145    int index = -1;
00146    for (unsigned int i = 0; i < this->GetNParameters(); i++)
00147       if (name == this -> GetParameter(i) -> GetName())
00148          index = i;
00149 
00150    if (index<0)
00151    {
00152       BCLog::OutWarning(
00153             Form(
00154                   "BCModel::GetParameter : Model %s has no parameter named '%s'",
00155                   (this -> GetName()).data(), name
00156                   )
00157             );
00158       return 0;
00159    }
00160 
00161    return this->GetParameter(index);
00162 }

BCParameter * BCModel::GetParameter ( int  index  ) 

Parameters:
index The index of the parameter in the parameter set.
Returns:
The parameter.

Definition at line 123 of file BCModel.cxx.

00124 {
00125    if (!fParameterSet)
00126       return 0;
00127 
00128    if (index < 0 || index >= (int)this -> GetNParameters())
00129    {
00130       BCLog::OutWarning(
00131             Form("BCModel::GetParameter. Parameter index %d not within range.", index));
00132       return 0;
00133    }
00134 
00135    return fParameterSet -> at(index);
00136 }

BCParameterSet* BCModel::GetParameterSet (  )  [inline]

Returns:
parameter set

Definition at line 169 of file BCModel.h.

00170          { return fParameterSet; };

double BCModel::GetPValue (  )  [inline]

Definition at line 518 of file BCModel.h.

00519          { return fPValue; };

double BCModel::GetPvalueFromChi2 ( std::vector< double >  par,
int  sigma_index 
)

Calculate p-value from Chi2 distribution for Gaussian problems

Parameters:
par Parameter set for the calculation of the likelihood
sigma_index Index of the sigma/uncertainty for the data points (for data in format "x y erry" the index would be 2)

Definition at line 1154 of file BCModel.cxx.

01155 {
01156    double ll = this -> LogLikelihood(par);
01157    int n = this -> GetNDataPoints();
01158 
01159    double sum_sigma=0;
01160    for (int i=0;i<n;i++)
01161       sum_sigma += log(this -> GetDataPoint(i) -> GetValue(sigma_index));
01162 
01163    double chi2 = -2.*(ll + (double)n/2. * log(2.*M_PI) + sum_sigma);
01164 
01165    fPValue = TMath::Prob(chi2,n);
01166 
01167    return fPValue;
01168 }

double BCModel::HessianMatrixElement ( BCParameter parameter1,
BCParameter parameter2,
std::vector< double >  point 
)

Calculates the matrix element of the Hessian matrix

Parameters:
parameter1 The parameter for the first derivative
parameter2 The parameter for the first derivative
Returns:
The matrix element of the Hessian matrix

Definition at line 1221 of file BCModel.cxx.

01222 {
01223    // check number of entries in vector
01224    if (point.size() != this -> GetNParameters())
01225    {
01226       BCLog::OutWarning("BCModel::HessianMatrixElement. Invalid number of entries in the vector.");
01227       return -1;
01228    }
01229 
01230    // define steps
01231    double nsteps = 1e5;
01232    double dx1 = par1 -> GetRangeWidth() / nsteps;
01233    double dx2 = par2 -> GetRangeWidth() / nsteps;
01234 
01235    // define points at which to evaluate
01236    std::vector<double> xpp = point;
01237    std::vector<double> xpm = point;
01238    std::vector<double> xmp = point;
01239    std::vector<double> xmm = point;
01240 
01241    int idx1 = par1 -> GetIndex();
01242    int idx2 = par2 -> GetIndex();
01243 
01244    xpp[idx1] += dx1;
01245    xpp[idx2] += dx2;
01246 
01247    xpm[idx1] += dx1;
01248    xpm[idx2] -= dx2;
01249 
01250    xmp[idx1] -= dx1;
01251    xmp[idx2] += dx2;
01252 
01253    xmm[idx1] -= dx1;
01254    xmm[idx2] -= dx2;
01255 
01256    // calculate probability at these points
01257    double ppp = this -> Likelihood(xpp);
01258    double ppm = this -> Likelihood(xpm);
01259    double pmp = this -> Likelihood(xmp);
01260    double pmm = this -> Likelihood(xmm);
01261 
01262    // return derivative
01263    return (ppp + pmm - ppm - pmp) / (4. * dx1 * dx2);
01264 }

double BCModel::Likelihood ( std::vector< double >  parameter  )  [inline]

Returns the likelihood

Parameters:
parameters A set of parameter values
Returns:
The likelihood

Definition at line 345 of file BCModel.h.

00346          { return exp( this->LogLikelihood(parameter) ); };

virtual double BCModel::LogAPrioriProbability ( std::vector< double >  parameters  )  [inline, virtual]

Returns natural logarithm of the prior probability. Method needs to be overloaded by the user.

Parameters:
parameters A set of parameter values
Returns:
The prior probability p(parameters)
See also:
GetPrior(std::vector <double> parameters)

Reimplemented in BCEfficiencyFitter, BCGoFTest, BCGraphFitter, and BCHistogramFitter.

Definition at line 338 of file BCModel.h.

00339          { return 0.; };

virtual double BCModel::LogConditionalProbabilityEntry ( BCDataPoint datapoint,
std::vector< double >  parameters 
) [inline, virtual]

Returns a natural logarithm of conditional probability. Method needs to be overloaded by the user.

Parameters:
datapoint A data point
parameters A set of parameter values
Returns:
The conditional probability p(datapoint|parameters)
See also:
GetConditionalEntry(BCDataPoint* datapoint, std::vector <double> parameters)

Definition at line 399 of file BCModel.h.

00400          { flag_ConditionalProbabilityEntry = false; return 0.; };

double BCModel::LogEval ( std::vector< double >  parameters  )  [virtual]

Overloaded function to evaluate integral.

Reimplemented from BCIntegrate.

Definition at line 517 of file BCModel.cxx.

00518 {
00519    return this -> LogProbabilityNN(parameters);
00520 }

double BCModel::LogLikelihood ( std::vector< double >  parameter  )  [virtual]

Calculates natural logarithm of the likelihood. Method needs to be overloaded by the user.

Parameters:
parameters A set of parameter values
Returns:
Natural logarithm of the likelihood

Reimplemented in BCEfficiencyFitter, BCGoFTest, BCGraphFitter, and BCHistogramFitter.

Definition at line 501 of file BCModel.cxx.

00502 {
00503    double logprob = 0.;
00504 
00505    // add log of conditional probabilities event-by-event
00506    for (unsigned int i=0;i<fDataSet -> GetNDataPoints();i++)
00507    {
00508       BCDataPoint * datapoint = this -> GetDataPoint(i);
00509       logprob += this -> LogConditionalProbabilityEntry(datapoint, parameters);
00510    }
00511 
00512    return logprob;
00513 }

double BCModel::LogProbability ( std::vector< double >  parameter  ) 

Returns natural logarithm of the a posteriori probability given a set of parameter values

Parameters:
parameters A set of parameter values
Returns:
The a posteriori probability

Definition at line 486 of file BCModel.cxx.

00487 {
00488    // check if normalized
00489    if (fNormalization<0. || fNormalization==0.)
00490    {
00491       BCLog::Out(BCLog::warning, BCLog::warning,
00492              "BCModel::LogProbability. Normalization not available or zero.");
00493       return 0.;
00494    }
00495 
00496    return this -> LogProbabilityNN(parameters) - log(fNormalization);
00497 }

double BCModel::LogProbabilityNN ( std::vector< double >  parameter  ) 

Returns the natural logarithm of likelihood times prior probability given a set of parameter values

Parameters:
parameters A set of parameter values
Returns:
The likelihood times prior probability

Definition at line 473 of file BCModel.cxx.

00474 {
00475    // add log of conditional probability
00476    double logprob = this -> LogLikelihood(parameters);
00477 
00478    // add log of prior probability
00479    logprob += this -> LogAPrioriProbability(parameters);
00480 
00481    return logprob;
00482 }

int BCModel::MarginalizeAll (  ) 

Marginalize all probabilities wrt. single parameters and all combinations of two parameters. The individual distributions can be retrieved using the GetMarginalized method.

Returns:
Total number of marginalized distributions

Definition at line 704 of file BCModel.cxx.

00705 {
00706    BCLog::OutSummary(Form("Running MCMC for model \'%s\'",this->GetName().data()));
00707 
00708    // prepare function fitting
00709    double dx = 0.0;
00710    double dy = 0.0;
00711 
00712    if (fFitFunctionIndexX >= 0)
00713    {
00714       dx = (fDataPointUpperBoundaries -> GetValue(fFitFunctionIndexX) - fDataPointLowerBoundaries -> GetValue(fFitFunctionIndexX))
00715             / (double)fErrorBandNbinsX;
00716 
00717       dy = (fDataPointUpperBoundaries -> GetValue(fFitFunctionIndexY) - fDataPointLowerBoundaries -> GetValue(fFitFunctionIndexY))
00718             / (double)fErrorBandNbinsY;
00719 
00720       fErrorBandXY = new TH2D(
00721             TString::Format("errorbandxy_%d",BCLog::GetHIndex()), "",
00722             fErrorBandNbinsX,
00723             fDataPointLowerBoundaries -> GetValue(fFitFunctionIndexX) - .5 * dx,
00724             fDataPointUpperBoundaries -> GetValue(fFitFunctionIndexX) + .5 * dx,
00725             fErrorBandNbinsY,
00726             fDataPointLowerBoundaries -> GetValue(fFitFunctionIndexY) - .5 * dy,
00727             fDataPointUpperBoundaries -> GetValue(fFitFunctionIndexY) + .5 * dy);
00728       fErrorBandXY -> SetStats(kFALSE);
00729 
00730       for (int ix = 1; ix <= fErrorBandNbinsX; ++ix)
00731          for (int iy = 1; iy <= fErrorBandNbinsX; ++iy)
00732             fErrorBandXY -> SetBinContent(ix, iy, 0.0);
00733    }
00734 
00735    this -> MCMCMetropolis();
00736    this -> FindModeMCMC();
00737 
00738    // this -> PrintResults(Form("%s.txt", this -> GetName().data()));
00739 
00740    return 1;
00741 }

double BCModel::Normalize (  ) 

Integrates over the un-normalized probability and updates fNormalization.

Definition at line 541 of file BCModel.cxx.

00542 {
00543    BCLog::OutSummary(Form("Model \'%s\': Normalizing probability",this->GetName().data()));
00544 
00545    unsigned int n = this -> GetNvar();
00546 
00547    // initialize BCIntegrate if not done already
00548    if (n == 0)
00549    {
00550       this->SetParameters(fParameterSet);
00551       n = this->GetNvar();
00552    }
00553 
00554    // integrate and get best fit parameters
00555    // maybe we have to remove the mode finding from here in the future
00556    fNormalization = this -> Integrate();
00557 
00558    BCLog::OutDetail(Form(" --> Normalization factor : %.6g", fNormalization));
00559 
00560    return fNormalization;
00561 }

int BCModel::PrintAllMarginalized ( const char *  file,
unsigned int  hdiv = 1,
unsigned int  ndiv = 1 
)

Definition at line 941 of file BCModel.cxx.

00942 {
00943    if(fMCMCH1Marginalized.size()==0 || (fMCMCH2Marginalized.size()==0 && this -> GetNParameters() > 1))
00944    {
00945       BCLog::OutError("BCModel::PrintAllMarginalized : Marginalized distributions not available.");
00946       return 0;
00947    }
00948 
00949    // if there's only one parameter, we just want to call Print()
00950    if (fMCMCH1Marginalized.size()==1 && fMCMCH2Marginalized.size()==0)
00951    {
00952       BCParameter * a = this->GetParameter(0);
00953       this -> GetMarginalized(a) -> Print(file);
00954       return 1;
00955    }
00956 
00957    int c_width=600; // default canvas width
00958    int c_height=600; // default canvas height
00959 
00960    int type = 112; // landscape
00961 
00962    if (hdiv > vdiv)
00963    {
00964       if(hdiv>3)
00965       {
00966          c_width=1000;
00967          c_height=700;
00968       }
00969       else
00970       {
00971          c_width=800;
00972          c_height=600;
00973       }
00974    }
00975    else if(hdiv < vdiv)
00976    {
00977       if(hdiv>3)
00978       {
00979          c_height=1000;
00980          c_width=700;
00981       }
00982       else
00983       {
00984          c_height=800;
00985          c_width=600;
00986       }
00987       type=111;
00988    }
00989 
00990    // get number of parameters of the model and calculate number of plots
00991    int npar = this -> GetNParameters();
00992    int nplots2d = npar * (npar-1)/2;
00993    int nplots = npar + nplots2d;
00994 
00995    // give out warning if too many plots
00996    BCLog::OutSummary(Form(
00997          "Printing all marginalized distributions (%d x 1D + %d x 2D = %d) into file %s",
00998          npar,nplots2d,nplots,file));
00999    if(nplots>100)
01000       BCLog::OutDetail("This can take a while...");
01001 
01002    // setup the canvas and postscript file
01003    TCanvas * c = new TCanvas( "c","canvas",c_width,c_height);
01004 
01005    TPostScript * ps = new TPostScript(file,type);
01006 
01007    if(type==112)
01008       ps->Range(24,16);
01009    else
01010       ps->Range(16,24);
01011 
01012    // draw all 1D distributions
01013    ps->NewPage();
01014    c->cd();
01015    c->Clear();
01016    c->Divide(hdiv,vdiv);
01017 
01018    int n=0;
01019    for(int i=0;i<npar;i++)
01020    {
01021       // if current page is full, swith to new page
01022       if(i!=0 && i%(hdiv*vdiv)==0)
01023       {
01024          c->Update();
01025          ps->NewPage();
01026          c->cd();
01027          c->Clear();
01028          c->Divide(hdiv,vdiv);
01029       }
01030 
01031       // go to next pad
01032       c->cd(i%(hdiv*vdiv)+1);
01033 
01034       BCParameter * a = this->GetParameter(i);
01035       this -> GetMarginalized(a) -> Draw();
01036       n++;
01037 
01038       if(n%100==0)
01039          BCLog::OutDetail(Form(" --> %d plots done",n));
01040    }
01041 
01042    c->Update();
01043 
01044    // draw all the 2D distributions
01045    ps->NewPage();
01046    c->cd();
01047    c->Clear();
01048    c->Divide(hdiv,vdiv);
01049 
01050    int k=0;
01051    for(int i=0;i<npar-1;i++)
01052    {
01053       for(int j=i+1;j<npar;j++)
01054       {
01055          // if current page is full, switch to new page
01056          if(k!=0 && k%(hdiv*vdiv)==0)
01057          {
01058             c->Update();
01059             ps->NewPage();
01060             c->cd();
01061             c->Clear();
01062             c->Divide(hdiv,vdiv);
01063          }
01064 
01065          // go to next pad
01066          c->cd(k%(hdiv*vdiv)+1);
01067 
01068          BCParameter * a = this->GetParameter(i);
01069          BCParameter * b = this->GetParameter(j);
01070 
01071          double meana = (a -> GetLowerLimit() + a -> GetUpperLimit()) / 2.0;
01072          double deltaa = (a -> GetUpperLimit() - a -> GetLowerLimit());
01073          if (deltaa <= 1e-7 * meana)
01074             continue;
01075 
01076          double meanb = (b -> GetLowerLimit() + b -> GetUpperLimit()) / 2.0;
01077          double deltab = (b -> GetUpperLimit() - b -> GetLowerLimit());
01078          if (deltab <= 1e-7 * meanb)
01079             continue;
01080 
01081          this -> GetMarginalized(a,b) -> Draw(52);
01082          k++;
01083 
01084          if((n+k)%100==0)
01085             BCLog::OutDetail(Form(" --> %d plots done",n+k));
01086       }
01087    }
01088 
01089    if( (n+k)>100 && (n+k)%100 != 0 )
01090       BCLog::OutDetail(Form(" --> %d plots done",n+k));
01091 
01092    c->Update();
01093    ps->Close();
01094 
01095    delete c;
01096    delete ps;
01097 
01098    // return total number of drawn histograms
01099    return n+k;
01100 }

int BCModel::PrintAllMarginalized1D ( const char *  filebase  ) 

Definition at line 884 of file BCModel.cxx.

00885 {
00886    if(fMCMCH1Marginalized.size()==0)
00887    {
00888       BCLog::OutError("BCModel::PrintAllMarginalized : Marginalized distributions not available.");
00889       return 0;
00890    }
00891 
00892    int n=this->GetNParameters();
00893    for(int i=0;i<n;i++)
00894    {
00895       BCParameter * a = this->GetParameter(i);
00896       this -> GetMarginalized(a) -> Print(Form("%s_1D_%s.ps",filebase,a->GetName().data()));
00897    }
00898 
00899    return n;
00900 }

int BCModel::PrintAllMarginalized2D ( const char *  filebase  ) 

Definition at line 904 of file BCModel.cxx.

00905 {
00906    if(fMCMCH2Marginalized.size()==0)
00907    {
00908       BCLog::OutError("BCModel::PrintAllMarginalized : Marginalized distributions not available.");
00909       return 0;
00910    }
00911 
00912    int k=0;
00913    int n=this->GetNParameters();
00914    for(int i=0;i<n-1;i++)
00915    {
00916       for(int j=i+1;j<n;j++)
00917       {
00918          BCParameter * a = this->GetParameter(i);
00919          BCParameter * b = this->GetParameter(j);
00920 
00921          double meana = (a -> GetLowerLimit() + a -> GetUpperLimit()) / 2.0;
00922          double deltaa = (a -> GetUpperLimit() - a -> GetLowerLimit());
00923          if (deltaa <= 1e-7 * meana)
00924             continue;
00925 
00926          double meanb = (b -> GetLowerLimit() + b -> GetUpperLimit()) / 2.0;
00927          double deltab = (b -> GetUpperLimit() - b -> GetLowerLimit());
00928          if (deltab <= 1e-7 * meanb)
00929             continue;
00930 
00931          this -> GetMarginalized(a,b) -> Print(Form("%s_2D_%s_%s.ps",filebase,a->GetName().data(),b->GetName().data()));
00932          k++;
00933       }
00934    }
00935 
00936    return k;
00937 }

void BCModel::PrintHessianMatrix ( std::vector< double >  parameters  ) 

Prints matrix elements of the Hessian matrix

Parameters:
parameters The parameter values at which point to evaluate the matrix

Definition at line 1533 of file BCModel.cxx.

01534 {
01535    // check number of entries in vector
01536    if (parameters.size() != this -> GetNParameters())
01537    {
01538       BCLog::OutError("BCModel::PrintHessianMatrix : Invalid number of entries in the vector");
01539       return;
01540    }
01541 
01542    // print to screen
01543    std::cout
01544       << std::endl
01545       << " Hessian matrix elements: " << std::endl
01546       << " Point: ";
01547 
01548    for (int i = 0; i < int(parameters.size()); i++)
01549       std::cout << parameters.at(i) << " ";
01550    std::cout << std::endl;
01551 
01552    // loop over all parameter pairs
01553    for (unsigned int i = 0; i < this -> GetNParameters(); i++)
01554       for (unsigned int j = 0; j < i; j++)
01555       {
01556          // calculate Hessian matrix element
01557          double hessianmatrixelement = this -> HessianMatrixElement(fParameterSet -> at(i),
01558                fParameterSet -> at(j), parameters);
01559 
01560          // print to screen
01561          std::cout << " " << i << " " << j << " : " << hessianmatrixelement << std::endl;
01562       }
01563 }

void BCModel::PrintResults ( const char *  file  ) 

Prints a summary of the Markov Chain Monte Carlo to a file.

Definition at line 1362 of file BCModel.cxx.

01363 {
01364    // print summary of Markov Chain Monte Carlo
01365 
01366    // open file
01367    ofstream ofi(file);
01368 
01369    // check if file is open
01370    if(!ofi.is_open())
01371    {
01372       std::cerr << "Couldn't open file " << file <<std::endl;
01373       return;
01374    }
01375 
01376    // number of parameters
01377    int npar = fParameterSet -> size();
01378 
01379    // check convergence
01380    bool flag_conv = ((this -> MCMCGetNIterationsConvergenceGlobal() > 0)?1:0);
01381 
01382    ofi
01383       << std::endl
01384       << " -----------------------------------------------------" << std::endl
01385       << " Summary of the Markov Chain Monte Carlo run" << std::endl
01386       << " -----------------------------------------------------" << std::endl
01387       << std::endl;
01388 
01389    if (!flag_conv)
01390    {
01391       ofi
01392          << " WARNING: the Markov Chain did not converge! Be" << std::endl
01393          << " cautious using the following results!" << std::endl
01394          << std::endl;
01395    }
01396 
01397    ofi
01398       << " Model summary" << std::endl
01399       << " =============" << std::endl
01400       << " Model: " << fName.data() << std::endl
01401       << " Number of parameters: " << npar << std::endl
01402       << " List of Parameters and ranges:" << std::endl;
01403    for (int i = 0; i < npar; ++i)
01404    {
01405       ofi
01406          << "  (" << i << ") Parameter \""
01407          << fParameterSet -> at(i) -> GetName().data() << "\""
01408          << ": " << fParameterSet -> at(i) -> GetLowerLimit()
01409          << " - "
01410          << fParameterSet -> at(i) -> GetUpperLimit() << std::endl;
01411    }
01412    ofi << std::endl;
01413 
01414    if (flag_conv)
01415    {
01416       ofi
01417          << " Results of the marginalization" << std::endl
01418          << " ==============================" << std::endl
01419          << " List of parameters and properties of the marginalized" << std::endl
01420          << " distributions:" << std::endl;
01421       for (int i = 0; i < npar; ++i)
01422       {
01423          BCH1D * bch1d = this -> GetMarginalized(fParameterSet -> at(i));
01424 
01425          ofi
01426             << "  (" << i << ") Parameter \""
01427                << fParameterSet -> at(i) -> GetName().data() << "\"" << std::endl
01428             << "      Mean +- RMS:         "
01429                << std::setprecision(4) << bch1d -> GetMean()
01430                << " +- "
01431                << std::setprecision(4) << bch1d -> GetRMS() << std::endl
01432             << "      Median +- sigma:     "
01433                << std::setprecision(4) << bch1d -> GetMedian()
01434                << " +  " << std::setprecision(4) << bch1d -> GetQuantile(0.84) - bch1d -> GetMedian()
01435                << " - " << std::setprecision(4) << bch1d -> GetMedian() - bch1d -> GetQuantile(0.16) << std::endl
01436             << "      (Marginalized) mode: " << bch1d -> GetMode() << std::endl
01437             << "      Smallest interval(s) containing 68% and local modes:" << std::endl;
01438 
01439          std::vector <double> v;
01440          v = bch1d -> GetSmallestIntervals(0.68);
01441          int ninter = int(v.size());
01442 
01443          for (int j = 0; j < ninter; j+=5)
01444             ofi << "       " << v.at(j) << " - " << v.at(j+1) << " (local mode at " << v.at(j+3) << " with rel. height " << v.at(j+2) << "; rel. area " << v.at(j+4) << ")" << std::endl;
01445       }
01446       ofi << std::endl;
01447    }
01448 
01449    ofi
01450       << " Results of the optimization" << std::endl
01451       << " ===========================" << std::endl
01452       << " Optimization algorithm used: ";
01453    switch(this -> GetOptimizationMethod())
01454    {
01455       case BCIntegrate::kOptSimulatedAnnealing:
01456          ofi << " Simulated Annealing" << std::endl;
01457          break;
01458       case BCIntegrate::kOptMinuit:
01459          ofi << " Minuit" << std::endl;
01460          break;
01461       case BCIntegrate::kOptMetropolis:
01462          ofi << " MCMC " << std::endl;
01463          break;
01464    }
01465 
01466    ofi << " List of parameters and global mode:" << std::endl;
01467    for (int i = 0; i < npar; ++i)
01468       ofi
01469          << "  (" << i << ") Parameter \""
01470          << fParameterSet -> at(i) -> GetName().data() << "\": "
01471          << fBestFitParameters.at(i) << std::endl;
01472    ofi << std::endl;
01473 
01474    if (fPValue >= 0.)
01475    {
01476       ofi
01477          << " Results of the model test" << std::endl
01478          << " =========================" << std::endl
01479          << " p-value at global mode: " << fPValue << std::endl
01480          << std::endl;
01481    }
01482 
01483    ofi
01484       << " Status of the MCMC" << std::endl
01485       << " ==================" << std::endl
01486       << " Convergence reached: " << ((flag_conv)?"yes":"no") << std::endl;
01487 
01488    if (flag_conv)
01489       ofi << " Number of iterations until convergence: " << this -> MCMCGetNIterationsConvergenceGlobal() << std::endl;
01490    else
01491       ofi
01492          << " WARNING: the Markov Chain did not converge! Be\n"
01493          << " cautious using the following results!" << std::endl
01494          << std::endl;
01495    ofi
01496       << " Number of chains:                       " << this -> MCMCGetNChains() << std::endl
01497       << " Number of iterations of each chain:     " << this -> MCMCGetNIterationsMax() << std::endl
01498       << std::endl;
01499 
01500    ofi
01501       << " -----------------------------------------------------" << std::endl
01502       << std::endl
01503       << " Notes" << std::endl
01504       << " =====" << std::endl
01505       << " (i) Median +- sigma denotes the median, m, of the" << std::endl
01506       << "     marginalized distribution and the intervals from" << std::endl
01507       << "     m to the 16% and 84% quantiles." << std::endl
01508       << " -----------------------------------------------------" << std::endl;
01509 
01510    // close file
01511 // ofi.close;
01512 }

void BCModel::PrintShortFitSummary (  ) 

Prints a short summary of the fit results on the screen.

Definition at line 1516 of file BCModel.cxx.

01517 {
01518    BCLog::OutSummary("---------------------------------------------------");
01519    BCLog::OutSummary(Form("Fit summary for model \'%s\':",this -> GetName().data()));
01520    BCLog::OutSummary(Form("   Number of parameters = %i", this -> GetNParameters()));
01521 
01522    BCLog::OutSummary("   Best fit parameters (global):");
01523    for (unsigned int i = 0; i < this -> GetNParameters(); ++i)
01524       BCLog::OutSummary(Form("      %s = %.2lf", this -> GetParameter(i) -> GetName().data(), this -> GetBestFitParameter(i)));
01525 
01526    BCLog::OutSummary("   Goodness-of-fit test:");
01527    BCLog::OutSummary(Form("      p-value = %.2lf", this -> GetPValue()));
01528    BCLog::OutSummary("---------------------------------------------------");
01529 }

void BCModel::PrintSummary (  ) 

Prints a summary on the screen.

Definition at line 1299 of file BCModel.cxx.

01300 {
01301    int nparameters = this -> GetNParameters();
01302 
01303    // model summary
01304    std::cout
01305       << std::endl
01306       << "   ---------------------------------" << std::endl
01307       << "    Model : " << fName.data() << std::endl
01308       << "   ---------------------------------"<< std::endl
01309       << "     Index                : " << fIndex << std::endl
01310       << "     Number of parameters : " << nparameters << std::endl
01311       << std::endl
01312       << "     - Parameters : " << std::endl
01313       << std::endl;
01314 
01315    // parameter summary
01316    for (int i=0; i<nparameters; i++)
01317       fParameterSet -> at(i) -> PrintSummary();
01318 
01319    // best fit parameters
01320    if (this -> GetBestFitParameters().size() > 0)
01321    {
01322       std::cout
01323          << std::endl
01324          << "     - Best fit parameters :" << std::endl
01325          << std::endl;
01326 
01327       for (int i=0; i<nparameters; i++)
01328       {
01329          std::cout
01330             << "       " << fParameterSet -> at(i) -> GetName().data()
01331             << " = " << this -> GetBestFitParameter(i)
01332             << " (overall)" << std::endl;
01333          if ((int)fBestFitParametersMarginalized.size() == nparameters)
01334             std::cout
01335                << "       " << fParameterSet -> at(i) -> GetName().data()
01336                << " = " << this -> GetBestFitParameterMarginalized(i)
01337                << " (marginalized)" << std::endl;
01338       }
01339    }
01340 
01341    std::cout << std::endl;
01342 
01343    // model testing
01344    if (fPValue >= 0)
01345    {
01346       double likelihood = this -> Likelihood(this -> GetBestFitParameters());
01347       std::cout
01348          << "   - Model testing:" << std::endl
01349          << std::endl
01350          << "       p(data|lambda*) = " << likelihood << std::endl
01351          << "       p-value         = " << fPValue << std::endl
01352          << std::endl;
01353    }
01354 
01355    // normalization
01356    if (fNormalization > 0)
01357       std::cout << "     Normalization        : " << fNormalization << std::endl;
01358 }

double BCModel::Probability ( std::vector< double >  parameter  )  [inline]

Returns the a posteriori probability given a set of parameter values

Parameters:
parameters A set of parameter values
Returns:
The a posteriori probability

Definition at line 373 of file BCModel.h.

00374          { return exp( this->LogProbability(parameter) ); };

double BCModel::ProbabilityNN ( std::vector< double >  parameter  )  [inline]

Returns the likelihood times prior probability given a set of parameter values

Parameters:
parameters A set of parameter values
Returns:
The likelihood times prior probability

Definition at line 359 of file BCModel.h.

00360          { return exp( this->LogProbabilityNN(parameter) ); };

int BCModel::ReadErrorBandFromFile ( const char *  file  ) 

Read

Definition at line 853 of file BCModel.cxx.

00854 {
00855    TFile * froot = new TFile(file);
00856    if(!froot->IsOpen())
00857    {
00858       BCLog::OutWarning(Form("BCModel::ReadErrorBandFromFile. Couldn't open file %s.",file));
00859       return 0;
00860    }
00861 
00862    int r=0;
00863 
00864    TH2D * h2 = (TH2D*)froot->Get("errorbandxy");
00865    if(h2)
00866    {
00867       h2->SetDirectory(0);
00868       h2->SetName(TString::Format("errorbandxy_%d",BCLog::GetHIndex()));
00869       this->SetErrorBandHisto(h2);
00870       r=1;
00871    }
00872    else
00873       BCLog::OutWarning(Form(
00874             "BCModel::ReadErrorBandFromFile : Couldn't read histogram \"errorbandxy\" from file %s.",
00875             file));
00876 
00877    froot->Close();
00878 
00879    return r;
00880 }

int BCModel::ReadMarginalizedFromFile ( const char *  file  ) 

Read

Definition at line 789 of file BCModel.cxx.

00790 {
00791    TFile * froot = new TFile(file);
00792    if(!froot->IsOpen())
00793    {
00794       BCLog::OutError(Form("BCModel::ReadMarginalizedFromFile : Couldn't open file %s.",file));
00795       return 0;
00796    }
00797 
00798    // We reset the MCMCEngine here for the moment.
00799    // In the future maybe we only want to do this if the engine
00800    // wans't initialized at all or when there were some changes
00801    // in the model.
00802    // But maybe we want reset everything since we're overwriting
00803    // the marginalized distributions anyway.
00804    this -> MCMCInitialize();
00805 
00806    int k=0;
00807    int n=this->GetNParameters();
00808    for(int i=0;i<n;i++)
00809    {
00810       BCParameter * a = this -> GetParameter(i);
00811       TKey * key = froot -> GetKey(TString::Format("hist_%s_%s", this -> GetName().data(), a -> GetName().data()));
00812       if(key)
00813       {
00814          TH1D * h1 = (TH1D*) key -> ReadObjectAny(TH1D::Class());
00815          h1->SetDirectory(0);
00816          if(this->SetMarginalized(i,h1))
00817             k++;
00818       }
00819       else
00820          BCLog::OutWarning(Form(
00821                "BCModel::ReadMarginalizedFromFile : Couldn't read histogram \"hist_%s_%s\" from file %s.",
00822                this -> GetName().data(), a -> GetName().data(), file));
00823    }
00824 
00825    for(int i=0;i<n-1;i++)
00826    {
00827       for(int j=i+1;j<n;j++)
00828       {
00829          BCParameter * a = this -> GetParameter(i);
00830          BCParameter * b = this -> GetParameter(j);
00831          TKey * key = froot -> GetKey(TString::Format("hist_%s_%s_%s", this -> GetName().data(), a -> GetName().data(), b -> GetName().data()));
00832          if(key)
00833          {
00834             TH2D * h2 = (TH2D*) key -> ReadObjectAny(TH2D::Class());
00835             h2->SetDirectory(0);
00836             if(this->SetMarginalized(i,j,h2))
00837                k++;
00838          }
00839          else
00840             BCLog::OutWarning(Form(
00841                   "BCModel::ReadMarginalizedFromFile : Couldn't read histogram \"hist_%s_%s_%s\" from file %s.",
00842                   this -> GetName().data(), a -> GetName().data(), b -> GetName().data(), file));
00843       }
00844    }
00845 
00846    froot->Close();
00847 
00848    return k;
00849 }

int BCModel::ReadMode ( const char *  file  ) 

Read mode from file created by WriteMode() call

Definition at line 664 of file BCModel.cxx.

00665 {
00666    ifstream ifi(file);
00667    if(!ifi.is_open())
00668    {
00669       BCLog::OutError(Form("BCModel::ReadMode : Couldn't open file %s.",file));
00670       return 0;
00671    }
00672 
00673    int npar = fParameterSet -> size();
00674    std::vector <double> mode;
00675 
00676    int i=0;
00677    while (i<npar && !ifi.eof())
00678    {
00679       double a;
00680       ifi>>a;
00681       mode.push_back(a);
00682       i++;
00683    }
00684 
00685    if(i<npar)
00686    {
00687       BCLog::OutError(Form("BCModel::ReadMode : Couldn't read mode from file %s.",file));
00688       BCLog::OutError(Form("BCModel::ReadMode : Expected %d parameters, found %d.",npar,i));
00689       return 0;
00690    }
00691 
00692    BCLog::OutSummary(Form("#  Read in best fit parameters (mode) for model \'%s\' from file %s:",fName.data(),file));
00693    this->SetMode(mode);
00694    for(int j=0 ; j<npar; j++)
00695       BCLog::OutSummary(Form("#    -> Parameter %d : %s = %e", j, fParameterSet->at(j)->GetName().data(), fBestFitParameters[j]));
00696 
00697    BCLog::OutWarning("#  ! Best fit values obtained before this call will be overwritten !");
00698 
00699    return npar;
00700 }

double BCModel::SamplingFunction ( std::vector< double >  parameters  )  [virtual]

Sampling function used for importance sampling. Method needs to be overloaded by the user.

Parameters:
parameters A set of parameter values
Returns:
The probability density at the parameter values

Definition at line 531 of file BCModel.cxx.

00532 {
00533    double probability = 1.0;
00534    for (std::vector<BCParameter*>::const_iterator it = fParameterSet -> begin(); it != fParameterSet -> end(); ++it)
00535       probability *= 1.0 / ((*it) -> GetUpperLimit() - (*it) -> GetLowerLimit());
00536    return probability;
00537 }

void BCModel::SetDataBoundaries ( unsigned int  index,
double  lowerboundary,
double  upperboundary,
bool  fixed = false 
)

Definition at line 371 of file BCModel.cxx.

00372 {
00373    // check if data set exists
00374    if (!fDataSet)
00375    {
00376       BCLog::OutError("BCModel::SetDataBoundaries : Need to define data set first.");
00377       return;
00378    }
00379 
00380    // check if index is within range
00381    if (index < 0 || index > fDataSet -> GetDataPoint(0) -> GetNValues())
00382    {
00383       BCLog::OutError("BCModel::SetDataBoundaries : Index out of range.");
00384       return;
00385    }
00386 
00387    // check if boundary data points exist
00388    if (!fDataPointLowerBoundaries)
00389       fDataPointLowerBoundaries = new BCDataPoint(fDataSet -> GetDataPoint(0) -> GetNValues());
00390 
00391    if (!fDataPointUpperBoundaries)
00392       fDataPointUpperBoundaries = new BCDataPoint(fDataSet -> GetDataPoint(0) -> GetNValues());
00393 
00394    if (fDataFixedValues.size() == 0)
00395       fDataFixedValues.assign(fDataSet -> GetDataPoint(0) -> GetNValues(), false);
00396 
00397    // set boundaries
00398    fDataPointLowerBoundaries -> SetValue(index, lowerboundary);
00399    fDataPointUpperBoundaries -> SetValue(index, upperboundary);
00400    fDataFixedValues[index] = fixed;
00401 }

void BCModel::SetDataSet ( BCDataSet dataset  )  [inline]

Sets the data set.

Parameters:
dataset A data set

Definition at line 273 of file BCModel.h.

00274          { fDataSet = dataset; fNormalization = -1.0; };

void BCModel::SetErrorBandContinuous ( bool  flag  ) 

Sets the error band flag to continuous function

Definition at line 405 of file BCModel.cxx.

00406 {
00407    fErrorBandContinuous = flag;
00408 
00409    if (flag)
00410       return;
00411 
00412    // clear x-values
00413    fErrorBandX.clear();
00414 
00415    // copy data x-values
00416    for (unsigned int i = 0; i < fDataSet -> GetNDataPoints(); ++i)
00417       fErrorBandX.push_back(fDataSet -> GetDataPoint(i) -> GetValue(fFitFunctionIndexX));
00418 }

void BCModel::SetGoFNChains ( int  n  )  [inline]

Definition at line 536 of file BCModel.h.

00537          { fGoFNChains=n; };

void BCModel::SetGoFNIterationsMax ( int  n  )  [inline]

Definition at line 524 of file BCModel.h.

00525          { fGoFNIterationsMax=n; };

void BCModel::SetGoFNIterationsRun ( int  n  )  [inline]

Definition at line 530 of file BCModel.h.

00531          { fGoFNIterationsRun=n; };

void BCModel::SetIndex ( int  index  )  [inline]

Sets the index of the model within the BCModelManager.

Parameters:
index The index of the model

Definition at line 240 of file BCModel.h.

00241          { fIndex = index; };

void BCModel::SetModelAPosterioriProbability ( double  probability  )  [inline]

Sets the a posteriori probability for a model.

Parameters:
model The model
probability The a posteriori probability

Definition at line 260 of file BCModel.h.

00261          { fModelAPosteriori = probability; };

void BCModel::SetModelAPrioriProbability ( double  probability  )  [inline]

Sets the a priori probability for a model.

Parameters:
model The model
probability The a priori probability

Definition at line 253 of file BCModel.h.

00254          { fModelAPriori = probability; };

void BCModel::SetName ( const char *  name  )  [inline]

Sets the name of the model.

Parameters:
name Name of the model

Definition at line 234 of file BCModel.h.

00235          { fName = name; };

void BCModel::SetNbins ( const char *  parname,
int  nbins 
)

Set the number of bins for the marginalized distribution of a parameter.

Parameters:
parname The name of the parameter in the parameter set
nbins Number of bins (default = 100)

Definition at line 166 of file BCModel.cxx.

00167 {
00168    BCParameter * p = this -> GetParameter(parname);
00169    if(!p)
00170    {
00171       BCLog::OutWarning(Form("BCModel::SetNbins : Parameter '%s' not found so Nbins not set",parname));
00172       return;
00173    }
00174 
00175    this -> BCIntegrate::SetNbins(nbins, p -> GetIndex());
00176 }

void BCModel::SetNDataPointsMaximum ( unsigned int  maximum  )  [inline]

Sets the maximum number of data points.

Definition at line 290 of file BCModel.h.

00291          { fNDataPointsMaximum = maximum; };

void BCModel::SetNDataPointsMinimum ( unsigned int  minimum  )  [inline]

Sets the minimum number of data points.

Definition at line 285 of file BCModel.h.

00286          { fNDataPointsMinimum = minimum; };

void BCModel::SetNormalization ( double  norm  )  [inline]

Sets the normalization of the likelihood. The normalization is the integral of likelihood over all parameters.

Parameters:
norm The normalization of the likelihood

Definition at line 267 of file BCModel.h.

00268          { fNormalization = norm; };

void BCModel::SetParameterSet ( BCParameterSet parset  )  [inline]

Set all parameters of the model using a BCParameterSet container.

pointer to parameter set

Definition at line 246 of file BCModel.h.

00247          { fParameterSet = parset; };

void BCModel::SetSingleDataPoint ( BCDataSet dataset,
unsigned int  index 
)

Definition at line 361 of file BCModel.cxx.

00362 {
00363    if (index < 0 || index > dataset -> GetNDataPoints())
00364       return;
00365 
00366    this -> SetSingleDataPoint(dataset -> GetDataPoint(index));
00367 }

void BCModel::SetSingleDataPoint ( BCDataPoint datapoint  ) 

Sets a single data point as data set.

Parameters:
datapoint A data point

Definition at line 346 of file BCModel.cxx.

00347 {
00348    // create new data set consisting of a single data point
00349    BCDataSet * dataset = new BCDataSet();
00350 
00351    // add the data point
00352    dataset -> AddDataPoint(datapoint);
00353 
00354    // set this new data set
00355    this -> SetDataSet(dataset);
00356 
00357 }

BCDataPoint * BCModel::VectorToDataPoint ( std::vector< double >  data  )  [private]

Converts a vector of doubles into a BCDataPoint

Definition at line 1567 of file BCModel.cxx.

01568 {
01569    int sizeofvector = int(data.size());
01570    BCDataPoint * datapoint = new BCDataPoint(sizeofvector);
01571    datapoint -> SetValues(data);
01572    return datapoint;
01573 }

void BCModel::WriteMode ( const char *  file  ) 

Write mode into file

Definition at line 623 of file BCModel.cxx.

00624 {
00625    ofstream ofi(file);
00626    if(!ofi.is_open())
00627    {
00628       std::cerr<<"Couldn't open file "<<file<<std::endl;
00629       return;
00630    }
00631 
00632    int npar = fParameterSet -> size();
00633    for (int i=0; i<npar; i++)
00634       ofi<<fBestFitParameters.at(i)<<std::endl;
00635 
00636    ofi<<std::endl;
00637    ofi<<"#######################################################################"<<std::endl;
00638    ofi<<"#"<<std::endl;
00639    ofi<<"#  This file was created automatically by BCModel::WriteMode() call."<<std::endl;
00640    ofi<<"#  It can be read in by call to BCModel::ReadMode()."<<std::endl;
00641    ofi<<"#  Do not modify it unless you know what you're doing."<<std::endl;
00642    ofi<<"#"<<std::endl;
00643    ofi<<"#######################################################################"<<std::endl;
00644    ofi<<"#"<<std::endl;
00645    ofi<<"#  Best fit parameters (mode) for model:"<<std::endl;
00646    ofi<<"#  \'"<<fName.data()<<"\'"<<std::endl;
00647    ofi<<"#"<<std::endl;
00648    ofi<<"#  Number of parameters: "<<npar<<std::endl;
00649    ofi<<"#  Parameters ordered as above:"<<std::endl;
00650 
00651    for (int i=0; i<npar; i++)
00652    {
00653       ofi<<"#     "<<i<<": ";
00654       ofi<<fParameterSet->at(i)->GetName().data()<<" = ";
00655       ofi<<fBestFitParameters.at(i)<<std::endl;
00656    }
00657 
00658    ofi<<"#"<<std::endl;
00659    ofi<<"########################################################################"<<std::endl;
00660 }


Member Data Documentation

A data set

Definition at line 591 of file BCModel.h.

int BCModel::fGoFNChains [protected]

Definition at line 622 of file BCModel.h.

int BCModel::fGoFNIterationsMax [protected]

Definition at line 612 of file BCModel.h.

int BCModel::fGoFNIterationsRun [protected]

Definition at line 617 of file BCModel.h.

int BCModel::fIndex [protected]

Index of the model.

Definition at line 571 of file BCModel.h.

A flag for overloading ConditionalProbabilityEntry

Definition at line 603 of file BCModel.h.

double BCModel::fModelAPosteriori [protected]

The model a posteriori probability.

Definition at line 583 of file BCModel.h.

double BCModel::fModelAPriori [protected]

The model prior probability.

Definition at line 579 of file BCModel.h.

std::string BCModel::fName [protected]

Name of the model.

Definition at line 575 of file BCModel.h.

unsigned int BCModel::fNDataPointsMaximum [protected]

Maximum number of data points

Definition at line 599 of file BCModel.h.

unsigned int BCModel::fNDataPointsMinimum [protected]

Minimum number of data points

Definition at line 595 of file BCModel.h.

double BCModel::fNormalization [private]

The Likelihood normalization.

Definition at line 636 of file BCModel.h.

A model parameter container.

Definition at line 587 of file BCModel.h.

double BCModel::fPValue [protected]

The p-value

Definition at line 607 of file BCModel.h.


Generated on Tue Apr 7 17:39:29 2009 for Bayesian Analysis Toolkit by  doxygen 1.5.8