Logo ROOT   6.18/05
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: interactive minimization with MINUIT
5///
6/// \macro_image
7/// \macro_output
8/// \macro_code
9/// \author 07/2008 - Wouter Verkerke
10
11#include "RooRealVar.h"
12#include "RooDataSet.h"
13#include "RooGaussian.h"
14#include "RooConstVar.h"
15#include "RooProdPdf.h"
16#include "RooAddPdf.h"
17#include "RooMinimizer.h"
18#include "RooFitResult.h"
19#include "RooPlot.h"
20#include "TCanvas.h"
21#include "TAxis.h"
22#include "TH1.h"
23using namespace RooFit;
24
25void rf601_intminuit()
26{
27 // S e t u p p d f a n d l i k e l i h o o d
28 // -----------------------------------------------
29
30 // Observable
31 RooRealVar x("x", "x", -20, 20);
32
33 // Model (intentional strong correlations)
34 RooRealVar mean("mean", "mean of g1 and g2", 0);
35 RooRealVar sigma_g1("sigma_g1", "width of g1", 3);
36 RooGaussian g1("g1", "g1", x, mean, sigma_g1);
37
38 RooRealVar sigma_g2("sigma_g2", "width of g2", 4, 3.0, 6.0);
39 RooGaussian g2("g2", "g2", x, mean, sigma_g2);
40
41 RooRealVar frac("frac", "frac", 0.5, 0.0, 1.0);
42 RooAddPdf model("model", "model", RooArgList(g1, g2), frac);
43
44 // Generate 1000 events
45 RooDataSet *data = model.generate(x, 1000);
46
47 // Construct unbinned likelihood of model w.r.t. data
48 RooAbsReal *nll = model.createNLL(*data);
49
50 // 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
51 // -------------------------------------------------------------------------------
52
53 // Create MINUIT interface object
54 RooMinimizer m(*nll);
55
56 // Activate verbose logging of MINUIT parameter space stepping
57 m.setVerbose(kTRUE);
58
59 // Call MIGRAD to minimize the likelihood
60 m.migrad();
61
62 // Print values of all parameters, that reflect values (and error estimates)
63 // that are back propagated from MINUIT
64 model.getParameters(x)->Print("s");
65
66 // Disable verbose logging
67 m.setVerbose(kFALSE);
68
69 // Run HESSE to calculate errors from d2L/dp2
70 m.hesse();
71
72 // Print value (and error) of sigma_g2 parameter, that reflects
73 // value and error back propagated from MINUIT
74 sigma_g2.Print();
75
76 // Run MINOS on sigma_g2 parameter only
77 m.minos(sigma_g2);
78
79 // Print value (and error) of sigma_g2 parameter, that reflects
80 // value and error back propagated from MINUIT
81 sigma_g2.Print();
82
83 // 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
84 // ---------------------------------------------------------
85
86 // Save a snapshot of the fit result. This object contains the initial
87 // fit parameters, the final fit parameters, the complete correlation
88 // matrix, the EDM, the minimized FCN , the last MINUIT status code and
89 // the number of times the RooFit function object has indicated evaluation
90 // problems (e.g. zero probabilities during likelihood evaluation)
91 RooFitResult *r = m.save();
92
93 // Make contour plot of mx vs sx at 1,2,3 sigma
94 RooPlot *frame = m.contour(frac, sigma_g2, 1, 2, 3);
95 frame->SetTitle("Minuit contour plot");
96
97 // Print the fit result snapshot
98 r->Print("v");
99
100 // 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
101 // -----------------------------------------------------------------
102
103 // At any moment you can manually change the value of a (constant)
104 // parameter
105 mean = 0.3;
106
107 // Rerun MIGRAD,HESSE
108 m.migrad();
109 m.hesse();
110 frac.Print();
111
112 // Now fix sigma_g2
113 sigma_g2.setConstant(kTRUE);
114
115 // Rerun MIGRAD,HESSE
116 m.migrad();
117 m.hesse();
118 frac.Print();
119
120 new TCanvas("rf601_intminuit", "rf601_intminuit", 600, 600);
121 gPad->SetLeftMargin(0.15);
122 frame->GetYaxis()->SetTitleOffset(1.4);
123 frame->Draw();
124}
ROOT::R::TRInterface & r
Definition: Object.C:4
const Bool_t kFALSE
Definition: RtypesCore.h:88
const Bool_t kTRUE
Definition: RtypesCore.h:87
#define gPad
Definition: TVirtualPad.h:286
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:53
RooAddPdf is an efficient implementation of a sum of PDFs of the form.
Definition: RooAddPdf.h:29
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgList.h:21
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:31
RooFitResult is a container class to hold the input and output of a PDF fit to a dataset.
Definition: RooFitResult.h:40
Plain Gaussian p.d.f.
Definition: RooGaussian.h:25
RooMinimizer is a wrapper class around ROOT::Fit:Fitter that provides a seamless interface between th...
Definition: RooMinimizer.h:38
A RooPlot is a plot frame and a container for graphics objects within that frame.
Definition: RooPlot.h:41
void SetTitle(const char *name)
Set the title of the RooPlot to 'title'.
Definition: RooPlot.cxx:1104
TAxis * GetYaxis() const
Definition: RooPlot.cxx:1123
virtual void Draw(Option_t *options=0)
Draw this plot and all of the elements it contains.
Definition: RooPlot.cxx:558
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:36
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:294
The Canvas class.
Definition: TCanvas.h:31
virtual void Print(Option_t *option="") const
This method must be overridden when a class wants to print itself.
Definition: TObject.cxx:550
Double_t x[n]
Definition: legend1.C:17
Template specialisation used in RooAbsArg:
auto * m
Definition: textangle.C:8