-//------------------------------------------------------------------------------
-// Implementation of AliComparisonEff class. It keeps information from
-// comparison of reconstructed and MC particle tracks. In addtion,
-// it keeps selection cuts used during comparison. The comparison
-// information is stored in the ROOT histograms. Analysis of these
-// histograms can be done by using Analyse() class function. The result of
-// the analysis (histograms) are stored in the output picture_eff.root file.
-//
-// Author: J.Otwinowski 04/02/2008
-//------------------------------------------------------------------------------
-
-/*
- // after running analysis, read the file, and get component
- gSystem->Load("libPWG1.so");
- TFile f("Output.root");
- AliComparisonEff * comp = (AliComparisonEff*)f.Get("AliComparisonEff");
-
-
- // analyse comparison data (output stored in pictures_eff.root)
- comp->Analyse();
-
- // paramtetrisation of the TPC track length (for information only)
- TF1 fl("fl","((min(250./(abs(x+0.000001)),250)-90))",0,2); // length function
- TF1 fl2("fl2","[0]/((min(250./(abs(x+0.000001)),250)-90))^[1]",0,2);
- fl2.SetParameter(1,1);
- fl2.SetParameter(0,1);
-*/
-
-
-#include <iostream>
-
-#include "TFile.h"
-#include "TCint.h"
-#include "TH3F.h"
-#include "TH2F.h"
-#include "TF1.h"
-#include "TProfile.h"
-#include "TProfile2D.h"
-#include "TGraph2D.h"
-#include "TCanvas.h"
-#include "TGraph.h"
-//
-#include "AliESDEvent.h"
-#include "AliESD.h"
-#include "AliESDfriend.h"
-#include "AliESDfriendTrack.h"
-#include "AliRecInfoCuts.h"
-#include "AliMCInfoCuts.h"
-#include "AliLog.h"
-//
-#include "AliMathBase.h"
-#include "AliTreeDraw.h"
-#include "AliMagFMaps.h"
-#include "AliESDVertex.h"
-#include "AliExternalTrackParam.h"
-#include "AliTracker.h"
-
-#include "AliMCInfo.h"
-#include "AliESDRecInfo.h"
-#include "AliComparisonEff.h"
-
-using namespace std;
-
-
-ClassImp(AliComparisonEff)
-
-//_____________________________________________________________________________
-AliComparisonEff::AliComparisonEff():
- TNamed("AliComparisonEff","AliComparisonEff"),
-
- // histograms
-
- fMCPt(0),
- fMCRecPt(0),
- fMCRecPrimPt(0),
- fMCRecSecPt(0),
-
- fEffTPCPt(0),
- fEffTPCPtMC(0),
- fEffTPCPtF(0),
- //
- fEffTPCPt_P(0),
- fEffTPCPtMC_P(0),
- fEffTPCPtF_P(0),
- //
- fEffTPCPt_Pi(0),
- fEffTPCPtMC_Pi(0),
- fEffTPCPtF_Pi(0),
- //
- fEffTPCPt_K(0),
- fEffTPCPtMC_K(0),
- fEffTPCPtF_K(0),
-
- fEffTPCTan(0),
- fEffTPCTanMC(0),
- fEffTPCTanF(0),
- //
- fEffTPCPtTan(0),
- fEffTPCPtTanMC(0),
- fEffTPCPtTanF(0),
-
- // Cuts
- fCutsRC(0),
- fCutsMC(0),
-
- fVertex(0)
-{
- // init vertex
- fVertex = new AliESDVertex();
- fVertex->SetXv(0.0); fVertex->SetYv(0.0); fVertex->SetZv(0.0);
-
- for(Int_t i=0; i<4; ++i)
- {
- fTPCPtDCASigmaIdeal[i]=0;
- fTPCPtDCASigmaFull[i]=0;
- fTPCPtDCASigmaDay0[i]=0;
-
- fTPCPtDCAXY[i]=0;
- fTPCPtDCAZ[i]=0;
-
- fTPCPtDCASigmaIdealPid[i]=0;
- fTPCPtDCASigmaFullPid[i]=0;
- fTPCPtDCASigmaDay0Pid[i]=0;
-
- fTPCPtDCAXYPid[i]=0;
- fTPCPtDCAZPid[i]=0;
- }
- InitHisto();
- InitCuts();
-}
-
-//_____________________________________________________________________________
-AliComparisonEff::~AliComparisonEff(){
-
- //
- if(fMCPt) delete fEffTPCPt; fEffTPCPt=0;
- if(fMCRecPt) delete fMCRecPt; fMCRecPt=0;
- if(fMCRecPrimPt) delete fMCRecPrimPt; fMCRecPrimPt=0;
- if(fMCRecSecPt) delete fMCRecSecPt; fMCRecSecPt=0;
-
- //
- if(fEffTPCPt) delete fEffTPCPt; fEffTPCPt=0;
- if(fEffTPCPtMC) delete fEffTPCPtMC; fEffTPCPtMC=0;
- if(fEffTPCPtF) delete fEffTPCPtF; fEffTPCPtF=0;
-
- //
- if(fEffTPCPt_P) delete fEffTPCPt_P; fEffTPCPt_P=0;
- if(fEffTPCPtMC_P) delete fEffTPCPtMC_P; fEffTPCPtMC_P=0;
- if(fEffTPCPtF_P) delete fEffTPCPtF_P; fEffTPCPtF_P=0;
-
- //
- if(fEffTPCPt_Pi) delete fEffTPCPt_Pi; fEffTPCPt_Pi=0;
- if(fEffTPCPtMC_Pi) delete fEffTPCPtMC_Pi; fEffTPCPtMC_Pi=0;
- if(fEffTPCPtF_Pi) delete fEffTPCPtF_Pi; fEffTPCPtF_Pi=0;
-
- //
- if(fEffTPCPt_K) delete fEffTPCPt_K; fEffTPCPt_K=0;
- if(fEffTPCPtMC_K) delete fEffTPCPtMC_K; fEffTPCPtMC_K=0;
- if(fEffTPCPtF_K) delete fEffTPCPtF_K; fEffTPCPtF_K=0;
-
- //
- if(fEffTPCTan) delete fEffTPCTan; fEffTPCTan=0;
- if(fEffTPCTanMC) delete fEffTPCTanMC; fEffTPCTanMC=0;
- if(fEffTPCTanF) delete fEffTPCTanF; fEffTPCTanF=0;
-
- //
- if(fEffTPCPtTan) delete fEffTPCPtTan; fEffTPCPtTan=0;
- if(fEffTPCPtTanMC) delete fEffTPCPtTanMC; fEffTPCPtTanMC=0;
- if(fEffTPCPtTanF) delete fEffTPCPtTanF; fEffTPCPtTanF=0;
-
- for(Int_t i=0; i<4; ++i)
- {
- if(fTPCPtDCASigmaIdeal[i]) delete fTPCPtDCASigmaIdeal[i]; fTPCPtDCASigmaIdeal[i]=0;
- if(fTPCPtDCASigmaFull[i]) delete fTPCPtDCASigmaFull[i]; fTPCPtDCASigmaFull[i]=0;
- if(fTPCPtDCASigmaDay0[i]) delete fTPCPtDCASigmaDay0[i]; fTPCPtDCASigmaDay0[i]=0;
-
- if(fTPCPtDCAXY[i]) delete fTPCPtDCAXY[i]; fTPCPtDCAXY[i]=0;
- if(fTPCPtDCAZ[i]) delete fTPCPtDCAZ[i]; fTPCPtDCAZ[i]=0;
-
- if(fTPCPtDCASigmaIdealPid[i]) delete fTPCPtDCASigmaIdealPid[i]; fTPCPtDCASigmaIdealPid[i]=0;
- if(fTPCPtDCASigmaFullPid[i]) delete fTPCPtDCASigmaFullPid[i]; fTPCPtDCASigmaFullPid[i]=0;
- if(fTPCPtDCASigmaDay0Pid[i]) delete fTPCPtDCASigmaDay0Pid[i]; fTPCPtDCASigmaDay0Pid[i]=0;
-
- if(fTPCPtDCAXYPid[i]) delete fTPCPtDCAXYPid[i]; fTPCPtDCAXYPid[i]=0;
- if(fTPCPtDCAZPid[i]) delete fTPCPtDCAZPid[i]; fTPCPtDCAZPid[i]=0;
- }
-}
-
-//_____________________________________________________________________________
-void AliComparisonEff::InitHisto(){
-
- // Init histograms
- //
- fMCPt = new TH1F("fMCPt","fMCPt",50,0.1,3);
- fMCPt->SetXTitle("p_{t}");
- fMCPt->SetYTitle("yield");
-
- fMCRecPt = new TH1F("fMCRecPt","fMCRecPt",50,0.1,3);
- fMCRecPt->SetXTitle("p_{t}");
- fMCRecPt->SetYTitle("yield");
-
- fMCRecPrimPt = new TH1F("fMCRecPrimPt","fMCRecPrimPt",50,0.1,3);
- fMCRecPrimPt->SetXTitle("p_{t}");
- fMCRecPrimPt->SetYTitle("yield");
-
- fMCRecSecPt = new TH1F("fMCRecSecPt","fMCRecSecPt",50,0.1,3);
- fMCRecSecPt->SetXTitle("p_{t}");
- fMCRecSecPt->SetYTitle("yield");
-
- // Efficiency as function of pt
- fEffTPCPt = new TProfile("Eff_pt","Eff_Pt",50,0.1,3);
- fEffTPCPt->SetXTitle("p_{t}");
- fEffTPCPt->SetYTitle("TPC Efficiency");
-
- fEffTPCPtMC = new TProfile("MC_Eff_pt","MC_Eff_Pt",50,0.1,3);
- fEffTPCPtMC->SetXTitle("p_{t}");
- fEffTPCPtMC->SetYTitle("MC TPC Efficiency");
-
- fEffTPCPtF = new TProfile("F_Eff_pt","F_Eff_Pt",50,0.1,3);
- fEffTPCPtF->SetXTitle("p_{t}");
- fEffTPCPtF->SetYTitle("TPC Findable Efficiency");
-
- // Efficiency as function of pt protons
- fEffTPCPt_P = new TProfile("Eff_pt_P","Eff_Pt_P",50,0.1,3);
- fEffTPCPt_P->SetXTitle("p_{t}");
- fEffTPCPt_P->SetYTitle("TPC Efficiency");
-
- fEffTPCPtMC_P = new TProfile("MC_Eff_pt_P","MC_Eff_Pt_P",50,0.1,3);
- fEffTPCPtMC_P->SetXTitle("p_{t}");
- fEffTPCPtMC_P->SetYTitle("MC TPC Efficiency");
-
- fEffTPCPtF_P = new TProfile("F_Eff_pt_P","F_Eff_Pt_P",50,0.1,3);
- fEffTPCPtF_P->SetXTitle("p_{t}");
- fEffTPCPtF_P->SetYTitle("TPC Findable Efficiency");
-
- // Efficiency as function of pt pions
- fEffTPCPt_Pi = new TProfile("Eff_pt_Pi","Eff_Pit_Pi",50,0.1,3);
- fEffTPCPt_Pi->SetXTitle("p_{t}");
- fEffTPCPt_Pi->SetYTitle("TPC Efficiency");
-
- fEffTPCPtMC_Pi = new TProfile("MC_Eff_pt_Pi","MC_Eff_Pit_Pi",50,0.1,3);
- fEffTPCPtMC_Pi->SetXTitle("p_{t}");
- fEffTPCPtMC_Pi->SetYTitle("MC TPC Efficiency");
-
- fEffTPCPtF_Pi = new TProfile("F_Eff_pt_Pi","F_Eff_Pit_Pi",50,0.1,3);
- fEffTPCPtF_Pi->SetXTitle("p_{t}");
- fEffTPCPtF_Pi->SetYTitle("TPC Findable Efficiency");
-
- // Efficiency as function of pt kaons
- fEffTPCPt_K = new TProfile("Eff_pt_K","Eff_Kt_K",50,0.1,3);
- fEffTPCPt_K->SetXTitle("p_{t}");
- fEffTPCPt_K->SetYTitle("TPC Efficiency");
-
- fEffTPCPtMC_K = new TProfile("MC_Eff_pt_K","MC_Eff_Kt_K",50,0.1,3);
- fEffTPCPtMC_K->SetXTitle("p_{t}");
- fEffTPCPtMC_K->SetYTitle("MC TPC Efficiency");
-
- fEffTPCPtF_K = new TProfile("F_Eff_pt_K","F_Eff_Kt_K",50,0.1,3);
- fEffTPCPtF_K->SetXTitle("p_{t}");
- fEffTPCPtF_K->SetYTitle("TPC Findable Efficiency");
-
- // Efficiency as function of tan(theta)
- fEffTPCTan = new TProfile("Eff_tan","Eff_tan",50,-2.5,2.5);
- fEffTPCTan->SetXTitle("tan(#theta)");
- fEffTPCTan->SetYTitle("TPC Efficiency");
-
- fEffTPCTanMC = new TProfile("MC_Eff_tan","MC_Eff_tan",50,-2.5,2.5);
- fEffTPCTanMC->SetXTitle("tan(#theta)");
- fEffTPCTanMC->SetYTitle("MC TPC Efficiency");
-
- fEffTPCTanF = new TProfile("F_Eff_tan","F_Eff_tan",50,-2.5,2.5);
- fEffTPCPtF->SetXTitle("tan(#theta)");
- fEffTPCPtF->SetYTitle("TPC Findable Efficiency");
-
- // Efficiency as function of pt and tan(theta)
- fEffTPCPtTan = new TProfile2D("Eff_pt_tan","Eff_Pt_tan",10,0.1,3,20,-2.,2.);
- fEffTPCPtTan->SetXTitle("tan(#theta)");
- fEffTPCPtTan->SetYTitle("p_{t}");
-
- fEffTPCPtTanMC = new TProfile2D("MC_Eff_pt_tan_MC","MC Eff Pt",10,0.1,3,20,-2.,2.);
- fEffTPCPtTanMC->SetXTitle("tan(#theta)");
- fEffTPCPtTanMC->SetYTitle("p_{t}");
-
- fEffTPCPtTanF = new TProfile2D("MC_Eff_pt_tan_F","MC Eff Pt",10,0.1,3,20,-2.,2.);
- fEffTPCPtTanF->SetXTitle("tan(#theta)");
- fEffTPCPtTanF->SetYTitle("p_{t}");
-
- char name[256];
- for(Int_t i=0; i<4; ++i)
- {
- sprintf(name, "fTPCPtDCASigmaIdeal_%d",i);
- fTPCPtDCASigmaIdeal[i] = new TH2F(name,name,50,0.1,3,100,0,100);
-
- sprintf(name, "fTPCPtDCASigmaFull_%d",i);
- fTPCPtDCASigmaFull[i] = new TH2F(name,name,50,0.1,3,100,0,100);
-
- sprintf(name, "fTPCPtDCASigmaDay0_%d",i);
- fTPCPtDCASigmaDay0[i] = new TH2F(name,name,50,0.1,3,100,0,100);
-
- sprintf(name, "fTPCPtDCAXY_%d",i);
- fTPCPtDCAXY[i]= new TH2F(name,name,50,0.1,3,100,0,100);
- sprintf(name, "fTPCPtDCAZ_%d",i);
- fTPCPtDCAZ[i]= new TH2F(name,name,50,0.1,3,100,0,100);
-
- sprintf(name, "fTPCPtDCASigmaIdealPid_%d",i);
- fTPCPtDCASigmaIdealPid[i] = new TH3F(name,name,50,0.1,3,100,0,100,5,0,5);
-
- sprintf(name, "fTPCPtDCASigmaFullPid_%d",i);
- fTPCPtDCASigmaFullPid[i]= new TH3F(name,name,50,0.1,3,100,0,100,5,0,5);
-
- sprintf(name, "fTPCPtDCASigmaDay0Pid_%d",i);
- fTPCPtDCASigmaDay0Pid[i]= new TH3F(name,name,50,0.1,3,100,0,100,5,0,5);
-
- sprintf(name, "fTPCPtDCAXYPid_%d",i);
- fTPCPtDCAXYPid[i]= new TH3F(name,name,50,0.1,3,100,0,100,5,0,5);
-
- sprintf(name, "fTPCPtDCAZPid_%d",i);
- fTPCPtDCAZPid[i]= new TH3F(name,name,50,0.1,3,100,0,100,5,0,5);
- }
-}
-
-//_____________________________________________________________________________
-void AliComparisonEff::InitCuts()
-{
-
- if(!fCutsMC)
- AliDebug(AliLog::kError, "ERROR: Cannot find AliMCInfoCuts object");
- if(!fCutsRC)
- AliDebug(AliLog::kError, "ERROR: Cannot find AliRecInfoCuts object");
-}
-
-//_____________________________________________________________________________
-void AliComparisonEff::Process(AliMCInfo* infoMC, AliESDRecInfo *infoRC)
-{
- // Fill efficiency comparison information
-
- AliExternalTrackParam *track = 0;
- Double_t kRadius = 3.0; // beam pipe radius
- Double_t kMaxStep = 5.0; // max step
- Double_t field = 0.4; // mag. field
- Double_t kMaxD = 123456.0; // max distance
-
- // systematics
- const Double_t kSigma2Full_xy = 0.25; // ExB full systematics [cm]
- const Double_t kSigma2Full_z = 5.0; // [cm]
-
- const Double_t kSigma2Day0_xy = 0.02; // ExB [cm]
- const Double_t kSigma2Day0_z = 0.1; // [cm] goofie
-
- //
- Double_t DCASigmaIdeal=0;
- Double_t DCASigmaFull=0;
- Double_t DCASigmaDay0=0;
-
- Double_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
- const Double_t* dv;
-
- // distance to Prim. vertex
- dv = infoMC->GetVDist();
-
- Float_t mcpt = infoMC->GetParticle().Pt();
- Float_t tantheta = TMath::Tan(infoMC->GetParticle().Theta()-TMath::Pi()*0.5);
- //Bool_t isPrim = infoMC->GetParticle().R()<fCutsMC->GetMaxR() && TMath::Abs(infoMC->GetParticle().Vz())<fCutsMC->GetMaxVz();
- Bool_t isPrim = TMath::Sqrt(dv[0]*dv[0] + dv[1]*dv[1])<fCutsMC->GetMaxR() && TMath::Abs(dv[2])<fCutsMC->GetMaxVz();
-
- // calculate and set prim. vertex
- fVertex->SetXv( infoMC->GetParticle().Vx() - dv[0] );
- fVertex->SetYv( infoMC->GetParticle().Vy() - dv[1] );
- fVertex->SetZv( infoMC->GetParticle().Vz() - dv[2] );
-
- // Check selection cuts
- if (fCutsMC->IsPdgParticle(TMath::Abs(infoMC->GetParticle().GetPdgCode())) == kFALSE) return;
-
- // transform Pdg to Pid
- Double_t pid = -1;
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetEP() ) pid = 0;
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetMuP() ) pid = 1;
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetPiP() ) pid = 2;
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetKP() ) pid = 3;
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetProt() ) pid = 4;
-
- //cout << "dv[0] " << dv[0] << " dv[1] " << dv[1] << " dv[2] " << dv[2] << endl;
- //cout << "v[0] " << fVertex->GetXv() << " v[1] " << fVertex->GetYv() << " v[2] " << fVertex->GetZv()<< endl;
- if (TMath::Abs(tantheta)<fCutsRC->GetMaxAbsTanTheta())
- {
-
- // MC pt
- fMCPt->Fill(mcpt);
-
-
- if (infoRC->GetESDtrack()->GetTPCInnerParam())
- {
- if ((track = new AliExternalTrackParam(*infoRC->GetESDtrack()->GetTPCInnerParam())) != 0 )
- {
- AliTracker::PropagateTrackTo(track,kRadius,infoMC->GetMass(),kMaxStep,kTRUE);
- track->PropagateToDCA(fVertex,field,kMaxD,dca,cov);
-
- //
- cov[2] = track->GetCovariance()[2];
-
- // Eff = infoRC->GetStatus(1)==3 && isPrim / isPrim;
- // Pt vs ( dca[0]^2/cov[0]^2 + dca[1]^2/cov[2]^2 )
- // Pt vs ( dca[0]^2/(cov[0]^2 + kSigma2Full_xy) + dca[1]^2/( cov[2]^2 + kSigma2Full_z )
- // Pt vs ( dca[0]^2/(cov[0]^2 + kSigma2_xy) + dca[1]^2/( cov[2]^2 + kSigma2_z )
-
- if(cov[0]>0.0 && cov[2]>0.0)
- {
- DCASigmaIdeal = TMath::Power(dca[0],2)/cov[0]
- + TMath::Power(dca[1],2)/cov[2];
-
- DCASigmaFull = TMath::Power(dca[0],2)/(cov[0]+kSigma2Full_xy)
- + TMath::Power(dca[1],2)/(cov[2]+kSigma2Full_z);
-
- DCASigmaDay0 = TMath::Power(dca[0],2)/(cov[0]+kSigma2Day0_xy)
- + TMath::Power(dca[1],2)/(cov[2]+kSigma2Day0_z);
-
- //cout << "dca[0] " << dca[0] << " dca[1] " << dca[1] << endl;
- //cout << "cov[0] " << cov[0] << " cov[2] " << cov[2] << endl;
- //cout << "DCASigmaIdeal " << DCASigmaIdeal << " DCASigmaFull " << DCASigmaFull << " DCASigmaDay0 " <<DCASigmaDay0 << endl;
- //cout << " -------------------------------------------------------- "<< endl;
- }
-
- if(infoRC->GetStatus(1)==3) fMCRecPt->Fill(mcpt);
- if(infoRC->GetStatus(1)==3 && isPrim) fMCRecPrimPt->Fill(mcpt);
- if(infoRC->GetStatus(1)==3 && !isPrim) fMCRecSecPt->Fill(mcpt);
-
- if(isPrim)
- {
- fTPCPtDCASigmaIdeal[0]->Fill(mcpt,DCASigmaIdeal);
- fTPCPtDCASigmaFull[0]->Fill(mcpt,DCASigmaFull);
- fTPCPtDCASigmaDay0[0]->Fill(mcpt,DCASigmaDay0);
-
- fTPCPtDCAXY[0]->Fill(mcpt,dca[0]);
- fTPCPtDCAZ[0]->Fill(mcpt,dca[1]);
-
- fTPCPtDCASigmaIdealPid[0]->Fill(mcpt,DCASigmaIdeal,pid);
- fTPCPtDCASigmaFullPid[0]->Fill(mcpt,DCASigmaFull,pid);
- fTPCPtDCASigmaDay0Pid[0]->Fill(mcpt,DCASigmaDay0,pid);
-
- fTPCPtDCAXYPid[0]->Fill(mcpt,dca[0],pid);
- fTPCPtDCAZPid[0]->Fill(mcpt,dca[1],pid);
-
- if(infoRC->GetStatus(1)==3)
- {
- fTPCPtDCASigmaIdeal[1]->Fill(mcpt,DCASigmaIdeal);
- fTPCPtDCASigmaFull[1]->Fill(mcpt,DCASigmaFull);
- fTPCPtDCASigmaDay0[1]->Fill(mcpt,DCASigmaDay0);
-
- fTPCPtDCAXY[1]->Fill(mcpt,dca[0]);
- fTPCPtDCAZ[1]->Fill(mcpt,dca[1]);
-
- fTPCPtDCASigmaIdealPid[1]->Fill(mcpt,DCASigmaIdeal,pid);
- fTPCPtDCASigmaFullPid[1]->Fill(mcpt,DCASigmaFull,pid);
- fTPCPtDCASigmaDay0Pid[1]->Fill(mcpt,DCASigmaDay0,pid);
-
- fTPCPtDCAXYPid[1]->Fill(mcpt,dca[0],pid);
- fTPCPtDCAZPid[1]->Fill(mcpt,dca[1],pid);
- }
- }
-
- // Cont = infoRC->GetStatus(1)==3 && !isPrim / infoRC->GetStatus(1)==3
- // Pt vs ( dz[0]^2/cov[0]^2 + dz[1]^2/cov[2]^2 )
- // Pt vs ( dz[0]^2/(cov[0]^2 + kSigma2Full_xy) + dz[1]^2/( cov[2]^2 + kSigma2Full_z )
- // Pt vs ( dz[0]^2/(cov[0]^2 + kSigma2_xy) + dz[1]^2/( cov[2]^2 + kSigma2_z )
-
- if(infoRC->GetStatus(1)==3)
- {
- fTPCPtDCASigmaIdeal[2]->Fill(mcpt,DCASigmaIdeal);
- fTPCPtDCASigmaFull[2]->Fill(mcpt,DCASigmaFull);
- fTPCPtDCASigmaDay0[2]->Fill(mcpt,DCASigmaDay0);
-
- fTPCPtDCAXY[2]->Fill(mcpt,dca[0]);
- fTPCPtDCAZ[2]->Fill(mcpt,dca[1]);
-
- fTPCPtDCASigmaIdealPid[2]->Fill(mcpt,DCASigmaIdeal,pid);
- fTPCPtDCASigmaFullPid[2]->Fill(mcpt,DCASigmaFull,pid);
- fTPCPtDCASigmaDay0Pid[2]->Fill(mcpt,DCASigmaDay0,pid);
-
- fTPCPtDCAXYPid[2]->Fill(mcpt,dca[0],pid);
- fTPCPtDCAZPid[2]->Fill(mcpt,dca[1],pid);
-
- if(isPrim==0)
- {
- fTPCPtDCASigmaIdeal[3]->Fill(mcpt,DCASigmaIdeal);
- fTPCPtDCASigmaFull[3]->Fill(mcpt,DCASigmaFull);
- fTPCPtDCASigmaDay0[3]->Fill(mcpt,DCASigmaDay0);
-
- fTPCPtDCAXY[3]->Fill(mcpt,dca[0]);
- fTPCPtDCAZ[3]->Fill(mcpt,dca[1]);
-
- fTPCPtDCASigmaIdealPid[3]->Fill(mcpt,DCASigmaIdeal,pid);
- fTPCPtDCASigmaFullPid[3]->Fill(mcpt,DCASigmaFull,pid);
- fTPCPtDCASigmaDay0Pid[3]->Fill(mcpt,DCASigmaDay0,pid);
-
- fTPCPtDCAXYPid[3]->Fill(mcpt,dca[0],pid);
- fTPCPtDCAZPid[3]->Fill(mcpt,dca[1],pid);
- }
- }
- delete track;
- }
- }
- else
- {
- if(isPrim)
- {
- fTPCPtDCASigmaIdeal[0]->Fill(mcpt,0.0);
- fTPCPtDCASigmaFull[0]->Fill(mcpt,0.0);
- fTPCPtDCASigmaDay0[0]->Fill(mcpt,0.0);
-
- fTPCPtDCAXY[0]->Fill(mcpt,0.0);
- fTPCPtDCAZ[0]->Fill(mcpt,0.0);
-
- fTPCPtDCASigmaIdealPid[0]->Fill(mcpt,0.0,pid);
- fTPCPtDCASigmaFullPid[0]->Fill(mcpt,0.0,pid);
- fTPCPtDCASigmaDay0Pid[0]->Fill(mcpt,0.0,pid);
-
- fTPCPtDCAXYPid[0]->Fill(mcpt,0.0,pid);
- fTPCPtDCAZPid[0]->Fill(mcpt,0.0,pid);
- }
- }
- }
-
- // only primary particles
- if (!isPrim) return;
-
- // pt
- if (TMath::Abs(tantheta)<fCutsRC->GetMaxAbsTanTheta()){
-
- fEffTPCPt->Fill(mcpt, infoRC->GetStatus(1)==3);
- fEffTPCPtMC->Fill(mcpt, infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits());
- if (infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits()){
- fEffTPCPtF->Fill(mcpt, infoRC->GetStatus(1)==3);
- }
-
- // protons
- if(TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetProt()) {
- fEffTPCPt_P->Fill(mcpt, infoRC->GetStatus(1)==3);
- fEffTPCPtMC_P->Fill(mcpt, infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits());
-
- if(infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits()) {
- fEffTPCPtF_P->Fill(mcpt, infoRC->GetStatus(1)==3);
- }
- }
-
- // pions
- if(TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetPiP()) {
- fEffTPCPt_Pi->Fill(mcpt, infoRC->GetStatus(1)==3);
- fEffTPCPtMC_Pi->Fill(mcpt, infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits());
-
- if(infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits()) {
- fEffTPCPtF_Pi->Fill(mcpt, infoRC->GetStatus(1)==3);
- }
- }
-
- // kaons
- if(TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetKP()) {
- fEffTPCPt_K->Fill(mcpt, infoRC->GetStatus(1)==3);
- fEffTPCPtMC_K->Fill(mcpt, infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits());
-
- if(infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits()) {
- fEffTPCPtF_K->Fill(mcpt, infoRC->GetStatus(1)==3);
- }
- }
- }
-
- // theta
- if (TMath::Abs(mcpt)>fCutsRC->GetPtMin()){
- fEffTPCTan->Fill(tantheta, infoRC->GetStatus(1)==3);
- fEffTPCTanMC->Fill(tantheta, infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits());
- if (infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits()){
- fEffTPCTanF->Fill(tantheta, infoRC->GetStatus(1)==3);
- }
- }
-
- // pt-theta
- fEffTPCPtTan->Fill(mcpt,tantheta,infoRC->GetStatus(1)==3);
- fEffTPCPtTanMC->Fill(mcpt,tantheta,infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits());
- if (infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits()){
- fEffTPCPtTanF->Fill(mcpt,tantheta,infoRC->GetStatus(1)==3);
- }
-}
-
-//_____________________________________________________________________________
-void AliComparisonEff::Exec(AliMCInfo* infoMC, AliESDRecInfo *infoRC)
-{
- // Process comparison information
- Process(infoMC,infoRC);
-}
-
-//_____________________________________________________________________________
-Long64_t AliComparisonEff::Merge(TCollection* list)
-{
- // Merge list of objects (needed by PROOF)
-
- if (!list)
- return 0;
-
- if (list->IsEmpty())
- return 1;
-
- TIterator* iter = list->MakeIterator();
- TObject* obj = 0;
-
- // collection of generated histograms
-
- Int_t count=0;
- while((obj = iter->Next()) != 0)
- {
- AliComparisonEff* entry = dynamic_cast<AliComparisonEff*>(obj);
- if (entry == 0)
- continue;
-
- fMCPt->Add(entry->fMCPt);
- fMCRecPt->Add(entry->fMCRecPt);
- fMCRecPrimPt->Add(entry->fMCRecPrimPt);
- fMCRecSecPt->Add(entry->fMCRecSecPt);
-
- fEffTPCPt->Add(entry->fEffTPCPt);
- fEffTPCPtMC->Add(entry->fEffTPCPtMC);
- fEffTPCPtF->Add(entry->fEffTPCPtF);
-
- fEffTPCPt_P->Add(entry->fEffTPCPt_P);
- fEffTPCPtMC_P->Add(entry->fEffTPCPtMC_P);
- fEffTPCPtF_P->Add(entry->fEffTPCPtF_P);
-
- fEffTPCPt_Pi->Add(entry->fEffTPCPt_Pi);
- fEffTPCPtMC_Pi->Add(entry->fEffTPCPtMC_Pi);
- fEffTPCPtF_Pi->Add(entry->fEffTPCPtF_Pi);
-
- fEffTPCPt_K->Add(entry->fEffTPCPt_K);
- fEffTPCPtMC_K->Add(entry->fEffTPCPtMC_K);
- fEffTPCPtF_K->Add(entry->fEffTPCPtF_K);
-
- fEffTPCTan->Add(entry->fEffTPCTan);
- fEffTPCTanMC->Add(entry->fEffTPCTanMC);
- fEffTPCTanF->Add(entry->fEffTPCTanF);
-
- fEffTPCPtTan->Add(entry->fEffTPCPtTan);
- fEffTPCPtTanMC->Add(entry->fEffTPCPtTanMC);
- fEffTPCPtTanF->Add(entry->fEffTPCPtTanF);
-
- for(Int_t i=0; i<4; ++i)
- {
- fTPCPtDCASigmaIdeal[i]->Add(entry->fTPCPtDCASigmaIdeal[i]);
- fTPCPtDCASigmaFull[i]->Add(entry->fTPCPtDCASigmaFull[i]);
- fTPCPtDCASigmaDay0[i]->Add(entry->fTPCPtDCASigmaDay0[i]);
-
- fTPCPtDCAXY[i]->Add(entry->fTPCPtDCAXY[i]);
- fTPCPtDCAZ[i]->Add(entry->fTPCPtDCAXY[i]);
-
- fTPCPtDCASigmaIdealPid[i]->Add(entry->fTPCPtDCASigmaIdealPid[i]);
- fTPCPtDCASigmaFullPid[i]->Add(entry->fTPCPtDCASigmaFullPid[i]);
- fTPCPtDCASigmaDay0Pid[i]->Add(entry->fTPCPtDCASigmaDay0Pid[i]);
-
- fTPCPtDCAXYPid[i]->Add(entry->fTPCPtDCAXYPid[i]);
- fTPCPtDCAZPid[i]->Add(entry->fTPCPtDCAXYPid[i]);
- }
-
- count++;
- }
-
-return count;
-}
-
-//_____________________________________________________________________________
-void AliComparisonEff::Analyse()
-{
- // Analyse output histograms
-
- AliComparisonEff * comp=this;
-
- // calculate efficiency and contamination (4 sigma)
-
- TH1 *h_sigmaidealpid[20];
- TH1 *h_sigmafullpid[20];
- TH1 *h_sigmaday0pid[20];
-
- //TH1 *h_sigmaday0pidclone[20];
-
- TH1 *h_sigmaidealpidtot[4];
- TH1 *h_sigmafullpidtot[4];
- TH1 *h_sigmaday0pidtot[4];
-
- //TH1 *h_sigmaday0pidtotclone[4];
-
- char name[256];
- char name1[256];
- Int_t idx;
-
- for(Int_t i=0; i<4; ++i)
- {
- //total
- comp->fTPCPtDCASigmaIdealPid[i]->GetYaxis()->SetRange(1,4);
- sprintf(name,"h_sigmaidealpidtot_%d",i);
- h_sigmaidealpidtot[i] = comp->fTPCPtDCASigmaIdealPid[i]->Project3D();
- h_sigmaidealpidtot[i]->SetName(name);
-
- comp->fTPCPtDCASigmaFullPid[i]->GetYaxis()->SetRange(1,4);
- sprintf(name,"h_sigmafullpidtot_%d",i);
- h_sigmafullpidtot[i] = comp->fTPCPtDCASigmaFullPid[i]->Project3D();
- h_sigmafullpidtot[i]->SetName(name);
-
- comp->fTPCPtDCASigmaDay0Pid[i]->GetYaxis()->SetRange(1,4);
- sprintf(name,"h_sigmaday0pidtot_%d",i);
- h_sigmaday0pidtot[i] = comp->fTPCPtDCASigmaDay0Pid[i]->Project3D();
- h_sigmaday0pidtot[i]->SetName(name);
-
- // pid wise
- for(Int_t j=0; j<5; ++j)
- {
- idx = i*5 + j;
-
- comp->fTPCPtDCASigmaIdealPid[i]->GetYaxis()->SetRange(1,4);
- comp->fTPCPtDCASigmaIdealPid[i]->GetZaxis()->SetRange(j+1,j+1);
-
- sprintf(name,"h_sigmaidealpid_%d",idx);
- h_sigmaidealpid[idx] = comp->fTPCPtDCASigmaIdealPid[i]->Project3D();
- h_sigmaidealpid[idx]->SetName(name);
-
-
- comp->fTPCPtDCASigmaFullPid[i]->GetYaxis()->SetRange(1,4);
- comp->fTPCPtDCASigmaFullPid[i]->GetZaxis()->SetRange(j+1,j+1);
-
- sprintf(name,"h_sigmafullpid_%d",idx);
- h_sigmafullpid[idx] = comp->fTPCPtDCASigmaFullPid[i]->Project3D();
- h_sigmafullpid[idx]->SetName(name);
-
- comp->fTPCPtDCASigmaDay0Pid[i]->GetYaxis()->SetRange(1,4);
- comp->fTPCPtDCASigmaDay0Pid[i]->GetZaxis()->SetRange(j+1,j+1);
-
- sprintf(name,"h_sigmaday0pid_%d",idx);
- h_sigmaday0pid[idx] = comp->fTPCPtDCASigmaDay0Pid[i]->Project3D();
- h_sigmaday0pid[idx]->SetName(name);
-
- }
- }
-
- // calculate efficiency and contamination (all pids)
- h_sigmaidealpidtot[0]->Sumw2();
- h_sigmaidealpidtot[1]->Divide(h_sigmaidealpidtot[0]);
- h_sigmaidealpidtot[2]->Sumw2();
- h_sigmaidealpidtot[3]->Divide(h_sigmaidealpidtot[2]);
-
- h_sigmafullpidtot[0]->Sumw2();
- h_sigmafullpidtot[1]->Divide(h_sigmafullpidtot[0]);
- h_sigmafullpidtot[2]->Sumw2();
- h_sigmafullpidtot[3]->Divide(h_sigmafullpidtot[2]);
-
- h_sigmaday0pidtot[0]->Sumw2();
- h_sigmaday0pidtot[1]->Divide(h_sigmaday0pidtot[0]);
- h_sigmaday0pidtot[2]->Sumw2();
- h_sigmaday0pidtot[3]->Divide(h_sigmaday0pidtot[2]);
-
- // calculate efficiency pid wise
- for(Int_t idx = 0; idx<5; idx++)
- {
- h_sigmaidealpid[idx]->Sumw2();
- h_sigmaidealpid[idx+5]->Divide(h_sigmaidealpid[idx]);
-
- h_sigmafullpid[idx]->Sumw2();
- h_sigmafullpid[idx+5]->Divide(h_sigmafullpid[idx]);
-
- h_sigmaday0pid[idx]->Sumw2();
- h_sigmaday0pid[idx+5]->Divide(h_sigmaday0pid[idx]);
- }
-
- // calculate cont. pid wise
- for(Int_t idx = 0; idx<5; idx++)
- {
- h_sigmaidealpid[idx+15]->Divide(h_sigmaidealpidtot[2]);
- h_sigmafullpid[idx+15]->Divide(h_sigmafullpidtot[2]);
- h_sigmaday0pid[idx+15]->Divide(h_sigmaday0pidtot[2]);
- }
-
- // write results
- TFile *fp = new TFile("pictures_eff.root","recreate");
- fp->cd();
-
- TCanvas * c = new TCanvas("Efficiency","Track efficiency");
- c->cd();
-
- fMCPt->Write();
- fMCRecPt->Write();
- fMCRecPrimPt->Write();
- fMCRecSecPt->Write();
-
- for(int i = 0; i<4;i++)
- {
- comp->fTPCPtDCASigmaIdealPid[i]->GetYaxis()->SetRange();
- comp->fTPCPtDCASigmaIdealPid[i]->GetZaxis()->SetRange();
- comp->fTPCPtDCASigmaFullPid[i]->GetYaxis()->SetRange();
- comp->fTPCPtDCASigmaFullPid[i]->GetZaxis()->SetRange();
- comp->fTPCPtDCASigmaDay0Pid[i]->GetYaxis()->SetRange();
- comp->fTPCPtDCASigmaDay0Pid[i]->GetZaxis()->SetRange();
-
- comp->fTPCPtDCASigmaIdealPid[i]->Write();
- comp->fTPCPtDCASigmaFullPid[i]->Write();
- comp->fTPCPtDCASigmaDay0Pid[i]->Write();
- }
- //
- comp->fEffTPCTanF->SetXTitle("Tan(#theta)");
- comp->fEffTPCTanF->SetYTitle("eff_{findable}");
- comp->fEffTPCTanF->Write("EffTanFindable");
- //
- comp->fEffTPCTan->SetXTitle("Tan(#theta)");
- comp->fEffTPCTan->SetYTitle("eff_{all}");
- comp->fEffTPCTan->Write("EffTanAll");
-
- h_sigmaidealpidtot[1]->Write("Eff_SigmaIdeal");
- h_sigmaidealpidtot[3]->Write("Cont_SigmaIdeal");
-
- h_sigmafullpidtot[1]->Write("Eff_SigmaFull");
- h_sigmafullpidtot[3]->Write("Cont_SigmaFull");
-
- h_sigmaday0pidtot[1]->Write("Eff_SigmaDay0");
- h_sigmaday0pidtot[3]->Write("Cont_SigmaDay0");
-
- for(Int_t idx = 0; idx<5; idx++)
- {
- sprintf(name,"Eff_SigmaIdeal_%d",idx);
- sprintf(name1,"Cont_SigmaIdeal_%d",idx);
-
- h_sigmaidealpid[idx+5]->Write(name);
- h_sigmaidealpid[idx+15]->Write(name1);
-
- sprintf(name,"Eff_SigmaFull_%d",idx);
- sprintf(name1,"Cont_SigmaFull_%d",idx);
-
- h_sigmafullpid[idx+5]->Write(name);
- h_sigmafullpid[idx+15]->Write(name1);
-
- sprintf(name,"Eff_SigmaDay0_%d",idx);
- sprintf(name1,"Cont_SigmaDay0_%d",idx);
-
- h_sigmaday0pid[idx+5]->Write(name);
- h_sigmaday0pid[idx+15]->Write(name1);
- }
-
- //
- fp->Close();
-}
+//------------------------------------------------------------------------------\r
+// Implementation of AliComparisonEff class. It keeps information from \r
+// comparison of reconstructed and MC particle tracks. In addtion, \r
+// it keeps selection cuts used during comparison. The comparison \r
+// information is stored in the ROOT histograms. Analysis of these \r
+// histograms can be done by using Analyse() class function. The result of \r
+// the analysis (histograms/graphs) are stored in the folder which is \r
+// a data member of AliComparisonEff.\r
+// \r
+// Author: J.Otwinowski 04/02/2008 \r
+//------------------------------------------------------------------------------\r
+\r
+/*\r
+ \r
+ // after running comparison task, read the file, and get component\r
+ gROOT->LoadMacro("$ALICE_ROOT/PWG1/Macros/LoadMyLibs.C");\r
+ LoadMyLibs();\r
+ TFile f("Output.root");\r
+ AliComparisonEff * compObj = (AliComparisonEff*)f.Get("AliComparisonEff");\r
+\r
+ // Analyse comparison data\r
+ compObj->Analyse();\r
+\r
+ // the output histograms/graphs will be stored in the folder "folderEff" \r
+ compObj->GetAnalysisFolder()->ls("*");\r
+\r
+ // user can save whole comparison object (or only folder with anlysed histograms) \r
+ // in the seperate output file (e.g.)\r
+ TFile fout("Analysed_Eff.root","recreate");\r
+ compObj->Write(); // compObj->GetAnalysisFolder()->Write();\r
+ fout.Close();\r
+\r
+*/\r
+\r
+\r
+#include <iostream>\r
+\r
+#include "TFile.h"\r
+#include "TCint.h"\r
+#include "TH3F.h"\r
+#include "TH2F.h"\r
+#include "TF1.h"\r
+#include "TProfile.h"\r
+#include "TProfile2D.h"\r
+#include "TGraph2D.h"\r
+#include "TCanvas.h"\r
+#include "TGraph.h"\r
+// \r
+#include "AliESDEvent.h" \r
+#include "AliESD.h"\r
+#include "AliESDfriend.h"\r
+#include "AliESDfriendTrack.h"\r
+#include "AliRecInfoCuts.h" \r
+#include "AliMCInfoCuts.h" \r
+#include "AliLog.h" \r
+//\r
+#include "AliMathBase.h"\r
+#include "AliTreeDraw.h" \r
+#include "AliMagFMaps.h" \r
+#include "AliESDVertex.h" \r
+#include "AliExternalTrackParam.h" \r
+#include "AliTracker.h" \r
+\r
+#include "AliMCInfo.h" \r
+#include "AliESDRecInfo.h" \r
+#include "AliComparisonEff.h" \r
+\r
+using namespace std;\r
+\r
+\r
+ClassImp(AliComparisonEff)\r
+\r
+//_____________________________________________________________________________\r
+AliComparisonEff::AliComparisonEff():\r
+ AliComparisonObject("AliComparisonEff"),\r
+\r
+ // histograms\r
+ \r
+ fMCPt(0),\r
+ fMCRecPt(0),\r
+ fMCRecPrimPt(0),\r
+ fMCRecSecPt(0),\r
+\r
+ fEffTPCPt(0),\r
+ fEffTPCPtMC(0),\r
+ fEffTPCPtF(0),\r
+ //\r
+ fEffTPCPt_P(0),\r
+ fEffTPCPtMC_P(0),\r
+ fEffTPCPtF_P(0),\r
+ //\r
+ fEffTPCPt_Pi(0),\r
+ fEffTPCPtMC_Pi(0),\r
+ fEffTPCPtF_Pi(0),\r
+ //\r
+ fEffTPCPt_K(0),\r
+ fEffTPCPtMC_K(0),\r
+ fEffTPCPtF_K(0),\r
+ \r
+ fEffTPCTan(0),\r
+ fEffTPCTanMC(0),\r
+ fEffTPCTanF(0),\r
+ //\r
+ fEffTPCPtTan(0),\r
+ fEffTPCPtTanMC(0),\r
+ fEffTPCPtTanF(0),\r
+\r
+ // Cuts \r
+ fCutsRC(0), \r
+ fCutsMC(0),\r
+\r
+ fVertex(0),\r
+\r
+ // histogram folder \r
+ fAnalysisFolder(0)\r
+{\r
+ // init vertex\r
+ fVertex = new AliESDVertex();\r
+ fVertex->SetXv(0.0); fVertex->SetYv(0.0); fVertex->SetZv(0.0); \r
+\r
+ for(Int_t i=0; i<4; ++i)\r
+ {\r
+ fTPCPtDCASigmaIdeal[i]=0;\r
+ fTPCPtDCASigmaFull[i]=0;\r
+ fTPCPtDCASigmaDay0[i]=0;\r
+\r
+ fTPCPtDCAXY[i]=0;\r
+ fTPCPtDCAZ[i]=0;\r
+\r
+ fTPCPtDCASigmaIdealPid[i]=0;\r
+ fTPCPtDCASigmaFullPid[i]=0;\r
+ fTPCPtDCASigmaDay0Pid[i]=0;\r
+\r
+ fTPCPtDCAXYPid[i]=0; \r
+ fTPCPtDCAZPid[i]=0; \r
+ }\r
+\r
+ Init();\r
+}\r
+\r
+//_____________________________________________________________________________\r
+AliComparisonEff::~AliComparisonEff(){\r
+\r
+ // \r
+ if(fMCPt) delete fMCPt; fMCPt=0;\r
+ if(fMCRecPt) delete fMCRecPt; fMCRecPt=0;\r
+ if(fMCRecPrimPt) delete fMCRecPrimPt; fMCRecPrimPt=0;\r
+ if(fMCRecSecPt) delete fMCRecSecPt; fMCRecSecPt=0;\r
+\r
+ // \r
+ if(fEffTPCPt) delete fEffTPCPt; fEffTPCPt=0;\r
+ if(fEffTPCPtMC) delete fEffTPCPtMC; fEffTPCPtMC=0;\r
+ if(fEffTPCPtF) delete fEffTPCPtF; fEffTPCPtF=0;\r
+\r
+ // \r
+ if(fEffTPCPt_P) delete fEffTPCPt_P; fEffTPCPt_P=0;\r
+ if(fEffTPCPtMC_P) delete fEffTPCPtMC_P; fEffTPCPtMC_P=0;\r
+ if(fEffTPCPtF_P) delete fEffTPCPtF_P; fEffTPCPtF_P=0;\r
+\r
+ // \r
+ if(fEffTPCPt_Pi) delete fEffTPCPt_Pi; fEffTPCPt_Pi=0;\r
+ if(fEffTPCPtMC_Pi) delete fEffTPCPtMC_Pi; fEffTPCPtMC_Pi=0;\r
+ if(fEffTPCPtF_Pi) delete fEffTPCPtF_Pi; fEffTPCPtF_Pi=0;\r
+\r
+ // \r
+ if(fEffTPCPt_K) delete fEffTPCPt_K; fEffTPCPt_K=0;\r
+ if(fEffTPCPtMC_K) delete fEffTPCPtMC_K; fEffTPCPtMC_K=0;\r
+ if(fEffTPCPtF_K) delete fEffTPCPtF_K; fEffTPCPtF_K=0;\r
+\r
+ // \r
+ if(fEffTPCTan) delete fEffTPCTan; fEffTPCTan=0;\r
+ if(fEffTPCTanMC) delete fEffTPCTanMC; fEffTPCTanMC=0;\r
+ if(fEffTPCTanF) delete fEffTPCTanF; fEffTPCTanF=0;\r
+\r
+ //\r
+ if(fEffTPCPtTan) delete fEffTPCPtTan; fEffTPCPtTan=0;\r
+ if(fEffTPCPtTanMC) delete fEffTPCPtTanMC; fEffTPCPtTanMC=0;\r
+ if(fEffTPCPtTanF) delete fEffTPCPtTanF; fEffTPCPtTanF=0;\r
+\r
+ for(Int_t i=0; i<4; ++i)\r
+ {\r
+ if(fTPCPtDCASigmaIdeal[i]) delete fTPCPtDCASigmaIdeal[i]; fTPCPtDCASigmaIdeal[i]=0;\r
+ if(fTPCPtDCASigmaFull[i]) delete fTPCPtDCASigmaFull[i]; fTPCPtDCASigmaFull[i]=0;\r
+ if(fTPCPtDCASigmaDay0[i]) delete fTPCPtDCASigmaDay0[i]; fTPCPtDCASigmaDay0[i]=0;\r
+\r
+ if(fTPCPtDCAXY[i]) delete fTPCPtDCAXY[i]; fTPCPtDCAXY[i]=0;\r
+ if(fTPCPtDCAZ[i]) delete fTPCPtDCAZ[i]; fTPCPtDCAZ[i]=0;\r
+\r
+ if(fTPCPtDCASigmaIdealPid[i]) delete fTPCPtDCASigmaIdealPid[i]; fTPCPtDCASigmaIdealPid[i]=0;\r
+ if(fTPCPtDCASigmaFullPid[i]) delete fTPCPtDCASigmaFullPid[i]; fTPCPtDCASigmaFullPid[i]=0;\r
+ if(fTPCPtDCASigmaDay0Pid[i]) delete fTPCPtDCASigmaDay0Pid[i]; fTPCPtDCASigmaDay0Pid[i]=0;\r
+\r
+ if(fTPCPtDCAXYPid[i]) delete fTPCPtDCAXYPid[i]; fTPCPtDCAXYPid[i]=0;\r
+ if(fTPCPtDCAZPid[i]) delete fTPCPtDCAZPid[i]; fTPCPtDCAZPid[i]=0;\r
+ }\r
+\r
+ if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0;\r
+}\r
+\r
+//_____________________________________________________________________________\r
+void AliComparisonEff::Init(){\r
+\r
+ // Init histograms\r
+ //\r
+ fMCPt = new TH1F("fMCPt","fMCPt",50,0.1,3); \r
+ fMCPt->SetXTitle("p_{t}");\r
+ fMCPt->SetYTitle("yield");\r
+\r
+ fMCRecPt = new TH1F("fMCRecPt","fMCRecPt",50,0.1,3); \r
+ fMCRecPt->SetXTitle("p_{t}");\r
+ fMCRecPt->SetYTitle("yield");\r
+\r
+ fMCRecPrimPt = new TH1F("fMCRecPrimPt","fMCRecPrimPt",50,0.1,3); \r
+ fMCRecPrimPt->SetXTitle("p_{t}");\r
+ fMCRecPrimPt->SetYTitle("yield");\r
+\r
+ fMCRecSecPt = new TH1F("fMCRecSecPt","fMCRecSecPt",50,0.1,3); \r
+ fMCRecSecPt->SetXTitle("p_{t}");\r
+ fMCRecSecPt->SetYTitle("yield");\r
+\r
+ // Efficiency as function of pt\r
+ fEffTPCPt = new TProfile("Eff_pt","Eff_Pt",50,0.1,3); \r
+ fEffTPCPt->SetXTitle("p_{t}");\r
+ fEffTPCPt->SetYTitle("TPC Efficiency");\r
+\r
+ fEffTPCPtMC = new TProfile("MC_Eff_pt","MC_Eff_Pt",50,0.1,3); \r
+ fEffTPCPtMC->SetXTitle("p_{t}");\r
+ fEffTPCPtMC->SetYTitle("MC TPC Efficiency");\r
+\r
+ fEffTPCPtF = new TProfile("F_Eff_pt","F_Eff_Pt",50,0.1,3); \r
+ fEffTPCPtF->SetXTitle("p_{t}");\r
+ fEffTPCPtF->SetYTitle("TPC Findable Efficiency");\r
+\r
+ // Efficiency as function of pt protons\r
+ fEffTPCPt_P = new TProfile("Eff_pt_P","Eff_Pt_P",50,0.1,3); \r
+ fEffTPCPt_P->SetXTitle("p_{t}");\r
+ fEffTPCPt_P->SetYTitle("TPC Efficiency");\r
+\r
+ fEffTPCPtMC_P = new TProfile("MC_Eff_pt_P","MC_Eff_Pt_P",50,0.1,3); \r
+ fEffTPCPtMC_P->SetXTitle("p_{t}");\r
+ fEffTPCPtMC_P->SetYTitle("MC TPC Efficiency");\r
+\r
+ fEffTPCPtF_P = new TProfile("F_Eff_pt_P","F_Eff_Pt_P",50,0.1,3); \r
+ fEffTPCPtF_P->SetXTitle("p_{t}");\r
+ fEffTPCPtF_P->SetYTitle("TPC Findable Efficiency");\r
+\r
+ // Efficiency as function of pt pions\r
+ fEffTPCPt_Pi = new TProfile("Eff_pt_Pi","Eff_Pit_Pi",50,0.1,3); \r
+ fEffTPCPt_Pi->SetXTitle("p_{t}");\r
+ fEffTPCPt_Pi->SetYTitle("TPC Efficiency");\r
+\r
+ fEffTPCPtMC_Pi = new TProfile("MC_Eff_pt_Pi","MC_Eff_Pit_Pi",50,0.1,3); \r
+ fEffTPCPtMC_Pi->SetXTitle("p_{t}");\r
+ fEffTPCPtMC_Pi->SetYTitle("MC TPC Efficiency");\r
+\r
+ fEffTPCPtF_Pi = new TProfile("F_Eff_pt_Pi","F_Eff_Pit_Pi",50,0.1,3); \r
+ fEffTPCPtF_Pi->SetXTitle("p_{t}");\r
+ fEffTPCPtF_Pi->SetYTitle("TPC Findable Efficiency");\r
+\r
+ // Efficiency as function of pt kaons\r
+ fEffTPCPt_K = new TProfile("Eff_pt_K","Eff_Kt_K",50,0.1,3); \r
+ fEffTPCPt_K->SetXTitle("p_{t}");\r
+ fEffTPCPt_K->SetYTitle("TPC Efficiency");\r
+\r
+ fEffTPCPtMC_K = new TProfile("MC_Eff_pt_K","MC_Eff_Kt_K",50,0.1,3); \r
+ fEffTPCPtMC_K->SetXTitle("p_{t}");\r
+ fEffTPCPtMC_K->SetYTitle("MC TPC Efficiency");\r
+\r
+ fEffTPCPtF_K = new TProfile("F_Eff_pt_K","F_Eff_Kt_K",50,0.1,3); \r
+ fEffTPCPtF_K->SetXTitle("p_{t}");\r
+ fEffTPCPtF_K->SetYTitle("TPC Findable Efficiency");\r
+\r
+ // Efficiency as function of tan(theta) \r
+ fEffTPCTan = new TProfile("Eff_tan","Eff_tan",50,-2.5,2.5); \r
+ fEffTPCTan->SetXTitle("tan(#theta)");\r
+ fEffTPCTan->SetYTitle("TPC Efficiency");\r
+\r
+ fEffTPCTanMC = new TProfile("MC_Eff_tan","MC_Eff_tan",50,-2.5,2.5); \r
+ fEffTPCTanMC->SetXTitle("tan(#theta)");\r
+ fEffTPCTanMC->SetYTitle("MC TPC Efficiency");\r
+\r
+ fEffTPCTanF = new TProfile("F_Eff_tan","F_Eff_tan",50,-2.5,2.5); \r
+ fEffTPCPtF->SetXTitle("tan(#theta)");\r
+ fEffTPCPtF->SetYTitle("TPC Findable Efficiency");\r
+\r
+ // Efficiency as function of pt and tan(theta) \r
+ fEffTPCPtTan = new TProfile2D("Eff_pt_tan","Eff_Pt_tan",10,0.1,3,20,-2.,2.);\r
+ fEffTPCPtTan->SetXTitle("tan(#theta)");\r
+ fEffTPCPtTan->SetYTitle("p_{t}");\r
+\r
+ fEffTPCPtTanMC = new TProfile2D("MC_Eff_pt_tan_MC","MC Eff Pt",10,0.1,3,20,-2.,2.);\r
+ fEffTPCPtTanMC->SetXTitle("tan(#theta)");\r
+ fEffTPCPtTanMC->SetYTitle("p_{t}");\r
+\r
+ fEffTPCPtTanF = new TProfile2D("MC_Eff_pt_tan_F","MC Eff Pt",10,0.1,3,20,-2.,2.);\r
+ fEffTPCPtTanF->SetXTitle("tan(#theta)");\r
+ fEffTPCPtTanF->SetYTitle("p_{t}");\r
+\r
+ char name[256];\r
+ for(Int_t i=0; i<4; ++i)\r
+ {\r
+ sprintf(name, "fTPCPtDCASigmaIdeal_%d",i);\r
+ fTPCPtDCASigmaIdeal[i] = new TH2F(name,name,50,0.1,3,100,0,100);\r
+\r
+ sprintf(name, "fTPCPtDCASigmaFull_%d",i);\r
+ fTPCPtDCASigmaFull[i] = new TH2F(name,name,50,0.1,3,100,0,100);\r
+\r
+ sprintf(name, "fTPCPtDCASigmaDay0_%d",i);\r
+ fTPCPtDCASigmaDay0[i] = new TH2F(name,name,50,0.1,3,100,0,100);\r
+\r
+ sprintf(name, "fTPCPtDCAXY_%d",i);\r
+ fTPCPtDCAXY[i]= new TH2F(name,name,50,0.1,3,100,0,100);\r
+ sprintf(name, "fTPCPtDCAZ_%d",i);\r
+ fTPCPtDCAZ[i]= new TH2F(name,name,50,0.1,3,100,0,100);\r
+\r
+ sprintf(name, "fTPCPtDCASigmaIdealPid_%d",i);\r
+ fTPCPtDCASigmaIdealPid[i] = new TH3F(name,name,50,0.1,3,100,0,100,5,0,5);\r
+\r
+ sprintf(name, "fTPCPtDCASigmaFullPid_%d",i);\r
+ fTPCPtDCASigmaFullPid[i]= new TH3F(name,name,50,0.1,3,100,0,100,5,0,5);\r
+\r
+ sprintf(name, "fTPCPtDCASigmaDay0Pid_%d",i);\r
+ fTPCPtDCASigmaDay0Pid[i]= new TH3F(name,name,50,0.1,3,100,0,100,5,0,5);\r
+\r
+ sprintf(name, "fTPCPtDCAXYPid_%d",i);\r
+ fTPCPtDCAXYPid[i]= new TH3F(name,name,50,0.1,3,100,0,100,5,0,5);\r
+\r
+ sprintf(name, "fTPCPtDCAZPid_%d",i);\r
+ fTPCPtDCAZPid[i]= new TH3F(name,name,50,0.1,3,100,0,100,5,0,5);\r
+ }\r
+\r
+ // init cuts\r
+ if(!fCutsMC) \r
+ AliDebug(AliLog::kError, "ERROR: Cannot find AliMCInfoCuts object");\r
+ if(!fCutsRC) \r
+ AliDebug(AliLog::kError, "ERROR: Cannot find AliRecInfoCuts object");\r
+\r
+ // init folder\r
+ fAnalysisFolder = CreateFolder("folderEff","Analysis Efficiency Folder");\r
+}\r
+\r
+//_____________________________________________________________________________\r
+void AliComparisonEff::Process(AliMCInfo* infoMC, AliESDRecInfo *infoRC)\r
+{\r
+ // Fill efficiency comparison information\r
+ \r
+ AliExternalTrackParam *track = 0;\r
+ Double_t field = AliTracker::GetBz(); // nominal Bz field [kG]\r
+ Double_t kMaxD = 123456.0; // max distance\r
+\r
+ // systematics\r
+ const Double_t kSigma2Full_xy = 0.25; // ExB full systematics [cm]\r
+ const Double_t kSigma2Full_z = 5.0; // drift velocity (goofie) [cm] \r
+\r
+ const Double_t kSigma2Day0_xy = 0.02; // ExB [cm]\r
+ const Double_t kSigma2Day0_z = 0.1; // drift velocity (goofie) [cm] \r
+\r
+ // \r
+ Double_t DCASigmaIdeal=0;\r
+ Double_t DCASigmaFull=0;\r
+ Double_t DCASigmaDay0=0;\r
+\r
+ Double_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z\r
+\r
+ // distance to Prim. vertex \r
+ const Double_t* dv = infoMC->GetVDist(); \r
+\r
+ Float_t mcpt = infoMC->GetParticle().Pt();\r
+ Float_t tantheta = TMath::Tan(infoMC->GetParticle().Theta()-TMath::Pi()*0.5);\r
+ Bool_t isPrim = TMath::Sqrt(dv[0]*dv[0] + dv[1]*dv[1])<fCutsMC->GetMaxR() && TMath::Abs(dv[2])<fCutsMC->GetMaxVz();\r
+ \r
+ // calculate and set prim. vertex\r
+ fVertex->SetXv( infoMC->GetParticle().Vx() - dv[0] );\r
+ fVertex->SetYv( infoMC->GetParticle().Vy() - dv[1] );\r
+ fVertex->SetZv( infoMC->GetParticle().Vz() - dv[2] );\r
+ \r
+ // Check selection cuts\r
+ if (fCutsMC->IsPdgParticle(TMath::Abs(infoMC->GetParticle().GetPdgCode())) == kFALSE) return; \r
+\r
+ // transform Pdg to Pid\r
+ // Pdg convension is different for hadrons and leptons \r
+ // (e.g. K+/K- = 321/-321; e+/e- = -11/11 ) \r
+ Double_t pid = -1;\r
+ if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetEM() ) pid = 0; \r
+ if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetMuM() ) pid = 1; \r
+ if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetPiP() ) pid = 2; \r
+ if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetKP() ) pid = 3; \r
+ if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetProt() ) pid = 4; \r
+\r
+ //cout << "dv[0] " << dv[0] << " dv[1] " << dv[1] << " dv[2] " << dv[2] << endl; \r
+ //cout << "v[0] " << fVertex->GetXv() << " v[1] " << fVertex->GetYv() << " v[2] " << fVertex->GetZv()<< endl; \r
+ if (TMath::Abs(tantheta)<fCutsRC->GetMaxAbsTanTheta())\r
+ {\r
+ if (infoRC->GetESDtrack()->GetTPCInnerParam()) \r
+ {\r
+ if ((track = new AliExternalTrackParam(*infoRC->GetESDtrack()->GetTPCInnerParam())) != 0 )\r
+ {\r
+ Bool_t bDCAStatus = track->PropagateToDCA(fVertex,field,kMaxD,dca,cov);\r
+ if(bDCAStatus) {\r
+ //\r
+ cov[2] = track->GetCovariance()[2];\r
+\r
+ // Eff = infoRC->GetStatus(1)==3 && isPrim / isPrim;\r
+ // Pt vs ( dca[0]^2/cov[0]^2 + dca[1]^2/cov[2]^2 ) \r
+ // Pt vs ( dca[0]^2/(cov[0]^2 + kSigma2Full_xy) + dca[1]^2/( cov[2]^2 + kSigma2Full_z ) \r
+ // Pt vs ( dca[0]^2/(cov[0]^2 + kSigma2_xy) + dca[1]^2/( cov[2]^2 + kSigma2_z ) \r
+\r
+ if(cov[0]>0.0 && cov[2]>0.0)\r
+ {\r
+ DCASigmaIdeal = TMath::Power(dca[0],2)/cov[0] \r
+ + TMath::Power(dca[1],2)/cov[2]; \r
+\r
+ DCASigmaFull = TMath::Power(dca[0],2)/(cov[0]+kSigma2Full_xy) \r
+ + TMath::Power(dca[1],2)/(cov[2]+kSigma2Full_z); \r
+\r
+ DCASigmaDay0 = TMath::Power(dca[0],2)/(cov[0]+kSigma2Day0_xy) \r
+ + TMath::Power(dca[1],2)/(cov[2]+kSigma2Day0_z); \r
+\r
+ //cout << "dca[0] " << dca[0] << " dca[1] " << dca[1] << endl; \r
+ //cout << "cov[0] " << cov[0] << " cov[2] " << cov[2] << endl; \r
+ //cout << "DCASigmaIdeal " << DCASigmaIdeal << " DCASigmaFull " << DCASigmaFull << " DCASigmaDay0 " <<DCASigmaDay0 << endl; \r
+ //cout << " -------------------------------------------------------- "<< endl; \r
+ }\r
+\r
+ // MC pt\r
+ fMCPt->Fill(mcpt); \r
+ if(infoRC->GetStatus(1)==3) fMCRecPt->Fill(mcpt); \r
+ if(infoRC->GetStatus(1)==3 && isPrim) fMCRecPrimPt->Fill(mcpt); \r
+ if(infoRC->GetStatus(1)==3 && !isPrim) fMCRecSecPt->Fill(mcpt); \r
+\r
+ if(isPrim)\r
+ {\r
+ fTPCPtDCASigmaIdeal[0]->Fill(mcpt,DCASigmaIdeal);\r
+ fTPCPtDCASigmaFull[0]->Fill(mcpt,DCASigmaFull);\r
+ fTPCPtDCASigmaDay0[0]->Fill(mcpt,DCASigmaDay0);\r
+\r
+ fTPCPtDCAXY[0]->Fill(mcpt,dca[0]);\r
+ fTPCPtDCAZ[0]->Fill(mcpt,dca[1]);\r
+\r
+ fTPCPtDCASigmaIdealPid[0]->Fill(mcpt,DCASigmaIdeal,pid);\r
+ fTPCPtDCASigmaFullPid[0]->Fill(mcpt,DCASigmaFull,pid);\r
+ fTPCPtDCASigmaDay0Pid[0]->Fill(mcpt,DCASigmaDay0,pid);\r
+\r
+ fTPCPtDCAXYPid[0]->Fill(mcpt,dca[0],pid);\r
+ fTPCPtDCAZPid[0]->Fill(mcpt,dca[1],pid);\r
+ \r
+ if(infoRC->GetStatus(1)==3)\r
+ {\r
+ fTPCPtDCASigmaIdeal[1]->Fill(mcpt,DCASigmaIdeal);\r
+ fTPCPtDCASigmaFull[1]->Fill(mcpt,DCASigmaFull);\r
+ fTPCPtDCASigmaDay0[1]->Fill(mcpt,DCASigmaDay0);\r
+\r
+ fTPCPtDCAXY[1]->Fill(mcpt,dca[0]);\r
+ fTPCPtDCAZ[1]->Fill(mcpt,dca[1]);\r
+\r
+ fTPCPtDCASigmaIdealPid[1]->Fill(mcpt,DCASigmaIdeal,pid);\r
+ fTPCPtDCASigmaFullPid[1]->Fill(mcpt,DCASigmaFull,pid);\r
+ fTPCPtDCASigmaDay0Pid[1]->Fill(mcpt,DCASigmaDay0,pid);\r
+\r
+ fTPCPtDCAXYPid[1]->Fill(mcpt,dca[0],pid);\r
+ fTPCPtDCAZPid[1]->Fill(mcpt,dca[1],pid);\r
+ }\r
+ }\r
+\r
+ // Cont = infoRC->GetStatus(1)==3 && !isPrim / infoRC->GetStatus(1)==3 \r
+ // Pt vs ( dz[0]^2/cov[0]^2 + dz[1]^2/cov[2]^2 ) \r
+ // Pt vs ( dz[0]^2/(cov[0]^2 + kSigma2Full_xy) + dz[1]^2/( cov[2]^2 + kSigma2Full_z ) \r
+ // Pt vs ( dz[0]^2/(cov[0]^2 + kSigma2_xy) + dz[1]^2/( cov[2]^2 + kSigma2_z ) \r
+\r
+ if(infoRC->GetStatus(1)==3) \r
+ {\r
+ fTPCPtDCASigmaIdeal[2]->Fill(mcpt,DCASigmaIdeal);\r
+ fTPCPtDCASigmaFull[2]->Fill(mcpt,DCASigmaFull);\r
+ fTPCPtDCASigmaDay0[2]->Fill(mcpt,DCASigmaDay0);\r
+\r
+ fTPCPtDCAXY[2]->Fill(mcpt,dca[0]);\r
+ fTPCPtDCAZ[2]->Fill(mcpt,dca[1]);\r
+\r
+ fTPCPtDCASigmaIdealPid[2]->Fill(mcpt,DCASigmaIdeal,pid);\r
+ fTPCPtDCASigmaFullPid[2]->Fill(mcpt,DCASigmaFull,pid);\r
+ fTPCPtDCASigmaDay0Pid[2]->Fill(mcpt,DCASigmaDay0,pid);\r
+\r
+ fTPCPtDCAXYPid[2]->Fill(mcpt,dca[0],pid);\r
+ fTPCPtDCAZPid[2]->Fill(mcpt,dca[1],pid);\r
+ \r
+ if(isPrim==0)\r
+ {\r
+ fTPCPtDCASigmaIdeal[3]->Fill(mcpt,DCASigmaIdeal);\r
+ fTPCPtDCASigmaFull[3]->Fill(mcpt,DCASigmaFull);\r
+ fTPCPtDCASigmaDay0[3]->Fill(mcpt,DCASigmaDay0);\r
+\r
+ fTPCPtDCAXY[3]->Fill(mcpt,dca[0]);\r
+ fTPCPtDCAZ[3]->Fill(mcpt,dca[1]);\r
+\r
+ fTPCPtDCASigmaIdealPid[3]->Fill(mcpt,DCASigmaIdeal,pid);\r
+ fTPCPtDCASigmaFullPid[3]->Fill(mcpt,DCASigmaFull,pid);\r
+ fTPCPtDCASigmaDay0Pid[3]->Fill(mcpt,DCASigmaDay0,pid);\r
+\r
+ fTPCPtDCAXYPid[3]->Fill(mcpt,dca[0],pid);\r
+ fTPCPtDCAZPid[3]->Fill(mcpt,dca[1],pid);\r
+ }\r
+ }\r
+ delete track;\r
+ }\r
+ }\r
+ } \r
+ else \r
+ {\r
+ if(isPrim)\r
+ {\r
+ fTPCPtDCASigmaIdeal[0]->Fill(mcpt,0.0);\r
+ fTPCPtDCASigmaFull[0]->Fill(mcpt,0.0);\r
+ fTPCPtDCASigmaDay0[0]->Fill(mcpt,0.0);\r
+\r
+ fTPCPtDCAXY[0]->Fill(mcpt,0.0);\r
+ fTPCPtDCAZ[0]->Fill(mcpt,0.0);\r
+\r
+ fTPCPtDCASigmaIdealPid[0]->Fill(mcpt,0.0,pid);\r
+ fTPCPtDCASigmaFullPid[0]->Fill(mcpt,0.0,pid);\r
+ fTPCPtDCASigmaDay0Pid[0]->Fill(mcpt,0.0,pid);\r
+\r
+ fTPCPtDCAXYPid[0]->Fill(mcpt,0.0,pid);\r
+ fTPCPtDCAZPid[0]->Fill(mcpt,0.0,pid);\r
+ }\r
+ }\r
+ }\r
+\r
+ // only primary particles\r
+ if (!isPrim) return;\r
+\r
+ // pt\r
+ if (TMath::Abs(tantheta)<fCutsRC->GetMaxAbsTanTheta()){\r
+\r
+ fEffTPCPt->Fill(mcpt, infoRC->GetStatus(1)==3);\r
+ fEffTPCPtMC->Fill(mcpt, infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits());\r
+ if (infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits()){\r
+ fEffTPCPtF->Fill(mcpt, infoRC->GetStatus(1)==3);\r
+ }\r
+\r
+ // protons\r
+ if(TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetProt()) { \r
+ fEffTPCPt_P->Fill(mcpt, infoRC->GetStatus(1)==3);\r
+ fEffTPCPtMC_P->Fill(mcpt, infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits());\r
+\r
+ if(infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits()) {\r
+ fEffTPCPtF_P->Fill(mcpt, infoRC->GetStatus(1)==3);\r
+ }\r
+ }\r
+\r
+ // pions\r
+ if(TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetPiP()) {\r
+ fEffTPCPt_Pi->Fill(mcpt, infoRC->GetStatus(1)==3);\r
+ fEffTPCPtMC_Pi->Fill(mcpt, infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits());\r
+\r
+ if(infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits()) {\r
+ fEffTPCPtF_Pi->Fill(mcpt, infoRC->GetStatus(1)==3);\r
+ }\r
+ }\r
+\r
+ // kaons\r
+ if(TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetKP()) {\r
+ fEffTPCPt_K->Fill(mcpt, infoRC->GetStatus(1)==3);\r
+ fEffTPCPtMC_K->Fill(mcpt, infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits());\r
+\r
+ if(infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits()) {\r
+ fEffTPCPtF_K->Fill(mcpt, infoRC->GetStatus(1)==3);\r
+ }\r
+ }\r
+ }\r
+\r
+ // theta\r
+ if (TMath::Abs(mcpt)>fCutsRC->GetPtMin()){\r
+ fEffTPCTan->Fill(tantheta, infoRC->GetStatus(1)==3);\r
+ fEffTPCTanMC->Fill(tantheta, infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits());\r
+ if (infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits()){\r
+ fEffTPCTanF->Fill(tantheta, infoRC->GetStatus(1)==3);\r
+ }\r
+ }\r
+\r
+ // pt-theta\r
+ fEffTPCPtTan->Fill(mcpt,tantheta,infoRC->GetStatus(1)==3);\r
+ fEffTPCPtTanMC->Fill(mcpt,tantheta,infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits()); \r
+ if (infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits()){\r
+ fEffTPCPtTanF->Fill(mcpt,tantheta,infoRC->GetStatus(1)==3); \r
+ }\r
+}\r
+\r
+//_____________________________________________________________________________\r
+void AliComparisonEff::Exec(AliMCInfo* infoMC, AliESDRecInfo *infoRC)\r
+{\r
+ // Process comparison information\r
+ Process(infoMC,infoRC);\r
+}\r
+\r
+//_____________________________________________________________________________\r
+Long64_t AliComparisonEff::Merge(TCollection* list) \r
+{\r
+ // Merge list of objects (needed by PROOF)\r
+\r
+ if (!list)\r
+ return 0;\r
+\r
+ if (list->IsEmpty())\r
+ return 1;\r
+\r
+ TIterator* iter = list->MakeIterator();\r
+ TObject* obj = 0;\r
+\r
+ // collection of generated histograms\r
+\r
+ Int_t count=0;\r
+ while((obj = iter->Next()) != 0) \r
+ {\r
+ AliComparisonEff* entry = dynamic_cast<AliComparisonEff*>(obj);\r
+ if (entry == 0) continue; \r
+ \r
+ fMCPt->Add(entry->fMCPt);\r
+ fMCRecPt->Add(entry->fMCRecPt);\r
+ fMCRecPrimPt->Add(entry->fMCRecPrimPt);\r
+ fMCRecSecPt->Add(entry->fMCRecSecPt);\r
+\r
+ fEffTPCPt->Add(entry->fEffTPCPt);\r
+ fEffTPCPtMC->Add(entry->fEffTPCPtMC);\r
+ fEffTPCPtF->Add(entry->fEffTPCPtF);\r
+\r
+ fEffTPCPt_P->Add(entry->fEffTPCPt_P);\r
+ fEffTPCPtMC_P->Add(entry->fEffTPCPtMC_P);\r
+ fEffTPCPtF_P->Add(entry->fEffTPCPtF_P);\r
+\r
+ fEffTPCPt_Pi->Add(entry->fEffTPCPt_Pi);\r
+ fEffTPCPtMC_Pi->Add(entry->fEffTPCPtMC_Pi);\r
+ fEffTPCPtF_Pi->Add(entry->fEffTPCPtF_Pi);\r
+\r
+ fEffTPCPt_K->Add(entry->fEffTPCPt_K);\r
+ fEffTPCPtMC_K->Add(entry->fEffTPCPtMC_K);\r
+ fEffTPCPtF_K->Add(entry->fEffTPCPtF_K);\r
+\r
+ fEffTPCTan->Add(entry->fEffTPCTan);\r
+ fEffTPCTanMC->Add(entry->fEffTPCTanMC);\r
+ fEffTPCTanF->Add(entry->fEffTPCTanF);\r
+ \r
+ fEffTPCPtTan->Add(entry->fEffTPCPtTan);\r
+ fEffTPCPtTanMC->Add(entry->fEffTPCPtTanMC);\r
+ fEffTPCPtTanF->Add(entry->fEffTPCPtTanF);\r
+ \r
+ for(Int_t i=0; i<4; ++i)\r
+ {\r
+ fTPCPtDCASigmaIdeal[i]->Add(entry->fTPCPtDCASigmaIdeal[i]);\r
+ fTPCPtDCASigmaFull[i]->Add(entry->fTPCPtDCASigmaFull[i]);\r
+ fTPCPtDCASigmaDay0[i]->Add(entry->fTPCPtDCASigmaDay0[i]);\r
+\r
+ fTPCPtDCAXY[i]->Add(entry->fTPCPtDCAXY[i]);\r
+ fTPCPtDCAZ[i]->Add(entry->fTPCPtDCAXY[i]);\r
+\r
+ fTPCPtDCASigmaIdealPid[i]->Add(entry->fTPCPtDCASigmaIdealPid[i]);\r
+ fTPCPtDCASigmaFullPid[i]->Add(entry->fTPCPtDCASigmaFullPid[i]);\r
+ fTPCPtDCASigmaDay0Pid[i]->Add(entry->fTPCPtDCASigmaDay0Pid[i]);\r
+\r
+ fTPCPtDCAXYPid[i]->Add(entry->fTPCPtDCAXYPid[i]);\r
+ fTPCPtDCAZPid[i]->Add(entry->fTPCPtDCAXYPid[i]);\r
+ }\r
+\r
+ count++;\r
+ }\r
+\r
+return count;\r
+}\r
+ \r
+//_____________________________________________________________________________\r
+void AliComparisonEff::Analyse() \r
+{\r
+ // Analyse comparison information and store output histograms\r
+ // in the folder "folderEff" \r
+ //\r
+ \r
+ TH1::AddDirectory(kFALSE);\r
+\r
+ AliComparisonEff * comp=this;\r
+ TObjArray *aFolderObj = new TObjArray;\r
+\r
+ // calculate efficiency and contamination (4 sigma) \r
+ TH1 *h_sigmaidealpid[20];\r
+ TH1 *h_sigmafullpid[20];\r
+ TH1 *h_sigmaday0pid[20];\r
+\r
+ TH1 *h_sigmaidealpidtot[4];\r
+ TH1 *h_sigmafullpidtot[4];\r
+ TH1 *h_sigmaday0pidtot[4];\r
+\r
+ char name[256];\r
+ char name1[256];\r
+ Int_t idx;\r
+\r
+ for(Int_t i=0; i<4; ++i)\r
+ {\r
+ //total\r
+ comp->fTPCPtDCASigmaIdealPid[i]->GetYaxis()->SetRange(1,4);\r
+ sprintf(name,"h_sigmaidealpidtot_%d",i);\r
+ h_sigmaidealpidtot[i] = comp->fTPCPtDCASigmaIdealPid[i]->Project3D();\r
+ h_sigmaidealpidtot[i]->SetName(name);\r
+\r
+ comp->fTPCPtDCASigmaFullPid[i]->GetYaxis()->SetRange(1,4);\r
+ sprintf(name,"h_sigmafullpidtot_%d",i);\r
+ h_sigmafullpidtot[i] = comp->fTPCPtDCASigmaFullPid[i]->Project3D();\r
+ h_sigmafullpidtot[i]->SetName(name);\r
+\r
+ comp->fTPCPtDCASigmaDay0Pid[i]->GetYaxis()->SetRange(1,4);\r
+ sprintf(name,"h_sigmaday0pidtot_%d",i);\r
+ h_sigmaday0pidtot[i] = comp->fTPCPtDCASigmaDay0Pid[i]->Project3D();\r
+ h_sigmaday0pidtot[i]->SetName(name);\r
+\r
+ // pid wise\r
+ for(Int_t j=0; j<5; ++j)\r
+ {\r
+ idx = i*5 + j;\r
+\r
+ comp->fTPCPtDCASigmaIdealPid[i]->GetYaxis()->SetRange(1,4);\r
+ comp->fTPCPtDCASigmaIdealPid[i]->GetZaxis()->SetRange(j+1,j+1);\r
+\r
+ sprintf(name,"h_sigmaidealpid_%d",idx);\r
+ h_sigmaidealpid[idx] = comp->fTPCPtDCASigmaIdealPid[i]->Project3D();\r
+ h_sigmaidealpid[idx]->SetName(name);\r
+ \r
+\r
+ comp->fTPCPtDCASigmaFullPid[i]->GetYaxis()->SetRange(1,4);\r
+ comp->fTPCPtDCASigmaFullPid[i]->GetZaxis()->SetRange(j+1,j+1);\r
+\r
+ sprintf(name,"h_sigmafullpid_%d",idx);\r
+ h_sigmafullpid[idx] = comp->fTPCPtDCASigmaFullPid[i]->Project3D();\r
+ h_sigmafullpid[idx]->SetName(name);\r
+ \r
+ comp->fTPCPtDCASigmaDay0Pid[i]->GetYaxis()->SetRange(1,4);\r
+ comp->fTPCPtDCASigmaDay0Pid[i]->GetZaxis()->SetRange(j+1,j+1);\r
+\r
+ sprintf(name,"h_sigmaday0pid_%d",idx);\r
+ h_sigmaday0pid[idx] = comp->fTPCPtDCASigmaDay0Pid[i]->Project3D();\r
+ h_sigmaday0pid[idx]->SetName(name);\r
+ } \r
+ }\r
+\r
+ // calculate efficiency and contamination (all pids)\r
+ h_sigmaidealpidtot[0]->Sumw2();\r
+ h_sigmaidealpidtot[1]->Divide(h_sigmaidealpidtot[0]);\r
+ h_sigmaidealpidtot[2]->Sumw2();\r
+ h_sigmaidealpidtot[3]->Divide(h_sigmaidealpidtot[2]);\r
+\r
+ h_sigmafullpidtot[0]->Sumw2();\r
+ h_sigmafullpidtot[1]->Divide(h_sigmafullpidtot[0]);\r
+ h_sigmafullpidtot[2]->Sumw2();\r
+ h_sigmafullpidtot[3]->Divide(h_sigmafullpidtot[2]);\r
+\r
+ h_sigmaday0pidtot[0]->Sumw2();\r
+ h_sigmaday0pidtot[1]->Divide(h_sigmaday0pidtot[0]);\r
+ h_sigmaday0pidtot[2]->Sumw2();\r
+ h_sigmaday0pidtot[3]->Divide(h_sigmaday0pidtot[2]);\r
+\r
+ // calculate efficiency pid wise\r
+ for(Int_t idx = 0; idx<5; idx++)\r
+ {\r
+ h_sigmaidealpid[idx]->Sumw2();\r
+ h_sigmaidealpid[idx+5]->Divide(h_sigmaidealpid[idx]);\r
+\r
+ h_sigmafullpid[idx]->Sumw2();\r
+ h_sigmafullpid[idx+5]->Divide(h_sigmafullpid[idx]);\r
+\r
+ h_sigmaday0pid[idx]->Sumw2();\r
+ h_sigmaday0pid[idx+5]->Divide(h_sigmaday0pid[idx]);\r
+ }\r
+\r
+ // calculate cont. pid wise\r
+ for(Int_t idx = 0; idx<5; idx++)\r
+ {\r
+ h_sigmaidealpid[idx+15]->Divide(h_sigmaidealpidtot[2]);\r
+ h_sigmafullpid[idx+15]->Divide(h_sigmafullpidtot[2]);\r
+ h_sigmaday0pid[idx+15]->Divide(h_sigmaday0pidtot[2]);\r
+ }\r
+\r
+ TCanvas * c = new TCanvas("Efficiency","Track efficiency");\r
+ c->cd();\r
+ c->Divide(1,2);\r
+\r
+ //\r
+ c->cd(1);\r
+ comp->fEffTPCTanF->SetXTitle("Tan(#theta)");\r
+ comp->fEffTPCTanF->SetYTitle("eff_{findable}");\r
+ comp->fEffTPCTanF->SetName("EffTanFindable");\r
+ comp->fEffTPCTanF->Draw();\r
+ //\r
+ c->cd(2);\r
+ comp->fEffTPCTan->SetXTitle("Tan(#theta)");\r
+ comp->fEffTPCTan->SetYTitle("eff_{all}");\r
+ comp->fEffTPCTan->SetName("EffTanAll");\r
+ comp->fEffTPCTan->Draw();\r
+\r
+ aFolderObj->Add(comp->fEffTPCTanF);\r
+ aFolderObj->Add(comp->fEffTPCTan);\r
+\r
+ h_sigmaidealpidtot[1]->SetXTitle("p_{t}");\r
+ h_sigmaidealpidtot[1]->SetYTitle("efficiency");\r
+ h_sigmaidealpidtot[1]->SetTitle("Eff_SigmaIdeal");\r
+ h_sigmaidealpidtot[1]->SetName("Eff_SigmaIdeal");\r
+\r
+ h_sigmaidealpidtot[3]->SetXTitle("p_{t}");\r
+ h_sigmaidealpidtot[3]->SetYTitle("contamination");\r
+ h_sigmaidealpidtot[3]->SetTitle("Cont_SigmaIdeal");\r
+ h_sigmaidealpidtot[3]->SetName("Cont_SigmaIdeal");\r
+\r
+ aFolderObj->Add(h_sigmaidealpidtot[1]);\r
+ aFolderObj->Add(h_sigmaidealpidtot[3]);\r
+\r
+ h_sigmafullpidtot[1]->SetXTitle("p_{t}");\r
+ h_sigmafullpidtot[1]->SetYTitle("efficiency");\r
+ h_sigmafullpidtot[1]->SetTitle("Eff_SigmaFull");\r
+ h_sigmafullpidtot[1]->SetName("Eff_SigmaFull");\r
+\r
+ h_sigmafullpidtot[3]->SetXTitle("p_{t}");\r
+ h_sigmafullpidtot[3]->SetYTitle("contamination");\r
+ h_sigmafullpidtot[3]->SetTitle("Cont_SigmaFull");\r
+ h_sigmafullpidtot[3]->SetName("Cont_SigmaFull");\r
+\r
+ aFolderObj->Add(h_sigmafullpidtot[1]);\r
+ aFolderObj->Add(h_sigmafullpidtot[3]);\r
+\r
+ h_sigmaday0pidtot[1]->SetXTitle("p_{t}");\r
+ h_sigmaday0pidtot[1]->SetYTitle("efficiency");\r
+ h_sigmaday0pidtot[1]->SetTitle("Eff_SigmaDay0");\r
+ h_sigmaday0pidtot[1]->SetName("Eff_SigmaDay0");\r
+\r
+ h_sigmaday0pidtot[3]->SetXTitle("p_{t}");\r
+ h_sigmaday0pidtot[3]->SetYTitle("contamination");\r
+ h_sigmaday0pidtot[3]->SetTitle("Cont_SigmaDay0");\r
+ h_sigmaday0pidtot[3]->SetName("Cont_SigmaDay0");\r
+\r
+ aFolderObj->Add(h_sigmaday0pidtot[1]);\r
+ aFolderObj->Add(h_sigmaday0pidtot[3]);\r
+\r
+ for(Int_t idx = 0; idx<5; idx++)\r
+ {\r
+ sprintf(name,"Eff_SigmaIdeal_%d",idx);\r
+ sprintf(name1,"Cont_SigmaIdeal_%d",idx);\r
+\r
+\r
+ h_sigmaidealpid[idx+5]->SetXTitle("p_{t}");\r
+ h_sigmaidealpid[idx+5]->SetYTitle("efficiency");\r
+ h_sigmaidealpid[idx+5]->SetTitle(name);\r
+ h_sigmaidealpid[idx+5]->SetName(name);\r
+\r
+ h_sigmaidealpid[idx+15]->SetXTitle("p_{t}");\r
+ h_sigmaidealpid[idx+15]->SetYTitle("contamination");\r
+ h_sigmaidealpid[idx+15]->SetTitle(name1);\r
+ h_sigmaidealpid[idx+15]->SetName(name1);\r
+\r
+ aFolderObj->Add(h_sigmaidealpid[idx+5]);\r
+ aFolderObj->Add(h_sigmaidealpid[idx+15]);\r
+\r
+ sprintf(name,"Eff_SigmaFull_%d",idx);\r
+ sprintf(name1,"Cont_SigmaFull_%d",idx);\r
+\r
+ h_sigmafullpid[idx+5]->SetXTitle("p_{t}");\r
+ h_sigmafullpid[idx+5]->SetYTitle("efficiency");\r
+ h_sigmafullpid[idx+5]->SetTitle(name);\r
+ h_sigmafullpid[idx+5]->SetName(name);\r
+\r
+ h_sigmafullpid[idx+15]->SetXTitle("p_{t}");\r
+ h_sigmafullpid[idx+15]->SetYTitle("contamination");\r
+ h_sigmafullpid[idx+15]->SetTitle(name1);\r
+ h_sigmafullpid[idx+15]->SetName(name1);\r
+\r
+ aFolderObj->Add(h_sigmafullpid[idx+5]);\r
+ aFolderObj->Add(h_sigmafullpid[idx+15]);\r
+\r
+ sprintf(name,"Eff_SigmaDay0_%d",idx);\r
+ sprintf(name1,"Cont_SigmaDay0_%d",idx);\r
+\r
+ h_sigmaday0pid[idx+5]->SetXTitle("p_{t}");\r
+ h_sigmaday0pid[idx+5]->SetYTitle("efficiency");\r
+ h_sigmaday0pid[idx+5]->SetTitle(name);\r
+ h_sigmaday0pid[idx+5]->SetName(name);\r
+\r
+ h_sigmaday0pid[idx+15]->SetXTitle("p_{t}");\r
+ h_sigmaday0pid[idx+15]->SetYTitle("contamination");\r
+ h_sigmaday0pid[idx+15]->SetTitle(name1);\r
+ h_sigmaday0pid[idx+15]->SetName(name1);\r
+\r
+ aFolderObj->Add(h_sigmaday0pid[idx+5]);\r
+ aFolderObj->Add(h_sigmaday0pid[idx+15]);\r
+ }\r
+\r
+ // export objects to analysis folder\r
+ fAnalysisFolder = ExportToFolder(aFolderObj);\r
+\r
+ // delete only TObjArray\r
+ if(aFolderObj) delete aFolderObj;\r
+}\r
+\r
+//_____________________________________________________________________________\r
+TFolder* AliComparisonEff::ExportToFolder(TObjArray * array) \r
+{\r
+ // recreate folder avery time and export objects to new one\r
+ //\r
+ AliComparisonEff * comp=this;\r
+ TFolder *folder = comp->GetAnalysisFolder();\r
+\r
+ TString name, title;\r
+ TFolder *newFolder = 0;\r
+ Int_t i = 0;\r
+ Int_t size = array->GetSize();\r
+\r
+ if(folder) { \r
+ // get name and title from old folder\r
+ name = folder->GetName(); \r
+ title = folder->GetTitle(); \r
+\r
+ // delete old one\r
+ delete folder;\r
+\r
+ // create new one\r
+ newFolder = CreateFolder(name.Data(),title.Data());\r
+ newFolder->SetOwner();\r
+\r
+ // add objects to folder\r
+ while(i < size) {\r
+ newFolder->Add(array->At(i));\r
+ i++;\r
+ }\r
+ }\r
+\r
+return newFolder;\r
+}\r
+\r
+\r
+//_____________________________________________________________________________\r
+TFolder* AliComparisonEff::CreateFolder(TString name,TString title) { \r
+// create folder for analysed histograms\r
+//\r
+TFolder *folder = 0;\r
+ folder = new TFolder(name.Data(),title.Data());\r
+\r
+ return folder;\r
+}\r
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff