including Andre Stahl's task
[u/mrichter/AliRoot.git] / PWGGA / EMCALTasks / AliAnalysisTaskEMCALPi0V2ShSh.cxx
CommitLineData
fc0f5568 1// $Id: AliAnalysisTaskEMCALPi0V2ShSh.cxx$
2
3#include "AliAnalysisTaskEMCALPi0V2ShSh.h"
4
5//Root include files
6//#include <Riostream.h>
7#include <TParticle.h>
8#include <TRefArray.h>
9#include <TCanvas.h>
10#include <TChain.h>
11#include <TFile.h>
12#include <TH1F.h>
13#include <TH1D.h>
14#include <TH2F.h>
15#include <THnSparse.h>
16#include <TList.h>
17#include <TMath.h>
18#include <TVirtualFFT.h>
19
20//AliRoot include files
21#include "AliAnalysisTaskSE.h"
22#include "AliRunLoader.h"
23#include "AliAnalysisManager.h"
24#include "AliAnalysisTask.h"
25#include "AliStack.h"
26#include "AliESDEvent.h"
27#include "AliESDVertex.h"
28#include "AliESDCaloCluster.h"
29#include "AliESDEvent.h"
30#include "AliESDHeader.h"
31#include "AliESDInputHandler.h"
32#include "AliESDtrack.h"
33#include "AliKFParticle.h"
34#include "AliVCluster.h"
35#include "AliCentrality.h"
36#include "AliEventplane.h"
37
38
39
40ClassImp(AliAnalysisTaskEMCALPi0V2ShSh)
41
42//________________________________________________________________________
43AliAnalysisTaskEMCALPi0V2ShSh::AliAnalysisTaskEMCALPi0V2ShSh() :
44 AliAnalysisTaskSE(),
45 fCaloClusters(0), fEventPlane(0), fCentralityV0M(99.), fESD(0), fOutputList(0),
46 fEPTPC(-999.), fEPTPCResolution(0.),
47 fEPV0(-999.), fEPV0A(-999.), fEPV0C(-999.), fEPV0Ar(-999.), fEPV0Cr(-999.), fEPV0r(-999.),
48 fEPV0A4r(-999.), fEPV0A5r(-999.), fEPV0A6r(-999.), fEPV0A7r(-999.), fEPV0C0r(-999.), fEPV0C1r(-999.), fEPV0C2r(-999.), fEPV0C3r(-999.),
49 fHistAllcentV0(0), fHistAllcentV0r(0), fHistAllcentV0A(0), fHistAllcentV0C(0), fHistAllcentTPC(0),
50 fHistEPTPC(0), fHistEPTPCResolution(0),
51 fHistEPV0(0), fHistEPV0A(0), fHistEPV0C(0), fHistEPV0Ar(0), fHistEPV0Cr(0), fHistEPV0r(0), fHistEPV0A4r(0), fHistEPV0A7r(0), fHistEPV0C0r(0), fHistEPV0C3r(0),
52 fHistdifV0A_V0C0r(0), fHistdifV0A_V0C3r(0), fHistdifV0C0r_V0C3r(0), fHistdifV0C_V0A4r(0), fHistdifV0C_V0A7r(0), fHistdifV0A4r_V0A7r(0), fHistdifV0Ar_V0Cr(0),
53 fHistClusterEta(0), fHistClusterPhi(0), fHistClusterE(0), fHistClusterEt(0), fHistClusterN(0), fHistClusterM02(0),
54 fHistClusterEN(0), fHistClusterEM02(0), fHistClusterPhiEta(0), fHistClusterEtN(0), fHistClusterEtM02(0), fHistClusterdphiV0(0),
55 fHistTrackPt(0), fHistTrackEta(0), fHistTrackPhi(0), fHistTrackPhiEta(0),
56 fClusterPbV0(0), fClusterPbV0A(0), fClusterPbV0C(0), fClusterPbTPC(0)
57{
58 // Default constructor.
59}
60
61//________________________________________________________________________
62AliAnalysisTaskEMCALPi0V2ShSh::AliAnalysisTaskEMCALPi0V2ShSh(const char *name) :
63 AliAnalysisTaskSE(name),
64 fCaloClusters(0), fEventPlane(0), fCentralityV0M(99.), fESD(0), fOutputList(0),
65 fEPTPC(-999.), fEPTPCResolution(0.),
66 fEPV0(-999.), fEPV0A(-999.), fEPV0C(-999.), fEPV0Ar(-999.), fEPV0Cr(-999.), fEPV0r(-999.),
67 fEPV0A4r(-999.), fEPV0A5r(-999.), fEPV0A6r(-999.), fEPV0A7r(-999.), fEPV0C0r(-999.), fEPV0C1r(-999.), fEPV0C2r(-999.), fEPV0C3r(-999.),
68 fHistAllcentV0(0), fHistAllcentV0r(0), fHistAllcentV0A(0), fHistAllcentV0C(0), fHistAllcentTPC(0),
69 fHistEPTPC(0), fHistEPTPCResolution(0),
70 fHistEPV0(0), fHistEPV0A(0), fHistEPV0C(0), fHistEPV0Ar(0), fHistEPV0Cr(0), fHistEPV0r(0), fHistEPV0A4r(0), fHistEPV0A7r(0), fHistEPV0C0r(0), fHistEPV0C3r(0),
71 fHistdifV0A_V0C0r(0), fHistdifV0A_V0C3r(0), fHistdifV0C0r_V0C3r(0), fHistdifV0C_V0A4r(0), fHistdifV0C_V0A7r(0), fHistdifV0A4r_V0A7r(0), fHistdifV0Ar_V0Cr(0),
72 fHistClusterEta(0), fHistClusterPhi(0), fHistClusterE(0), fHistClusterEt(0), fHistClusterN(0), fHistClusterM02(0),
73 fHistClusterEN(0), fHistClusterEM02(0), fHistClusterPhiEta(0), fHistClusterEtN(0), fHistClusterEtM02(0), fHistClusterdphiV0(0),
74 fHistTrackPt(0), fHistTrackEta(0), fHistTrackPhi(0), fHistTrackPhiEta(0),
75 fClusterPbV0(0), fClusterPbV0A(0), fClusterPbV0C(0), fClusterPbTPC(0)
76{
77 // Constructor
78 // Define input and output slots here
79 // Input slot #0 works with a TChain
80 DefineInput(0, TChain::Class());
81 // Output slot #0 id reserved by the base class for AOD
82 // Output slot #1 writes into a TH1 container
83 DefineOutput(1, TList::Class());
84}
85
86//________________________________________________________________________
87void AliAnalysisTaskEMCALPi0V2ShSh::UserCreateOutputObjects()
88{
89 // Create histograms, called once.
90
91 fCaloClusters = new TRefArray();
92
93 fOutputList = new TList();
94 fOutputList->SetOwner();// Container cleans up all histos (avoids leaks in merging)
95
96 fHistTrackPt = new TH1F("fHistTrackPt","Track Transverse Momentum Distribution of Pb+Pb",100,0.0,30.0);
97 fHistTrackPt->GetYaxis()->SetTitle("Entries"); fHistTrackPt->GetXaxis()->SetTitle("P_{t} [GeV/c]");
98 fOutputList->Add(fHistTrackPt);
99
100 fHistTrackEta = new TH1F("fHistTrackEta","Track Pseudorapidity Distribution of Pb+Pb",100,-1.0,1.0);
101 fHistTrackEta->GetYaxis()->SetTitle("Entries"); fHistTrackEta->GetXaxis()->SetTitle("#eta");
102 fOutputList->Add(fHistTrackEta);
103
104 fHistTrackPhi = new TH1F("fHistTrackPhi","Track #phi Distribution of Pb+Pb",100,0.0,6.29);
105 fHistTrackPhi->GetYaxis()->SetTitle("Entries"); fHistTrackPhi->GetXaxis()->SetTitle("#phi [rad]");
106 fOutputList->Add(fHistTrackPhi);
107
108 fHistTrackPhiEta = new TH2F("fHistTrackPhiEta","Track Pseudorapidity vs #phi of Pb+Pb",100,-1.0,1.0,100,0.0,6.29);
109 fHistTrackPhiEta->GetXaxis()->SetTitle("#eta"); fHistTrackPhiEta->GetYaxis()->SetTitle("#phi [rad]");
110 fOutputList->Add(fHistTrackPhiEta);
111
112 fHistClusterEta = new TH1F("fHistClusterEta","Cluster Pseudorapidity Distribution of Pb+Pb",100,-1.0,1.0);
113 fHistClusterEta->GetYaxis()->SetTitle("Entries"); fHistClusterEta->GetXaxis()->SetTitle("#eta");
114 fOutputList->Add(fHistClusterEta);
115
116 fHistClusterPhi = new TH1F("fHistClusterPhi","Cluster #phi Distribution of Pb+Pb",100,0.0,6.29);
117 fHistClusterPhi->GetYaxis()->SetTitle("Entries"); fHistClusterPhi->GetXaxis()->SetTitle("#phi [rad]");
118 fOutputList->Add(fHistClusterPhi);
119
120 fHistClusterPhiEta = new TH2F("fHistClusterPhiEta","Cluster Pseudorapidity vs #phi of Pb+Pb",100,-1.0,1.0,100,0.0,6.29);
121 fHistClusterPhiEta->GetXaxis()->SetTitle("#eta"); fHistClusterPhiEta->GetYaxis()->SetTitle("#phi [rad]");
122 fOutputList->Add(fHistClusterPhiEta);
123
124 fHistClusterM02 = new TH1F("fHistClusterM02","Cluster M02 Distribution of Pb+Pb",100,0.0,3.0);
125 fHistClusterM02->GetYaxis()->SetTitle("Entries"); fHistClusterM02->GetXaxis()->SetTitle("M02");
126 fOutputList->Add(fHistClusterM02);
127
128 fHistClusterE = new TH1F("fHistClusterE","Cluster Energy Distribution of Pb+Pb",100,0.0,20.0);
129 fHistClusterE->GetYaxis()->SetTitle("Entries"); fHistClusterE->GetXaxis()->SetTitle("Energy [GeV]");
130 fOutputList->Add(fHistClusterE);
131
132 fHistClusterEt = new TH1F("fHistClusterEt","Cluster Transverse Energy Distribution of Pb+Pb",100,0.0,20.0);
133 fHistClusterEt->GetYaxis()->SetTitle("Entries"); fHistClusterEt->GetXaxis()->SetTitle("Transverse Energy [GeV]");
134 fOutputList->Add(fHistClusterEt);
135
136 fHistClusterEM02 = new TH2F("fHistClusterEM02","Cluster Energy vs M02 of Pb+Pb",100,0.0,20.0,100,0.0,3.0);
137 fHistClusterEM02->GetYaxis()->SetTitle("M02"); fHistClusterEM02->GetXaxis()->SetTitle("Energy [GeV]");
138 fOutputList->Add(fHistClusterEM02);
139
140 fHistClusterEtM02 = new TH2F("fHistClusterEtM02","Cluster Transverse Energy vs M02 of Pb+Pb",100,0.0,20.0,100,0.0,3.0);
141 fHistClusterEtM02->GetYaxis()->SetTitle("M02"); fHistClusterEtM02->GetXaxis()->SetTitle("Transverse Energy [GeV]");
142 fOutputList->Add(fHistClusterEtM02);
143
144 fHistClusterN = new TH1F("fHistClusterN","Cluster N Distribution of Pb+Pb",30,0.0,30.0);
145 fHistClusterN->GetYaxis()->SetTitle("Entries"); fHistClusterN->GetXaxis()->SetTitle("N");
146 fOutputList->Add(fHistClusterN);
147
148 fHistClusterEN = new TH2F("fHistClusterEN","N vs Cluster Energy of Pb+Pb",100,0.0,20.0,30,0.0,30.0);
149 fHistClusterEN->GetYaxis()->SetTitle("N"); fHistClusterEN->GetXaxis()->SetTitle("Energy [GeV]");
150 fOutputList->Add(fHistClusterEN);
151
152 fHistClusterEtN = new TH2F("fHistClusterEtN","N vs Cluster Transverse Energy of Pb+Pb",100,0.0,20.0,30,0.0,30.0);
153 fHistClusterEtN->GetYaxis()->SetTitle("N"); fHistClusterEtN->GetXaxis()->SetTitle("Transverse Energy [GeV]");
154 fOutputList->Add(fHistClusterEtN);
155
156 fHistClusterdphiV0 = new TH1D("fHistClusterdphiV0","Cluster dphiV0 Distribution of Pb+Pb",100,0.0,TMath::Pi());
157 fHistClusterdphiV0->GetYaxis()->SetTitle("Entries"); fHistClusterdphiV0->GetXaxis()->SetTitle("dphiV0 [rad]");
158 fOutputList->Add(fHistClusterdphiV0);
159
160 fHistEPV0 = new TH2F("fHistEPV0","V0 Event Plane Angle vs Centrality of Pb+Pb",100,0,100,100,0.0,TMath::Pi());
161 fHistEPV0->GetYaxis()->SetTitle("V0 Event Plane Angle [rad]"); fHistEPV0->GetXaxis()->SetTitle("V0M Centrality");
162 fOutputList->Add(fHistEPV0);
163
164 fHistEPV0A = new TH2F("fHistEPV0A","V0A Event Plane Angle vs Centrality of Pb+Pb",100,0,100,100,0.0,TMath::Pi());
165 fHistEPV0A->GetYaxis()->SetTitle("V0A Event Plane Angle [rad]"); fHistEPV0A->GetXaxis()->SetTitle("V0M Centrality");
166 fOutputList->Add(fHistEPV0A);
167
168 fHistEPV0C = new TH2F("fHistEPV0","V0C Event Plane Angle vs Centrality of Pb+Pb",100,0,100,100,0.0,TMath::Pi());
169 fHistEPV0C->GetYaxis()->SetTitle("V0C Event Plane Angle [rad]"); fHistEPV0C->GetXaxis()->SetTitle("V0M Centrality");
170 fOutputList->Add(fHistEPV0C);
171
172 fHistEPV0r = new TH2F("fHistEPV0r","V0r Event Plane Angle vs Centrality of Pb+Pb",100,0,100,100,0.0,TMath::Pi());
173 fHistEPV0r->GetYaxis()->SetTitle("V0r Event Plane Angle [rad]"); fHistEPV0r->GetXaxis()->SetTitle("V0M Centrality");
174 fOutputList->Add(fHistEPV0r);
175
176 fHistEPV0Ar = new TH2F("fHistEPV0","V0Ar Event Plane Angle vs Centrality of Pb+Pb",100,0,100,100,0.0,TMath::Pi());
177 fHistEPV0Ar->GetYaxis()->SetTitle("V0Ar Event Plane Angle [rad]"); fHistEPV0Ar->GetXaxis()->SetTitle("V0M Centrality");
178 fOutputList->Add(fHistEPV0Ar);
179
180 fHistEPV0Cr = new TH2F("fHistEPV0Cr","V0Cr Event Plane Angle vs Centrality of Pb+Pb",100,0,100,100,0.0,TMath::Pi());
181 fHistEPV0Cr->GetYaxis()->SetTitle("V0Cr Event Plane Angle [rad]"); fHistEPV0Cr->GetXaxis()->SetTitle("V0M Centrality");
182 fOutputList->Add(fHistEPV0Cr);
183
184 fHistEPV0A4r = new TH2F("fHistEPV0Cr","V0A4r Event Plane Angle vs Centrality of Pb+Pb",100,0,100,100,0.0,TMath::Pi());
185 fHistEPV0A4r->GetYaxis()->SetTitle("V0A4r Event Plane Angle [rad]"); fHistEPV0A4r->GetXaxis()->SetTitle("V0M Centrality");
186 fOutputList->Add(fHistEPV0A4r);
187
188 fHistEPV0A7r = new TH2F("fHistEPV0Cr","V0A7r Event Plane Angle vs Centrality of Pb+Pb",100,0,100,100,0.0,TMath::Pi());
189 fHistEPV0A7r->GetYaxis()->SetTitle("V0A7r Event Plane Angle [rad]"); fHistEPV0A7r->GetXaxis()->SetTitle("V0M Centrality");
190 fOutputList->Add(fHistEPV0A7r);
191
192 fHistEPV0C0r = new TH2F("fHistEPV0Cr","V0C0r Event Plane Angle vs Centrality of Pb+Pb",100,0,100,100,0.0,TMath::Pi());
193 fHistEPV0C0r->GetYaxis()->SetTitle("V0C0r Event Plane Angle [rad]"); fHistEPV0C0r->GetXaxis()->SetTitle("V0M Centrality");
194 fOutputList->Add(fHistEPV0C0r);
195
196 fHistEPV0C3r = new TH2F("fHistEPV0Cr","V0C3r Event Plane Angle vs Centrality of Pb+Pb",100,0,100,100,0.0,TMath::Pi());
197 fHistEPV0C3r->GetYaxis()->SetTitle("V0C3r Event Plane Angle [rad]"); fHistEPV0C3r->GetXaxis()->SetTitle("V0M Centrality");
198 fOutputList->Add(fHistEPV0C3r);
199
200 fHistdifV0A_V0C0r = new TH2F("fHistdifV0A_V0C0r","(V0A - V0C0r) vs Centrality of Pb+Pb",100,0,100,100,-1.0,1.0);
201 fHistdifV0A_V0C0r->GetYaxis()->SetTitle("Cos[2*(V0A - V0C0r)]"); fHistdifV0A_V0C0r->GetXaxis()->SetTitle("V0M Centrality");
202 fOutputList->Add(fHistdifV0A_V0C0r);
203
204 fHistdifV0A_V0C3r = new TH2F("fHistdifV0A_V0C3r","(V0A - V0C3r) vs Centrality of Pb+Pb",100,0,100,100,-1.0,1.0);
205 fHistdifV0A_V0C3r->GetYaxis()->SetTitle("Cos[2*(V0A - V0C3r)]"); fHistdifV0A_V0C3r->GetXaxis()->SetTitle("V0M Centrality");
206 fOutputList->Add(fHistdifV0A_V0C3r);
207
208 fHistdifV0C0r_V0C3r = new TH2F("fHistdifV0C0r_V0C3r","(V0C0r - V0C3r) vs Centrality of Pb+Pb",100,0,100,100,-1.0,1.0);
209 fHistdifV0C0r_V0C3r->GetYaxis()->SetTitle("Cos[2*(V0C0r - V0C3r)]"); fHistdifV0C0r_V0C3r->GetXaxis()->SetTitle("V0M Centrality");
210 fOutputList->Add(fHistdifV0C0r_V0C3r);
211
212 fHistdifV0C_V0A4r = new TH2F("fHistdifV0C_V0A4r","(V0C - V0A4r) vs Centrality of Pb+Pb",100,0,100,100,-1.0,1.0);
213 fHistdifV0C_V0A4r->GetYaxis()->SetTitle("Cos[2*(V0C - V0C4r)]"); fHistdifV0C_V0A4r->GetXaxis()->SetTitle("V0M Centrality");
214 fOutputList->Add(fHistdifV0C_V0A4r);
215
216 fHistdifV0C_V0A7r = new TH2F("fHistdifV0C_V0A3r","(V0C - V0A7r) vs Centrality of Pb+Pb",100,0,100,100,-1.0,1.0);
217 fHistdifV0C_V0A7r->GetYaxis()->SetTitle("Cos[2*(V0C - V0A7r)]"); fHistdifV0C_V0A7r->GetXaxis()->SetTitle("V0M Centrality");
218 fOutputList->Add(fHistdifV0C_V0A7r);
219
220 fHistdifV0A4r_V0A7r = new TH2F("fHistdifV0A4r_V0A7r","(V0A4r - V0A7r) vs Centrality of Pb+Pb",100,0,100,100,-1.0,1.0);
221 fHistdifV0A4r_V0A7r->GetYaxis()->SetTitle("Cos[2*(V0A4r - V0A7r)]"); fHistdifV0A4r_V0A7r->GetXaxis()->SetTitle("V0M Centrality");
222 fOutputList->Add(fHistdifV0A4r_V0A7r);
223
224 fHistdifV0Ar_V0Cr = new TH2F("fHistdifV0Ar_V0Cr","(V0Ar - V0Cr) vs Centrality of Pb+Pb",100,0,100,100,-1.0,1.0);
225 fHistdifV0Ar_V0Cr->GetYaxis()->SetTitle("Cos[2*(V0Ar - V0Cr)]"); fHistdifV0Ar_V0Cr->GetXaxis()->SetTitle("V0M Centrality");
226 fOutputList->Add(fHistdifV0Ar_V0Cr);
227
228 fHistAllcentV0 = new TH1F("fHistAllcentV0","V0 Event Plane Angle of Pb+Pb",100,0.0,TMath::Pi());
229 fHistAllcentV0->GetXaxis()->SetTitle("V0 Event Plane Angle [rad]"); fHistAllcentV0->GetYaxis()->SetTitle("");
230 fOutputList->Add(fHistAllcentV0);
231
232 fHistAllcentV0r = new TH1F("fHistAllcentV0r","V0r Event Plane Angle of Pb+Pb",100,0.0,TMath::Pi());
233 fHistAllcentV0r->GetXaxis()->SetTitle("V0r Event Plane Angle [rad]"); fHistAllcentV0r->GetYaxis()->SetTitle("");
234 fOutputList->Add(fHistAllcentV0r);
235
236 fHistAllcentV0A = new TH1F("fHistAllcentV0A","V0A Event Plane Angle of Pb+Pb",100,0.0,TMath::Pi());
237 fHistAllcentV0A->GetXaxis()->SetTitle("V0A Event Plane Angle [rad]"); fHistAllcentV0A->GetYaxis()->SetTitle("");
238 fOutputList->Add(fHistAllcentV0A);
239
240 fHistAllcentV0C = new TH1F("fHistAllcentV0C","V0C Event Plane Angle of Pb+Pb",100,0.0,TMath::Pi());
241 fHistAllcentV0C->GetXaxis()->SetTitle("V0C Event Plane Angle [rad]"); fHistAllcentV0C->GetYaxis()->SetTitle("");
242 fOutputList->Add(fHistAllcentV0C);
243
244 fHistAllcentTPC = new TH1F("fHistAllcentTPC","TPC Event Plane Angle of Pb+Pb",100,0.0,TMath::Pi());
245 fHistAllcentTPC->GetXaxis()->SetTitle("TPC Event Plane Angle [rad]"); fHistAllcentTPC->GetYaxis()->SetTitle("");
246 fOutputList->Add(fHistAllcentTPC);
247
248 fHistEPTPC = new TH2F("fHistEPTPC","TPC Event Plane Angle vs Centrality of Pb+Pb",100,0,100,100,0.0,TMath::Pi());
249 fHistEPTPC->GetYaxis()->SetTitle("Event Plane Angle [rad]"); fHistEPTPC->GetXaxis()->SetTitle("V0M Centrality");
250 fOutputList->Add(fHistEPTPC);
251
252 fHistEPTPCResolution = new TH2F("fHistEPTPCResolution","TPC Resolution vs Centrality of Pb+Pb",100,0,100,100,0.0,1.0);
253 fHistEPTPCResolution->GetYaxis()->SetTitle("TPC Resolution"); fHistEPTPCResolution->GetXaxis()->SetTitle("V0M Centrality");
254 fOutputList->Add(fHistEPTPCResolution);
255
256
257 // Et M02 V0Mcent DeltaPhi Cos[2*DeltaPhi]
258 Int_t bins[5] = { 500, 250, 100, 100, 100 }; // binning
259 Double_t min[5] = { 0.0, 0.0, 0, 0.0, -1.0}; // min x
260 Double_t max[5] = { 50.0, 2.5, 100, TMath::Pi(), 1.0}; // max x
261
262 fClusterPbV0 = new THnSparseF("fClusterPbV0","",5,bins,min,max);
263 fClusterPbV0->GetAxis(0)->SetTitle("Transverse Energy [GeV]"); fClusterPbV0->GetAxis(1)->SetTitle("M02"); fClusterPbV0->GetAxis(2)->SetTitle("V0M Centrality");
264 fClusterPbV0->GetAxis(3)->SetTitle("Delta(#phi) [rad]"); fClusterPbV0->GetAxis(4)->SetTitle("Cos[2*Delta(#phi)]");
265 fOutputList->Add(fClusterPbV0);
266
267 fClusterPbV0A = new THnSparseF("fClusterPbV0A","",5,bins,min,max);
268 fClusterPbV0A->GetAxis(0)->SetTitle("Transverse Energy [GeV]"); fClusterPbV0A->GetAxis(1)->SetTitle("M02"); fClusterPbV0A->GetAxis(2)->SetTitle("V0M Centrality");
269 fClusterPbV0A->GetAxis(3)->SetTitle("Delta(#phi) [rad]"); fClusterPbV0A->GetAxis(4)->SetTitle("Cos[2*Delta(#phi)]");
270 fOutputList->Add(fClusterPbV0A);
271
272 fClusterPbV0C = new THnSparseF("fClusterPbV0C","",5,bins,min,max);
273 fClusterPbV0C->GetAxis(0)->SetTitle("Transverse Energy [GeV]"); fClusterPbV0C->GetAxis(1)->SetTitle("M02"); fClusterPbV0C->GetAxis(2)->SetTitle("V0M Centrality");
274 fClusterPbV0C->GetAxis(3)->SetTitle("Delta(#phi) [rad]"); fClusterPbV0C->GetAxis(4)->SetTitle("Cos[2*Delta(#phi)]");
275 fOutputList->Add(fClusterPbV0C);
276
277 fClusterPbTPC = new THnSparseF("fClusterPbTPC","",5,bins,min,max);
278 fClusterPbTPC->GetAxis(0)->SetTitle("Transverse Energy [GeV]"); fClusterPbTPC->GetAxis(1)->SetTitle("M02"); fClusterPbTPC->GetAxis(2)->SetTitle("V0M Centrality");
279 fClusterPbTPC->GetAxis(3)->SetTitle("Delta(#phi) [rad]"); fClusterPbTPC->GetAxis(4)->SetTitle("Cos[2*Delta(#phi)]");
280 fOutputList->Add(fClusterPbTPC);
281
282
283 PostData(1, fOutputList);
284}
285
286//________________________________________________________________________
287void AliAnalysisTaskEMCALPi0V2ShSh::UserExec(Option_t *)
288{
289 // Main loop, called for each event.
290
291 // Create pointer to reconstructed event
292 AliVEvent *event = InputEvent();
293 if (!event) { Printf("ERROR: Could not retrieve event\n"); return; }
294
295 fESD = dynamic_cast<AliESDEvent*>(event);
296 if (!fESD) {
297 printf("ERROR: Could not retrieve the ESD event\n");
298 return;
299 }
300
301 //Bool_t isSelected =0;
302 //isSelected = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & (AliVEvent::kSemiCentral));
303 //if(!isSelected) { return; }
304
305 if(fESD->GetCentrality()){
306 fCentralityV0M = fESD->GetCentrality()->GetCentralityPercentile("V0M");
307 } else {
308 printf("ERROR: Could not retrieve the Centrality\n");
309 return;
310 }
311
312 fEventPlane = fESD->GetEventplane();
313 if (fEventPlane) {
314 VZEROEventPlane();
315 } else {
316 printf("ERROR: Could not retrieve the Centrality\n");
317 return;
318 }
319
320 if(fESD->GetEMCALClusters(fCaloClusters)) { //select EMCAL clusters only
321 FillClusterHists();
322 } else{
323 printf("ERROR: Could not get EMCAL Cluster\n");
324 return;
325 }
326
327 FillTrackHists();
328
329 fCaloClusters->Clear();
330 PostData(1, fOutputList);
331}
332
333//________________________________________________________________________
334void AliAnalysisTaskEMCALPi0V2ShSh::VZEROEventPlane()
335{// Calculate the V0 Event Plane
336
337 if (fEventPlane->GetQVector()) {
338 fEPTPC = TVector2::Phi_0_2pi(fEventPlane->GetQVector()->Phi())/2.0; //if(fEPTPC>TMath::Pi()) {fEPTPC-=TMath::Pi();}
339 } else { fEPTPC = -999.; }
340
341 if (fEventPlane->GetQsub1()&&fEventPlane->GetQsub2()) {
342 fEPTPCResolution = TMath::Cos(2.0*(fEventPlane->GetQsub1()->Phi()/2.0-fEventPlane->GetQsub2()->Phi()/2.0)); }
343 else { fEPTPCResolution = -1; }
344
345 fEPV0 = TVector2::Phi_0_2pi(fEventPlane->GetEventplane("V0", fESD)); if(fEPV0>TMath::Pi()) {fEPV0-=TMath::Pi();}
346 fEPV0A = TVector2::Phi_0_2pi(fEventPlane->GetEventplane("V0A", fESD)); if(fEPV0A>TMath::Pi()) {fEPV0A-=TMath::Pi();}
347 fEPV0C = TVector2::Phi_0_2pi(fEventPlane->GetEventplane("V0C", fESD)); if(fEPV0C>TMath::Pi()) {fEPV0C-=TMath::Pi();}
348
349 Double_t qx=0, qy=0, qxr=0, qyr=0;
350 fEPV0Ar = TVector2::Phi_0_2pi(fEventPlane->CalculateVZEROEventPlane(fESD, 4, 5, 2, qxr, qyr)); if(fEPV0Ar>TMath::Pi()) {fEPV0Ar-=TMath::Pi();}
351 fEPV0Cr = TVector2::Phi_0_2pi(fEventPlane->CalculateVZEROEventPlane(fESD, 2, 3, 2, qx, qy)); if(fEPV0Cr>TMath::Pi()) {fEPV0Cr-=TMath::Pi();}
352 qxr += qx; qyr += qy;
353 fEPV0r = TVector2::Phi_0_2pi(TMath::ATan2(qyr,qxr))/2.0;
354
355 fEPV0A4r = TVector2::Phi_0_2pi(fEventPlane->CalculateVZEROEventPlane(fESD, 4, 2, qx, qy)); if(fEPV0A4r>TMath::Pi()) {fEPV0A4r-=TMath::Pi();}
356 fEPV0A5r = TVector2::Phi_0_2pi(fEventPlane->CalculateVZEROEventPlane(fESD, 5, 2, qx, qy)); if(fEPV0A5r>TMath::Pi()) {fEPV0A5r-=TMath::Pi();}
357 fEPV0A6r = TVector2::Phi_0_2pi(fEventPlane->CalculateVZEROEventPlane(fESD, 6, 2, qx, qy)); if(fEPV0A6r>TMath::Pi()) {fEPV0A6r-=TMath::Pi();}
358 fEPV0A7r = TVector2::Phi_0_2pi(fEventPlane->CalculateVZEROEventPlane(fESD, 7, 2, qx, qy)); if(fEPV0A7r>TMath::Pi()) {fEPV0A7r-=TMath::Pi();}
359 fEPV0C0r = TVector2::Phi_0_2pi(fEventPlane->CalculateVZEROEventPlane(fESD, 0, 2, qx, qy)); if(fEPV0C0r>TMath::Pi()) {fEPV0C0r-=TMath::Pi();}
360 fEPV0C1r = TVector2::Phi_0_2pi(fEventPlane->CalculateVZEROEventPlane(fESD, 1, 2, qx, qy)); if(fEPV0C1r>TMath::Pi()) {fEPV0C1r-=TMath::Pi();}
361 fEPV0C2r = TVector2::Phi_0_2pi(fEventPlane->CalculateVZEROEventPlane(fESD, 2, 2, qx, qy)); if(fEPV0C2r>TMath::Pi()) {fEPV0C2r-=TMath::Pi();}
362 fEPV0C3r = TVector2::Phi_0_2pi(fEventPlane->CalculateVZEROEventPlane(fESD, 3, 2, qx, qy)); if(fEPV0C3r>TMath::Pi()) {fEPV0C3r-=TMath::Pi();}
363
364 fHistEPTPC->Fill(fCentralityV0M, fEPTPC);
365 if(fEPTPCResolution!=-1) { fHistEPTPCResolution->Fill(fCentralityV0M, fEPTPCResolution); }
366 fHistEPV0->Fill(fCentralityV0M, fEPV0);
367 fHistEPV0A->Fill(fCentralityV0M, fEPV0A);
368 fHistEPV0C->Fill(fCentralityV0M, fEPV0C);
369 fHistEPV0Ar->Fill(fCentralityV0M, fEPV0Ar);
370 fHistEPV0Cr->Fill(fCentralityV0M, fEPV0Cr);
371 fHistEPV0r->Fill(fCentralityV0M, fEPV0r);
372 fHistEPV0A4r->Fill(fCentralityV0M, fEPV0A4r);
373 fHistEPV0A7r->Fill(fCentralityV0M, fEPV0A7r);
374 fHistEPV0C0r->Fill(fCentralityV0M, fEPV0C0r);
375 fHistEPV0C3r->Fill(fCentralityV0M, fEPV0C3r);
376
377 fHistAllcentV0->Fill(fEPV0);
378 fHistAllcentV0r->Fill(fEPV0r);
379 fHistAllcentV0A->Fill(fEPV0A);
380 fHistAllcentV0C->Fill(fEPV0C);
381 fHistAllcentTPC->Fill(fEPTPC);
382
383 fHistdifV0A_V0C0r->Fill(fCentralityV0M, TMath::Cos(2.0*(fEPV0A - fEPV0C0r)));
384 fHistdifV0A_V0C3r->Fill(fCentralityV0M, TMath::Cos(2.0*(fEPV0A - fEPV0C3r)));
385 fHistdifV0C0r_V0C3r->Fill(fCentralityV0M, TMath::Cos(2.0*(fEPV0C0r - fEPV0C3r)));
386 fHistdifV0C_V0A4r->Fill(fCentralityV0M, TMath::Cos(2.0*(fEPV0C - fEPV0A4r)));
387 fHistdifV0C_V0A7r->Fill(fCentralityV0M, TMath::Cos(2.0*(fEPV0C - fEPV0A7r)));
388 fHistdifV0A4r_V0A7r->Fill(fCentralityV0M, TMath::Cos(2.0*(fEPV0A4r - fEPV0A7r)));
389 fHistdifV0Ar_V0Cr->Fill(fCentralityV0M, TMath::Cos(2.0*(fEPV0Ar - fEPV0Cr)));
390
391}
392
393//________________________________________________________________________
394void AliAnalysisTaskEMCALPi0V2ShSh::FillClusterHists()
395{// Fill cluster histograms.
396
397 Float_t pos[3] ;
398
399 if(!fCaloClusters)
400 return;
401
402 const Int_t nclus = fCaloClusters->GetEntries();
403
404 if(nclus==0)
405 return;
406
407 for(Int_t iclus=0;iclus<nclus;iclus++){
408 AliESDCaloCluster *clus = (AliESDCaloCluster *) fCaloClusters->At(iclus); // retrieve cluster from esd
409 if(!clus)
410 continue;
411 if(!clus->IsEMCAL()){
412 printf("ERROR: Current Cluster is not an EMCAL Cluster\n");
413 continue;
414 }
415 clus->GetPosition(pos);
416 TVector3 vpos(pos[0],pos[1],pos[2]);
417 Double_t Transverse_Energy = ((clus->E())/ (TMath::CosH(vpos.Eta())));
418 Double_t dphiV0 = TVector2::Phi_0_2pi(vpos.Phi()-fEPV0r); if(dphiV0>TMath::Pi()) {dphiV0-=TMath::Pi();}
419 Double_t dphiV0A = TVector2::Phi_0_2pi(vpos.Phi()-fEPV0A); if(dphiV0A>TMath::Pi()) {dphiV0A-=TMath::Pi();}
420 Double_t dphiV0C = TVector2::Phi_0_2pi(vpos.Phi()-fEPV0C); if(dphiV0C>TMath::Pi()) {dphiV0C-=TMath::Pi();}
421 Double_t dphiTPC = TVector2::Phi_0_2pi(vpos.Phi()-fEPTPC); if(dphiTPC>TMath::Pi()) {dphiTPC-=TMath::Pi();}
422
423 Double_t DataV0[5];
424 DataV0[0] = Transverse_Energy;
425 DataV0[1] = clus->GetM02();
426 DataV0[2] = fCentralityV0M;
427 DataV0[3] = dphiV0;
428 DataV0[4] = TMath::Cos(2.0*(dphiV0));
429 fClusterPbV0->Fill(DataV0);
430
431 Double_t DataV0A[5];
432 DataV0A[0] = Transverse_Energy;
433 DataV0A[1] = clus->GetM02();
434 DataV0A[2] = fCentralityV0M;
435 DataV0A[3] = dphiV0A;
436 DataV0A[4] = TMath::Cos(2.0*(dphiV0A));
437 fClusterPbV0A->Fill(DataV0A);
438
439 Double_t DataV0C[5];
440 DataV0C[0] = Transverse_Energy;
441 DataV0C[1] = clus->GetM02();
442 DataV0C[2] = fCentralityV0M;
443 DataV0C[3] = dphiV0C;
444 DataV0C[4] = TMath::Cos(2.0*(dphiV0C));
445 fClusterPbV0C->Fill(DataV0C);
446
447 Double_t DataTPC[5];
448 DataTPC[0] = Transverse_Energy;
449 DataTPC[1] = clus->GetM02();
450 DataTPC[2] = fCentralityV0M;
451 DataTPC[3] = dphiTPC;
452 DataTPC[4] = TMath::Cos(2.0*(dphiTPC));
453 fClusterPbTPC->Fill(DataTPC);
454
455 fHistClusterE->Fill(clus->E());
456 fHistClusterEt->Fill(Transverse_Energy);
457 fHistClusterM02->Fill(clus->GetM02());
458 fHistClusterN->Fill(clus->GetNCells());
459 fHistClusterPhi->Fill(vpos.Phi());
460 fHistClusterEta->Fill(vpos.Eta());
461 fHistClusterPhiEta->Fill(vpos.Eta(),vpos.Phi());
462 fHistClusterEN->Fill(clus->E(),clus->GetNCells());
463 fHistClusterEM02->Fill(clus->E(),clus->GetM02());
464 fHistClusterEtN->Fill(Transverse_Energy,clus->GetNCells());
465 fHistClusterEtM02->Fill(Transverse_Energy,clus->GetM02());
466 fHistClusterdphiV0->Fill(dphiV0);
467
468 }
469}
470
471//________________________________________________________________________
472void AliAnalysisTaskEMCALPi0V2ShSh::FillTrackHists()
473{// Fill track histograms.
474
475 for (Int_t iTracks = 0; iTracks < fESD->GetNumberOfTracks(); iTracks++) {
476 AliESDtrack* track = (AliESDtrack*)fESD->GetTrack(iTracks);
477 if (!track){
478 printf("ERROR: Could not retreive esdtrack\n");
479 continue;
480 }
481 fHistTrackPt->Fill(track->Pt());
482 fHistTrackEta->Fill(track->Eta());
483 fHistTrackPhi->Fill(track->Phi());
484 fHistTrackPhiEta->Fill(track->Eta(),track->Phi());
485
486 }
487}
488
489//________________________________________________________________________
490void AliAnalysisTaskEMCALPi0V2ShSh::Terminate(Option_t *)
491{
492 // Draw result to screen, or perform fitting, normalizations
493 // Called once at the end of the query
494
495 /* fOutputList = dynamic_cast<TList*> (GetOutputData(1));
496 if(!fOutputList) { Printf("ERROR: could not retrieve TList fOutputList"); return; }
497
498 fHistTrackPt = dynamic_cast<TH1F*> (fOutputList->FindObject("fHistTrackPt"));
499 if (!fHistTrackPt) { Printf("ERROR: could not retrieve fHistTrackPt"); return;}
500 fHistTrackEta = dynamic_cast<TH1F*> (fOutputList->FindObject("fHistTrackEta"));
501 if (!fHistTrackEta) { Printf("ERROR: could not retrieve fHistTrackEta"); return;}
502 fHistTrackPhi = dynamic_cast<TH1F*> (fOutputList->FindObject("fHistTrackPhi"));
503 if (!fHistTrackPhi) { Printf("ERROR: could not retrieve fHistTrackPhi"); return;}
504 fHistClusterEta = dynamic_cast<TH1F*> (fOutputList->FindObject("fHistClusterEta"));
505 if (!fHistClusterEta) { Printf("ERROR: could not retrieve fHistClusterEta"); return;}
506 fHistClusterPhi = dynamic_cast<TH1F*> (fOutputList->FindObject("fHistClusterPhi"));
507 if (!fHistClusterPhi) { Printf("ERROR: could not retrieve fHistClusterPhi"); return;}
508 fHistClusterE = dynamic_cast<TH1F*> (fOutputList->FindObject("fHistClusterE"));
509 if (!fHistClusterE) { Printf("ERROR: could not retrieve fHistClusterE"); return;}
510 fHistClusterN = dynamic_cast<TH1F*> (fOutputList->FindObject("fHistClusterN"));
511 if (!fHistClusterN) { Printf("ERROR: could not retrieve fHistClusterN"); return;}
512 fHistClusterM02 = dynamic_cast<TH1F*> (fOutputList->FindObject("fHistClusterM02"));
513 if (!fHistClusterM02) { Printf("ERROR: could not retrieve fHistClusterM02"); return;}
514
515 // Get the physics selection histograms with the selection statistics
516 //AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
517 //AliESDInputHandler *inputH = dynamic_cast<AliESDInputHandler*>(mgr->GetInputEventHandler());
518 //TH2F *histStat = (TH2F*)inputH->GetStatistics();
519
520
521 TCanvas *canvas1 = new TCanvas("canvas1","Track P_{T} & #eta & #phi",10,10,1020,510);
522 canvas1->Divide(3,1);
523 canvas1->cd(1)->SetLogy();
524 fHistTrackPt->DrawCopy("");
525 canvas1->cd(2);
526 fHistTrackEta->DrawCopy("");
527 canvas1->cd(3);
528 fHistTrackPhi->DrawCopy("");
529
530 TCanvas *canvas2 = new TCanvas("canvas2","Cluster #eta & #phi",10,10,1020,510);
531 canvas2->Divide(2,1);
532 canvas2->cd(1);
533 fHistClusterEta->DrawCopy("");
534 canvas2->cd(2);
535 fHistClusterPhi->DrawCopy("");
536
537 TCanvas *canvas3 = new TCanvas("canvas3","Cluster E & N & M02",10,10,1020,510);
538 canvas3->Divide(3,1);
539 canvas3->cd(1)->SetLogy();
540 fHistClusterE->DrawCopy("");
541 canvas3->cd(2)->SetLogy();
542 fHistClusterN->DrawCopy("");
543 canvas3->cd(3)->SetLogy();
544 fHistClusterM02->DrawCopy("");
545
546 TCanvas *canvas4 = new TCanvas("canvas4","Track #phi vs #eta & Cluster #phi vs #eta",10,10,1020,510);
547 canvas4->Divide(2,1);
548 canvas4->cd(1)->SetLogz();
549 fHistTrackPhiEta->DrawCopy("colorz");
550 canvas4->cd(2)->SetLogz();
551 fHistClusterPhiEta->DrawCopy("colorz");
552
553 TCanvas *canvas5 = new TCanvas("canvas5","Cluster E vs N & E vs M02",10,10,1020,510);
554 canvas5->Divide(2,1);
555 canvas5->cd(1)->SetLogz();
556 fHistClusterEN->DrawCopy("colorz");
557 canvas5->cd(2)->SetLogz();
558 fHistClusterEM02->DrawCopy("colorz");
559
560 TCanvas *canvas6 = new TCanvas("canvas6","Cluster Et & Et vs N & Et vs M02",10,10,1020,510);
561 canvas6->cd(1)->SetLogy();
562 fHistClusterEt->DrawCopy("");
563 canvas6->cd(2)->SetLogz();
564 fHistClusterEtN->DrawCopy("colorz");
565 canvas6->cd(3)->SetLogz();
566 fHistClusterEtM02->DrawCopy("colorz");
567
568 canvas1->SaveAs("lhc11h_2_Track_PT_Eta_Phi.jpg");
569 canvas2->SaveAs("lhc11h_2_Cluster_Eta_Phi.jpg");
570 canvas3->SaveAs("lhc11h_2_Cluster_E_N_M02.jpg");
571 canvas4->SaveAs("lhc11h_2_Track_PhivsEta_Cluster_PhivsEta.jpg");
572 canvas5->SaveAs("lhc11h_2_Cluster_EvsN_EvsM02.jpg");
573 canvas6->SaveAs("lhc11h_2_Cluster_Et_EtvsN_EtvsM02.jpg"); */
574}