Logo ROOT  
Reference Guide
RooRealSumPdf.cxx
Go to the documentation of this file.
1/*****************************************************************************
2 * Project: RooFit *
3 * Package: RooFitCore *
4 * @(#)root/roofitcore:$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/** \class RooRealSumPdf
19 \ingroup Roofitcore
20
21
22The class RooRealSumPdf implements a PDF constructed from a sum of functions:
23\f[
24 \mathrm{PDF}(x) = \frac{ \sum_{i=1}^{n-1} \mathrm{coef}_i * \mathrm{func}_i(x) + \left[ 1 - \sum_{i=1}^{n-1} \mathrm{coef}_i \right] * \mathrm{func}_n(x) }
25 {\sum_{i=1}^{n-1} \mathrm{coef}_i * \int \mathrm{func}_i(x)dx + \left[ 1 - \sum_{i=1}^{n-1} \mathrm{coef}_i \right] * \int \mathrm{func}_n(x) dx }
26\f]
27
28where \f$\mathrm{coef}_i\f$ and \f$\mathrm{func}_i\f$ are RooAbsReal objects, and \f$ x \f$ is the collection of dependents.
29In the present version \f$\mathrm{coef}_i\f$ may not depend on \f$ x \f$, but this limitation could be removed should the need arise.
30
31If the number of coefficients is one less than the number of functions, the PDF is assumed to be normalised. Due to this additional constraint,
32\f$\mathrm{coef}_n\f$ is computed from the other coefficients.
33
34### Extending the PDF
35If an \f$ n^\mathrm{th} \f$ coefficient is provided, the PDF **can** be used as an extended PDF, *i.e.* the total number of events will be measured in addition
36to the fractions of the various functions. **This requires setting the last argument of the constructor to `true`.**
37\note For the RooAddPdf, the extension happens automatically.
38
39### Difference to RooAddPdf / RooRealSumFunc
40- RooAddPdf is a PDF of PDFs, *i.e.* its components need to be normalised and non-negative.
41- RooRealSumPdf is a PDF of functions, *i.e.*, its components can be negative, but their sum cannot be. The normalisation
42 is computed automatically, unless the PDF is extended (see above).
43- RooRealSumFunc is a sum of functions. It is neither normalised, nor need it be positive.
44
45*/
46
47#include "RooRealSumPdf.h"
48
49#include "RooNormalizedPdf.h"
50#include "RooRealIntegral.h"
51#include "RooRealProxy.h"
52#include "RooRealVar.h"
53#include "RooMsgService.h"
54#include "RooNaNPacker.h"
55#include "RunContext.h"
56
57#include <TError.h>
58
59#include <algorithm>
60#include <memory>
61#include <stdexcept>
62
63using namespace std;
64
66
68
69////////////////////////////////////////////////////////////////////////////////
70/// Default constructor
71/// coverity[UNINIT_CTOR]
72
73RooRealSumPdf::RooRealSumPdf() : _normIntMgr(this,10)
74{
75
76}
77
78
79
80////////////////////////////////////////////////////////////////////////////////
81/// Constructor with name and title
82
83RooRealSumPdf::RooRealSumPdf(const char *name, const char *title) :
84 RooAbsPdf(name,title),
85 _normIntMgr(this,10),
86 _funcList("!funcList","List of functions",this),
87 _coefList("!coefList","List of coefficients",this),
88 _extended(false),
89 _doFloor(false)
90{
91
92}
93
94
95
96////////////////////////////////////////////////////////////////////////////////
97/// Construct p.d.f consisting of \f$ \mathrm{coef}_1 * \mathrm{func}_1 + (1-\mathrm{coef}_1) * \mathrm{func}_2 \f$.
98/// The input coefficients and functions are allowed to be negative
99/// but the resulting sum is not, which is enforced at runtime.
100
101RooRealSumPdf::RooRealSumPdf(const char *name, const char *title,
102 RooAbsReal& func1, RooAbsReal& func2, RooAbsReal& coef1) :
103 RooRealSumPdf(name, title)
104{
105 // Special constructor with two functions and one coefficient
106
107 _funcList.add(func1) ;
108 _funcList.add(func2) ;
109 _coefList.add(coef1) ;
110
111}
112
113
114////////////////////////////////////////////////////////////////////////////////
115/// Constructor for a PDF from a list of functions and coefficients.
116/// It implements
117/// \f[
118/// \sum_i \mathrm{coef}_i \cdot \mathrm{func}_i,
119/// \f]
120/// if \f$ N_\mathrm{coef} = N_\mathrm{func} \f$. With `extended=true`, the coefficients can take any values. With `extended=false`,
121/// there is the danger of getting a degenerate minimisation problem because a PDF has to be normalised, which needs one degree
122/// of freedom less.
123///
124/// A plain (normalised) PDF can therefore be implemented with one less coefficient. RooFit then computes
125/// \f[
126/// \sum_i^{N-1} \mathrm{coef}_i \cdot \mathrm{func}_i + (1 - \sum_i \mathrm{coef}_i ) \cdot \mathrm{func}_N,
127/// \f]
128/// if \f$ N_\mathrm{coef} = N_\mathrm{func} - 1 \f$.
129///
130/// All coefficients and functions are allowed to be negative
131/// but the sum (*i.e.* the PDF) is not, which is enforced at runtime.
132///
133/// \param name Name of the PDF
134/// \param title Title (for plotting)
135/// \param inFuncList List of functions to sum
136/// \param inCoefList List of coefficients
137/// \param extended Interpret as extended PDF (requires equal number of functions and coefficients)
138
139RooRealSumPdf::RooRealSumPdf(const char *name, const char *title,
140 const RooArgList& inFuncList, const RooArgList& inCoefList, bool extended) :
141 RooRealSumPdf(name, title)
142{
143 _extended = extended;
144 RooRealSumPdf::initializeFuncsAndCoefs(*this, inFuncList, inCoefList, _funcList, _coefList);
145}
146
147
149 const RooArgList& inFuncList, const RooArgList& inCoefList,
150 RooArgList& funcList, RooArgList& coefList)
151{
152 const std::string className = caller.ClassName();
153 const std::string constructorName = className + "::" + className;
154
155 if (!(inFuncList.size()==inCoefList.size()+1 || inFuncList.size()==inCoefList.size())) {
156 oocoutE(&caller, InputArguments) << constructorName << "(" << caller.GetName()
157 << ") number of pdfs and coefficients inconsistent, must have Nfunc=Ncoef or Nfunc=Ncoef+1" << std::endl;
158 throw std::invalid_argument(className + ": Number of PDFs and coefficients is inconsistent.");
159 }
160
161 // Constructor with N functions and N or N-1 coefs
162 for (unsigned int i = 0; i < inCoefList.size(); ++i) {
163 const auto& func = inFuncList[i];
164 const auto& coef = inCoefList[i];
165
166 if (!dynamic_cast<const RooAbsReal*>(&coef)) {
167 oocoutW(&caller, InputArguments) << constructorName << "(" << caller.GetName() << ") coefficient " << coef.GetName() << " is not of type RooAbsReal, ignored" << std::endl ;
168 continue ;
169 }
170 if (!dynamic_cast<const RooAbsReal*>(&func)) {
171 oocoutW(&caller, InputArguments) << constructorName << "(" << caller.GetName() << ") func " << func.GetName() << " is not of type RooAbsReal, ignored" << std::endl ;
172 continue ;
173 }
174 funcList.add(func) ;
175 coefList.add(coef) ;
176 }
177
178 if (inFuncList.size() == inCoefList.size() + 1) {
179 const auto& func = inFuncList[inFuncList.size()-1];
180 if (!dynamic_cast<const RooAbsReal*>(&func)) {
181 oocoutE(&caller, InputArguments) << constructorName << "(" << caller.GetName() << ") last func " << func.GetName() << " is not of type RooAbsReal, fatal error" << std::endl ;
182 throw std::invalid_argument(className + ": Function passed as is not of type RooAbsReal.");
183 }
184 funcList.add(func);
185 }
186
187}
188
189
190
191
192////////////////////////////////////////////////////////////////////////////////
193/// Copy constructor
194
196 RooAbsPdf(other,name),
197 _normIntMgr(other._normIntMgr,this),
198 _funcList("!funcList",this,other._funcList),
199 _coefList("!coefList",this,other._coefList),
200 _extended(other._extended),
201 _doFloor(other._doFloor)
202{
203
204}
205
206
207
208////////////////////////////////////////////////////////////////////////////////
209/// Destructor
210
212{
213
214}
215
216
217
218
219
220////////////////////////////////////////////////////////////////////////////////
221
223{
225}
226
227
229 RooArgList const& funcList,
230 RooArgList const& coefList,
231 bool doFloor,
232 bool & hasWarnedBefore)
233{
234 // Do running sum of coef/func pairs, calculate lastCoef.
235 double value = 0;
236 double sumCoeff = 0.;
237 for (unsigned int i = 0; i < funcList.size(); ++i) {
238 const auto func = static_cast<RooAbsReal*>(&funcList[i]);
239 const auto coef = static_cast<RooAbsReal*>(i < coefList.size() ? &coefList[i] : nullptr);
240 const double coefVal = coef != nullptr ? coef->getVal() : (1. - sumCoeff);
241
242 // Warn about degeneration of last coefficient
243 if (coef == nullptr && (coefVal < 0 || coefVal > 1.)) {
244 if (!hasWarnedBefore) {
245 oocoutW(&caller, Eval) << caller.ClassName() << "::evaluate(" << caller.GetName()
246 << ") WARNING: sum of FUNC coefficients not in range [0-1], value="
247 << sumCoeff << ". This means that the PDF is not properly normalised."
248 << " If the PDF was meant to be extended, provide as many coefficients as functions." << std::endl;
249 hasWarnedBefore = true;
250 }
251 // Signal that we are in an undefined region:
252 value = RooNaNPacker::packFloatIntoNaN(100.f * (coefVal < 0. ? -coefVal : coefVal - 1.));
253 }
254
255 if (func->isSelectedComp()) {
256 value += func->getVal() * coefVal;
257 }
258
259 sumCoeff += coefVal;
260 }
261
262 // Introduce floor if so requested
263 return value < 0 && doFloor ? 0.0 : value;
264}
265
266
267////////////////////////////////////////////////////////////////////////////////
268/// Calculate the current value
269
271{
273}
274
275
276void RooRealSumPdf::computeBatch(cudaStream_t* /*stream*/, double* output, size_t nEvents, RooFit::Detail::DataMap const& dataMap) const {
277 // Do running sum of coef/func pairs, calculate lastCoef.
278 for (unsigned int j = 0; j < nEvents; ++j) {
279 output[j] = 0.0;
280 }
281
282 double sumCoeff = 0.;
283 for (unsigned int i = 0; i < _funcList.size(); ++i) {
284 const auto func = static_cast<RooAbsReal*>(&_funcList[i]);
285 const auto coef = static_cast<RooAbsReal*>(i < _coefList.size() ? &_coefList[i] : nullptr);
286 const double coefVal = coef != nullptr ? dataMap.at(coef)[0] : (1. - sumCoeff);
287
288 if (func->isSelectedComp()) {
289 auto funcValues = dataMap.at(func);
290 if(funcValues.size() == 1) {
291 for (unsigned int j = 0; j < nEvents; ++j) {
292 output[j] += funcValues[0] * coefVal;
293 }
294 } else {
295 for (unsigned int j = 0; j < nEvents; ++j) {
296 output[j] += funcValues[j] * coefVal;
297 }
298 }
299 }
300
301 // Warn about degeneration of last coefficient
302 if (coef == nullptr && (coefVal < 0 || coefVal > 1.)) {
303 if (!_haveWarned) {
304 coutW(Eval) << "RooRealSumPdf::evaluateSpan(" << GetName()
305 << ") WARNING: sum of FUNC coefficients not in range [0-1], value="
306 << sumCoeff << ". This means that the PDF is not properly normalised. If the PDF was meant to be extended, provide as many coefficients as functions." << endl ;
307 _haveWarned = true;
308 }
309 // Signal that we are in an undefined region by handing back one NaN.
310 output[0] = RooNaNPacker::packFloatIntoNaN(100.f * (coefVal < 0. ? -coefVal : coefVal - 1.));
311 }
312
313 sumCoeff += coefVal;
314 }
315
316 // Introduce floor if so requested
317 if (_doFloor || _doFloorGlobal) {
318 for (unsigned int j = 0; j < nEvents; ++j) {
319 output[j] += std::max(0., output[j]);
320 }
321 }
322}
323
324
326 RooArgList const& funcList, RooArgList const& coefList)
327{
328 bool ret(false) ;
329
330 for (unsigned int i=0; i < coefList.size(); ++i) {
331 const auto& coef = coefList[i];
332 const auto& func = funcList[i];
333
334 if (func.observableOverlaps(nset, coef)) {
335 oocoutE(&caller, InputArguments) << caller.ClassName() << "::checkObservables(" << caller.GetName()
336 << "): ERROR: coefficient " << coef.GetName()
337 << " and FUNC " << func.GetName() << " have one or more observables in common" << std::endl;
338 ret = true ;
339 }
340 if (coef.dependsOn(*nset)) {
341 oocoutE(&caller, InputArguments) << caller.ClassName() << "::checkObservables(" << caller.GetName()
342 << "): ERROR coefficient " << coef.GetName()
343 << " depends on one or more of the following observables" ; nset->Print("1") ;
344 ret = true ;
345 }
346 }
347
348 return ret ;
349}
350
351
352////////////////////////////////////////////////////////////////////////////////
353/// Check if FUNC is valid for given normalization set.
354/// Coefficient and FUNC must be non-overlapping, but func-coefficient
355/// pairs may overlap each other.
356///
357/// In the present implementation, coefficients may not be observables or derive
358/// from observables.
359
361{
362 return checkObservables(*this, nset, _funcList, _coefList);
363}
364
365
366////////////////////////////////////////////////////////////////////////////////
367/// Advertise that all integrals can be handled internally.
368
370 const RooArgSet* normSet2, const char* rangeName) const
371{
372 return getAnalyticalIntegralWN(*this, _normIntMgr, _funcList, _coefList, allVars, analVars, normSet2, rangeName);
373}
374
375
377 RooArgList const& funcList, RooArgList const& /*coefList*/,
378 RooArgSet& allVars, RooArgSet& analVars,
379 const RooArgSet* normSet2, const char* rangeName)
380{
381 // Handle trivial no-integration scenario
382 if (allVars.empty()) return 0 ;
383 if (caller.getForceNumInt()) return 0 ;
384
385 // Select subset of allVars that are actual dependents
386 analVars.add(allVars) ;
387 std::unique_ptr<RooArgSet> normSet;
388 if(normSet2) {
389 normSet = std::make_unique<RooArgSet>();
390 caller.getObservables(normSet2, *normSet);
391 }
392
393
394 // Check if this configuration was created before
395 Int_t sterileIdx(-1) ;
396 auto* cache = static_cast<CacheElem*>(normIntMgr.getObj(normSet.get(),&analVars,&sterileIdx,RooNameReg::ptr(rangeName)));
397 if (cache) {
398 //cout << "RooRealSumPdf("<<this<<")::getAnalyticalIntegralWN:"<<GetName()<<"("<<allVars<<","<<analVars<<","<<(normSet2?*normSet2:RooArgSet())<<","<<(rangeName?rangeName:"<none>") << " -> " << _normIntMgr.lastIndex()+1 << " (cached)" << endl;
399 return normIntMgr.lastIndex()+1 ;
400 }
401
402 // Create new cache element
403 cache = new CacheElem ;
404
405 // Make list of function projection and normalization integrals
406 for (const auto elm : funcList) {
407 const auto func = static_cast<RooAbsReal*>(elm);
408
409 RooAbsReal* funcInt = func->createIntegral(analVars,rangeName) ;
410 if(funcInt->InheritsFrom(RooRealIntegral::Class())) ((RooRealIntegral*)funcInt)->setAllowComponentSelection(true);
411 cache->_funcIntList.addOwned(*funcInt) ;
412 if (normSet && !normSet->empty()) {
413 RooAbsReal* funcNorm = func->createIntegral(*normSet) ;
414 cache->_funcNormList.addOwned(*funcNorm) ;
415 }
416 }
417
418 // Store cache element
419 Int_t code = normIntMgr.setObj(normSet.get(),&analVars,(RooAbsCacheElement*)cache,RooNameReg::ptr(rangeName)) ;
420
421 //cout << "RooRealSumPdf("<<this<<")::getAnalyticalIntegralWN:"<<GetName()<<"("<<allVars<<","<<analVars<<","<<(normSet2?*normSet2:RooArgSet())<<","<<(rangeName?rangeName:"<none>") << " -> " << code+1 << endl;
422 return code+1 ;
423}
424
425
426
427
428////////////////////////////////////////////////////////////////////////////////
429/// Implement analytical integrations by deferring integration of component
430/// functions to integrators of components.
431
432double RooRealSumPdf::analyticalIntegralWN(Int_t code, const RooArgSet* normSet2, const char* rangeName) const
433{
434 return analyticalIntegralWN(*this, _normIntMgr, _funcList, _coefList, code, normSet2, rangeName, _haveWarned);
435}
436
438 RooArgList const& funcList, RooArgList const& coefList,
439 Int_t code, const RooArgSet* normSet2, const char* rangeName,
440 bool hasWarnedBefore)
441{
442 // Handle trivial passthrough scenario
443 if (code==0) return caller.getVal(normSet2) ;
444
445
446 // WVE needs adaptation for rangeName feature
447 auto* cache = static_cast<CacheElem*>(normIntMgr.getObjByIndex(code-1));
448 if (cache==0) { // revive the (sterilized) cache
449 //cout << "RooRealSumPdf("<<this<<")::analyticalIntegralWN:"<<GetName()<<"("<<code<<","<<(normSet2?*normSet2:RooArgSet())<<","<<(rangeName?rangeName:"<none>") << ": reviving cache "<< endl;
450 RooArgSet vars;
451 caller.getParameters(nullptr, vars);
452 RooArgSet iset = normIntMgr.selectFromSet2(vars, code-1);
453 RooArgSet nset = normIntMgr.selectFromSet1(vars, code-1);
454 RooArgSet dummy;
455 Int_t code2 = caller.getAnalyticalIntegralWN(iset,dummy,&nset,rangeName);
456 R__ASSERT(code==code2); // must have revived the right (sterilized) slot...
457 cache = (CacheElem*) normIntMgr.getObjByIndex(code-1) ;
458 R__ASSERT(cache!=0);
459 }
460
461 double value(0) ;
462
463 // N funcs, N-1 coefficients
464 double lastCoef(1) ;
465 auto funcIt = funcList.begin();
466 auto funcIntIt = cache->_funcIntList.begin();
467 for (const auto coefArg : coefList) {
468 assert(funcIt != funcList.end());
469 const auto coef = static_cast<const RooAbsReal*>(coefArg);
470 const auto func = static_cast<const RooAbsReal*>(*funcIt++);
471 const auto funcInt = static_cast<RooAbsReal*>(*funcIntIt++);
472
473 double coefVal = coef->getVal(normSet2) ;
474 if (coefVal) {
475 assert(func);
476 if (normSet2 ==0 || func->isSelectedComp()) {
477 assert(funcInt);
478 value += funcInt->getVal()*coefVal ;
479 }
480 lastCoef -= coef->getVal(normSet2) ;
481 }
482 }
483
484 const bool haveLastCoef = funcList.size() == coefList.size();
485
486 if (!haveLastCoef) {
487 // Add last func with correct coefficient
488 const auto func = static_cast<const RooAbsReal*>(*funcIt);
489 const auto funcInt = static_cast<RooAbsReal*>(*funcIntIt);
490 assert(func);
491
492 if (normSet2 ==0 || func->isSelectedComp()) {
493 assert(funcInt);
494 value += funcInt->getVal()*lastCoef ;
495 }
496
497 // Warn about coefficient degeneration
498 if (!hasWarnedBefore && (lastCoef<0 || lastCoef>1)) {
499 oocoutW(&caller, Eval) << caller.ClassName() << "::evaluate(" << caller.GetName()
500 << " WARNING: sum of FUNC coefficients not in range [0-1], value="
501 << 1-lastCoef << endl ;
502 }
503 }
504
505 double normVal(1) ;
506 if (normSet2 && normSet2->getSize()>0) {
507 normVal = 0 ;
508
509 // N funcs, N-1 coefficients
510 auto funcNormIter = cache->_funcNormList.begin();
511 for (const auto coefAsArg : coefList) {
512 auto coef = static_cast<RooAbsReal*>(coefAsArg);
513 auto funcNorm = static_cast<RooAbsReal*>(*funcNormIter++);
514
515 double coefVal = coef->getVal(normSet2);
516 if (coefVal) {
517 assert(funcNorm);
518 normVal += funcNorm->getVal()*coefVal ;
519 }
520 }
521
522 // Add last func with correct coefficient
523 if (!haveLastCoef) {
524 auto funcNorm = static_cast<RooAbsReal*>(*funcNormIter);
525 assert(funcNorm);
526
527 normVal += funcNorm->getVal()*lastCoef;
528 }
529 }
530
531 return value / normVal;
532}
533
534
535////////////////////////////////////////////////////////////////////////////////
536
538{
539 double n = getNorm(nset) ;
540 if (n<0) {
541 logEvalError("Expected number of events is negative") ;
542 }
543 return n ;
544}
545
546
547////////////////////////////////////////////////////////////////////////////////
548
549std::list<double>* RooRealSumPdf::binBoundaries(RooAbsRealLValue& obs, double xlo, double xhi) const
550{
551 return binBoundaries(_funcList, obs, xlo, xhi);
552}
553
554
555std::list<double>* RooRealSumPdf::binBoundaries(RooArgList const& funcList, RooAbsRealLValue& obs, double xlo, double xhi)
556{
557 std::list<double>* sumBinB = nullptr;
558 bool needClean(false) ;
559
560 // Loop over components pdf
561 for (auto * func : static_range_cast<RooAbsReal*>(funcList)) {
562
563 list<double>* funcBinB = func->binBoundaries(obs,xlo,xhi) ;
564
565 // Process hint
566 if (funcBinB) {
567 if (!sumBinB) {
568 // If this is the first hint, then just save it
569 sumBinB = funcBinB ;
570 } else {
571
572 std::list<double>* newSumBinB = new list<double>(sumBinB->size()+funcBinB->size()) ;
573
574 // Merge hints into temporary array
575 merge(funcBinB->begin(),funcBinB->end(),sumBinB->begin(),sumBinB->end(),newSumBinB->begin()) ;
576
577 // Copy merged array without duplicates to new sumBinBArrau
578 delete sumBinB ;
579 delete funcBinB ;
580 sumBinB = newSumBinB ;
581 needClean = true ;
582 }
583 }
584 }
585
586 // Remove consecutive duplicates
587 if (needClean) {
588 list<double>::iterator new_end = unique(sumBinB->begin(),sumBinB->end()) ;
589 sumBinB->erase(new_end,sumBinB->end()) ;
590 }
591
592 return sumBinB ;
593}
594
595
596
597/// Check if all components that depend on `obs` are binned.
599{
600 return isBinnedDistribution(_funcList, obs);
601}
602
603
605{
606 for (auto* func : static_range_cast<RooAbsReal*>(funcList)) {
607
608 if (func->dependsOn(obs) && !func->isBinnedDistribution(obs)) {
609 return false ;
610 }
611 }
612
613 return true ;
614}
615
616
617
618////////////////////////////////////////////////////////////////////////////////
619
620std::list<double>* RooRealSumPdf::plotSamplingHint(RooAbsRealLValue& obs, double xlo, double xhi) const
621{
622 return plotSamplingHint(_funcList, obs, xlo, xhi);
623}
624
625std::list<double>* RooRealSumPdf::plotSamplingHint(RooArgList const& funcList, RooAbsRealLValue& obs, double xlo, double xhi)
626{
627 std::list<double>* sumHint = nullptr;
628 bool needClean(false) ;
629
630 // Loop over components pdf
631 for (const auto elm : funcList) {
632 auto func = static_cast<RooAbsReal*>(elm);
633
634 list<double>* funcHint = func->plotSamplingHint(obs,xlo,xhi) ;
635
636 // Process hint
637 if (funcHint) {
638 if (!sumHint) {
639
640 // If this is the first hint, then just save it
641 sumHint = funcHint ;
642
643 } else {
644
645 auto* newSumHint = new std::list<double>(sumHint->size()+funcHint->size()) ;
646
647 // Merge hints into temporary array
648 merge(funcHint->begin(),funcHint->end(),sumHint->begin(),sumHint->end(),newSumHint->begin()) ;
649
650 // Copy merged array without duplicates to new sumHintArrau
651 delete sumHint ;
652 sumHint = newSumHint ;
653 needClean = true ;
654 }
655 }
656 }
657
658 // Remove consecutive duplicates
659 if (needClean) {
660 sumHint->erase(std::unique(sumHint->begin(),sumHint->end()), sumHint->end()) ;
661 }
662
663 return sumHint ;
664}
665
666
667
668
669////////////////////////////////////////////////////////////////////////////////
670/// Label OK'ed components of a RooRealSumPdf with cache-and-track
671
673{
674 setCacheAndTrackHints(_funcList, trackNodes);
675}
676
677
679{
680 for (const auto sarg : funcList) {
681 if (sarg->canNodeBeCached()==Always) {
682 trackNodes.add(*sarg) ;
683 //cout << "tracking node RealSumPdf component " << sarg->ClassName() << "::" << sarg->GetName() << endl ;
684 }
685 }
686}
687
688
689////////////////////////////////////////////////////////////////////////////////
690/// Customized printing of arguments of a RooRealSumPdf to more intuitively reflect the contents of the
691/// product operator construction
692
693void RooRealSumPdf::printMetaArgs(ostream& os) const
694{
696}
697
698
699void RooRealSumPdf::printMetaArgs(RooArgList const& funcList, RooArgList const& coefList, ostream& os)
700{
701
702 bool first(true) ;
703
704 if (!coefList.empty()) {
705 auto funcIter = funcList.begin();
706
707 for (const auto coef : coefList) {
708 if (!first) {
709 os << " + " ;
710 } else {
711 first = false ;
712 }
713 const auto func = *(funcIter++);
714 os << coef->GetName() << " * " << func->GetName();
715 }
716
717 if (funcIter != funcList.end()) {
718 os << " + [%] * " << (*funcIter)->GetName() ;
719 }
720 } else {
721
722 for (const auto func : funcList) {
723 if (!first) {
724 os << " + " ;
725 } else {
726 first = false ;
727 }
728 os << func->GetName() ;
729 }
730 }
731
732 os << " " ;
733}
734
735std::unique_ptr<RooAbsArg> RooRealSumPdf::compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext & ctx) const
736{
737 if(normSet.empty() || selfNormalized()) {
738 return RooAbsPdf::compileForNormSet({}, ctx);
739 }
740 std::unique_ptr<RooAbsPdf> pdfClone(static_cast<RooAbsPdf*>(this->Clone()));
741 ctx.compileServers(*pdfClone, {});
742
743 auto depList = new RooArgSet; // INTENTIONAL LEAK FOR NOW!
744 pdfClone->getObservables(&normSet, *depList);
745
746 auto newArg = std::make_unique<RooNormalizedPdf>(*pdfClone, *depList);
747
748 // The direct servers are this pdf and the normalization integral, which
749 // don't need to be compiled further.
750 for(RooAbsArg * server : newArg->servers()) {
751 server->setAttribute("_COMPILED");
752 }
753 newArg->setAttribute("_COMPILED");
754 newArg->addOwnedComponents(std::move(pdfClone));
755 return newArg;
756}
#define oocoutW(o, a)
Definition: RooMsgService.h:51
#define coutW(a)
Definition: RooMsgService.h:36
#define oocoutE(o, a)
Definition: RooMsgService.h:52
#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 value
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:72
RooArgSet * getObservables(const RooArgSet &set, bool valueOnly=true) const
Given a set of possible observables, return the observables that this PDF depends on.
Definition: RooAbsArg.h:296
const RefCountList_t & servers() const
List of all servers of this object.
Definition: RooAbsArg.h:199
TObject * Clone(const char *newname=nullptr) const override
Make a clone of an object using the Streamer facility.
Definition: RooAbsArg.h:84
RooArgSet * getParameters(const RooAbsData *data, bool stripDisconnected=true) const
Create a list of leaf nodes in the arg tree starting with ourself as top node that don't match any of...
Definition: RooAbsArg.cxx:563
RooAbsCacheElement is the abstract base class for objects to be stored in RooAbsCache cache manager o...
bool empty() const
Int_t getSize() const
Return the number of elements in the collection.
virtual bool add(const RooAbsArg &var, bool silent=false)
Add the specified argument to list.
const_iterator end() const
Storage_t::size_type size() const
const_iterator begin() const
void Print(Option_t *options=nullptr) const override
This method must be overridden when a class wants to print itself.
double getNorm(const RooArgSet &nset) const
Get normalisation term needed to normalise the raw values returned by getVal().
Definition: RooAbsPdf.h:255
std::unique_ptr< RooAbsArg > compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext &ctx) const override
Definition: RooAbsPdf.cxx:3640
@ CanBeExtended
Definition: RooAbsPdf.h:272
@ CanNotBeExtended
Definition: RooAbsPdf.h:272
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:60
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:104
virtual Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &analVars, const RooArgSet *normSet, const char *rangeName=nullptr) const
Variant of getAnalyticalIntegral that is also passed the normalization set that should be applied to ...
Definition: RooAbsReal.cxx:361
RooAbsReal * createIntegral(const RooArgSet &iset, const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
Create an object that represents the integral of the function over one or more observables std::liste...
Definition: RooAbsReal.cxx:523
bool getForceNumInt() const
Definition: RooAbsReal.h:191
void logEvalError(const char *message, const char *serverValueString=nullptr) const
Log evaluation error message.
virtual std::list< double > * plotSamplingHint(RooAbsRealLValue &obs, double xlo, double xhi) const
Interface for returning an optional hint for initial sampling points when constructing a curve projec...
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
Int_t setObj(const RooArgSet *nset, T *obj, const TNamed *isetRangeName=nullptr)
Setter function without integration set.
RooArgSet selectFromSet1(RooArgSet const &argSet, int index) const
Create RooArgSet contatining the objects that are both in the cached set 1.
T * getObjByIndex(Int_t index) const
Retrieve payload object by slot index.
RooArgSet selectFromSet2(RooArgSet const &argSet, int index) const
Create RooArgSet contatining the objects that are both in the cached set 2.
Int_t lastIndex() const
Return index of slot used in last get or set operation.
T * getObj(const RooArgSet *nset, Int_t *sterileIndex=nullptr, const TNamed *isetRangeName=nullptr)
Getter function without integration set.
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...
void compileServers(RooAbsArg &arg, RooArgSet const &normSet)
auto & at(RooAbsArg const *arg, RooAbsArg const *=nullptr)
Definition: DataMap.h:88
static const TNamed * ptr(const char *stringPtr)
Return a unique TNamed pointer for given C++ string.
Definition: RooNameReg.cxx:81
Class RooObjCacheManager is an implementation of class RooCacheManager<RooAbsCacheElement> and specia...
RooRealIntegral performs hybrid numerical/analytical integrals of RooAbsReal objects.
static TClass * Class()
The class RooRealSumPdf implements a PDF constructed from a sum of functions:
Definition: RooRealSumPdf.h:24
void setCacheAndTrackHints(RooArgSet &) override
Label OK'ed components of a RooRealSumPdf with cache-and-track.
static bool _doFloorGlobal
Global flag for introducing floor at zero in pdf.
Definition: RooRealSumPdf.h:92
bool _doFloor
Introduce floor at zero in pdf.
Definition: RooRealSumPdf.h:90
bool checkObservables(const RooArgSet *nset) const override
Check if FUNC is valid for given normalization set.
const RooArgList & funcList() const
Definition: RooRealSumPdf.h:45
~RooRealSumPdf() override
Destructor.
RooObjCacheManager _normIntMgr
! The integration cache manager
Definition: RooRealSumPdf.h:83
void computeBatch(cudaStream_t *, double *output, size_t size, RooFit::Detail::DataMap const &) const override
Base function for computing multiple values of a RooAbsReal.
std::unique_ptr< RooAbsArg > compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext &ctx) const override
std::list< double > * plotSamplingHint(RooAbsRealLValue &, double, double) const override
Interface for returning an optional hint for initial sampling points when constructing a curve projec...
double evaluate() const override
Calculate the current value.
ExtendMode extendMode() const override
Returns ability of PDF to provide extended likelihood terms.
RooListProxy _coefList
List of coefficients.
Definition: RooRealSumPdf.h:87
bool selfNormalized() const override
Shows if a PDF is self-normalized, which means that no attempt is made to add a normalization term.
Definition: RooRealSumPdf.h:54
RooListProxy _funcList
List of component FUNCs.
Definition: RooRealSumPdf.h:86
std::list< double > * binBoundaries(RooAbsRealLValue &, double, double) const override
Retrieve bin boundaries if this distribution is binned in obs.
RooRealSumPdf()
Default constructor coverity[UNINIT_CTOR].
void printMetaArgs(std::ostream &os) const override
Customized printing of arguments of a RooRealSumPdf to more intuitively reflect the contents of the p...
Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &numVars, const RooArgSet *normSet, const char *rangeName=nullptr) const override
Advertise that all integrals can be handled internally.
double expectedEvents(const RooArgSet *nset) const override
Return expected number of events for extended likelihood calculation, which is the sum of all coeffic...
double analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=nullptr) const override
Implement analytical integrations by deferring integration of component functions to integrators of c...
static void initializeFuncsAndCoefs(RooAbsReal const &caller, const RooArgList &inFuncList, const RooArgList &inCoefList, RooArgList &funcList, RooArgList &coefList)
const RooArgList & coefList() const
Definition: RooRealSumPdf.h:46
bool _extended
Allow use as extended p.d.f.
Definition: RooRealSumPdf.h:88
bool isBinnedDistribution(const RooArgSet &obs) const override
Check if all components that depend on obs are binned.
const char * GetName() const override
Returns name of object.
Definition: TNamed.h:47
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition: TObject.cxx:207
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition: TObject.cxx:525
const Int_t n
Definition: legend1.C:16
@ InputArguments
Definition: RooGlobalFunc.h:63
Definition: first.py:1
static double packFloatIntoNaN(float payload)
Pack float into mantissa of a NaN.
Definition: RooNaNPacker.h:109
static void output()