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