]> git.uio.no Git - u/mrichter/AliRoot.git/blob - PWGGA/EMCALTasks/AliAnalysisTaskPi0V2.cxx
changes from fzhou
[u/mrichter/AliRoot.git] / PWGGA / EMCALTasks / AliAnalysisTaskPi0V2.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 /* $Id: AliAnalysisTaskPi0V2.cxx 55404 2012-03-29 10:10:19Z fca $ */
17
18 /* AliAnalysisTaskPi0V2.cxx
19  *
20  * Template task producing a P_t spectrum and pseudorapidity distribution.
21  * Includes explanations of physics and primary track selections
22  *
23  * Instructions for adding histograms can be found below, starting with NEW HISTO
24  *
25  * Based on tutorial example from offline pages
26  * Edited by Arvinder Palaha
27  */
28 #include "AliAnalysisTaskPi0V2.h"
29
30 #include "Riostream.h"
31 #include "TChain.h"
32 #include "TTree.h"
33 #include "TH1F.h"
34 #include "TH2F.h"
35 #include "TH3F.h"
36 #include "TCanvas.h"
37 #include "TList.h"
38
39 #include "AliAnalysisTaskSE.h"
40 #include "AliAnalysisManager.h"
41 #include "AliStack.h"
42 #include "AliESDtrackCuts.h"
43 #include "AliESDEvent.h"
44 #include "AliESDInputHandler.h"
45 #include "AliAODEvent.h"
46 #include "AliMCEvent.h"
47 #include "AliVCluster.h"
48
49 #include "AliEventplane.h"
50 #include "AliEMCALGeometry.h"
51 #include "THnSparse.h"
52 #include "TClonesArray.h"
53 #include "TString.h"
54
55 using std::cout;
56 using std::endl;
57
58 ClassImp(AliAnalysisTaskPi0V2)
59
60 //________________________________________________________________________
61 AliAnalysisTaskPi0V2::AliAnalysisTaskPi0V2(const char *name) // All data members should be initialised here
62    :AliAnalysisTaskSE(name),
63     fOutput(0),
64     fESD(0),
65     fTracksName("PicoTrack"), fV1ClusName("CaloCluster"), fV2ClusName("CaloCluster"),
66     fTrigClass("CVLN_|CSEMI_|CCENT|CVHN"),
67     fTracks(0), fV1Clus(0), fV2Clus(0),
68     fRunNumber(-999.),
69     fEvtSelect(1),
70     fVtxCut(15.),
71     fNcellCut(2), fECut(1), fEtaCut(0.65), fM02Cut(0.5), fPi0AsyCut(0), isV1Clus(1),
72     fCentrality(99.),
73     fEPTPC(-999.),
74     fEPTPCreso(0.), 
75     fEPV0(-999.), fEPV0A(-999.), fEPV0C(-999.), fEPV0Ar(-999.), fEPV0Cr(-999.), fEPV0r(-999.),
76     fEPV0AR4(-999.), fEPV0AR5(-999.), fEPV0AR6(-999.), fEPV0AR7(-999.), fEPV0CR0(-999.), fEPV0CR1(-999.), fEPV0CR2(-999.), fEPV0CR3(-999.),
77     hEvtCount(0), hAllcentV0(0), hAllcentV0r(0), hAllcentV0A(0), hAllcentV0C(0), hAllcentTPC(0),
78     h2DcosV0r(0), h2DsinV0r(0), h2DcosV0A(0), h2DsinV0A(0), h2DcosV0C(0), h2DsinV0C(0), h2DcosTPC(0), h2DsinTPC(0), 
79     hEPTPC(0), hresoTPC(0),
80     hEPV0(0), hEPV0A(0), hEPV0C(0), hEPV0Ar(0), hEPV0Cr(0), hEPV0r(0), hEPV0AR4(0), hEPV0AR7(0), hEPV0CR0(0), hEPV0CR3(0),
81     hdifV0Ar_V0Cr(0), hdifV0A_V0CR0(0), hdifV0A_V0CR3(0), hdifV0ACR0_V0CR3(0), hdifV0C_V0AR4(0), hdifV0C_V0AR7(0), hdifV0AR4_V0AR7(0),
82     hdifV0A_V0C(0), hdifV0A_TPC(0), hdifV0C_TPC(0), hdifV0C_V0A(0), 
83     hM02vsPtA(0), hM02vsPtB(0), hClusDxDZA(0), hClusDxDZB(0),
84     hdifEMC_EPV0(0), hdifEMC_EPV0A(0), hdifEMC_EPV0C(0), hdifful_EPV0(0), hdifful_EPV0A(0), hdifful_EPV0C(0), 
85     hdifout_EPV0(0), hdifout_EPV0A(0), hdifout_EPV0C(0), hdifEMC_EPTPC(0), hdifful_EPTPC(0), hdifout_EPTPC(0),
86     hdifClus_EPV0(0), hdifClus_EPV0A(0), hdifClus_EPV0C(0), hdifClus_EPTPC(0),
87     fHEPV0r(0), fHEPV0A(0), fHEPV0C(0), fHEPTPC(0)
88
89 {
90     // Dummy constructor ALWAYS needed for I/O.
91     DefineInput(0, TChain::Class());
92     DefineOutput(1, TList::Class());                                            // for output list
93 }
94
95 //________________________________________________________________________
96 AliAnalysisTaskPi0V2::AliAnalysisTaskPi0V2() // All data members should be initialised here
97    :AliAnalysisTaskSE("default_name"),
98     fOutput(0),
99     fESD(0),
100     fTracksName("PicoTrack"), fV1ClusName("CaloCluster"), fV2ClusName("CaloCluster"),
101     fTrigClass("CVLN_|CSEMI_|CCENT|CVHN"),
102     fTracks(0), fV1Clus(0), fV2Clus(0),
103     fRunNumber(-999.),
104     fEvtSelect(1),
105     fVtxCut(15.),
106     fNcellCut(2), fECut(1), fEtaCut(0.65), fM02Cut(0.5), fPi0AsyCut(0), isV1Clus(1),
107     fCentrality(99.),
108     fEPTPC(-999.),
109     fEPTPCreso(0.),
110     fEPV0(-999.), fEPV0A(-999.), fEPV0C(-999.), fEPV0Ar(-999.), fEPV0Cr(-999.), fEPV0r(-999.),
111     fEPV0AR4(-999.), fEPV0AR5(-999.), fEPV0AR6(-999.), fEPV0AR7(-999.), fEPV0CR0(-999.), fEPV0CR1(-999.), fEPV0CR2(-999.), fEPV0CR3(-999.),
112     hEvtCount(0), hAllcentV0(0), hAllcentV0r(0), hAllcentV0A(0), hAllcentV0C(0), hAllcentTPC(0),
113     h2DcosV0r(0), h2DsinV0r(0), h2DcosV0A(0), h2DsinV0A(0), h2DcosV0C(0), h2DsinV0C(0), h2DcosTPC(0), h2DsinTPC(0),
114     hEPTPC(0), hresoTPC(0),
115     hEPV0(0), hEPV0A(0), hEPV0C(0), hEPV0Ar(0), hEPV0Cr(0), hEPV0r(0), hEPV0AR4(0), hEPV0AR7(0), hEPV0CR0(0), hEPV0CR3(0),
116     hdifV0Ar_V0Cr(0), hdifV0A_V0CR0(0), hdifV0A_V0CR3(0), hdifV0ACR0_V0CR3(0), hdifV0C_V0AR4(0), hdifV0C_V0AR7(0), hdifV0AR4_V0AR7(0),
117     hdifV0A_V0C(0), hdifV0A_TPC(0), hdifV0C_TPC(0), hdifV0C_V0A(0),
118     hM02vsPtA(0), hM02vsPtB(0), hClusDxDZA(0), hClusDxDZB(0),
119     hdifEMC_EPV0(0), hdifEMC_EPV0A(0), hdifEMC_EPV0C(0), hdifful_EPV0(0), hdifful_EPV0A(0), hdifful_EPV0C(0),
120     hdifout_EPV0(0), hdifout_EPV0A(0), hdifout_EPV0C(0), hdifEMC_EPTPC(0), hdifful_EPTPC(0), hdifout_EPTPC(0),
121     hdifClus_EPV0(0), hdifClus_EPV0A(0), hdifClus_EPV0C(0), hdifClus_EPTPC(0),
122     fHEPV0r(0), fHEPV0A(0), fHEPV0C(0), fHEPTPC(0)
123 {
124     // Constructor
125     // Define input and output slots here (never in the dummy constructor)
126     // Input slot #0 works with a TChain - it is connected to the default input container
127     // Output slot #1 writes into a TH1 container
128     DefineInput(0, TChain::Class());
129     DefineOutput(1, TList::Class());                                            // for output list
130 }
131
132 //________________________________________________________________________
133 AliAnalysisTaskPi0V2::~AliAnalysisTaskPi0V2()
134 {
135     // Destructor. Clean-up the output list, but not the histograms that are put inside
136     // (the list is owner and will clean-up these histograms). Protect in PROOF case.
137      delete fOutput;
138 }
139 //_____________________________________________________________________
140 Double_t AliAnalysisTaskPi0V2::GetMaxCellEnergy(const AliVCluster *cluster, Short_t &id) const
141 {
142   // Get maximum energy of attached cell.
143
144   id = -1;
145
146   AliVCaloCells *cells = 0;
147   if (fESD)
148     cells = fESD->GetEMCALCells();
149   if (!cells)
150     return 0;
151
152   Double_t maxe = 0;
153   const Int_t ncells = cluster->GetNCells();
154   for (Int_t i=0; i<ncells; i++) {
155     Double_t e = cells->GetCellAmplitude(TMath::Abs(cluster->GetCellAbsId(i)));
156     if (e>maxe) {
157       maxe = e;
158       id   = cluster->GetCellAbsId(i);
159     }
160   }
161   return maxe;
162 }
163 //_____________________________________________________________________
164 Double_t AliAnalysisTaskPi0V2::GetCrossEnergy(const AliVCluster *cluster, Short_t &idmax) const
165 {
166   // Calculate the energy of cross cells around the leading cell.
167
168   AliVCaloCells *cells = 0;
169   if (fESD)
170     cells = fESD->GetEMCALCells();
171   if (!cells)
172     return 0;
173
174   AliEMCALGeometry *geom = AliEMCALGeometry::GetInstance();
175   if (!geom)
176     return 0;
177
178   Int_t iSupMod = -1;
179   Int_t iTower  = -1;
180   Int_t iIphi   = -1;
181   Int_t iIeta   = -1;
182   Int_t iphi    = -1;
183   Int_t ieta    = -1;
184   Int_t iphis   = -1;
185   Int_t ietas   = -1;
186
187   Double_t crossEnergy = 0.;
188
189   geom->GetCellIndex(idmax,iSupMod,iTower,iIphi,iIeta);
190   geom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi, iIeta,iphis,ietas);
191
192   Int_t ncells = cluster->GetNCells();
193   for (Int_t i=0; i<ncells; i++) {
194     Int_t cellAbsId = cluster->GetCellAbsId(i);
195     geom->GetCellIndex(cellAbsId,iSupMod,iTower,iIphi,iIeta);
196     geom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi, iIeta,iphi,ieta);
197     Int_t aphidiff = TMath::Abs(iphi-iphis);
198     if (aphidiff>1)
199       continue;
200     Int_t aetadiff = TMath::Abs(ieta-ietas);
201     if (aetadiff>1)
202       continue;
203     if ( (aphidiff==1 && aetadiff==0) ||
204         (aphidiff==0 && aetadiff==1) ) {
205       crossEnergy += cells->GetCellAmplitude(cellAbsId);
206     }
207   }
208
209   return crossEnergy;
210 }
211 //_____________________________________________________________________
212 Bool_t AliAnalysisTaskPi0V2::IsWithinFiducialVolume(Short_t id) const
213 {
214   // Check if cell is within given fiducial volume.
215
216   Double_t fNFiducial = 1;
217
218   Int_t iSupMod = -1;
219   Int_t iTower  = -1;
220   Int_t iIphi   = -1;
221   Int_t iIeta   = -1;
222   Int_t iphi    = -1;
223   Int_t ieta    = -1;
224
225   Bool_t okrow = kFALSE;
226   Bool_t okcol = kFALSE;
227
228   AliEMCALGeometry *geom = AliEMCALGeometry::GetInstance();
229   if (!geom)
230     return kFALSE;
231
232   Int_t cellAbsId = id;
233   geom->GetCellIndex(cellAbsId,iSupMod,iTower,iIphi,iIeta);
234   geom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi, iIeta,iphi,ieta);
235
236   // Check rows/phi
237   if (iSupMod < 10) {
238     if (iphi >= fNFiducial && iphi < 24-fNFiducial)
239       okrow = kTRUE;
240   } else {
241     if (iphi >= fNFiducial && iphi < 12-fNFiducial)
242       okrow = kTRUE;
243   }
244   // Check columns/eta
245   Bool_t noEMCALBorderAtEta0 = kTRUE;
246   if (!noEMCALBorderAtEta0) {
247     if (ieta > fNFiducial && ieta < 48-fNFiducial)
248       okcol = kTRUE;
249   } else {
250     if (iSupMod%2==0) {
251       if (ieta >= fNFiducial)
252         okcol = kTRUE;
253     } else {
254       if (ieta < 48-fNFiducial)
255         okcol = kTRUE;
256     }
257   }
258   if (okrow && okcol)
259      return kTRUE;
260
261   return kFALSE;
262 }
263 //______________________________________________________________________
264 Bool_t AliAnalysisTaskPi0V2::IsGoodCluster(const AliVCluster *c) const
265 {
266
267   if(!c)
268     return kFALSE;
269
270   if(c->GetNCells() < fNcellCut)
271    return kFALSE;
272
273   if(c->E() < fECut)
274    return kFALSE;
275
276   Short_t id = -1;
277   Double_t maxE = GetMaxCellEnergy(c, id); 
278      if((1. - double(GetCrossEnergy(c,id))/maxE) > 0.97)
279     return kFALSE;
280
281
282   Float_t pos1[3];
283   c->GetPosition(pos1);
284   TVector3 clsPos(pos1);
285   Double_t eta = clsPos.Eta();
286
287   if(TMath::Abs(eta) > fEtaCut)
288     return kFALSE;  
289
290   if (!IsWithinFiducialVolume(id))
291     return kFALSE;
292
293   if(c->GetM02() >fM02Cut)
294     return kFALSE;
295
296
297   return kTRUE;
298
299 }
300 //________________________________________________________________________________________________
301 Bool_t AliAnalysisTaskPi0V2::IsGoodClusterV1(const AliVCluster *c) const
302 {
303
304   if(!c)
305     return kFALSE;
306
307   if(c->GetNCells() < fNcellCut)
308    return kFALSE;
309
310   if(c->E() < fECut)
311    return kFALSE;
312
313   Short_t id = -1;
314   Double_t maxE = GetMaxCellEnergy(c, id);
315      if((1. - double(GetCrossEnergy(c,id))/maxE) > 0.97)
316     return kFALSE;
317
318
319   Float_t pos1[3];
320   c->GetPosition(pos1);
321   TVector3 clsPos(pos1);
322   Double_t eta = clsPos.Eta();
323
324   if(TMath::Abs(eta) > fEtaCut)
325     return kFALSE;
326
327   if (!IsWithinFiducialVolume(id))
328     return kFALSE;
329
330   if(c->GetM02() <0.5)
331     return kFALSE;
332
333   Double_t dr = TMath::Sqrt(c->GetTrackDx()*c->GetTrackDx() + c->GetTrackDz()*c->GetTrackDz());
334   if(dr>0.025)
335     return kFALSE;
336
337   return kTRUE;
338
339 }
340 //_____________________________________________________________________
341 Bool_t AliAnalysisTaskPi0V2::IsGoodPion(const TLorentzVector &p1, const TLorentzVector &p2) const
342 {
343   // Is good pion?
344
345   if(fPi0AsyCut){
346     Double_t asym = TMath::Abs(p1.E()-p2.E())/(p1.E()+p2.E());
347     if (asym>0.7)
348       return kFALSE;
349   }
350   TLorentzVector pion;
351   pion = p1 + p2;
352   Double_t eta = pion.Eta();
353   if(TMath::Abs(eta) > fEtaCut)
354     return kFALSE;
355
356   return kTRUE;
357 }
358 //_______________________________________________________________________
359 void AliAnalysisTaskPi0V2::FillPion(const TLorentzVector& p1, const TLorentzVector& p2, Double_t EPV0r, Double_t EPV0A, Double_t EPV0C, Double_t EPTPC)
360 {
361   // Fill histogram.
362
363   if (!IsGoodPion(p1,p2))
364     return;
365   TLorentzVector pion;
366   pion = p1 + p2;
367
368   Double_t mass = pion.M();
369   Double_t pt   = pion.Pt();
370   Double_t phi  = pion.Phi();
371
372   Double_t dphiV0   = phi-EPV0r;
373   Double_t dphiV0A  = phi-EPV0A;
374   Double_t dphiV0C  = phi-EPV0C;
375   Double_t dphiTPC  = phi-EPTPC;
376
377   Double_t cos2phiV0  = TMath::Cos(2.*(dphiV0));
378   Double_t cos2phiV0A = TMath::Cos(2.*(dphiV0A));
379   Double_t cos2phiV0C = TMath::Cos(2.*(dphiV0C));
380   Double_t cos2phiTPC = TMath::Cos(2.*(dphiTPC));
381
382   dphiV0  = TVector2::Phi_0_2pi(dphiV0);  if(dphiV0  >TMath::Pi())  dphiV0 -= TMath::Pi();
383   dphiV0A = TVector2::Phi_0_2pi(dphiV0A); if(dphiV0A >TMath::Pi())  dphiV0A -= TMath::Pi();
384   dphiV0C = TVector2::Phi_0_2pi(dphiV0C); if(dphiV0C >TMath::Pi())  dphiV0C -= TMath::Pi();
385   dphiTPC = TVector2::Phi_0_2pi(dphiTPC); if(dphiTPC >TMath::Pi())  dphiTPC -= TMath::Pi();
386
387   Double_t xV0[5]; // Match ndims in fH  V0 EP
388   xV0[0]       = mass;
389   xV0[1]       = pt;
390   xV0[2]       = fCentrality;
391   xV0[3]       = dphiV0;
392   xV0[4]       = cos2phiV0;
393   fHEPV0r->Fill(xV0);
394
395   Double_t xV0A[5]; // Match ndims in fH V0A EP
396   xV0A[0]       = mass;
397   xV0A[1]       = pt;
398   xV0A[2]       = fCentrality;
399   xV0A[3]       = dphiV0A;
400   xV0A[4]       = cos2phiV0A;
401   fHEPV0A->Fill(xV0A);
402
403   Double_t xV0C[5]; // Match ndims in fH V0C EP
404   xV0C[0]       = mass;
405   xV0C[1]       = pt;
406   xV0C[2]       = fCentrality;
407   xV0C[3]       = dphiV0C;
408   xV0C[4]       = cos2phiV0C;
409   fHEPV0C->Fill(xV0C);
410
411   Double_t xTPC[5]; // Match ndims in fH TPC EP
412   xTPC[0]       = mass;
413   xTPC[1]       = pt;
414   xTPC[2]       = fCentrality;
415   xTPC[3]       = dphiTPC;
416   xTPC[4]       = cos2phiTPC;
417   fHEPTPC->Fill(xTPC);
418
419
420 }
421 //________________________________________________________________________________________________________________________________
422 void AliAnalysisTaskPi0V2::FillCluster(const TLorentzVector& p1, Double_t EPV0r, Double_t EPV0A, Double_t EPV0C, Double_t EPTPC)
423 {
424   //cluster(photon) v2 method
425   Double_t Pt   = p1.Pt();
426   Double_t Phi  = p1.Phi();
427
428   Double_t difClusV0 = TVector2::Phi_0_2pi(Phi-EPV0r);   if(difClusV0 >TMath::Pi()) difClusV0  -= TMath::Pi();
429   Double_t difClusV0A = TVector2::Phi_0_2pi(Phi-EPV0A);  if(difClusV0A >TMath::Pi()) difClusV0A -= TMath::Pi();
430   Double_t difClusV0C = TVector2::Phi_0_2pi(Phi-EPV0C);  if(difClusV0C >TMath::Pi()) difClusV0C -= TMath::Pi();
431   Double_t difClusTPC = TVector2::Phi_0_2pi(Phi-EPTPC);  if(difClusTPC >TMath::Pi()) difClusTPC -= TMath::Pi();
432
433   hdifClus_EPV0->Fill(fCentrality,  difClusV0, Pt);
434   hdifClus_EPV0A->Fill(fCentrality, difClusV0A, Pt);
435   hdifClus_EPV0C->Fill(fCentrality, difClusV0C, Pt);
436   hdifClus_EPTPC->Fill(fCentrality, difClusTPC, Pt);
437
438
439 }
440 //_________________________________________________________________________________________________
441 void AliAnalysisTaskPi0V2::GetMom(TLorentzVector& p, const AliVCluster *c, Double_t *vertex)
442 {
443   // Calculate momentum.
444   Float_t posMom[3];
445   c->GetPosition(posMom);
446   TVector3 clsPos2(posMom);
447
448   Double_t e   = c->E();
449   Double_t r   = clsPos2.Perp();
450   Double_t eta = clsPos2.Eta();
451   Double_t phi = clsPos2.Phi();
452
453   TVector3 pos;
454   pos.SetPtEtaPhi(r,eta,phi);
455
456   if (vertex) { //calculate direction relative to vertex
457     pos -= vertex;
458   }
459
460   Double_t rad = pos.Mag();
461   p.SetPxPyPzE(e*pos.x()/rad, e*pos.y()/rad, e*pos.z()/rad, e);
462
463 }
464 //________________________________________________________________________
465 void AliAnalysisTaskPi0V2::UserCreateOutputObjects()
466 {
467     // Create histograms
468     // Called once (on the worker node)
469         
470     fOutput = new TList();
471     fOutput->SetOwner();  // IMPORTANT!
472
473     hEvtCount = new TH1F("hEvtCount", " Event Plane", 10, 0.5, 10.5);
474     hEvtCount->GetXaxis()->SetBinLabel(1,"All");
475     hEvtCount->GetXaxis()->SetBinLabel(2,"Evt Cut");
476     hEvtCount->GetXaxis()->SetBinLabel(3,"Trg Class");
477     hEvtCount->GetXaxis()->SetBinLabel(4,"Vtx");
478     hEvtCount->GetXaxis()->SetBinLabel(5,"Cent");
479     hEvtCount->GetXaxis()->SetBinLabel(5,"EPtask");
480     hEvtCount->GetXaxis()->SetBinLabel(7,"EPvalue");
481     hEvtCount->GetXaxis()->SetBinLabel(8,"Pass");
482     fOutput->Add(hEvtCount);
483     
484     hEPTPC   = new TH2F("hEPTPC",   "EPTPC     vs cent", 100, 0., 100., 100, 0., TMath::Pi());
485     hresoTPC = new TH2F("hresoTPC", "TPc reso  vs cent", 100, 0., 100., 100, 0., 1.);
486     hEPV0    = new TH2F("hEPV0",    "EPV0      vs cent", 100, 0., 100., 100, 0., TMath::Pi());
487     hEPV0A   = new TH2F("hEPV0A",   "EPV0A     vs cent", 100, 0., 100., 100, 0., TMath::Pi());
488     hEPV0C   = new TH2F("hEPV0C",   "EPV0C     vs cent", 100, 0., 100., 100, 0., TMath::Pi());
489     hEPV0Ar  = new TH2F("hEPV0Ar",  "EPV0Ar    vs cent", 100, 0., 100., 100, 0., TMath::Pi());
490     hEPV0Cr  = new TH2F("hEPV0Cr",  "EPV0Cr    vs cent", 100, 0., 100., 100, 0., TMath::Pi());
491     hEPV0r   = new TH2F("hEPV0r",   "EPV0r     vs cent", 100, 0., 100., 100, 0., TMath::Pi());
492     hEPV0AR4 = new TH2F("hEPV0AR4", "EPV0AR4   vs cent", 100, 0., 100., 100, 0., TMath::Pi());
493     hEPV0AR7 = new TH2F("hEPV0AR7", "EPV0AR7   vs cent", 100, 0., 100., 100, 0., TMath::Pi());
494     hEPV0CR0 = new TH2F("hEPV0CR0", "EPV0CR0   vs cent", 100, 0., 100., 100, 0., TMath::Pi());
495     hEPV0CR3 = new TH2F("hEPV0CR3", "EPV0CR3   vs cent", 100, 0., 100., 100, 0., TMath::Pi());
496     fOutput->Add(hEPTPC);
497     fOutput->Add(hresoTPC);
498     fOutput->Add(hEPV0);
499     fOutput->Add(hEPV0A);
500     fOutput->Add(hEPV0C);
501     fOutput->Add(hEPV0Ar);
502     fOutput->Add(hEPV0Cr);
503     fOutput->Add(hEPV0r);
504     fOutput->Add(hEPV0AR4);
505     fOutput->Add(hEPV0AR7);
506     fOutput->Add(hEPV0CR0);
507     fOutput->Add(hEPV0CR3);
508
509     hdifV0Ar_V0Cr    = new TH2F("hdifV0Ar_V0Cr",    "EP Ar-Cr ", 100, 0., 100., 100, -1., 1.);    
510     hdifV0A_V0CR0    = new TH2F("hdifV0A_V0CR0",    "EP A-R0 ",  100, 0., 100., 100, -1., 1.);    
511     hdifV0A_V0CR3    = new TH2F("hdifV0A_V0CR3",    "EP A-R3 ",  100, 0., 100., 100, -1., 1.);    
512     hdifV0ACR0_V0CR3 = new TH2F("hdifV0ACR0_V0CR3", "EP R0-R3 ", 100, 0., 100., 100, -1., 1.);    
513     hdifV0C_V0AR4    = new TH2F("hdifV0C_V0AR4",    "EP C-R4 ",  100, 0., 100., 100, -1., 1.);    
514     hdifV0C_V0AR7    = new TH2F("hdifV0C_V0AR7",    "EP C-R7 ",  100, 0., 100., 100, -1., 1.);    
515     hdifV0AR4_V0AR7  = new TH2F("hdifV0AR4_V0AR7",  "EP R4-R7 ", 100, 0., 100., 100, -1., 1.);    
516     fOutput->Add(hdifV0Ar_V0Cr);
517     fOutput->Add(hdifV0A_V0CR0);
518     fOutput->Add(hdifV0A_V0CR3);
519     fOutput->Add(hdifV0ACR0_V0CR3);
520     fOutput->Add(hdifV0C_V0AR4);
521     fOutput->Add(hdifV0C_V0AR7);
522     fOutput->Add(hdifV0AR4_V0AR7);
523
524     hdifV0A_V0C = new TH2F("hdifV0A_V0C", "EP A-C  ", 100, 0., 100., 100, -1., 1.);
525     hdifV0A_TPC = new TH2F("hdifV0A_TPC", "EP A-TPC", 100, 0., 100., 100, -1., 1.);
526     hdifV0C_TPC = new TH2F("hdifV0C_TPC", "EP C-TPC", 100, 0., 100., 100, -1., 1.);
527     hdifV0C_V0A = new TH2F("hdifV0C_V0A", "EP C-A  ", 100, 0., 100., 100, -1., 1.);
528     fOutput->Add(hdifV0A_V0C);
529     fOutput->Add(hdifV0A_TPC);
530     fOutput->Add(hdifV0C_TPC);
531     fOutput->Add(hdifV0C_V0A);
532
533
534
535     hdifEMC_EPV0  = new TH3F("hdifEMC_EPV0",  "dif phi in EMC with EP",  100, 0., 100., 100, 0., TMath::Pi(), 15, 0., 15.);
536     hdifEMC_EPV0A = new TH3F("hdifEMC_EPV0A", "dif phi in EMC with EP",  100, 0., 100., 100, 0., TMath::Pi(), 15, 0., 15.);
537     hdifEMC_EPV0C = new TH3F("hdifEMC_EPV0C", "dif phi in EMC with EP",  100, 0., 100., 100, 0., TMath::Pi(), 15, 0., 15.);
538     fOutput->Add(hdifEMC_EPV0);
539     fOutput->Add(hdifEMC_EPV0A);
540     fOutput->Add(hdifEMC_EPV0C);
541
542     hdifful_EPV0 = new TH3F("hdifful_EPV0",    "dif phi in full with EP", 100, 0., 100., 100, 0., TMath::Pi(), 15, 0., 15.);
543     hdifful_EPV0A = new TH3F("hdifful_EPV0A",  "dif phi in full with EP", 100, 0., 100., 100, 0., TMath::Pi(), 15, 0., 15.);
544     hdifful_EPV0C = new TH3F("hdifful_EPV0C",  "dif phi in full with EP", 100, 0., 100., 100, 0., TMath::Pi(), 15, 0., 15.);
545     fOutput->Add(hdifful_EPV0);
546     fOutput->Add(hdifful_EPV0A);
547     fOutput->Add(hdifful_EPV0C);
548
549     hdifout_EPV0  = new TH3F("hdifout_EPV0",  "dif phi NOT in EMC with EP", 100, 0., 100., 100, 0., TMath::Pi(), 15, 0., 15.);
550     hdifout_EPV0A = new TH3F("hdifout_EPV0A", "dif phi NOT in EMC with EP", 100, 0., 100., 100, 0., TMath::Pi(), 15, 0., 15.);
551     hdifout_EPV0C = new TH3F("hdifout_EPV0C", "dif phi NOT in EMC with EP", 100, 0., 100., 100, 0., TMath::Pi(), 15, 0., 15.);
552     fOutput->Add(hdifout_EPV0);
553     fOutput->Add(hdifout_EPV0A);
554     fOutput->Add(hdifout_EPV0C);
555
556     hdifEMC_EPTPC = new TH3F("hdifEMC_EPTPC", "dif phi in EMC with EP",  100, 0., 100., 100, 0., TMath::Pi(), 15, 0., 15.);
557     hdifful_EPTPC = new TH3F("hdifful_EPTPC",  "dif phi in full with EP", 100, 0., 100., 100, 0., TMath::Pi(), 15, 0., 15.);
558     hdifout_EPTPC = new TH3F("hdifout_EPTPC", "dif phi NOT in EMC with EP", 100, 0., 100., 100, 0., TMath::Pi(), 15, 0., 15.);
559     fOutput->Add(hdifEMC_EPTPC);
560     fOutput->Add(hdifful_EPTPC);
561     fOutput->Add(hdifout_EPTPC);
562
563     hdifClus_EPV0 = new TH3F("hdifClus_EPV0",   "dif phi in EMC Clus with EP", 100, 0., 100., 100, 0., TMath::Pi(), 90, 5., 50.);
564     hdifClus_EPV0A = new TH3F("hdifClus_EPV0A", "dif phi in EMC Clus with EP", 100, 0., 100., 100, 0., TMath::Pi(), 90, 5., 50.);
565     hdifClus_EPV0C = new TH3F("hdifClus_EPV0C", "dif phi in EMC Clus with EP", 100, 0., 100., 100, 0., TMath::Pi(), 90, 5., 50.);
566     hdifClus_EPTPC = new TH3F("hdifClus_EPTPC", "dif phi in EMC Clus with EP", 100, 0., 100., 100, 0., TMath::Pi(), 90, 5., 50.);
567     fOutput->Add(hdifClus_EPV0);
568     fOutput->Add(hdifClus_EPV0A);
569     fOutput->Add(hdifClus_EPV0C);
570     fOutput->Add(hdifClus_EPTPC);
571
572     hAllcentV0  = new TH1F("hAllcentV0",  "All cent EP V0",  100, 0., TMath::Pi());
573     hAllcentV0r = new TH1F("hAllcentV0r", "All cent EP V0r", 100, 0., TMath::Pi());
574     hAllcentV0A = new TH1F("hAllcentV0A", "All cent EP V0A", 100, 0., TMath::Pi());
575     hAllcentV0C = new TH1F("hAllcentV0C", "All cent EP V0C", 100, 0., TMath::Pi());
576     hAllcentTPC = new TH1F("hAllcentTPC", "All cent EP TPC", 100, 0., TMath::Pi());
577     fOutput->Add(hAllcentV0);
578     fOutput->Add(hAllcentV0r);
579     fOutput->Add(hAllcentV0A);
580     fOutput->Add(hAllcentV0C);
581     fOutput->Add(hAllcentTPC);
582
583     h2DcosV0r = new TH2F("h2DcosV0r", "cos(Phi) V0r vs Run NUmber", 200, 0, 200, 100, -1, 1);
584     h2DsinV0r = new TH2F("h2DsinV0r", "sin(Phi) V0r vs Run NUmber", 200, 0, 200, 100, -1, 1);
585     h2DcosV0A = new TH2F("h2DcosV0A", "cos(Phi) V0r vs Run NUmber", 200, 0, 200, 100, -1, 1);
586     h2DsinV0A = new TH2F("h2DsinV0A", "sin(Phi) V0r vs Run NUmber", 200, 0, 200, 100, -1, 1);
587     h2DcosV0C = new TH2F("h2DcosV0C", "cos(Phi) V0r vs Run NUmber", 200, 0, 200, 100, -1, 1);
588     h2DsinV0C = new TH2F("h2DsinV0C", "sin(Phi) V0r vs Run NUmber", 200, 0, 200, 100, -1, 1);
589     h2DcosTPC = new TH2F("h2DcosTPC", "cos(Phi) V0r vs Run NUmber", 200, 0, 200, 100, -1, 1);
590     h2DsinTPC = new TH2F("h2DsinTPC", "sin(Phi) V0r vs Run NUmber", 200, 0, 200, 100, -1, 1);
591     fOutput->Add(h2DcosV0r);
592     fOutput->Add(h2DsinV0r);
593     fOutput->Add(h2DcosV0A);
594     fOutput->Add(h2DsinV0A);
595     fOutput->Add(h2DcosV0C);
596     fOutput->Add(h2DsinV0C);
597     fOutput->Add(h2DcosTPC);
598     fOutput->Add(h2DsinTPC);
599
600     hM02vsPtA = new TH2F("hM02vsPtA", "M02 vs Et before cut", 5000, 0, 50, 400, 0, 4.);
601     hM02vsPtB = new TH2F("hM02vsPtB", "M02 vs Et before cut", 5000, 0, 50, 400, 0, 4.);
602     fOutput->Add(hM02vsPtA);
603     fOutput->Add(hM02vsPtB);
604
605     hClusDxDZA = new TH2F("hClusDxDZA", "clus Dx vs Dz", 1000, -1., 1., 1000, -1., 1);  
606     hClusDxDZB = new TH2F("hClusDxDZB", "clus Dx vs Dz", 1000, -1., 1., 1000, -1., 1);
607     fOutput->Add(hClusDxDZA);
608     fOutput->Add(hClusDxDZB);
609
610     const Int_t ndims = 5;
611     Int_t nMgg=500, nPt=40, nCent=20, nDeltaPhi=315,  ncos2phi=500;
612     Int_t bins[ndims] = {nMgg, nPt, nCent, nDeltaPhi, ncos2phi};
613     Double_t xmin[ndims] = { 0,   0.,  0,   0.,     -1.};
614     Double_t xmax[ndims] = { 0.5, 20., 100, 3.15,   1.};
615     fHEPV0r  = new THnSparseF("fHEPV0r",  "Flow histogram EPV0",  ndims, bins, xmin, xmax);
616     fHEPV0A = new THnSparseF("fHEPV0A",   "Flow histogram EPV0A", ndims, bins, xmin, xmax);
617     fHEPV0C = new THnSparseF("fHEPV0C",   "Flow histogram EPV0C", ndims, bins, xmin, xmax);
618     fHEPTPC = new THnSparseF("fHEPTPC",   "Flow histogram EPTPC", ndims, bins, xmin, xmax);
619     fOutput->Add(fHEPV0r);
620     fOutput->Add(fHEPV0A);
621     fOutput->Add(fHEPV0C);
622     fOutput->Add(fHEPTPC);
623
624     
625
626     PostData(1, fOutput); // Post data for ALL output slots >0 here, to get at least an empty histogram
627 }
628
629 //________________________________________________________________________
630 void AliAnalysisTaskPi0V2::UserExec(Option_t *) 
631 {
632     // Main loop
633     // Called for each event
634         
635     // Create pointer to reconstructed event
636    AliVEvent *event = InputEvent();
637    if (!event) { Printf("ERROR: Could not retrieve event"); return; }
638    // create pointer to event
639    fESD = dynamic_cast<AliESDEvent*>(event);
640    if (!fESD) {
641      AliError("Cannot get the ESD event");
642      return;
643    }
644    hEvtCount->Fill(1);
645     
646   Int_t AbsRunNumber = fESD->GetRunNumber();
647   fRunNumber = ConvertToInternalRunNumber(AbsRunNumber);
648
649   Bool_t isSelected =0;      
650   if(fEvtSelect == 1){  //MB+SemiCentral
651     isSelected = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & (AliVEvent::kAnyINT | AliVEvent::kSemiCentral));
652   } else if (fEvtSelect == 2){  //MB+Central+SemiCentral
653     isSelected = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & (AliVEvent::kAnyINT | AliVEvent::kSemiCentral | AliVEvent::kCentral));
654   } else if(fEvtSelect == 3){  //MB
655  isSelected = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & (AliVEvent::kAnyINT ));
656   }
657
658   if(!isSelected )
659         return; 
660
661   hEvtCount->Fill(2);
662   if(!fTrigClass.IsNull()){
663     TString fired;
664     fired = fESD->GetFiredTriggerClasses();
665     if (!fired.Contains("-B-"))
666       return;
667     TObjArray *arr = fTrigClass.Tokenize("|");
668     if (!arr)
669       return;
670     Bool_t match = 0;
671     for (Int_t i=0;i<arr->GetEntriesFast();++i) {
672       TObject *obj = arr->At(i);
673       if (!obj)
674         continue;
675       if (fired.Contains(obj->GetName())) {
676         match = 1;
677         break;
678       }
679     }
680     delete arr;
681     if (
682         !match && //select by Trigger classes in KCentral and KSemiCentral
683         !(((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & (AliVEvent::kAnyINT ))              // always accept MB
684         ) 
685       return; //Not match skip this event
686   }
687
688     hEvtCount->Fill(3);
689     const AliESDVertex* fvertex = fESD->GetPrimaryVertex();
690     if(TMath::Abs(fvertex->GetZ())>fVtxCut)
691       return;
692     Double_t vertex[3] = {fvertex->GetX(), fvertex->GetY(), fvertex->GetZ()};
693
694     hEvtCount->Fill(4);
695
696     if(fESD->GetCentrality()) {
697       fCentrality = 
698         fESD->GetCentrality()->GetCentralityPercentile("CL1"); //spd vertex
699     } else{
700            return;
701     }
702
703     hEvtCount->Fill(5);
704     AliEventplane *ep = fESD->GetEventplane();
705     if (ep) {
706       if (ep->GetQVector())
707         fEPTPC    = ep->GetQVector()->Phi()/2. ;
708       else
709         fEPTPC = -999.;
710       if (ep->GetQsub1()&&ep->GetQsub2())
711         fEPTPCreso  = TMath::Cos(2.*(ep->GetQsub1()->Phi()/2.-ep->GetQsub2()->Phi()/2.));
712       else
713         fEPTPCreso = -1;
714
715       fEPV0    = ep->GetEventplane("V0",  fESD);
716       fEPV0A   = ep->GetEventplane("V0A", fESD);
717       fEPV0C   = ep->GetEventplane("V0C", fESD);
718       Double_t qx=0, qy=0;
719       Double_t qxr=0, qyr=0;
720       fEPV0Ar  = ep->CalculateVZEROEventPlane(fESD, 4, 5, 2, qxr, qyr);
721       fEPV0Cr  = ep->CalculateVZEROEventPlane(fESD, 2, 3, 2, qx,  qy);
722       qxr += qx;
723       qyr += qy;
724       fEPV0r   = TMath::ATan2(qyr,qxr)/2.;
725       fEPV0AR4 = ep->CalculateVZEROEventPlane(fESD, 4, 2, qx, qy);
726       fEPV0AR5 = ep->CalculateVZEROEventPlane(fESD, 5, 2, qx, qy);
727       fEPV0AR6 = ep->CalculateVZEROEventPlane(fESD, 6, 2, qx, qy);
728       fEPV0AR7 = ep->CalculateVZEROEventPlane(fESD, 7, 2, qx, qy);
729       fEPV0CR0 = ep->CalculateVZEROEventPlane(fESD, 0, 2, qx, qy);
730       fEPV0CR1 = ep->CalculateVZEROEventPlane(fESD, 1, 2, qx, qy);
731       fEPV0CR2 = ep->CalculateVZEROEventPlane(fESD, 2, 2, qx, qy);
732       fEPV0CR3 = ep->CalculateVZEROEventPlane(fESD, 3, 2, qx, qy);
733     }
734     FillEPQA(); //Fill the EP QA
735
736     hEvtCount->Fill(6);
737
738     if( fEPV0A<-2. || fEPV0C<-2. || fEPTPC<-2. || fEPV0r<-2.) 
739       return;
740
741     hEvtCount->Fill(7);
742
743     fEPV0   = TVector2::Phi_0_2pi(fEPV0);    if(fEPV0>TMath::Pi())   fEPV0  = fEPV0  - TMath::Pi();
744     fEPV0r  = TVector2::Phi_0_2pi(fEPV0r);   if(fEPV0r>TMath::Pi())  fEPV0r = fEPV0r - TMath::Pi();
745     fEPV0A  = TVector2::Phi_0_2pi(fEPV0A);   if(fEPV0A>TMath::Pi())  fEPV0A = fEPV0A - TMath::Pi();
746     fEPV0C  = TVector2::Phi_0_2pi(fEPV0C);   if(fEPV0C>TMath::Pi())  fEPV0C = fEPV0C - TMath::Pi();
747     fEPV0Ar = TVector2::Phi_0_2pi(fEPV0Ar);  if(fEPV0Ar>TMath::Pi()) fEPV0Ar = fEPV0Ar - TMath::Pi();
748     fEPV0Cr = TVector2::Phi_0_2pi(fEPV0Cr);  if(fEPV0Cr>TMath::Pi()) fEPV0Cr = fEPV0Cr - TMath::Pi();
749     fEPV0AR4   = TVector2::Phi_0_2pi(fEPV0AR4);    if(fEPV0AR4>TMath::Pi())   fEPV0AR4  = fEPV0AR4 - TMath::Pi();
750     fEPV0AR7   = TVector2::Phi_0_2pi(fEPV0AR7);    if(fEPV0AR7>TMath::Pi())   fEPV0AR7  = fEPV0AR7 - TMath::Pi();
751     fEPV0CR0   = TVector2::Phi_0_2pi(fEPV0CR0);    if(fEPV0CR0>TMath::Pi())   fEPV0CR0  = fEPV0CR0 - TMath::Pi();
752     fEPV0CR3   = TVector2::Phi_0_2pi(fEPV0CR3);    if(fEPV0CR3>TMath::Pi())   fEPV0CR3  = fEPV0CR3 - TMath::Pi();
753
754    if(fEPTPC != -999.)
755    hEPTPC->Fill(fCentrality,  fEPTPC); 
756    if(fEPTPCreso!=-1) hresoTPC->Fill(fCentrality, fEPTPCreso);
757    hEPV0->Fill(fCentrality,   fEPV0);
758    hEPV0A->Fill(fCentrality,  fEPV0A);
759    hEPV0C->Fill(fCentrality,  fEPV0C);
760    hEPV0Ar->Fill(fCentrality, fEPV0Ar);
761    hEPV0Cr->Fill(fCentrality, fEPV0Cr);
762    hEPV0r->Fill(fCentrality,  fEPV0r);
763    hEPV0AR4->Fill(fCentrality, fEPV0AR4);
764    hEPV0AR7->Fill(fCentrality, fEPV0AR7);
765    hEPV0CR0->Fill(fCentrality, fEPV0CR0);
766    hEPV0CR3->Fill(fCentrality, fEPV0CR3);
767
768    hAllcentV0->Fill(fEPV0);
769    hAllcentV0r->Fill(fEPV0r);
770    hAllcentV0A->Fill(fEPV0A);
771    hAllcentV0C->Fill(fEPV0C);  
772    hAllcentTPC->Fill(fEPTPC);
773
774    hdifV0Ar_V0Cr->Fill(fCentrality, TMath::Cos(2.*(fEPV0Ar - fEPV0Cr)));
775    hdifV0A_V0CR0->Fill(fCentrality, TMath::Cos(2.*(fEPV0A - fEPV0CR0)));
776    hdifV0A_V0CR3->Fill(fCentrality, TMath::Cos(2.*(fEPV0A - fEPV0CR3)));
777    hdifV0ACR0_V0CR3->Fill(fCentrality, TMath::Cos(2*(fEPV0CR0 - fEPV0CR3)));
778    hdifV0C_V0AR4->Fill(fCentrality, TMath::Cos(2*(fEPV0C - fEPV0AR4)));
779    hdifV0C_V0AR7->Fill(fCentrality, TMath::Cos(2*(fEPV0C - fEPV0AR7)));
780    hdifV0AR4_V0AR7->Fill(fCentrality, TMath::Cos(2*(fEPV0AR4 - fEPV0AR7)));
781         
782    hdifV0A_V0C->Fill(fCentrality, TMath::Cos(2*(fEPV0A - fEPV0C)));
783    hdifV0A_TPC->Fill(fCentrality, TMath::Cos(2*(fEPV0A - fEPTPC)));
784    hdifV0C_TPC->Fill(fCentrality, TMath::Cos(2*(fEPV0C - fEPTPC)));
785    hdifV0C_V0A->Fill(fCentrality, TMath::Cos(2*(fEPV0C - fEPV0A)));
786     // Cluster loop for reconstructed event
787
788 //================ for v2 clusterize analysis==============================================
789   if(!isV1Clus){
790     if (!fV2ClusName.IsNull() && !fV2Clus) {
791       fV2Clus = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fV2ClusName));
792       if (!fV2Clus) {
793         AliError(Form("%s: Could not retrieve tracks %s!", GetName(), fV2ClusName.Data()));
794         return;
795       }
796     }
797     Int_t nCluster =  fV2Clus->GetEntries(); 
798     for(Int_t i=0; i<nCluster; ++i){
799       AliVCluster *c1 = static_cast<AliVCluster*>(fV2Clus->At(i));      
800       hClusDxDZA->Fill(c1->GetTrackDz(), c1->GetTrackDx());
801       if(!c1->IsEMCAL()) continue;
802       if(!IsGoodCluster(c1)) continue;
803       hClusDxDZB->Fill(c1->GetTrackDz(), c1->GetTrackDx());
804       TLorentzVector p1;
805       GetMom(p1, c1, vertex);
806       for(Int_t j=i+1; j<nCluster; ++j){
807         AliVCluster *c2 = static_cast<AliVCluster*>(fV2Clus->At(i));      
808         if(!c2->IsEMCAL()) continue;
809         if(!IsGoodCluster(c2)) continue;
810         TLorentzVector p2;
811         GetMom(p2, c2, vertex);
812         FillPion(p1, p2, fEPV0r, fEPV0A, fEPV0C, fEPTPC);
813       }
814     }
815   }
816   //================ for v1 clusterize analysis==============================================
817   if(isV1Clus){
818     if (!fV2ClusName.IsNull() && !fV1Clus) {
819       fV1Clus = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fV1ClusName));
820       if (!fV1Clus) {
821         AliError(Form("%s: Could not retrieve tracks %s!", GetName(), fV1ClusName.Data()));
822         return;
823       }
824     }
825     Int_t nClusterV1 = fV1Clus->GetEntries();
826     for(Int_t i=0; i<nClusterV1; ++i){
827       AliVCluster *c3 = static_cast<AliVCluster*>(fV1Clus->At(i));      
828       hClusDxDZA->Fill(c3->GetTrackDz(), c3->GetTrackDx());
829       Float_t clsPosEt[3] = {0,0,0};
830       c3->GetPosition(clsPosEt);
831       TVector3 clsVec(clsPosEt);
832       Double_t Et = c3->E()*TMath::Sin(clsVec.Theta());
833       hM02vsPtA->Fill(Et, c3->GetM02());
834       if(!c3->IsEMCAL()) continue;
835       if(!IsGoodClusterV1(c3)) continue;
836       hM02vsPtB->Fill(Et, c3->GetM02());
837       hClusDxDZB->Fill(c3->GetTrackDz(), c3->GetTrackDx());
838       TLorentzVector p3;
839       GetMom(p3, c3, vertex);
840       FillCluster(p3, fEPV0r, fEPV0A, fEPV0C, fEPTPC);
841     }
842   }
843
844
845    if (!fTracksName.IsNull() && !fTracks) {
846      fTracks = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fTracksName));
847      if (!fTracks) {
848        AliError(Form("%s: Could not retrieve tracks %s!", GetName(), fTracksName.Data())); 
849        return;
850      }
851   }
852
853    Int_t ntracks = fTracks->GetEntries();
854    for(Int_t i=0; i<ntracks; ++i){
855      AliVTrack *track = static_cast<AliVTrack*>(fTracks->At(i));
856      Double_t tPhi = track->Phi();
857      Double_t tPt  = track->Pt();
858
859      Double_t difTrackV0  = TVector2::Phi_0_2pi(tPhi-fEPV0);   if(difTrackV0  >TMath::Pi()) difTrackV0  -= TMath::Pi();
860      Double_t difTrackV0A = TVector2::Phi_0_2pi(tPhi-fEPV0A);  if(difTrackV0A >TMath::Pi()) difTrackV0A -= TMath::Pi();
861      Double_t difTrackV0C = TVector2::Phi_0_2pi(tPhi-fEPV0C);  if(difTrackV0C >TMath::Pi()) difTrackV0C -= TMath::Pi();
862      Double_t difTrackTPC = TVector2::Phi_0_2pi(tPhi-fEPTPC);  if(difTrackTPC >TMath::Pi()) difTrackTPC -= TMath::Pi();
863      if(track->IsEMCAL()){      
864        hdifEMC_EPV0->Fill(fCentrality, difTrackV0, tPt);
865        hdifEMC_EPV0A->Fill(fCentrality, difTrackV0A, tPt);
866        hdifEMC_EPV0C->Fill(fCentrality, difTrackV0C, tPt);
867        hdifEMC_EPTPC->Fill(fCentrality, difTrackTPC, tPt);
868      }else{
869        hdifout_EPV0->Fill(fCentrality, difTrackV0, tPt);
870        hdifout_EPV0A->Fill(fCentrality, difTrackV0A, tPt);
871        hdifout_EPV0C->Fill(fCentrality, difTrackV0C, tPt);
872        hdifout_EPTPC->Fill(fCentrality, difTrackTPC, tPt);
873      }
874      hdifful_EPV0->Fill(fCentrality,    difTrackV0, tPt);
875      hdifful_EPV0A->Fill(fCentrality,   difTrackV0A, tPt);
876      hdifful_EPV0C->Fill(fCentrality,   difTrackV0C, tPt);
877      hdifful_EPTPC->Fill(fCentrality,   difTrackTPC, tPt);
878     }
879     hEvtCount->Fill(8);
880
881     // NEW HISTO should be filled before this point, as PostData puts the
882     // information for this iteration of the UserExec in the container
883     PostData(1, fOutput);
884 }
885 //____________________________________________________________________
886 Int_t AliAnalysisTaskPi0V2::ConvertToInternalRunNumber(Int_t n)
887 {
888     switch(n){
889     case  170593 : return 179 ;
890     case  170572 : return 178 ;
891     case  170556 : return 177 ;
892     case  170552 : return 176 ;
893     case  170546 : return 175 ;
894     case  170390 : return 174 ;
895     case  170389 : return 173 ;
896     case  170388 : return 172 ;
897     case  170387 : return 171 ;
898     case  170315 : return 170 ;
899     case  170313 : return 169 ;
900     case  170312 : return 168 ;
901     case  170311 : return 167 ;
902     case  170309 : return 166 ;
903     case  170308 : return 165 ;
904     case  170306 : return 164 ;
905     case  170270 : return 163 ;
906     case  170269 : return 162 ;
907     case  170268 : return 161 ;
908     case  170267 : return 160 ;
909     case  170264 : return 159 ;
910     case  170230 : return 158 ;
911     case  170228 : return 157 ;
912     case  170208 : return 156 ;
913     case  170207 : return 155 ;
914     case  170205 : return 154 ;
915     case  170204 : return 153 ;
916     case  170203 : return 152 ;
917     case  170195 : return 151 ;
918     case  170193 : return 150 ;
919     case  170163 : return 149 ;
920     case  170162 : return 148 ;
921     case  170159 : return 147 ;
922     case  170155 : return 146 ;
923     case  170152 : return 145 ;
924     case  170091 : return 144 ;
925     case  170089 : return 143 ;
926     case  170088 : return 142 ;
927     case  170085 : return 141 ;
928     case  170084 : return 140 ;
929     case  170083 : return 139 ;
930     case  170081 : return 138 ;
931     case  170040 : return 137 ;
932     case  170038 : return 136 ;
933     case  170036 : return 135 ;
934     case  170027 : return 134 ;
935     case  169981 : return 133 ;
936     case  169975 : return 132 ;
937     case  169969 : return 131 ;
938     case  169965 : return 130 ;
939     case  169961 : return 129 ;
940     case  169956 : return 128 ;
941     case  169926 : return 127 ;
942     case  169924 : return 126 ;
943     case  169923 : return 125 ;
944     case  169922 : return 124 ;
945     case  169919 : return 123 ;
946     case  169918 : return 122 ;
947     case  169914 : return 121 ;
948     case  169859 : return 120 ;
949     case  169858 : return 119 ;
950     case  169855 : return 118 ;
951     case  169846 : return 117 ;
952     case  169838 : return 116 ;
953     case  169837 : return 115 ;
954     case  169835 : return 114 ;
955     case  169683 : return 113 ;
956     case  169628 : return 112 ;
957     case  169591 : return 111 ;
958     case  169590 : return 110 ;
959     case  169588 : return 109 ;
960     case  169587 : return 108 ;
961     case  169586 : return 107 ;
962     case  169584 : return 106 ;
963     case  169557 : return 105 ;
964     case  169555 : return 104 ;
965     case  169554 : return 103 ;
966     case  169553 : return 102 ;
967     case  169550 : return 101 ;
968     case  169515 : return 100 ;
969     case  169512 : return 99 ;
970     case  169506 : return 98 ;
971     case  169504 : return 97 ;
972     case  169498 : return 96 ;
973     case  169475 : return 95 ;
974     case  169420 : return 94 ;
975     case  169419 : return 93 ;
976     case  169418 : return 92 ;
977     case  169417 : return 91 ;
978     case  169415 : return 90 ;
979     case  169411 : return 89 ;
980     case  169238 : return 88 ;
981     case  169236 : return 87 ;
982     case  169167 : return 86 ;
983     case  169160 : return 85 ;
984     case  169156 : return 84 ;
985     case  169148 : return 83 ;
986     case  169145 : return 82 ;
987     case  169144 : return 81 ;
988     case  169143 : return 80 ;
989     case  169138 : return 79 ;
990     case  169099 : return 78 ;
991     case  169094 : return 77 ;
992     case  169091 : return 76 ;
993     case  169045 : return 75 ;
994     case  169044 : return 74 ;
995     case  169040 : return 73 ;
996     case  169035 : return 72 ;
997     case  168992 : return 71 ;
998     case  168988 : return 70 ;
999     case  168984 : return 69 ;
1000     case  168826 : return 68 ;
1001     case  168777 : return 67 ;
1002     case  168514 : return 66 ;
1003     case  168512 : return 65 ;
1004     case  168511 : return 64 ;
1005     case  168467 : return 63 ;
1006     case  168464 : return 62 ;
1007     case  168461 : return 61 ;
1008     case  168460 : return 60 ;
1009     case  168458 : return 59 ;
1010     case  168362 : return 58 ;
1011     case  168361 : return 57 ;
1012     case  168356 : return 56 ;
1013     case  168342 : return 55 ;
1014     case  168341 : return 54 ;
1015     case  168325 : return 53 ;
1016     case  168322 : return 52 ;
1017     case  168318 : return 51 ;
1018     case  168311 : return 50 ;
1019     case  168310 : return 49 ;
1020     case  168213 : return 48 ;
1021     case  168212 : return 47 ;
1022     case  168208 : return 46 ;
1023     case  168207 : return 45 ;
1024     case  168206 : return 44 ;
1025     case  168205 : return 43 ;
1026     case  168204 : return 42 ;
1027     case  168203 : return 41 ;
1028     case  168181 : return 40 ;
1029     case  168177 : return 39 ;
1030     case  168175 : return 38 ;
1031     case  168173 : return 37 ;
1032     case  168172 : return 36 ;
1033     case  168171 : return 35 ;
1034     case  168115 : return 34 ;
1035     case  168108 : return 33 ;
1036     case  168107 : return 32 ;
1037     case  168105 : return 31 ;
1038     case  168104 : return 30 ;
1039     case  168103 : return 29 ;
1040     case  168076 : return 28 ;
1041     case  168069 : return 27 ;
1042     case  168068 : return 26 ;
1043     case  168066 : return 25 ;
1044     case  167988 : return 24 ;
1045     case  167987 : return 23 ;
1046     case  167986 : return 22 ;
1047     case  167985 : return 21 ;
1048     case  167921 : return 20 ;
1049     case  167920 : return 19 ;
1050     case  167915 : return 18 ;
1051     case  167909 : return 17 ;
1052     case  167903 : return 16 ;
1053     case  167902 : return 15 ;
1054     case  167818 : return 14 ;
1055     case  167814 : return 13 ;
1056     case  167813 : return 12 ;
1057     case  167808 : return 11 ;
1058     case  167807 : return 10 ;
1059     case  167806 : return 9 ;
1060     case  167713 : return 8 ;
1061     case  167712 : return 7 ;
1062     case  167711 : return 6 ;
1063     case  167706 : return 5 ;
1064     case  167693 : return 4 ;
1065     case  166532 : return 3 ;
1066     case  166530 : return 2 ;
1067     case  166529 : return 1 ;
1068
1069     default : return 199;
1070     }
1071 }
1072 //_______________________________________________________________________
1073 void AliAnalysisTaskPi0V2::FillEPQA()
1074 {
1075   
1076   h2DcosV0r->Fill(fRunNumber, TMath::Cos(fEPV0r));
1077   h2DsinV0r->Fill(fRunNumber, TMath::Sin(fEPV0r));
1078   h2DcosV0A->Fill(fRunNumber, TMath::Cos(fEPV0A));
1079   h2DsinV0A->Fill(fRunNumber, TMath::Sin(fEPV0A));
1080   h2DcosV0C->Fill(fRunNumber, TMath::Cos(fEPV0C));
1081   h2DsinV0C->Fill(fRunNumber, TMath::Sin(fEPV0C));
1082   h2DcosTPC->Fill(fRunNumber, TMath::Cos(fEPTPC));
1083   h2DsinTPC->Fill(fRunNumber, TMath::Sin(fEPTPC));
1084
1085
1086 }
1087 //________________________________________________________________________
1088 void AliAnalysisTaskPi0V2::Terminate(Option_t *) 
1089 {
1090     // Draw result to screen, or perform fitting, normalizations
1091     // Called once at the end of the query
1092 //    fOutput = dynamic_cast<TList*> (GetOutputData(1));
1093  //   if(!fOutput) { Printf("ERROR: could not retrieve TList fOutput"); return; }
1094         
1095     // Get the physics selection histograms with the selection statistics
1096     //AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
1097     //AliESDInputHandler *inputH = dynamic_cast<AliESDInputHandler*>(mgr->GetInputEventHandler());
1098     //TH2F *histStat = (TH2F*)inputH->GetStatistics();
1099    
1100    
1101     // NEW HISTO should be retrieved from the TList container in the above way,
1102     // so it is available to draw on a canvas such as below
1103
1104 }