Logo ROOT  
Reference Guide
RooGenContext.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 RooGenContext.cxx
19\class RooGenContext
20\ingroup Roofitcore
21
22Class RooGenContext implement a universal generator context for all
23RooAbsPdf classes that do not have or need a specialized generator
24context. This generator context queries the input p.d.f which observables
25it can generate internally and delegates generation of those observables
26to the p.d.f if it deems that safe. The other observables are generated
27use a RooAcceptReject sampling technique.
28**/
29
30#include "RooMsgService.h"
31#include "Riostream.h"
32
33#include "RooGenContext.h"
34#include "RooAbsPdf.h"
35#include "RooDataSet.h"
36#include "RooRealIntegral.h"
37#include "RooAcceptReject.h"
38#include "RooRealVar.h"
39#include "RooDataHist.h"
40#include "RooErrorHandler.h"
41#include "RooNumGenConfig.h"
42#include "RooNumGenFactory.h"
43
44#include "TString.h"
45#include "TIterator.h"
46#include "TClass.h"
47
48
49
50using namespace std;
51
53 ;
54
55
56
57////////////////////////////////////////////////////////////////////////////////
58/// Initialize a new context for generating events with the specified
59/// variables, using the specified PDF model. A prototype dataset (if provided)
60/// is not cloned and still belongs to the caller. The contents and shape
61/// of this dataset can be changed between calls to generate() as long as the
62/// expected columns to be copied to the generated dataset are present.
63/// Any argument supplied in the forceDirect RooArgSet are always offered
64/// for internal generation to the p.d.f., even if this is deemed unsafe by
65/// the logic of RooGenContext.
66
68 const RooDataSet *prototype, const RooArgSet* auxProto,
69 Bool_t verbose, const RooArgSet* forceDirect) :
70 RooAbsGenContext(model,vars,prototype,auxProto,verbose),
71 _pdfClone(0), _acceptRejectFunc(0), _generator(0),
72 _maxVar(0), _uniIter(0), _updateFMaxPerEvent(0)
73{
74 cxcoutI(Generation) << "RooGenContext::ctor() setting up event generator context for p.d.f. " << model.GetName()
75 << " for generation of observable(s) " << vars ;
76 if (prototype) ccxcoutI(Generation) << " with prototype data for " << *prototype->get() ;
77 if (auxProto && auxProto->getSize()>0) ccxcoutI(Generation) << " with auxiliary prototypes " << *auxProto ;
78 if (forceDirect && forceDirect->getSize()>0) ccxcoutI(Generation) << " with internal generation forced for observables " << *forceDirect ;
79 ccxcoutI(Generation) << endl ;
80
81
82 // Clone the model and all nodes that it depends on so that this context
83 // is independent of any existing objects.
84 RooArgSet nodes(model,model.GetName());
85 if (nodes.snapshot(_cloneSet, true)) {
86 coutE(Generation) << "RooGenContext::RooGenContext(" << GetName() << ") Couldn't deep-clone PDF, abort," << endl ;
88 }
89
90 // Find the clone in the snapshot list
91 _pdfClone = static_cast<RooAbsPdf*>(_cloneSet.find(model.GetName()));
93
94 // Optionally fix RooAddPdf normalizations
95 if (prototype&&_pdfClone->dependsOn(*prototype->get())) {
96 RooArgSet fullNormSet(vars) ;
97 fullNormSet.add(*prototype->get()) ;
99 }
100
101 // Analyze the list of variables to generate...
103 TIterator *iterator= vars.createIterator();
105 const RooAbsArg *tmp = 0;
106 const RooAbsArg *arg = 0;
107 while((_isValid && (tmp= (const RooAbsArg*)iterator->Next()))) {
108 // is this argument derived?
109 if(!tmp->isFundamental()) {
110 coutE(Generation) << "RooGenContext::ctor(): cannot generate values for derived \"" << tmp->GetName() << "\"" << endl;
112 continue;
113 }
114 // lookup this argument in the cloned set of PDF dependents
115 arg= static_cast<RooAbsArg const*>(_cloneSet.find(tmp->GetName()));
116 if(0 == arg) {
117 coutI(Generation) << "RooGenContext::ctor() WARNING model does not depend on \"" << tmp->GetName()
118 << "\" which will have uniform distribution" << endl;
119 _uniformVars.add(*tmp);
120 }
121 else {
122
123 // does the model depend on this variable directly, ie, like "x" in
124 // f(x) or f(x,g(x,y)) or even f(x,x) ?
125 const RooAbsArg *direct= arg ;
126 if (forceDirect==0 || !forceDirect->find(direct->GetName())) {
127 if (!_pdfClone->isDirectGenSafe(*arg)) {
128 cxcoutD(Generation) << "RooGenContext::ctor() observable " << arg->GetName() << " has been determined to be unsafe for internal generation" << endl;
129 direct=0 ;
130 }
131 }
132
133 // does the model depend indirectly on this variable through an lvalue chain?
134 // otherwise, this variable will have to be generated with accept/reject
135 if(direct) {
136 _directVars.add(*arg);
137 } else {
138 _otherVars.add(*arg);
139 }
140 }
141 }
142 delete servers;
143 delete iterator;
144 if(!isValid()) {
145 coutE(Generation) << "RooGenContext::ctor() constructor failed with errors" << endl;
146 return;
147 }
148
149 // At this point directVars are all variables that are safe to be generated directly
150 // otherVars are all variables that are _not_ safe to be generated directly
151 if (_directVars.getSize()>0) {
152 cxcoutD(Generation) << "RooGenContext::ctor() observables " << _directVars << " are safe for internal generator (if supported by p.d.f)" << endl ;
153 }
154 if (_otherVars.getSize()>0) {
155 cxcoutD(Generation) << "RooGenContext::ctor() observables " << _otherVars << " are NOT safe for internal generator (if supported by p.d.f)" << endl ;
156 }
157
158 // If PDF depends on prototype data, direct generator cannot use static initialization
159 // in initGenerator()
160 Bool_t staticInitOK = !_pdfClone->dependsOn(_protoVars) ;
161 if (!staticInitOK) {
162 cxcoutD(Generation) << "RooGenContext::ctor() Model depends on supplied protodata observables, static initialization of internal generator will not be allowed" << endl ;
163 }
164
165 // Can the model generate any of the direct variables itself?
166 RooArgSet generatedVars;
167 if (_directVars.getSize()>0) {
168 _code= _pdfClone->getGenerator(_directVars,generatedVars,staticInitOK);
169
170 cxcoutD(Generation) << "RooGenContext::ctor() Model indicates that it can internally generate observables "
171 << generatedVars << " with configuration identifier " << _code << endl ;
172 } else {
173 _code = 0 ;
174 }
175
176 // Move variables which cannot be generated into the list to be generated with accept/reject
177 _directVars.remove(generatedVars);
179
180 // Update _directVars to only include variables that will actually be directly generated
182 _directVars.add(generatedVars);
183
184 cxcoutI(Generation) << "RooGenContext::ctor() Context will" ;
185 if (_directVars.getSize()>0) ccxcoutI(Generation) << " let model internally generate observables " << _directVars ;
186 if (_directVars.getSize()>0 && _otherVars.getSize()>0) ccxcoutI(Generation) << " and" ;
187 if (_otherVars.getSize()>0) ccxcoutI(Generation) << " generate variables " << _otherVars << " with accept/reject sampling" ;
188 if (_uniformVars.getSize()>0) ccxcoutI(Generation) << ", non-observable variables " << _uniformVars << " will be generated with flat distribution" ;
189 ccxcoutI(Generation) << endl ;
190
191 // initialize the accept-reject generator
192 RooArgSet depList;
194 depList.remove(_otherVars);
195
196 TString nname(_pdfClone->GetName()) ;
197 nname.Append("_AccRej") ;
198 TString ntitle(_pdfClone->GetTitle()) ;
199 ntitle.Append(" (Accept/Reject)") ;
200
201
202 RooArgSet protoDeps;
203 model.getObservables(&_protoVars, protoDeps);
204
205
206 if (_protoVars.getSize()==0) {
207
208 // No prototype variables
209
210 if(depList.empty()) {
211 // All variable are generated with accept-reject
212
213 // Check if PDF supports maximum finding
214 Int_t maxFindCode = _pdfClone->getMaxVal(_otherVars) ;
215 if (maxFindCode != 0) {
216 coutI(Generation) << "RooGenContext::ctor() no prototype data provided, all observables are generated with numerically and "
217 << "model supports analytical maximum findin:, can provide analytical pdf maximum to numeric generator" << endl ;
218 Double_t maxVal = _pdfClone->maxVal(maxFindCode) / _pdfClone->getNorm(&_theEvent) ;
219 _maxVar = new RooRealVar("funcMax","function maximum",maxVal) ;
220 cxcoutD(Generation) << "RooGenContext::ctor() maximum value returned by RooAbsPdf::maxVal() is " << maxVal << endl ;
221 }
222 }
223
224 if (_otherVars.getSize()>0) {
225 _pdfClone->getVal(&vars) ; // WVE debug
226 _acceptRejectFunc= new RooRealIntegral(nname,ntitle,*_pdfClone,depList,&vars);
227 cxcoutI(Generation) << "RooGenContext::ctor() accept/reject sampling function is " << _acceptRejectFunc->GetName() << endl ;
228 } else {
230 }
231
232 } else {
233
234 // Generation _with_ prototype variable
235 depList.remove(_protoVars,true,true);
237 cxcoutI(Generation) << "RooGenContext::ctor() accept/reject sampling function is " << _acceptRejectFunc->GetName() << endl ;
238
239 // Check if PDF supports maximum finding for the entire phase space
240 RooArgSet allVars(_otherVars);
241 allVars.add(_directVars);
242 Int_t maxFindCode = _pdfClone->getMaxVal(allVars) ;
243 if (maxFindCode != 0) {
244 // Special case: PDF supports max-finding in otherVars, no need to scan other+proto space for maximum
245 coutI(Generation) << "RooGenContext::ctor() prototype data provided, and "
246 << "model supports analytical maximum finding in the full phase space: "
247 << "can provide analytical pdf maximum to numeric generator" << endl ;
248 _maxVar = new RooRealVar("funcMax","function maximum",_pdfClone->maxVal(maxFindCode)) ;
249 } else {
250 maxFindCode = _pdfClone->getMaxVal(_otherVars) ;
251 if (maxFindCode != 0) {
252 _updateFMaxPerEvent = maxFindCode ;
253 coutI(Generation) << "RooGenContext::ctor() prototype data provided, and "
254 << "model supports analytical maximum finding in the variables which are not"
255 << " internally generated. Can provide analytical pdf maximum to numeric "
256 << "generator" << endl;
257 cxcoutD(Generation) << "RooGenContext::ctor() maximum value must be reevaluated for each "
258 << "event with configuration code " << maxFindCode << endl ;
259 _maxVar = new RooRealVar("funcMax","function maximum",1) ;
260 }
261 }
262
263 if (!_maxVar) {
264
265 // Regular case: First find maximum in other+proto space
266 RooArgSet otherAndProto(_otherVars) ;
267
268 otherAndProto.add(protoDeps) ;
269
270 if (_otherVars.getSize()>0) {
271
272 cxcoutD(Generation) << "RooGenContext::ctor() prototype data provided, observables are generated numericaly no "
273 << "analytical estimate of maximum function value provided by model, must determine maximum value through initial sampling space "
274 << "of accept/reject observables plus prototype observables: " << otherAndProto << endl ;
275
276 // Calculate maximum in other+proto space if there are any accept/reject generated observables
278 *model.getGeneratorConfig(),_verbose) ;
279// RooAcceptReject maxFinder(*_acceptRejectFunc,otherAndProto,RooNumGenConfig::defaultConfig(),_verbose) ;
280 Double_t max = maxFinder->getFuncMax() ;
281 _maxVar = new RooRealVar("funcMax","function maximum",max) ;
282
283 if (max==0) {
284 coutE(Generation) << "RooGenContext::ctor(" << model.GetName()
285 << ") ERROR: generating conditional p.d.f. which requires prior knowledge of function maximum, "
286 << "but chosen numeric generator (" << maxFinder->IsA()->GetName() << ") does not support maximum finding" << endl ;
287 delete maxFinder ;
288 throw string("RooGenContext::ctor()") ;
289 }
290 delete maxFinder ;
291
292 cxcoutD(Generation) << "RooGenContext::ctor() maximum function value found through initial sampling is " << max << endl ;
293 }
294 }
295
296 }
297
300 cxcoutD(Generation) << "RooGenContext::ctor() creating MC sampling generator " << _generator->IsA()->GetName() << " from function for observables " << _otherVars << endl ;
301 //_generator= new RooAcceptReject(*_acceptRejectFunc,_otherVars,RooNumGenConfig::defaultConfig(),_verbose,_maxVar);
302 } else {
303 _generator = 0 ;
304 }
305
307}
308
309
310////////////////////////////////////////////////////////////////////////////////
311/// Destructor.
312
314{
315 // Clean up our accept/reject generator
316 if (_generator) delete _generator;
318 if (_maxVar) delete _maxVar ;
319 if (_uniIter) delete _uniIter ;
320}
321
322
323
324////////////////////////////////////////////////////////////////////////////////
325/// Attach the cloned model to the event buffer we will be filling.
326
328{
330 if (_acceptRejectFunc) {
332 }
333
334 // Attach the RooAcceptReject generator the event buffer
335 if (_generator) {
337 }
338
339}
340
341
342
343////////////////////////////////////////////////////////////////////////////////
344/// Perform one-time initialization of the generator context
345
347{
348 RooFIter iter = theEvent.fwdIterator() ;
349 RooAbsArg* arg ;
350 while((arg=iter.next())) {
352 }
353
354 attach(theEvent) ;
355
356 // Reset the cloned model's error counters.
358
359 // Initialize the PDFs internal generator
360 if (_directVars.getSize()>0) {
361 cxcoutD(Generation) << "RooGenContext::initGenerator() initializing internal generator of model with code " << _code << endl ;
363 }
364
365 // Create iterator for uniform vars (if any)
366 if (_uniformVars.getSize()>0) {
367 delete _uniIter;
369 }
370}
371
372
373////////////////////////////////////////////////////////////////////////////////
374/// Generate one event. The 'remaining' integer is not used other than
375/// for printing messages
376
378{
379 if(_otherVars.getSize() > 0) {
380 // call the accept-reject generator to generate its variables
381
382 if (_updateFMaxPerEvent!=0) {
384 cxcoutD(Generation) << "RooGenContext::initGenerator() reevaluation of maximum function value is required for each event, new value is " << max << endl ;
385 _maxVar->setVal(max) ;
386 }
387
388 if (_generator) {
389 Double_t resampleRatio(1) ;
390 const RooArgSet *subEvent= _generator->generateEvent(remaining,resampleRatio);
391 if (resampleRatio<1) {
392 coutI(Generation) << "RooGenContext::generateEvent INFO: accept/reject generator requests resampling of previously produced events by factor "
393 << resampleRatio << " due to increased maximum weight" << endl ;
394 resampleData(resampleRatio) ;
395 }
396 if(0 == subEvent) {
397 coutE(Generation) << "RooGenContext::generateEvent ERROR accept/reject generator failed" << endl ;
398 return;
399 }
400 theEvent.assignValueOnly(*subEvent) ;
401 //theEvent= *subEvent;
402
403 }
404 }
405
406 // Use the model's optimized generator, if one is available.
407 // The generator writes directly into our local 'event' since we attached it above.
408 if(_directVars.getSize() > 0) {
410 }
411
412 // Generate uniform variables (non-dependents)
413 if (_uniIter) {
414 _uniIter->Reset() ;
415 RooAbsArg* uniVar ;
416 while((uniVar=(RooAbsArg*)_uniIter->Next())) {
417 RooAbsLValue* arglv = dynamic_cast<RooAbsLValue*>(uniVar) ;
418 if (!arglv) {
419 coutE(Generation) << "RooGenContext::generateEvent(" << GetName() << ") ERROR: uniform variable " << uniVar->GetName() << " is not an lvalue" << endl ;
421 }
422 arglv->randomize() ;
423 }
424 theEvent.assign(_uniformVars) ;
425 }
426
427}
428
429
430////////////////////////////////////////////////////////////////////////////////
431/// Printing interface
432
434{
436 os << indent << " --- RooGenContext --- " << endl ;
437 os << indent << "Using PDF ";
439
440 if(verbose) {
441 os << indent << "Use PDF generator for " << _directVars << endl ;
442 os << indent << "Use MC sampling generator " << (_generator ? _generator->IsA()->GetName() : "<none>") << " for " << _otherVars << endl ;
443 if (_protoVars.getSize()>0) {
444 os << indent << "Prototype observables are " << _protoVars << endl ;
445 }
446 }
447}
#define coutI(a)
Definition: RooMsgService.h:30
#define cxcoutI(a)
Definition: RooMsgService.h:85
#define cxcoutD(a)
Definition: RooMsgService.h:81
#define coutE(a)
Definition: RooMsgService.h:33
#define ccxcoutI(a)
Definition: RooMsgService.h:86
const Bool_t kFALSE
Definition: RtypesCore.h:101
const Bool_t kTRUE
Definition: RtypesCore.h:100
#define ClassImp(name)
Definition: Rtypes.h:375
static void indent(ostringstream &buf, int indent_level)
winID h direct
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition: RooAbsArg.h:78
RooArgSet * getObservables(const RooArgSet &set, Bool_t valueOnly=kTRUE) const
Given a set of possible observables, return the observables that this PDF depends on.
Definition: RooAbsArg.h:318
Bool_t dependsOn(const RooAbsCollection &serverList, const RooAbsArg *ignoreArg=0, Bool_t valueOnly=kFALSE) const
Test whether we depend on (ie, are served by) any object in the specified collection.
Definition: RooAbsArg.cxx:835
virtual Bool_t isFundamental() const
Is this object a fundamental type that can be added to a dataset? Fundamental-type subclasses overrid...
Definition: RooAbsArg.h:248
void setOperMode(OperMode mode, Bool_t recurseADirty=kTRUE)
Set the operation mode of this node.
Definition: RooAbsArg.cxx:1876
TIterator * serverIterator() const
Definition: RooAbsArg.h:146
Bool_t recursiveRedirectServers(const RooAbsCollection &newServerList, Bool_t mustReplaceAll=kFALSE, Bool_t nameChange=kFALSE, Bool_t recurseInNewSet=kTRUE)
Recursively replace all servers with the new servers in newSet.
Definition: RooAbsArg.cxx:1197
virtual void removeAll()
Remove all arguments from our set, deleting them if we own them.
RooAbsCollection & assignValueOnly(const RooAbsCollection &other, bool forceIfSizeOne=false)
Sets the value of any argument in our set that also appears in the other set.
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
bool empty() const
Int_t getSize() const
Return the number of elements in the collection.
RooFIter fwdIterator() const
One-time forward iterator.
void assign(const RooAbsCollection &other) const
Sets the value, cache and constant attribute of any argument in our set that also appears in the othe...
virtual Bool_t remove(const RooAbsArg &var, Bool_t silent=kFALSE, Bool_t matchByNameOnly=kFALSE)
Remove the specified argument from our list.
TIterator * createIterator(Bool_t dir=kIterForward) const
TIterator-style iteration over contained elements.
RooAbsArg * find(const char *name) const
Find object with given name in list.
RooAbsGenContext is the abstract base class for generator contexts of RooAbsPdf objects.
void printMultiline(std::ostream &os, Int_t contents, Bool_t verbose=kFALSE, TString indent="") const override
Interface for multi-line printing.
Bool_t isValid() const
RooArgSet _theEvent
Pointer to observable event being generated.
Bool_t _isValid
Is context in valid state?
void resampleData(Double_t &ratio)
Rescale existing output buffer with given ratio.
RooArgSet _protoVars
Prototype observables.
Bool_t _verbose
Verbose messaging?
Abstract base class for objects that are lvalues, i.e.
Definition: RooAbsLValue.h:26
virtual void randomize(const char *rangeName=0)=0
Class RooAbsNumGenerator is the abstract base class for MC event generator implementations like RooAc...
TClass * IsA() const override
virtual Double_t getFuncMax()
void attachParameters(const RooArgSet &vars)
Reattach original parameters to function clone.
virtual const RooArgSet * generateEvent(UInt_t remaining, Double_t &resampleRatio)=0
Double_t getNorm(const RooArgSet &nset) const
Get normalisation term needed to normalise the raw values returned by getVal().
Definition: RooAbsPdf.h:240
virtual void generateEvent(Int_t code)
Interface for generation of an event using the algorithm corresponding to the specified code.
Definition: RooAbsPdf.cxx:2324
virtual void resetErrorCounters(Int_t resetValue=10)
Reset error counter to given value, limiting the number of future error messages for this pdf to 'res...
Definition: RooAbsPdf.cxx:604
virtual Int_t getGenerator(const RooArgSet &directVars, RooArgSet &generateVars, Bool_t staticInitOK=kTRUE) const
Load generatedVars with the subset of directVars that we can generate events for, and return a code t...
Definition: RooAbsPdf.cxx:2302
virtual Bool_t isDirectGenSafe(const RooAbsArg &arg) const
Check if given observable can be safely generated using the pdfs internal generator mechanism (if tha...
Definition: RooAbsPdf.cxx:2337
const RooNumGenConfig * getGeneratorConfig() const
Return the numeric MC generator configuration used for this object.
Definition: RooAbsPdf.cxx:3470
virtual void initGenerator(Int_t code)
Interface for one-time initialization to setup the generator for the specified code.
Definition: RooAbsPdf.cxx:2312
RooAbsReal * createIntegral(const RooArgSet &iset, const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
Create an object that represents the integral of the function over one or more observables listed in ...
Definition: RooAbsReal.cxx:553
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:93
virtual Double_t maxVal(Int_t code) const
Return maximum value for set of observables identified by code assigned in getMaxVal.
virtual void fixAddCoefNormalization(const RooArgSet &addNormSet=RooArgSet(), Bool_t force=kTRUE)
Fix the interpretation of the coefficient of any RooAddPdf component in the expression tree headed by...
virtual Int_t getMaxVal(const RooArgSet &vars) const
Advertise capability to determine maximum value of function for given set of observables.
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:57
RooArgSet * snapshot(bool deepCopy=true) const
Use RooAbsCollection::snapshot(), but return as RooArgSet.
Definition: RooArgSet.h:180
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:55
const RooArgSet * get(Int_t index) const override
Return RooArgSet with coordinates of event 'index'.
static void softAbort()
Soft abort function that interrupts macro execution but doesn't kill ROOT.
A one-time forward iterator working on RooLinkedList or RooAbsCollection.
RooAbsArg * next()
Return next element or nullptr if at end.
Class RooGenContext implement a universal generator context for all RooAbsPdf classes that do not hav...
Definition: RooGenContext.h:30
void generateEvent(RooArgSet &theEvent, Int_t remaining) override
Generate one event.
RooArgSet _otherVars
List of observables generated internally, randomly, and by accept/reject sampling.
Definition: RooGenContext.h:47
void printMultiline(std::ostream &os, Int_t content, Bool_t verbose=kFALSE, TString indent="") const override
Printing interface.
Int_t _code
Internal generation code.
Definition: RooGenContext.h:48
RooArgSet _uniformVars
Definition: RooGenContext.h:47
void initGenerator(const RooArgSet &theEvent) override
Perform one-time initialization of the generator context.
~RooGenContext() override
Destructor.
RooArgSet _directVars
Definition: RooGenContext.h:47
TIterator * _uniIter
Iterator over uniform observables.
Definition: RooGenContext.h:53
RooAbsPdf * _pdfClone
Clone of input p.d.f.
Definition: RooGenContext.h:46
RooArgSet _cloneSet
Clone of all nodes of input p.d.f.
Definition: RooGenContext.h:45
RooGenContext(const RooAbsPdf &model, const RooArgSet &vars, const RooDataSet *prototype=0, const RooArgSet *auxProto=0, Bool_t verbose=kFALSE, const RooArgSet *forceDirect=0)
Initialize a new context for generating events with the specified variables, using the specified PDF ...
Int_t _updateFMaxPerEvent
If true, maximum p.d.f value needs to be recalculated for each event.
Definition: RooGenContext.h:54
RooRealIntegral * _acceptRejectFunc
Projection function to be passed to accept/reject sampler.
Definition: RooGenContext.h:50
void attach(const RooArgSet &params) override
Attach the cloned model to the event buffer we will be filling.
RooAbsNumGenerator * _generator
MC sampling generation engine.
Definition: RooGenContext.h:51
RooRealVar * _maxVar
Variable holding maximum value of p.d.f.
Definition: RooGenContext.h:52
static RooNumGenFactory & instance()
Static method returning reference to singleton instance of factory.
RooAbsNumGenerator * createSampler(RooAbsReal &func, const RooArgSet &genVars, const RooArgSet &condVars, const RooNumGenConfig &config, Bool_t verbose=kFALSE, RooAbsReal *maxFuncVal=0)
Construct a numeric integrator instance that operates on function 'func' and is configured with 'conf...
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,...
RooRealIntegral performs hybrid numerical/analytical integrals of RooAbsReal objects.
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:40
void setVal(Double_t value) override
Set value of variable to 'value'.
Definition: RooRealVar.cxx:281
Iterator abstract base class.
Definition: TIterator.h:30
virtual void Reset()=0
virtual TObject * Next()=0
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
TString & Append(const char *cs)
Definition: TString.h:564
@ Generation
Definition: RooGlobalFunc.h:63