Loading [MathJax]/extensions/tex2jax.js
Logo ROOT  
Reference Guide
 
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
Loading...
Searching...
No Matches
minuit2FitBench.C File Reference

Detailed Description

View in nbviewer Open in SWAN Demonstrate performance and usage of Minuit2 and Fumili2 for monodimensional fits.

*********************************************************************************
Minuit
*********************************************************************************
pass : 0
................... FCN=205.276 FROM MINOS STATUS=SUCCESSFUL 44 CALLS 429 TOTAL
EDM=3.83288e-10 STRATEGY= 1 ERROR MATRIX ACCURATE
EXT PARAMETER STEP FIRST
NO. NAME VALUE ERROR SIZE DERIVATIVE
1 p0 5.13639e+01 2.01329e+00 -2.79418e-04 -2.05471e-06
2 p1 5.57813e+01 4.80582e+00 3.09127e-03 -9.98919e-07
3 p2 7.42112e+01 1.87041e+00 -1.20311e-03 -1.93173e-07
4 p3 4.27344e+02 2.93232e+00 -1.66243e-02 -7.80957e-07
5 p4 3.58604e-02 3.47005e-04 1.74159e-07 9.80777e-02
6 p5 1.00001e+00 1.64203e-04 1.64203e-04 3.19213e-02
Minuit, npass=20 : RT= 0.175 s, Cpu= 0.180 s
*********************************************************************************
Fumili
*********************************************************************************
pass : 0
...................
****************************************
Minimizer is Fumili
Chi2 = 206.284
NDf = 194
NCalls = 4
p0 = 51.4325 +/- 2.01397
p1 = 55.5412 +/- 4.81253
p2 = 74.2976 +/- 1.87298
p3 = 427.425 +/- 2.93868
p4 = 0.0358559 +/- 0.000357243
p5 = 1.00001 +/- 0.00016009
Fumili, npass=20 : RT= 0.092 s, Cpu= 0.090 s
*********************************************************************************
Minuit2
*********************************************************************************
pass : 0
...................
****************************************
Minimizer is Minuit2 / Migrad
Chi2 = 205.34
NDf = 194
Edm = 1.91395e-10
NCalls = 85
p0 = 51.3576 +/- 2.0133 -2.01329 +2.0133 (Minos)
p1 = 55.8172 +/- 4.80582 -4.80605 +4.80561 (Minos)
p2 = 74.1981 +/- 1.8704 -1.87032 +1.8705 (Minos)
p3 = 427.317 +/- 2.93218 -2.93197 +2.93246 (Minos)
p4 = 0.0358574 +/- 0.000346961 -0.000345228 +0.000348733 (Minos)
p5 = 1.00001 +/- 0.000164204 -0.000164437 +0.000163973 (Minos)
Minuit2, npass=20 : RT= 0.085 s, Cpu= 0.090 s
*********************************************************************************
Fumili2
*********************************************************************************
pass : 0
...................
****************************************
Minimizer is Minuit2 / Fumili
Chi2 = 207.495
NDf = 194
Edm = 4.34894e-10
NCalls = 92
p0 = 51.4278 +/- 2.01381
p1 = 55.5389 +/- 4.8062
p2 = 74.3005 +/- 1.8705
p3 = 427.406 +/- 2.93236
p4 = 0.0358584 +/- 0.000346902
p5 = 1.00001 +/- 0.000164198
Fumili2, npass=20 : RT= 0.063 s, Cpu= 0.060 s
(int) 0
#include "TH1.h"
#include "TF1.h"
#include "TCanvas.h"
#include "TStopwatch.h"
#include "TSystem.h"
#include "TRandom3.h"
#include "TPaveLabel.h"
#include "TStyle.h"
#include "TMath.h"
#include "TROOT.h"
#include "TFrame.h"
/*#include "Fit/FitConfig.h"*/
TF1 *fitFcn;
TH1 *histo;
// Quadratic background function
Double_t background(Double_t *x, Double_t *par) {
return par[0] + par[1]*x[0] + par[2]*x[0]*x[0];
}
// Lorenzian Peak function
Double_t lorentzianPeak(Double_t *x, Double_t *par) {
return (0.5*par[0]*par[1]/TMath::Pi()) /
TMath::Max( 1.e-10,(x[0]-par[2])*(x[0]-par[2]) + .25*par[1]*par[1]);
}
// Sum of background and peak function
Double_t fitFunction(Double_t *x, Double_t *par) {
return background(x,par) + lorentzianPeak(x,&par[3]);
}
bool DoFit(const char* fitter, TVirtualPad *pad, Int_t npass) {
printf("\n*********************************************************************************\n");
printf("\t %s \n",fitter);
printf("*********************************************************************************\n");
gRandom = new TRandom3();
TStopwatch timer;
// timer.Start();
pad->SetGrid();
pad->SetLogy();
fitFcn->SetParameters(1,1,1,6,.03,1);
fitFcn->Update();
std::string title = std::string(fitter) + " fit bench";
histo = new TH1D(fitter,title.c_str(),200,0,3);
TString fitterType(fitter);
timer.Start();
bool ok = true;
// fill histogram many times
// every time increase its statistics and re-use previous fitted
// parameter values as starting point
for (Int_t pass=0;pass<npass;pass++) {
if (pass%100 == 0) printf("pass : %d\n",pass);
else printf(".");
if (pass == 0)fitFcn->SetParameters(1,1,1,6,.03,1);
for (Int_t i=0;i<5000;i++) {
histo->Fill(fitFcn->GetRandom());
}
int iret = histo->Fit(fitFcn,"Q0");
ok &= (iret == 0);
if (iret!=0) Error("DoFit","Fit pass %d failed !",pass);
}
// do last fit computing Minos Errors (except for Fumili)
if (!fitterType.Contains("Fumili")) // Fumili does not implement Error options (MINOS)
histo->Fit(fitFcn,"E");
else
histo->Fit(fitFcn,"");
timer.Stop();
(histo->GetFunction("fitFcn"))->SetLineColor(kRed+3);
gPad->SetFillColor(kYellow-10);
Double_t cputime = timer.CpuTime();
printf("%s, npass=%d : RT=%7.3f s, Cpu=%7.3f s\n",fitter,npass,timer.RealTime(),cputime);
TPaveLabel *p = new TPaveLabel(0.45,0.7,0.88,0.8,Form("%s CPU= %g s",fitter,cputime),"brNDC");
p->Draw();
pad->Update();
return ok;
}
int minuit2FitBench(Int_t npass=20) {
TCanvas *c1 = new TCanvas("FitBench","Fitting Demo",10,10,900,900);
c1->Divide(2,2);
c1->SetFillColor(kYellow-9);
// create a TF1 with the range from 0 to 3 and 6 parameters
fitFcn = new TF1("fitFcn",fitFunction,0,3,6);
fitFcn->SetNpx(200);
bool ok = true;
//with Minuit
c1->cd(1);
ok &= DoFit("Minuit",gPad,npass);
//with Fumili
c1->cd(2);
ok &= DoFit("Fumili",gPad,npass);
//with Minuit2
c1->cd(3);
ok &= DoFit("Minuit2",gPad,npass);
//with Fumili2
c1->cd(4);
ok &= DoFit("Fumili2",gPad,npass);
c1->SaveAs("FitBench.root");
return (ok) ? 0 : 1;
}
int Int_t
Definition RtypesCore.h:45
const Bool_t kFALSE
Definition RtypesCore.h:101
double Double_t
Definition RtypesCore.h:59
@ kRed
Definition Rtypes.h:66
@ kYellow
Definition Rtypes.h:66
void Error(const char *location, const char *msgfmt,...)
Use this function in case an error occurred.
Definition TError.cxx:187
R__EXTERN TRandom * gRandom
Definition TRandom.h:62
char * Form(const char *fmt,...)
R__EXTERN TStyle * gStyle
Definition TStyle.h:413
#define gPad
static void SetDefaultMinimizer(const char *type, const char *algo=0)
virtual void SetFillColor(Color_t fcolor)
Set the fill area color.
Definition TAttFill.h:37
virtual void SetTextColor(Color_t tcolor=1)
Set the text color.
Definition TAttText.h:44
The Canvas class.
Definition TCanvas.h:23
1-Dim function class
Definition TF1.h:213
virtual void Update()
Called by functions such as SetRange, SetNpx, SetParameters to force the deletion of the associated h...
Definition TF1.cxx:3630
virtual Double_t GetRandom(TRandom *rng=nullptr, Option_t *opt=nullptr)
Return a random number following this function shape.
Definition TF1.cxx:2186
virtual void SetNpx(Int_t npx=100)
Set the number of points used to draw the function.
Definition TF1.cxx:3453
virtual void SetParameters(const Double_t *params)
Definition TF1.h:644
1-D histogram with a double per channel (see TH1 documentation)}
Definition TH1.h:618
TH1 is the base class of all histogram classes in ROOT.
Definition TH1.h:58
static void AddDirectory(Bool_t add=kTRUE)
Sets the flag controlling the automatic add of histograms in memory.
Definition TH1.cxx:1283
virtual TFitResultPtr Fit(const char *formula, Option_t *option="", Option_t *goption="", Double_t xmin=0, Double_t xmax=0)
Fit histogram with function fname.
Definition TH1.cxx:3893
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition TH1.cxx:3351
virtual TF1 * GetFunction(const char *name) const
Return pointer to function with name.
Definition TH1.cxx:8881
A Pave (see TPave) with a text centered in the Pave.
Definition TPaveLabel.h:20
virtual void Draw(Option_t *option="")
Draw this pavelabel with its current attributes.
Random number generator class based on M.
Definition TRandom3.h:27
Stopwatch class.
Definition TStopwatch.h:28
Double_t RealTime()
Stop the stopwatch (if it is running) and return the realtime (in seconds) passed between the start a...
void Start(Bool_t reset=kTRUE)
Start the stopwatch.
Double_t CpuTime()
Stop the stopwatch (if it is running) and return the cputime (in seconds) passed between the start an...
void Stop()
Stop the stopwatch.
Basic string class.
Definition TString.h:136
void SetStatY(Float_t y=0)
Definition TStyle.h:381
void SetOptFit(Int_t fit=1)
The type of information about fit parameters printed in the histogram statistics box can be selected ...
Definition TStyle.cxx:1541
TVirtualPad is an abstract base class for the Pad and Canvas classes.
Definition TVirtualPad.h:51
virtual void SetGrid(Int_t valuex=1, Int_t valuey=1)=0
virtual void Update()=0
virtual void SetLogy(Int_t value=1)=0
return c1
Definition legend1.C:41
Double_t x[n]
Definition legend1.C:17
Short_t Max(Short_t a, Short_t b)
Definition TMathBase.h:208
constexpr Double_t Pi()
Definition TMath.h:37
lv SetLineColor(kBlue)
Author
Lorenzo Moneta

Definition in file minuit2FitBench.C.