54 #include "RConfigure.h"
65 #ifdef R__HAS_MATHMORE
74 struct LikelihoodFunction {
76 fFunc(f), fPrior(prior),
77 fOffset(offset), fMaxL(0) {
78 fFunc.binding().resetNumCall();
81 void SetPrior(
RooFunctor * prior) { fPrior = prior; }
84 double nll = fFunc(x) - fOffset;
87 if (fPrior) likelihood *= (*fPrior)(
x);
89 int nCalls = fFunc.binding().numCall();
90 if (nCalls > 0 && nCalls % 1000 == 0) {
92 <<
" x0 " << x[0] <<
" nll = " << nll+fOffset;
95 <<
" max Likelihood " << fMaxL << std::endl;
98 if (likelihood > fMaxL ) {
100 if ( likelihood > 1.E10) {
101 ooccoutW((
TObject*)0,
Eval) <<
"LikelihoodFunction::() WARNING - Huge likelihood value found for parameters ";
102 for (
int i = 0; i < fFunc.nObs(); ++i)
114 assert(fFunc.nObs() == 1);
116 return (*
this)(&tmp);
122 mutable double fMaxL;
137 fLikelihood(fFunctor, 0, nllMinimum),
139 fXmin(bindParams.getSize() ),
140 fXmax(bindParams.getSize() ),
141 fNorm(1.0), fNormErr(0.0), fOffset(0), fMaxPOI(0),
142 fHasNorm(
false), fUseOldValues(true), fError(
false)
147 fLikelihood.SetPrior(fPriorFunc.get() );
150 fIntegrator.SetFunction(fLikelihood, bindParams.
getSize() );
153 <<
" nllMinimum is " << nllMinimum << std::endl;
155 std::vector<double>
par(bindParams.
getSize());
156 for (
unsigned int i = 0; i < fXmin.size(); ++i) {
162 <<
" in interval [ " << fXmin[i] <<
" , " << fXmax[i] <<
" ] " << std::endl;
171 fNorm = (*this)( fMaxPOI );
174 fNormCdfValues.insert(std::make_pair(fXmin[0], 0) );
175 fNormCdfValues.insert(std::make_pair(fXmax[0], 1.0) );
182 PosteriorCdfFunction(
const PosteriorCdfFunction & rhs) :
184 fFunctor(rhs.fFunctor),
186 fPriorFunc(rhs.fPriorFunc),
187 fLikelihood(fFunctor, fPriorFunc.get(), rhs.fLikelihood.fOffset),
192 fNormErr( rhs.fNormErr),
193 fOffset(rhs.fOffset),
194 fMaxPOI(rhs.fMaxPOI),
195 fHasNorm(rhs.fHasNorm),
196 fUseOldValues(rhs.fUseOldValues),
198 fNormCdfValues(rhs.fNormCdfValues)
200 fIntegrator.SetFunction(fLikelihood, fXmin.size() );
209 bool HasError()
const {
return fError; }
214 return new PosteriorCdfFunction(*
this);
218 void SetOffset(
double offset) { fOffset = offset; }
223 PosteriorCdfFunction&
operator=(
const PosteriorCdfFunction &) {
227 double DoEval (
double x)
const {
231 if (x <= fXmin[0] )
return -fOffset;
233 if (x >= fMaxPOI && fHasNorm)
return 1. - fOffset;
237 if (fHasNorm && fUseOldValues) {
239 std::map<double,double>::iterator itr = fNormCdfValues.upper_bound(x);
241 if (itr != fNormCdfValues.end() ) {
242 fXmin[0] = itr->first;
243 normcdf0 = itr->second;
249 fFunctor.binding().resetNumCall();
251 double cdf = fIntegrator.Integral(&fXmin[0],&fXmax[0]);
252 double error = fIntegrator.Error();
253 double normcdf = cdf/fNorm;
256 << fXmax[0] <<
"] integral = " << cdf <<
" +/- " << error
257 <<
" norm-integ = " << normcdf <<
" cdf(x) = " << normcdf+normcdf0
258 <<
" ncalls = " << fFunctor.binding().numCall() << std::endl;
266 if (cdf != 0 && error/cdf > 0.2 )
268 <<
" x = " << x <<
" cdf(x) = " << cdf <<
" +/- " << error << std::endl;
272 << cdf <<
" +/- " << error << std::endl;
281 fNormCdfValues.insert(std::make_pair(x, normcdf) );
284 double errnorm =
sqrt( error*error + normcdf*normcdf * fNormErr * fNormErr )/fNorm;
285 if (normcdf > 1. + 3 * errnorm) {
287 <<
" x = " << x <<
" normcdf(x) = " << normcdf <<
" +/- " << error/fNorm << std::endl;
290 return normcdf - fOffset;
294 mutable std::shared_ptr<RooFunctor> fPriorFunc;
295 LikelihoodFunction fLikelihood;
297 mutable std::vector<double> fXmin;
298 mutable std::vector<double> fXmax;
300 mutable double fNormErr;
306 mutable std::map<double,double> fNormCdfValues;
320 norm = 1.0,
double nllOffset = 0,
int niter = 0) :
323 fLikelihood(fFunctor, 0, nllOffset),
325 fXmin(nuisParams.getSize() ),
326 fXmax(nuisParams.getSize() ),
333 fLikelihood.SetPrior(fPriorFunc.get() );
337 for (
unsigned int i = 0; i < fXmin.size(); ++i) {
342 <<
" in interval [" << fXmin[i] <<
" , " << fXmax[i] <<
" ] " << std::endl;
344 if (fXmin.size() == 1) {
347 fIntegratorOneDim->SetFunction(fLikelihood);
352 else if (fXmin.size() > 1) {
354 fIntegratorMultiDim->SetFunction(fLikelihood, fXmin.size());
358 fIntegratorMultiDim->SetOptions(opt);
372 double Error()
const {
return fError;}
376 double DoEval (
double x)
const {
381 fFunctor.binding().resetNumCall();
385 if (fXmin.size() == 1) {
386 f = fIntegratorOneDim->Integral(fXmin[0],fXmax[0]);
387 error = fIntegratorOneDim->Error();
389 else if (fXmin.size() > 1) {
390 f = fIntegratorMultiDim->Integral(&fXmin[0],&fXmax[0]);
391 error = fIntegratorMultiDim->Error();
399 << x <<
"\tf(x) = " << f <<
" +/- " << error
400 <<
" norm-f(x) = " << f/fNorm
401 <<
" ncalls = " << fFunctor.binding().numCall() << std::endl;
406 if (f != 0 && error/f > 0.2 )
408 << fXmin.size() <<
" Dim is larger than 20 % "
409 <<
"x = " << x <<
" p(x) = " << f <<
" +/- " << error << std::endl;
411 fError = error / fNorm;
416 mutable std::shared_ptr<RooFunctor> fPriorFunc;
417 LikelihoodFunction fLikelihood;
419 std::unique_ptr<ROOT::Math::Integrator> fIntegratorOneDim;
420 std::unique_ptr<ROOT::Math::IntegratorMultiDim> fIntegratorMultiDim;
421 std::vector<double> fXmin;
422 std::vector<double> fXmax;
424 mutable double fError;
435 nllOffset = 0,
int niter = 0,
bool redoToys =
true ) :
438 fLikelihood(fFunctor, 0, nllOffset),
441 fNuisParams(nuisParams),
443 fNumIterations(niter),
447 if (niter == 0) fNumIterations = 100;
451 fLikelihood.SetPrior(fPriorFunc.get() );
454 ooccoutI((
TObject*)0,
InputArguments) <<
"PosteriorFunctionFromToyMC::Evaluate the posterior function by randomizing the nuisances: niter " << fNumIterations << std::endl;
456 ooccoutI((
TObject*)0,
InputArguments) <<
"PosteriorFunctionFromToyMC::Pdf used for randomizing the nuisance is " << fPdf->GetName() << std::endl;
459 for (
int i = 0; i < fNuisParams.getSize(); ++i) {
460 if (!vars->
find( fNuisParams[i].GetName() ) ) {
462 <<
" is not part of sampling pdf. "
463 <<
"they will be trated as constant " << std::endl;
469 ooccoutI((
TObject*)0,
InputArguments) <<
"PosteriorFunctionFromToyMC::Generate nuisance toys only one time (for all POI points)" << std::endl;
474 virtual ~PosteriorFunctionFromToyMC() {
if (fGenParams)
delete fGenParams; }
477 void GenerateToys()
const {
478 if (fGenParams)
delete fGenParams;
479 fGenParams = fPdf->generate(fNuisParams, fNumIterations);
485 double Error()
const {
return fError;}
497 double DoEval(
double x)
const {
499 int npar = fNuisParams.getSize();
504 if (fRedoToys) GenerateToys();
505 if (!fGenParams)
return 0;
518 std::vector<double> p(npar);
519 for (
int i = 0; i < npar; ++i) {
529 double fval = fLikelihood( &p.front() );
535 double nuisPdfVal = fPdf->getVal(&arg);
541 <<
"Likelihood evaluates to infinity " << std::endl;
544 for (
int i = 0; i < npar; ++i)
553 <<
"Likelihood is a NaN " << std::endl;
556 for (
int i = 0; i < npar; ++i)
570 double val = sum/
double(fNumIterations);
571 double dval2 =
std::max( sum2/
double(fNumIterations) - val*val, 0.0);
572 fError =
std::sqrt( dval2 / fNumIterations);
576 << x <<
"\tp(x) = " << val <<
" +/- " << fError << std::endl;
579 if (val != 0 && fError/val > 0.2 ) {
581 <<
" - Error in estimating posterior is larger than 20% ! "
582 <<
"x = " << x <<
" p(x) = " << val <<
" +/- " << fError << std::endl;
590 mutable std::shared_ptr<RooFunctor> fPriorFunc;
591 LikelihoodFunction fLikelihood;
597 mutable double fError;
606 BayesianCalculator::BayesianCalculator() :
611 fProductPdf (0), fLogLike(0), fLikelihood (0), fIntegratedLikelihood (0), fPosteriorPdf(0),
612 fPosteriorFunction(0), fApproxPosterior(0),
613 fLower(0), fUpper(0),
615 fSize(0.05), fLeftSideFraction(0.5),
616 fBrfPrecision(0.00005),
619 fValidInterval(
false)
633 fPriorPdf(&priorPdf),
635 fProductPdf (0), fLogLike(0), fLikelihood (0), fIntegratedLikelihood (0), fPosteriorPdf(0),
636 fPosteriorFunction(0), fApproxPosterior(0),
637 fLower(0), fUpper(0),
639 fSize(0.05), fLeftSideFraction(0.5),
640 fBrfPrecision(0.00005),
643 fValidInterval(
false)
657 fPdf(model.GetPdf()),
658 fPriorPdf( model.GetPriorPdf()),
660 fProductPdf (0), fLogLike(0), fLikelihood (0), fIntegratedLikelihood (0), fPosteriorPdf(0),
661 fPosteriorFunction(0), fApproxPosterior(0),
662 fLower(0), fUpper(0),
664 fSize(0.05), fLeftSideFraction(0.5),
665 fBrfPrecision(0.00005),
668 fValidInterval(
false)
741 coutE(
InputArguments) <<
"BayesianCalculator::GetPosteriorPdf - missing pdf model" << std::endl;
745 coutE(
InputArguments) <<
"BayesianCalculator::GetPosteriorPdf - missing parameter of interest" << std::endl;
749 coutE(
InputArguments) <<
"BayesianCalculator::GetPosteriorPdf - current implementation works only on 1D intervals" << std::endl;
765 ccoutD(
Eval) <<
"BayesianCalculator::GetPosteriorFunction : "
775 coutE(
Eval) <<
"BayesianCalculator::GetPosteriorFunction : "
776 <<
" Negative log likelihood evaluates to infinity " << std::endl
777 <<
" Non-const Parameter values : ";
779 for (
int i = 0; i < p.
getSize(); ++i) {
784 ccoutE(
Eval) <<
"-- Perform a full likelihood fit of the model before or set more reasanable parameter values"
786 coutE(
Eval) <<
"BayesianCalculator::GetPosteriorFunction : " <<
" cannot compute posterior function " << std::endl;
806 coutI(
Eval) <<
"BayesianCalculator::GetPosteriorFunction : "
807 <<
" nll value " << nllVal <<
" poi value = " << poi->
getVal() << std::endl;
816 coutI(
Eval) <<
"BayesianCalculator::GetPosteriorFunction : minimum of NLL vs POI for POI = "
821 delete constrainedParams;
826 ccoutD(
Eval) <<
"BayesianCalculator::GetPosteriorFunction : use ROOFIT integration "
829 #ifdef DOLATER // (not clear why this does not work)
884 bool doToysEveryIteration =
true;
890 ccoutI(
Eval) <<
"BayesianCalculator::GetPosteriorFunction : no nuisance pdf is provided, try using global pdf (this will be slower)"
894 fNumIterations, doToysEveryIteration );
936 if (!plike)
return 0;
966 if (!posterior)
return 0;
1045 coutW(
Eval) <<
"BayesianCalculator::GetInterval - recomputing interval for the same CL and same model" << std::endl;
1049 coutE(
Eval) <<
"BayesianCalculator::GetInterval - no parameter of interest is set " << std::endl;
1088 coutW(
Eval) <<
"BayesianCalculator::GetInterval - computing integral from cdf failed - do a scan in "
1106 coutE(
Eval) <<
"BayesianCalculator::GetInterval - cannot compute a valid interval - return a dummy [1,0] interval"
1110 coutI(
Eval) <<
"BayesianCalculator::GetInterval - found a valid interval : [" <<
fLower <<
" , "
1111 <<
fUpper <<
" ]" << std::endl;
1116 interval->
SetTitle(
"SimpleInterval from BayesianCalculator");
1137 coutI(
Eval) <<
"BayesianCalculator: Compute interval using RooFit: posteriorPdf + createCdf + RooBrentRootFinder " << std::endl;
1150 if (!cdf_bind)
return;
1155 double tmpVal = poi->
getVal();
1158 if (lowerCutOff > 0) {
1159 double y = lowerCutOff;
1165 if (upperCutOff < 1.0) {
1166 double y=upperCutOff;
1171 if (!ret)
coutE(
Eval) <<
"BayesianCalculator::GetInterval "
1172 <<
"Error returned from Root finder, estimated interval is not fully correct"
1189 coutI(
InputArguments) <<
"BayesianCalculator:GetInterval Compute the interval from the posterior cdf " << std::endl;
1194 coutE(
InputArguments) <<
"BayesianCalculator::GetInterval() cannot make posterior Function " << std::endl;
1208 if( cdf.HasError() ) {
1209 coutE(
Eval) <<
"BayesianCalculator: Numerical error computing CDF integral - try a different method " << std::endl;
1217 ccoutD(
Eval) <<
"BayesianCalculator::GetInterval - finding roots of posterior using RF " << rf.
Name()
1220 if (lowerCutOff > 0) {
1221 cdf.SetOffset(lowerCutOff);
1222 ccoutD(
NumIntegration) <<
"Integrating posterior to get cdf and search lower limit at p =" << lowerCutOff << std::endl;
1224 if( cdf.HasError() )
1225 coutW(
Eval) <<
"BayesianCalculator: Numerical error integrating the CDF " << std::endl;
1227 coutE(
NumIntegration) <<
"BayesianCalculator::GetInterval - Error from root finder when searching lower limit !" << std::endl;
1235 if (upperCutOff < 1.0) {
1236 cdf.SetOffset(upperCutOff);
1237 ccoutD(
NumIntegration) <<
"Integrating posterior to get cdf and search upper interval limit at p =" << upperCutOff << std::endl;
1239 if( cdf.HasError() )
1240 coutW(
Eval) <<
"BayesianCalculator: Numerical error integrating the CDF " << std::endl;
1242 coutE(
NumIntegration) <<
"BayesianCalculator::GetInterval - Error from root finder when searching upper limit !" << std::endl;
1271 if (!posterior)
return;
1279 coutI(
Eval) <<
"BayesianCalculator - scan posterior function in nbins = " << tmp->
GetNpx() << std::endl;
1304 ccoutD(
Eval) <<
"BayesianCalculator: Compute interval from the approximate posterior " << std::endl;
1310 double limits[2] = {0,0};
1311 prob[0] = lowerCutOff;
1312 prob[1] = upperCutOff;
1321 coutI(
Eval) <<
"BayesianCalculator - computing shortest interval with CL = " << 1.-
fSize << std::endl;
1333 std::vector<int> index(n);
1337 double actualCL = 0;
1342 for (
int i = 0; i <
n; ++i) {
1344 double p = bins[ idx] /
norm;
1346 if (sum > 1.-
fSize ) {
1357 ccoutD(
Eval) <<
"BayesianCalculator::ComputeShortestInterval - actual interval CL = "
1358 << actualCL <<
" difference from requested is " << (actualCL-(1.-
fSize))/
fSize*100. <<
"% "
1359 <<
" limits are [ " << lower <<
" , " <<
" upper ] " << std::endl;
1362 if (lower < upper) {
1369 coutW(
Eval) <<
"BayesianCalculator::ComputeShortestInterval - actual interval CL = "
1370 << actualCL <<
" differs more than 10% from desired CL value - must increase nbins "
1371 << n <<
" to an higher value " << std::endl;
1374 coutE(
Eval) <<
"BayesianCalculator::ComputeShortestInterval " << n <<
" bins are not sufficient " << std::endl;
User Class to find the Root of one dimensional functions.
virtual RooAbsReal * createNLL(RooAbsData &data, const RooLinkedList &cmdList)
Construct representation of -log(L) of PDFwith given dataset.
RooPlot * GetPosteriorPlot(bool norm=false, double precision=0.01) const
get the plot with option to get it normalized
virtual const char * GetTitle() const
Returns title of object.
virtual Int_t GetQuantiles(Int_t nprobSum, Double_t *q, const Double_t *probSum)
Compute Quantiles for density distribution of this function.
virtual RooPlot * plotOn(RooPlot *frame, const RooCmdArg &arg1=RooCmdArg(), const RooCmdArg &arg2=RooCmdArg(), const RooCmdArg &arg3=RooCmdArg(), const RooCmdArg &arg4=RooCmdArg(), const RooCmdArg &arg5=RooCmdArg(), const RooCmdArg &arg6=RooCmdArg(), const RooCmdArg &arg7=RooCmdArg(), const RooCmdArg &arg8=RooCmdArg(), const RooCmdArg &arg9=RooCmdArg(), const RooCmdArg &arg10=RooCmdArg()) const
Plot (project) PDF on specified frame.
static IntegrationMultiDim::Type GetType(const char *name)
static function to get the enumeration from a string
RooAbsPdf * GetPosteriorPdf() const
return posterior pdf (object is managed by the BayesianCalculator class)
ROOT::Math::IGenFunction * fPosteriorFunction
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
RooCmdArg DrawOption(const char *opt)
ModelConfig is a simple class that holds configuration information specifying how a model should be u...
RooAbsPdf * GetPdf() const
get model PDF (return NULL if pdf has not been specified or does not exist)
Interface (abstract class) for generic functions objects of one-dimension Provides a method to evalua...
virtual void SetNpx(Int_t npx=100)
Set the number of points used to draw the function.
Namespace for new ROOT classes and functions.
virtual Int_t GetMaximumBin() const
Return location of bin with maximum value in the range.
const RooArgSet * GetConditionalObservables() const
get RooArgSet for conditional observables (return NULL if not existing)
void SetIntegrationType(const char *type)
set the integration type (possible type are) :
void ToUpper()
Change string to upper case.
static void clearEvalErrorLog()
Clear the stack of evaluation error messages.
virtual Double_t getMin(const char *name=0) const
void SetNCalls(unsigned int calls)
set maximum number of function calls
void SetTitle(const char *name)
Set the title of the RooPlot to 'title'.
Type GetType(const std::string &Name)
const TKDTreeBinning * bins
double Root() const
Return the current and latest estimate of the Root.
ClassImp(TIterator) Bool_t TIterator return false
Compare two iterator objects.
const RooArgSet * GetNuisanceParameters() const
get RooArgSet containing the nuisance parameters (return NULL if not existing)
virtual TH1 * GetHistogram() const
Return a pointer to the histogram used to vusualize the function.
VecExpr< UnaryOp< Sqrt< T >, VecExpr< A, T, D >, T >, T, D > sqrt(const VecExpr< A, T, D > &rhs)
double GetMode() const
return the mode (most probable value of the posterior function)
virtual SimpleInterval * GetInterval() const
compute the interval.
void ComputeShortestInterval() const
RooAbsReal * createCdf(const RooArgSet &iset, const RooArgSet &nset=RooArgSet())
Create a cumulative distribution function of this p.d.f in terms of the observables listed in iset...
RooAbsArg * first() const
RooAbsReal * GetPosteriorFunction() const
return posterior function (object is managed by the BayesianCalculator class)
RooCmdArg Range(const char *rangeName, Bool_t adjustNorm=kTRUE)
void ComputeIntervalUsingRooFit(double c1, double c2) const
static void setEvalErrorLoggingMode(ErrorLoggingMode m)
Set evaluation error logging mode.
virtual Double_t GetBinLowEdge(Int_t bin) const
return bin lower edge for 1D historam Better to use h1.GetXaxis().GetBinLowEdge(bin) ...
virtual void removeAll()
Remove all arguments from our set, deleting them if we own them.
double cdf(double *x, double *p)
void Print(std::ostream &os=std::cout) const
print all the options
std::map< std::string, std::string >::const_iterator iter
static Vc_ALWAYS_INLINE Vector< T > abs(const Vector< T > &x)
TF1 * asTF(const RooArgList &obs, const RooArgList &pars=RooArgList(), const RooArgSet &nset=RooArgSet()) const
Return a ROOT TF1,2,3 object bound to this RooAbsReal with given definition of observables and parame...
void Sort(Index n, const Element *a, Index *index, Bool_t down=kTRUE)
Double_t getVal(const RooArgSet *set=0) const
virtual TObject * Clone(const char *newname="") const
Make a clone of an object using the Streamer facility.
ClassImp(RooStats::BayesianCalculator) using namespace std
BayesianCalculator implementation.
User class for performing function minimization.
RooAbsReal * fIntegratedLikelihood
void ComputeIntervalFromApproxPosterior(double c1, double c2) const
RooArgSet * getParameters(const RooAbsData *data, Bool_t stripDisconnected=kTRUE) const
Create a list of leaf nodes in the arg tree starting with ourself as top node that don't match any of...
virtual void setVal(Double_t value)
Set value of variable to 'value'.
void ComputeIntervalFromCdf(double c1, double c2) const
virtual Double_t GetBinCenter(Int_t bin) const
return bin center for 1D historam Better to use h1.GetXaxis().GetBinCenter(bin)
Numerical multi dimensional integration options.
RooAbsArg * find(const char *name) const
Find object with given name in list.
static Int_t numEvalErrors()
Return the number of logged evaluation errors since the last clearing.
IBaseFunctionOneDim IGenFunction
virtual Double_t GetBinUpEdge(Int_t bin) const
Return up edge of bin.
static IntegrationOneDim::Type GetType(const char *name)
static function to get the enumeration from a string
void SetFunction(const ROOT::Math::IGenFunction &f, double xlow, double xup)
Sets function to be minimized.
void Form(const char *fmt,...)
Formats a string using a printf style format descriptor.
User Class for performing numerical integration of a function in one dimension.
virtual Double_t GetSumOfWeights() const
Return the sum of weights excluding under/overflows.
virtual const char * GetName() const
Returns name of object.
void setTol(Double_t tol)
virtual Double_t ConfidenceLevel() const
Get the Confidence level for the test.
RooPlot * frame(const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
1-D histogram with a double per channel (see TH1 documentation)}
const Double_t * GetArray() const
virtual void SetModel(const ModelConfig &model)
set the model via the ModelConfig
virtual void SetName(const char *name)
Change the name of this histogram.
virtual Int_t GetNpx() const
TRObject operator()(const T1 &t1) const
const char * Name() const
Return the current and latest estimate of the lower value of the Root-finding interval (for bracketin...
RooCmdArg Precision(Double_t prec)
RooCmdArg FillColor(Color_t color)
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
RooArgSet fConditionalObs
virtual double DoEval(double x) const =0
implementation of the evaluation function. Must be implemented by derived classes ...
static Vc_ALWAYS_INLINE int_v max(const int_v &x, const int_v &y)
Namespace for new Math classes and functions.
RooArgSet fNuisanceParameters
Binding & operator=(OUT(*fun)(void))
Mother of all ROOT objects.
virtual Double_t getMax(const char *name=0) const
void RemoveConstantParameters(RooArgSet *set)
virtual Bool_t findRoot(Double_t &result, Double_t xlo, Double_t xhi, Double_t value=0) const
Do the root finding using the Brent-Decker method.
RooAbsPdf * fPosteriorPdf
BayesianCalculator()
constructor
virtual double FValMinimum() const
Return function value at current estimate of the minimum.
RooAbsPdf is the abstract interface for all probability density functions The class provides hybrid a...
User class for performing multidimensional integration.
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
virtual IBaseFunctionOneDim * Clone() const =0
Clone a function.
virtual ~BayesianCalculator()
destructor
const RooArgSet * GetParametersOfInterest() const
get RooArgSet containing the parameter of interest (return NULL if not existing)
virtual bool Minimize(int maxIter, double absTol=1.E-8, double relTol=1.E-10)
Find minimum position iterating until convergence specified by the absolute and relative tolerance or...
Functor1D class for one-dimensional functions.
RooAbsArg is the common abstract base class for objects that represent a value (of arbitrary type) an...
RooAbsPdf * GetPriorPdf() const
get parameters prior pdf (return NULL if not existing)
RooFunctor * functor(const RooArgList &obs, const RooArgList &pars=RooArgList(), const RooArgSet &nset=RooArgSet()) const
Return a RooFunctor object bound to this RooAbsReal with given definition of observables and paramete...
void ApproximatePosterior() const
RooCmdArg ConditionalObservables(const RooArgSet &set)
virtual void SetTitle(const char *title="")
Change (i.e. set) the title of the TNamed.
double norm(double *x, double *p)
bool Solve(Function &f, Derivative &d, double start, int maxIter=100, double absTol=1E-8, double relTol=1E-10)
RooAbsFunc * bindVars(const RooArgSet &vars, const RooArgSet *nset=0, Bool_t clipInvalid=kFALSE) const
Create an interface adaptor f(vars) that binds us to the specified variables (in arbitrary order)...
BayesianCalculator is a concrete implementation of IntervalCalculator, providing the computation of a...
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add element to non-owning set.
RooCmdArg Constrain(const RooArgSet ¶ms)
void Error(ErrorHandler_t func, int code, const char *va_(fmt),...)
Write error message and call a handler, if required.
RooAbsReal * createIntegral(const RooArgSet &iset, const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
Create an object that represents the integral of the function over one or more observables listed in ...