Logo ROOT  
Reference Guide
RooAddition.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\file RooAddition.cxx
19\class RooAddition
20\ingroup Roofitcore
21
22RooAddition calculates the sum of a set of RooAbsReal terms, or
23when constructed with two sets, it sums the product of the terms
24in the two sets. This class does not (yet) do any smart handling of integrals,
25i.e. all integrals of the product are handled numerically.
26**/
27
28
29#include "RooFit.h"
30
31#include "Riostream.h"
32#include "RooAddition.h"
33#include "RooProduct.h"
34#include "RooAbsReal.h"
35#include "RooErrorHandler.h"
36#include "RooArgSet.h"
37#include "RooNameReg.h"
38#include "RooNLLVar.h"
39#include "RooChi2Var.h"
40#include "RooMsgService.h"
41
42#include <algorithm>
43#include <cmath>
44
45using namespace std;
46
48;
49
50
51////////////////////////////////////////////////////////////////////////////////
52/// Empty constructor
54{
55}
56
57
58
59////////////////////////////////////////////////////////////////////////////////
60/// Constructor with a single set consisting of RooAbsReal.
61/// \param[in] name Name of the PDF
62/// \param[in] title Title
63/// \param[in] sumSet The value of the function will be the sum of the values in this set
64/// \param[in] takeOwnership If true, the RooAddition object will take ownership of the arguments in `sumSet`
65
66RooAddition::RooAddition(const char* name, const char* title, const RooArgList& sumSet, Bool_t takeOwnership)
67 : RooAbsReal(name, title)
68 , _set("!set","set of components",this)
69 , _cacheMgr(this,10)
70{
71 for (const auto comp : sumSet) {
72 if (!dynamic_cast<RooAbsReal*>(comp)) {
73 coutE(InputArguments) << "RooAddition::ctor(" << GetName() << ") ERROR: component " << comp->GetName()
74 << " is not of type RooAbsReal" << endl ;
76 }
77 _set.add(*comp) ;
78 if (takeOwnership) _ownedList.addOwned(*comp) ;
79 }
80
81}
82
83
84
85////////////////////////////////////////////////////////////////////////////////
86/// Constructor with two sets of RooAbsReals.
87///
88/// The sum of pair-wise products of elements in the sets will be computed:
89/// \f[
90/// A = \sum_i \mathrm{Set1}[i] * \mathrm{Set2}[i]
91/// \f]
92///
93/// \param[in] name Name of the PDF
94/// \param[in] title Title
95/// \param[in] sumSet1 Left-hand element of the pair-wise products
96/// \param[in] sumSet2 Right-hand element of the pair-wise products
97/// \param[in] takeOwnership If true, the RooAddition object will take ownership of the arguments in the `sumSets`
98///
99RooAddition::RooAddition(const char* name, const char* title, const RooArgList& sumSet1, const RooArgList& sumSet2, Bool_t takeOwnership)
100 : RooAbsReal(name, title)
101 , _set("!set","set of components",this)
102 , _cacheMgr(this,10)
103{
104 if (sumSet1.getSize() != sumSet2.getSize()) {
105 coutE(InputArguments) << "RooAddition::ctor(" << GetName() << ") ERROR: input lists should be of equal length" << endl ;
107 }
108
109 for (unsigned int i = 0; i < sumSet1.size(); ++i) {
110 const auto comp1 = &sumSet1[i];
111 const auto comp2 = &sumSet2[i];
112
113 if (!dynamic_cast<RooAbsReal*>(comp1)) {
114 coutE(InputArguments) << "RooAddition::ctor(" << GetName() << ") ERROR: component " << comp1->GetName()
115 << " in first list is not of type RooAbsReal" << endl ;
117 }
118
119 if (!dynamic_cast<RooAbsReal*>(comp2)) {
120 coutE(InputArguments) << "RooAddition::ctor(" << GetName() << ") ERROR: component " << comp2->GetName()
121 << " in first list is not of type RooAbsReal" << endl ;
123 }
124 // TODO: add flag to RooProduct c'tor to make it assume ownership...
125 TString _name(name);
126 _name.Append( "_[");
127 _name.Append(comp1->GetName());
128 _name.Append( "_x_");
129 _name.Append(comp2->GetName());
130 _name.Append( "]");
131 RooProduct *prod = new RooProduct( _name, _name , RooArgSet(*comp1, *comp2) /*, takeOwnership */ ) ;
132 _set.add(*prod);
133 _ownedList.addOwned(*prod) ;
134 if (takeOwnership) {
135 _ownedList.addOwned(*comp1) ;
136 _ownedList.addOwned(*comp2) ;
137 }
138 }
139}
140
141
142
143////////////////////////////////////////////////////////////////////////////////
144/// Copy constructor
145
146RooAddition::RooAddition(const RooAddition& other, const char* name)
147 : RooAbsReal(other, name)
148 , _set("!set",this,other._set)
149 , _cacheMgr(other._cacheMgr,this)
150{
151 // Member _ownedList is intentionally not copy-constructed -- ownership is not transferred
152}
153
154
155////////////////////////////////////////////////////////////////////////////////
156
158{ // Destructor
159
160}
161
162////////////////////////////////////////////////////////////////////////////////
163/// Calculate and return current value of self
164
166{
167 Double_t sum(0);
168 const RooArgSet* nset = _set.nset() ;
169
170// cout << "RooAddition::eval sum = " ;
171
172 for (const auto arg : _set) {
173 const auto comp = static_cast<RooAbsReal*>(arg);
174 const Double_t tmp = comp->getVal(nset);
175// cout << tmp << " " ;
176 sum += tmp ;
177 }
178// cout << " = " << sum << endl ;
179 return sum ;
180}
181
182
183////////////////////////////////////////////////////////////////////////////////
184/// Return the default error level for MINUIT error analysis
185/// If the addition contains one or more RooNLLVars and
186/// no RooChi2Vars, return the defaultErrorLevel() of
187/// RooNLLVar. If the addition contains one ore more RooChi2Vars
188/// and no RooNLLVars, return the defaultErrorLevel() of
189/// RooChi2Var. If the addition contains neither or both
190/// issue a warning message and return a value of 1
191
193{
194 RooAbsReal* nllArg(0) ;
195 RooAbsReal* chi2Arg(0) ;
196
197 RooAbsArg* arg ;
198
199 RooArgSet* comps = getComponents() ;
200 TIterator* iter = comps->createIterator() ;
201 while((arg=(RooAbsArg*)iter->Next())) {
202 if (dynamic_cast<RooNLLVar*>(arg)) {
203 nllArg = (RooAbsReal*)arg ;
204 }
205 if (dynamic_cast<RooChi2Var*>(arg)) {
206 chi2Arg = (RooAbsReal*)arg ;
207 }
208 }
209 delete iter ;
210 delete comps ;
211
212 if (nllArg && !chi2Arg) {
213 coutI(Fitting) << "RooAddition::defaultErrorLevel(" << GetName()
214 << ") Summation contains a RooNLLVar, using its error level" << endl ;
215 return nllArg->defaultErrorLevel() ;
216 } else if (chi2Arg && !nllArg) {
217 coutI(Fitting) << "RooAddition::defaultErrorLevel(" << GetName()
218 << ") Summation contains a RooChi2Var, using its error level" << endl ;
219 return chi2Arg->defaultErrorLevel() ;
220 } else if (!nllArg && !chi2Arg) {
221 coutI(Fitting) << "RooAddition::defaultErrorLevel(" << GetName() << ") WARNING: "
222 << "Summation contains neither RooNLLVar nor RooChi2Var server, using default level of 1.0" << endl ;
223 } else {
224 coutI(Fitting) << "RooAddition::defaultErrorLevel(" << GetName() << ") WARNING: "
225 << "Summation contains BOTH RooNLLVar and RooChi2Var server, using default level of 1.0" << endl ;
226 }
227
228 return 1.0 ;
229}
230
231
232////////////////////////////////////////////////////////////////////////////////
233
235{
236 for (auto arg : _set) {
237 static_cast<RooAbsReal*>(arg)->enableOffsetting(flag) ;
238 }
239}
240
241
242
243////////////////////////////////////////////////////////////////////////////////
244
246{
247 for (const auto arg : _set) {
248 static_cast<RooAbsReal*>(arg)->setData(data,cloneData) ;
249 }
250 return kTRUE ;
251}
252
253
254
255////////////////////////////////////////////////////////////////////////////////
256
257void RooAddition::printMetaArgs(ostream& os) const
258{
260 for (const auto arg : _set) {
261 if (!first) {
262 os << " + " ;
263 } else {
264 first = kFALSE ;
265 }
266 os << arg->GetName() ;
267 }
268 os << " " ;
269}
270
271////////////////////////////////////////////////////////////////////////////////
272
273Int_t RooAddition::getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* rangeName) const
274{
275 // we always do things ourselves -- actually, always delegate further down the line ;-)
276 analVars.add(allVars);
277
278 // check if we already have integrals for this combination of factors
279 Int_t sterileIndex(-1);
280 CacheElem* cache = (CacheElem*) _cacheMgr.getObj(&analVars,&analVars,&sterileIndex,RooNameReg::ptr(rangeName));
281 if (cache!=0) {
282 Int_t code = _cacheMgr.lastIndex();
283 return code+1;
284 }
285
286 // we don't, so we make it right here....
287 cache = new CacheElem;
288 for (const auto arg : _set) {// checked in c'tor that this will work...
289 RooAbsReal *I = static_cast<const RooAbsReal*>(arg)->createIntegral(analVars,rangeName);
290 cache->_I.addOwned(*I);
291 }
292
293 Int_t code = _cacheMgr.setObj(&analVars,&analVars,(RooAbsCacheElement*)cache,RooNameReg::ptr(rangeName));
294 return 1+code;
295}
296
297////////////////////////////////////////////////////////////////////////////////
298/// Calculate integral internally from appropriate integral cache
299
300Double_t RooAddition::analyticalIntegral(Int_t code, const char* rangeName) const
301{
302 // note: rangeName implicit encoded in code: see _cacheMgr.setObj in getPartIntList...
303 CacheElem *cache = (CacheElem*) _cacheMgr.getObjByIndex(code-1);
304 if (cache==0) {
305 // cache got sterilized, trigger repopulation of this slot, then try again...
306 std::unique_ptr<RooArgSet> vars( getParameters(RooArgSet()) );
307 std::unique_ptr<RooArgSet> iset( _cacheMgr.nameSet2ByIndex(code-1)->select(*vars) );
309 Int_t code2 = getAnalyticalIntegral(*iset,dummy,rangeName);
310 assert(code==code2); // must have revived the right (sterilized) slot...
311 return analyticalIntegral(code2,rangeName);
312 }
313 assert(cache!=0);
314
315 // loop over cache, and sum...
316 double result(0);
317 for (auto I : cache->_I) {
318 result += static_cast<const RooAbsReal*>(I)->getVal();
319 }
320 return result;
321
322}
323
324
325
326////////////////////////////////////////////////////////////////////////////////
327
328std::list<Double_t>* RooAddition::binBoundaries(RooAbsRealLValue& obs, Double_t xlo, Double_t xhi) const
329{
330 std::list<Double_t>* sumBinB = 0 ;
331 Bool_t needClean(kFALSE) ;
332
333 RooFIter iter = _set.fwdIterator() ;
334 RooAbsReal* func ;
335 // Loop over components pdf
336 while((func=(RooAbsReal*)iter.next())) {
337
338 std::list<Double_t>* funcBinB = func->binBoundaries(obs,xlo,xhi) ;
339
340 // Process hint
341 if (funcBinB) {
342 if (!sumBinB) {
343 // If this is the first hint, then just save it
344 sumBinB = funcBinB ;
345 } else {
346
347 std::list<Double_t>* newSumBinB = new std::list<Double_t>(sumBinB->size()+funcBinB->size()) ;
348
349 // Merge hints into temporary array
350 merge(funcBinB->begin(),funcBinB->end(),sumBinB->begin(),sumBinB->end(),newSumBinB->begin()) ;
351
352 // Copy merged array without duplicates to new sumBinBArrau
353 delete sumBinB ;
354 delete funcBinB ;
355 sumBinB = newSumBinB ;
356 needClean = kTRUE ;
357 }
358 }
359 }
360
361 // Remove consecutive duplicates
362 if (needClean) {
363 std::list<Double_t>::iterator new_end = unique(sumBinB->begin(),sumBinB->end()) ;
364 sumBinB->erase(new_end,sumBinB->end()) ;
365 }
366
367 return sumBinB ;
368}
369
370
371//_____________________________________________________________________________B
373{
374 // If all components that depend on obs are binned that so is the product
375
376 RooFIter iter = _set.fwdIterator() ;
377 RooAbsReal* func ;
378 while((func=(RooAbsReal*)iter.next())) {
379 if (func->dependsOn(obs) && !func->isBinnedDistribution(obs)) {
380 return kFALSE ;
381 }
382 }
383
384 return kTRUE ;
385}
386
387
388
389
390////////////////////////////////////////////////////////////////////////////////
391
392std::list<Double_t>* RooAddition::plotSamplingHint(RooAbsRealLValue& obs, Double_t xlo, Double_t xhi) const
393{
394 std::list<Double_t>* sumHint = 0 ;
395 Bool_t needClean(kFALSE) ;
396
397 RooFIter iter = _set.fwdIterator() ;
398 RooAbsReal* func ;
399 // Loop over components pdf
400 while((func=(RooAbsReal*)iter.next())) {
401
402 std::list<Double_t>* funcHint = func->plotSamplingHint(obs,xlo,xhi) ;
403
404 // Process hint
405 if (funcHint) {
406 if (!sumHint) {
407
408 // If this is the first hint, then just save it
409 sumHint = funcHint ;
410
411 } else {
412
413 std::list<Double_t>* newSumHint = new std::list<Double_t>(sumHint->size()+funcHint->size()) ;
414
415 // Merge hints into temporary array
416 merge(funcHint->begin(),funcHint->end(),sumHint->begin(),sumHint->end(),newSumHint->begin()) ;
417
418 // Copy merged array without duplicates to new sumHintArrau
419 delete sumHint ;
420 sumHint = newSumHint ;
421 needClean = kTRUE ;
422 }
423 }
424 }
425
426 // Remove consecutive duplicates
427 if (needClean) {
428 std::list<Double_t>::iterator new_end = unique(sumHint->begin(),sumHint->end()) ;
429 sumHint->erase(new_end,sumHint->end()) ;
430 }
431
432 return sumHint ;
433}
434
435
436
437////////////////////////////////////////////////////////////////////////////////
438/// Return list of all RooAbsArgs in cache element
439
441{
442 RooArgList ret(_I) ;
443 return ret ;
444}
445
447{
448 // Destructor
449}
450
451
static RooMathCoreReg dummy
#define coutI(a)
Definition: RooMsgService.h:31
#define coutE(a)
Definition: RooMsgService.h:34
int Int_t
Definition: RtypesCore.h:41
const Bool_t kFALSE
Definition: RtypesCore.h:88
bool Bool_t
Definition: RtypesCore.h:59
double Double_t
Definition: RtypesCore.h:55
const Bool_t kTRUE
Definition: RtypesCore.h:87
#define ClassImp(name)
Definition: Rtypes.h:365
char name[80]
Definition: TGX11.cxx:109
RooAbsArg is the common abstract base class for objects that represent a value (of arbitrary type) an...
Definition: RooAbsArg.h:71
Bool_t dependsOn(const RooAbsCollection &serverList, const RooAbsArg *ignoreArg=0, Bool_t valueOnly=kFALSE) const
Test whether we depend on (ie, are served by) any object in the specified collection.
Definition: RooAbsArg.cxx:734
friend class RooArgSet
Definition: RooAbsArg.h:551
RooArgSet * getParameters(const RooAbsData *data, Bool_t stripDisconnected=kTRUE) 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:548
RooArgSet * getComponents() const
Definition: RooAbsArg.cxx:678
RooAbsCacheElement is the abstract base class for objects to be stored in RooAbsCache cache manager o...
RooFIter fwdIterator() const R__SUGGEST_ALTERNATIVE("begin()
One-time forward iterator.
Int_t getSize() const
virtual Bool_t addOwned(RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
Storage_t::size_type size() const
TIterator * createIterator(Bool_t dir=kIterForward) const R__SUGGEST_ALTERNATIVE("begin()
TIterator-style iteration over contained elements.
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition: RooAbsData.h:39
const RooArgSet * nset() const
Definition: RooAbsProxy.h:46
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:59
virtual Double_t defaultErrorLevel() const
Definition: RooAbsReal.h:221
virtual std::list< Double_t > * plotSamplingHint(RooAbsRealLValue &, Double_t, Double_t) const
Definition: RooAbsReal.h:311
virtual std::list< Double_t > * binBoundaries(RooAbsRealLValue &, Double_t, Double_t) const
Definition: RooAbsReal.h:310
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:87
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 listed in ...
Definition: RooAbsReal.cxx:562
virtual Bool_t isBinnedDistribution(const RooArgSet &) const
Definition: RooAbsReal.h:309
virtual RooArgList containedArgs(Action)
Return list of all RooAbsArgs in cache element.
RooAddition calculates the sum of a set of RooAbsReal terms, or when constructed with two sets,...
Definition: RooAddition.h:26
void printMetaArgs(std::ostream &os) const
virtual std::list< Double_t > * plotSamplingHint(RooAbsRealLValue &, Double_t, Double_t) const
Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &numVars, const char *rangeName=0) const
Interface function getAnalyticalIntergral advertises the analytical integrals that are supported.
Bool_t isBinnedDistribution(const RooArgSet &obs) const
RooArgList _ownedList
Definition: RooAddition.h:61
virtual void enableOffsetting(Bool_t)
Double_t analyticalIntegral(Int_t code, const char *rangeName=0) const
Calculate integral internally from appropriate integral cache.
RooListProxy _set
Definition: RooAddition.h:62
Double_t evaluate() const
Calculate and return current value of self.
virtual ~RooAddition()
RooAddition()
Empty constructor.
Definition: RooAddition.cxx:53
RooObjCacheManager _cacheMgr
Definition: RooAddition.h:71
virtual std::list< Double_t > * binBoundaries(RooAbsRealLValue &, Double_t, Double_t) const
Bool_t setData(RooAbsData &data, Bool_t cloneData=kTRUE)
virtual Double_t defaultErrorLevel() const
Return the default error level for MINUIT error analysis If the addition contains one or more RooNLLV...
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgList.h:21
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:28
virtual Bool_t add(const RooAbsCollection &col, Bool_t silent=kFALSE)
Add a collection of arguments to this collection by calling add() for each element in the source coll...
Definition: RooArgSet.h:88
T * getObj(const RooArgSet *nset, Int_t *sterileIndex=0, const TNamed *isetRangeName=0)
const RooNameSet * nameSet2ByIndex(Int_t index) const
Retrieve RooNameSet associated with slot at given index.
T * getObjByIndex(Int_t index) const
Retrieve payload object by slot index.
Int_t lastIndex() const
Int_t setObj(const RooArgSet *nset, T *obj, const TNamed *isetRangeName=0)
Class RooChi2Var implements a simple chi^2 calculation from a binned dataset and a PDF.
Definition: RooChi2Var.h:25
static void softAbort()
A one-time forward iterator working on RooLinkedList or RooAbsCollection.
RooAbsArg * next()
Return next element or nullptr if at end.
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Reimplementation of standard RooArgList::add()
Class RooNLLVar implements a -log(likelihood) calculation from a dataset and a PDF.
Definition: RooNLLVar.h:26
static const TNamed * ptr(const char *stringPtr)
Return a unique TNamed pointer for given C++ string.
Definition: RooNameReg.cxx:93
RooArgSet * select(const RooArgSet &list) const
Construct a RooArgSet of objects in input 'list' whose names match to those in the internal name list...
Definition: RooNameSet.cxx:187
A RooProduct represents the product of a given set of RooAbsReal objects.
Definition: RooProduct.h:32
Iterator abstract base class.
Definition: TIterator.h:30
virtual TObject * Next()=0
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
Basic string class.
Definition: TString.h:131
TString & Append(const char *cs)
Definition: TString.h:559
#define I(x, y, z)
@ InputArguments
Definition: RooGlobalFunc.h:68
Definition: first.py:1
static long int sum(long int i)
Definition: Factory.cxx:2276