Logo ROOT   6.14/05
Reference Guide
rf511_wsfactory_basic.C
Go to the documentation of this file.
1 /// \file
2 /// \ingroup tutorial_roofit
3 /// \notebook -nodraw
4 /// 'ORGANIZATION AND SIMULTANEOUS FITS' RooFit tutorial macro #511
5 ///
6 /// Basic use of the 'object factory' associated with a workspace
7 /// to rapidly build p.d.f.s functions and their parameter components
8 ///
9 /// \macro_output
10 /// \macro_code
11 /// \author 04/2009 - Wouter Verkerke
12 
13 
14 
15 #include "RooRealVar.h"
16 #include "RooDataSet.h"
17 #include "RooGaussian.h"
18 #include "RooConstVar.h"
19 #include "RooChebychev.h"
20 #include "RooAddPdf.h"
21 #include "RooWorkspace.h"
22 #include "RooPlot.h"
23 #include "TCanvas.h"
24 #include "TAxis.h"
25 using namespace RooFit ;
26 
27 
28 void rf511_wsfactory_basic(Bool_t compact=kFALSE)
29 {
30  RooWorkspace* w = new RooWorkspace("w") ;
31 
32  // C r e a t i n g a n d a d d i n g b a s i c p . d . f . s
33  // ----------------------------------------------------------------
34 
35 
36  // Remake example p.d.f. of tutorial rs502_wspacewrite.C:
37  //
38  // Basic p.d.f. construction: ClassName::ObjectName(constructor arguments)
39  // Variable construction : VarName[x,xlo,xhi], VarName[xlo,xhi], VarName[x]
40  // P.d.f. addition : SUM::ObjectName(coef1*pdf1,...coefM*pdfM,pdfN)
41  //
42 
43  if (!compact) {
44 
45  // Use object factory to build p.d.f. of tutorial rs502_wspacewrite
46  w->factory("Gaussian::sig1(x[-10,10],mean[5,0,10],0.5)") ;
47  w->factory("Gaussian::sig2(x,mean,1)") ;
48  w->factory("Chebychev::bkg(x,{a0[0.5,0.,1],a1[0.2,0.,1.]})") ;
49  w->factory("SUM::sig(sig1frac[0.8,0.,1.]*sig1,sig2)") ;
50  w->factory("SUM::model(bkgfrac[0.5,0.,1.]*bkg,sig)") ;
51 
52  } else {
53 
54  // Use object factory to build p.d.f. of tutorial rs502_wspacewrite but
55  // - Contracted to a single line recursive expression,
56  // - Omitting explicit names for components that are not referred to explicitly later
57 
58  w->factory("SUM::model(bkgfrac[0.5,0.,1.]*Chebychev::bkg(x[-10,10],{a0[0.5,0.,1],a1[0.2,0.,1.]}),"
59  "SUM(sig1frac[0.8,0.,1.]*Gaussian(x,mean[5,0,10],0.5), Gaussian(x,mean,1)))") ;
60  }
61 
62 
63  // A d v a n c e d p . d . f . c o n s t r u c t o r a r g u m e n t s
64  // ----------------------------------------------------------------
65  //
66  // P.d.f. constructor arguments may by any type of RooAbsArg, but also
67  //
68  // Double_t --> converted to RooConst(...)
69  // {a,b,c} --> converted to RooArgSet() or RooArgList() depending on required ctor arg
70  // dataset name --> converted to RooAbsData reference for any dataset residing in the workspace
71  // enum --> Any enum label that belongs to an enum defined in the (base) class
72 
73 
74  // Make a dummy dataset p.d.f. 'model' and import it in the workspace
75  RooDataSet* data = w->pdf("model")->generate(*w->var("x"),1000) ;
76  w->import(*data,Rename("data")) ;
77 
78  // Construct a KEYS p.d.f. passing a dataset name and an enum type defining the
79  // mirroring strategy
80  w->factory("KeysPdf::k(x,data,NoMirror,0.2)") ;
81 
82 
83  // Print workspace contents
84  w->Print() ;
85 
86  // Make workspace visible on command line
87  gDirectory->Add(w) ;
88 
89 }
90 
91 
bool Bool_t
Definition: RtypesCore.h:59
RooCmdArg Rename(const char *suffix)
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:29
const Bool_t kFALSE
Definition: RtypesCore.h:88
RooAbsPdf * pdf(const char *name) const
Retrieve p.d.f (RooAbsPdf) with given name. A null pointer is returned if not found.
RooRealVar * var(const char *name) const
Retrieve real-valued variable (RooRealVar) with given name. A null pointer is returned if not found...
RooFactoryWSTool & factory()
Return instance to factory tool.
Bool_t import(const RooAbsArg &arg, const RooCmdArg &arg1=RooCmdArg(), const RooCmdArg &arg2=RooCmdArg(), const RooCmdArg &arg3=RooCmdArg(), const RooCmdArg &arg4=RooCmdArg(), const RooCmdArg &arg5=RooCmdArg(), const RooCmdArg &arg6=RooCmdArg(), const RooCmdArg &arg7=RooCmdArg(), const RooCmdArg &arg8=RooCmdArg(), const RooCmdArg &arg9=RooCmdArg())
Import a RooAbsArg object, e.g.
RooDataSet * generate(const RooArgSet &whatVars, Int_t nEvents, const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none())
Generate a new dataset containing the specified variables with events sampled from our distribution...
Definition: RooAbsPdf.cxx:1725
#define gDirectory
Definition: TDirectory.h:213
void Print(Option_t *opts=0) const
Print contents of the workspace.
The RooWorkspace is a persistable container for RooFit projects.
Definition: RooWorkspace.h:43