PWG0/esdTrackCuts/AliESDtrackCuts.cxx

   1 #include "AliESDtrackCuts.h"
   2
   3
   4 #include <AliESDtrack.h>
   5 #include <AliESD.h>
   6 #include <AliLog.h>
   7
   8 //____________________________________________________________________
   9 ClassImp(AliESDtrackCuts)
  10
  11 // Cut names
  12 const Char_t* AliESDtrackCuts::fgkCutNames[kNCuts] = {
  13  "require TPC refit",
  14  "require ITS refit",
  15  "n clusters TPC",
  16  "n clusters ITS",
  17  "#Chi^{2}/clusters TPC",
  18  "#Chi^{2}/clusters ITS",
  19  "cov 11",
  20  "cov 22",
  21  "cov 33",
  22  "cov 44",
  23  "cov 55",
  24  "trk-to-vtx",
  25  "trk-to-vtx failed",
  26  "kink daughters",
  27  "p",
  28  "p_{T}",
  29  "p_{x}",
  30  "p_{y}",
  31  "p_{z}",
  32  "y",
  33  "eta"
  34 };
  35
  36 //____________________________________________________________________
  37 AliESDtrackCuts::AliESDtrackCuts() : TObject(),
  38   fCutMinNClusterTPC(0),
  39   fCutMinNClusterITS(0),
  40   fCutMaxChi2PerClusterTPC(0),
  41   fCutMaxChi2PerClusterITS(0),
  42   fCutMaxC11(0),
  43   fCutMaxC22(0),
  44   fCutMaxC33(0),
  45   fCutMaxC44(0),
  46   fCutMaxC55(0),
  47   fCutAcceptKinkDaughters(0),
  48   fCutRequireTPCRefit(0),
  49   fCutRequireITSRefit(0),
  50   fCutNsigmaToVertex(0),
  51   fCutSigmaToVertexRequired(0),
  52   fPMin(0),
  53   fPMax(0),
  54   fPtMin(0),
  55   fPtMax(0),
  56   fPxMin(0),
  57   fPxMax(0),
  58   fPyMin(0),
  59   fPyMax(0),
  60   fPzMin(0),
  61   fPzMax(0),
  62   fEtaMin(0),
  63   fEtaMax(0),
  64   fRapMin(0),
  65   fRapMax(0),
  66   fHistogramsOn(0),
  67   fhCutStatistics(0),
  68   fhCutCorrelation(0)
  69 {
  70   //
  71   // constructor
  72   //
  73
  74   Init();
  75
  76   //##############################################################################
  77   // setting default cuts
  78   SetMinNClustersTPC();
  79   SetMinNClustersITS();
  80   SetMaxChi2PerClusterTPC();
  81   SetMaxChi2PerClusterITS();
  82   SetMaxCovDiagonalElements();
  83   SetRequireTPCRefit();
  84   SetRequireITSRefit();
  85   SetAcceptKingDaughters();
  86   SetMinNsigmaToVertex();
  87   SetRequireSigmaToVertex();
  88   SetPRange();
  89   SetPtRange();
  90   SetPxRange();
  91   SetPyRange();
  92   SetPzRange();
  93   SetEtaRange();
  94   SetRapRange();
  95
  96   SetHistogramsOn();
  97 }
  98
  99 //_____________________________________________________________________________
 100 AliESDtrackCuts::AliESDtrackCuts(const AliESDtrackCuts &c) : TObject(c),
 101   fCutMinNClusterTPC(0),
 102   fCutMinNClusterITS(0),
 103   fCutMaxChi2PerClusterTPC(0),
 104   fCutMaxChi2PerClusterITS(0),
 105   fCutMaxC11(0),
 106   fCutMaxC22(0),
 107   fCutMaxC33(0),
 108   fCutMaxC44(0),
 109   fCutMaxC55(0),
 110   fCutAcceptKinkDaughters(0),
 111   fCutRequireTPCRefit(0),
 112   fCutRequireITSRefit(0),
 113   fCutNsigmaToVertex(0),
 114   fCutSigmaToVertexRequired(0),
 115   fPMin(0),
 116   fPMax(0),
 117   fPtMin(0),
 118   fPtMax(0),
 119   fPxMin(0),
 120   fPxMax(0),
 121   fPyMin(0),
 122   fPyMax(0),
 123   fPzMin(0),
 124   fPzMax(0),
 125   fEtaMin(0),
 126   fEtaMax(0),
 127   fRapMin(0),
 128   fRapMax(0),
 129   fHistogramsOn(0),
 130   fhCutStatistics(0),
 131   fhCutCorrelation(0)
 132 {
 133   //
 134   // copy constructor
 135   //
 136
 137   ((AliESDtrackCuts &) c).Copy(*this);
 138 }
 139
 140 AliESDtrackCuts::~AliESDtrackCuts()
 141 {
 142   //
 143   // destructor
 144   //
 145
 146   // ## TODO to be implemented
 147 }
 148
 149 void AliESDtrackCuts::Init()
 150 {
 151   //
 152   // sets everything to zero
 153   //
 154
 155   fCutMinNClusterTPC = 0;
 156   fCutMinNClusterITS = 0;
 157
 158   fCutMaxChi2PerClusterTPC = 0;
 159   fCutMaxChi2PerClusterITS = 0;
 160
 161   fCutMaxC11 = 0;
 162   fCutMaxC22 = 0;
 163   fCutMaxC33 = 0;
 164   fCutMaxC44 = 0;
 165   fCutMaxC55 = 0;
 166
 167   fCutAcceptKinkDaughters = 0;
 168   fCutRequireTPCRefit = 0;
 169   fCutRequireITSRefit = 0;
 170
 171   fCutNsigmaToVertex = 0;
 172   fCutSigmaToVertexRequired = 0;
 173
 174   fPMin = 0;
 175   fPMax = 0;
 176   fPtMin = 0;
 177   fPtMax = 0;
 178   fPxMin = 0;
 179   fPxMax = 0;
 180   fPyMin = 0;
 181   fPyMax = 0;
 182   fPzMin = 0;
 183   fPzMax = 0;
 184   fEtaMin = 0;
 185   fEtaMax = 0;
 186   fRapMin = 0;
 187   fRapMax = 0;
 188
 189   fHistogramsOn = kFALSE;
 190
 191   for (Int_t i=0; i<2; ++i)
 192   {
 193     fhNClustersITS[i] = 0;
 194     fhNClustersTPC[i] = 0;
 195
 196     fhChi2PerClusterITS[i] = 0;
 197     fhChi2PerClusterTPC[i] = 0;
 198
 199     fhC11[i] = 0;
 200     fhC22[i] = 0;
 201     fhC33[i] = 0;
 202     fhC44[i] = 0;
 203     fhC55[i] = 0;
 204
 205     fhDXY[i] = 0;
 206     fhDZ[i] = 0;
 207     fhDXYvsDZ[i] = 0;
 208
 209     fhDXYNormalized[i] = 0;
 210     fhDZNormalized[i] = 0;
 211     fhDXYvsDZNormalized[i] = 0;
 212   }
 213
 214   fhCutStatistics = 0;
 215   fhCutCorrelation = 0;
 216 }
 217
 218 //_____________________________________________________________________________
 219 AliESDtrackCuts &AliESDtrackCuts::operator=(const AliESDtrackCuts &c)
 220 {
 221   //
 222   // Assignment operator
 223   //
 224
 225   if (this != &c) ((AliESDtrackCuts &) c).Copy(*this);
 226   return *this;
 227 }
 228
 229 //_____________________________________________________________________________
 230 void AliESDtrackCuts::Copy(TObject &c) const
 231 {
 232   //
 233   // Copy function
 234   //
 235
 236   AliESDtrackCuts& target = (AliESDtrackCuts &) c;
 237
 238   target.Init();
 239
 240   target.fCutMinNClusterTPC = fCutMinNClusterTPC;
 241   target.fCutMinNClusterITS = fCutMinNClusterITS;
 242
 243   target.fCutMaxChi2PerClusterTPC = fCutMaxChi2PerClusterTPC;
 244   target.fCutMaxChi2PerClusterITS = fCutMaxChi2PerClusterITS;
 245
 246   target.fCutMaxC11 = fCutMaxC11;
 247   target.fCutMaxC22 = fCutMaxC22;
 248   target.fCutMaxC33 = fCutMaxC33;
 249   target.fCutMaxC44 = fCutMaxC44;
 250   target.fCutMaxC55 = fCutMaxC55;
 251
 252   target.fCutAcceptKinkDaughters = fCutAcceptKinkDaughters;
 253   target.fCutRequireTPCRefit = fCutRequireTPCRefit;
 254   target.fCutRequireITSRefit = fCutRequireITSRefit;
 255
 256   target.fCutNsigmaToVertex = fCutNsigmaToVertex;
 257   target.fCutSigmaToVertexRequired = fCutSigmaToVertexRequired;
 258
 259   target.fPMin = fPMin;
 260   target.fPMax = fPMax;
 261   target.fPtMin = fPtMin;
 262   target.fPtMax = fPtMax;
 263   target.fPxMin = fPxMin;
 264   target.fPxMax = fPxMax;
 265   target.fPyMin = fPyMin;
 266   target.fPyMax = fPyMax;
 267   target.fPzMin = fPzMin;
 268   target.fPzMax = fPzMax;
 269   target.fEtaMin = fEtaMin;
 270   target.fEtaMax = fEtaMax;
 271   target.fRapMin = fRapMin;
 272   target.fRapMax = fRapMax;
 273
 274   target.fHistogramsOn = fHistogramsOn;
 275
 276   for (Int_t i=0; i<2; ++i)
 277   {
 278     if (fhNClustersITS[i]) target.fhNClustersITS[i] = (TH1F*) fhNClustersITS[i]->Clone();
 279     if (fhNClustersTPC[i]) target.fhNClustersTPC[i] = (TH1F*) fhNClustersTPC[i]->Clone();
 280
 281     if (fhChi2PerClusterITS[i]) target.fhChi2PerClusterITS[i] = (TH1F*) fhChi2PerClusterITS[i]->Clone();
 282     if (fhChi2PerClusterTPC[i]) target.fhChi2PerClusterTPC[i] = (TH1F*) fhChi2PerClusterTPC[i]->Clone();
 283
 284     if (fhC11[i]) target.fhC11[i] = (TH1F*) fhC11[i]->Clone();
 285     if (fhC22[i]) target.fhC22[i] = (TH1F*) fhC22[i]->Clone();
 286     if (fhC33[i]) target.fhC33[i] = (TH1F*) fhC33[i]->Clone();
 287     if (fhC44[i]) target.fhC44[i] = (TH1F*) fhC44[i]->Clone();
 288     if (fhC55[i]) target.fhC55[i] = (TH1F*) fhC55[i]->Clone();
 289
 290     if (fhDXY[i]) target.fhDXY[i] = (TH1F*) fhDXY[i]->Clone();
 291     if (fhDZ[i]) target.fhDZ[i] = (TH1F*) fhDZ[i]->Clone();
 292     if (fhDXYvsDZ[i]) target.fhDXYvsDZ[i] = (TH2F*) fhDXYvsDZ[i]->Clone();
 293
 294     if (fhDXYNormalized[i]) target.fhDXYNormalized[i] = (TH1F*) fhDXYNormalized[i]->Clone();
 295     if (fhDZNormalized[i]) target.fhDZNormalized[i] = (TH1F*) fhDZNormalized[i]->Clone();
 296     if (fhDXYvsDZNormalized[i]) target.fhDXYvsDZNormalized[i] = (TH2F*) fhDXYvsDZNormalized[i]->Clone();
 297   }
 298
 299   if (fhCutStatistics) target.fhCutStatistics = (TH1F*) fhCutStatistics->Clone();
 300   if (fhCutCorrelation) target.fhCutCorrelation = (TH2F*) fhCutCorrelation->Clone();
 301
 302   TObject::Copy(c);
 303 }
 304
 305 //____________________________________________________________________
 306 Float_t AliESDtrackCuts::GetSigmaToVertex(AliESDtrack* esdTrack)
 307 {
 308   //
 309
 310   Float_t b[2];
 311   Float_t bRes[2];
 312   Float_t bCov[3];
 313   esdTrack->GetImpactParameters(b,bCov);
 314   if (bCov[0]<=0 || bCov[2]<=0) {
 315     AliDebug(1, "Estimated b resolution lower or equal zero!");
 316     bCov[0]=0; bCov[2]=0;
 317   }
 318   bRes[0] = TMath::Sqrt(bCov[0]);
 319   bRes[1] = TMath::Sqrt(bCov[2]);
 320
 321   // -----------------------------------
 322   // How to get to a n-sigma cut?
 323   //
 324   // The accumulated statistics from 0 to d is
 325   //
 326   // ->  Erf(d/Sqrt(2)) for a 1-dim gauss (d = n_sigma)
 327   // ->  1 - Exp(-d**2) for a 2-dim gauss (d*d = dx*dx + dy*dy != n_sigma)
 328   //
 329   // It means that for a 2-dim gauss: n_sigma(d) = Sqrt(2)*ErfInv(1 - Exp((-x**2)/2)
 330   // Can this be expressed in a different way?
 331   //
 332   //
 333   // FIX: I don't think this is correct!!! Keeping d as n_sigma for now...
 334
 335   if (bRes[0] == 0 || bRes[1] ==0)
 336     return -1;
 337
 338   Float_t d = TMath::Sqrt(TMath::Power(b[0]/bRes[0],2) + TMath::Power(b[1]/bRes[1],2));
 339
 340   // stupid rounding problem screws up everything:
 341   // if d is too big, TMath::Exp(...) gets 0, and TMath::ErfInverse(1) that should be infinite, gets 0 :(
 342   if (TMath::Exp(-d * d / 2) < 1e-10)
 343     return 1000;
 344
 345   d = TMath::ErfInverse(1 - TMath::Exp(-d * d / 2)) * TMath::Sqrt(2);
 346   return d;
 347 }
 348
 349 //____________________________________________________________________
 350 Bool_t
 351 AliESDtrackCuts::AcceptTrack(AliESDtrack* esdTrack) {
 352   //
 353   // figure out if the tracks survives all the track cuts defined
 354   //
 355   // the different quality parameter and kinematic values are first
 356   // retrieved from the track. then it is found out what cuts the
 357   // track did not survive and finally the cuts are imposed.
 358
 359
 360
 361   UInt_t status = esdTrack->GetStatus();
 362
 363   // dummy array
 364   Int_t  fIdxInt[200];
 365
 366   // getting quality parameters from the ESD track
 367   Int_t nClustersITS = esdTrack->GetITSclusters(fIdxInt);
 368   Int_t nClustersTPC = esdTrack->GetTPCclusters(fIdxInt);
 369
 370
 371
 372   Float_t chi2PerClusterITS = -1;
 373   Float_t chi2PerClusterTPC = -1;
 374   if (nClustersITS!=0)
 375     chi2PerClusterITS = esdTrack->GetITSchi2()/Float_t(nClustersITS);
 376   if (nClustersTPC!=0)
 377     chi2PerClusterTPC = esdTrack->GetTPCchi2()/Float_t(nClustersTPC);
 378
 379   Double_t extCov[15];
 380   esdTrack->GetExternalCovariance(extCov);
 381
 382   // getting the track to vertex parameters
 383   Float_t nSigmaToVertex = GetSigmaToVertex(esdTrack);
 384
 385   // getting the kinematic variables of the track
 386   // (assuming the mass is known)
 387   Double_t p[3];
 388   esdTrack->GetPxPyPz(p);
 389   Float_t momentum = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2) + TMath::Power(p[2],2));
 390   Float_t pt       = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2));
 391   Float_t energy   = TMath::Sqrt(TMath::Power(esdTrack->GetMass(),2) + TMath::Power(momentum,2));
 392
 393
 394   //y-eta related calculations
 395   Float_t eta = -100.;
 396   Float_t y   = -100.;
 397   if((momentum != TMath::Abs(p[2]))&&(momentum != 0))
 398     eta = 0.5*TMath::Log((momentum + p[2])/(momentum - p[2]));
 399   if((energy != TMath::Abs(p[2]))&&(momentum != 0))
 400     y = 0.5*TMath::Log((energy + p[2])/(energy - p[2]));
 401
 402
 403   //########################################################################
 404   // cut the track?
 405
 406   Bool_t cuts[kNCuts];
 407   for (Int_t i=0; i<kNCuts; i++) cuts[i]=kFALSE;
 408
 409   // track quality cuts
 410   if (fCutRequireTPCRefit && (status&AliESDtrack::kTPCrefit)==0)
 411     cuts[0]=kTRUE;
 412   if (fCutRequireITSRefit && (status&AliESDtrack::kITSrefit)==0)
 413     cuts[1]=kTRUE;
 414   if (nClustersTPC<fCutMinNClusterTPC)
 415     cuts[2]=kTRUE;
 416   if (nClustersITS<fCutMinNClusterITS)
 417     cuts[3]=kTRUE;
 418   if (chi2PerClusterTPC>fCutMaxChi2PerClusterTPC)
 419     cuts[4]=kTRUE;
 420   if (chi2PerClusterITS>fCutMaxChi2PerClusterITS)
 421     cuts[5]=kTRUE;
 422   if (extCov[0]  > fCutMaxC11)
 423     cuts[6]=kTRUE;
 424   if (extCov[2]  > fCutMaxC22)
 425     cuts[7]=kTRUE;
 426   if (extCov[5]  > fCutMaxC33)
 427     cuts[8]=kTRUE;
 428   if (extCov[9]  > fCutMaxC44)
 429     cuts[9]=kTRUE;
 430   if (extCov[14]  > fCutMaxC55)
 431     cuts[10]=kTRUE;
 432   if (nSigmaToVertex > fCutNsigmaToVertex)
 433     cuts[11] = kTRUE;
 434   // if n sigma could not be calculated
 435   if (nSigmaToVertex<0 && fCutSigmaToVertexRequired)
 436     cuts[12]=kTRUE;
 437   if (!fCutAcceptKinkDaughters && esdTrack->GetKinkIndex(0)>0)
 438     cuts[13]=kTRUE;
 439   // track kinematics cut
 440   if((momentum < fPMin) || (momentum > fPMax))
 441     cuts[14]=kTRUE;
 442   if((pt < fPtMin) || (pt > fPtMax))
 443     cuts[15] = kTRUE;
 444   if((p[0] < fPxMin) || (p[0] > fPxMax))
 445     cuts[16] = kTRUE;
 446   if((p[1] < fPyMin) || (p[1] > fPyMax))
 447     cuts[17] = kTRUE;
 448   if((p[2] < fPzMin) || (p[2] > fPzMax))
 449     cuts[18] = kTRUE;
 450   if((eta < fEtaMin) || (eta > fEtaMax))
 451     cuts[19] = kTRUE;
 452   if((y < fRapMin) || (y > fRapMax))
 453     cuts[20] = kTRUE;
 454
 455   Bool_t cut=kFALSE;
 456   for (Int_t i=0; i<kNCuts; i++)
 457     if (cuts[i]) cut = kTRUE;
 458
 459   //########################################################################
 460   // filling histograms
 461   if (fHistogramsOn) {
 462     fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n tracks")));
 463
 464     if (cut)
 465       fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n cut tracks")));
 466
 467     for (Int_t i=0; i<kNCuts; i++) {
 468       if (cuts[i])
 469         fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin(fgkCutNames[i])));
 470
 471       for (Int_t j=i; j<kNCuts; j++) {
 472         if (cuts[i] && cuts[j]) {
 473           Float_t x = fhCutCorrelation->GetXaxis()->GetBinCenter(fhCutCorrelation->GetXaxis()->FindBin(fgkCutNames[i]));
 474           Float_t y = fhCutCorrelation->GetYaxis()->GetBinCenter(fhCutCorrelation->GetYaxis()->FindBin(fgkCutNames[j]));
 475           fhCutCorrelation->Fill(x,y);
 476         }
 477       }
 478     }
 479
 480
 481     fhNClustersITS[0]->Fill(nClustersITS);
 482     fhNClustersTPC[0]->Fill(nClustersTPC);
 483     fhChi2PerClusterITS[0]->Fill(chi2PerClusterITS);
 484     fhChi2PerClusterTPC[0]->Fill(chi2PerClusterTPC);
 485
 486     fhC11[0]->Fill(extCov[0]);
 487     fhC22[0]->Fill(extCov[2]);
 488     fhC33[0]->Fill(extCov[5]);
 489     fhC44[0]->Fill(extCov[9]);
 490     fhC55[0]->Fill(extCov[14]);
 491
 492     Float_t b[2];
 493     Float_t bRes[2];
 494     Float_t bCov[3];
 495     esdTrack->GetImpactParameters(b,bCov);
 496     if (bCov[0]<=0 || bCov[2]<=0) {
 497       AliDebug(1, "Estimated b resolution lower or equal zero!");
 498       bCov[0]=0; bCov[2]=0;
 499     }
 500     bRes[0] = TMath::Sqrt(bCov[0]);
 501     bRes[1] = TMath::Sqrt(bCov[2]);
 502
 503     fhDZ[0]->Fill(b[1]);
 504     fhDXY[0]->Fill(b[0]);
 505     fhDXYvsDZ[0]->Fill(b[1],b[0]);
 506
 507     if (bRes[0]!=0 && bRes[1]!=0) {
 508       fhDZNormalized[0]->Fill(b[1]/bRes[1]);
 509       fhDXYNormalized[0]->Fill(b[0]/bRes[0]);
 510       fhDXYvsDZNormalized[0]->Fill(b[1]/bRes[1], b[0]/bRes[0]);
 511     }
 512   }
 513
 514   //########################################################################
 515   // cut the track!
 516   if (cut) return kFALSE;
 517
 518   //########################################################################
 519   // filling histograms after cut
 520   if (fHistogramsOn) {
 521     fhNClustersITS[1]->Fill(nClustersITS);
 522     fhNClustersTPC[1]->Fill(nClustersTPC);
 523     fhChi2PerClusterITS[1]->Fill(chi2PerClusterITS);
 524     fhChi2PerClusterTPC[1]->Fill(chi2PerClusterTPC);
 525
 526     fhC11[1]->Fill(extCov[0]);
 527     fhC22[1]->Fill(extCov[2]);
 528     fhC33[1]->Fill(extCov[5]);
 529     fhC44[1]->Fill(extCov[9]);
 530     fhC55[1]->Fill(extCov[14]);
 531
 532     Float_t b[2];
 533     Float_t bRes[2];
 534     Float_t bCov[3];
 535     esdTrack->GetImpactParameters(b,bCov);
 536     if (bCov[0]<=0 || bCov[2]<=0) {
 537       AliDebug(1, "Estimated b resolution lower or equal zero!");
 538       bCov[0]=0; bCov[2]=0;
 539     }
 540     bRes[0] = TMath::Sqrt(bCov[0]);
 541     bRes[1] = TMath::Sqrt(bCov[2]);
 542
 543     fhDZ[1]->Fill(b[1]);
 544     fhDXY[1]->Fill(b[0]);
 545     fhDXYvsDZ[1]->Fill(b[1],b[0]);
 546
 547     if (bRes[0]!=0 && bRes[1]!=0)
 548     {
 549       fhDZNormalized[1]->Fill(b[1]/bRes[1]);
 550       fhDXYNormalized[1]->Fill(b[0]/bRes[0]);
 551       fhDXYvsDZNormalized[1]->Fill(b[1]/bRes[1], b[0]/bRes[0]);
 552     }
 553   }
 554
 555   return kTRUE;
 556 }
 557
 558 //____________________________________________________________________
 559 TObjArray*
 560 AliESDtrackCuts::GetAcceptedTracks(AliESD* esd)
 561 {
 562   //
 563   // returns an array of all tracks that pass the cuts
 564   //
 565
 566   TObjArray* acceptedTracks = new TObjArray();
 567
 568   // loop over esd tracks
 569   for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
 570     AliESDtrack* track = esd->GetTrack(iTrack);
 571
 572     if (AcceptTrack(track))
 573       acceptedTracks->Add(track);
 574   }
 575
 576   return acceptedTracks;
 577 }
 578
 579 //____________________________________________________________________
 580 Int_t
 581 AliESDtrackCuts::CountAcceptedTracks(AliESD* esd)
 582 {
 583   //
 584   // returns an the number of tracks that pass the cuts
 585   //
 586
 587   Int_t count = 0;
 588
 589   // loop over esd tracks
 590   for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
 591     AliESDtrack* track = esd->GetTrack(iTrack);
 592
 593     if (AcceptTrack(track))
 594       count++;
 595   }
 596
 597   return count;
 598 }
 599
 600 //____________________________________________________________________
 601  void AliESDtrackCuts::DefineHistograms(Int_t color) {
 602    //
 603    // diagnostics histograms are defined
 604    //
 605
 606    fHistogramsOn=kTRUE;
 607
 608    //###################################################################################
 609    // defining histograms
 610
 611    fhCutStatistics = new TH1F("cut_statistics","cut statistics",kNCuts+4,-0.5,kNCuts+3.5);
 612
 613    fhCutStatistics->GetXaxis()->SetBinLabel(1,"n tracks");
 614    fhCutStatistics->GetXaxis()->SetBinLabel(2,"n cut tracks");
 615
 616    fhCutCorrelation = new TH2F("cut_correlation","cut correlation",kNCuts,-0.5,kNCuts-0.5,kNCuts,-0.5,kNCuts-0.5);;
 617
 618    for (Int_t i=0; i<kNCuts; i++) {
 619      fhCutStatistics->GetXaxis()->SetBinLabel(i+4,fgkCutNames[i]);
 620      fhCutCorrelation->GetXaxis()->SetBinLabel(i+1,fgkCutNames[i]);
 621      fhCutCorrelation->GetYaxis()->SetBinLabel(i+1,fgkCutNames[i]);
 622    }
 623
 624   fhCutStatistics  ->SetLineColor(color);
 625   fhCutCorrelation ->SetLineColor(color);
 626   fhCutStatistics  ->SetLineWidth(2);
 627   fhCutCorrelation ->SetLineWidth(2);
 628
 629   Char_t str[256];
 630   for (Int_t i=0; i<2; i++) {
 631     if (i==0) sprintf(str," ");
 632     else sprintf(str,"_cut");
 633
 634     fhNClustersITS[i]        = new TH1F(Form("nClustersITS%s",str),"",8,-0.5,7.5);
 635     fhNClustersTPC[i]        = new TH1F(Form("nClustersTPC%s",str),"",165,-0.5,164.5);
 636     fhChi2PerClusterITS[i]   = new TH1F(Form("chi2PerClusterITS%s",str),"",500,0,10);
 637     fhChi2PerClusterTPC[i]   = new TH1F(Form("chi2PerClusterTPC%s",str),"",500,0,10);
 638
 639     fhC11[i]                 = new TH1F(Form("covMatrixDiagonal11%s",str),"",1000,0,5);
 640     fhC22[i]                 = new  TH1F(Form("covMatrixDiagonal22%s",str),"",1000,0,5);
 641     fhC33[i]                 = new  TH1F(Form("covMatrixDiagonal33%s",str),"",1000,0,0.5);
 642     fhC44[i]                 = new  TH1F(Form("covMatrixDiagonal44%s",str),"",1000,0,5);
 643     fhC55[i]                 = new  TH1F(Form("covMatrixDiagonal55%s",str),"",1000,0,5);
 644
 645     fhDXY[i]                 = new  TH1F(Form("dXY%s",str),"",500,-10,10);
 646     fhDZ[i]                  = new  TH1F(Form("dZ%s",str),"",500,-10,10);
 647     fhDXYvsDZ[i]             = new  TH2F(Form("dXYvsDZ%s",str),"",200,-10,10,200,-10,10);
 648
 649     fhDXYNormalized[i]       = new  TH1F(Form("dXYNormalized%s",str),"",500,-10,10);
 650     fhDZNormalized[i]        = new  TH1F(Form("dZNormalized%s",str),"",500,-10,10);
 651     fhDXYvsDZNormalized[i]   = new  TH2F(Form("dXYvsDZNormalized%s",str),"",200,-10,10,200,-10,10);
 652
 653
 654     fhNClustersITS[i]->SetXTitle("n ITS clusters");
 655     fhNClustersTPC[i]->SetXTitle("n TPC clusters");
 656     fhChi2PerClusterITS[i]->SetXTitle("#Chi^{2} per ITS cluster");
 657     fhChi2PerClusterTPC[i]->SetXTitle("#Chi^{2} per TPC cluster");
 658
 659     fhC11[i]->SetXTitle("cov 11 : #sigma_{y}^{2} [cm^{2}]");
 660     fhC22[i]->SetXTitle("cov 22 : #sigma_{z}^{2} [cm^{2}]");
 661     fhC33[i]->SetXTitle("cov 33 : #sigma_{sin(#phi)}^{2}");
 662     fhC44[i]->SetXTitle("cov 44 : #sigma_{tan(#theta_{dip})}^{2}");
 663     fhC55[i]->SetXTitle("cov 55 : #sigma_{1/p_{T}}^{2} [(c/GeV)^2]");
 664
 665     fhDXY[i]->SetXTitle("transverse impact parameter");
 666     fhDZ[i]->SetXTitle("longitudinal impact parameter");
 667     fhDXYvsDZ[i]->SetXTitle("longitudinal impact parameter");
 668     fhDXYvsDZ[i]->SetYTitle("transverse impact parameter");
 669
 670     fhDXYNormalized[i]->SetXTitle("normalized trans impact par");
 671     fhDZNormalized[i]->SetXTitle("normalized long impact par");
 672     fhDXYvsDZNormalized[i]->SetXTitle("normalized long impact par");
 673     fhDXYvsDZNormalized[i]->SetYTitle("normalized trans impact par");
 674
 675     fhNClustersITS[i]->SetLineColor(color);   fhNClustersITS[i]->SetLineWidth(2);
 676     fhNClustersTPC[i]->SetLineColor(color);   fhNClustersTPC[i]->SetLineWidth(2);
 677     fhChi2PerClusterITS[i]->SetLineColor(color);   fhChi2PerClusterITS[i]->SetLineWidth(2);
 678     fhChi2PerClusterTPC[i]->SetLineColor(color);   fhChi2PerClusterTPC[i]->SetLineWidth(2);
 679
 680     fhC11[i]->SetLineColor(color);   fhC11[i]->SetLineWidth(2);
 681     fhC22[i]->SetLineColor(color);   fhC22[i]->SetLineWidth(2);
 682     fhC33[i]->SetLineColor(color);   fhC33[i]->SetLineWidth(2);
 683     fhC44[i]->SetLineColor(color);   fhC44[i]->SetLineWidth(2);
 684     fhC55[i]->SetLineColor(color);   fhC55[i]->SetLineWidth(2);
 685
 686     fhDXY[i]->SetLineColor(color);   fhDXY[i]->SetLineWidth(2);
 687     fhDZ[i]->SetLineColor(color);   fhDZ[i]->SetLineWidth(2);
 688
 689     fhDXYNormalized[i]->SetLineColor(color);   fhDXYNormalized[i]->SetLineWidth(2);
 690     fhDZNormalized[i]->SetLineColor(color);   fhDZNormalized[i]->SetLineWidth(2);
 691   }
 692 }
 693
 694 //____________________________________________________________________
 695 void
 696 AliESDtrackCuts::Print(const Option_t*) const {
 697   //
 698   // print method - still to be implemented
 699   //
 700
 701   AliInfo("AliESDtrackCuts...");
 702 }
 703
 704
 705 //____________________________________________________________________
 706 void AliESDtrackCuts::SaveHistograms(Char_t* dir) {
 707   //
 708   // saves the histograms in a directory (dir)
 709   //
 710
 711
 712   if (!fHistogramsOn) {
 713     AliDebug(0, "Histograms not on - cannot save histograms!!!");
 714     return;
 715   }
 716
 717   gDirectory->mkdir(dir);
 718   gDirectory->cd(dir);
 719
 720   gDirectory->mkdir("before_cuts");
 721   gDirectory->mkdir("after_cuts");
 722
 723   fhCutStatistics->Write();
 724   fhCutCorrelation->Write();
 725
 726   for (Int_t i=0; i<2; i++) {
 727     if (i==0)
 728       gDirectory->cd("before_cuts");
 729     else
 730       gDirectory->cd("after_cuts");
 731
 732     fhNClustersITS[i]        ->Write();
 733     fhNClustersTPC[i]        ->Write();
 734     fhChi2PerClusterITS[i]   ->Write();
 735     fhChi2PerClusterTPC[i]   ->Write();
 736
 737     fhC11[i]                 ->Write();
 738     fhC22[i]                 ->Write();
 739     fhC33[i]                 ->Write();
 740     fhC44[i]                 ->Write();
 741     fhC55[i]                 ->Write();
 742
 743     fhDXY[i]                 ->Write();
 744     fhDZ[i]                  ->Write();
 745     fhDXYvsDZ[i]             ->Write();
 746
 747     fhDXYNormalized[i]       ->Write();
 748     fhDZNormalized[i]        ->Write();
 749     fhDXYvsDZNormalized[i]   ->Write();
 750
 751     gDirectory->cd("../");
 752   }
 753
 754   gDirectory->cd("../");
 755 }
 756
 757
 758