+++ /dev/null
-#if !defined(__CINT__) || defined(__MAKECINT__)
-#include <Riostream.h>
-#include <TLatex.h>
-#include <TImage.h>
-#include <TSystem.h>
-#include <TPaveText.h>
-#include <TH1F.h>
-#include <TF1.h>
-#include <TH1D.h>
-#include <TH2F.h>
-#include <TMath.h>
-#include <TNtuple.h>
-#include <TGraphErrors.h>
-#include <TList.h>
-#include <TLegend.h>
-#include <TLegendEntry.h>
-#include <TCanvas.h>
-#include <TFile.h>
-#include <TStyle.h>
-#include <Rtypes.h>
-#include "AliITSsadEdxFitter.h"
-#endif
-
-
-Double_t BetheBloch(Double_t *mom, Double_t *mass);
-void Logo();
-//
-
-//______________________________________________________________________
-void MakeRawITSsaSpectraMultiBin(Bool_t optMC=kTRUE, Int_t multibin=0){
-
-
- AliITSsadEdxFitter *ITSsa=new AliITSsadEdxFitter(optMC);
-
- TString filename, dirName, ps0, ps1, ps2;
- Char_t listname[50];
- if(optMC) dirName=(Form("../gridmultiplicitybins/LHC10d1_1.5sigma_7DCA_negmag/Spectra_MC_negmag_MultBin%d",multibin));
- else dirName=(Form("../gridmultiplicitybins/data_1.5sigma_7DCA_negmag/Spectra_data_negmag_MultBin%d",multibin));
- switch(multibin){
- case 0:
- sprintf(listname,"clistITSsaMult0to9999");
- break;
- case 1:
- sprintf(listname,"clistITSsaMult0to5");
- break;
- case 2:
- sprintf(listname,"clistITSsaMult6to9");
- break;
- case 3:
- sprintf(listname,"clistITSsaMult10to14");
- break;
- case 4:
- sprintf(listname,"clistITSsaMult15to22");
- break;
- case 5:
- sprintf(listname,"clistITSsaMult23to9999");
- break;
- }
- if(optMC){
- filename=Form("%s/AnalysisResults.root",dirName.Data());
- ps0 = Form("%s/outSIM.ps[",dirName.Data());
- ps1 = Form("%s/outSIM.ps",dirName.Data());
- ps2 = Form("%s/outSIM.ps]",dirName.Data());
- }
- else{
- filename=Form("%s/AnalysisResults.root",dirName.Data());
- ps0 = Form("%s/outDATA.ps[",dirName.Data());
- ps1 = Form("%s/outDATA.ps",dirName.Data());
- ps2 = Form("%s/outDATA.ps]",dirName.Data());
- }
- TString openfilename="./";
- openfilename+=filename;
- cout<<openfilename<<endl;
- TFile *fi=new TFile(openfilename.Data());
- if(!fi){
- cout<<"TFile loading failed"<<endl;
- return;
- }
- TDirectoryFile *di=(TDirectoryFile*) fi->Get("PWG2SpectraITSsa");
- if(!di){
- cout<<"TDirectory loading failed!"<<endl;
- return;
- }
- TList *li=(TList*)di->Get(listname);
- if(!li){
- cout<<"TList loading failed"<<endl;
- return;
- }
- cout<<"File loaded"<<endl;
-
- TCanvas *cdummy=new TCanvas("dummy","IntegralMethod",1000,800);
- cdummy->Print(ps0.Data());
-
- //binning
- const Int_t nbins = 22;
- Double_t xbins[nbins+1]={0.08,0.10,0.12,0.14,0.16,0.18,0.20,0.25,0.30,0.35,0.40,0.45,0.50,0.55,0.60,0.65,0.70,0.75,0.80,0.85,0.90,0.95,1.0};
-
- //histograms
- TH1F *fHistMCPosPi[nbins];
- TH1F *fHistMCPosK[nbins];
- TH1F *fHistMCPosP[nbins];
- TH1F *fHistMCNegPi[nbins];
- TH1F *fHistMCNegK[nbins];
- TH1F *fHistMCNegP[nbins];
- TH1F *fHistPosPi[nbins];
- TH1F *fHistPosK[nbins];
- TH1F *fHistPosP[nbins];
- TH1F *fHistNegPi[nbins];
- TH1F *fHistNegK[nbins];
- TH1F *fHistNegP[nbins];
- for(Int_t m=0;m<nbins;m++){
- fHistMCNegPi[m] = (TH1F*)li->FindObject(Form("fHistMCNegPi%d",m));
- fHistMCNegK[m] = (TH1F*)li->FindObject(Form("fHistMCNegK%d",m));
- fHistMCNegP[m] = (TH1F*)li->FindObject(Form("fHistMCNegP%d",m));
- fHistMCPosPi[m] = (TH1F*)li->FindObject(Form("fHistMCPosPi%d",m));
- fHistMCPosK[m] = (TH1F*)li->FindObject(Form("fHistMCPosK%d",m));
- fHistMCPosP[m] = (TH1F*)li->FindObject(Form("fHistMCPosP%d",m));
- fHistPosPi[m] = (TH1F*)li->FindObject(Form("fHistPosPi%d",m));
- fHistPosK[m] = (TH1F*)li->FindObject(Form("fHistPosK%d",m));
- fHistPosP[m] = (TH1F*)li->FindObject(Form("fHistPosP%d",m));
- fHistNegPi[m] = (TH1F*)li->FindObject(Form("fHistNegPi%d",m));
- fHistNegK[m] = (TH1F*)li->FindObject(Form("fHistNegK%d",m));
- fHistNegP[m] = (TH1F*)li->FindObject(Form("fHistNegP%d",m));
- }
-
- TH1F *fHistNEvents = (TH1F*)li->FindObject("fHistNEvents");
- TH2F *fHistDEDX = (TH2F*)li->FindObject("fHistDEDX");
- TH2F *fHistDEDXdouble = (TH2F*)li->FindObject("fHistDEDXdouble");
- TH1F *fHistCharge[4];
- for(Int_t j=0;j<4;j++) fHistCharge[j] = (TH1F*)li->FindObject((Form("fHistChargeLay%d",j)));
-
- TH1F *hEffPos[3];
- TH1F *hEffNeg[3];
- TH1F *hCorrFacPos[3];
- TH1F *hCorrFacNeg[3];
- TH1F *hEffMCPIDPos[3];
- TH1F *hEffMCPIDNeg[3];
- TH1F *hCorrFacMCPIDNeg[3];
- TH1F *hCorrFacMCPIDPos[3];
- TH1F *hSpectraPrimPosMC[3];
- TH1F *hSpectraPrimNegMC[3];
- TH1F *hSpectraPrimPosMCBefEvSel[3];
- TH1F *hSpectraPrimNegMCBefEvSel[3];
- TH1F *hSpectraPos[3];
- TH1F *hSpectraNeg[3];
- TH1F *hSpectraMCPIDPos[3];
- TH1F *hSpectraMCPIDNeg[3];
- TH1F* hMeanPos[3];
- TH1F* hMeanNeg[3];
- TH1F* hSigmaPos[3];
- TH1F* hSigmaNeg[3];
- TGraph *gres[6][22];
-
- for(Int_t i=0; i<3; i++){
- hSpectraPrimPosMC[i]=(TH1F*)li->FindObject(Form("fHistPrimMCpos%d",i));
- hSpectraPrimNegMC[i]=(TH1F*)li->FindObject(Form("fHistPrimMCneg%d",i));
- hSpectraPrimPosMCBefEvSel[i]=(TH1F*)li->FindObject(Form("fHistPrimMCposBefEvSel%d",i));
- hSpectraPrimNegMCBefEvSel[i]=(TH1F*)li->FindObject(Form("fHistPrimMCnegBefEvSel%d",i));
- hSpectraMCPIDPos[i]=new TH1F(Form("hSpectraMCPIDPos%d",i),Form("hSpectraMCPIDPos%d",i),nbins,xbins);
- hSpectraMCPIDNeg[i]=new TH1F(Form("hSpectraMCPIDNeg%d",i),Form("hSpectraMCPIDNeg%d",i),nbins,xbins);
- hSpectraPos[i]=new TH1F(Form("hSpectraPos%d",i),Form("hSpectraPos%d",i),nbins,xbins);
- hSpectraNeg[i]=new TH1F(Form("hSpectraNeg%d",i),Form("hSpectraNeg%d",i),nbins,xbins);
- hMeanPos[i]=new TH1F(Form("hMeanPos%d",i),Form("hMeanPos%d",i),nbins,xbins);
- hMeanNeg[i]=new TH1F(Form("hMeanNeg%d",i),Form("hMeanNeg%d",i),nbins,xbins);
- hSigmaPos[i]=new TH1F(Form("hSigmaPos%d",i),Form("hSigmaPos%d",i),nbins,xbins);
- hSigmaNeg[i]=new TH1F(Form("hSigmaNeg%d",i),Form("hSigmaNeg%d",i),nbins,xbins);
- }
-
- //division for DeltaPt (for MC spectra)
- for(Int_t ipt=0;ipt<3;ipt++){
- for(Int_t bin=1; bin <= hSpectraPrimPosMC[ipt]->GetNbinsX(); bin++){
- Float_t binSize=hSpectraPrimPosMC[ipt]->GetBinLowEdge(bin+1) - hSpectraPrimPosMC[ipt]->GetBinLowEdge(bin);
- hSpectraPrimPosMC[ipt]->SetBinContent(bin, hSpectraPrimPosMC[ipt]->GetBinContent(bin) / binSize);
- hSpectraPrimPosMC[ipt]->SetBinError(bin, 0);
- hSpectraPrimNegMC[ipt]->SetBinContent(bin, hSpectraPrimNegMC[ipt]->GetBinContent(bin) / binSize);
- hSpectraPrimNegMC[ipt]->SetBinError(bin, 0);
- if(hSpectraPrimPosMCBefEvSel[ipt]){
- hSpectraPrimPosMCBefEvSel[ipt]->SetBinContent(bin, hSpectraPrimPosMCBefEvSel[ipt]->GetBinContent(bin) / binSize);
- hSpectraPrimPosMCBefEvSel[ipt]->SetBinError(bin, 0);
- }
- if(hSpectraPrimNegMCBefEvSel[ipt]){
- hSpectraPrimNegMCBefEvSel[ipt]->SetBinContent(bin, hSpectraPrimNegMCBefEvSel[ipt]->GetBinContent(bin) / binSize);
- hSpectraPrimNegMCBefEvSel[ipt]->SetBinError(bin, 0);
- }
- }
- }
- cout<<"All plots loaded"<<endl;
-
- //open output file to store the dedx distribution and the spectra
- TString savename=filename.Data();
- savename.ReplaceAll("AnalysisResults","SpectraReco");
- TFile *fout=new TFile(savename.Data(),"recreate");
- fout->cd();
-
- gStyle->SetOptStat(0);
- gStyle->SetOptFit(111);
- gStyle->SetPalette(1);
-
- //propaganda plot
- Float_t pdgmass[5]={0.13957,0.493677,0.938272,1.875612762,0.00596}; //mass for pi, K, P, d (Gev/c^2)
- TF1 *funPos[5];
- TF1 *funNeg[5];
- for(Int_t m=0;m<5;m++){
- funPos[m] = new TF1(Form("funPos%d",m),BetheBloch,0.02,5,1);
- funPos[m]->SetParameter(0,pdgmass[m]);
- funPos[m]->SetLineWidth(2);
- funNeg[m] = new TF1(Form("funNeg%d",m),BetheBloch,-5,0.02,1);
- funNeg[m]->SetParameter(0,-pdgmass[m]);
- funNeg[m]->SetLineWidth(2);
- }
-
- TCanvas *cEvents=new TCanvas("cEvents","cEvents",500,400);
- cEvents->cd();
- fHistNEvents->Draw("text");
- fHistNEvents->Write();
-
- TCanvas *cDEDX=new TCanvas("cDEDX","DEDX",1000,800);
- cDEDX->cd();
- gPad->SetLogx();
- gPad->SetLogz();
- fHistDEDX->GetXaxis()->SetRangeUser(0.08,5);
- fHistDEDX->GetYaxis()->SetRangeUser(0.,700);
- fHistDEDX->GetXaxis()->SetTitle("momentum [GeV/c]");
- fHistDEDX->GetYaxis()->SetTitle("dE [keV/300#mum]");
- fHistDEDX->Draw("colz");
- fHistDEDX->Write();
- for(Int_t m=0;m<3;m++) funPos[m]->Draw("same");
- Logo();
-
- TCanvas *cDEDXdouble=new TCanvas("cDEDXdouble","DEDXdouble",1000,800);
- cDEDXdouble->cd();
- gPad->SetLogz();
- fHistDEDXdouble->GetXaxis()->SetRangeUser(-3,3);
- fHistDEDXdouble->GetXaxis()->SetTitle("momentum * sign [GeV/c]");
- fHistDEDXdouble->GetYaxis()->SetTitle("dE [keV/300#mum]");
- fHistDEDXdouble->Draw("colz");
- fHistDEDXdouble->Write();
- for(Int_t m=0;m<3;m++) {
- funPos[m]->Draw("same");
- funNeg[m]->Draw("same");
- }
- Logo();
-
- //calibration check histo
- TCanvas *cs=new TCanvas("cs","cs",1000,800);
- cs->Divide(2,2);
- for(Int_t j=0;j<4;j++){
- cs->cd(j+1);
- fHistCharge[j]->GetXaxis()->SetTitle("Charge [keV]");
- if(j==0) fHistCharge[j]->SetTitle("Drift inner layer");
- if(j==1) fHistCharge[j]->SetTitle("Drift outer layer");
- if(j==2) fHistCharge[j]->SetTitle("Strip inner layer");
- if(j==3) fHistCharge[j]->SetTitle("Strip inner layer");
- fHistCharge[j]->Draw();
- fHistCharge[j]->SetFillColor(7);
- fHistCharge[j]->Fit("gaus","QR","",70,100);
- fHistCharge[j]->Write();
- }
-
- //canvas MC spectra
- if(optMC){
- TCanvas *cspectraMC=new TCanvas("cspectraMC","cspectraMC",1000,800);
- cspectraMC->Divide(2,1);
- cspectraMC->cd(1);
- gPad->SetLogy();
- gPad->SetGridy();
- for(Int_t i=0;i<3;i++){
- if(i==0){
- hSpectraPrimPosMC[i]->Draw("");
- hSpectraPrimPosMC[i]->SetTitle("ITS EvMC truth positive");
- hSpectraPrimPosMC[i]->SetMinimum(100);
- hSpectraPrimPosMC[i]->GetYaxis()->SetTitle("d^{2}N/dp_{t}dy");
- hSpectraPrimPosMC[i]->GetXaxis()->SetTitle("p_{t} [GeV/c]");
- }
- hSpectraPrimPosMC[i]->Draw("same");
- hSpectraPrimPosMC[i]->SetLineColor(i+2);
- hSpectraPrimPosMC[i]->SetMarkerColor(i+2);
- hSpectraPrimPosMC[i]->SetMarkerStyle(21+i);
- hSpectraPrimPosMC[i]->Write();
- }
- cspectraMC->cd(2);
- gPad->SetLogy();
- gPad->SetGridy();
- for(Int_t i=0;i<3;i++){
- if(i==0){
- hSpectraPrimNegMC[i]->Draw("");
- hSpectraPrimNegMC[i]->SetTitle("ITS MC truth negative");
- hSpectraPrimNegMC[i]->SetMinimum(100);
- hSpectraPrimNegMC[i]->GetYaxis()->SetTitle("d^{2}N/dp_{t}dy");
- hSpectraPrimNegMC[i]->GetXaxis()->SetTitle("p_{t} [GeV/c]");
- }
- hSpectraPrimNegMC[i]->Draw("same");
- hSpectraPrimNegMC[i]->SetLineColor(i+2);
- hSpectraPrimNegMC[i]->SetMarkerColor(i+2);
- hSpectraPrimNegMC[i]->SetMarkerStyle(21+i);
- hSpectraPrimNegMC[i]->Write();
- }
- cspectraMC->Print(ps0.Data());
- }
-
- //dedx distribution to be fitted
- TCanvas *cgausPipos=new TCanvas("cgausPipos","PIONS pos",1000,800);
- cgausPipos->Divide(5,4,0.001,0.001);
- TCanvas *cgausKpos=new TCanvas("cgausKpos","KAONS pos",1000,800);
- cgausKpos->Divide(5,4,0.001,0.001);
- TCanvas *cgausPpos=new TCanvas("cgausPpos","PROTONS pos",1000,800);
- cgausPpos->Divide(5,4,0.001,0.001);
- TCanvas *cgausPineg=new TCanvas("cgausPineg","PIONS neg",1000,800);
- cgausPineg->Divide(5,4,0.001,0.001);
- TCanvas *cgausKneg=new TCanvas("cgausKneg","KAONS neg",1000,800);
- cgausKneg->Divide(5,4,0.001,0.001);
- TCanvas *cgausPneg=new TCanvas("cgausPneg","PROTONS neg",1000,800);
- cgausPneg->Divide(5,4,0.001,0.001);
-
- //binsize for one histo only because are all equal
- Float_t binsize = fHistPosPi[0]->GetBinWidth(1);
- Double_t fpar[5],efpar[5];
- for(Int_t i=0; i<nbins-2; i++){
- //------------------- positive particles
- //pions
- cgausPipos->cd(i+1);
- gPad->SetLogy();
- fHistPosPi[i]->Write();
- gres[0][i]=new TGraph();
- ITSsa->DoFit(fHistPosPi[i],i,211,gres[0][i]);
- ITSsa->FillHisto(hSpectraPos[0],i+1,binsize,211);
- if(optMC) ITSsa->FillHistoMC(hSpectraMCPIDPos[0],i+1,211,fHistMCPosPi[i]);
- if(ITSsa->IsGoodBin(i,211)){
- ITSsa->GetFitPar(fpar,efpar);
- hMeanPos[0]->SetBinContent(i+1,fpar[1]);
- hMeanPos[0]->SetBinError(i+1,efpar[1]);
- hSigmaPos[0]->SetBinContent(i+1,fpar[2]);
- hSigmaPos[0]->SetBinError(i+1,efpar[2]);
- }
- cgausPipos->Update();
- //kaons
- cgausKpos->cd(i+1);
- gPad->SetLogy();
- fHistPosK[i]->Write();
- gres[1][i]=new TGraph();
- ITSsa->DoFit(fHistPosK[i],i,321,gres[1][i]);
- ITSsa->FillHisto(hSpectraPos[1],i+1,binsize,321);
- if(optMC) ITSsa->FillHistoMC(hSpectraMCPIDPos[1],i+1,321,fHistMCPosK[i]);
- if(ITSsa->IsGoodBin(i,321)){
- ITSsa->GetFitPar(fpar,efpar);
- hMeanPos[1]->SetBinContent(i+1,fpar[1]);
- hMeanPos[1]->SetBinError(i+1,efpar[1]);
- hSigmaPos[1]->SetBinContent(i+1,fpar[2]);
- hSigmaPos[1]->SetBinError(i+1,efpar[2]);
- }
- cgausKpos->Update();
- //protons
- cgausPpos->cd(i+1);
- gPad->SetLogy();
- fHistPosP[i]->Write();
- fHistPosP[i]->SetFillColor(16);
- gres[2][i]=new TGraph();
- ITSsa->DoFitProton(fHistPosP[i],i,2212,gres[2][i]);
- ITSsa->FillHisto(hSpectraPos[2],i+1,binsize,2212);
- if(optMC) ITSsa->FillHistoMC(hSpectraMCPIDPos[2],i+1,2212,fHistMCPosP[i]);
- if(ITSsa->IsGoodBin(i,2212)){
- ITSsa->GetFitPar(fpar,efpar);
- hMeanPos[2]->SetBinContent(i+1,fpar[1]);
- hMeanPos[2]->SetBinError(i+1,efpar[1]);
- hSigmaPos[2]->SetBinContent(i+1,fpar[2]);
- hSigmaPos[2]->SetBinError(i+1,efpar[2]);
- }
- cgausPpos->Update();
-
- //------------------- negative particles
- //pions
- cgausPineg->cd(i+1);
- gPad->SetLogy();
- fHistNegPi[i]->Write();
- gres[3][i]=new TGraph();
- ITSsa->DoFit(fHistNegPi[i],i,211,gres[3][i]);
- ITSsa->FillHisto(hSpectraNeg[0],i+1,binsize,211);
- if(optMC) ITSsa->FillHistoMC(hSpectraMCPIDNeg[0],i+1,211,fHistMCNegPi[i]);
- if(ITSsa->IsGoodBin(i,211)){
- ITSsa->GetFitPar(fpar,efpar);
- hMeanNeg[0]->SetBinContent(i+1,fpar[1]);
- hMeanNeg[0]->SetBinError(i+1,efpar[1]);
- hSigmaNeg[0]->SetBinContent(i+1,fpar[2]);
- hSigmaNeg[0]->SetBinError(i+1,efpar[2]);
- }
- cgausPineg->Update();
- //kaons
- cgausKneg->cd(i+1);
- gPad->SetLogy();
- fHistNegK[i]->Write();
- gres[4][i]=new TGraph();
- ITSsa->DoFit(fHistNegK[i],i,321,gres[4][i]);
- ITSsa->FillHisto(hSpectraNeg[1],i+1,binsize,321);
- if(optMC) ITSsa->FillHistoMC(hSpectraMCPIDNeg[1],i+1,321,fHistMCNegK[i]);
- if(ITSsa->IsGoodBin(i,321)){
- ITSsa->GetFitPar(fpar,efpar);
- hMeanNeg[1]->SetBinContent(i+1,fpar[1]);
- hMeanNeg[1]->SetBinError(i+1,efpar[1]);
- hSigmaNeg[1]->SetBinContent(i+1,fpar[2]);
- hSigmaNeg[1]->SetBinError(i+1,efpar[2]);
- }
- cgausKneg->Update();
- //protons
- cgausPneg->cd(i+1);
- gPad->SetLogy();
- fHistNegP[i]->Write();
- gres[5][i]=new TGraph();
- ITSsa->DoFitProton(fHistNegP[i],i,2212,gres[5][i]);
- ITSsa->FillHisto(hSpectraNeg[2],i+1,binsize,2212);
- if(optMC) ITSsa->FillHistoMC(hSpectraMCPIDNeg[2],i+1,2212,fHistMCNegP[i]);
- if(ITSsa->IsGoodBin(i,2212)){
- ITSsa->GetFitPar(fpar,efpar);
- hMeanNeg[2]->SetBinContent(i+1,fpar[1]);
- hMeanNeg[2]->SetBinError(i+1,efpar[1]);
- hSigmaNeg[2]->SetBinContent(i+1,fpar[2]);
- hSigmaNeg[2]->SetBinError(i+1,efpar[2]);
- }
- cgausPneg->Update();
- }
-
- //save histograms in the ps file
- cgausPipos->Print(ps1.Data());
- cgausKpos->Print(ps1.Data());
- cgausPpos->Print(ps1.Data());
- cgausPineg->Print(ps1.Data());
- cgausKneg->Print(ps1.Data());
- cgausPneg->Print(ps1.Data());
-
- //spectra REC
- TCanvas *cspecREC=new TCanvas("cspecREC","SPECTRA rec",1000,800);
- cspecREC->Divide(2,1);
- cspecREC->cd(1);
- gPad->SetLogy();
- gPad->SetGridy();
- for(Int_t i=0;i<3;i++){
- if(i==0){
- hSpectraPos[i]->Draw("text");
- hSpectraPos[i]->SetTitle("ITSsa RAW positive");
- hSpectraPos[i]->SetMinimum(10);
- hSpectraPos[i]->GetYaxis()->SetTitle("d^{2}N/dp_{t}dy");
- hSpectraPos[i]->GetXaxis()->SetTitle("p_{t} [GeV/c]");
- }
- hSpectraPos[i]->Draw("esame");
- hSpectraPos[i]->SetLineColor(i+2);
- hSpectraPos[i]->SetMarkerColor(i+2);
- hSpectraPos[i]->SetMarkerStyle(21+i);
- hSpectraPos[i]->Write();
- }
- TLegend *leg=new TLegend(0.51,0.11,0.84,0.35,NULL,"brNDC");
- leg->SetFillColor(0);
- leg->SetBorderSize(0);
- TLegendEntry *entry0=leg->AddEntry(hSpectraPos[0],"pions","p");
- entry0->SetTextColor(2);
- TLegendEntry *entry2=leg->AddEntry(hSpectraPos[1],"kaons","p");
- entry2->SetTextColor(3);
- TLegendEntry *entry4=leg->AddEntry(hSpectraPos[2],"protons","p");
- entry4->SetTextColor(4);
- leg->Draw("same");
-
- cspecREC->cd(2);
- gPad->SetLogy();
- gPad->SetGridy();
- for(Int_t i=0;i<3;i++){
- if(i==0){
- hSpectraNeg[i]->Draw("e");
- hSpectraNeg[i]->SetTitle("ITSsa RAW negative");
- hSpectraNeg[i]->SetMinimum(10);
- hSpectraNeg[i]->GetYaxis()->SetTitle("d^{2}N/dp_{t}dy");
- hSpectraNeg[i]->GetXaxis()->SetTitle("p_{t} [GeV/c]");
- }
- hSpectraNeg[i]->Draw("esame");
- hSpectraNeg[i]->SetLineColor(i+2);
- hSpectraNeg[i]->SetMarkerColor(i+2);
- hSpectraNeg[i]->SetMarkerStyle(21+i);
- hSpectraNeg[i]->Write();
- }
-
-
- //FitParameters
- TCanvas *cFitPar=new TCanvas("cFitPar","FitParameters",1000,800);
- cFitPar->Divide(2,3);
- TLatex** tplus=new TLatex*[3];
- TLatex** tminus=new TLatex*[3];
- tplus[0]=new TLatex(0.7,0.75,"#pi^{+}");
- tplus[1]=new TLatex(0.7,0.75,"K^{+}");
- tplus[2]=new TLatex(0.7,0.75,"p");
- tminus[0]=new TLatex(0.7,0.65,"#pi^{-}");
- tminus[1]=new TLatex(0.7,0.65,"K^{-}");
- tminus[2]=new TLatex(0.7,0.65,"#bar{p}");
-
- for(Int_t ipart=0; ipart<3; ipart++){
- tplus[ipart]->SetNDC();
- tplus[ipart]->SetTextColor(1);
- tplus[ipart]->SetTextSize(0.05);
- tminus[ipart]->SetNDC();
- tminus[ipart]->SetTextColor(4);
- tminus[ipart]->SetTextSize(0.05);
- cFitPar->cd(1+2*ipart);
- gPad->SetGridy();
- hMeanPos[ipart]->SetMarkerStyle(20);
- hMeanPos[ipart]->Draw("e");
- hMeanPos[ipart]->GetYaxis()->SetRangeUser(-1,1);
- hMeanNeg[ipart]->SetMarkerStyle(24);
- hMeanNeg[ipart]->SetMarkerColor(4);
- hMeanNeg[ipart]->SetLineColor(4);
- hMeanNeg[ipart]->Draw("esame");
- tplus[ipart]->Draw();
- tminus[ipart]->Draw();
- cFitPar->cd(2+2*ipart);
- gPad->SetGridy();
- hSigmaPos[ipart]->SetMarkerStyle(20);
- hSigmaPos[ipart]->Draw("e");
- hSigmaPos[ipart]->GetYaxis()->SetRangeUser(0.1,0.3);
- hSigmaNeg[ipart]->SetMarkerStyle(24);
- hSigmaNeg[ipart]->SetMarkerColor(4);
- hSigmaNeg[ipart]->SetLineColor(4);
- hSigmaNeg[ipart]->Draw("esame");
- hSigmaPos[ipart]->Write();
- hSigmaNeg[ipart]->Write();
- tplus[ipart]->Draw();
- tminus[ipart]->Draw();
- }
-
- //plus/minus ratio plots
- TH1F* hRatioPions=new TH1F(*hSpectraPos[0]);
- hRatioPions->Divide(hSpectraNeg[0]);
- TH1F* hRatioKaons=new TH1F(*hSpectraPos[1]);
- hRatioKaons->Divide(hSpectraNeg[1]);
- TH1F* hRatioProtons=new TH1F(*hSpectraPos[2]);
- hRatioProtons->Divide(hSpectraNeg[2]);
-
- TCanvas *cratios=new TCanvas("cratios","Ratios +/-",1000,800);
- cratios->Divide(1,3);
- cratios->cd(1);
- hRatioPions->SetMinimum(0.7);
- hRatioPions->SetMaximum(1.3);
- hRatioPions->GetXaxis()->SetTitle("p_{t} [GeV/c]");
- hRatioPions->GetYaxis()->SetTitle("#pi^{+}/#pi^{-}");
- hRatioPions->Draw("e");
- cratios->cd(2);
- hRatioKaons->SetMinimum(0.7);
- hRatioKaons->SetMaximum(1.3);
- hRatioKaons->GetXaxis()->SetTitle("p_{t} [GeV/c]");
- hRatioKaons->GetYaxis()->SetTitle("K^{+}/K^{-}");
- hRatioKaons->Draw("e");
- cratios->cd(3);
- hRatioProtons->SetMinimum(0.7);
- hRatioProtons->SetMaximum(1.3);
- hRatioProtons->GetXaxis()->SetTitle("p_{t} [GeV/c]");
- hRatioProtons->GetYaxis()->SetTitle("p/#bar{p}");
- hRatioProtons->Draw("e");
-
- //efficiency and correction factor histograms
- if(optMC){
- for(Int_t i=0;i<3;i++){
- hEffPos[i] = (TH1F*)hSpectraPos[i]->Clone(Form("hEffPos%d",i));
- hEffPos[i]->SetTitle(Form("hEffPos%d",i));
- hEffPos[i]->Divide(hEffPos[i], hSpectraPrimPosMC[i], 1.0, 1.0, "B");//binomial errors
- hEffPos[i]->SetLineColor(i+2);
- hEffPos[i]->SetMarkerColor(i+2);
- hEffPos[i]->SetMarkerStyle(21+i);
- hEffPos[i]->Write();
-
- if(hSpectraPrimPosMCBefEvSel[i]){
- hCorrFacPos[i] = (TH1F*)hSpectraPos[i]->Clone(Form("hCorrFacPos%d",i));
- hCorrFacPos[i]->SetTitle(Form("hCorrFacPos%d",i));
- hCorrFacPos[i]->Divide(hCorrFacPos[i], hSpectraPrimPosMCBefEvSel[i], 1.0, 1.0, "B");//binomial errors
- hCorrFacPos[i]->SetLineColor(i+2);
- hCorrFacPos[i]->SetMarkerColor(i+2);
- hCorrFacPos[i]->SetMarkerStyle(25+i);
- hCorrFacPos[i]->Write();
- }
- hEffMCPIDPos[i] = (TH1F*)hSpectraMCPIDPos[i]->Clone(Form("hEffMCPIDPos%d",i));
- hEffMCPIDPos[i]->SetTitle(Form("hEffMCPIDPos%d",i));
- hEffMCPIDPos[i]->Divide(hEffMCPIDPos[i], hSpectraPrimPosMC[i], 1.0, 1.0, "B");//binomial errors
- hEffMCPIDPos[i]->SetLineColor(i+2);
- hEffMCPIDPos[i]->SetMarkerColor(i+2);
- hEffMCPIDPos[i]->SetMarkerStyle(25+i);
- hEffMCPIDPos[i]->Write();
-
- if(hSpectraPrimPosMCBefEvSel[i]){
- hCorrFacMCPIDPos[i] = (TH1F*)hSpectraMCPIDPos[i]->Clone(Form("hCorrFacMCPIDPos%d",i));
- hCorrFacMCPIDPos[i]->SetTitle(Form("hCorrFacMCPIDPos%d",i));
- hCorrFacMCPIDPos[i]->Divide(hCorrFacMCPIDPos[i], hSpectraPrimPosMCBefEvSel[i], 1.0, 1.0, "B");//binomial errors
- hCorrFacMCPIDPos[i]->SetLineColor(i+2);
- hCorrFacMCPIDPos[i]->SetMarkerColor(i+2);
- hCorrFacMCPIDPos[i]->SetMarkerStyle(25+i);
- hCorrFacMCPIDPos[i]->Write();
- }
-
- hEffNeg[i] = (TH1F*)hSpectraNeg[i]->Clone(Form("hEffNeg%d",i));
- hEffNeg[i]->SetTitle(Form("hEffNeg%d",i));
- hEffNeg[i]->Divide(hEffNeg[i], hSpectraPrimNegMC[i], 1.0, 1.0, "B");//binomial errors
- hEffNeg[i]->SetLineColor(i+2);
- hEffNeg[i]->SetMarkerColor(i+2);
- hEffNeg[i]->SetMarkerStyle(21+i);
- hEffNeg[i]->Write();
-
- if(hSpectraPrimNegMCBefEvSel[i]){
- hCorrFacNeg[i] = (TH1F*)hSpectraNeg[i]->Clone(Form("hCorrFacNeg%d",i));
- hCorrFacNeg[i]->SetTitle(Form("hCorrFacNeg%d",i));
- hCorrFacNeg[i]->Divide(hCorrFacNeg[i], hSpectraPrimNegMCBefEvSel[i], 1.0, 1.0, "B");//binomial errors
- hCorrFacNeg[i]->SetLineColor(i+2);
- hCorrFacNeg[i]->SetMarkerColor(i+2);
- hCorrFacNeg[i]->SetMarkerStyle(25+i);
- hCorrFacNeg[i]->Write();
- }
-
- hEffMCPIDNeg[i] = (TH1F*)hSpectraMCPIDNeg[i]->Clone(Form("hEffMCPIDNeg%d",i));
- hEffMCPIDNeg[i]->SetTitle(Form("hEffMCPIDNeg%d",i));
- hEffMCPIDNeg[i]->Divide(hEffMCPIDNeg[i], hSpectraPrimNegMC[i], 1.0, 1.0, "B");//binomial errors
- hEffMCPIDNeg[i]->SetLineColor(i+2);
- hEffMCPIDNeg[i]->SetMarkerColor(i+2);
- hEffMCPIDNeg[i]->SetMarkerStyle(25+i);
- hEffMCPIDNeg[i]->Write();
-
- if(hSpectraPrimPosMCBefEvSel[i]){
- hCorrFacMCPIDNeg[i] = (TH1F*)hSpectraMCPIDPos[i]->Clone(Form("hCorrFacMCPIDNeg%d",i));
- hCorrFacMCPIDNeg[i]->SetTitle(Form("hCorrFacMCPIDNeg%d",i));
- hCorrFacMCPIDNeg[i]->Divide(hCorrFacMCPIDNeg[i], hSpectraPrimNegMCBefEvSel[i], 1.0, 1.0, "B");//binomial errors
- hCorrFacMCPIDNeg[i]->SetLineColor(i+2);
- hCorrFacMCPIDNeg[i]->SetMarkerColor(i+2);
- hCorrFacMCPIDNeg[i]->SetMarkerStyle(25+i);
- hCorrFacMCPIDNeg[i]->Write();
- }
-
- hSpectraMCPIDPos[i]->SetLineColor(i+2);
- hSpectraMCPIDPos[i]->SetMarkerColor(i+2);
- hSpectraMCPIDPos[i]->SetMarkerStyle(25+i);
- hSpectraMCPIDPos[i]->Write();
-
- hSpectraMCPIDNeg[i]->SetLineColor(i+2);
- hSpectraMCPIDNeg[i]->SetMarkerColor(i+2);
- hSpectraMCPIDNeg[i]->SetMarkerStyle(25+i);
- hSpectraMCPIDNeg[i]->Write();
- }
-
- //spectra REC Ideal PID
- TCanvas *cspecREC2=new TCanvas("cspecIdPid","SPECTRA rec Ideal PID",1000,800);
- cspecREC2->Divide(2,1);
- cspecREC2->cd(1);
- gPad->SetLogy();
- gPad->SetGridy();
- for(Int_t i=0;i<3;i++){
- if(i==0){
- hSpectraMCPIDPos[i]->Draw("text");
- hSpectraMCPIDPos[i]->SetTitle("ITSsa RAW positive - Ideal PID");
- hSpectraMCPIDPos[i]->SetMinimum(10);
- hSpectraMCPIDPos[i]->GetYaxis()->SetTitle("d^{2}N/dp_{t}dy");
- hSpectraMCPIDPos[i]->GetXaxis()->SetTitle("p_{t} [GeV/c]");
- }
- hSpectraMCPIDPos[i]->Draw("esame");
- }
- leg->Draw("same");
- cspecREC2->cd(2);
- gPad->SetLogy();
- gPad->SetGridy();
- for(Int_t i=0;i<3;i++){
- if(i==0){
- hSpectraMCPIDNeg[i]->Draw("e");
- hSpectraMCPIDNeg[i]->SetTitle("ITSsa RAW negative -Ideal PID");
- hSpectraMCPIDNeg[i]->SetMinimum(10);
- hSpectraMCPIDNeg[i]->GetYaxis()->SetTitle("d^{2}N/dp_{t}dy");
- hSpectraMCPIDNeg[i]->GetXaxis()->SetTitle("p_{t} [GeV/c]");
- }
- hSpectraMCPIDNeg[i]->Draw("esame");
- }
-
- // Efficiency
- TCanvas *ceff=new TCanvas("ceff","ceff",1000,800);
- ceff->Divide(2,1);
- ceff->cd(1);
- gPad->SetGridy();
- for(Int_t i=0;i<3;i++){
- if(i==0){
- hEffPos[i]->SetTitle("ITSsa efficiency (positive)");
- hEffPos[i]->Draw("e");
- hEffPos[i]->GetYaxis()->SetRangeUser(0.1,0.9);
- hEffPos[i]->GetYaxis()->SetTitle("#epsilon");
- hEffPos[i]->GetXaxis()->SetTitle("p_{t} [GeV/c]");
- }
- hEffPos[i]->Draw("esame");
- hEffMCPIDPos[i]->Draw("esame");
- }
- TLegend *legeff=new TLegend(0.51,0.11,0.89,0.5,NULL,"brNDC");
- legeff->SetBorderSize(0);
- legeff->SetFillColor(0);
- TLegendEntry *entryeff0=legeff->AddEntry(hEffPos[0],"pions from fit","p");
- entryeff0->SetTextColor(2);
- TLegendEntry *entryeff2=legeff->AddEntry(hEffPos[1],"kaons from fit","p");
- entryeff2->SetTextColor(3);
- TLegendEntry *entryeff4=legeff->AddEntry(hEffPos[2],"protons from fit","p");
- entryeff4->SetTextColor(4);
- TLegendEntry *entryeff1=legeff->AddEntry(hEffMCPIDPos[0],"pions ideal PID","p");
- entryeff1->SetTextColor(2);
- TLegendEntry *entryeff3=legeff->AddEntry(hEffMCPIDPos[1],"kaons ideal PID","p");
- entryeff3->SetTextColor(3);
- TLegendEntry *entryeff5=legeff->AddEntry(hEffMCPIDPos[2],"protons ideal PID","p");
- entryeff5->SetTextColor(4);
- legeff->Draw("same");
-
- ceff->cd(2);
- gPad->SetGridy();
- for(Int_t i=0;i<3;i++){
- if(i==0){
- hEffNeg[i]->SetTitle("ITSsa efficiency (negative)");
- hEffNeg[i]->Draw("e");
- hEffNeg[i]->GetYaxis()->SetRangeUser(0.1,0.9);
- hEffNeg[i]->GetYaxis()->SetTitle("#epsilon");
- hEffNeg[i]->GetXaxis()->SetTitle("p_{t} [GeV/c]");
- }
- hEffNeg[i]->Draw("esame");
- hEffMCPIDNeg[i]->Draw("esame");
- }
- ceff->Update();
-
- // Correction factor
- TCanvas *ccf=new TCanvas("ccf","Corr Factor",1000,800);
- ccf->Divide(2,1);
- ccf->cd(1);
- gPad->SetGridy();
- for(Int_t i=0;i<3;i++){
- if(i==0){
- hEffPos[i]->SetTitle("ITSsa efficiency (positive)");
- hEffPos[i]->Draw("e");
- hEffPos[i]->GetYaxis()->SetRangeUser(0.1,0.9);
- hEffPos[i]->GetYaxis()->SetTitle("#epsilon");
- hEffPos[i]->GetXaxis()->SetTitle("p_{t} [GeV/c]");
- }
- hEffPos[i]->Draw("esame");
- if(hSpectraPrimPosMCBefEvSel[i]) hCorrFacPos[i]->Draw("esame");
- }
- TLegend *legcf=new TLegend(0.51,0.11,0.89,0.5,NULL,"brNDC");
- legcf->SetBorderSize(0);
- legcf->SetFillColor(0);
- TLegendEntry *entrycf0=legcf->AddEntry(hEffPos[0],"pions - eff.","p");
- entrycf0->SetTextColor(2);
- TLegendEntry *entrycf2=legcf->AddEntry(hEffPos[1],"kaons - eff.","p");
- entrycf2->SetTextColor(3);
- TLegendEntry *entrycf4=legcf->AddEntry(hEffPos[2],"protons - eff.","p");
- entrycf4->SetTextColor(4);
- TLegendEntry *entrycf1=legcf->AddEntry(hCorrFacPos[0],"pions - Corr. Fac.","p");
- entrycf1->SetTextColor(2);
- TLegendEntry *entrycf3=legcf->AddEntry(hCorrFacPos[1],"kaons - Corr. Fac.","p");
- entrycf3->SetTextColor(3);
- TLegendEntry *entrycf5=legcf->AddEntry(hCorrFacPos[2],"protons - Corr. Fac.","p");
- entrycf5->SetTextColor(4);
- legcf->Draw("same");
-
- ccf->cd(2);
- gPad->SetGridy();
- for(Int_t i=0;i<3;i++){
- if(i==0){
- hEffNeg[i]->SetTitle("ITSsa efficiency (negative)");
- hEffNeg[i]->Draw("e");
- hEffNeg[i]->GetYaxis()->SetRangeUser(0.1,0.9);
- hEffNeg[i]->GetYaxis()->SetTitle("#epsilon");
- hEffNeg[i]->GetXaxis()->SetTitle("p_{t} [GeV/c]");
- }
- hEffNeg[i]->Draw("esame");
- if(hSpectraPrimNegMCBefEvSel[i]) hCorrFacNeg[i]->Draw("esame");
- }
- ccf->Update();
-
- cspecREC2->Print(ps1.Data());
- ceff->Print(ps1.Data());
- ccf->Print(ps1.Data());
- }
-
- //save canvas in the ps file
- cratios->Print(ps1.Data());
- cspecREC->Print(ps1.Data());
- cFitPar->Print(ps1.Data());
- cs->Print(ps1.Data());
- cDEDX->Print(ps1.Data());
- cDEDXdouble->Print(ps1.Data());
- cdummy->Print(ps2.Data());
- // fout->Close();
- return;
-
-}//end of the main
-
-
-//______________________________________________________________________
-Double_t BetheBloch(Double_t *mom, Double_t *mass) {
- // BB PHobos parametrization fine tuned by G. Ortona
- // on data and MC in June 2010
- Double_t bg=mom[0]/mass[0];
- Double_t par[5];
- Bool_t MC=kFALSE;
- if(MC){
- par[0]=1.39126e+02;//tuned on LHC10d (sim pass4)
- par[1]=2.33589e+01;
- par[2]=6.05219e-02;
- par[3]=2.04336e-01;
- par[4]=-4.99854e-04;
- }
- else{
- par[0]=5.33458e+04;//tuned on data pass6
- par[1]=1.65303e+01;
- par[2]=2.60065e-03;
- par[3]=3.59533e-04;
- par[4]=7.51168e-05;
- }
- Double_t beta = bg/TMath::Sqrt(1.+ bg*bg);
- Double_t gamma=bg/beta;
- Double_t eff=1.0;
- if(bg<par[2]) eff=(bg-par[3])*(bg-par[3])+par[4];
- else eff=(par[2]-par[3])*(par[2]-par[3])+par[4];
- return (par[1]+2.0*TMath::Log(gamma)-beta*beta)*(par[0]/(beta*beta))*eff;
-}
-
-//______________________________________________________________________
-void Logo(){
- //method to plot the ALICE logo in the propaganda plot
- TPaveText *tpave=new TPaveText(0.3,0.7,0.7,0.89,"brNDC");
- tpave->SetBorderSize(0);
- tpave->SetFillStyle(0);
- tpave->SetFillColor(0);
- TText *txt1=tpave->AddText("ALICE performance");
- TText *txt2=tpave->AddText("ITS stand-alone tracks");
- TText *txt3=tpave->AddText("pp @ #sqrt{s} = 7 TeV");
- txt1->SetTextFont(62);
- txt2->SetTextFont(62);
- txt3->SetTextFont(62);
- tpave->Draw();
- TPad *pad1=new TPad("pad1", "pad1", 0.75, 0.7, 0.89, 0.89);
- pad1->Draw();
- pad1->cd();
- TImage *img1 = TImage::Open("ALICElogo.gif");
- img1->FillRectangle(0);
- img1->Draw();
-}
-
-//EOF
-
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff