Logo ROOT   6.10/09
Reference Guide
rf204_extrangefit.C
Go to the documentation of this file.
1 /// \file
2 /// \ingroup tutorial_roofit
3 /// \notebook -nodraw
4 /// 'ADDITION AND CONVOLUTION' RooFit tutorial macro #204
5 ///
6 /// Extended maximum likelihood fit with alternate range definition
7 /// for observed number of events.
8 ///
9 /// \macro_output
10 /// \macro_code
11 /// \author 07/2008 - Wouter Verkerke
12 
13 
14 #include "RooRealVar.h"
15 #include "RooDataSet.h"
16 #include "RooGaussian.h"
17 #include "RooChebychev.h"
18 #include "RooAddPdf.h"
19 #include "RooExtendPdf.h"
20 #include "RooFitResult.h"
21 #include "TCanvas.h"
22 #include "TAxis.h"
23 #include "RooPlot.h"
24 using namespace RooFit ;
25 
26 
27 void rf204_extrangefit()
28 {
29 
30 
31  // S e t u p c o m p o n e n t p d f s
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) ;
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 
55  // C o n s t r u c t e x t e n d e d c o m p s wi t h r a n g e s p e c
56  // ------------------------------------------------------------------------------
57 
58  // Define signal range in which events counts are to be defined
59  x.setRange("signalRange",4,6) ;
60 
61  // Associated nsig/nbkg as expected number of events with sig/bkg _in_the_range_ "signalRange"
62  RooRealVar nsig("nsig","number of signal events in signalRange",500,0.,10000) ;
63  RooRealVar nbkg("nbkg","number of background events in signalRange",500,0,10000) ;
64  RooExtendPdf esig("esig","extended signal p.d.f",sig,nsig,"signalRange") ;
65  RooExtendPdf ebkg("ebkg","extended background p.d.f",bkg,nbkg,"signalRange") ;
66 
67 
68  // S u m e x t e n d e d c o m p o n e n t s
69  // ---------------------------------------------
70 
71  // Construct sum of two extended p.d.f. (no coefficients required)
72  RooAddPdf model("model","(g1+g2)+a",RooArgList(ebkg,esig)) ;
73 
74 
75  // S a m p l e d a t a , f i t m o d e l
76  // -------------------------------------------
77 
78  // Generate 1000 events from model so that nsig,nbkg come out to numbers <<500 in fit
79  RooDataSet *data = model.generate(x,1000) ;
80 
81 
82  // Perform unbinned extended ML fit to data
83  RooFitResult* r = model.fitTo(*data,Extended(kTRUE),Save()) ;
84  r->Print() ;
85 
86 }
RooAddPdf is an efficient implementation of a sum of PDFs of the form.
Definition: RooAddPdf.h:29
RooCmdArg Extended(Bool_t flag=kTRUE)
virtual void Print(Option_t *options=0) const
Print TNamed name and title.
Definition: RooFitResult.h:64
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
TRandom2 r(17)
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:29
RooCmdArg Save(Bool_t flag=kTRUE)
Chebychev polynomial p.d.f.
Definition: RooChebychev.h:25
const Bool_t kTRUE
Definition: RtypesCore.h:91