]>
Commit | Line | Data |
---|---|---|
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 | ||
40 | ClassImp(AliAnalysisTaskEMCALPi0V2ShSh) | |
41 | ||
42 | //________________________________________________________________________ | |
43 | AliAnalysisTaskEMCALPi0V2ShSh::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 | //________________________________________________________________________ | |
62 | AliAnalysisTaskEMCALPi0V2ShSh::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 | //________________________________________________________________________ | |
87 | void 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 | //________________________________________________________________________ | |
287 | void 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 | //________________________________________________________________________ | |
334 | void 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 | //________________________________________________________________________ | |
394 | void 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 | //________________________________________________________________________ | |
472 | void 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 | //________________________________________________________________________ | |
490 | void 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 | } |