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 
22 RooHistPdf implements a probablity density function sampled from a
23 multidimensional histogram. The histogram distribution is explicitly
24 normalized by RooHistPdf and can have an arbitrary number of real or
25 discrete dimensions.
26 **/
27 
28 #include "RooFit.h"
29 #include "Riostream.h"
30 
31 #include "RooHistPdf.h"
32 #include "RooDataHist.h"
33 #include "RooMsgService.h"
34 #include "RooRealVar.h"
35 #include "RooCategory.h"
36 #include "RooWorkspace.h"
37 #include "RooGlobalFunc.h"
38 
39 #include "TError.h"
40 #include "TBuffer.h"
41 
42 using namespace std;
43 
45 
46 
47 
48 
49 ////////////////////////////////////////////////////////////////////////////////
50 /// Default constructor
51 /// coverity[UNINIT_CTOR]
52 
53 RooHistPdf::RooHistPdf() : _dataHist(0), _totVolume(0), _unitNorm(kFALSE)
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 
65 RooHistPdf::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(kFALSE),
73  _totVolume(0),
74  _unitNorm(kFALSE)
75 {
76  _histObsList.addClone(vars) ;
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 
117 RooHistPdf::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(kFALSE),
125  _totVolume(0),
126  _unitNorm(kFALSE)
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 
169 RooHistPdf::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,kTRUE) ;
213  if (!harg->inRange(0)) {
214  return 0 ;
215  }
216  }
217  }
218 
220 // cout << "RooHistPdf::evaluate(" << GetName() << ") ret = " << ret << " ";
221 // cout << _histObsList[0] << " ";
222 // _histObsList[0]->Print("");
223 // _dataHist->Print("V");
224 // _dataHist->dump2();
225 
226  if (ret<0) {
227  ret=0 ;
228  }
229  return ret ;
230 }
231 
232 
233 ////////////////////////////////////////////////////////////////////////////////
234 /// Return the total volume spanned by the observables of the RooHistPdf
235 
237 {
238  // Return previously calculated value, if any
239  if (_totVolume>0) {
240  return _totVolume ;
241  }
242  _totVolume = 1. ;
243 
244  for (const auto arg : _histObsList) {
245  RooRealVar* real = dynamic_cast<RooRealVar*>(arg) ;
246  if (real) {
247  _totVolume *= (real->getMax()-real->getMin()) ;
248  } else {
249  RooCategory* cat = dynamic_cast<RooCategory*>(arg) ;
250  if (cat) {
251  _totVolume *= cat->numTypes() ;
252  }
253  }
254  }
255 
256  return _totVolume ;
257 }
258 
259 namespace {
260 bool fullRange(const RooAbsArg& x, const RooAbsArg& y ,const char* range)
261 {
262  const RooAbsRealLValue *_x = dynamic_cast<const RooAbsRealLValue*>(&x);
263  const RooAbsRealLValue *_y = dynamic_cast<const RooAbsRealLValue*>(&y);
264  if (!_x || !_y) return false;
265  if (!range || !strlen(range) || !_x->hasRange(range) ||
266  _x->getBinningPtr(range)->isParameterized()) {
267  // parameterized ranges may be full range now, but that might change,
268  // so return false
269  if (range && strlen(range) && _x->getBinningPtr(range)->isParameterized())
270  return false;
271  return (_x->getMin() == _y->getMin() && _x->getMax() == _y->getMax());
272  }
273  return (_x->getMin(range) == _y->getMin() && _x->getMax(range) == _y->getMax());
274 }
275 }
276 
277 
278 ////////////////////////////////////////////////////////////////////////////////
279 /// Determine integration scenario. If no interpolation is used,
280 /// RooHistPdf can perform all integrals over its dependents
281 /// analytically via partial or complete summation of the input
282 /// histogram. If interpolation is used on the integral over
283 /// all histogram observables is supported
284 
285 Int_t RooHistPdf::getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* rangeName) const
286 {
287  // First make list of pdf observables to histogram observables
288  // and select only those for which the integral is over the full range
289 
290  Int_t code = 0, frcode = 0;
291  for (unsigned int n=0; n < _pdfObsList.size() && n < _histObsList.size(); ++n) {
292  const auto pa = _pdfObsList[n];
293  const auto ha = _histObsList[n];
294 
295  if (allVars.find(*pa)) {
296  code |= 2 << n;
297  analVars.add(*pa);
298  if (fullRange(*pa, *ha, rangeName)) {
299  frcode |= 2 << n;
300  }
301  }
302  }
303 
304  if (code == frcode) {
305  // integrate over full range of all observables - use bit 0 to indicate
306  // full range integration over all observables
307  code |= 1;
308  }
309  // Disable partial analytical integrals if interpolation is used, and we
310  // integrate over sub-ranges, but leave them enabled when we integrate over
311  // the full range of one or several variables
312  if (_intOrder > 1 && !(code & 1)) {
313  analVars.removeAll();
314  return 0;
315  }
316  return (code >= 2) ? code : 0;
317 }
318 
319 
320 
321 ////////////////////////////////////////////////////////////////////////////////
322 /// Return integral identified by 'code'. The actual integration
323 /// is deferred to RooDataHist::sum() which implements partial
324 /// or complete summation over the histograms contents.
325 
326 Double_t RooHistPdf::analyticalIntegral(Int_t code, const char* rangeName) const
327 {
328  // Simplest scenario, full-range integration over all dependents
329  if (((2 << _histObsList.getSize()) - 1) == code) {
330  return _dataHist->sum(kFALSE);
331  }
332 
333  // Partial integration scenario, retrieve set of variables, calculate partial
334  // sum, figure out integration ranges (if needed)
335  RooArgSet intSet;
336  std::map<const RooAbsArg*, std::pair<Double_t, Double_t> > ranges;
337  for (unsigned int n=0; n < _pdfObsList.size() && n < _histObsList.size(); ++n) {
338  const auto pa = _pdfObsList[n];
339  const auto ha = _histObsList[n];
340 
341  if (code & (2 << n)) {
342  intSet.add(*ha);
343  }
344  if (!(code & 1)) {
345  RooAbsRealLValue* rlv = dynamic_cast<RooAbsRealLValue*>(pa);
346  if (rlv) {
347  const RooAbsBinning* binning = rlv->getBinningPtr(rangeName);
348  if (rangeName && rlv->hasRange(rangeName)) {
349  ranges[ha] = std::make_pair(
350  rlv->getMin(rangeName), rlv->getMax(rangeName));
351  } else if (binning) {
352  if (!binning->isParameterized()) {
353  ranges[ha] = std::make_pair(
354  binning->lowBound(), binning->highBound());
355  } else {
356  ranges[ha] = std::make_pair(
357  binning->lowBoundFunc()->getVal(), binning->highBoundFunc()->getVal());
358  }
359  }
360  }
361  }
362  // WVE must sync hist slice list values to pdf slice list
363  // Transfer values from
364  if (ha != pa) {
365  pa->syncCache();
366  ha->copyCache(pa,kTRUE);
367  }
368  }
369 
370  Double_t ret = (code & 1) ?
372  _dataHist->sum(intSet,_histObsList,kFALSE,kTRUE, ranges);
373 
374  // cout << "intSet = " << intSet << endl ;
375  // cout << "slice position = " << endl ;
376  // _histObsList.Print("v") ;
377  // cout << "RooHistPdf::ai(" << GetName() << ") code = " << code << " ret = " << ret << endl ;
378 
379  return ret ;
380 }
381 
382 
383 
384 ////////////////////////////////////////////////////////////////////////////////
385 /// Return sampling hint for making curves of (projections) of this function
386 /// as the recursive division strategy of RooCurve cannot deal efficiently
387 /// with the vertical lines that occur in a non-interpolated histogram
388 
389 list<Double_t>* RooHistPdf::plotSamplingHint(RooAbsRealLValue& obs, Double_t xlo, Double_t xhi) const
390 {
391  // No hints are required when interpolation is used
392  if (_intOrder>0) {
393  return 0 ;
394  }
395 
396  // Check that observable is in dataset, if not no hint is generated
397  RooAbsArg* dhObs = nullptr;
398  for (unsigned int i=0; i < _pdfObsList.size(); ++i) {
399  RooAbsArg* histObs = _histObsList[i];
400  RooAbsArg* pdfObs = _pdfObsList[i];
401  if (TString(obs.GetName())==pdfObs->GetName()) {
402  dhObs = _dataHist->get()->find(histObs->GetName()) ;
403  break;
404  }
405  }
406 
407  if (!dhObs) {
408  return 0 ;
409  }
410  RooAbsLValue* lval = dynamic_cast<RooAbsLValue*>(dhObs) ;
411  if (!lval) {
412  return 0 ;
413  }
414 
415  // Retrieve position of all bin boundaries
416 
417  const RooAbsBinning* binning = lval->getBinningPtr(0) ;
418  Double_t* boundaries = binning->array() ;
419 
420  list<Double_t>* hint = new list<Double_t> ;
421 
422  // Widen range slighty
423  xlo = xlo - 0.01*(xhi-xlo) ;
424  xhi = xhi + 0.01*(xhi-xlo) ;
425 
426  Double_t delta = (xhi-xlo)*1e-8 ;
427 
428  // Construct array with pairs of points positioned epsilon to the left and
429  // right of the bin boundaries
430  for (Int_t i=0 ; i<binning->numBoundaries() ; i++) {
431  if (boundaries[i]>=xlo && boundaries[i]<=xhi) {
432  hint->push_back(boundaries[i]-delta) ;
433  hint->push_back(boundaries[i]+delta) ;
434  }
435  }
436 
437  return hint ;
438 }
439 
440 
441 
442 ////////////////////////////////////////////////////////////////////////////////
443 /// Return sampling hint for making curves of (projections) of this function
444 /// as the recursive division strategy of RooCurve cannot deal efficiently
445 /// with the vertical lines that occur in a non-interpolated histogram
446 
447 std::list<Double_t>* RooHistPdf::binBoundaries(RooAbsRealLValue& obs, Double_t xlo, Double_t xhi) const
448 {
449  // No hints are required when interpolation is used
450  if (_intOrder>0) {
451  return 0 ;
452  }
453 
454  // Check that observable is in dataset, if not no hint is generated
455  RooAbsLValue* lvarg = dynamic_cast<RooAbsLValue*>(_dataHist->get()->find(obs.GetName())) ;
456  if (!lvarg) {
457  return 0 ;
458  }
459 
460  // Retrieve position of all bin boundaries
461  const RooAbsBinning* binning = lvarg->getBinningPtr(0) ;
462  Double_t* boundaries = binning->array() ;
463 
464  list<Double_t>* hint = new list<Double_t> ;
465 
466  // Construct array with pairs of points positioned epsilon to the left and
467  // right of the bin boundaries
468  for (Int_t i=0 ; i<binning->numBoundaries() ; i++) {
469  if (boundaries[i]>=xlo && boundaries[i]<=xhi) {
470  hint->push_back(boundaries[i]) ;
471  }
472  }
473 
474  return hint ;
475 }
476 
477 
478 
479 
480 ////////////////////////////////////////////////////////////////////////////////
481 /// Only handle case of maximum in all variables
482 
484 {
485  RooAbsCollection* common = _pdfObsList.selectCommon(vars) ;
486  if (common->getSize()==_pdfObsList.getSize()) {
487  delete common ;
488  return 1;
489  }
490  delete common ;
491  return 0 ;
492 }
493 
494 
495 ////////////////////////////////////////////////////////////////////////////////
496 
498 {
499  R__ASSERT(code==1) ;
500 
501  Double_t max(-1) ;
502  for (Int_t i=0 ; i<_dataHist->numEntries() ; i++) {
503  _dataHist->get(i) ;
504  Double_t wgt = _dataHist->weight() ;
505  if (wgt>max) max=wgt ;
506  }
507 
508  return max*1.05 ;
509 }
510 
511 
512 
513 
514 ////////////////////////////////////////////////////////////////////////////////
515 
517 {
518  if (fabs(dh1.sumEntries()-dh2.sumEntries())>1e-8) return kFALSE ;
519  if (dh1.numEntries() != dh2.numEntries()) return kFALSE ;
520  for (int i=0 ; i < dh1.numEntries() ; i++) {
521  dh1.get(i) ;
522  dh2.get(i) ;
523  if (fabs(dh1.weight()-dh2.weight())>1e-8) return kFALSE ;
524  }
525  return kTRUE ;
526 }
527 
528 
529 
530 ////////////////////////////////////////////////////////////////////////////////
531 /// Check if our datahist is already in the workspace
532 
534 {
535  std::list<RooAbsData*> allData = ws.allData() ;
536  std::list<RooAbsData*>::const_iterator iter ;
537  for (iter = allData.begin() ; iter != allData.end() ; ++iter) {
538  // If your dataset is already in this workspace nothing needs to be done
539  if (*iter == _dataHist) {
540  return kFALSE ;
541  }
542  }
543 
544  // Check if dataset with given name already exists
545  RooAbsData* wsdata = ws.embeddedData(_dataHist->GetName()) ;
546 
547  if (wsdata) {
548 
549  // Yes it exists - now check if it is identical to our internal histogram
550  if (wsdata->InheritsFrom(RooDataHist::Class())) {
551 
552  // Check if histograms are identical
553  if (areIdentical((RooDataHist&)*wsdata,*_dataHist)) {
554 
555  // Exists and is of correct type, and identical -- adjust internal pointer to WS copy
556  _dataHist = (RooDataHist*) wsdata ;
557  } else {
558 
559  // not identical, clone rename and import
560  TString uniqueName = Form("%s_%s",_dataHist->GetName(),GetName()) ;
561  Bool_t flag = ws.import(*_dataHist,RooFit::Rename(uniqueName.Data()),RooFit::Embedded()) ;
562  if (flag) {
563  coutE(ObjectHandling) << " RooHistPdf::importWorkspaceHook(" << GetName() << ") unable to import clone of underlying RooDataHist with unique name " << uniqueName << ", abort" << endl ;
564  return kTRUE ;
565  }
566  _dataHist = (RooDataHist*) ws.embeddedData(uniqueName.Data()) ;
567  }
568 
569  } else {
570 
571  // Exists and is NOT of correct type: clone rename and import
572  TString uniqueName = Form("%s_%s",_dataHist->GetName(),GetName()) ;
573  Bool_t flag = ws.import(*_dataHist,RooFit::Rename(uniqueName.Data()),RooFit::Embedded()) ;
574  if (flag) {
575  coutE(ObjectHandling) << " RooHistPdf::importWorkspaceHook(" << GetName() << ") unable to import clone of underlying RooDataHist with unique name " << uniqueName << ", abort" << endl ;
576  return kTRUE ;
577  }
578  _dataHist = (RooDataHist*) ws.embeddedData(uniqueName.Data()) ;
579 
580  }
581  return kFALSE ;
582  }
583 
584  // We need to import our datahist into the workspace
585  ws.import(*_dataHist,RooFit::Embedded()) ;
586 
587  // Redirect our internal pointer to the copy in the workspace
588  _dataHist = (RooDataHist*) ws.embeddedData(_dataHist->GetName()) ;
589  return kFALSE ;
590 }
591 
592 
593 ////////////////////////////////////////////////////////////////////////////////
594 /// Stream an object of class RooHistPdf.
595 
596 void RooHistPdf::Streamer(TBuffer &R__b)
597 {
598  if (R__b.IsReading()) {
599  R__b.ReadClassBuffer(RooHistPdf::Class(),this);
600  // WVE - interim solution - fix proxies here
601  //_proxyList.Clear() ;
602  //registerProxy(_pdfObsList) ;
603  } else {
604  R__b.WriteClassBuffer(RooHistPdf::Class(),this);
605  }
606 }
607 
RooHistPdf::_cdfBoundaries
Bool_t _cdfBoundaries
Definition: RooHistPdf.h:103
n
const Int_t n
Definition: legend1.C:16
RooWorkspace.h
RooHistPdf::maxVal
virtual Double_t maxVal(Int_t code) const
Return maximum value for set of observables identified by code assigned in getMaxVal.
Definition: RooHistPdf.cxx:497
kTRUE
const Bool_t kTRUE
Definition: RtypesCore.h:91
e
#define e(i)
Definition: RSha256.hxx:103
RooArgSet::addClone
RooAbsArg * addClone(const RooAbsArg &var, Bool_t silent=kFALSE) override
Add clone of specified element to an owning set.
Definition: RooArgSet.cxx:288
RooMsgService.h
RooHistPdf::totVolume
Double_t totVolume() const
Return the total volume spanned by the observables of the RooHistPdf.
Definition: RooHistPdf.cxx:236
RooAbsData
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition: RooAbsData.h:49
RooAbsBinning::highBoundFunc
virtual RooAbsReal * highBoundFunc() const
Return pointer to RooAbsReal parameterized upper bound, if any.
Definition: RooAbsBinning.h:88
RooAbsRealLValue::getMax
virtual Double_t getMax(const char *name=0) const
Get maximum of currently defined range.
Definition: RooAbsRealLValue.h:89
RooFit.h
RooFit::InputArguments
@ InputArguments
Definition: RooGlobalFunc.h:68
RooAbsBinning::lowBoundFunc
virtual RooAbsReal * lowBoundFunc() const
Return pointer to RooAbsReal parameterized lower bound, if any.
Definition: RooAbsBinning.h:84
RooHistPdf::getAnalyticalIntegral
Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const
Determine integration scenario.
Definition: RooHistPdf.cxx:285
RooAbsBinning::highBound
virtual Double_t highBound() const =0
RooHistPdf::RooHistPdf
RooHistPdf()
Default constructor coverity[UNINIT_CTOR].
Definition: RooHistPdf.cxx:53
TString::Data
const char * Data() const
Definition: TString.h:369
ClassImp
#define ClassImp(name)
Definition: Rtypes.h:364
Form
char * Form(const char *fmt,...)
TBuffer::ReadClassBuffer
virtual Int_t ReadClassBuffer(const TClass *cl, void *pointer, const TClass *onfile_class=0)=0
coutE
#define coutE(a)
Definition: RooMsgService.h:33
RooArgList
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgList.h:21
RooAbsReal::getVal
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:91
RooHistPdf::analyticalIntegral
Double_t analyticalIntegral(Int_t code, const char *rangeName=0) const
Return integral identified by 'code'.
Definition: RooHistPdf.cxx:326
RooDataHist::get
virtual const RooArgSet * get() const
Definition: RooDataHist.h:78
RooAbsLValue::getBinningPtr
virtual const RooAbsBinning * getBinningPtr(const char *rangeName) const =0
RooAbsCollection::find
RooAbsArg * find(const char *name) const
Find object with given name in list.
Definition: RooAbsCollection.cxx:813
RooDataHist::weight
virtual Double_t weight() const
Definition: RooDataHist.h:100
x
Double_t x[n]
Definition: legend1.C:17
TBuffer
Buffer base class used for serializing objects.
Definition: TBuffer.h:43
RooHistPdf::_unitNorm
Bool_t _unitNorm
Total volume of space (product of ranges of observables)
Definition: RooHistPdf.h:105
TString
Basic string class.
Definition: TString.h:136
RooRealVar::setRange
void setRange(const char *name, Double_t min, Double_t max)
Set a fit or plotting range.
Definition: RooRealVar.cxx:526
TObject::InheritsFrom
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition: TObject.cxx:445
RooHistPdf::_histObsList
RooArgSet _histObsList
Definition: RooHistPdf.h:98
RooAbsCategory::numTypes
Int_t numTypes(const char *=0) const
Return number of types defined (in range named rangeName if rangeName!=0)
Definition: RooAbsCategory.h:128
bool
RooAbsRealLValue::getMin
virtual Double_t getMin(const char *name=0) const
Get miniminum of currently defined range.
Definition: RooAbsRealLValue.h:86
RooAbsBinning::array
virtual Double_t * array() const =0
RooAbsRealLValue::getBinningPtr
virtual const RooAbsBinning * getBinningPtr(const char *rangeName) const
Definition: RooAbsRealLValue.h:59
ROOT::Math::fabs
VecExpr< UnaryOp< Fabs< T >, VecExpr< A, T, D >, T >, T, D > fabs(const VecExpr< A, T, D > &rhs)
Definition: UnaryOperators.h:131
RooFit::Rename
RooCmdArg Rename(const char *suffix)
Definition: RooGlobalFunc.cxx:321
RooDataHist
The RooDataHist is a container class to hold N-dimensional binned data.
Definition: RooDataHist.h:39
ws
void ws()
Definition: ws.C:66
TBuffer.h
RooHistPdf::_intOrder
Int_t _intOrder
Auxiliary class keeping tracking of analytical integration code.
Definition: RooHistPdf.h:102
RooAbsArg::syncCache
virtual void syncCache(const RooArgSet *nset=0)=0
RooDataHist::sumEntries
virtual Double_t sumEntries() const
Return effective number of entries in dataset, i.e., sum all weights.
Definition: RooDataHist.cxx:1704
TBuffer::WriteClassBuffer
virtual Int_t WriteClassBuffer(const TClass *cl, void *pointer)=0
kFALSE
const Bool_t kFALSE
Definition: RtypesCore.h:92
RooDataHist.h
RooArgSet::add
Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE) override
Add element to non-owning set.
Definition: RooArgSet.cxx:261
RooAbsCollection::selectCommon
RooAbsCollection * selectCommon(const RooAbsCollection &refColl) const
Create a subset of the current collection, consisting only of those elements that are contained as we...
Definition: RooAbsCollection.cxx:700
RooAbsBinning
RooAbsBinning is the abstract base class for RooRealVar binning definitions.
Definition: RooAbsBinning.h:26
RooAbsCollection
RooAbsCollection is an abstract container object that can hold multiple RooAbsArg objects.
Definition: RooAbsCollection.h:31
RooFit::Embedded
RooCmdArg Embedded(Bool_t flag=kTRUE)
Definition: RooGlobalFunc.cxx:322
RooHistPdf::_dataHist
RooDataHist * _dataHist
Definition: RooHistPdf.h:100
RooAbsCollection::size
Storage_t::size_type size() const
Definition: RooAbsCollection.h:155
RooCategory.h
y
Double_t y[n]
Definition: legend1.C:17
RooFit::ObjectHandling
@ ObjectHandling
Definition: RooGlobalFunc.h:68
RooHistPdf::_pdfObsList
RooSetProxy _pdfObsList
Definition: RooHistPdf.h:99
RooHistPdf::importWorkspaceHook
Bool_t importWorkspaceHook(RooWorkspace &ws)
Check if our datahist is already in the workspace.
Definition: RooHistPdf.cxx:533
RooRealVar.h
RooHistPdf
RooHistPdf implements a probablity density function sampled from a multidimensional histogram.
Definition: RooHistPdf.h:29
RooGlobalFunc.h
RooAbsArg::copyCache
virtual void copyCache(const RooAbsArg *source, Bool_t valueOnly=kFALSE, Bool_t setValDirty=kTRUE)=0
RooHistPdf.h
TBuffer::IsReading
Bool_t IsReading() const
Definition: TBuffer.h:86
RooWorkspace
The RooWorkspace is a persistable container for RooFit projects.
Definition: RooWorkspace.h:43
RooDataHist::sum
Double_t sum(Bool_t correctForBinSize, Bool_t inverseCorr=kFALSE) const
Return the sum of the weights of all hist bins.
Definition: RooDataHist.cxx:1431
RooHistPdf::evaluate
Double_t evaluate() const
Return the current value: The value of the bin enclosing the current coordinates of the observables,...
Definition: RooHistPdf.cxx:203
RooAbsRealLValue::hasRange
virtual Bool_t hasRange(const char *name) const
Check if variable has a binning with given name.
Definition: RooAbsRealLValue.h:102
Double_t
double Double_t
Definition: RtypesCore.h:59
R__ASSERT
#define R__ASSERT(e)
Definition: TError.h:120
RooHistPdf::binBoundaries
virtual std::list< Double_t > * binBoundaries(RooAbsRealLValue &, Double_t, Double_t) const
Return sampling hint for making curves of (projections) of this function as the recursive division st...
Definition: RooHistPdf.cxx:447
RooCategory
RooCategory is an object to represent discrete states.
Definition: RooCategory.h:27
RooDataHist::numEntries
virtual Int_t numEntries() const
Return the number of bins.
Definition: RooDataHist.cxx:1695
RooHistPdf::areIdentical
Bool_t areIdentical(const RooDataHist &dh1, const RooDataHist &dh2)
Definition: RooHistPdf.cxx:516
RooAbsBinning::numBoundaries
virtual Int_t numBoundaries() const =0
name
char name[80]
Definition: TGX11.cxx:110
RooAbsBinning::isParameterized
virtual Bool_t isParameterized() const
Interface function.
Definition: RooAbsBinning.h:80
RooHistPdf::_totVolume
Double_t _totVolume
Definition: RooHistPdf.h:104
RooAbsArg
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition: RooAbsArg.h:72
RooHistPdf::plotSamplingHint
virtual std::list< Double_t > * plotSamplingHint(RooAbsRealLValue &obs, Double_t xlo, Double_t xhi) const
Return sampling hint for making curves of (projections) of this function as the recursive division st...
Definition: RooHistPdf.cxx:389
RooAbsPdf
Definition: RooAbsPdf.h:43
RooAbsLValue
Abstract base class for objects that are lvalues, i.e.
Definition: RooAbsLValue.h:26
TNamed::GetName
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
RooHistPdf::getMaxVal
virtual Int_t getMaxVal(const RooArgSet &vars) const
Only handle case of maximum in all variables.
Definition: RooHistPdf.cxx:483
Class
void Class()
Definition: Class.C:29
RooRealVar
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:39
RooAbsCollection::removeAll
virtual void removeAll()
Remove all arguments from our set, deleting them if we own them.
Definition: RooAbsCollection.cxx:643
Riostream.h
RooAbsRealLValue
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
Definition: RooAbsRealLValue.h:31
RooAbsBinning::lowBound
virtual Double_t lowBound() const =0
RooAbsArg::inRange
virtual Bool_t inRange(const char *) const
Definition: RooAbsArg.h:376
RooAbsCollection::getSize
Int_t getSize() const
Definition: RooAbsCollection.h:172
RooSetProxy::add
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE) override
Overloaded RooArgSet::add() method inserts 'var' into set and registers 'var' as server to owner with...
Definition: RooSetProxy.cxx:165
RooArgSet
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:29
RooHistPdf::~RooHistPdf
virtual ~RooHistPdf()
Destructor.
Definition: RooHistPdf.cxx:189
int
TError.h