line3Dfit.C: Fitting of a TGraph2D with a 3D straight line | Fitting tutorials | minuit2FitBench2D.C: Quadratic background function |
//+ Fitting 1-D histograms with minuit2 // @(#)root/minuit2:$Id: minuit2FitBench.C 26949 2008-12-16 11:38:17Z brun $ // Author: L. Moneta 10/2005 /********************************************************************** * * * Copyright (c) 2005 ROOT Foundation, CERN/PH-SFT * * * **********************************************************************/ #include "TH1.h" #include "TF1.h" #include "TCanvas.h" #include "TStopwatch.h" #include "TSystem.h" #include "TRandom3.h" #include "TVirtualFitter.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]); } void 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(); TVirtualFitter::SetDefaultFitter(fitter); //ROOT::Fit::FitConfig::SetDefaultMinimizer(fitter); 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); timer.Start(); for (Int_t pass=0;pass<npass;pass++) { if (pass%100 == 0) printf("pass : %d\n",pass); fitFcn->SetParameters(1,1,1,6,.03,1); for (Int_t i=0;i<5000;i++) { histo->Fill(fitFcn->GetRandom()); } histo->Fit(fitFcn,"Q0"); } histo->Fit(fitFcn,"EV"); 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(); p->SetTextColor(kRed+3); p->SetFillColor(kYellow-8); pad->Update(); } void minuit2FitBench(Int_t npass=20) { TH1::AddDirectory(kFALSE); 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); gStyle->SetOptFit(); gStyle->SetStatY(0.6); //with Minuit c1->cd(1); DoFit("Minuit",gPad,npass); //with Fumili c1->cd(2); DoFit("Fumili",gPad,npass); //with Minuit2 c1->cd(3); DoFit("Minuit2",gPad,npass); //with Fumili2 c1->cd(4); DoFit("Fumili2",gPad,npass); c1->SaveAs("FitBench.root"); } minuit2FitBench.C:1 minuit2FitBench.C:2 minuit2FitBench.C:3 minuit2FitBench.C:4 minuit2FitBench.C:5 minuit2FitBench.C:6 minuit2FitBench.C:7 minuit2FitBench.C:8 minuit2FitBench.C:9 minuit2FitBench.C:10 minuit2FitBench.C:11 minuit2FitBench.C:12 minuit2FitBench.C:13 minuit2FitBench.C:14 minuit2FitBench.C:15 minuit2FitBench.C:16 minuit2FitBench.C:17 minuit2FitBench.C:18 minuit2FitBench.C:19 minuit2FitBench.C:20 minuit2FitBench.C:21 minuit2FitBench.C:22 minuit2FitBench.C:23 minuit2FitBench.C:24 minuit2FitBench.C:25 minuit2FitBench.C:26 minuit2FitBench.C:27 minuit2FitBench.C:28 minuit2FitBench.C:29 minuit2FitBench.C:30 minuit2FitBench.C:31 minuit2FitBench.C:32 minuit2FitBench.C:33 minuit2FitBench.C:34 minuit2FitBench.C:35 minuit2FitBench.C:36 minuit2FitBench.C:37 minuit2FitBench.C:38 minuit2FitBench.C:39 minuit2FitBench.C:40 minuit2FitBench.C:41 minuit2FitBench.C:42 minuit2FitBench.C:43 minuit2FitBench.C:44 minuit2FitBench.C:45 minuit2FitBench.C:46 minuit2FitBench.C:47 minuit2FitBench.C:48 minuit2FitBench.C:49 minuit2FitBench.C:50 minuit2FitBench.C:51 minuit2FitBench.C:52 minuit2FitBench.C:53 minuit2FitBench.C:54 minuit2FitBench.C:55 minuit2FitBench.C:56 minuit2FitBench.C:57 minuit2FitBench.C:58 minuit2FitBench.C:59 minuit2FitBench.C:60 minuit2FitBench.C:61 minuit2FitBench.C:62 minuit2FitBench.C:63 minuit2FitBench.C:64 minuit2FitBench.C:65 minuit2FitBench.C:66 minuit2FitBench.C:67 minuit2FitBench.C:68 minuit2FitBench.C:69 minuit2FitBench.C:70 minuit2FitBench.C:71 minuit2FitBench.C:72 minuit2FitBench.C:73 minuit2FitBench.C:74 minuit2FitBench.C:75 minuit2FitBench.C:76 minuit2FitBench.C:77 minuit2FitBench.C:78 minuit2FitBench.C:79 minuit2FitBench.C:80 minuit2FitBench.C:81 minuit2FitBench.C:82 minuit2FitBench.C:83 minuit2FitBench.C:84 minuit2FitBench.C:85 minuit2FitBench.C:86 minuit2FitBench.C:87 minuit2FitBench.C:88 minuit2FitBench.C:89 minuit2FitBench.C:90 minuit2FitBench.C:91 minuit2FitBench.C:92 minuit2FitBench.C:93 minuit2FitBench.C:94 minuit2FitBench.C:95 minuit2FitBench.C:96 minuit2FitBench.C:97 minuit2FitBench.C:98 minuit2FitBench.C:99 minuit2FitBench.C:100 minuit2FitBench.C:101 minuit2FitBench.C:102 minuit2FitBench.C:103 minuit2FitBench.C:104 minuit2FitBench.C:105 minuit2FitBench.C:106 minuit2FitBench.C:107 minuit2FitBench.C:108 minuit2FitBench.C:109 minuit2FitBench.C:110 minuit2FitBench.C:111 minuit2FitBench.C:112 minuit2FitBench.C:113 minuit2FitBench.C:114 minuit2FitBench.C:115 minuit2FitBench.C:116 minuit2FitBench.C:117 minuit2FitBench.C:118 |
|