]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskSOH.cxx
Use pass-aware container
[u/mrichter/AliRoot.git] / PWGJE / EMCALJetTasks / UserTasks / AliAnalysisTaskSOH.cxx
index 78fe51156ff2e6e74f5f23e321aadd49752ac66d..f72c2747ec767fac33b8edc645974e2329c7c066 100644 (file)
@@ -11,6 +11,7 @@
 #include <TH3F.h>
 #include <THnSparse.h>
 
+#include "TChain.h"
 #include "AliAnalysisManager.h"
 #include "AliAnalysisTask.h"
 #include "AliEMCALRecoUtils.h"
@@ -34,34 +35,19 @@ AliAnalysisTaskSOH::AliAnalysisTaskSOH() :
   AliAnalysisTaskSE(), 
   fESD(0), 
   fMC(0), 
+  fZVtxMax(10),
   fEsdTrackCuts(0x0),
   fHybridTrackCuts1(0x0),
   fHybridTrackCuts2(0x0),
   fTrackIndices(0x0),
   fClusterIndices(0x0),
   fClusterArray(0x0),
+  fMcProcess(kTRUE),
+  fTrackProcess(kTRUE),
+  fSFProcess(kFALSE),
+  fClusterProcess(kFALSE),
   fOutputList(0x0),        
-  fHEventStat(0),      
-  fHTrkEffParGenPtEtaPhi(0), 
-  fHTrkEffDetGenPtEtaPhi(0), 
-  fHTrkEffDetRecPtEtaPhi(0),
-  fHTrkEffDetRecFakePtEtaPhi(0), 
-  fHTrkEffParGenPt_rmInj(0),
-  fHTrkEffDetGenPt_rmInj(0),
-  fHTrkEffParGenPt_PiRmInj(0),
-  fHTrkEffDetGenPt_PiRmInj(0),
-  fHTrkEffParGenPt_all(0),
-  fHTrkEffDetGenPt_all(0),
-  fHTrkEffParGenPt_Dch(0),
-  fHTrkEffDetGenPt_Dch(0),
-  fHTrkEffParGenPt_Dn(0),
-  fHTrkEffDetGenPt_Dn(0),
-  fHTrkEffParGenPt_Ds(0),
-  fHTrkEffDetGenPt_Ds(0),
-  fHTrkEffParGenPt_cL(0),
-  fHTrkEffDetGenPt_cL(0),
-  fHTrkEffParGenPt_JPsi(0),
-  fHTrkEffDetGenPt_JPsi(0),
+  fHEventStat(0),
   fHScaleFactor(0),
   fHScaleFactor100HC(0),
   fHEOverPVsPt(0x0),
@@ -88,10 +74,30 @@ AliAnalysisTaskSOH::AliAnalysisTaskSOH() :
   fHClsEoverMcE_All(0x0),
   fHClsEoverMcE_Photon(0x0),
   fHClsEoverMcE_Elec(0x0),
-  fHClsEoverMcE_Pion(0x0)
+  fHClsEoverMcE_Pion(0x0),
+  fHParGenPion_p(0x0),
+  fHParGenPion_m(0x0),
+  fHParGenPion_rmInj_p(0x0),
+  fHParGenPion_rmInj_m(0x0),
+  fHDetGenFakePion(0x0),
+  fHDetRecFakePion(0x0),  
+  fHDetGenSecPion(0x0),
+  fHDetRecSecPion(0x0)
 {
-  // Constructor
-
+  for(Int_t i=0; i<3; i++)
+  {
+    fHDetGenPion_p[i]        = 0x0;   
+    fHDetRecPion_p[i]        = 0x0;  
+    fHDetGenPion_m[i]        = 0x0;  
+    fHDetRecPion_m[i]        = 0x0;  
+    fHDetGenPion_rmInj_p[i]  = 0x0;  
+    fHDetRecPion_rmInj_p[i]  = 0x0;
+    fHDetGenPion_rmInj_m[i]  = 0x0;
+    fHDetRecPion_rmInj_m[i]  = 0x0;
+  }
+  DefineInput (0, TChain::Class());
+  DefineOutput(1, TList::Class());
+  
 }
 
 
@@ -100,34 +106,19 @@ AliAnalysisTaskSOH::AliAnalysisTaskSOH(const char *name) :
   AliAnalysisTaskSE(name), 
   fESD(0), 
   fMC(0), 
+  fZVtxMax(10),
   fEsdTrackCuts(0x0),
   fHybridTrackCuts1(0x0),
   fHybridTrackCuts2(0x0),
   fTrackIndices(0x0),
   fClusterIndices(0x0),
   fClusterArray(0x0),
+  fMcProcess(kTRUE),
+  fTrackProcess(kTRUE),
+  fSFProcess(kFALSE),
+  fClusterProcess(kFALSE),
   fOutputList(0x0),        
-  fHEventStat(0),      
-  fHTrkEffParGenPtEtaPhi(0), 
-  fHTrkEffDetGenPtEtaPhi(0), 
-  fHTrkEffDetRecPtEtaPhi(0),
-  fHTrkEffDetRecFakePtEtaPhi(0), 
-  fHTrkEffParGenPt_rmInj(0),
-  fHTrkEffDetGenPt_rmInj(0),
-  fHTrkEffParGenPt_PiRmInj(0),
-  fHTrkEffDetGenPt_PiRmInj(0),
-  fHTrkEffParGenPt_all(0),
-  fHTrkEffDetGenPt_all(0),
-  fHTrkEffParGenPt_Dch(0),
-  fHTrkEffDetGenPt_Dch(0),
-  fHTrkEffParGenPt_Dn(0),
-  fHTrkEffDetGenPt_Dn(0),
-  fHTrkEffParGenPt_Ds(0),
-  fHTrkEffDetGenPt_Ds(0),
-  fHTrkEffParGenPt_cL(0),
-  fHTrkEffDetGenPt_cL(0),
-  fHTrkEffParGenPt_JPsi(0),
-  fHTrkEffDetGenPt_JPsi(0),
+  fHEventStat(0), 
   fHScaleFactor(0),
   fHScaleFactor100HC(0),
   fHEOverPVsPt(0x0),
@@ -154,10 +145,29 @@ AliAnalysisTaskSOH::AliAnalysisTaskSOH(const char *name) :
   fHClsEoverMcE_All(0x0),
   fHClsEoverMcE_Photon(0x0),
   fHClsEoverMcE_Elec(0x0),
-  fHClsEoverMcE_Pion(0x0)
+  fHClsEoverMcE_Pion(0x0),
+  fHParGenPion_p(0x0),
+  fHParGenPion_m(0x0),
+  fHParGenPion_rmInj_p(0x0),
+  fHParGenPion_rmInj_m(0x0),
+  fHDetGenFakePion(0x0),
+  fHDetRecFakePion(0x0),  
+  fHDetGenSecPion(0x0),
+  fHDetRecSecPion(0x0)
 {
-  // Constructor
+  for(Int_t i=0; i<3; i++)
+  {
+    fHDetGenPion_p[i]        = 0x0;   
+    fHDetRecPion_p[i]        = 0x0;  
+    fHDetGenPion_m[i]        = 0x0;  
+    fHDetRecPion_m[i]        = 0x0;  
+    fHDetGenPion_rmInj_p[i]  = 0x0;  
+    fHDetRecPion_rmInj_p[i]  = 0x0;
+    fHDetGenPion_rmInj_m[i]  = 0x0;
+    fHDetRecPion_rmInj_m[i]  = 0x0;
+  }
 
+  // Constructor
   // Output slot #1 writes into a TH1 container
   DefineOutput(1, TList::Class());
 }
@@ -168,12 +178,12 @@ AliAnalysisTaskSOH::~AliAnalysisTaskSOH()
 {
   // Destructor.
 
-  delete fEsdTrackCuts;
-  delete fHybridTrackCuts1;
-  delete fHybridTrackCuts2;
-  delete fTrackIndices;
-  delete fClusterIndices;
-  delete fClusterArray;
+  if(fEsdTrackCuts) delete fEsdTrackCuts;
+  if(fHybridTrackCuts1) delete fHybridTrackCuts1;
+  if(fHybridTrackCuts2) delete fHybridTrackCuts2;
+  if(fTrackIndices) delete fTrackIndices;
+  if(fClusterIndices) delete fClusterIndices;
+  if(fClusterArray) delete fClusterArray;
 }
 
 
@@ -198,169 +208,215 @@ void AliAnalysisTaskSOH::UserCreateOutputObjects()
   fHEventStat->GetXaxis()->SetBinLabel(8,"cls/>3-truth");
   fOutputList->Add(fHEventStat);
 
-  fHTrkEffParGenPtEtaPhi = new TH3F("fHTrkEffParGenPtEtaPhi","Particle level truth Phi-Eta-p_{T} distribution of primary charged pions", 500, 0, 100, 60, -1.5, 1.5, 128, 0, 6.4);
-  fHTrkEffParGenPtEtaPhi->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
-  fHTrkEffParGenPtEtaPhi->GetYaxis()->SetTitle("#eta");
-  fHTrkEffParGenPtEtaPhi->GetZaxis()->SetTitle("#phi");
-  fOutputList->Add(fHTrkEffParGenPtEtaPhi);
-  
-  fHTrkEffDetGenPtEtaPhi = new TH3F("fHTrkEffDetGenPtEtaPhi","Detector level truth Phi-Eta-p_{T} distribution of primary charged pions", 500, 0, 100, 60, -1.5, 1.5, 128, 0, 6.4);
-  fHTrkEffDetGenPtEtaPhi->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
-  fHTrkEffDetGenPtEtaPhi->GetYaxis()->SetTitle("#eta");
-  fHTrkEffDetGenPtEtaPhi->GetZaxis()->SetTitle("#phi");
-  fOutputList->Add(fHTrkEffDetGenPtEtaPhi);
-
-  fHTrkEffDetRecPtEtaPhi = new TH3F("fHTrkEffDetRecPtEtaPhi","Reconstructed track Phi-Eta-p_{T} distribution of primary charged pions", 500, 0, 100, 60, -1.5, 1.5, 128, 0, 6.4);
-  fHTrkEffDetRecPtEtaPhi->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
-  fHTrkEffDetRecPtEtaPhi->GetYaxis()->SetTitle("#eta");
-  fHTrkEffDetRecPtEtaPhi->GetZaxis()->SetTitle("#phi");
-  fOutputList->Add(fHTrkEffDetRecPtEtaPhi);
-
-  fHTrkEffDetRecFakePtEtaPhi = new TH3F("fHTrkEffDetRecFakePtEtaPhi","Reconstructed fake track Phi-Eta-p_{T} distribution of pions", 500, 0, 100, 60, -1.5, 1.5, 128, 0, 6.4);
-  fHTrkEffDetRecFakePtEtaPhi->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
-  fHTrkEffDetRecFakePtEtaPhi->GetYaxis()->SetTitle("#eta");
-  fHTrkEffDetRecFakePtEtaPhi->GetZaxis()->SetTitle("#phi");
-  fOutputList->Add(fHTrkEffDetRecFakePtEtaPhi);
-
-  fHTrkEffParGenPt_rmInj = new TH1F("fHTrkEffParGenPt_rmInj", "Truth w/o injected signal p_{T} distribution;p_{T}^{rec} (GeV/c)",500,0.0,100);
-  fOutputList->Add(fHTrkEffParGenPt_rmInj);
-
-  fHTrkEffDetGenPt_rmInj = new TH1F("fHTrkEffDetGenPt_rmInj", "generated detector level w/o injected signal track p_{T} distribution;p_{T}^{rec} (GeV/c)",500,0.0,100);
-  fOutputList->Add(fHTrkEffDetGenPt_rmInj);
-
-  fHTrkEffParGenPt_PiRmInj = new TH1F("fHTrkEffParGenPt_PiRmInj", "Truth pion w/o injected signal p_{T} distribution;p_{T}^{rec} (GeV/c)",500,0.0,100);
-  fOutputList->Add(fHTrkEffParGenPt_PiRmInj);
-
-  fHTrkEffDetGenPt_PiRmInj = new TH1F("fHTrkEffDetGenPt_PiRmInj", "generated detector level pion w/o injected signal track p_{T} distribution;p_{T}^{rec} (GeV/c)",500,0.0,100);
-  fOutputList->Add(fHTrkEffDetGenPt_PiRmInj);
-
-  fHTrkEffParGenPt_all = new TH1F("fHTrkEffParGenPt_all", "Truth p_{T} distribution;p_{T}^{rec} (GeV/c)",500,0.0,100);
-  fOutputList->Add(fHTrkEffParGenPt_all);
-
-  fHTrkEffDetGenPt_all = new TH1F("fHTrkEffDetGenPt_all", "generated detector level track p_{T} distribution;p_{T}^{rec} (GeV/c)",500,0.0,100);
-  fOutputList->Add(fHTrkEffDetGenPt_all);
-
-  fHTrkEffParGenPt_Dch = new TH1F("fHTrkEffParGenPt_Dch", "Charged D meson truth p_{T} distribution;p_{T}^{rec} (GeV/c)",500,0.0,100);
-  fOutputList->Add(fHTrkEffParGenPt_Dch);
-
-  fHTrkEffDetGenPt_Dch = new TH1F("fHTrkEffDetGenPt_Dch", "Charged D meson generated detector level track p_{T} distribution;p_{T}^{rec} (GeV/c)",500,0.0,100);
-  fOutputList->Add(fHTrkEffDetGenPt_Dch);
-
-  fHTrkEffParGenPt_Dn = new TH1F("fHTrkEffParGenPt_Dn", "Neutral D meson truth p_{T} distribution;p_{T}^{rec} (GeV/c)",500,0.0,100);
-  fOutputList->Add(fHTrkEffParGenPt_Dn);
-
-  fHTrkEffDetGenPt_Dn = new TH1F("fHTrkEffDetGenPt_Dn", "Neutral D meson generated detector level track p_{T} distribution;p_{T}^{rec} (GeV/c)",500,0.0,100);
-  fOutputList->Add(fHTrkEffDetGenPt_Dn);
-
-  fHTrkEffParGenPt_Ds = new TH1F("fHTrkEffParGenPt_Ds", "D_{s} truth p_{T} distribution;p_{T}^{rec} (GeV/c)",500,0.0,100);
-  fOutputList->Add(fHTrkEffParGenPt_Ds);
-
-  fHTrkEffDetGenPt_Ds = new TH1F("fHTrkEffDetGenPt_Ds", "D_{s} generated detector level track p_{T} distribution;p_{T}^{rec} (GeV/c)",500,0.0,100);
-  fOutputList->Add(fHTrkEffDetGenPt_Ds);
-
-  fHTrkEffParGenPt_cL = new TH1F("fHTrkEffParGenPt_cL", "#Lambda_{c} truth p_{T} distribution;p_{T}^{rec} (GeV/c)",500,0.0,100);
-  fOutputList->Add(fHTrkEffParGenPt_cL);
-
-  fHTrkEffDetGenPt_cL = new TH1F("fHTrkEffDetGenPt_cL", "#Lambda_{c} generated detector level track p_{T} distribution;p_{T}^{rec} (GeV/c)",500,0.0,100);
-  fOutputList->Add(fHTrkEffDetGenPt_cL);
-
-  fHTrkEffParGenPt_JPsi = new TH1F("fHTrkEffParGenPt_JPsi", "J/#Psi truth p_{T} distribution;p_{T}^{rec} (GeV/c)",500,0.0,100);
-  fOutputList->Add(fHTrkEffParGenPt_JPsi);
-
-  fHTrkEffDetGenPt_JPsi = new TH1F("fHTrkEffDetGenPt_JPsi", "J/#Psi generated detector level track p_{T} distribution;p_{T}^{rec} (GeV/c)",500,0.0,100);
-  fOutputList->Add(fHTrkEffDetGenPt_JPsi);
-
-  fHScaleFactor = new TH1F("fHScaleFactor", "Scale factor distribution without hadronic correction;Scale factor",100,0,10);
-  fOutputList->Add(fHScaleFactor);
-
-  fHScaleFactor100HC = new TH1F("fHScaleFactor100HC", "Scale factor distribution with 100% hadronic correction;Scale factor",100,0,10);
-  fOutputList->Add(fHScaleFactor100HC);
-
-  fHEOverPVsPt = new TH2F("fHEOverPVsPt", "E/P vs track p_{T}; p_{T} (GeV/c); E/P", 200 , 0, 4, 200, 0, 3.2);
-  fOutputList->Add(fHEOverPVsPt);
-
-  fHEMCalResponsePion = new TH2F("fHEMCalResponsePion", "Pion E/P vs track p_{T}; p_{T} (GeV/c); E/P", 100 , 0, 4, 100, 0, 3.2);
-  fOutputList->Add(fHEMCalResponsePion);
-
-  fHEMCalResponseElec = new TH2F("fHEMCalResponseElec", "Electron E/P vs track p_{T};  p_{T} (GeV/c); E/P", 100 , 0, 4, 100, 0, 3.2);
-  fOutputList->Add(fHEMCalResponseElec);
-
-  fHEMCalResponseProton = new TH2F("fHEMCalResponseProton", "Proton E/P vs track p_{T};  p_{T} (GeV/c); E/P", 100 , 0, 4, 100, 0, 3.2);
-  fOutputList->Add(fHEMCalResponseProton);
-
-  fHEMCalRecdPhidEta = new TH2F("fHEMCalRecdPhidEta","EMCAL Cluster-Track #Delta#phi-#Delta#eta; #Delta#eta; #Delta#phi",1000,-0.1,0.1,1000,-0.5,0.5);
-  fOutputList->Add(fHEMCalRecdPhidEta);
-
-  fHEMCalRecdPhidEtaP = new TH2F("fHEMCalRecdPhidEtaP","EMCAL Charge+ Cluster-Track #Delta#phi-#Delta#eta; #Delta#eta; #Delta#phi",1000,-0.1,0.1,1000,-0.5,0.5);
-  fOutputList->Add(fHEMCalRecdPhidEtaP);
-
-  fHEMCalRecdPhidEtaM = new TH2F("fHEMCalRecdPhidEtaM","EMCAL Charge- Cluster-Track #Delta#phi-#Delta#eta; #Delta#eta; #Delta#phi",1000,-0.1,0.1,1000,-0.5,0.5);
-  fOutputList->Add(fHEMCalRecdPhidEtaM);
-
-  fHEMCalRecdPhidEta_Truth = new TH2F("fHEMCalRecdPhidEta_Truth","EMCAL Cluster-Track(Truth matched) #Delta#phi-#Delta#eta; #Delta#eta; #Delta#phi",1000,-0.1,0.1,1000,-0.5,0.5);
-  fOutputList->Add(fHEMCalRecdPhidEta_Truth);
-
-  fHEMCalRecdPhidEtaP_Truth = new TH2F("fHEMCalRecdPhidEtaP_Truth","EMCAL Charge+ Cluster-Track(Truth matched) #Delta#phi-#Delta#eta; #Delta#eta; #Delta#phi",1000,-0.1,0.1,1000,-0.5,0.5);
-  fOutputList->Add(fHEMCalRecdPhidEtaP_Truth);
-
-  fHEMCalRecdPhidEtaM_Truth = new TH2F("fHEMCalRecdPhidEtaM_Truth","EMCAL Charge- Cluster-Track(Truth matched) #Delta#phi-#Delta#eta; #Delta#eta; #Delta#phi",1000,-0.1,0.1,1000,-0.5,0.5);
-  fOutputList->Add(fHEMCalRecdPhidEtaM_Truth);
-
-  fHEMCalRecdPhidEtaposEta = new TH2F("fHEMCalRecdPhidEtaposEta","(+eta track) EMCAL Cluster-Track #Delta#phi-#Delta#eta; #Delta#eta; #Delta#phi",1000,-0.1,0.1,1000,-0.5,0.5);
-  fOutputList->Add(fHEMCalRecdPhidEtaposEta);
-
-  fHEMCalRecdPhidEtanegEta = new TH2F("fHEMCalRecdPhidEtanegEta","(-eta track) EMCAL Cluster-Track #Delta#phi-#Delta#eta; #Delta#eta; #Delta#phi",1000,-0.1,0.1,1000,-0.5,0.5);
-  fOutputList->Add(fHEMCalRecdPhidEtanegEta);
-
-  fHPhotonEdiff100HC = new TH2F("fHPhotonEdiff100HC","Photon (E_{Truth}- E_{calc,100% HC})/E_{Truth} vs. E_{Truth}; E_{Truth} (GeV); (E_{Truth}- E_{calc,100% HC})/E_{Truth}",1000,0,10,600,-4.9,1.1);
-  fOutputList->Add(fHPhotonEdiff100HC);
-
-  fHPhotonEdiff70HC = new TH2F("fHPhotonEdiff70HC","Photon (E_{Truth}- E_{calc,70% HC})/E_{Truth} vs. E_{Truth}; E_{Truth} (GeV); (E_{Truth}- E_{calc,30% HC})/E_{Truth}",1000,0,10,600,-4.9,1.1);
-  fOutputList->Add(fHPhotonEdiff70HC);
-
-  fHPhotonEdiff30HC = new TH2F("fHPhotonEdiff30HC","Photon (E_{Truth}- E_{calc,30% HC})/E_{Truth} vs. E_{Truth}; E_{Truth} (GeV); (E_{Truth}- E_{calc,30% HC})/E_{Truth}",1000,0,10,600,-4.9,1.1);
-  fOutputList->Add(fHPhotonEdiff30HC);
-
-  fHPhotonEdiff0HC = new TH2F("fHPhotonEdiff0HC","Photon (E_{Truth}- E_{calc,0% HC})/E_{Truth} vs. E_{Truth}; E_{Truth} (GeV); (E_{Truth}- E_{cls})/E_{Truth}",1000,0,10,600,-4.9,1.1);
-  fOutputList->Add(fHPhotonEdiff0HC);
-
-  fHPhotonEVsClsE = new TH2F("fHPhotonEVsClsE","Cluster E vs. photon E_{Truth}; photon E_{Truth} (GeV); Cluster E (GeV)",500,0,5,500,0,5);
-  fOutputList->Add(fHPhotonEVsClsE);
-
-  fHistEsub1Pch =new  TH2F("fHistEsub1Pch", "(subtracted E in 100% HC) vs. total track P, clusters with 1 matching track; total track P (GeV/c); E_{sub}(GeV)" , 1000, 0., 10, 1000, 0., 10.);
-  fOutputList->Add(fHistEsub1Pch);
-
-  fHistEsub2Pch =new  TH2F("fHistEsub2Pch", "(subtracted E in 100% HC) vs. total track P, clusters with 2 matching tracks; total track P (GeV/c); E_{sub}(GeV)" , 1000, 0., 10, 1000, 0., 10.);
-  fOutputList->Add(fHistEsub2Pch);
-
-  fHistEsub1PchRat =new  TH2F("fHistEsub1PchRat", "(subtracted E in 100% HC)/total track P vs. total track P, clusters with 1 matching track; total track P (GeV/c); E_{sub}/P_{tot}" , 1000, 0., 10, 1100, 0., 1.1);
-  fOutputList->Add(fHistEsub1PchRat);
-
-  fHistEsub2PchRat =new  TH2F("fHistEsub2PchRat", "(subtracted E in 100% HC)/total track P vs. total track P, clusters with 2 matching tracks; total track P (GeV/c); E_{sub}/P_{tot}" , 1000, 0., 10, 1100, 0., 1.1);
-  fOutputList->Add(fHistEsub2PchRat);
-
-  Int_t bins[4] = {150, 150, 100, 200};
-  Double_t xmin[4] = {1.3, -0.8, 0, 0};
-  Double_t xmax[4] = {3.2, 0.8, 10, 2};
-
-  fHClsEoverMcE_All = new THnSparseF("fHClsEoverMcE_All", "Cluster E/MC E, clusters with 1 matching particle; #phi; #eta; E (GeV); ClsE/McE", 4, bins, xmin, xmax);
-  fOutputList->Add(fHClsEoverMcE_All);
+  if(fSFProcess)
+  {
+    fHScaleFactor = new TH1F("fHScaleFactor", "Scale factor distribution without hadronic correction;Scale factor",100,0,10);
+    fOutputList->Add(fHScaleFactor);
+    
+    fHScaleFactor100HC = new TH1F("fHScaleFactor100HC", "Scale factor distribution with 100% hadronic correction;Scale factor",100,0,10);
+    fOutputList->Add(fHScaleFactor100HC);
+  }
 
-  fHClsEoverMcE_Photon = new THnSparseF("fHClsEoverMcE_Photon", "Cluster E/MC E, clusters with 1 matching photon; #phi; #eta; E (GeV); ClsE/McE", 4, bins, xmin, xmax);
-  fOutputList->Add(fHClsEoverMcE_Photon);
+  if(fClusterProcess)
+  {
+    fHEOverPVsPt = new TH2F("fHEOverPVsPt", "E/P vs track p_{T}; p_{T} (GeV/c); E/P", 200 , 0, 4, 200, 0, 3.2);
+    fOutputList->Add(fHEOverPVsPt);
+    
+    fHEMCalResponsePion = new TH2F("fHEMCalResponsePion", "Pion E/P vs track p_{T}; p_{T} (GeV/c); E/P", 100 , 0, 4, 100, 0, 3.2);
+    fOutputList->Add(fHEMCalResponsePion);
+    
+    fHEMCalResponseElec = new TH2F("fHEMCalResponseElec", "Electron E/P vs track p_{T};  p_{T} (GeV/c); E/P", 100 , 0, 4, 100, 0, 3.2);
+    fOutputList->Add(fHEMCalResponseElec);
+    
+    fHEMCalResponseProton = new TH2F("fHEMCalResponseProton", "Proton E/P vs track p_{T};  p_{T} (GeV/c); E/P", 100 , 0, 4, 100, 0, 3.2);
+    fOutputList->Add(fHEMCalResponseProton);
+    
+    fHEMCalRecdPhidEta = new TH2F("fHEMCalRecdPhidEta","EMCAL Cluster-Track #Delta#phi-#Delta#eta; #Delta#eta; #Delta#phi",1000,-0.1,0.1,1000,-0.5,0.5);
+    fOutputList->Add(fHEMCalRecdPhidEta);
+    
+    fHEMCalRecdPhidEtaP = new TH2F("fHEMCalRecdPhidEtaP","EMCAL Charge+ Cluster-Track #Delta#phi-#Delta#eta; #Delta#eta; #Delta#phi",1000,-0.1,0.1,1000,-0.5,0.5);
+    fOutputList->Add(fHEMCalRecdPhidEtaP);
+    
+    fHEMCalRecdPhidEtaM = new TH2F("fHEMCalRecdPhidEtaM","EMCAL Charge- Cluster-Track #Delta#phi-#Delta#eta; #Delta#eta; #Delta#phi",1000,-0.1,0.1,1000,-0.5,0.5);
+    fOutputList->Add(fHEMCalRecdPhidEtaM);
+    
+    fHEMCalRecdPhidEta_Truth = new TH2F("fHEMCalRecdPhidEta_Truth","EMCAL Cluster-Track(Truth matched) #Delta#phi-#Delta#eta; #Delta#eta; #Delta#phi",1000,-0.1,0.1,1000,-0.5,0.5);
+    fOutputList->Add(fHEMCalRecdPhidEta_Truth);
+    
+    fHEMCalRecdPhidEtaP_Truth = new TH2F("fHEMCalRecdPhidEtaP_Truth","EMCAL Charge+ Cluster-Track(Truth matched) #Delta#phi-#Delta#eta; #Delta#eta; #Delta#phi",1000,-0.1,0.1,1000,-0.5,0.5);
+    fOutputList->Add(fHEMCalRecdPhidEtaP_Truth);
+    
+    fHEMCalRecdPhidEtaM_Truth = new TH2F("fHEMCalRecdPhidEtaM_Truth","EMCAL Charge- Cluster-Track(Truth matched) #Delta#phi-#Delta#eta; #Delta#eta; #Delta#phi",1000,-0.1,0.1,1000,-0.5,0.5);
+    fOutputList->Add(fHEMCalRecdPhidEtaM_Truth);
+    
+    fHEMCalRecdPhidEtaposEta = new TH2F("fHEMCalRecdPhidEtaposEta","(+eta track) EMCAL Cluster-Track #Delta#phi-#Delta#eta; #Delta#eta; #Delta#phi",1000,-0.1,0.1,1000,-0.5,0.5);
+    fOutputList->Add(fHEMCalRecdPhidEtaposEta);
+    
+    fHEMCalRecdPhidEtanegEta = new TH2F("fHEMCalRecdPhidEtanegEta","(-eta track) EMCAL Cluster-Track #Delta#phi-#Delta#eta; #Delta#eta; #Delta#phi",1000,-0.1,0.1,1000,-0.5,0.5);
+    fOutputList->Add(fHEMCalRecdPhidEtanegEta);
+    
+    fHPhotonEdiff100HC = new TH2F("fHPhotonEdiff100HC","Photon (E_{Truth}- E_{calc,100% HC})/E_{Truth} vs. E_{Truth}; E_{Truth} (GeV); (E_{Truth}- E_{calc,100% HC})/E_{Truth}",1000,0,10,600,-4.9,1.1);
+    fOutputList->Add(fHPhotonEdiff100HC);
+    
+    fHPhotonEdiff70HC = new TH2F("fHPhotonEdiff70HC","Photon (E_{Truth}- E_{calc,70% HC})/E_{Truth} vs. E_{Truth}; E_{Truth} (GeV); (E_{Truth}- E_{calc,30% HC})/E_{Truth}",1000,0,10,600,-4.9,1.1);
+    fOutputList->Add(fHPhotonEdiff70HC);
+    
+    fHPhotonEdiff30HC = new TH2F("fHPhotonEdiff30HC","Photon (E_{Truth}- E_{calc,30% HC})/E_{Truth} vs. E_{Truth}; E_{Truth} (GeV); (E_{Truth}- E_{calc,30% HC})/E_{Truth}",1000,0,10,600,-4.9,1.1);
+    fOutputList->Add(fHPhotonEdiff30HC);
+    
+    fHPhotonEdiff0HC = new TH2F("fHPhotonEdiff0HC","Photon (E_{Truth}- E_{calc,0% HC})/E_{Truth} vs. E_{Truth}; E_{Truth} (GeV); (E_{Truth}- E_{cls})/E_{Truth}",1000,0,10,600,-4.9,1.1);
+    fOutputList->Add(fHPhotonEdiff0HC);
+    
+    fHPhotonEVsClsE = new TH2F("fHPhotonEVsClsE","Cluster E vs. photon E_{Truth}; photon E_{Truth} (GeV); Cluster E (GeV)",500,0,5,500,0,5);
+    fOutputList->Add(fHPhotonEVsClsE);
+    
+    fHistEsub1Pch =new  TH2F("fHistEsub1Pch", "(subtracted E in 100% HC) vs. total track P, clusters with 1 matching track; total track P (GeV/c); E_{sub}(GeV)" , 1000, 0., 10, 1000, 0., 10.);
+    fOutputList->Add(fHistEsub1Pch);
+    
+    fHistEsub2Pch =new  TH2F("fHistEsub2Pch", "(subtracted E in 100% HC) vs. total track P, clusters with 2 matching tracks; total track P (GeV/c); E_{sub}(GeV)" , 1000, 0., 10, 1000, 0., 10.);
+    fOutputList->Add(fHistEsub2Pch);
+    
+    fHistEsub1PchRat =new  TH2F("fHistEsub1PchRat", "(subtracted E in 100% HC)/total track P vs. total track P, clusters with 1 matching track; total track P (GeV/c); E_{sub}/P_{tot}" , 1000, 0., 10, 1100, 0., 1.1);
+    fOutputList->Add(fHistEsub1PchRat);
+    
+    fHistEsub2PchRat =new  TH2F("fHistEsub2PchRat", "(subtracted E in 100% HC)/total track P vs. total track P, clusters with 2 matching tracks; total track P (GeV/c); E_{sub}/P_{tot}" , 1000, 0., 10, 1100, 0., 1.1);
+    fOutputList->Add(fHistEsub2PchRat);
+    
+    Int_t bins[4] = {150, 150, 100, 200};
+    Double_t xmin[4] = {1.3, -0.8, 0, 0};
+    Double_t xmax[4] = {3.2, 0.8, 10, 2};
+    
+    fHClsEoverMcE_All = new THnSparseF("fHClsEoverMcE_All", "Cluster E/MC E, clusters with 1 matching particle; #phi; #eta; E (GeV); ClsE/McE", 4, bins, xmin, xmax);
+    fOutputList->Add(fHClsEoverMcE_All);
+    
+    fHClsEoverMcE_Photon = new THnSparseF("fHClsEoverMcE_Photon", "Cluster E/MC E, clusters with 1 matching photon; #phi; #eta; E (GeV); ClsE/McE", 4, bins, xmin, xmax);
+    fOutputList->Add(fHClsEoverMcE_Photon);
+    
+    fHClsEoverMcE_Elec = new THnSparseF("fHClsEoverMcE_Elec", "Cluster E/MC E, clusters with 1 matching electron; #phi; #eta; E (GeV); ClsE/McE", 4, bins, xmin, xmax);
+    fOutputList->Add(fHClsEoverMcE_Elec);
+    
+    fHClsEoverMcE_Pion = new THnSparseF("fHClsEoverMcE_Pion", "Cluster E/MC E, clusters with 1 matching pion; #phi; #eta; E (GeV); ClsE/McE",4, bins, xmin, xmax);
+    fOutputList->Add(fHClsEoverMcE_Pion);
+  }
+  
+  fHParGenPion_p = new TH3F("fHParGenPion_p","Particle level truth Phi-Eta-p_{T} distribution of #pi+",  500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
+  fHParGenPion_p->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
+  fHParGenPion_p->GetYaxis()->SetTitle("#eta");
+  fHParGenPion_p->GetZaxis()->SetTitle("#phi");
+  fOutputList->Add(fHParGenPion_p);
+  
+  fHParGenPion_m = new TH3F("fHParGenPion_m", "Particle level truth Phi-Eta-p_{T} distribution of all #pi-", 500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
+  fHParGenPion_m->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
+  fHParGenPion_m->GetYaxis()->SetTitle("#eta");
+  fHParGenPion_m->GetZaxis()->SetTitle("#phi");
+  fOutputList->Add(fHParGenPion_m);
+  
+  fHParGenPion_rmInj_p = new TH3F("fHParGenPion_rmInj_p","Particle level truth Phi-Eta-p_{T} distribution of all #pi+ without injected signal",  500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
+  fHParGenPion_rmInj_p->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
+  fHParGenPion_rmInj_p->GetYaxis()->SetTitle("#eta");
+  fHParGenPion_rmInj_p->GetZaxis()->SetTitle("#phi");
+  fOutputList->Add(fHParGenPion_rmInj_p);
+  
+  fHParGenPion_rmInj_m = new TH3F("fHParGenPion_rmInj_m","Particle level truth Phi-Eta-p_{T} distribution of #pi- without injected signal", 500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
+  fHParGenPion_rmInj_m->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
+  fHParGenPion_rmInj_m->GetYaxis()->SetTitle("#eta");
+  fHParGenPion_rmInj_m->GetZaxis()->SetTitle("#phi");
+  fOutputList->Add(fHParGenPion_rmInj_m);
 
-  fHClsEoverMcE_Elec = new THnSparseF("fHClsEoverMcE_Elec", "Cluster E/MC E, clusters with 1 matching electron; #phi; #eta; E (GeV); ClsE/McE", 4, bins, xmin, xmax);
-  fOutputList->Add(fHClsEoverMcE_Elec);
 
-  fHClsEoverMcE_Pion = new THnSparseF("fHClsEoverMcE_Pion", "Cluster E/MC E, clusters with 1 matching pion; #phi; #eta; E (GeV); ClsE/McE",4, bins, xmin, xmax);
-  fOutputList->Add(fHClsEoverMcE_Pion);
+  
+  const char* trackCut[3] = {"cut1","cut2", "cut3"};
+  if(fMcProcess && fTrackProcess)
+  {    
+    //Fake
+    fHDetGenFakePion = new TH3F("fHDetGenFakePion", "fake charged pion track Phi-Eta-p_{T} distribution",500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
+    fHDetGenFakePion->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
+    fHDetGenFakePion->GetYaxis()->SetTitle("#eta");
+    fHDetGenFakePion->GetZaxis()->SetTitle("#phi");
+    fOutputList->Add(fHDetGenFakePion);
+
+    fHDetRecFakePion = new TH3F("fHDetRecFakePion", "fake charged pion track Phi-Eta-p_{T} distribution", 500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
+    fHDetRecFakePion->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
+    fHDetRecFakePion->GetYaxis()->SetTitle("#eta");
+    fHDetRecFakePion->GetZaxis()->SetTitle("#phi");
+    fOutputList->Add(fHDetRecFakePion);
+    
+    //Secondary
+    fHDetGenSecPion = new TH3F("fHDetGenSecPion", "secondary charged pion charged pion track Phi-Eta-p_{T} distribution", 500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
+    fHDetGenSecPion->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
+    fHDetGenSecPion->GetYaxis()->SetTitle("#eta");
+    fHDetGenSecPion->GetZaxis()->SetTitle("#phi");
+    fOutputList->Add(fHDetGenSecPion);
+
+    fHDetRecSecPion = new TH3F("fHDetRecSecPion", "secondary charged pion charged pion track Phi-Eta-p_{T} distribution", 500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
+    fHDetRecSecPion->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
+    fHDetRecSecPion->GetYaxis()->SetTitle("#eta");
+    fHDetRecSecPion->GetZaxis()->SetTitle("#phi");
+    fOutputList->Add(fHDetRecSecPion);
+    for(Int_t i=0; i<3; i++)
+    {
+      // pi+
+      fHDetGenPion_p[i] = new TH3F(Form("fHDetGenPion_p_%s", trackCut[i]), Form("%s: Detector level truth Phi-Eta-p_{T} distribution of #pi+", trackCut[i]),  500, 0, 100, 60, -1.2, 1.2, 128, 0, 6.4);
+      fHDetGenPion_p[i]->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
+      fHDetGenPion_p[i]->GetYaxis()->SetTitle("#eta");
+      fHDetGenPion_p[i]->GetZaxis()->SetTitle("#phi");
+      fOutputList->Add(fHDetGenPion_p[i]);
+      
+      fHDetRecPion_p[i] = new TH3F(Form("fHDetRecPion_p_%s", trackCut[i]), Form("%s: Reconstructed track Phi-Eta-p_{T} distribution of all #pi+", trackCut[i]),  500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
+      fHDetRecPion_p[i]->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
+      fHDetRecPion_p[i]->GetYaxis()->SetTitle("#eta");
+      fHDetRecPion_p[i]->GetZaxis()->SetTitle("#phi");
+      fOutputList->Add(fHDetRecPion_p[i]);
+
+      // pi-
+      fHDetGenPion_m[i] = new TH3F(Form("fHDetGenPion_m_%s", trackCut[i]), Form("%s: Detector level truth Phi-Eta-p_{T} distribution of #pi-", trackCut[i]),  500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
+      fHDetGenPion_m[i]->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
+      fHDetGenPion_m[i]->GetYaxis()->SetTitle("#eta");
+      fHDetGenPion_m[i]->GetZaxis()->SetTitle("#phi");
+      fOutputList->Add(fHDetGenPion_m[i]);
+      
+      fHDetRecPion_m[i] = new TH3F(Form("fHDetRecPion_m_%s", trackCut[i]), Form("%s: Reconstructed track Phi-Eta-p_{T} distribution of #pi-", trackCut[i]),  500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
+      fHDetRecPion_m[i]->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
+      fHDetRecPion_m[i]->GetYaxis()->SetTitle("#eta");
+      fHDetRecPion_m[i]->GetZaxis()->SetTitle("#phi");
+      fOutputList->Add(fHDetRecPion_m[i]);
+
+        //pi+ without injected signal
+      fHDetGenPion_rmInj_p[i] = new TH3F(Form("fHDetGenPion_rmInj_p_%s", trackCut[i]), Form("%s: Detector level truth Phi-Eta-p_{T} distribution of #pi+ without injected signal", trackCut[i]), 500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
+      fHDetGenPion_rmInj_p[i]->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
+      fHDetGenPion_rmInj_p[i]->GetYaxis()->SetTitle("#eta");
+      fHDetGenPion_rmInj_p[i]->GetZaxis()->SetTitle("#phi");
+      fOutputList->Add(fHDetGenPion_rmInj_p[i]);
+      
+      fHDetRecPion_rmInj_p[i] = new TH3F(Form("fHDetRecPion_rmInj_p_%s", trackCut[i]), Form("%s: Reconstructed track Phi-Eta-p_{T} distribution of #pi+ without injected signal", trackCut[i]), 500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
+      fHDetRecPion_rmInj_p[i]->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
+      fHDetRecPion_rmInj_p[i]->GetYaxis()->SetTitle("#eta");
+      fHDetRecPion_rmInj_p[i]->GetZaxis()->SetTitle("#phi");
+      fOutputList->Add(fHDetRecPion_rmInj_p[i]);
+
+      //pi- charged particle without injected signal
+      fHDetGenPion_rmInj_m[i] = new TH3F(Form("fHDetGenPion_rmInj_m_%s", trackCut[i]), Form("%s: Detector level truth Phi-Eta-p_{T} distribution of #pi- without injected signal", trackCut[i]), 500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
+      fHDetGenPion_rmInj_m[i]->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
+      fHDetGenPion_rmInj_m[i]->GetYaxis()->SetTitle("#eta");
+      fHDetGenPion_rmInj_m[i]->GetZaxis()->SetTitle("#phi");
+      fOutputList->Add(fHDetGenPion_rmInj_m[i]);
+      
+      fHDetRecPion_rmInj_m[i] = new TH3F(Form("fHDetRecPion_rmInj_m_%s", trackCut[i]), Form("%s: Reconstructed track Phi-Eta-p_{T} distribution of #pi- without injected signal", trackCut[i]), 500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
+      fHDetRecPion_rmInj_m[i]->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
+      fHDetRecPion_rmInj_m[i]->GetYaxis()->SetTitle("#eta");
+      fHDetRecPion_rmInj_m[i]->GetZaxis()->SetTitle("#phi");
+      fOutputList->Add(fHDetRecPion_rmInj_m[i]);
+    }
+  }    
   
   fTrackIndices = new TArrayI();
   fClusterIndices = new TArrayI();
-
+  
   fClusterArray = new TObjArray();
   fClusterArray->SetOwner(1);
-
+  
   PostData(1, fOutputList);
 }
 
@@ -375,6 +431,7 @@ void AliAnalysisTaskSOH::UserExec(Option_t *)
     printf("ERROR: fESD not available\n");
     return;
   }
+  if (!EsdVertexOk())  return;  // Vetex cut
 
   fMC = MCEvent();
   if (!fMC) {
@@ -391,10 +448,14 @@ void AliAnalysisTaskSOH::UserExec(Option_t *)
   if(fClusterArray)
     fClusterArray->Delete();
 
-  ProcessTrack();
-  ProcessCluster();
-  ProcessMc();
-  ProcessScaleFactor();
+  if(fTrackProcess) 
+    ProcessTrack();
+  if(fClusterProcess)
+    ProcessCluster();
+  if(fMcProcess)
+    ProcessMc();
+  if(fSFProcess)
+    ProcessScaleFactor();
   
   PostData(1, fOutputList);
 }   
@@ -403,7 +464,7 @@ void AliAnalysisTaskSOH::UserExec(Option_t *)
 void  AliAnalysisTaskSOH::ProcessTrack()
 {
   // Process track.
-
+  
   fTrackIndices->Set(fESD->GetNumberOfTracks());
   AliDebug(3,Form("%s:%d Selecting tracks",(char*)__FILE__,__LINE__));
 
@@ -419,79 +480,81 @@ void  AliAnalysisTaskSOH::ProcessTrack()
     if(!esdtrack)continue;
     AliESDtrack *newTrack = GetAcceptTrack(esdtrack);
     if(!newTrack) continue;
-  
-    Double_t clsE = -1;
-    Int_t clsIndex = newTrack->GetEMCALcluster();
-    if(newTrack->GetEMCALcluster()>-1)
-    {
-      AliESDCaloCluster *cluster = fESD->GetCaloCluster(clsIndex);
-      if(IsGoodCluster(cluster))
+    if(newTrack->Pt()<0.15 || TMath::Abs(newTrack->Eta())>0.9) {delete newTrack; continue;}
+    
+    if(fClusterProcess)
+    {  
+      Double_t clsE = -1;
+      Int_t clsIndex = newTrack->GetEMCALcluster();
+      if(newTrack->GetEMCALcluster()>-1)
       {
-        isMth=1;
-             
-        cluster->GetPosition(ClsPos);
-        TVector3 VClsPos(ClsPos[0], ClsPos[1], ClsPos[2]);
-             
-        AliEMCALTrack EMCTrk(*newTrack);
-        if(!EMCTrk.PropagateToGlobal(ClsPos[0], ClsPos[1], ClsPos[2], 0.0, 0.0)) {
-          continue;
-        }
-        EMCTrk.GetXYZ(emcTrkpos);
-        TVector3 VemcTrkPos(emcTrkpos[0],emcTrkpos[1],emcTrkpos[2]);
-             
-        Double_t dPhi = VClsPos.Phi() - VemcTrkPos.Phi();
-        if (dPhi < -1*TMath::Pi()) dPhi += (2*TMath::Pi());
-        else if (dPhi > TMath::Pi()) dPhi -= (2*TMath::Pi());
-             
-        Double_t dEta = VClsPos.Eta() - VemcTrkPos.Eta();
-             
-        fHEMCalRecdPhidEta->Fill(dEta, dPhi);  
-
-        if((newTrack->GetLabel())>-1 && (newTrack->GetLabel()) < fMC->GetNumberOfTracks())
+       AliESDCaloCluster *cluster = fESD->GetCaloCluster(clsIndex);
+       if(IsGoodCluster(cluster))
         {
-          AliVParticle *vParticle = fMC->GetTrack(newTrack->GetLabel());
-          if(IsGoodMcParticle(vParticle, newTrack->GetLabel()))
-          {
-            fHEMCalRecdPhidEta_Truth->Fill(dEta, dPhi);
-            if(vParticle->Charge() > 0) fHEMCalRecdPhidEtaP_Truth->Fill(dEta, dPhi);
-            if(vParticle->Charge() < 0) fHEMCalRecdPhidEtaM_Truth->Fill(dEta, dPhi);
-          }
-        }
-
-        if(esdtrack->Charge() > 0) {fHEMCalRecdPhidEtaP->Fill(dEta, dPhi);}
-        if(esdtrack->Charge() < 0) {fHEMCalRecdPhidEtaM->Fill(dEta, dPhi);}
-
-        if(VemcTrkPos.Eta() > 0) fHEMCalRecdPhidEtaposEta->Fill(dEta, dPhi);
-        if(VemcTrkPos.Eta() < 0) fHEMCalRecdPhidEtanegEta->Fill(dEta, dPhi);
-             
-        clsE = cluster->E();
-        if(newTrack->P()>0) fHEOverPVsPt->Fill(newTrack->Pt(),clsE/newTrack->P());
-      }
+         isMth=1;
          
-      Int_t ipart = newTrack->GetLabel();
-      if(ipart>-1 && ipart<fMC->GetNumberOfTracks())
-      {
-        AliVParticle* vParticle = fMC->GetTrack(ipart);
-        if(isMth && vParticle)
-        {
-          if(TMath::Abs(vParticle->PdgCode())==211)
-          {
-            fHEMCalResponsePion->Fill(newTrack->Pt(),clsE/newTrack->P());
-          }
-          if(TMath::Abs(vParticle->PdgCode())==11)
-          {
-            fHEMCalResponseElec->Fill(newTrack->Pt(),clsE/newTrack->P());
-          }
-          if(TMath::Abs(vParticle->PdgCode())==2212)
-          {
-            fHEMCalResponseProton->Fill(newTrack->Pt(),clsE/newTrack->P());
-          }
-        }
+         cluster->GetPosition(ClsPos);
+         TVector3 VClsPos(ClsPos[0], ClsPos[1], ClsPos[2]);
+         
+         AliEMCALTrack EMCTrk(*newTrack);
+         if(!EMCTrk.PropagateToGlobal(ClsPos[0], ClsPos[1], ClsPos[2], 0.0, 0.0)) {continue;}
+         EMCTrk.GetXYZ(emcTrkpos);
+         TVector3 VemcTrkPos(emcTrkpos[0],emcTrkpos[1],emcTrkpos[2]);
+         
+         Double_t dPhi = VClsPos.Phi() - VemcTrkPos.Phi();
+         if (dPhi < -1*TMath::Pi()) dPhi += (2*TMath::Pi());
+         else if (dPhi > TMath::Pi()) dPhi -= (2*TMath::Pi());
+         
+         Double_t dEta = VClsPos.Eta() - VemcTrkPos.Eta();
+         
+         fHEMCalRecdPhidEta->Fill(dEta, dPhi); 
+         
+         if((newTrack->GetLabel())>-1 && (newTrack->GetLabel()) < fMC->GetNumberOfTracks())
+         {
+           AliVParticle *vParticle = fMC->GetTrack(newTrack->GetLabel());
+           if(IsGoodMcParticle(vParticle, newTrack->GetLabel()))
+           {
+             fHEMCalRecdPhidEta_Truth->Fill(dEta, dPhi);
+             if(vParticle->Charge() > 0) fHEMCalRecdPhidEtaP_Truth->Fill(dEta, dPhi);
+             if(vParticle->Charge() < 0) fHEMCalRecdPhidEtaM_Truth->Fill(dEta, dPhi);
+           }
+         }
+         
+         if(esdtrack->Charge() > 0) {fHEMCalRecdPhidEtaP->Fill(dEta, dPhi);}
+         if(esdtrack->Charge() < 0) {fHEMCalRecdPhidEtaM->Fill(dEta, dPhi);}
+         
+         if(VemcTrkPos.Eta() > 0) fHEMCalRecdPhidEtaposEta->Fill(dEta, dPhi);
+         if(VemcTrkPos.Eta() < 0) fHEMCalRecdPhidEtanegEta->Fill(dEta, dPhi);
+         
+         clsE = cluster->E();
+         if(newTrack->P()>0) fHEOverPVsPt->Fill(newTrack->Pt(),clsE/newTrack->P());
+       }
+       
+       Int_t ipart = newTrack->GetLabel();
+       if(ipart>-1 && ipart<fMC->GetNumberOfTracks())
+       {
+         AliVParticle* vParticle = fMC->GetTrack(ipart);
+         if(isMth && vParticle)
+         {
+           if(TMath::Abs(vParticle->PdgCode())==211)
+           {
+             fHEMCalResponsePion->Fill(newTrack->Pt(),clsE/newTrack->P());
+           }
+           if(TMath::Abs(vParticle->PdgCode())==11)
+           {
+             fHEMCalResponseElec->Fill(newTrack->Pt(),clsE/newTrack->P());
+           }
+           if(TMath::Abs(vParticle->PdgCode())==2212)
+           {
+             fHEMCalResponseProton->Fill(newTrack->Pt(),clsE/newTrack->P());
+           }
+         }
+       }
       }
     }
-    if(newTrack->GetLabel()<0 && newTrack->GetPID()==2) fHTrkEffDetRecFakePtEtaPhi->Fill(esdtrack->Pt(), esdtrack->Eta(), esdtrack->Phi());
+    if(newTrack) delete newTrack;
 
- // Track Indices
   // Track Indices
     fTrackIndices->AddAt(itr,nTracks);
     nTracks++;
   }
@@ -528,7 +591,7 @@ void AliAnalysisTaskSOH::ProcessCluster()
     
     if(TrackLabels->GetSize() == 1)
     {
-      AliESDtrack *esdtrack = fESD->GetTrack(TrackLabels->GetAt(0));
+      AliESDtrack *esdtrack = fESD->GetTrack(TrackLabels->operator[](0));
       AliESDtrack *newTrack = GetAcceptTrack(esdtrack);
       if(newTrack && TMath::Abs(newTrack->Eta())<0.7)
       {
@@ -541,8 +604,8 @@ void AliAnalysisTaskSOH::ProcessCluster()
 
     if(TrackLabels->GetSize() == 2)
     {
-      AliESDtrack *esdtrack1 = fESD->GetTrack(TrackLabels->GetAt(0));
-      AliESDtrack *esdtrack2 = fESD->GetTrack(TrackLabels->GetAt(1));
+      AliESDtrack *esdtrack1 = fESD->GetTrack(TrackLabels->operator[](0));
+      AliESDtrack *esdtrack2 = fESD->GetTrack(TrackLabels->operator[](1));
       AliESDtrack *newTrack1 = GetAcceptTrack(esdtrack1);
       AliESDtrack *newTrack2 = GetAcceptTrack(esdtrack2);
       if(newTrack1 && newTrack2 && TMath::Abs(newTrack1->Eta())<0.7  && TMath::Abs(newTrack2->Eta())<0.7)
@@ -575,8 +638,8 @@ void AliAnalysisTaskSOH::ProcessCluster()
     if(MCLabels->GetSize() == 1) 
     {
       fHEventStat->Fill(4.5);
-      AliVParticle* vParticle1 = fMC->GetTrack(MCLabels->GetAt(0));
-      if(IsGoodMcParticle(vParticle1, MCLabels->GetAt(0)))
+      AliVParticle* vParticle1 = fMC->GetTrack(MCLabels->operator[](0));
+      if(IsGoodMcParticle(vParticle1, MCLabels->operator[](0)))
       {
        Double_t Entries[4] = {VClsPos.Phi(), VClsPos.Eta(), vParticle1->E(), cluster->E()/vParticle1->E()}; 
        fHClsEoverMcE_All->Fill(Entries);
@@ -597,9 +660,9 @@ void AliAnalysisTaskSOH::ProcessCluster()
     if(MCLabels->GetSize() == 2) 
     {
       fHEventStat->Fill(5.5);
-      AliVParticle* vParticle1 = fMC->GetTrack(MCLabels->GetAt(0));
-      AliVParticle* vParticle2 = fMC->GetTrack(MCLabels->GetAt(1));
-      if(IsGoodMcParticle(vParticle1, MCLabels->GetAt(0)) && IsGoodMcParticle(vParticle2, MCLabels->GetAt(1))) 
+      AliVParticle* vParticle1 = fMC->GetTrack(MCLabels->operator[](0));
+      AliVParticle* vParticle2 = fMC->GetTrack(MCLabels->operator[](1));
+      if(IsGoodMcParticle(vParticle1, MCLabels->operator[](0)) && IsGoodMcParticle(vParticle2, MCLabels->operator[](1))) 
       {
        fHEventStat->Fill(6.5);
        if((vParticle1->PdgCode()==22) && (vParticle2->PdgCode()==22)) {;}
@@ -646,14 +709,33 @@ void AliAnalysisTaskSOH::ProcessCluster()
 void AliAnalysisTaskSOH::ProcessMc()
 {
   // Process MC.
-
-  for(Int_t i=0; i<fTrackIndices->GetSize(); i++)
+   for(Int_t i=0; i<fTrackIndices->GetSize(); i++)
   {
     AliESDtrack *esdtrack = fESD->GetTrack(fTrackIndices->At(i));
+    if(!esdtrack)continue;
+    AliESDtrack *newTrack = GetAcceptTrack(esdtrack);
+    if(!newTrack) continue;
+    if(newTrack->Pt()<0.15 || TMath::Abs(newTrack->Eta())>0.9) {delete newTrack; continue;}
+
     Int_t index = esdtrack->GetLabel();
-    if(index < 0 || index > fESD->GetNumberOfTracks()) continue;
-    AliVParticle* vParticle = fMC->GetTrack(index);
-    if(!IsGoodMcParticle(vParticle, index) && esdtrack->GetPID() == 2) fHTrkEffDetRecFakePtEtaPhi->Fill(esdtrack->Pt(), esdtrack->Eta(), esdtrack->Phi());
+    if(index < 0) 
+    {
+      AliVParticle *vParticle1  = (AliVParticle*)fMC->GetTrack(-1*index);
+      if((TMath::Abs(vParticle1->PdgCode())==211) && IsGoodMcParticle(vParticle1, -1*index)) 
+      {
+       fHDetRecFakePion->Fill(newTrack->Pt(), newTrack->Eta(), newTrack->Phi());
+       fHDetGenFakePion->Fill(vParticle1->Pt(), vParticle1->Eta(), vParticle1->Phi());
+      }
+    }
+    
+    AliVParticle* vParticle2 = fMC->GetTrack(TMath::Abs(index));
+    if(!IsGoodMcParticle(vParticle2, TMath::Abs(index)) && (TMath::Abs(vParticle2->PdgCode())==211)) 
+    {
+      fHDetRecSecPion->Fill(newTrack->Pt(), newTrack->Eta(), newTrack->Phi());
+      fHDetGenSecPion->Fill(vParticle2->Pt(), vParticle2->Eta(), vParticle2->Phi());
+    }
+
+    if(newTrack) delete newTrack;
   }
 
   //tracking effciency
@@ -677,102 +759,107 @@ void AliAnalysisTaskSOH::ProcessMc()
     Int_t pdgCode = vParticle->PdgCode();
     AliMCParticle *McParticle  = (AliMCParticle*)fMC->GetTrack(ipart);
     if(!IsGoodMcParticle(vParticle, ipart)) continue;  
-
-    if((McParticle->GetMother() > eventHeader->NProduced()-1) && TMath::Abs(vParticle->Eta())<0.9 && TMath::Abs(pdgCode)==211)
-    {
-      AliMCParticle *McMother = (AliMCParticle*)fMC->GetTrack(McParticle->GetMother());
-      if(TMath::Abs(McMother->PdgCode()) == 411) fHTrkEffParGenPt_Dch->Fill(vParticle->Pt());
-      else if(TMath::Abs(McMother->PdgCode()) == 421) fHTrkEffParGenPt_Dn->Fill(vParticle->Pt());
-      else if(TMath::Abs(McMother->PdgCode()) == 431) fHTrkEffParGenPt_Ds->Fill(vParticle->Pt());       
-      else if(TMath::Abs(McMother->PdgCode()) == 4122) fHTrkEffParGenPt_cL->Fill(vParticle->Pt()); 
-      else if(TMath::Abs(McMother->PdgCode()) == 443) fHTrkEffParGenPt_JPsi->Fill(vParticle->Pt()); 
-      else {;}
-    }
-
-    if((McParticle->GetMother() < eventHeader->NProduced()) && TMath::Abs(vParticle->Eta())<0.9) 
+    if(TMath::Abs(vParticle->Eta())<0.9)
     {
-      fHTrkEffParGenPt_rmInj->Fill(vParticle->Pt());
-      if(TMath::Abs(pdgCode)==211) fHTrkEffParGenPt_PiRmInj->Fill(vParticle->Pt());
+      if(pdgCode==211) 
+      {
+       fHParGenPion_p->Fill(vParticle->Pt(), vParticle->Eta(), vParticle->Phi());
+       if(McParticle->GetMother() < eventHeader->NProduced()) fHParGenPion_rmInj_p->Fill(vParticle->Pt(), vParticle->Eta(), vParticle->Phi()); 
+      }
+    
+      else if(pdgCode==-211) 
+      {
+       fHParGenPion_m->Fill(vParticle->Pt(), vParticle->Eta(), vParticle->Phi());
+       if(McParticle->GetMother() < eventHeader->NProduced()) fHParGenPion_rmInj_m->Fill(vParticle->Pt(), vParticle->Eta(), vParticle->Phi()); 
+      }
     }
-
-    if(TMath::Abs(vParticle->Eta())<0.9)
+      
+    for(Int_t j=0; j<fTrackIndices->GetSize(); j++)
     {
-      fHTrkEffParGenPt_all->Fill(vParticle->Pt());
-      if(TMath::Abs(pdgCode)==211) fHTrkEffParGenPtEtaPhi->Fill(vParticle->Pt(), vParticle->Eta(), vParticle->Phi());
+      AliESDtrack *esdtrack = fESD->GetTrack(fTrackIndices->At(j));
+      if(!esdtrack)continue;
+      AliESDtrack *newTrack = GetAcceptTrack(esdtrack);
+      if(!newTrack) continue;
+      if(newTrack->Pt()<0.15 || TMath::Abs(newTrack->Eta())>0.9) {delete newTrack; continue;}
       
-      for(Int_t j=0; j<fTrackIndices->GetSize(); j++)
+      Int_t cutType = (Int_t)newTrack->GetTRDQuality();
+      
+      if(newTrack->GetLabel()==ipart && TMath::Abs(vParticle->Eta())<0.9)
       {
-        AliESDtrack *esdtrack = fESD->GetTrack(fTrackIndices->At(j));
-        if(esdtrack && (TMath::Abs(esdtrack->GetLabel())==ipart))
-        {
-         fHTrkEffDetGenPt_all->Fill(vParticle->Pt());
-         if(TMath::Abs(pdgCode)==211)
-         {
-           fHTrkEffDetGenPtEtaPhi->Fill(vParticle->Pt(), vParticle->Eta(), vParticle->Phi());
-           fHTrkEffDetRecPtEtaPhi->Fill(esdtrack->Pt(), esdtrack->Eta(), esdtrack->Phi());
-         }
-         if((McParticle->GetMother() > eventHeader->NProduced()-1) && TMath::Abs(pdgCode) == 211)
+       if(pdgCode==211) 
+       {
+         fHDetGenPion_p[cutType]->Fill(vParticle->Pt(), vParticle->Eta(), vParticle->Phi());
+         fHDetRecPion_p[cutType]->Fill(newTrack->Pt(), newTrack->Eta(), newTrack->Phi());
+         if(McParticle->GetMother() < eventHeader->NProduced())
          {
-           AliMCParticle *McMother = (AliMCParticle*)fMC->GetTrack(McParticle->GetMother());
-           if(TMath::Abs(McMother->PdgCode()) == 411) fHTrkEffDetGenPt_Dch->Fill(vParticle->Pt());
-           else if(TMath::Abs(McMother->PdgCode()) == 421) fHTrkEffDetGenPt_Dn->Fill(vParticle->Pt());
-           else if(TMath::Abs(McMother->PdgCode()) == 431) fHTrkEffDetGenPt_Ds->Fill(vParticle->Pt());  
-           else if(TMath::Abs(McMother->PdgCode()) == 4122) fHTrkEffDetGenPt_cL->Fill(vParticle->Pt()); 
-           else if(TMath::Abs(McMother->PdgCode()) == 443) fHTrkEffDetGenPt_JPsi->Fill(vParticle->Pt()); 
-           else {;}
-         }
+           fHDetGenPion_rmInj_p[cutType]->Fill(vParticle->Pt(), vParticle->Eta(), vParticle->Phi());
+           fHDetRecPion_rmInj_p[cutType]->Fill(newTrack->Pt(), newTrack->Eta(), newTrack->Phi());
+         }       
+       }
+       else if(pdgCode==-211) 
+       {
+         fHDetGenPion_m[cutType]->Fill(vParticle->Pt(), vParticle->Eta(), vParticle->Phi());
+         fHDetRecPion_m[cutType]->Fill(newTrack->Pt(), newTrack->Eta(), newTrack->Phi());
          if(McParticle->GetMother() < eventHeader->NProduced())
          {
-           fHTrkEffDetGenPt_rmInj->Fill(vParticle->Pt());
-           if(TMath::Abs(pdgCode)==211) fHTrkEffDetGenPt_PiRmInj->Fill(vParticle->Pt());
-         }
+           fHDetGenPion_rmInj_m[cutType]->Fill(vParticle->Pt(), vParticle->Eta(), vParticle->Phi());
+           fHDetRecPion_rmInj_m[cutType]->Fill(newTrack->Pt(), newTrack->Eta(), newTrack->Phi());
+         }       
+       }
+         
 
     //cluster E vs. truth photon energy
-         for(Int_t k=0; k<fClusterIndices->GetSize(); k++)
-         {
-           AliESDCaloCluster *cluster = fESD->GetCaloCluster(fClusterIndices->At(k));
-           Double_t clsE = cluster->E();
-           TArrayI *MCLabels = cluster->GetLabelsArray();
-           AliVParticle* vParticle1 = fMC->GetTrack(MCLabels->GetAt(0));
-           AliVParticle* vParticle2 = fMC->GetTrack(MCLabels->GetAt(1));
-           
-           if(vParticle1->PdgCode()==22 && vParticle2 == vParticle)
-           {
-             fHPhotonEdiff0HC->Fill(vParticle1->E(), (vParticle1->E() - clsE)/vParticle1->E());
-             fHPhotonEVsClsE->Fill(vParticle1->E(), clsE);
-
-             if((clsE - 0.3*esdtrack->E())<0) fHPhotonEdiff30HC->Fill(vParticle1->E(), 1);
-             else  fHPhotonEdiff30HC->Fill(vParticle1->E(), (vParticle1->E() + 0.3*esdtrack->E() - clsE)/vParticle1->E());
-
-             if((clsE - 0.7*esdtrack->E())<0) fHPhotonEdiff70HC->Fill(vParticle1->E(), 1);
-             else  fHPhotonEdiff70HC->Fill(vParticle1->E(), (vParticle1->E() + 0.7*esdtrack->E() - clsE)/vParticle1->E());
-
-             if((clsE - esdtrack->E())<0) fHPhotonEdiff100HC->Fill(vParticle1->E(), 1);
-             else  fHPhotonEdiff100HC->Fill(vParticle1->E(), (vParticle1->E() + esdtrack->E() - clsE)/vParticle1->E());
-             continue;
-           }
-           if(vParticle2->PdgCode()==22 && vParticle1 == vParticle)
-           {
-             fHPhotonEdiff0HC->Fill(vParticle2->E(), (vParticle2->E() - clsE)/vParticle2->E());
-             fHPhotonEVsClsE->Fill(vParticle2->E(), clsE);
-
-             if((clsE - 0.3*esdtrack->E())<0) fHPhotonEdiff30HC->Fill(vParticle2->E(), 1);
-             else fHPhotonEdiff30HC->Fill(vParticle2->E(), (vParticle2->E() + 0.3*esdtrack->E() - clsE)/vParticle2->E());
-
-             if((clsE - 0.7*esdtrack->E())<0) fHPhotonEdiff70HC->Fill(vParticle2->E(), 1);
-             else fHPhotonEdiff70HC->Fill(vParticle2->E(), (vParticle2->E() + 0.7*esdtrack->E() - clsE)/vParticle2->E());
-
-             if((clsE-esdtrack->E())<0) fHPhotonEdiff100HC->Fill(vParticle2->E(), 1);
-             else fHPhotonEdiff100HC->Fill(vParticle2->E(), (vParticle2->E() + esdtrack->E() - clsE)/vParticle2->E());
-           }
-         }
-          break;
-        }
+       if(fClusterProcess)
+        {
+          for(Int_t k=0; k<fClusterIndices->GetSize(); k++)
+          {
+            AliESDCaloCluster *cluster = fESD->GetCaloCluster(fClusterIndices->At(k));
+            Double_t clsE = cluster->E();
+            TArrayI *MCLabels = cluster->GetLabelsArray();
+            AliVParticle* vParticle1 = fMC->GetTrack(MCLabels->operator[](0));
+            AliVParticle* vParticle2 = fMC->GetTrack(MCLabels->operator[](1));
+            
+            if(vParticle1->PdgCode()==22 && vParticle2 == vParticle)
+            {
+              fHPhotonEdiff0HC->Fill(vParticle1->E(), (vParticle1->E() - clsE)/vParticle1->E());
+              fHPhotonEVsClsE->Fill(vParticle1->E(), clsE);
+              
+              if((clsE - 0.3*esdtrack->E())<0) fHPhotonEdiff30HC->Fill(vParticle1->E(), 1);
+              else  fHPhotonEdiff30HC->Fill(vParticle1->E(), (vParticle1->E() + 0.3*esdtrack->E() - clsE)/vParticle1->E());
+              
+              if((clsE - 0.7*esdtrack->E())<0) fHPhotonEdiff70HC->Fill(vParticle1->E(), 1);
+              else  fHPhotonEdiff70HC->Fill(vParticle1->E(), (vParticle1->E() + 0.7*esdtrack->E() - clsE)/vParticle1->E());
+              
+              if((clsE - esdtrack->E())<0) fHPhotonEdiff100HC->Fill(vParticle1->E(), 1);
+              else  fHPhotonEdiff100HC->Fill(vParticle1->E(), (vParticle1->E() + esdtrack->E() - clsE)/vParticle1->E());
+              continue;
+            }
+            if(vParticle2->PdgCode()==22 && vParticle1 == vParticle)
+            {
+              fHPhotonEdiff0HC->Fill(vParticle2->E(), (vParticle2->E() - clsE)/vParticle2->E());
+              fHPhotonEVsClsE->Fill(vParticle2->E(), clsE);
+              
+              if((clsE - 0.3*esdtrack->E())<0) fHPhotonEdiff30HC->Fill(vParticle2->E(), 1);
+              else fHPhotonEdiff30HC->Fill(vParticle2->E(), (vParticle2->E() + 0.3*esdtrack->E() - clsE)/vParticle2->E());
+              
+              if((clsE - 0.7*esdtrack->E())<0) fHPhotonEdiff70HC->Fill(vParticle2->E(), 1);
+              else fHPhotonEdiff70HC->Fill(vParticle2->E(), (vParticle2->E() + 0.7*esdtrack->E() - clsE)/vParticle2->E());
+              
+              if((clsE-esdtrack->E())<0) fHPhotonEdiff100HC->Fill(vParticle2->E(), 1);
+              else fHPhotonEdiff100HC->Fill(vParticle2->E(), (vParticle2->E() + esdtrack->E() - clsE)/vParticle2->E());
+            }
+          }
+        }
+       if(newTrack) delete newTrack;
+       break;
       }
+      if(newTrack) delete newTrack;
     }
-  } 
+  }
 }
 
+
 //________________________________________________________________________
 void AliAnalysisTaskSOH::ProcessScaleFactor()
 {
@@ -821,21 +908,20 @@ void AliAnalysisTaskSOH::ProcessScaleFactor()
 AliESDtrack *AliAnalysisTaskSOH::GetAcceptTrack(AliESDtrack *esdtrack)
 {
   // Get accepted track.
-
-  static AliESDtrack newTrack;
+  AliESDtrack *newTrack = 0x0;
   if(fEsdTrackCuts->AcceptTrack(esdtrack))
   {
-    esdtrack->Copy(newTrack);
-    newTrack.SetTRDQuality(0);
+    newTrack = new AliESDtrack(*esdtrack);
+    newTrack->SetTRDQuality(0);
   }
   else if(fHybridTrackCuts1->AcceptTrack(esdtrack))
   {
     if(esdtrack->GetConstrainedParam())
     {
-      esdtrack->Copy(newTrack);
+      newTrack = new AliESDtrack(*esdtrack);
       const AliExternalTrackParam* constrainParam = esdtrack->GetConstrainedParam();
-      newTrack.Set(constrainParam->GetX(),constrainParam->GetAlpha(),constrainParam->GetParameter(),constrainParam->GetCovariance());
-      newTrack.SetTRDQuality(1);               
+      newTrack->Set(constrainParam->GetX(),constrainParam->GetAlpha(),constrainParam->GetParameter(),constrainParam->GetCovariance());
+      newTrack->SetTRDQuality(1);              
     }
     else 
       return 0x0;
@@ -844,10 +930,10 @@ AliESDtrack *AliAnalysisTaskSOH::GetAcceptTrack(AliESDtrack *esdtrack)
   {
     if(esdtrack->GetConstrainedParam())
     {
-      esdtrack->Copy(newTrack);
+      newTrack = new AliESDtrack(*esdtrack);
       const AliExternalTrackParam* constrainParam = esdtrack->GetConstrainedParam();
-      newTrack.Set(constrainParam->GetX(),constrainParam->GetAlpha(),constrainParam->GetParameter(),constrainParam->GetCovariance());
-      newTrack.SetTRDQuality(2);               
+      newTrack->Set(constrainParam->GetX(),constrainParam->GetAlpha(),constrainParam->GetParameter(),constrainParam->GetCovariance());
+      newTrack->SetTRDQuality(2);              
     }
     else 
       return 0x0;
@@ -857,7 +943,7 @@ AliESDtrack *AliAnalysisTaskSOH::GetAcceptTrack(AliESDtrack *esdtrack)
     return 0x0;
   }
 
-  return &newTrack;
+  return newTrack;
 }
 
 //________________________________________________________________________
@@ -887,3 +973,16 @@ void AliAnalysisTaskSOH::Terminate(Option_t *)
 {
   // Terminate analysis.
 }
+
+//________________________________________________________________________
+Bool_t AliAnalysisTaskSOH::EsdVertexOk() const
+{
+  // Modified from AliAnalyseLeadingTrackUE::VertexSelection()  
+
+  const AliESDVertex* vtx = fESD->GetPrimaryVertex();
+  if (!vtx) return kFALSE;
+  Int_t nContributors = vtx->GetNContributors();
+  Double_t zVertex    = vtx->GetZ();
+  if( nContributors < 1 || TMath::Abs(zVertex) > fZVtxMax ) return kFALSE;
+  return kTRUE;
+}