+++ /dev/null
-
-#include <TFile.h>
-#include <TH1.h>
-#include <TProfile.h>
-#include <TF1.h>
-#include <TGraphErrors.h>
-#include <TCanvas.h>
-#include <TAxis.h>
-#include <TLegend.h>
-#include <TStyle.h>
-#include <TGraphAsymmErrors.h>
-#include <TList.h>
-#include <TMath.h>
-#include <TString.h>
-#include <TSystem.h>
-#include <TLatex.h>
-#include <TF1.h>
-#include <TF2.h>
-#include <TFitResultPtr.h>
-#include <TFitResult.h>
-
-#include "AliHighPtDeDxCalib.h"
-
-#include "my_tools.C"
-#include "my_functions.C"
-
-#include <iostream>
-#include <fstream>
-#include <string>
-
-using namespace std;
-
-/*
- Ideas to improve:
- =================
-
- Use the real mean p -> Might give some effect
- Effect: Push the <dE/dx> down (let'a see if it is still needed in the fits)
-
- Use the real <nCl>
- Effect: Increase sigma for p<15. For p>15 seems to saturate.
-
- To use:
-
- =======
- root is enough
-
- gSystem->AddIncludePath("-I$ALICE_ROOT/TPC")
- gSystem->AddIncludePath("-I../lib")
- gSystem->AddIncludePath("-I../grid")
- gSystem->AddIncludePath("-I../macros")
- gROOT->SetMacroPath(".:../macros:../grid:../lib/")
- .L libMyDeDxAnalysis.so
- .L my_functions.C+
- .L my_tools.C+
- .L DebugClasses.C+
- .L calibrate_de_dx.C+
-
- // Step 2: test calibrations done in step 1
- DrawStep2("calib_eta/7tev_b.root", kFALSE);
- DrawStep2("calib_eta/7tev_b.root", kTRUE);
- TestResolutionVsDeDx("calib_eta/7tev_b.root");
- TestResolutionVsEta("calib_eta/7tev_b.root", kTRUE);
-
- // Here we want to see that the data is symmetric in eta and roughly symmetric for high vs low eta
- CompareYields("calib_eta/7tev_b.root", 0, 3.0, 30.0, "eta-80", "eta08")
- CompareYields("calib_eta/7tev_b.root", 0, 3.0, 30.0, "etaabs04", "etaabs48")
-
- // Step 3: extract
- FitDeDxVsP("calib_eta/7tev_b.root", 3.0, 30.0, 4, 2, kTRUE)
-
-
- // Step 4: inspect the results in results/calibdedx and results/calibplots
- // using e.g. gwenview. Event the fit did not converge it is still ok if it
- // describes the data!
-
-
-
- //
- // Test
- //
- // Step 2: test calibrations done in step 1
- DrawStep2("calib_eta/7tev_b_test.root", kFALSE);
- DrawStep2("calib_eta/7tev_b_test.root", kTRUE);
- TestResolutionVsDeDx("calib_eta/7tev_b_test.root");
- TestResolutionVsEta("calib_eta/7tev_b_test.root", kTRUE);
- // Step 3: extract
- FitDeDxVsP("calib_eta/7tev_b_test.root", 2.0, 10.0, 4, 2, kTRUE)
-
-
- */
-
-
-//____________________________________________________________________________
-void FitDeDxVsP(const Char_t* calibFileName,
- Double_t pStart, Double_t pStop,
- Int_t optionDeDx, Int_t optionSigma,
- Bool_t performFit = kFALSE,
- Int_t run = 0,
- Int_t filter = 1,
- Bool_t usePhiCut = kTRUE,
- const Char_t* v0FileName=0,
- Bool_t fixAllPar=kFALSE)
-{
- gStyle->SetOptStat(0);
-
- TString baseName(gSystem->BaseName(calibFileName));
- baseName.ReplaceAll(".root", "");
-
- TFile* calibFile = FindFileFresh(calibFileName);
- if(!calibFile)
- return;
- AliHighPtDeDxCalib* calib = (AliHighPtDeDxCalib*)GetObject(calibFile, filter, usePhiCut, run);
- calib->Print();
-
- fixMIP = calib->GetHistDeDx(kTRUE, 0)->GetMean();
- fixPlateau = calib->GetHistDeDx(kFALSE, 0)->GetMean();
- hDeDxVsP = calib->GetHistDeDxVsP(0);
- hMeanP = calib->GetHistMeanP();
-
- AliHighPtDeDxCalib* lambdaCalib = 0;
- AliHighPtDeDxCalib* kaonCalib = 0;
- TH2D* hDeDxVsPV0Lambda = 0;
- TH2D* hDeDxVsPV0Kaon = 0;
- if(v0FileName) {
- TFile* v0File = FindFileFresh(v0FileName);
- if(!v0File)
- return;
- lambdaCalib = (AliHighPtDeDxCalib*)GetObject(v0File, 0, usePhiCut, run, "lambda");
- lambdaCalib->Print();
- hDeDxVsPV0Lambda = lambdaCalib->GetHistDeDxVsP(0);
- kaonCalib = (AliHighPtDeDxCalib*)GetObject(v0File, 0, usePhiCut, run, "kaon");
- kaonCalib->Print();
- hDeDxVsPV0Kaon = kaonCalib->GetHistDeDxVsP(0);
- }
-
- CreateDir("old");
- gSystem->Exec("mv results/calibdedx/* old/");
- if(calib->IsMc())
- gSystem->Exec("mv debugmc/* old/");
-
- if(hDeDxVsPV0Lambda) {
- gSystem->Exec("mv debuglambda/* old/");
- gSystem->Exec("mv debugkaon/* old/");
- }
-
- TH2D* hDeDxVsPPi = 0;
- TH2D* hDeDxVsPK = 0;
- TH2D* hDeDxVsPP = 0;
- TH2D* hDeDxVsPMu = 0;
-
- if(calib->IsMc()) {
-
- hDeDxVsPPi = calib->GetHistDeDxVsP(1);
- hDeDxVsPK = calib->GetHistDeDxVsP(2);
- hDeDxVsPP = calib->GetHistDeDxVsP(3);
- hDeDxVsPMu = calib->GetHistDeDxVsP(5);
- }
-
- TCanvas* cDeDxVsP = new TCanvas("cDeDxVsP", "dE/dx vs p", 400, 300);
- cDeDxVsP->Clear();
- cDeDxVsP->cd();
- cDeDxVsP->SetLogz();
- hDeDxVsP->SetTitle(0);
- hDeDxVsP->Draw("COLZ");
-
- TCanvas* cDeDxVsPLogX = new TCanvas("cDeDxVsPLogX", "dE/dx vs p", 400, 300);
- cDeDxVsPLogX->Clear();
- cDeDxVsPLogX->cd();
- cDeDxVsPLogX->SetLogz();
- cDeDxVsPLogX->SetLogx();
- hDeDxVsP->Draw("COLZ");
-
- // Root is a bit stupid with finidng bins so we have to add and subtract a
- // little to be sure we get the right bin as we typically put edges as
- // limits
- const Int_t binStart = hDeDxVsP->GetXaxis()->FindBin(pStart+0.01);
- pStart = hDeDxVsP->GetXaxis()->GetBinLowEdge(binStart);
- const Int_t binStop = hDeDxVsP->GetXaxis()->FindBin(pStop-0.01);
- pStop = hDeDxVsP->GetXaxis()->GetBinUpEdge(binStop);
- const Int_t nBins = binStop - binStart + 1;
-
- cout << "Doing 2d fit from pTlow = " << pStart << " (bin: " << binStart
- << ") to pThigh = " << pStop << " (bin: " << binStop << ")" << endl;
-
- // Double_t binSize = (histdEdxvsP->GetXaxis()->GetXmax() - histdEdxvsP->GetXaxis()->GetXmin())/ histdEdxvsP->GetXaxis()->GetNbins();
-
- Double_t parDeDx[3] = {0, 0, 0};
- Double_t parSigma[3] = {0, 0, 0};
-
- const Int_t nLocalParams = 3; // pi, K, p yields
- Int_t nDeDxPar = 0;
- Int_t nSigmaPar = 0;
-
- switch(optionDeDx) {
-
- case 1:
- nDeDxPar = 2;
- parDeDx[0] = 39.7;
- parDeDx[1] = 6.3;
- break;
- case 2:
- nDeDxPar = 1;
- parDeDx[0] = 7.3;
- break;
- case 3:
- nDeDxPar = 2;
- parDeDx[0] = 6.85097;
- parDeDx[1] = 29.5933;
- break;
- // case 4:
- // nDeDxPar = 3;
- // parDeDx[0] = 35.5471;
- // parDeDx[1] = 6.85097;
- // parDeDx[2] = 29.5933;
- // break;
- case 4:
- nDeDxPar = 3;
- parDeDx[0] = 35.6694;
- parDeDx[1] = 6.80835;
- parDeDx[2] = 7.6542;
- break;
- default:
-
- cout << "optionDeDx does not support option: " << optionDeDx << endl;
- return;
- break;
- }
-
- switch(optionSigma) {
-
- case 1:
- nSigmaPar = 1;
- parSigma[0] = 3.0;
- break;
- case 2:
- nSigmaPar = 1;
- parSigma[0] = 0.0655688;
- break;
- case 3:
- nSigmaPar = 2;
- parSigma[0] = 0.06;
- parSigma[1] = -0.001;
- case 4:
- nSigmaPar = 2;
- parSigma[0] = 0.06;
- parSigma[1] = 1.0;
- break;
- case 5:
- nSigmaPar = 2;
- parSigma[0] = 0.06;
- parSigma[1] = 1.0;
- break;
- default:
-
- cout << "optionSigma does not support option: " << optionSigma << endl;
- return;
- break;
- }
-
- const Int_t nGlobalParams = 2 // binStart, binStop,
- + 2 + nDeDxPar // optionDeDx, nDeDxPar, dedxpar0 ....
- + 2 + nSigmaPar; // nSigmaPar, sigmapar0 .....
-
- const Int_t nParams = nBins*nLocalParams + nGlobalParams;
-
-
- TF2* fitAll = new TF2("fitAll", fitf3G, pStart, pStop, 30, 90, nParams);
- Double_t parametersIn[nParams];
-
- parametersIn[0] = binStart;
- fitAll->SetParName(0,"binStart");
- fitAll->FixParameter(0, parametersIn[0]);
-
- parametersIn[1] = binStop;
- fitAll->SetParName(1,"binStop");
- fitAll->FixParameter(1, parametersIn[1]);
-
- // dE/dx parameters
- parametersIn[2] = nDeDxPar;
- fitAll->SetParName(2,"nDeDxPar");
- fitAll->FixParameter(2, parametersIn[2]);
-
- parametersIn[3] = optionDeDx;
- fitAll->SetParName(3,"optionDeDx");
- fitAll->FixParameter(3, parametersIn[3]);
-
- for(Int_t i = 0; i < nDeDxPar; i++) {
-
- parametersIn[4+i] = parDeDx[i];
- fitAll->SetParName(4+i,Form("dEdXpar%d", i));
- // if(optionDeDx == 4 && i <2) {
-
- // fitAll->FixParameter(4+i, parametersIn[4+i]);
- // }
- if(fixAllPar) {
-
- fitAll->FixParameter(4+i, parametersIn[4+i]);
- }
- }
-
- // sigma parameters
- parametersIn[4+nDeDxPar] = nSigmaPar;
- fitAll->SetParName(4+nDeDxPar,"nSigmaPar");
- fitAll->FixParameter(4+nDeDxPar, parametersIn[4+nDeDxPar]);
-
- parametersIn[5+nDeDxPar] = optionSigma;
- fitAll->SetParName(5+nDeDxPar,"optionSigma");
- fitAll->FixParameter(5+nDeDxPar, parametersIn[5+nDeDxPar]);
-
- for(Int_t i = 0; i < nSigmaPar; i++) {
-
- parametersIn[6+nDeDxPar+i] = parSigma[i];
- fitAll->SetParName(6+nDeDxPar+i,Form("sigmapar%d", i));
- if(fixAllPar) {
-
- fitAll->FixParameter(6+nDeDxPar+i, parametersIn[6+nDeDxPar+i]);
- }
- }
-
- // Initial yields
-
- for(Int_t bin = binStart; bin <= binStop; bin++) {
-
- TH1D* hDeDxVsPProj =(TH1D*)hDeDxVsP->ProjectionY("hDeDxVsPProj", bin, bin);
-
- const Int_t offset = nGlobalParams + (bin - binStart)*nLocalParams;
- const Double_t all = hDeDxVsPProj->Integral();
- parametersIn[offset + 0] = (all)*0.6;
- parametersIn[offset + 1] = (all)*0.3;
- parametersIn[offset + 2] = (all)*0.1;
-
- fitAll->SetParName(offset + 0, Form("piYield%d", bin));
- fitAll->SetParLimits(offset + 0, 0, 10*all);
- fitAll->SetParName(offset + 1, Form("kYield%d", bin));
- fitAll->SetParLimits(offset + 1, 0, 10*all);
- fitAll->SetParName(offset + 2, Form("pYield%d", bin));
- fitAll->SetParLimits(offset + 2, 0, 10*all);
- // fitAll->SetParLimits(offset + 0, 0., histdEdxvsPpy->GetEntries());
- // fitAll->SetParLimits(offset + 1, 0., histdEdxvsPpy->GetEntries());
- // fitAll->SetParLimits(offset + 2, 0., histdEdxvsPpy->GetEntries());
- }
-
- fitAll->SetParameters(parametersIn);
- fitAll->Print();
-
- Bool_t converged = kFALSE;
-
- if(performFit) {
- for(Int_t i = 0; i < 4; i++) {
- TFitResultPtr fitResultPtr = hDeDxVsP->Fit(fitAll, "0S", "", pStart+0.01, pStop-0.01);
- if(!fitResultPtr->Status()) {
- // if(!fitResultPtr->Status() && !fitResultPtr->CovMatrixStatus()) {
-
- converged = kTRUE;
- break;
- }
- }
- }
- // else we just draw how the results looks with the input parameters
-
- Double_t parametersOut[nParams];
- fitAll->GetParameters(parametersOut);
- const Double_t* parameterErrorsOut = fitAll->GetParErrors();
-
- // overlay the fitfunction
-
-
- TF1* fDeDxPi = new TF1("fDeDxPi", FitFunc, 0, 50, nDeDxPar+1); // +1 because of option!
- fDeDxPi->SetParameters(¶metersOut[3]);
- fDeDxPi->SetLineWidth(2);
- cDeDxVsP->cd();
- fDeDxPi->Draw("SAME");
- cDeDxVsPLogX->cd();
- fDeDxPi->Draw("SAME");
-
- TF1* fDeDxK = new TF1("fDeDxK", FitFunc, 0, 50, nDeDxPar+1);
- fDeDxK->SetParameters(¶metersOut[3]);
- fDeDxK->SetParameter(0, fDeDxK->GetParameter(0)+10);
- fDeDxK->SetLineWidth(2);
- cDeDxVsP->cd();
- fDeDxK->Draw("SAME");
- cDeDxVsPLogX->cd();
- fDeDxK->Draw("SAME");
-
- TF1* fDeDxP = new TF1("fDeDxP", FitFunc, 0, 50, nDeDxPar+1);
- fDeDxP->SetParameters(¶metersOut[3]);
- fDeDxP->SetParameter(0, fDeDxP->GetParameter(0)+20);
- fDeDxP->SetLineWidth(2);
- cDeDxVsP->cd();
- fDeDxP->Draw("SAME");
- cDeDxVsPLogX->cd();
- fDeDxP->Draw("SAME");
-
- TF1* fDeDxE = new TF1("fDeDxE", FitFunc, 0, 50, nDeDxPar+1);
- fDeDxE->SetParameters(¶metersOut[3]);
- fDeDxE->SetParameter(0, fDeDxE->GetParameter(0)+30);
- fDeDxE->SetLineWidth(2);
- cDeDxVsP->cd();
- fDeDxE->Draw("SAME");
- cDeDxVsPLogX->cd();
- fDeDxE->Draw("SAME");
-
- TF1* fSigma = new TF1("fSigma", SigmaFunc, 0, 50, nSigmaPar+1);
- fSigma->SetParameters(¶metersOut[5 + nDeDxPar]);
-
- // fitAll->Draw("same cont3");
-
- CreateDir("results/calibdedx");
-
- cDeDxVsP->cd();
- cDeDxVsP->Modified();
- cDeDxVsP->Update();
- gROOT->ProcessLine(".x drawText.C");
- cDeDxVsP->SaveAs("results/calibdedx/dedx_vs_p.gif");
- cDeDxVsP->SaveAs("results/calibdedx/dedx_vs_p.pdf");
-
- cDeDxVsPLogX->cd();
- cDeDxVsPLogX->Modified();
- cDeDxVsPLogX->Update();
- gROOT->ProcessLine(".x drawText.C");
- cDeDxVsPLogX->SaveAs("results/calibdedx/dedx_vs_p_logx.gif");
- cDeDxVsPLogX->SaveAs("results/calibdedx/dedx_vs_p_logx.pdf");
-
- //cross check
- TCanvas* cFits = new TCanvas("cFits", "Fit comparison to data", 1200, 800);
- cFits->Clear();
- cFits->Divide(7, 5);
-
- TF1* pion = new TF1("pion", "gausn", 30, 90);
- pion->SetLineWidth(2);
- pion->SetLineColor(kRed);
- TF1* kaon = new TF1("kaon", "gausn", 30, 90);
- kaon->SetLineWidth(2);
- kaon->SetLineColor(kGreen);
- TF1* proton = new TF1("proton", "gausn", 30, 90);
- proton->SetLineWidth(2);
- proton->SetLineColor(kBlue);
- TF1* total = new TF1("total", "gausn(0)+gausn(3)+gausn(6)", 30, 90);
- total->SetLineColor(kBlack);
- total->SetLineWidth(2);
- total->SetLineStyle(2);
-
- TLegend* legend = new TLegend(0.11, 0.6, 0.35, 0.85);
- legend->SetBorderSize(0);
- legend->SetFillColor(0);
- legend->AddEntry(total, "3-Gauss fit", "L");
- legend->AddEntry(pion, "#pi", "L");
- legend->AddEntry(kaon, "K", "L");
- legend->AddEntry(proton, "p", "L");
-
-
- TCanvas* cSingleFit = new TCanvas("cSingleFit", "single fit");
-
- TH1D* hPionRatio =(TH1D*)hDeDxVsP->ProjectionX("hPionRatio", 1, 1);
- hPionRatio->Reset();
- hPionRatio->GetXaxis()->SetRangeUser(pStart+0.001, pStop-0.001);
- hPionRatio->GetYaxis()->SetRangeUser(0.0, 1.0);
- hPionRatio->SetTitle("particle fractions; p [GeV/c]; particle fraction");
- TH1D* hKaonRatio = (TH1D*)hPionRatio->Clone("hKaonRatio");
- TH1D* hProtonRatio = (TH1D*)hPionRatio->Clone("hProtonRatio");
-
- TH1D* hPionRatioMc = (TH1D*)hPionRatio->Clone("hPionRatioMc");
- TH1D* hKaonRatioMc = (TH1D*)hPionRatio->Clone("hKaonRatioMc");
- TH1D* hProtonRatioMc = (TH1D*)hPionRatio->Clone("hProtonRatioMc");
- TH1D* hMuonRatioMc = (TH1D*)hPionRatio->Clone("hMuonRatioMc");
-
- for(Int_t bin = binStart; bin <= binStop; bin++){
-
- cout << "Making projection for bin: " << bin << endl;
-
- const Int_t j = bin-binStart;
- cFits->cd();
- cFits->cd(j + 1);
-
- TH1D* hDeDxVsPProj =(TH1D*)hDeDxVsP->ProjectionY(Form("hDeDxVsPProj%d", bin), bin, bin);
- hDeDxVsPProj->GetXaxis()->SetRangeUser(40, 90);
- hDeDxVsPProj->SetTitle(Form("%.1f<p<%.1f GeV/c",
- hDeDxVsP->GetXaxis()->GetBinLowEdge(bin),
- hDeDxVsP->GetXaxis()->GetBinUpEdge(bin)));
- hDeDxVsPProj->Draw();
-
- const Int_t offset = nGlobalParams + j*nLocalParams;
- const Double_t p = hDeDxVsP->GetXaxis()->GetBinCenter(bin);
- const Double_t pKeff = p*piMass/kMass; // corresponding p of a pion with same dE/dx
- const Double_t pPeff = p*piMass/pMass; // corresponding p of a pion with same dE/dx
- const Double_t meanDeDxPi = fDeDxPi->Eval(p);
- const Double_t meanDeDxK = fDeDxPi->Eval(pKeff);
- const Double_t meanDeDxP = fDeDxPi->Eval(pPeff);
- Double_t gausParams[9] = {
- parametersOut[offset + 0],
- meanDeDxPi,
- fSigma->Eval(meanDeDxPi) ,
- parametersOut[offset + 1],
- meanDeDxK,
- fSigma->Eval(meanDeDxK) ,
- parametersOut[offset + 2],
- meanDeDxP,
- fSigma->Eval(meanDeDxP) ,
- };
-
- for(Int_t i = 0; i < 9; i++)
- cout << gausParams[i] << ", ";
-
- cout << endl;
-
- pion->SetParameters(&gausParams[0]);
- pion->DrawCopy("same");
- Double_t all = hDeDxVsPProj->Integral();
- hPionRatio->SetBinContent(bin, parametersOut[offset + 0]/all);
- hPionRatio->SetBinError(bin, parameterErrorsOut[offset + 0]/all);
-
- kaon->SetParameters(&gausParams[3]);
- kaon->DrawCopy("same");
- hKaonRatio->SetBinContent(bin, parametersOut[offset + 1]/all);
- hKaonRatio->SetBinError(bin, parameterErrorsOut[offset + 1]/all);
-
- proton->SetParameters(&gausParams[6]);
- proton->DrawCopy("same");
- hProtonRatio->SetBinContent(bin, parametersOut[offset + 2]/all);
- hProtonRatio->SetBinError(bin, parameterErrorsOut[offset + 2]/all);
-
- total->SetParameters(gausParams);
- total->DrawCopy("same");
-
- cSingleFit->cd();
- cSingleFit->Clear();
- // cSingleFit->SetLogy();
- hDeDxVsPProj->Draw();
- pion->DrawCopy("same");
- kaon->DrawCopy("same");
- proton->DrawCopy("same");
- total->DrawCopy("same");
-
- gROOT->ProcessLine(".x drawText.C(2)");
- cSingleFit->SaveAs(Form("results/calibdedx/ptspectrum_bin%d.gif", bin));
- cSingleFit->SaveAs(Form("results/calibdedx/ptspectrum_bin%d.pdf", bin));
- // legend->Draw();
-
- if(calib->IsMc()) {
-
- cSingleFit->cd();
- cSingleFit->Clear();
- TH1D* hDeDxVsPPiProj =(TH1D*)hDeDxVsPPi->ProjectionY(Form("hDeDxVsPPiProj%d", bin), bin, bin);
- hDeDxVsPPiProj->SetMarkerStyle(20);
- hDeDxVsPPiProj->SetMarkerColor(2);
- hDeDxVsPPiProj->GetXaxis()->SetRangeUser(40, 90);
- hDeDxVsPPiProj->SetTitle(Form("%.1f<p<%.1f GeV/c",
- hDeDxVsP->GetXaxis()->GetBinLowEdge(bin),
- hDeDxVsP->GetXaxis()->GetBinUpEdge(bin)));
- hDeDxVsPPiProj->Draw("P");
- hPionRatioMc->SetBinContent(bin, hDeDxVsPPiProj->Integral()/all);
- TH1D* hDeDxVsPKProj =(TH1D*)hDeDxVsPK->ProjectionY(Form("hDeDxVsPKProj%d", bin), bin, bin);
- hDeDxVsPKProj->SetMarkerStyle(21);
- hDeDxVsPKProj->SetMarkerColor(3);
- hDeDxVsPKProj->Draw("SAME P");
- hKaonRatioMc->SetBinContent(bin, hDeDxVsPKProj->Integral()/all);
- TH1D* hDeDxVsPPProj =(TH1D*)hDeDxVsPP->ProjectionY(Form("hDeDxVsPPProj%d", bin), bin, bin);
- hDeDxVsPPProj->SetMarkerStyle(22);
- hDeDxVsPPProj->SetMarkerColor(4);
- hDeDxVsPPProj->Draw("SAME P");
- hProtonRatioMc->SetBinContent(bin, hDeDxVsPPProj->Integral()/all);
- TH1D* hDeDxVsPMuProj =(TH1D*)hDeDxVsPMu->ProjectionY(Form("hDeDxVsPMuProj%d", bin), bin, bin);
- hDeDxVsPMuProj->SetMarkerStyle(22);
- hDeDxVsPMuProj->SetMarkerColor(6);
- // hDeDxVsPMuProj->Draw("SAME P");
- hMuonRatioMc->SetBinContent(bin, hDeDxVsPMuProj->Integral()/all);
- // cSingleFit->SetLogy();
- pion->DrawCopy("same");
- kaon->DrawCopy("same");
- proton->DrawCopy("same");
- CreateDir("results/calibdedx/debugmc");
- cSingleFit->SaveAs(Form("results/calibdedx/debugmc/ptspectrum_bin%d.gif", bin));
- }
-
- if(hDeDxVsPV0Lambda) {
- cSingleFit->cd();
- cSingleFit->Clear();
- TH1D* hDeDxVsPV0LambdaProj =(TH1D*)hDeDxVsPV0Lambda->ProjectionY(Form("hDeDxVsPV0LambdaProj%d", bin), bin, bin);
- hDeDxVsPV0LambdaProj->SetMarkerStyle(20);
- hDeDxVsPV0LambdaProj->SetMarkerColor(4);
- hDeDxVsPV0LambdaProj->GetXaxis()->SetRangeUser(40, 90);
- hDeDxVsPV0LambdaProj->SetTitle(Form("%.1f<p<%.1f GeV/c",
- hDeDxVsP->GetXaxis()->GetBinLowEdge(bin),
- hDeDxVsP->GetXaxis()->GetBinUpEdge(bin)));
- hDeDxVsPV0LambdaProj->Draw("P");
- proton->SetParameter(0, hDeDxVsPV0LambdaProj->Integral());
- proton->DrawCopy("same");
-
- CreateDir("results/calibdedx/debuglambda");
- cSingleFit->SaveAs(Form("results/calibdedx/debuglambda/ptspectrum_bin%d.gif", bin));
- }
-
- if(hDeDxVsPV0Kaon) {
- cSingleFit->cd();
- cSingleFit->Clear();
- TH1D* hDeDxVsPV0KaonProj =(TH1D*)hDeDxVsPV0Kaon->ProjectionY(Form("hDeDxVsPV0KaonProj%d", bin), bin, bin);
- hDeDxVsPV0KaonProj->SetMarkerStyle(20);
- hDeDxVsPV0KaonProj->SetMarkerColor(2);
- hDeDxVsPV0KaonProj->GetXaxis()->SetRangeUser(40, 90);
- hDeDxVsPV0KaonProj->SetTitle(Form("%.1f<p<%.1f GeV/c",
- hDeDxVsP->GetXaxis()->GetBinLowEdge(bin),
- hDeDxVsP->GetXaxis()->GetBinUpEdge(bin)));
- hDeDxVsPV0KaonProj->Draw("P");
- pion->SetParameter(0, hDeDxVsPV0KaonProj->Integral());
- pion->DrawCopy("same");
- CreateDir("results/calibdedx/debugkaon");
- cSingleFit->SaveAs(Form("results/calibdedx/debugkaon/ptspectrum_bin%d.gif", bin));
- }
- }
-
- TCanvas* cRatio = new TCanvas("cRatio", "ratios/all vs p", 600, 400);
- cRatio->Clear();
- hPionRatio->SetMaximum(0.8);
- hPionRatio->SetMarkerStyle(20);
- hPionRatio->SetMarkerColor(2);
- hPionRatio->Draw("P E");
-
- hKaonRatio->SetMarkerStyle(20);
- hKaonRatio->SetMarkerColor(3);
- hKaonRatio->Draw("SAME P E");
-
- hProtonRatio->SetMarkerStyle(20);
- hProtonRatio->SetMarkerColor(4);
- hProtonRatio->Draw("SAME P E");
- gROOT->ProcessLine(".x drawText.C(2)");
- cRatio->SaveAs("results/calibdedx/particle_ratios.gif");
- cRatio->SaveAs("results/calibdedx/particle_ratios.pdf");
-
- if(calib->IsMc()) {
-
- hPionRatioMc->SetMarkerStyle(24);
- hPionRatioMc->SetMarkerColor(2);
- hPionRatioMc->Draw("SAME P");
-
- hKaonRatioMc->SetMarkerStyle(24);
- hKaonRatioMc->SetMarkerColor(3);
- hKaonRatioMc->Draw("SAME P");
-
- hProtonRatioMc->SetMarkerStyle(24);
- hProtonRatioMc->SetMarkerColor(4);
- hProtonRatioMc->Draw("SAME P");
-
- hMuonRatioMc->SetMarkerStyle(24);
- hMuonRatioMc->SetMarkerColor(6);
- // hMuonRatioMc->Draw("SAME P");
- cRatio->SaveAs("results/calibdedx/debugmc/particle_ratios_mc.gif");
- }
-
-
- cout << "MIP was fixed to: " << fixMIP << endl
- << "Plateau was fixed to: " << fixPlateau << endl;
-
- //
- // Store the <dE/dx> parameters in a file that we can get them back to use for the Delta-pi!
- //
- DeDxFitInfo* fitInfo = new DeDxFitInfo();
- fitInfo->MIP = fixMIP;
- fitInfo->plateau = fixPlateau;
- fitInfo->optionDeDx = optionDeDx;
- fitInfo->nDeDxPar = nDeDxPar;
- for(Int_t i = 0; i < nDeDxPar; i++) {
- fitInfo->parDeDx[i] = fDeDxPi->GetParameter(i+1); // 1st parameter is option
- }
- fitInfo->optionSigma = optionSigma;
- fitInfo->nSigmaPar = nSigmaPar;
- for(Int_t i = 0; i < nSigmaPar; i++) {
- fitInfo->parSigma[i] = fSigma->GetParameter(i+1); // 1st parameter is option
- }
- if(!strstr(calibFileName, "no_calib_eta"))
- fitInfo->calibFileName = calibFileName;
-
- fitInfo->Print();
-
- CreateDir("fitparameters");
- TFile* outFile = new TFile(Form("fitparameters/%s", gSystem->BaseName(calibFileName)), "RECREATE");
- fitInfo->Write("fitInfo");
- outFile->Close();
-
- if(converged) {
-
- cout << "Fit converged and error matrix was ok" << endl;
- } else {
-
- cout << "WARNING!!!!!!!!!!!!!!! Fit did not converge, or error matrix was not ok!" << endl;
- }
-}
-
-
-
-//____________________________________________________________________________
-void DrawStep2(const Char_t* calibFileName,
- Bool_t forMIP = kTRUE,
- Int_t run = 0,
- Int_t filter = 1,
- Bool_t usePhiCut = kTRUE)
-{
- // option:
- // 0 - compare to ALICE INEL
- // 1 - compare to MC input spectrum
- // 2 - compare to MC correct spectrum
- gStyle->SetOptStat(0);
-
- CreateDir("results/calibplots");
- TString baseName(gSystem->BaseName(calibFileName));
- baseName.ReplaceAll(".root", "");
-
- TFile* calibFile = FindFileFresh(calibFileName);
- if(!calibFile)
- return;
- AliHighPtDeDxCalib* calib = (AliHighPtDeDxCalib*)GetObject(calibFile, filter, usePhiCut, run);
- calib->Print();
-
- TCanvas* cPhiCut = calib->DrawPhiCutHistograms();
- gROOT->ProcessLine(".x drawText.C(2)");
- cPhiCut->SaveAs(Form("results/calibplots/phicut_%s.gif", baseName.Data()));
- cPhiCut->SaveAs(Form("results/calibplots/phicut_%s.pdf", baseName.Data()));
- TCanvas* cEta = calib->DrawEta(forMIP);
- gROOT->ProcessLine(".x drawText.C(2)");
- TCanvas* cEtaCal = calib->DrawEtaCalibrated(forMIP);
- gROOT->ProcessLine(".x drawText.C(2)");
- if(forMIP) {
- cEta->cd();
- cEta->SaveAs(Form("results/calibplots/MIP_dedx_vs_eta_nocal_%s.gif", baseName.Data()));
- cEta->SaveAs(Form("results/calibplots/MIP_dedx_vs_eta_nocal_%s.pdf", baseName.Data()));
- cEtaCal->cd();
- cEtaCal->SaveAs(Form("results/calibplots/MIP_dedx_vs_eta_final_%s.gif", baseName.Data()));
- cEtaCal->SaveAs(Form("results/calibplots/MIP_dedx_vs_eta_final_%s.pdf", baseName.Data()));
- } else {
- cEta->cd();
- cEta->SaveAs(Form("results/calibplots/Electrons_dedx_vs_eta_nocal_%s.gif", baseName.Data()));
- cEta->SaveAs(Form("results/calibplots/Electrons_dedx_vs_eta_nocal_%s.pdf", baseName.Data()));
- cEtaCal->cd();
- cEtaCal->SaveAs(Form("results/calibplots/Electrons_dedx_vs_eta_final_%s.gif", baseName.Data()));
- cEtaCal->SaveAs(Form("results/calibplots/Electrons_dedx_vs_eta_final_%s.pdf", baseName.Data()));
- }
- TCanvas* cEtaSel = calib->DrawSelectionHistograms();
- gROOT->ProcessLine(".x drawText.C(2)");
- cEtaSel->SaveAs(Form("results/calibplots/selection_%s.gif", baseName.Data()));
- cEtaSel->SaveAs(Form("results/calibplots/selection_%s.pdf", baseName.Data()));
-
- // TCanvas* cNcl = calib->DrawNclCal();
- // gROOT->ProcessLine(".x drawText.C");
- // cNcl->SaveAs(Form("results/calibplots/dedx_vs_ncl_calib__%s.gif", baseName.Data()));
- // cNcl->SaveAs(Form("results/calibplots/dedx_vs_ncl_calib__%s.pdf", baseName.Data()));
-}
-
-//____________________________________________________________________________
-void TestResolutionVsDeDx(const Char_t* calibFileName,
- Int_t run = 0,
- Int_t filter = 1,
- Bool_t usePhiCut = kTRUE)
-{
- gStyle->SetOptStat(0);
-
- TFile* calibFile = FindFileFresh(calibFileName);
- if(!calibFile)
- return;
- AliHighPtDeDxCalib* calib = (AliHighPtDeDxCalib*)GetObject(calibFile, filter, usePhiCut, run);
- calib->Print();
-
- TH2D* hDeDxVsNcl = calib->GetHistDeDxVsNcl(kTRUE, 0);
- TH2D* hDeDxVsNclE = calib->GetHistDeDxVsNcl(kFALSE, 0);
- TH1D* hDeDx = calib->GetHistDeDx(kTRUE, 0);
- TH1D* hDeDxE = calib->GetHistDeDx(kFALSE, 0);
-
- Double_t mdedx = hDeDx->GetMean();
- Double_t mdedxE = hDeDxE->GetMean();
-
- TObjArray* helpArray = new TObjArray();
-
- hDeDxVsNcl->FitSlicesY(0, 0, -1, 0, "QNR", helpArray);
- TH1D* hMeanVsNcl = (TH1D*)helpArray->At(1);
- hMeanVsNcl->SetName("hMeanVsNcl");
- hMeanVsNcl->SetMarkerStyle(29);
- TH1D* hSigmaVsNcl = (TH1D*)helpArray->At(2);
- hSigmaVsNcl->SetName("hSigmaVsNcl");
- hSigmaVsNcl->SetTitle("#sigma vs Ncl; Ncl; #sigma");
- hSigmaVsNcl->SetMarkerStyle(29);
- hSigmaVsNcl->SetMarkerColor(2);
- hSigmaVsNcl->Scale(1.0/mdedx);
-
- hDeDxVsNclE->FitSlicesY(0, 0, -1, 0, "QNR", helpArray);
- TH1D* hMeanVsNclE = (TH1D*)helpArray->At(1);
- hMeanVsNclE->SetName("hSigmaVsNclE");
- hMeanVsNclE->SetMarkerStyle(29);
- TH1D* hSigmaVsNclE = (TH1D*)helpArray->At(2);
- hSigmaVsNclE->SetName("hSigmaVsNclE");
- hSigmaVsNclE->SetTitle("#sigma vs Ncl; Ncl; #sigma");
- hSigmaVsNclE->SetMarkerStyle(29);
- hSigmaVsNclE->SetMarkerColor(kMagenta);
- hSigmaVsNclE->Scale(1.0/mdedxE);
-
- TLegend* legend = new TLegend(0.54, 0.67, 0.84, 0.87);
- legend->SetFillColor(0);
- legend->SetBorderSize(0);
- legend->SetTextSize(0.05);
- legend->AddEntry(hSigmaVsNcl, "MIP pions", "P");
- legend->AddEntry(hSigmaVsNclE, "electrons", "P");
-
- TCanvas* cTestDeDx = new TCanvas("cTestDeDx", "Comparing resolution for MIPs and electrons",
- 1200, 400);
- cTestDeDx->Clear();
- cTestDeDx->Divide(3, 1);
-
- cTestDeDx->cd(1);
- hDeDxVsNcl->Draw("COL");
- hMeanVsNcl->Draw("SAME");
-
- cTestDeDx->cd(2);
- hDeDxVsNclE->Draw("COL");
- hMeanVsNclE->Draw("SAME");
-
- cTestDeDx->cd(3);
- TF1* fSDedxVsNclSqrt = new TF1("fSDedxVsNclSqrt", "[0]*sqrt(159.0/x)", 0, 160);
- fSDedxVsNclSqrt->SetParameter(0, 0.05);
- hSigmaVsNcl->SetMinimum(0.0);
- hSigmaVsNcl->SetMaximum(0.15);
- hSigmaVsNcl->SetTitle("#sigma/<dE/dx> vs Ncl; Ncl; rel. #sigma");
- // hSigmaVsNcl->Fit(fSDedxVsNclSqrt);
- hSigmaVsNcl->Draw();
- hSigmaVsNclE->Draw("SAME");
- legend->Draw();
- gROOT->ProcessLine(".x drawText.C(2)");
- TString baseName(gSystem->BaseName(calibFileName));
- baseName.ReplaceAll(".root", "");
- cTestDeDx->SaveAs(Form("results/calibplots/test_resolution_vs_dedx_%s.gif", baseName.Data()));
- cTestDeDx->SaveAs(Form("results/calibplots/test_resolution_vs_dedx_%s.pdf", baseName.Data()));
-}
-
-//___________________________________________________________________________________________
-void TestResolutionVsEta(const Char_t* calibFileName,
- Bool_t forMIP,
- Int_t run = 0,
- Int_t filter = 1,
- Bool_t usePhiCut = kTRUE)
-{
- gStyle->SetOptStat(0);
-
- TFile* calibFile = FindFileFresh(calibFileName);
- if(!calibFile)
- return;
- AliHighPtDeDxCalib* calib = (AliHighPtDeDxCalib*)GetObject(calibFile, filter, usePhiCut, run);
- calib->Print();
-
- const Int_t nEta = 4;
- Int_t color[nEta] = {kBlack, kBlue, kGreen+2, kRed};
- TH1D* hDeDx[nEta] = {0, 0, 0, 0};
- TH1D* hMeanVsNcl[nEta] = {0, 0, 0, 0};
- TH1D* hSigmaVsNcl[nEta] = {0, 0, 0, 0};
-
- const Double_t mdedx = (calib->GetHistDeDx(forMIP, 0))->GetMean();
-
- TObjArray* helpArray = new TObjArray();
-
- for(Int_t bin = 1; bin <= nEta; bin++) {
-
- const Int_t i = bin - 1;
- hDeDx[i] = calib->GetHistDeDx(forMIP, bin);
- hDeDx[i]->Scale(1.0/hDeDx[i]->GetEntries());
- hDeDx[i]->SetMarkerStyle(24+i);
- hDeDx[i]->SetMarkerColor(color[i]);
-
- TH2D* hDeDxVsNcl = calib->GetHistDeDxVsNcl(forMIP, bin);
- hDeDxVsNcl->FitSlicesY(0, 0, -1, 0, "QNR", helpArray);
-
- hMeanVsNcl[i] = (TH1D*)helpArray->At(1);
- hMeanVsNcl[i]->SetTitle("<dE/dx> vs Ncl; Ncl; #sigma");
- hMeanVsNcl[i]->SetName(Form("hMeanVsNcl%d", bin));
- hMeanVsNcl[i]->SetMarkerStyle(24+i);
- hMeanVsNcl[i]->SetMarkerColor(color[i]);
- hMeanVsNcl[i]->SetMinimum(mdedx*0.9);
- hMeanVsNcl[i]->SetMaximum(mdedx*1.1);
-
- hSigmaVsNcl[i] = (TH1D*)helpArray->At(2);
- hSigmaVsNcl[i]->SetName("hSigmaVsNcl");
- hSigmaVsNcl[i]->SetTitle("#sigma/<dE/dx> vs Ncl; Ncl; rel. #sigma");
- hSigmaVsNcl[i]->SetMarkerStyle(24+i);
- hSigmaVsNcl[i]->SetMarkerColor(color[i]);
- hSigmaVsNcl[i]->Scale(1.0/mdedx);
- hSigmaVsNcl[i]->SetMinimum(0.0);
- hSigmaVsNcl[i]->SetMaximum(0.15);
- }
-
- TLegend* legend = new TLegend(0.54, 0.67, 0.84, 0.87);
- legend->SetFillColor(0);
- legend->SetBorderSize(0);
- legend->SetTextSize(0.05);
-
- legend->AddEntry(hDeDx[0], "|#eta|<0.2", "P");
- legend->AddEntry(hDeDx[1], "0.2<|#eta|<0.4", "P");
- legend->AddEntry(hDeDx[2], "0.4<|#eta|<0.6", "P");
- legend->AddEntry(hDeDx[3], "0.6<|#eta|<0.8", "P");
-
- TCanvas* cTestEta = new TCanvas("cTestEta", "Comparing resolution for MIPs and electrons",
- 1200, 400);
- cTestEta->Clear();
- cTestEta->Divide(3, 1);
-
- for(Int_t i = 0; i < nEta; i++) {
-
- cTestEta->cd(1);
- if(i==0)
- hDeDx[i]->Draw("P");
- else
- hDeDx[i]->Draw("SAME P");
-
- cTestEta->cd(2);
- if(i==0)
- hMeanVsNcl[i]->Draw("P");
- else
- hMeanVsNcl[i]->Draw("SAME P");
-
- cTestEta->cd(3);
- if(i==0) {
- hSigmaVsNcl[i]->Draw("P");
- legend->Draw();
- } else
- hSigmaVsNcl[i]->Draw("SAME P");
- }
-
- gROOT->ProcessLine(".x drawText.C(2)");
- TString baseName(gSystem->BaseName(calibFileName));
- baseName.ReplaceAll(".root", "");
- cTestEta->SaveAs(Form("results/calibplots/test_resolution_vs_eta_%s.gif", baseName.Data()));
- cTestEta->SaveAs(Form("results/calibplots/test_resolution_vs_eta_%s.pdf", baseName.Data()));
-}
-
- // hMeanVsNcl = (TH1D*)helpArray->At(1);
- // hMeanVsNcl->SetName("hMeanVsNcl");
- // hMeanVsNcl->SetTitle("Mean vs Ncl; Ncl; Mean");
- // hMeanVsNcl->SetMarkerStyle(29);
- // hMeanVsNcl->SetMinimum(45);
- // hMeanVsNcl->SetMaximum(55);
-
- // if(isMC) {
- // hDeDxVsNclElectrons->FitSlicesY(0, 0, -1, 0, "QNR", helpArray);
- // hSigmaVsNclElectrons = (TH1D*)helpArray->At(2);
- // hSigmaVsNclElectrons->SetName("hSigmaVsNcl");
- // hSigmaVsNclElectrons->SetTitle("#sigma vs Ncl; Ncl; #sigma");
- // }
-
- // TCanvas* cNcl2 = new TCanvas("cNcl2", "sigma dE/dx vs Ncl", 600, 400);
- // cNcl2->Clear();
- // cNcl2->cd();
-
- // hSigmaVsNcl->DrawCopy();
- // hSigmaVsNcl->Fit(fSDeDxVsNcl, "0");
- // hSigmaVsNcl->Fit(fSDeDxVsNclSqrt, "0");
-
- // fSDeDxVsNcl->DrawCopy("SAME");
- // fSDeDxVsNclSqrt->SetLineStyle(2);
- // fSDeDxVsNclSqrt->DrawCopy("SAME");
-
- // TLegend* legFit = new TLegend(0.64, 0.67, 0.84, 0.87);
- // legFit->SetFillColor(0);
- // legFit->SetBorderSize(0);
- // legFit->SetTextSize(0.05);
- // legFit->AddEntry(fSDeDxVsNcl, "Linear fit", "L");
- // legFit->AddEntry(fSDeDxVsNclSqrt, "1/Sqrt fit", "L");
- // legFit->Draw();
-
- // cNcl2->SaveAs(Form("%s/sigma_vs_ncl.gif", baseName.Data()));
-
- // if(isMC) {
- // TCanvas* cNcl2Electrons = new TCanvas("cNcl2Electrons", "sigma dE/dx vs Ncl (comparison)", 600, 400);
- // cNcl2Electrons->Clear();
- // cNcl2Electrons->cd();
-
- // hSigmaVsNcl->SetMarkerStyle(29);
- // hSigmaVsNcl->DrawCopy();
- // hSigmaVsNclElectrons->Scale(50.0 / hDeDxVsEtaCalibratedElectrons->GetMean(2));
- // hSigmaVsNclElectrons->SetMarkerStyle(29);
- // hSigmaVsNclElectrons->SetMarkerColor(6);
- // hSigmaVsNclElectrons->DrawCopy("SAME");
- // fSDeDxVsNcl->DrawCopy("SAME");
- // fSDeDxVsNclSqrt->DrawCopy("SAME");
-
- // TLegend* legSigma = new TLegend(0.64, 0.67, 0.84, 0.87);
- // legSigma->SetFillColor(0);
- // legSigma->SetBorderSize(0);
- // legSigma->SetTextSize(0.05);
- // legSigma->AddEntry(hSigmaVsNcl, "MIP pions", "P");
- // legSigma->AddEntry(fSDeDxVsNcl, "Linear fit", "L");
- // legSigma->AddEntry(fSDeDxVsNclSqrt, "1/Sqrt fit", "L");
- // legSigma->AddEntry(hSigmaVsNclElectrons, "Scaled electrons", "P");
- // legSigma->Draw();
-
- // cNcl2Electrons->SaveAs(Form("%s/electrons_sigma_vs_ncl_comparison.gif", baseName.Data()));
- // }
-
- // TCanvas* cDeDx = new TCanvas("cDeDx", "dE/dx", 600, 400);
- // cDeDx->Clear();
- // hDeDx->Scale(1.0/hDeDx->GetEntries());
- // hDeDx->SetMarkerStyle(29);
- // hDeDx->DrawCopy();
- // TCanvas* cNcl3 = new TCanvas("cNcl3", "#sigma vs Ncl", 600, 400);
- // cNcl3->Clear();
- // hSigmaVsNcl->DrawCopy();
- // TCanvas* cNcl3Mean = new TCanvas("cNcl3Mean", "Mean vs Ncl", 600, 400);
- // cNcl3Mean->Clear();
- // hMeanVsNcl->DrawCopy();
-
- // TLegend* legDeDx = new TLegend(0.7, 0.67, 0.9, 0.87);
- // legDeDx->SetFillColor(0);
- // legDeDx->SetBorderSize(0);
- // legDeDx->SetTextSize(0.05);
- // legDeDx->AddEntry(hDeDx, "No #eta cut", "P");
-
- // TLegend* legScan = new TLegend(0.44, 0.67, 0.64, 0.87);
- // legScan->SetFillColor(0);
- // legScan->SetBorderSize(0);
- // legScan->SetTextSize(0.05);
- // legScan->AddEntry(hSigmaVsNcl, "No #eta cut", "P");
-
- // TH2D* hArray[4] = {hDeDxVsNcl0, hDeDxVsNcl1, hDeDxVsNcl2, hDeDxVsNcl3};
- // TH1D* hArrayDeDx[4] = {hDeDx0, hDeDx1, hDeDx3, hDeDx4};
-
- // Int_t color[4] = {2, kGreen+1, 4, 6};
- // Double_t x[4];
- // Double_t x_err[4];
- // Double_t slope[4];
- // Double_t slope_err[4];
- // Double_t offset[4];
- // Double_t offset_err[4];
- // for(Int_t i = 0; i < 4; i++) {
-
- // hArrayDeDx[i]->SetMarkerStyle(20);
- // hArrayDeDx[i]->SetMarkerColor(color[i]);
- // hArrayDeDx[i]->Scale(1.0/hArrayDeDx[i]->GetEntries());
- // cDeDx->cd();
- // hArrayDeDx[i]->DrawCopy("SAME");
-
- // x[i] = hMeanEta->GetBinContent(i+1);
- // x_err[i] = hMeanEta->GetBinError(i+1);
-
- // hArray[i]->FitSlicesY(0, 0, -1, 0, "QNR", helpArray);
- // hSigmaVsNclHelp = (TH1D*)helpArray->At(2);
- // hSigmaVsNclHelp->SetName(Form("hSigmaVsNcl%d", i));
- // hSigmaVsNclHelp->SetTitle("#sigma vs Ncl; Ncl; #sigma");
- // hSigmaVsNclHelp->SetMarkerStyle(20);
- // hSigmaVsNclHelp->SetMarkerColor(color[i]);
- // hSigmaVsNclHelp->Fit(fSDeDxVsNclHelp, "0");
- // slope[i] = fSDeDxVsNclHelp->GetParameter(0);
- // slope_err[i] = fSDeDxVsNclHelp->GetParError(0);
- // offset[i] = fSDeDxVsNclHelp->GetParameter(1);
- // offset_err[i] = fSDeDxVsNclHelp->GetParError(1);
- // cNcl3->cd();
- // hSigmaVsNclHelp->DrawCopy("SAME");
- // // fSDeDxVsNclHelp->DrawCopy("SAME");
-
- // if(i==0) {
- // legDeDx->AddEntry(hArrayDeDx[i], "|#eta|<0.2", "P");
- // legScan->AddEntry(hSigmaVsNclHelp, "|#eta|<0.2", "P");
- // } else if(i==1) {
- // legDeDx->AddEntry(hArrayDeDx[i], "0.2<|#eta|<0.4", "P");
- // legScan->AddEntry(hSigmaVsNclHelp, "0.2<|#eta|<0.4", "P");
- // } else if(i==2) {
- // legDeDx->AddEntry(hArrayDeDx[i], "0.4<|#eta|<0.6", "P");
- // legScan->AddEntry(hSigmaVsNclHelp, "0.4<|#eta|<0.6", "P");
- // } else if(i==3) {
- // legDeDx->AddEntry(hArrayDeDx[i], "0.6<|#eta|<0.8", "P");
- // legScan->AddEntry(hSigmaVsNclHelp, "0.6<|#eta|<0.8", "P");
- // }
-
- // cNcl3Mean->cd();
- // hMeanVsNclHelp = (TH1D*)helpArray->At(1);
- // hMeanVsNclHelp->SetName(Form("hMeanVsNcl%d", i));
- // hMeanVsNclHelp->SetTitle("Mean vs Ncl; Ncl; Mean");
- // hMeanVsNclHelp->SetMarkerStyle(20);
- // hMeanVsNclHelp->SetMarkerColor(color[i]);
- // hMeanVsNclHelp->DrawCopy("SAME");
- // }
-
- // cDeDx->cd();
- // hDeDx->DrawCopy("SAME");
- // legDeDx->Draw();
- // cDeDx->SaveAs(Form("%s/dedx_resolution_vs_eta.gif", baseName.Data()));
-
- // cNcl3->cd();
- // hSigmaVsNcl->DrawCopy("SAME");
- // legScan->Draw();
- // cNcl3->SaveAs(Form("%s/sigma_vs_ncl_vs_eta.gif", baseName.Data()));
-
- // cNcl3Mean->cd();
- // hMeanVsNcl->DrawCopy("SAME");
- // legScan->Draw();
- // cNcl3Mean->SaveAs(Form("%s/mean_vs_ncl_vs_eta.gif", baseName.Data()));
-
- // TGraphErrors* graphSlope = new TGraphErrors(4, x, slope, x_err, slope_err);
- // graphSlope->SetTitle("slope vs <#eta>; <#eta>; slope");
- // graphSlope->SetMarkerStyle(29);
- // TCanvas* cNcl4 = new TCanvas("cNcl4", "Slope vs <eta>", 600, 400);
- // cNcl4->Clear();
- // graphSlope->Draw("AP");
- // cNcl4->SaveAs(Form("%s/fit_slope_vs_eta.gif", baseName.Data()));
-
- // TGraphErrors* graphOffset = new TGraphErrors(4, x, offset, x_err, offset_err);
- // graphOffset->SetMarkerStyle(29);
- // graphOffset->SetTitle("offset vs <#eta>; <#eta>; offset");
- // TCanvas* cNcl5 = new TCanvas("cNcl5", "Offset vs <eta>", 600, 400);
- // cNcl5->Clear();
- // graphOffset->Draw("AP");
- // cNcl5->SaveAs(Form("%s/fit_offset_vs_eta.gif", baseName.Data()));
-
- // }
-
-//____________________________________________________________________________
-void CompareYields(const Char_t* dataFileName1,
- const Char_t* dataFileName2,
- Double_t ptStart, Double_t ptStop,
- const Char_t* endName1=0,
- const Char_t* endName2=0,
- Int_t run = 0,
- Int_t filter = 1,
- Bool_t usePhiCut = kTRUE)
-{
- gStyle->SetOptStat(0);
-
-
- TFile* dataFile1 = FindFileFresh(dataFileName1);
- if(!dataFile1)
- return;
- AliHighPtDeDxCalib* data1 = (AliHighPtDeDxCalib*)GetObject(dataFile1, filter, usePhiCut, run, "filter", endName1);
- data1->Print();
-
- // gSystem->Exec("mv debug/* olddebug/");
-
-
- TH2D* hDeltaPiVsPt1 = data1->GetHistDeDxVsP(0);
-
- AliHighPtDeDxCalib* data2 = data1;
- if(dataFileName2) {
- TFile* dataFile2 = FindFileFresh(dataFileName2);
- if(!dataFile2)
- return;
- data2 = (AliHighPtDeDxCalib*)GetObject(dataFile2, filter, usePhiCut, run, "filter", endName2);
- data2->Print();
- } else if (endName2) {
-
- data2 = (AliHighPtDeDxCalib*)GetObject(dataFile1, filter, usePhiCut, run, "filter", endName2);
- }
-
- TH2D* hDeltaPiVsPt2 = data2->GetHistDeDxVsP(0);
-
- // Root is a bit stupid with finidng bins so we have to add and subtract a
- // little to be sure we get the right bin as we typically put edges as
- // limits
- const Int_t binStart = hDeltaPiVsPt1->GetXaxis()->FindBin(ptStart+0.01);
- ptStart = hDeltaPiVsPt1->GetXaxis()->GetBinLowEdge(binStart);
- const Int_t binStop = hDeltaPiVsPt1->GetXaxis()->FindBin(ptStop-0.01);
- ptStop = hDeltaPiVsPt1->GetXaxis()->GetBinUpEdge(binStop);
- // const Int_t nBins = binStop - binStart + 1;
-
- cout << "Doing fits from pTlow = " << ptStart << " (bin: " << binStart
- << ") to pThigh = " << ptStop << " (bin: " << binStop << ")" << endl;
-
-
- //cross check
- TCanvas* cFits = new TCanvas("cFits", "Fit comparison to data", 1200, 800);
- cFits->Clear();
- cFits->Divide(7, 4);
-
- TLegend* legend = new TLegend(0.11, 0.6, 0.35, 0.85);
- legend->SetBorderSize(0);
- legend->SetFillColor(0);
-
- TCanvas* cSingleFit = new TCanvas("cSingleFit", "single fit");
-
- for(Int_t bin = binStart; bin <= binStop; bin++){
-
- cout << "Making projection for bin: " << bin << endl;
-
- const Int_t j = bin-binStart;
-
- TH1D* hDeltaPiVsPtProj1 =(TH1D*)hDeltaPiVsPt1->ProjectionY(Form("hDeltaPiVsPtProj1%d", bin), bin, bin);
- // hDeltaPiVsPtProj1->GetXaxis()->SetRangeUser(-20, 20);
- hDeltaPiVsPtProj1->SetTitle(Form("%.1f<p<%.1f GeV/c",
- hDeltaPiVsPt1->GetXaxis()->GetBinLowEdge(bin),
- hDeltaPiVsPt1->GetXaxis()->GetBinUpEdge(bin)));
-
- TH1D* hDeltaPiVsPtProj2 =(TH1D*)hDeltaPiVsPt2->ProjectionY(Form("hDeltaPiVsPtProj2%d", bin), bin, bin);
-
- hDeltaPiVsPtProj1->SetLineColor(4);
- hDeltaPiVsPtProj2->SetLineColor(2);
- hDeltaPiVsPtProj1->SetNormFactor();
- hDeltaPiVsPtProj2->SetNormFactor();
-
- cFits->cd();
- cFits->cd(j + 1);
- hDeltaPiVsPtProj1->DrawCopy();
- hDeltaPiVsPtProj2->DrawCopy("SAME");
-
- cSingleFit->cd();
- cSingleFit->Clear();
- hDeltaPiVsPtProj1->DrawCopy();
- hDeltaPiVsPtProj2->DrawCopy("SAME");
-
- CreateDir("results/debug");
- cSingleFit->SaveAs(Form("results/debug/ptspectrum_bin%d.gif", bin));
- }
-}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff