Adding pt and eta control plots, increasing ClassDef version
[u/mrichter/AliRoot.git] / PWG0 / esdTrackCuts / AliESDtrackCuts.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 /* $Id$ */
17
18 #include "AliESDtrackCuts.h"
19
20 #include <AliESDtrack.h>
21 #include <AliESD.h>
22 #include <AliESDEvent.h>
23 #include <AliLog.h>
24
25 #include <TTree.h>
26 #include <TCanvas.h>
27 #include <TDirectory.h>
28
29 //____________________________________________________________________
30 ClassImp(AliESDtrackCuts)
31
32 // Cut names
33 const Char_t* AliESDtrackCuts::fgkCutNames[kNCuts] = {
34  "require TPC refit",
35  "require ITS refit",
36  "n clusters TPC",
37  "n clusters ITS",
38  "#Chi^{2}/clusters TPC",
39  "#Chi^{2}/clusters ITS",
40  "cov 11",
41  "cov 22",
42  "cov 33",
43  "cov 44",
44  "cov 55",
45  "trk-to-vtx",
46  "trk-to-vtx failed",
47  "kink daughters",
48  "p",
49  "p_{T}",
50  "p_{x}",
51  "p_{y}",
52  "p_{z}",
53  "y",
54  "eta"
55 };
56
57 //____________________________________________________________________
58 AliESDtrackCuts::AliESDtrackCuts(const Char_t* name, const Char_t* title) : AliAnalysisCuts(name,title),
59   fCutMinNClusterTPC(0),
60   fCutMinNClusterITS(0),
61   fCutMaxChi2PerClusterTPC(0),
62   fCutMaxChi2PerClusterITS(0),
63   fCutMaxC11(0),
64   fCutMaxC22(0),
65   fCutMaxC33(0),
66   fCutMaxC44(0),
67   fCutMaxC55(0),
68   fCutAcceptKinkDaughters(0),
69   fCutRequireTPCRefit(0),
70   fCutRequireITSRefit(0),
71   fCutNsigmaToVertex(0),
72   fCutSigmaToVertexRequired(0),
73   fPMin(0),
74   fPMax(0),
75   fPtMin(0),
76   fPtMax(0),
77   fPxMin(0),
78   fPxMax(0),
79   fPyMin(0),
80   fPyMax(0),
81   fPzMin(0),
82   fPzMax(0),
83   fEtaMin(0),
84   fEtaMax(0),
85   fRapMin(0),
86   fRapMax(0),
87   fHistogramsOn(0),
88   fhCutStatistics(0),         
89   fhCutCorrelation(0)
90 {
91   //
92   // constructor
93   //
94
95   Init();
96
97   //##############################################################################
98   // setting default cuts
99   SetMinNClustersTPC();
100   SetMinNClustersITS();
101   SetMaxChi2PerClusterTPC();
102   SetMaxChi2PerClusterITS();                                
103   SetMaxCovDiagonalElements();                                      
104   SetRequireTPCRefit();
105   SetRequireITSRefit();
106   SetAcceptKingDaughters();
107   SetMinNsigmaToVertex();
108   SetRequireSigmaToVertex();
109   SetPRange();
110   SetPtRange();
111   SetPxRange();
112   SetPyRange();
113   SetPzRange();
114   SetEtaRange();
115   SetRapRange();
116
117   SetHistogramsOn();
118 }
119
120 //_____________________________________________________________________________
121 AliESDtrackCuts::AliESDtrackCuts(const AliESDtrackCuts &c) : AliAnalysisCuts(c),
122   fCutMinNClusterTPC(0),
123   fCutMinNClusterITS(0),
124   fCutMaxChi2PerClusterTPC(0),
125   fCutMaxChi2PerClusterITS(0),
126   fCutMaxC11(0),
127   fCutMaxC22(0),
128   fCutMaxC33(0),
129   fCutMaxC44(0),
130   fCutMaxC55(0),
131   fCutAcceptKinkDaughters(0),
132   fCutRequireTPCRefit(0),
133   fCutRequireITSRefit(0),
134   fCutNsigmaToVertex(0),
135   fCutSigmaToVertexRequired(0),
136   fPMin(0),
137   fPMax(0),
138   fPtMin(0),
139   fPtMax(0),
140   fPxMin(0),
141   fPxMax(0),
142   fPyMin(0),
143   fPyMax(0),
144   fPzMin(0),
145   fPzMax(0),
146   fEtaMin(0),
147   fEtaMax(0),
148   fRapMin(0),
149   fRapMax(0),
150   fHistogramsOn(0),
151   ffDTheoretical(0),                                 
152   fhCutStatistics(0),         
153   fhCutCorrelation(0)
154 {
155   //
156   // copy constructor
157   //
158
159   ((AliESDtrackCuts &) c).Copy(*this);
160 }
161
162 AliESDtrackCuts::~AliESDtrackCuts()
163 {
164   //
165   // destructor
166   //
167
168   for (Int_t i=0; i<2; i++) {
169     
170     if (fhNClustersITS[i])
171       delete fhNClustersITS[i];            
172     if (fhNClustersTPC[i])
173       delete fhNClustersTPC[i];                
174     if (fhChi2PerClusterITS[i])
175       delete fhChi2PerClusterITS[i];       
176     if (fhChi2PerClusterTPC[i])
177       delete fhChi2PerClusterTPC[i];       
178     if (fhC11[i])
179       delete fhC11[i];                     
180     if (fhC22[i])
181       delete fhC22[i];                     
182     if (fhC33[i])
183       delete fhC33[i];                     
184     if (fhC44[i])
185       delete fhC44[i];                     
186     if (fhC55[i])
187     delete fhC55[i];                     
188     
189     if (fhDXY[i])
190       delete fhDXY[i];                     
191     if (fhDZ[i])
192       delete fhDZ[i];                      
193     if (fhDXYvsDZ[i])
194       delete fhDXYvsDZ[i];                 
195     
196     if (fhDXYNormalized[i])
197       delete fhDXYNormalized[i];           
198     if (fhDZNormalized[i])
199       delete fhDZNormalized[i];
200     if (fhDXYvsDZNormalized[i])
201       delete fhDXYvsDZNormalized[i];       
202     if (fhNSigmaToVertex[i])
203       delete fhNSigmaToVertex[i];
204     if (fhPt[i])
205       delete fhPt[i];
206     if (fhEta[i])
207       delete fhEta[i];
208   }
209
210   if (ffDTheoretical)
211     delete ffDTheoretical;
212
213   if (fhCutStatistics)
214     delete fhCutStatistics;             
215   if (fhCutCorrelation)
216     delete fhCutCorrelation;            
217 }
218
219 void AliESDtrackCuts::Init()
220 {
221   //
222   // sets everything to zero
223   //
224
225   fCutMinNClusterTPC = 0;
226   fCutMinNClusterITS = 0;
227
228   fCutMaxChi2PerClusterTPC = 0;
229   fCutMaxChi2PerClusterITS = 0;
230
231   fCutMaxC11 = 0;
232   fCutMaxC22 = 0;
233   fCutMaxC33 = 0;
234   fCutMaxC44 = 0;
235   fCutMaxC55 = 0;
236
237   fCutAcceptKinkDaughters = 0;
238   fCutRequireTPCRefit = 0;
239   fCutRequireITSRefit = 0;
240
241   fCutNsigmaToVertex = 0;
242   fCutSigmaToVertexRequired = 0;
243
244   fPMin = 0;
245   fPMax = 0;
246   fPtMin = 0;
247   fPtMax = 0;
248   fPxMin = 0;
249   fPxMax = 0;
250   fPyMin = 0;
251   fPyMax = 0;
252   fPzMin = 0;
253   fPzMax = 0;
254   fEtaMin = 0;
255   fEtaMax = 0;
256   fRapMin = 0;
257   fRapMax = 0;
258
259   fHistogramsOn = kFALSE;
260
261   for (Int_t i=0; i<2; ++i)
262   {
263     fhNClustersITS[i] = 0;
264     fhNClustersTPC[i] = 0;
265
266     fhChi2PerClusterITS[i] = 0;
267     fhChi2PerClusterTPC[i] = 0;
268
269     fhC11[i] = 0;
270     fhC22[i] = 0;
271     fhC33[i] = 0;
272     fhC44[i] = 0;
273     fhC55[i] = 0;
274
275     fhDXY[i] = 0;
276     fhDZ[i] = 0;
277     fhDXYvsDZ[i] = 0;
278
279     fhDXYNormalized[i] = 0;
280     fhDZNormalized[i] = 0;
281     fhDXYvsDZNormalized[i] = 0;
282     fhNSigmaToVertex[i] = 0;
283     
284     fhPt[i] = 0;
285     fhEta[i] = 0;
286   }
287   ffDTheoretical = 0;
288
289   fhCutStatistics = 0;
290   fhCutCorrelation = 0;
291 }
292
293 //_____________________________________________________________________________
294 AliESDtrackCuts &AliESDtrackCuts::operator=(const AliESDtrackCuts &c)
295 {
296   //
297   // Assignment operator
298   //
299
300   if (this != &c) ((AliESDtrackCuts &) c).Copy(*this);
301   return *this;
302 }
303
304 //_____________________________________________________________________________
305 void AliESDtrackCuts::Copy(TObject &c) const
306 {
307   //
308   // Copy function
309   //
310
311   AliESDtrackCuts& target = (AliESDtrackCuts &) c;
312
313   target.Init();
314
315   target.fCutMinNClusterTPC = fCutMinNClusterTPC;
316   target.fCutMinNClusterITS = fCutMinNClusterITS;
317
318   target.fCutMaxChi2PerClusterTPC = fCutMaxChi2PerClusterTPC;
319   target.fCutMaxChi2PerClusterITS = fCutMaxChi2PerClusterITS;
320
321   target.fCutMaxC11 = fCutMaxC11;
322   target.fCutMaxC22 = fCutMaxC22;
323   target.fCutMaxC33 = fCutMaxC33;
324   target.fCutMaxC44 = fCutMaxC44;
325   target.fCutMaxC55 = fCutMaxC55;
326
327   target.fCutAcceptKinkDaughters = fCutAcceptKinkDaughters;
328   target.fCutRequireTPCRefit = fCutRequireTPCRefit;
329   target.fCutRequireITSRefit = fCutRequireITSRefit;
330
331   target.fCutNsigmaToVertex = fCutNsigmaToVertex;
332   target.fCutSigmaToVertexRequired = fCutSigmaToVertexRequired;
333
334   target.fPMin = fPMin;
335   target.fPMax = fPMax;
336   target.fPtMin = fPtMin;
337   target.fPtMax = fPtMax;
338   target.fPxMin = fPxMin;
339   target.fPxMax = fPxMax;
340   target.fPyMin = fPyMin;
341   target.fPyMax = fPyMax;
342   target.fPzMin = fPzMin;
343   target.fPzMax = fPzMax;
344   target.fEtaMin = fEtaMin;
345   target.fEtaMax = fEtaMax;
346   target.fRapMin = fRapMin;
347   target.fRapMax = fRapMax;
348
349   target.fHistogramsOn = fHistogramsOn;
350
351   for (Int_t i=0; i<2; ++i)
352   {
353     if (fhNClustersITS[i]) target.fhNClustersITS[i] = (TH1F*) fhNClustersITS[i]->Clone();
354     if (fhNClustersTPC[i]) target.fhNClustersTPC[i] = (TH1F*) fhNClustersTPC[i]->Clone();
355
356     if (fhChi2PerClusterITS[i]) target.fhChi2PerClusterITS[i] = (TH1F*) fhChi2PerClusterITS[i]->Clone();
357     if (fhChi2PerClusterTPC[i]) target.fhChi2PerClusterTPC[i] = (TH1F*) fhChi2PerClusterTPC[i]->Clone();
358
359     if (fhC11[i]) target.fhC11[i] = (TH1F*) fhC11[i]->Clone();
360     if (fhC22[i]) target.fhC22[i] = (TH1F*) fhC22[i]->Clone();
361     if (fhC33[i]) target.fhC33[i] = (TH1F*) fhC33[i]->Clone();
362     if (fhC44[i]) target.fhC44[i] = (TH1F*) fhC44[i]->Clone();
363     if (fhC55[i]) target.fhC55[i] = (TH1F*) fhC55[i]->Clone();
364
365     if (fhDXY[i]) target.fhDXY[i] = (TH1F*) fhDXY[i]->Clone();
366     if (fhDZ[i]) target.fhDZ[i] = (TH1F*) fhDZ[i]->Clone();
367     if (fhDXYvsDZ[i]) target.fhDXYvsDZ[i] = (TH2F*) fhDXYvsDZ[i]->Clone();
368
369     if (fhDXYNormalized[i]) target.fhDXYNormalized[i] = (TH1F*) fhDXYNormalized[i]->Clone();
370     if (fhDZNormalized[i]) target.fhDZNormalized[i] = (TH1F*) fhDZNormalized[i]->Clone();
371     if (fhDXYvsDZNormalized[i]) target.fhDXYvsDZNormalized[i] = (TH2F*) fhDXYvsDZNormalized[i]->Clone();
372     if (fhNSigmaToVertex[i]) target.fhNSigmaToVertex[i] = (TH1F*) fhNSigmaToVertex[i]->Clone();
373     
374     if (fhPt[i]) target.fhPt[i] = (TH1F*) fhPt[i]->Clone();
375     if (fhEta[i]) target.fhEta[i] = (TH1F*) fhEta[i]->Clone();
376   }
377   if (ffDTheoretical) target.ffDTheoretical = (TF1*) ffDTheoretical->Clone();
378
379   if (fhCutStatistics) target.fhCutStatistics = (TH1F*) fhCutStatistics->Clone();
380   if (fhCutCorrelation) target.fhCutCorrelation = (TH2F*) fhCutCorrelation->Clone();
381
382   TNamed::Copy(c);
383 }
384
385 //_____________________________________________________________________________
386 Long64_t AliESDtrackCuts::Merge(TCollection* list) {
387   // Merge a list of AliESDtrackCuts objects with this (needed for PROOF)
388   // Returns the number of merged objects (including this)
389
390   if (!list)
391     return 0;
392   
393   if (list->IsEmpty())
394     return 1;
395
396   if (!fHistogramsOn)
397     return 0;
398
399   TIterator* iter = list->MakeIterator();
400   TObject* obj;
401
402
403   // collection of measured and generated histograms
404   Int_t count = 0;
405   while ((obj = iter->Next())) {
406
407     AliESDtrackCuts* entry = dynamic_cast<AliESDtrackCuts*>(obj);
408     if (entry == 0)
409       continue;
410
411     if (!entry->fHistogramsOn)
412       continue;
413     
414     for (Int_t i=0; i<2; i++) {
415       
416       fhNClustersITS[i]      ->Add(entry->fhNClustersITS[i]     );      
417       fhNClustersTPC[i]      ->Add(entry->fhNClustersTPC[i]     ); 
418                                                                     
419       fhChi2PerClusterITS[i] ->Add(entry->fhChi2PerClusterITS[i]); 
420       fhChi2PerClusterTPC[i] ->Add(entry->fhChi2PerClusterTPC[i]); 
421                                                                     
422       fhC11[i]               ->Add(entry->fhC11[i]              ); 
423       fhC22[i]               ->Add(entry->fhC22[i]              ); 
424       fhC33[i]               ->Add(entry->fhC33[i]              ); 
425       fhC44[i]               ->Add(entry->fhC44[i]              ); 
426       fhC55[i]               ->Add(entry->fhC55[i]              ); 
427                                                                     
428       fhDXY[i]               ->Add(entry->fhDXY[i]              ); 
429       fhDZ[i]                ->Add(entry->fhDZ[i]               ); 
430       fhDXYvsDZ[i]           ->Add(entry->fhDXYvsDZ[i]          ); 
431                                                                     
432       fhDXYNormalized[i]     ->Add(entry->fhDXYNormalized[i]    ); 
433       fhDZNormalized[i]      ->Add(entry->fhDZNormalized[i]     );
434       fhDXYvsDZNormalized[i] ->Add(entry->fhDXYvsDZNormalized[i]); 
435       fhNSigmaToVertex[i]    ->Add(entry->fhNSigmaToVertex[i]); 
436
437       fhPt[i]                ->Add(entry->fhPt[i]); 
438       fhEta[i]               ->Add(entry->fhEta[i]); 
439     }      
440
441     fhCutStatistics  ->Add(entry->fhCutStatistics);        
442     fhCutCorrelation ->Add(entry->fhCutCorrelation);      
443
444     count++;
445   }
446
447   return count+1;
448 }
449
450
451 //____________________________________________________________________
452 Float_t AliESDtrackCuts::GetSigmaToVertex(AliESDtrack* esdTrack)
453 {
454   // Calculates the number of sigma to the vertex.
455
456   Float_t b[2];
457   Float_t bRes[2];
458   Float_t bCov[3];
459   esdTrack->GetImpactParameters(b,bCov);
460   if (bCov[0]<=0 || bCov[2]<=0) {
461     AliDebug(1, "Estimated b resolution lower or equal zero!");
462     bCov[0]=0; bCov[2]=0;
463   }
464   bRes[0] = TMath::Sqrt(bCov[0]);
465   bRes[1] = TMath::Sqrt(bCov[2]);
466
467   // -----------------------------------
468   // How to get to a n-sigma cut?
469   //
470   // The accumulated statistics from 0 to d is
471   //
472   // ->  Erf(d/Sqrt(2)) for a 1-dim gauss (d = n_sigma)
473   // ->  1 - Exp(-d**2) for a 2-dim gauss (d*d = dx*dx + dy*dy != n_sigma)
474   //
475   // It means that for a 2-dim gauss: n_sigma(d) = Sqrt(2)*ErfInv(1 - Exp((-x**2)/2)
476   // Can this be expressed in a different way?
477
478   if (bRes[0] == 0 || bRes[1] ==0)
479     return -1;
480
481   Float_t d = TMath::Sqrt(TMath::Power(b[0]/bRes[0],2) + TMath::Power(b[1]/bRes[1],2));
482
483   // stupid rounding problem screws up everything:
484   // if d is too big, TMath::Exp(...) gets 0, and TMath::ErfInverse(1) that should be infinite, gets 0 :(
485   if (TMath::Exp(-d * d / 2) < 1e-10)
486     return 1000;
487
488   d = TMath::ErfInverse(1 - TMath::Exp(-d * d / 2)) * TMath::Sqrt(2);
489   return d;
490 }
491
492 void AliESDtrackCuts::EnableNeededBranches(TTree* tree)
493 {
494   // enables the branches needed by AcceptTrack, for a list see comment of AcceptTrack
495
496   tree->SetBranchStatus("fTracks.fFlags", 1);
497   tree->SetBranchStatus("fTracks.fITSncls", 1);
498   tree->SetBranchStatus("fTracks.fTPCncls", 1);
499   tree->SetBranchStatus("fTracks.fITSchi2", 1);
500   tree->SetBranchStatus("fTracks.fTPCchi2", 1);
501   tree->SetBranchStatus("fTracks.fC*", 1);
502   tree->SetBranchStatus("fTracks.fD", 1);
503   tree->SetBranchStatus("fTracks.fZ", 1);
504   tree->SetBranchStatus("fTracks.fCdd", 1);
505   tree->SetBranchStatus("fTracks.fCdz", 1);
506   tree->SetBranchStatus("fTracks.fCzz", 1);
507   tree->SetBranchStatus("fTracks.fP*", 1);
508   tree->SetBranchStatus("fTracks.fR*", 1);
509   tree->SetBranchStatus("fTracks.fKinkIndexes*", 1);
510 }
511
512 //____________________________________________________________________
513 Bool_t
514 AliESDtrackCuts::AcceptTrack(AliESDtrack* esdTrack) {
515   // 
516   // figure out if the tracks survives all the track cuts defined
517   //
518   // the different quality parameter and kinematic values are first
519   // retrieved from the track. then it is found out what cuts the
520   // track did not survive and finally the cuts are imposed.
521
522   // this function needs the following branches:
523   // fTracks.fFlags
524   // fTracks.fITSncls
525   // fTracks.fTPCncls
526   // fTracks.fITSchi2
527   // fTracks.fTPCchi2
528   // fTracks.fC   //GetExternalCovariance
529   // fTracks.fD   //GetImpactParameters
530   // fTracks.fZ   //GetImpactParameters
531   // fTracks.fCdd //GetImpactParameters
532   // fTracks.fCdz //GetImpactParameters
533   // fTracks.fCzz //GetImpactParameters
534   // fTracks.fP   //GetPxPyPz
535   // fTracks.fR   //GetMass
536   // fTracks.fP   //GetMass
537   // fTracks.fKinkIndexes
538
539   UInt_t status = esdTrack->GetStatus();
540
541   // dummy array
542   Int_t  fIdxInt[200];
543
544   // getting quality parameters from the ESD track
545   Int_t nClustersITS = esdTrack->GetITSclusters(fIdxInt);
546   Int_t nClustersTPC = esdTrack->GetTPCclusters(fIdxInt);
547   
548
549
550   Float_t chi2PerClusterITS = -1;
551   Float_t chi2PerClusterTPC = -1;
552   if (nClustersITS!=0)
553     chi2PerClusterITS = esdTrack->GetITSchi2()/Float_t(nClustersITS);
554   if (nClustersTPC!=0)
555     chi2PerClusterTPC = esdTrack->GetTPCchi2()/Float_t(nClustersTPC);
556
557   Double_t extCov[15];
558   esdTrack->GetExternalCovariance(extCov);
559
560   // getting the track to vertex parameters
561   Float_t nSigmaToVertex = GetSigmaToVertex(esdTrack);
562
563   // getting the kinematic variables of the track
564   // (assuming the mass is known)
565   Double_t p[3];
566   esdTrack->GetPxPyPz(p);
567   Float_t momentum = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2) + TMath::Power(p[2],2));
568   Float_t pt       = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2));
569   Float_t energy   = TMath::Sqrt(TMath::Power(esdTrack->GetMass(),2) + TMath::Power(momentum,2));
570
571
572   //y-eta related calculations
573   Float_t eta = -100.;
574   Float_t y   = -100.;
575   if((momentum != TMath::Abs(p[2]))&&(momentum != 0))
576     eta = 0.5*TMath::Log((momentum + p[2])/(momentum - p[2]));
577   if((energy != TMath::Abs(p[2]))&&(momentum != 0))
578     y = 0.5*TMath::Log((energy + p[2])/(energy - p[2]));
579
580   
581   //########################################################################
582   // cut the track?
583   
584   Bool_t cuts[kNCuts];
585   for (Int_t i=0; i<kNCuts; i++) cuts[i]=kFALSE;
586   
587   // track quality cuts
588   if (fCutRequireTPCRefit && (status&AliESDtrack::kTPCrefit)==0)
589     cuts[0]=kTRUE;
590   if (fCutRequireITSRefit && (status&AliESDtrack::kITSrefit)==0)
591     cuts[1]=kTRUE;
592   if (nClustersTPC<fCutMinNClusterTPC)
593     cuts[2]=kTRUE;
594   if (nClustersITS<fCutMinNClusterITS) 
595     cuts[3]=kTRUE;
596   if (chi2PerClusterTPC>fCutMaxChi2PerClusterTPC) 
597     cuts[4]=kTRUE; 
598   if (chi2PerClusterITS>fCutMaxChi2PerClusterITS) 
599     cuts[5]=kTRUE;
600   if (extCov[0]  > fCutMaxC11) 
601     cuts[6]=kTRUE;  
602   if (extCov[2]  > fCutMaxC22) 
603     cuts[7]=kTRUE;  
604   if (extCov[5]  > fCutMaxC33) 
605     cuts[8]=kTRUE;  
606   if (extCov[9]  > fCutMaxC44) 
607     cuts[9]=kTRUE;  
608   if (extCov[14]  > fCutMaxC55) 
609     cuts[10]=kTRUE;  
610   if (nSigmaToVertex > fCutNsigmaToVertex && fCutSigmaToVertexRequired)
611     cuts[11] = kTRUE;
612   // if n sigma could not be calculated
613   if (nSigmaToVertex<0 && fCutSigmaToVertexRequired)
614     cuts[12]=kTRUE;
615   if (!fCutAcceptKinkDaughters && esdTrack->GetKinkIndex(0)>0)
616     cuts[13]=kTRUE;
617   // track kinematics cut
618   if((momentum < fPMin) || (momentum > fPMax)) 
619     cuts[14]=kTRUE;
620   if((pt < fPtMin) || (pt > fPtMax)) 
621     cuts[15] = kTRUE;
622   if((p[0] < fPxMin) || (p[0] > fPxMax)) 
623     cuts[16] = kTRUE;
624   if((p[1] < fPyMin) || (p[1] > fPyMax)) 
625     cuts[17] = kTRUE;
626   if((p[2] < fPzMin) || (p[2] > fPzMax))
627     cuts[18] = kTRUE;
628   if((eta < fEtaMin) || (eta > fEtaMax)) 
629     cuts[19] = kTRUE;
630   if((y < fRapMin) || (y > fRapMax)) 
631     cuts[20] = kTRUE;
632
633   Bool_t cut=kFALSE;
634   for (Int_t i=0; i<kNCuts; i++) 
635     if (cuts[i]) cut = kTRUE;
636   
637   //########################################################################
638   // filling histograms
639   if (fHistogramsOn) {
640     fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n tracks")));
641     
642     if (cut)
643       fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n cut tracks")));
644     
645     for (Int_t i=0; i<kNCuts; i++) {
646       if (cuts[i])
647         fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin(fgkCutNames[i])));
648       
649       for (Int_t j=i; j<kNCuts; j++) {
650         if (cuts[i] && cuts[j]) {
651           Float_t x = fhCutCorrelation->GetXaxis()->GetBinCenter(fhCutCorrelation->GetXaxis()->FindBin(fgkCutNames[i]));
652           Float_t y = fhCutCorrelation->GetYaxis()->GetBinCenter(fhCutCorrelation->GetYaxis()->FindBin(fgkCutNames[j]));
653           fhCutCorrelation->Fill(x,y);
654         }
655       }
656     }
657     
658     fhNClustersITS[0]->Fill(nClustersITS);
659     fhNClustersTPC[0]->Fill(nClustersTPC);
660     fhChi2PerClusterITS[0]->Fill(chi2PerClusterITS);
661     fhChi2PerClusterTPC[0]->Fill(chi2PerClusterTPC);
662
663     fhC11[0]->Fill(extCov[0]);
664     fhC22[0]->Fill(extCov[2]);
665     fhC33[0]->Fill(extCov[5]);
666     fhC44[0]->Fill(extCov[9]);
667     fhC55[0]->Fill(extCov[14]);
668     
669     fhPt[0]->Fill(pt);
670     fhEta[0]->Fill(eta);
671
672     Float_t b[2];
673     Float_t bRes[2];
674     Float_t bCov[3];
675     esdTrack->GetImpactParameters(b,bCov);
676     if (bCov[0]<=0 || bCov[2]<=0) {
677       AliDebug(1, "Estimated b resolution lower or equal zero!");
678       bCov[0]=0; bCov[2]=0;
679     }
680     bRes[0] = TMath::Sqrt(bCov[0]);
681     bRes[1] = TMath::Sqrt(bCov[2]);
682
683     fhDZ[0]->Fill(b[1]);
684     fhDXY[0]->Fill(b[0]);
685     fhDXYvsDZ[0]->Fill(b[1],b[0]);
686
687     if (bRes[0]!=0 && bRes[1]!=0) {
688       fhDZNormalized[0]->Fill(b[1]/bRes[1]);
689       fhDXYNormalized[0]->Fill(b[0]/bRes[0]);
690       fhDXYvsDZNormalized[0]->Fill(b[1]/bRes[1], b[0]/bRes[0]);
691       fhNSigmaToVertex[0]->Fill(nSigmaToVertex);
692     }
693   }
694
695   //########################################################################
696   // cut the track!
697   if (cut) return kFALSE;
698
699   //########################################################################
700   // filling histograms after cut
701   if (fHistogramsOn) {
702     fhNClustersITS[1]->Fill(nClustersITS);
703     fhNClustersTPC[1]->Fill(nClustersTPC);
704     fhChi2PerClusterITS[1]->Fill(chi2PerClusterITS);
705     fhChi2PerClusterTPC[1]->Fill(chi2PerClusterTPC);
706
707     fhC11[1]->Fill(extCov[0]);
708     fhC22[1]->Fill(extCov[2]);
709     fhC33[1]->Fill(extCov[5]);
710     fhC44[1]->Fill(extCov[9]);
711     fhC55[1]->Fill(extCov[14]);
712
713     fhPt[1]->Fill(pt);
714     fhEta[1]->Fill(eta);
715     
716     Float_t b[2];
717     Float_t bRes[2];
718     Float_t bCov[3];
719     esdTrack->GetImpactParameters(b,bCov);
720     if (bCov[0]<=0 || bCov[2]<=0) {
721       AliDebug(1, "Estimated b resolution lower or equal zero!");
722       bCov[0]=0; bCov[2]=0;
723     }
724     bRes[0] = TMath::Sqrt(bCov[0]);
725     bRes[1] = TMath::Sqrt(bCov[2]);
726
727     fhDZ[1]->Fill(b[1]);
728     fhDXY[1]->Fill(b[0]);
729     fhDXYvsDZ[1]->Fill(b[1],b[0]);
730
731     if (bRes[0]!=0 && bRes[1]!=0)
732     {
733       fhDZNormalized[1]->Fill(b[1]/bRes[1]);
734       fhDXYNormalized[1]->Fill(b[0]/bRes[0]);
735       fhDXYvsDZNormalized[1]->Fill(b[1]/bRes[1], b[0]/bRes[0]);
736       fhNSigmaToVertex[1]->Fill(nSigmaToVertex);
737     }
738   }
739
740   return kTRUE;
741 }
742
743 //____________________________________________________________________
744 TObjArray* AliESDtrackCuts::GetAcceptedTracks(AliESD* esd)
745 {
746   //
747   // returns an array of all tracks that pass the cuts
748   //
749
750   TObjArray* acceptedTracks = new TObjArray();
751
752   // loop over esd tracks
753   for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
754     AliESDtrack* track = esd->GetTrack(iTrack);
755
756     if (AcceptTrack(track))
757       acceptedTracks->Add(track);
758   }
759
760   return acceptedTracks;
761 }
762
763 //____________________________________________________________________
764 Int_t AliESDtrackCuts::CountAcceptedTracks(AliESD* esd)
765 {
766   //
767   // returns an the number of tracks that pass the cuts
768   //
769
770   Int_t count = 0;
771
772   // loop over esd tracks
773   for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
774     AliESDtrack* track = esd->GetTrack(iTrack);
775
776     if (AcceptTrack(track))
777       count++;
778   }
779
780   return count;
781 }
782
783 //____________________________________________________________________
784 TObjArray* AliESDtrackCuts::GetAcceptedTracks(AliESDEvent* esd)
785 {
786   //
787   // returns an array of all tracks that pass the cuts
788   //
789
790   TObjArray* acceptedTracks = new TObjArray();
791
792   // loop over esd tracks
793   for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
794     AliESDtrack* track = esd->GetTrack(iTrack);
795
796     if (AcceptTrack(track))
797       acceptedTracks->Add(track);
798   }
799
800   return acceptedTracks;
801 }
802
803 //____________________________________________________________________
804 Int_t AliESDtrackCuts::CountAcceptedTracks(AliESDEvent* esd)
805 {
806   //
807   // returns an the number of tracks that pass the cuts
808   //
809
810   Int_t count = 0;
811
812   // loop over esd tracks
813   for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
814     AliESDtrack* track = esd->GetTrack(iTrack);
815
816     if (AcceptTrack(track))
817       count++;
818   }
819
820   return count;
821 }
822
823 //____________________________________________________________________
824  void AliESDtrackCuts::DefineHistograms(Int_t color) {
825    // 
826    // diagnostics histograms are defined
827    // 
828
829    fHistogramsOn=kTRUE;
830    
831    Bool_t oldStatus = TH1::AddDirectoryStatus();
832    TH1::AddDirectory(kFALSE);
833    
834    //###################################################################################
835    // defining histograms
836
837    fhCutStatistics = new TH1F("cut_statistics","cut statistics",kNCuts+4,-0.5,kNCuts+3.5);
838
839    fhCutStatistics->GetXaxis()->SetBinLabel(1,"n tracks");
840    fhCutStatistics->GetXaxis()->SetBinLabel(2,"n cut tracks");
841
842    fhCutCorrelation = new TH2F("cut_correlation","cut correlation",kNCuts,-0.5,kNCuts-0.5,kNCuts,-0.5,kNCuts-0.5);;
843   
844    for (Int_t i=0; i<kNCuts; i++) {
845      fhCutStatistics->GetXaxis()->SetBinLabel(i+4,fgkCutNames[i]);
846      fhCutCorrelation->GetXaxis()->SetBinLabel(i+1,fgkCutNames[i]);
847      fhCutCorrelation->GetYaxis()->SetBinLabel(i+1,fgkCutNames[i]);
848    } 
849
850   fhCutStatistics  ->SetLineColor(color);
851   fhCutCorrelation ->SetLineColor(color);
852   fhCutStatistics  ->SetLineWidth(2);
853   fhCutCorrelation ->SetLineWidth(2);
854
855   Char_t str[256];
856   for (Int_t i=0; i<2; i++) {
857     if (i==0) sprintf(str," ");
858     else sprintf(str,"_cut");
859
860     fhNClustersITS[i]        = new TH1F(Form("nClustersITS%s",str)     ,"",8,-0.5,7.5);
861     fhNClustersTPC[i]        = new TH1F(Form("nClustersTPC%s",str)     ,"",165,-0.5,164.5);
862     fhChi2PerClusterITS[i]   = new TH1F(Form("chi2PerClusterITS%s",str),"",500,0,10);
863     fhChi2PerClusterTPC[i]   = new TH1F(Form("chi2PerClusterTPC%s",str),"",500,0,10);
864
865     fhC11[i]                 = new  TH1F(Form("covMatrixDiagonal11%s",str),"",2000,0,20);
866     fhC22[i]                 = new  TH1F(Form("covMatrixDiagonal22%s",str),"",2000,0,20);
867     fhC33[i]                 = new  TH1F(Form("covMatrixDiagonal33%s",str),"",1000,0,1);
868     fhC44[i]                 = new  TH1F(Form("covMatrixDiagonal44%s",str),"",1000,0,5);
869     fhC55[i]                 = new  TH1F(Form("covMatrixDiagonal55%s",str),"",1000,0,5);
870
871     fhDXY[i]                 = new  TH1F(Form("dXY%s",str)    ,"",500,-10,10);
872     fhDZ[i]                  = new  TH1F(Form("dZ%s",str)     ,"",500,-10,10);
873     fhDXYvsDZ[i]             = new  TH2F(Form("dXYvsDZ%s",str),"",200,-10,10,200,-10,10);
874
875     fhDXYNormalized[i]       = new  TH1F(Form("dXYNormalized%s",str)    ,"",500,-10,10);
876     fhDZNormalized[i]        = new  TH1F(Form("dZNormalized%s",str)     ,"",500,-10,10);
877     fhDXYvsDZNormalized[i]   = new  TH2F(Form("dXYvsDZNormalized%s",str),"",200,-10,10,200,-10,10);
878
879     fhNSigmaToVertex[i]      = new  TH1F(Form("nSigmaToVertex%s",str),"",500,0,50);
880
881     fhPt[i]                  = new TH1F(Form("pt%s",str)     ,"p_{T} distribution;p_{T} (GeV/c)",500,0.0,100.0);
882     fhEta[i]                 = new TH1F(Form("eta%s",str)     ,"#eta distribution;#eta",40,-2.0,2.0);
883     
884     fhNClustersITS[i]->SetTitle("n ITS clusters");
885     fhNClustersTPC[i]->SetTitle("n TPC clusters");
886     fhChi2PerClusterITS[i]->SetTitle("#Chi^{2} per ITS cluster");
887     fhChi2PerClusterTPC[i]->SetTitle("#Chi^{2} per TPC cluster");
888
889     fhC11[i]->SetTitle("cov 11 : #sigma_{y}^{2} [cm^{2}]");
890     fhC22[i]->SetTitle("cov 22 : #sigma_{z}^{2} [cm^{2}]");
891     fhC33[i]->SetTitle("cov 33 : #sigma_{sin(#phi)}^{2}");
892     fhC44[i]->SetTitle("cov 44 : #sigma_{tan(#theta_{dip})}^{2}");
893     fhC55[i]->SetTitle("cov 55 : #sigma_{1/p_{T}}^{2} [(c/GeV)^2]");
894
895     fhDXY[i]->SetTitle("transverse impact parameter");
896     fhDZ[i]->SetTitle("longitudinal impact parameter");
897     fhDXYvsDZ[i]->SetTitle("longitudinal impact parameter");
898     fhDXYvsDZ[i]->SetYTitle("transverse impact parameter");
899
900     fhDXYNormalized[i]->SetTitle("normalized trans impact par");
901     fhDZNormalized[i]->SetTitle("normalized long impact par");
902     fhDXYvsDZNormalized[i]->SetTitle("normalized long impact par");
903     fhDXYvsDZNormalized[i]->SetYTitle("normalized trans impact par");
904     fhNSigmaToVertex[i]->SetTitle("n #sigma to vertex");
905
906     fhNClustersITS[i]->SetLineColor(color);   fhNClustersITS[i]->SetLineWidth(2);
907     fhNClustersTPC[i]->SetLineColor(color);   fhNClustersTPC[i]->SetLineWidth(2);
908     fhChi2PerClusterITS[i]->SetLineColor(color);   fhChi2PerClusterITS[i]->SetLineWidth(2);
909     fhChi2PerClusterTPC[i]->SetLineColor(color);   fhChi2PerClusterTPC[i]->SetLineWidth(2);
910
911     fhC11[i]->SetLineColor(color);   fhC11[i]->SetLineWidth(2);
912     fhC22[i]->SetLineColor(color);   fhC22[i]->SetLineWidth(2);
913     fhC33[i]->SetLineColor(color);   fhC33[i]->SetLineWidth(2);
914     fhC44[i]->SetLineColor(color);   fhC44[i]->SetLineWidth(2);
915     fhC55[i]->SetLineColor(color);   fhC55[i]->SetLineWidth(2);
916
917     fhDXY[i]->SetLineColor(color);   fhDXY[i]->SetLineWidth(2);
918     fhDZ[i]->SetLineColor(color);   fhDZ[i]->SetLineWidth(2);
919
920     fhDXYNormalized[i]->SetLineColor(color);   fhDXYNormalized[i]->SetLineWidth(2);
921     fhDZNormalized[i]->SetLineColor(color);    fhDZNormalized[i]->SetLineWidth(2);
922     fhNSigmaToVertex[i]->SetLineColor(color);  fhNSigmaToVertex[i]->SetLineWidth(2); 
923   }
924
925   // The number of sigmas to the vertex is per definition gaussian
926   ffDTheoretical = new TF1("nSigmaToVertexTheoretical","([0]/2.506628274)*exp(-(x**2)/2)",0,50);
927   ffDTheoretical->SetParameter(0,1);
928   
929   TH1::AddDirectory(oldStatus);
930 }
931
932 //____________________________________________________________________
933 Bool_t AliESDtrackCuts::LoadHistograms(const Char_t* dir)
934 {
935   //
936   // loads the histograms from a file
937   // if dir is empty a directory with the name of this object is taken (like in SaveHistogram)
938   //
939
940   if (!dir)
941     dir = GetName();
942
943   if (!gDirectory->cd(dir))
944     return kFALSE;
945
946   ffDTheoretical = dynamic_cast<TF1*> (gDirectory->Get("nSigmaToVertexTheory"));
947
948   fhCutStatistics = dynamic_cast<TH1F*> (gDirectory->Get("cut_statistics"));
949   fhCutCorrelation = dynamic_cast<TH2F*> (gDirectory->Get("cut_correlation"));
950
951   Char_t str[5];
952   for (Int_t i=0; i<2; i++) {
953     if (i==0)
954     {
955       gDirectory->cd("before_cuts");
956       str[0] = 0;
957     }
958     else
959     {
960       gDirectory->cd("after_cuts");
961       sprintf(str,"_cut");
962     }
963
964     fhNClustersITS[i]      = dynamic_cast<TH1F*> (gDirectory->Get(Form("nClustersITS%s",str)     ));
965     fhNClustersTPC[i]      = dynamic_cast<TH1F*> (gDirectory->Get(Form("nClustersTPC%s",str)     ));
966     fhChi2PerClusterITS[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("chi2PerClusterITS%s",str)));
967     fhChi2PerClusterTPC[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("chi2PerClusterTPC%s",str)));
968
969     fhC11[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("covMatrixDiagonal11%s",str)));
970     fhC22[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("covMatrixDiagonal22%s",str)));
971     fhC33[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("covMatrixDiagonal33%s",str)));
972     fhC44[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("covMatrixDiagonal44%s",str)));
973     fhC55[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("covMatrixDiagonal55%s",str)));
974
975     fhDXY[i] =     dynamic_cast<TH1F*> (gDirectory->Get(Form("dXY%s",str)    ));
976     fhDZ[i] =      dynamic_cast<TH1F*> (gDirectory->Get(Form("dZ%s",str)     ));
977     fhDXYvsDZ[i] = dynamic_cast<TH2F*> (gDirectory->Get(Form("dXYvsDZ%s",str)));
978
979     fhDXYNormalized[i] =     dynamic_cast<TH1F*> (gDirectory->Get(Form("dXYNormalized%s",str)    ));
980     fhDZNormalized[i] =      dynamic_cast<TH1F*> (gDirectory->Get(Form("dZNormalized%s",str)     ));
981     fhDXYvsDZNormalized[i] = dynamic_cast<TH2F*> (gDirectory->Get(Form("dXYvsDZNormalized%s",str)));
982     fhNSigmaToVertex[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("nSigmaToVertex%s",str)));
983
984     fhPt[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("pt%s",str)));
985     fhEta[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("eta%s",str)));
986
987     gDirectory->cd("../");
988   }
989
990   gDirectory->cd("..");
991
992   return kTRUE;
993 }
994
995 //____________________________________________________________________
996 void AliESDtrackCuts::SaveHistograms(const Char_t* dir) {
997   //
998   // saves the histograms in a directory (dir)
999   //
1000
1001   if (!fHistogramsOn) {
1002     AliDebug(0, "Histograms not on - cannot save histograms!!!");
1003     return;
1004   }
1005
1006   if (!dir)
1007     dir = GetName();
1008
1009   gDirectory->mkdir(dir);
1010   gDirectory->cd(dir);
1011
1012   gDirectory->mkdir("before_cuts");
1013   gDirectory->mkdir("after_cuts");
1014  
1015   // a factor of 2 is needed since n sigma is positive
1016   ffDTheoretical->SetParameter(0,2*fhNSigmaToVertex[0]->Integral("width"));
1017   ffDTheoretical->Write("nSigmaToVertexTheory");
1018
1019   fhCutStatistics->Write();
1020   fhCutCorrelation->Write();
1021
1022   for (Int_t i=0; i<2; i++) {
1023     if (i==0)
1024       gDirectory->cd("before_cuts");
1025     else
1026       gDirectory->cd("after_cuts");
1027
1028     fhNClustersITS[i]        ->Write();
1029     fhNClustersTPC[i]        ->Write();
1030     fhChi2PerClusterITS[i]   ->Write();
1031     fhChi2PerClusterTPC[i]   ->Write();
1032
1033     fhC11[i]                 ->Write();
1034     fhC22[i]                 ->Write();
1035     fhC33[i]                 ->Write();
1036     fhC44[i]                 ->Write();
1037     fhC55[i]                 ->Write();
1038
1039     fhDXY[i]                 ->Write();
1040     fhDZ[i]                  ->Write();
1041     fhDXYvsDZ[i]             ->Write();
1042
1043     fhDXYNormalized[i]       ->Write();
1044     fhDZNormalized[i]        ->Write();
1045     fhDXYvsDZNormalized[i]   ->Write();
1046     fhNSigmaToVertex[i]      ->Write();
1047
1048     fhPt[i]                  ->Write();
1049     fhEta[i]                 ->Write();
1050     
1051     gDirectory->cd("../");
1052   }
1053
1054   gDirectory->cd("../");
1055 }
1056
1057 //____________________________________________________________________
1058 void AliESDtrackCuts::DrawHistograms()
1059 {
1060   // draws some histograms
1061
1062   TCanvas* canvas1 = new TCanvas(Form("%s_1", GetName()), "Track Quality Results1", 800, 800);
1063   canvas1->Divide(2, 2);
1064
1065   canvas1->cd(1);
1066   fhNClustersTPC[0]->SetStats(kFALSE);
1067   fhNClustersTPC[0]->Draw();
1068
1069   canvas1->cd(2);
1070   fhChi2PerClusterTPC[0]->SetStats(kFALSE);
1071   fhChi2PerClusterTPC[0]->Draw();
1072
1073   canvas1->cd(3);
1074   fhNSigmaToVertex[0]->SetStats(kFALSE);
1075   fhNSigmaToVertex[0]->GetXaxis()->SetRangeUser(0, 10);
1076   fhNSigmaToVertex[0]->Draw();
1077
1078   canvas1->SaveAs(Form("%s_%s.gif", GetName(), canvas1->GetName()));
1079
1080   TCanvas* canvas2 = new TCanvas(Form("%s_2", GetName()), "Track Quality Results2", 1200, 800);
1081   canvas2->Divide(3, 2);
1082
1083   canvas2->cd(1);
1084   fhC11[0]->SetStats(kFALSE);
1085   gPad->SetLogy();
1086   fhC11[0]->Draw();
1087
1088   canvas2->cd(2);
1089   fhC22[0]->SetStats(kFALSE);
1090   gPad->SetLogy();
1091   fhC22[0]->Draw();
1092
1093   canvas2->cd(3);
1094   fhC33[0]->SetStats(kFALSE);
1095   gPad->SetLogy();
1096   fhC33[0]->Draw();
1097
1098   canvas2->cd(4);
1099   fhC44[0]->SetStats(kFALSE);
1100   gPad->SetLogy();
1101   fhC44[0]->Draw();
1102
1103   canvas2->cd(5);
1104   fhC55[0]->SetStats(kFALSE);
1105   gPad->SetLogy();
1106   fhC55[0]->Draw();
1107
1108   canvas2->SaveAs(Form("%s_%s.gif", GetName(), canvas2->GetName()));
1109
1110   TCanvas* canvas3 = new TCanvas(Form("%s_3", GetName()), "Track Quality Results3", 1200, 800);
1111   canvas3->Divide(3, 2);
1112
1113   canvas3->cd(1);
1114   fhDXY[0]->SetStats(kFALSE);
1115   gPad->SetLogy();
1116   fhDXY[0]->Draw();
1117
1118   canvas3->cd(2);
1119   fhDZ[0]->SetStats(kFALSE);
1120   gPad->SetLogy();
1121   fhDZ[0]->Draw();
1122
1123   canvas3->cd(3);
1124   fhDXYvsDZ[0]->SetStats(kFALSE);
1125   gPad->SetLogz();
1126   gPad->SetRightMargin(0.15);
1127   fhDXYvsDZ[0]->Draw("COLZ");
1128
1129   canvas3->cd(4);
1130   fhDXYNormalized[0]->SetStats(kFALSE);
1131   gPad->SetLogy();
1132   fhDXYNormalized[0]->Draw();
1133
1134   canvas3->cd(5);
1135   fhDZNormalized[0]->SetStats(kFALSE);
1136   gPad->SetLogy();
1137   fhDZNormalized[0]->Draw();
1138
1139   canvas3->cd(6);
1140   fhDXYvsDZNormalized[0]->SetStats(kFALSE);
1141   gPad->SetLogz();
1142   gPad->SetRightMargin(0.15);
1143   fhDXYvsDZNormalized[0]->Draw("COLZ");
1144
1145   canvas3->SaveAs(Form("%s_%s.gif", GetName(), canvas3->GetName()));
1146
1147   TCanvas* canvas4 = new TCanvas(Form("%s_4", GetName()), "Track Quality Results4", 800, 500);
1148   canvas4->Divide(2, 1);
1149
1150   canvas4->cd(1);
1151   fhCutStatistics->SetStats(kFALSE);
1152   fhCutStatistics->LabelsOption("v");
1153   gPad->SetBottomMargin(0.3);
1154   fhCutStatistics->Draw();
1155
1156   canvas4->cd(2);
1157   fhCutCorrelation->SetStats(kFALSE);
1158   fhCutCorrelation->LabelsOption("v");
1159   gPad->SetBottomMargin(0.3);
1160   gPad->SetLeftMargin(0.3);
1161   fhCutCorrelation->Draw("COLZ");
1162
1163   canvas4->SaveAs(Form("%s_%s.gif", GetName(), canvas4->GetName()));
1164
1165   /*canvas->cd(1);
1166   fhDXYvsDZNormalized[0]->SetStats(kFALSE);
1167   fhDXYvsDZNormalized[0]->DrawCopy("COLZ");
1168
1169   canvas->cd(2);
1170   fhNClustersTPC[0]->SetStats(kFALSE);
1171   fhNClustersTPC[0]->DrawCopy();
1172
1173   canvas->cd(3);
1174   fhChi2PerClusterITS[0]->SetStats(kFALSE);
1175   fhChi2PerClusterITS[0]->DrawCopy();
1176   fhChi2PerClusterITS[1]->SetLineColor(2);
1177   fhChi2PerClusterITS[1]->DrawCopy("SAME");
1178
1179   canvas->cd(4);
1180   fhChi2PerClusterTPC[0]->SetStats(kFALSE);
1181   fhChi2PerClusterTPC[0]->DrawCopy();
1182   fhChi2PerClusterTPC[1]->SetLineColor(2);
1183   fhChi2PerClusterTPC[1]->DrawCopy("SAME");*/
1184 }
1185