ROOT logo

From $ROOTSYS/tutorials/roofit/rf601_intminuit.C

/////////////////////////////////////////////////////////////////////////
//
// 'LIKELIHOOD AND MINIMIZATION' RooFit tutorial macro #601
// 
// Interactive minimization with MINUIT
//
//
// 07/2008 - Wouter Verkerke 
//
/////////////////////////////////////////////////////////////////////////

#ifndef __CINT__
#include "RooGlobalFunc.h"
#endif
#include "RooRealVar.h"
#include "RooDataSet.h"
#include "RooGaussian.h"
#include "RooConstVar.h"
#include "RooProdPdf.h"
#include "RooAddPdf.h"
#include "RooMinuit.h"
#include "RooFitResult.h"
#include "RooPlot.h"
#include "TCanvas.h"
#include "TAxis.h"
#include "TH1.h"
using namespace RooFit ;


void rf601_intminuit()
{
  // S e t u p   p d f   a n d   l i k e l i h o o d 
  // -----------------------------------------------

  // Observable
  RooRealVar x("x","x",-20,20) ;

  // Model (intentional strong correlations)
  RooRealVar mean("mean","mean of g1 and g2",0) ;
  RooRealVar sigma_g1("sigma_g1","width of g1",3) ;
  RooGaussian g1("g1","g1",x,mean,sigma_g1) ;

  RooRealVar sigma_g2("sigma_g2","width of g2",4,3.0,6.0) ;
  RooGaussian g2("g2","g2",x,mean,sigma_g2) ;

  RooRealVar frac("frac","frac",0.5,0.0,1.0) ;
  RooAddPdf model("model","model",RooArgList(g1,g2),frac) ;

  // Generate 1000 events
  RooDataSet* data = model.generate(x,1000) ;
  
  // Construct unbinned likelihood of model w.r.t. data
  RooAbsReal* nll = model.createNLL(*data) ;

  // 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
  // -------------------------------------------------------------------------------

  // Create MINUIT interface object
  RooMinuit m(*nll) ;

  // Activate verbose logging of MINUIT parameter space stepping
  m.setVerbose(kTRUE) ;

  // Call MIGRAD to minimize the likelihood
  m.migrad() ;

  // Print values of all paramaters, that reflect values (and error estimates)
  // that are back propagated from MINUIT
  model.getParameters(x)->Print("s") ;

  // Disable verbose logging
  m.setVerbose(kFALSE) ;

  // Run HESSE to calculate errors from d2L/dp2
  m.hesse() ;

  // Print value (and error) of sigma_g2 parameter, that reflects
  // value and error back propagated from MINUIT
  sigma_g2.Print() ;

  // Run MINOS on sigma_g2 parameter only
  m.minos(sigma_g2) ;

  // Print value (and error) of sigma_g2 parameter, that reflects
  // value and error back propagated from MINUIT
  sigma_g2.Print() ;



  // 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 
  // ---------------------------------------------------------

  // Save a snapshot of the fit result. This object contains the initial
  // fit parameters, the final fit parameters, the complete correlation
  // matrix, the EDM, the minimized FCN , the last MINUIT status code and
  // the number of times the RooFit function object has indicated evaluation
  // problems (e.g. zero probabilities during likelihood evaluation)
  RooFitResult* r = m.save() ;

  // Make contour plot of mx vs sx at 1,2,3 sigma
  RooPlot* frame = m.contour(frac,sigma_g2,1,2,3) ;
  frame->SetTitle("RooMinuit contour plot") ;

  // Print the fit result snapshot
  r->Print("v") ;



  // 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
  // -----------------------------------------------------------------

  // At any moment you can manually change the value of a (constant)
  // parameter
  mean = 0.3 ;
  
  // Rerun MIGRAD,HESSE
  m.migrad() ;
  m.hesse() ;
  frac.Print() ;

  // Now fix sigma_g2
  sigma_g2.setConstant(kTRUE) ;

  // Rerun MIGRAD,HESSE
  m.migrad() ;
  m.hesse() ;
  frac.Print() ;



  new TCanvas("rf601_intminuit","rf601_intminuit",600,600) ;
  gPad->SetLeftMargin(0.15) ; frame->GetYaxis()->SetTitleOffset(1.4) ; frame->Draw() ;

}
 rf601_intminuit.C:1
 rf601_intminuit.C:2
 rf601_intminuit.C:3
 rf601_intminuit.C:4
 rf601_intminuit.C:5
 rf601_intminuit.C:6
 rf601_intminuit.C:7
 rf601_intminuit.C:8
 rf601_intminuit.C:9
 rf601_intminuit.C:10
 rf601_intminuit.C:11
 rf601_intminuit.C:12
 rf601_intminuit.C:13
 rf601_intminuit.C:14
 rf601_intminuit.C:15
 rf601_intminuit.C:16
 rf601_intminuit.C:17
 rf601_intminuit.C:18
 rf601_intminuit.C:19
 rf601_intminuit.C:20
 rf601_intminuit.C:21
 rf601_intminuit.C:22
 rf601_intminuit.C:23
 rf601_intminuit.C:24
 rf601_intminuit.C:25
 rf601_intminuit.C:26
 rf601_intminuit.C:27
 rf601_intminuit.C:28
 rf601_intminuit.C:29
 rf601_intminuit.C:30
 rf601_intminuit.C:31
 rf601_intminuit.C:32
 rf601_intminuit.C:33
 rf601_intminuit.C:34
 rf601_intminuit.C:35
 rf601_intminuit.C:36
 rf601_intminuit.C:37
 rf601_intminuit.C:38
 rf601_intminuit.C:39
 rf601_intminuit.C:40
 rf601_intminuit.C:41
 rf601_intminuit.C:42
 rf601_intminuit.C:43
 rf601_intminuit.C:44
 rf601_intminuit.C:45
 rf601_intminuit.C:46
 rf601_intminuit.C:47
 rf601_intminuit.C:48
 rf601_intminuit.C:49
 rf601_intminuit.C:50
 rf601_intminuit.C:51
 rf601_intminuit.C:52
 rf601_intminuit.C:53
 rf601_intminuit.C:54
 rf601_intminuit.C:55
 rf601_intminuit.C:56
 rf601_intminuit.C:57
 rf601_intminuit.C:58
 rf601_intminuit.C:59
 rf601_intminuit.C:60
 rf601_intminuit.C:61
 rf601_intminuit.C:62
 rf601_intminuit.C:63
 rf601_intminuit.C:64
 rf601_intminuit.C:65
 rf601_intminuit.C:66
 rf601_intminuit.C:67
 rf601_intminuit.C:68
 rf601_intminuit.C:69
 rf601_intminuit.C:70
 rf601_intminuit.C:71
 rf601_intminuit.C:72
 rf601_intminuit.C:73
 rf601_intminuit.C:74
 rf601_intminuit.C:75
 rf601_intminuit.C:76
 rf601_intminuit.C:77
 rf601_intminuit.C:78
 rf601_intminuit.C:79
 rf601_intminuit.C:80
 rf601_intminuit.C:81
 rf601_intminuit.C:82
 rf601_intminuit.C:83
 rf601_intminuit.C:84
 rf601_intminuit.C:85
 rf601_intminuit.C:86
 rf601_intminuit.C:87
 rf601_intminuit.C:88
 rf601_intminuit.C:89
 rf601_intminuit.C:90
 rf601_intminuit.C:91
 rf601_intminuit.C:92
 rf601_intminuit.C:93
 rf601_intminuit.C:94
 rf601_intminuit.C:95
 rf601_intminuit.C:96
 rf601_intminuit.C:97
 rf601_intminuit.C:98
 rf601_intminuit.C:99
 rf601_intminuit.C:100
 rf601_intminuit.C:101
 rf601_intminuit.C:102
 rf601_intminuit.C:103
 rf601_intminuit.C:104
 rf601_intminuit.C:105
 rf601_intminuit.C:106
 rf601_intminuit.C:107
 rf601_intminuit.C:108
 rf601_intminuit.C:109
 rf601_intminuit.C:110
 rf601_intminuit.C:111
 rf601_intminuit.C:112
 rf601_intminuit.C:113
 rf601_intminuit.C:114
 rf601_intminuit.C:115
 rf601_intminuit.C:116
 rf601_intminuit.C:117
 rf601_intminuit.C:118
 rf601_intminuit.C:119
 rf601_intminuit.C:120
 rf601_intminuit.C:121
 rf601_intminuit.C:122
 rf601_intminuit.C:123
 rf601_intminuit.C:124
 rf601_intminuit.C:125
 rf601_intminuit.C:126
 rf601_intminuit.C:127
 rf601_intminuit.C:128
 rf601_intminuit.C:129
 rf601_intminuit.C:130
 rf601_intminuit.C:131
 rf601_intminuit.C:132
 rf601_intminuit.C:133
 rf601_intminuit.C:134
 rf601_intminuit.C:135