Adding comments
[u/mrichter/AliRoot.git] / PWG0 / esdTrackCuts / AliESDtrackCuts.cxx
CommitLineData
0b75bef2 1#include "AliESDtrackCuts.h"
2
aa93c798 3
4#include <AliESDtrack.h>
5#include <AliESD.h>
6#include <AliLog.h>
5c495d37 7#include <TTree.h>
aa93c798 8
0b75bef2 9//____________________________________________________________________
f58f1a93 10ClassImp(AliESDtrackCuts)
0b75bef2 11
ebd60d8e 12// Cut names
ab2c1f0d 13const Char_t* AliESDtrackCuts::fgkCutNames[kNCuts] = {
ebd60d8e 14 "require TPC refit",
15 "require ITS refit",
16 "n clusters TPC",
17 "n clusters ITS",
18 "#Chi^{2}/clusters TPC",
19 "#Chi^{2}/clusters ITS",
20 "cov 11",
21 "cov 22",
22 "cov 33",
23 "cov 44",
24 "cov 55",
25 "trk-to-vtx",
26 "trk-to-vtx failed",
27 "kink daughters",
ebd60d8e 28 "p",
29 "p_{T}",
30 "p_{x}",
31 "p_{y}",
32 "p_{z}",
33 "y",
34 "eta"
35};
36
0b75bef2 37//____________________________________________________________________
d8594d03 38AliESDtrackCuts::AliESDtrackCuts() : TNamed(),
39 fCutMinNClusterTPC(0),
40 fCutMinNClusterITS(0),
41 fCutMaxChi2PerClusterTPC(0),
42 fCutMaxChi2PerClusterITS(0),
43 fCutMaxC11(0),
44 fCutMaxC22(0),
45 fCutMaxC33(0),
46 fCutMaxC44(0),
47 fCutMaxC55(0),
48 fCutAcceptKinkDaughters(0),
49 fCutRequireTPCRefit(0),
50 fCutRequireITSRefit(0),
51 fCutNsigmaToVertex(0),
52 fCutSigmaToVertexRequired(0),
53 fPMin(0),
54 fPMax(0),
55 fPtMin(0),
56 fPtMax(0),
57 fPxMin(0),
58 fPxMax(0),
59 fPyMin(0),
60 fPyMax(0),
61 fPzMin(0),
62 fPzMax(0),
63 fEtaMin(0),
64 fEtaMax(0),
65 fRapMin(0),
66 fRapMax(0),
67 fHistogramsOn(0),
68 fhCutStatistics(0),
69 fhCutCorrelation(0)
70{
71 //
72 // default constructor
73 //
0b8a2421 74
75 Init();
d8594d03 76}
77
78//____________________________________________________________________
79AliESDtrackCuts::AliESDtrackCuts(Char_t* name, Char_t* title) : TNamed(name,title),
aa617804 80 fCutMinNClusterTPC(0),
81 fCutMinNClusterITS(0),
82 fCutMaxChi2PerClusterTPC(0),
83 fCutMaxChi2PerClusterITS(0),
84 fCutMaxC11(0),
85 fCutMaxC22(0),
86 fCutMaxC33(0),
87 fCutMaxC44(0),
88 fCutMaxC55(0),
89 fCutAcceptKinkDaughters(0),
90 fCutRequireTPCRefit(0),
91 fCutRequireITSRefit(0),
92 fCutNsigmaToVertex(0),
93 fCutSigmaToVertexRequired(0),
94 fPMin(0),
95 fPMax(0),
96 fPtMin(0),
97 fPtMax(0),
98 fPxMin(0),
99 fPxMax(0),
100 fPyMin(0),
101 fPyMax(0),
102 fPzMin(0),
103 fPzMax(0),
104 fEtaMin(0),
105 fEtaMax(0),
106 fRapMin(0),
107 fRapMax(0),
108 fHistogramsOn(0),
109 fhCutStatistics(0),
110 fhCutCorrelation(0)
ebd60d8e 111{
112 //
113 // constructor
114 //
115
116 Init();
0b75bef2 117
118 //##############################################################################
119 // setting default cuts
0b75bef2 120 SetMinNClustersTPC();
ebd60d8e 121 SetMinNClustersITS();
0b75bef2 122 SetMaxChi2PerClusterTPC();
123 SetMaxChi2PerClusterITS();
124 SetMaxCovDiagonalElements();
125 SetRequireTPCRefit();
126 SetRequireITSRefit();
127 SetAcceptKingDaughters();
128 SetMinNsigmaToVertex();
129 SetRequireSigmaToVertex();
130 SetPRange();
131 SetPtRange();
132 SetPxRange();
133 SetPyRange();
134 SetPzRange();
135 SetEtaRange();
136 SetRapRange();
137
138 SetHistogramsOn();
ebd60d8e 139}
140
141//_____________________________________________________________________________
3633478d 142AliESDtrackCuts::AliESDtrackCuts(const AliESDtrackCuts &c) : TNamed(c),
aa617804 143 fCutMinNClusterTPC(0),
144 fCutMinNClusterITS(0),
145 fCutMaxChi2PerClusterTPC(0),
146 fCutMaxChi2PerClusterITS(0),
147 fCutMaxC11(0),
148 fCutMaxC22(0),
149 fCutMaxC33(0),
150 fCutMaxC44(0),
151 fCutMaxC55(0),
152 fCutAcceptKinkDaughters(0),
153 fCutRequireTPCRefit(0),
154 fCutRequireITSRefit(0),
155 fCutNsigmaToVertex(0),
156 fCutSigmaToVertexRequired(0),
157 fPMin(0),
158 fPMax(0),
159 fPtMin(0),
160 fPtMax(0),
161 fPxMin(0),
162 fPxMax(0),
163 fPyMin(0),
164 fPyMax(0),
165 fPzMin(0),
166 fPzMax(0),
167 fEtaMin(0),
168 fEtaMax(0),
169 fRapMin(0),
170 fRapMax(0),
171 fHistogramsOn(0),
172 fhCutStatistics(0),
173 fhCutCorrelation(0)
ebd60d8e 174{
175 //
176 // copy constructor
177 //
178
179 ((AliESDtrackCuts &) c).Copy(*this);
180}
181
182AliESDtrackCuts::~AliESDtrackCuts()
183{
184 //
185 // destructor
186 //
187
927c71e5 188 for (Int_t i=0; i<2; i++) {
189
190 if (fhNClustersITS[i])
191 delete fhNClustersITS[i];
192 if (fhNClustersTPC[i])
193 delete fhNClustersTPC[i];
194 if (fhChi2PerClusterITS[i])
195 delete fhChi2PerClusterITS[i];
196 if (fhChi2PerClusterTPC[i])
197 delete fhChi2PerClusterTPC[i];
198 if (fhC11[i])
199 delete fhC11[i];
200 if (fhC22[i])
201 delete fhC22[i];
202 if (fhC33[i])
203 delete fhC33[i];
204 if (fhC44[i])
205 delete fhC44[i];
206 if (fhC55[i])
207 delete fhC55[i];
208
209 if (fhDXY[i])
210 delete fhDXY[i];
211 if (fhDZ[i])
212 delete fhDZ[i];
213 if (fhDXYvsDZ[i])
214 delete fhDXYvsDZ[i];
215
216 if (fhDXYNormalized[i])
217 delete fhDXYNormalized[i];
218 if (fhDZNormalized[i])
219 delete fhDZNormalized[i];
220 if (fhDXYvsDZNormalized[i])
221 delete fhDXYvsDZNormalized[i];
222 }
223
224 if (fhCutStatistics)
225 delete fhCutStatistics;
226 if (fhCutCorrelation)
227 delete fhCutCorrelation;
ebd60d8e 228}
229
230void AliESDtrackCuts::Init()
231{
232 //
233 // sets everything to zero
234 //
235
ab2c1f0d 236 fCutMinNClusterTPC = 0;
237 fCutMinNClusterITS = 0;
ebd60d8e 238
ab2c1f0d 239 fCutMaxChi2PerClusterTPC = 0;
240 fCutMaxChi2PerClusterITS = 0;
ebd60d8e 241
ab2c1f0d 242 fCutMaxC11 = 0;
243 fCutMaxC22 = 0;
244 fCutMaxC33 = 0;
245 fCutMaxC44 = 0;
246 fCutMaxC55 = 0;
ebd60d8e 247
ab2c1f0d 248 fCutAcceptKinkDaughters = 0;
249 fCutRequireTPCRefit = 0;
250 fCutRequireITSRefit = 0;
ebd60d8e 251
ab2c1f0d 252 fCutNsigmaToVertex = 0;
253 fCutSigmaToVertexRequired = 0;
ebd60d8e 254
255 fPMin = 0;
256 fPMax = 0;
257 fPtMin = 0;
258 fPtMax = 0;
259 fPxMin = 0;
260 fPxMax = 0;
261 fPyMin = 0;
262 fPyMax = 0;
263 fPzMin = 0;
264 fPzMax = 0;
265 fEtaMin = 0;
266 fEtaMax = 0;
267 fRapMin = 0;
268 fRapMax = 0;
269
270 fHistogramsOn = kFALSE;
271
272 for (Int_t i=0; i<2; ++i)
273 {
274 fhNClustersITS[i] = 0;
275 fhNClustersTPC[i] = 0;
276
277 fhChi2PerClusterITS[i] = 0;
278 fhChi2PerClusterTPC[i] = 0;
279
280 fhC11[i] = 0;
281 fhC22[i] = 0;
282 fhC33[i] = 0;
283 fhC44[i] = 0;
284 fhC55[i] = 0;
285
286 fhDXY[i] = 0;
287 fhDZ[i] = 0;
288 fhDXYvsDZ[i] = 0;
289
290 fhDXYNormalized[i] = 0;
291 fhDZNormalized[i] = 0;
292 fhDXYvsDZNormalized[i] = 0;
293 }
294
295 fhCutStatistics = 0;
296 fhCutCorrelation = 0;
297}
298
299//_____________________________________________________________________________
300AliESDtrackCuts &AliESDtrackCuts::operator=(const AliESDtrackCuts &c)
301{
302 //
303 // Assignment operator
304 //
0b75bef2 305
ebd60d8e 306 if (this != &c) ((AliESDtrackCuts &) c).Copy(*this);
307 return *this;
308}
309
310//_____________________________________________________________________________
311void AliESDtrackCuts::Copy(TObject &c) const
312{
313 //
314 // Copy function
315 //
316
317 AliESDtrackCuts& target = (AliESDtrackCuts &) c;
318
319 target.Init();
320
ab2c1f0d 321 target.fCutMinNClusterTPC = fCutMinNClusterTPC;
322 target.fCutMinNClusterITS = fCutMinNClusterITS;
ebd60d8e 323
ab2c1f0d 324 target.fCutMaxChi2PerClusterTPC = fCutMaxChi2PerClusterTPC;
325 target.fCutMaxChi2PerClusterITS = fCutMaxChi2PerClusterITS;
ebd60d8e 326
ab2c1f0d 327 target.fCutMaxC11 = fCutMaxC11;
328 target.fCutMaxC22 = fCutMaxC22;
329 target.fCutMaxC33 = fCutMaxC33;
330 target.fCutMaxC44 = fCutMaxC44;
331 target.fCutMaxC55 = fCutMaxC55;
ebd60d8e 332
ab2c1f0d 333 target.fCutAcceptKinkDaughters = fCutAcceptKinkDaughters;
334 target.fCutRequireTPCRefit = fCutRequireTPCRefit;
335 target.fCutRequireITSRefit = fCutRequireITSRefit;
ebd60d8e 336
ab2c1f0d 337 target.fCutNsigmaToVertex = fCutNsigmaToVertex;
338 target.fCutSigmaToVertexRequired = fCutSigmaToVertexRequired;
ebd60d8e 339
340 target.fPMin = fPMin;
341 target.fPMax = fPMax;
342 target.fPtMin = fPtMin;
343 target.fPtMax = fPtMax;
344 target.fPxMin = fPxMin;
345 target.fPxMax = fPxMax;
346 target.fPyMin = fPyMin;
347 target.fPyMax = fPyMax;
348 target.fPzMin = fPzMin;
349 target.fPzMax = fPzMax;
350 target.fEtaMin = fEtaMin;
351 target.fEtaMax = fEtaMax;
352 target.fRapMin = fRapMin;
353 target.fRapMax = fRapMax;
354
355 target.fHistogramsOn = fHistogramsOn;
356
357 for (Int_t i=0; i<2; ++i)
358 {
359 if (fhNClustersITS[i]) target.fhNClustersITS[i] = (TH1F*) fhNClustersITS[i]->Clone();
360 if (fhNClustersTPC[i]) target.fhNClustersTPC[i] = (TH1F*) fhNClustersTPC[i]->Clone();
361
362 if (fhChi2PerClusterITS[i]) target.fhChi2PerClusterITS[i] = (TH1F*) fhChi2PerClusterITS[i]->Clone();
363 if (fhChi2PerClusterTPC[i]) target.fhChi2PerClusterTPC[i] = (TH1F*) fhChi2PerClusterTPC[i]->Clone();
364
365 if (fhC11[i]) target.fhC11[i] = (TH1F*) fhC11[i]->Clone();
366 if (fhC22[i]) target.fhC22[i] = (TH1F*) fhC22[i]->Clone();
367 if (fhC33[i]) target.fhC33[i] = (TH1F*) fhC33[i]->Clone();
368 if (fhC44[i]) target.fhC44[i] = (TH1F*) fhC44[i]->Clone();
369 if (fhC55[i]) target.fhC55[i] = (TH1F*) fhC55[i]->Clone();
370
371 if (fhDXY[i]) target.fhDXY[i] = (TH1F*) fhDXY[i]->Clone();
372 if (fhDZ[i]) target.fhDZ[i] = (TH1F*) fhDZ[i]->Clone();
373 if (fhDXYvsDZ[i]) target.fhDXYvsDZ[i] = (TH2F*) fhDXYvsDZ[i]->Clone();
374
375 if (fhDXYNormalized[i]) target.fhDXYNormalized[i] = (TH1F*) fhDXYNormalized[i]->Clone();
376 if (fhDZNormalized[i]) target.fhDZNormalized[i] = (TH1F*) fhDZNormalized[i]->Clone();
377 if (fhDXYvsDZNormalized[i]) target.fhDXYvsDZNormalized[i] = (TH2F*) fhDXYvsDZNormalized[i]->Clone();
378 }
379
380 if (fhCutStatistics) target.fhCutStatistics = (TH1F*) fhCutStatistics->Clone();
381 if (fhCutCorrelation) target.fhCutCorrelation = (TH2F*) fhCutCorrelation->Clone();
0b75bef2 382
0b8a2421 383 TNamed::Copy(c);
0b75bef2 384}
385
927c71e5 386//_____________________________________________________________________________
387Long64_t AliESDtrackCuts::Merge(TCollection* list) {
388 // Merge a list of AliESDtrackCuts objects with this (needed for PROOF)
389 // Returns the number of merged objects (including this)
390
391 if (!list)
392 return 0;
393
394 if (list->IsEmpty())
395 return 1;
396
397 if (!fHistogramsOn)
398 return 0;
399
400 TIterator* iter = list->MakeIterator();
401 TObject* obj;
402
403
404 // collection of measured and generated histograms
405 Int_t count = 0;
406 while ((obj = iter->Next())) {
407
408 AliESDtrackCuts* entry = dynamic_cast<AliESDtrackCuts*>(obj);
409 if (entry == 0)
410 continue;
411
412 if (!entry->fHistogramsOn)
413 continue;
414
415 for (Int_t i=0; i<2; i++) {
416
417 fhNClustersITS[i] ->Add(entry->fhNClustersITS[i] );
418 fhNClustersTPC[i] ->Add(entry->fhNClustersTPC[i] );
419
420 fhChi2PerClusterITS[i] ->Add(entry->fhChi2PerClusterITS[i]);
421 fhChi2PerClusterTPC[i] ->Add(entry->fhChi2PerClusterTPC[i]);
422
423 fhC11[i] ->Add(entry->fhC11[i] );
424 fhC22[i] ->Add(entry->fhC22[i] );
425 fhC33[i] ->Add(entry->fhC33[i] );
426 fhC44[i] ->Add(entry->fhC44[i] );
427 fhC55[i] ->Add(entry->fhC55[i] );
428
429 fhDXY[i] ->Add(entry->fhDXY[i] );
430 fhDZ[i] ->Add(entry->fhDZ[i] );
431 fhDXYvsDZ[i] ->Add(entry->fhDXYvsDZ[i] );
432
433 fhDXYNormalized[i] ->Add(entry->fhDXYNormalized[i] );
434 fhDZNormalized[i] ->Add(entry->fhDZNormalized[i] );
435 fhDXYvsDZNormalized[i] ->Add(entry->fhDXYvsDZNormalized[i]);
436 }
437
438 fhCutStatistics ->Add(entry->fhCutStatistics);
439 fhCutCorrelation ->Add(entry->fhCutCorrelation);
440
441 count++;
442 }
443
444 return count+1;
445}
446
447
0b75bef2 448//____________________________________________________________________
2b9019e4 449Float_t AliESDtrackCuts::GetSigmaToVertex(AliESDtrack* esdTrack)
450{
451 //
452
453 Float_t b[2];
454 Float_t bRes[2];
455 Float_t bCov[3];
456 esdTrack->GetImpactParameters(b,bCov);
457 if (bCov[0]<=0 || bCov[2]<=0) {
458 AliDebug(1, "Estimated b resolution lower or equal zero!");
459 bCov[0]=0; bCov[2]=0;
460 }
461 bRes[0] = TMath::Sqrt(bCov[0]);
462 bRes[1] = TMath::Sqrt(bCov[2]);
463
464 // -----------------------------------
465 // How to get to a n-sigma cut?
466 //
467 // The accumulated statistics from 0 to d is
468 //
469 // -> Erf(d/Sqrt(2)) for a 1-dim gauss (d = n_sigma)
470 // -> 1 - Exp(-d**2) for a 2-dim gauss (d*d = dx*dx + dy*dy != n_sigma)
471 //
472 // It means that for a 2-dim gauss: n_sigma(d) = Sqrt(2)*ErfInv(1 - Exp((-x**2)/2)
473 // Can this be expressed in a different way?
2b9019e4 474
475 if (bRes[0] == 0 || bRes[1] ==0)
476 return -1;
477
478 Float_t d = TMath::Sqrt(TMath::Power(b[0]/bRes[0],2) + TMath::Power(b[1]/bRes[1],2));
479
480 // stupid rounding problem screws up everything:
481 // if d is too big, TMath::Exp(...) gets 0, and TMath::ErfInverse(1) that should be infinite, gets 0 :(
482 if (TMath::Exp(-d * d / 2) < 1e-10)
483 return 1000;
484
485 d = TMath::ErfInverse(1 - TMath::Exp(-d * d / 2)) * TMath::Sqrt(2);
486 return d;
487}
488
5c495d37 489void AliESDtrackCuts::EnableNeededBranches(TTree* tree)
490{
491 // enables the branches needed by AcceptTrack, for a list see comment of AcceptTrack
492
493 tree->SetBranchStatus("fTracks.fFlags", 1);
494 tree->SetBranchStatus("fTracks.fITSncls", 1);
495 tree->SetBranchStatus("fTracks.fTPCncls", 1);
496 tree->SetBranchStatus("fTracks.fITSchi2", 1);
497 tree->SetBranchStatus("fTracks.fTPCchi2", 1);
498 tree->SetBranchStatus("fTracks.fC*", 1);
499 tree->SetBranchStatus("fTracks.fD", 1);
500 tree->SetBranchStatus("fTracks.fZ", 1);
501 tree->SetBranchStatus("fTracks.fCdd", 1);
502 tree->SetBranchStatus("fTracks.fCdz", 1);
503 tree->SetBranchStatus("fTracks.fCzz", 1);
504 tree->SetBranchStatus("fTracks.fP*", 1);
505 tree->SetBranchStatus("fTracks.fR*", 1);
506 tree->SetBranchStatus("fTracks.fKinkIndexes*", 1);
507}
508
2b9019e4 509//____________________________________________________________________
510Bool_t
0b75bef2 511AliESDtrackCuts::AcceptTrack(AliESDtrack* esdTrack) {
512 //
513 // figure out if the tracks survives all the track cuts defined
514 //
ab2c1f0d 515 // the different quality parameter and kinematic values are first
516 // retrieved from the track. then it is found out what cuts the
517 // track did not survive and finally the cuts are imposed.
0b75bef2 518
5c495d37 519 // this function needs the following branches:
520 // fTracks.fFlags
521 // fTracks.fITSncls
522 // fTracks.fTPCncls
523 // fTracks.fITSchi2
524 // fTracks.fTPCchi2
525 // fTracks.fC //GetExternalCovariance
526 // fTracks.fD //GetImpactParameters
527 // fTracks.fZ //GetImpactParameters
528 // fTracks.fCdd //GetImpactParameters
529 // fTracks.fCdz //GetImpactParameters
530 // fTracks.fCzz //GetImpactParameters
531 // fTracks.fP //GetPxPyPz
532 // fTracks.fR //GetMass
533 // fTracks.fP //GetMass
534 // fTracks.fKinkIndexes
1ae6c515 535
0b75bef2 536 UInt_t status = esdTrack->GetStatus();
ebd60d8e 537
538 // dummy array
539 Int_t fIdxInt[200];
540
0b75bef2 541 // getting quality parameters from the ESD track
542 Int_t nClustersITS = esdTrack->GetITSclusters(fIdxInt);
543 Int_t nClustersTPC = esdTrack->GetTPCclusters(fIdxInt);
544
1ae6c515 545
546
0b75bef2 547 Float_t chi2PerClusterITS = -1;
548 Float_t chi2PerClusterTPC = -1;
549 if (nClustersITS!=0)
550 chi2PerClusterITS = esdTrack->GetITSchi2()/Float_t(nClustersITS);
551 if (nClustersTPC!=0)
5c495d37 552 chi2PerClusterTPC = esdTrack->GetTPCchi2()/Float_t(nClustersTPC);
0b75bef2 553
554 Double_t extCov[15];
ebd60d8e 555 esdTrack->GetExternalCovariance(extCov);
0b75bef2 556
557 // getting the track to vertex parameters
2b9019e4 558 Float_t nSigmaToVertex = GetSigmaToVertex(esdTrack);
0b75bef2 559
2b9019e4 560 // getting the kinematic variables of the track
0b75bef2 561 // (assuming the mass is known)
562 Double_t p[3];
563 esdTrack->GetPxPyPz(p);
564 Float_t momentum = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2) + TMath::Power(p[2],2));
565 Float_t pt = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2));
566 Float_t energy = TMath::Sqrt(TMath::Power(esdTrack->GetMass(),2) + TMath::Power(momentum,2));
567
1ae6c515 568
0b75bef2 569 //y-eta related calculations
570 Float_t eta = -100.;
571 Float_t y = -100.;
572 if((momentum != TMath::Abs(p[2]))&&(momentum != 0))
573 eta = 0.5*TMath::Log((momentum + p[2])/(momentum - p[2]));
574 if((energy != TMath::Abs(p[2]))&&(momentum != 0))
575 y = 0.5*TMath::Log((energy + p[2])/(energy - p[2]));
576
577
578 //########################################################################
579 // cut the track?
580
ebd60d8e 581 Bool_t cuts[kNCuts];
582 for (Int_t i=0; i<kNCuts; i++) cuts[i]=kFALSE;
0b75bef2 583
584 // track quality cuts
ab2c1f0d 585 if (fCutRequireTPCRefit && (status&AliESDtrack::kTPCrefit)==0)
0b75bef2 586 cuts[0]=kTRUE;
ab2c1f0d 587 if (fCutRequireITSRefit && (status&AliESDtrack::kITSrefit)==0)
0b75bef2 588 cuts[1]=kTRUE;
2b9019e4 589 if (nClustersTPC<fCutMinNClusterTPC)
0b75bef2 590 cuts[2]=kTRUE;
ab2c1f0d 591 if (nClustersITS<fCutMinNClusterITS)
0b75bef2 592 cuts[3]=kTRUE;
ab2c1f0d 593 if (chi2PerClusterTPC>fCutMaxChi2PerClusterTPC)
0b75bef2 594 cuts[4]=kTRUE;
ab2c1f0d 595 if (chi2PerClusterITS>fCutMaxChi2PerClusterITS)
0b75bef2 596 cuts[5]=kTRUE;
ab2c1f0d 597 if (extCov[0] > fCutMaxC11)
0b75bef2 598 cuts[6]=kTRUE;
ab2c1f0d 599 if (extCov[2] > fCutMaxC22)
0b75bef2 600 cuts[7]=kTRUE;
ab2c1f0d 601 if (extCov[5] > fCutMaxC33)
0b75bef2 602 cuts[8]=kTRUE;
ab2c1f0d 603 if (extCov[9] > fCutMaxC44)
0b75bef2 604 cuts[9]=kTRUE;
ab2c1f0d 605 if (extCov[14] > fCutMaxC55)
0b75bef2 606 cuts[10]=kTRUE;
2b9019e4 607 if (nSigmaToVertex > fCutNsigmaToVertex)
0b75bef2 608 cuts[11] = kTRUE;
609 // if n sigma could not be calculated
10ebe68d 610 if (nSigmaToVertex<0 && fCutSigmaToVertexRequired)
0b75bef2 611 cuts[12]=kTRUE;
5c495d37 612 if (!fCutAcceptKinkDaughters && esdTrack->GetKinkIndex(0)>0)
0b75bef2 613 cuts[13]=kTRUE;
614 // track kinematics cut
615 if((momentum < fPMin) || (momentum > fPMax))
616 cuts[14]=kTRUE;
617 if((pt < fPtMin) || (pt > fPtMax))
618 cuts[15] = kTRUE;
619 if((p[0] < fPxMin) || (p[0] > fPxMax))
620 cuts[16] = kTRUE;
621 if((p[1] < fPyMin) || (p[1] > fPyMax))
622 cuts[17] = kTRUE;
2b9019e4 623 if((p[2] < fPzMin) || (p[2] > fPzMax))
0b75bef2 624 cuts[18] = kTRUE;
625 if((eta < fEtaMin) || (eta > fEtaMax))
626 cuts[19] = kTRUE;
627 if((y < fRapMin) || (y > fRapMax))
628 cuts[20] = kTRUE;
629
630 Bool_t cut=kFALSE;
ebd60d8e 631 for (Int_t i=0; i<kNCuts; i++)
0b75bef2 632 if (cuts[i]) cut = kTRUE;
633
634 //########################################################################
635 // filling histograms
636 if (fHistogramsOn) {
f58f1a93 637 fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n tracks")));
0b75bef2 638
639 if (cut)
f58f1a93 640 fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n cut tracks")));
0b75bef2 641
ebd60d8e 642 for (Int_t i=0; i<kNCuts; i++) {
0b75bef2 643 if (cuts[i])
ab2c1f0d 644 fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin(fgkCutNames[i])));
0b75bef2 645
ebd60d8e 646 for (Int_t j=i; j<kNCuts; j++) {
0b75bef2 647 if (cuts[i] && cuts[j]) {
ab2c1f0d 648 Float_t x = fhCutCorrelation->GetXaxis()->GetBinCenter(fhCutCorrelation->GetXaxis()->FindBin(fgkCutNames[i]));
649 Float_t y = fhCutCorrelation->GetYaxis()->GetBinCenter(fhCutCorrelation->GetYaxis()->FindBin(fgkCutNames[j]));
f58f1a93 650 fhCutCorrelation->Fill(x,y);
0b75bef2 651 }
652 }
653 }
654
655
2b9019e4 656 fhNClustersITS[0]->Fill(nClustersITS);
657 fhNClustersTPC[0]->Fill(nClustersTPC);
f58f1a93 658 fhChi2PerClusterITS[0]->Fill(chi2PerClusterITS);
2b9019e4 659 fhChi2PerClusterTPC[0]->Fill(chi2PerClusterTPC);
660
661 fhC11[0]->Fill(extCov[0]);
662 fhC22[0]->Fill(extCov[2]);
663 fhC33[0]->Fill(extCov[5]);
664 fhC44[0]->Fill(extCov[9]);
665 fhC55[0]->Fill(extCov[14]);
666
667 Float_t b[2];
668 Float_t bRes[2];
669 Float_t bCov[3];
670 esdTrack->GetImpactParameters(b,bCov);
671 if (bCov[0]<=0 || bCov[2]<=0) {
672 AliDebug(1, "Estimated b resolution lower or equal zero!");
673 bCov[0]=0; bCov[2]=0;
674 }
675 bRes[0] = TMath::Sqrt(bCov[0]);
676 bRes[1] = TMath::Sqrt(bCov[2]);
677
678 fhDZ[0]->Fill(b[1]);
679 fhDXY[0]->Fill(b[0]);
f58f1a93 680 fhDXYvsDZ[0]->Fill(b[1],b[0]);
0b75bef2 681
682 if (bRes[0]!=0 && bRes[1]!=0) {
10ebe68d 683 fhDZNormalized[0]->Fill(b[1]/bRes[1]);
2b9019e4 684 fhDXYNormalized[0]->Fill(b[0]/bRes[0]);
f58f1a93 685 fhDXYvsDZNormalized[0]->Fill(b[1]/bRes[1], b[0]/bRes[0]);
0b75bef2 686 }
687 }
688
2b9019e4 689 //########################################################################
0b75bef2 690 // cut the track!
691 if (cut) return kFALSE;
692
2b9019e4 693 //########################################################################
0b75bef2 694 // filling histograms after cut
695 if (fHistogramsOn) {
2b9019e4 696 fhNClustersITS[1]->Fill(nClustersITS);
697 fhNClustersTPC[1]->Fill(nClustersTPC);
f58f1a93 698 fhChi2PerClusterITS[1]->Fill(chi2PerClusterITS);
2b9019e4 699 fhChi2PerClusterTPC[1]->Fill(chi2PerClusterTPC);
700
701 fhC11[1]->Fill(extCov[0]);
702 fhC22[1]->Fill(extCov[2]);
10ebe68d 703 fhC33[1]->Fill(extCov[5]);
2b9019e4 704 fhC44[1]->Fill(extCov[9]);
705 fhC55[1]->Fill(extCov[14]);
706
707 Float_t b[2];
708 Float_t bRes[2];
709 Float_t bCov[3];
710 esdTrack->GetImpactParameters(b,bCov);
711 if (bCov[0]<=0 || bCov[2]<=0) {
712 AliDebug(1, "Estimated b resolution lower or equal zero!");
713 bCov[0]=0; bCov[2]=0;
714 }
715 bRes[0] = TMath::Sqrt(bCov[0]);
716 bRes[1] = TMath::Sqrt(bCov[2]);
717
718 fhDZ[1]->Fill(b[1]);
719 fhDXY[1]->Fill(b[0]);
f58f1a93 720 fhDXYvsDZ[1]->Fill(b[1],b[0]);
0b75bef2 721
10ebe68d 722 if (bRes[0]!=0 && bRes[1]!=0)
723 {
724 fhDZNormalized[1]->Fill(b[1]/bRes[1]);
725 fhDXYNormalized[1]->Fill(b[0]/bRes[0]);
726 fhDXYvsDZNormalized[1]->Fill(b[1]/bRes[1], b[0]/bRes[0]);
727 }
0b75bef2 728 }
2b9019e4 729
0b75bef2 730 return kTRUE;
731}
732
733//____________________________________________________________________
734TObjArray*
ebd60d8e 735AliESDtrackCuts::GetAcceptedTracks(AliESD* esd)
736{
737 //
0b75bef2 738 // returns an array of all tracks that pass the cuts
ebd60d8e 739 //
740
741 TObjArray* acceptedTracks = new TObjArray();
c5858666 742
0b75bef2 743 // loop over esd tracks
744 for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
745 AliESDtrack* track = esd->GetTrack(iTrack);
c5858666 746
ebd60d8e 747 if (AcceptTrack(track))
748 acceptedTracks->Add(track);
0b75bef2 749 }
750
ebd60d8e 751 return acceptedTracks;
0b75bef2 752}
753
c5858666 754//____________________________________________________________________
755Int_t
756AliESDtrackCuts::CountAcceptedTracks(AliESD* esd)
757{
758 //
759 // returns an the number of tracks that pass the cuts
760 //
761
762 Int_t count = 0;
763
764 // loop over esd tracks
765 for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
766 AliESDtrack* track = esd->GetTrack(iTrack);
767
768 if (AcceptTrack(track))
769 count++;
770 }
771
772 return count;
773}
774
0b75bef2 775//____________________________________________________________________
f58f1a93 776 void AliESDtrackCuts::DefineHistograms(Int_t color) {
ab2c1f0d 777 //
778 // diagnostics histograms are defined
779 //
0b75bef2 780
f58f1a93 781 fHistogramsOn=kTRUE;
0b75bef2 782
ab2c1f0d 783 //###################################################################################
f58f1a93 784 // defining histograms
0b75bef2 785
ebd60d8e 786 fhCutStatistics = new TH1F("cut_statistics","cut statistics",kNCuts+4,-0.5,kNCuts+3.5);
0b75bef2 787
f58f1a93 788 fhCutStatistics->GetXaxis()->SetBinLabel(1,"n tracks");
789 fhCutStatistics->GetXaxis()->SetBinLabel(2,"n cut tracks");
0b75bef2 790
ebd60d8e 791 fhCutCorrelation = new TH2F("cut_correlation","cut correlation",kNCuts,-0.5,kNCuts-0.5,kNCuts,-0.5,kNCuts-0.5);;
0b75bef2 792
ebd60d8e 793 for (Int_t i=0; i<kNCuts; i++) {
ab2c1f0d 794 fhCutStatistics->GetXaxis()->SetBinLabel(i+4,fgkCutNames[i]);
795 fhCutCorrelation->GetXaxis()->SetBinLabel(i+1,fgkCutNames[i]);
796 fhCutCorrelation->GetYaxis()->SetBinLabel(i+1,fgkCutNames[i]);
f58f1a93 797 }
798
799 fhCutStatistics ->SetLineColor(color);
800 fhCutCorrelation ->SetLineColor(color);
801 fhCutStatistics ->SetLineWidth(2);
802 fhCutCorrelation ->SetLineWidth(2);
803
0b75bef2 804 Char_t str[256];
805 for (Int_t i=0; i<2; i++) {
806 if (i==0) sprintf(str," ");
807 else sprintf(str,"_cut");
808
ebd60d8e 809 fhNClustersITS[i] = new TH1F(Form("nClustersITS%s",str),"",8,-0.5,7.5);
810 fhNClustersTPC[i] = new TH1F(Form("nClustersTPC%s",str),"",165,-0.5,164.5);
811 fhChi2PerClusterITS[i] = new TH1F(Form("chi2PerClusterITS%s",str),"",500,0,10);
812 fhChi2PerClusterTPC[i] = new TH1F(Form("chi2PerClusterTPC%s",str),"",500,0,10);
813
d8594d03 814 fhC11[i] = new TH1F(Form("covMatrixDiagonal11%s",str),"",2000,0,20);
815 fhC22[i] = new TH1F(Form("covMatrixDiagonal22%s",str),"",2000,0,20);
816 fhC33[i] = new TH1F(Form("covMatrixDiagonal33%s",str),"",1000,0,1);
ebd60d8e 817 fhC44[i] = new TH1F(Form("covMatrixDiagonal44%s",str),"",1000,0,5);
818 fhC55[i] = new TH1F(Form("covMatrixDiagonal55%s",str),"",1000,0,5);
819
820 fhDXY[i] = new TH1F(Form("dXY%s",str),"",500,-10,10);
821 fhDZ[i] = new TH1F(Form("dZ%s",str),"",500,-10,10);
822 fhDXYvsDZ[i] = new TH2F(Form("dXYvsDZ%s",str),"",200,-10,10,200,-10,10);
823
824 fhDXYNormalized[i] = new TH1F(Form("dXYNormalized%s",str),"",500,-10,10);
825 fhDZNormalized[i] = new TH1F(Form("dZNormalized%s",str),"",500,-10,10);
826 fhDXYvsDZNormalized[i] = new TH2F(Form("dXYvsDZNormalized%s",str),"",200,-10,10,200,-10,10);
827
828
829 fhNClustersITS[i]->SetXTitle("n ITS clusters");
830 fhNClustersTPC[i]->SetXTitle("n TPC clusters");
831 fhChi2PerClusterITS[i]->SetXTitle("#Chi^{2} per ITS cluster");
832 fhChi2PerClusterTPC[i]->SetXTitle("#Chi^{2} per TPC cluster");
833
834 fhC11[i]->SetXTitle("cov 11 : #sigma_{y}^{2} [cm^{2}]");
835 fhC22[i]->SetXTitle("cov 22 : #sigma_{z}^{2} [cm^{2}]");
836 fhC33[i]->SetXTitle("cov 33 : #sigma_{sin(#phi)}^{2}");
837 fhC44[i]->SetXTitle("cov 44 : #sigma_{tan(#theta_{dip})}^{2}");
838 fhC55[i]->SetXTitle("cov 55 : #sigma_{1/p_{T}}^{2} [(c/GeV)^2]");
839
840 fhDXY[i]->SetXTitle("transverse impact parameter");
841 fhDZ[i]->SetXTitle("longitudinal impact parameter");
842 fhDXYvsDZ[i]->SetXTitle("longitudinal impact parameter");
843 fhDXYvsDZ[i]->SetYTitle("transverse impact parameter");
844
845 fhDXYNormalized[i]->SetXTitle("normalized trans impact par");
846 fhDZNormalized[i]->SetXTitle("normalized long impact par");
847 fhDXYvsDZNormalized[i]->SetXTitle("normalized long impact par");
848 fhDXYvsDZNormalized[i]->SetYTitle("normalized trans impact par");
849
850 fhNClustersITS[i]->SetLineColor(color); fhNClustersITS[i]->SetLineWidth(2);
851 fhNClustersTPC[i]->SetLineColor(color); fhNClustersTPC[i]->SetLineWidth(2);
852 fhChi2PerClusterITS[i]->SetLineColor(color); fhChi2PerClusterITS[i]->SetLineWidth(2);
853 fhChi2PerClusterTPC[i]->SetLineColor(color); fhChi2PerClusterTPC[i]->SetLineWidth(2);
854
855 fhC11[i]->SetLineColor(color); fhC11[i]->SetLineWidth(2);
856 fhC22[i]->SetLineColor(color); fhC22[i]->SetLineWidth(2);
857 fhC33[i]->SetLineColor(color); fhC33[i]->SetLineWidth(2);
858 fhC44[i]->SetLineColor(color); fhC44[i]->SetLineWidth(2);
859 fhC55[i]->SetLineColor(color); fhC55[i]->SetLineWidth(2);
860
861 fhDXY[i]->SetLineColor(color); fhDXY[i]->SetLineWidth(2);
862 fhDZ[i]->SetLineColor(color); fhDZ[i]->SetLineWidth(2);
863
864 fhDXYNormalized[i]->SetLineColor(color); fhDXYNormalized[i]->SetLineWidth(2);
865 fhDZNormalized[i]->SetLineColor(color); fhDZNormalized[i]->SetLineWidth(2);
0b75bef2 866 }
867}
868
d8594d03 869
870
0b75bef2 871//____________________________________________________________________
872void
49dc84d9 873AliESDtrackCuts::Print(const Option_t*) const {
ab2c1f0d 874 //
875 // print method - still to be implemented
876 //
0b75bef2 877
878 AliInfo("AliESDtrackCuts...");
879}
880
881
882//____________________________________________________________________
f58f1a93 883void AliESDtrackCuts::SaveHistograms(Char_t* dir) {
ab2c1f0d 884 //
885 // saves the histograms in a directory (dir)
886 //
887
0b75bef2 888
889 if (!fHistogramsOn) {
890 AliDebug(0, "Histograms not on - cannot save histograms!!!");
891 return;
892 }
893
894 gDirectory->mkdir(dir);
895 gDirectory->cd(dir);
896
897 gDirectory->mkdir("before_cuts");
898 gDirectory->mkdir("after_cuts");
899
f58f1a93 900 fhCutStatistics->Write();
901 fhCutCorrelation->Write();
0b75bef2 902
903 for (Int_t i=0; i<2; i++) {
904 if (i==0)
905 gDirectory->cd("before_cuts");
906 else
907 gDirectory->cd("after_cuts");
908
f58f1a93 909 fhNClustersITS[i] ->Write();
910 fhNClustersTPC[i] ->Write();
911 fhChi2PerClusterITS[i] ->Write();
912 fhChi2PerClusterTPC[i] ->Write();
0b75bef2 913
f58f1a93 914 fhC11[i] ->Write();
915 fhC22[i] ->Write();
916 fhC33[i] ->Write();
917 fhC44[i] ->Write();
918 fhC55[i] ->Write();
919
920 fhDXY[i] ->Write();
921 fhDZ[i] ->Write();
922 fhDXYvsDZ[i] ->Write();
0b75bef2 923
f58f1a93 924 fhDXYNormalized[i] ->Write();
925 fhDZNormalized[i] ->Write();
926 fhDXYvsDZNormalized[i] ->Write();
0b75bef2 927
928 gDirectory->cd("../");
929 }
930
931 gDirectory->cd("../");
932}
933
934
935