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"
71 //____________________________________________________________________
72 AliESDtrackCuts::AliESDtrackCuts(const Char_t* name, const Char_t* title) : AliAnalysisCuts(name,title),
73 fCutMinNClusterTPC(0),
74 fCutMinNClusterITS(0),
75 fCutMaxChi2PerClusterTPC(0),
76 fCutMaxChi2PerClusterITS(0),
82 fCutMaxRel1PtUncertainty(0),
83 fCutAcceptKinkDaughters(0),
84 fCutAcceptSharedTPCClusters(0),
85 fCutMaxFractionSharedTPCClusters(0),
86 fCutRequireTPCRefit(0),
87 fCutRequireTPCStandAlone(0),
88 fCutRequireITSRefit(0),
89 fCutRequireITSStandAlone(0),
90 fCutRejectITSpureSA(0),
91 fCutNsigmaToVertex(0),
92 fCutSigmaToVertexRequired(0),
93 fCutMaxDCAToVertexXY(0),
94 fCutMaxDCAToVertexZ(0),
95 fCutMinDCAToVertexXY(0),
96 fCutMinDCAToVertexZ(0),
97 fCutMaxDCAToVertexXYPtDep(""),
98 fCutMaxDCAToVertexZPtDep(""),
99 fCutMinDCAToVertexXYPtDep(""),
100 fCutMinDCAToVertexZPtDep(""),
101 f1CutMaxDCAToVertexXYPtDep(0x0),
102 f1CutMaxDCAToVertexZPtDep(0x0),
103 f1CutMinDCAToVertexXYPtDep(0x0),
104 f1CutMinDCAToVertexZPtDep(0x0),
105 fCutDCAToVertex2D(0),
131 //##############################################################################
132 // setting default cuts
133 SetMinNClustersTPC();
134 SetMinNClustersITS();
135 SetMaxChi2PerClusterTPC();
136 SetMaxChi2PerClusterITS();
137 SetMaxCovDiagonalElements();
138 SetMaxRel1PtUncertainty();
139 SetRequireTPCRefit();
140 SetRequireTPCStandAlone();
141 SetRequireITSRefit();
142 SetRequireITSStandAlone(kFALSE);
143 SetAcceptKinkDaughters();
144 SetMaxNsigmaToVertex();
145 SetMaxDCAToVertexXY();
146 SetMaxDCAToVertexZ();
148 SetMinDCAToVertexXY();
149 SetMinDCAToVertexZ();
157 SetClusterRequirementITS(kSPD);
158 SetClusterRequirementITS(kSDD);
159 SetClusterRequirementITS(kSSD);
164 //_____________________________________________________________________________
165 AliESDtrackCuts::AliESDtrackCuts(const AliESDtrackCuts &c) : AliAnalysisCuts(c),
166 fCutMinNClusterTPC(0),
167 fCutMinNClusterITS(0),
168 fCutMaxChi2PerClusterTPC(0),
169 fCutMaxChi2PerClusterITS(0),
175 fCutMaxRel1PtUncertainty(0),
176 fCutAcceptKinkDaughters(0),
177 fCutAcceptSharedTPCClusters(0),
178 fCutMaxFractionSharedTPCClusters(0),
179 fCutRequireTPCRefit(0),
180 fCutRequireTPCStandAlone(0),
181 fCutRequireITSRefit(0),
182 fCutRequireITSStandAlone(0),
183 fCutRejectITSpureSA(0),
184 fCutNsigmaToVertex(0),
185 fCutSigmaToVertexRequired(0),
186 fCutMaxDCAToVertexXY(0),
187 fCutMaxDCAToVertexZ(0),
188 fCutMinDCAToVertexXY(0),
189 fCutMinDCAToVertexZ(0),
190 fCutMaxDCAToVertexXYPtDep(""),
191 fCutMaxDCAToVertexZPtDep(""),
192 fCutMinDCAToVertexXYPtDep(""),
193 fCutMinDCAToVertexZPtDep(""),
194 f1CutMaxDCAToVertexXYPtDep(0x0),
195 f1CutMaxDCAToVertexZPtDep(0x0),
196 f1CutMinDCAToVertexXYPtDep(0x0),
197 f1CutMinDCAToVertexZPtDep(0x0),
198 fCutDCAToVertex2D(0),
222 ((AliESDtrackCuts &) c).Copy(*this);
225 AliESDtrackCuts::~AliESDtrackCuts()
231 for (Int_t i=0; i<2; i++) {
233 if (fhNClustersITS[i])
234 delete fhNClustersITS[i];
235 if (fhNClustersTPC[i])
236 delete fhNClustersTPC[i];
237 if (fhChi2PerClusterITS[i])
238 delete fhChi2PerClusterITS[i];
239 if (fhChi2PerClusterTPC[i])
240 delete fhChi2PerClusterTPC[i];
252 if (fhRel1PtUncertainty[i])
253 delete fhRel1PtUncertainty[i];
264 if (fhDXYNormalized[i])
265 delete fhDXYNormalized[i];
266 if (fhDZNormalized[i])
267 delete fhDZNormalized[i];
268 if (fhDXYvsDZNormalized[i])
269 delete fhDXYvsDZNormalized[i];
270 if (fhNSigmaToVertex[i])
271 delete fhNSigmaToVertex[i];
278 if(f1CutMaxDCAToVertexXYPtDep)delete f1CutMaxDCAToVertexXYPtDep;
279 f1CutMaxDCAToVertexXYPtDep = 0;
280 if( f1CutMaxDCAToVertexZPtDep) delete f1CutMaxDCAToVertexZPtDep;
281 f1CutMaxDCAToVertexZPtDep = 0;
282 if( f1CutMinDCAToVertexXYPtDep)delete f1CutMinDCAToVertexXYPtDep;
283 f1CutMinDCAToVertexXYPtDep = 0;
284 if(f1CutMinDCAToVertexZPtDep)delete f1CutMinDCAToVertexZPtDep;
285 f1CutMinDCAToVertexZPtDep = 0;
289 delete ffDTheoretical;
292 delete fhCutStatistics;
293 if (fhCutCorrelation)
294 delete fhCutCorrelation;
297 void AliESDtrackCuts::Init()
300 // sets everything to zero
303 fCutMinNClusterTPC = 0;
304 fCutMinNClusterITS = 0;
306 fCutMaxChi2PerClusterTPC = 0;
307 fCutMaxChi2PerClusterITS = 0;
309 for (Int_t i = 0; i < 3; i++)
310 fCutClusterRequirementITS[i] = kOff;
318 fCutMaxRel1PtUncertainty = 0;
320 fCutAcceptKinkDaughters = 0;
321 fCutAcceptSharedTPCClusters = 0;
322 fCutMaxFractionSharedTPCClusters = 0;
323 fCutRequireTPCRefit = 0;
324 fCutRequireTPCStandAlone = 0;
325 fCutRequireITSRefit = 0;
326 fCutRequireITSStandAlone = 0;
327 fCutRejectITSpureSA = 0;
329 fCutNsigmaToVertex = 0;
330 fCutSigmaToVertexRequired = 0;
331 fCutMaxDCAToVertexXY = 0;
332 fCutMaxDCAToVertexZ = 0;
333 fCutDCAToVertex2D = 0;
334 fCutMinDCAToVertexXY = 0;
335 fCutMinDCAToVertexZ = 0;
336 fCutMaxDCAToVertexXYPtDep = "";
337 fCutMaxDCAToVertexZPtDep = "";
338 fCutMinDCAToVertexXYPtDep = "";
339 fCutMinDCAToVertexZPtDep = "";
341 if(f1CutMaxDCAToVertexXYPtDep)delete f1CutMaxDCAToVertexXYPtDep;
342 f1CutMaxDCAToVertexXYPtDep = 0;
343 if( f1CutMaxDCAToVertexXYPtDep) delete f1CutMaxDCAToVertexXYPtDep;
344 f1CutMaxDCAToVertexXYPtDep = 0;
345 if( f1CutMaxDCAToVertexZPtDep) delete f1CutMaxDCAToVertexZPtDep;
346 f1CutMaxDCAToVertexZPtDep = 0;
347 if( f1CutMinDCAToVertexXYPtDep)delete f1CutMinDCAToVertexXYPtDep;
348 f1CutMinDCAToVertexXYPtDep = 0;
349 if(f1CutMinDCAToVertexZPtDep)delete f1CutMinDCAToVertexZPtDep;
350 f1CutMinDCAToVertexZPtDep = 0;
368 fHistogramsOn = kFALSE;
370 for (Int_t i=0; i<2; ++i)
372 fhNClustersITS[i] = 0;
373 fhNClustersTPC[i] = 0;
375 fhChi2PerClusterITS[i] = 0;
376 fhChi2PerClusterTPC[i] = 0;
384 fhRel1PtUncertainty[i] = 0;
391 fhDXYNormalized[i] = 0;
392 fhDZNormalized[i] = 0;
393 fhDXYvsDZNormalized[i] = 0;
394 fhNSigmaToVertex[i] = 0;
402 fhCutCorrelation = 0;
405 //_____________________________________________________________________________
406 AliESDtrackCuts &AliESDtrackCuts::operator=(const AliESDtrackCuts &c)
409 // Assignment operator
412 if (this != &c) ((AliESDtrackCuts &) c).Copy(*this);
416 //_____________________________________________________________________________
417 void AliESDtrackCuts::Copy(TObject &c) const
423 AliESDtrackCuts& target = (AliESDtrackCuts &) c;
427 target.fCutMinNClusterTPC = fCutMinNClusterTPC;
428 target.fCutMinNClusterITS = fCutMinNClusterITS;
430 target.fCutMaxChi2PerClusterTPC = fCutMaxChi2PerClusterTPC;
431 target.fCutMaxChi2PerClusterITS = fCutMaxChi2PerClusterITS;
433 for (Int_t i = 0; i < 3; i++)
434 target.fCutClusterRequirementITS[i] = fCutClusterRequirementITS[i];
436 target.fCutMaxC11 = fCutMaxC11;
437 target.fCutMaxC22 = fCutMaxC22;
438 target.fCutMaxC33 = fCutMaxC33;
439 target.fCutMaxC44 = fCutMaxC44;
440 target.fCutMaxC55 = fCutMaxC55;
442 target.fCutMaxRel1PtUncertainty = fCutMaxRel1PtUncertainty;
444 target.fCutAcceptKinkDaughters = fCutAcceptKinkDaughters;
445 target.fCutAcceptSharedTPCClusters = fCutAcceptSharedTPCClusters;
446 target.fCutMaxFractionSharedTPCClusters = fCutMaxFractionSharedTPCClusters;
447 target.fCutRequireTPCRefit = fCutRequireTPCRefit;
448 target.fCutRequireTPCStandAlone = fCutRequireTPCStandAlone;
449 target.fCutRequireITSRefit = fCutRequireITSRefit;
450 target.fCutRequireITSStandAlone = fCutRequireITSStandAlone;
451 target.fCutRejectITSpureSA = fCutRejectITSpureSA;
453 target.fCutNsigmaToVertex = fCutNsigmaToVertex;
454 target.fCutSigmaToVertexRequired = fCutSigmaToVertexRequired;
455 target.fCutMaxDCAToVertexXY = fCutMaxDCAToVertexXY;
456 target.fCutMaxDCAToVertexZ = fCutMaxDCAToVertexZ;
457 target.fCutDCAToVertex2D = fCutDCAToVertex2D;
458 target.fCutMinDCAToVertexXY = fCutMinDCAToVertexXY;
459 target.fCutMinDCAToVertexZ = fCutMinDCAToVertexZ;
461 target.fCutMaxDCAToVertexXYPtDep = fCutMaxDCAToVertexXYPtDep;
462 target.SetMaxDCAToVertexXYPtDep(fCutMaxDCAToVertexXYPtDep.Data());
464 target.fCutMaxDCAToVertexZPtDep = fCutMaxDCAToVertexZPtDep;
465 target.SetMaxDCAToVertexZPtDep(fCutMaxDCAToVertexZPtDep.Data());
467 target.fCutMinDCAToVertexXYPtDep = fCutMinDCAToVertexXYPtDep;
468 target.SetMinDCAToVertexXYPtDep(fCutMinDCAToVertexXYPtDep.Data());
470 target.fCutMinDCAToVertexZPtDep = fCutMinDCAToVertexZPtDep;
471 target.SetMinDCAToVertexZPtDep(fCutMinDCAToVertexZPtDep.Data());
473 target.fPMin = fPMin;
474 target.fPMax = fPMax;
475 target.fPtMin = fPtMin;
476 target.fPtMax = fPtMax;
477 target.fPxMin = fPxMin;
478 target.fPxMax = fPxMax;
479 target.fPyMin = fPyMin;
480 target.fPyMax = fPyMax;
481 target.fPzMin = fPzMin;
482 target.fPzMax = fPzMax;
483 target.fEtaMin = fEtaMin;
484 target.fEtaMax = fEtaMax;
485 target.fRapMin = fRapMin;
486 target.fRapMax = fRapMax;
488 target.fHistogramsOn = fHistogramsOn;
490 for (Int_t i=0; i<2; ++i)
492 if (fhNClustersITS[i]) target.fhNClustersITS[i] = (TH1F*) fhNClustersITS[i]->Clone();
493 if (fhNClustersTPC[i]) target.fhNClustersTPC[i] = (TH1F*) fhNClustersTPC[i]->Clone();
495 if (fhChi2PerClusterITS[i]) target.fhChi2PerClusterITS[i] = (TH1F*) fhChi2PerClusterITS[i]->Clone();
496 if (fhChi2PerClusterTPC[i]) target.fhChi2PerClusterTPC[i] = (TH1F*) fhChi2PerClusterTPC[i]->Clone();
498 if (fhC11[i]) target.fhC11[i] = (TH1F*) fhC11[i]->Clone();
499 if (fhC22[i]) target.fhC22[i] = (TH1F*) fhC22[i]->Clone();
500 if (fhC33[i]) target.fhC33[i] = (TH1F*) fhC33[i]->Clone();
501 if (fhC44[i]) target.fhC44[i] = (TH1F*) fhC44[i]->Clone();
502 if (fhC55[i]) target.fhC55[i] = (TH1F*) fhC55[i]->Clone();
504 if (fhRel1PtUncertainty[i]) target.fhRel1PtUncertainty[i] = (TH1F*) fhRel1PtUncertainty[i]->Clone();
506 if (fhDXY[i]) target.fhDXY[i] = (TH1F*) fhDXY[i]->Clone();
507 if (fhDZ[i]) target.fhDZ[i] = (TH1F*) fhDZ[i]->Clone();
508 if (fhDXYDZ[i]) target.fhDXYDZ[i] = (TH1F*) fhDXYDZ[i]->Clone();
509 if (fhDXYvsDZ[i]) target.fhDXYvsDZ[i] = (TH2F*) fhDXYvsDZ[i]->Clone();
511 if (fhDXYNormalized[i]) target.fhDXYNormalized[i] = (TH1F*) fhDXYNormalized[i]->Clone();
512 if (fhDZNormalized[i]) target.fhDZNormalized[i] = (TH1F*) fhDZNormalized[i]->Clone();
513 if (fhDXYvsDZNormalized[i]) target.fhDXYvsDZNormalized[i] = (TH2F*) fhDXYvsDZNormalized[i]->Clone();
514 if (fhNSigmaToVertex[i]) target.fhNSigmaToVertex[i] = (TH1F*) fhNSigmaToVertex[i]->Clone();
516 if (fhPt[i]) target.fhPt[i] = (TH1F*) fhPt[i]->Clone();
517 if (fhEta[i]) target.fhEta[i] = (TH1F*) fhEta[i]->Clone();
519 if (ffDTheoretical) target.ffDTheoretical = (TF1*) ffDTheoretical->Clone();
521 if (fhCutStatistics) target.fhCutStatistics = (TH1F*) fhCutStatistics->Clone();
522 if (fhCutCorrelation) target.fhCutCorrelation = (TH2F*) fhCutCorrelation->Clone();
527 //_____________________________________________________________________________
528 Long64_t AliESDtrackCuts::Merge(TCollection* list) {
529 // Merge a list of AliESDtrackCuts objects with this (needed for PROOF)
530 // Returns the number of merged objects (including this)
537 TIterator* iter = list->MakeIterator();
540 // collection of measured and generated histograms
542 while ((obj = iter->Next())) {
544 AliESDtrackCuts* entry = dynamic_cast<AliESDtrackCuts*>(obj);
548 if (!entry->fHistogramsOn)
551 for (Int_t i=0; i<2; i++) {
553 fhNClustersITS[i] ->Add(entry->fhNClustersITS[i] );
554 fhNClustersTPC[i] ->Add(entry->fhNClustersTPC[i] );
556 fhChi2PerClusterITS[i] ->Add(entry->fhChi2PerClusterITS[i]);
557 fhChi2PerClusterTPC[i] ->Add(entry->fhChi2PerClusterTPC[i]);
559 fhC11[i] ->Add(entry->fhC11[i] );
560 fhC22[i] ->Add(entry->fhC22[i] );
561 fhC33[i] ->Add(entry->fhC33[i] );
562 fhC44[i] ->Add(entry->fhC44[i] );
563 fhC55[i] ->Add(entry->fhC55[i] );
565 fhRel1PtUncertainty[i] ->Add(entry->fhRel1PtUncertainty[i]);
567 fhDXY[i] ->Add(entry->fhDXY[i] );
568 fhDZ[i] ->Add(entry->fhDZ[i] );
569 fhDXYDZ[i] ->Add(entry->fhDXYDZ[i] );
570 fhDXYvsDZ[i] ->Add(entry->fhDXYvsDZ[i] );
572 fhDXYNormalized[i] ->Add(entry->fhDXYNormalized[i] );
573 fhDZNormalized[i] ->Add(entry->fhDZNormalized[i] );
574 fhDXYvsDZNormalized[i] ->Add(entry->fhDXYvsDZNormalized[i]);
575 fhNSigmaToVertex[i] ->Add(entry->fhNSigmaToVertex[i]);
577 fhPt[i] ->Add(entry->fhPt[i]);
578 fhEta[i] ->Add(entry->fhEta[i]);
581 fhCutStatistics ->Add(entry->fhCutStatistics);
582 fhCutCorrelation ->Add(entry->fhCutCorrelation);
589 //____________________________________________________________________
590 AliESDtrackCuts* AliESDtrackCuts::GetStandardTPCOnlyTrackCuts()
592 // creates an AliESDtrackCuts object and fills it with standard (pre data-taking) values for TPC-only cuts
594 Printf("AliESDtrackCuts::GetStandardTPCOnlyTrackCuts: Creating track cuts for TPC-only.");
596 AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
598 esdTrackCuts->SetMinNClustersTPC(50);
599 esdTrackCuts->SetMaxChi2PerClusterTPC(4);
600 esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
602 esdTrackCuts->SetMaxDCAToVertexZ(3.2);
603 esdTrackCuts->SetMaxDCAToVertexXY(2.4);
604 esdTrackCuts->SetDCAToVertex2D(kTRUE);
609 //____________________________________________________________________
610 AliESDtrackCuts* AliESDtrackCuts::GetStandardITSTPCTrackCuts2009(Bool_t selPrimaries)
612 // creates an AliESDtrackCuts object and fills it with standard values for ITS-TPC cuts for pp 2009 data
614 Printf("AliESDtrackCuts::GetStandardITSTPCTrackCuts: Creating track cuts for ITS+TPC.");
616 AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
619 esdTrackCuts->SetRequireTPCStandAlone(kTRUE); // to get chi2 and ncls of kTPCin
620 esdTrackCuts->SetMinNClustersTPC(70);
621 esdTrackCuts->SetMaxChi2PerClusterTPC(4);
622 esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
623 esdTrackCuts->SetRequireTPCRefit(kTRUE);
625 esdTrackCuts->SetRequireITSRefit(kTRUE);
626 esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
627 AliESDtrackCuts::kAny);
629 // 7*(0.0050+0.0060/pt^0.9)
630 esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0350+0.0420/pt^0.9");
632 esdTrackCuts->SetMaxDCAToVertexZ(1.e6);
633 esdTrackCuts->SetDCAToVertex2D(kFALSE);
634 esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
635 //esdTrackCuts->SetEtaRange(-0.8,+0.8);
640 //____________________________________________________________________
641 Int_t AliESDtrackCuts::GetReferenceMultiplicity(AliESDEvent* esd, Bool_t tpcOnly)
643 // Gets reference multiplicity following the standard cuts and a defined fiducial volume
644 // tpcOnly = kTRUE -> consider TPC-only tracks
645 // = kFALSE -> consider global tracks
649 Printf("AliESDtrackCuts::GetReferenceMultiplicity: Not implemented for global tracks!");
653 AliESDtrackCuts* esdTrackCuts = GetStandardTPCOnlyTrackCuts();
654 esdTrackCuts->SetEtaRange(-0.8, 0.8);
655 esdTrackCuts->SetPtRange(0.15);
657 Int_t nTracks = esdTrackCuts->CountAcceptedTracks(esd);
665 //____________________________________________________________________
666 Float_t AliESDtrackCuts::GetSigmaToVertex(AliESDtrack* esdTrack)
668 // Calculates the number of sigma to the vertex.
673 esdTrack->GetImpactParameters(b,bCov);
675 if (bCov[0]<=0 || bCov[2]<=0) {
676 AliDebugClass(1, "Estimated b resolution lower or equal zero!");
677 bCov[0]=0; bCov[2]=0;
679 bRes[0] = TMath::Sqrt(bCov[0]);
680 bRes[1] = TMath::Sqrt(bCov[2]);
682 // -----------------------------------
683 // How to get to a n-sigma cut?
685 // The accumulated statistics from 0 to d is
687 // -> Erf(d/Sqrt(2)) for a 1-dim gauss (d = n_sigma)
688 // -> 1 - Exp(-d**2) for a 2-dim gauss (d*d = dx*dx + dy*dy != n_sigma)
690 // It means that for a 2-dim gauss: n_sigma(d) = Sqrt(2)*ErfInv(1 - Exp((-d**2)/2)
691 // Can this be expressed in a different way?
693 if (bRes[0] == 0 || bRes[1] ==0)
696 Float_t d = TMath::Sqrt(TMath::Power(b[0]/bRes[0],2) + TMath::Power(b[1]/bRes[1],2));
698 // work around precision problem
699 // if d is too big, TMath::Exp(...) gets 0, and TMath::ErfInverse(1) that should be infinite, gets 0 :(
700 // 1e-15 corresponds to nsigma ~ 7.7
701 if (TMath::Exp(-d * d / 2) < 1e-15)
704 Float_t nSigma = TMath::ErfInverse(1 - TMath::Exp(-d * d / 2)) * TMath::Sqrt(2);
708 void AliESDtrackCuts::EnableNeededBranches(TTree* tree)
710 // enables the branches needed by AcceptTrack, for a list see comment of AcceptTrack
712 tree->SetBranchStatus("fTracks.fFlags", 1);
713 tree->SetBranchStatus("fTracks.fITSncls", 1);
714 tree->SetBranchStatus("fTracks.fTPCncls", 1);
715 tree->SetBranchStatus("fTracks.fITSchi2", 1);
716 tree->SetBranchStatus("fTracks.fTPCchi2", 1);
717 tree->SetBranchStatus("fTracks.fC*", 1);
718 tree->SetBranchStatus("fTracks.fD", 1);
719 tree->SetBranchStatus("fTracks.fZ", 1);
720 tree->SetBranchStatus("fTracks.fCdd", 1);
721 tree->SetBranchStatus("fTracks.fCdz", 1);
722 tree->SetBranchStatus("fTracks.fCzz", 1);
723 tree->SetBranchStatus("fTracks.fP*", 1);
724 tree->SetBranchStatus("fTracks.fR*", 1);
725 tree->SetBranchStatus("fTracks.fKinkIndexes*", 1);
728 //____________________________________________________________________
729 Bool_t AliESDtrackCuts::AcceptTrack(AliESDtrack* esdTrack)
732 // figure out if the tracks survives all the track cuts defined
734 // the different quality parameter and kinematic values are first
735 // retrieved from the track. then it is found out what cuts the
736 // track did not survive and finally the cuts are imposed.
738 // this function needs the following branches:
744 // fTracks.fC //GetExternalCovariance
745 // fTracks.fD //GetImpactParameters
746 // fTracks.fZ //GetImpactParameters
747 // fTracks.fCdd //GetImpactParameters
748 // fTracks.fCdz //GetImpactParameters
749 // fTracks.fCzz //GetImpactParameters
750 // fTracks.fP //GetPxPyPz
751 // fTracks.fR //GetMass
752 // fTracks.fP //GetMass
753 // fTracks.fKinkIndexes
755 UInt_t status = esdTrack->GetStatus();
757 // getting quality parameters from the ESD track
758 Int_t nClustersITS = esdTrack->GetITSclusters(0);
759 Int_t nClustersTPC = -1;
760 if(fCutRequireTPCStandAlone) {
761 nClustersTPC = esdTrack->GetTPCNclsIter1();
764 nClustersTPC = esdTrack->GetTPCclusters(0);
767 Int_t nClustersTPCShared = esdTrack->GetTPCnclsS();
768 Float_t fracClustersTPCShared = -1.;
770 Float_t chi2PerClusterITS = -1;
771 Float_t chi2PerClusterTPC = -1;
773 chi2PerClusterITS = esdTrack->GetITSchi2()/Float_t(nClustersITS);
774 if (nClustersTPC!=0) {
775 if(fCutRequireTPCStandAlone) {
776 chi2PerClusterTPC = esdTrack->GetTPCchi2Iter1()/Float_t(nClustersTPC);
778 chi2PerClusterTPC = esdTrack->GetTPCchi2()/Float_t(nClustersTPC);
780 fracClustersTPCShared = Float_t(nClustersTPCShared)/Float_t(nClustersTPC);
784 esdTrack->GetExternalCovariance(extCov);
786 // getting the track to vertex parameters
787 Float_t nSigmaToVertex = GetSigmaToVertex(esdTrack);
791 esdTrack->GetImpactParameters(b,bCov);
792 if (bCov[0]<=0 || bCov[2]<=0) {
793 AliDebug(1, "Estimated b resolution lower or equal zero!");
794 bCov[0]=0; bCov[2]=0;
798 // set pt-dependent DCA cuts, if requested
799 SetPtDepDCACuts(esdTrack->Pt());
802 Float_t dcaToVertexXY = b[0];
803 Float_t dcaToVertexZ = b[1];
805 Float_t dcaToVertex = -1;
807 if (fCutDCAToVertex2D)
809 dcaToVertex = TMath::Sqrt(dcaToVertexXY*dcaToVertexXY/fCutMaxDCAToVertexXY/fCutMaxDCAToVertexXY + dcaToVertexZ*dcaToVertexZ/fCutMaxDCAToVertexZ/fCutMaxDCAToVertexZ);
812 dcaToVertex = TMath::Sqrt(dcaToVertexXY*dcaToVertexXY + dcaToVertexZ*dcaToVertexZ);
814 // getting the kinematic variables of the track
815 // (assuming the mass is known)
817 esdTrack->GetPxPyPz(p);
819 Float_t momentum = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2) + TMath::Power(p[2],2));
820 Float_t pt = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2));
821 Float_t energy = TMath::Sqrt(TMath::Power(esdTrack->GetMass(),2) + TMath::Power(momentum,2));
823 //y-eta related calculations
826 if((momentum != TMath::Abs(p[2]))&&(momentum != 0))
827 eta = 0.5*TMath::Log((momentum + p[2])/(momentum - p[2]));
828 if((energy != TMath::Abs(p[2]))&&(momentum != 0))
829 y = 0.5*TMath::Log((energy + p[2])/(energy - p[2]));
833 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]);
836 Float_t relUncertainty1Pt = TMath::Sqrt(extCov[14])*pt;
838 //########################################################################
842 for (Int_t i=0; i<kNCuts; i++) cuts[i]=kFALSE;
844 // track quality cuts
845 if (fCutRequireTPCRefit && (status&AliESDtrack::kTPCrefit)==0)
847 if (fCutRequireTPCStandAlone && (status&AliESDtrack::kTPCin)==0)
849 if (fCutRequireITSRefit && (status&AliESDtrack::kITSrefit)==0)
851 if (nClustersTPC<fCutMinNClusterTPC)
853 if (nClustersITS<fCutMinNClusterITS)
855 if (chi2PerClusterTPC>fCutMaxChi2PerClusterTPC)
857 if (chi2PerClusterITS>fCutMaxChi2PerClusterITS)
859 if (extCov[0] > fCutMaxC11)
861 if (extCov[2] > fCutMaxC22)
863 if (extCov[5] > fCutMaxC33)
865 if (extCov[9] > fCutMaxC44)
867 if (extCov[14] > fCutMaxC55)
869 if (nSigmaToVertex > fCutNsigmaToVertex && fCutSigmaToVertexRequired)
871 // if n sigma could not be calculated
872 if (nSigmaToVertex<0 && fCutSigmaToVertexRequired)
874 if (!fCutAcceptKinkDaughters && esdTrack->GetKinkIndex(0)>0)
876 // track kinematics cut
877 if((momentum < fPMin) || (momentum > fPMax))
879 if((pt < fPtMin) || (pt > fPtMax))
881 if((p[0] < fPxMin) || (p[0] > fPxMax))
883 if((p[1] < fPyMin) || (p[1] > fPyMax))
885 if((p[2] < fPzMin) || (p[2] > fPzMax))
887 if((eta < fEtaMin) || (eta > fEtaMax))
889 if((y < fRapMin) || (y > fRapMax))
891 if (fCutDCAToVertex2D && dcaToVertex > 1)
893 if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexXY) > fCutMaxDCAToVertexXY)
895 if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexZ) > fCutMaxDCAToVertexZ)
897 if (fCutDCAToVertex2D && fCutMinDCAToVertexXY > 0 && fCutMinDCAToVertexZ > 0 && dcaToVertexXY*dcaToVertexXY/fCutMinDCAToVertexXY/fCutMinDCAToVertexXY + dcaToVertexZ*dcaToVertexZ/fCutMinDCAToVertexZ/fCutMinDCAToVertexZ < 1)
899 if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexXY) < fCutMinDCAToVertexXY)
901 if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexZ) < fCutMinDCAToVertexZ)
904 for (Int_t i = 0; i < 3; i++)
905 cuts[28+i] = !CheckITSClusterRequirement(fCutClusterRequirementITS[i], esdTrack->HasPointOnITSLayer(i*2), esdTrack->HasPointOnITSLayer(i*2+1));
907 if (fCutRequireITSStandAlone && ((status & AliESDtrack::kITSin) == 0 || (status & AliESDtrack::kTPCin)))
910 if (fCutRejectITSpureSA && (status & AliESDtrack::kITSpureSA))
914 if (relUncertainty1Pt > fCutMaxRel1PtUncertainty)
917 if (!fCutAcceptSharedTPCClusters && nClustersTPCShared!=0)
920 if (fracClustersTPCShared > fCutMaxFractionSharedTPCClusters)
924 for (Int_t i=0; i<kNCuts; i++)
925 if (cuts[i]) {cut = kTRUE;}
928 //########################################################################
929 // filling histograms
931 fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n tracks")));
933 fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n cut tracks")));
935 for (Int_t i=0; i<kNCuts; i++) {
937 fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin(fgkCutNames[i])));
939 for (Int_t j=i; j<kNCuts; j++) {
940 if (cuts[i] && cuts[j]) {
941 Float_t xC = fhCutCorrelation->GetXaxis()->GetBinCenter(fhCutCorrelation->GetXaxis()->FindBin(fgkCutNames[i]));
942 Float_t yC = fhCutCorrelation->GetYaxis()->GetBinCenter(fhCutCorrelation->GetYaxis()->FindBin(fgkCutNames[j]));
943 fhCutCorrelation->Fill(xC, yC);
949 // now we loop over the filling of the histograms twice: once "before" the cut, once "after"
950 // the code is not in a function due to too many local variables that would need to be passed
952 for (Int_t id = 0; id < 2; id++)
954 // id = 0 --> before cut
955 // id = 1 --> after cut
959 fhNClustersITS[id]->Fill(nClustersITS);
960 fhNClustersTPC[id]->Fill(nClustersTPC);
961 fhChi2PerClusterITS[id]->Fill(chi2PerClusterITS);
962 fhChi2PerClusterTPC[id]->Fill(chi2PerClusterTPC);
964 fhC11[id]->Fill(extCov[0]);
965 fhC22[id]->Fill(extCov[2]);
966 fhC33[id]->Fill(extCov[5]);
967 fhC44[id]->Fill(extCov[9]);
968 fhC55[id]->Fill(extCov[14]);
970 fhRel1PtUncertainty[id]->Fill(relUncertainty1Pt);
973 fhEta[id]->Fill(eta);
976 bRes[0] = TMath::Sqrt(bCov[0]);
977 bRes[1] = TMath::Sqrt(bCov[2]);
979 fhDZ[id]->Fill(b[1]);
980 fhDXY[id]->Fill(b[0]);
981 fhDXYDZ[id]->Fill(dcaToVertex);
982 fhDXYvsDZ[id]->Fill(b[1],b[0]);
984 if (bRes[0]!=0 && bRes[1]!=0) {
985 fhDZNormalized[id]->Fill(b[1]/bRes[1]);
986 fhDXYNormalized[id]->Fill(b[0]/bRes[0]);
987 fhDXYvsDZNormalized[id]->Fill(b[1]/bRes[1], b[0]/bRes[0]);
988 fhNSigmaToVertex[id]->Fill(nSigmaToVertex);
1000 //____________________________________________________________________
1001 Bool_t AliESDtrackCuts::CheckITSClusterRequirement(ITSClusterRequirement req, Bool_t clusterL1, Bool_t clusterL2)
1003 // checks if the cluster requirement is fullfilled (in this case: return kTRUE)
1007 case kOff: return kTRUE;
1008 case kNone: return !clusterL1 && !clusterL2;
1009 case kAny: return clusterL1 || clusterL2;
1010 case kFirst: return clusterL1;
1011 case kOnlyFirst: return clusterL1 && !clusterL2;
1012 case kSecond: return clusterL2;
1013 case kOnlySecond: return clusterL2 && !clusterL1;
1014 case kBoth: return clusterL1 && clusterL2;
1020 //____________________________________________________________________
1021 AliESDtrack* AliESDtrackCuts::GetTPCOnlyTrack(AliESDEvent* esd, Int_t iTrack)
1024 // Utility function to
1025 // create a TPC only track from the given esd track
1027 // IMPORTANT: The track has to be deleted by the user
1029 // NB. most of the functionality to get a TPC only track from an ESD track is in AliESDtrack, where it should be
1030 // there are only missing propagations here that are needed for old data
1031 // this function will therefore become obsolete
1033 // adapted from code provided by CKB
1035 if (!esd->GetPrimaryVertexTPC())
1036 return 0; // No TPC vertex no TPC tracks
1038 if(!esd->GetPrimaryVertexTPC()->GetStatus())
1039 return 0; // TPC Vertex is created by default in AliESDEvent, do not use in this case
1041 AliESDtrack* track = esd->GetTrack(iTrack);
1045 AliESDtrack *tpcTrack = new AliESDtrack();
1047 // only true if we have a tpc track
1048 if (!track->FillTPCOnlyTrack(*tpcTrack))
1054 // propagate to Vertex
1055 // not needed for normal reconstructed ESDs...
1056 // Double_t pTPC[2],covTPC[3];
1057 // tpcTrack->PropagateToDCA(esd->GetPrimaryVertexTPC(), esd->GetMagneticField(), 10000, pTPC, covTPC);
1062 //____________________________________________________________________
1063 TObjArray* AliESDtrackCuts::GetAcceptedTracks(AliESDEvent* esd,Bool_t bTPC)
1066 // returns an array of all tracks that pass the cuts
1067 // or an array of TPC only tracks (propagated to the TPC vertex during reco)
1068 // tracks that pass the cut
1070 TObjArray* acceptedTracks = new TObjArray();
1072 // loop over esd tracks
1073 for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
1075 if(!esd->GetPrimaryVertexTPC())return acceptedTracks; // No TPC vertex no TPC tracks
1076 if(!esd->GetPrimaryVertexTPC()->GetStatus())return acceptedTracks; // No proper TPC vertex, only the default
1078 AliESDtrack *tpcTrack = GetTPCOnlyTrack(esd, iTrack);
1082 if (AcceptTrack(tpcTrack)) {
1083 acceptedTracks->Add(tpcTrack);
1090 AliESDtrack* track = esd->GetTrack(iTrack);
1091 if(AcceptTrack(track))
1092 acceptedTracks->Add(track);
1095 if(bTPC)acceptedTracks->SetOwner(kTRUE);
1096 return acceptedTracks;
1099 //____________________________________________________________________
1100 Int_t AliESDtrackCuts::CountAcceptedTracks(AliESDEvent* esd)
1103 // returns an the number of tracks that pass the cuts
1108 // loop over esd tracks
1109 for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
1110 AliESDtrack* track = esd->GetTrack(iTrack);
1111 if (AcceptTrack(track))
1118 //____________________________________________________________________
1119 void AliESDtrackCuts::DefineHistograms(Int_t color) {
1121 // diagnostics histograms are defined
1124 fHistogramsOn=kTRUE;
1126 Bool_t oldStatus = TH1::AddDirectoryStatus();
1127 TH1::AddDirectory(kFALSE);
1129 //###################################################################################
1130 // defining histograms
1132 fhCutStatistics = new TH1F("cut_statistics","cut statistics",kNCuts+4,-0.5,kNCuts+3.5);
1134 fhCutStatistics->GetXaxis()->SetBinLabel(1,"n tracks");
1135 fhCutStatistics->GetXaxis()->SetBinLabel(2,"n cut tracks");
1137 fhCutCorrelation = new TH2F("cut_correlation","cut correlation",kNCuts,-0.5,kNCuts-0.5,kNCuts,-0.5,kNCuts-0.5);;
1139 for (Int_t i=0; i<kNCuts; i++) {
1140 fhCutStatistics->GetXaxis()->SetBinLabel(i+4,fgkCutNames[i]);
1141 fhCutCorrelation->GetXaxis()->SetBinLabel(i+1,fgkCutNames[i]);
1142 fhCutCorrelation->GetYaxis()->SetBinLabel(i+1,fgkCutNames[i]);
1145 fhCutStatistics ->SetLineColor(color);
1146 fhCutCorrelation ->SetLineColor(color);
1147 fhCutStatistics ->SetLineWidth(2);
1148 fhCutCorrelation ->SetLineWidth(2);
1150 for (Int_t i=0; i<2; i++) {
1151 fhNClustersITS[i] = new TH1F("nClustersITS" ,"",8,-0.5,7.5);
1152 fhNClustersTPC[i] = new TH1F("nClustersTPC" ,"",165,-0.5,164.5);
1153 fhChi2PerClusterITS[i] = new TH1F("chi2PerClusterITS","",500,0,10);
1154 fhChi2PerClusterTPC[i] = new TH1F("chi2PerClusterTPC","",500,0,10);
1156 fhC11[i] = new TH1F("covMatrixDiagonal11","",2000,0,20);
1157 fhC22[i] = new TH1F("covMatrixDiagonal22","",2000,0,20);
1158 fhC33[i] = new TH1F("covMatrixDiagonal33","",1000,0,0.1);
1159 fhC44[i] = new TH1F("covMatrixDiagonal44","",1000,0,0.1);
1160 fhC55[i] = new TH1F("covMatrixDiagonal55","",1000,0,5);
1162 fhRel1PtUncertainty[i] = new TH1F("rel1PtUncertainty","",1000,0,5);
1164 fhDXY[i] = new TH1F("dXY" ,"",500,-10,10);
1165 fhDZ[i] = new TH1F("dZ" ,"",500,-10,10);
1166 fhDXYDZ[i] = new TH1F("dXYDZ" ,"",500,0,10);
1167 fhDXYvsDZ[i] = new TH2F("dXYvsDZ","",200,-10,10,200,-10,10);
1169 fhDXYNormalized[i] = new TH1F("dXYNormalized" ,"",500,-10,10);
1170 fhDZNormalized[i] = new TH1F("dZNormalized" ,"",500,-10,10);
1171 fhDXYvsDZNormalized[i] = new TH2F("dXYvsDZNormalized","",200,-10,10,200,-10,10);
1173 fhNSigmaToVertex[i] = new TH1F("nSigmaToVertex","",500,0,10);
1175 fhPt[i] = new TH1F("pt" ,"p_{T} distribution;p_{T} (GeV/c)", 800, 0.0, 10.0);
1176 fhEta[i] = new TH1F("eta" ,"#eta distribution;#eta",40,-2.0,2.0);
1178 fhNClustersITS[i]->SetTitle("n ITS clusters");
1179 fhNClustersTPC[i]->SetTitle("n TPC clusters");
1180 fhChi2PerClusterITS[i]->SetTitle("#Chi^{2} per ITS cluster");
1181 fhChi2PerClusterTPC[i]->SetTitle("#Chi^{2} per TPC cluster");
1183 fhC11[i]->SetTitle("cov 11 : #sigma_{y}^{2} [cm^{2}]");
1184 fhC22[i]->SetTitle("cov 22 : #sigma_{z}^{2} [cm^{2}]");
1185 fhC33[i]->SetTitle("cov 33 : #sigma_{sin(#phi)}^{2}");
1186 fhC44[i]->SetTitle("cov 44 : #sigma_{tan(#theta_{dip})}^{2}");
1187 fhC55[i]->SetTitle("cov 55 : #sigma_{1/p_{T}}^{2} [(c/GeV)^2]");
1189 fhRel1PtUncertainty[i]->SetTitle("rel. uncertainty of 1/p_{T}");
1191 fhDXY[i]->SetXTitle("transverse impact parameter (cm)");
1192 fhDZ[i]->SetXTitle("longitudinal impact parameter (cm)");
1193 fhDXYDZ[i]->SetTitle("absolute impact parameter;sqrt(dXY**2 + dZ**2) (cm)");
1194 fhDXYvsDZ[i]->SetXTitle("longitudinal impact parameter (cm)");
1195 fhDXYvsDZ[i]->SetYTitle("transverse impact parameter (cm)");
1197 fhDXYNormalized[i]->SetTitle("normalized trans impact par (n#sigma)");
1198 fhDZNormalized[i]->SetTitle("normalized long impact par (n#sigma)");
1199 fhDXYvsDZNormalized[i]->SetTitle("normalized long impact par (n#sigma)");
1200 fhDXYvsDZNormalized[i]->SetYTitle("normalized trans impact par (n#sigma)");
1201 fhNSigmaToVertex[i]->SetTitle("n #sigma to vertex");
1203 fhNClustersITS[i]->SetLineColor(color); fhNClustersITS[i]->SetLineWidth(2);
1204 fhNClustersTPC[i]->SetLineColor(color); fhNClustersTPC[i]->SetLineWidth(2);
1205 fhChi2PerClusterITS[i]->SetLineColor(color); fhChi2PerClusterITS[i]->SetLineWidth(2);
1206 fhChi2PerClusterTPC[i]->SetLineColor(color); fhChi2PerClusterTPC[i]->SetLineWidth(2);
1208 fhC11[i]->SetLineColor(color); fhC11[i]->SetLineWidth(2);
1209 fhC22[i]->SetLineColor(color); fhC22[i]->SetLineWidth(2);
1210 fhC33[i]->SetLineColor(color); fhC33[i]->SetLineWidth(2);
1211 fhC44[i]->SetLineColor(color); fhC44[i]->SetLineWidth(2);
1212 fhC55[i]->SetLineColor(color); fhC55[i]->SetLineWidth(2);
1214 fhRel1PtUncertainty[i]->SetLineColor(color); fhRel1PtUncertainty[i]->SetLineWidth(2);
1216 fhDXY[i]->SetLineColor(color); fhDXY[i]->SetLineWidth(2);
1217 fhDZ[i]->SetLineColor(color); fhDZ[i]->SetLineWidth(2);
1218 fhDXYDZ[i]->SetLineColor(color); fhDXYDZ[i]->SetLineWidth(2);
1220 fhDXYNormalized[i]->SetLineColor(color); fhDXYNormalized[i]->SetLineWidth(2);
1221 fhDZNormalized[i]->SetLineColor(color); fhDZNormalized[i]->SetLineWidth(2);
1222 fhNSigmaToVertex[i]->SetLineColor(color); fhNSigmaToVertex[i]->SetLineWidth(2);
1225 // The number of sigmas to the vertex is per definition gaussian
1226 ffDTheoretical = new TF1("nSigmaToVertexTheoretical","([0]/2.506628274)*exp(-(x**2)/2)",0,50);
1227 ffDTheoretical->SetParameter(0,1);
1229 TH1::AddDirectory(oldStatus);
1232 //____________________________________________________________________
1233 Bool_t AliESDtrackCuts::LoadHistograms(const Char_t* dir)
1236 // loads the histograms from a file
1237 // if dir is empty a directory with the name of this object is taken (like in SaveHistogram)
1243 if (!gDirectory->cd(dir))
1246 ffDTheoretical = dynamic_cast<TF1*> (gDirectory->Get("nSigmaToVertexTheory"));
1248 fhCutStatistics = dynamic_cast<TH1F*> (gDirectory->Get("cut_statistics"));
1249 fhCutCorrelation = dynamic_cast<TH2F*> (gDirectory->Get("cut_correlation"));
1251 for (Int_t i=0; i<2; i++) {
1254 gDirectory->cd("before_cuts");
1257 gDirectory->cd("after_cuts");
1259 fhNClustersITS[i] = dynamic_cast<TH1F*> (gDirectory->Get("nClustersITS" ));
1260 fhNClustersTPC[i] = dynamic_cast<TH1F*> (gDirectory->Get("nClustersTPC" ));
1261 fhChi2PerClusterITS[i] = dynamic_cast<TH1F*> (gDirectory->Get("chi2PerClusterITS"));
1262 fhChi2PerClusterTPC[i] = dynamic_cast<TH1F*> (gDirectory->Get("chi2PerClusterTPC"));
1264 fhC11[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal11"));
1265 fhC22[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal22"));
1266 fhC33[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal33"));
1267 fhC44[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal44"));
1268 fhC55[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal55"));
1270 fhRel1PtUncertainty[i] = dynamic_cast<TH1F*> (gDirectory->Get("rel1PtUncertainty"));
1272 fhDXY[i] = dynamic_cast<TH1F*> (gDirectory->Get("dXY" ));
1273 fhDZ[i] = dynamic_cast<TH1F*> (gDirectory->Get("dZ" ));
1274 fhDXYDZ[i] = dynamic_cast<TH1F*> (gDirectory->Get("dXYDZ"));
1275 fhDXYvsDZ[i] = dynamic_cast<TH2F*> (gDirectory->Get("dXYvsDZ"));
1277 fhDXYNormalized[i] = dynamic_cast<TH1F*> (gDirectory->Get("dXYNormalized" ));
1278 fhDZNormalized[i] = dynamic_cast<TH1F*> (gDirectory->Get("dZNormalized" ));
1279 fhDXYvsDZNormalized[i] = dynamic_cast<TH2F*> (gDirectory->Get("dXYvsDZNormalized"));
1280 fhNSigmaToVertex[i] = dynamic_cast<TH1F*> (gDirectory->Get("nSigmaToVertex"));
1282 fhPt[i] = dynamic_cast<TH1F*> (gDirectory->Get("pt"));
1283 fhEta[i] = dynamic_cast<TH1F*> (gDirectory->Get("eta"));
1285 gDirectory->cd("../");
1288 gDirectory->cd("..");
1293 //____________________________________________________________________
1294 void AliESDtrackCuts::SaveHistograms(const Char_t* dir) {
1296 // saves the histograms in a directory (dir)
1299 if (!fHistogramsOn) {
1300 AliDebug(0, "Histograms not on - cannot save histograms!!!");
1307 gDirectory->mkdir(dir);
1308 gDirectory->cd(dir);
1310 gDirectory->mkdir("before_cuts");
1311 gDirectory->mkdir("after_cuts");
1313 // a factor of 2 is needed since n sigma is positive
1314 ffDTheoretical->SetParameter(0,2*fhNSigmaToVertex[0]->Integral("width"));
1315 ffDTheoretical->Write("nSigmaToVertexTheory");
1317 fhCutStatistics->Write();
1318 fhCutCorrelation->Write();
1320 for (Int_t i=0; i<2; i++) {
1322 gDirectory->cd("before_cuts");
1324 gDirectory->cd("after_cuts");
1326 fhNClustersITS[i] ->Write();
1327 fhNClustersTPC[i] ->Write();
1328 fhChi2PerClusterITS[i] ->Write();
1329 fhChi2PerClusterTPC[i] ->Write();
1337 fhRel1PtUncertainty[i] ->Write();
1341 fhDXYDZ[i] ->Write();
1342 fhDXYvsDZ[i] ->Write();
1344 fhDXYNormalized[i] ->Write();
1345 fhDZNormalized[i] ->Write();
1346 fhDXYvsDZNormalized[i] ->Write();
1347 fhNSigmaToVertex[i] ->Write();
1352 gDirectory->cd("../");
1355 gDirectory->cd("../");
1358 //____________________________________________________________________
1359 void AliESDtrackCuts::DrawHistograms()
1361 // draws some histograms
1363 TCanvas* canvas1 = new TCanvas(Form("%s_1", GetName()), "Track Quality Results1", 800, 800);
1364 canvas1->Divide(2, 2);
1367 fhNClustersTPC[0]->SetStats(kFALSE);
1368 fhNClustersTPC[0]->Draw();
1371 fhChi2PerClusterTPC[0]->SetStats(kFALSE);
1372 fhChi2PerClusterTPC[0]->Draw();
1375 fhNSigmaToVertex[0]->SetStats(kFALSE);
1376 fhNSigmaToVertex[0]->GetXaxis()->SetRangeUser(0, 10);
1377 fhNSigmaToVertex[0]->Draw();
1379 canvas1->SaveAs(Form("%s_%s.gif", GetName(), canvas1->GetName()));
1381 TCanvas* canvas2 = new TCanvas(Form("%s_2", GetName()), "Track Quality Results2", 1200, 800);
1382 canvas2->Divide(3, 2);
1385 fhC11[0]->SetStats(kFALSE);
1390 fhC22[0]->SetStats(kFALSE);
1395 fhC33[0]->SetStats(kFALSE);
1400 fhC44[0]->SetStats(kFALSE);
1405 fhC55[0]->SetStats(kFALSE);
1410 fhRel1PtUncertainty[0]->SetStats(kFALSE);
1412 fhRel1PtUncertainty[0]->Draw();
1414 canvas2->SaveAs(Form("%s_%s.gif", GetName(), canvas2->GetName()));
1416 TCanvas* canvas3 = new TCanvas(Form("%s_3", GetName()), "Track Quality Results3", 1200, 800);
1417 canvas3->Divide(3, 2);
1420 fhDXY[0]->SetStats(kFALSE);
1425 fhDZ[0]->SetStats(kFALSE);
1430 fhDXYvsDZ[0]->SetStats(kFALSE);
1432 gPad->SetRightMargin(0.15);
1433 fhDXYvsDZ[0]->Draw("COLZ");
1436 fhDXYNormalized[0]->SetStats(kFALSE);
1438 fhDXYNormalized[0]->Draw();
1441 fhDZNormalized[0]->SetStats(kFALSE);
1443 fhDZNormalized[0]->Draw();
1446 fhDXYvsDZNormalized[0]->SetStats(kFALSE);
1448 gPad->SetRightMargin(0.15);
1449 fhDXYvsDZNormalized[0]->Draw("COLZ");
1451 canvas3->SaveAs(Form("%s_%s.gif", GetName(), canvas3->GetName()));
1453 TCanvas* canvas4 = new TCanvas(Form("%s_4", GetName()), "Track Quality Results4", 800, 500);
1454 canvas4->Divide(2, 1);
1457 fhCutStatistics->SetStats(kFALSE);
1458 fhCutStatistics->LabelsOption("v");
1459 gPad->SetBottomMargin(0.3);
1460 fhCutStatistics->Draw();
1463 fhCutCorrelation->SetStats(kFALSE);
1464 fhCutCorrelation->LabelsOption("v");
1465 gPad->SetBottomMargin(0.3);
1466 gPad->SetLeftMargin(0.3);
1467 fhCutCorrelation->Draw("COLZ");
1469 canvas4->SaveAs(Form("%s_%s.gif", GetName(), canvas4->GetName()));
1472 fhDXYvsDZNormalized[0]->SetStats(kFALSE);
1473 fhDXYvsDZNormalized[0]->DrawCopy("COLZ");
1476 fhNClustersTPC[0]->SetStats(kFALSE);
1477 fhNClustersTPC[0]->DrawCopy();
1480 fhChi2PerClusterITS[0]->SetStats(kFALSE);
1481 fhChi2PerClusterITS[0]->DrawCopy();
1482 fhChi2PerClusterITS[1]->SetLineColor(2);
1483 fhChi2PerClusterITS[1]->DrawCopy("SAME");
1486 fhChi2PerClusterTPC[0]->SetStats(kFALSE);
1487 fhChi2PerClusterTPC[0]->DrawCopy();
1488 fhChi2PerClusterTPC[1]->SetLineColor(2);
1489 fhChi2PerClusterTPC[1]->DrawCopy("SAME");*/
1491 //--------------------------------------------------------------------------
1492 void AliESDtrackCuts::SetPtDepDCACuts(Double_t pt) {
1494 // set the pt-dependent DCA cuts
1497 if(f1CutMaxDCAToVertexXYPtDep) {
1498 fCutMaxDCAToVertexXY=f1CutMaxDCAToVertexXYPtDep->Eval(pt);
1501 if(f1CutMaxDCAToVertexZPtDep) {
1502 fCutMaxDCAToVertexZ=f1CutMaxDCAToVertexZPtDep->Eval(pt);
1505 if(f1CutMinDCAToVertexXYPtDep) {
1506 fCutMinDCAToVertexXY=f1CutMinDCAToVertexXYPtDep->Eval(pt);
1509 if(f1CutMinDCAToVertexZPtDep) {
1510 fCutMinDCAToVertexZ=f1CutMinDCAToVertexZPtDep->Eval(pt);
1519 //--------------------------------------------------------------------------
1520 Bool_t AliESDtrackCuts::CheckPtDepDCA(TString dist,Bool_t print) const {
1522 // Check the correctness of the string syntax
1524 Bool_t retval=kTRUE;
1526 if(!dist.Contains("pt")) {
1527 if(print) printf("AliESDtrackCuts::CheckPtDepDCA(): string must contain \"pt\"\n");
1533 void AliESDtrackCuts::SetMaxDCAToVertexXYPtDep(const char *dist){
1534 if(!CheckPtDepDCA(dist,kTRUE)) return;
1535 if(f1CutMaxDCAToVertexXYPtDep)delete f1CutMaxDCAToVertexXYPtDep;
1536 fCutMaxDCAToVertexXYPtDep = dist;
1538 tmp.ReplaceAll("pt","x");
1539 f1CutMaxDCAToVertexXYPtDep = new TFormula("f1CutMaxDCAToVertexXYPtDep",tmp.Data());
1543 void AliESDtrackCuts::SetMaxDCAToVertexZPtDep(const char *dist){
1544 if(!CheckPtDepDCA(dist,kTRUE)) return;
1545 if(f1CutMaxDCAToVertexZPtDep)delete f1CutMaxDCAToVertexZPtDep;
1546 fCutMaxDCAToVertexZPtDep = dist;
1548 tmp.ReplaceAll("pt","x");
1549 f1CutMaxDCAToVertexZPtDep = new TFormula("f1CutMaxDCAToVertexZPtDep",tmp.Data());
1555 void AliESDtrackCuts::SetMinDCAToVertexXYPtDep(const char *dist){
1556 if(!CheckPtDepDCA(dist,kTRUE)) return;
1557 if(f1CutMinDCAToVertexXYPtDep)delete f1CutMinDCAToVertexXYPtDep;
1558 fCutMinDCAToVertexXYPtDep = dist;
1560 tmp.ReplaceAll("pt","x");
1561 f1CutMinDCAToVertexXYPtDep = new TFormula("f1CutMinDCAToVertexXYPtDep",tmp.Data());
1566 void AliESDtrackCuts::SetMinDCAToVertexZPtDep(const char *dist){
1567 if(!CheckPtDepDCA(dist,kTRUE)) return;
1568 if(f1CutMinDCAToVertexZPtDep)delete f1CutMinDCAToVertexZPtDep;
1569 fCutMinDCAToVertexZPtDep = dist;
1571 tmp.ReplaceAll("pt","x");
1572 f1CutMinDCAToVertexZPtDep = new TFormula("f1CutMinDCAToVertexZPtDep",tmp.Data());