1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 /* $Id: AliESDtrackCuts.cxx 24534 2008-03-16 22:22:11Z fca $ */
18 #include "AliESDtrackCuts.h"
20 #include <AliESDtrack.h>
21 #include <AliESDVertex.h>
22 #include <AliESDEvent.h>
27 #include <TDirectory.h>
31 //____________________________________________________________________
32 ClassImp(AliESDtrackCuts)
35 const Char_t* AliESDtrackCuts::fgkCutNames[kNCuts] = {
37 "require TPC standalone",
41 "#Chi^{2}/cluster TPC",
42 "#Chi^{2}/cluster ITS",
58 "trk-to-vtx max dca 2D absolute",
59 "trk-to-vtx max dca xy absolute",
60 "trk-to-vtx max dca z absolute",
61 "trk-to-vtx min dca 2D absolute",
62 "trk-to-vtx min dca xy absolute",
63 "trk-to-vtx min dca z absolute",
64 "SPD cluster requirement",
65 "SDD cluster requirement",
66 "SSD cluster requirement",
67 "require ITS stand-alone",
68 "rel 1/pt uncertainty",
72 //____________________________________________________________________
73 AliESDtrackCuts::AliESDtrackCuts(const Char_t* name, const Char_t* title) : AliAnalysisCuts(name,title),
74 fCutMinNClusterTPC(0),
75 fCutMinNClusterITS(0),
76 fCutMaxChi2PerClusterTPC(0),
77 fCutMaxChi2PerClusterITS(0),
83 fCutMaxRel1PtUncertainty(0),
84 fCutAcceptKinkDaughters(0),
85 fCutAcceptSharedTPCClusters(0),
86 fCutMaxFractionSharedTPCClusters(0),
87 fCutRequireTPCRefit(0),
88 fCutRequireTPCStandAlone(0),
89 fCutRequireITSRefit(0),
91 fCutRequireITSStandAlone(0),
92 fCutRequireITSpureSA(0),
93 fCutNsigmaToVertex(0),
94 fCutSigmaToVertexRequired(0),
95 fCutMaxDCAToVertexXY(0),
96 fCutMaxDCAToVertexZ(0),
97 fCutMinDCAToVertexXY(0),
98 fCutMinDCAToVertexZ(0),
99 fCutMaxDCAToVertexXYPtDep(""),
100 fCutMaxDCAToVertexZPtDep(""),
101 fCutMinDCAToVertexXYPtDep(""),
102 fCutMinDCAToVertexZPtDep(""),
103 f1CutMaxDCAToVertexXYPtDep(0x0),
104 f1CutMaxDCAToVertexZPtDep(0x0),
105 f1CutMinDCAToVertexXYPtDep(0x0),
106 f1CutMinDCAToVertexZPtDep(0x0),
107 fCutDCAToVertex2D(0),
133 //##############################################################################
134 // setting default cuts
135 SetMinNClustersTPC();
136 SetMinNClustersITS();
137 SetMaxChi2PerClusterTPC();
138 SetMaxChi2PerClusterITS();
139 SetMaxCovDiagonalElements();
140 SetMaxRel1PtUncertainty();
141 SetRequireTPCRefit();
142 SetRequireTPCStandAlone();
143 SetRequireITSRefit();
144 SetRequireITSPid(kFALSE);
145 SetRequireITSStandAlone(kFALSE);
146 SetRequireITSPureStandAlone(kFALSE);
147 SetAcceptKinkDaughters();
148 SetMaxNsigmaToVertex();
149 SetMaxDCAToVertexXY();
150 SetMaxDCAToVertexZ();
152 SetMinDCAToVertexXY();
153 SetMinDCAToVertexZ();
161 SetClusterRequirementITS(kSPD);
162 SetClusterRequirementITS(kSDD);
163 SetClusterRequirementITS(kSSD);
168 //_____________________________________________________________________________
169 AliESDtrackCuts::AliESDtrackCuts(const AliESDtrackCuts &c) : AliAnalysisCuts(c),
170 fCutMinNClusterTPC(0),
171 fCutMinNClusterITS(0),
172 fCutMaxChi2PerClusterTPC(0),
173 fCutMaxChi2PerClusterITS(0),
179 fCutMaxRel1PtUncertainty(0),
180 fCutAcceptKinkDaughters(0),
181 fCutAcceptSharedTPCClusters(0),
182 fCutMaxFractionSharedTPCClusters(0),
183 fCutRequireTPCRefit(0),
184 fCutRequireTPCStandAlone(0),
185 fCutRequireITSRefit(0),
186 fCutRequireITSPid(0),
187 fCutRequireITSStandAlone(0),
188 fCutRequireITSpureSA(0),
189 fCutNsigmaToVertex(0),
190 fCutSigmaToVertexRequired(0),
191 fCutMaxDCAToVertexXY(0),
192 fCutMaxDCAToVertexZ(0),
193 fCutMinDCAToVertexXY(0),
194 fCutMinDCAToVertexZ(0),
195 fCutMaxDCAToVertexXYPtDep(""),
196 fCutMaxDCAToVertexZPtDep(""),
197 fCutMinDCAToVertexXYPtDep(""),
198 fCutMinDCAToVertexZPtDep(""),
199 f1CutMaxDCAToVertexXYPtDep(0x0),
200 f1CutMaxDCAToVertexZPtDep(0x0),
201 f1CutMinDCAToVertexXYPtDep(0x0),
202 f1CutMinDCAToVertexZPtDep(0x0),
203 fCutDCAToVertex2D(0),
227 ((AliESDtrackCuts &) c).Copy(*this);
230 AliESDtrackCuts::~AliESDtrackCuts()
236 for (Int_t i=0; i<2; i++) {
238 if (fhNClustersITS[i])
239 delete fhNClustersITS[i];
240 if (fhNClustersTPC[i])
241 delete fhNClustersTPC[i];
242 if (fhChi2PerClusterITS[i])
243 delete fhChi2PerClusterITS[i];
244 if (fhChi2PerClusterTPC[i])
245 delete fhChi2PerClusterTPC[i];
257 if (fhRel1PtUncertainty[i])
258 delete fhRel1PtUncertainty[i];
269 if (fhDXYNormalized[i])
270 delete fhDXYNormalized[i];
271 if (fhDZNormalized[i])
272 delete fhDZNormalized[i];
273 if (fhDXYvsDZNormalized[i])
274 delete fhDXYvsDZNormalized[i];
275 if (fhNSigmaToVertex[i])
276 delete fhNSigmaToVertex[i];
283 if(f1CutMaxDCAToVertexXYPtDep)delete f1CutMaxDCAToVertexXYPtDep;
284 f1CutMaxDCAToVertexXYPtDep = 0;
285 if( f1CutMaxDCAToVertexZPtDep) delete f1CutMaxDCAToVertexZPtDep;
286 f1CutMaxDCAToVertexZPtDep = 0;
287 if( f1CutMinDCAToVertexXYPtDep)delete f1CutMinDCAToVertexXYPtDep;
288 f1CutMinDCAToVertexXYPtDep = 0;
289 if(f1CutMinDCAToVertexZPtDep)delete f1CutMinDCAToVertexZPtDep;
290 f1CutMinDCAToVertexZPtDep = 0;
294 delete ffDTheoretical;
297 delete fhCutStatistics;
298 if (fhCutCorrelation)
299 delete fhCutCorrelation;
302 void AliESDtrackCuts::Init()
305 // sets everything to zero
308 fCutMinNClusterTPC = 0;
309 fCutMinNClusterITS = 0;
311 fCutMaxChi2PerClusterTPC = 0;
312 fCutMaxChi2PerClusterITS = 0;
314 for (Int_t i = 0; i < 3; i++)
315 fCutClusterRequirementITS[i] = kOff;
323 fCutMaxRel1PtUncertainty = 0;
325 fCutAcceptKinkDaughters = 0;
326 fCutAcceptSharedTPCClusters = 0;
327 fCutMaxFractionSharedTPCClusters = 0;
328 fCutRequireTPCRefit = 0;
329 fCutRequireTPCStandAlone = 0;
330 fCutRequireITSRefit = 0;
331 fCutRequireITSPid = 0;
332 fCutRequireITSStandAlone = 0;
333 fCutRequireITSpureSA = 0;
335 fCutNsigmaToVertex = 0;
336 fCutSigmaToVertexRequired = 0;
337 fCutMaxDCAToVertexXY = 0;
338 fCutMaxDCAToVertexZ = 0;
339 fCutDCAToVertex2D = 0;
340 fCutMinDCAToVertexXY = 0;
341 fCutMinDCAToVertexZ = 0;
342 fCutMaxDCAToVertexXYPtDep = "";
343 fCutMaxDCAToVertexZPtDep = "";
344 fCutMinDCAToVertexXYPtDep = "";
345 fCutMinDCAToVertexZPtDep = "";
347 if(f1CutMaxDCAToVertexXYPtDep)delete f1CutMaxDCAToVertexXYPtDep;
348 f1CutMaxDCAToVertexXYPtDep = 0;
349 if( f1CutMaxDCAToVertexXYPtDep) delete f1CutMaxDCAToVertexXYPtDep;
350 f1CutMaxDCAToVertexXYPtDep = 0;
351 if( f1CutMaxDCAToVertexZPtDep) delete f1CutMaxDCAToVertexZPtDep;
352 f1CutMaxDCAToVertexZPtDep = 0;
353 if( f1CutMinDCAToVertexXYPtDep)delete f1CutMinDCAToVertexXYPtDep;
354 f1CutMinDCAToVertexXYPtDep = 0;
355 if(f1CutMinDCAToVertexZPtDep)delete f1CutMinDCAToVertexZPtDep;
356 f1CutMinDCAToVertexZPtDep = 0;
374 fHistogramsOn = kFALSE;
376 for (Int_t i=0; i<2; ++i)
378 fhNClustersITS[i] = 0;
379 fhNClustersTPC[i] = 0;
381 fhChi2PerClusterITS[i] = 0;
382 fhChi2PerClusterTPC[i] = 0;
390 fhRel1PtUncertainty[i] = 0;
397 fhDXYNormalized[i] = 0;
398 fhDZNormalized[i] = 0;
399 fhDXYvsDZNormalized[i] = 0;
400 fhNSigmaToVertex[i] = 0;
408 fhCutCorrelation = 0;
411 //_____________________________________________________________________________
412 AliESDtrackCuts &AliESDtrackCuts::operator=(const AliESDtrackCuts &c)
415 // Assignment operator
418 if (this != &c) ((AliESDtrackCuts &) c).Copy(*this);
422 //_____________________________________________________________________________
423 void AliESDtrackCuts::Copy(TObject &c) const
429 AliESDtrackCuts& target = (AliESDtrackCuts &) c;
433 target.fCutMinNClusterTPC = fCutMinNClusterTPC;
434 target.fCutMinNClusterITS = fCutMinNClusterITS;
436 target.fCutMaxChi2PerClusterTPC = fCutMaxChi2PerClusterTPC;
437 target.fCutMaxChi2PerClusterITS = fCutMaxChi2PerClusterITS;
439 for (Int_t i = 0; i < 3; i++)
440 target.fCutClusterRequirementITS[i] = fCutClusterRequirementITS[i];
442 target.fCutMaxC11 = fCutMaxC11;
443 target.fCutMaxC22 = fCutMaxC22;
444 target.fCutMaxC33 = fCutMaxC33;
445 target.fCutMaxC44 = fCutMaxC44;
446 target.fCutMaxC55 = fCutMaxC55;
448 target.fCutMaxRel1PtUncertainty = fCutMaxRel1PtUncertainty;
450 target.fCutAcceptKinkDaughters = fCutAcceptKinkDaughters;
451 target.fCutAcceptSharedTPCClusters = fCutAcceptSharedTPCClusters;
452 target.fCutMaxFractionSharedTPCClusters = fCutMaxFractionSharedTPCClusters;
453 target.fCutRequireTPCRefit = fCutRequireTPCRefit;
454 target.fCutRequireTPCStandAlone = fCutRequireTPCStandAlone;
455 target.fCutRequireITSRefit = fCutRequireITSRefit;
456 target.fCutRequireITSPid = fCutRequireITSPid;
457 target.fCutRequireITSStandAlone = fCutRequireITSStandAlone;
458 target.fCutRequireITSpureSA = fCutRequireITSpureSA;
460 target.fCutNsigmaToVertex = fCutNsigmaToVertex;
461 target.fCutSigmaToVertexRequired = fCutSigmaToVertexRequired;
462 target.fCutMaxDCAToVertexXY = fCutMaxDCAToVertexXY;
463 target.fCutMaxDCAToVertexZ = fCutMaxDCAToVertexZ;
464 target.fCutDCAToVertex2D = fCutDCAToVertex2D;
465 target.fCutMinDCAToVertexXY = fCutMinDCAToVertexXY;
466 target.fCutMinDCAToVertexZ = fCutMinDCAToVertexZ;
468 target.fCutMaxDCAToVertexXYPtDep = fCutMaxDCAToVertexXYPtDep;
469 target.SetMaxDCAToVertexXYPtDep(fCutMaxDCAToVertexXYPtDep.Data());
471 target.fCutMaxDCAToVertexZPtDep = fCutMaxDCAToVertexZPtDep;
472 target.SetMaxDCAToVertexZPtDep(fCutMaxDCAToVertexZPtDep.Data());
474 target.fCutMinDCAToVertexXYPtDep = fCutMinDCAToVertexXYPtDep;
475 target.SetMinDCAToVertexXYPtDep(fCutMinDCAToVertexXYPtDep.Data());
477 target.fCutMinDCAToVertexZPtDep = fCutMinDCAToVertexZPtDep;
478 target.SetMinDCAToVertexZPtDep(fCutMinDCAToVertexZPtDep.Data());
480 target.fPMin = fPMin;
481 target.fPMax = fPMax;
482 target.fPtMin = fPtMin;
483 target.fPtMax = fPtMax;
484 target.fPxMin = fPxMin;
485 target.fPxMax = fPxMax;
486 target.fPyMin = fPyMin;
487 target.fPyMax = fPyMax;
488 target.fPzMin = fPzMin;
489 target.fPzMax = fPzMax;
490 target.fEtaMin = fEtaMin;
491 target.fEtaMax = fEtaMax;
492 target.fRapMin = fRapMin;
493 target.fRapMax = fRapMax;
495 target.fHistogramsOn = fHistogramsOn;
497 for (Int_t i=0; i<2; ++i)
499 if (fhNClustersITS[i]) target.fhNClustersITS[i] = (TH1F*) fhNClustersITS[i]->Clone();
500 if (fhNClustersTPC[i]) target.fhNClustersTPC[i] = (TH1F*) fhNClustersTPC[i]->Clone();
502 if (fhChi2PerClusterITS[i]) target.fhChi2PerClusterITS[i] = (TH1F*) fhChi2PerClusterITS[i]->Clone();
503 if (fhChi2PerClusterTPC[i]) target.fhChi2PerClusterTPC[i] = (TH1F*) fhChi2PerClusterTPC[i]->Clone();
505 if (fhC11[i]) target.fhC11[i] = (TH1F*) fhC11[i]->Clone();
506 if (fhC22[i]) target.fhC22[i] = (TH1F*) fhC22[i]->Clone();
507 if (fhC33[i]) target.fhC33[i] = (TH1F*) fhC33[i]->Clone();
508 if (fhC44[i]) target.fhC44[i] = (TH1F*) fhC44[i]->Clone();
509 if (fhC55[i]) target.fhC55[i] = (TH1F*) fhC55[i]->Clone();
511 if (fhRel1PtUncertainty[i]) target.fhRel1PtUncertainty[i] = (TH1F*) fhRel1PtUncertainty[i]->Clone();
513 if (fhDXY[i]) target.fhDXY[i] = (TH1F*) fhDXY[i]->Clone();
514 if (fhDZ[i]) target.fhDZ[i] = (TH1F*) fhDZ[i]->Clone();
515 if (fhDXYDZ[i]) target.fhDXYDZ[i] = (TH1F*) fhDXYDZ[i]->Clone();
516 if (fhDXYvsDZ[i]) target.fhDXYvsDZ[i] = (TH2F*) fhDXYvsDZ[i]->Clone();
518 if (fhDXYNormalized[i]) target.fhDXYNormalized[i] = (TH1F*) fhDXYNormalized[i]->Clone();
519 if (fhDZNormalized[i]) target.fhDZNormalized[i] = (TH1F*) fhDZNormalized[i]->Clone();
520 if (fhDXYvsDZNormalized[i]) target.fhDXYvsDZNormalized[i] = (TH2F*) fhDXYvsDZNormalized[i]->Clone();
521 if (fhNSigmaToVertex[i]) target.fhNSigmaToVertex[i] = (TH1F*) fhNSigmaToVertex[i]->Clone();
523 if (fhPt[i]) target.fhPt[i] = (TH1F*) fhPt[i]->Clone();
524 if (fhEta[i]) target.fhEta[i] = (TH1F*) fhEta[i]->Clone();
526 if (ffDTheoretical) target.ffDTheoretical = (TF1*) ffDTheoretical->Clone();
528 if (fhCutStatistics) target.fhCutStatistics = (TH1F*) fhCutStatistics->Clone();
529 if (fhCutCorrelation) target.fhCutCorrelation = (TH2F*) fhCutCorrelation->Clone();
534 //_____________________________________________________________________________
535 Long64_t AliESDtrackCuts::Merge(TCollection* list) {
536 // Merge a list of AliESDtrackCuts objects with this (needed for PROOF)
537 // Returns the number of merged objects (including this)
544 TIterator* iter = list->MakeIterator();
547 // collection of measured and generated histograms
549 while ((obj = iter->Next())) {
551 AliESDtrackCuts* entry = dynamic_cast<AliESDtrackCuts*>(obj);
555 if (!entry->fHistogramsOn)
558 for (Int_t i=0; i<2; i++) {
560 fhNClustersITS[i] ->Add(entry->fhNClustersITS[i] );
561 fhNClustersTPC[i] ->Add(entry->fhNClustersTPC[i] );
563 fhChi2PerClusterITS[i] ->Add(entry->fhChi2PerClusterITS[i]);
564 fhChi2PerClusterTPC[i] ->Add(entry->fhChi2PerClusterTPC[i]);
566 fhC11[i] ->Add(entry->fhC11[i] );
567 fhC22[i] ->Add(entry->fhC22[i] );
568 fhC33[i] ->Add(entry->fhC33[i] );
569 fhC44[i] ->Add(entry->fhC44[i] );
570 fhC55[i] ->Add(entry->fhC55[i] );
572 fhRel1PtUncertainty[i] ->Add(entry->fhRel1PtUncertainty[i]);
574 fhDXY[i] ->Add(entry->fhDXY[i] );
575 fhDZ[i] ->Add(entry->fhDZ[i] );
576 fhDXYDZ[i] ->Add(entry->fhDXYDZ[i] );
577 fhDXYvsDZ[i] ->Add(entry->fhDXYvsDZ[i] );
579 fhDXYNormalized[i] ->Add(entry->fhDXYNormalized[i] );
580 fhDZNormalized[i] ->Add(entry->fhDZNormalized[i] );
581 fhDXYvsDZNormalized[i] ->Add(entry->fhDXYvsDZNormalized[i]);
582 fhNSigmaToVertex[i] ->Add(entry->fhNSigmaToVertex[i]);
584 fhPt[i] ->Add(entry->fhPt[i]);
585 fhEta[i] ->Add(entry->fhEta[i]);
588 fhCutStatistics ->Add(entry->fhCutStatistics);
589 fhCutCorrelation ->Add(entry->fhCutCorrelation);
596 //____________________________________________________________________
597 AliESDtrackCuts* AliESDtrackCuts::GetStandardTPCOnlyTrackCuts()
599 // creates an AliESDtrackCuts object and fills it with standard (pre data-taking) values for TPC-only cuts
601 Printf("AliESDtrackCuts::GetStandardTPCOnlyTrackCuts: Creating track cuts for TPC-only.");
603 AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
605 esdTrackCuts->SetMinNClustersTPC(50);
606 esdTrackCuts->SetMaxChi2PerClusterTPC(4);
607 esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
609 esdTrackCuts->SetMaxDCAToVertexZ(3.2);
610 esdTrackCuts->SetMaxDCAToVertexXY(2.4);
611 esdTrackCuts->SetDCAToVertex2D(kTRUE);
616 //____________________________________________________________________
617 AliESDtrackCuts* AliESDtrackCuts::GetStandardITSTPCTrackCuts2009(Bool_t selPrimaries)
619 // creates an AliESDtrackCuts object and fills it with standard values for ITS-TPC cuts for pp 2009 data
621 Printf("AliESDtrackCuts::GetStandardITSTPCTrackCuts: Creating track cuts for ITS+TPC.");
623 AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
626 esdTrackCuts->SetRequireTPCStandAlone(kTRUE); // to get chi2 and ncls of kTPCin
627 esdTrackCuts->SetMinNClustersTPC(70);
628 esdTrackCuts->SetMaxChi2PerClusterTPC(4);
629 esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
630 esdTrackCuts->SetRequireTPCRefit(kTRUE);
632 esdTrackCuts->SetRequireITSRefit(kTRUE);
633 esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
634 AliESDtrackCuts::kAny);
636 // 7*(0.0050+0.0060/pt^0.9)
637 esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0350+0.0420/pt^0.9");
639 esdTrackCuts->SetMaxDCAToVertexZ(1.e6);
640 esdTrackCuts->SetDCAToVertex2D(kFALSE);
641 esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
642 //esdTrackCuts->SetEtaRange(-0.8,+0.8);
646 //____________________________________________________________________
647 AliESDtrackCuts* AliESDtrackCuts::GetStandardITSPureSATrackCuts2009(Bool_t selPrimaries, Bool_t useForPid)
649 // creates an AliESDtrackCuts object and fills it with standard values for ITS pure SA tracks
651 AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
652 esdTrackCuts->SetRequireITSPureStandAlone(kTRUE);
653 esdTrackCuts->SetRequireITSRefit(kTRUE);
654 esdTrackCuts->SetMinNClustersITS(4);
655 esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
656 AliESDtrackCuts::kAny);
657 esdTrackCuts->SetMaxChi2PerClusterITS(1.);
660 // 7*(0.0085+0.0026/pt^1.55)
661 esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0595+0.0182/pt^1.55");
664 esdTrackCuts->SetRequireITSPid(kTRUE);
668 //____________________________________________________________________
669 AliESDtrackCuts* AliESDtrackCuts::GetStandardITSSATrackCuts2009(Bool_t selPrimaries, Bool_t useForPid)
671 // creates an AliESDtrackCuts object and fills it with standard values for ITS pure SA tracks
673 AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
674 esdTrackCuts->SetRequireITSStandAlone(kTRUE);
675 esdTrackCuts->SetRequireITSPureStandAlone(kFALSE);
676 esdTrackCuts->SetRequireITSRefit(kTRUE);
677 esdTrackCuts->SetMinNClustersITS(4);
678 esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
679 AliESDtrackCuts::kAny);
680 esdTrackCuts->SetMaxChi2PerClusterITS(1.);
683 // 7*(0.0085+0.0026/pt^1.55)
684 esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0595+0.0182/pt^1.55");
687 esdTrackCuts->SetRequireITSPid(kTRUE);
693 //____________________________________________________________________
694 Int_t AliESDtrackCuts::GetReferenceMultiplicity(AliESDEvent* esd, Bool_t tpcOnly)
696 // Gets reference multiplicity following the standard cuts and a defined fiducial volume
697 // tpcOnly = kTRUE -> consider TPC-only tracks
698 // = kFALSE -> consider global tracks
702 Printf("AliESDtrackCuts::GetReferenceMultiplicity: Not implemented for global tracks!");
706 AliESDtrackCuts* esdTrackCuts = GetStandardTPCOnlyTrackCuts();
707 esdTrackCuts->SetEtaRange(-0.8, 0.8);
708 esdTrackCuts->SetPtRange(0.15);
710 Int_t nTracks = esdTrackCuts->CountAcceptedTracks(esd);
718 //____________________________________________________________________
719 Float_t AliESDtrackCuts::GetSigmaToVertex(AliESDtrack* esdTrack)
721 // Calculates the number of sigma to the vertex.
726 esdTrack->GetImpactParameters(b,bCov);
728 if (bCov[0]<=0 || bCov[2]<=0) {
729 AliDebugClass(1, "Estimated b resolution lower or equal zero!");
730 bCov[0]=0; bCov[2]=0;
732 bRes[0] = TMath::Sqrt(bCov[0]);
733 bRes[1] = TMath::Sqrt(bCov[2]);
735 // -----------------------------------
736 // How to get to a n-sigma cut?
738 // The accumulated statistics from 0 to d is
740 // -> Erf(d/Sqrt(2)) for a 1-dim gauss (d = n_sigma)
741 // -> 1 - Exp(-d**2) for a 2-dim gauss (d*d = dx*dx + dy*dy != n_sigma)
743 // It means that for a 2-dim gauss: n_sigma(d) = Sqrt(2)*ErfInv(1 - Exp((-d**2)/2)
744 // Can this be expressed in a different way?
746 if (bRes[0] == 0 || bRes[1] ==0)
749 Float_t d = TMath::Sqrt(TMath::Power(b[0]/bRes[0],2) + TMath::Power(b[1]/bRes[1],2));
751 // work around precision problem
752 // if d is too big, TMath::Exp(...) gets 0, and TMath::ErfInverse(1) that should be infinite, gets 0 :(
753 // 1e-15 corresponds to nsigma ~ 7.7
754 if (TMath::Exp(-d * d / 2) < 1e-15)
757 Float_t nSigma = TMath::ErfInverse(1 - TMath::Exp(-d * d / 2)) * TMath::Sqrt(2);
761 void AliESDtrackCuts::EnableNeededBranches(TTree* tree)
763 // enables the branches needed by AcceptTrack, for a list see comment of AcceptTrack
765 tree->SetBranchStatus("fTracks.fFlags", 1);
766 tree->SetBranchStatus("fTracks.fITSncls", 1);
767 tree->SetBranchStatus("fTracks.fTPCncls", 1);
768 tree->SetBranchStatus("fTracks.fITSchi2", 1);
769 tree->SetBranchStatus("fTracks.fTPCchi2", 1);
770 tree->SetBranchStatus("fTracks.fC*", 1);
771 tree->SetBranchStatus("fTracks.fD", 1);
772 tree->SetBranchStatus("fTracks.fZ", 1);
773 tree->SetBranchStatus("fTracks.fCdd", 1);
774 tree->SetBranchStatus("fTracks.fCdz", 1);
775 tree->SetBranchStatus("fTracks.fCzz", 1);
776 tree->SetBranchStatus("fTracks.fP*", 1);
777 tree->SetBranchStatus("fTracks.fR*", 1);
778 tree->SetBranchStatus("fTracks.fKinkIndexes*", 1);
781 //____________________________________________________________________
782 Bool_t AliESDtrackCuts::AcceptTrack(AliESDtrack* esdTrack)
785 // figure out if the tracks survives all the track cuts defined
787 // the different quality parameter and kinematic values are first
788 // retrieved from the track. then it is found out what cuts the
789 // track did not survive and finally the cuts are imposed.
791 // this function needs the following branches:
797 // fTracks.fC //GetExternalCovariance
798 // fTracks.fD //GetImpactParameters
799 // fTracks.fZ //GetImpactParameters
800 // fTracks.fCdd //GetImpactParameters
801 // fTracks.fCdz //GetImpactParameters
802 // fTracks.fCzz //GetImpactParameters
803 // fTracks.fP //GetPxPyPz
804 // fTracks.fR //GetMass
805 // fTracks.fP //GetMass
806 // fTracks.fKinkIndexes
808 UInt_t status = esdTrack->GetStatus();
810 // getting quality parameters from the ESD track
811 Int_t nClustersITS = esdTrack->GetITSclusters(0);
812 Int_t nClustersTPC = -1;
813 if(fCutRequireTPCStandAlone) {
814 nClustersTPC = esdTrack->GetTPCNclsIter1();
817 nClustersTPC = esdTrack->GetTPCclusters(0);
820 Int_t nClustersTPCShared = esdTrack->GetTPCnclsS();
821 Float_t fracClustersTPCShared = -1.;
823 Float_t chi2PerClusterITS = -1;
824 Float_t chi2PerClusterTPC = -1;
826 chi2PerClusterITS = esdTrack->GetITSchi2()/Float_t(nClustersITS);
827 if (nClustersTPC!=0) {
828 if(fCutRequireTPCStandAlone) {
829 chi2PerClusterTPC = esdTrack->GetTPCchi2Iter1()/Float_t(nClustersTPC);
831 chi2PerClusterTPC = esdTrack->GetTPCchi2()/Float_t(nClustersTPC);
833 fracClustersTPCShared = Float_t(nClustersTPCShared)/Float_t(nClustersTPC);
837 esdTrack->GetExternalCovariance(extCov);
839 // getting the track to vertex parameters
840 Float_t nSigmaToVertex = GetSigmaToVertex(esdTrack);
844 esdTrack->GetImpactParameters(b,bCov);
845 if (bCov[0]<=0 || bCov[2]<=0) {
846 AliDebug(1, "Estimated b resolution lower or equal zero!");
847 bCov[0]=0; bCov[2]=0;
851 // set pt-dependent DCA cuts, if requested
852 SetPtDepDCACuts(esdTrack->Pt());
855 Float_t dcaToVertexXY = b[0];
856 Float_t dcaToVertexZ = b[1];
858 Float_t dcaToVertex = -1;
860 if (fCutDCAToVertex2D)
862 dcaToVertex = TMath::Sqrt(dcaToVertexXY*dcaToVertexXY/fCutMaxDCAToVertexXY/fCutMaxDCAToVertexXY + dcaToVertexZ*dcaToVertexZ/fCutMaxDCAToVertexZ/fCutMaxDCAToVertexZ);
865 dcaToVertex = TMath::Sqrt(dcaToVertexXY*dcaToVertexXY + dcaToVertexZ*dcaToVertexZ);
867 // getting the kinematic variables of the track
868 // (assuming the mass is known)
870 esdTrack->GetPxPyPz(p);
872 Float_t momentum = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2) + TMath::Power(p[2],2));
873 Float_t pt = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2));
874 Float_t energy = TMath::Sqrt(TMath::Power(esdTrack->GetMass(),2) + TMath::Power(momentum,2));
876 //y-eta related calculations
879 if((momentum != TMath::Abs(p[2]))&&(momentum != 0))
880 eta = 0.5*TMath::Log((momentum + p[2])/(momentum - p[2]));
881 if((energy != TMath::Abs(p[2]))&&(momentum != 0))
882 y = 0.5*TMath::Log((energy + p[2])/(energy - p[2]));
886 Printf("AliESDtrackCuts::AcceptTrack: WARNING: GetSigma1Pt2() returns negative value for external covariance matrix element fC[14]: %f. Corrupted track information, track will not be accepted!", extCov[14]);
889 Float_t relUncertainty1Pt = TMath::Sqrt(extCov[14])*pt;
891 //########################################################################
895 for (Int_t i=0; i<kNCuts; i++) cuts[i]=kFALSE;
897 // track quality cuts
898 if (fCutRequireTPCRefit && (status&AliESDtrack::kTPCrefit)==0)
900 if (fCutRequireTPCStandAlone && (status&AliESDtrack::kTPCin)==0)
902 if (fCutRequireITSRefit && (status&AliESDtrack::kITSrefit)==0)
904 if (nClustersTPC<fCutMinNClusterTPC)
906 if (nClustersITS<fCutMinNClusterITS)
908 if (chi2PerClusterTPC>fCutMaxChi2PerClusterTPC)
910 if (chi2PerClusterITS>fCutMaxChi2PerClusterITS)
912 if (extCov[0] > fCutMaxC11)
914 if (extCov[2] > fCutMaxC22)
916 if (extCov[5] > fCutMaxC33)
918 if (extCov[9] > fCutMaxC44)
920 if (extCov[14] > fCutMaxC55)
922 if (nSigmaToVertex > fCutNsigmaToVertex && fCutSigmaToVertexRequired)
924 // if n sigma could not be calculated
925 if (nSigmaToVertex<0 && fCutSigmaToVertexRequired)
927 if (!fCutAcceptKinkDaughters && esdTrack->GetKinkIndex(0)>0)
929 // track kinematics cut
930 if((momentum < fPMin) || (momentum > fPMax))
932 if((pt < fPtMin) || (pt > fPtMax))
934 if((p[0] < fPxMin) || (p[0] > fPxMax))
936 if((p[1] < fPyMin) || (p[1] > fPyMax))
938 if((p[2] < fPzMin) || (p[2] > fPzMax))
940 if((eta < fEtaMin) || (eta > fEtaMax))
942 if((y < fRapMin) || (y > fRapMax))
944 if (fCutDCAToVertex2D && dcaToVertex > 1)
946 if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexXY) > fCutMaxDCAToVertexXY)
948 if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexZ) > fCutMaxDCAToVertexZ)
950 if (fCutDCAToVertex2D && fCutMinDCAToVertexXY > 0 && fCutMinDCAToVertexZ > 0 && dcaToVertexXY*dcaToVertexXY/fCutMinDCAToVertexXY/fCutMinDCAToVertexXY + dcaToVertexZ*dcaToVertexZ/fCutMinDCAToVertexZ/fCutMinDCAToVertexZ < 1)
952 if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexXY) < fCutMinDCAToVertexXY)
954 if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexZ) < fCutMinDCAToVertexZ)
957 for (Int_t i = 0; i < 3; i++)
958 cuts[28+i] = !CheckITSClusterRequirement(fCutClusterRequirementITS[i], esdTrack->HasPointOnITSLayer(i*2), esdTrack->HasPointOnITSLayer(i*2+1));
960 if(fCutRequireITSStandAlone || fCutRequireITSpureSA){
961 if ((status & AliESDtrack::kITSin) == 0 || (status & AliESDtrack::kTPCin)){
965 // ITS standalone tracks
966 if(fCutRequireITSStandAlone && !fCutRequireITSpureSA){
967 if(status & AliESDtrack::kITSpureSA) cuts[31] = kTRUE;
968 }else if(fCutRequireITSpureSA){
969 if(!(status & AliESDtrack::kITSpureSA)) cuts[31] = kTRUE;
974 if (relUncertainty1Pt > fCutMaxRel1PtUncertainty)
977 if (!fCutAcceptSharedTPCClusters && nClustersTPCShared!=0)
980 if (fracClustersTPCShared > fCutMaxFractionSharedTPCClusters)
983 if(fCutRequireITSPid){
984 UChar_t clumap=esdTrack->GetITSClusterMap();
985 Int_t nPointsForPid=0;
986 for(Int_t i=2; i<6; i++){
987 if(clumap&(1<<i)) ++nPointsForPid;
989 if(nPointsForPid<3) cuts[35] = kTRUE;
993 for (Int_t i=0; i<kNCuts; i++)
994 if (cuts[i]) {cut = kTRUE;}
997 //########################################################################
998 // filling histograms
1000 fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n tracks")));
1002 fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n cut tracks")));
1004 for (Int_t i=0; i<kNCuts; i++) {
1006 fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin(fgkCutNames[i])));
1008 for (Int_t j=i; j<kNCuts; j++) {
1009 if (cuts[i] && cuts[j]) {
1010 Float_t xC = fhCutCorrelation->GetXaxis()->GetBinCenter(fhCutCorrelation->GetXaxis()->FindBin(fgkCutNames[i]));
1011 Float_t yC = fhCutCorrelation->GetYaxis()->GetBinCenter(fhCutCorrelation->GetYaxis()->FindBin(fgkCutNames[j]));
1012 fhCutCorrelation->Fill(xC, yC);
1018 // now we loop over the filling of the histograms twice: once "before" the cut, once "after"
1019 // the code is not in a function due to too many local variables that would need to be passed
1021 for (Int_t id = 0; id < 2; id++)
1023 // id = 0 --> before cut
1024 // id = 1 --> after cut
1028 fhNClustersITS[id]->Fill(nClustersITS);
1029 fhNClustersTPC[id]->Fill(nClustersTPC);
1030 fhChi2PerClusterITS[id]->Fill(chi2PerClusterITS);
1031 fhChi2PerClusterTPC[id]->Fill(chi2PerClusterTPC);
1033 fhC11[id]->Fill(extCov[0]);
1034 fhC22[id]->Fill(extCov[2]);
1035 fhC33[id]->Fill(extCov[5]);
1036 fhC44[id]->Fill(extCov[9]);
1037 fhC55[id]->Fill(extCov[14]);
1039 fhRel1PtUncertainty[id]->Fill(relUncertainty1Pt);
1042 fhEta[id]->Fill(eta);
1045 bRes[0] = TMath::Sqrt(bCov[0]);
1046 bRes[1] = TMath::Sqrt(bCov[2]);
1048 fhDZ[id]->Fill(b[1]);
1049 fhDXY[id]->Fill(b[0]);
1050 fhDXYDZ[id]->Fill(dcaToVertex);
1051 fhDXYvsDZ[id]->Fill(b[1],b[0]);
1053 if (bRes[0]!=0 && bRes[1]!=0) {
1054 fhDZNormalized[id]->Fill(b[1]/bRes[1]);
1055 fhDXYNormalized[id]->Fill(b[0]/bRes[0]);
1056 fhDXYvsDZNormalized[id]->Fill(b[1]/bRes[1], b[0]/bRes[0]);
1057 fhNSigmaToVertex[id]->Fill(nSigmaToVertex);
1069 //____________________________________________________________________
1070 Bool_t AliESDtrackCuts::CheckITSClusterRequirement(ITSClusterRequirement req, Bool_t clusterL1, Bool_t clusterL2)
1072 // checks if the cluster requirement is fullfilled (in this case: return kTRUE)
1076 case kOff: return kTRUE;
1077 case kNone: return !clusterL1 && !clusterL2;
1078 case kAny: return clusterL1 || clusterL2;
1079 case kFirst: return clusterL1;
1080 case kOnlyFirst: return clusterL1 && !clusterL2;
1081 case kSecond: return clusterL2;
1082 case kOnlySecond: return clusterL2 && !clusterL1;
1083 case kBoth: return clusterL1 && clusterL2;
1089 //____________________________________________________________________
1090 AliESDtrack* AliESDtrackCuts::GetTPCOnlyTrack(AliESDEvent* esd, Int_t iTrack)
1093 // Utility function to
1094 // create a TPC only track from the given esd track
1096 // IMPORTANT: The track has to be deleted by the user
1098 // NB. most of the functionality to get a TPC only track from an ESD track is in AliESDtrack, where it should be
1099 // there are only missing propagations here that are needed for old data
1100 // this function will therefore become obsolete
1102 // adapted from code provided by CKB
1104 if (!esd->GetPrimaryVertexTPC())
1105 return 0; // No TPC vertex no TPC tracks
1107 if(!esd->GetPrimaryVertexTPC()->GetStatus())
1108 return 0; // TPC Vertex is created by default in AliESDEvent, do not use in this case
1110 AliESDtrack* track = esd->GetTrack(iTrack);
1114 AliESDtrack *tpcTrack = new AliESDtrack();
1116 // only true if we have a tpc track
1117 if (!track->FillTPCOnlyTrack(*tpcTrack))
1123 // propagate to Vertex
1124 // not needed for normal reconstructed ESDs...
1125 // Double_t pTPC[2],covTPC[3];
1126 // tpcTrack->PropagateToDCA(esd->GetPrimaryVertexTPC(), esd->GetMagneticField(), 10000, pTPC, covTPC);
1131 //____________________________________________________________________
1132 TObjArray* AliESDtrackCuts::GetAcceptedTracks(AliESDEvent* esd,Bool_t bTPC)
1135 // returns an array of all tracks that pass the cuts
1136 // or an array of TPC only tracks (propagated to the TPC vertex during reco)
1137 // tracks that pass the cut
1139 TObjArray* acceptedTracks = new TObjArray();
1141 // loop over esd tracks
1142 for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
1144 if(!esd->GetPrimaryVertexTPC())return acceptedTracks; // No TPC vertex no TPC tracks
1145 if(!esd->GetPrimaryVertexTPC()->GetStatus())return acceptedTracks; // No proper TPC vertex, only the default
1147 AliESDtrack *tpcTrack = GetTPCOnlyTrack(esd, iTrack);
1151 if (AcceptTrack(tpcTrack)) {
1152 acceptedTracks->Add(tpcTrack);
1159 AliESDtrack* track = esd->GetTrack(iTrack);
1160 if(AcceptTrack(track))
1161 acceptedTracks->Add(track);
1164 if(bTPC)acceptedTracks->SetOwner(kTRUE);
1165 return acceptedTracks;
1168 //____________________________________________________________________
1169 Int_t AliESDtrackCuts::CountAcceptedTracks(AliESDEvent* esd)
1172 // returns an the number of tracks that pass the cuts
1177 // loop over esd tracks
1178 for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
1179 AliESDtrack* track = esd->GetTrack(iTrack);
1180 if (AcceptTrack(track))
1187 //____________________________________________________________________
1188 void AliESDtrackCuts::DefineHistograms(Int_t color) {
1190 // diagnostics histograms are defined
1193 fHistogramsOn=kTRUE;
1195 Bool_t oldStatus = TH1::AddDirectoryStatus();
1196 TH1::AddDirectory(kFALSE);
1198 //###################################################################################
1199 // defining histograms
1201 fhCutStatistics = new TH1F("cut_statistics","cut statistics",kNCuts+4,-0.5,kNCuts+3.5);
1203 fhCutStatistics->GetXaxis()->SetBinLabel(1,"n tracks");
1204 fhCutStatistics->GetXaxis()->SetBinLabel(2,"n cut tracks");
1206 fhCutCorrelation = new TH2F("cut_correlation","cut correlation",kNCuts,-0.5,kNCuts-0.5,kNCuts,-0.5,kNCuts-0.5);;
1208 for (Int_t i=0; i<kNCuts; i++) {
1209 fhCutStatistics->GetXaxis()->SetBinLabel(i+4,fgkCutNames[i]);
1210 fhCutCorrelation->GetXaxis()->SetBinLabel(i+1,fgkCutNames[i]);
1211 fhCutCorrelation->GetYaxis()->SetBinLabel(i+1,fgkCutNames[i]);
1214 fhCutStatistics ->SetLineColor(color);
1215 fhCutCorrelation ->SetLineColor(color);
1216 fhCutStatistics ->SetLineWidth(2);
1217 fhCutCorrelation ->SetLineWidth(2);
1219 for (Int_t i=0; i<2; i++) {
1220 fhNClustersITS[i] = new TH1F("nClustersITS" ,"",8,-0.5,7.5);
1221 fhNClustersTPC[i] = new TH1F("nClustersTPC" ,"",165,-0.5,164.5);
1222 fhChi2PerClusterITS[i] = new TH1F("chi2PerClusterITS","",500,0,10);
1223 fhChi2PerClusterTPC[i] = new TH1F("chi2PerClusterTPC","",500,0,10);
1225 fhC11[i] = new TH1F("covMatrixDiagonal11","",2000,0,20);
1226 fhC22[i] = new TH1F("covMatrixDiagonal22","",2000,0,20);
1227 fhC33[i] = new TH1F("covMatrixDiagonal33","",1000,0,0.1);
1228 fhC44[i] = new TH1F("covMatrixDiagonal44","",1000,0,0.1);
1229 fhC55[i] = new TH1F("covMatrixDiagonal55","",1000,0,5);
1231 fhRel1PtUncertainty[i] = new TH1F("rel1PtUncertainty","",1000,0,5);
1233 fhDXY[i] = new TH1F("dXY" ,"",500,-10,10);
1234 fhDZ[i] = new TH1F("dZ" ,"",500,-10,10);
1235 fhDXYDZ[i] = new TH1F("dXYDZ" ,"",500,0,10);
1236 fhDXYvsDZ[i] = new TH2F("dXYvsDZ","",200,-10,10,200,-10,10);
1238 fhDXYNormalized[i] = new TH1F("dXYNormalized" ,"",500,-10,10);
1239 fhDZNormalized[i] = new TH1F("dZNormalized" ,"",500,-10,10);
1240 fhDXYvsDZNormalized[i] = new TH2F("dXYvsDZNormalized","",200,-10,10,200,-10,10);
1242 fhNSigmaToVertex[i] = new TH1F("nSigmaToVertex","",500,0,10);
1244 fhPt[i] = new TH1F("pt" ,"p_{T} distribution;p_{T} (GeV/c)", 800, 0.0, 10.0);
1245 fhEta[i] = new TH1F("eta" ,"#eta distribution;#eta",40,-2.0,2.0);
1247 fhNClustersITS[i]->SetTitle("n ITS clusters");
1248 fhNClustersTPC[i]->SetTitle("n TPC clusters");
1249 fhChi2PerClusterITS[i]->SetTitle("#Chi^{2} per ITS cluster");
1250 fhChi2PerClusterTPC[i]->SetTitle("#Chi^{2} per TPC cluster");
1252 fhC11[i]->SetTitle("cov 11 : #sigma_{y}^{2} [cm^{2}]");
1253 fhC22[i]->SetTitle("cov 22 : #sigma_{z}^{2} [cm^{2}]");
1254 fhC33[i]->SetTitle("cov 33 : #sigma_{sin(#phi)}^{2}");
1255 fhC44[i]->SetTitle("cov 44 : #sigma_{tan(#theta_{dip})}^{2}");
1256 fhC55[i]->SetTitle("cov 55 : #sigma_{1/p_{T}}^{2} [(c/GeV)^2]");
1258 fhRel1PtUncertainty[i]->SetTitle("rel. uncertainty of 1/p_{T}");
1260 fhDXY[i]->SetXTitle("transverse impact parameter (cm)");
1261 fhDZ[i]->SetXTitle("longitudinal impact parameter (cm)");
1262 fhDXYDZ[i]->SetTitle("absolute impact parameter;sqrt(dXY**2 + dZ**2) (cm)");
1263 fhDXYvsDZ[i]->SetXTitle("longitudinal impact parameter (cm)");
1264 fhDXYvsDZ[i]->SetYTitle("transverse impact parameter (cm)");
1266 fhDXYNormalized[i]->SetTitle("normalized trans impact par (n#sigma)");
1267 fhDZNormalized[i]->SetTitle("normalized long impact par (n#sigma)");
1268 fhDXYvsDZNormalized[i]->SetTitle("normalized long impact par (n#sigma)");
1269 fhDXYvsDZNormalized[i]->SetYTitle("normalized trans impact par (n#sigma)");
1270 fhNSigmaToVertex[i]->SetTitle("n #sigma to vertex");
1272 fhNClustersITS[i]->SetLineColor(color); fhNClustersITS[i]->SetLineWidth(2);
1273 fhNClustersTPC[i]->SetLineColor(color); fhNClustersTPC[i]->SetLineWidth(2);
1274 fhChi2PerClusterITS[i]->SetLineColor(color); fhChi2PerClusterITS[i]->SetLineWidth(2);
1275 fhChi2PerClusterTPC[i]->SetLineColor(color); fhChi2PerClusterTPC[i]->SetLineWidth(2);
1277 fhC11[i]->SetLineColor(color); fhC11[i]->SetLineWidth(2);
1278 fhC22[i]->SetLineColor(color); fhC22[i]->SetLineWidth(2);
1279 fhC33[i]->SetLineColor(color); fhC33[i]->SetLineWidth(2);
1280 fhC44[i]->SetLineColor(color); fhC44[i]->SetLineWidth(2);
1281 fhC55[i]->SetLineColor(color); fhC55[i]->SetLineWidth(2);
1283 fhRel1PtUncertainty[i]->SetLineColor(color); fhRel1PtUncertainty[i]->SetLineWidth(2);
1285 fhDXY[i]->SetLineColor(color); fhDXY[i]->SetLineWidth(2);
1286 fhDZ[i]->SetLineColor(color); fhDZ[i]->SetLineWidth(2);
1287 fhDXYDZ[i]->SetLineColor(color); fhDXYDZ[i]->SetLineWidth(2);
1289 fhDXYNormalized[i]->SetLineColor(color); fhDXYNormalized[i]->SetLineWidth(2);
1290 fhDZNormalized[i]->SetLineColor(color); fhDZNormalized[i]->SetLineWidth(2);
1291 fhNSigmaToVertex[i]->SetLineColor(color); fhNSigmaToVertex[i]->SetLineWidth(2);
1294 // The number of sigmas to the vertex is per definition gaussian
1295 ffDTheoretical = new TF1("nSigmaToVertexTheoretical","([0]/2.506628274)*exp(-(x**2)/2)",0,50);
1296 ffDTheoretical->SetParameter(0,1);
1298 TH1::AddDirectory(oldStatus);
1301 //____________________________________________________________________
1302 Bool_t AliESDtrackCuts::LoadHistograms(const Char_t* dir)
1305 // loads the histograms from a file
1306 // if dir is empty a directory with the name of this object is taken (like in SaveHistogram)
1312 if (!gDirectory->cd(dir))
1315 ffDTheoretical = dynamic_cast<TF1*> (gDirectory->Get("nSigmaToVertexTheory"));
1317 fhCutStatistics = dynamic_cast<TH1F*> (gDirectory->Get("cut_statistics"));
1318 fhCutCorrelation = dynamic_cast<TH2F*> (gDirectory->Get("cut_correlation"));
1320 for (Int_t i=0; i<2; i++) {
1323 gDirectory->cd("before_cuts");
1326 gDirectory->cd("after_cuts");
1328 fhNClustersITS[i] = dynamic_cast<TH1F*> (gDirectory->Get("nClustersITS" ));
1329 fhNClustersTPC[i] = dynamic_cast<TH1F*> (gDirectory->Get("nClustersTPC" ));
1330 fhChi2PerClusterITS[i] = dynamic_cast<TH1F*> (gDirectory->Get("chi2PerClusterITS"));
1331 fhChi2PerClusterTPC[i] = dynamic_cast<TH1F*> (gDirectory->Get("chi2PerClusterTPC"));
1333 fhC11[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal11"));
1334 fhC22[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal22"));
1335 fhC33[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal33"));
1336 fhC44[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal44"));
1337 fhC55[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal55"));
1339 fhRel1PtUncertainty[i] = dynamic_cast<TH1F*> (gDirectory->Get("rel1PtUncertainty"));
1341 fhDXY[i] = dynamic_cast<TH1F*> (gDirectory->Get("dXY" ));
1342 fhDZ[i] = dynamic_cast<TH1F*> (gDirectory->Get("dZ" ));
1343 fhDXYDZ[i] = dynamic_cast<TH1F*> (gDirectory->Get("dXYDZ"));
1344 fhDXYvsDZ[i] = dynamic_cast<TH2F*> (gDirectory->Get("dXYvsDZ"));
1346 fhDXYNormalized[i] = dynamic_cast<TH1F*> (gDirectory->Get("dXYNormalized" ));
1347 fhDZNormalized[i] = dynamic_cast<TH1F*> (gDirectory->Get("dZNormalized" ));
1348 fhDXYvsDZNormalized[i] = dynamic_cast<TH2F*> (gDirectory->Get("dXYvsDZNormalized"));
1349 fhNSigmaToVertex[i] = dynamic_cast<TH1F*> (gDirectory->Get("nSigmaToVertex"));
1351 fhPt[i] = dynamic_cast<TH1F*> (gDirectory->Get("pt"));
1352 fhEta[i] = dynamic_cast<TH1F*> (gDirectory->Get("eta"));
1354 gDirectory->cd("../");
1357 gDirectory->cd("..");
1362 //____________________________________________________________________
1363 void AliESDtrackCuts::SaveHistograms(const Char_t* dir) {
1365 // saves the histograms in a directory (dir)
1368 if (!fHistogramsOn) {
1369 AliDebug(0, "Histograms not on - cannot save histograms!!!");
1376 gDirectory->mkdir(dir);
1377 gDirectory->cd(dir);
1379 gDirectory->mkdir("before_cuts");
1380 gDirectory->mkdir("after_cuts");
1382 // a factor of 2 is needed since n sigma is positive
1383 ffDTheoretical->SetParameter(0,2*fhNSigmaToVertex[0]->Integral("width"));
1384 ffDTheoretical->Write("nSigmaToVertexTheory");
1386 fhCutStatistics->Write();
1387 fhCutCorrelation->Write();
1389 for (Int_t i=0; i<2; i++) {
1391 gDirectory->cd("before_cuts");
1393 gDirectory->cd("after_cuts");
1395 fhNClustersITS[i] ->Write();
1396 fhNClustersTPC[i] ->Write();
1397 fhChi2PerClusterITS[i] ->Write();
1398 fhChi2PerClusterTPC[i] ->Write();
1406 fhRel1PtUncertainty[i] ->Write();
1410 fhDXYDZ[i] ->Write();
1411 fhDXYvsDZ[i] ->Write();
1413 fhDXYNormalized[i] ->Write();
1414 fhDZNormalized[i] ->Write();
1415 fhDXYvsDZNormalized[i] ->Write();
1416 fhNSigmaToVertex[i] ->Write();
1421 gDirectory->cd("../");
1424 gDirectory->cd("../");
1427 //____________________________________________________________________
1428 void AliESDtrackCuts::DrawHistograms()
1430 // draws some histograms
1432 TCanvas* canvas1 = new TCanvas(Form("%s_1", GetName()), "Track Quality Results1", 800, 800);
1433 canvas1->Divide(2, 2);
1436 fhNClustersTPC[0]->SetStats(kFALSE);
1437 fhNClustersTPC[0]->Draw();
1440 fhChi2PerClusterTPC[0]->SetStats(kFALSE);
1441 fhChi2PerClusterTPC[0]->Draw();
1444 fhNSigmaToVertex[0]->SetStats(kFALSE);
1445 fhNSigmaToVertex[0]->GetXaxis()->SetRangeUser(0, 10);
1446 fhNSigmaToVertex[0]->Draw();
1448 canvas1->SaveAs(Form("%s_%s.gif", GetName(), canvas1->GetName()));
1450 TCanvas* canvas2 = new TCanvas(Form("%s_2", GetName()), "Track Quality Results2", 1200, 800);
1451 canvas2->Divide(3, 2);
1454 fhC11[0]->SetStats(kFALSE);
1459 fhC22[0]->SetStats(kFALSE);
1464 fhC33[0]->SetStats(kFALSE);
1469 fhC44[0]->SetStats(kFALSE);
1474 fhC55[0]->SetStats(kFALSE);
1479 fhRel1PtUncertainty[0]->SetStats(kFALSE);
1481 fhRel1PtUncertainty[0]->Draw();
1483 canvas2->SaveAs(Form("%s_%s.gif", GetName(), canvas2->GetName()));
1485 TCanvas* canvas3 = new TCanvas(Form("%s_3", GetName()), "Track Quality Results3", 1200, 800);
1486 canvas3->Divide(3, 2);
1489 fhDXY[0]->SetStats(kFALSE);
1494 fhDZ[0]->SetStats(kFALSE);
1499 fhDXYvsDZ[0]->SetStats(kFALSE);
1501 gPad->SetRightMargin(0.15);
1502 fhDXYvsDZ[0]->Draw("COLZ");
1505 fhDXYNormalized[0]->SetStats(kFALSE);
1507 fhDXYNormalized[0]->Draw();
1510 fhDZNormalized[0]->SetStats(kFALSE);
1512 fhDZNormalized[0]->Draw();
1515 fhDXYvsDZNormalized[0]->SetStats(kFALSE);
1517 gPad->SetRightMargin(0.15);
1518 fhDXYvsDZNormalized[0]->Draw("COLZ");
1520 canvas3->SaveAs(Form("%s_%s.gif", GetName(), canvas3->GetName()));
1522 TCanvas* canvas4 = new TCanvas(Form("%s_4", GetName()), "Track Quality Results4", 800, 500);
1523 canvas4->Divide(2, 1);
1526 fhCutStatistics->SetStats(kFALSE);
1527 fhCutStatistics->LabelsOption("v");
1528 gPad->SetBottomMargin(0.3);
1529 fhCutStatistics->Draw();
1532 fhCutCorrelation->SetStats(kFALSE);
1533 fhCutCorrelation->LabelsOption("v");
1534 gPad->SetBottomMargin(0.3);
1535 gPad->SetLeftMargin(0.3);
1536 fhCutCorrelation->Draw("COLZ");
1538 canvas4->SaveAs(Form("%s_%s.gif", GetName(), canvas4->GetName()));
1541 fhDXYvsDZNormalized[0]->SetStats(kFALSE);
1542 fhDXYvsDZNormalized[0]->DrawCopy("COLZ");
1545 fhNClustersTPC[0]->SetStats(kFALSE);
1546 fhNClustersTPC[0]->DrawCopy();
1549 fhChi2PerClusterITS[0]->SetStats(kFALSE);
1550 fhChi2PerClusterITS[0]->DrawCopy();
1551 fhChi2PerClusterITS[1]->SetLineColor(2);
1552 fhChi2PerClusterITS[1]->DrawCopy("SAME");
1555 fhChi2PerClusterTPC[0]->SetStats(kFALSE);
1556 fhChi2PerClusterTPC[0]->DrawCopy();
1557 fhChi2PerClusterTPC[1]->SetLineColor(2);
1558 fhChi2PerClusterTPC[1]->DrawCopy("SAME");*/
1560 //--------------------------------------------------------------------------
1561 void AliESDtrackCuts::SetPtDepDCACuts(Double_t pt) {
1563 // set the pt-dependent DCA cuts
1566 if(f1CutMaxDCAToVertexXYPtDep) {
1567 fCutMaxDCAToVertexXY=f1CutMaxDCAToVertexXYPtDep->Eval(pt);
1570 if(f1CutMaxDCAToVertexZPtDep) {
1571 fCutMaxDCAToVertexZ=f1CutMaxDCAToVertexZPtDep->Eval(pt);
1574 if(f1CutMinDCAToVertexXYPtDep) {
1575 fCutMinDCAToVertexXY=f1CutMinDCAToVertexXYPtDep->Eval(pt);
1578 if(f1CutMinDCAToVertexZPtDep) {
1579 fCutMinDCAToVertexZ=f1CutMinDCAToVertexZPtDep->Eval(pt);
1588 //--------------------------------------------------------------------------
1589 Bool_t AliESDtrackCuts::CheckPtDepDCA(TString dist,Bool_t print) const {
1591 // Check the correctness of the string syntax
1593 Bool_t retval=kTRUE;
1595 if(!dist.Contains("pt")) {
1596 if(print) printf("AliESDtrackCuts::CheckPtDepDCA(): string must contain \"pt\"\n");
1602 void AliESDtrackCuts::SetMaxDCAToVertexXYPtDep(const char *dist){
1603 if(!CheckPtDepDCA(dist,kTRUE)) return;
1604 if(f1CutMaxDCAToVertexXYPtDep)delete f1CutMaxDCAToVertexXYPtDep;
1605 fCutMaxDCAToVertexXYPtDep = dist;
1607 tmp.ReplaceAll("pt","x");
1608 f1CutMaxDCAToVertexXYPtDep = new TFormula("f1CutMaxDCAToVertexXYPtDep",tmp.Data());
1612 void AliESDtrackCuts::SetMaxDCAToVertexZPtDep(const char *dist){
1613 if(!CheckPtDepDCA(dist,kTRUE)) return;
1614 if(f1CutMaxDCAToVertexZPtDep)delete f1CutMaxDCAToVertexZPtDep;
1615 fCutMaxDCAToVertexZPtDep = dist;
1617 tmp.ReplaceAll("pt","x");
1618 f1CutMaxDCAToVertexZPtDep = new TFormula("f1CutMaxDCAToVertexZPtDep",tmp.Data());
1624 void AliESDtrackCuts::SetMinDCAToVertexXYPtDep(const char *dist){
1625 if(!CheckPtDepDCA(dist,kTRUE)) return;
1626 if(f1CutMinDCAToVertexXYPtDep)delete f1CutMinDCAToVertexXYPtDep;
1627 fCutMinDCAToVertexXYPtDep = dist;
1629 tmp.ReplaceAll("pt","x");
1630 f1CutMinDCAToVertexXYPtDep = new TFormula("f1CutMinDCAToVertexXYPtDep",tmp.Data());
1635 void AliESDtrackCuts::SetMinDCAToVertexZPtDep(const char *dist){
1636 if(!CheckPtDepDCA(dist,kTRUE)) return;
1637 if(f1CutMinDCAToVertexZPtDep)delete f1CutMinDCAToVertexZPtDep;
1638 fCutMinDCAToVertexZPtDep = dist;
1640 tmp.ReplaceAll("pt","x");
1641 f1CutMinDCAToVertexZPtDep = new TFormula("f1CutMinDCAToVertexZPtDep",tmp.Data());