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