Logo ROOT   6.08/07
Reference Guide
RooRandomizeParamMCSModule.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 RooRandomizeParamMCSModule.cxx
19 \class RooRandomizeParamMCSModule
20 \ingroup Roofitcore
21 
22 RooRandomizeParamMCSModule is an add-on modules to RooMCStudy that
23 allows you to randomize input generation parameters. Randomized generation
24 parameters can be sampled from a uniform or Gaussian distribution.
25 For every randomized parameter, an extra variable is added to
26 RooMCStudy::fitParDataSet() named <tt><parname>_gen</tt> that indicates the actual
27 value used for generation for each trial.
28 You can also choose to randomize the sum of N parameters, rather
29 than a single parameter. In that case common multiplicative scale
30 factor is applied to each component to bring the sum to the desired
31 target value taken from either uniform or Gaussian sampling. This
32 latter option is for example useful if you want to change the total
33 number of expected events of an extended p.d.f
34 **/
35 
36 
37 #include "Riostream.h"
38 #include "RooDataSet.h"
39 #include "RooRealVar.h"
40 #include "RooRandom.h"
41 #include "TString.h"
42 #include "RooFit.h"
43 #include "RooFitResult.h"
44 #include "RooAddition.h"
45 #include "RooMsgService.h"
47 
48 using namespace std ;
49 
51  ;
52 
53 
54 
55 ////////////////////////////////////////////////////////////////////////////////
56 /// Constructor
57 
59  RooAbsMCStudyModule("RooRandomizeParamMCSModule","RooRandomizeParamMCSModule"), _data(0)
60 {
61 }
62 
63 
64 
65 ////////////////////////////////////////////////////////////////////////////////
66 /// Copy constructor
67 
69  RooAbsMCStudyModule(other),
70  _unifParams(other._unifParams),
71  _gausParams(other._gausParams),
72  _data(0)
73 {
74 }
75 
76 
77 
78 ////////////////////////////////////////////////////////////////////////////////
79 /// Destructor
80 
82 {
83  if (_data) {
84  delete _data ;
85  }
86 }
87 
88 
89 
90 ////////////////////////////////////////////////////////////////////////////////
91 /// Request uniform smearing of param in range [lo,hi] in RooMCStudy
92 /// generation cycle
93 
95 {
96  // If we're already attached to a RooMCStudy, check that given param is actual generator model parameter
97  // If not attached, this check is repeated at the attachment moment
98  if (genParams()) {
99  RooRealVar* actualPar = static_cast<RooRealVar*>(genParams()->find(param.GetName())) ;
100  if (!actualPar) {
101  oocoutW((TObject*)0,InputArguments) << "RooRandomizeParamMCSModule::initializeInstance: variable " << param.GetName() << " is not a parameter of RooMCStudy model and is ignored!" << endl ;
102  return ;
103  }
104  }
105 
106  _unifParams.push_back(UniParam(&param,lo,hi)) ;
107 }
108 
109 
110 
111 ////////////////////////////////////////////////////////////////////////////////
112 /// Request Gaussian smearing of param in with mean 'mean' and width
113 /// 'sigma' in RooMCStudy generation cycle
114 
116 {
117  // If we're already attached to a RooMCStudy, check that given param is actual generator model parameter
118  // If not attached, this check is repeated at the attachment moment
119  if (genParams()) {
120  RooRealVar* actualPar = static_cast<RooRealVar*>(genParams()->find(param.GetName())) ;
121  if (!actualPar) {
122  oocoutW((TObject*)0,InputArguments) << "RooRandomizeParamMCSModule::initializeInstance: variable " << param.GetName() << " is not a parameter of RooMCStudy model and is ignored!" << endl ;
123  return ;
124  }
125  }
126 
127  _gausParams.push_back(GausParam(&param,mean,sigma)) ;
128 }
129 
130 
131 
132 
133 ////////////////////////////////////////////////////////////////////////////////
134 /// Request uniform smearing of sum of parameters in paramSet uniform
135 /// smearing in range [lo,hi] in RooMCStudy generation cycle. This
136 /// option applies a common multiplicative factor to each parameter
137 /// in paramSet to make the sum of the parameters add up to the
138 /// sampled value in the range [lo,hi]
139 
141 {
142  // Check that all args are RooRealVars
143  RooArgSet okset ;
144  TIterator* iter = paramSet.createIterator() ;
145  RooAbsArg* arg ;
146  while((arg=(RooAbsArg*)iter->Next())) {
147  // Check that arg is a RooRealVar
148  RooRealVar* rrv = dynamic_cast<RooRealVar*>(arg) ;
149  if (!rrv) {
150  oocoutW((TObject*)0,InputArguments) << "RooRandomizeParamMCSModule::sampleSumUniform() ERROR: input parameter " << arg->GetName() << " is not a RooRealVar and is ignored" << endl ;
151  continue;
152  }
153  okset.add(*rrv) ;
154  }
155  delete iter ;
156 
157  // If we're already attached to a RooMCStudy, check that given param is actual generator model parameter
158  // If not attached, this check is repeated at the attachment moment
159  RooArgSet okset2 ;
160  if (genParams()) {
161  TIterator* psiter = okset.createIterator() ;
162  RooAbsArg* arg2 ;
163  while ((arg2=(RooAbsArg*)psiter->Next())) {
164  RooRealVar* actualVar= static_cast<RooRealVar*>(genParams()->find(arg2->GetName())) ;
165  if (!actualVar) {
166  oocoutW((TObject*)0,InputArguments) << "RooRandomizeParamMCSModule::sampleSumUniform: variable " << arg2->GetName() << " is not a parameter of RooMCStudy model and is ignored!" << endl ;
167  } else {
168  okset2.add(*actualVar) ;
169  }
170  }
171  delete psiter ;
172  } else {
173 
174  // If genParams() are not available, skip this check for now
175  okset2.add(okset) ;
176 
177  }
178 
179 
180  _unifParamSets.push_back(UniParamSet(okset2,lo,hi)) ;
181 
182 }
183 
184 
185 
186 
187 ////////////////////////////////////////////////////////////////////////////////
188 /// Request gaussian smearing of sum of parameters in paramSet
189 /// uniform smearing with mean 'mean' and width 'sigma' in RooMCStudy
190 /// generation cycle. This option applies a common multiplicative
191 /// factor to each parameter in paramSet to make the sum of the
192 /// parameters add up to the sampled value from the
193 /// gaussian(mean,sigma)
194 
196 {
197  // Check that all args are RooRealVars
198  RooArgSet okset ;
199  TIterator* iter = paramSet.createIterator() ;
200  RooAbsArg* arg ;
201  while((arg=(RooAbsArg*)iter->Next())) {
202  // Check that arg is a RooRealVar
203  RooRealVar* rrv = dynamic_cast<RooRealVar*>(arg) ;
204  if (!rrv) {
205  oocoutW((TObject*)0,InputArguments) << "RooRandomizeParamMCSModule::sampleSumGauss() ERROR: input parameter " << arg->GetName() << " is not a RooRealVar and is ignored" << endl ;
206  continue;
207  }
208  okset.add(*rrv) ;
209  }
210 
211  // If we're already attached to a RooMCStudy, check that given param is actual generator model parameter
212  // If not attached, this check is repeated at the attachment moment
213  RooArgSet okset2 ;
214  if (genParams()) {
215  TIterator* psiter = okset.createIterator() ;
216  RooAbsArg* arg2 ;
217  while ((arg2=(RooAbsArg*)psiter->Next())) {
218  RooRealVar* actualVar= static_cast<RooRealVar*>(genParams()->find(arg2->GetName())) ;
219  if (!actualVar) {
220  oocoutW((TObject*)0,InputArguments) << "RooRandomizeParamMCSModule::sampleSumUniform: variable " << arg2->GetName() << " is not a parameter of RooMCStudy model and is ignored!" << endl ;
221  } else {
222  okset2.add(*actualVar) ;
223  }
224  }
225  delete psiter ;
226  } else {
227 
228  // If genParams() are not available, skip this check for now
229  okset2.add(okset) ;
230 
231  }
232 
233  _gausParamSets.push_back(GausParamSet(okset,mean,sigma)) ;
234 
235 }
236 
237 
238 
239 
240 ////////////////////////////////////////////////////////////////////////////////
241 /// Initialize module after attachment to RooMCStudy object
242 
244 {
245  // Loop over all uniform smearing parameters
246  std::list<UniParam>::iterator uiter ;
247  for (uiter= _unifParams.begin() ; uiter!= _unifParams.end() ; ++uiter) {
248 
249  // Check that listed variable is actual generator model parameter
250  RooRealVar* actualPar = static_cast<RooRealVar*>(genParams()->find(uiter->_param->GetName())) ;
251  if (!actualPar) {
252  oocoutW((TObject*)0,InputArguments) << "RooRandomizeParamMCSModule::initializeInstance: variable " << uiter->_param->GetName() << " is not a parameter of RooMCStudy model and is ignored!" << endl ;
253  uiter = _unifParams.erase(uiter) ;
254  continue ;
255  }
256  uiter->_param = actualPar ;
257 
258  // Add variable to summary dataset to hold generator value
259  TString parName = Form("%s_gen",uiter->_param->GetName()) ;
260  TString parTitle = Form("%s as generated",uiter->_param->GetTitle()) ;
261  RooRealVar* par_gen = new RooRealVar(parName.Data(),parTitle.Data(),0) ;
262  _genParSet.addOwned(*par_gen) ;
263  }
264 
265  // Loop over all gaussian smearing parameters
266  std::list<UniParam>::iterator giter ;
267  for (giter= _unifParams.begin() ; giter!= _unifParams.end() ; ++giter) {
268 
269  // Check that listed variable is actual generator model parameter
270  RooRealVar* actualPar = static_cast<RooRealVar*>(genParams()->find(giter->_param->GetName())) ;
271  if (!actualPar) {
272  oocoutW((TObject*)0,InputArguments) << "RooRandomizeParamMCSModule::initializeInstance: variable " << giter->_param->GetName() << " is not a parameter of RooMCStudy model and is ignored!" << endl ;
273  giter = _unifParams.erase(giter) ;
274  continue ;
275  }
276  giter->_param = actualPar ;
277 
278  // Add variable to summary dataset to hold generator value
279  TString parName = Form("%s_gen",giter->_param->GetName()) ;
280  TString parTitle = Form("%s as generated",giter->_param->GetTitle()) ;
281  RooRealVar* par_gen = new RooRealVar(parName.Data(),parTitle.Data(),0) ;
282  _genParSet.addOwned(*par_gen) ;
283  }
284 
285 
286  // Loop over all uniform smearing set of parameters
287  std::list<UniParamSet>::iterator usiter ;
288  for (usiter= _unifParamSets.begin() ; usiter!= _unifParamSets.end() ; ++usiter) {
289 
290  // Check that all listed variables are actual generator model parameters
291  RooArgSet actualPSet ;
292  TIterator* psiter = usiter->_pset.createIterator() ;
293  RooAbsArg* arg ;
294  while ((arg=(RooAbsArg*)psiter->Next())) {
295  RooRealVar* actualVar= static_cast<RooRealVar*>(genParams()->find(arg->GetName())) ;
296  if (!actualVar) {
297  oocoutW((TObject*)0,InputArguments) << "RooRandomizeParamMCSModule::initializeInstance: variable " << arg->GetName() << " is not a parameter of RooMCStudy model and is ignored!" << endl ;
298  } else {
299  actualPSet.add(*actualVar) ;
300  }
301  }
302  delete psiter ;
303  usiter->_pset.removeAll() ;
304  usiter->_pset.add(actualPSet) ;
305 
306  // Add variables to summary dataset to hold generator values
307  TIterator* iter = usiter->_pset.createIterator() ;
308  RooRealVar* param ;
309  while((param=(RooRealVar*)iter->Next())) {
310  TString parName = Form("%s_gen",param->GetName()) ;
311  TString parTitle = Form("%s as generated",param->GetTitle()) ;
312  RooRealVar* par_gen = new RooRealVar(parName.Data(),parTitle.Data(),0) ;
313  _genParSet.addOwned(*par_gen) ;
314  }
315  delete iter ;
316  }
317 
318  // Loop over all gaussian smearing set of parameters
319  std::list<GausParamSet>::iterator ugiter ;
320  for (ugiter= _gausParamSets.begin() ; ugiter!= _gausParamSets.end() ; ++ugiter) {
321 
322  // Check that all listed variables are actual generator model parameters
323  RooArgSet actualPSet ;
324  TIterator* psiter = ugiter->_pset.createIterator() ;
325  RooAbsArg* arg ;
326  while ((arg=(RooAbsArg*)psiter->Next())) {
327  RooRealVar* actualVar= static_cast<RooRealVar*>(genParams()->find(arg->GetName())) ;
328  if (!actualVar) {
329  oocoutW((TObject*)0,InputArguments) << "RooRandomizeParamMCSModule::initializeInstance: variable " << arg->GetName() << " is not a parameter of RooMCStudy model and is ignored!" << endl ;
330  } else {
331  actualPSet.add(*actualVar) ;
332  }
333  }
334  ugiter->_pset.removeAll() ;
335  ugiter->_pset.add(actualPSet) ;
336 
337  // Add variables to summary dataset to hold generator values
338  TIterator* iter = ugiter->_pset.createIterator() ;
339  RooRealVar* param ;
340  while((param=(RooRealVar*)iter->Next())) {
341  TString parName = Form("%s_gen",param->GetName()) ;
342  TString parTitle = Form("%s as generated",param->GetTitle()) ;
343  RooRealVar* par_gen = new RooRealVar(parName.Data(),parTitle.Data(),0) ;
344  _genParSet.addOwned(*par_gen) ;
345  }
346  }
347 
348  // Create new dataset to be merged with RooMCStudy::fitParDataSet
349  _data = new RooDataSet("DeltaLLSigData","Additional data for Delta(-log(L)) study",_genParSet) ;
350 
351  return kTRUE ;
352 }
353 
354 
355 
356 ////////////////////////////////////////////////////////////////////////////////
357 /// Initialize module at beginning of RooCMStudy run
358 
360 {
361  // Clear dataset at beginning of run
362  _data->reset() ;
363  return kTRUE ;
364 }
365 
366 
367 
368 ////////////////////////////////////////////////////////////////////////////////
369 /// Apply all smearings to generator parameters
370 
372 {
373  // Apply uniform smearing to all generator parameters for which it is requested
374  std::list<UniParam>::iterator uiter ;
375  for (uiter= _unifParams.begin() ; uiter!= _unifParams.end() ; ++uiter) {
376  Double_t newVal = RooRandom::randomGenerator()->Uniform(uiter->_lo,uiter->_hi) ;
377  oocoutE((TObject*)0,Generation) << "RooRandomizeParamMCSModule::processBeforeGen: applying uniform smearing to generator parameter "
378  << uiter->_param->GetName() << " in range [" << uiter->_lo << "," << uiter->_hi << "], chosen value for this sample is " << newVal << endl ;
379  uiter->_param->setVal(newVal) ;
380 
381  RooRealVar* genpar = static_cast<RooRealVar*>(_genParSet.find(Form("%s_gen",uiter->_param->GetName()))) ;
382  genpar->setVal(newVal) ;
383  }
384 
385  // Apply gaussian smearing to all generator parameters for which it is requested
386  std::list<GausParam>::iterator giter ;
387  for (giter= _gausParams.begin() ; giter!= _gausParams.end() ; ++giter) {
388  Double_t newVal = RooRandom::randomGenerator()->Gaus(giter->_mean,giter->_sigma) ;
389  oocoutI((TObject*)0,Generation) << "RooRandomizeParamMCSModule::processBeforeGen: applying gaussian smearing to generator parameter "
390  << giter->_param->GetName() << " with a mean of " << giter->_mean << " and a width of " << giter->_sigma << ", chosen value for this sample is " << newVal << endl ;
391  giter->_param->setVal(newVal) ;
392 
393  RooRealVar* genpar = static_cast<RooRealVar*>(_genParSet.find(Form("%s_gen",giter->_param->GetName()))) ;
394  genpar->setVal(newVal) ;
395  }
396 
397  // Apply uniform smearing to all sets of generator parameters for which it is requested
398  std::list<UniParamSet>::iterator usiter ;
399  for (usiter= _unifParamSets.begin() ; usiter!= _unifParamSets.end() ; ++usiter) {
400 
401  // Calculate new value for sum
402  Double_t newVal = RooRandom::randomGenerator()->Uniform(usiter->_lo,usiter->_hi) ;
403  oocoutI((TObject*)0,Generation) << "RooRandomizeParamMCSModule::processBeforeGen: applying uniform smearing to sum of set of generator parameters "
404  << usiter->_pset
405  << " in range [" << usiter->_lo << "," << usiter->_hi << "], chosen sum value for this sample is " << newVal << endl ;
406 
407  // Determine original value of sum and calculate per-component scale factor to obtain new valye for sum
408  RooAddition sumVal("sumVal","sumVal",usiter->_pset) ;
409  Double_t compScaleFactor = newVal/sumVal.getVal() ;
410 
411  // Apply multiplicative correction to each term of the sum
412  TIterator* iter = usiter->_pset.createIterator() ;
413  RooRealVar* param ;
414  while((param=(RooRealVar*)iter->Next())) {
415  param->setVal(param->getVal()*compScaleFactor) ;
416  RooRealVar* genpar = static_cast<RooRealVar*>(_genParSet.find(Form("%s_gen",param->GetName()))) ;
417  genpar->setVal(param->getVal()) ;
418  }
419  delete iter ;
420  }
421 
422  // Apply gaussian smearing to all sets of generator parameters for which it is requested
423  std::list<GausParamSet>::iterator gsiter ;
424  for (gsiter= _gausParamSets.begin() ; gsiter!= _gausParamSets.end() ; ++gsiter) {
425 
426  // Calculate new value for sum
427  Double_t newVal = RooRandom::randomGenerator()->Gaus(gsiter->_mean,gsiter->_sigma) ;
428  oocoutI((TObject*)0,Generation) << "RooRandomizeParamMCSModule::processBeforeGen: applying gaussian smearing to sum of set of generator parameters "
429  << gsiter->_pset
430  << " with a mean of " << gsiter->_mean << " and a width of " << gsiter->_sigma
431  << ", chosen value for this sample is " << newVal << endl ;
432 
433  // Determine original value of sum and calculate per-component scale factor to obtain new valye for sum
434  RooAddition sumVal("sumVal","sumVal",gsiter->_pset) ;
435  Double_t compScaleFactor = newVal/sumVal.getVal() ;
436 
437  // Apply multiplicative correction to each term of the sum
438  TIterator* iter = gsiter->_pset.createIterator() ;
439  RooRealVar* param ;
440  while((param=(RooRealVar*)iter->Next())) {
441  param->setVal(param->getVal()*compScaleFactor) ;
442  RooRealVar* genpar = static_cast<RooRealVar*>(_genParSet.find(Form("%s_gen",param->GetName()))) ;
443  genpar->setVal(param->getVal()) ;
444  }
445  }
446 
447  // Store generator values for all modified parameters
448  _data->add(_genParSet) ;
449 
450  return kTRUE ;
451 }
452 
453 
454 
455 ////////////////////////////////////////////////////////////////////////////////
456 /// Return auxiliary data of this module so that it is merged with
457 /// RooMCStudy::fitParDataSet()
458 
460 {
461  return _data ;
462 }
463 
464 
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:51
void sampleSumUniform(const RooArgSet &paramSet, Double_t lo, Double_t hi)
Request uniform smearing of sum of parameters in paramSet uniform smearing in range [lo...
TIterator * createIterator(Bool_t dir=kIterForward) const
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:86
RooAbsMCStudyModule is a base class for add-on modules to RooMCStudy that can perform additional calc...
Bool_t initializeInstance()
Initialize module after attachment to RooMCStudy object.
void sampleGaussian(RooRealVar &param, Double_t mean, Double_t sigma)
Request Gaussian smearing of param in with mean &#39;mean&#39; and width &#39;sigma&#39; in RooMCStudy generation cyc...
virtual Double_t Gaus(Double_t mean=0, Double_t sigma=1)
Samples a random number from the standard Normal (Gaussian) Distribution with the given mean and sigm...
Definition: TRandom.cxx:235
Double_t getVal(const RooArgSet *set=0) const
Definition: RooAbsReal.h:64
#define oocoutI(o, a)
Definition: RooMsgService.h:45
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
STL namespace.
Iterator abstract base class.
Definition: TIterator.h:32
Bool_t initializeRun(Int_t)
Initialize module at beginning of RooCMStudy run.
virtual void reset()
Definition: RooAbsData.cxx:276
#define oocoutE(o, a)
Definition: RooMsgService.h:48
if on multiple lines(like in C++). **The " * configuration fragment. * * The "import myobject continue
Parses the configuration file.
Definition: HLFactory.cxx:368
static TRandom * randomGenerator()
Return a pointer to a singleton random-number generator implementation.
Definition: RooRandom.cxx:54
Bool_t processBeforeGen(Int_t)
Apply all smearings to generator parameters.
const Double_t sigma
std::list< UniParamSet > _unifParamSets
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:37
std::list< GausParamSet > _gausParamSets
virtual void setVal(Double_t value)
Set value of variable to &#39;value&#39;.
Definition: RooRealVar.cxx:205
virtual ~RooRandomizeParamMCSModule()
Destructor.
RooRandomizeParamMCSModule is an add-on modules to RooMCStudy that allows you to randomize input gene...
char * Form(const char *fmt,...)
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:90
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:29
virtual void add(const RooArgSet &row, Double_t weight=1.0, Double_t weightError=0)
Add a data point, with its coordinates specified in the &#39;data&#39; argset, to the data set...
RooDataSet * finalizeRun()
Return auxiliary data of this module so that it is merged with RooMCStudy::fitParDataSet() ...
void sampleUniform(RooRealVar &param, Double_t lo, Double_t hi)
Request uniform smearing of param in range [lo,hi] in RooMCStudy generation cycle.
#define ClassImp(name)
Definition: Rtypes.h:279
RooAbsArg * find(const char *name) const
Find object with given name in list.
double Double_t
Definition: RtypesCore.h:55
virtual Double_t Uniform(Double_t x1=1)
Returns a uniform deviate on the interval (0, x1).
Definition: TRandom.cxx:606
#define oocoutW(o, a)
Definition: RooMsgService.h:47
Mother of all ROOT objects.
Definition: TObject.h:37
virtual TObject * Next()=0
float type_of_call hi(const int &, const int &)
RooAddition calculates the sum of a set of RooAbsReal terms, or when constructed with two sets...
Definition: RooAddition.h:26
RooAbsArg is the common abstract base class for objects that represent a value (of arbitrary type) an...
Definition: RooAbsArg.h:66
void sampleSumGauss(const RooArgSet &paramSet, Double_t lo, Double_t hi)
Request gaussian smearing of sum of parameters in paramSet uniform smearing with mean &#39;mean&#39; and widt...
const Bool_t kTRUE
Definition: Rtypes.h:91
return
Definition: HLFactory.cxx:514
virtual const char * GetTitle() const
Returns title of object.
Definition: TNamed.h:52