]>
Commit | Line | Data |
---|---|---|
11476c41 | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: The ALICE Off-line Project. * | |
5 | * Contributors are mentioned in the code where appropriate. * | |
6 | * * | |
7 | * Permission to use, copy, modify and distribute this software and its * | |
8 | * documentation strictly for non-commercial purposes is hereby granted * | |
9 | * without fee, provided that the above copyright notice appears in all * | |
10 | * copies and that both the copyright notice and this permission notice * | |
11 | * appear in the supporting documentation. The authors make no claims * | |
12 | * about the suitability of this software for any purpose. It is * | |
13 | * provided "as is" without express or implied warranty. * | |
14 | **************************************************************************/ | |
15 | // | |
16 | // | |
17 | // Jet fragmentation transverse momentum (j_T) analysis task | |
18 | // | |
19 | // Author: T.Snellman | |
20 | ||
21 | #include <TClonesArray.h> | |
22 | #include <TH1F.h> | |
23 | #include <TH2F.h> | |
24 | #include <TH3F.h> | |
25 | #include <TList.h> | |
26 | #include <TProfile.h> | |
27 | #include <TLorentzVector.h> | |
28 | #include <TVector.h> | |
29 | #include <TGraphErrors.h> | |
30 | #include <TGrid.h> | |
31 | #include <TSystem.h> | |
32 | #include <TFile.h> | |
33 | ||
34 | #include "AliCentrality.h" | |
35 | ||
36 | ||
37 | ||
38 | #include "AliVCluster.h" | |
39 | #include "AliAODCaloCluster.h" | |
40 | #include "AliESDCaloCluster.h" | |
41 | #include "AliVTrack.h" | |
42 | #include "AliEmcalJet.h" | |
43 | #include "AliRhoParameter.h" | |
44 | #include "AliLog.h" | |
45 | #include "AliJetContainer.h" | |
46 | #include "AliParticleContainer.h" | |
47 | #include "AliClusterContainer.h" | |
48 | #include "AliPicoTrack.h" | |
49 | ||
50 | #include "AliAnalysisTaskJetJTJT.h" | |
51 | ||
52 | ||
53 | ClassImp(AliAnalysisTaskJetJTJT) | |
54 | ||
55 | //________________________________________________________________________ | |
56 | AliAnalysisTaskJetJTJT::AliAnalysisTaskJetJTJT() : | |
57 | AliAnalysisTaskEmcalJet("AliAnalysisTaskJetJTJT", kTRUE), | |
58 | fHistTracksPt(0), | |
59 | fHistTracksJt(0), | |
60 | fHistClustersPt(0), | |
61 | fHistLeadingJetPt(0), | |
62 | fHistJetsPt(0), | |
63 | fHistBackgroundDone(0), | |
64 | fHistJTPta(0), | |
65 | fHistLogJTPta(0), | |
66 | fHistJTPta_all(0), | |
67 | fHistJTBg(0), | |
68 | fHistLogJTBg(0), | |
69 | fHistBgMulti(0), | |
70 | fHistBgPt(0), | |
71 | fHistJetEta(0), | |
72 | fHistJetMulti(0), | |
73 | fHistJetTracksPt(0), | |
74 | fhTrackingEfficiency(0), | |
75 | fNpttBins(1), | |
76 | fNptaBins(1), | |
77 | fJetsCont(0), | |
78 | //fJetsConts(0), | |
79 | //nJetsConts(0), | |
80 | fTracksCont(0), | |
81 | fCaloClustersCont(0), | |
82 | fTracks(0), | |
83 | fTrackArrayName("nonejk"), | |
84 | runPeriod(""), | |
85 | fEfficiency(0), | |
86 | debug(0) | |
87 | ||
88 | ||
89 | { | |
90 | // Default constructor. | |
91 | ||
92 | fHistTracksPt = new TH1*[fNcentBins]; | |
93 | fHistTracksJt = new TH1*[fNcentBins]; | |
94 | fHistClustersPt = new TH1*[fNcentBins]; | |
95 | fHistLeadingJetPt = new TH1*[fNcentBins]; | |
96 | fHistJetsPt = new TH1**[fNcentBins]; | |
97 | fHistBackgroundDone = new TH1**[fNcentBins]; | |
98 | fHistJTPta = new TH1***[fNcentBins]; | |
99 | fHistLogJTPta = new TH1***[fNcentBins]; | |
100 | fHistJTPta_all = new TH1***[fNcentBins]; | |
101 | fHistJTBg = new TH1***[fNcentBins]; | |
102 | fHistLogJTBg = new TH1***[fNcentBins]; | |
103 | fHistBgMulti = new TH1**[fNcentBins]; | |
104 | fHistBgPt = new TH1**[fNcentBins]; | |
105 | fHistJetEta = new TH1**[fNcentBins]; | |
106 | fHistJetMulti = new TH1**[fNcentBins]; | |
107 | fHistJetTracksPt = new TH1**[fNcentBins]; | |
108 | fhTrackingEfficiency = new TProfile*[fNcentBins]; | |
109 | //CentBinBorders = new Double_t[10]; | |
110 | ||
111 | ||
112 | for (Int_t i = 0; i < fNcentBins; i++) { | |
113 | fHistJTPta[i] = 0; | |
114 | fHistLogJTPta[i] = 0; | |
115 | fHistJTPta_all[i] = 0; | |
116 | fHistJTBg[i] = 0; | |
117 | fHistLogJTBg[i] = 0; | |
118 | fHistBackgroundDone[i] = 0; | |
119 | fHistTracksPt[i] = 0; | |
120 | fHistTracksJt[i] = 0; | |
121 | fHistClustersPt[i] = 0; | |
122 | fHistLeadingJetPt[i] = 0; | |
123 | fHistJetsPt[i] = 0; | |
124 | fHistBgMulti[i] = 0; | |
125 | fHistBgPt[i] = 0; | |
126 | fHistJetEta[i] = 0; | |
127 | fHistJetMulti[i] = 0; | |
128 | fHistJetTracksPt[i] = 0; | |
129 | fhTrackingEfficiency[i] = 0; | |
130 | } | |
131 | ||
132 | /*for(Int_t i = 0; i < nJetsConts; i++){ | |
133 | fJetsConts[i] = 0; | |
134 | }*/ | |
135 | SetMakeGeneralHistograms(kTRUE); | |
136 | } | |
137 | ||
138 | //________________________________________________________________________ | |
139 | AliAnalysisTaskJetJTJT::AliAnalysisTaskJetJTJT(const char *name) : | |
140 | AliAnalysisTaskEmcalJet(name, kTRUE), | |
141 | fHistTracksPt(0), | |
142 | fHistTracksJt(0), | |
143 | fHistClustersPt(0), | |
144 | fHistLeadingJetPt(0), | |
145 | fHistJetsPt(0), | |
146 | fHistBackgroundDone(0), | |
147 | fHistJTPta(0), | |
148 | fHistLogJTPta(0), | |
149 | fHistJTPta_all(0), | |
150 | fHistJTBg(0), | |
151 | fHistLogJTBg(0), | |
152 | fHistBgMulti(0), | |
153 | fHistBgPt(0), | |
154 | fHistJetEta(0), | |
155 | fHistJetMulti(0), | |
156 | fHistJetTracksPt(0), | |
157 | fhTrackingEfficiency(0), | |
158 | fNpttBins(1), | |
159 | fNptaBins(1), | |
160 | fJetsCont(0), | |
161 | //fJetsConts(0), | |
162 | //nJetsConts(0), | |
163 | fTracksCont(0), | |
164 | fCaloClustersCont(0), | |
165 | fTracks(0), | |
166 | fTrackArrayName("nonejk"), | |
167 | runPeriod(""), | |
168 | fEfficiency(0), | |
169 | debug(0) | |
170 | { | |
171 | // Standard constructor. | |
172 | fHistTracksPt = new TH1*[fNcentBins]; | |
173 | fHistTracksJt = new TH1*[fNcentBins]; | |
174 | fHistClustersPt = new TH1*[fNcentBins]; | |
175 | fHistLeadingJetPt = new TH1*[fNcentBins]; | |
176 | fHistJetsPt = new TH1**[fNcentBins]; | |
177 | fHistBackgroundDone = new TH1**[fNcentBins]; | |
178 | fHistJTPta = new TH1***[fNcentBins]; | |
179 | fHistLogJTPta = new TH1***[fNcentBins]; | |
180 | fHistJTPta_all = new TH1***[fNcentBins]; | |
181 | fHistJTBg = new TH1***[fNcentBins]; | |
182 | fHistLogJTBg = new TH1***[fNcentBins]; | |
183 | fHistBgMulti = new TH1**[fNcentBins]; | |
184 | fHistBgPt = new TH1**[fNcentBins]; | |
185 | fHistJetEta = new TH1**[fNcentBins]; | |
186 | fHistJetMulti = new TH1**[fNcentBins]; | |
187 | fHistJetTracksPt = new TH1**[fNcentBins]; | |
188 | fhTrackingEfficiency = new TProfile*[fNcentBins]; | |
189 | //CentBinBorders = new Double_t[10]; | |
190 | ||
191 | ||
192 | for (Int_t i = 0; i < fNcentBins; i++) { | |
193 | fHistJTPta[i] = 0; | |
194 | fHistLogJTPta[i] = 0; | |
195 | fHistJTPta_all[i] = 0; | |
196 | fHistJTBg[i] = 0; | |
197 | fHistLogJTBg[i] = 0; | |
198 | fHistBackgroundDone[i] = 0; | |
199 | fHistTracksPt[i] = 0; | |
200 | fHistTracksJt[i] = 0; | |
201 | fHistClustersPt[i] = 0; | |
202 | fHistLeadingJetPt[i] = 0; | |
203 | fHistJetsPt[i] = 0; | |
204 | fHistBgMulti[i] = 0; | |
205 | fHistBgPt[i] = 0; | |
206 | fHistJetEta[i] = 0; | |
207 | fHistJetMulti[i] = 0; | |
208 | fHistJetTracksPt[i] = 0; | |
209 | fhTrackingEfficiency[i] = 0; | |
210 | } | |
211 | ||
212 | /*for(Int_t i = 0; i < nJetsConts; i++){ | |
213 | fJetsConts[i] = 0; | |
214 | }*/ | |
215 | SetMakeGeneralHistograms(kTRUE); | |
216 | } | |
217 | ||
218 | //________________________________________________________________________ | |
219 | AliAnalysisTaskJetJTJT::~AliAnalysisTaskJetJTJT() | |
220 | { | |
221 | // Destructor. | |
222 | } | |
223 | ||
224 | ||
225 | void AliAnalysisTaskJetJTJT::setCentBinBorders( int n, Double_t *c){ | |
226 | fNcentBins=n; | |
227 | if(debug > 0){ | |
228 | cout << "AliAnalysisTaskJetJTJT::setCentBinBorders: " << endl; | |
229 | } | |
230 | for(int i= 0 ; i < fNcentBins; i++){ | |
231 | CentBinBorders[i]= c[i]; | |
232 | if(debug > 0) | |
233 | cout << CentBinBorders[i] << endl; | |
234 | } | |
235 | } | |
236 | ||
237 | void AliAnalysisTaskJetJTJT::setTriggPtBorders( int n, Double_t *c){ | |
238 | fNpttBins=n; | |
239 | if(debug > 0) | |
240 | cout << "AliAnalysisTaskJetJTJT::setTriggPtBorders: " << endl; | |
241 | for(int i= 0 ; i < fNpttBins; i++){ | |
242 | TriggPtBorders[i]= c[i]; | |
243 | if(debug > 0) | |
244 | cout << TriggPtBorders[i] << endl; | |
245 | } | |
246 | } | |
247 | ||
248 | void AliAnalysisTaskJetJTJT::setAssocPtBorders( int n, Double_t *c){ | |
249 | fNptaBins=n; | |
250 | if(debug > 0) | |
251 | cout << "AliAnalysisTaskJetJTJT::setAssocPtBorders: " << endl; | |
252 | for(int i= 0 ; i < fNptaBins; i++){ | |
253 | AssocPtBorders[i]= c[i]; | |
254 | if(debug > 0) | |
255 | cout << AssocPtBorders[i] << endl; | |
256 | } | |
257 | } | |
258 | ||
259 | ||
260 | //________________________________________________________________________ | |
261 | void AliAnalysisTaskJetJTJT::UserCreateOutputObjects() | |
262 | { | |
263 | // Create user output. | |
264 | ||
265 | AliAnalysisTaskEmcalJet::UserCreateOutputObjects(); | |
266 | if(debug > 0) | |
267 | cout << "Creating Histograms" << endl; | |
268 | ||
269 | fJetsCont = GetJetContainer(0); | |
270 | /*for(int i = 0; i < nJetsConts ; i++){ | |
271 | fJetsConts[0] = GetJetContainer(i); | |
272 | }*/ | |
273 | /*if(fJetsCont) { //get particles and clusters connected to jets | |
274 | fTracksCont = fJetsCont->GetParticleContainer(); | |
275 | fCaloClustersCont = fJetsCont->GetClusterContainer(); | |
276 | } else { //no jets, just analysis tracks and clusters | |
277 | fTracksCont = GetParticleContainer(0); | |
278 | fCaloClustersCont = GetClusterContainer(0); | |
279 | }*/ | |
280 | fTracksCont = GetParticleContainer(0); | |
281 | fCaloClustersCont = GetClusterContainer(0); | |
282 | fTracksCont->SetClassName("AliVTrack"); | |
283 | fCaloClustersCont->SetClassName("AliAODCaloCluster"); | |
284 | ||
285 | TString histname; | |
286 | //Int_t fMinBinJt = 0; | |
287 | //Int_t fMaxBinJt = 5; | |
288 | ||
289 | Int_t NBINSJt=150; | |
290 | double LogBinsJt[NBINSJt+1], LimLJt=0.01, LimHJt=10; | |
291 | double logBWJt = (TMath::Log(LimHJt)-TMath::Log(LimLJt))/(NBINSJt-1); | |
292 | LogBinsJt[0] = 0; | |
293 | for(int ij=1;ij<=NBINSJt;ij++) LogBinsJt[ij]=LimLJt*exp(ij*logBWJt); | |
294 | ||
295 | int NBINSJtW=150; | |
296 | double LimLJtW=TMath::Log(0.01), LimHJtW=TMath::Log(10); | |
297 | ||
298 | //==== Efficiency ==== | |
299 | if(debug > 0) | |
300 | cout << "AliAnalysisTaskJetJTJT::UserCreateOutputObjects: Creating efficiency" << endl; | |
301 | fEfficiency = new JTJTEfficiency; | |
302 | fEfficiency->SetMode( 1 ); // 0:NoEff, 1:Period 2:RunNum 3:Auto | |
303 | fEfficiency->SetDataPath("alien:///alice/cern.ch/user/d/djkim/legotrain/efficieny/data"); // Efficiency root file location local or alien | |
304 | ||
305 | ||
306 | ||
307 | for (Int_t ic = 0; ic < fNcentBins; ic++) { | |
308 | if (fParticleCollArray.GetEntriesFast()>0) { | |
309 | histname = "fHistTracksPt_"; | |
310 | histname += ic; | |
311 | fHistTracksPt[ic] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt / 4); | |
312 | fHistTracksPt[ic]->GetXaxis()->SetTitle("p_{T,track} (GeV/c)"); | |
313 | fHistTracksPt[ic]->GetYaxis()->SetTitle("counts"); | |
314 | fOutput->Add(fHistTracksPt[ic]); | |
315 | } | |
316 | ||
317 | if (fParticleCollArray.GetEntriesFast()>0) { | |
318 | histname = "fHistTracksJt_"; | |
319 | histname += ic; | |
320 | fHistTracksJt[ic] = new TH1F(histname.Data(), histname.Data(), NBINSJt, LogBinsJt); | |
321 | fHistTracksJt[ic]->GetXaxis()->SetTitle("J_{T,track} (GeV/c)"); | |
322 | fHistTracksJt[ic]->GetYaxis()->SetTitle("counts"); | |
323 | fOutput->Add(fHistTracksJt[ic]); | |
324 | } | |
325 | ||
326 | if (fParticleCollArray.GetEntriesFast()>0) { | |
327 | histname = "fhTrackingEfficiency_"; | |
328 | histname += ic; | |
329 | fhTrackingEfficiency[ic] = new TProfile(histname.Data(), histname.Data(), fNbins / 2, fMinBinPt, fMaxBinPt / 2); | |
330 | fhTrackingEfficiency[ic]->GetXaxis()->SetTitle("p_{T,track} (GeV/c)"); | |
331 | fhTrackingEfficiency[ic]->GetYaxis()->SetTitle("counts"); | |
332 | fOutput->Add(fhTrackingEfficiency[ic]); | |
333 | } | |
334 | fHistJTPta[ic] = new TH1**[fNpttBins]; | |
335 | fHistLogJTPta[ic] = new TH1**[fNpttBins]; | |
336 | fHistJTPta_all[ic] = new TH1**[fNpttBins]; | |
337 | fHistJTBg[ic] = new TH1**[fNpttBins]; | |
338 | fHistLogJTBg[ic] = new TH1**[fNpttBins]; | |
339 | fHistJetsPt[ic] = new TH1*[fNpttBins]; | |
340 | fHistBackgroundDone[ic] = new TH1*[fNpttBins]; | |
341 | fHistBgMulti[ic] = new TH1*[fNpttBins]; | |
342 | fHistBgPt[ic] = new TH1*[fNpttBins]; | |
343 | fHistJetEta[ic] = new TH1*[fNpttBins]; | |
344 | fHistJetMulti[ic] = new TH1*[fNpttBins]; | |
345 | fHistJetTracksPt[ic] = new TH1*[fNpttBins]; | |
346 | for(Int_t j=0; j < fNpttBins; j++){ | |
347 | fHistJTPta[ic][j] = new TH1*[fNptaBins]; | |
348 | fHistLogJTPta[ic][j] = new TH1*[fNptaBins]; | |
349 | fHistJTPta_all[ic][j] = new TH1*[fNptaBins]; | |
350 | fHistJTBg[ic][j] = new TH1*[fNptaBins]; | |
351 | fHistLogJTBg[ic][j] = new TH1*[fNptaBins]; | |
352 | for(Int_t k=0; k < fNptaBins; k++){ | |
353 | fHistJTPta[ic][j][k] = 0; | |
354 | fHistLogJTPta[ic][j][k] = 0; | |
355 | fHistJTPta_all[ic][j][k] = 0; | |
356 | fHistJTBg[ic][j][k] = 0; | |
357 | fHistLogJTBg[ic][j][k] = 0; | |
358 | } | |
359 | fHistJetsPt[ic][j] = 0; | |
360 | fHistBackgroundDone[ic][j] = 0; | |
361 | fHistBgMulti[ic][j] = 0; | |
362 | fHistBgPt[ic][j] = 0; | |
363 | fHistJetEta[ic][j] = 0; | |
364 | fHistJetMulti[ic][j] =0; | |
365 | fHistJetTracksPt[ic][j] = 0; | |
366 | } | |
367 | ||
368 | ||
369 | if (fParticleCollArray.GetEntriesFast()>0) { | |
370 | for(Int_t iptt = 0 ; iptt < fNpttBins; iptt++){ | |
371 | for(Int_t ipta = 0 ; ipta < fNptaBins; ipta++){ | |
372 | histname = "hJTPtaD00C"; | |
373 | //histname += ic; | |
374 | histname += Form("C%02dT%02dA%02d", ic, iptt, ipta); | |
375 | if(debug > 1) | |
376 | cout << histname << endl; | |
377 | fHistJTPta[ic][iptt][ipta] = new TH1F(histname.Data(), histname.Data(),NBINSJt, LogBinsJt); | |
378 | fHistJTPta[ic][iptt][ipta]->GetXaxis()->SetTitle("J_{T,track} (GeV/c)"); | |
379 | fHistJTPta[ic][iptt][ipta]->GetYaxis()->SetTitle("counts"); | |
380 | fOutput->Add(fHistJTPta[ic][iptt][ipta]); | |
381 | ||
382 | histname = "hLogJTPtaD00C"; | |
383 | //histname += ic; | |
384 | histname += Form("C%02dT%02dA%02d", ic, iptt, ipta); | |
385 | if(debug > 1) | |
386 | cout << histname << endl; | |
387 | fHistLogJTPta[ic][iptt][ipta] = new TH1F(histname.Data(), histname.Data(), NBINSJtW, LimLJtW, LimHJtW); | |
388 | fHistLogJTPta[ic][iptt][ipta]->GetXaxis()->SetTitle("ln(J_{T,track}/ (GeV/c))"); | |
389 | fHistLogJTPta[ic][iptt][ipta]->GetYaxis()->SetTitle("counts"); | |
390 | fOutput->Add(fHistLogJTPta[ic][iptt][ipta]); | |
391 | ||
392 | histname = "hJTPta_allD00C"; | |
393 | //histname += ic; | |
394 | histname += Form("C%02dT%02dA%02d", ic, iptt, ipta); | |
395 | if(debug > 1) | |
396 | cout << histname << endl; | |
397 | fHistJTPta_all[ic][iptt][ipta] = new TH1F(histname.Data(), histname.Data(), NBINSJt, LogBinsJt); | |
398 | fHistJTPta_all[ic][iptt][ipta]->GetXaxis()->SetTitle("J_{T,track} (GeV/c)"); | |
399 | fHistJTPta_all[ic][iptt][ipta]->GetYaxis()->SetTitle("counts"); | |
400 | fOutput->Add(fHistJTPta_all[ic][iptt][ipta]); | |
401 | ||
402 | histname = "hJTBg"; | |
403 | //histname += ic; | |
404 | histname += Form("C%02dT%02dA%02d", ic, iptt, ipta); | |
405 | if(debug > 1) | |
406 | cout << histname << endl; | |
407 | fHistJTBg[ic][iptt][ipta] = new TH1F(histname.Data(), histname.Data(), NBINSJt, LogBinsJt); | |
408 | fHistJTBg[ic][iptt][ipta]->GetXaxis()->SetTitle("J_{T,track} (GeV/c)"); | |
409 | fHistJTBg[ic][iptt][ipta]->GetYaxis()->SetTitle("counts"); | |
410 | fOutput->Add(fHistJTBg[ic][iptt][ipta]); | |
411 | ||
412 | histname = "hLogJTBg"; | |
413 | //histname += ic; | |
414 | histname += Form("C%02dT%02dA%02d", ic, iptt, ipta); | |
415 | if(debug > 1) | |
416 | cout << histname << endl; | |
417 | fHistLogJTBg[ic][iptt][ipta] = new TH1F(histname.Data(), histname.Data(), NBINSJtW, LimLJtW, LimHJtW); | |
418 | fHistLogJTBg[ic][iptt][ipta]->GetXaxis()->SetTitle("ln(J_{T,track}/ (GeV/c))"); | |
419 | fHistLogJTBg[ic][iptt][ipta]->GetYaxis()->SetTitle("counts"); | |
420 | fOutput->Add(fHistLogJTBg[ic][iptt][ipta]); | |
421 | ||
422 | ||
423 | } | |
424 | } | |
425 | } | |
426 | ||
427 | if (fClusterCollArray.GetEntriesFast()>0) { | |
428 | histname = "fHistClustersPt_"; | |
429 | histname += ic; | |
430 | fHistClustersPt[ic] = new TH1F(histname.Data(), histname.Data(), fNbins / 2, fMinBinPt, fMaxBinPt / 2); | |
431 | fHistClustersPt[ic]->GetXaxis()->SetTitle("p_{T,clus} (GeV/c)"); | |
432 | fHistClustersPt[ic]->GetYaxis()->SetTitle("counts"); | |
433 | fOutput->Add(fHistClustersPt[ic]); | |
434 | } | |
435 | ||
436 | if (fJetCollArray.GetEntriesFast()>0) { | |
437 | histname = "fHistLeadingJetPt_"; | |
438 | histname += ic; | |
439 | fHistLeadingJetPt[ic] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt); | |
440 | fHistLeadingJetPt[ic]->GetXaxis()->SetTitle("p_{T}^{raw} (GeV/c)"); | |
441 | fHistLeadingJetPt[ic]->GetYaxis()->SetTitle("counts"); | |
442 | fOutput->Add(fHistLeadingJetPt[ic]); | |
443 | ||
444 | for(Int_t iptt = 0 ; iptt < fNpttBins; iptt++){ | |
445 | ||
446 | histname = "fHistJetsPt_"; | |
447 | histname += Form("C%02dT%02d", ic, iptt); | |
448 | fHistJetsPt[ic][iptt] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt); | |
449 | fHistJetsPt[ic][iptt]->GetXaxis()->SetTitle("p_{T}^{raw} (GeV/c)"); | |
450 | fHistJetsPt[ic][iptt]->GetYaxis()->SetTitle("counts"); | |
451 | fOutput->Add(fHistJetsPt[ic][iptt]); | |
452 | ||
453 | histname = "fHistBackgroundDone_"; | |
454 | histname += Form("C%02dT%02d", ic, iptt);; | |
455 | fHistBackgroundDone[ic][iptt] = new TH1F(histname.Data(), histname.Data(), 2, -1, 2); | |
456 | fHistBackgroundDone[ic][iptt]->GetXaxis()->SetTitle("Number of jets"); | |
457 | fHistBackgroundDone[ic][iptt]->GetYaxis()->SetTitle("0 = not done, 1 = done"); | |
458 | fOutput->Add(fHistBackgroundDone[ic][iptt]); | |
459 | ||
460 | histname = "fHistJetEta_"; | |
461 | histname += Form("C%02dT%02d", ic, iptt); | |
462 | fHistJetEta[ic][iptt] = new TH1F(histname.Data(), histname.Data(), fNbins, -5, 5); | |
463 | fHistJetEta[ic][iptt]->GetXaxis()->SetTitle("#eta"); | |
464 | fHistJetEta[ic][iptt]->GetYaxis()->SetTitle("jets"); | |
465 | fOutput->Add(fHistJetEta[ic][iptt]); | |
466 | ||
467 | histname = "fHistJetMulti_"; | |
468 | histname += Form("C%02dT%02d", ic, iptt); | |
469 | fHistJetMulti[ic][iptt] = new TH1F(histname.Data(), histname.Data(), 200, 0, 200); | |
470 | fHistJetMulti[ic][iptt]->GetXaxis()->SetTitle("Multiplicity"); | |
471 | fHistJetMulti[ic][iptt]->GetYaxis()->SetTitle("jets"); | |
472 | fOutput->Add(fHistJetMulti[ic][iptt]); | |
473 | ||
474 | histname = "fHistBgMulti_"; | |
475 | histname += Form("C%02dT%02d", ic, iptt); | |
476 | fHistBgMulti[ic][iptt] = new TH1F(histname.Data(), histname.Data(), 200, 0, 200); | |
477 | fHistBgMulti[ic][iptt]->GetXaxis()->SetTitle("Multiplicity"); | |
478 | fHistBgMulti[ic][iptt]->GetYaxis()->SetTitle("Events"); | |
479 | fOutput->Add(fHistBgMulti[ic][iptt]); | |
480 | ||
481 | histname = "fHistBgPt_"; | |
482 | histname += Form("C%02dT%02d", ic, iptt); | |
483 | fHistBgPt[ic][iptt] = new TH1F(histname.Data(), histname.Data(), fNbins,fMinBinPt, fMaxBinPt/4); | |
484 | fHistBgPt[ic][iptt]->GetXaxis()->SetTitle("Multiplicity"); | |
485 | fHistBgPt[ic][iptt]->GetYaxis()->SetTitle("jets"); | |
486 | fOutput->Add(fHistBgPt[ic][iptt]); | |
487 | ||
488 | histname = "fHistJetTracksPt_"; | |
489 | //histname += ic; | |
490 | histname += Form("C%02dT%02d", ic, iptt); | |
491 | if(debug > 1) | |
492 | cout << histname << endl; | |
493 | fHistJetTracksPt[ic][iptt] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt/4); | |
494 | fHistJetTracksPt[ic][iptt]->GetXaxis()->SetTitle("p_{T,track} (GeV/c)"); | |
495 | fHistJetTracksPt[ic][iptt]->GetYaxis()->SetTitle("counts"); | |
496 | fOutput->Add(fHistJetTracksPt[ic][iptt]); | |
497 | } | |
498 | ||
499 | /* | |
500 | if (!(GetJetContainer()->GetRhoName().IsNull())) { | |
501 | histname = "fHistJetsCorrPtArea_"; | |
502 | histname += i; | |
503 | fHistJetsCorrPtArea[i] = new TH2F(histname.Data(), histname.Data(), fNbins*2, -fMaxBinPt, fMaxBinPt, 30, 0, 3); | |
504 | fHistJetsCorrPtArea[i]->GetXaxis()->SetTitle("p_{T}^{corr} [GeV/c]"); | |
505 | fHistJetsCorrPtArea[i]->GetYaxis()->SetTitle("area"); | |
506 | fOutput->Add(fHistJetsCorrPtArea[i]); | |
507 | } | |
508 | */ | |
509 | } | |
510 | } | |
511 | ||
512 | PostData(1, fOutput); // Post data for ALL output slots > 0 here. | |
513 | } | |
514 | ||
515 | //________________________________________________________________________ | |
516 | Bool_t AliAnalysisTaskJetJTJT::FillHistograms() | |
517 | { | |
518 | // Fill histograms. | |
519 | ||
520 | AliCentrality *aliCent = InputEvent()->GetCentrality(); | |
521 | fCentBin = 0; | |
522 | if (aliCent) { | |
523 | //fCent = aliCent->GetCentralityPercentile(fCentEst.Data()); | |
524 | fCent = aliCent->GetCentralityPercentile("V0M"); | |
525 | /*if(debug > 0){ | |
526 | cout << "Centrality " << fCent << endl; | |
527 | }*/ | |
528 | for(int ic = 0; ic < fNcentBins; ic++){ | |
529 | /*if(debug > 0){ | |
530 | cout << "ic: " << ic << " / " << fNcentBins << endl; | |
531 | cout << "Centrality bin " << fCentBin << endl; | |
532 | cout << "Border: " << CentBinBorders[ic] << endl; | |
533 | } */ | |
534 | if(fCent > CentBinBorders[ic]){ | |
535 | fCentBin = ic; | |
536 | } | |
537 | } | |
538 | //cout << "Centrality bin: " << fCentBin << endl; | |
539 | } else { | |
540 | AliWarning(Form("%s: Could not retrieve centrality information! Assuming 99", GetName())); | |
541 | fCentBin = 3; | |
542 | } | |
543 | int fHadronSelectionCut = 5; //5=Hybrid cut | |
544 | if (fTracksCont) { | |
545 | AliVTrack *track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle(0)); | |
546 | while(track) { | |
547 | double ptt = track->Pt(); | |
548 | ||
549 | //<<<<<<<<<<<< Efficiency >>>>>>>>>>> | |
550 | //AliJBaseTrack *triggTr = (AliJBaseTrack*)fchargedHadronList->At(i); | |
551 | ||
552 | double effCorr = 1./fEfficiency->GetCorrection(ptt, fHadronSelectionCut, fCent); // here you generate warning if ptt>30 | |
553 | //double effCorr = 1.; | |
554 | fhTrackingEfficiency[fCentBin]->Fill( ptt, 1./effCorr ); | |
555 | //triggTr->SetTrackEff( 1./effCorr ); | |
556 | //<<<<<<<<<<<< Efficiency >>>>>>>>>>> | |
557 | ||
558 | if(ptt > 0 && 1.0/ptt > 0){ | |
559 | fHistTracksPt[fCentBin]->Fill(ptt,1./ptt*effCorr); | |
560 | } | |
561 | ||
562 | ||
563 | track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle()); | |
564 | } | |
565 | } | |
566 | ||
567 | if (fCaloClustersCont) { | |
568 | AliVCluster *cluster = fCaloClustersCont->GetNextAcceptCluster(0); | |
569 | while(cluster) { | |
570 | TLorentzVector nPart; | |
571 | cluster->GetMomentum(nPart, fVertex); | |
572 | fHistClustersPt[fCentBin]->Fill(nPart.Pt()); | |
573 | ||
574 | cluster = fCaloClustersCont->GetNextAcceptCluster(); | |
575 | } | |
576 | } | |
577 | ||
578 | Int_t fPttBin, fPtaBin; | |
579 | fPtaBin = 0; | |
580 | ||
581 | ||
582 | if (fJetsCont) { | |
583 | //Int_t Njets = fJetsCont->GetNJets(); | |
584 | Int_t Njets = 0; | |
585 | if(debug > 1){ | |
586 | cout << "Number of Jets: " << Njets << endl; | |
587 | } | |
588 | ||
589 | //Make a array to hold jets to be tested in background jT | |
590 | Float_t jetPhis[200] = {}; | |
591 | Float_t jetEtas[200] = {}; | |
592 | AliEmcalJet *jet = fJetsCont->GetNextAcceptJet(0); | |
593 | Int_t ij = 0; | |
594 | while(jet) { | |
595 | //cout << "Jet found " << ij << " pt: " << jet->Pt() << endl; | |
596 | if(jet->Pt() > 5){ //Only consider jets with pT > 5 GeV | |
597 | jetPhis[ij] = jet->Phi(); | |
598 | jetEtas[ij] = jet->Eta(); | |
599 | ij++; | |
600 | Njets++; | |
601 | if(debug > 1) | |
602 | cout << "i: " << ij << " jetPhi: " << jetPhis[ij] << " jetEta: " << jetEtas[ij] << endl; | |
603 | }else{ | |
604 | //jetPhis[ij] = 100; | |
605 | //jetEtas[ij] = 100; | |
606 | //Njets--; | |
607 | if(debug > 1) | |
608 | cout << "jetPt: " << jet->Pt() << " jetPhi: " << jet->Phi() << " jetEta: " << jet->Eta() << endl; | |
609 | } | |
610 | //i++; | |
611 | jet = fJetsCont->GetNextAcceptJet(); | |
612 | } | |
613 | ||
614 | ||
615 | //fTracks =dynamic_cast <TClonesArray*>(InputEvent()->FindListObject( fTrackArrayName.Data() )); | |
616 | jet = fJetsCont->GetNextAcceptJet(0); | |
617 | while(jet) { | |
618 | if(jet->Pt() > 5){ | |
619 | if(jet->Eta() < -0.4 || jet->Eta() > 0.4){ | |
620 | if(debug > 0) | |
621 | cout << "Jet outside eta range, Eta: " << jet->Eta() << endl; | |
622 | jet = fJetsCont->GetNextAcceptJet(); | |
623 | continue; | |
624 | } | |
625 | //cout << "Jet found " << ij << " pt: " << jet->Pt() << endl; | |
626 | //Get the trigger pT bin | |
627 | fPttBin = 0; | |
628 | for(int iptt = 0 ; iptt < fNpttBins; iptt++){ | |
629 | if(jet->Pt() > TriggPtBorders[iptt]){ | |
630 | fPttBin = iptt; | |
631 | } | |
632 | ||
633 | } | |
634 | fHistJetEta[fCentBin][fPttBin]->Fill(jet->Eta()); | |
635 | if(jet->Pt() > 0 && 1.0/jet->Pt() > 0){ | |
636 | fHistJetsPt[fCentBin][fPttBin]->Fill(jet->Pt(),1.0/jet->Pt()); //Fill jet dN/(pT dpT) | |
637 | } | |
638 | ||
639 | /*if (fHistJetsCorrPtArea[fCentBin]) { | |
640 | Float_t corrPt = jet->Pt() - fJetsCont->GetRhoVal() * jet->Area(); | |
641 | fHistJetsCorrPtArea[fCentBin]->Fill(corrPt, jet->Area()); | |
642 | }*/ | |
643 | //Float_t jetp = sqrt(jet->Px()*jet->Px()+jet->Py()*jet->Py()+jet->Pz()*jet->Pz()); //Jet pT norm | |
644 | ||
645 | ||
646 | Int_t nTrack = jet->GetNumberOfTracks(); | |
647 | if (debug > 0) | |
648 | cout << "Number of tracks " << nTrack << " Jet Pt: " << jet->Pt() << endl; | |
649 | fHistJetMulti[fCentBin][fPttBin]->Fill(nTrack); | |
650 | for(Int_t it = 0; it < nTrack; it++ ){ | |
651 | AliVParticle *track = (AliVParticle*)jet->TrackAt( it, fTracks ); | |
652 | if( !track ){ | |
653 | cout << "No Track found" << endl; | |
654 | continue; | |
655 | } | |
656 | fPtaBin = 0; //Get the associated pT bin | |
657 | for(int ipta = 0 ; ipta < fNptaBins; ipta++){ | |
658 | if(track->Pt() > AssocPtBorders[ipta]){ | |
659 | fPtaBin = ipta; | |
660 | } | |
661 | } | |
662 | fHistJetTracksPt[fCentBin][fPttBin]->Fill(track->Pt()); | |
663 | if(debug > 2) | |
664 | cout << "Filling fHistJetTracksPt C" << fCentBin << " T" << fPttBin << endl; | |
665 | /*Float_t dotproduct = track->Px()*jet->Px()+track->Py()*jet->Py()+track->Pz()*jet->Pz(); // p_track dot p_jet | |
666 | Float_t constp = sqrt(track->Px()*track->Px()+track->Py()*track->Py()+track->Pz()*track->Pz()); // track pT norm | |
667 | Float_t normproduct = constp*jetp; // jet pT norm times track pT norm | |
668 | Float_t costheta2 = dotproduct/normproduct; | |
669 | //Float_t sintheta = sqrt(1-costheta2*costheta2); | |
670 | Float_t jt = constp*sqrt(1-costheta2*costheta2);*/ | |
671 | Float_t jt = getJt(track,jet,0); | |
672 | double effCorr = 1./fEfficiency->GetCorrection(track->Pt(), fHadronSelectionCut, fCent); // here you generate warning if ptt>30 | |
673 | //double effCorr = 1.; | |
674 | if(jt > 0 && 1.0/jt > 0){ | |
675 | fHistTracksJt[fCentBin]->Fill(jt,1.0/jt*effCorr); //Fill dN/(djT jT) | |
676 | fHistJTPta[fCentBin][fPttBin][fPtaBin]->Fill(jt,1.0/jt*effCorr); //Fill dN/(djT jT) | |
677 | fHistLogJTPta[fCentBin][fPttBin][fPtaBin]->Fill(TMath::Log(jt),1.0/jt*effCorr); //Fill logarithmic dN/(djT jT) | |
678 | } | |
679 | if(debug > 1) | |
680 | cout << "Filling JT C" << fCentBin << "T" << fPttBin << "A" << fPtaBin << " jt:" << jt << " with " << 1.0/jt*effCorr << endl; | |
681 | } | |
682 | ||
683 | //Get Jet azimuth and rapidity of jet | |
684 | Float_t jetAngle = jet->Phi(); | |
685 | Float_t jetRap = jet->Eta(); | |
686 | ||
687 | //Rotate jet angle for background cone | |
688 | Float_t rotatedAngle = jetAngle+TMath::Pi()/2; | |
689 | if(rotatedAngle > TMath::Pi()*2){ | |
690 | rotatedAngle = rotatedAngle- TMath::Pi()*2; | |
691 | } | |
692 | Float_t jetArea = jet->Area(); | |
693 | Float_t testRadius = TMath::Sqrt(jetArea/TMath::Pi()); | |
694 | ||
695 | Bool_t doBg = 1; | |
696 | ||
697 | //Test if there are jets in the background test cone | |
698 | for(int i_j = 0; i_j < Njets; i_j++){ | |
699 | Float_t diffR = TMath::Sqrt(TMath::Power(jetPhis[i_j]-rotatedAngle,2)+TMath::Power(jetEtas[i_j]-jetRap,2)); | |
700 | if(debug > 1){ | |
701 | cout << "i_j: " << i_j << " JetPhi: " << jetPhis[i_j] << " jetEta: " << jetEtas[i_j] << endl; | |
702 | cout << "DiffR: " << diffR << " doBG: " << doBg <<endl; | |
703 | } | |
704 | if(diffR < testRadius *2){ //Jets muts be at least 2*cone radius away from the background cone axis | |
705 | doBg =0; | |
706 | break; | |
707 | } | |
708 | ||
709 | } | |
710 | ||
711 | // Do jT for all particles in respect to jet axis | |
712 | if (fTracksCont) { | |
713 | int counter = 0; | |
714 | AliVTrack *track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle(0)); | |
715 | while(track) { | |
716 | /*Float_t dotproduct = track->Px()*jet->Px()+track->Py()*jet->Py()+track->Pz()*jet->Pz(); | |
717 | Float_t constp = sqrt(track->Px()*track->Px()+track->Py()*track->Py()+track->Pz()*track->Pz()); | |
718 | Float_t normproduct = constp*jetp; | |
719 | Float_t costheta2 = dotproduct/normproduct; | |
720 | //Float_t sintheta = sqrt(1-costheta2*costheta2); | |
721 | Float_t jt = constp*sqrt(1-costheta2*costheta2);*/ | |
722 | Double_t jt = getJt(track,jet,0); | |
723 | double effCorr = 1./fEfficiency->GetCorrection(track->Pt(), fHadronSelectionCut, fCent); // here you generate warning if ptt>30 | |
724 | if(jt > 0 && 1.0/jt > 0){ | |
725 | fHistJTPta_all[fCentBin][fPttBin][fPtaBin]->Fill(jt,1.0/jt*effCorr); | |
726 | } | |
727 | for(int ipta = 0 ; ipta < fNptaBins; ipta++){ | |
728 | if(track->Pt() > AssocPtBorders[ipta]){ | |
729 | fPtaBin = ipta; | |
730 | } | |
731 | } | |
732 | ||
733 | //If background is to be filled | |
734 | if(doBg){ | |
735 | Float_t diffR = TMath::Sqrt(TMath::Power(track->Phi()-rotatedAngle,2)+TMath::Power(track->Eta()-jetRap,2)); | |
736 | //Particles in the rotated cone | |
737 | if(diffR < testRadius){ | |
738 | counter++; | |
739 | fHistBgPt[fCentBin][fPttBin]->Fill(track->Pt()); | |
740 | /*dotproduct = -track->Px()*jet->Py()+track->Py()*jet->Px()+track->Pz()*jet->Pz(); | |
741 | constp = sqrt(track->Px()*track->Px()+track->Py()*track->Py()+track->Pz()*track->Pz()); | |
742 | normproduct = constp*jetp; | |
743 | costheta2 = dotproduct/normproduct; | |
744 | //sintheta = sqrt(1-costheta2*costheta2); | |
745 | jt = constp*sqrt(1-costheta2*costheta2);*/ | |
746 | jt = getJt(track,jet,1); | |
747 | if(jt > 0 && 1.0/jt > 0){ | |
748 | fHistJTBg[fCentBin][fPttBin][fPtaBin]->Fill(jt,1.0/jt*effCorr); | |
749 | fHistLogJTBg[fCentBin][fPttBin][fPtaBin]->Fill(TMath::Log(jt),1.0/jt*effCorr); | |
750 | } | |
751 | if(debug > 1) | |
752 | cout << "Filling Background C" << fCentBin << "T" << fPttBin << "A" << fPtaBin << " jt:" << jt << " with " << 1.0/jt*effCorr << endl; | |
753 | //Fill background jT | |
754 | ||
755 | } | |
756 | } | |
757 | track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle()); | |
758 | } | |
759 | if(doBg){ | |
760 | fHistBgMulti[fCentBin][fPttBin]->Fill(counter); | |
761 | } | |
762 | } | |
763 | if(doBg){ | |
764 | fHistBackgroundDone[fCentBin][fPttBin]->Fill(1); | |
765 | }else{ | |
766 | fHistBackgroundDone[fCentBin][fPttBin]->Fill(0); | |
767 | } | |
768 | ||
769 | } | |
770 | jet = fJetsCont->GetNextAcceptJet(); | |
771 | ||
772 | } | |
773 | jet = fJetsCont->GetLeadingJet(); | |
774 | if(jet){ | |
775 | if(jet->Pt() > 0 && 1.0/jet->Pt() > 0){ | |
776 | fHistLeadingJetPt[fCentBin]->Fill(jet->Pt(),1.0/jet->Pt()); | |
777 | } | |
778 | } | |
779 | } | |
780 | ||
781 | CheckClusTrackMatching(); | |
782 | ||
783 | return kTRUE; | |
784 | } | |
785 | ||
786 | ||
787 | //----------------------------------------------------------------------- | |
788 | Double_t AliAnalysisTaskJetJTJT::getJt(AliVTrack *track, AliEmcalJet *jet,int reverse){ | |
789 | Float_t dotproduct = 0; | |
790 | Float_t jetp = sqrt(jet->Px()*jet->Px()+jet->Py()*jet->Py()+jet->Pz()*jet->Pz()); //Jet pT norm | |
791 | if(reverse){ | |
792 | dotproduct = -track->Px()*jet->Py()+track->Py()*jet->Px()+track->Pz()*jet->Pz(); | |
793 | } else{ | |
794 | dotproduct = track->Px()*jet->Px()+track->Py()*jet->Py()+track->Pz()*jet->Pz(); | |
795 | } | |
796 | Float_t constp = sqrt(track->Px()*track->Px()+track->Py()*track->Py()+track->Pz()*track->Pz()); | |
797 | Float_t normproduct = constp*jetp; | |
798 | Float_t costheta2 = dotproduct/normproduct; | |
799 | //Float_t sintheta = sqrt(1-costheta2*costheta2); | |
800 | Float_t jt = constp*sqrt(1-costheta2*costheta2); | |
801 | return jt; | |
802 | } | |
803 | ||
804 | Double_t AliAnalysisTaskJetJTJT::getJt(AliVParticle *track, AliEmcalJet *jet,int reverse){ | |
805 | Float_t dotproduct = 0; | |
806 | Float_t jetp = sqrt(jet->Px()*jet->Px()+jet->Py()*jet->Py()+jet->Pz()*jet->Pz()); //Jet pT norm | |
807 | if(reverse){ | |
808 | dotproduct = -track->Px()*jet->Py()+track->Py()*jet->Px()+track->Pz()*jet->Pz(); | |
809 | } else{ | |
810 | dotproduct = track->Px()*jet->Px()+track->Py()*jet->Py()+track->Pz()*jet->Pz(); | |
811 | } | |
812 | Float_t constp = sqrt(track->Px()*track->Px()+track->Py()*track->Py()+track->Pz()*track->Pz()); | |
813 | Float_t normproduct = constp*jetp; | |
814 | Float_t costheta2 = dotproduct/normproduct; | |
815 | //Float_t sintheta = sqrt(1-costheta2*costheta2); | |
816 | Float_t jt = constp*sqrt(1-costheta2*costheta2); | |
817 | return jt; | |
818 | } | |
819 | ||
820 | //________________________________________________________________________ | |
821 | void AliAnalysisTaskJetJTJT::CheckClusTrackMatching() | |
822 | { | |
823 | ||
824 | if(!fTracksCont || !fCaloClustersCont) | |
825 | return; | |
826 | ||
827 | Double_t deta = 999; | |
828 | Double_t dphi = 999; | |
829 | ||
830 | //Get closest cluster to track | |
831 | AliVTrack *track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle(0)); | |
832 | while(track) { | |
833 | //Get matched cluster | |
834 | Int_t emc1 = track->GetEMCALcluster(); | |
835 | if(fCaloClustersCont && emc1>=0) { | |
836 | AliVCluster *clusMatch = fCaloClustersCont->GetCluster(emc1); | |
837 | if(clusMatch) { | |
838 | AliPicoTrack::GetEtaPhiDiff(track, clusMatch, dphi, deta); | |
839 | //fHistPtDEtaDPhiTrackClus->Fill(track->Pt(),deta,dphi); | |
840 | } | |
841 | } | |
842 | track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle()); | |
843 | } | |
844 | ||
845 | //Get closest track to cluster | |
846 | AliVCluster *cluster = fCaloClustersCont->GetNextAcceptCluster(0); | |
847 | while(cluster) { | |
848 | TLorentzVector nPart; | |
849 | cluster->GetMomentum(nPart, fVertex); | |
850 | fHistClustersPt[fCentBin]->Fill(nPart.Pt()); | |
851 | ||
852 | //Get matched track | |
853 | AliVTrack *mt = NULL; | |
854 | AliAODCaloCluster *acl = dynamic_cast<AliAODCaloCluster*>(cluster); | |
855 | if(acl) { | |
856 | if(acl->GetNTracksMatched()>1) | |
857 | mt = static_cast<AliVTrack*>(acl->GetTrackMatched(0)); | |
858 | } | |
859 | else { | |
860 | AliESDCaloCluster *ecl = dynamic_cast<AliESDCaloCluster*>(cluster); | |
861 | Int_t im = ecl->GetTrackMatchedIndex(); | |
862 | if(fTracksCont && im>=0) { | |
863 | mt = static_cast<AliVTrack*>(fTracksCont->GetParticle(im)); | |
864 | } | |
865 | } | |
866 | if(mt) { | |
867 | AliPicoTrack::GetEtaPhiDiff(mt, cluster, dphi, deta); | |
868 | //fHistPtDEtaDPhiClusTrack->Fill(nPart.Pt(),deta,dphi); | |
869 | ||
870 | //debugging | |
871 | /* | |
872 | if(mt->IsEMCAL()) { | |
873 | Int_t emc1 = mt->GetEMCALcluster(); | |
874 | Printf("current id: %d emc1: %d",fCaloClustersCont->GetCurrentID(),emc1); | |
875 | AliVCluster *clm = fCaloClustersCont->GetCluster(emc1); | |
876 | AliPicoTrack::GetEtaPhiDiff(mt, clm, dphi, deta); | |
877 | Printf("deta: %f dphi: %f",deta,dphi); | |
878 | } | |
879 | */ | |
880 | } | |
881 | cluster = fCaloClustersCont->GetNextAcceptCluster(); | |
882 | } | |
883 | } | |
884 | ||
885 | //________________________________________________________________________ | |
886 | /*AliJetContainer* AliAnalysisTaskJetJTJT::AddJetContainer(const char *n, TString defaultCutType, Float_t jetRadius) { | |
887 | ||
888 | AliAnalysisTaskEmcalJet::ExecOnce(); | |
889 | nJetsConts++; | |
890 | AliJetContainer *cont = 0x0; | |
891 | cont = AliAnalysisTaskEmcalJet::AddJetContainer(n,defaultCutType,jetRadius); | |
892 | return cont; | |
893 | }*/ | |
894 | ||
895 | ||
896 | //________________________________________________________________________ | |
897 | void AliAnalysisTaskJetJTJT::ExecOnce() { | |
898 | ||
899 | if(debug > 0){ | |
900 | cout << "AliAnalysisTaskJetJTJT::ExecOnce(): " << endl; | |
901 | cout << "Get fTracks from " << fTrackArrayName.Data() << endl; | |
902 | } | |
903 | fTracks =dynamic_cast <TClonesArray*>(InputEvent()->FindListObject( fTrackArrayName.Data() )); | |
904 | ||
905 | AliAnalysisTaskEmcalJet::ExecOnce(); | |
906 | if(debug > 1) | |
907 | cout << "Efficiency: Set Run Period Name " << runPeriod << endl; | |
908 | fEfficiency->SetPeriodName(runPeriod); | |
909 | if(debug > 1) | |
910 | cout << "Efficiency: Set Run number " << InputEvent()->GetRunNumber() << endl; | |
911 | fEfficiency->SetRunNumber( InputEvent()->GetRunNumber() ); //TODO Get run Number | |
912 | if(debug > 1) | |
913 | cout << "Efficiency: Load()" << endl; | |
914 | fEfficiency->Load(); | |
915 | if(debug > 1) | |
916 | cout << "fEfficiency loaded" << endl; | |
917 | ||
918 | ||
919 | ||
920 | if (fJetsCont && fJetsCont->GetArray() == 0) fJetsCont = 0; | |
921 | if (fTracksCont && fTracksCont->GetArray() == 0) fTracksCont = 0; | |
922 | if (fCaloClustersCont && fCaloClustersCont->GetArray() == 0) fCaloClustersCont = 0; | |
923 | ||
924 | } | |
925 | ||
926 | //________________________________________________________________________ | |
927 | Bool_t AliAnalysisTaskJetJTJT::Run() | |
928 | { | |
929 | // Run analysis code here, if needed. It will be executed before FillHistograms(). | |
930 | ||
931 | return kTRUE; // If return kFALSE FillHistogram() will NOT be executed. | |
932 | } | |
933 | ||
934 | //________________________________________________________________________ | |
935 | void AliAnalysisTaskJetJTJT::Terminate(Option_t *) | |
936 | { | |
937 | // Called once at the end of the analysis. | |
938 | } | |
939 | ||
940 | ||
941 | ||
942 | //________________________________________________________________________ | |
943 | JTJTEfficiency::JTJTEfficiency(): | |
944 | fMode(kAuto), | |
945 | fPeriod(-1), | |
946 | fDataPath(""), | |
947 | fName(""), | |
948 | fPeriodStr(""), | |
949 | fMCPeriodStr(""), | |
950 | fRunNumber(0), | |
951 | fTag(""), | |
952 | fInputRootName(""), | |
953 | fInputRoot(NULL), | |
954 | fCentBinAxis(0x0) | |
955 | { | |
66b055a5 | 956 | for (Int_t i=0; i < 3; i++) fEffDir[i] = 0; |
11476c41 | 957 | |
958 | } | |
959 | ||
960 | JTJTEfficiency::JTJTEfficiency(const JTJTEfficiency& obj) : | |
961 | fMode(obj.fMode), | |
962 | fPeriod(obj.fPeriod), | |
963 | fDataPath(obj.fDataPath), | |
964 | fName(obj.fName), | |
965 | fPeriodStr(obj.fPeriodStr), | |
966 | fMCPeriodStr(obj.fMCPeriodStr), | |
967 | fRunNumber(obj.fRunNumber), | |
968 | fTag(obj.fTag), | |
969 | fInputRootName(obj.fInputRootName), | |
970 | fInputRoot(obj.fInputRoot), | |
971 | fCentBinAxis(obj.fCentBinAxis) | |
972 | { | |
973 | // copy constructor TODO: handling of pointer members | |
974 | JUNUSED(obj); | |
66b055a5 | 975 | for (Int_t i=0; i < 3; i++) fEffDir[i] = obj.fEffDir[i]; |
11476c41 | 976 | } |
977 | ||
978 | JTJTEfficiency& JTJTEfficiency::operator=(const JTJTEfficiency& obj){ | |
979 | // equal sign operator TODO: content | |
980 | JUNUSED(obj); | |
981 | return *this; | |
982 | } | |
983 | ||
984 | ||
985 | //________________________________________________________________________ | |
986 | double JTJTEfficiency::GetCorrection( double pt, int icut , double cent ) const { | |
987 | if( fMode == kNotUse ) return 1; | |
988 | int icent = fCentBinAxis->FindBin( cent ) -1 ; | |
989 | if( icent < 0 || icent > fCentBinAxis->GetNbins()-1 ) { | |
990 | cout<<"J_WARNING : Centrality "<<cent<<" is out of CentBinBorder"<<endl; | |
991 | return 1; | |
992 | } | |
993 | // TODO error check for icent; | |
994 | int ivtx = 0; | |
995 | if( ! fCorrection[ivtx][icent][icut] ) { | |
996 | cout<<"J_WARNING : No Eff Info "<<pt<<"\t"<<icut<<"\t"<<cent<<"\t"<<icent<<endl; | |
997 | return 1; | |
998 | } | |
999 | TGraphErrors * gr = fCorrection[ivtx][icent][icut]; | |
1000 | //=== TEMPERORY SETTING. IT will be removed soon. | |
1001 | if( pt > 30 ) pt = 30; // Getting eff of 30GeV for lager pt | |
1002 | double cor = gr->Eval(pt); | |
1003 | if ( cor < 0.2 ) cor = 0.2; | |
1004 | return cor; | |
1005 | } | |
1006 | ||
1007 | ||
1008 | TString JTJTEfficiency::GetEffName() { | |
1009 | /* | |
1010 | 1. kNotUse : no Load, efficiency is 1 always | |
1011 | 2. has fInputRootName : Load that or crash | |
1012 | 3. has fName : Load fName [+runnumber] or crash | |
1013 | 4. has runnumber : Find Good MC period from AliJRunTable, or crash | |
1014 | 3. has period : Find Good MC period from AliJRunTable, or crash | |
1015 | ||
1016 | } | |
1017 | */ | |
1018 | if(fPeriodStr == "LHC10b"){ | |
1019 | fInputRootName = "Eff--LHC10b-LHC10d1-0-.root"; | |
1020 | } | |
1021 | if(fPeriodStr == "LHC10c"){ | |
1022 | fInputRootName = "Eff--LHC10c-LHC10d4-0-.root"; | |
1023 | } | |
1024 | if(fPeriodStr == "LHC10d"){ | |
1025 | fInputRootName = "Eff--LHC10d-LHC10f6a-0-.root"; | |
1026 | } | |
1027 | if(fPeriodStr == "LHC10e"){ | |
1028 | fInputRootName = "Eff--LHC10e-LHC10e20-0-.root"; | |
1029 | } | |
1030 | if(fPeriodStr == "LHC10h"){ | |
1031 | fInputRootName = "Eff--LHC10h-LHC11a10a_bis-0-.root"; | |
1032 | } | |
1033 | if(fPeriodStr == "LHC11a"){ | |
1034 | fInputRootName = "Eff--LHC11a-LHC11b10a-0-.root"; | |
1035 | } | |
1036 | if(fPeriodStr == "LHC13b"){ | |
1037 | fInputRootName = "Eff--LHC13b-LHC13b2-efix_p1-0-.root"; | |
1038 | } | |
1039 | ||
1040 | if(fPeriodStr == "LHC13c"){ | |
1041 | fInputRootName = "Eff--LHC13c-LHC13b2-efix_p1-0-.root"; | |
1042 | } | |
1043 | if(fPeriodStr == "LHC13d"){ | |
1044 | fInputRootName = "Eff--LHC13d-LHC13b2-efix_p1-0-.root"; | |
1045 | } | |
1046 | if(fPeriodStr == "LHC13e"){ | |
1047 | fInputRootName = "Eff--LHC13e-LHC13b2-efix_p1-0-.root"; | |
1048 | } | |
1049 | ||
1050 | return fInputRootName; | |
1051 | } | |
1052 | ||
1053 | TString JTJTEfficiency::GetEffFullName() { | |
1054 | GetEffName(); | |
1055 | fInputRootName = fDataPath + "/" + fInputRootName; | |
1056 | return fInputRootName; | |
1057 | } | |
1058 | ||
1059 | ||
1060 | //________________________________________________________________________ | |
1061 | bool JTJTEfficiency::Load(){ | |
1062 | // Load Efficiency File based on fMode | |
1063 | if( fMode == kNotUse ) { | |
1064 | cout<<"J_WARNING : Eff Mode is \"NOTUSE\". eff is 1 !!!"<<endl; | |
1065 | return true; | |
1066 | } | |
1067 | GetEffFullName(); | |
1068 | if (TString(fInputRootName).BeginsWith("alien:")) TGrid::Connect("alien:"); | |
1069 | fInputRoot = TFile::Open( fInputRootName); | |
1070 | //fInputRoot = new TFile( fInputRootName,"READ"); | |
1071 | if( !fInputRoot ) { | |
66b055a5 | 1072 | cout << "J_ERROR : %s does not exist" << fInputRootName << endl; |
1073 | return false; | |
11476c41 | 1074 | } |
1075 | ||
1076 | //fEffDir[0] = (TDirectory*)fInputRoot->Get("EffRE"); | |
1077 | ///fEffDir[1] = (TDirectory*)fInputRoot->Get("EffMC"); | |
1078 | fEffDir[2] = (TDirectory*)fInputRoot->Get("Efficiency"); | |
1079 | //iif( fEffDir[0] && fEffDir[1] && fEffDir[2] ) | |
1080 | if( !fEffDir[2] ) | |
1081 | { | |
66b055a5 | 1082 | cout << "J_ERROR : Directory EFF is not exist"<<endl; |
1083 | return false; | |
11476c41 | 1084 | } |
1085 | ||
1086 | fCentBinAxis = (TAxis*)fEffDir[2]->Get("CentralityBin"); | |
1087 | if( !fCentBinAxis ){ | |
66b055a5 | 1088 | cout << "J_ERROR : No CentralityBin in directory" << endl; |
81fe377a | 1089 | return false; |
11476c41 | 1090 | } |
1091 | ||
1092 | ||
1093 | int nVtx = 1; | |
1094 | int nCentBin = fCentBinAxis->GetNbins(); | |
1095 | for( int ivtx=0;ivtx<nVtx;ivtx++ ){ | |
1096 | for( int icent=0;icent<nCentBin;icent++ ){ | |
1097 | for( int icut=0;icut<kJNTrackCuts;icut++ ){ | |
1098 | fCorrection[ivtx][icent][icut] | |
1099 | = (TGraphErrors*) fEffDir[2]->Get(Form("gCor%02d%02d%02d", ivtx,icent,icut)); | |
1100 | //cout<<"J_LOG : Eff graph - "<<Form("gCor%02d%02d%02d", ivtx,icent,icut)<<" - "<<g<<endl; | |
1101 | } | |
1102 | } | |
1103 | } | |
1104 | cout<<"J_LOG : Eff file is "<<fInputRootName<<endl; | |
1105 | cout<<"J_LOG : Eff Cent Bins are "; | |
1106 | for( int i=0;i<=nCentBin;i++ ){ | |
1107 | //cout<<fCentBinAxis->GetXbins()->At(i)<<" "; | |
1108 | } | |
1109 | //cout<<endl; | |
1110 | return true; | |
1111 | } | |
1112 | ||
1113 | ||
1114 |