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 #include <TH2F.h>
  29 #include <TF1.h>
  30
  31 //____________________________________________________________________
  32 ClassImp(AliESDtrackCuts)
  33
  34 // Cut names
  35 const Char_t* AliESDtrackCuts::fgkCutNames[kNCuts] = {
  36  "require TPC refit",
  37  "require TPC standalone",
  38  "require ITS refit",
  39  "n clusters TPC",
  40  "n clusters ITS",
  41  "#Chi^{2}/cluster TPC",
  42  "#Chi^{2}/cluster ITS",
  43  "cov 11",
  44  "cov 22",
  45  "cov 33",
  46  "cov 44",
  47  "cov 55",
  48  "trk-to-vtx",
  49  "trk-to-vtx failed",
  50  "kink daughters",
  51  "p",
  52  "p_{T}",
  53  "p_{x}",
  54  "p_{y}",
  55  "p_{z}",
  56  "eta",
  57  "y",
  58  "trk-to-vtx max dca 2D absolute",
  59  "trk-to-vtx max dca xy absolute",
  60  "trk-to-vtx max dca z absolute",
  61  "trk-to-vtx min dca 2D absolute",
  62  "trk-to-vtx min dca xy absolute",
  63  "trk-to-vtx min dca z absolute",
  64  "SPD cluster requirement",
  65  "SDD cluster requirement",
  66  "SSD cluster requirement",
  67  "require ITS stand-alone",
  68  "rel 1/pt uncertainty",
  69  "TPC n shared clusters",
  70  "TPC rel shared clusters",
  71  "require ITS Pid"
  72 };
  73
  74 //____________________________________________________________________
  75 AliESDtrackCuts::AliESDtrackCuts(const Char_t* name, const Char_t* title) : AliAnalysisCuts(name,title),
  76   fCutMinNClusterTPC(0),
  77   fCutMinNClusterITS(0),
  78   fCutMaxChi2PerClusterTPC(0),
  79   fCutMaxChi2PerClusterITS(0),
  80   fCutMaxC11(0),
  81   fCutMaxC22(0),
  82   fCutMaxC33(0),
  83   fCutMaxC44(0),
  84   fCutMaxC55(0),
  85   fCutMaxRel1PtUncertainty(0),
  86   fCutAcceptKinkDaughters(0),
  87   fCutAcceptSharedTPCClusters(0),
  88   fCutMaxFractionSharedTPCClusters(0),
  89   fCutRequireTPCRefit(0),
  90   fCutRequireTPCStandAlone(0),
  91   fCutRequireITSRefit(0),
  92   fCutRequireITSPid(0),
  93   fCutRequireITSStandAlone(0),
  94   fCutRequireITSpureSA(0),
  95   fCutNsigmaToVertex(0),
  96   fCutSigmaToVertexRequired(0),
  97   fCutMaxDCAToVertexXY(0),
  98   fCutMaxDCAToVertexZ(0),
  99   fCutMinDCAToVertexXY(0),
 100   fCutMinDCAToVertexZ(0),
 101   fCutMaxDCAToVertexXYPtDep(""),
 102   fCutMaxDCAToVertexZPtDep(""),
 103   fCutMinDCAToVertexXYPtDep(""),
 104   fCutMinDCAToVertexZPtDep(""),
 105   f1CutMaxDCAToVertexXYPtDep(0x0),
 106   f1CutMaxDCAToVertexZPtDep(0x0),
 107   f1CutMinDCAToVertexXYPtDep(0x0),
 108   f1CutMinDCAToVertexZPtDep(0x0),
 109   fCutDCAToVertex2D(0),
 110   fPMin(0),
 111   fPMax(0),
 112   fPtMin(0),
 113   fPtMax(0),
 114   fPxMin(0),
 115   fPxMax(0),
 116   fPyMin(0),
 117   fPyMax(0),
 118   fPzMin(0),
 119   fPzMax(0),
 120   fEtaMin(0),
 121   fEtaMax(0),
 122   fRapMin(0),
 123   fRapMax(0),
 124   fHistogramsOn(0),
 125   ffDTheoretical(0),
 126   fhCutStatistics(0),
 127   fhCutCorrelation(0)
 128 {
 129   //
 130   // constructor
 131   //
 132
 133   Init();
 134
 135   //##############################################################################
 136   // setting default cuts
 137   SetMinNClustersTPC();
 138   SetMinNClustersITS();
 139   SetMaxChi2PerClusterTPC();
 140   SetMaxChi2PerClusterITS();
 141   SetMaxCovDiagonalElements();
 142   SetMaxRel1PtUncertainty();
 143   SetRequireTPCRefit();
 144   SetRequireTPCStandAlone();
 145   SetRequireITSRefit();
 146   SetRequireITSPid(kFALSE);
 147   SetRequireITSStandAlone(kFALSE);
 148   SetRequireITSPureStandAlone(kFALSE);
 149   SetAcceptKinkDaughters();
 150   SetAcceptSharedTPCClusters();
 151   SetMaxFractionSharedTPCClusters();
 152   SetMaxNsigmaToVertex();
 153   SetMaxDCAToVertexXY();
 154   SetMaxDCAToVertexZ();
 155   SetDCAToVertex2D();
 156   SetMinDCAToVertexXY();
 157   SetMinDCAToVertexZ();
 158   SetPRange();
 159   SetPtRange();
 160   SetPxRange();
 161   SetPyRange();
 162   SetPzRange();
 163   SetEtaRange();
 164   SetRapRange();
 165   SetClusterRequirementITS(kSPD);
 166   SetClusterRequirementITS(kSDD);
 167   SetClusterRequirementITS(kSSD);
 168
 169   SetHistogramsOn();
 170 }
 171
 172 //_____________________________________________________________________________
 173 AliESDtrackCuts::AliESDtrackCuts(const AliESDtrackCuts &c) : AliAnalysisCuts(c),
 174   fCutMinNClusterTPC(0),
 175   fCutMinNClusterITS(0),
 176   fCutMaxChi2PerClusterTPC(0),
 177   fCutMaxChi2PerClusterITS(0),
 178   fCutMaxC11(0),
 179   fCutMaxC22(0),
 180   fCutMaxC33(0),
 181   fCutMaxC44(0),
 182   fCutMaxC55(0),
 183   fCutMaxRel1PtUncertainty(0),
 184   fCutAcceptKinkDaughters(0),
 185   fCutAcceptSharedTPCClusters(0),
 186   fCutMaxFractionSharedTPCClusters(0),
 187   fCutRequireTPCRefit(0),
 188   fCutRequireTPCStandAlone(0),
 189   fCutRequireITSRefit(0),
 190   fCutRequireITSPid(0),
 191   fCutRequireITSStandAlone(0),
 192   fCutRequireITSpureSA(0),
 193   fCutNsigmaToVertex(0),
 194   fCutSigmaToVertexRequired(0),
 195   fCutMaxDCAToVertexXY(0),
 196   fCutMaxDCAToVertexZ(0),
 197   fCutMinDCAToVertexXY(0),
 198   fCutMinDCAToVertexZ(0),
 199   fCutMaxDCAToVertexXYPtDep(""),
 200   fCutMaxDCAToVertexZPtDep(""),
 201   fCutMinDCAToVertexXYPtDep(""),
 202   fCutMinDCAToVertexZPtDep(""),
 203   f1CutMaxDCAToVertexXYPtDep(0x0),
 204   f1CutMaxDCAToVertexZPtDep(0x0),
 205   f1CutMinDCAToVertexXYPtDep(0x0),
 206   f1CutMinDCAToVertexZPtDep(0x0),
 207   fCutDCAToVertex2D(0),
 208   fPMin(0),
 209   fPMax(0),
 210   fPtMin(0),
 211   fPtMax(0),
 212   fPxMin(0),
 213   fPxMax(0),
 214   fPyMin(0),
 215   fPyMax(0),
 216   fPzMin(0),
 217   fPzMax(0),
 218   fEtaMin(0),
 219   fEtaMax(0),
 220   fRapMin(0),
 221   fRapMax(0),
 222   fHistogramsOn(0),
 223   ffDTheoretical(0),
 224   fhCutStatistics(0),
 225   fhCutCorrelation(0)
 226 {
 227   //
 228   // copy constructor
 229   //
 230
 231   ((AliESDtrackCuts &) c).Copy(*this);
 232 }
 233
 234 AliESDtrackCuts::~AliESDtrackCuts()
 235 {
 236   //
 237   // destructor
 238   //
 239
 240   for (Int_t i=0; i<2; i++) {
 241
 242     if (fhNClustersITS[i])
 243       delete fhNClustersITS[i];
 244     if (fhNClustersTPC[i])
 245       delete fhNClustersTPC[i];
 246     if (fhChi2PerClusterITS[i])
 247       delete fhChi2PerClusterITS[i];
 248     if (fhChi2PerClusterTPC[i])
 249       delete fhChi2PerClusterTPC[i];
 250     if (fhC11[i])
 251       delete fhC11[i];
 252     if (fhC22[i])
 253       delete fhC22[i];
 254     if (fhC33[i])
 255       delete fhC33[i];
 256     if (fhC44[i])
 257       delete fhC44[i];
 258     if (fhC55[i])
 259       delete fhC55[i];
 260
 261     if (fhRel1PtUncertainty[i])
 262       delete fhRel1PtUncertainty[i];
 263
 264     if (fhDXY[i])
 265       delete fhDXY[i];
 266     if (fhDZ[i])
 267       delete fhDZ[i];
 268     if (fhDXYDZ[i])
 269       delete fhDXYDZ[i];
 270     if (fhDXYvsDZ[i])
 271       delete fhDXYvsDZ[i];
 272
 273     if (fhDXYNormalized[i])
 274       delete fhDXYNormalized[i];
 275     if (fhDZNormalized[i])
 276       delete fhDZNormalized[i];
 277     if (fhDXYvsDZNormalized[i])
 278       delete fhDXYvsDZNormalized[i];
 279     if (fhNSigmaToVertex[i])
 280       delete fhNSigmaToVertex[i];
 281     if (fhPt[i])
 282       delete fhPt[i];
 283     if (fhEta[i])
 284       delete fhEta[i];
 285   }
 286
 287   if(f1CutMaxDCAToVertexXYPtDep)delete f1CutMaxDCAToVertexXYPtDep;
 288   f1CutMaxDCAToVertexXYPtDep = 0;
 289   if( f1CutMaxDCAToVertexZPtDep) delete  f1CutMaxDCAToVertexZPtDep;
 290   f1CutMaxDCAToVertexZPtDep = 0;
 291   if( f1CutMinDCAToVertexXYPtDep)delete  f1CutMinDCAToVertexXYPtDep;
 292   f1CutMinDCAToVertexXYPtDep = 0;
 293   if(f1CutMinDCAToVertexZPtDep)delete  f1CutMinDCAToVertexZPtDep;
 294   f1CutMinDCAToVertexZPtDep = 0;
 295
 296
 297   if (ffDTheoretical)
 298     delete ffDTheoretical;
 299
 300   if (fhCutStatistics)
 301     delete fhCutStatistics;
 302   if (fhCutCorrelation)
 303     delete fhCutCorrelation;
 304 }
 305
 306 void AliESDtrackCuts::Init()
 307 {
 308   //
 309   // sets everything to zero
 310   //
 311
 312   fCutMinNClusterTPC = 0;
 313   fCutMinNClusterITS = 0;
 314
 315   fCutMaxChi2PerClusterTPC = 0;
 316   fCutMaxChi2PerClusterITS = 0;
 317
 318   for (Int_t i = 0; i < 3; i++)
 319         fCutClusterRequirementITS[i] = kOff;
 320
 321   fCutMaxC11 = 0;
 322   fCutMaxC22 = 0;
 323   fCutMaxC33 = 0;
 324   fCutMaxC44 = 0;
 325   fCutMaxC55 = 0;
 326
 327   fCutMaxRel1PtUncertainty = 0;
 328
 329   fCutAcceptKinkDaughters = 0;
 330   fCutAcceptSharedTPCClusters = 0;
 331   fCutMaxFractionSharedTPCClusters = 0;
 332   fCutRequireTPCRefit = 0;
 333   fCutRequireTPCStandAlone = 0;
 334   fCutRequireITSRefit = 0;
 335   fCutRequireITSPid = 0;
 336   fCutRequireITSStandAlone = 0;
 337   fCutRequireITSpureSA = 0;
 338
 339   fCutNsigmaToVertex = 0;
 340   fCutSigmaToVertexRequired = 0;
 341   fCutMaxDCAToVertexXY = 0;
 342   fCutMaxDCAToVertexZ = 0;
 343   fCutDCAToVertex2D = 0;
 344   fCutMinDCAToVertexXY = 0;
 345   fCutMinDCAToVertexZ = 0;
 346   fCutMaxDCAToVertexXYPtDep = "";
 347   fCutMaxDCAToVertexZPtDep = "";
 348   fCutMinDCAToVertexXYPtDep = "";
 349   fCutMinDCAToVertexZPtDep = "";
 350
 351   if(f1CutMaxDCAToVertexXYPtDep)delete f1CutMaxDCAToVertexXYPtDep;
 352   f1CutMaxDCAToVertexXYPtDep = 0;
 353   if( f1CutMaxDCAToVertexXYPtDep) delete  f1CutMaxDCAToVertexXYPtDep;
 354   f1CutMaxDCAToVertexXYPtDep = 0;
 355   if( f1CutMaxDCAToVertexZPtDep) delete  f1CutMaxDCAToVertexZPtDep;
 356   f1CutMaxDCAToVertexZPtDep = 0;
 357   if( f1CutMinDCAToVertexXYPtDep)delete  f1CutMinDCAToVertexXYPtDep;
 358   f1CutMinDCAToVertexXYPtDep = 0;
 359   if(f1CutMinDCAToVertexZPtDep)delete f1CutMinDCAToVertexZPtDep;
 360   f1CutMinDCAToVertexZPtDep = 0;
 361
 362
 363   fPMin = 0;
 364   fPMax = 0;
 365   fPtMin = 0;
 366   fPtMax = 0;
 367   fPxMin = 0;
 368   fPxMax = 0;
 369   fPyMin = 0;
 370   fPyMax = 0;
 371   fPzMin = 0;
 372   fPzMax = 0;
 373   fEtaMin = 0;
 374   fEtaMax = 0;
 375   fRapMin = 0;
 376   fRapMax = 0;
 377
 378   fHistogramsOn = kFALSE;
 379
 380   for (Int_t i=0; i<2; ++i)
 381   {
 382     fhNClustersITS[i] = 0;
 383     fhNClustersTPC[i] = 0;
 384
 385     fhChi2PerClusterITS[i] = 0;
 386     fhChi2PerClusterTPC[i] = 0;
 387
 388     fhC11[i] = 0;
 389     fhC22[i] = 0;
 390     fhC33[i] = 0;
 391     fhC44[i] = 0;
 392     fhC55[i] = 0;
 393
 394     fhRel1PtUncertainty[i] = 0;
 395
 396     fhDXY[i] = 0;
 397     fhDZ[i] = 0;
 398     fhDXYDZ[i] = 0;
 399     fhDXYvsDZ[i] = 0;
 400
 401     fhDXYNormalized[i] = 0;
 402     fhDZNormalized[i] = 0;
 403     fhDXYvsDZNormalized[i] = 0;
 404     fhNSigmaToVertex[i] = 0;
 405
 406     fhPt[i] = 0;
 407     fhEta[i] = 0;
 408   }
 409   ffDTheoretical = 0;
 410
 411   fhCutStatistics = 0;
 412   fhCutCorrelation = 0;
 413 }
 414
 415 //_____________________________________________________________________________
 416 AliESDtrackCuts &AliESDtrackCuts::operator=(const AliESDtrackCuts &c)
 417 {
 418   //
 419   // Assignment operator
 420   //
 421
 422   if (this != &c) ((AliESDtrackCuts &) c).Copy(*this);
 423   return *this;
 424 }
 425
 426 //_____________________________________________________________________________
 427 void AliESDtrackCuts::Copy(TObject &c) const
 428 {
 429   //
 430   // Copy function
 431   //
 432
 433   AliESDtrackCuts& target = (AliESDtrackCuts &) c;
 434
 435   target.Init();
 436
 437   target.fCutMinNClusterTPC = fCutMinNClusterTPC;
 438   target.fCutMinNClusterITS = fCutMinNClusterITS;
 439
 440   target.fCutMaxChi2PerClusterTPC = fCutMaxChi2PerClusterTPC;
 441   target.fCutMaxChi2PerClusterITS = fCutMaxChi2PerClusterITS;
 442
 443   for (Int_t i = 0; i < 3; i++)
 444     target.fCutClusterRequirementITS[i] = fCutClusterRequirementITS[i];
 445
 446   target.fCutMaxC11 = fCutMaxC11;
 447   target.fCutMaxC22 = fCutMaxC22;
 448   target.fCutMaxC33 = fCutMaxC33;
 449   target.fCutMaxC44 = fCutMaxC44;
 450   target.fCutMaxC55 = fCutMaxC55;
 451
 452   target.fCutMaxRel1PtUncertainty = fCutMaxRel1PtUncertainty;
 453
 454   target.fCutAcceptKinkDaughters = fCutAcceptKinkDaughters;
 455   target.fCutAcceptSharedTPCClusters = fCutAcceptSharedTPCClusters;
 456   target.fCutMaxFractionSharedTPCClusters = fCutMaxFractionSharedTPCClusters;
 457   target.fCutRequireTPCRefit = fCutRequireTPCRefit;
 458   target.fCutRequireTPCStandAlone = fCutRequireTPCStandAlone;
 459   target.fCutRequireITSRefit = fCutRequireITSRefit;
 460   target.fCutRequireITSPid = fCutRequireITSPid;
 461   target.fCutRequireITSStandAlone = fCutRequireITSStandAlone;
 462   target.fCutRequireITSpureSA = fCutRequireITSpureSA;
 463
 464   target.fCutNsigmaToVertex = fCutNsigmaToVertex;
 465   target.fCutSigmaToVertexRequired = fCutSigmaToVertexRequired;
 466   target.fCutMaxDCAToVertexXY = fCutMaxDCAToVertexXY;
 467   target.fCutMaxDCAToVertexZ = fCutMaxDCAToVertexZ;
 468   target.fCutDCAToVertex2D = fCutDCAToVertex2D;
 469   target.fCutMinDCAToVertexXY = fCutMinDCAToVertexXY;
 470   target.fCutMinDCAToVertexZ = fCutMinDCAToVertexZ;
 471
 472   target.fCutMaxDCAToVertexXYPtDep = fCutMaxDCAToVertexXYPtDep;
 473   target.SetMaxDCAToVertexXYPtDep(fCutMaxDCAToVertexXYPtDep.Data());
 474
 475   target.fCutMaxDCAToVertexZPtDep = fCutMaxDCAToVertexZPtDep;
 476   target.SetMaxDCAToVertexZPtDep(fCutMaxDCAToVertexZPtDep.Data());
 477
 478   target.fCutMinDCAToVertexXYPtDep = fCutMinDCAToVertexXYPtDep;
 479   target.SetMinDCAToVertexXYPtDep(fCutMinDCAToVertexXYPtDep.Data());
 480
 481   target.fCutMinDCAToVertexZPtDep = fCutMinDCAToVertexZPtDep;
 482   target.SetMinDCAToVertexZPtDep(fCutMinDCAToVertexZPtDep.Data());
 483
 484   target.fPMin = fPMin;
 485   target.fPMax = fPMax;
 486   target.fPtMin = fPtMin;
 487   target.fPtMax = fPtMax;
 488   target.fPxMin = fPxMin;
 489   target.fPxMax = fPxMax;
 490   target.fPyMin = fPyMin;
 491   target.fPyMax = fPyMax;
 492   target.fPzMin = fPzMin;
 493   target.fPzMax = fPzMax;
 494   target.fEtaMin = fEtaMin;
 495   target.fEtaMax = fEtaMax;
 496   target.fRapMin = fRapMin;
 497   target.fRapMax = fRapMax;
 498
 499   target.fHistogramsOn = fHistogramsOn;
 500
 501   for (Int_t i=0; i<2; ++i)
 502   {
 503     if (fhNClustersITS[i]) target.fhNClustersITS[i] = (TH1F*) fhNClustersITS[i]->Clone();
 504     if (fhNClustersTPC[i]) target.fhNClustersTPC[i] = (TH1F*) fhNClustersTPC[i]->Clone();
 505
 506     if (fhChi2PerClusterITS[i]) target.fhChi2PerClusterITS[i] = (TH1F*) fhChi2PerClusterITS[i]->Clone();
 507     if (fhChi2PerClusterTPC[i]) target.fhChi2PerClusterTPC[i] = (TH1F*) fhChi2PerClusterTPC[i]->Clone();
 508
 509     if (fhC11[i]) target.fhC11[i] = (TH1F*) fhC11[i]->Clone();
 510     if (fhC22[i]) target.fhC22[i] = (TH1F*) fhC22[i]->Clone();
 511     if (fhC33[i]) target.fhC33[i] = (TH1F*) fhC33[i]->Clone();
 512     if (fhC44[i]) target.fhC44[i] = (TH1F*) fhC44[i]->Clone();
 513     if (fhC55[i]) target.fhC55[i] = (TH1F*) fhC55[i]->Clone();
 514
 515     if (fhRel1PtUncertainty[i]) target.fhRel1PtUncertainty[i] = (TH1F*) fhRel1PtUncertainty[i]->Clone();
 516
 517     if (fhDXY[i]) target.fhDXY[i] = (TH1F*) fhDXY[i]->Clone();
 518     if (fhDZ[i]) target.fhDZ[i] = (TH1F*) fhDZ[i]->Clone();
 519     if (fhDXYDZ[i]) target.fhDXYDZ[i] = (TH1F*) fhDXYDZ[i]->Clone();
 520     if (fhDXYvsDZ[i]) target.fhDXYvsDZ[i] = (TH2F*) fhDXYvsDZ[i]->Clone();
 521
 522     if (fhDXYNormalized[i]) target.fhDXYNormalized[i] = (TH1F*) fhDXYNormalized[i]->Clone();
 523     if (fhDZNormalized[i]) target.fhDZNormalized[i] = (TH1F*) fhDZNormalized[i]->Clone();
 524     if (fhDXYvsDZNormalized[i]) target.fhDXYvsDZNormalized[i] = (TH2F*) fhDXYvsDZNormalized[i]->Clone();
 525     if (fhNSigmaToVertex[i]) target.fhNSigmaToVertex[i] = (TH1F*) fhNSigmaToVertex[i]->Clone();
 526
 527     if (fhPt[i]) target.fhPt[i] = (TH1F*) fhPt[i]->Clone();
 528     if (fhEta[i]) target.fhEta[i] = (TH1F*) fhEta[i]->Clone();
 529   }
 530   if (ffDTheoretical) target.ffDTheoretical = (TF1*) ffDTheoretical->Clone();
 531
 532   if (fhCutStatistics) target.fhCutStatistics = (TH1F*) fhCutStatistics->Clone();
 533   if (fhCutCorrelation) target.fhCutCorrelation = (TH2F*) fhCutCorrelation->Clone();
 534
 535   TNamed::Copy(c);
 536 }
 537
 538 //_____________________________________________________________________________
 539 Long64_t AliESDtrackCuts::Merge(TCollection* list) {
 540   // Merge a list of AliESDtrackCuts objects with this (needed for PROOF)
 541   // Returns the number of merged objects (including this)
 542   if (!list)
 543     return 0;
 544   if (list->IsEmpty())
 545     return 1;
 546   if (!fHistogramsOn)
 547     return 0;
 548   TIterator* iter = list->MakeIterator();
 549   TObject* obj;
 550
 551   // collection of measured and generated histograms
 552   Int_t count = 0;
 553   while ((obj = iter->Next())) {
 554
 555     AliESDtrackCuts* entry = dynamic_cast<AliESDtrackCuts*>(obj);
 556     if (entry == 0)
 557       continue;
 558
 559     if (!entry->fHistogramsOn)
 560       continue;
 561
 562     for (Int_t i=0; i<2; i++) {
 563
 564       fhNClustersITS[i]      ->Add(entry->fhNClustersITS[i]     );
 565       fhNClustersTPC[i]      ->Add(entry->fhNClustersTPC[i]     );
 566
 567       fhChi2PerClusterITS[i] ->Add(entry->fhChi2PerClusterITS[i]);
 568       fhChi2PerClusterTPC[i] ->Add(entry->fhChi2PerClusterTPC[i]);
 569
 570       fhC11[i]               ->Add(entry->fhC11[i]              );
 571       fhC22[i]               ->Add(entry->fhC22[i]              );
 572       fhC33[i]               ->Add(entry->fhC33[i]              );
 573       fhC44[i]               ->Add(entry->fhC44[i]              );
 574       fhC55[i]               ->Add(entry->fhC55[i]              );
 575
 576       fhRel1PtUncertainty[i] ->Add(entry->fhRel1PtUncertainty[i]);
 577
 578       fhDXY[i]               ->Add(entry->fhDXY[i]              );
 579       fhDZ[i]                ->Add(entry->fhDZ[i]               );
 580       fhDXYDZ[i]             ->Add(entry->fhDXYDZ[i]          );
 581       fhDXYvsDZ[i]           ->Add(entry->fhDXYvsDZ[i]          );
 582
 583       fhDXYNormalized[i]     ->Add(entry->fhDXYNormalized[i]    );
 584       fhDZNormalized[i]      ->Add(entry->fhDZNormalized[i]     );
 585       fhDXYvsDZNormalized[i] ->Add(entry->fhDXYvsDZNormalized[i]);
 586       fhNSigmaToVertex[i]    ->Add(entry->fhNSigmaToVertex[i]);
 587
 588       fhPt[i]                ->Add(entry->fhPt[i]);
 589       fhEta[i]               ->Add(entry->fhEta[i]);
 590     }
 591
 592     fhCutStatistics  ->Add(entry->fhCutStatistics);
 593     fhCutCorrelation ->Add(entry->fhCutCorrelation);
 594
 595     count++;
 596   }
 597   return count+1;
 598 }
 599
 600 //____________________________________________________________________
 601 AliESDtrackCuts* AliESDtrackCuts::GetStandardTPCOnlyTrackCuts()
 602 {
 603   // creates an AliESDtrackCuts object and fills it with standard (pre data-taking) values for TPC-only cuts
 604
 605   Printf("AliESDtrackCuts::GetStandardTPCOnlyTrackCuts: Creating track cuts for TPC-only.");
 606
 607   AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
 608
 609   esdTrackCuts->SetMinNClustersTPC(50);
 610   esdTrackCuts->SetMaxChi2PerClusterTPC(4);
 611   esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
 612
 613   esdTrackCuts->SetMaxDCAToVertexZ(3.2);
 614   esdTrackCuts->SetMaxDCAToVertexXY(2.4);
 615   esdTrackCuts->SetDCAToVertex2D(kTRUE);
 616
 617   return esdTrackCuts;
 618 }
 619
 620 //____________________________________________________________________
 621 AliESDtrackCuts* AliESDtrackCuts::GetStandardITSTPCTrackCuts2009(Bool_t selPrimaries)
 622 {
 623   // creates an AliESDtrackCuts object and fills it with standard values for ITS-TPC cuts for pp 2009 data
 624
 625   Printf("AliESDtrackCuts::GetStandardITSTPCTrackCuts: Creating track cuts for ITS+TPC.");
 626
 627   AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
 628
 629   // TPC
 630   esdTrackCuts->SetRequireTPCStandAlone(kTRUE); // to get chi2 and ncls of kTPCin
 631   esdTrackCuts->SetMinNClustersTPC(70);
 632   esdTrackCuts->SetMaxChi2PerClusterTPC(4);
 633   esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
 634   esdTrackCuts->SetRequireTPCRefit(kTRUE);
 635   // ITS
 636   esdTrackCuts->SetRequireITSRefit(kTRUE);
 637   esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
 638                                          AliESDtrackCuts::kAny);
 639   if(selPrimaries) {
 640     // 7*(0.0050+0.0060/pt^0.9)
 641     esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0350+0.0420/pt^0.9");
 642   }
 643   esdTrackCuts->SetMaxDCAToVertexZ(1.e6);
 644   esdTrackCuts->SetDCAToVertex2D(kFALSE);
 645   esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
 646   //esdTrackCuts->SetEtaRange(-0.8,+0.8);
 647
 648   return esdTrackCuts;
 649 }
 650
 651 //____________________________________________________________________
 652 AliESDtrackCuts* AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(Bool_t selPrimaries)
 653 {
 654   // creates an AliESDtrackCuts object and fills it with standard values for ITS-TPC cuts for pp 2010 data
 655
 656   Printf("AliESDtrackCuts::GetStandardITSTPCTrackCuts: Creating track cuts for ITS+TPC.");
 657
 658   AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
 659
 660   // TPC
 661   esdTrackCuts->SetMinNClustersTPC(70);
 662   esdTrackCuts->SetMaxChi2PerClusterTPC(4);
 663   esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
 664   esdTrackCuts->SetRequireTPCRefit(kTRUE);
 665   // ITS
 666   esdTrackCuts->SetRequireITSRefit(kTRUE);
 667   esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
 668                                          AliESDtrackCuts::kAny);
 669   if(selPrimaries) {
 670     // 7*(0.0026+0.0050/pt^1.01)
 671     esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0182+0.0350/pt^1.01");
 672   }
 673   esdTrackCuts->SetMaxDCAToVertexZ(2);
 674   esdTrackCuts->SetDCAToVertex2D(kFALSE);
 675   esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
 676
 677   return esdTrackCuts;
 678 }
 679
 680 //____________________________________________________________________
 681 AliESDtrackCuts* AliESDtrackCuts::GetStandardITSPureSATrackCuts2009(Bool_t selPrimaries, Bool_t useForPid)
 682 {
 683   // creates an AliESDtrackCuts object and fills it with standard values for ITS pure SA tracks
 684
 685   AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
 686   esdTrackCuts->SetRequireITSPureStandAlone(kTRUE);
 687   esdTrackCuts->SetRequireITSRefit(kTRUE);
 688   esdTrackCuts->SetMinNClustersITS(4);
 689   esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
 690                                          AliESDtrackCuts::kAny);
 691   esdTrackCuts->SetMaxChi2PerClusterITS(1.);
 692
 693   if(selPrimaries) {
 694     // 7*(0.0085+0.0026/pt^1.55)
 695     esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0595+0.0182/pt^1.55");
 696   }
 697   if(useForPid){
 698     esdTrackCuts->SetRequireITSPid(kTRUE);
 699   }
 700   return esdTrackCuts;
 701 }
 702
 703 //____________________________________________________________________
 704 AliESDtrackCuts* AliESDtrackCuts::GetStandardITSPureSATrackCuts2010(Bool_t selPrimaries, Bool_t useForPid)
 705 {
 706   // creates an AliESDtrackCuts object and fills it with standard values for ITS pure SA tracks - pp 2010
 707
 708   AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
 709   esdTrackCuts->SetRequireITSPureStandAlone(kTRUE);
 710   esdTrackCuts->SetRequireITSRefit(kTRUE);
 711   esdTrackCuts->SetMinNClustersITS(4);
 712   esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
 713                                          AliESDtrackCuts::kAny);
 714   esdTrackCuts->SetMaxChi2PerClusterITS(2.5);
 715
 716   if(selPrimaries) {
 717     // 7*(0.0033+0.0045/pt^1.3)
 718     esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0231+0.0315/pt^1.3");
 719   }
 720   if(useForPid){
 721     esdTrackCuts->SetRequireITSPid(kTRUE);
 722   }
 723   return esdTrackCuts;
 724 }
 725
 726 //____________________________________________________________________
 727 AliESDtrackCuts* AliESDtrackCuts::GetStandardITSSATrackCuts2009(Bool_t selPrimaries, Bool_t useForPid)
 728 {
 729   // creates an AliESDtrackCuts object and fills it with standard values for ITS pure SA tracks
 730
 731   AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
 732   esdTrackCuts->SetRequireITSStandAlone(kTRUE);
 733   esdTrackCuts->SetRequireITSPureStandAlone(kFALSE);
 734   esdTrackCuts->SetRequireITSRefit(kTRUE);
 735   esdTrackCuts->SetMinNClustersITS(4);
 736   esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
 737                                          AliESDtrackCuts::kAny);
 738   esdTrackCuts->SetMaxChi2PerClusterITS(1.);
 739
 740   if(selPrimaries) {
 741     // 7*(0.0085+0.0026/pt^1.55)
 742     esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0595+0.0182/pt^1.55");
 743   }
 744   if(useForPid){
 745     esdTrackCuts->SetRequireITSPid(kTRUE);
 746   }
 747   return esdTrackCuts;
 748 }
 749
 750 //____________________________________________________________________
 751 AliESDtrackCuts* AliESDtrackCuts::GetStandardITSSATrackCuts2010(Bool_t selPrimaries, Bool_t useForPid)
 752 {
 753   // creates an AliESDtrackCuts object and fills it with standard values for ITS pure SA tracks --pp 2010
 754
 755   AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
 756   esdTrackCuts->SetRequireITSStandAlone(kTRUE);
 757   esdTrackCuts->SetRequireITSPureStandAlone(kFALSE);
 758   esdTrackCuts->SetRequireITSRefit(kTRUE);
 759   esdTrackCuts->SetMinNClustersITS(4);
 760   esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
 761                                          AliESDtrackCuts::kAny);
 762   esdTrackCuts->SetMaxChi2PerClusterITS(2.5);
 763
 764   if(selPrimaries) {
 765     // 7*(0.0033+0.0045/pt^1.3)
 766     esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0231+0.0315/pt^1.3");
 767   }
 768   if(useForPid){
 769     esdTrackCuts->SetRequireITSPid(kTRUE);
 770   }
 771   return esdTrackCuts;
 772 }
 773
 774 //____________________________________________________________________
 775 Int_t AliESDtrackCuts::GetReferenceMultiplicity(AliESDEvent* esd, Bool_t tpcOnly)
 776 {
 777   // Gets reference multiplicity following the standard cuts and a defined fiducial volume
 778   // tpcOnly = kTRUE -> consider TPC-only tracks
 779   //         = kFALSE -> consider global tracks
 780
 781   if (!tpcOnly)
 782   {
 783     Printf("AliESDtrackCuts::GetReferenceMultiplicity: Not implemented for global tracks!");
 784     return -1;
 785   }
 786
 787   static AliESDtrackCuts* esdTrackCuts = 0;
 788   if (!esdTrackCuts)
 789   {
 790     esdTrackCuts = GetStandardTPCOnlyTrackCuts();
 791     esdTrackCuts->SetEtaRange(-0.8, 0.8);
 792     esdTrackCuts->SetPtRange(0.15);
 793   }
 794
 795   Int_t nTracks = esdTrackCuts->CountAcceptedTracks(esd);
 796
 797   return nTracks;
 798 }
 799
 800 //____________________________________________________________________
 801 Float_t AliESDtrackCuts::GetSigmaToVertex(AliESDtrack* const esdTrack)
 802 {
 803   // Calculates the number of sigma to the vertex.
 804
 805   Float_t b[2];
 806   Float_t bRes[2];
 807   Float_t bCov[3];
 808   esdTrack->GetImpactParameters(b,bCov);
 809
 810   if (bCov[0]<=0 || bCov[2]<=0) {
 811     AliDebugClass(1, "Estimated b resolution lower or equal zero!");
 812     bCov[0]=0; bCov[2]=0;
 813   }
 814   bRes[0] = TMath::Sqrt(bCov[0]);
 815   bRes[1] = TMath::Sqrt(bCov[2]);
 816
 817   // -----------------------------------
 818   // How to get to a n-sigma cut?
 819   //
 820   // The accumulated statistics from 0 to d is
 821   //
 822   // ->  Erf(d/Sqrt(2)) for a 1-dim gauss (d = n_sigma)
 823   // ->  1 - Exp(-d**2) for a 2-dim gauss (d*d = dx*dx + dy*dy != n_sigma)
 824   //
 825   // It means that for a 2-dim gauss: n_sigma(d) = Sqrt(2)*ErfInv(1 - Exp((-d**2)/2)
 826   // Can this be expressed in a different way?
 827
 828   if (bRes[0] == 0 || bRes[1] ==0)
 829     return -1;
 830
 831   Float_t d = TMath::Sqrt(TMath::Power(b[0]/bRes[0],2) + TMath::Power(b[1]/bRes[1],2));
 832
 833   // work around precision problem
 834   // if d is too big, TMath::Exp(...) gets 0, and TMath::ErfInverse(1) that should be infinite, gets 0 :(
 835   // 1e-15 corresponds to nsigma ~ 7.7
 836   if (TMath::Exp(-d * d / 2) < 1e-15)
 837     return 1000;
 838
 839   Float_t nSigma = TMath::ErfInverse(1 - TMath::Exp(-d * d / 2)) * TMath::Sqrt(2);
 840   return nSigma;
 841 }
 842
 843 void AliESDtrackCuts::EnableNeededBranches(TTree* tree)
 844 {
 845   // enables the branches needed by AcceptTrack, for a list see comment of AcceptTrack
 846
 847   tree->SetBranchStatus("fTracks.fFlags", 1);
 848   tree->SetBranchStatus("fTracks.fITSncls", 1);
 849   tree->SetBranchStatus("fTracks.fTPCncls", 1);
 850   tree->SetBranchStatus("fTracks.fITSchi2", 1);
 851   tree->SetBranchStatus("fTracks.fTPCchi2", 1);
 852   tree->SetBranchStatus("fTracks.fC*", 1);
 853   tree->SetBranchStatus("fTracks.fD", 1);
 854   tree->SetBranchStatus("fTracks.fZ", 1);
 855   tree->SetBranchStatus("fTracks.fCdd", 1);
 856   tree->SetBranchStatus("fTracks.fCdz", 1);
 857   tree->SetBranchStatus("fTracks.fCzz", 1);
 858   tree->SetBranchStatus("fTracks.fP*", 1);
 859   tree->SetBranchStatus("fTracks.fR*", 1);
 860   tree->SetBranchStatus("fTracks.fKinkIndexes*", 1);
 861 }
 862
 863 //____________________________________________________________________
 864 Bool_t AliESDtrackCuts::AcceptTrack(AliESDtrack* esdTrack)
 865 {
 866   //
 867   // figure out if the tracks survives all the track cuts defined
 868   //
 869   // the different quality parameter and kinematic values are first
 870   // retrieved from the track. then it is found out what cuts the
 871   // track did not survive and finally the cuts are imposed.
 872
 873   // this function needs the following branches:
 874   // fTracks.fFlags
 875   // fTracks.fITSncls
 876   // fTracks.fTPCncls
 877   // fTracks.fITSchi2
 878   // fTracks.fTPCchi2
 879   // fTracks.fC   //GetExternalCovariance
 880   // fTracks.fD   //GetImpactParameters
 881   // fTracks.fZ   //GetImpactParameters
 882   // fTracks.fCdd //GetImpactParameters
 883   // fTracks.fCdz //GetImpactParameters
 884   // fTracks.fCzz //GetImpactParameters
 885   // fTracks.fP   //GetPxPyPz
 886   // fTracks.fR   //GetMass
 887   // fTracks.fP   //GetMass
 888   // fTracks.fKinkIndexes
 889
 890   UInt_t status = esdTrack->GetStatus();
 891
 892   // getting quality parameters from the ESD track
 893   Int_t nClustersITS = esdTrack->GetITSclusters(0);
 894   Int_t nClustersTPC = -1;
 895   if(fCutRequireTPCStandAlone) {
 896     nClustersTPC = esdTrack->GetTPCNclsIter1();
 897   }
 898   else {
 899     nClustersTPC = esdTrack->GetTPCclusters(0);
 900   }
 901
 902   Int_t nClustersTPCShared = esdTrack->GetTPCnclsS();
 903   Float_t fracClustersTPCShared = -1.;
 904
 905   Float_t chi2PerClusterITS = -1;
 906   Float_t chi2PerClusterTPC = -1;
 907   if (nClustersITS!=0)
 908     chi2PerClusterITS = esdTrack->GetITSchi2()/Float_t(nClustersITS);
 909   if (nClustersTPC!=0) {
 910     if(fCutRequireTPCStandAlone) {
 911       chi2PerClusterTPC = esdTrack->GetTPCchi2Iter1()/Float_t(nClustersTPC);
 912     } else {
 913       chi2PerClusterTPC = esdTrack->GetTPCchi2()/Float_t(nClustersTPC);
 914     }
 915     fracClustersTPCShared = Float_t(nClustersTPCShared)/Float_t(nClustersTPC);
 916   }
 917
 918   Double_t extCov[15];
 919   esdTrack->GetExternalCovariance(extCov);
 920
 921   // getting the track to vertex parameters
 922   Float_t nSigmaToVertex = GetSigmaToVertex(esdTrack);
 923
 924   Float_t b[2];
 925   Float_t bCov[3];
 926   esdTrack->GetImpactParameters(b,bCov);
 927   if (bCov[0]<=0 || bCov[2]<=0) {
 928     AliDebug(1, "Estimated b resolution lower or equal zero!");
 929     bCov[0]=0; bCov[2]=0;
 930   }
 931
 932
 933   // set pt-dependent DCA cuts, if requested
 934   SetPtDepDCACuts(esdTrack->Pt());
 935
 936
 937   Float_t dcaToVertexXY = b[0];
 938   Float_t dcaToVertexZ = b[1];
 939
 940   Float_t dcaToVertex = -1;
 941
 942   if (fCutDCAToVertex2D)
 943   {
 944     dcaToVertex = TMath::Sqrt(dcaToVertexXY*dcaToVertexXY/fCutMaxDCAToVertexXY/fCutMaxDCAToVertexXY + dcaToVertexZ*dcaToVertexZ/fCutMaxDCAToVertexZ/fCutMaxDCAToVertexZ);
 945   }
 946   else
 947     dcaToVertex = TMath::Sqrt(dcaToVertexXY*dcaToVertexXY + dcaToVertexZ*dcaToVertexZ);
 948
 949   // getting the kinematic variables of the track
 950   // (assuming the mass is known)
 951   Double_t p[3];
 952   esdTrack->GetPxPyPz(p);
 953
 954   Float_t momentum = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2) + TMath::Power(p[2],2));
 955   Float_t pt       = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2));
 956   Float_t energy   = TMath::Sqrt(TMath::Power(esdTrack->GetMass(),2) + TMath::Power(momentum,2));
 957
 958   //y-eta related calculations
 959   Float_t eta = -100.;
 960   Float_t y   = -100.;
 961   if((momentum != TMath::Abs(p[2]))&&(momentum != 0))
 962     eta = 0.5*TMath::Log((momentum + p[2])/(momentum - p[2]));
 963   if((energy != TMath::Abs(p[2]))&&(momentum != 0))
 964     y = 0.5*TMath::Log((energy + p[2])/(energy - p[2]));
 965
 966   if (extCov[14] < 0)
 967   {
 968     Printf("AliESDtrackCuts::AcceptTrack: WARNING: GetSigma1Pt2() returns negative value for external covariance matrix element fC[14]: %f. Corrupted track information, track will not be accepted!", extCov[14]);
 969     return kFALSE;
 970   }
 971   Float_t relUncertainty1Pt = TMath::Sqrt(extCov[14])*pt;
 972
 973   //########################################################################
 974   // cut the track?
 975
 976   Bool_t cuts[kNCuts];
 977   for (Int_t i=0; i<kNCuts; i++) cuts[i]=kFALSE;
 978
 979   // track quality cuts
 980   if (fCutRequireTPCRefit && (status&AliESDtrack::kTPCrefit)==0)
 981     cuts[0]=kTRUE;
 982   if (fCutRequireTPCStandAlone && (status&AliESDtrack::kTPCin)==0)
 983     cuts[1]=kTRUE;
 984   if (fCutRequireITSRefit && (status&AliESDtrack::kITSrefit)==0)
 985     cuts[2]=kTRUE;
 986   if (nClustersTPC<fCutMinNClusterTPC)
 987     cuts[3]=kTRUE;
 988   if (nClustersITS<fCutMinNClusterITS)
 989     cuts[4]=kTRUE;
 990   if (chi2PerClusterTPC>fCutMaxChi2PerClusterTPC)
 991     cuts[5]=kTRUE;
 992   if (chi2PerClusterITS>fCutMaxChi2PerClusterITS)
 993     cuts[6]=kTRUE;
 994   if (extCov[0]  > fCutMaxC11)
 995     cuts[7]=kTRUE;
 996   if (extCov[2]  > fCutMaxC22)
 997     cuts[8]=kTRUE;
 998   if (extCov[5]  > fCutMaxC33)
 999     cuts[9]=kTRUE;
1000   if (extCov[9]  > fCutMaxC44)
1001     cuts[10]=kTRUE;
1002   if (extCov[14]  > fCutMaxC55)
1003     cuts[11]=kTRUE;
1004   if (nSigmaToVertex > fCutNsigmaToVertex && fCutSigmaToVertexRequired)
1005     cuts[12] = kTRUE;
1006   // if n sigma could not be calculated
1007   if (nSigmaToVertex<0 && fCutSigmaToVertexRequired)
1008     cuts[13]=kTRUE;
1009   if (!fCutAcceptKinkDaughters && esdTrack->GetKinkIndex(0)>0)
1010     cuts[14]=kTRUE;
1011   // track kinematics cut
1012   if((momentum < fPMin) || (momentum > fPMax))
1013     cuts[15]=kTRUE;
1014   if((pt < fPtMin) || (pt > fPtMax))
1015     cuts[16] = kTRUE;
1016   if((p[0] < fPxMin) || (p[0] > fPxMax))
1017     cuts[17] = kTRUE;
1018   if((p[1] < fPyMin) || (p[1] > fPyMax))
1019     cuts[18] = kTRUE;
1020   if((p[2] < fPzMin) || (p[2] > fPzMax))
1021     cuts[19] = kTRUE;
1022   if((eta < fEtaMin) || (eta > fEtaMax))
1023     cuts[20] = kTRUE;
1024   if((y < fRapMin) || (y > fRapMax))
1025     cuts[21] = kTRUE;
1026   if (fCutDCAToVertex2D && dcaToVertex > 1)
1027     cuts[22] = kTRUE;
1028   if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexXY) > fCutMaxDCAToVertexXY)
1029     cuts[23] = kTRUE;
1030   if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexZ) > fCutMaxDCAToVertexZ)
1031     cuts[24] = kTRUE;
1032   if (fCutDCAToVertex2D && fCutMinDCAToVertexXY > 0 && fCutMinDCAToVertexZ > 0 && dcaToVertexXY*dcaToVertexXY/fCutMinDCAToVertexXY/fCutMinDCAToVertexXY + dcaToVertexZ*dcaToVertexZ/fCutMinDCAToVertexZ/fCutMinDCAToVertexZ < 1)
1033     cuts[25] = kTRUE;
1034   if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexXY) < fCutMinDCAToVertexXY)
1035     cuts[26] = kTRUE;
1036   if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexZ) < fCutMinDCAToVertexZ)
1037     cuts[27] = kTRUE;
1038
1039   for (Int_t i = 0; i < 3; i++)
1040     cuts[28+i] = !CheckITSClusterRequirement(fCutClusterRequirementITS[i], esdTrack->HasPointOnITSLayer(i*2), esdTrack->HasPointOnITSLayer(i*2+1));
1041
1042   if(fCutRequireITSStandAlone || fCutRequireITSpureSA){
1043     if ((status & AliESDtrack::kITSin) == 0 || (status & AliESDtrack::kTPCin)){
1044       // TPC tracks
1045       cuts[31] = kTRUE;
1046     }else{
1047       // ITS standalone tracks
1048       if(fCutRequireITSStandAlone && !fCutRequireITSpureSA){
1049         if(status & AliESDtrack::kITSpureSA) cuts[31] = kTRUE;
1050       }else if(fCutRequireITSpureSA){
1051         if(!(status & AliESDtrack::kITSpureSA)) cuts[31] = kTRUE;
1052       }
1053     }
1054   }
1055
1056   if (relUncertainty1Pt > fCutMaxRel1PtUncertainty)
1057      cuts[32] = kTRUE;
1058
1059   if (!fCutAcceptSharedTPCClusters && nClustersTPCShared!=0)
1060     cuts[33] = kTRUE;
1061
1062   if (fracClustersTPCShared > fCutMaxFractionSharedTPCClusters)
1063     cuts[34] = kTRUE;
1064
1065   if(fCutRequireITSPid){
1066     UChar_t clumap=esdTrack->GetITSClusterMap();
1067     Int_t nPointsForPid=0;
1068     for(Int_t i=2; i<6; i++){
1069       if(clumap&(1<<i)) ++nPointsForPid;
1070     }
1071     if(nPointsForPid<3) cuts[35] = kTRUE;
1072   }
1073
1074   Bool_t cut=kFALSE;
1075   for (Int_t i=0; i<kNCuts; i++)
1076     if (cuts[i]) {cut = kTRUE;}
1077
1078
1079   //########################################################################
1080   // filling histograms
1081   if (fHistogramsOn) {
1082     fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n tracks")));
1083     if (cut)
1084       fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n cut tracks")));
1085
1086     for (Int_t i=0; i<kNCuts; i++) {
1087       if (fhCutStatistics->GetXaxis()->FindBin(fgkCutNames[i]) < 1)
1088         AliFatal(Form("Inconsistency! Cut %d with name %s not found", i, fgkCutNames[i]));
1089
1090       if (cuts[i])
1091         fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin(fgkCutNames[i])));
1092
1093       for (Int_t j=i; j<kNCuts; j++) {
1094         if (cuts[i] && cuts[j]) {
1095           Float_t xC = fhCutCorrelation->GetXaxis()->GetBinCenter(fhCutCorrelation->GetXaxis()->FindBin(fgkCutNames[i]));
1096           Float_t yC = fhCutCorrelation->GetYaxis()->GetBinCenter(fhCutCorrelation->GetYaxis()->FindBin(fgkCutNames[j]));
1097           fhCutCorrelation->Fill(xC, yC);
1098         }
1099       }
1100     }
1101   }
1102
1103   // now we loop over the filling of the histograms twice: once "before" the cut, once "after"
1104   // the code is not in a function due to too many local variables that would need to be passed
1105
1106   for (Int_t id = 0; id < 2; id++)
1107   {
1108     // id = 0 --> before cut
1109     // id = 1 --> after cut
1110
1111     if (fHistogramsOn)
1112     {
1113       fhNClustersITS[id]->Fill(nClustersITS);
1114       fhNClustersTPC[id]->Fill(nClustersTPC);
1115       fhChi2PerClusterITS[id]->Fill(chi2PerClusterITS);
1116       fhChi2PerClusterTPC[id]->Fill(chi2PerClusterTPC);
1117
1118       fhC11[id]->Fill(extCov[0]);
1119       fhC22[id]->Fill(extCov[2]);
1120       fhC33[id]->Fill(extCov[5]);
1121       fhC44[id]->Fill(extCov[9]);
1122       fhC55[id]->Fill(extCov[14]);
1123
1124       fhRel1PtUncertainty[id]->Fill(relUncertainty1Pt);
1125
1126       fhPt[id]->Fill(pt);
1127       fhEta[id]->Fill(eta);
1128
1129       Float_t bRes[2];
1130       bRes[0] = TMath::Sqrt(bCov[0]);
1131       bRes[1] = TMath::Sqrt(bCov[2]);
1132
1133       fhDZ[id]->Fill(b[1]);
1134       fhDXY[id]->Fill(b[0]);
1135       fhDXYDZ[id]->Fill(dcaToVertex);
1136       fhDXYvsDZ[id]->Fill(b[1],b[0]);
1137
1138       if (bRes[0]!=0 && bRes[1]!=0) {
1139         fhDZNormalized[id]->Fill(b[1]/bRes[1]);
1140         fhDXYNormalized[id]->Fill(b[0]/bRes[0]);
1141         fhDXYvsDZNormalized[id]->Fill(b[1]/bRes[1], b[0]/bRes[0]);
1142         fhNSigmaToVertex[id]->Fill(nSigmaToVertex);
1143       }
1144     }
1145
1146     // cut the track
1147     if (cut)
1148       return kFALSE;
1149   }
1150
1151   return kTRUE;
1152 }
1153
1154 //____________________________________________________________________
1155 Bool_t AliESDtrackCuts::CheckITSClusterRequirement(ITSClusterRequirement req, Bool_t clusterL1, Bool_t clusterL2)
1156 {
1157   // checks if the cluster requirement is fullfilled (in this case: return kTRUE)
1158
1159   switch (req)
1160   {
1161         case kOff:        return kTRUE;
1162         case kNone:       return !clusterL1 && !clusterL2;
1163         case kAny:        return clusterL1 || clusterL2;
1164         case kFirst:      return clusterL1;
1165         case kOnlyFirst:  return clusterL1 && !clusterL2;
1166         case kSecond:     return clusterL2;
1167         case kOnlySecond: return clusterL2 && !clusterL1;
1168         case kBoth:       return clusterL1 && clusterL2;
1169   }
1170
1171   return kFALSE;
1172 }
1173
1174 //____________________________________________________________________
1175 AliESDtrack* AliESDtrackCuts::GetTPCOnlyTrack(AliESDEvent* esd, Int_t iTrack)
1176 {
1177
1178   // Utility function to
1179   // create a TPC only track from the given esd track
1180   //
1181   // IMPORTANT: The track has to be deleted by the user
1182   //
1183   // NB. most of the functionality to get a TPC only track from an ESD track is in AliESDtrack, where it should be
1184   // there are only missing propagations here that are needed for old data
1185   // this function will therefore become obsolete
1186   //
1187   // adapted from code provided by CKB
1188
1189   if (!esd->GetPrimaryVertexTPC())
1190     return 0; // No TPC vertex no TPC tracks
1191
1192   if(!esd->GetPrimaryVertexTPC()->GetStatus())
1193     return 0; // TPC Vertex is created by default in AliESDEvent, do not use in this case
1194
1195   AliESDtrack* track = esd->GetTrack(iTrack);
1196   if (!track)
1197     return 0;
1198
1199   AliESDtrack *tpcTrack = new AliESDtrack();
1200
1201   // only true if we have a tpc track
1202   if (!track->FillTPCOnlyTrack(*tpcTrack))
1203   {
1204     delete tpcTrack;
1205     return 0;
1206   }
1207
1208   // propagate to Vertex
1209   // not needed for normal reconstructed ESDs...
1210   // Double_t pTPC[2],covTPC[3];
1211   // tpcTrack->PropagateToDCA(esd->GetPrimaryVertexTPC(), esd->GetMagneticField(), 10000,  pTPC, covTPC);
1212
1213   return tpcTrack;
1214 }
1215
1216 //____________________________________________________________________
1217 TObjArray* AliESDtrackCuts::GetAcceptedTracks(AliESDEvent* esd,Bool_t bTPC)
1218 {
1219   //
1220   // returns an array of all tracks that pass the cuts
1221   // or an array of TPC only tracks (propagated to the TPC vertex during reco)
1222   // tracks that pass the cut
1223   //
1224   // NOTE: List has to be deleted by the user
1225
1226   TObjArray* acceptedTracks = new TObjArray();
1227
1228   // loop over esd tracks
1229   for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
1230     if(bTPC){
1231       if(!esd->GetPrimaryVertexTPC())return acceptedTracks; // No TPC vertex no TPC tracks
1232       if(!esd->GetPrimaryVertexTPC()->GetStatus())return acceptedTracks; // No proper TPC vertex, only the default
1233
1234       AliESDtrack *tpcTrack = GetTPCOnlyTrack(esd, iTrack);
1235       if (!tpcTrack)
1236         continue;
1237
1238       if (AcceptTrack(tpcTrack)) {
1239         acceptedTracks->Add(tpcTrack);
1240       }
1241       else
1242         delete tpcTrack;
1243     }
1244     else
1245     {
1246       AliESDtrack* track = esd->GetTrack(iTrack);
1247       if(AcceptTrack(track))
1248         acceptedTracks->Add(track);
1249     }
1250   }
1251   if(bTPC)acceptedTracks->SetOwner(kTRUE);
1252   return acceptedTracks;
1253 }
1254
1255 //____________________________________________________________________
1256 Int_t AliESDtrackCuts::CountAcceptedTracks(AliESDEvent* const esd)
1257 {
1258   //
1259   // returns an the number of tracks that pass the cuts
1260   //
1261
1262   Int_t count = 0;
1263
1264   // loop over esd tracks
1265   for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
1266     AliESDtrack* track = esd->GetTrack(iTrack);
1267     if (AcceptTrack(track))
1268       count++;
1269   }
1270
1271   return count;
1272 }
1273
1274 //____________________________________________________________________
1275  void AliESDtrackCuts::DefineHistograms(Int_t color) {
1276    //
1277    // diagnostics histograms are defined
1278    //
1279
1280    fHistogramsOn=kTRUE;
1281
1282    Bool_t oldStatus = TH1::AddDirectoryStatus();
1283    TH1::AddDirectory(kFALSE);
1284
1285    //###################################################################################
1286    // defining histograms
1287
1288    fhCutStatistics = new TH1F("cut_statistics","cut statistics",kNCuts+4,-0.5,kNCuts+3.5);
1289
1290    fhCutStatistics->GetXaxis()->SetBinLabel(1,"n tracks");
1291    fhCutStatistics->GetXaxis()->SetBinLabel(2,"n cut tracks");
1292
1293    fhCutCorrelation = new TH2F("cut_correlation","cut correlation",kNCuts,-0.5,kNCuts-0.5,kNCuts,-0.5,kNCuts-0.5);;
1294
1295    for (Int_t i=0; i<kNCuts; i++) {
1296      fhCutStatistics->GetXaxis()->SetBinLabel(i+4,fgkCutNames[i]);
1297      fhCutCorrelation->GetXaxis()->SetBinLabel(i+1,fgkCutNames[i]);
1298      fhCutCorrelation->GetYaxis()->SetBinLabel(i+1,fgkCutNames[i]);
1299    }
1300
1301   fhCutStatistics  ->SetLineColor(color);
1302   fhCutCorrelation ->SetLineColor(color);
1303   fhCutStatistics  ->SetLineWidth(2);
1304   fhCutCorrelation ->SetLineWidth(2);
1305
1306   for (Int_t i=0; i<2; i++) {
1307     fhNClustersITS[i]        = new TH1F("nClustersITS"    ,"",8,-0.5,7.5);
1308     fhNClustersTPC[i]        = new TH1F("nClustersTPC"     ,"",165,-0.5,164.5);
1309     fhChi2PerClusterITS[i]   = new TH1F("chi2PerClusterITS","",500,0,10);
1310     fhChi2PerClusterTPC[i]   = new TH1F("chi2PerClusterTPC","",500,0,10);
1311
1312     fhC11[i]                 = new TH1F("covMatrixDiagonal11","",2000,0,20);
1313     fhC22[i]                 = new TH1F("covMatrixDiagonal22","",2000,0,20);
1314     fhC33[i]                 = new TH1F("covMatrixDiagonal33","",1000,0,0.1);
1315     fhC44[i]                 = new TH1F("covMatrixDiagonal44","",1000,0,0.1);
1316     fhC55[i]                 = new TH1F("covMatrixDiagonal55","",1000,0,5);
1317
1318     fhRel1PtUncertainty[i]   = new TH1F("rel1PtUncertainty","",1000,0,5);
1319
1320     fhDXY[i]                 = new TH1F("dXY"    ,"",500,-10,10);
1321     fhDZ[i]                  = new TH1F("dZ"     ,"",500,-10,10);
1322     fhDXYDZ[i]               = new TH1F("dXYDZ"  ,"",500,0,10);
1323     fhDXYvsDZ[i]             = new TH2F("dXYvsDZ","",200,-10,10,200,-10,10);
1324
1325     fhDXYNormalized[i]       = new TH1F("dXYNormalized"    ,"",500,-10,10);
1326     fhDZNormalized[i]        = new TH1F("dZNormalized"     ,"",500,-10,10);
1327     fhDXYvsDZNormalized[i]   = new TH2F("dXYvsDZNormalized","",200,-10,10,200,-10,10);
1328
1329     fhNSigmaToVertex[i]      = new TH1F("nSigmaToVertex","",500,0,10);
1330
1331     fhPt[i]                  = new TH1F("pt"     ,"p_{T} distribution;p_{T} (GeV/c)", 800, 0.0, 10.0);
1332     fhEta[i]                 = new TH1F("eta"     ,"#eta distribution;#eta",40,-2.0,2.0);
1333
1334     fhNClustersITS[i]->SetTitle("n ITS clusters");
1335     fhNClustersTPC[i]->SetTitle("n TPC clusters");
1336     fhChi2PerClusterITS[i]->SetTitle("#Chi^{2} per ITS cluster");
1337     fhChi2PerClusterTPC[i]->SetTitle("#Chi^{2} per TPC cluster");
1338
1339     fhC11[i]->SetTitle("cov 11 : #sigma_{y}^{2} [cm^{2}]");
1340     fhC22[i]->SetTitle("cov 22 : #sigma_{z}^{2} [cm^{2}]");
1341     fhC33[i]->SetTitle("cov 33 : #sigma_{sin(#phi)}^{2}");
1342     fhC44[i]->SetTitle("cov 44 : #sigma_{tan(#theta_{dip})}^{2}");
1343     fhC55[i]->SetTitle("cov 55 : #sigma_{1/p_{T}}^{2} [(c/GeV)^2]");
1344
1345     fhRel1PtUncertainty[i]->SetTitle("rel. uncertainty of 1/p_{T}");
1346
1347     fhDXY[i]->SetXTitle("transverse impact parameter (cm)");
1348     fhDZ[i]->SetXTitle("longitudinal impact parameter (cm)");
1349     fhDXYDZ[i]->SetTitle("absolute impact parameter;sqrt(dXY**2 + dZ**2)  (cm)");
1350     fhDXYvsDZ[i]->SetXTitle("longitudinal impact parameter (cm)");
1351     fhDXYvsDZ[i]->SetYTitle("transverse impact parameter (cm)");
1352
1353     fhDXYNormalized[i]->SetTitle("normalized trans impact par (n#sigma)");
1354     fhDZNormalized[i]->SetTitle("normalized long impact par (n#sigma)");
1355     fhDXYvsDZNormalized[i]->SetTitle("normalized long impact par (n#sigma)");
1356     fhDXYvsDZNormalized[i]->SetYTitle("normalized trans impact par (n#sigma)");
1357     fhNSigmaToVertex[i]->SetTitle("n #sigma to vertex");
1358
1359     fhNClustersITS[i]->SetLineColor(color);   fhNClustersITS[i]->SetLineWidth(2);
1360     fhNClustersTPC[i]->SetLineColor(color);   fhNClustersTPC[i]->SetLineWidth(2);
1361     fhChi2PerClusterITS[i]->SetLineColor(color);   fhChi2PerClusterITS[i]->SetLineWidth(2);
1362     fhChi2PerClusterTPC[i]->SetLineColor(color);   fhChi2PerClusterTPC[i]->SetLineWidth(2);
1363
1364     fhC11[i]->SetLineColor(color);   fhC11[i]->SetLineWidth(2);
1365     fhC22[i]->SetLineColor(color);   fhC22[i]->SetLineWidth(2);
1366     fhC33[i]->SetLineColor(color);   fhC33[i]->SetLineWidth(2);
1367     fhC44[i]->SetLineColor(color);   fhC44[i]->SetLineWidth(2);
1368     fhC55[i]->SetLineColor(color);   fhC55[i]->SetLineWidth(2);
1369
1370     fhRel1PtUncertainty[i]->SetLineColor(color); fhRel1PtUncertainty[i]->SetLineWidth(2);
1371
1372     fhDXY[i]->SetLineColor(color);   fhDXY[i]->SetLineWidth(2);
1373     fhDZ[i]->SetLineColor(color);    fhDZ[i]->SetLineWidth(2);
1374     fhDXYDZ[i]->SetLineColor(color); fhDXYDZ[i]->SetLineWidth(2);
1375
1376     fhDXYNormalized[i]->SetLineColor(color);   fhDXYNormalized[i]->SetLineWidth(2);
1377     fhDZNormalized[i]->SetLineColor(color);    fhDZNormalized[i]->SetLineWidth(2);
1378     fhNSigmaToVertex[i]->SetLineColor(color);  fhNSigmaToVertex[i]->SetLineWidth(2);
1379   }
1380
1381   // The number of sigmas to the vertex is per definition gaussian
1382   ffDTheoretical = new TF1("nSigmaToVertexTheoretical","([0]/2.506628274)*exp(-(x**2)/2)",0,50);
1383   ffDTheoretical->SetParameter(0,1);
1384
1385   TH1::AddDirectory(oldStatus);
1386 }
1387
1388 //____________________________________________________________________
1389 Bool_t AliESDtrackCuts::LoadHistograms(const Char_t* dir)
1390 {
1391   //
1392   // loads the histograms from a file
1393   // if dir is empty a directory with the name of this object is taken (like in SaveHistogram)
1394   //
1395
1396   if (!dir)
1397     dir = GetName();
1398
1399   if (!gDirectory->cd(dir))
1400     return kFALSE;
1401
1402   ffDTheoretical = dynamic_cast<TF1*> (gDirectory->Get("nSigmaToVertexTheory"));
1403
1404   fhCutStatistics = dynamic_cast<TH1F*> (gDirectory->Get("cut_statistics"));
1405   fhCutCorrelation = dynamic_cast<TH2F*> (gDirectory->Get("cut_correlation"));
1406
1407   for (Int_t i=0; i<2; i++) {
1408     if (i==0)
1409     {
1410       gDirectory->cd("before_cuts");
1411     }
1412     else
1413       gDirectory->cd("after_cuts");
1414
1415     fhNClustersITS[i]      = dynamic_cast<TH1F*> (gDirectory->Get("nClustersITS"     ));
1416     fhNClustersTPC[i]      = dynamic_cast<TH1F*> (gDirectory->Get("nClustersTPC"     ));
1417     fhChi2PerClusterITS[i] = dynamic_cast<TH1F*> (gDirectory->Get("chi2PerClusterITS"));
1418     fhChi2PerClusterTPC[i] = dynamic_cast<TH1F*> (gDirectory->Get("chi2PerClusterTPC"));
1419
1420     fhC11[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal11"));
1421     fhC22[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal22"));
1422     fhC33[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal33"));
1423     fhC44[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal44"));
1424     fhC55[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal55"));
1425
1426     fhRel1PtUncertainty[i] = dynamic_cast<TH1F*> (gDirectory->Get("rel1PtUncertainty"));
1427
1428     fhDXY[i] =     dynamic_cast<TH1F*> (gDirectory->Get("dXY"    ));
1429     fhDZ[i] =      dynamic_cast<TH1F*> (gDirectory->Get("dZ"     ));
1430     fhDXYDZ[i] =   dynamic_cast<TH1F*> (gDirectory->Get("dXYDZ"));
1431     fhDXYvsDZ[i] = dynamic_cast<TH2F*> (gDirectory->Get("dXYvsDZ"));
1432
1433     fhDXYNormalized[i] =     dynamic_cast<TH1F*> (gDirectory->Get("dXYNormalized"    ));
1434     fhDZNormalized[i] =      dynamic_cast<TH1F*> (gDirectory->Get("dZNormalized"     ));
1435     fhDXYvsDZNormalized[i] = dynamic_cast<TH2F*> (gDirectory->Get("dXYvsDZNormalized"));
1436     fhNSigmaToVertex[i] = dynamic_cast<TH1F*> (gDirectory->Get("nSigmaToVertex"));
1437
1438     fhPt[i] = dynamic_cast<TH1F*> (gDirectory->Get("pt"));
1439     fhEta[i] = dynamic_cast<TH1F*> (gDirectory->Get("eta"));
1440
1441     gDirectory->cd("../");
1442   }
1443
1444   gDirectory->cd("..");
1445
1446   return kTRUE;
1447 }
1448
1449 //____________________________________________________________________
1450 void AliESDtrackCuts::SaveHistograms(const Char_t* dir) {
1451   //
1452   // saves the histograms in a directory (dir)
1453   //
1454
1455   if (!fHistogramsOn) {
1456     AliDebug(0, "Histograms not on - cannot save histograms!!!");
1457     return;
1458   }
1459
1460   if (!dir)
1461     dir = GetName();
1462
1463   gDirectory->mkdir(dir);
1464   gDirectory->cd(dir);
1465
1466   gDirectory->mkdir("before_cuts");
1467   gDirectory->mkdir("after_cuts");
1468
1469   // a factor of 2 is needed since n sigma is positive
1470   ffDTheoretical->SetParameter(0,2*fhNSigmaToVertex[0]->Integral("width"));
1471   ffDTheoretical->Write("nSigmaToVertexTheory");
1472
1473   fhCutStatistics->Write();
1474   fhCutCorrelation->Write();
1475
1476   for (Int_t i=0; i<2; i++) {
1477     if (i==0)
1478       gDirectory->cd("before_cuts");
1479     else
1480       gDirectory->cd("after_cuts");
1481
1482     fhNClustersITS[i]        ->Write();
1483     fhNClustersTPC[i]        ->Write();
1484     fhChi2PerClusterITS[i]   ->Write();
1485     fhChi2PerClusterTPC[i]   ->Write();
1486
1487     fhC11[i]                 ->Write();
1488     fhC22[i]                 ->Write();
1489     fhC33[i]                 ->Write();
1490     fhC44[i]                 ->Write();
1491     fhC55[i]                 ->Write();
1492
1493     fhRel1PtUncertainty[i]   ->Write();
1494
1495     fhDXY[i]                 ->Write();
1496     fhDZ[i]                  ->Write();
1497     fhDXYDZ[i]               ->Write();
1498     fhDXYvsDZ[i]             ->Write();
1499
1500     fhDXYNormalized[i]       ->Write();
1501     fhDZNormalized[i]        ->Write();
1502     fhDXYvsDZNormalized[i]   ->Write();
1503     fhNSigmaToVertex[i]      ->Write();
1504
1505     fhPt[i]                  ->Write();
1506     fhEta[i]                 ->Write();
1507
1508     gDirectory->cd("../");
1509   }
1510
1511   gDirectory->cd("../");
1512 }
1513
1514 //____________________________________________________________________
1515 void AliESDtrackCuts::DrawHistograms()
1516 {
1517   // draws some histograms
1518
1519   TCanvas* canvas1 = new TCanvas(Form("%s_1", GetName()), "Track Quality Results1", 800, 800);
1520   canvas1->Divide(2, 2);
1521
1522   canvas1->cd(1);
1523   fhNClustersTPC[0]->SetStats(kFALSE);
1524   fhNClustersTPC[0]->Draw();
1525
1526   canvas1->cd(2);
1527   fhChi2PerClusterTPC[0]->SetStats(kFALSE);
1528   fhChi2PerClusterTPC[0]->Draw();
1529
1530   canvas1->cd(3);
1531   fhNSigmaToVertex[0]->SetStats(kFALSE);
1532   fhNSigmaToVertex[0]->GetXaxis()->SetRangeUser(0, 10);
1533   fhNSigmaToVertex[0]->Draw();
1534
1535   canvas1->SaveAs(Form("%s_%s.gif", GetName(), canvas1->GetName()));
1536
1537   TCanvas* canvas2 = new TCanvas(Form("%s_2", GetName()), "Track Quality Results2", 1200, 800);
1538   canvas2->Divide(3, 2);
1539
1540   canvas2->cd(1);
1541   fhC11[0]->SetStats(kFALSE);
1542   gPad->SetLogy();
1543   fhC11[0]->Draw();
1544
1545   canvas2->cd(2);
1546   fhC22[0]->SetStats(kFALSE);
1547   gPad->SetLogy();
1548   fhC22[0]->Draw();
1549
1550   canvas2->cd(3);
1551   fhC33[0]->SetStats(kFALSE);
1552   gPad->SetLogy();
1553   fhC33[0]->Draw();
1554
1555   canvas2->cd(4);
1556   fhC44[0]->SetStats(kFALSE);
1557   gPad->SetLogy();
1558   fhC44[0]->Draw();
1559
1560   canvas2->cd(5);
1561   fhC55[0]->SetStats(kFALSE);
1562   gPad->SetLogy();
1563   fhC55[0]->Draw();
1564
1565   canvas2->cd(6);
1566   fhRel1PtUncertainty[0]->SetStats(kFALSE);
1567   gPad->SetLogy();
1568   fhRel1PtUncertainty[0]->Draw();
1569
1570   canvas2->SaveAs(Form("%s_%s.gif", GetName(), canvas2->GetName()));
1571
1572   TCanvas* canvas3 = new TCanvas(Form("%s_3", GetName()), "Track Quality Results3", 1200, 800);
1573   canvas3->Divide(3, 2);
1574
1575   canvas3->cd(1);
1576   fhDXY[0]->SetStats(kFALSE);
1577   gPad->SetLogy();
1578   fhDXY[0]->Draw();
1579
1580   canvas3->cd(2);
1581   fhDZ[0]->SetStats(kFALSE);
1582   gPad->SetLogy();
1583   fhDZ[0]->Draw();
1584
1585   canvas3->cd(3);
1586   fhDXYvsDZ[0]->SetStats(kFALSE);
1587   gPad->SetLogz();
1588   gPad->SetRightMargin(0.15);
1589   fhDXYvsDZ[0]->Draw("COLZ");
1590
1591   canvas3->cd(4);
1592   fhDXYNormalized[0]->SetStats(kFALSE);
1593   gPad->SetLogy();
1594   fhDXYNormalized[0]->Draw();
1595
1596   canvas3->cd(5);
1597   fhDZNormalized[0]->SetStats(kFALSE);
1598   gPad->SetLogy();
1599   fhDZNormalized[0]->Draw();
1600
1601   canvas3->cd(6);
1602   fhDXYvsDZNormalized[0]->SetStats(kFALSE);
1603   gPad->SetLogz();
1604   gPad->SetRightMargin(0.15);
1605   fhDXYvsDZNormalized[0]->Draw("COLZ");
1606
1607   canvas3->SaveAs(Form("%s_%s.gif", GetName(), canvas3->GetName()));
1608
1609   TCanvas* canvas4 = new TCanvas(Form("%s_4", GetName()), "Track Quality Results4", 800, 500);
1610   canvas4->Divide(2, 1);
1611
1612   canvas4->cd(1);
1613   fhCutStatistics->SetStats(kFALSE);
1614   fhCutStatistics->LabelsOption("v");
1615   gPad->SetBottomMargin(0.3);
1616   fhCutStatistics->Draw();
1617
1618   canvas4->cd(2);
1619   fhCutCorrelation->SetStats(kFALSE);
1620   fhCutCorrelation->LabelsOption("v");
1621   gPad->SetBottomMargin(0.3);
1622   gPad->SetLeftMargin(0.3);
1623   fhCutCorrelation->Draw("COLZ");
1624
1625   canvas4->SaveAs(Form("%s_%s.gif", GetName(), canvas4->GetName()));
1626
1627   /*canvas->cd(1);
1628   fhDXYvsDZNormalized[0]->SetStats(kFALSE);
1629   fhDXYvsDZNormalized[0]->DrawCopy("COLZ");
1630
1631   canvas->cd(2);
1632   fhNClustersTPC[0]->SetStats(kFALSE);
1633   fhNClustersTPC[0]->DrawCopy();
1634
1635   canvas->cd(3);
1636   fhChi2PerClusterITS[0]->SetStats(kFALSE);
1637   fhChi2PerClusterITS[0]->DrawCopy();
1638   fhChi2PerClusterITS[1]->SetLineColor(2);
1639   fhChi2PerClusterITS[1]->DrawCopy("SAME");
1640
1641   canvas->cd(4);
1642   fhChi2PerClusterTPC[0]->SetStats(kFALSE);
1643   fhChi2PerClusterTPC[0]->DrawCopy();
1644   fhChi2PerClusterTPC[1]->SetLineColor(2);
1645   fhChi2PerClusterTPC[1]->DrawCopy("SAME");*/
1646 }
1647 //--------------------------------------------------------------------------
1648 void AliESDtrackCuts::SetPtDepDCACuts(Double_t pt) {
1649   //
1650   // set the pt-dependent DCA cuts
1651   //
1652
1653   if(f1CutMaxDCAToVertexXYPtDep) {
1654      fCutMaxDCAToVertexXY=f1CutMaxDCAToVertexXYPtDep->Eval(pt);
1655   }
1656
1657   if(f1CutMaxDCAToVertexZPtDep) {
1658     fCutMaxDCAToVertexZ=f1CutMaxDCAToVertexZPtDep->Eval(pt);
1659   }
1660
1661   if(f1CutMinDCAToVertexXYPtDep) {
1662     fCutMinDCAToVertexXY=f1CutMinDCAToVertexXYPtDep->Eval(pt);
1663   }
1664
1665   if(f1CutMinDCAToVertexZPtDep) {
1666     fCutMinDCAToVertexZ=f1CutMinDCAToVertexZPtDep->Eval(pt);
1667   }
1668
1669
1670   return;
1671 }
1672
1673
1674
1675 //--------------------------------------------------------------------------
1676 Bool_t AliESDtrackCuts::CheckPtDepDCA(TString dist,Bool_t print) const {
1677   //
1678   // Check the correctness of the string syntax
1679   //
1680   Bool_t retval=kTRUE;
1681
1682   if(!dist.Contains("pt")) {
1683     if(print) printf("AliESDtrackCuts::CheckPtDepDCA(): string must contain \"pt\"\n");
1684     retval= kFALSE;
1685   }
1686   return retval;
1687 }
1688
1689  void AliESDtrackCuts::SetMaxDCAToVertexXYPtDep(const char *dist){
1690
1691    if(f1CutMaxDCAToVertexXYPtDep){
1692      delete f1CutMaxDCAToVertexXYPtDep;
1693      // resetiing both
1694      f1CutMaxDCAToVertexXYPtDep = 0;
1695      fCutMaxDCAToVertexXYPtDep = "";
1696    }
1697    if(!CheckPtDepDCA(dist,kTRUE)){
1698      return;
1699    }
1700    fCutMaxDCAToVertexXYPtDep = dist;
1701    TString tmp(dist);
1702    tmp.ReplaceAll("pt","x");
1703    f1CutMaxDCAToVertexXYPtDep = new TFormula("f1CutMaxDCAToVertexXYPtDep",tmp.Data());
1704
1705 }
1706
1707  void AliESDtrackCuts::SetMaxDCAToVertexZPtDep(const char *dist){
1708
1709
1710    if(f1CutMaxDCAToVertexZPtDep){
1711      delete f1CutMaxDCAToVertexZPtDep;
1712      // resetiing both
1713      f1CutMaxDCAToVertexZPtDep = 0;
1714      fCutMaxDCAToVertexZPtDep = "";
1715    }
1716    if(!CheckPtDepDCA(dist,kTRUE))return;
1717
1718    fCutMaxDCAToVertexZPtDep = dist;
1719    TString tmp(dist);
1720    tmp.ReplaceAll("pt","x");
1721    f1CutMaxDCAToVertexZPtDep = new TFormula("f1CutMaxDCAToVertexZPtDep",tmp.Data());
1722
1723
1724 }
1725
1726
1727  void AliESDtrackCuts::SetMinDCAToVertexXYPtDep(const char *dist){
1728
1729
1730    if(f1CutMinDCAToVertexXYPtDep){
1731      delete f1CutMinDCAToVertexXYPtDep;
1732      // resetiing both
1733      f1CutMinDCAToVertexXYPtDep = 0;
1734      fCutMinDCAToVertexXYPtDep = "";
1735    }
1736    if(!CheckPtDepDCA(dist,kTRUE))return;
1737
1738    fCutMinDCAToVertexXYPtDep = dist;
1739    TString tmp(dist);
1740    tmp.ReplaceAll("pt","x");
1741    f1CutMinDCAToVertexXYPtDep = new TFormula("f1CutMinDCAToVertexXYPtDep",tmp.Data());
1742
1743 }
1744
1745
1746  void AliESDtrackCuts::SetMinDCAToVertexZPtDep(const char *dist){
1747
1748
1749
1750    if(f1CutMinDCAToVertexZPtDep){
1751      delete f1CutMinDCAToVertexZPtDep;
1752      // resetiing both
1753      f1CutMinDCAToVertexZPtDep = 0;
1754      fCutMinDCAToVertexZPtDep = "";
1755    }
1756    if(!CheckPtDepDCA(dist,kTRUE))return;
1757    fCutMinDCAToVertexZPtDep = dist;
1758    TString tmp(dist);
1759    tmp.ReplaceAll("pt","x");
1760    f1CutMinDCAToVertexZPtDep = new TFormula("f1CutMinDCAToVertexZPtDep",tmp.Data());
1761 }