Logo ROOT   6.10/09
Reference Guide
rf316_llratioplot.cxx
Go to the documentation of this file.
1 //////////////////////////////////////////////////////////////////////////
2 //
3 // 'MULTIDIMENSIONAL MODELS' RooFit tutorial macro #316
4 //
5 // Using the likelihood ratio techique to construct a signal enhanced
6 // one-dimensional projection of a multi-dimensional p.d.f.
7 //
8 //
9 //
10 // 07/2008 - Wouter Verkerke
11 //
12 /////////////////////////////////////////////////////////////////////////
13 
14 #ifndef __CINT__
15 #include "RooGlobalFunc.h"
16 #endif
17 #include "RooRealVar.h"
18 #include "RooDataSet.h"
19 #include "RooGaussian.h"
20 #include "RooPolynomial.h"
21 #include "RooAddPdf.h"
22 #include "RooProdPdf.h"
23 #include "TCanvas.h"
24 #include "RooPlot.h"
25 using namespace RooFit ;
26 
27 
28 class TestBasic316 : public RooFitTestUnit
29 {
30 public:
31  TestBasic316(TFile* refFile, Bool_t writeRef, Int_t verbose) : RooFitTestUnit("Likelihood ratio projection plot",refFile,writeRef,verbose) {} ;
32  Bool_t testCode() {
33 
34  // C r e a t e 3 D p d f a n d d a t a
35  // -------------------------------------------
36 
37  // Create observables
38  RooRealVar x("x","x",-5,5) ;
39  RooRealVar y("y","y",-5,5) ;
40  RooRealVar z("z","z",-5,5) ;
41 
42  // Create signal pdf gauss(x)*gauss(y)*gauss(z)
43  RooGaussian gx("gx","gx",x,RooConst(0),RooConst(1)) ;
44  RooGaussian gy("gy","gy",y,RooConst(0),RooConst(1)) ;
45  RooGaussian gz("gz","gz",z,RooConst(0),RooConst(1)) ;
46  RooProdPdf sig("sig","sig",RooArgSet(gx,gy,gz)) ;
47 
48  // Create background pdf poly(x)*poly(y)*poly(z)
49  RooPolynomial px("px","px",x,RooArgSet(RooConst(-0.1),RooConst(0.004))) ;
50  RooPolynomial py("py","py",y,RooArgSet(RooConst(0.1),RooConst(-0.004))) ;
51  RooPolynomial pz("pz","pz",z) ;
52  RooProdPdf bkg("bkg","bkg",RooArgSet(px,py,pz)) ;
53 
54  // Create composite pdf sig+bkg
55  RooRealVar fsig("fsig","signal fraction",0.1,0.,1.) ;
56  RooAddPdf model("model","model",RooArgList(sig,bkg),fsig) ;
57 
58  RooDataSet* data = model.generate(RooArgSet(x,y,z),20000) ;
59 
60 
61 
62  // P r o j e c t p d f a n d d a t a o n x
63  // -------------------------------------------------
64 
65  // Make plain projection of data and pdf on x observable
66  RooPlot* frame = x.frame(Title("Projection of 3D data and pdf on X"),Bins(40)) ;
67  data->plotOn(frame) ;
68  model.plotOn(frame) ;
69 
70 
71 
72  // D e f i n e p r o j e c t e d s i g n a l l i k e l i h o o d r a t i o
73  // ----------------------------------------------------------------------------------
74 
75  // Calculate projection of signal and total likelihood on (y,z) observables
76  // i.e. integrate signal and composite model over x
77  RooAbsPdf* sigyz = sig.createProjection(x) ;
78  RooAbsPdf* totyz = model.createProjection(x) ;
79 
80  // Construct the log of the signal / signal+background probability
81  RooFormulaVar llratio_func("llratio","log10(@0)-log10(@1)",RooArgList(*sigyz,*totyz)) ;
82 
83 
84 
85  // P l o t d a t a w i t h a L L r a t i o c u t
86  // -------------------------------------------------------
87 
88  // Calculate the llratio value for each event in the dataset
89  data->addColumn(llratio_func) ;
90 
91  // Extract the subset of data with large signal likelihood
92  RooDataSet* dataSel = (RooDataSet*) data->reduce(Cut("llratio>0.7")) ;
93 
94  // Make plot frame
95  RooPlot* frame2 = x.frame(Title("Same projection on X with LLratio(y,z)>0.7"),Bins(40)) ;
96 
97  // Plot select data on frame
98  dataSel->plotOn(frame2) ;
99 
100 
101 
102  // M a k e M C p r o j e c t i o n o f p d f w i t h s a m e L L r a t i o c u t
103  // ---------------------------------------------------------------------------------------------
104 
105  // Generate large number of events for MC integration of pdf projection
106  RooDataSet* mcprojData = model.generate(RooArgSet(x,y,z),10000) ;
107 
108  // Calculate LL ratio for each generated event and select MC events with llratio)0.7
109  mcprojData->addColumn(llratio_func) ;
110  RooDataSet* mcprojDataSel = (RooDataSet*) mcprojData->reduce(Cut("llratio>0.7")) ;
111 
112  // Project model on x, integrating projected observables (y,z) with Monte Carlo technique
113  // on set of events with the same llratio cut as was applied to data
114  model.plotOn(frame2,ProjWData(*mcprojDataSel)) ;
115 
116 
117  regPlot(frame,"rf316_plot1") ;
118  regPlot(frame2,"rf316_plot2") ;
119 
120  delete data ;
121  delete dataSel ;
122  delete mcprojData ;
123  delete mcprojDataSel ;
124  delete sigyz ;
125  delete totyz ;
126 
127  return kTRUE ;
128  }
129 } ;
RooAddPdf is an efficient implementation of a sum of PDFs of the form.
Definition: RooAddPdf.h:29
virtual RooPlot * plotOn(RooPlot *frame, const RooCmdArg &arg1=RooCmdArg::none(), const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
Plot dataset on specified frame.
Definition: RooAbsData.cxx:552
RooCmdArg Cut(const char *cutSpec)
RooProdPdf is an efficient implementation of a product of PDFs of the form.
Definition: RooProdPdf.h:31
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
RooAbsData * reduce(const RooCmdArg &arg1, 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())
Create a reduced copy of this dataset.
Definition: RooAbsData.cxx:343
RooCmdArg ProjWData(const RooAbsData &projData, Bool_t binData=kFALSE)
RooCmdArg Title(const char *name)
Double_t x[n]
Definition: legend1.C:17
virtual RooAbsArg * addColumn(RooAbsArg &var, Bool_t adjustRange=kTRUE)
Add a column with the values of the given (function) argument to this dataset.
Plain Gaussian p.d.f.
Definition: RooGaussian.h:25
virtual RooAbsPdf * createProjection(const RooArgSet &iset)
Return a p.d.f that represent a projection of this p.d.f integrated over given observables.
Definition: RooAbsPdf.cxx:2989
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:36
bool verbose
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:29
A RooPlot is a plot frame and a container for graphics objects within that frame. ...
Definition: RooPlot.h:41
Double_t y[n]
Definition: legend1.C:17
you should not use this method at all Int_t Int_t z
Definition: TRolke.cxx:630
RooConstVar & RooConst(Double_t val)
RooAbsPdf is the abstract interface for all probability density functions The class provides hybrid a...
Definition: RooAbsPdf.h:41
RooCmdArg Bins(Int_t nbin)
RooPolynomial implements a polynomial p.d.f of the form By default coefficient a_0 is chosen to be 1...
Definition: RooPolynomial.h:28
const Bool_t kTRUE
Definition: RtypesCore.h:91