Logo ROOT  
Reference Guide
RooJeffreysPrior.cxx
Go to the documentation of this file.
1/** \class RooJeffreysPrior
2\ingroup Roofit
3
4Implementation of Jeffrey's prior. This class estimates the fisher information matrix by generating
5a binned Asimov dataset from the supplied PDFs, fitting it, retrieving the covariance matrix and inverting
6it. It returns the square root of the determinant of this matrix.
7Numerical integration is used to normalise. Since each integration step requires fits to be run,
8evaluating complicated PDFs may take long.
9
10Check the tutorial rs302_JeffreysPriorDemo.C for a demonstration with a simple PDF.
11**/
12
13#include "RooJeffreysPrior.h"
14
15#include "RooAbsReal.h"
16#include "RooAbsPdf.h"
17#include "RooErrorHandler.h"
18#include "RooArgSet.h"
19#include "RooMsgService.h"
20#include "RooFitResult.h"
21#include "TMatrixDSym.h"
22#include "RooDataHist.h"
23#include "RooNumIntConfig.h"
24#include "RooRealVar.h"
25#include "RooHelpers.h"
26
27using namespace std;
28
30
31using namespace RooFit;
32
33////////////////////////////////////////////////////////////////////////////////
34/// Construct a new JeffreysPrior.
35/// \param[in] name Name of this object.
36/// \param[in] title Title (for plotting)
37/// \param[in] nominal The PDF to base Jeffrey's prior on.
38/// \param[in] paramSet Parameters of the PDF.
39/// \param[in] obsSet Observables of the PDF.
40
41RooJeffreysPrior::RooJeffreysPrior(const char* name, const char* title,
42 RooAbsPdf& nominal,
43 const RooArgList& paramSet,
44 const RooArgList& obsSet) :
45 RooAbsPdf(name, title),
46 _nominal("nominal","nominal",this, nominal, false, false),
47 _obsSet("!obsSet","Observables",this, false, false),
48 _paramSet("!paramSet","Parameters",this),
49 _cacheMgr(this, 1, true, false)
50{
51 for (const auto comp : obsSet) {
52 if (!dynamic_cast<RooAbsReal*>(comp)) {
53 coutE(InputArguments) << "RooJeffreysPrior::ctor(" << GetName() << ") ERROR: component " << comp->GetName()
54 << " in observable list is not of type RooAbsReal" << endl ;
56 }
57 _obsSet.add(*comp) ;
58 }
59
60 for (const auto comp : paramSet) {
61 if (!dynamic_cast<RooAbsReal*>(comp)) {
62 coutE(InputArguments) << "RooJeffreysPrior::ctor(" << GetName() << ") ERROR: component " << comp->GetName()
63 << " in parameter list is not of type RooAbsReal" << endl ;
65 }
66 _paramSet.add(*comp) ;
67 }
68
69 // use a different integrator by default.
70 if(paramSet.getSize()==1)
71 this->specialIntegratorConfig(true)->method1D().setLabel("RooAdaptiveGaussKronrodIntegrator1D") ;
72}
73
74////////////////////////////////////////////////////////////////////////////////
75/// Copy constructor
76
78 RooAbsPdf(other, name),
79 _nominal("!nominal",this,other._nominal),
80 _obsSet("!obsSet",this,other._obsSet),
81 _paramSet("!paramSet",this,other._paramSet),
82 _cacheMgr(this, 1, true, false)
83{
84
85}
86
87////////////////////////////////////////////////////////////////////////////////
88/// Destructor
89
91{
92
93}
94
95////////////////////////////////////////////////////////////////////////////////
96/// Calculate and return current value of self
97
99{
101
102
103 CacheElem* cacheElm = (CacheElem*) _cacheMgr.getObj(nullptr);
104 if (!cacheElm) {
105 //Internally, we have to enlarge the range of fit parameters to make
106 //fits converge even if we are close to the limit of a parameter. Therefore, we clone the pdf and its
107 //observables here. If something happens to the external PDF, the cache is wiped,
108 //and we start to clone again.
109 auto& pdf = _nominal.arg();
110 RooAbsPdf* clonePdf = static_cast<RooAbsPdf*>(pdf.cloneTree());
111 auto vars = clonePdf->getParameters(_obsSet);
112 for (auto varTmp : *vars) {
113 auto& var = static_cast<RooRealVar&>(*varTmp);
114 auto range = var.getRange();
115 double span = range.second - range.first;
116 var.setRange(range.first - 0.1*span, range.second + 0.1 * span);
117 }
118
119 cacheElm = new CacheElem;
120 cacheElm->_pdf.reset(clonePdf);
121 cacheElm->_pdfVariables.reset(vars);
122
123 _cacheMgr.setObj(nullptr, cacheElm);
124 }
125
126 auto& cachedPdf = *cacheElm->_pdf;
127 auto& pdfVars = *cacheElm->_pdfVariables;
128 pdfVars.assign(_paramSet);
129
130 std::unique_ptr<RooDataHist> data( cachedPdf.generateBinned(_obsSet,ExpectedData()) );
131 std::unique_ptr<RooFitResult> res( cachedPdf.fitTo(*data, Save(),PrintLevel(-1),Minos(false),SumW2Error(false)) );
132 TMatrixDSym cov = res->covarianceMatrix();
133 cov.Invert();
134
135 return sqrt(cov.Determinant());
136}
137
#define coutE(a)
Definition: RooMsgService.h:37
#define ClassImp(name)
Definition: Rtypes.h:375
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
char name[80]
Definition: TGX11.cxx:110
RooArgSet * getParameters(const RooAbsData *data, bool stripDisconnected=true) const
Create a list of leaf nodes in the arg tree starting with ourself as top node that don't match any of...
Definition: RooAbsArg.cxx:565
Int_t getSize() const
Return the number of elements in the collection.
friend class CacheElem
The cache manager.
Definition: RooAbsPdf.h:367
std::pair< double, double > getRange(const char *name=0) const
Get low and high bound of the variable.
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:64
RooNumIntConfig * specialIntegratorConfig() const
Returns the specialized integrator configuration for this RooAbsReal.
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgList.h:22
T * getObj(const RooArgSet *nset, Int_t *sterileIndex=0, const TNamed *isetRangeName=0)
Getter function without integration set.
Int_t setObj(const RooArgSet *nset, T *obj, const TNamed *isetRangeName=0)
Setter function without integration set.
bool setLabel(const char *label, bool printError=true) override
Set value by specifying the name of the desired state.
bool add(const RooAbsArg &var, bool valueServer, bool shapeServer, bool silent)
Overloaded RooCollection_t::add() method insert object into set and registers object as server to own...
static void softAbort()
Soft abort function that interrupts macro execution but doesn't kill ROOT.
Switches the message service to a different level while the instance is alive.
Definition: RooHelpers.h:42
Implementation of Jeffrey's prior.
RooTemplateProxy< RooAbsPdf > _nominal
~RooJeffreysPrior() override
Destructor.
double evaluate() const override
Calculate and return current value of self.
RooObjCacheManager _cacheMgr
RooListProxy _obsSet
RooListProxy _paramSet
RooCategory & method1D()
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:40
const T & arg() const
Return reference to object held in proxy.
TMatrixTSym< Element > & Invert(Double_t *det=0)
Invert the matrix and calculate its determinant Notice that the LU decomposition is used instead of B...
Double_t Determinant() const override
const char * GetName() const override
Returns name of object.
Definition: TNamed.h:47
RooCmdArg Save(bool flag=true)
RooCmdArg SumW2Error(bool flag)
RooCmdArg Minos(bool flag=true)
RooCmdArg PrintLevel(Int_t code)
RooCmdArg ExpectedData(bool flag=true)
VecExpr< UnaryOp< Sqrt< T >, VecExpr< A, T, D >, T >, T, D > sqrt(const VecExpr< A, T, D > &rhs)
The namespace RooFit contains mostly switches that change the behaviour of functions of PDFs (or othe...
Definition: Common.h:18
@ InputArguments
Definition: RooGlobalFunc.h:64
std::unique_ptr< RooAbsPdf > _pdf
std::unique_ptr< RooArgSet > _pdfVariables