Logo ROOT  
Reference Guide
RooMultiVarGaussian.cxx
Go to the documentation of this file.
1/*****************************************************************************
2 * Project: RooFit *
3 * Package: RooFitModels *
4 * @(#)root/roofit:$Id$
5 * Authors: *
6 * WV, Wouter Verkerke, UC Santa Barbara, verkerke@slac.stanford.edu *
7 * DK, David Kirkby, UC Irvine, dkirkby@uci.edu *
8 * *
9 * Copyright (c) 2000-2005, Regents of the University of California *
10 * and Stanford University. All rights reserved. *
11 * *
12 * Redistribution and use in source and binary forms, *
13 * with or without modification, are permitted according to the terms *
14 * listed in LICENSE (http://roofit.sourceforge.net/license.txt) *
15 *****************************************************************************/
16
17/**
18\file RooMultiVarGaussian.cxx
19\class RooMultiVarGaussian
20\ingroup Roofitcore
21
22Multivariate Gaussian p.d.f. with correlations
23**/
24
25#include "Riostream.h"
26#include <math.h>
27
28#include "RooMultiVarGaussian.h"
29#include "RooAbsReal.h"
30#include "RooRealVar.h"
31#include "RooRandom.h"
32#include "RooMath.h"
33#include "RooGlobalFunc.h"
34#include "RooConstVar.h"
35#include "TDecompChol.h"
36#include "RooFitResult.h"
37
38using namespace std;
39
41 ;
42
43////////////////////////////////////////////////////////////////////////////////
44
45RooMultiVarGaussian::RooMultiVarGaussian(const char *name, const char *title,
46 const RooArgList& xvec, const RooArgList& mu, const TMatrixDSym& cov) :
47 RooAbsPdf(name,title),
48 _x("x","Observables",this,kTRUE,kFALSE),
49 _mu("mu","Offset vector",this,kTRUE,kFALSE),
50 _cov(cov),
51 _covI(cov),
52 _z(4)
53{
54 _x.add(xvec) ;
55
56 _mu.add(mu) ;
57
59
60 // Invert covariance matrix
61 _covI.Invert() ;
62}
63
64
65////////////////////////////////////////////////////////////////////////////////
66
67RooMultiVarGaussian::RooMultiVarGaussian(const char *name, const char *title,
68 const RooArgList& xvec, const RooFitResult& fr, Bool_t reduceToConditional) :
69 RooAbsPdf(name,title),
70 _x("x","Observables",this,kTRUE,kFALSE),
71 _mu("mu","Offset vector",this,kTRUE,kFALSE),
72 _cov(reduceToConditional ? fr.conditionalCovarianceMatrix(xvec) : fr.reducedCovarianceMatrix(xvec)),
73 _covI(_cov),
74 _z(4)
75{
77
78 // Fill mu vector with constant RooRealVars
79 list<string> munames ;
80 const RooArgList& fpf = fr.floatParsFinal() ;
81 for (Int_t i=0 ; i<fpf.getSize() ; i++) {
82 if (xvec.find(fpf.at(i)->GetName())) {
83 RooRealVar* parclone = (RooRealVar*) fpf.at(i)->Clone(Form("%s_centralvalue",fpf.at(i)->GetName())) ;
84 parclone->setConstant(kTRUE) ;
85 _mu.addOwned(*parclone) ;
86 munames.push_back(fpf.at(i)->GetName()) ;
87 }
88 }
89
90 // Fill X vector in same order as mu vector
91 for (list<string>::iterator iter=munames.begin() ; iter!=munames.end() ; ++iter) {
92 RooRealVar* xvar = (RooRealVar*) xvec.find(iter->c_str()) ;
93 _x.add(*xvar) ;
94 }
95
96 // Invert covariance matrix
97 _covI.Invert() ;
98
99}
100
101
102////////////////////////////////////////////////////////////////////////////////
103
104RooMultiVarGaussian::RooMultiVarGaussian(const char *name, const char *title,
105 const RooArgList& xvec, const TVectorD& mu, const TMatrixDSym& cov) :
106 RooAbsPdf(name,title),
107 _x("x","Observables",this,kTRUE,kFALSE),
108 _mu("mu","Offset vector",this,kTRUE,kFALSE),
109 _cov(cov),
110 _covI(cov),
111 _z(4)
112{
113 _x.add(xvec) ;
114
115 for (Int_t i=0 ; i<mu.GetNrows() ; i++) {
116 _mu.add(RooFit::RooConst(mu(i))) ;
117 }
118
119 _det = _cov.Determinant() ;
120
121 // Invert covariance matrix
122 _covI.Invert() ;
123}
124
125////////////////////////////////////////////////////////////////////////////////
126
127RooMultiVarGaussian::RooMultiVarGaussian(const char *name, const char *title,
128 const RooArgList& xvec, const TMatrixDSym& cov) :
129 RooAbsPdf(name,title),
130 _x("x","Observables",this,kTRUE,kFALSE),
131 _mu("mu","Offset vector",this,kTRUE,kFALSE),
132 _cov(cov),
133 _covI(cov),
134 _z(4)
135{
136 _x.add(xvec) ;
137
138 for (Int_t i=0 ; i<xvec.getSize() ; i++) {
140 }
141
142 _det = _cov.Determinant() ;
143
144 // Invert covariance matrix
145 _covI.Invert() ;
146}
147
148
149
150////////////////////////////////////////////////////////////////////////////////
151
153 RooAbsPdf(other,name), _aicMap(other._aicMap), _x("x",this,other._x), _mu("mu",this,other._mu),
154 _cov(other._cov), _covI(other._covI), _det(other._det), _z(other._z)
155{
156}
157
158
159
160////////////////////////////////////////////////////////////////////////////////
161
163{
165 for (Int_t i=0 ; i<_mu.getSize() ; i++) {
166 _muVec[i] = ((RooAbsReal*)_mu.at(i))->getVal() ;
167 }
168}
169
170
171////////////////////////////////////////////////////////////////////////////////
172/// Represent observables as vector
173
175{
176 TVectorD x(_x.getSize()) ;
177 for (int i=0 ; i<_x.getSize() ; i++) {
178 x[i] = ((RooAbsReal*)_x.at(i))->getVal() ;
179 }
180
181 // Calculate return value
182 syncMuVec() ;
183 TVectorD x_min_mu = x - _muVec ;
184
185 Double_t alpha = x_min_mu * (_covI * x_min_mu) ;
186 return exp(-0.5*alpha) ;
187}
188
189
190
191////////////////////////////////////////////////////////////////////////////////
192
193Int_t RooMultiVarGaussian::getAnalyticalIntegral(RooArgSet& allVarsIn, RooArgSet& analVars, const char* rangeName) const
194{
195 RooArgSet allVars(allVarsIn) ;
196
197 // If allVars contains x_i it cannot contain mu_i
198 for (Int_t i=0 ; i<_x.getSize() ; i++) {
199 if (allVars.contains(*_x.at(i))) {
200 allVars.remove(*_mu.at(i),kTRUE,kTRUE) ;
201 }
202 }
203
204
205 // Analytical integral known over all observables
206 if (allVars.getSize()==_x.getSize() && !rangeName) {
207 analVars.add(allVars) ;
208 return -1 ;
209 }
210
211 Int_t code(0) ;
212
213 Int_t nx = _x.getSize() ;
214 if (nx>127) {
215 // Warn that analytical integration is only provided for the first 127 observables
216 coutW(Integration) << "RooMultiVarGaussian::getAnalyticalIntegral(" << GetName() << ") WARNING: p.d.f. has " << _x.getSize()
217 << " observables, analytical integration is only implemented for the first 127 observables" << endl ;
218 nx=127 ;
219 }
220
221 // Advertise partial analytical integral over all observables for which is wide enough to
222 // use asymptotic integral calculation
223 BitBlock bits ;
224 Bool_t anyBits(kFALSE) ;
225 syncMuVec() ;
226 for (int i=0 ; i<_x.getSize() ; i++) {
227
228 // Check if integration over observable #i is requested
229 if (allVars.find(_x.at(i)->GetName())) {
230 // Check if range is wider than Z sigma
231 RooRealVar* xi = (RooRealVar*)_x.at(i) ;
232 if (xi->getMin(rangeName)<_muVec(i)-_z*sqrt(_cov(i,i)) && xi->getMax(rangeName) > _muVec(i)+_z*sqrt(_cov(i,i))) {
233 cxcoutD(Integration) << "RooMultiVarGaussian::getAnalyticalIntegral(" << GetName()
234 << ") Advertising analytical integral over " << xi->GetName() << " as range is >" << _z << " sigma" << endl ;
235 bits.setBit(i) ;
236 anyBits = kTRUE ;
237 analVars.add(*allVars.find(_x.at(i)->GetName())) ;
238 } else {
239 cxcoutD(Integration) << "RooMultiVarGaussian::getAnalyticalIntegral(" << GetName() << ") Range of " << xi->GetName() << " is <"
240 << _z << " sigma, relying on numeric integral" << endl ;
241 }
242 }
243
244 // Check if integration over parameter #i is requested
245 if (allVars.find(_mu.at(i)->GetName())) {
246 // Check if range is wider than Z sigma
247 RooRealVar* pi = (RooRealVar*)_mu.at(i) ;
248 if (pi->getMin(rangeName)<_muVec(i)-_z*sqrt(_cov(i,i)) && pi->getMax(rangeName) > _muVec(i)+_z*sqrt(_cov(i,i))) {
249 cxcoutD(Integration) << "RooMultiVarGaussian::getAnalyticalIntegral(" << GetName()
250 << ") Advertising analytical integral over " << pi->GetName() << " as range is >" << _z << " sigma" << endl ;
251 bits.setBit(i) ;
252 anyBits = kTRUE ;
253 analVars.add(*allVars.find(_mu.at(i)->GetName())) ;
254 } else {
255 cxcoutD(Integration) << "RooMultiVarGaussian::getAnalyticalIntegral(" << GetName() << ") Range of " << pi->GetName() << " is <"
256 << _z << " sigma, relying on numeric integral" << endl ;
257 }
258 }
259
260
261 }
262
263 // Full numeric integration over requested observables maps always to code zero
264 if (!anyBits) {
265 return 0 ;
266 }
267
268 // Map BitBlock into return code
269 for (UInt_t i=0 ; i<_aicMap.size() ; i++) {
270 if (_aicMap[i]==bits) {
271 code = i+1 ;
272 }
273 }
274 if (code==0) {
275 _aicMap.push_back(bits) ;
276 code = _aicMap.size() ;
277 }
278
279 return code ;
280}
281
282
283
284////////////////////////////////////////////////////////////////////////////////
285/// Handle full integral here
286
287Double_t RooMultiVarGaussian::analyticalIntegral(Int_t code, const char* /*rangeName*/) const
288{
289 if (code==-1) {
290 return pow(2*3.14159268,_x.getSize()/2.)*sqrt(fabs(_det)) ;
291 }
292
293 // Handle partial integrals here
294
295 // Retrieve |S22|, S22bar from cache
296 AnaIntData& aid = anaIntData(code) ;
297
298 // Fill position vector for non-integrated observables
299 syncMuVec() ;
300 TVectorD u(aid.pmap.size()) ;
301 for (UInt_t i=0 ; i<aid.pmap.size() ; i++) {
302 u(i) = ((RooAbsReal*)_x.at(aid.pmap[i]))->getVal() - _muVec(aid.pmap[i]) ;
303 }
304
305 // Calculate partial integral
306 Double_t ret = pow(2*3.14159268,aid.nint/2.)/sqrt(fabs(aid.S22det))*exp(-0.5*u*(aid.S22bar*u)) ;
307
308 return ret ;
309}
310
311
312
313////////////////////////////////////////////////////////////////////////////////
314/// Check if cache entry was previously created
315
317{
318 map<int,AnaIntData>::iterator iter = _anaIntCache.find(code) ;
319 if (iter != _anaIntCache.end()) {
320 return iter->second ;
321 }
322
323 // Calculate cache contents
324
325 // Decode integration code
326 vector<int> map1,map2 ;
327 decodeCode(code,map1,map2) ;
328
329 // Rearrage observables so that all non-integrated observables
330 // go first (preserving relative order) and all integrated observables
331 // go last (preserving relative order)
332 TMatrixDSym S11, S22 ;
333 TMatrixD S12, S21 ;
334 blockDecompose(_covI,map1,map2,S11,S12,S21,S22) ;
335
336 // Begin calculation of partial integrals
337 // ___
338 // sqrt(2pi)^(#intObs) (-0.5 * u1T S22 u1 )
339 // I = ------------------- * e
340 // sqrt(|det(S22)|)
341 // ___
342 // Where S22 is the sub-matrix of covI for the integrated observables and S22
343 // is the Schur complement of S22
344 // ___ -1
345 // S22 = S11 - S12 * S22 * S21
346 //
347 // and u1 is the vector of non-integrated observables
348
349 // Calculate Schur complement S22bar
350 TMatrixD S22inv(S22) ;
351 S22inv.Invert() ;
352 TMatrixD S22bar = S11 - S12*S22inv*S21 ;
353
354 // Create new cache entry
355 AnaIntData& cacheData = _anaIntCache[code] ;
356 cacheData.S22bar.ResizeTo(S22bar) ;
357 cacheData.S22bar=S22bar ;
358 cacheData.S22det= S22.Determinant() ;
359 cacheData.pmap = map1 ;
360 cacheData.nint = map2.size() ;
361
362 return cacheData ;
363}
364
365
366
367////////////////////////////////////////////////////////////////////////////////
368/// Special case: generate all observables
369
370Int_t RooMultiVarGaussian::getGenerator(const RooArgSet& directVars, RooArgSet &generateVars, Bool_t /*staticInitOK*/) const
371{
372 if (directVars.getSize()==_x.getSize()) {
373 generateVars.add(directVars) ;
374 return -1 ;
375 }
376
377 Int_t nx = _x.getSize() ;
378 if (nx>127) {
379 // Warn that analytical integration is only provided for the first 127 observables
380 coutW(Integration) << "RooMultiVarGaussian::getGenerator(" << GetName() << ") WARNING: p.d.f. has " << _x.getSize()
381 << " observables, partial internal generation is only implemented for the first 127 observables" << endl ;
382 nx=127 ;
383 }
384
385 // Advertise partial generation over all permutations of observables
386 Int_t code(0) ;
387 BitBlock bits ;
388 for (int i=0 ; i<_x.getSize() ; i++) {
389 RooAbsArg* arg = directVars.find(_x.at(i)->GetName()) ;
390 if (arg) {
391 bits.setBit(i) ;
392// code |= (1<<i) ;
393 generateVars.add(*arg) ;
394 }
395 }
396
397 // Map BitBlock into return code
398 for (UInt_t i=0 ; i<_aicMap.size() ; i++) {
399 if (_aicMap[i]==bits) {
400 code = i+1 ;
401 }
402 }
403 if (code==0) {
404 _aicMap.push_back(bits) ;
405 code = _aicMap.size() ;
406 }
407
408
409 return code ;
410}
411
412
413
414////////////////////////////////////////////////////////////////////////////////
415/// Clear the GenData cache as its content is not invariant under changes in
416/// the mu vector.
417
419{
420 _genCache.clear() ;
421
422}
423
424
425
426
427////////////////////////////////////////////////////////////////////////////////
428/// Retrieve generator config from cache
429
431{
432 GenData& gd = genData(code) ;
433 TMatrixD& TU = gd.UT ;
434 Int_t nobs = TU.GetNcols() ;
435 vector<int>& omap = gd.omap ;
436
437 while(1) {
438
439 // Create unit Gaussian vector
440 TVectorD xgen(nobs);
441 for(Int_t k= 0; k <nobs; k++) {
442 xgen(k)= RooRandom::gaussian();
443 }
444
445 // Apply transformation matrix
446 xgen *= TU ;
447
448 // Apply shift
449 if (code == -1) {
450
451 // Simple shift if we generate all observables
452 xgen += gd.mu1 ;
453
454 } else {
455
456 // Non-generated observable dependent shift for partial generations
457
458 // mubar = mu1 + S12 S22Inv ( x2 - mu2)
459 TVectorD mubar(gd.mu1) ;
460 TVectorD x2(gd.pmap.size()) ;
461 for (UInt_t i=0 ; i<gd.pmap.size() ; i++) {
462 x2(i) = ((RooAbsReal*)_x.at(gd.pmap[i]))->getVal() ;
463 }
464 mubar += gd.S12S22I * (x2 - gd.mu2) ;
465
466 xgen += mubar ;
467
468 }
469
470 // Transfer values and check if values are in range
471 Bool_t ok(kTRUE) ;
472 for (int i=0 ; i<nobs ; i++) {
473 RooRealVar* xi = (RooRealVar*)_x.at(omap[i]) ;
474 if (xgen(i)<xi->getMin() || xgen(i)>xi->getMax()) {
475 ok = kFALSE ;
476 break ;
477 } else {
478 xi->setVal(xgen(i)) ;
479 }
480 }
481
482 // If all values are in range, accept event and return
483 // otherwise retry
484 if (ok) {
485 break ;
486 }
487 }
488
489 return;
490}
491
492
493
494////////////////////////////////////////////////////////////////////////////////
495/// WVE -- CHECK THAT GENDATA IS VALID GIVEN CURRENT VALUES OF _MU
496
498{
499 // Check if cache entry was previously created
500 map<int,GenData>::iterator iter = _genCache.find(code) ;
501 if (iter != _genCache.end()) {
502 return iter->second ;
503 }
504
505 // Create new entry
506 GenData& cacheData = _genCache[code] ;
507
508 if (code==-1) {
509
510 // Do eigen value decomposition
511 TDecompChol tdc(_cov) ;
512 tdc.Decompose() ;
513 TMatrixD U = tdc.GetU() ;
515
516 // Fill cache data
517 cacheData.UT.ResizeTo(TU) ;
518 cacheData.UT = TU ;
519 cacheData.omap.resize(_x.getSize()) ;
520 for (int i=0 ; i<_x.getSize() ; i++) {
521 cacheData.omap[i] = i ;
522 }
523 syncMuVec() ;
524 cacheData.mu1.ResizeTo(_muVec) ;
525 cacheData.mu1 = _muVec ;
526
527 } else {
528
529 // Construct observables: map1 = generated, map2 = given
530 vector<int> map1, map2 ;
531 decodeCode(code,map2,map1) ;
532
533 // Do block decomposition of covariance matrix
534 TMatrixDSym S11, S22 ;
535 TMatrixD S12, S21 ;
536 blockDecompose(_cov,map1,map2,S11,S12,S21,S22) ;
537
538 // Constructed conditional matrix form
539 // -1
540 // F(X1|X2) --> CovI --> S22bar = S11 - S12 S22 S21
541 // -1
542 // --> mu --> mubar = mu1 + S12 S22 ( x2 - mu2)
543
544 // Do eigenvalue decomposition
545 TMatrixD S22Inv(TMatrixD::kInverted,S22) ;
546 TMatrixD S22bar = S11 - S12 * (S22Inv * S21) ;
547
548 // Do eigen value decomposition of S22bar
549 TDecompChol tdc(S22bar) ;
550 tdc.Decompose() ;
551 TMatrixD U = tdc.GetU() ;
553
554 // Split mu vector into mu1 and mu2
555 TVectorD mu1(map1.size()),mu2(map2.size()) ;
556 syncMuVec() ;
557 for (UInt_t i=0 ; i<map1.size() ; i++) {
558 mu1(i) = _muVec(map1[i]) ;
559 }
560 for (UInt_t i=0 ; i<map2.size() ; i++) {
561 mu2(i) = _muVec(map2[i]) ;
562 }
563
564 // Calculate rotation matrix for mu vector
565 TMatrixD S12S22Inv = S12 * S22Inv ;
566
567 // Fill cache data
568 cacheData.UT.ResizeTo(TU) ;
569 cacheData.UT = TU ;
570 cacheData.omap = map1 ;
571 cacheData.pmap = map2 ;
572 cacheData.mu1.ResizeTo(mu1) ;
573 cacheData.mu2.ResizeTo(mu2) ;
574 cacheData.mu1 = mu1 ;
575 cacheData.mu2 = mu2 ;
576 cacheData.S12S22I.ResizeTo(S12S22Inv) ;
577 cacheData.S12S22I = S12S22Inv ;
578
579 }
580
581
582 return cacheData ;
583}
584
585
586
587
588////////////////////////////////////////////////////////////////////////////////
589/// Decode analytical integration/generation code into index map of integrated/generated (map2)
590/// and non-integrated/generated observables (map1)
591
592void RooMultiVarGaussian::decodeCode(Int_t code, vector<int>& map1, vector<int>& map2) const
593{
594 if (code<0 || code> (Int_t)_aicMap.size()) {
595 cout << "RooMultiVarGaussian::decodeCode(" << GetName() << ") ERROR don't have bit pattern for code " << code << endl ;
596 throw string("RooMultiVarGaussian::decodeCode() ERROR don't have bit pattern for code") ;
597 }
598
599 BitBlock b = _aicMap[code-1] ;
600 map1.clear() ;
601 map2.clear() ;
602 for (int i=0 ; i<_x.getSize() ; i++) {
603 if (b.getBit(i)) {
604 map2.push_back(i) ;
605 } else {
606 map1.push_back(i) ;
607 }
608 }
609}
610
611
612////////////////////////////////////////////////////////////////////////////////
613/// Block decomposition of covI according to given maps of observables
614
615void RooMultiVarGaussian::blockDecompose(const TMatrixD& input, const vector<int>& map1, const vector<int>& map2, TMatrixDSym& S11, TMatrixD& S12, TMatrixD& S21, TMatrixDSym& S22)
616{
617 // Allocate and fill reordered covI matrix in 2x2 block structure
618
619 S11.ResizeTo(map1.size(),map1.size()) ;
620 S12.ResizeTo(map1.size(),map2.size()) ;
621 S21.ResizeTo(map2.size(),map1.size()) ;
622 S22.ResizeTo(map2.size(),map2.size()) ;
623
624 for (UInt_t i=0 ; i<map1.size() ; i++) {
625 for (UInt_t j=0 ; j<map1.size() ; j++)
626 S11(i,j) = input(map1[i],map1[j]) ;
627 for (UInt_t j=0 ; j<map2.size() ; j++)
628 S12(i,j) = input(map1[i],map2[j]) ;
629 }
630 for (UInt_t i=0 ; i<map2.size() ; i++) {
631 for (UInt_t j=0 ; j<map1.size() ; j++)
632 S21(i,j) = input(map2[i],map1[j]) ;
633 for (UInt_t j=0 ; j<map2.size() ; j++)
634 S22(i,j) = input(map2[i],map2[j]) ;
635 }
636
637}
638
639
641{
642 if (ibit<32) { b0 |= (1<<ibit) ; return ; }
643 if (ibit<64) { b1 |= (1<<(ibit-32)) ; return ; }
644 if (ibit<96) { b2 |= (1<<(ibit-64)) ; return ; }
645 if (ibit<128) { b3 |= (1<<(ibit-96)) ; return ; }
646}
647
649{
650 if (ibit<32) return (b0 & (1<<ibit)) ;
651 if (ibit<64) return (b1 & (1<<(ibit-32))) ;
652 if (ibit<96) return (b2 & (1<<(ibit-64))) ;
653 if (ibit<128) return (b3 & (1<<(ibit-96))) ;
654 return kFALSE ;
655}
656
658{
659 if (b0 != other.b0) return kFALSE ;
660 if (b1 != other.b1) return kFALSE ;
661 if (b2 != other.b2) return kFALSE ;
662 if (b3 != other.b3) return kFALSE ;
663 return kTRUE ;
664}
665
666
667
668
669
#define cxcoutD(a)
Definition: RooMsgService.h:81
#define coutW(a)
Definition: RooMsgService.h:32
const Bool_t kFALSE
Definition: RtypesCore.h:101
const Bool_t kTRUE
Definition: RtypesCore.h:100
#define ClassImp(name)
Definition: Rtypes.h:375
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void input
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t b
Option_t Option_t TPoint TPoint const char x2
char name[80]
Definition: TGX11.cxx:110
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
Definition: TString.cxx:2447
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition: RooAbsArg.h:78
TObject * Clone(const char *newname=0) const override
Make a clone of an object using the Streamer facility.
Definition: RooAbsArg.h:90
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
Int_t getSize() const
Return the number of elements in the collection.
Bool_t contains(const RooAbsArg &var) const
Check if collection contains an argument with the same name as var.
virtual Bool_t remove(const RooAbsArg &var, Bool_t silent=kFALSE, Bool_t matchByNameOnly=kFALSE)
Remove the specified argument from our list.
RooAbsArg * find(const char *name) const
Find object with given name in list.
virtual Double_t getMax(const char *name=0) const
Get maximum of currently defined range.
void setConstant(Bool_t value=kTRUE)
virtual Double_t getMin(const char *name=0) const
Get minimum of currently defined range.
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:63
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:93
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
bool addOwned(RooAbsArg &var, bool silent=false) override
Overloaded RooCollection_t::addOwned() method insert object into owning set and registers object as s...
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...
RooFitResult is a container class to hold the input and output of a PDF fit to a dataset.
Definition: RooFitResult.h:40
const RooArgList & floatParsFinal() const
Return list of floating parameters after fit.
Definition: RooFitResult.h:111
Bool_t operator==(const BitBlock &other)
Multivariate Gaussian p.d.f.
void decodeCode(Int_t code, std::vector< int > &map1, std::vector< int > &map2) const
Decode analytical integration/generation code into index map of integrated/generated (map2) and non-i...
std::vector< BitBlock > _aicMap
!
static void blockDecompose(const TMatrixD &input, const std::vector< int > &map1, const std::vector< int > &map2, TMatrixDSym &S11, TMatrixD &S12, TMatrixD &S21, TMatrixDSym &S22)
Block decomposition of covI according to given maps of observables.
AnaIntData & anaIntData(Int_t code) const
Check if cache entry was previously created.
std::map< int, GenData > _genCache
!
Double_t evaluate() const override
Do not persist.
Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const override
Interface function getAnalyticalIntergral advertises the analytical integrals that are supported.
void initGenerator(Int_t code) override
Clear the GenData cache as its content is not invariant under changes in the mu vector.
Double_t analyticalIntegral(Int_t code, const char *rangeName=0) const override
Handle full integral here.
GenData & genData(Int_t code) const
WVE – CHECK THAT GENDATA IS VALID GIVEN CURRENT VALUES OF _MU.
Int_t getGenerator(const RooArgSet &directVars, RooArgSet &generateVars, Bool_t staticInitOK=kTRUE) const override
Special case: generate all observables.
std::map< int, AnaIntData > _anaIntCache
!
void generateEvent(Int_t code) override
Retrieve generator config from cache.
static Double_t gaussian(TRandom *generator=randomGenerator())
Return a Gaussian random variable with mean 0 and variance 1.
Definition: RooRandom.cxx:108
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:40
void setVal(Double_t value) override
Set value of variable to 'value'.
Definition: RooRealVar.cxx:281
Cholesky Decomposition class.
Definition: TDecompChol.h:25
Bool_t Decompose() override
Matrix A is decomposed in component U so that A = U^T * U If the decomposition succeeds,...
const TMatrixD & GetU() const
Definition: TDecompChol.h:45
Int_t GetNcols() const
Definition: TMatrixTBase.h:126
TMatrixTBase< Element > & ResizeTo(Int_t nrows, Int_t ncols, Int_t=-1) override
Set size of the matrix to nrows x ncols New dynamic elements are created, the overlapping part of the...
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
TMatrixT< Element > & Invert(Double_t *det=0)
Invert the matrix and calculate its determinant.
Definition: TMatrixT.cxx:1397
TMatrixTBase< Element > & ResizeTo(Int_t nrows, Int_t ncols, Int_t=-1) override
Set size of the matrix to nrows x ncols New dynamic elements are created, the overlapping part of the...
Definition: TMatrixT.cxx:1211
const char * GetName() const override
Returns name of object.
Definition: TNamed.h:47
TVectorT< Element > & ResizeTo(Int_t lwb, Int_t upb)
Resize the vector to [lwb:upb] .
Definition: TVectorT.cxx:294
Int_t GetNrows() const
Definition: TVectorT.h:75
RooConstVar & RooConst(Double_t val)
RVec< PromoteTypes< T0, T1 > > pow(const T0 &x, const RVec< T1 > &v)
Definition: RVec.hxx:1753
RVec< PromoteType< T > > exp(const RVec< T > &v)
Definition: RVec.hxx:1744
Double_t x[n]
Definition: legend1.C:17
VecExpr< UnaryOp< Sqrt< T >, VecExpr< A, T, D >, T >, T, D > sqrt(const VecExpr< A, T, D > &rhs)
VecExpr< UnaryOp< Fabs< T >, VecExpr< A, T, D >, T >, T, D > fabs(const VecExpr< A, T, D > &rhs)
@ Integration
Definition: RooGlobalFunc.h:63
static constexpr double pi