PWGHF/hfe/AliAnalysisTaskEMCalHFEpA.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16         ////////////////////////////////////////////////////////////////////////
  17         //                                                                    //
  18         //      Task for Heavy-flavour electron analysis in pPb collisions    //
  19         //      (+ Electron-Hadron Jetlike Azimuthal Correlation)             //
  20         //                                                                                                                                        //
  21         //              version: September 12th, 2013.                                                            //
  22         //                                                                    //
  23         //          Authors                                                                               //
  24         //              Elienos Pereira de Oliveira Filho (epereira@cern.ch)          //
  25         //              Cristiane Jahnke                (cristiane.jahnke@cern.ch)                    //
  26         //                                                                    //
  27         ////////////////////////////////////////////////////////////////////////
  28
  29 #include "TChain.h"
  30 #include "TTree.h"
  31 #include "TNtuple.h"
  32 #include "TH1F.h"
  33 #include "TH2F.h"
  34 #include "TCanvas.h"
  35 #include "AliAnalysisTask.h"
  36 #include "AliAnalysisManager.h"
  37 #include "AliESDEvent.h"
  38 #include "AliAODEvent.h"
  39 #include "AliVEvent.h"
  40 #include "AliESDInputHandler.h"
  41 #include "AliESDtrackCuts.h"
  42 #include "AliESDCaloCluster.h"
  43 #include "AliESDCaloCells.h"
  44 #include "AliEMCALTrack.h"
  45 #include "AliExternalTrackParam.h"
  46 #include "AliPhysicsSelection.h"
  47 #include "TGeoGlobalMagField.h"
  48 #include "AliMagF.h"
  49 #include "AliLog.h"
  50 #include "AliStack.h"
  51 #include "AliCentrality.h"
  52 #include "AliAODMCParticle.h"
  53 #include "AliAODMCHeader.h"
  54 #include "AliPID.h"
  55 #include "AliPIDResponse.h"
  56 #include "AliHFEcontainer.h"
  57 #include "AliHFEcuts.h"
  58 #include "AliHFEpid.h"
  59 #include "AliHFEpidBase.h"
  60 #include "AliHFEpidQAmanager.h"
  61 #include "AliHFEtools.h"
  62 #include "AliCFContainer.h"
  63 #include "AliCFManager.h"
  64 #include "AliSelectNonHFE.h"
  65 #include "AliHFEpidTPC.h"
  66 #include "AliAnalysisTaskEMCalHFEpA.h"
  67 #include "TMath.h"
  68 #include "THnSparse.h"
  69 #include "TLorentzVector.h"
  70 #include "TString.h"
  71 #include "TFile.h"
  72 #include "AliESDHandler.h"
  73 #include "AliMCEventHandler.h"
  74 #include "AliMCEvent.h"
  75 #include "AliStack.h"
  76 #include "TParticle.h"
  77 #include "AliLog.h"
  78 #include "AliAnalysisTaskSE.h"
  79 #include "TRefArray.h"
  80 #include "TVector.h"
  81 #include "stdio.h"
  82 #include "TGeoManager.h"
  83 #include "iostream"
  84 #include "fstream"
  85 #include "AliKFParticle.h"
  86 #include "AliKFVertex.h"
  87 #include "AliVParticle.h"
  88 #include "AliVTrack.h"
  89 #include "AliEventPoolManager.h"
  90 #include "TObjArray.h"
  91         //______________________________________________________________________
  92
  93         //______________________________________________________________________
  94 ClassImp(AliAnalysisTaskEMCalHFEpA)
  95
  96         //______________________________________________________________________
  97 AliAnalysisTaskEMCalHFEpA::AliAnalysisTaskEMCalHFEpA(const char *name)
  98 : AliAnalysisTaskSE(name)
  99 ,fCorrelationFlag(0)
 100 ,fIsMC(0)
 101 ,fUseEMCal(kFALSE)
 102 ,fUseShowerShapeCut(kFALSE)
 103 ,fFillBackground(kFALSE)
 104
 105
 106 ,fEMCEG1(kFALSE)
 107 ,fEMCEG2(kFALSE)
 108 ,fIsHFE1(kFALSE)
 109 ,fIsHFE2(kFALSE)
 110 ,fIsNonHFE(kFALSE)
 111 ,fIsFromD(kFALSE)
 112 ,fIsFromB(kFALSE)
 113 ,fIsFromPi0(kFALSE)
 114 ,fIsFromEta(kFALSE)
 115 ,fIsFromGamma(kFALSE)
 116 ,fESD(0)
 117 ,fAOD(0)
 118 ,fVevent(0)
 119 ,fPartnerCuts(new AliESDtrackCuts())
 120 ,fOutputList(0)
 121 ,fPidResponse(0)
 122 ,fNonHFE(new AliSelectNonHFE())
 123 ,fIsAOD(kFALSE)
 124 ,fCentrality(0)
 125 ,fCentralityMin(0)
 126 ,fCentralityMax(100)
 127 ,fHasCentralitySelection(kFALSE)
 128 ,fCentralityHist(0)
 129 ,fCentralityHistPass(0)
 130 ,fZvtx(0)
 131 ,fEstimator(0)
 132 ,fClus(0)
 133         //,fClusESD(0)
 134 ,fNevent(0)
 135 ,fPtElec_Inc(0)
 136
 137 ,fCharge_n(0)
 138 ,fCharge_p(0)
 139
 140 ,fTime(0)
 141 ,fTime2(0)
 142 ,ftimingEle(0)
 143 ,ftimingEle2(0)
 144
 145
 146 ,fPtElec_ULS(0)
 147 ,fPtElec_LS(0)
 148 ,fPtElec_ULS2(0)
 149 ,fPtElec_LS2(0)
 150 ,fpid(0)
 151 ,fEoverP_pt(0)
 152 ,fEoverP_tpc(0)
 153 ,fTPC_pt(0)
 154 ,fTPC_p(0)
 155 ,fTPCnsigma_pt(0)
 156 ,fTPCnsigma_p(0)
 157 ,fTPCnsigma_pt_2D(0)
 158 ,fShowerShapeCut(0)
 159
 160 ,fShowerShapeM02_EoverP(0)
 161 ,fShowerShapeM20_EoverP(0)
 162
 163 ,fShowerShape_ha(0)
 164 ,fShowerShape_ele(0)
 165
 166 ,fTPCnsigma_eta(0)
 167 ,fTPCnsigma_phi(0)
 168 ,fECluster(0)
 169 ,fEtaPhi(0)
 170 ,fVtxZ(0)
 171 ,fNTracks(0)
 172 ,fNClusters(0)
 173 ,fTPCNcls_EoverP(0)
 174 ,fEta(0)
 175 ,fPhi(0)
 176 ,fR(0)
 177 ,fR_EoverP(0)
 178 ,fNcells(0)
 179 ,fNcells_EoverP(0)
 180 ,fNcells_electrons(0)
 181 ,fNcells_hadrons(0)
 182 ,fECluster_ptbins(0)
 183 ,fEoverP_ptbins(0)
 184 ,fEoverP_wSSCut(0)
 185 ,fM02_EoverP(0)
 186 ,fM20_EoverP(0)
 187 ,fTPCnsigma_eta_electrons(0)
 188 ,fTPCnsigma_eta_hadrons(0)
 189 ,fEoverP_pt_pions(0)
 190 ,ftpc_p_EoverPcut(0)
 191 ,fnsigma_p_EoverPcut(0)
 192 ,fEoverP_pt_pions2(0)
 193 ,fNcells_pt(0)
 194 ,fEoverP_pt_hadrons(0)
 195 ,fCEtaPhi_Inc(0)
 196 ,fCEtaPhi_ULS(0)
 197 ,fCEtaPhi_LS(0)
 198 ,fCEtaPhi_ULS_NoP(0)
 199 ,fCEtaPhi_LS_NoP(0)
 200 ,fCEtaPhi_ULS_Weight(0)
 201 ,fCEtaPhi_LS_Weight(0)
 202 ,fCEtaPhi_ULS_NoP_Weight(0)
 203 ,fCEtaPhi_LS_NoP_Weight(0)
 204
 205 ,fInvMass(0)
 206 ,fInvMassBack(0)
 207 ,fDCA(0)
 208 ,fDCABack(0)
 209 ,fOpAngle(0)
 210 ,fOpAngleBack(0)
 211
 212 ,fInvMass2(0)
 213 ,fInvMassBack2(0)
 214 ,fDCA2(0)
 215 ,fDCABack2(0)
 216 ,fOpAngle2(0)
 217 ,fOpAngleBack2(0)
 218
 219 ,fMassCut(0.1)
 220 ,fEtaCutMin(-0.9)
 221 ,fEtaCutMax(0.9)
 222
 223 ,fdPhiCut(0.05)
 224 ,fdEtaCut(0.05)
 225
 226
 227 ,fEoverPCutMin(0.8)
 228 ,fEoverPCutMax(1.2)
 229
 230 ,fM20CutMin(0.0)
 231 ,fM20CutMax(10)
 232 ,fM02CutMin(0.0)
 233 ,fM02CutMax(10)
 234
 235 ,fAngleCut(999)
 236 ,fChi2Cut(3.5)
 237 ,fDCAcut(999)
 238 ,fMassCutFlag(kTRUE)
 239 ,fAngleCutFlag(kFALSE)
 240 ,fChi2CutFlag(kFALSE)
 241 ,fDCAcutFlag(kFALSE)
 242 ,fPtBackgroundBeforeReco(0)
 243 ,fPtBackgroundBeforeReco2(0)
 244 ,fPtBackgroundAfterReco(0)
 245
 246 ,fPtMinAsso(0.3)
 247 ,fTpcNclsAsso(80)
 248
 249 ,fPtMCparticleAll(0)
 250 ,fPtMCparticleAll_nonPrimary(0)
 251 ,fPtMCparticleAlle_nonPrimary(0)
 252 ,fPtMCparticleReco(0)
 253 ,fPtMCparticleReco_nonPrimary(0)
 254 ,fPtMCparticleAllHfe1(0)
 255 ,fPtMCparticleRecoHfe1(0)
 256 ,fPtMCparticleAllHfe2(0)
 257 ,fPtMCparticleRecoHfe2(0)
 258 ,fPtMCelectronAfterAll(0)
 259 ,fPtMCelectronAfterAll_nonPrimary(0)
 260 ,fPtMCpi0(0)
 261 ,fPtMCeta(0)
 262 ,fPtMC_EMCal_All(0)
 263 ,fPtMC_EMCal_Selected(0)
 264 ,fPtMC_TPC_All(0)
 265 ,fPtMC_TPC_Selected(0)
 266 ,fPtMCWithLabel(0)
 267 ,fPtMCWithoutLabel(0)
 268 ,fPtIsPhysicaPrimary(0)
 269 ,fCuts(0)
 270 ,fCFM(0)
 271 ,fPID(new AliHFEpid("hfePid"))
 272 ,fPIDqa(0)
 273 ,fMCstack(0)
 274 ,fRejectKinkMother(kFALSE)
 275 ,fMCtrack(0)
 276 ,fMCtrackMother(0)
 277 ,fMCtrackGMother(0)
 278 ,fMCtrackGGMother(0)
 279 ,fMCtrackGGGMother(0)
 280 ,fMCarray(0)
 281 ,fMCheader(0)
 282 ,fMCparticle(0)
 283 ,fMCparticleMother(0)
 284 ,fMCparticleGMother(0)
 285 ,fMCparticleGGMother(0)
 286 ,fMCparticleGGGMother(0)
 287 ,fEventHandler(0)
 288 ,fMCevent(0)
 289 ,fPoolMgr(0)
 290 ,fPool(0)
 291 ,fTracksClone(0)
 292 ,fTracks(0)
 293 ,fCEtaPhi_Inc_EM(0)
 294 ,fCEtaPhi_ULS_EM(0)
 295 ,fCEtaPhi_LS_EM(0)
 296 ,fCEtaPhi_ULS_Weight_EM(0)
 297 ,fCEtaPhi_LS_Weight_EM(0)
 298 ,fPoolNevents(0)
 299 ,fEventMixingFlag(0)
 300 ,fCEtaPhi_Inc_DiHadron(0)
 301 ,fPtTrigger_Inc(0)
 302 {
 303                 //Named constructor
 304                 // Define input and output slots here
 305                 // Input slot #0 works with a TChain
 306         DefineInput(0, TChain::Class());
 307                 // Output slot #0 id reserved by the base class for AOD
 308                 // Output slot #1 writes into a TH1 container
 309                 // DefineOutput(1, TH1I::Class());
 310         DefineOutput(1, TList::Class());
 311                 //  DefineOutput(3, TTree::Class());
 312 }
 313
 314         //________________________________________________________________________
 315 AliAnalysisTaskEMCalHFEpA::AliAnalysisTaskEMCalHFEpA()
 316 : AliAnalysisTaskSE("DefaultAnalysis_AliAnalysisTaskEMCalHFEpA")
 317 ,fCorrelationFlag(0)
 318 ,fIsMC(0)
 319 ,fUseEMCal(kFALSE)
 320 ,fUseShowerShapeCut(kFALSE)
 321 ,fFillBackground(kFALSE)
 322
 323 ,fEMCEG1(kFALSE)
 324 ,fEMCEG2(kFALSE)
 325 ,fIsHFE1(kFALSE)
 326 ,fIsHFE2(kFALSE)
 327 ,fIsNonHFE(kFALSE)
 328 ,fIsFromD(kFALSE)
 329 ,fIsFromB(kFALSE)
 330 ,fIsFromPi0(kFALSE)
 331 ,fIsFromEta(kFALSE)
 332 ,fIsFromGamma(kFALSE)
 333 ,fESD(0)
 334 ,fAOD(0)
 335 ,fVevent(0)
 336 ,fPartnerCuts(new AliESDtrackCuts())
 337 ,fOutputList(0)
 338 ,fPidResponse(0)
 339 ,fNonHFE(new AliSelectNonHFE())
 340 ,fIsAOD(kFALSE)
 341 ,fCentrality(0)
 342 ,fCentralityMin(0)
 343 ,fCentralityMax(100)
 344 ,fHasCentralitySelection(kFALSE)
 345 ,fCentralityHist(0)
 346 ,fCentralityHistPass(0)
 347 ,fZvtx(0)
 348 ,fEstimator(0)
 349 ,fClus(0)
 350         //,fClusESD(0)
 351 ,fNevent(0)
 352 ,fPtElec_Inc(0)
 353
 354 ,fCharge_n(0)
 355 ,fCharge_p(0)
 356
 357 ,fTime(0)
 358 ,fTime2(0)
 359 ,ftimingEle(0)
 360 ,ftimingEle2(0)
 361
 362
 363
 364 ,fPtElec_ULS(0)
 365 ,fPtElec_LS(0)
 366 ,fPtElec_ULS2(0)
 367 ,fPtElec_LS2(0)
 368 ,fpid(0)
 369 ,fEoverP_pt(0)
 370 ,fEoverP_tpc(0)
 371 ,fTPC_pt(0)
 372 ,fTPC_p(0)
 373 ,fTPCnsigma_pt(0)
 374 ,fTPCnsigma_p(0)
 375 ,fTPCnsigma_pt_2D(0)
 376 ,fShowerShapeCut(0)
 377 ,fShowerShapeM02_EoverP(0)
 378 ,fShowerShapeM20_EoverP(0)
 379
 380 ,fShowerShape_ha(0)
 381 ,fShowerShape_ele(0)
 382
 383 ,fTPCnsigma_eta(0)
 384 ,fTPCnsigma_phi(0)
 385 ,fECluster(0)
 386 ,fEtaPhi(0)
 387 ,fVtxZ(0)
 388 ,fNTracks(0)
 389 ,fNClusters(0)
 390 ,fTPCNcls_EoverP(0)
 391 ,fEta(0)
 392 ,fPhi(0)
 393 ,fR(0)
 394 ,fR_EoverP(0)
 395 ,fNcells(0)
 396 ,fNcells_EoverP(0)
 397 ,fNcells_electrons(0)
 398 ,fNcells_hadrons(0)
 399 ,fECluster_ptbins(0)
 400 ,fEoverP_ptbins(0)
 401 ,fEoverP_wSSCut(0)
 402 ,fM02_EoverP(0)
 403 ,fM20_EoverP(0)
 404 ,fTPCnsigma_eta_electrons(0)
 405 ,fTPCnsigma_eta_hadrons(0)
 406 ,fEoverP_pt_pions(0)
 407 ,ftpc_p_EoverPcut(0)
 408 ,fnsigma_p_EoverPcut(0)
 409 ,fEoverP_pt_pions2(0)
 410 ,fNcells_pt(0)
 411 ,fEoverP_pt_hadrons(0)
 412 ,fCEtaPhi_Inc(0)
 413 ,fCEtaPhi_ULS(0)
 414 ,fCEtaPhi_LS(0)
 415 ,fCEtaPhi_ULS_NoP(0)
 416 ,fCEtaPhi_LS_NoP(0)
 417 ,fCEtaPhi_ULS_Weight(0)
 418 ,fCEtaPhi_LS_Weight(0)
 419 ,fCEtaPhi_ULS_NoP_Weight(0)
 420 ,fCEtaPhi_LS_NoP_Weight(0)
 421
 422 ,fInvMass(0)
 423 ,fInvMassBack(0)
 424 ,fDCA(0)
 425 ,fDCABack(0)
 426 ,fOpAngle(0)
 427 ,fOpAngleBack(0)
 428
 429 ,fInvMass2(0)
 430 ,fInvMassBack2(0)
 431 ,fDCA2(0)
 432 ,fDCABack2(0)
 433 ,fOpAngle2(0)
 434 ,fOpAngleBack2(0)
 435
 436 ,fMassCut(0.1)
 437 ,fEtaCutMin(-0.9)
 438 ,fEtaCutMax(0.9)
 439
 440 ,fdPhiCut(0.05)
 441 ,fdEtaCut(0.05)
 442
 443 ,fEoverPCutMin(0.8)
 444 ,fEoverPCutMax(1.2)
 445
 446 ,fM20CutMin(0)
 447 ,fM20CutMax(10)
 448 ,fM02CutMin(0)
 449 ,fM02CutMax(10)
 450
 451 ,fAngleCut(999)
 452 ,fChi2Cut(3.5)
 453 ,fDCAcut(999)
 454 ,fMassCutFlag(kTRUE)
 455 ,fAngleCutFlag(kFALSE)
 456 ,fChi2CutFlag(kFALSE)
 457 ,fDCAcutFlag(kFALSE)
 458 ,fPtBackgroundBeforeReco(0)
 459 ,fPtBackgroundBeforeReco2(0)
 460 ,fPtBackgroundAfterReco(0)
 461
 462 ,fPtMinAsso(0.3)
 463 ,fTpcNclsAsso(80)
 464
 465
 466 ,fPtMCparticleAll(0)
 467 ,fPtMCparticleAll_nonPrimary(0)
 468 ,fPtMCparticleAlle_nonPrimary(0)
 469 ,fPtMCparticleReco(0)
 470 ,fPtMCparticleReco_nonPrimary(0)
 471
 472 ,fPtMCparticleAllHfe1(0)
 473 ,fPtMCparticleRecoHfe1(0)
 474 ,fPtMCparticleAllHfe2(0)
 475 ,fPtMCparticleRecoHfe2(0)
 476 ,fPtMCelectronAfterAll(0)
 477 ,fPtMCelectronAfterAll_nonPrimary(0)
 478
 479 ,fPtMCpi0(0)
 480 ,fPtMCeta(0)
 481 ,fPtMC_EMCal_All(0)
 482 ,fPtMC_EMCal_Selected(0)
 483 ,fPtMC_TPC_All(0)
 484 ,fPtMC_TPC_Selected(0)
 485 ,fPtMCWithLabel(0)
 486 ,fPtMCWithoutLabel(0)
 487 ,fPtIsPhysicaPrimary(0)
 488 ,fCuts(0)
 489 ,fCFM(0)
 490 ,fPID(new AliHFEpid("hfePid"))
 491 ,fPIDqa(0)
 492 ,fMCstack(0)
 493 ,fRejectKinkMother(kFALSE)
 494 ,fMCtrack(0)
 495 ,fMCtrackMother(0)
 496 ,fMCtrackGMother(0)
 497 ,fMCtrackGGMother(0)
 498 ,fMCtrackGGGMother(0)
 499 ,fMCarray(0)
 500 ,fMCheader(0)
 501 ,fMCparticle(0)
 502 ,fMCparticleMother(0)
 503 ,fMCparticleGMother(0)
 504 ,fMCparticleGGMother(0)
 505 ,fMCparticleGGGMother(0)
 506 ,fEventHandler(0)
 507 ,fMCevent(0)
 508 ,fPoolMgr(0)
 509 ,fPool(0)
 510 ,fTracksClone(0)
 511 ,fTracks(0)
 512 ,fCEtaPhi_Inc_EM(0)
 513 ,fCEtaPhi_ULS_EM(0)
 514 ,fCEtaPhi_LS_EM(0)
 515 ,fCEtaPhi_ULS_Weight_EM(0)
 516 ,fCEtaPhi_LS_Weight_EM(0)
 517 ,fPoolNevents(0)
 518 ,fEventMixingFlag(0)
 519 ,fCEtaPhi_Inc_DiHadron(0)
 520 ,fPtTrigger_Inc(0)
 521 {
 522                 // Constructor
 523                 // Define input and output slots here
 524                 // Input slot #0 works with a TChain
 525         DefineInput(0, TChain::Class());
 526                 // Output slot #0 id reserved by the base class for AOD
 527                 // Output slot #1 writes into a TH1 container
 528                 // DefineOutput(1, TH1I::Class());
 529         DefineOutput(1, TList::Class());
 530                 //DefineOutput(3, TTree::Class());
 531 }
 532
 533         //______________________________________________________________________
 534 AliAnalysisTaskEMCalHFEpA::~AliAnalysisTaskEMCalHFEpA()
 535 {
 536                 //Destructor
 537         delete fOutputList;
 538         delete fPID;
 539         delete fCFM;
 540         delete fPIDqa;
 541 }
 542
 543         //______________________________________________________________________
 544         //Create Output Objects
 545         //Here we can define the histograms and others output files
 546         //Called once
 547 void AliAnalysisTaskEMCalHFEpA::UserCreateOutputObjects()
 548 {
 549                 //______________________________________________________________________
 550                 //Initialize PID
 551         if(!fPID->GetNumberOfPIDdetectors())
 552     {
 553                 fPID->AddDetector("TPC", 0);
 554     }
 555
 556         fPID->SortDetectors();
 557
 558         fPIDqa = new AliHFEpidQAmanager();
 559         fPIDqa->Initialize(fPID);
 560                 //______________________________________________________________________
 561
 562                 //______________________________________________________________________
 563                 //Initialize correction Framework and Cuts
 564         fCFM = new AliCFManager;
 565         const Int_t kNcutSteps = AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kNcutStepsRecTrack + AliHFEcuts::kNcutStepsDETrack;
 566         fCFM->SetNStepParticle(kNcutSteps);
 567         for(Int_t istep = 0; istep < kNcutSteps; istep++) fCFM->SetParticleCutsList(istep, NULL);
 568
 569         if(!fCuts)
 570         {
 571                 AliWarning("Cuts not available. Default cuts will be used");
 572                 fCuts = new AliHFEcuts;
 573                 fCuts->CreateStandardCuts();
 574         }
 575
 576         fCuts->Initialize(fCFM);
 577                 //______________________________________________________________________
 578
 579                 ///______________________________________________________________________
 580                 ///Output Tlist
 581                 //Create TList
 582         fOutputList = new TList();
 583         fOutputList->SetOwner();
 584
 585                 //PIDqa
 586         fOutputList->Add(fPIDqa->MakeList("PIDQA"));
 587
 588                 //Store the number of events
 589                 //Define the histo
 590         fNevent = new TH1F("fNevent","Number of Events",5,-0.5,4.5);
 591                 //And then, add to the output list
 592         fOutputList->Add(fNevent);
 593
 594         fpid = new TH1F("fpid","PID flag",5,0,5);
 595         fOutputList->Add(fpid);
 596
 597                 //pt Distribution
 598         fPtElec_Inc = new TH1F("fPtElec_Inc","Inclusive Electrons; p_{T} (GeV/c); Count",300,0,30);
 599
 600         fPtElec_ULS = new TH1F("fPtElec_ULS","Inclusive Electrons; p_{T} (GeV/c); Count",300,0,30);
 601         fPtElec_LS = new TH1F("fPtElec_LS","Inclusive Electrons; p_{T} (GeV/c); Count",300,0,30);
 602
 603         if(fFillBackground){
 604                 fPtElec_ULS2 = new TH1F("fPtElec_ULS2","Inclusive Electrons; p_{T} (GeV/c); Count",300,0,30);
 605                 fPtElec_LS2 = new TH1F("fPtElec_LS2","Inclusive Electrons; p_{T} (GeV/c); Count",300,0,30);
 606         }
 607
 608         fPtTrigger_Inc = new TH1F("fPtTrigger_Inc","pT dist for Hadron Contamination; p_{t} (GeV/c); Count",300,0,30);
 609         fTPCnsigma_pt_2D = new TH2F("fTPCnsigma_pt_2D",";pt (GeV/c);TPC Electron N#sigma",1000,0.3,30,1000,-15,10);
 610         fShowerShapeCut = new TH2F("fShowerShapeCut","Shower Shape;M02;M20",500,0,1.8,500,0,1.8);
 611
 612
 613
 614         fCharge_n = new TH1F("fCharge_n","Inclusive Electrons (Negative Charge); p_{t} (GeV/c); Count",200,0,30);
 615         fCharge_p = new TH1F("fCharge_p","Inclusive Positrons (Positive Charge); p_{t} (GeV/c); Count",200,0,30);
 616
 617         if(fUseEMCal){
 618
 619                 if(!fIsAOD){
 620                         fTime = new TH2D("fTime","Cells Cluster Time; p_{T} (GeV/c); Time (s)",300,0,30,1000,1e-8,1e-5);
 621                         fTime2 = new TH2D("fTime2","Cells Cluster Time;  p_{T} (GeV/c); Time (s)",300,0,30,1000,1e-8,1e-5);
 622                 }
 623
 624                 ftimingEle = new TH2D("ftimingEle","Cluster Time;  p_{T} (GeV/c); Time (s)",300,0,30,1000,1e-8,1e-5);
 625                 ftimingEle2 = new TH2D("ftimingEle2","Cluster Time;  p_{T} (GeV/c); Time (s)",300,0,30,1000,1e-8,1e-5);
 626
 627                 fShowerShape_ha = new TH2F("fShowerShape_ha","Shower Shape hadrons;M02;M20",500,0,1.8,500,0,1.8);
 628                 fShowerShape_ele = new TH2F("fShowerShape_ele","Shower Shape electrons;M02;M20",500,0,1.8,500,0,1.8);
 629
 630                 fShowerShapeM02_EoverP = new TH2F("fShowerShapeM02_EoverP","Shower Shape;M02;E/p",500,0,1.8,500,0,1.8);
 631                 fShowerShapeM20_EoverP = new TH2F("fShowerShapeM20_EoverP","Shower Shape;M20;E/p",500,0,1.8,500,0,1.8);
 632
 633         }
 634
 635         fOutputList->Add(fPtElec_Inc);
 636         fOutputList->Add(fPtElec_ULS);
 637         fOutputList->Add(fPtElec_LS);
 638
 639         if(fFillBackground){
 640                 fOutputList->Add(fPtElec_ULS2);
 641                 fOutputList->Add(fPtElec_LS2);
 642         }
 643
 644
 645         fOutputList->Add(fPtTrigger_Inc);
 646         fOutputList->Add(fTPCnsigma_pt_2D);
 647         fOutputList->Add(fShowerShapeCut);
 648
 649         fOutputList->Add(fCharge_n);
 650         fOutputList->Add(fCharge_p);
 651
 652         if(fUseEMCal){
 653
 654                 if(!fIsAOD){
 655                         fOutputList->Add(fTime);
 656                         fOutputList->Add(fTime2);
 657
 658                 }
 659
 660                 fOutputList->Add(ftimingEle);
 661                 fOutputList->Add(ftimingEle2);
 662
 663                 fOutputList->Add(fShowerShape_ha);
 664                 fOutputList->Add(fShowerShape_ele);
 665
 666                 fOutputList->Add(fShowerShapeM02_EoverP);
 667                 fOutputList->Add(fShowerShapeM20_EoverP);
 668
 669
 670
 671         }
 672
 673                 //General Histograms
 674
 675                 //Steps
 676                 //Step 1: Before Track cuts
 677                 //Step 2: Before PID
 678                 //Step 3: After PID
 679
 680         fEoverP_pt = new TH2F *[3];
 681         fTPC_p = new TH2F *[3];
 682         fTPCnsigma_p = new TH2F *[3];
 683
 684         fECluster= new TH1F *[3];
 685         fEtaPhi= new TH2F *[3];
 686         fVtxZ= new  TH1F *[3];
 687         fNTracks= new  TH1F *[3];
 688         fNClusters= new TH1F *[3];
 689         fTPCNcls_EoverP= new TH2F *[3];
 690
 691         for(Int_t i = 0; i < 3; i++)
 692         {
 693                 fEoverP_pt[i] = new TH2F(Form("fEoverP_pt%d",i),";p_{t} (GeV/c);E / p ",1000,0,30,500,0,2);
 694                 fTPC_p[i] = new TH2F(Form("fTPC_p%d",i),";pt (GeV/c);TPC dE/dx (a. u.)",1000,0.3,15,1000,-20,200);
 695                 fTPCnsigma_p[i] = new TH2F(Form("fTPCnsigma_p%d",i),";p (GeV/c);TPC Electron N#sigma",1000,0.3,15,1000,-15,10);
 696
 697
 698                 fECluster[i]= new TH1F(Form("fECluster%d",i), ";ECluster",2000, 0,100);
 699                 fEtaPhi[i]= new TH2F(Form("fEtaPhi%d",i),"#eta x #phi Clusters;#phi;#eta",200,0.,5,50,-1.,1.);
 700                 fVtxZ[i]= new  TH1F(Form("fVtxZ%d",i),"VtxZ",1000, -50,50);
 701                 fNTracks[i]= new  TH1F(Form("fNTracks%d",i),"NTracks",1000, 0,1000);
 702                 fNClusters[i]= new TH1F(Form("fNClusters%d",i),"fNClusters0",200, 0,100);
 703                 fTPCNcls_EoverP[i]= new TH2F(Form("fTPCNcls_EoverP%d",i),"TPCNcls_EoverP",1000,0,200,200,0,2);
 704
 705
 706                 fOutputList->Add(fEoverP_pt[i]);
 707                 fOutputList->Add(fTPC_p[i]);
 708                 fOutputList->Add(fTPCnsigma_p[i]);
 709
 710
 711                 fOutputList->Add(fECluster[i]);
 712                 fOutputList->Add(fEtaPhi[i]);
 713                 fOutputList->Add(fVtxZ[i]);
 714                 fOutputList->Add(fNTracks[i]);
 715                 fOutputList->Add(fNClusters[i]);
 716                 fOutputList->Add(fTPCNcls_EoverP[i]);
 717         }
 718
 719                 //pt bin
 720         Int_t fPtBin[6] = {2,4,6,8,10,15};
 721
 722         fEoverP_tpc = new TH2F *[5];
 723         fTPC_pt = new TH1F *[5];
 724         fTPCnsigma_pt = new TH1F *[5];
 725
 726         fEta=new TH1F *[5];
 727         fPhi=new TH1F *[5];
 728         fR=new TH1F *[5];
 729         fR_EoverP=new TH2F *[5];
 730         fNcells=new TH1F *[5];
 731         fNcells_EoverP=new TH2F *[5];
 732         fM02_EoverP= new TH2F *[5];
 733         fM20_EoverP= new TH2F *[5];
 734         fEoverP_ptbins=new TH1F *[5];
 735         fECluster_ptbins=new TH1F *[5];
 736         fEoverP_wSSCut=new TH1F *[5];
 737         fNcells_electrons=new TH1F *[5];
 738         fNcells_hadrons=new TH1F *[5];
 739         fTPCnsigma_eta_electrons=new TH2F *[5];
 740         fTPCnsigma_eta_hadrons=new TH2F *[5];
 741
 742         if(fCorrelationFlag)
 743         {
 744                 fCEtaPhi_Inc = new TH2F *[5];
 745                 fCEtaPhi_Inc_DiHadron = new TH2F *[5];
 746
 747                 fCEtaPhi_ULS = new TH2F *[5];
 748                 fCEtaPhi_LS = new TH2F *[5];
 749                 fCEtaPhi_ULS_NoP = new TH2F *[5];
 750                 fCEtaPhi_LS_NoP = new TH2F *[5];
 751
 752                 fCEtaPhi_ULS_Weight = new TH2F *[5];
 753                 fCEtaPhi_LS_Weight = new TH2F *[5];
 754                 fCEtaPhi_ULS_NoP_Weight = new TH2F *[5];
 755                 fCEtaPhi_LS_NoP_Weight = new TH2F *[5];
 756
 757                 fCEtaPhi_Inc_EM = new TH2F *[5];
 758
 759                 fCEtaPhi_ULS_EM = new TH2F *[5];
 760                 fCEtaPhi_LS_EM = new TH2F *[5];
 761
 762                 fCEtaPhi_ULS_Weight_EM = new TH2F *[5];
 763                 fCEtaPhi_LS_Weight_EM = new TH2F *[5];
 764
 765                 fInvMass = new TH1F("fInvMass","",200,0,0.3);
 766                 fInvMassBack = new TH1F("fInvMassBack","",200,0,0.3);
 767                 fDCA = new TH1F("fDCA","",200,0,1);
 768                 fDCABack = new TH1F("fDCABack","",200,0,1);
 769                 fOpAngle = new TH1F("fOpAngle","",200,0,0.5);
 770                 fOpAngleBack = new TH1F("fOpAngleBack","",200,0,0.5);
 771
 772                 fOutputList->Add(fInvMass);
 773                 fOutputList->Add(fInvMassBack);
 774                 fOutputList->Add(fDCA);
 775                 fOutputList->Add(fDCABack);
 776                 fOutputList->Add(fOpAngle);
 777                 fOutputList->Add(fOpAngleBack);
 778         }
 779
 780         if(fFillBackground){
 781
 782                 fInvMass = new TH1F("fInvMass","",200,0,0.3);
 783                 fInvMassBack = new TH1F("fInvMassBack","",200,0,0.3);
 784                 fDCA = new TH1F("fDCA","",200,0,1);
 785                 fDCABack = new TH1F("fDCABack","",200,0,1);
 786                 fOpAngle = new TH1F("fOpAngle","",200,0,0.5);
 787                 fOpAngleBack = new TH1F("fOpAngleBack","",200,0,0.5);
 788
 789                 fOutputList->Add(fInvMass);
 790                 fOutputList->Add(fInvMassBack);
 791                 fOutputList->Add(fDCA);
 792                 fOutputList->Add(fDCABack);
 793                 fOutputList->Add(fOpAngle);
 794                 fOutputList->Add(fOpAngleBack);
 795
 796                         //histos for TPC-only
 797                 fInvMass2 = new TH1F("fInvMass2","",200,0,0.3);
 798                 fInvMassBack2 = new TH1F("fInvMassBack2","",200,0,0.3);
 799                 fDCA2 = new TH1F("fDCA2","",200,0,1);
 800                 fDCABack2 = new TH1F("fDCABack2","",200,0,1);
 801                 fOpAngle2 = new TH1F("fOpAngle2","",200,0,0.5);
 802                 fOpAngleBack2 = new TH1F("fOpAngleBack2","",200,0,0.5);
 803
 804                 fOutputList->Add(fInvMass2);
 805                 fOutputList->Add(fInvMassBack2);
 806                 fOutputList->Add(fDCA2);
 807                 fOutputList->Add(fDCABack2);
 808                 fOutputList->Add(fOpAngle2);
 809                 fOutputList->Add(fOpAngleBack2);
 810
 811         }
 812
 813         for(Int_t i = 0; i < 5; i++)
 814         {
 815                 fEoverP_tpc[i] = new TH2F(Form("fEoverP_tpc%d",i),Form("%d < p_{t} < %d GeV/c;TPC Electron N#sigma;E / p ",fPtBin[i],fPtBin[i+1]),1000,-15,15,100,0,2);
 816                 fTPC_pt[i] = new TH1F(Form("fTPC_pt%d",i),Form("%d < p_{t} < %d GeV/c;TPC Electron N#sigma;Count",fPtBin[i],fPtBin[i+1]),200,20,200);
 817                 fTPCnsigma_pt[i] = new TH1F(Form("fTPCnsigma_pt%d",i),Form("%d < p_{t} < %d GeV/c;TPC Electron N#sigma;Count",fPtBin[i],fPtBin[i+1]),200,-15,10);
 818
 819                 fEta[i]=new TH1F(Form("fEta%d",i), Form("%d < p_{t} < %d GeV/c;#eta; counts",fPtBin[i],fPtBin[i+1]),100, -0.1,0.1);
 820                 fPhi[i]=new TH1F(Form("fPhi%d",i),Form("%d < p_{t} < %d GeV/c;#phi; counts )",fPtBin[i],fPtBin[i+1]), 100, -0.1,0.1);
 821                 fR[i]=new TH1F(Form("fR%d",i),Form("%d < p_{t} < %d GeV/c;R;counts )",fPtBin[i],fPtBin[i+1]), 100, -0.1,0.1);
 822                 fR_EoverP[i]=new TH2F(Form("fR_EoverP%d",i),Form("%d < p_{t} < %d GeV/c;R;E / p ",fPtBin[i],fPtBin[i+1]),100, 0,0.1,1000,0,10);
 823                 fNcells[i]=new TH1F(Form("fNcells%d",i), Form("%d < p_{t} < %d GeV/c;ncells;counts ",fPtBin[i],fPtBin[i+1]),100, 0, 30);
 824                 fNcells_electrons[i]=new TH1F(Form("fNcells_electrons%d",i), Form("%d < p_{t} < %d GeV/c;ncells;counts ",fPtBin[i],fPtBin[i+1]),100, 0, 30);
 825                 fNcells_hadrons[i]=new TH1F(Form("fNcells_hadrons%d",i), Form("%d < p_{t} < %d GeV/c;ncells;counts ",fPtBin[i],fPtBin[i+1]),100, 0, 30);
 826                 fNcells_EoverP[i]=new TH2F(Form("fNcells_EoverP%d",i),Form("%d < p_{t} < %d GeV/c; Ncells; E / p ",fPtBin[i],fPtBin[i+1]),1000, 0,20,100,0,30);
 827                 fM02_EoverP[i]= new TH2F(Form("fM02_EoverP%d",i),Form("%d < p_{t} < %d GeV/c; M02; E / p ",fPtBin[i],fPtBin[i+1]),1000,0,100,100,0,2);
 828                 fM20_EoverP[i]= new TH2F(Form("fM20_EoverP%d",i),Form("%d < p_{t} < %d GeV/c; M20; E / p ",fPtBin[i],fPtBin[i+1]),1000,0,100,100,0,2);
 829                 fEoverP_ptbins[i] = new TH1F(Form("fEoverP_ptbins%d",i),Form("%d < p_{t} < %d GeV/c;E / p ",fPtBin[i],fPtBin[i+1]),500,0,2);
 830                 fECluster_ptbins[i]= new TH1F(Form("fECluster_ptbins%d",i), Form("%d < p_{t} < %d GeV/c;ECluster; Counts ",fPtBin[i],fPtBin[i+1]),2000, 0,100);
 831                 fEoverP_wSSCut[i]=new TH1F(Form("fEoverP_wSSCut%d",i),Form("%d < p_{t} < %d GeV/c;E / p ; Counts",fPtBin[i],fPtBin[i+1]),500,0,2);
 832                 fTPCnsigma_eta_electrons[i]=new TH2F(Form("fTPCnsigma_eta_electrons%d",i),Form("%d < p_{t} < %d GeV/c;TPC Electron N#sigma;Eta ",fPtBin[i],fPtBin[i+1]),1000,-15,15,100,-1,1);
 833                 fTPCnsigma_eta_hadrons[i]=new TH2F(Form("fTPCnsigma_eta_hadrons%d",i),Form("%d < p_{t} < %d GeV/c;TPC Electron N#sigma;Eta ",fPtBin[i],fPtBin[i+1]),1000,-15,15,100,-1,1);
 834
 835                 fOutputList->Add(fEoverP_tpc[i]);
 836                 fOutputList->Add(fTPC_pt[i]);
 837                 fOutputList->Add(fTPCnsigma_pt[i]);
 838
 839                 fOutputList->Add(fEta[i]);
 840                 fOutputList->Add(fPhi[i]);
 841                 fOutputList->Add(fR[i]);
 842                 fOutputList->Add(fR_EoverP[i]);
 843                 fOutputList->Add(fNcells[i]);
 844                 fOutputList->Add(fNcells_electrons[i]);
 845                 fOutputList->Add(fNcells_hadrons[i]);
 846                 fOutputList->Add(fNcells_EoverP[i]);
 847                 fOutputList->Add(fECluster_ptbins[i]);
 848                 fOutputList->Add(fEoverP_ptbins[i]);
 849                 fOutputList->Add(fEoverP_wSSCut[i]);
 850                 fOutputList->Add(fM02_EoverP[i]);
 851                 fOutputList->Add(fM20_EoverP[i]);
 852                 fOutputList->Add(fTPCnsigma_eta_electrons[i]);
 853                 fOutputList->Add(fTPCnsigma_eta_hadrons[i]);
 854
 855
 856                 if(fCorrelationFlag)
 857                 {
 858                         fCEtaPhi_Inc[i] = new TH2F(Form("fCEtaPhi_Inc%d",i),Form("%d < p_{t} < %d GeV/c;#DeltaPhi (rad);#Delta#eta",fPtBin[i],fPtBin[i+1]),200,-0.5*TMath::Pi(),1.5*TMath::Pi(),200,-2,2);
 859                         fCEtaPhi_Inc_DiHadron[i] = new TH2F(Form("fCEtaPhi_Inc_DiHadron%d",i),Form("%d < p_{t} < %d GeV/c;#DeltaPhi (rad);#Delta#eta",fPtBin[i],fPtBin[i+1]),200,-0.5*TMath::Pi(),1.5*TMath::Pi(),200,-2,2);
 860
 861                         fCEtaPhi_ULS[i] = new TH2F(Form("fCEtaPhi_ULS%d",i),Form("%d < p_{t} < %d GeV/c;#DeltaPhi (rad);#Delta#eta",fPtBin[i],fPtBin[i+1]),200,-0.5*TMath::Pi(),1.5*TMath::Pi(),200,-2,2);
 862                         fCEtaPhi_LS[i] = new TH2F(Form("fCEtaPhi_LS%d",i),Form("%d < p_{t} < %d GeV/c;#DeltaPhi (rad);#Delta#eta",fPtBin[i],fPtBin[i+1]),200,-0.5*TMath::Pi(),1.5*TMath::Pi(),200,-2,2);
 863                         fCEtaPhi_ULS_NoP[i] = new TH2F(Form("fCEtaPhi_ULS_NoP%d",i),Form("%d < p_{t} < %d GeV/c;#DeltaPhi (rad);#Delta#eta",fPtBin[i],fPtBin[i+1]),200,-0.5*TMath::Pi(),1.5*TMath::Pi(),200,-2,2);
 864                         fCEtaPhi_LS_NoP[i] = new TH2F(Form("fCEtaPhi_LS_NoP%d",i),Form("%d < p_{t} < %d GeV/c;#DeltaPhi (rad);#Delta#eta",fPtBin[i],fPtBin[i+1]),200,-0.5*TMath::Pi(),1.5*TMath::Pi(),200,-2,2);
 865
 866                         fCEtaPhi_ULS_Weight[i] = new TH2F(Form("fCEtaPhi_ULS_Weight%d",i),Form("%d < p_{t} < %d GeV/c;#DeltaPhi (rad);#Delta#eta",fPtBin[i],fPtBin[i+1]),200,-0.5*TMath::Pi(),1.5*TMath::Pi(),200,-2,2);
 867                         fCEtaPhi_LS_Weight[i] = new TH2F(Form("fCEtaPhi_LS_Weight%d",i),Form("%d < p_{t} < %d GeV/c;#DeltaPhi (rad);#Delta#eta",fPtBin[i],fPtBin[i+1]),200,-0.5*TMath::Pi(),1.5*TMath::Pi(),200,-2,2);
 868                         fCEtaPhi_ULS_NoP_Weight[i] = new TH2F(Form("fCEtaPhi_ULS_NoP_Weight%d",i),Form("%d < p_{t} < %d GeV/c;#DeltaPhi (rad);#Delta#eta",fPtBin[i],fPtBin[i+1]),200,-0.5*TMath::Pi(),1.5*TMath::Pi(),200,-2,2);
 869                         fCEtaPhi_LS_NoP_Weight[i] = new TH2F(Form("fCEtaPhi_LS_NoP_Weight%d",i),Form("%d < p_{t} < %d GeV/c;#DeltaPhi (rad);#Delta#eta",fPtBin[i],fPtBin[i+1]),200,-0.5*TMath::Pi(),1.5*TMath::Pi(),200,-2,2);
 870
 871                         fOutputList->Add(fCEtaPhi_Inc[i]);
 872                         fOutputList->Add(fCEtaPhi_Inc_DiHadron[i]);
 873
 874                         fOutputList->Add(fCEtaPhi_ULS[i]);
 875                         fOutputList->Add(fCEtaPhi_LS[i]);
 876                         fOutputList->Add(fCEtaPhi_ULS_NoP[i]);
 877                         fOutputList->Add(fCEtaPhi_LS_NoP[i]);
 878
 879                         fOutputList->Add(fCEtaPhi_ULS_Weight[i]);
 880                         fOutputList->Add(fCEtaPhi_LS_Weight[i]);
 881                         fOutputList->Add(fCEtaPhi_ULS_NoP_Weight[i]);
 882                         fOutputList->Add(fCEtaPhi_LS_NoP_Weight[i]);
 883
 884                         if(fEventMixingFlag)
 885                         {
 886                                 fCEtaPhi_Inc_EM[i] = new TH2F(Form("fCEtaPhi_Inc_EM%d",i),Form("%d < p_{t} < %d GeV/c;#DeltaPhi (rad);#Delta#eta",fPtBin[i],fPtBin[i+1]),200,-0.5*TMath::Pi(),1.5*TMath::Pi(),200,-2,2);
 887
 888                                 fCEtaPhi_ULS_EM[i] = new TH2F(Form("fCEtaPhi_ULS_EM%d",i),Form("%d < p_{t} < %d GeV/c;#DeltaPhi (rad);#Delta#eta",fPtBin[i],fPtBin[i+1]),200,-0.5*TMath::Pi(),1.5*TMath::Pi(),200,-2,2);
 889                                 fCEtaPhi_LS_EM[i] = new TH2F(Form("fCEtaPhi_LS_EM%d",i),Form("%d < p_{t} < %d GeV/c;#DeltaPhi (rad);#Delta#eta",fPtBin[i],fPtBin[i+1]),200,-0.5*TMath::Pi(),1.5*TMath::Pi(),200,-2,2);
 890
 891                                 fCEtaPhi_ULS_Weight_EM[i] = new TH2F(Form("fCEtaPhi_ULS_Weight_EM%d",i),Form("%d < p_{t} < %d GeV/c;#DeltaPhi (rad);#Delta#eta",fPtBin[i],fPtBin[i+1]),200,-0.5*TMath::Pi(),1.5*TMath::Pi(),200,-2,2);
 892                                 fCEtaPhi_LS_Weight_EM[i] = new TH2F(Form("fCEtaPhi_LS_Weight_EM%d",i),Form("%d < p_{t} < %d GeV/c;#DeltaPhi (rad);#Delta#eta",fPtBin[i],fPtBin[i+1]),200,-0.5*TMath::Pi(),1.5*TMath::Pi(),200,-2,2);
 893
 894                                 fOutputList->Add(fCEtaPhi_Inc_EM[i]);
 895
 896                                 fOutputList->Add(fCEtaPhi_ULS_EM[i]);
 897                                 fOutputList->Add(fCEtaPhi_LS_EM[i]);
 898
 899                                 fOutputList->Add(fCEtaPhi_ULS_Weight_EM[i]);
 900                                 fOutputList->Add(fCEtaPhi_LS_Weight_EM[i]);
 901                         }
 902                 }
 903         }
 904
 905                 //pt integrated
 906         fTPCnsigma_eta = new TH2F("fTPCnsigma_eta",";Pseudorapidity #eta; TPC signal - <TPC signal>_{elec} (#sigma)",200,-0.9,0.9,200,-15,15);
 907         fTPCnsigma_phi = new TH2F("fTPCnsigma_phi",";Azimuthal Angle #phi; TPC signal - <TPC signal>_{elec} (#sigma)",200,0,2*TMath::Pi(),200,-15,15);
 908
 909
 910         fNcells_pt=new TH2F("fNcells_pt","fNcells_pt",1000, 0,20,100,0,30);
 911         fEoverP_pt_pions= new TH2F("fEoverP_pt_pions","fEoverP_pt_pions",1000,0,30,500,0,2);
 912
 913         ftpc_p_EoverPcut= new TH2F("ftpc_p_EoverPcut","ftpc_p_EoverPcut",1000,0,30,200,20,200);
 914         fnsigma_p_EoverPcut= new TH2F("fnsigma_p_EoverPcut","fnsigma_p_EoverPcut",1000,0,30,500,-15,15);
 915
 916         fEoverP_pt_pions2= new TH2F("fEoverP_pt_pions2","fEoverP_pt_pions2",1000,0,30,500,0,2);
 917         fEoverP_pt_hadrons= new TH2F("fEoverP_pt_hadrons","fEoverP_pt_hadrons",1000,0,30,500,0,2);
 918
 919
 920         fOutputList->Add(fTPCnsigma_eta);
 921         fOutputList->Add(fTPCnsigma_phi);
 922
 923         fOutputList->Add(fNcells_pt);
 924         fOutputList->Add(fEoverP_pt_pions);
 925
 926         fOutputList->Add(ftpc_p_EoverPcut);
 927         fOutputList->Add(fnsigma_p_EoverPcut);
 928
 929         fOutputList->Add(fEoverP_pt_pions2);
 930         fOutputList->Add(fEoverP_pt_hadrons);
 931
 932                 //__________________________________________________________________
 933                 //Efficiency studies
 934         if(fIsMC)
 935         {
 936                 fPtBackgroundBeforeReco = new TH1F("fPtBackgroundBeforeReco",";p_{T} (GeV/c);Count",300,0,30);
 937                 if(fFillBackground)fPtBackgroundBeforeReco2 = new TH1F("fPtBackgroundBeforeReco2",";p_{T} (GeV/c);Count",300,0,30);
 938                 fPtBackgroundAfterReco = new TH1F("fPtBackgroundAfterReco",";p_{T} (GeV/c);Count",300,0,30);
 939                 fPtMCparticleAll = new TH1F("fPtMCparticleAll",";p_{T} (GeV/c);Count",200,0,40);
 940                 fPtMCparticleReco = new TH1F("fPtMCparticleReco",";p_{T} (GeV/c);Count",200,0,40);
 941
 942                 fPtMCparticleAll_nonPrimary = new TH1F("fPtMCparticleAll_nonPrimary",";p_{T} (GeV/c);Count",200,0,40);
 943                 fPtMCparticleAlle_nonPrimary = new TH1F("fPtMCparticleAlle_nonPrimary",";p_{T} (GeV/c);Count",200,0,40);
 944                 fPtMCparticleReco_nonPrimary = new TH1F("fPtMCparticleReco_nonPrimary",";p_{T} (GeV/c);Count",200,0,40);
 945
 946                 fPtMCparticleAllHfe1 = new TH1F("fPtMCparticleAllHfe1",";p_{t} (GeV/c);Count",200,0,40);
 947                 fPtMCparticleRecoHfe1 = new TH1F("fPtMCparticleRecoHfe1",";p_{t} (GeV/c);Count",200,0,40);
 948                 fPtMCparticleAllHfe2 = new TH1F("fPtMCparticleAllHfe2",";p_{t} (GeV/c);Count",200,0,40);
 949                 fPtMCparticleRecoHfe2 = new TH1F("fPtMCparticleRecoHfe2",";p_{t} (GeV/c);Count",200,0,40);
 950
 951                 fPtMCelectronAfterAll = new TH1F("fPtMCelectronAfterAll",";p_{T} (GeV/c);Count",200,0,40);
 952                 fPtMCelectronAfterAll_nonPrimary = new TH1F("fPtMCelectronAfterAll_nonPrimary",";p_{T} (GeV/c);Count",200,0,40);
 953
 954
 955
 956                 fPtMCpi0 = new TH1F("fPtMCpi0",";p_{t} (GeV/c);Count",200,0,30);
 957                 fPtMCeta = new TH1F("fPtMCeta",";p_{T} (GeV/c);Count",200,0,30);
 958                 fPtMC_EMCal_All= new TH1F("fPtMC_EMCal_All",";p_{t} (GeV/c);Count",200,0,40);
 959                 fPtMC_EMCal_Selected= new TH1F("fPtMC_EMCal_Selected",";p_{t} (GeV/c);Count",200,0,40);
 960                 fPtMC_TPC_All= new TH1F("fPtMC_TPC_All",";p_{t} (GeV/c);Count",200,0,40);
 961                 fPtMC_TPC_Selected = new TH1F("fPtMC_TPC_Selected",";p_{t} (GeV/c);Count",200,0,40);
 962                 fPtMCWithLabel = new TH1F("fPtMCWithLabel",";p_{t} (GeV/c);Count",200,0,40);
 963                 fPtMCWithoutLabel = new TH1F("fPtMCWithoutLabel",";p_{t} (GeV/c);Count",200,0,40);
 964                 fPtIsPhysicaPrimary = new TH1F("fPtIsPhysicaPrimary",";p_{t} (GeV/c);Count",200,0,40);
 965
 966                 fOutputList->Add(fPtBackgroundBeforeReco);
 967                 if(fFillBackground) fOutputList->Add(fPtBackgroundBeforeReco2);
 968                 fOutputList->Add(fPtBackgroundAfterReco);
 969                 fOutputList->Add(fPtMCparticleAll);
 970                 fOutputList->Add(fPtMCparticleReco);
 971
 972                 fOutputList->Add(fPtMCparticleAll_nonPrimary);
 973                 fOutputList->Add(fPtMCparticleAlle_nonPrimary);
 974                 fOutputList->Add(fPtMCparticleReco_nonPrimary);
 975
 976                 fOutputList->Add(fPtMCparticleAllHfe1);
 977                 fOutputList->Add(fPtMCparticleRecoHfe1);
 978                 fOutputList->Add(fPtMCparticleAllHfe2);
 979                 fOutputList->Add(fPtMCparticleRecoHfe2);
 980                 fOutputList->Add(fPtMCelectronAfterAll);
 981
 982                 fOutputList->Add(fPtMCelectronAfterAll_nonPrimary);
 983
 984
 985                 fOutputList->Add(fPtMCpi0);
 986                 fOutputList->Add(fPtMCeta);
 987                 fOutputList->Add(fPtMC_EMCal_All);
 988                 fOutputList->Add(fPtMC_EMCal_Selected);
 989                 fOutputList->Add(fPtMC_TPC_All);
 990                 fOutputList->Add(fPtMC_TPC_Selected);
 991                 fOutputList->Add(fPtMCWithLabel);
 992                 fOutputList->Add(fPtMCWithoutLabel);
 993                 fOutputList->Add(fPtIsPhysicaPrimary);
 994         }
 995
 996         fCentralityHist = new TH1F("fCentralityHist",";Centrality (%); Count",1000000,0,100);
 997         fCentralityHistPass = new TH1F("fCentralityHistPass",";Centrality (%); Count",1000000,0,100);
 998         fOutputList->Add(fCentralityHist);
 999         fOutputList->Add(fCentralityHistPass);
1000
1001                 //______________________________________________________________________
1002                 //Mixed event analysis
1003         if(fEventMixingFlag)
1004         {
1005                 fPoolNevents = new TH1F("fPoolNevents","Event Mixing Statistics; Number of events; Count",1000,0,1000);
1006                 fOutputList->Add(fPoolNevents);
1007
1008                 Int_t trackDepth = 2000; // number of objects (tracks) kept per event buffer bin. Once the number of stored objects (tracks) is above that limit, the oldest ones are removed.
1009                 Int_t poolsize   = 1000;  // Maximum number of events, ignored in the present implemented of AliEventPoolManager
1010
1011                 Int_t nCentralityBins  = 15;
1012                 Double_t centralityBins[] = { 0, 1, 2, 3, 4, 5, 10, 20, 30, 40, 50, 60, 70, 80, 90, 100.1 };
1013
1014                 Int_t nZvtxBins  = 9;
1015                 Double_t vertexBins[] = {-10, -7, -5, -3, -1, 1, 3, 5, 7, 10};
1016
1017                 fPoolMgr = new AliEventPoolManager(poolsize, trackDepth, nCentralityBins, (Double_t*) centralityBins, nZvtxBins, (Double_t*) vertexBins);
1018         }
1019                 //______________________________________________________________________
1020
1021         PostData(1, fOutputList);
1022
1023                 ///______________________________________________________________________
1024 }
1025
1026         //______________________________________________________________________
1027         //Main loop
1028         //Called for each event
1029 void AliAnalysisTaskEMCalHFEpA::UserExec(Option_t *)
1030 {
1031                 //Check Event
1032         fESD = dynamic_cast<AliESDEvent*>(InputEvent());
1033         fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
1034
1035         if(!(fESD || fAOD))
1036         {
1037                 printf("ERROR: fESD & fAOD not available\n");
1038                 return;
1039         }
1040
1041         fVevent = dynamic_cast<AliVEvent*>(InputEvent());
1042
1043         if(!fVevent)
1044         {
1045                 printf("ERROR: fVEvent not available\n");
1046                 return;
1047         }
1048
1049                 //Check Cuts
1050         if(!fCuts)
1051         {
1052                 AliError("HFE cuts not available");
1053                 return;
1054         }
1055                 //Check PID
1056         if(!fPID->IsInitialized())
1057         {
1058                         // Initialize PID with the given run number
1059                 AliWarning("PID not initialised, get from Run no");
1060
1061                 if(fIsAOD)
1062                 {
1063                         fPID->InitializePID(fAOD->GetRunNumber());
1064                 }
1065                 else
1066                 {
1067                         fPID->InitializePID(fESD->GetRunNumber());
1068                 }
1069         }
1070
1071                 //PID response
1072         fPidResponse = fInputHandler->GetPIDResponse();
1073
1074
1075                 //Check PID response
1076         if(!fPidResponse)
1077         {
1078                 AliDebug(1, "Using default PID Response");
1079                 fPidResponse = AliHFEtools::GetDefaultPID(kFALSE, fInputEvent->IsA() == AliAODEvent::Class());
1080         }
1081
1082         fPID->SetPIDResponse(fPidResponse);
1083
1084         fCFM->SetRecEventInfo(fVevent);
1085
1086         Double_t *fListOfmotherkink = 0;
1087         Int_t fNumberOfVertices = 0;
1088         Int_t fNumberOfMotherkink = 0;
1089
1090                 //______________________________________________________________________
1091                 //Vertex Selection
1092         if(fIsAOD)
1093         {
1094                 const AliAODVertex* trkVtx = fAOD->GetPrimaryVertex();
1095                 if(!trkVtx || trkVtx->GetNContributors()<=0) return;
1096                 TString vtxTtl = trkVtx->GetTitle();
1097                 if(!vtxTtl.Contains("VertexerTracks")) return;
1098                 Float_t zvtx = trkVtx->GetZ();
1099                 fZvtx = zvtx;
1100                 const AliAODVertex* spdVtx = fAOD->GetPrimaryVertexSPD();
1101                 if(spdVtx->GetNContributors()<=0) return;
1102                 TString vtxTyp = spdVtx->GetTitle();
1103                 Double_t cov[6]={0};
1104                 spdVtx->GetCovarianceMatrix(cov);
1105                 Double_t zRes = TMath::Sqrt(cov[5]);
1106                 if(vtxTyp.Contains("vertexer:Z") && (zRes>0.25)) return;
1107                 if(TMath::Abs(spdVtx->GetZ() - trkVtx->GetZ())>0.5) return;
1108                 if(TMath::Abs(zvtx) > 10) return;
1109
1110                         //Look for kink mother for AOD
1111
1112                 fNumberOfVertices = 0;
1113                 fNumberOfMotherkink = 0;
1114
1115                 if(fIsAOD)
1116                 {
1117                         fNumberOfVertices = fAOD->GetNumberOfVertices();
1118
1119                         fListOfmotherkink = new Double_t[fNumberOfVertices];
1120
1121                         for(Int_t ivertex=0; ivertex < fNumberOfVertices; ivertex++)
1122                         {
1123                                 AliAODVertex *aodvertex = fAOD->GetVertex(ivertex);
1124                                 if(!aodvertex) continue;
1125                                 if(aodvertex->GetType()==AliAODVertex::kKink)
1126                                 {
1127                                         AliAODTrack *mother1 = (AliAODTrack *) aodvertex->GetParent();
1128                                         if(!mother1) continue;
1129                                         Int_t idmother = mother1->GetID();
1130                                         fListOfmotherkink[fNumberOfMotherkink] = idmother;
1131                                         fNumberOfMotherkink++;
1132                                 }
1133                         }
1134                 }
1135         }
1136         else
1137         {
1138
1139
1140
1141                         /// ESD
1142                 const AliESDVertex* trkVtx = fESD->GetPrimaryVertex();
1143                 if(!trkVtx || trkVtx->GetNContributors()<=0) return;
1144                 TString vtxTtl = trkVtx->GetTitle();
1145                 if(!vtxTtl.Contains("VertexerTracks")) return;
1146                 Float_t zvtx = trkVtx->GetZ();
1147
1148                 const AliESDVertex* spdVtx = fESD->GetPrimaryVertexSPD();
1149                 if(spdVtx->GetNContributors()<=0) return;
1150                 TString vtxTyp = spdVtx->GetTitle();
1151                 Double_t cov[6]={0};
1152                 spdVtx->GetCovarianceMatrix(cov);
1153                 Double_t zRes = TMath::Sqrt(cov[5]);
1154                 if(vtxTyp.Contains("vertexer:Z") && (zRes>0.25)) return;
1155                 if(TMath::Abs(spdVtx->GetZ() - trkVtx->GetZ())>0.5) return;
1156                 if(TMath::Abs(zvtx) > 10) return;
1157         }
1158
1159                 //______________________________________________________________________
1160
1161                 //Only events with at least 2 tracks are accepted
1162         Int_t fNOtrks =  fVevent->GetNumberOfTracks();
1163         if(fNOtrks<2) return;
1164
1165                 //______________________________________________________________________
1166                 //Centrality Selection
1167         if(fHasCentralitySelection)
1168         {
1169                 Float_t centrality = -1;
1170
1171                 if(fIsAOD)
1172                 {
1173                         fCentrality = fAOD->GetHeader()->GetCentralityP();
1174                 }
1175                 else
1176                 {
1177                         fCentrality = fESD->GetCentrality();
1178                 }
1179
1180                 if(fEstimator==1) centrality = fCentrality->GetCentralityPercentile("ZDC");
1181                 else centrality = fCentrality->GetCentralityPercentile("V0A");
1182
1183                         //cout << "Centrality = " << centrality << " %" << endl;
1184
1185                 fCentralityHist->Fill(centrality);
1186
1187                 if(centrality<fCentralityMin || centrality>fCentralityMax) return;
1188
1189                 fCentralityHistPass->Fill(centrality);
1190         }
1191                 //______________________________________________________________________
1192
1193                 //______________________________________________________________________
1194
1195         if(fIsMC)
1196         {
1197                 if(fIsAOD)
1198                 {
1199                         fMCarray = dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
1200
1201                         if(!fMCarray)
1202                         {
1203                                 AliError("Array of MC particles not found");
1204                                 return;
1205                         }
1206
1207                         fMCheader = dynamic_cast<AliAODMCHeader*>(fAOD->GetList()->FindObject(AliAODMCHeader::StdBranchName()));
1208
1209                         if(!fMCheader)
1210                         {
1211                                 AliError("Could not find MC Header in AOD");
1212                                 return;
1213                         }
1214
1215                         for(Int_t iMC = 0; iMC < fMCarray->GetEntries(); iMC++)
1216                         {
1217                                 fMCparticle = (AliAODMCParticle*) fMCarray->At(iMC);
1218
1219                                 Int_t pdg = fMCparticle->GetPdgCode();
1220                                 if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax && fMCparticle->Charge()!=0)
1221                                 {
1222                                         if (TMath::Abs(pdg) == 11) fPtMCparticleAlle_nonPrimary->Fill(fMCparticle->Pt()); //denominator for total efficiency for all electrons, and not primary
1223
1224                                         if( TMath::Abs(pdg) == 211 || TMath::Abs(pdg) == 2212 || TMath::Abs(pdg) == 321 || TMath::Abs(pdg) == 11 || TMath::Abs(pdg) == 13 )
1225                                         {
1226
1227                                                 fPtMCparticleAll_nonPrimary->Fill(fMCparticle->Pt()); //denominator for total efficiency for all particles, and not primary
1228                                                 if(fMCparticle->IsPhysicalPrimary())
1229                                                 {
1230                                                         fPtMCparticleAll->Fill(fMCparticle->Pt());
1231
1232                                                         Bool_t MotherFound = FindMother(iMC);
1233                                                         if(MotherFound)
1234                                                         {
1235                                                                 if(fIsHFE1) fPtMCparticleAllHfe1->Fill(fMCparticle->Pt()); //denominator for total efficiency
1236                                                                 if(fIsHFE2) fPtMCparticleAllHfe2->Fill(fMCparticle->Pt());
1237                                                         }
1238                                                 }
1239                                         }
1240                                 }
1241                                 if(TMath::Abs(pdg)==111) fPtMCpi0->Fill(fMCparticle->Pt());
1242                                 if(TMath::Abs(pdg)==211) fPtMCeta->Fill(fMCparticle->Pt());
1243                         }
1244                 }
1245                 else
1246                 {
1247                         fEventHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
1248                         if (!fEventHandler) {
1249                                 Printf("ERROR: Could not retrieve MC event handler");
1250                                 return;
1251                         }
1252
1253                         fMCevent = fEventHandler->MCEvent();
1254                         if (!fMCevent) {
1255                                 Printf("ERROR: Could not retrieve MC event");
1256                                 return;
1257                         }
1258
1259                         fMCstack = fMCevent->Stack();
1260
1261                 for(Int_t iMC = 0; iMC < fMCstack->GetNtrack(); iMC++)
1262                 {
1263
1264                                 fMCtrack = fMCstack->Particle(iMC);
1265                                 Int_t pdg = fMCtrack->GetPdgCode();
1266
1267                                 if(TMath::Abs(pdg)==111) fPtMCpi0->Fill(fMCtrack->Pt());
1268                                 if(TMath::Abs(pdg)==211) fPtMCeta->Fill(fMCtrack->Pt());
1269
1270
1271
1272
1273                                 if(fMCtrack->Eta()>=fEtaCutMin && fMCtrack->Eta()<=fEtaCutMax)
1274                                 {
1275
1276                                         if (TMath::Abs(pdg) == 11)  fPtMCparticleAlle_nonPrimary->Fill(fMCtrack->Pt());//denominator for total efficiency for all electrons, and not primary
1277
1278
1279                                         if( TMath::Abs(pdg) == 211 || TMath::Abs(pdg) == 2212 || TMath::Abs(pdg) == 321 || TMath::Abs(pdg) == 11 || TMath::Abs(pdg) == 13 )
1280                                         {
1281                                                 fPtMCparticleAll_nonPrimary->Fill(fMCtrack->Pt());//denominator for total efficiency for all particle, non Primary track
1282
1283                                                 if(!fMCstack->IsPhysicalPrimary(iMC)) continue;
1284                                                 fPtMCparticleAll->Fill(fMCtrack->Pt());
1285
1286                                                 Bool_t MotherFound = FindMother(iMC);
1287                                                 if(MotherFound)
1288                                                 {
1289                                                         if(fIsHFE1) fPtMCparticleAllHfe1->Fill(fMCtrack->Pt());//denominator for total efficiency
1290                                                         if(fIsHFE2) fPtMCparticleAllHfe2->Fill(fMCtrack->Pt());
1291                                                 }
1292                                         }
1293                                 }
1294                 }
1295                 }
1296         }
1297
1298                 //______________________________________________________________________
1299                 //EMCal Trigger Selection (Threshould selection)
1300         TString firedTrigger;
1301         TString TriggerEG1("EG1");
1302         TString TriggerEG2("EG2");
1303
1304         if(fAOD) firedTrigger = fAOD->GetFiredTriggerClasses();
1305         else if(fESD) firedTrigger = fESD->GetFiredTriggerClasses();
1306
1307         fNevent->Fill(0);
1308         if(firedTrigger.Contains(TriggerEG1)) fNevent->Fill(1);
1309         if(firedTrigger.Contains(TriggerEG2)) fNevent->Fill(2);
1310
1311                 //EG1
1312         if(firedTrigger.Contains(TriggerEG1))
1313         {
1314                 fNevent->Fill(3);
1315         }
1316         else
1317         {
1318                 if(fEMCEG1) return;
1319         }
1320
1321                 //EG2
1322         if(firedTrigger.Contains(TriggerEG2))
1323         {
1324                 fNevent->Fill(4);
1325         }
1326         else
1327         {
1328                 if(fEMCEG2) return;
1329         }
1330
1331                 //______________________________________________________________________
1332
1333         Int_t ClsNo = -999;
1334         if(!fIsAOD) ClsNo = fESD->GetNumberOfCaloClusters();
1335         else ClsNo = fAOD->GetNumberOfCaloClusters();
1336
1337                 //______________________________________________________________________
1338
1339                 ///______________________________________________________________________
1340                 ///Track loop
1341         for(Int_t iTracks = 0; iTracks < fVevent->GetNumberOfTracks(); iTracks++)
1342         {
1343                 AliVParticle* Vtrack = fVevent->GetTrack(iTracks);
1344                 if (!Vtrack)
1345                 {
1346                         printf("ERROR: Could not receive track %d\n", iTracks);
1347                         continue;
1348                 }
1349
1350                 AliVTrack *track = dynamic_cast<AliVTrack*>(Vtrack);
1351                 AliESDtrack *etrack = dynamic_cast<AliESDtrack*>(Vtrack);
1352                 AliAODTrack *atrack = dynamic_cast<AliAODTrack*>(Vtrack);
1353
1354                 Double_t fTPCnSigma = -999;
1355                 Double_t fTPCnSigma_pion = -999;
1356                 Double_t fTPCnSigma_proton = -999;
1357                 Double_t fTPCnSigma_kaon = -999;
1358                 Double_t fTPCsignal = -999;
1359                 Double_t fPt = -999;
1360                 Double_t fP = -999;
1361
1362                         ///_____________________________________________________________________________
1363                         ///Fill QA plots without track selection
1364                 fPt = track->Pt();
1365                 fP = TMath::Sqrt((track->Pt())*(track->Pt()) + (track->Pz())*(track->Pz()));
1366
1367                 fTPCsignal = track->GetTPCsignal();
1368                 fTPCnSigma = fPidResponse->NumberOfSigmasTPC(track, AliPID::kElectron);
1369                 fTPCnSigma_pion = fPidResponse->NumberOfSigmasTPC(track, AliPID::kPion);
1370                 fTPCnSigma_proton = fPidResponse->NumberOfSigmasTPC(track, AliPID::kProton);
1371                 fTPCnSigma_kaon = fPidResponse->NumberOfSigmasTPC(track, AliPID::kKaon);
1372
1373                 fTPC_p[0]->Fill(fPt,fTPCsignal);
1374                 fTPCnsigma_p[0]->Fill(fP,fTPCnSigma);
1375
1376
1377                 Float_t TPCNcls = track->GetTPCNcls();
1378                 Float_t pos[3]={0,0,0};
1379
1380                 Double_t fEMCflag = kFALSE;
1381                 if(track->GetEMCALcluster()>0)
1382                 {
1383                         fClus = fVevent->GetCaloCluster(track->GetEMCALcluster());
1384                         if(fClus->IsEMCAL())
1385                         {
1386                                 if(TMath::Abs(fClus->GetTrackDx())<=fdPhiCut && TMath::Abs(fClus->GetTrackDz())<=fdEtaCut)
1387                                 {
1388                                         fEMCflag = kTRUE;
1389                                         fEoverP_pt[0]->Fill(fPt,(fClus->E() / fP));
1390
1391
1392                                         Float_t Energy  = fClus->E();
1393                                         Float_t EoverP  = Energy/track->P();
1394                                                 //Float_t M02   = fClus->GetM02();
1395                                                 //Float_t M20   = fClus->GetM20();
1396
1397                                                 /////////////// for Eta Phi distribution
1398                                         fClus->GetPosition(pos);
1399                                         TVector3 vpos(pos[0],pos[1],pos[2]);
1400                                         Double_t cphi = vpos.Phi();
1401                                         Double_t ceta = vpos.Eta();
1402                                         fEtaPhi[0]->Fill(cphi,ceta);
1403
1404                                         fECluster[0]->Fill(Energy);
1405                                         fTPCNcls_EoverP[0]->Fill(TPCNcls, EoverP);
1406                                 }
1407                         }
1408                 }
1409
1410                         //______________________________________________________________
1411                         // Vertex
1412
1413                 fVtxZ[0]->Fill(fZvtx);
1414                 fNTracks[0]->Fill(fNOtrks);
1415                 fNClusters[0]->Fill(ClsNo);
1416
1417                         //______________________________________________________________
1418
1419                         ///Fill QA plots without track selection
1420                         ///_____________________________________________________________________________
1421                         //______________________________________________________________________________________
1422                         //Track Selection Cuts
1423
1424                         //AOD (Test Filter Bit)
1425                 if(fIsAOD)
1426                 {
1427                                 // standard cuts with very loose DCA - BIT(4)
1428                                 // Description:
1429                         /*
1430                          GetStandardITSTPCTrackCuts2011(kFALSE)
1431                          SetMaxChi2PerClusterTPC(4);
1432                          SetAcceptKinkDaughters(kFALSE);
1433                          SetRequireTPCRefit(kTRUE);
1434                          SetRequireITSRefit(kTRUE);
1435                          SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kAny);
1436                          SetMaxDCAToVertexZ(2);
1437                          SetMaxDCAToVertex2D(kFALSE);
1438                          SetRequireSigmaToVertex(kFALSE);
1439                          SetMaxChi2PerClusterITS(36);
1440                          SetMaxDCAToVertexXY(2.4)
1441                          SetMaxDCAToVertexZ(3.2)
1442                          SetDCaToVertex2D(kTRUE)
1443                          */
1444
1445                         if(!atrack->TestFilterMask(AliAODTrack::kTrkGlobalNoDCA)) continue;
1446                 }
1447
1448                         //RecKine: ITSTPC cuts
1449                 if(!ProcessCutStep(AliHFEcuts::kStepRecKineITSTPC, track)) continue;
1450                         //RecKink
1451                 if(fRejectKinkMother)
1452                 {
1453                         if(fIsAOD)
1454                         {
1455                                 Bool_t kinkmotherpass = kTRUE;
1456                                 for(Int_t kinkmother = 0; kinkmother < fNumberOfMotherkink; kinkmother++)
1457                                 {
1458                                         if(track->GetID() == fListOfmotherkink[kinkmother])
1459                                         {
1460                                                 kinkmotherpass = kFALSE;
1461                                                 continue;
1462                                         }
1463                                 }
1464                                 if(!kinkmotherpass) continue;
1465                         }
1466                         else
1467                         {
1468                                 if(etrack->GetKinkIndex(0) != 0) continue;
1469                         }
1470                 }
1471
1472                         //RecPrim
1473                 if(!fIsAOD)
1474                 {
1475                         if(!ProcessCutStep(AliHFEcuts::kStepRecPrim, track)) continue;
1476                 }
1477
1478                         //HFEcuts: ITS layers cuts
1479                 if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsITS, track)) continue;
1480
1481                         //HFE cuts: TPC PID cleanup
1482                 if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTPC, track)) continue;
1483                         //______________________________________________________________________________________
1484
1485                         ///_____________________________________________________________
1486                         ///QA plots after track selection
1487                 if(fIsMC)
1488                 {
1489                         if(track->GetLabel()>=0) fPtMCWithLabel->Fill(fPt);
1490                         else fPtMCWithoutLabel->Fill(fPt);
1491                 }
1492
1493                 if(fIsMC && fIsAOD && track->GetLabel()>=0)
1494                 {
1495                         fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
1496
1497                         if(fMCparticle->IsPhysicalPrimary()) fPtIsPhysicaPrimary->Fill(fPt);
1498
1499                         Int_t pdg = fMCparticle->GetPdgCode();
1500                         if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax && fMCparticle->Charge()!=0)
1501                         {
1502
1503                                 fPtMCparticleReco_nonPrimary->Fill(fMCparticle->Pt()); //not Primary track
1504
1505                                 if( TMath::Abs(pdg) == 211 || TMath::Abs(pdg) == 2212 || TMath::Abs(pdg) == 321 || TMath::Abs(pdg) == 11 || TMath::Abs(pdg) == 13 )
1506                                 {
1507                                         fPtMCparticleReco_nonPrimary->Fill(fMCparticle->Pt()); //not Primary track
1508
1509                                         if(fMCparticle->IsPhysicalPrimary())
1510                                         {
1511                                                 fPtMCparticleReco->Fill(fMCparticle->Pt());
1512
1513                                                 Bool_t MotherFound = FindMother(track->GetLabel());
1514                                                 if(MotherFound)
1515                                                 {
1516                                                         if(fIsHFE1) fPtMCparticleRecoHfe1->Fill(fMCparticle->Pt());
1517                                                         if(fIsHFE2) fPtMCparticleRecoHfe2->Fill(fMCparticle->Pt());
1518                                                 }
1519                                         }
1520                                 }
1521                         }
1522                 }
1523                 else if(fIsMC && track->GetLabel()>=0)
1524                 {
1525                         if(fMCtrack->Eta()>=fEtaCutMin && fMCtrack->Eta()<=fEtaCutMax)
1526                         {
1527
1528                                 fMCtrack = fMCstack->Particle(track->GetLabel());
1529                                 Int_t pdg = fMCtrack->GetPdgCode();
1530
1531                                 if( TMath::Abs(pdg) == 211 || TMath::Abs(pdg) == 2212 || TMath::Abs(pdg) == 321 || TMath::Abs(pdg) == 11 || TMath::Abs(pdg) == 13 )
1532                                 {
1533                                         fPtMCparticleReco_nonPrimary->Fill(fMCtrack->Pt());//not Primary track
1534                                 }
1535
1536
1537                                 if(fMCstack->IsPhysicalPrimary(track->GetLabel()))
1538                                 {
1539                                                 fPtIsPhysicaPrimary->Fill(fPt);
1540
1541
1542
1543
1544                                                 if( TMath::Abs(pdg) == 211 || TMath::Abs(pdg) == 2212 || TMath::Abs(pdg) == 321 || TMath::Abs(pdg) == 11 || TMath::Abs(pdg) == 13 )
1545                                                 {
1546                                                         fPtMCparticleReco->Fill(fMCtrack->Pt());
1547
1548                                                         Bool_t MotherFound = FindMother(track->GetLabel());
1549                                                         if(MotherFound)
1550                                                         {
1551                                                                 if(fIsHFE1) fPtMCparticleRecoHfe1->Fill(fMCtrack->Pt());
1552                                                                 if(fIsHFE2) fPtMCparticleRecoHfe2->Fill(fMCtrack->Pt());
1553                                                         }
1554                                                 }
1555                                 }
1556                         }
1557                 }
1558
1559                 fTPC_p[1]->Fill(fPt,fTPCsignal);
1560                 fTPCnsigma_p[1]->Fill(fP,fTPCnSigma);
1561                 Double_t fPtBin[6] = {2,4,6,8,10,15};
1562
1563                 TPCNcls = track->GetTPCNcls();
1564                 Float_t pos2[3]={0,0,0};
1565
1566                 if(track->GetEMCALcluster()>0)
1567                 {
1568                         fClus = fVevent->GetCaloCluster(track->GetEMCALcluster());
1569                         if(fClus->IsEMCAL())
1570                         {
1571                                 if(TMath::Abs(fClus->GetTrackDx())<=fdPhiCut && TMath::Abs(fClus->GetTrackDz())<=fdEtaCut)
1572                                 {
1573                                         fEoverP_pt[1]->Fill(fPt,(fClus->E() / fP));
1574
1575                                         Float_t Energy  = fClus->E();
1576                                         Float_t EoverP  = Energy/track->P();
1577                                         Float_t M02     = fClus->GetM02();
1578                                         Float_t M20     = fClus->GetM20();
1579                                         Float_t ncells  = fClus->GetNCells();
1580
1581                                                 /////////////// for Eta Phi distribution
1582                                         fClus->GetPosition(pos2);
1583                                         TVector3 vpos(pos2[0],pos2[1],pos2[2]);
1584                                         Double_t cphi = vpos.Phi();
1585                                         Double_t ceta = vpos.Eta();
1586                                         fEtaPhi[1]->Fill(cphi,ceta);
1587
1588                                         fECluster[1]->Fill(Energy);
1589                                         fTPCNcls_EoverP[1]->Fill(TPCNcls, EoverP);
1590
1591
1592                                                 //for EMCal trigger performance
1593                                         if(EoverP > 0.9){
1594                                                 ftpc_p_EoverPcut->Fill(track->P(), fTPCsignal);
1595                                                 fnsigma_p_EoverPcut->Fill(track->P(), fTPCnSigma);
1596
1597                                         }
1598
1599
1600                                                 //for hadron contamination calculations
1601
1602
1603                                                 // EtaCut -> dados
1604                                         if(track->Eta()>=fEtaCutMin && track->Eta()<=fEtaCutMax ){
1605                                                         //main
1606                                                 if(TMath::Abs(fTPCnSigma_pion)<3 || TMath::Abs(fTPCnSigma_proton)<3 || TMath::Abs(fTPCnSigma_kaon)<3 ){
1607
1608                                                         if(fTPCnSigma<-3.5){
1609                                                                 fEoverP_pt_hadrons->Fill(fPt,EoverP);
1610                                                                 if(fUseEMCal) fShowerShape_ha->Fill(M02,M20);
1611                                                         }
1612                                                 }
1613
1614                                                 if(fTPCnSigma < -3.5){
1615                                                         fEoverP_pt_pions->Fill(fPt, EoverP);
1616
1617                                                 }
1618
1619                                                 if(fTPCnSigma < -3.5 && fTPCnSigma > -10){
1620                                                         fEoverP_pt_pions2->Fill(fPt, EoverP);
1621
1622                                                 }
1623
1624
1625                                         }
1626
1627
1628
1629
1630                                         for(Int_t i = 0; i < 5; i++)
1631                                         {
1632                                                 if(fPt>=fPtBin[i] && fPt<fPtBin[i+1])
1633                                                 {
1634
1635                                                         if(fTPCnSigma < -5 && fTPCnSigma > -10){
1636                                                                 fNcells_hadrons[i]->Fill(ncells);
1637                                                         }
1638                                                                 //hadrons selection using E/p
1639                                                         if((fClus->E() / fP) >= 0.2 && (fClus->E() / fP) <= 0.7){
1640                                                                 fTPCnsigma_eta_hadrons[i]->Fill(fTPCnSigma, ceta);
1641                                                         }
1642                                                                 //electrons selection using E/p
1643                                                         if((fClus->E() / fP) >= fEoverPCutMin && (fClus->E() / fP) <= fEoverPCutMax) {
1644                                                                 fTPCnsigma_eta_electrons[i]->Fill(fTPCnSigma, ceta);
1645                                                         }
1646                                                 }
1647                                         }
1648
1649                                         if((fClus->E() / fP) >= fEoverPCutMin && (fClus->E() / fP) <= fEoverPCutMax)
1650                                         {
1651                                                 fTPCnsigma_eta->Fill(track->Eta(),fTPCnSigma);
1652                                                 fTPCnsigma_phi->Fill(track->Phi(),fTPCnSigma);
1653
1654                                                 if(fUseEMCal)
1655                                                 {
1656                                                         if(fTPCnSigma < 3.5 && fCorrelationFlag)
1657                                                         {
1658                                                                 fPtTrigger_Inc->Fill(fPt);
1659                                                                 DiHadronCorrelation(track, iTracks);
1660                                                         }
1661                                                 }
1662                                         }
1663
1664                                 }
1665                         }
1666                 }
1667
1668                         //______________________________________________________________
1669                         // Vertex
1670
1671                 fVtxZ[1]->Fill(fZvtx);
1672                 fNTracks[1]->Fill(fNOtrks);
1673                 fNClusters[1]->Fill(ClsNo);
1674                         //______________________________________________________________
1675
1676                         ///______________________________________________________________________
1677                         ///Histograms for PID Studies
1678                         //Double_t fPtBin[6] = {2,4,6,8,10,15};
1679
1680                 for(Int_t i = 0; i < 5; i++)
1681                 {
1682                         if(fPt>=fPtBin[i] && fPt<fPtBin[i+1])
1683                         {
1684                                 if(fEMCflag) fEoverP_tpc[i]->Fill(fTPCnSigma,(fClus->E() / fP));
1685
1686
1687                                 fTPC_pt[i]->Fill(fTPCsignal);
1688                                 fTPCnsigma_pt[i]->Fill(fTPCnSigma);
1689
1690                         }
1691                 }
1692
1693                         ///QA plots after track selection
1694                         ///_____________________________________________________________
1695
1696                         //_______________________________________________________
1697                         //Correlation Analysis - DiHadron
1698                 if(!fUseEMCal)
1699                 {
1700                         if(fTPCnSigma < 3.5 && fCorrelationFlag)
1701                         {
1702                                 fPtTrigger_Inc->Fill(fPt);
1703                                 DiHadronCorrelation(track, iTracks);
1704                         }
1705                 }
1706                         //_______________________________________________________
1707
1708
1709                         ///////////////////////////////////////////////////////////////////
1710                         ///TPC - efficiency calculation //
1711
1712                         /// changing start here
1713                 if(fIsMC && fIsAOD && track->GetLabel()>=0)
1714                 {
1715                         fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
1716                         Int_t pdg = fMCparticle->GetPdgCode();
1717
1718                                 //
1719                         if(fMCparticle->IsPhysicalPrimary()){
1720
1721
1722                                 if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax ){
1723
1724                                         Bool_t MotherFound = FindMother(track->GetLabel());
1725                                         if(MotherFound){
1726                                                 fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
1727                                         if( TMath::Abs(pdg) == 11 && fMCparticleMother->GetPdgCode()!=22 ){
1728                                                         if(fIsHFE1) fPtMC_TPC_All->Fill(fMCparticle->Pt());
1729                                         }
1730                                         }
1731                                 }
1732                         }
1733                 }///until here
1734
1735                 else if(fIsMC && track->GetLabel()>=0)//ESD
1736                 {
1737
1738                         if(fMCstack->IsPhysicalPrimary(track->GetLabel())){
1739                                 fMCtrack = fMCstack->Particle(track->GetLabel());
1740
1741                                 Int_t pdg = fMCtrack->GetPdgCode();
1742                                 if(fMCtrack->Eta()>=fEtaCutMin && fMCtrack->Eta()<=fEtaCutMax ){
1743
1744                                         if(fMCtrack->GetFirstMother()>0){
1745                                             fMCtrackMother = fMCstack->Particle(fMCtrack->GetFirstMother());
1746                                         if( TMath::Abs(pdg) == 11 && fMCtrackMother->GetPdgCode()!=22 ){
1747                                                         if(fIsHFE1) fPtMC_TPC_All->Fill(fMCtrack->Pt());
1748                                         }
1749                                         }
1750                                 }
1751                         }
1752                 }
1753
1754
1755                         ///________________________________________________________________________
1756                         ///PID
1757                         ///Here the PID cuts defined in the file "ConfigEMCalHFEpA.C" is applied
1758             Int_t pidpassed = 1;
1759             AliHFEpidObject hfetrack;
1760             hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);
1761             hfetrack.SetRecTrack(track);
1762             hfetrack.SetPP();   //proton-proton analysis
1763             if(!fPID->IsSelected(&hfetrack, NULL, "", fPIDqa)) pidpassed = 0;
1764             fpid->Fill(pidpassed);
1765
1766             if(pidpassed==0) continue;
1767                         ///________________________________________________________________________
1768
1769
1770                         ////////////////////////////////////////////////////////////////////
1771                         ///TPC efficiency calculations
1772
1773                         /// changing start here
1774                 if(fIsMC && fIsAOD && track->GetLabel()>=0)
1775                 {
1776                         fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
1777                         Int_t pdg = fMCparticle->GetPdgCode();
1778
1779                                 //
1780                         if(fMCparticle->IsPhysicalPrimary()){
1781
1782
1783                                 if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax ){
1784
1785                                         Bool_t MotherFound = FindMother(track->GetLabel());
1786                                         if(MotherFound){
1787                                                 fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
1788                                         if( TMath::Abs(pdg) == 11 && fMCparticleMother->GetPdgCode()!=22 ){
1789                                                         if(fIsHFE1) fPtMC_TPC_Selected->Fill(fMCparticle->Pt());
1790                                         }
1791                                         }
1792                                 }
1793                         }
1794                 }///until here
1795
1796                 else if(fIsMC && track->GetLabel()>=0)//ESD
1797                 {
1798
1799                         if(fMCstack->IsPhysicalPrimary(track->GetLabel())){
1800                                 fMCtrack = fMCstack->Particle(track->GetLabel());
1801
1802                                 Int_t pdg = fMCtrack->GetPdgCode();
1803                                 if(fMCtrack->Eta()>=fEtaCutMin && fMCtrack->Eta()<=fEtaCutMax ){
1804
1805                                         if(fMCtrack->GetFirstMother()>0){
1806                                             fMCtrackMother = fMCstack->Particle(fMCtrack->GetFirstMother());
1807                                         if( TMath::Abs(pdg) == 11 && fMCtrackMother->GetPdgCode()!=22 ){
1808                                                         if(fIsHFE1) fPtMC_TPC_Selected->Fill(fMCtrack->Pt());
1809                                         }
1810                                         }
1811                                 }
1812                         }
1813                 }
1814
1815                         //Eta Cut for TPC only
1816                 if(track->Eta()>=fEtaCutMin && track->Eta()<=fEtaCutMax ){
1817                         fTPCnsigma_pt_2D->Fill(fPt,fTPCnSigma);
1818                 }
1819
1820                         //Background for TPC only
1821                 if(fFillBackground){
1822                         Background(track, iTracks, Vtrack, kTRUE); //IsTPConly=kTRUE
1823                 }
1824
1825
1826                 fTPCnsigma_p[2]->Fill(fP,fTPCnSigma);
1827                 fTPC_p[2]->Fill(fP,fTPCsignal);
1828                 TPCNcls = track->GetTPCNcls();
1829                 Float_t pos3[3]={0,0,0};
1830
1831                 if(track->GetEMCALcluster()>0)
1832                 {
1833                         fClus = fVevent->GetCaloCluster(track->GetEMCALcluster());
1834                         if(fClus->IsEMCAL())
1835                         {
1836
1837                                         //________________________________________________________________________
1838
1839
1840                                         //Cluster timing distribution
1841                                 if(fUseEMCal && !fIsAOD ){
1842                                         AliESDCaloCells &cells=*(fESD->GetEMCALCells());
1843                                                 //      Int_t nTotalCells = cells.GetNumberOfCells() ;
1844                                                 //Int_t type        = cells.GetType();
1845                                                 //for (Int_t icell=  0; icell <  nTotalCells; icell++) {
1846                                                 //fTime->Fill(cells.GetTime(icell));
1847                                                 //}
1848
1849                                         TRefArray* caloClusters = new TRefArray();
1850                                         fESD->GetEMCALClusters(caloClusters);
1851
1852                                         Int_t nclus = caloClusters->GetEntries();
1853
1854                                                 //fClusESD = fESD->GetCaloCluster(track->GetEMCALcluster());
1855
1856                                         for (Int_t icl = 0; icl < nclus; icl++) {
1857
1858                                                 AliESDCaloCluster* clus = (AliESDCaloCluster*)caloClusters->At(icl);
1859
1860                                                 if(fClus->IsEMCAL()){
1861                                                         Float_t ncells  = fClus->GetNCells();
1862                                                         UShort_t * index = clus->GetCellsAbsId() ;
1863                                                         UShort_t * index2 = fClus->GetCellsAbsId() ;
1864
1865
1866                                                         for(Int_t i = 0; i < ncells ; i++){
1867
1868                                                                 Int_t absId =   index[i];
1869                                                                 fTime->Fill(fPt,cells.GetCellTime(absId));
1870
1871                                                                 Int_t absId2 =   index2[i];
1872                                                                 fTime2->Fill(fPt,cells.GetCellTime(absId2));
1873                                                         }
1874
1875                                                 }
1876                                         }
1877
1878
1879
1880
1881                                 }
1882
1883                                 if(fUseEMCal){
1884                                         double emctof = fClus->GetTOF();
1885                                         ftimingEle->Fill(fPt,emctof);
1886                                 }
1887                                         //________________________________________________________________________
1888
1889
1890
1891
1892                                         /////////////// Residuals
1893                                 Double_t Dx = fClus->GetTrackDx();
1894                                 Double_t Dz = fClus->GetTrackDz();
1895                                 Double_t R=TMath::Sqrt(Dx*Dx+Dz*Dz);
1896                                 for(Int_t i = 0; i < 5; i++)
1897                                 {
1898                                         if(fPt>=fPtBin[i] && fPt<fPtBin[i+1])
1899                                         {
1900
1901                                                 fEta[i]->Fill(Dz);
1902                                                 fPhi[i]->Fill(Dx);
1903                                                 fR[i]->Fill(R);
1904                                         }
1905                                 }
1906
1907                                 if(TMath::Abs(fClus->GetTrackDx())<=fdPhiCut && TMath::Abs(fClus->GetTrackDz())<=fdEtaCut)
1908                                 {
1909                                         Float_t Energy  = fClus->E();
1910                                         Float_t EoverP  = Energy/track->P();
1911                                         Float_t M02     = fClus->GetM02();
1912                                         Float_t M20     = fClus->GetM20();
1913                                         Float_t ncells  = fClus->GetNCells();
1914                                                 //----------------------------------------------------------------------------------------
1915                                                 // EtaCut electrons histogram
1916                                                 //Shower Shape Cut
1917                                         if(track->Eta()>=fEtaCutMin && track->Eta()<=fEtaCutMax ){
1918
1919                                                 if(fUseShowerShapeCut){
1920                                                         if(M02 >= fM02CutMin && M02<=fM02CutMax && M20>=fM20CutMin && M20<=fM20CutMax){
1921                                                                 fEoverP_pt[2]->Fill(fPt,(fClus->E() / fP));
1922                                                                 fShowerShapeCut->Fill(M02,M20);
1923
1924                                                         }
1925
1926                                                 }
1927                                                 if(!fUseShowerShapeCut){
1928                                                         fEoverP_pt[2]->Fill(fPt,(fClus->E() / fP));
1929                                                         fShowerShapeCut->Fill(M02,M20);
1930
1931                                                 }
1932                                                 if(fUseEMCal) fShowerShape_ele->Fill(M02,M20);
1933
1934                                                         //for shower shape cut studies - now with TPC PID Cut
1935                                                 if(fUseEMCal){
1936                                                         fShowerShapeM02_EoverP->Fill(M02,EoverP);
1937                                                         fShowerShapeM20_EoverP->Fill(M20,EoverP);
1938                                                 }
1939
1940                                         }
1941
1942                                                 //----------------------------------------------------------------------------------------
1943
1944
1945
1946                                                 /////////////// for Eta Phi distribution
1947                                         fClus->GetPosition(pos3);
1948                                         TVector3 vpos(pos3[0],pos3[1],pos3[2]);
1949                                         Double_t cphi = vpos.Phi();
1950                                         Double_t ceta = vpos.Eta();
1951                                         fEtaPhi[2]->Fill(cphi,ceta);
1952
1953
1954
1955                                         fTPCNcls_EoverP[2]->Fill(TPCNcls, EoverP);
1956
1957                                         for(Int_t i = 0; i < 5; i++)
1958                                         {
1959                                                 if(fPt>=fPtBin[i] && fPt<fPtBin[i+1])
1960                                                 {
1961
1962                                                         fR_EoverP[i]->Fill(R, EoverP);
1963                                                         fNcells[i]->Fill(ncells);
1964                                                         fNcells_EoverP[i]->Fill(EoverP, ncells);
1965                                                         fM02_EoverP[i]->Fill(M02,EoverP);
1966                                                         fM20_EoverP[i]->Fill(M20,EoverP);
1967                                                         fECluster_ptbins[i]->Fill(Energy);
1968                                                         fEoverP_ptbins[i]->Fill(EoverP);
1969
1970                                                         if((fClus->E() / fP) >= fEoverPCutMin && (fClus->E() / fP) <= fEoverPCutMax) {
1971                                                                 fNcells_electrons[i]->Fill(ncells);
1972                                                         }
1973
1974                                                         if(M02<0.5 && M20<0.3) {
1975                                                                 fEoverP_wSSCut[i]->Fill(EoverP);
1976                                                         }
1977                                                 }
1978                                         }
1979
1980                                         fNcells_pt->Fill(fPt, ncells);
1981
1982
1983                                                 ////////////////////////////////////////////////////////////////////
1984                                                 ///EMCal - Efficiency calculations
1985
1986                                                 /// changing start here
1987                                         if(fIsMC && fIsAOD && track->GetLabel()>=0)
1988                                         {
1989                                                 fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
1990                                                 Int_t pdg = fMCparticle->GetPdgCode();
1991
1992                                                         //
1993                                                 if(fMCparticle->IsPhysicalPrimary()){
1994
1995
1996                                                         if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax ){
1997
1998                                                                 Bool_t MotherFound = FindMother(track->GetLabel());
1999                                                                 if(MotherFound){
2000                                                                         fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
2001                                                                         if( TMath::Abs(pdg) == 11 && fMCparticleMother->GetPdgCode()!=22 ){
2002                                                                                 fPtMC_EMCal_All->Fill(fMCparticle->Pt());
2003                                                                         }
2004                                                                 }
2005                                                         }
2006                                                 }
2007                                         }///until here
2008
2009                                         else if(fIsMC && track->GetLabel()>=0)//ESD
2010                                         {
2011
2012                                                 if(fMCstack->IsPhysicalPrimary(track->GetLabel()))
2013                                                 {
2014
2015                                                         fMCtrack = fMCstack->Particle(track->GetLabel());
2016
2017                                                         Int_t pdg = fMCtrack->GetPdgCode();
2018                                                         if(fMCtrack->Eta()>=fEtaCutMin && fMCtrack->Eta()<=fEtaCutMax && fMCtrack->Phi()>=(TMath::Pi()*80/180) && fMCtrack->Phi()<=TMath::Pi())
2019                                                         {
2020                                                                 if(fMCtrack->GetFirstMother()>0){
2021                                                                         fMCtrackMother = fMCstack->Particle(fMCtrack->GetFirstMother());
2022                                                                         if( TMath::Abs(pdg) == 11 && fMCtrackMother->GetPdgCode()!=22 ){
2023
2024                                                                                 fPtMC_EMCal_All->Fill(fMCtrack->Pt());
2025                                                                         }
2026                                                                 }
2027                                                         }
2028                                                 }
2029                                         }
2030
2031
2032                                         if((fClus->E() / fP) >= fEoverPCutMin && (fClus->E() / fP) <= fEoverPCutMax)
2033                                         {
2034
2035                                             fECluster[2]->Fill(Energy);
2036                                                         //_______________________________________________________
2037                                                         //Correlation Analysis
2038                                                 if(fUseEMCal)
2039                                                 {
2040                                                         fPtElec_Inc->Fill(fPt);
2041                                                                 //new function to fill non-HFE histos
2042                                                         if(fFillBackground){
2043                                                                 Background(track, iTracks, Vtrack, kFALSE);
2044                                                         }
2045
2046                                                         double emctof2 = fClus->GetTOF();
2047                                                         ftimingEle2->Fill(fPt,emctof2);
2048
2049                                                         if(fCorrelationFlag)
2050                                                         {
2051                                                                 ElectronHadronCorrelation(track, iTracks, Vtrack);
2052                                                         }
2053                                                 }
2054                                                         //_______________________________________________________
2055
2056                                                         ////////////////////////////////////////////////////////////////////
2057                                                         ///EMCal - efficiency calculations
2058
2059                                                 if(track->Charge()<0)  fCharge_n->Fill(fPt);
2060                                                 if(track->Charge()>0)  fCharge_p->Fill(fPt);
2061
2062
2063                                                         /// changing start here
2064                                                 if(fIsMC && fIsAOD && track->GetLabel()>=0)
2065                                                 {
2066                                                         if(track->Charge()<0)  fCharge_n->Fill(fPt);
2067                                                         if(track->Charge()>0)  fCharge_p->Fill(fPt);
2068
2069                                                         fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
2070                                                         Int_t pdg = fMCparticle->GetPdgCode();
2071
2072                                                         if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax ){
2073
2074                                                                 if( TMath::Abs(pdg) == 11) fPtMCelectronAfterAll_nonPrimary->Fill(fMCparticle->Pt()); //numerator for the total efficiency, non Primary track
2075                                                         }
2076
2077
2078                                                         //
2079                                                         if(fMCparticle->IsPhysicalPrimary()){
2080
2081
2082                                                                 if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax ){
2083                                                                                 //check
2084                                                                         if(fIsHFE1) fPtMCelectronAfterAll->Fill(fMCparticle->Pt()); //numerator for the total efficiency
2085
2086
2087                                                                         Bool_t MotherFound = FindMother(track->GetLabel());
2088                                                                         if(MotherFound){
2089                                                                                 fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
2090                                                                                 if( TMath::Abs(pdg) == 11 && fMCparticleMother->GetPdgCode()!=22 ){
2091                                                                                         fPtMC_EMCal_Selected->Fill(fMCparticle->Pt());
2092                                                                                 }
2093                                                                         }
2094                                                                 }
2095                                                         }
2096                                                 }///until here
2097
2098                                                 else if(fIsMC && track->GetLabel()>=0)//ESD
2099                                                 {
2100                                                         if(track->Charge()<0)  fCharge_n->Fill(fPt);
2101                                                         if(track->Charge()>0)  fCharge_p->Fill(fPt);
2102
2103                                                         fMCtrack = fMCstack->Particle(track->GetLabel());
2104                                                         Int_t pdg = fMCtrack->GetPdgCode();
2105
2106                                                         if(fMCtrack->Eta()>=fEtaCutMin && fMCtrack->Eta()<=fEtaCutMax){
2107                                                                 if( TMath::Abs(pdg) == 11) fPtMCelectronAfterAll_nonPrimary->Fill(fMCtrack->Pt()); //numerator for the total efficiency, non Primary track
2108                                                         }
2109
2110                                                         if(fMCstack->IsPhysicalPrimary(track->GetLabel()))
2111                                                         {
2112                                                                 Bool_t MotherFound = FindMother(track->GetLabel());
2113
2114                                                                 if(MotherFound)
2115                                                                 {
2116                                                                         if(fMCtrack->Eta()>=fEtaCutMin && fMCtrack->Eta()<=fEtaCutMax){
2117                                                                                 if(fIsHFE1) fPtMCelectronAfterAll->Fill(fMCtrack->Pt()); //numerator for the total efficiency
2118                                                                         }
2119                                                                 }
2120
2121
2122
2123
2124                                                                 if(fMCtrack->Eta()>=fEtaCutMin && fMCtrack->Eta()<=fEtaCutMax && fMCtrack->Phi()>=(TMath::Pi()*80/180) && fMCtrack->Phi()<=TMath::Pi())
2125                                                                 {
2126                                                                         if(fMCtrack->GetFirstMother()>0){
2127                                                                                 fMCtrackMother = fMCstack->Particle(fMCtrack->GetFirstMother());
2128                                                                                 if( TMath::Abs(pdg) == 11 && fMCtrackMother->GetPdgCode()!=22 ){
2129
2130                                                                                         fPtMC_EMCal_Selected->Fill(fMCtrack->Pt());
2131                                                                                 }
2132                                                                         }
2133                                                                 }
2134                                                         }
2135                                                 }
2136                                                         ///////////////////////////////////////////////////////////////////
2137
2138
2139                                         }
2140                                 }
2141                         }
2142                 }
2143
2144                         //______________________________________________________________
2145                         // Vertex
2146
2147                 fVtxZ[2]->Fill(fZvtx);
2148                 fNTracks[2]->Fill(fNOtrks);
2149                 fNClusters[2]->Fill(ClsNo);
2150
2151                         //______________________________________________________________
2152
2153                         //_______________________________________________________
2154                         //Correlation Analysis
2155                 if(!fUseEMCal)
2156                 {
2157                         fPtElec_Inc->Fill(fPt);
2158
2159                         if(fCorrelationFlag)
2160                         {
2161                                 ElectronHadronCorrelation(track, iTracks, Vtrack);
2162                         }
2163                 }
2164                         //_______________________________________________________
2165
2166                         ///________________________________________________________________________
2167         }
2168
2169                 //__________________________________________________________________
2170                 //Event Mixing Analysis
2171                 //Filling pool
2172         if(fEventMixingFlag)
2173         {
2174                 fPool = fPoolMgr->GetEventPool(fCentrality->GetCentralityPercentile("V0A"), fZvtx); // Get the buffer associated with the current centrality and z-vtx
2175
2176                 if(!fPool) AliFatal(Form("No pool found for centrality = %f, zVtx = %f", fCentrality->GetCentralityPercentile("V0A"), fZvtx));
2177
2178                 fPool->UpdatePool(SelectedHadrons()); // fill the tracks in the event buffer. The ownership is handed over to the event mixing class. We are not allowed to delete tracksClone anymore!
2179
2180
2181         }
2182
2183                 //__________________________________________________________________
2184
2185         delete fListOfmotherkink;
2186         PostData(1, fOutputList);
2187 }
2188
2189         //______________________________________________________________________
2190 void AliAnalysisTaskEMCalHFEpA::Terminate(Option_t *)
2191 {
2192                 //Draw result to the screen
2193                 //Called once at the end of the query
2194
2195         fOutputList = dynamic_cast<TList*> (GetOutputData(1));
2196
2197         if(!fOutputList)
2198         {
2199                 printf("ERROR: Output list not available\n");
2200                 return;
2201         }
2202 }
2203
2204         //______________________________________________________________________
2205 Bool_t AliAnalysisTaskEMCalHFEpA::ProcessCutStep(Int_t cutStep, AliVParticle *track)
2206 {
2207                 //Check single track cuts for a given cut step
2208                 //Note this function is called inside the UserExec function
2209         const Int_t kMCOffset = AliHFEcuts::kNcutStepsMCTrack;
2210         if(!fCFM->CheckParticleCuts(cutStep + kMCOffset, track)) return kFALSE;
2211         return kTRUE;
2212 }
2213
2214
2215         //______________________________________________________________________
2216
2217
2218 void AliAnalysisTaskEMCalHFEpA::Background(AliVTrack *track, Int_t trackIndex, AliVParticle *vtrack, Bool_t IsTPConly)
2219 {
2220                 ///_________________________________________________________________
2221                 ///MC analysis
2222         if(fIsMC)
2223         {
2224                 if(track->GetLabel() < 0)
2225         {
2226                         AliWarning(Form("The track %d does not have a valid MC label",trackIndex));
2227                         return;
2228         }
2229
2230                 if(fIsAOD)
2231                 {
2232                         fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
2233
2234                         if(fMCparticle->GetMother()<0) return;
2235
2236                 fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
2237
2238                 if(TMath::Abs(fMCparticle->GetPdgCode())==11 && (TMath::Abs(fMCparticleMother->GetPdgCode())==22 || TMath::Abs(fMCparticleMother->GetPdgCode())==111 || TMath::Abs(fMCparticleMother->GetPdgCode())==221))
2239                 {
2240                                         //Is Background
2241                                 if(!IsTPConly)fPtBackgroundBeforeReco->Fill(track->Pt());
2242                                 if(IsTPConly)fPtBackgroundBeforeReco2->Fill(track->Pt());
2243                         }
2244                 }
2245                 else
2246                 {
2247                 fMCtrack = fMCstack->Particle(track->GetLabel());
2248
2249                 if(fMCtrack->GetFirstMother()<0) return;
2250
2251                 fMCtrackMother = fMCstack->Particle(fMCtrack->GetFirstMother());
2252
2253                 if(TMath::Abs(fMCtrack->GetPdgCode())==11 && (TMath::Abs(fMCtrackMother->GetPdgCode())==22 || TMath::Abs(fMCtrackMother->GetPdgCode())==111 || TMath::Abs(fMCtrackMother->GetPdgCode())==221))
2254                 {
2255                                         //Is Background
2256                                 if(!IsTPConly)fPtBackgroundBeforeReco->Fill(track->Pt());
2257                                 if(IsTPConly)fPtBackgroundBeforeReco2->Fill(track->Pt());
2258                         }
2259                 }
2260         }
2261                 ///_________________________________________________________________
2262
2263                 //________________________________________________
2264                 //Associated particle cut
2265         fPartnerCuts->SetAcceptKinkDaughters(kFALSE);
2266     fPartnerCuts->SetRequireITSRefit(kTRUE);
2267     fPartnerCuts->SetRequireTPCRefit(kTRUE);
2268     fPartnerCuts->SetEtaRange(-0.9,0.9);
2269     fPartnerCuts->SetMaxChi2PerClusterTPC(4.0);
2270     fPartnerCuts->SetMinNClustersTPC(80);
2271     fPartnerCuts->SetPtRange(0.3,1e10);
2272                 //fPartnerCuts->SetRequireSigmaToVertex(kTRUE);
2273                 //fPartnerCuts->SetMaxDCAToVertexXY(1);
2274                 //fPartnerCuts->SetMaxDCAToVertexZ(3);
2275                 //_________________________________________________
2276
2277                 ///#################################################################
2278                 //Non-HFE reconstruction
2279         fNonHFE = new AliSelectNonHFE();
2280         fNonHFE->SetAODanalysis(fIsAOD);
2281         if(fMassCutFlag) fNonHFE->SetInvariantMassCut(fMassCut);
2282         if(fAngleCutFlag) fNonHFE->SetOpeningAngleCut(fAngleCut);
2283         if(fChi2CutFlag) fNonHFE->SetChi2OverNDFCut(fChi2Cut);
2284         if(fDCAcutFlag) fNonHFE->SetDCACut(fDCAcut);
2285         fNonHFE->SetAlgorithm("DCA"); //KF
2286         fNonHFE->SetPIDresponse(fPidResponse);
2287         fNonHFE->SetTrackCuts(-3.5,3.5,fPartnerCuts);
2288         fNonHFE->SetAdditionalCuts(fPtMinAsso,fTpcNclsAsso);
2289
2290         if(!IsTPConly){
2291                 fNonHFE->SetHistAngleBack(fOpAngleBack);
2292                 fNonHFE->SetHistAngle(fOpAngle);
2293                 fNonHFE->SetHistDCABack(fDCABack);
2294                 fNonHFE->SetHistDCA(fDCA);
2295                 fNonHFE->SetHistMassBack(fInvMassBack);
2296                 fNonHFE->SetHistMass(fInvMass);
2297         }
2298         if(IsTPConly){
2299                 fNonHFE->SetHistAngleBack(fOpAngleBack2);
2300                 fNonHFE->SetHistAngle(fOpAngle2);
2301                 fNonHFE->SetHistDCABack(fDCABack2);
2302                 fNonHFE->SetHistDCA(fDCA2);
2303                 fNonHFE->SetHistMassBack(fInvMassBack2);
2304                 fNonHFE->SetHistMass(fInvMass2);
2305         }
2306
2307         fNonHFE->FindNonHFE(trackIndex,vtrack,fVevent);
2308
2309
2310
2311                 //Electron Information
2312         Double_t fPhiE = -999;
2313         Double_t fEtaE = -999;
2314         Double_t fPtE = -999;
2315         fPhiE = track->Phi();
2316         fEtaE = track->Eta();
2317         fPtE = track->Pt();
2318
2319                 ///_________________________________________________________________
2320                 ///MC analysis
2321         if(fIsMC)
2322         {
2323                 if(fIsAOD)
2324                 {
2325                         if(TMath::Abs(fMCparticle->GetPdgCode())==11 && (TMath::Abs(fMCparticleMother->GetPdgCode())==22 || TMath::Abs(fMCparticleMother->GetPdgCode())==111 || TMath::Abs(fMCparticleMother->GetPdgCode())==221))
2326                 {
2327                                 if(!IsTPConly){
2328                                         if(fNonHFE->IsULS()) fPtElec_ULS->Fill(fPtE,fNonHFE->GetNULS());
2329                                         if(fNonHFE->IsLS()) fPtElec_LS->Fill(fPtE,fNonHFE->GetNLS());
2330                                 }
2331
2332                                 if(IsTPConly){
2333                                         if(fNonHFE->IsULS()) fPtElec_ULS2->Fill(fPtE,fNonHFE->GetNULS());
2334                                         if(fNonHFE->IsLS()) fPtElec_LS2->Fill(fPtE,fNonHFE->GetNLS());
2335                                 }
2336                         }
2337                 }
2338                 else
2339                 {
2340                         if(TMath::Abs(fMCtrack->GetPdgCode())==11 && (TMath::Abs(fMCtrackMother->GetPdgCode())==22 || TMath::Abs(fMCtrackMother->GetPdgCode())==111 || TMath::Abs(fMCtrackMother->GetPdgCode())==221))
2341                         {
2342                                 if(!IsTPConly){
2343                                         if(fNonHFE->IsULS()) fPtElec_ULS->Fill(fPtE,fNonHFE->GetNULS());
2344                                         if(fNonHFE->IsLS()) fPtElec_LS->Fill(fPtE,fNonHFE->GetNLS());
2345                                 }
2346
2347                                 if(IsTPConly){
2348                                         if(fNonHFE->IsULS()) fPtElec_ULS2->Fill(fPtE,fNonHFE->GetNULS());
2349                                         if(fNonHFE->IsLS()) fPtElec_LS2->Fill(fPtE,fNonHFE->GetNLS());
2350                                 }
2351                         }
2352                 }
2353         }
2354                 ///_________________________________________________________________
2355         else
2356         {
2357                 if(!IsTPConly){
2358                         if(fNonHFE->IsULS()) fPtElec_ULS->Fill(fPtE,fNonHFE->GetNULS());
2359                         if(fNonHFE->IsLS()) fPtElec_LS->Fill(fPtE,fNonHFE->GetNLS());
2360                 }
2361
2362                 if(IsTPConly){
2363                         if(fNonHFE->IsULS()) fPtElec_ULS2->Fill(fPtE,fNonHFE->GetNULS());
2364                         if(fNonHFE->IsLS()) fPtElec_LS2->Fill(fPtE,fNonHFE->GetNLS());
2365                 }
2366         }
2367
2368
2369
2370 }
2371
2372
2373         //______________________________________________________________________
2374 void AliAnalysisTaskEMCalHFEpA::ElectronHadronCorrelation(AliVTrack *track, Int_t trackIndex, AliVParticle *vtrack)
2375 {
2376
2377                 ///_________________________________________________________________
2378                 ///MC analysis
2379         if(fIsMC)
2380         {
2381                 if(track->GetLabel() < 0)
2382         {
2383                         AliWarning(Form("The track %d does not have a valid MC label",trackIndex));
2384                         return;
2385         }
2386
2387                 if(fIsAOD)
2388                 {
2389                         fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
2390
2391                         if(fMCparticle->GetMother()<0) return;
2392
2393                 fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
2394
2395                 if(TMath::Abs(fMCparticle->GetPdgCode())==11 && (TMath::Abs(fMCparticleMother->GetPdgCode())==22 || TMath::Abs(fMCparticleMother->GetPdgCode())==111 || TMath::Abs(fMCparticleMother->GetPdgCode())==221))
2396                 {
2397                                         //Is Background
2398                                 fPtBackgroundBeforeReco->Fill(track->Pt());
2399                         }
2400                 }
2401                 else
2402                 {
2403                 fMCtrack = fMCstack->Particle(track->GetLabel());
2404
2405                 if(fMCtrack->GetFirstMother()<0) return;
2406
2407                 fMCtrackMother = fMCstack->Particle(fMCtrack->GetFirstMother());
2408
2409                 if(TMath::Abs(fMCtrack->GetPdgCode())==11 && (TMath::Abs(fMCtrackMother->GetPdgCode())==22 || TMath::Abs(fMCtrackMother->GetPdgCode())==111 || TMath::Abs(fMCtrackMother->GetPdgCode())==221))
2410                 {
2411                                         //Is Background
2412                                 fPtBackgroundBeforeReco->Fill(track->Pt());
2413                         }
2414                 }
2415         }
2416                 ///_________________________________________________________________
2417
2418                 //________________________________________________
2419                 //Associated particle cut
2420         fPartnerCuts->SetAcceptKinkDaughters(kFALSE);
2421     fPartnerCuts->SetRequireITSRefit(kTRUE);
2422     fPartnerCuts->SetRequireTPCRefit(kTRUE);
2423     fPartnerCuts->SetEtaRange(-0.9,0.9);
2424     fPartnerCuts->SetMaxChi2PerClusterTPC(4.0);
2425     fPartnerCuts->SetMinNClustersTPC(80);
2426     fPartnerCuts->SetPtRange(0.3,1e10);
2427                 //fPartnerCuts->SetRequireSigmaToVertex(kTRUE);
2428                 //fPartnerCuts->SetMaxDCAToVertexXY(1);
2429                 //fPartnerCuts->SetMaxDCAToVertexZ(3);
2430                 //_________________________________________________
2431
2432                 ///#################################################################
2433                 //Non-HFE reconstruction
2434         fNonHFE = new AliSelectNonHFE();
2435         fNonHFE->SetAODanalysis(fIsAOD);
2436         if(fMassCutFlag) fNonHFE->SetInvariantMassCut(fMassCut);
2437         if(fAngleCutFlag) fNonHFE->SetOpeningAngleCut(fAngleCut);
2438         if(fChi2CutFlag) fNonHFE->SetChi2OverNDFCut(fChi2Cut);
2439         if(fDCAcutFlag) fNonHFE->SetDCACut(fDCAcut);
2440         fNonHFE->SetAlgorithm("DCA"); //KF
2441         fNonHFE->SetPIDresponse(fPidResponse);
2442         fNonHFE->SetTrackCuts(-3.5,3.5,fPartnerCuts);
2443         fNonHFE->SetAdditionalCuts(fPtMinAsso,fTpcNclsAsso);
2444
2445
2446         fNonHFE->SetHistAngleBack(fOpAngleBack);
2447         fNonHFE->SetHistAngle(fOpAngle);
2448         fNonHFE->SetHistDCABack(fDCABack);
2449         fNonHFE->SetHistDCA(fDCA);
2450         fNonHFE->SetHistMassBack(fInvMassBack);
2451         fNonHFE->SetHistMass(fInvMass);
2452
2453         fNonHFE->FindNonHFE(trackIndex,vtrack,fVevent);
2454
2455         Int_t *fUlsPartner = fNonHFE->GetPartnersULS();
2456         Int_t *fLsPartner = fNonHFE->GetPartnersLS();
2457         Bool_t fUlsIsPartner = kFALSE;
2458         Bool_t fLsIsPartner = kFALSE;
2459                 ///#################################################################
2460
2461
2462                 //Electron Information
2463         Double_t fPhiE = -999;
2464         Double_t fEtaE = -999;
2465         Double_t fPhiH = -999;
2466         Double_t fEtaH = -999;
2467         Double_t fDphi = -999;
2468         Double_t fDeta = -999;
2469         Double_t fPtE = -999;
2470         Double_t fPtH = -999;
2471
2472         Double_t pi = TMath::Pi();
2473
2474         fPhiE = track->Phi();
2475         fEtaE = track->Eta();
2476         fPtE = track->Pt();
2477
2478
2479                 ///_________________________________________________________________
2480                 ///MC analysis
2481         if(fIsMC)
2482         {
2483                 if(fIsAOD)
2484                 {
2485                         if(TMath::Abs(fMCparticle->GetPdgCode())==11 && (TMath::Abs(fMCparticleMother->GetPdgCode())==22 || TMath::Abs(fMCparticleMother->GetPdgCode())==111 || TMath::Abs(fMCparticleMother->GetPdgCode())==221))
2486                 {
2487                                 if(fNonHFE->IsULS()) fPtElec_ULS->Fill(fPtE,fNonHFE->GetNULS());
2488                                 if(fNonHFE->IsLS()) fPtElec_LS->Fill(fPtE,fNonHFE->GetNLS());
2489                         }
2490                 }
2491                 else
2492                 {
2493                         if(TMath::Abs(fMCtrack->GetPdgCode())==11 && (TMath::Abs(fMCtrackMother->GetPdgCode())==22 || TMath::Abs(fMCtrackMother->GetPdgCode())==111 || TMath::Abs(fMCtrackMother->GetPdgCode())==221))
2494                         {
2495                                 if(fNonHFE->IsULS()) fPtElec_ULS->Fill(fPtE,fNonHFE->GetNULS());
2496                                 if(fNonHFE->IsLS()) fPtElec_LS->Fill(fPtE,fNonHFE->GetNLS());
2497                         }
2498                 }
2499         }
2500                 ///_________________________________________________________________
2501         else
2502         {
2503                 if(fNonHFE->IsULS()) fPtElec_ULS->Fill(fPtE,fNonHFE->GetNULS());
2504                 if(fNonHFE->IsLS()) fPtElec_LS->Fill(fPtE,fNonHFE->GetNLS());
2505         }
2506
2507
2508
2509
2510                 //__________________________________________________________________
2511                 //Event Mixing Analysis - Hadron Loop
2512                 //Retrieve
2513         if(fEventMixingFlag)
2514         {
2515                 fPool = fPoolMgr->GetEventPool(fCentrality->GetCentralityPercentile("V0A"), fZvtx); // Get the buffer associated with the current centrality and z-vtx
2516
2517                 if(!fPool) AliFatal(Form("No pool found for centrality = %f, zVtx = %f",fCentrality->GetCentralityPercentile("V0A"), fZvtx));
2518
2519                 if(fPool->GetCurrentNEvents() >= 5) // start mixing when 5 events are in the buffer
2520                 {
2521                         fPoolNevents->Fill(fPool->GetCurrentNEvents());
2522
2523                         for (Int_t jMix = 0; jMix < fPool->GetCurrentNEvents(); jMix++)  // mix with each event in the buffer
2524                         {
2525                                 TObjArray* bgTracks = fPool->GetEvent(jMix);
2526
2527                                 for (Int_t kMix = 0; kMix < bgTracks->GetEntriesFast(); kMix++)  // mix with each track in the event
2528                                 {
2529                                         const AliEHCParticle* MixedTrack(dynamic_cast<AliEHCParticle*>(bgTracks->At(kMix)));
2530                                         if (NULL == MixedTrack) continue;
2531
2532                                         fPhiH = MixedTrack->Phi();
2533                                         fEtaH = MixedTrack->Eta();
2534                                         fPtH = MixedTrack->Pt();
2535
2536                                         if(fPtH>fPtE || fPtH<1) continue;
2537
2538                                         fDphi = fPhiE - fPhiH;
2539
2540                                         if (fDphi > 3*pi/2) fDphi = fDphi - 2*pi;
2541                                         if (fDphi < -pi/2)  fDphi = fDphi + 2*pi;
2542
2543                                         fDeta = fEtaE - fEtaH;
2544
2545                                         Double_t fPtBin[6] = {2,4,6,8,10,15};
2546
2547                                         for(Int_t i = 0; i < 5; i++)
2548                                         {
2549                                             if(fPtE>=fPtBin[i] && fPtE<fPtBin[i+1])
2550                                             {
2551                                                         fCEtaPhi_Inc_EM[i]->Fill(fDphi,fDeta);
2552
2553                                                         if(fNonHFE->IsULS()) fCEtaPhi_ULS_EM[i]->Fill(fDphi,fDeta);
2554                                                         if(fNonHFE->IsLS()) fCEtaPhi_LS_EM[i]->Fill(fDphi,fDeta);
2555
2556                                                         if(fNonHFE->IsULS()) fCEtaPhi_ULS_Weight_EM[i]->Fill(fDphi,fDeta,fNonHFE->GetNULS());
2557                                                         if(fNonHFE->IsLS()) fCEtaPhi_LS_Weight_EM[i]->Fill(fDphi,fDeta,fNonHFE->GetNLS());
2558                                             }
2559                                         }
2560
2561                                                 // TODO your code: do event mixing with current event and bgTracks
2562                                                 // note that usually the content filled now is weighted by 1 / pool->GetCurrentNEvents()
2563                                 }
2564                         }
2565                 }
2566         }
2567                 //__________________________________________________________________
2568
2569                 //__________________________________________________________________
2570                 //Same Event Analysis - Hadron Loop
2571         for(Int_t iTracks = 0; iTracks < fVevent->GetNumberOfTracks(); iTracks++)
2572         {
2573                 if(trackIndex==iTracks) continue;
2574
2575                 AliVParticle* Vtrack2 = fVevent->GetTrack(iTracks);
2576                 if (!Vtrack2)
2577                 {
2578                         printf("ERROR: Could not receive track %d\n", iTracks);
2579                         continue;
2580                 }
2581
2582                 AliVTrack *track2 = dynamic_cast<AliVTrack*>(Vtrack2);
2583
2584                 if(fIsAOD)
2585                 {
2586                         AliAODTrack *atrack2 = dynamic_cast<AliAODTrack*>(Vtrack2);
2587                         if(!atrack2->TestFilterMask(AliAODTrack::kTrkTPCOnly)) continue;
2588                         if((!(atrack2->GetStatus()&AliESDtrack::kITSrefit)|| (!(atrack2->GetStatus()&AliESDtrack::kTPCrefit)))) continue;
2589                         if(atrack2->GetTPCNcls() < 80) continue;
2590                 }
2591                 else
2592                 {
2593                         AliESDtrack *etrack2 = dynamic_cast<AliESDtrack*>(Vtrack2);
2594                         if(!fPartnerCuts->AcceptTrack(etrack2)) continue;
2595                 }
2596
2597                 fPhiH = track2->Phi();
2598                 fEtaH = track2->Eta();
2599                 fPtH = track2->Pt();
2600
2601                 if(fPtH>fPtE || fPtH<1) continue;
2602
2603                 fDphi = fPhiE - fPhiH;
2604
2605                 if (fDphi > 3*pi/2) fDphi = fDphi - 2*pi;
2606                 if (fDphi < -pi/2)  fDphi = fDphi + 2*pi;
2607
2608                 fDeta = fEtaE - fEtaH;
2609
2610                 Double_t fPtBin[6] = {2,4,6,8,10,15};
2611
2612                         //______________________________________________________________
2613                         //Check if this track is a Non-HFE partner
2614                 for(Int_t i = 0; i < fNonHFE->GetNULS(); i++)
2615                 {
2616                         if(fUlsPartner[i]==iTracks) fUlsIsPartner=kTRUE;
2617                 }
2618                 for(Int_t i = 0; i < fNonHFE->GetNLS(); i++)
2619                 {
2620                         if(fLsPartner[i]==iTracks) fLsIsPartner=kTRUE;
2621                 }
2622                         //______________________________________________________________
2623
2624                 for(Int_t i = 0; i < 5; i++)
2625                 {
2626                     if(fPtE>=fPtBin[i] && fPtE<fPtBin[i+1])
2627                     {
2628                                 fCEtaPhi_Inc[i]->Fill(fDphi,fDeta);
2629
2630                                 if(fNonHFE->IsULS()) fCEtaPhi_ULS[i]->Fill(fDphi,fDeta);
2631                                 if(fNonHFE->IsLS()) fCEtaPhi_LS[i]->Fill(fDphi,fDeta);
2632                                 if(fNonHFE->IsULS() && !fUlsIsPartner) fCEtaPhi_ULS_NoP[i]->Fill(fDphi,fDeta);
2633                                 if(fNonHFE->IsLS() && !fLsIsPartner) fCEtaPhi_LS_NoP[i]->Fill(fDphi,fDeta);
2634
2635                                 if(fNonHFE->IsULS()) fCEtaPhi_ULS_Weight[i]->Fill(fDphi,fDeta,fNonHFE->GetNULS());
2636                                 if(fNonHFE->IsLS()) fCEtaPhi_LS_Weight[i]->Fill(fDphi,fDeta,fNonHFE->GetNLS());
2637                                 if(fNonHFE->IsULS() && !fUlsIsPartner) fCEtaPhi_ULS_NoP_Weight[i]->Fill(fDphi,fDeta,fNonHFE->GetNULS());
2638                                 if(fNonHFE->IsLS() && !fLsIsPartner) fCEtaPhi_LS_NoP_Weight[i]->Fill(fDphi,fDeta,fNonHFE->GetNLS());
2639                     }
2640                 }
2641         }
2642 }
2643
2644         //____________________________________________________________________________________________________________
2645         //Create a TObjArray with selected hadrons, for the mixed event analysis
2646 TObjArray* AliAnalysisTaskEMCalHFEpA::SelectedHadrons()
2647 {
2648         fTracksClone = new TObjArray;
2649         fTracksClone->SetOwner(kTRUE);
2650
2651         for(Int_t iTracks = 0; iTracks < fVevent->GetNumberOfTracks(); iTracks++)
2652         {
2653                 AliVParticle* Vtrack2 = fVevent->GetTrack(iTracks);
2654                 if (!Vtrack2)
2655                 {
2656                         printf("ERROR: Could not receive track %d\n", iTracks);
2657                         continue;
2658                 }
2659
2660                 AliVTrack *track2 = dynamic_cast<AliVTrack*>(Vtrack2);
2661
2662                 if(fIsAOD)
2663                 {
2664                         AliAODTrack *atrack2 = dynamic_cast<AliAODTrack*>(Vtrack2);
2665                         if(!atrack2->TestFilterMask(AliAODTrack::kTrkTPCOnly)) continue;
2666                         if((!(atrack2->GetStatus()&AliESDtrack::kITSrefit)|| (!(atrack2->GetStatus()&AliESDtrack::kTPCrefit)))) continue;
2667                         if(atrack2->GetTPCNcls() < 80) continue;
2668                 }
2669                 else
2670                 {
2671                         AliESDtrack *etrack2 = dynamic_cast<AliESDtrack*>(Vtrack2);
2672                         if(!fPartnerCuts->AcceptTrack(etrack2)) continue;
2673                 }
2674
2675                 fTracksClone->Add(new AliEHCParticle(track2->Eta(), track2->Phi(), track2->Pt()));
2676         }
2677         return fTracksClone;
2678 }
2679         //____________________________________________________________________________________________________________
2680
2681         //______________________________________________________________________
2682 void AliAnalysisTaskEMCalHFEpA::DiHadronCorrelation(AliVTrack *track, Int_t trackIndex)
2683 {
2684                 //________________________________________________
2685                 //Associated particle cut
2686         fPartnerCuts->SetAcceptKinkDaughters(kFALSE);
2687     fPartnerCuts->SetRequireITSRefit(kTRUE);
2688     fPartnerCuts->SetRequireTPCRefit(kTRUE);
2689     fPartnerCuts->SetEtaRange(-0.9,0.9);
2690     fPartnerCuts->SetMaxChi2PerClusterTPC(4.0);
2691     fPartnerCuts->SetMinNClustersTPC(80);
2692     fPartnerCuts->SetPtRange(0.3,1e10);
2693                 //fPartnerCuts->SetRequireSigmaToVertex(kTRUE);
2694                 //fPartnerCuts->SetMaxDCAToVertexXY(1);
2695                 //fPartnerCuts->SetMaxDCAToVertexZ(3);
2696                 //_________________________________________________
2697
2698                 //Electron Information
2699         Double_t fPhiE = -999;
2700         Double_t fEtaE = -999;
2701         Double_t fPhiH = -999;
2702         Double_t fEtaH = -999;
2703         Double_t fDphi = -999;
2704         Double_t fDeta = -999;
2705         Double_t fPtE = -999;
2706         Double_t fPtH = -999;
2707
2708         Double_t pi = TMath::Pi();
2709
2710         fPhiE = track->Phi();
2711         fEtaE = track->Eta();
2712         fPtE = track->Pt();
2713
2714                 //__________________________________________________________________
2715                 //Same Event Analysis - Hadron Loop
2716         for(Int_t iTracks = 0; iTracks < fVevent->GetNumberOfTracks(); iTracks++)
2717         {
2718                 if(trackIndex==iTracks) continue;
2719
2720                 AliVParticle* Vtrack2 = fVevent->GetTrack(iTracks);
2721                 if (!Vtrack2)
2722                 {
2723                         printf("ERROR: Could not receive track %d\n", iTracks);
2724                         continue;
2725                 }
2726
2727                 AliVTrack *track2 = dynamic_cast<AliVTrack*>(Vtrack2);
2728
2729                 if(fIsAOD)
2730                 {
2731                         AliAODTrack *atrack2 = dynamic_cast<AliAODTrack*>(Vtrack2);
2732                         if(!atrack2->TestFilterMask(AliAODTrack::kTrkTPCOnly)) continue;
2733                         if((!(atrack2->GetStatus()&AliESDtrack::kITSrefit)|| (!(atrack2->GetStatus()&AliESDtrack::kTPCrefit)))) continue;
2734                         if(atrack2->GetTPCNcls() < 80) continue;
2735                 }
2736                 else
2737                 {
2738                         AliESDtrack *etrack2 = dynamic_cast<AliESDtrack*>(Vtrack2);
2739                         if(!fPartnerCuts->AcceptTrack(etrack2)) continue;
2740                 }
2741
2742                 fPhiH = track2->Phi();
2743                 fEtaH = track2->Eta();
2744                 fPtH = track2->Pt();
2745
2746                 if(fPtH>fPtE || fPtH<1) continue;
2747
2748                 fDphi = fPhiE - fPhiH;
2749
2750                 if (fDphi > 3*pi/2) fDphi = fDphi - 2*pi;
2751                 if (fDphi < -pi/2)  fDphi = fDphi + 2*pi;
2752
2753                 fDeta = fEtaE - fEtaH;
2754
2755                 Double_t fPtBin[6] = {2,4,6,8,10,15};
2756
2757                 for(Int_t i = 0; i < 5; i++)
2758                 {
2759                     if(fPtE>=fPtBin[i] && fPtE<fPtBin[i+1])
2760                     {
2761                                 fCEtaPhi_Inc_DiHadron[i]->Fill(fDphi,fDeta);
2762                     }
2763                 }
2764         }
2765 }
2766         //____________________________________________________________________________________________________________
2767
2768         //______________________________________________________________________
2769 Bool_t AliAnalysisTaskEMCalHFEpA::FindMother(Int_t mcIndex)
2770 {
2771         fIsHFE1 = kFALSE;
2772         fIsHFE2 = kFALSE;
2773         fIsNonHFE = kFALSE;
2774         fIsFromD = kFALSE;
2775         fIsFromB = kFALSE;
2776         fIsFromPi0 = kFALSE;
2777         fIsFromEta = kFALSE;
2778         fIsFromGamma = kFALSE;
2779
2780         if(mcIndex < 0 || !fIsMC)
2781         {
2782                 return kFALSE;
2783         }
2784
2785         Int_t pdg = -99999;
2786         Int_t mpdg = -99999;
2787         Int_t gmpdg = -99999;
2788         Int_t ggmpdg = -99999;
2789         Int_t gggmpdg = -99999;
2790
2791         if(fIsAOD)
2792         {
2793                 fMCparticle = (AliAODMCParticle*) fMCarray->At(mcIndex);
2794
2795                 pdg = TMath::Abs(fMCparticle->GetPdgCode());
2796
2797
2798                 if(pdg!=11)
2799                 {
2800                         fIsHFE1 = kFALSE;
2801                         fIsHFE2 = kFALSE;
2802                         fIsNonHFE = kFALSE;
2803                         fIsFromD = kFALSE;
2804                         fIsFromB = kFALSE;
2805                         fIsFromPi0 = kFALSE;
2806                         fIsFromEta = kFALSE;
2807                         fIsFromGamma = kFALSE;
2808                         return kFALSE;
2809                 }
2810
2811                 if(fMCparticle->GetMother()<0)
2812                 {
2813                         fIsHFE1 = kFALSE;
2814                         fIsHFE2 = kFALSE;
2815                         fIsNonHFE = kFALSE;
2816                         fIsFromD = kFALSE;
2817                         fIsFromB = kFALSE;
2818                         fIsFromPi0 = kFALSE;
2819                         fIsFromEta = kFALSE;
2820                         fIsFromGamma = kFALSE;
2821                         return kFALSE;
2822                 }
2823
2824                 fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
2825                 mpdg = TMath::Abs(fMCparticleMother->GetPdgCode());
2826
2827                 if(fMCparticleMother->GetMother()<0)
2828                 {
2829                         gmpdg = 0;
2830                         ggmpdg = 0;
2831                         gggmpdg = 0;
2832                 }
2833                 else
2834                 {
2835                         fMCparticleGMother = (AliAODMCParticle*) fMCarray->At(fMCparticleMother->GetMother());
2836                         gmpdg = TMath::Abs(fMCparticleGMother->GetPdgCode());
2837                         if(fMCparticleGMother->GetMother()<0)
2838                         {
2839                                 ggmpdg = 0;
2840                                 gggmpdg = 0;
2841                         }
2842                         else
2843                         {
2844                                 fMCparticleGGMother = (AliAODMCParticle*) fMCarray->At(fMCparticleGMother->GetMother());
2845                                 ggmpdg = TMath::Abs(fMCparticleGGMother->GetPdgCode());
2846                                 if(fMCparticleGGMother->GetMother()<0)
2847                                 {
2848                                         gggmpdg = 0;
2849                                 }
2850                                 else
2851                                 {
2852                                         fMCparticleGGGMother = (AliAODMCParticle*) fMCarray->At(fMCparticleGGMother->GetMother());
2853                                         gggmpdg = TMath::Abs(fMCparticleGGGMother->GetPdgCode());
2854                                 }
2855                         }
2856                 }
2857         }
2858         else
2859         {
2860                 fMCtrack = fMCstack->Particle(mcIndex);
2861
2862                 pdg = TMath::Abs(fMCtrack->GetPdgCode());
2863
2864                 if(pdg!=11)
2865                 {
2866                         fIsHFE1 = kFALSE;
2867                         fIsHFE2 = kFALSE;
2868                         fIsNonHFE = kFALSE;
2869                         fIsFromD = kFALSE;
2870                         fIsFromB = kFALSE;
2871                         fIsFromPi0 = kFALSE;
2872                         fIsFromEta = kFALSE;
2873                         fIsFromGamma = kFALSE;
2874                         return kFALSE;
2875                 }
2876
2877                 if(fMCtrack->GetFirstMother()<0)
2878                 {
2879                         fIsHFE1 = kFALSE;
2880                         fIsHFE2 = kFALSE;
2881                         fIsNonHFE = kFALSE;
2882                         fIsFromD = kFALSE;
2883                         fIsFromB = kFALSE;
2884                         fIsFromPi0 = kFALSE;
2885                         fIsFromEta = kFALSE;
2886                         fIsFromGamma = kFALSE;
2887                         return kFALSE;
2888                 }
2889
2890                 fMCtrackMother = fMCstack->Particle(fMCtrack->GetFirstMother());
2891                 mpdg = TMath::Abs(fMCtrackMother->GetPdgCode());
2892
2893                 if(fMCtrackMother->GetFirstMother()<0)
2894                 {
2895                         gmpdg = 0;
2896                         ggmpdg = 0;
2897                         gggmpdg = 0;
2898                 }
2899                 else
2900                 {
2901                         fMCtrackGMother = fMCstack->Particle(fMCtrackMother->GetFirstMother());
2902                         gmpdg = TMath::Abs(fMCtrackGMother->GetPdgCode());
2903
2904                         if(fMCtrackGMother->GetFirstMother()<0)
2905                         {
2906                                 ggmpdg = 0;
2907                                 gggmpdg = 0;
2908                         }
2909                         else
2910                         {
2911                                 fMCtrackGGMother = fMCstack->Particle(fMCtrackGMother->GetFirstMother());
2912                                 ggmpdg = TMath::Abs(fMCtrackGGMother->GetPdgCode());
2913
2914                                 if(fMCtrackGGMother->GetFirstMother()<0)
2915                                 {
2916                                         gggmpdg = 0;
2917                                 }
2918                                 else
2919                                 {
2920                                         fMCtrackGGGMother = fMCstack->Particle(fMCtrackGGMother->GetFirstMother());
2921                                         gggmpdg = TMath::Abs(fMCtrackGGGMother->GetPdgCode());
2922                                 }
2923                         }
2924                 }
2925         }
2926
2927         //Tag Electron Source
2928         if(mpdg==111 || mpdg==221 || mpdg==22)
2929         {
2930                 fIsHFE1 = kFALSE;
2931                 fIsHFE2 = kFALSE;
2932                 fIsNonHFE = kTRUE;
2933                 fIsFromD = kFALSE;
2934                 fIsFromB = kFALSE;
2935
2936                 fIsFromPi0 = kFALSE;
2937                 fIsFromEta = kFALSE;
2938                 fIsFromGamma = kFALSE;
2939
2940                 if(mpdg==111) fIsFromPi0 = kFALSE;
2941                 if(mpdg==221)fIsFromEta = kFALSE;
2942                 if(mpdg==22) fIsFromGamma = kFALSE;
2943
2944                 return kTRUE;
2945         }
2946         else
2947         {
2948                 fIsHFE1 = kFALSE;
2949                 fIsHFE2 = kTRUE;
2950
2951                 fIsFromPi0 = kFALSE;
2952                 fIsFromEta = kFALSE;
2953                 fIsFromGamma = kFALSE;
2954
2955                 fIsNonHFE = kFALSE;
2956
2957                 fIsFromD = kFALSE;
2958                 fIsFromB = kFALSE;
2959
2960                 if(mpdg>400 && mpdg<500)
2961                 {
2962                         if((gmpdg>500 && gmpdg<600) || (ggmpdg>500 && ggmpdg<600) || (gggmpdg>500 && gggmpdg<600))
2963                         {
2964                                 fIsHFE1 = kTRUE;
2965                                 fIsFromD = kFALSE;
2966                                 fIsFromB = kTRUE;
2967                                 return kTRUE;
2968                         }
2969                         else
2970                         {
2971                                 fIsHFE1 = kTRUE;
2972                                 fIsFromD = kTRUE;
2973                                 fIsFromB = kFALSE;
2974                                 return kTRUE;
2975                         }
2976                 }
2977                 else if(mpdg>500 && mpdg<600)
2978                 {
2979                         fIsHFE1 = kTRUE;
2980                         fIsFromD = kFALSE;
2981                         fIsFromB = kTRUE;
2982                         return kTRUE;
2983                 }
2984                 else
2985                 {
2986                         fIsHFE1 = kFALSE;
2987                         fIsFromD = kFALSE;
2988                         fIsFromB = kFALSE;
2989                         return kFALSE;
2990                 }
2991         }
2992 }