Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
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
22A propability 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 "RooCategory.h"
31#include "RooCurve.h"
32#include "RooDataHist.h"
33#include "RooFitImplHelpers.h"
34#include "RooGlobalFunc.h"
35#include "RooHistPdf.h"
36#include "RooMsgService.h"
37#include "RooRealVar.h"
38#include "RooUniformBinning.h"
39#include "RooWorkspace.h"
40
41#include "TError.h"
42#include "TBuffer.h"
43
44
46
47
48////////////////////////////////////////////////////////////////////////////////
49/// Constructor from a RooDataHist. RooDataHist dimensions
50/// can be either real or discrete. See RooDataHist::RooDataHist for details on the binning.
51/// RooHistPdf neither owns or clone 'dhist' and the user must ensure the input histogram exists
52/// for the entire life span of this PDF.
53
54RooHistPdf::RooHistPdf(const char *name, const char *title, const RooArgSet& vars,
55 const RooDataHist& dhist, Int_t intOrder) :
56 RooAbsPdf(name,title),
57 _pdfObsList("pdfObs","List of p.d.f. observables",this),
58 _dataHist(const_cast<RooDataHist*>(&dhist)),
59 _codeReg(10),
60 _intOrder(intOrder)
61{
63 _pdfObsList.add(vars) ;
64
65 // Verify that vars and dhist.get() have identical contents
66 const RooArgSet* dvars = dhist.get() ;
67 if (vars.size()!=dvars->size()) {
68 coutE(InputArguments) << "RooHistPdf::ctor(" << GetName()
69 << ") ERROR variable list and RooDataHist must contain the same variables." << std::endl ;
70 assert(0) ;
71 }
72 for (const auto arg : vars) {
73 if (!dvars->find(arg->GetName())) {
74 coutE(InputArguments) << "RooHistPdf::ctor(" << GetName()
75 << ") ERROR variable list and RooDataHist must contain the same variables." << std::endl ;
76 assert(0) ;
77 }
78 }
79
80
81 // Adjust ranges of _histObsList to those of _dataHist
82 for (const auto hobs : _histObsList) {
83 // Guaranteed to succeed, since checked above in constructor
84 RooAbsArg* dhobs = dhist.get()->find(hobs->GetName()) ;
85 RooRealVar* dhreal = dynamic_cast<RooRealVar*>(dhobs) ;
86 if (dhreal){
87 (static_cast<RooRealVar*>(hobs))->setRange(dhreal->getMin(),dhreal->getMax()) ;
88 }
89 }
90
91}
92
93
94
95
96////////////////////////////////////////////////////////////////////////////////
97/// Constructor from a RooDataHist. The first list of observables are the p.d.f.
98/// observables, which may any RooAbsReal (function or variable). The second list
99/// are the corresponding observables in the RooDataHist which must be of type
100/// RooRealVar or RooCategory This constructor thus allows to apply a coordinate transformation
101/// on the histogram data to be applied.
102
103RooHistPdf::RooHistPdf(const char *name, const char *title, const RooArgList& pdfObs,
104 const RooArgList& histObs, const RooDataHist& dhist, Int_t intOrder) :
105 RooAbsPdf(name,title),
106 _pdfObsList("pdfObs","List of p.d.f. observables",this),
107 _dataHist(const_cast<RooDataHist*>(&dhist)),
108 _codeReg(10),
109 _intOrder(intOrder)
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 constructor
139 RooAbsArg* dhobs = dhist.get()->find(hobs->GetName()) ;
140 RooRealVar* dhreal = dynamic_cast<RooRealVar*>(dhobs) ;
141 if (dhreal){
142 (static_cast<RooRealVar*>(hobs))->setRange(dhreal->getMin(),dhreal->getMax()) ;
143 }
144 }
145}
146
147RooHistPdf::RooHistPdf(const char *name, const char *title, const RooArgSet &vars, std::unique_ptr<RooDataHist> dhist,
148 int intOrder)
149 : RooHistPdf{name, title, vars, *dhist, intOrder}
150{
151 initializeOwnedDataHist(std::move(dhist));
152}
153RooHistPdf::RooHistPdf(const char *name, const char *title, const RooArgList &pdfObs, const RooArgList &histObs,
154 std::unique_ptr<RooDataHist> dhist, int intOrder)
155 : RooHistPdf{name, title, pdfObs, histObs, *dhist, intOrder}
156{
157 initializeOwnedDataHist(std::move(dhist));
158}
159
160
161////////////////////////////////////////////////////////////////////////////////
162/// Copy constructor
163
164RooHistPdf::RooHistPdf(const RooHistPdf& other, const char* name) :
165 RooAbsPdf(other,name),
166 _pdfObsList("pdfObs",this,other._pdfObsList),
167 _dataHist(other._dataHist),
168 _codeReg(other._codeReg),
169 _intOrder(other._intOrder),
170 _cdfBoundaries(other._cdfBoundaries),
171 _totVolume(other._totVolume),
172 _unitNorm(other._unitNorm)
173{
175}
176
178 if (_ownedDataHist) return _ownedDataHist.get();
179 _ownedDataHist.reset(static_cast<RooDataHist*>(_dataHist->Clone(newname)));
181 return _dataHist;
182}
183
185{
186 std::span<double> output = ctx.output();
187
188 // For interpolation and histograms of higher dimension, use base function
189 if (_pdfObsList.size() > 1) {
191 return;
192 }
193
194 auto xVals = ctx.at(_pdfObsList[0]);
195 _dataHist->weights(output.data(), xVals, _intOrder, true, _cdfBoundaries);
196}
197
198
199////////////////////////////////////////////////////////////////////////////////
200/// Return the current value: The value of the bin enclosing the current coordinates
201/// of the observables, normalized by the histograms contents. Interpolation
202/// is applied if the RooHistPdf is configured to do that.
203
205{
206 // Transfer values from
207 for (unsigned int i=0; i < _pdfObsList.size(); ++i) {
208 RooAbsArg* harg = _histObsList[i];
209 RooAbsArg* parg = _pdfObsList[i];
210
211 if (harg != parg) {
212 parg->syncCache() ;
213 harg->copyCache(parg,true) ;
214 if (!harg->inRange(nullptr)) {
215 return 0 ;
216 }
217 }
218 }
219
221
222 return std::max(ret, 0.0);
223}
224
226 RooDataHist const *dataHist, const RooArgSet &obs, bool correctForBinSize)
227{
228 if (intOrder != 0) {
229 ooccoutE(klass, InputArguments) << "RooHistPdf::weight(" << klass->GetName()
230 << ") ERROR: Code Squashing currently only supports non-interpolation cases."
231 << std::endl;
232 return;
233 }
234
235 std::string const &idxName = dataHist->calculateTreeIndexForCodeSquash(klass, ctx, obs);
236 std::string const &weightName = dataHist->declWeightArrayForCodeSquash(klass, ctx, correctForBinSize);
237 std::string res = weightName;
238 if (weightName[0] == '_')
239 res += "[" + idxName + "]";
240 ctx.addResult(klass, res);
241}
242
244{
246}
247
248////////////////////////////////////////////////////////////////////////////////
249/// Return the total volume spanned by the observables of the RooHistPdf
250
252{
253 // Return previously calculated value, if any
254 if (_totVolume>0) {
255 return _totVolume ;
256 }
257 _totVolume = 1. ;
258
259 for (const auto arg : _histObsList) {
260 RooRealVar* real = dynamic_cast<RooRealVar*>(arg) ;
261 if (real) {
262 _totVolume *= (real->getMax()-real->getMin()) ;
263 } else {
264 RooCategory* cat = dynamic_cast<RooCategory*>(arg) ;
265 if (cat) {
266 _totVolume *= cat->numTypes() ;
267 }
268 }
269 }
270
271 return _totVolume ;
272}
273
274namespace {
275
276bool fullRange(const RooAbsArg& x, const RooAbsArg& y ,const char* range)
277{
278 const RooAbsRealLValue *_x = dynamic_cast<const RooAbsRealLValue*>(&x);
279 const RooAbsRealLValue *_y = dynamic_cast<const RooAbsRealLValue*>(&y);
280 if (!_x || !_y) return false;
281 if (!range || !strlen(range) || !_x->hasRange(range) ||
282 _x->getBinningPtr(range)->isParameterized()) {
283 // parameterized ranges may be full range now, but that might change,
284 // so return false
285 if (range && strlen(range) && _x->getBinningPtr(range)->isParameterized())
286 return false;
287 return (_x->getMin() == _y->getMin() && _x->getMax() == _y->getMax());
288 }
289 return (_x->getMin(range) == _y->getMin() && _x->getMax(range) == _y->getMax());
290}
291
292bool okayForAnalytical(RooAbsArg const& obs, RooArgSet const& allVars)
293{
294 auto lobs = dynamic_cast<RooAbsRealLValue const*>(&obs);
295 if(lobs == nullptr) return false;
296
297 bool isOkayForAnalyticalInt = false;
298
299 for(RooAbsArg *var : allVars) {
300 if(obs.dependsOn(*var)) {
301 if(!lobs->isJacobianOK(*var)) return false;
302 isOkayForAnalyticalInt = true;
303 }
304 }
305
306 return isOkayForAnalyticalInt;
307}
308
309} // namespace
310
311
313 RooArgSet& analVars,
314 const char* rangeName,
315 RooArgSet const& histObsList,
316 RooArgSet const& pdfObsList,
317 Int_t intOrder)
318{
319 // First make list of pdf observables to histogram observables
320 // and select only those for which the integral is over the full range
321
322 Int_t code = 0;
323 Int_t frcode = 0;
324 for (unsigned int n=0; n < pdfObsList.size() && n < histObsList.size(); ++n) {
325 const auto pa = pdfObsList[n];
326 const auto ha = histObsList[n];
327
328 if (okayForAnalytical(*pa, allVars)) {
329 code |= 2 << n;
330 analVars.add(*pa);
331 if (fullRange(*pa, *ha, rangeName)) {
332 frcode |= 2 << n;
333 }
334 }
335 }
336
337 if (code == frcode) {
338 // integrate over full range of all observables - use bit 0 to indicate
339 // full range integration over all observables
340 code |= 1;
341 }
342
343 // Disable partial analytical integrals if interpolation is used, and we
344 // integrate over sub-ranges, but leave them enabled when we integrate over
345 // the full range of one or several variables
346 if (intOrder > 1 && !(code & 1)) {
347 analVars.removeAll();
348 return 0;
349 }
350 return (code >= 2) ? code : 0;
351}
352
353
355 const char* rangeName,
356 RooArgSet const& histObsList,
357 RooArgSet const& pdfObsList,
358 RooDataHist& dataHist,
359 bool histFuncMode) {
360 // Simplest scenario, full-range integration over all dependents
361 if (((2 << histObsList.size()) - 1) == code) {
362 return dataHist.sum(histFuncMode);
363 }
364
365 // Partial integration scenario, retrieve set of variables, calculate partial
366 // sum, figure out integration ranges (if needed)
367 RooArgSet intSet;
368 std::map<const RooAbsArg*, std::pair<double, double> > ranges;
369 for (unsigned int n=0; n < pdfObsList.size() && n < histObsList.size(); ++n) {
370 const auto pa = pdfObsList[n];
371 const auto ha = histObsList[n];
372
373 if (code & (2 << n)) {
374 intSet.add(*ha);
375 }
376 if (!(code & 1)) {
377 ranges[ha] = RooHelpers::getRangeOrBinningInterval(pa, rangeName);
378 }
379 // WVE must sync hist slice list values to pdf slice list
380 // Transfer values from
381 if (ha != pa) {
382 pa->syncCache();
383 ha->copyCache(pa,true);
384 }
385 }
386
387 double ret = (code & 1) ? dataHist.sum(intSet,histObsList,true,!histFuncMode) :
388 dataHist.sum(intSet,histObsList,true,!histFuncMode, ranges);
389
390 return ret ;
391}
392
393std::string RooHistPdf::rooHistIntegralTranslateImpl(int code, RooAbsArg const *klass, RooDataHist const *dataHist,
394 const RooArgSet &obs, bool histFuncMode)
395{
396 if (((2 << obs.size()) - 1) != code) {
397 oocoutE(klass, InputArguments)
398 << "RooHistPdf::integral(" << klass->GetName()
399 << ") ERROR: AD currently only supports integrating over all histogram observables." << std::endl;
400 return "";
401 }
402 return std::to_string(dataHist->sum(histFuncMode));
403}
404
405std::string RooHistPdf::buildCallToAnalyticIntegral(int code, const char * /*rangeName */,
406 RooFit::Detail::CodeSquashContext & /* ctx */) const
407{
408 return rooHistIntegralTranslateImpl(code, this, _dataHist, _pdfObsList, false);
409}
410
411////////////////////////////////////////////////////////////////////////////////
412/// Determine integration scenario. If no interpolation is used,
413/// RooHistPdf can perform all integrals over its dependents
414/// analytically via partial or complete summation of the input
415/// histogram. If interpolation is used on the integral over
416/// all histogram observables is supported
417
418Int_t RooHistPdf::getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* rangeName) const
419{
420 return getAnalyticalIntegral(allVars, analVars, rangeName, _histObsList, _pdfObsList, _intOrder);
421}
422
423
424////////////////////////////////////////////////////////////////////////////////
425/// Return integral identified by 'code'. The actual integration
426/// is deferred to RooDataHist::sum() which implements partial
427/// or complete summation over the histograms contents.
428
429double RooHistPdf::analyticalIntegral(Int_t code, const char* rangeName) const
430{
431 return analyticalIntegral(code, rangeName, _histObsList, _pdfObsList, *_dataHist, false);
432}
433
434
435bool RooHistPdf::forceAnalyticalInt(RooArgSet const& pdfObsList, const RooAbsArg& dep)
436{
437 bool isOkayForAnalyticalInt = false;
438
439 for (RooAbsArg * obs : pdfObsList) {
440 if(obs->dependsOn(dep)) {
441 // If the observable doesn't depend linearly on the integration
442 // variable we will not do analytical integration.
443 auto lvalue = dynamic_cast<RooAbsRealLValue const*>(obs);
444 if(!(lvalue && lvalue->isJacobianOK(dep))) return false;
445 isOkayForAnalyticalInt = true;
446 }
447 }
448
449 return isOkayForAnalyticalInt;
450}
451
452
454{
455 return forceAnalyticalInt(_pdfObsList, dep);
456}
457
458
459////////////////////////////////////////////////////////////////////////////////
460/// Return sampling hint for making curves of (projections) of this function
461/// as the recursive division strategy of RooCurve cannot deal efficiently
462/// with the vertical lines that occur in a non-interpolated histogram
463
464std::list<double>* RooHistPdf::plotSamplingHint(RooAbsRealLValue& obs, double xlo, double xhi) const
465{
467}
468
469
470std::list<double>* RooHistPdf::plotSamplingHint(RooDataHist const& dataHist,
471 RooArgSet const& pdfObsList,
472 RooArgSet const& histObsList,
473 int intOrder,
474 RooAbsRealLValue& obs,
475 double xlo,
476 double xhi)
477{
478 // No hints are required when interpolation is used
479 if (intOrder>0) {
480 return nullptr;
481 }
482
483 // Check that observable is in dataset, if not no hint is generated
484 RooAbsArg* dhObs = nullptr;
485 for (unsigned int i=0; i < pdfObsList.size(); ++i) {
486 RooAbsArg* histObs = histObsList[i];
487 RooAbsArg* pdfObs = pdfObsList[i];
488 if (std::string(obs.GetName())==pdfObs->GetName()) {
489 dhObs = dataHist.get()->find(histObs->GetName()) ;
490 break;
491 }
492 }
493
494 if (!dhObs) {
495 return nullptr;
496 }
497 RooAbsLValue* lval = dynamic_cast<RooAbsLValue*>(dhObs) ;
498 if (!lval) {
499 return nullptr;
500 }
501
502 // Retrieve position of all bin boundaries
503
504 const RooAbsBinning* binning = lval->getBinningPtr(nullptr);
505 std::span<const double> boundaries{binning->array(), static_cast<std::size_t>(binning->numBoundaries())};
506
507 // Use the helper function from RooCurve to make sure to get sampling hints
508 // that work with the RooFitPlotting.
509 return RooCurve::plotSamplingHintForBinBoundaries(boundaries, xlo, xhi);
510}
511
512
513////////////////////////////////////////////////////////////////////////////////
514/// Return sampling hint for making curves of (projections) of this function
515/// as the recursive division strategy of RooCurve cannot deal efficiently
516/// with the vertical lines that occur in a non-interpolated histogram
517
518std::list<double>* RooHistPdf::binBoundaries(RooAbsRealLValue& obs, double xlo, double xhi) const
519{
520 // No hints are required when interpolation is used
521 if (_intOrder>0) {
522 return nullptr;
523 }
524
525 // Check that observable is in dataset, if not no hint is generated
526 RooAbsLValue* lvarg = dynamic_cast<RooAbsLValue*>(_dataHist->get()->find(obs.GetName())) ;
527 if (!lvarg) {
528 return nullptr ;
529 }
530
531 // Retrieve position of all bin boundaries
532 const RooAbsBinning* binning = lvarg->getBinningPtr(nullptr);
533 double* boundaries = binning->array() ;
534
535 auto hint = new std::list<double> ;
536
537 // Construct array with pairs of points positioned epsilon to the left and
538 // right of the bin boundaries
539 for (Int_t i=0 ; i<binning->numBoundaries() ; i++) {
540 if (boundaries[i]>=xlo && boundaries[i]<=xhi) {
541 hint->push_back(boundaries[i]) ;
542 }
543 }
544
545 return hint ;
546}
547
548
549
550
551////////////////////////////////////////////////////////////////////////////////
552/// Only handle case of maximum in all variables
553
555{
556 std::unique_ptr<RooAbsCollection> common{_pdfObsList.selectCommon(vars)};
557 if (common->size()==_pdfObsList.size()) {
558 return 1;
559 }
560 return 0 ;
561}
562
563
564////////////////////////////////////////////////////////////////////////////////
565
566double RooHistPdf::maxVal(Int_t code) const
567{
568 R__ASSERT(code==1) ;
569
570 double max(-1) ;
571 for (Int_t i=0 ; i<_dataHist->numEntries() ; i++) {
572 _dataHist->get(i) ;
573 double wgt = _dataHist->weight() ;
574 if (wgt>max) max=wgt ;
575 }
576
577 return max*1.05 ;
578}
579
580
581
582
583////////////////////////////////////////////////////////////////////////////////
584
586{
587 if (std::abs(dh1.sumEntries()-dh2.sumEntries())>1e-8) return false ;
588 if (dh1.numEntries() != dh2.numEntries()) return false ;
589 for (int i=0 ; i < dh1.numEntries() ; i++) {
590 dh1.get(i) ;
591 dh2.get(i) ;
592 if (std::abs(dh1.weight()-dh2.weight())>1e-8) return false ;
593 }
594 return true ;
595}
596
597
598
599////////////////////////////////////////////////////////////////////////////////
600/// Check if our datahist is already in the workspace
601
603{
604 for(auto const& data : ws.allData()) {
605 // If your dataset is already in this workspace nothing needs to be done
606 if (data == _dataHist) {
607 return false ;
608 }
609 }
610
611 // Check if dataset with given name already exists
612 if (RooAbsData* wsdata = ws.embeddedData(_dataHist->GetName())) {
613
614 // Yes it exists - now check if it is identical to our internal histogram
615 if (wsdata->InheritsFrom(RooDataHist::Class())) {
616
617 // Check if histograms are identical
618 if (areIdentical(static_cast<RooDataHist&>(*wsdata),*_dataHist)) {
619
620 // Exists and is of correct type, and identical -- adjust internal pointer to WS copy
621 _dataHist = static_cast<RooDataHist*>(wsdata) ;
622 } else {
623
624 // not identical, clone rename and import
625 auto uniqueName = std::string(_dataHist->GetName()) + "_" + GetName();
626 bool flag = ws.import(*_dataHist,RooFit::Rename(uniqueName.c_str()),RooFit::Embedded()) ;
627 if (flag) {
628 coutE(ObjectHandling) << " RooHistPdf::importWorkspaceHook(" << GetName() << ") unable to import clone of underlying RooDataHist with unique name " << uniqueName << ", abort" << std::endl ;
629 return true ;
630 }
631 _dataHist = static_cast<RooDataHist*>(ws.embeddedData(uniqueName)) ;
632 }
633
634 } else {
635
636 // Exists and is NOT of correct type: clone rename and import
637 auto uniqueName = std::string(_dataHist->GetName()) + "_" + GetName();
638 bool flag = ws.import(*_dataHist,RooFit::Rename(uniqueName.c_str()),RooFit::Embedded()) ;
639 if (flag) {
640 coutE(ObjectHandling) << " RooHistPdf::importWorkspaceHook(" << GetName() << ") unable to import clone of underlying RooDataHist with unique name " << uniqueName << ", abort" << std::endl ;
641 return true ;
642 }
643 _dataHist = static_cast<RooDataHist*>(ws.embeddedData(uniqueName));
644
645 }
646 return false ;
647 }
648
649 // We need to import our datahist into the workspace
651
652 // Redirect our internal pointer to the copy in the workspace
653 _dataHist = static_cast<RooDataHist*>(ws.embeddedData(_dataHist->GetName())) ;
654 return false ;
655}
656
657
658////////////////////////////////////////////////////////////////////////////////
659/// Stream an object of class RooHistPdf.
660
662{
663 if (R__b.IsReading()) {
665 // WVE - interim solution - fix proxies here
666 //_proxyList.Clear() ;
667 //registerProxy(_pdfObsList) ;
668 } else {
670 }
671}
#define e(i)
Definition RSha256.hxx:103
#define oocoutE(o, a)
#define coutE(a)
#define ooccoutE(o, a)
#define ClassImp(name)
Definition Rtypes.h:377
#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
Common abstract base class for objects that represent a value and a "shape" in RooFit.
Definition RooAbsArg.h:77
virtual void copyCache(const RooAbsArg *source, bool valueOnly=false, bool setValDirty=true)=0
bool dependsOn(const RooAbsCollection &serverList, const RooAbsArg *ignoreArg=nullptr, bool valueOnly=false) const
Test whether we depend on (ie, are served by) any object in the specified collection.
virtual void syncCache(const RooArgSet *nset=nullptr)=0
friend void RooRefArray::Streamer(TBuffer &)
virtual bool inRange(const char *) const
Definition RooAbsArg.h:376
Abstract base class for RooRealVar binning definitions.
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.
Abstract base class for binned and unbinned datasets.
Definition RooAbsData.h:57
virtual Int_t numEntries() const
Return number of entries in dataset, i.e., count unweighted entries.
Abstract base class for objects that are lvalues, i.e.
virtual const RooAbsBinning * getBinningPtr(const char *rangeName) const =0
Abstract interface for all probability density functions.
Definition RooAbsPdf.h:40
Abstract base class for objects that represent a real value that may appear on the left hand side of ...
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 doEval(RooFit::EvalContext &) 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:55
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...
static std::list< double > * plotSamplingHintForBinBoundaries(std::span< const double > boundaries, double xlo, double xhi)
Returns sampling hints for a histogram with given boundaries.
Definition RooCurve.cxx:921
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.
void weights(double *output, std::span< double const > xVals, int intOrder, bool correctForBinSize, bool cdfBoundaries)
A vectorized version of RooDataHist::weight() for one dimensional histograms with up to one dimension...
static TClass * Class()
TObject * Clone(const char *newname="") const override
Make a clone of an object using the Streamer facility.
Definition RooDataHist.h:55
std::string calculateTreeIndexForCodeSquash(RooAbsArg const *klass, RooFit::Detail::CodeSquashContext &ctx, const RooAbsCollection &coords, bool reverse=false) const
double weight(std::size_t i) const
Return weight of i-th bin.
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...
std::string declWeightArrayForCodeSquash(RooAbsArg const *klass, RooFit::Detail::CodeSquashContext &ctx, bool correctForBinSize) const
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.
A class to maintain the context for squashing of RooFit models into code.
void addResult(RooAbsArg const *key, std::string const &value)
A function to save an expression that includes/depends on the result of the input node.
std::span< const double > at(RooAbsArg const *arg, RooAbsArg const *caller=nullptr)
std::span< double > output()
A propability density function sampled from a multidimensional histogram.
Definition RooHistPdf.h:30
static void rooHistTranslateImpl(RooAbsArg const *klass, RooFit::Detail::CodeSquashContext &ctx, int intOrder, RooDataHist const *dataHist, const RooArgSet &obs, bool correctForBinSize)
RooArgSet _histObsList
List of observables defining dimensions of histogram.
Definition RooHistPdf.h:112
Int_t _intOrder
Interpolation order.
Definition RooHistPdf.h:117
bool forceAnalyticalInt(const RooAbsArg &dep) const override
bool areIdentical(const RooDataHist &dh1, const RooDataHist &dh2)
RooDataHist * _dataHist
Unowned pointer to underlying histogram.
Definition RooHistPdf.h:114
bool _cdfBoundaries
Use boundary conditions for CDFs.
Definition RooHistPdf.h:118
static std::string rooHistIntegralTranslateImpl(int code, RooAbsArg const *klass, RooDataHist const *dataHist, const RooArgSet &obs, bool histFuncMode)
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...
double totVolume() const
Return the total volume spanned by the observables of the RooHistPdf.
void initializeOwnedDataHist(std::unique_ptr< RooDataHist > &&dataHist)
Definition RooHistPdf.h:157
bool importWorkspaceHook(RooWorkspace &ws) override
Check if our datahist is already in the workspace.
static TClass * Class()
std::string buildCallToAnalyticIntegral(int code, const char *rangeName, RooFit::Detail::CodeSquashContext &ctx) const override
This function defines the analytical integral translation for the class.
RooSetProxy _pdfObsList
List of observables mapped onto histogram observables.
Definition RooHistPdf.h:113
double maxVal(Int_t code) const override
Return maximum value for set of observables identified by code assigned in getMaxVal.
double analyticalIntegral(Int_t code, const char *rangeName=nullptr) const override
Return integral identified by 'code'.
void translate(RooFit::Detail::CodeSquashContext &ctx) const override
This function defines a translation for each RooAbsReal based object that can be used to express the ...
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...
RooDataHist & dataHist()
Definition RooHistPdf.h:42
Int_t getMaxVal(const RooArgSet &vars) const override
Only handle case of maximum in all variables.
double _totVolume
! Total volume of space (product of ranges of observables)
Definition RooHistPdf.h:119
RooDataHist * cloneAndOwnDataHist(const char *newname="")
Replaces underlying RooDataHist with a clone, which is now owned, and returns the clone.
std::unique_ptr< RooDataHist > _ownedDataHist
! Owned pointer to underlying histogram
Definition RooHistPdf.h:115
void doEval(RooFit::EvalContext &) const override
Base function for computing multiple values of a RooAbsReal.
Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=nullptr) const override
Determine integration scenario.
double evaluate() const override
Return the current value: The value of the bin enclosing the current coordinates of the observables,...
bool _unitNorm
Assume contents is unit normalized (for use as pdf cache)
Definition RooHistPdf.h:120
Variable that can be changed from the outside.
Definition RooRealVar.h:37
void setRange(const char *name, double min, double max)
Set a fit or plotting range.
Persistable container for RooFit projects.
RooAbsData * embeddedData(RooStringView name) const
Retrieve dataset (binned or unbinned) with given name. A null pointer is returned if not found.
std::list< RooAbsData * > allData() const
Return list of all dataset in the workspace.
bool import(const RooAbsArg &arg, const RooCmdArg &arg1={}, const RooCmdArg &arg2={}, const RooCmdArg &arg3={}, const RooCmdArg &arg4={}, const RooCmdArg &arg5={}, const RooCmdArg &arg6={}, const RooCmdArg &arg7={}, const RooCmdArg &arg8={}, const RooCmdArg &arg9={})
Import a RooAbsArg object, e.g.
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
std::pair< double, double > getRangeOrBinningInterval(RooAbsArg const *arg, const char *rangeName)
static void output()