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_t 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 "RooFit.h"
66 
67 #include "TClass.h"
68 #include "TMath.h"
69 #include "Riostream.h"
70 #include "RooResolutionModel.h"
71 #include "RooMsgService.h"
72 
73 using namespace std;
74 
76 
77 
78 ////////////////////////////////////////////////////////////////////////////////
79 /// Constructor with convolution variable 'x'.
80 /// The convolution variable needs to be convertable to real values, and be able
81 /// to give information about its range. This is supported by e.g. RooRealVar or RooLinearVar, which
82 /// accepts offsetting and scaling an observable.
83 RooResolutionModel::RooResolutionModel(const char *name, const char *title, RooAbsRealLValue& _x) :
84  RooAbsPdf(name,title),
85  x("x","Dependent or convolution variable",this,_x),
86  _basisCode(0), _basis(0),
87  _ownBasis(kFALSE)
88 {
89 
90 }
91 
92 
93 
94 ////////////////////////////////////////////////////////////////////////////////
95 /// Copy constructor
96 
98  RooAbsPdf(other,name),
99  x("x",this,other.x),
100  _basisCode(other._basisCode), _basis(0),
101  _ownBasis(kFALSE)
102 {
103  if (other._basis) {
104  _basis = (RooFormulaVar*) other._basis->Clone() ;
105  _ownBasis = kTRUE ;
106  //_basis = other._basis ;
107  }
108 
109  if (_basis) {
110  TIterator* bsIter = _basis->serverIterator() ;
111  RooAbsArg* basisServer ;
112  while((basisServer = (RooAbsArg*)bsIter->Next())) {
113  addServer(*basisServer,kTRUE,kFALSE) ;
114  }
115  delete bsIter ;
116  }
117 }
118 
119 
120 
121 ////////////////////////////////////////////////////////////////////////////////
122 /// Destructor
123 
125 {
126  if (_ownBasis && _basis) {
127  delete _basis ;
128  }
129 }
130 
131 
132 
133 ////////////////////////////////////////////////////////////////////////////////
134 /// Return identity formula pointer
135 
137 {
138  static RooFormulaVar identity("identity","1",RooArgSet(""));
139  return &identity;
140 }
141 
142 
143 
144 ////////////////////////////////////////////////////////////////////////////////
145 /// Instantiate a clone of this resolution model representing a convolution with given
146 /// basis function. The owners object name is incorporated in the clones name
147 /// to avoid multiple convolution objects with the same name in complex PDF structures.
148 ///
149 /// Note: The 'inBasis' formula expression must be a RooFormulaVar that encodes the formula
150 /// in the title of the object and this expression must be an exact match against the
151 /// implemented basis function strings (see derived class implementation of method basisCode()
152 /// for those strings
153 
155 {
156  // Check that primary variable of basis functions is our convolution variable
157  if (inBasis->getParameter(0) != x.absArg()) {
158  coutE(InputArguments) << "RooResolutionModel::convolution(" << GetName() << "," << this
159  << ") convolution parameter of basis function and PDF don't match" << endl
160  << "basis->findServer(0) = " << inBasis->findServer(0) << endl
161  << "x.absArg() = " << x.absArg() << endl ;
162  return 0 ;
163  }
164 
165  if (basisCode(inBasis->GetTitle())==0) {
166  coutE(InputArguments) << "RooResolutionModel::convolution(" << GetName() << "," << this
167  << ") basis function '" << inBasis->GetTitle() << "' is not supported." << endl ;
168  return 0 ;
169  }
170 
171  TString newName(GetName()) ;
172  newName.Append("_conv_") ;
173  newName.Append(inBasis->GetName()) ;
174  newName.Append("_[") ;
175  newName.Append(owner->GetName()) ;
176  newName.Append("]") ;
177 
178  RooResolutionModel* conv = (RooResolutionModel*) clone(newName) ;
179 
180  TString newTitle(conv->GetTitle()) ;
181  newTitle.Append(" convoluted with basis function ") ;
182  newTitle.Append(inBasis->GetName()) ;
183  conv->SetTitle(newTitle.Data()) ;
184 
185  conv->changeBasis(inBasis) ;
186 
187  return conv ;
188 }
189 
190 
191 
192 ////////////////////////////////////////////////////////////////////////////////
193 /// Change the basis function we convolute with.
194 /// For one-time use by convolution() only.
195 
197 {
198  // Remove client-server link to old basis
199  if (_basis) {
200  TIterator* bsIter = _basis->serverIterator() ;
201  RooAbsArg* basisServer ;
202  while((basisServer = (RooAbsArg*)bsIter->Next())) {
203  removeServer(*basisServer) ;
204  }
205  delete bsIter ;
206 
207  if (_ownBasis) {
208  delete _basis ;
209  }
210  }
211  _ownBasis = kFALSE ;
212 
213  // Change basis pointer and update client-server link
214  _basis = inBasis ;
215  if (_basis) {
216  TIterator* bsIter = _basis->serverIterator() ;
217  RooAbsArg* basisServer ;
218  while((basisServer = (RooAbsArg*)bsIter->Next())) {
219  addServer(*basisServer,kTRUE,kFALSE) ;
220  }
221  delete bsIter ;
222  }
223 
224  _basisCode = inBasis?basisCode(inBasis->GetTitle()):0 ;
225 }
226 
227 
228 
229 ////////////////////////////////////////////////////////////////////////////////
230 /// Return the convolution variable of the selection basis function.
231 /// This is, by definition, the first parameter of the basis function
232 
234 {
235  // Convolution variable is by definition first server of basis function
236  TIterator* sIter = basis().serverIterator() ;
237  RooRealVar* var = (RooRealVar*) sIter->Next() ;
238  delete sIter ;
239 
240  return *var ;
241 }
242 
243 
244 ////////////////////////////////////////////////////////////////////////////////
245 /// Modified version of RooAbsPdf::getValF(). If used as regular PDF,
246 /// call RooAbsPdf::getValF(), otherwise return unnormalized value
247 /// regardless of specified normalization set
248 
250 {
251  if (!_basis) return RooAbsPdf::getValV(nset) ;
252 
253  // Return value of object. Calculated if dirty, otherwise cached value is returned.
254  if (isValueDirty()) {
255  _value = evaluate() ;
256 
257  // WVE insert traceEval traceEval
258  if (_verboseDirty) cxcoutD(Tracing) << "RooResolutionModel(" << GetName() << ") value = " << _value << endl ;
259 
260  clearValueDirty() ;
261  clearShapeDirty() ;
262  }
263 
264  return _value ;
265 }
266 
267 
268 
269 ////////////////////////////////////////////////////////////////////////////////
270 /// Forward redirectServers call to our basis function, which is not connected to either resolution
271 /// model or the physics model.
272 
273 Bool_t RooResolutionModel::redirectServersHook(const RooAbsCollection& newServerList, Bool_t mustReplaceAll, Bool_t nameChange, Bool_t /*isRecursive*/)
274 {
275  if (!_basis) {
276  _norm = 0 ;
277  return kFALSE ;
278  }
279 
280  RooFormulaVar* newBasis = (RooFormulaVar*) newServerList.find(_basis->GetName()) ;
281  if (newBasis) {
282 
283  if (_ownBasis) {
284  delete _basis ;
285  }
286 
287  _basis = newBasis ;
288  _ownBasis = kFALSE ;
289  }
290 
291  _basis->redirectServers(newServerList,mustReplaceAll,nameChange) ;
292 
293  return (mustReplaceAll && !newBasis) ;
294 }
295 
296 
297 
298 ////////////////////////////////////////////////////////////////////////////////
299 /// Floating point error checking and tracing for given float value
300 
301 //Bool_t RooResolutionModel::traceEvalHook(Double_t value) const
302 //{
303 // // check for a math error or negative value
304 // return TMath::IsNaN(value) ;
305 //}
306 
307 
308 
309 ////////////////////////////////////////////////////////////////////////////////
310 /// Return the list of servers used by our normalization integral
311 
313 {
314  _norm->leafNodeServerList(&list) ;
315 }
316 
317 
318 
319 ////////////////////////////////////////////////////////////////////////////////
320 /// Return the integral of this PDF over all elements of 'nset'.
321 
323 {
324  if (!nset) {
325  return getVal() ;
326  }
327 
328  syncNormalization(nset,kFALSE) ;
329  if (_verboseEval>1) cxcoutD(Tracing) << IsA()->GetName() << "::getNorm(" << GetName()
330  << "): norm(" << _norm << ") = " << _norm->getVal() << endl ;
331 
332  Double_t ret = _norm->getVal() ;
333  return ret ;
334 }
335 
336 
337 
338 ////////////////////////////////////////////////////////////////////////////////
339 /// Print info about this object to the specified stream. In addition to the info
340 /// from RooAbsArg::printStream() we add:
341 ///
342 /// Shape : value, units, plot range
343 /// Verbose : default binning and print label
344 
346 {
348 
349  if(verbose) {
350  os << indent << "--- RooResolutionModel ---" << endl;
351  os << indent << "basis function = " ;
352  if (_basis) {
354  } else {
355  os << "<none>" << endl ;
356  }
357  }
358 }
359 
RooAbsArg::Clone
virtual TObject * Clone(const char *newname=0) const
Make a clone of an object using the Streamer facility.
Definition: RooAbsArg.h:85
RooAbsArg::serverIterator
TIterator * serverIterator() const
Definition: RooAbsArg.h:141
kTRUE
const Bool_t kTRUE
Definition: RtypesCore.h:91
RooMsgService.h
RooResolutionModel::_ownBasis
Bool_t _ownBasis
Definition: RooResolutionModel.h:75
RooFit.h
RooFit::InputArguments
@ InputArguments
Definition: RooGlobalFunc.h:68
RooResolutionModel::_basis
RooFormulaVar * _basis
Definition: RooResolutionModel.h:74
RooResolutionModel::getValV
Double_t getValV(const RooArgSet *nset=0) const
Modified version of RooAbsPdf::getValF().
Definition: RooResolutionModel.cxx:249
TString::Data
const char * Data() const
Definition: TString.h:369
ClassImp
#define ClassImp(name)
Definition: Rtypes.h:364
TNamed::GetTitle
virtual const char * GetTitle() const
Returns title of object.
Definition: TNamed.h:54
coutE
#define coutE(a)
Definition: RooMsgService.h:33
RooAbsReal::getVal
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:91
RooResolutionModel.h
RooAbsPdf::getValV
virtual Double_t getValV(const RooArgSet *set=0) const
Return current value, normalized by integrating over the observables in nset.
Definition: RooAbsPdf.cxx:278
RooAbsPdf::printMultiline
virtual void printMultiline(std::ostream &os, Int_t contents, Bool_t verbose=kFALSE, TString indent="") const
Print multi line detailed information of this RooAbsPdf.
Definition: RooAbsPdf.cxx:1978
RooAbsCollection::find
RooAbsArg * find(const char *name) const
Find object with given name in list.
Definition: RooAbsCollection.cxx:813
x
Double_t x[n]
Definition: legend1.C:17
TClass.h
indent
static void indent(ostringstream &buf, int indent_level)
Definition: TClingCallFunc.cxx:87
RooPrintable::kAddress
@ kAddress
Definition: RooPrintable.h:33
RooResolutionModel::basisConvVar
const RooRealVar & basisConvVar() const
Return the convolution variable of the selection basis function.
Definition: RooResolutionModel.cxx:233
TString
Definition: TString.h:136
bool
TIterator
Definition: TIterator.h:30
RooResolutionModel::normLeafServerList
virtual void normLeafServerList(RooArgSet &list) const
Floating point error checking and tracing for given float value.
Definition: RooResolutionModel.cxx:312
RooResolutionModel::redirectServersHook
virtual Bool_t redirectServersHook(const RooAbsCollection &newServerList, Bool_t mustReplaceAll, Bool_t nameChange, Bool_t isRecursive)
Forward redirectServers call to our basis function, which is not connected to either resolution model...
Definition: RooResolutionModel.cxx:273
RooResolutionModel::convolution
virtual RooResolutionModel * convolution(RooFormulaVar *basis, RooAbsArg *owner) const
Instantiate a clone of this resolution model representing a convolution with given basis function.
Definition: RooResolutionModel.cxx:154
RooResolutionModel::printMultiline
virtual void printMultiline(std::ostream &os, Int_t content, Bool_t verbose=kFALSE, TString indent="") const
Print info about this object to the specified stream.
Definition: RooResolutionModel.cxx:345
RooResolutionModel::basisCode
virtual Int_t basisCode(const char *name) const =0
RooPrintable::kName
@ kName
Definition: RooPrintable.h:33
RooResolutionModel::~RooResolutionModel
virtual ~RooResolutionModel()
Destructor.
Definition: RooResolutionModel.cxx:124
RooFormulaVar
Definition: RooFormulaVar.h:30
RooResolutionModel::x
RooTemplateProxy< RooAbsRealLValue > x
Definition: RooResolutionModel.h:65
RooResolutionModel
Definition: RooResolutionModel.h:26
RooFormulaVar::getParameter
RooAbsArg * getParameter(const char *name) const
Definition: RooFormulaVar.h:41
RooAbsArg::removeServer
void removeServer(RooAbsArg &server, Bool_t force=kFALSE)
Unregister another RooAbsArg as a server to us, ie, declare that we no longer depend on its value and...
Definition: RooAbsArg.cxx:407
RooAbsArg::isValueDirty
Bool_t isValueDirty() const
Definition: RooAbsArg.h:429
RooArgProxy::absArg
RooAbsArg * absArg() const
Definition: RooArgProxy.h:51
kFALSE
const Bool_t kFALSE
Definition: RtypesCore.h:92
TString::Append
TString & Append(const char *cs)
Definition: TString.h:564
RooAbsReal::evaluate
virtual Double_t evaluate() const =0
Evaluate this PDF / function / constant. Needs to be overridden by all derived classes.
RooAbsCollection
Definition: RooAbsCollection.h:30
RooAbsArg::RooArgSet
friend class RooArgSet
Definition: RooAbsArg.h:579
RooAbsReal::_value
Double_t _value
Definition: RooAbsReal.h:450
TNamed::SetTitle
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
Definition: TNamed.cxx:164
TIterator::Next
virtual TObject * Next()=0
RooResolutionModel::RooResolutionModel
RooResolutionModel()
Definition: RooResolutionModel.h:30
RooAbsArg::addServer
void addServer(RooAbsArg &server, Bool_t valueProp=kTRUE, Bool_t shapeProp=kFALSE, std::size_t refCount=1)
Register another RooAbsArg as a server to us, ie, declare that we depend on it.
Definition: RooAbsArg.cxx:356
RooAbsArg::clearShapeDirty
void clearShapeDirty() const
Definition: RooAbsArg.h:557
RooFit::Tracing
@ Tracing
Definition: RooGlobalFunc.h:68
RooResolutionModel::_basisCode
Int_t _basisCode
Definition: RooResolutionModel.h:73
Double_t
double Double_t
Definition: RtypesCore.h:59
RooAbsArg::leafNodeServerList
void leafNodeServerList(RooAbsCollection *list, const RooAbsArg *arg=0, Bool_t recurseNonDerived=kFALSE) const
Fill supplied list with all leaf nodes of the arg tree, starting with ourself as top node.
Definition: RooAbsArg.cxx:477
RooResolutionModel::clone
virtual TObject * clone(const char *newname) const =0
RooAbsPdf::_norm
RooAbsReal * _norm
Definition: RooAbsPdf.h:320
RooResolutionModel::identity
static RooFormulaVar * identity()
Return identity formula pointer.
Definition: RooResolutionModel.cxx:136
RooResolutionModel::getNorm
Double_t getNorm(const RooArgSet *nset=0) const
Return the integral of this PDF over all elements of 'nset'.
Definition: RooResolutionModel.cxx:322
RooAbsPdf::_verboseEval
static Int_t _verboseEval
Definition: RooAbsPdf.h:314
name
char name[80]
Definition: TGX11.cxx:110
RooAbsArg::redirectServers
Bool_t redirectServers(const RooAbsCollection &newServerList, Bool_t mustReplaceAll=kFALSE, Bool_t nameChange=kFALSE, Bool_t isRecursionStep=kFALSE)
Substitute our servers with those listed in newSet.
Definition: RooAbsArg.cxx:911
RooPrintable::kTitle
@ kTitle
Definition: RooPrintable.h:33
RooPrintable::printStream
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,...
Definition: RooPrintable.cxx:75
genreflex::verbose
bool verbose
Definition: rootcling_impl.cxx:133
RooPrintable::kSingleLine
@ kSingleLine
Definition: RooPrintable.h:34
RooAbsArg
Definition: RooAbsArg.h:73
RooResolutionModel::changeBasis
virtual void changeBasis(RooFormulaVar *basis)
Change the basis function we convolute with.
Definition: RooResolutionModel.cxx:196
RooAbsPdf
Definition: RooAbsPdf.h:40
TNamed::GetName
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:53
RooAbsArg::findServer
RooAbsArg * findServer(const char *name) const
Return server of this with name name. Returns nullptr if not found.
Definition: RooAbsArg.h:203
RooAbsPdf::syncNormalization
virtual Bool_t syncNormalization(const RooArgSet *dset, Bool_t adjustProxies=kTRUE) const
Verify that the normalization integral cached with this PDF is valid for given set of normalization o...
Definition: RooAbsPdf.cxx:522
RooRealVar
Definition: RooRealVar.h:35
RooAbsArg::_verboseDirty
static Bool_t _verboseDirty
Definition: RooAbsArg.h:647
Riostream.h
RooAbsRealLValue
Definition: RooAbsRealLValue.h:31
RooAbsArg::clearValueDirty
void clearValueDirty() const
Definition: RooAbsArg.h:554
TMath.h
RooArgSet
Definition: RooArgSet.h:28
cxcoutD
#define cxcoutD(a)
Definition: RooMsgService.h:81
RooResolutionModel::basis
const RooFormulaVar & basis() const
Definition: RooResolutionModel.h:52
int