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