Logo ROOT  
Reference Guide
rf601_intminuit.C
Go to the documentation of this file.
1 /// \file
2 /// \ingroup tutorial_roofit
3 /// \notebook -js
4 ///
5 /// Likelihood and minimization: interactive minimization with MINUIT
6 ///
7 /// \macro_image
8 /// \macro_output
9 /// \macro_code
10 ///
11 /// \date 07/2008
12 /// \author Wouter Verkerke
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 void rf601_intminuit()
29 {
30  // S e t u p p d f a n d l i k e l i h o o d
31  // -----------------------------------------------
32 
33  // Observable
34  RooRealVar x("x", "x", -20, 20);
35 
36  // Model (intentional strong correlations)
37  RooRealVar mean("mean", "mean of g1 and g2", 0);
38  RooRealVar sigma_g1("sigma_g1", "width of g1", 3);
39  RooGaussian g1("g1", "g1", x, mean, sigma_g1);
40 
41  RooRealVar sigma_g2("sigma_g2", "width of g2", 4, 3.0, 6.0);
42  RooGaussian g2("g2", "g2", x, mean, sigma_g2);
43 
44  RooRealVar frac("frac", "frac", 0.5, 0.0, 1.0);
45  RooAddPdf model("model", "model", RooArgList(g1, g2), frac);
46 
47  // Generate 1000 events
48  RooDataSet *data = model.generate(x, 1000);
49 
50  // Construct unbinned likelihood of model w.r.t. data
51  RooAbsReal *nll = model.createNLL(*data);
52 
53  // 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
54  // -------------------------------------------------------------------------------
55 
56  // Create MINUIT interface object
57  RooMinimizer m(*nll);
58 
59  // Activate verbose logging of MINUIT parameter space stepping
60  m.setVerbose(kTRUE);
61 
62  // Call MIGRAD to minimize the likelihood
63  m.migrad();
64 
65  // Print values of all parameters, that reflect values (and error estimates)
66  // that are back propagated from MINUIT
67  model.getParameters(x)->Print("s");
68 
69  // Disable verbose logging
70  m.setVerbose(kFALSE);
71 
72  // Run HESSE to calculate errors from d2L/dp2
73  m.hesse();
74 
75  // Print value (and error) of sigma_g2 parameter, that reflects
76  // value and error back propagated from MINUIT
77  sigma_g2.Print();
78 
79  // Run MINOS on sigma_g2 parameter only
80  m.minos(sigma_g2);
81 
82  // Print value (and error) of sigma_g2 parameter, that reflects
83  // value and error back propagated from MINUIT
84  sigma_g2.Print();
85 
86  // 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
87  // ---------------------------------------------------------
88 
89  // Save a snapshot of the fit result. This object contains the initial
90  // fit parameters, the final fit parameters, the complete correlation
91  // matrix, the EDM, the minimized FCN , the last MINUIT status code and
92  // the number of times the RooFit function object has indicated evaluation
93  // problems (e.g. zero probabilities during likelihood evaluation)
94  RooFitResult *r = m.save();
95 
96  // Make contour plot of mx vs sx at 1,2,3 sigma
97  RooPlot *frame = m.contour(frac, sigma_g2, 1, 2, 3);
98  frame->SetTitle("Minuit contour plot");
99 
100  // Print the fit result snapshot
101  r->Print("v");
102 
103  // 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
104  // -----------------------------------------------------------------
105 
106  // At any moment you can manually change the value of a (constant)
107  // parameter
108  mean = 0.3;
109 
110  // Rerun MIGRAD,HESSE
111  m.migrad();
112  m.hesse();
113  frac.Print();
114 
115  // Now fix sigma_g2
116  sigma_g2.setConstant(kTRUE);
117 
118  // Rerun MIGRAD,HESSE
119  m.migrad();
120  m.hesse();
121  frac.Print();
122 
123  new TCanvas("rf601_intminuit", "rf601_intminuit", 600, 600);
124  gPad->SetLeftMargin(0.15);
125  frame->GetYaxis()->SetTitleOffset(1.4);
126  frame->Draw();
127 }
m
auto * m
Definition: textangle.C:8
RooPlot::Draw
virtual void Draw(Option_t *options=0)
Draw this plot and all of the elements it contains.
Definition: RooPlot.cxx:691
kTRUE
const Bool_t kTRUE
Definition: RtypesCore.h:91
RooAddPdf
Definition: RooAddPdf.h:32
TObject::Print
virtual void Print(Option_t *option="") const
This method must be overridden when a class wants to print itself.
Definition: TObject.cxx:552
r
ROOT::R::TRInterface & r
Definition: Object.C:4
RooArgList
Definition: RooArgList.h:21
RooGaussian.h
x
Double_t x[n]
Definition: legend1.C:17
RooGaussian
Definition: RooGaussian.h:25
RooAddPdf.h
RooAbsReal
Definition: RooAbsReal.h:61
TCanvas.h
RooDataSet.h
RooFitResult
Definition: RooFitResult.h:40
RooProdPdf.h
RooFit
Definition: RooCFunction1Binding.h:29
kFALSE
const Bool_t kFALSE
Definition: RtypesCore.h:92
RooPlot.h
RooPlot::GetYaxis
TAxis * GetYaxis() const
Definition: RooPlot.cxx:1256
RooPlot
Definition: RooPlot.h:44
RooRealVar.h
RooFitResult.h
RooConstVar.h
RooPlot::SetTitle
void SetTitle(const char *name)
Set the title of the RooPlot to 'title'.
Definition: RooPlot.cxx:1237
TCanvas
Definition: TCanvas.h:23
RooMinimizer.h
TAxis.h
gPad
#define gPad
Definition: TVirtualPad.h:287
RooDataSet
Definition: RooDataSet.h:33
make_cnn_model.model
model
Definition: make_cnn_model.py:6
TAttAxis::SetTitleOffset
virtual void SetTitleOffset(Float_t offset=1)
Set distance between the axis and the axis title.
Definition: TAttAxis.cxx:293
RooRealVar
Definition: RooRealVar.h:35
TH1.h
RooMinimizer
Definition: RooMinimizer.h:40