Logo ROOT  
Reference Guide
RooHistPdf.cxx
Go to the documentation of this file.
1/*****************************************************************************
2 * Project: RooFit *
3 * Package: RooFitCore *
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 RooHistPdf.cxx
19\class RooHistPdf
20\ingroup Roofitcore
21
22RooHistPdf implements a probablity density function sampled from a
23multidimensional histogram. The histogram distribution is explicitly
24normalized by RooHistPdf and can have an arbitrary number of real or
25discrete dimensions.
26**/
27
28#include "Riostream.h"
29
30#include "RooHistPdf.h"
31#include "RooDataHist.h"
32#include "RooMsgService.h"
33#include "RooRealVar.h"
34#include "RooCategory.h"
35#include "RooWorkspace.h"
36#include "RooGlobalFunc.h"
37#include "RooHelpers.h"
38
39#include "TError.h"
40#include "TBuffer.h"
41
42using namespace std;
43
45
46
47
48
49////////////////////////////////////////////////////////////////////////////////
50/// Default constructor
51/// coverity[UNINIT_CTOR]
52
53RooHistPdf::RooHistPdf() : _dataHist(0), _totVolume(0), _unitNorm(false)
54{
55
56}
57
58
59////////////////////////////////////////////////////////////////////////////////
60/// Constructor from a RooDataHist. RooDataHist dimensions
61/// can be either real or discrete. See RooDataHist::RooDataHist for details on the binning.
62/// RooHistPdf neither owns or clone 'dhist' and the user must ensure the input histogram exists
63/// for the entire life span of this PDF.
64
65RooHistPdf::RooHistPdf(const char *name, const char *title, const RooArgSet& vars,
66 const RooDataHist& dhist, Int_t intOrder) :
67 RooAbsPdf(name,title),
68 _pdfObsList("pdfObs","List of p.d.f. observables",this),
69 _dataHist((RooDataHist*)&dhist),
70 _codeReg(10),
71 _intOrder(intOrder),
72 _cdfBoundaries(false),
73 _totVolume(0),
74 _unitNorm(false)
75{
77 _pdfObsList.add(vars) ;
78
79 // Verify that vars and dhist.get() have identical contents
80 const RooArgSet* dvars = dhist.get() ;
81 if (vars.getSize()!=dvars->getSize()) {
82 coutE(InputArguments) << "RooHistPdf::ctor(" << GetName()
83 << ") ERROR variable list and RooDataHist must contain the same variables." << endl ;
84 assert(0) ;
85 }
86 for (const auto arg : vars) {
87 if (!dvars->find(arg->GetName())) {
88 coutE(InputArguments) << "RooHistPdf::ctor(" << GetName()
89 << ") ERROR variable list and RooDataHist must contain the same variables." << endl ;
90 assert(0) ;
91 }
92 }
93
94
95 // Adjust ranges of _histObsList to those of _dataHist
96 for (const auto hobs : _histObsList) {
97 // Guaranteed to succeed, since checked above in ctor
98 RooAbsArg* dhobs = dhist.get()->find(hobs->GetName()) ;
99 RooRealVar* dhreal = dynamic_cast<RooRealVar*>(dhobs) ;
100 if (dhreal){
101 ((RooRealVar*)hobs)->setRange(dhreal->getMin(),dhreal->getMax()) ;
102 }
103 }
104
105}
106
107
108
109
110////////////////////////////////////////////////////////////////////////////////
111/// Constructor from a RooDataHist. The first list of observables are the p.d.f.
112/// observables, which may any RooAbsReal (function or variable). The second list
113/// are the corresponding observables in the RooDataHist which must be of type
114/// RooRealVar or RooCategory This constructor thus allows to apply a coordinate transformation
115/// on the histogram data to be applied.
116
117RooHistPdf::RooHistPdf(const char *name, const char *title, const RooArgList& pdfObs,
118 const RooArgList& histObs, const RooDataHist& dhist, Int_t intOrder) :
119 RooAbsPdf(name,title),
120 _pdfObsList("pdfObs","List of p.d.f. observables",this),
121 _dataHist((RooDataHist*)&dhist),
122 _codeReg(10),
123 _intOrder(intOrder),
124 _cdfBoundaries(false),
125 _totVolume(0),
126 _unitNorm(false)
127{
128 _histObsList.addClone(histObs) ;
129 _pdfObsList.add(pdfObs) ;
130
131 // Verify that vars and dhist.get() have identical contents
132 const RooArgSet* dvars = dhist.get() ;
133 if (histObs.getSize()!=dvars->getSize()) {
134 coutE(InputArguments) << "RooHistPdf::ctor(" << GetName()
135 << ") ERROR histogram variable list and RooDataHist must contain the same variables." << endl ;
136 throw(string("RooHistPdf::ctor() ERROR: histogram variable list and RooDataHist must contain the same variables")) ;
137 }
138
139 for (const auto arg : histObs) {
140 if (!dvars->find(arg->GetName())) {
141 coutE(InputArguments) << "RooHistPdf::ctor(" << GetName()
142 << ") ERROR variable list and RooDataHist must contain the same variables." << endl ;
143 throw(string("RooHistPdf::ctor() ERROR: histogram variable list and RooDataHist must contain the same variables")) ;
144 }
145 if (!arg->isFundamental()) {
146 coutE(InputArguments) << "RooHistPdf::ctor(" << GetName()
147 << ") ERROR all elements of histogram observables set must be of type RooRealVar or RooCategory." << endl ;
148 throw(string("RooHistPdf::ctor() ERROR all elements of histogram observables set must be of type RooRealVar or RooCategory.")) ;
149 }
150 }
151
152
153 // Adjust ranges of _histObsList to those of _dataHist
154 for (const auto hobs : _histObsList) {
155 // Guaranteed to succeed, since checked above in ctor
156 RooAbsArg* dhobs = dhist.get()->find(hobs->GetName()) ;
157 RooRealVar* dhreal = dynamic_cast<RooRealVar*>(dhobs) ;
158 if (dhreal){
159 ((RooRealVar*)hobs)->setRange(dhreal->getMin(),dhreal->getMax()) ;
160 }
161 }
162}
163
164
165
166////////////////////////////////////////////////////////////////////////////////
167/// Copy constructor
168
169RooHistPdf::RooHistPdf(const RooHistPdf& other, const char* name) :
170 RooAbsPdf(other,name),
171 _pdfObsList("pdfObs",this,other._pdfObsList),
172 _dataHist(other._dataHist),
173 _codeReg(other._codeReg),
174 _intOrder(other._intOrder),
175 _cdfBoundaries(other._cdfBoundaries),
176 _totVolume(other._totVolume),
177 _unitNorm(other._unitNorm)
178{
180
181}
182
183
184
185
186////////////////////////////////////////////////////////////////////////////////
187/// Destructor
188
190{
191
192}
193
194
195
196
197
198////////////////////////////////////////////////////////////////////////////////
199/// Return the current value: The value of the bin enclosing the current coordinates
200/// of the observables, normalized by the histograms contents. Interpolation
201/// is applied if the RooHistPdf is configured to do that.
202
204{
205 // Transfer values from
206 for (unsigned int i=0; i < _pdfObsList.size(); ++i) {
207 RooAbsArg* harg = _histObsList[i];
208 RooAbsArg* parg = _pdfObsList[i];
209
210 if (harg != parg) {
211 parg->syncCache() ;
212 harg->copyCache(parg,true) ;
213 if (!harg->inRange(0)) {
214 return 0 ;
215 }
216 }
217 }
218
220
221 return std::max(ret, 0.0);
222}
223
224
225////////////////////////////////////////////////////////////////////////////////
226/// Return the total volume spanned by the observables of the RooHistPdf
227
229{
230 // Return previously calculated value, if any
231 if (_totVolume>0) {
232 return _totVolume ;
233 }
234 _totVolume = 1. ;
235
236 for (const auto arg : _histObsList) {
237 RooRealVar* real = dynamic_cast<RooRealVar*>(arg) ;
238 if (real) {
239 _totVolume *= (real->getMax()-real->getMin()) ;
240 } else {
241 RooCategory* cat = dynamic_cast<RooCategory*>(arg) ;
242 if (cat) {
243 _totVolume *= cat->numTypes() ;
244 }
245 }
246 }
247
248 return _totVolume ;
249}
250
251namespace {
252bool fullRange(const RooAbsArg& x, const RooAbsArg& y ,const char* range)
253{
254 const RooAbsRealLValue *_x = dynamic_cast<const RooAbsRealLValue*>(&x);
255 const RooAbsRealLValue *_y = dynamic_cast<const RooAbsRealLValue*>(&y);
256 if (!_x || !_y) return false;
257 if (!range || !strlen(range) || !_x->hasRange(range) ||
258 _x->getBinningPtr(range)->isParameterized()) {
259 // parameterized ranges may be full range now, but that might change,
260 // so return false
261 if (range && strlen(range) && _x->getBinningPtr(range)->isParameterized())
262 return false;
263 return (_x->getMin() == _y->getMin() && _x->getMax() == _y->getMax());
264 }
265 return (_x->getMin(range) == _y->getMin() && _x->getMax(range) == _y->getMax());
266}
267}
268
269
271 RooArgSet& analVars,
272 const char* rangeName,
273 RooArgSet const& histObsList,
274 RooSetProxy const& pdfObsList,
275 Int_t intOrder) {
276 // First make list of pdf observables to histogram observables
277 // and select only those for which the integral is over the full range
278
279 Int_t code = 0;
280 Int_t frcode = 0;
281 for (unsigned int n=0; n < pdfObsList.size() && n < histObsList.size(); ++n) {
282 const auto pa = pdfObsList[n];
283 const auto ha = histObsList[n];
284
285 if (allVars.find(*pa)) {
286 code |= 2 << n;
287 analVars.add(*pa);
288 if (fullRange(*pa, *ha, rangeName)) {
289 frcode |= 2 << n;
290 }
291 }
292 }
293
294 if (code == frcode) {
295 // integrate over full range of all observables - use bit 0 to indicate
296 // full range integration over all observables
297 code |= 1;
298 }
299
300 // Disable partial analytical integrals if interpolation is used, and we
301 // integrate over sub-ranges, but leave them enabled when we integrate over
302 // the full range of one or several variables
303 if (intOrder > 1 && !(code & 1)) {
304 analVars.removeAll();
305 return 0;
306 }
307 return (code >= 2) ? code : 0;
308}
309
310
312 const char* rangeName,
313 RooArgSet const& histObsList,
314 RooSetProxy const& pdfObsList,
315 RooDataHist& dataHist,
316 bool histFuncMode) {
317 // Simplest scenario, full-range integration over all dependents
318 if (((2 << histObsList.getSize()) - 1) == code) {
319 return dataHist.sum(histFuncMode);
320 }
321
322 // Partial integration scenario, retrieve set of variables, calculate partial
323 // sum, figure out integration ranges (if needed)
324 RooArgSet intSet;
325 std::map<const RooAbsArg*, std::pair<double, double> > ranges;
326 for (unsigned int n=0; n < pdfObsList.size() && n < histObsList.size(); ++n) {
327 const auto pa = pdfObsList[n];
328 const auto ha = histObsList[n];
329
330 if (code & (2 << n)) {
331 intSet.add(*ha);
332 }
333 if (!(code & 1)) {
334 ranges[ha] = RooHelpers::getRangeOrBinningInterval(pa, rangeName);
335 }
336 // WVE must sync hist slice list values to pdf slice list
337 // Transfer values from
338 if (ha != pa) {
339 pa->syncCache();
340 ha->copyCache(pa,true);
341 }
342 }
343
344 double ret = (code & 1) ? dataHist.sum(intSet,histObsList,true,!histFuncMode) :
345 dataHist.sum(intSet,histObsList,true,!histFuncMode, ranges);
346
347 return ret ;
348}
349
350////////////////////////////////////////////////////////////////////////////////
351/// Determine integration scenario. If no interpolation is used,
352/// RooHistPdf can perform all integrals over its dependents
353/// analytically via partial or complete summation of the input
354/// histogram. If interpolation is used on the integral over
355/// all histogram observables is supported
356
357Int_t RooHistPdf::getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* rangeName) const
358{
359 return getAnalyticalIntegral(allVars, analVars, rangeName, _histObsList, _pdfObsList, _intOrder);
360}
361
362
363////////////////////////////////////////////////////////////////////////////////
364/// Return integral identified by 'code'. The actual integration
365/// is deferred to RooDataHist::sum() which implements partial
366/// or complete summation over the histograms contents.
367
368double RooHistPdf::analyticalIntegral(Int_t code, const char* rangeName) const
369{
370 return analyticalIntegral(code, rangeName, _histObsList, _pdfObsList, *_dataHist, false);
371}
372
373
374////////////////////////////////////////////////////////////////////////////////
375/// Return sampling hint for making curves of (projections) of this function
376/// as the recursive division strategy of RooCurve cannot deal efficiently
377/// with the vertical lines that occur in a non-interpolated histogram
378
379list<double>* RooHistPdf::plotSamplingHint(RooAbsRealLValue& obs, double xlo, double xhi) const
380{
381 // No hints are required when interpolation is used
382 if (_intOrder>0) {
383 return 0 ;
384 }
385
386 // Check that observable is in dataset, if not no hint is generated
387 RooAbsArg* dhObs = nullptr;
388 for (unsigned int i=0; i < _pdfObsList.size(); ++i) {
389 RooAbsArg* histObs = _histObsList[i];
390 RooAbsArg* pdfObs = _pdfObsList[i];
391 if (TString(obs.GetName())==pdfObs->GetName()) {
392 dhObs = _dataHist->get()->find(histObs->GetName()) ;
393 break;
394 }
395 }
396
397 if (!dhObs) {
398 return 0 ;
399 }
400 RooAbsLValue* lval = dynamic_cast<RooAbsLValue*>(dhObs) ;
401 if (!lval) {
402 return 0 ;
403 }
404
405 // Retrieve position of all bin boundaries
406
407 const RooAbsBinning* binning = lval->getBinningPtr(0) ;
408 double* boundaries = binning->array() ;
409
410 list<double>* hint = new list<double> ;
411
412 // Widen range slighty
413 xlo = xlo - 0.01*(xhi-xlo) ;
414 xhi = xhi + 0.01*(xhi-xlo) ;
415
416 double delta = (xhi-xlo)*1e-8 ;
417
418 // Construct array with pairs of points positioned epsilon to the left and
419 // right of the bin boundaries
420 for (Int_t i=0 ; i<binning->numBoundaries() ; i++) {
421 if (boundaries[i]>=xlo && boundaries[i]<=xhi) {
422 hint->push_back(boundaries[i]-delta) ;
423 hint->push_back(boundaries[i]+delta) ;
424 }
425 }
426
427 return hint ;
428}
429
430
431
432////////////////////////////////////////////////////////////////////////////////
433/// Return sampling hint for making curves of (projections) of this function
434/// as the recursive division strategy of RooCurve cannot deal efficiently
435/// with the vertical lines that occur in a non-interpolated histogram
436
437std::list<double>* RooHistPdf::binBoundaries(RooAbsRealLValue& obs, double xlo, double xhi) const
438{
439 // No hints are required when interpolation is used
440 if (_intOrder>0) {
441 return 0 ;
442 }
443
444 // Check that observable is in dataset, if not no hint is generated
445 RooAbsLValue* lvarg = dynamic_cast<RooAbsLValue*>(_dataHist->get()->find(obs.GetName())) ;
446 if (!lvarg) {
447 return 0 ;
448 }
449
450 // Retrieve position of all bin boundaries
451 const RooAbsBinning* binning = lvarg->getBinningPtr(0) ;
452 double* boundaries = binning->array() ;
453
454 list<double>* hint = new list<double> ;
455
456 // Construct array with pairs of points positioned epsilon to the left and
457 // right of the bin boundaries
458 for (Int_t i=0 ; i<binning->numBoundaries() ; i++) {
459 if (boundaries[i]>=xlo && boundaries[i]<=xhi) {
460 hint->push_back(boundaries[i]) ;
461 }
462 }
463
464 return hint ;
465}
466
467
468
469
470////////////////////////////////////////////////////////////////////////////////
471/// Only handle case of maximum in all variables
472
474{
476 if (common->getSize()==_pdfObsList.getSize()) {
477 delete common ;
478 return 1;
479 }
480 delete common ;
481 return 0 ;
482}
483
484
485////////////////////////////////////////////////////////////////////////////////
486
487double RooHistPdf::maxVal(Int_t code) const
488{
489 R__ASSERT(code==1) ;
490
491 double max(-1) ;
492 for (Int_t i=0 ; i<_dataHist->numEntries() ; i++) {
493 _dataHist->get(i) ;
494 double wgt = _dataHist->weight() ;
495 if (wgt>max) max=wgt ;
496 }
497
498 return max*1.05 ;
499}
500
501
502
503
504////////////////////////////////////////////////////////////////////////////////
505
507{
508 if (fabs(dh1.sumEntries()-dh2.sumEntries())>1e-8) return false ;
509 if (dh1.numEntries() != dh2.numEntries()) return false ;
510 for (int i=0 ; i < dh1.numEntries() ; i++) {
511 dh1.get(i) ;
512 dh2.get(i) ;
513 if (fabs(dh1.weight()-dh2.weight())>1e-8) return false ;
514 }
515 return true ;
516}
517
518
519
520////////////////////////////////////////////////////////////////////////////////
521/// Check if our datahist is already in the workspace
522
524{
525 std::list<RooAbsData*> allData = ws.allData() ;
526 std::list<RooAbsData*>::const_iterator iter ;
527 for (iter = allData.begin() ; iter != allData.end() ; ++iter) {
528 // If your dataset is already in this workspace nothing needs to be done
529 if (*iter == _dataHist) {
530 return false ;
531 }
532 }
533
534 // Check if dataset with given name already exists
535 RooAbsData* wsdata = ws.embeddedData(_dataHist->GetName()) ;
536
537 if (wsdata) {
538
539 // Yes it exists - now check if it is identical to our internal histogram
540 if (wsdata->InheritsFrom(RooDataHist::Class())) {
541
542 // Check if histograms are identical
543 if (areIdentical((RooDataHist&)*wsdata,*_dataHist)) {
544
545 // Exists and is of correct type, and identical -- adjust internal pointer to WS copy
546 _dataHist = (RooDataHist*) wsdata ;
547 } else {
548
549 // not identical, clone rename and import
550 TString uniqueName = Form("%s_%s",_dataHist->GetName(),GetName()) ;
551 bool flag = ws.import(*_dataHist,RooFit::Rename(uniqueName.Data()),RooFit::Embedded()) ;
552 if (flag) {
553 coutE(ObjectHandling) << " RooHistPdf::importWorkspaceHook(" << GetName() << ") unable to import clone of underlying RooDataHist with unique name " << uniqueName << ", abort" << endl ;
554 return true ;
555 }
556 _dataHist = (RooDataHist*) ws.embeddedData(uniqueName.Data()) ;
557 }
558
559 } else {
560
561 // Exists and is NOT of correct type: clone rename and import
562 TString uniqueName = Form("%s_%s",_dataHist->GetName(),GetName()) ;
563 bool flag = ws.import(*_dataHist,RooFit::Rename(uniqueName.Data()),RooFit::Embedded()) ;
564 if (flag) {
565 coutE(ObjectHandling) << " RooHistPdf::importWorkspaceHook(" << GetName() << ") unable to import clone of underlying RooDataHist with unique name " << uniqueName << ", abort" << endl ;
566 return true ;
567 }
568 _dataHist = (RooDataHist*) ws.embeddedData(uniqueName.Data()) ;
569
570 }
571 return false ;
572 }
573
574 // We need to import our datahist into the workspace
575 ws.import(*_dataHist,RooFit::Embedded()) ;
576
577 // Redirect our internal pointer to the copy in the workspace
578 _dataHist = (RooDataHist*) ws.embeddedData(_dataHist->GetName()) ;
579 return false ;
580}
581
582
583////////////////////////////////////////////////////////////////////////////////
584/// Stream an object of class RooHistPdf.
585
587{
588 if (R__b.IsReading()) {
590 // WVE - interim solution - fix proxies here
591 //_proxyList.Clear() ;
592 //registerProxy(_pdfObsList) ;
593 } else {
595 }
596}
597
#define e(i)
Definition: RSha256.hxx:103
#define coutE(a)
Definition: RooMsgService.h:37
#define ClassImp(name)
Definition: Rtypes.h:375
#define R__ASSERT(e)
Definition: TError.h:118
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:72
virtual void copyCache(const RooAbsArg *source, bool valueOnly=false, bool setValDirty=true)=0
virtual bool inRange(const char *) const
Definition: RooAbsArg.h:396
virtual void syncCache(const RooArgSet *nset=0)=0
RooAbsBinning is the abstract base class for RooRealVar binning definitions.
Definition: RooAbsBinning.h:26
virtual bool isParameterized() const
Interface function.
Definition: RooAbsBinning.h:80
virtual Int_t numBoundaries() const =0
virtual double * array() const =0
Int_t numTypes(const char *=0) const
Return number of types defined (in range named rangeName if rangeName!=0)
RooAbsCollection is an abstract container object that can hold multiple RooAbsArg objects.
virtual void removeAll()
Remove all arguments from our set, deleting them if we own them.
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.
virtual RooAbsArg * addClone(const RooAbsArg &var, bool silent=false)
Add a clone of the specified argument to list.
Storage_t::size_type size() const
bool selectCommon(const RooAbsCollection &refColl, RooAbsCollection &outColl) const
Create a subset of the current collection, consisting only of those elements that are contained as we...
RooAbsArg * find(const char *name) const
Find object with given name in list.
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition: RooAbsData.h:62
Abstract base class for objects that are lvalues, i.e.
Definition: RooAbsLValue.h:26
virtual const RooAbsBinning * getBinningPtr(const char *rangeName) const =0
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
virtual double getMax(const char *name=0) const
Get maximum of currently defined range.
virtual double getMin(const char *name=0) const
Get minimum of currently defined range.
const RooAbsBinning * getBinningPtr(const char *rangeName) const override
bool hasRange(const char *name) const override
Check if variable has a binning with given name.
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgList.h:22
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:57
RooCategory is an object to represent discrete states.
Definition: RooCategory.h:28
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...
The RooDataHist is a container class to hold N-dimensional binned data.
Definition: RooDataHist.h:45
double sum(bool correctForBinSize, bool inverseCorr=false) const
Return the sum of the weights of all bins in the histogram.
static TClass * Class()
double weight(std::size_t i) const
Return weight of i-th bin.
Definition: RooDataHist.h:110
double weightFast(const RooArgSet &bin, int intOrder, bool correctForBinSize, bool cdfBoundaries)
A faster version of RooDataHist::weight that assumes the passed arguments are aligned with the histog...
Int_t numEntries() const override
Return the number of bins.
const RooArgSet * get() const override
Get bin centre of current bin.
Definition: RooDataHist.h:82
double sumEntries() const override
Sum the weights of all bins.
RooHistPdf implements a probablity density function sampled from a multidimensional histogram.
Definition: RooHistPdf.h:29
void Streamer(TBuffer &) override
Stream an object of class RooHistPdf.
Definition: RooHistPdf.cxx:586
RooArgSet _histObsList
List of observables defining dimensions of histogram.
Definition: RooHistPdf.h:112
Int_t _intOrder
Interpolation order.
Definition: RooHistPdf.h:116
RooHistPdf()
Default constructor coverity[UNINIT_CTOR].
Definition: RooHistPdf.cxx:53
bool areIdentical(const RooDataHist &dh1, const RooDataHist &dh2)
Definition: RooHistPdf.cxx:506
RooDataHist * _dataHist
Unowned pointer to underlying histogram.
Definition: RooHistPdf.h:114
bool _cdfBoundaries
Use boundary conditions for CDFs.
Definition: RooHistPdf.h:117
std::list< double > * binBoundaries(RooAbsRealLValue &, double, double) const override
Return sampling hint for making curves of (projections) of this function as the recursive division st...
Definition: RooHistPdf.cxx:437
double totVolume() const
Return the total volume spanned by the observables of the RooHistPdf.
Definition: RooHistPdf.cxx:228
bool importWorkspaceHook(RooWorkspace &ws) override
Check if our datahist is already in the workspace.
Definition: RooHistPdf.cxx:523
static TClass * Class()
std::list< double > * plotSamplingHint(RooAbsRealLValue &obs, double xlo, double xhi) const override
Return sampling hint for making curves of (projections) of this function as the recursive division st...
Definition: RooHistPdf.cxx:379
RooSetProxy _pdfObsList
List of observables mapped onto histogram observables.
Definition: RooHistPdf.h:113
static double analyticalIntegral(Int_t code, const char *rangeName, RooArgSet const &histObsList, RooSetProxy const &pdfObsList, RooDataHist &dataHist, bool histFuncMode)
Definition: RooHistPdf.cxx:311
double maxVal(Int_t code) const override
Return maximum value for set of observables identified by code assigned in getMaxVal.
Definition: RooHistPdf.cxx:487
RooDataHist & dataHist()
Definition: RooHistPdf.h:38
Int_t getMaxVal(const RooArgSet &vars) const override
Only handle case of maximum in all variables.
Definition: RooHistPdf.cxx:473
double _totVolume
! Total volume of space (product of ranges of observables)
Definition: RooHistPdf.h:118
static Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName, RooArgSet const &histObsList, RooSetProxy const &pdfObsList, Int_t intOrder)
Definition: RooHistPdf.cxx:270
~RooHistPdf() override
Destructor.
Definition: RooHistPdf.cxx:189
double evaluate() const override
Return the current value: The value of the bin enclosing the current coordinates of the observables,...
Definition: RooHistPdf.cxx:203
bool _unitNorm
Assume contents is unit normalized (for use as pdf cache)
Definition: RooHistPdf.h:119
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:40
void setRange(const char *name, double min, double max)
Set a fit or plotting range.
Definition: RooRealVar.cxx:552
The RooWorkspace is a persistable container for RooFit projects.
Definition: RooWorkspace.h:43
Buffer base class used for serializing objects.
Definition: TBuffer.h:43
virtual Int_t ReadClassBuffer(const TClass *cl, void *pointer, const TClass *onfile_class=0)=0
Bool_t IsReading() const
Definition: TBuffer.h:86
virtual Int_t WriteClassBuffer(const TClass *cl, void *pointer)=0
const char * GetName() const override
Returns name of object.
Definition: TNamed.h:47
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition: TObject.cxx:445
Basic string class.
Definition: TString.h:136
const char * Data() const
Definition: TString.h:369
RooCmdArg Rename(const char *suffix)
RooCmdArg Embedded(bool flag=true)
Double_t y[n]
Definition: legend1.C:17
Double_t x[n]
Definition: legend1.C:17
const Int_t n
Definition: legend1.C:16
VecExpr< UnaryOp< Fabs< T >, VecExpr< A, T, D >, T >, T, D > fabs(const VecExpr< A, T, D > &rhs)
@ InputArguments
Definition: RooGlobalFunc.h:64
@ ObjectHandling
Definition: RooGlobalFunc.h:64
std::pair< double, double > getRangeOrBinningInterval(RooAbsArg const *arg, const char *rangeName)
Get the lower and upper bound of parameter range if arg can be casted to RooAbsRealLValue.
Definition: RooHelpers.cxx:168
void ws()
Definition: ws.C:66