BayesianAnalysisToolkit  0.9.3
Public Member Functions | Private Attributes | List of all members
BCMVCDataModel Class Reference

#include <BCMVCDataModel.h>

Inheritance diagram for BCMVCDataModel:
Inheritance graph
[legend]
Collaboration diagram for BCMVCDataModel:
Collaboration graph
[legend]

Public Member Functions

 BCMVCDataModel (BCMVCombination *mvc)
 
 ~BCMVCDataModel ()
 
void SetNMeasurements (int n, double min, double max)
 
void SetVectorMeasurements (TVectorD measurements)
 
void SetParameters (std::vector< double > parameters)
 
void SetMeasurementRanges (const std::vector< double > &min, const std::vector< double > &max)
 
void SetMeasurementRanges (double min, double max)
 
void SetVectorObservable (std::vector< int > vec)
 
void SetCovarianceMatrix (TMatrixD matrix)
 
void SetHistChi2 (TH1D *hist)
 
void PrintToys (std::string filename)
 
void PrintSummary ()
 
double Chi2 (TVectorD observables, TVectorD measurements)
 
double LogAPrioriProbability (const std::vector< double > &parameters)
 
double LogLikelihood (const std::vector< double > &parameters)
 
void MCMCUserIterationInterface ()
 
- Public Member Functions inherited from BCModel
 BCModel (const char *name="model")
 
 BCModel (const BCModel &bcmodel)
 
virtual ~BCModel ()
 
BCModeloperator= (const BCModel &bcmodel)
 
const std::string & GetName () const
 
double GetModelAPrioriProbability () const
 
double GetModelAPosterioriProbability () const
 
BCDataSetGetDataSet () const
 
BCDataPointGetDataPointLowerBoundaries () const
 
BCDataPointGetDataPointUpperBoundaries () const
 
double GetDataPointLowerBoundary (unsigned int index) const
 
double GetDataPointUpperBoundary (unsigned int index) const
 
bool GetFlagBoundaries () const
 
unsigned GetNDataPoints () const
 
BCDataPointGetDataPoint (unsigned int index) const
 
void SetName (const char *name)
 
void SetModelAPrioriProbability (double probability)
 
void SetModelAPosterioriProbability (double probability)
 
virtual int AddParameter (BCParameter *parameter)
 
void SetDataSet (BCDataSet *dataset)
 
void SetSingleDataPoint (BCDataPoint *datapoint)
 
void SetSingleDataPoint (BCDataSet *dataset, unsigned int index)
 
void SetDataBoundaries (unsigned int index, double lowerboundary, double upperboundary, bool fixed=false)
 
void SetDataPointLowerBoundaries (BCDataPoint *datasetlowerboundaries)
 
void SetDataPointUpperBoundaries (BCDataPoint *datasetupperboundaries)
 
void SetDataPointLowerBoundary (int index, double lowerboundary)
 
void SetDataPointUpperBoundary (int index, double upperboundary)
 
int SetPrior (int index, TF1 *f)
 
int SetPrior (const char *name, TF1 *f)
 
int SetPriorDelta (int index, double value)
 
int SetPriorDelta (const char *name, double value)
 
int SetPriorGauss (int index, double mean, double sigma)
 
int SetPriorGauss (const char *name, double mean, double sigma)
 
int SetPriorGauss (int index, double mean, double sigmadown, double sigmaup)
 
int SetPriorGauss (const char *name, double mean, double sigmadown, double sigmaup)
 
int SetPrior (int index, TH1 *h, bool flag=false)
 
int SetPrior (const char *name, TH1 *h, bool flag=false)
 
int SetPriorConstant (int index)
 
int SetPriorConstant (const char *name)
 
int SetPriorConstantAll ()
 
void Copy (const BCModel &bcmodel)
 
double APrioriProbability (const std::vector< double > &parameters)
 
virtual double Likelihood (const std::vector< double > &params)
 
double ProbabilityNN (const std::vector< double > &params)
 
double LogProbabilityNN (const std::vector< double > &parameters)
 
double Probability (const std::vector< double > &parameter)
 
double LogProbability (const std::vector< double > &parameter)
 
virtual double SamplingFunction (const std::vector< double > &parameters)
 
double Eval (const std::vector< double > &parameters)
 
virtual double LogEval (const std::vector< double > &parameters)
 
virtual void CorrelateDataPointValues (std::vector< double > &x)
 
double GetPvalueFromChi2 (const std::vector< double > &par, int sigma_index)
 
double GetPvalueFromKolmogorov (const std::vector< double > &par, int index)
 
double GetPvalueFromChi2NDoF (std::vector< double > par, int sigma_index)
 
BCH1DCalculatePValue (std::vector< double > par, bool flag_histogram=false)
 
double GetPValue ()
 
double GetPValueNDoF ()
 
double GetChi2NDoF ()
 
std::vector< double > GetChi2Runs (int dataIndex, int sigmaIndex)
 
void SetGoFNIterationsMax (int n)
 
void SetGoFNIterationsRun (int n)
 
void SetGoFNChains (int n)
 
double HessianMatrixElement (const BCParameter *parameter1, const BCParameter *parameter2, std::vector< double > point)
 
void PrintSummary ()
 
void PrintResults (const char *file)
 
void PrintShortFitSummary (int chi2flag=0)
 
void PrintHessianMatrix (std::vector< double > parameters)
 
virtual double CDF (const std::vector< double > &, int, bool)
 
- Public Member Functions inherited from BCIntegrate
 BCIntegrate ()
 
 BCIntegrate (const BCIntegrate &bcintegrate)
 
virtual ~BCIntegrate ()
 
BCIntegrateoperator= (const BCIntegrate &bcintegrate)
 
int ReadMarginalizedFromFile (const char *file)
 
BCH1DGetMarginalized (const BCParameter *parameter)
 
BCH1DGetMarginalized (const char *name)
 
BCH1DGetMarginalized (unsigned index)
 
BCH2DGetMarginalized (const BCParameter *parameter1, const BCParameter *parameter2)
 
BCH2DGetMarginalized (const char *name1, const char *name2)
 
BCH2DGetMarginalized (unsigned index1, unsigned index2)
 
int PrintAllMarginalized1D (const char *filebase)
 
int PrintAllMarginalized2D (const char *filebase)
 
int PrintAllMarginalized (const char *file, std::string options1d="BTsiB3CS1D0pdf0Lmeanmode", std::string options2d="BTfB3CS1meangmode", unsigned int hdiv=1, unsigned int ndiv=1)
 
double GetIntegral () const
 
BCIntegrate::BCOptimizationMethod GetOptimizationMethod () const
 
BCIntegrate::BCIntegrationMethod GetIntegrationMethod () const
 
BCIntegrate::BCMarginalizationMethod GetMarginalizationMethod () const
 
BCIntegrate::BCSASchedule GetSASchedule () const
 
void GetRandomVectorUnitHypercube (std::vector< double > &x) const
 
void GetRandomVectorInParameterSpace (std::vector< double > &x) const
 
double GetRandomPoint (std::vector< double > &x)
 
int GetNIterationsMin () const
 
int GetNIterationsMax () const
 
int GetNIterationsPrecisionCheck () const
 
int GetNIterationsOutput () const
 
int GetNIterations () const
 
double GetRelativePrecision () const
 
double GetAbsolutePrecision () const
 
BCCubaMethod GetCubaIntegrationMethod () const
 
const BCCubaOptions::VegasGetCubaVegasOptions () const
 
const BCCubaOptions::SuaveGetCubaSuaveOptions () const
 
const BCCubaOptions::DivonneGetCubaDivonneOptions () const
 
const BCCubaOptions::CuhreGetCubaCuhreOptions () const
 
BCH1DGetSlice (const BCParameter *parameter, const std::vector< double > parameters=std::vector< double >(0), int bins=0, bool flag_norm=true)
 
BCH1DGetSlice (const char *name, const std::vector< double > parameters=std::vector< double >(0), int nbins=0, bool flag_norm=true)
 
BCH2DGetSlice (const BCParameter *parameter1, const BCParameter *parameter2, const std::vector< double > parameters=std::vector< double >(0), int bins=0, bool flag_norm=true)
 
BCH2DGetSlice (const char *name1, const char *name2, const std::vector< double > parameters=std::vector< double >(0), int nbins=0, bool flag_norm=true)
 
BCH2DGetSlice (unsigned index1, unsigned index2, const std::vector< double > parameters=std::vector< double >(0), int nbins=0, bool flag_norm=true)
 
double GetError () const
 
TMinuit * GetMinuit ()
 
int GetMinuitErrorFlag () const
 
double GetSAT0 () const
 
double GetSATmin () const
 
double GetBestFitParameter (unsigned index) const
 
double GetBestFitParameterError (unsigned index) const
 
double GetLogMaximum ()
 
const std::vector< double > & GetBestFitParameters () const
 
const std::vector< double > & GetBestFitParameterErrors () const
 
void SetMinuitArlist (double *arglist)
 
void SetFlagIgnorePrevOptimization (bool flag)
 
void SetOptimizationMethod (BCIntegrate::BCOptimizationMethod method)
 
void SetIntegrationMethod (BCIntegrate::BCIntegrationMethod method)
 
void SetMarginalizationMethod (BCIntegrate::BCMarginalizationMethod method)
 
void SetSASchedule (BCIntegrate::BCSASchedule schedule)
 
void SetNIterationsMin (int niterations)
 
void SetNIterationsMax (int niterations)
 
void SetNIterationsPrecisionCheck (int niterations)
 
void SetNIterationsOutput (int niterations)
 
void SetRelativePrecision (double relprecision)
 
void SetAbsolutePrecision (double absprecision)
 
void SetCubaIntegrationMethod (BCCubaMethod type)
 
void SetCubaOptions (const BCCubaOptions::Vegas &options)
 
void SetCubaOptions (const BCCubaOptions::Suave &options)
 
void SetCubaOptions (const BCCubaOptions::Divonne &options)
 
void SetCubaOptions (const BCCubaOptions::Cuhre &options)
 
void SetSAT0 (double T0)
 
void SetSATmin (double Tmin)
 
void SetFlagWriteSAToFile (bool flag)
 
TTree * GetSATree ()
 
void InitializeSATree ()
 
double Normalize ()
 
double Integrate (BCIntegrationMethod intmethod)
 
double Integrate ()
 
double Integrate (BCIntegrationMethod type, tRandomizer randomizer, tEvaluator evaluator, tIntegralUpdater updater, std::vector< double > &sums)
 
double EvaluatorMC (std::vector< double > &sums, const std::vector< double > &point, bool &accepted)
 
int MarginalizeAll ()
 
int MarginalizeAll (BCMarginalizationMethod margmethod)
 
virtual void MarginalizePreprocess ()
 
virtual void MarginalizePostprocess ()
 
void SAInitialize ()
 
std::vector< double > FindMode (std::vector< double > start=std::vector< double >())
 
std::vector< double > FindMode (BCIntegrate::BCOptimizationMethod optmethod, std::vector< double > start=std::vector< double >())
 
double SATemperature (double t)
 
double SATemperatureBoltzmann (double t)
 
double SATemperatureCauchy (double t)
 
virtual double SATemperatureCustom (double t)
 
std::vector< double > GetProposalPointSA (const std::vector< double > &x, int t)
 
std::vector< double > GetProposalPointSABoltzmann (const std::vector< double > &x, int t)
 
std::vector< double > GetProposalPointSACauchy (const std::vector< double > &x, int t)
 
virtual std::vector< double > GetProposalPointSACustom (const std::vector< double > &x, int t)
 
std::vector< double > SAHelperGetRandomPointOnHypersphere ()
 
double SAHelperGetRadialCauchy ()
 
double SAHelperSinusToNIntegral (int dim, double theta)
 
virtual void ResetResults ()
 
std::string DumpIntegrationMethod (BCIntegrationMethod type)
 
std::string DumpCurrentIntegrationMethod ()
 
std::string DumpUsedIntegrationMethod ()
 
std::string DumpMarginalizationMethod (BCMarginalizationMethod type)
 
std::string DumpCurrentMarginalizationMethod ()
 
std::string DumpUsedMarginalizationMethod ()
 
std::string DumpOptimizationMethod (BCOptimizationMethod type)
 
std::string DumpCurrentOptimizationMethod ()
 
std::string DumpUsedOptimizationMethod ()
 
std::string DumpCubaIntegrationMethod (BCCubaMethod type)
 
std::string DumpCubaIntegrationMethod ()
 
void SetBestFitParameters (const std::vector< double > &x)
 
void SetBestFitParameters (const std::vector< double > &x, const double &new_value, double &old_value)
 
unsigned GetNIntegrationVariables ()
 
double CalculateIntegrationVolume ()
 
bool CheckMarginalizationAvailability (BCMarginalizationMethod type)
 
bool CheckMarginalizationIndices (TH1 *hist, const std::vector< unsigned > &index)
 
- Public Member Functions inherited from BCEngineMCMC
void WriteMarkovChain (bool flag)
 
 BCEngineMCMC ()
 
 BCEngineMCMC (const BCEngineMCMC &enginemcmc)
 
virtual ~BCEngineMCMC ()
 
BCEngineMCMCoperator= (const BCEngineMCMC &engineMCMC)
 
unsigned MCMCGetNChains () const
 
unsigned MCMCGetNLag () const
 
const std::vector< unsigned > & MCMCGetNIterations () const
 
unsigned MCMCGetCurrentIteration () const
 
unsigned MCMCGetCurrentChain () const
 
unsigned MCMCGetNIterationsConvergenceGlobal () const
 
bool MCMCGetFlagConvergenceGlobal () const
 
unsigned MCMCGetNIterationsMax () const
 
unsigned MCMCGetNIterationsRun () const
 
unsigned MCMCGetNIterationsPreRunMin () const
 
unsigned MCMCGetNIterationsUpdate () const
 
unsigned MCMCGetNIterationsUpdateMax () const
 
std::vector< int > MCMCGetNTrialsTrue () const
 
int MCMCGetNTrials () const
 
const std::vector< double > & MCMCGetprobMean () const
 
const std::vector< double > & MCMCGetVariance () const
 
const std::vector< double > & MCMCGetTrialFunctionScaleFactor () const
 
std::vector< double > MCMCGetTrialFunctionScaleFactor (unsigned ichain) const
 
double MCMCGetTrialFunctionScaleFactor (unsigned ichain, unsigned ipar)
 
const std::vector< double > & MCMCGetx () const
 
std::vector< double > MCMCGetx (unsigned ichain)
 
double MCMCGetx (unsigned ichain, unsigned ipar) const
 
const std::vector< double > & MCMCGetLogProbx () const
 
double MCMCGetLogProbx (unsigned ichain)
 
int MCMCGetPhase () const
 
const std::vector< double > & MCMCGetMaximumPoints () const
 
std::vector< double > MCMCGetMaximumPoint (unsigned i) const
 
const std::vector< double > & MCMCGetMaximumLogProb () const
 
int MCMCGetFlagInitialPosition () const
 
double MCMCGetRValueCriterion () const
 
double MCMCGetRValueParametersCriterion () const
 
double MCMCGetRValue () const
 
double MCMCGetRValueParameters (unsigned i)
 
bool MCMCGetRValueStrict () const
 
bool MCMCGetFlagRun () const
 
TTree * MCMCGetMarkovChainTree (unsigned i)
 
BCH1DMCMCGetH1Marginalized (unsigned i)
 
BCH2DMCMCGetH2Marginalized (unsigned i, unsigned j)
 
BCParameterGetParameter (int index) const
 
BCParameterGetParameter (const char *name) const
 
unsigned int GetNParameters () const
 
unsigned int GetNFixedParameters ()
 
unsigned int GetNFreeParameters ()
 
const std::vector< double > & GetBestFitParametersMarginalized () const
 
void MCMCSetTrialFunctionScaleFactor (std::vector< double > scale)
 
void MCMCSetNChains (unsigned n)
 
void MCMCSetNLag (unsigned n)
 
void MCMCSetNIterationsMax (unsigned n)
 
void MCMCSetNIterationsRun (unsigned n)
 
void MCMCSetNIterationsPreRunMin (unsigned n)
 
void MCMCSetNIterationsUpdate (unsigned n)
 
void MCMCSetNIterationsUpdateMax (unsigned n)
 
void MCMCSetMinimumEfficiency (double efficiency)
 
void MCMCSetMaximumEfficiency (double efficiency)
 
void MCMCSetRandomSeed (unsigned seed)
 
void MCMCSetWriteChainToFile (bool flag)
 
void MCMCSetWritePreRunToFile (bool flag)
 
void MCMCSetInitialPositions (const std::vector< double > &x0s)
 
void MCMCSetInitialPositions (std::vector< std::vector< double > > x0s)
 
void MCMCSetFlagInitialPosition (int flag)
 
void MCMCSetFlagOrderParameters (bool flag)
 
void MCMCSetFlagFillHistograms (bool flag)
 
void MCMCSetFlagPreRun (bool flag)
 
void MCMCSetRValueCriterion (double r)
 
void MCMCSetRValueParametersCriterion (double r)
 
void MCMCSetRValueStrict (bool strict=true)
 
void MCMCSetMarkovChainTrees (const std::vector< TTree * > &trees)
 
void MCMCInitializeMarkovChainTrees ()
 
int SetMarginalized (unsigned index, TH1D *h)
 
int SetMarginalized (unsigned index1, unsigned index2, TH2D *h)
 
void MCMCSetValuesDefault ()
 
void MCMCSetValuesQuick ()
 
void MCMCSetValuesDetail ()
 
void MCMCSetPrecision (BCEngineMCMC::Precision precision)
 
void SetNbins (unsigned int nbins)
 
void Copy (const BCEngineMCMC &enginemcmc)
 
virtual int AddParameter (const char *name, double min, double max, const char *latexname="")
 
virtual void MCMCTrialFunction (unsigned ichain, std::vector< double > &x)
 
virtual double MCMCTrialFunctionSingle (unsigned ichain, unsigned ipar)
 
bool MCMCGetProposalPointMetropolis (unsigned chain, std::vector< double > &x)
 
bool MCMCGetProposalPointMetropolis (unsigned chain, unsigned parameter, std::vector< double > &x)
 
bool MCMCGetNewPointMetropolis (unsigned chain=0)
 
bool MCMCGetNewPointMetropolis (unsigned chain, unsigned parameter)
 
void MCMCInChainCheckMaximum ()
 
void MCMCInChainUpdateStatistics ()
 
void MCMCInChainFillHistograms ()
 
void MCMCInChainTestConvergenceAllChains ()
 
void MCMCInChainWriteChains ()
 
int MCMCMetropolis ()
 
int MCMCMetropolisPreRun ()
 
void MCMCResetRunStatistics ()
 
void MCMCInitializeMarkovChains ()
 
int MCMCInitialize ()
 
void ClearParameters (bool hard=false)
 
virtual void MCMCIterationInterface ()
 
virtual void MCMCCurrentPointInterface (std::vector< double > &, int, bool)
 

Private Attributes

TVectorD fVectorMeasurements
 
TMatrixD fCovarianceMatrix
 
TMatrixD fInvCovarianceMatrix
 
double fDetCovariance
 
std::vector< double > fPars
 
TVectorD fVectorObservables
 
std::vector< int > fVectorObservable
 
TH1D * fHistChi2
 

Additional Inherited Members

- Public Types inherited from BCIntegrate
enum  BCOptimizationMethod {
  kOptEmpty, kOptSimAnn, kOptMetropolis, kOptMinuit,
  kOptDefault, NOptMethods
}
 
enum  BCIntegrationMethod {
  kIntEmpty, kIntMonteCarlo, kIntCuba, kIntGrid,
  kIntDefault, NIntMethods
}
 
enum  BCMarginalizationMethod {
  kMargEmpty, kMargMetropolis, kMargMonteCarlo, kMargGrid,
  kMargDefault, NMargMethods
}
 
enum  BCSASchedule { kSACauchy, kSABoltzmann, kSACustom, NSAMethods }
 
enum  BCCubaMethod {
  kCubaVegas, kCubaSuave, kCubaDivonne, kCubaCuhre,
  NCubaMethods
}
 
typedef void(BCIntegrate::* tRandomizer )(std::vector< double > &) const
 
typedef double(BCIntegrate::* tEvaluator )(std::vector< double > &, const std::vector< double > &, bool &)
 
typedef void(* tIntegralUpdater )(const std::vector< double > &, const int &, double &, double &)
 
- Static Public Member Functions inherited from BCIntegrate
static void IntegralUpdaterMC (const std::vector< double > &sums, const int &nIterations, double &integral, double &absprecision)
 
static int CubaIntegrand (const int *ndim, const double xx[], const int *ncomp, double ff[], void *userdata)
 
static void FCNLikelihood (int &npar, double *grad, double &fval, double *par, int flag)
 
- Protected Member Functions inherited from BCIntegrate
unsigned IntegrationOutputFrequency () const
 
void LogOutputAtStartOfIntegration (BCIntegrationMethod type, BCCubaMethod cubatype)
 
void LogOutputAtIntegrationStatusUpdate (BCIntegrationMethod type, double integral, double absprecision, int nIterations)
 
void LogOutputAtEndOfIntegration (double integral, double absprecision, double relprecision, int nIterations)
 
void Copy (const BCIntegrate &bcintegrate)
 
- Protected Attributes inherited from BCModel
std::string fName
 
double fModelAPriori
 
double fModelAPosteriori
 
BCDataSetfDataSet
 
BCDataPointfDataPointLowerBoundaries
 
BCDataPointfDataPointUpperBoundaries
 
std::vector< bool > fDataFixedValues
 
double fPValue
 
double fChi2NDoF
 
double fPValueNDoF
 
bool flag_discrete
 
int fGoFNIterationsMax
 
int fGoFNIterationsRun
 
int fGoFNChains
 
std::vector< TNamed * > fPriorContainer
 
bool fPriorConstantAll
 
std::vector< bool > fPriorContainerConstant
 
std::vector< bool > fPriorContainerInterpolate
 

Detailed Description

Definition at line 24 of file BCMVCDataModel.h.

Constructor & Destructor Documentation

BCMVCDataModel::BCMVCDataModel ( BCMVCombination mvc)

Definition at line 26 of file BCMVCDataModel.cxx.

26  : BCModel("BCMVCDataModel")
27 {
29  mvc->GetParameter(0)->GetLowerLimit(),
30  mvc->GetParameter(0)->GetUpperLimit());
34 }
BCMVCDataModel::~BCMVCDataModel ( )

Definition at line 37 of file BCMVCDataModel.cxx.

38 {
39 }

Member Function Documentation

double BCMVCDataModel::Chi2 ( TVectorD  observables,
TVectorD  measurements 
)

Definition at line 263 of file BCMVCDataModel.cxx.

264 {
265  TVectorD prod1 = observables - measurements;
266  TVectorD prod2 = fInvCovarianceMatrix * prod1;
267  double chi2 = prod1 * prod2;
268 
269  return chi2;
270 }
double BCMVCDataModel::LogAPrioriProbability ( const std::vector< double > &  parameters)
virtual

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

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

Reimplemented from BCModel.

Definition at line 128 of file BCMVCDataModel.cxx.

129 {
130  double logprob = 0.;
131 
132  (void) parameters; // suppress compiler warning about unused parameters
133 
134  return logprob;
135 }
double BCMVCDataModel::LogLikelihood ( const std::vector< double > &  params)
virtual

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

Parameters
paramsA set of parameter values
Returns
Natural logarithm of the likelihood

Implements BCModel.

Definition at line 103 of file BCMVCDataModel.cxx.

104 {
105  double logprob = 0.;
106 
107  int nmeas = GetNParameters();
108 
109  // copy parameters into a vector
110  TVectorD observables(nmeas);
111  TVectorD measurements(nmeas);
112 
113  for (int i = 0; i < nmeas; ++i) {
114  observables[i] = fPars[fVectorObservable[i]];
115  measurements[i] = parameters[i];
116  }
117 
118  TVectorD prod1 = observables - measurements;
119  TVectorD prod2 = fInvCovarianceMatrix * prod1;
120  double prod = prod1 * prod2;
121 
122  logprob = -0.5 * prod - log(TMath::Power(2*TMath::Pi(), nmeas/2.) * sqrt(fDetCovariance));
123 
124  return logprob;
125 }
void BCMVCDataModel::MCMCUserIterationInterface ( )
virtual

Method executed for every iteration of the MCMC. User's code should be provided via overloading in the derived class

Reimplemented from BCEngineMCMC.

Definition at line 138 of file BCMVCDataModel.cxx.

139 {
140  // get number of chains
141  int nchains = MCMCGetNChains();
142 
143  // get number of parameters
144  int npar = GetNParameters();
145 
146  // loop over all chains and fill histogram
147  for (int i = 0; i < nchains; ++i) {
148  // get the current values of the parameters x and y. These are
149 
150  // copy parameters into a vector
151  TVectorD observables(npar);
152  TVectorD measurements(npar);
153 
154  for (int j = 0; j < npar; ++j) {
155  observables[j] = fPars[fVectorObservable[j]];
156  measurements[j] = fMCMCx.at(i * npar + j);
157  }
158 
159  double chi2 = Chi2(observables, measurements);
160 
161  // fill the histogram
162  if (fHistChi2)
163  fHistChi2->Fill(chi2);
164  }
165 }
void BCMVCDataModel::PrintSummary ( )

Definition at line 242 of file BCMVCDataModel.cxx.

242  {
243 
244  std::cout << " Goodness-of-fit test:" << std::endl << std::endl;
245 
246  // calculate observed chi2
248 
249  // calculate expected chi2 distribution
250  BCH1D* hist_chi2 = new BCH1D(fHistChi2);
251  hist_chi2->GetHistogram()->Scale(1.0/hist_chi2->GetHistogram()->Integral("width"));
252 
253  // calculate p-value
254  double pvalue = hist_chi2->GetHistogram()->Integral(hist_chi2->GetHistogram()->FindBin(chi2), hist_chi2->GetHistogram()->GetNbinsX(), "width");
255 
256  std::cout << " chi2 : " << chi2 << std::endl;
257  std::cout << " p-value : " << pvalue << std::endl;
258  std::cout << std::endl;
259 
260 }
void BCMVCDataModel::PrintToys ( std::string  filename)

Definition at line 168 of file BCMVCDataModel.cxx.

169 {
170  // check if file extension is pdf
171  if ( (filename.find_last_of(".") != std::string::npos) &&
172  (filename.substr(filename.find_last_of(".")+1) == "pdf") ) {
173  ; // it's a PDF file
174 
175  }
176  else if ( (filename.find_last_of(".") != std::string::npos) &&
177  (filename.substr(filename.find_last_of(".")+1) == "ps") ) {
178  ; // it's a PS file
179  }
180  else {
181  ; // make it a PDF file
182  filename += ".pdf";
183  }
184 
185  int npars = GetNParameters();
186 
187  TCanvas* c1 = new TCanvas("");
188  c1->cd();
189 
190  // calculate observed chi2
192 
193  // calculate expected chi2 distribution
194  BCH1D* hist_chi2 = new BCH1D(fHistChi2);
195  if (hist_chi2->GetHistogram()->Integral("width") > 0)
196  hist_chi2->GetHistogram()->Scale(1.0/hist_chi2->GetHistogram()->Integral("width"));
197  else {
198  BCLog::OutWarning("BCMVCDataModel::PrintToys. Chi2 histogram not normalized. Abort printing.");
199  return;
200  }
201 
202  // calculate p-value
203  double pvalue = hist_chi2->GetHistogram()->Integral(hist_chi2->GetHistogram()->FindBin(chi2), hist_chi2->GetHistogram()->GetNbinsX(), "width");
204 
205  // draw expected chi2 distribution
206  hist_chi2->Draw("BTllB1CS1L", 1-pvalue);
207 
208  c1->Print(std::string(filename+"(").c_str());
209 
210  // draw all 1D histograms indicating observed value
211  for (int i = 0; i < npars; ++i) {
212  double obs = fVectorMeasurements[i];
213  const BCParameter* par = GetParameter(i);
214  BCH1D* hist_par = GetMarginalized(par);
215  double p = hist_par->GetHistogram()->Integral(hist_par->GetHistogram()->FindBin(obs), hist_par->GetHistogram()->GetNbinsX(), "width");
216  hist_par->Draw("BTllB1CS1L", 1-p);
217  c1->Print(filename.c_str());
218  }
219 
220  // draw all 2D histograms indicating observed values
221  for (int i = 0; i < npars; ++i) {
222  for (int j = 0; j < i; ++j) {
223  double obs_i = fVectorMeasurements[i];
224  double obs_j = fVectorMeasurements[j];
225  const BCParameter* par_i = GetParameter(i);
226  const BCParameter* par_j = GetParameter(j);
227  BCH2D* hist_par = GetMarginalized(par_i, par_j);
228  hist_par->Draw("BTfB3CS1");
229  TMarker* m = new TMarker();
230  m->SetMarkerStyle(21);
231  m->DrawMarker(obs_j, obs_i);
232  if (i == npars - 1 && j == i-1)
233  c1->Print(std::string(filename+")").c_str());
234  else
235  c1->Print(filename.c_str());
236  }
237  }
238 
239 }
void BCMVCDataModel::SetCovarianceMatrix ( TMatrixD  matrix)

Definition at line 49 of file BCMVCDataModel.cxx.

50 {
51  fCovarianceMatrix.Clear();
52  fCovarianceMatrix.ResizeTo(matrix);
53  fCovarianceMatrix = matrix;
54 
55  fInvCovarianceMatrix.Clear();
58  fInvCovarianceMatrix.Invert();
59 
60  fDetCovariance = fCovarianceMatrix.Determinant();
61 }
void BCMVCDataModel::SetHistChi2 ( TH1D *  hist)
inline

Definition at line 63 of file BCMVCDataModel.h.

64  { fHistChi2 = hist; };
void BCMVCDataModel::SetMeasurementRanges ( const std::vector< double > &  min,
const std::vector< double > &  max 
)

Definition at line 78 of file BCMVCDataModel.cxx.

79 {
80  unsigned npar = GetNParameters();
81 
82  if ( (min.size() != npar) || (max.size() != npar) ) {
83  BCLog::OutWarning("BCMVCDataModel::SetnMeasurementRanges. Size of ranges does not fit the number of measurements.");
84  return;
85  }
86 
87  // set the parameter ranges
88  for (unsigned i = 0; i < npar; ++i) {
89  GetParameter(i)->SetLimits(min.at(i), max.at(i));
90  }
91 
92 }
void BCMVCDataModel::SetMeasurementRanges ( double  min,
double  max 
)

Definition at line 95 of file BCMVCDataModel.cxx.

96 {
97  std::vector<double> min_vec(GetNParameters(), min);
98  std::vector<double> max_vec(GetNParameters(), max);
99  SetMeasurementRanges(min_vec, max_vec);
100 }
void BCMVCDataModel::SetNMeasurements ( int  n,
double  min,
double  max 
)

Definition at line 42 of file BCMVCDataModel.cxx.

43 {
44  for (int i = 0; i < n; ++i)
45  AddParameter(Form("measurement_%i", i), min, max);
46 }
void BCMVCDataModel::SetParameters ( std::vector< double >  parameters)

Definition at line 64 of file BCMVCDataModel.cxx.

65 {
66  fPars = parameters;
67  int nmeas = GetNParameters();
68 
69  fVectorObservables.Clear();
70  fVectorObservables.ResizeTo(nmeas);
71 
72  for (int i = 0; i < nmeas; ++i) {
74  }
75 }
void BCMVCDataModel::SetVectorMeasurements ( TVectorD  measurements)
inline

Definition at line 40 of file BCMVCDataModel.h.

41  { fVectorMeasurements.Clear();
42  fVectorMeasurements.ResizeTo(measurements);
43  fVectorMeasurements=measurements; };
void BCMVCDataModel::SetVectorObservable ( std::vector< int >  vec)
inline

Definition at line 56 of file BCMVCDataModel.h.

57  { fVectorObservable = vec; };

Member Data Documentation

TMatrixD BCMVCDataModel::fCovarianceMatrix
private

Definition at line 93 of file BCMVCDataModel.h.

double BCMVCDataModel::fDetCovariance
private

Definition at line 99 of file BCMVCDataModel.h.

TH1D* BCMVCDataModel::fHistChi2
private

Definition at line 111 of file BCMVCDataModel.h.

TMatrixD BCMVCDataModel::fInvCovarianceMatrix
private

Definition at line 96 of file BCMVCDataModel.h.

std::vector<double> BCMVCDataModel::fPars
private

Definition at line 102 of file BCMVCDataModel.h.

TVectorD BCMVCDataModel::fVectorMeasurements
private

Definition at line 90 of file BCMVCDataModel.h.

std::vector<int> BCMVCDataModel::fVectorObservable
private

Definition at line 108 of file BCMVCDataModel.h.

TVectorD BCMVCDataModel::fVectorObservables
private

Definition at line 105 of file BCMVCDataModel.h.


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