]> git.uio.no Git - u/mrichter/AliRoot.git/blame_incremental - ANALYSIS/AliESDtrackCuts.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame_incremental
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (non-incremental) | history | HEAD
Fixed bin size bug
[u/mrichter/AliRoot.git] / ANALYSIS / AliESDtrackCuts.cxx
... / ...
CommitLineData
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/* $Id: AliESDtrackCuts.cxx 24534 2008-03-16 22:22:11Z fca $ */
17
18#include "AliESDtrackCuts.h"
19
20#include <AliESDtrack.h>
21#include <AliESDVertex.h>
22#include <AliESDEvent.h>
23#include <AliMultiplicity.h>
24#include <AliLog.h>
25
26#include <TTree.h>
27#include <TCanvas.h>
28#include <TDirectory.h>
29#include <TH2F.h>
30#include <TF1.h>
31#include <TBits.h>
32
33//____________________________________________________________________
34ClassImp(AliESDtrackCuts)
35
36// Cut names
37const Char_t* AliESDtrackCuts::fgkCutNames[kNCuts] = {
38 "require TPC refit",
39 "require TPC standalone",
40 "require ITS refit",
41 "n clusters TPC",
42 "n clusters ITS",
43 "#Chi^{2}/cluster TPC",
44 "#Chi^{2}/cluster ITS",
45 "cov 11",
46 "cov 22",
47 "cov 33",
48 "cov 44",
49 "cov 55",
50 "trk-to-vtx",
51 "trk-to-vtx failed",
52 "kink daughters",
53 "p",
54 "p_{T}",
55 "p_{x}",
56 "p_{y}",
57 "p_{z}",
58 "eta",
59 "y",
60 "trk-to-vtx max dca 2D absolute",
61 "trk-to-vtx max dca xy absolute",
62 "trk-to-vtx max dca z absolute",
63 "trk-to-vtx min dca 2D absolute",
64 "trk-to-vtx min dca xy absolute",
65 "trk-to-vtx min dca z absolute",
66 "SPD cluster requirement",
67 "SDD cluster requirement",
68 "SSD cluster requirement",
69 "require ITS stand-alone",
70 "rel 1/pt uncertainty",
71 "TPC n shared clusters",
72 "TPC rel shared clusters",
73 "require ITS Pid",
74 "n crossed rows TPC",
75 "n crossed rows / n findable clusters",
76 "missing ITS points",
77 "#Chi^{2} TPC constrained vs. global",
78 "require TOF out",
79 "TOF Distance cut",
80 "min length in active volume TPC"
81};
82
83AliESDtrackCuts* AliESDtrackCuts::fgMultEstTrackCuts[AliESDtrackCuts::kNMultEstTrackCuts] = { 0, 0, 0, 0 };
84Char_t AliESDtrackCuts::fgBeamTypeFlag = -1;
85
86//____________________________________________________________________
87AliESDtrackCuts::AliESDtrackCuts(const Char_t* name, const Char_t* title) : AliAnalysisCuts(name,title),
88 fCutMinNClusterTPC(0),
89 fCutMinNClusterITS(0),
90 fCutMinNCrossedRowsTPC(0),
91 fCutMinRatioCrossedRowsOverFindableClustersTPC(0),
92 f1CutMinNClustersTPCPtDep(0x0),
93 fCutMaxPtDepNClustersTPC(0),
94 fCutMinLengthActiveVolumeTPC(0),
95 fCutMaxChi2PerClusterTPC(0),
96 fCutMaxChi2PerClusterITS(0),
97 fCutMaxChi2TPCConstrainedVsGlobal(0),
98 fCutMaxChi2TPCConstrainedVsGlobalVertexType(kVertexTracks | kVertexSPD),
99 fCutMaxMissingITSPoints(0),
100 fCutMaxC11(0),
101 fCutMaxC22(0),
102 fCutMaxC33(0),
103 fCutMaxC44(0),
104 fCutMaxC55(0),
105 fCutMaxRel1PtUncertainty(0),
106 fCutAcceptKinkDaughters(0),
107 fCutAcceptSharedTPCClusters(0),
108 fCutMaxFractionSharedTPCClusters(0),
109 fCutRequireTPCRefit(0),
110 fCutRequireTPCStandAlone(0),
111 fCutRequireITSRefit(0),
112 fCutRequireITSPid(0),
113 fCutRequireITSStandAlone(0),
114 fCutRequireITSpureSA(0),
115 fCutNsigmaToVertex(0),
116 fCutSigmaToVertexRequired(0),
117 fCutMaxDCAToVertexXY(0),
118 fCutMaxDCAToVertexZ(0),
119 fCutMinDCAToVertexXY(0),
120 fCutMinDCAToVertexZ(0),
121 fCutMaxDCAToVertexXYPtDep(""),
122 fCutMaxDCAToVertexZPtDep(""),
123 fCutMinDCAToVertexXYPtDep(""),
124 fCutMinDCAToVertexZPtDep(""),
125 f1CutMaxDCAToVertexXYPtDep(0x0),
126 f1CutMaxDCAToVertexZPtDep(0x0),
127 f1CutMinDCAToVertexXYPtDep(0x0),
128 f1CutMinDCAToVertexZPtDep(0x0),
129 fCutDCAToVertex2D(0),
130 fPMin(0),
131 fPMax(0),
132 fPtMin(0),
133 fPtMax(0),
134 fPxMin(0),
135 fPxMax(0),
136 fPyMin(0),
137 fPyMax(0),
138 fPzMin(0),
139 fPzMax(0),
140 fEtaMin(0),
141 fEtaMax(0),
142 fRapMin(0),
143 fRapMax(0),
144 fCutRequireTOFout(kFALSE),
145 fFlagCutTOFdistance(kFALSE),
146 fCutTOFdistance(3.),
147 fHistogramsOn(0),
148 ffDTheoretical(0),
149 fhCutStatistics(0),
150 fhCutCorrelation(0)
151{
152 //
153 // constructor
154 //
155
156 Init();
157
158 //##############################################################################
159 // setting default cuts
160 SetMinNClustersTPC();
161 SetMinNClustersITS();
162 SetMinNCrossedRowsTPC();
163 SetMinRatioCrossedRowsOverFindableClustersTPC();
164 SetMaxChi2PerClusterTPC();
165 SetMaxChi2PerClusterITS();
166 SetMaxChi2TPCConstrainedGlobal();
167 SetMaxChi2TPCConstrainedGlobalVertexType();
168 SetMaxNOfMissingITSPoints();
169 SetMaxCovDiagonalElements();
170 SetMaxRel1PtUncertainty();
171 SetRequireTPCRefit();
172 SetRequireTPCStandAlone();
173 SetRequireITSRefit();
174 SetRequireITSPid(kFALSE);
175 SetRequireITSStandAlone(kFALSE);
176 SetRequireITSPureStandAlone(kFALSE);
177 SetAcceptKinkDaughters();
178 SetAcceptSharedTPCClusters();
179 SetMaxFractionSharedTPCClusters();
180 SetMaxNsigmaToVertex();
181 SetMaxDCAToVertexXY();
182 SetMaxDCAToVertexZ();
183 SetDCAToVertex2D();
184 SetMinDCAToVertexXY();
185 SetMinDCAToVertexZ();
186 SetPRange();
187 SetPtRange();
188 SetPxRange();
189 SetPyRange();
190 SetPzRange();
191 SetEtaRange();
192 SetRapRange();
193 SetClusterRequirementITS(kSPD);
194 SetClusterRequirementITS(kSDD);
195 SetClusterRequirementITS(kSSD);
196
197 SetHistogramsOn();
198}
199
200//_____________________________________________________________________________
201AliESDtrackCuts::AliESDtrackCuts(const AliESDtrackCuts &c) : AliAnalysisCuts(c),
202 fCutMinNClusterTPC(0),
203 fCutMinNClusterITS(0),
204 fCutMinNCrossedRowsTPC(0),
205 fCutMinRatioCrossedRowsOverFindableClustersTPC(0),
206 f1CutMinNClustersTPCPtDep(0x0),
207 fCutMaxPtDepNClustersTPC(0),
208 fCutMinLengthActiveVolumeTPC(0),
209 fCutMaxChi2PerClusterTPC(0),
210 fCutMaxChi2PerClusterITS(0),
211 fCutMaxChi2TPCConstrainedVsGlobal(0),
212 fCutMaxChi2TPCConstrainedVsGlobalVertexType(kVertexTracks | kVertexSPD),
213 fCutMaxMissingITSPoints(0),
214 fCutMaxC11(0),
215 fCutMaxC22(0),
216 fCutMaxC33(0),
217 fCutMaxC44(0),
218 fCutMaxC55(0),
219 fCutMaxRel1PtUncertainty(0),
220 fCutAcceptKinkDaughters(0),
221 fCutAcceptSharedTPCClusters(0),
222 fCutMaxFractionSharedTPCClusters(0),
223 fCutRequireTPCRefit(0),
224 fCutRequireTPCStandAlone(0),
225 fCutRequireITSRefit(0),
226 fCutRequireITSPid(0),
227 fCutRequireITSStandAlone(0),
228 fCutRequireITSpureSA(0),
229 fCutNsigmaToVertex(0),
230 fCutSigmaToVertexRequired(0),
231 fCutMaxDCAToVertexXY(0),
232 fCutMaxDCAToVertexZ(0),
233 fCutMinDCAToVertexXY(0),
234 fCutMinDCAToVertexZ(0),
235 fCutMaxDCAToVertexXYPtDep(""),
236 fCutMaxDCAToVertexZPtDep(""),
237 fCutMinDCAToVertexXYPtDep(""),
238 fCutMinDCAToVertexZPtDep(""),
239 f1CutMaxDCAToVertexXYPtDep(0x0),
240 f1CutMaxDCAToVertexZPtDep(0x0),
241 f1CutMinDCAToVertexXYPtDep(0x0),
242 f1CutMinDCAToVertexZPtDep(0x0),
243 fCutDCAToVertex2D(0),
244 fPMin(0),
245 fPMax(0),
246 fPtMin(0),
247 fPtMax(0),
248 fPxMin(0),
249 fPxMax(0),
250 fPyMin(0),
251 fPyMax(0),
252 fPzMin(0),
253 fPzMax(0),
254 fEtaMin(0),
255 fEtaMax(0),
256 fRapMin(0),
257 fRapMax(0),
258 fCutRequireTOFout(kFALSE),
259 fFlagCutTOFdistance(kFALSE),
260 fCutTOFdistance(3.),
261 fHistogramsOn(0),
262 ffDTheoretical(0),
263 fhCutStatistics(0),
264 fhCutCorrelation(0)
265{
266 //
267 // copy constructor
268 //
269
270 ((AliESDtrackCuts &) c).Copy(*this);
271}
272
273AliESDtrackCuts::~AliESDtrackCuts()
274{
275 //
276 // destructor
277 //
278
279 for (Int_t i=0; i<2; i++) {
280
281 if (fhNClustersITS[i])
282 delete fhNClustersITS[i];
283 if (fhNClustersTPC[i])
284 delete fhNClustersTPC[i];
285 if (fhNSharedClustersTPC[i])
286 delete fhNSharedClustersTPC[i];
287 if (fhNCrossedRowsTPC[i])
288 delete fhNCrossedRowsTPC[i];
289 if (fhRatioCrossedRowsOverFindableClustersTPC[i])
290 delete fhRatioCrossedRowsOverFindableClustersTPC[i];
291 if (fhChi2PerClusterITS[i])
292 delete fhChi2PerClusterITS[i];
293 if (fhChi2PerClusterTPC[i])
294 delete fhChi2PerClusterTPC[i];
295 if (fhChi2TPCConstrainedVsGlobal[i])
296 delete fhChi2TPCConstrainedVsGlobal[i];
297 if(fhNClustersForITSPID[i])
298 delete fhNClustersForITSPID[i];
299 if(fhNMissingITSPoints[i])
300 delete fhNMissingITSPoints[i];
301 if (fhC11[i])
302 delete fhC11[i];
303 if (fhC22[i])
304 delete fhC22[i];
305 if (fhC33[i])
306 delete fhC33[i];
307 if (fhC44[i])
308 delete fhC44[i];
309 if (fhC55[i])
310 delete fhC55[i];
311
312 if (fhRel1PtUncertainty[i])
313 delete fhRel1PtUncertainty[i];
314
315 if (fhDXY[i])
316 delete fhDXY[i];
317 if (fhDZ[i])
318 delete fhDZ[i];
319 if (fhDXYDZ[i])
320 delete fhDXYDZ[i];
321 if (fhDXYvsDZ[i])
322 delete fhDXYvsDZ[i];
323
324 if (fhDXYNormalized[i])
325 delete fhDXYNormalized[i];
326 if (fhDZNormalized[i])
327 delete fhDZNormalized[i];
328 if (fhDXYvsDZNormalized[i])
329 delete fhDXYvsDZNormalized[i];
330 if (fhNSigmaToVertex[i])
331 delete fhNSigmaToVertex[i];
332 if (fhPt[i])
333 delete fhPt[i];
334 if (fhEta[i])
335 delete fhEta[i];
336 if (fhTOFdistance[i])
337 delete fhTOFdistance[i];
338 }
339
340 if(f1CutMaxDCAToVertexXYPtDep)delete f1CutMaxDCAToVertexXYPtDep;
341 f1CutMaxDCAToVertexXYPtDep = 0;
342 if( f1CutMaxDCAToVertexZPtDep) delete f1CutMaxDCAToVertexZPtDep;
343 f1CutMaxDCAToVertexZPtDep = 0;
344 if( f1CutMinDCAToVertexXYPtDep)delete f1CutMinDCAToVertexXYPtDep;
345 f1CutMinDCAToVertexXYPtDep = 0;
346 if(f1CutMinDCAToVertexZPtDep)delete f1CutMinDCAToVertexZPtDep;
347 f1CutMinDCAToVertexZPtDep = 0;
348
349
350 if (ffDTheoretical)
351 delete ffDTheoretical;
352
353 if (fhCutStatistics)
354 delete fhCutStatistics;
355 if (fhCutCorrelation)
356 delete fhCutCorrelation;
357
358 if(f1CutMinNClustersTPCPtDep)
359 delete f1CutMinNClustersTPCPtDep;
360
361}
362
363void AliESDtrackCuts::Init()
364{
365 //
366 // sets everything to zero
367 //
368
369 fCutMinNClusterTPC = 0;
370 fCutMinNClusterITS = 0;
371
372 fCutMaxChi2PerClusterTPC = 0;
373 fCutMaxChi2PerClusterITS = 0;
374 fCutMaxChi2TPCConstrainedVsGlobal = 0;
375 fCutMaxChi2TPCConstrainedVsGlobalVertexType = kVertexTracks | kVertexSPD;
376 fCutMaxMissingITSPoints = 0;
377
378 for (Int_t i = 0; i < 3; i++)
379 fCutClusterRequirementITS[i] = kOff;
380
381 fCutMaxC11 = 0;
382 fCutMaxC22 = 0;
383 fCutMaxC33 = 0;
384 fCutMaxC44 = 0;
385 fCutMaxC55 = 0;
386
387 fCutMaxRel1PtUncertainty = 0;
388
389 fCutAcceptKinkDaughters = 0;
390 fCutAcceptSharedTPCClusters = 0;
391 fCutMaxFractionSharedTPCClusters = 0;
392 fCutRequireTPCRefit = 0;
393 fCutRequireTPCStandAlone = 0;
394 fCutRequireITSRefit = 0;
395 fCutRequireITSPid = 0;
396 fCutRequireITSStandAlone = 0;
397 fCutRequireITSpureSA = 0;
398
399 fCutNsigmaToVertex = 0;
400 fCutSigmaToVertexRequired = 0;
401 fCutMaxDCAToVertexXY = 0;
402 fCutMaxDCAToVertexZ = 0;
403 fCutDCAToVertex2D = 0;
404 fCutMinDCAToVertexXY = 0;
405 fCutMinDCAToVertexZ = 0;
406 fCutMaxDCAToVertexXYPtDep = "";
407 fCutMaxDCAToVertexZPtDep = "";
408 fCutMinDCAToVertexXYPtDep = "";
409 fCutMinDCAToVertexZPtDep = "";
410
411 if(f1CutMaxDCAToVertexXYPtDep)delete f1CutMaxDCAToVertexXYPtDep;
412 f1CutMaxDCAToVertexXYPtDep = 0;
413 if( f1CutMaxDCAToVertexXYPtDep) delete f1CutMaxDCAToVertexXYPtDep;
414 f1CutMaxDCAToVertexXYPtDep = 0;
415 if( f1CutMaxDCAToVertexZPtDep) delete f1CutMaxDCAToVertexZPtDep;
416 f1CutMaxDCAToVertexZPtDep = 0;
417 if( f1CutMinDCAToVertexXYPtDep)delete f1CutMinDCAToVertexXYPtDep;
418 f1CutMinDCAToVertexXYPtDep = 0;
419 if(f1CutMinDCAToVertexZPtDep)delete f1CutMinDCAToVertexZPtDep;
420 f1CutMinDCAToVertexZPtDep = 0;
421
422
423 fPMin = 0;
424 fPMax = 0;
425 fPtMin = 0;
426 fPtMax = 0;
427 fPxMin = 0;
428 fPxMax = 0;
429 fPyMin = 0;
430 fPyMax = 0;
431 fPzMin = 0;
432 fPzMax = 0;
433 fEtaMin = 0;
434 fEtaMax = 0;
435 fRapMin = 0;
436 fRapMax = 0;
437
438 fHistogramsOn = kFALSE;
439
440 for (Int_t i=0; i<2; ++i)
441 {
442 fhNClustersITS[i] = 0;
443 fhNClustersTPC[i] = 0;
444 fhNSharedClustersTPC[i] = 0;
445 fhNCrossedRowsTPC[i] = 0;
446 fhRatioCrossedRowsOverFindableClustersTPC[i] = 0;
447
448 fhChi2PerClusterITS[i] = 0;
449 fhChi2PerClusterTPC[i] = 0;
450 fhChi2TPCConstrainedVsGlobal[i] = 0;
451 fhNClustersForITSPID[i] = 0;
452 fhNMissingITSPoints[i] = 0;
453
454 fhC11[i] = 0;
455 fhC22[i] = 0;
456 fhC33[i] = 0;
457 fhC44[i] = 0;
458 fhC55[i] = 0;
459
460 fhRel1PtUncertainty[i] = 0;
461
462 fhDXY[i] = 0;
463 fhDZ[i] = 0;
464 fhDXYDZ[i] = 0;
465 fhDXYvsDZ[i] = 0;
466
467 fhDXYNormalized[i] = 0;
468 fhDZNormalized[i] = 0;
469 fhDXYvsDZNormalized[i] = 0;
470 fhNSigmaToVertex[i] = 0;
471
472 fhPt[i] = 0;
473 fhEta[i] = 0;
474 fhTOFdistance[i] = 0;
475 }
476 ffDTheoretical = 0;
477
478 fhCutStatistics = 0;
479 fhCutCorrelation = 0;
480}
481
482//_____________________________________________________________________________
483AliESDtrackCuts &AliESDtrackCuts::operator=(const AliESDtrackCuts &c)
484{
485 //
486 // Assignment operator
487 //
488
489 if (this != &c) ((AliESDtrackCuts &) c).Copy(*this);
490 return *this;
491}
492
493//_____________________________________________________________________________
494void AliESDtrackCuts::Copy(TObject &c) const
495{
496 //
497 // Copy function
498 //
499
500 AliESDtrackCuts& target = (AliESDtrackCuts &) c;
501
502 target.Init();
503
504 target.fCutMinNClusterTPC = fCutMinNClusterTPC;
505 target.fCutMinNClusterITS = fCutMinNClusterITS;
506 target.fCutMinNCrossedRowsTPC = fCutMinNCrossedRowsTPC;
507 target.fCutMinRatioCrossedRowsOverFindableClustersTPC = fCutMinRatioCrossedRowsOverFindableClustersTPC;
508 if(f1CutMinNClustersTPCPtDep){
509 target.f1CutMinNClustersTPCPtDep = (TFormula*) f1CutMinNClustersTPCPtDep->Clone("f1CutMinNClustersTPCPtDep");
510 }
511 target.fCutMaxPtDepNClustersTPC = fCutMaxPtDepNClustersTPC;
512 target.fCutMinLengthActiveVolumeTPC = fCutMinLengthActiveVolumeTPC;
513
514 target.fCutMaxChi2PerClusterTPC = fCutMaxChi2PerClusterTPC;
515 target.fCutMaxChi2PerClusterITS = fCutMaxChi2PerClusterITS;
516 target.fCutMaxChi2TPCConstrainedVsGlobal = fCutMaxChi2TPCConstrainedVsGlobal;
517 target.fCutMaxChi2TPCConstrainedVsGlobalVertexType = fCutMaxChi2TPCConstrainedVsGlobalVertexType;
518 target.fCutMaxMissingITSPoints = fCutMaxMissingITSPoints;
519
520 for (Int_t i = 0; i < 3; i++)
521 target.fCutClusterRequirementITS[i] = fCutClusterRequirementITS[i];
522
523 target.fCutMaxC11 = fCutMaxC11;
524 target.fCutMaxC22 = fCutMaxC22;
525 target.fCutMaxC33 = fCutMaxC33;
526 target.fCutMaxC44 = fCutMaxC44;
527 target.fCutMaxC55 = fCutMaxC55;
528
529 target.fCutMaxRel1PtUncertainty = fCutMaxRel1PtUncertainty;
530
531 target.fCutAcceptKinkDaughters = fCutAcceptKinkDaughters;
532 target.fCutAcceptSharedTPCClusters = fCutAcceptSharedTPCClusters;
533 target.fCutMaxFractionSharedTPCClusters = fCutMaxFractionSharedTPCClusters;
534 target.fCutRequireTPCRefit = fCutRequireTPCRefit;
535 target.fCutRequireTPCStandAlone = fCutRequireTPCStandAlone;
536 target.fCutRequireITSRefit = fCutRequireITSRefit;
537 target.fCutRequireITSPid = fCutRequireITSPid;
538 target.fCutRequireITSStandAlone = fCutRequireITSStandAlone;
539 target.fCutRequireITSpureSA = fCutRequireITSpureSA;
540
541 target.fCutNsigmaToVertex = fCutNsigmaToVertex;
542 target.fCutSigmaToVertexRequired = fCutSigmaToVertexRequired;
543 target.fCutMaxDCAToVertexXY = fCutMaxDCAToVertexXY;
544 target.fCutMaxDCAToVertexZ = fCutMaxDCAToVertexZ;
545 target.fCutDCAToVertex2D = fCutDCAToVertex2D;
546 target.fCutMinDCAToVertexXY = fCutMinDCAToVertexXY;
547 target.fCutMinDCAToVertexZ = fCutMinDCAToVertexZ;
548
549 target.fCutMaxDCAToVertexXYPtDep = fCutMaxDCAToVertexXYPtDep;
550 if(fCutMaxDCAToVertexXYPtDep.Length()>0)target.SetMaxDCAToVertexXYPtDep(fCutMaxDCAToVertexXYPtDep.Data());
551
552 target.fCutMaxDCAToVertexZPtDep = fCutMaxDCAToVertexZPtDep;
553 if(fCutMaxDCAToVertexZPtDep.Length()>0)target.SetMaxDCAToVertexZPtDep(fCutMaxDCAToVertexZPtDep.Data());
554
555 target.fCutMinDCAToVertexXYPtDep = fCutMinDCAToVertexXYPtDep;
556 if(fCutMinDCAToVertexXYPtDep.Length()>0)target.SetMinDCAToVertexXYPtDep(fCutMinDCAToVertexXYPtDep.Data());
557
558 target.fCutMinDCAToVertexZPtDep = fCutMinDCAToVertexZPtDep;
559 if(fCutMinDCAToVertexZPtDep.Length()>0)target.SetMinDCAToVertexZPtDep(fCutMinDCAToVertexZPtDep.Data());
560
561 target.fPMin = fPMin;
562 target.fPMax = fPMax;
563 target.fPtMin = fPtMin;
564 target.fPtMax = fPtMax;
565 target.fPxMin = fPxMin;
566 target.fPxMax = fPxMax;
567 target.fPyMin = fPyMin;
568 target.fPyMax = fPyMax;
569 target.fPzMin = fPzMin;
570 target.fPzMax = fPzMax;
571 target.fEtaMin = fEtaMin;
572 target.fEtaMax = fEtaMax;
573 target.fRapMin = fRapMin;
574 target.fRapMax = fRapMax;
575
576 target.fFlagCutTOFdistance = fFlagCutTOFdistance;
577 target.fCutTOFdistance = fCutTOFdistance;
578 target.fCutRequireTOFout = fCutRequireTOFout;
579
580 target.fHistogramsOn = fHistogramsOn;
581
582 for (Int_t i=0; i<2; ++i)
583 {
584 if (fhNClustersITS[i]) target.fhNClustersITS[i] = (TH1F*) fhNClustersITS[i]->Clone();
585 if (fhNClustersTPC[i]) target.fhNClustersTPC[i] = (TH1F*) fhNClustersTPC[i]->Clone();
586 if (fhNSharedClustersTPC[i]) target.fhNSharedClustersTPC[i] = (TH1F*) fhNSharedClustersTPC[i]->Clone();
587 if (fhNCrossedRowsTPC[i]) target.fhNCrossedRowsTPC[i] = (TH1F*) fhNCrossedRowsTPC[i]->Clone();
588 if (fhRatioCrossedRowsOverFindableClustersTPC[i]) target.fhRatioCrossedRowsOverFindableClustersTPC[i] = (TH1F*) fhRatioCrossedRowsOverFindableClustersTPC[i]->Clone();
589
590 if (fhChi2PerClusterITS[i]) target.fhChi2PerClusterITS[i] = (TH1F*) fhChi2PerClusterITS[i]->Clone();
591 if (fhChi2PerClusterTPC[i]) target.fhChi2PerClusterTPC[i] = (TH1F*) fhChi2PerClusterTPC[i]->Clone();
592 if (fhChi2TPCConstrainedVsGlobal[i]) target.fhChi2TPCConstrainedVsGlobal[i] = (TH1F*) fhChi2TPCConstrainedVsGlobal[i]->Clone();
593 if (fhNClustersForITSPID[i]) target.fhNClustersForITSPID[i] = (TH1F*) fhNClustersForITSPID[i]->Clone();
594 if (fhNMissingITSPoints[i]) target.fhNMissingITSPoints[i] = (TH1F*) fhNMissingITSPoints[i]->Clone();
595
596 if (fhC11[i]) target.fhC11[i] = (TH1F*) fhC11[i]->Clone();
597 if (fhC22[i]) target.fhC22[i] = (TH1F*) fhC22[i]->Clone();
598 if (fhC33[i]) target.fhC33[i] = (TH1F*) fhC33[i]->Clone();
599 if (fhC44[i]) target.fhC44[i] = (TH1F*) fhC44[i]->Clone();
600 if (fhC55[i]) target.fhC55[i] = (TH1F*) fhC55[i]->Clone();
601
602 if (fhRel1PtUncertainty[i]) target.fhRel1PtUncertainty[i] = (TH1F*) fhRel1PtUncertainty[i]->Clone();
603
604 if (fhDXY[i]) target.fhDXY[i] = (TH1F*) fhDXY[i]->Clone();
605 if (fhDZ[i]) target.fhDZ[i] = (TH1F*) fhDZ[i]->Clone();
606 if (fhDXYDZ[i]) target.fhDXYDZ[i] = (TH1F*) fhDXYDZ[i]->Clone();
607 if (fhDXYvsDZ[i]) target.fhDXYvsDZ[i] = (TH2F*) fhDXYvsDZ[i]->Clone();
608
609 if (fhDXYNormalized[i]) target.fhDXYNormalized[i] = (TH1F*) fhDXYNormalized[i]->Clone();
610 if (fhDZNormalized[i]) target.fhDZNormalized[i] = (TH1F*) fhDZNormalized[i]->Clone();
611 if (fhDXYvsDZNormalized[i]) target.fhDXYvsDZNormalized[i] = (TH2F*) fhDXYvsDZNormalized[i]->Clone();
612 if (fhNSigmaToVertex[i]) target.fhNSigmaToVertex[i] = (TH1F*) fhNSigmaToVertex[i]->Clone();
613
614 if (fhPt[i]) target.fhPt[i] = (TH1F*) fhPt[i]->Clone();
615 if (fhEta[i]) target.fhEta[i] = (TH1F*) fhEta[i]->Clone();
616 if (fhTOFdistance[i]) target.fhTOFdistance[i] = (TH2F*) fhTOFdistance[i]->Clone();
617 }
618 if (ffDTheoretical) target.ffDTheoretical = (TF1*) ffDTheoretical->Clone();
619
620 if (fhCutStatistics) target.fhCutStatistics = (TH1F*) fhCutStatistics->Clone();
621 if (fhCutCorrelation) target.fhCutCorrelation = (TH2F*) fhCutCorrelation->Clone();
622
623 TNamed::Copy(c);
624}
625
626//_____________________________________________________________________________
627Long64_t AliESDtrackCuts::Merge(TCollection* list) {
628 // Merge a list of AliESDtrackCuts objects with this (needed for PROOF)
629 // Returns the number of merged objects (including this)
630 if (!list)
631 return 0;
632 if (list->IsEmpty())
633 return 1;
634 if (!fHistogramsOn)
635 return 0;
636 TIterator* iter = list->MakeIterator();
637 TObject* obj;
638
639 // collection of measured and generated histograms
640 Int_t count = 0;
641 while ((obj = iter->Next())) {
642
643 AliESDtrackCuts* entry = dynamic_cast<AliESDtrackCuts*>(obj);
644 if (entry == 0)
645 continue;
646
647 if (!entry->fHistogramsOn)
648 continue;
649
650 for (Int_t i=0; i<2; i++) {
651
652 fhNClustersITS[i] ->Add(entry->fhNClustersITS[i] );
653 fhNClustersTPC[i] ->Add(entry->fhNClustersTPC[i] );
654 if (fhNSharedClustersTPC[i])
655 fhNSharedClustersTPC[i] ->Add(entry->fhNSharedClustersTPC[i] );
656 if (fhNCrossedRowsTPC[i])
657 fhNCrossedRowsTPC[i] ->Add(entry->fhNCrossedRowsTPC[i] );
658 if (fhRatioCrossedRowsOverFindableClustersTPC[i])
659 fhRatioCrossedRowsOverFindableClustersTPC[i] ->Add(entry->fhRatioCrossedRowsOverFindableClustersTPC[i] );
660
661 fhChi2PerClusterITS[i] ->Add(entry->fhChi2PerClusterITS[i]);
662 fhChi2PerClusterTPC[i] ->Add(entry->fhChi2PerClusterTPC[i]);
663 if (fhChi2TPCConstrainedVsGlobal[i])
664 fhChi2TPCConstrainedVsGlobal[i]->Add(entry->fhChi2TPCConstrainedVsGlobal[i]);
665 if (fhNClustersForITSPID[i])
666 fhNClustersForITSPID[i]->Add(entry->fhNClustersForITSPID[i]);
667 if (fhNMissingITSPoints[i])
668 fhNMissingITSPoints[i] ->Add(entry->fhNMissingITSPoints[i]);
669
670 fhC11[i] ->Add(entry->fhC11[i] );
671 fhC22[i] ->Add(entry->fhC22[i] );
672 fhC33[i] ->Add(entry->fhC33[i] );
673 fhC44[i] ->Add(entry->fhC44[i] );
674 fhC55[i] ->Add(entry->fhC55[i] );
675
676 fhRel1PtUncertainty[i] ->Add(entry->fhRel1PtUncertainty[i]);
677
678 fhDXY[i] ->Add(entry->fhDXY[i] );
679 fhDZ[i] ->Add(entry->fhDZ[i] );
680 fhDXYDZ[i] ->Add(entry->fhDXYDZ[i] );
681 fhDXYvsDZ[i] ->Add(entry->fhDXYvsDZ[i] );
682
683 fhDXYNormalized[i] ->Add(entry->fhDXYNormalized[i] );
684 fhDZNormalized[i] ->Add(entry->fhDZNormalized[i] );
685 fhDXYvsDZNormalized[i] ->Add(entry->fhDXYvsDZNormalized[i]);
686 fhNSigmaToVertex[i] ->Add(entry->fhNSigmaToVertex[i]);
687
688 fhPt[i] ->Add(entry->fhPt[i]);
689 fhEta[i] ->Add(entry->fhEta[i]);
690 fhTOFdistance[i] ->Add(entry->fhTOFdistance[i]);
691 }
692
693 fhCutStatistics ->Add(entry->fhCutStatistics);
694 fhCutCorrelation ->Add(entry->fhCutCorrelation);
695
696 count++;
697 }
698 return count+1;
699}
700
701void AliESDtrackCuts::SetMinNClustersTPCPtDep(TFormula *f1, Float_t ptmax)
702{
703 //
704 // Sets the pT dependent NClustersTPC cut
705 //
706
707 if(f1){
708 delete f1CutMinNClustersTPCPtDep;
709 f1CutMinNClustersTPCPtDep = (TFormula*)f1->Clone("f1CutMinNClustersTPCPtDep");
710 }
711 fCutMaxPtDepNClustersTPC=ptmax;
712}
713
714//____________________________________________________________________
715AliESDtrackCuts* AliESDtrackCuts::GetStandardTPCOnlyTrackCuts()
716{
717 // creates an AliESDtrackCuts object and fills it with standard (pre data-taking) values for TPC-only cuts
718
719 AliInfoClass("Creating track cuts for TPC-only.");
720
721 AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
722
723 esdTrackCuts->SetMinNClustersTPC(50);
724 esdTrackCuts->SetMaxChi2PerClusterTPC(4);
725 esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
726
727 esdTrackCuts->SetMaxDCAToVertexZ(3.2);
728 esdTrackCuts->SetMaxDCAToVertexXY(2.4);
729 esdTrackCuts->SetDCAToVertex2D(kTRUE);
730
731 return esdTrackCuts;
732}
733
734//____________________________________________________________________
735AliESDtrackCuts* AliESDtrackCuts::GetStandardITSTPCTrackCuts2009(Bool_t selPrimaries)
736{
737 // creates an AliESDtrackCuts object and fills it with standard values for ITS-TPC cuts for pp 2009 data
738
739 AliInfoClass("Creating track cuts for ITS+TPC (2009 definition).");
740
741 AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
742
743 // TPC
744 esdTrackCuts->SetRequireTPCStandAlone(kTRUE); // to get chi2 and ncls of kTPCin
745 esdTrackCuts->SetMinNClustersTPC(70);
746 esdTrackCuts->SetMaxChi2PerClusterTPC(4);
747 esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
748 esdTrackCuts->SetRequireTPCRefit(kTRUE);
749 // ITS
750 esdTrackCuts->SetRequireITSRefit(kTRUE);
751 esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
752 AliESDtrackCuts::kAny);
753 if(selPrimaries) {
754 // 7*(0.0050+0.0060/pt^0.9)
755 esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0350+0.0420/pt^0.9");
756 esdTrackCuts->SetMaxChi2TPCConstrainedGlobal(36);
757 }
758 esdTrackCuts->SetMaxDCAToVertexZ(1.e6);
759 esdTrackCuts->SetDCAToVertex2D(kFALSE);
760 esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
761 //esdTrackCuts->SetEtaRange(-0.8,+0.8);
762
763 esdTrackCuts->SetMaxChi2PerClusterITS(36);
764
765 return esdTrackCuts;
766}
767
768//____________________________________________________________________
769AliESDtrackCuts* AliESDtrackCuts::GetStandardITSTPCTrackCuts2011(Bool_t selPrimaries, Int_t clusterCut)
770{
771 // creates an AliESDtrackCuts object and fills it with standard values for ITS-TPC cuts for pp 2011 data
772 // if clusterCut = 1, the cut on the number of clusters is replaced by
773 // a cut on the number of crossed rows and on the ration crossed
774 // rows/findable clusters
775
776 AliInfoClass("Creating track cuts for ITS+TPC (2011 definition).");
777
778 AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
779
780 // TPC
781 if(clusterCut == 0) esdTrackCuts->SetMinNClustersTPC(50);
782 else if (clusterCut == 1) {
783 esdTrackCuts->SetMinNCrossedRowsTPC(70);
784 esdTrackCuts->SetMinRatioCrossedRowsOverFindableClustersTPC(0.8);
785 }
786 else {
787 AliWarningClass(Form("Wrong value of the clusterCut parameter (%d), using cut on Nclusters",clusterCut));
788 esdTrackCuts->SetMinNClustersTPC(50);
789 }
790 esdTrackCuts->SetMaxChi2PerClusterTPC(4);
791 esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
792 esdTrackCuts->SetRequireTPCRefit(kTRUE);
793 // ITS
794 esdTrackCuts->SetRequireITSRefit(kTRUE);
795 esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
796 AliESDtrackCuts::kAny);
797 if(selPrimaries) {
798 // 7*(0.0015+0.0050/pt^1.1)
799 esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0105+0.0350/pt^1.1");
800 esdTrackCuts->SetMaxChi2TPCConstrainedGlobal(36);
801 }
802 esdTrackCuts->SetMaxDCAToVertexZ(2);
803 esdTrackCuts->SetDCAToVertex2D(kFALSE);
804 esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
805
806 esdTrackCuts->SetMaxChi2PerClusterITS(36);
807
808 return esdTrackCuts;
809}
810
811//____________________________________________________________________
812AliESDtrackCuts* AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(Bool_t selPrimaries,Int_t clusterCut)
813{
814 // creates an AliESDtrackCuts object and fills it with standard values for ITS-TPC cuts for pp 2010 data
815 // if clusterCut = 1, the cut on the number of clusters is replaced by
816 // a cut on the number of crossed rows and on the ration crossed
817 // rows/findable clusters
818
819 AliInfoClass("Creating track cuts for ITS+TPC (2010 definition).");
820
821 AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
822
823 // TPC
824 if(clusterCut == 0) esdTrackCuts->SetMinNClustersTPC(70);
825 else if (clusterCut == 1) {
826 esdTrackCuts->SetMinNCrossedRowsTPC(70);
827 esdTrackCuts->SetMinRatioCrossedRowsOverFindableClustersTPC(0.8);
828 }
829 else {
830 AliWarningClass(Form("Wrong value of the clusterCut parameter (%d), using cut on Nclusters",clusterCut));
831 esdTrackCuts->SetMinNClustersTPC(70);
832 }
833 esdTrackCuts->SetMaxChi2PerClusterTPC(4);
834 esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
835 esdTrackCuts->SetRequireTPCRefit(kTRUE);
836 // ITS
837 esdTrackCuts->SetRequireITSRefit(kTRUE);
838 esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
839 AliESDtrackCuts::kAny);
840 if(selPrimaries) {
841 // 7*(0.0026+0.0050/pt^1.01)
842 esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0182+0.0350/pt^1.01");
843 esdTrackCuts->SetMaxChi2TPCConstrainedGlobal(36);
844 }
845 esdTrackCuts->SetMaxDCAToVertexZ(2);
846 esdTrackCuts->SetDCAToVertex2D(kFALSE);
847 esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
848
849 esdTrackCuts->SetMaxChi2PerClusterITS(36);
850
851 return esdTrackCuts;
852}
853
854//____________________________________________________________________
855AliESDtrackCuts* AliESDtrackCuts::GetStandardITSPureSATrackCuts2009(Bool_t selPrimaries, Bool_t useForPid)
856{
857 // creates an AliESDtrackCuts object and fills it with standard values for ITS pure SA tracks
858
859 AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
860 esdTrackCuts->SetRequireITSPureStandAlone(kTRUE);
861 esdTrackCuts->SetRequireITSRefit(kTRUE);
862 esdTrackCuts->SetMinNClustersITS(4);
863 esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
864 AliESDtrackCuts::kAny);
865 esdTrackCuts->SetMaxChi2PerClusterITS(1.);
866
867 if(selPrimaries) {
868 // 7*(0.0085+0.0026/pt^1.55)
869 esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0595+0.0182/pt^1.55");
870 }
871 if(useForPid){
872 esdTrackCuts->SetRequireITSPid(kTRUE);
873 }
874 return esdTrackCuts;
875}
876
877//____________________________________________________________________
878AliESDtrackCuts* AliESDtrackCuts::GetStandardITSPureSATrackCuts2010(Bool_t selPrimaries, Bool_t useForPid)
879{
880 // creates an AliESDtrackCuts object and fills it with standard values for ITS pure SA tracks - pp 2010
881
882 AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
883 esdTrackCuts->SetRequireITSPureStandAlone(kTRUE);
884 esdTrackCuts->SetRequireITSRefit(kTRUE);
885 esdTrackCuts->SetMinNClustersITS(4);
886 esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
887 AliESDtrackCuts::kAny);
888 esdTrackCuts->SetMaxChi2PerClusterITS(2.5);
889
890 if(selPrimaries) {
891 // 7*(0.0033+0.0045/pt^1.3)
892 esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0231+0.0315/pt^1.3");
893 }
894 if(useForPid){
895 esdTrackCuts->SetRequireITSPid(kTRUE);
896 }
897 return esdTrackCuts;
898}
899
900//____________________________________________________________________
901AliESDtrackCuts* AliESDtrackCuts::GetStandardITSSATrackCuts2009(Bool_t selPrimaries, Bool_t useForPid)
902{
903 // creates an AliESDtrackCuts object and fills it with standard values for ITS pure SA tracks
904
905 AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
906 esdTrackCuts->SetRequireITSStandAlone(kTRUE);
907 esdTrackCuts->SetRequireITSPureStandAlone(kFALSE);
908 esdTrackCuts->SetRequireITSRefit(kTRUE);
909 esdTrackCuts->SetMinNClustersITS(4);
910 esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
911 AliESDtrackCuts::kAny);
912 esdTrackCuts->SetMaxChi2PerClusterITS(1.);
913
914 if(selPrimaries) {
915 // 7*(0.0085+0.0026/pt^1.55)
916 esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0595+0.0182/pt^1.55");
917 }
918 if(useForPid){
919 esdTrackCuts->SetRequireITSPid(kTRUE);
920 }
921 return esdTrackCuts;
922}
923
924//____________________________________________________________________
925AliESDtrackCuts* AliESDtrackCuts::GetStandardITSSATrackCuts2010(Bool_t selPrimaries, Bool_t useForPid)
926{
927 // creates an AliESDtrackCuts object and fills it with standard values for ITS pure SA tracks --pp 2010
928
929 AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
930 esdTrackCuts->SetRequireITSStandAlone(kTRUE);
931 esdTrackCuts->SetRequireITSPureStandAlone(kFALSE);
932 esdTrackCuts->SetRequireITSRefit(kTRUE);
933 esdTrackCuts->SetMinNClustersITS(4);
934 esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
935 AliESDtrackCuts::kAny);
936 esdTrackCuts->SetMaxChi2PerClusterITS(2.5);
937
938 if(selPrimaries) {
939 // 7*(0.0033+0.0045/pt^1.3)
940 esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0231+0.0315/pt^1.3");
941 }
942 if(useForPid){
943 esdTrackCuts->SetRequireITSPid(kTRUE);
944 }
945 return esdTrackCuts;
946}
947
948//____________________________________________________________________
949AliESDtrackCuts* AliESDtrackCuts::GetStandardITSSATrackCutsPbPb2010(Bool_t selPrimaries, Bool_t useForPid)
950{
951 // creates an AliESDtrackCuts object and fills it with standard values for ITS pure SA tracks -- PbPb 2010
952
953 AliESDtrackCuts* esdTrackCuts = GetStandardITSSATrackCuts2010(selPrimaries, useForPid);
954 esdTrackCuts->SetMaxNOfMissingITSPoints(1);
955
956 return esdTrackCuts;
957}
958//____________________________________________________________________
959
960AliESDtrackCuts* AliESDtrackCuts::GetStandardV0DaughterCuts()
961{
962 // creates a AliESDtrackCuts object and fills it with standard cuts for V0 daughters
963 AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
964 esdTrackCuts->SetRequireTPCRefit(kTRUE);
965 esdTrackCuts->SetMinNClustersTPC(70);
966 esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
967 return esdTrackCuts;
968}
969
970//____________________________________________________________________
971Int_t AliESDtrackCuts::GetReferenceMultiplicity(const AliESDEvent* esd, Bool_t tpcOnly)
972{
973 // Gets reference multiplicity following the standard cuts and a defined fiducial volume
974 // tpcOnly = kTRUE -> consider TPC-only tracks
975 // = kFALSE -> consider global tracks
976 //
977 // DEPRECATED Use GetReferenceMultiplicity with the enum as second argument instead
978
979 if (!tpcOnly)
980 {
981 AliErrorClass("Not implemented for global tracks!");
982 return -1;
983 }
984
985 static AliESDtrackCuts* esdTrackCuts = 0;
986 if (!esdTrackCuts)
987 {
988 esdTrackCuts = GetStandardTPCOnlyTrackCuts();
989 esdTrackCuts->SetEtaRange(-0.8, 0.8);
990 esdTrackCuts->SetPtRange(0.15);
991 }
992
993 Int_t nTracks = esdTrackCuts->CountAcceptedTracks(esd);
994
995 return nTracks;
996}
997
998//____________________________________________________________________
999Float_t AliESDtrackCuts::GetSigmaToVertex(const AliESDtrack* const esdTrack)
1000{
1001 // Calculates the number of sigma to the vertex.
1002
1003 Float_t b[2];
1004 Float_t bRes[2];
1005 Float_t bCov[3];
1006 esdTrack->GetImpactParameters(b,bCov);
1007
1008 if (bCov[0]<=0 || bCov[2]<=0) {
1009 AliDebugClass(1, "Estimated b resolution lower or equal zero!");
1010 bCov[0]=0; bCov[2]=0;
1011 }
1012 bRes[0] = TMath::Sqrt(bCov[0]);
1013 bRes[1] = TMath::Sqrt(bCov[2]);
1014
1015 // -----------------------------------
1016 // How to get to a n-sigma cut?
1017 //
1018 // The accumulated statistics from 0 to d is
1019 //
1020 // -> Erf(d/Sqrt(2)) for a 1-dim gauss (d = n_sigma)
1021 // -> 1 - Exp(-d**2) for a 2-dim gauss (d*d = dx*dx + dy*dy != n_sigma)
1022 //
1023 // It means that for a 2-dim gauss: n_sigma(d) = Sqrt(2)*ErfInv(1 - Exp((-d**2)/2)
1024 // Can this be expressed in a different way?
1025
1026 if (bRes[0] == 0 || bRes[1] ==0)
1027 return -1;
1028
1029 Float_t d = TMath::Sqrt(TMath::Power(b[0]/bRes[0],2) + TMath::Power(b[1]/bRes[1],2));
1030
1031 // work around precision problem
1032 // if d is too big, TMath::Exp(...) gets 0, and TMath::ErfInverse(1) that should be infinite, gets 0 :(
1033 // 1e-15 corresponds to nsigma ~ 7.7
1034 if (TMath::Exp(-d * d / 2) < 1e-15)
1035 return 1000;
1036
1037 Float_t nSigma = TMath::ErfInverse(1 - TMath::Exp(-d * d / 2)) * TMath::Sqrt(2);
1038 return nSigma;
1039}
1040
1041void AliESDtrackCuts::EnableNeededBranches(TTree* tree)
1042{
1043 // enables the branches needed by AcceptTrack, for a list see comment of AcceptTrack
1044
1045 tree->SetBranchStatus("fTracks.fFlags", 1);
1046 tree->SetBranchStatus("fTracks.fITSncls", 1);
1047 tree->SetBranchStatus("fTracks.fTPCncls", 1);
1048 tree->SetBranchStatus("fTracks.fITSchi2", 1);
1049 tree->SetBranchStatus("fTracks.fTPCchi2", 1);
1050 tree->SetBranchStatus("fTracks.fC*", 1);
1051 tree->SetBranchStatus("fTracks.fD", 1);
1052 tree->SetBranchStatus("fTracks.fZ", 1);
1053 tree->SetBranchStatus("fTracks.fCdd", 1);
1054 tree->SetBranchStatus("fTracks.fCdz", 1);
1055 tree->SetBranchStatus("fTracks.fCzz", 1);
1056 tree->SetBranchStatus("fTracks.fP*", 1);
1057 tree->SetBranchStatus("fTracks.fR*", 1);
1058 tree->SetBranchStatus("fTracks.fKinkIndexes*", 1);
1059}
1060
1061//____________________________________________________________________
1062Bool_t AliESDtrackCuts::AcceptTrack(const AliESDtrack* esdTrack)
1063{
1064 //
1065 // figure out if the tracks survives all the track cuts defined
1066 //
1067 // the different quality parameter and kinematic values are first
1068 // retrieved from the track. then it is found out what cuts the
1069 // track did not survive and finally the cuts are imposed.
1070
1071 // this function needs the following branches:
1072 // fTracks.fFlags
1073 // fTracks.fITSncls
1074 // fTracks.fTPCncls
1075 // fTracks.fITSchi2
1076 // fTracks.fTPCchi2
1077 // fTracks.fC //GetExternalCovariance
1078 // fTracks.fD //GetImpactParameters
1079 // fTracks.fZ //GetImpactParameters
1080 // fTracks.fCdd //GetImpactParameters
1081 // fTracks.fCdz //GetImpactParameters
1082 // fTracks.fCzz //GetImpactParameters
1083 // fTracks.fP //GetPxPyPz
1084 // fTracks.fR //GetMass
1085 // fTracks.fP //GetMass
1086 // fTracks.fKinkIndexes
1087 //
1088 // esdEvent is only required for the MaxChi2TPCConstrainedVsGlobal
1089
1090 UInt_t status = esdTrack->GetStatus();
1091
1092 // getting quality parameters from the ESD track
1093 Int_t nClustersITS = esdTrack->GetITSclusters(0);
1094 Int_t nClustersTPC = -1;
1095 if(fCutRequireTPCStandAlone) {
1096 nClustersTPC = esdTrack->GetTPCNclsIter1();
1097 }
1098 else {
1099 nClustersTPC = esdTrack->GetTPCclusters(0);
1100 }
1101
1102 //Pt dependent NClusters Cut
1103 if(f1CutMinNClustersTPCPtDep) {
1104 if(esdTrack->Pt()<fCutMaxPtDepNClustersTPC)
1105 fCutMinNClusterTPC = (Int_t)(f1CutMinNClustersTPCPtDep->Eval(esdTrack->Pt()));
1106 else
1107 fCutMinNClusterTPC = (Int_t)(f1CutMinNClustersTPCPtDep->Eval(fCutMaxPtDepNClustersTPC));
1108 }
1109
1110 Float_t nCrossedRowsTPC = esdTrack->GetTPCCrossedRows();
1111 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
1112 if (esdTrack->GetTPCNclsF()>0) {
1113 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC / esdTrack->GetTPCNclsF();
1114 }
1115
1116 Int_t nClustersTPCShared = esdTrack->GetTPCnclsS();
1117 Float_t fracClustersTPCShared = -1.;
1118
1119 Float_t chi2PerClusterITS = -1;
1120 Float_t chi2PerClusterTPC = -1;
1121 if (nClustersITS!=0)
1122 chi2PerClusterITS = esdTrack->GetITSchi2()/Float_t(nClustersITS);
1123 if (nClustersTPC!=0) {
1124 if(fCutRequireTPCStandAlone) {
1125 chi2PerClusterTPC = esdTrack->GetTPCchi2Iter1()/Float_t(nClustersTPC);
1126 } else {
1127 chi2PerClusterTPC = esdTrack->GetTPCchi2()/Float_t(nClustersTPC);
1128 }
1129 fracClustersTPCShared = Float_t(nClustersTPCShared)/Float_t(nClustersTPC);
1130 }
1131
1132 Double_t extCov[15];
1133 esdTrack->GetExternalCovariance(extCov);
1134
1135 Float_t b[2];
1136 Float_t bCov[3];
1137 esdTrack->GetImpactParameters(b,bCov);
1138 if (bCov[0]<=0 || bCov[2]<=0) {
1139 AliDebug(1, "Estimated b resolution lower or equal zero!");
1140 bCov[0]=0; bCov[2]=0;
1141 }
1142
1143
1144 // set pt-dependent DCA cuts, if requested
1145 SetPtDepDCACuts(esdTrack->Pt());
1146
1147
1148 Float_t dcaToVertexXY = b[0];
1149 Float_t dcaToVertexZ = b[1];
1150
1151 Float_t dcaToVertex = -1;
1152
1153 if (fCutDCAToVertex2D)
1154 {
1155 dcaToVertex = TMath::Sqrt(dcaToVertexXY*dcaToVertexXY/fCutMaxDCAToVertexXY/fCutMaxDCAToVertexXY + dcaToVertexZ*dcaToVertexZ/fCutMaxDCAToVertexZ/fCutMaxDCAToVertexZ);
1156 }
1157 else
1158 dcaToVertex = TMath::Sqrt(dcaToVertexXY*dcaToVertexXY + dcaToVertexZ*dcaToVertexZ);
1159
1160 // getting the kinematic variables of the track
1161 // (assuming the mass is known)
1162 Double_t p[3];
1163 esdTrack->GetPxPyPz(p);
1164
1165 // Changed from float to double to prevent rounding errors leading to negative
1166 // log arguments (M.G.)
1167 Double_t momentum = TMath::Sqrt(p[0]*p[0] + p[1]*p[1] + p[2]*p[2]);
1168 Double_t pt = TMath::Sqrt(p[0]*p[0] + p[1]*p[1]);
1169 Double_t mass = esdTrack->GetMass();
1170 Double_t energy = TMath::Sqrt(mass*mass + momentum*momentum);
1171
1172 //y-eta related calculations
1173 Float_t eta = -100.;
1174 Float_t y = -100.;
1175 if((momentum != TMath::Abs(p[2]))&&(momentum != 0))
1176 eta = 0.5*TMath::Log((momentum + p[2])/(momentum - p[2]));
1177 if((energy != TMath::Abs(p[2]))&&(momentum != 0))
1178 y = 0.5*TMath::Log((energy + p[2])/(energy - p[2]));
1179
1180 if (extCov[14] < 0)
1181 {
1182 AliWarning(Form("GetSigma1Pt2() returns negative value for external covariance matrix element fC[14]: %f. Corrupted track information, track will not be accepted!", extCov[14]));
1183 return kFALSE;
1184 }
1185 Float_t relUncertainty1Pt = TMath::Sqrt(extCov[14])*pt;
1186
1187 //########################################################################
1188 // cut the track?
1189
1190 Bool_t cuts[kNCuts];
1191 for (Int_t i=0; i<kNCuts; i++) cuts[i]=kFALSE;
1192
1193 // track quality cuts
1194 if (fCutRequireTPCRefit && (status&AliESDtrack::kTPCrefit)==0)
1195 cuts[0]=kTRUE;
1196 if (fCutRequireTPCStandAlone && (status&AliESDtrack::kTPCin)==0)
1197 cuts[1]=kTRUE;
1198 if (fCutRequireITSRefit && (status&AliESDtrack::kITSrefit)==0)
1199 cuts[2]=kTRUE;
1200 if (nClustersTPC<fCutMinNClusterTPC)
1201 cuts[3]=kTRUE;
1202 if (nClustersITS<fCutMinNClusterITS)
1203 cuts[4]=kTRUE;
1204 if (chi2PerClusterTPC>fCutMaxChi2PerClusterTPC)
1205 cuts[5]=kTRUE;
1206 if (chi2PerClusterITS>fCutMaxChi2PerClusterITS)
1207 cuts[6]=kTRUE;
1208 if (extCov[0] > fCutMaxC11)
1209 cuts[7]=kTRUE;
1210 if (extCov[2] > fCutMaxC22)
1211 cuts[8]=kTRUE;
1212 if (extCov[5] > fCutMaxC33)
1213 cuts[9]=kTRUE;
1214 if (extCov[9] > fCutMaxC44)
1215 cuts[10]=kTRUE;
1216 if (extCov[14] > fCutMaxC55)
1217 cuts[11]=kTRUE;
1218
1219 // cut 12 and 13 see below
1220
1221 if (!fCutAcceptKinkDaughters && esdTrack->GetKinkIndex(0)>0)
1222 cuts[14]=kTRUE;
1223 // track kinematics cut
1224 if((momentum < fPMin) || (momentum > fPMax))
1225 cuts[15]=kTRUE;
1226 if((pt < fPtMin) || (pt > fPtMax))
1227 cuts[16] = kTRUE;
1228 if((p[0] < fPxMin) || (p[0] > fPxMax))
1229 cuts[17] = kTRUE;
1230 if((p[1] < fPyMin) || (p[1] > fPyMax))
1231 cuts[18] = kTRUE;
1232 if((p[2] < fPzMin) || (p[2] > fPzMax))
1233 cuts[19] = kTRUE;
1234 if((eta < fEtaMin) || (eta > fEtaMax))
1235 cuts[20] = kTRUE;
1236 if((y < fRapMin) || (y > fRapMax))
1237 cuts[21] = kTRUE;
1238 if (fCutDCAToVertex2D && dcaToVertex > 1)
1239 cuts[22] = kTRUE;
1240 if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexXY) > fCutMaxDCAToVertexXY)
1241 cuts[23] = kTRUE;
1242 if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexZ) > fCutMaxDCAToVertexZ)
1243 cuts[24] = kTRUE;
1244 if (fCutDCAToVertex2D && fCutMinDCAToVertexXY > 0 && fCutMinDCAToVertexZ > 0 && dcaToVertexXY*dcaToVertexXY/fCutMinDCAToVertexXY/fCutMinDCAToVertexXY + dcaToVertexZ*dcaToVertexZ/fCutMinDCAToVertexZ/fCutMinDCAToVertexZ < 1)
1245 cuts[25] = kTRUE;
1246 if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexXY) < fCutMinDCAToVertexXY)
1247 cuts[26] = kTRUE;
1248 if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexZ) < fCutMinDCAToVertexZ)
1249 cuts[27] = kTRUE;
1250
1251 for (Int_t i = 0; i < 3; i++) {
1252 if(!(esdTrack->GetStatus()&AliESDtrack::kITSupg)) { // current ITS
1253 cuts[28+i] = !CheckITSClusterRequirement(fCutClusterRequirementITS[i], esdTrack->HasPointOnITSLayer(i*2), esdTrack->HasPointOnITSLayer(i*2+1));
1254 } else { // upgraded ITS (7 layers)
1255 // at the moment, for L012 the layers 12 are considered together
1256 if(i==0) { // L012
1257 cuts[28+i] = !CheckITSClusterRequirement(fCutClusterRequirementITS[i], esdTrack->HasPointOnITSLayer(0), (esdTrack->HasPointOnITSLayer(1))&(esdTrack->HasPointOnITSLayer(2)));
1258 } else { // L34 or L56
1259 cuts[28+i] = !CheckITSClusterRequirement(fCutClusterRequirementITS[i], esdTrack->HasPointOnITSLayer(i*2+1), esdTrack->HasPointOnITSLayer(i*2+2));
1260 }
1261 }
1262 }
1263
1264 if(fCutRequireITSStandAlone || fCutRequireITSpureSA){
1265 if ((status & AliESDtrack::kITSin) == 0 || (status & AliESDtrack::kTPCin)){
1266 // TPC tracks
1267 cuts[31] = kTRUE;
1268 }else{
1269 // ITS standalone tracks
1270 if(fCutRequireITSStandAlone && !fCutRequireITSpureSA){
1271 if(status & AliESDtrack::kITSpureSA) cuts[31] = kTRUE;
1272 }else if(fCutRequireITSpureSA){
1273 if(!(status & AliESDtrack::kITSpureSA)) cuts[31] = kTRUE;
1274 }
1275 }
1276 }
1277
1278 if (relUncertainty1Pt > fCutMaxRel1PtUncertainty)
1279 cuts[32] = kTRUE;
1280
1281 if (!fCutAcceptSharedTPCClusters && nClustersTPCShared!=0)
1282 cuts[33] = kTRUE;
1283
1284 if (fracClustersTPCShared > fCutMaxFractionSharedTPCClusters)
1285 cuts[34] = kTRUE;
1286
1287 Int_t nITSPointsForPid=0;
1288 UChar_t clumap=esdTrack->GetITSClusterMap();
1289 for(Int_t i=2; i<6; i++){
1290 if(clumap&(1<<i)) ++nITSPointsForPid;
1291 }
1292 if(fCutRequireITSPid && nITSPointsForPid<3) cuts[35] = kTRUE;
1293
1294
1295 if (nCrossedRowsTPC<fCutMinNCrossedRowsTPC)
1296 cuts[36]=kTRUE;
1297 if (ratioCrossedRowsOverFindableClustersTPC<fCutMinRatioCrossedRowsOverFindableClustersTPC)
1298 cuts[37]=kTRUE;
1299
1300 Int_t nMissITSpts=0;
1301 Int_t idet,statusLay;
1302 Float_t xloc,zloc;
1303 for(Int_t iLay=0; iLay<6; iLay++){
1304 Bool_t retc=esdTrack->GetITSModuleIndexInfo(iLay,idet,statusLay,xloc,zloc);
1305 if(retc && statusLay==5) ++nMissITSpts;
1306 }
1307 if(nMissITSpts>fCutMaxMissingITSPoints) cuts[38] = kTRUE;
1308
1309 //kTOFout
1310 if (fCutRequireTOFout && (status&AliESDtrack::kTOFout)==0)
1311 cuts[40]=kTRUE;
1312
1313 // TOF signal Dz cut
1314 Float_t dxTOF = esdTrack->GetTOFsignalDx();
1315 Float_t dzTOF = esdTrack->GetTOFsignalDz();
1316 if (fFlagCutTOFdistance && (esdTrack->GetStatus() & AliESDtrack::kTOFout) == AliESDtrack::kTOFout){ // applying the TOF distance cut only if requested, and only on tracks that reached the TOF and where associated with a TOF hit
1317 if (fgBeamTypeFlag < 0) { // the check on the beam type was not done yet
1318 const AliESDEvent* event = esdTrack->GetESDEvent();
1319 if (event){
1320 TString beamTypeESD = event->GetBeamType();
1321 AliDebug(2,Form("Beam type from ESD event = %s",beamTypeESD.Data()));
1322 if (beamTypeESD.CompareTo("A-A",TString::kIgnoreCase) == 0){ // we are in PbPb collisions --> fgBeamTypeFlag will be set to 1, to apply the cut on TOF signal Dz
1323 fgBeamTypeFlag = 1;
1324 }
1325 else { // we are NOT in PbPb collisions --> fgBeamTypeFlag will be set to 0, to NOT apply the cu6 on TOF signal Dz
1326 fgBeamTypeFlag = 0;
1327 }
1328 }
1329 else{
1330 AliFatal("Beam type not available, but it is needed to apply the TOF cut!");
1331 }
1332 }
1333
1334 if (fgBeamTypeFlag == 1){ // we are in PbPb collisions --> apply the cut on TOF signal Dz
1335 Float_t radiusTOF = TMath::Sqrt(dxTOF*dxTOF + dzTOF*dzTOF);
1336 AliDebug(3,Form("TOF check (with fCutTOFdistance = %f) --> dx = %f, dz = %f, radius = %f", fCutTOFdistance, dxTOF, dzTOF, radiusTOF));
1337 if (radiusTOF > fCutTOFdistance){
1338 AliDebug(2, Form("************* the radius is outside the range! %f > %f, the track will be skipped", radiusTOF, fCutTOFdistance));
1339 cuts[41] = kTRUE;
1340 }
1341 }
1342 }
1343
1344 Bool_t cut=kFALSE;
1345 for (Int_t i=0; i<kNCuts; i++)
1346 if (cuts[i]) {cut = kTRUE;}
1347
1348 // for performance evaluate the CPU intensive cuts only when the others have passed, and when they are requested
1349 Double_t chi2TPCConstrainedVsGlobal = -2;
1350 Float_t nSigmaToVertex = -2;
1351 if (!cut)
1352 {
1353 // getting the track to vertex parameters
1354 if (fCutSigmaToVertexRequired)
1355 {
1356 nSigmaToVertex = GetSigmaToVertex(esdTrack);
1357 if (nSigmaToVertex > fCutNsigmaToVertex && fCutSigmaToVertexRequired)
1358 {
1359 cuts[12] = kTRUE;
1360 cut = kTRUE;
1361 }
1362 // if n sigma could not be calculated
1363 if (nSigmaToVertex<0 && fCutSigmaToVertexRequired)
1364 {
1365 cuts[13] = kTRUE;
1366 cut = kTRUE;
1367 }
1368 }
1369
1370 // max chi2 TPC constrained vs global track only if track passed the other cut
1371 if (fCutMaxChi2TPCConstrainedVsGlobal < 1e9)
1372 {
1373 const AliESDEvent* esdEvent = esdTrack->GetESDEvent();
1374
1375 if (!esdEvent)
1376 AliFatal("fCutMaxChi2TPCConstrainedVsGlobal set but ESD event not set in AliESDTrack. Use AliESDTrack::SetESDEvent before calling AliESDtrackCuts.");
1377
1378 // get vertex
1379 const AliESDVertex* vertex = 0;
1380 if (fCutMaxChi2TPCConstrainedVsGlobalVertexType & kVertexTracks)
1381 vertex = esdEvent->GetPrimaryVertexTracks();
1382
1383 if ((!vertex || !vertex->GetStatus()) && fCutMaxChi2TPCConstrainedVsGlobalVertexType & kVertexSPD)
1384 vertex = esdEvent->GetPrimaryVertexSPD();
1385
1386 if ((!vertex || !vertex->GetStatus()) && fCutMaxChi2TPCConstrainedVsGlobalVertexType & kVertexTPC)
1387 vertex = esdEvent->GetPrimaryVertexTPC();
1388
1389 if (vertex->GetStatus())
1390 chi2TPCConstrainedVsGlobal = esdTrack->GetChi2TPCConstrainedVsGlobal(vertex);
1391
1392 if (chi2TPCConstrainedVsGlobal < 0 || chi2TPCConstrainedVsGlobal > fCutMaxChi2TPCConstrainedVsGlobal)
1393 {
1394 cuts[39] = kTRUE;
1395 cut = kTRUE;
1396 }
1397 }
1398
1399 // max length in active volume
1400 Float_t lengthInActiveZoneTPC = -1;
1401 if (fCutMinLengthActiveVolumeTPC > 1.) { // do the calculation only if needed to save cpu-time
1402 if (esdTrack->GetESDEvent()) {
1403 if (esdTrack->GetInnerParam()) lengthInActiveZoneTPC = esdTrack->GetLengthInActiveZone(1, 1.8, 220, esdTrack->GetESDEvent()->GetMagneticField());
1404 //
1405 if (lengthInActiveZoneTPC < fCutMinLengthActiveVolumeTPC ) {
1406 cuts[42] = kTRUE;
1407 cut = kTRUE;
1408 }
1409 }
1410 }
1411
1412
1413 }
1414
1415 //########################################################################
1416 // filling histograms
1417 if (fHistogramsOn) {
1418 fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n tracks")));
1419 if (cut)
1420 fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n cut tracks")));
1421
1422 for (Int_t i=0; i<kNCuts; i++) {
1423 if (fhCutStatistics->GetXaxis()->FindBin(fgkCutNames[i]) < 1)
1424 AliFatal(Form("Inconsistency! Cut %d with name %s not found", i, fgkCutNames[i]));
1425
1426 if (cuts[i])
1427 fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin(fgkCutNames[i])));
1428
1429 for (Int_t j=i; j<kNCuts; j++) {
1430 if (cuts[i] && cuts[j]) {
1431 Float_t xC = fhCutCorrelation->GetXaxis()->GetBinCenter(fhCutCorrelation->GetXaxis()->FindBin(fgkCutNames[i]));
1432 Float_t yC = fhCutCorrelation->GetYaxis()->GetBinCenter(fhCutCorrelation->GetYaxis()->FindBin(fgkCutNames[j]));
1433 fhCutCorrelation->Fill(xC, yC);
1434 }
1435 }
1436 }
1437 }
1438
1439 // now we loop over the filling of the histograms twice: once "before" the cut, once "after"
1440 // the code is not in a function due to too many local variables that would need to be passed
1441
1442 for (Int_t id = 0; id < 2; id++)
1443 {
1444 // id = 0 --> before cut
1445 // id = 1 --> after cut
1446
1447 if (fHistogramsOn)
1448 {
1449 fhNClustersITS[id]->Fill(nClustersITS);
1450 fhNClustersTPC[id]->Fill(nClustersTPC);
1451 fhNSharedClustersTPC[id]->Fill(nClustersTPCShared);
1452 fhNCrossedRowsTPC[id]->Fill(nCrossedRowsTPC);
1453 fhRatioCrossedRowsOverFindableClustersTPC[id]->Fill(ratioCrossedRowsOverFindableClustersTPC);
1454 fhChi2PerClusterITS[id]->Fill(chi2PerClusterITS);
1455 fhChi2PerClusterTPC[id]->Fill(chi2PerClusterTPC);
1456 fhChi2TPCConstrainedVsGlobal[id]->Fill(chi2TPCConstrainedVsGlobal);
1457 fhNClustersForITSPID[id]->Fill(nITSPointsForPid);
1458 fhNMissingITSPoints[id]->Fill(nMissITSpts);
1459
1460 fhC11[id]->Fill(extCov[0]);
1461 fhC22[id]->Fill(extCov[2]);
1462 fhC33[id]->Fill(extCov[5]);
1463 fhC44[id]->Fill(extCov[9]);
1464 fhC55[id]->Fill(extCov[14]);
1465
1466 fhRel1PtUncertainty[id]->Fill(relUncertainty1Pt);
1467
1468 fhPt[id]->Fill(pt);
1469 fhEta[id]->Fill(eta);
1470 fhTOFdistance[id]->Fill(dxTOF, dzTOF);
1471
1472 Float_t bRes[2];
1473 bRes[0] = TMath::Sqrt(bCov[0]);
1474 bRes[1] = TMath::Sqrt(bCov[2]);
1475
1476 fhDZ[id]->Fill(b[1]);
1477 fhDXY[id]->Fill(b[0]);
1478 fhDXYDZ[id]->Fill(dcaToVertex);
1479 fhDXYvsDZ[id]->Fill(b[1],b[0]);
1480
1481 if (bRes[0]!=0 && bRes[1]!=0) {
1482 fhDZNormalized[id]->Fill(b[1]/bRes[1]);
1483 fhDXYNormalized[id]->Fill(b[0]/bRes[0]);
1484 fhDXYvsDZNormalized[id]->Fill(b[1]/bRes[1], b[0]/bRes[0]);
1485 fhNSigmaToVertex[id]->Fill(nSigmaToVertex);
1486 }
1487 }
1488
1489 // cut the track
1490 if (cut)
1491 return kFALSE;
1492 }
1493
1494 return kTRUE;
1495}
1496
1497//____________________________________________________________________
1498Bool_t AliESDtrackCuts::CheckITSClusterRequirement(ITSClusterRequirement req, Bool_t clusterL1, Bool_t clusterL2)
1499{
1500 // checks if the cluster requirement is fullfilled (in this case: return kTRUE)
1501
1502 switch (req)
1503 {
1504 case kOff: return kTRUE;
1505 case kNone: return !clusterL1 && !clusterL2;
1506 case kAny: return clusterL1 || clusterL2;
1507 case kFirst: return clusterL1;
1508 case kOnlyFirst: return clusterL1 && !clusterL2;
1509 case kSecond: return clusterL2;
1510 case kOnlySecond: return clusterL2 && !clusterL1;
1511 case kBoth: return clusterL1 && clusterL2;
1512 }
1513
1514 return kFALSE;
1515}
1516
1517//____________________________________________________________________
1518AliESDtrack* AliESDtrackCuts::GetTPCOnlyTrack(const AliESDEvent* esd, Int_t iTrack)
1519{
1520 // Utility function to create a TPC only track from the given esd track
1521 //
1522 // IMPORTANT: The track has to be deleted by the user
1523 //
1524 // NB. most of the functionality to get a TPC only track from an ESD track is in AliESDtrack, where it should be
1525 // there are only missing propagations here that are needed for old data
1526 // this function will therefore become obsolete
1527 //
1528 // adapted from code provided by CKB
1529
1530 if (!esd->GetPrimaryVertexTPC())
1531 return 0; // No TPC vertex no TPC tracks
1532
1533 if(!esd->GetPrimaryVertexTPC()->GetStatus())
1534 return 0; // TPC Vertex is created by default in AliESDEvent, do not use in this case
1535
1536 AliESDtrack* track = esd->GetTrack(iTrack);
1537 if (!track)
1538 return 0;
1539
1540 AliESDtrack *tpcTrack = new AliESDtrack();
1541
1542 // only true if we have a tpc track
1543 if (!track->FillTPCOnlyTrack(*tpcTrack))
1544 {
1545 delete tpcTrack;
1546 return 0;
1547 }
1548
1549 return tpcTrack;
1550}
1551
1552//____________________________________________________________________
1553TObjArray* AliESDtrackCuts::GetAcceptedTracks(const AliESDEvent* esd, Bool_t bTPC)
1554{
1555 //
1556 // returns an array of all tracks that pass the cuts
1557 // or an array of TPC only tracks (propagated to the TPC vertex during reco)
1558 // tracks that pass the cut
1559 //
1560 // NOTE: List has to be deleted by the user
1561
1562 TObjArray* acceptedTracks = new TObjArray();
1563
1564 // loop over esd tracks
1565 for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
1566 if(bTPC){
1567 if(!esd->GetPrimaryVertexTPC())return acceptedTracks; // No TPC vertex no TPC tracks
1568 if(!esd->GetPrimaryVertexTPC()->GetStatus())return acceptedTracks; // No proper TPC vertex, only the default
1569
1570 AliESDtrack *tpcTrack = GetTPCOnlyTrack(esd, iTrack);
1571 if (!tpcTrack)
1572 continue;
1573
1574 if (AcceptTrack(tpcTrack)) {
1575 acceptedTracks->Add(tpcTrack);
1576 }
1577 else
1578 delete tpcTrack;
1579 }
1580 else
1581 {
1582 AliESDtrack* track = esd->GetTrack(iTrack);
1583 if(AcceptTrack(track))
1584 acceptedTracks->Add(track);
1585 }
1586 }
1587 if(bTPC)acceptedTracks->SetOwner(kTRUE);
1588 return acceptedTracks;
1589}
1590
1591//____________________________________________________________________
1592Int_t AliESDtrackCuts::CountAcceptedTracks(const AliESDEvent* const esd)
1593{
1594 //
1595 // returns an the number of tracks that pass the cuts
1596 //
1597
1598 Int_t count = 0;
1599
1600 // loop over esd tracks
1601 for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
1602 AliESDtrack* track = esd->GetTrack(iTrack);
1603 if (AcceptTrack(track))
1604 count++;
1605 }
1606
1607 return count;
1608}
1609
1610//____________________________________________________________________
1611 void AliESDtrackCuts::DefineHistograms(Int_t color) {
1612 //
1613 // diagnostics histograms are defined
1614 //
1615
1616 fHistogramsOn=kTRUE;
1617
1618 Bool_t oldStatus = TH1::AddDirectoryStatus();
1619 TH1::AddDirectory(kFALSE);
1620
1621 //###################################################################################
1622 // defining histograms
1623
1624 fhCutStatistics = new TH1F("cut_statistics","cut statistics",kNCuts+4,-0.5,kNCuts+3.5);
1625
1626 fhCutStatistics->GetXaxis()->SetBinLabel(1,"n tracks");
1627 fhCutStatistics->GetXaxis()->SetBinLabel(2,"n cut tracks");
1628
1629 fhCutCorrelation = new TH2F("cut_correlation","cut correlation",kNCuts,-0.5,kNCuts-0.5,kNCuts,-0.5,kNCuts-0.5);;
1630
1631 for (Int_t i=0; i<kNCuts; i++) {
1632 fhCutStatistics->GetXaxis()->SetBinLabel(i+4,fgkCutNames[i]);
1633 fhCutCorrelation->GetXaxis()->SetBinLabel(i+1,fgkCutNames[i]);
1634 fhCutCorrelation->GetYaxis()->SetBinLabel(i+1,fgkCutNames[i]);
1635 }
1636
1637 fhCutStatistics ->SetLineColor(color);
1638 fhCutCorrelation ->SetLineColor(color);
1639 fhCutStatistics ->SetLineWidth(2);
1640 fhCutCorrelation ->SetLineWidth(2);
1641
1642 for (Int_t i=0; i<2; i++) {
1643 fhNClustersITS[i] = new TH1F("nClustersITS" ,"",8,-0.5,7.5);
1644 fhNClustersTPC[i] = new TH1F("nClustersTPC" ,"",165,-0.5,164.5);
1645 fhNSharedClustersTPC[i] = new TH1F("nSharedClustersTPC" ,"",165,-0.5,164.5);
1646 fhNCrossedRowsTPC[i] = new TH1F("nCrossedRowsTPC" ,"",165,-0.5,164.5);
1647 fhRatioCrossedRowsOverFindableClustersTPC[i] = new TH1F("ratioCrossedRowsOverFindableClustersTPC" ,"",60,0,1.5);
1648 fhChi2PerClusterITS[i] = new TH1F("chi2PerClusterITS","",500,0,10);
1649 fhChi2PerClusterTPC[i] = new TH1F("chi2PerClusterTPC","",500,0,10);
1650 fhChi2TPCConstrainedVsGlobal[i] = new TH1F("chi2TPCConstrainedVsGlobal","",600,-2,50);
1651 fhNClustersForITSPID[i] = new TH1F("nPointsForITSpid","",5,-0.5,4.5);
1652 fhNMissingITSPoints[i] = new TH1F("nMissingITSClusters","",7,-0.5,6.5);
1653
1654 fhC11[i] = new TH1F("covMatrixDiagonal11","",2000,0,20);
1655 fhC22[i] = new TH1F("covMatrixDiagonal22","",2000,0,20);
1656 fhC33[i] = new TH1F("covMatrixDiagonal33","",1000,0,0.1);
1657 fhC44[i] = new TH1F("covMatrixDiagonal44","",1000,0,0.1);
1658 fhC55[i] = new TH1F("covMatrixDiagonal55","",1000,0,5);
1659
1660 fhRel1PtUncertainty[i] = new TH1F("rel1PtUncertainty","",1000,0,5);
1661
1662 fhDXY[i] = new TH1F("dXY" ,"",500,-10,10);
1663 fhDZ[i] = new TH1F("dZ" ,"",500,-10,10);
1664 fhDXYDZ[i] = new TH1F("dXYDZ" ,"",500,0,10);
1665 fhDXYvsDZ[i] = new TH2F("dXYvsDZ","",200,-10,10,200,-10,10);
1666
1667 fhDXYNormalized[i] = new TH1F("dXYNormalized" ,"",500,-10,10);
1668 fhDZNormalized[i] = new TH1F("dZNormalized" ,"",500,-10,10);
1669 fhDXYvsDZNormalized[i] = new TH2F("dXYvsDZNormalized","",200,-10,10,200,-10,10);
1670
1671 fhNSigmaToVertex[i] = new TH1F("nSigmaToVertex","",500,0,10);
1672
1673 fhPt[i] = new TH1F("pt" ,"p_{T} distribution;p_{T} (GeV/c)", 800, 0.0, 10.0);
1674 fhEta[i] = new TH1F("eta" ,"#eta distribution;#eta",40,-2.0,2.0);
1675 fhTOFdistance[i] = new TH2F("TOFdistance" ,"TOF distance;dx (cm};dz (cm)", 150, -15, 15, 150, -15, 15);
1676
1677 fhNClustersITS[i]->SetTitle("n ITS clusters");
1678 fhNClustersTPC[i]->SetTitle("n TPC clusters");
1679 fhNSharedClustersTPC[i]->SetTitle("n TPC shared clusters");
1680 fhChi2PerClusterITS[i]->SetTitle("#Chi^{2} per ITS cluster");
1681 fhChi2PerClusterTPC[i]->SetTitle("#Chi^{2} per TPC cluster");
1682 fhChi2TPCConstrainedVsGlobal[i]->SetTitle("#Chi^{2} TPC constrained track vs global track");
1683 fhNClustersForITSPID[i]->SetTitle("n ITS points for PID");
1684 fhNMissingITSPoints[i]->SetTitle("n ITS layers with missing cluster");
1685
1686 fhC11[i]->SetTitle("cov 11 : #sigma_{y}^{2} [cm^{2}]");
1687 fhC22[i]->SetTitle("cov 22 : #sigma_{z}^{2} [cm^{2}]");
1688 fhC33[i]->SetTitle("cov 33 : #sigma_{sin(#phi)}^{2}");
1689 fhC44[i]->SetTitle("cov 44 : #sigma_{tan(#theta_{dip})}^{2}");
1690 fhC55[i]->SetTitle("cov 55 : #sigma_{1/p_{T}}^{2} [(c/GeV)^2]");
1691
1692 fhRel1PtUncertainty[i]->SetTitle("rel. uncertainty of 1/p_{T}");
1693
1694 fhDXY[i]->SetXTitle("transverse impact parameter (cm)");
1695 fhDZ[i]->SetXTitle("longitudinal impact parameter (cm)");
1696 fhDXYDZ[i]->SetTitle("absolute impact parameter;sqrt(dXY**2 + dZ**2) (cm)");
1697 fhDXYvsDZ[i]->SetXTitle("longitudinal impact parameter (cm)");
1698 fhDXYvsDZ[i]->SetYTitle("transverse impact parameter (cm)");
1699
1700 fhDXYNormalized[i]->SetTitle("normalized trans impact par (n#sigma)");
1701 fhDZNormalized[i]->SetTitle("normalized long impact par (n#sigma)");
1702 fhDXYvsDZNormalized[i]->SetTitle("normalized long impact par (n#sigma)");
1703 fhDXYvsDZNormalized[i]->SetYTitle("normalized trans impact par (n#sigma)");
1704 fhNSigmaToVertex[i]->SetTitle("n #sigma to vertex");
1705
1706 fhNClustersITS[i]->SetLineColor(color); fhNClustersITS[i]->SetLineWidth(2);
1707 fhNClustersTPC[i]->SetLineColor(color); fhNClustersTPC[i]->SetLineWidth(2);
1708 fhNSharedClustersTPC[i]->SetLineColor(color); fhNSharedClustersTPC[i]->SetLineWidth(2);
1709 fhChi2PerClusterITS[i]->SetLineColor(color); fhChi2PerClusterITS[i]->SetLineWidth(2);
1710 fhChi2PerClusterTPC[i]->SetLineColor(color); fhChi2PerClusterTPC[i]->SetLineWidth(2);
1711 fhChi2TPCConstrainedVsGlobal[i]->SetLineColor(color); fhChi2TPCConstrainedVsGlobal[i]->SetLineWidth(2);
1712 fhNClustersForITSPID[i]->SetLineColor(color); fhNClustersForITSPID[i]->SetLineWidth(2);
1713 fhNMissingITSPoints[i]->SetLineColor(color); fhNMissingITSPoints[i]->SetLineWidth(2);
1714
1715 fhC11[i]->SetLineColor(color); fhC11[i]->SetLineWidth(2);
1716 fhC22[i]->SetLineColor(color); fhC22[i]->SetLineWidth(2);
1717 fhC33[i]->SetLineColor(color); fhC33[i]->SetLineWidth(2);
1718 fhC44[i]->SetLineColor(color); fhC44[i]->SetLineWidth(2);
1719 fhC55[i]->SetLineColor(color); fhC55[i]->SetLineWidth(2);
1720
1721 fhRel1PtUncertainty[i]->SetLineColor(color); fhRel1PtUncertainty[i]->SetLineWidth(2);
1722
1723 fhDXY[i]->SetLineColor(color); fhDXY[i]->SetLineWidth(2);
1724 fhDZ[i]->SetLineColor(color); fhDZ[i]->SetLineWidth(2);
1725 fhDXYDZ[i]->SetLineColor(color); fhDXYDZ[i]->SetLineWidth(2);
1726
1727 fhDXYNormalized[i]->SetLineColor(color); fhDXYNormalized[i]->SetLineWidth(2);
1728 fhDZNormalized[i]->SetLineColor(color); fhDZNormalized[i]->SetLineWidth(2);
1729 fhNSigmaToVertex[i]->SetLineColor(color); fhNSigmaToVertex[i]->SetLineWidth(2);
1730 }
1731
1732 // The number of sigmas to the vertex is per definition gaussian
1733 ffDTheoretical = new TF1("nSigmaToVertexTheoretical","([0]/2.506628274)*exp(-(x**2)/2)",0,50);
1734 ffDTheoretical->SetParameter(0,1);
1735
1736 TH1::AddDirectory(oldStatus);
1737}
1738
1739//____________________________________________________________________
1740Bool_t AliESDtrackCuts::LoadHistograms(const Char_t* dir)
1741{
1742 //
1743 // loads the histograms from a file
1744 // if dir is empty a directory with the name of this object is taken (like in SaveHistogram)
1745 //
1746
1747 if (!dir)
1748 dir = GetName();
1749
1750 if (!gDirectory->cd(dir))
1751 return kFALSE;
1752
1753 ffDTheoretical = dynamic_cast<TF1*> (gDirectory->Get("nSigmaToVertexTheory"));
1754
1755 fhCutStatistics = dynamic_cast<TH1F*> (gDirectory->Get("cut_statistics"));
1756 fhCutCorrelation = dynamic_cast<TH2F*> (gDirectory->Get("cut_correlation"));
1757
1758 for (Int_t i=0; i<2; i++) {
1759 if (i==0)
1760 {
1761 gDirectory->cd("before_cuts");
1762 }
1763 else
1764 gDirectory->cd("after_cuts");
1765
1766 fhNClustersITS[i] = dynamic_cast<TH1F*> (gDirectory->Get("nClustersITS" ));
1767 fhNClustersTPC[i] = dynamic_cast<TH1F*> (gDirectory->Get("nClustersTPC" ));
1768 fhNSharedClustersTPC[i] = dynamic_cast<TH1F*> (gDirectory->Get("nSharedClustersTPC" ));
1769 fhNCrossedRowsTPC[i] = dynamic_cast<TH1F*> (gDirectory->Get("nCrossedRowsTPC" ));
1770 fhRatioCrossedRowsOverFindableClustersTPC[i] = dynamic_cast<TH1F*> (gDirectory->Get("ratioCrossedRowsOverFindableClustersTPC" ));
1771 fhChi2PerClusterITS[i] = dynamic_cast<TH1F*> (gDirectory->Get("chi2PerClusterITS"));
1772 fhChi2PerClusterTPC[i] = dynamic_cast<TH1F*> (gDirectory->Get("chi2PerClusterTPC"));
1773 fhChi2TPCConstrainedVsGlobal[i] = dynamic_cast<TH1F*> (gDirectory->Get("fhChi2TPCConstrainedVsGlobal"));
1774 fhNClustersForITSPID[i] = dynamic_cast<TH1F*> (gDirectory->Get("nPointsForITSpid"));
1775 fhNMissingITSPoints[i] = dynamic_cast<TH1F*> (gDirectory->Get("nMissingITSClusters"));
1776
1777 fhC11[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal11"));
1778 fhC22[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal22"));
1779 fhC33[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal33"));
1780 fhC44[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal44"));
1781 fhC55[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal55"));
1782
1783 fhRel1PtUncertainty[i] = dynamic_cast<TH1F*> (gDirectory->Get("rel1PtUncertainty"));
1784
1785 fhDXY[i] = dynamic_cast<TH1F*> (gDirectory->Get("dXY" ));
1786 fhDZ[i] = dynamic_cast<TH1F*> (gDirectory->Get("dZ" ));
1787 fhDXYDZ[i] = dynamic_cast<TH1F*> (gDirectory->Get("dXYDZ"));
1788 fhDXYvsDZ[i] = dynamic_cast<TH2F*> (gDirectory->Get("dXYvsDZ"));
1789
1790 fhDXYNormalized[i] = dynamic_cast<TH1F*> (gDirectory->Get("dXYNormalized" ));
1791 fhDZNormalized[i] = dynamic_cast<TH1F*> (gDirectory->Get("dZNormalized" ));
1792 fhDXYvsDZNormalized[i] = dynamic_cast<TH2F*> (gDirectory->Get("dXYvsDZNormalized"));
1793 fhNSigmaToVertex[i] = dynamic_cast<TH1F*> (gDirectory->Get("nSigmaToVertex"));
1794
1795 fhPt[i] = dynamic_cast<TH1F*> (gDirectory->Get("pt"));
1796 fhEta[i] = dynamic_cast<TH1F*> (gDirectory->Get("eta"));
1797 fhTOFdistance[i] = dynamic_cast<TH2F*> (gDirectory->Get("TOFdistance"));
1798
1799 gDirectory->cd("../");
1800 }
1801
1802 gDirectory->cd("..");
1803
1804 return kTRUE;
1805}
1806
1807//____________________________________________________________________
1808void AliESDtrackCuts::SaveHistograms(const Char_t* dir) {
1809 //
1810 // saves the histograms in a directory (dir)
1811 //
1812
1813 if (!fHistogramsOn) {
1814 AliDebug(0, "Histograms not on - cannot save histograms!!!");
1815 return;
1816 }
1817
1818 if (!dir)
1819 dir = GetName();
1820
1821 gDirectory->mkdir(dir);
1822 gDirectory->cd(dir);
1823
1824 gDirectory->mkdir("before_cuts");
1825 gDirectory->mkdir("after_cuts");
1826
1827 // a factor of 2 is needed since n sigma is positive
1828 ffDTheoretical->SetParameter(0,2*fhNSigmaToVertex[0]->Integral("width"));
1829 ffDTheoretical->Write("nSigmaToVertexTheory");
1830
1831 fhCutStatistics->Write();
1832 fhCutCorrelation->Write();
1833
1834 for (Int_t i=0; i<2; i++) {
1835 if (i==0)
1836 gDirectory->cd("before_cuts");
1837 else
1838 gDirectory->cd("after_cuts");
1839
1840 fhNClustersITS[i] ->Write();
1841 fhNClustersTPC[i] ->Write();
1842 fhNSharedClustersTPC[i] ->Write();
1843 fhNCrossedRowsTPC[i] ->Write();
1844 fhRatioCrossedRowsOverFindableClustersTPC[i] ->Write();
1845 fhChi2PerClusterITS[i] ->Write();
1846 fhChi2PerClusterTPC[i] ->Write();
1847 fhChi2TPCConstrainedVsGlobal[i] ->Write();
1848 fhNClustersForITSPID[i] ->Write();
1849 fhNMissingITSPoints[i] ->Write();
1850
1851 fhC11[i] ->Write();
1852 fhC22[i] ->Write();
1853 fhC33[i] ->Write();
1854 fhC44[i] ->Write();
1855 fhC55[i] ->Write();
1856
1857 fhRel1PtUncertainty[i] ->Write();
1858
1859 fhDXY[i] ->Write();
1860 fhDZ[i] ->Write();
1861 fhDXYDZ[i] ->Write();
1862 fhDXYvsDZ[i] ->Write();
1863
1864 fhDXYNormalized[i] ->Write();
1865 fhDZNormalized[i] ->Write();
1866 fhDXYvsDZNormalized[i] ->Write();
1867 fhNSigmaToVertex[i] ->Write();
1868
1869 fhPt[i] ->Write();
1870 fhEta[i] ->Write();
1871 fhTOFdistance[i] ->Write();
1872
1873 gDirectory->cd("../");
1874 }
1875
1876 gDirectory->cd("../");
1877}
1878
1879//____________________________________________________________________
1880void AliESDtrackCuts::DrawHistograms()
1881{
1882 // draws some histograms
1883
1884 TCanvas* canvas1 = new TCanvas(Form("%s_1", GetName()), "Track Quality Results1", 800, 800);
1885 canvas1->Divide(2, 2);
1886
1887 canvas1->cd(1);
1888 fhNClustersTPC[0]->SetStats(kFALSE);
1889 fhNClustersTPC[0]->Draw();
1890
1891 canvas1->cd(2);
1892 fhChi2PerClusterTPC[0]->SetStats(kFALSE);
1893 fhChi2PerClusterTPC[0]->Draw();
1894
1895 canvas1->cd(3);
1896 fhNSigmaToVertex[0]->SetStats(kFALSE);
1897 fhNSigmaToVertex[0]->GetXaxis()->SetRangeUser(0, 10);
1898 fhNSigmaToVertex[0]->Draw();
1899
1900 canvas1->SaveAs(Form("%s_%s.gif", GetName(), canvas1->GetName()));
1901
1902 TCanvas* canvas2 = new TCanvas(Form("%s_2", GetName()), "Track Quality Results2", 1200, 800);
1903 canvas2->Divide(3, 2);
1904
1905 canvas2->cd(1);
1906 fhC11[0]->SetStats(kFALSE);
1907 gPad->SetLogy();
1908 fhC11[0]->Draw();
1909
1910 canvas2->cd(2);
1911 fhC22[0]->SetStats(kFALSE);
1912 gPad->SetLogy();
1913 fhC22[0]->Draw();
1914
1915 canvas2->cd(3);
1916 fhC33[0]->SetStats(kFALSE);
1917 gPad->SetLogy();
1918 fhC33[0]->Draw();
1919
1920 canvas2->cd(4);
1921 fhC44[0]->SetStats(kFALSE);
1922 gPad->SetLogy();
1923 fhC44[0]->Draw();
1924
1925 canvas2->cd(5);
1926 fhC55[0]->SetStats(kFALSE);
1927 gPad->SetLogy();
1928 fhC55[0]->Draw();
1929
1930 canvas2->cd(6);
1931 fhRel1PtUncertainty[0]->SetStats(kFALSE);
1932 gPad->SetLogy();
1933 fhRel1PtUncertainty[0]->Draw();
1934
1935 canvas2->SaveAs(Form("%s_%s.gif", GetName(), canvas2->GetName()));
1936
1937 TCanvas* canvas3 = new TCanvas(Form("%s_3", GetName()), "Track Quality Results3", 1200, 800);
1938 canvas3->Divide(3, 2);
1939
1940 canvas3->cd(1);
1941 fhDXY[0]->SetStats(kFALSE);
1942 gPad->SetLogy();
1943 fhDXY[0]->Draw();
1944
1945 canvas3->cd(2);
1946 fhDZ[0]->SetStats(kFALSE);
1947 gPad->SetLogy();
1948 fhDZ[0]->Draw();
1949
1950 canvas3->cd(3);
1951 fhDXYvsDZ[0]->SetStats(kFALSE);
1952 gPad->SetLogz();
1953 gPad->SetRightMargin(0.15);
1954 fhDXYvsDZ[0]->Draw("COLZ");
1955
1956 canvas3->cd(4);
1957 fhDXYNormalized[0]->SetStats(kFALSE);
1958 gPad->SetLogy();
1959 fhDXYNormalized[0]->Draw();
1960
1961 canvas3->cd(5);
1962 fhDZNormalized[0]->SetStats(kFALSE);
1963 gPad->SetLogy();
1964 fhDZNormalized[0]->Draw();
1965
1966 canvas3->cd(6);
1967 fhDXYvsDZNormalized[0]->SetStats(kFALSE);
1968 gPad->SetLogz();
1969 gPad->SetRightMargin(0.15);
1970 fhDXYvsDZNormalized[0]->Draw("COLZ");
1971
1972 canvas3->SaveAs(Form("%s_%s.gif", GetName(), canvas3->GetName()));
1973
1974 TCanvas* canvas4 = new TCanvas(Form("%s_4", GetName()), "Track Quality Results4", 800, 500);
1975 canvas4->Divide(2, 1);
1976
1977 canvas4->cd(1);
1978 fhCutStatistics->SetStats(kFALSE);
1979 fhCutStatistics->LabelsOption("v");
1980 gPad->SetBottomMargin(0.3);
1981 fhCutStatistics->Draw();
1982
1983 canvas4->cd(2);
1984 fhCutCorrelation->SetStats(kFALSE);
1985 fhCutCorrelation->LabelsOption("v");
1986 gPad->SetBottomMargin(0.3);
1987 gPad->SetLeftMargin(0.3);
1988 fhCutCorrelation->Draw("COLZ");
1989
1990 canvas4->SaveAs(Form("%s_%s.gif", GetName(), canvas4->GetName()));
1991
1992 /*canvas->cd(1);
1993 fhDXYvsDZNormalized[0]->SetStats(kFALSE);
1994 fhDXYvsDZNormalized[0]->DrawCopy("COLZ");
1995
1996 canvas->cd(2);
1997 fhNClustersTPC[0]->SetStats(kFALSE);
1998 fhNClustersTPC[0]->DrawCopy();
1999
2000 canvas->cd(3);
2001 fhChi2PerClusterITS[0]->SetStats(kFALSE);
2002 fhChi2PerClusterITS[0]->DrawCopy();
2003 fhChi2PerClusterITS[1]->SetLineColor(2);
2004 fhChi2PerClusterITS[1]->DrawCopy("SAME");
2005
2006 canvas->cd(4);
2007 fhChi2PerClusterTPC[0]->SetStats(kFALSE);
2008 fhChi2PerClusterTPC[0]->DrawCopy();
2009 fhChi2PerClusterTPC[1]->SetLineColor(2);
2010 fhChi2PerClusterTPC[1]->DrawCopy("SAME");*/
2011}
2012//--------------------------------------------------------------------------
2013void AliESDtrackCuts::SetPtDepDCACuts(Double_t pt) {
2014 //
2015 // set the pt-dependent DCA cuts
2016 //
2017
2018 if(f1CutMaxDCAToVertexXYPtDep) {
2019 fCutMaxDCAToVertexXY=f1CutMaxDCAToVertexXYPtDep->Eval(pt);
2020 }
2021
2022 if(f1CutMaxDCAToVertexZPtDep) {
2023 fCutMaxDCAToVertexZ=f1CutMaxDCAToVertexZPtDep->Eval(pt);
2024 }
2025
2026 if(f1CutMinDCAToVertexXYPtDep) {
2027 fCutMinDCAToVertexXY=f1CutMinDCAToVertexXYPtDep->Eval(pt);
2028 }
2029
2030 if(f1CutMinDCAToVertexZPtDep) {
2031 fCutMinDCAToVertexZ=f1CutMinDCAToVertexZPtDep->Eval(pt);
2032 }
2033
2034
2035 return;
2036}
2037
2038
2039
2040//--------------------------------------------------------------------------
2041Bool_t AliESDtrackCuts::CheckPtDepDCA(TString dist,Bool_t print) const {
2042 //
2043 // Check the correctness of the string syntax
2044 //
2045 Bool_t retval=kTRUE;
2046
2047 if(!dist.Contains("pt")) {
2048 if(print) AliError("string must contain \"pt\"");
2049 retval= kFALSE;
2050 }
2051 return retval;
2052}
2053
2054 void AliESDtrackCuts::SetMaxDCAToVertexXYPtDep(const char *dist){
2055
2056 if(f1CutMaxDCAToVertexXYPtDep){
2057 delete f1CutMaxDCAToVertexXYPtDep;
2058 // resetiing both
2059 f1CutMaxDCAToVertexXYPtDep = 0;
2060 fCutMaxDCAToVertexXYPtDep = "";
2061 }
2062 if(!CheckPtDepDCA(dist,kTRUE)){
2063 return;
2064 }
2065 fCutMaxDCAToVertexXYPtDep = dist;
2066 TString tmp(dist);
2067 tmp.ReplaceAll("pt","x");
2068 f1CutMaxDCAToVertexXYPtDep = new TFormula("f1CutMaxDCAToVertexXYPtDep",tmp.Data());
2069
2070}
2071
2072 void AliESDtrackCuts::SetMaxDCAToVertexZPtDep(const char *dist){
2073
2074
2075 if(f1CutMaxDCAToVertexZPtDep){
2076 delete f1CutMaxDCAToVertexZPtDep;
2077 // resetiing both
2078 f1CutMaxDCAToVertexZPtDep = 0;
2079 fCutMaxDCAToVertexZPtDep = "";
2080 }
2081 if(!CheckPtDepDCA(dist,kTRUE))return;
2082
2083 fCutMaxDCAToVertexZPtDep = dist;
2084 TString tmp(dist);
2085 tmp.ReplaceAll("pt","x");
2086 f1CutMaxDCAToVertexZPtDep = new TFormula("f1CutMaxDCAToVertexZPtDep",tmp.Data());
2087
2088
2089}
2090
2091
2092 void AliESDtrackCuts::SetMinDCAToVertexXYPtDep(const char *dist){
2093
2094
2095 if(f1CutMinDCAToVertexXYPtDep){
2096 delete f1CutMinDCAToVertexXYPtDep;
2097 // resetiing both
2098 f1CutMinDCAToVertexXYPtDep = 0;
2099 fCutMinDCAToVertexXYPtDep = "";
2100 }
2101 if(!CheckPtDepDCA(dist,kTRUE))return;
2102
2103 fCutMinDCAToVertexXYPtDep = dist;
2104 TString tmp(dist);
2105 tmp.ReplaceAll("pt","x");
2106 f1CutMinDCAToVertexXYPtDep = new TFormula("f1CutMinDCAToVertexXYPtDep",tmp.Data());
2107
2108}
2109
2110
2111 void AliESDtrackCuts::SetMinDCAToVertexZPtDep(const char *dist){
2112
2113
2114
2115 if(f1CutMinDCAToVertexZPtDep){
2116 delete f1CutMinDCAToVertexZPtDep;
2117 // resetiing both
2118 f1CutMinDCAToVertexZPtDep = 0;
2119 fCutMinDCAToVertexZPtDep = "";
2120 }
2121 if(!CheckPtDepDCA(dist,kTRUE))return;
2122 fCutMinDCAToVertexZPtDep = dist;
2123 TString tmp(dist);
2124 tmp.ReplaceAll("pt","x");
2125 f1CutMinDCAToVertexZPtDep = new TFormula("f1CutMinDCAToVertexZPtDep",tmp.Data());
2126}
2127
2128AliESDtrackCuts* AliESDtrackCuts::GetMultEstTrackCuts(MultEstTrackCuts cut)
2129{
2130 // returns the multiplicity estimator track cuts objects to allow for user configuration
2131 // upon first call the objects are created
2132 //
2133 // the cut defined here correspond to GetStandardITSTPCTrackCuts2010 (apart from the one for without SPD)
2134
2135 if (!fgMultEstTrackCuts[kMultEstTrackCutGlobal])
2136 {
2137 // quality cut on ITS+TPC tracks
2138 fgMultEstTrackCuts[kMultEstTrackCutGlobal] = new AliESDtrackCuts();
2139 fgMultEstTrackCuts[kMultEstTrackCutGlobal]->SetMinNClustersTPC(70);
2140 fgMultEstTrackCuts[kMultEstTrackCutGlobal]->SetMaxChi2PerClusterTPC(4);
2141 fgMultEstTrackCuts[kMultEstTrackCutGlobal]->SetAcceptKinkDaughters(kFALSE);
2142 fgMultEstTrackCuts[kMultEstTrackCutGlobal]->SetRequireTPCRefit(kTRUE);
2143 fgMultEstTrackCuts[kMultEstTrackCutGlobal]->SetRequireITSRefit(kTRUE);
2144 //multiplicity underestimate if we use global tracks with |eta| > 0.9
2145 fgMultEstTrackCuts[kMultEstTrackCutGlobal]->SetEtaRange(-0.9, 0.9);
2146
2147 // quality cut on ITS_SA tracks (complementary to ITS+TPC)
2148 fgMultEstTrackCuts[kMultEstTrackCutITSSA] = new AliESDtrackCuts();
2149 fgMultEstTrackCuts[kMultEstTrackCutITSSA]->SetRequireITSRefit(kTRUE);
2150
2151 // primary selection for tracks with SPD hits
2152 fgMultEstTrackCuts[kMultEstTrackCutDCAwSPD] = new AliESDtrackCuts();
2153 fgMultEstTrackCuts[kMultEstTrackCutDCAwSPD]->SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kAny);
2154 fgMultEstTrackCuts[kMultEstTrackCutDCAwSPD]->SetMaxDCAToVertexXYPtDep("0.0182+0.0350/pt^1.01");
2155 fgMultEstTrackCuts[kMultEstTrackCutDCAwSPD]->SetMaxDCAToVertexZ(2);
2156
2157 // primary selection for tracks w/o SPD hits
2158 fgMultEstTrackCuts[kMultEstTrackCutDCAwoSPD] = new AliESDtrackCuts();
2159 fgMultEstTrackCuts[kMultEstTrackCutDCAwoSPD]->SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kNone);
2160 fgMultEstTrackCuts[kMultEstTrackCutDCAwoSPD]->SetMaxDCAToVertexXYPtDep("1.5*(0.0182+0.0350/pt^1.01)");
2161 fgMultEstTrackCuts[kMultEstTrackCutDCAwoSPD]->SetMaxDCAToVertexZ(2);
2162 }
2163
2164 return fgMultEstTrackCuts[cut];
2165}
2166
2167Int_t AliESDtrackCuts::GetReferenceMultiplicity(const AliESDEvent* esd, MultEstTrackType trackType, Float_t etaRange)
2168{
2169 // Get multiplicity estimate based on TPC/ITS tracks and tracklets
2170 // Adapted for AliESDtrackCuts from a version developed by: Ruben Shahoyan, Anton Alkin, Arvinder Palaha
2171 //
2172 // Returns a negative value if no reliable estimate can be provided:
2173 // -1 SPD vertex missing
2174 // -2 SPD VertexerZ dispersion too large
2175 // -3 Track vertex missing (not checked for kTracklets)
2176 // -4 Distance between SPD and track vertices too large (not checked for kTracklets)
2177 //
2178 // WARNING This functions does not cut on the z vtx. Depending on the eta range requested, you need to restrict your z vertex range!
2179 //
2180 // Strategy for combined estimators
2181 // 1. Count global tracks and flag them
2182 // 2. Count ITSSA as complementaries for ITSTPC+ or as main tracks
2183 // 3. Count the complementary SPD tracklets
2184
2185 const AliESDVertex* vertices[2];
2186 vertices[0] = esd->GetPrimaryVertexSPD();
2187 vertices[1] = esd->GetPrimaryVertexTracks();
2188
2189 if (!vertices[0]->GetStatus())
2190 {
2191 AliDebugClass(AliLog::kDebug, "No SPD vertex. Not able to make a reliable multiplicity estimate.");
2192 return -1;
2193 }
2194
2195 if (vertices[0]->IsFromVertexerZ() && vertices[0]->GetDispersion() > 0.02)
2196 {
2197 AliDebugClass(AliLog::kDebug, "Vertexer z dispersion > 0.02. Not able to make a reliable multiplicity estimate.");
2198 return -2;
2199 }
2200
2201 Int_t multiplicityEstimate = 0;
2202
2203 // SPD tracklet-only estimate
2204 if (trackType == kTracklets)
2205 {
2206 const AliMultiplicity* spdmult = esd->GetMultiplicity(); // spd multiplicity object
2207 for (Int_t i=0; i<spdmult->GetNumberOfTracklets(); ++i)
2208 {
2209 if (TMath::Abs(spdmult->GetEta(i)) > etaRange)
2210 continue; // eta selection for tracklets
2211 multiplicityEstimate++;
2212 }
2213 return multiplicityEstimate;
2214 }
2215
2216 if (!vertices[1]->GetStatus())
2217 {
2218 AliDebugClass(AliLog::kDebug, "No track vertex. Not able to make a reliable multiplicity estimate.");
2219 return -3;
2220 }
2221
2222 // TODO value of displacement to be studied
2223 const Float_t maxDisplacement = 0.5;
2224 //check for displaced vertices
2225 Double_t displacement = TMath::Abs(vertices[0]->GetZ() - vertices[1]->GetZ());
2226 if (displacement > maxDisplacement)
2227 {
2228 AliDebugClass(AliLog::kDebug, Form("Displaced vertices %f > %f",displacement,maxDisplacement));
2229 return -4;
2230 }
2231
2232 // update eta range in track cuts
2233 GetMultEstTrackCuts(kMultEstTrackCutITSSA)->SetEtaRange(-etaRange, etaRange);
2234 GetMultEstTrackCuts(kMultEstTrackCutDCAwSPD)->SetEtaRange(-etaRange, etaRange);
2235 GetMultEstTrackCuts(kMultEstTrackCutDCAwoSPD)->SetEtaRange(-etaRange, etaRange);
2236
2237 //*******************************************************************************************************
2238 //set counters to initial zeros
2239 Int_t tracksITSTPC = 0; //number of global tracks for a given event
2240 Int_t tracksITSSA = 0; //number of ITS standalone tracks for a given event
2241 Int_t tracksITSTPCSA_complementary = 0; //number of ITS standalone tracks complementary to TPC for a given event
2242 Int_t trackletsITSTPC_complementary = 0;//number of SPD tracklets complementary to global/ITSSA tracks for a given events
2243 Int_t trackletsITSSA_complementary = 0; //number of SPD tracklets complementary to ITSSA tracks for a given event
2244
2245 const Int_t nESDTracks = esd->GetNumberOfTracks();
2246
2247 // flags for secondary and rejected tracks
2248 const Int_t kRejBit = BIT(15); // set this bit in global tracks if it is rejected by a cut
2249 const Int_t kSecBit = BIT(16); // set this bit in global tracks if it is secondary according to a cut
2250
2251 for(Int_t itracks=0; itracks < nESDTracks; itracks++) {
2252 esd->GetTrack(itracks)->ResetBit(kSecBit|kRejBit); //reset bits used for flagging secondaries and rejected tracks in case they were changed before this analysis
2253 }
2254 const Int_t maxid = nESDTracks; // used to define bool array for check multiple associations of tracklets to one track. array starts at 0.
2255
2256 // bit mask for esd tracks, to check if multiple tracklets are associated to it
2257 TBits globalBits(maxid), pureITSBits(maxid);
2258 // why labels are used with the data? RS
2259 // Bool_t globalBits[maxid], pureITSBits[maxid];
2260 // for(Int_t i=0; i<maxid; i++){ // set all bools to false
2261 // globalBits[i]=kFALSE;
2262 // pureITSBits[i]=kFALSE;
2263 // }
2264
2265 //*******************************************************************************************************
2266 // get multiplicity from global tracks
2267 for(Int_t itracks = 0; itracks < nESDTracks; itracks++) { // flag the tracks
2268 AliESDtrack* track = esd->GetTrack(itracks);
2269
2270 // if track is a secondary from a V0, flag as a secondary
2271 if (track->IsOn(AliESDtrack::kMultInV0)) {
2272 track->SetBit(kSecBit);
2273 continue;
2274 }
2275 /* done via proper DCA cut
2276 //secondary?
2277 if (track->IsOn(AliESDtrack::kMultSec)) {
2278 track->SetBit(kSecBit);
2279 continue;
2280 }
2281 */
2282 // check tracks with ITS part
2283 //*******************************************************************************************************
2284 if (track->IsOn(AliESDtrack::kITSin) && !track->IsOn(AliESDtrack::kITSpureSA) && trackType == kTrackletsITSTPC) { // track has ITS part but is not an ITS_SA
2285 //*******************************************************************************************************
2286 // TPC+ITS
2287 if (track->IsOn(AliESDtrack::kTPCin)) { // Global track, has ITS and TPC contributions
2288 if (fgMultEstTrackCuts[kMultEstTrackCutGlobal]->AcceptTrack(track)) { // good ITSTPC track
2289 if (fgMultEstTrackCuts[kMultEstTrackCutDCAwSPD]->AcceptTrack(track) || fgMultEstTrackCuts[kMultEstTrackCutDCAwoSPD]->AcceptTrack(track)) {
2290 tracksITSTPC++; //global track counted
2291 globalBits.SetBitNumber(itracks);
2292 }
2293 else track->SetBit(kSecBit); // large DCA -> secondary, don't count either track not associated tracklet
2294 }
2295 else track->SetBit(kRejBit); // bad quality, don't count the track, but may count tracklet if associated
2296 }
2297 //*******************************************************************************************************
2298 // ITS complementary
2299 else if (fgMultEstTrackCuts[kMultEstTrackCutITSSA]->AcceptTrack(track)) { // good ITS complementary track
2300 if (fgMultEstTrackCuts[kMultEstTrackCutDCAwSPD]->AcceptTrack(track) || fgMultEstTrackCuts[kMultEstTrackCutDCAwoSPD]->AcceptTrack(track)) {
2301 tracksITSTPCSA_complementary++;
2302 globalBits.SetBitNumber(itracks);
2303 }
2304 else track->SetBit(kSecBit); // large DCA -> secondary, don't count either track not associated tracklet
2305 }
2306 else track->SetBit(kRejBit); // bad quality, don't count the track, but may count tracklet if associated
2307 }
2308 //*******************************************************************************************************
2309 // check tracks from ITS_SA_PURE
2310 if (track->IsOn(AliESDtrack::kITSin) && track->IsOn(AliESDtrack::kITSpureSA) && trackType == kTrackletsITSSA){
2311 if (fgMultEstTrackCuts[kMultEstTrackCutITSSA]->AcceptTrack(track)) { // good ITSSA track
2312 if (fgMultEstTrackCuts[kMultEstTrackCutDCAwSPD]->AcceptTrack(track) || fgMultEstTrackCuts[kMultEstTrackCutDCAwoSPD]->AcceptTrack(track)) {
2313 tracksITSSA++;
2314 pureITSBits.SetBitNumber(itracks);
2315 }
2316 else track->SetBit(kSecBit);
2317 }
2318 else track->SetBit(kRejBit);
2319 }
2320 }//ESD tracks counted
2321
2322 //*******************************************************************************************************
2323 // get multiplicity from ITS tracklets to complement TPC+ITS, and ITSpureSA
2324 const AliMultiplicity* spdmult = esd->GetMultiplicity(); // spd multiplicity object
2325 for (Int_t i=0; i<spdmult->GetNumberOfTracklets(); ++i) {
2326 if (TMath::Abs(spdmult->GetEta(i)) > etaRange) continue; // eta selection for tracklets
2327
2328 // if counting tracks+tracklets, check if clusters were already used in tracks
2329 Int_t id1, id2, id3, id4;
2330 spdmult->GetTrackletTrackIDs ( i, 0, id1, id2 ); // references for eventual Global/ITS_SA tracks
2331 spdmult->GetTrackletTrackIDs ( i, 1, id3, id4 ); // references for eventual ITS_SA_pure tracks
2332
2333 // are both clusters from the same tracks? If not, skip the tracklet (shouldn't change things much)
2334 if ( ( id1 != id2 && id1 >= 0 && id2 >= 0 ) || ( id3 != id4 && id3 >= 0 && id4 >= 0 ) ) continue;
2335
2336 //referenced track
2337 //at this point we either have id1 = id2 (id3 = id4) or only one of the ids pair is -1
2338 // id1>=0, id2>=0, id1=id2 : tracklet has associated track
2339 // id1>=0, id2 = -1 : 1st layer cluster has associated track
2340 // id1=-1, id2>=0 : 2nd layer cluster has associated track
2341 // id1=-1, id2=-1 : tracklet has no associated track
2342 //
2343 Int_t bUsedInGlobal(-1);
2344 if ( id1 != -1 ) bUsedInGlobal = globalBits.TestBitNumber(id1) ? id1 : -1;
2345 else if ( id2 != -1) bUsedInGlobal = globalBits.TestBitNumber(id2) ? id2 : -1;
2346 Int_t bUsedInPureITS(-1);
2347 if ( id3 != -1 ) bUsedInPureITS = pureITSBits.TestBitNumber(id3) ? id3 : -1;
2348 else if ( id4 != -1) bUsedInPureITS = pureITSBits.TestBitNumber(id4) ? id4 : -1;
2349 //
2350 AliESDtrack* tr_global = bUsedInGlobal >= 0 ? esd->GetTrack ( bUsedInGlobal ) : 0;
2351 AliESDtrack* tr_itssa = bUsedInPureITS >= 0 ? esd->GetTrack ( bUsedInPureITS ) : 0;
2352 //
2353 // has associated pure ITS track been associated to a previous tracklet?
2354 //*******************************************************************************************************
2355 if (trackType == kTrackletsITSTPC) {
2356 //*******************************************************************************************************
2357 // count tracklets towards global+complimentary tracks
2358 if ( ( tr_global && !tr_global->TestBit ( kSecBit ) ) && ( tr_global && tr_global->TestBit ( kRejBit ) ) ) { // count tracklet as bad quality track
2359 globalBits.SetBitNumber( bUsedInGlobal ); // mark global track linked to this tracklet as used
2360 ++trackletsITSTPC_complementary;
2361 }
2362
2363 if ( bUsedInGlobal < 0 ) ++trackletsITSTPC_complementary; //no associated track, just count the tracklet
2364 } else {
2365 //*******************************************************************************************************
2366 // count tracklets towards ITS_SA_pure tracks
2367 if ( ( tr_itssa && !tr_itssa->TestBit ( kSecBit ) ) && ( tr_itssa && tr_itssa->TestBit ( kRejBit ) ) ) { // count tracklet as bad quality track
2368 pureITSBits.SetBitNumber( bUsedInPureITS ); // mark ITS pure SA track linked to this tracklet as used
2369 ++trackletsITSSA_complementary;
2370 }
2371
2372 if ( bUsedInPureITS < 0 ) ++trackletsITSSA_complementary; //no associated track, just count the tracklet
2373 }
2374 }
2375
2376 //*******************************************************************************************************
2377 if (trackType == kTrackletsITSTPC)
2378 multiplicityEstimate = tracksITSTPC + tracksITSTPCSA_complementary + trackletsITSTPC_complementary;
2379 else
2380 multiplicityEstimate = tracksITSSA + trackletsITSSA_complementary;
2381
2382 return multiplicityEstimate;
2383}
2384
2385//____________________________________________________________________
2386void AliESDtrackCuts::SetRequireStandardTOFmatchCuts(){
2387
2388 // setting the TOF cuts flags (kTOFout = TOF matching distance) to true, to include the selection on the standard TOF matching
2389
2390 SetRequireTOFout(kTRUE);
2391 SetFlagCutTOFdistance(kTRUE);
2392 SetCutTOFdistance(3.);
2393
2394}
2395