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