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 <AliMultiplicity.h>
  24 #include <AliLog.h>
  25
  26 #include <TTree.h>
  27 #include <TCanvas.h>
  28 #include <TDirectory.h>
  29 #include <TH2F.h>
  30 #include <TF1.h>
  31 #include <TBits.h>
  32
  33 //____________________________________________________________________
  34 ClassImp(AliESDtrackCuts)
  35
  36 // Cut names
  37 const Char_t* AliESDtrackCuts::fgkCutNames[kNCuts] = {
  38  "require TPC refit",
  39  "require TPC standalone",
  40  "require ITS refit",
  41  "n clusters TPC",
  42  "n clusters ITS",
  43  "#Chi^{2}/cluster TPC",
  44  "#Chi^{2}/cluster ITS",
  45  "cov 11",
  46  "cov 22",
  47  "cov 33",
  48  "cov 44",
  49  "cov 55",
  50  "trk-to-vtx",
  51  "trk-to-vtx failed",
  52  "kink daughters",
  53  "p",
  54  "p_{T}",
  55  "p_{x}",
  56  "p_{y}",
  57  "p_{z}",
  58  "eta",
  59  "y",
  60  "trk-to-vtx max dca 2D absolute",
  61  "trk-to-vtx max dca xy absolute",
  62  "trk-to-vtx max dca z absolute",
  63  "trk-to-vtx min dca 2D absolute",
  64  "trk-to-vtx min dca xy absolute",
  65  "trk-to-vtx min dca z absolute",
  66  "SPD cluster requirement",
  67  "SDD cluster requirement",
  68  "SSD cluster requirement",
  69  "require ITS stand-alone",
  70  "rel 1/pt uncertainty",
  71  "TPC n shared clusters",
  72  "TPC rel shared clusters",
  73  "require ITS Pid",
  74  "n crossed rows TPC",
  75  "n crossed rows / n findable clusters",
  76  "missing ITS points",
  77  "#Chi^{2} TPC constrained vs. global",
  78  "require TOF out",
  79  "TOF Distance cut",
  80  "min length in active volume TPC"
  81 };
  82
  83 AliESDtrackCuts* AliESDtrackCuts::fgMultEstTrackCuts[AliESDtrackCuts::kNMultEstTrackCuts] = { 0, 0, 0, 0 };
  84 Char_t AliESDtrackCuts::fgBeamTypeFlag = -1;
  85
  86 //____________________________________________________________________
  87 AliESDtrackCuts::AliESDtrackCuts(const Char_t* name, const Char_t* title) : AliAnalysisCuts(name,title),
  88   fCutMinNClusterTPC(0),
  89   fCutMinNClusterITS(0),
  90   fCutMinNCrossedRowsTPC(0),
  91   fCutMinRatioCrossedRowsOverFindableClustersTPC(0),
  92   f1CutMinNClustersTPCPtDep(0x0),
  93   fCutMaxPtDepNClustersTPC(0),
  94   fCutMinLengthActiveVolumeTPC(0),
  95   fCutMaxChi2PerClusterTPC(0),
  96   fCutMaxChi2PerClusterITS(0),
  97   fCutMaxChi2TPCConstrainedVsGlobal(0),
  98   fCutMaxChi2TPCConstrainedVsGlobalVertexType(kVertexTracks | kVertexSPD),
  99   fCutMaxMissingITSPoints(0),
 100   fCutMaxC11(0),
 101   fCutMaxC22(0),
 102   fCutMaxC33(0),
 103   fCutMaxC44(0),
 104   fCutMaxC55(0),
 105   fCutMaxRel1PtUncertainty(0),
 106   fCutAcceptKinkDaughters(0),
 107   fCutAcceptSharedTPCClusters(0),
 108   fCutMaxFractionSharedTPCClusters(0),
 109   fCutRequireTPCRefit(0),
 110   fCutRequireTPCStandAlone(0),
 111   fCutRequireITSRefit(0),
 112   fCutRequireITSPid(0),
 113   fCutRequireITSStandAlone(0),
 114   fCutRequireITSpureSA(0),
 115   fCutNsigmaToVertex(0),
 116   fCutSigmaToVertexRequired(0),
 117   fCutMaxDCAToVertexXY(0),
 118   fCutMaxDCAToVertexZ(0),
 119   fCutMinDCAToVertexXY(0),
 120   fCutMinDCAToVertexZ(0),
 121   fCutMaxDCAToVertexXYPtDep(""),
 122   fCutMaxDCAToVertexZPtDep(""),
 123   fCutMinDCAToVertexXYPtDep(""),
 124   fCutMinDCAToVertexZPtDep(""),
 125   f1CutMaxDCAToVertexXYPtDep(0x0),
 126   f1CutMaxDCAToVertexZPtDep(0x0),
 127   f1CutMinDCAToVertexXYPtDep(0x0),
 128   f1CutMinDCAToVertexZPtDep(0x0),
 129   fCutDCAToVertex2D(0),
 130   fPMin(0),
 131   fPMax(0),
 132   fPtMin(0),
 133   fPtMax(0),
 134   fPxMin(0),
 135   fPxMax(0),
 136   fPyMin(0),
 137   fPyMax(0),
 138   fPzMin(0),
 139   fPzMax(0),
 140   fEtaMin(0),
 141   fEtaMax(0),
 142   fRapMin(0),
 143   fRapMax(0),
 144   fCutRequireTOFout(kFALSE),
 145   fFlagCutTOFdistance(kFALSE),
 146   fCutTOFdistance(3.),
 147   fHistogramsOn(0),
 148   ffDTheoretical(0),
 149   fhCutStatistics(0),
 150   fhCutCorrelation(0)
 151 {
 152   //
 153   // constructor
 154   //
 155
 156   Init();
 157
 158   //##############################################################################
 159   // setting default cuts
 160   SetMinNClustersTPC();
 161   SetMinNClustersITS();
 162   SetMinNCrossedRowsTPC();
 163   SetMinRatioCrossedRowsOverFindableClustersTPC();
 164   SetMaxChi2PerClusterTPC();
 165   SetMaxChi2PerClusterITS();
 166   SetMaxChi2TPCConstrainedGlobal();
 167   SetMaxChi2TPCConstrainedGlobalVertexType();
 168   SetMaxNOfMissingITSPoints();
 169   SetMaxCovDiagonalElements();
 170   SetMaxRel1PtUncertainty();
 171   SetRequireTPCRefit();
 172   SetRequireTPCStandAlone();
 173   SetRequireITSRefit();
 174   SetRequireITSPid(kFALSE);
 175   SetRequireITSStandAlone(kFALSE);
 176   SetRequireITSPureStandAlone(kFALSE);
 177   SetAcceptKinkDaughters();
 178   SetAcceptSharedTPCClusters();
 179   SetMaxFractionSharedTPCClusters();
 180   SetMaxNsigmaToVertex();
 181   SetMaxDCAToVertexXY();
 182   SetMaxDCAToVertexZ();
 183   SetDCAToVertex2D();
 184   SetMinDCAToVertexXY();
 185   SetMinDCAToVertexZ();
 186   SetPRange();
 187   SetPtRange();
 188   SetPxRange();
 189   SetPyRange();
 190   SetPzRange();
 191   SetEtaRange();
 192   SetRapRange();
 193   SetClusterRequirementITS(kSPD);
 194   SetClusterRequirementITS(kSDD);
 195   SetClusterRequirementITS(kSSD);
 196
 197   SetHistogramsOn();
 198 }
 199
 200 //_____________________________________________________________________________
 201 AliESDtrackCuts::AliESDtrackCuts(const AliESDtrackCuts &c) : AliAnalysisCuts(c),
 202   fCutMinNClusterTPC(0),
 203   fCutMinNClusterITS(0),
 204   fCutMinNCrossedRowsTPC(0),
 205   fCutMinRatioCrossedRowsOverFindableClustersTPC(0),
 206   f1CutMinNClustersTPCPtDep(0x0),
 207   fCutMaxPtDepNClustersTPC(0),
 208   fCutMinLengthActiveVolumeTPC(0),
 209   fCutMaxChi2PerClusterTPC(0),
 210   fCutMaxChi2PerClusterITS(0),
 211   fCutMaxChi2TPCConstrainedVsGlobal(0),
 212   fCutMaxChi2TPCConstrainedVsGlobalVertexType(kVertexTracks | kVertexSPD),
 213   fCutMaxMissingITSPoints(0),
 214   fCutMaxC11(0),
 215   fCutMaxC22(0),
 216   fCutMaxC33(0),
 217   fCutMaxC44(0),
 218   fCutMaxC55(0),
 219   fCutMaxRel1PtUncertainty(0),
 220   fCutAcceptKinkDaughters(0),
 221   fCutAcceptSharedTPCClusters(0),
 222   fCutMaxFractionSharedTPCClusters(0),
 223   fCutRequireTPCRefit(0),
 224   fCutRequireTPCStandAlone(0),
 225   fCutRequireITSRefit(0),
 226   fCutRequireITSPid(0),
 227   fCutRequireITSStandAlone(0),
 228   fCutRequireITSpureSA(0),
 229   fCutNsigmaToVertex(0),
 230   fCutSigmaToVertexRequired(0),
 231   fCutMaxDCAToVertexXY(0),
 232   fCutMaxDCAToVertexZ(0),
 233   fCutMinDCAToVertexXY(0),
 234   fCutMinDCAToVertexZ(0),
 235   fCutMaxDCAToVertexXYPtDep(""),
 236   fCutMaxDCAToVertexZPtDep(""),
 237   fCutMinDCAToVertexXYPtDep(""),
 238   fCutMinDCAToVertexZPtDep(""),
 239   f1CutMaxDCAToVertexXYPtDep(0x0),
 240   f1CutMaxDCAToVertexZPtDep(0x0),
 241   f1CutMinDCAToVertexXYPtDep(0x0),
 242   f1CutMinDCAToVertexZPtDep(0x0),
 243   fCutDCAToVertex2D(0),
 244   fPMin(0),
 245   fPMax(0),
 246   fPtMin(0),
 247   fPtMax(0),
 248   fPxMin(0),
 249   fPxMax(0),
 250   fPyMin(0),
 251   fPyMax(0),
 252   fPzMin(0),
 253   fPzMax(0),
 254   fEtaMin(0),
 255   fEtaMax(0),
 256   fRapMin(0),
 257   fRapMax(0),
 258   fCutRequireTOFout(kFALSE),
 259   fFlagCutTOFdistance(kFALSE),
 260   fCutTOFdistance(3.),
 261   fHistogramsOn(0),
 262   ffDTheoretical(0),
 263   fhCutStatistics(0),
 264   fhCutCorrelation(0)
 265 {
 266   //
 267   // copy constructor
 268   //
 269
 270   ((AliESDtrackCuts &) c).Copy(*this);
 271 }
 272
 273 AliESDtrackCuts::~AliESDtrackCuts()
 274 {
 275   //
 276   // destructor
 277   //
 278
 279   for (Int_t i=0; i<2; i++) {
 280
 281     if (fhNClustersITS[i])
 282       delete fhNClustersITS[i];
 283     if (fhNClustersTPC[i])
 284       delete fhNClustersTPC[i];
 285     if (fhNSharedClustersTPC[i])
 286       delete fhNSharedClustersTPC[i];
 287     if (fhNCrossedRowsTPC[i])
 288       delete fhNCrossedRowsTPC[i];
 289     if (fhRatioCrossedRowsOverFindableClustersTPC[i])
 290       delete fhRatioCrossedRowsOverFindableClustersTPC[i];
 291     if (fhChi2PerClusterITS[i])
 292       delete fhChi2PerClusterITS[i];
 293     if (fhChi2PerClusterTPC[i])
 294       delete fhChi2PerClusterTPC[i];
 295     if (fhChi2TPCConstrainedVsGlobal[i])
 296       delete fhChi2TPCConstrainedVsGlobal[i];
 297     if(fhNClustersForITSPID[i])
 298       delete fhNClustersForITSPID[i];
 299     if(fhNMissingITSPoints[i])
 300       delete fhNMissingITSPoints[i];
 301     if (fhC11[i])
 302       delete fhC11[i];
 303     if (fhC22[i])
 304       delete fhC22[i];
 305     if (fhC33[i])
 306       delete fhC33[i];
 307     if (fhC44[i])
 308       delete fhC44[i];
 309     if (fhC55[i])
 310       delete fhC55[i];
 311
 312     if (fhRel1PtUncertainty[i])
 313       delete fhRel1PtUncertainty[i];
 314
 315     if (fhDXY[i])
 316       delete fhDXY[i];
 317     if (fhDZ[i])
 318       delete fhDZ[i];
 319     if (fhDXYDZ[i])
 320       delete fhDXYDZ[i];
 321     if (fhDXYvsDZ[i])
 322       delete fhDXYvsDZ[i];
 323
 324     if (fhDXYNormalized[i])
 325       delete fhDXYNormalized[i];
 326     if (fhDZNormalized[i])
 327       delete fhDZNormalized[i];
 328     if (fhDXYvsDZNormalized[i])
 329       delete fhDXYvsDZNormalized[i];
 330     if (fhNSigmaToVertex[i])
 331       delete fhNSigmaToVertex[i];
 332     if (fhPt[i])
 333       delete fhPt[i];
 334     if (fhEta[i])
 335       delete fhEta[i];
 336     if (fhTOFdistance[i])
 337       delete fhTOFdistance[i];
 338   }
 339
 340   if(f1CutMaxDCAToVertexXYPtDep)delete f1CutMaxDCAToVertexXYPtDep;
 341   f1CutMaxDCAToVertexXYPtDep = 0;
 342   if( f1CutMaxDCAToVertexZPtDep) delete  f1CutMaxDCAToVertexZPtDep;
 343   f1CutMaxDCAToVertexZPtDep = 0;
 344   if( f1CutMinDCAToVertexXYPtDep)delete  f1CutMinDCAToVertexXYPtDep;
 345   f1CutMinDCAToVertexXYPtDep = 0;
 346   if(f1CutMinDCAToVertexZPtDep)delete  f1CutMinDCAToVertexZPtDep;
 347   f1CutMinDCAToVertexZPtDep = 0;
 348
 349
 350   if (ffDTheoretical)
 351     delete ffDTheoretical;
 352
 353   if (fhCutStatistics)
 354     delete fhCutStatistics;
 355   if (fhCutCorrelation)
 356     delete fhCutCorrelation;
 357
 358   if(f1CutMinNClustersTPCPtDep)
 359     delete f1CutMinNClustersTPCPtDep;
 360
 361 }
 362
 363 void AliESDtrackCuts::Init()
 364 {
 365   //
 366   // sets everything to zero
 367   //
 368
 369   fCutMinNClusterTPC = 0;
 370   fCutMinNClusterITS = 0;
 371
 372   fCutMaxChi2PerClusterTPC = 0;
 373   fCutMaxChi2PerClusterITS = 0;
 374   fCutMaxChi2TPCConstrainedVsGlobal = 0;
 375   fCutMaxChi2TPCConstrainedVsGlobalVertexType = kVertexTracks | kVertexSPD;
 376   fCutMaxMissingITSPoints  = 0;
 377
 378   for (Int_t i = 0; i < 3; i++)
 379         fCutClusterRequirementITS[i] = kOff;
 380
 381   fCutMaxC11 = 0;
 382   fCutMaxC22 = 0;
 383   fCutMaxC33 = 0;
 384   fCutMaxC44 = 0;
 385   fCutMaxC55 = 0;
 386
 387   fCutMaxRel1PtUncertainty = 0;
 388
 389   fCutAcceptKinkDaughters = 0;
 390   fCutAcceptSharedTPCClusters = 0;
 391   fCutMaxFractionSharedTPCClusters = 0;
 392   fCutRequireTPCRefit = 0;
 393   fCutRequireTPCStandAlone = 0;
 394   fCutRequireITSRefit = 0;
 395   fCutRequireITSPid = 0;
 396   fCutRequireITSStandAlone = 0;
 397   fCutRequireITSpureSA = 0;
 398
 399   fCutNsigmaToVertex = 0;
 400   fCutSigmaToVertexRequired = 0;
 401   fCutMaxDCAToVertexXY = 0;
 402   fCutMaxDCAToVertexZ = 0;
 403   fCutDCAToVertex2D = 0;
 404   fCutMinDCAToVertexXY = 0;
 405   fCutMinDCAToVertexZ = 0;
 406   fCutMaxDCAToVertexXYPtDep = "";
 407   fCutMaxDCAToVertexZPtDep = "";
 408   fCutMinDCAToVertexXYPtDep = "";
 409   fCutMinDCAToVertexZPtDep = "";
 410
 411   if(f1CutMaxDCAToVertexXYPtDep)delete f1CutMaxDCAToVertexXYPtDep;
 412   f1CutMaxDCAToVertexXYPtDep = 0;
 413   if( f1CutMaxDCAToVertexXYPtDep) delete  f1CutMaxDCAToVertexXYPtDep;
 414   f1CutMaxDCAToVertexXYPtDep = 0;
 415   if( f1CutMaxDCAToVertexZPtDep) delete  f1CutMaxDCAToVertexZPtDep;
 416   f1CutMaxDCAToVertexZPtDep = 0;
 417   if( f1CutMinDCAToVertexXYPtDep)delete  f1CutMinDCAToVertexXYPtDep;
 418   f1CutMinDCAToVertexXYPtDep = 0;
 419   if(f1CutMinDCAToVertexZPtDep)delete f1CutMinDCAToVertexZPtDep;
 420   f1CutMinDCAToVertexZPtDep = 0;
 421
 422
 423   fPMin = 0;
 424   fPMax = 0;
 425   fPtMin = 0;
 426   fPtMax = 0;
 427   fPxMin = 0;
 428   fPxMax = 0;
 429   fPyMin = 0;
 430   fPyMax = 0;
 431   fPzMin = 0;
 432   fPzMax = 0;
 433   fEtaMin = 0;
 434   fEtaMax = 0;
 435   fRapMin = 0;
 436   fRapMax = 0;
 437
 438   fHistogramsOn = kFALSE;
 439
 440   for (Int_t i=0; i<2; ++i)
 441   {
 442     fhNClustersITS[i] = 0;
 443     fhNClustersTPC[i] = 0;
 444     fhNSharedClustersTPC[i] = 0;
 445     fhNCrossedRowsTPC[i] = 0;
 446     fhRatioCrossedRowsOverFindableClustersTPC[i] = 0;
 447
 448     fhChi2PerClusterITS[i] = 0;
 449     fhChi2PerClusterTPC[i] = 0;
 450     fhChi2TPCConstrainedVsGlobal[i] = 0;
 451     fhNClustersForITSPID[i] = 0;
 452     fhNMissingITSPoints[i] = 0;
 453
 454     fhC11[i] = 0;
 455     fhC22[i] = 0;
 456     fhC33[i] = 0;
 457     fhC44[i] = 0;
 458     fhC55[i] = 0;
 459
 460     fhRel1PtUncertainty[i] = 0;
 461
 462     fhDXY[i] = 0;
 463     fhDZ[i] = 0;
 464     fhDXYDZ[i] = 0;
 465     fhDXYvsDZ[i] = 0;
 466
 467     fhDXYNormalized[i] = 0;
 468     fhDZNormalized[i] = 0;
 469     fhDXYvsDZNormalized[i] = 0;
 470     fhNSigmaToVertex[i] = 0;
 471
 472     fhPt[i] = 0;
 473     fhEta[i] = 0;
 474     fhTOFdistance[i] = 0;
 475   }
 476   ffDTheoretical = 0;
 477
 478   fhCutStatistics = 0;
 479   fhCutCorrelation = 0;
 480 }
 481
 482 //_____________________________________________________________________________
 483 AliESDtrackCuts &AliESDtrackCuts::operator=(const AliESDtrackCuts &c)
 484 {
 485   //
 486   // Assignment operator
 487   //
 488
 489   if (this != &c) ((AliESDtrackCuts &) c).Copy(*this);
 490   return *this;
 491 }
 492
 493 //_____________________________________________________________________________
 494 void AliESDtrackCuts::Copy(TObject &c) const
 495 {
 496   //
 497   // Copy function
 498   //
 499
 500   AliESDtrackCuts& target = (AliESDtrackCuts &) c;
 501
 502   target.Init();
 503
 504   target.fCutMinNClusterTPC = fCutMinNClusterTPC;
 505   target.fCutMinNClusterITS = fCutMinNClusterITS;
 506   target.fCutMinNCrossedRowsTPC = fCutMinNCrossedRowsTPC;
 507   target.fCutMinRatioCrossedRowsOverFindableClustersTPC = fCutMinRatioCrossedRowsOverFindableClustersTPC;
 508   if(f1CutMinNClustersTPCPtDep){
 509     target.f1CutMinNClustersTPCPtDep = (TFormula*) f1CutMinNClustersTPCPtDep->Clone("f1CutMinNClustersTPCPtDep");
 510   }
 511   target.fCutMaxPtDepNClustersTPC =   fCutMaxPtDepNClustersTPC;
 512   target.fCutMinLengthActiveVolumeTPC = fCutMinLengthActiveVolumeTPC;
 513
 514   target.fCutMaxChi2PerClusterTPC = fCutMaxChi2PerClusterTPC;
 515   target.fCutMaxChi2PerClusterITS = fCutMaxChi2PerClusterITS;
 516   target.fCutMaxChi2TPCConstrainedVsGlobal = fCutMaxChi2TPCConstrainedVsGlobal;
 517   target.fCutMaxChi2TPCConstrainedVsGlobalVertexType = fCutMaxChi2TPCConstrainedVsGlobalVertexType;
 518   target.fCutMaxMissingITSPoints = fCutMaxMissingITSPoints;
 519
 520   for (Int_t i = 0; i < 3; i++)
 521     target.fCutClusterRequirementITS[i] = fCutClusterRequirementITS[i];
 522
 523   target.fCutMaxC11 = fCutMaxC11;
 524   target.fCutMaxC22 = fCutMaxC22;
 525   target.fCutMaxC33 = fCutMaxC33;
 526   target.fCutMaxC44 = fCutMaxC44;
 527   target.fCutMaxC55 = fCutMaxC55;
 528
 529   target.fCutMaxRel1PtUncertainty = fCutMaxRel1PtUncertainty;
 530
 531   target.fCutAcceptKinkDaughters = fCutAcceptKinkDaughters;
 532   target.fCutAcceptSharedTPCClusters = fCutAcceptSharedTPCClusters;
 533   target.fCutMaxFractionSharedTPCClusters = fCutMaxFractionSharedTPCClusters;
 534   target.fCutRequireTPCRefit = fCutRequireTPCRefit;
 535   target.fCutRequireTPCStandAlone = fCutRequireTPCStandAlone;
 536   target.fCutRequireITSRefit = fCutRequireITSRefit;
 537   target.fCutRequireITSPid = fCutRequireITSPid;
 538   target.fCutRequireITSStandAlone = fCutRequireITSStandAlone;
 539   target.fCutRequireITSpureSA = fCutRequireITSpureSA;
 540
 541   target.fCutNsigmaToVertex = fCutNsigmaToVertex;
 542   target.fCutSigmaToVertexRequired = fCutSigmaToVertexRequired;
 543   target.fCutMaxDCAToVertexXY = fCutMaxDCAToVertexXY;
 544   target.fCutMaxDCAToVertexZ = fCutMaxDCAToVertexZ;
 545   target.fCutDCAToVertex2D = fCutDCAToVertex2D;
 546   target.fCutMinDCAToVertexXY = fCutMinDCAToVertexXY;
 547   target.fCutMinDCAToVertexZ = fCutMinDCAToVertexZ;
 548
 549   target.fCutMaxDCAToVertexXYPtDep = fCutMaxDCAToVertexXYPtDep;
 550   if(fCutMaxDCAToVertexXYPtDep.Length()>0)target.SetMaxDCAToVertexXYPtDep(fCutMaxDCAToVertexXYPtDep.Data());
 551
 552   target.fCutMaxDCAToVertexZPtDep = fCutMaxDCAToVertexZPtDep;
 553   if(fCutMaxDCAToVertexZPtDep.Length()>0)target.SetMaxDCAToVertexZPtDep(fCutMaxDCAToVertexZPtDep.Data());
 554
 555   target.fCutMinDCAToVertexXYPtDep = fCutMinDCAToVertexXYPtDep;
 556   if(fCutMinDCAToVertexXYPtDep.Length()>0)target.SetMinDCAToVertexXYPtDep(fCutMinDCAToVertexXYPtDep.Data());
 557
 558   target.fCutMinDCAToVertexZPtDep = fCutMinDCAToVertexZPtDep;
 559   if(fCutMinDCAToVertexZPtDep.Length()>0)target.SetMinDCAToVertexZPtDep(fCutMinDCAToVertexZPtDep.Data());
 560
 561   target.fPMin = fPMin;
 562   target.fPMax = fPMax;
 563   target.fPtMin = fPtMin;
 564   target.fPtMax = fPtMax;
 565   target.fPxMin = fPxMin;
 566   target.fPxMax = fPxMax;
 567   target.fPyMin = fPyMin;
 568   target.fPyMax = fPyMax;
 569   target.fPzMin = fPzMin;
 570   target.fPzMax = fPzMax;
 571   target.fEtaMin = fEtaMin;
 572   target.fEtaMax = fEtaMax;
 573   target.fRapMin = fRapMin;
 574   target.fRapMax = fRapMax;
 575
 576   target.fFlagCutTOFdistance = fFlagCutTOFdistance;
 577   target.fCutTOFdistance = fCutTOFdistance;
 578   target.fCutRequireTOFout = fCutRequireTOFout;
 579
 580   target.fHistogramsOn = fHistogramsOn;
 581
 582   for (Int_t i=0; i<2; ++i)
 583   {
 584     if (fhNClustersITS[i]) target.fhNClustersITS[i] = (TH1F*) fhNClustersITS[i]->Clone();
 585     if (fhNClustersTPC[i]) target.fhNClustersTPC[i] = (TH1F*) fhNClustersTPC[i]->Clone();
 586     if (fhNSharedClustersTPC[i]) target.fhNSharedClustersTPC[i] = (TH1F*) fhNSharedClustersTPC[i]->Clone();
 587     if (fhNCrossedRowsTPC[i]) target.fhNCrossedRowsTPC[i] = (TH1F*) fhNCrossedRowsTPC[i]->Clone();
 588     if (fhRatioCrossedRowsOverFindableClustersTPC[i]) target.fhRatioCrossedRowsOverFindableClustersTPC[i] = (TH1F*) fhRatioCrossedRowsOverFindableClustersTPC[i]->Clone();
 589
 590     if (fhChi2PerClusterITS[i]) target.fhChi2PerClusterITS[i] = (TH1F*) fhChi2PerClusterITS[i]->Clone();
 591     if (fhChi2PerClusterTPC[i]) target.fhChi2PerClusterTPC[i] = (TH1F*) fhChi2PerClusterTPC[i]->Clone();
 592     if (fhChi2TPCConstrainedVsGlobal[i]) target.fhChi2TPCConstrainedVsGlobal[i] = (TH1F*) fhChi2TPCConstrainedVsGlobal[i]->Clone();
 593     if (fhNClustersForITSPID[i]) target.fhNClustersForITSPID[i] = (TH1F*) fhNClustersForITSPID[i]->Clone();
 594     if (fhNMissingITSPoints[i]) target.fhNMissingITSPoints[i] = (TH1F*) fhNMissingITSPoints[i]->Clone();
 595
 596     if (fhC11[i]) target.fhC11[i] = (TH1F*) fhC11[i]->Clone();
 597     if (fhC22[i]) target.fhC22[i] = (TH1F*) fhC22[i]->Clone();
 598     if (fhC33[i]) target.fhC33[i] = (TH1F*) fhC33[i]->Clone();
 599     if (fhC44[i]) target.fhC44[i] = (TH1F*) fhC44[i]->Clone();
 600     if (fhC55[i]) target.fhC55[i] = (TH1F*) fhC55[i]->Clone();
 601
 602     if (fhRel1PtUncertainty[i]) target.fhRel1PtUncertainty[i] = (TH1F*) fhRel1PtUncertainty[i]->Clone();
 603
 604     if (fhDXY[i]) target.fhDXY[i] = (TH1F*) fhDXY[i]->Clone();
 605     if (fhDZ[i]) target.fhDZ[i] = (TH1F*) fhDZ[i]->Clone();
 606     if (fhDXYDZ[i]) target.fhDXYDZ[i] = (TH1F*) fhDXYDZ[i]->Clone();
 607     if (fhDXYvsDZ[i]) target.fhDXYvsDZ[i] = (TH2F*) fhDXYvsDZ[i]->Clone();
 608
 609     if (fhDXYNormalized[i]) target.fhDXYNormalized[i] = (TH1F*) fhDXYNormalized[i]->Clone();
 610     if (fhDZNormalized[i]) target.fhDZNormalized[i] = (TH1F*) fhDZNormalized[i]->Clone();
 611     if (fhDXYvsDZNormalized[i]) target.fhDXYvsDZNormalized[i] = (TH2F*) fhDXYvsDZNormalized[i]->Clone();
 612     if (fhNSigmaToVertex[i]) target.fhNSigmaToVertex[i] = (TH1F*) fhNSigmaToVertex[i]->Clone();
 613
 614     if (fhPt[i]) target.fhPt[i] = (TH1F*) fhPt[i]->Clone();
 615     if (fhEta[i]) target.fhEta[i] = (TH1F*) fhEta[i]->Clone();
 616     if (fhTOFdistance[i]) target.fhTOFdistance[i] = (TH2F*) fhTOFdistance[i]->Clone();
 617   }
 618   if (ffDTheoretical) target.ffDTheoretical = (TF1*) ffDTheoretical->Clone();
 619
 620   if (fhCutStatistics) target.fhCutStatistics = (TH1F*) fhCutStatistics->Clone();
 621   if (fhCutCorrelation) target.fhCutCorrelation = (TH2F*) fhCutCorrelation->Clone();
 622
 623   TNamed::Copy(c);
 624 }
 625
 626 //_____________________________________________________________________________
 627 Long64_t AliESDtrackCuts::Merge(TCollection* list) {
 628   // Merge a list of AliESDtrackCuts objects with this (needed for PROOF)
 629   // Returns the number of merged objects (including this)
 630   if (!list)
 631     return 0;
 632   if (list->IsEmpty())
 633     return 1;
 634   if (!fHistogramsOn)
 635     return 0;
 636   TIterator* iter = list->MakeIterator();
 637   TObject* obj;
 638
 639   // collection of measured and generated histograms
 640   Int_t count = 0;
 641   while ((obj = iter->Next())) {
 642
 643     AliESDtrackCuts* entry = dynamic_cast<AliESDtrackCuts*>(obj);
 644     if (entry == 0)
 645       continue;
 646
 647     if (!entry->fHistogramsOn)
 648       continue;
 649
 650     for (Int_t i=0; i<2; i++) {
 651
 652       fhNClustersITS[i]      ->Add(entry->fhNClustersITS[i]     );
 653       fhNClustersTPC[i]      ->Add(entry->fhNClustersTPC[i]     );
 654       if (fhNSharedClustersTPC[i])
 655         fhNSharedClustersTPC[i]      ->Add(entry->fhNSharedClustersTPC[i]     );
 656       if (fhNCrossedRowsTPC[i])
 657         fhNCrossedRowsTPC[i]   ->Add(entry->fhNCrossedRowsTPC[i]     );
 658       if (fhRatioCrossedRowsOverFindableClustersTPC[i])
 659         fhRatioCrossedRowsOverFindableClustersTPC[i]      ->Add(entry->fhRatioCrossedRowsOverFindableClustersTPC[i]     );
 660
 661       fhChi2PerClusterITS[i] ->Add(entry->fhChi2PerClusterITS[i]);
 662       fhChi2PerClusterTPC[i] ->Add(entry->fhChi2PerClusterTPC[i]);
 663       if (fhChi2TPCConstrainedVsGlobal[i])
 664         fhChi2TPCConstrainedVsGlobal[i]->Add(entry->fhChi2TPCConstrainedVsGlobal[i]);
 665       if (fhNClustersForITSPID[i])
 666         fhNClustersForITSPID[i]->Add(entry->fhNClustersForITSPID[i]);
 667       if (fhNMissingITSPoints[i])
 668         fhNMissingITSPoints[i] ->Add(entry->fhNMissingITSPoints[i]);
 669
 670       fhC11[i]               ->Add(entry->fhC11[i]              );
 671       fhC22[i]               ->Add(entry->fhC22[i]              );
 672       fhC33[i]               ->Add(entry->fhC33[i]              );
 673       fhC44[i]               ->Add(entry->fhC44[i]              );
 674       fhC55[i]               ->Add(entry->fhC55[i]              );
 675
 676       fhRel1PtUncertainty[i] ->Add(entry->fhRel1PtUncertainty[i]);
 677
 678       fhDXY[i]               ->Add(entry->fhDXY[i]              );
 679       fhDZ[i]                ->Add(entry->fhDZ[i]               );
 680       fhDXYDZ[i]             ->Add(entry->fhDXYDZ[i]          );
 681       fhDXYvsDZ[i]           ->Add(entry->fhDXYvsDZ[i]          );
 682
 683       fhDXYNormalized[i]     ->Add(entry->fhDXYNormalized[i]    );
 684       fhDZNormalized[i]      ->Add(entry->fhDZNormalized[i]     );
 685       fhDXYvsDZNormalized[i] ->Add(entry->fhDXYvsDZNormalized[i]);
 686       fhNSigmaToVertex[i]    ->Add(entry->fhNSigmaToVertex[i]);
 687
 688       fhPt[i]                ->Add(entry->fhPt[i]);
 689       fhEta[i]               ->Add(entry->fhEta[i]);
 690       fhTOFdistance[i]       ->Add(entry->fhTOFdistance[i]);
 691     }
 692
 693     fhCutStatistics  ->Add(entry->fhCutStatistics);
 694     fhCutCorrelation ->Add(entry->fhCutCorrelation);
 695
 696     count++;
 697   }
 698   return count+1;
 699 }
 700
 701 void AliESDtrackCuts::SetMinNClustersTPCPtDep(TFormula *f1, Float_t ptmax)
 702 {
 703   //
 704   // Sets the pT dependent NClustersTPC cut
 705   //
 706
 707   if(f1){
 708     delete f1CutMinNClustersTPCPtDep;
 709     f1CutMinNClustersTPCPtDep = (TFormula*)f1->Clone("f1CutMinNClustersTPCPtDep");
 710   }
 711   fCutMaxPtDepNClustersTPC=ptmax;
 712 }
 713
 714 //____________________________________________________________________
 715 AliESDtrackCuts* AliESDtrackCuts::GetStandardTPCOnlyTrackCuts()
 716 {
 717   // creates an AliESDtrackCuts object and fills it with standard (pre data-taking) values for TPC-only cuts
 718
 719   AliInfoClass("Creating track cuts for TPC-only.");
 720
 721   AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
 722
 723   esdTrackCuts->SetMinNClustersTPC(50);
 724   esdTrackCuts->SetMaxChi2PerClusterTPC(4);
 725   esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
 726
 727   esdTrackCuts->SetMaxDCAToVertexZ(3.2);
 728   esdTrackCuts->SetMaxDCAToVertexXY(2.4);
 729   esdTrackCuts->SetDCAToVertex2D(kTRUE);
 730
 731   return esdTrackCuts;
 732 }
 733
 734 //____________________________________________________________________
 735 AliESDtrackCuts* AliESDtrackCuts::GetStandardITSTPCTrackCuts2009(Bool_t selPrimaries)
 736 {
 737   // creates an AliESDtrackCuts object and fills it with standard values for ITS-TPC cuts for pp 2009 data
 738
 739   AliInfoClass("Creating track cuts for ITS+TPC (2009 definition).");
 740
 741   AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
 742
 743   // TPC
 744   esdTrackCuts->SetRequireTPCStandAlone(kTRUE); // to get chi2 and ncls of kTPCin
 745   esdTrackCuts->SetMinNClustersTPC(70);
 746   esdTrackCuts->SetMaxChi2PerClusterTPC(4);
 747   esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
 748   esdTrackCuts->SetRequireTPCRefit(kTRUE);
 749   // ITS
 750   esdTrackCuts->SetRequireITSRefit(kTRUE);
 751   esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
 752                                          AliESDtrackCuts::kAny);
 753   if(selPrimaries) {
 754     // 7*(0.0050+0.0060/pt^0.9)
 755     esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0350+0.0420/pt^0.9");
 756     esdTrackCuts->SetMaxChi2TPCConstrainedGlobal(36);
 757   }
 758   esdTrackCuts->SetMaxDCAToVertexZ(1.e6);
 759   esdTrackCuts->SetDCAToVertex2D(kFALSE);
 760   esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
 761   //esdTrackCuts->SetEtaRange(-0.8,+0.8);
 762
 763   esdTrackCuts->SetMaxChi2PerClusterITS(36);
 764
 765   return esdTrackCuts;
 766 }
 767
 768 //____________________________________________________________________
 769 AliESDtrackCuts* AliESDtrackCuts::GetStandardITSTPCTrackCuts2011(Bool_t selPrimaries, Int_t clusterCut)
 770 {
 771   // creates an AliESDtrackCuts object and fills it with standard values for ITS-TPC cuts for pp 2011 data
 772   // if clusterCut = 1, the cut on the number of clusters is replaced by
 773   // a cut on the number of crossed rows and on the ration crossed
 774   // rows/findable clusters
 775
 776   AliInfoClass("Creating track cuts for ITS+TPC (2011 definition).");
 777
 778   AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
 779
 780   // TPC
 781   if(clusterCut == 0)  esdTrackCuts->SetMinNClustersTPC(50);
 782   else if (clusterCut == 1) {
 783     esdTrackCuts->SetMinNCrossedRowsTPC(70);
 784     esdTrackCuts->SetMinRatioCrossedRowsOverFindableClustersTPC(0.8);
 785   }
 786   else {
 787     AliWarningClass(Form("Wrong value of the clusterCut parameter (%d), using cut on Nclusters",clusterCut));
 788     esdTrackCuts->SetMinNClustersTPC(50);
 789   }
 790   esdTrackCuts->SetMaxChi2PerClusterTPC(4);
 791   esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
 792   esdTrackCuts->SetRequireTPCRefit(kTRUE);
 793   // ITS
 794   esdTrackCuts->SetRequireITSRefit(kTRUE);
 795   esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
 796                                          AliESDtrackCuts::kAny);
 797   if(selPrimaries) {
 798     // 7*(0.0015+0.0050/pt^1.1)
 799     esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0105+0.0350/pt^1.1");
 800     esdTrackCuts->SetMaxChi2TPCConstrainedGlobal(36);
 801   }
 802   esdTrackCuts->SetMaxDCAToVertexZ(2);
 803   esdTrackCuts->SetDCAToVertex2D(kFALSE);
 804   esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
 805
 806   esdTrackCuts->SetMaxChi2PerClusterITS(36);
 807
 808   return esdTrackCuts;
 809 }
 810
 811 //____________________________________________________________________
 812 AliESDtrackCuts* AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(Bool_t selPrimaries,Int_t clusterCut)
 813 {
 814   // creates an AliESDtrackCuts object and fills it with standard values for ITS-TPC cuts for pp 2010 data
 815   // if clusterCut = 1, the cut on the number of clusters is replaced by
 816   // a cut on the number of crossed rows and on the ration crossed
 817   // rows/findable clusters
 818
 819   AliInfoClass("Creating track cuts for ITS+TPC (2010 definition).");
 820
 821   AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
 822
 823   // TPC
 824   if(clusterCut == 0)  esdTrackCuts->SetMinNClustersTPC(70);
 825   else if (clusterCut == 1) {
 826     esdTrackCuts->SetMinNCrossedRowsTPC(70);
 827     esdTrackCuts->SetMinRatioCrossedRowsOverFindableClustersTPC(0.8);
 828   }
 829   else {
 830     AliWarningClass(Form("Wrong value of the clusterCut parameter (%d), using cut on Nclusters",clusterCut));
 831     esdTrackCuts->SetMinNClustersTPC(70);
 832   }
 833   esdTrackCuts->SetMaxChi2PerClusterTPC(4);
 834   esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
 835   esdTrackCuts->SetRequireTPCRefit(kTRUE);
 836   // ITS
 837   esdTrackCuts->SetRequireITSRefit(kTRUE);
 838   esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
 839                                          AliESDtrackCuts::kAny);
 840   if(selPrimaries) {
 841     // 7*(0.0026+0.0050/pt^1.01)
 842     esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0182+0.0350/pt^1.01");
 843     esdTrackCuts->SetMaxChi2TPCConstrainedGlobal(36);
 844   }
 845   esdTrackCuts->SetMaxDCAToVertexZ(2);
 846   esdTrackCuts->SetDCAToVertex2D(kFALSE);
 847   esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
 848
 849   esdTrackCuts->SetMaxChi2PerClusterITS(36);
 850
 851   return esdTrackCuts;
 852 }
 853
 854 //____________________________________________________________________
 855 AliESDtrackCuts* AliESDtrackCuts::GetStandardITSPureSATrackCuts2009(Bool_t selPrimaries, Bool_t useForPid)
 856 {
 857   // creates an AliESDtrackCuts object and fills it with standard values for ITS pure SA tracks
 858
 859   AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
 860   esdTrackCuts->SetRequireITSPureStandAlone(kTRUE);
 861   esdTrackCuts->SetRequireITSRefit(kTRUE);
 862   esdTrackCuts->SetMinNClustersITS(4);
 863   esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
 864                                          AliESDtrackCuts::kAny);
 865   esdTrackCuts->SetMaxChi2PerClusterITS(1.);
 866
 867   if(selPrimaries) {
 868     // 7*(0.0085+0.0026/pt^1.55)
 869     esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0595+0.0182/pt^1.55");
 870   }
 871   if(useForPid){
 872     esdTrackCuts->SetRequireITSPid(kTRUE);
 873   }
 874   return esdTrackCuts;
 875 }
 876
 877 //____________________________________________________________________
 878 AliESDtrackCuts* AliESDtrackCuts::GetStandardITSPureSATrackCuts2010(Bool_t selPrimaries, Bool_t useForPid)
 879 {
 880   // creates an AliESDtrackCuts object and fills it with standard values for ITS pure SA tracks - pp 2010
 881
 882   AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
 883   esdTrackCuts->SetRequireITSPureStandAlone(kTRUE);
 884   esdTrackCuts->SetRequireITSRefit(kTRUE);
 885   esdTrackCuts->SetMinNClustersITS(4);
 886   esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
 887                                          AliESDtrackCuts::kAny);
 888   esdTrackCuts->SetMaxChi2PerClusterITS(2.5);
 889
 890   if(selPrimaries) {
 891     // 7*(0.0033+0.0045/pt^1.3)
 892     esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0231+0.0315/pt^1.3");
 893   }
 894   if(useForPid){
 895     esdTrackCuts->SetRequireITSPid(kTRUE);
 896   }
 897   return esdTrackCuts;
 898 }
 899
 900 //____________________________________________________________________
 901 AliESDtrackCuts* AliESDtrackCuts::GetStandardITSSATrackCuts2009(Bool_t selPrimaries, Bool_t useForPid)
 902 {
 903   // creates an AliESDtrackCuts object and fills it with standard values for ITS pure SA tracks
 904
 905   AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
 906   esdTrackCuts->SetRequireITSStandAlone(kTRUE);
 907   esdTrackCuts->SetRequireITSPureStandAlone(kFALSE);
 908   esdTrackCuts->SetRequireITSRefit(kTRUE);
 909   esdTrackCuts->SetMinNClustersITS(4);
 910   esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
 911                                          AliESDtrackCuts::kAny);
 912   esdTrackCuts->SetMaxChi2PerClusterITS(1.);
 913
 914   if(selPrimaries) {
 915     // 7*(0.0085+0.0026/pt^1.55)
 916     esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0595+0.0182/pt^1.55");
 917   }
 918   if(useForPid){
 919     esdTrackCuts->SetRequireITSPid(kTRUE);
 920   }
 921   return esdTrackCuts;
 922 }
 923
 924 //____________________________________________________________________
 925 AliESDtrackCuts* AliESDtrackCuts::GetStandardITSSATrackCuts2010(Bool_t selPrimaries, Bool_t useForPid)
 926 {
 927   // creates an AliESDtrackCuts object and fills it with standard values for ITS pure SA tracks --pp 2010
 928
 929   AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
 930   esdTrackCuts->SetRequireITSStandAlone(kTRUE);
 931   esdTrackCuts->SetRequireITSPureStandAlone(kFALSE);
 932   esdTrackCuts->SetRequireITSRefit(kTRUE);
 933   esdTrackCuts->SetMinNClustersITS(4);
 934   esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
 935                                          AliESDtrackCuts::kAny);
 936   esdTrackCuts->SetMaxChi2PerClusterITS(2.5);
 937
 938   if(selPrimaries) {
 939     // 7*(0.0033+0.0045/pt^1.3)
 940     esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0231+0.0315/pt^1.3");
 941   }
 942   if(useForPid){
 943     esdTrackCuts->SetRequireITSPid(kTRUE);
 944   }
 945   return esdTrackCuts;
 946 }
 947
 948 //____________________________________________________________________
 949 AliESDtrackCuts* AliESDtrackCuts::GetStandardITSSATrackCutsPbPb2010(Bool_t selPrimaries, Bool_t useForPid)
 950 {
 951   // creates an AliESDtrackCuts object and fills it with standard values for ITS pure SA tracks -- PbPb 2010
 952
 953   AliESDtrackCuts* esdTrackCuts = GetStandardITSSATrackCuts2010(selPrimaries, useForPid);
 954   esdTrackCuts->SetMaxNOfMissingITSPoints(1);
 955
 956   return esdTrackCuts;
 957 }
 958 //____________________________________________________________________
 959
 960 AliESDtrackCuts* AliESDtrackCuts::GetStandardV0DaughterCuts()
 961 {
 962   // creates a AliESDtrackCuts object and fills it with standard cuts for V0 daughters
 963   AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
 964   esdTrackCuts->SetRequireTPCRefit(kTRUE);
 965   esdTrackCuts->SetMinNClustersTPC(70);
 966   esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
 967   return esdTrackCuts;
 968 }
 969
 970 //____________________________________________________________________
 971 Int_t AliESDtrackCuts::GetReferenceMultiplicity(const AliESDEvent* esd, Bool_t tpcOnly)
 972 {
 973   // Gets reference multiplicity following the standard cuts and a defined fiducial volume
 974   // tpcOnly = kTRUE -> consider TPC-only tracks
 975   //         = kFALSE -> consider global tracks
 976   //
 977   // DEPRECATED Use GetReferenceMultiplicity with the enum as second argument instead
 978
 979   if (!tpcOnly)
 980   {
 981     AliErrorClass("Not implemented for global tracks!");
 982     return -1;
 983   }
 984
 985   static AliESDtrackCuts* esdTrackCuts = 0;
 986   if (!esdTrackCuts)
 987   {
 988     esdTrackCuts = GetStandardTPCOnlyTrackCuts();
 989     esdTrackCuts->SetEtaRange(-0.8, 0.8);
 990     esdTrackCuts->SetPtRange(0.15);
 991   }
 992
 993   Int_t nTracks = esdTrackCuts->CountAcceptedTracks(esd);
 994
 995   return nTracks;
 996 }
 997
 998 //____________________________________________________________________
 999 Float_t AliESDtrackCuts::GetSigmaToVertex(const AliESDtrack* const esdTrack)
1000 {
1001   // Calculates the number of sigma to the vertex.
1002
1003   Float_t b[2];
1004   Float_t bRes[2];
1005   Float_t bCov[3];
1006   esdTrack->GetImpactParameters(b,bCov);
1007
1008   if (bCov[0]<=0 || bCov[2]<=0) {
1009     AliDebugClass(1, "Estimated b resolution lower or equal zero!");
1010     bCov[0]=0; bCov[2]=0;
1011   }
1012   bRes[0] = TMath::Sqrt(bCov[0]);
1013   bRes[1] = TMath::Sqrt(bCov[2]);
1014
1015   // -----------------------------------
1016   // How to get to a n-sigma cut?
1017   //
1018   // The accumulated statistics from 0 to d is
1019   //
1020   // ->  Erf(d/Sqrt(2)) for a 1-dim gauss (d = n_sigma)
1021   // ->  1 - Exp(-d**2) for a 2-dim gauss (d*d = dx*dx + dy*dy != n_sigma)
1022   //
1023   // It means that for a 2-dim gauss: n_sigma(d) = Sqrt(2)*ErfInv(1 - Exp((-d**2)/2)
1024   // Can this be expressed in a different way?
1025
1026   if (bRes[0] == 0 || bRes[1] ==0)
1027     return -1;
1028
1029   Float_t d = TMath::Sqrt(TMath::Power(b[0]/bRes[0],2) + TMath::Power(b[1]/bRes[1],2));
1030
1031   // work around precision problem
1032   // if d is too big, TMath::Exp(...) gets 0, and TMath::ErfInverse(1) that should be infinite, gets 0 :(
1033   // 1e-15 corresponds to nsigma ~ 7.7
1034   if (TMath::Exp(-d * d / 2) < 1e-15)
1035     return 1000;
1036
1037   Float_t nSigma = TMath::ErfInverse(1 - TMath::Exp(-d * d / 2)) * TMath::Sqrt(2);
1038   return nSigma;
1039 }
1040
1041 void AliESDtrackCuts::EnableNeededBranches(TTree* tree)
1042 {
1043   // enables the branches needed by AcceptTrack, for a list see comment of AcceptTrack
1044
1045   tree->SetBranchStatus("fTracks.fFlags", 1);
1046   tree->SetBranchStatus("fTracks.fITSncls", 1);
1047   tree->SetBranchStatus("fTracks.fTPCncls", 1);
1048   tree->SetBranchStatus("fTracks.fITSchi2", 1);
1049   tree->SetBranchStatus("fTracks.fTPCchi2", 1);
1050   tree->SetBranchStatus("fTracks.fC*", 1);
1051   tree->SetBranchStatus("fTracks.fD", 1);
1052   tree->SetBranchStatus("fTracks.fZ", 1);
1053   tree->SetBranchStatus("fTracks.fCdd", 1);
1054   tree->SetBranchStatus("fTracks.fCdz", 1);
1055   tree->SetBranchStatus("fTracks.fCzz", 1);
1056   tree->SetBranchStatus("fTracks.fP*", 1);
1057   tree->SetBranchStatus("fTracks.fR*", 1);
1058   tree->SetBranchStatus("fTracks.fKinkIndexes*", 1);
1059 }
1060
1061 //____________________________________________________________________
1062 Bool_t AliESDtrackCuts::AcceptTrack(const AliESDtrack* esdTrack)
1063 {
1064   //
1065   // figure out if the tracks survives all the track cuts defined
1066   //
1067   // the different quality parameter and kinematic values are first
1068   // retrieved from the track. then it is found out what cuts the
1069   // track did not survive and finally the cuts are imposed.
1070
1071   // this function needs the following branches:
1072   // fTracks.fFlags
1073   // fTracks.fITSncls
1074   // fTracks.fTPCncls
1075   // fTracks.fITSchi2
1076   // fTracks.fTPCchi2
1077   // fTracks.fC   //GetExternalCovariance
1078   // fTracks.fD   //GetImpactParameters
1079   // fTracks.fZ   //GetImpactParameters
1080   // fTracks.fCdd //GetImpactParameters
1081   // fTracks.fCdz //GetImpactParameters
1082   // fTracks.fCzz //GetImpactParameters
1083   // fTracks.fP   //GetPxPyPz
1084   // fTracks.fR   //GetMass
1085   // fTracks.fP   //GetMass
1086   // fTracks.fKinkIndexes
1087   //
1088   // esdEvent is only required for the MaxChi2TPCConstrainedVsGlobal
1089
1090   UInt_t status = esdTrack->GetStatus();
1091
1092   // getting quality parameters from the ESD track
1093   Int_t nClustersITS = esdTrack->GetITSclusters(0);
1094   Int_t nClustersTPC = -1;
1095   if(fCutRequireTPCStandAlone) {
1096     nClustersTPC = esdTrack->GetTPCNclsIter1();
1097   }
1098   else {
1099     nClustersTPC = esdTrack->GetTPCclusters(0);
1100   }
1101
1102   //Pt dependent NClusters Cut
1103   if(f1CutMinNClustersTPCPtDep) {
1104     if(esdTrack->Pt()<fCutMaxPtDepNClustersTPC)
1105       fCutMinNClusterTPC = (Int_t)(f1CutMinNClustersTPCPtDep->Eval(esdTrack->Pt()));
1106     else
1107       fCutMinNClusterTPC = (Int_t)(f1CutMinNClustersTPCPtDep->Eval(fCutMaxPtDepNClustersTPC));
1108   }
1109
1110   Float_t nCrossedRowsTPC = esdTrack->GetTPCCrossedRows();
1111   Float_t  ratioCrossedRowsOverFindableClustersTPC = 1.0;
1112   if (esdTrack->GetTPCNclsF()>0) {
1113     ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC / esdTrack->GetTPCNclsF();
1114   }
1115
1116   Int_t nClustersTPCShared = esdTrack->GetTPCnclsS();
1117   Float_t fracClustersTPCShared = -1.;
1118
1119   Float_t chi2PerClusterITS = -1;
1120   Float_t chi2PerClusterTPC = -1;
1121   if (nClustersITS!=0)
1122     chi2PerClusterITS = esdTrack->GetITSchi2()/Float_t(nClustersITS);
1123   if (nClustersTPC!=0) {
1124     if(fCutRequireTPCStandAlone) {
1125       chi2PerClusterTPC = esdTrack->GetTPCchi2Iter1()/Float_t(nClustersTPC);
1126     } else {
1127       chi2PerClusterTPC = esdTrack->GetTPCchi2()/Float_t(nClustersTPC);
1128     }
1129     fracClustersTPCShared = Float_t(nClustersTPCShared)/Float_t(nClustersTPC);
1130   }
1131
1132   Double_t extCov[15];
1133   esdTrack->GetExternalCovariance(extCov);
1134
1135   Float_t b[2];
1136   Float_t bCov[3];
1137   esdTrack->GetImpactParameters(b,bCov);
1138   if (bCov[0]<=0 || bCov[2]<=0) {
1139     AliDebug(1, "Estimated b resolution lower or equal zero!");
1140     bCov[0]=0; bCov[2]=0;
1141   }
1142
1143
1144   // set pt-dependent DCA cuts, if requested
1145   SetPtDepDCACuts(esdTrack->Pt());
1146
1147
1148   Float_t dcaToVertexXY = b[0];
1149   Float_t dcaToVertexZ = b[1];
1150
1151   Float_t dcaToVertex = -1;
1152
1153   if (fCutDCAToVertex2D)
1154   {
1155     dcaToVertex = TMath::Sqrt(dcaToVertexXY*dcaToVertexXY/fCutMaxDCAToVertexXY/fCutMaxDCAToVertexXY + dcaToVertexZ*dcaToVertexZ/fCutMaxDCAToVertexZ/fCutMaxDCAToVertexZ);
1156   }
1157   else
1158     dcaToVertex = TMath::Sqrt(dcaToVertexXY*dcaToVertexXY + dcaToVertexZ*dcaToVertexZ);
1159
1160   // getting the kinematic variables of the track
1161   // (assuming the mass is known)
1162   Double_t p[3];
1163   esdTrack->GetPxPyPz(p);
1164
1165   // Changed from float to double to prevent rounding errors leading to negative
1166   // log arguments (M.G.)
1167   Double_t momentum = TMath::Sqrt(p[0]*p[0] + p[1]*p[1] + p[2]*p[2]);
1168   Double_t pt       = TMath::Sqrt(p[0]*p[0] + p[1]*p[1]);
1169   Double_t mass     = esdTrack->GetMass();
1170   Double_t energy   = TMath::Sqrt(mass*mass + momentum*momentum);
1171
1172   //y-eta related calculations
1173   Float_t eta = -100.;
1174   Float_t y   = -100.;
1175   if((momentum != TMath::Abs(p[2]))&&(momentum != 0))
1176     eta = 0.5*TMath::Log((momentum + p[2])/(momentum - p[2]));
1177   if((energy != TMath::Abs(p[2]))&&(momentum != 0))
1178     y = 0.5*TMath::Log((energy + p[2])/(energy - p[2]));
1179
1180   if (extCov[14] < 0)
1181   {
1182     AliWarning(Form("GetSigma1Pt2() returns negative value for external covariance matrix element fC[14]: %f. Corrupted track information, track will not be accepted!", extCov[14]));
1183     return kFALSE;
1184   }
1185   Float_t relUncertainty1Pt = TMath::Sqrt(extCov[14])*pt;
1186
1187   //########################################################################
1188   // cut the track?
1189
1190   Bool_t cuts[kNCuts];
1191   for (Int_t i=0; i<kNCuts; i++) cuts[i]=kFALSE;
1192
1193   // track quality cuts
1194   if (fCutRequireTPCRefit && (status&AliESDtrack::kTPCrefit)==0)
1195     cuts[0]=kTRUE;
1196   if (fCutRequireTPCStandAlone && (status&AliESDtrack::kTPCin)==0)
1197     cuts[1]=kTRUE;
1198   if (fCutRequireITSRefit && (status&AliESDtrack::kITSrefit)==0)
1199     cuts[2]=kTRUE;
1200   if (nClustersTPC<fCutMinNClusterTPC)
1201     cuts[3]=kTRUE;
1202   if (nClustersITS<fCutMinNClusterITS)
1203     cuts[4]=kTRUE;
1204   if (chi2PerClusterTPC>fCutMaxChi2PerClusterTPC)
1205     cuts[5]=kTRUE;
1206   if (chi2PerClusterITS>fCutMaxChi2PerClusterITS)
1207     cuts[6]=kTRUE;
1208   if (extCov[0]  > fCutMaxC11)
1209     cuts[7]=kTRUE;
1210   if (extCov[2]  > fCutMaxC22)
1211     cuts[8]=kTRUE;
1212   if (extCov[5]  > fCutMaxC33)
1213     cuts[9]=kTRUE;
1214   if (extCov[9]  > fCutMaxC44)
1215     cuts[10]=kTRUE;
1216   if (extCov[14]  > fCutMaxC55)
1217     cuts[11]=kTRUE;
1218
1219   // cut 12 and 13 see below
1220
1221   if (!fCutAcceptKinkDaughters && esdTrack->GetKinkIndex(0)>0)
1222     cuts[14]=kTRUE;
1223   // track kinematics cut
1224   if((momentum < fPMin) || (momentum > fPMax))
1225     cuts[15]=kTRUE;
1226   if((pt < fPtMin) || (pt > fPtMax))
1227     cuts[16] = kTRUE;
1228   if((p[0] < fPxMin) || (p[0] > fPxMax))
1229     cuts[17] = kTRUE;
1230   if((p[1] < fPyMin) || (p[1] > fPyMax))
1231     cuts[18] = kTRUE;
1232   if((p[2] < fPzMin) || (p[2] > fPzMax))
1233     cuts[19] = kTRUE;
1234   if((eta < fEtaMin) || (eta > fEtaMax))
1235     cuts[20] = kTRUE;
1236   if((y < fRapMin) || (y > fRapMax))
1237     cuts[21] = kTRUE;
1238   if (fCutDCAToVertex2D && dcaToVertex > 1)
1239     cuts[22] = kTRUE;
1240   if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexXY) > fCutMaxDCAToVertexXY)
1241     cuts[23] = kTRUE;
1242   if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexZ) > fCutMaxDCAToVertexZ)
1243     cuts[24] = kTRUE;
1244   if (fCutDCAToVertex2D && fCutMinDCAToVertexXY > 0 && fCutMinDCAToVertexZ > 0 && dcaToVertexXY*dcaToVertexXY/fCutMinDCAToVertexXY/fCutMinDCAToVertexXY + dcaToVertexZ*dcaToVertexZ/fCutMinDCAToVertexZ/fCutMinDCAToVertexZ < 1)
1245     cuts[25] = kTRUE;
1246   if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexXY) < fCutMinDCAToVertexXY)
1247     cuts[26] = kTRUE;
1248   if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexZ) < fCutMinDCAToVertexZ)
1249     cuts[27] = kTRUE;
1250
1251   for (Int_t i = 0; i < 3; i++) {
1252     if(!(esdTrack->GetStatus()&AliESDtrack::kITSupg)) { // current ITS
1253       cuts[28+i] = !CheckITSClusterRequirement(fCutClusterRequirementITS[i], esdTrack->HasPointOnITSLayer(i*2), esdTrack->HasPointOnITSLayer(i*2+1));
1254     } else { // upgraded ITS (7 layers)
1255       // at the moment, for L012 the layers 12 are considered together
1256       if(i==0) { // L012
1257         cuts[28+i] = !CheckITSClusterRequirement(fCutClusterRequirementITS[i], esdTrack->HasPointOnITSLayer(0), (esdTrack->HasPointOnITSLayer(1))&(esdTrack->HasPointOnITSLayer(2)));
1258       } else { // L34 or L56
1259         cuts[28+i] = !CheckITSClusterRequirement(fCutClusterRequirementITS[i], esdTrack->HasPointOnITSLayer(i*2+1), esdTrack->HasPointOnITSLayer(i*2+2));
1260       }
1261     }
1262   }
1263
1264   if(fCutRequireITSStandAlone || fCutRequireITSpureSA){
1265     if ((status & AliESDtrack::kITSin) == 0 || (status & AliESDtrack::kTPCin)){
1266       // TPC tracks
1267       cuts[31] = kTRUE;
1268     }else{
1269       // ITS standalone tracks
1270       if(fCutRequireITSStandAlone && !fCutRequireITSpureSA){
1271         if(status & AliESDtrack::kITSpureSA) cuts[31] = kTRUE;
1272       }else if(fCutRequireITSpureSA){
1273         if(!(status & AliESDtrack::kITSpureSA)) cuts[31] = kTRUE;
1274       }
1275     }
1276   }
1277
1278   if (relUncertainty1Pt > fCutMaxRel1PtUncertainty)
1279      cuts[32] = kTRUE;
1280
1281   if (!fCutAcceptSharedTPCClusters && nClustersTPCShared!=0)
1282     cuts[33] = kTRUE;
1283
1284   if (fracClustersTPCShared > fCutMaxFractionSharedTPCClusters)
1285     cuts[34] = kTRUE;
1286
1287   Int_t nITSPointsForPid=0;
1288   UChar_t clumap=esdTrack->GetITSClusterMap();
1289   for(Int_t i=2; i<6; i++){
1290     if(clumap&(1<<i)) ++nITSPointsForPid;
1291   }
1292   if(fCutRequireITSPid && nITSPointsForPid<3) cuts[35] = kTRUE;
1293
1294
1295   if (nCrossedRowsTPC<fCutMinNCrossedRowsTPC)
1296     cuts[36]=kTRUE;
1297   if (ratioCrossedRowsOverFindableClustersTPC<fCutMinRatioCrossedRowsOverFindableClustersTPC)
1298     cuts[37]=kTRUE;
1299
1300   Int_t nMissITSpts=0;
1301   Int_t idet,statusLay;
1302   Float_t xloc,zloc;
1303   for(Int_t iLay=0; iLay<6; iLay++){
1304     Bool_t retc=esdTrack->GetITSModuleIndexInfo(iLay,idet,statusLay,xloc,zloc);
1305     if(retc && statusLay==5) ++nMissITSpts;
1306   }
1307   if(nMissITSpts>fCutMaxMissingITSPoints) cuts[38] = kTRUE;
1308
1309   //kTOFout
1310   if (fCutRequireTOFout && (status&AliESDtrack::kTOFout)==0)
1311     cuts[40]=kTRUE;
1312
1313   // TOF signal Dz cut
1314   Float_t dxTOF = esdTrack->GetTOFsignalDx();
1315   Float_t dzTOF = esdTrack->GetTOFsignalDz();
1316   if (fFlagCutTOFdistance && (esdTrack->GetStatus() & AliESDtrack::kTOFout) == AliESDtrack::kTOFout){ // applying the TOF distance cut only if requested, and only on tracks that reached the TOF and where associated with a TOF hit
1317           if (fgBeamTypeFlag < 0) {  // the check on the beam type was not done yet
1318                   const AliESDEvent* event = esdTrack->GetESDEvent();
1319                   if (event){
1320                           TString beamTypeESD = event->GetBeamType();
1321                           AliDebug(2,Form("Beam type from ESD event = %s",beamTypeESD.Data()));
1322                           if (beamTypeESD.CompareTo("A-A",TString::kIgnoreCase) == 0){ // we are in PbPb collisions --> fgBeamTypeFlag will be set to 1, to apply the cut on TOF signal Dz
1323                                   fgBeamTypeFlag = 1;
1324                           }
1325                           else { // we are NOT in PbPb collisions --> fgBeamTypeFlag will be set to 0, to NOT apply the cu6 on TOF signal Dz
1326                                   fgBeamTypeFlag = 0;
1327                           }
1328                   }
1329                   else{
1330                           AliFatal("Beam type not available, but it is needed to apply the TOF cut!");
1331                   }
1332           }
1333
1334           if (fgBeamTypeFlag == 1){ // we are in PbPb collisions --> apply the cut on TOF signal Dz
1335                   Float_t radiusTOF = TMath::Sqrt(dxTOF*dxTOF + dzTOF*dzTOF);
1336                   AliDebug(3,Form("TOF check (with fCutTOFdistance = %f) --> dx = %f, dz = %f, radius = %f", fCutTOFdistance, dxTOF, dzTOF, radiusTOF));
1337                   if (radiusTOF > fCutTOFdistance){
1338                           AliDebug(2, Form("************* the radius is outside the range! %f > %f, the track will be skipped", radiusTOF, fCutTOFdistance));
1339                           cuts[41] = kTRUE;
1340                   }
1341           }
1342   }
1343
1344   Bool_t cut=kFALSE;
1345   for (Int_t i=0; i<kNCuts; i++)
1346     if (cuts[i]) {cut = kTRUE;}
1347
1348   // for performance evaluate the CPU intensive cuts only when the others have passed, and when they are requested
1349   Double_t chi2TPCConstrainedVsGlobal = -2;
1350   Float_t nSigmaToVertex = -2;
1351   if (!cut)
1352   {
1353     // getting the track to vertex parameters
1354     if (fCutSigmaToVertexRequired)
1355     {
1356       nSigmaToVertex = GetSigmaToVertex(esdTrack);
1357       if (nSigmaToVertex > fCutNsigmaToVertex && fCutSigmaToVertexRequired)
1358       {
1359         cuts[12] = kTRUE;
1360         cut = kTRUE;
1361       }
1362       // if n sigma could not be calculated
1363       if (nSigmaToVertex<0 && fCutSigmaToVertexRequired)
1364       {
1365         cuts[13] = kTRUE;
1366         cut = kTRUE;
1367       }
1368     }
1369
1370     // max chi2 TPC constrained vs global track only if track passed the other cut
1371     if (fCutMaxChi2TPCConstrainedVsGlobal < 1e9)
1372     {
1373       const AliESDEvent* esdEvent = esdTrack->GetESDEvent();
1374
1375       if (!esdEvent)
1376         AliFatal("fCutMaxChi2TPCConstrainedVsGlobal set but ESD event not set in AliESDTrack. Use AliESDTrack::SetESDEvent before calling AliESDtrackCuts.");
1377
1378       // get vertex
1379       const AliESDVertex* vertex = 0;
1380       if (fCutMaxChi2TPCConstrainedVsGlobalVertexType & kVertexTracks)
1381         vertex = esdEvent->GetPrimaryVertexTracks();
1382
1383       if ((!vertex || !vertex->GetStatus()) && fCutMaxChi2TPCConstrainedVsGlobalVertexType & kVertexSPD)
1384         vertex = esdEvent->GetPrimaryVertexSPD();
1385
1386       if ((!vertex || !vertex->GetStatus()) && fCutMaxChi2TPCConstrainedVsGlobalVertexType & kVertexTPC)
1387         vertex = esdEvent->GetPrimaryVertexTPC();
1388
1389       if (vertex->GetStatus())
1390         chi2TPCConstrainedVsGlobal = esdTrack->GetChi2TPCConstrainedVsGlobal(vertex);
1391
1392       if (chi2TPCConstrainedVsGlobal < 0 || chi2TPCConstrainedVsGlobal > fCutMaxChi2TPCConstrainedVsGlobal)
1393       {
1394         cuts[39] = kTRUE;
1395         cut = kTRUE;
1396       }
1397     }
1398
1399     // max length in active volume
1400     Float_t lengthInActiveZoneTPC = -1;
1401     if (fCutMinLengthActiveVolumeTPC > 1.) { // do the calculation only if needed to save cpu-time
1402       if (esdTrack->GetESDEvent()) {
1403         if (esdTrack->GetInnerParam()) lengthInActiveZoneTPC = esdTrack->GetLengthInActiveZone(1, 1.8, 220, esdTrack->GetESDEvent()->GetMagneticField());
1404         //
1405         if (lengthInActiveZoneTPC < fCutMinLengthActiveVolumeTPC ) {
1406           cuts[42] = kTRUE;
1407           cut = kTRUE;
1408         }
1409       }
1410     }
1411
1412
1413   }
1414
1415   //########################################################################
1416   // filling histograms
1417   if (fHistogramsOn) {
1418     fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n tracks")));
1419     if (cut)
1420       fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n cut tracks")));
1421
1422     for (Int_t i=0; i<kNCuts; i++) {
1423       if (fhCutStatistics->GetXaxis()->FindBin(fgkCutNames[i]) < 1)
1424         AliFatal(Form("Inconsistency! Cut %d with name %s not found", i, fgkCutNames[i]));
1425
1426       if (cuts[i])
1427         fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin(fgkCutNames[i])));
1428
1429       for (Int_t j=i; j<kNCuts; j++) {
1430         if (cuts[i] && cuts[j]) {
1431           Float_t xC = fhCutCorrelation->GetXaxis()->GetBinCenter(fhCutCorrelation->GetXaxis()->FindBin(fgkCutNames[i]));
1432           Float_t yC = fhCutCorrelation->GetYaxis()->GetBinCenter(fhCutCorrelation->GetYaxis()->FindBin(fgkCutNames[j]));
1433           fhCutCorrelation->Fill(xC, yC);
1434         }
1435       }
1436     }
1437   }
1438
1439   // now we loop over the filling of the histograms twice: once "before" the cut, once "after"
1440   // the code is not in a function due to too many local variables that would need to be passed
1441
1442   for (Int_t id = 0; id < 2; id++)
1443   {
1444     // id = 0 --> before cut
1445     // id = 1 --> after cut
1446
1447     if (fHistogramsOn)
1448     {
1449       fhNClustersITS[id]->Fill(nClustersITS);
1450       fhNClustersTPC[id]->Fill(nClustersTPC);
1451       fhNSharedClustersTPC[id]->Fill(nClustersTPCShared);
1452       fhNCrossedRowsTPC[id]->Fill(nCrossedRowsTPC);
1453       fhRatioCrossedRowsOverFindableClustersTPC[id]->Fill(ratioCrossedRowsOverFindableClustersTPC);
1454       fhChi2PerClusterITS[id]->Fill(chi2PerClusterITS);
1455       fhChi2PerClusterTPC[id]->Fill(chi2PerClusterTPC);
1456       fhChi2TPCConstrainedVsGlobal[id]->Fill(chi2TPCConstrainedVsGlobal);
1457       fhNClustersForITSPID[id]->Fill(nITSPointsForPid);
1458       fhNMissingITSPoints[id]->Fill(nMissITSpts);
1459
1460       fhC11[id]->Fill(extCov[0]);
1461       fhC22[id]->Fill(extCov[2]);
1462       fhC33[id]->Fill(extCov[5]);
1463       fhC44[id]->Fill(extCov[9]);
1464       fhC55[id]->Fill(extCov[14]);
1465
1466       fhRel1PtUncertainty[id]->Fill(relUncertainty1Pt);
1467
1468       fhPt[id]->Fill(pt);
1469       fhEta[id]->Fill(eta);
1470       fhTOFdistance[id]->Fill(dxTOF, dzTOF);
1471
1472       Float_t bRes[2];
1473       bRes[0] = TMath::Sqrt(bCov[0]);
1474       bRes[1] = TMath::Sqrt(bCov[2]);
1475
1476       fhDZ[id]->Fill(b[1]);
1477       fhDXY[id]->Fill(b[0]);
1478       fhDXYDZ[id]->Fill(dcaToVertex);
1479       fhDXYvsDZ[id]->Fill(b[1],b[0]);
1480
1481       if (bRes[0]!=0 && bRes[1]!=0) {
1482         fhDZNormalized[id]->Fill(b[1]/bRes[1]);
1483         fhDXYNormalized[id]->Fill(b[0]/bRes[0]);
1484         fhDXYvsDZNormalized[id]->Fill(b[1]/bRes[1], b[0]/bRes[0]);
1485         fhNSigmaToVertex[id]->Fill(nSigmaToVertex);
1486       }
1487     }
1488
1489     // cut the track
1490     if (cut)
1491       return kFALSE;
1492   }
1493
1494   return kTRUE;
1495 }
1496
1497 //____________________________________________________________________
1498 Bool_t AliESDtrackCuts::CheckITSClusterRequirement(ITSClusterRequirement req, Bool_t clusterL1, Bool_t clusterL2)
1499 {
1500   // checks if the cluster requirement is fullfilled (in this case: return kTRUE)
1501
1502   switch (req)
1503   {
1504         case kOff:        return kTRUE;
1505         case kNone:       return !clusterL1 && !clusterL2;
1506         case kAny:        return clusterL1 || clusterL2;
1507         case kFirst:      return clusterL1;
1508         case kOnlyFirst:  return clusterL1 && !clusterL2;
1509         case kSecond:     return clusterL2;
1510         case kOnlySecond: return clusterL2 && !clusterL1;
1511         case kBoth:       return clusterL1 && clusterL2;
1512   }
1513
1514   return kFALSE;
1515 }
1516
1517 //____________________________________________________________________
1518 AliESDtrack* AliESDtrackCuts::GetTPCOnlyTrack(const AliESDEvent* esd, Int_t iTrack)
1519 {
1520   // Utility function to create a TPC only track from the given esd track
1521   //
1522   // IMPORTANT: The track has to be deleted by the user
1523   //
1524   // NB. most of the functionality to get a TPC only track from an ESD track is in AliESDtrack, where it should be
1525   // there are only missing propagations here that are needed for old data
1526   // this function will therefore become obsolete
1527   //
1528   // adapted from code provided by CKB
1529
1530   if (!esd->GetPrimaryVertexTPC())
1531     return 0; // No TPC vertex no TPC tracks
1532
1533   if(!esd->GetPrimaryVertexTPC()->GetStatus())
1534     return 0; // TPC Vertex is created by default in AliESDEvent, do not use in this case
1535
1536   AliESDtrack* track = esd->GetTrack(iTrack);
1537   if (!track)
1538     return 0;
1539
1540   AliESDtrack *tpcTrack = new AliESDtrack();
1541
1542   // only true if we have a tpc track
1543   if (!track->FillTPCOnlyTrack(*tpcTrack))
1544   {
1545     delete tpcTrack;
1546     return 0;
1547   }
1548
1549   return tpcTrack;
1550 }
1551
1552 //____________________________________________________________________
1553 TObjArray* AliESDtrackCuts::GetAcceptedTracks(const AliESDEvent* esd, Bool_t bTPC)
1554 {
1555   //
1556   // returns an array of all tracks that pass the cuts
1557   // or an array of TPC only tracks (propagated to the TPC vertex during reco)
1558   // tracks that pass the cut
1559   //
1560   // NOTE: List has to be deleted by the user
1561
1562   TObjArray* acceptedTracks = new TObjArray();
1563
1564   // loop over esd tracks
1565   for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
1566     if(bTPC){
1567       if(!esd->GetPrimaryVertexTPC())return acceptedTracks; // No TPC vertex no TPC tracks
1568       if(!esd->GetPrimaryVertexTPC()->GetStatus())return acceptedTracks; // No proper TPC vertex, only the default
1569
1570       AliESDtrack *tpcTrack = GetTPCOnlyTrack(esd, iTrack);
1571       if (!tpcTrack)
1572         continue;
1573
1574       if (AcceptTrack(tpcTrack)) {
1575         acceptedTracks->Add(tpcTrack);
1576       }
1577       else
1578         delete tpcTrack;
1579     }
1580     else
1581     {
1582       AliESDtrack* track = esd->GetTrack(iTrack);
1583       if(AcceptTrack(track))
1584         acceptedTracks->Add(track);
1585     }
1586   }
1587   if(bTPC)acceptedTracks->SetOwner(kTRUE);
1588   return acceptedTracks;
1589 }
1590
1591 //____________________________________________________________________
1592 Int_t AliESDtrackCuts::CountAcceptedTracks(const AliESDEvent* const esd)
1593 {
1594   //
1595   // returns an the number of tracks that pass the cuts
1596   //
1597
1598   Int_t count = 0;
1599
1600   // loop over esd tracks
1601   for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
1602     AliESDtrack* track = esd->GetTrack(iTrack);
1603     if (AcceptTrack(track))
1604       count++;
1605   }
1606
1607   return count;
1608 }
1609
1610 //____________________________________________________________________
1611  void AliESDtrackCuts::DefineHistograms(Int_t color) {
1612    //
1613    // diagnostics histograms are defined
1614    //
1615
1616    fHistogramsOn=kTRUE;
1617
1618    Bool_t oldStatus = TH1::AddDirectoryStatus();
1619    TH1::AddDirectory(kFALSE);
1620
1621    //###################################################################################
1622    // defining histograms
1623
1624    fhCutStatistics = new TH1F("cut_statistics","cut statistics",kNCuts+4,-0.5,kNCuts+3.5);
1625
1626    fhCutStatistics->GetXaxis()->SetBinLabel(1,"n tracks");
1627    fhCutStatistics->GetXaxis()->SetBinLabel(2,"n cut tracks");
1628
1629    fhCutCorrelation = new TH2F("cut_correlation","cut correlation",kNCuts,-0.5,kNCuts-0.5,kNCuts,-0.5,kNCuts-0.5);;
1630
1631    for (Int_t i=0; i<kNCuts; i++) {
1632      fhCutStatistics->GetXaxis()->SetBinLabel(i+4,fgkCutNames[i]);
1633      fhCutCorrelation->GetXaxis()->SetBinLabel(i+1,fgkCutNames[i]);
1634      fhCutCorrelation->GetYaxis()->SetBinLabel(i+1,fgkCutNames[i]);
1635    }
1636
1637   fhCutStatistics  ->SetLineColor(color);
1638   fhCutCorrelation ->SetLineColor(color);
1639   fhCutStatistics  ->SetLineWidth(2);
1640   fhCutCorrelation ->SetLineWidth(2);
1641
1642   for (Int_t i=0; i<2; i++) {
1643     fhNClustersITS[i]        = new TH1F("nClustersITS"    ,"",8,-0.5,7.5);
1644     fhNClustersTPC[i]        = new TH1F("nClustersTPC"     ,"",165,-0.5,164.5);
1645     fhNSharedClustersTPC[i]  = new TH1F("nSharedClustersTPC"     ,"",165,-0.5,164.5);
1646     fhNCrossedRowsTPC[i]     = new TH1F("nCrossedRowsTPC"  ,"",165,-0.5,164.5);
1647     fhRatioCrossedRowsOverFindableClustersTPC[i]     = new TH1F("ratioCrossedRowsOverFindableClustersTPC"  ,"",60,0,1.5);
1648     fhChi2PerClusterITS[i]   = new TH1F("chi2PerClusterITS","",500,0,10);
1649     fhChi2PerClusterTPC[i]   = new TH1F("chi2PerClusterTPC","",500,0,10);
1650     fhChi2TPCConstrainedVsGlobal[i]   = new TH1F("chi2TPCConstrainedVsGlobal","",600,-2,50);
1651     fhNClustersForITSPID[i]  = new TH1F("nPointsForITSpid","",5,-0.5,4.5);
1652     fhNMissingITSPoints[i]   = new TH1F("nMissingITSClusters","",7,-0.5,6.5);
1653
1654     fhC11[i]                 = new TH1F("covMatrixDiagonal11","",2000,0,20);
1655     fhC22[i]                 = new TH1F("covMatrixDiagonal22","",2000,0,20);
1656     fhC33[i]                 = new TH1F("covMatrixDiagonal33","",1000,0,0.1);
1657     fhC44[i]                 = new TH1F("covMatrixDiagonal44","",1000,0,0.1);
1658     fhC55[i]                 = new TH1F("covMatrixDiagonal55","",1000,0,5);
1659
1660     fhRel1PtUncertainty[i]   = new TH1F("rel1PtUncertainty","",1000,0,5);
1661
1662     fhDXY[i]                 = new TH1F("dXY"    ,"",500,-10,10);
1663     fhDZ[i]                  = new TH1F("dZ"     ,"",500,-10,10);
1664     fhDXYDZ[i]               = new TH1F("dXYDZ"  ,"",500,0,10);
1665     fhDXYvsDZ[i]             = new TH2F("dXYvsDZ","",200,-10,10,200,-10,10);
1666
1667     fhDXYNormalized[i]       = new TH1F("dXYNormalized"    ,"",500,-10,10);
1668     fhDZNormalized[i]        = new TH1F("dZNormalized"     ,"",500,-10,10);
1669     fhDXYvsDZNormalized[i]   = new TH2F("dXYvsDZNormalized","",200,-10,10,200,-10,10);
1670
1671     fhNSigmaToVertex[i]      = new TH1F("nSigmaToVertex","",500,0,10);
1672
1673     fhPt[i]                  = new TH1F("pt"     ,"p_{T} distribution;p_{T} (GeV/c)", 800, 0.0, 10.0);
1674     fhEta[i]                 = new TH1F("eta"     ,"#eta distribution;#eta",40,-2.0,2.0);
1675     fhTOFdistance[i]         = new TH2F("TOFdistance"     ,"TOF distance;dx (cm};dz (cm)", 150, -15, 15, 150, -15, 15);
1676
1677     fhNClustersITS[i]->SetTitle("n ITS clusters");
1678     fhNClustersTPC[i]->SetTitle("n TPC clusters");
1679     fhNSharedClustersTPC[i]->SetTitle("n TPC shared clusters");
1680     fhChi2PerClusterITS[i]->SetTitle("#Chi^{2} per ITS cluster");
1681     fhChi2PerClusterTPC[i]->SetTitle("#Chi^{2} per TPC cluster");
1682     fhChi2TPCConstrainedVsGlobal[i]->SetTitle("#Chi^{2} TPC constrained track vs global track");
1683     fhNClustersForITSPID[i]->SetTitle("n ITS points for PID");
1684     fhNMissingITSPoints[i]->SetTitle("n ITS layers with missing cluster");
1685
1686     fhC11[i]->SetTitle("cov 11 : #sigma_{y}^{2} [cm^{2}]");
1687     fhC22[i]->SetTitle("cov 22 : #sigma_{z}^{2} [cm^{2}]");
1688     fhC33[i]->SetTitle("cov 33 : #sigma_{sin(#phi)}^{2}");
1689     fhC44[i]->SetTitle("cov 44 : #sigma_{tan(#theta_{dip})}^{2}");
1690     fhC55[i]->SetTitle("cov 55 : #sigma_{1/p_{T}}^{2} [(c/GeV)^2]");
1691
1692     fhRel1PtUncertainty[i]->SetTitle("rel. uncertainty of 1/p_{T}");
1693
1694     fhDXY[i]->SetXTitle("transverse impact parameter (cm)");
1695     fhDZ[i]->SetXTitle("longitudinal impact parameter (cm)");
1696     fhDXYDZ[i]->SetTitle("absolute impact parameter;sqrt(dXY**2 + dZ**2)  (cm)");
1697     fhDXYvsDZ[i]->SetXTitle("longitudinal impact parameter (cm)");
1698     fhDXYvsDZ[i]->SetYTitle("transverse impact parameter (cm)");
1699
1700     fhDXYNormalized[i]->SetTitle("normalized trans impact par (n#sigma)");
1701     fhDZNormalized[i]->SetTitle("normalized long impact par (n#sigma)");
1702     fhDXYvsDZNormalized[i]->SetTitle("normalized long impact par (n#sigma)");
1703     fhDXYvsDZNormalized[i]->SetYTitle("normalized trans impact par (n#sigma)");
1704     fhNSigmaToVertex[i]->SetTitle("n #sigma to vertex");
1705
1706     fhNClustersITS[i]->SetLineColor(color);   fhNClustersITS[i]->SetLineWidth(2);
1707     fhNClustersTPC[i]->SetLineColor(color);   fhNClustersTPC[i]->SetLineWidth(2);
1708     fhNSharedClustersTPC[i]->SetLineColor(color);  fhNSharedClustersTPC[i]->SetLineWidth(2);
1709     fhChi2PerClusterITS[i]->SetLineColor(color);   fhChi2PerClusterITS[i]->SetLineWidth(2);
1710     fhChi2PerClusterTPC[i]->SetLineColor(color);   fhChi2PerClusterTPC[i]->SetLineWidth(2);
1711     fhChi2TPCConstrainedVsGlobal[i]->SetLineColor(color);   fhChi2TPCConstrainedVsGlobal[i]->SetLineWidth(2);
1712     fhNClustersForITSPID[i]->SetLineColor(color);  fhNClustersForITSPID[i]->SetLineWidth(2);
1713     fhNMissingITSPoints[i]->SetLineColor(color);   fhNMissingITSPoints[i]->SetLineWidth(2);
1714
1715     fhC11[i]->SetLineColor(color);   fhC11[i]->SetLineWidth(2);
1716     fhC22[i]->SetLineColor(color);   fhC22[i]->SetLineWidth(2);
1717     fhC33[i]->SetLineColor(color);   fhC33[i]->SetLineWidth(2);
1718     fhC44[i]->SetLineColor(color);   fhC44[i]->SetLineWidth(2);
1719     fhC55[i]->SetLineColor(color);   fhC55[i]->SetLineWidth(2);
1720
1721     fhRel1PtUncertainty[i]->SetLineColor(color); fhRel1PtUncertainty[i]->SetLineWidth(2);
1722
1723     fhDXY[i]->SetLineColor(color);   fhDXY[i]->SetLineWidth(2);
1724     fhDZ[i]->SetLineColor(color);    fhDZ[i]->SetLineWidth(2);
1725     fhDXYDZ[i]->SetLineColor(color); fhDXYDZ[i]->SetLineWidth(2);
1726
1727     fhDXYNormalized[i]->SetLineColor(color);   fhDXYNormalized[i]->SetLineWidth(2);
1728     fhDZNormalized[i]->SetLineColor(color);    fhDZNormalized[i]->SetLineWidth(2);
1729     fhNSigmaToVertex[i]->SetLineColor(color);  fhNSigmaToVertex[i]->SetLineWidth(2);
1730   }
1731
1732   // The number of sigmas to the vertex is per definition gaussian
1733   ffDTheoretical = new TF1("nSigmaToVertexTheoretical","([0]/2.506628274)*exp(-(x**2)/2)",0,50);
1734   ffDTheoretical->SetParameter(0,1);
1735
1736   TH1::AddDirectory(oldStatus);
1737 }
1738
1739 //____________________________________________________________________
1740 Bool_t AliESDtrackCuts::LoadHistograms(const Char_t* dir)
1741 {
1742   //
1743   // loads the histograms from a file
1744   // if dir is empty a directory with the name of this object is taken (like in SaveHistogram)
1745   //
1746
1747   if (!dir)
1748     dir = GetName();
1749
1750   if (!gDirectory->cd(dir))
1751     return kFALSE;
1752
1753   ffDTheoretical = dynamic_cast<TF1*> (gDirectory->Get("nSigmaToVertexTheory"));
1754
1755   fhCutStatistics = dynamic_cast<TH1F*> (gDirectory->Get("cut_statistics"));
1756   fhCutCorrelation = dynamic_cast<TH2F*> (gDirectory->Get("cut_correlation"));
1757
1758   for (Int_t i=0; i<2; i++) {
1759     if (i==0)
1760     {
1761       gDirectory->cd("before_cuts");
1762     }
1763     else
1764       gDirectory->cd("after_cuts");
1765
1766     fhNClustersITS[i]      = dynamic_cast<TH1F*> (gDirectory->Get("nClustersITS"     ));
1767     fhNClustersTPC[i]      = dynamic_cast<TH1F*> (gDirectory->Get("nClustersTPC"     ));
1768     fhNSharedClustersTPC[i] = dynamic_cast<TH1F*> (gDirectory->Get("nSharedClustersTPC"     ));
1769     fhNCrossedRowsTPC[i]   = dynamic_cast<TH1F*> (gDirectory->Get("nCrossedRowsTPC"  ));
1770     fhRatioCrossedRowsOverFindableClustersTPC[i]   = dynamic_cast<TH1F*> (gDirectory->Get("ratioCrossedRowsOverFindableClustersTPC"  ));
1771     fhChi2PerClusterITS[i] = dynamic_cast<TH1F*> (gDirectory->Get("chi2PerClusterITS"));
1772     fhChi2PerClusterTPC[i] = dynamic_cast<TH1F*> (gDirectory->Get("chi2PerClusterTPC"));
1773     fhChi2TPCConstrainedVsGlobal[i] = dynamic_cast<TH1F*> (gDirectory->Get("fhChi2TPCConstrainedVsGlobal"));
1774     fhNClustersForITSPID[i] = dynamic_cast<TH1F*> (gDirectory->Get("nPointsForITSpid"));
1775     fhNMissingITSPoints[i] = dynamic_cast<TH1F*> (gDirectory->Get("nMissingITSClusters"));
1776
1777     fhC11[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal11"));
1778     fhC22[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal22"));
1779     fhC33[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal33"));
1780     fhC44[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal44"));
1781     fhC55[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal55"));
1782
1783     fhRel1PtUncertainty[i] = dynamic_cast<TH1F*> (gDirectory->Get("rel1PtUncertainty"));
1784
1785     fhDXY[i] =     dynamic_cast<TH1F*> (gDirectory->Get("dXY"    ));
1786     fhDZ[i] =      dynamic_cast<TH1F*> (gDirectory->Get("dZ"     ));
1787     fhDXYDZ[i] =   dynamic_cast<TH1F*> (gDirectory->Get("dXYDZ"));
1788     fhDXYvsDZ[i] = dynamic_cast<TH2F*> (gDirectory->Get("dXYvsDZ"));
1789
1790     fhDXYNormalized[i] =     dynamic_cast<TH1F*> (gDirectory->Get("dXYNormalized"    ));
1791     fhDZNormalized[i] =      dynamic_cast<TH1F*> (gDirectory->Get("dZNormalized"     ));
1792     fhDXYvsDZNormalized[i] = dynamic_cast<TH2F*> (gDirectory->Get("dXYvsDZNormalized"));
1793     fhNSigmaToVertex[i] = dynamic_cast<TH1F*> (gDirectory->Get("nSigmaToVertex"));
1794
1795     fhPt[i] = dynamic_cast<TH1F*> (gDirectory->Get("pt"));
1796     fhEta[i] = dynamic_cast<TH1F*> (gDirectory->Get("eta"));
1797     fhTOFdistance[i] = dynamic_cast<TH2F*> (gDirectory->Get("TOFdistance"));
1798
1799     gDirectory->cd("../");
1800   }
1801
1802   gDirectory->cd("..");
1803
1804   return kTRUE;
1805 }
1806
1807 //____________________________________________________________________
1808 void AliESDtrackCuts::SaveHistograms(const Char_t* dir) {
1809   //
1810   // saves the histograms in a directory (dir)
1811   //
1812
1813   if (!fHistogramsOn) {
1814     AliDebug(0, "Histograms not on - cannot save histograms!!!");
1815     return;
1816   }
1817
1818   if (!dir)
1819     dir = GetName();
1820
1821   gDirectory->mkdir(dir);
1822   gDirectory->cd(dir);
1823
1824   gDirectory->mkdir("before_cuts");
1825   gDirectory->mkdir("after_cuts");
1826
1827   // a factor of 2 is needed since n sigma is positive
1828   ffDTheoretical->SetParameter(0,2*fhNSigmaToVertex[0]->Integral("width"));
1829   ffDTheoretical->Write("nSigmaToVertexTheory");
1830
1831   fhCutStatistics->Write();
1832   fhCutCorrelation->Write();
1833
1834   for (Int_t i=0; i<2; i++) {
1835     if (i==0)
1836       gDirectory->cd("before_cuts");
1837     else
1838       gDirectory->cd("after_cuts");
1839
1840     fhNClustersITS[i]        ->Write();
1841     fhNClustersTPC[i]        ->Write();
1842     fhNSharedClustersTPC[i]  ->Write();
1843     fhNCrossedRowsTPC[i]     ->Write();
1844     fhRatioCrossedRowsOverFindableClustersTPC[i]     ->Write();
1845     fhChi2PerClusterITS[i]   ->Write();
1846     fhChi2PerClusterTPC[i]   ->Write();
1847     fhChi2TPCConstrainedVsGlobal[i]   ->Write();
1848     fhNClustersForITSPID[i]  ->Write();
1849     fhNMissingITSPoints[i]   ->Write();
1850
1851     fhC11[i]                 ->Write();
1852     fhC22[i]                 ->Write();
1853     fhC33[i]                 ->Write();
1854     fhC44[i]                 ->Write();
1855     fhC55[i]                 ->Write();
1856
1857     fhRel1PtUncertainty[i]   ->Write();
1858
1859     fhDXY[i]                 ->Write();
1860     fhDZ[i]                  ->Write();
1861     fhDXYDZ[i]               ->Write();
1862     fhDXYvsDZ[i]             ->Write();
1863
1864     fhDXYNormalized[i]       ->Write();
1865     fhDZNormalized[i]        ->Write();
1866     fhDXYvsDZNormalized[i]   ->Write();
1867     fhNSigmaToVertex[i]      ->Write();
1868
1869     fhPt[i]                  ->Write();
1870     fhEta[i]                 ->Write();
1871     fhTOFdistance[i]         ->Write();
1872
1873     gDirectory->cd("../");
1874   }
1875
1876   gDirectory->cd("../");
1877 }
1878
1879 //____________________________________________________________________
1880 void AliESDtrackCuts::DrawHistograms()
1881 {
1882   // draws some histograms
1883
1884   TCanvas* canvas1 = new TCanvas(Form("%s_1", GetName()), "Track Quality Results1", 800, 800);
1885   canvas1->Divide(2, 2);
1886
1887   canvas1->cd(1);
1888   fhNClustersTPC[0]->SetStats(kFALSE);
1889   fhNClustersTPC[0]->Draw();
1890
1891   canvas1->cd(2);
1892   fhChi2PerClusterTPC[0]->SetStats(kFALSE);
1893   fhChi2PerClusterTPC[0]->Draw();
1894
1895   canvas1->cd(3);
1896   fhNSigmaToVertex[0]->SetStats(kFALSE);
1897   fhNSigmaToVertex[0]->GetXaxis()->SetRangeUser(0, 10);
1898   fhNSigmaToVertex[0]->Draw();
1899
1900   canvas1->SaveAs(Form("%s_%s.gif", GetName(), canvas1->GetName()));
1901
1902   TCanvas* canvas2 = new TCanvas(Form("%s_2", GetName()), "Track Quality Results2", 1200, 800);
1903   canvas2->Divide(3, 2);
1904
1905   canvas2->cd(1);
1906   fhC11[0]->SetStats(kFALSE);
1907   gPad->SetLogy();
1908   fhC11[0]->Draw();
1909
1910   canvas2->cd(2);
1911   fhC22[0]->SetStats(kFALSE);
1912   gPad->SetLogy();
1913   fhC22[0]->Draw();
1914
1915   canvas2->cd(3);
1916   fhC33[0]->SetStats(kFALSE);
1917   gPad->SetLogy();
1918   fhC33[0]->Draw();
1919
1920   canvas2->cd(4);
1921   fhC44[0]->SetStats(kFALSE);
1922   gPad->SetLogy();
1923   fhC44[0]->Draw();
1924
1925   canvas2->cd(5);
1926   fhC55[0]->SetStats(kFALSE);
1927   gPad->SetLogy();
1928   fhC55[0]->Draw();
1929
1930   canvas2->cd(6);
1931   fhRel1PtUncertainty[0]->SetStats(kFALSE);
1932   gPad->SetLogy();
1933   fhRel1PtUncertainty[0]->Draw();
1934
1935   canvas2->SaveAs(Form("%s_%s.gif", GetName(), canvas2->GetName()));
1936
1937   TCanvas* canvas3 = new TCanvas(Form("%s_3", GetName()), "Track Quality Results3", 1200, 800);
1938   canvas3->Divide(3, 2);
1939
1940   canvas3->cd(1);
1941   fhDXY[0]->SetStats(kFALSE);
1942   gPad->SetLogy();
1943   fhDXY[0]->Draw();
1944
1945   canvas3->cd(2);
1946   fhDZ[0]->SetStats(kFALSE);
1947   gPad->SetLogy();
1948   fhDZ[0]->Draw();
1949
1950   canvas3->cd(3);
1951   fhDXYvsDZ[0]->SetStats(kFALSE);
1952   gPad->SetLogz();
1953   gPad->SetRightMargin(0.15);
1954   fhDXYvsDZ[0]->Draw("COLZ");
1955
1956   canvas3->cd(4);
1957   fhDXYNormalized[0]->SetStats(kFALSE);
1958   gPad->SetLogy();
1959   fhDXYNormalized[0]->Draw();
1960
1961   canvas3->cd(5);
1962   fhDZNormalized[0]->SetStats(kFALSE);
1963   gPad->SetLogy();
1964   fhDZNormalized[0]->Draw();
1965
1966   canvas3->cd(6);
1967   fhDXYvsDZNormalized[0]->SetStats(kFALSE);
1968   gPad->SetLogz();
1969   gPad->SetRightMargin(0.15);
1970   fhDXYvsDZNormalized[0]->Draw("COLZ");
1971
1972   canvas3->SaveAs(Form("%s_%s.gif", GetName(), canvas3->GetName()));
1973
1974   TCanvas* canvas4 = new TCanvas(Form("%s_4", GetName()), "Track Quality Results4", 800, 500);
1975   canvas4->Divide(2, 1);
1976
1977   canvas4->cd(1);
1978   fhCutStatistics->SetStats(kFALSE);
1979   fhCutStatistics->LabelsOption("v");
1980   gPad->SetBottomMargin(0.3);
1981   fhCutStatistics->Draw();
1982
1983   canvas4->cd(2);
1984   fhCutCorrelation->SetStats(kFALSE);
1985   fhCutCorrelation->LabelsOption("v");
1986   gPad->SetBottomMargin(0.3);
1987   gPad->SetLeftMargin(0.3);
1988   fhCutCorrelation->Draw("COLZ");
1989
1990   canvas4->SaveAs(Form("%s_%s.gif", GetName(), canvas4->GetName()));
1991
1992   /*canvas->cd(1);
1993   fhDXYvsDZNormalized[0]->SetStats(kFALSE);
1994   fhDXYvsDZNormalized[0]->DrawCopy("COLZ");
1995
1996   canvas->cd(2);
1997   fhNClustersTPC[0]->SetStats(kFALSE);
1998   fhNClustersTPC[0]->DrawCopy();
1999
2000   canvas->cd(3);
2001   fhChi2PerClusterITS[0]->SetStats(kFALSE);
2002   fhChi2PerClusterITS[0]->DrawCopy();
2003   fhChi2PerClusterITS[1]->SetLineColor(2);
2004   fhChi2PerClusterITS[1]->DrawCopy("SAME");
2005
2006   canvas->cd(4);
2007   fhChi2PerClusterTPC[0]->SetStats(kFALSE);
2008   fhChi2PerClusterTPC[0]->DrawCopy();
2009   fhChi2PerClusterTPC[1]->SetLineColor(2);
2010   fhChi2PerClusterTPC[1]->DrawCopy("SAME");*/
2011 }
2012 //--------------------------------------------------------------------------
2013 void AliESDtrackCuts::SetPtDepDCACuts(Double_t pt) {
2014   //
2015   // set the pt-dependent DCA cuts
2016   //
2017
2018   if(f1CutMaxDCAToVertexXYPtDep) {
2019      fCutMaxDCAToVertexXY=f1CutMaxDCAToVertexXYPtDep->Eval(pt);
2020   }
2021
2022   if(f1CutMaxDCAToVertexZPtDep) {
2023     fCutMaxDCAToVertexZ=f1CutMaxDCAToVertexZPtDep->Eval(pt);
2024   }
2025
2026   if(f1CutMinDCAToVertexXYPtDep) {
2027     fCutMinDCAToVertexXY=f1CutMinDCAToVertexXYPtDep->Eval(pt);
2028   }
2029
2030   if(f1CutMinDCAToVertexZPtDep) {
2031     fCutMinDCAToVertexZ=f1CutMinDCAToVertexZPtDep->Eval(pt);
2032   }
2033
2034
2035   return;
2036 }
2037
2038
2039
2040 //--------------------------------------------------------------------------
2041 Bool_t AliESDtrackCuts::CheckPtDepDCA(TString dist,Bool_t print) const {
2042   //
2043   // Check the correctness of the string syntax
2044   //
2045   Bool_t retval=kTRUE;
2046
2047   if(!dist.Contains("pt")) {
2048     if(print) AliError("string must contain \"pt\"");
2049     retval= kFALSE;
2050   }
2051   return retval;
2052 }
2053
2054  void AliESDtrackCuts::SetMaxDCAToVertexXYPtDep(const char *dist){
2055
2056    if(f1CutMaxDCAToVertexXYPtDep){
2057      delete f1CutMaxDCAToVertexXYPtDep;
2058      // resetiing both
2059      f1CutMaxDCAToVertexXYPtDep = 0;
2060      fCutMaxDCAToVertexXYPtDep = "";
2061    }
2062    if(!CheckPtDepDCA(dist,kTRUE)){
2063      return;
2064    }
2065    fCutMaxDCAToVertexXYPtDep = dist;
2066    TString tmp(dist);
2067    tmp.ReplaceAll("pt","x");
2068    f1CutMaxDCAToVertexXYPtDep = new TFormula("f1CutMaxDCAToVertexXYPtDep",tmp.Data());
2069
2070 }
2071
2072  void AliESDtrackCuts::SetMaxDCAToVertexZPtDep(const char *dist){
2073
2074
2075    if(f1CutMaxDCAToVertexZPtDep){
2076      delete f1CutMaxDCAToVertexZPtDep;
2077      // resetiing both
2078      f1CutMaxDCAToVertexZPtDep = 0;
2079      fCutMaxDCAToVertexZPtDep = "";
2080    }
2081    if(!CheckPtDepDCA(dist,kTRUE))return;
2082
2083    fCutMaxDCAToVertexZPtDep = dist;
2084    TString tmp(dist);
2085    tmp.ReplaceAll("pt","x");
2086    f1CutMaxDCAToVertexZPtDep = new TFormula("f1CutMaxDCAToVertexZPtDep",tmp.Data());
2087
2088
2089 }
2090
2091
2092  void AliESDtrackCuts::SetMinDCAToVertexXYPtDep(const char *dist){
2093
2094
2095    if(f1CutMinDCAToVertexXYPtDep){
2096      delete f1CutMinDCAToVertexXYPtDep;
2097      // resetiing both
2098      f1CutMinDCAToVertexXYPtDep = 0;
2099      fCutMinDCAToVertexXYPtDep = "";
2100    }
2101    if(!CheckPtDepDCA(dist,kTRUE))return;
2102
2103    fCutMinDCAToVertexXYPtDep = dist;
2104    TString tmp(dist);
2105    tmp.ReplaceAll("pt","x");
2106    f1CutMinDCAToVertexXYPtDep = new TFormula("f1CutMinDCAToVertexXYPtDep",tmp.Data());
2107
2108 }
2109
2110
2111  void AliESDtrackCuts::SetMinDCAToVertexZPtDep(const char *dist){
2112
2113
2114
2115    if(f1CutMinDCAToVertexZPtDep){
2116      delete f1CutMinDCAToVertexZPtDep;
2117      // resetiing both
2118      f1CutMinDCAToVertexZPtDep = 0;
2119      fCutMinDCAToVertexZPtDep = "";
2120    }
2121    if(!CheckPtDepDCA(dist,kTRUE))return;
2122    fCutMinDCAToVertexZPtDep = dist;
2123    TString tmp(dist);
2124    tmp.ReplaceAll("pt","x");
2125    f1CutMinDCAToVertexZPtDep = new TFormula("f1CutMinDCAToVertexZPtDep",tmp.Data());
2126 }
2127
2128 AliESDtrackCuts* AliESDtrackCuts::GetMultEstTrackCuts(MultEstTrackCuts cut)
2129 {
2130   // returns the multiplicity estimator track cuts objects to allow for user configuration
2131   // upon first call the objects are created
2132   //
2133   // the cut defined here correspond to GetStandardITSTPCTrackCuts2010 (apart from the one for without SPD)
2134
2135   if (!fgMultEstTrackCuts[kMultEstTrackCutGlobal])
2136   {
2137     // quality cut on ITS+TPC tracks
2138     fgMultEstTrackCuts[kMultEstTrackCutGlobal] = new AliESDtrackCuts();
2139     fgMultEstTrackCuts[kMultEstTrackCutGlobal]->SetMinNClustersTPC(70);
2140     fgMultEstTrackCuts[kMultEstTrackCutGlobal]->SetMaxChi2PerClusterTPC(4);
2141     fgMultEstTrackCuts[kMultEstTrackCutGlobal]->SetAcceptKinkDaughters(kFALSE);
2142     fgMultEstTrackCuts[kMultEstTrackCutGlobal]->SetRequireTPCRefit(kTRUE);
2143     fgMultEstTrackCuts[kMultEstTrackCutGlobal]->SetRequireITSRefit(kTRUE);
2144     //multiplicity underestimate if we use global tracks with |eta| > 0.9
2145     fgMultEstTrackCuts[kMultEstTrackCutGlobal]->SetEtaRange(-0.9, 0.9);
2146
2147     // quality cut on ITS_SA tracks (complementary to ITS+TPC)
2148     fgMultEstTrackCuts[kMultEstTrackCutITSSA] = new AliESDtrackCuts();
2149     fgMultEstTrackCuts[kMultEstTrackCutITSSA]->SetRequireITSRefit(kTRUE);
2150
2151     // primary selection for tracks with SPD hits
2152     fgMultEstTrackCuts[kMultEstTrackCutDCAwSPD] = new AliESDtrackCuts();
2153     fgMultEstTrackCuts[kMultEstTrackCutDCAwSPD]->SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kAny);
2154     fgMultEstTrackCuts[kMultEstTrackCutDCAwSPD]->SetMaxDCAToVertexXYPtDep("0.0182+0.0350/pt^1.01");
2155     fgMultEstTrackCuts[kMultEstTrackCutDCAwSPD]->SetMaxDCAToVertexZ(2);
2156
2157     // primary selection for tracks w/o SPD hits
2158     fgMultEstTrackCuts[kMultEstTrackCutDCAwoSPD] = new AliESDtrackCuts();
2159     fgMultEstTrackCuts[kMultEstTrackCutDCAwoSPD]->SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kNone);
2160     fgMultEstTrackCuts[kMultEstTrackCutDCAwoSPD]->SetMaxDCAToVertexXYPtDep("1.5*(0.0182+0.0350/pt^1.01)");
2161     fgMultEstTrackCuts[kMultEstTrackCutDCAwoSPD]->SetMaxDCAToVertexZ(2);
2162   }
2163
2164   return fgMultEstTrackCuts[cut];
2165 }
2166
2167 Int_t AliESDtrackCuts::GetReferenceMultiplicity(const AliESDEvent* esd, MultEstTrackType trackType, Float_t etaRange)
2168 {
2169   // Get multiplicity estimate based on TPC/ITS tracks and tracklets
2170   // Adapted for AliESDtrackCuts from a version developed by: Ruben Shahoyan, Anton Alkin, Arvinder Palaha
2171   //
2172   // Returns a negative value if no reliable estimate can be provided:
2173   //   -1 SPD vertex missing
2174   //   -2 SPD VertexerZ dispersion too large
2175   //   -3 Track vertex missing (not checked for kTracklets)
2176   //   -4 Distance between SPD and track vertices too large (not checked for kTracklets)
2177   //
2178   // WARNING This functions does not cut on the z vtx. Depending on the eta range requested, you need to restrict your z vertex range!
2179   //
2180   // Strategy for combined estimators
2181   //   1. Count global tracks and flag them
2182   //   2. Count ITSSA as complementaries for ITSTPC+ or as main tracks
2183   //   3. Count the complementary SPD tracklets
2184
2185   const AliESDVertex* vertices[2];
2186   vertices[0] = esd->GetPrimaryVertexSPD();
2187   vertices[1] = esd->GetPrimaryVertexTracks();
2188
2189   if (!vertices[0]->GetStatus())
2190   {
2191     AliDebugClass(AliLog::kDebug, "No SPD vertex. Not able to make a reliable multiplicity estimate.");
2192     return -1;
2193   }
2194
2195   if (vertices[0]->IsFromVertexerZ() && vertices[0]->GetDispersion() > 0.02)
2196   {
2197     AliDebugClass(AliLog::kDebug, "Vertexer z dispersion > 0.02. Not able to make a reliable multiplicity estimate.");
2198     return -2;
2199   }
2200
2201   Int_t multiplicityEstimate = 0;
2202
2203   // SPD tracklet-only estimate
2204   if (trackType == kTracklets)
2205   {
2206     const AliMultiplicity* spdmult = esd->GetMultiplicity();    // spd multiplicity object
2207     for (Int_t i=0; i<spdmult->GetNumberOfTracklets(); ++i)
2208     {
2209       if (TMath::Abs(spdmult->GetEta(i)) > etaRange)
2210         continue; // eta selection for tracklets
2211       multiplicityEstimate++;
2212     }
2213     return multiplicityEstimate;
2214   }
2215
2216   if (!vertices[1]->GetStatus())
2217   {
2218     AliDebugClass(AliLog::kDebug, "No track vertex. Not able to make a reliable multiplicity estimate.");
2219     return -3;
2220   }
2221
2222   // TODO value of displacement to be studied
2223   const Float_t maxDisplacement = 0.5;
2224   //check for displaced vertices
2225   Double_t displacement = TMath::Abs(vertices[0]->GetZ() - vertices[1]->GetZ());
2226   if (displacement > maxDisplacement)
2227   {
2228     AliDebugClass(AliLog::kDebug, Form("Displaced vertices %f > %f",displacement,maxDisplacement));
2229     return -4;
2230   }
2231
2232   // update eta range in track cuts
2233   GetMultEstTrackCuts(kMultEstTrackCutITSSA)->SetEtaRange(-etaRange, etaRange);
2234   GetMultEstTrackCuts(kMultEstTrackCutDCAwSPD)->SetEtaRange(-etaRange, etaRange);
2235   GetMultEstTrackCuts(kMultEstTrackCutDCAwoSPD)->SetEtaRange(-etaRange, etaRange);
2236
2237   //*******************************************************************************************************
2238   //set counters to initial zeros
2239   Int_t tracksITSTPC = 0;     //number of global tracks for a given event
2240   Int_t tracksITSSA = 0;      //number of ITS standalone tracks for a given event
2241   Int_t tracksITSTPCSA_complementary = 0; //number of ITS standalone tracks complementary to TPC for a given event
2242   Int_t trackletsITSTPC_complementary = 0;//number of SPD tracklets complementary to global/ITSSA tracks for a given events
2243   Int_t trackletsITSSA_complementary = 0; //number of SPD tracklets complementary to ITSSA tracks for a given event
2244
2245   const Int_t nESDTracks = esd->GetNumberOfTracks();
2246
2247   // flags for secondary and rejected tracks
2248   const Int_t kRejBit = BIT(15); // set this bit in global tracks if it is rejected by a cut
2249   const Int_t kSecBit = BIT(16); // set this bit in global tracks if it is secondary according to a cut
2250
2251   for(Int_t itracks=0; itracks < nESDTracks; itracks++) {
2252     esd->GetTrack(itracks)->ResetBit(kSecBit|kRejBit); //reset bits used for flagging secondaries and rejected tracks in case they were changed before this analysis
2253   }
2254   const Int_t maxid = nESDTracks; // used to define bool array for check multiple associations of tracklets to one track. array starts at 0.
2255
2256   // bit mask for esd tracks, to check if multiple tracklets are associated to it
2257   TBits globalBits(maxid), pureITSBits(maxid);
2258   // why labels are used with the data? RS
2259   //  Bool_t globalBits[maxid], pureITSBits[maxid];
2260   //  for(Int_t i=0; i<maxid; i++){ // set all bools to false
2261   //      globalBits[i]=kFALSE;
2262   //      pureITSBits[i]=kFALSE;
2263   //  }
2264
2265   //*******************************************************************************************************
2266   // get multiplicity from global tracks
2267   for(Int_t itracks = 0; itracks < nESDTracks; itracks++) { // flag the tracks
2268     AliESDtrack* track = esd->GetTrack(itracks);
2269
2270     // if track is a secondary from a V0, flag as a secondary
2271     if (track->IsOn(AliESDtrack::kMultInV0)) {
2272       track->SetBit(kSecBit);
2273       continue;
2274     }
2275     /* done via proper DCA cut
2276     //secondary?
2277     if (track->IsOn(AliESDtrack::kMultSec)) {
2278       track->SetBit(kSecBit);
2279       continue;
2280     }
2281     */
2282     // check tracks with ITS part
2283     //*******************************************************************************************************
2284     if (track->IsOn(AliESDtrack::kITSin) && !track->IsOn(AliESDtrack::kITSpureSA) && trackType == kTrackletsITSTPC) { // track has ITS part but is not an ITS_SA
2285       //*******************************************************************************************************
2286       // TPC+ITS
2287       if (track->IsOn(AliESDtrack::kTPCin)) {  // Global track, has ITS and TPC contributions
2288           if (fgMultEstTrackCuts[kMultEstTrackCutGlobal]->AcceptTrack(track)) { // good ITSTPC track
2289             if (fgMultEstTrackCuts[kMultEstTrackCutDCAwSPD]->AcceptTrack(track) || fgMultEstTrackCuts[kMultEstTrackCutDCAwoSPD]->AcceptTrack(track)) {
2290               tracksITSTPC++; //global track counted
2291               globalBits.SetBitNumber(itracks);
2292             }
2293             else track->SetBit(kSecBit); // large DCA -> secondary, don't count either track not associated tracklet
2294           }
2295           else track->SetBit(kRejBit); // bad quality, don't count the track, but may count tracklet if associated
2296       }
2297       //*******************************************************************************************************
2298       // ITS complementary
2299       else if (fgMultEstTrackCuts[kMultEstTrackCutITSSA]->AcceptTrack(track)) { // good ITS complementary track
2300         if (fgMultEstTrackCuts[kMultEstTrackCutDCAwSPD]->AcceptTrack(track) || fgMultEstTrackCuts[kMultEstTrackCutDCAwoSPD]->AcceptTrack(track)) {
2301           tracksITSTPCSA_complementary++;
2302           globalBits.SetBitNumber(itracks);
2303         }
2304         else track->SetBit(kSecBit); // large DCA -> secondary, don't count either track not associated tracklet
2305       }
2306       else track->SetBit(kRejBit); // bad quality, don't count the track, but may count tracklet if associated
2307     }
2308     //*******************************************************************************************************
2309     // check tracks from ITS_SA_PURE
2310     if (track->IsOn(AliESDtrack::kITSin) && track->IsOn(AliESDtrack::kITSpureSA) && trackType == kTrackletsITSSA){
2311       if (fgMultEstTrackCuts[kMultEstTrackCutITSSA]->AcceptTrack(track)) { // good ITSSA track
2312         if (fgMultEstTrackCuts[kMultEstTrackCutDCAwSPD]->AcceptTrack(track) || fgMultEstTrackCuts[kMultEstTrackCutDCAwoSPD]->AcceptTrack(track)) {
2313           tracksITSSA++;
2314           pureITSBits.SetBitNumber(itracks);
2315         }
2316         else track->SetBit(kSecBit);
2317       }
2318       else track->SetBit(kRejBit);
2319     }
2320   }//ESD tracks counted
2321
2322   //*******************************************************************************************************
2323   // get multiplicity from ITS tracklets to complement TPC+ITS, and ITSpureSA
2324   const AliMultiplicity* spdmult = esd->GetMultiplicity();    // spd multiplicity object
2325   for (Int_t i=0; i<spdmult->GetNumberOfTracklets(); ++i) {
2326     if (TMath::Abs(spdmult->GetEta(i)) > etaRange) continue; // eta selection for tracklets
2327
2328     // if counting tracks+tracklets, check if clusters were already used in tracks
2329     Int_t id1, id2, id3, id4;
2330     spdmult->GetTrackletTrackIDs ( i, 0, id1, id2 ); // references for eventual Global/ITS_SA tracks
2331     spdmult->GetTrackletTrackIDs ( i, 1, id3, id4 ); // references for eventual ITS_SA_pure tracks
2332
2333     // are both clusters from the same tracks? If not, skip the tracklet (shouldn't change things much)
2334     if ( ( id1 != id2 && id1 >= 0 && id2 >= 0 ) || ( id3 != id4 && id3 >= 0 && id4 >= 0 ) ) continue;
2335
2336     //referenced track
2337     //at this point we either have id1 = id2 (id3 = id4) or only one of the ids pair is -1
2338     // id1>=0, id2>=0, id1=id2  : tracklet has associated track
2339     // id1>=0, id2 = -1         : 1st layer cluster has associated track
2340     // id1=-1, id2>=0           : 2nd layer cluster has associated track
2341     // id1=-1, id2=-1           : tracklet has no associated track
2342     //
2343     Int_t bUsedInGlobal(-1);
2344     if ( id1 != -1 ) bUsedInGlobal = globalBits.TestBitNumber(id1) ? id1 : -1;
2345     else if ( id2 != -1) bUsedInGlobal = globalBits.TestBitNumber(id2) ? id2 : -1;
2346     Int_t bUsedInPureITS(-1);
2347     if ( id3 != -1 ) bUsedInPureITS = pureITSBits.TestBitNumber(id3) ? id3 : -1;
2348     else if ( id4 != -1) bUsedInPureITS = pureITSBits.TestBitNumber(id4) ? id4 : -1;
2349     //
2350     AliESDtrack* tr_global = bUsedInGlobal >= 0 ? esd->GetTrack ( bUsedInGlobal ) : 0;
2351     AliESDtrack* tr_itssa = bUsedInPureITS >= 0 ? esd->GetTrack ( bUsedInPureITS ) : 0;
2352     //
2353     // has associated pure ITS track been associated to a previous tracklet?
2354     //*******************************************************************************************************
2355     if (trackType == kTrackletsITSTPC) {
2356             //*******************************************************************************************************
2357             // count tracklets towards global+complimentary tracks
2358             if ( ( tr_global && !tr_global->TestBit ( kSecBit ) ) && ( tr_global &&  tr_global->TestBit ( kRejBit ) ) ) {  // count tracklet as bad quality track
2359                     globalBits.SetBitNumber( bUsedInGlobal ); // mark global track linked to this tracklet as used
2360                     ++trackletsITSTPC_complementary;
2361             }
2362
2363             if ( bUsedInGlobal < 0 ) ++trackletsITSTPC_complementary; //no associated track, just count the tracklet
2364     } else {
2365             //*******************************************************************************************************
2366             // count tracklets towards ITS_SA_pure tracks
2367             if ( ( tr_itssa && !tr_itssa->TestBit ( kSecBit ) ) && ( tr_itssa &&  tr_itssa->TestBit ( kRejBit ) ) ) {  // count tracklet as bad quality track
2368                     pureITSBits.SetBitNumber( bUsedInPureITS ); // mark ITS pure SA track linked to this tracklet as used
2369                     ++trackletsITSSA_complementary;
2370             }
2371
2372             if ( bUsedInPureITS < 0 ) ++trackletsITSSA_complementary; //no associated track, just count the tracklet
2373     }
2374   }
2375
2376   //*******************************************************************************************************
2377   if (trackType == kTrackletsITSTPC)
2378     multiplicityEstimate = tracksITSTPC + tracksITSTPCSA_complementary + trackletsITSTPC_complementary;
2379   else
2380     multiplicityEstimate = tracksITSSA + trackletsITSSA_complementary;
2381
2382   return multiplicityEstimate;
2383 }
2384
2385 //____________________________________________________________________
2386 void AliESDtrackCuts::SetRequireStandardTOFmatchCuts(){
2387
2388         // setting the TOF cuts flags (kTOFout = TOF matching distance) to true, to include the selection on the standard TOF matching
2389
2390         SetRequireTOFout(kTRUE);
2391         SetFlagCutTOFdistance(kTRUE);
2392         SetCutTOFdistance(3.);
2393
2394 }
2395