Double_t langaufun(Double_t *x, Double_t *par) { //Fit parameters: //par[0]=Width (scale) parameter of Landau density //par[1]=Most Probable (MP, location) parameter of Landau density //par[2]=Total area (integral -inf to inf, normalization constant) //par[3]=Width (sigma) of convoluted Gaussian function // //In the Landau distribution (represented by the CERNLIB approximation), //the maximum is located at x=-0.22278298 with the location parameter=0. //This shift is corrected within this function, so that the actual //maximum is identical to the MP parameter. // Numeric constants Double_t invsq2pi = 0.3989422804014; // (2 pi)^(-1/2) Double_t mpshift = -0.22278298; // Landau maximum location // Control constants Double_t np = 100.0; // number of convolution steps Double_t sc = 5.0; // convolution extends to +-sc Gaussian sigmas // Variables Double_t xx; Double_t mpc; Double_t fland; Double_t sum = 0.0; Double_t xlow,xupp; Double_t step; Double_t i; // MP shift correction mpc = par[1] - mpshift * par[0]; // Range of convolution integral xlow = x[0] - sc * par[3]; xupp = x[0] + sc * par[3]; step = (xupp-xlow) / np; // Convolution integral of Landau and Gaussian by sum for(i=1.0; i<=np/2; i++) { xx = xlow + (i-.5) * step; fland = TMath::Landau(xx,mpc,par[0]) / par[0]; sum += fland * TMath::Gaus(x[0],xx,par[3]); xx = xupp - (i-.5) * step; fland = TMath::Landau(xx,mpc,par[0]) / par[0]; sum += fland * TMath::Gaus(x[0],xx,par[3]); } return (par[2] * step * sum * invsq2pi / par[3]); } void mip(TString region="WHA") { SetStyle(); TFile* f = TFile::Open("mip_p10.root"); TH1D* p5 = (TH1D*)f->Get("h_EH_"+region+"_p5CMX"); TH1D* p6 = (TH1D*)f->Get("h_EH_"+region+"_p6CMX"); TH1D* p7 = (TH1D*)f->Get("h_EH_"+region+"_p7CMX"); TH1D* p8 = (TH1D*)f->Get("h_EH_"+region+"_p8CMX"); TH1D* p9 = (TH1D*)f->Get("h_EH_"+region+"_p9CMX"); TH1D* p10 = (TH1D*)f->Get("h_EH_"+region+"_p10CMX"); p5->Add((TH1D*)f->Get("h_EH_"+region+"_p5CMIOCES")); p6->Add((TH1D*)f->Get("h_EH_"+region+"_p6CMIOCES")); p7->Add((TH1D*)f->Get("h_EH_"+region+"_p7CMIOCES")); p8->Add((TH1D*)f->Get("h_EH_"+region+"_p8CMIOCES")); p9->Add((TH1D*)f->Get("h_EH_"+region+"_p9CMIOCES")); p10->Add((TH1D*)f->Get("h_EH_"+region+"_p10CMIOCES")); p5->Add((TH1D*)f->Get("h_EH_"+region+"_p5CMUP")); p6->Add((TH1D*)f->Get("h_EH_"+region+"_p6CMUP")); p7->Add((TH1D*)f->Get("h_EH_"+region+"_p7CMUP")); p8->Add((TH1D*)f->Get("h_EH_"+region+"_p8CMUP")); p9->Add((TH1D*)f->Get("h_EH_"+region+"_p9CMUP")); p10->Add((TH1D*)f->Get("h_EH_"+region+"_p10CMUP")); p5->Sumw2(); p6->Sumw2(); p7->Sumw2(); p8->Sumw2(); p9->Sumw2(); p10->Sumw2(); TH1D* peaks = new TH1D("peaks","; Data Period; Peak Position [GeV];",6,4.5,10.5); TH1D* gwidth = new TH1D("gwidth","; Data Period; Guassian Width [GeV];",6,4.5,10.5); TH1D* lwidth = new TH1D("lwidth","; Data Period; Landau Width [GeV];",6,4.5,10.5); TCanvas* c = new TCanvas("c","c",10,10,800,600); p5->Rebin(2); p6->Rebin(2); p7->Rebin(2); p8->Rebin(2); p9->Rebin(2); p10->Rebin(2); double lowbound =0.8; double highbound =3.5; int lowbin = p5->FindBin(lowbound); int highbin = p5->FindBin(highbound); p5->Scale(1.0/p5->Integral(lowbin,highbin)); p6->Scale(1.0/p6->Integral(lowbin,highbin)); p7->Scale(1.0/p7->Integral(lowbin,highbin)); p8->Scale(1.0/p8->Integral(lowbin,highbin)); p9->Scale(1.0/p9->Integral(lowbin,highbin)); p10->Scale(1.0/p10->Integral(lowbin,highbin)); SetHStyle(p5,0.045); SetHStyle(p6,0.045); SetHStyle(p7,0.045); SetHStyle(p8,0.045); SetHStyle(p9,0.045); SetHStyle(p10,0.045); // p5->Draw("e"); // p6->SetLineColor(3); // p6->Draw("esame"); // p7->SetLineColor(4); // p7->Draw("esame"); // p8->SetLineColor(5); // p8->Draw("esame"); // p9->SetLineColor(2); // p9->Draw("esame"); // TF1* func = new TF1("func","landau(0)",3); TF1* func = new TF1("func",langaufun,lowbound,highbound,4); func->SetParameter(0,0.17); func->SetParameter(1,1.65); func->SetParameter(2,0.12); func->SetParameter(3,0.5); func->SetParameter(4,1.8); func->SetParameter(5,2.0); // func->FixParameter(0,0.0); // func->FixParameter(3,0.01); // func->FixParameter(4,-0.1); p5->GetXaxis()->SetRangeUser(0.2,4.2); p5->Fit("func","","e",lowbound,highbound); p5->GetFunction("func")->SetLineWidth(3); lwidth->SetBinContent(1,func->GetParameter(0)); peaks->SetBinContent (1,func->GetParameter(1)); gwidth->SetBinContent(1,func->GetParameter(3)); lwidth->SetBinError (1,func->GetParError(0)); peaks->SetBinError (1,func->GetParError(1)); gwidth->SetBinError (1,func->GetParError(3)); p6->SetLineColor(3); p6->SetMarkerColor(3); p6->Fit("func","","samee",lowbound,highbound); p6->GetFunction("func")->SetLineColor(3); p6->GetFunction("func")->SetLineWidth(3); lwidth->SetBinContent(2,func->GetParameter(0)); peaks->SetBinContent (2,func->GetParameter(1)); gwidth->SetBinContent(2,func->GetParameter(3)); lwidth->SetBinError (2,func->GetParError(0)); peaks->SetBinError (2,func->GetParError(1)); gwidth->SetBinError (2,func->GetParError(3)); p7->SetLineColor(4); p7->SetMarkerColor(4); p7->Fit("func","","samee",lowbound,highbound); p7->GetFunction("func")->SetLineColor(4); p7->GetFunction("func")->SetLineWidth(3); lwidth->SetBinContent(3,func->GetParameter(0)); peaks->SetBinContent (3,func->GetParameter(1)); gwidth->SetBinContent(3,func->GetParameter(3)); lwidth->SetBinError (3,func->GetParError(0)); peaks->SetBinError (3,func->GetParError(1)); gwidth->SetBinError (3,func->GetParError(3)); p8->SetLineColor(5); p8->SetMarkerColor(5); p8->Fit("func","","samee",lowbound,highbound); p8->GetFunction("func")->SetLineColor(5); p8->GetFunction("func")->SetLineWidth(3); lwidth->SetBinContent(4,func->GetParameter(0)); peaks->SetBinContent (4,func->GetParameter(1)); gwidth->SetBinContent(4,func->GetParameter(3)); lwidth->SetBinError (4,func->GetParError(0)); peaks->SetBinError (4,func->GetParError(1)); gwidth->SetBinError (4,func->GetParError(3)); p9->SetLineColor(2); p9->SetMarkerColor(2); p9->Fit("func","","samee",lowbound,highbound); p9->GetFunction("func")->SetLineColor(2); p9->GetFunction("func")->SetLineWidth(3); lwidth->SetBinContent(5,func->GetParameter(0)); peaks->SetBinContent (5,func->GetParameter(1)); gwidth->SetBinContent(5,func->GetParameter(3)); lwidth->SetBinError (5,func->GetParError(0)); peaks->SetBinError (5,func->GetParError(1)); gwidth->SetBinError (5,func->GetParError(3)); p10->SetLineColor(6); p10->SetMarkerColor(2); p10->Fit("func","","samee",lowbound,highbound); p10->GetFunction("func")->SetLineColor(6); p10->GetFunction("func")->SetLineWidth(3); lwidth->SetBinContent(6,func->GetParameter(0)); peaks->SetBinContent (6,func->GetParameter(1)); gwidth->SetBinContent(6,func->GetParameter(3)); lwidth->SetBinError (6,func->GetParError(0)); peaks->SetBinError (6,func->GetParError(1)); gwidth->SetBinError (6,func->GetParError(3)); TLegend* l = new TLegend(0.7,0.60,0.90,0.90,0,"NDC"); l->AddEntry(p5,"Period 5","l"); l->AddEntry(p6,"Period 6","l"); l->AddEntry(p7,"Period 7","l"); l->AddEntry(p8,"Period 8","l"); l->AddEntry(p9,"Period 9","l"); l->AddEntry(p10,"Period 10","l"); DrawLegend(l); c->SaveAs("EH_"+region+"_fits.eps"); c->SaveAs("EH_"+region+"_fits.gif"); SetStyle(); SetHStyle(peaks,0.045); peaks->Draw(); peaks->SetMinimum(1.5); peaks->SetMaximum(1.8); c->SaveAs("EH_"+region+"_peaks.eps"); c->SaveAs("EH_"+region+"_peaks.gif"); SetHStyle(lwidth,0.045); lwidth->SetMinimum(0.0); // lwidth->SetMaximum(0.0); lwidth->Draw(); c->SaveAs("EH_"+region+"_lwidth.eps"); c->SaveAs("EH_"+region+"_lwidth.gif"); SetHStyle(gwidth,0.045); gwidth->SetMinimum(0.0); // gwidth->SetMaximum(0.75); gwidth->Draw(); c->SaveAs("EH_"+region+"_gwidth.eps"); c->SaveAs("EH_"+region+"_gwidth.gif"); }