Logo ROOT   6.18/05
Reference Guide
RooProdGenContext.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 RooProdGenContext.cxx
19\class RooProdGenContext
20\ingroup Roofitcore
21
22RooProdGenContext is an efficient implementation of the generator context
23specific for RooProdPdf PDFs. The sim-context owns a list of
24component generator contexts that are used to generate the dependents
25for each component PDF sequentially.
26**/
27
28#include "RooFit.h"
29#include "Riostream.h"
30#include "RooMsgService.h"
31
32#include "RooProdGenContext.h"
33#include "RooProdGenContext.h"
34#include "RooProdPdf.h"
35#include "RooDataSet.h"
36#include "RooRealVar.h"
37#include "RooGlobalFunc.h"
38
39
40
41using namespace std;
42
44;
45
46
47////////////////////////////////////////////////////////////////////////////////
48
50 const RooDataSet *prototype, const RooArgSet* auxProto, Bool_t verbose) :
51 RooAbsGenContext(model,vars,prototype,auxProto,verbose), _uniIter(0), _pdf(&model)
52{
53 // Constructor of optimization generator context for RooProdPdf objects
54
55 //Build an array of generator contexts for each product component PDF
56 cxcoutI(Generation) << "RooProdGenContext::ctor() setting up event special generator context for product p.d.f. " << model.GetName()
57 << " for generation of observable(s) " << vars ;
58 if (prototype) ccxcoutI(Generation) << " with prototype data for " << *prototype->get() ;
59 if (auxProto && auxProto->getSize()>0) ccxcoutI(Generation) << " with auxiliary prototypes " << *auxProto ;
60 ccxcoutI(Generation) << endl ;
61
62 // Make full list of dependents (generated & proto)
63 RooArgSet deps(vars) ;
64 if (prototype) {
65 RooArgSet* protoDeps = model.getObservables(*prototype->get()) ;
66 deps.remove(*protoDeps,kTRUE,kTRUE) ;
67 delete protoDeps ;
68 }
69
70 // Factorize product in irreducible terms
71 RooLinkedList termList,depsList,impDepList,crossDepList,intList ;
72 model.factorizeProduct(deps,RooArgSet(),termList,depsList,impDepList,crossDepList,intList) ;
73 TIterator* termIter = termList.MakeIterator() ;
74 TIterator* normIter = depsList.MakeIterator() ;
75 TIterator* impIter = impDepList.MakeIterator() ;
76
77 if (dologD(Generation)) {
78 cxcoutD(Generation) << "RooProdGenContext::ctor() factorizing product expression in irriducible terms " ;
79 while(RooArgSet* t=(RooArgSet*)termIter->Next()) {
80 ccxcoutD(Generation) << *t ;
81 }
82 ccxcoutD(Generation) << endl ;
83 }
84
85 RooArgSet genDeps ;
86 // First add terms that do not import observables
87
88 Bool_t anyAction = kTRUE ;
89 Bool_t go=kTRUE ;
90 while(go) {
91
92 RooAbsPdf* pdf ;
93 RooArgSet* term ;
94 RooArgSet* impDeps ;
95 RooArgSet* termDeps ;
96
97 termIter->Reset() ;
98 impIter->Reset() ;
99 normIter->Reset() ;
100
101 Bool_t anyPrevAction=anyAction ;
102 anyAction=kFALSE ;
103
104 if (termList.GetSize()==0) {
105 break ;
106 }
107
108 while((term=(RooArgSet*)termIter->Next())) {
109
110 impDeps = (RooArgSet*)impIter->Next() ;
111 termDeps = (RooArgSet*)normIter->Next() ;
112 if (impDeps==0 || termDeps==0) {
113 break ;
114 }
115
116 cxcoutD(Generation) << "RooProdGenContext::ctor() analyzing product term " << *term << " with observable(s) " << *termDeps ;
117 if (impDeps->getSize()>0) {
118 ccxcoutD(Generation) << " which has dependence of external observable(s) " << *impDeps << " that to be generated first by other terms" ;
119 }
120 ccxcoutD(Generation) << endl ;
121
122 // Add this term if we have no imported dependents, or imported dependents are already generated
123 RooArgSet neededDeps(*impDeps) ;
124 neededDeps.remove(genDeps,kTRUE,kTRUE) ;
125
126 if (neededDeps.getSize()>0) {
127 if (!anyPrevAction) {
128 cxcoutD(Generation) << "RooProdGenContext::ctor() no convergence in single term analysis loop, terminating loop and process remainder of terms as single unit " << endl ;
129 go=kFALSE ;
130 break ;
131 }
132 cxcoutD(Generation) << "RooProdGenContext::ctor() skipping this term for now because it needs imported dependents that are not generated yet" << endl ;
133 continue ;
134 }
135
136 // Check if this component has any dependents that need to be generated
137 // e.g. it can happen that there are none if all dependents of this component are prototyped
138 if (termDeps->getSize()==0) {
139 cxcoutD(Generation) << "RooProdGenContext::ctor() term has no observables requested to be generated, removing it" << endl ;
140 termList.Remove(term) ;
141 depsList.Remove(termDeps) ;
142 impDepList.Remove(impDeps) ;
143 delete term ;
144 delete termDeps ;
145 delete impDeps ;
146 anyAction=kTRUE ;
147 continue ;
148 }
149
150 TIterator* pdfIter = term->createIterator() ;
151 if (term->getSize()==1) {
152 // Simple term
153
154 pdf = (RooAbsPdf*) pdfIter->Next() ;
155 RooArgSet* pdfDep = pdf->getObservables(termDeps) ;
156 if (pdfDep->getSize()>0) {
157 coutI(Generation) << "RooProdGenContext::ctor() creating subcontext for generation of observables " << *pdfDep << " from model " << pdf->GetName() << endl ;
158 RooArgSet* auxProto2 = pdf->getObservables(impDeps) ;
159 RooAbsGenContext* cx = pdf->genContext(*pdfDep,prototype,auxProto2,verbose) ;
160 delete auxProto2 ;
161 _gcList.push_back(cx) ;
162 }
163
164// cout << "adding following dependents to list of generated observables: " ; pdfDep->Print("1") ;
165 genDeps.add(*pdfDep) ;
166
167 delete pdfDep ;
168
169 } else {
170
171 // Composite term
172 if (termDeps->getSize()>0) {
173 const std::string name = model.makeRGPPName("PRODGEN_",*term,RooArgSet(),RooArgSet(),0) ;
174
175 // Construct auxiliary PDF expressing product of composite terms,
176 // following Conditional component specification of input model
177 RooLinkedList cmdList ;
178 RooLinkedList pdfSetList ;
179 pdfIter->Reset() ;
180 RooArgSet fullPdfSet ;
181 while((pdf=(RooAbsPdf*)pdfIter->Next())) {
182
183 RooArgSet* pdfnset = model.findPdfNSet(*pdf) ;
184 RooArgSet* pdfSet = new RooArgSet(*pdf) ;
185 pdfSetList.Add(pdfSet) ;
186
187 if (pdfnset && pdfnset->getSize()>0) {
188 // This PDF requires a Conditional() construction
189 cmdList.Add(RooFit::Conditional(*pdfSet,*pdfnset).Clone()) ;
190// cout << "Conditional " << pdf->GetName() << " " ; pdfnset->Print("1") ;
191 } else {
192 fullPdfSet.add(*pdfSet) ;
193 }
194
195 }
196 RooProdPdf* multiPdf = new RooProdPdf(name.c_str(),name.c_str(),fullPdfSet,cmdList) ;
197 cmdList.Delete() ;
198 pdfSetList.Delete() ;
199
201 multiPdf->useDefaultGen(kTRUE) ;
202 _ownedMultiProds.addOwned(*multiPdf) ;
203
204 coutI(Generation) << "RooProdGenContext()::ctor creating subcontext for generation of observables " << *termDeps
205 << "for irriducuble composite term using sub-product object " << multiPdf->GetName() ;
206 RooAbsGenContext* cx = multiPdf->genContext(*termDeps,prototype,auxProto,verbose) ;
207 _gcList.push_back(cx) ;
208
209 genDeps.add(*termDeps) ;
210
211 }
212 }
213
214 delete pdfIter ;
215
216// cout << "added generator for this term, removing from list" << endl ;
217
218 termList.Remove(term) ;
219 depsList.Remove(termDeps) ;
220 impDepList.Remove(impDeps) ;
221 delete term ;
222 delete termDeps ;
223 delete impDeps ;
224 anyAction=kTRUE ;
225 }
226 }
227
228 // Check if there are any left over terms that cannot be generated
229 // separately due to cross dependency of observables
230 if (termList.GetSize()>0) {
231
232 cxcoutD(Generation) << "RooProdGenContext::ctor() there are left-over terms that need to be generated separately" << endl ;
233
234 RooAbsPdf* pdf ;
235 RooArgSet* term ;
236
237 // Concatenate remaining terms
238 termIter->Reset() ;
239 normIter->Reset() ;
240 RooArgSet trailerTerm ;
241 RooArgSet trailerTermDeps ;
242 while((term=(RooArgSet*)termIter->Next())) {
243 RooArgSet* termDeps = (RooArgSet*)normIter->Next() ;
244 trailerTerm.add(*term) ;
245 trailerTermDeps.add(*termDeps) ;
246 }
247
248 const std::string name = model.makeRGPPName("PRODGEN_",trailerTerm,RooArgSet(),RooArgSet(),0) ;
249
250 // Construct auxiliary PDF expressing product of composite terms,
251 // following Partial/Full component specification of input model
252 RooLinkedList cmdList ;
253 RooLinkedList pdfSetList ;
254 RooArgSet fullPdfSet ;
255
256 TIterator* pdfIter = trailerTerm.createIterator() ;
257 while((pdf=(RooAbsPdf*)pdfIter->Next())) {
258
259 RooArgSet* pdfnset = model.findPdfNSet(*pdf) ;
260 RooArgSet* pdfSet = new RooArgSet(*pdf) ;
261 pdfSetList.Add(pdfSet) ;
262
263 if (pdfnset && pdfnset->getSize()>0) {
264 // This PDF requires a Conditional() construction
265 cmdList.Add(RooFit::Conditional(*pdfSet,*pdfnset).Clone()) ;
266 } else {
267 fullPdfSet.add(*pdfSet) ;
268 }
269
270 }
271// cmdList.Print("v") ;
272 RooProdPdf* multiPdf = new RooProdPdf(name.c_str(),name.c_str(),fullPdfSet,cmdList) ;
273 cmdList.Delete() ;
274 pdfSetList.Delete() ;
275
277 multiPdf->useDefaultGen(kTRUE) ;
278 _ownedMultiProds.addOwned(*multiPdf) ;
279
280 cxcoutD(Generation) << "RooProdGenContext(" << model.GetName() << "): creating context for irreducible composite trailer term "
281 << multiPdf->GetName() << " that generates observables " << trailerTermDeps << endl ;
282 RooAbsGenContext* cx = multiPdf->genContext(trailerTermDeps,prototype,auxProto,verbose) ;
283 _gcList.push_back(cx) ;
284 }
285
286 // Now check if the are observables in vars that are not generated by any of the above p.d.f.s
287 // If not, generate uniform distributions for these using a special context
288 _uniObs.add(vars) ;
289 _uniObs.remove(genDeps,kTRUE,kTRUE) ;
290 if (_uniObs.getSize()>0) {
292 coutI(Generation) << "RooProdGenContext(" << model.GetName() << "): generating uniform distribution for non-dependent observable(s) " << _uniObs << endl ;
293 }
294
295
296 delete termIter ;
297 delete impIter ;
298 delete normIter ;
299
300
301 // We own contents of lists filled by factorizeProduct()
302 termList.Delete() ;
303 depsList.Delete() ;
304 impDepList.Delete() ;
305 crossDepList.Delete() ;
306 intList.Delete() ;
307
308}
309
310
311
312////////////////////////////////////////////////////////////////////////////////
313/// Destructor. Delete all owned subgenerator contexts
314
316{
317 delete _uniIter ;
318 for (list<RooAbsGenContext*>::iterator iter=_gcList.begin() ; iter!=_gcList.end() ; ++iter) {
319 delete (*iter) ;
320 }
321}
322
323
324////////////////////////////////////////////////////////////////////////////////
325/// Attach generator to given event buffer
326
328{
329 //Forward initGenerator call to all components
330 for (list<RooAbsGenContext*>::iterator iter=_gcList.begin() ; iter!=_gcList.end() ; ++iter) {
331 (*iter)->attach(args) ;
332 }
333}
334
335
336////////////////////////////////////////////////////////////////////////////////
337/// One-time initialization of generator context, forward to component generators
338
340{
341 // Forward initGenerator call to all components
342 for (list<RooAbsGenContext*>::iterator iter=_gcList.begin() ; iter!=_gcList.end() ; ++iter) {
343 (*iter)->initGenerator(theEvent) ;
344 }
345}
346
347
348
349////////////////////////////////////////////////////////////////////////////////
350/// Generate a single event of the product by generating the components
351/// of the products sequentially. The subcontext have been order such
352/// that all conditional dependencies are correctly taken into account
353/// when processed in sequential order
354
356{
357 // Loop over the component generators
358
359 for (list<RooAbsGenContext*>::iterator iter=_gcList.begin() ; iter!=_gcList.end() ; ++iter) {
360 (*iter)->generateEvent(theEvent,remaining) ;
361 }
362
363 // Generate uniform variables (non-dependents)
364 if (_uniIter) {
365 _uniIter->Reset() ;
366 RooAbsArg* uniVar ;
367 while((uniVar=(RooAbsArg*)_uniIter->Next())) {
368 RooAbsLValue* arglv = dynamic_cast<RooAbsLValue*>(uniVar) ;
369 if (arglv) {
370 arglv->randomize() ;
371 }
372 }
373 theEvent = _uniObs ;
374 }
375
376}
377
378
379
380////////////////////////////////////////////////////////////////////////////////
381/// Set the traversal order of the prototype dataset by the
382/// given lookup table
383
385{
386 // Forward call to component generators
388
389 for (list<RooAbsGenContext*>::iterator iter=_gcList.begin() ; iter!=_gcList.end() ; ++iter) {
390 (*iter)->setProtoDataOrder(lut) ;
391 }
392
393}
394
395
396
397////////////////////////////////////////////////////////////////////////////////
398/// Detailed printing interface
399
401{
403 os << indent << "--- RooProdGenContext ---" << endl ;
404 os << indent << "Using PDF ";
406 os << indent << "List of component generators" << endl ;
407
408 TString indent2(indent) ;
409 indent2.Append(" ") ;
410
411 for (list<RooAbsGenContext*>::const_iterator iter=_gcList.begin() ; iter!=_gcList.end() ; ++iter) {
412 (*iter)->printMultiline(os,content,verbose,indent2) ;
413 }
414}
#define coutI(a)
Definition: RooMsgService.h:31
#define cxcoutI(a)
Definition: RooMsgService.h:83
#define cxcoutD(a)
Definition: RooMsgService.h:79
#define dologD(a)
Definition: RooMsgService.h:63
#define ccxcoutD(a)
Definition: RooMsgService.h:80
#define ccxcoutI(a)
Definition: RooMsgService.h:84
int Int_t
Definition: RtypesCore.h:41
const Bool_t kFALSE
Definition: RtypesCore.h:88
bool Bool_t
Definition: RtypesCore.h:59
const Bool_t kTRUE
Definition: RtypesCore.h:87
#define ClassImp(name)
Definition: Rtypes.h:365
static void indent(ostringstream &buf, int indent_level)
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:70
RooArgSet * getObservables(const RooArgSet &set, Bool_t valueOnly=kTRUE) const
Return the observables of this pdf given a set of observables.
Definition: RooAbsArg.h:240
void setOperMode(OperMode mode, Bool_t recurseADirty=kTRUE)
Change cache operation mode to given mode.
Definition: RooAbsArg.cxx:1719
Int_t getSize() const
TIterator * createIterator(Bool_t dir=kIterForward) const R__SUGGEST_ALTERNATIVE("begin()
TIterator-style iteration over contained elements.
virtual Bool_t remove(const RooAbsArg &var, Bool_t silent=kFALSE, Bool_t matchByNameOnly=kFALSE)
Remove the specified argument from our list.
RooAbsGenContext is the abstract base class for generator contexts of RooAbsPdf objects.
virtual void printMultiline(std::ostream &os, Int_t contents, Bool_t verbose=kFALSE, TString indent="") const
Interface for multi-line printing.
virtual void setProtoDataOrder(Int_t *lut)
Set the traversal order of prototype data to that in the lookup tables passed as argument.
Abstract base class for objects that are lvalues, i.e.
Definition: RooAbsLValue.h:26
virtual void randomize(const char *rangeName=0)=0
virtual RooAbsGenContext * genContext(const RooArgSet &vars, const RooDataSet *prototype=0, const RooArgSet *auxProto=0, Bool_t verbose=kFALSE) const
Interface function to create a generator context from a p.d.f.
Definition: RooAbsPdf.cxx:1709
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:28
virtual Bool_t addOwned(const RooAbsCollection &col, Bool_t silent=kFALSE)
Add a collection of arguments to this collection by calling addOwned() for each element in the source...
Definition: RooArgSet.h:92
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
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:31
virtual const RooArgSet * get(Int_t index) const
Return RooArgSet with coordinates of event 'index'.
RooLinkedList is an collection class for internal use, storing a collection of RooAbsArg pointers in ...
Definition: RooLinkedList.h:36
Int_t GetSize() const
Definition: RooLinkedList.h:61
TIterator * MakeIterator(Bool_t forward=kTRUE) const
Create a TIterator for this list.
void Delete(Option_t *o=0)
Remove all elements in collection and delete all elements NB: Collection does not own elements,...
virtual void Add(TObject *arg)
Definition: RooLinkedList.h:63
virtual Bool_t Remove(TObject *arg)
Remove object from collection.
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,...
RooProdGenContext is an efficient implementation of the generator context specific for RooProdPdf PDF...
RooArgSet _ownedMultiProds
std::list< RooAbsGenContext * > _gcList
virtual void attach(const RooArgSet &params)
Attach generator to given event buffer.
const RooProdPdf * _pdf
RooProdGenContext(const RooProdPdf &model, const RooArgSet &vars, const RooDataSet *prototype=0, const RooArgSet *auxProto=0, Bool_t _verbose=kFALSE)
virtual void generateEvent(RooArgSet &theEvent, Int_t remaining)
Generate a single event of the product by generating the components of the products sequentially.
virtual void setProtoDataOrder(Int_t *lut)
Set the traversal order of the prototype dataset by the given lookup table.
virtual void initGenerator(const RooArgSet &theEvent)
One-time initialization of generator context, forward to component generators.
virtual void printMultiline(std::ostream &os, Int_t content, Bool_t verbose=kFALSE, TString indent="") const
Detailed printing interface.
virtual ~RooProdGenContext()
Destructor. Delete all owned subgenerator contexts.
RooProdPdf is an efficient implementation of a product of PDFs of the form.
Definition: RooProdPdf.h:31
void factorizeProduct(const RooArgSet &normSet, const RooArgSet &intSet, RooLinkedList &termList, RooLinkedList &normList, RooLinkedList &impDepList, RooLinkedList &crossDepList, RooLinkedList &intList) const
Factorize product in irreducible terms for given choice of integration/normalization.
Definition: RooProdPdf.cxx:542
void useDefaultGen(Bool_t flag=kTRUE)
Definition: RooProdPdf.h:170
std::string makeRGPPName(const char *pfx, const RooArgSet &term, const RooArgSet &iset, const RooArgSet &nset, const char *isetRangeName) const
Make an appropriate automatic name for a RooGenProdProj object in getPartIntList()
RooArgSet * findPdfNSet(RooAbsPdf &pdf) const
Look up user specified normalization set for given input PDF component.
virtual RooAbsGenContext * genContext(const RooArgSet &vars, const RooDataSet *prototype=0, const RooArgSet *auxProto=0, Bool_t verbose=kFALSE) const
Return generator context optimized for generating events from product p.d.f.s.
Iterator abstract base class.
Definition: TIterator.h:30
virtual void Reset()=0
virtual TObject * Next()=0
virtual TObject * Clone(const char *newname="") const
Make a clone of an object using the Streamer facility.
Definition: TNamed.cxx:74
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
@ Generation
Definition: RooGlobalFunc.h:57
RooCmdArg Conditional(const RooArgSet &pdfSet, const RooArgSet &depSet, Bool_t depsAreCond=kFALSE)