Coverity fixes
[u/mrichter/AliRoot.git] / PWGJE / EMCALJetTasks / UserTasks / AliAnalysisTaskJetJTJT.cxx
CommitLineData
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
53ClassImp(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//________________________________________________________________________
139AliAnalysisTaskJetJTJT::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//________________________________________________________________________
219AliAnalysisTaskJetJTJT::~AliAnalysisTaskJetJTJT()
220{
221 // Destructor.
222}
223
224
225void 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
237void 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
248void 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//________________________________________________________________________
261void 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//________________________________________________________________________
516Bool_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//-----------------------------------------------------------------------
788Double_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
804Double_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//________________________________________________________________________
821void 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//________________________________________________________________________
897void 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//________________________________________________________________________
927Bool_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//________________________________________________________________________
935void AliAnalysisTaskJetJTJT::Terminate(Option_t *)
936{
937 // Called once at the end of the analysis.
938}
939
940
941
942//________________________________________________________________________
943JTJTEfficiency::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
960JTJTEfficiency::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
978JTJTEfficiency& JTJTEfficiency::operator=(const JTJTEfficiency& obj){
979 // equal sign operator TODO: content
980 JUNUSED(obj);
981 return *this;
982}
983
984
985//________________________________________________________________________
986double 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
1008TString 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 */
1018if(fPeriodStr == "LHC10b"){
1019 fInputRootName = "Eff--LHC10b-LHC10d1-0-.root";
1020}
1021if(fPeriodStr == "LHC10c"){
1022 fInputRootName = "Eff--LHC10c-LHC10d4-0-.root";
1023}
1024if(fPeriodStr == "LHC10d"){
1025 fInputRootName = "Eff--LHC10d-LHC10f6a-0-.root";
1026}
1027if(fPeriodStr == "LHC10e"){
1028 fInputRootName = "Eff--LHC10e-LHC10e20-0-.root";
1029}
1030if(fPeriodStr == "LHC10h"){
1031 fInputRootName = "Eff--LHC10h-LHC11a10a_bis-0-.root";
1032}
1033if(fPeriodStr == "LHC11a"){
1034 fInputRootName = "Eff--LHC11a-LHC11b10a-0-.root";
1035}
1036if(fPeriodStr == "LHC13b"){
1037 fInputRootName = "Eff--LHC13b-LHC13b2-efix_p1-0-.root";
1038}
1039
1040if(fPeriodStr == "LHC13c"){
1041 fInputRootName = "Eff--LHC13c-LHC13b2-efix_p1-0-.root";
1042}
1043if(fPeriodStr == "LHC13d"){
1044 fInputRootName = "Eff--LHC13d-LHC13b2-efix_p1-0-.root";
1045}
1046if(fPeriodStr == "LHC13e"){
1047 fInputRootName = "Eff--LHC13e-LHC13b2-efix_p1-0-.root";
1048}
1049
1050return fInputRootName;
1051}
1052
1053TString JTJTEfficiency::GetEffFullName() {
1054 GetEffName();
1055 fInputRootName = fDataPath + "/" + fInputRootName;
1056 return fInputRootName;
1057}
1058
1059
1060//________________________________________________________________________
1061bool 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;
11476c41 1089 gSystem->Exit(1);
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