ANALYSIS/AliESDtrackCuts.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 /* $Id: AliESDtrackCuts.cxx 24534 2008-03-16 22:22:11Z fca $ */
  17
  18 #include "AliESDtrackCuts.h"
  19
  20 #include <AliESDtrack.h>
  21 #include <AliESDVertex.h>
  22 #include <AliESDEvent.h>
  23 #include <AliLog.h>
  24
  25 #include <TTree.h>
  26 #include <TCanvas.h>
  27 #include <TDirectory.h>
  28
  29 //____________________________________________________________________
  30 ClassImp(AliESDtrackCuts)
  31
  32 // Cut names
  33 const Char_t* AliESDtrackCuts::fgkCutNames[kNCuts] = {
  34  "require TPC refit",
  35  "require ITS refit",
  36  "n clusters TPC",
  37  "n clusters ITS",
  38  "#Chi^{2}/cluster TPC",
  39  "#Chi^{2}/cluster ITS",
  40  "cov 11",
  41  "cov 22",
  42  "cov 33",
  43  "cov 44",
  44  "cov 55",
  45  "trk-to-vtx",
  46  "trk-to-vtx failed",
  47  "kink daughters",
  48  "p",
  49  "p_{T}",
  50  "p_{x}",
  51  "p_{y}",
  52  "p_{z}",
  53  "eta",
  54  "y",
  55  "trk-to-vtx max dca 2D absolute",
  56  "trk-to-vtx max dca xy absolute",
  57  "trk-to-vtx max dca z absolute",
  58  "trk-to-vtx min dca 2D absolute",
  59  "trk-to-vtx min dca xy absolute",
  60  "trk-to-vtx min dca z absolute",
  61  "SPD cluster requirement",
  62  "SDD cluster requirement",
  63  "SSD cluster requirement"
  64 };
  65
  66 //____________________________________________________________________
  67 AliESDtrackCuts::AliESDtrackCuts(const Char_t* name, const Char_t* title) : AliAnalysisCuts(name,title),
  68   fCutMinNClusterTPC(0),
  69   fCutMinNClusterITS(0),
  70   fCutMaxChi2PerClusterTPC(0),
  71   fCutMaxChi2PerClusterITS(0),
  72   fCutMaxC11(0),
  73   fCutMaxC22(0),
  74   fCutMaxC33(0),
  75   fCutMaxC44(0),
  76   fCutMaxC55(0),
  77   fCutAcceptKinkDaughters(0),
  78   fCutRequireTPCRefit(0),
  79   fCutRequireITSRefit(0),
  80   fCutNsigmaToVertex(0),
  81   fCutSigmaToVertexRequired(0),
  82   fCutMaxDCAToVertexXY(0),
  83   fCutMaxDCAToVertexZ(0),
  84   fCutMinDCAToVertexXY(0),
  85   fCutMinDCAToVertexZ(0),
  86   fCutDCAToVertex2D(0),
  87   fPMin(0),
  88   fPMax(0),
  89   fPtMin(0),
  90   fPtMax(0),
  91   fPxMin(0),
  92   fPxMax(0),
  93   fPyMin(0),
  94   fPyMax(0),
  95   fPzMin(0),
  96   fPzMax(0),
  97   fEtaMin(0),
  98   fEtaMax(0),
  99   fRapMin(0),
 100   fRapMax(0),
 101   fHistogramsOn(0),
 102   ffDTheoretical(0),
 103   fhCutStatistics(0),
 104   fhCutCorrelation(0)
 105 {
 106   //
 107   // constructor
 108   //
 109
 110   Init();
 111
 112   //##############################################################################
 113   // setting default cuts
 114   SetMinNClustersTPC();
 115   SetMinNClustersITS();
 116   SetMaxChi2PerClusterTPC();
 117   SetMaxChi2PerClusterITS();
 118   SetMaxCovDiagonalElements();
 119   SetRequireTPCRefit();
 120   SetRequireITSRefit();
 121   SetAcceptKingDaughters();
 122   SetMaxNsigmaToVertex();
 123   SetMaxDCAToVertexXY();
 124   SetMaxDCAToVertexZ();
 125   SetDCAToVertex2D();
 126   SetMinDCAToVertexXY();
 127   SetMinDCAToVertexZ();
 128   SetPRange();
 129   SetPtRange();
 130   SetPxRange();
 131   SetPyRange();
 132   SetPzRange();
 133   SetEtaRange();
 134   SetRapRange();
 135   SetClusterRequirementITS(kSPD);
 136   SetClusterRequirementITS(kSDD);
 137   SetClusterRequirementITS(kSSD);
 138
 139   SetHistogramsOn();
 140 }
 141
 142 //_____________________________________________________________________________
 143 AliESDtrackCuts::AliESDtrackCuts(const AliESDtrackCuts &c) : AliAnalysisCuts(c),
 144   fCutMinNClusterTPC(0),
 145   fCutMinNClusterITS(0),
 146   fCutMaxChi2PerClusterTPC(0),
 147   fCutMaxChi2PerClusterITS(0),
 148   fCutMaxC11(0),
 149   fCutMaxC22(0),
 150   fCutMaxC33(0),
 151   fCutMaxC44(0),
 152   fCutMaxC55(0),
 153   fCutAcceptKinkDaughters(0),
 154   fCutRequireTPCRefit(0),
 155   fCutRequireITSRefit(0),
 156   fCutNsigmaToVertex(0),
 157   fCutSigmaToVertexRequired(0),
 158   fCutMaxDCAToVertexXY(0),
 159   fCutMaxDCAToVertexZ(0),
 160   fCutMinDCAToVertexXY(0),
 161   fCutMinDCAToVertexZ(0),
 162   fCutDCAToVertex2D(0),
 163   fPMin(0),
 164   fPMax(0),
 165   fPtMin(0),
 166   fPtMax(0),
 167   fPxMin(0),
 168   fPxMax(0),
 169   fPyMin(0),
 170   fPyMax(0),
 171   fPzMin(0),
 172   fPzMax(0),
 173   fEtaMin(0),
 174   fEtaMax(0),
 175   fRapMin(0),
 176   fRapMax(0),
 177   fHistogramsOn(0),
 178   ffDTheoretical(0),
 179   fhCutStatistics(0),
 180   fhCutCorrelation(0)
 181 {
 182   //
 183   // copy constructor
 184   //
 185
 186   ((AliESDtrackCuts &) c).Copy(*this);
 187 }
 188
 189 AliESDtrackCuts::~AliESDtrackCuts()
 190 {
 191   //
 192   // destructor
 193   //
 194
 195   for (Int_t i=0; i<2; i++) {
 196
 197     if (fhNClustersITS[i])
 198       delete fhNClustersITS[i];
 199     if (fhNClustersTPC[i])
 200       delete fhNClustersTPC[i];
 201     if (fhChi2PerClusterITS[i])
 202       delete fhChi2PerClusterITS[i];
 203     if (fhChi2PerClusterTPC[i])
 204       delete fhChi2PerClusterTPC[i];
 205     if (fhC11[i])
 206       delete fhC11[i];
 207     if (fhC22[i])
 208       delete fhC22[i];
 209     if (fhC33[i])
 210       delete fhC33[i];
 211     if (fhC44[i])
 212       delete fhC44[i];
 213     if (fhC55[i])
 214     delete fhC55[i];
 215
 216     if (fhDXY[i])
 217       delete fhDXY[i];
 218     if (fhDZ[i])
 219       delete fhDZ[i];
 220     if (fhDXYDZ[i])
 221       delete fhDXYDZ[i];
 222     if (fhDXYvsDZ[i])
 223       delete fhDXYvsDZ[i];
 224
 225     if (fhDXYNormalized[i])
 226       delete fhDXYNormalized[i];
 227     if (fhDZNormalized[i])
 228       delete fhDZNormalized[i];
 229     if (fhDXYvsDZNormalized[i])
 230       delete fhDXYvsDZNormalized[i];
 231     if (fhNSigmaToVertex[i])
 232       delete fhNSigmaToVertex[i];
 233     if (fhPt[i])
 234       delete fhPt[i];
 235     if (fhEta[i])
 236       delete fhEta[i];
 237   }
 238
 239   if (ffDTheoretical)
 240     delete ffDTheoretical;
 241
 242   if (fhCutStatistics)
 243     delete fhCutStatistics;
 244   if (fhCutCorrelation)
 245     delete fhCutCorrelation;
 246 }
 247
 248 void AliESDtrackCuts::Init()
 249 {
 250   //
 251   // sets everything to zero
 252   //
 253
 254   fCutMinNClusterTPC = 0;
 255   fCutMinNClusterITS = 0;
 256
 257   fCutMaxChi2PerClusterTPC = 0;
 258   fCutMaxChi2PerClusterITS = 0;
 259
 260   for (Int_t i = 0; i < 3; i++)
 261         fCutClusterRequirementITS[i] = kOff;
 262
 263   fCutMaxC11 = 0;
 264   fCutMaxC22 = 0;
 265   fCutMaxC33 = 0;
 266   fCutMaxC44 = 0;
 267   fCutMaxC55 = 0;
 268
 269   fCutAcceptKinkDaughters = 0;
 270   fCutRequireTPCRefit = 0;
 271   fCutRequireITSRefit = 0;
 272
 273   fCutNsigmaToVertex = 0;
 274   fCutSigmaToVertexRequired = 0;
 275   fCutMaxDCAToVertexXY = 0;
 276   fCutMaxDCAToVertexZ = 0;
 277   fCutDCAToVertex2D = 0;
 278   fCutMinDCAToVertexXY = 0;
 279   fCutMinDCAToVertexZ = 0;
 280
 281
 282   fPMin = 0;
 283   fPMax = 0;
 284   fPtMin = 0;
 285   fPtMax = 0;
 286   fPxMin = 0;
 287   fPxMax = 0;
 288   fPyMin = 0;
 289   fPyMax = 0;
 290   fPzMin = 0;
 291   fPzMax = 0;
 292   fEtaMin = 0;
 293   fEtaMax = 0;
 294   fRapMin = 0;
 295   fRapMax = 0;
 296
 297   fHistogramsOn = kFALSE;
 298
 299   for (Int_t i=0; i<2; ++i)
 300   {
 301     fhNClustersITS[i] = 0;
 302     fhNClustersTPC[i] = 0;
 303
 304     fhChi2PerClusterITS[i] = 0;
 305     fhChi2PerClusterTPC[i] = 0;
 306
 307     fhC11[i] = 0;
 308     fhC22[i] = 0;
 309     fhC33[i] = 0;
 310     fhC44[i] = 0;
 311     fhC55[i] = 0;
 312
 313     fhDXY[i] = 0;
 314     fhDZ[i] = 0;
 315     fhDXYDZ[i] = 0;
 316     fhDXYvsDZ[i] = 0;
 317
 318     fhDXYNormalized[i] = 0;
 319     fhDZNormalized[i] = 0;
 320     fhDXYvsDZNormalized[i] = 0;
 321     fhNSigmaToVertex[i] = 0;
 322
 323     fhPt[i] = 0;
 324     fhEta[i] = 0;
 325   }
 326   ffDTheoretical = 0;
 327
 328   fhCutStatistics = 0;
 329   fhCutCorrelation = 0;
 330 }
 331
 332 //_____________________________________________________________________________
 333 AliESDtrackCuts &AliESDtrackCuts::operator=(const AliESDtrackCuts &c)
 334 {
 335   //
 336   // Assignment operator
 337   //
 338
 339   if (this != &c) ((AliESDtrackCuts &) c).Copy(*this);
 340   return *this;
 341 }
 342
 343 //_____________________________________________________________________________
 344 void AliESDtrackCuts::Copy(TObject &c) const
 345 {
 346   //
 347   // Copy function
 348   //
 349
 350   AliESDtrackCuts& target = (AliESDtrackCuts &) c;
 351
 352   target.Init();
 353
 354   target.fCutMinNClusterTPC = fCutMinNClusterTPC;
 355   target.fCutMinNClusterITS = fCutMinNClusterITS;
 356
 357   target.fCutMaxChi2PerClusterTPC = fCutMaxChi2PerClusterTPC;
 358   target.fCutMaxChi2PerClusterITS = fCutMaxChi2PerClusterITS;
 359
 360   for (Int_t i = 0; i < 3; i++)
 361     target.fCutClusterRequirementITS[i] = fCutClusterRequirementITS[i];
 362
 363   target.fCutMaxC11 = fCutMaxC11;
 364   target.fCutMaxC22 = fCutMaxC22;
 365   target.fCutMaxC33 = fCutMaxC33;
 366   target.fCutMaxC44 = fCutMaxC44;
 367   target.fCutMaxC55 = fCutMaxC55;
 368
 369   target.fCutAcceptKinkDaughters = fCutAcceptKinkDaughters;
 370   target.fCutRequireTPCRefit = fCutRequireTPCRefit;
 371   target.fCutRequireITSRefit = fCutRequireITSRefit;
 372
 373   target.fCutNsigmaToVertex = fCutNsigmaToVertex;
 374   target.fCutSigmaToVertexRequired = fCutSigmaToVertexRequired;
 375   target.fCutMaxDCAToVertexXY = fCutMaxDCAToVertexXY;
 376   target.fCutMaxDCAToVertexZ = fCutMaxDCAToVertexZ;
 377   target.fCutDCAToVertex2D = fCutDCAToVertex2D;
 378   target.fCutMinDCAToVertexXY = fCutMinDCAToVertexXY;
 379   target.fCutMinDCAToVertexZ = fCutMinDCAToVertexZ;
 380
 381   target.fPMin = fPMin;
 382   target.fPMax = fPMax;
 383   target.fPtMin = fPtMin;
 384   target.fPtMax = fPtMax;
 385   target.fPxMin = fPxMin;
 386   target.fPxMax = fPxMax;
 387   target.fPyMin = fPyMin;
 388   target.fPyMax = fPyMax;
 389   target.fPzMin = fPzMin;
 390   target.fPzMax = fPzMax;
 391   target.fEtaMin = fEtaMin;
 392   target.fEtaMax = fEtaMax;
 393   target.fRapMin = fRapMin;
 394   target.fRapMax = fRapMax;
 395
 396   target.fHistogramsOn = fHistogramsOn;
 397
 398   for (Int_t i=0; i<2; ++i)
 399   {
 400     if (fhNClustersITS[i]) target.fhNClustersITS[i] = (TH1F*) fhNClustersITS[i]->Clone();
 401     if (fhNClustersTPC[i]) target.fhNClustersTPC[i] = (TH1F*) fhNClustersTPC[i]->Clone();
 402
 403     if (fhChi2PerClusterITS[i]) target.fhChi2PerClusterITS[i] = (TH1F*) fhChi2PerClusterITS[i]->Clone();
 404     if (fhChi2PerClusterTPC[i]) target.fhChi2PerClusterTPC[i] = (TH1F*) fhChi2PerClusterTPC[i]->Clone();
 405
 406     if (fhC11[i]) target.fhC11[i] = (TH1F*) fhC11[i]->Clone();
 407     if (fhC22[i]) target.fhC22[i] = (TH1F*) fhC22[i]->Clone();
 408     if (fhC33[i]) target.fhC33[i] = (TH1F*) fhC33[i]->Clone();
 409     if (fhC44[i]) target.fhC44[i] = (TH1F*) fhC44[i]->Clone();
 410     if (fhC55[i]) target.fhC55[i] = (TH1F*) fhC55[i]->Clone();
 411
 412     if (fhDXY[i]) target.fhDXY[i] = (TH1F*) fhDXY[i]->Clone();
 413     if (fhDZ[i]) target.fhDZ[i] = (TH1F*) fhDZ[i]->Clone();
 414     if (fhDXYDZ[i]) target.fhDXYDZ[i] = (TH1F*) fhDXYDZ[i]->Clone();
 415     if (fhDXYvsDZ[i]) target.fhDXYvsDZ[i] = (TH2F*) fhDXYvsDZ[i]->Clone();
 416
 417     if (fhDXYNormalized[i]) target.fhDXYNormalized[i] = (TH1F*) fhDXYNormalized[i]->Clone();
 418     if (fhDZNormalized[i]) target.fhDZNormalized[i] = (TH1F*) fhDZNormalized[i]->Clone();
 419     if (fhDXYvsDZNormalized[i]) target.fhDXYvsDZNormalized[i] = (TH2F*) fhDXYvsDZNormalized[i]->Clone();
 420     if (fhNSigmaToVertex[i]) target.fhNSigmaToVertex[i] = (TH1F*) fhNSigmaToVertex[i]->Clone();
 421
 422     if (fhPt[i]) target.fhPt[i] = (TH1F*) fhPt[i]->Clone();
 423     if (fhEta[i]) target.fhEta[i] = (TH1F*) fhEta[i]->Clone();
 424   }
 425   if (ffDTheoretical) target.ffDTheoretical = (TF1*) ffDTheoretical->Clone();
 426
 427   if (fhCutStatistics) target.fhCutStatistics = (TH1F*) fhCutStatistics->Clone();
 428   if (fhCutCorrelation) target.fhCutCorrelation = (TH2F*) fhCutCorrelation->Clone();
 429
 430   TNamed::Copy(c);
 431 }
 432
 433 //_____________________________________________________________________________
 434 Long64_t AliESDtrackCuts::Merge(TCollection* list) {
 435   // Merge a list of AliESDtrackCuts objects with this (needed for PROOF)
 436   // Returns the number of merged objects (including this)
 437   if (!list)
 438     return 0;
 439   if (list->IsEmpty())
 440     return 1;
 441   if (!fHistogramsOn)
 442     return 0;
 443   TIterator* iter = list->MakeIterator();
 444   TObject* obj;
 445
 446   // collection of measured and generated histograms
 447   Int_t count = 0;
 448   while ((obj = iter->Next())) {
 449
 450     AliESDtrackCuts* entry = dynamic_cast<AliESDtrackCuts*>(obj);
 451     if (entry == 0)
 452       continue;
 453
 454     if (!entry->fHistogramsOn)
 455       continue;
 456
 457     for (Int_t i=0; i<2; i++) {
 458
 459       fhNClustersITS[i]      ->Add(entry->fhNClustersITS[i]     );
 460       fhNClustersTPC[i]      ->Add(entry->fhNClustersTPC[i]     );
 461
 462       fhChi2PerClusterITS[i] ->Add(entry->fhChi2PerClusterITS[i]);
 463       fhChi2PerClusterTPC[i] ->Add(entry->fhChi2PerClusterTPC[i]);
 464
 465       fhC11[i]               ->Add(entry->fhC11[i]              );
 466       fhC22[i]               ->Add(entry->fhC22[i]              );
 467       fhC33[i]               ->Add(entry->fhC33[i]              );
 468       fhC44[i]               ->Add(entry->fhC44[i]              );
 469       fhC55[i]               ->Add(entry->fhC55[i]              );
 470
 471       fhDXY[i]               ->Add(entry->fhDXY[i]              );
 472       fhDZ[i]                ->Add(entry->fhDZ[i]               );
 473       fhDXYDZ[i]             ->Add(entry->fhDXYDZ[i]          );
 474       fhDXYvsDZ[i]           ->Add(entry->fhDXYvsDZ[i]          );
 475
 476       fhDXYNormalized[i]     ->Add(entry->fhDXYNormalized[i]    );
 477       fhDZNormalized[i]      ->Add(entry->fhDZNormalized[i]     );
 478       fhDXYvsDZNormalized[i] ->Add(entry->fhDXYvsDZNormalized[i]);
 479       fhNSigmaToVertex[i]    ->Add(entry->fhNSigmaToVertex[i]);
 480
 481       fhPt[i]                ->Add(entry->fhPt[i]);
 482       fhEta[i]               ->Add(entry->fhEta[i]);
 483     }
 484
 485     fhCutStatistics  ->Add(entry->fhCutStatistics);
 486     fhCutCorrelation ->Add(entry->fhCutCorrelation);
 487
 488     count++;
 489   }
 490   return count+1;
 491 }
 492
 493
 494 //____________________________________________________________________
 495 Float_t AliESDtrackCuts::GetSigmaToVertex(AliESDtrack* esdTrack)
 496 {
 497   // Calculates the number of sigma to the vertex.
 498
 499   Float_t b[2];
 500   Float_t bRes[2];
 501   Float_t bCov[3];
 502   esdTrack->GetImpactParameters(b,bCov);
 503
 504   if (bCov[0]<=0 || bCov[2]<=0) {
 505     AliDebugClass(1, "Estimated b resolution lower or equal zero!");
 506     bCov[0]=0; bCov[2]=0;
 507   }
 508   bRes[0] = TMath::Sqrt(bCov[0]);
 509   bRes[1] = TMath::Sqrt(bCov[2]);
 510
 511   // -----------------------------------
 512   // How to get to a n-sigma cut?
 513   //
 514   // The accumulated statistics from 0 to d is
 515   //
 516   // ->  Erf(d/Sqrt(2)) for a 1-dim gauss (d = n_sigma)
 517   // ->  1 - Exp(-d**2) for a 2-dim gauss (d*d = dx*dx + dy*dy != n_sigma)
 518   //
 519   // It means that for a 2-dim gauss: n_sigma(d) = Sqrt(2)*ErfInv(1 - Exp((-x**2)/2)
 520   // Can this be expressed in a different way?
 521
 522   if (bRes[0] == 0 || bRes[1] ==0)
 523     return -1;
 524
 525   Float_t d = TMath::Sqrt(TMath::Power(b[0]/bRes[0],2) + TMath::Power(b[1]/bRes[1],2));
 526
 527   // work around precision problem
 528   // if d is too big, TMath::Exp(...) gets 0, and TMath::ErfInverse(1) that should be infinite, gets 0 :(
 529   // 1e-15 corresponds to nsigma ~ 7.7
 530   if (TMath::Exp(-d * d / 2) < 1e-15)
 531     return 1000;
 532
 533   Float_t nSigma = TMath::ErfInverse(1 - TMath::Exp(-d * d / 2)) * TMath::Sqrt(2);
 534   return nSigma;
 535 }
 536
 537 void AliESDtrackCuts::EnableNeededBranches(TTree* tree)
 538 {
 539   // enables the branches needed by AcceptTrack, for a list see comment of AcceptTrack
 540
 541   tree->SetBranchStatus("fTracks.fFlags", 1);
 542   tree->SetBranchStatus("fTracks.fITSncls", 1);
 543   tree->SetBranchStatus("fTracks.fTPCncls", 1);
 544   tree->SetBranchStatus("fTracks.fITSchi2", 1);
 545   tree->SetBranchStatus("fTracks.fTPCchi2", 1);
 546   tree->SetBranchStatus("fTracks.fC*", 1);
 547   tree->SetBranchStatus("fTracks.fD", 1);
 548   tree->SetBranchStatus("fTracks.fZ", 1);
 549   tree->SetBranchStatus("fTracks.fCdd", 1);
 550   tree->SetBranchStatus("fTracks.fCdz", 1);
 551   tree->SetBranchStatus("fTracks.fCzz", 1);
 552   tree->SetBranchStatus("fTracks.fP*", 1);
 553   tree->SetBranchStatus("fTracks.fR*", 1);
 554   tree->SetBranchStatus("fTracks.fKinkIndexes*", 1);
 555 }
 556
 557 //____________________________________________________________________
 558 Bool_t
 559 AliESDtrackCuts::AcceptTrack(AliESDtrack* esdTrack) {
 560   //
 561   // figure out if the tracks survives all the track cuts defined
 562   //
 563   // the different quality parameter and kinematic values are first
 564   // retrieved from the track. then it is found out what cuts the
 565   // track did not survive and finally the cuts are imposed.
 566
 567   // this function needs the following branches:
 568   // fTracks.fFlags
 569   // fTracks.fITSncls
 570   // fTracks.fTPCncls
 571   // fTracks.fITSchi2
 572   // fTracks.fTPCchi2
 573   // fTracks.fC   //GetExternalCovariance
 574   // fTracks.fD   //GetImpactParameters
 575   // fTracks.fZ   //GetImpactParameters
 576   // fTracks.fCdd //GetImpactParameters
 577   // fTracks.fCdz //GetImpactParameters
 578   // fTracks.fCzz //GetImpactParameters
 579   // fTracks.fP   //GetPxPyPz
 580   // fTracks.fR   //GetMass
 581   // fTracks.fP   //GetMass
 582   // fTracks.fKinkIndexes
 583
 584   UInt_t status = esdTrack->GetStatus();
 585
 586   // getting quality parameters from the ESD track
 587   Int_t nClustersITS = esdTrack->GetITSclusters(0);
 588   Int_t nClustersTPC = esdTrack->GetTPCclusters(0);
 589
 590   Float_t chi2PerClusterITS = -1;
 591   Float_t chi2PerClusterTPC = -1;
 592   if (nClustersITS!=0)
 593     chi2PerClusterITS = esdTrack->GetITSchi2()/Float_t(nClustersITS);
 594   if (nClustersTPC!=0)
 595     chi2PerClusterTPC = esdTrack->GetTPCchi2()/Float_t(nClustersTPC);
 596   Double_t extCov[15];
 597   esdTrack->GetExternalCovariance(extCov);
 598
 599   // getting the track to vertex parameters
 600   Float_t nSigmaToVertex = GetSigmaToVertex(esdTrack);
 601
 602   Float_t b[2];
 603   Float_t bCov[3];
 604   esdTrack->GetImpactParameters(b,bCov);
 605   if (bCov[0]<=0 || bCov[2]<=0) {
 606     AliDebug(1, "Estimated b resolution lower or equal zero!");
 607     bCov[0]=0; bCov[2]=0;
 608   }
 609
 610   Float_t dcaToVertexXY = b[0];
 611   Float_t dcaToVertexZ = b[1];
 612
 613   Float_t dcaToVertex = -1;
 614
 615   if (fCutDCAToVertex2D)
 616   {
 617     dcaToVertex = TMath::Sqrt(dcaToVertexXY*dcaToVertexXY/fCutMaxDCAToVertexXY/fCutMaxDCAToVertexXY + dcaToVertexZ*dcaToVertexZ/fCutMaxDCAToVertexZ/fCutMaxDCAToVertexZ);
 618   }
 619   else
 620     dcaToVertex = TMath::Sqrt(dcaToVertexXY*dcaToVertexXY + dcaToVertexZ*dcaToVertexZ);
 621
 622   // getting the kinematic variables of the track
 623   // (assuming the mass is known)
 624   Double_t p[3];
 625   esdTrack->GetPxPyPz(p);
 626
 627   Float_t momentum = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2) + TMath::Power(p[2],2));
 628   Float_t pt       = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2));
 629   Float_t energy   = TMath::Sqrt(TMath::Power(esdTrack->GetMass(),2) + TMath::Power(momentum,2));
 630
 631
 632   //y-eta related calculations
 633   Float_t eta = -100.;
 634   Float_t y   = -100.;
 635   if((momentum != TMath::Abs(p[2]))&&(momentum != 0))
 636     eta = 0.5*TMath::Log((momentum + p[2])/(momentum - p[2]));
 637   if((energy != TMath::Abs(p[2]))&&(momentum != 0))
 638     y = 0.5*TMath::Log((energy + p[2])/(energy - p[2]));
 639
 640
 641   //########################################################################
 642   // cut the track?
 643
 644   Bool_t cuts[kNCuts];
 645   for (Int_t i=0; i<kNCuts; i++) cuts[i]=kFALSE;
 646
 647   // track quality cuts
 648   if (fCutRequireTPCRefit && (status&AliESDtrack::kTPCrefit)==0)
 649     cuts[0]=kTRUE;
 650   if (fCutRequireITSRefit && (status&AliESDtrack::kITSrefit)==0)
 651     cuts[1]=kTRUE;
 652   if (nClustersTPC<fCutMinNClusterTPC)
 653     cuts[2]=kTRUE;
 654   if (nClustersITS<fCutMinNClusterITS)
 655     cuts[3]=kTRUE;
 656   if (chi2PerClusterTPC>fCutMaxChi2PerClusterTPC)
 657     cuts[4]=kTRUE;
 658   if (chi2PerClusterITS>fCutMaxChi2PerClusterITS)
 659     cuts[5]=kTRUE;
 660   if (extCov[0]  > fCutMaxC11)
 661     cuts[6]=kTRUE;
 662   if (extCov[2]  > fCutMaxC22)
 663     cuts[7]=kTRUE;
 664   if (extCov[5]  > fCutMaxC33)
 665     cuts[8]=kTRUE;
 666   if (extCov[9]  > fCutMaxC44)
 667     cuts[9]=kTRUE;
 668   if (extCov[14]  > fCutMaxC55)
 669     cuts[10]=kTRUE;
 670   if (nSigmaToVertex > fCutNsigmaToVertex && fCutSigmaToVertexRequired)
 671     cuts[11] = kTRUE;
 672   // if n sigma could not be calculated
 673   if (nSigmaToVertex<0 && fCutSigmaToVertexRequired)
 674     cuts[12]=kTRUE;
 675   if (!fCutAcceptKinkDaughters && esdTrack->GetKinkIndex(0)>0)
 676     cuts[13]=kTRUE;
 677   // track kinematics cut
 678   if((momentum < fPMin) || (momentum > fPMax))
 679     cuts[14]=kTRUE;
 680   if((pt < fPtMin) || (pt > fPtMax))
 681     cuts[15] = kTRUE;
 682   if((p[0] < fPxMin) || (p[0] > fPxMax))
 683     cuts[16] = kTRUE;
 684   if((p[1] < fPyMin) || (p[1] > fPyMax))
 685     cuts[17] = kTRUE;
 686   if((p[2] < fPzMin) || (p[2] > fPzMax))
 687     cuts[18] = kTRUE;
 688   if((eta < fEtaMin) || (eta > fEtaMax))
 689     cuts[19] = kTRUE;
 690   if((y < fRapMin) || (y > fRapMax))
 691     cuts[20] = kTRUE;
 692   if (fCutDCAToVertex2D && dcaToVertex > 1)
 693     cuts[21] = kTRUE;
 694   if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexXY) > fCutMaxDCAToVertexXY)
 695     cuts[22] = kTRUE;
 696   if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexZ) > fCutMaxDCAToVertexZ)
 697     cuts[23] = kTRUE;
 698   if (fCutDCAToVertex2D && fCutMinDCAToVertexXY > 0 && fCutMinDCAToVertexZ > 0 && dcaToVertexXY*dcaToVertexXY/fCutMinDCAToVertexXY/fCutMinDCAToVertexXY + dcaToVertexZ*dcaToVertexZ/fCutMinDCAToVertexZ/fCutMinDCAToVertexZ < 1)
 699     cuts[24] = kTRUE;
 700   if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexXY) < fCutMinDCAToVertexXY)
 701     cuts[25] = kTRUE;
 702   if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexZ) < fCutMinDCAToVertexZ)
 703     cuts[26] = kTRUE;
 704
 705   for (Int_t i = 0; i < 3; i++)
 706     cuts[27+i] = !CheckITSClusterRequirement(fCutClusterRequirementITS[i], esdTrack->HasPointOnITSLayer(i*2), esdTrack->HasPointOnITSLayer(i*2+1));
 707
 708   Bool_t cut=kFALSE;
 709   for (Int_t i=0; i<kNCuts; i++)
 710     if (cuts[i]) {cut = kTRUE;}
 711
 712
 713
 714   //########################################################################
 715   // filling histograms
 716   if (fHistogramsOn) {
 717     fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n tracks")));
 718     if (cut)
 719       fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n cut tracks")));
 720
 721     for (Int_t i=0; i<kNCuts; i++) {
 722       if (cuts[i])
 723         fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin(fgkCutNames[i])));
 724
 725       for (Int_t j=i; j<kNCuts; j++) {
 726         if (cuts[i] && cuts[j]) {
 727           Float_t xC = fhCutCorrelation->GetXaxis()->GetBinCenter(fhCutCorrelation->GetXaxis()->FindBin(fgkCutNames[i]));
 728           Float_t yC = fhCutCorrelation->GetYaxis()->GetBinCenter(fhCutCorrelation->GetYaxis()->FindBin(fgkCutNames[j]));
 729           fhCutCorrelation->Fill(xC, yC);
 730         }
 731       }
 732     }
 733   }
 734
 735   // now we loop over the filling of the histograms twice: once "before" the cut, once "after"
 736   // the code is not in a function due to too many local variables that would need to be passed
 737
 738   for (Int_t id = 0; id < 2; id++)
 739   {
 740     // id = 0 --> before cut
 741     // id = 1 --> after cut
 742
 743     if (fHistogramsOn)
 744     {
 745       fhNClustersITS[id]->Fill(nClustersITS);
 746       fhNClustersTPC[id]->Fill(nClustersTPC);
 747       fhChi2PerClusterITS[id]->Fill(chi2PerClusterITS);
 748       fhChi2PerClusterTPC[id]->Fill(chi2PerClusterTPC);
 749
 750       fhC11[id]->Fill(extCov[0]);
 751       fhC22[id]->Fill(extCov[2]);
 752       fhC33[id]->Fill(extCov[5]);
 753       fhC44[id]->Fill(extCov[9]);
 754       fhC55[id]->Fill(extCov[14]);
 755
 756       fhPt[id]->Fill(pt);
 757       fhEta[id]->Fill(eta);
 758
 759       Float_t bRes[2];
 760       bRes[0] = TMath::Sqrt(bCov[0]);
 761       bRes[1] = TMath::Sqrt(bCov[2]);
 762
 763       fhDZ[id]->Fill(b[1]);
 764       fhDXY[id]->Fill(b[0]);
 765       fhDXYDZ[id]->Fill(dcaToVertex);
 766       fhDXYvsDZ[id]->Fill(b[1],b[0]);
 767
 768       if (bRes[0]!=0 && bRes[1]!=0) {
 769         fhDZNormalized[id]->Fill(b[1]/bRes[1]);
 770         fhDXYNormalized[id]->Fill(b[0]/bRes[0]);
 771         fhDXYvsDZNormalized[id]->Fill(b[1]/bRes[1], b[0]/bRes[0]);
 772         fhNSigmaToVertex[id]->Fill(nSigmaToVertex);
 773       }
 774     }
 775
 776     // cut the track
 777     if (cut)
 778       return kFALSE;
 779   }
 780
 781   return kTRUE;
 782 }
 783
 784 //____________________________________________________________________
 785 Bool_t AliESDtrackCuts::CheckITSClusterRequirement(ITSClusterRequirement req, Bool_t clusterL1, Bool_t clusterL2)
 786 {
 787   // checks if the cluster requirement is fullfilled (in this case: return kTRUE)
 788
 789   switch (req)
 790   {
 791         case kOff:        return kTRUE;
 792         case kNone:       return !clusterL1 && !clusterL2;
 793         case kAny:        return clusterL1 || clusterL2;
 794         case kFirst:      return clusterL1;
 795         case kOnlyFirst:  return clusterL1 && !clusterL2;
 796         case kSecond:     return clusterL2;
 797         case kOnlySecond: return clusterL2 && !clusterL1;
 798         case kBoth:       return clusterL1 && clusterL2;
 799   }
 800
 801   return kFALSE;
 802 }
 803
 804 //____________________________________________________________________
 805 AliESDtrack* AliESDtrackCuts::GetTPCOnlyTrack(AliESDEvent* esd, Int_t iTrack)
 806 {
 807   // creates a TPC only track from the given esd track
 808   // the track has to be deleted by the user
 809   //
 810   // NB. most of the functionality to get a TPC only track from an ESD track is in AliESDtrack, where it should be
 811   // there are only missing propagations here that are needed for old data
 812   // this function will therefore become obsolete
 813   //
 814   // adapted from code provided by CKB
 815
 816   if (!esd->GetPrimaryVertexTPC())
 817     return 0; // No TPC vertex no TPC tracks
 818
 819   if(!esd->GetPrimaryVertexTPC()->GetStatus())
 820     return 0; // TPC Vertex is created by default in AliESDEvent, do not use in this case
 821
 822   AliESDtrack* track = esd->GetTrack(iTrack);
 823   if (!track)
 824     return 0;
 825
 826   AliESDtrack *tpcTrack = new AliESDtrack();
 827
 828   // This should have been done during the reconstruction
 829   // fixed by Juri in r26675
 830   // but recalculate for older data CKB
 831   Float_t p[2],cov[3];
 832   track->GetImpactParametersTPC(p,cov);
 833   if(p[0]==0&&p[1]==0)
 834     track->RelateToVertexTPC(esd->GetPrimaryVertexTPC(),esd->GetMagneticField(),kVeryBig);
 835   // BKC
 836
 837   // only true if we have a tpc track
 838   if (!track->FillTPCOnlyTrack(*tpcTrack))
 839   {
 840     delete tpcTrack;
 841     return 0;
 842   }
 843
 844   // propagate to Vertex
 845   // not needed for normal reconstructed ESDs...
 846   // Double_t pTPC[2],covTPC[3];
 847   // tpcTrack->PropagateToDCA(esd->GetPrimaryVertexTPC(), esd->GetMagneticField(), 10000,  pTPC, covTPC);
 848
 849   return tpcTrack;
 850 }
 851
 852 //____________________________________________________________________
 853 TObjArray* AliESDtrackCuts::GetAcceptedTracks(AliESDEvent* esd,Bool_t bTPC)
 854 {
 855   //
 856   // returns an array of all tracks that pass the cuts
 857   // or an array of TPC only tracks (propagated to the TPC vertex during reco)
 858   // tracks that pass the cut
 859
 860   TObjArray* acceptedTracks = new TObjArray();
 861
 862   // loop over esd tracks
 863   for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
 864     if(bTPC){
 865       if(!esd->GetPrimaryVertexTPC())return acceptedTracks; // No TPC vertex no TPC tracks
 866       if(!esd->GetPrimaryVertexTPC()->GetStatus())return acceptedTracks; // No proper TPC vertex, only the default
 867
 868       AliESDtrack *tpcTrack = GetTPCOnlyTrack(esd, iTrack);
 869       if (!tpcTrack)
 870         continue;
 871
 872       if (AcceptTrack(tpcTrack)) {
 873         acceptedTracks->Add(tpcTrack);
 874       }
 875       else
 876         delete tpcTrack;
 877     }
 878     else
 879     {
 880       AliESDtrack* track = esd->GetTrack(iTrack);
 881       if(AcceptTrack(track))
 882         acceptedTracks->Add(track);
 883     }
 884   }
 885   if(bTPC)acceptedTracks->SetOwner(kTRUE);
 886   return acceptedTracks;
 887 }
 888
 889 //____________________________________________________________________
 890 Int_t AliESDtrackCuts::CountAcceptedTracks(AliESDEvent* esd)
 891 {
 892   //
 893   // returns an the number of tracks that pass the cuts
 894   //
 895
 896   Int_t count = 0;
 897
 898   // loop over esd tracks
 899   for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
 900     AliESDtrack* track = esd->GetTrack(iTrack);
 901     if (AcceptTrack(track))
 902       count++;
 903   }
 904
 905   return count;
 906 }
 907
 908 //____________________________________________________________________
 909  void AliESDtrackCuts::DefineHistograms(Int_t color) {
 910    //
 911    // diagnostics histograms are defined
 912    //
 913
 914    fHistogramsOn=kTRUE;
 915
 916    Bool_t oldStatus = TH1::AddDirectoryStatus();
 917    TH1::AddDirectory(kFALSE);
 918
 919    //###################################################################################
 920    // defining histograms
 921
 922    fhCutStatistics = new TH1F("cut_statistics","cut statistics",kNCuts+4,-0.5,kNCuts+3.5);
 923
 924    fhCutStatistics->GetXaxis()->SetBinLabel(1,"n tracks");
 925    fhCutStatistics->GetXaxis()->SetBinLabel(2,"n cut tracks");
 926
 927    fhCutCorrelation = new TH2F("cut_correlation","cut correlation",kNCuts,-0.5,kNCuts-0.5,kNCuts,-0.5,kNCuts-0.5);;
 928
 929    for (Int_t i=0; i<kNCuts; i++) {
 930      fhCutStatistics->GetXaxis()->SetBinLabel(i+4,fgkCutNames[i]);
 931      fhCutCorrelation->GetXaxis()->SetBinLabel(i+1,fgkCutNames[i]);
 932      fhCutCorrelation->GetYaxis()->SetBinLabel(i+1,fgkCutNames[i]);
 933    }
 934
 935   fhCutStatistics  ->SetLineColor(color);
 936   fhCutCorrelation ->SetLineColor(color);
 937   fhCutStatistics  ->SetLineWidth(2);
 938   fhCutCorrelation ->SetLineWidth(2);
 939
 940   for (Int_t i=0; i<2; i++) {
 941     fhNClustersITS[i]        = new TH1F("nClustersITS"    ,"",8,-0.5,7.5);
 942     fhNClustersTPC[i]        = new TH1F("nClustersTPC"     ,"",165,-0.5,164.5);
 943     fhChi2PerClusterITS[i]   = new TH1F("chi2PerClusterITS","",500,0,10);
 944     fhChi2PerClusterTPC[i]   = new TH1F("chi2PerClusterTPC","",500,0,10);
 945
 946     fhC11[i]                 = new TH1F("covMatrixDiagonal11","",2000,0,20);
 947     fhC22[i]                 = new TH1F("covMatrixDiagonal22","",2000,0,20);
 948     fhC33[i]                 = new TH1F("covMatrixDiagonal33","",1000,0,0.1);
 949     fhC44[i]                 = new TH1F("covMatrixDiagonal44","",1000,0,0.1);
 950     fhC55[i]                 = new TH1F("covMatrixDiagonal55","",1000,0,5);
 951
 952     fhDXY[i]                 = new TH1F("dXY"    ,"",500,-10,10);
 953     fhDZ[i]                  = new TH1F("dZ"     ,"",500,-10,10);
 954     fhDXYDZ[i]               = new TH1F("dXYDZ"  ,"",500,0,10);
 955     fhDXYvsDZ[i]             = new TH2F("dXYvsDZ","",200,-10,10,200,-10,10);
 956
 957     fhDXYNormalized[i]       = new TH1F("dXYNormalized"    ,"",500,-10,10);
 958     fhDZNormalized[i]        = new TH1F("dZNormalized"     ,"",500,-10,10);
 959     fhDXYvsDZNormalized[i]   = new TH2F("dXYvsDZNormalized","",200,-10,10,200,-10,10);
 960
 961     fhNSigmaToVertex[i]      = new TH1F("nSigmaToVertex","",500,0,10);
 962
 963     fhPt[i]                  = new TH1F("pt"     ,"p_{T} distribution;p_{T} (GeV/c)", 800, 0.0, 10.0);
 964     fhEta[i]                 = new TH1F("eta"     ,"#eta distribution;#eta",40,-2.0,2.0);
 965
 966     fhNClustersITS[i]->SetTitle("n ITS clusters");
 967     fhNClustersTPC[i]->SetTitle("n TPC clusters");
 968     fhChi2PerClusterITS[i]->SetTitle("#Chi^{2} per ITS cluster");
 969     fhChi2PerClusterTPC[i]->SetTitle("#Chi^{2} per TPC cluster");
 970
 971     fhC11[i]->SetTitle("cov 11 : #sigma_{y}^{2} [cm^{2}]");
 972     fhC22[i]->SetTitle("cov 22 : #sigma_{z}^{2} [cm^{2}]");
 973     fhC33[i]->SetTitle("cov 33 : #sigma_{sin(#phi)}^{2}");
 974     fhC44[i]->SetTitle("cov 44 : #sigma_{tan(#theta_{dip})}^{2}");
 975     fhC55[i]->SetTitle("cov 55 : #sigma_{1/p_{T}}^{2} [(c/GeV)^2]");
 976
 977     fhDXY[i]->SetXTitle("transverse impact parameter (cm)");
 978     fhDZ[i]->SetXTitle("longitudinal impact parameter (cm)");
 979     fhDXYDZ[i]->SetTitle("absolute impact parameter;sqrt(dXY**2 + dZ**2)  (cm)");
 980     fhDXYvsDZ[i]->SetXTitle("longitudinal impact parameter (cm)");
 981     fhDXYvsDZ[i]->SetYTitle("transverse impact parameter (cm)");
 982
 983     fhDXYNormalized[i]->SetTitle("normalized trans impact par (n#sigma)");
 984     fhDZNormalized[i]->SetTitle("normalized long impact par (n#sigma)");
 985     fhDXYvsDZNormalized[i]->SetTitle("normalized long impact par (n#sigma)");
 986     fhDXYvsDZNormalized[i]->SetYTitle("normalized trans impact par (n#sigma)");
 987     fhNSigmaToVertex[i]->SetTitle("n #sigma to vertex");
 988
 989     fhNClustersITS[i]->SetLineColor(color);   fhNClustersITS[i]->SetLineWidth(2);
 990     fhNClustersTPC[i]->SetLineColor(color);   fhNClustersTPC[i]->SetLineWidth(2);
 991     fhChi2PerClusterITS[i]->SetLineColor(color);   fhChi2PerClusterITS[i]->SetLineWidth(2);
 992     fhChi2PerClusterTPC[i]->SetLineColor(color);   fhChi2PerClusterTPC[i]->SetLineWidth(2);
 993
 994     fhC11[i]->SetLineColor(color);   fhC11[i]->SetLineWidth(2);
 995     fhC22[i]->SetLineColor(color);   fhC22[i]->SetLineWidth(2);
 996     fhC33[i]->SetLineColor(color);   fhC33[i]->SetLineWidth(2);
 997     fhC44[i]->SetLineColor(color);   fhC44[i]->SetLineWidth(2);
 998     fhC55[i]->SetLineColor(color);   fhC55[i]->SetLineWidth(2);
 999
1000     fhDXY[i]->SetLineColor(color);   fhDXY[i]->SetLineWidth(2);
1001     fhDZ[i]->SetLineColor(color);    fhDZ[i]->SetLineWidth(2);
1002     fhDXYDZ[i]->SetLineColor(color); fhDXYDZ[i]->SetLineWidth(2);
1003
1004     fhDXYNormalized[i]->SetLineColor(color);   fhDXYNormalized[i]->SetLineWidth(2);
1005     fhDZNormalized[i]->SetLineColor(color);    fhDZNormalized[i]->SetLineWidth(2);
1006     fhNSigmaToVertex[i]->SetLineColor(color);  fhNSigmaToVertex[i]->SetLineWidth(2);
1007   }
1008
1009   // The number of sigmas to the vertex is per definition gaussian
1010   ffDTheoretical = new TF1("nSigmaToVertexTheoretical","([0]/2.506628274)*exp(-(x**2)/2)",0,50);
1011   ffDTheoretical->SetParameter(0,1);
1012
1013   TH1::AddDirectory(oldStatus);
1014 }
1015
1016 //____________________________________________________________________
1017 Bool_t AliESDtrackCuts::LoadHistograms(const Char_t* dir)
1018 {
1019   //
1020   // loads the histograms from a file
1021   // if dir is empty a directory with the name of this object is taken (like in SaveHistogram)
1022   //
1023
1024   if (!dir)
1025     dir = GetName();
1026
1027   if (!gDirectory->cd(dir))
1028     return kFALSE;
1029
1030   ffDTheoretical = dynamic_cast<TF1*> (gDirectory->Get("nSigmaToVertexTheory"));
1031
1032   fhCutStatistics = dynamic_cast<TH1F*> (gDirectory->Get("cut_statistics"));
1033   fhCutCorrelation = dynamic_cast<TH2F*> (gDirectory->Get("cut_correlation"));
1034
1035   for (Int_t i=0; i<2; i++) {
1036     if (i==0)
1037     {
1038       gDirectory->cd("before_cuts");
1039     }
1040     else
1041       gDirectory->cd("after_cuts");
1042
1043     fhNClustersITS[i]      = dynamic_cast<TH1F*> (gDirectory->Get("nClustersITS"     ));
1044     fhNClustersTPC[i]      = dynamic_cast<TH1F*> (gDirectory->Get("nClustersTPC"     ));
1045     fhChi2PerClusterITS[i] = dynamic_cast<TH1F*> (gDirectory->Get("chi2PerClusterITS"));
1046     fhChi2PerClusterTPC[i] = dynamic_cast<TH1F*> (gDirectory->Get("chi2PerClusterTPC"));
1047
1048     fhC11[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal11"));
1049     fhC22[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal22"));
1050     fhC33[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal33"));
1051     fhC44[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal44"));
1052     fhC55[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal55"));
1053
1054     fhDXY[i] =     dynamic_cast<TH1F*> (gDirectory->Get("dXY"    ));
1055     fhDZ[i] =      dynamic_cast<TH1F*> (gDirectory->Get("dZ"     ));
1056     fhDXYDZ[i] =   dynamic_cast<TH1F*> (gDirectory->Get("dXYDZ"));
1057     fhDXYvsDZ[i] = dynamic_cast<TH2F*> (gDirectory->Get("dXYvsDZ"));
1058
1059     fhDXYNormalized[i] =     dynamic_cast<TH1F*> (gDirectory->Get("dXYNormalized"    ));
1060     fhDZNormalized[i] =      dynamic_cast<TH1F*> (gDirectory->Get("dZNormalized"     ));
1061     fhDXYvsDZNormalized[i] = dynamic_cast<TH2F*> (gDirectory->Get("dXYvsDZNormalized"));
1062     fhNSigmaToVertex[i] = dynamic_cast<TH1F*> (gDirectory->Get("nSigmaToVertex"));
1063
1064     fhPt[i] = dynamic_cast<TH1F*> (gDirectory->Get("pt"));
1065     fhEta[i] = dynamic_cast<TH1F*> (gDirectory->Get("eta"));
1066
1067     gDirectory->cd("../");
1068   }
1069
1070   gDirectory->cd("..");
1071
1072   return kTRUE;
1073 }
1074
1075 //____________________________________________________________________
1076 void AliESDtrackCuts::SaveHistograms(const Char_t* dir) {
1077   //
1078   // saves the histograms in a directory (dir)
1079   //
1080
1081   if (!fHistogramsOn) {
1082     AliDebug(0, "Histograms not on - cannot save histograms!!!");
1083     return;
1084   }
1085
1086   if (!dir)
1087     dir = GetName();
1088
1089   gDirectory->mkdir(dir);
1090   gDirectory->cd(dir);
1091
1092   gDirectory->mkdir("before_cuts");
1093   gDirectory->mkdir("after_cuts");
1094
1095   // a factor of 2 is needed since n sigma is positive
1096   ffDTheoretical->SetParameter(0,2*fhNSigmaToVertex[0]->Integral("width"));
1097   ffDTheoretical->Write("nSigmaToVertexTheory");
1098
1099   fhCutStatistics->Write();
1100   fhCutCorrelation->Write();
1101
1102   for (Int_t i=0; i<2; i++) {
1103     if (i==0)
1104       gDirectory->cd("before_cuts");
1105     else
1106       gDirectory->cd("after_cuts");
1107
1108     fhNClustersITS[i]        ->Write();
1109     fhNClustersTPC[i]        ->Write();
1110     fhChi2PerClusterITS[i]   ->Write();
1111     fhChi2PerClusterTPC[i]   ->Write();
1112
1113     fhC11[i]                 ->Write();
1114     fhC22[i]                 ->Write();
1115     fhC33[i]                 ->Write();
1116     fhC44[i]                 ->Write();
1117     fhC55[i]                 ->Write();
1118
1119     fhDXY[i]                 ->Write();
1120     fhDZ[i]                  ->Write();
1121     fhDXYDZ[i]               ->Write();
1122     fhDXYvsDZ[i]             ->Write();
1123
1124     fhDXYNormalized[i]       ->Write();
1125     fhDZNormalized[i]        ->Write();
1126     fhDXYvsDZNormalized[i]   ->Write();
1127     fhNSigmaToVertex[i]      ->Write();
1128
1129     fhPt[i]                  ->Write();
1130     fhEta[i]                 ->Write();
1131
1132     gDirectory->cd("../");
1133   }
1134
1135   gDirectory->cd("../");
1136 }
1137
1138 //____________________________________________________________________
1139 void AliESDtrackCuts::DrawHistograms()
1140 {
1141   // draws some histograms
1142
1143   TCanvas* canvas1 = new TCanvas(Form("%s_1", GetName()), "Track Quality Results1", 800, 800);
1144   canvas1->Divide(2, 2);
1145
1146   canvas1->cd(1);
1147   fhNClustersTPC[0]->SetStats(kFALSE);
1148   fhNClustersTPC[0]->Draw();
1149
1150   canvas1->cd(2);
1151   fhChi2PerClusterTPC[0]->SetStats(kFALSE);
1152   fhChi2PerClusterTPC[0]->Draw();
1153
1154   canvas1->cd(3);
1155   fhNSigmaToVertex[0]->SetStats(kFALSE);
1156   fhNSigmaToVertex[0]->GetXaxis()->SetRangeUser(0, 10);
1157   fhNSigmaToVertex[0]->Draw();
1158
1159   canvas1->SaveAs(Form("%s_%s.gif", GetName(), canvas1->GetName()));
1160
1161   TCanvas* canvas2 = new TCanvas(Form("%s_2", GetName()), "Track Quality Results2", 1200, 800);
1162   canvas2->Divide(3, 2);
1163
1164   canvas2->cd(1);
1165   fhC11[0]->SetStats(kFALSE);
1166   gPad->SetLogy();
1167   fhC11[0]->Draw();
1168
1169   canvas2->cd(2);
1170   fhC22[0]->SetStats(kFALSE);
1171   gPad->SetLogy();
1172   fhC22[0]->Draw();
1173
1174   canvas2->cd(3);
1175   fhC33[0]->SetStats(kFALSE);
1176   gPad->SetLogy();
1177   fhC33[0]->Draw();
1178
1179   canvas2->cd(4);
1180   fhC44[0]->SetStats(kFALSE);
1181   gPad->SetLogy();
1182   fhC44[0]->Draw();
1183
1184   canvas2->cd(5);
1185   fhC55[0]->SetStats(kFALSE);
1186   gPad->SetLogy();
1187   fhC55[0]->Draw();
1188
1189   canvas2->SaveAs(Form("%s_%s.gif", GetName(), canvas2->GetName()));
1190
1191   TCanvas* canvas3 = new TCanvas(Form("%s_3", GetName()), "Track Quality Results3", 1200, 800);
1192   canvas3->Divide(3, 2);
1193
1194   canvas3->cd(1);
1195   fhDXY[0]->SetStats(kFALSE);
1196   gPad->SetLogy();
1197   fhDXY[0]->Draw();
1198
1199   canvas3->cd(2);
1200   fhDZ[0]->SetStats(kFALSE);
1201   gPad->SetLogy();
1202   fhDZ[0]->Draw();
1203
1204   canvas3->cd(3);
1205   fhDXYvsDZ[0]->SetStats(kFALSE);
1206   gPad->SetLogz();
1207   gPad->SetRightMargin(0.15);
1208   fhDXYvsDZ[0]->Draw("COLZ");
1209
1210   canvas3->cd(4);
1211   fhDXYNormalized[0]->SetStats(kFALSE);
1212   gPad->SetLogy();
1213   fhDXYNormalized[0]->Draw();
1214
1215   canvas3->cd(5);
1216   fhDZNormalized[0]->SetStats(kFALSE);
1217   gPad->SetLogy();
1218   fhDZNormalized[0]->Draw();
1219
1220   canvas3->cd(6);
1221   fhDXYvsDZNormalized[0]->SetStats(kFALSE);
1222   gPad->SetLogz();
1223   gPad->SetRightMargin(0.15);
1224   fhDXYvsDZNormalized[0]->Draw("COLZ");
1225
1226   canvas3->SaveAs(Form("%s_%s.gif", GetName(), canvas3->GetName()));
1227
1228   TCanvas* canvas4 = new TCanvas(Form("%s_4", GetName()), "Track Quality Results4", 800, 500);
1229   canvas4->Divide(2, 1);
1230
1231   canvas4->cd(1);
1232   fhCutStatistics->SetStats(kFALSE);
1233   fhCutStatistics->LabelsOption("v");
1234   gPad->SetBottomMargin(0.3);
1235   fhCutStatistics->Draw();
1236
1237   canvas4->cd(2);
1238   fhCutCorrelation->SetStats(kFALSE);
1239   fhCutCorrelation->LabelsOption("v");
1240   gPad->SetBottomMargin(0.3);
1241   gPad->SetLeftMargin(0.3);
1242   fhCutCorrelation->Draw("COLZ");
1243
1244   canvas4->SaveAs(Form("%s_%s.gif", GetName(), canvas4->GetName()));
1245
1246   /*canvas->cd(1);
1247   fhDXYvsDZNormalized[0]->SetStats(kFALSE);
1248   fhDXYvsDZNormalized[0]->DrawCopy("COLZ");
1249
1250   canvas->cd(2);
1251   fhNClustersTPC[0]->SetStats(kFALSE);
1252   fhNClustersTPC[0]->DrawCopy();
1253
1254   canvas->cd(3);
1255   fhChi2PerClusterITS[0]->SetStats(kFALSE);
1256   fhChi2PerClusterITS[0]->DrawCopy();
1257   fhChi2PerClusterITS[1]->SetLineColor(2);
1258   fhChi2PerClusterITS[1]->DrawCopy("SAME");
1259
1260   canvas->cd(4);
1261   fhChi2PerClusterTPC[0]->SetStats(kFALSE);
1262   fhChi2PerClusterTPC[0]->DrawCopy();
1263   fhChi2PerClusterTPC[1]->SetLineColor(2);
1264   fhChi2PerClusterTPC[1]->DrawCopy("SAME");*/
1265 }
1266
1267 Float_t AliESDtrackCuts::GetMinNsigmaToVertex() const
1268 {
1269   // deprecated, please use GetMaxNsigmaToVertex
1270
1271   Printf("WARNING: AliESDtrackCuts::GetMinNsigmaToVertex is DEPRECATED and will be removed in the next release. Please use GetMaxNsigmaToVertex instead. Renaming was done to improve code readability.");
1272
1273   return GetMaxNsigmaToVertex();
1274 }
1275
1276 void AliESDtrackCuts::SetMinNsigmaToVertex(Float_t sigma)
1277 {
1278   // deprecated, will be removed in next release
1279
1280   Printf("WARNING: AliESDtrackCuts::SetMinNsigmaToVertex is DEPRECATED and will be removed in the next release. Please use SetMaxNsigmaToVertex instead. Renaming was done to improve code readability.");
1281
1282   SetMaxNsigmaToVertex(sigma);
1283 }
1284
1285 void AliESDtrackCuts::SetDCAToVertex(Float_t dist)
1286 {
1287   // deprecated, will be removed in next release
1288
1289   SetMaxDCAToVertex(dist);
1290 }
1291
1292 void AliESDtrackCuts::SetMaxDCAToVertex(Float_t dist)
1293 {
1294   // deprecated, will be removed in next release
1295
1296   Printf("WARNING: AliESDtrackCuts::SetMaxDCAToVertex is DEPRECATED and will be removed in the next release. Please use SetDCAToVertexXY(dist) and SetDCAToVertexZ(dist) and SetDCAToVertex2D(kTRUE)");
1297
1298   SetDCAToVertexXY(dist);
1299   SetDCAToVertexZ(dist);
1300   SetDCAToVertex2D(kTRUE);
1301 }
1302
1303 void AliESDtrackCuts::SetDCAToVertexXY(Float_t dist)
1304 {
1305   // deprecated, will be removed in next release
1306
1307   Printf("WARNING: AliESDtrackCuts::SetDCAToVertexXY is DEPRECATED and will be removed in the next release. Please use SetMaxDCAToVertexXY instead. Renaming was done to improve code readability.");
1308
1309   SetMaxDCAToVertexXY(dist);
1310 }
1311
1312 void AliESDtrackCuts::SetDCAToVertexZ(Float_t dist)
1313 {
1314   // deprecated, will be removed in next release
1315
1316   Printf("WARNING: AliESDtrackCuts::SetDCAToVertexZ is DEPRECATED and will be removed in the next release. Please use SetMaxDCAToVertexZ instead. Renaming was done to improve code readability.");
1317
1318   SetMaxDCAToVertexZ(dist);
1319 }