Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
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 "Riostream.h"
29#include "RooAddition.h"
30#include "RooRealSumFunc.h"
31#include "RooRealSumPdf.h"
32#include "RooProduct.h"
33#include "RooErrorHandler.h"
34#include "RooArgSet.h"
35#include "RooNameReg.h"
36#include "RooNLLVar.h"
37#include "RooNLLVarNew.h"
38#include "RooChi2Var.h"
39#include "RooMsgService.h"
40#include "RooBatchCompute.h"
41
42#include <algorithm>
43#include <cmath>
44
46
47
48////////////////////////////////////////////////////////////////////////////////
49/// Constructor with a single set consisting of RooAbsReal.
50/// \param[in] name Name of the PDF
51/// \param[in] title Title
52/// \param[in] sumSet The value of the function will be the sum of the values in this set
53/// \param[in] takeOwnership If true, the RooAddition object will take ownership of the arguments in `sumSet`
54
55RooAddition::RooAddition(const char* name, const char* title, const RooArgList& sumSet
56#ifndef ROOFIT_MEMORY_SAFE_INTERFACES
57 , bool takeOwnership
58#endif
59 )
60 : RooAbsReal(name, title)
61 , _set("!set","set of components",this)
62 , _cacheMgr(this,10)
63{
64 for (RooAbsArg *comp : sumSet) {
65 if (!dynamic_cast<RooAbsReal*>(comp)) {
66 coutE(InputArguments) << "RooAddition::ctor(" << GetName() << ") ERROR: component " << comp->GetName()
67 << " is not of type RooAbsReal" << std::endl;
69 }
70 _set.add(*comp) ;
71#ifndef ROOFIT_MEMORY_SAFE_INTERFACES
72 if (takeOwnership) _ownedList.addOwned(std::unique_ptr<RooAbsArg>{comp});
73#endif
74 }
75
76}
77
78
79
80////////////////////////////////////////////////////////////////////////////////
81/// Constructor with two sets of RooAbsReals.
82///
83/// The sum of pair-wise products of elements in the sets will be computed:
84/// \f[
85/// A = \sum_i \mathrm{Set1}[i] * \mathrm{Set2}[i]
86/// \f]
87///
88/// \param[in] name Name of the PDF
89/// \param[in] title Title
90/// \param[in] sumSet1 Left-hand element of the pair-wise products
91/// \param[in] sumSet2 Right-hand element of the pair-wise products
92/// \param[in] takeOwnership If true, the RooAddition object will take ownership of the arguments in the `sumSets`
93///
94RooAddition::RooAddition(const char* name, const char* title, const RooArgList& sumSet1, const RooArgList& sumSet2
95#ifndef ROOFIT_MEMORY_SAFE_INTERFACES
96 , bool takeOwnership
97#endif
98 )
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" << std::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" << std::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" << std::endl;
122 }
123 TString _name(name);
124 _name.Append( "_[");
125 _name.Append(comp1->GetName());
126 _name.Append( "_x_");
127 _name.Append(comp2->GetName());
128 _name.Append( "]");
129 auto prod = std::make_unique<RooProduct>( _name, _name , RooArgSet(*comp1, *comp2));
130 _set.add(*prod);
131 _ownedList.addOwned(std::move(prod));
132#ifndef ROOFIT_MEMORY_SAFE_INTERFACES
133 if (takeOwnership) {
134 _ownedList.addOwned(std::unique_ptr<RooAbsArg>{comp1});
135 _ownedList.addOwned(std::unique_ptr<RooAbsArg>{comp2});
136 }
137#endif
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/// Calculate and return current value of self
156
158{
159 double sum(0);
160 const RooArgSet* nset = _set.nset() ;
161
162 for (auto* comp : static_range_cast<RooAbsReal*>(_set)) {
163 const double tmp = comp->getVal(nset);
164 sum += tmp ;
165 }
166 return sum ;
167}
168
169
170////////////////////////////////////////////////////////////////////////////////
171/// Compute addition of PDFs in batches.
172void RooAddition::computeBatch(double* output, size_t nEvents, RooFit::Detail::DataMap const& dataMap) const
173{
176 pdfs.reserve(_set.size());
177 coefs.reserve(_set.size());
178 for (const auto arg : _set)
179 {
180 pdfs.push_back(dataMap.at(arg));
181 coefs.push_back(1.0);
182 }
183 RooBatchCompute::compute(dataMap.config(this), RooBatchCompute::AddPdf, output, nEvents, pdfs, coefs);
184}
185
186////////////////////////////////////////////////////////////////////////////////
187
189{
190 // If the number of elements to sum is less than 3, just build a sum expression.
191 // else build a loop to sum over the values.
192 unsigned int eleSize = _set.size();
193 std::string result;
194 if (eleSize > 3) {
195 std::string className = GetName();
196 std::string varName = "elements" + className;
197 std::string sumName = "sum" + className;
198 std::string code = "";
199 std::string decl = "double " + varName + "[" + std::to_string(eleSize) + "]{";
200 int idx = 0;
201 for (RooAbsArg *it : _set) {
202 decl += ctx.getResult(*it) + ",";
203 ctx.addResult(it, varName + "[" + std::to_string(idx) + "]");
204 idx++;
205 }
206 decl.back() = '}';
207 code += decl + ";\n";
208
209 ctx.addToGlobalScope("double " + sumName + " = 0;\n");
210 std::string iterator = "i_" + className;
211 code += "for(int " + iterator + " = 0; " + iterator + " < " + std::to_string(eleSize) + "; " + iterator +
212 "++) {\n" + sumName + " += " + varName + "[" + iterator + "];\n}\n";
213 result = sumName;
214 ctx.addResult(this, result);
215
216 ctx.addToCodeBody(this, code);
217 }
218
219 result = "(";
220 for (RooAbsArg *it : _set) {
221 result += ctx.getResult(*it) + '+';
222 }
223 result.back() = ')';
224 ctx.addResult(this, result);
225}
226
227////////////////////////////////////////////////////////////////////////////////
228/// Return the default error level for MINUIT error analysis
229/// If the addition contains one or more RooNLLVars and
230/// no RooChi2Vars, return the defaultErrorLevel() of
231/// RooNLLVar. If the addition contains one ore more RooChi2Vars
232/// and no RooNLLVars, return the defaultErrorLevel() of
233/// RooChi2Var. If the addition contains neither or both
234/// issue a warning message and return a value of 1
235
237{
238 RooAbsReal* nllArg(nullptr) ;
239 RooAbsReal* chi2Arg(nullptr) ;
240
241 std::unique_ptr<RooArgSet> comps{getComponents()};
242 for(RooAbsArg * arg : *comps) {
243 if (dynamic_cast<RooNLLVar*>(arg) || dynamic_cast<RooNLLVarNew*>(arg)) {
244 nllArg = (RooAbsReal*)arg ;
245 }
246 if (dynamic_cast<RooChi2Var*>(arg)) {
247 chi2Arg = (RooAbsReal*)arg ;
248 }
249 }
250
251 if (nllArg && !chi2Arg) {
252 coutI(Fitting) << "RooAddition::defaultErrorLevel(" << GetName()
253 << ") Summation contains a RooNLLVar, using its error level" << std::endl;
254 return nllArg->defaultErrorLevel() ;
255 } else if (chi2Arg && !nllArg) {
256 coutI(Fitting) << "RooAddition::defaultErrorLevel(" << GetName()
257 << ") Summation contains a RooChi2Var, using its error level" << std::endl;
258 return chi2Arg->defaultErrorLevel() ;
259 } else if (!nllArg && !chi2Arg) {
260 coutI(Fitting) << "RooAddition::defaultErrorLevel(" << GetName() << ") WARNING: "
261 << "Summation contains neither RooNLLVar nor RooChi2Var server, using default level of 1.0" << std::endl;
262 } else {
263 coutI(Fitting) << "RooAddition::defaultErrorLevel(" << GetName() << ") WARNING: "
264 << "Summation contains BOTH RooNLLVar and RooChi2Var server, using default level of 1.0" << std::endl;
265 }
266
267 return 1.0 ;
268}
269
270
271
272////////////////////////////////////////////////////////////////////////////////
273
275{
276 for (const auto arg : _set) {
277 static_cast<RooAbsReal*>(arg)->setData(data,cloneData) ;
278 }
279 return true ;
280}
281
282
283
284////////////////////////////////////////////////////////////////////////////////
285
286void RooAddition::printMetaArgs(std::ostream& os) const
287{
288 // We can use the implementation of RooRealSumPdf with an empty coefficient list.
289 static const RooArgList coefs{};
291}
292
293////////////////////////////////////////////////////////////////////////////////
294
295Int_t RooAddition::getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* rangeName) const
296{
297 // we always do things ourselves -- actually, always delegate further down the line ;-)
298 analVars.add(allVars);
299
300 // check if we already have integrals for this combination of factors
301 Int_t sterileIndex(-1);
302 CacheElem* cache = (CacheElem*) _cacheMgr.getObj(&analVars,&analVars,&sterileIndex,RooNameReg::ptr(rangeName));
303 if (cache!=nullptr) {
304 Int_t code = _cacheMgr.lastIndex();
305 return code+1;
306 }
307
308 // we don't, so we make it right here....
309 cache = new CacheElem;
310 for (auto *arg : static_range_cast<RooAbsReal const*>(_set)) {// checked in c'tor that this will work...
311 cache->_I.addOwned(std::unique_ptr<RooAbsReal>{arg->createIntegral(analVars,rangeName)});
312 }
313
314 Int_t code = _cacheMgr.setObj(&analVars,&analVars,(RooAbsCacheElement*)cache,RooNameReg::ptr(rangeName));
315 return 1+code;
316}
317
318////////////////////////////////////////////////////////////////////////////////
319/// Calculate integral internally from appropriate integral cache
320
321double RooAddition::analyticalIntegral(Int_t code, const char* rangeName) const
322{
323 // note: rangeName implicit encoded in code: see _cacheMgr.setObj in getPartIntList...
324 CacheElem *cache = (CacheElem*) _cacheMgr.getObjByIndex(code-1);
325 if (cache==nullptr) {
326 // cache got sterilized, trigger repopulation of this slot, then try again...
327 std::unique_ptr<RooArgSet> vars( getParameters(RooArgSet()) );
328 RooArgSet iset = _cacheMgr.selectFromSet2(*vars, code-1);
329 RooArgSet dummy;
330 Int_t code2 = getAnalyticalIntegral(iset,dummy,rangeName);
331 assert(code==code2); // must have revived the right (sterilized) slot...
332 return analyticalIntegral(code2,rangeName);
333 }
334 assert(cache!=nullptr);
335
336 // loop over cache, and sum...
337 double result(0);
338 for (auto I : cache->_I) {
339 result += static_cast<const RooAbsReal*>(I)->getVal();
340 }
341 return result;
342
343}
344
345
346////////////////////////////////////////////////////////////////////////////////
347
348std::list<double>* RooAddition::binBoundaries(RooAbsRealLValue& obs, double xlo, double xhi) const
349{
350 return RooRealSumPdf::binBoundaries(_set, obs, xlo, xhi);
351}
352
353
355{
357}
358
359
360////////////////////////////////////////////////////////////////////////////////
361
362std::list<double>* RooAddition::plotSamplingHint(RooAbsRealLValue& obs, double xlo, double xhi) const
363{
364 return RooRealSumPdf::plotSamplingHint(_set, obs, xlo, xhi);
365}
#define coutI(a)
#define coutE(a)
#define ClassImp(name)
Definition Rtypes.h:377
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t result
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:80
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 > getComponents() const
Create a RooArgSet with all components (branch nodes) of the expression tree headed by this object.
RooAbsCacheElement is the abstract base class for objects to be stored in RooAbsCache cache manager o...
Int_t getSize() const
Return the number of elements in the collection.
const char * GetName() const override
Returns name of object.
virtual bool add(const RooAbsArg &var, bool silent=false)
Add the specified argument to list.
Storage_t::size_type size() const
virtual bool addOwned(RooAbsArg &var, bool silent=false)
Add an argument and transfer the ownership to the collection.
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition RooAbsData.h:57
const RooArgSet * nset() const
Definition RooAbsProxy.h:52
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
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition RooAbsReal.h:103
virtual double defaultErrorLevel() const
Definition RooAbsReal.h:248
RooAddition calculates the sum of a set of RooAbsReal terms, or when constructed with two sets,...
Definition RooAddition.h:27
RooArgList _ownedList
List of owned components.
Definition RooAddition.h:71
Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &numVars, const char *rangeName=nullptr) const override
Interface function getAnalyticalIntergral advertises the analytical integrals that are supported.
void computeBatch(double *output, size_t nEvents, RooFit::Detail::DataMap const &) const override
Compute addition of PDFs in batches.
std::list< double > * binBoundaries(RooAbsRealLValue &, double, double) const override
Retrieve bin boundaries if this distribution is binned in obs.
RooListProxy _set
set of terms to be summed
Definition RooAddition.h:72
void printMetaArgs(std::ostream &os) const override
bool setData(RooAbsData &data, bool cloneData=true) override
double analyticalIntegral(Int_t code, const char *rangeName=nullptr) const override
Calculate integral internally from appropriate integral cache.
RooObjCacheManager _cacheMgr
! The cache manager
Definition RooAddition.h:80
double defaultErrorLevel() const override
Return the default error level for MINUIT error analysis If the addition contains one or more RooNLLV...
double evaluate() const override
Calculate and return current value of self.
std::list< double > * plotSamplingHint(RooAbsRealLValue &, double, double) const override
Interface for returning an optional hint for initial sampling points when constructing a curve projec...
bool isBinnedDistribution(const RooArgSet &obs) const override
Tests if the distribution is binned. Unless overridden by derived classes, this always returns false.
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 ...
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.
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.
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.
RooChi2Var implements a simple calculation from a binned dataset and a PDF.
Definition RooChi2Var.h:25
bool add(const RooAbsArg &var, bool valueServer, bool shapeServer, bool silent)
Overloaded RooCollection_t::add() method insert object into set and registers object as server to own...
static void softAbort()
Soft abort function that interrupts macro execution but doesn't kill ROOT.
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 const & getResult(RooAbsArg const &arg)
Gets the result for the given node using the node name.
void addToGlobalScope(std::string const &str)
Adds the given string to the string block that will be emitted at the top of the squashed function.
RooBatchCompute::Config config(RooAbsArg const *arg) const
Definition DataMap.cxx:40
std::span< const double > at(RooAbsArg const *arg, RooAbsArg const *caller=nullptr)
Definition DataMap.cxx:22
This is a simple class designed to produce the nll values needed by the fitter.
Class RooNLLVar implements a -log(likelihood) calculation from a dataset and a PDF.
Definition RooNLLVar.h:25
static const TNamed * ptr(const char *stringPtr)
Return a unique TNamed pointer for given C++ string.
std::list< double > * plotSamplingHint(RooAbsRealLValue &, double, double) const override
Interface for returning an optional hint for initial sampling points when constructing a curve projec...
std::list< double > * binBoundaries(RooAbsRealLValue &, double, double) const override
Retrieve bin boundaries if this distribution is binned in obs.
void printMetaArgs(std::ostream &os) const override
Customized printing of arguments of a RooRealSumPdf to more intuitively reflect the contents of the p...
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
Basic string class.
Definition TString.h:139
TString & Append(const char *cs)
Definition TString.h:576
#define I(x, y, z)
std::vector< std::span< const double > > VarVector
std::vector< double > ArgVector
void compute(Config cfg, Computer comp, RestrictArr output, size_t size, const VarVector &vars, ArgVector &extraArgs)
static uint64_t sum(uint64_t i)
Definition Factory.cxx:2345
static void output()