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