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