PWGHF/hfe/AliAnalysisTaskEMCalHFEpA.cxx

   1
   2 /**************************************************************************
   3  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   4  *                                                                        *
   5  * Author: The ALICE Off-line Project.                                    *
   6  * Contributors are mentioned in the code where appropriate.              *
   7  *                                                                        *
   8  * Permission to use, copy, modify and distribute this software and its   *
   9  * documentation strictly for non-commercial purposes is hereby granted   *
  10  * without fee, provided that the above copyright notice appears in all   *
  11  * copies and that both the copyright notice and this permission notice   *
  12  * appear in the supporting documentation. The authors make no claims     *
  13  * about the suitability of this software for any purpose. It is          *
  14  * provided "as is" without express or implied warranty.                  *
  15  **************************************************************************/
  16
  17         ////////////////////////////////////////////////////////////////////////
  18         //                                                                    //
  19         //      Task for Heavy-flavour electron analysis in pPb collisions    //
  20         //      (+ Electron-Hadron Jetlike Azimuthal Correlation)             //
  21         //                                                                                                                                        //
  22         //              version: September 16, 2014.                                                                  //
  23         //                                                                    //
  24         //          Authors                                                                               //
  25         //              Elienos Pereira de Oliveira Filho (epereira@cern.ch)          //
  26         //              Cristiane Jahnke                (cristiane.jahnke@cern.ch)                    //
  27         //                                                                    //
  28         ////////////////////////////////////////////////////////////////////////
  29
  30 #include "TChain.h"
  31 #include "TTree.h"
  32 #include "TNtuple.h"
  33 #include "TH1F.h"
  34 #include "TH2F.h"
  35 #include "TCanvas.h"
  36 #include "AliAnalysisTask.h"
  37 #include "AliAnalysisManager.h"
  38 #include "AliESDEvent.h"
  39 #include "AliAODEvent.h"
  40 #include "AliVEvent.h"
  41 #include "AliESDInputHandler.h"
  42 #include "AliESDtrackCuts.h"
  43 #include "AliESDCaloCluster.h"
  44 #include "AliESDCaloCells.h"
  45 #include "AliEMCALTrack.h"
  46 #include "AliExternalTrackParam.h"
  47 #include "AliPhysicsSelection.h"
  48 #include "TGeoGlobalMagField.h"
  49 #include "AliMagF.h"
  50 #include "AliLog.h"
  51 #include "AliStack.h"
  52 #include "AliCentrality.h"
  53 #include "AliAODMCParticle.h"
  54 #include "AliAODMCHeader.h"
  55 #include "AliPID.h"
  56 #include "AliPIDResponse.h"
  57 #include "AliHFEcontainer.h"
  58 #include "AliHFEcuts.h"
  59 #include "AliHFEpid.h"
  60 #include "AliHFEpidBase.h"
  61 #include "AliHFEpidQAmanager.h"
  62 #include "AliHFEtools.h"
  63 #include "AliCFContainer.h"
  64 #include "AliCFManager.h"
  65 #include "AliSelectNonHFE.h"
  66 #include "AliHFEpidTPC.h"
  67 #include "AliAnalysisTaskEMCalHFEpA.h"
  68 #include "TMath.h"
  69 #include "THnSparse.h"
  70 #include "TLorentzVector.h"
  71 #include "TString.h"
  72 #include "TFile.h"
  73 #include "AliESDHandler.h"
  74 #include "AliMCEventHandler.h"
  75 #include "AliMCEvent.h"
  76 #include "AliStack.h"
  77 #include "TParticle.h"
  78 #include "AliLog.h"
  79 #include "AliAnalysisTaskSE.h"
  80 #include "TRefArray.h"
  81 #include "TVector.h"
  82 #include "stdio.h"
  83 #include "TGeoManager.h"
  84 #include "iostream"
  85 #include "fstream"
  86 #include "AliKFParticle.h"
  87 #include "AliKFVertex.h"
  88 #include "AliVParticle.h"
  89 #include "AliVTrack.h"
  90 #include "AliEventPoolManager.h"
  91 #include "TObjArray.h"
  92         //include to use reader as Lucile does
  93 #include "AliCaloTrackAODReader.h"
  94 #include "AliCaloTrackReader.h"
  95 #include "AliEMCALRecoUtils.h" //to remove exotics
  96 #include "AliAODHeader.h"
  97 #include "AliEMCALGeometry.h"
  98
  99
 100
 101         // --- ANALYSIS system ---
 102 #include "AliCalorimeterUtils.h"
 103 #include "AliESDEvent.h"
 104 #include "AliMCEvent.h"
 105 #include "AliStack.h"
 106 #include "AliAODPWG4Particle.h"
 107 #include "AliVCluster.h"
 108 #include "AliVCaloCells.h"
 109 #include "AliMixedEvent.h"
 110 #include "AliAODCaloCluster.h"
 111 #include "AliOADBContainer.h"
 112 #include "AliAnalysisManager.h"
 113
 114         // --- Detector ---
 115 #include "AliEMCALGeometry.h"
 116 #include "AliPHOSGeoUtils.h"
 117
 118         //______________________________________________________________________
 119
 120         //______________________________________________________________________
 121 ClassImp(AliAnalysisTaskEMCalHFEpA)
 122
 123         //______________________________________________________________________
 124 AliAnalysisTaskEMCalHFEpA::AliAnalysisTaskEMCalHFEpA(const char *name)
 125 : AliAnalysisTaskSE(name)
 126 ,fCorrelationFlag(0)
 127 ,fIsMC(0)
 128 ,fUseEMCal(kFALSE)
 129 ,fUseTrigger(kFALSE)
 130 ,fUseTender(kFALSE)
 131 ,fUseShowerShapeCut(kFALSE)
 132 ,fFillBackground(kFALSE)
 133 ,fEoverPnsigma(kFALSE)
 134 ,fAssocWithSPD(kFALSE)
 135 ,fEMCEG1(kFALSE)
 136 ,fEMCEG2(kFALSE)
 137 ,fIsHFE1(kFALSE)
 138 ,fIsHFE2(kFALSE)
 139 ,fIsNonHFE(kFALSE)
 140 ,fIsFromD(kFALSE)
 141 ,fIsFromB(kFALSE)
 142 ,fIsFromPi0(kFALSE)
 143 ,fIsFromEta(kFALSE)
 144 ,fIsFromGamma(kFALSE)
 145 ,fESD(0)
 146 ,fAOD(0)
 147 ,fVevent(0)
 148 ,fPartnerCuts(new AliESDtrackCuts())
 149 ,fOutputList(0)
 150 ,fPidResponse(0)
 151 ,fNonHFE(new AliSelectNonHFE())
 152 ,fIsAOD(kFALSE)
 153 ,fCentrality(0)
 154 ,fCentralityMin(0)
 155 ,fCentralityMax(100)
 156 ,fHasCentralitySelection(kFALSE)
 157 ,fCentralityHist(0)
 158 ,fCentralityHistPass(0)
 159 ,fZvtx(0)
 160 ,fEstimator(0)
 161 ,fClus(0)
 162 //,fClusESD(0)
 163 ,fNevent(0)
 164 ,fNevent2(0)
 165 ,fPtElec_Inc(0)
 166 ,fPtPrim(0)
 167 ,fPtSec(0)
 168 ,fPtPrim2(0)
 169 ,fPtSec2(0)
 170 ,fCharge_n(0)
 171 ,fCharge_p(0)
 172 ,fTime(0)
 173 ,fTime2(0)
 174 ,ftimingEle(0)
 175 ,ftimingEle2(0)
 176 ,fPtElec_ULS(0)
 177 ,fPtElec_LS(0)
 178
 179 ,fPtElec_ULS_NoPid(0)
 180 ,fPtElec_LS_NoPid(0)
 181 ,fPtElec_ULS_MC(0)
 182 ,fPtElec_ULS_MC_weight(0)
 183 ,fPtElec_ULS2(0)
 184 ,fPtElec_LS2(0)
 185
 186 ,fPtElec_ULS_mc_closure(0)
 187 ,fPtElec_LS_mc_closure(0)
 188 ,fPtElec_ULS2_mc_closure(0)
 189 ,fPtElec_LS2_mc_closure(0)
 190
 191 ,fPtElec_ULS_weight(0)
 192 ,fPtElec_LS_weight(0)
 193 ,fPtElec_ULS2_weight(0)
 194 ,fPtElec_LS2_weight(0)
 195 ,fTOF01(0)
 196 ,fTOF02(0)
 197 ,fTOF03(0)
 198 ,fpid(0)
 199 ,fEoverP_pt_true_electrons(0)
 200 ,fEoverP_pt_true_hadrons(0)
 201 ,fEoverP_pt_true_electrons0(0)
 202 ,fEoverP_pt_true_hadrons0(0)
 203 ,fEoverP_pt(0)
 204 ,fEoverP_tpc(0)
 205 ,fEoverP_tpc_p_trigger(0)
 206 ,fEoverP_tpc_pt_trigger(0)
 207 ,fTPC_pt(0)
 208 ,fTPC_p(0)
 209 ,fTPCnsigma_pt(0)
 210 ,fTPCnsigma_p(0)
 211 ,fTPCnsigma_p_TPC(0)
 212 ,fTPCnsigma_p_TPC_on_EMCal_acc(0)
 213 ,fTPCnsigma_p_TPC_EoverP_cut(0)
 214 ,fTPCnsigma_pt_2D(0)
 215 ,fShowerShapeCut(0)
 216 ,fShowerShapeM02_EoverP(0)
 217 ,fShowerShapeM20_EoverP(0)
 218 ,fShowerShape_ha(0)
 219 ,fShowerShape_ele(0)
 220 ,fTPCnsigma_eta(0)
 221 ,fTPCnsigma_phi(0)
 222 ,fECluster(0)
 223 ,fECluster_pure(0)
 224 ,fECluster_not_exotic(0)
 225 ,fECluster_not_exotic1(0)
 226 ,fECluster_not_exotic2(0)
 227 ,fECluster_exotic(0)
 228 ,fNCluster_pure(0)
 229 ,fNCluster_pure_aod(0)
 230 ,fNCluster_ECluster(0)
 231 ,fNcells_energy(0)
 232 ,fNcells_energy_elec_selected(0)
 233 ,fNcells_energy_not_exotic(0)
 234
 235 ,fEtaPhi(0)
 236 ,fEtaPhi_num(0)
 237 ,fEtaPhi_den(0)
 238 ,fEtaPhi_data(0)
 239 ,fpt_reco_pt_MC_num(0)
 240 ,fpt_reco_pt_MC_den(0)
 241 ,fVtxZ(0)
 242
 243 ,fVtxZ_new1(0)
 244 ,fVtxZ_new2(0)
 245 ,fVtxZ_new3(0)
 246 ,fVtxZ_new4(0)
 247
 248 ,fzRes1(0)
 249 ,fzRes2(0)
 250 ,fSPD_track_vtx1(0)
 251 ,fSPD_track_vtx2(0)
 252
 253
 254 ,fEtad(0)
 255 ,fNTracks(0)
 256 ,fTrack_Multi(0)
 257 ,fNTracks_pt(0)
 258 ,fNTracks_eta(0)
 259 ,fNTracks_phi(0)
 260 ,fNClusters(0)
 261 ,fTPCNcls_EoverP(0)
 262 ,fTPCNcls_pid(0)
 263 ,fEta(0)
 264 ,fPhi(0)
 265 ,fR(0)
 266 ,fR_EoverP(0)
 267 ,fNcells(0)
 268 ,fNcells_EoverP(0)
 269 ,fNcells_electrons(0)
 270 ,fNcells_hadrons(0)
 271 ,fECluster_ptbins(0)
 272 ,fEoverP_ptbins(0)
 273 ,fEoverP_wSSCut(0)
 274 ,fM02_EoverP(0)
 275 ,fM20_EoverP(0)
 276 ,fTPCnsigma_eta_electrons(0)
 277 ,fTPCnsigma_eta_hadrons(0)
 278 ,fEoverP_pt_pions(0)
 279 ,ftpc_p_EoverPcut(0)
 280 ,fnsigma_p_EoverPcut(0)
 281 ,fEoverP_pt_pions2(0)
 282 ,fNcells_pt(0)
 283 ,fEoverP_pt_hadrons(0)
 284 ,fCEtaPhi_Inc(0)
 285 ,fCEtaPhi_ULS(0)
 286 ,fCEtaPhi_LS(0)
 287 ,fCEtaPhi_ULS_NoP(0)
 288 ,fCEtaPhi_LS_NoP(0)
 289 ,fCEtaPhi_ULS_Weight(0)
 290 ,fCEtaPhi_LS_Weight(0)
 291 ,fCEtaPhi_ULS_NoP_Weight(0)
 292 ,fCEtaPhi_LS_NoP_Weight(0)
 293 ,fInvMass(0)
 294 ,fInvMassBack(0)
 295 ,fDCA(0)
 296 ,fDCABack(0)
 297 ,fOpAngle(0)
 298 ,fOpAngleBack(0)
 299 ,fInvMass2(0)
 300 ,fInvMassBack2(0)
 301 ,fDCA2(0)
 302 ,fDCABack2(0)
 303 ,fOpAngle2(0)
 304 ,fOpAngleBack2(0)
 305 ,fMassCut(0.1)
 306 ,fEtaCutMin(-0.9)
 307 ,fEtaCutMax(0.9)
 308 ,fdPhiCut(0.05)
 309 ,fdEtaCut(0.05)
 310 ,fEoverPCutMin(0.8)
 311 ,fEoverPCutMax(1.2)
 312 ,fM20CutMin(0.0)
 313 ,fM20CutMax(10)
 314 ,fM02CutMin(0.0)
 315 ,fM02CutMax(10)
 316 ,fAngleCut(999)
 317 ,fChi2Cut(3.5)
 318 ,fDCAcut(999)//dca between two tracks
 319 ,fDCAcutr(999)//dca to vertex
 320 ,fDCAcutz(999)//dca to vertex
 321 ,fMassCutFlag(kTRUE)
 322 ,fAngleCutFlag(kFALSE)
 323 ,fChi2CutFlag(kFALSE)
 324 ,fDCAcutFlag(kFALSE)
 325 ,fAssHadronPtMin(0.5)
 326 ,fAssHadronPtMax(2.0)
 327 ,fPtBackgroundBeforeReco(0)
 328 ,fPtBackgroundBeforeReco2(0)
 329 ,fPtBackgroundBeforeReco_weight(0)
 330 ,fPtBackgroundBeforeReco2_weight(0)
 331 ,fpT_m_electron(0)
 332 ,fpT_gm_electron(0)
 333 ,fPtBackgroundAfterReco(0)
 334 ,fPtMinAsso(0.3)
 335 ,fTpcNclsAsso(80)
 336 ,fPtMCparticleAll(0)
 337 ,fPtMCparticleAll_nonPrimary(0)
 338 ,fPtMCparticleAlle_nonPrimary(0)
 339 ,fPtMCparticleAlle_Primary(0)
 340 ,fPtMCparticleReco(0)
 341 ,fPtMCparticleReco_nonPrimary(0)
 342 ,fPtMCparticleAllHfe1(0)
 343 ,fPtMCparticleRecoHfe1(0)
 344 ,fPtMCparticleAllHfe2(0)
 345 ,fPtMCparticleRecoHfe2(0)
 346 ,fPtMCelectronAfterAll(0)
 347 ,fPtMCelectronAfterAll_unfolding(0)
 348 ,fPtMCelectronAfterAll_nonPrimary(0)
 349 ,fPtMCelectronAfterAll_Primary(0)
 350 ,fPtMCpi0(0)
 351 ,fPtMCeta(0)
 352 ,fPtMCpi02(0)
 353 ,fPtMCeta2(0)
 354 ,fPtMCpi03(0)
 355 ,fPtMCeta3(0)
 356 ,fPtMC_EMCal_All(0)
 357 ,fPtMC_EMCal_Selected(0)
 358 ,fPtMC_TPC_All(0)
 359 ,fPtMC_TPC_Selected(0)
 360 ,fPt_track_match_den(0)
 361 ,fPt_track_match_num(0)
 362 ,fPtMCWithLabel(0)
 363 ,fPtMCWithoutLabel(0)
 364 ,fPtIsPhysicaPrimary(0)
 365 ,fCuts(0)
 366         //,reader(0)
 367 ,fCFM(0)
 368 ,fPID(new AliHFEpid("hfePid"))
 369         //Lucile
 370 ,fPIDqa(0)
 371 ,fMCstack(0)
 372 ,fRejectKinkMother(kFALSE)
 373 ,fMCtrack(0)
 374 ,fMCtrackMother(0)
 375 ,fMCtrackGMother(0)
 376 ,fMCtrackGGMother(0)
 377 ,fMCtrackGGGMother(0)
 378 ,fMCarray(0)
 379 ,fMCheader(0)
 380 ,fMCparticle(0)
 381 ,fMCparticle2(0)
 382 ,fMCparticleMother(0)
 383 ,fMCparticleGMother(0)
 384 ,fMCparticleGGMother(0)
 385 ,fMCparticleGGGMother(0)
 386 ,fEventHandler(0)
 387 ,fMCevent(0)
 388 ,fPoolMgr(0)
 389 ,fPool(0)
 390 ,fTracksClone(0)
 391 ,fTracks(0)
 392 ,fCEtaPhi_Inc_EM(0)
 393 ,fCEtaPhi_ULS_EM(0)
 394 ,fCEtaPhi_LS_EM(0)
 395 ,fCEtaPhi_ULS_Weight_EM(0)
 396 ,fCEtaPhi_LS_Weight_EM(0)
 397 ,fPoolNevents(0)
 398 ,fEventMixingFlag(0)
 399 ,fCEtaPhi_Inc_DiHadron(0)
 400 ,fPtTrigger_Inc(0)
 401         //,fEMCALRecoUtils(new AliEMCALRecoUtils)
 402         //,fEMCALGeo(0x0)
 403         //,fCaloUtils(0x0)
 404
 405 ,fBitEGA(0)
 406 //,fEMCALRecoUtils(0)//exotic
 407
 408 {
 409                 //Named constructor
 410                 // Define input and output slots here
 411                 // Input slot #0 works with a TChain
 412
 413                 //exotic
 414                 //fEMCALRecoUtils  = new AliEMCALRecoUtils();
 415
 416         DefineInput(0, TChain::Class());
 417                 // Output slot #0 id reserved by the base class for AOD
 418                 // Output slot #1 writes into a TH1 container
 419                 // DefineOutput(1, TH1I::Class());
 420         DefineOutput(1, TList::Class());
 421                 //  DefineOutput(3, TTree::Class());
 422 }
 423
 424         //________________________________________________________________________
 425 AliAnalysisTaskEMCalHFEpA::AliAnalysisTaskEMCalHFEpA()
 426 : AliAnalysisTaskSE("DefaultAnalysis_AliAnalysisTaskEMCalHFEpA")
 427 ,fCorrelationFlag(0)
 428 ,fIsMC(0)
 429 ,fUseEMCal(kFALSE)
 430 ,fUseTrigger(kFALSE)
 431 ,fUseTender(kFALSE)
 432 ,fUseShowerShapeCut(kFALSE)
 433 ,fFillBackground(kFALSE)
 434 ,fEoverPnsigma(kFALSE)
 435 ,fAssocWithSPD(kFALSE)
 436 ,fEMCEG1(kFALSE)
 437 ,fEMCEG2(kFALSE)
 438 ,fIsHFE1(kFALSE)
 439 ,fIsHFE2(kFALSE)
 440 ,fIsNonHFE(kFALSE)
 441 ,fIsFromD(kFALSE)
 442 ,fIsFromB(kFALSE)
 443 ,fIsFromPi0(kFALSE)
 444 ,fIsFromEta(kFALSE)
 445 ,fIsFromGamma(kFALSE)
 446 ,fESD(0)
 447 ,fAOD(0)
 448 ,fVevent(0)
 449 ,fPartnerCuts(new AliESDtrackCuts())
 450 ,fOutputList(0)
 451 ,fPidResponse(0)
 452 ,fNonHFE(new AliSelectNonHFE())
 453 ,fIsAOD(kFALSE)
 454 ,fCentrality(0)
 455 ,fCentralityMin(0)
 456 ,fCentralityMax(100)
 457 ,fHasCentralitySelection(kFALSE)
 458 ,fCentralityHist(0)
 459 ,fCentralityHistPass(0)
 460 ,fZvtx(0)
 461 ,fEstimator(0)
 462 ,fClus(0)
 463         //,fClusESD(0)
 464 ,fNevent(0)
 465 ,fNevent2(0)
 466 ,fPtElec_Inc(0)
 467 ,fPtPrim(0)
 468 ,fPtSec(0)
 469 ,fPtPrim2(0)
 470 ,fPtSec2(0)
 471 ,fCharge_n(0)
 472 ,fCharge_p(0)
 473 ,fTime(0)
 474 ,fTime2(0)
 475 ,ftimingEle(0)
 476 ,ftimingEle2(0)
 477 ,fPtElec_ULS(0)
 478 ,fPtElec_LS(0)
 479 ,fPtElec_ULS_NoPid(0)
 480 ,fPtElec_LS_NoPid(0)
 481 ,fPtElec_ULS_MC(0)
 482 ,fPtElec_ULS_MC_weight(0)
 483 ,fPtElec_ULS2(0)
 484 ,fPtElec_LS2(0)
 485
 486 ,fPtElec_ULS_mc_closure(0)
 487 ,fPtElec_LS_mc_closure(0)
 488 ,fPtElec_ULS2_mc_closure(0)
 489 ,fPtElec_LS2_mc_closure(0)
 490
 491 ,fPtElec_ULS_weight(0)
 492 ,fPtElec_LS_weight(0)
 493 ,fPtElec_ULS2_weight(0)
 494 ,fPtElec_LS2_weight(0)
 495 ,fTOF01(0)
 496 ,fTOF02(0)
 497 ,fTOF03(0)
 498 ,fpid(0)
 499 ,fEoverP_pt_true_electrons(0)
 500 ,fEoverP_pt_true_hadrons(0)
 501 ,fEoverP_pt_true_electrons0(0)
 502 ,fEoverP_pt_true_hadrons0(0)
 503 ,fEoverP_pt(0)
 504 ,fEoverP_tpc(0)
 505 ,fEoverP_tpc_p_trigger(0)
 506 ,fEoverP_tpc_pt_trigger(0)
 507 ,fTPC_pt(0)
 508 ,fTPC_p(0)
 509 ,fTPCnsigma_pt(0)
 510 ,fTPCnsigma_p(0)
 511 ,fTPCnsigma_p_TPC(0)
 512 ,fTPCnsigma_p_TPC_on_EMCal_acc(0)
 513 ,fTPCnsigma_p_TPC_EoverP_cut(0)
 514 ,fTPCnsigma_pt_2D(0)
 515 ,fShowerShapeCut(0)
 516 ,fShowerShapeM02_EoverP(0)
 517 ,fShowerShapeM20_EoverP(0)
 518 ,fShowerShape_ha(0)
 519 ,fShowerShape_ele(0)
 520 ,fTPCnsigma_eta(0)
 521 ,fTPCnsigma_phi(0)
 522 ,fECluster(0)
 523 ,fECluster_pure(0)
 524 ,fECluster_not_exotic(0)
 525 ,fECluster_not_exotic1(0)
 526 ,fECluster_not_exotic2(0)
 527 ,fECluster_exotic(0)
 528 ,fNCluster_pure(0)
 529 ,fNCluster_pure_aod(0)
 530 ,fNCluster_ECluster(0)
 531 ,fNcells_energy(0)
 532 ,fNcells_energy_elec_selected(0)
 533 ,fNcells_energy_not_exotic(0)
 534 ,fEtaPhi(0)
 535 ,fEtaPhi_num(0)
 536 ,fEtaPhi_den(0)
 537 ,fEtaPhi_data(0)
 538 ,fpt_reco_pt_MC_num(0)
 539 ,fpt_reco_pt_MC_den(0)
 540 ,fVtxZ(0)
 541
 542 ,fVtxZ_new1(0)
 543 ,fVtxZ_new2(0)
 544 ,fVtxZ_new3(0)
 545 ,fVtxZ_new4(0)
 546
 547 ,fzRes1(0)
 548 ,fzRes2(0)
 549 ,fSPD_track_vtx1(0)
 550 ,fSPD_track_vtx2(0)
 551
 552
 553
 554 ,fEtad(0)
 555 ,fNTracks(0)
 556 ,fTrack_Multi(0)
 557 ,fNTracks_pt(0)
 558 ,fNTracks_eta(0)
 559 ,fNTracks_phi(0)
 560 ,fNClusters(0)
 561 ,fTPCNcls_EoverP(0)
 562 ,fTPCNcls_pid(0)
 563 ,fEta(0)
 564 ,fPhi(0)
 565 ,fR(0)
 566 ,fR_EoverP(0)
 567 ,fNcells(0)
 568 ,fNcells_EoverP(0)
 569 ,fNcells_electrons(0)
 570 ,fNcells_hadrons(0)
 571 ,fECluster_ptbins(0)
 572 ,fEoverP_ptbins(0)
 573 ,fEoverP_wSSCut(0)
 574 ,fM02_EoverP(0)
 575 ,fM20_EoverP(0)
 576 ,fTPCnsigma_eta_electrons(0)
 577 ,fTPCnsigma_eta_hadrons(0)
 578 ,fEoverP_pt_pions(0)
 579 ,ftpc_p_EoverPcut(0)
 580 ,fnsigma_p_EoverPcut(0)
 581 ,fEoverP_pt_pions2(0)
 582 ,fNcells_pt(0)
 583 ,fEoverP_pt_hadrons(0)
 584 ,fCEtaPhi_Inc(0)
 585 ,fCEtaPhi_ULS(0)
 586 ,fCEtaPhi_LS(0)
 587 ,fCEtaPhi_ULS_NoP(0)
 588 ,fCEtaPhi_LS_NoP(0)
 589 ,fCEtaPhi_ULS_Weight(0)
 590 ,fCEtaPhi_LS_Weight(0)
 591 ,fCEtaPhi_ULS_NoP_Weight(0)
 592 ,fCEtaPhi_LS_NoP_Weight(0)
 593 ,fInvMass(0)
 594 ,fInvMassBack(0)
 595 ,fDCA(0)
 596 ,fDCABack(0)
 597 ,fOpAngle(0)
 598 ,fOpAngleBack(0)
 599 ,fInvMass2(0)
 600 ,fInvMassBack2(0)
 601 ,fDCA2(0)
 602 ,fDCABack2(0)
 603 ,fOpAngle2(0)
 604 ,fOpAngleBack2(0)
 605 ,fMassCut(0.1)
 606 ,fEtaCutMin(-0.9)
 607 ,fEtaCutMax(0.9)
 608 ,fdPhiCut(0.05)
 609 ,fdEtaCut(0.05)
 610 ,fEoverPCutMin(0.8)
 611 ,fEoverPCutMax(1.2)
 612 ,fM20CutMin(0)
 613 ,fM20CutMax(10)
 614 ,fM02CutMin(0)
 615 ,fM02CutMax(10)
 616 ,fAngleCut(999)
 617 ,fChi2Cut(3.5)
 618 ,fDCAcut(999)//dca between two tracks
 619 ,fDCAcutr(999)//dca to vertex
 620 ,fDCAcutz(999)//dca to vertex
 621 ,fMassCutFlag(kTRUE)
 622 ,fAngleCutFlag(kFALSE)
 623 ,fChi2CutFlag(kFALSE)
 624 ,fDCAcutFlag(kFALSE)
 625 ,fAssHadronPtMin(0.5)
 626 ,fAssHadronPtMax(2.0)
 627 ,fPtBackgroundBeforeReco(0)
 628 ,fPtBackgroundBeforeReco2(0)
 629 ,fPtBackgroundBeforeReco_weight(0)
 630 ,fPtBackgroundBeforeReco2_weight(0)
 631 ,fpT_m_electron(0)
 632 ,fpT_gm_electron(0)
 633 ,fPtBackgroundAfterReco(0)
 634 ,fPtMinAsso(0.3)
 635 ,fTpcNclsAsso(80)
 636 ,fPtMCparticleAll(0)
 637 ,fPtMCparticleAll_nonPrimary(0)
 638 ,fPtMCparticleAlle_nonPrimary(0)
 639 ,fPtMCparticleAlle_Primary(0)
 640 ,fPtMCparticleReco(0)
 641 ,fPtMCparticleReco_nonPrimary(0)
 642 ,fPtMCparticleAllHfe1(0)
 643 ,fPtMCparticleRecoHfe1(0)
 644 ,fPtMCparticleAllHfe2(0)
 645 ,fPtMCparticleRecoHfe2(0)
 646 ,fPtMCelectronAfterAll(0)
 647 ,fPtMCelectronAfterAll_unfolding(0)
 648 ,fPtMCelectronAfterAll_nonPrimary(0)
 649 ,fPtMCelectronAfterAll_Primary(0)
 650 ,fPtMCpi0(0)
 651 ,fPtMCeta(0)
 652 ,fPtMCpi02(0)
 653 ,fPtMCeta2(0)
 654 ,fPtMCpi03(0)
 655 ,fPtMCeta3(0)
 656 ,fPtMC_EMCal_All(0)
 657 ,fPtMC_EMCal_Selected(0)
 658 ,fPtMC_TPC_All(0)
 659 ,fPtMC_TPC_Selected(0)
 660 ,fPt_track_match_den(0)
 661 ,fPt_track_match_num(0)
 662 ,fPtMCWithLabel(0)
 663 ,fPtMCWithoutLabel(0)
 664 ,fPtIsPhysicaPrimary(0)
 665 ,fCuts(0)
 666         //Lucile
 667         //,reader(0)
 668 ,fCFM(0)
 669 ,fPID(new AliHFEpid("hfePid"))
 670 ,fPIDqa(0)
 671 ,fMCstack(0)
 672 ,fRejectKinkMother(kFALSE)
 673 ,fMCtrack(0)
 674 ,fMCtrackMother(0)
 675 ,fMCtrackGMother(0)
 676 ,fMCtrackGGMother(0)
 677 ,fMCtrackGGGMother(0)
 678 ,fMCarray(0)
 679 ,fMCheader(0)
 680 ,fMCparticle(0)
 681 ,fMCparticle2(0)
 682 ,fMCparticleMother(0)
 683 ,fMCparticleGMother(0)
 684 ,fMCparticleGGMother(0)
 685 ,fMCparticleGGGMother(0)
 686 ,fEventHandler(0)
 687 ,fMCevent(0)
 688 ,fPoolMgr(0)
 689 ,fPool(0)
 690 ,fTracksClone(0)
 691 ,fTracks(0)
 692 ,fCEtaPhi_Inc_EM(0)
 693 ,fCEtaPhi_ULS_EM(0)
 694 ,fCEtaPhi_LS_EM(0)
 695 ,fCEtaPhi_ULS_Weight_EM(0)
 696 ,fCEtaPhi_LS_Weight_EM(0)
 697 ,fPoolNevents(0)
 698 ,fEventMixingFlag(0)
 699 ,fCEtaPhi_Inc_DiHadron(0)
 700 ,fPtTrigger_Inc(0)
 701         //,fEMCALRecoUtils(new AliEMCALRecoUtils)
 702         //,fEMCALGeo(0x0)
 703         //,fCaloUtils(0x0)
 704 ,fBitEGA(0)
 705         //,fEMCALRecoUtils(0)//exotic
 706 {
 707                 // Constructor
 708                 // Define input and output slots here
 709                 // Input slot #0 works with a TChain
 710
 711                 //exotic
 712                 // fEMCALRecoUtils  = new AliEMCALRecoUtils();
 713
 714         DefineInput(0, TChain::Class());
 715                 // Output slot #0 id reserved by the base class for AOD
 716                 // Output slot #1 writes into a TH1 container
 717                 // DefineOutput(1, TH1I::Class());
 718         DefineOutput(1, TList::Class());
 719                 //DefineOutput(3, TTree::Class());
 720 }
 721
 722         //______________________________________________________________________
 723 AliAnalysisTaskEMCalHFEpA::~AliAnalysisTaskEMCalHFEpA()
 724 {
 725         //Destructor
 726         delete fOutputList;
 727         delete fPID;
 728         delete fCFM;
 729         delete fPIDqa;
 730         //Lucile
 731         //delete reader;
 732         //if(fEMCALRecoUtils)   delete fEMCALRecoUtils ;
 733 }
 734
 735         //______________________________________________________________________
 736         //Create Output Objects
 737         //Here we can define the histograms and others output files
 738         //Called once
 739 void AliAnalysisTaskEMCalHFEpA::UserCreateOutputObjects()
 740 {
 741                 //______________________________________________________________________
 742                 //Initialize PID
 743         if(!fPID->GetNumberOfPIDdetectors())
 744     {
 745                 fPID->AddDetector("TPC", 0);
 746     }
 747
 748         fPID->SortDetectors();
 749
 750         fPIDqa = new AliHFEpidQAmanager();
 751         fPIDqa->Initialize(fPID);
 752                 //______________________________________________________________________
 753
 754                 //______________________________________________________________________
 755                 //Initialize correction Framework and Cuts
 756         fCFM = new AliCFManager;
 757         const Int_t kNcutSteps = AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kNcutStepsRecTrack + AliHFEcuts::kNcutStepsDETrack;
 758         fCFM->SetNStepParticle(kNcutSteps);
 759         for(Int_t istep = 0; istep < kNcutSteps; istep++) fCFM->SetParticleCutsList(istep, NULL);
 760
 761         if(!fCuts)
 762         {
 763                 AliWarning("Cuts not available. Default cuts will be used");
 764                 fCuts = new AliHFEcuts;
 765                 fCuts->CreateStandardCuts();
 766         }
 767
 768         fCuts->Initialize(fCFM);
 769                 //______________________________________________________________________
 770
 771                 ///___________________//Lucile
 772
 773                 //if(fIsAOD) {
 774                 //      reader = new AliCaloTrackAODReader();
 775                 //}
 776                         //___________________________________________________
 777                 ///Output Tlist
 778                 //Create TList
 779         fOutputList = new TList();
 780         fOutputList->SetOwner();
 781
 782                 //PIDqa
 783         fOutputList->Add(fPIDqa->MakeList("PIDQA"));
 784
 785                 //Store the number of events
 786                 //Define the histo
 787         fNevent = new TH1F("fNevent","Number of Events",30,0,30);
 788         fNevent2 = new TH1F("fNevent2","Number of Events 2",30,0,30);
 789                 //And then, add to the output list
 790         fOutputList->Add(fNevent);
 791         fOutputList->Add(fNevent2);
 792
 793         fpid = new TH1F("fpid","PID flag",5,0,5);
 794         fOutputList->Add(fpid);
 795
 796         //pt Distribution
 797         fPtElec_Inc = new TH1F("fPtElec_Inc","Inclusive Electrons; p_{T} (GeV/c); Count",300,0,30);
 798
 799         fPtPrim = new TH1F("fPtPrim","Primary Electrons aod track; p_{T} (GeV/c); Count",300,0,30);
 800         fPtSec = new TH1F("fPtSec","Secundary Electrons aod track; p_{T} (GeV/c); Count",300,0,30);
 801         fPtPrim2 = new TH1F("fPtPrim2","Primary Electrons vtrack; p_{T} (GeV/c); Count",300,0,30);
 802         fPtSec2 = new TH1F("fPtSec2","Secundary Electrons vtrack; p_{T} (GeV/c); Count",300,0,30);
 803
 804
 805         fPtElec_ULS = new TH1F("fPtElec_ULS","ULS; p_{T} (GeV/c); Count",300,0,30);
 806         fPtElec_LS = new TH1F("fPtElec_LS","LS; p_{T} (GeV/c); Count",300,0,30);
 807
 808         fPtElec_ULS_NoPid = new TH1F("fPtElec_ULS_NoPid","ULS; p_{T} (GeV/c); Count",300,0,30);
 809         fPtElec_LS_NoPid = new TH1F("fPtElec_LS_NoPid","LS; p_{T} (GeV/c); Count",300,0,30);
 810
 811         fPtElec_ULS_MC = new TH1F("fPtElec_ULS_MC","ULS; p_{T} (GeV/c); Count",300,0,30);
 812         fPtElec_ULS_MC_weight = new TH1F("fPtElec_ULS_MC_weight","ULS; p_{T} (GeV/c); Count",300,0,30);
 813
 814         fPtElec_ULS_weight = new TH1F("fPtElec_ULS_weight","ULS; p_{T} (GeV/c); Count",300,0,30);
 815         fPtElec_LS_weight = new TH1F("fPtElec_LS_weight","LS; p_{T} (GeV/c); Count",300,0,30);
 816
 817         fTOF01 = new TH2F("fTOF01","",200,-20,20,200,-20,20);
 818         fTOF02 = new TH2F("fTOF02","",200,-20,20,200,-20,20);
 819         fTOF03 = new TH2F("fTOF03","",200,-20,20,200,-20,20);
 820
 821         if(fFillBackground){
 822                 fPtElec_ULS2 = new TH1F("fPtElec_ULS2","ULS; p_{T} (GeV/c); Count",300,0,30);
 823                 fPtElec_LS2 = new TH1F("fPtElec_LS2","LS; p_{T} (GeV/c); Count",300,0,30);
 824
 825                 fPtElec_ULS2_weight = new TH1F("fPtElec_ULS2_weight","ULS; p_{T} (GeV/c); Count",300,0,30);
 826                 fPtElec_LS2_weight = new TH1F("fPtElec_LS2_weight","LS; p_{T} (GeV/c); Count",300,0,30);
 827
 828                         //mc closure
 829                 fPtElec_ULS_mc_closure = new TH1F("fPtElec_ULS_mc_closure","ULS; p_{T} (GeV/c); Count",300,0,30);
 830                 fPtElec_LS_mc_closure = new TH1F("fPtElec_LS_mc_closure","LS; p_{T} (GeV/c); Count",300,0,30);
 831                 fPtElec_ULS2_mc_closure = new TH1F("fPtElec_ULS2_mc_closure","ULS; p_{T} (GeV/c); Count",300,0,30);
 832                 fPtElec_LS2_mc_closure = new TH1F("fPtElec_LS2_mc_closure","LS; p_{T} (GeV/c); Count",300,0,30);
 833
 834
 835
 836         }
 837
 838         fPtTrigger_Inc = new TH1F("fPtTrigger_Inc","pT dist for Hadron Contamination; p_{t} (GeV/c); Count",300,0,30);
 839         fTPCnsigma_pt_2D = new TH2F("fTPCnsigma_pt_2D",";pt (GeV/c);TPC Electron N#sigma",1000,0.3,30,1000,-15,10);
 840
 841         //new histos for TPC signal -> Can be used for any p range
 842         fTPCnsigma_p_TPC = new TH2F("fTPCnsigma_p_TPC",";p (GeV/c);TPC Electron N#sigma",3000,0,30,1000,-15,10);
 843         fTPCnsigma_p_TPC_on_EMCal_acc = new TH2F("fTPCnsigma_p_TPC_on_EMCal_acc",";p (GeV/c);TPC Electron N#sigma",3000,0,30,1000,-15,10);
 844         fTPCnsigma_p_TPC_EoverP_cut = new TH2F("fTPCnsigma_p_TPC_EoverP_cut",";p (GeV/c);TPC Electron N#sigma",3000,0,30,1000,-15,10);
 845
 846
 847         fShowerShapeCut = new TH2F("fShowerShapeCut","Shower Shape;M02;M20",500,0,1.8,500,0,1.8);
 848         fEtaPhi_num=new TH2F("fEtaPhi_num","#eta x #phi track;#phi;#eta",200,0.,5,50,-1.,1.);
 849         fEtaPhi_den=new TH2F("fEtaPhi_den","#eta x #phi track;#phi;#eta",200,0.,5,50,-1.,1.);
 850         fEtaPhi_data=new TH2F("fEtaPhi_data","#eta x #phi track;#phi;#eta",200,0.,5,50,-1.,1.);
 851
 852         fpt_reco_pt_MC_num=new TH2F("fpt_reco_pt_MC_num","pt reco x pt MC;pt reco; pt MC",300,0.,30,300,0.,30);
 853         fpt_reco_pt_MC_den=new TH2F("fpt_reco_pt_MC_den","pt reco x pt MC;pt reco; pt MC",300,0.,30,300,0.,30);
 854
 855
 856         fCharge_n = new TH1F("fCharge_n","Inclusive Electrons (Negative Charge); p_{t} (GeV/c); Count",200,0,30);
 857         fCharge_p = new TH1F("fCharge_p","Inclusive Positrons (Positive Charge); p_{t} (GeV/c); Count",200,0,30);
 858
 859         fECluster_pure= new TH1F("fECluster_pure", ";ECluster pure",2000,0,100);
 860         fECluster_not_exotic= new TH1F("fECluster_not_exotic", ";ECluster not exotic - function ",2000,0,100);
 861
 862         fECluster_not_exotic1= new TH1F("fECluster_not_exotic1", ";ECluster not exotic Ncells > E/3+1",2000,0,100);
 863
 864         fECluster_not_exotic2= new TH1F("fECluster_not_exotic2", ";ECluster not exotic 2",2000,0,100);
 865         fECluster_exotic= new TH1F("fECluster_exotic", ";ECluster exotic",2000,0,100);
 866
 867         //not associated with tracks
 868         fNCluster_pure= new TH1F("fNCluster_pure", ";Number of clusters - pure",2000,-1,1999);
 869         fNCluster_pure_aod= new TH1F("fNCluster_pure_aod", ";Number of clusters - pure -aod",2000,-1,1999);
 870         fNCluster_ECluster= new TH2F("fNCluster_ECluster", ";Number of clusters vs. Energy of Cluster",2000,-1,1999, 4000, -1, 1999);
 871
 872         fNcells_energy= new TH2F("fNcells_energy", "all clusters;Number of cells;Energy of Cluster",100,0,100, 2000, 0, 100);
 873         fNcells_energy_elec_selected= new TH2F("fNcells_energy_elec_selected", "clusters for electrons on TPC;Number of cells;Energy of Cluster",100,0,100, 2000, 0, 100);
 874         fNcells_energy_not_exotic= new TH2F("fNcells_energy_not_exotic", "not exotic cluster;Number of cells;Energy of Cluster ",100,0,100, 2000, 0, 100);
 875
 876         if(fUseEMCal){
 877
 878                 if(!fIsAOD){
 879                         fTime = new TH2D("fTime","Cells Cluster Time; p_{T} (GeV/c); Time (s)",300,0,30,1000,1e-8,1e-5);
 880                         fTime2 = new TH2D("fTime2","Cells Cluster Time;  p_{T} (GeV/c); Time (s)",300,0,30,1000,1e-8,1e-5);
 881                 }
 882
 883                 ftimingEle = new TH2D("ftimingEle","Cluster Time;  p_{T} (GeV/c); Time (s)",300,0,30,1000,1e-8,1e-5);
 884                 ftimingEle2 = new TH2D("ftimingEle2","Cluster Time;  p_{T} (GeV/c); Time (s)",300,0,30,1000,1e-8,1e-5);
 885
 886                 fShowerShape_ha = new TH2F("fShowerShape_ha","Shower Shape hadrons;M02;M20",500,0,1.8,500,0,1.8);
 887                 fShowerShape_ele = new TH2F("fShowerShape_ele","Shower Shape electrons;M02;M20",500,0,1.8,500,0,1.8);
 888
 889                 fShowerShapeM02_EoverP = new TH2F("fShowerShapeM02_EoverP","Shower Shape;M02;E/p",500,0,1.8,500,0,1.8);
 890                 fShowerShapeM20_EoverP = new TH2F("fShowerShapeM20_EoverP","Shower Shape;M20;E/p",500,0,1.8,500,0,1.8);
 891
 892         }
 893
 894         fOutputList->Add(fTOF01);
 895         fOutputList->Add(fTOF02);
 896         fOutputList->Add(fTOF03);
 897
 898         fOutputList->Add(fEtaPhi_num);
 899         fOutputList->Add(fEtaPhi_den);
 900         fOutputList->Add(fEtaPhi_data);
 901
 902         fOutputList->Add(fpt_reco_pt_MC_num);
 903         fOutputList->Add(fpt_reco_pt_MC_den);
 904
 905
 906         fOutputList->Add(fPtElec_Inc);
 907         fOutputList->Add(fPtElec_ULS);
 908         fOutputList->Add(fPtElec_LS);
 909         fOutputList->Add(fPtElec_ULS_NoPid);
 910         fOutputList->Add(fPtElec_LS_NoPid);
 911         fOutputList->Add(fPtElec_ULS_MC);
 912         fOutputList->Add(fPtElec_ULS_MC_weight);
 913
 914         fOutputList->Add(fPtPrim);
 915         fOutputList->Add(fPtSec);
 916         fOutputList->Add(fPtPrim2);
 917         fOutputList->Add(fPtSec2);
 918
 919
 920
 921         fOutputList->Add(fPtElec_ULS_weight);
 922         fOutputList->Add(fPtElec_LS_weight);
 923
 924         if(fFillBackground){
 925                 fOutputList->Add(fPtElec_ULS2);
 926                 fOutputList->Add(fPtElec_LS2);
 927                 fOutputList->Add(fPtElec_ULS2_weight);
 928                 fOutputList->Add(fPtElec_LS2_weight);
 929
 930                 fOutputList->Add(fPtElec_ULS_mc_closure);
 931                 fOutputList->Add(fPtElec_LS_mc_closure);
 932                 fOutputList->Add(fPtElec_ULS2_mc_closure);
 933                 fOutputList->Add(fPtElec_LS2_mc_closure);
 934         }
 935
 936
 937         fOutputList->Add(fPtTrigger_Inc);
 938         fOutputList->Add(fTPCnsigma_pt_2D);
 939
 940         fOutputList->Add(fTPCnsigma_p_TPC);
 941         fOutputList->Add(fTPCnsigma_p_TPC_on_EMCal_acc);
 942         fOutputList->Add(fTPCnsigma_p_TPC_EoverP_cut);
 943
 944
 945
 946         fOutputList->Add(fShowerShapeCut);
 947
 948         fOutputList->Add(fCharge_n);
 949         fOutputList->Add(fCharge_p);
 950
 951         fOutputList->Add(fECluster_pure);
 952         fOutputList->Add(fECluster_not_exotic);
 953         fOutputList->Add(fECluster_not_exotic1);
 954         fOutputList->Add(fECluster_not_exotic2);
 955         fOutputList->Add(fECluster_exotic);
 956
 957         fOutputList->Add(fNCluster_pure);
 958         fOutputList->Add(fNCluster_pure_aod);
 959
 960         fOutputList->Add(fNCluster_ECluster);
 961         fOutputList->Add(fNcells_energy);
 962         fOutputList->Add(fNcells_energy_elec_selected);
 963         fOutputList->Add(fNcells_energy_not_exotic);
 964
 965
 966         if(fUseEMCal){
 967
 968                 if(!fIsAOD){
 969                         fOutputList->Add(fTime);
 970                         fOutputList->Add(fTime2);
 971
 972                 }
 973
 974                 fOutputList->Add(ftimingEle);
 975                 fOutputList->Add(ftimingEle2);
 976
 977                 fOutputList->Add(fShowerShape_ha);
 978                 fOutputList->Add(fShowerShape_ele);
 979
 980                 fOutputList->Add(fShowerShapeM02_EoverP);
 981                 fOutputList->Add(fShowerShapeM20_EoverP);
 982
 983
 984
 985         }
 986
 987         fVtxZ_new1= new  TH1F("fVtxZ_new1","fVtxZ_new1",4000, -50,50);
 988         fVtxZ_new2= new  TH1F("fVtxZ_new2","fVtxZ_new2",4000, -50,50);
 989         fVtxZ_new3= new  TH1F("fVtxZ_new3","fVtxZ_new3",4000, -50,50);
 990         fVtxZ_new4= new  TH1F("fVtxZ_new4","fVtxZ_new4",4000, -50,50);
 991
 992         fzRes1= new  TH1F("fzRes1","fzRes1",4000, 0,1);
 993         fzRes2= new  TH1F("fzRes2","fzRes2",4000, 0,1);
 994         fSPD_track_vtx1= new  TH1F("fSPD_track_vtx1","fSPD_track_vtx1",4000, -5,5);
 995         fSPD_track_vtx2= new  TH1F("fSPD_track_vtx2","fSPD_track_vtx2",4000, -5,5);
 996
 997
 998                 //General Histograms
 999
1000                 //Steps
1001                 //Step 1: Before Track cuts
1002                 //Step 2: Before PID
1003                 //Step 3: After PID
1004
1005         fEoverP_pt = new TH2F *[3];
1006         fTPC_p = new TH2F *[3];
1007         fTPCnsigma_p = new TH2F *[3];
1008
1009         fECluster= new TH1F *[3];
1010         fEtaPhi= new TH2F *[3];
1011         fVtxZ= new  TH1F *[3];
1012         fEtad= new  TH1F *[3];
1013         fNTracks= new  TH1F *[3];
1014
1015         fNTracks_pt= new  TH2F *[3];
1016         fNTracks_eta= new  TH2F *[3];
1017         fNTracks_phi= new  TH2F *[3];
1018
1019         fNClusters= new TH1F *[3];
1020         fTPCNcls_EoverP= new TH2F *[3];
1021         fTPCNcls_pid=new TH2F *[4];
1022
1023         fEoverP_pt_true_electrons = new TH2F("fEoverP_pt_true_electrons",";p_{T} (GeV/c);E/p ",1000,0,30,2000,0,2);
1024         fOutputList->Add(fEoverP_pt_true_electrons);
1025         fEoverP_pt_true_electrons0 = new TH2F("fEoverP_pt_true_electrons0",";p_{T} (GeV/c);E/p ",1000,0,30,2000,0,2);
1026         fOutputList->Add(fEoverP_pt_true_electrons0);
1027
1028
1029         fEoverP_pt_true_hadrons = new TH2F("fEoverP_pt_true_hadrons",";p_{T} (GeV/c);E/p ",1000,0,30,2000,0,2);
1030         fOutputList->Add(fEoverP_pt_true_hadrons);
1031         fEoverP_pt_true_hadrons0 = new TH2F("fEoverP_pt_true_hadrons0",";p_{T} (GeV/c);E/p ",1000,0,30,2000,0,2);
1032         fOutputList->Add(fEoverP_pt_true_hadrons0);
1033
1034
1035
1036
1037         for(Int_t i = 0; i < 3; i++)
1038         {
1039                 fEoverP_pt[i] = new TH2F(Form("fEoverP_pt%d",i),";p_{t} (GeV/c);E / p ",1000,0,30,2000,0,2);
1040                 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);
1041                 fTPCnsigma_p[i] = new TH2F(Form("fTPCnsigma_p%d",i),";p (GeV/c);TPC Electron N#sigma",1000,0.3,15,1000,-15,10);
1042
1043
1044                 fECluster[i]= new TH1F(Form("fECluster%d",i), ";ECluster",2000, 0,100);
1045                 fEtaPhi[i]= new TH2F(Form("fEtaPhi%d",i),"#eta x #phi Clusters;#phi;#eta",200,0.,5,50,-1.,1.);
1046                 fVtxZ[i]= new  TH1F(Form("fVtxZ%d",i),"VtxZ",1000, -50,50);
1047                 fEtad[i]= new  TH1F(Form("fEtad%d",i),"Eta distribution",200, -1.2,1.2);
1048                 fNTracks[i]= new  TH1F(Form("fNTracks%d",i),"NTracks",1000, 0,5000);
1049
1050                 fNTracks_pt[i]= new  TH2F(Form("fNTracks_pt%d",i),"NTracks vs. pt",1000, 0,5000, 1000, 0, 100);
1051                 fNTracks_eta[i]= new  TH2F(Form("fNTracks_eta%d",i),"NTracks vs. pt",1000, 0,5000, 500, -1.0, 1.0);
1052                 fNTracks_phi[i]= new  TH2F(Form("fNTracks_phi%d",i),"NTracks vs. pt",1000, 0,5000, 500, 0, 5.0);
1053
1054
1055
1056                 fNClusters[i]= new TH1F(Form("fNClusters%d",i),"fNClusters0",200, 0,100);
1057                 fTPCNcls_EoverP[i]= new TH2F(Form("fTPCNcls_EoverP%d",i),"TPCNcls_EoverP",1000,0,200,200,0,2);
1058
1059
1060                         fOutputList->Add(fEoverP_pt[i]);
1061                 fOutputList->Add(fTPC_p[i]);
1062                 fOutputList->Add(fTPCnsigma_p[i]);
1063
1064
1065                 fOutputList->Add(fECluster[i]);
1066                 fOutputList->Add(fEtaPhi[i]);
1067                 fOutputList->Add(fVtxZ[i]);
1068                 fOutputList->Add(fEtad[i]);
1069                 fOutputList->Add(fNTracks[i]);
1070
1071                 fOutputList->Add(fNTracks_pt[i]);
1072                 fOutputList->Add(fNTracks_eta[i]);
1073                 fOutputList->Add(fNTracks_phi[i]);
1074
1075                 fOutputList->Add(fNClusters[i]);
1076                 fOutputList->Add(fTPCNcls_EoverP[i]);
1077         }
1078
1079
1080
1081
1082         fTrack_Multi= new  TH1F("fTrack_Multi","fTrack_Multi",1000, 0,1000);
1083
1084         for(Int_t i = 0; i < 4; i++)
1085         {
1086                 fTPCNcls_pid[i]= new TH2F(Form("fTPCNcls_pid%d",i),"fTPCNcls_pid;NCls;NCls for PID",200,0,200,200,0,200);
1087                 fOutputList->Add(fTPCNcls_pid[i]);
1088         }
1089
1090                 //pt bin
1091         Int_t fPtBin[7] = {1,2,4,6,8,10,15};
1092
1093
1094         Int_t fPtBin_trigger[11] = {1,2,4,6,8,10,12,14,16,18,20};
1095         fEoverP_tpc_p_trigger = new TH2F *[10];
1096         fEoverP_tpc_pt_trigger = new TH2F *[10];
1097
1098
1099         fEoverP_tpc = new TH2F *[6];
1100         fTPC_pt = new TH1F *[6];
1101         fTPCnsigma_pt = new TH1F *[6];
1102
1103         fEta=new TH1F *[6];
1104         fPhi=new TH1F *[6];
1105         fR=new TH1F *[6];
1106         fR_EoverP=new TH2F *[6];
1107         fNcells=new TH1F *[6];
1108         fNcells_EoverP=new TH2F *[6];
1109         fM02_EoverP= new TH2F *[6];
1110         fM20_EoverP= new TH2F *[6];
1111         fEoverP_ptbins=new TH1F *[6];
1112         fECluster_ptbins=new TH1F *[6];
1113         fEoverP_wSSCut=new TH1F *[6];
1114         fNcells_electrons=new TH1F *[6];
1115         fNcells_hadrons=new TH1F *[6];
1116         fTPCnsigma_eta_electrons=new TH2F *[6];
1117         fTPCnsigma_eta_hadrons=new TH2F *[6];
1118
1119         if(fCorrelationFlag)
1120         {
1121                 fCEtaPhi_Inc = new TH2F *[6];
1122                 fCEtaPhi_Inc_DiHadron = new TH2F *[6];
1123
1124                 fCEtaPhi_ULS = new TH2F *[6];
1125                 fCEtaPhi_LS = new TH2F *[6];
1126                 fCEtaPhi_ULS_NoP = new TH2F *[6];
1127                 fCEtaPhi_LS_NoP = new TH2F *[6];
1128
1129                 fCEtaPhi_ULS_Weight = new TH2F *[6];
1130                 fCEtaPhi_LS_Weight = new TH2F *[6];
1131                 fCEtaPhi_ULS_NoP_Weight = new TH2F *[6];
1132                 fCEtaPhi_LS_NoP_Weight = new TH2F *[6];
1133
1134                 fCEtaPhi_Inc_EM = new TH2F *[6];
1135
1136                 fCEtaPhi_ULS_EM = new TH2F *[6];
1137                 fCEtaPhi_LS_EM = new TH2F *[6];
1138
1139                 fCEtaPhi_ULS_Weight_EM = new TH2F *[6];
1140                 fCEtaPhi_LS_Weight_EM = new TH2F *[6];
1141
1142                 fInvMass = new TH1F("fInvMass","",200,0,0.3);
1143                 fInvMassBack = new TH1F("fInvMassBack","",200,0,0.3);
1144                 fDCA = new TH1F("fDCA","",200,0,1);
1145                 fDCABack = new TH1F("fDCABack","",200,0,1);
1146                 fOpAngle = new TH1F("fOpAngle","",200,0,0.5);
1147                 fOpAngleBack = new TH1F("fOpAngleBack","",200,0,0.5);
1148
1149                 fOutputList->Add(fInvMass);
1150                 fOutputList->Add(fInvMassBack);
1151                 fOutputList->Add(fDCA);
1152                 fOutputList->Add(fDCABack);
1153                 fOutputList->Add(fOpAngle);
1154                 fOutputList->Add(fOpAngleBack);
1155         }
1156
1157         if(fFillBackground){
1158
1159                 fInvMass = new TH1F("fInvMass","",200,0,0.3);
1160                 fInvMassBack = new TH1F("fInvMassBack","",200,0,0.3);
1161                 fDCA = new TH1F("fDCA","",200,0,1);
1162                 fDCABack = new TH1F("fDCABack","",200,0,1);
1163                 fOpAngle = new TH1F("fOpAngle","",200,0,0.5);
1164                 fOpAngleBack = new TH1F("fOpAngleBack","",200,0,0.5);
1165
1166                 fOutputList->Add(fInvMass);
1167                 fOutputList->Add(fInvMassBack);
1168                 fOutputList->Add(fDCA);
1169                 fOutputList->Add(fDCABack);
1170                 fOutputList->Add(fOpAngle);
1171                 fOutputList->Add(fOpAngleBack);
1172
1173                         //histos for TPC-only
1174                 fInvMass2 = new TH1F("fInvMass2","",200,0,0.3);
1175                 fInvMassBack2 = new TH1F("fInvMassBack2","",200,0,0.3);
1176                 fDCA2 = new TH1F("fDCA2","",200,0,1);
1177                 fDCABack2 = new TH1F("fDCABack2","",200,0,1);
1178                 fOpAngle2 = new TH1F("fOpAngle2","",200,0,0.5);
1179                 fOpAngleBack2 = new TH1F("fOpAngleBack2","",200,0,0.5);
1180
1181                 fOutputList->Add(fInvMass2);
1182                 fOutputList->Add(fInvMassBack2);
1183                 fOutputList->Add(fDCA2);
1184                 fOutputList->Add(fDCABack2);
1185                 fOutputList->Add(fOpAngle2);
1186                 fOutputList->Add(fOpAngleBack2);
1187
1188         }
1189
1190                 //new histo for trigger data
1191
1192         for(Int_t i = 0; i < 10; i++)
1193         {
1194                         fEoverP_tpc_pt_trigger[i] = new TH2F(Form("fEoverP_tpc_pt_trigger%d",i),Form("%d < p_{t} < %d GeV/c;TPC Electron N#sigma; E/p ",fPtBin_trigger[i],fPtBin_trigger[i+1]),1000,-15,15,100,0,2);
1195                         fEoverP_tpc_p_trigger[i] = new TH2F(Form("fEoverP_tpc_p_trigger%d",i),Form("%d < p < %d GeV/c;TPC Electron N#sigma; E/p ",fPtBin_trigger[i],fPtBin_trigger[i+1]),1000,-15,15,100,0,2);
1196                         fOutputList->Add(fEoverP_tpc_pt_trigger[i]);
1197                         fOutputList->Add(fEoverP_tpc_p_trigger[i]);
1198
1199         }
1200
1201
1202         for(Int_t i = 0; i < 6; i++)
1203         {
1204                 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);
1205                 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);
1206                 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);
1207
1208                 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);
1209                 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);
1210                 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);
1211                 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);
1212                 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);
1213                 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);
1214                 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);
1215                 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);
1216                 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);
1217                 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);
1218                 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);
1219                 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);
1220                 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);
1221                 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);
1222                 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);
1223
1224                 fOutputList->Add(fEoverP_tpc[i]);
1225                 fOutputList->Add(fTPC_pt[i]);
1226                 fOutputList->Add(fTPCnsigma_pt[i]);
1227
1228                 fOutputList->Add(fEta[i]);
1229                 fOutputList->Add(fPhi[i]);
1230                 fOutputList->Add(fR[i]);
1231                 fOutputList->Add(fR_EoverP[i]);
1232                 fOutputList->Add(fNcells[i]);
1233                 fOutputList->Add(fNcells_electrons[i]);
1234                 fOutputList->Add(fNcells_hadrons[i]);
1235                 fOutputList->Add(fNcells_EoverP[i]);
1236                 fOutputList->Add(fECluster_ptbins[i]);
1237                 fOutputList->Add(fEoverP_ptbins[i]);
1238                 fOutputList->Add(fEoverP_wSSCut[i]);
1239                 fOutputList->Add(fM02_EoverP[i]);
1240                 fOutputList->Add(fM20_EoverP[i]);
1241                 fOutputList->Add(fTPCnsigma_eta_electrons[i]);
1242                 fOutputList->Add(fTPCnsigma_eta_hadrons[i]);
1243
1244
1245                 if(fCorrelationFlag)
1246                 {
1247                         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);
1248                         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);
1249
1250                         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);
1251                         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);
1252                         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);
1253                         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);
1254
1255                         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);
1256                         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);
1257                         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);
1258                         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);
1259
1260                         fOutputList->Add(fCEtaPhi_Inc[i]);
1261                         fOutputList->Add(fCEtaPhi_Inc_DiHadron[i]);
1262
1263                         fOutputList->Add(fCEtaPhi_ULS[i]);
1264                         fOutputList->Add(fCEtaPhi_LS[i]);
1265                         fOutputList->Add(fCEtaPhi_ULS_NoP[i]);
1266                         fOutputList->Add(fCEtaPhi_LS_NoP[i]);
1267
1268                         fOutputList->Add(fCEtaPhi_ULS_Weight[i]);
1269                         fOutputList->Add(fCEtaPhi_LS_Weight[i]);
1270                         fOutputList->Add(fCEtaPhi_ULS_NoP_Weight[i]);
1271                         fOutputList->Add(fCEtaPhi_LS_NoP_Weight[i]);
1272
1273                         if(fEventMixingFlag)
1274                         {
1275                                 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);
1276
1277                                 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);
1278                                 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);
1279
1280                                 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);
1281                                 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);
1282
1283                                 fOutputList->Add(fCEtaPhi_Inc_EM[i]);
1284
1285                                 fOutputList->Add(fCEtaPhi_ULS_EM[i]);
1286                                 fOutputList->Add(fCEtaPhi_LS_EM[i]);
1287
1288                                 fOutputList->Add(fCEtaPhi_ULS_Weight_EM[i]);
1289                                 fOutputList->Add(fCEtaPhi_LS_Weight_EM[i]);
1290                         }
1291                 }
1292         }
1293
1294                 //pt integrated
1295         fTPCnsigma_eta = new TH2F("fTPCnsigma_eta",";Pseudorapidity #eta; TPC signal - <TPC signal>_{elec} (#sigma)",200,-0.9,0.9,200,-15,15);
1296         fTPCnsigma_phi = new TH2F("fTPCnsigma_phi",";Azimuthal Angle #phi; TPC signal - <TPC signal>_{elec} (#sigma)",200,0,2*TMath::Pi(),200,-15,15);
1297
1298
1299
1300         fNcells_pt=new TH2F("fNcells_pt","fNcells_pt",1000, 0,20,100,0,30);
1301         fEoverP_pt_pions= new TH2F("fEoverP_pt_pions","fEoverP_pt_pions",1000,0,30,2000,0,2);
1302
1303         ftpc_p_EoverPcut= new TH2F("ftpc_p_EoverPcut","ftpc_p_EoverPcut",1000,0,30,200,20,200);
1304         fnsigma_p_EoverPcut= new TH2F("fnsigma_p_EoverPcut","fnsigma_p_EoverPcut",1000,0,30,500,-15,15);
1305
1306         fEoverP_pt_pions2= new TH2F("fEoverP_pt_pions2","fEoverP_pt_pions2",1000,0,30,2000,0,2);
1307         fEoverP_pt_hadrons= new TH2F("fEoverP_pt_hadrons","fEoverP_pt_hadrons",1000,0,30,2000,0,2);
1308
1309
1310         fOutputList->Add(fTPCnsigma_eta);
1311         fOutputList->Add(fTPCnsigma_phi);
1312
1313         fOutputList->Add(fNcells_pt);
1314         fOutputList->Add(fEoverP_pt_pions);
1315
1316         fOutputList->Add(ftpc_p_EoverPcut);
1317         fOutputList->Add(fnsigma_p_EoverPcut);
1318
1319         fOutputList->Add(fEoverP_pt_pions2);
1320         fOutputList->Add(fEoverP_pt_hadrons);
1321
1322         fOutputList->Add(fVtxZ_new1);
1323         fOutputList->Add(fVtxZ_new2);
1324         fOutputList->Add(fVtxZ_new3);
1325         fOutputList->Add(fVtxZ_new4);
1326
1327         fOutputList->Add(fzRes1);
1328         fOutputList->Add(fzRes2);
1329         fOutputList->Add(fSPD_track_vtx1);
1330         fOutputList->Add(fSPD_track_vtx2);
1331
1332
1333
1334                 //__________________________________________________________________
1335                 //Efficiency studies
1336         if(fIsMC)
1337         {
1338                 fPtBackgroundBeforeReco = new TH1F("fPtBackgroundBeforeReco",";p_{T} (GeV/c);Count",300,0,30);
1339                 fPtBackgroundBeforeReco_weight = new TH1F("fPtBackgroundBeforeReco_weight",";p_{T} (GeV/c);Count",300,0,30);
1340                 if(fFillBackground)fPtBackgroundBeforeReco2 = new TH1F("fPtBackgroundBeforeReco2",";p_{T} (GeV/c);Count",300,0,30);
1341                 if(fFillBackground)fPtBackgroundBeforeReco2_weight = new TH1F("fPtBackgroundBeforeReco2_weight",";p_{T} (GeV/c);Count",300,0,30);
1342                 fpT_m_electron= new TH2F("fpT_m_electron","fpT_m_electron",300,0,30,300,0,30);
1343                 fpT_gm_electron= new TH2F("fpT_gm_electron","fpT_gm_electron",300,0,30,300,0,30);
1344
1345                 fPtBackgroundAfterReco = new TH1F("fPtBackgroundAfterReco",";p_{T} (GeV/c);Count",300,0,30);
1346                 fPtMCparticleAll = new TH1F("fPtMCparticleAll",";p_{T} (GeV/c);Count",200,0,40);
1347                 fPtMCparticleReco = new TH1F("fPtMCparticleReco",";p_{T} (GeV/c);Count",200,0,40);
1348
1349                 fPtMCparticleAll_nonPrimary = new TH1F("fPtMCparticleAll_nonPrimary",";p_{T} (GeV/c);Count",200,0,40);
1350                 fPtMCparticleAlle_nonPrimary = new TH1F("fPtMCparticleAlle_nonPrimary",";p_{T} (GeV/c);Count",200,0,40);
1351                 fPtMCparticleAlle_Primary = new TH1F("fPtMCparticleAlle_Primary",";p_{T} (GeV/c);Count",200,0,40);
1352
1353                 fPtMCparticleReco_nonPrimary = new TH1F("fPtMCparticleReco_nonPrimary",";p_{T} (GeV/c);Count",200,0,40);
1354
1355                 fPtMCparticleAllHfe1 = new TH1F("fPtMCparticleAllHfe1",";p_{t} (GeV/c);Count",200,0,40);
1356                 fPtMCparticleRecoHfe1 = new TH1F("fPtMCparticleRecoHfe1",";p_{t} (GeV/c);Count",200,0,40);
1357                 fPtMCparticleAllHfe2 = new TH1F("fPtMCparticleAllHfe2",";p_{t} (GeV/c);Count",200,0,40);
1358                 fPtMCparticleRecoHfe2 = new TH1F("fPtMCparticleRecoHfe2",";p_{t} (GeV/c);Count",200,0,40);
1359
1360                 fPtMCelectronAfterAll = new TH1F("fPtMCelectronAfterAll",";p_{T} (GeV/c);Count",200,0,40);
1361                 fPtMCelectronAfterAll_unfolding = new TH1F("fPtMCelectronAfterAll_unfolding",";p_{T} (GeV/c);Count",200,0,40);
1362                 fPtMCelectronAfterAll_nonPrimary = new TH1F("fPtMCelectronAfterAll_nonPrimary",";p_{T} (GeV/c);Count",200,0,40);
1363                 fPtMCelectronAfterAll_Primary = new TH1F("fPtMCelectronAfterAll_Primary",";p_{T} (GeV/c);Count",200,0,40);
1364
1365
1366
1367                 fPtMCpi0 = new TH1F("fPtMCpi0",";p_{t} (GeV/c);Count",200,0,30);
1368                 fPtMCeta = new TH1F("fPtMCeta",";p_{T} (GeV/c);Count",200,0,30);
1369                 fPtMCpi02 = new TH1F("fPtMCpi02",";p_{t} (GeV/c);Count",200,0,30);
1370                 fPtMCeta2 = new TH1F("fPtMCeta2",";p_{T} (GeV/c);Count",200,0,30);
1371                 fPtMCpi03 = new TH1F("fPtMCpi03",";p_{t} (GeV/c);Count",200,0,30);
1372                 fPtMCeta3 = new TH1F("fPtMCeta3",";p_{T} (GeV/c);Count",200,0,30);
1373
1374                 fPtMC_EMCal_All= new TH1F("fPtMC_EMCal_All",";p_{t} (GeV/c);Count",200,0,40);
1375                 fPtMC_EMCal_Selected= new TH1F("fPtMC_EMCal_Selected",";p_{t} (GeV/c);Count",200,0,40);
1376                 fPtMC_TPC_All= new TH1F("fPtMC_TPC_All",";p_{T} (GeV/c);Count",200,0,40);
1377                 fPtMC_TPC_Selected = new TH1F("fPtMC_TPC_Selected",";p_{T} (GeV/c);Count",200,0,40);
1378
1379                 fPt_track_match_den = new TH1F("fPt_track_match_den",";p_{T} (GeV/c);Count",200,0,40);
1380                 fPt_track_match_num = new TH1F("fPt_track_match_num",";p_{T} (GeV/c);Count",200,0,40);
1381
1382                 fPtMCWithLabel = new TH1F("fPtMCWithLabel",";p_{t} (GeV/c);Count",200,0,40);
1383                 fPtMCWithoutLabel = new TH1F("fPtMCWithoutLabel",";p_{t} (GeV/c);Count",200,0,40);
1384                 fPtIsPhysicaPrimary = new TH1F("fPtIsPhysicaPrimary",";p_{t} (GeV/c);Count",200,0,40);
1385
1386                 fOutputList->Add(fPtBackgroundBeforeReco);
1387                 fOutputList->Add(fPtBackgroundBeforeReco_weight);
1388
1389                 if(fFillBackground) fOutputList->Add(fPtBackgroundBeforeReco2);
1390                 if(fFillBackground) fOutputList->Add(fPtBackgroundBeforeReco2_weight);
1391
1392                 fOutputList->Add(fpT_m_electron);
1393                 fOutputList->Add(fpT_gm_electron);
1394
1395                 fOutputList->Add(fPtBackgroundAfterReco);
1396                 fOutputList->Add(fPtMCparticleAll);
1397                 fOutputList->Add(fPtMCparticleReco);
1398
1399                 fOutputList->Add(fPtMCparticleAll_nonPrimary);
1400                 fOutputList->Add(fPtMCparticleAlle_nonPrimary);
1401
1402                 fOutputList->Add(fPtMCparticleAlle_Primary);
1403                 fOutputList->Add(fPtMCparticleReco_nonPrimary);
1404
1405                 fOutputList->Add(fPtMCparticleAllHfe1);
1406                 fOutputList->Add(fPtMCparticleRecoHfe1);
1407                 fOutputList->Add(fPtMCparticleAllHfe2);
1408                 fOutputList->Add(fPtMCparticleRecoHfe2);
1409                 fOutputList->Add(fPtMCelectronAfterAll);
1410                 fOutputList->Add(fPtMCelectronAfterAll_unfolding);
1411
1412                 fOutputList->Add(fPtMCelectronAfterAll_nonPrimary);
1413                 fOutputList->Add(fPtMCelectronAfterAll_Primary);
1414
1415
1416
1417                 fOutputList->Add(fPtMCpi0);
1418                 fOutputList->Add(fPtMCeta);
1419                 fOutputList->Add(fPtMCpi02);
1420                 fOutputList->Add(fPtMCeta2);
1421                 fOutputList->Add(fPtMCpi03);
1422                 fOutputList->Add(fPtMCeta3);
1423                 fOutputList->Add(fPtMC_EMCal_All);
1424                 fOutputList->Add(fPtMC_EMCal_Selected);
1425                 fOutputList->Add(fPtMC_TPC_All);
1426                 fOutputList->Add(fPtMC_TPC_Selected);
1427
1428                 fOutputList->Add(fPt_track_match_den);
1429                 fOutputList->Add(fPt_track_match_num);
1430
1431                 fOutputList->Add(fPtMCWithLabel);
1432                 fOutputList->Add(fPtMCWithoutLabel);
1433                 fOutputList->Add(fPtIsPhysicaPrimary);
1434         }
1435
1436         fCentralityHist = new TH1F("fCentralityHist",";Centrality (%); Count",1000000,0,100);
1437         fCentralityHistPass = new TH1F("fCentralityHistPass",";Centrality (%); Count",1000000,0,100);
1438         fOutputList->Add(fCentralityHist);
1439         fOutputList->Add(fCentralityHistPass);
1440
1441         //______________________________________________________________________
1442         //Mixed event analysis
1443         if(fEventMixingFlag)
1444         {
1445                 fPoolNevents = new TH1F("fPoolNevents","Event Mixing Statistics; Number of events; Count",1000,0,1000);
1446                 fOutputList->Add(fPoolNevents);
1447
1448                 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.
1449                 Int_t poolsize   = 1000;  // Maximum number of events, ignored in the present implemented of AliEventPoolManager
1450
1451                 Int_t nCentralityBins  = 15;
1452                 Double_t centralityBins[] = { 0, 1, 2, 3, 4, 5, 10, 20, 30, 40, 50, 60, 70, 80, 90, 100.1 };
1453
1454                 Int_t nZvtxBins  = 9;
1455                 Double_t vertexBins[] = {-10, -7, -5, -3, -1, 1, 3, 5, 7, 10};
1456
1457                 fPoolMgr = new AliEventPoolManager(poolsize, trackDepth, nCentralityBins, (Double_t*) centralityBins, nZvtxBins, (Double_t*) vertexBins);
1458         }
1459                 //______________________________________________________________________
1460
1461         PostData(1, fOutputList);
1462
1463                 ///______________________________________________________________________
1464 }
1465
1466 //______________________________________________________________________
1467 //Main loop
1468 //Called for each event
1469 void AliAnalysisTaskEMCalHFEpA::UserExec(Option_t *)
1470 {
1471                 //Check Event
1472         fESD = dynamic_cast<AliESDEvent*>(InputEvent());
1473         fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
1474
1475         if(!(fESD || fAOD))
1476         {
1477                 printf("ERROR: fESD & fAOD not available\n");
1478                 return;
1479         }
1480
1481         fVevent = dynamic_cast<AliVEvent*>(InputEvent());
1482
1483         if(!fVevent)
1484         {
1485                 printf("ERROR: fVEvent not available\n");
1486                 return;
1487         }
1488
1489                 //Check Cuts
1490         if(!fCuts)
1491         {
1492                 AliError("HFE cuts not available");
1493                 return;
1494         }
1495                 //Check PID
1496         if(!fPID->IsInitialized())
1497         {
1498                         // Initialize PID with the given run number
1499                 AliWarning("PID not initialised, get from Run no");
1500
1501                 if(fIsAOD)
1502                 {
1503                         fPID->InitializePID(fAOD->GetRunNumber());
1504                 }
1505                 else
1506                 {
1507                         fPID->InitializePID(fESD->GetRunNumber());
1508                 }
1509         }
1510
1511                 //PID response
1512         fPidResponse = fInputHandler->GetPIDResponse();
1513
1514
1515                 //Check PID response
1516         if(!fPidResponse)
1517         {
1518                 AliDebug(1, "Using default PID Response");
1519                 fPidResponse = AliHFEtools::GetDefaultPID(kFALSE, fInputEvent->IsA() == AliAODEvent::Class());
1520         }
1521
1522         fPID->SetPIDResponse(fPidResponse);
1523
1524         fCFM->SetRecEventInfo(fVevent);
1525
1526         Double_t *fListOfmotherkink = 0;
1527         Int_t fNumberOfVertices = 0;
1528         Int_t fNumberOfMotherkink = 0;
1529
1530
1531         //total event before event selection
1532         fNevent->Fill(1);
1533
1534         //______________________________________________________________________
1535         //Vertex Selection
1536         if(fIsAOD)
1537         {
1538                 const AliAODVertex* trkVtx = fAOD->GetPrimaryVertex();
1539                 if(!trkVtx || trkVtx->GetNContributors()<=0) return;
1540                 TString vtxTtl = trkVtx->GetTitle();
1541                 if(!vtxTtl.Contains("VertexerTracks")) return;
1542                         //Float_t zvtx = trkVtx->GetZ();
1543                 Float_t zvtx = -100;
1544                 zvtx=trkVtx->GetZ();
1545                 fZvtx = zvtx;
1546
1547                 fVtxZ_new1->Fill(fZvtx);
1548
1549                 const AliAODVertex* spdVtx = fAOD->GetPrimaryVertexSPD();
1550                 if(spdVtx->GetNContributors()<=0) return;
1551                 TString vtxTyp = spdVtx->GetTitle();
1552                 Double_t cov[6]={0};
1553                 spdVtx->GetCovarianceMatrix(cov);
1554                 Double_t zRes = TMath::Sqrt(cov[5]);
1555
1556                 fzRes1->Fill(zRes);
1557                 if(vtxTyp.Contains("vertexer:Z") && (zRes>0.25)) return;
1558                 fzRes2->Fill(zRes);
1559
1560                 fSPD_track_vtx1->Fill(spdVtx->GetZ() - trkVtx->GetZ());
1561                 if(TMath::Abs(spdVtx->GetZ() - trkVtx->GetZ())>0.5) return;
1562                 fSPD_track_vtx2->Fill(spdVtx->GetZ() - trkVtx->GetZ());
1563
1564
1565                 if(TMath::Abs(zvtx) > 10) return;
1566                 fVtxZ_new2->Fill(fZvtx);
1567
1568                 if(fabs(zvtx>10.0))return;
1569                 fVtxZ_new3->Fill(fZvtx);
1570
1571
1572                 //Look for kink mother for AOD
1573
1574                 fNumberOfVertices = 0;
1575                 fNumberOfMotherkink = 0;
1576
1577                 if(fIsAOD)
1578                 {
1579                         fNumberOfVertices = fAOD->GetNumberOfVertices();
1580
1581                         fListOfmotherkink = new Double_t[fNumberOfVertices];
1582
1583                         for(Int_t ivertex=0; ivertex < fNumberOfVertices; ivertex++)
1584                         {
1585                                 AliAODVertex *aodvertex = fAOD->GetVertex(ivertex);
1586                                 if(!aodvertex) continue;
1587                                 if(aodvertex->GetType()==AliAODVertex::kKink)
1588                                 {
1589                                         AliAODTrack *mother1 = (AliAODTrack *) aodvertex->GetParent();
1590                                         if(!mother1) continue;
1591                                         Int_t idmother = mother1->GetID();
1592                                         fListOfmotherkink[fNumberOfMotherkink] = idmother;
1593                                         fNumberOfMotherkink++;
1594                                 }
1595                         }
1596                 }
1597         }
1598         else
1599         {
1600
1601
1602
1603                 /// ESD
1604                 const AliESDVertex* trkVtx = fESD->GetPrimaryVertex();
1605                 if(!trkVtx || trkVtx->GetNContributors()<=0) return;
1606                 TString vtxTtl = trkVtx->GetTitle();
1607                 if(!vtxTtl.Contains("VertexerTracks")) return;
1608                 Float_t zvtx = -100;
1609                 zvtx=trkVtx->GetZ();
1610
1611
1612                 const AliESDVertex* spdVtx = fESD->GetPrimaryVertexSPD();
1613                 if(spdVtx->GetNContributors()<=0) return;
1614                 TString vtxTyp = spdVtx->GetTitle();
1615                 Double_t cov[6]={0};
1616                 spdVtx->GetCovarianceMatrix(cov);
1617                 Double_t zRes = TMath::Sqrt(cov[5]);
1618                 if(vtxTyp.Contains("vertexer:Z") && (zRes>0.25)) return;
1619                 if(TMath::Abs(spdVtx->GetZ() - trkVtx->GetZ())>0.5) return;
1620                 if(TMath::Abs(zvtx) > 10) return;
1621         }
1622
1623         //______________________________________________________________________
1624         //after vertex selection
1625         fNevent->Fill(10);
1626
1627         //______________________________________________________________________
1628         //EMCal Trigger Selection (Threshold selection)
1629
1630         TString firedTrigger;
1631         TString TriggerEG1("EG1"); //takes trigger with name with EG1, ex: CEMC7EG1-B-NOPF-CENTNOTRD
1632         TString TriggerEG2("EG2");
1633         //Jan 17, 2014
1634         //TString TriggerEJE("EJE");
1635
1636         if(fAOD) firedTrigger = fAOD->GetFiredTriggerClasses();
1637         else if(fESD) firedTrigger = fESD->GetFiredTriggerClasses();
1638
1639                 //Bool_t IsEventEMCALL0=kTRUE;
1640         Bool_t IsEventEMCALL1=kFALSE;
1641
1642         if(firedTrigger.Contains(TriggerEG1)){
1643                 fNevent->Fill(2);
1644                 IsEventEMCALL1=kTRUE;
1645         }
1646         if(firedTrigger.Contains(TriggerEG2)){
1647                 fNevent->Fill(3);
1648                 IsEventEMCALL1=kTRUE;
1649         }
1650
1651         //if the flag is for a given threshold and it was not fired, return.
1652
1653         if(fEMCEG1){
1654                 if(!firedTrigger.Contains(TriggerEG1))return;
1655                 if(firedTrigger.Contains(TriggerEG2)){
1656                         fNevent->Fill(4);
1657
1658                 }
1659
1660         }
1661
1662
1663         if(fEMCEG2){
1664                 if(!firedTrigger.Contains(TriggerEG2))return;
1665                 if(firedTrigger.Contains(TriggerEG1)){
1666                         fNevent->Fill(5);
1667                 }
1668
1669         }
1670
1671
1672
1673         //______________________________________________________________________
1674         //Testing if there is an overlap EGA and EJE
1675         //none
1676         /*
1677         if(!(firedTrigger.Contains(TriggerEG1) && firedTrigger.Contains(TriggerEG2) ) && !firedTrigger.Contains(TriggerEJE))
1678         {
1679                 fNevent->Fill(6);
1680         }
1681                 //only GA
1682         if((firedTrigger.Contains(TriggerEG1) || firedTrigger.Contains(TriggerEG2)) && !firedTrigger.Contains(TriggerEJE))
1683         {
1684                 fNevent->Fill(7);
1685         }
1686                 //only JE
1687         if(!(firedTrigger.Contains(TriggerEG1) && firedTrigger.Contains(TriggerEG2)) && firedTrigger.Contains(TriggerEJE))
1688         {
1689                 fNevent->Fill(8);
1690         }
1691                 //both
1692         if((firedTrigger.Contains(TriggerEG1) || firedTrigger.Contains(TriggerEG2)) && firedTrigger.Contains(TriggerEJE))
1693         {
1694                 fNevent->Fill(9);
1695         }
1696         */
1697
1698
1699
1700
1701         //Only events with at least 2 tracks are accepted
1702         Int_t fNOtrks =  fVevent->GetNumberOfTracks();
1703         if(fNOtrks<2) return;
1704
1705         fNevent->Fill(11);
1706
1707         if(fIsAOD){
1708                 Int_t fNOtrks2 =  fAOD->GetNumberOfTracks();
1709                 if(fNOtrks2<2) return;
1710         }
1711         fNevent->Fill(12);
1712
1713         //______________________________________________________________________
1714         //new track loop to select events
1715         //track pt cut (at least 2)
1716
1717         if(fUseTrigger){
1718                 if(fIsAOD){
1719                         double fTrackMulti=0;
1720                         for(Int_t iTracks = 0; iTracks < fVevent->GetNumberOfTracks(); iTracks++)
1721                         {
1722                                 AliVParticle* Vtrack = fVevent->GetTrack(iTracks);
1723                                 if (!Vtrack) continue;
1724
1725                                 AliVTrack *track = dynamic_cast<AliVTrack*>(Vtrack);
1726                                         //AliAODTrack *atrack = dynamic_cast<AliAODTrack*>(Vtrack);
1727
1728                                 if((track->Pt())<0.2 || (track->Pt())>1000.0) continue;
1729                                 else fTrackMulti=fTrackMulti+1;
1730
1731                         }
1732                                 //Only take event if track multiplicity is bigger than 2.
1733                         if(fTrackMulti<2) return;
1734                 }
1735         }
1736         fNevent->Fill(13);
1737         //______________________________________________________________________
1738         //Using more cuts than I have beeing using
1739         //eta cut and primary (at least 2)
1740         if(fUseTrigger){
1741                 if(fIsAOD){
1742                         double fTrackMulti2=0;
1743                         for(Int_t i = 0; i < fVevent->GetNumberOfTracks(); i++)
1744                         {
1745                                 AliVParticle* Vtrack2 = fVevent->GetTrack(i);
1746                                 if (!Vtrack2) continue;
1747
1748
1749                                 AliVTrack *track_new = dynamic_cast<AliVTrack*>(Vtrack2);
1750                                 AliAODTrack *aodtrack = dynamic_cast<AliAODTrack*>(Vtrack2);
1751
1752
1753                                 if(aodtrack)
1754                                 {
1755
1756
1757                                     if(TMath::Abs(track_new->Eta())> 0.9) continue;
1758                                         if (aodtrack->GetType()!= AliAODTrack::kPrimary) continue ;
1759                                     else fTrackMulti2=fTrackMulti2+1;
1760                                 }
1761                         }
1762                                 //Only take event if track multiplicity is bigger than 2.
1763                         if(fTrackMulti2<2) return;
1764
1765
1766                 }
1767         }
1768         fNevent->Fill(14);
1769 //______________________________________________________________________
1770 //Using more cuts than I have beeing using
1771 //hybrid (at least2)
1772         if(fUseTrigger){
1773                 if(fIsAOD){
1774                         double fTrackMulti3=0;
1775                         for(Int_t i = 0; i < fVevent->GetNumberOfTracks(); i++)
1776                         {
1777                                 AliVParticle* Vtrack3 = fVevent->GetTrack(i);
1778                                 if (!Vtrack3) continue;
1779
1780                                         //AliVTrack *track_new = dynamic_cast<AliVTrack*>(Vtrack3);
1781                                 AliAODTrack *aodtrack = dynamic_cast<AliAODTrack*>(Vtrack3);
1782
1783
1784                                 if(aodtrack)
1785                                 {
1786
1787                                         if (!aodtrack->IsHybridGlobalConstrainedGlobal()) continue ;
1788                                                 //another option if I don't want to use hybrid
1789                                                 //if ( aodtrack->TestFilterBit(128)==kFALSE) continue ;
1790                                     else fTrackMulti3=fTrackMulti3+1;
1791                                 }
1792                         }
1793                         //Only take event if track multiplicity is bigger than 2.
1794                         if(fTrackMulti3<2) return;
1795
1796                 }
1797         }
1798         fNevent->Fill(15);
1799 //______________________________________________________________________
1800
1801
1802         if(fUseTrigger){
1803                 if(fIsAOD){
1804                         double fTrackMulti4=0;
1805                         for(Int_t iTracks = 0; iTracks < fVevent->GetNumberOfTracks(); iTracks++)
1806                         {
1807                                 AliVParticle* Vtrack4 = fVevent->GetTrack(iTracks);
1808                                 if (!Vtrack4) continue;
1809
1810
1811                                         //AliVTrack *track = dynamic_cast<AliVTrack*>(Vtrack4);
1812                                 AliAODTrack *atrack = dynamic_cast<AliAODTrack*>(Vtrack4);
1813
1814                                 if(!atrack->TestFilterBit(768)) continue;
1815                                 if(!atrack->IsHybridGlobalConstrainedGlobal()) continue ;
1816
1817
1818                                 else fTrackMulti4=fTrackMulti4+1;
1819
1820                         }
1821                         //Only take event if track multiplicity is bigger than 2.
1822                         if(fTrackMulti4<2) return;
1823                         fTrack_Multi->Fill(fTrackMulti4);
1824                 }
1825         }
1826         fNevent->Fill(16);
1827 //______________________________________________________________________
1828
1829
1830 //______________________________________________________________________
1831 //Centrality Selection
1832         if(fHasCentralitySelection)
1833         {
1834                 Float_t centrality = -1;
1835
1836                 if(fIsAOD)
1837                 {
1838                         fCentrality = ((AliVAODHeader*)fAOD->GetHeader())->GetCentralityP();
1839                 }
1840                 else
1841                 {
1842                         fCentrality = fESD->GetCentrality();
1843                 }
1844
1845                 if(fEstimator==1) centrality = fCentrality->GetCentralityPercentile("ZDC");
1846                 else centrality = fCentrality->GetCentralityPercentile("V0A");
1847
1848                         //cout << "Centrality = " << centrality << " %" << endl;
1849
1850                 fCentralityHist->Fill(centrality);
1851
1852                 if(centrality<fCentralityMin || centrality>fCentralityMax) return;
1853
1854                 fCentralityHistPass->Fill(centrality);
1855         }
1856         //______________________________________________________________________
1857
1858
1859         fNevent->Fill(17);
1860
1861         //______________________________________________________________________
1862
1863         if(fIsMC)
1864         {
1865                 if(fIsAOD)
1866                 {
1867                         fMCarray = dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
1868
1869                         if(!fMCarray)
1870                         {
1871                                 AliError("Array of MC particles not found");
1872                                 return;
1873                         }
1874
1875                         fMCheader = dynamic_cast<AliAODMCHeader*>(fAOD->GetList()->FindObject(AliAODMCHeader::StdBranchName()));
1876
1877                         if(!fMCheader)
1878                         {
1879                                 AliError("Could not find MC Header in AOD");
1880                                 return;
1881                         }
1882
1883                         for(Int_t iMC = 0; iMC < fMCarray->GetEntries(); iMC++)
1884                         {
1885                                 fMCparticle = (AliAODMCParticle*) fMCarray->At(iMC);
1886                                 if(fMCparticle->GetMother()>0) fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
1887
1888                                 Int_t pdg = fMCparticle->GetPdgCode();
1889
1890                                 //====================================================================
1891                                 //trying take pions spectra 27/May/2014
1892                                 //IsPrimary only take events from pythia
1893                                 //IsPhysicalPrimariee: (all prompt particles, including strong decay products plus weak decay product from heavy-flavor).
1894                                 //eta cut same as MinJung
1895
1896                                 if(fMCparticle->Eta()>=-0.8 && fMCparticle->Eta()<=0.8)
1897                                 {
1898                                         if(fMCparticle->IsPrimary()){
1899
1900                                                 if(TMath::Abs(pdg)==111) fPtMCpi0->Fill(fMCparticle->Pt());
1901                                                 if(TMath::Abs(pdg)==221) fPtMCeta->Fill(fMCparticle->Pt());
1902                                                 //eta cut same as MinJung
1903                                         }
1904
1905                                         if(fMCparticle->IsPhysicalPrimary()){
1906                                                 if(TMath::Abs(pdg)==111) fPtMCpi02->Fill(fMCparticle->Pt());
1907                                                 if(TMath::Abs(pdg)==221) fPtMCeta2->Fill(fMCparticle->Pt());
1908
1909                                         }
1910
1911                                         if(TMath::Abs(pdg)==111) fPtMCpi03->Fill(fMCparticle->Pt());
1912                                         if(TMath::Abs(pdg)==221) fPtMCeta3->Fill(fMCparticle->Pt());
1913                                 }
1914                                 //====================================================================
1915
1916                                 double proX = fMCparticle->Xv();
1917                                 double proY = fMCparticle->Yv();
1918                                 double proR = sqrt(pow(proX,2)+pow(proY,2));
1919
1920
1921                                 if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax && fMCparticle->Charge()!=0)
1922                                 {
1923                                                 //to correct background
1924                                         if (TMath::Abs(pdg) == 11 && fMCparticle->GetMother()>0){
1925                                                 Int_t mpdg = fMCparticleMother->GetPdgCode();
1926
1927                                                 if(TMath::Abs(mpdg) == 221 || TMath::Abs(mpdg) == 22 || TMath::Abs(mpdg) == 111){
1928
1929                                                         if(proR<7){
1930                                                                 fPtMCparticleAlle_nonPrimary->Fill(fMCparticle->Pt()); //denominator for total efficiency for all electrons, and not primary
1931
1932                                                         }
1933                                                 }
1934                                         }
1935
1936                                         if (TMath::Abs(pdg) == 11 && fMCparticle->IsPhysicalPrimary()){
1937                                                 fPtMCparticleAlle_Primary->Fill(fMCparticle->Pt()); //denominator for total efficiency for all electrons primary
1938                                         }
1939
1940                                         if( TMath::Abs(pdg) == 211 || TMath::Abs(pdg) == 2212 || TMath::Abs(pdg) == 321 || TMath::Abs(pdg) == 11 || TMath::Abs(pdg) == 13 )
1941                                         {
1942
1943                                                 fPtMCparticleAll_nonPrimary->Fill(fMCparticle->Pt()); //denominator for total efficiency for all particles, and not primary
1944                                                 if(fMCparticle->IsPhysicalPrimary())
1945                                                 {
1946                                                         fPtMCparticleAll->Fill(fMCparticle->Pt());
1947
1948                                                         Bool_t MotherFound = FindMother(iMC);
1949                                                         //Bool_t MotherFound = FindMother(track->GetLabel());
1950                                                         if(MotherFound)
1951                                                         {
1952                                                                 if(fIsHFE1){
1953                                                                         //denominator for total efficiency and tracking
1954                                                                         //unfolding: denominator is pt_MC and numerator is pt_reco
1955                                                                         fPtMCparticleAllHfe1->Fill(fMCparticle->Pt());
1956                                                                         fEtaPhi_den->Fill(fMCparticle->Phi(),fMCparticle->Eta());
1957
1958
1959                                                                 } //denominator for total efficiency and tracking
1960                                                                 if(fIsHFE2){
1961                                                                         fPtMCparticleAllHfe2->Fill(fMCparticle->Pt());
1962                                                                 }
1963                                                         }
1964                                                 }
1965                                         }
1966                                 }//eta cut
1967
1968
1969
1970                         }//loop tracks
1971
1972
1973
1974                         //second loop over track, but only the primaries ones
1975                         //only primary pions --> how to take the primaries ones in AOD?
1976                         /*
1977                         for(Int_t iMC = 0; iMC < fMCarray->GetNPrimary(); iMC++){
1978                                 fMCparticle = (AliAODMCParticle*) fMCarray->At(iMC);
1979                                 pdg = fMCparticle->GetPdgCode();
1980
1981                                 if(TMath::Abs(pdg)==111) fPtMCpi0->Fill(fMCparticle->Pt());
1982                                 if(TMath::Abs(pdg)==221) fPtMCeta->Fill(fMCparticle->Pt());
1983
1984                                 if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax)
1985                                 {
1986
1987                                         if(TMath::Abs(pdg)==111) fPtMCpi02->Fill(fMCparticle->Pt());
1988                                         if(TMath::Abs(pdg)==221) fPtMCeta2->Fill(fMCparticle->Pt());
1989
1990                                 }
1991                         }
1992                          */
1993
1994
1995                 }//AOD
1996                 else
1997                 {
1998                         fEventHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
1999                         if (!fEventHandler) {
2000                                 Printf("ERROR: Could not retrieve MC event handler");
2001                                 return;
2002                         }
2003
2004                         fMCevent = fEventHandler->MCEvent();
2005                         if (!fMCevent) {
2006                                 Printf("ERROR: Could not retrieve MC event");
2007                                 return;
2008                         }
2009
2010                         fMCstack = fMCevent->Stack();
2011
2012                         //pion and eta spectrum
2013                         //MinJung code
2014
2015                         //----------------------------------------------------------------------------------------------------
2016                         AliVParticle *mctrack2 = NULL;
2017                         AliMCParticle *mctrack0 = NULL;
2018
2019
2020                         for(Int_t imc = 0; imc <fMCEvent->GetNumberOfPrimaries(); imc++){
2021                                 if(!(mctrack2 = fMCEvent->GetTrack(imc))) continue;
2022                                 TParticle* mcpart0 = fMCEvent->Stack()->Particle(imc);
2023                                 if(!mcpart0) continue;
2024                                 mctrack0 = dynamic_cast<AliMCParticle *>(mctrack2);
2025                                 if(!mctrack0) continue;
2026
2027                                 if(TMath::Abs(AliHFEtools::GetRapidity(mcpart0))<0.8){
2028
2029                                         if(TMath::Abs(mctrack0->PdgCode()) == 111) // pi0
2030                                         {
2031                                                 fPtMCpi0->Fill(mctrack0->Pt());
2032                                         }
2033
2034                                         if(TMath::Abs(mctrack0->PdgCode()) == 221) // eta
2035                                         {
2036                                                 fPtMCeta->Fill(mctrack0->Pt());
2037                                         }
2038
2039                                 }
2040
2041                         }
2042                         // end of MinJung
2043                         //----------------------------------------------------------------------------------------------------
2044
2045
2046                 for(Int_t iMC = 0; iMC < fMCstack->GetNtrack(); iMC++)
2047                 {
2048
2049                                 fMCtrack = fMCstack->Particle(iMC);
2050                                 if(fMCtrack->GetFirstMother()>0) fMCtrackMother = fMCstack->Particle(fMCtrack->GetFirstMother());
2051                                 TParticle *particle=fMCstack->Particle(iMC);
2052
2053                                 Int_t pdg = fMCtrack->GetPdgCode();
2054
2055
2056                                 if(TMath::Abs(pdg)==111) fPtMCpi0->Fill(fMCtrack->Pt());
2057                                 if(TMath::Abs(pdg)==221) fPtMCeta->Fill(fMCtrack->Pt());
2058
2059
2060                                 if(fMCtrack->Eta()>=fEtaCutMin && fMCtrack->Eta()<=fEtaCutMax)
2061                                 {
2062
2063                                                 //to correct background
2064                                         if (TMath::Abs(pdg) == 11 && fMCtrack->GetFirstMother()>0){
2065                                                 Int_t mpdg = fMCtrackMother->GetPdgCode();
2066                                                 if(TMath::Abs(mpdg) == 221 || TMath::Abs(mpdg) == 22 || TMath::Abs(mpdg) == 111){
2067                                                                 Double_t proR=particle->R();
2068                                                                 if(proR<7){
2069                                                                         fPtMCparticleAlle_nonPrimary->Fill(fMCtrack->Pt()); //denominator for total efficiency for all electrons, and not primary
2070                                                                 }
2071                                                 }
2072                                         }
2073
2074                                         if(TMath::Abs(pdg) == 11 && fMCstack->IsPhysicalPrimary(iMC)){
2075
2076                                                 fPtMCparticleAlle_Primary->Fill(fMCtrack->Pt());
2077                                         }
2078
2079                                         if( TMath::Abs(pdg) == 211 || TMath::Abs(pdg) == 2212 || TMath::Abs(pdg) == 321 || TMath::Abs(pdg) == 11 || TMath::Abs(pdg) == 13 )
2080                                         {
2081                                                 fPtMCparticleAll_nonPrimary->Fill(fMCtrack->Pt());//denominator for total efficiency for all particle, non Primary track
2082
2083                                                 if(fMCstack->IsPhysicalPrimary(iMC))
2084                                                 {
2085                                                         fPtMCparticleAll->Fill(fMCtrack->Pt());
2086
2087                                                         Bool_t MotherFound = FindMother(iMC);
2088                                                         //Bool_t MotherFound = FindMother(track->GetLabel());
2089                                                         if(MotherFound)
2090                                                         {
2091                                                                 if(fIsHFE1){
2092                                                                         fPtMCparticleAllHfe1->Fill(fMCtrack->Pt());//denominator for total efficiency and tracking
2093                                                                         fEtaPhi_den->Fill(fMCtrack->Phi(),fMCtrack->Eta());
2094                                                                 }
2095                                                                 if(fIsHFE2){
2096                                                                         fPtMCparticleAllHfe2->Fill(fMCtrack->Pt());
2097                                                                 }
2098                                                         }
2099                                                 }//Is Physical primary
2100                                         }
2101                                 }//eta cut
2102                 }//loop tracks
2103                 }//ESD
2104         }//Is MC
2105
2106 //______________________________________________________________________
2107 //threshold selection was here
2108 //______________________________________________________________________
2109 //all events selected
2110
2111         fNevent->Fill(0);
2112
2113
2114 //______________________________________________________________________
2115 //events in the threshold
2116
2117         if(firedTrigger.Contains(TriggerEG1))
2118         {
2119                 if(fEMCEG1){
2120                         fNevent->Fill(18);
2121                     if(!firedTrigger.Contains(TriggerEG2)) fNevent->Fill(19);
2122                                 //if(firedTrigger.Contains(TriggerEG2)) return;
2123                 }
2124         }
2125
2126
2127         //EG2
2128         if(firedTrigger.Contains(TriggerEG2))
2129         {
2130                 if(fEMCEG2){
2131                         fNevent->Fill(20);
2132                         if(!firedTrigger.Contains(TriggerEG1)) fNevent->Fill(21);
2133                                 //if(firedTrigger.Contains(TriggerEG1)) return;
2134                 }
2135         }
2136
2137
2138
2139         //New cluster information
2140         //after trigger threshold selection
2141         Int_t ClsNo2 = -999;
2142         ClsNo2 = fVevent->GetNumberOfCaloClusters();
2143         fNCluster_pure->Fill(ClsNo2);
2144
2145
2146
2147         if(ClsNo2<=0){
2148                 fNevent->Fill(22); //events with no cluster
2149                 return;
2150         }
2151
2152         //in order to include time cut
2153         //fEMCALCells = fAOD->GetEMCALCells();
2154         //Double_t tof = clus->GetTOF();
2155         //clus->GetNCells()
2156         //if ( clus->E() < minE ) continue ;
2157
2158
2159
2160         if(fIsAOD){
2161
2162                         //AliAODHeader * aodh = dynamic_cast<AliAODHeader*>(fAOD->GetHeader());
2163           //  if(!aodh) AliFatal("Not a standard AOD");
2164                         //Int_t bc= aodh->GetBunchCrossNumber();
2165
2166
2167                 Int_t ClsNo_aod = -999;
2168                 ClsNo_aod = fAOD->GetNumberOfCaloClusters();
2169                 fNCluster_pure_aod->Fill(ClsNo_aod);
2170                         //Bool_t exotic=kTRUE;
2171
2172                 for (Int_t i=0; i< ClsNo_aod; i++ ){
2173
2174                         //fClus = fVevent->GetCaloCluster(i);
2175                         //to be compatible with Shingo
2176                     AliVCluster *clust = 0x0;
2177                         clust = (AliVCluster*) fAOD->GetCaloCluster(i);
2178
2179                         if(clust && clust->IsEMCAL())
2180                         {
2181                                 //pure cluster information
2182                                 fECluster_pure->Fill(clust->E());
2183
2184                                 fNcells_energy->Fill(clust->GetNCells(),clust->E());
2185                                 fNCluster_ECluster->Fill(ClsNo_aod,clust->E());
2186
2187                                 if(clust->E()>1000) fNevent->Fill(23);
2188
2189                                 //exotic
2190                                 /*
2191                                 exotic   = fEMCALRecoUtils->IsExoticCluster(clust, (AliVCaloCells*)fAOD->GetEMCALCells(), bc);
2192                                 if(exotic == kFALSE){
2193                                         fECluster_not_exotic->Fill(clust->E());
2194                                         fNcells_energy_not_exotic->Fill(clust->GetNCells(),clust->E());
2195                                 }
2196                                 */
2197
2198                                 //approximation to remove exotics
2199                                 if(clust->GetNCells()<5 && clust->E()>15.0){
2200                                         fECluster_exotic->Fill(clust->E());
2201                                 }
2202                                 //Marcel cut (another approximation to remove exotics)
2203                                 else if((clust->GetNCells())> ((clust->E())/3+1)){
2204                                         fECluster_not_exotic1->Fill(clust->E());
2205                                 }
2206                                 else{
2207                                         fECluster_not_exotic2->Fill(clust->E());
2208                                 }
2209
2210
2211                         }
2212                         /*
2213                         //______________________________________________________________________
2214                         //Trying to remove events with bad cells and find patches
2215                         //First, I will try to count them
2216                         //______________________________________________________________________
2217
2218                         if(clust && clust->IsEMCAL())
2219                         {
2220                                 Bool_t badchannel = ContainsBadChannel("EMCAL", clust->GetCellsAbsId(),clust->GetNCells() );
2221                                 printf("Contém bad channel? %d ", badchannel);
2222                                 if(badchannel)fNevent2->Fill(0);
2223
2224                                 //trying to find patches
2225                                 TArrayI patches_found=GetTriggerPatches(IsEventEMCALL0, IsEventEMCALL1);
2226                                 printf("N patches %d, first %d, last %d\n",patches_found.GetSize(),  patches_found.At(0), patches_found.At(patches_found.GetSize()-1));
2227
2228                         }
2229
2230                         //end of bad cells
2231                         //______________________________________________________________________
2232                         */
2233
2234                 }
2235         }
2236
2237
2238         fNevent->Fill(24); //events with cluster
2239
2240
2241         fVtxZ_new4->Fill(fZvtx);
2242
2243
2244         Int_t ClsNo = -999;
2245         if(!fIsAOD) ClsNo = fESD->GetNumberOfCaloClusters();
2246         else ClsNo = fAOD->GetNumberOfCaloClusters();
2247
2248         //______________________________________________________________________
2249
2250         ///_____________________________________________________________________
2251         ///Track loop
2252         Int_t NTracks=fVevent->GetNumberOfTracks();
2253
2254         ///////////////////////////////////////////////////////////////////////////////////////////////////////////
2255         //To use tender
2256         if(fUseTender){
2257                 TClonesArray  *fTracks_tender = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("AODFilterTracks"));
2258                 NTracks = fTracks_tender->GetEntries();
2259                 TClonesArray  *fCaloClusters_tender = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("EmcCaloClusters"));
2260                 ClsNo = fCaloClusters_tender->GetEntries();
2261
2262         }
2263         ///////////////////////////////////////////////////////////////////////////////////////////////////////////
2264
2265         for(Int_t iTracks = 0; iTracks < NTracks; iTracks++)
2266         {
2267                 AliVParticle* Vtrack = fVevent->GetTrack(iTracks);
2268                 if (!Vtrack)
2269                 {
2270                         printf("ERROR: Could not receive track %d\n", iTracks);
2271                         continue;
2272                 }
2273
2274                         //main place where track is defined
2275
2276                 AliVTrack *track = dynamic_cast<AliVTrack*>(Vtrack);
2277                 AliESDtrack *etrack = dynamic_cast<AliESDtrack*>(Vtrack);
2278                 AliAODTrack *atrack = dynamic_cast<AliAODTrack*>(Vtrack);
2279
2280                 ///////////////////////////////////////////////////////////////////////////////////////////////////////////
2281                 //To use tender
2282                 if(fUseTender){
2283                         TClonesArray  *fTracks_tender = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("AODFilterTracks"));
2284                         track = static_cast<AliAODTrack*>(fTracks_tender->At(iTracks));
2285                         if (!track) {
2286                                 printf("ERROR: Could not receive track calibrated %d\n", iTracks);
2287                                 continue;
2288                         }
2289                 }
2290                 ///////////////////////////////////////////////////////////////////////////////////////////////////////////
2291
2292                 Double_t fTPCnSigma = -999;
2293                 Double_t fTOFnSigma = -999;
2294                 Double_t fTPCnSigma_pion = -999;
2295                 Double_t fTPCnSigma_proton = -999;
2296                 Double_t fTPCnSigma_kaon = -999;
2297                 Double_t fTPCsignal = -999;
2298                 Double_t fPt = -999;
2299                 Double_t fP = -999;
2300
2301                 //December 9th 2013
2302                 //Etacut test on the begging
2303                 fEtad[0]->Fill(track->Eta());
2304                         //if(track->Eta()<fEtaCutMin || track->Eta()>fEtaCutMax) continue;
2305                 fEtad[1]->Fill(track->Eta());
2306
2307
2308
2309
2310                 ///_____________________________________________________________________________
2311                 ///Fill QA plots without track selection
2312                 fPt = track->Pt();
2313                 fP = TMath::Sqrt((track->Pt())*(track->Pt()) + (track->Pz())*(track->Pz()));
2314
2315                 fTPCsignal = track->GetTPCsignal();
2316                 fTPCnSigma = fPidResponse->NumberOfSigmasTPC(track, AliPID::kElectron);
2317                 fTOFnSigma = fPidResponse->NumberOfSigmasTOF(track, AliPID::kElectron);
2318                 fTPCnSigma_pion = fPidResponse->NumberOfSigmasTPC(track, AliPID::kPion);
2319                 fTPCnSigma_proton = fPidResponse->NumberOfSigmasTPC(track, AliPID::kProton);
2320                 fTPCnSigma_kaon = fPidResponse->NumberOfSigmasTPC(track, AliPID::kKaon);
2321
2322                 fTPC_p[0]->Fill(fPt,fTPCsignal);
2323                 fTPCnsigma_p[0]->Fill(fP,fTPCnSigma);
2324
2325
2326                 Float_t TPCNcls = track->GetTPCNcls();
2327                 //TPC Ncls for pid
2328                 Float_t TPCNcls_pid = track->GetTPCsignalN();
2329
2330
2331
2332                 Float_t pos[3]={0,0,0};
2333
2334                 Double_t fEMCflag = kFALSE;
2335                 if(track->GetEMCALcluster()>0)
2336                 {
2337                         fClus = fVevent->GetCaloCluster(track->GetEMCALcluster());
2338                         ///////////////////////////////////////////////////////////////////////////////////////////////////////////
2339                         //to use tender
2340                         if(fUseTender){
2341                                 int EMCalIndex = -1;
2342                                 EMCalIndex = track->GetEMCALcluster();
2343                                 TClonesArray  *fCaloClusters_tender = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("EmcCaloClusters"));
2344                                 fClus = static_cast<AliVCluster*>(fCaloClusters_tender->At(EMCalIndex));
2345                         }
2346             ///////////////////////////////////////////////////////////////////////////////////////////////////////////
2347
2348                         if(fClus->IsEMCAL())
2349                         {
2350
2351                                 //only for charged tracks
2352                                 fECluster[0]->Fill(fClus->E());
2353
2354
2355                                 if(TMath::Abs(fClus->GetTrackDx())<=fdPhiCut && TMath::Abs(fClus->GetTrackDz())<=fdEtaCut)
2356                                 {
2357                                         fEMCflag = kTRUE;
2358                                         fEoverP_pt[0]->Fill(fPt,(fClus->E() / fP));
2359
2360
2361                                         Float_t Energy  = fClus->E();
2362                                         Float_t EoverP  = Energy/track->P();
2363                                                 //Float_t M02   = fClus->GetM02();
2364                                                 //Float_t M20   = fClus->GetM20();
2365
2366                                                 /////////////// for Eta Phi distribution
2367                                         fClus->GetPosition(pos);
2368                                         TVector3 vpos(pos[0],pos[1],pos[2]);
2369                                         Double_t cphi = vpos.Phi();
2370                                         Double_t ceta = vpos.Eta();
2371                                         fEtaPhi[0]->Fill(cphi,ceta);
2372
2373
2374                                         fTPCNcls_EoverP[0]->Fill(TPCNcls, EoverP);
2375                                 }
2376                         }
2377                 }
2378
2379                         //______________________________________________________________
2380                         // Vertex
2381
2382                 fVtxZ[0]->Fill(fZvtx);
2383                 if(iTracks == 0)fNTracks[0]->Fill(fNOtrks);
2384                 fNTracks_pt[0]->Fill(fNOtrks, fPt);
2385                 fNTracks_eta[0]->Fill(fNOtrks, track->Eta());
2386                 fNTracks_phi[0]->Fill(fNOtrks, track->Phi());
2387
2388
2389                 fNClusters[0]->Fill(ClsNo);
2390                 fTPCNcls_pid[0]->Fill(TPCNcls, TPCNcls_pid);
2391                         //______________________________________________________________
2392
2393 ///Fill QA plots without track selection
2394 ///_____________________________________________________________________________
2395 //______________________________________________________________________________________
2396 //Track Selection Cuts
2397
2398 //AOD (Test Filter Bit)
2399                 if(fIsAOD)
2400                 {
2401                                 // standard cuts with very loose DCA - BIT(4)
2402                                 // Description:
2403                         /*
2404                          GetStandardITSTPCTrackCuts2011(kFALSE)
2405                          SetMaxChi2PerClusterTPC(4);
2406                          SetAcceptKinkDaughters(kFALSE);
2407                          SetRequireTPCRefit(kTRUE);
2408                          SetRequireITSRefit(kTRUE);
2409                          SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kAny);
2410                          SetMaxDCAToVertexZ(2);
2411                          SetMaxDCAToVertex2D(kFALSE);
2412                          SetRequireSigmaToVertex(kFALSE);
2413                          SetMaxChi2PerClusterITS(36);
2414                          SetMaxDCAToVertexXY(2.4)
2415                          SetMaxDCAToVertexZ(3.2)
2416                          SetDCaToVertex2D(kTRUE)
2417                          */
2418
2419                         if(!atrack->TestFilterMask(AliAODTrack::kTrkGlobalNoDCA)) continue;
2420                 }
2421
2422 //RecKine: ITSTPC cuts
2423                 if(!ProcessCutStep(AliHFEcuts::kStepRecKineITSTPC, track)) continue;
2424 //RecKink
2425                 if(fRejectKinkMother)
2426                 {
2427                         if(fIsAOD)
2428                         {
2429                                 Bool_t kinkmotherpass = kTRUE;
2430                                 for(Int_t kinkmother = 0; kinkmother < fNumberOfMotherkink; kinkmother++)
2431                                 {
2432                                         if(track->GetID() == fListOfmotherkink[kinkmother])
2433                                         {
2434                                                 kinkmotherpass = kFALSE;
2435                                                 continue;
2436                                         }
2437                                 }
2438                                 if(!kinkmotherpass) continue;
2439                         }
2440                         else
2441                         {
2442                                 if(etrack->GetKinkIndex(0) != 0) continue;
2443                         }
2444                 }
2445
2446 //RecPrim
2447
2448                         //it was not working on aod... testing again
2449                         //July 29th, 2014: aparently working again
2450
2451 ///if(!fIsAOD)
2452 //{
2453
2454                 if(!ProcessCutStep(AliHFEcuts::kStepRecPrim, track)) continue;
2455 //}
2456
2457 //HFEcuts: ITS layers cuts
2458                 if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsITS, track)) continue;
2459
2460 //HFE cuts: TPC PID cleanup
2461                 if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTPC, track)) continue;
2462
2463 //DCA cut done by hand -- July 29th, 2014  -> fix it... has to be in absolute values!!!!
2464                 if(fIsAOD){
2465                         Double_t d0z0[2], cov[3];
2466                         AliAODVertex *prim_vtx = fAOD->GetPrimaryVertex();
2467                         track->PropagateToDCA(prim_vtx, fAOD->GetMagneticField(), 20., d0z0, cov);
2468                         Double_t DCAxy = d0z0[0];
2469                         Double_t DCAz = d0z0[1];
2470
2471                         if(DCAxy >= fDCAcutr || DCAz>=fDCAcutz) continue;
2472                 }
2473
2474
2475
2476
2477
2478 //______________________________________________________________________________________
2479
2480                 if(fIsAOD){
2481                                 //AOD test -- Francesco suggestion
2482                                 //aod test -- Francesco suggestion
2483                         AliAODTrack *aod_track=dynamic_cast<AliAODTrack*>(fAOD->GetTrack(iTracks));
2484                         if(!aod_track) AliFatal("Not a standard AOD");
2485
2486                         Int_t type=aod_track->GetType();
2487                         if(type==AliAODTrack::kPrimary) fPtPrim->Fill(aod_track->Pt());
2488                         if(type==AliAODTrack::kFromDecayVtx) fPtSec->Fill(aod_track->Pt());
2489
2490                                 //Int_t type2=track->GetType();
2491                         if(type==AliAODTrack::kPrimary) fPtPrim2->Fill(track->Pt());
2492                         if(type==AliAODTrack::kFromDecayVtx) fPtSec2->Fill(track->Pt());
2493                 }
2494
2495
2496 ///_____________________________________________________________
2497 ///QA plots after track selection
2498                 if(fIsMC)
2499                 {
2500                         if(track->GetLabel()>=0) fPtMCWithLabel->Fill(fPt);
2501                         else fPtMCWithoutLabel->Fill(fPt);
2502                 }
2503
2504                 if(fIsMC  && track->GetLabel()>=0)
2505                 {
2506                         if(fIsAOD)
2507                         {
2508                                 fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
2509
2510
2511                                 if(fMCparticle->IsPhysicalPrimary())
2512                                 {
2513                                         fPtIsPhysicaPrimary->Fill(fPt);
2514                                 }
2515
2516                                 Int_t pdg = fMCparticle->GetPdgCode();
2517                                 if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax && fMCparticle->Charge()!=0)
2518                                 {
2519
2520
2521                                         if( TMath::Abs(pdg) == 211 || TMath::Abs(pdg) == 2212 || TMath::Abs(pdg) == 321 || TMath::Abs(pdg) == 11 || TMath::Abs(pdg) == 13 )
2522                                         {
2523                                                 fPtMCparticleReco_nonPrimary->Fill(fMCparticle->Pt()); //not Primary track
2524
2525                                                 if(fMCparticle->IsPhysicalPrimary())
2526                                                 {
2527                                                         fPtMCparticleReco->Fill(fMCparticle->Pt());
2528
2529                                                         Bool_t MotherFound = FindMother(track->GetLabel());
2530                                                         if(MotherFound)
2531                                                         {
2532                                                                 if(fIsHFE1)
2533                                                                 {
2534                                                                         fPtMCparticleRecoHfe1->Fill(fMCparticle->Pt());//numerator tracking
2535                                                                         //unfolding
2536                                                                         fpt_reco_pt_MC_den->Fill(track->Pt(),fMCparticle->Pt());
2537
2538                                                                 }
2539                                                                 if(fIsHFE2){
2540                                                                         fPtMCparticleRecoHfe2->Fill(fMCparticle->Pt());
2541                                                                 }
2542                                                         }
2543                                                 }
2544                                         }
2545                                 }
2546
2547
2548
2549                         }//close AOD
2550                          //ESD
2551                         else
2552                         {
2553
2554
2555                                 if(fMCtrack->Eta()>=fEtaCutMin && fMCtrack->Eta()<=fEtaCutMax)
2556                                 {
2557
2558                                         fMCtrack = fMCstack->Particle(track->GetLabel());
2559                                         Int_t pdg = fMCtrack->GetPdgCode();
2560
2561                                         if( TMath::Abs(pdg) == 211 || TMath::Abs(pdg) == 2212 || TMath::Abs(pdg) == 321 || TMath::Abs(pdg) == 11 || TMath::Abs(pdg) == 13 )
2562                                         {
2563                                                 fPtMCparticleReco_nonPrimary->Fill(fMCtrack->Pt());//not Primary track
2564                                         }
2565
2566
2567                                         if(fMCstack->IsPhysicalPrimary(track->GetLabel()))
2568                                         {
2569                                                 fPtIsPhysicaPrimary->Fill(fPt);
2570
2571
2572                                                 if( TMath::Abs(pdg) == 211 || TMath::Abs(pdg) == 2212 || TMath::Abs(pdg) == 321 || TMath::Abs(pdg) == 11 || TMath::Abs(pdg) == 13 )
2573                                                 {
2574                                                         fPtMCparticleReco->Fill(fMCtrack->Pt());
2575
2576                                                         Bool_t MotherFound = FindMother(track->GetLabel());
2577                                                         if(MotherFound)
2578                                                         {
2579                                                                 if(fIsHFE1) fPtMCparticleRecoHfe1->Fill(fMCtrack->Pt());//numerator tracking
2580                                                                 if(fIsHFE2) fPtMCparticleRecoHfe2->Fill(fMCtrack->Pt());
2581                                                         }
2582                                                 }
2583                                         }
2584                                 }
2585                         }//close ESD
2586                 }//close IsMC
2587
2588                 fTPC_p[1]->Fill(fPt,fTPCsignal);
2589                 fTPCnsigma_p[1]->Fill(fP,fTPCnSigma);
2590                 Double_t fPtBin[7] = {1,2,4,6,8,10,15};
2591                 Double_t fPtBin_trigger[11] = {1,2,4,6,8,10,12,14,16,18,20};
2592
2593                 TPCNcls = track->GetTPCNcls();
2594                 Float_t pos2[3]={0,0,0};
2595
2596                 if(track->GetEMCALcluster()>0)
2597                 {
2598                         fClus = fVevent->GetCaloCluster(track->GetEMCALcluster());
2599                         ///////////////////////////////////////////////////////////////////////////////////////////////////////////
2600                         //to use tender
2601                         if(fUseTender){
2602                                 int EMCalIndex = -1;
2603                                 EMCalIndex = track->GetEMCALcluster();
2604                                 TClonesArray  *fCaloClusters_tender = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("EmcCaloClusters"));
2605                                 fClus = static_cast<AliVCluster*>(fCaloClusters_tender->At(EMCalIndex));
2606                         }
2607                         ///////////////////////////////////////////////////////////////////////////////////////////////////////////
2608
2609                         if(fClus->IsEMCAL())
2610                         {
2611                                 if(TMath::Abs(fClus->GetTrackDx())<=fdPhiCut && TMath::Abs(fClus->GetTrackDz())<=fdEtaCut)
2612                                 {
2613                                         fEoverP_pt[1]->Fill(fPt,(fClus->E() / fP));
2614
2615                                         Float_t Energy  = fClus->E();
2616                                         Float_t EoverP  = Energy/track->P();
2617                                         Float_t M02     = fClus->GetM02();
2618                                         Float_t M20     = fClus->GetM20();
2619                                         Float_t ncells  = fClus->GetNCells();
2620
2621                                                 /////////////// for Eta Phi distribution
2622                                         fClus->GetPosition(pos2);
2623                                         TVector3 vpos(pos2[0],pos2[1],pos2[2]);
2624                                         Double_t cphi = vpos.Phi();
2625                                         Double_t ceta = vpos.Eta();
2626                                         fEtaPhi[1]->Fill(cphi,ceta);
2627
2628                                         fECluster[1]->Fill(Energy);
2629                                         fTPCNcls_EoverP[1]->Fill(TPCNcls, EoverP);
2630
2631
2632                                         //for EMCal trigger performance
2633                                         if(EoverP > 0.9){
2634                                                 ftpc_p_EoverPcut->Fill(track->P(), fTPCsignal);
2635                                                 fnsigma_p_EoverPcut->Fill(track->P(), fTPCnSigma);
2636
2637                                         }
2638
2639
2640                                         //for hadron contamination calculations
2641
2642
2643                                         // EtaCut -> dados
2644                                         if(track->Eta()>=fEtaCutMin && track->Eta()<=fEtaCutMax ){
2645
2646                                                 //-------------------------------------------------------------------
2647                                                 //true hadrons E/p shape before cuts
2648                                                 if(fIsMC && fIsAOD && track->GetLabel()>=0){
2649                                                         fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
2650                                                         Int_t pdg = fMCparticle->GetPdgCode();
2651
2652                                                         if( TMath::Abs(pdg) != 11){
2653                                                                 fEoverP_pt_true_hadrons0->Fill(fPt,(fClus->E() / fP));
2654                                                         }
2655                                                 }
2656
2657
2658                                                 //true electrons E/p shape before cuts
2659                                                 if(fIsMC && fIsAOD && track->GetLabel()>=0){
2660                                                         fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
2661                                                         Int_t pdg = fMCparticle->GetPdgCode();
2662
2663                                                         if( TMath::Abs(pdg) == 11){
2664                                                                 fEoverP_pt_true_electrons0->Fill(fPt,(fClus->E() / fP));
2665                                                         }
2666                                                 }
2667                                                 //-------------------------------------------------------------------
2668
2669
2670                                                 //main
2671                                                 if(TMath::Abs(fTPCnSigma_pion)<3 || TMath::Abs(fTPCnSigma_proton)<3 || TMath::Abs(fTPCnSigma_kaon)<3 ){
2672
2673                                                         if(fTPCnSigma<-3.5){
2674                                                                 fEoverP_pt_hadrons->Fill(fPt,EoverP);
2675                                                                 if(fUseEMCal) fShowerShape_ha->Fill(M02,M20);
2676
2677
2678                                                                 //after cut -> Are they really hadrons? E/p shape
2679                                                                 if(fIsMC && fIsAOD && track->GetLabel()>=0){
2680                                                                         fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
2681                                                                         Int_t pdg = fMCparticle->GetPdgCode();
2682
2683                                                                         if( TMath::Abs(pdg) != 11){
2684                                                                                 fEoverP_pt_true_hadrons->Fill(fPt,(fClus->E() / fP));
2685                                                                         }
2686                                                                 }
2687
2688
2689                                                         }
2690                                                 }
2691
2692                                                 //for systematic studies of hadron contamination
2693                                                 if(fTPCnSigma < -3){
2694                                                         fEoverP_pt_pions->Fill(fPt, EoverP);
2695
2696                                                 }
2697
2698                                                 if(fTPCnSigma < -3.5){
2699                                                         fEoverP_pt_pions2->Fill(fPt, EoverP);
2700
2701                                                 }
2702
2703
2704                                         }
2705
2706
2707
2708
2709                                         for(Int_t i = 0; i < 6; i++)
2710                                         {
2711                                                 if(fPt>=fPtBin[i] && fPt<fPtBin[i+1])
2712                                                 {
2713
2714                                                         if(fTPCnSigma < -5 && fTPCnSigma > -10){
2715                                                                 fNcells_hadrons[i]->Fill(ncells);
2716                                                         }
2717                                                                 //hadrons selection using E/p
2718                                                         if((fClus->E() / fP) >= 0.2 && (fClus->E() / fP) <= 0.7){
2719                                                                 fTPCnsigma_eta_hadrons[i]->Fill(fTPCnSigma, ceta);
2720                                                         }
2721                                                                 //electrons selection using E/p
2722                                                         if((fClus->E() / fP) >= fEoverPCutMin && (fClus->E() / fP) <= fEoverPCutMax) {
2723                                                                 fTPCnsigma_eta_electrons[i]->Fill(fTPCnSigma, ceta);
2724                                                         }
2725                                                 }
2726                                         }
2727
2728                                         if((fClus->E() / fP) >= fEoverPCutMin && (fClus->E() / fP) <= fEoverPCutMax)
2729                                         {
2730                                                 fTPCnsigma_eta->Fill(track->Eta(),fTPCnSigma);
2731                                                 fTPCnsigma_phi->Fill(track->Phi(),fTPCnSigma);
2732
2733                                                 if(fUseEMCal)
2734                                                 {
2735                                                         if(fTPCnSigma < 3.5 && fCorrelationFlag)
2736                                                         {
2737                                                                 fPtTrigger_Inc->Fill(fPt);
2738                                                                 DiHadronCorrelation(track, iTracks);
2739                                                         }
2740                                                 }
2741                                         }
2742
2743                                 }
2744                         }
2745                 }
2746
2747                 //______________________________________________________________
2748                 // Vertex
2749
2750                 fVtxZ[1]->Fill(fZvtx);
2751                 if(iTracks == 0)fNTracks[1]->Fill(fNOtrks);
2752                 fNTracks_pt[1]->Fill(fNOtrks, fPt);
2753                 fNTracks_eta[1]->Fill(fNOtrks, track->Eta());
2754                 fNTracks_phi[1]->Fill(fNOtrks, track->Phi());
2755                 fNClusters[1]->Fill(ClsNo);
2756                 fTPCNcls_pid[1]->Fill(TPCNcls, TPCNcls_pid);
2757
2758                 //______________________________________________________________
2759
2760                 ///______________________________________________________________________
2761                 ///Histograms for PID Studies
2762                 //Double_t fPtBin[6] = {2,4,6,8,10,15};
2763
2764                 for(Int_t i = 0; i < 6; i++)
2765                 {
2766                         if(fPt>=fPtBin[i] && fPt<fPtBin[i+1])
2767                         {
2768                                 if(fEMCflag) fEoverP_tpc[i]->Fill(fTPCnSigma,(fClus->E() / fP));
2769                                 fTPC_pt[i]->Fill(fTPCsignal);
2770                                 fTPCnsigma_pt[i]->Fill(fTPCnSigma);
2771
2772                         }
2773                 }
2774
2775                 if(track->Eta()>=fEtaCutMin && track->Eta()<=fEtaCutMax ){
2776
2777                 //new pt bins for trigger data
2778
2779                         for(Int_t i = 0; i < 10; i++)
2780                         {
2781                                 if(fP>=fPtBin_trigger[i] && fP<fPtBin_trigger[i+1])
2782                                 {
2783                                         if(fEMCflag)fEoverP_tpc_p_trigger[i]->Fill(fTPCnSigma,(fClus->E() / fP));
2784
2785                                 }
2786
2787                                 if(fPt>=fPtBin_trigger[i] && fPt<fPtBin_trigger[i+1])
2788                                 {
2789                                         if(fEMCflag)fEoverP_tpc_pt_trigger[i]->Fill(fTPCnSigma,(fClus->E() / fP));
2790
2791                                 }
2792                         }
2793
2794
2795                         //new way to calculate TPCnsigma distribution: TPCnsigma in function of p, with/without E/p cut
2796                         fTPCnsigma_p_TPC->Fill(fP, fTPCnSigma);
2797                         if(fEMCflag){
2798
2799                                 fTPCnsigma_p_TPC_on_EMCal_acc->Fill(fP, fTPCnSigma);
2800
2801                                 if((fClus->E() / fP) >= fEoverPCutMin && (fClus->E() / fP) <= fEoverPCutMax){
2802                                         fTPCnsigma_p_TPC_EoverP_cut->Fill(fP, fTPCnSigma);
2803                                 }
2804
2805                         }//close EMCflag
2806
2807                 }//close eta cut
2808
2809                 ///QA plots after track selection
2810                 ///_____________________________________________________________
2811
2812                 //_______________________________________________________
2813                 //Correlation Analysis - DiHadron
2814                 if(!fUseEMCal)
2815                 {
2816                         if(fTPCnSigma < 3.5 && fCorrelationFlag)
2817                         {
2818                                 fPtTrigger_Inc->Fill(fPt);
2819                                 DiHadronCorrelation(track, iTracks);
2820                         }
2821                 }
2822                         //_______________________________________________________
2823
2824
2825                         ///////////////////////////////////////////////////////////////////
2826                         ///TPC - efficiency calculation //
2827
2828                         /// changing start here
2829                 if(fIsMC && fIsAOD && track->GetLabel()>=0)
2830                 {
2831                         fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
2832                         Int_t pdg = fMCparticle->GetPdgCode();
2833
2834                                 //
2835                         if(fMCparticle->IsPhysicalPrimary()){
2836
2837
2838                                 if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax ){
2839
2840                                         Bool_t MotherFound = FindMother(track->GetLabel());
2841                                         if(MotherFound){
2842                                                 fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
2843                                         if( TMath::Abs(pdg) == 11 && fMCparticleMother->GetPdgCode()!=22 ){
2844                                                         if(fIsHFE1) fPtMC_TPC_All->Fill(fMCparticle->Pt());
2845                                         }
2846                                         }
2847                                 }
2848                         }
2849                 }///until here
2850
2851                 else if(fIsMC && track->GetLabel()>=0)//ESD
2852                 {
2853
2854                         if(fMCstack->IsPhysicalPrimary(track->GetLabel())){
2855                                 fMCtrack = fMCstack->Particle(track->GetLabel());
2856
2857                                 Int_t pdg = fMCtrack->GetPdgCode();
2858                                 if(fMCtrack->Eta()>=fEtaCutMin && fMCtrack->Eta()<=fEtaCutMax ){
2859
2860                                         if(fMCtrack->GetFirstMother()>0){
2861                                             fMCtrackMother = fMCstack->Particle(fMCtrack->GetFirstMother());
2862                                                 Bool_t MotherFound = FindMother(track->GetLabel());
2863                                                 if(MotherFound){
2864                                                         if( TMath::Abs(pdg) == 11 && fMCtrackMother->GetPdgCode()!=22 ){
2865                                                                 if(fIsHFE1) fPtMC_TPC_All->Fill(fMCtrack->Pt());
2866                                                         }
2867                                                 }
2868                                         }
2869                                 }
2870                         }
2871                 }
2872
2873
2874                         if(fPt>1 && fPt<2) fTOF01->Fill(fTOFnSigma,fTPCnSigma);
2875                         if(fPt>2 && fPt<4) fTOF02->Fill(fTOFnSigma,fTPCnSigma);
2876                         if(fPt>4 && fPt<6) fTOF03->Fill(fTOFnSigma,fTPCnSigma);
2877
2878 ///________________________________________________________________________
2879 ///PID
2880 ///Here the PID cuts defined in the file "ConfigEMCalHFEpA.C" is applied
2881             Int_t pidpassed = 1;
2882             AliHFEpidObject hfetrack;
2883             hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);
2884             hfetrack.SetRecTrack(track);
2885             hfetrack.SetPP();   //proton-proton analysis
2886             if(!fPID->IsSelected(&hfetrack, NULL, "", fPIDqa)) pidpassed = 0;
2887             fpid->Fill(pidpassed);
2888
2889             if(pidpassed==0) continue;
2890 ///________________________________________________________________________
2891
2892
2893 ////////////////////////////////////////////////////////////////////
2894 ///TPC efficiency calculations
2895
2896
2897                 if(fIsMC && fIsAOD && track->GetLabel()>=0)
2898                 {
2899                         fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
2900                         Int_t pdg = fMCparticle->GetPdgCode();
2901
2902                                 //
2903                         if(fMCparticle->IsPhysicalPrimary()){
2904
2905
2906                                 if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax ){
2907
2908                                         Bool_t MotherFound = FindMother(track->GetLabel());
2909                                         if(MotherFound){
2910                                                 fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
2911                                         if( TMath::Abs(pdg) == 11 && fMCparticleMother->GetPdgCode()!=22 ){
2912                                                         if(fIsHFE1) fPtMC_TPC_Selected->Fill(fMCparticle->Pt());
2913                                         }
2914                                         }
2915                                 }
2916                         }
2917                 }///until here
2918
2919                 else if(fIsMC && track->GetLabel()>=0)//ESD
2920                 {
2921
2922                         if(fMCstack->IsPhysicalPrimary(track->GetLabel())){
2923                                 fMCtrack = fMCstack->Particle(track->GetLabel());
2924
2925                                 Int_t pdg = fMCtrack->GetPdgCode();
2926                                 if(fMCtrack->Eta()>=fEtaCutMin && fMCtrack->Eta()<=fEtaCutMax ){
2927
2928                                         if(fMCtrack->GetFirstMother()>0){
2929                                             fMCtrackMother = fMCstack->Particle(fMCtrack->GetFirstMother());
2930                                                 Bool_t MotherFound = FindMother(track->GetLabel());
2931                                                 if(MotherFound){
2932                                                         if( TMath::Abs(pdg) == 11 && fMCtrackMother->GetPdgCode()!=22 ){
2933                                                                 if(fIsHFE1) fPtMC_TPC_Selected->Fill(fMCtrack->Pt());
2934                                                         }
2935                                                 }
2936                                         }
2937                                 }
2938                         }
2939                 }
2940
2941                 //Eta Cut for TPC only
2942                 if(track->Eta()>=fEtaCutMin && track->Eta()<=fEtaCutMax ){
2943                         fTPCnsigma_pt_2D->Fill(fPt,fTPCnSigma);
2944                 }
2945
2946                 //Background for TPC only
2947                 if(fFillBackground){
2948
2949                         //efficiency without SS cut for TPC only
2950                         if(track->Eta()>=fEtaCutMin && track->Eta()<=fEtaCutMax){
2951                                 Background(track, iTracks, Vtrack, kTRUE); //IsTPConly=kTRUE
2952                         } //Eta cut to be consistent with other efficiency
2953                 }
2954
2955
2956                 fTPCnsigma_p[2]->Fill(fP,fTPCnSigma);
2957                 fTPC_p[2]->Fill(fP,fTPCsignal);
2958                 TPCNcls = track->GetTPCNcls();
2959                 Float_t pos3[3]={0,0,0};
2960
2961
2962                 //here denominator for track-matching efficiency
2963                 if(fIsMC && fIsAOD && track->GetLabel()>=0)
2964                 {
2965                         fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
2966                         Int_t pdg = fMCparticle->GetPdgCode();
2967
2968                                 //
2969                         if(fMCparticle->IsPhysicalPrimary()){
2970
2971
2972                                 if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax ){
2973
2974                                         Bool_t MotherFound = FindMother(track->GetLabel());
2975                                         if(MotherFound){
2976                                                 fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
2977                                         if( TMath::Abs(pdg) == 11 && fMCparticleMother->GetPdgCode()!=22 ){
2978                                                         if(fIsHFE1) fPt_track_match_den->Fill(fMCparticle->Pt());
2979                                         }
2980                                         }
2981                                 }
2982                         }
2983                 }///until here
2984
2985
2986                 if(track->GetEMCALcluster()>0)
2987                 {
2988                         fClus = fVevent->GetCaloCluster(track->GetEMCALcluster());
2989                         ///////////////////////////////////////////////////////////////////////////////////////////////////////////
2990                         //to use tender
2991                         if(fUseTender){
2992                                 int EMCalIndex = -1;
2993                                 EMCalIndex = track->GetEMCALcluster();
2994                                 TClonesArray  *fCaloClusters_tender = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("EmcCaloClusters"));
2995                                 fClus = static_cast<AliVCluster*>(fCaloClusters_tender->At(EMCalIndex));
2996                         }
2997                         ///////////////////////////////////////////////////////////////////////////////////////////////////////////
2998
2999                         if(fClus->IsEMCAL())
3000                         {
3001
3002                 //________________________________________________________________________
3003
3004
3005                                 //Cluster timing distribution -- for ESD
3006                                 if(fUseEMCal && !fIsAOD){
3007
3008                                         AliESDCaloCells &cells_esd=*(fESD->GetEMCALCells());
3009                                         TRefArray* caloClusters_esd = new TRefArray();
3010                                     fESD->GetEMCALClusters(caloClusters_esd);
3011                                         Int_t nclus_esd = caloClusters_esd->GetEntries();
3012
3013
3014                                         for (Int_t icl = 0; icl < nclus_esd; icl++) {
3015
3016                                                 AliESDCaloCluster* clus_esd = (AliESDCaloCluster*)caloClusters_esd->At(icl);
3017
3018                                                 if(clus_esd->IsEMCAL()){
3019                                                         Float_t ncells_esd      = fClus->GetNCells();
3020                                                         UShort_t * index_esd = clus_esd->GetCellsAbsId() ;
3021                                                         UShort_t * index2_esd = fClus->GetCellsAbsId() ;
3022
3023
3024                                                         for(Int_t i = 0; i < ncells_esd ; i++){
3025
3026                                                                 Int_t absId_esd =   index_esd[i];
3027                                                                 fTime->Fill(fPt,cells_esd.GetCellTime(absId_esd));
3028
3029                                                                 Int_t absId2_esd =   index2_esd[i];
3030                                                                 fTime2->Fill(fPt,cells_esd.GetCellTime(absId2_esd));
3031                                                         }
3032
3033                                                 }
3034                                         }
3035
3036                                 }
3037                                 /* not working!
3038                                 //Cluster timing distribution -- for AOD
3039                                 if(fUseEMCal && fIsAOD){
3040
3041                                         AliAODCaloCells &cells_aod=*(fAOD->GetEMCALCells());
3042
3043                                         TRefArray* caloClusters_aod = new TRefArray();
3044                                         fAOD->GetEMCALClusters(caloClusters_aod);
3045
3046                                         Int_t nclus_aod = caloClusters_aod->GetEntries();
3047
3048                                         for (Int_t icl = 0; icl < nclus_aod; icl++) {
3049
3050                                                 AliAODCaloCluster* clus_aod = (AliAODCaloCluster*)caloClusters_aod->At(icl);
3051
3052
3053                                                 if(clus_aod->IsEMCAL()){
3054                                                         Float_t ncells_aod      = fClus->GetNCells();
3055                                                         UShort_t * index_aod = clus_aod->GetCellsAbsId() ;
3056                                                         UShort_t * index2_aod = fClus->GetCellsAbsId() ;
3057
3058
3059                                                         for(Int_t i = 0; i < ncells_aod ; i++){
3060
3061                                                                 Int_t absId_aod =   index_aod[i];
3062                                                                 fTime->Fill(fPt,cells_aod.GetCellTime(absId_aod));
3063
3064                                                                 Int_t absId2_aod =   index2_aod[i];
3065                                                                 fTime2->Fill(fPt,cells_aod.GetCellTime(absId2_aod));
3066                                                         }
3067
3068                                                 }
3069                                         }
3070
3071                                 }
3072                                  */
3073
3074
3075                                 if(fUseEMCal){
3076                                         double emctof = fClus->GetTOF();
3077                                         ftimingEle->Fill(fPt,emctof);
3078                                 }
3079                                         //________________________________________________________________________
3080
3081
3082
3083
3084                                 // Residuals
3085                                 Double_t Dx = fClus->GetTrackDx();
3086                                 Double_t Dz = fClus->GetTrackDz();
3087                                 Double_t R=TMath::Sqrt(Dx*Dx+Dz*Dz);
3088
3089                                 for(Int_t i = 0; i < 6; i++)
3090                                 {
3091                                         if(fPt>=fPtBin[i] && fPt<fPtBin[i+1])
3092                                         {
3093
3094                                                 fEta[i]->Fill(Dz);
3095                                                 fPhi[i]->Fill(Dx);
3096                                                 fR[i]->Fill(R);
3097                                         }
3098                                 }
3099
3100                                 if(TMath::Abs(fClus->GetTrackDx())<=fdPhiCut && TMath::Abs(fClus->GetTrackDz())<=fdEtaCut)
3101                                 {
3102                                         Float_t Energy  = fClus->E();
3103                                         Float_t EoverP  = Energy/track->P();
3104                                         Float_t M02     = fClus->GetM02();
3105                                         Float_t M20     = fClus->GetM20();
3106                                         Float_t ncells  = fClus->GetNCells();
3107
3108                                         //here numerator for track-matching efficiency
3109                                         if(fIsMC && fIsAOD && track->GetLabel()>=0)
3110                                         {
3111                                                 fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
3112                                                 Int_t pdg = fMCparticle->GetPdgCode();
3113
3114                                                         //
3115                                                 if(fMCparticle->IsPhysicalPrimary()){
3116
3117
3118                                                         if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax ){
3119
3120                                                                 Bool_t MotherFound = FindMother(track->GetLabel());
3121                                                                 if(MotherFound){
3122                                                                         fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
3123                                                                         if( TMath::Abs(pdg) == 11 && fMCparticleMother->GetPdgCode()!=22 ){
3124                                                                                 if(fIsHFE1) fPt_track_match_num->Fill(fMCparticle->Pt());
3125                                                                         }
3126                                                                 }
3127                                                         }
3128                                                 }
3129                                         }///until here
3130
3131
3132
3133                                         //----------------------------------------------------------------------------------------
3134                                         //
3135                                         //EtaCut electrons histogram
3136                                         //Shower Shape Cut
3137                                         if(track->Eta()>=fEtaCutMin && track->Eta()<=fEtaCutMax ){
3138
3139                                                 if(fUseShowerShapeCut){
3140                                                         if(M02 >= fM02CutMin && M02<=fM02CutMax && M20>=fM20CutMin && M20<=fM20CutMax){
3141                                                                 fEoverP_pt[2]->Fill(fPt,(fClus->E() / fP));
3142                                                                 fShowerShapeCut->Fill(M02,M20);
3143                                                                 //in order to check if there are exotic cluster in this selected cluster (27 may 2014)
3144                                                                 fNcells_energy_elec_selected->Fill(ncells,Energy);
3145
3146                                                                         //true electrons E/p shape
3147                                                                 if(fIsMC && fIsAOD && track->GetLabel()>=0){
3148                                                                         fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
3149                                                                         Int_t pdg = fMCparticle->GetPdgCode();
3150
3151                                                                         if( TMath::Abs(pdg) == 11){
3152                                                                                 fEoverP_pt_true_electrons->Fill(fPt,(fClus->E() / fP));
3153                                                                         }
3154                                                                 }
3155
3156
3157
3158                                                         }
3159
3160                                                 }
3161                                                 if(!fUseShowerShapeCut){
3162                                                         fEoverP_pt[2]->Fill(fPt,(fClus->E() / fP));
3163                                                         fShowerShapeCut->Fill(M02,M20);
3164                                                         fNcells_energy_elec_selected->Fill(ncells,Energy);
3165
3166                                                         //after cut -> Are they really electrons ? E/p shape
3167                                                         if(fIsMC && fIsAOD && track->GetLabel()>=0){
3168                                                                 fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
3169                                                                 Int_t pdg = fMCparticle->GetPdgCode();
3170
3171                                                                 if( TMath::Abs(pdg) == 11){
3172                                                                         fEoverP_pt_true_electrons->Fill(fPt,(fClus->E() / fP));
3173                                                                 }
3174                                                         }
3175
3176
3177                                                 }
3178                                                 if(fUseEMCal) fShowerShape_ele->Fill(M02,M20);
3179
3180                                                 //for shower shape cut studies - now with TPC PID Cut
3181                                                 if(fUseEMCal){
3182                                                         fShowerShapeM02_EoverP->Fill(M02,EoverP);
3183                                                         fShowerShapeM20_EoverP->Fill(M20,EoverP);
3184                                                 }
3185
3186                                         }
3187
3188                                         //----------------------------------------------------------------------------------------
3189
3190
3191
3192                                         // for Eta Phi distribution
3193                                         fClus->GetPosition(pos3);
3194                                         TVector3 vpos(pos3[0],pos3[1],pos3[2]);
3195                                         Double_t cphi = vpos.Phi();
3196                                         Double_t ceta = vpos.Eta();
3197                                         fEtaPhi[2]->Fill(cphi,ceta);
3198
3199
3200
3201                                         fTPCNcls_EoverP[2]->Fill(TPCNcls, EoverP);
3202
3203                                         for(Int_t i = 0; i < 6; i++)
3204                                         {
3205                                                 if(fPt>=fPtBin[i] && fPt<fPtBin[i+1])
3206                                                 {
3207
3208                                                         fR_EoverP[i]->Fill(R, EoverP);
3209                                                         fNcells[i]->Fill(ncells);
3210                                                         fNcells_EoverP[i]->Fill(EoverP, ncells);
3211                                                         fM02_EoverP[i]->Fill(M02,EoverP);
3212                                                         fM20_EoverP[i]->Fill(M20,EoverP);
3213                                                         fECluster_ptbins[i]->Fill(Energy);
3214                                                         fEoverP_ptbins[i]->Fill(EoverP);
3215
3216                                                         if((fClus->E() / fP) >= fEoverPCutMin && (fClus->E() / fP) <= fEoverPCutMax) {
3217                                                                 fNcells_electrons[i]->Fill(ncells);
3218                                                         }
3219
3220                                                         if(M02<0.5 && M20<0.3) {
3221                                                                 fEoverP_wSSCut[i]->Fill(EoverP);
3222                                                         }
3223                                                 }
3224                                         }
3225
3226                                         fNcells_pt->Fill(fPt, ncells);
3227
3228
3229                                         ////////////////////////////////////////////////////////////////////
3230                                         ///EMCal - Efficiency calculations
3231
3232
3233                                         if(fIsMC && fIsAOD && track->GetLabel()>=0)
3234                                         {
3235                                                 fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
3236                                                 Int_t pdg = fMCparticle->GetPdgCode();
3237
3238                                                 //
3239                                                 if(fMCparticle->IsPhysicalPrimary()){
3240
3241                                                         //Phi cut && fMCparticle->Phi()>=(TMath::Pi()*80/180) && fMCparticle->Phi()<=TMath::Pi()
3242                                                         if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax  ){
3243
3244                                                                 Bool_t MotherFound = FindMother(track->GetLabel());
3245                                                                 if(MotherFound){
3246                                                                         fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
3247                                                                         if( TMath::Abs(pdg) == 11 && fMCparticleMother->GetPdgCode()!=22 ){
3248                                                                                 if(fIsHFE1)fPtMC_EMCal_All->Fill(fMCparticle->Pt());
3249                                                                         }
3250                                                                 }
3251                                                         }
3252                                                 }
3253                                         }
3254
3255                                         else if(fIsMC && track->GetLabel()>=0)//ESD
3256                                         {
3257
3258                                                 if(fMCstack->IsPhysicalPrimary(track->GetLabel()))
3259                                                 {
3260
3261                                                         fMCtrack = fMCstack->Particle(track->GetLabel());
3262
3263                                                         Int_t pdg = fMCtrack->GetPdgCode();
3264                                                         //Phi cut && fMCtrack->Phi()>=(TMath::Pi()*80/180) && fMCtrack->Phi()<=TMath::Pi()
3265                                                         if(fMCtrack->Eta()>=fEtaCutMin && fMCtrack->Eta()<=fEtaCutMax )
3266                                                         {
3267                                                                 Bool_t MotherFound = FindMother(track->GetLabel());
3268                                                                         //MotherFound included
3269                                                                 if(MotherFound){
3270                                                                         if(fMCtrack->GetFirstMother()>0){
3271                                                                                 fMCtrackMother = fMCstack->Particle(fMCtrack->GetFirstMother());
3272                                                                                 if( TMath::Abs(pdg) == 11 && fMCtrackMother->GetPdgCode()!=22 ){
3273                                                                                         if(fIsHFE1)fPtMC_EMCal_All->Fill(fMCtrack->Pt());
3274                                                                                 }
3275                                                                         }
3276                                                                 }
3277                                                         }
3278                                                 }
3279                                         }
3280
3281                                         //_______________________________________________________
3282                                         //PID using EMCal
3283                                         if(fEoverPnsigma)SetEoverPCutPtDependentMC(track->Pt());
3284
3285                                         if((fClus->E() / fP) >= fEoverPCutMin && (fClus->E() / fP) <= fEoverPCutMax)
3286                                         {
3287
3288
3289                                                         fECluster[2]->Fill(Energy);
3290                                                         fTPCNcls_pid[3]->Fill(TPCNcls, TPCNcls_pid);
3291
3292
3293                                                 if(fUseEMCal)
3294                                                 {
3295                                                         if(track->Eta()>=fEtaCutMin && track->Eta()<=fEtaCutMax ){
3296                                                                 fPtElec_Inc->Fill(fPt);
3297                                                                 //eta phi distribution for data
3298                                                                 fEtaPhi_data->Fill(track->Phi(),track->Eta());
3299                                                         }
3300
3301                                                         //Eta cut for background
3302                                                         if(fFillBackground){
3303                                                                 fEtad[2]->Fill(track->Eta());
3304
3305                                                                 //background for triggered data: trigger electron must have same cuts on shower shape  06/Jan/2014
3306                                                                 if(fUseShowerShapeCut){
3307                                                                         if(M02 >= fM02CutMin && M02<=fM02CutMax && M20>=fM20CutMin && M20<=fM20CutMax){
3308                                                                                 if(track->Eta()>=fEtaCutMin && track->Eta()<=fEtaCutMax){
3309                                                                                         Background(track, iTracks, Vtrack, kFALSE);
3310                                                                                 }
3311                                                                         }
3312                                                                 }
3313                                                                 else{
3314                                                                         if(track->Eta()>=fEtaCutMin && track->Eta()<=fEtaCutMax){
3315                                                                                 Background(track, iTracks, Vtrack, kFALSE);
3316                                                                         }
3317                                                                 }
3318
3319                                                         }
3320
3321                                                         double emctof2 = fClus->GetTOF();
3322                                                         ftimingEle2->Fill(fPt,emctof2);
3323                                                         //Correlation Analysis
3324                                                         if(fCorrelationFlag)
3325                                                         {
3326                                                                 ElectronHadronCorrelation(track, iTracks, Vtrack);
3327                                                         }
3328                                                 }
3329                                                 //_______________________________________________________
3330
3331                                                 ////////////////////////////////////////////////////////////////////
3332                                                 ///EMCal - efficiency calculations
3333
3334                                                 if(track->Charge()<0)  fCharge_n->Fill(fPt);
3335                                                 if(track->Charge()>0)  fCharge_p->Fill(fPt);
3336
3337
3338
3339                                                 if(fIsMC && fIsAOD && track->GetLabel()>=0)//AOD
3340                                                 {
3341                                                         if(track->Charge()<0)  fCharge_n->Fill(fPt);
3342                                                         if(track->Charge()>0)  fCharge_p->Fill(fPt);
3343
3344                                                         fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
3345                                                         if(fMCparticle->GetMother()>0) fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
3346                                                         Int_t pdg = fMCparticle->GetPdgCode();
3347
3348                                                         double proX = fMCparticle->Xv();
3349                                                         double proY = fMCparticle->Yv();
3350                                                         double proR = sqrt(pow(proX,2)+pow(proY,2));
3351
3352
3353                                                         if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax && fMCparticle->Phi()>=(TMath::Pi()*80/180) && fMCparticle->Phi()<=TMath::Pi()  ){
3354
3355                                                                 if( TMath::Abs(pdg) == 11 && fMCparticle->GetMother()>0 ){
3356                                                                         Int_t mpdg = fMCparticleMother->GetPdgCode();
3357                                                                         if(TMath::Abs(mpdg) == 221 || TMath::Abs(mpdg) == 22 || TMath::Abs(mpdg) == 111){
3358                                                                                 if(proR<7)fPtMCelectronAfterAll_nonPrimary->Fill(fMCparticle->Pt()); //numerator for the total efficiency, non Primary track
3359                                                                         }
3360                                                                 }
3361                                                                 if( TMath::Abs(pdg) == 11 && fMCparticle->IsPhysicalPrimary()) fPtMCelectronAfterAll_Primary->Fill(fMCparticle->Pt());
3362                                                         }
3363
3364
3365                                                         //
3366                                                         if(fMCparticle->IsPhysicalPrimary()){
3367
3368
3369                                                                 if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax ){
3370
3371                                                                         Bool_t MotherFound = FindMother(track->GetLabel());
3372                                                                         if(MotherFound){
3373
3374                                                                                 if(!fUseShowerShapeCut){
3375                                                                                         if(fIsHFE1){
3376                                                                                                 //Unfolding   pt_reco/pt_MC  in the efficiency
3377                                                                                                 fPtMCelectronAfterAll->Fill(fMCparticle->Pt());
3378                                                                                                 fPtMCelectronAfterAll_unfolding->Fill(track->Pt());
3379                                                                                                 fEtaPhi_num->Fill(fMCparticle->Phi(),fMCparticle->Eta());
3380
3381                                                                                                 //new histo to estimate how different is pt reco from pt MC
3382                                                                                                 fpt_reco_pt_MC_num->Fill(track->Pt(),fMCparticle->Pt());
3383                                                                                         }//numerator for the total efficiency AOD
3384                                                                                 }
3385                                                                                         //November 11 for efficiency of triggered data
3386                                                                                 if(fUseShowerShapeCut){
3387                                                                                         if(M02 >= fM02CutMin && M02<=fM02CutMax && M20>=fM20CutMin && M20<=fM20CutMax){
3388                                                                                                 if(fIsHFE1){
3389                                                                                                         //Unfolding   pt_reco/pt_MC  in the efficiency
3390                                                                                                         fPtMCelectronAfterAll->Fill(fMCparticle->Pt());
3391                                                                                                         fPtMCelectronAfterAll_unfolding->Fill(track->Pt());
3392                                                                                                         fEtaPhi_num->Fill(fMCparticle->Phi(),fMCparticle->Eta());
3393
3394                                                                                                         //new histo to estimate how different is pt reco from pt MC
3395                                                                                                         fpt_reco_pt_MC_num->Fill(track->Pt(),fMCparticle->Pt());
3396                                                                                                 }//numerator for the total efficiency AOD
3397                                                                                         }
3398                                                                                 }
3399
3400
3401
3402                                                                                         fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
3403                                                                                         if( TMath::Abs(pdg) == 11 && fMCparticleMother->GetPdgCode()!=22 ){
3404                                                                                                 if(fIsHFE1)fPtMC_EMCal_Selected->Fill(fMCparticle->Pt());
3405                                                                                         }
3406                                                                         }//if MotherFound
3407                                                                 }//eta cut
3408                                                         }
3409                                                 }///close AOD
3410
3411                                                 else if(fIsMC && track->GetLabel()>=0)//ESD
3412                                                 {
3413                                                         if(track->Charge()<0)  fCharge_n->Fill(fPt);
3414                                                         if(track->Charge()>0)  fCharge_p->Fill(fPt);
3415
3416                                                         fMCtrack = fMCstack->Particle(track->GetLabel());
3417                                                         if(fMCtrack->GetFirstMother()>0)
3418                                                         {
3419                                                                 fMCtrackMother = fMCstack->Particle(fMCtrack->GetFirstMother());
3420                                                         }
3421                                                         TParticle *particle=fMCstack->Particle(track->GetLabel());
3422
3423                                                         Int_t pdg = fMCtrack->GetPdgCode();
3424
3425
3426                                                         if(fMCtrack->Eta()>=fEtaCutMin && fMCtrack->Eta()<=fEtaCutMax)
3427                                                         {
3428                                                                 if( TMath::Abs(pdg) == 11 && fMCtrack->GetFirstMother()>0 )
3429                                                                 {
3430                                                                         Int_t mpdg = fMCtrackMother->GetPdgCode();
3431                                                                         if(TMath::Abs(mpdg) == 221 || TMath::Abs(mpdg) == 22 || TMath::Abs(mpdg) == 111)
3432                                                                         {
3433                                                                                 Double_t proR=particle->R();
3434                                                                                 if(proR<7)
3435                                                                                 {
3436                                                                                   fPtMCelectronAfterAll_nonPrimary->Fill(fMCtrack->Pt()); //numerator for the total efficiency, non Primary track
3437                                                                                 }
3438                                                                         }
3439                                                                 }
3440                                                                 if( TMath::Abs(pdg) == 11 && fMCstack->IsPhysicalPrimary(track->GetLabel()))
3441                                                                 {
3442                                                                         fPtMCelectronAfterAll_Primary->Fill(fMCtrack->Pt());
3443                                                                 }
3444                                                         }
3445
3446                                                         if(fMCstack->IsPhysicalPrimary(track->GetLabel()))
3447                                                         {
3448                                                                 Bool_t MotherFound = FindMother(track->GetLabel());
3449
3450                                                                 if(MotherFound)
3451                                                                 {
3452                                                                         if(fMCtrack->Eta()>=fEtaCutMin && fMCtrack->Eta()<=fEtaCutMax)
3453                                                                         {
3454                                                                                 if(!fUseShowerShapeCut){
3455                                                                                         if(fIsHFE1){
3456                                                                                                 fPtMCelectronAfterAll->Fill(fMCtrack->Pt()); //numerator for the total efficiency ESD
3457                                                                                                 fEtaPhi_num->Fill(fMCtrack->Phi(),fMCtrack->Eta());
3458                                                                                         }
3459                                                                                 }
3460                                                                                 //November 11 for efficiency of triggered data
3461                                                                                 if(fUseShowerShapeCut){
3462                                                                                         if(M02 >= fM02CutMin && M02<=fM02CutMax && M20>=fM20CutMin && M20<=fM20CutMax){
3463                                                                                                 if(fIsHFE1){
3464                                                                                                         fPtMCelectronAfterAll->Fill(fMCtrack->Pt()); //numerator for the total efficiency ESD
3465                                                                                                         fEtaPhi_num->Fill(fMCtrack->Phi(),fMCtrack->Eta());
3466                                                                                                 }
3467                                                                                         }
3468                                                                                 }
3469
3470
3471
3472                                                                         }
3473                                                                 }
3474
3475
3476
3477                                                                 // Phi cut && fMCtrack->Phi()>=(TMath::Pi()*80/180) && fMCtrack->Phi()<=TMath::Pi()
3478                                                                 if(fMCtrack->Eta()>=fEtaCutMin && fMCtrack->Eta()<=fEtaCutMax)
3479                                                                 {
3480                                                                         //included MotherFound
3481
3482                                                                         if(MotherFound)
3483                                                                         {
3484                                                                                 if(fMCtrack->GetFirstMother()>0){
3485                                                                                         fMCtrackMother = fMCstack->Particle(fMCtrack->GetFirstMother());
3486                                                                                         if( TMath::Abs(pdg) == 11 && fMCtrackMother->GetPdgCode()!=22 ){
3487
3488                                                                                                 if(fIsHFE1){fPtMC_EMCal_Selected->Fill(fMCtrack->Pt());}
3489                                                                                         }
3490                                                                                 }
3491                                                                         }
3492                                                                 }
3493                                                         }
3494                                                 }//close ESD
3495                                                 ///////////////////////////////////////////////////////////////////
3496
3497
3498                                         }
3499                                 }
3500                         }
3501                 }
3502
3503                         //______________________________________________________________
3504                         // Vertex
3505
3506                 fVtxZ[2]->Fill(fZvtx);
3507                 if(iTracks == 0)fNTracks[2]->Fill(fNOtrks);
3508                 fNTracks_pt[2]->Fill(fNOtrks, fPt);
3509                 fNTracks_eta[2]->Fill(fNOtrks, track->Eta());
3510                 fNTracks_phi[2]->Fill(fNOtrks, track->Phi());
3511                 fNClusters[2]->Fill(ClsNo);
3512                 fTPCNcls_pid[2]->Fill(TPCNcls, TPCNcls_pid);
3513
3514                         //______________________________________________________________
3515
3516                         //_______________________________________________________
3517                         //Correlation Analysis
3518                 if(!fUseEMCal)
3519                 {
3520                         fPtElec_Inc->Fill(fPt);
3521
3522                         if(fCorrelationFlag)
3523                         {
3524                                 ElectronHadronCorrelation(track, iTracks, Vtrack);
3525                         }
3526                 }
3527                         //_______________________________________________________
3528
3529                         ///________________________________________________________________________
3530         }
3531
3532                 //__________________________________________________________________
3533                 //Event Mixing Analysis
3534                 //Filling pool
3535         if(fEventMixingFlag)
3536         {
3537                 fPool = fPoolMgr->GetEventPool(fCentrality->GetCentralityPercentile("V0A"), fZvtx); // Get the buffer associated with the current centrality and z-vtx
3538
3539                 if(!fPool) AliFatal(Form("No pool found for centrality = %f, zVtx = %f", fCentrality->GetCentralityPercentile("V0A"), fZvtx));
3540
3541                 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!
3542
3543
3544         }
3545
3546                 //__________________________________________________________________
3547
3548         delete fListOfmotherkink;
3549         PostData(1, fOutputList);
3550 }
3551
3552         //______________________________________________________________________
3553 void AliAnalysisTaskEMCalHFEpA::Terminate(Option_t *)
3554 {
3555                 //Draw result to the screen
3556                 //Called once at the end of the query
3557
3558         fOutputList = dynamic_cast<TList*> (GetOutputData(1));
3559
3560         if(!fOutputList)
3561         {
3562                 printf("ERROR: Output list not available\n");
3563                 return;
3564         }
3565 }
3566
3567         //______________________________________________________________________
3568 Bool_t AliAnalysisTaskEMCalHFEpA::ProcessCutStep(Int_t cutStep, AliVParticle *track)
3569 {
3570                 //Check single track cuts for a given cut step
3571                 //Note this function is called inside the UserExec function
3572         const Int_t kMCOffset = AliHFEcuts::kNcutStepsMCTrack;
3573         if(!fCFM->CheckParticleCuts(cutStep + kMCOffset, track)) return kFALSE;
3574         return kTRUE;
3575 }
3576
3577
3578         //______________________________________________________________________
3579
3580
3581 void AliAnalysisTaskEMCalHFEpA::Background(AliVTrack *track, Int_t trackIndex, AliVParticle *vtrack, Bool_t IsTPConly)
3582 {
3583         ///_________________________________________________________________
3584         ///MC analysis
3585         //Bool_t IsMCefix=kFALSE; //to make correction on efix, use kTRUE (do not change the efficiency, so I will keep the correction only for d3)
3586
3587                 if(fIsMC)
3588                 {
3589                         if(track->GetLabel() < 0)
3590                         {
3591                                 AliWarning(Form("The track %d does not have a valid MC label",trackIndex));
3592                                 return;
3593                         }
3594
3595                         if(fIsAOD)
3596                         {
3597                                 fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
3598
3599                                 if(fMCparticle->GetMother()<0) return;
3600
3601                                 fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
3602                                 if(fMCparticleMother->GetMother()>0)fMCparticleGMother = (AliAODMCParticle*) fMCarray->At(fMCparticleMother->GetMother());
3603
3604                                 if(TMath::Abs(fMCparticle->GetPdgCode())==11 && (TMath::Abs(fMCparticleMother->GetPdgCode())==22 || TMath::Abs(fMCparticleMother->GetPdgCode())==111 || TMath::Abs(fMCparticleMother->GetPdgCode())==221))
3605                                 {
3606                                                 //Is Background
3607                                         if(!IsTPConly)fPtBackgroundBeforeReco->Fill(track->Pt());
3608                                         if(IsTPConly)fPtBackgroundBeforeReco2->Fill(track->Pt());
3609
3610
3611                                                 //October 08th weighted histos
3612                                         if(TMath::Abs(fMCparticleMother->GetPdgCode())==111 || TMath::Abs(fMCparticleMother->GetPdgCode())==221 ){
3613
3614                                                 Double_t mPt=fMCparticleMother->Pt();
3615                                                 Double_t mweight=1;
3616
3617
3618                                                 //________________________________________________________________
3619                                                 //correction for d3 based on data
3620
3621                                                         if(TMath::Abs(fMCparticleMother->GetPdgCode())==111){
3622                                                                 Double_t x=mPt;
3623                                                                 mweight=CalculateWeight(111, x);
3624
3625                                                         }
3626                                                     if(TMath::Abs(fMCparticleMother->GetPdgCode())==221){
3627                                                             Double_t x=mPt;
3628                                                             mweight=CalculateWeight(221, x);
3629
3630                                                     }
3631
3632                                                 //________________________________________________________________
3633
3634                                                 //Histo pT mother versus pT electron
3635                                                 fpT_m_electron->Fill(mPt, track->Pt());
3636
3637                                                 if(!IsTPConly)fPtBackgroundBeforeReco_weight->Fill(track->Pt(), 1./mweight);
3638                                                 if(IsTPConly)fPtBackgroundBeforeReco2_weight->Fill(track->Pt(), 1./mweight);
3639                                         }
3640                                         else if(fMCparticleMother->GetMother()>0 && (TMath::Abs(fMCparticleGMother->GetPdgCode())==111 || TMath::Abs(fMCparticleGMother->GetPdgCode())==221 )){
3641
3642                                                 Double_t gmPt=fMCparticleGMother->Pt();
3643                                                 Double_t gmweight=1;
3644
3645                                                 //________________________________________________________________
3646                                                 //correction for d3 based on data
3647
3648                                                         if(TMath::Abs(fMCparticleGMother->GetPdgCode())==111){
3649                                                                 Double_t x=gmPt;
3650                                                                 gmweight=CalculateWeight(111, x);
3651                                                         }
3652                                                         if(TMath::Abs(fMCparticleGMother->GetPdgCode())==221){
3653                                                                 Double_t x=gmPt;
3654                                                                 gmweight=CalculateWeight(221, x);
3655                                                         }
3656
3657
3658                                                 //________________________________________________________________
3659                                                 //Histo pT gmother versus pT electron
3660
3661                                                 fpT_gm_electron->Fill(gmPt, track->Pt());
3662
3663                                                 if(!IsTPConly)fPtBackgroundBeforeReco_weight->Fill(track->Pt(), 1./gmweight);
3664                                                 if(IsTPConly)fPtBackgroundBeforeReco2_weight->Fill(track->Pt(), 1./gmweight);
3665                                         }
3666                                         else{
3667                                                 if(!IsTPConly)fPtBackgroundBeforeReco_weight->Fill(track->Pt());
3668                                                 if(IsTPConly)fPtBackgroundBeforeReco2_weight->Fill(track->Pt());
3669                                         }
3670                                 }//particle kind
3671                         }//IsAOD
3672                          //ESD
3673                         else
3674                         {
3675                                 fMCtrack = fMCstack->Particle(track->GetLabel());
3676
3677                                 if(fMCtrack->GetFirstMother()<0) return;
3678
3679                                 fMCtrackMother = fMCstack->Particle(fMCtrack->GetFirstMother());
3680
3681                                 if(TMath::Abs(fMCtrack->GetPdgCode())==11 && (TMath::Abs(fMCtrackMother->GetPdgCode())==22 || TMath::Abs(fMCtrackMother->GetPdgCode())==111 || TMath::Abs(fMCtrackMother->GetPdgCode())==221))
3682                                 {
3683                                                 //Is Background
3684                                         if(!IsTPConly)fPtBackgroundBeforeReco->Fill(track->Pt());
3685                                         if(IsTPConly)fPtBackgroundBeforeReco2->Fill(track->Pt());
3686                                 }
3687                         }
3688                 }//IsMC
3689
3690                 ///_________________________________________________________________
3691
3692                 //________________________________________________
3693                 //Associated particle cut
3694                 fPartnerCuts->SetAcceptKinkDaughters(kFALSE);
3695                 fPartnerCuts->SetRequireITSRefit(kTRUE);
3696                 fPartnerCuts->SetRequireTPCRefit(kTRUE);
3697                 fPartnerCuts->SetEtaRange(-0.9,0.9);
3698                 fPartnerCuts->SetMaxChi2PerClusterTPC(4.0);
3699                 fPartnerCuts->SetMinNClustersTPC(80);
3700                 fPartnerCuts->SetPtRange(0,1e10);
3701                 //fPartnerCuts->SetRequireSigmaToVertex(kTRUE);
3702                 //fPartnerCuts->SetMaxDCAToVertexXY(1);
3703                 //fPartnerCuts->SetMaxDCAToVertexZ(3);
3704                 //_________________________________________________
3705
3706                 ///#################################################################
3707                 //Non-HFE reconstruction
3708                 fNonHFE = new AliSelectNonHFE();
3709                 fNonHFE->SetAODanalysis(fIsAOD);
3710                 if(fMassCutFlag) fNonHFE->SetInvariantMassCut(fMassCut);
3711                 if(fAngleCutFlag) fNonHFE->SetOpeningAngleCut(fAngleCut);
3712                 if(fChi2CutFlag) fNonHFE->SetChi2OverNDFCut(fChi2Cut);
3713                 if(fDCAcutFlag) fNonHFE->SetDCACut(fDCAcut);
3714                 fNonHFE->SetAlgorithm("DCA"); //KF
3715                 fNonHFE->SetPIDresponse(fPidResponse);
3716                 fNonHFE->SetTrackCuts(-3.5,3.5,fPartnerCuts);
3717                 fNonHFE->SetAdditionalCuts(fPtMinAsso,fTpcNclsAsso);
3718
3719                 if(!IsTPConly){
3720                         fNonHFE->SetHistAngleBack(fOpAngleBack);
3721                         fNonHFE->SetHistAngle(fOpAngle);
3722                         fNonHFE->SetHistDCABack(fDCABack);
3723                         fNonHFE->SetHistDCA(fDCA);
3724                         fNonHFE->SetHistMassBack(fInvMassBack);
3725                         fNonHFE->SetHistMass(fInvMass);
3726                 }
3727                 if(IsTPConly){
3728                         fNonHFE->SetHistAngleBack(fOpAngleBack2);
3729                         fNonHFE->SetHistAngle(fOpAngle2);
3730                         fNonHFE->SetHistDCABack(fDCABack2);
3731                         fNonHFE->SetHistDCA(fDCA2);
3732                         fNonHFE->SetHistMassBack(fInvMassBack2);
3733                         fNonHFE->SetHistMass(fInvMass2);
3734                 }
3735
3736                 fNonHFE->FindNonHFE(trackIndex,vtrack,fVevent);
3737
3738                         //index of track selected as partner
3739                 Int_t *fUlsPartner = fNonHFE->GetPartnersULS();
3740
3741                         //Electron Information
3742                 Double_t fPhiE = -999;
3743                 Double_t fEtaE = -999;
3744                 Double_t fPtE = -999;
3745                 fPhiE = track->Phi();
3746                 fEtaE = track->Eta();
3747                 fPtE = track->Pt();
3748
3749                         ///_________________________________________________________________
3750                         ///MC analysis
3751                 if(fIsMC)
3752                 {
3753                         if(fIsAOD)
3754                         {
3755                                 if(TMath::Abs(fMCparticle->GetPdgCode())==11 && (TMath::Abs(fMCparticleMother->GetPdgCode())==22 || TMath::Abs(fMCparticleMother->GetPdgCode())==111 || TMath::Abs(fMCparticleMother->GetPdgCode())==221))
3756                                 {
3757
3758                                         Double_t weight=1;
3759
3760                                         if(!IsTPConly){
3761                                                 if(fNonHFE->IsULS()) fPtElec_ULS->Fill(fPtE,fNonHFE->GetNULS());
3762                                                 if(fNonHFE->IsLS()) fPtElec_LS->Fill(fPtE,fNonHFE->GetNLS());
3763
3764
3765
3766
3767                                                 if(TMath::Abs(fMCparticleMother->GetPdgCode())==111 || TMath::Abs(fMCparticleMother->GetPdgCode())==221){
3768                                                         Double_t mPt=fMCparticleMother->Pt();
3769                                                                 Double_t mweight1=1;
3770                                                                 Double_t mweight2=1;
3771                                                                 //Double_t weight=1;
3772
3773
3774                                                          //----------------------------------------------------------------------------
3775                                                          //correction based on data only
3776                                                                 if(TMath::Abs(fMCparticleMother->GetPdgCode())==111){
3777                                                                         Double_t x=mPt;
3778                                                                         weight=CalculateWeight(111, x);
3779                                                                 }
3780                                                                 if(TMath::Abs(fMCparticleMother->GetPdgCode())==221){
3781                                                                         Double_t x=mPt;
3782                                                                         weight=CalculateWeight(221, x);
3783                                                                 }
3784
3785
3786                                                                 //----------------------------------------------------------------------------
3787
3788                                                                 //check this
3789                                                         if(fNonHFE->IsULS()) mweight1=(fNonHFE->GetNULS())/weight;
3790                                                         if(fNonHFE->IsLS())  mweight2=(fNonHFE->GetNLS())/weight;
3791
3792                                                                 //fill histos
3793                                                         if(fNonHFE->IsULS())fPtElec_ULS_weight->Fill(fPtE, mweight1);
3794                                                         if(fNonHFE->IsLS())fPtElec_LS_weight->Fill(fPtE, mweight2);
3795                                                 }
3796                                                 else if(fMCparticleMother->GetMother()>0 && (TMath::Abs(fMCparticleGMother->GetPdgCode())==111 || TMath::Abs(fMCparticleGMother->GetPdgCode())==221 )){
3797                                                         Double_t gmPt=fMCparticleGMother->Pt();
3798                                                                 Double_t gmweight1=1;
3799                                                                 Double_t gmweight2=1;
3800                                                                 //Double_t weight=1;
3801
3802                                                         //----------------------------------------------------------------------------
3803
3804                                                         //----------------------------------------------------------------------------
3805
3806                                                         //correction based on data only for pi0
3807
3808                                                                 if(TMath::Abs(fMCparticleGMother->GetPdgCode())==111){
3809                                                                         Double_t x=gmPt;
3810                                                                         weight=CalculateWeight(111, x);
3811                                                                 }
3812                                                                 if(TMath::Abs(fMCparticleGMother->GetPdgCode())==221){
3813                                                                         Double_t x=gmPt;
3814                                                                         weight=CalculateWeight(221, x);
3815                                                         }
3816
3817
3818
3819
3820
3821                                                                 //check this
3822                                                         if(fNonHFE->IsULS()) gmweight1=(fNonHFE->GetNULS())/weight;
3823                                                         if(fNonHFE->IsLS())  gmweight2=(fNonHFE->GetNLS())/weight;
3824
3825                                                                 //fill histos
3826                                                         if(fNonHFE->IsULS())fPtElec_ULS_weight->Fill(fPtE, gmweight1);
3827                                                         if(fNonHFE->IsLS())fPtElec_LS_weight->Fill(fPtE, gmweight2);
3828                                                 }
3829                                                 else{
3830                                                         if(fNonHFE->IsULS()) fPtElec_ULS_weight->Fill(fPtE,fNonHFE->GetNULS());
3831                                                         if(fNonHFE->IsLS()) fPtElec_LS_weight->Fill(fPtE,fNonHFE->GetNLS());
3832                                                 }
3833
3834
3835                                         }//!IsTPConly
3836
3837                                         if(IsTPConly){
3838                                                 if(fNonHFE->IsULS()) fPtElec_ULS2->Fill(fPtE,fNonHFE->GetNULS());
3839                                                 if(fNonHFE->IsLS()) fPtElec_LS2->Fill(fPtE,fNonHFE->GetNLS());
3840
3841                                                         //new 08 October        //weighted histograms
3842                                                 if(TMath::Abs(fMCparticleMother->GetPdgCode())==111 || TMath::Abs(fMCparticleMother->GetPdgCode())==221){
3843                                                         Double_t mPt=fMCparticleMother->Pt();
3844
3845                                                                 Double_t mweight1=1;
3846                                                                 Double_t mweight2=1;
3847                                                                 //Double_t weight=1;
3848
3849
3850                                                          //----------------------------------------------------------------------------
3851
3852                                                         //correction based on data for d3
3853
3854                                                                 if(TMath::Abs(fMCparticleMother->GetPdgCode())==111){
3855                                                                         Double_t x=mPt;
3856                                                                         weight=CalculateWeight(111, x);
3857
3858                                                                 }
3859                                                                 if(TMath::Abs(fMCparticleMother->GetPdgCode())==221){
3860                                                                         Double_t x=mPt;
3861                                                                         weight=CalculateWeight(221, x);
3862
3863                                                                 }
3864
3865
3866                                                                 //check this
3867                                                         if(fNonHFE->IsULS()) mweight1=(fNonHFE->GetNULS())/weight;
3868                                                         if(fNonHFE->IsLS())  mweight2=(fNonHFE->GetNLS())/weight;
3869
3870                                                                 //fill histos
3871                                                         if(fNonHFE->IsULS())fPtElec_ULS2_weight->Fill(fPtE, mweight1);
3872                                                         if(fNonHFE->IsLS())fPtElec_LS2_weight->Fill(fPtE, mweight2);
3873                                                 }
3874                                                 else if(fMCparticleMother->GetMother()>0 && (TMath::Abs(fMCparticleGMother->GetPdgCode())==111 || TMath::Abs(fMCparticleGMother->GetPdgCode())==221 )){
3875                                                         Double_t gmPt=fMCparticleGMother->Pt();
3876                                                         Double_t gmweight1=1;
3877                                                         Double_t gmweight2=1;
3878                                                                 //Double_t weight=1;
3879
3880
3881                                                         //----------------------------------------------------------------------------
3882                                                          //correction based on data only for pi0
3883
3884                                                                 if(TMath::Abs(fMCparticleGMother->GetPdgCode())==111){
3885                                                                         Double_t x=gmPt;
3886                                                                         weight=CalculateWeight(111, x);
3887                                                                 }
3888
3889
3890                                                                 if(TMath::Abs(fMCparticleGMother->GetPdgCode())==221){
3891                                                                         Double_t x=gmPt;
3892                                                                         weight=CalculateWeight(221, x);
3893                                                                 }
3894
3895                                                         //----------------------------------------------------------------------------
3896
3897
3898
3899                                                                 //check this
3900                                                         if(fNonHFE->IsULS()) gmweight1=(fNonHFE->GetNULS())/weight;
3901                                                         if(fNonHFE->IsLS())  gmweight2=(fNonHFE->GetNLS())/weight;
3902
3903                                                                 //fill histos
3904                                                         if(fNonHFE->IsULS())fPtElec_ULS2_weight->Fill(fPtE, gmweight1);
3905                                                         if(fNonHFE->IsLS())fPtElec_LS2_weight->Fill(fPtE, gmweight2);
3906                                                 }
3907                                                 else{
3908                                                         if(fNonHFE->IsULS()) fPtElec_ULS2_weight->Fill(fPtE,fNonHFE->GetNULS());
3909                                                         if(fNonHFE->IsLS()) fPtElec_LS2_weight->Fill(fPtE,fNonHFE->GetNLS());
3910                                                 }
3911
3912
3913                                                 //----------------------------------------------------------------------------
3914                                                 //to check other way to calculate efficiency
3915                                                 //ULS with no weight from ULS-LS original
3916                                                 //we have to know if track2 comes from same mother!!!
3917                                                 //----------------------------------------------------------------------------
3918                                                 if(fNonHFE->IsULS()){
3919
3920                                                         for(Int_t iTracks = 0; iTracks < fVevent->GetNumberOfTracks(); iTracks++)
3921                                                         {
3922
3923                                                                 AliVParticle* Vtrack2 = fVevent->GetTrack(iTracks);
3924                                                                 if (!Vtrack2)
3925                                                                 {
3926                                                                         printf("ERROR: Could not receive track %d\n", iTracks);
3927                                                                         continue;
3928                                                                 }
3929                                                                 AliVTrack *track2 = dynamic_cast<AliVTrack*>(Vtrack2);
3930                                                                 if(track2->GetLabel()<0) continue;
3931                                                                 fMCparticle2 = (AliAODMCParticle*) fMCarray->At(track2->GetLabel());
3932                                                                 if(fMCparticle2->GetMother()<0) continue;
3933
3934                                                                 for(Int_t i = 0; i < fNonHFE->GetNULS(); i++)
3935                                                                 {
3936                                                                         if(fUlsPartner[i]==iTracks){
3937                                                                                         //only fill if it has same mother
3938                                                                                         //with weight to take into account the number of partners
3939                                                                                 if(fMCparticle2->GetMother()==fMCparticle->GetMother()) fPtElec_ULS_MC->Fill(fPtE, fNonHFE->GetNULS());
3940
3941                                                                                 //-----------------------------------------------------------------------------------------------------------
3942                                                                                 //weight for mother
3943                                                                                 //Double_t weight2=1;
3944                                                                                 Double_t mPt=fMCparticleMother->Pt();
3945
3946
3947                                                                                 if(TMath::Abs(fMCparticleMother->GetPdgCode())==111){
3948                                                                                         Double_t x=mPt;
3949                                                                                         weight=CalculateWeight(111, x);
3950
3951
3952                                                                                 }
3953                                                                                 if(TMath::Abs(fMCparticleMother->GetPdgCode())==221){
3954                                                                                         Double_t x=mPt;
3955                                                                                         weight=CalculateWeight(221, x);
3956
3957
3958                                                                                 }
3959
3960                                                                                 //weight for grandmother
3961
3962                                                                                 if((fMCparticleMother->GetMother()>0) && TMath::Abs(((fMCparticleGMother->GetPdgCode())==111))){
3963                                                                                         Double_t gmPt=fMCparticleGMother->Pt();
3964                                                                                         Double_t x=gmPt;
3965                                                                                         weight=CalculateWeight(111, x);
3966
3967                                                                                 }
3968                                                                                 if((fMCparticleMother->GetMother()>0) && TMath::Abs(((fMCparticleGMother->GetPdgCode())==221))){
3969                                                                                         Double_t gmPt=fMCparticleGMother->Pt();
3970                                                                                         Double_t x=gmPt;
3971                                                                                         weight=CalculateWeight(221, x);
3972
3973                                                                                 }
3974
3975
3976                                                                                 if(fMCparticle2->GetMother()==fMCparticle->GetMother()) fPtElec_ULS_MC_weight->Fill(fPtE, (fNonHFE->GetNULS())*1./weight);
3977
3978                                                                                 //-----------------------------------------------------------------------------------------------------------
3979                                                                                 //end of weight
3980
3981                                                                         }//partner found same as track
3982                                                                 }//loop in all partner
3983
3984                                                         }//track
3985                                                 }//is ULS
3986                                                 //----------------------------------------------------------------------------
3987                                                 //end of part to check other way to calculate efficiency
3988                                                 //----------------------------------------------------------------------------
3989
3990                                         }//IsTPConly
3991
3992                                 }//particle kind
3993
3994                                 if(IsTPConly){
3995                                                 //ULS-LS with no pid AOD
3996                                         if(fNonHFE->IsULS()) fPtElec_ULS_NoPid->Fill(fPtE,fNonHFE->GetNULS());
3997                                         if(fNonHFE->IsLS()) fPtElec_LS_NoPid->Fill(fPtE,fNonHFE->GetNLS());
3998                                 }
3999
4000                         }//close IsAOD
4001                          //It is ESD
4002                         else
4003                         {
4004                                 if(TMath::Abs(fMCtrack->GetPdgCode())==11 && (TMath::Abs(fMCtrackMother->GetPdgCode())==22 || TMath::Abs(fMCtrackMother->GetPdgCode())==111 || TMath::Abs(fMCtrackMother->GetPdgCode())==221))
4005                                 {
4006                                         if(!IsTPConly){
4007                                                 if(fNonHFE->IsULS()) fPtElec_ULS->Fill(fPtE,fNonHFE->GetNULS());
4008                                                 if(fNonHFE->IsLS()) fPtElec_LS->Fill(fPtE,fNonHFE->GetNLS());
4009                                         }
4010
4011                                         if(IsTPConly){
4012                                                 if(fNonHFE->IsULS()) fPtElec_ULS2->Fill(fPtE,fNonHFE->GetNULS());
4013                                                 if(fNonHFE->IsLS()) fPtElec_LS2->Fill(fPtE,fNonHFE->GetNLS());
4014                                         }
4015                                 }
4016
4017
4018                         }
4019                 }//close IsMC
4020                  ///_________________________________________________________________
4021                  //not MC
4022                 else
4023                 {
4024                         if(!IsTPConly){
4025                                 if(fNonHFE->IsULS()) fPtElec_ULS->Fill(fPtE,fNonHFE->GetNULS());
4026                                 if(fNonHFE->IsLS()) fPtElec_LS->Fill(fPtE,fNonHFE->GetNLS());
4027                         }
4028
4029                         if(IsTPConly){
4030                                 if(fNonHFE->IsULS()) fPtElec_ULS2->Fill(fPtE,fNonHFE->GetNULS());
4031                                 if(fNonHFE->IsLS()) fPtElec_LS2->Fill(fPtE,fNonHFE->GetNLS());
4032                         }
4033                 }
4034
4035         //for MC closure test
4036         // to be used as data, with MC input
4037         if(!IsTPConly){
4038                 if(fNonHFE->IsULS()) fPtElec_ULS_mc_closure->Fill(fPtE,fNonHFE->GetNULS());
4039                 if(fNonHFE->IsLS()) fPtElec_LS_mc_closure->Fill(fPtE,fNonHFE->GetNLS());
4040         }
4041
4042         if(IsTPConly){
4043                 if(fNonHFE->IsULS()) fPtElec_ULS2_mc_closure->Fill(fPtE,fNonHFE->GetNULS());
4044                 if(fNonHFE->IsLS()) fPtElec_LS2_mc_closure->Fill(fPtE,fNonHFE->GetNLS());
4045         }
4046
4047
4048
4049 }//end of Background function
4050
4051         //______________________________________________________________________
4052 void AliAnalysisTaskEMCalHFEpA::ElectronHadronCorrelation(AliVTrack *track, Int_t trackIndex, AliVParticle *vtrack)
4053 {
4054
4055         ///_________________________________________________________________
4056         ///MC analysis
4057         if(fIsMC)
4058         {
4059                 if(track->GetLabel() < 0)
4060         {
4061                         AliWarning(Form("The track %d does not have a valid MC label",trackIndex));
4062                         return;
4063         }
4064
4065                 if(fIsAOD)
4066                 {
4067                         fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
4068
4069                         if(fMCparticle->GetMother()<0) return;
4070
4071                 fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
4072
4073                 if(TMath::Abs(fMCparticle->GetPdgCode())==11 && (TMath::Abs(fMCparticleMother->GetPdgCode())==22 || TMath::Abs(fMCparticleMother->GetPdgCode())==111 || TMath::Abs(fMCparticleMother->GetPdgCode())==221))
4074                 {
4075                                         //Is Background
4076                                 fPtBackgroundBeforeReco->Fill(track->Pt());
4077                         }
4078                 }
4079                 else
4080                 {
4081                 fMCtrack = fMCstack->Particle(track->GetLabel());
4082
4083                 if(fMCtrack->GetFirstMother()<0) return;
4084
4085                 fMCtrackMother = fMCstack->Particle(fMCtrack->GetFirstMother());
4086
4087                 if(TMath::Abs(fMCtrack->GetPdgCode())==11 && (TMath::Abs(fMCtrackMother->GetPdgCode())==22 || TMath::Abs(fMCtrackMother->GetPdgCode())==111 || TMath::Abs(fMCtrackMother->GetPdgCode())==221))
4088                 {
4089                                         //Is Background
4090                                 fPtBackgroundBeforeReco->Fill(track->Pt());
4091                         }
4092                 }
4093         }
4094                 ///_________________________________________________________________
4095
4096                 //________________________________________________
4097                 //Associated particle cut
4098         fPartnerCuts->SetAcceptKinkDaughters(kFALSE);
4099     fPartnerCuts->SetRequireITSRefit(kTRUE);
4100     fPartnerCuts->SetRequireTPCRefit(kTRUE);
4101     fPartnerCuts->SetEtaRange(-0.9,0.9);
4102     fPartnerCuts->SetMaxChi2PerClusterTPC(4.0);
4103     fPartnerCuts->SetMinNClustersTPC(80);
4104     fPartnerCuts->SetPtRange(0.3,1e10);
4105                 //fPartnerCuts->SetRequireSigmaToVertex(kTRUE);
4106                 //fPartnerCuts->SetMaxDCAToVertexXY(1);
4107                 //fPartnerCuts->SetMaxDCAToVertexZ(3);
4108                 //_________________________________________________
4109
4110                 ///#################################################################
4111                 //Non-HFE reconstruction
4112         fNonHFE = new AliSelectNonHFE();
4113         fNonHFE->SetAODanalysis(fIsAOD);
4114         if(fMassCutFlag) fNonHFE->SetInvariantMassCut(fMassCut);
4115         if(fAngleCutFlag) fNonHFE->SetOpeningAngleCut(fAngleCut);
4116         if(fChi2CutFlag) fNonHFE->SetChi2OverNDFCut(fChi2Cut);
4117         if(fDCAcutFlag) fNonHFE->SetDCACut(fDCAcut);
4118         fNonHFE->SetAlgorithm("DCA"); //KF
4119         fNonHFE->SetPIDresponse(fPidResponse);
4120         fNonHFE->SetTrackCuts(-3.5,3.5,fPartnerCuts);
4121         fNonHFE->SetAdditionalCuts(fPtMinAsso,fTpcNclsAsso);
4122
4123
4124         fNonHFE->SetHistAngleBack(fOpAngleBack);
4125         fNonHFE->SetHistAngle(fOpAngle);
4126         fNonHFE->SetHistDCABack(fDCABack);
4127         fNonHFE->SetHistDCA(fDCA);
4128         fNonHFE->SetHistMassBack(fInvMassBack);
4129         fNonHFE->SetHistMass(fInvMass);
4130
4131         fNonHFE->FindNonHFE(trackIndex,vtrack,fVevent);
4132
4133         Int_t *fUlsPartner = fNonHFE->GetPartnersULS();
4134         Int_t *fLsPartner = fNonHFE->GetPartnersLS();
4135         Bool_t fUlsIsPartner = kFALSE;
4136         Bool_t fLsIsPartner = kFALSE;
4137                 ///#################################################################
4138
4139
4140                 //Electron Information
4141         Double_t fPhiE = -999;
4142         Double_t fEtaE = -999;
4143         Double_t fPhiH = -999;
4144         Double_t fEtaH = -999;
4145         Double_t fDphi = -999;
4146         Double_t fDeta = -999;
4147         Double_t fPtE = -999;
4148         Double_t fPtH = -999;
4149
4150         Double_t pi = TMath::Pi();
4151
4152         fPhiE = track->Phi();
4153         fEtaE = track->Eta();
4154         fPtE = track->Pt();
4155
4156
4157                 ///_________________________________________________________________
4158                 ///MC analysis
4159         if(fIsMC)
4160         {
4161                 if(fIsAOD)
4162                 {
4163                         if(TMath::Abs(fMCparticle->GetPdgCode())==11 && (TMath::Abs(fMCparticleMother->GetPdgCode())==22 || TMath::Abs(fMCparticleMother->GetPdgCode())==111 || TMath::Abs(fMCparticleMother->GetPdgCode())==221))
4164                 {
4165                                 if(fNonHFE->IsULS()) fPtElec_ULS->Fill(fPtE,fNonHFE->GetNULS());
4166                                 if(fNonHFE->IsLS()) fPtElec_LS->Fill(fPtE,fNonHFE->GetNLS());
4167
4168
4169                         }
4170
4171                         if(fNonHFE->IsULS()) fPtElec_ULS_NoPid->Fill(fPtE,fNonHFE->GetNULS());
4172                         if(fNonHFE->IsLS()) fPtElec_LS_NoPid->Fill(fPtE,fNonHFE->GetNLS());
4173                 }
4174                 else
4175                 {
4176                         if(TMath::Abs(fMCtrack->GetPdgCode())==11 && (TMath::Abs(fMCtrackMother->GetPdgCode())==22 || TMath::Abs(fMCtrackMother->GetPdgCode())==111 || TMath::Abs(fMCtrackMother->GetPdgCode())==221))
4177                         {
4178                                 if(fNonHFE->IsULS()) fPtElec_ULS->Fill(fPtE,fNonHFE->GetNULS());
4179                                 if(fNonHFE->IsLS()) fPtElec_LS->Fill(fPtE,fNonHFE->GetNLS());
4180
4181                         }
4182
4183
4184                 }
4185         }
4186                 ///_________________________________________________________________
4187         else
4188         {
4189                 if(fNonHFE->IsULS()) fPtElec_ULS->Fill(fPtE,fNonHFE->GetNULS());
4190                 if(fNonHFE->IsLS()) fPtElec_LS->Fill(fPtE,fNonHFE->GetNLS());
4191         }
4192
4193
4194
4195
4196                 //__________________________________________________________________
4197                 //Event Mixing Analysis - Hadron Loop
4198                 //Retrieve
4199         if(fEventMixingFlag)
4200         {
4201                 fPool = fPoolMgr->GetEventPool(fCentrality->GetCentralityPercentile("V0A"), fZvtx); // Get the buffer associated with the current centrality and z-vtx
4202
4203                 if(!fPool) AliFatal(Form("No pool found for centrality = %f, zVtx = %f",fCentrality->GetCentralityPercentile("V0A"), fZvtx));
4204
4205                 if(fPool->GetCurrentNEvents() >= 5) // start mixing when 5 events are in the buffer
4206                 {
4207                         fPoolNevents->Fill(fPool->GetCurrentNEvents());
4208
4209                         for (Int_t jMix = 0; jMix < fPool->GetCurrentNEvents(); jMix++)  // mix with each event in the buffer
4210                         {
4211                                 TObjArray* bgTracks = fPool->GetEvent(jMix);
4212
4213                                 for (Int_t kMix = 0; kMix < bgTracks->GetEntriesFast(); kMix++)  // mix with each track in the event
4214                                 {
4215                                         const AliEHCParticle* MixedTrack(dynamic_cast<AliEHCParticle*>(bgTracks->At(kMix)));
4216                                         if (NULL == MixedTrack) continue;
4217
4218                                         fPhiH = MixedTrack->Phi();
4219                                         fEtaH = MixedTrack->Eta();
4220                                         fPtH = MixedTrack->Pt();
4221
4222                                         if(fPtH<fAssHadronPtMin || fPtH>fAssHadronPtMax) continue;
4223
4224                                         fDphi = fPhiE - fPhiH;
4225
4226                                         if (fDphi > 3*pi/2) fDphi = fDphi - 2*pi;
4227                                         if (fDphi < -pi/2)  fDphi = fDphi + 2*pi;
4228
4229                                         fDeta = fEtaE - fEtaH;
4230
4231                                         Double_t fPtBin[7] = {1,2,4,6,8,10,15};
4232
4233                                         for(Int_t i = 0; i < 6; i++)
4234                                         {
4235                                             if(fPtE>=fPtBin[i] && fPtE<fPtBin[i+1])
4236                                             {
4237                                                         fCEtaPhi_Inc_EM[i]->Fill(fDphi,fDeta);
4238
4239                                                         if(fNonHFE->IsULS()) fCEtaPhi_ULS_EM[i]->Fill(fDphi,fDeta);
4240                                                         if(fNonHFE->IsLS()) fCEtaPhi_LS_EM[i]->Fill(fDphi,fDeta);
4241
4242                                                         if(fNonHFE->IsULS()) fCEtaPhi_ULS_Weight_EM[i]->Fill(fDphi,fDeta,fNonHFE->GetNULS());
4243                                                         if(fNonHFE->IsLS()) fCEtaPhi_LS_Weight_EM[i]->Fill(fDphi,fDeta,fNonHFE->GetNLS());
4244                                             }
4245                                         }
4246
4247                                                 // TODO your code: do event mixing with current event and bgTracks
4248                                                 // note that usually the content filled now is weighted by 1 / pool->GetCurrentNEvents()
4249                                 }
4250                         }
4251                 }
4252         }
4253                 //__________________________________________________________________
4254
4255                 //__________________________________________________________________
4256                 //Same Event Analysis - Hadron Loop
4257         for(Int_t iTracks = 0; iTracks < fVevent->GetNumberOfTracks(); iTracks++)
4258         {
4259                 if(trackIndex==iTracks) continue;
4260
4261                 AliVParticle* Vtrack2 = fVevent->GetTrack(iTracks);
4262                 if (!Vtrack2)
4263                 {
4264                         printf("ERROR: Could not receive track %d\n", iTracks);
4265                         continue;
4266                 }
4267
4268                 AliVTrack *track2 = dynamic_cast<AliVTrack*>(Vtrack2);
4269
4270                 if(track2->Eta()<fEtaCutMin || track2->Eta()>fEtaCutMax ) continue;
4271
4272                 if(fIsAOD)
4273                 {
4274                         AliAODTrack *atrack2 = dynamic_cast<AliAODTrack*>(Vtrack2);
4275                         if(!atrack2->TestFilterMask(AliAODTrack::kTrkTPCOnly)) continue;
4276                         if((!(atrack2->GetStatus()&AliESDtrack::kITSrefit)|| (!(atrack2->GetStatus()&AliESDtrack::kTPCrefit)))) continue;
4277                         if(atrack2->GetTPCNcls() < 80) continue;
4278                         if(fAssocWithSPD && ((!(atrack2->HasPointOnITSLayer(0))) && (!(atrack2->HasPointOnITSLayer(1))))) continue;
4279                 }
4280                 else
4281                 {
4282                         AliESDtrack *etrack2 = dynamic_cast<AliESDtrack*>(Vtrack2);
4283                         if(!fPartnerCuts->AcceptTrack(etrack2)) continue;
4284                 }
4285
4286                 fPhiH = track2->Phi();
4287                 fEtaH = track2->Eta();
4288                 fPtH = track2->Pt();
4289
4290                 if(fPtH<fAssHadronPtMin || fPtH>fAssHadronPtMax) continue;
4291
4292                 fDphi = fPhiE - fPhiH;
4293
4294                 if (fDphi > 3*pi/2) fDphi = fDphi - 2*pi;
4295                 if (fDphi < -pi/2)  fDphi = fDphi + 2*pi;
4296
4297                 fDeta = fEtaE - fEtaH;
4298
4299                 Double_t fPtBin[7] = {1,2,4,6,8,10,15};
4300
4301                 //______________________________________________________________
4302                 //Check if this track is a Non-HFE partner
4303                 fUlsIsPartner = kFALSE;
4304                 fLsIsPartner = kFALSE;
4305                 for(Int_t i = 0; i < fNonHFE->GetNULS(); i++)
4306                 {
4307                         if(fUlsPartner[i]==iTracks) fUlsIsPartner=kTRUE;
4308                 }
4309                 for(Int_t i = 0; i < fNonHFE->GetNLS(); i++)
4310                 {
4311                         if(fLsPartner[i]==iTracks) fLsIsPartner=kTRUE;
4312                 }
4313                 //______________________________________________________________
4314
4315                 for(Int_t i = 0; i < 6; i++)
4316                 {
4317                     if(fPtE>=fPtBin[i] && fPtE<fPtBin[i+1])
4318                     {
4319                                 fCEtaPhi_Inc[i]->Fill(fDphi,fDeta);
4320
4321                                 if(fNonHFE->IsULS()) fCEtaPhi_ULS[i]->Fill(fDphi,fDeta);
4322                                 if(fNonHFE->IsLS()) fCEtaPhi_LS[i]->Fill(fDphi,fDeta);
4323                                 if(fNonHFE->IsULS() && !fUlsIsPartner) fCEtaPhi_ULS_NoP[i]->Fill(fDphi,fDeta);
4324                                 if(fNonHFE->IsLS() && !fLsIsPartner) fCEtaPhi_LS_NoP[i]->Fill(fDphi,fDeta);
4325
4326                                 if(fNonHFE->IsULS()) fCEtaPhi_ULS_Weight[i]->Fill(fDphi,fDeta,fNonHFE->GetNULS());
4327                                 if(fNonHFE->IsLS()) fCEtaPhi_LS_Weight[i]->Fill(fDphi,fDeta,fNonHFE->GetNLS());
4328                                 if(fNonHFE->IsULS() && !fUlsIsPartner) fCEtaPhi_ULS_NoP_Weight[i]->Fill(fDphi,fDeta,fNonHFE->GetNULS());
4329                                 if(fNonHFE->IsLS() && !fLsIsPartner) fCEtaPhi_LS_NoP_Weight[i]->Fill(fDphi,fDeta,fNonHFE->GetNLS());
4330                     }
4331                 }
4332         }
4333 }
4334
4335         //____________________________________________________________________________________________________________
4336         //Create a TObjArray with selected hadrons, for the mixed event analysis
4337 TObjArray* AliAnalysisTaskEMCalHFEpA::SelectedHadrons()
4338 {
4339         fTracksClone = new TObjArray;
4340         fTracksClone->SetOwner(kTRUE);
4341
4342         for(Int_t iTracks = 0; iTracks < fVevent->GetNumberOfTracks(); iTracks++)
4343         {
4344                 AliVParticle* Vtrack2 = fVevent->GetTrack(iTracks);
4345                 if (!Vtrack2)
4346                 {
4347                         printf("ERROR: Could not receive track %d\n", iTracks);
4348                         continue;
4349                 }
4350
4351                 AliVTrack *track2 = dynamic_cast<AliVTrack*>(Vtrack2);
4352
4353                 if(track2->Eta()<fEtaCutMin || track2->Eta()>fEtaCutMax ) continue;
4354
4355                 if(fIsAOD)
4356                 {
4357                         AliAODTrack *atrack2 = dynamic_cast<AliAODTrack*>(Vtrack2);
4358                         if(!atrack2->TestFilterMask(AliAODTrack::kTrkTPCOnly)) continue;
4359                         if((!(atrack2->GetStatus()&AliESDtrack::kITSrefit)|| (!(atrack2->GetStatus()&AliESDtrack::kTPCrefit)))) continue;
4360                         if(atrack2->GetTPCNcls() < 80) continue;
4361                         if(fAssocWithSPD && ((!(atrack2->HasPointOnITSLayer(0))) && (!(atrack2->HasPointOnITSLayer(1))))) continue;
4362                 }
4363                 else
4364                 {
4365                         AliESDtrack *etrack2 = dynamic_cast<AliESDtrack*>(Vtrack2);
4366                         if(!fPartnerCuts->AcceptTrack(etrack2)) continue;
4367                 }
4368
4369                 fTracksClone->Add(new AliEHCParticle(track2->Eta(), track2->Phi(), track2->Pt()));
4370         }
4371         return fTracksClone;
4372 }
4373         //____________________________________________________________________________________________________________
4374
4375         //______________________________________________________________________
4376 void AliAnalysisTaskEMCalHFEpA::DiHadronCorrelation(AliVTrack *track, Int_t trackIndex)
4377 {
4378                 //________________________________________________
4379                 //Associated particle cut
4380         fPartnerCuts->SetAcceptKinkDaughters(kFALSE);
4381     fPartnerCuts->SetRequireITSRefit(kTRUE);
4382     fPartnerCuts->SetRequireTPCRefit(kTRUE);
4383     fPartnerCuts->SetEtaRange(-0.9,0.9);
4384     fPartnerCuts->SetMaxChi2PerClusterTPC(4.0);
4385     fPartnerCuts->SetMinNClustersTPC(80);
4386     fPartnerCuts->SetPtRange(0.3,1e10);
4387                 //fPartnerCuts->SetRequireSigmaToVertex(kTRUE);
4388                 //fPartnerCuts->SetMaxDCAToVertexXY(1);
4389                 //fPartnerCuts->SetMaxDCAToVertexZ(3);
4390                 //_________________________________________________
4391
4392                 //Electron Information
4393         Double_t fPhiE = -999;
4394         Double_t fEtaE = -999;
4395         Double_t fPhiH = -999;
4396         Double_t fEtaH = -999;
4397         Double_t fDphi = -999;
4398         Double_t fDeta = -999;
4399         Double_t fPtE = -999;
4400         Double_t fPtH = -999;
4401
4402         Double_t pi = TMath::Pi();
4403
4404         fPhiE = track->Phi();
4405         fEtaE = track->Eta();
4406         fPtE = track->Pt();
4407
4408                 //__________________________________________________________________
4409                 //Same Event Analysis - Hadron Loop
4410         for(Int_t iTracks = 0; iTracks < fVevent->GetNumberOfTracks(); iTracks++)
4411         {
4412                 if(trackIndex==iTracks) continue;
4413
4414                 AliVParticle* Vtrack2 = fVevent->GetTrack(iTracks);
4415                 if (!Vtrack2)
4416                 {
4417                         printf("ERROR: Could not receive track %d\n", iTracks);
4418                         continue;
4419                 }
4420
4421                 AliVTrack *track2 = dynamic_cast<AliVTrack*>(Vtrack2);
4422                 if(track2->Eta()<fEtaCutMin || track2->Eta()>fEtaCutMax ) continue;
4423
4424                 if(fIsAOD)
4425                 {
4426                         AliAODTrack *atrack2 = dynamic_cast<AliAODTrack*>(Vtrack2);
4427                         if(!atrack2->TestFilterMask(AliAODTrack::kTrkTPCOnly)) continue;
4428                         if((!(atrack2->GetStatus()&AliESDtrack::kITSrefit)|| (!(atrack2->GetStatus()&AliESDtrack::kTPCrefit)))) continue;
4429                         if(atrack2->GetTPCNcls() < 80) continue;
4430                         if(fAssocWithSPD && ((!(atrack2->HasPointOnITSLayer(0))) && (!(atrack2->HasPointOnITSLayer(1))))) continue;
4431                 }
4432                 else
4433                 {
4434                         AliESDtrack *etrack2 = dynamic_cast<AliESDtrack*>(Vtrack2);
4435                         if(!fPartnerCuts->AcceptTrack(etrack2)) continue;
4436                 }
4437
4438                 fPhiH = track2->Phi();
4439                 fEtaH = track2->Eta();
4440                 fPtH = track2->Pt();
4441
4442                 if(fPtH<fAssHadronPtMin || fPtH>fAssHadronPtMax) continue;
4443
4444                 fDphi = fPhiE - fPhiH;
4445
4446                 if (fDphi > 3*pi/2) fDphi = fDphi - 2*pi;
4447                 if (fDphi < -pi/2)  fDphi = fDphi + 2*pi;
4448
4449                 fDeta = fEtaE - fEtaH;
4450
4451                 Double_t fPtBin[7] = {1,2,4,6,8,10,15};
4452
4453                 for(Int_t i = 0; i < 6; i++)
4454                 {
4455                     if(fPtE>=fPtBin[i] && fPtE<fPtBin[i+1])
4456                     {
4457                                 fCEtaPhi_Inc_DiHadron[i]->Fill(fDphi,fDeta);
4458                     }
4459                 }
4460         }
4461 }
4462         //____________________________________________________________________________________________________________
4463
4464         //______________________________________________________________________
4465 Bool_t AliAnalysisTaskEMCalHFEpA::FindMother(Int_t mcIndex)
4466 {
4467         fIsHFE1 = kFALSE;
4468         fIsHFE2 = kFALSE;
4469         fIsNonHFE = kFALSE;
4470         fIsFromD = kFALSE;
4471         fIsFromB = kFALSE;
4472         fIsFromPi0 = kFALSE;
4473         fIsFromEta = kFALSE;
4474         fIsFromGamma = kFALSE;
4475
4476         if(mcIndex < 0 || !fIsMC)
4477         {
4478                 return kFALSE;
4479         }
4480
4481         Int_t pdg = -99999;
4482         Int_t mpdg = -99999;
4483         Int_t gmpdg = -99999;
4484         Int_t ggmpdg = -99999;
4485         Int_t gggmpdg = -99999;
4486
4487         if(fIsAOD)
4488         {
4489                 fMCparticle = (AliAODMCParticle*) fMCarray->At(mcIndex);
4490
4491                 pdg = TMath::Abs(fMCparticle->GetPdgCode());
4492
4493
4494                 if(pdg!=11)
4495                 {
4496                         fIsHFE1 = kFALSE;
4497                         fIsHFE2 = kFALSE;
4498                         fIsNonHFE = kFALSE;
4499                         fIsFromD = kFALSE;
4500                         fIsFromB = kFALSE;
4501                         fIsFromPi0 = kFALSE;
4502                         fIsFromEta = kFALSE;
4503                         fIsFromGamma = kFALSE;
4504                         return kFALSE;
4505                 }
4506
4507                 if(fMCparticle->GetMother()<0)
4508                 {
4509                         fIsHFE1 = kFALSE;
4510                         fIsHFE2 = kFALSE;
4511                         fIsNonHFE = kFALSE;
4512                         fIsFromD = kFALSE;
4513                         fIsFromB = kFALSE;
4514                         fIsFromPi0 = kFALSE;
4515                         fIsFromEta = kFALSE;
4516                         fIsFromGamma = kFALSE;
4517                         return kFALSE;
4518                 }
4519
4520                 fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
4521                 mpdg = TMath::Abs(fMCparticleMother->GetPdgCode());
4522
4523                 if(fMCparticleMother->GetMother()<0)
4524                 {
4525                         gmpdg = 0;
4526                         ggmpdg = 0;
4527                         gggmpdg = 0;
4528                 }
4529                 else
4530                 {
4531                         fMCparticleGMother = (AliAODMCParticle*) fMCarray->At(fMCparticleMother->GetMother());
4532                         gmpdg = TMath::Abs(fMCparticleGMother->GetPdgCode());
4533                         if(fMCparticleGMother->GetMother()<0)
4534                         {
4535                                 ggmpdg = 0;
4536                                 gggmpdg = 0;
4537                         }
4538                         else
4539                         {
4540                                 fMCparticleGGMother = (AliAODMCParticle*) fMCarray->At(fMCparticleGMother->GetMother());
4541                                 ggmpdg = TMath::Abs(fMCparticleGGMother->GetPdgCode());
4542                                 if(fMCparticleGGMother->GetMother()<0)
4543                                 {
4544                                         gggmpdg = 0;
4545                                 }
4546                                 else
4547                                 {
4548                                         fMCparticleGGGMother = (AliAODMCParticle*) fMCarray->At(fMCparticleGGMother->GetMother());
4549                                         gggmpdg = TMath::Abs(fMCparticleGGGMother->GetPdgCode());
4550                                 }
4551                         }
4552                 }
4553         }
4554         else
4555         {
4556                 fMCtrack = fMCstack->Particle(mcIndex);
4557
4558                 pdg = TMath::Abs(fMCtrack->GetPdgCode());
4559
4560                 if(pdg!=11)
4561                 {
4562                         fIsHFE1 = kFALSE;
4563                         fIsHFE2 = kFALSE;
4564                         fIsNonHFE = kFALSE;
4565                         fIsFromD = kFALSE;
4566                         fIsFromB = kFALSE;
4567                         fIsFromPi0 = kFALSE;
4568                         fIsFromEta = kFALSE;
4569                         fIsFromGamma = kFALSE;
4570                         return kFALSE;
4571                 }
4572
4573                 if(fMCtrack->GetFirstMother()<0)
4574                 {
4575                         fIsHFE1 = kFALSE;
4576                         fIsHFE2 = kFALSE;
4577                         fIsNonHFE = kFALSE;
4578                         fIsFromD = kFALSE;
4579                         fIsFromB = kFALSE;
4580                         fIsFromPi0 = kFALSE;
4581                         fIsFromEta = kFALSE;
4582                         fIsFromGamma = kFALSE;
4583                         return kFALSE;
4584                 }
4585
4586                 fMCtrackMother = fMCstack->Particle(fMCtrack->GetFirstMother());
4587                 mpdg = TMath::Abs(fMCtrackMother->GetPdgCode());
4588
4589                 if(fMCtrackMother->GetFirstMother()<0)
4590                 {
4591                         gmpdg = 0;
4592                         ggmpdg = 0;
4593                         gggmpdg = 0;
4594                 }
4595                 else
4596                 {
4597                         fMCtrackGMother = fMCstack->Particle(fMCtrackMother->GetFirstMother());
4598                         gmpdg = TMath::Abs(fMCtrackGMother->GetPdgCode());
4599
4600                         if(fMCtrackGMother->GetFirstMother()<0)
4601                         {
4602                                 ggmpdg = 0;
4603                                 gggmpdg = 0;
4604                         }
4605                         else
4606                         {
4607                                 fMCtrackGGMother = fMCstack->Particle(fMCtrackGMother->GetFirstMother());
4608                                 ggmpdg = TMath::Abs(fMCtrackGGMother->GetPdgCode());
4609
4610                                 if(fMCtrackGGMother->GetFirstMother()<0)
4611                                 {
4612                                         gggmpdg = 0;
4613                                 }
4614                                 else
4615                                 {
4616                                         fMCtrackGGGMother = fMCstack->Particle(fMCtrackGGMother->GetFirstMother());
4617                                         gggmpdg = TMath::Abs(fMCtrackGGGMother->GetPdgCode());
4618                                 }
4619                         }
4620                 }
4621         }
4622
4623         //Tag Electron Source
4624         if(mpdg==111 || mpdg==221 || mpdg==22)
4625         {
4626                 fIsHFE1 = kFALSE;
4627                 fIsHFE2 = kFALSE;
4628                 fIsNonHFE = kTRUE;
4629                 fIsFromD = kFALSE;
4630                 fIsFromB = kFALSE;
4631
4632                 fIsFromPi0 = kFALSE;
4633                 fIsFromEta = kFALSE;
4634                 fIsFromGamma = kFALSE;
4635
4636                 if(mpdg==111) fIsFromPi0 = kFALSE;
4637                 if(mpdg==221)fIsFromEta = kFALSE;
4638                 if(mpdg==22) fIsFromGamma = kFALSE;
4639
4640                 return kTRUE;
4641         }
4642         else
4643         {
4644                 fIsHFE1 = kFALSE;
4645                 fIsHFE2 = kTRUE;
4646
4647                 fIsFromPi0 = kFALSE;
4648                 fIsFromEta = kFALSE;
4649                 fIsFromGamma = kFALSE;
4650
4651                 fIsNonHFE = kFALSE;
4652
4653                 fIsFromD = kFALSE;
4654                 fIsFromB = kFALSE;
4655
4656                 if(mpdg>400 && mpdg<500)
4657                 {
4658                         if((gmpdg>500 && gmpdg<600) || (ggmpdg>500 && ggmpdg<600) || (gggmpdg>500 && gggmpdg<600))
4659                         {
4660                                 fIsHFE1 = kTRUE;
4661                                 fIsFromD = kFALSE;
4662                                 fIsFromB = kTRUE;
4663                                 return kTRUE;
4664                         }
4665                         else
4666                         {
4667                                 fIsHFE1 = kTRUE;
4668                                 fIsFromD = kTRUE;
4669                                 fIsFromB = kFALSE;
4670                                 return kTRUE;
4671                         }
4672                 }
4673                 else if(mpdg>500 && mpdg<600)
4674                 {
4675                         fIsHFE1 = kTRUE;
4676                         fIsFromD = kFALSE;
4677                         fIsFromB = kTRUE;
4678                         return kTRUE;
4679                 }
4680                 else
4681                 {
4682                         fIsHFE1 = kFALSE;
4683                         fIsFromD = kFALSE;
4684                         fIsFromB = kFALSE;
4685                         return kFALSE;
4686                 }
4687         }
4688 }
4689 /*
4690 Bool_t AliAnalysisTaskEMCalHFEpA::ContainsBadChannel(TString calorimeter,UShort_t* cellList, Int_t nCells)
4691 {
4692                 // Check that in the cluster cells, there is no bad channel of those stored
4693                 // in fEMCALBadChannelMap
4694                 // adapted from AliCalorimeterUtils
4695
4696                 //if (!fRemoveBadChannels) return kFALSE;
4697                 //printf("fEMCALBadChannelMap %p, fPHOSBadChannelMap %p \n",fEMCALBadChannelMap,fPHOSBadChannelMap);
4698         if( !fEMCALRecoUtils->GetEMCALChannelStatusMap(0)) return kFALSE;
4699
4700                 //Int_t icol = -1;
4701                 //Int_t irow = -1;
4702                 //Int_t imod = -1;
4703         for(Int_t iCell = 0; iCell<nCells; iCell++){
4704
4705                         //Get the column and row
4706                 if(calorimeter == "EMCAL"){
4707                         return fEMCALRecoUtils->ClusterContainsBadChannel((AliEMCALGeometry*)fEMCALGeo,cellList,nCells);
4708                 }
4709                 else return kFALSE;
4710
4711         }// cell cluster loop
4712
4713         return kFALSE;
4714
4715 }
4716 */
4717 /*
4718
4719 //________________________________________________________________________________
4720 TArrayI AliAnalysisTaskEMCalHFEpA::GetTriggerPatches(Bool_t IsEventEMCALL0, Bool_t IsEventEMCALL1)
4721 {
4722         // Select the patches that triggered
4723         // Depend on L0 or L1
4724
4725         // some variables
4726         //Int_t  trigtimes[30], globCol, globRow,ntimes, i;
4727         Int_t   globCol, globRow;
4728
4729         Int_t  absId  = -1; //[100];
4730         Int_t  nPatch = 0;
4731
4732         TArrayI patches(0);
4733
4734                 // get object pointer
4735         AliVCaloTrigger *caloTrigger = InputEvent()->GetCaloTrigger( "EMCAL" );
4736
4737                 // class is not empty
4738         if( caloTrigger->GetEntries() > 0 )
4739         {
4740                         // must reset before usage, or the class will fail
4741                 caloTrigger->Reset();
4742
4743                         // go throuth the trigger channels
4744                 while( caloTrigger->Next() )
4745                 {
4746                                 // get position in global 2x2 tower coordinates
4747                         caloTrigger->GetPosition( globCol, globRow );
4748
4749                         //L0
4750                         if(IsEventEMCALL0)
4751                         {
4752                             //not implemented
4753                         }
4754
4755
4756                         else if(IsEventEMCALL1) // L1
4757                         {
4758                                         Int_t bit = 0;
4759                                         caloTrigger->GetTriggerBits(bit);
4760
4761                                         Bool_t isEGA = ((bit >> fBitEGA) & 0x1);
4762                                         //Bool_t isEJE = ((bit >> fBitEJE) & 0x1) && IsEventEMCALL1Jet  () ;
4763
4764                                         if(!isEGA) continue;
4765
4766                                         Int_t patchsize = -1;
4767                                         if(isEGA) patchsize =  2;
4768                                         //else if (isEJE) patchsize = 16;
4769
4770                                         // add 2x2 (EGA) or 16x16 (EJE) patches
4771                                         for(Int_t irow=0; irow < patchsize; irow++)
4772                                         {
4773                                                 for(Int_t icol=0; icol < patchsize; icol++)
4774                                                 {
4775                                                         GetCaloUtils()->GetEMCALGeometry()->GetAbsFastORIndexFromPositionInEMCAL(globCol+icol,globRow+irow, absId);
4776                                                         //printf("pass the time cut globCol %d, globRow %d absId %d\n",globCol,globRow, absIDTrig[nTrig]);
4777                                                         patches.Set(nPatch+1); // Set size of this array to nPatch+1 ints.
4778                                                         patches.AddAt(absId,nPatch++); //Add Int_t absId at position nPatch++
4779                                                 }
4780                                         }
4781
4782                         } // L1
4783
4784                 } // trigger iterator
4785         } // go thorough triggers
4786
4787         printf("N patches %d, test %d,first %d, last %d\n",patches.GetSize(), nPatch, patches.At(0), patches.At(patches.GetSize()-1));
4788
4789         return patches;
4790 }
4791  */
4792 Double_t AliAnalysisTaskEMCalHFEpA::CalculateWeight(Int_t pdg_particle, Double_t x)
4793 {
4794                 //weight for d3 based on MinJung parametrization //sent by Jan
4795             Double_t weight=1;
4796                 if(pdg_particle==111){
4797                         if(x>= 0.100000 &&  x < 0.112797 ) weight=1.262120;
4798                         if(x>= 0.112797 &&  x < 0.127231 ) weight=1.277765;
4799                         if(x>= 0.127231 &&  x < 0.143512 ) weight=1.295605;
4800                         if(x>= 0.143512 &&  x < 0.161877 ) weight=1.318155;
4801                         if(x>= 0.161877 &&  x < 0.182592 ) weight=1.348693;
4802                         if(x>= 0.182592 &&  x < 0.205957 ) weight=1.388636;
4803                         if(x>= 0.205957 &&  x < 0.232313 ) weight=1.439122;
4804                         if(x>= 0.232313 &&  x < 0.262041 ) weight=1.497452;
4805                         if(x>= 0.262041 &&  x < 0.295573 ) weight=1.559409;
4806                         if(x>= 0.295573 &&  x < 0.333397 ) weight=1.615169;
4807                         if(x>= 0.333397 &&  x < 0.376060 ) weight=1.654954;
4808                         if(x>= 0.376060 &&  x < 0.424183 ) weight=1.668753;
4809                         if(x>= 0.424183 &&  x < 0.478465 ) weight=1.652225;
4810                         if(x>= 0.478465 &&  x < 0.539692 ) weight=1.603119;
4811                         if(x>= 0.539692 &&  x < 0.608754 ) weight=1.526049;
4812                         if(x>= 0.608754 &&  x < 0.686654 ) weight=1.426724;
4813                         if(x>= 0.686654 &&  x < 0.774523 ) weight=1.312684;
4814                         if(x>= 0.774523 &&  x < 0.873636 ) weight=1.195395;
4815                         if(x>= 0.873636 &&  x < 0.985432 ) weight=1.086264;
4816                         if(x>= 0.985432 &&  x < 1.111534 ) weight=0.993666;
4817                         if(x>= 1.111534 &&  x < 1.253773 ) weight=0.922587;
4818                         if(x>= 1.253773 &&  x < 1.414214 ) weight=0.875739;
4819                         if(x>= 1.414214 &&  x < 1.595185 ) weight=0.852181;
4820                         if(x>= 1.595185 &&  x < 1.799315 ) weight=0.847828;
4821                         if(x>= 1.799315 &&  x < 2.029567 ) weight=0.863875;
4822                         if(x>= 2.029567 &&  x < 2.289283 ) weight=0.899112;
4823                         if(x>= 2.289283 &&  x < 2.582235 ) weight=0.955194;
4824                         if(x>= 2.582235 &&  x < 2.912674 ) weight=1.033824;
4825                         if(x>= 2.912674 &&  x < 3.285398 ) weight=1.133714;
4826                         if(x>= 3.285398 &&  x < 3.705818 ) weight=1.259471;
4827                         if(x>= 3.705818 &&  x < 4.180038 ) weight=1.406883;
4828                         if(x>= 4.180038 &&  x < 4.714942 ) weight=1.578923;
4829                         if(x>= 4.714942 &&  x < 5.318296 ) weight=1.778513;
4830                         if(x>= 5.318296 &&  x < 5.998859 ) weight=2.001171;
4831                         if(x>= 5.998859 &&  x < 6.766511 ) weight=2.223161;
4832                         if(x>= 6.766511 &&  x < 7.632396 ) weight=2.449445;
4833                         if(x>= 7.632396 &&  x < 8.609086 ) weight=2.661734;
4834                         if(x>= 8.609086 &&  x < 9.710759 ) weight=2.851935;
4835                         if(x>= 9.710759 &&  x < 10.953409 ) weight=2.974319;
4836                         if(x>= 10.953409 &&  x < 12.355077 ) weight=3.106314;
4837                         if(x>= 12.355077 &&  x < 13.936111 ) weight=3.126815;
4838                         if(x>= 13.936111 &&  x < 15.719464 ) weight=3.150053;
4839                         if(x>= 15.719464 &&  x < 17.731026 ) weight=3.218509;
4840                         if(x>= 17.731026 &&  x < 20.000000 ) weight=3.252141;
4841
4842                 }
4843             else if(pdg_particle==221){
4844                         if(x>= 0.100000 &&  x < 0.112797 ) weight=2.159358;
4845                         if(x>= 0.112797 &&  x < 0.127231 ) weight=2.145546;
4846                         if(x>= 0.127231 &&  x < 0.143512 ) weight=2.132390;
4847                         if(x>= 0.143512 &&  x < 0.161877 ) weight=2.109918;
4848                         if(x>= 0.161877 &&  x < 0.182592 ) weight=2.084920;
4849                         if(x>= 0.182592 &&  x < 0.205957 ) weight=2.054302;
4850                         if(x>= 0.205957 &&  x < 0.232313 ) weight=2.015202;
4851                         if(x>= 0.232313 &&  x < 0.262041 ) weight=1.966068;
4852                         if(x>= 0.262041 &&  x < 0.295573 ) weight=1.912255;
4853                         if(x>= 0.295573 &&  x < 0.333397 ) weight=1.844087;
4854                         if(x>= 0.333397 &&  x < 0.376060 ) weight=1.767913;
4855                         if(x>= 0.376060 &&  x < 0.424183 ) weight=1.680366;
4856                         if(x>= 0.424183 &&  x < 0.478465 ) weight=1.583468;
4857                         if(x>= 0.478465 &&  x < 0.539692 ) weight=1.475131;
4858                         if(x>= 0.539692 &&  x < 0.608754 ) weight=1.361436;
4859                         if(x>= 0.608754 &&  x < 0.686654 ) weight=1.244388;
4860                         if(x>= 0.686654 &&  x < 0.774523 ) weight=1.125973;
4861                         if(x>= 0.774523 &&  x < 0.873636 ) weight=1.015769;
4862                         if(x>= 0.873636 &&  x < 0.985432 ) weight=0.914579;
4863                         if(x>= 0.985432 &&  x < 1.111534 ) weight=0.830529;
4864                         if(x>= 1.111534 &&  x < 1.253773 ) weight=0.766397;
4865                         if(x>= 1.253773 &&  x < 1.414214 ) weight=0.723663;
4866                         if(x>= 1.414214 &&  x < 1.595185 ) weight=0.701236;
4867                         if(x>= 1.595185 &&  x < 1.799315 ) weight=0.695605;
4868                         if(x>= 1.799315 &&  x < 2.029567 ) weight=0.707578;
4869                         if(x>= 2.029567 &&  x < 2.289283 ) weight=0.735194;
4870                         if(x>= 2.289283 &&  x < 2.582235 ) weight=0.781052;
4871                         if(x>= 2.582235 &&  x < 2.912674 ) weight=0.842350;
4872                         if(x>= 2.912674 &&  x < 3.285398 ) weight=0.923676;
4873                         if(x>= 3.285398 &&  x < 3.705818 ) weight=1.028317;
4874                         if(x>= 3.705818 &&  x < 4.180038 ) weight=1.154029;
4875                         if(x>= 4.180038 &&  x < 4.714942 ) weight=1.296915;
4876                         if(x>= 4.714942 &&  x < 5.318296 ) weight=1.463674;
4877                         if(x>= 5.318296 &&  x < 5.998859 ) weight=1.651985;
4878                         if(x>= 5.998859 &&  x < 6.766511 ) weight=1.847912;
4879                         if(x>= 6.766511 &&  x < 7.632396 ) weight=2.066284;
4880                         if(x>= 7.632396 &&  x < 8.609086 ) weight=2.202231;
4881                         if(x>= 8.609086 &&  x < 9.710759 ) weight=2.399942;
4882                         if(x>= 9.710759 &&  x < 10.953409 ) weight=2.555106;
4883                         if(x>= 10.953409 &&  x < 12.355077 ) weight=2.632377;
4884                         if(x>= 12.355077 &&  x < 13.936111 ) weight=2.609660;
4885                         if(x>= 13.936111 &&  x < 15.719464 ) weight=2.656343;
4886                         if(x>= 15.719464 &&  x < 17.731026 ) weight=2.615438;
4887                         if(x>= 17.731026 &&  x < 20.000000 ) weight=2.576269;
4888
4889                 }
4890             else weight=1;
4891
4892             return weight;
4893
4894 }
4895 Double_t AliAnalysisTaskEMCalHFEpA::SetEoverPCutPtDependentMC(Double_t pt)
4896 {
4897
4898         fEoverPCutMin=0.8;
4899         fEoverPCutMax=1.2;
4900
4901         /*
4902         ////======================================================================================
4903         ////======================================================================================
4904         ////============================= data 2 sigma ===========================================
4905         if(!fIsMC){
4906                         //data 2 sigma
4907                 if(pt>= 2.00 &&  pt < 2.50){
4908                         fEoverPCutMin=0.7719;
4909                         fEoverPCutMax=1.1023;
4910                 }
4911                 if(pt>= 2.50 &&  pt < 3.00){
4912                         fEoverPCutMin=0.7966;
4913                         fEoverPCutMax=1.1088;
4914                 }
4915                 if(pt>= 3.00 &&  pt < 4.00){
4916                         fEoverPCutMin=0.8175;
4917                         fEoverPCutMax=1.1101;
4918                 }
4919                 if(pt>= 4.00 &&  pt < 5.00){
4920                         fEoverPCutMin=0.8302;
4921                         fEoverPCutMax=1.1194;
4922                 }
4923                 if(pt>= 5.00 &&  pt < 6.00){
4924                         fEoverPCutMin=0.8517;
4925                         fEoverPCutMax=1.1177;
4926                 }
4927                 if(pt>= 6.00 &&  pt < 8.00){
4928                         fEoverPCutMin=0.8901;
4929                         fEoverPCutMax=1.1139;
4930                 }
4931                 if(pt>= 8.00 &&  pt < 10.00){
4932                         fEoverPCutMin=0.8703;
4933                         fEoverPCutMax=1.1377;
4934                 }
4935                 if(pt>= 10.00 &&  pt < 12.00){
4936                         fEoverPCutMin=0.9043;
4937                         fEoverPCutMax=1.1977;
4938                 }
4939         }
4940
4941         ////=========================================== MC 2 sigma ===========================================
4942
4943
4944         if(fIsMC){
4945
4946                 if(pt>= 2.00 &&  pt < 2.50){
4947                         fEoverPCutMin=0.7712;
4948                         fEoverPCutMax=1.0746;
4949                 }
4950                 if(pt>= 2.50 &&  pt < 3.00){
4951                         fEoverPCutMin=0.7946;
4952                         fEoverPCutMax=1.0708;
4953                 }
4954                 if(pt>= 3.00 &&  pt < 4.00){
4955                         fEoverPCutMin=0.8196;
4956                         fEoverPCutMax=1.0678;
4957                 }
4958                 if(pt>= 4.00 &&  pt < 5.00){
4959                         fEoverPCutMin=0.8396;
4960                         fEoverPCutMax=1.0646;
4961                 }
4962                 if(pt>= 5.00 &&  pt < 6.00){
4963                         fEoverPCutMin=0.8527;
4964                         fEoverPCutMax=1.0647;
4965                 }
4966                 if(pt>= 6.00 &&  pt < 8.00){
4967                         fEoverPCutMin=0.8652;
4968                         fEoverPCutMax=1.0670;
4969                 }
4970                 if(pt>= 8.00 &&  pt < 10.00){
4971                         fEoverPCutMin=0.8748;
4972                         fEoverPCutMax=1.0804;
4973                 }
4974                 if(pt>= 10.00 &&  pt < 12.00){
4975                         fEoverPCutMin=0.8814;
4976                         fEoverPCutMax=1.0842;
4977                 }
4978
4979         }
4980         */
4981                 ////======================================================================================
4982                 ////======================================================================================
4983                 ////=========================== data 1.5 sigma ===========================================
4984         /*
4985         if(!fIsMC){
4986                 //data 1.5 sigmas
4987                 if(pt>= 2.00 &&  pt < 2.50){
4988                         fEoverPCutMin=0.8132;
4989                         fEoverPCutMax=1.0610;
4990                 }
4991                 if(pt>= 2.50 &&  pt < 3.00){
4992                         fEoverPCutMin=0.8356;
4993                         fEoverPCutMax=1.0698;
4994                 }
4995                 if(pt>= 3.00 &&  pt < 4.00){
4996                         fEoverPCutMin=0.8541;
4997                         fEoverPCutMax=1.0735;
4998                 }
4999                 if(pt>= 4.00 &&  pt < 5.00){
5000                         fEoverPCutMin=0.8664;
5001                         fEoverPCutMax=1.0832;
5002                 }
5003                 if(pt>= 5.00 &&  pt < 6.00){
5004                         fEoverPCutMin=0.8849;
5005                         fEoverPCutMax=1.0845;
5006                 }
5007                 if(pt>= 6.00 &&  pt < 8.00){
5008                         fEoverPCutMin=0.9180;
5009                         fEoverPCutMax=1.0860;
5010                 }
5011                 if(pt>= 8.00 &&  pt < 10.00){
5012                         fEoverPCutMin=0.9037;
5013                         fEoverPCutMax=1.1043;
5014                 }
5015                 if(pt>= 10.00 &&  pt < 12.00){
5016                         fEoverPCutMin=0.9409;
5017                         fEoverPCutMax=1.1611;
5018                 }
5019         }
5020
5021         ////======================================================================================
5022         ////======================================================================================
5023         ////=========================== MC 1.5 sigma ===========================================
5024
5025         if(fIsMC){
5026                 //MC 1.5 sigmas
5027                 if(pt>= 2.00 &&  pt < 2.50){
5028                         fEoverPCutMin=0.8091;
5029                         fEoverPCutMax=1.0367;
5030                 }
5031                 if(pt>= 2.50 &&  pt < 3.00){
5032                 fEoverPCutMin=0.8292;
5033                 fEoverPCutMax=1.0362;
5034                 }
5035                 if(pt>= 3.00 &&  pt < 4.00){
5036                         fEoverPCutMin=0.8506;
5037                         fEoverPCutMax=1.0368;
5038                 }
5039                 if(pt>= 4.00 &&  pt < 5.00){
5040                         fEoverPCutMin=0.8677;
5041                         fEoverPCutMax=1.0365;
5042                 }
5043                 if(pt>= 5.00 &&  pt < 6.00){
5044                         fEoverPCutMin=0.8792;
5045                         fEoverPCutMax=1.0382;
5046                 }
5047                 if(pt>= 6.00 &&  pt < 8.00){
5048                         fEoverPCutMin=0.8904;
5049                         fEoverPCutMax=1.0418;
5050                 }
5051                 if(pt>= 8.00 &&  pt < 10.00){
5052                         fEoverPCutMin=0.9005;
5053                         fEoverPCutMax=1.0547;
5054                 }
5055                 if(pt>= 10.00 &&  pt < 12.00){
5056                         fEoverPCutMin=0.9067;
5057                         fEoverPCutMax=1.0589;
5058                 }
5059
5060         }
5061         */
5062                 ////======================================================================================
5063                 ////======================================================================================
5064                 ////=========================== data 2.5 sigma ===========================================
5065
5066         if(!fIsMC){
5067                         //data 2.5 sigmas
5068                 if(pt>= 2.00 &&  pt < 2.50){
5069                         fEoverPCutMin=0.7306;
5070                         fEoverPCutMax=1.1436;
5071                 }
5072                 if(pt>= 2.50 &&  pt < 3.00){
5073                         fEoverPCutMin=0.7575;
5074                         fEoverPCutMax=1.1479;
5075                 }
5076                 if(pt>= 3.00 &&  pt < 4.00){
5077                         fEoverPCutMin=0.7809;
5078                         fEoverPCutMax=1.1467;
5079                 }
5080                 if(pt>= 4.00 &&  pt < 5.00){
5081                         fEoverPCutMin=0.7941;
5082                         fEoverPCutMax=1.1555;
5083                 }
5084                 if(pt>= 5.00 &&  pt < 6.00){
5085                         fEoverPCutMin=0.8184;
5086                         fEoverPCutMax=1.1510;
5087                 }
5088                 if(pt>= 6.00 &&  pt < 8.00){
5089                         fEoverPCutMin=0.8621;
5090                         fEoverPCutMax=1.1419;
5091                 }
5092                 if(pt>= 8.00 &&  pt < 10.00){
5093                         fEoverPCutMin=0.8368;
5094                         fEoverPCutMax=1.1712;
5095                 }
5096                 if(pt>= 10.00 &&  pt < 12.00){
5097                         fEoverPCutMin=0.8676;
5098                         fEoverPCutMax=1.2344;
5099                 }
5100
5101 }
5102
5103                 ////======================================================================================
5104                 ////======================================================================================
5105                 ////=========================== MC 2.5 sigma ===========================================
5106
5107         if(fIsMC){
5108                         //MC 2.5 sigmas
5109                 if(pt>= 2.00 &&  pt < 2.50){
5110                         fEoverPCutMin=0.7333;
5111                         fEoverPCutMax=1.1125;
5112                 }
5113                 if(pt>= 2.50 &&  pt < 3.00){
5114                         fEoverPCutMin=0.7601;
5115                         fEoverPCutMax=1.1053;
5116                 }
5117                 if(pt>= 3.00 &&  pt < 4.00){
5118                         fEoverPCutMin=0.7885;
5119                         fEoverPCutMax=1.0989;
5120                 }
5121                 if(pt>= 4.00 &&  pt < 5.00){
5122                         fEoverPCutMin=0.8115;
5123                         fEoverPCutMax=1.0927;
5124                 }
5125                 if(pt>= 5.00 &&  pt < 6.00){
5126                         fEoverPCutMin=0.8262;
5127                         fEoverPCutMax=1.0912;
5128                 }
5129                 if(pt>= 6.00 &&  pt < 8.00){
5130                         fEoverPCutMin=0.8400;
5131                         fEoverPCutMax=1.0922;
5132                 }
5133                 if(pt>= 8.00 &&  pt < 10.00){
5134                         fEoverPCutMin=0.8492;
5135                         fEoverPCutMax=1.1060;
5136                 }
5137                 if(pt>= 10.00 &&  pt < 12.00){
5138                         fEoverPCutMin=0.8560;
5139                         fEoverPCutMax=1.1096;
5140                 }
5141         }
5142
5143
5144         return fEoverPCutMin;
5145         return fEoverPCutMax;
5146
5147 }
5148
5149
5150