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 bool cdfBoundaries)
228{
229 if (intOrder != 0 && !(!cdfBoundaries && !correctForBinSize && intOrder == 1 && obs.size() == 1)) {
230 ooccoutE(klass, InputArguments) << "RooHistPdf::weight(" << klass->GetName()
231 << ") ERROR: Code Squashing currently only supports non-interpolation cases."
232 << std::endl;
233 return;
234 }
235
236 if (intOrder == 1) {
237 RooAbsBinning const &binning = *dataHist->getBinnings()[0];
238 std::string weightArr = dataHist->declWeightArrayForCodeSquash(ctx, correctForBinSize);
239 ctx.addResult(klass, ctx.buildCall("RooFit::Detail::MathFuncs::interpolate1d", binning.lowBound(),
240 binning.highBound(), *obs[0], binning.numBins(), weightArr));
241 return;
242 }
243 std::string const &offset = dataHist->calculateTreeIndexForCodeSquash(klass, ctx, obs);
244 std::string weightArr = dataHist->declWeightArrayForCodeSquash(ctx, correctForBinSize);
245 ctx.addResult(klass, "*(" + weightArr + " + " + offset + ")");
246}
247
249{
251}
252
253////////////////////////////////////////////////////////////////////////////////
254/// Return the total volume spanned by the observables of the RooHistPdf
255
257{
258 // Return previously calculated value, if any
259 if (_totVolume>0) {
260 return _totVolume ;
261 }
262 _totVolume = 1. ;
263
264 for (const auto arg : _histObsList) {
265 RooRealVar* real = dynamic_cast<RooRealVar*>(arg) ;
266 if (real) {
267 _totVolume *= (real->getMax()-real->getMin()) ;
268 } else {
269 RooCategory* cat = dynamic_cast<RooCategory*>(arg) ;
270 if (cat) {
271 _totVolume *= cat->numTypes() ;
272 }
273 }
274 }
275
276 return _totVolume ;
277}
278
279namespace {
280
281bool fullRange(const RooAbsArg& x, const RooAbsArg& y ,const char* range)
282{
283 const RooAbsRealLValue *_x = dynamic_cast<const RooAbsRealLValue*>(&x);
284 const RooAbsRealLValue *_y = dynamic_cast<const RooAbsRealLValue*>(&y);
285 if (!_x || !_y) return false;
286 if (!range || !strlen(range) || !_x->hasRange(range) ||
287 _x->getBinningPtr(range)->isParameterized()) {
288 // parameterized ranges may be full range now, but that might change,
289 // so return false
290 if (range && strlen(range) && _x->getBinningPtr(range)->isParameterized())
291 return false;
292 return (_x->getMin() == _y->getMin() && _x->getMax() == _y->getMax());
293 }
294 return (_x->getMin(range) == _y->getMin() && _x->getMax(range) == _y->getMax());
295}
296
297bool okayForAnalytical(RooAbsArg const& obs, RooArgSet const& allVars)
298{
299 auto lobs = dynamic_cast<RooAbsRealLValue const*>(&obs);
300 if(lobs == nullptr) return false;
301
302 bool isOkayForAnalyticalInt = false;
303
304 for(RooAbsArg *var : allVars) {
305 if(obs.dependsOn(*var)) {
306 if(!lobs->isJacobianOK(*var)) return false;
307 isOkayForAnalyticalInt = true;
308 }
309 }
310
311 return isOkayForAnalyticalInt;
312}
313
314} // namespace
315
316
318 RooArgSet& analVars,
319 const char* rangeName,
320 RooArgSet const& histObsList,
321 RooArgSet const& pdfObsList,
322 Int_t intOrder)
323{
324 // First make list of pdf observables to histogram observables
325 // and select only those for which the integral is over the full range
326
327 Int_t code = 0;
328 Int_t frcode = 0;
329 for (unsigned int n=0; n < pdfObsList.size() && n < histObsList.size(); ++n) {
330 const auto pa = pdfObsList[n];
331 const auto ha = histObsList[n];
332
333 if (okayForAnalytical(*pa, allVars)) {
334 code |= 2 << n;
335 analVars.add(*pa);
336 if (fullRange(*pa, *ha, rangeName)) {
337 frcode |= 2 << n;
338 }
339 }
340 }
341
342 if (code == frcode) {
343 // integrate over full range of all observables - use bit 0 to indicate
344 // full range integration over all observables
345 code |= 1;
346 }
347
348 // Disable partial analytical integrals if interpolation is used, and we
349 // integrate over sub-ranges, but leave them enabled when we integrate over
350 // the full range of one or several variables
351 if (intOrder > 1 && !(code & 1)) {
352 analVars.removeAll();
353 return 0;
354 }
355 return (code >= 2) ? code : 0;
356}
357
358
360 const char* rangeName,
361 RooArgSet const& histObsList,
362 RooArgSet const& pdfObsList,
363 RooDataHist& dataHist,
364 bool histFuncMode) {
365 // Simplest scenario, full-range integration over all dependents
366 if (((2 << histObsList.size()) - 1) == code) {
367 return dataHist.sum(histFuncMode);
368 }
369
370 // Partial integration scenario, retrieve set of variables, calculate partial
371 // sum, figure out integration ranges (if needed)
372 RooArgSet intSet;
373 std::map<const RooAbsArg*, std::pair<double, double> > ranges;
374 for (unsigned int n=0; n < pdfObsList.size() && n < histObsList.size(); ++n) {
375 const auto pa = pdfObsList[n];
376 const auto ha = histObsList[n];
377
378 if (code & (2 << n)) {
379 intSet.add(*ha);
380 }
381 if (!(code & 1)) {
382 ranges[ha] = RooHelpers::getRangeOrBinningInterval(pa, rangeName);
383 }
384 // WVE must sync hist slice list values to pdf slice list
385 // Transfer values from
386 if (ha != pa) {
387 pa->syncCache();
388 ha->copyCache(pa,true);
389 }
390 }
391
392 double ret = (code & 1) ? dataHist.sum(intSet,histObsList,true,!histFuncMode) :
393 dataHist.sum(intSet,histObsList,true,!histFuncMode, ranges);
394
395 return ret ;
396}
397
398std::string RooHistPdf::rooHistIntegralTranslateImpl(int code, RooAbsArg const *klass, RooDataHist const *dataHist,
399 const RooArgSet &obs, bool histFuncMode)
400{
401 if (((2 << obs.size()) - 1) != code) {
402 oocoutE(klass, InputArguments)
403 << "RooHistPdf::integral(" << klass->GetName()
404 << ") ERROR: AD currently only supports integrating over all histogram observables." << std::endl;
405 return "";
406 }
407 return std::to_string(dataHist->sum(histFuncMode));
408}
409
410std::string RooHistPdf::buildCallToAnalyticIntegral(int code, const char * /*rangeName */,
411 RooFit::Detail::CodeSquashContext & /* ctx */) const
412{
413 return rooHistIntegralTranslateImpl(code, this, _dataHist, _pdfObsList, false);
414}
415
416////////////////////////////////////////////////////////////////////////////////
417/// Determine integration scenario. If no interpolation is used,
418/// RooHistPdf can perform all integrals over its dependents
419/// analytically via partial or complete summation of the input
420/// histogram. If interpolation is used on the integral over
421/// all histogram observables is supported
422
423Int_t RooHistPdf::getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* rangeName) const
424{
425 return getAnalyticalIntegral(allVars, analVars, rangeName, _histObsList, _pdfObsList, _intOrder);
426}
427
428
429////////////////////////////////////////////////////////////////////////////////
430/// Return integral identified by 'code'. The actual integration
431/// is deferred to RooDataHist::sum() which implements partial
432/// or complete summation over the histograms contents.
433
434double RooHistPdf::analyticalIntegral(Int_t code, const char* rangeName) const
435{
436 return analyticalIntegral(code, rangeName, _histObsList, _pdfObsList, *_dataHist, false);
437}
438
439
440bool RooHistPdf::forceAnalyticalInt(RooArgSet const& pdfObsList, const RooAbsArg& dep)
441{
442 bool isOkayForAnalyticalInt = false;
443
444 for (RooAbsArg * obs : pdfObsList) {
445 if(obs->dependsOn(dep)) {
446 // If the observable doesn't depend linearly on the integration
447 // variable we will not do analytical integration.
448 auto lvalue = dynamic_cast<RooAbsRealLValue const*>(obs);
449 if(!(lvalue && lvalue->isJacobianOK(dep))) return false;
450 isOkayForAnalyticalInt = true;
451 }
452 }
453
454 return isOkayForAnalyticalInt;
455}
456
457
459{
460 return forceAnalyticalInt(_pdfObsList, dep);
461}
462
463
464////////////////////////////////////////////////////////////////////////////////
465/// Return sampling hint for making curves of (projections) of this function
466/// as the recursive division strategy of RooCurve cannot deal efficiently
467/// with the vertical lines that occur in a non-interpolated histogram
468
469std::list<double>* RooHistPdf::plotSamplingHint(RooAbsRealLValue& obs, double xlo, double xhi) const
470{
472}
473
474
475std::list<double>* RooHistPdf::plotSamplingHint(RooDataHist const& dataHist,
476 RooArgSet const& pdfObsList,
477 RooArgSet const& histObsList,
478 int intOrder,
479 RooAbsRealLValue& obs,
480 double xlo,
481 double xhi)
482{
483 // No hints are required when interpolation is used
484 if (intOrder>0) {
485 return nullptr;
486 }
487
488 // Check that observable is in dataset, if not no hint is generated
489 RooAbsArg* dhObs = nullptr;
490 for (unsigned int i=0; i < pdfObsList.size(); ++i) {
491 RooAbsArg* histObs = histObsList[i];
492 RooAbsArg* pdfObs = pdfObsList[i];
493 if (std::string(obs.GetName())==pdfObs->GetName()) {
494 dhObs = dataHist.get()->find(histObs->GetName()) ;
495 break;
496 }
497 }
498
499 if (!dhObs) {
500 return nullptr;
501 }
502 RooAbsLValue* lval = dynamic_cast<RooAbsLValue*>(dhObs) ;
503 if (!lval) {
504 return nullptr;
505 }
506
507 // Retrieve position of all bin boundaries
508
509 const RooAbsBinning* binning = lval->getBinningPtr(nullptr);
510 std::span<const double> boundaries{binning->array(), static_cast<std::size_t>(binning->numBoundaries())};
511
512 // Use the helper function from RooCurve to make sure to get sampling hints
513 // that work with the RooFitPlotting.
514 return RooCurve::plotSamplingHintForBinBoundaries(boundaries, xlo, xhi);
515}
516
517
518////////////////////////////////////////////////////////////////////////////////
519/// Return sampling hint for making curves of (projections) of this function
520/// as the recursive division strategy of RooCurve cannot deal efficiently
521/// with the vertical lines that occur in a non-interpolated histogram
522
523std::list<double>* RooHistPdf::binBoundaries(RooAbsRealLValue& obs, double xlo, double xhi) const
524{
525 // No hints are required when interpolation is used
526 if (_intOrder>0) {
527 return nullptr;
528 }
529
530 // Check that observable is in dataset, if not no hint is generated
531 RooAbsLValue* lvarg = dynamic_cast<RooAbsLValue*>(_dataHist->get()->find(obs.GetName())) ;
532 if (!lvarg) {
533 return nullptr ;
534 }
535
536 // Retrieve position of all bin boundaries
537 const RooAbsBinning* binning = lvarg->getBinningPtr(nullptr);
538 double* boundaries = binning->array() ;
539
540 auto hint = new std::list<double> ;
541
542 // Construct array with pairs of points positioned epsilon to the left and
543 // right of the bin boundaries
544 for (Int_t i=0 ; i<binning->numBoundaries() ; i++) {
545 if (boundaries[i]>=xlo && boundaries[i]<=xhi) {
546 hint->push_back(boundaries[i]) ;
547 }
548 }
549
550 return hint ;
551}
552
553
554
555
556////////////////////////////////////////////////////////////////////////////////
557/// Only handle case of maximum in all variables
558
560{
561 std::unique_ptr<RooAbsCollection> common{_pdfObsList.selectCommon(vars)};
562 if (common->size()==_pdfObsList.size()) {
563 return 1;
564 }
565 return 0 ;
566}
567
568
569////////////////////////////////////////////////////////////////////////////////
570
571double RooHistPdf::maxVal(Int_t code) const
572{
573 R__ASSERT(code==1) ;
574
575 double max(-1) ;
576 for (Int_t i=0 ; i<_dataHist->numEntries() ; i++) {
577 _dataHist->get(i) ;
578 double wgt = _dataHist->weight() ;
579 if (wgt>max) max=wgt ;
580 }
581
582 return max*1.05 ;
583}
584
585
586
587
588////////////////////////////////////////////////////////////////////////////////
589
591{
592 if (std::abs(dh1.sumEntries()-dh2.sumEntries())>1e-8) return false ;
593 if (dh1.numEntries() != dh2.numEntries()) return false ;
594 for (int i=0 ; i < dh1.numEntries() ; i++) {
595 dh1.get(i) ;
596 dh2.get(i) ;
597 if (std::abs(dh1.weight()-dh2.weight())>1e-8) return false ;
598 }
599 return true ;
600}
601
602
603
604////////////////////////////////////////////////////////////////////////////////
605/// Check if our datahist is already in the workspace
606
608{
609 for(auto const& data : ws.allData()) {
610 // If your dataset is already in this workspace nothing needs to be done
611 if (data == _dataHist) {
612 return false ;
613 }
614 }
615
616 // Check if dataset with given name already exists
617 if (RooAbsData* wsdata = ws.embeddedData(_dataHist->GetName())) {
618
619 // Yes it exists - now check if it is identical to our internal histogram
620 if (wsdata->InheritsFrom(RooDataHist::Class())) {
621
622 // Check if histograms are identical
623 if (areIdentical(static_cast<RooDataHist&>(*wsdata),*_dataHist)) {
624
625 // Exists and is of correct type, and identical -- adjust internal pointer to WS copy
626 _dataHist = static_cast<RooDataHist*>(wsdata) ;
627 } else {
628
629 // not identical, clone rename and import
630 auto uniqueName = std::string(_dataHist->GetName()) + "_" + GetName();
631 bool flag = ws.import(*_dataHist,RooFit::Rename(uniqueName.c_str()),RooFit::Embedded()) ;
632 if (flag) {
633 coutE(ObjectHandling) << " RooHistPdf::importWorkspaceHook(" << GetName() << ") unable to import clone of underlying RooDataHist with unique name " << uniqueName << ", abort" << std::endl ;
634 return true ;
635 }
636 _dataHist = static_cast<RooDataHist*>(ws.embeddedData(uniqueName)) ;
637 }
638
639 } else {
640
641 // Exists and is NOT of correct type: clone rename and import
642 auto uniqueName = std::string(_dataHist->GetName()) + "_" + GetName();
643 bool flag = ws.import(*_dataHist,RooFit::Rename(uniqueName.c_str()),RooFit::Embedded()) ;
644 if (flag) {
645 coutE(ObjectHandling) << " RooHistPdf::importWorkspaceHook(" << GetName() << ") unable to import clone of underlying RooDataHist with unique name " << uniqueName << ", abort" << std::endl ;
646 return true ;
647 }
648 _dataHist = static_cast<RooDataHist*>(ws.embeddedData(uniqueName));
649
650 }
651 return false ;
652 }
653
654 // We need to import our datahist into the workspace
656
657 // Redirect our internal pointer to the copy in the workspace
658 _dataHist = static_cast<RooDataHist*>(ws.embeddedData(_dataHist->GetName())) ;
659 return false ;
660}
661
662
663////////////////////////////////////////////////////////////////////////////////
664/// Stream an object of class RooHistPdf.
665
667{
668 if (R__b.IsReading()) {
670 // WVE - interim solution - fix proxies here
671 //_proxyList.Clear() ;
672 //registerProxy(_pdfObsList) ;
673 } else {
675 }
676}
#define e(i)
Definition RSha256.hxx:103
#define oocoutE(o, a)
#define coutE(a)
#define ooccoutE(o, a)
#define ClassImp(name)
Definition Rtypes.h:382
#define R__ASSERT(e)
Checks condition e and reports a fatal error if it's false.
Definition TError.h:125
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:79
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.
R__DEPRECATED(6, 36, "Use getObservables().") RooFit R__DEPRECATED(6, 36, "Use getObservables().") RooFit R__DEPRECATED(6, 36, "Use getObservables().") RooFit const RooAbsArg &testArg const
Definition RooAbsArg.h:145
virtual void syncCache(const RooArgSet *nset=nullptr)=0
friend void RooRefArray::Streamer(TBuffer &)
virtual bool inRange(const char *) const
Definition RooAbsArg.h:350
Abstract base class for RooRealVar binning definitions.
Int_t numBins() const
Return number of bins.
virtual bool isParameterized() const
Interface function.
virtual double highBound() const =0
virtual double lowBound() const =0
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.
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.
virtual double offset() const
Definition RooAbsReal.h:371
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:24
RooArgSet * selectCommon(const RooAbsCollection &refColl) const
Use RooAbsCollection::selecCommon(), but return as RooArgSet.
Definition RooArgSet.h:149
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:893
Container class to hold N-dimensional binned data.
Definition RooDataHist.h:40
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...
std::vector< std::unique_ptr< const RooAbsBinning > > const & getBinnings() const
static TClass * Class()
TObject * Clone(const char *newname="") const override
Make a clone of an object using the Streamer facility.
Definition RooDataHist.h:61
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(RooFit::Detail::CodeSquashContext &ctx, bool correctForBinSize) const
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.
A class to maintain the context for squashing of RooFit models into code.
std::string buildCall(std::string const &funcname, Args_t const &...args)
Build the code to call the function with name funcname, passing some arguments.
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
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.
static void rooHistTranslateImpl(RooAbsArg const *klass, RooFit::Detail::CodeSquashContext &ctx, int intOrder, RooDataHist const *dataHist, const RooArgSet &obs, bool correctForBinSize, bool cdfBoundaries)
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()