Perform fits with different configurations using Minuit2
Author: Lorenzo Moneta
This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Monday, March 27, 2023 at 09:47 AM.
Definition of a helper function:
%%cpp -d
#include "TH1.h"
#include "TF1.h"
#include "TCanvas.h"
#include "TRandom3.h"
#include "TVirtualFitter.h"
#include "TPaveLabel.h"
#include "TStyle.h"
#include <iostream>
#include <string>
void testGausFit( std::string type = "Minuit2", int n = 1000) {
gRandom = new TRandom3();
TVirtualFitter::SetDefaultFitter(type.c_str() );
std::string name;
name = "h1_" + type;
TH1D * h1 = new TH1D(name.c_str(),"Chi2 Fit",100, -5, 5. );
name = "h2_" + type;
TH1D * h2 = new TH1D(name.c_str(),"Chi2 Fit with Minos Error",100, -5, 5. );
name = "h3_" + type;
TH1D * h3 = new TH1D(name.c_str(),"Chi2 Fit with Integral and Minos",100, -5, 5. );
name = "h4_" + type;
TH1D * h4 = new TH1D(name.c_str(),"Likelihood Fit with Minos Error",100, -5, 5. );
gStyle->SetOptStat(1111111);
gStyle->SetOptFit(1111111);
for (int i = 0; i < n; ++i) {
double x = gRandom->Gaus(0,1);
h1->Fill( x );
h2->Fill( x );
h3->Fill( x );
h4->Fill( x );
}
std::string cname = type + "Canvas" ;
std::string ctitle = type + " Gaussian Fit" ;
TCanvas *c1 = new TCanvas(cname.c_str(),cname.c_str(),10,10,900,900);
c1->Divide(2,2);
c1->cd(1);
std::cout << "\nDo Fit 1\n";
h1->Fit("gaus","Q");
h1->Draw();
c1->cd(2);
std::cout << "\nDo Fit 2\n";
h2->Fit("gaus","E");
h2->Draw();
c1->cd(3);
std::cout << "\nDo Fit 3\n";
h3->Fit("gaus","IGE");
h3->Draw();
c1->cd(4);
std::cout << "\nDo Fit 4\n";
h4->Fit("gaus","LE");
h4->Draw();
}
int n = 1000;
testGausFit("Minuit2",n);
testGausFit("Fumili2",n);
Do Fit 1 Do Fit 2 **************************************** Minimizer is Minuit2 / Migrad Chi2 = 65.1586 NDf = 56 Edm = 1.93774e-09 NCalls = 69 Constant = 36.3132 +/- 1.52625 -1.51651 +1.53547 (Minos) Mean = 0.013082 +/- 0.0347499 -0.0347674 +0.0347613 (Minos) Sigma = 1.03413 +/- 0.0288039 -0.0286274 +0.0290102 (Minos) (limited) Do Fit 3 **************************************** Minimizer is Minuit2 / Migrad Chi2 = 65.1586 NDf = 56 Edm = 5.21419e-08 NCalls = 71 Constant = 36.3276 +/- 2 -1.5174 +1.53672 (Minos) Mean = 0.0130818 +/- 2 Sigma = 1.03372 +/- 6.72117 (limited) Do Fit 4 **************************************** Minimizer is Minuit2 / Migrad MinFCN = 43.3935 Chi2 = 86.7869 NDf = 97 Edm = 9.97216e-08 NCalls = 62 Constant = 38.427 +/- 1.48837 -1.46667 +1.51031 (Minos) Mean = 0.027601 +/- 0.032831 -0.0328395 +0.0328395 (Minos) Sigma = 1.03819 +/- 0.0232194 -0.0227841 +0.0236699 (Minos) (limited) Do Fit 1 Do Fit 2 **************************************** Minimizer is Minuit2 / Fumili Chi2 = 65.1586 NDf = 56 Edm = 8.05711e-09 NCalls = 45 Constant = 36.3131 +/- 1.52625 -1.51642 +1.53556 (Minos) Mean = 0.0130818 +/- 0.0347499 -0.0347671 +0.0347615 (Minos) Sigma = 1.03413 +/- 0.0288039 -0.0286291 +0.0290085 (Minos) (limited) Do Fit 3 **************************************** Minimizer is Minuit2 / Fumili Chi2 = 65.1586 NDf = 56 Edm = 1.52369e-08 NCalls = 45 Constant = 36.3272 +/- 1.52734 -1.51745 +1.53671 (Minos) Mean = 0.0130818 +/- 0.0347499 -0.0347671 +0.0347615 (Minos) Sigma = 1.03373 +/- 0.0288151 -0.0286415 +0.0290186 (Minos) (limited) Do Fit 4 **************************************** Minimizer is Minuit2 / Fumili MinFCN = 43.3935 Chi2 = 86.7869 NDf = 97 Edm = 3.1863e-08 NCalls = 45 Constant = 38.4264 +/- 1.48835 -1.46601 +1.51097 (Minos) Mean = 0.0275931 +/- 0.0328313 -0.0328316 +0.0328474 (Minos) Sigma = 1.0382 +/- 0.0232197 -0.0227928 +0.0236612 (Minos) (limited)
Draw all canvases
%jsroot on
gROOT->GetListOfCanvases()->Draw()