-#include "TChain.h"\r
-#include "TTree.h"\r
-#include "TH1F.h"\r
-#include "TH2F.h"\r
-#include "TProfile.h"\r
-\r
-#include "TCanvas.h"\r
-#include "TList.h"\r
-#include "TParticle.h"\r
-#include "TParticlePDG.h"\r
-#include "TProfile.h"\r
-\r
-#include "AliAnalysisTask.h"\r
-#include "AliAnalysisManager.h"\r
-\r
-#include "AliESDEvent.h"\r
-#include "AliStack.h"\r
-#include "AliESDVertex.h"\r
-#include "AliESDInputHandler.h"\r
-#include "AliESDtrackCuts.h"\r
-#include "AliMultiplicity.h"\r
-\r
-#include "AliMCParticle.h"\r
-#include "AliMCEvent.h"\r
-#include "AliAnalysisTaskEfficiency.h"\r
-#include "AliExternalTrackParam.h"\r
-#include "AliTrackReference.h"\r
-#include "AliStack.h"\r
-#include "AliHeader.h"\r
-#include "AliGenEventHeader.h"\r
-#include "AliGenCocktailEventHeader.h"\r
-#include "AliGenPythiaEventHeader.h"\r
-#include "AliGenDPMjetEventHeader.h"\r
-\r
-// Analysis Task for basic QA on the ESD\r
-// Authors: Jan Fiete Grosse-Oetringhaus, Christian Klein-Boesing, Veronica Canoa, AM, Eva Sicking\r
-\r
-ClassImp(AliAnalysisTaskEfficiency)\r
-\r
-//________________________________________________________________________\r
- AliAnalysisTaskEfficiency::AliAnalysisTaskEfficiency(const char *name) \r
- : AliAnalysisTaskSE(name) \r
- ,fFieldOn(kTRUE)\r
- ,fHists(0)\r
- ,fHistRECpt(0)\r
- ,fHistMCpt(0)\r
- ,fHistMCNRpt(0)\r
- ,fHistFAKEpt(0)\r
- ,fHistMULTpt(0)\r
- ,fHistRECeta(0)\r
- ,fHistMCeta(0)\r
- ,fHistMCNReta(0)\r
- ,fHistFAKEeta(0)\r
- ,fHistMULTeta(0)\r
- ,fHistRECphi(0)\r
- ,fHistMCphi(0)\r
- ,fHistMCNRphi(0)\r
- ,fHistFAKEphi(0)\r
- ,fHistMULTphi(0)\r
- ,fHistRECHPTeta(0)\r
- ,fHistMCHPTeta(0)\r
- ,fHistMCNRHPTeta(0)\r
- ,fHistFAKEHPTeta(0)\r
- ,fHistRECHPTphi(0)\r
- ,fHistMCHPTphi(0)\r
- ,fHistMCNRHPTphi(0)\r
- ,fHistFAKEHPTphi(0)\r
- ,fHistRecMult(0)\r
- ,fHistNCluster(0) \r
- ,fh1VtxEff(0) \r
- ,fh2VtxTpcSpd(0)\r
- ,fh2EtaCorrelation(0) \r
- ,fh2EtaCorrelationShift(0) \r
- ,fh2PhiCorrelation(0) \r
- ,fh2PhiCorrelationShift(0) \r
- ,fh2PtCorrelation(0) \r
- ,fh2PtCorrelationShift(0)\r
- ,fHistDeltaphiprimaries(0)\r
- ,fHistDeltaphisecondaries(0) \r
- ,fHistDeltaphireject(0) \r
- ,fHistTPCITSdeltax(0)\r
- ,fHistTPCITSdeltay(0)\r
- ,fHistTPCITSdeltaz(0)\r
- ,fHistTPCITSdeltar(0)\r
- ,fHistTPCTRDdeltax(0)\r
- ,fHistTPCTRDdeltay(0)\r
- ,fHistTPCTRDdeltaz(0)\r
- ,fHistTPCTRDdeltar(0)\r
- ,fHistTPCTRDdeltaphi(0)\r
- ,fPtPullsPos(0)\r
- ,fPtPullsNeg(0)\r
- ,fPtShiftPos(0)\r
- ,fPtShiftNeg(0)\r
- ,fTrackType(0)\r
- ,fCuts(0)\r
- ,fVertexRvsZ(0)\r
- ,fVertexRvsZC(0)\r
- ,fEtaMultiFluc(0)\r
- ,fEtaMulti(0)\r
- ,fEtaMultiH(0)\r
- ,fPhiMultiFluc(0)\r
- ,fPhiMulti(0)\r
- ,fPhiMultiH(0)\r
-{\r
- //\r
- // Constructor\r
- //\r
-\r
- for(int i = 0;i< 2;++i){\r
- fh2MultSpdChips[i] = 0;\r
- }\r
-\r
- for(int i = 0;i < kTotalAC;++i){\r
- fh2PhiPadRow[i] = fh2PhiLayer[i] = 0;\r
- }\r
-\r
- for(int i = 0; i < kTotalVtx; ++i){ \r
- fh2VertexCorrelation[i] = \r
- fh2VertexCorrelationShift[i] = 0;\r
- fh1VertexShift[i] = \r
- fh1VertexShiftNorm[i] = 0; \r
- }\r
-\r
- for(int i = 0;i< 8;++i){\r
- fHistRECptCharge[i]=0;\r
- fHistMCptCharge[i]=0;\r
- fHistMCNRptCharge[i]=0;\r
- fHistFAKEptCharge[i]=0;\r
-\r
- fHistRECetaCharge[i]=0;\r
- fHistMCetaCharge[i]=0;\r
- fHistMCNRetaCharge[i]=0;\r
- fHistFAKEetaCharge[i]=0;\r
-\r
- fHistRECphiCharge[i]=0;\r
- fHistMCphiCharge[i]=0;\r
- fHistMCNRphiCharge[i]=0;\r
- fHistFAKEphiCharge[i]=0;\r
-\r
- }\r
-\r
- DefineOutput(1, TList::Class());\r
-}\r
-\r
-\r
-//________________________________________________________________________\r
-void AliAnalysisTaskEfficiency::UserCreateOutputObjects()\r
-{\r
- // Create histograms\r
-\r
-\r
- fHists = new TList();\r
- fHistRECpt = new TH1F("fHistRECpt", " p_{T} distribution: all reconstructed", 100, 0., 20.);\r
- fHistMCpt = new TH1F("fHistMCpt", " p_{T} distribution: all MC", 100, 0., 20.);\r
- fHistMCNRpt = new TH1F("fHistMCNRpt", " p_{T} distribution: all not-reconstructed MC", 100, 0., 20.);\r
- fHistFAKEpt = new TH1F("fHistFAKEpt", " p_{T} distribution: all Fake", 100, 0., 20.);\r
- fHistMULTpt = new TH1F("fHistMULTpt", " p_{T} distribution: multiply rec.", 100, 0., 20.);\r
- \r
- fHistRECeta = new TH1F("fHistRECeta", " #eta distribution: all reconstructed", 100,-1.0,1.0);\r
- fHistMCeta = new TH1F("fHistMCeta", " #eta distribution: all MC", 100,-1.0,1.0);\r
- fHistMCNReta = new TH1F("fHistMCNReta", " #eta distribution: all not-reconstructed MC", 100,-1.0,1.0);\r
- fHistFAKEeta = new TH1F("fHistFAKEeta", " #eta distribution: all Fake", 100,-1.0,1.0);\r
- fHistMULTeta = new TH1F("fHistMULTeta", " #eta distribution: multiply rec.", 100,-1.0,1.0);\r
-\r
- fHistRECphi = new TH1F("fHistRECphi", " #phi distribution: all reconstructed", 314, 0., 6.28);\r
- fHistMCphi = new TH1F("fHistMCphi", " #phi distribution: all MC", 314, 0., 6.28);\r
- fHistMCNRphi = new TH1F("fHistMCNRphi", " #phi distribution: all not-reconstructed MC", 314, 0., 6.28);\r
- fHistFAKEphi = new TH1F("fHistFAKEphi", " #phi distribution: all Fake", 314, 0., 6.28);\r
- fHistMULTphi = new TH1F("fHistMULTphi", " #phi distribution: multipli rec.", 314, 0., 6.28);\r
-\r
- fHistRECHPTeta = new TH1F("fHistRECHPTeta", " #eta distribution: all reconstructed", 100,-1.0,1.0);\r
- fHistMCHPTeta = new TH1F("fHistMCHPTeta", " #eta distribution: all MC", 100,-1.0,1.0);\r
- fHistMCNRHPTeta = new TH1F("fHistMCNRHPTeta", " #eta distribution: all not-reconstructed MC", 100,-1.0,1.0);\r
- fHistFAKEHPTeta = new TH1F("fHistFAKEHPTeta", " #eta distribution: all Fake", 100,-1.0,1.0);\r
-\r
- fHistRECHPTphi = new TH1F("fHistRECHPTphi", " #phi distribution: all reconstructed", 314, 0., 6.28);\r
- fHistMCHPTphi = new TH1F("fHistMCHPTphi", " #phi distribution: all MC", 314, 0., 6.28);\r
- fHistMCNRHPTphi = new TH1F("fHistMCNRHPTphi", " #phi distribution: all not-reconstructed MC", 314, 0., 6.28);\r
- fHistFAKEHPTphi = new TH1F("fHistFAKEHPTphi", " #phi distribution: all Fake", 314, 0., 6.28);\r
-\r
- fHistRecMult = new TH1F("fHistRecMult", "Multiple reconstructed tracks", 50, 0., 50.);\r
- fHistNCluster = new TH1F("fHistNCluster", "Number of clusters for suspicious tracks", 300, 0., 300.);\r
-\r
- // CKB\r
-\r
- fh1VtxEff = new TH1F("fh1VtxEff","VtxEff TPC",4,-0.5,3.5);\r
- fh1VtxEff->GetXaxis()->SetBinLabel(1,"NO TPC VTX");\r
- fh1VtxEff->GetXaxis()->SetBinLabel(2,"TPC VTX");\r
- fh1VtxEff->GetXaxis()->SetBinLabel(3,"NO SPD VTX");\r
- fh1VtxEff->GetXaxis()->SetBinLabel(4,"SPD VTX");\r
-\r
- TString labels[kTotalAC];\r
- labels[kNegA]="NegA";\r
- labels[kPosA]="PosA";\r
- labels[kNegC]="NegC";\r
- labels[kPosC]="PosC";\r
-\r
- for(int i = 0;i< kTotalAC;++i){\r
- fh2PhiPadRow[i] = new TH2F(Form("fh2PhiPadRow_%d",i),Form("Padrow vs phi AC %s;phi;padrow",labels[i].Data()),360,0.,360.,159,-0.5,158.5);\r
- fh2PhiLayer[i] = new TH2F(Form("fh2PhiLayer_%d",i),Form("layer vs phi AC %s;phi;layer",labels[i].Data()),360,0.,360.,6,-0.5,5.5);\r
- }\r
- for(int i = 0;i<2;++i){\r
- fh2MultSpdChips[i] = new TH2F(Form("fh2ChipsSpdMult_%d",i),"mult SPD vs chips;chips;Mult",201,-1.5,199.5,201,-1.5,199.5);\r
- }\r
- fh2VtxTpcSpd = new TH2F("fh2VtxTpcSpd","SPD vtx vs TPC;tpc-z;spd-z",200,-20.,20.,200,-20,20.);\r
-\r
- TString labelsVtx[kTotalVtx];\r
- labelsVtx[kTPC] = "TPC";\r
- labelsVtx[kSPD] = "SPD";\r
- for(int i = 0; i < kTotalVtx; ++i){ \r
- fh2VertexCorrelation[i] = new TH2F(Form("fh2VertexCorrelation_%d",i),Form("VertexCorrelation %s;MC z-vtx;ESD z-vtx",labelsVtx[i].Data()), 120, -30, 30, 120, -30, 30);\r
- fh2VertexCorrelationShift[i] = new TH2F(Form("fh2VertexCorrelationShift_%d",i), Form("VertexCorrelationShift %s;MC z-vtx;MC z-vtx - ESD z-vtx",labelsVtx[i].Data()), 120, -30, 30, 100, -1, 1); \r
- fh1VertexShift[i] = new TH1F(Form("fh1VertexShift_%d",i), Form("VertexShift %s;(MC z-vtx - ESD z-vtx);Entries",labelsVtx[i].Data()), 201, -2, 2); \r
- fh1VertexShiftNorm[i] = new TH1F(Form("fh1VertexShiftNorm_%d",i), Form("VertexShiftNorm %s;(MC z-vtx - ESD z-vtx) / #sigma_{ESD z-vtx};Entries",labelsVtx[i].Data()), 200, -100, 100); \r
- }\r
-\r
- \r
- fh2EtaCorrelation = new TH2F("fh2EtaCorrelation", "EtaCorrelation;MC #eta;ESD #eta", 120, -3, 3, 120, -3, 3);\r
- fh2EtaCorrelationShift = new TH2F("fh2EtaCorrelationShift", "EtaCorrelationShift;ESD #eta;MC #eta - ESD #eta", 120, -3, 3, 100, -0.1, 0.1);\r
- fh2PhiCorrelation = new TH2F("fh2PhiCorrelation", "PhiCorrelation;MC #phi;ESD #phi", 120, 0, 360, 120, 0, 360);\r
- fh2PhiCorrelationShift = new TH2F("fh2PhiCorrelationShift", "PhiCorrelationShift;ESD #phi;MC #phi - ESD #phi", 120, 0, 360, 100, -10, 10);\r
- fh2PtCorrelation = new TH2F("fh2PtCorrelation", "PtCorrelation;MC p_T;ESD p_T", 200, 0., 200., 200, 0, 200);\r
- fh2PtCorrelationShift = new TH2F("fh2PtCorrelationShift", "PtCorrelationShift;ESD p_T;MC p_T - ESD #p_T", 120, 0, 10, 100, -2, 2);\r
- \r
- fHistDeltaphiprimaries =new TH1F("fHistDeltaphiprimaries", " #Delta#phi distribution: primaries", 314, -0.2,0.2); \r
- fHistDeltaphisecondaries=new TH1F("fHistDeltaphisecondaries", "#Delta#phi distribution: secondaries", 314, -0.2,0.2); \r
- fHistDeltaphireject =new TH1F("fHistDeltaphireject", " #Delta#phi distribution: reject trackled", 314, -0.2,0.2); \r
- \r
- fHistTPCITSdeltax =new TH1F("fHistTPCITSdeltax", "TPC-ITS matching:#Delta x", 100,-20,20); \r
- fHistTPCITSdeltay =new TH1F("fHistTPCITSdeltay", "TPC-ITS matching:#Delta y", 100,-20,20); \r
- fHistTPCITSdeltaz =new TH1F("fHistTPCITSdeltaz", "TPC-ITS matching:#Delta z", 100,-20,20); \r
- fHistTPCITSdeltar =new TH1F("fHistTPCITSdeltar", "TPC-ITS matching:#Delta r", 100,0,20); \r
- \r
- fHistTPCTRDdeltax = new TH1F("fHistTPCTRDdeltax", "TPC-TRD matching:#Delta x", 400,-400, 400); \r
- fHistTPCTRDdeltay = new TH1F("fHistTPCTRDdeltay", "TPC-TRD matching:#Delta y", 400,-400, 400); \r
- fHistTPCTRDdeltaz = new TH1F("fHistTPCTRDdeltaz", "TPC-TRD matching:#Delta z", 400,-400, 400); \r
- fHistTPCTRDdeltar = new TH1F("fHistTPCTRDdeltar", "TPC-TRD matching:#Delta r", 100,0,20); \r
- fHistTPCTRDdeltaphi = new TH1F("fHistTPCTRDdeltarphi", "TPC-TRD matching:#Delta #phi", 128,-3.14 , 3.14); \r
- // Pulls\r
- fPtPullsPos = new TProfile("fPtPullsPos", "Pt Pulls for pos. primaries", 50 , 0., 10., -5., 5.,"S");\r
- fPtPullsNeg = new TProfile("fPtPullsNeg", "Pt Pulls for neg. primaries", 50 , 0., 10., -5., 5.,"S");\r
- fPtShiftPos = new TProfile("fPtShiftPos", "Pt Shift for pos. primaries", 50 , 0., 10., -0.2, 0.2,"S");\r
- fPtShiftNeg = new TProfile("fPtShiftNeg", "Pt Shift for neg. primaries", 50 , 0., 10., -0.2, 0.2,"S");\r
-\r
- fEtaMultiFluc = new TProfile("fEtaMultiFluc", "eta multiplicity fluctuations", 10, -2., 2., 0., 600., "S");\r
- fEtaMulti = new TH1F("fEtaMulti", "eta multiplicity fluctuations", 10, -2., 2.);\r
- fEtaMultiH = new TH1F("fEtaMultiH", "eta multiplicity fluctuations", 600, 0., 600.);\r
-\r
- fPhiMultiFluc = new TProfile("fPhiMultiFluc", "phi multiplicity fluctuations", 64, 0., 6.4, 0., 100., "S");\r
- fPhiMulti = new TH1F("fPhiMulti", "phi multiplicity fluctuations", 64, 0., 6.4);\r
- fPhiMultiH = new TH1F("fPhiMultiH", "phi multiplicity fluctuations", 100, 0., 100.);\r
-\r
- // SPD Vertex\r
- fVertexRvsZ = new TH2F("fVertexRvsZ", "SPD Vertex R vs Z", 200, -1., 1., 200, -1., 1.);\r
- fVertexRvsZC = new TH2F("fVertexRvsZC", "SPD Vertex R vs Z", 200, -1., 1., 200, -1., 1.);\r
- \r
- TString charge[8]={"Deu", "AntiDeu", "Tri", "AntiTri", "He3", "AntiHe3", "He4", "AntiHe4"};\r
- for(Int_t i=0;i<8;i++){\r
- fHistRECptCharge[i] = new TH1F(Form("fHistRECptCharge%s",charge[i].Data()), \r
- "p_{T} distribution: all reconstructed",\r
- 100, 0., 20.);\r
- fHistMCptCharge[i] = new TH1F(Form("fHistMCptCharge%s",charge[i].Data()), \r
- "p_{T} distribution: all MC",\r
- 100, 0., 20.);\r
- fHistMCNRptCharge[i] = new TH1F(Form("fHistMCNRptCharge%s",charge[i].Data()), \r
- "p_{T} distribution: all not-reconstructed MC", \r
- 100, 0., 20.);\r
- fHistFAKEptCharge[i] = new TH1F( Form("fHistFAKEptCharge%s",charge[i].Data()), \r
- "p_{T} distribution: all Fake", \r
- 100, 0., 20.);\r
-\r
- fHistRECetaCharge[i] = new TH1F(Form("fHistRECetaCharge%s",charge[i].Data()), \r
- "p_{T} distribution: all reconstructed",\r
- 100, 0., 20.);\r
- fHistMCetaCharge[i] = new TH1F(Form("fHistMCetaCharge%s",charge[i].Data()), \r
- "p_{T} distribution: all MC",\r
- 100, 0., 20.);\r
- fHistMCNRetaCharge[i] = new TH1F(Form("fHistMCNRetaCharge%s",charge[i].Data()), \r
- "p_{T} distribution: all not-reconstructed MC", \r
- 100, 0., 20.);\r
- fHistFAKEetaCharge[i] = new TH1F( Form("fHistFAKEetaCharge%s",charge[i].Data()), \r
- "p_{T} distribution: all Fake", \r
- 100, 0., 20.);\r
-\r
- fHistRECphiCharge[i] = new TH1F(Form("fHistRECphiCharge%s",charge[i].Data()), \r
- "p_{T} distribution: all reconstructed",\r
- 100, 0., 20.);\r
- fHistMCphiCharge[i] = new TH1F(Form("fHistMCphiCharge%s",charge[i].Data()), \r
- "p_{T} distribution: all MC",\r
- 100, 0., 20.);\r
- fHistMCNRphiCharge[i] = new TH1F(Form("fHistMCNRphiCharge%s",charge[i].Data()), \r
- "p_{T} distribution: all not-reconstructed MC", \r
- 100, 0., 20.);\r
- fHistFAKEphiCharge[i] = new TH1F( Form("fHistFAKEphiCharge%s",charge[i].Data()), \r
- "p_{T} distribution: all Fake", \r
- 100, 0., 20.);\r
-\r
- }\r
- // BKC\r
-\r
-\r
- fHists->SetOwner();\r
-\r
- fHists->Add(fHistRECpt);\r
- fHists->Add(fHistMCpt);\r
- fHists->Add(fHistMCNRpt);\r
- fHists->Add(fHistFAKEpt);\r
- fHists->Add(fHistMULTpt);\r
-\r
- fHists->Add(fHistRECeta);\r
- fHists->Add(fHistMCeta);\r
- fHists->Add(fHistMCNReta);\r
- fHists->Add(fHistFAKEeta);\r
- fHists->Add(fHistMULTeta);\r
-\r
- fHists->Add(fHistRECphi);\r
- fHists->Add(fHistMCphi);\r
- fHists->Add(fHistMCNRphi);\r
- fHists->Add(fHistFAKEphi);\r
- fHists->Add(fHistMULTphi);\r
-\r
- fHists->Add(fHistRECHPTeta);\r
- fHists->Add(fHistMCHPTeta);\r
- fHists->Add(fHistMCNRHPTeta);\r
- fHists->Add(fHistFAKEHPTeta);\r
-\r
- fHists->Add(fHistRECHPTphi);\r
- fHists->Add(fHistMCHPTphi);\r
- fHists->Add(fHistMCNRHPTphi);\r
- fHists->Add(fHistFAKEHPTphi);\r
-\r
- // CKB\r
-\r
- fHists->Add(fh1VtxEff);\r
- for(int i = 0;i < kTotalAC;++i){\r
- fHists->Add(fh2PhiPadRow[i]);\r
- fHists->Add(fh2PhiLayer[i]);\r
- }\r
- for(int i = 0;i<2;++i){\r
- fHists->Add(fh2MultSpdChips[i]);\r
- }\r
- fHists->Add(fh2VtxTpcSpd);\r
- \r
- for(int i = 0; i < kTotalVtx; ++i){ \r
- fHists->Add(fh2VertexCorrelation[i]);\r
- fHists->Add(fh2VertexCorrelationShift[i]);\r
- fHists->Add(fh1VertexShift[i]); \r
- fHists->Add(fh1VertexShiftNorm[i]);\r
- }\r
-\r
-\r
- fHists->Add(fh2EtaCorrelation);\r
- fHists->Add(fh2EtaCorrelationShift);\r
- fHists->Add(fh2PhiCorrelation);\r
- fHists->Add(fh2PhiCorrelationShift);\r
- fHists->Add(fh2PtCorrelation);\r
- fHists->Add(fh2PtCorrelationShift);\r
-\r
- fHists->Add(fHistDeltaphiprimaries); \r
- fHists->Add(fHistDeltaphisecondaries);\r
- fHists->Add(fHistDeltaphireject); \r
-\r
- fHists->Add(fHistTPCITSdeltax); \r
- fHists->Add(fHistTPCITSdeltay); \r
- fHists->Add(fHistTPCITSdeltaz); \r
- fHists->Add(fHistTPCITSdeltar); \r
- \r
- fHists->Add(fHistTPCTRDdeltax); \r
- fHists->Add(fHistTPCTRDdeltay); \r
- fHists->Add(fHistTPCTRDdeltaz); \r
- fHists->Add(fHistTPCTRDdeltar); \r
- fHists->Add(fHistTPCTRDdeltaphi); \r
-\r
- fHists->Add(fHistRecMult);\r
- fHists->Add(fHistNCluster);\r
- \r
-\r
- fHists->Add(fPtPullsPos);\r
- fHists->Add(fPtPullsNeg);\r
- fHists->Add(fPtShiftPos);\r
- fHists->Add(fPtShiftNeg);\r
-\r
- fHists->Add(fVertexRvsZ);\r
- fHists->Add(fVertexRvsZC);\r
- fHists->Add(fEtaMultiFluc);\r
- fHists->Add(fEtaMulti);\r
- fHists->Add(fEtaMultiH);\r
- fHists->Add(fPhiMultiFluc);\r
- fHists->Add(fPhiMulti);\r
- fHists->Add(fPhiMultiH);\r
-\r
- for(Int_t i=0;i<8;i++){\r
- fHists->Add(fHistRECptCharge[i]);\r
- fHists->Add(fHistMCptCharge[i]);\r
- fHists->Add(fHistMCNRptCharge[i]);\r
- fHists->Add(fHistFAKEptCharge[i]);\r
-\r
- fHists->Add(fHistRECetaCharge[i]);\r
- fHists->Add(fHistMCetaCharge[i]);\r
- fHists->Add(fHistMCNRetaCharge[i]);\r
- fHists->Add(fHistFAKEetaCharge[i]);\r
-\r
- fHists->Add(fHistRECphiCharge[i]);\r
- fHists->Add(fHistMCphiCharge[i]);\r
- fHists->Add(fHistMCNRphiCharge[i]);\r
- fHists->Add(fHistFAKEphiCharge[i]);\r
- }\r
-\r
- for (Int_t i=0; i<fHists->GetEntries(); ++i) {\r
- TH1 *h1 = dynamic_cast<TH1*>(fHists->At(i));\r
- if (h1){\r
- h1->Sumw2();\r
- }\r
- }\r
- // BKC\r
-\r
-\r
- // Post output data.\r
- PostData(1, fHists);\r
-\r
-}\r
-\r
-//________________________________________________________________________\r
-void AliAnalysisTaskEfficiency::UserExec(Option_t *) \r
-{\r
- // Event loop\r
- // Analysis of MC true particles and reconstructed tracks\r
- // Different track types (Global, TPC, ITS) can be selected\r
-\r
- if (!fInputEvent) {\r
- Printf("ERROR: fESD not available");\r
- return;\r
- }\r
- static Int_t nfc = 0;\r
-\r
- //AliESDInputHandler* esdI = (AliESDInputHandler*) \r
- //(AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler();\r
- //AliESDEvent* hltEvent = esdI->GetHLTEvent();\r
-\r
- AliStack* stack = MCEvent()->Stack();\r
- AliESDEvent* esdE = (AliESDEvent*) fInputEvent;\r
- \r
- // CKB\r
- // Fetch the SPD vertex:\r
- Double_t vtxSPD[3];\r
- const AliESDVertex* spdVertex = esdE->GetPrimaryVertexSPD();\r
- if (spdVertex) {\r
- if(spdVertex->GetNContributors()<=0){\r
- spdVertex = 0;\r
- }\r
- else{\r
- spdVertex->GetXYZ(vtxSPD);\r
- }\r
- }\r
- \r
- // Fetch the TPC vertex\r
- Double_t vtxTPC[3];\r
- const AliESDVertex* tpcVertex = esdE->GetPrimaryVertexTPC();\r
- if (tpcVertex) {\r
- if(tpcVertex->GetNContributors()<=0){\r
- tpcVertex = 0;\r
- }\r
- else{\r
- tpcVertex->GetXYZ(vtxTPC);\r
- }\r
- }\r
- // SPD Vertex\r
- if (spdVertex) {\r
- Double_t x,y,z;\r
- x = spdVertex->GetX() + 0.07;\r
- y = spdVertex->GetY() - 0.25;\r
- z = spdVertex->GetZ();\r
- if (TMath::Abs(z) < 10.) {\r
- fVertexRvsZ->Fill(x,y);\r
- if (TMath::Sqrt(x*x + y*y) > 5. * 0.028) \r
- fVertexRvsZC->Fill(x,y);\r
- }\r
- }\r
- // BKC\r
- //Printf("%s:%d %5d",(char*)__FILE__,__LINE__, esdE->GetNumberOfTracks());\r
- TArrayI labels(esdE->GetNumberOfTracks());\r
- Int_t igood = 0;\r
- // Track loop to fill a pT spectrum\r
- //printf("ESD track loop \n");\r
-\r
- AliESDtrack *tpcP = 0x0;\r
-\r
- for (Int_t iTracks = 0; iTracks < esdE->GetNumberOfTracks(); iTracks++) {\r
-\r
- //prevent mem leak for TPConly track\r
- if(fTrackType==2&&tpcP){\r
- delete tpcP;\r
- tpcP = 0;\r
- }\r
-\r
- AliVParticle *track = esdE->GetTrack(iTracks);\r
- AliESDtrack *esdtrack = dynamic_cast<AliESDtrack*>(track);\r
- esdtrack->PropagateToDCA(esdE->GetPrimaryVertex(),\r
- esdE->GetMagneticField(), 10000.);\r
-\r
- if (!track) {\r
- Printf("ERROR: Could not receive track %d", iTracks);\r
- continue;\r
- }\r
- //__________\r
- // run Task for global tracks or ITS tracks or TPC tracks\r
- const AliExternalTrackParam *tpcPin = 0x0;\r
- //Double_t phiIn=0.;\r
- Float_t phi = 0.;\r
- Float_t eta = 0.;\r
- Float_t pt = 0.;\r
- Int_t indexAC = GetIndexAC(esdtrack); \r
-\r
- //select in the steering macro, which track type should be analysed. \r
- // Four track types are selectable:\r
- // 0. Global tracks\r
- // 1. ITS tracks (SA or Pure SA)\r
- // 2. TPC only tracks\r
- // Track selection copied from PWGPP/AliAnalysisTaskQAsym\r
-\r
- if(fTrackType==0){\r
- //Fill all histograms with global tracks\r
- tpcP = esdtrack;\r
- if (!tpcP) continue;\r
- if (!fCuts->AcceptTrack(tpcP)) continue;\r
- phi = tpcP->Phi();\r
- eta = tpcP->Eta();\r
- pt = tpcP->Pt();\r
- }\r
- else if(fTrackType==1){\r
- //Fill all histograms with ITS tracks\r
- tpcP = esdtrack;\r
- if (!tpcP) continue;\r
- if (!fCuts->AcceptTrack(tpcP)) continue;\r
- phi = tpcP->Phi();\r
- eta = tpcP->Eta();\r
- pt = tpcP->Pt();\r
- }\r
- else if(fTrackType==2){ \r
- //Fill all histograms with TPC track information\r
- tpcPin = esdtrack->GetInnerParam();\r
- if (!tpcPin) continue;\r
- tpcP = AliESDtrackCuts::\r
- GetTPCOnlyTrack(dynamic_cast<AliESDEvent*>(esdE),\r
- esdtrack->GetID());\r
- if (!tpcP) continue;\r
- if (!fCuts->AcceptTrack(tpcP)) continue;\r
- if(tpcP->GetNcls(1)>160)continue;//jacek's track cut\r
- if(tpcP->GetConstrainedChi2TPC()<0)continue; // jacek's track cut\r
- phi=tpcPin->Phi();\r
- eta=tpcPin->Eta();\r
- pt=tpcPin->Pt();\r
- }\r
- else{\r
- Printf("ERROR: wrong track type \n");\r
- continue;\r
- }\r
- //___________\r
- //\r
-\r
-\r
- labels.AddAt(-1, iTracks);\r
- \r
-\r
- // Selected\r
- if (TMath::Abs(eta) > 0.9) {\r
- } else {\r
-\r
- //Int_t charge=track->Charge();\r
- // Within acceptance\r
- Int_t ind = TMath::Abs(esdtrack->GetLabel());\r
- AliMCParticle* mcPtest = (AliMCParticle*) MCEvent()->GetTrack(ind);\r
- if (stack->IsPhysicalPrimary(mcPtest->Label())){\r
- if(TMath::Abs(mcPtest->PdgCode())>=1000020030){//all helium (+-1000020030,+-1000020040)\r
- pt*=2;// reconstruction takes charge=1 -> for particles with charge=2, pT,rec need to be corrected\r
- }\r
- \r
- Int_t index=ConvertePDG(mcPtest->PdgCode());\r
- if(fDebug>1)Printf("PDG=%d, index=%d", mcPtest->PdgCode(), index);\r
- //fill tracks comming from d,t,3He, 4He\r
- if(index<8){\r
- fHistRECptCharge [index] ->Fill(pt);\r
- fHistRECetaCharge[index]->Fill(eta);\r
- fHistRECphiCharge[index]->Fill(phi);\r
- }\r
-\r
- }\r
- fHistRECpt ->Fill(pt);\r
- fHistRECeta->Fill(eta);\r
- fHistRECphi->Fill(phi);\r
-\r
- if (pt > 2.) {\r
- fHistRECHPTeta->Fill(eta);\r
- fHistRECHPTphi->Fill(phi);\r
- }\r
- \r
- // Int_t ind = TMath::Abs(esdtrack->GetLabel());\r
- AliMCParticle* mcP = (AliMCParticle*) MCEvent()->GetTrack(ind);\r
- if (!(stack->IsPhysicalPrimary(mcP->Label()))) {\r
- fHistFAKEpt ->Fill(pt);\r
- fHistFAKEeta->Fill(eta); \r
- fHistFAKEphi->Fill(phi);\r
- \r
- if (pt > 2.) {\r
- fHistFAKEHPTeta->Fill(eta);\r
- fHistFAKEHPTphi->Fill(phi);\r
- }\r
- \r
- }\r
- \r
- labels.AddAt(TMath::Abs(esdtrack->GetLabel()), iTracks);\r
- igood++;\r
- \r
- \r
- Float_t phiDeg = TMath::RadToDeg()*track->Phi();\r
- \r
- //TPC-ITS matching\r
- //TPC-TRD matching \r
- \r
- if (tpcP){\r
- Int_t labeltpcits = esdtrack->GetLabel();\r
- AliExternalTrackParam * tpcPCopy = new AliExternalTrackParam(*tpcP);\r
- Double_t xk = 43.6; // cm\r
- Double_t bz = 5.0; // kG\r
- if(tpcPCopy->PropagateTo(xk,bz)){\r
- Double_t alpha=tpcPCopy->GetAlpha();\r
- // if(tpcPCopy->Rotate(0.)){ \r
- Float_t xtpc=tpcPCopy->GetX();\r
- Float_t ytpc=tpcPCopy->GetY();\r
- Float_t ztpc=tpcPCopy->GetZ();\r
- Float_t xpr,ypr ; \r
- xpr=xtpc*TMath::Cos(alpha)-ytpc*TMath::Sin(alpha);\r
- ypr=xtpc*TMath::Sin(alpha)+ytpc*TMath::Cos(alpha);\r
- \r
- AliMCParticle* mcPart = (AliMCParticle*) MCEvent()->GetTrack(abs(labeltpcits));\r
- Int_t ntref = mcPart->GetNumberOfTrackReferences();\r
- \r
- for (Int_t k = 0; k < ntref; k++) {\r
- AliTrackReference* ref = mcPart->GetTrackReference(k);\r
- Float_t xhits = ref->X();\r
- Float_t yhits = ref->Y();\r
- Float_t radio = TMath::Sqrt(xhits*xhits+yhits*yhits);\r
- \r
- if(ref->DetectorId() == 0 && (radio > 42.8 && radio < 43.7)) {\r
- \r
- Float_t xits = ref->X();\r
- Float_t yits = ref->Y();\r
- Float_t zits = ref->Z();\r
- \r
- Float_t difx=(xits-xpr);\r
- fHistTPCITSdeltax->Fill(difx);\r
- Float_t dify=(yits-ypr);\r
- fHistTPCITSdeltay->Fill(dify);\r
- Float_t difz=(zits-ztpc);\r
- fHistTPCITSdeltaz->Fill(difz);\r
- Float_t deltar = TMath::Sqrt(difx * difx + dify * dify + difz * difz);\r
- fHistTPCITSdeltar->Fill(deltar);\r
- break;\r
- }\r
- // }\r
- } // trackrefs\r
- } \r
- } //ITS-TPC maching\r
- //TPC-TRD maching \r
- \r
- const AliExternalTrackParam *trd = esdtrack->GetOuterParam();\r
- \r
- if (trd){Int_t labeltpctrd = track->GetLabel();\r
- \r
- AliExternalTrackParam * trdCopy = new AliExternalTrackParam(*trd);\r
- Double_t xktrd=296.0;\r
- Double_t bztrd=5.0;\r
- if(trdCopy->PropagateTo(xktrd,bztrd)){\r
- Float_t xtpc2=trdCopy->GetX();\r
- Float_t ytpc2=trdCopy->GetY();\r
- Float_t ztpc2=trdCopy->GetZ();\r
- Double_t alpha=trdCopy->GetAlpha();\r
- Float_t xpr,ypr ; \r
- xpr=xtpc2*TMath::Cos(alpha)-ytpc2*TMath::Sin(alpha);\r
- ypr=xtpc2*TMath::Sin(alpha)+ytpc2*TMath::Cos(alpha);\r
- AliMCParticle* mcPart = (AliMCParticle*) MCEvent()->GetTrack(abs(labeltpctrd));\r
- Int_t ntreftrd = mcPart->GetNumberOfTrackReferences(); \r
- \r
- for (Int_t k = 0; k < ntreftrd; k++) {\r
- AliTrackReference* ref = mcPart->GetTrackReference(k);\r
- if(ref->DetectorId() == 3 && ref->Label() == labeltpctrd){\r
- \r
- Float_t xtrd = ref->X();\r
- Float_t ytrd = ref->Y();\r
- Float_t ztrd = ref->Z();\r
- \r
- Float_t difx=(xtrd-xpr);\r
- fHistTPCTRDdeltax->Fill(difx);\r
- Float_t dify=(ytrd-ypr);\r
- fHistTPCTRDdeltay->Fill(dify);\r
- Float_t difz=(ztrd-ztpc2);\r
- fHistTPCTRDdeltaz->Fill(difz);\r
- Float_t deltar = TMath::Sqrt(difx * difx + dify * dify + difz * difz);\r
- fHistTPCTRDdeltar->Fill(deltar);\r
- Float_t phi_tpc = TMath::ATan2(ypr, xpr);\r
- Float_t phi_trd = TMath::ATan2(ytrd, xtrd);\r
- Float_t dphi = phi_trd - phi_tpc;\r
- if (dphi > TMath::Pi()) dphi -= 2. * TMath::Pi();\r
- if (dphi < - TMath::Pi()) dphi += 2. * TMath::Pi();\r
-\r
- fHistTPCTRDdeltaphi->Fill(dphi);\r
- break;\r
- }\r
- \r
- }\r
- }\r
-\r
- }//TRD-TPC maching\r
-\r
- // CKB\r
- if(pt>2.&&fFieldOn){\r
- TBits bits(esdtrack->GetTPCClusterMap());\r
- for(unsigned int ip = 0;ip < bits.GetNbits();++ip){\r
- if(bits[ip]){\r
- fh2PhiPadRow[indexAC]->Fill(phiDeg,ip);\r
- }\r
- }\r
- for(int i = 0;i < 6;++i){\r
- if(esdtrack->HasPointOnITSLayer(i))fh2PhiLayer[indexAC]->Fill(phiDeg,i);\r
- }\r
- }\r
- else if(!fFieldOn){ // field not on all momenta are set to MPV \r
- TBits bits(esdtrack->GetTPCClusterMap());\r
- for(unsigned int ip = 0;ip < bits.GetNbits();++ip){\r
- if(bits[ip]){\r
- fh2PhiPadRow[indexAC]->Fill(phiDeg,ip);\r
- }\r
- for(int i = 0;i < 6;++i){\r
- if(esdtrack->HasPointOnITSLayer(i))fh2PhiLayer[indexAC]->Fill(phiDeg,i);\r
- }\r
- } \r
-\r
- }\r
- // Fill the correlation\r
- if(MCEvent()){\r
- TParticle *part = MCEvent()->Stack()->Particle(TMath::Abs(esdtrack->GetLabel()));\r
- if(part){\r
- Float_t mcEta = part->Eta();\r
- Float_t mcPhi = TMath::RadToDeg()*part->Phi();\r
- Float_t mcPt = part->Pt();\r
- // if (pt - mcPt > 20.) {\r
- fh2EtaCorrelation->Fill(mcEta,eta);\r
- fh2EtaCorrelationShift->Fill(eta,mcEta-eta);\r
- fh2PhiCorrelation->Fill(mcPhi,phiDeg);\r
- fh2PhiCorrelationShift->Fill(phiDeg,mcPhi-phiDeg);\r
- fh2PtCorrelation->Fill(mcPt,pt);\r
- fh2PtCorrelationShift->Fill(pt,mcPt-pt);\r
- //}\r
- Double_t sigmaPt = TMath::Sqrt(esdtrack->GetSigma1Pt2());\r
- if (track->Charge() > 0.) {\r
- fPtPullsPos->Fill(mcPt, (1./pt - 1./mcPt) / sigmaPt); \r
- fPtShiftPos->Fill(mcPt, (1./pt - 1./mcPt));\r
- } else {\r
- fPtPullsNeg->Fill(mcPt, (1./pt - 1./mcPt) / sigmaPt); \r
- fPtShiftNeg->Fill(mcPt, (1./pt - 1./mcPt)); \r
- }\r
- }\r
- // BKC\r
- }\r
- }\r
- \r
- } //track loop \r
-\r
- //prevent mem leak for TPConly track\r
- if(fTrackType==2&&tpcP){\r
- delete tpcP;\r
- tpcP = 0;\r
- }\r
-\r
- //Printf("%s:%d",(char*)__FILE__,__LINE__);\r
- // CKB\r
- // Vertex eff\r
- if(tpcVertex){\r
- fh1VtxEff->Fill("TPC VTX",1);\r
- }\r
- else{\r
- fh1VtxEff->Fill("NO TPC VTX",1);\r
- }\r
-\r
- if(spdVertex){\r
- fh1VtxEff->Fill("SPD VTX",1);\r
- }\r
- else{\r
- fh1VtxEff->Fill("NO SPD VTX",1);\r
- }\r
-\r
-\r
- // Vertex correlation SPD vs. TPC\r
- if(tpcVertex&&spdVertex){\r
- fh2VtxTpcSpd->Fill(vtxTPC[2],vtxSPD[2]);\r
- }\r
- // Printf("%s:%d",(char*)__FILE__,__LINE__);\r
- // Multiplicity checks in the SPD\r
- const AliMultiplicity *spdMult = esdE->GetMultiplicity();\r
- // Multiplicity Analysis\r
- Int_t nt = spdMult->GetNumberOfTracklets();\r
- nfc++;\r
- if (nfc == 520) {\r
- nfc = 0;\r
- for (Int_t ib = 1; ib <= 10; ib++) {\r
- Int_t ic = Int_t(fEtaMulti->GetBinContent(ib));\r
- Float_t eta = -1.8 + Float_t(ib-1) * 0.4;\r
- fEtaMultiFluc->Fill(eta, Float_t(ic));\r
- if (ib == 5) fEtaMultiH->Fill(Float_t(ic));\r
- }\r
-\r
- for (Int_t ib = 1; ib <= 64; ib++) {\r
- Int_t ic = Int_t(fPhiMulti->GetBinContent(ib));\r
- Float_t phi = 0.05 + Float_t(ib-1) * 0.1;\r
- fPhiMultiFluc->Fill(phi, Float_t(ic));\r
- if (ib == 2) fPhiMultiH->Fill(Float_t(ic));\r
- }\r
-\r
- fEtaMulti->Reset();\r
- fPhiMulti->Reset(); \r
- }\r
- \r
- for (Int_t j = 0; j < nt; j++) {\r
- fEtaMulti->Fill(spdMult->GetEta(j));\r
- fPhiMulti->Fill(spdMult->GetPhi(j));\r
- }\r
-\r
- Short_t chips0 = spdMult->GetNumberOfFiredChips(0);\r
- Short_t chips1 = spdMult->GetNumberOfFiredChips(1);\r
- //Printf("%s:%d",(char*)__FILE__,__LINE__);\r
- Int_t inputCount = 0;\r
- \r
- if(spdVertex){\r
- // get multiplicity from ITS tracklets\r
- for (Int_t i=0; i<spdMult->GetNumberOfTracklets(); ++i){\r
- Float_t deltaPhi = spdMult->GetDeltaPhi(i);\r
- // prevent values to be shifted by 2 Pi()\r
- if (deltaPhi < -TMath::Pi())\r
- deltaPhi += TMath::Pi() * 2;\r
- if (deltaPhi > TMath::Pi())\r
- deltaPhi -= TMath::Pi() * 2; \r
- if (TMath::Abs(deltaPhi) > 1)\r
- printf("WARNING: Very high Delta Phi: %d %f %f %f\n", i, spdMult->GetTheta(i), spdMult->GetPhi(i), deltaPhi);\r
- //trackled Deltaphi\r
- Int_t label1=spdMult->GetLabel(i,0);\r
- Int_t label2=spdMult->GetLabel(i,1);\r
- if ( label1==label2){\r
- if(stack->IsPhysicalPrimary(label1) == kTRUE)\r
- fHistDeltaphiprimaries->Fill(deltaPhi);\r
- else\r
- fHistDeltaphisecondaries->Fill(deltaPhi);\r
- }\r
- else{\r
- fHistDeltaphireject->Fill(deltaPhi);\r
- }\r
- ++inputCount;\r
- }\r
- // Printf("%s:%d",(char*)__FILE__,__LINE__);\r
- fh2MultSpdChips[0]->Fill(chips0,inputCount);\r
- fh2MultSpdChips[1]->Fill(chips1,inputCount);\r
- //Printf("%s:%d",(char*)__FILE__,__LINE__);\r
- }\r
- // BKC\r
-\r
- // MC analysis\r
- Int_t nmc = MCEvent()->GetNumberOfTracks();\r
- AliGenEventHeader* header = MCEvent()->GenEventHeader();\r
- // some test\r
- AliGenCocktailEventHeader* cH = dynamic_cast<AliGenCocktailEventHeader*> (header);\r
- AliGenPythiaEventHeader* pH;\r
- if (cH == 0) \r
- {\r
- header->Dump();\r
- pH = dynamic_cast<AliGenPythiaEventHeader*>(header);\r
- } else {\r
- TObject* entry = cH->GetHeaders()->FindObject("Pythia");\r
- pH = dynamic_cast<AliGenPythiaEventHeader*>(entry);\r
- }\r
- Int_t iproc = -1;\r
-\r
- if (pH) iproc = pH->ProcessType();\r
-\r
- TArrayF mcV;\r
- header->PrimaryVertex(mcV);\r
- Float_t vzMC = mcV[2];\r
-\r
- // CKB\r
- // Fill the correlation with MC\r
- if(spdVertex&&MCEvent()){\r
- Float_t diff = vzMC-vtxSPD[2];\r
- fh2VertexCorrelation[kSPD]->Fill(vzMC,vtxSPD[2]);\r
- fh2VertexCorrelationShift[kSPD]->Fill(vzMC,diff);\r
- fh1VertexShift[kSPD]->Fill(diff);\r
- if(spdVertex->GetZRes()>0)fh1VertexShiftNorm[kSPD]->Fill(diff/spdVertex->GetZRes());\r
- }\r
- if(tpcVertex&&MCEvent()){\r
- Float_t diff = vzMC-vtxTPC[2];\r
- fh2VertexCorrelation[kTPC]->Fill(vzMC,vtxTPC[2]);\r
- fh2VertexCorrelationShift[kTPC]->Fill(vzMC,diff);\r
- fh1VertexShift[kTPC]->Fill(diff);\r
- if(tpcVertex->GetZRes()>0)fh1VertexShiftNorm[kTPC]->Fill(diff/tpcVertex->GetZRes());\r
- }\r
- // BKC\r
-\r
-\r
- // MC event loop\r
- //printf("MC particle loop \n");\r
- for (Int_t i = 0; i < nmc; i++) {\r
- AliMCParticle* mcP = (AliMCParticle*) MCEvent()->GetTrack(i);\r
- //printf("MC particle loop %5d \n", i);\r
- // Primaries only\r
- if (!(stack->IsPhysicalPrimary(mcP->Label()))) continue;\r
- if (mcP->Particle()->GetPDG()->Charge() == 0) continue;\r
- // Int_t charge=Int_t(mcP->Particle()->GetPDG()->Charge());\r
- Float_t phi = mcP->Phi();\r
- Float_t eta = mcP->Eta();\r
- Float_t pt = mcP->Pt();\r
- if (TMath::Abs(eta) > 0.9) continue;\r
- Int_t ntr = mcP->GetNumberOfTrackReferences();\r
- Int_t nITS = 0;\r
- Int_t nTPC = 0;\r
- Int_t nFRA = 0;\r
- Float_t x = 0.;\r
- Float_t y = 0.;\r
- \r
- for (Int_t j = 0; j < ntr; j++) {\r
- \r
- AliTrackReference* ref = mcP->GetTrackReference(j);\r
- \r
- if (ref->DetectorId() == 0) nITS++;\r
- if (ref->DetectorId() == 1) nTPC++;\r
- if (ref->DetectorId() == 2) nFRA++;\r
- if (nTPC == 1) {\r
- x = ref->X();\r
- y = ref->Y();\r
- break;\r
- }\r
- }\r
-\r
- fHistMCpt ->Fill(pt);\r
- fHistMCeta->Fill(eta);\r
- fHistMCphi->Fill(phi);\r
- Int_t index=ConvertePDG(mcP->PdgCode());\r
- if(index<8){\r
- fHistMCptCharge [index] ->Fill(pt);\r
- fHistMCetaCharge[index]->Fill(eta);\r
- fHistMCphiCharge[index]->Fill(phi);\r
- }\r
-\r
-\r
- if (pt > 2.) {\r
- fHistMCHPTeta->Fill(eta);\r
- fHistMCHPTphi->Fill(phi);\r
- }\r
-\r
- Bool_t reco = kFALSE;\r
- Int_t multR = 0;\r
- Int_t jold = -1;\r
- \r
- for (Int_t j = 0; j < esdE->GetNumberOfTracks(); j++) {\r
- if (i == labels.At(j)) {\r
- reco = kTRUE;\r
- multR++;\r
- //AliESDtrack* test = esdE->GetTrack(j);\r
- if (multR > 1) {\r
- Int_t nclus = 0;\r
- AliESDtrack* track = esdE->GetTrack(jold);\r
- nclus = track->GetTPCclusters(0); \r
- fHistNCluster->Fill(nclus);\r
- \r
- \r
- track = esdE->GetTrack(j);\r
- nclus = track->GetTPCclusters(0);\r
- fHistNCluster->Fill(nclus);\r
- \r
-\r
- }\r
- jold = j;\r
- }\r
- }\r
-\r
- fHistRecMult->Fill(Float_t(multR), 1.);\r
- if (multR == 0) {\r
- fHistMCNRpt ->Fill(pt);\r
- fHistMCNReta->Fill(eta);\r
- fHistMCNRphi->Fill(phi);\r
- if(index<8){\r
- fHistMCNRptCharge [index] ->Fill(pt);\r
- fHistMCNRetaCharge[index]->Fill(eta);\r
- fHistMCNRphiCharge[index]->Fill(phi);\r
- }\r
- \r
-\r
- if (pt > 2.) {\r
- fHistMCNRHPTeta->Fill(eta);\r
- fHistMCNRHPTphi->Fill(phi);\r
- }\r
- }\r
- else if (multR > 1)\r
- {\r
- fHistMULTpt ->Fill(pt);\r
- fHistMULTeta->Fill(eta);\r
- fHistMULTphi->Fill(phi);\r
- } \r
- } // MC particle loop\r
- //printf("End of MC particle loop \n");\r
- \r
- \r
-} \r
-\r
-//________________________________________________________________________\r
-void AliAnalysisTaskEfficiency::Terminate(Option_t *) \r
-{\r
-\r
-} \r
-\r
-\r
-Bool_t AliAnalysisTaskEfficiency::SelectJouri(AliESDtrack* track) \r
-{\r
-\r
- Bool_t selected = kTRUE;\r
- AliESDEvent* esdE = (AliESDEvent*) fInputEvent; \r
- // > 50 Clusters\r
- if (track->GetTPCclusters(0) < 50) selected = kFALSE;\r
-\r
- const AliESDVertex *vtx=esdE->GetPrimaryVertexSPD();\r
- if (!vtx->GetStatus()) selected = kFALSE;\r
- \r
- Double_t zv=vtx->GetZv();\r
-\r
- \r
- const AliExternalTrackParam *ct=track->GetTPCInnerParam();\r
- if (!ct) selected = kFALSE;\r
- \r
- Double_t xyz[3];\r
- ct->GetXYZ(xyz);\r
- Float_t rv = TMath::Sqrt(xyz[0] * xyz[0] + xyz[1] * xyz[1]);\r
- if (rv > 3.0) selected = kFALSE;\r
- if (TMath::Abs(xyz[2] - zv)>2.5) selected = kFALSE;\r
- return (selected);\r
- \r
-}\r
-\r
-Int_t AliAnalysisTaskEfficiency::GetIndexAC(AliESDtrack *track){\r
- if(!track)return -1;\r
-\r
- // Crude selection for A and C side\r
- // just with eta\r
- if (track->Eta() > 0) { \r
- if (track->Charge() > 0) \r
- return kPosA;\r
- else\r
- return kNegA; \r
- }\r
- else { // C side\r
- if (track->Charge() > 0) \r
- return kPosC;\r
- else\r
- return kNegC; \r
- }\r
- \r
- return -1;\r
-\r
-}\r
-\r
-\r
-Int_t AliAnalysisTaskEfficiency::ConvertePDG(Int_t pdg)\r
-{\r
-\r
- // function converts the pdg values for nuclei (d, t, 3He, \r
- // 4He and their antiparticles) into indices for histo-array\r
- // filled in UserExec\r
- \r
- if ( pdg == 1000010020 )return 0;\r
- else if ( pdg == -1000010020 )return 1;\r
- else if ( pdg == 1000010030 )return 2;\r
- else if ( pdg == -1000010030 )return 3;\r
- else if ( pdg == 1000020030 )return 4;\r
- else if ( pdg == -1000020030 )return 5;\r
- else if ( pdg == 1000020040 )return 6;\r
- else if ( pdg == -1000020040 )return 7;\r
- else return 9;\r
- \r
- \r
-}\r
+#include "TChain.h"
+#include "TTree.h"
+#include "TH1F.h"
+#include "TH2F.h"
+#include "TProfile.h"
+
+#include "TCanvas.h"
+#include "TList.h"
+#include "TParticle.h"
+#include "TParticlePDG.h"
+#include "TProfile.h"
+
+#include "AliAnalysisTask.h"
+#include "AliAnalysisManager.h"
+
+#include "AliESDEvent.h"
+#include "AliStack.h"
+#include "AliESDVertex.h"
+#include "AliESDInputHandler.h"
+#include "AliESDtrackCuts.h"
+#include "AliMultiplicity.h"
+
+#include "AliMCParticle.h"
+#include "AliMCEvent.h"
+#include "AliAnalysisTaskEfficiency.h"
+#include "AliExternalTrackParam.h"
+#include "AliTrackReference.h"
+#include "AliStack.h"
+#include "AliHeader.h"
+#include "AliGenEventHeader.h"
+#include "AliGenCocktailEventHeader.h"
+#include "AliGenPythiaEventHeader.h"
+#include "AliGenDPMjetEventHeader.h"
+
+// Analysis Task for basic QA on the ESD
+// Authors: Jan Fiete Grosse-Oetringhaus, Christian Klein-Boesing, Veronica Canoa, AM, Eva Sicking
+
+ClassImp(AliAnalysisTaskEfficiency)
+
+//________________________________________________________________________
+ AliAnalysisTaskEfficiency::AliAnalysisTaskEfficiency(const char *name)
+ : AliAnalysisTaskSE(name)
+ ,fFieldOn(kTRUE)
+ ,fHists(0)
+ ,fHistRECpt(0)
+ ,fHistMCpt(0)
+ ,fHistMCNRpt(0)
+ ,fHistFAKEpt(0)
+ ,fHistMULTpt(0)
+ ,fHistRECeta(0)
+ ,fHistMCeta(0)
+ ,fHistMCNReta(0)
+ ,fHistFAKEeta(0)
+ ,fHistMULTeta(0)
+ ,fHistRECphi(0)
+ ,fHistMCphi(0)
+ ,fHistMCNRphi(0)
+ ,fHistFAKEphi(0)
+ ,fHistMULTphi(0)
+ ,fHistRECHPTeta(0)
+ ,fHistMCHPTeta(0)
+ ,fHistMCNRHPTeta(0)
+ ,fHistFAKEHPTeta(0)
+ ,fHistRECHPTphi(0)
+ ,fHistMCHPTphi(0)
+ ,fHistMCNRHPTphi(0)
+ ,fHistFAKEHPTphi(0)
+ ,fHistRecMult(0)
+ ,fHistNCluster(0)
+ ,fh1VtxEff(0)
+ ,fh2VtxTpcSpd(0)
+ ,fh2EtaCorrelation(0)
+ ,fh2EtaCorrelationShift(0)
+ ,fh2PhiCorrelation(0)
+ ,fh2PhiCorrelationShift(0)
+ ,fh2PtCorrelation(0)
+ ,fh2PtCorrelationShift(0)
+ ,fHistDeltaphiprimaries(0)
+ ,fHistDeltaphisecondaries(0)
+ ,fHistDeltaphireject(0)
+ ,fHistTPCITSdeltax(0)
+ ,fHistTPCITSdeltay(0)
+ ,fHistTPCITSdeltaz(0)
+ ,fHistTPCITSdeltar(0)
+ ,fHistTPCTRDdeltax(0)
+ ,fHistTPCTRDdeltay(0)
+ ,fHistTPCTRDdeltaz(0)
+ ,fHistTPCTRDdeltar(0)
+ ,fHistTPCTRDdeltaphi(0)
+ ,fPtPullsPos(0)
+ ,fPtPullsNeg(0)
+ ,fPtShiftPos(0)
+ ,fPtShiftNeg(0)
+ ,fTrackType(0)
+ ,fCuts(0)
+ ,fVertexRvsZ(0)
+ ,fVertexRvsZC(0)
+ ,fEtaMultiFluc(0)
+ ,fEtaMulti(0)
+ ,fEtaMultiH(0)
+ ,fPhiMultiFluc(0)
+ ,fPhiMulti(0)
+ ,fPhiMultiH(0)
+{
+ //
+ // Constructor
+ //
+
+ for(int i = 0;i< 2;++i){
+ fh2MultSpdChips[i] = 0;
+ }
+
+ for(int i = 0;i < kTotalAC;++i){
+ fh2PhiPadRow[i] = fh2PhiLayer[i] = 0;
+ }
+
+ for(int i = 0; i < kTotalVtx; ++i){
+ fh2VertexCorrelation[i] =
+ fh2VertexCorrelationShift[i] = 0;
+ fh1VertexShift[i] =
+ fh1VertexShiftNorm[i] = 0;
+ }
+
+ for(int i = 0;i< 8;++i){
+ fHistRECptCharge[i]=0;
+ fHistMCptCharge[i]=0;
+ fHistMCNRptCharge[i]=0;
+ fHistFAKEptCharge[i]=0;
+
+ fHistRECetaCharge[i]=0;
+ fHistMCetaCharge[i]=0;
+ fHistMCNRetaCharge[i]=0;
+ fHistFAKEetaCharge[i]=0;
+
+ fHistRECphiCharge[i]=0;
+ fHistMCphiCharge[i]=0;
+ fHistMCNRphiCharge[i]=0;
+ fHistFAKEphiCharge[i]=0;
+
+ }
+
+ DefineOutput(1, TList::Class());
+}
+
+
+//________________________________________________________________________
+void AliAnalysisTaskEfficiency::UserCreateOutputObjects()
+{
+ // Create histograms
+
+
+ fHists = new TList();
+ fHistRECpt = new TH1F("fHistRECpt", " p_{T} distribution: all reconstructed", 100, 0., 20.);
+ fHistMCpt = new TH1F("fHistMCpt", " p_{T} distribution: all MC", 100, 0., 20.);
+ fHistMCNRpt = new TH1F("fHistMCNRpt", " p_{T} distribution: all not-reconstructed MC", 100, 0., 20.);
+ fHistFAKEpt = new TH1F("fHistFAKEpt", " p_{T} distribution: all Fake", 100, 0., 20.);
+ fHistMULTpt = new TH1F("fHistMULTpt", " p_{T} distribution: multiply rec.", 100, 0., 20.);
+
+ fHistRECeta = new TH1F("fHistRECeta", " #eta distribution: all reconstructed", 100,-1.0,1.0);
+ fHistMCeta = new TH1F("fHistMCeta", " #eta distribution: all MC", 100,-1.0,1.0);
+ fHistMCNReta = new TH1F("fHistMCNReta", " #eta distribution: all not-reconstructed MC", 100,-1.0,1.0);
+ fHistFAKEeta = new TH1F("fHistFAKEeta", " #eta distribution: all Fake", 100,-1.0,1.0);
+ fHistMULTeta = new TH1F("fHistMULTeta", " #eta distribution: multiply rec.", 100,-1.0,1.0);
+
+ fHistRECphi = new TH1F("fHistRECphi", " #phi distribution: all reconstructed", 314, 0., 6.28);
+ fHistMCphi = new TH1F("fHistMCphi", " #phi distribution: all MC", 314, 0., 6.28);
+ fHistMCNRphi = new TH1F("fHistMCNRphi", " #phi distribution: all not-reconstructed MC", 314, 0., 6.28);
+ fHistFAKEphi = new TH1F("fHistFAKEphi", " #phi distribution: all Fake", 314, 0., 6.28);
+ fHistMULTphi = new TH1F("fHistMULTphi", " #phi distribution: multipli rec.", 314, 0., 6.28);
+
+ fHistRECHPTeta = new TH1F("fHistRECHPTeta", " #eta distribution: all reconstructed", 100,-1.0,1.0);
+ fHistMCHPTeta = new TH1F("fHistMCHPTeta", " #eta distribution: all MC", 100,-1.0,1.0);
+ fHistMCNRHPTeta = new TH1F("fHistMCNRHPTeta", " #eta distribution: all not-reconstructed MC", 100,-1.0,1.0);
+ fHistFAKEHPTeta = new TH1F("fHistFAKEHPTeta", " #eta distribution: all Fake", 100,-1.0,1.0);
+
+ fHistRECHPTphi = new TH1F("fHistRECHPTphi", " #phi distribution: all reconstructed", 314, 0., 6.28);
+ fHistMCHPTphi = new TH1F("fHistMCHPTphi", " #phi distribution: all MC", 314, 0., 6.28);
+ fHistMCNRHPTphi = new TH1F("fHistMCNRHPTphi", " #phi distribution: all not-reconstructed MC", 314, 0., 6.28);
+ fHistFAKEHPTphi = new TH1F("fHistFAKEHPTphi", " #phi distribution: all Fake", 314, 0., 6.28);
+
+ fHistRecMult = new TH1F("fHistRecMult", "Multiple reconstructed tracks", 50, 0., 50.);
+ fHistNCluster = new TH1F("fHistNCluster", "Number of clusters for suspicious tracks", 300, 0., 300.);
+
+ // CKB
+
+ fh1VtxEff = new TH1F("fh1VtxEff","VtxEff TPC",4,-0.5,3.5);
+ fh1VtxEff->GetXaxis()->SetBinLabel(1,"NO TPC VTX");
+ fh1VtxEff->GetXaxis()->SetBinLabel(2,"TPC VTX");
+ fh1VtxEff->GetXaxis()->SetBinLabel(3,"NO SPD VTX");
+ fh1VtxEff->GetXaxis()->SetBinLabel(4,"SPD VTX");
+
+ TString labels[kTotalAC];
+ labels[kNegA]="NegA";
+ labels[kPosA]="PosA";
+ labels[kNegC]="NegC";
+ labels[kPosC]="PosC";
+
+ for(int i = 0;i< kTotalAC;++i){
+ fh2PhiPadRow[i] = new TH2F(Form("fh2PhiPadRow_%d",i),Form("Padrow vs phi AC %s;phi;padrow",labels[i].Data()),360,0.,360.,159,-0.5,158.5);
+ fh2PhiLayer[i] = new TH2F(Form("fh2PhiLayer_%d",i),Form("layer vs phi AC %s;phi;layer",labels[i].Data()),360,0.,360.,6,-0.5,5.5);
+ }
+ for(int i = 0;i<2;++i){
+ fh2MultSpdChips[i] = new TH2F(Form("fh2ChipsSpdMult_%d",i),"mult SPD vs chips;chips;Mult",201,-1.5,199.5,201,-1.5,199.5);
+ }
+ fh2VtxTpcSpd = new TH2F("fh2VtxTpcSpd","SPD vtx vs TPC;tpc-z;spd-z",200,-20.,20.,200,-20,20.);
+
+ TString labelsVtx[kTotalVtx];
+ labelsVtx[kTPC] = "TPC";
+ labelsVtx[kSPD] = "SPD";
+ for(int i = 0; i < kTotalVtx; ++i){
+ fh2VertexCorrelation[i] = new TH2F(Form("fh2VertexCorrelation_%d",i),Form("VertexCorrelation %s;MC z-vtx;ESD z-vtx",labelsVtx[i].Data()), 120, -30, 30, 120, -30, 30);
+ fh2VertexCorrelationShift[i] = new TH2F(Form("fh2VertexCorrelationShift_%d",i), Form("VertexCorrelationShift %s;MC z-vtx;MC z-vtx - ESD z-vtx",labelsVtx[i].Data()), 120, -30, 30, 100, -1, 1);
+ fh1VertexShift[i] = new TH1F(Form("fh1VertexShift_%d",i), Form("VertexShift %s;(MC z-vtx - ESD z-vtx);Entries",labelsVtx[i].Data()), 201, -2, 2);
+ fh1VertexShiftNorm[i] = new TH1F(Form("fh1VertexShiftNorm_%d",i), Form("VertexShiftNorm %s;(MC z-vtx - ESD z-vtx) / #sigma_{ESD z-vtx};Entries",labelsVtx[i].Data()), 200, -100, 100);
+ }
+
+
+ fh2EtaCorrelation = new TH2F("fh2EtaCorrelation", "EtaCorrelation;MC #eta;ESD #eta", 120, -3, 3, 120, -3, 3);
+ fh2EtaCorrelationShift = new TH2F("fh2EtaCorrelationShift", "EtaCorrelationShift;ESD #eta;MC #eta - ESD #eta", 120, -3, 3, 100, -0.1, 0.1);
+ fh2PhiCorrelation = new TH2F("fh2PhiCorrelation", "PhiCorrelation;MC #phi;ESD #phi", 120, 0, 360, 120, 0, 360);
+ fh2PhiCorrelationShift = new TH2F("fh2PhiCorrelationShift", "PhiCorrelationShift;ESD #phi;MC #phi - ESD #phi", 120, 0, 360, 100, -10, 10);
+ fh2PtCorrelation = new TH2F("fh2PtCorrelation", "PtCorrelation;MC p_T;ESD p_T", 200, 0., 200., 200, 0, 200);
+ fh2PtCorrelationShift = new TH2F("fh2PtCorrelationShift", "PtCorrelationShift;ESD p_T;MC p_T - ESD #p_T", 120, 0, 10, 100, -2, 2);
+
+ fHistDeltaphiprimaries =new TH1F("fHistDeltaphiprimaries", " #Delta#phi distribution: primaries", 314, -0.2,0.2);
+ fHistDeltaphisecondaries=new TH1F("fHistDeltaphisecondaries", "#Delta#phi distribution: secondaries", 314, -0.2,0.2);
+ fHistDeltaphireject =new TH1F("fHistDeltaphireject", " #Delta#phi distribution: reject trackled", 314, -0.2,0.2);
+
+ fHistTPCITSdeltax =new TH1F("fHistTPCITSdeltax", "TPC-ITS matching:#Delta x", 100,-20,20);
+ fHistTPCITSdeltay =new TH1F("fHistTPCITSdeltay", "TPC-ITS matching:#Delta y", 100,-20,20);
+ fHistTPCITSdeltaz =new TH1F("fHistTPCITSdeltaz", "TPC-ITS matching:#Delta z", 100,-20,20);
+ fHistTPCITSdeltar =new TH1F("fHistTPCITSdeltar", "TPC-ITS matching:#Delta r", 100,0,20);
+
+ fHistTPCTRDdeltax = new TH1F("fHistTPCTRDdeltax", "TPC-TRD matching:#Delta x", 400,-400, 400);
+ fHistTPCTRDdeltay = new TH1F("fHistTPCTRDdeltay", "TPC-TRD matching:#Delta y", 400,-400, 400);
+ fHistTPCTRDdeltaz = new TH1F("fHistTPCTRDdeltaz", "TPC-TRD matching:#Delta z", 400,-400, 400);
+ fHistTPCTRDdeltar = new TH1F("fHistTPCTRDdeltar", "TPC-TRD matching:#Delta r", 100,0,20);
+ fHistTPCTRDdeltaphi = new TH1F("fHistTPCTRDdeltarphi", "TPC-TRD matching:#Delta #phi", 128,-3.14 , 3.14);
+ // Pulls
+ fPtPullsPos = new TProfile("fPtPullsPos", "Pt Pulls for pos. primaries", 50 , 0., 10., -5., 5.,"S");
+ fPtPullsNeg = new TProfile("fPtPullsNeg", "Pt Pulls for neg. primaries", 50 , 0., 10., -5., 5.,"S");
+ fPtShiftPos = new TProfile("fPtShiftPos", "Pt Shift for pos. primaries", 50 , 0., 10., -0.2, 0.2,"S");
+ fPtShiftNeg = new TProfile("fPtShiftNeg", "Pt Shift for neg. primaries", 50 , 0., 10., -0.2, 0.2,"S");
+
+ fEtaMultiFluc = new TProfile("fEtaMultiFluc", "eta multiplicity fluctuations", 10, -2., 2., 0., 600., "S");
+ fEtaMulti = new TH1F("fEtaMulti", "eta multiplicity fluctuations", 10, -2., 2.);
+ fEtaMultiH = new TH1F("fEtaMultiH", "eta multiplicity fluctuations", 600, 0., 600.);
+
+ fPhiMultiFluc = new TProfile("fPhiMultiFluc", "phi multiplicity fluctuations", 64, 0., 6.4, 0., 100., "S");
+ fPhiMulti = new TH1F("fPhiMulti", "phi multiplicity fluctuations", 64, 0., 6.4);
+ fPhiMultiH = new TH1F("fPhiMultiH", "phi multiplicity fluctuations", 100, 0., 100.);
+
+ // SPD Vertex
+ fVertexRvsZ = new TH2F("fVertexRvsZ", "SPD Vertex R vs Z", 200, -1., 1., 200, -1., 1.);
+ fVertexRvsZC = new TH2F("fVertexRvsZC", "SPD Vertex R vs Z", 200, -1., 1., 200, -1., 1.);
+
+ TString charge[8]={"Deu", "AntiDeu", "Tri", "AntiTri", "He3", "AntiHe3", "He4", "AntiHe4"};
+ for(Int_t i=0;i<8;i++){
+ fHistRECptCharge[i] = new TH1F(Form("fHistRECptCharge%s",charge[i].Data()),
+ "p_{T} distribution: all reconstructed",
+ 100, 0., 20.);
+ fHistMCptCharge[i] = new TH1F(Form("fHistMCptCharge%s",charge[i].Data()),
+ "p_{T} distribution: all MC",
+ 100, 0., 20.);
+ fHistMCNRptCharge[i] = new TH1F(Form("fHistMCNRptCharge%s",charge[i].Data()),
+ "p_{T} distribution: all not-reconstructed MC",
+ 100, 0., 20.);
+ fHistFAKEptCharge[i] = new TH1F( Form("fHistFAKEptCharge%s",charge[i].Data()),
+ "p_{T} distribution: all Fake",
+ 100, 0., 20.);
+
+ fHistRECetaCharge[i] = new TH1F(Form("fHistRECetaCharge%s",charge[i].Data()),
+ "p_{T} distribution: all reconstructed",
+ 100, 0., 20.);
+ fHistMCetaCharge[i] = new TH1F(Form("fHistMCetaCharge%s",charge[i].Data()),
+ "p_{T} distribution: all MC",
+ 100, 0., 20.);
+ fHistMCNRetaCharge[i] = new TH1F(Form("fHistMCNRetaCharge%s",charge[i].Data()),
+ "p_{T} distribution: all not-reconstructed MC",
+ 100, 0., 20.);
+ fHistFAKEetaCharge[i] = new TH1F( Form("fHistFAKEetaCharge%s",charge[i].Data()),
+ "p_{T} distribution: all Fake",
+ 100, 0., 20.);
+
+ fHistRECphiCharge[i] = new TH1F(Form("fHistRECphiCharge%s",charge[i].Data()),
+ "p_{T} distribution: all reconstructed",
+ 100, 0., 20.);
+ fHistMCphiCharge[i] = new TH1F(Form("fHistMCphiCharge%s",charge[i].Data()),
+ "p_{T} distribution: all MC",
+ 100, 0., 20.);
+ fHistMCNRphiCharge[i] = new TH1F(Form("fHistMCNRphiCharge%s",charge[i].Data()),
+ "p_{T} distribution: all not-reconstructed MC",
+ 100, 0., 20.);
+ fHistFAKEphiCharge[i] = new TH1F( Form("fHistFAKEphiCharge%s",charge[i].Data()),
+ "p_{T} distribution: all Fake",
+ 100, 0., 20.);
+
+ }
+ // BKC
+
+
+ fHists->SetOwner();
+
+ fHists->Add(fHistRECpt);
+ fHists->Add(fHistMCpt);
+ fHists->Add(fHistMCNRpt);
+ fHists->Add(fHistFAKEpt);
+ fHists->Add(fHistMULTpt);
+
+ fHists->Add(fHistRECeta);
+ fHists->Add(fHistMCeta);
+ fHists->Add(fHistMCNReta);
+ fHists->Add(fHistFAKEeta);
+ fHists->Add(fHistMULTeta);
+
+ fHists->Add(fHistRECphi);
+ fHists->Add(fHistMCphi);
+ fHists->Add(fHistMCNRphi);
+ fHists->Add(fHistFAKEphi);
+ fHists->Add(fHistMULTphi);
+
+ fHists->Add(fHistRECHPTeta);
+ fHists->Add(fHistMCHPTeta);
+ fHists->Add(fHistMCNRHPTeta);
+ fHists->Add(fHistFAKEHPTeta);
+
+ fHists->Add(fHistRECHPTphi);
+ fHists->Add(fHistMCHPTphi);
+ fHists->Add(fHistMCNRHPTphi);
+ fHists->Add(fHistFAKEHPTphi);
+
+ // CKB
+
+ fHists->Add(fh1VtxEff);
+ for(int i = 0;i < kTotalAC;++i){
+ fHists->Add(fh2PhiPadRow[i]);
+ fHists->Add(fh2PhiLayer[i]);
+ }
+ for(int i = 0;i<2;++i){
+ fHists->Add(fh2MultSpdChips[i]);
+ }
+ fHists->Add(fh2VtxTpcSpd);
+
+ for(int i = 0; i < kTotalVtx; ++i){
+ fHists->Add(fh2VertexCorrelation[i]);
+ fHists->Add(fh2VertexCorrelationShift[i]);
+ fHists->Add(fh1VertexShift[i]);
+ fHists->Add(fh1VertexShiftNorm[i]);
+ }
+
+
+ fHists->Add(fh2EtaCorrelation);
+ fHists->Add(fh2EtaCorrelationShift);
+ fHists->Add(fh2PhiCorrelation);
+ fHists->Add(fh2PhiCorrelationShift);
+ fHists->Add(fh2PtCorrelation);
+ fHists->Add(fh2PtCorrelationShift);
+
+ fHists->Add(fHistDeltaphiprimaries);
+ fHists->Add(fHistDeltaphisecondaries);
+ fHists->Add(fHistDeltaphireject);
+
+ fHists->Add(fHistTPCITSdeltax);
+ fHists->Add(fHistTPCITSdeltay);
+ fHists->Add(fHistTPCITSdeltaz);
+ fHists->Add(fHistTPCITSdeltar);
+
+ fHists->Add(fHistTPCTRDdeltax);
+ fHists->Add(fHistTPCTRDdeltay);
+ fHists->Add(fHistTPCTRDdeltaz);
+ fHists->Add(fHistTPCTRDdeltar);
+ fHists->Add(fHistTPCTRDdeltaphi);
+
+ fHists->Add(fHistRecMult);
+ fHists->Add(fHistNCluster);
+
+
+ fHists->Add(fPtPullsPos);
+ fHists->Add(fPtPullsNeg);
+ fHists->Add(fPtShiftPos);
+ fHists->Add(fPtShiftNeg);
+
+ fHists->Add(fVertexRvsZ);
+ fHists->Add(fVertexRvsZC);
+ fHists->Add(fEtaMultiFluc);
+ fHists->Add(fEtaMulti);
+ fHists->Add(fEtaMultiH);
+ fHists->Add(fPhiMultiFluc);
+ fHists->Add(fPhiMulti);
+ fHists->Add(fPhiMultiH);
+
+ for(Int_t i=0;i<8;i++){
+ fHists->Add(fHistRECptCharge[i]);
+ fHists->Add(fHistMCptCharge[i]);
+ fHists->Add(fHistMCNRptCharge[i]);
+ fHists->Add(fHistFAKEptCharge[i]);
+
+ fHists->Add(fHistRECetaCharge[i]);
+ fHists->Add(fHistMCetaCharge[i]);
+ fHists->Add(fHistMCNRetaCharge[i]);
+ fHists->Add(fHistFAKEetaCharge[i]);
+
+ fHists->Add(fHistRECphiCharge[i]);
+ fHists->Add(fHistMCphiCharge[i]);
+ fHists->Add(fHistMCNRphiCharge[i]);
+ fHists->Add(fHistFAKEphiCharge[i]);
+ }
+
+ for (Int_t i=0; i<fHists->GetEntries(); ++i) {
+ TH1 *h1 = dynamic_cast<TH1*>(fHists->At(i));
+ if (h1){
+ h1->Sumw2();
+ }
+ }
+ // BKC
+
+
+ // Post output data.
+ PostData(1, fHists);
+
+}
+
+//________________________________________________________________________
+void AliAnalysisTaskEfficiency::UserExec(Option_t *)
+{
+ // Event loop
+ // Analysis of MC true particles and reconstructed tracks
+ // Different track types (Global, TPC, ITS) can be selected
+
+ if (!fInputEvent) {
+ Printf("ERROR: fESD not available");
+ return;
+ }
+ static Int_t nfc = 0;
+
+ //AliESDInputHandler* esdI = (AliESDInputHandler*)
+ //(AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler();
+ //AliESDEvent* hltEvent = esdI->GetHLTEvent();
+
+ AliStack* stack = MCEvent()->Stack();
+ AliESDEvent* esdE = (AliESDEvent*) fInputEvent;
+
+ // CKB
+ // Fetch the SPD vertex:
+ Double_t vtxSPD[3];
+ const AliESDVertex* spdVertex = esdE->GetPrimaryVertexSPD();
+ if (spdVertex) {
+ if(spdVertex->GetNContributors()<=0){
+ spdVertex = 0;
+ }
+ else{
+ spdVertex->GetXYZ(vtxSPD);
+ }
+ }
+
+ // Fetch the TPC vertex
+ Double_t vtxTPC[3];
+ const AliESDVertex* tpcVertex = esdE->GetPrimaryVertexTPC();
+ if (tpcVertex) {
+ if(tpcVertex->GetNContributors()<=0){
+ tpcVertex = 0;
+ }
+ else{
+ tpcVertex->GetXYZ(vtxTPC);
+ }
+ }
+ // SPD Vertex
+ if (spdVertex) {
+ Double_t x,y,z;
+ x = spdVertex->GetX() + 0.07;
+ y = spdVertex->GetY() - 0.25;
+ z = spdVertex->GetZ();
+ if (TMath::Abs(z) < 10.) {
+ fVertexRvsZ->Fill(x,y);
+ if (TMath::Sqrt(x*x + y*y) > 5. * 0.028)
+ fVertexRvsZC->Fill(x,y);
+ }
+ }
+ // BKC
+ //Printf("%s:%d %5d",(char*)__FILE__,__LINE__, esdE->GetNumberOfTracks());
+ TArrayI labels(esdE->GetNumberOfTracks());
+ Int_t igood = 0;
+ // Track loop to fill a pT spectrum
+ //printf("ESD track loop \n");
+
+ AliESDtrack *tpcP = 0x0;
+
+ for (Int_t iTracks = 0; iTracks < esdE->GetNumberOfTracks(); iTracks++) {
+
+ //prevent mem leak for TPConly track
+ if(fTrackType==2&&tpcP){
+ delete tpcP;
+ tpcP = 0;
+ }
+
+ AliVParticle *track = esdE->GetTrack(iTracks);
+ AliESDtrack *esdtrack = dynamic_cast<AliESDtrack*>(track);
+ esdtrack->PropagateToDCA(esdE->GetPrimaryVertex(),
+ esdE->GetMagneticField(), 10000.);
+
+ if (!track) {
+ Printf("ERROR: Could not receive track %d", iTracks);
+ continue;
+ }
+ //__________
+ // run Task for global tracks or ITS tracks or TPC tracks
+ const AliExternalTrackParam *tpcPin = 0x0;
+ //Double_t phiIn=0.;
+ Float_t phi = 0.;
+ Float_t eta = 0.;
+ Float_t pt = 0.;
+ Int_t indexAC = GetIndexAC(esdtrack);
+
+ //select in the steering macro, which track type should be analysed.
+ // Four track types are selectable:
+ // 0. Global tracks
+ // 1. ITS tracks (SA or Pure SA)
+ // 2. TPC only tracks
+ // Track selection copied from PWGPP/AliAnalysisTaskQAsym
+
+ if(fTrackType==0){
+ //Fill all histograms with global tracks
+ tpcP = esdtrack;
+ if (!tpcP) continue;
+ if (!fCuts->AcceptTrack(tpcP)) continue;
+ phi = tpcP->Phi();
+ eta = tpcP->Eta();
+ pt = tpcP->Pt();
+ }
+ else if(fTrackType==1){
+ //Fill all histograms with ITS tracks
+ tpcP = esdtrack;
+ if (!tpcP) continue;
+ if (!fCuts->AcceptTrack(tpcP)) continue;
+ phi = tpcP->Phi();
+ eta = tpcP->Eta();
+ pt = tpcP->Pt();
+ }
+ else if(fTrackType==2){
+ //Fill all histograms with TPC track information
+ tpcPin = esdtrack->GetInnerParam();
+ if (!tpcPin) continue;
+ tpcP = AliESDtrackCuts::
+ GetTPCOnlyTrack(dynamic_cast<AliESDEvent*>(esdE),
+ esdtrack->GetID());
+ if (!tpcP) continue;
+ if (!fCuts->AcceptTrack(tpcP)) continue;
+ if(tpcP->GetNcls(1)>160)continue;//jacek's track cut
+ if(tpcP->GetConstrainedChi2TPC()<0)continue; // jacek's track cut
+ phi=tpcPin->Phi();
+ eta=tpcPin->Eta();
+ pt=tpcPin->Pt();
+ }
+ else{
+ Printf("ERROR: wrong track type \n");
+ continue;
+ }
+ //___________
+ //
+
+
+ labels.AddAt(-1, iTracks);
+
+
+ // Selected
+ if (TMath::Abs(eta) > 0.9) {
+ } else {
+
+ //Int_t charge=track->Charge();
+ // Within acceptance
+ Int_t ind = TMath::Abs(esdtrack->GetLabel());
+ AliMCParticle* mcPtest = (AliMCParticle*) MCEvent()->GetTrack(ind);
+ if (stack->IsPhysicalPrimary(mcPtest->Label())){
+ if(TMath::Abs(mcPtest->PdgCode())>=1000020030){//all helium (+-1000020030,+-1000020040)
+ pt*=2;// reconstruction takes charge=1 -> for particles with charge=2, pT,rec need to be corrected
+ }
+
+ Int_t index=ConvertePDG(mcPtest->PdgCode());
+ if(fDebug>1)Printf("PDG=%d, index=%d", mcPtest->PdgCode(), index);
+ //fill tracks comming from d,t,3He, 4He
+ if(index<8){
+ fHistRECptCharge [index] ->Fill(pt);
+ fHistRECetaCharge[index]->Fill(eta);
+ fHistRECphiCharge[index]->Fill(phi);
+ }
+
+ }
+ fHistRECpt ->Fill(pt);
+ fHistRECeta->Fill(eta);
+ fHistRECphi->Fill(phi);
+
+ if (pt > 2.) {
+ fHistRECHPTeta->Fill(eta);
+ fHistRECHPTphi->Fill(phi);
+ }
+
+ // Int_t ind = TMath::Abs(esdtrack->GetLabel());
+ AliMCParticle* mcP = (AliMCParticle*) MCEvent()->GetTrack(ind);
+ if (!(stack->IsPhysicalPrimary(mcP->Label()))) {
+ fHistFAKEpt ->Fill(pt);
+ fHistFAKEeta->Fill(eta);
+ fHistFAKEphi->Fill(phi);
+
+ if (pt > 2.) {
+ fHistFAKEHPTeta->Fill(eta);
+ fHistFAKEHPTphi->Fill(phi);
+ }
+
+ }
+
+ labels.AddAt(TMath::Abs(esdtrack->GetLabel()), iTracks);
+ igood++;
+
+
+ Float_t phiDeg = TMath::RadToDeg()*track->Phi();
+
+ //TPC-ITS matching
+ //TPC-TRD matching
+
+ if (tpcP){
+ Int_t labeltpcits = esdtrack->GetLabel();
+ AliExternalTrackParam * tpcPCopy = new AliExternalTrackParam(*tpcP);
+ Double_t xk = 43.6; // cm
+ Double_t bz = 5.0; // kG
+ if(tpcPCopy->PropagateTo(xk,bz)){
+ Double_t alpha=tpcPCopy->GetAlpha();
+ // if(tpcPCopy->Rotate(0.)){
+ Float_t xtpc=tpcPCopy->GetX();
+ Float_t ytpc=tpcPCopy->GetY();
+ Float_t ztpc=tpcPCopy->GetZ();
+ Float_t xpr,ypr ;
+ xpr=xtpc*TMath::Cos(alpha)-ytpc*TMath::Sin(alpha);
+ ypr=xtpc*TMath::Sin(alpha)+ytpc*TMath::Cos(alpha);
+
+ AliMCParticle* mcPart = (AliMCParticle*) MCEvent()->GetTrack(abs(labeltpcits));
+ Int_t ntref = mcPart->GetNumberOfTrackReferences();
+
+ for (Int_t k = 0; k < ntref; k++) {
+ AliTrackReference* ref = mcPart->GetTrackReference(k);
+ Float_t xhits = ref->X();
+ Float_t yhits = ref->Y();
+ Float_t radio = TMath::Sqrt(xhits*xhits+yhits*yhits);
+
+ if(ref->DetectorId() == 0 && (radio > 42.8 && radio < 43.7)) {
+
+ Float_t xits = ref->X();
+ Float_t yits = ref->Y();
+ Float_t zits = ref->Z();
+
+ Float_t difx=(xits-xpr);
+ fHistTPCITSdeltax->Fill(difx);
+ Float_t dify=(yits-ypr);
+ fHistTPCITSdeltay->Fill(dify);
+ Float_t difz=(zits-ztpc);
+ fHistTPCITSdeltaz->Fill(difz);
+ Float_t deltar = TMath::Sqrt(difx * difx + dify * dify + difz * difz);
+ fHistTPCITSdeltar->Fill(deltar);
+ break;
+ }
+ // }
+ } // trackrefs
+ }
+ } //ITS-TPC maching
+ //TPC-TRD maching
+
+ const AliExternalTrackParam *trd = esdtrack->GetOuterParam();
+
+ if (trd){Int_t labeltpctrd = track->GetLabel();
+
+ AliExternalTrackParam * trdCopy = new AliExternalTrackParam(*trd);
+ Double_t xktrd=296.0;
+ Double_t bztrd=5.0;
+ if(trdCopy->PropagateTo(xktrd,bztrd)){
+ Float_t xtpc2=trdCopy->GetX();
+ Float_t ytpc2=trdCopy->GetY();
+ Float_t ztpc2=trdCopy->GetZ();
+ Double_t alpha=trdCopy->GetAlpha();
+ Float_t xpr,ypr ;
+ xpr=xtpc2*TMath::Cos(alpha)-ytpc2*TMath::Sin(alpha);
+ ypr=xtpc2*TMath::Sin(alpha)+ytpc2*TMath::Cos(alpha);
+ AliMCParticle* mcPart = (AliMCParticle*) MCEvent()->GetTrack(abs(labeltpctrd));
+ Int_t ntreftrd = mcPart->GetNumberOfTrackReferences();
+
+ for (Int_t k = 0; k < ntreftrd; k++) {
+ AliTrackReference* ref = mcPart->GetTrackReference(k);
+ if(ref->DetectorId() == 3 && ref->Label() == labeltpctrd){
+
+ Float_t xtrd = ref->X();
+ Float_t ytrd = ref->Y();
+ Float_t ztrd = ref->Z();
+
+ Float_t difx=(xtrd-xpr);
+ fHistTPCTRDdeltax->Fill(difx);
+ Float_t dify=(ytrd-ypr);
+ fHistTPCTRDdeltay->Fill(dify);
+ Float_t difz=(ztrd-ztpc2);
+ fHistTPCTRDdeltaz->Fill(difz);
+ Float_t deltar = TMath::Sqrt(difx * difx + dify * dify + difz * difz);
+ fHistTPCTRDdeltar->Fill(deltar);
+ Float_t phi_tpc = TMath::ATan2(ypr, xpr);
+ Float_t phi_trd = TMath::ATan2(ytrd, xtrd);
+ Float_t dphi = phi_trd - phi_tpc;
+ if (dphi > TMath::Pi()) dphi -= 2. * TMath::Pi();
+ if (dphi < - TMath::Pi()) dphi += 2. * TMath::Pi();
+
+ fHistTPCTRDdeltaphi->Fill(dphi);
+ break;
+ }
+
+ }
+ }
+
+ }//TRD-TPC maching
+
+ // CKB
+ if(pt>2.&&fFieldOn){
+ TBits bits(esdtrack->GetTPCClusterMap());
+ for(unsigned int ip = 0;ip < bits.GetNbits();++ip){
+ if(bits[ip]){
+ fh2PhiPadRow[indexAC]->Fill(phiDeg,ip);
+ }
+ }
+ for(int i = 0;i < 6;++i){
+ if(esdtrack->HasPointOnITSLayer(i))fh2PhiLayer[indexAC]->Fill(phiDeg,i);
+ }
+ }
+ else if(!fFieldOn){ // field not on all momenta are set to MPV
+ TBits bits(esdtrack->GetTPCClusterMap());
+ for(unsigned int ip = 0;ip < bits.GetNbits();++ip){
+ if(bits[ip]){
+ fh2PhiPadRow[indexAC]->Fill(phiDeg,ip);
+ }
+ for(int i = 0;i < 6;++i){
+ if(esdtrack->HasPointOnITSLayer(i))fh2PhiLayer[indexAC]->Fill(phiDeg,i);
+ }
+ }
+
+ }
+ // Fill the correlation
+ if(MCEvent()){
+ TParticle *part = MCEvent()->Stack()->Particle(TMath::Abs(esdtrack->GetLabel()));
+ if(part){
+ Float_t mcEta = part->Eta();
+ Float_t mcPhi = TMath::RadToDeg()*part->Phi();
+ Float_t mcPt = part->Pt();
+ // if (pt - mcPt > 20.) {
+ fh2EtaCorrelation->Fill(mcEta,eta);
+ fh2EtaCorrelationShift->Fill(eta,mcEta-eta);
+ fh2PhiCorrelation->Fill(mcPhi,phiDeg);
+ fh2PhiCorrelationShift->Fill(phiDeg,mcPhi-phiDeg);
+ fh2PtCorrelation->Fill(mcPt,pt);
+ fh2PtCorrelationShift->Fill(pt,mcPt-pt);
+ //}
+ Double_t sigmaPt = TMath::Sqrt(esdtrack->GetSigma1Pt2());
+ if (track->Charge() > 0.) {
+ fPtPullsPos->Fill(mcPt, (1./pt - 1./mcPt) / sigmaPt);
+ fPtShiftPos->Fill(mcPt, (1./pt - 1./mcPt));
+ } else {
+ fPtPullsNeg->Fill(mcPt, (1./pt - 1./mcPt) / sigmaPt);
+ fPtShiftNeg->Fill(mcPt, (1./pt - 1./mcPt));
+ }
+ }
+ // BKC
+ }
+ }
+
+ } //track loop
+
+ //prevent mem leak for TPConly track
+ if(fTrackType==2&&tpcP){
+ delete tpcP;
+ tpcP = 0;
+ }
+
+ //Printf("%s:%d",(char*)__FILE__,__LINE__);
+ // CKB
+ // Vertex eff
+ if(tpcVertex){
+ fh1VtxEff->Fill("TPC VTX",1);
+ }
+ else{
+ fh1VtxEff->Fill("NO TPC VTX",1);
+ }
+
+ if(spdVertex){
+ fh1VtxEff->Fill("SPD VTX",1);
+ }
+ else{
+ fh1VtxEff->Fill("NO SPD VTX",1);
+ }
+
+
+ // Vertex correlation SPD vs. TPC
+ if(tpcVertex&&spdVertex){
+ fh2VtxTpcSpd->Fill(vtxTPC[2],vtxSPD[2]);
+ }
+ // Printf("%s:%d",(char*)__FILE__,__LINE__);
+ // Multiplicity checks in the SPD
+ const AliMultiplicity *spdMult = esdE->GetMultiplicity();
+ // Multiplicity Analysis
+ Int_t nt = spdMult->GetNumberOfTracklets();
+ nfc++;
+ if (nfc == 520) {
+ nfc = 0;
+ for (Int_t ib = 1; ib <= 10; ib++) {
+ Int_t ic = Int_t(fEtaMulti->GetBinContent(ib));
+ Float_t eta = -1.8 + Float_t(ib-1) * 0.4;
+ fEtaMultiFluc->Fill(eta, Float_t(ic));
+ if (ib == 5) fEtaMultiH->Fill(Float_t(ic));
+ }
+
+ for (Int_t ib = 1; ib <= 64; ib++) {
+ Int_t ic = Int_t(fPhiMulti->GetBinContent(ib));
+ Float_t phi = 0.05 + Float_t(ib-1) * 0.1;
+ fPhiMultiFluc->Fill(phi, Float_t(ic));
+ if (ib == 2) fPhiMultiH->Fill(Float_t(ic));
+ }
+
+ fEtaMulti->Reset();
+ fPhiMulti->Reset();
+ }
+
+ for (Int_t j = 0; j < nt; j++) {
+ fEtaMulti->Fill(spdMult->GetEta(j));
+ fPhiMulti->Fill(spdMult->GetPhi(j));
+ }
+
+ Short_t chips0 = spdMult->GetNumberOfFiredChips(0);
+ Short_t chips1 = spdMult->GetNumberOfFiredChips(1);
+ //Printf("%s:%d",(char*)__FILE__,__LINE__);
+ Int_t inputCount = 0;
+
+ if(spdVertex){
+ // get multiplicity from ITS tracklets
+ for (Int_t i=0; i<spdMult->GetNumberOfTracklets(); ++i){
+ Float_t deltaPhi = spdMult->GetDeltaPhi(i);
+ // prevent values to be shifted by 2 Pi()
+ if (deltaPhi < -TMath::Pi())
+ deltaPhi += TMath::Pi() * 2;
+ if (deltaPhi > TMath::Pi())
+ deltaPhi -= TMath::Pi() * 2;
+ if (TMath::Abs(deltaPhi) > 1)
+ printf("WARNING: Very high Delta Phi: %d %f %f %f\n", i, spdMult->GetTheta(i), spdMult->GetPhi(i), deltaPhi);
+ //trackled Deltaphi
+ Int_t label1=spdMult->GetLabel(i,0);
+ Int_t label2=spdMult->GetLabel(i,1);
+ if ( label1==label2){
+ if(stack->IsPhysicalPrimary(label1) == kTRUE)
+ fHistDeltaphiprimaries->Fill(deltaPhi);
+ else
+ fHistDeltaphisecondaries->Fill(deltaPhi);
+ }
+ else{
+ fHistDeltaphireject->Fill(deltaPhi);
+ }
+ ++inputCount;
+ }
+ // Printf("%s:%d",(char*)__FILE__,__LINE__);
+ fh2MultSpdChips[0]->Fill(chips0,inputCount);
+ fh2MultSpdChips[1]->Fill(chips1,inputCount);
+ //Printf("%s:%d",(char*)__FILE__,__LINE__);
+ }
+ // BKC
+
+ // MC analysis
+ Int_t nmc = MCEvent()->GetNumberOfTracks();
+ AliGenEventHeader* header = MCEvent()->GenEventHeader();
+ // some test
+ AliGenCocktailEventHeader* cH = dynamic_cast<AliGenCocktailEventHeader*> (header);
+ AliGenPythiaEventHeader* pH;
+ if (cH == 0)
+ {
+ header->Dump();
+ pH = dynamic_cast<AliGenPythiaEventHeader*>(header);
+ } else {
+ TObject* entry = cH->GetHeaders()->FindObject("Pythia");
+ pH = dynamic_cast<AliGenPythiaEventHeader*>(entry);
+ }
+ Int_t iproc = -1;
+
+ if (pH) iproc = pH->ProcessType();
+
+ TArrayF mcV;
+ header->PrimaryVertex(mcV);
+ Float_t vzMC = mcV[2];
+
+ // CKB
+ // Fill the correlation with MC
+ if(spdVertex&&MCEvent()){
+ Float_t diff = vzMC-vtxSPD[2];
+ fh2VertexCorrelation[kSPD]->Fill(vzMC,vtxSPD[2]);
+ fh2VertexCorrelationShift[kSPD]->Fill(vzMC,diff);
+ fh1VertexShift[kSPD]->Fill(diff);
+ if(spdVertex->GetZRes()>0)fh1VertexShiftNorm[kSPD]->Fill(diff/spdVertex->GetZRes());
+ }
+ if(tpcVertex&&MCEvent()){
+ Float_t diff = vzMC-vtxTPC[2];
+ fh2VertexCorrelation[kTPC]->Fill(vzMC,vtxTPC[2]);
+ fh2VertexCorrelationShift[kTPC]->Fill(vzMC,diff);
+ fh1VertexShift[kTPC]->Fill(diff);
+ if(tpcVertex->GetZRes()>0)fh1VertexShiftNorm[kTPC]->Fill(diff/tpcVertex->GetZRes());
+ }
+ // BKC
+
+
+ // MC event loop
+ //printf("MC particle loop \n");
+ for (Int_t i = 0; i < nmc; i++) {
+ AliMCParticle* mcP = (AliMCParticle*) MCEvent()->GetTrack(i);
+ //printf("MC particle loop %5d \n", i);
+ // Primaries only
+ if (!(stack->IsPhysicalPrimary(mcP->Label()))) continue;
+ if (mcP->Particle()->GetPDG()->Charge() == 0) continue;
+ // Int_t charge=Int_t(mcP->Particle()->GetPDG()->Charge());
+ Float_t phi = mcP->Phi();
+ Float_t eta = mcP->Eta();
+ Float_t pt = mcP->Pt();
+ if (TMath::Abs(eta) > 0.9) continue;
+ Int_t ntr = mcP->GetNumberOfTrackReferences();
+ Int_t nITS = 0;
+ Int_t nTPC = 0;
+ Int_t nFRA = 0;
+ Float_t x = 0.;
+ Float_t y = 0.;
+
+ for (Int_t j = 0; j < ntr; j++) {
+
+ AliTrackReference* ref = mcP->GetTrackReference(j);
+
+ if (ref->DetectorId() == 0) nITS++;
+ if (ref->DetectorId() == 1) nTPC++;
+ if (ref->DetectorId() == 2) nFRA++;
+ if (nTPC == 1) {
+ x = ref->X();
+ y = ref->Y();
+ break;
+ }
+ }
+
+ fHistMCpt ->Fill(pt);
+ fHistMCeta->Fill(eta);
+ fHistMCphi->Fill(phi);
+ Int_t index=ConvertePDG(mcP->PdgCode());
+ if(index<8){
+ fHistMCptCharge [index] ->Fill(pt);
+ fHistMCetaCharge[index]->Fill(eta);
+ fHistMCphiCharge[index]->Fill(phi);
+ }
+
+
+ if (pt > 2.) {
+ fHistMCHPTeta->Fill(eta);
+ fHistMCHPTphi->Fill(phi);
+ }
+
+ Bool_t reco = kFALSE;
+ Int_t multR = 0;
+ Int_t jold = -1;
+
+ for (Int_t j = 0; j < esdE->GetNumberOfTracks(); j++) {
+ if (i == labels.At(j)) {
+ reco = kTRUE;
+ multR++;
+ //AliESDtrack* test = esdE->GetTrack(j);
+ if (multR > 1) {
+ Int_t nclus = 0;
+ AliESDtrack* track = esdE->GetTrack(jold);
+ nclus = track->GetTPCclusters(0);
+ fHistNCluster->Fill(nclus);
+
+
+ track = esdE->GetTrack(j);
+ nclus = track->GetTPCclusters(0);
+ fHistNCluster->Fill(nclus);
+
+
+ }
+ jold = j;
+ }
+ }
+
+ fHistRecMult->Fill(Float_t(multR), 1.);
+ if (multR == 0) {
+ fHistMCNRpt ->Fill(pt);
+ fHistMCNReta->Fill(eta);
+ fHistMCNRphi->Fill(phi);
+ if(index<8){
+ fHistMCNRptCharge [index] ->Fill(pt);
+ fHistMCNRetaCharge[index]->Fill(eta);
+ fHistMCNRphiCharge[index]->Fill(phi);
+ }
+
+
+ if (pt > 2.) {
+ fHistMCNRHPTeta->Fill(eta);
+ fHistMCNRHPTphi->Fill(phi);
+ }
+ }
+ else if (multR > 1)
+ {
+ fHistMULTpt ->Fill(pt);
+ fHistMULTeta->Fill(eta);
+ fHistMULTphi->Fill(phi);
+ }
+ } // MC particle loop
+ //printf("End of MC particle loop \n");
+
+
+}
+
+//________________________________________________________________________
+void AliAnalysisTaskEfficiency::Terminate(Option_t *)
+{
+
+}
+
+
+Bool_t AliAnalysisTaskEfficiency::SelectJouri(AliESDtrack* track)
+{
+
+ Bool_t selected = kTRUE;
+ AliESDEvent* esdE = (AliESDEvent*) fInputEvent;
+ // > 50 Clusters
+ if (track->GetTPCclusters(0) < 50) selected = kFALSE;
+
+ const AliESDVertex *vtx=esdE->GetPrimaryVertexSPD();
+ if (!vtx->GetStatus()) selected = kFALSE;
+
+ Double_t zv=vtx->GetZv();
+
+
+ const AliExternalTrackParam *ct=track->GetTPCInnerParam();
+ if (!ct) selected = kFALSE;
+
+ Double_t xyz[3];
+ ct->GetXYZ(xyz);
+ Float_t rv = TMath::Sqrt(xyz[0] * xyz[0] + xyz[1] * xyz[1]);
+ if (rv > 3.0) selected = kFALSE;
+ if (TMath::Abs(xyz[2] - zv)>2.5) selected = kFALSE;
+ return (selected);
+
+}
+
+Int_t AliAnalysisTaskEfficiency::GetIndexAC(AliESDtrack *track){
+ if(!track)return -1;
+
+ // Crude selection for A and C side
+ // just with eta
+ if (track->Eta() > 0) {
+ if (track->Charge() > 0)
+ return kPosA;
+ else
+ return kNegA;
+ }
+ else { // C side
+ if (track->Charge() > 0)
+ return kPosC;
+ else
+ return kNegC;
+ }
+
+ return -1;
+
+}
+
+
+Int_t AliAnalysisTaskEfficiency::ConvertePDG(Int_t pdg)
+{
+
+ // function converts the pdg values for nuclei (d, t, 3He,
+ // 4He and their antiparticles) into indices for histo-array
+ // filled in UserExec
+
+ if ( pdg == 1000010020 )return 0;
+ else if ( pdg == -1000010020 )return 1;
+ else if ( pdg == 1000010030 )return 2;
+ else if ( pdg == -1000010030 )return 3;
+ else if ( pdg == 1000020030 )return 4;
+ else if ( pdg == -1000020030 )return 5;
+ else if ( pdg == 1000020040 )return 6;
+ else if ( pdg == -1000020040 )return 7;
+ else return 9;
+
+
+}