Logo ROOT  
Reference Guide
RooResolutionModel.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/**
19 * \class RooResolutionModel
20 * RooResolutionModel is the base class for PDFs that represent a
21 * resolution model that can be convoluted with a physics model of the form
22 * \f[
23 * \mathrm{Phys}(x,a,b) = \sum_k \mathrm{coef}_k(a) * \mathrm{basis}_k(x,b)
24 * \f]
25 * where basis_k are a limited number of functions in terms of the variable
26 * to be convoluted and coef_k are coefficients independent of the convolution
27 * variable.
28 *
29 * Classes derived from RooResolutionModel implement
30 * \f[
31 * R_k(x,\bar{b},\bar{c}) = \int \mathrm{basis}_k(x',\bar{b}) * \mathrm{resModel}(x-x',\bar{c}) \; \mathrm{d} x',
32 * \f]
33 * which RooAbsAnaConvPdf uses to construct the pdf for [ Phys (x) R ] :
34 * \f[
35 * \mathrm{PDF}(x,\bar a, \bar b, \bar c) = \sum_k \mathrm{coef}_k(\bar a) * R_k(x, \bar b, \bar c)
36 * \f]
37 *
38 * A minimal implementation of a RooResolutionModel consists of a
39 * ```
40 * Int_t basisCode(const char* name)
41 * ```
42 * function indicating which basis functions this resolution model supports, and
43 * ```
44 * double evaluate(),
45 * ```
46 * which should implement the resolution model (optionally convoluted with one of the
47 * supported basis functions). RooResolutionModel objects can be used as regular
48 * PDFs (They inherit from RooAbsPdf), or as resolution model convoluted with
49 * a basis function. The implementation of evaluate() can identify the requested
50 * mode using basisCode(). If zero, the regular PDF value
51 * should be calculated. If non-zero, the model's value convoluted with the
52 * basis function identified by the code should be calculated.
53 *
54 * Optionally, analytical integrals can be advertised and implemented, in the
55 * same way as done for regular PDFS (see RooAbsPdf for further details).
56 * Also in getAnalyticalIntegral() / analyticalIntegral(), the implementation
57 * should use basisCode() to determine for which scenario the integral is
58 * requested.
59 *
60 * The choice of basis returned by basisCode() is guaranteed not to change
61 * during the lifetime of a RooResolutionModel object.
62 *
63 */
64
65#include "TClass.h"
66#include "TMath.h"
67#include "Riostream.h"
68#include "RooResolutionModel.h"
69#include "RooMsgService.h"
70
71using namespace std;
72
74
75
76////////////////////////////////////////////////////////////////////////////////
77/// Constructor with convolution variable 'x'.
78/// The convolution variable needs to be convertable to real values, and be able
79/// to give information about its range. This is supported by e.g. RooRealVar or RooLinearVar, which
80/// accepts offsetting and scaling an observable.
81RooResolutionModel::RooResolutionModel(const char *name, const char *title, RooAbsRealLValue& _x) :
82 RooAbsPdf(name,title),
83 x("x","Dependent or convolution variable",this,_x),
84 _basisCode(0), _basis(0),
85 _ownBasis(false)
86{
87
88}
89
90
91
92////////////////////////////////////////////////////////////////////////////////
93/// Copy constructor
94
96 RooAbsPdf(other,name),
97 x("x",this,other.x),
98 _basisCode(other._basisCode), _basis(0),
99 _ownBasis(false)
100{
101 if (other._basis) {
102 _basis = (RooFormulaVar*) other._basis->Clone() ;
103 _ownBasis = true ;
104 //_basis = other._basis ;
105 }
106
107 if (_basis) {
108 TIterator* bsIter = _basis->serverIterator() ;
109 RooAbsArg* basisServer ;
110 while((basisServer = (RooAbsArg*)bsIter->Next())) {
111 addServer(*basisServer,true,false) ;
112 }
113 delete bsIter ;
114 }
115}
116
117
118
119////////////////////////////////////////////////////////////////////////////////
120/// Destructor
121
123{
124 if (_ownBasis && _basis) {
125 delete _basis ;
126 }
127}
128
129
130
131////////////////////////////////////////////////////////////////////////////////
132/// Return identity formula pointer
133
135{
136 static RooFormulaVar identity("identity","1",RooArgSet(""));
137 return &identity;
138}
139
140
141
142////////////////////////////////////////////////////////////////////////////////
143/// Instantiate a clone of this resolution model representing a convolution with given
144/// basis function. The owners object name is incorporated in the clones name
145/// to avoid multiple convolution objects with the same name in complex PDF structures.
146///
147/// Note: The 'inBasis' formula expression must be a RooFormulaVar that encodes the formula
148/// in the title of the object and this expression must be an exact match against the
149/// implemented basis function strings (see derived class implementation of method basisCode()
150/// for those strings
151
153{
154 // Check that primary variable of basis functions is our convolution variable
155 if (inBasis->getParameter(0) != x.absArg()) {
156 coutE(InputArguments) << "RooResolutionModel::convolution(" << GetName() << "," << this
157 << ") convolution parameter of basis function and PDF don't match" << endl
158 << "basis->findServer(0) = " << inBasis->findServer(0) << endl
159 << "x.absArg() = " << x.absArg() << endl ;
160 return 0 ;
161 }
162
163 if (basisCode(inBasis->GetTitle())==0) {
164 coutE(InputArguments) << "RooResolutionModel::convolution(" << GetName() << "," << this
165 << ") basis function '" << inBasis->GetTitle() << "' is not supported." << endl ;
166 return 0 ;
167 }
168
169 TString newName(GetName()) ;
170 newName.Append("_conv_") ;
171 newName.Append(inBasis->GetName()) ;
172 newName.Append("_[") ;
173 newName.Append(owner->GetName()) ;
174 newName.Append("]") ;
175
176 RooResolutionModel* conv = (RooResolutionModel*) clone(newName) ;
177
178 TString newTitle(conv->GetTitle()) ;
179 newTitle.Append(" convoluted with basis function ") ;
180 newTitle.Append(inBasis->GetName()) ;
181 conv->SetTitle(newTitle.Data()) ;
182
183 conv->changeBasis(inBasis) ;
184
185 return conv ;
186}
187
188
189
190////////////////////////////////////////////////////////////////////////////////
191/// Change the basis function we convolute with.
192/// For one-time use by convolution() only.
193
195{
196 // Remove client-server link to old basis
197 if (_basis) {
198 TIterator* bsIter = _basis->serverIterator() ;
199 RooAbsArg* basisServer ;
200 while((basisServer = (RooAbsArg*)bsIter->Next())) {
201 removeServer(*basisServer) ;
202 }
203 delete bsIter ;
204
205 if (_ownBasis) {
206 delete _basis ;
207 }
208 }
209 _ownBasis = false ;
210
211 // Change basis pointer and update client-server link
212 _basis = inBasis ;
213 if (_basis) {
214 TIterator* bsIter = _basis->serverIterator() ;
215 RooAbsArg* basisServer ;
216 while((basisServer = (RooAbsArg*)bsIter->Next())) {
217 addServer(*basisServer,true,false) ;
218 }
219 delete bsIter ;
220 }
221
222 _basisCode = inBasis?basisCode(inBasis->GetTitle()):0 ;
223}
224
225
226
227////////////////////////////////////////////////////////////////////////////////
228/// Return the convolution variable of the selection basis function.
229/// This is, by definition, the first parameter of the basis function
230
232{
233 // Convolution variable is by definition first server of basis function
234 TIterator* sIter = basis().serverIterator() ;
235 RooRealVar* var = (RooRealVar*) sIter->Next() ;
236 delete sIter ;
237
238 return *var ;
239}
240
241
242////////////////////////////////////////////////////////////////////////////////
243/// Modified version of RooAbsPdf::getValF(). If used as regular PDF,
244/// call RooAbsPdf::getValF(), otherwise return unnormalized value
245/// regardless of specified normalization set
246
247double RooResolutionModel::getValV(const RooArgSet* nset) const
248{
249 if (!_basis) return RooAbsPdf::getValV(nset) ;
250
251 // Return value of object. Calculated if dirty, otherwise cached value is returned.
252 if (isValueDirty()) {
253 _value = evaluate() ;
254
255 // WVE insert traceEval traceEval
256 if (_verboseDirty) cxcoutD(Tracing) << "RooResolutionModel(" << GetName() << ") value = " << _value << endl ;
257
260 }
261
262 return _value ;
263}
264
265
266
267////////////////////////////////////////////////////////////////////////////////
268/// Forward redirectServers call to our basis function, which is not connected to either resolution
269/// model or the physics model.
270
271bool RooResolutionModel::redirectServersHook(const RooAbsCollection& newServerList, bool mustReplaceAll, bool nameChange, bool /*isRecursive*/)
272{
273 if (!_basis) {
274 _norm = 0 ;
275 return false ;
276 }
277
278 RooFormulaVar* newBasis = (RooFormulaVar*) newServerList.find(_basis->GetName()) ;
279 if (newBasis) {
280
281 if (_ownBasis) {
282 delete _basis ;
283 }
284
285 _basis = newBasis ;
286 _ownBasis = false ;
287 }
288
289 _basis->redirectServers(newServerList,mustReplaceAll,nameChange) ;
290
291 return (mustReplaceAll && !newBasis) ;
292}
293
294
295
296////////////////////////////////////////////////////////////////////////////////
297/// Floating point error checking and tracing for given float value
298
299//bool RooResolutionModel::traceEvalHook(double value) const
300//{
301// // check for a math error or negative value
302// return TMath::IsNaN(value) ;
303//}
304
305
306
307////////////////////////////////////////////////////////////////////////////////
308/// Return the list of servers used by our normalization integral
309
311{
312 _norm->leafNodeServerList(&list) ;
313}
314
315
316
317////////////////////////////////////////////////////////////////////////////////
318/// Return the integral of this PDF over all elements of 'nset'.
319
320double RooResolutionModel::getNorm(const RooArgSet* nset) const
321{
322 if (!nset) {
323 return getVal() ;
324 }
325
326 syncNormalization(nset,false) ;
327 if (_verboseEval>1) cxcoutD(Tracing) << IsA()->GetName() << "::getNorm(" << GetName()
328 << "): norm(" << _norm << ") = " << _norm->getVal() << endl ;
329
330 double ret = _norm->getVal() ;
331 return ret ;
332}
333
334
335
336////////////////////////////////////////////////////////////////////////////////
337/// Print info about this object to the specified stream. In addition to the info
338/// from RooAbsArg::printStream() we add:
339///
340/// Shape : value, units, plot range
341/// Verbose : default binning and print label
342
343void RooResolutionModel::printMultiline(ostream& os, Int_t content, bool verbose, TString indent) const
344{
346
347 if(verbose) {
348 os << indent << "--- RooResolutionModel ---" << endl;
349 os << indent << "basis function = " ;
350 if (_basis) {
352 } else {
353 os << "<none>" << endl ;
354 }
355 }
356}
357
#define cxcoutD(a)
Definition: RooMsgService.h:85
#define coutE(a)
Definition: RooMsgService.h:37
#define ClassImp(name)
Definition: Rtypes.h:375
static void indent(ostringstream &buf, int indent_level)
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:77
void leafNodeServerList(RooAbsCollection *list, const RooAbsArg *arg=0, bool recurseNonDerived=false) const
Fill supplied list with all leaf nodes of the arg tree, starting with ourself as top node.
Definition: RooAbsArg.cxx:500
TObject * Clone(const char *newname=0) const override
Make a clone of an object using the Streamer facility.
Definition: RooAbsArg.h:89
void removeServer(RooAbsArg &server, bool force=false)
Unregister another RooAbsArg as a server to us, ie, declare that we no longer depend on its value and...
Definition: RooAbsArg.cxx:430
friend class RooArgSet
Definition: RooAbsArg.h:648
static bool _verboseDirty
cache of the list of proxies. Avoids type casting.
Definition: RooAbsArg.h:719
bool redirectServers(const RooAbsCollection &newServerList, bool mustReplaceAll=false, bool nameChange=false, bool isRecursionStep=false)
Replace all direct servers of this object with the new servers in newServerList.
Definition: RooAbsArg.cxx:1031
void clearValueDirty() const
Definition: RooAbsArg.h:622
void addServer(RooAbsArg &server, bool valueProp=true, bool shapeProp=false, std::size_t refCount=1)
Register another RooAbsArg as a server to us, ie, declare that we depend on it.
Definition: RooAbsArg.cxx:379
bool isValueDirty() const
Definition: RooAbsArg.h:444
void clearShapeDirty() const
Definition: RooAbsArg.h:625
TIterator * serverIterator() const
Definition: RooAbsArg.h:145
RooAbsArg * findServer(const char *name) const
Return server of this with name name. Returns nullptr if not found.
Definition: RooAbsArg.h:208
RooAbsCollection is an abstract container object that can hold multiple RooAbsArg objects.
RooAbsArg * find(const char *name) const
Find object with given name in list.
virtual bool syncNormalization(const RooArgSet *dset, bool adjustProxies=true) const
Verify that the normalization integral cached with this PDF is valid for given set of normalization o...
Definition: RooAbsPdf.cxx:511
void printMultiline(std::ostream &os, Int_t contents, bool verbose=false, TString indent="") const override
Print multi line detailed information of this RooAbsPdf.
Definition: RooAbsPdf.cxx:1876
RooAbsReal * _norm
Definition: RooAbsPdf.h:353
double getValV(const RooArgSet *set=0) const override
Return current value, normalized by integrating over the observables in nset.
Definition: RooAbsPdf.cxx:326
static Int_t _verboseEval
Definition: RooAbsPdf.h:347
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:94
double _value
Cache for current value of object.
Definition: RooAbsReal.h:484
virtual double evaluate() const =0
Evaluate this PDF / function / constant. Needs to be overridden by all derived classes.
RooAbsArg * absArg() const
Return pointer to contained argument.
Definition: RooArgProxy.h:40
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:57
A RooFormulaVar is a generic implementation of a real-valued object, which takes a RooArgList of serv...
Definition: RooFormulaVar.h:30
RooAbsArg * getParameter(const char *name) const
Return pointer to parameter with given name.
Definition: RooFormulaVar.h:44
virtual void printStream(std::ostream &os, Int_t contents, StyleOption style, TString indent="") const
Print description of object on ostream, printing contents set by contents integer,...
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:40
RooResolutionModel is the base class for PDFs that represent a resolution model that can be convolute...
~RooResolutionModel() override
Destructor.
TObject * clone(const char *newname) const override=0
virtual void changeBasis(RooFormulaVar *basis)
Change the basis function we convolute with.
virtual Int_t basisCode(const char *name) const =0
virtual RooResolutionModel * convolution(RooFormulaVar *basis, RooAbsArg *owner) const
Instantiate a clone of this resolution model representing a convolution with given basis function.
double getValV(const RooArgSet *nset=0) const override
Modified version of RooAbsPdf::getValF().
static RooFormulaVar * identity()
Return identity formula pointer.
bool redirectServersHook(const RooAbsCollection &newServerList, bool mustReplaceAll, bool nameChange, bool isRecursive) override
Forward redirectServers call to our basis function, which is not connected to either resolution model...
const RooRealVar & basisConvVar() const
Return the convolution variable of the selection basis function.
bool _ownBasis
Flag indicating ownership of _basis.
double getNorm(const RooArgSet *nset=0) const override
Return the integral of this PDF over all elements of 'nset'.
TClass * IsA() const override
virtual void normLeafServerList(RooArgSet &list) const
Floating point error checking and tracing for given float value.
void printMultiline(std::ostream &os, Int_t content, bool verbose=false, TString indent="") const override
Print info about this object to the specified stream.
Int_t _basisCode
Identifier code for selected basis function.
RooFormulaVar * _basis
Basis function convolved with this resolution model.
const RooFormulaVar & basis() const
RooTemplateProxy< RooAbsRealLValue > x
Dependent/convolution variable.
Iterator abstract base class.
Definition: TIterator.h:30
virtual TObject * Next()=0
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
Definition: TNamed.cxx:164
const char * GetName() const override
Returns name of object.
Definition: TNamed.h:47
const char * GetTitle() const override
Returns title of object.
Definition: TNamed.h:48
Basic string class.
Definition: TString.h:136
const char * Data() const
Definition: TString.h:369
TString & Append(const char *cs)
Definition: TString.h:564
Double_t x[n]
Definition: legend1.C:17
@ InputArguments
Definition: RooGlobalFunc.h:64