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