Logo ROOT   6.08/07
Reference Guide
rf312_multirangefit.C
Go to the documentation of this file.
1 /// \file
2 /// \ingroup tutorial_roofit
3 /// \notebook -nodraw
4 /// 'MULTIDIMENSIONAL MODELS' RooFit tutorial macro #312
5 ///
6 /// Performing fits in multiple (disjoint) ranges in one or more dimensions
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 "RooProdPdf.h"
18 #include "RooAddPdf.h"
19 #include "RooPolynomial.h"
20 #include "TCanvas.h"
21 #include "TAxis.h"
22 #include "RooPlot.h"
23 #include "RooFitResult.h"
24 using namespace RooFit ;
25 
26 
27 void rf312_multirangefit()
28 {
29 
30  // C r e a t e 2 D p d f a n d d a t a
31  // -------------------------------------------
32 
33  // Define observables x,y
34  RooRealVar x("x","x",-10,10) ;
35  RooRealVar y("y","y",-10,10) ;
36 
37  // Construct the signal pdf gauss(x)*gauss(y)
38  RooRealVar mx("mx","mx",1,-10,10) ;
39  RooRealVar my("my","my",1,-10,10) ;
40 
41  RooGaussian gx("gx","gx",x,mx,RooConst(1)) ;
42  RooGaussian gy("gy","gy",y,my,RooConst(1)) ;
43 
44  RooProdPdf sig("sig","sig",gx,gy) ;
45 
46  // Construct the background pdf (flat in x,y)
47  RooPolynomial px("px","px",x) ;
48  RooPolynomial py("py","py",y) ;
49  RooProdPdf bkg("bkg","bkg",px,py) ;
50 
51  // Construct the composite model sig+bkg
52  RooRealVar f("f","f",0.,1.) ;
53  RooAddPdf model("model","model",RooArgList(sig,bkg),f) ;
54 
55  // Sample 10000 events in (x,y) from the model
56  RooDataSet* modelData = model.generate(RooArgSet(x,y),10000) ;
57 
58 
59 
60  // D e f i n e s i g n a l a n d s i d e b a n d r e g i o n s
61  // -------------------------------------------------------------------
62 
63  // Construct the SideBand1,SideBand2,Signal regions
64  //
65  // |
66  // +-------------+-----------+
67  // | | |
68  // | Side | Sig |
69  // | Band1 | nal |
70  // | | |
71  // --+-------------+-----------+--
72  // | |
73  // | Side |
74  // | Band2 |
75  // | |
76  // +-------------+-----------+
77  // |
78 
79  x.setRange("SB1",-10,+10) ;
80  y.setRange("SB1",-10,0) ;
81 
82  x.setRange("SB2",-10,0) ;
83  y.setRange("SB2",0,+10) ;
84 
85  x.setRange("SIG",0,+10) ;
86  y.setRange("SIG",0,+10) ;
87 
88  x.setRange("FULL",-10,+10) ;
89  y.setRange("FULL",-10,+10) ;
90 
91 
92  // P e r f o r m f i t s i n i n d i v i d u a l s i d e b a n d r e g i o n s
93  // -------------------------------------------------------------------------------------
94 
95  // Perform fit in SideBand1 region (RooAddPdf coefficients will be interpreted in full range)
96  RooFitResult* r_sb1 = model.fitTo(*modelData,Range("SB1"),Save()) ;
97 
98  // Perform fit in SideBand2 region (RooAddPdf coefficients will be interpreted in full range)
99  RooFitResult* r_sb2 = model.fitTo(*modelData,Range("SB2"),Save()) ;
100 
101 
102 
103  // P e r f o r m f i t s i n j o i n t s i d e b a n d r e g i o n s
104  // -----------------------------------------------------------------------------
105 
106  // Now perform fit to joint 'L-shaped' sideband region 'SB1|SB2'
107  // (RooAddPdf coefficients will be interpreted in full range)
108  RooFitResult* r_sb12 = model.fitTo(*modelData,Range("SB1,SB2"),Save()) ;
109 
110 
111  // Print results for comparison
112  r_sb1->Print() ;
113  r_sb2->Print() ;
114  r_sb12->Print() ;
115 
116 
117 }
RooAddPdf is an efficient implementation of a sum of PDFs of the form.
Definition: RooAddPdf.h:29
RooProdPdf is an efficient implementation of a product of PDFs of the form.
Definition: RooProdPdf.h:31
virtual void Print(Option_t *options=0) const
Print TNamed name and title.
Definition: RooFitResult.h:65
RooCmdArg Range(const char *rangeName, Bool_t adjustNorm=kTRUE)
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:37
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:29
double f(double x)
Double_t y[n]
Definition: legend1.C:17
RooCmdArg Save(Bool_t flag=kTRUE)
RooConstVar & RooConst(Double_t val)
RooPolynomial implements a polynomial p.d.f of the form By default coefficient a_0 is chosen to be 1...
Definition: RooPolynomial.h:28