CORRFW/AliCFTrackIsPrimaryCuts.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 // The class AliCFTrackIsPrimaryCut is designed to select reconstructed tracks
  17 // with a small impact parameter and tracks which are (not) daughters of kink
  18 // decays and to provide corresponding QA histograms.
  19 // This class inherits from the Analysis' Framework abstract base class
  20 // AliAnalysisCuts and is a part of the Correction Framework.
  21 // This class acts on single, reconstructed tracks, it is applicable on
  22 // ESD and AOD data.
  23 // It mainly consists of a IsSelected function that returns a boolean.
  24 // This function checks whether the considered track passes a set of cuts:
  25 // - distance to main vertex in units of sigma (resolution)
  26 // - require that the dca calculation doesn't fail
  27 // - accept or not accept daughter tracks of kink decays
  28 //
  29 // The cut values for these cuts are set with the corresponding set functions.
  30 // All cut classes provided by the correction framework are supposed to be
  31 // added in the Analysis Framwork's class AliAnalysisFilter and applied by
  32 // the filter via a loop.
  33 //
  34 // author: I. Kraus (Ingrid.Kraus@cern.ch)
  35 // idea taken form
  36 // AliESDtrackCuts writte by Jan Fiete Grosse-Oetringhaus and
  37 // AliRsnDaughterCut class written by A. Pulvirenti.
  38
  39 #include <TCanvas.h>
  40 #include <TDirectory.h>
  41 #include <TH2.h>
  42 #include <TBits.h>
  43
  44 #include <AliESDtrack.h>
  45 #include <AliLog.h>
  46 #include "AliCFTrackIsPrimaryCuts.h"
  47
  48 ClassImp(AliCFTrackIsPrimaryCuts)
  49
  50 //__________________________________________________________________________________
  51 AliCFTrackIsPrimaryCuts::AliCFTrackIsPrimaryCuts() :
  52   AliCFCutBase(),
  53   fNSigmaToVertex(0),
  54   fNSigmaToVertexMax(0),
  55   fRequireSigmaToVertex(0),
  56   fAODType(AliAODTrack::kUndef),
  57   fAcceptKinkDaughters(0),
  58   fhCutStatistics(0),
  59   fhCutCorrelation(0),
  60   fBitmap(0x0),
  61   fhNBinsNSigma(0),
  62   fhNBinsRequireSigma(0),
  63   fhNBinsAcceptKink(0),
  64   fhNBinsDcaXY(0),
  65   fhNBinsDcaZ(0),
  66   fhNBinsDcaXYnorm(0),
  67   fhNBinsDcaZnorm(0),
  68   fhBinLimNSigma(0x0),
  69   fhBinLimRequireSigma(0x0),
  70   fhBinLimAcceptKink(0x0),
  71   fhBinLimDcaXY(0x0),
  72   fhBinLimDcaZ(0x0),
  73   fhBinLimDcaXYnorm(0x0),
  74   fhBinLimDcaZnorm(0x0)
  75 {
  76   //
  77   // Default constructor
  78   //
  79   Initialise();
  80 }
  81 //__________________________________________________________________________________
  82 AliCFTrackIsPrimaryCuts::AliCFTrackIsPrimaryCuts(Char_t* name, Char_t* title) :
  83   AliCFCutBase(name,title),
  84   fNSigmaToVertex(0),
  85   fNSigmaToVertexMax(0),
  86   fRequireSigmaToVertex(0),
  87   fAODType(AliAODTrack::kUndef),
  88   fAcceptKinkDaughters(0),
  89   fhCutStatistics(0),
  90   fhCutCorrelation(0),
  91   fBitmap(0x0),
  92   fhNBinsNSigma(0),
  93   fhNBinsRequireSigma(0),
  94   fhNBinsAcceptKink(0),
  95   fhNBinsDcaXY(0),
  96   fhNBinsDcaZ(0),
  97   fhNBinsDcaXYnorm(0),
  98   fhNBinsDcaZnorm(0),
  99   fhBinLimNSigma(0x0),
 100   fhBinLimRequireSigma(0x0),
 101   fhBinLimAcceptKink(0x0),
 102   fhBinLimDcaXY(0x0),
 103   fhBinLimDcaZ(0x0),
 104   fhBinLimDcaXYnorm(0x0),
 105   fhBinLimDcaZnorm(0x0)
 106 {
 107   //
 108   // Constructor
 109   //
 110   Initialise();
 111 }
 112 //__________________________________________________________________________________
 113 AliCFTrackIsPrimaryCuts::AliCFTrackIsPrimaryCuts(const AliCFTrackIsPrimaryCuts& c) :
 114   AliCFCutBase(c),
 115   fNSigmaToVertex(c.fNSigmaToVertex),
 116   fNSigmaToVertexMax(c.fNSigmaToVertexMax),
 117   fRequireSigmaToVertex(c.fRequireSigmaToVertex),
 118   fAODType(c.fAODType),
 119   fAcceptKinkDaughters(c.fAcceptKinkDaughters),
 120   fhCutStatistics(c.fhCutStatistics),
 121   fhCutCorrelation(c.fhCutCorrelation),
 122   fBitmap(c.fBitmap),
 123   fhNBinsNSigma(c.fhNBinsNSigma),
 124   fhNBinsRequireSigma(c.fhNBinsRequireSigma),
 125   fhNBinsAcceptKink(c.fhNBinsAcceptKink),
 126   fhNBinsDcaXY(c.fhNBinsDcaXY),
 127   fhNBinsDcaZ(c.fhNBinsDcaZ),
 128   fhNBinsDcaXYnorm(c.fhNBinsDcaXYnorm),
 129   fhNBinsDcaZnorm(c.fhNBinsDcaZnorm),
 130   fhBinLimNSigma(c.fhBinLimNSigma),
 131   fhBinLimRequireSigma(c.fhBinLimRequireSigma),
 132   fhBinLimAcceptKink(c.fhBinLimAcceptKink),
 133   fhBinLimDcaXY(c.fhBinLimDcaXY),
 134   fhBinLimDcaZ(c.fhBinLimDcaZ),
 135   fhBinLimDcaXYnorm(c.fhBinLimDcaXYnorm),
 136   fhBinLimDcaZnorm(c.fhBinLimDcaZnorm)
 137 {
 138   //
 139   // copy constructor
 140   //
 141   ((AliCFTrackIsPrimaryCuts &) c).Copy(*this);
 142 }
 143 //__________________________________________________________________________________
 144 AliCFTrackIsPrimaryCuts& AliCFTrackIsPrimaryCuts::operator=(const AliCFTrackIsPrimaryCuts& c)
 145 {
 146   //
 147   // Assignment operator
 148   //
 149   if (this != &c) {
 150     AliCFCutBase::operator=(c) ;
 151     fNSigmaToVertex = c.fNSigmaToVertex ;
 152     fNSigmaToVertexMax = c.fNSigmaToVertexMax ;
 153     fRequireSigmaToVertex = c.fRequireSigmaToVertex ;
 154     fAODType = c.fAODType ;
 155     fAcceptKinkDaughters = c.fAcceptKinkDaughters ;
 156     fhCutStatistics = c.fhCutStatistics ;
 157     fhCutCorrelation = c.fhCutCorrelation ;
 158     fBitmap =  c.fBitmap;
 159     fhNBinsNSigma = c.fhNBinsNSigma;
 160     fhNBinsRequireSigma = c.fhNBinsRequireSigma;
 161     fhNBinsAcceptKink = c.fhNBinsAcceptKink;
 162     fhNBinsDcaXY = c.fhNBinsDcaXY;
 163     fhNBinsDcaZ = c.fhNBinsDcaZ;
 164     fhNBinsDcaXYnorm = c.fhNBinsDcaXYnorm;
 165     fhNBinsDcaZnorm = c.fhNBinsDcaZnorm;
 166     fhBinLimNSigma = c.fhBinLimNSigma;
 167     fhBinLimRequireSigma = c.fhBinLimRequireSigma;
 168     fhBinLimAcceptKink = c.fhBinLimAcceptKink;
 169     fhBinLimDcaXY = c.fhBinLimDcaXY;
 170     fhBinLimDcaZ = c.fhBinLimDcaZ;
 171     fhBinLimDcaXYnorm = c.fhBinLimDcaXYnorm;
 172     fhBinLimDcaZnorm = c.fhBinLimDcaZnorm;
 173
 174     for (Int_t j=0; j<c.kNStepQA; j++){
 175       if(c.fhDcaXYvsDcaZ[j]) fhDcaXYvsDcaZ[j] = (TH2F*)c.fhDcaXYvsDcaZ[j]->Clone();
 176       if(c.fhDcaXYvsDcaZnorm[j]) fhDcaXYvsDcaZnorm[j] = (TH2F*)c.fhDcaXYvsDcaZnorm[j]->Clone();
 177       for (Int_t i=0; i<c.kNHist; i++){
 178         if(c.fhQA[i][j]) fhQA[i][j] = (TH1F*)c.fhQA[i][j]->Clone();
 179       }
 180     }
 181     ((AliCFTrackIsPrimaryCuts &) c).Copy(*this);
 182  }
 183   return *this;
 184 }
 185 //__________________________________________________________________________________
 186 AliCFTrackIsPrimaryCuts::~AliCFTrackIsPrimaryCuts()
 187 {
 188   //
 189   // destructor
 190   //
 191   if (fhCutStatistics)                  delete fhCutStatistics;
 192   if (fhCutCorrelation)                 delete fhCutCorrelation;
 193
 194   for (Int_t j=0; j<kNStepQA; j++){
 195     if(fhDcaXYvsDcaZ[j])        delete fhDcaXYvsDcaZ[j];
 196     if(fhDcaXYvsDcaZnorm[j])    delete fhDcaXYvsDcaZnorm[j];
 197     for (Int_t i=0; i<kNHist; i++)
 198       if(fhQA[i][j])            delete fhQA[i][j];
 199   }
 200   if(fBitmap)                   delete fBitmap;
 201   if(fhBinLimNSigma)            delete fhBinLimNSigma;
 202   if(fhBinLimRequireSigma)      delete fhBinLimRequireSigma;
 203   if(fhBinLimAcceptKink)        delete fhBinLimAcceptKink;
 204   if(fhBinLimDcaXY)             delete fhBinLimDcaXY;
 205   if(fhBinLimDcaZ)              delete fhBinLimDcaZ;
 206   if(fhBinLimDcaXYnorm)         delete fhBinLimDcaXYnorm;
 207   if(fhBinLimDcaZnorm)          delete fhBinLimDcaZnorm;
 208 }
 209 //__________________________________________________________________________________
 210 void AliCFTrackIsPrimaryCuts::Initialise()
 211 {
 212   //
 213   // sets everything to zero
 214   //
 215   fNSigmaToVertex = 0;
 216   fNSigmaToVertexMax = 0;
 217   fRequireSigmaToVertex = 0;
 218   fAcceptKinkDaughters = 0;
 219   fAODType = AliAODTrack::kUndef;
 220
 221   SetMaxNSigmaToVertex();
 222   SetRequireSigmaToVertex();
 223   SetAcceptKinkDaughters();
 224   SetAODType();
 225
 226   for (Int_t j=0; j<kNStepQA; j++)  {
 227     fhDcaXYvsDcaZ[j] = 0x0;
 228     fhDcaXYvsDcaZnorm[j] = 0x0;
 229     for (Int_t i=0; i<kNHist; i++)
 230       fhQA[i][j] = 0x0;
 231   }
 232   fhCutStatistics = 0;
 233   fhCutCorrelation = 0;
 234   fBitmap=new TBits(0);
 235
 236   //set default bining for QA histograms
 237   SetHistogramBins(kCutNSigmaToVertex,500,0.,50.);
 238   SetHistogramBins(kCutRequireSigmaToVertex,5,-0.75,1.75);
 239   SetHistogramBins(kCutAcceptKinkDaughters,5,-0.75,1.75);
 240   SetHistogramBins(kDcaXY,500,-10.,10.);
 241   SetHistogramBins(kDcaZ,500,-10.,10.);
 242   SetHistogramBins(kDcaXYnorm,500,-10.,10.);
 243   SetHistogramBins(kDcaZnorm,500,-10.,10.);
 244 }
 245 //__________________________________________________________________________________
 246 void AliCFTrackIsPrimaryCuts::Copy(TObject &c) const
 247 {
 248   //
 249   // Copy function
 250   //
 251   AliCFTrackIsPrimaryCuts& target = (AliCFTrackIsPrimaryCuts &) c;
 252
 253   target.Initialise();
 254
 255   if (fhCutStatistics)  target.fhCutStatistics = (TH1F*) fhCutStatistics->Clone();
 256   if (fhCutCorrelation) target.fhCutCorrelation = (TH2F*) fhCutCorrelation->Clone();
 257
 258   for (Int_t j=0; j<kNStepQA; j++){
 259     if(fhDcaXYvsDcaZ[j]) target.fhDcaXYvsDcaZ[j] = (TH2F*)fhDcaXYvsDcaZ[j]->Clone();
 260     if(fhDcaXYvsDcaZnorm[j]) target.fhDcaXYvsDcaZnorm[j] = (TH2F*)fhDcaXYvsDcaZnorm[j]->Clone();
 261     for (Int_t i=0; i<kNHist; i++)
 262       if(fhQA[i][j]) target.fhQA[i][j] = (TH1F*)fhQA[i][j]->Clone();
 263   }
 264   TNamed::Copy(c);
 265 }
 266 //____________________________________________________________________
 267 void AliCFTrackIsPrimaryCuts::GetSigmaToVertex(AliESDtrack* esdTrack)
 268 {
 269   //
 270   // Calculates the number of sigma to the vertex.
 271   // Currently applicable for ESD tracks only
 272   //
 273   Float_t b[2];
 274   Float_t bRes[2];
 275   Float_t bCov[3];
 276   esdTrack->GetImpactParameters(b,bCov);
 277   if (bCov[0]<=0 || bCov[2]<=0) {
 278     AliDebug(1, "Estimated b resolution lower or equal zero!");
 279     bCov[0]=0; bCov[2]=0;
 280   }
 281   bRes[0] = TMath::Sqrt(bCov[0]);
 282   bRes[1] = TMath::Sqrt(bCov[2]);
 283
 284   // -----------------------------------
 285   // How to get to a n-sigma cut?
 286   //
 287   // The accumulated statistics from 0 to d is
 288   //
 289   // ->  Erf(d/Sqrt(2)) for a 1-dim gauss (d = n_sigma)
 290   // ->  1 - Exp(-d**2) for a 2-dim gauss (d*d = dx*dx + dy*dy != n_sigma)
 291   //
 292   // It means that for a 2-dim gauss: n_sigma(d) = Sqrt(2)*ErfInv(1 - Exp((-x**2)/2)
 293   // Can this be expressed in a different way?
 294
 295   if (bRes[0] == 0 || bRes[1] ==0) {
 296     fNSigmaToVertex = -1;
 297     return;
 298   }
 299
 300   Float_t d = TMath::Sqrt(TMath::Power(b[0]/bRes[0],2) + TMath::Power(b[1]/bRes[1],2));
 301
 302   // stupid rounding problem screws up everything:
 303   // if d is too big, TMath::Exp(...) gets 0, and TMath::ErfInverse(1) that should be infinite, gets 0 :(
 304   if (TMath::Exp(-d * d / 2) < 1e-10) {
 305     fNSigmaToVertex = 1000;
 306     return;
 307   }
 308
 309   d = TMath::ErfInverse(1 - TMath::Exp(-d * d / 2)) * TMath::Sqrt(2);
 310   fNSigmaToVertex = d;
 311   return;
 312 }
 313 //__________________________________________________________________________________
 314 void AliCFTrackIsPrimaryCuts::SelectionBitMap(TObject* obj)
 315 {
 316   //
 317   // test if the track passes the single cuts
 318   // and store the information in a bitmap
 319   //
 320
 321   // bitmap stores the decision of each single cut
 322   for(Int_t i=0; i<kNCuts; i++)fBitmap->SetBitNumber(i,kFALSE);
 323
 324   // check TObject and cast into ESDtrack
 325   if (!obj) return;
 326   if (!obj->InheritsFrom("AliVParticle")) {
 327     AliError("object must derived from AliVParticle !");
 328     return;
 329   }
 330
 331   Bool_t isESDTrack = strcmp(obj->ClassName(),"AliESDtrack") == 0 ? kTRUE : kFALSE ;
 332   Bool_t isAODTrack = strcmp(obj->ClassName(),"AliAODTrack") == 0 ? kTRUE : kFALSE ;
 333
 334   AliESDtrack * esdTrack = 0x0 ;
 335   AliAODTrack * aodTrack = 0x0 ;
 336   if (isESDTrack) esdTrack = dynamic_cast<AliESDtrack*>(obj);
 337   if (isAODTrack) aodTrack = dynamic_cast<AliAODTrack*>(obj);
 338
 339   // get the track to vertex parameter for ESD track
 340   if (isESDTrack) GetSigmaToVertex(esdTrack);
 341
 342   // fill the bitmap
 343   Int_t iCutBit = 0;
 344
 345   if (isESDTrack) {
 346     if (fNSigmaToVertex <= fNSigmaToVertexMax) {
 347       fBitmap->SetBitNumber(iCutBit,kTRUE);
 348     }
 349   }
 350   else fBitmap->SetBitNumber(iCutBit,kTRUE);
 351
 352   iCutBit++;
 353
 354   if (isESDTrack) {
 355     if (!fRequireSigmaToVertex || (fNSigmaToVertex>=0 && fRequireSigmaToVertex)) {
 356       fBitmap->SetBitNumber(iCutBit,kTRUE);
 357     }
 358   }
 359   else fBitmap->SetBitNumber(iCutBit,kTRUE);
 360
 361   iCutBit++;
 362
 363   if (esdTrack) {
 364     if (fAcceptKinkDaughters || (!fAcceptKinkDaughters && esdTrack->GetKinkIndex(0)<=0)) {
 365       fBitmap->SetBitNumber(iCutBit,kTRUE);
 366     }
 367   }
 368   else fBitmap->SetBitNumber(iCutBit,kTRUE);
 369
 370   iCutBit++;
 371
 372   if (isAODTrack) {
 373     if (fAODType==AliAODTrack::kUndef || fAODType == aodTrack->GetType()) {
 374       fBitmap->SetBitNumber(iCutBit,kTRUE);
 375     }
 376   }
 377   else fBitmap->SetBitNumber(iCutBit,kTRUE);
 378
 379   return;
 380 }
 381 //__________________________________________________________________________________
 382 Bool_t AliCFTrackIsPrimaryCuts::IsSelected(TObject* obj) {
 383   //
 384   // loops over decisions of single cuts and returns if the track is accepted
 385   //
 386   SelectionBitMap(obj);
 387
 388   if (fIsQAOn) FillHistograms(obj,0);
 389   Bool_t isSelected = kTRUE;
 390
 391   for (UInt_t icut=0; icut<fBitmap->GetNbits();icut++) {
 392     if(!fBitmap->TestBitNumber(icut)) {
 393         isSelected = kFALSE;
 394         break;
 395     }
 396   }
 397   if (!isSelected) return kFALSE ;
 398   if (fIsQAOn) FillHistograms(obj,1);
 399   return kTRUE;
 400 }
 401 //__________________________________________________________________________________
 402 void AliCFTrackIsPrimaryCuts::SetHistogramBins(Int_t index, Int_t nbins, Double_t *bins)
 403 {
 404   //
 405   // variable bin size
 406   //
 407   if(!fIsQAOn) return;
 408
 409   switch(index){
 410   case kCutNSigmaToVertex:
 411     fhNBinsNSigma=nbins+1;
 412     fhBinLimNSigma=new Double_t[nbins+1];
 413     for(Int_t i=0;i<nbins+1;i++)fhBinLimNSigma[i]=bins[i];
 414     break;
 415
 416   case kCutRequireSigmaToVertex:
 417     fhNBinsRequireSigma=nbins+1;
 418     fhBinLimRequireSigma=new Double_t[nbins+1];
 419     for(Int_t i=0;i<nbins+1;i++)fhBinLimRequireSigma[i]=bins[i];
 420     break;
 421
 422   case kCutAcceptKinkDaughters:
 423     fhNBinsAcceptKink=nbins+1;
 424     fhBinLimAcceptKink=new Double_t[nbins+1];
 425     for(Int_t i=0;i<nbins+1;i++)fhBinLimAcceptKink[i]=bins[i];
 426     break;
 427
 428   case kDcaXY:
 429     fhNBinsDcaXY=nbins+1;
 430     fhBinLimDcaXY=new Double_t[nbins+1];
 431     for(Int_t i=0;i<nbins+1;i++)fhBinLimDcaXY[i]=bins[i];
 432     break;
 433
 434   case kDcaZ:
 435     fhNBinsDcaZ=nbins+1;
 436     fhBinLimDcaZ=new Double_t[nbins+1];
 437     for(Int_t i=0;i<nbins+1;i++)fhBinLimDcaZ[i]=bins[i];
 438     break;
 439
 440   case kDcaXYnorm:
 441     fhNBinsDcaXYnorm=nbins+1;
 442     fhBinLimDcaXYnorm=new Double_t[nbins+1];
 443     for(Int_t i=0;i<nbins+1;i++)fhBinLimDcaXYnorm[i]=bins[i];
 444     break;
 445
 446   case kDcaZnorm:
 447     fhNBinsDcaZnorm=nbins+1;
 448     fhBinLimDcaZnorm=new Double_t[nbins+1];
 449     for(Int_t i=0;i<nbins+1;i++)fhBinLimDcaZnorm[i]=bins[i];
 450     break;
 451   }
 452 }
 453 //__________________________________________________________________________________
 454 void AliCFTrackIsPrimaryCuts::SetHistogramBins(Int_t index, Int_t nbins, Double_t xmin, Double_t xmax)
 455 {
 456   //
 457   // fixed bin size
 458   //
 459   switch(index){
 460   case kCutNSigmaToVertex:
 461     fhNBinsNSigma=nbins+1;
 462     fhBinLimNSigma=new Double_t[nbins+1];
 463     for(Int_t i=0;i<nbins+1;i++)fhBinLimNSigma[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
 464     break;
 465
 466   case kCutRequireSigmaToVertex:
 467     fhNBinsRequireSigma=nbins+1;
 468     fhBinLimRequireSigma=new Double_t[nbins+1];
 469     for(Int_t i=0;i<nbins+1;i++)fhBinLimRequireSigma[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
 470     break;
 471
 472   case kCutAcceptKinkDaughters:
 473     fhNBinsAcceptKink=nbins+1;
 474     fhBinLimAcceptKink=new Double_t[nbins+1];
 475     for(Int_t i=0;i<nbins+1;i++)fhBinLimAcceptKink[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
 476     break;
 477
 478   case kDcaXY:
 479     fhNBinsDcaXY=nbins+1;
 480     fhBinLimDcaXY=new Double_t[nbins+1];
 481     for(Int_t i=0;i<nbins+1;i++)fhBinLimDcaXY[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
 482     break;
 483
 484   case kDcaZ:
 485     fhNBinsDcaZ=nbins+1;
 486     fhBinLimDcaZ=new Double_t[nbins+1];
 487     for(Int_t i=0;i<nbins+1;i++)fhBinLimDcaZ[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
 488     break;
 489
 490   case kDcaXYnorm:
 491     fhNBinsDcaXYnorm=nbins+1;
 492     fhBinLimDcaXYnorm=new Double_t[nbins+1];
 493     for(Int_t i=0;i<nbins+1;i++)fhBinLimDcaXYnorm[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
 494     break;
 495
 496   case kDcaZnorm:
 497     fhNBinsDcaZnorm=nbins+1;
 498     fhBinLimDcaZnorm=new Double_t[nbins+1];
 499     for(Int_t i=0;i<nbins+1;i++)fhBinLimDcaZnorm[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
 500     break;
 501   }
 502 }
 503 //__________________________________________________________________________________
 504  void AliCFTrackIsPrimaryCuts::DefineHistograms() {
 505   //
 506   // histograms for cut variables, cut statistics and cut correlations
 507   //
 508   Int_t color = 2;
 509
 510   // book cut statistics and cut correlation histograms
 511   fhCutStatistics = new TH1F(Form("%s_cut_statistics",GetName()), Form("%s cut statistics",GetName()), kNCuts,0.5,kNCuts+0.5);
 512   fhCutStatistics->SetLineWidth(2);
 513   fhCutStatistics->GetXaxis()->SetBinLabel(1,"n dca");
 514   fhCutStatistics->GetXaxis()->SetBinLabel(2,"require dca");
 515   fhCutStatistics->GetXaxis()->SetBinLabel(3,"kink daughter");
 516   fhCutStatistics->GetXaxis()->SetBinLabel(4,"AOD type");
 517
 518   fhCutCorrelation = new TH2F(Form("%s_cut_correlation",GetName()), Form("%s cut  correlation",GetName()), kNCuts,0.5,kNCuts+0.5,kNCuts,0.5,kNCuts+0.5);
 519   fhCutCorrelation->SetLineWidth(2);
 520   fhCutCorrelation->GetXaxis()->SetBinLabel(1,"n dca");
 521   fhCutCorrelation->GetXaxis()->SetBinLabel(2,"require dca");
 522   fhCutCorrelation->GetXaxis()->SetBinLabel(3,"kink daughter");
 523   fhCutCorrelation->GetXaxis()->SetBinLabel(4,"AOD type");
 524
 525   fhCutCorrelation->GetYaxis()->SetBinLabel(1,"n dca");
 526   fhCutCorrelation->GetYaxis()->SetBinLabel(2,"require dca");
 527   fhCutCorrelation->GetYaxis()->SetBinLabel(3,"kink daughter");
 528   fhCutCorrelation->GetYaxis()->SetBinLabel(4,"AOD type");
 529
 530   // book QA histograms
 531   Char_t str[256];
 532   for (Int_t i=0; i<kNStepQA; i++) {
 533     if (i==0) sprintf(str," ");
 534     else sprintf(str,"_cut");
 535
 536     fhDcaXYvsDcaZ[i]            = new  TH2F(Form("%s_dcaXYvsDcaZ%s",GetName(),str),"",200,-10,10,200,-10,10);
 537     fhDcaXYvsDcaZnorm[i]        = new  TH2F(Form("%s_dcaXYvsDcaZnorm%s",GetName(),str),"",200,-10,10,200,-10,10);
 538     fhQA[kCutNSigmaToVertex][i] = new TH1F(Form("%s_nSigmaToVertex%s",GetName(),str),"",fhNBinsNSigma-1,fhBinLimNSigma);
 539     fhQA[kCutRequireSigmaToVertex][i] = new TH1F(Form("%s_requireSigmaToVertex%s",GetName(),str),"",fhNBinsRequireSigma-1,fhBinLimRequireSigma);
 540     fhQA[kCutAcceptKinkDaughters][i] = new TH1F(Form("%s_acceptKinkDaughters%s",GetName(),str),"",fhNBinsAcceptKink-1,fhBinLimAcceptKink);
 541     fhQA[kDcaXY][i]             = new TH1F(Form("%s_dcaXY%s",GetName(),str),"",fhNBinsDcaXY-1,fhBinLimDcaXY);
 542     fhQA[kDcaZ][i]              = new TH1F(Form("%s_dcaZ%s",GetName(),str),"",fhNBinsDcaZ-1,fhBinLimDcaZ);
 543     fhQA[kDcaXYnorm][i]         = new TH1F(Form("%s_dcaXYnorm%s",GetName(),str),"",fhNBinsDcaXYnorm-1,fhBinLimDcaXYnorm);
 544     fhQA[kDcaZnorm][i]          = new TH1F(Form("%s_dcaZnorm%s",GetName(),str),"",fhNBinsDcaZnorm-1,fhBinLimDcaZnorm);
 545
 546     fhDcaXYvsDcaZ[i]->SetXTitle("impact par. d_{z}");
 547     fhDcaXYvsDcaZ[i]->SetYTitle("impact par. d_{xy}");
 548     fhDcaXYvsDcaZnorm[i]->SetXTitle("norm. impact par. d_{z} / #sigma_{z}");
 549     fhDcaXYvsDcaZnorm[i]->SetYTitle("norm. impact par. d_{xy} / #sigma_{xy}");
 550
 551     fhQA[kCutNSigmaToVertex][i]->SetXTitle("n #sigma to vertex");
 552     fhQA[kCutRequireSigmaToVertex][i]->SetXTitle("require #sigma to vertex");
 553     fhQA[kCutAcceptKinkDaughters][i]->SetXTitle("accept kink daughters");
 554     fhQA[kDcaXY][i]->SetXTitle("impact par. d_{xy}");
 555     fhQA[kDcaZ][i]->SetXTitle("impact par. d_{z}");
 556     fhQA[kDcaXYnorm][i]->SetXTitle("norm. impact par. d_{xy} / #sigma_{xy}");
 557     fhQA[kDcaZnorm][i]->SetXTitle("norm. impact par. d_{z} / #sigma_{z}");
 558   }
 559   for(Int_t i=0; i<kNHist; i++) fhQA[i][1]->SetLineColor(color);
 560 }
 561 //__________________________________________________________________________________
 562 void AliCFTrackIsPrimaryCuts::FillHistograms(TObject* obj, Bool_t f)
 563 {
 564   //
 565   // fill the QA histograms
 566   //
 567
 568   if (!obj) return;
 569
 570   Bool_t isESDTrack = strcmp(obj->ClassName(),"AliESDtrack") == 0 ? kTRUE : kFALSE ;
 571   Bool_t isAODTrack = strcmp(obj->ClassName(),"AliAODTrack") == 0 ? kTRUE : kFALSE ;
 572
 573   AliESDtrack * esdTrack = 0x0 ;
 574   AliAODTrack * aodTrack = 0x0 ;
 575   if (isESDTrack) esdTrack = dynamic_cast<AliESDtrack*>(obj);
 576   if (isAODTrack) aodTrack = dynamic_cast<AliAODTrack*>(obj);
 577
 578   // f = 0: fill histograms before cuts
 579   // f = 1: fill histograms after cuts
 580
 581   if (esdTrack) {
 582     Float_t b[2];
 583     Float_t bRes[2];
 584     Float_t bCov[3];
 585     esdTrack->GetImpactParameters(b,bCov);
 586     if (bCov[0]<=0 || bCov[2]<=0) {
 587       AliDebug(1, "Estimated b resolution lower or equal zero!");
 588       bCov[0]=0; bCov[2]=0;
 589     }
 590     bRes[0] = TMath::Sqrt(bCov[0]);
 591     bRes[1] = TMath::Sqrt(bCov[2]);
 592
 593     fhQA[kDcaZ][f]->Fill(b[1]);
 594     fhQA[kDcaXY][f]->Fill(b[0]);
 595     fhDcaXYvsDcaZ[f]->Fill(b[1],b[0]);
 596
 597     if (bRes[0]!=0 && bRes[1]!=0) {
 598       fhQA[kDcaZnorm][f]->Fill(b[1]/bRes[1]);
 599       fhQA[kDcaXYnorm][f]->Fill(b[0]/bRes[0]);
 600       fhDcaXYvsDcaZnorm[f]->Fill(b[1]/bRes[1], b[0]/bRes[0]);
 601     }
 602
 603     fhQA[kCutNSigmaToVertex][f]->Fill(fNSigmaToVertex);
 604     if (fNSigmaToVertex<0 && fRequireSigmaToVertex) fhQA[kCutRequireSigmaToVertex][f]->Fill(0.);
 605     if (!(fNSigmaToVertex<0 && fRequireSigmaToVertex)) fhQA[kCutRequireSigmaToVertex][f]->Fill(1.);
 606
 607     if (!fAcceptKinkDaughters && esdTrack->GetKinkIndex(0)>0) fhQA[kCutAcceptKinkDaughters][f]->Fill(0.);
 608     if (!(!fAcceptKinkDaughters && esdTrack->GetKinkIndex(0)>0)) fhQA[kCutAcceptKinkDaughters][f]->Fill(0.);
 609   }
 610
 611   // fill cut statistics and cut correlation histograms with information from the bitmap
 612   if (f) return;
 613
 614   // Get the bitmap of the single cuts
 615   if ( !obj ) return;
 616   SelectionBitMap(obj);
 617
 618   // Number of single cuts in this class
 619   UInt_t ncuts = fBitmap->GetNbits();
 620   for(UInt_t bit=0; bit<ncuts;bit++) {
 621     if (!fBitmap->TestBitNumber(bit)) {
 622         fhCutStatistics->Fill(bit+1);
 623         for (UInt_t bit2=bit; bit2<ncuts;bit2++) {
 624           if (!fBitmap->TestBitNumber(bit2))
 625             fhCutCorrelation->Fill(bit+1,bit2+1);
 626         }
 627     }
 628   }
 629 }
 630 //__________________________________________________________________________________
 631 void AliCFTrackIsPrimaryCuts::SaveHistograms(const Char_t* dir) {
 632   //
 633   // saves the histograms in a directory (dir)
 634   //
 635   if(!fIsQAOn) return;
 636
 637   if (!dir)
 638     dir = GetName();
 639
 640   gDirectory->mkdir(dir);
 641   gDirectory->cd(dir);
 642
 643   gDirectory->mkdir("before_cuts");
 644   gDirectory->mkdir("after_cuts");
 645
 646   fhCutStatistics->Write();
 647   fhCutCorrelation->Write();
 648
 649   for (Int_t j=0; j<kNStepQA; j++) {
 650     if (j==0)
 651       gDirectory->cd("before_cuts");
 652     else
 653       gDirectory->cd("after_cuts");
 654
 655     fhDcaXYvsDcaZ[j]    ->Write();
 656     fhDcaXYvsDcaZnorm[j]->Write();
 657
 658     for(Int_t i=0; i<kNHist; i++) fhQA[i][j]->Write();
 659
 660     gDirectory->cd("../");
 661   }
 662   gDirectory->cd("../");
 663 }
 664 //__________________________________________________________________________________
 665 void AliCFTrackIsPrimaryCuts::DrawHistograms()
 666 {
 667   //
 668   // draws some histograms
 669   //
 670   if(!fIsQAOn) return;
 671
 672   // pad margins
 673   Float_t right = 0.03;
 674   Float_t left = 0.175;
 675   Float_t top = 0.03;
 676   Float_t bottom = 0.175;
 677
 678   TCanvas* canvas1 = new TCanvas("Track_QA_Primary_1", "Track QA Primary 1", 800, 500);
 679   canvas1->Divide(2, 1);
 680
 681   canvas1->cd(1);
 682   fhCutStatistics->SetStats(kFALSE);
 683   fhCutStatistics->LabelsOption("v");
 684   gPad->SetLeftMargin(left);
 685   gPad->SetBottomMargin(0.25);
 686   gPad->SetRightMargin(right);
 687   gPad->SetTopMargin(0.1);
 688   fhCutStatistics->Draw();
 689
 690   canvas1->cd(2);
 691   fhCutCorrelation->SetStats(kFALSE);
 692   fhCutCorrelation->LabelsOption("v");
 693   gPad->SetLeftMargin(0.30);
 694   gPad->SetRightMargin(bottom);
 695   gPad->SetTopMargin(0.1);
 696   gPad->SetBottomMargin(0.25);
 697   fhCutCorrelation->Draw("COLZ");
 698
 699   canvas1->SaveAs(Form("%s.eps", canvas1->GetName()));
 700   canvas1->SaveAs(Form("%s.ps", canvas1->GetName()));
 701
 702   // -----
 703
 704   TCanvas* canvas2 = new TCanvas("Track_QA_Primary_2", "Track QA Primary 2", 800, 800);
 705   canvas2->Divide(2, 2);
 706
 707   canvas2->cd(1);
 708   gPad->SetRightMargin(right);
 709   gPad->SetLeftMargin(left);
 710   gPad->SetTopMargin(top);
 711   gPad->SetBottomMargin(bottom);
 712   fhQA[kDcaXY][0]->SetStats(kFALSE);
 713   fhQA[kDcaXY][0]->Draw();
 714   fhQA[kDcaXY][1]->Draw("same");
 715
 716   canvas2->cd(2);
 717   gPad->SetRightMargin(right);
 718   gPad->SetLeftMargin(left);
 719   gPad->SetTopMargin(top);
 720   gPad->SetBottomMargin(bottom);
 721   fhQA[kDcaZ][0]->SetStats(kFALSE);
 722   fhQA[kDcaZ][0]->Draw();
 723   fhQA[kDcaZ][1]->Draw("same");
 724
 725   canvas2->cd(3);
 726 //   fhDXYvsDZ[0]->SetStats(kFALSE);
 727   gPad->SetLogz();
 728   gPad->SetLeftMargin(bottom);
 729   gPad->SetTopMargin(0.1);
 730   gPad->SetBottomMargin(bottom);
 731   gPad->SetRightMargin(0.2);
 732   fhDcaXYvsDcaZ[0]->Draw("COLZ");
 733
 734   canvas2->cd(4);
 735 //   fhDXYvsDZ[1]->SetStats(kFALSE);
 736   gPad->SetLogz();
 737   gPad->SetLeftMargin(bottom);
 738   gPad->SetTopMargin(0.1);
 739   gPad->SetBottomMargin(bottom);
 740   gPad->SetRightMargin(0.2);
 741   fhDcaXYvsDcaZ[1]->Draw("COLZ");
 742
 743   canvas2->SaveAs(Form("%s.eps", canvas2->GetName()));
 744   canvas2->SaveAs(Form("%s.ps", canvas2->GetName()));
 745
 746   // -----
 747
 748   TCanvas* canvas3 = new TCanvas("Track_QA_Primary_3", "Track QA Primary 3", 800, 800);
 749   canvas3->Divide(2, 2);
 750
 751   canvas3->cd(1);
 752   gPad->SetRightMargin(right);
 753   gPad->SetLeftMargin(left);
 754   gPad->SetTopMargin(top);
 755   gPad->SetBottomMargin(bottom);
 756   fhQA[kDcaXYnorm][0]->SetStats(kFALSE);
 757   fhQA[kDcaXYnorm][0]->Draw();
 758   fhQA[kDcaXYnorm][1]->Draw("same");
 759
 760   canvas3->cd(2);
 761   gPad->SetRightMargin(right);
 762   gPad->SetLeftMargin(left);
 763   gPad->SetTopMargin(top);
 764   gPad->SetBottomMargin(bottom);
 765   fhQA[kDcaZnorm][0]->SetStats(kFALSE);
 766   fhQA[kDcaZnorm][0]->Draw();
 767   fhQA[kDcaZnorm][1]->Draw("same");
 768
 769   canvas3->cd(3);
 770 //   fhDXYvsDZ[0]->SetStats(kFALSE);
 771   gPad->SetLogz();
 772   gPad->SetLeftMargin(bottom);
 773   gPad->SetTopMargin(0.1);
 774   gPad->SetBottomMargin(bottom);
 775   gPad->SetRightMargin(0.2);
 776   fhDcaXYvsDcaZnorm[0]->Draw("COLZ");
 777
 778   canvas3->cd(4);
 779 //   fhDXYvsDZ[1]->SetStats(kFALSE);
 780   gPad->SetLogz();
 781   gPad->SetLeftMargin(bottom);
 782   gPad->SetTopMargin(0.1);
 783   gPad->SetBottomMargin(bottom);
 784   gPad->SetRightMargin(0.2);
 785   fhDcaXYvsDcaZnorm[1]->Draw("COLZ");
 786
 787   canvas3->SaveAs(Form("%s.eps", canvas3->GetName()));
 788   canvas3->SaveAs(Form("%s.ps", canvas3->GetName()));
 789
 790   // -----
 791
 792   TCanvas* canvas4 = new TCanvas("Track_QA_Primary_4", "Track QA Primary 4", 1200, 500);
 793   canvas4->Divide(3, 1);
 794
 795   canvas4->cd(1);
 796   gPad->SetRightMargin(right);
 797   gPad->SetLeftMargin(left);
 798   gPad->SetTopMargin(top);
 799   gPad->SetBottomMargin(bottom);
 800   fhQA[kCutNSigmaToVertex][0]->SetStats(kFALSE);
 801   fhQA[kCutNSigmaToVertex][0]->Draw();
 802   fhQA[kCutNSigmaToVertex][1]->Draw("same");
 803
 804   canvas4->cd(2);
 805   gPad->SetRightMargin(right);
 806   gPad->SetLeftMargin(left);
 807   gPad->SetTopMargin(top);
 808   gPad->SetBottomMargin(bottom);
 809   fhQA[kCutRequireSigmaToVertex][0]->SetStats(kFALSE);
 810   fhQA[kCutRequireSigmaToVertex][0]->Draw();
 811   fhQA[kCutRequireSigmaToVertex][1]->Draw("same");
 812
 813   canvas4->cd(3);
 814   gPad->SetRightMargin(right);
 815   gPad->SetLeftMargin(left);
 816   gPad->SetTopMargin(top);
 817   gPad->SetBottomMargin(bottom);
 818   fhQA[kCutAcceptKinkDaughters][0]->SetStats(kFALSE);
 819   fhQA[kCutAcceptKinkDaughters][0]->Draw();
 820   fhQA[kCutAcceptKinkDaughters][1]->Draw("same");
 821
 822   canvas4->SaveAs(Form("%s.eps", canvas4->GetName()));
 823   canvas4->SaveAs(Form("%s.ps", canvas4->GetName()));
 824 }
 825 //__________________________________________________________________________________
 826 void AliCFTrackIsPrimaryCuts::AddQAHistograms(TList *qaList) {
 827   //
 828   // saves the histograms in a TList
 829   //
 830   DefineHistograms();
 831
 832   qaList->Add(fhCutStatistics);
 833   qaList->Add(fhCutCorrelation);
 834
 835   for (Int_t j=0; j<kNStepQA; j++) {
 836     qaList->Add(fhDcaXYvsDcaZ[j]);
 837     qaList->Add(fhDcaXYvsDcaZnorm[j]);
 838     for(Int_t i=0; i<kNHist; i++)
 839         qaList->Add(fhQA[i][j]);
 840   }
 841 }