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
38 : _cacheMgr(this,10,true,true), _curNormSet(0), _mref(0), _M(0), _useHorizMorph(true)
39{
40}
41
42////////////////////////////////////////////////////////////////////////////////
43/// CTOR
44
45RooMomentMorph::RooMomentMorph(const char *name, const char *title,
46 RooAbsReal& _m,
47 const RooArgList& varList,
48 const RooArgList& pdfList,
49 const TVectorD& mrefpoints,
50 Setting setting) :
51 RooAbsPdf(name,title),
52 _cacheMgr(this,10,true,true),
53 m("m","m",this,_m),
54 _varList("varList","List of variables",this),
55 _pdfList("pdfList","List of pdfs",this),
56 _setting(setting),
57 _useHorizMorph(true)
58{
59 // observables
60 for (auto *var : varList) {
61 if (!dynamic_cast<RooAbsReal*>(var)) {
62 coutE(InputArguments) << "RooMomentMorph::ctor(" << GetName() << ") ERROR: variable " << var->GetName() << " is not of type RooAbsReal" << endl ;
63 throw string("RooPolyMorh::ctor() ERROR variable is not of type RooAbsReal") ;
64 }
65 _varList.add(*var) ;
66 }
67
68 // reference p.d.f.s
69
70 for (auto const *pdf : pdfList) {
71 if (!dynamic_cast<RooAbsPdf const*>(pdf)) {
72 coutE(InputArguments) << "RooMomentMorph::ctor(" << GetName() << ") ERROR: pdf " << pdf->GetName() << " is not of type RooAbsPdf" << endl ;
73 throw string("RooPolyMorh::ctor() ERROR pdf is not of type RooAbsPdf") ;
74 }
75 _pdfList.add(*pdf) ;
76 }
77
78 _mref = new TVectorD(mrefpoints);
79
80 // initialization
81 initialize();
82}
83
84////////////////////////////////////////////////////////////////////////////////
85/// CTOR
86
87RooMomentMorph::RooMomentMorph(const char *name, const char *title,
88 RooAbsReal& _m,
89 const RooArgList& varList,
90 const RooArgList& pdfList,
91 const RooArgList& mrefList,
92 Setting setting) :
93 RooAbsPdf(name,title),
94 _cacheMgr(this,10,true,true),
95 m("m","m",this,_m),
96 _varList("varList","List of variables",this),
97 _pdfList("pdfList","List of pdfs",this),
98 _setting(setting),
99 _useHorizMorph(true)
100{
101 // observables
102 for (auto *var : varList) {
103 if (!dynamic_cast<RooAbsReal*>(var)) {
104 coutE(InputArguments) << "RooMomentMorph::ctor(" << GetName() << ") ERROR: variable " << var->GetName() << " is not of type RooAbsReal" << endl ;
105 throw string("RooPolyMorh::ctor() ERROR variable is not of type RooAbsReal") ;
106 }
107 _varList.add(*var) ;
108 }
109
110 // reference p.d.f.s
111 for (auto const *pdf : pdfList) {
112 if (!dynamic_cast<RooAbsPdf const*>(pdf)) {
113 coutE(InputArguments) << "RooMomentMorph::ctor(" << GetName() << ") ERROR: pdf " << pdf->GetName() << " is not of type RooAbsPdf" << endl ;
114 throw string("RooPolyMorh::ctor() ERROR pdf is not of type RooAbsPdf") ;
115 }
116 _pdfList.add(*pdf) ;
117 }
118
119 // reference points in m
120 _mref = new TVectorD(mrefList.getSize());
121 Int_t i = 0;
122 for (auto *mref : mrefList) {
123 if (!dynamic_cast<RooAbsReal*>(mref)) {
124 coutE(InputArguments) << "RooMomentMorph::ctor(" << GetName() << ") ERROR: mref " << mref->GetName() << " is not of type RooAbsReal" << endl ;
125 throw string("RooPolyMorh::ctor() ERROR mref is not of type RooAbsReal") ;
126 }
127 if (!dynamic_cast<RooConstVar*>(mref)) {
128 coutW(InputArguments) << "RooMomentMorph::ctor(" << GetName() << ") WARNING mref point " << i << " is not a constant, taking a snapshot of its value" << endl ;
129 }
130 (*_mref)[i] = static_cast<RooAbsReal*>(mref)->getVal() ;
131 i++;
132 }
133
134 // initialization
135 initialize();
136}
137
138////////////////////////////////////////////////////////////////////////////////
139
141 RooAbsPdf(other,name),
142 _cacheMgr(other._cacheMgr,this),
143 _curNormSet(0),
144 m("m",this,other.m),
145 _varList("varList",this,other._varList),
146 _pdfList("pdfList",this,other._pdfList),
147 _setting(other._setting),
148 _useHorizMorph(other._useHorizMorph)
149{
150 _mref = new TVectorD(*other._mref) ;
151
152 // initialization
153 initialize();
154}
155
156////////////////////////////////////////////////////////////////////////////////
157
159{
160 if (_mref) delete _mref;
161 if (_varItr) delete _varItr;
162 if (_pdfItr) delete _pdfItr;
163 if (_M) delete _M;
164}
165
166////////////////////////////////////////////////////////////////////////////////
167
169{
170 Int_t nPdf = _pdfList.getSize();
171
172 // other quantities needed
173 if (nPdf!=_mref->GetNrows()) {
174 coutE(InputArguments) << "RooMomentMorph::initialize(" << GetName() << ") ERROR: nPdf != nRefPoints" << endl ;
175 assert(0) ;
176 }
177
178 TVectorD* dm = new TVectorD(nPdf);
179 _M = new TMatrixD(nPdf,nPdf);
180
181 // transformation matrix for non-linear extrapolation, needed in evaluate()
182 TMatrixD M(nPdf,nPdf);
183 for (Int_t i=0; i<_mref->GetNrows(); ++i) {
184 (*dm)[i] = (*_mref)[i]-(*_mref)[0];
185 M(i,0) = 1.;
186 if (i>0) M(0,i) = 0.;
187 }
188 for (Int_t i=1; i<_mref->GetNrows(); ++i) {
189 for (Int_t j=1; j<_mref->GetNrows(); ++j) {
190 M(i,j) = TMath::Power((*dm)[i],(double)j);
191 }
192 }
193 (*_M) = M.Invert();
194
195 delete dm ;
196}
197
198////////////////////////////////////////////////////////////////////////////////
199
201{
202 CacheElem* cache = (CacheElem*) _cacheMgr.getObj(0,(RooArgSet*)0) ;
203 if (cache) {
204 return cache ;
205 }
206 Int_t nVar = _varList.getSize();
207 Int_t nPdf = _pdfList.getSize();
208
209 RooAbsReal* null = 0 ;
210 vector<RooAbsReal*> meanrv(nPdf*nVar,null);
211 vector<RooAbsReal*> sigmarv(nPdf*nVar,null);
212 vector<RooAbsReal*> myrms(nVar,null);
213 vector<RooAbsReal*> mypos(nVar,null);
214 vector<RooAbsReal*> slope(nPdf*nVar,null);
215 vector<RooAbsReal*> offs(nPdf*nVar,null);
216 vector<RooAbsReal*> transVar(nPdf*nVar,null);
217 vector<RooAbsReal*> transPdf(nPdf,null);
218
219 RooArgSet ownedComps ;
220
221 RooArgList fracl ;
222
223 // fraction parameters
224 RooArgList coefList("coefList");
225 RooArgList coefList2("coefList2");
226 for (Int_t i=0; i<2*nPdf; ++i) {
227 std::string fracName = Form("frac_%d",i);
228
229 RooRealVar* frac = new RooRealVar(fracName.c_str(),fracName.c_str(),1.) ;
230
231 fracl.add(*frac); // to be set later
232 if (i<nPdf) coefList.add(*(RooRealVar*)(fracl.at(i))) ;
233 else coefList2.add(*(RooRealVar*)(fracl.at(i))) ;
234 ownedComps.add(*(RooRealVar*)(fracl.at(i))) ;
235 }
236
237 RooAddPdf* theSumPdf = 0;
238 std::string sumpdfName = Form("%s_sumpdf",GetName());
239
240 if (_useHorizMorph) {
241 // mean and sigma
242 RooArgList varList(_varList) ;
243 for (Int_t i=0; i<nPdf; ++i) {
244 for (Int_t j=0; j<nVar; ++j) {
245
246 std::string meanName = Form("%s_mean_%d_%d",GetName(),i,j);
247 std::string sigmaName = Form("%s_sigma_%d_%d",GetName(),i,j);
248
249 RooAbsMoment* mom = nVar==1 ?
250 ((RooAbsPdf*)_pdfList.at(i))->sigma((RooRealVar&)*varList.at(j)) :
251 ((RooAbsPdf*)_pdfList.at(i))->sigma((RooRealVar&)*varList.at(j),varList) ;
252
253 mom->setLocalNoDirtyInhibit(true) ;
254 mom->mean()->setLocalNoDirtyInhibit(true) ;
255
256 sigmarv[ij(i,j)] = mom ;
257 meanrv[ij(i,j)] = mom->mean() ;
258
259 ownedComps.add(*sigmarv[ij(i,j)]) ;
260 }
261 }
262
263 // slope and offset (to be set later, depend on m)
264 for (Int_t j=0; j<nVar; ++j) {
265 RooArgList meanList("meanList");
266 RooArgList rmsList("rmsList");
267 for (Int_t i=0; i<nPdf; ++i) {
268 meanList.add(*meanrv[ij(i,j)]);
269 rmsList.add(*sigmarv[ij(i,j)]);
270 }
271 std::string myrmsName = Form("%s_rms_%d",GetName(),j);
272 std::string myposName = Form("%s_pos_%d",GetName(),j);
273 myrms[j] = new RooAddition(myrmsName.c_str(),myrmsName.c_str(),rmsList,coefList2);
274 mypos[j] = new RooAddition(myposName.c_str(),myposName.c_str(),meanList,coefList2);
275 ownedComps.add(RooArgSet(*myrms[j],*mypos[j])) ;
276 }
277
278 // construction of unit pdfs
279 _pdfItr->Reset();
280 RooAbsPdf* pdf;
281 RooArgList transPdfList;
282
283 for (Int_t i=0; i<nPdf; ++i) {
284 _varItr->Reset() ;
285 RooRealVar* var ;
286
287 pdf = (RooAbsPdf*)_pdfItr->Next();
288 std::string pdfName = Form("pdf_%d",i);
289 RooCustomizer cust(*pdf,pdfName.c_str());
290
291 for (Int_t j=0; j<nVar; ++j) {
292 // slope and offset formulas
293 std::string slopeName = Form("%s_slope_%d_%d",GetName(),i,j);
294 std::string offsetName = Form("%s_offset_%d_%d",GetName(),i,j);
295 slope[ij(i,j)] = new RooFormulaVar(slopeName.c_str(),"@0/@1",RooArgList(*sigmarv[ij(i,j)],*myrms[j]));
296 offs[ij(i,j)] = new RooFormulaVar(offsetName.c_str(),"@0-(@1*@2)",RooArgList(*meanrv[ij(i,j)],*mypos[j],*slope[ij(i,j)]));
297 ownedComps.add(RooArgSet(*slope[ij(i,j)],*offs[ij(i,j)])) ;
298 // linear transformations, so pdf can be renormalized
299 var = (RooRealVar*)(_varItr->Next());
300 std::string transVarName = Form("%s_transVar_%d_%d",GetName(),i,j);
301 //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)]));
302
303 transVar[ij(i,j)] = new RooLinearVar(transVarName.c_str(),transVarName.c_str(),*var,*slope[ij(i,j)],*offs[ij(i,j)]);
304
305 // *** WVE this is important *** this declares that frac effectively depends on the morphing parameters
306 // This will prevent the likelihood optimizers from erroneously declaring terms constant
307 transVar[ij(i,j)]->addServer((RooAbsArg&)m.arg()) ;
308
309 ownedComps.add(*transVar[ij(i,j)]) ;
310 cust.replaceArg(*var,*transVar[ij(i,j)]);
311 }
312 transPdf[i] = (RooAbsPdf*) cust.build() ;
313 transPdfList.add(*transPdf[i]);
314 ownedComps.add(*transPdf[i]) ;
315 }
316 // sum pdf
317 theSumPdf = new RooAddPdf(sumpdfName.c_str(),sumpdfName.c_str(),transPdfList,coefList);
318 }
319 else {
320 theSumPdf = new RooAddPdf(sumpdfName.c_str(),sumpdfName.c_str(),_pdfList,coefList);
321 }
322
323 // *** WVE this is important *** this declares that frac effectively depends on the morphing parameters
324 // This will prevent the likelihood optimizers from erroneously declaring terms constant
325 theSumPdf->addServer((RooAbsArg&)m.arg()) ;
326 theSumPdf->addOwnedComponents(ownedComps) ;
327
328 // change tracker for fraction parameters
329 std::string trackerName = Form("%s_frac_tracker",GetName()) ;
330 RooChangeTracker* tracker = new RooChangeTracker(trackerName.c_str(),trackerName.c_str(),m.arg(),true) ;
331
332 // Store it in the cache
333 cache = new CacheElem(*theSumPdf,*tracker,fracl) ;
334 _cacheMgr.setObj(0,0,cache,0) ;
335
336 cache->calculateFractions(*this, false);
337 return cache ;
338}
339
340////////////////////////////////////////////////////////////////////////////////
341
343{
344 return RooArgList(*_sumPdf,*_tracker) ;
345}
346
347////////////////////////////////////////////////////////////////////////////////
348
350{
351 delete _sumPdf ;
352 delete _tracker ;
353}
354
355////////////////////////////////////////////////////////////////////////////////
356/// Special version of getVal() overrides RooAbsReal::getVal() to save value of current normalization set
357
358double RooMomentMorph::getVal(const RooArgSet* set) const
359{
360 _curNormSet = set ? (RooArgSet*)set : (RooArgSet*)&_varList ;
361 return RooAbsPdf::getVal(set) ;
362}
363
364////////////////////////////////////////////////////////////////////////////////
365
367{
368 CacheElem* cache = getCache(nset ? nset : _curNormSet) ;
369
370 if (cache->_tracker->hasChanged(true)) {
371 cache->calculateFractions(*this,false); // verbose turned off
372 }
373
374 return cache->_sumPdf ;
375}
376
377////////////////////////////////////////////////////////////////////////////////
378
380{
381 CacheElem* cache = getCache(_curNormSet) ;
382
383 if (cache->_tracker->hasChanged(true)) {
384 cache->calculateFractions(*this,false); // verbose turned off
385 }
386
387 double ret = cache->_sumPdf->getVal(_pdfList.nset());
388 return ret ;
389}
390
391////////////////////////////////////////////////////////////////////////////////
392
394{
395 return (RooRealVar*)(_frac.at(i)) ;
396}
397
398////////////////////////////////////////////////////////////////////////////////
399
401{
402 return (RooRealVar*)(_frac.at(i)) ;
403}
404
405////////////////////////////////////////////////////////////////////////////////
406
408{
409 Int_t nPdf = self._pdfList.getSize();
410
411 double dm = self.m - (*self._mref)[0];
412
413 // fully non-linear
414 double sumposfrac=0.;
415 for (Int_t i=0; i<nPdf; ++i) {
416 double ffrac=0.;
417 for (Int_t j=0; j<nPdf; ++j) { ffrac += (*self._M)(j,i) * (j==0?1.:TMath::Power(dm,(double)j)); }
418 if (ffrac>=0) sumposfrac+=ffrac;
419 // fractions for pdf
420 ((RooRealVar*)frac(i))->setVal(ffrac);
421 // fractions for rms and mean
422 ((RooRealVar*)frac(nPdf+i))->setVal(ffrac);
423 if (verbose) { cout << ffrac << endl; }
424 }
425
426 // various mode settings
427 int imin = self.idxmin(self.m);
428 int imax = self.idxmax(self.m);
429 double mfrac = (self.m-(*self._mref)[imin])/((*self._mref)[imax]-(*self._mref)[imin]);
430 switch (self._setting) {
431 case NonLinear:
432 // default already set above
433 break;
434
435 case SineLinear:
436 mfrac = TMath::Sin( TMath::PiOver2()*mfrac ); // this gives a continuous differentiable transition between grid points.
437
438 // now fall through to Linear case
439
440 case Linear:
441 for (Int_t i=0; i<2*nPdf; ++i)
442 ((RooRealVar*)frac(i))->setVal(0.);
443 if (imax>imin) { // m in between mmin and mmax
444 ((RooRealVar*)frac(imin))->setVal(1.-mfrac);
445 ((RooRealVar*)frac(nPdf+imin))->setVal(1.-mfrac);
446 ((RooRealVar*)frac(imax))->setVal(mfrac);
447 ((RooRealVar*)frac(nPdf+imax))->setVal(mfrac);
448 } else if (imax==imin) { // m outside mmin and mmax
449 ((RooRealVar*)frac(imin))->setVal(1.);
450 ((RooRealVar*)frac(nPdf+imin))->setVal(1.);
451 }
452 break;
454 for (Int_t i=0; i<nPdf; ++i)
455 ((RooRealVar*)frac(i))->setVal(0.);
456 if (imax>imin) { // m in between mmin and mmax
457 ((RooRealVar*)frac(imin))->setVal(1.-mfrac);
458 ((RooRealVar*)frac(imax))->setVal(mfrac);
459 } else if (imax==imin) { // m outside mmin and mmax
460 ((RooRealVar*)frac(imin))->setVal(1.);
461 }
462 break;
464 for (Int_t i=0; i<nPdf; ++i) {
465 if (((RooRealVar*)frac(i))->getVal()<0) ((RooRealVar*)frac(i))->setVal(0.);
466 ((RooRealVar*)frac(i))->setVal(((RooRealVar*)frac(i))->getVal()/sumposfrac);
467 }
468 break;
469 }
470
471}
472
473////////////////////////////////////////////////////////////////////////////////
474
475int RooMomentMorph::idxmin(const double& mval) const
476{
477 int imin(0);
478 Int_t nPdf = _pdfList.getSize();
479 double mmin=-DBL_MAX;
480 for (Int_t i=0; i<nPdf; ++i)
481 if ( (*_mref)[i]>mmin && (*_mref)[i]<=mval ) { mmin=(*_mref)[i]; imin=i; }
482 return imin;
483}
484
485
486////////////////////////////////////////////////////////////////////////////////
487
488int RooMomentMorph::idxmax(const double& mval) const
489{
490 int imax(0);
491 Int_t nPdf = _pdfList.getSize();
492 double mmax=DBL_MAX;
493 for (Int_t i=0; i<nPdf; ++i)
494 if ( (*_mref)[i]<mmax && (*_mref)[i]>=mval ) { mmax=(*_mref)[i]; imax=i; }
495 return imax;
496}
#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:2189
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
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:64
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:94
friend class RooAddPdf
Definition: RooAbsReal.h:558
RooAbsMoment * sigma(RooRealVar &obs)
Definition: RooAbsReal.h:363
RooAddPdf is an efficient implementation of a sum of PDFs of the form.
Definition: RooAddPdf.h:32
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
RooArgList containedArgs(Action) override
RooRealVar * frac(Int_t i)
void calculateFractions(const RooMomentMorph &self, bool verbose=true) const
RooChangeTracker * _tracker
RooObjCacheManager _cacheMgr
friend class CacheElem
Current normalization set.
~RooMomentMorph() override
TIterator * _pdfItr
do not persist
RooSetProxy _varList
double evaluate() const override
Evaluate this PDF / function / constant. Needs to be overridden by all derived classes.
RooMomentMorph()
coverity[UNINIT_CTOR]
RooRealProxy m
RooAbsPdf * sumPdf(const RooArgSet *nset)
Int_t ij(const Int_t &i, const Int_t &j) const
virtual double getVal(const RooArgSet *set=0) const
Special version of getVal() overrides RooAbsReal::getVal() to save value of current normalization set...
RooArgSet * _curNormSet
The cache manager.
RooListProxy _pdfList
TIterator * _varItr
TMatrixD * _M
int idxmin(const double &m) const
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: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
TMarker m
Definition: textangle.C:8