Collision event selection.
[u/mrichter/AliRoot.git] / PWG1 / AliAnalysisTaskQASym.cxx
CommitLineData
9f9631f3 1#include "TChain.h"
2#include "TTree.h"
3#include "TH1F.h"
4#include "TH2F.h"
5#include "TCanvas.h"
6#include "TList.h"
7#include "TParticle.h"
8#include "TParticlePDG.h"
9#include "TProfile.h"
10#include "TNtuple.h"
11#include "TFile.h"
12
13#include "AliAnalysisTask.h"
14#include "AliAnalysisManager.h"
15
16#include "AliESDEvent.h"
17
18#include "AliESDVertex.h"
19#include "AliESDInputHandler.h"
20#include "AliESDtrackCuts.h"
21#include "AliMultiplicity.h"
22
23
24#include "AliAnalysisTaskQASym.h"
25#include "AliExternalTrackParam.h"
26#include "AliTrackReference.h"
27
28#include "AliHeader.h"
29#include "AliGenEventHeader.h"
30#include "AliGenDPMjetEventHeader.h"
31
32// Analysis Task for basic QA on the ESD
33
34// Authors: Jan Fiete Grosse-Oetringhaus, Christian Klein-Boesing,
35// Andreas Morsch, Eva Sicking
36
37ClassImp(AliAnalysisTaskQASym)
38
39//________________________________________________________________________
40AliAnalysisTaskQASym::AliAnalysisTaskQASym(const char *name)
41 : AliAnalysisTaskSE(name)
42 ,fFieldOn(kTRUE)
43 ,fHists(0)
44 ,fHistRECpt(0)
45 ,fEta(0)
46 ,fEtaPt(0)
47 ,fQPt(0)
48 ,fDca(0)
49 ,fqRec(0)
50 ,fsigmaPt(0)
51
52 ,fRecPtPos(0)
53 ,fRecPtNeg(0)
54 ,fRecPhiPos(0)
55 ,fRecPhiNeg(0)
56 ,fRecEtaPos(0)
57 ,fRecEtaNeg(0)
58 ,fRecEtaPtPos(0)
59 ,fRecEtaPtNeg(0)
60 ,fRecDcaPos(0)
61 ,fRecDcaNeg(0)
62 ,fRecDcaNegInv(0)
63 ,fRecDPos(0)
64 ,fRecDNeg(0)
65
66 ,fRecQPtPosEta(0)
67 ,fRecQPtNegEta(0)
68 ,fRecPtPosEta(0)
69 ,fRecPtNegEta(0)
70 ,fRecPhiPosEta(0)
71 ,fRecPhiNegEta(0)
72 ,fRecDcaPosEta(0)
73 ,fRecDcaNegEta(0)
74 ,fRecDPosEta(0)
75 ,fRecDNegEta(0)
76
77 ,fRecPtPosVz(0)
78 ,fRecPtNegVz(0)
79 ,fRecEtaPosVz(0)
80 ,fRecEtaNegVz(0)
81 ,fRecPhiPosVz(0)
82 ,fRecPhiNegVz(0)
83 ,fSignedDcaPosVz(0)
84 ,fSignedDcaNegVz(0)
85 ,fRecQPtPosEtaVz(0)
86 ,fRecQPtNegEtaVz(0)
87 ,fRecEtaPtPosVz(0)
88 ,fRecEtaPtNegVz(0)
89
90
91 ,fDeltaPhiAll(0)
92 ,fDeltaPhiLeading(0)
93 ,fDiffDcaD(0)
94 ,fPhiRec(0)
95 ,fThetaRec(0)
96 ,fNumber(0)
97 ,fVx(0)
98 ,fVy(0)
99 ,fVz(0)
100 ,test(0)
101
102
103 ,fCuts(0)
104
105
106
107{
108 // Constructor
109 for(Int_t i = 0;i<18;++i){
110 fRecPtTpcSector[i] = 0;
111 fRecEtaTpcSector[i] = 0;
112 fSignedDcaTpcSector[i] = 0;
113 fRecQPtTpcSector[i] = 0;
114 fRecEtaPtTpcSector[i] = 0;
115 }
116
117 for(Int_t i = 0;i< 7;++i){
118 fRecPtPosLadder[i] = 0;
119 fRecPtNegLadder[i] = 0;
120 fRecPhiPosLadder[i] = 0;
121 fRecPhiNegLadder[i] = 0;
122 fRecEtaPosLadder[i] = 0;
123 fRecEtaNegLadder[i] = 0;
124 fSignDcaPos[i] = 0;
125 fSignDcaNeg[i] = 0;
126 fSignDcaNegInv[i] = 0;
127 fPtSigmaPos[i] =0;
128 fPtSigmaNeg[i] =0;
129 fqPtRec[i] =0;
130 fDcaSigmaPos[i] =0;
131 fDcaSigmaNeg[i] =0;
132 }
133
134 DefineOutput(1, TList::Class());
135
136
137
138}
139
140
141//________________________________________________________________________
142void AliAnalysisTaskQASym::UserCreateOutputObjects()
143{
144 // Create histograms
145 // Called once
146
147 Double_t range = 300.;
148 Double_t pt = 20.;
149
150 fHists = new TList();
151 test = new TNtuple("test","test",
152 "pt:phi:theta:x:y:z:charge");
153 fHistRECpt = new TH1F("fHistRECpt",
154 " p_{T}",
155 100, 0., pt);
156 fEta = new TH1F("fEta",
157 " #eta",
158 200, -2., 2.);
159 fThetaRec = new TH1F("fThetaRec",
160 " #theta",
161 180, 0., TMath::Pi());
162 fPhiRec = new TH1F("fPhiRec",
163 " #phi",
164 180, 0., 2*TMath::Pi());
165 fNumber = new TH1F("fNumber",
166 "number of tracks per event",
167 50, 0.5, 49.5);
168 fVx = new TH1F("fVx",
169 "X of vertex",
170 500, -500., 500.);
171 fVy = new TH1F("fVy",
172 "Y of vertex",
173 500, -500., 500.);
174 fVz = new TH1F("fVz",
175 "Z of vertex",
176 500, -500., 500.);
177
178 fEtaPt = new TH1F("fEtaPt",
179 " #eta/p_{T} ",
180 100, -1., 1.);
181
182 fQPt = new TH1F("fQPt",
183 " charge/p_{T} ",
184 100, -1., 1.);
185
186 fDca = new TH1F("fDca",
187 " dca ",
188 200, -range, range);
189
190
191 fqRec = new TH1F("fqRec",
192 " charge all reconstructed particle",
193 21, -9.5, 10.5);
194
195 fsigmaPt = new TH1F("fsigmaPt",
196 "Log_{10}(#sigma_{p_{T}})",
197 200, -2., 8.);
198
199
200
201
202 //------------
203 for(Int_t ITSlayer_case=0;ITSlayer_case<7;ITSlayer_case++){
204
205 fSignDcaPos[ITSlayer_case] = new TH1F(Form("fSignDcaPos%d", ITSlayer_case),
206 " Signed dca",
207 200, -range, range);
208 fSignDcaPos[ITSlayer_case]->GetXaxis()->SetTitle("dca");
209 fSignDcaPos[ITSlayer_case]->GetYaxis()->SetTitle("");
210
211
212 fSignDcaNeg[ITSlayer_case] = new TH1F(Form("fSignDcaNeg%d", ITSlayer_case),
213 " Signed dcas",
214 200, -range, range);
215 fSignDcaNeg[ITSlayer_case]->GetXaxis()->SetTitle("dca");
216 fSignDcaNeg[ITSlayer_case]->GetYaxis()->SetTitle("");
217
218 fSignDcaNegInv[ITSlayer_case] = new TH1F(Form("fSignDcaNegInv%d", ITSlayer_case),
219 " inverse Signed dca ",
220 200, -range, range);
221 fSignDcaNegInv[ITSlayer_case]->GetXaxis()->SetTitle("-dca");
222 fSignDcaNegInv[ITSlayer_case]->GetYaxis()->SetTitle("");
223
224
225
226
227 fPtSigmaPos[ITSlayer_case] = new TH1F(Form("fPtSigmaPos%d", ITSlayer_case),
228 " #sigma_{pT} ",
229 208, -2., 8.);
230 fPtSigmaPos[ITSlayer_case]->GetXaxis()->SetTitle("Log_{10}(#sigma_{pT})");
231 fPtSigmaPos[ITSlayer_case]->GetYaxis()->SetTitle("");
232
233
234 fPtSigmaNeg[ITSlayer_case] = new TH1F(Form("fPtSigmaNeg%d",ITSlayer_case),
235 " #sigma_{pT}",
236 208, -2., 8.);
237 fPtSigmaNeg[ITSlayer_case]->GetXaxis()->SetTitle("Log_{10}(#sigma_{pT})");
238 fPtSigmaNeg[ITSlayer_case]->GetYaxis()->SetTitle("");
239
240
241
242
243
244 fqPtRec[ITSlayer_case] = new TH1F(Form("fqPtRec%d",ITSlayer_case),
245 "q/ p_{T}",
246 200, -100., 100.);
247 fqPtRec[ITSlayer_case]->GetXaxis()->SetTitle("q_{tr}/p_{T, tr} (GeV/c)");
248 fqPtRec[ITSlayer_case]->GetYaxis()->SetTitle("");
249
250
251
252
253
254 fDcaSigmaPos[ITSlayer_case] = new TH2F(Form("fDcaSigmaPos%d", ITSlayer_case),
255 " p_{T} shift vs #sigma_{pT} ",
256 200, -range, range,200, -4., 4. );
257 fDcaSigmaPos[ITSlayer_case]->GetXaxis()->SetTitle("signed DCA)");
258 fDcaSigmaPos[ITSlayer_case]->GetYaxis()->SetTitle("log_{10}(#sigma_{pT})");
259
260
261 fDcaSigmaNeg[ITSlayer_case] = new TH2F(Form("fDcaSigmaNeg%d", ITSlayer_case),
262 " p_{T} shift vs #sigma_{pT} ",
263 200, -range, range,200, -4., 4. );
264 fDcaSigmaNeg[ITSlayer_case]->GetXaxis()->SetTitle("signed DCA");
265 fDcaSigmaNeg[ITSlayer_case]->GetYaxis()->SetTitle("log_{10}(#sigma_{pT})");
266
267
268
269
270
271
272 // YIELDs---------- positive and negative particles
273
274 fRecPtPos = new TH1F("fRecPtPos",
275 " p_{T}",
276 100, 0.,pt);
277 fRecPtPos->GetXaxis()->SetTitle("p_{T} (GeV/c)");
278 fRecPtNeg = new TH1F("fRecPtNeg",
279 " p_{T} ",
280 100, 0., pt);
281 fRecPtNeg->GetXaxis()->SetTitle("p_{T} (GeV/c)");
282
283
284 fRecPhiPos = new TH1F("fRecPhiPos",
285 "#phi",
286 361, 0., 360.);
287 fRecPhiPos->GetXaxis()->SetTitle("#phi (deg)");
288
289 fRecPhiNeg = new TH1F("fRecPhiNeg",
290 "#phi ",
291 361, 0., 360.);
292 fRecPhiNeg->GetXaxis()->SetTitle("#phi (deg)");
293
294 fRecEtaPos = new TH1F("fRecEtaPos",
295 "#eta",
296 200, -2., 2.);
297 fRecEtaPos->GetXaxis()->SetTitle("#eta");
298
299 fRecEtaNeg = new TH1F("fRecEtaNeg",
300 "#eta",
301 200, -2., 2.);
302 fRecEtaNeg->GetXaxis()->SetTitle("#eta");
303
304 fRecEtaPtPos = new TH1F("fRecEtaPtPos",
305 "#eta/p_{T}",
306 200, -0.1, .1);
307 fRecEtaPtPos->GetXaxis()->SetTitle("#eta/p_{T}");
308
309 fRecEtaPtNeg = new TH1F("fRecEtaPtNeg",
310 "#eta/p_{T}",
311 200, -.1, .1);
312 fRecEtaPtNeg->GetXaxis()->SetTitle("#eta/p_{T}");
313
314 fRecDcaPos = new TH1F("fRecDcaPos",
315 " dca",
316 100, -range, range);
317 fRecDcaPos->GetXaxis()->SetTitle("dca (cm)");
318 fRecDcaNeg = new TH1F("fRecDcaNeg",
319 " dca",
320 100, -range, range);
321 fRecDcaNeg->GetXaxis()->SetTitle("dca (cm)");
322
323 fRecDcaNegInv = new TH1F("fRecDcaNegInv",
324 " dca",
325 100, -range, range);
326 fRecDcaNegInv->GetXaxis()->SetTitle("dca (cm)");
327
328
329 fRecDPos = new TH1F("fRecDPos",
330 " d",
331 100, -range, range);
332 fRecDPos->GetXaxis()->SetTitle("d (cm)");
333 fRecDNeg = new TH1F("fRecDNeg",
334 "d",
335 100, -range, range);
336 fRecDNeg->GetXaxis()->SetTitle("d (cm)");
337
338
339 // YIELDs ---------------- positive and negative eta
340
341
342 fRecQPtPosEta = new TH1F("fRecQPtPosEta",
343 "q/p_{T}",
344 200, -0.5, 0.5);
345 fRecQPtPosEta->GetXaxis()->SetTitle("q/p_{T} ");
346
347 fRecQPtNegEta = new TH1F("fRecQPtNegEta",
348 "q/p_{T}",
349 200, -0.5, 0.5);
350 fRecQPtNegEta->GetXaxis()->SetTitle("q/p_{T}");
351
352 fRecPtPosEta = new TH1F("fRecPtPosEta",
353 " p_{T} ",
354 100, 0., pt);
355 fRecPtPosEta->GetXaxis()->SetTitle("p_{T} (GeV/c)");
356
357 fRecPtNegEta = new TH1F("fRecPtNegEta",
358 " p_{T}",
359 100, 0., pt);
360 fRecPtNegEta->GetXaxis()->SetTitle("p_{T} (GeV/c)");
361
362 fRecPhiPosEta = new TH1F("fRecPhiPosEta",
363 "#phi",
364 361, 0., 360);
365 fRecPhiPosEta->GetXaxis()->SetTitle("#phi (deg)");
366
367 fRecPhiNegEta = new TH1F("fRecPhiNegEta",
368 "#phi ",
369 361, 0, 360);
370 fRecPhiNegEta->GetXaxis()->SetTitle("#phi (deg)");
371
372 fRecDcaPosEta = new TH1F("fRecDcaPosEta",
373 " dca ",
374 100, -range, range);
375 fRecDcaPosEta->GetXaxis()->SetTitle("dca (cm)");
376 fRecDcaNegEta = new TH1F("fRecDcaNegEta",
377 " dca",
378 100, -range, range);
379 fRecDcaNegEta->GetXaxis()->SetTitle("dca (cm)");
380
381 fRecDPosEta = new TH1F("fRecDPosEta",
382 " d",
383 100, -range, range);
384 fRecDPosEta->GetXaxis()->SetTitle("d (cm)");
385 fRecDNegEta = new TH1F("fRecDNegEta",
386 "d",
387 100, -5., 5.);
388 fRecDNegEta->GetXaxis()->SetTitle("d (cm)");
389
390
391
392 // YIELDs ---------------- for TPC sectors
393 for(Int_t sector=0; sector<18;sector++){
394
395
396 fRecPtTpcSector[sector] = new TH1F(Form("fRecPtTpcSector%02d",sector),
397 Form("p_{T} distribution: TPC sector %d",
398 sector),100, 0., pt);
399 fRecPtTpcSector[sector]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
400
401 fRecEtaTpcSector[sector] = new TH1F(Form("fRecEtaTpcSector%02d",sector),
402 Form("#eta distribution: TPC sector %d",
403 sector),200, -2., 2.);
404 fRecEtaTpcSector[sector]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
405
406
407 fSignedDcaTpcSector[sector] = new TH1F(Form("fSignedDcaTpcSector%02d",sector),
408 Form("dca distribution: TPC sector %d",
409 sector),200, -range, range );
410 fSignedDcaTpcSector[sector]->GetXaxis()->SetTitle("dca");
411
412 fRecQPtTpcSector[sector] = new TH1F(Form("fRecQPtTpcSector%02d",sector),
413 Form("Q/ p_{T} distribution: TPC sector %d",
414 sector),100, -1., 1.);
415 fRecQPtTpcSector[sector]->GetXaxis()->SetTitle("Q/p_{T} (GeV/c)");
416
417 fRecEtaPtTpcSector[sector] = new TH1F(Form("fRecEtaPtTpcSector%02d",sector),
418 Form("#eta/ p_{T} distribution: TPC sector %d",
419 sector),100, -1., 1.);
420 fRecEtaPtTpcSector[sector]->GetXaxis()->SetTitle("#eta/p_{T} (GeV/c)");
421
422 }
423 // YIELDS ITS ladder
424 for(Int_t i=0;i<7;i++){
425 fRecPtPosLadder[i] = new TH1F(Form("fRecPtPosLadder%d", i),
426 " p_{T} distribution",
427 100, 0., pt);
428 fRecPtPosLadder[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
429 fRecPtNegLadder[i] = new TH1F(Form("fRecPtNegLadder%d",i),
430 " p_{T} distribution ",
431 100, 0., pt);
432 fRecPtNegLadder[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
433
434
435 fRecPhiPosLadder[i] = new TH1F(Form("fRecPhiPosLadder%d",i),
436 "#phi distribution: all pos eta",
437 361, 0., 360);
438 fRecPhiPosLadder[i]->GetXaxis()->SetTitle("#phi (deg)");
439
440 fRecPhiNegLadder[i] = new TH1F(Form("fRecPhiNegLadder%d", i),
441 "#phi distribution: all neg eta",
442 361, 0, 360);
443 fRecPhiNegLadder[i]->GetXaxis()->SetTitle("#phi (deg)");
444
445
446
447 fRecEtaPosLadder[i] = new TH1F(Form("fRecEtaPosLadder%d",i),
448 "#eta distribution",
449 200, -2., 2.);
450 fRecEtaPosLadder[i]->GetXaxis()->SetTitle("#eta)");
451
452 fRecEtaNegLadder[i] = new TH1F(Form("fRecEtaNegLadder%d", i),
453 "#eta distribution",
454 200, -2., 2.);
455 fRecEtaNegLadder[i]->GetXaxis()->SetTitle("#eta");
456 }
457
458 Double_t vzmax = 30;
459
460 fRecPtPosVz = new TH2F("fRecPtPosVz",
461 "p_{T} distribution vs Vz()",
462 100, -1., 2., 200,-vzmax,vzmax);
463 fRecPtPosVz->GetXaxis()->SetTitle("log_{10}(p_{T})");
464
465 fRecPtNegVz = new TH2F("fRecPtNegVz",
466 "p_{T} distribution vs Vz()",
467 100, -1., 2.,200,-vzmax,vzmax);
468 fRecPtNegVz->GetXaxis()->SetTitle("Log_{10}(p_{T})");
469
470
471 fRecEtaPosVz= new TH2F("fRecEtaPosVz",
472 "#eta distribution vs Vz()",
473 100, -2., 2., 200,-vzmax,vzmax);
474 fRecEtaPosVz->GetXaxis()->SetTitle("#eta");
475 fRecEtaNegVz = new TH2F("fRecEtaNegVz",
476 "#eta distribution vs Vz()",
477 100, -2., 2.,200,-vzmax,vzmax);
478 fRecEtaNegVz->GetXaxis()->SetTitle("#eta");
479
480 fRecPhiPosVz= new TH2F("fRecPhiPosVz",
481 "#eta distribution vs Vz()",
482 361, 0., 360., 200,-vzmax,vzmax);
483 fRecPhiPosVz->GetXaxis()->SetTitle("#phi (deg)");
484 fRecPhiNegVz = new TH2F("fRecPhiNegVz",
485 "dca vs Vz()",
486 361, 0., 360.,200,-vzmax,vzmax);
487 fRecPhiNegVz->GetXaxis()->SetTitle("#phi (deg)");
488
489 fSignedDcaPosVz= new TH2F("fSignedDcaPosVz",
490 "#eta distribution vs Vz()",
491 200, -range, range, 200,-vzmax,vzmax);
492 fSignedDcaPosVz->GetXaxis()->SetTitle("dca (cm)");
493 fSignedDcaNegVz = new TH2F("fSignedDcaNegVz",
494 "dca vs Vz()",
495 200, -range, range,200,-vzmax,vzmax);
496 fSignedDcaNegVz->GetXaxis()->SetTitle("dca (cm)");
497
498 fRecQPtPosEtaVz= new TH2F("fRecQPtPosEtaVz",
499 " Q/p_{T} distribution vs Vz()",
500 100, -1., 1., 200,-vzmax,vzmax);
501 fRecQPtPosEtaVz->GetXaxis()->SetTitle("Q/p_{T}");
502 fRecQPtNegEtaVz = new TH2F("fRecQPtNegEtaVz",
503 " Q/p_{T} distribution vs Vz()",
504 100, -1., 1.,200,-vzmax,vzmax);
505 fRecQPtNegEtaVz->GetXaxis()->SetTitle("Q/p_{T}");
506
507
508 fRecEtaPtPosVz= new TH2F("fRecEtaPtPosVz",
509 " #eta/p_{T} distribution vs Vz()",
510 100, -1., 1., 200,-vzmax,vzmax);
511 fRecEtaPtPosVz->GetXaxis()->SetTitle("#eta/p_{T");
512 fRecEtaPtNegVz = new TH2F("fRecEtaPtNegVz",
513 " #eta/p_{T} distribution vs Vz()",
514 100, -1., 1.,200,-vzmax,vzmax);
515 fRecEtaPtNegVz->GetXaxis()->SetTitle("#eta/p_{T}");
516
517 //new
518 fDeltaPhiAll = new TH1F("fDeltaPhiAll",
519 " #Delta #phi",200,-360,360);
520 fDeltaPhiAll->GetXaxis()->SetTitle("#Delta #phi");
521
522
523 fDeltaPhiLeading = new TH2F("fDeltaPhiLeading",
524 " #Delta #phi",361,-360,360, 361,0, 360);
525 fDeltaPhiLeading->GetXaxis()->SetTitle("#Delta #phi (deg.)");
526 fDeltaPhiLeading->GetYaxis()->SetTitle("#phi_{leading particle} (deg.)");
527
528 fDiffDcaD = new TH1F("fDiffDcaD",
529 "dca-d",
530 200, -5., 5.);
531
532 }
533
534 fHists->SetOwner();
535
536 fHists->Add(test);
537 fHists->Add(fHistRECpt);
538 fHists->Add(fEta);
539 fHists->Add(fThetaRec);
540 fHists->Add(fPhiRec);
541 fHists->Add(fNumber);
542 fHists->Add(fVx);
543 fHists->Add(fVy);
544 fHists->Add(fVz);
545
546 fHists->Add(fEtaPt);
547 fHists->Add(fQPt);
548 fHists->Add(fDca);
549
550 fHists->Add(fDeltaPhiAll);
551 fHists->Add(fDeltaPhiLeading);
552 fHists->Add(fDiffDcaD);
553
554 fHists->Add(fqRec);
555 fHists->Add(fsigmaPt);
556
557 fHists->Add(fRecPtPos);
558 fHists->Add(fRecPtNeg);
559 fHists->Add(fRecPhiPos);
560 fHists->Add(fRecPhiNeg);
561 fHists->Add(fRecEtaPos);
562 fHists->Add(fRecEtaNeg);
563 fHists->Add(fRecEtaPtPos);
564 fHists->Add(fRecEtaPtNeg);
565 fHists->Add(fRecDcaPos);
566 fHists->Add(fRecDcaNeg);
567 fHists->Add(fRecDcaNegInv);
568 fHists->Add(fRecDPos);
569 fHists->Add(fRecDNeg);
570
571
572 fHists->Add(fRecQPtPosEta);
573 fHists->Add(fRecQPtNegEta);
574 fHists->Add(fRecPtPosEta);
575 fHists->Add(fRecPtNegEta);
576 fHists->Add(fRecPhiPosEta);
577 fHists->Add(fRecPhiNegEta);
578 fHists->Add(fRecDcaPosEta);
579 fHists->Add(fRecDcaNegEta);
580 fHists->Add(fRecDPosEta);
581 fHists->Add(fRecDNegEta);
582
583
584 for(Int_t i=0;i<18;i++){
585 fHists->Add(fRecPtTpcSector[i]);
586 fHists->Add(fRecEtaTpcSector[i]);
587 fHists->Add(fSignedDcaTpcSector[i]);
588 fHists->Add(fRecQPtTpcSector[i]);
589 fHists->Add(fRecEtaPtTpcSector[i]);
590 }
591
592 for(Int_t i=0;i<7;i++){
593 fHists->Add(fRecPtPosLadder[i]);
594 fHists->Add(fRecPtNegLadder[i]);
595 fHists->Add(fRecPhiPosLadder[i]);
596 fHists->Add(fRecPhiNegLadder[i]);
597 fHists->Add(fRecEtaPosLadder[i]);
598 fHists->Add(fRecEtaNegLadder[i]);
599 }
600
601 fHists->Add(fRecPtPosVz);
602 fHists->Add(fRecPtNegVz);
603 fHists->Add(fRecEtaPosVz);
604 fHists->Add(fRecEtaNegVz);
605 fHists->Add(fRecPhiPosVz);
606 fHists->Add(fRecPhiNegVz);
607 fHists->Add(fSignedDcaPosVz);
608 fHists->Add(fSignedDcaNegVz);
609 fHists->Add(fRecQPtPosEtaVz);
610 fHists->Add(fRecQPtNegEtaVz);
611 fHists->Add(fRecEtaPtPosVz);
612 fHists->Add(fRecEtaPtNegVz);
613
614
615 for(Int_t i=0;i<7;i++){
616 fHists->Add(fSignDcaPos[i]);
617 fHists->Add(fSignDcaNeg[i]);
618 fHists->Add(fSignDcaNegInv[i]);
619
620 fHists->Add(fPtSigmaPos[i]);
621 fHists->Add(fPtSigmaNeg[i]);
622 fHists->Add(fqPtRec[i]);
623
624 fHists->Add(fDcaSigmaPos[i]);
625 fHists->Add(fDcaSigmaNeg[i]);
626
627
628 }
629
630
631
632 for (Int_t i=0; i<fHists->GetEntries(); ++i) {
633 TH1 *h1 = dynamic_cast<TH1*>(fHists->At(i));
634 if (h1){
635 // Printf("%s ",h1->GetName());
636 h1->Sumw2();
637 }
638 }
639 // BKC
640
641
642}
643
644//__________________________________________________________
645
646void AliAnalysisTaskQASym::UserExec(Option_t *)
647{
648 AliVEvent *event = InputEvent();
649 if (!event) {
650 Printf("ERROR: Could not retrieve event");
651 return;
652 }
653
654
655 if(Entry()==0){
656 AliESDEvent* esd = dynamic_cast<AliESDEvent*>(event);
657 if(esd){
658 Printf("We are reading from ESD");
659 }
660
661 }
662
663 Printf("There are %d tracks in this event", event->GetNumberOfTracks());
664
665
666 Int_t LeadingTrack = 0;
667 Float_t LeadingEnergy = -20.;
668 Float_t LeadingPhi = 0;//TMath::Pi();
669
670
671 if(event->GetNumberOfTracks()!=0) fNumber->Fill(event->GetNumberOfTracks());
672
673 for (Int_t iTrack = 0; iTrack < event->GetNumberOfTracks(); iTrack++) {
674
675
676 AliVParticle *track = event->GetTrack(iTrack);
677 AliESDtrack *esdtrack = dynamic_cast<AliESDtrack*>(track);
678 if (!track) {
679 Printf("ERROR: Could not receive track %d", iTrack);
680 continue;
681 }
682
683 if (!fCuts->AcceptTrack(esdtrack)) continue;
684 const AliExternalTrackParam * tpcP = esdtrack->GetTPCInnerParam();
685 if (!tpcP) continue;
686
687
688
689 if(tpcP->E()>LeadingEnergy){
690 LeadingTrack=iTrack;
691 LeadingEnergy=tpcP->E();
692 LeadingPhi=tpcP->Phi();
693 }
694
695
696
697
698 //propagate to dca
699 esdtrack->PropagateToDCA(event->GetPrimaryVertex(),
700 event->GetMagneticField(), 10000.);
701
702
703 // if(tpcP->Pt()<2.)continue;
704
705
706 fqRec->Fill(tpcP->Charge());
707
708
709 Double_t sigmapt = tpcP->GetSigma1Pt2();
710 sigmapt= sqrt(sigmapt);
711 sigmapt= sigmapt *(tpcP->Pt()*tpcP->Pt());
712
713 if(sigmapt == 0.)continue;
714 fsigmaPt->Fill(TMath::Log10(sigmapt));
715
716
717
718
719
720
721
722 // hits in ITS layer
723 Int_t cas=-1;
724 if(esdtrack->HasPointOnITSLayer(0))
725 cas=0;
726 else if(!esdtrack->HasPointOnITSLayer(0)
727 && esdtrack->HasPointOnITSLayer(1))
728 cas=1;
729 else if(!esdtrack->HasPointOnITSLayer(0)
730 && !esdtrack->HasPointOnITSLayer(1)
731 && esdtrack->HasPointOnITSLayer(2))
732 cas=2;
733 else if(!esdtrack->HasPointOnITSLayer(0)
734 && !esdtrack->HasPointOnITSLayer(1)
735 && !esdtrack->HasPointOnITSLayer(2)
736 && esdtrack->HasPointOnITSLayer(3))
737 cas=3;
738 else if(!esdtrack->HasPointOnITSLayer(0)
739 && !esdtrack->HasPointOnITSLayer(1)
740 && !esdtrack->HasPointOnITSLayer(2)
741 && !esdtrack->HasPointOnITSLayer(3)
742 && esdtrack->HasPointOnITSLayer(4))
743 cas=4;
744 else if( !esdtrack->HasPointOnITSLayer(0)
745 && !esdtrack->HasPointOnITSLayer(1)
746 && !esdtrack->HasPointOnITSLayer(2)
747 && !esdtrack->HasPointOnITSLayer(3)
748 && !esdtrack->HasPointOnITSLayer(4)
749 && esdtrack->HasPointOnITSLayer(5))
750 cas=5;
751 else
752 cas=6;
753
754
755
756 //-------------------
757
758 sdca_tr = (tpcP->Py()*tpcP->Xv()
759 - tpcP->Px()*tpcP->Yv())/tpcP->Pt();
760
761
762
763 fqPtRec[cas]->Fill(tpcP->Charge()/tpcP->Pt());
764
765
766 test->Fill(tpcP->Pt(),tpcP->Phi(), tpcP->Theta(),
767 tpcP->Xv(),tpcP->Yv(), tpcP->Zv(), tpcP->Charge());
768 fHistRECpt->Fill(tpcP->Pt());
769 fEta->Fill(tpcP->Eta());
770 fThetaRec->Fill(tpcP->Theta());
771 fPhiRec->Fill(tpcP->Phi());
772 fVx->Fill(tpcP->Xv());
773 fVy->Fill(tpcP->Yv());
774 fVz->Fill(tpcP->Zv());
775
776
777 fEtaPt->Fill(tpcP->Eta()/tpcP->Pt());
778 fQPt->Fill(tpcP->Charge()/tpcP->Pt());
779 fDca->Fill(sdca_tr);
780
781
782
783
784 esdtrack->GetImpactParameters(xy,z);
785 fDiffDcaD->Fill(sdca_tr+xy);
786
787
788 //for positive particles
789 if(tpcP->Charge()>0){
790 fRecPtPos->Fill(tpcP->Pt());
791 fRecPtPosLadder[cas]->Fill(tpcP->Pt());
792 fRecPtPosVz->Fill(TMath::Log10(tpcP->Pt()),tpcP->Zv());
793 fRecPhiPos->Fill(TMath::RadToDeg()*tpcP->Phi());
794
795
796 fRecPhiPosLadder[cas]->Fill(TMath::RadToDeg()*tpcP->Phi());
797 fRecPhiPosVz->Fill(TMath::RadToDeg()*tpcP->Phi(),tpcP->Zv());
798 fSignedDcaPosVz->Fill(sdca_tr,tpcP->Zv());
799
800 fRecEtaPos->Fill(tpcP->Eta());
801 fRecEtaPosLadder[cas]->Fill(tpcP->Eta());
802 fRecEtaPtPos->Fill(tpcP->Eta()/tpcP->Pt());
803 fRecEtaPosVz->Fill(tpcP->Eta(),tpcP->Zv());
804 fRecEtaPtPosVz->Fill(tpcP->Eta()/tpcP->Pt(),tpcP->Zv());
805
806 fRecDcaPos->Fill(sdca_tr);
807 fRecDPos->Fill(xy);
808 fSignDcaPos[cas]->Fill(sdca_tr);
809
810
811 fDcaSigmaPos[cas]->Fill(sdca_tr, TMath::Log10(sigmapt));
812
813 fPtSigmaPos[cas]->Fill(TMath::Log10(sigmapt));
814 }
815 //and negative particles
816 else {
817 fRecPtNeg->Fill(tpcP->Pt());
818 fRecPtNegLadder[cas]->Fill(tpcP->Pt());
819 fRecPtNegVz->Fill(TMath::Log10(tpcP->Pt()),tpcP->Zv());
820
821 fRecPhiNeg->Fill(TMath::RadToDeg()*tpcP->Phi());
822 fRecPhiNegLadder[cas]->Fill(TMath::RadToDeg()*tpcP->Phi());
823 fRecPhiNegVz->Fill(TMath::RadToDeg()*tpcP->Phi(),tpcP->Zv());
824 fSignedDcaNegVz->Fill(sdca_tr,tpcP->Zv());
825 fRecEtaPtNegVz->Fill(tpcP->Eta()/tpcP->Pt(),tpcP->Zv());
826
827 fRecEtaNeg->Fill(tpcP->Eta());
828 fRecEtaNegLadder[cas]->Fill(tpcP->Eta());
829 fRecEtaPtNeg->Fill(tpcP->Eta()/tpcP->Pt());
830 fRecEtaNegVz->Fill(tpcP->Eta(),tpcP->Zv());
831
832 fRecDcaNeg->Fill(sdca_tr);
833 fRecDcaNegInv->Fill(-sdca_tr);
834 fRecDNeg->Fill(xy);
835 fSignDcaNeg[cas]->Fill(sdca_tr);
836 fSignDcaNegInv[cas]->Fill(-sdca_tr);
837
838
839 fDcaSigmaNeg[cas]->Fill(sdca_tr,TMath::Log10(sigmapt));
840
841 fPtSigmaNeg[cas]->Fill(TMath::Log10(sigmapt));
842 }
843
844
845
846 //all particles with positive eta
847 if(tpcP->Eta()>0){
848 fRecQPtPosEta->Fill(tpcP->Charge()/tpcP->Pt());
849 fRecPtPosEta->Fill(tpcP->Pt());
850 fRecPhiPosEta->Fill(TMath::RadToDeg()*tpcP->Phi());
851 fRecQPtPosEtaVz->Fill(tpcP->Charge()/tpcP->Pt(),tpcP->Zv());
852 fRecDcaPosEta->Fill(sdca_tr);
853 fRecDPosEta->Fill(xy);
854 }
855 //all particles with negative eta (and eta==0)
856 else{
857 fRecQPtNegEta->Fill(tpcP->Charge()/tpcP->Pt());
858 fRecPtNegEta->Fill(tpcP->Pt());
859 fRecPhiNegEta->Fill(TMath::RadToDeg()*tpcP->Phi());
860 fRecQPtNegEtaVz->Fill(tpcP->Charge()/tpcP->Pt(),tpcP->Zv());
861 fRecDcaNegEta->Fill(sdca_tr);
862 fRecDNegEta->Fill(xy);
863
864 }
865
866
867 //spectren detected by TPC sectors
868 //pt cut on 1 GeV/c ?!
869 // if(tpcP->Pt()<1.) continue;
870 fRecPtTpcSector[Int_t(tpcP->Phi()*
871 TMath::RadToDeg()/20)]->Fill(tpcP->Pt());
872 fRecEtaTpcSector[Int_t(tpcP->Phi()*
873 TMath::RadToDeg()/20)]->Fill(tpcP->Eta());
874 fSignedDcaTpcSector[Int_t(tpcP->Phi()*
875 TMath::RadToDeg()/20)]->Fill(sdca_tr);
876 fRecQPtTpcSector[Int_t(tpcP->Phi()*
877 TMath::RadToDeg()/20)]->Fill(tpcP->Charge()/tpcP->Pt());
878 fRecEtaPtTpcSector[Int_t(tpcP->Phi()*
879 TMath::RadToDeg()/20)]->Fill(tpcP->Eta()/tpcP->Pt());
880
881
882
883
884
885
886 // another track loop
887 for (Int_t iTrack2 = 0; iTrack2 < event->GetNumberOfTracks(); iTrack2++) {
888
889 if(LeadingTrack==iTrack2) continue;
890
891 AliVParticle *track2 = event->GetTrack(iTrack2);
892 AliESDtrack* esdtrack2 = dynamic_cast<AliESDtrack*>(track2);
893 if (!track2) {
894 Printf("ERROR: Could not receive track %d", iTrack);
895 continue;
896 }
897 if (!fCuts->AcceptTrack(esdtrack2)) continue;
898 //propagate to dca
899 esdtrack2->PropagateToDCA(event->GetPrimaryVertex(),
900 event->GetMagneticField(), 10000.);
901
902 fDeltaPhiLeading->Fill((LeadingPhi-esdtrack2->Phi())*TMath::RadToDeg(),
903 LeadingPhi*TMath::RadToDeg() );
904
905
906
907 }//second track loop
908 }//first track loop
909
910
911
912
913 // Post output data.
914 // PostData(1, fHistPt);
915 PostData(1, fHists);
916}
917
918
919
920
921
922//________________________________________________________________________
923void AliAnalysisTaskQASym::Terminate(Option_t *)
924{
925
926
927}
928
929
930
931
932