ROOT  6.06/09
Reference Guide
RooSimSplitGenContext.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 // BEGIN_HTML
20 // RooSimSplitGenContext is an efficient implementation of the generator context
21 // specific for RooSimultaneous PDFs when generating more than one of the
22 // component pdfs.
23 // END_HTML
24 //
25 
26 #include "RooFit.h"
27 #include "Riostream.h"
28 
29 #include "RooSimSplitGenContext.h"
30 #include "RooSimultaneous.h"
31 #include "RooRealProxy.h"
32 #include "RooDataSet.h"
33 #include "Roo1DTable.h"
34 #include "RooCategory.h"
35 #include "RooMsgService.h"
36 #include "RooRandom.h"
37 #include "RooGlobalFunc.h"
38 
39 using namespace RooFit ;
40 
41 #include <string>
42 
43 using namespace std;
44 
46 ;
47 
48 
49 ////////////////////////////////////////////////////////////////////////////////
50 /// Constructor of specialized generator context for RooSimultaneous p.d.f.s. This
51 /// context creates a dedicated context for each component p.d.f.s and delegates
52 /// generation of events to the appropriate component generator context
53 
54 RooSimSplitGenContext::RooSimSplitGenContext(const RooSimultaneous &model, const RooArgSet &vars, Bool_t verbose, Bool_t autoBinned, const char* binnedTag) :
55  RooAbsGenContext(model,vars,0,0,verbose), _pdf(&model)
56 {
57  // Determine if we are requested to generate the index category
58  RooAbsCategory *idxCat = (RooAbsCategory*) model._indexCat.absArg() ;
59  RooArgSet pdfVars(vars) ;
60 
61  RooArgSet allPdfVars(pdfVars) ;
62 
63  if (!idxCat->isDerived()) {
64  pdfVars.remove(*idxCat,kTRUE,kTRUE) ;
65  Bool_t doGenIdx = allPdfVars.find(idxCat->GetName())?kTRUE:kFALSE ;
66 
67  if (!doGenIdx) {
68  oocoutE(_pdf,Generation) << "RooSimSplitGenContext::ctor(" << GetName() << ") ERROR: This context must"
69  << " generate the index category" << endl ;
70  _isValid = kFALSE ;
71  _numPdf = 0 ;
72  // coverity[UNINIT_CTOR]
73  return ;
74  }
75  } else {
76  TIterator* sIter = idxCat->serverIterator() ;
77  RooAbsArg* server ;
78  Bool_t anyServer(kFALSE), allServers(kTRUE) ;
79  while((server=(RooAbsArg*)sIter->Next())) {
80  if (vars.find(server->GetName())) {
81  anyServer=kTRUE ;
82  pdfVars.remove(*server,kTRUE,kTRUE) ;
83  } else {
84  allServers=kFALSE ;
85  }
86  }
87  delete sIter ;
88 
89  if (anyServer && !allServers) {
90  oocoutE(_pdf,Generation) << "RooSimSplitGenContext::ctor(" << GetName() << ") ERROR: This context must"
91  << " generate all components of a derived index category" << endl ;
92  _isValid = kFALSE ;
93  _numPdf = 0 ;
94  // coverity[UNINIT_CTOR]
95  return ;
96  }
97  }
98 
99  // We must extended likelihood to determine the relative fractions of the components
100  _idxCatName = idxCat->GetName() ;
101  if (!model.canBeExtended()) {
102  oocoutE(_pdf,Generation) << "RooSimSplitGenContext::ctor(" << GetName() << ") ERROR: Need either extended mode"
103  << " to calculate number of events per category" << endl ;
104  _isValid = kFALSE ;
105  _numPdf = 0 ;
106  // coverity[UNINIT_CTOR]
107  return ;
108  }
109 
110  // Initialize fraction threshold array (used only in extended mode)
111  _numPdf = model._pdfProxyList.GetSize() ;
112  _fracThresh = new Double_t[_numPdf+1] ;
113  _fracThresh[0] = 0 ;
114 
115  // Generate index category and all registered PDFS
117  _allVarsPdf.add(allPdfVars) ;
118  RooRealProxy* proxy ;
119  RooAbsPdf* pdf ;
120  Int_t i(1) ;
121  while((proxy=(RooRealProxy*)_proxyIter->Next())) {
122  pdf=(RooAbsPdf*)proxy->absArg() ;
123 
124  // Create generator context for this PDF
125  RooArgSet* compVars = pdf->getObservables(pdfVars) ;
126  RooAbsGenContext* cx = pdf->autoGenContext(*compVars,0,0,verbose,autoBinned,binnedTag) ;
127  delete compVars ;
128 
129  const RooCatType* state = idxCat->lookupType(proxy->name()) ;
130 
131  cx->SetName(proxy->name()) ;
132  _gcList.push_back(cx) ;
133  _gcIndex.push_back(state->getVal()) ;
134 
135  // Fill fraction threshold array
136  _fracThresh[i] = _fracThresh[i-1] + pdf->expectedEvents(&allPdfVars) ;
137  i++ ;
138  }
139 
140  for(i=0 ; i<_numPdf ; i++) {
142  }
143 
144  // Clone the index category
145  _idxCatSet = (RooArgSet*) RooArgSet(model._indexCat.arg()).snapshot(kTRUE) ;
146  if (!_idxCatSet) {
147  oocoutE(_pdf,Generation) << "RooSimSplitGenContext::RooSimSplitGenContext(" << GetName() << ") Couldn't deep-clone index category, abort," << endl ;
148  throw std::string("RooSimSplitGenContext::RooSimSplitGenContext() Couldn't deep-clone index category, abort") ;
149  }
150 
152 }
153 
154 
155 
156 ////////////////////////////////////////////////////////////////////////////////
157 /// Destructor. Delete all owned subgenerator contexts
158 
160 {
161  delete[] _fracThresh ;
162  delete _idxCatSet ;
163  for (vector<RooAbsGenContext*>::iterator iter = _gcList.begin() ; iter!=_gcList.end() ; ++iter) {
164  delete (*iter) ;
165  }
166  delete _proxyIter ;
167 }
168 
169 
170 
171 ////////////////////////////////////////////////////////////////////////////////
172 /// Attach the index category clone to the given event buffer
173 
175 {
176  if (_idxCat->isDerived()) {
178  }
179 
180  // Forward initGenerator call to all components
181  for (vector<RooAbsGenContext*>::iterator iter = _gcList.begin() ; iter!=_gcList.end() ; ++iter) {
182  (*iter)->attach(args) ;
183  }
184 
185 }
186 
187 
188 ////////////////////////////////////////////////////////////////////////////////
189 /// Perform one-time initialization of generator context
190 
192 {
193  // Attach the index category clone to the event
194  if (_idxCat->isDerived()) {
196  } else {
197  _idxCat = (RooAbsCategoryLValue*) theEvent.find(_idxCat->GetName()) ;
198  }
199 
200  // Forward initGenerator call to all components
201  for (vector<RooAbsGenContext*>::iterator iter = _gcList.begin() ; iter!=_gcList.end() ; ++iter) {
202  (*iter)->initGenerator(theEvent) ;
203  }
204 
205 }
206 
207 
208 
209 ////////////////////////////////////////////////////////////////////////////////
210 
212 {
213  if(!isValid()) {
214  coutE(Generation) << ClassName() << "::" << GetName() << ": context is not valid" << endl;
215  return 0;
216  }
217 
218 
219  // Calculate the expected number of events if necessary
220  if(nEvents <= 0) {
221  nEvents= _expectedEvents;
222  }
223  coutI(Generation) << ClassName() << "::" << GetName() << ":generate: will generate "
224  << nEvents << " events" << endl;
225 
226  if (_verbose) Print("v") ;
227 
228  // Perform any subclass implementation-specific initialization
229  // Can be skipped if this is a rerun with an identical configuration
230  if (!skipInit) {
232  }
233 
234  // Generate lookup table from expected event counts
235  vector<Double_t> nGen(_numPdf) ;
236  if (extendedMode ) {
237  _proxyIter->Reset() ;
238  RooRealProxy* proxy ;
239  Int_t i(0) ;
240  while((proxy=(RooRealProxy*)_proxyIter->Next())) {
241  RooAbsPdf* pdf=(RooAbsPdf*)proxy->absArg() ;
242  //nGen[i] = Int_t(pdf->expectedEvents(&_allVarsPdf)+0.5) ;
243  nGen[i] = pdf->expectedEvents(&_allVarsPdf) ;
244  i++ ;
245  }
246 
247  } else {
248  _proxyIter->Reset() ;
249  RooRealProxy* proxy ;
250  Int_t i(1) ;
251  _fracThresh[0] = 0 ;
252  while((proxy=(RooRealProxy*)_proxyIter->Next())) {
253  RooAbsPdf* pdf=(RooAbsPdf*)proxy->absArg() ;
255  i++ ;
256  }
257  for(i=0 ; i<_numPdf ; i++) {
259  }
260 
261  // Determine from that total number of events to be generated for each component
262  Double_t nGenSoFar(0) ;
263  while (nGenSoFar<nEvents) {
264  Double_t rand = RooRandom::uniform() ;
265  i=0 ;
266  for (i=0 ; i<_numPdf ; i++) {
267  if (rand>_fracThresh[i] && rand<_fracThresh[i+1]) {
268  nGen[i]++ ;
269  nGenSoFar++ ;
270  break ;
271  }
272  }
273  }
274  }
275 
276 
277 
278  // Now loop over states
279  _proxyIter->Reset() ;
280  map<string,RooAbsData*> dataMap ;
281  Int_t icomp(0) ;
282  RooRealProxy* proxy ;
283  while((proxy=(RooRealProxy*)_proxyIter->Next())) {
284 
285  // Calculate number of events to generate for this state
286  if (_gcList[icomp]) {
287  dataMap[proxy->GetName()] = _gcList[icomp]->generate(nGen[icomp],skipInit,extendedMode) ;
288  }
289 
290  icomp++ ;
291  }
292 
293  // Put all datasets together in a composite-store RooDataSet that links and owns the component datasets
294  RooDataSet* hmaster = new RooDataSet("hmaster","hmaster",_allVarsPdf,RooFit::Index((RooCategory&)*_idxCat),RooFit::Link(dataMap),RooFit::OwnLinked()) ;
295  return hmaster ;
296 }
297 
298 
299 
300 ////////////////////////////////////////////////////////////////////////////////
301 /// Forward to components
302 
304 {
305  for (vector<RooAbsGenContext*>::iterator iter=_gcList.begin() ; iter!=_gcList.end() ; ++iter) {
306  (*iter)->setExpectedData(flag) ;
307  }
308 }
309 
310 
311 
312 ////////////////////////////////////////////////////////////////////////////////
313 /// this method is empty because it is not used by this context
314 
315 RooDataSet* RooSimSplitGenContext::createDataSet(const char* , const char* , const RooArgSet& )
316 {
317  return 0 ;
318 }
319 
320 
321 
322 ////////////////////////////////////////////////////////////////////////////////
323 /// this method is empty because it is not used in this type of context
324 
326 {
327  assert(0) ;
328 }
329 
330 
331 
332 
333 ////////////////////////////////////////////////////////////////////////////////
334 /// this method is empty because proto datasets are not supported by this context
335 
337 {
338  assert(0) ;
339 }
340 
341 
342 ////////////////////////////////////////////////////////////////////////////////
343 /// Detailed printing interface
344 
346 {
347  RooAbsGenContext::printMultiline(os,content,verbose,indent) ;
348  os << indent << "--- RooSimSplitGenContext ---" << endl ;
349  os << indent << "Using PDF ";
351 }
#define coutE(a)
Definition: RooMsgService.h:35
RooCategoryProxy _indexCat
virtual Bool_t isDerived() const
Definition: RooAbsArg.h:81
virtual void generateEvent(RooArgSet &theEvent, Int_t remaining)
this method is empty because it is not used in this type of context
TIterator * serverIterator() const
Definition: RooAbsArg.h:112
virtual void Reset()=0
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, which is interpreted as an OR of 'enum ContentsOptions' values and in the style given by 'enum StyleOption'.
#define coutI(a)
Definition: RooMsgService.h:32
virtual void setExpectedData(Bool_t)
Forward to components.
ClassImp(RooSimSplitGenContext)
#define assert(cond)
Definition: unittest.h:542
virtual void SetName(const char *name)
Change (i.e.
Definition: TNamed.cxx:128
RooArgSet * getObservables(const RooArgSet &set, Bool_t valueOnly=kTRUE) const
Definition: RooAbsArg.h:194
RooCmdArg Link(const char *state, RooAbsData &data)
virtual const char * name() const
Definition: RooArgProxy.h:42
Bool_t isValid() const
Basic string class.
Definition: TString.h:137
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
const Bool_t kFALSE
Definition: Rtypes.h:92
STL namespace.
virtual Double_t expectedEvents(const RooArgSet *nset) const
Return expected number of events from this p.d.f for use in extended likelihood calculations.
Definition: RooAbsPdf.cxx:2930
Iterator abstract base class.
Definition: TIterator.h:32
RooDataSet * createDataSet(const char *name, const char *title, const RooArgSet &obs)
this method is empty because it is not used by this context
#define oocoutE(o, a)
Definition: RooMsgService.h:48
virtual void attach(const RooArgSet &params)
Attach the index category clone to the given event buffer.
std::map< std::string, std::string >::const_iterator iter
Definition: TAlienJob.cxx:54
virtual void initGenerator(const RooArgSet &theEvent)
Perform one-time initialization of generator context.
virtual RooAbsGenContext * autoGenContext(const RooArgSet &vars, const RooDataSet *prototype=0, const RooArgSet *auxProto=0, Bool_t verbose=kFALSE, Bool_t autoBinned=kTRUE, const char *binnedTag="") const
Definition: RooAbsPdf.cxx:1647
RooCmdArg OwnLinked()
const int nEvents
Definition: testRooFit.cxx:42
RooAbsArg * find(const char *name) const
Find object with given name in list.
return
Definition: TBase64.cxx:62
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition: TObject.cxx:187
bool verbose
RooAbsArg * absArg() const
Definition: RooArgProxy.h:37
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:51
static void indent(ostringstream &buf, int indent_level)
virtual void Print(Option_t *options=0) const
Print TNamed name and title.
const RooSimultaneous * _pdf
RooCmdArg Index(RooCategory &icat)
RooArgSet * _theEvent
virtual Int_t GetSize() const
Definition: TCollection.h:95
static Double_t uniform(TRandom *generator=randomGenerator())
Return a number uniformly distributed from (0,1)
Definition: RooRandom.cxx:83
double Double_t
Definition: RtypesCore.h:55
Bool_t canBeExtended() const
Definition: RooAbsPdf.h:216
virtual void printMultiline(std::ostream &os, Int_t contents, Bool_t verbose=kFALSE, TString indent="") const
Interface for multi-line printing.
std::vector< int > _gcIndex
virtual void setProtoDataOrder(Int_t *lut)
this method is empty because proto datasets are not supported by this context
Int_t getVal() const
Definition: RooCatType.h:80
virtual Bool_t remove(const RooAbsArg &var, Bool_t silent=kFALSE, Bool_t matchByNameOnly=kFALSE)
Remove the specified argument from our list.
RooAbsPdf is the abstract interface for all probability density functions The class provides hybrid a...
Definition: RooAbsPdf.h:41
virtual TObject * Next()=0
virtual ~RooSimSplitGenContext()
Destructor. Delete all owned subgenerator contexts.
virtual TIterator * MakeIterator(Bool_t dir=kIterForward) const
Return a list iterator.
Definition: TList.cxx:603
const RooAbsCategory & arg() const
RooAbsCategoryLValue * _idxCat
std::vector< RooAbsGenContext * > _gcList
virtual RooDataSet * generate(Double_t nEvents=0, Bool_t skipInit=kFALSE, Bool_t extendedMode=kFALSE)
Generate the specified number of events with nEvents>0 and and return a dataset containing the genera...
RooAbsArg is the common abstract base class for objects that represent a value (of arbitrary type) an...
Definition: RooAbsArg.h:66
virtual void printMultiline(std::ostream &os, Int_t content, Bool_t verbose=kFALSE, TString indent="") const
Detailed printing interface.
Bool_t recursiveRedirectServers(const RooAbsCollection &newServerList, Bool_t mustReplaceAll=kFALSE, Bool_t nameChange=kFALSE, Bool_t recurseInNewSet=kTRUE)
Definition: RooAbsArg.cxx:1087
const Bool_t kTRUE
Definition: Rtypes.h:91
RooSimSplitGenContext(const RooSimultaneous &model, const RooArgSet &vars, Bool_t _verbose=kFALSE, Bool_t autoBinned=kTRUE, const char *binnedTag="")
Constructor of specialized generator context for RooSimultaneous p.d.f.s.
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add element to non-owning set.
Definition: RooArgSet.cxx:448
const RooCatType * lookupType(Int_t index, Bool_t printError=kFALSE) const
Find our type corresponding to the specified index, or return 0 for no match.