Logo ROOT  
Reference Guide
RooMomentMorph.cxx
Go to the documentation of this file.
1/*****************************************************************************
2 * Project: RooFit *
3 * *
4 * This code was autogenerated by RooClassFactory *
5 *****************************************************************************/
6
7/** \class RooMomentMorph
8 \ingroup Roofit
9
10**/
11
12#include "Riostream.h"
13
14#include "RooMomentMorph.h"
15#include "RooAbsCategory.h"
16#include "RooRealIntegral.h"
17#include "RooRealConstant.h"
18#include "RooRealVar.h"
19#include "RooFormulaVar.h"
20#include "RooCustomizer.h"
21#include "RooAddPdf.h"
22#include "RooAddition.h"
23#include "RooMoment.h"
24#include "RooLinearVar.h"
25#include "RooChangeTracker.h"
26
27#include "TMath.h"
28#include "TH1.h"
29
30using namespace std;
31
33
34////////////////////////////////////////////////////////////////////////////////
35/// coverity[UNINIT_CTOR]
36
37RooMomentMorph::RooMomentMorph() : _curNormSet(0), _mref(0), _M(0), _useHorizMorph(true)
38{
41}
42
43////////////////////////////////////////////////////////////////////////////////
44/// CTOR
45
46RooMomentMorph::RooMomentMorph(const char *name, const char *title,
47 RooAbsReal& _m,
48 const RooArgList& varList,
49 const RooArgList& pdfList,
50 const TVectorD& mrefpoints,
51 Setting setting) :
52 RooAbsPdf(name,title),
53 _cacheMgr(this,10,kTRUE,kTRUE),
54 m("m","m",this,_m),
55 _varList("varList","List of variables",this),
56 _pdfList("pdfList","List of pdfs",this),
57 _setting(setting),
58 _useHorizMorph(true)
59{
60 // observables
61 TIterator* varItr = varList.createIterator() ;
62 RooAbsArg* var ;
63 for (Int_t i=0; (var = (RooAbsArg*)varItr->Next()); ++i) {
64 if (!dynamic_cast<RooAbsReal*>(var)) {
65 coutE(InputArguments) << "RooMomentMorph::ctor(" << GetName() << ") ERROR: variable " << var->GetName() << " is not of type RooAbsReal" << endl ;
66 throw string("RooPolyMorh::ctor() ERROR variable is not of type RooAbsReal") ;
67 }
68 _varList.add(*var) ;
69 }
70 delete varItr ;
71
72 // reference p.d.f.s
73 TIterator* pdfItr = pdfList.createIterator() ;
74 RooAbsPdf* pdf ;
75 for (Int_t i=0; (pdf = dynamic_cast<RooAbsPdf*>(pdfItr->Next())); ++i) {
76 if (!pdf) {
77 coutE(InputArguments) << "RooMomentMorph::ctor(" << GetName() << ") ERROR: pdf " << pdf->GetName() << " is not of type RooAbsPdf" << endl ;
78 throw string("RooPolyMorh::ctor() ERROR pdf is not of type RooAbsPdf") ;
79 }
80 _pdfList.add(*pdf) ;
81 }
82 delete pdfItr ;
83
84 _mref = new TVectorD(mrefpoints);
87
88 // initialization
89 initialize();
90}
91
92////////////////////////////////////////////////////////////////////////////////
93/// CTOR
94
95RooMomentMorph::RooMomentMorph(const char *name, const char *title,
96 RooAbsReal& _m,
97 const RooArgList& varList,
98 const RooArgList& pdfList,
99 const RooArgList& mrefList,
100 Setting setting) :
101 RooAbsPdf(name,title),
102 _cacheMgr(this,10,kTRUE,kTRUE),
103 m("m","m",this,_m),
104 _varList("varList","List of variables",this),
105 _pdfList("pdfList","List of pdfs",this),
106 _setting(setting),
107 _useHorizMorph(true)
108{
109 // observables
110 TIterator* varItr = varList.createIterator() ;
111 RooAbsArg* var ;
112 for (Int_t i=0; (var = (RooAbsArg*)varItr->Next()); ++i) {
113 if (!dynamic_cast<RooAbsReal*>(var)) {
114 coutE(InputArguments) << "RooMomentMorph::ctor(" << GetName() << ") ERROR: variable " << var->GetName() << " is not of type RooAbsReal" << endl ;
115 throw string("RooPolyMorh::ctor() ERROR variable is not of type RooAbsReal") ;
116 }
117 _varList.add(*var) ;
118 }
119 delete varItr ;
120
121 // reference p.d.f.s
122 TIterator* pdfItr = pdfList.createIterator() ;
123 RooAbsPdf* pdf ;
124 for (Int_t i=0; (pdf = dynamic_cast<RooAbsPdf*>(pdfItr->Next())); ++i) {
125 if (!pdf) {
126 coutE(InputArguments) << "RooMomentMorph::ctor(" << GetName() << ") ERROR: pdf " << pdf->GetName() << " is not of type RooAbsPdf" << endl ;
127 throw string("RooPolyMorh::ctor() ERROR pdf is not of type RooAbsPdf") ;
128 }
129 _pdfList.add(*pdf) ;
130 }
131 delete pdfItr ;
132
133 // reference points in m
134 _mref = new TVectorD(mrefList.getSize());
135 TIterator* mrefItr = mrefList.createIterator() ;
136 RooAbsReal* mref ;
137 for (Int_t i=0; (mref = dynamic_cast<RooAbsReal*>(mrefItr->Next())); ++i) {
138 if (!mref) {
139 coutE(InputArguments) << "RooMomentMorph::ctor(" << GetName() << ") ERROR: mref " << mref->GetName() << " is not of type RooAbsReal" << endl ;
140 throw string("RooPolyMorh::ctor() ERROR mref is not of type RooAbsReal") ;
141 }
142 if (!dynamic_cast<RooConstVar*>(mref)) {
143 coutW(InputArguments) << "RooMomentMorph::ctor(" << GetName() << ") WARNING mref point " << i << " is not a constant, taking a snapshot of its value" << endl ;
144 }
145 (*_mref)[i] = mref->getVal() ;
146 }
147 delete mrefItr ;
148
151
152 // initialization
153 initialize();
154}
155
156////////////////////////////////////////////////////////////////////////////////
157
159 RooAbsPdf(other,name),
160 _cacheMgr(other._cacheMgr,this),
161 _curNormSet(0),
162 m("m",this,other.m),
163 _varList("varList",this,other._varList),
164 _pdfList("pdfList",this,other._pdfList),
165 _setting(other._setting),
166 _useHorizMorph(other._useHorizMorph)
167{
168 _mref = new TVectorD(*other._mref) ;
171
172 // initialization
173 initialize();
174}
175
176////////////////////////////////////////////////////////////////////////////////
177
179{
180 if (_mref) delete _mref;
181 if (_varItr) delete _varItr;
182 if (_pdfItr) delete _pdfItr;
183 if (_M) delete _M;
184}
185
186////////////////////////////////////////////////////////////////////////////////
187
189{
190 Int_t nPdf = _pdfList.getSize();
191
192 // other quantities needed
193 if (nPdf!=_mref->GetNrows()) {
194 coutE(InputArguments) << "RooMomentMorph::initialize(" << GetName() << ") ERROR: nPdf != nRefPoints" << endl ;
195 assert(0) ;
196 }
197
198 TVectorD* dm = new TVectorD(nPdf);
199 _M = new TMatrixD(nPdf,nPdf);
200
201 // transformation matrix for non-linear extrapolation, needed in evaluate()
202 TMatrixD M(nPdf,nPdf);
203 for (Int_t i=0; i<_mref->GetNrows(); ++i) {
204 (*dm)[i] = (*_mref)[i]-(*_mref)[0];
205 M(i,0) = 1.;
206 if (i>0) M(0,i) = 0.;
207 }
208 for (Int_t i=1; i<_mref->GetNrows(); ++i) {
209 for (Int_t j=1; j<_mref->GetNrows(); ++j) {
210 M(i,j) = TMath::Power((*dm)[i],(double)j);
211 }
212 }
213 (*_M) = M.Invert();
214
215 delete dm ;
216}
217
218////////////////////////////////////////////////////////////////////////////////
219
221{
222 CacheElem* cache = (CacheElem*) _cacheMgr.getObj(0,(RooArgSet*)0) ;
223 if (cache) {
224 return cache ;
225 }
226 Int_t nVar = _varList.getSize();
227 Int_t nPdf = _pdfList.getSize();
228
229 RooAbsReal* null = 0 ;
230 vector<RooAbsReal*> meanrv(nPdf*nVar,null);
231 vector<RooAbsReal*> sigmarv(nPdf*nVar,null);
232 vector<RooAbsReal*> myrms(nVar,null);
233 vector<RooAbsReal*> mypos(nVar,null);
234 vector<RooAbsReal*> slope(nPdf*nVar,null);
235 vector<RooAbsReal*> offs(nPdf*nVar,null);
236 vector<RooAbsReal*> transVar(nPdf*nVar,null);
237 vector<RooAbsReal*> transPdf(nPdf,null);
238
239 RooArgSet ownedComps ;
240
241 RooArgList fracl ;
242
243 // fraction parameters
244 RooArgList coefList("coefList");
245 RooArgList coefList2("coefList2");
246 for (Int_t i=0; i<2*nPdf; ++i) {
247 std::string fracName = Form("frac_%d",i);
248
249 RooRealVar* frac = new RooRealVar(fracName.c_str(),fracName.c_str(),1.) ;
250
251 fracl.add(*frac); // to be set later
252 if (i<nPdf) coefList.add(*(RooRealVar*)(fracl.at(i))) ;
253 else coefList2.add(*(RooRealVar*)(fracl.at(i))) ;
254 ownedComps.add(*(RooRealVar*)(fracl.at(i))) ;
255 }
256
257 RooAddPdf* theSumPdf = 0;
258 std::string sumpdfName = Form("%s_sumpdf",GetName());
259
260 if (_useHorizMorph) {
261 // mean and sigma
262 RooArgList varList(_varList) ;
263 for (Int_t i=0; i<nPdf; ++i) {
264 for (Int_t j=0; j<nVar; ++j) {
265
266 std::string meanName = Form("%s_mean_%d_%d",GetName(),i,j);
267 std::string sigmaName = Form("%s_sigma_%d_%d",GetName(),i,j);
268
269 RooAbsMoment* mom = nVar==1 ?
270 ((RooAbsPdf*)_pdfList.at(i))->sigma((RooRealVar&)*varList.at(j)) :
271 ((RooAbsPdf*)_pdfList.at(i))->sigma((RooRealVar&)*varList.at(j),varList) ;
272
275
276 sigmarv[ij(i,j)] = mom ;
277 meanrv[ij(i,j)] = mom->mean() ;
278
279 ownedComps.add(*sigmarv[ij(i,j)]) ;
280 }
281 }
282
283 // slope and offset (to be set later, depend on m)
284 for (Int_t j=0; j<nVar; ++j) {
285 RooArgList meanList("meanList");
286 RooArgList rmsList("rmsList");
287 for (Int_t i=0; i<nPdf; ++i) {
288 meanList.add(*meanrv[ij(i,j)]);
289 rmsList.add(*sigmarv[ij(i,j)]);
290 }
291 std::string myrmsName = Form("%s_rms_%d",GetName(),j);
292 std::string myposName = Form("%s_pos_%d",GetName(),j);
293 myrms[j] = new RooAddition(myrmsName.c_str(),myrmsName.c_str(),rmsList,coefList2);
294 mypos[j] = new RooAddition(myposName.c_str(),myposName.c_str(),meanList,coefList2);
295 ownedComps.add(RooArgSet(*myrms[j],*mypos[j])) ;
296 }
297
298 // construction of unit pdfs
299 _pdfItr->Reset();
300 RooAbsPdf* pdf;
301 RooArgList transPdfList;
302
303 for (Int_t i=0; i<nPdf; ++i) {
304 _varItr->Reset() ;
305 RooRealVar* var ;
306
307 pdf = (RooAbsPdf*)_pdfItr->Next();
308 std::string pdfName = Form("pdf_%d",i);
309 RooCustomizer cust(*pdf,pdfName.c_str());
310
311 for (Int_t j=0; j<nVar; ++j) {
312 // slope and offset formulas
313 std::string slopeName = Form("%s_slope_%d_%d",GetName(),i,j);
314 std::string offsetName = Form("%s_offset_%d_%d",GetName(),i,j);
315 slope[ij(i,j)] = new RooFormulaVar(slopeName.c_str(),"@0/@1",RooArgList(*sigmarv[ij(i,j)],*myrms[j]));
316 offs[ij(i,j)] = new RooFormulaVar(offsetName.c_str(),"@0-(@1*@2)",RooArgList(*meanrv[ij(i,j)],*mypos[j],*slope[ij(i,j)]));
317 ownedComps.add(RooArgSet(*slope[ij(i,j)],*offs[ij(i,j)])) ;
318 // linear transformations, so pdf can be renormalized
319 var = (RooRealVar*)(_varItr->Next());
320 std::string transVarName = Form("%s_transVar_%d_%d",GetName(),i,j);
321 //transVar[ij(i,j)] = new RooFormulaVar(transVarName.c_str(),transVarName.c_str(),"@0*@1+@2",RooArgList(*var,*slope[ij(i,j)],*offs[ij(i,j)]));
322
323 transVar[ij(i,j)] = new RooLinearVar(transVarName.c_str(),transVarName.c_str(),*var,*slope[ij(i,j)],*offs[ij(i,j)]);
324
325 // *** WVE this is important *** this declares that frac effectively depends on the morphing parameters
326 // This will prevent the likelihood optimizers from erroneously declaring terms constant
327 transVar[ij(i,j)]->addServer((RooAbsArg&)m.arg()) ;
328
329 ownedComps.add(*transVar[ij(i,j)]) ;
330 cust.replaceArg(*var,*transVar[ij(i,j)]);
331 }
332 transPdf[i] = (RooAbsPdf*) cust.build() ;
333 transPdfList.add(*transPdf[i]);
334 ownedComps.add(*transPdf[i]) ;
335 }
336 // sum pdf
337 theSumPdf = new RooAddPdf(sumpdfName.c_str(),sumpdfName.c_str(),transPdfList,coefList);
338 }
339 else {
340 theSumPdf = new RooAddPdf(sumpdfName.c_str(),sumpdfName.c_str(),_pdfList,coefList);
341 }
342
343 // *** WVE this is important *** this declares that frac effectively depends on the morphing parameters
344 // This will prevent the likelihood optimizers from erroneously declaring terms constant
345 theSumPdf->addServer((RooAbsArg&)m.arg()) ;
346 theSumPdf->addOwnedComponents(ownedComps) ;
347
348 // change tracker for fraction parameters
349 std::string trackerName = Form("%s_frac_tracker",GetName()) ;
350 RooChangeTracker* tracker = new RooChangeTracker(trackerName.c_str(),trackerName.c_str(),m.arg(),kTRUE) ;
351
352 // Store it in the cache
353 cache = new CacheElem(*theSumPdf,*tracker,fracl) ;
354 _cacheMgr.setObj(0,0,cache,0) ;
355
356 cache->calculateFractions(*this, kFALSE);
357 return cache ;
358}
359
360////////////////////////////////////////////////////////////////////////////////
361
363{
364 return RooArgList(*_sumPdf,*_tracker) ;
365}
366
367////////////////////////////////////////////////////////////////////////////////
368
370{
371 delete _sumPdf ;
372 delete _tracker ;
373}
374
375////////////////////////////////////////////////////////////////////////////////
376/// Special version of getVal() overrides RooAbsReal::getVal() to save value of current normalization set
377
379{
380 _curNormSet = set ? (RooArgSet*)set : (RooArgSet*)&_varList ;
381 return RooAbsPdf::getVal(set) ;
382}
383
384////////////////////////////////////////////////////////////////////////////////
385
387{
388 CacheElem* cache = getCache(nset ? nset : _curNormSet) ;
389
390 if (cache->_tracker->hasChanged(kTRUE)) {
391 cache->calculateFractions(*this,kFALSE); // verbose turned off
392 }
393
394 return cache->_sumPdf ;
395}
396
397////////////////////////////////////////////////////////////////////////////////
398
400{
401 CacheElem* cache = getCache(_curNormSet) ;
402
403 if (cache->_tracker->hasChanged(kTRUE)) {
404 cache->calculateFractions(*this,kFALSE); // verbose turned off
405 }
406
407 Double_t ret = cache->_sumPdf->getVal(_pdfList.nset());
408 return ret ;
409}
410
411////////////////////////////////////////////////////////////////////////////////
412
414{
415 return (RooRealVar*)(_frac.at(i)) ;
416}
417
418////////////////////////////////////////////////////////////////////////////////
419
421{
422 return (RooRealVar*)(_frac.at(i)) ;
423}
424
425////////////////////////////////////////////////////////////////////////////////
426
428{
429 Int_t nPdf = self._pdfList.getSize();
430
431 Double_t dm = self.m - (*self._mref)[0];
432
433 // fully non-linear
434 double sumposfrac=0.;
435 for (Int_t i=0; i<nPdf; ++i) {
436 double ffrac=0.;
437 for (Int_t j=0; j<nPdf; ++j) { ffrac += (*self._M)(j,i) * (j==0?1.:TMath::Power(dm,(double)j)); }
438 if (ffrac>=0) sumposfrac+=ffrac;
439 // fractions for pdf
440 ((RooRealVar*)frac(i))->setVal(ffrac);
441 // fractions for rms and mean
442 ((RooRealVar*)frac(nPdf+i))->setVal(ffrac);
443 if (verbose) { cout << ffrac << endl; }
444 }
445
446 // various mode settings
447 int imin = self.idxmin(self.m);
448 int imax = self.idxmax(self.m);
449 double mfrac = (self.m-(*self._mref)[imin])/((*self._mref)[imax]-(*self._mref)[imin]);
450 switch (self._setting) {
451 case NonLinear:
452 // default already set above
453 break;
454
455 case SineLinear:
456 mfrac = TMath::Sin( TMath::PiOver2()*mfrac ); // this gives a continuous differentiable transition between grid points.
457
458 // now fall through to Linear case
459
460 case Linear:
461 for (Int_t i=0; i<2*nPdf; ++i)
462 ((RooRealVar*)frac(i))->setVal(0.);
463 if (imax>imin) { // m in between mmin and mmax
464 ((RooRealVar*)frac(imin))->setVal(1.-mfrac);
465 ((RooRealVar*)frac(nPdf+imin))->setVal(1.-mfrac);
466 ((RooRealVar*)frac(imax))->setVal(mfrac);
467 ((RooRealVar*)frac(nPdf+imax))->setVal(mfrac);
468 } else if (imax==imin) { // m outside mmin and mmax
469 ((RooRealVar*)frac(imin))->setVal(1.);
470 ((RooRealVar*)frac(nPdf+imin))->setVal(1.);
471 }
472 break;
474 for (Int_t i=0; i<nPdf; ++i)
475 ((RooRealVar*)frac(i))->setVal(0.);
476 if (imax>imin) { // m in between mmin and mmax
477 ((RooRealVar*)frac(imin))->setVal(1.-mfrac);
478 ((RooRealVar*)frac(imax))->setVal(mfrac);
479 } else if (imax==imin) { // m outside mmin and mmax
480 ((RooRealVar*)frac(imin))->setVal(1.);
481 }
482 break;
484 for (Int_t i=0; i<nPdf; ++i) {
485 if (((RooRealVar*)frac(i))->getVal()<0) ((RooRealVar*)frac(i))->setVal(0.);
486 ((RooRealVar*)frac(i))->setVal(((RooRealVar*)frac(i))->getVal()/sumposfrac);
487 }
488 break;
489 }
490
491}
492
493////////////////////////////////////////////////////////////////////////////////
494
495int RooMomentMorph::idxmin(const double& mval) const
496{
497 int imin(0);
498 Int_t nPdf = _pdfList.getSize();
499 double mmin=-DBL_MAX;
500 for (Int_t i=0; i<nPdf; ++i)
501 if ( (*_mref)[i]>mmin && (*_mref)[i]<=mval ) { mmin=(*_mref)[i]; imin=i; }
502 return imin;
503}
504
505
506////////////////////////////////////////////////////////////////////////////////
507
508int RooMomentMorph::idxmax(const double& mval) const
509{
510 int imax(0);
511 Int_t nPdf = _pdfList.getSize();
512 double mmax=DBL_MAX;
513 for (Int_t i=0; i<nPdf; ++i)
514 if ( (*_mref)[i]<mmax && (*_mref)[i]>=mval ) { mmax=(*_mref)[i]; imax=i; }
515 return imax;
516}
#define coutW(a)
Definition: RooMsgService.h:32
#define coutE(a)
Definition: RooMsgService.h:33
const Bool_t kFALSE
Definition: RtypesCore.h:90
const Bool_t kTRUE
Definition: RtypesCore.h:89
#define ClassImp(name)
Definition: Rtypes.h:361
char name[80]
Definition: TGX11.cxx:109
TMatrixT< Double_t > TMatrixD
Definition: TMatrixDfwd.h:22
char * Form(const char *fmt,...)
TVectorT< Double_t > TVectorD
Definition: TVectorDfwd.h:22
RooAbsArg is the common abstract base class for objects that represent a value (of arbitrary type) an...
Definition: RooAbsArg.h:73
friend class RooArgSet
Definition: RooAbsArg.h:572
void addServer(RooAbsArg &server, Bool_t valueProp=kTRUE, Bool_t shapeProp=kFALSE, std::size_t refCount=1)
Register another RooAbsArg as a server to us, ie, declare that we depend on it.
Definition: RooAbsArg.cxx:354
Bool_t addOwnedComponents(const RooArgSet &comps)
Take ownership of the contents of 'comps'.
Definition: RooAbsArg.cxx:2107
void setLocalNoDirtyInhibit(Bool_t flag) const
Definition: RooAbsArg.h:647
Int_t getSize() const
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
TIterator * createIterator(Bool_t dir=kIterForward) const
TIterator-style iteration over contained elements.
RooAbsMoment represents the first, second, or third order derivative of any RooAbsReal as calculated ...
Definition: RooAbsMoment.h:27
RooAbsReal * mean()
Definition: RooAbsMoment.h:37
const RooArgSet * nset() const
Definition: RooAbsProxy.h:46
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:60
friend class RooAddPdf
Definition: RooAbsReal.h:524
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:90
RooAbsMoment * sigma(RooRealVar &obs)
Definition: RooAbsReal.h:334
RooAddPdf is an efficient implementation of a sum of PDFs of the form.
Definition: RooAddPdf.h:29
RooAddition calculates the sum of a set of RooAbsReal terms, or when constructed with two sets,...
Definition: RooAddition.h:26
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgList.h:21
RooAbsArg * at(Int_t idx) const
Return object at given index, or nullptr if index is out of range.
Definition: RooArgList.h:74
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:28
virtual Bool_t add(const RooAbsCollection &col, Bool_t silent=kFALSE)
Add a collection of arguments to this collection by calling add() for each element in the source coll...
Definition: RooArgSet.h:88
T * getObj(const RooArgSet *nset, Int_t *sterileIndex=0, const TNamed *isetRangeName=0)
Int_t setObj(const RooArgSet *nset, T *obj, const TNamed *isetRangeName=0)
RooChangeTracker is a meta object that tracks value changes in a given set of RooAbsArgs by registeri...
Bool_t hasChanged(Bool_t clearState)
Returns true if state has changed since last call with clearState=kTRUE.
RooConstVar represent a constant real-valued object.
Definition: RooConstVar.h:25
RooCustomizer is a factory class to produce clones of a prototype composite PDF object with the same ...
Definition: RooCustomizer.h:32
void replaceArg(const RooAbsArg &orig, const RooAbsArg &subst)
Replace any occurence of arg 'orig' with arg 'subst'.
RooAbsArg * build(const char *masterCatState, Bool_t verbose=kFALSE)
Build a clone of the prototype executing all registered 'replace' rules and 'split' rules for the mas...
A RooFormulaVar is a generic implementation of a real-valued object, which takes a RooArgList of serv...
Definition: RooFormulaVar.h:29
RooLinearVar is the most general form of a derived real-valued object that can be used by RooRealInte...
Definition: RooLinearVar.h:29
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Reimplementation of standard RooArgList::add()
void calculateFractions(const RooMomentMorph &self, Bool_t verbose=kTRUE) const
RooRealVar * frac(Int_t i)
virtual RooArgList containedArgs(Action)
RooChangeTracker * _tracker
RooObjCacheManager _cacheMgr
friend class CacheElem
Current normalization set.
TIterator * _pdfItr
do not persist
RooSetProxy _varList
Double_t evaluate() const
Evaluate this PDF / function / constant. Needs to be overridden by all derived classes.
virtual ~RooMomentMorph()
RooMomentMorph()
coverity[UNINIT_CTOR]
RooRealProxy m
RooAbsPdf * sumPdf(const RooArgSet *nset)
Int_t ij(const Int_t &i, const Int_t &j) const
RooArgSet * _curNormSet
The cache manager.
RooListProxy _pdfList
TIterator * _varItr
TMatrixD * _M
int idxmin(const double &m) const
virtual Double_t getVal(const RooArgSet *set=0) const
Special version of getVal() overrides RooAbsReal::getVal() to save value of current normalization set...
int idxmax(const double &m) const
CacheElem * getCache(const RooArgSet *nset) const
TVectorD * _mref
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:35
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Overloaded RooArgSet::add() method inserts 'var' into set and registers 'var' as server to owner with...
const T & arg() const
Return reference to object held in proxy.
Iterator abstract base class.
Definition: TIterator.h:30
virtual void Reset()=0
virtual TObject * Next()=0
TMatrixT< Element > & Invert(Double_t *det=0)
Invert the matrix and calculate its determinant.
Definition: TMatrixT.cxx:1399
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
Int_t GetNrows() const
Definition: TVectorT.h:75
@ InputArguments
Definition: RooGlobalFunc.h:68
null_t< F > null()
constexpr Double_t PiOver2()
Definition: TMath.h:52
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Definition: TMath.h:725
Double_t Sin(Double_t)
Definition: TMath.h:627
auto * m
Definition: textangle.C:8