Logo ROOT   6.10/09
Reference Guide
rf502_wspacewrite.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 #502
5 ///
6 /// Creating and writing a workspace
7 ///
8 /// \macro_output
9 /// \macro_code
10 /// \author 07/2008 - Wouter Verkerke
11 
12 
13 #include "RooRealVar.h"
14 #include "RooDataSet.h"
15 #include "RooGaussian.h"
16 #include "RooConstVar.h"
17 #include "RooChebychev.h"
18 #include "RooAddPdf.h"
19 #include "RooWorkspace.h"
20 #include "RooPlot.h"
21 #include "TCanvas.h"
22 #include "TAxis.h"
23 #include "TFile.h"
24 #include "TH1.h"
25 using namespace RooFit ;
26 
27 
28 
29 void rf502_wspacewrite()
30 {
31  // C r e a t e m o d e l a n d d a t a s e t
32  // -----------------------------------------------
33 
34  // Declare observable x
35  RooRealVar x("x","x",0,10) ;
36 
37  // Create two Gaussian PDFs g1(x,mean1,sigma) anf g2(x,mean2,sigma) and their parameters
38  RooRealVar mean("mean","mean of gaussians",5,0,10) ;
39  RooRealVar sigma1("sigma1","width of gaussians",0.5) ;
40  RooRealVar sigma2("sigma2","width of gaussians",1) ;
41 
42  RooGaussian sig1("sig1","Signal component 1",x,mean,sigma1) ;
43  RooGaussian sig2("sig2","Signal component 2",x,mean,sigma2) ;
44 
45  // Build Chebychev polynomial p.d.f.
46  RooRealVar a0("a0","a0",0.5,0.,1.) ;
47  RooRealVar a1("a1","a1",0.2,0,1.) ;
48  RooChebychev bkg("bkg","Background",x,RooArgSet(a0,a1)) ;
49 
50  // Sum the signal components into a composite signal p.d.f.
51  RooRealVar sig1frac("sig1frac","fraction of component 1 in signal",0.8,0.,1.) ;
52  RooAddPdf sig("sig","Signal",RooArgList(sig1,sig2),sig1frac) ;
53 
54  // Sum the composite signal and background
55  RooRealVar bkgfrac("bkgfrac","fraction of background",0.5,0.,1.) ;
56  RooAddPdf model("model","g1+g2+a",RooArgList(bkg,sig),bkgfrac) ;
57 
58  // Generate a data sample of 1000 events in x from model
59  RooDataSet *data = model.generate(x,1000) ;
60 
61 
62 
63  // C r e a t e w o r k s p a c e , i m p o r t d a t a a n d m o d e l
64  // -----------------------------------------------------------------------------
65 
66  // Create a new empty workspace
67  RooWorkspace *w = new RooWorkspace("w","workspace") ;
68 
69  // Import model and all its components into the workspace
70  w->import(model) ;
71 
72  // Import data into the workspace
73  w->import(*data) ;
74 
75  // Print workspace contents
76  w->Print() ;
77 
78 
79 
80  // S a v e w o r k s p a c e i n f i l e
81  // -------------------------------------------
82 
83  // Save the workspace into a ROOT file
84  w->writeToFile("rf502_workspace.root") ;
85 
86 
87  // Workspace will remain in memory after macro finishes
88  gDirectory->Add(w) ;
89 
90 }
RooAddPdf is an efficient implementation of a sum of PDFs of the form.
Definition: RooAddPdf.h:29
Double_t x[n]
Definition: legend1.C:17
Plain Gaussian p.d.f.
Definition: RooGaussian.h:25
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:36
Bool_t writeToFile(const char *fileName, Bool_t recreate=kTRUE)
Save this current workspace into given file.
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:29
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.
#define gDirectory
Definition: TDirectory.h:211
Chebychev polynomial p.d.f.
Definition: RooChebychev.h:25
void Print(Option_t *opts=0) const
Print contents of the workspace.
The RooWorkspace is a persistable container for RooFit projects.
Definition: RooWorkspace.h:42