minuit2FitBench2D.C: Quadratic background function | Fitting tutorials | multidimfit.C: Multi-Dimensional Parametrisation and Fitting |
// @(#)root/minuit2:$Id: minuit2GausFit.C 20881 2007-11-19 11:23:56Z rdm $ // Author: L. Moneta 10/2005 /********************************************************************** * * * Copyright (c) 2005 ROOT Foundation, CERN/PH-SFT * * * **********************************************************************/ #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); cout << "\nDo Fit 1\n"; h1->Fit("gaus","Q"); h1->Draw(); c1->cd(2); cout << "\nDo Fit 2\n"; h2->Fit("gaus","VE"); h2->Draw(); c1->cd(3); cout << "\nDo Fit 3\n"; h3->Fit("gaus","IE"); h3->Draw(); c1->cd(4); cout << "\nDo Fit 4\n"; h4->Fit("gaus","VLE"); h4->Draw(); } void minuit2GausFit() { int n = 1000; testGausFit("Minuit2",n); testGausFit("Fumili2",n); } minuit2GausFit.C:1 minuit2GausFit.C:2 minuit2GausFit.C:3 minuit2GausFit.C:4 minuit2GausFit.C:5 minuit2GausFit.C:6 minuit2GausFit.C:7 minuit2GausFit.C:8 minuit2GausFit.C:9 minuit2GausFit.C:10 minuit2GausFit.C:11 minuit2GausFit.C:12 minuit2GausFit.C:13 minuit2GausFit.C:14 minuit2GausFit.C:15 minuit2GausFit.C:16 minuit2GausFit.C:17 minuit2GausFit.C:18 minuit2GausFit.C:19 minuit2GausFit.C:20 minuit2GausFit.C:21 minuit2GausFit.C:22 minuit2GausFit.C:23 minuit2GausFit.C:24 minuit2GausFit.C:25 minuit2GausFit.C:26 minuit2GausFit.C:27 minuit2GausFit.C:28 minuit2GausFit.C:29 minuit2GausFit.C:30 minuit2GausFit.C:31 minuit2GausFit.C:32 minuit2GausFit.C:33 minuit2GausFit.C:34 minuit2GausFit.C:35 minuit2GausFit.C:36 minuit2GausFit.C:37 minuit2GausFit.C:38 minuit2GausFit.C:39 minuit2GausFit.C:40 minuit2GausFit.C:41 minuit2GausFit.C:42 minuit2GausFit.C:43 minuit2GausFit.C:44 minuit2GausFit.C:45 minuit2GausFit.C:46 minuit2GausFit.C:47 minuit2GausFit.C:48 minuit2GausFit.C:49 minuit2GausFit.C:50 minuit2GausFit.C:51 minuit2GausFit.C:52 minuit2GausFit.C:53 minuit2GausFit.C:54 minuit2GausFit.C:55 minuit2GausFit.C:56 minuit2GausFit.C:57 minuit2GausFit.C:58 minuit2GausFit.C:59 minuit2GausFit.C:60 minuit2GausFit.C:61 minuit2GausFit.C:62 minuit2GausFit.C:63 minuit2GausFit.C:64 minuit2GausFit.C:65 minuit2GausFit.C:66 minuit2GausFit.C:67 minuit2GausFit.C:68 minuit2GausFit.C:69 minuit2GausFit.C:70 minuit2GausFit.C:71 minuit2GausFit.C:72 minuit2GausFit.C:73 minuit2GausFit.C:74 minuit2GausFit.C:75 minuit2GausFit.C:76 minuit2GausFit.C:77 minuit2GausFit.C:78 minuit2GausFit.C:79 minuit2GausFit.C:80 minuit2GausFit.C:81 minuit2GausFit.C:82 minuit2GausFit.C:83 minuit2GausFit.C:84 |
|