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