Logo ROOT   6.10/09
Reference Guide
RooJeffreysPrior.cxx
Go to the documentation of this file.
1 /** \class RooJeffreysPrior
2 \ingroup Roofit
3 
4 RooJeffreysPrior
5 **/
6 
7 
8 #include "RooFit.h"
9 
10 #include "Riostream.h"
11 #include "Riostream.h"
12 #include <math.h>
13 
14 #include "RooJeffreysPrior.h"
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 "RooFitResult.h"
24 #include "RooNumIntConfig.h"
25 #include "RooRealVar.h"
26 
27 #include "TError.h"
28 
29 using namespace std;
30 
32 
33 using namespace RooFit;
34 
35 ////////////////////////////////////////////////////////////////////////////////
36 ///_obsSet("!obsSet","obs-side variation",this),
37 
38 RooJeffreysPrior::RooJeffreysPrior(const char* name, const char* title,
39  RooAbsPdf& nominal,
40  const RooArgList& paramSet,
41  const RooArgList& obsSet) :
42  RooAbsPdf(name, title),
43  _nominal("nominal","nominal",this,nominal,kFALSE,kFALSE),
44  _obsSet("!obsSet","obs-side variation",this,kFALSE,kFALSE),
45  _paramSet("!paramSet","high-side variation",this)
46 {
49 
50 
51  TIterator* inputIter1 = obsSet.createIterator() ;
52  RooAbsArg* comp ;
53  while((comp = (RooAbsArg*)inputIter1->Next())) {
54  if (!dynamic_cast<RooAbsReal*>(comp)) {
55  coutE(InputArguments) << "RooJeffreysPrior::ctor(" << GetName() << ") ERROR: component " << comp->GetName()
56  << " in first list is not of type RooAbsReal" << endl ;
58  }
59  _obsSet.add(*comp) ;
60  // if (takeOwnership) {
61  // _ownedList.addOwned(*comp) ;
62  // }
63  }
64  delete inputIter1 ;
65 
66  TIterator* inputIter3 = paramSet.createIterator() ;
67  while((comp = (RooAbsArg*)inputIter3->Next())) {
68  if (!dynamic_cast<RooAbsReal*>(comp)) {
69  coutE(InputArguments) << "RooJeffreysPrior::ctor(" << GetName() << ") ERROR: component " << comp->GetName()
70  << " in first list is not of type RooAbsReal" << endl ;
72  }
73  _paramSet.add(*comp) ;
74  // if (takeOwnership) {
75  // _ownedList.addOwned(*comp) ;
76  // }
77  }
78  delete inputIter3 ;
79 
80  // use a different integrator by default.
81  if(paramSet.getSize()==1)
82  this->specialIntegratorConfig(kTRUE)->method1D().setLabel("RooAdaptiveGaussKronrodIntegrator1D") ;
83 
84 }
85 
86 ////////////////////////////////////////////////////////////////////////////////
87 /// Copy constructor
88 
90  RooAbsPdf(other, name),
91  _nominal("!nominal",this,other._nominal),
92  _obsSet("!obsSet",this,other._obsSet),
93  _paramSet("!paramSet",this,other._paramSet)
94 {
97 
98  // Member _ownedList is intentionally not copy-constructed -- ownership is not transferred
99 }
100 
101 ////////////////////////////////////////////////////////////////////////////////
102 /// Default constructor
103 
105 {
106  _obsIter = NULL;
107  _paramIter = NULL;
108 
109 }
110 
111 ////////////////////////////////////////////////////////////////////////////////
112 /// Destructor
113 
115 {
116  if (_obsIter) delete _obsIter ;
117  if (_paramIter) delete _paramIter ;
118 }
119 
120 ////////////////////////////////////////////////////////////////////////////////
121 /// Calculate and return current value of self
122 
124 {
127  // create Asimov dataset
128  // _paramSet.Print("v");
129  // return sqrt(_paramSet.getRealValue("mu"));
131  /*
132  cout << "_______________"<<endl;
133  _paramSet.Print("v");
134  _nominal->getVariables()->Print("v");
135  cout << "_______________"<<endl;
136  */
137  RooDataHist* data = ((RooAbsPdf&)(_nominal.arg())).generateBinned(_obsSet,ExpectedData());
138  // data->Print("v");
139  //RooFitResult* res = _nominal->fitTo(*data, Save(),PrintLevel(-1),Minos(kFALSE),SumW2Error(kTRUE));
140  RooFitResult* res = ((RooAbsPdf&)(_nominal.arg())).fitTo(*data, Save(),PrintLevel(-1),Minos(kFALSE),SumW2Error(kFALSE));
141  TMatrixDSym cov = res->covarianceMatrix();
142  cov.Invert();
143  double ret = sqrt(cov.Determinant());
144 
145  /*
146  // for 1 parameter can avoid making TMatrix etc.
147  // but number of params may be > 1 with others held constant
148  if(_paramSet.getSize()==1){
149  RooRealVar* var = (RooRealVar*) _paramSet.first();
150  // also, the _paramSet proxy one does not pick up a different value
151  cout << "eval at "<< ret << " " << 1/(var->getError()) << endl;
152  // need to get the actual variable instance out of the pdf like below
153  var = (RooRealVar*) _nominal->getVariables()->find(var->GetName());
154  cout << "eval at "<< ret << " " << 1/(var->getError()) << endl;
155  }
156  */
157 
158  // res->Print();
159  delete data;
160  delete res;
162 
163  // cout << "eval at "<< ret << endl;
164  // _paramSet.Print("v");
165  return ret;
166 
167 }
168 
169 ////////////////////////////////////////////////////////////////////////////////
170 /// if (matchArgs(allVars,analVars,x)) return 1 ;
171 /// if (matchArgs(allVars,analVars,mean)) return 2 ;
172 /// return 1;
173 
174 Int_t RooJeffreysPrior::getAnalyticalIntegral(RooArgSet& /*allVars*/, RooArgSet& /*analVars*/, const char* /*rangeName*/) const
175 {
176  return 0 ;
177 }
178 
179 ////////////////////////////////////////////////////////////////////////////////
180 
181 Double_t RooJeffreysPrior::analyticalIntegral(Int_t code, const char* /*rangeName*/) const
182 {
183  R__ASSERT(code==1 );
184  //cout << "evaluating analytic integral" << endl;
185  return 1.;
186 }
RooArgSet * getVariables(Bool_t stripDisconnected=kTRUE) const
Return RooArgSet with all variables (tree leaf nodes of expresssion tree)
Definition: RooAbsArg.cxx:2082
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
TIterator * createIterator(Bool_t dir=kIterForward) const
#define coutE(a)
Definition: RooMsgService.h:34
RooListProxy _obsSet
static void softAbort()
RooCmdArg PrintLevel(Int_t code)
RooCmdArg Minos(Bool_t flag=kTRUE)
const TMatrixDSym & covarianceMatrix() const
Return covariance matrix.
RooJeffreysPrior()
Default constructor.
RooFit::MsgLevel globalKillBelow() const
#define R__ASSERT(e)
Definition: TError.h:96
int Int_t
Definition: RtypesCore.h:41
static RooMsgService & instance()
Return reference to singleton instance.
STL namespace.
#define NULL
Definition: RtypesCore.h:88
Double_t analyticalIntegral(Int_t code, const char *rangeName=0) const
Implements the actual analytical integral(s) advertised by getAnalyticalIntegral. ...
Iterator abstract base class.
Definition: TIterator.h:30
virtual Bool_t setLabel(const char *label, Bool_t printError=kTRUE)
Set value by specifying the name of the desired state If printError is set, a message will be printed...
double sqrt(double)
RooCategory & method1D()
RooDataSet is a container class to hold N-dimensional binned data.
Definition: RooDataHist.h:40
virtual ~RooJeffreysPrior()
Destructor.
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Reimplementation of standard RooArgList::add()
RooCmdArg ExpectedData(Bool_t flag=kTRUE)
RooNumIntConfig * specialIntegratorConfig() const
Returns the specialized integrator configuration for this RooAbsReal.
Int_t getSize() const
virtual Double_t Determinant() const
void setGlobalKillBelow(RooFit::MsgLevel level)
TIterator * _paramIter
RooCmdArg SumW2Error(Bool_t flag)
const Bool_t kFALSE
Definition: RtypesCore.h:92
#define ClassImp(name)
Definition: Rtypes.h:336
double Double_t
Definition: RtypesCore.h:55
Double_t evaluate() const
Iterator over lowSet.
TMatrixTSym< Element > & Invert(Double_t *det=0)
Invert the matrix and calculate its determinant Notice that the LU decomposition is used instead of B...
RooListProxy _paramSet
RooCmdArg Save(Bool_t flag=kTRUE)
RooJeffreysPrior.
RooAbsPdf is the abstract interface for all probability density functions The class provides hybrid a...
Definition: RooAbsPdf.h:41
Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const
if (matchArgs(allVars,analVars,x)) return 1 ; if (matchArgs(allVars,analVars,mean)) return 2 ; return...
virtual TObject * Next()=0
RooAbsArg is the common abstract base class for objects that represent a value (of arbitrary type) an...
Definition: RooAbsArg.h:66
const RooAbsReal & arg() const
Definition: RooRealProxy.h:43
const Bool_t kTRUE
Definition: RtypesCore.h:91
RooRealProxy _nominal
TIterator * _obsIter
Iterator over paramSet.