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 = 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 = fAOD->GetHeader();
2163                         //Int_t bc= aodh->GetBunchCrossNumber();
2164
2165
2166                 Int_t ClsNo_aod = -999;
2167                 ClsNo_aod = fAOD->GetNumberOfCaloClusters();
2168                 fNCluster_pure_aod->Fill(ClsNo_aod);
2169                         //Bool_t exotic=kTRUE;
2170
2171                 for (Int_t i=0; i< ClsNo_aod; i++ ){
2172
2173                         //fClus = fVevent->GetCaloCluster(i);
2174                         //to be compatible with Shingo
2175                     AliVCluster *clust = 0x0;
2176                         clust = (AliVCluster*) fAOD->GetCaloCluster(i);
2177
2178                         if(clust && clust->IsEMCAL())
2179                         {
2180                                 //pure cluster information
2181                                 fECluster_pure->Fill(clust->E());
2182
2183                                 fNcells_energy->Fill(clust->GetNCells(),clust->E());
2184                                 fNCluster_ECluster->Fill(ClsNo_aod,clust->E());
2185
2186                                 if(clust->E()>1000) fNevent->Fill(23);
2187
2188                                 //exotic
2189                                 /*
2190                                 exotic   = fEMCALRecoUtils->IsExoticCluster(clust, (AliVCaloCells*)fAOD->GetEMCALCells(), bc);
2191                                 if(exotic == kFALSE){
2192                                         fECluster_not_exotic->Fill(clust->E());
2193                                         fNcells_energy_not_exotic->Fill(clust->GetNCells(),clust->E());
2194                                 }
2195                                 */
2196
2197                                 //approximation to remove exotics
2198                                 if(clust->GetNCells()<5 && clust->E()>15.0){
2199                                         fECluster_exotic->Fill(clust->E());
2200                                 }
2201                                 //Marcel cut (another approximation to remove exotics)
2202                                 else if((clust->GetNCells())> ((clust->E())/3+1)){
2203                                         fECluster_not_exotic1->Fill(clust->E());
2204                                 }
2205                                 else{
2206                                         fECluster_not_exotic2->Fill(clust->E());
2207                                 }
2208
2209
2210                         }
2211                         /*
2212                         //______________________________________________________________________
2213                         //Trying to remove events with bad cells and find patches
2214                         //First, I will try to count them
2215                         //______________________________________________________________________
2216
2217                         if(clust && clust->IsEMCAL())
2218                         {
2219                                 Bool_t badchannel = ContainsBadChannel("EMCAL", clust->GetCellsAbsId(),clust->GetNCells() );
2220                                 printf("Contém bad channel? %d ", badchannel);
2221                                 if(badchannel)fNevent2->Fill(0);
2222
2223                                 //trying to find patches
2224                                 TArrayI patches_found=GetTriggerPatches(IsEventEMCALL0, IsEventEMCALL1);
2225                                 printf("N patches %d, first %d, last %d\n",patches_found.GetSize(),  patches_found.At(0), patches_found.At(patches_found.GetSize()-1));
2226
2227                         }
2228
2229                         //end of bad cells
2230                         //______________________________________________________________________
2231                         */
2232
2233                 }
2234         }
2235
2236
2237         fNevent->Fill(24); //events with cluster
2238
2239
2240         fVtxZ_new4->Fill(fZvtx);
2241
2242
2243         Int_t ClsNo = -999;
2244         if(!fIsAOD) ClsNo = fESD->GetNumberOfCaloClusters();
2245         else ClsNo = fAOD->GetNumberOfCaloClusters();
2246
2247         //______________________________________________________________________
2248
2249         ///_____________________________________________________________________
2250         ///Track loop
2251         Int_t NTracks=fVevent->GetNumberOfTracks();
2252
2253         ///////////////////////////////////////////////////////////////////////////////////////////////////////////
2254         //To use tender
2255         if(fUseTender){
2256                 TClonesArray  *fTracks_tender = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("AODFilterTracks"));
2257                 NTracks = fTracks_tender->GetEntries();
2258                 TClonesArray  *fCaloClusters_tender = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("EmcCaloClusters"));
2259                 ClsNo = fCaloClusters_tender->GetEntries();
2260
2261         }
2262         ///////////////////////////////////////////////////////////////////////////////////////////////////////////
2263
2264         for(Int_t iTracks = 0; iTracks < NTracks; iTracks++)
2265         {
2266                 AliVParticle* Vtrack = fVevent->GetTrack(iTracks);
2267                 if (!Vtrack)
2268                 {
2269                         printf("ERROR: Could not receive track %d\n", iTracks);
2270                         continue;
2271                 }
2272
2273                         //main place where track is defined
2274
2275                 AliVTrack *track = dynamic_cast<AliVTrack*>(Vtrack);
2276                 AliESDtrack *etrack = dynamic_cast<AliESDtrack*>(Vtrack);
2277                 AliAODTrack *atrack = dynamic_cast<AliAODTrack*>(Vtrack);
2278
2279                 ///////////////////////////////////////////////////////////////////////////////////////////////////////////
2280                 //To use tender
2281                 if(fUseTender){
2282                         TClonesArray  *fTracks_tender = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("AODFilterTracks"));
2283                         track = static_cast<AliAODTrack*>(fTracks_tender->At(iTracks));
2284                         if (!track) {
2285                                 printf("ERROR: Could not receive track calibrated %d\n", iTracks);
2286                                 continue;
2287                         }
2288                 }
2289                 ///////////////////////////////////////////////////////////////////////////////////////////////////////////
2290
2291                 Double_t fTPCnSigma = -999;
2292                 Double_t fTOFnSigma = -999;
2293                 Double_t fTPCnSigma_pion = -999;
2294                 Double_t fTPCnSigma_proton = -999;
2295                 Double_t fTPCnSigma_kaon = -999;
2296                 Double_t fTPCsignal = -999;
2297                 Double_t fPt = -999;
2298                 Double_t fP = -999;
2299
2300                 //December 9th 2013
2301                 //Etacut test on the begging
2302                 fEtad[0]->Fill(track->Eta());
2303                         //if(track->Eta()<fEtaCutMin || track->Eta()>fEtaCutMax) continue;
2304                 fEtad[1]->Fill(track->Eta());
2305
2306
2307
2308
2309                 ///_____________________________________________________________________________
2310                 ///Fill QA plots without track selection
2311                 fPt = track->Pt();
2312                 fP = TMath::Sqrt((track->Pt())*(track->Pt()) + (track->Pz())*(track->Pz()));
2313
2314                 fTPCsignal = track->GetTPCsignal();
2315                 fTPCnSigma = fPidResponse->NumberOfSigmasTPC(track, AliPID::kElectron);
2316                 fTOFnSigma = fPidResponse->NumberOfSigmasTOF(track, AliPID::kElectron);
2317                 fTPCnSigma_pion = fPidResponse->NumberOfSigmasTPC(track, AliPID::kPion);
2318                 fTPCnSigma_proton = fPidResponse->NumberOfSigmasTPC(track, AliPID::kProton);
2319                 fTPCnSigma_kaon = fPidResponse->NumberOfSigmasTPC(track, AliPID::kKaon);
2320
2321                 fTPC_p[0]->Fill(fPt,fTPCsignal);
2322                 fTPCnsigma_p[0]->Fill(fP,fTPCnSigma);
2323
2324
2325                 Float_t TPCNcls = track->GetTPCNcls();
2326                 //TPC Ncls for pid
2327                 Float_t TPCNcls_pid = track->GetTPCsignalN();
2328
2329
2330
2331                 Float_t pos[3]={0,0,0};
2332
2333                 Double_t fEMCflag = kFALSE;
2334                 if(track->GetEMCALcluster()>0)
2335                 {
2336                         fClus = fVevent->GetCaloCluster(track->GetEMCALcluster());
2337                         ///////////////////////////////////////////////////////////////////////////////////////////////////////////
2338                         //to use tender
2339                         if(fUseTender){
2340                                 int EMCalIndex = -1;
2341                                 EMCalIndex = track->GetEMCALcluster();
2342                                 TClonesArray  *fCaloClusters_tender = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("EmcCaloClusters"));
2343                                 fClus = static_cast<AliVCluster*>(fCaloClusters_tender->At(EMCalIndex));
2344                         }
2345             ///////////////////////////////////////////////////////////////////////////////////////////////////////////
2346
2347                         if(fClus->IsEMCAL())
2348                         {
2349
2350                                 //only for charged tracks
2351                                 fECluster[0]->Fill(fClus->E());
2352
2353
2354                                 if(TMath::Abs(fClus->GetTrackDx())<=fdPhiCut && TMath::Abs(fClus->GetTrackDz())<=fdEtaCut)
2355                                 {
2356                                         fEMCflag = kTRUE;
2357                                         fEoverP_pt[0]->Fill(fPt,(fClus->E() / fP));
2358
2359
2360                                         Float_t Energy  = fClus->E();
2361                                         Float_t EoverP  = Energy/track->P();
2362                                                 //Float_t M02   = fClus->GetM02();
2363                                                 //Float_t M20   = fClus->GetM20();
2364
2365                                                 /////////////// for Eta Phi distribution
2366                                         fClus->GetPosition(pos);
2367                                         TVector3 vpos(pos[0],pos[1],pos[2]);
2368                                         Double_t cphi = vpos.Phi();
2369                                         Double_t ceta = vpos.Eta();
2370                                         fEtaPhi[0]->Fill(cphi,ceta);
2371
2372
2373                                         fTPCNcls_EoverP[0]->Fill(TPCNcls, EoverP);
2374                                 }
2375                         }
2376                 }
2377
2378                         //______________________________________________________________
2379                         // Vertex
2380
2381                 fVtxZ[0]->Fill(fZvtx);
2382                 if(iTracks == 0)fNTracks[0]->Fill(fNOtrks);
2383                 fNTracks_pt[0]->Fill(fNOtrks, fPt);
2384                 fNTracks_eta[0]->Fill(fNOtrks, track->Eta());
2385                 fNTracks_phi[0]->Fill(fNOtrks, track->Phi());
2386
2387
2388                 fNClusters[0]->Fill(ClsNo);
2389                 fTPCNcls_pid[0]->Fill(TPCNcls, TPCNcls_pid);
2390                         //______________________________________________________________
2391
2392 ///Fill QA plots without track selection
2393 ///_____________________________________________________________________________
2394 //______________________________________________________________________________________
2395 //Track Selection Cuts
2396
2397 //AOD (Test Filter Bit)
2398                 if(fIsAOD)
2399                 {
2400                                 // standard cuts with very loose DCA - BIT(4)
2401                                 // Description:
2402                         /*
2403                          GetStandardITSTPCTrackCuts2011(kFALSE)
2404                          SetMaxChi2PerClusterTPC(4);
2405                          SetAcceptKinkDaughters(kFALSE);
2406                          SetRequireTPCRefit(kTRUE);
2407                          SetRequireITSRefit(kTRUE);
2408                          SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kAny);
2409                          SetMaxDCAToVertexZ(2);
2410                          SetMaxDCAToVertex2D(kFALSE);
2411                          SetRequireSigmaToVertex(kFALSE);
2412                          SetMaxChi2PerClusterITS(36);
2413                          SetMaxDCAToVertexXY(2.4)
2414                          SetMaxDCAToVertexZ(3.2)
2415                          SetDCaToVertex2D(kTRUE)
2416                          */
2417
2418                         if(!atrack->TestFilterMask(AliAODTrack::kTrkGlobalNoDCA)) continue;
2419                 }
2420
2421 //RecKine: ITSTPC cuts
2422                 if(!ProcessCutStep(AliHFEcuts::kStepRecKineITSTPC, track)) continue;
2423 //RecKink
2424                 if(fRejectKinkMother)
2425                 {
2426                         if(fIsAOD)
2427                         {
2428                                 Bool_t kinkmotherpass = kTRUE;
2429                                 for(Int_t kinkmother = 0; kinkmother < fNumberOfMotherkink; kinkmother++)
2430                                 {
2431                                         if(track->GetID() == fListOfmotherkink[kinkmother])
2432                                         {
2433                                                 kinkmotherpass = kFALSE;
2434                                                 continue;
2435                                         }
2436                                 }
2437                                 if(!kinkmotherpass) continue;
2438                         }
2439                         else
2440                         {
2441                                 if(etrack->GetKinkIndex(0) != 0) continue;
2442                         }
2443                 }
2444
2445 //RecPrim
2446
2447                         //it was not working on aod... testing again
2448                         //July 29th, 2014: aparently working again
2449
2450 ///if(!fIsAOD)
2451 //{
2452
2453                 if(!ProcessCutStep(AliHFEcuts::kStepRecPrim, track)) continue;
2454 //}
2455
2456 //HFEcuts: ITS layers cuts
2457                 if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsITS, track)) continue;
2458
2459 //HFE cuts: TPC PID cleanup
2460                 if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTPC, track)) continue;
2461
2462 //DCA cut done by hand -- July 29th, 2014  -> fix it... has to be in absolute values!!!!
2463                 if(fIsAOD){
2464                         Double_t d0z0[2], cov[3];
2465                         AliAODVertex *prim_vtx = fAOD->GetPrimaryVertex();
2466                         track->PropagateToDCA(prim_vtx, fAOD->GetMagneticField(), 20., d0z0, cov);
2467                         Double_t DCAxy = d0z0[0];
2468                         Double_t DCAz = d0z0[1];
2469
2470                         if(DCAxy >= fDCAcutr || DCAz>=fDCAcutz) continue;
2471                 }
2472
2473
2474
2475
2476
2477 //______________________________________________________________________________________
2478
2479                 if(fIsAOD){
2480                                 //AOD test -- Francesco suggestion
2481                                 //aod test -- Francesco suggestion
2482                         AliAODTrack *aod_track=fAOD->GetTrack(iTracks);
2483
2484                         Int_t type=aod_track->GetType();
2485                         if(type==AliAODTrack::kPrimary) fPtPrim->Fill(aod_track->Pt());
2486                         if(type==AliAODTrack::kSecondary) fPtSec->Fill(aod_track->Pt());
2487
2488                                 //Int_t type2=track->GetType();
2489                         if(type==AliAODTrack::kPrimary) fPtPrim2->Fill(track->Pt());
2490                         if(type==AliAODTrack::kSecondary) fPtSec2->Fill(track->Pt());
2491                 }
2492
2493
2494 ///_____________________________________________________________
2495 ///QA plots after track selection
2496                 if(fIsMC)
2497                 {
2498                         if(track->GetLabel()>=0) fPtMCWithLabel->Fill(fPt);
2499                         else fPtMCWithoutLabel->Fill(fPt);
2500                 }
2501
2502                 if(fIsMC  && track->GetLabel()>=0)
2503                 {
2504                         if(fIsAOD)
2505                         {
2506                                 fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
2507
2508
2509                                 if(fMCparticle->IsPhysicalPrimary())
2510                                 {
2511                                         fPtIsPhysicaPrimary->Fill(fPt);
2512                                 }
2513
2514                                 Int_t pdg = fMCparticle->GetPdgCode();
2515                                 if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax && fMCparticle->Charge()!=0)
2516                                 {
2517
2518
2519                                         if( TMath::Abs(pdg) == 211 || TMath::Abs(pdg) == 2212 || TMath::Abs(pdg) == 321 || TMath::Abs(pdg) == 11 || TMath::Abs(pdg) == 13 )
2520                                         {
2521                                                 fPtMCparticleReco_nonPrimary->Fill(fMCparticle->Pt()); //not Primary track
2522
2523                                                 if(fMCparticle->IsPhysicalPrimary())
2524                                                 {
2525                                                         fPtMCparticleReco->Fill(fMCparticle->Pt());
2526
2527                                                         Bool_t MotherFound = FindMother(track->GetLabel());
2528                                                         if(MotherFound)
2529                                                         {
2530                                                                 if(fIsHFE1)
2531                                                                 {
2532                                                                         fPtMCparticleRecoHfe1->Fill(fMCparticle->Pt());//numerator tracking
2533                                                                         //unfolding
2534                                                                         fpt_reco_pt_MC_den->Fill(track->Pt(),fMCparticle->Pt());
2535
2536                                                                 }
2537                                                                 if(fIsHFE2){
2538                                                                         fPtMCparticleRecoHfe2->Fill(fMCparticle->Pt());
2539                                                                 }
2540                                                         }
2541                                                 }
2542                                         }
2543                                 }
2544
2545
2546
2547                         }//close AOD
2548                          //ESD
2549                         else
2550                         {
2551
2552
2553                                 if(fMCtrack->Eta()>=fEtaCutMin && fMCtrack->Eta()<=fEtaCutMax)
2554                                 {
2555
2556                                         fMCtrack = fMCstack->Particle(track->GetLabel());
2557                                         Int_t pdg = fMCtrack->GetPdgCode();
2558
2559                                         if( TMath::Abs(pdg) == 211 || TMath::Abs(pdg) == 2212 || TMath::Abs(pdg) == 321 || TMath::Abs(pdg) == 11 || TMath::Abs(pdg) == 13 )
2560                                         {
2561                                                 fPtMCparticleReco_nonPrimary->Fill(fMCtrack->Pt());//not Primary track
2562                                         }
2563
2564
2565                                         if(fMCstack->IsPhysicalPrimary(track->GetLabel()))
2566                                         {
2567                                                 fPtIsPhysicaPrimary->Fill(fPt);
2568
2569
2570                                                 if( TMath::Abs(pdg) == 211 || TMath::Abs(pdg) == 2212 || TMath::Abs(pdg) == 321 || TMath::Abs(pdg) == 11 || TMath::Abs(pdg) == 13 )
2571                                                 {
2572                                                         fPtMCparticleReco->Fill(fMCtrack->Pt());
2573
2574                                                         Bool_t MotherFound = FindMother(track->GetLabel());
2575                                                         if(MotherFound)
2576                                                         {
2577                                                                 if(fIsHFE1) fPtMCparticleRecoHfe1->Fill(fMCtrack->Pt());//numerator tracking
2578                                                                 if(fIsHFE2) fPtMCparticleRecoHfe2->Fill(fMCtrack->Pt());
2579                                                         }
2580                                                 }
2581                                         }
2582                                 }
2583                         }//close ESD
2584                 }//close IsMC
2585
2586                 fTPC_p[1]->Fill(fPt,fTPCsignal);
2587                 fTPCnsigma_p[1]->Fill(fP,fTPCnSigma);
2588                 Double_t fPtBin[7] = {1,2,4,6,8,10,15};
2589                 Double_t fPtBin_trigger[11] = {1,2,4,6,8,10,12,14,16,18,20};
2590
2591                 TPCNcls = track->GetTPCNcls();
2592                 Float_t pos2[3]={0,0,0};
2593
2594                 if(track->GetEMCALcluster()>0)
2595                 {
2596                         fClus = fVevent->GetCaloCluster(track->GetEMCALcluster());
2597                         ///////////////////////////////////////////////////////////////////////////////////////////////////////////
2598                         //to use tender
2599                         if(fUseTender){
2600                                 int EMCalIndex = -1;
2601                                 EMCalIndex = track->GetEMCALcluster();
2602                                 TClonesArray  *fCaloClusters_tender = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("EmcCaloClusters"));
2603                                 fClus = static_cast<AliVCluster*>(fCaloClusters_tender->At(EMCalIndex));
2604                         }
2605                         ///////////////////////////////////////////////////////////////////////////////////////////////////////////
2606
2607                         if(fClus->IsEMCAL())
2608                         {
2609                                 if(TMath::Abs(fClus->GetTrackDx())<=fdPhiCut && TMath::Abs(fClus->GetTrackDz())<=fdEtaCut)
2610                                 {
2611                                         fEoverP_pt[1]->Fill(fPt,(fClus->E() / fP));
2612
2613                                         Float_t Energy  = fClus->E();
2614                                         Float_t EoverP  = Energy/track->P();
2615                                         Float_t M02     = fClus->GetM02();
2616                                         Float_t M20     = fClus->GetM20();
2617                                         Float_t ncells  = fClus->GetNCells();
2618
2619                                                 /////////////// for Eta Phi distribution
2620                                         fClus->GetPosition(pos2);
2621                                         TVector3 vpos(pos2[0],pos2[1],pos2[2]);
2622                                         Double_t cphi = vpos.Phi();
2623                                         Double_t ceta = vpos.Eta();
2624                                         fEtaPhi[1]->Fill(cphi,ceta);
2625
2626                                         fECluster[1]->Fill(Energy);
2627                                         fTPCNcls_EoverP[1]->Fill(TPCNcls, EoverP);
2628
2629
2630                                         //for EMCal trigger performance
2631                                         if(EoverP > 0.9){
2632                                                 ftpc_p_EoverPcut->Fill(track->P(), fTPCsignal);
2633                                                 fnsigma_p_EoverPcut->Fill(track->P(), fTPCnSigma);
2634
2635                                         }
2636
2637
2638                                         //for hadron contamination calculations
2639
2640
2641                                         // EtaCut -> dados
2642                                         if(track->Eta()>=fEtaCutMin && track->Eta()<=fEtaCutMax ){
2643
2644                                                 //-------------------------------------------------------------------
2645                                                 //true hadrons E/p shape before cuts
2646                                                 if(fIsMC && fIsAOD && track->GetLabel()>=0){
2647                                                         fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
2648                                                         Int_t pdg = fMCparticle->GetPdgCode();
2649
2650                                                         if( TMath::Abs(pdg) != 11){
2651                                                                 fEoverP_pt_true_hadrons0->Fill(fPt,(fClus->E() / fP));
2652                                                         }
2653                                                 }
2654
2655
2656                                                 //true electrons E/p shape before cuts
2657                                                 if(fIsMC && fIsAOD && track->GetLabel()>=0){
2658                                                         fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
2659                                                         Int_t pdg = fMCparticle->GetPdgCode();
2660
2661                                                         if( TMath::Abs(pdg) == 11){
2662                                                                 fEoverP_pt_true_electrons0->Fill(fPt,(fClus->E() / fP));
2663                                                         }
2664                                                 }
2665                                                 //-------------------------------------------------------------------
2666
2667
2668                                                 //main
2669                                                 if(TMath::Abs(fTPCnSigma_pion)<3 || TMath::Abs(fTPCnSigma_proton)<3 || TMath::Abs(fTPCnSigma_kaon)<3 ){
2670
2671                                                         if(fTPCnSigma<-3.5){
2672                                                                 fEoverP_pt_hadrons->Fill(fPt,EoverP);
2673                                                                 if(fUseEMCal) fShowerShape_ha->Fill(M02,M20);
2674
2675
2676                                                                 //after cut -> Are they really hadrons? E/p shape
2677                                                                 if(fIsMC && fIsAOD && track->GetLabel()>=0){
2678                                                                         fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
2679                                                                         Int_t pdg = fMCparticle->GetPdgCode();
2680
2681                                                                         if( TMath::Abs(pdg) != 11){
2682                                                                                 fEoverP_pt_true_hadrons->Fill(fPt,(fClus->E() / fP));
2683                                                                         }
2684                                                                 }
2685
2686
2687                                                         }
2688                                                 }
2689
2690                                                 //for systematic studies of hadron contamination
2691                                                 if(fTPCnSigma < -3){
2692                                                         fEoverP_pt_pions->Fill(fPt, EoverP);
2693
2694                                                 }
2695
2696                                                 if(fTPCnSigma < -3.5){
2697                                                         fEoverP_pt_pions2->Fill(fPt, EoverP);
2698
2699                                                 }
2700
2701
2702                                         }
2703
2704
2705
2706
2707                                         for(Int_t i = 0; i < 6; i++)
2708                                         {
2709                                                 if(fPt>=fPtBin[i] && fPt<fPtBin[i+1])
2710                                                 {
2711
2712                                                         if(fTPCnSigma < -5 && fTPCnSigma > -10){
2713                                                                 fNcells_hadrons[i]->Fill(ncells);
2714                                                         }
2715                                                                 //hadrons selection using E/p
2716                                                         if((fClus->E() / fP) >= 0.2 && (fClus->E() / fP) <= 0.7){
2717                                                                 fTPCnsigma_eta_hadrons[i]->Fill(fTPCnSigma, ceta);
2718                                                         }
2719                                                                 //electrons selection using E/p
2720                                                         if((fClus->E() / fP) >= fEoverPCutMin && (fClus->E() / fP) <= fEoverPCutMax) {
2721                                                                 fTPCnsigma_eta_electrons[i]->Fill(fTPCnSigma, ceta);
2722                                                         }
2723                                                 }
2724                                         }
2725
2726                                         if((fClus->E() / fP) >= fEoverPCutMin && (fClus->E() / fP) <= fEoverPCutMax)
2727                                         {
2728                                                 fTPCnsigma_eta->Fill(track->Eta(),fTPCnSigma);
2729                                                 fTPCnsigma_phi->Fill(track->Phi(),fTPCnSigma);
2730
2731                                                 if(fUseEMCal)
2732                                                 {
2733                                                         if(fTPCnSigma < 3.5 && fCorrelationFlag)
2734                                                         {
2735                                                                 fPtTrigger_Inc->Fill(fPt);
2736                                                                 DiHadronCorrelation(track, iTracks);
2737                                                         }
2738                                                 }
2739                                         }
2740
2741                                 }
2742                         }
2743                 }
2744
2745                 //______________________________________________________________
2746                 // Vertex
2747
2748                 fVtxZ[1]->Fill(fZvtx);
2749                 if(iTracks == 0)fNTracks[1]->Fill(fNOtrks);
2750                 fNTracks_pt[1]->Fill(fNOtrks, fPt);
2751                 fNTracks_eta[1]->Fill(fNOtrks, track->Eta());
2752                 fNTracks_phi[1]->Fill(fNOtrks, track->Phi());
2753                 fNClusters[1]->Fill(ClsNo);
2754                 fTPCNcls_pid[1]->Fill(TPCNcls, TPCNcls_pid);
2755
2756                 //______________________________________________________________
2757
2758                 ///______________________________________________________________________
2759                 ///Histograms for PID Studies
2760                 //Double_t fPtBin[6] = {2,4,6,8,10,15};
2761
2762                 for(Int_t i = 0; i < 6; i++)
2763                 {
2764                         if(fPt>=fPtBin[i] && fPt<fPtBin[i+1])
2765                         {
2766                                 if(fEMCflag) fEoverP_tpc[i]->Fill(fTPCnSigma,(fClus->E() / fP));
2767                                 fTPC_pt[i]->Fill(fTPCsignal);
2768                                 fTPCnsigma_pt[i]->Fill(fTPCnSigma);
2769
2770                         }
2771                 }
2772
2773                 if(track->Eta()>=fEtaCutMin && track->Eta()<=fEtaCutMax ){
2774
2775                 //new pt bins for trigger data
2776
2777                         for(Int_t i = 0; i < 10; i++)
2778                         {
2779                                 if(fP>=fPtBin_trigger[i] && fP<fPtBin_trigger[i+1])
2780                                 {
2781                                         if(fEMCflag)fEoverP_tpc_p_trigger[i]->Fill(fTPCnSigma,(fClus->E() / fP));
2782
2783                                 }
2784
2785                                 if(fPt>=fPtBin_trigger[i] && fPt<fPtBin_trigger[i+1])
2786                                 {
2787                                         if(fEMCflag)fEoverP_tpc_pt_trigger[i]->Fill(fTPCnSigma,(fClus->E() / fP));
2788
2789                                 }
2790                         }
2791
2792
2793                         //new way to calculate TPCnsigma distribution: TPCnsigma in function of p, with/without E/p cut
2794                         fTPCnsigma_p_TPC->Fill(fP, fTPCnSigma);
2795                         if(fEMCflag){
2796
2797                                 fTPCnsigma_p_TPC_on_EMCal_acc->Fill(fP, fTPCnSigma);
2798
2799                                 if((fClus->E() / fP) >= fEoverPCutMin && (fClus->E() / fP) <= fEoverPCutMax){
2800                                         fTPCnsigma_p_TPC_EoverP_cut->Fill(fP, fTPCnSigma);
2801                                 }
2802
2803                         }//close EMCflag
2804
2805                 }//close eta cut
2806
2807                 ///QA plots after track selection
2808                 ///_____________________________________________________________
2809
2810                 //_______________________________________________________
2811                 //Correlation Analysis - DiHadron
2812                 if(!fUseEMCal)
2813                 {
2814                         if(fTPCnSigma < 3.5 && fCorrelationFlag)
2815                         {
2816                                 fPtTrigger_Inc->Fill(fPt);
2817                                 DiHadronCorrelation(track, iTracks);
2818                         }
2819                 }
2820                         //_______________________________________________________
2821
2822
2823                         ///////////////////////////////////////////////////////////////////
2824                         ///TPC - efficiency calculation //
2825
2826                         /// changing start here
2827                 if(fIsMC && fIsAOD && track->GetLabel()>=0)
2828                 {
2829                         fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
2830                         Int_t pdg = fMCparticle->GetPdgCode();
2831
2832                                 //
2833                         if(fMCparticle->IsPhysicalPrimary()){
2834
2835
2836                                 if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax ){
2837
2838                                         Bool_t MotherFound = FindMother(track->GetLabel());
2839                                         if(MotherFound){
2840                                                 fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
2841                                         if( TMath::Abs(pdg) == 11 && fMCparticleMother->GetPdgCode()!=22 ){
2842                                                         if(fIsHFE1) fPtMC_TPC_All->Fill(fMCparticle->Pt());
2843                                         }
2844                                         }
2845                                 }
2846                         }
2847                 }///until here
2848
2849                 else if(fIsMC && track->GetLabel()>=0)//ESD
2850                 {
2851
2852                         if(fMCstack->IsPhysicalPrimary(track->GetLabel())){
2853                                 fMCtrack = fMCstack->Particle(track->GetLabel());
2854
2855                                 Int_t pdg = fMCtrack->GetPdgCode();
2856                                 if(fMCtrack->Eta()>=fEtaCutMin && fMCtrack->Eta()<=fEtaCutMax ){
2857
2858                                         if(fMCtrack->GetFirstMother()>0){
2859                                             fMCtrackMother = fMCstack->Particle(fMCtrack->GetFirstMother());
2860                                                 Bool_t MotherFound = FindMother(track->GetLabel());
2861                                                 if(MotherFound){
2862                                                         if( TMath::Abs(pdg) == 11 && fMCtrackMother->GetPdgCode()!=22 ){
2863                                                                 if(fIsHFE1) fPtMC_TPC_All->Fill(fMCtrack->Pt());
2864                                                         }
2865                                                 }
2866                                         }
2867                                 }
2868                         }
2869                 }
2870
2871
2872                         if(fPt>1 && fPt<2) fTOF01->Fill(fTOFnSigma,fTPCnSigma);
2873                         if(fPt>2 && fPt<4) fTOF02->Fill(fTOFnSigma,fTPCnSigma);
2874                         if(fPt>4 && fPt<6) fTOF03->Fill(fTOFnSigma,fTPCnSigma);
2875
2876 ///________________________________________________________________________
2877 ///PID
2878 ///Here the PID cuts defined in the file "ConfigEMCalHFEpA.C" is applied
2879             Int_t pidpassed = 1;
2880             AliHFEpidObject hfetrack;
2881             hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);
2882             hfetrack.SetRecTrack(track);
2883             hfetrack.SetPP();   //proton-proton analysis
2884             if(!fPID->IsSelected(&hfetrack, NULL, "", fPIDqa)) pidpassed = 0;
2885             fpid->Fill(pidpassed);
2886
2887             if(pidpassed==0) continue;
2888 ///________________________________________________________________________
2889
2890
2891 ////////////////////////////////////////////////////////////////////
2892 ///TPC efficiency calculations
2893
2894
2895                 if(fIsMC && fIsAOD && track->GetLabel()>=0)
2896                 {
2897                         fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
2898                         Int_t pdg = fMCparticle->GetPdgCode();
2899
2900                                 //
2901                         if(fMCparticle->IsPhysicalPrimary()){
2902
2903
2904                                 if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax ){
2905
2906                                         Bool_t MotherFound = FindMother(track->GetLabel());
2907                                         if(MotherFound){
2908                                                 fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
2909                                         if( TMath::Abs(pdg) == 11 && fMCparticleMother->GetPdgCode()!=22 ){
2910                                                         if(fIsHFE1) fPtMC_TPC_Selected->Fill(fMCparticle->Pt());
2911                                         }
2912                                         }
2913                                 }
2914                         }
2915                 }///until here
2916
2917                 else if(fIsMC && track->GetLabel()>=0)//ESD
2918                 {
2919
2920                         if(fMCstack->IsPhysicalPrimary(track->GetLabel())){
2921                                 fMCtrack = fMCstack->Particle(track->GetLabel());
2922
2923                                 Int_t pdg = fMCtrack->GetPdgCode();
2924                                 if(fMCtrack->Eta()>=fEtaCutMin && fMCtrack->Eta()<=fEtaCutMax ){
2925
2926                                         if(fMCtrack->GetFirstMother()>0){
2927                                             fMCtrackMother = fMCstack->Particle(fMCtrack->GetFirstMother());
2928                                                 Bool_t MotherFound = FindMother(track->GetLabel());
2929                                                 if(MotherFound){
2930                                                         if( TMath::Abs(pdg) == 11 && fMCtrackMother->GetPdgCode()!=22 ){
2931                                                                 if(fIsHFE1) fPtMC_TPC_Selected->Fill(fMCtrack->Pt());
2932                                                         }
2933                                                 }
2934                                         }
2935                                 }
2936                         }
2937                 }
2938
2939                 //Eta Cut for TPC only
2940                 if(track->Eta()>=fEtaCutMin && track->Eta()<=fEtaCutMax ){
2941                         fTPCnsigma_pt_2D->Fill(fPt,fTPCnSigma);
2942                 }
2943
2944                 //Background for TPC only
2945                 if(fFillBackground){
2946
2947                         //efficiency without SS cut for TPC only
2948                         if(track->Eta()>=fEtaCutMin && track->Eta()<=fEtaCutMax){
2949                                 Background(track, iTracks, Vtrack, kTRUE); //IsTPConly=kTRUE
2950                         } //Eta cut to be consistent with other efficiency
2951                 }
2952
2953
2954                 fTPCnsigma_p[2]->Fill(fP,fTPCnSigma);
2955                 fTPC_p[2]->Fill(fP,fTPCsignal);
2956                 TPCNcls = track->GetTPCNcls();
2957                 Float_t pos3[3]={0,0,0};
2958
2959
2960                 //here denominator for track-matching efficiency
2961                 if(fIsMC && fIsAOD && track->GetLabel()>=0)
2962                 {
2963                         fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
2964                         Int_t pdg = fMCparticle->GetPdgCode();
2965
2966                                 //
2967                         if(fMCparticle->IsPhysicalPrimary()){
2968
2969
2970                                 if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax ){
2971
2972                                         Bool_t MotherFound = FindMother(track->GetLabel());
2973                                         if(MotherFound){
2974                                                 fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
2975                                         if( TMath::Abs(pdg) == 11 && fMCparticleMother->GetPdgCode()!=22 ){
2976                                                         if(fIsHFE1) fPt_track_match_den->Fill(fMCparticle->Pt());
2977                                         }
2978                                         }
2979                                 }
2980                         }
2981                 }///until here
2982
2983
2984                 if(track->GetEMCALcluster()>0)
2985                 {
2986                         fClus = fVevent->GetCaloCluster(track->GetEMCALcluster());
2987                         ///////////////////////////////////////////////////////////////////////////////////////////////////////////
2988                         //to use tender
2989                         if(fUseTender){
2990                                 int EMCalIndex = -1;
2991                                 EMCalIndex = track->GetEMCALcluster();
2992                                 TClonesArray  *fCaloClusters_tender = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("EmcCaloClusters"));
2993                                 fClus = static_cast<AliVCluster*>(fCaloClusters_tender->At(EMCalIndex));
2994                         }
2995                         ///////////////////////////////////////////////////////////////////////////////////////////////////////////
2996
2997                         if(fClus->IsEMCAL())
2998                         {
2999
3000                 //________________________________________________________________________
3001
3002
3003                                 //Cluster timing distribution -- for ESD
3004                                 if(fUseEMCal && !fIsAOD){
3005
3006                                         AliESDCaloCells &cells_esd=*(fESD->GetEMCALCells());
3007                                         TRefArray* caloClusters_esd = new TRefArray();
3008                                     fESD->GetEMCALClusters(caloClusters_esd);
3009                                         Int_t nclus_esd = caloClusters_esd->GetEntries();
3010
3011
3012                                         for (Int_t icl = 0; icl < nclus_esd; icl++) {
3013
3014                                                 AliESDCaloCluster* clus_esd = (AliESDCaloCluster*)caloClusters_esd->At(icl);
3015
3016                                                 if(clus_esd->IsEMCAL()){
3017                                                         Float_t ncells_esd      = fClus->GetNCells();
3018                                                         UShort_t * index_esd = clus_esd->GetCellsAbsId() ;
3019                                                         UShort_t * index2_esd = fClus->GetCellsAbsId() ;
3020
3021
3022                                                         for(Int_t i = 0; i < ncells_esd ; i++){
3023
3024                                                                 Int_t absId_esd =   index_esd[i];
3025                                                                 fTime->Fill(fPt,cells_esd.GetCellTime(absId_esd));
3026
3027                                                                 Int_t absId2_esd =   index2_esd[i];
3028                                                                 fTime2->Fill(fPt,cells_esd.GetCellTime(absId2_esd));
3029                                                         }
3030
3031                                                 }
3032                                         }
3033
3034                                 }
3035                                 /* not working!
3036                                 //Cluster timing distribution -- for AOD
3037                                 if(fUseEMCal && fIsAOD){
3038
3039                                         AliAODCaloCells &cells_aod=*(fAOD->GetEMCALCells());
3040
3041                                         TRefArray* caloClusters_aod = new TRefArray();
3042                                         fAOD->GetEMCALClusters(caloClusters_aod);
3043
3044                                         Int_t nclus_aod = caloClusters_aod->GetEntries();
3045
3046                                         for (Int_t icl = 0; icl < nclus_aod; icl++) {
3047
3048                                                 AliAODCaloCluster* clus_aod = (AliAODCaloCluster*)caloClusters_aod->At(icl);
3049
3050
3051                                                 if(clus_aod->IsEMCAL()){
3052                                                         Float_t ncells_aod      = fClus->GetNCells();
3053                                                         UShort_t * index_aod = clus_aod->GetCellsAbsId() ;
3054                                                         UShort_t * index2_aod = fClus->GetCellsAbsId() ;
3055
3056
3057                                                         for(Int_t i = 0; i < ncells_aod ; i++){
3058
3059                                                                 Int_t absId_aod =   index_aod[i];
3060                                                                 fTime->Fill(fPt,cells_aod.GetCellTime(absId_aod));
3061
3062                                                                 Int_t absId2_aod =   index2_aod[i];
3063                                                                 fTime2->Fill(fPt,cells_aod.GetCellTime(absId2_aod));
3064                                                         }
3065
3066                                                 }
3067                                         }
3068
3069                                 }
3070                                  */
3071
3072
3073                                 if(fUseEMCal){
3074                                         double emctof = fClus->GetTOF();
3075                                         ftimingEle->Fill(fPt,emctof);
3076                                 }
3077                                         //________________________________________________________________________
3078
3079
3080
3081
3082                                 // Residuals
3083                                 Double_t Dx = fClus->GetTrackDx();
3084                                 Double_t Dz = fClus->GetTrackDz();
3085                                 Double_t R=TMath::Sqrt(Dx*Dx+Dz*Dz);
3086
3087                                 for(Int_t i = 0; i < 6; i++)
3088                                 {
3089                                         if(fPt>=fPtBin[i] && fPt<fPtBin[i+1])
3090                                         {
3091
3092                                                 fEta[i]->Fill(Dz);
3093                                                 fPhi[i]->Fill(Dx);
3094                                                 fR[i]->Fill(R);
3095                                         }
3096                                 }
3097
3098                                 if(TMath::Abs(fClus->GetTrackDx())<=fdPhiCut && TMath::Abs(fClus->GetTrackDz())<=fdEtaCut)
3099                                 {
3100                                         Float_t Energy  = fClus->E();
3101                                         Float_t EoverP  = Energy/track->P();
3102                                         Float_t M02     = fClus->GetM02();
3103                                         Float_t M20     = fClus->GetM20();
3104                                         Float_t ncells  = fClus->GetNCells();
3105
3106                                         //here numerator for track-matching efficiency
3107                                         if(fIsMC && fIsAOD && track->GetLabel()>=0)
3108                                         {
3109                                                 fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
3110                                                 Int_t pdg = fMCparticle->GetPdgCode();
3111
3112                                                         //
3113                                                 if(fMCparticle->IsPhysicalPrimary()){
3114
3115
3116                                                         if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax ){
3117
3118                                                                 Bool_t MotherFound = FindMother(track->GetLabel());
3119                                                                 if(MotherFound){
3120                                                                         fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
3121                                                                         if( TMath::Abs(pdg) == 11 && fMCparticleMother->GetPdgCode()!=22 ){
3122                                                                                 if(fIsHFE1) fPt_track_match_num->Fill(fMCparticle->Pt());
3123                                                                         }
3124                                                                 }
3125                                                         }
3126                                                 }
3127                                         }///until here
3128
3129
3130
3131                                         //----------------------------------------------------------------------------------------
3132                                         //
3133                                         //EtaCut electrons histogram
3134                                         //Shower Shape Cut
3135                                         if(track->Eta()>=fEtaCutMin && track->Eta()<=fEtaCutMax ){
3136
3137                                                 if(fUseShowerShapeCut){
3138                                                         if(M02 >= fM02CutMin && M02<=fM02CutMax && M20>=fM20CutMin && M20<=fM20CutMax){
3139                                                                 fEoverP_pt[2]->Fill(fPt,(fClus->E() / fP));
3140                                                                 fShowerShapeCut->Fill(M02,M20);
3141                                                                 //in order to check if there are exotic cluster in this selected cluster (27 may 2014)
3142                                                                 fNcells_energy_elec_selected->Fill(ncells,Energy);
3143
3144                                                                         //true electrons E/p shape
3145                                                                 if(fIsMC && fIsAOD && track->GetLabel()>=0){
3146                                                                         fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
3147                                                                         Int_t pdg = fMCparticle->GetPdgCode();
3148
3149                                                                         if( TMath::Abs(pdg) == 11){
3150                                                                                 fEoverP_pt_true_electrons->Fill(fPt,(fClus->E() / fP));
3151                                                                         }
3152                                                                 }
3153
3154
3155
3156                                                         }
3157
3158                                                 }
3159                                                 if(!fUseShowerShapeCut){
3160                                                         fEoverP_pt[2]->Fill(fPt,(fClus->E() / fP));
3161                                                         fShowerShapeCut->Fill(M02,M20);
3162                                                         fNcells_energy_elec_selected->Fill(ncells,Energy);
3163
3164                                                         //after cut -> Are they really electrons ? E/p shape
3165                                                         if(fIsMC && fIsAOD && track->GetLabel()>=0){
3166                                                                 fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
3167                                                                 Int_t pdg = fMCparticle->GetPdgCode();
3168
3169                                                                 if( TMath::Abs(pdg) == 11){
3170                                                                         fEoverP_pt_true_electrons->Fill(fPt,(fClus->E() / fP));
3171                                                                 }
3172                                                         }
3173
3174
3175                                                 }
3176                                                 if(fUseEMCal) fShowerShape_ele->Fill(M02,M20);
3177
3178                                                 //for shower shape cut studies - now with TPC PID Cut
3179                                                 if(fUseEMCal){
3180                                                         fShowerShapeM02_EoverP->Fill(M02,EoverP);
3181                                                         fShowerShapeM20_EoverP->Fill(M20,EoverP);
3182                                                 }
3183
3184                                         }
3185
3186                                         //----------------------------------------------------------------------------------------
3187
3188
3189
3190                                         // for Eta Phi distribution
3191                                         fClus->GetPosition(pos3);
3192                                         TVector3 vpos(pos3[0],pos3[1],pos3[2]);
3193                                         Double_t cphi = vpos.Phi();
3194                                         Double_t ceta = vpos.Eta();
3195                                         fEtaPhi[2]->Fill(cphi,ceta);
3196
3197
3198
3199                                         fTPCNcls_EoverP[2]->Fill(TPCNcls, EoverP);
3200
3201                                         for(Int_t i = 0; i < 6; i++)
3202                                         {
3203                                                 if(fPt>=fPtBin[i] && fPt<fPtBin[i+1])
3204                                                 {
3205
3206                                                         fR_EoverP[i]->Fill(R, EoverP);
3207                                                         fNcells[i]->Fill(ncells);
3208                                                         fNcells_EoverP[i]->Fill(EoverP, ncells);
3209                                                         fM02_EoverP[i]->Fill(M02,EoverP);
3210                                                         fM20_EoverP[i]->Fill(M20,EoverP);
3211                                                         fECluster_ptbins[i]->Fill(Energy);
3212                                                         fEoverP_ptbins[i]->Fill(EoverP);
3213
3214                                                         if((fClus->E() / fP) >= fEoverPCutMin && (fClus->E() / fP) <= fEoverPCutMax) {
3215                                                                 fNcells_electrons[i]->Fill(ncells);
3216                                                         }
3217
3218                                                         if(M02<0.5 && M20<0.3) {
3219                                                                 fEoverP_wSSCut[i]->Fill(EoverP);
3220                                                         }
3221                                                 }
3222                                         }
3223
3224                                         fNcells_pt->Fill(fPt, ncells);
3225
3226
3227                                         ////////////////////////////////////////////////////////////////////
3228                                         ///EMCal - Efficiency calculations
3229
3230
3231                                         if(fIsMC && fIsAOD && track->GetLabel()>=0)
3232                                         {
3233                                                 fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
3234                                                 Int_t pdg = fMCparticle->GetPdgCode();
3235
3236                                                 //
3237                                                 if(fMCparticle->IsPhysicalPrimary()){
3238
3239                                                         //Phi cut && fMCparticle->Phi()>=(TMath::Pi()*80/180) && fMCparticle->Phi()<=TMath::Pi()
3240                                                         if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax  ){
3241
3242                                                                 Bool_t MotherFound = FindMother(track->GetLabel());
3243                                                                 if(MotherFound){
3244                                                                         fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
3245                                                                         if( TMath::Abs(pdg) == 11 && fMCparticleMother->GetPdgCode()!=22 ){
3246                                                                                 if(fIsHFE1)fPtMC_EMCal_All->Fill(fMCparticle->Pt());
3247                                                                         }
3248                                                                 }
3249                                                         }
3250                                                 }
3251                                         }
3252
3253                                         else if(fIsMC && track->GetLabel()>=0)//ESD
3254                                         {
3255
3256                                                 if(fMCstack->IsPhysicalPrimary(track->GetLabel()))
3257                                                 {
3258
3259                                                         fMCtrack = fMCstack->Particle(track->GetLabel());
3260
3261                                                         Int_t pdg = fMCtrack->GetPdgCode();
3262                                                         //Phi cut && fMCtrack->Phi()>=(TMath::Pi()*80/180) && fMCtrack->Phi()<=TMath::Pi()
3263                                                         if(fMCtrack->Eta()>=fEtaCutMin && fMCtrack->Eta()<=fEtaCutMax )
3264                                                         {
3265                                                                 Bool_t MotherFound = FindMother(track->GetLabel());
3266                                                                         //MotherFound included
3267                                                                 if(MotherFound){
3268                                                                         if(fMCtrack->GetFirstMother()>0){
3269                                                                                 fMCtrackMother = fMCstack->Particle(fMCtrack->GetFirstMother());
3270                                                                                 if( TMath::Abs(pdg) == 11 && fMCtrackMother->GetPdgCode()!=22 ){
3271                                                                                         if(fIsHFE1)fPtMC_EMCal_All->Fill(fMCtrack->Pt());
3272                                                                                 }
3273                                                                         }
3274                                                                 }
3275                                                         }
3276                                                 }
3277                                         }
3278
3279                                         //_______________________________________________________
3280                                         //PID using EMCal
3281                                         if(fEoverPnsigma)SetEoverPCutPtDependentMC(track->Pt());
3282
3283                                         if((fClus->E() / fP) >= fEoverPCutMin && (fClus->E() / fP) <= fEoverPCutMax)
3284                                         {
3285
3286
3287                                                         fECluster[2]->Fill(Energy);
3288                                                         fTPCNcls_pid[3]->Fill(TPCNcls, TPCNcls_pid);
3289
3290
3291                                                 if(fUseEMCal)
3292                                                 {
3293                                                         if(track->Eta()>=fEtaCutMin && track->Eta()<=fEtaCutMax ){
3294                                                                 fPtElec_Inc->Fill(fPt);
3295                                                                 //eta phi distribution for data
3296                                                                 fEtaPhi_data->Fill(track->Phi(),track->Eta());
3297                                                         }
3298
3299                                                         //Eta cut for background
3300                                                         if(fFillBackground){
3301                                                                 fEtad[2]->Fill(track->Eta());
3302
3303                                                                 //background for triggered data: trigger electron must have same cuts on shower shape  06/Jan/2014
3304                                                                 if(fUseShowerShapeCut){
3305                                                                         if(M02 >= fM02CutMin && M02<=fM02CutMax && M20>=fM20CutMin && M20<=fM20CutMax){
3306                                                                                 if(track->Eta()>=fEtaCutMin && track->Eta()<=fEtaCutMax){
3307                                                                                         Background(track, iTracks, Vtrack, kFALSE);
3308                                                                                 }
3309                                                                         }
3310                                                                 }
3311                                                                 else{
3312                                                                         if(track->Eta()>=fEtaCutMin && track->Eta()<=fEtaCutMax){
3313                                                                                 Background(track, iTracks, Vtrack, kFALSE);
3314                                                                         }
3315                                                                 }
3316
3317                                                         }
3318
3319                                                         double emctof2 = fClus->GetTOF();
3320                                                         ftimingEle2->Fill(fPt,emctof2);
3321                                                         //Correlation Analysis
3322                                                         if(fCorrelationFlag)
3323                                                         {
3324                                                                 ElectronHadronCorrelation(track, iTracks, Vtrack);
3325                                                         }
3326                                                 }
3327                                                 //_______________________________________________________
3328
3329                                                 ////////////////////////////////////////////////////////////////////
3330                                                 ///EMCal - efficiency calculations
3331
3332                                                 if(track->Charge()<0)  fCharge_n->Fill(fPt);
3333                                                 if(track->Charge()>0)  fCharge_p->Fill(fPt);
3334
3335
3336
3337                                                 if(fIsMC && fIsAOD && track->GetLabel()>=0)//AOD
3338                                                 {
3339                                                         if(track->Charge()<0)  fCharge_n->Fill(fPt);
3340                                                         if(track->Charge()>0)  fCharge_p->Fill(fPt);
3341
3342                                                         fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
3343                                                         if(fMCparticle->GetMother()>0) fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
3344                                                         Int_t pdg = fMCparticle->GetPdgCode();
3345
3346                                                         double proX = fMCparticle->Xv();
3347                                                         double proY = fMCparticle->Yv();
3348                                                         double proR = sqrt(pow(proX,2)+pow(proY,2));
3349
3350
3351                                                         if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax && fMCparticle->Phi()>=(TMath::Pi()*80/180) && fMCparticle->Phi()<=TMath::Pi()  ){
3352
3353                                                                 if( TMath::Abs(pdg) == 11 && fMCparticle->GetMother()>0 ){
3354                                                                         Int_t mpdg = fMCparticleMother->GetPdgCode();
3355                                                                         if(TMath::Abs(mpdg) == 221 || TMath::Abs(mpdg) == 22 || TMath::Abs(mpdg) == 111){
3356                                                                                 if(proR<7)fPtMCelectronAfterAll_nonPrimary->Fill(fMCparticle->Pt()); //numerator for the total efficiency, non Primary track
3357                                                                         }
3358                                                                 }
3359                                                                 if( TMath::Abs(pdg) == 11 && fMCparticle->IsPhysicalPrimary()) fPtMCelectronAfterAll_Primary->Fill(fMCparticle->Pt());
3360                                                         }
3361
3362
3363                                                         //
3364                                                         if(fMCparticle->IsPhysicalPrimary()){
3365
3366
3367                                                                 if(fMCparticle->Eta()>=fEtaCutMin && fMCparticle->Eta()<=fEtaCutMax ){
3368
3369                                                                         Bool_t MotherFound = FindMother(track->GetLabel());
3370                                                                         if(MotherFound){
3371
3372                                                                                 if(!fUseShowerShapeCut){
3373                                                                                         if(fIsHFE1){
3374                                                                                                 //Unfolding   pt_reco/pt_MC  in the efficiency
3375                                                                                                 fPtMCelectronAfterAll->Fill(fMCparticle->Pt());
3376                                                                                                 fPtMCelectronAfterAll_unfolding->Fill(track->Pt());
3377                                                                                                 fEtaPhi_num->Fill(fMCparticle->Phi(),fMCparticle->Eta());
3378
3379                                                                                                 //new histo to estimate how different is pt reco from pt MC
3380                                                                                                 fpt_reco_pt_MC_num->Fill(track->Pt(),fMCparticle->Pt());
3381                                                                                         }//numerator for the total efficiency AOD
3382                                                                                 }
3383                                                                                         //November 11 for efficiency of triggered data
3384                                                                                 if(fUseShowerShapeCut){
3385                                                                                         if(M02 >= fM02CutMin && M02<=fM02CutMax && M20>=fM20CutMin && M20<=fM20CutMax){
3386                                                                                                 if(fIsHFE1){
3387                                                                                                         //Unfolding   pt_reco/pt_MC  in the efficiency
3388                                                                                                         fPtMCelectronAfterAll->Fill(fMCparticle->Pt());
3389                                                                                                         fPtMCelectronAfterAll_unfolding->Fill(track->Pt());
3390                                                                                                         fEtaPhi_num->Fill(fMCparticle->Phi(),fMCparticle->Eta());
3391
3392                                                                                                         //new histo to estimate how different is pt reco from pt MC
3393                                                                                                         fpt_reco_pt_MC_num->Fill(track->Pt(),fMCparticle->Pt());
3394                                                                                                 }//numerator for the total efficiency AOD
3395                                                                                         }
3396                                                                                 }
3397
3398
3399
3400                                                                                         fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
3401                                                                                         if( TMath::Abs(pdg) == 11 && fMCparticleMother->GetPdgCode()!=22 ){
3402                                                                                                 if(fIsHFE1)fPtMC_EMCal_Selected->Fill(fMCparticle->Pt());
3403                                                                                         }
3404                                                                         }//if MotherFound
3405                                                                 }//eta cut
3406                                                         }
3407                                                 }///close AOD
3408
3409                                                 else if(fIsMC && track->GetLabel()>=0)//ESD
3410                                                 {
3411                                                         if(track->Charge()<0)  fCharge_n->Fill(fPt);
3412                                                         if(track->Charge()>0)  fCharge_p->Fill(fPt);
3413
3414                                                         fMCtrack = fMCstack->Particle(track->GetLabel());
3415                                                         if(fMCtrack->GetFirstMother()>0)
3416                                                         {
3417                                                                 fMCtrackMother = fMCstack->Particle(fMCtrack->GetFirstMother());
3418                                                         }
3419                                                         TParticle *particle=fMCstack->Particle(track->GetLabel());
3420
3421                                                         Int_t pdg = fMCtrack->GetPdgCode();
3422
3423
3424                                                         if(fMCtrack->Eta()>=fEtaCutMin && fMCtrack->Eta()<=fEtaCutMax)
3425                                                         {
3426                                                                 if( TMath::Abs(pdg) == 11 && fMCtrack->GetFirstMother()>0 )
3427                                                                 {
3428                                                                         Int_t mpdg = fMCtrackMother->GetPdgCode();
3429                                                                         if(TMath::Abs(mpdg) == 221 || TMath::Abs(mpdg) == 22 || TMath::Abs(mpdg) == 111)
3430                                                                         {
3431                                                                                 Double_t proR=particle->R();
3432                                                                                 if(proR<7)
3433                                                                                 {
3434                                                                                   fPtMCelectronAfterAll_nonPrimary->Fill(fMCtrack->Pt()); //numerator for the total efficiency, non Primary track
3435                                                                                 }
3436                                                                         }
3437                                                                 }
3438                                                                 if( TMath::Abs(pdg) == 11 && fMCstack->IsPhysicalPrimary(track->GetLabel()))
3439                                                                 {
3440                                                                         fPtMCelectronAfterAll_Primary->Fill(fMCtrack->Pt());
3441                                                                 }
3442                                                         }
3443
3444                                                         if(fMCstack->IsPhysicalPrimary(track->GetLabel()))
3445                                                         {
3446                                                                 Bool_t MotherFound = FindMother(track->GetLabel());
3447
3448                                                                 if(MotherFound)
3449                                                                 {
3450                                                                         if(fMCtrack->Eta()>=fEtaCutMin && fMCtrack->Eta()<=fEtaCutMax)
3451                                                                         {
3452                                                                                 if(!fUseShowerShapeCut){
3453                                                                                         if(fIsHFE1){
3454                                                                                                 fPtMCelectronAfterAll->Fill(fMCtrack->Pt()); //numerator for the total efficiency ESD
3455                                                                                                 fEtaPhi_num->Fill(fMCtrack->Phi(),fMCtrack->Eta());
3456                                                                                         }
3457                                                                                 }
3458                                                                                 //November 11 for efficiency of triggered data
3459                                                                                 if(fUseShowerShapeCut){
3460                                                                                         if(M02 >= fM02CutMin && M02<=fM02CutMax && M20>=fM20CutMin && M20<=fM20CutMax){
3461                                                                                                 if(fIsHFE1){
3462                                                                                                         fPtMCelectronAfterAll->Fill(fMCtrack->Pt()); //numerator for the total efficiency ESD
3463                                                                                                         fEtaPhi_num->Fill(fMCtrack->Phi(),fMCtrack->Eta());
3464                                                                                                 }
3465                                                                                         }
3466                                                                                 }
3467
3468
3469
3470                                                                         }
3471                                                                 }
3472
3473
3474
3475                                                                 // Phi cut && fMCtrack->Phi()>=(TMath::Pi()*80/180) && fMCtrack->Phi()<=TMath::Pi()
3476                                                                 if(fMCtrack->Eta()>=fEtaCutMin && fMCtrack->Eta()<=fEtaCutMax)
3477                                                                 {
3478                                                                         //included MotherFound
3479
3480                                                                         if(MotherFound)
3481                                                                         {
3482                                                                                 if(fMCtrack->GetFirstMother()>0){
3483                                                                                         fMCtrackMother = fMCstack->Particle(fMCtrack->GetFirstMother());
3484                                                                                         if( TMath::Abs(pdg) == 11 && fMCtrackMother->GetPdgCode()!=22 ){
3485
3486                                                                                                 if(fIsHFE1){fPtMC_EMCal_Selected->Fill(fMCtrack->Pt());}
3487                                                                                         }
3488                                                                                 }
3489                                                                         }
3490                                                                 }
3491                                                         }
3492                                                 }//close ESD
3493                                                 ///////////////////////////////////////////////////////////////////
3494
3495
3496                                         }
3497                                 }
3498                         }
3499                 }
3500
3501                         //______________________________________________________________
3502                         // Vertex
3503
3504                 fVtxZ[2]->Fill(fZvtx);
3505                 if(iTracks == 0)fNTracks[2]->Fill(fNOtrks);
3506                 fNTracks_pt[2]->Fill(fNOtrks, fPt);
3507                 fNTracks_eta[2]->Fill(fNOtrks, track->Eta());
3508                 fNTracks_phi[2]->Fill(fNOtrks, track->Phi());
3509                 fNClusters[2]->Fill(ClsNo);
3510                 fTPCNcls_pid[2]->Fill(TPCNcls, TPCNcls_pid);
3511
3512                         //______________________________________________________________
3513
3514                         //_______________________________________________________
3515                         //Correlation Analysis
3516                 if(!fUseEMCal)
3517                 {
3518                         fPtElec_Inc->Fill(fPt);
3519
3520                         if(fCorrelationFlag)
3521                         {
3522                                 ElectronHadronCorrelation(track, iTracks, Vtrack);
3523                         }
3524                 }
3525                         //_______________________________________________________
3526
3527                         ///________________________________________________________________________
3528         }
3529
3530                 //__________________________________________________________________
3531                 //Event Mixing Analysis
3532                 //Filling pool
3533         if(fEventMixingFlag)
3534         {
3535                 fPool = fPoolMgr->GetEventPool(fCentrality->GetCentralityPercentile("V0A"), fZvtx); // Get the buffer associated with the current centrality and z-vtx
3536
3537                 if(!fPool) AliFatal(Form("No pool found for centrality = %f, zVtx = %f", fCentrality->GetCentralityPercentile("V0A"), fZvtx));
3538
3539                 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!
3540
3541
3542         }
3543
3544                 //__________________________________________________________________
3545
3546         delete fListOfmotherkink;
3547         PostData(1, fOutputList);
3548 }
3549
3550         //______________________________________________________________________
3551 void AliAnalysisTaskEMCalHFEpA::Terminate(Option_t *)
3552 {
3553                 //Draw result to the screen
3554                 //Called once at the end of the query
3555
3556         fOutputList = dynamic_cast<TList*> (GetOutputData(1));
3557
3558         if(!fOutputList)
3559         {
3560                 printf("ERROR: Output list not available\n");
3561                 return;
3562         }
3563 }
3564
3565         //______________________________________________________________________
3566 Bool_t AliAnalysisTaskEMCalHFEpA::ProcessCutStep(Int_t cutStep, AliVParticle *track)
3567 {
3568                 //Check single track cuts for a given cut step
3569                 //Note this function is called inside the UserExec function
3570         const Int_t kMCOffset = AliHFEcuts::kNcutStepsMCTrack;
3571         if(!fCFM->CheckParticleCuts(cutStep + kMCOffset, track)) return kFALSE;
3572         return kTRUE;
3573 }
3574
3575
3576         //______________________________________________________________________
3577
3578
3579 void AliAnalysisTaskEMCalHFEpA::Background(AliVTrack *track, Int_t trackIndex, AliVParticle *vtrack, Bool_t IsTPConly)
3580 {
3581         ///_________________________________________________________________
3582         ///MC analysis
3583         //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)
3584
3585                 if(fIsMC)
3586                 {
3587                         if(track->GetLabel() < 0)
3588                         {
3589                                 AliWarning(Form("The track %d does not have a valid MC label",trackIndex));
3590                                 return;
3591                         }
3592
3593                         if(fIsAOD)
3594                         {
3595                                 fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
3596
3597                                 if(fMCparticle->GetMother()<0) return;
3598
3599                                 fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
3600                                 if(fMCparticleMother->GetMother()>0)fMCparticleGMother = (AliAODMCParticle*) fMCarray->At(fMCparticleMother->GetMother());
3601
3602                                 if(TMath::Abs(fMCparticle->GetPdgCode())==11 && (TMath::Abs(fMCparticleMother->GetPdgCode())==22 || TMath::Abs(fMCparticleMother->GetPdgCode())==111 || TMath::Abs(fMCparticleMother->GetPdgCode())==221))
3603                                 {
3604                                                 //Is Background
3605                                         if(!IsTPConly)fPtBackgroundBeforeReco->Fill(track->Pt());
3606                                         if(IsTPConly)fPtBackgroundBeforeReco2->Fill(track->Pt());
3607
3608
3609                                                 //October 08th weighted histos
3610                                         if(TMath::Abs(fMCparticleMother->GetPdgCode())==111 || TMath::Abs(fMCparticleMother->GetPdgCode())==221 ){
3611
3612                                                 Double_t mPt=fMCparticleMother->Pt();
3613                                                 Double_t mweight=1;
3614
3615
3616                                                 //________________________________________________________________
3617                                                 //correction for d3 based on data
3618
3619                                                         if(TMath::Abs(fMCparticleMother->GetPdgCode())==111){
3620                                                                 Double_t x=mPt;
3621                                                                 mweight=CalculateWeight(111, x);
3622
3623                                                         }
3624                                                     if(TMath::Abs(fMCparticleMother->GetPdgCode())==221){
3625                                                             Double_t x=mPt;
3626                                                             mweight=CalculateWeight(221, x);
3627
3628                                                     }
3629
3630                                                 //________________________________________________________________
3631
3632                                                 //Histo pT mother versus pT electron
3633                                                 fpT_m_electron->Fill(mPt, track->Pt());
3634
3635                                                 if(!IsTPConly)fPtBackgroundBeforeReco_weight->Fill(track->Pt(), 1./mweight);
3636                                                 if(IsTPConly)fPtBackgroundBeforeReco2_weight->Fill(track->Pt(), 1./mweight);
3637                                         }
3638                                         else if(fMCparticleMother->GetMother()>0 && (TMath::Abs(fMCparticleGMother->GetPdgCode())==111 || TMath::Abs(fMCparticleGMother->GetPdgCode())==221 )){
3639
3640                                                 Double_t gmPt=fMCparticleGMother->Pt();
3641                                                 Double_t gmweight=1;
3642
3643                                                 //________________________________________________________________
3644                                                 //correction for d3 based on data
3645
3646                                                         if(TMath::Abs(fMCparticleGMother->GetPdgCode())==111){
3647                                                                 Double_t x=gmPt;
3648                                                                 gmweight=CalculateWeight(111, x);
3649                                                         }
3650                                                         if(TMath::Abs(fMCparticleGMother->GetPdgCode())==221){
3651                                                                 Double_t x=gmPt;
3652                                                                 gmweight=CalculateWeight(221, x);
3653                                                         }
3654
3655
3656                                                 //________________________________________________________________
3657                                                 //Histo pT gmother versus pT electron
3658
3659                                                 fpT_gm_electron->Fill(gmPt, track->Pt());
3660
3661                                                 if(!IsTPConly)fPtBackgroundBeforeReco_weight->Fill(track->Pt(), 1./gmweight);
3662                                                 if(IsTPConly)fPtBackgroundBeforeReco2_weight->Fill(track->Pt(), 1./gmweight);
3663                                         }
3664                                         else{
3665                                                 if(!IsTPConly)fPtBackgroundBeforeReco_weight->Fill(track->Pt());
3666                                                 if(IsTPConly)fPtBackgroundBeforeReco2_weight->Fill(track->Pt());
3667                                         }
3668                                 }//particle kind
3669                         }//IsAOD
3670                          //ESD
3671                         else
3672                         {
3673                                 fMCtrack = fMCstack->Particle(track->GetLabel());
3674
3675                                 if(fMCtrack->GetFirstMother()<0) return;
3676
3677                                 fMCtrackMother = fMCstack->Particle(fMCtrack->GetFirstMother());
3678
3679                                 if(TMath::Abs(fMCtrack->GetPdgCode())==11 && (TMath::Abs(fMCtrackMother->GetPdgCode())==22 || TMath::Abs(fMCtrackMother->GetPdgCode())==111 || TMath::Abs(fMCtrackMother->GetPdgCode())==221))
3680                                 {
3681                                                 //Is Background
3682                                         if(!IsTPConly)fPtBackgroundBeforeReco->Fill(track->Pt());
3683                                         if(IsTPConly)fPtBackgroundBeforeReco2->Fill(track->Pt());
3684                                 }
3685                         }
3686                 }//IsMC
3687
3688                 ///_________________________________________________________________
3689
3690                 //________________________________________________
3691                 //Associated particle cut
3692                 fPartnerCuts->SetAcceptKinkDaughters(kFALSE);
3693                 fPartnerCuts->SetRequireITSRefit(kTRUE);
3694                 fPartnerCuts->SetRequireTPCRefit(kTRUE);
3695                 fPartnerCuts->SetEtaRange(-0.9,0.9);
3696                 fPartnerCuts->SetMaxChi2PerClusterTPC(4.0);
3697                 fPartnerCuts->SetMinNClustersTPC(80);
3698                 fPartnerCuts->SetPtRange(0,1e10);
3699                 //fPartnerCuts->SetRequireSigmaToVertex(kTRUE);
3700                 //fPartnerCuts->SetMaxDCAToVertexXY(1);
3701                 //fPartnerCuts->SetMaxDCAToVertexZ(3);
3702                 //_________________________________________________
3703
3704                 ///#################################################################
3705                 //Non-HFE reconstruction
3706                 fNonHFE = new AliSelectNonHFE();
3707                 fNonHFE->SetAODanalysis(fIsAOD);
3708                 if(fMassCutFlag) fNonHFE->SetInvariantMassCut(fMassCut);
3709                 if(fAngleCutFlag) fNonHFE->SetOpeningAngleCut(fAngleCut);
3710                 if(fChi2CutFlag) fNonHFE->SetChi2OverNDFCut(fChi2Cut);
3711                 if(fDCAcutFlag) fNonHFE->SetDCACut(fDCAcut);
3712                 fNonHFE->SetAlgorithm("DCA"); //KF
3713                 fNonHFE->SetPIDresponse(fPidResponse);
3714                 fNonHFE->SetTrackCuts(-3.5,3.5,fPartnerCuts);
3715                 fNonHFE->SetAdditionalCuts(fPtMinAsso,fTpcNclsAsso);
3716
3717                 if(!IsTPConly){
3718                         fNonHFE->SetHistAngleBack(fOpAngleBack);
3719                         fNonHFE->SetHistAngle(fOpAngle);
3720                         fNonHFE->SetHistDCABack(fDCABack);
3721                         fNonHFE->SetHistDCA(fDCA);
3722                         fNonHFE->SetHistMassBack(fInvMassBack);
3723                         fNonHFE->SetHistMass(fInvMass);
3724                 }
3725                 if(IsTPConly){
3726                         fNonHFE->SetHistAngleBack(fOpAngleBack2);
3727                         fNonHFE->SetHistAngle(fOpAngle2);
3728                         fNonHFE->SetHistDCABack(fDCABack2);
3729                         fNonHFE->SetHistDCA(fDCA2);
3730                         fNonHFE->SetHistMassBack(fInvMassBack2);
3731                         fNonHFE->SetHistMass(fInvMass2);
3732                 }
3733
3734                 fNonHFE->FindNonHFE(trackIndex,vtrack,fVevent);
3735
3736                         //index of track selected as partner
3737                 Int_t *fUlsPartner = fNonHFE->GetPartnersULS();
3738
3739                         //Electron Information
3740                 Double_t fPhiE = -999;
3741                 Double_t fEtaE = -999;
3742                 Double_t fPtE = -999;
3743                 fPhiE = track->Phi();
3744                 fEtaE = track->Eta();
3745                 fPtE = track->Pt();
3746
3747                         ///_________________________________________________________________
3748                         ///MC analysis
3749                 if(fIsMC)
3750                 {
3751                         if(fIsAOD)
3752                         {
3753                                 if(TMath::Abs(fMCparticle->GetPdgCode())==11 && (TMath::Abs(fMCparticleMother->GetPdgCode())==22 || TMath::Abs(fMCparticleMother->GetPdgCode())==111 || TMath::Abs(fMCparticleMother->GetPdgCode())==221))
3754                                 {
3755
3756                                         Double_t weight=1;
3757
3758                                         if(!IsTPConly){
3759                                                 if(fNonHFE->IsULS()) fPtElec_ULS->Fill(fPtE,fNonHFE->GetNULS());
3760                                                 if(fNonHFE->IsLS()) fPtElec_LS->Fill(fPtE,fNonHFE->GetNLS());
3761
3762
3763
3764
3765                                                 if(TMath::Abs(fMCparticleMother->GetPdgCode())==111 || TMath::Abs(fMCparticleMother->GetPdgCode())==221){
3766                                                         Double_t mPt=fMCparticleMother->Pt();
3767                                                                 Double_t mweight1=1;
3768                                                                 Double_t mweight2=1;
3769                                                                 //Double_t weight=1;
3770
3771
3772                                                          //----------------------------------------------------------------------------
3773                                                          //correction based on data only
3774                                                                 if(TMath::Abs(fMCparticleMother->GetPdgCode())==111){
3775                                                                         Double_t x=mPt;
3776                                                                         weight=CalculateWeight(111, x);
3777                                                                 }
3778                                                                 if(TMath::Abs(fMCparticleMother->GetPdgCode())==221){
3779                                                                         Double_t x=mPt;
3780                                                                         weight=CalculateWeight(221, x);
3781                                                                 }
3782
3783
3784                                                                 //----------------------------------------------------------------------------
3785
3786                                                                 //check this
3787                                                         if(fNonHFE->IsULS()) mweight1=(fNonHFE->GetNULS())/weight;
3788                                                         if(fNonHFE->IsLS())  mweight2=(fNonHFE->GetNLS())/weight;
3789
3790                                                                 //fill histos
3791                                                         if(fNonHFE->IsULS())fPtElec_ULS_weight->Fill(fPtE, mweight1);
3792                                                         if(fNonHFE->IsLS())fPtElec_LS_weight->Fill(fPtE, mweight2);
3793                                                 }
3794                                                 else if(fMCparticleMother->GetMother()>0 && (TMath::Abs(fMCparticleGMother->GetPdgCode())==111 || TMath::Abs(fMCparticleGMother->GetPdgCode())==221 )){
3795                                                         Double_t gmPt=fMCparticleGMother->Pt();
3796                                                                 Double_t gmweight1=1;
3797                                                                 Double_t gmweight2=1;
3798                                                                 //Double_t weight=1;
3799
3800                                                         //----------------------------------------------------------------------------
3801
3802                                                         //----------------------------------------------------------------------------
3803
3804                                                         //correction based on data only for pi0
3805
3806                                                                 if(TMath::Abs(fMCparticleGMother->GetPdgCode())==111){
3807                                                                         Double_t x=gmPt;
3808                                                                         weight=CalculateWeight(111, x);
3809                                                                 }
3810                                                                 if(TMath::Abs(fMCparticleGMother->GetPdgCode())==221){
3811                                                                         Double_t x=gmPt;
3812                                                                         weight=CalculateWeight(221, x);
3813                                                         }
3814
3815
3816
3817
3818
3819                                                                 //check this
3820                                                         if(fNonHFE->IsULS()) gmweight1=(fNonHFE->GetNULS())/weight;
3821                                                         if(fNonHFE->IsLS())  gmweight2=(fNonHFE->GetNLS())/weight;
3822
3823                                                                 //fill histos
3824                                                         if(fNonHFE->IsULS())fPtElec_ULS_weight->Fill(fPtE, gmweight1);
3825                                                         if(fNonHFE->IsLS())fPtElec_LS_weight->Fill(fPtE, gmweight2);
3826                                                 }
3827                                                 else{
3828                                                         if(fNonHFE->IsULS()) fPtElec_ULS_weight->Fill(fPtE,fNonHFE->GetNULS());
3829                                                         if(fNonHFE->IsLS()) fPtElec_LS_weight->Fill(fPtE,fNonHFE->GetNLS());
3830                                                 }
3831
3832
3833                                         }//!IsTPConly
3834
3835                                         if(IsTPConly){
3836                                                 if(fNonHFE->IsULS()) fPtElec_ULS2->Fill(fPtE,fNonHFE->GetNULS());
3837                                                 if(fNonHFE->IsLS()) fPtElec_LS2->Fill(fPtE,fNonHFE->GetNLS());
3838
3839                                                         //new 08 October        //weighted histograms
3840                                                 if(TMath::Abs(fMCparticleMother->GetPdgCode())==111 || TMath::Abs(fMCparticleMother->GetPdgCode())==221){
3841                                                         Double_t mPt=fMCparticleMother->Pt();
3842
3843                                                                 Double_t mweight1=1;
3844                                                                 Double_t mweight2=1;
3845                                                                 //Double_t weight=1;
3846
3847
3848                                                          //----------------------------------------------------------------------------
3849
3850                                                         //correction based on data for d3
3851
3852                                                                 if(TMath::Abs(fMCparticleMother->GetPdgCode())==111){
3853                                                                         Double_t x=mPt;
3854                                                                         weight=CalculateWeight(111, x);
3855
3856                                                                 }
3857                                                                 if(TMath::Abs(fMCparticleMother->GetPdgCode())==221){
3858                                                                         Double_t x=mPt;
3859                                                                         weight=CalculateWeight(221, x);
3860
3861                                                                 }
3862
3863
3864                                                                 //check this
3865                                                         if(fNonHFE->IsULS()) mweight1=(fNonHFE->GetNULS())/weight;
3866                                                         if(fNonHFE->IsLS())  mweight2=(fNonHFE->GetNLS())/weight;
3867
3868                                                                 //fill histos
3869                                                         if(fNonHFE->IsULS())fPtElec_ULS2_weight->Fill(fPtE, mweight1);
3870                                                         if(fNonHFE->IsLS())fPtElec_LS2_weight->Fill(fPtE, mweight2);
3871                                                 }
3872                                                 else if(fMCparticleMother->GetMother()>0 && (TMath::Abs(fMCparticleGMother->GetPdgCode())==111 || TMath::Abs(fMCparticleGMother->GetPdgCode())==221 )){
3873                                                         Double_t gmPt=fMCparticleGMother->Pt();
3874                                                         Double_t gmweight1=1;
3875                                                         Double_t gmweight2=1;
3876                                                                 //Double_t weight=1;
3877
3878
3879                                                         //----------------------------------------------------------------------------
3880                                                          //correction based on data only for pi0
3881
3882                                                                 if(TMath::Abs(fMCparticleGMother->GetPdgCode())==111){
3883                                                                         Double_t x=gmPt;
3884                                                                         weight=CalculateWeight(111, x);
3885                                                                 }
3886
3887
3888                                                                 if(TMath::Abs(fMCparticleGMother->GetPdgCode())==221){
3889                                                                         Double_t x=gmPt;
3890                                                                         weight=CalculateWeight(221, x);
3891                                                                 }
3892
3893                                                         //----------------------------------------------------------------------------
3894
3895
3896
3897                                                                 //check this
3898                                                         if(fNonHFE->IsULS()) gmweight1=(fNonHFE->GetNULS())/weight;
3899                                                         if(fNonHFE->IsLS())  gmweight2=(fNonHFE->GetNLS())/weight;
3900
3901                                                                 //fill histos
3902                                                         if(fNonHFE->IsULS())fPtElec_ULS2_weight->Fill(fPtE, gmweight1);
3903                                                         if(fNonHFE->IsLS())fPtElec_LS2_weight->Fill(fPtE, gmweight2);
3904                                                 }
3905                                                 else{
3906                                                         if(fNonHFE->IsULS()) fPtElec_ULS2_weight->Fill(fPtE,fNonHFE->GetNULS());
3907                                                         if(fNonHFE->IsLS()) fPtElec_LS2_weight->Fill(fPtE,fNonHFE->GetNLS());
3908                                                 }
3909
3910
3911                                                 //----------------------------------------------------------------------------
3912                                                 //to check other way to calculate efficiency
3913                                                 //ULS with no weight from ULS-LS original
3914                                                 //we have to know if track2 comes from same mother!!!
3915                                                 //----------------------------------------------------------------------------
3916                                                 if(fNonHFE->IsULS()){
3917
3918                                                         for(Int_t iTracks = 0; iTracks < fVevent->GetNumberOfTracks(); iTracks++)
3919                                                         {
3920
3921                                                                 AliVParticle* Vtrack2 = fVevent->GetTrack(iTracks);
3922                                                                 if (!Vtrack2)
3923                                                                 {
3924                                                                         printf("ERROR: Could not receive track %d\n", iTracks);
3925                                                                         continue;
3926                                                                 }
3927                                                                 AliVTrack *track2 = dynamic_cast<AliVTrack*>(Vtrack2);
3928                                                                 if(track2->GetLabel()<0) continue;
3929                                                                 fMCparticle2 = (AliAODMCParticle*) fMCarray->At(track2->GetLabel());
3930                                                                 if(fMCparticle2->GetMother()<0) continue;
3931
3932                                                                 for(Int_t i = 0; i < fNonHFE->GetNULS(); i++)
3933                                                                 {
3934                                                                         if(fUlsPartner[i]==iTracks){
3935                                                                                         //only fill if it has same mother
3936                                                                                         //with weight to take into account the number of partners
3937                                                                                 if(fMCparticle2->GetMother()==fMCparticle->GetMother()) fPtElec_ULS_MC->Fill(fPtE, fNonHFE->GetNULS());
3938
3939                                                                                 //-----------------------------------------------------------------------------------------------------------
3940                                                                                 //weight for mother
3941                                                                                 //Double_t weight2=1;
3942                                                                                 Double_t mPt=fMCparticleMother->Pt();
3943
3944
3945                                                                                 if(TMath::Abs(fMCparticleMother->GetPdgCode())==111){
3946                                                                                         Double_t x=mPt;
3947                                                                                         weight=CalculateWeight(111, x);
3948
3949
3950                                                                                 }
3951                                                                                 if(TMath::Abs(fMCparticleMother->GetPdgCode())==221){
3952                                                                                         Double_t x=mPt;
3953                                                                                         weight=CalculateWeight(221, x);
3954
3955
3956                                                                                 }
3957
3958                                                                                 //weight for grandmother
3959
3960                                                                                 if((fMCparticleMother->GetMother()>0) && TMath::Abs(((fMCparticleGMother->GetPdgCode())==111))){
3961                                                                                         Double_t gmPt=fMCparticleGMother->Pt();
3962                                                                                         Double_t x=gmPt;
3963                                                                                         weight=CalculateWeight(111, x);
3964
3965                                                                                 }
3966                                                                                 if((fMCparticleMother->GetMother()>0) && TMath::Abs(((fMCparticleGMother->GetPdgCode())==221))){
3967                                                                                         Double_t gmPt=fMCparticleGMother->Pt();
3968                                                                                         Double_t x=gmPt;
3969                                                                                         weight=CalculateWeight(221, x);
3970
3971                                                                                 }
3972
3973
3974                                                                                 if(fMCparticle2->GetMother()==fMCparticle->GetMother()) fPtElec_ULS_MC_weight->Fill(fPtE, (fNonHFE->GetNULS())*1./weight);
3975
3976                                                                                 //-----------------------------------------------------------------------------------------------------------
3977                                                                                 //end of weight
3978
3979                                                                         }//partner found same as track
3980                                                                 }//loop in all partner
3981
3982                                                         }//track
3983                                                 }//is ULS
3984                                                 //----------------------------------------------------------------------------
3985                                                 //end of part to check other way to calculate efficiency
3986                                                 //----------------------------------------------------------------------------
3987
3988                                         }//IsTPConly
3989
3990                                 }//particle kind
3991
3992                                 if(IsTPConly){
3993                                                 //ULS-LS with no pid AOD
3994                                         if(fNonHFE->IsULS()) fPtElec_ULS_NoPid->Fill(fPtE,fNonHFE->GetNULS());
3995                                         if(fNonHFE->IsLS()) fPtElec_LS_NoPid->Fill(fPtE,fNonHFE->GetNLS());
3996                                 }
3997
3998                         }//close IsAOD
3999                          //It is ESD
4000                         else
4001                         {
4002                                 if(TMath::Abs(fMCtrack->GetPdgCode())==11 && (TMath::Abs(fMCtrackMother->GetPdgCode())==22 || TMath::Abs(fMCtrackMother->GetPdgCode())==111 || TMath::Abs(fMCtrackMother->GetPdgCode())==221))
4003                                 {
4004                                         if(!IsTPConly){
4005                                                 if(fNonHFE->IsULS()) fPtElec_ULS->Fill(fPtE,fNonHFE->GetNULS());
4006                                                 if(fNonHFE->IsLS()) fPtElec_LS->Fill(fPtE,fNonHFE->GetNLS());
4007                                         }
4008
4009                                         if(IsTPConly){
4010                                                 if(fNonHFE->IsULS()) fPtElec_ULS2->Fill(fPtE,fNonHFE->GetNULS());
4011                                                 if(fNonHFE->IsLS()) fPtElec_LS2->Fill(fPtE,fNonHFE->GetNLS());
4012                                         }
4013                                 }
4014
4015
4016                         }
4017                 }//close IsMC
4018                  ///_________________________________________________________________
4019                  //not MC
4020                 else
4021                 {
4022                         if(!IsTPConly){
4023                                 if(fNonHFE->IsULS()) fPtElec_ULS->Fill(fPtE,fNonHFE->GetNULS());
4024                                 if(fNonHFE->IsLS()) fPtElec_LS->Fill(fPtE,fNonHFE->GetNLS());
4025                         }
4026
4027                         if(IsTPConly){
4028                                 if(fNonHFE->IsULS()) fPtElec_ULS2->Fill(fPtE,fNonHFE->GetNULS());
4029                                 if(fNonHFE->IsLS()) fPtElec_LS2->Fill(fPtE,fNonHFE->GetNLS());
4030                         }
4031                 }
4032
4033         //for MC closure test
4034         // to be used as data, with MC input
4035         if(!IsTPConly){
4036                 if(fNonHFE->IsULS()) fPtElec_ULS_mc_closure->Fill(fPtE,fNonHFE->GetNULS());
4037                 if(fNonHFE->IsLS()) fPtElec_LS_mc_closure->Fill(fPtE,fNonHFE->GetNLS());
4038         }
4039
4040         if(IsTPConly){
4041                 if(fNonHFE->IsULS()) fPtElec_ULS2_mc_closure->Fill(fPtE,fNonHFE->GetNULS());
4042                 if(fNonHFE->IsLS()) fPtElec_LS2_mc_closure->Fill(fPtE,fNonHFE->GetNLS());
4043         }
4044
4045
4046
4047 }//end of Background function
4048
4049         //______________________________________________________________________
4050 void AliAnalysisTaskEMCalHFEpA::ElectronHadronCorrelation(AliVTrack *track, Int_t trackIndex, AliVParticle *vtrack)
4051 {
4052
4053         ///_________________________________________________________________
4054         ///MC analysis
4055         if(fIsMC)
4056         {
4057                 if(track->GetLabel() < 0)
4058         {
4059                         AliWarning(Form("The track %d does not have a valid MC label",trackIndex));
4060                         return;
4061         }
4062
4063                 if(fIsAOD)
4064                 {
4065                         fMCparticle = (AliAODMCParticle*) fMCarray->At(track->GetLabel());
4066
4067                         if(fMCparticle->GetMother()<0) return;
4068
4069                 fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
4070
4071                 if(TMath::Abs(fMCparticle->GetPdgCode())==11 && (TMath::Abs(fMCparticleMother->GetPdgCode())==22 || TMath::Abs(fMCparticleMother->GetPdgCode())==111 || TMath::Abs(fMCparticleMother->GetPdgCode())==221))
4072                 {
4073                                         //Is Background
4074                                 fPtBackgroundBeforeReco->Fill(track->Pt());
4075                         }
4076                 }
4077                 else
4078                 {
4079                 fMCtrack = fMCstack->Particle(track->GetLabel());
4080
4081                 if(fMCtrack->GetFirstMother()<0) return;
4082
4083                 fMCtrackMother = fMCstack->Particle(fMCtrack->GetFirstMother());
4084
4085                 if(TMath::Abs(fMCtrack->GetPdgCode())==11 && (TMath::Abs(fMCtrackMother->GetPdgCode())==22 || TMath::Abs(fMCtrackMother->GetPdgCode())==111 || TMath::Abs(fMCtrackMother->GetPdgCode())==221))
4086                 {
4087                                         //Is Background
4088                                 fPtBackgroundBeforeReco->Fill(track->Pt());
4089                         }
4090                 }
4091         }
4092                 ///_________________________________________________________________
4093
4094                 //________________________________________________
4095                 //Associated particle cut
4096         fPartnerCuts->SetAcceptKinkDaughters(kFALSE);
4097     fPartnerCuts->SetRequireITSRefit(kTRUE);
4098     fPartnerCuts->SetRequireTPCRefit(kTRUE);
4099     fPartnerCuts->SetEtaRange(-0.9,0.9);
4100     fPartnerCuts->SetMaxChi2PerClusterTPC(4.0);
4101     fPartnerCuts->SetMinNClustersTPC(80);
4102     fPartnerCuts->SetPtRange(0.3,1e10);
4103                 //fPartnerCuts->SetRequireSigmaToVertex(kTRUE);
4104                 //fPartnerCuts->SetMaxDCAToVertexXY(1);
4105                 //fPartnerCuts->SetMaxDCAToVertexZ(3);
4106                 //_________________________________________________
4107
4108                 ///#################################################################
4109                 //Non-HFE reconstruction
4110         fNonHFE = new AliSelectNonHFE();
4111         fNonHFE->SetAODanalysis(fIsAOD);
4112         if(fMassCutFlag) fNonHFE->SetInvariantMassCut(fMassCut);
4113         if(fAngleCutFlag) fNonHFE->SetOpeningAngleCut(fAngleCut);
4114         if(fChi2CutFlag) fNonHFE->SetChi2OverNDFCut(fChi2Cut);
4115         if(fDCAcutFlag) fNonHFE->SetDCACut(fDCAcut);
4116         fNonHFE->SetAlgorithm("DCA"); //KF
4117         fNonHFE->SetPIDresponse(fPidResponse);
4118         fNonHFE->SetTrackCuts(-3.5,3.5,fPartnerCuts);
4119         fNonHFE->SetAdditionalCuts(fPtMinAsso,fTpcNclsAsso);
4120
4121
4122         fNonHFE->SetHistAngleBack(fOpAngleBack);
4123         fNonHFE->SetHistAngle(fOpAngle);
4124         fNonHFE->SetHistDCABack(fDCABack);
4125         fNonHFE->SetHistDCA(fDCA);
4126         fNonHFE->SetHistMassBack(fInvMassBack);
4127         fNonHFE->SetHistMass(fInvMass);
4128
4129         fNonHFE->FindNonHFE(trackIndex,vtrack,fVevent);
4130
4131         Int_t *fUlsPartner = fNonHFE->GetPartnersULS();
4132         Int_t *fLsPartner = fNonHFE->GetPartnersLS();
4133         Bool_t fUlsIsPartner = kFALSE;
4134         Bool_t fLsIsPartner = kFALSE;
4135                 ///#################################################################
4136
4137
4138                 //Electron Information
4139         Double_t fPhiE = -999;
4140         Double_t fEtaE = -999;
4141         Double_t fPhiH = -999;
4142         Double_t fEtaH = -999;
4143         Double_t fDphi = -999;
4144         Double_t fDeta = -999;
4145         Double_t fPtE = -999;
4146         Double_t fPtH = -999;
4147
4148         Double_t pi = TMath::Pi();
4149
4150         fPhiE = track->Phi();
4151         fEtaE = track->Eta();
4152         fPtE = track->Pt();
4153
4154
4155                 ///_________________________________________________________________
4156                 ///MC analysis
4157         if(fIsMC)
4158         {
4159                 if(fIsAOD)
4160                 {
4161                         if(TMath::Abs(fMCparticle->GetPdgCode())==11 && (TMath::Abs(fMCparticleMother->GetPdgCode())==22 || TMath::Abs(fMCparticleMother->GetPdgCode())==111 || TMath::Abs(fMCparticleMother->GetPdgCode())==221))
4162                 {
4163                                 if(fNonHFE->IsULS()) fPtElec_ULS->Fill(fPtE,fNonHFE->GetNULS());
4164                                 if(fNonHFE->IsLS()) fPtElec_LS->Fill(fPtE,fNonHFE->GetNLS());
4165
4166
4167                         }
4168
4169                         if(fNonHFE->IsULS()) fPtElec_ULS_NoPid->Fill(fPtE,fNonHFE->GetNULS());
4170                         if(fNonHFE->IsLS()) fPtElec_LS_NoPid->Fill(fPtE,fNonHFE->GetNLS());
4171                 }
4172                 else
4173                 {
4174                         if(TMath::Abs(fMCtrack->GetPdgCode())==11 && (TMath::Abs(fMCtrackMother->GetPdgCode())==22 || TMath::Abs(fMCtrackMother->GetPdgCode())==111 || TMath::Abs(fMCtrackMother->GetPdgCode())==221))
4175                         {
4176                                 if(fNonHFE->IsULS()) fPtElec_ULS->Fill(fPtE,fNonHFE->GetNULS());
4177                                 if(fNonHFE->IsLS()) fPtElec_LS->Fill(fPtE,fNonHFE->GetNLS());
4178
4179                         }
4180
4181
4182                 }
4183         }
4184                 ///_________________________________________________________________
4185         else
4186         {
4187                 if(fNonHFE->IsULS()) fPtElec_ULS->Fill(fPtE,fNonHFE->GetNULS());
4188                 if(fNonHFE->IsLS()) fPtElec_LS->Fill(fPtE,fNonHFE->GetNLS());
4189         }
4190
4191
4192
4193
4194                 //__________________________________________________________________
4195                 //Event Mixing Analysis - Hadron Loop
4196                 //Retrieve
4197         if(fEventMixingFlag)
4198         {
4199                 fPool = fPoolMgr->GetEventPool(fCentrality->GetCentralityPercentile("V0A"), fZvtx); // Get the buffer associated with the current centrality and z-vtx
4200
4201                 if(!fPool) AliFatal(Form("No pool found for centrality = %f, zVtx = %f",fCentrality->GetCentralityPercentile("V0A"), fZvtx));
4202
4203                 if(fPool->GetCurrentNEvents() >= 5) // start mixing when 5 events are in the buffer
4204                 {
4205                         fPoolNevents->Fill(fPool->GetCurrentNEvents());
4206
4207                         for (Int_t jMix = 0; jMix < fPool->GetCurrentNEvents(); jMix++)  // mix with each event in the buffer
4208                         {
4209                                 TObjArray* bgTracks = fPool->GetEvent(jMix);
4210
4211                                 for (Int_t kMix = 0; kMix < bgTracks->GetEntriesFast(); kMix++)  // mix with each track in the event
4212                                 {
4213                                         const AliEHCParticle* MixedTrack(dynamic_cast<AliEHCParticle*>(bgTracks->At(kMix)));
4214                                         if (NULL == MixedTrack) continue;
4215
4216                                         fPhiH = MixedTrack->Phi();
4217                                         fEtaH = MixedTrack->Eta();
4218                                         fPtH = MixedTrack->Pt();
4219
4220                                         if(fPtH<fAssHadronPtMin || fPtH>fAssHadronPtMax) continue;
4221
4222                                         fDphi = fPhiE - fPhiH;
4223
4224                                         if (fDphi > 3*pi/2) fDphi = fDphi - 2*pi;
4225                                         if (fDphi < -pi/2)  fDphi = fDphi + 2*pi;
4226
4227                                         fDeta = fEtaE - fEtaH;
4228
4229                                         Double_t fPtBin[7] = {1,2,4,6,8,10,15};
4230
4231                                         for(Int_t i = 0; i < 6; i++)
4232                                         {
4233                                             if(fPtE>=fPtBin[i] && fPtE<fPtBin[i+1])
4234                                             {
4235                                                         fCEtaPhi_Inc_EM[i]->Fill(fDphi,fDeta);
4236
4237                                                         if(fNonHFE->IsULS()) fCEtaPhi_ULS_EM[i]->Fill(fDphi,fDeta);
4238                                                         if(fNonHFE->IsLS()) fCEtaPhi_LS_EM[i]->Fill(fDphi,fDeta);
4239
4240                                                         if(fNonHFE->IsULS()) fCEtaPhi_ULS_Weight_EM[i]->Fill(fDphi,fDeta,fNonHFE->GetNULS());
4241                                                         if(fNonHFE->IsLS()) fCEtaPhi_LS_Weight_EM[i]->Fill(fDphi,fDeta,fNonHFE->GetNLS());
4242                                             }
4243                                         }
4244
4245                                                 // TODO your code: do event mixing with current event and bgTracks
4246                                                 // note that usually the content filled now is weighted by 1 / pool->GetCurrentNEvents()
4247                                 }
4248                         }
4249                 }
4250         }
4251                 //__________________________________________________________________
4252
4253                 //__________________________________________________________________
4254                 //Same Event Analysis - Hadron Loop
4255         for(Int_t iTracks = 0; iTracks < fVevent->GetNumberOfTracks(); iTracks++)
4256         {
4257                 if(trackIndex==iTracks) continue;
4258
4259                 AliVParticle* Vtrack2 = fVevent->GetTrack(iTracks);
4260                 if (!Vtrack2)
4261                 {
4262                         printf("ERROR: Could not receive track %d\n", iTracks);
4263                         continue;
4264                 }
4265
4266                 AliVTrack *track2 = dynamic_cast<AliVTrack*>(Vtrack2);
4267
4268                 if(track2->Eta()<fEtaCutMin || track2->Eta()>fEtaCutMax ) continue;
4269
4270                 if(fIsAOD)
4271                 {
4272                         AliAODTrack *atrack2 = dynamic_cast<AliAODTrack*>(Vtrack2);
4273                         if(!atrack2->TestFilterMask(AliAODTrack::kTrkTPCOnly)) continue;
4274                         if((!(atrack2->GetStatus()&AliESDtrack::kITSrefit)|| (!(atrack2->GetStatus()&AliESDtrack::kTPCrefit)))) continue;
4275                         if(atrack2->GetTPCNcls() < 80) continue;
4276                         if(fAssocWithSPD && ((!(atrack2->HasPointOnITSLayer(0))) && (!(atrack2->HasPointOnITSLayer(1))))) continue;
4277                 }
4278                 else
4279                 {
4280                         AliESDtrack *etrack2 = dynamic_cast<AliESDtrack*>(Vtrack2);
4281                         if(!fPartnerCuts->AcceptTrack(etrack2)) continue;
4282                 }
4283
4284                 fPhiH = track2->Phi();
4285                 fEtaH = track2->Eta();
4286                 fPtH = track2->Pt();
4287
4288                 if(fPtH<fAssHadronPtMin || fPtH>fAssHadronPtMax) continue;
4289
4290                 fDphi = fPhiE - fPhiH;
4291
4292                 if (fDphi > 3*pi/2) fDphi = fDphi - 2*pi;
4293                 if (fDphi < -pi/2)  fDphi = fDphi + 2*pi;
4294
4295                 fDeta = fEtaE - fEtaH;
4296
4297                 Double_t fPtBin[7] = {1,2,4,6,8,10,15};
4298
4299                 //______________________________________________________________
4300                 //Check if this track is a Non-HFE partner
4301                 fUlsIsPartner = kFALSE;
4302                 fLsIsPartner = kFALSE;
4303                 for(Int_t i = 0; i < fNonHFE->GetNULS(); i++)
4304                 {
4305                         if(fUlsPartner[i]==iTracks) fUlsIsPartner=kTRUE;
4306                 }
4307                 for(Int_t i = 0; i < fNonHFE->GetNLS(); i++)
4308                 {
4309                         if(fLsPartner[i]==iTracks) fLsIsPartner=kTRUE;
4310                 }
4311                 //______________________________________________________________
4312
4313                 for(Int_t i = 0; i < 6; i++)
4314                 {
4315                     if(fPtE>=fPtBin[i] && fPtE<fPtBin[i+1])
4316                     {
4317                                 fCEtaPhi_Inc[i]->Fill(fDphi,fDeta);
4318
4319                                 if(fNonHFE->IsULS()) fCEtaPhi_ULS[i]->Fill(fDphi,fDeta);
4320                                 if(fNonHFE->IsLS()) fCEtaPhi_LS[i]->Fill(fDphi,fDeta);
4321                                 if(fNonHFE->IsULS() && !fUlsIsPartner) fCEtaPhi_ULS_NoP[i]->Fill(fDphi,fDeta);
4322                                 if(fNonHFE->IsLS() && !fLsIsPartner) fCEtaPhi_LS_NoP[i]->Fill(fDphi,fDeta);
4323
4324                                 if(fNonHFE->IsULS()) fCEtaPhi_ULS_Weight[i]->Fill(fDphi,fDeta,fNonHFE->GetNULS());
4325                                 if(fNonHFE->IsLS()) fCEtaPhi_LS_Weight[i]->Fill(fDphi,fDeta,fNonHFE->GetNLS());
4326                                 if(fNonHFE->IsULS() && !fUlsIsPartner) fCEtaPhi_ULS_NoP_Weight[i]->Fill(fDphi,fDeta,fNonHFE->GetNULS());
4327                                 if(fNonHFE->IsLS() && !fLsIsPartner) fCEtaPhi_LS_NoP_Weight[i]->Fill(fDphi,fDeta,fNonHFE->GetNLS());
4328                     }
4329                 }
4330         }
4331 }
4332
4333         //____________________________________________________________________________________________________________
4334         //Create a TObjArray with selected hadrons, for the mixed event analysis
4335 TObjArray* AliAnalysisTaskEMCalHFEpA::SelectedHadrons()
4336 {
4337         fTracksClone = new TObjArray;
4338         fTracksClone->SetOwner(kTRUE);
4339
4340         for(Int_t iTracks = 0; iTracks < fVevent->GetNumberOfTracks(); iTracks++)
4341         {
4342                 AliVParticle* Vtrack2 = fVevent->GetTrack(iTracks);
4343                 if (!Vtrack2)
4344                 {
4345                         printf("ERROR: Could not receive track %d\n", iTracks);
4346                         continue;
4347                 }
4348
4349                 AliVTrack *track2 = dynamic_cast<AliVTrack*>(Vtrack2);
4350
4351                 if(track2->Eta()<fEtaCutMin || track2->Eta()>fEtaCutMax ) continue;
4352
4353                 if(fIsAOD)
4354                 {
4355                         AliAODTrack *atrack2 = dynamic_cast<AliAODTrack*>(Vtrack2);
4356                         if(!atrack2->TestFilterMask(AliAODTrack::kTrkTPCOnly)) continue;
4357                         if((!(atrack2->GetStatus()&AliESDtrack::kITSrefit)|| (!(atrack2->GetStatus()&AliESDtrack::kTPCrefit)))) continue;
4358                         if(atrack2->GetTPCNcls() < 80) continue;
4359                         if(fAssocWithSPD && ((!(atrack2->HasPointOnITSLayer(0))) && (!(atrack2->HasPointOnITSLayer(1))))) continue;
4360                 }
4361                 else
4362                 {
4363                         AliESDtrack *etrack2 = dynamic_cast<AliESDtrack*>(Vtrack2);
4364                         if(!fPartnerCuts->AcceptTrack(etrack2)) continue;
4365                 }
4366
4367                 fTracksClone->Add(new AliEHCParticle(track2->Eta(), track2->Phi(), track2->Pt()));
4368         }
4369         return fTracksClone;
4370 }
4371         //____________________________________________________________________________________________________________
4372
4373         //______________________________________________________________________
4374 void AliAnalysisTaskEMCalHFEpA::DiHadronCorrelation(AliVTrack *track, Int_t trackIndex)
4375 {
4376                 //________________________________________________
4377                 //Associated particle cut
4378         fPartnerCuts->SetAcceptKinkDaughters(kFALSE);
4379     fPartnerCuts->SetRequireITSRefit(kTRUE);
4380     fPartnerCuts->SetRequireTPCRefit(kTRUE);
4381     fPartnerCuts->SetEtaRange(-0.9,0.9);
4382     fPartnerCuts->SetMaxChi2PerClusterTPC(4.0);
4383     fPartnerCuts->SetMinNClustersTPC(80);
4384     fPartnerCuts->SetPtRange(0.3,1e10);
4385                 //fPartnerCuts->SetRequireSigmaToVertex(kTRUE);
4386                 //fPartnerCuts->SetMaxDCAToVertexXY(1);
4387                 //fPartnerCuts->SetMaxDCAToVertexZ(3);
4388                 //_________________________________________________
4389
4390                 //Electron Information
4391         Double_t fPhiE = -999;
4392         Double_t fEtaE = -999;
4393         Double_t fPhiH = -999;
4394         Double_t fEtaH = -999;
4395         Double_t fDphi = -999;
4396         Double_t fDeta = -999;
4397         Double_t fPtE = -999;
4398         Double_t fPtH = -999;
4399
4400         Double_t pi = TMath::Pi();
4401
4402         fPhiE = track->Phi();
4403         fEtaE = track->Eta();
4404         fPtE = track->Pt();
4405
4406                 //__________________________________________________________________
4407                 //Same Event Analysis - Hadron Loop
4408         for(Int_t iTracks = 0; iTracks < fVevent->GetNumberOfTracks(); iTracks++)
4409         {
4410                 if(trackIndex==iTracks) continue;
4411
4412                 AliVParticle* Vtrack2 = fVevent->GetTrack(iTracks);
4413                 if (!Vtrack2)
4414                 {
4415                         printf("ERROR: Could not receive track %d\n", iTracks);
4416                         continue;
4417                 }
4418
4419                 AliVTrack *track2 = dynamic_cast<AliVTrack*>(Vtrack2);
4420                 if(track2->Eta()<fEtaCutMin || track2->Eta()>fEtaCutMax ) continue;
4421
4422                 if(fIsAOD)
4423                 {
4424                         AliAODTrack *atrack2 = dynamic_cast<AliAODTrack*>(Vtrack2);
4425                         if(!atrack2->TestFilterMask(AliAODTrack::kTrkTPCOnly)) continue;
4426                         if((!(atrack2->GetStatus()&AliESDtrack::kITSrefit)|| (!(atrack2->GetStatus()&AliESDtrack::kTPCrefit)))) continue;
4427                         if(atrack2->GetTPCNcls() < 80) continue;
4428                         if(fAssocWithSPD && ((!(atrack2->HasPointOnITSLayer(0))) && (!(atrack2->HasPointOnITSLayer(1))))) continue;
4429                 }
4430                 else
4431                 {
4432                         AliESDtrack *etrack2 = dynamic_cast<AliESDtrack*>(Vtrack2);
4433                         if(!fPartnerCuts->AcceptTrack(etrack2)) continue;
4434                 }
4435
4436                 fPhiH = track2->Phi();
4437                 fEtaH = track2->Eta();
4438                 fPtH = track2->Pt();
4439
4440                 if(fPtH<fAssHadronPtMin || fPtH>fAssHadronPtMax) continue;
4441
4442                 fDphi = fPhiE - fPhiH;
4443
4444                 if (fDphi > 3*pi/2) fDphi = fDphi - 2*pi;
4445                 if (fDphi < -pi/2)  fDphi = fDphi + 2*pi;
4446
4447                 fDeta = fEtaE - fEtaH;
4448
4449                 Double_t fPtBin[7] = {1,2,4,6,8,10,15};
4450
4451                 for(Int_t i = 0; i < 6; i++)
4452                 {
4453                     if(fPtE>=fPtBin[i] && fPtE<fPtBin[i+1])
4454                     {
4455                                 fCEtaPhi_Inc_DiHadron[i]->Fill(fDphi,fDeta);
4456                     }
4457                 }
4458         }
4459 }
4460         //____________________________________________________________________________________________________________
4461
4462         //______________________________________________________________________
4463 Bool_t AliAnalysisTaskEMCalHFEpA::FindMother(Int_t mcIndex)
4464 {
4465         fIsHFE1 = kFALSE;
4466         fIsHFE2 = kFALSE;
4467         fIsNonHFE = kFALSE;
4468         fIsFromD = kFALSE;
4469         fIsFromB = kFALSE;
4470         fIsFromPi0 = kFALSE;
4471         fIsFromEta = kFALSE;
4472         fIsFromGamma = kFALSE;
4473
4474         if(mcIndex < 0 || !fIsMC)
4475         {
4476                 return kFALSE;
4477         }
4478
4479         Int_t pdg = -99999;
4480         Int_t mpdg = -99999;
4481         Int_t gmpdg = -99999;
4482         Int_t ggmpdg = -99999;
4483         Int_t gggmpdg = -99999;
4484
4485         if(fIsAOD)
4486         {
4487                 fMCparticle = (AliAODMCParticle*) fMCarray->At(mcIndex);
4488
4489                 pdg = TMath::Abs(fMCparticle->GetPdgCode());
4490
4491
4492                 if(pdg!=11)
4493                 {
4494                         fIsHFE1 = kFALSE;
4495                         fIsHFE2 = kFALSE;
4496                         fIsNonHFE = kFALSE;
4497                         fIsFromD = kFALSE;
4498                         fIsFromB = kFALSE;
4499                         fIsFromPi0 = kFALSE;
4500                         fIsFromEta = kFALSE;
4501                         fIsFromGamma = kFALSE;
4502                         return kFALSE;
4503                 }
4504
4505                 if(fMCparticle->GetMother()<0)
4506                 {
4507                         fIsHFE1 = kFALSE;
4508                         fIsHFE2 = kFALSE;
4509                         fIsNonHFE = kFALSE;
4510                         fIsFromD = kFALSE;
4511                         fIsFromB = kFALSE;
4512                         fIsFromPi0 = kFALSE;
4513                         fIsFromEta = kFALSE;
4514                         fIsFromGamma = kFALSE;
4515                         return kFALSE;
4516                 }
4517
4518                 fMCparticleMother = (AliAODMCParticle*) fMCarray->At(fMCparticle->GetMother());
4519                 mpdg = TMath::Abs(fMCparticleMother->GetPdgCode());
4520
4521                 if(fMCparticleMother->GetMother()<0)
4522                 {
4523                         gmpdg = 0;
4524                         ggmpdg = 0;
4525                         gggmpdg = 0;
4526                 }
4527                 else
4528                 {
4529                         fMCparticleGMother = (AliAODMCParticle*) fMCarray->At(fMCparticleMother->GetMother());
4530                         gmpdg = TMath::Abs(fMCparticleGMother->GetPdgCode());
4531                         if(fMCparticleGMother->GetMother()<0)
4532                         {
4533                                 ggmpdg = 0;
4534                                 gggmpdg = 0;
4535                         }
4536                         else
4537                         {
4538                                 fMCparticleGGMother = (AliAODMCParticle*) fMCarray->At(fMCparticleGMother->GetMother());
4539                                 ggmpdg = TMath::Abs(fMCparticleGGMother->GetPdgCode());
4540                                 if(fMCparticleGGMother->GetMother()<0)
4541                                 {
4542                                         gggmpdg = 0;
4543                                 }
4544                                 else
4545                                 {
4546                                         fMCparticleGGGMother = (AliAODMCParticle*) fMCarray->At(fMCparticleGGMother->GetMother());
4547                                         gggmpdg = TMath::Abs(fMCparticleGGGMother->GetPdgCode());
4548                                 }
4549                         }
4550                 }
4551         }
4552         else
4553         {
4554                 fMCtrack = fMCstack->Particle(mcIndex);
4555
4556                 pdg = TMath::Abs(fMCtrack->GetPdgCode());
4557
4558                 if(pdg!=11)
4559                 {
4560                         fIsHFE1 = kFALSE;
4561                         fIsHFE2 = kFALSE;
4562                         fIsNonHFE = kFALSE;
4563                         fIsFromD = kFALSE;
4564                         fIsFromB = kFALSE;
4565                         fIsFromPi0 = kFALSE;
4566                         fIsFromEta = kFALSE;
4567                         fIsFromGamma = kFALSE;
4568                         return kFALSE;
4569                 }
4570
4571                 if(fMCtrack->GetFirstMother()<0)
4572                 {
4573                         fIsHFE1 = kFALSE;
4574                         fIsHFE2 = kFALSE;
4575                         fIsNonHFE = kFALSE;
4576                         fIsFromD = kFALSE;
4577                         fIsFromB = kFALSE;
4578                         fIsFromPi0 = kFALSE;
4579                         fIsFromEta = kFALSE;
4580                         fIsFromGamma = kFALSE;
4581                         return kFALSE;
4582                 }
4583
4584                 fMCtrackMother = fMCstack->Particle(fMCtrack->GetFirstMother());
4585                 mpdg = TMath::Abs(fMCtrackMother->GetPdgCode());
4586
4587                 if(fMCtrackMother->GetFirstMother()<0)
4588                 {
4589                         gmpdg = 0;
4590                         ggmpdg = 0;
4591                         gggmpdg = 0;
4592                 }
4593                 else
4594                 {
4595                         fMCtrackGMother = fMCstack->Particle(fMCtrackMother->GetFirstMother());
4596                         gmpdg = TMath::Abs(fMCtrackGMother->GetPdgCode());
4597
4598                         if(fMCtrackGMother->GetFirstMother()<0)
4599                         {
4600                                 ggmpdg = 0;
4601                                 gggmpdg = 0;
4602                         }
4603                         else
4604                         {
4605                                 fMCtrackGGMother = fMCstack->Particle(fMCtrackGMother->GetFirstMother());
4606                                 ggmpdg = TMath::Abs(fMCtrackGGMother->GetPdgCode());
4607
4608                                 if(fMCtrackGGMother->GetFirstMother()<0)
4609                                 {
4610                                         gggmpdg = 0;
4611                                 }
4612                                 else
4613                                 {
4614                                         fMCtrackGGGMother = fMCstack->Particle(fMCtrackGGMother->GetFirstMother());
4615                                         gggmpdg = TMath::Abs(fMCtrackGGGMother->GetPdgCode());
4616                                 }
4617                         }
4618                 }
4619         }
4620
4621         //Tag Electron Source
4622         if(mpdg==111 || mpdg==221 || mpdg==22)
4623         {
4624                 fIsHFE1 = kFALSE;
4625                 fIsHFE2 = kFALSE;
4626                 fIsNonHFE = kTRUE;
4627                 fIsFromD = kFALSE;
4628                 fIsFromB = kFALSE;
4629
4630                 fIsFromPi0 = kFALSE;
4631                 fIsFromEta = kFALSE;
4632                 fIsFromGamma = kFALSE;
4633
4634                 if(mpdg==111) fIsFromPi0 = kFALSE;
4635                 if(mpdg==221)fIsFromEta = kFALSE;
4636                 if(mpdg==22) fIsFromGamma = kFALSE;
4637
4638                 return kTRUE;
4639         }
4640         else
4641         {
4642                 fIsHFE1 = kFALSE;
4643                 fIsHFE2 = kTRUE;
4644
4645                 fIsFromPi0 = kFALSE;
4646                 fIsFromEta = kFALSE;
4647                 fIsFromGamma = kFALSE;
4648
4649                 fIsNonHFE = kFALSE;
4650
4651                 fIsFromD = kFALSE;
4652                 fIsFromB = kFALSE;
4653
4654                 if(mpdg>400 && mpdg<500)
4655                 {
4656                         if((gmpdg>500 && gmpdg<600) || (ggmpdg>500 && ggmpdg<600) || (gggmpdg>500 && gggmpdg<600))
4657                         {
4658                                 fIsHFE1 = kTRUE;
4659                                 fIsFromD = kFALSE;
4660                                 fIsFromB = kTRUE;
4661                                 return kTRUE;
4662                         }
4663                         else
4664                         {
4665                                 fIsHFE1 = kTRUE;
4666                                 fIsFromD = kTRUE;
4667                                 fIsFromB = kFALSE;
4668                                 return kTRUE;
4669                         }
4670                 }
4671                 else if(mpdg>500 && mpdg<600)
4672                 {
4673                         fIsHFE1 = kTRUE;
4674                         fIsFromD = kFALSE;
4675                         fIsFromB = kTRUE;
4676                         return kTRUE;
4677                 }
4678                 else
4679                 {
4680                         fIsHFE1 = kFALSE;
4681                         fIsFromD = kFALSE;
4682                         fIsFromB = kFALSE;
4683                         return kFALSE;
4684                 }
4685         }
4686 }
4687 /*
4688 Bool_t AliAnalysisTaskEMCalHFEpA::ContainsBadChannel(TString calorimeter,UShort_t* cellList, Int_t nCells)
4689 {
4690                 // Check that in the cluster cells, there is no bad channel of those stored
4691                 // in fEMCALBadChannelMap
4692                 // adapted from AliCalorimeterUtils
4693
4694                 //if (!fRemoveBadChannels) return kFALSE;
4695                 //printf("fEMCALBadChannelMap %p, fPHOSBadChannelMap %p \n",fEMCALBadChannelMap,fPHOSBadChannelMap);
4696         if( !fEMCALRecoUtils->GetEMCALChannelStatusMap(0)) return kFALSE;
4697
4698                 //Int_t icol = -1;
4699                 //Int_t irow = -1;
4700                 //Int_t imod = -1;
4701         for(Int_t iCell = 0; iCell<nCells; iCell++){
4702
4703                         //Get the column and row
4704                 if(calorimeter == "EMCAL"){
4705                         return fEMCALRecoUtils->ClusterContainsBadChannel((AliEMCALGeometry*)fEMCALGeo,cellList,nCells);
4706                 }
4707                 else return kFALSE;
4708
4709         }// cell cluster loop
4710
4711         return kFALSE;
4712
4713 }
4714 */
4715 /*
4716
4717 //________________________________________________________________________________
4718 TArrayI AliAnalysisTaskEMCalHFEpA::GetTriggerPatches(Bool_t IsEventEMCALL0, Bool_t IsEventEMCALL1)
4719 {
4720         // Select the patches that triggered
4721         // Depend on L0 or L1
4722
4723         // some variables
4724         //Int_t  trigtimes[30], globCol, globRow,ntimes, i;
4725         Int_t   globCol, globRow;
4726
4727         Int_t  absId  = -1; //[100];
4728         Int_t  nPatch = 0;
4729
4730         TArrayI patches(0);
4731
4732                 // get object pointer
4733         AliVCaloTrigger *caloTrigger = InputEvent()->GetCaloTrigger( "EMCAL" );
4734
4735                 // class is not empty
4736         if( caloTrigger->GetEntries() > 0 )
4737         {
4738                         // must reset before usage, or the class will fail
4739                 caloTrigger->Reset();
4740
4741                         // go throuth the trigger channels
4742                 while( caloTrigger->Next() )
4743                 {
4744                                 // get position in global 2x2 tower coordinates
4745                         caloTrigger->GetPosition( globCol, globRow );
4746
4747                         //L0
4748                         if(IsEventEMCALL0)
4749                         {
4750                             //not implemented
4751                         }
4752
4753
4754                         else if(IsEventEMCALL1) // L1
4755                         {
4756                                         Int_t bit = 0;
4757                                         caloTrigger->GetTriggerBits(bit);
4758
4759                                         Bool_t isEGA = ((bit >> fBitEGA) & 0x1);
4760                                         //Bool_t isEJE = ((bit >> fBitEJE) & 0x1) && IsEventEMCALL1Jet  () ;
4761
4762                                         if(!isEGA) continue;
4763
4764                                         Int_t patchsize = -1;
4765                                         if(isEGA) patchsize =  2;
4766                                         //else if (isEJE) patchsize = 16;
4767
4768                                         // add 2x2 (EGA) or 16x16 (EJE) patches
4769                                         for(Int_t irow=0; irow < patchsize; irow++)
4770                                         {
4771                                                 for(Int_t icol=0; icol < patchsize; icol++)
4772                                                 {
4773                                                         GetCaloUtils()->GetEMCALGeometry()->GetAbsFastORIndexFromPositionInEMCAL(globCol+icol,globRow+irow, absId);
4774                                                         //printf("pass the time cut globCol %d, globRow %d absId %d\n",globCol,globRow, absIDTrig[nTrig]);
4775                                                         patches.Set(nPatch+1); // Set size of this array to nPatch+1 ints.
4776                                                         patches.AddAt(absId,nPatch++); //Add Int_t absId at position nPatch++
4777                                                 }
4778                                         }
4779
4780                         } // L1
4781
4782                 } // trigger iterator
4783         } // go thorough triggers
4784
4785         printf("N patches %d, test %d,first %d, last %d\n",patches.GetSize(), nPatch, patches.At(0), patches.At(patches.GetSize()-1));
4786
4787         return patches;
4788 }
4789  */
4790 Double_t AliAnalysisTaskEMCalHFEpA::CalculateWeight(Int_t pdg_particle, Double_t x)
4791 {
4792                 //weight for d3 based on MinJung parametrization //sent by Jan
4793             Double_t weight=1;
4794                 if(pdg_particle==111){
4795                         if(x>= 0.100000 &&  x < 0.112797 ) weight=1.262120;
4796                         if(x>= 0.112797 &&  x < 0.127231 ) weight=1.277765;
4797                         if(x>= 0.127231 &&  x < 0.143512 ) weight=1.295605;
4798                         if(x>= 0.143512 &&  x < 0.161877 ) weight=1.318155;
4799                         if(x>= 0.161877 &&  x < 0.182592 ) weight=1.348693;
4800                         if(x>= 0.182592 &&  x < 0.205957 ) weight=1.388636;
4801                         if(x>= 0.205957 &&  x < 0.232313 ) weight=1.439122;
4802                         if(x>= 0.232313 &&  x < 0.262041 ) weight=1.497452;
4803                         if(x>= 0.262041 &&  x < 0.295573 ) weight=1.559409;
4804                         if(x>= 0.295573 &&  x < 0.333397 ) weight=1.615169;
4805                         if(x>= 0.333397 &&  x < 0.376060 ) weight=1.654954;
4806                         if(x>= 0.376060 &&  x < 0.424183 ) weight=1.668753;
4807                         if(x>= 0.424183 &&  x < 0.478465 ) weight=1.652225;
4808                         if(x>= 0.478465 &&  x < 0.539692 ) weight=1.603119;
4809                         if(x>= 0.539692 &&  x < 0.608754 ) weight=1.526049;
4810                         if(x>= 0.608754 &&  x < 0.686654 ) weight=1.426724;
4811                         if(x>= 0.686654 &&  x < 0.774523 ) weight=1.312684;
4812                         if(x>= 0.774523 &&  x < 0.873636 ) weight=1.195395;
4813                         if(x>= 0.873636 &&  x < 0.985432 ) weight=1.086264;
4814                         if(x>= 0.985432 &&  x < 1.111534 ) weight=0.993666;
4815                         if(x>= 1.111534 &&  x < 1.253773 ) weight=0.922587;
4816                         if(x>= 1.253773 &&  x < 1.414214 ) weight=0.875739;
4817                         if(x>= 1.414214 &&  x < 1.595185 ) weight=0.852181;
4818                         if(x>= 1.595185 &&  x < 1.799315 ) weight=0.847828;
4819                         if(x>= 1.799315 &&  x < 2.029567 ) weight=0.863875;
4820                         if(x>= 2.029567 &&  x < 2.289283 ) weight=0.899112;
4821                         if(x>= 2.289283 &&  x < 2.582235 ) weight=0.955194;
4822                         if(x>= 2.582235 &&  x < 2.912674 ) weight=1.033824;
4823                         if(x>= 2.912674 &&  x < 3.285398 ) weight=1.133714;
4824                         if(x>= 3.285398 &&  x < 3.705818 ) weight=1.259471;
4825                         if(x>= 3.705818 &&  x < 4.180038 ) weight=1.406883;
4826                         if(x>= 4.180038 &&  x < 4.714942 ) weight=1.578923;
4827                         if(x>= 4.714942 &&  x < 5.318296 ) weight=1.778513;
4828                         if(x>= 5.318296 &&  x < 5.998859 ) weight=2.001171;
4829                         if(x>= 5.998859 &&  x < 6.766511 ) weight=2.223161;
4830                         if(x>= 6.766511 &&  x < 7.632396 ) weight=2.449445;
4831                         if(x>= 7.632396 &&  x < 8.609086 ) weight=2.661734;
4832                         if(x>= 8.609086 &&  x < 9.710759 ) weight=2.851935;
4833                         if(x>= 9.710759 &&  x < 10.953409 ) weight=2.974319;
4834                         if(x>= 10.953409 &&  x < 12.355077 ) weight=3.106314;
4835                         if(x>= 12.355077 &&  x < 13.936111 ) weight=3.126815;
4836                         if(x>= 13.936111 &&  x < 15.719464 ) weight=3.150053;
4837                         if(x>= 15.719464 &&  x < 17.731026 ) weight=3.218509;
4838                         if(x>= 17.731026 &&  x < 20.000000 ) weight=3.252141;
4839
4840                 }
4841             else if(pdg_particle==221){
4842                         if(x>= 0.100000 &&  x < 0.112797 ) weight=2.159358;
4843                         if(x>= 0.112797 &&  x < 0.127231 ) weight=2.145546;
4844                         if(x>= 0.127231 &&  x < 0.143512 ) weight=2.132390;
4845                         if(x>= 0.143512 &&  x < 0.161877 ) weight=2.109918;
4846                         if(x>= 0.161877 &&  x < 0.182592 ) weight=2.084920;
4847                         if(x>= 0.182592 &&  x < 0.205957 ) weight=2.054302;
4848                         if(x>= 0.205957 &&  x < 0.232313 ) weight=2.015202;
4849                         if(x>= 0.232313 &&  x < 0.262041 ) weight=1.966068;
4850                         if(x>= 0.262041 &&  x < 0.295573 ) weight=1.912255;
4851                         if(x>= 0.295573 &&  x < 0.333397 ) weight=1.844087;
4852                         if(x>= 0.333397 &&  x < 0.376060 ) weight=1.767913;
4853                         if(x>= 0.376060 &&  x < 0.424183 ) weight=1.680366;
4854                         if(x>= 0.424183 &&  x < 0.478465 ) weight=1.583468;
4855                         if(x>= 0.478465 &&  x < 0.539692 ) weight=1.475131;
4856                         if(x>= 0.539692 &&  x < 0.608754 ) weight=1.361436;
4857                         if(x>= 0.608754 &&  x < 0.686654 ) weight=1.244388;
4858                         if(x>= 0.686654 &&  x < 0.774523 ) weight=1.125973;
4859                         if(x>= 0.774523 &&  x < 0.873636 ) weight=1.015769;
4860                         if(x>= 0.873636 &&  x < 0.985432 ) weight=0.914579;
4861                         if(x>= 0.985432 &&  x < 1.111534 ) weight=0.830529;
4862                         if(x>= 1.111534 &&  x < 1.253773 ) weight=0.766397;
4863                         if(x>= 1.253773 &&  x < 1.414214 ) weight=0.723663;
4864                         if(x>= 1.414214 &&  x < 1.595185 ) weight=0.701236;
4865                         if(x>= 1.595185 &&  x < 1.799315 ) weight=0.695605;
4866                         if(x>= 1.799315 &&  x < 2.029567 ) weight=0.707578;
4867                         if(x>= 2.029567 &&  x < 2.289283 ) weight=0.735194;
4868                         if(x>= 2.289283 &&  x < 2.582235 ) weight=0.781052;
4869                         if(x>= 2.582235 &&  x < 2.912674 ) weight=0.842350;
4870                         if(x>= 2.912674 &&  x < 3.285398 ) weight=0.923676;
4871                         if(x>= 3.285398 &&  x < 3.705818 ) weight=1.028317;
4872                         if(x>= 3.705818 &&  x < 4.180038 ) weight=1.154029;
4873                         if(x>= 4.180038 &&  x < 4.714942 ) weight=1.296915;
4874                         if(x>= 4.714942 &&  x < 5.318296 ) weight=1.463674;
4875                         if(x>= 5.318296 &&  x < 5.998859 ) weight=1.651985;
4876                         if(x>= 5.998859 &&  x < 6.766511 ) weight=1.847912;
4877                         if(x>= 6.766511 &&  x < 7.632396 ) weight=2.066284;
4878                         if(x>= 7.632396 &&  x < 8.609086 ) weight=2.202231;
4879                         if(x>= 8.609086 &&  x < 9.710759 ) weight=2.399942;
4880                         if(x>= 9.710759 &&  x < 10.953409 ) weight=2.555106;
4881                         if(x>= 10.953409 &&  x < 12.355077 ) weight=2.632377;
4882                         if(x>= 12.355077 &&  x < 13.936111 ) weight=2.609660;
4883                         if(x>= 13.936111 &&  x < 15.719464 ) weight=2.656343;
4884                         if(x>= 15.719464 &&  x < 17.731026 ) weight=2.615438;
4885                         if(x>= 17.731026 &&  x < 20.000000 ) weight=2.576269;
4886
4887                 }
4888             else weight=1;
4889
4890             return weight;
4891
4892 }
4893 Double_t AliAnalysisTaskEMCalHFEpA::SetEoverPCutPtDependentMC(Double_t pt)
4894 {
4895
4896         fEoverPCutMin=0.8;
4897         fEoverPCutMax=1.2;
4898
4899         /*
4900         ////======================================================================================
4901         ////======================================================================================
4902         ////============================= data 2 sigma ===========================================
4903         if(!fIsMC){
4904                         //data 2 sigma
4905                 if(pt>= 2.00 &&  pt < 2.50){
4906                         fEoverPCutMin=0.7719;
4907                         fEoverPCutMax=1.1023;
4908                 }
4909                 if(pt>= 2.50 &&  pt < 3.00){
4910                         fEoverPCutMin=0.7966;
4911                         fEoverPCutMax=1.1088;
4912                 }
4913                 if(pt>= 3.00 &&  pt < 4.00){
4914                         fEoverPCutMin=0.8175;
4915                         fEoverPCutMax=1.1101;
4916                 }
4917                 if(pt>= 4.00 &&  pt < 5.00){
4918                         fEoverPCutMin=0.8302;
4919                         fEoverPCutMax=1.1194;
4920                 }
4921                 if(pt>= 5.00 &&  pt < 6.00){
4922                         fEoverPCutMin=0.8517;
4923                         fEoverPCutMax=1.1177;
4924                 }
4925                 if(pt>= 6.00 &&  pt < 8.00){
4926                         fEoverPCutMin=0.8901;
4927                         fEoverPCutMax=1.1139;
4928                 }
4929                 if(pt>= 8.00 &&  pt < 10.00){
4930                         fEoverPCutMin=0.8703;
4931                         fEoverPCutMax=1.1377;
4932                 }
4933                 if(pt>= 10.00 &&  pt < 12.00){
4934                         fEoverPCutMin=0.9043;
4935                         fEoverPCutMax=1.1977;
4936                 }
4937         }
4938
4939         ////=========================================== MC 2 sigma ===========================================
4940
4941
4942         if(fIsMC){
4943
4944                 if(pt>= 2.00 &&  pt < 2.50){
4945                         fEoverPCutMin=0.7712;
4946                         fEoverPCutMax=1.0746;
4947                 }
4948                 if(pt>= 2.50 &&  pt < 3.00){
4949                         fEoverPCutMin=0.7946;
4950                         fEoverPCutMax=1.0708;
4951                 }
4952                 if(pt>= 3.00 &&  pt < 4.00){
4953                         fEoverPCutMin=0.8196;
4954                         fEoverPCutMax=1.0678;
4955                 }
4956                 if(pt>= 4.00 &&  pt < 5.00){
4957                         fEoverPCutMin=0.8396;
4958                         fEoverPCutMax=1.0646;
4959                 }
4960                 if(pt>= 5.00 &&  pt < 6.00){
4961                         fEoverPCutMin=0.8527;
4962                         fEoverPCutMax=1.0647;
4963                 }
4964                 if(pt>= 6.00 &&  pt < 8.00){
4965                         fEoverPCutMin=0.8652;
4966                         fEoverPCutMax=1.0670;
4967                 }
4968                 if(pt>= 8.00 &&  pt < 10.00){
4969                         fEoverPCutMin=0.8748;
4970                         fEoverPCutMax=1.0804;
4971                 }
4972                 if(pt>= 10.00 &&  pt < 12.00){
4973                         fEoverPCutMin=0.8814;
4974                         fEoverPCutMax=1.0842;
4975                 }
4976
4977         }
4978         */
4979                 ////======================================================================================
4980                 ////======================================================================================
4981                 ////=========================== data 1.5 sigma ===========================================
4982         /*
4983         if(!fIsMC){
4984                 //data 1.5 sigmas
4985                 if(pt>= 2.00 &&  pt < 2.50){
4986                         fEoverPCutMin=0.8132;
4987                         fEoverPCutMax=1.0610;
4988                 }
4989                 if(pt>= 2.50 &&  pt < 3.00){
4990                         fEoverPCutMin=0.8356;
4991                         fEoverPCutMax=1.0698;
4992                 }
4993                 if(pt>= 3.00 &&  pt < 4.00){
4994                         fEoverPCutMin=0.8541;
4995                         fEoverPCutMax=1.0735;
4996                 }
4997                 if(pt>= 4.00 &&  pt < 5.00){
4998                         fEoverPCutMin=0.8664;
4999                         fEoverPCutMax=1.0832;
5000                 }
5001                 if(pt>= 5.00 &&  pt < 6.00){
5002                         fEoverPCutMin=0.8849;
5003                         fEoverPCutMax=1.0845;
5004                 }
5005                 if(pt>= 6.00 &&  pt < 8.00){
5006                         fEoverPCutMin=0.9180;
5007                         fEoverPCutMax=1.0860;
5008                 }
5009                 if(pt>= 8.00 &&  pt < 10.00){
5010                         fEoverPCutMin=0.9037;
5011                         fEoverPCutMax=1.1043;
5012                 }
5013                 if(pt>= 10.00 &&  pt < 12.00){
5014                         fEoverPCutMin=0.9409;
5015                         fEoverPCutMax=1.1611;
5016                 }
5017         }
5018
5019         ////======================================================================================
5020         ////======================================================================================
5021         ////=========================== MC 1.5 sigma ===========================================
5022
5023         if(fIsMC){
5024                 //MC 1.5 sigmas
5025                 if(pt>= 2.00 &&  pt < 2.50){
5026                         fEoverPCutMin=0.8091;
5027                         fEoverPCutMax=1.0367;
5028                 }
5029                 if(pt>= 2.50 &&  pt < 3.00){
5030                 fEoverPCutMin=0.8292;
5031                 fEoverPCutMax=1.0362;
5032                 }
5033                 if(pt>= 3.00 &&  pt < 4.00){
5034                         fEoverPCutMin=0.8506;
5035                         fEoverPCutMax=1.0368;
5036                 }
5037                 if(pt>= 4.00 &&  pt < 5.00){
5038                         fEoverPCutMin=0.8677;
5039                         fEoverPCutMax=1.0365;
5040                 }
5041                 if(pt>= 5.00 &&  pt < 6.00){
5042                         fEoverPCutMin=0.8792;
5043                         fEoverPCutMax=1.0382;
5044                 }
5045                 if(pt>= 6.00 &&  pt < 8.00){
5046                         fEoverPCutMin=0.8904;
5047                         fEoverPCutMax=1.0418;
5048                 }
5049                 if(pt>= 8.00 &&  pt < 10.00){
5050                         fEoverPCutMin=0.9005;
5051                         fEoverPCutMax=1.0547;
5052                 }
5053                 if(pt>= 10.00 &&  pt < 12.00){
5054                         fEoverPCutMin=0.9067;
5055                         fEoverPCutMax=1.0589;
5056                 }
5057
5058         }
5059         */
5060                 ////======================================================================================
5061                 ////======================================================================================
5062                 ////=========================== data 2.5 sigma ===========================================
5063
5064         if(!fIsMC){
5065                         //data 2.5 sigmas
5066                 if(pt>= 2.00 &&  pt < 2.50){
5067                         fEoverPCutMin=0.7306;
5068                         fEoverPCutMax=1.1436;
5069                 }
5070                 if(pt>= 2.50 &&  pt < 3.00){
5071                         fEoverPCutMin=0.7575;
5072                         fEoverPCutMax=1.1479;
5073                 }
5074                 if(pt>= 3.00 &&  pt < 4.00){
5075                         fEoverPCutMin=0.7809;
5076                         fEoverPCutMax=1.1467;
5077                 }
5078                 if(pt>= 4.00 &&  pt < 5.00){
5079                         fEoverPCutMin=0.7941;
5080                         fEoverPCutMax=1.1555;
5081                 }
5082                 if(pt>= 5.00 &&  pt < 6.00){
5083                         fEoverPCutMin=0.8184;
5084                         fEoverPCutMax=1.1510;
5085                 }
5086                 if(pt>= 6.00 &&  pt < 8.00){
5087                         fEoverPCutMin=0.8621;
5088                         fEoverPCutMax=1.1419;
5089                 }
5090                 if(pt>= 8.00 &&  pt < 10.00){
5091                         fEoverPCutMin=0.8368;
5092                         fEoverPCutMax=1.1712;
5093                 }
5094                 if(pt>= 10.00 &&  pt < 12.00){
5095                         fEoverPCutMin=0.8676;
5096                         fEoverPCutMax=1.2344;
5097                 }
5098
5099 }
5100
5101                 ////======================================================================================
5102                 ////======================================================================================
5103                 ////=========================== MC 2.5 sigma ===========================================
5104
5105         if(fIsMC){
5106                         //MC 2.5 sigmas
5107                 if(pt>= 2.00 &&  pt < 2.50){
5108                         fEoverPCutMin=0.7333;
5109                         fEoverPCutMax=1.1125;
5110                 }
5111                 if(pt>= 2.50 &&  pt < 3.00){
5112                         fEoverPCutMin=0.7601;
5113                         fEoverPCutMax=1.1053;
5114                 }
5115                 if(pt>= 3.00 &&  pt < 4.00){
5116                         fEoverPCutMin=0.7885;
5117                         fEoverPCutMax=1.0989;
5118                 }
5119                 if(pt>= 4.00 &&  pt < 5.00){
5120                         fEoverPCutMin=0.8115;
5121                         fEoverPCutMax=1.0927;
5122                 }
5123                 if(pt>= 5.00 &&  pt < 6.00){
5124                         fEoverPCutMin=0.8262;
5125                         fEoverPCutMax=1.0912;
5126                 }
5127                 if(pt>= 6.00 &&  pt < 8.00){
5128                         fEoverPCutMin=0.8400;
5129                         fEoverPCutMax=1.0922;
5130                 }
5131                 if(pt>= 8.00 &&  pt < 10.00){
5132                         fEoverPCutMin=0.8492;
5133                         fEoverPCutMax=1.1060;
5134                 }
5135                 if(pt>= 10.00 &&  pt < 12.00){
5136                         fEoverPCutMin=0.8560;
5137                         fEoverPCutMax=1.1096;
5138                 }
5139         }
5140
5141
5142         return fEoverPCutMin;
5143         return fEoverPCutMax;
5144
5145 }
5146
5147
5148