Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
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
22Implements 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
56#include <TError.h>
57
58#include <algorithm>
59#include <memory>
60#include <stdexcept>
61
62using std::list, std::endl, std::ostream;
63
65
67
68////////////////////////////////////////////////////////////////////////////////
69/// Default constructor
70/// coverity[UNINIT_CTOR]
71
72RooRealSumPdf::RooRealSumPdf() : _normIntMgr(this, 10) {}
73
74////////////////////////////////////////////////////////////////////////////////
75/// Constructor with name and title
76
77RooRealSumPdf::RooRealSumPdf(const char *name, const char *title) :
78 RooAbsPdf(name,title),
79 _normIntMgr(this,10),
80 _funcList("!funcList","List of functions",this),
81 _coefList("!coefList","List of coefficients",this),
82 _extended(false)
83{
84}
85
86////////////////////////////////////////////////////////////////////////////////
87/// Construct p.d.f consisting of \f$ \mathrm{coef}_1 * \mathrm{func}_1 + (1-\mathrm{coef}_1) * \mathrm{func}_2 \f$.
88/// The input coefficients and functions are allowed to be negative
89/// but the resulting sum is not, which is enforced at runtime.
90
91RooRealSumPdf::RooRealSumPdf(const char *name, const char *title,
92 RooAbsReal& func1, RooAbsReal& func2, RooAbsReal& coef1) :
93 RooRealSumPdf(name, title)
94{
95 // Special constructor with two functions and one coefficient
96
97 _funcList.add(func1) ;
98 _funcList.add(func2) ;
99 _coefList.add(coef1) ;
100
101}
102
103
104////////////////////////////////////////////////////////////////////////////////
105/// Constructor for a PDF from a list of functions and coefficients.
106/// It implements
107/// \f[
108/// \sum_i \mathrm{coef}_i \cdot \mathrm{func}_i,
109/// \f]
110/// if \f$ N_\mathrm{coef} = N_\mathrm{func} \f$. With `extended=true`, the coefficients can take any values. With `extended=false`,
111/// there is the danger of getting a degenerate minimisation problem because a PDF has to be normalised, which needs one degree
112/// of freedom less.
113///
114/// A plain (normalised) PDF can therefore be implemented with one less coefficient. RooFit then computes
115/// \f[
116/// \sum_i^{N-1} \mathrm{coef}_i \cdot \mathrm{func}_i + (1 - \sum_i \mathrm{coef}_i ) \cdot \mathrm{func}_N,
117/// \f]
118/// if \f$ N_\mathrm{coef} = N_\mathrm{func} - 1 \f$.
119///
120/// All coefficients and functions are allowed to be negative
121/// but the sum (*i.e.* the PDF) is not, which is enforced at runtime.
122///
123/// \param name Name of the PDF
124/// \param title Title (for plotting)
125/// \param inFuncList List of functions to sum
126/// \param inCoefList List of coefficients
127/// \param extended Interpret as extended PDF (requires equal number of functions and coefficients)
128
129RooRealSumPdf::RooRealSumPdf(const char *name, const char *title, const RooArgList &inFuncList,
130 const RooArgList &inCoefList, bool extended)
131 : RooRealSumPdf(name, title)
132{
133 setExtended(extended);
134
135 RooRealSumPdf::initializeFuncsAndCoefs(*this, inFuncList, inCoefList, _funcList, _coefList);
136}
137
138
140 const RooArgList& inFuncList, const RooArgList& inCoefList,
141 RooArgList& funcList, RooArgList& coefList)
142{
143 const std::string className = caller.ClassName();
144 const std::string constructorName = className + "::" + className;
145
146 if (!(inFuncList.size()==inCoefList.size()+1 || inFuncList.size()==inCoefList.size())) {
147 oocoutE(&caller, InputArguments) << constructorName << "(" << caller.GetName()
148 << ") number of pdfs and coefficients inconsistent, must have Nfunc=Ncoef or Nfunc=Ncoef+1" << std::endl;
149 throw std::invalid_argument(className + ": Number of PDFs and coefficients is inconsistent.");
150 }
151
152 // Constructor with N functions and N or N-1 coefs
153 for (unsigned int i = 0; i < inCoefList.size(); ++i) {
154 const auto& func = inFuncList[i];
155 const auto& coef = inCoefList[i];
156
157 if (!dynamic_cast<const RooAbsReal*>(&coef)) {
158 oocoutW(&caller, InputArguments) << constructorName << "(" << caller.GetName() << ") coefficient " << coef.GetName() << " is not of type RooAbsReal, ignored" << std::endl ;
159 continue ;
160 }
161 if (!dynamic_cast<const RooAbsReal*>(&func)) {
162 oocoutW(&caller, InputArguments) << constructorName << "(" << caller.GetName() << ") func " << func.GetName() << " is not of type RooAbsReal, ignored" << std::endl ;
163 continue ;
164 }
165 funcList.add(func) ;
166 coefList.add(coef) ;
167 }
168
169 if (inFuncList.size() == inCoefList.size() + 1) {
170 const auto& func = inFuncList[inFuncList.size()-1];
171 if (!dynamic_cast<const RooAbsReal*>(&func)) {
172 oocoutE(&caller, InputArguments) << constructorName << "(" << caller.GetName() << ") last func " << func.GetName() << " is not of type RooAbsReal, fatal error" << std::endl ;
173 throw std::invalid_argument(className + ": Function passed as is not of type RooAbsReal.");
174 }
175 funcList.add(func);
176 }
177
178}
179
180
181
182
183////////////////////////////////////////////////////////////////////////////////
184/// Copy constructor
185
187 RooAbsPdf(other,name),
188 _normIntMgr(other._normIntMgr,this),
189 _funcList("!funcList",this,other._funcList),
190 _coefList("!coefList",this,other._coefList),
191 _extended(other._extended),
192 _doFloor(other._doFloor)
193{
194
195}
196
197////////////////////////////////////////////////////////////////////////////////
198
200{
202}
203
204
206 RooArgList const& funcList,
207 RooArgList const& coefList,
208 bool doFloor,
209 bool & hasWarnedBefore)
210{
211 // Do running sum of coef/func pairs, calculate lastCoef.
212 double value = 0;
213 double sumCoeff = 0.;
214 for (unsigned int i = 0; i < funcList.size(); ++i) {
215 const auto func = static_cast<RooAbsReal*>(&funcList[i]);
216 const auto coef = static_cast<RooAbsReal*>(i < coefList.size() ? &coefList[i] : nullptr);
217 const double coefVal = coef != nullptr ? coef->getVal() : (1. - sumCoeff);
218
219 // Warn about degeneration of last coefficient
220 if (coef == nullptr && (coefVal < 0 || coefVal > 1.)) {
221 if (!hasWarnedBefore) {
222 oocoutW(&caller, Eval) << caller.ClassName() << "::evaluate(" << caller.GetName()
223 << ") WARNING: sum of FUNC coefficients not in range [0-1], value="
224 << sumCoeff << ". This means that the PDF is not properly normalised."
225 << " If the PDF was meant to be extended, provide as many coefficients as functions." << std::endl;
226 hasWarnedBefore = true;
227 }
228 // Signal that we are in an undefined region:
229 value = RooNaNPacker::packFloatIntoNaN(100.f * (coefVal < 0. ? -coefVal : coefVal - 1.));
230 }
231
232 if (func->isSelectedComp()) {
233 value += func->getVal() * coefVal;
234 }
235
236 sumCoeff += coefVal;
237 }
238
239 // Introduce floor if so requested
240 return value < 0 && doFloor ? 0.0 : value;
241}
242
243////////////////////////////////////////////////////////////////////////////////
244/// Calculate the current value
245
247{
249}
250
251
253{
254 std::span<double> output = ctx.output();
255 std::size_t nEvents = output.size();
256
257 // Do running sum of coef/func pairs, calculate lastCoef.
258 for (unsigned int j = 0; j < nEvents; ++j) {
259 output[j] = 0.0;
260 }
261
262 double sumCoeff = 0.;
263 for (unsigned int i = 0; i < _funcList.size(); ++i) {
264 const auto func = static_cast<RooAbsReal*>(&_funcList[i]);
265 const auto coef = static_cast<RooAbsReal*>(i < _coefList.size() ? &_coefList[i] : nullptr);
266 const double coefVal = coef != nullptr ? ctx.at(coef)[0] : (1. - sumCoeff);
267
268 if (func->isSelectedComp()) {
269 auto funcValues = ctx.at(func);
270 if(funcValues.size() == 1) {
271 for (unsigned int j = 0; j < nEvents; ++j) {
272 output[j] += funcValues[0] * coefVal;
273 }
274 } else {
275 for (unsigned int j = 0; j < nEvents; ++j) {
276 output[j] += funcValues[j] * coefVal;
277 }
278 }
279 }
280
281 // Warn about degeneration of last coefficient
282 if (coef == nullptr && (coefVal < 0 || coefVal > 1.)) {
283 if (!_haveWarned) {
284 coutW(Eval) << "RooRealSumPdf::doEval(" << GetName()
285 << ") WARNING: sum of FUNC coefficients not in range [0-1], value="
286 << 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 ;
287 _haveWarned = true;
288 }
289 // Signal that we are in an undefined region by handing back one NaN.
290 output[0] = RooNaNPacker::packFloatIntoNaN(100.f * (coefVal < 0. ? -coefVal : coefVal - 1.));
291 }
292
293 sumCoeff += coefVal;
294 }
295
296 // Introduce floor if so requested
297 if (_doFloor || _doFloorGlobal) {
298 for (unsigned int j = 0; j < nEvents; ++j) {
299 output[j] += std::max(0., output[j]);
300 }
301 }
302}
303
305 RooArgList const &funcList, RooArgList const &coefList)
306{
307 bool noLastCoeff = funcList.size() != coefList.size();
308 std::string sum;
309 std::string lastCoeff;
310
311 if (funcList.size() < 3) {
312 lastCoeff = "(1";
313
314 std::size_t i;
315 for (i = 0; i < coefList.size(); ++i) {
316 auto coeff = ctx.getResult(coefList[i]);
317 sum += "(" + ctx.getResult(funcList[i]) + " * " + coeff + ") +";
318 if (noLastCoeff)
319 lastCoeff += " - " + coeff;
320 }
321 lastCoeff += ")";
322
323 if (noLastCoeff) {
324 sum += "(" + ctx.getResult(funcList[i]) + " * " + lastCoeff + ")";
325 } else {
326 sum.pop_back();
327 }
328
329 } else {
330 std::string const &funcName = ctx.buildArg(funcList);
331 std::string const &coeffName = ctx.buildArg(coefList);
332 std::string const &coeffSize = std::to_string(coefList.size());
333
334 sum = ctx.getTmpVarName();
335 lastCoeff = ctx.getTmpVarName();
336 ctx.addToCodeBody(klass, "double " + sum + " = 0, " + lastCoeff + "= 0;\n");
337
338 std::string iterator = "i_" + ctx.getTmpVarName();
339 std::string subscriptExpr = "[" + iterator + "]";
340
341 std::string code = "for(int " + iterator + " = 0; " + iterator + " < " + coeffSize + "; " + iterator + "++) {\n" +
342 sum + " += " + funcName + subscriptExpr + " * " + coeffName + subscriptExpr + ";\n";
343 if (noLastCoeff)
344 code += lastCoeff + " += " + coeffName + subscriptExpr + ";\n";
345 code += "}\n";
346
347 if (noLastCoeff)
348 code += sum + " += " + funcName + "[" + coeffSize + "]" + " * (1 - " + lastCoeff + ");\n";
349 ctx.addToCodeBody(klass, code);
350 }
351
352 ctx.addResult(klass, sum);
353}
354
356{
357 translateImpl(ctx, this, _funcList, _coefList);
358}
359
361 RooArgList const& funcList, RooArgList const& coefList)
362{
363 bool ret(false) ;
364
365 for (unsigned int i=0; i < coefList.size(); ++i) {
366 const auto& coef = coefList[i];
367 const auto& func = funcList[i];
368
369 if (func.observableOverlaps(nset, coef)) {
370 oocoutE(&caller, InputArguments) << caller.ClassName() << "::checkObservables(" << caller.GetName()
371 << "): ERROR: coefficient " << coef.GetName()
372 << " and FUNC " << func.GetName() << " have one or more observables in common" << std::endl;
373 ret = true ;
374 }
375 if (coef.dependsOn(*nset)) {
376 oocoutE(&caller, InputArguments) << caller.ClassName() << "::checkObservables(" << caller.GetName()
377 << "): ERROR coefficient " << coef.GetName()
378 << " depends on one or more of the following observables" ; nset->Print("1") ;
379 ret = true ;
380 }
381 }
382
383 return ret ;
384}
385
386
387////////////////////////////////////////////////////////////////////////////////
388/// Check if FUNC is valid for given normalization set.
389/// Coefficient and FUNC must be non-overlapping, but func-coefficient
390/// pairs may overlap each other.
391///
392/// In the present implementation, coefficients may not be observables or derive
393/// from observables.
394
396{
397 return checkObservables(*this, nset, _funcList, _coefList);
398}
399
400
401////////////////////////////////////////////////////////////////////////////////
402/// Advertise that all integrals can be handled internally.
403
405 const RooArgSet* normSet2, const char* rangeName) const
406{
407 return getAnalyticalIntegralWN(*this, _normIntMgr, _funcList, _coefList, allVars, analVars, normSet2, rangeName);
408}
409
410
412 RooArgList const& funcList, RooArgList const& /*coefList*/,
413 RooArgSet& allVars, RooArgSet& analVars,
414 const RooArgSet* normSet2, const char* rangeName)
415{
416 // Handle trivial no-integration scenario
417 if (allVars.empty()) return 0 ;
418 if (caller.getForceNumInt()) return 0 ;
419
420 // Select subset of allVars that are actual dependents
421 analVars.add(allVars) ;
422 std::unique_ptr<RooArgSet> normSet;
423 if(normSet2) {
424 normSet = std::make_unique<RooArgSet>();
425 caller.getObservables(normSet2, *normSet);
426 }
427
428
429 // Check if this configuration was created before
430 Int_t sterileIdx(-1) ;
431 auto* cache = static_cast<CacheElem*>(normIntMgr.getObj(normSet.get(),&analVars,&sterileIdx,RooNameReg::ptr(rangeName)));
432 if (cache) {
433 //cout << "RooRealSumPdf("<<this<<")::getAnalyticalIntegralWN:"<<GetName()<<"("<<allVars<<","<<analVars<<","<<(normSet2?*normSet2:RooArgSet())<<","<<(rangeName?rangeName:"<none>") << " -> " << _normIntMgr.lastIndex()+1 << " (cached)" << endl;
434 return normIntMgr.lastIndex()+1 ;
435 }
436
437 // Create new cache element
438 cache = new CacheElem ;
439
440 // Make list of function projection and normalization integrals
441 for (const auto elm : funcList) {
442 const auto func = static_cast<RooAbsReal*>(elm);
443
444 std::unique_ptr<RooAbsReal> funcInt{func->createIntegral(analVars,rangeName)};
445 if(auto funcRealInt = dynamic_cast<RooRealIntegral*>(funcInt.get())) funcRealInt->setAllowComponentSelection(true);
446 cache->_funcIntList.addOwned(std::move(funcInt));
447 if (normSet && !normSet->empty()) {
448 cache->_funcNormList.addOwned(std::unique_ptr<RooAbsReal>{func->createIntegral(*normSet)});
449 }
450 }
451
452 // Store cache element
453 Int_t code = normIntMgr.setObj(normSet.get(),&analVars,(RooAbsCacheElement*)cache,RooNameReg::ptr(rangeName)) ;
454
455 //cout << "RooRealSumPdf("<<this<<")::getAnalyticalIntegralWN:"<<GetName()<<"("<<allVars<<","<<analVars<<","<<(normSet2?*normSet2:RooArgSet())<<","<<(rangeName?rangeName:"<none>") << " -> " << code+1 << endl;
456 return code+1 ;
457}
458
459
460
461
462////////////////////////////////////////////////////////////////////////////////
463/// Implement analytical integrations by deferring integration of component
464/// functions to integrators of components.
465
466double RooRealSumPdf::analyticalIntegralWN(Int_t code, const RooArgSet* normSet2, const char* rangeName) const
467{
468 return analyticalIntegralWN(*this, _normIntMgr, _funcList, _coefList, code, normSet2, rangeName, _haveWarned);
469}
470
472 RooArgList const& funcList, RooArgList const& coefList,
473 Int_t code, const RooArgSet* normSet2, const char* rangeName,
474 bool hasWarnedBefore)
475{
476 // Handle trivial passthrough scenario
477 if (code==0) return caller.getVal(normSet2) ;
478
479
480 // WVE needs adaptation for rangeName feature
481 auto* cache = static_cast<CacheElem*>(normIntMgr.getObjByIndex(code-1));
482 if (cache==nullptr) { // revive the (sterilized) cache
483 //cout << "RooRealSumPdf("<<this<<")::analyticalIntegralWN:"<<GetName()<<"("<<code<<","<<(normSet2?*normSet2:RooArgSet())<<","<<(rangeName?rangeName:"<none>") << ": reviving cache "<< endl;
484 RooArgSet vars;
485 caller.getParameters(nullptr, vars);
486 RooArgSet iset = normIntMgr.selectFromSet2(vars, code-1);
487 RooArgSet nset = normIntMgr.selectFromSet1(vars, code-1);
488 RooArgSet dummy;
489 Int_t code2 = caller.getAnalyticalIntegralWN(iset,dummy,&nset,rangeName);
490 R__ASSERT(code==code2); // must have revived the right (sterilized) slot...
491 cache = static_cast<CacheElem*>(normIntMgr.getObjByIndex(code-1)) ;
492 R__ASSERT(cache!=nullptr);
493 }
494
495 double value(0) ;
496
497 // N funcs, N-1 coefficients
498 double lastCoef(1) ;
499 auto funcIt = funcList.begin();
500 auto funcIntIt = cache->_funcIntList.begin();
501 for (const auto coefArg : coefList) {
502 assert(funcIt != funcList.end());
503 const auto coef = static_cast<const RooAbsReal*>(coefArg);
504 const auto func = static_cast<const RooAbsReal*>(*funcIt++);
505 const auto funcInt = static_cast<RooAbsReal*>(*funcIntIt++);
506
507 double coefVal = coef->getVal(normSet2) ;
508 if (coefVal) {
509 assert(func);
510 if (normSet2 ==nullptr || func->isSelectedComp()) {
511 assert(funcInt);
512 value += funcInt->getVal()*coefVal ;
513 }
514 lastCoef -= coef->getVal(normSet2) ;
515 }
516 }
517
518 const bool haveLastCoef = funcList.size() == coefList.size();
519
520 if (!haveLastCoef) {
521 // Add last func with correct coefficient
522 const auto func = static_cast<const RooAbsReal*>(*funcIt);
523 const auto funcInt = static_cast<RooAbsReal*>(*funcIntIt);
524 assert(func);
525
526 if (normSet2 ==nullptr || func->isSelectedComp()) {
527 assert(funcInt);
528 value += funcInt->getVal()*lastCoef ;
529 }
530
531 // Warn about coefficient degeneration
532 if (!hasWarnedBefore && (lastCoef<0 || lastCoef>1)) {
533 oocoutW(&caller, Eval) << caller.ClassName() << "::evaluate(" << caller.GetName()
534 << " WARNING: sum of FUNC coefficients not in range [0-1], value="
535 << 1-lastCoef << endl ;
536 }
537 }
538
539 double normVal(1) ;
540 if (normSet2 && !normSet2->empty()) {
541 normVal = 0 ;
542
543 // N funcs, N-1 coefficients
544 auto funcNormIter = cache->_funcNormList.begin();
545 for (const auto coefAsArg : coefList) {
546 auto coef = static_cast<RooAbsReal*>(coefAsArg);
547 auto funcNorm = static_cast<RooAbsReal*>(*funcNormIter++);
548
549 double coefVal = coef->getVal(normSet2);
550 if (coefVal) {
551 assert(funcNorm);
552 normVal += funcNorm->getVal()*coefVal ;
553 }
554 }
555
556 // Add last func with correct coefficient
557 if (!haveLastCoef) {
558 auto funcNorm = static_cast<RooAbsReal*>(*funcNormIter);
559 assert(funcNorm);
560
561 normVal += funcNorm->getVal()*lastCoef;
562 }
563 }
564
565 return value / normVal;
566}
567
568
569////////////////////////////////////////////////////////////////////////////////
570
572{
573 double n = getNorm(nset) ;
574 if (n<0) {
575 logEvalError("Expected number of events is negative") ;
576 }
577 return n ;
578}
579
580
581////////////////////////////////////////////////////////////////////////////////
582
583std::list<double>* RooRealSumPdf::binBoundaries(RooAbsRealLValue& obs, double xlo, double xhi) const
584{
585 return binBoundaries(_funcList, obs, xlo, xhi);
586}
587
588
589std::list<double>* RooRealSumPdf::binBoundaries(RooArgList const& funcList, RooAbsRealLValue& obs, double xlo, double xhi)
590{
591 std::list<double>* sumBinB = nullptr;
592 bool needClean(false) ;
593
594 // Loop over components pdf
595 for (auto * func : static_range_cast<RooAbsReal*>(funcList)) {
596
597 list<double>* funcBinB = func->binBoundaries(obs,xlo,xhi) ;
598
599 // Process hint
600 if (funcBinB) {
601 if (!sumBinB) {
602 // If this is the first hint, then just save it
603 sumBinB = funcBinB ;
604 } else {
605
606 std::list<double>* newSumBinB = new list<double>(sumBinB->size()+funcBinB->size()) ;
607
608 // Merge hints into temporary array
609 merge(funcBinB->begin(),funcBinB->end(),sumBinB->begin(),sumBinB->end(),newSumBinB->begin()) ;
610
611 // Copy merged array without duplicates to new sumBinBArrau
612 delete sumBinB ;
613 delete funcBinB ;
614 sumBinB = newSumBinB ;
615 needClean = true ;
616 }
617 }
618 }
619
620 // Remove consecutive duplicates
621 if (needClean) {
622 list<double>::iterator new_end = unique(sumBinB->begin(),sumBinB->end()) ;
623 sumBinB->erase(new_end,sumBinB->end()) ;
624 }
625
626 return sumBinB ;
627}
628
629
630
631/// Check if all components that depend on `obs` are binned.
633{
634 return isBinnedDistribution(_funcList, obs);
635}
636
637
639{
640 for (auto* func : static_range_cast<RooAbsReal*>(funcList)) {
641
642 if (func->dependsOn(obs) && !func->isBinnedDistribution(obs)) {
643 return false ;
644 }
645 }
646
647 return true ;
648}
649
650
651
652////////////////////////////////////////////////////////////////////////////////
653
654std::list<double>* RooRealSumPdf::plotSamplingHint(RooAbsRealLValue& obs, double xlo, double xhi) const
655{
656 return plotSamplingHint(_funcList, obs, xlo, xhi);
657}
658
659std::list<double>* RooRealSumPdf::plotSamplingHint(RooArgList const& funcList, RooAbsRealLValue& obs, double xlo, double xhi)
660{
661 std::list<double>* sumHint = nullptr;
662 bool needClean(false) ;
663
664 // Loop over components pdf
665 for (const auto elm : funcList) {
666 auto func = static_cast<RooAbsReal*>(elm);
667
668 list<double>* funcHint = func->plotSamplingHint(obs,xlo,xhi) ;
669
670 // Process hint
671 if (funcHint) {
672 if (!sumHint) {
673
674 // If this is the first hint, then just save it
675 sumHint = funcHint ;
676
677 } else {
678
679 auto* newSumHint = new std::list<double>(sumHint->size()+funcHint->size()) ;
680
681 // the lists must be sorted before merging them
682 funcHint->sort();
683 sumHint->sort();
684 // Merge hints into temporary array
685 merge(funcHint->begin(),funcHint->end(),sumHint->begin(),sumHint->end(),newSumHint->begin()) ;
686
687 // Copy merged array without duplicates to new sumHintArrau
688 delete sumHint ;
689 sumHint = newSumHint ;
690 needClean = true ;
691 }
692 }
693 }
694
695 // Remove consecutive duplicates
696 if (needClean) {
697 sumHint->erase(std::unique(sumHint->begin(),sumHint->end()), sumHint->end()) ;
698 }
699
700 return sumHint ;
701}
702
703
704
705
706////////////////////////////////////////////////////////////////////////////////
707/// Label OK'ed components of a RooRealSumPdf with cache-and-track
708
710{
711 setCacheAndTrackHints(_funcList, trackNodes);
712}
713
714
716{
717 for (const auto sarg : funcList) {
718 if (sarg->canNodeBeCached()==Always) {
719 trackNodes.add(*sarg) ;
720 //cout << "tracking node RealSumPdf component " << sarg->ClassName() << "::" << sarg->GetName() << endl ;
721 }
722 }
723}
724
725
726////////////////////////////////////////////////////////////////////////////////
727/// Customized printing of arguments of a RooRealSumPdf to more intuitively reflect the contents of the
728/// product operator construction
729
730void RooRealSumPdf::printMetaArgs(ostream& os) const
731{
733}
734
735
736void RooRealSumPdf::printMetaArgs(RooArgList const& funcList, RooArgList const& coefList, ostream& os)
737{
738
739 bool first(true) ;
740
741 if (!coefList.empty()) {
742 auto funcIter = funcList.begin();
743
744 for (const auto coef : coefList) {
745 if (!first) {
746 os << " + " ;
747 } else {
748 first = false ;
749 }
750 const auto func = *(funcIter++);
751 os << coef->GetName() << " * " << func->GetName();
752 }
753
754 if (funcIter != funcList.end()) {
755 os << " + [%] * " << (*funcIter)->GetName() ;
756 }
757 } else {
758
759 for (const auto func : funcList) {
760 if (!first) {
761 os << " + " ;
762 } else {
763 first = false ;
764 }
765 os << func->GetName() ;
766 }
767 }
768
769 os << " " ;
770}
771
772std::unique_ptr<RooAbsArg>
774{
775 if (normSet.empty() || selfNormalized()) {
776 return RooAbsPdf::compileForNormSet({}, ctx);
777 }
778 std::unique_ptr<RooAbsPdf> pdfClone(static_cast<RooAbsPdf *>(this->Clone()));
779
780 if (ctx.likelihoodMode() && pdfClone->getAttribute("BinnedLikelihood")) {
781
782 // This has to be done before compiling the servers, such that the
783 // RooBinWidthFunctions know to disable themselves.
784 ctx.setBinnedLikelihoodMode(true);
785
786 ctx.markAsCompiled(*pdfClone);
787 ctx.compileServers(*pdfClone, {});
788
789 pdfClone->setAttribute("BinnedLikelihoodActive");
790 // If this is a binned likelihood, we're flagging it in the context.
791 // Then, the RooBinWidthFunctions know that they should not put
792 // themselves in the computation graph. Like this, the pdf values can
793 // directly be interpreted as yields, without multiplying them with the
794 // bin widths again in the NLL. However, the NLL class has to be careful
795 // to only skip the bin with multiplication when there actually were
796 // RooBinWidthFunctions! This is not the case for old workspace before
797 // ROOT 6.26. Therefore, we use the "BinnedLikelihoodActiveYields"
798 // attribute to let the NLL know what it should do.
799 if (ctx.binWidthFuncFlag()) {
800 pdfClone->setAttribute("BinnedLikelihoodActiveYields");
801 }
802 return pdfClone;
803 }
804
805 ctx.compileServers(*pdfClone, {});
806
807 RooArgSet depList;
808 pdfClone->getObservables(&normSet, depList);
809
810 auto newArg = std::make_unique<RooNormalizedPdf>(*pdfClone, depList);
811
812 // The direct servers are this pdf and the normalization integral, which
813 // don't need to be compiled further.
814 for (RooAbsArg *server : newArg->servers()) {
815 ctx.markAsCompiled(*server);
816 }
817 ctx.markAsCompiled(*newArg);
818 newArg->addOwnedComponents(std::move(pdfClone));
819 return newArg;
820}
821
822std::unique_ptr<RooAbsReal> RooRealSumPdf::createExpectedEventsFunc(const RooArgSet *nset) const
823{
824 if (nset == nullptr)
825 return nullptr;
826 return std::unique_ptr<RooAbsReal>{createIntegral(*nset, *getIntegratorConfig(), normRange())};
827}
#define oocoutW(o, a)
#define coutW(a)
#define oocoutE(o, a)
bool _extended
Definition RooNLLVar.h:43
#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 value
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
RooFit::OwningPtr< 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...
RooFit::OwningPtr< RooArgSet > getObservables(const RooArgSet &set, bool valueOnly=true) const
Given a set of possible observables, return the observables that this PDF depends on.
TObject * Clone(const char *newname=nullptr) const override
Make a clone of an object using the Streamer facility.
Definition RooAbsArg.h:89
Abstract base class for objects to be stored in RooAbsCache cache manager objects.
virtual bool add(const RooAbsArg &var, bool silent=false)
Add the specified argument to list.
Storage_t::size_type size() const
TIterator begin()
TIterator end() and range-based for loops.")
void Print(Option_t *options=nullptr) const override
This method must be overridden when a class wants to print itself.
Abstract interface for all probability density functions.
Definition RooAbsPdf.h:40
double getNorm(const RooArgSet &nset) const
Get normalisation term needed to normalise the raw values returned by getVal().
Definition RooAbsPdf.h:196
std::unique_ptr< RooAbsArg > compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext &ctx) const override
@ CanBeExtended
Definition RooAbsPdf.h:213
@ CanNotBeExtended
Definition RooAbsPdf.h:213
const char * normRange() const
Definition RooAbsPdf.h:251
Abstract base class for objects that represent a real value that may appear on the left hand side of ...
Abstract base class for objects that represent a real value and implements functionality common to al...
Definition RooAbsReal.h:59
RooFit::OwningPtr< RooAbsReal > createIntegral(const RooArgSet &iset, 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
Create an object that represents the integral of the function over one or more observables listed in ...
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition RooAbsReal.h:103
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 ...
bool getForceNumInt() const
Definition RooAbsReal.h:174
void logEvalError(const char *message, const char *serverValueString=nullptr) const
Log evaluation error message.
const RooNumIntConfig * getIntegratorConfig() const
Return the numeric integration configuration used for this object.
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:55
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 containing the objects that are both in the cached set 1 with a given index and an i...
T * getObjByIndex(Int_t index) const
Retrieve payload object by slot index.
RooArgSet selectFromSet2(RooArgSet const &argSet, int index) const
Create RooArgSet containing the objects that are both in the cached set 2 with a given index and an i...
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...
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.
void addToCodeBody(RooAbsArg const *klass, std::string const &in)
Adds the input string to the squashed code body.
std::string getTmpVarName() const
Get a unique variable name to be used in the generated code.
std::string const & getResult(RooAbsArg const &arg)
Gets the result for the given node using the node name.
std::string buildArg(RooAbsCollection const &x)
Function to save a RooListProxy as an array in the squashed code.
void markAsCompiled(RooAbsArg &arg) const
void compileServers(RooAbsArg &arg, RooArgSet const &normSet)
std::span< const double > at(RooAbsArg const *arg, RooAbsArg const *caller=nullptr)
std::span< double > output()
static const TNamed * ptr(const char *stringPtr)
Return a unique TNamed pointer for given C++ string.
Implementation of a RooCacheManager<RooAbsCacheElement> that specializes in the storage of cache elem...
Performs hybrid numerical/analytical integrals of RooAbsReal objects.
Implements a PDF constructed from a sum of functions:
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.
bool _doFloor
Introduce floor at zero in pdf.
std::unique_ptr< RooAbsReal > createExpectedEventsFunc(const RooArgSet *nset) const override
Returns an object that represents the expected number of events for a given normalization set,...
bool checkObservables(const RooArgSet *nset) const override
Check if FUNC is valid for given normalization set.
void setExtended(bool extended)
const RooArgList & funcList() const
RooObjCacheManager _normIntMgr
! The integration cache manager
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.
bool selfNormalized() const override
Shows if a PDF is self-normalized, which means that no attempt is made to add a normalization term.
RooListProxy _funcList
List of component FUNCs.
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...
void doEval(RooFit::EvalContext &) const override
Base function for computing multiple values of a RooAbsReal.
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 translateImpl(RooFit::Detail::CodeSquashContext &ctx, RooAbsArg const *klass, RooArgList const &funcList, RooArgList const &coefList)
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 ...
static void initializeFuncsAndCoefs(RooAbsReal const &caller, const RooArgList &inFuncList, const RooArgList &inCoefList, RooArgList &funcList, RooArgList &coefList)
const RooArgList & coefList() const
bool _extended
Allow use as extended p.d.f.
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
const Int_t n
Definition legend1.C:16
__roodevice__ static __roohost__ double packFloatIntoNaN(float payload)
Pack float into mantissa of a NaN.
static uint64_t sum(uint64_t i)
Definition Factory.cxx:2345
static void output()