PWG4/PartCorrDep/AliAnaElectron.cxx

   1 /**************************************************************************\r
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
   3  *                                                                        *\r
   4  * Author: The ALICE Off-line Project.                                    *\r
   5  * Contributors are mentioned in the code where appropriate.              *\r
   6  *                                                                        *\r
   7  * Permission to use, copy, modify and distribute this software and its   *\r
   8  * documentation strictly for non-commercial purposes hereby granted      *\r
   9  * without fee, provided that the above copyright notice appears in all   *\r
  10  * copies and that both the copyright notice and this permission notice   *\r
  11  * appear in the supporting documentation. The authors make no claims     *\r
  12  * about the suitability of this software for any purpose. It is          *\r
  13  * provided "as is" without express or implied warranty.                  *\r
  14  **************************************************************************/\r
  15 /* $Id: $ */\r
  16 \r
  17 //_________________________________________________________________________\r
  18 //\r
  19 // Class for the electron identification.\r
  20 // Clusters from EMCAL matched to tracks\r
  21 // and kept in the AOD. Few histograms produced.\r
  22 //\r
  23 // -- Author: J.L. Klay (Cal Poly), M. Heinz (Yale)\r
  24 //////////////////////////////////////////////////////////////////////////////\r
  25   \r
  26 // --- ROOT system --- \r
  27 #include <TH2F.h>\r
  28 #include <TParticle.h>\r
  29 #include <TNtuple.h>\r
  30 #include <TClonesArray.h>\r
  31 #include <TObjString.h>\r
  32 //#include <Riostream.h>\r
  33 \r
  34 // --- Analysis system --- \r
  35 #include "AliAnaElectron.h" \r
  36 #include "AliCaloTrackReader.h"\r
  37 #include "AliMCAnalysisUtils.h"\r
  38 #include "AliAODCaloCluster.h"\r
  39 #include "AliFidutialCut.h"\r
  40 #include "AliAODTrack.h"\r
  41 #include "AliAODPid.h"\r
  42 #include "AliCaloPID.h"\r
  43 #include "AliAODMCParticle.h"\r
  44 #include "AliStack.h"\r
  45 #include "AliExternalTrackParam.h"\r
  46 #include "AliESDv0.h"\r
  47 #include "AliESDtrack.h"\r
  48 #include "AliAODJet.h"\r
  49 #include "AliAODEvent.h"\r
  50 \r
  51 ClassImp(AliAnaElectron)\r
  52   \r
  53 //____________________________________________________________________________\r
  54 AliAnaElectron::AliAnaElectron() \r
  55 : AliAnaPartCorrBaseClass(),fCalorimeter(""),\r
  56   fpOverEmin(0.),fpOverEmax(0.),fResidualCut(0.),\r
  57   fDrCut(0.),fPairDcaCut(0.),fDecayLenCut(0.),fImpactCut(0.),\r
  58   fAssocPtCut(0.),fMassCut(0.),fSdcaCut(0.),fITSCut(0),\r
  59   fNTagTrkCut(0),fIPSigCut(0.),\r
  60   fWriteNtuple(kFALSE),\r
  61   //matching checks\r
  62   fEleNtuple(0),\r
  63   fh1pOverE(0),fh1dR(0),fh2EledEdx(0),fh2MatchdEdx(0),fh2dEtadPhi(0),\r
  64   fh2dEtadPhiMatched(0),fh2dEtadPhiUnmatched(0),\r
  65   fh2TrackPVsClusterE(0),fh2TrackPtVsClusterE(0),fh2TrackPhiVsClusterPhi(0),fh2TrackEtaVsClusterEta(0),\r
  66   //Photonic electron checks\r
  67   fh1OpeningAngle(0),fh1MinvPhoton(0),\r
  68   //reco\r
  69   fhPtElectron(0),fhPhiElectron(0),fhEtaElectron(0),\r
  70   fhPtNPE(0),fhPhiNPE(0),fhEtaNPE(0),\r
  71   fhPtPE(0),fhPhiPE(0),fhEtaPE(0),\r
  72   fhPtConversion(0),fhPhiConversion(0),fhEtaConversion(0),\r
  73   fhPtBottom(0),fhPhiBottom(0),fhEtaBottom(0),\r
  74   fhPtCharm(0),fhPhiCharm(0),fhEtaCharm(0),\r
  75   fhPtCFromB(0),fhPhiCFromB(0),fhEtaCFromB(0),\r
  76   fhPtDalitz(0),fhPhiDalitz(0),fhEtaDalitz(0),\r
  77   fhPtWDecay(0),fhPhiWDecay(0),fhEtaWDecay(0),\r
  78   fhPtZDecay(0),fhPhiZDecay(0),fhEtaZDecay(0),\r
  79   fhPtAll(0),fhPhiAll(0),fhEtaAll(0),\r
  80   fhPtUnknown(0),fhPhiUnknown(0),fhEtaUnknown(0),\r
  81   fhPtMisidentified(0),fhPhiMisidentified(0),fhEtaMisidentified(0),\r
  82   fhPtHadron(0),fhPtEleTrkDet(0),\r
  83   //event QA\r
  84   fhImpactXY(0),fhRefMult(0),fhRefMult2(0),\r
  85   //B-tagging\r
  86   fhDVMBtagCut1(0),fhDVMBtagCut2(0),fhDVMBtagCut3(0),fhDVMBtagQA1(0),fhDVMBtagQA2(0),fhDVMBtagQA3(0),\r
  87   fhDVMBtagQA4(0),fhDVMBtagQA5(0),fhIPSigBtagQA1(0),fhIPSigBtagQA2(0),\r
  88   //B-jets\r
  89   fhJetType(0),fhBJetXsiFF(0),fhBJetPtFF(0),fhBJetEtaPhi(0),fhNonBJetXsiFF(0),fhNonBJetPtFF(0),fhNonBJetEtaPhi(0),\r
  90   //MC\r
  91   fMCEleNtuple(0),fhPtMCHadron(0),fhPtMCBottom(0),fhPtMCCharm(0),fhPtMCCFromB(0),fhPtMCConversion(0),\r
  92   fhPtMCDalitz(0),fhPtMCWDecay(0),fhPtMCZDecay(0),fhPtMCUnknown(0)\r
  93 {\r
  94   //default ctor\r
  95   \r
  96   //Initialize parameters\r
  97   InitParameters();\r
  98 \r
  99 }\r
 100 \r
 101 //____________________________________________________________________________\r
 102 AliAnaElectron::AliAnaElectron(const AliAnaElectron & g) \r
 103   : AliAnaPartCorrBaseClass(g),fCalorimeter(g.fCalorimeter),\r
 104     fpOverEmin(g.fpOverEmin),fpOverEmax(g.fpOverEmax),fResidualCut(g.fResidualCut),\r
 105     fDrCut(g.fDrCut),fPairDcaCut(g.fPairDcaCut),fDecayLenCut(g.fDecayLenCut),fImpactCut(g.fImpactCut),\r
 106     fAssocPtCut(g.fAssocPtCut),fMassCut(g.fMassCut),fSdcaCut(g.fSdcaCut),fITSCut(g.fITSCut),\r
 107     fNTagTrkCut(g.fNTagTrkCut),fIPSigCut(g.fIPSigCut),\r
 108     fWriteNtuple(g.fWriteNtuple),\r
 109     //matching checks\r
 110     fEleNtuple(g.fEleNtuple),\r
 111     fh1pOverE(g.fh1pOverE),fh1dR(g.fh1dR),\r
 112     fh2EledEdx(g.fh2EledEdx),fh2MatchdEdx(g.fh2MatchdEdx),fh2dEtadPhi(g.fh2dEtadPhi),\r
 113     fh2dEtadPhiMatched(g.fh2dEtadPhiMatched),fh2dEtadPhiUnmatched(g.fh2dEtadPhiUnmatched),\r
 114     fh2TrackPVsClusterE(g.fh2TrackPVsClusterE),fh2TrackPtVsClusterE(g.fh2TrackPtVsClusterE),\r
 115     fh2TrackPhiVsClusterPhi(g.fh2TrackPhiVsClusterPhi),fh2TrackEtaVsClusterEta(g.fh2TrackEtaVsClusterEta),\r
 116     //Photonic electron checks\r
 117     fh1OpeningAngle(g.fh1OpeningAngle),fh1MinvPhoton(g.fh1MinvPhoton),\r
 118     //reco\r
 119     fhPtElectron(g.fhPtElectron),fhPhiElectron(g.fhPhiElectron),fhEtaElectron(g.fhEtaElectron),\r
 120     fhPtNPE(g.fhPtNPE),fhPhiNPE(g.fhPhiNPE),fhEtaNPE(g.fhEtaNPE),\r
 121     fhPtPE(g.fhPtPE),fhPhiPE(g.fhPhiPE),fhEtaPE(g.fhEtaPE),\r
 122     fhPtConversion(g.fhPtConversion),fhPhiConversion(g.fhPhiConversion),fhEtaConversion(g.fhEtaConversion),\r
 123     fhPtBottom(g.fhPtBottom),fhPhiBottom(g.fhPhiBottom),fhEtaBottom(g.fhEtaBottom),\r
 124     fhPtCharm(g.fhPtCharm),fhPhiCharm(g.fhPhiCharm),fhEtaCharm(g.fhEtaCharm),\r
 125     fhPtCFromB(g.fhPtCFromB),fhPhiCFromB(g.fhPhiCFromB),fhEtaCFromB(g.fhEtaCFromB),\r
 126     fhPtDalitz(g.fhPtDalitz),fhPhiDalitz(g.fhPhiDalitz),fhEtaDalitz(g.fhEtaDalitz),\r
 127     fhPtWDecay(g.fhPtWDecay),fhPhiWDecay(g.fhPhiWDecay),fhEtaWDecay(g.fhEtaWDecay),\r
 128     fhPtZDecay(g.fhPtZDecay),fhPhiZDecay(g.fhPhiZDecay),fhEtaZDecay(g.fhEtaZDecay),\r
 129     fhPtAll(g.fhPtAll),fhPhiAll(g.fhPhiAll),fhEtaAll(g.fhEtaAll),\r
 130     fhPtUnknown(g.fhPtUnknown),fhPhiUnknown(g.fhPhiUnknown),fhEtaUnknown(g.fhEtaUnknown),\r
 131     fhPtMisidentified(g.fhPtMisidentified),fhPhiMisidentified(g.fhPhiMisidentified),fhEtaMisidentified(g.fhEtaMisidentified),\r
 132     fhPtHadron(g.fhPtHadron),fhPtEleTrkDet(g.fhPtEleTrkDet),\r
 133     //event QA\r
 134     fhImpactXY(g.fhImpactXY),fhRefMult(g.fhRefMult),fhRefMult2(g.fhRefMult2),\r
 135     //B-tagging\r
 136     fhDVMBtagCut1(g.fhDVMBtagCut1),fhDVMBtagCut2(g.fhDVMBtagCut2),fhDVMBtagCut3(g.fhDVMBtagCut3),\r
 137     fhDVMBtagQA1(g.fhDVMBtagQA1),fhDVMBtagQA2(g.fhDVMBtagQA2),fhDVMBtagQA3(g.fhDVMBtagQA3),\r
 138     fhDVMBtagQA4(g.fhDVMBtagQA4),fhDVMBtagQA5(g.fhDVMBtagQA5),fhIPSigBtagQA1(g.fhIPSigBtagQA1),\r
 139     fhIPSigBtagQA2(g.fhIPSigBtagQA2),\r
 140     //B-jets\r
 141     fhJetType(g.fhJetType),fhBJetXsiFF(g.fhBJetXsiFF),fhBJetPtFF(g.fhBJetPtFF),fhBJetEtaPhi(g.fhBJetEtaPhi),\r
 142     fhNonBJetXsiFF(g.fhNonBJetXsiFF),fhNonBJetPtFF(g.fhNonBJetPtFF),fhNonBJetEtaPhi(g.fhNonBJetEtaPhi),\r
 143     //MC\r
 144     fMCEleNtuple(g.fMCEleNtuple),fhPtMCHadron(g.fhPtMCHadron),fhPtMCBottom(g.fhPtMCBottom),\r
 145     fhPtMCCharm(g.fhPtMCCharm),fhPtMCCFromB(g.fhPtMCCFromB),fhPtMCConversion(g.fhPtMCConversion),\r
 146     fhPtMCDalitz(g.fhPtMCDalitz),fhPtMCWDecay(g.fhPtMCWDecay),\r
 147     fhPtMCZDecay(g.fhPtMCZDecay),fhPtMCUnknown(g.fhPtMCUnknown)\r
 148 {\r
 149   // cpy ctor\r
 150   \r
 151 }\r
 152 \r
 153 //_________________________________________________________________________\r
 154 AliAnaElectron & AliAnaElectron::operator = (const AliAnaElectron & g)\r
 155 {\r
 156   // assignment operator\r
 157   \r
 158   if(&g == this) return *this;\r
 159   fCalorimeter = g.fCalorimeter;\r
 160   fpOverEmin = g.fpOverEmin;\r
 161   fpOverEmax = g.fpOverEmax;\r
 162   fResidualCut = g.fResidualCut;\r
 163   fDrCut = g.fDrCut;\r
 164   fPairDcaCut = g.fPairDcaCut;\r
 165   fDecayLenCut = g.fDecayLenCut;\r
 166   fImpactCut = g.fImpactCut;\r
 167   fAssocPtCut = g.fAssocPtCut;\r
 168   fMassCut = g.fMassCut;\r
 169   fSdcaCut = g.fSdcaCut;\r
 170   fITSCut = g.fITSCut;\r
 171   fNTagTrkCut = g.fNTagTrkCut;\r
 172   fIPSigCut = g.fIPSigCut;\r
 173   fWriteNtuple = g.fWriteNtuple;\r
 174   fEleNtuple = g.fEleNtuple;\r
 175   fh1pOverE = g.fh1pOverE;\r
 176   fh1dR = g.fh1dR;\r
 177   fh2EledEdx = g.fh2EledEdx;\r
 178   fh2MatchdEdx = g.fh2MatchdEdx;\r
 179   fh2dEtadPhi = g.fh2dEtadPhi;\r
 180   fh2dEtadPhiMatched = g.fh2dEtadPhiMatched;\r
 181   fh2dEtadPhiUnmatched = g.fh2dEtadPhiUnmatched;\r
 182   fh2TrackPVsClusterE = g.fh2TrackPVsClusterE;\r
 183   fh2TrackPtVsClusterE = g.fh2TrackPtVsClusterE;\r
 184   fh2TrackPhiVsClusterPhi = g.fh2TrackPhiVsClusterPhi;\r
 185   fh2TrackEtaVsClusterEta = g.fh2TrackEtaVsClusterEta;   \r
 186   fh1OpeningAngle = g.fh1OpeningAngle;\r
 187   fh1MinvPhoton = g.fh1MinvPhoton;\r
 188   fhPtElectron = g.fhPtElectron;\r
 189   fhPhiElectron = g.fhPhiElectron;\r
 190   fhEtaElectron = g.fhEtaElectron;\r
 191   fhPtNPE = g.fhPtNPE;\r
 192   fhPhiNPE = g.fhPhiNPE;\r
 193   fhEtaNPE = g.fhEtaNPE;\r
 194   fhPtPE = g.fhPtPE;\r
 195   fhPhiPE = g.fhPhiPE;\r
 196   fhEtaPE = g.fhEtaPE;\r
 197   fhPtConversion = g.fhPtConversion;\r
 198   fhPhiConversion = g.fhPhiConversion;\r
 199   fhEtaConversion = g.fhEtaConversion;\r
 200   fhPtBottom = g.fhPtBottom;\r
 201   fhPhiBottom = g.fhPhiBottom;\r
 202   fhEtaBottom = g.fhEtaBottom;\r
 203   fhPtCharm = g.fhPtCharm;\r
 204   fhPhiCharm = g.fhPhiCharm;\r
 205   fhEtaCharm = g.fhEtaCharm;\r
 206   fhPtCFromB = g.fhPtCFromB;\r
 207   fhPhiCFromB = g.fhPhiCFromB;\r
 208   fhEtaCFromB = g.fhEtaCFromB;\r
 209   fhPtDalitz = g.fhPtDalitz;\r
 210   fhPhiDalitz = g.fhPhiDalitz;\r
 211   fhEtaDalitz = g.fhEtaDalitz;\r
 212   fhPtWDecay = g.fhPtWDecay;\r
 213   fhPhiWDecay = g.fhPhiWDecay;\r
 214   fhEtaWDecay = g.fhEtaWDecay;\r
 215   fhPtZDecay = g.fhPtZDecay;\r
 216   fhPhiZDecay = g.fhPhiZDecay;\r
 217   fhEtaZDecay = g.fhEtaZDecay;\r
 218   fhPtAll = g.fhPtAll;\r
 219   fhPhiAll = g.fhPhiAll;\r
 220   fhEtaAll = g.fhEtaAll;\r
 221   fhPtUnknown = g.fhPtUnknown;\r
 222   fhPhiUnknown = g.fhPhiUnknown;\r
 223   fhEtaUnknown = g.fhEtaUnknown;\r
 224   fhPtMisidentified = g.fhPtMisidentified;\r
 225   fhPhiMisidentified = g.fhPhiMisidentified;\r
 226   fhEtaMisidentified = g.fhEtaMisidentified;\r
 227 \r
 228   fhPtHadron = g.fhPtHadron;\r
 229   fhPtEleTrkDet = g.fhPtEleTrkDet;\r
 230 \r
 231   //event QA\r
 232   fhImpactXY = g.fhImpactXY;\r
 233   fhRefMult = g.fhRefMult;\r
 234   fhRefMult2 = g.fhRefMult2;\r
 235 \r
 236   //B-tagging\r
 237   fhDVMBtagCut1 = g.fhDVMBtagCut1;\r
 238   fhDVMBtagCut2 = g.fhDVMBtagCut2;\r
 239   fhDVMBtagCut3 = g.fhDVMBtagCut3;\r
 240   fhDVMBtagQA1 = g.fhDVMBtagQA1;\r
 241   fhDVMBtagQA2 = g.fhDVMBtagQA2;\r
 242   fhDVMBtagQA3 = g.fhDVMBtagQA3;\r
 243   fhDVMBtagQA4 = g.fhDVMBtagQA4;\r
 244   fhDVMBtagQA5 = g.fhDVMBtagQA5;\r
 245   fhIPSigBtagQA1 = g.fhIPSigBtagQA1;\r
 246   fhIPSigBtagQA2 = g.fhIPSigBtagQA2;\r
 247 \r
 248   fhJetType = g.fhJetType;\r
 249   fhBJetXsiFF = g.fhBJetXsiFF;\r
 250   fhBJetPtFF = g.fhBJetPtFF;\r
 251   fhBJetEtaPhi = g.fhBJetEtaPhi;\r
 252   fhNonBJetXsiFF = g.fhNonBJetXsiFF;\r
 253   fhNonBJetPtFF = g.fhNonBJetPtFF;\r
 254   fhNonBJetEtaPhi = g.fhNonBJetEtaPhi;\r
 255 \r
 256   fMCEleNtuple = g.fMCEleNtuple;\r
 257   fhPtMCHadron = g.fhPtMCHadron;\r
 258   fhPtMCBottom = g.fhPtMCBottom;\r
 259   fhPtMCCharm = g.fhPtMCCharm;\r
 260   fhPtMCCFromB = g.fhPtMCCFromB;\r
 261   fhPtMCConversion = g.fhPtMCConversion;\r
 262   fhPtMCDalitz = g.fhPtMCDalitz;\r
 263   fhPtMCWDecay = g.fhPtMCWDecay;\r
 264   fhPtMCZDecay = g.fhPtMCZDecay;\r
 265   fhPtMCUnknown = g.fhPtMCUnknown;\r
 266 \r
 267   return *this;\r
 268   \r
 269 }\r
 270 \r
 271 //____________________________________________________________________________\r
 272 AliAnaElectron::~AliAnaElectron() \r
 273 {\r
 274   //dtor\r
 275 \r
 276 }\r
 277 \r
 278 \r
 279 //________________________________________________________________________\r
 280 TList *  AliAnaElectron::GetCreateOutputObjects()\r
 281 {  \r
 282   // Create histograms to be saved in output file and \r
 283   // store them in outputContainer\r
 284   TList * outputContainer = new TList() ; \r
 285   outputContainer->SetName("ElectronHistos") ; \r
 286   \r
 287   //created ele ntuple for further analysis\r
 288   if(fWriteNtuple) {\r
 289       fEleNtuple = new TNtuple("EleNtuple","Electron Ntuple","tmctag:cmctag:pt:phi:eta:p:E:deta:dphi:nCells:dEdx:pidProb:impXY:impZ");\r
 290     outputContainer->Add(fEleNtuple) ;\r
 291   }\r
 292 \r
 293   Int_t nptbins  = GetHistoNPtBins();\r
 294   Int_t nphibins = GetHistoNPhiBins();\r
 295   Int_t netabins = GetHistoNEtaBins();\r
 296   Float_t ptmax  = GetHistoPtMax();\r
 297   Float_t phimax = GetHistoPhiMax();\r
 298   Float_t etamax = GetHistoEtaMax();\r
 299   Float_t ptmin  = GetHistoPtMin();\r
 300   Float_t phimin = GetHistoPhiMin();\r
 301   Float_t etamin = GetHistoEtaMin();    \r
 302 \r
 303   fh1pOverE = new TH1F("h1pOverE","EMCAL-TRACK matches p/E",100,0.,10.);\r
 304   fh1dR = new TH1F("h1dR","EMCAL-TRACK matches dR",300, 0.,TMath::Pi());\r
 305   fh2EledEdx = new TH2F("h2EledEdx","dE/dx vs. p for electrons",200,0.,50.,200,0.,400.);\r
 306   fh2MatchdEdx = new TH2F("h2MatchdEdx","dE/dx vs. p for all matches",200,0.,50.,200,0.,400.);\r
 307   fh2dEtadPhi = new TH2F("h2dEtadPhi","#Delta#eta vs. #Delta#phi for all track-cluster pairs",200,0.,1.4,300,0.,TMath::Pi());\r
 308   fh2dEtadPhiMatched = new TH2F("h2dEtadPhiMatched","#Delta#eta vs. #Delta#phi for matched track-cluster pairs",200,0.,1.4,300,0.,TMath::Pi());\r
 309   fh2dEtadPhiUnmatched = new TH2F("h2dEtadPhiUnmatched","#Delta#eta vs. #Delta#phi for unmatched track-cluster pairs",200,0.,1.4,300,0.,TMath::Pi());\r
 310 \r
 311   fh2TrackPVsClusterE = new TH2F("h2TrackPVsClusterE","h2TrackPVsClusterE",nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);\r
 312   fh2TrackPtVsClusterE = new TH2F("h2TrackPtVsClusterE","h2TrackPtVsClusterE",nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);\r
 313   fh2TrackPhiVsClusterPhi = new TH2F("h2TrackPhiVsClusterPhi","h2TrackPhiVsClusterPhi",nphibins,phimin,phimax,nphibins,phimin,phimax);\r
 314   fh2TrackEtaVsClusterEta = new TH2F("h2TrackEtaVsClusterEta","h2TrackEtaVsClusterEta",netabins,etamin,etamax,netabins,etamin,etamax);\r
 315 \r
 316   outputContainer->Add(fh1pOverE) ; \r
 317   outputContainer->Add(fh1dR) ; \r
 318   outputContainer->Add(fh2EledEdx) ;\r
 319   outputContainer->Add(fh2MatchdEdx) ;\r
 320   outputContainer->Add(fh2dEtadPhi) ;\r
 321   outputContainer->Add(fh2dEtadPhiMatched) ;\r
 322   outputContainer->Add(fh2dEtadPhiUnmatched) ;\r
 323   outputContainer->Add(fh2TrackPVsClusterE) ;\r
 324   outputContainer->Add(fh2TrackPtVsClusterE) ;\r
 325   outputContainer->Add(fh2TrackPhiVsClusterPhi) ;\r
 326   outputContainer->Add(fh2TrackEtaVsClusterEta) ;\r
 327   \r
 328   //photonic electron checks\r
 329   fh1OpeningAngle = new TH1F("hOpeningAngle","Opening angle between e+e- pairs",100,0.,TMath::Pi());\r
 330   fh1MinvPhoton = new TH1F("hMinvPhoton","Invariant mass of e+e- pairs",200,0.,2.);\r
 331 \r
 332   outputContainer->Add(fh1OpeningAngle);\r
 333   outputContainer->Add(fh1MinvPhoton);\r
 334 \r
 335   fhPtElectron = new TH1F("hPtElectron","Electron pT",nptbins,ptmin,ptmax);\r
 336   fhPhiElectron = new TH2F("hPhiElectron","Electron phi vs pT",nptbins,ptmin,ptmax,nphibins,phimin,phimax);\r
 337   fhEtaElectron = new TH2F("hEtaElectron","Electron eta vs. eta",nptbins,ptmin,ptmax,netabins,etamin,etamax);\r
 338   fhPtNPE = new TH1F("hPtNPE","Non-photonic Electron pT",nptbins,ptmin,ptmax);\r
 339   fhPhiNPE = new TH2F("hPhiNPE","Non-photonic Electron phi vs pT",nptbins,ptmin,ptmax,nphibins,phimin,phimax);\r
 340   fhEtaNPE = new TH2F("hEtaNPE","Non-photonic Electron eta vs. eta",nptbins,ptmin,ptmax,netabins,etamin,etamax);\r
 341   fhPtPE = new TH1F("hPtPE","Photonic Electron pT",nptbins,ptmin,ptmax);\r
 342   fhPhiPE = new TH2F("hPhiPE","Photonic Electron phi vs pT",nptbins,ptmin,ptmax,nphibins,phimin,phimax);\r
 343   fhEtaPE = new TH2F("hEtaPE","Photonic Electron eta vs. eta",nptbins,ptmin,ptmax,netabins,etamin,etamax);\r
 344 \r
 345   outputContainer->Add(fhPtElectron) ; \r
 346   outputContainer->Add(fhPhiElectron) ; \r
 347   outputContainer->Add(fhEtaElectron) ;\r
 348   outputContainer->Add(fhPtNPE) ; \r
 349   outputContainer->Add(fhPhiNPE) ; \r
 350   outputContainer->Add(fhEtaNPE) ;\r
 351   outputContainer->Add(fhPtPE) ; \r
 352   outputContainer->Add(fhPhiPE) ; \r
 353   outputContainer->Add(fhEtaPE) ;\r
 354 \r
 355   fhPtHadron = new TH1F("hPtHadron","Charged hadrons w/in EMCAL acceptance",nptbins,ptmin,ptmax);\r
 356   fhPtEleTrkDet = new TH1F("hPtEleTrkDet","Electrons identified by tracking detectors w/in EMCAL acceptance",nptbins,ptmin,ptmax);\r
 357 \r
 358   outputContainer->Add(fhPtHadron);\r
 359   outputContainer->Add(fhPtEleTrkDet);\r
 360 \r
 361   //event QA\r
 362   fhImpactXY = new TH1F("hImpactXY","Impact parameter for all tracks",200,-10,10.);\r
 363   fhRefMult = new TH1F("hRefMult" ,"refmult QA: " ,100,0,200);\r
 364   fhRefMult2  = new TH1F("hRefMult2" ,"refmult2 QA: " ,100,0,200);\r
 365 \r
 366   outputContainer->Add(fhImpactXY);\r
 367   outputContainer->Add(fhRefMult);\r
 368   outputContainer->Add(fhRefMult2);\r
 369 \r
 370   //B-tagging\r
 371   fhDVMBtagCut1 = new TH2F("hdvmbtag_cut1","DVM B-tag result cut1", 10,0,10 ,nptbins,ptmin,ptmax);\r
 372   fhDVMBtagCut2 = new TH2F("hdvmbtag_cut2","DVM B-tag result cut2", 10,0,10 ,nptbins,ptmin,ptmax);\r
 373   fhDVMBtagCut3 = new TH2F("hdvmbtag_cut3","DVM B-tag result cut3", 10,0,10 ,nptbins,ptmin,ptmax);\r
 374   fhDVMBtagQA1  = new TH2F("hdvmbtag_qa1" ,"DVM B-tag QA: pairDCA vs length", 100,0,0.2 ,100,0,1.0);\r
 375   fhDVMBtagQA2  = new TH2F("hdvmbtag_qa2" ,"DVM B-tag QA: signDCA vs mass"  , 200,-0.5,0.5 ,100,0,10);\r
 376   fhDVMBtagQA3  = new TH1F("hdvmbtag_qa3" ,"DVM B-tag QA: ITS-Hits electron" ,7,0,7);\r
 377   fhDVMBtagQA4  = new TH1F("hdvmbtag_qa4" ,"DVM B-tag QA: IP d electron" ,200,-3,3);\r
 378   fhDVMBtagQA5  = new TH1F("hdvmbtag_qa5" ,"DVM B-tag QA: IP z electron" ,200,-3,3);\r
 379   fhIPSigBtagQA1  = new TH1F("hipsigbtag_qa1" ,"IPSig B-tag QA: # tag tracks", 20,0,20);\r
 380   fhIPSigBtagQA2  = new TH1F("hipsigbtag_qa2" ,"IPSig B-tag QA: IP significance", 200,-10.,10.);\r
 381 \r
 382   outputContainer->Add(fhDVMBtagCut1) ;\r
 383   outputContainer->Add(fhDVMBtagCut2) ;\r
 384   outputContainer->Add(fhDVMBtagCut3) ;\r
 385   outputContainer->Add(fhDVMBtagQA1) ;\r
 386   outputContainer->Add(fhDVMBtagQA2) ;\r
 387   outputContainer->Add(fhDVMBtagQA3) ;\r
 388   outputContainer->Add(fhDVMBtagQA4) ;\r
 389   outputContainer->Add(fhDVMBtagQA5) ;\r
 390   outputContainer->Add(fhIPSigBtagQA1) ;\r
 391   outputContainer->Add(fhIPSigBtagQA2) ;\r
 392 \r
 393   fhJetType = new TH2F("hJetType","# jets passing each tag method vs jet pt",10,0,10,300,0.,300.);\r
 394   fhBJetXsiFF = new TH2F("hBJetXsiFF","B-jet #Xsi Frag. Fn.",100,0.,10.,300,0.,300.);\r
 395   fhBJetPtFF = new TH2F("hBJetPtFF","B-jet p_{T} Frag. Fn.",nptbins,ptmin,ptmax,300,0.,300.);\r
 396   fhBJetEtaPhi = new TH2F("hBJetEtaPhi","B-jet eta-phi distribution",netabins,etamin,etamax,nphibins,phimin,phimax);\r
 397   fhNonBJetXsiFF = new TH2F("hNonBJetXsiFF","Non B-jet #Xsi Frag. Fn.",100,0.,10.,300,0.,300.);\r
 398   fhNonBJetPtFF = new TH2F("hNonBJetPtFF","Non B-jet p_{T} Frag. Fn.",nptbins,ptmin,ptmax,300,0.,300.);\r
 399   fhNonBJetEtaPhi = new TH2F("hNonBJetEtaPhi","Non B-jet eta-phi distribution",netabins,etamin,etamax,nphibins,phimin,phimax);\r
 400 \r
 401   outputContainer->Add(fhJetType);\r
 402   outputContainer->Add(fhBJetXsiFF);\r
 403   outputContainer->Add(fhBJetPtFF);\r
 404   outputContainer->Add(fhBJetEtaPhi);\r
 405   outputContainer->Add(fhNonBJetXsiFF);\r
 406   outputContainer->Add(fhNonBJetPtFF);\r
 407   outputContainer->Add(fhNonBJetEtaPhi);\r
 408 \r
 409   if(IsDataMC()){\r
 410     \r
 411     fhPtConversion = new TH1F("hPtConversion","Conversion electron pT",nptbins,ptmin,ptmax);\r
 412     fhPhiConversion = new TH2F("hPhiConversion","Conversion Electron phi vs pT",nptbins,ptmin,ptmax,nphibins,phimin,phimax);\r
 413     fhEtaConversion = new TH2F("hEtaConversion","Conversion Electron eta vs. eta",nptbins,ptmin,ptmax,netabins,etamin,etamax);\r
 414     fhPtBottom = new TH1F("hPtBottom","Bottom electron pT",nptbins,ptmin,ptmax);\r
 415     fhPhiBottom = new TH2F("hPhiBottom","Bottom Electron phi vs pT",nptbins,ptmin,ptmax,nphibins,phimin,phimax);\r
 416     fhEtaBottom = new TH2F("hEtaBottom","Bottom Electron eta vs. eta",nptbins,ptmin,ptmax,netabins,etamin,etamax);\r
 417     fhPtCharm = new TH1F("hPtCharm","Charm electron pT",nptbins,ptmin,ptmax);\r
 418     fhPhiCharm = new TH2F("hPhiCharm","Charm Electron phi vs pT",nptbins,ptmin,ptmax,nphibins,phimin,phimax);\r
 419     fhEtaCharm = new TH2F("hEtaCharm","Charm Electron eta vs. eta",nptbins,ptmin,ptmax,netabins,etamin,etamax);\r
 420     fhPtCFromB = new TH1F("hPtCFromB","Charm from Bottom electron pT",nptbins,ptmin,ptmax);\r
 421     fhPhiCFromB = new TH2F("hPhiCFromB","Charm from Bottom Electron phi vs pT",nptbins,ptmin,ptmax,nphibins,phimin,phimax);\r
 422     fhEtaCFromB = new TH2F("hEtaCFromB","Charm from Bottom Electron eta vs. eta",nptbins,ptmin,ptmax,netabins,etamin,etamax);\r
 423     fhPtDalitz = new TH1F("hPtDalitz","Dalitz electron pT",nptbins,ptmin,ptmax);\r
 424     fhPhiDalitz = new TH2F("hPhiDalitz","Dalitz Electron phi vs pT",nptbins,ptmin,ptmax,nphibins,phimin,phimax);\r
 425     fhEtaDalitz = new TH2F("hEtaDalitz","Dalitz Electron eta vs. eta",nptbins,ptmin,ptmax,netabins,etamin,etamax);\r
 426     fhPtWDecay = new TH1F("hPtWDecay","W-boson Electron pT",nptbins,ptmin,ptmax);\r
 427     fhPhiWDecay = new TH2F("hPhiWDecay","W-boson electron phi vs pT",nptbins,ptmin,ptmax,nphibins,phimin,phimax);\r
 428     fhEtaWDecay = new TH2F("hEtaWDecay","W-boson Electron eta vs. eta",nptbins,ptmin,ptmax,netabins,etamin,etamax);\r
 429     fhPtZDecay = new TH1F("hPtZDecay","Z-boson electron pT",nptbins,ptmin,ptmax);\r
 430     fhPhiZDecay = new TH2F("hPhiZDecay","Z-boson Electron phi vs pT",nptbins,ptmin,ptmax,nphibins,phimin,phimax);\r
 431     fhEtaZDecay = new TH2F("hEtaZDecay","Z-boson Electron eta vs. eta",nptbins,ptmin,ptmax,netabins,etamin,etamax);\r
 432     fhPtAll = new TH1F("hPtAll","All electron pT",nptbins,ptmin,ptmax);\r
 433     fhPhiAll = new TH2F("hPhiAll","All Electron phi vs pT",nptbins,ptmin,ptmax,nphibins,phimin,phimax);\r
 434     fhEtaAll = new TH2F("hEtaAll","All Electron eta vs. eta",nptbins,ptmin,ptmax,netabins,etamin,etamax);\r
 435     fhPtUnknown = new TH1F("hPtUnknown","Unknown electron pT",nptbins,ptmin,ptmax);\r
 436     fhPhiUnknown = new TH2F("hPhiUnknown","Unknown Electron phi vs pT",nptbins,ptmin,ptmax,nphibins,phimin,phimax);\r
 437     fhEtaUnknown = new TH2F("hEtaUnknown","Unknown Electron eta vs. eta",nptbins,ptmin,ptmax,netabins,etamin,etamax);\r
 438     fhPtMisidentified = new TH1F("hPtMisidentified","Misidentified electron pT",nptbins,ptmin,ptmax);\r
 439     fhPhiMisidentified = new TH2F("hPhiMisidentified","Misidentified Electron phi vs pT",nptbins,ptmin,ptmax,nphibins,phimin,phimax);\r
 440     fhEtaMisidentified = new TH2F("hEtaMisidentified","Misidentified Electron eta vs. eta",nptbins,ptmin,ptmax,netabins,etamin,etamax);\r
 441 \r
 442     outputContainer->Add(fhPtConversion);\r
 443     outputContainer->Add(fhPhiConversion);\r
 444     outputContainer->Add(fhEtaConversion);\r
 445     outputContainer->Add(fhPtBottom);\r
 446     outputContainer->Add(fhPhiBottom);\r
 447     outputContainer->Add(fhEtaBottom);\r
 448     outputContainer->Add(fhPtCharm);\r
 449     outputContainer->Add(fhPhiCharm);\r
 450     outputContainer->Add(fhEtaCharm);\r
 451     outputContainer->Add(fhPtCFromB);\r
 452     outputContainer->Add(fhPhiCFromB);\r
 453     outputContainer->Add(fhEtaCFromB);\r
 454     outputContainer->Add(fhPtDalitz);\r
 455     outputContainer->Add(fhPhiDalitz);\r
 456     outputContainer->Add(fhEtaDalitz);\r
 457     outputContainer->Add(fhPtWDecay);\r
 458     outputContainer->Add(fhPhiWDecay);\r
 459     outputContainer->Add(fhEtaWDecay);\r
 460     outputContainer->Add(fhPtZDecay);\r
 461     outputContainer->Add(fhPhiZDecay);\r
 462     outputContainer->Add(fhEtaZDecay);\r
 463     outputContainer->Add(fhPtAll);\r
 464     outputContainer->Add(fhPhiAll);\r
 465     outputContainer->Add(fhEtaAll);\r
 466     outputContainer->Add(fhPtUnknown);\r
 467     outputContainer->Add(fhPhiUnknown);\r
 468     outputContainer->Add(fhEtaUnknown);\r
 469     outputContainer->Add(fhPtMisidentified);\r
 470     outputContainer->Add(fhPhiMisidentified);\r
 471     outputContainer->Add(fhEtaMisidentified);\r
 472     \r
 473     //created ele ntuple for further analysis\r
 474     if(fWriteNtuple) {\r
 475       fMCEleNtuple = new TNtuple("MCEleNtuple","MC Electron Ntuple","mctag:pt:phi:eta:x:y:z");\r
 476       outputContainer->Add(fMCEleNtuple) ;\r
 477     }\r
 478 \r
 479     fhPtMCHadron = new TH1F("hPtMCHadron","MC Charged hadrons w/in EMCAL acceptance",nptbins,ptmin,ptmax);\r
 480     fhPtMCBottom = new TH1F("hPtMCBottom","MC Bottom electrons w/in EMCAL acceptance",nptbins,ptmin,ptmax);\r
 481     fhPtMCCharm = new TH1F("hPtMCCharm","MC Charm electrons w/in EMCAL acceptance",nptbins,ptmin,ptmax);\r
 482     fhPtMCCFromB = new TH1F("hPtMCCFromB","MC Charm from Bottom electrons w/in EMCAL acceptance",nptbins,ptmin,ptmax);\r
 483     fhPtMCConversion = new TH1F("hPtMCConversion","MC Conversion electrons w/in EMCAL acceptance",nptbins,ptmin,ptmax);\r
 484     fhPtMCDalitz  = new TH1F("hPtMCDalitz","MC Dalitz electrons w/in EMCAL acceptance",nptbins,ptmin,ptmax);\r
 485     fhPtMCWDecay = new TH1F("hPtMCWDecay","MC W Decay electrons w/in EMCAL acceptance",nptbins,ptmin,ptmax);\r
 486     fhPtMCZDecay = new TH1F("hPtMCZDecay","MC Z Decay electrons w/in EMCAL acceptance",nptbins,ptmin,ptmax);\r
 487     fhPtMCUnknown = new TH1F("hPtMCUnknown","MC Unknown electrons w/in EMCAL acceptance",nptbins,ptmin,ptmax);\r
 488 \r
 489     outputContainer->Add(fhPtMCHadron);\r
 490     outputContainer->Add(fhPtMCBottom);\r
 491     outputContainer->Add(fhPtMCCharm);\r
 492     outputContainer->Add(fhPtMCCFromB);\r
 493     outputContainer->Add(fhPtMCConversion);\r
 494     outputContainer->Add(fhPtMCDalitz);\r
 495     outputContainer->Add(fhPtMCWDecay);\r
 496     outputContainer->Add(fhPtMCZDecay);\r
 497     outputContainer->Add(fhPtMCUnknown);\r
 498 \r
 499   }//Histos with MC\r
 500   \r
 501   //Save parameters used for analysis\r
 502   TString parList ; //this will be list of parameters used for this analysis.\r
 503   char onePar[500] ;\r
 504   \r
 505   sprintf(onePar,"--- AliAnaElectron ---\n") ;\r
 506   parList+=onePar ;     \r
 507   sprintf(onePar,"fCalorimeter: %s\n",fCalorimeter.Data()) ;\r
 508   parList+=onePar ;  \r
 509   sprintf(onePar,"fpOverEmin: %f\n",fpOverEmin) ;\r
 510   parList+=onePar ;  \r
 511   sprintf(onePar,"fpOverEmax: %f\n",fpOverEmax) ;\r
 512   parList+=onePar ;  \r
 513   sprintf(onePar,"fResidualCut: %f\n",fResidualCut) ;\r
 514   parList+=onePar ;  \r
 515   sprintf(onePar,"---DVM Btagging\n");\r
 516   parList+=onePar ;\r
 517   sprintf(onePar,"max IP-cut (e,h): %f\n",fImpactCut);\r
 518   parList+=onePar ;\r
 519   sprintf(onePar,"min ITS-hits: %d\n",fITSCut);\r
 520   parList+=onePar ;\r
 521   sprintf(onePar,"max dR (e,h): %f\n",fDrCut);\r
 522   parList+=onePar ;\r
 523   sprintf(onePar,"max pairDCA: %f\n",fPairDcaCut);\r
 524   parList+=onePar ;\r
 525   sprintf(onePar,"max decaylength: %f\n",fDecayLenCut);\r
 526   parList+=onePar ;\r
 527   sprintf(onePar,"min Associated Pt: %f\n",fAssocPtCut);\r
 528   parList+=onePar ;\r
 529   sprintf(onePar,"---IPSig Btagging\n");\r
 530   parList+=onePar ;\r
 531   sprintf(onePar,"min tag track: %d\n",fNTagTrkCut);\r
 532   parList+=onePar ;\r
 533   sprintf(onePar,"min IP significance: %f\n",fIPSigCut);\r
 534   parList+=onePar ;\r
 535 \r
 536   //Get parameters set in base class.\r
 537   parList += GetBaseParametersList() ;\r
 538   \r
 539   //Get parameters set in FidutialCut class (not available yet)\r
 540   //parlist += GetFidCut()->GetFidCutParametersList() \r
 541   \r
 542   TObjString *oString= new TObjString(parList) ;\r
 543   outputContainer->Add(oString);\r
 544   \r
 545   return outputContainer ;\r
 546   \r
 547 }\r
 548 \r
 549 //____________________________________________________________________________\r
 550 void AliAnaElectron::Init()\r
 551 {\r
 552 \r
 553   //do some initialization\r
 554   if(fCalorimeter == "PHOS") {\r
 555     printf("AliAnaElectron::Init() - !!STOP: You want to use PHOS in analysis but this is not (yet) supported!!\n!!Check the configuration file!!\n");\r
 556     fCalorimeter = "EMCAL";\r
 557   }\r
 558   if(fCalorimeter == "EMCAL" && !GetReader()->IsEMCALSwitchedOn()){\r
 559     printf("AliAnaElectron::Init() - !!STOP: You want to use EMCAL in analysis but it is not read!!\n!!Check the configuration file!!\n");\r
 560     abort();\r
 561   }\r
 562 \r
 563 }\r
 564 \r
 565 \r
 566 //____________________________________________________________________________\r
 567 void AliAnaElectron::InitParameters()\r
 568 {\r
 569   \r
 570   //Initialize the parameters of the analysis.\r
 571   SetOutputAODClassName("AliAODPWG4Particle");\r
 572   SetOutputAODName("PWG4Particle");\r
 573 \r
 574   AddToHistogramsName("AnaElectron_");\r
 575 \r
 576   fCalorimeter = "EMCAL" ;\r
 577   fpOverEmin = 0.5;\r
 578   fpOverEmax = 1.5;\r
 579   fResidualCut = 0.02;\r
 580   //DVM B-tagging\r
 581   fDrCut       = 1.0; \r
 582   fPairDcaCut  = 0.02;\r
 583   fDecayLenCut = 1.0;\r
 584   fImpactCut   = 0.5;\r
 585   fAssocPtCut  = 1.0;\r
 586   fMassCut     = 1.5;\r
 587   fSdcaCut     = 0.1;\r
 588   fITSCut      = 4;\r
 589   //IPSig B-tagging\r
 590   fNTagTrkCut  = 2;\r
 591   fIPSigCut    = 3.0;\r
 592 }\r
 593 \r
 594 //__________________________________________________________________\r
 595 void  AliAnaElectron::MakeAnalysisFillAOD() \r
 596 {\r
 597   //\r
 598   // Do analysis and fill aods with electron candidates\r
 599   // These AODs will be used to do subsequent histogram filling\r
 600   //\r
 601   // Also fill some QA histograms\r
 602   //\r
 603 \r
 604   TObjArray *cl = new TObjArray();\r
 605 \r
 606   Double_t bfield = 0.;\r
 607   if(GetReader()->GetDataType() != AliCaloTrackReader::kMC) bfield = GetReader()->GetBField();\r
 608 \r
 609   //Select the calorimeter of the electron\r
 610   if(fCalorimeter != "EMCAL") {\r
 611     printf("This class not yet implemented for PHOS\n");\r
 612     abort();\r
 613   }\r
 614   cl = GetAODEMCAL();\r
 615   \r
 616   ////////////////////////////////////////////////\r
 617   //Start from tracks and get associated clusters \r
 618   ////////////////////////////////////////////////\r
 619   if(!GetAODCTS() || GetAODCTS()->GetEntriesFast() == 0) return ;\r
 620   Int_t ntracks = GetAODCTS()->GetEntriesFast();\r
 621   Int_t refmult = 0; Int_t refmult2 = 0;\r
 622   if(GetDebug() > 0)\r
 623     printf("AliAnaElectron::MakeAnalysisFillAOD() - In CTS aod entries %d\n", ntracks);\r
 624 \r
 625   //Unfortunately, AliAODTracks don't have associated EMCAL clusters.\r
 626   //we have to redo track-matching, I guess\r
 627   Int_t iCluster = -999;\r
 628   Int_t bt = 0; //counter for event b-tags\r
 629 \r
 630   for (Int_t itrk =  0; itrk <  ntracks; itrk++) {////////////// track loop\r
 631     iCluster = -999; //start with no match\r
 632     AliAODTrack * track = (AliAODTrack*) (GetAODCTS()->At(itrk)) ;\r
 633     if (TMath::Abs(track->Eta())< 0.5) refmult++;\r
 634     Double_t imp[2] = {-999.,-999.}; Double_t cov[3] = {-999.,-999.,-999.};\r
 635     Bool_t dcaOkay = GetDCA(track,imp,cov);  //homegrown dca calculation until AOD is fixed\r
 636     if(!dcaOkay) printf("AliAnaElectron::Problem computing DCA to primary vertex for track %d.  Skipping it...\n",itrk);\r
 637     if(TMath::Abs(track->Eta())< 0.5 && TMath::Abs(imp[0])<1.0 && TMath::Abs(imp[1])<1.0) refmult2++;\r
 638     fhImpactXY->Fill(imp[0]);\r
 639 \r
 640     AliAODPid* pid = (AliAODPid*) track->GetDetPid();\r
 641     if(pid == 0) {\r
 642       if(GetDebug() > 0) printf("AliAnaElectron::MakeAnalysisFillAOD() - No PID object - skipping track %d",itrk);\r
 643       continue;\r
 644     } else {\r
 645       Double_t emcpos[3];\r
 646       pid->GetEMCALPosition(emcpos);\r
 647       Double_t emcmom[3];\r
 648       pid->GetEMCALMomentum(emcmom);\r
 649       \r
 650       TVector3 pos(emcpos[0],emcpos[1],emcpos[2]);\r
 651       TVector3 mom(emcmom[0],emcmom[1],emcmom[2]);\r
 652       Double_t tphi = pos.Phi();\r
 653       Double_t teta = pos.Eta();\r
 654       Double_t tmom = mom.Mag();\r
 655       \r
 656       TLorentzVector mom2(mom,0.);\r
 657       Bool_t in =  GetFidutialCut()->IsInFidutialCut(mom2,fCalorimeter) ;\r
 658       if(GetDebug() > 1) printf("AliAnaElectron::MakeAnalysisFillAOD() - Track pt %2.2f, phi %2.2f, eta %2.2f in fidutial cut %d\n",track->Pt(), track->Phi(), track->Eta(), in);\r
 659       if(mom.Pt() > GetMinPt() && in) {\r
 660         \r
 661         Double_t dEdx = pid->GetTPCsignal();\r
 662                 \r
 663         Int_t ntot = cl->GetEntriesFast();\r
 664         Double_t res = 999.;\r
 665         Double_t pOverE = -999.;\r
 666         \r
 667         Int_t pidProb = track->GetMostProbablePID();\r
 668         if(pidProb == AliAODTrack::kPion || pidProb == AliAODTrack::kKaon || pidProb == AliAODTrack::kProton) fhPtHadron->Fill(track->Pt());\r
 669         if(pidProb == AliAODTrack::kElectron) fhPtEleTrkDet->Fill(track->Pt());\r
 670         \r
 671         Bool_t isElectron = kFALSE;      \r
 672         //For tracks in EMCAL acceptance, pair them with all clusters\r
 673         //and fill the dEta vs dPhi for these pairs:\r
 674         for(Int_t iclus = 0; iclus < ntot; iclus++) {\r
 675           AliAODCaloCluster * clus = (AliAODCaloCluster*) (cl->At(iclus));\r
 676           if(!clus) continue;\r
 677           \r
 678           Double_t x[3];\r
 679           clus->GetPosition(x);\r
 680           TVector3 cluspos(x[0],x[1],x[2]);\r
 681           Double_t deta = teta - cluspos.Eta();\r
 682           Double_t dphi = tphi - cluspos.Phi();\r
 683           if(dphi > TMath::Pi()) dphi -= 2*TMath::Pi();\r
 684           if(dphi < -TMath::Pi()) dphi += 2*TMath::Pi();\r
 685           fh2dEtadPhi->Fill(deta,dphi);\r
 686           fh2TrackPVsClusterE->Fill(clus->E(),track->P());\r
 687           fh2TrackPtVsClusterE->Fill(clus->E(),track->Pt());\r
 688           fh2TrackPhiVsClusterPhi->Fill(cluspos.Phi(),mom.Phi());\r
 689           fh2TrackEtaVsClusterEta->Fill(cluspos.Eta(),mom.Eta());\r
 690           \r
 691           res = sqrt(dphi*dphi + deta*deta);\r
 692           fh1dR->Fill(res);\r
 693           \r
 694           /////////////////////////////////\r
 695           //Perform electron cut analysis//\r
 696           /////////////////////////////////\r
 697           //Good match\r
 698           if(res < fResidualCut) {\r
 699             fh2dEtadPhiMatched->Fill(deta,dphi);\r
 700             iCluster = iclus;\r
 701             \r
 702             Int_t tmctag = -1;\r
 703             Int_t cmctag = -1;\r
 704             \r
 705             if(IsDataMC()) {\r
 706               //Input from second AOD?\r
 707               Int_t input = 0;\r
 708               if(GetReader()->GetAODCTSNormalInputEntries() <= itrk) input = 1;\r
 709               tmctag = GetMCAnalysisUtils()->CheckOrigin(track->GetLabel(),GetReader(),input);\r
 710               \r
 711               //Do you want the cluster or the track label?\r
 712               input = 0;\r
 713               if(GetReader()->GetAODEMCALNormalInputEntries() <= iclus) input = 1;\r
 714               cmctag = GetMCAnalysisUtils()->CheckOrigin(clus->GetLabel(0),GetReader(),input);\r
 715             }\r
 716             \r
 717             if(fWriteNtuple) {\r
 718               fEleNtuple->Fill(tmctag,cmctag,track->Pt(),track->Phi(),track->Eta(),track->P(),clus->E(),deta,dphi,clus->GetNCells(),dEdx,pidProb,imp[0],imp[1]);\r
 719             }\r
 720             \r
 721             fh2MatchdEdx->Fill(track->P(),dEdx);\r
 722             \r
 723             Double_t energy = clus->E(); \r
 724             if(energy > 0) pOverE = tmom/energy;\r
 725             fh1pOverE->Fill(pOverE);\r
 726             \r
 727             Int_t mult = clus->GetNCells();\r
 728             if(mult < 2 &&  GetDebug() > 0) printf("Single digit cluster.\n");\r
 729             \r
 730             //////////////////////////////\r
 731             //Electron cuts happen here!//\r
 732             //////////////////////////////\r
 733             if(pOverE > fpOverEmin && pOverE < fpOverEmax) isElectron = kTRUE;\r
 734           } else {\r
 735             fh2dEtadPhiUnmatched->Fill(deta,dphi);\r
 736           }\r
 737           \r
 738         } //calocluster loop\r
 739         \r
 740         ///////////////////////////\r
 741         //Fill AOD with electrons//\r
 742         ///////////////////////////\r
 743         if(isElectron) {\r
 744           \r
 745           //B-tagging\r
 746           if(GetDebug() > 1) printf("Found Electron - do b-tagging\n");\r
 747           Int_t dvmbtag = GetDVMBtag(track); bt += dvmbtag;\r
 748 \r
 749           fh2EledEdx->Fill(track->P(),dEdx);\r
 750           \r
 751           Double_t eMass = 0.511/1000; //mass in GeV\r
 752           Double_t eleE = sqrt(track->P()*track->P() + eMass*eMass);\r
 753           AliAODPWG4Particle tr = AliAODPWG4Particle(track->Px(),track->Py(),track->Pz(),eleE);\r
 754           tr.SetLabel(track->GetLabel());\r
 755           tr.SetCaloLabel(iCluster,-1); //sets the indices of the original caloclusters\r
 756           tr.SetTrackLabel(itrk,-1); //sets the indices of the original tracks\r
 757           tr.SetDetector(fCalorimeter);\r
 758           if(GetReader()->GetAODCTSNormalInputEntries() <= itrk) tr.SetInputFileIndex(1);\r
 759           //Make this preserve sign of particle\r
 760           if(track->Charge() < 0) tr.SetPdg(11); //electron is 11\r
 761           else  tr.SetPdg(-11); //positron is -11\r
 762           Int_t btag = -1;\r
 763           if(dvmbtag > 0) tr.SetBTagBit(btag,tr.kDVMTag0);\r
 764           if(dvmbtag > 1) tr.SetBTagBit(btag,tr.kDVMTag1);\r
 765           if(dvmbtag > 2) tr.SetBTagBit(btag,tr.kDVMTag2);\r
 766           tr.SetBtag(btag);\r
 767           \r
 768           //Play with the MC stack if available\r
 769           //Check origin of the candidates\r
 770           if(IsDataMC()){\r
 771             \r
 772             //FIXME:  Need to re-think this for track-oriented analysis\r
 773             //JLK DO WE WANT TRACK TAG OR CLUSTER TAG?\r
 774             tr.SetTag(GetMCAnalysisUtils()->CheckOrigin(tr.GetLabel(),GetReader(),tr.GetInputFileIndex()));\r
 775             \r
 776             if(GetDebug() > 0) printf("AliAnaElectron::MakeAnalysisFillAOD() - Origin of candidate %d\n",tr.GetTag());\r
 777           }//Work with stack also   \r
 778           \r
 779           AddAODParticle(tr);\r
 780           \r
 781           if(GetDebug() > 1) printf("AliAnaElectron::MakeAnalysisFillAOD() - Electron selection cuts passed: pT %3.2f, pdg %d\n",tr.Pt(),tr.GetPdg());  \r
 782         }//electron\r
 783       }//pt, fiducial selection\r
 784     }//pid check\r
 785   }//track loop                         \r
 786   \r
 787   //FIXME:  Should we also check from the calocluster side, just in\r
 788   //case?\r
 789   fhRefMult->Fill(refmult);\r
 790   fhRefMult2->Fill(refmult2);\r
 791 \r
 792   if(GetDebug() > 1 && bt > 0) printf("AliAnaElectron::MakeAnalysisFillAOD() *** Event Btagged *** \n");\r
 793   if(GetDebug() > 1) printf("AliAnaElectron::MakeAnalysisFillAOD()  End fill AODs \n");  \r
 794   \r
 795 }\r
 796 \r
 797 //__________________________________________________________________\r
 798 void  AliAnaElectron::MakeAnalysisFillHistograms() \r
 799 {\r
 800   //Do analysis and fill histograms\r
 801 \r
 802   AliStack * stack = 0x0;\r
 803   TParticle * primary = 0x0;\r
 804   TClonesArray * mcparticles0 = 0x0;\r
 805   TClonesArray * mcparticles1 = 0x0;\r
 806   AliAODMCParticle * aodprimary = 0x0;\r
 807 \r
 808   Int_t ph1 = 0;  //photonic 1 count\r
 809   Int_t ph2 = 0;  //photonic 2 count\r
 810   Int_t phB = 0;  //both count\r
 811 \r
 812   if(IsDataMC()) {\r
 813     if(GetReader()->ReadStack()){\r
 814       stack =  GetMCStack() ;\r
 815       \r
 816       if(!stack)\r
 817         printf("AliAnaElectron::MakeAnalysisFillHistograms() *** no stack ***: \n");\r
 818       \r
 819     }\r
 820     else if(GetReader()->ReadAODMCParticles()){\r
 821       //Get the list of MC particles\r
 822       mcparticles0 = GetReader()->GetAODMCParticles(0);\r
 823       if(!mcparticles0 && GetDebug() > 0)     {\r
 824         printf("AliAnaElectron::MakeAnalysisFillHistograms() -  Standard MCParticles not available!\n");\r
 825       }\r
 826       if(GetReader()->GetSecondInputAODTree()){\r
 827         mcparticles1 = GetReader()->GetAODMCParticles(1);\r
 828         if(!mcparticles1 && GetDebug() > 0)     {\r
 829           printf("AliAnaElectron::MakeAnalysisFillHistograms() -  Second input MCParticles not available!\n");\r
 830         }\r
 831       }\r
 832       \r
 833     }\r
 834   }// is data and MC\r
 835 \r
 836   ////////////////////////////////////\r
 837   //Loop over jets and check for b-tag\r
 838   ////////////////////////////////////\r
 839   Int_t njets = (GetReader()->GetOutputEvent())->GetNJets();\r
 840   if(njets > 0) {\r
 841     if(GetDebug() > 0) printf("AliAnaElectron::MakeAnalysisFillHistograms() - Jet AOD branch has %d jets.  Performing b-jet tag analysis",njets);\r
 842     \r
 843     for(Int_t ijet = 0; ijet < njets ; ijet++) {\r
 844       AliAODJet * jet = (AliAODJet*)(GetReader()->GetOutputEvent())->GetJet(ijet) ;\r
 845       if(GetDebug() > 3) {\r
 846         printf("AliAODJet ijet = %d\n",ijet);\r
 847         jet->Print("");\r
 848       }\r
 849       //To "tag" the jet, we will look for it to pass our various criteria\r
 850       //For e jet tag, we just look to see which ones have NPEs\r
 851       //For DVM jet tag, we will look for DVM electrons\r
 852       //For IPSig, we compute the IPSig for all tracks and if the\r
 853       //number passing is above the cut, it passes\r
 854       Bool_t eJet = kFALSE;  \r
 855       Bool_t dvmJet = kFALSE;  \r
 856       Bool_t ipsigJet = kFALSE;\r
 857       TRefArray* rt = jet->GetRefTracks();\r
 858       Int_t ntrk = rt->GetEntries();\r
 859       Int_t trackCounter = 0; //for ipsig\r
 860       for(Int_t itrk = 0; itrk < ntrk; itrk++) {\r
 861         AliAODTrack* jetTrack = (AliAODTrack*)jet->GetTrack(itrk);\r
 862         if( GetIPSignificance(jetTrack, jet->Phi()) > fIPSigCut) trackCounter++;\r
 863         Bool_t isNPE = CheckTrack(jetTrack,"NPE");\r
 864         if(isNPE) eJet = kTRUE;\r
 865         Bool_t isDVM = CheckTrack(jetTrack,"DVM");\r
 866         if(isDVM) dvmJet = kTRUE;\r
 867       }\r
 868       fhIPSigBtagQA1->Fill(trackCounter);\r
 869       if(trackCounter > fNTagTrkCut) ipsigJet = kTRUE;\r
 870 \r
 871       //Fill bjet histograms here\r
 872       if(!(eJet || ipsigJet || dvmJet)) fhJetType->Fill(0.,jet->Pt()); //none\r
 873       if(eJet && !(ipsigJet || dvmJet)) fhJetType->Fill(1.,jet->Pt()); //only ejet\r
 874       if(dvmJet && !(eJet || ipsigJet)) fhJetType->Fill(2.,jet->Pt()); //only dvm\r
 875       if(ipsigJet && !(eJet || dvmJet)) fhJetType->Fill(3.,jet->Pt()); //only ipsig\r
 876       if(eJet && dvmJet && !ipsigJet)   fhJetType->Fill(4.,jet->Pt()); //ejet & dvm\r
 877       if(eJet && ipsigJet && !dvmJet)   fhJetType->Fill(5.,jet->Pt()); //ejet & ipsig\r
 878       if(dvmJet && ipsigJet && !eJet)   fhJetType->Fill(6.,jet->Pt()); //dvm & ipsig\r
 879       if(dvmJet && ipsigJet && eJet)    fhJetType->Fill(7.,jet->Pt()); //all\r
 880 \r
 881       for(Int_t itrk = 0; itrk < ntrk; itrk++) {\r
 882         AliAODTrack* jetTrack = (AliAODTrack*)jet->GetTrack(itrk);\r
 883         Double_t xsi = TMath::Log(jet->Pt()/jetTrack->Pt());\r
 884         if(eJet || ipsigJet || dvmJet) {\r
 885           if(GetDebug() > 0) printf("AliAnaElectron::MakeAnalysisFillHistograms - We have a bjet!\n");\r
 886 \r
 887           fhBJetXsiFF->Fill(xsi,jet->Pt());\r
 888           fhBJetPtFF->Fill(jetTrack->Pt(),jet->Pt());\r
 889           fhBJetEtaPhi->Fill(jet->Eta(),jet->Phi());\r
 890         } else {\r
 891           //Fill non-bjet histograms here\r
 892           fhNonBJetXsiFF->Fill(xsi,jet->Pt());\r
 893           fhNonBJetPtFF->Fill(jetTrack->Pt(),jet->Pt());\r
 894           fhNonBJetEtaPhi->Fill(jet->Eta(),jet->Phi());\r
 895         }\r
 896       }\r
 897     } //jet loop\r
 898   } //jets exist\r
 899   \r
 900   //////////////////////////////\r
 901   //Loop on stored AOD electrons\r
 902   //////////////////////////////\r
 903   Int_t naod = GetOutputAODBranch()->GetEntriesFast();\r
 904   if(GetDebug() > 0) printf("AliAnaElectron::MakeAnalysisFillHistograms() - aod branch entries %d\n", naod);\r
 905   \r
 906   for(Int_t iaod = 0; iaod < naod ; iaod++){\r
 907     AliAODPWG4Particle* ele =  (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));\r
 908     Int_t pdg = ele->GetPdg();\r
 909     \r
 910     if(GetDebug() > 3) \r
 911       printf("AliAnaElectron::MakeAnalysisFillHistograms() - PDG %d, MC TAG %d, Calorimeter %s\n", ele->GetPdg(),ele->GetTag(), (ele->GetDetector()).Data()) ;\r
 912     \r
 913     if(TMath::Abs(pdg) != AliCaloPID::kElectron) continue; \r
 914     if(ele->GetDetector() != fCalorimeter) continue;\r
 915     \r
 916     if(GetDebug() > 1) \r
 917       printf("AliAnaElectron::MakeAnalysisFillHistograms() - ID Electron: pt %f, phi %f, eta %f\n", ele->Pt(),ele->Phi(),ele->Eta()) ;\r
 918 \r
 919     //Filter for photonic electrons based on opening angle and Minv\r
 920     //cuts, also fill histograms\r
 921     Bool_t photonic = kFALSE;\r
 922     Bool_t photonic1 = kFALSE;\r
 923     photonic1 = IsItPhotonic(ele); //check against primaries\r
 924     if(photonic1) ph1++;\r
 925     Bool_t photonic2 = kFALSE;\r
 926     photonic2 = IsItPhotonic2(ele); //check against V0s\r
 927     if(photonic2) ph2++;\r
 928     if(photonic1 && photonic2) phB++;\r
 929     if(photonic1 || photonic2) photonic = kTRUE;\r
 930 \r
 931     //Fill electron histograms \r
 932     Float_t ptele = ele->Pt();\r
 933     Float_t phiele = ele->Phi();\r
 934     Float_t etaele = ele->Eta();\r
 935     \r
 936     fhPtElectron  ->Fill(ptele);\r
 937     fhPhiElectron ->Fill(ptele,phiele);\r
 938     fhEtaElectron ->Fill(ptele,etaele);\r
 939 \r
 940     if(photonic) {\r
 941       fhPtPE->Fill(ptele);\r
 942       fhPhiPE->Fill(ptele,phiele);\r
 943       fhEtaPE->Fill(ptele,etaele);\r
 944     } else {\r
 945       fhPtNPE->Fill(ptele);\r
 946       fhPhiNPE->Fill(ptele,phiele);\r
 947       fhEtaNPE->Fill(ptele,etaele);\r
 948     }\r
 949 \r
 950     if(IsDataMC()){\r
 951       Int_t tag = ele->GetTag();\r
 952       if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron)) {\r
 953         fhPtAll  ->Fill(ptele);\r
 954         fhPhiAll ->Fill(ptele,phiele);\r
 955         fhEtaAll ->Fill(ptele,etaele);    \r
 956 \r
 957         if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion)){\r
 958           fhPtConversion  ->Fill(ptele);\r
 959           fhPhiConversion ->Fill(ptele,phiele);\r
 960           fhEtaConversion ->Fill(ptele,etaele);\r
 961         }\r
 962         else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEFromB)){\r
 963           fhPtBottom  ->Fill(ptele);\r
 964           fhPhiBottom ->Fill(ptele,phiele);\r
 965           fhEtaBottom ->Fill(ptele,etaele);\r
 966         }\r
 967         else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEFromC) && !GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEFromB)){\r
 968           fhPtCharm  ->Fill(ptele);\r
 969           fhPhiCharm ->Fill(ptele,phiele);\r
 970           fhEtaCharm ->Fill(ptele,etaele);\r
 971         }\r
 972         else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEFromCFromB)){\r
 973           fhPtCFromB  ->Fill(ptele);\r
 974           fhPhiCFromB ->Fill(ptele,phiele);\r
 975           fhEtaCFromB ->Fill(ptele,etaele);\r
 976         }\r
 977         else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay) || GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay) || GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay)) {\r
 978           fhPtDalitz  ->Fill(ptele);\r
 979           fhPhiDalitz ->Fill(ptele,phiele);\r
 980           fhEtaDalitz ->Fill(ptele,etaele);\r
 981         }\r
 982         else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCWDecay)) {\r
 983           fhPtWDecay  ->Fill(ptele);\r
 984           fhPhiWDecay ->Fill(ptele,phiele);\r
 985           fhEtaWDecay ->Fill(ptele,etaele);\r
 986         }\r
 987         else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCZDecay)) {\r
 988           fhPtZDecay  ->Fill(ptele);\r
 989           fhPhiZDecay ->Fill(ptele,phiele);\r
 990           fhEtaZDecay ->Fill(ptele,etaele);\r
 991         } \r
 992         else {\r
 993           fhPtUnknown  ->Fill(ptele);\r
 994           fhPhiUnknown ->Fill(ptele,phiele);\r
 995           fhEtaUnknown ->Fill(ptele,etaele);\r
 996         }\r
 997         //not electron:\r
 998       } else {\r
 999         fhPtMisidentified  ->Fill(ptele);\r
1000         fhPhiMisidentified ->Fill(ptele,phiele);\r
1001         fhEtaMisidentified ->Fill(ptele,etaele);\r
1002       }\r
1003     }//Histograms with MC\r
1004     \r
1005   }// aod loop\r
1006 \r
1007   ////////////////////////////////////////////////////////\r
1008   //Fill histograms of pure MC kinematics from the stack//\r
1009   ////////////////////////////////////////////////////////\r
1010   if(IsDataMC()) {\r
1011     if(GetReader()->ReadStack()) {\r
1012       for(Int_t ipart = 0; ipart < stack->GetNtrack(); ipart++) {\r
1013         primary = stack->Particle(ipart);\r
1014         TLorentzVector mom;\r
1015         primary->Momentum(mom);\r
1016         Bool_t in = GetFidutialCut()->IsInFidutialCut(mom,fCalorimeter);\r
1017         if(primary->Pt() < GetMinPt()) continue;\r
1018         if(!in) continue;\r
1019 \r
1020         Int_t pdgcode = primary->GetPdgCode();\r
1021         if(TMath::Abs(pdgcode) == 211 || TMath::Abs(pdgcode) == 321 || TMath::Abs(pdgcode) == 2212)\r
1022           fhPtMCHadron->Fill(primary->Pt());\r
1023 \r
1024         //we only care about electrons\r
1025         if(TMath::Abs(pdgcode) != 11) continue;\r
1026         //we only want TRACKABLE electrons (TPC 85-250cm)\r
1027         if(primary->R() > 200.) continue;\r
1028         //Ignore low pt electrons\r
1029         if(primary->Pt() < 0.2) continue;\r
1030         \r
1031         //find out what the ancestry of this electron is\r
1032         Int_t mctag = -1;\r
1033         Int_t input = 0;\r
1034         mctag = GetMCAnalysisUtils()->CheckOrigin(ipart,GetReader(),input);\r
1035 \r
1036         //fill ntuple\r
1037         if(fWriteNtuple) {\r
1038           fMCEleNtuple->Fill(mctag,primary->Pt(),primary->Phi(),primary->Eta(),primary->Vx(),primary->Vy(),primary->Vz());\r
1039         }\r
1040         if(!GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCConversion)) {\r
1041           if(GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCEFromB))\r
1042             fhPtMCBottom->Fill(primary->Pt());\r
1043           else if(GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCEFromC))\r
1044             fhPtMCCharm->Fill(primary->Pt());\r
1045           else if(GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCEFromCFromB))\r
1046             fhPtMCCFromB->Fill(primary->Pt());\r
1047           else if(GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCWDecay))\r
1048             fhPtMCWDecay->Fill(primary->Pt());\r
1049           else if(GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCZDecay))\r
1050             fhPtMCZDecay->Fill(primary->Pt());\r
1051           else if(GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCPi0Decay) \r
1052                   || GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCEtaDecay)\r
1053                   || GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCOtherDecay))\r
1054             fhPtMCDalitz->Fill(primary->Pt());\r
1055           else \r
1056             fhPtMCUnknown->Fill(primary->Pt());\r
1057         } else if(GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCConversion)) {\r
1058           fhPtMCConversion->Fill(primary->Pt());\r
1059         }\r
1060       }      \r
1061     } else if(GetReader()->ReadAODMCParticles()) {\r
1062       Int_t npart0 = mcparticles0->GetEntriesFast();\r
1063       Int_t npart1 = 0;\r
1064       if(mcparticles1) npart1 = mcparticles1->GetEntriesFast();\r
1065       Int_t npart = npart0+npart1;\r
1066       for(Int_t ipart = 0; ipart < npart; ipart++) {\r
1067         if(ipart < npart0) aodprimary = (AliAODMCParticle*)mcparticles0->At(ipart);\r
1068         else aodprimary = (AliAODMCParticle*)mcparticles1->At(ipart-npart0);\r
1069         if(!aodprimary) {\r
1070           printf("AliAnaElectron::MakeAnalysisFillHistograms() *** no primary ***:  label %d \n", ipart);\r
1071           continue;\r
1072         }\r
1073 \r
1074         Double_t mom[3] = {0.,0.,0.};\r
1075         aodprimary->PxPyPz(mom);\r
1076         TLorentzVector mom2(mom,0.);    \r
1077         Bool_t in = GetFidutialCut()->IsInFidutialCut(mom2,fCalorimeter);\r
1078         if(aodprimary->Pt() < GetMinPt()) continue;\r
1079         if(!in) continue;\r
1080 \r
1081         Int_t pdgcode = aodprimary->GetPdgCode();\r
1082         if(TMath::Abs(pdgcode) == 211 || TMath::Abs(pdgcode) == 321 || TMath::Abs(pdgcode) == 2212)\r
1083           fhPtMCHadron->Fill(aodprimary->Pt());\r
1084 \r
1085         //we only care about electrons\r
1086         if(TMath::Abs(pdgcode) != 11) continue;\r
1087         //we only want TRACKABLE electrons (TPC 85-250cm)\r
1088         Double_t radius = TMath::Sqrt(aodprimary->Xv()*aodprimary->Xv() + aodprimary->Yv()*aodprimary->Yv());\r
1089         if(radius > 200.) continue;\r
1090 \r
1091         //find out what the ancestry of this electron is\r
1092         Int_t mctag = -1;\r
1093         Int_t input = 0;\r
1094         Int_t ival = ipart;\r
1095         if(ipart > npart0) { ival -= npart0; input = 1;}\r
1096         mctag = GetMCAnalysisUtils()->CheckOrigin(ival,GetReader(),input);\r
1097         \r
1098         //fill ntuple\r
1099         if(fWriteNtuple) {\r
1100           fMCEleNtuple->Fill(mctag,aodprimary->Pt(),aodprimary->Phi(),aodprimary->Eta(),aodprimary->Xv(),aodprimary->Yv(),aodprimary->Zv());\r
1101         }\r
1102         if(!GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCConversion)) {\r
1103           if(GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCEFromB))\r
1104             fhPtMCBottom->Fill(aodprimary->Pt());\r
1105           else if(GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCEFromC))\r
1106             fhPtMCCharm->Fill(aodprimary->Pt());\r
1107           else if(GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCEFromCFromB))\r
1108             fhPtMCCFromB->Fill(aodprimary->Pt());\r
1109           else if(GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCWDecay))\r
1110             fhPtMCWDecay->Fill(aodprimary->Pt());\r
1111           else if(GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCZDecay))\r
1112             fhPtMCZDecay->Fill(aodprimary->Pt());\r
1113           else if(GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCPi0Decay) \r
1114                   || GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCEtaDecay)\r
1115                   || GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCOtherDecay))\r
1116             fhPtMCDalitz->Fill(aodprimary->Pt());\r
1117           else \r
1118             fhPtMCUnknown->Fill(aodprimary->Pt());\r
1119         } else if(GetMCAnalysisUtils()->CheckTagBit(mctag,AliMCAnalysisUtils::kMCConversion)) {\r
1120           fhPtMCConversion->Fill(aodprimary->Pt());\r
1121         }       \r
1122       }\r
1123     }\r
1124   } //pure MC kine histos\r
1125     \r
1126   //if(GetDebug() > 0) \r
1127     printf("\tAliAnaElectron::Photonic electron counts: ph1 %d, ph2 %d, Both %d\n",ph1,ph2,phB);\r
1128 }\r
1129 \r
1130 //__________________________________________________________________\r
1131 Int_t AliAnaElectron::GetDVMBtag(AliAODTrack * tr )\r
1132 {\r
1133   //This method uses the Displaced Vertex between electron-hadron\r
1134   //pairs and the primary vertex to determine whether an electron is\r
1135   //likely from a B hadron.\r
1136 \r
1137   Int_t ncls1 = 0;\r
1138   for(Int_t l = 0; l < 6; l++) if(TESTBIT(tr->GetITSClusterMap(),l)) ncls1++;\r
1139 \r
1140   fhDVMBtagQA3->Fill(ncls1);\r
1141   if (ncls1 < fITSCut) return 0;\r
1142 \r
1143   Double_t imp[2] = {-999.,-999.}; Double_t cov[3] = {-999.,-999.,-999.};\r
1144   Bool_t dcaOkay = GetDCA(tr,imp,cov);  //homegrown dca calculation until AOD is fixed                  \r
1145   if(!dcaOkay) {\r
1146     printf("AliAnaElectron::Problem computing DCA to primary vertex for track %d",tr->GetID());\r
1147     return 0;\r
1148   }\r
1149 \r
1150   fhDVMBtagQA4->Fill(imp[0]);\r
1151   if (TMath::Abs(imp[0])   > fImpactCut ) return 0;\r
1152   fhDVMBtagQA5->Fill(imp[1]);\r
1153   if (TMath::Abs(imp[1])   > fImpactCut ) return 0;\r
1154 \r
1155   Int_t nvtx1 = 0;\r
1156   Int_t nvtx2 = 0;\r
1157   Int_t nvtx3 = 0;\r
1158 \r
1159   for (Int_t k2 =0; k2 < GetAODCTS()->GetEntriesFast() ; k2++) {\r
1160     //loop over assoc\r
1161     AliAODTrack* track2 = (AliAODTrack*) (GetAODCTS()->At(k2));\r
1162     Int_t id1 = tr->GetID();\r
1163     Int_t id2 = track2->GetID();\r
1164     if(id1 == id2) continue;\r
1165 \r
1166     Int_t ncls2 = 0;\r
1167     for(Int_t l = 0; l < 6; l++) if(TESTBIT(track2->GetITSClusterMap(),l)) ncls2++;\r
1168     if (ncls2 < fITSCut) continue;\r
1169 \r
1170     if(track2->Pt() < fAssocPtCut) continue;\r
1171 \r
1172     Double_t dphi = tr->Phi() - track2->Phi();\r
1173     if(dphi > TMath::Pi()) dphi -= 2*TMath::Pi();\r
1174     if(dphi < -TMath::Pi()) dphi += 2*TMath::Pi();\r
1175     Double_t deta = tr->Eta() - track2->Eta();\r
1176     Double_t dr = sqrt(deta*deta + dphi*dphi);\r
1177 \r
1178     if(dr > fDrCut) continue;\r
1179     \r
1180     Double_t sDca1 = ComputeSignDca(tr, track2, 1.0);\r
1181     if (sDca1 > fSdcaCut) nvtx1++;\r
1182     Double_t sDca2 = ComputeSignDca(tr, track2, 1.5);\r
1183     if (sDca2 > fSdcaCut) nvtx2++;\r
1184     Double_t sDca3 = ComputeSignDca(tr, track2, 1.8);\r
1185     if (sDca3 > fSdcaCut) nvtx3++;\r
1186 \r
1187   } //loop over hadrons\r
1188 \r
1189   if(GetDebug() > 0) {\r
1190     if (nvtx1>0) printf("result1 of btagging: %d \n",nvtx1);\r
1191     if (nvtx2>0) printf("result2 of btagging: %d \n",nvtx2);\r
1192     if (nvtx3>0) printf("result3 of btagging: %d \n",nvtx3);\r
1193   }\r
1194 \r
1195   //fill QA histograms\r
1196   fhDVMBtagCut1->Fill(nvtx1,tr->Pt());\r
1197   fhDVMBtagCut2->Fill(nvtx2,tr->Pt());\r
1198   fhDVMBtagCut3->Fill(nvtx3,tr->Pt());\r
1199 \r
1200   return nvtx2;\r
1201 \r
1202 }\r
1203 \r
1204 //__________________________________________________________________\r
1205 Double_t AliAnaElectron::ComputeSignDca(AliAODTrack *tr, AliAODTrack *tr2 , float masscut)\r
1206 {\r
1207   //Compute the signed dca between two tracks\r
1208   //and return the result\r
1209 \r
1210   Double_t signDca=-999.;\r
1211   if(GetDebug() > 2 ) printf(">>ComputeSdca:: track1 %d, track2 %d, masscut %f \n", tr->GetLabel(), tr2->GetLabel(), masscut);\r
1212 \r
1213   //=====Now calculate DCA between both tracks=======  \r
1214   Double_t massE = 0.000511;\r
1215   Double_t massK = 0.493677;\r
1216 \r
1217   Double_t bfield = 5.; //kG\r
1218   if(GetReader()->GetDataType() != AliCaloTrackReader::kMC) bfield = GetReader()->GetBField();\r
1219 \r
1220   Double_t vertex[3] = {-999.,-999.,-999}; //vertex\r
1221   if(GetReader()->GetDataType() != AliCaloTrackReader::kMC) {\r
1222     GetReader()->GetVertex(vertex); //If only one file, get the vertex from there\r
1223     //FIXME:  Add a check for whether file 2 is PYTHIA or HIJING\r
1224     //If PYTHIA, then set the vertex from file 2, if not, use the\r
1225     //vertex from file 1\r
1226     if(GetReader()->GetSecondInputAODTree()) GetReader()->GetSecondInputAODVertex(vertex);\r
1227   }\r
1228   \r
1229   TVector3 primV(vertex[0],vertex[1],vertex[2]) ;\r
1230 \r
1231   if(GetDebug() > 5) printf(">>ComputeSdca:: primary vertex = %2.2f,%2.2f,%2.2f \n",vertex[0],vertex[1],vertex[2]) ;\r
1232 \r
1233   AliExternalTrackParam *param1 = new AliExternalTrackParam(tr);\r
1234   AliExternalTrackParam *param2 = new AliExternalTrackParam(tr2);\r
1235 \r
1236   Double_t xplane1 = 0.; Double_t xplane2 = 0.;\r
1237   Double_t pairdca = param1->GetDCA(param2,bfield,xplane1,xplane2);\r
1238 \r
1239   Int_t id1 = 0, id2 = 0;\r
1240   AliESDv0 bvertex(*param1,id1,*param2,id2);\r
1241   Double_t vx,vy,vz;\r
1242   bvertex.GetXYZ(vx,vy,vz);\r
1243 \r
1244   Double_t emom[3];\r
1245   Double_t hmom[3];\r
1246   param1->PxPyPz(emom);\r
1247   param2->PxPyPz(hmom);\r
1248   TVector3 emomAtB(emom[0],emom[1],emom[2]);\r
1249   TVector3 hmomAtB(hmom[0],hmom[1],hmom[2]);\r
1250   TVector3 secvtxpt(vx,vy,vz);\r
1251   TVector3 decayvector(0,0,0);\r
1252   decayvector = secvtxpt - primV; //decay vector from PrimVtx\r
1253   Double_t decaylength = decayvector.Mag();\r
1254 \r
1255   if(GetDebug() > 0) {\r
1256     printf(">>ComputeSdca:: mom1=%f, mom2=%f \n", emomAtB.Perp(), hmomAtB.Perp() );\r
1257     printf(">>ComputeSdca:: pairDCA=%f, length=%f \n", pairdca,decaylength );\r
1258   }\r
1259 \r
1260   if (masscut<1.1) fhDVMBtagQA1->Fill(pairdca,decaylength);\r
1261 \r
1262   if (emomAtB.Mag()>0 && pairdca < fPairDcaCut && decaylength < fDecayLenCut ) {\r
1263     TVector3 sumMom = emomAtB+hmomAtB;\r
1264     Double_t ener1 = sqrt(pow(emomAtB.Mag(),2) + massE*massE);\r
1265     Double_t ener2 = sqrt(pow(hmomAtB.Mag(),2) + massK*massK);\r
1266     Double_t ener3 = sqrt(pow(hmomAtB.Mag(),2) + massE*massE);\r
1267     Double_t mass = sqrt(pow((ener1+ener2),2) - pow(sumMom.Mag(),2));\r
1268     Double_t massPhot = sqrt(pow((ener1+ener3),2) - pow(sumMom.Mag(),2));\r
1269     Double_t sDca = decayvector.Dot(emomAtB)/emomAtB.Mag();\r
1270 \r
1271     if (masscut<1.1) fhDVMBtagQA2->Fill(sDca, mass);\r
1272 \r
1273     if (mass > masscut && massPhot > 0.1) signDca = sDca;\r
1274     \r
1275     if(GetDebug() > 0) printf("\t>>ComputeSdca:: mass=%f \n", mass);\r
1276     if(GetDebug() > 0) printf("\t>>ComputeSdca:: sec vtx-signdca :%f\n",signDca);\r
1277   }\r
1278 \r
1279   //clean up\r
1280   delete param1;\r
1281   delete param2;\r
1282 \r
1283   return signDca;\r
1284 }\r
1285 \r
1286 //__________________________________________________________________\r
1287 Double_t AliAnaElectron::GetIPSignificance(AliAODTrack *tr, Double_t jetPhi)\r
1288 {\r
1289   //get signed impact parameter significance of the given AOD track\r
1290   //for the given jet\r
1291   Int_t trackIndex = 0;\r
1292   for (int k2 =0; k2 < GetAODCTS()->GetEntriesFast() ; k2++) {\r
1293     //loop over assoc\r
1294     AliAODTrack* track2 = (AliAODTrack*) (GetAODCTS()->At(k2));\r
1295     int id1 = tr->GetID();\r
1296     int id2 = track2->GetID();\r
1297     if(id1 == id2) {\r
1298       trackIndex = k2;//FIXME: check if GetAODCTS stores tracks in the\r
1299                       //same order of the event\r
1300       break;\r
1301     }\r
1302   }\r
1303   \r
1304   Double_t significance=0;\r
1305   Double_t magField = 0;\r
1306   Double_t maxD = 10000.;\r
1307   Double_t impPar[] = {0,0};\r
1308   Double_t ipCov[]={0,0,0};\r
1309   Double_t ipVec2D[] = {0,0};\r
1310 \r
1311   AliVEvent* vEvent = (AliVEvent*)GetReader()->GetInputEvent();\r
1312   if(!vEvent) return -97;\r
1313   AliVVertex* vv = (AliVVertex*)vEvent->GetPrimaryVertex();\r
1314   if(!vv) return -98;\r
1315   AliVTrack* vTrack = (AliVTrack*)vEvent->GetTrack(trackIndex);\r
1316   if(!vTrack) return -99;\r
1317   AliESDtrack esdTrack(vTrack);\r
1318   if(!esdTrack.PropagateToDCA(vv, magField, maxD, impPar, ipCov)) return -100;\r
1319   if(ipCov[0]<0) return -101;\r
1320 \r
1321   Double_t Pxy[] = {esdTrack.Px(), esdTrack.Py()};\r
1322   Double_t Txy[] = {esdTrack.Xv(), esdTrack.Yv()};\r
1323   Double_t Vxy[] = {vv->GetX(),  vv->GetY()};\r
1324   GetImpactParamVect(Pxy, Txy, Vxy, ipVec2D);\r
1325   Double_t phiIP = TMath::ATan2(ipVec2D[1], ipVec2D[0]) + (fabs(ipVec2D[1])-ipVec2D[1])/fabs(ipVec2D[1])*TMath::Pi();\r
1326   Double_t cosTheta = TMath::Cos(jetPhi - phiIP);\r
1327   Double_t sign = cosTheta/TMath::Abs(cosTheta);\r
1328   significance = TMath::Abs(impPar[0])/TMath::Sqrt(ipCov[0])*sign;\r
1329   //ip = fabs(impPar[0]);\r
1330   fhIPSigBtagQA2->Fill(significance);\r
1331   return significance;\r
1332 }\r
1333 \r
1334 //__________________________________________________________________\r
1335 void AliAnaElectron::GetImpactParamVect(Double_t Pxy[2], Double_t Txy[2], Double_t Vxy[2], Double_t IPxy[2])\r
1336 {\r
1337   //px,py: momentum components at the origin of the track; tx, ty:\r
1338   //origin (x,y) of track; vx, vy: coordinates of primary vertex\r
1339   // analytical geometry auxiliary variables\r
1340   Double_t mr = Pxy[1]/Pxy[0]; //angular coeficient of the straight\r
1341                               //line that lies on top of track\r
1342                               //momentum\r
1343   Double_t br = Txy[1] - mr*Txy[0]; //linear coeficient of the straight\r
1344                                    //line that lies on top of track\r
1345                                    //momentum\r
1346   Double_t ms = -1./mr; //angular coeficient of the straight line that\r
1347                        //lies on top of the impact parameter line\r
1348   //  Double_t bs = Vxy[1] - ms*Vxy[0]; //linear coeficient of the straight\r
1349                                    //line that lies on top of the\r
1350                                    //impact parameter line \r
1351   Double_t xIntersection = (mr*Txy[0] - ms*Vxy[0] + Vxy[1] - Txy[1])/(mr - ms);\r
1352   Double_t yIntersection = mr*xIntersection + br;\r
1353   //if(ceil(10000*yIntersection) - ceil(10000*(ms*xIntersection + bs))\r
1354   //!= 0 )cout<<yIntersection<<", "<<ms*xIntersection + bs<<endl;\r
1355   IPxy[0] = xIntersection - Vxy[0];\r
1356   IPxy[1] = yIntersection - Vxy[1];\r
1357   return;\r
1358 }\r
1359 \r
1360 //__________________________________________________________________\r
1361 Bool_t AliAnaElectron::IsItPhotonic(const AliAODPWG4Particle* part) \r
1362 {\r
1363   //This method checks the opening angle and invariant mass of\r
1364   //electron pairs within the AliAODPWG4Particle list to see if \r
1365   //they are likely to be photonic electrons\r
1366 \r
1367   Bool_t itIS = kFALSE;\r
1368 \r
1369   Double_t massE = 0.000511;\r
1370   Double_t bfield = 5.; //kG\r
1371   if(GetReader()->GetDataType() != AliCaloTrackReader::kMC) bfield = GetReader()->GetBField();\r
1372 \r
1373   Int_t pdg1 = part->GetPdg();\r
1374   Int_t trackId = part->GetTrackLabel(0);\r
1375   AliAODTrack* track = (AliAODTrack*)GetAODCTS()->At(trackId);\r
1376   if(!track) {\r
1377     if(GetDebug() > 0) printf("AliAnaElectron::IsItPhotonic - can't get the AOD Track from the particle!  Skipping the photonic check");\r
1378     return kFALSE; //Don't proceed because we can't get the track\r
1379   }\r
1380 \r
1381   AliExternalTrackParam *param1 = new AliExternalTrackParam(track);\r
1382 \r
1383   //Loop on stored AOD electrons and compute the angle differences and Minv\r
1384   for (Int_t k2 =0; k2 < GetOutputAODBranch()->GetEntriesFast() ; k2++) {\r
1385     AliAODPWG4Particle* part2 = (AliAODPWG4Particle*) GetOutputAODBranch()->At(k2);\r
1386     Int_t track2Id = part2->GetTrackLabel(0);\r
1387     if(trackId == track2Id) continue;\r
1388     Int_t pdg2 = part2->GetPdg();\r
1389     if(TMath::Abs(pdg2) != AliCaloPID::kElectron) continue;\r
1390     if(part2->GetDetector() != fCalorimeter) continue;\r
1391 \r
1392     //JLK: Check opp. sign pairs only\r
1393     if(pdg1*pdg2 > 0) continue; //skip same-sign pairs\r
1394 \r
1395     //propagate to common vertex and check opening angle\r
1396     AliAODTrack* track2 = (AliAODTrack*)GetAODCTS()->At(track2Id);\r
1397     if(!track2) {\r
1398       if(GetDebug() >0) printf("AliAnaElectron::IsItPhotonic - problem getting the partner track.  Continuing on to the next one");\r
1399       continue;\r
1400     }\r
1401     AliExternalTrackParam *param2 = new AliExternalTrackParam(track2);\r
1402     Int_t id1 = 0, id2 = 0;\r
1403     AliESDv0 photonVtx(*param1,id1,*param2,id2);\r
1404     Double_t vx,vy,vz;\r
1405     photonVtx.GetXYZ(vx,vy,vz);\r
1406 \r
1407     Double_t p1mom[3];\r
1408     Double_t p2mom[3];\r
1409     param1->PxPyPz(p1mom);\r
1410     param2->PxPyPz(p2mom);\r
1411 \r
1412     TVector3 p1momAtB(p1mom[0],p1mom[1],p1mom[2]);\r
1413     TVector3 p2momAtB(p2mom[0],p2mom[1],p2mom[2]);\r
1414     TVector3 sumMom = p1momAtB+p2momAtB;\r
1415 \r
1416     Double_t ener1 = sqrt(pow(p1momAtB.Mag(),2) + massE*massE);\r
1417     Double_t ener2 = sqrt(pow(p2momAtB.Mag(),2) + massE*massE);\r
1418     Double_t mass = sqrt(pow((ener1+ener2),2) - pow(sumMom.Mag(),2));\r
1419 \r
1420     Double_t dphi = p1momAtB.DeltaPhi(p2momAtB);\r
1421     fh1OpeningAngle->Fill(dphi);\r
1422     fh1MinvPhoton->Fill(mass);\r
1423 \r
1424     if(mass < 0.1) {\r
1425       if(GetDebug() > 0) printf("######PROBABLY A PHOTON\n");\r
1426       itIS = kTRUE;\r
1427     }\r
1428 \r
1429     //clean up\r
1430     delete param2;\r
1431 \r
1432   }\r
1433 \r
1434   delete param1;\r
1435   return itIS;\r
1436 \r
1437 }\r
1438 \r
1439 //__________________________________________________________________\r
1440 Bool_t AliAnaElectron::IsItPhotonic2(const AliAODPWG4Particle* part) \r
1441 {\r
1442   //This method checks to see whether a track that has been flagged as\r
1443   //an electron was determined to match to a V0 candidate with\r
1444   //invariant mass consistent with photon conversion\r
1445 \r
1446   Bool_t itIS = kFALSE;\r
1447   Int_t id = part->GetTrackLabel(0);\r
1448   \r
1449   //---Get V0s---\r
1450   AliAODEvent *aod = (AliAODEvent*) GetReader()->GetInputEvent();\r
1451   int nv0s = aod->GetNumberOfV0s();\r
1452   for (Int_t iV0 = 0; iV0 < nv0s; iV0++) {\r
1453     AliAODv0 *v0 = aod->GetV0(iV0);\r
1454     if (!v0) continue;\r
1455     double radius = v0->RadiusV0();\r
1456     double mass = v0->InvMass2Prongs(0,1,11,11);\r
1457     if(GetDebug() > 0) {\r
1458       printf("## IsItPhotonic2() :: v0: %d, radius: %f \n", iV0 , radius );\r
1459       printf("## IsItPhotonic2() :: neg-id: %d, pos-id: %d \n", v0->GetNegID(), v0->GetPosID() );\r
1460       printf("## IsItPhotonic2() :: Minv(e,e): %f \n", v0->InvMass2Prongs(0,1,11,11) );\r
1461     }\r
1462     if (mass < 0.100) {\r
1463       if ( id == v0->GetNegID() || id == v0->GetPosID()) {\r
1464         itIS=kTRUE;\r
1465         printf("## IsItPhotonic2() :: It's a conversion electron!!! \n" );\r
1466       }\r
1467     } }\r
1468   return itIS;\r
1469 \r
1470 }\r
1471 \r
1472 //__________________________________________________________________\r
1473 Bool_t AliAnaElectron::GetDCA(const AliAODTrack* track,Double_t impPar[2], Double_t cov[3]) \r
1474 {\r
1475   //Use the Event vertex and AOD track information to get\r
1476   //a real impact parameter for the track\r
1477   //Once alice-off gets its act together and fixes the AOD, this\r
1478   //should become obsolete.\r
1479 \r
1480   Double_t bfield = 5.; //kG\r
1481   Double_t maxD = 100000.; //max transverse IP\r
1482   if(GetReader()->GetDataType() != AliCaloTrackReader::kMC) {\r
1483     bfield = GetReader()->GetBField();\r
1484     AliVEvent* ve = (AliVEvent*)GetReader()->GetInputEvent();\r
1485     AliVVertex *vv = (AliVVertex*)ve->GetPrimaryVertex();\r
1486     AliESDtrack esdTrack(track);\r
1487     Bool_t gotit = esdTrack.PropagateToDCA(vv,bfield,maxD,impPar,cov);\r
1488     return gotit;\r
1489   }\r
1490 \r
1491   return kFALSE;\r
1492 \r
1493 }\r
1494 \r
1495 \r
1496 //__________________________________________________________________\r
1497 Bool_t AliAnaElectron::CheckTrack(const AliAODTrack* track, const char* type) \r
1498 {\r
1499   //Check this track to see if it is also tagged as an electron in the\r
1500   //AliAODPWG4Particle list and if it is non-photonic\r
1501 \r
1502   Bool_t pass = kFALSE;\r
1503 \r
1504   Int_t trackId = track->GetID(); //get the index in the reader\r
1505 \r
1506   Int_t naod = GetOutputAODBranch()->GetEntriesFast();\r
1507   if(GetDebug() > 3) printf("AliAnaElectron::CheckTrack() - aod branch entries %d\n", naod);\r
1508   for(Int_t iaod = 0; iaod < naod ; iaod++){\r
1509     AliAODPWG4Particle* ele =  (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));\r
1510     Int_t label = ele->GetTrackLabel(0);\r
1511     if(label != trackId) continue;  //skip to the next one if they don't match\r
1512 \r
1513     if(type=="DVM") { \r
1514       if(ele->CheckBTagBit(ele->GetBtag(),AliAODPWG4Particle::kDVMTag0) ||\r
1515          ele->CheckBTagBit(ele->GetBtag(),AliAODPWG4Particle::kDVMTag1) ||\r
1516          ele->CheckBTagBit(ele->GetBtag(),AliAODPWG4Particle::kDVMTag2))\r
1517         pass = kTRUE;\r
1518 \r
1519     } else if (type=="NPE") {\r
1520 \r
1521       Bool_t photonic = kFALSE;\r
1522       Bool_t photonic1 = kFALSE;\r
1523       photonic1 = IsItPhotonic(ele); //check against primaries\r
1524       Bool_t photonic2 = kFALSE;\r
1525       photonic2 = IsItPhotonic2(ele); //check against V0s\r
1526       if(photonic1 || photonic2) photonic = kTRUE;\r
1527       \r
1528       if(!photonic) pass = kTRUE;\r
1529 \r
1530     } else {\r
1531       return kFALSE;\r
1532     }\r
1533   }\r
1534 \r
1535   return pass;\r
1536 \r
1537 }\r
1538 \r
1539 //__________________________________________________________________\r
1540 void AliAnaElectron::Print(const Option_t * opt) const\r
1541 {\r
1542   //Print some relevant parameters set for the analysis\r
1543   \r
1544   if(! opt)\r
1545     return;\r
1546   \r
1547   printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;\r
1548   AliAnaPartCorrBaseClass::Print(" ");\r
1549 \r
1550   printf("Calorimeter            =     %s\n", fCalorimeter.Data()) ;\r
1551   printf("pOverE range           =     %f - %f\n",fpOverEmin,fpOverEmax);\r
1552   printf("residual cut           =     %f\n",fResidualCut);\r
1553   printf("---DVM Btagging\n");\r
1554   printf("max IP-cut (e,h)       =     %f\n",fImpactCut);\r
1555   printf("min ITS-hits           =     %d\n",fITSCut);\r
1556   printf("max dR (e,h)           =     %f\n",fDrCut);\r
1557   printf("max pairDCA            =     %f\n",fPairDcaCut);\r
1558   printf("max decaylength        =     %f\n",fDecayLenCut);\r
1559   printf("min Associated Pt      =     %f\n",fAssocPtCut);\r
1560   printf("---IPSig Btagging\n");\r
1561   printf("min tag track          =     %d\n",fNTagTrkCut);\r
1562   printf("min IP significance    =     %f\n",fIPSigCut);\r
1563   printf("    \n") ;\r
1564         \r
1565 } \r
1566 \r
1567 //________________________________________________________________________\r
1568 void AliAnaElectron::ReadHistograms(TList* outputList)\r
1569 {\r
1570   // Needed when Terminate is executed in distributed environment                             \r
1571   // Refill analysis histograms of this class with corresponding\r
1572   // histograms in output list.   \r
1573 \r
1574   // Histograms of this analsys are kept in the same list as other\r
1575   // analysis, recover the position of\r
1576   // the first one and then add the next                                                      \r
1577   Int_t index = outputList->IndexOf(outputList->FindObject(GetAddedHistogramsStringToName()+"fh1pOverE"));\r
1578 \r
1579   //Read histograms, must be in the same order as in\r
1580   //GetCreateOutputObject.                   \r
1581   fh1pOverE     = (TH1F *) outputList->At(index);\r
1582   fh1dR         = (TH1F *) outputList->At(index++);\r
1583   fh2EledEdx    = (TH2F *) outputList->At(index++);\r
1584   fh2MatchdEdx  = (TH2F *) outputList->At(index++);\r
1585   \r
1586 }\r
1587 \r
1588 //__________________________________________________________________\r
1589 void  AliAnaElectron::Terminate(TList* outputList)\r
1590 {\r
1591 \r
1592   //Do some plots to end\r
1593   //Recover histograms from output histograms list, needed for\r
1594   //distributed analysis.                \r
1595   //ReadHistograms(outputList);\r
1596 \r
1597   printf(" AliAnaElectron::Terminate()  *** %s Report: %d outputs\n", GetName(), outputList->GetEntries()) ;\r
1598 \r
1599 }\r
1600 \r