Logo ROOT  
Reference Guide
RooConvGenContext.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 RooConvGenContext.cxx
19\class RooConvGenContext
20\ingroup Roofitcore
21
22RooConvGenContext is an efficient implementation of the generator context
23specific for RooAbsAnaConvPdf objects. The physics model is generated
24with a truth resolution model and the requested resolution model is generated
25separately as a PDF. The convolution variable of the physics model is
26subsequently explicitly smeared with the resolution model distribution.
27**/
28
29#include "RooMsgService.h"
30#include "RooErrorHandler.h"
31#include "RooConvGenContext.h"
32#include "RooAbsAnaConvPdf.h"
33#include "RooNumConvPdf.h"
34#include "RooFFTConvPdf.h"
35#include "RooProdPdf.h"
36#include "RooDataSet.h"
37#include "RooArgSet.h"
38#include "RooTruthModel.h"
39#include "Riostream.h"
40
41
42using namespace std;
43
45;
46
47
48////////////////////////////////////////////////////////////////////////////////
49/// Constructor for specialized generator context for analytical convolutions.
50///
51/// Builds a generator for the physics PDF convoluted with the truth model
52/// and a generator for the resolution model as PDF. Events are generated
53/// by sampling events from the p.d.f and smearings from the resolution model
54/// and adding these to obtain a distribution of events consistent with the
55/// convolution of these two. The advantage of this procedure is so that
56/// both p.d.f and resolution model can take advantage of any internal
57/// generators that may be defined.
58
60 const RooDataSet *prototype, const RooArgSet* auxProto, bool verbose) :
61 RooAbsGenContext(model,vars,prototype,auxProto,verbose)
62{
63 cxcoutI(Generation) << "RooConvGenContext::ctor() setting up special generator context for analytical convolution p.d.f. " << model.GetName()
64 << " for generation of observable(s) " << vars << endl ;
65
66 // Clone PDF and change model to internal truth model
67 _pdfCloneSet.reset(RooArgSet(model).snapshot(true));
68 if (!_pdfCloneSet) {
69 coutE(Generation) << "RooConvGenContext::RooConvGenContext(" << GetName() << ") Couldn't deep-clone PDF, abort," << endl ;
71 }
72
73 RooAbsAnaConvPdf* pdfClone = (RooAbsAnaConvPdf*) _pdfCloneSet->find(model.GetName()) ;
74 RooTruthModel truthModel("truthModel","Truth resolution model",*pdfClone->convVar()) ;
75 pdfClone->changeModel(truthModel) ;
76 auto convV = dynamic_cast<RooRealVar*>(pdfClone->convVar());
77 if (!convV) {
78 throw std::runtime_error("RooConvGenContext only works with convolution variables of type RooRealVar.");
79 }
80 convV->removeRange();
81
82 // Create generator for physics X truth model
83 _pdfVars.reset(pdfClone->getObservables(&vars));
84 _pdfGen.reset(pdfClone->genContext(*_pdfVars,prototype,auxProto,verbose));
85
86 // Clone resolution model and use as normal PDF
87 _modelCloneSet.reset(RooArgSet(*model._convSet.at(0)).snapshot(true));
88 if (!_modelCloneSet) {
89 coutE(Generation) << "RooConvGenContext::RooConvGenContext(" << GetName() << ") Couldn't deep-clone resolution model, abort," << std::endl;
91 }
92 auto modelClone = static_cast<RooResolutionModel*>(_modelCloneSet->find(model._convSet.at(0)->GetName())->Clone("smearing"));
93 _modelCloneSet->addOwned(*modelClone) ;
94 modelClone->changeBasis(0) ;
95 convV = dynamic_cast<RooRealVar*>(&modelClone->convVar());
96 if (!convV) {
97 throw std::runtime_error("RooConvGenContext only works with convolution variables of type RooRealVar.");
98 }
99 convV->removeRange();
100
101 // Create generator for resolution model as PDF
102 _modelVars.reset(modelClone->getObservables(&vars));
103
104 _modelVars->add(modelClone->convVar()) ;
105 _convVarName = modelClone->convVar().GetName() ;
106 _modelGen.reset(modelClone->genContext(*_modelVars,prototype,auxProto,verbose));
107
108 if (prototype) {
109 _pdfVars->add(*prototype->get()) ;
110 _modelVars->add(*prototype->get()) ;
111 }
112
113 // WVE ADD FOR DEBUGGING
114 if (auxProto) {
115 _pdfVars->add(*auxProto) ;
116 _modelVars->add(*auxProto) ;
117 }
118}
119
120
121
122////////////////////////////////////////////////////////////////////////////////
123/// Constructor for specialized generator context for numerical convolutions.
124///
125/// Builds a generator for the physics PDF convoluted with the truth model
126/// and a generator for the resolution model as PDF. Events are generated
127/// by sampling events from the p.d.f and smearings from the resolution model
128/// and adding these to obtain a distribution of events consistent with the
129/// convolution of these two. The advantage of this procedure is so that
130/// both p.d.f and resolution model can take advantage of any internal
131/// generators that may be defined.
132
134 const RooDataSet *prototype, const RooArgSet* auxProto, bool verbose) :
135 RooAbsGenContext(model,vars,prototype,auxProto,verbose)
136{
137 cxcoutI(Generation) << "RooConvGenContext::ctor() setting up special generator context for numeric convolution p.d.f. " << model.GetName()
138 << " for generation of observable(s) " << vars << endl ;
139
140 // Create generator for physics X truth model
141 {
142 RooArgSet clonedPdfObservables;
143 model.conv().clonePdf().getObservables(&vars, clonedPdfObservables);
144 _pdfVarsOwned.reset(clonedPdfObservables.snapshot(true));
145 }
146 _pdfVars = std::make_unique<RooArgSet>(*_pdfVarsOwned) ;
147 _pdfGen.reset(static_cast<RooAbsPdf&>(model.conv().clonePdf()).genContext(*_pdfVars,prototype,auxProto,verbose));
148
149 // Create generator for resolution model as PDF
150 {
151 RooArgSet clonedModelObservables;
152 model.conv().cloneModel().getObservables(&vars, clonedModelObservables);
153 _modelVarsOwned.reset(clonedModelObservables.snapshot(true));
154 }
155 _modelVars = std::make_unique<RooArgSet>(*_modelVarsOwned) ;
156 _convVarName = model.conv().cloneVar().GetName() ;
157 _modelGen.reset(static_cast<RooAbsPdf&>(model.conv().cloneModel()).genContext(*_modelVars,prototype,auxProto,verbose));
158 _modelCloneSet = std::make_unique<RooArgSet>();
159 _modelCloneSet->add(model.conv().cloneModel()) ;
160
161 if (prototype) {
162 _pdfVars->add(*prototype->get()) ;
163 _modelVars->add(*prototype->get()) ;
164 }
165}
166
167
168
169////////////////////////////////////////////////////////////////////////////////
170/// Constructor for specialized generator context for FFT numerical convolutions.
171///
172/// Builds a generator for the physics PDF convoluted with the truth model
173/// and a generator for the resolution model as PDF. Events are generated
174/// by sampling events from the p.d.f and smearings from the resolution model
175/// and adding these to obtain a distribution of events consistent with the
176/// convolution of these two. The advantage of this procedure is so that
177/// both p.d.f and resolution model can take advantage of any internal
178/// generators that may be defined.
179
181 const RooDataSet *prototype, const RooArgSet* auxProto, bool verbose) :
182 RooAbsGenContext(model,vars,prototype,auxProto,verbose)
183{
184 cxcoutI(Generation) << "RooConvGenContext::ctor() setting up special generator context for fft convolution p.d.f. " << model.GetName()
185 << " for generation of observable(s) " << vars << endl ;
186
187 _convVarName = model._x.arg().GetName() ;
188
189 // Create generator for physics model
190 _pdfCloneSet.reset(RooArgSet(model._pdf1.arg()).snapshot(true));
191 RooAbsPdf* pdfClone = (RooAbsPdf*) _pdfCloneSet->find(model._pdf1.arg().GetName()) ;
192 RooRealVar* cvPdf = (RooRealVar*) _pdfCloneSet->find(model._x.arg().GetName()) ;
193 cvPdf->removeRange() ;
194 RooArgSet tmp1;
195 pdfClone->getObservables(&vars, tmp1) ;
196 _pdfVarsOwned.reset(tmp1.snapshot(true));
197 _pdfVars = std::make_unique<RooArgSet>(*_pdfVarsOwned) ;
198 _pdfGen.reset(pdfClone->genContext(*_pdfVars,prototype,auxProto,verbose));
199
200 // Create generator for resolution model
201 _modelCloneSet.reset(RooArgSet(model._pdf2.arg()).snapshot(true));
202 RooAbsPdf* modelClone = (RooAbsPdf*) _modelCloneSet->find(model._pdf2.arg().GetName()) ;
203 RooRealVar* cvModel = (RooRealVar*) _modelCloneSet->find(model._x.arg().GetName()) ;
204 cvModel->removeRange() ;
205 RooArgSet tmp2;
206 modelClone->getObservables(&vars, tmp2) ;
207 _modelVarsOwned.reset(tmp2.snapshot(true));
208 _modelVars = std::make_unique<RooArgSet>(*_modelVarsOwned) ;
209 _modelGen.reset(modelClone->genContext(*_pdfVars,prototype,auxProto,verbose));
210
211 if (prototype) {
212 _pdfVars->add(*prototype->get()) ;
213 _modelVars->add(*prototype->get()) ;
214 }
215}
216
217
218////////////////////////////////////////////////////////////////////////////////
219/// Attach given set of arguments to internal clones of
220/// pdf and resolution model
221
223{
224 // Find convolution variable in input and output sets
225 auto* cvModel = static_cast<RooRealVar*>(_modelVars->find(_convVarName));
226 auto* cvPdf = static_cast<RooRealVar*>(_pdfVars->find(_convVarName));
227
228 // Replace all servers in _pdfVars and _modelVars with those in theEvent, except for the convolution variable
229 std::unique_ptr<RooArgSet> pdfCommon{static_cast<RooArgSet*>(args.selectCommon(*_pdfVars))};
230 pdfCommon->remove(*cvPdf,true,true) ;
231
232 std::unique_ptr<RooArgSet> modelCommon{static_cast<RooArgSet*>(args.selectCommon(*_modelVars))};
233 modelCommon->remove(*cvModel,true,true) ;
234
235 _pdfGen->attach(*pdfCommon) ;
236 _modelGen->attach(*modelCommon) ;
237}
238
239
240////////////////////////////////////////////////////////////////////////////////
241/// One-time initialization of generator context, attaches
242/// the context to the supplied event container
243
245{
246 // Find convolution variable in input and output sets
249 _cvOut = (RooRealVar*) theEvent.find(_convVarName) ;
250
251 // Replace all servers in _pdfVars and _modelVars with those in theEvent, except for the convolution variable
252 std::unique_ptr<RooArgSet> pdfCommon{static_cast<RooArgSet*>(theEvent.selectCommon(*_pdfVars))};
253 pdfCommon->remove(*_cvPdf,true,true) ;
254 _pdfVars->replace(*pdfCommon) ;
255
256 std::unique_ptr<RooArgSet> modelCommon{static_cast<RooArgSet*>(theEvent.selectCommon(*_modelVars))};
257 modelCommon->remove(*_cvModel,true,true) ;
258 _modelVars->replace(*modelCommon) ;
259
260 // Initialize component generators
261 _pdfGen->initGenerator(*_pdfVars) ;
262 _modelGen->initGenerator(*_modelVars) ;
263}
264
265
266
267////////////////////////////////////////////////////////////////////////////////
268/// Generate a single event
269
271{
272 while(1) {
273
274 // Generate pdf and model data
275 _modelGen->generateEvent(*_modelVars,remaining) ;
276 _pdfGen->generateEvent(*_pdfVars,remaining) ;
277
278 // Construct smeared convolution variable
279 double convValSmeared = _cvPdf->getVal() + _cvModel->getVal() ;
280 if (_cvOut->isValidReal(convValSmeared)) {
281 // Smeared value in acceptance range, transfer values to output set
282 theEvent.assign(*_modelVars) ;
283 theEvent.assign(*_pdfVars) ;
284 _cvOut->setVal(convValSmeared) ;
285 return ;
286 }
287 }
288}
289
290
291
292////////////////////////////////////////////////////////////////////////////////
293/// Set the traversal order for events in the prototype dataset
294/// The argument is a array of integers with a size identical
295/// to the number of events in the prototype dataset. Each element
296/// should contain an integer in the range 1-N.
297
299{
301 _modelGen->setProtoDataOrder(lut) ;
302 _pdfGen->setProtoDataOrder(lut) ;
303}
304
305
306////////////////////////////////////////////////////////////////////////////////
307/// Print the details of this generator context
308
309void RooConvGenContext::printMultiline(ostream &os, Int_t content, bool verbose, TString indent) const
310{
312 os << indent << "--- RooConvGenContext ---" << endl ;
313 os << indent << "List of component generators" << endl ;
314
315 TString indent2(indent) ;
316 indent2.Append(" ") ;
317
318 _modelGen->printMultiline(os,content,verbose,indent2);
319 _pdfGen->printMultiline(os,content,verbose,indent2);
320}
#define cxcoutI(a)
Definition: RooMsgService.h:89
#define coutE(a)
Definition: RooMsgService.h:37
#define ClassImp(name)
Definition: Rtypes.h:375
static void indent(ostringstream &buf, int indent_level)
RooAbsAnaConvPdf is the base class for PDFs that represent a physics model that can be analytically c...
virtual bool changeModel(const RooResolutionModel &newModel)
Change the current resolution model to newModel.
RooAbsGenContext * genContext(const RooArgSet &vars, const RooDataSet *prototype=nullptr, const RooArgSet *auxProto=nullptr, bool verbose=false) const override
Create a generator context for this p.d.f.
RooAbsRealLValue * convVar()
Retrieve the convolution variable.
RooListProxy _convSet
Set of (resModel (x) basisFunc) convolution objects.
RooArgSet * getObservables(const RooArgSet &set, bool valueOnly=true) const
Given a set of possible observables, return the observables that this PDF depends on.
Definition: RooAbsArg.h:293
virtual bool remove(const RooAbsArg &var, bool silent=false, bool matchByNameOnly=false)
Remove the specified argument from our list.
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...
bool selectCommon(const RooAbsCollection &refColl, RooAbsCollection &outColl) const
Create a subset of the current collection, consisting only of those elements that are contained as we...
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 verbose=false, TString indent="") const override
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.
virtual RooAbsGenContext * genContext(const RooArgSet &vars, const RooDataSet *prototype=nullptr, const RooArgSet *auxProto=nullptr, bool verbose=false) const
Interface function to create a generator context from a p.d.f.
Definition: RooAbsPdf.cxx:1864
bool isValidReal(double value, bool printError=false) const override
Check if given value is valid.
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:91
RooAbsArg * at(Int_t idx) const
Return object at given index, or nullptr if index is out of range.
Definition: RooArgList.h:110
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:56
RooArgSet * snapshot(bool deepCopy=true) const
Use RooAbsCollection::snapshot(), but return as RooArgSet.
Definition: RooArgSet.h:179
RooConvGenContext is an efficient implementation of the generator context specific for RooAbsAnaConvP...
std::unique_ptr< RooAbsGenContext > _modelGen
Resolution model generator context.
RooConvGenContext(const RooFFTConvPdf &model, const RooArgSet &vars, const RooDataSet *prototype=nullptr, const RooArgSet *auxProto=nullptr, bool _verbose=false)
Constructor for specialized generator context for FFT numerical convolutions.
std::unique_ptr< RooArgSet > _pdfCloneSet
Owner of PDF clone.
void setProtoDataOrder(Int_t *lut) override
Set the traversal order for events in the prototype dataset The argument is a array of integers with ...
std::unique_ptr< RooArgSet > _pdfVars
Holder of PDF x truth event.
std::unique_ptr< RooArgSet > _modelVars
Holder of resModel event.
RooRealVar * _cvPdf
Convolution variable in PDFxTruth event.
std::unique_ptr< RooArgSet > _modelVarsOwned
Owning version of modelVars ;.
void attach(const RooArgSet &params) override
Attach given set of arguments to internal clones of pdf and resolution model.
std::unique_ptr< RooArgSet > _modelCloneSet
Owner of resModel clone.
RooRealVar * _cvOut
Convolution variable in output event.
RooRealVar * _cvModel
Convolution variable in resModel event.
void generateEvent(RooArgSet &theEvent, Int_t remaining) override
Generate a single event.
void initGenerator(const RooArgSet &theEvent) override
One-time initialization of generator context, attaches the context to the supplied event container.
std::unique_ptr< RooArgSet > _pdfVarsOwned
Owning version of pdfVars ;.
TString _convVarName
Name of convolution variable.
void printMultiline(std::ostream &os, Int_t content, bool verbose=false, TString indent="") const override
Print the details of this generator context.
std::unique_ptr< RooAbsGenContext > _pdfGen
Physics model generator context.
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.
PDF for the numerical (FFT) convolution of two PDFs.
Definition: RooFFTConvPdf.h:25
RooRealProxy _pdf1
First input p.d.f.
Definition: RooFFTConvPdf.h:69
RooRealProxy _x
Convolution observable.
Definition: RooFFTConvPdf.h:67
RooRealProxy _pdf2
Second input p.d.f.
Definition: RooFFTConvPdf.h:70
Numeric 1-dimensional convolution operator PDF.
Definition: RooNumConvPdf.h:26
RooNumConvolution & conv() const
Definition: RooNumConvPdf.h:63
RooAbsReal & cloneModel() const
RooRealVar & cloneVar() const
RooAbsReal & clonePdf() const
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:40
void setVal(double value) override
Set value of variable to 'value'.
Definition: RooRealVar.cxx:254
void removeRange(const char *name=nullptr)
Remove range limits for binning with given name. Empty name means default range.
Definition: RooRealVar.cxx:398
RooResolutionModel is the base class for PDFs that represent a resolution model that can be convolute...
const T & arg() const
Return reference to object held in proxy.
RooTruthModel is an implementation of RooResolution model that provides a delta-function resolution m...
Definition: RooTruthModel.h:21
const char * GetName() const override
Returns name of object.
Definition: TNamed.h:47
Basic string class.
Definition: TString.h:136
TString & Append(const char *cs)
Definition: TString.h:564
@ Generation
Definition: RooGlobalFunc.h:61