Logo ROOT   6.08/07
Reference Guide
rf601_intminuit.C
Go to the documentation of this file.
1 /// \file
2 /// \ingroup tutorial_roofit
3 /// \notebook -js
4 /// 'LIKELIHOOD AND MINIMIZATION' RooFit tutorial macro #601
5 ///
6 /// Interactive minimization with MINUIT
7 ///
8 /// \macro_image
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 "RooConstVar.h"
18 #include "RooProdPdf.h"
19 #include "RooAddPdf.h"
20 #include "RooMinimizer.h"
21 #include "RooFitResult.h"
22 #include "RooPlot.h"
23 #include "TCanvas.h"
24 #include "TAxis.h"
25 #include "TH1.h"
26 using namespace RooFit ;
27 
28 
29 void rf601_intminuit()
30 {
31  // S e t u p p d f a n d l i k e l i h o o d
32  // -----------------------------------------------
33 
34  // Observable
35  RooRealVar x("x","x",-20,20) ;
36 
37  // Model (intentional strong correlations)
38  RooRealVar mean("mean","mean of g1 and g2",0) ;
39  RooRealVar sigma_g1("sigma_g1","width of g1",3) ;
40  RooGaussian g1("g1","g1",x,mean,sigma_g1) ;
41 
42  RooRealVar sigma_g2("sigma_g2","width of g2",4,3.0,6.0) ;
43  RooGaussian g2("g2","g2",x,mean,sigma_g2) ;
44 
45  RooRealVar frac("frac","frac",0.5,0.0,1.0) ;
46  RooAddPdf model("model","model",RooArgList(g1,g2),frac) ;
47 
48  // Generate 1000 events
49  RooDataSet* data = model.generate(x,1000) ;
50 
51  // Construct unbinned likelihood of model w.r.t. data
52  RooAbsReal* nll = model.createNLL(*data) ;
53 
54  // I n t e r a c t i v e m i n i m i z a t i o n , e r r o r a n a l y s i s
55  // -------------------------------------------------------------------------------
56 
57  // Create MINUIT interface object
58  RooMinimizer m(*nll) ;
59 
60  // Activate verbose logging of MINUIT parameter space stepping
61  m.setVerbose(kTRUE) ;
62 
63  // Call MIGRAD to minimize the likelihood
64  m.migrad() ;
65 
66  // Print values of all parameters, that reflect values (and error estimates)
67  // that are back propagated from MINUIT
68  model.getParameters(x)->Print("s") ;
69 
70  // Disable verbose logging
71  m.setVerbose(kFALSE) ;
72 
73  // Run HESSE to calculate errors from d2L/dp2
74  m.hesse() ;
75 
76  // Print value (and error) of sigma_g2 parameter, that reflects
77  // value and error back propagated from MINUIT
78  sigma_g2.Print() ;
79 
80  // Run MINOS on sigma_g2 parameter only
81  m.minos(sigma_g2) ;
82 
83  // Print value (and error) of sigma_g2 parameter, that reflects
84  // value and error back propagated from MINUIT
85  sigma_g2.Print() ;
86 
87 
88 
89  // S a v i n g r e s u l t s , c o n t o u r p l o t s
90  // ---------------------------------------------------------
91 
92  // Save a snapshot of the fit result. This object contains the initial
93  // fit parameters, the final fit parameters, the complete correlation
94  // matrix, the EDM, the minimized FCN , the last MINUIT status code and
95  // the number of times the RooFit function object has indicated evaluation
96  // problems (e.g. zero probabilities during likelihood evaluation)
97  RooFitResult* r = m.save() ;
98 
99  // Make contour plot of mx vs sx at 1,2,3 sigma
100  RooPlot* frame = m.contour(frac,sigma_g2,1,2,3) ;
101  frame->SetTitle("Minuit contour plot") ;
102 
103  // Print the fit result snapshot
104  r->Print("v") ;
105 
106 
107 
108  // C h a n g e p a r a m e t e r v a l u e s , f l o a t i n g
109  // -----------------------------------------------------------------
110 
111  // At any moment you can manually change the value of a (constant)
112  // parameter
113  mean = 0.3 ;
114 
115  // Rerun MIGRAD,HESSE
116  m.migrad() ;
117  m.hesse() ;
118  frac.Print() ;
119 
120  // Now fix sigma_g2
121  sigma_g2.setConstant(kTRUE) ;
122 
123  // Rerun MIGRAD,HESSE
124  m.migrad() ;
125  m.hesse() ;
126  frac.Print() ;
127 
128 
129 
130  new TCanvas("rf601_intminuit","rf601_intminuit",600,600) ;
131  gPad->SetLeftMargin(0.15) ; frame->GetYaxis()->SetTitleOffset(1.4) ; frame->Draw() ;
132 
133 }
virtual void SetTitleOffset(Float_t offset=1)
Set distance between the axis and the axis title Offset is a correction factor with respect to the "s...
Definition: TAttAxis.cxx:262
RooAddPdf is an efficient implementation of a sum of PDFs of the form.
Definition: RooAddPdf.h:29
TAxis * GetYaxis() const
Definition: RooPlot.cxx:1118
virtual void Print(Option_t *option="") const
Dump this marker with its attributes.
Definition: TMarker.cxx:277
const Bool_t kFALSE
Definition: Rtypes.h:92
void SetTitle(const char *name)
Set the title of the RooPlot to 'title'.
Definition: RooPlot.cxx:1099
virtual void Print(Option_t *options=0) const
Print TNamed name and title.
Definition: RooFitResult.h:65
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
TRandom2 r(17)
TMarker * m
Definition: textangle.C:8
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
The Canvas class.
Definition: TCanvas.h:41
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:53
RooMinimizer is a wrapper class around ROOT::Fit:Fitter that provides a seamless interface between th...
Definition: RooMinimizer.h:38
#define gPad
Definition: TVirtualPad.h:289
const Bool_t kTRUE
Definition: Rtypes.h:91
virtual void Draw(Option_t *options=0)
Draw this plot and all of the elements it contains.
Definition: RooPlot.cxx:559