Logo ROOT  
Reference Guide
RooBinnedGenContext.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 RooBinnedGenContext.cxx
19\class RooBinnedGenContext
20\ingroup Roofitcore
21
22RooBinnedGenContext is an efficient implementation of the
23generator context specific for binned pdfs.
24**/
25
26
27#include "RooFit.h"
28
29#include "Riostream.h"
30
31
32#include "RooMsgService.h"
33#include "RooBinnedGenContext.h"
34#include "RooAbsPdf.h"
35#include "RooRealVar.h"
36#include "RooDataHist.h"
37#include "RooDataSet.h"
38#include "RooRandom.h"
39
40using namespace std;
41
43;
44
45
46////////////////////////////////////////////////////////////////////////////////
47/// Constructor
48
50 const RooDataSet *prototype, const RooArgSet* auxProto,
52 RooAbsGenContext(model,vars,prototype,auxProto,verbose)
53{
54 cxcoutI(Generation) << "RooBinnedGenContext::ctor() setting up event special generator context for sum p.d.f. " << model.GetName()
55 << " for generation of observable(s) " << vars ;
56 if (prototype) ccxcoutI(Generation) << " with prototype data for " << *prototype->get() ;
57 if (auxProto && auxProto->getSize()>0) ccxcoutI(Generation) << " with auxiliary prototypes " << *auxProto ;
58 ccxcoutI(Generation) << endl ;
59
60 // Constructor. Build an array of generator contexts for each product component PDF
62 _pdf = (RooAbsPdf*) _pdfSet->find(model.GetName()) ;
64
65 // Fix normalization set of this RooAddPdf
66 if (prototype)
67 {
68 RooArgSet coefNSet(vars) ;
69 coefNSet.add(*prototype->get()) ;
70 _pdf->fixAddCoefNormalization(coefNSet) ;
71 }
72
74 _vars = _pdf->getObservables(vars) ;
75
76 // If pdf has boundary definitions, follow those for the binning
77 RooFIter viter = _vars->fwdIterator() ;
78 RooAbsArg* var ;
79 while((var=viter.next())) {
80 RooRealVar* rvar = dynamic_cast<RooRealVar*>(var) ;
81 if (rvar) {
82 list<Double_t>* binb = model.binBoundaries(*rvar,rvar->getMin(),rvar->getMax()) ;
83 delete binb ;
84 }
85 }
86
87
88 // Create empty RooDataHist
89 _hist = new RooDataHist("genData","genData",*_vars) ;
90
92}
93
94
95////////////////////////////////////////////////////////////////////////////////
96/// Destructor. Delete all owned subgenerator contexts
97
99{
100 delete _vars ;
101 delete _pdfSet ;
102 delete _hist ;
103}
104
105
106
107////////////////////////////////////////////////////////////////////////////////
108/// Attach given set of variables to internal p.d.f. clone
109
111{
113}
114
115
116
117////////////////////////////////////////////////////////////////////////////////
118/// One-time initialization of generator contex. Attach theEvent
119/// to internal p.d.f clone and forward initialization call to
120/// the component generators
121
123{
124 _pdf->recursiveRedirectServers(theEvent) ;
125
126}
127
128
129////////////////////////////////////////////////////////////////////////////////
130
132{
133 _expectedData = flag ;
134}
135
136
137////////////////////////////////////////////////////////////////////////////////
138
140{
141 // Scale to number of events and introduce Poisson fluctuations
142 _hist->reset() ;
143
144 Double_t nEvents = nEvt ;
145
146 if (nEvents<=0) {
147 if (!_pdf->canBeExtended()) {
148 coutE(InputArguments) << "RooAbsPdf::generateBinned(" << GetName()
149 << ") ERROR: No event count provided and p.d.f does not provide expected number of events" << endl ;
150 return 0 ;
151 } else {
152 // Don't round in expectedData mode
153 if (_expectedData || extended) {
154 nEvents = _pdf->expectedEvents(_vars) ;
155 } else {
156 nEvents = Int_t(_pdf->expectedEvents(_vars)+0.5) ;
157 }
158 }
159 }
160
161 // Sample p.d.f. distribution
163
164 // Output container
165 RooRealVar weight("weight","weight",0,1e9) ;
166 RooArgSet tmp(*_vars) ;
167 tmp.add(weight) ;
168 RooDataSet* wudata = new RooDataSet("wu","wu",tmp,RooFit::WeightVar("weight")) ;
169
170 vector<int> histOut(_hist->numEntries()) ;
171 Double_t histMax(-1) ;
172 Int_t histOutSum(0) ;
173 for (int i=0 ; i<_hist->numEntries() ; i++) {
174 _hist->get(i) ;
175 if (_expectedData) {
176
177 // Expected data, multiply p.d.f by nEvents
178 Double_t w=_hist->weight()*nEvents ;
179 wudata->add(*_hist->get(),w) ;
180
181 } else if (extended) {
182
183 // Extended mode, set contents to Poisson(pdf*nEvents)
185 wudata->add(*_hist->get(),w) ;
186
187 } else {
188
189 // Regular mode, fill array of weights with Poisson(pdf*nEvents), but to not fill
190 // histogram yet.
191 if (_hist->weight()>histMax) {
192 histMax = _hist->weight() ;
193 }
194 histOut[i] = RooRandom::randomGenerator()->Poisson(_hist->weight()*nEvents) ;
195 histOutSum += histOut[i] ;
196 }
197 }
198
199
200 if (!_expectedData && !extended) {
201
202 // Second pass for regular mode - Trim/Extend dataset to exact number of entries
203
204 // Calculate difference between what is generated so far and what is requested
205 Int_t nEvtExtra = abs(Int_t(nEvents)-histOutSum) ;
206 Int_t wgt = (histOutSum>nEvents) ? -1 : 1 ;
207
208 // Perform simple binned accept/reject procedure to get to exact event count
209 while(nEvtExtra>0) {
210
212 _hist->get(ibinRand) ;
213 Double_t ranY = RooRandom::randomGenerator()->Uniform(histMax) ;
214
215 if (ranY<_hist->weight()) {
216 if (wgt==1) {
217 histOut[ibinRand]++ ;
218 } else {
219 // If weight is negative, prior bin content must be at least 1
220 if (histOut[ibinRand]>0) {
221 histOut[ibinRand]-- ;
222 } else {
223 continue ;
224 }
225 }
226 nEvtExtra-- ;
227 }
228 }
229
230 // Transfer working array to histogram
231 for (int i=0 ; i<_hist->numEntries() ; i++) {
232 _hist->get(i) ;
233 wudata->add(*_hist->get(),histOut[i]) ;
234 }
235
236 }
237
238 return wudata ;
239
240}
241
242
243
244////////////////////////////////////////////////////////////////////////////////
245/// this method is not implemented for this context
246
248{
249 assert(0) ;
250}
251
252
253
254////////////////////////////////////////////////////////////////////////////////
255/// Print the details of the context
256
258{
260 os << indent << "--- RooBinnedGenContext ---" << endl ;
261 os << indent << "Using PDF ";
263}
#define cxcoutI(a)
Definition: RooMsgService.h:85
#define coutE(a)
Definition: RooMsgService.h:33
#define ccxcoutI(a)
Definition: RooMsgService.h:86
int Int_t
Definition: RtypesCore.h:43
const Bool_t kFALSE
Definition: RtypesCore.h:90
bool Bool_t
Definition: RtypesCore.h:61
double Double_t
Definition: RtypesCore.h:57
const Bool_t kTRUE
Definition: RtypesCore.h:89
#define ClassImp(name)
Definition: Rtypes.h:361
static void indent(ostringstream &buf, int indent_level)
RooAbsArg is the common abstract base class for objects that represent a value (of arbitrary type) an...
Definition: RooAbsArg.h:73
RooArgSet * getObservables(const RooArgSet &set, Bool_t valueOnly=kTRUE) const
Return the observables of this pdf given a set of observables.
Definition: RooAbsArg.h:276
void setOperMode(OperMode mode, Bool_t recurseADirty=kTRUE)
Change cache operation mode to given mode.
Definition: RooAbsArg.cxx:1718
Bool_t recursiveRedirectServers(const RooAbsCollection &newServerList, Bool_t mustReplaceAll=kFALSE, Bool_t nameChange=kFALSE, Bool_t recurseInNewSet=kTRUE)
Definition: RooAbsArg.cxx:1064
Int_t getSize() const
RooFIter fwdIterator() const
One-time forward iterator.
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.
RooArgSet * _theEvent
virtual void printMultiline(std::ostream &os, Int_t contents, Bool_t verbose=kFALSE, TString indent="") const
Interface for multi-line printing.
Bool_t canBeExtended() const
Definition: RooAbsPdf.h:236
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:3283
virtual Double_t getMax(const char *name=0) const
Get maximum of currently defined range.
virtual Double_t getMin(const char *name=0) const
Get miniminum of currently defined range.
virtual std::list< Double_t > * binBoundaries(RooAbsRealLValue &, Double_t, Double_t) const
Definition: RooAbsReal.h:313
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...
RooDataHist * fillDataHist(RooDataHist *hist, const RooArgSet *nset, Double_t scaleFactor, Bool_t correctForBinVolume=kFALSE, Bool_t showProgress=kFALSE) const
Fill a RooDataHist with values sampled from this function at the bin centers.
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:28
RooArgSet * snapshot(bool deepCopy=true) const
Use RooAbsCollection::snapshot(), but return as RooArgSet.
Definition: RooArgSet.h:126
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
RooBinnedGenContext is an efficient implementation of the generator context specific for binned pdfs.
RooBinnedGenContext(const RooAbsPdf &model, const RooArgSet &vars, const RooDataSet *prototype=0, const RooArgSet *auxProto=0, Bool_t _verbose=kFALSE)
Constructor.
virtual void initGenerator(const RooArgSet &theEvent)
One-time initialization of generator contex.
const RooArgSet * _vars
virtual void attach(const RooArgSet &params)
Attach given set of variables to internal p.d.f. clone.
virtual void generateEvent(RooArgSet &theEvent, Int_t remaining)
this method is not implemented for this context
virtual ~RooBinnedGenContext()
Destructor. Delete all owned subgenerator contexts.
virtual void setExpectedData(Bool_t)
virtual void printMultiline(std::ostream &os, Int_t content, Bool_t verbose=kFALSE, TString indent="") const
Print the details of the context.
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...
The RooDataHist is a container class to hold N-dimensional binned data.
Definition: RooDataHist.h:40
virtual void reset()
Reset all bin weights to zero.
virtual Double_t weight() const
Definition: RooDataHist.h:106
virtual Int_t numEntries() const
Return the number of bins.
virtual const RooArgSet * get() const
Definition: RooDataHist.h:79
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:33
virtual const RooArgSet * get(Int_t index) const override
Return RooArgSet with coordinates of event 'index'.
virtual void add(const RooArgSet &row, Double_t weight=1.0, Double_t weightError=0) override
Add a data point, with its coordinates specified in the 'data' argset, to the data set.
A one-time forward iterator working on RooLinkedList or RooAbsCollection.
RooAbsArg * next()
Return next element or nullptr if at end.
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,...
static TRandom * randomGenerator()
Return a pointer to a singleton random-number generator implementation.
Definition: RooRandom.cxx:53
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:35
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
virtual Int_t Poisson(Double_t mean)
Generates a random integer N according to a Poisson law.
Definition: TRandom.cxx:391
virtual Double_t Uniform(Double_t x1=1)
Returns a uniform deviate on the interval (0, x1).
Definition: TRandom.cxx:635
virtual UInt_t Integer(UInt_t imax)
Returns a random integer uniformly distributed on the interval [ 0, imax-1 ].
Definition: TRandom.cxx:349
Basic string class.
Definition: TString.h:131
RooCmdArg WeightVar(const char *name, Bool_t reinterpretAsWeight=kFALSE)
@ Generation
Definition: RooGlobalFunc.h:67
@ InputArguments
Definition: RooGlobalFunc.h:68