include check for centrality == 0
[u/mrichter/AliRoot.git] / PWGGA / EMCALTasks / AliAnalysisTaskPi0V2.cxx
CommitLineData
3c40321c 1/* $Id: AliAnalysisTaskPi0V2.cxx 55404 2012-03-29 10:10:19Z fca $ */
2
3c40321c 3#include "AliAnalysisTaskPi0V2.h"
4
70d53162 5#include <Riostream.h>
6#include <TCanvas.h>
7#include <TChain.h>
8#include <TH1F.h>
9#include <TH2F.h>
10#include <THnSparse.h>
11#include <TList.h>
12#include <TTree.h>
3c40321c 13
70d53162 14#include "AliAODEvent.h"
3c40321c 15#include "AliAnalysisManager.h"
70d53162 16#include "AliAnalysisTaskSE.h"
17#include "AliEMCALGeometry.h"
3c40321c 18#include "AliESDEvent.h"
19#include "AliESDInputHandler.h"
70d53162 20#include "AliESDtrackCuts.h"
3c40321c 21#include "AliEventplane.h"
70d53162 22#include "AliMCEvent.h"
23#include "AliStack.h"
3c40321c 24
25ClassImp(AliAnalysisTaskPi0V2)
26
27//________________________________________________________________________
70d53162 28AliAnalysisTaskPi0V2::AliAnalysisTaskPi0V2() :
29 AliAnalysisTaskSE(),
30 fOutput(0),
31 fTrackCuts(0),
32 fESD(0),
33 fcheckEP2sub(1),
34 fCentrality(99.),
35 fEPTPC(-999.),
36 fEPTPCreso(0.),
37 fEPV0(-999.), fEPV0A(-999.), fEPV0C(-999.), fEPV0Ar(-999.), fEPV0Cr(-999.), fEPV0r(-999.),
38 fEPV0AR4(-999.), fEPV0AR5(-999.), fEPV0AR6(-999.), fEPV0AR7(-999.), fEPV0CR0(-999.), fEPV0CR1(-999.), fEPV0CR2(-999.), fEPV0CR3(-999.),
39 hEPTPC(0), hresoTPC(0),
40 hEPV0(0), hEPV0A(0), hEPV0C(0), hEPV0Ar(0), hEPV0Cr(0), hEPV0r(0), hEPV0AR4(0), hEPV0AR7(0), hEPV0CR0(0), hEPV0CR3(0),
41 hdifV0A_V0CR0(0), hdifV0A_V0CR3(0), hdifV0ACR0_V0CR3(0), hdifV0C_V0AR4(0), hdifV0C_V0AR7(0), hdifV0AR4_V0AR7(0),
42 fHEPV0r(0), fHEPV0A(0), fHEPV0C(0), fHEPTPC(0)
3c40321c 43{
70d53162 44 // Dummy constructor.
3c40321c 45}
46
47//________________________________________________________________________
70d53162 48AliAnalysisTaskPi0V2::AliAnalysisTaskPi0V2(const char *name) :
49 AliAnalysisTaskSE(name),
50 fOutput(0),
51 fTrackCuts(0),
52 fESD(0),
53 fcheckEP2sub(1),
54 fCentrality(99.),
55 fEPTPC(-999.),
56 fEPTPCreso(0.),
57 fEPV0(-999.), fEPV0A(-999.), fEPV0C(-999.), fEPV0Ar(-999.), fEPV0Cr(-999.), fEPV0r(-999.),
58 fEPV0AR4(-999.), fEPV0AR5(-999.), fEPV0AR6(-999.), fEPV0AR7(-999.), fEPV0CR0(-999.), fEPV0CR1(-999.), fEPV0CR2(-999.), fEPV0CR3(-999.),
59 hEPTPC(0), hresoTPC(0),
60 hEPV0(0), hEPV0A(0), hEPV0C(0), hEPV0Ar(0), hEPV0Cr(0), hEPV0r(0), hEPV0AR4(0), hEPV0AR7(0), hEPV0CR0(0), hEPV0CR3(0),
61 hdifV0A_V0CR0(0), hdifV0A_V0CR3(0), hdifV0ACR0_V0CR3(0), hdifV0C_V0AR4(0), hdifV0C_V0AR7(0), hdifV0AR4_V0AR7(0),
62 fHEPV0r(0), fHEPV0A(0), fHEPV0C(0), fHEPTPC(0)
3c40321c 63{
70d53162 64 // Constructor
65 // Define input and output slots here (never in the dummy constructor)
66 // Input slot #0 works with a TChain - it is connected to the default input container
67 // Output slot #1 writes into a TH1 container
68 DefineInput(0, TChain::Class());
69 DefineOutput(1, TList::Class());
3c40321c 70}
71
72//________________________________________________________________________
73AliAnalysisTaskPi0V2::~AliAnalysisTaskPi0V2()
74{
70d53162 75 // Destructor. Clean-up the output list, but not the histograms that are put inside
76 // (the list is owner and will clean-up these histograms). Protect in PROOF case.
77 if (fOutput && !AliAnalysisManager::GetAnalysisManager()->IsProofMode()) {
78 delete fOutput;
79 }
80 delete fTrackCuts;
3c40321c 81}
70d53162 82
3c40321c 83//_____________________________________________________________________
84Double_t AliAnalysisTaskPi0V2::GetMaxCellEnergy(const AliVCluster *cluster, Short_t &id) const
85{
86 // Get maximum energy of attached cell.
87
88 id = -1;
89
90 AliVCaloCells *cells = 0;
91 if (fESD)
92 cells = fESD->GetEMCALCells();
93// else
94// cells = fAOD->GetEMCALCells();
95 if (!cells)
96 return 0;
97
98 Double_t maxe = 0;
99 const Int_t ncells = cluster->GetNCells();
100 for (Int_t i=0; i<ncells; i++) {
101 Double_t e = cells->GetCellAmplitude(TMath::Abs(cluster->GetCellAbsId(i)));
102 if (e>maxe) {
103 maxe = e;
104 id = cluster->GetCellAbsId(i);
105 }
106 }
107 return maxe;
108}
70d53162 109
3c40321c 110//_____________________________________________________________________
111Double_t AliAnalysisTaskPi0V2::GetCrossEnergy(const AliVCluster *cluster, Short_t &idmax) const
112{
113 // Calculate the energy of cross cells around the leading cell.
114
115 AliVCaloCells *cells = 0;
116 if (fESD)
117 cells = fESD->GetEMCALCells();
118// else
119// cells = fAOD->GetEMCALCells();
120 if (!cells)
121 return 0;
122
123 AliEMCALGeometry *geom = AliEMCALGeometry::GetInstance();
124 if (!geom)
125 return 0;
126
127 Int_t iSupMod = -1;
128 Int_t iTower = -1;
129 Int_t iIphi = -1;
130 Int_t iIeta = -1;
131 Int_t iphi = -1;
132 Int_t ieta = -1;
133 Int_t iphis = -1;
134 Int_t ietas = -1;
135
136 Double_t crossEnergy = 0;
137
138 geom->GetCellIndex(idmax,iSupMod,iTower,iIphi,iIeta);
139 geom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi, iIeta,iphis,ietas);
140
141 Int_t ncells = cluster->GetNCells();
142 for (Int_t i=0; i<ncells; i++) {
143 Int_t cellAbsId = cluster->GetCellAbsId(i);
144 geom->GetCellIndex(cellAbsId,iSupMod,iTower,iIphi,iIeta);
145 geom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi, iIeta,iphi,ieta);
146 Int_t aphidiff = TMath::Abs(iphi-iphis);
147 if (aphidiff>1)
148 continue;
149 Int_t aetadiff = TMath::Abs(ieta-ietas);
150 if (aetadiff>1)
151 continue;
152 if ( (aphidiff==1 && aetadiff==0) ||
153 (aphidiff==0 && aetadiff==1) ) {
154 crossEnergy += cells->GetCellAmplitude(cellAbsId);
155 }
156 }
157
158 return crossEnergy;
159}
70d53162 160
3c40321c 161//_____________________________________________________________________
162Bool_t AliAnalysisTaskPi0V2::IsWithinFiducialVolume(Short_t id) const
163{
164 // Check if cell is within given fiducial volume.
165
166 Double_t fNFiducial = 1;
167
168 Int_t iSupMod = -1;
169 Int_t iTower = -1;
170 Int_t iIphi = -1;
171 Int_t iIeta = -1;
172 Int_t iphi = -1;
173 Int_t ieta = -1;
174
175 Bool_t okrow = kFALSE;
176 Bool_t okcol = kFALSE;
177
178 AliEMCALGeometry *geom = AliEMCALGeometry::GetInstance("EMCAL_COMPLETEV1");
179 if (!geom)
180 return kFALSE;
181
182 Int_t cellAbsId = id;
183 geom->GetCellIndex(cellAbsId,iSupMod,iTower,iIphi,iIeta);
184 geom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi, iIeta,iphi,ieta);
185
186 // Check rows/phi
187 if (iSupMod < 10) {
188 if (iphi >= fNFiducial && iphi < 24-fNFiducial)
189 okrow = kTRUE;
190 } else {
191 if (iphi >= fNFiducial && iphi < 12-fNFiducial)
192 okrow = kTRUE;
193 }
194 // Check columns/eta
195 Bool_t noEMCALBorderAtEta0 = kTRUE;
196 if (!noEMCALBorderAtEta0) {
197 if (ieta > fNFiducial && ieta < 48-fNFiducial)
198 okcol = kTRUE;
199 } else {
200 if (iSupMod%2==0) {
201 if (ieta >= fNFiducial)
202 okcol = kTRUE;
203 } else {
204 if (ieta < 48-fNFiducial)
205 okcol = kTRUE;
206 }
207 }
208 if (okrow && okcol)
209 return kTRUE;
210
211 return kFALSE;
212}
70d53162 213
3c40321c 214//______________________________________________________________________
215Bool_t AliAnalysisTaskPi0V2::IsGoodCluster(const AliESDCaloCluster *c) const
216{
70d53162 217 // Check if cluster is good cluster.
3c40321c 218
219 if(!c)
220 return kFALSE;
221
222 if(c->GetNCells() < 2)
223 return kFALSE;
224
225 if(c->E() < 1.)
226 return kFALSE;
227
228 Short_t id = -1;
229 Double_t maxE = GetMaxCellEnergy(c, id);
230 if((1. - GetCrossEnergy(c,id) / maxE) > 0.97)
231 return kFALSE;
232
233 Float_t pos1[3];
234 c->GetPosition(pos1);
235 TVector3 clsPos(pos1);
236 Double_t eta = clsPos.Eta();
237
238 if(eta > 0.65 || eta < -0.65)
239 return kFALSE;
240
241 if (!IsWithinFiducialVolume(id))
242 return kFALSE;
243
244 if(c->GetM02() >0.5 )
245 return kFALSE;
246
3c40321c 247 return kTRUE;
3c40321c 248}
70d53162 249
3c40321c 250//_____________________________________________________________________
251Bool_t AliAnalysisTaskPi0V2::IsGoodPion(const TLorentzVector &p1, const TLorentzVector &p2) const
252{
253 // Is good pion?
254
3c40321c 255 Double_t asym = TMath::Abs(p1.E()-p2.E())/(p1.E()+p2.E());
256 if (asym>0.7)
257 return kFALSE;
258
259// if (TMath::Abs(p1.Eta()-p2.Eta())>0.5)
260// return kFALSE;
261
262 TLorentzVector pion;
263 pion = p1 + p2;
264 Double_t eta = pion.Eta();
265 if (eta<-0.65)
266 return kFALSE;
267 if (eta>0.65)
268 return kFALSE;
269
270 return kTRUE;
271}
70d53162 272
3c40321c 273//_______________________________________________________________________
274void AliAnalysisTaskPi0V2::FillPion(const TLorentzVector& p1, const TLorentzVector& p2, Double_t EPV0r, Double_t EPV0A, Double_t EPV0C, Double_t EPTPC)
275{
276 // Fill histogram.
277
278 if (!IsGoodPion(p1,p2))
279 return;
280 TLorentzVector pion;
281 pion = p1 + p2;
282
283 Double_t mass = pion.M();
284 Double_t pt = pion.Pt();
285 Double_t phi = pion.Phi();
286
287 Double_t dphiV0 = phi-EPV0r;
288 Double_t dphiV0A = phi-EPV0A;
289 Double_t dphiV0C = phi-EPV0C;
290 Double_t dphiTPC = phi-EPTPC;
291
292 Double_t cos2phiV0 = TMath::Cos(2.*(dphiV0));
293 Double_t cos2phiV0A = TMath::Cos(2.*(dphiV0A));
294 Double_t cos2phiV0C = TMath::Cos(2.*(dphiV0C));
295 Double_t cos2phiTPC = TMath::Cos(2.*(dphiTPC));
296
297 dphiV0 = TVector2::Phi_0_2pi(dphiV0);
298 dphiV0A = TVector2::Phi_0_2pi(dphiV0A);
299 dphiV0C = TVector2::Phi_0_2pi(dphiV0C);
300 dphiTPC = TVector2::Phi_0_2pi(dphiTPC);
301
302 if(dphiV0 > TMath::Pi()) dphiV0 -= TMath::Pi();
303 if(dphiV0A > TMath::Pi()) dphiV0A -= TMath::Pi();
304 if(dphiV0C > TMath::Pi()) dphiV0C -= TMath::Pi();
305 if(dphiTPC > TMath::Pi()) dphiTPC -= TMath::Pi();
306
307 Double_t xV0[5]; // Match ndims in fH V0 EP
308 xV0[0] = mass;
309 xV0[1] = pt;
310 xV0[2] = fCentrality;
311 xV0[3] = dphiV0;
312 xV0[4] = cos2phiV0;
313 fHEPV0r->Fill(xV0);
314
315 Double_t xV0A[5]; // Match ndims in fH V0A EP
316 xV0A[0] = mass;
317 xV0A[1] = pt;
318 xV0A[2] = fCentrality;
319 xV0A[3] = dphiV0A;
320 xV0A[4] = cos2phiV0A;
321 fHEPV0A->Fill(xV0A);
322
323 Double_t xV0C[5]; // Match ndims in fH V0C EP
324 xV0C[0] = mass;
325 xV0C[1] = pt;
326 xV0C[2] = fCentrality;
327 xV0C[3] = dphiV0C;
328 xV0C[4] = cos2phiV0C;
329 fHEPV0C->Fill(xV0C);
330
331 Double_t xTPC[5]; // Match ndims in fH TPC EP
332 xTPC[0] = mass;
333 xTPC[1] = pt;
334 xTPC[2] = fCentrality;
335 xTPC[3] = dphiTPC;
336 xTPC[4] = cos2phiTPC;
337 fHEPTPC->Fill(xTPC);
3c40321c 338}
70d53162 339
3c40321c 340//_________________________________________________________________________________________________
341void AliAnalysisTaskPi0V2::GetMom(TLorentzVector& p, const AliESDCaloCluster *c, Double_t *vertex)
342{
343 // Calculate momentum.
344 Float_t posMom[3];
345 c->GetPosition(posMom);
346 TVector3 clsPos2(posMom);
347
348 Double_t e = c->E();
349 Double_t r = clsPos2.Perp();
350 Double_t eta = clsPos2.Eta();
351 Double_t phi = clsPos2.Phi();
352
353 TVector3 pos;
354 pos.SetPtEtaPhi(r,eta,phi);
355
356 if (vertex) { //calculate direction relative to vertex
357 pos -= vertex;
358 }
359
360 Double_t rad = pos.Mag();
361 p.SetPxPyPzE(e*pos.x()/rad, e*pos.y()/rad, e*pos.z()/rad, e);
3c40321c 362}
70d53162 363
3c40321c 364//________________________________________________________________________
365void AliAnalysisTaskPi0V2::UserCreateOutputObjects()
366{
70d53162 367 // Create histograms.
3c40321c 368
70d53162 369 fOutput = new TList();
370 fOutput->SetOwner();
3c40321c 371
70d53162 372 fTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(kTRUE);
373
374 hEPTPC = new TH2F("hEPTPC", "EPTPC vs cent", 100, 0., 100., 152, 0., 3.2);
375 hresoTPC = new TH2F("hresoTPC", "TPc reso vs cent", 100, 0., 100., 152, 0., 3.2);
376 hEPV0 = new TH2F("hEPV0", "EPV0 vs cent", 100, 0., 100., 152, 0., 3.2);
377 hEPV0A = new TH2F("hEPV0A", "EPV0A vs cent", 100, 0., 100., 152, 0., 3.2);
378 hEPV0C = new TH2F("hEPV0C", "EPV0C vs cent", 100, 0., 100., 152, 0., 3.2);
379 hEPV0Ar = new TH2F("hEPV0Ar", "EPV0Ar vs cent", 100, 0., 100., 152, 0., 3.2);
380 hEPV0Cr = new TH2F("hEPV0Cr", "EPV0Cr vs cent", 100, 0., 100., 152, 0., 3.2);
381 hEPV0r = new TH2F("hEPV0r", "EPV0r vs cent", 100, 0., 100., 152, 0., 3.2);
382 hEPV0AR4 = new TH2F("hEPV0AR4", "EPV0AR4 vs cent", 100, 0., 100., 152, 0., 3.2);
383 hEPV0AR7 = new TH2F("hEPV0AR7", "EPV0AR7 vs cent", 100, 0., 100., 152, 0., 3.2);
384 hEPV0CR0 = new TH2F("hEPV0CR0", "EPV0CR0 vs cent", 100, 0., 100., 152, 0., 3.2);
385 hEPV0CR3 = new TH2F("hEPV0CR3", "EPV0CR3 vs cent", 100, 0., 100., 152, 0., 3.2);
386 fOutput->Add(hEPTPC);
387 fOutput->Add(hresoTPC);
388 fOutput->Add(hEPV0);
389 fOutput->Add(hEPV0A);
390 fOutput->Add(hEPV0C);
391 fOutput->Add(hEPV0Ar);
392 fOutput->Add(hEPV0Cr);
393 fOutput->Add(hEPV0r);
394 fOutput->Add(hEPV0AR4);
395 fOutput->Add(hEPV0AR7);
396 fOutput->Add(hEPV0CR0);
397 fOutput->Add(hEPV0CR3);
398
399 hdifV0A_V0CR0 = new TH2F("hdifV0A_V0CR0", "EP A-R0 ", 100, 0., 100., 152, 0., 3.2);
400 hdifV0A_V0CR3 = new TH2F("hdifV0A_V0CR3", "EP A-R3 ", 100, 0., 100., 152, 0., 3.2);
401 hdifV0ACR0_V0CR3 = new TH2F("hdifV0ACR0_V0CR3", "EP R0-R3 ", 100, 0., 100., 152, 0., 3.2);
402 hdifV0C_V0AR4 = new TH2F("hdifV0C_V0AR4", "EP C-R4 ", 100, 0., 100., 152, 0., 3.2);
403 hdifV0C_V0AR7 = new TH2F("hdifV0C_V0AR7", "EP C-R7 ", 100, 0., 100., 152, 0., 3.2);
404 hdifV0AR4_V0AR7 = new TH2F("hdifV0AR4_V0AR7", "EP R4-R7 ", 100, 0., 100., 152, 0., 3.2);
405 fOutput->Add(hdifV0A_V0CR0);
406 fOutput->Add(hdifV0A_V0CR3);
407 fOutput->Add(hdifV0ACR0_V0CR3);
408 fOutput->Add(hdifV0C_V0AR4);
409 fOutput->Add(hdifV0C_V0AR7);
410 fOutput->Add(hdifV0AR4_V0AR7);
411
412 const Int_t ndims = 5;
413 Int_t nMgg=500, nPt=40, nCent=20, nDeltaPhi=315, ncos2phi=500;
414 Int_t bins[ndims] = {nMgg, nPt, nCent, nDeltaPhi, ncos2phi};
415 Double_t xmin[ndims] = { 0, 0., 0, 0., -1.};
416 Double_t xmax[ndims] = { 0.5, 20., 100, 3.15, 1.};
417 fHEPV0r = new THnSparseF("fHEPV0r", "Flow histogram EPV0", ndims, bins, xmin, xmax);
418 fHEPV0A = new THnSparseF("fHEPV0A", "Flow histogram EPV0A", ndims, bins, xmin, xmax);
419 fHEPV0C = new THnSparseF("fHEPV0C", "Flow histogram EPV0C", ndims, bins, xmin, xmax);
420 fHEPTPC = new THnSparseF("fHEPTPC", "Flow histogram EPTPC", ndims, bins, xmin, xmax);
421 fOutput->Add(fHEPV0r);
422 fOutput->Add(fHEPV0A);
423 fOutput->Add(fHEPV0C);
424 fOutput->Add(fHEPTPC);
425
426 PostData(1, fOutput); // Post data for ALL output slots >0 here, to get at least an empty histogram
3c40321c 427}
428
429//________________________________________________________________________
430void AliAnalysisTaskPi0V2::UserExec(Option_t *)
431{
70d53162 432 // Main loop.
433
434 // Create pointer to reconstructed event
435 AliVEvent *event = InputEvent();
436 if (!event) { Printf("ERROR: Could not retrieve event"); return; }
437
438 // create pointer to event
439 fESD = dynamic_cast<AliESDEvent*>(event);
440 if (!fESD) {
441 AliError("Cannot get the ESD event");
442 return;
443 }
444 const AliESDVertex* fvertex = fESD->GetPrimaryVertex();
3c40321c 445
70d53162 446 if(!(fvertex->GetStatus())) return;
447 // if vertex is from spd vertexZ, require more stringent cut
448 if (fvertex->IsFromVertexerZ()) {
449 if (fvertex->GetDispersion()>0.02 || fvertex->GetZRes()>0.25 ) return; // bad vertex from VertexerZ
450 }
451 Double_t vertex[3] = {fvertex->GetX(), fvertex->GetY(), fvertex->GetZ()};
452
453 if(fESD->GetCentrality()) {
454 fCentrality =
455 fESD->GetCentrality()->GetCentralityPercentile("V0M");
456 }
3c40321c 457
70d53162 458 AliEventplane *ep = fESD->GetEventplane();
459 if (ep) {
460 if (ep->GetQVector())
461 fEPTPC = ep->GetQVector()->Phi()/2. ;
462 else
463 fEPTPC = -999.;
464 if (ep->GetQsub1()&&ep->GetQsub2())
465 fEPTPCreso = TMath::Cos(2.*(ep->GetQsub1()->Phi()/2.-ep->GetQsub2()->Phi()/2.));
466 else
467 fEPTPCreso = -1;
468
469 fEPV0 = ep->GetEventplane("V0", fESD);
470 fEPV0A = ep->GetEventplane("V0A", fESD);
471 fEPV0C = ep->GetEventplane("V0C", fESD);
472 Double_t qx=0, qy=0;
473 Double_t qxr=0, qyr=0;
474 fEPV0Ar = ep->CalculateVZEROEventPlane(fESD, 4, 5, 2, qxr, qyr);
475 fEPV0Cr = ep->CalculateVZEROEventPlane(fESD, 2, 3, 2, qx, qy);
476 qxr += qx;
477 qyr += qy;
478 fEPV0r = TMath::ATan2(qyr,qxr)/2.;
479 fEPV0AR4 = ep->CalculateVZEROEventPlane(fESD, 4, 2, qx, qy);
480 fEPV0AR5 = ep->CalculateVZEROEventPlane(fESD, 5, 2, qx, qy);
481 fEPV0AR6 = ep->CalculateVZEROEventPlane(fESD, 6, 2, qx, qy);
482 fEPV0AR7 = ep->CalculateVZEROEventPlane(fESD, 7, 2, qx, qy);
483 fEPV0CR0 = ep->CalculateVZEROEventPlane(fESD, 0, 2, qx, qy);
484 fEPV0CR1 = ep->CalculateVZEROEventPlane(fESD, 1, 2, qx, qy);
485 fEPV0CR2 = ep->CalculateVZEROEventPlane(fESD, 2, 2, qx, qy);
486 fEPV0CR3 = ep->CalculateVZEROEventPlane(fESD, 3, 2, qx, qy);
3c40321c 487
70d53162 488 }
3c40321c 489//cout<<" fEPV0:"<<fEPV0<<" fEPV0A:"<<fEPV0A<<" fEPV0C:"<<fEPV0C<<" fEPV0Ar:"<<fEPV0Ar<<" fEPV0Cr:"<<fEPV0Cr<<" fEPV0r:"<<fEPV0AR4<<" fEPV0AR7:"<<fEPV0AR7<<" fEPV0CR0:"<<fEPV0CR0<<" fEPV0CR3:"<<fEPV0CR3<<"--------------------------------------------"<<endl;
70d53162 490 if(fcheckEP2sub){
491 if(fEPV0r<-2. || fEPV0Ar<-2. || fEPV0Cr<-2.) return;
492 }
493 if(fcheckEP2sub){
494 if(fEPV0A<-2. || fEPV0C<-2. || fEPV0AR4<-2. || fEPV0AR7<-2. || fEPV0CR0<-2. || fEPV0CR3<-2.) return;
495 }
3c40321c 496
70d53162 497 if(fEPV0 < TMath::Pi()) fEPV0 += TMath::Pi();
498 if(fEPV0A < TMath::Pi()) fEPV0A += TMath::Pi();
499 if(fEPV0C < TMath::Pi()) fEPV0C += TMath::Pi();
500 if(fEPV0Ar < TMath::Pi()) fEPV0A += TMath::Pi();
501 if(fEPV0Cr < TMath::Pi()) fEPV0C += TMath::Pi();
502 if(fEPV0r < TMath::Pi()) fEPV0A += TMath::Pi();
503 if(fEPV0AR4 < TMath::Pi()) fEPV0AR4 += TMath::Pi();
504 if(fEPV0AR7 < TMath::Pi()) fEPV0AR7 += TMath::Pi();
505 if(fEPV0CR0 < TMath::Pi()) fEPV0CR0 += TMath::Pi();
506 if(fEPV0CR3 < TMath::Pi()) fEPV0CR3 += TMath::Pi();
3c40321c 507
508//cout<<" fEPV0:"<<fEPV0<<" fEPV0A:"<<fEPV0A<<" fEPV0C:"<<fEPV0C<<" fEPV0Ar:"<<fEPV0Ar<<" fEPV0Cr:"<<fEPV0Cr<<" fEPV0r:"<<fEPV0AR4<<" fEPV0AR7:"<<fEPV0AR7<<" fEPV0CR0:"<<fEPV0CR0<<" fEPV0CR3:"<<fEPV0CR3<<"--------------------------------------------"<<endl;
70d53162 509 hEPTPC->Fill(fCentrality, fEPTPC);
510 if(fEPTPCreso!=-1) hresoTPC->Fill(fCentrality, fEPTPCreso);
511 hEPV0->Fill(fCentrality, fEPV0);
512 hEPV0A->Fill(fCentrality, fEPV0A);
513 hEPV0C->Fill(fCentrality, fEPV0C);
514 hEPV0Ar->Fill(fCentrality, fEPV0Ar);
515 hEPV0Cr->Fill(fCentrality, fEPV0Cr);
516 hEPV0r->Fill(fCentrality, fEPV0r);
517
518 hdifV0A_V0CR0->Fill(fCentrality, TMath::Cos(2.*(fEPV0A - fEPV0CR0)));
519 hdifV0A_V0CR3->Fill(fCentrality, TMath::Cos(2.*(fEPV0A - fEPV0CR3)));
520 hdifV0ACR0_V0CR3->Fill(fCentrality, TMath::Cos(2*(fEPV0CR0 - fEPV0CR3)));
521 hdifV0C_V0AR4->Fill(fCentrality, TMath::Cos(2*(fEPV0C - fEPV0AR4)));
522 hdifV0C_V0AR7->Fill(fCentrality, TMath::Cos(2*(fEPV0C - fEPV0AR7)));
523 hdifV0AR4_V0AR7->Fill(fCentrality, TMath::Cos(2*(fEPV0AR4 - fEPV0AR7)));
3c40321c 524
70d53162 525 // Cluster loop for reconstructed event
3c40321c 526
70d53162 527 Int_t nCluster = fESD->GetNumberOfCaloClusters();
528 for(Int_t i=0; i<nCluster; ++i){
529 AliESDCaloCluster *c1 = fESD->GetCaloCluster(i);
530 if(!IsGoodCluster(c1)) continue;
531 for(Int_t j=i+1; j<nCluster; ++j){
532 AliESDCaloCluster *c2 = fESD->GetCaloCluster(j);
533 if(!IsGoodCluster(c2)) continue;
534 TLorentzVector p1;
535 GetMom(p1, c1, vertex);
536 TLorentzVector p2;
537 GetMom(p2, c2, vertex);
538 FillPion(p1, p2, fEPV0r, fEPV0A, fEPV0C, fEPTPC);
539 }
540 }
3c40321c 541
70d53162 542 // NEW HISTO should be filled before this point, as PostData puts the
543 // information for this iteration of the UserExec in the container
544 PostData(1, fOutput);
3c40321c 545}
546
3c40321c 547//________________________________________________________________________
548void AliAnalysisTaskPi0V2::Terminate(Option_t *)
549{
70d53162 550 // Called once at the end of the query.
3c40321c 551}