]> git.uio.no Git - u/mrichter/AliRoot.git/blob - ANALYSIS/AliESDtrackCuts.cxx
More user-friendly selection for ITS standalone tracks and allows to select tracks...
[u/mrichter/AliRoot.git] / ANALYSIS / 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: AliESDtrackCuts.cxx 24534 2008-03-16 22:22:11Z fca $ */
17
18 #include "AliESDtrackCuts.h"
19
20 #include <AliESDtrack.h>
21 #include <AliESDVertex.h>
22 #include <AliESDEvent.h>
23 #include <AliLog.h>
24
25 #include <TTree.h>
26 #include <TCanvas.h>
27 #include <TDirectory.h>
28 #include <TH2F.h>
29 #include <TF1.h>
30
31 //____________________________________________________________________
32 ClassImp(AliESDtrackCuts)
33
34 // Cut names
35 const Char_t* AliESDtrackCuts::fgkCutNames[kNCuts] = {
36  "require TPC refit",
37  "require TPC standalone",
38  "require ITS refit",
39  "n clusters TPC",
40  "n clusters ITS",
41  "#Chi^{2}/cluster TPC",
42  "#Chi^{2}/cluster ITS",
43  "cov 11",
44  "cov 22",
45  "cov 33",
46  "cov 44",
47  "cov 55",
48  "trk-to-vtx",
49  "trk-to-vtx failed",
50  "kink daughters",
51  "p",
52  "p_{T}",
53  "p_{x}",
54  "p_{y}",
55  "p_{z}",
56  "eta",
57  "y",
58  "trk-to-vtx max dca 2D absolute",
59  "trk-to-vtx max dca xy absolute",
60  "trk-to-vtx max dca z absolute",
61  "trk-to-vtx min dca 2D absolute",
62  "trk-to-vtx min dca xy absolute",
63  "trk-to-vtx min dca z absolute",
64  "SPD cluster requirement",
65  "SDD cluster requirement",
66  "SSD cluster requirement",
67  "require ITS stand-alone",
68  "rel 1/pt uncertainty",
69  "require ITS Pid"
70 };
71
72 //____________________________________________________________________
73 AliESDtrackCuts::AliESDtrackCuts(const Char_t* name, const Char_t* title) : AliAnalysisCuts(name,title),
74   fCutMinNClusterTPC(0),
75   fCutMinNClusterITS(0),
76   fCutMaxChi2PerClusterTPC(0),
77   fCutMaxChi2PerClusterITS(0),
78   fCutMaxC11(0),
79   fCutMaxC22(0),
80   fCutMaxC33(0),
81   fCutMaxC44(0),
82   fCutMaxC55(0),
83   fCutMaxRel1PtUncertainty(0),
84   fCutAcceptKinkDaughters(0),
85   fCutAcceptSharedTPCClusters(0),
86   fCutMaxFractionSharedTPCClusters(0),
87   fCutRequireTPCRefit(0),
88   fCutRequireTPCStandAlone(0),
89   fCutRequireITSRefit(0), 
90   fCutRequireITSPid(0),
91   fCutRequireITSStandAlone(0),
92   fCutRequireITSpureSA(0),
93   fCutNsigmaToVertex(0),
94   fCutSigmaToVertexRequired(0),
95   fCutMaxDCAToVertexXY(0),
96   fCutMaxDCAToVertexZ(0),
97   fCutMinDCAToVertexXY(0),
98   fCutMinDCAToVertexZ(0),
99   fCutMaxDCAToVertexXYPtDep(""),
100   fCutMaxDCAToVertexZPtDep(""),
101   fCutMinDCAToVertexXYPtDep(""),
102   fCutMinDCAToVertexZPtDep(""),
103   f1CutMaxDCAToVertexXYPtDep(0x0),
104   f1CutMaxDCAToVertexZPtDep(0x0),
105   f1CutMinDCAToVertexXYPtDep(0x0),
106   f1CutMinDCAToVertexZPtDep(0x0),
107   fCutDCAToVertex2D(0),
108   fPMin(0),
109   fPMax(0),
110   fPtMin(0),
111   fPtMax(0),
112   fPxMin(0),
113   fPxMax(0),
114   fPyMin(0),
115   fPyMax(0),
116   fPzMin(0),
117   fPzMax(0),
118   fEtaMin(0),
119   fEtaMax(0),
120   fRapMin(0),
121   fRapMax(0),
122   fHistogramsOn(0),
123   ffDTheoretical(0),
124   fhCutStatistics(0),         
125   fhCutCorrelation(0)
126 {
127   //
128   // constructor
129   //
130
131   Init();
132
133   //##############################################################################
134   // setting default cuts
135   SetMinNClustersTPC();
136   SetMinNClustersITS();
137   SetMaxChi2PerClusterTPC();
138   SetMaxChi2PerClusterITS();                                
139   SetMaxCovDiagonalElements();
140   SetMaxRel1PtUncertainty();
141   SetRequireTPCRefit();
142   SetRequireTPCStandAlone();
143   SetRequireITSRefit();
144   SetRequireITSPid(kFALSE);
145   SetRequireITSStandAlone(kFALSE);
146   SetRequireITSPureStandAlone(kFALSE);
147   SetAcceptKinkDaughters();
148   SetMaxNsigmaToVertex();
149   SetMaxDCAToVertexXY();
150   SetMaxDCAToVertexZ();
151   SetDCAToVertex2D();
152   SetMinDCAToVertexXY();
153   SetMinDCAToVertexZ();
154   SetPRange();
155   SetPtRange();
156   SetPxRange();
157   SetPyRange();
158   SetPzRange();
159   SetEtaRange();
160   SetRapRange();
161   SetClusterRequirementITS(kSPD);
162   SetClusterRequirementITS(kSDD);
163   SetClusterRequirementITS(kSSD);
164
165   SetHistogramsOn();
166 }
167
168 //_____________________________________________________________________________
169 AliESDtrackCuts::AliESDtrackCuts(const AliESDtrackCuts &c) : AliAnalysisCuts(c),
170   fCutMinNClusterTPC(0),
171   fCutMinNClusterITS(0),
172   fCutMaxChi2PerClusterTPC(0),
173   fCutMaxChi2PerClusterITS(0),
174   fCutMaxC11(0),
175   fCutMaxC22(0),
176   fCutMaxC33(0),
177   fCutMaxC44(0),
178   fCutMaxC55(0),
179   fCutMaxRel1PtUncertainty(0),
180   fCutAcceptKinkDaughters(0),
181   fCutAcceptSharedTPCClusters(0),
182   fCutMaxFractionSharedTPCClusters(0),
183   fCutRequireTPCRefit(0),
184   fCutRequireTPCStandAlone(0),
185   fCutRequireITSRefit(0),
186   fCutRequireITSPid(0),
187   fCutRequireITSStandAlone(0),
188   fCutRequireITSpureSA(0),
189   fCutNsigmaToVertex(0),
190   fCutSigmaToVertexRequired(0),
191   fCutMaxDCAToVertexXY(0),
192   fCutMaxDCAToVertexZ(0),
193   fCutMinDCAToVertexXY(0),
194   fCutMinDCAToVertexZ(0),
195   fCutMaxDCAToVertexXYPtDep(""),
196   fCutMaxDCAToVertexZPtDep(""),
197   fCutMinDCAToVertexXYPtDep(""),
198   fCutMinDCAToVertexZPtDep(""),
199   f1CutMaxDCAToVertexXYPtDep(0x0),
200   f1CutMaxDCAToVertexZPtDep(0x0),
201   f1CutMinDCAToVertexXYPtDep(0x0),
202   f1CutMinDCAToVertexZPtDep(0x0),
203   fCutDCAToVertex2D(0),
204   fPMin(0),
205   fPMax(0),
206   fPtMin(0),
207   fPtMax(0),
208   fPxMin(0),
209   fPxMax(0),
210   fPyMin(0),
211   fPyMax(0),
212   fPzMin(0),
213   fPzMax(0),
214   fEtaMin(0),
215   fEtaMax(0),
216   fRapMin(0),
217   fRapMax(0),
218   fHistogramsOn(0),
219   ffDTheoretical(0),                                 
220   fhCutStatistics(0),         
221   fhCutCorrelation(0)
222 {
223   //
224   // copy constructor
225   //
226
227   ((AliESDtrackCuts &) c).Copy(*this);
228 }
229
230 AliESDtrackCuts::~AliESDtrackCuts()
231 {
232   //
233   // destructor
234   //
235
236   for (Int_t i=0; i<2; i++) {
237     
238     if (fhNClustersITS[i])
239       delete fhNClustersITS[i];            
240     if (fhNClustersTPC[i])
241       delete fhNClustersTPC[i];                
242     if (fhChi2PerClusterITS[i])
243       delete fhChi2PerClusterITS[i];       
244     if (fhChi2PerClusterTPC[i])
245       delete fhChi2PerClusterTPC[i];       
246     if (fhC11[i])
247       delete fhC11[i];                     
248     if (fhC22[i])
249       delete fhC22[i];                     
250     if (fhC33[i])
251       delete fhC33[i];                     
252     if (fhC44[i])
253       delete fhC44[i];                     
254     if (fhC55[i])
255       delete fhC55[i];
256
257     if (fhRel1PtUncertainty[i])
258       delete fhRel1PtUncertainty[i];
259     
260     if (fhDXY[i])
261       delete fhDXY[i];                     
262     if (fhDZ[i])
263       delete fhDZ[i];
264     if (fhDXYDZ[i])
265       delete fhDXYDZ[i];
266     if (fhDXYvsDZ[i])
267       delete fhDXYvsDZ[i];
268
269     if (fhDXYNormalized[i])
270       delete fhDXYNormalized[i];           
271     if (fhDZNormalized[i])
272       delete fhDZNormalized[i];
273     if (fhDXYvsDZNormalized[i])
274       delete fhDXYvsDZNormalized[i];
275     if (fhNSigmaToVertex[i])
276       delete fhNSigmaToVertex[i];
277     if (fhPt[i])
278       delete fhPt[i];
279     if (fhEta[i])
280       delete fhEta[i];
281   }
282
283   if(f1CutMaxDCAToVertexXYPtDep)delete f1CutMaxDCAToVertexXYPtDep;
284   f1CutMaxDCAToVertexXYPtDep = 0;
285   if( f1CutMaxDCAToVertexZPtDep) delete  f1CutMaxDCAToVertexZPtDep;
286   f1CutMaxDCAToVertexZPtDep = 0;
287   if( f1CutMinDCAToVertexXYPtDep)delete  f1CutMinDCAToVertexXYPtDep;
288   f1CutMinDCAToVertexXYPtDep = 0;
289   if(f1CutMinDCAToVertexZPtDep)delete  f1CutMinDCAToVertexZPtDep; 
290   f1CutMinDCAToVertexZPtDep = 0;
291
292
293   if (ffDTheoretical)
294     delete ffDTheoretical;
295
296   if (fhCutStatistics)
297     delete fhCutStatistics;             
298   if (fhCutCorrelation)
299     delete fhCutCorrelation;            
300 }
301
302 void AliESDtrackCuts::Init()
303 {
304   //
305   // sets everything to zero
306   //
307
308   fCutMinNClusterTPC = 0;
309   fCutMinNClusterITS = 0;
310
311   fCutMaxChi2PerClusterTPC = 0;
312   fCutMaxChi2PerClusterITS = 0;
313   
314   for (Int_t i = 0; i < 3; i++)
315         fCutClusterRequirementITS[i] = kOff;
316
317   fCutMaxC11 = 0;
318   fCutMaxC22 = 0;
319   fCutMaxC33 = 0;
320   fCutMaxC44 = 0;
321   fCutMaxC55 = 0;
322   
323   fCutMaxRel1PtUncertainty = 0;
324
325   fCutAcceptKinkDaughters = 0;
326   fCutAcceptSharedTPCClusters = 0;
327   fCutMaxFractionSharedTPCClusters = 0;
328   fCutRequireTPCRefit = 0;
329   fCutRequireTPCStandAlone = 0;
330   fCutRequireITSRefit = 0;
331   fCutRequireITSPid = 0;
332   fCutRequireITSStandAlone = 0;
333   fCutRequireITSpureSA = 0;
334
335   fCutNsigmaToVertex = 0;
336   fCutSigmaToVertexRequired = 0;
337   fCutMaxDCAToVertexXY = 0;
338   fCutMaxDCAToVertexZ = 0;
339   fCutDCAToVertex2D = 0;
340   fCutMinDCAToVertexXY = 0;
341   fCutMinDCAToVertexZ = 0;
342   fCutMaxDCAToVertexXYPtDep = "";
343   fCutMaxDCAToVertexZPtDep = "";
344   fCutMinDCAToVertexXYPtDep = "";
345   fCutMinDCAToVertexZPtDep = "";
346
347   if(f1CutMaxDCAToVertexXYPtDep)delete f1CutMaxDCAToVertexXYPtDep;
348   f1CutMaxDCAToVertexXYPtDep = 0;
349   if( f1CutMaxDCAToVertexXYPtDep) delete  f1CutMaxDCAToVertexXYPtDep;
350   f1CutMaxDCAToVertexXYPtDep = 0;
351   if( f1CutMaxDCAToVertexZPtDep) delete  f1CutMaxDCAToVertexZPtDep;
352   f1CutMaxDCAToVertexZPtDep = 0;
353   if( f1CutMinDCAToVertexXYPtDep)delete  f1CutMinDCAToVertexXYPtDep;
354   f1CutMinDCAToVertexXYPtDep = 0;
355   if(f1CutMinDCAToVertexZPtDep)delete f1CutMinDCAToVertexZPtDep;
356   f1CutMinDCAToVertexZPtDep = 0;
357
358   
359   fPMin = 0;
360   fPMax = 0;
361   fPtMin = 0;
362   fPtMax = 0;
363   fPxMin = 0;
364   fPxMax = 0;
365   fPyMin = 0;
366   fPyMax = 0;
367   fPzMin = 0;
368   fPzMax = 0;
369   fEtaMin = 0;
370   fEtaMax = 0;
371   fRapMin = 0;
372   fRapMax = 0;
373
374   fHistogramsOn = kFALSE;
375
376   for (Int_t i=0; i<2; ++i)
377   {
378     fhNClustersITS[i] = 0;
379     fhNClustersTPC[i] = 0;
380
381     fhChi2PerClusterITS[i] = 0;
382     fhChi2PerClusterTPC[i] = 0;
383
384     fhC11[i] = 0;
385     fhC22[i] = 0;
386     fhC33[i] = 0;
387     fhC44[i] = 0;
388     fhC55[i] = 0;
389
390     fhRel1PtUncertainty[i] = 0;
391
392     fhDXY[i] = 0;
393     fhDZ[i] = 0;
394     fhDXYDZ[i] = 0;
395     fhDXYvsDZ[i] = 0;
396
397     fhDXYNormalized[i] = 0;
398     fhDZNormalized[i] = 0;
399     fhDXYvsDZNormalized[i] = 0;
400     fhNSigmaToVertex[i] = 0;
401     
402     fhPt[i] = 0;
403     fhEta[i] = 0;
404   }
405   ffDTheoretical = 0;
406
407   fhCutStatistics = 0;
408   fhCutCorrelation = 0;
409 }
410
411 //_____________________________________________________________________________
412 AliESDtrackCuts &AliESDtrackCuts::operator=(const AliESDtrackCuts &c)
413 {
414   //
415   // Assignment operator
416   //
417
418   if (this != &c) ((AliESDtrackCuts &) c).Copy(*this);
419   return *this;
420 }
421
422 //_____________________________________________________________________________
423 void AliESDtrackCuts::Copy(TObject &c) const
424 {
425   //
426   // Copy function
427   //
428
429   AliESDtrackCuts& target = (AliESDtrackCuts &) c;
430
431   target.Init();
432
433   target.fCutMinNClusterTPC = fCutMinNClusterTPC;
434   target.fCutMinNClusterITS = fCutMinNClusterITS;
435
436   target.fCutMaxChi2PerClusterTPC = fCutMaxChi2PerClusterTPC;
437   target.fCutMaxChi2PerClusterITS = fCutMaxChi2PerClusterITS;
438
439   for (Int_t i = 0; i < 3; i++)
440     target.fCutClusterRequirementITS[i] = fCutClusterRequirementITS[i];
441
442   target.fCutMaxC11 = fCutMaxC11;
443   target.fCutMaxC22 = fCutMaxC22;
444   target.fCutMaxC33 = fCutMaxC33;
445   target.fCutMaxC44 = fCutMaxC44;
446   target.fCutMaxC55 = fCutMaxC55;
447
448   target.fCutMaxRel1PtUncertainty = fCutMaxRel1PtUncertainty;
449
450   target.fCutAcceptKinkDaughters = fCutAcceptKinkDaughters;
451   target.fCutAcceptSharedTPCClusters = fCutAcceptSharedTPCClusters;
452   target.fCutMaxFractionSharedTPCClusters = fCutMaxFractionSharedTPCClusters;
453   target.fCutRequireTPCRefit = fCutRequireTPCRefit;
454   target.fCutRequireTPCStandAlone = fCutRequireTPCStandAlone;
455   target.fCutRequireITSRefit = fCutRequireITSRefit;
456   target.fCutRequireITSPid = fCutRequireITSPid;
457   target.fCutRequireITSStandAlone = fCutRequireITSStandAlone;
458   target.fCutRequireITSpureSA = fCutRequireITSpureSA;
459
460   target.fCutNsigmaToVertex = fCutNsigmaToVertex;
461   target.fCutSigmaToVertexRequired = fCutSigmaToVertexRequired;
462   target.fCutMaxDCAToVertexXY = fCutMaxDCAToVertexXY;
463   target.fCutMaxDCAToVertexZ = fCutMaxDCAToVertexZ;
464   target.fCutDCAToVertex2D = fCutDCAToVertex2D;
465   target.fCutMinDCAToVertexXY = fCutMinDCAToVertexXY;
466   target.fCutMinDCAToVertexZ = fCutMinDCAToVertexZ;
467
468   target.fCutMaxDCAToVertexXYPtDep = fCutMaxDCAToVertexXYPtDep;
469   target.SetMaxDCAToVertexXYPtDep(fCutMaxDCAToVertexXYPtDep.Data());
470
471   target.fCutMaxDCAToVertexZPtDep = fCutMaxDCAToVertexZPtDep;
472   target.SetMaxDCAToVertexZPtDep(fCutMaxDCAToVertexZPtDep.Data());
473
474   target.fCutMinDCAToVertexXYPtDep = fCutMinDCAToVertexXYPtDep;
475   target.SetMinDCAToVertexXYPtDep(fCutMinDCAToVertexXYPtDep.Data());
476
477   target.fCutMinDCAToVertexZPtDep = fCutMinDCAToVertexZPtDep;
478   target.SetMinDCAToVertexZPtDep(fCutMinDCAToVertexZPtDep.Data());
479
480   target.fPMin = fPMin;
481   target.fPMax = fPMax;
482   target.fPtMin = fPtMin;
483   target.fPtMax = fPtMax;
484   target.fPxMin = fPxMin;
485   target.fPxMax = fPxMax;
486   target.fPyMin = fPyMin;
487   target.fPyMax = fPyMax;
488   target.fPzMin = fPzMin;
489   target.fPzMax = fPzMax;
490   target.fEtaMin = fEtaMin;
491   target.fEtaMax = fEtaMax;
492   target.fRapMin = fRapMin;
493   target.fRapMax = fRapMax;
494
495   target.fHistogramsOn = fHistogramsOn;
496
497   for (Int_t i=0; i<2; ++i)
498   {
499     if (fhNClustersITS[i]) target.fhNClustersITS[i] = (TH1F*) fhNClustersITS[i]->Clone();
500     if (fhNClustersTPC[i]) target.fhNClustersTPC[i] = (TH1F*) fhNClustersTPC[i]->Clone();
501
502     if (fhChi2PerClusterITS[i]) target.fhChi2PerClusterITS[i] = (TH1F*) fhChi2PerClusterITS[i]->Clone();
503     if (fhChi2PerClusterTPC[i]) target.fhChi2PerClusterTPC[i] = (TH1F*) fhChi2PerClusterTPC[i]->Clone();
504
505     if (fhC11[i]) target.fhC11[i] = (TH1F*) fhC11[i]->Clone();
506     if (fhC22[i]) target.fhC22[i] = (TH1F*) fhC22[i]->Clone();
507     if (fhC33[i]) target.fhC33[i] = (TH1F*) fhC33[i]->Clone();
508     if (fhC44[i]) target.fhC44[i] = (TH1F*) fhC44[i]->Clone();
509     if (fhC55[i]) target.fhC55[i] = (TH1F*) fhC55[i]->Clone();
510
511     if (fhRel1PtUncertainty[i]) target.fhRel1PtUncertainty[i] = (TH1F*) fhRel1PtUncertainty[i]->Clone();
512
513     if (fhDXY[i]) target.fhDXY[i] = (TH1F*) fhDXY[i]->Clone();
514     if (fhDZ[i]) target.fhDZ[i] = (TH1F*) fhDZ[i]->Clone();
515     if (fhDXYDZ[i]) target.fhDXYDZ[i] = (TH1F*) fhDXYDZ[i]->Clone();
516     if (fhDXYvsDZ[i]) target.fhDXYvsDZ[i] = (TH2F*) fhDXYvsDZ[i]->Clone();
517
518     if (fhDXYNormalized[i]) target.fhDXYNormalized[i] = (TH1F*) fhDXYNormalized[i]->Clone();
519     if (fhDZNormalized[i]) target.fhDZNormalized[i] = (TH1F*) fhDZNormalized[i]->Clone();
520     if (fhDXYvsDZNormalized[i]) target.fhDXYvsDZNormalized[i] = (TH2F*) fhDXYvsDZNormalized[i]->Clone();
521     if (fhNSigmaToVertex[i]) target.fhNSigmaToVertex[i] = (TH1F*) fhNSigmaToVertex[i]->Clone();
522     
523     if (fhPt[i]) target.fhPt[i] = (TH1F*) fhPt[i]->Clone();
524     if (fhEta[i]) target.fhEta[i] = (TH1F*) fhEta[i]->Clone();
525   }
526   if (ffDTheoretical) target.ffDTheoretical = (TF1*) ffDTheoretical->Clone();
527
528   if (fhCutStatistics) target.fhCutStatistics = (TH1F*) fhCutStatistics->Clone();
529   if (fhCutCorrelation) target.fhCutCorrelation = (TH2F*) fhCutCorrelation->Clone();
530
531   TNamed::Copy(c);
532 }
533
534 //_____________________________________________________________________________
535 Long64_t AliESDtrackCuts::Merge(TCollection* list) {
536   // Merge a list of AliESDtrackCuts objects with this (needed for PROOF)
537   // Returns the number of merged objects (including this)
538   if (!list)
539     return 0;
540   if (list->IsEmpty())
541     return 1;
542   if (!fHistogramsOn)
543     return 0;
544   TIterator* iter = list->MakeIterator();
545   TObject* obj;
546
547   // collection of measured and generated histograms
548   Int_t count = 0;
549   while ((obj = iter->Next())) {
550
551     AliESDtrackCuts* entry = dynamic_cast<AliESDtrackCuts*>(obj);
552     if (entry == 0)
553       continue;
554
555     if (!entry->fHistogramsOn)
556       continue;
557
558     for (Int_t i=0; i<2; i++) {
559       
560       fhNClustersITS[i]      ->Add(entry->fhNClustersITS[i]     );      
561       fhNClustersTPC[i]      ->Add(entry->fhNClustersTPC[i]     ); 
562                                                                     
563       fhChi2PerClusterITS[i] ->Add(entry->fhChi2PerClusterITS[i]); 
564       fhChi2PerClusterTPC[i] ->Add(entry->fhChi2PerClusterTPC[i]); 
565                                                                     
566       fhC11[i]               ->Add(entry->fhC11[i]              ); 
567       fhC22[i]               ->Add(entry->fhC22[i]              ); 
568       fhC33[i]               ->Add(entry->fhC33[i]              ); 
569       fhC44[i]               ->Add(entry->fhC44[i]              ); 
570       fhC55[i]               ->Add(entry->fhC55[i]              );
571
572       fhRel1PtUncertainty[i] ->Add(entry->fhRel1PtUncertainty[i]);
573                                                                     
574       fhDXY[i]               ->Add(entry->fhDXY[i]              ); 
575       fhDZ[i]                ->Add(entry->fhDZ[i]               ); 
576       fhDXYDZ[i]             ->Add(entry->fhDXYDZ[i]          );
577       fhDXYvsDZ[i]           ->Add(entry->fhDXYvsDZ[i]          );
578
579       fhDXYNormalized[i]     ->Add(entry->fhDXYNormalized[i]    );
580       fhDZNormalized[i]      ->Add(entry->fhDZNormalized[i]     );
581       fhDXYvsDZNormalized[i] ->Add(entry->fhDXYvsDZNormalized[i]);
582       fhNSigmaToVertex[i]    ->Add(entry->fhNSigmaToVertex[i]); 
583
584       fhPt[i]                ->Add(entry->fhPt[i]); 
585       fhEta[i]               ->Add(entry->fhEta[i]); 
586     }      
587
588     fhCutStatistics  ->Add(entry->fhCutStatistics);        
589     fhCutCorrelation ->Add(entry->fhCutCorrelation);      
590
591     count++;
592   }
593   return count+1;
594 }
595
596 //____________________________________________________________________
597 AliESDtrackCuts* AliESDtrackCuts::GetStandardTPCOnlyTrackCuts()
598 {
599   // creates an AliESDtrackCuts object and fills it with standard (pre data-taking) values for TPC-only cuts
600   
601   Printf("AliESDtrackCuts::GetStandardTPCOnlyTrackCuts: Creating track cuts for TPC-only.");
602   
603   AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
604   
605   esdTrackCuts->SetMinNClustersTPC(50);
606   esdTrackCuts->SetMaxChi2PerClusterTPC(4);
607   esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
608   
609   esdTrackCuts->SetMaxDCAToVertexZ(3.2);
610   esdTrackCuts->SetMaxDCAToVertexXY(2.4);
611   esdTrackCuts->SetDCAToVertex2D(kTRUE);
612   
613   return esdTrackCuts;
614 }
615
616 //____________________________________________________________________
617 AliESDtrackCuts* AliESDtrackCuts::GetStandardITSTPCTrackCuts2009(Bool_t selPrimaries)
618 {
619   // creates an AliESDtrackCuts object and fills it with standard values for ITS-TPC cuts for pp 2009 data
620   
621   Printf("AliESDtrackCuts::GetStandardITSTPCTrackCuts: Creating track cuts for ITS+TPC.");
622   
623   AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
624
625   // TPC  
626   esdTrackCuts->SetRequireTPCStandAlone(kTRUE); // to get chi2 and ncls of kTPCin
627   esdTrackCuts->SetMinNClustersTPC(70);
628   esdTrackCuts->SetMaxChi2PerClusterTPC(4);
629   esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
630   esdTrackCuts->SetRequireTPCRefit(kTRUE);
631   // ITS
632   esdTrackCuts->SetRequireITSRefit(kTRUE);
633   esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
634                                          AliESDtrackCuts::kAny);
635   if(selPrimaries) {
636     // 7*(0.0050+0.0060/pt^0.9)
637     esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0350+0.0420/pt^0.9");
638   }
639   esdTrackCuts->SetMaxDCAToVertexZ(1.e6);
640   esdTrackCuts->SetDCAToVertex2D(kFALSE);
641   esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
642   //esdTrackCuts->SetEtaRange(-0.8,+0.8);
643   
644   return esdTrackCuts;
645 }
646 //____________________________________________________________________
647 AliESDtrackCuts* AliESDtrackCuts::GetStandardITSPureSATrackCuts2009(Bool_t selPrimaries, Bool_t useForPid)
648 {
649   // creates an AliESDtrackCuts object and fills it with standard values for ITS pure SA tracks
650   
651   AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
652   esdTrackCuts->SetRequireITSPureStandAlone(kTRUE);
653   esdTrackCuts->SetRequireITSRefit(kTRUE); 
654   esdTrackCuts->SetMinNClustersITS(4);
655   esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
656                                          AliESDtrackCuts::kAny);
657   esdTrackCuts->SetMaxChi2PerClusterITS(1.);
658
659   if(selPrimaries) {
660     // 7*(0.0085+0.0026/pt^1.55)
661     esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0595+0.0182/pt^1.55");
662   }
663   if(useForPid){
664     esdTrackCuts->SetRequireITSPid(kTRUE);
665   }
666   return esdTrackCuts;
667 }
668 //____________________________________________________________________
669 AliESDtrackCuts* AliESDtrackCuts::GetStandardITSSATrackCuts2009(Bool_t selPrimaries, Bool_t useForPid)
670 {
671   // creates an AliESDtrackCuts object and fills it with standard values for ITS pure SA tracks
672   
673   AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts;
674   esdTrackCuts->SetRequireITSStandAlone(kTRUE);
675   esdTrackCuts->SetRequireITSPureStandAlone(kFALSE);
676   esdTrackCuts->SetRequireITSRefit(kTRUE); 
677   esdTrackCuts->SetMinNClustersITS(4);
678   esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
679                                          AliESDtrackCuts::kAny);
680   esdTrackCuts->SetMaxChi2PerClusterITS(1.);
681
682   if(selPrimaries) {
683     // 7*(0.0085+0.0026/pt^1.55)
684     esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0595+0.0182/pt^1.55");
685   }
686   if(useForPid){
687     esdTrackCuts->SetRequireITSPid(kTRUE);
688   }
689   return esdTrackCuts;
690 }
691
692
693 //____________________________________________________________________
694 Int_t AliESDtrackCuts::GetReferenceMultiplicity(AliESDEvent* esd, Bool_t tpcOnly)
695 {
696   // Gets reference multiplicity following the standard cuts and a defined fiducial volume
697   // tpcOnly = kTRUE -> consider TPC-only tracks
698   //         = kFALSE -> consider global tracks
699   
700   if (!tpcOnly)
701   {
702     Printf("AliESDtrackCuts::GetReferenceMultiplicity: Not implemented for global tracks!");
703     return -1;
704   }
705   
706   AliESDtrackCuts* esdTrackCuts = GetStandardTPCOnlyTrackCuts();
707   esdTrackCuts->SetEtaRange(-0.8, 0.8);
708   esdTrackCuts->SetPtRange(0.15);
709   
710   Int_t nTracks = esdTrackCuts->CountAcceptedTracks(esd);
711   
712   delete esdTrackCuts;
713   esdTrackCuts = 0;
714   
715   return nTracks;
716 }
717
718 //____________________________________________________________________
719 Float_t AliESDtrackCuts::GetSigmaToVertex(AliESDtrack* esdTrack)
720 {
721   // Calculates the number of sigma to the vertex.
722
723   Float_t b[2];
724   Float_t bRes[2];
725   Float_t bCov[3];
726   esdTrack->GetImpactParameters(b,bCov);
727   
728   if (bCov[0]<=0 || bCov[2]<=0) {
729     AliDebugClass(1, "Estimated b resolution lower or equal zero!");
730     bCov[0]=0; bCov[2]=0;
731   }
732   bRes[0] = TMath::Sqrt(bCov[0]);
733   bRes[1] = TMath::Sqrt(bCov[2]);
734
735   // -----------------------------------
736   // How to get to a n-sigma cut?
737   //
738   // The accumulated statistics from 0 to d is
739   //
740   // ->  Erf(d/Sqrt(2)) for a 1-dim gauss (d = n_sigma)
741   // ->  1 - Exp(-d**2) for a 2-dim gauss (d*d = dx*dx + dy*dy != n_sigma)
742   //
743   // It means that for a 2-dim gauss: n_sigma(d) = Sqrt(2)*ErfInv(1 - Exp((-d**2)/2)
744   // Can this be expressed in a different way?
745
746   if (bRes[0] == 0 || bRes[1] ==0)
747     return -1;
748
749   Float_t d = TMath::Sqrt(TMath::Power(b[0]/bRes[0],2) + TMath::Power(b[1]/bRes[1],2));
750
751   // work around precision problem
752   // if d is too big, TMath::Exp(...) gets 0, and TMath::ErfInverse(1) that should be infinite, gets 0 :(
753   // 1e-15 corresponds to nsigma ~ 7.7
754   if (TMath::Exp(-d * d / 2) < 1e-15)
755     return 1000;
756
757   Float_t nSigma = TMath::ErfInverse(1 - TMath::Exp(-d * d / 2)) * TMath::Sqrt(2);
758   return nSigma;
759 }
760
761 void AliESDtrackCuts::EnableNeededBranches(TTree* tree)
762 {
763   // enables the branches needed by AcceptTrack, for a list see comment of AcceptTrack
764
765   tree->SetBranchStatus("fTracks.fFlags", 1);
766   tree->SetBranchStatus("fTracks.fITSncls", 1);
767   tree->SetBranchStatus("fTracks.fTPCncls", 1);
768   tree->SetBranchStatus("fTracks.fITSchi2", 1);
769   tree->SetBranchStatus("fTracks.fTPCchi2", 1);
770   tree->SetBranchStatus("fTracks.fC*", 1);
771   tree->SetBranchStatus("fTracks.fD", 1);
772   tree->SetBranchStatus("fTracks.fZ", 1);
773   tree->SetBranchStatus("fTracks.fCdd", 1);
774   tree->SetBranchStatus("fTracks.fCdz", 1);
775   tree->SetBranchStatus("fTracks.fCzz", 1);
776   tree->SetBranchStatus("fTracks.fP*", 1);
777   tree->SetBranchStatus("fTracks.fR*", 1);
778   tree->SetBranchStatus("fTracks.fKinkIndexes*", 1);
779 }
780
781 //____________________________________________________________________
782 Bool_t AliESDtrackCuts::AcceptTrack(AliESDtrack* esdTrack) 
783 {
784   // 
785   // figure out if the tracks survives all the track cuts defined
786   //
787   // the different quality parameter and kinematic values are first
788   // retrieved from the track. then it is found out what cuts the
789   // track did not survive and finally the cuts are imposed.
790
791   // this function needs the following branches:
792   // fTracks.fFlags
793   // fTracks.fITSncls
794   // fTracks.fTPCncls
795   // fTracks.fITSchi2
796   // fTracks.fTPCchi2
797   // fTracks.fC   //GetExternalCovariance
798   // fTracks.fD   //GetImpactParameters
799   // fTracks.fZ   //GetImpactParameters
800   // fTracks.fCdd //GetImpactParameters
801   // fTracks.fCdz //GetImpactParameters
802   // fTracks.fCzz //GetImpactParameters
803   // fTracks.fP   //GetPxPyPz
804   // fTracks.fR   //GetMass
805   // fTracks.fP   //GetMass
806   // fTracks.fKinkIndexes
807
808   UInt_t status = esdTrack->GetStatus();
809
810   // getting quality parameters from the ESD track
811   Int_t nClustersITS = esdTrack->GetITSclusters(0);
812   Int_t nClustersTPC = -1;
813   if(fCutRequireTPCStandAlone) {
814     nClustersTPC = esdTrack->GetTPCNclsIter1();
815   }
816   else {
817     nClustersTPC = esdTrack->GetTPCclusters(0);
818   }
819
820   Int_t nClustersTPCShared = esdTrack->GetTPCnclsS();
821   Float_t fracClustersTPCShared = -1.;
822
823   Float_t chi2PerClusterITS = -1;
824   Float_t chi2PerClusterTPC = -1;
825   if (nClustersITS!=0)
826     chi2PerClusterITS = esdTrack->GetITSchi2()/Float_t(nClustersITS);
827   if (nClustersTPC!=0) {
828     if(fCutRequireTPCStandAlone) {
829       chi2PerClusterTPC = esdTrack->GetTPCchi2Iter1()/Float_t(nClustersTPC);
830     } else {
831       chi2PerClusterTPC = esdTrack->GetTPCchi2()/Float_t(nClustersTPC);
832     }
833     fracClustersTPCShared = Float_t(nClustersTPCShared)/Float_t(nClustersTPC);
834   }
835
836   Double_t extCov[15];
837   esdTrack->GetExternalCovariance(extCov);
838
839   // getting the track to vertex parameters
840   Float_t nSigmaToVertex = GetSigmaToVertex(esdTrack);
841       
842   Float_t b[2];
843   Float_t bCov[3];
844   esdTrack->GetImpactParameters(b,bCov);
845   if (bCov[0]<=0 || bCov[2]<=0) {
846     AliDebug(1, "Estimated b resolution lower or equal zero!");
847     bCov[0]=0; bCov[2]=0;
848   }
849
850
851   // set pt-dependent DCA cuts, if requested
852   SetPtDepDCACuts(esdTrack->Pt());
853
854
855   Float_t dcaToVertexXY = b[0];
856   Float_t dcaToVertexZ = b[1];
857
858   Float_t dcaToVertex = -1;
859   
860   if (fCutDCAToVertex2D)
861   {
862     dcaToVertex = TMath::Sqrt(dcaToVertexXY*dcaToVertexXY/fCutMaxDCAToVertexXY/fCutMaxDCAToVertexXY + dcaToVertexZ*dcaToVertexZ/fCutMaxDCAToVertexZ/fCutMaxDCAToVertexZ);
863   }
864   else
865     dcaToVertex = TMath::Sqrt(dcaToVertexXY*dcaToVertexXY + dcaToVertexZ*dcaToVertexZ);
866     
867   // getting the kinematic variables of the track
868   // (assuming the mass is known)
869   Double_t p[3];
870   esdTrack->GetPxPyPz(p);
871
872   Float_t momentum = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2) + TMath::Power(p[2],2));
873   Float_t pt       = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2));
874   Float_t energy   = TMath::Sqrt(TMath::Power(esdTrack->GetMass(),2) + TMath::Power(momentum,2));
875
876   //y-eta related calculations
877   Float_t eta = -100.;
878   Float_t y   = -100.;
879   if((momentum != TMath::Abs(p[2]))&&(momentum != 0))
880     eta = 0.5*TMath::Log((momentum + p[2])/(momentum - p[2]));
881   if((energy != TMath::Abs(p[2]))&&(momentum != 0))
882     y = 0.5*TMath::Log((energy + p[2])/(energy - p[2]));
883     
884   if (extCov[14] < 0) 
885   {
886     Printf("AliESDtrackCuts::AcceptTrack: WARNING: GetSigma1Pt2() returns negative value for external covariance matrix element fC[14]: %f. Corrupted track information, track will not be accepted!", extCov[14]);
887     return kFALSE;
888   }
889   Float_t relUncertainty1Pt = TMath::Sqrt(extCov[14])*pt;
890   
891   //########################################################################
892   // cut the track?
893   
894   Bool_t cuts[kNCuts];
895   for (Int_t i=0; i<kNCuts; i++) cuts[i]=kFALSE;
896   
897   // track quality cuts
898   if (fCutRequireTPCRefit && (status&AliESDtrack::kTPCrefit)==0)
899     cuts[0]=kTRUE;
900   if (fCutRequireTPCStandAlone && (status&AliESDtrack::kTPCin)==0)
901     cuts[1]=kTRUE;
902   if (fCutRequireITSRefit && (status&AliESDtrack::kITSrefit)==0)
903     cuts[2]=kTRUE;
904   if (nClustersTPC<fCutMinNClusterTPC)
905     cuts[3]=kTRUE;
906   if (nClustersITS<fCutMinNClusterITS) 
907     cuts[4]=kTRUE;
908   if (chi2PerClusterTPC>fCutMaxChi2PerClusterTPC) 
909     cuts[5]=kTRUE; 
910   if (chi2PerClusterITS>fCutMaxChi2PerClusterITS) 
911     cuts[6]=kTRUE;
912   if (extCov[0]  > fCutMaxC11) 
913     cuts[7]=kTRUE;  
914   if (extCov[2]  > fCutMaxC22) 
915     cuts[8]=kTRUE;  
916   if (extCov[5]  > fCutMaxC33) 
917     cuts[9]=kTRUE;  
918   if (extCov[9]  > fCutMaxC44) 
919     cuts[10]=kTRUE;  
920   if (extCov[14]  > fCutMaxC55) 
921     cuts[11]=kTRUE;  
922   if (nSigmaToVertex > fCutNsigmaToVertex && fCutSigmaToVertexRequired)
923     cuts[12] = kTRUE;
924   // if n sigma could not be calculated
925   if (nSigmaToVertex<0 && fCutSigmaToVertexRequired)
926     cuts[13]=kTRUE;
927   if (!fCutAcceptKinkDaughters && esdTrack->GetKinkIndex(0)>0)
928     cuts[14]=kTRUE;
929   // track kinematics cut
930   if((momentum < fPMin) || (momentum > fPMax)) 
931     cuts[15]=kTRUE;
932   if((pt < fPtMin) || (pt > fPtMax)) 
933     cuts[16] = kTRUE;
934   if((p[0] < fPxMin) || (p[0] > fPxMax)) 
935     cuts[17] = kTRUE;
936   if((p[1] < fPyMin) || (p[1] > fPyMax)) 
937     cuts[18] = kTRUE;
938   if((p[2] < fPzMin) || (p[2] > fPzMax))
939     cuts[19] = kTRUE;
940   if((eta < fEtaMin) || (eta > fEtaMax))
941     cuts[20] = kTRUE;
942   if((y < fRapMin) || (y > fRapMax)) 
943     cuts[21] = kTRUE;
944   if (fCutDCAToVertex2D && dcaToVertex > 1)
945     cuts[22] = kTRUE;
946   if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexXY) > fCutMaxDCAToVertexXY)
947     cuts[23] = kTRUE;
948   if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexZ) > fCutMaxDCAToVertexZ)
949     cuts[24] = kTRUE;
950   if (fCutDCAToVertex2D && fCutMinDCAToVertexXY > 0 && fCutMinDCAToVertexZ > 0 && dcaToVertexXY*dcaToVertexXY/fCutMinDCAToVertexXY/fCutMinDCAToVertexXY + dcaToVertexZ*dcaToVertexZ/fCutMinDCAToVertexZ/fCutMinDCAToVertexZ < 1)
951     cuts[25] = kTRUE;
952   if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexXY) < fCutMinDCAToVertexXY)
953     cuts[26] = kTRUE;
954   if (!fCutDCAToVertex2D && TMath::Abs(dcaToVertexZ) < fCutMinDCAToVertexZ)
955     cuts[27] = kTRUE;
956   
957   for (Int_t i = 0; i < 3; i++)
958     cuts[28+i] = !CheckITSClusterRequirement(fCutClusterRequirementITS[i], esdTrack->HasPointOnITSLayer(i*2), esdTrack->HasPointOnITSLayer(i*2+1));
959   
960   if(fCutRequireITSStandAlone || fCutRequireITSpureSA){
961     if ((status & AliESDtrack::kITSin) == 0 || (status & AliESDtrack::kTPCin)){
962       // TPC tracks
963       cuts[31] = kTRUE; 
964     }else{
965       // ITS standalone tracks
966       if(fCutRequireITSStandAlone && !fCutRequireITSpureSA){
967         if(status & AliESDtrack::kITSpureSA) cuts[31] = kTRUE;
968       }else if(fCutRequireITSpureSA){
969         if(!(status & AliESDtrack::kITSpureSA)) cuts[31] = kTRUE;
970       }
971     }
972   }
973
974   if (relUncertainty1Pt > fCutMaxRel1PtUncertainty)
975      cuts[32] = kTRUE;
976
977   if (!fCutAcceptSharedTPCClusters && nClustersTPCShared!=0)
978     cuts[33] = kTRUE;
979
980   if (fracClustersTPCShared > fCutMaxFractionSharedTPCClusters)
981     cuts[34] = kTRUE;  
982
983   if(fCutRequireITSPid){
984     UChar_t clumap=esdTrack->GetITSClusterMap();
985     Int_t nPointsForPid=0;
986     for(Int_t i=2; i<6; i++){
987       if(clumap&(1<<i)) ++nPointsForPid;
988     }
989     if(nPointsForPid<3) cuts[35] = kTRUE;
990   }
991
992   Bool_t cut=kFALSE;
993   for (Int_t i=0; i<kNCuts; i++) 
994     if (cuts[i]) {cut = kTRUE;}
995
996
997   //########################################################################
998   // filling histograms
999   if (fHistogramsOn) {
1000     fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n tracks")));
1001     if (cut)
1002       fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n cut tracks")));
1003     
1004     for (Int_t i=0; i<kNCuts; i++) {
1005       if (cuts[i])
1006         fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin(fgkCutNames[i])));
1007
1008       for (Int_t j=i; j<kNCuts; j++) {
1009         if (cuts[i] && cuts[j]) {
1010           Float_t xC = fhCutCorrelation->GetXaxis()->GetBinCenter(fhCutCorrelation->GetXaxis()->FindBin(fgkCutNames[i]));
1011           Float_t yC = fhCutCorrelation->GetYaxis()->GetBinCenter(fhCutCorrelation->GetYaxis()->FindBin(fgkCutNames[j]));
1012           fhCutCorrelation->Fill(xC, yC);
1013         }
1014       }
1015     }
1016   }
1017
1018   // now we loop over the filling of the histograms twice: once "before" the cut, once "after"
1019   // the code is not in a function due to too many local variables that would need to be passed
1020
1021   for (Int_t id = 0; id < 2; id++)
1022   {
1023     // id = 0 --> before cut
1024     // id = 1 --> after cut
1025
1026     if (fHistogramsOn)
1027     {
1028       fhNClustersITS[id]->Fill(nClustersITS);
1029       fhNClustersTPC[id]->Fill(nClustersTPC);
1030       fhChi2PerClusterITS[id]->Fill(chi2PerClusterITS);
1031       fhChi2PerClusterTPC[id]->Fill(chi2PerClusterTPC);
1032
1033       fhC11[id]->Fill(extCov[0]);
1034       fhC22[id]->Fill(extCov[2]);
1035       fhC33[id]->Fill(extCov[5]);
1036       fhC44[id]->Fill(extCov[9]);
1037       fhC55[id]->Fill(extCov[14]);
1038
1039       fhRel1PtUncertainty[id]->Fill(relUncertainty1Pt);
1040
1041       fhPt[id]->Fill(pt);
1042       fhEta[id]->Fill(eta);
1043
1044       Float_t bRes[2];
1045       bRes[0] = TMath::Sqrt(bCov[0]);
1046       bRes[1] = TMath::Sqrt(bCov[2]);
1047
1048       fhDZ[id]->Fill(b[1]);
1049       fhDXY[id]->Fill(b[0]);
1050       fhDXYDZ[id]->Fill(dcaToVertex);
1051       fhDXYvsDZ[id]->Fill(b[1],b[0]);
1052
1053       if (bRes[0]!=0 && bRes[1]!=0) {
1054         fhDZNormalized[id]->Fill(b[1]/bRes[1]);
1055         fhDXYNormalized[id]->Fill(b[0]/bRes[0]);
1056         fhDXYvsDZNormalized[id]->Fill(b[1]/bRes[1], b[0]/bRes[0]);
1057         fhNSigmaToVertex[id]->Fill(nSigmaToVertex);
1058       }
1059     }
1060
1061     // cut the track
1062     if (cut)
1063       return kFALSE;
1064   }
1065
1066   return kTRUE;
1067 }
1068
1069 //____________________________________________________________________
1070 Bool_t AliESDtrackCuts::CheckITSClusterRequirement(ITSClusterRequirement req, Bool_t clusterL1, Bool_t clusterL2)
1071 {
1072   // checks if the cluster requirement is fullfilled (in this case: return kTRUE)
1073   
1074   switch (req)
1075   {
1076         case kOff:        return kTRUE;
1077         case kNone:       return !clusterL1 && !clusterL2;
1078         case kAny:        return clusterL1 || clusterL2;
1079         case kFirst:      return clusterL1;
1080         case kOnlyFirst:  return clusterL1 && !clusterL2;
1081         case kSecond:     return clusterL2;
1082         case kOnlySecond: return clusterL2 && !clusterL1;
1083         case kBoth:       return clusterL1 && clusterL2;
1084   }
1085   
1086   return kFALSE;
1087 }
1088
1089 //____________________________________________________________________
1090 AliESDtrack* AliESDtrackCuts::GetTPCOnlyTrack(AliESDEvent* esd, Int_t iTrack)
1091 {
1092   
1093   // Utility function to 
1094   // create a TPC only track from the given esd track
1095   // 
1096   // IMPORTANT: The track has to be deleted by the user
1097   //
1098   // NB. most of the functionality to get a TPC only track from an ESD track is in AliESDtrack, where it should be
1099   // there are only missing propagations here that are needed for old data
1100   // this function will therefore become obsolete
1101   //
1102   // adapted from code provided by CKB
1103
1104   if (!esd->GetPrimaryVertexTPC())
1105     return 0; // No TPC vertex no TPC tracks
1106
1107   if(!esd->GetPrimaryVertexTPC()->GetStatus())
1108     return 0; // TPC Vertex is created by default in AliESDEvent, do not use in this case
1109
1110   AliESDtrack* track = esd->GetTrack(iTrack);
1111   if (!track)
1112     return 0;
1113
1114   AliESDtrack *tpcTrack = new AliESDtrack();
1115
1116   // only true if we have a tpc track
1117   if (!track->FillTPCOnlyTrack(*tpcTrack))
1118   {
1119     delete tpcTrack;
1120     return 0;
1121   }
1122
1123   // propagate to Vertex
1124   // not needed for normal reconstructed ESDs...
1125   // Double_t pTPC[2],covTPC[3];
1126   // tpcTrack->PropagateToDCA(esd->GetPrimaryVertexTPC(), esd->GetMagneticField(), 10000,  pTPC, covTPC);
1127
1128   return tpcTrack;
1129 }
1130
1131 //____________________________________________________________________
1132 TObjArray* AliESDtrackCuts::GetAcceptedTracks(AliESDEvent* esd,Bool_t bTPC)
1133 {
1134   //
1135   // returns an array of all tracks that pass the cuts
1136   // or an array of TPC only tracks (propagated to the TPC vertex during reco)
1137   // tracks that pass the cut
1138
1139   TObjArray* acceptedTracks = new TObjArray();
1140
1141   // loop over esd tracks
1142   for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
1143     if(bTPC){
1144       if(!esd->GetPrimaryVertexTPC())return acceptedTracks; // No TPC vertex no TPC tracks
1145       if(!esd->GetPrimaryVertexTPC()->GetStatus())return acceptedTracks; // No proper TPC vertex, only the default
1146
1147       AliESDtrack *tpcTrack = GetTPCOnlyTrack(esd, iTrack);
1148       if (!tpcTrack)
1149         continue;
1150
1151       if (AcceptTrack(tpcTrack)) {
1152         acceptedTracks->Add(tpcTrack);
1153       }
1154       else
1155         delete tpcTrack;
1156     }
1157     else
1158     {
1159       AliESDtrack* track = esd->GetTrack(iTrack);
1160       if(AcceptTrack(track))
1161         acceptedTracks->Add(track);
1162     }
1163   } 
1164   if(bTPC)acceptedTracks->SetOwner(kTRUE);
1165   return acceptedTracks;
1166 }
1167
1168 //____________________________________________________________________
1169 Int_t AliESDtrackCuts::CountAcceptedTracks(AliESDEvent* esd)
1170 {
1171   //
1172   // returns an the number of tracks that pass the cuts
1173   //
1174
1175   Int_t count = 0;
1176
1177   // loop over esd tracks
1178   for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
1179     AliESDtrack* track = esd->GetTrack(iTrack);
1180     if (AcceptTrack(track))
1181       count++;
1182   }
1183
1184   return count;
1185 }
1186
1187 //____________________________________________________________________
1188  void AliESDtrackCuts::DefineHistograms(Int_t color) {
1189    // 
1190    // diagnostics histograms are defined
1191    // 
1192
1193    fHistogramsOn=kTRUE;
1194    
1195    Bool_t oldStatus = TH1::AddDirectoryStatus();
1196    TH1::AddDirectory(kFALSE);
1197    
1198    //###################################################################################
1199    // defining histograms
1200
1201    fhCutStatistics = new TH1F("cut_statistics","cut statistics",kNCuts+4,-0.5,kNCuts+3.5);
1202
1203    fhCutStatistics->GetXaxis()->SetBinLabel(1,"n tracks");
1204    fhCutStatistics->GetXaxis()->SetBinLabel(2,"n cut tracks");
1205
1206    fhCutCorrelation = new TH2F("cut_correlation","cut correlation",kNCuts,-0.5,kNCuts-0.5,kNCuts,-0.5,kNCuts-0.5);;
1207   
1208    for (Int_t i=0; i<kNCuts; i++) {
1209      fhCutStatistics->GetXaxis()->SetBinLabel(i+4,fgkCutNames[i]);
1210      fhCutCorrelation->GetXaxis()->SetBinLabel(i+1,fgkCutNames[i]);
1211      fhCutCorrelation->GetYaxis()->SetBinLabel(i+1,fgkCutNames[i]);
1212    }
1213
1214   fhCutStatistics  ->SetLineColor(color);
1215   fhCutCorrelation ->SetLineColor(color);
1216   fhCutStatistics  ->SetLineWidth(2);
1217   fhCutCorrelation ->SetLineWidth(2);
1218
1219   for (Int_t i=0; i<2; i++) {
1220     fhNClustersITS[i]        = new TH1F("nClustersITS"    ,"",8,-0.5,7.5);
1221     fhNClustersTPC[i]        = new TH1F("nClustersTPC"     ,"",165,-0.5,164.5);
1222     fhChi2PerClusterITS[i]   = new TH1F("chi2PerClusterITS","",500,0,10);
1223     fhChi2PerClusterTPC[i]   = new TH1F("chi2PerClusterTPC","",500,0,10);
1224
1225     fhC11[i]                 = new TH1F("covMatrixDiagonal11","",2000,0,20);
1226     fhC22[i]                 = new TH1F("covMatrixDiagonal22","",2000,0,20);
1227     fhC33[i]                 = new TH1F("covMatrixDiagonal33","",1000,0,0.1);
1228     fhC44[i]                 = new TH1F("covMatrixDiagonal44","",1000,0,0.1);
1229     fhC55[i]                 = new TH1F("covMatrixDiagonal55","",1000,0,5);
1230
1231     fhRel1PtUncertainty[i]   = new TH1F("rel1PtUncertainty","",1000,0,5);
1232
1233     fhDXY[i]                 = new TH1F("dXY"    ,"",500,-10,10);
1234     fhDZ[i]                  = new TH1F("dZ"     ,"",500,-10,10);
1235     fhDXYDZ[i]               = new TH1F("dXYDZ"  ,"",500,0,10);
1236     fhDXYvsDZ[i]             = new TH2F("dXYvsDZ","",200,-10,10,200,-10,10);
1237
1238     fhDXYNormalized[i]       = new TH1F("dXYNormalized"    ,"",500,-10,10);
1239     fhDZNormalized[i]        = new TH1F("dZNormalized"     ,"",500,-10,10);
1240     fhDXYvsDZNormalized[i]   = new TH2F("dXYvsDZNormalized","",200,-10,10,200,-10,10);
1241
1242     fhNSigmaToVertex[i]      = new TH1F("nSigmaToVertex","",500,0,10);
1243
1244     fhPt[i]                  = new TH1F("pt"     ,"p_{T} distribution;p_{T} (GeV/c)", 800, 0.0, 10.0);
1245     fhEta[i]                 = new TH1F("eta"     ,"#eta distribution;#eta",40,-2.0,2.0);
1246     
1247     fhNClustersITS[i]->SetTitle("n ITS clusters");
1248     fhNClustersTPC[i]->SetTitle("n TPC clusters");
1249     fhChi2PerClusterITS[i]->SetTitle("#Chi^{2} per ITS cluster");
1250     fhChi2PerClusterTPC[i]->SetTitle("#Chi^{2} per TPC cluster");
1251
1252     fhC11[i]->SetTitle("cov 11 : #sigma_{y}^{2} [cm^{2}]");
1253     fhC22[i]->SetTitle("cov 22 : #sigma_{z}^{2} [cm^{2}]");
1254     fhC33[i]->SetTitle("cov 33 : #sigma_{sin(#phi)}^{2}");
1255     fhC44[i]->SetTitle("cov 44 : #sigma_{tan(#theta_{dip})}^{2}");
1256     fhC55[i]->SetTitle("cov 55 : #sigma_{1/p_{T}}^{2} [(c/GeV)^2]");
1257
1258     fhRel1PtUncertainty[i]->SetTitle("rel. uncertainty of 1/p_{T}");
1259
1260     fhDXY[i]->SetXTitle("transverse impact parameter (cm)");
1261     fhDZ[i]->SetXTitle("longitudinal impact parameter (cm)");
1262     fhDXYDZ[i]->SetTitle("absolute impact parameter;sqrt(dXY**2 + dZ**2)  (cm)");
1263     fhDXYvsDZ[i]->SetXTitle("longitudinal impact parameter (cm)");
1264     fhDXYvsDZ[i]->SetYTitle("transverse impact parameter (cm)");
1265
1266     fhDXYNormalized[i]->SetTitle("normalized trans impact par (n#sigma)");
1267     fhDZNormalized[i]->SetTitle("normalized long impact par (n#sigma)");
1268     fhDXYvsDZNormalized[i]->SetTitle("normalized long impact par (n#sigma)");
1269     fhDXYvsDZNormalized[i]->SetYTitle("normalized trans impact par (n#sigma)");
1270     fhNSigmaToVertex[i]->SetTitle("n #sigma to vertex");
1271
1272     fhNClustersITS[i]->SetLineColor(color);   fhNClustersITS[i]->SetLineWidth(2);
1273     fhNClustersTPC[i]->SetLineColor(color);   fhNClustersTPC[i]->SetLineWidth(2);
1274     fhChi2PerClusterITS[i]->SetLineColor(color);   fhChi2PerClusterITS[i]->SetLineWidth(2);
1275     fhChi2PerClusterTPC[i]->SetLineColor(color);   fhChi2PerClusterTPC[i]->SetLineWidth(2);
1276
1277     fhC11[i]->SetLineColor(color);   fhC11[i]->SetLineWidth(2);
1278     fhC22[i]->SetLineColor(color);   fhC22[i]->SetLineWidth(2);
1279     fhC33[i]->SetLineColor(color);   fhC33[i]->SetLineWidth(2);
1280     fhC44[i]->SetLineColor(color);   fhC44[i]->SetLineWidth(2);
1281     fhC55[i]->SetLineColor(color);   fhC55[i]->SetLineWidth(2);
1282
1283     fhRel1PtUncertainty[i]->SetLineColor(color); fhRel1PtUncertainty[i]->SetLineWidth(2);
1284
1285     fhDXY[i]->SetLineColor(color);   fhDXY[i]->SetLineWidth(2);
1286     fhDZ[i]->SetLineColor(color);    fhDZ[i]->SetLineWidth(2);
1287     fhDXYDZ[i]->SetLineColor(color); fhDXYDZ[i]->SetLineWidth(2);
1288
1289     fhDXYNormalized[i]->SetLineColor(color);   fhDXYNormalized[i]->SetLineWidth(2);
1290     fhDZNormalized[i]->SetLineColor(color);    fhDZNormalized[i]->SetLineWidth(2);
1291     fhNSigmaToVertex[i]->SetLineColor(color);  fhNSigmaToVertex[i]->SetLineWidth(2); 
1292   }
1293
1294   // The number of sigmas to the vertex is per definition gaussian
1295   ffDTheoretical = new TF1("nSigmaToVertexTheoretical","([0]/2.506628274)*exp(-(x**2)/2)",0,50);
1296   ffDTheoretical->SetParameter(0,1);
1297
1298   TH1::AddDirectory(oldStatus);
1299 }
1300
1301 //____________________________________________________________________
1302 Bool_t AliESDtrackCuts::LoadHistograms(const Char_t* dir)
1303 {
1304   //
1305   // loads the histograms from a file
1306   // if dir is empty a directory with the name of this object is taken (like in SaveHistogram)
1307   //
1308
1309   if (!dir)
1310     dir = GetName();
1311
1312   if (!gDirectory->cd(dir))
1313     return kFALSE;
1314
1315   ffDTheoretical = dynamic_cast<TF1*> (gDirectory->Get("nSigmaToVertexTheory"));
1316
1317   fhCutStatistics = dynamic_cast<TH1F*> (gDirectory->Get("cut_statistics"));
1318   fhCutCorrelation = dynamic_cast<TH2F*> (gDirectory->Get("cut_correlation"));
1319
1320   for (Int_t i=0; i<2; i++) {
1321     if (i==0)
1322     {
1323       gDirectory->cd("before_cuts");
1324     }
1325     else
1326       gDirectory->cd("after_cuts");
1327
1328     fhNClustersITS[i]      = dynamic_cast<TH1F*> (gDirectory->Get("nClustersITS"     ));
1329     fhNClustersTPC[i]      = dynamic_cast<TH1F*> (gDirectory->Get("nClustersTPC"     ));
1330     fhChi2PerClusterITS[i] = dynamic_cast<TH1F*> (gDirectory->Get("chi2PerClusterITS"));
1331     fhChi2PerClusterTPC[i] = dynamic_cast<TH1F*> (gDirectory->Get("chi2PerClusterTPC"));
1332
1333     fhC11[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal11"));
1334     fhC22[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal22"));
1335     fhC33[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal33"));
1336     fhC44[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal44"));
1337     fhC55[i] = dynamic_cast<TH1F*> (gDirectory->Get("covMatrixDiagonal55"));
1338
1339     fhRel1PtUncertainty[i] = dynamic_cast<TH1F*> (gDirectory->Get("rel1PtUncertainty"));
1340
1341     fhDXY[i] =     dynamic_cast<TH1F*> (gDirectory->Get("dXY"    ));
1342     fhDZ[i] =      dynamic_cast<TH1F*> (gDirectory->Get("dZ"     ));
1343     fhDXYDZ[i] =   dynamic_cast<TH1F*> (gDirectory->Get("dXYDZ"));
1344     fhDXYvsDZ[i] = dynamic_cast<TH2F*> (gDirectory->Get("dXYvsDZ"));
1345
1346     fhDXYNormalized[i] =     dynamic_cast<TH1F*> (gDirectory->Get("dXYNormalized"    ));
1347     fhDZNormalized[i] =      dynamic_cast<TH1F*> (gDirectory->Get("dZNormalized"     ));
1348     fhDXYvsDZNormalized[i] = dynamic_cast<TH2F*> (gDirectory->Get("dXYvsDZNormalized"));
1349     fhNSigmaToVertex[i] = dynamic_cast<TH1F*> (gDirectory->Get("nSigmaToVertex"));
1350
1351     fhPt[i] = dynamic_cast<TH1F*> (gDirectory->Get("pt"));
1352     fhEta[i] = dynamic_cast<TH1F*> (gDirectory->Get("eta"));
1353
1354     gDirectory->cd("../");
1355   }
1356
1357   gDirectory->cd("..");
1358
1359   return kTRUE;
1360 }
1361
1362 //____________________________________________________________________
1363 void AliESDtrackCuts::SaveHistograms(const Char_t* dir) {
1364   //
1365   // saves the histograms in a directory (dir)
1366   //
1367
1368   if (!fHistogramsOn) {
1369     AliDebug(0, "Histograms not on - cannot save histograms!!!");
1370     return;
1371   }
1372
1373   if (!dir)
1374     dir = GetName();
1375
1376   gDirectory->mkdir(dir);
1377   gDirectory->cd(dir);
1378
1379   gDirectory->mkdir("before_cuts");
1380   gDirectory->mkdir("after_cuts");
1381  
1382   // a factor of 2 is needed since n sigma is positive
1383   ffDTheoretical->SetParameter(0,2*fhNSigmaToVertex[0]->Integral("width"));
1384   ffDTheoretical->Write("nSigmaToVertexTheory");
1385
1386   fhCutStatistics->Write();
1387   fhCutCorrelation->Write();
1388
1389   for (Int_t i=0; i<2; i++) {
1390     if (i==0)
1391       gDirectory->cd("before_cuts");
1392     else
1393       gDirectory->cd("after_cuts");
1394
1395     fhNClustersITS[i]        ->Write();
1396     fhNClustersTPC[i]        ->Write();
1397     fhChi2PerClusterITS[i]   ->Write();
1398     fhChi2PerClusterTPC[i]   ->Write();
1399
1400     fhC11[i]                 ->Write();
1401     fhC22[i]                 ->Write();
1402     fhC33[i]                 ->Write();
1403     fhC44[i]                 ->Write();
1404     fhC55[i]                 ->Write();
1405
1406     fhRel1PtUncertainty[i]   ->Write();
1407
1408     fhDXY[i]                 ->Write();
1409     fhDZ[i]                  ->Write();
1410     fhDXYDZ[i]               ->Write();
1411     fhDXYvsDZ[i]             ->Write();
1412
1413     fhDXYNormalized[i]       ->Write();
1414     fhDZNormalized[i]        ->Write();
1415     fhDXYvsDZNormalized[i]   ->Write();
1416     fhNSigmaToVertex[i]      ->Write();
1417
1418     fhPt[i]                  ->Write();
1419     fhEta[i]                 ->Write();
1420     
1421     gDirectory->cd("../");
1422   }
1423
1424   gDirectory->cd("../");
1425 }
1426
1427 //____________________________________________________________________
1428 void AliESDtrackCuts::DrawHistograms()
1429 {
1430   // draws some histograms
1431
1432   TCanvas* canvas1 = new TCanvas(Form("%s_1", GetName()), "Track Quality Results1", 800, 800);
1433   canvas1->Divide(2, 2);
1434
1435   canvas1->cd(1);
1436   fhNClustersTPC[0]->SetStats(kFALSE);
1437   fhNClustersTPC[0]->Draw();
1438
1439   canvas1->cd(2);
1440   fhChi2PerClusterTPC[0]->SetStats(kFALSE);
1441   fhChi2PerClusterTPC[0]->Draw();
1442
1443   canvas1->cd(3);
1444   fhNSigmaToVertex[0]->SetStats(kFALSE);
1445   fhNSigmaToVertex[0]->GetXaxis()->SetRangeUser(0, 10);
1446   fhNSigmaToVertex[0]->Draw();
1447
1448   canvas1->SaveAs(Form("%s_%s.gif", GetName(), canvas1->GetName()));
1449
1450   TCanvas* canvas2 = new TCanvas(Form("%s_2", GetName()), "Track Quality Results2", 1200, 800);
1451   canvas2->Divide(3, 2);
1452
1453   canvas2->cd(1);
1454   fhC11[0]->SetStats(kFALSE);
1455   gPad->SetLogy();
1456   fhC11[0]->Draw();
1457
1458   canvas2->cd(2);
1459   fhC22[0]->SetStats(kFALSE);
1460   gPad->SetLogy();
1461   fhC22[0]->Draw();
1462
1463   canvas2->cd(3);
1464   fhC33[0]->SetStats(kFALSE);
1465   gPad->SetLogy();
1466   fhC33[0]->Draw();
1467
1468   canvas2->cd(4);
1469   fhC44[0]->SetStats(kFALSE);
1470   gPad->SetLogy();
1471   fhC44[0]->Draw();
1472
1473   canvas2->cd(5);
1474   fhC55[0]->SetStats(kFALSE);
1475   gPad->SetLogy();
1476   fhC55[0]->Draw();
1477
1478   canvas2->cd(6);
1479   fhRel1PtUncertainty[0]->SetStats(kFALSE);
1480   gPad->SetLogy();
1481   fhRel1PtUncertainty[0]->Draw();
1482
1483   canvas2->SaveAs(Form("%s_%s.gif", GetName(), canvas2->GetName()));
1484
1485   TCanvas* canvas3 = new TCanvas(Form("%s_3", GetName()), "Track Quality Results3", 1200, 800);
1486   canvas3->Divide(3, 2);
1487
1488   canvas3->cd(1);
1489   fhDXY[0]->SetStats(kFALSE);
1490   gPad->SetLogy();
1491   fhDXY[0]->Draw();
1492
1493   canvas3->cd(2);
1494   fhDZ[0]->SetStats(kFALSE);
1495   gPad->SetLogy();
1496   fhDZ[0]->Draw();
1497
1498   canvas3->cd(3);
1499   fhDXYvsDZ[0]->SetStats(kFALSE);
1500   gPad->SetLogz();
1501   gPad->SetRightMargin(0.15);
1502   fhDXYvsDZ[0]->Draw("COLZ");
1503
1504   canvas3->cd(4);
1505   fhDXYNormalized[0]->SetStats(kFALSE);
1506   gPad->SetLogy();
1507   fhDXYNormalized[0]->Draw();
1508
1509   canvas3->cd(5);
1510   fhDZNormalized[0]->SetStats(kFALSE);
1511   gPad->SetLogy();
1512   fhDZNormalized[0]->Draw();
1513
1514   canvas3->cd(6);
1515   fhDXYvsDZNormalized[0]->SetStats(kFALSE);
1516   gPad->SetLogz();
1517   gPad->SetRightMargin(0.15);
1518   fhDXYvsDZNormalized[0]->Draw("COLZ");
1519
1520   canvas3->SaveAs(Form("%s_%s.gif", GetName(), canvas3->GetName()));
1521
1522   TCanvas* canvas4 = new TCanvas(Form("%s_4", GetName()), "Track Quality Results4", 800, 500);
1523   canvas4->Divide(2, 1);
1524
1525   canvas4->cd(1);
1526   fhCutStatistics->SetStats(kFALSE);
1527   fhCutStatistics->LabelsOption("v");
1528   gPad->SetBottomMargin(0.3);
1529   fhCutStatistics->Draw();
1530
1531   canvas4->cd(2);
1532   fhCutCorrelation->SetStats(kFALSE);
1533   fhCutCorrelation->LabelsOption("v");
1534   gPad->SetBottomMargin(0.3);
1535   gPad->SetLeftMargin(0.3);
1536   fhCutCorrelation->Draw("COLZ");
1537
1538   canvas4->SaveAs(Form("%s_%s.gif", GetName(), canvas4->GetName()));
1539
1540   /*canvas->cd(1);
1541   fhDXYvsDZNormalized[0]->SetStats(kFALSE);
1542   fhDXYvsDZNormalized[0]->DrawCopy("COLZ");
1543
1544   canvas->cd(2);
1545   fhNClustersTPC[0]->SetStats(kFALSE);
1546   fhNClustersTPC[0]->DrawCopy();
1547
1548   canvas->cd(3);
1549   fhChi2PerClusterITS[0]->SetStats(kFALSE);
1550   fhChi2PerClusterITS[0]->DrawCopy();
1551   fhChi2PerClusterITS[1]->SetLineColor(2);
1552   fhChi2PerClusterITS[1]->DrawCopy("SAME");
1553
1554   canvas->cd(4);
1555   fhChi2PerClusterTPC[0]->SetStats(kFALSE);
1556   fhChi2PerClusterTPC[0]->DrawCopy();
1557   fhChi2PerClusterTPC[1]->SetLineColor(2);
1558   fhChi2PerClusterTPC[1]->DrawCopy("SAME");*/
1559 }
1560 //--------------------------------------------------------------------------
1561 void AliESDtrackCuts::SetPtDepDCACuts(Double_t pt) {
1562   //
1563   // set the pt-dependent DCA cuts
1564   //
1565
1566   if(f1CutMaxDCAToVertexXYPtDep) {
1567      fCutMaxDCAToVertexXY=f1CutMaxDCAToVertexXYPtDep->Eval(pt);
1568   }
1569
1570   if(f1CutMaxDCAToVertexZPtDep) {
1571     fCutMaxDCAToVertexZ=f1CutMaxDCAToVertexZPtDep->Eval(pt);
1572   }
1573
1574   if(f1CutMinDCAToVertexXYPtDep) {
1575     fCutMinDCAToVertexXY=f1CutMinDCAToVertexXYPtDep->Eval(pt);
1576   }
1577
1578   if(f1CutMinDCAToVertexZPtDep) {
1579     fCutMinDCAToVertexZ=f1CutMinDCAToVertexZPtDep->Eval(pt);
1580   }
1581
1582
1583   return;
1584 }
1585
1586
1587
1588 //--------------------------------------------------------------------------
1589 Bool_t AliESDtrackCuts::CheckPtDepDCA(TString dist,Bool_t print) const {
1590   //
1591   // Check the correctness of the string syntax
1592   //
1593   Bool_t retval=kTRUE;
1594
1595   if(!dist.Contains("pt")) {
1596     if(print) printf("AliESDtrackCuts::CheckPtDepDCA(): string must contain \"pt\"\n");
1597     retval= kFALSE;
1598   } 
1599   return retval;
1600 }
1601
1602  void AliESDtrackCuts::SetMaxDCAToVertexXYPtDep(const char *dist){
1603    if(!CheckPtDepDCA(dist,kTRUE)) return;
1604    if(f1CutMaxDCAToVertexXYPtDep)delete f1CutMaxDCAToVertexXYPtDep;
1605    fCutMaxDCAToVertexXYPtDep = dist;
1606    TString tmp(dist);
1607    tmp.ReplaceAll("pt","x");
1608    f1CutMaxDCAToVertexXYPtDep = new TFormula("f1CutMaxDCAToVertexXYPtDep",tmp.Data());
1609  
1610 }
1611
1612  void AliESDtrackCuts::SetMaxDCAToVertexZPtDep(const char *dist){
1613    if(!CheckPtDepDCA(dist,kTRUE)) return;
1614    if(f1CutMaxDCAToVertexZPtDep)delete f1CutMaxDCAToVertexZPtDep;
1615    fCutMaxDCAToVertexZPtDep = dist;
1616    TString tmp(dist);
1617    tmp.ReplaceAll("pt","x");
1618    f1CutMaxDCAToVertexZPtDep = new TFormula("f1CutMaxDCAToVertexZPtDep",tmp.Data());
1619
1620    
1621 }
1622
1623
1624  void AliESDtrackCuts::SetMinDCAToVertexXYPtDep(const char *dist){
1625    if(!CheckPtDepDCA(dist,kTRUE)) return;
1626    if(f1CutMinDCAToVertexXYPtDep)delete f1CutMinDCAToVertexXYPtDep;
1627    fCutMinDCAToVertexXYPtDep = dist;
1628    TString tmp(dist);
1629    tmp.ReplaceAll("pt","x");
1630    f1CutMinDCAToVertexXYPtDep = new TFormula("f1CutMinDCAToVertexXYPtDep",tmp.Data());
1631
1632 }
1633
1634
1635  void AliESDtrackCuts::SetMinDCAToVertexZPtDep(const char *dist){
1636    if(!CheckPtDepDCA(dist,kTRUE)) return;
1637    if(f1CutMinDCAToVertexZPtDep)delete f1CutMinDCAToVertexZPtDep;
1638    fCutMinDCAToVertexZPtDep = dist;
1639    TString tmp(dist);
1640    tmp.ReplaceAll("pt","x");
1641    f1CutMinDCAToVertexZPtDep = new TFormula("f1CutMinDCAToVertexZPtDep",tmp.Data());
1642 }