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