Cut class for the filtering of ESD V0s (Boris Hippolyte <hippolyt@in2p3.fr>)
[u/mrichter/AliRoot.git] / ANALYSIS / AliESDv0Cuts.cxx
CommitLineData
94b50f5f 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: AliESDv0Cuts.cxx 24534 2008-03-16 22:22:11Z fca $ */
17
18#include "AliESDv0Cuts.h"
19
20#include <AliESDVertex.h>
21#include <AliESDtrack.h>
22#include <AliESDv0.h>
23#include <AliESD.h>
24#include <AliESDEvent.h>
25#include <AliLog.h>
26
27#include <TTree.h>
28#include <TCanvas.h>
29#include <TDirectory.h>
30
31//____________________________________________________________________
32ClassImp(AliESDv0Cuts)
33
34// Cut names
35const Char_t* AliESDv0Cuts::fgkCutNames[kNCuts] = {
36 "dca positive to pvtx",
37 "dca negative to pvtx",
38 "#Chi^{2}",
39 "dca v0 daughters",
40 "min decay radius",
41 "max decay radius",
42 "cosine pointing angle",
43 "on-the-fly status",
44 "dca v0 to pvtx",
45 "p",
46 "p_{T}",
47 "p_{x}",
48 "p_{y}",
49 "p_{z}"
50};
51
52//____________________________________________________________________
53AliESDv0Cuts::AliESDv0Cuts(const Char_t* name, const Char_t* title) : AliAnalysisCuts(name,title),
54 fCutMinDcaPosToVertex(0),
55 fCutMinDcaNegToVertex(0),
56 fCutMaxChi2(40),
57 fCutMaxDcaV0Daughters(0),
58 fCutMinRadius(0),
59 fCutMaxRadius(0),
60 fCutMinCosinePointingAngle(0),
61 fCutRequireOnFlyStatus(0),
62 fCutMaxDcaV0ToVertex(0),
63 fPMin(0),
64 fPMax(0),
65 fPtMin(0),
66 fPtMax(0),
67 fPxMin(0),
68 fPxMax(0),
69 fPyMin(0),
70 fPyMax(0),
71 fPzMin(0),
72 fPzMax(0),
73 fHistogramsOn(0),
74 fhCutStatistics(0),
75 fhCutCorrelation(0)
76{
77 //
78 // constructor
79 //
80
81 Init();
82
83 //##############################################################################
84 // setting default cuts
85 SetMinDcaPosToVertex();
86 SetMinDcaNegToVertex();
87 SetMaxChi2();
88 SetMaxDcaV0Daughters();
89 SetMinRadius();
90 SetMaxRadius();
91 SetMinCosinePointingAngle();
92 SetRequireOnFlyStatus();
93 SetMaxDcaV0ToVertex();
94 SetPRange();
95 SetPtRange();
96 SetPxRange();
97 SetPyRange();
98 SetPzRange();
99
100 SetHistogramsOn();
101}
102
103//_____________________________________________________________________________
104AliESDv0Cuts::AliESDv0Cuts(const AliESDv0Cuts &c) : AliAnalysisCuts(c),
105 fCutMinDcaPosToVertex(0),
106 fCutMinDcaNegToVertex(0),
107 fCutMaxChi2(0),
108 fCutMaxDcaV0Daughters(0),
109 fCutMinRadius(0),
110 fCutMaxRadius(0),
111 fCutMinCosinePointingAngle(0),
112 fCutRequireOnFlyStatus(0),
113 fCutMaxDcaV0ToVertex(0),
114 fPMin(0),
115 fPMax(0),
116 fPtMin(0),
117 fPtMax(0),
118 fPxMin(0),
119 fPxMax(0),
120 fPyMin(0),
121 fPyMax(0),
122 fPzMin(0),
123 fPzMax(0),
124 fHistogramsOn(0),
125 fhCutStatistics(0),
126 fhCutCorrelation(0)
127{
128 //
129 // copy constructor
130 //
131
132 ((AliESDv0Cuts &) c).Copy(*this);
133}
134
135AliESDv0Cuts::~AliESDv0Cuts()
136{
137 //
138 // destructor
139 //
140
141 for (Int_t i=0; i<2; i++) {
142
143 if (fhDcaPosToVertex[i])
144 delete fhDcaPosToVertex[i];
145 if (fhDcaNegToVertex[i])
146 delete fhDcaNegToVertex[i];
147 if (fhChi2[i])
148 delete fhChi2[i];
149 if (fhDcaV0Daughters[i])
150 delete fhDcaV0Daughters[i];
151 if (fhRadius[i])
152 delete fhRadius[i];
153 if (fhCosinePointingAngle[i])
154 delete fhCosinePointingAngle[i];
155 if (fhOnFlyStatus[i])
156 delete fhOnFlyStatus[i];
157 if (fhDcaV0ToVertex[i])
158 delete fhDcaV0ToVertex[i];
159 if (fhPt[i])
160 delete fhPt[i];
161 }
162
163 if (fhCutStatistics)
164 delete fhCutStatistics;
165 if (fhCutCorrelation)
166 delete fhCutCorrelation;
167}
168
169void AliESDv0Cuts::Init()
170{
171 //
172 // sets everything to zero
173 //
174 fCutMinDcaPosToVertex = 0;
175 fCutMinDcaNegToVertex = 0;
176 fCutMaxChi2 = 0;
177 fCutMaxDcaV0Daughters = 0;
178 fCutMinRadius = 0;
179 fCutMaxRadius = 0;
180 fCutMinCosinePointingAngle = 0;
181 fCutRequireOnFlyStatus = 0;
182 fCutMaxDcaV0ToVertex = 0;
183
184 fPMin = 0;
185 fPMax = 0;
186 fPtMin = 0;
187 fPtMax = 0;
188 fPxMin = 0;
189 fPxMax = 0;
190 fPyMin = 0;
191 fPyMax = 0;
192 fPzMin = 0;
193 fPzMax = 0;
194
195 fHistogramsOn = kFALSE;
196
197 for (Int_t i=0; i<2; ++i)
198 {
199 fhDcaPosToVertex[i] = 0;
200 fhDcaNegToVertex[i] = 0;
201 fhChi2[i] = 0;
202 fhDcaV0Daughters[i] = 0;
203 fhRadius[i] = 0;
204 fhCosinePointingAngle[i] = 0;
205 fhOnFlyStatus[i] = 0;
206 fhDcaV0ToVertex[i] = 0;
207
208 fhPt[i] = 0;
209 }
210 fhCutStatistics = 0;
211 fhCutCorrelation = 0;
212}
213
214//_____________________________________________________________________________
215AliESDv0Cuts &AliESDv0Cuts::operator=(const AliESDv0Cuts &c)
216{
217 //
218 // Assignment operator
219 //
220
221 if (this != &c) ((AliESDv0Cuts &) c).Copy(*this);
222 return *this;
223}
224
225//_____________________________________________________________________________
226void AliESDv0Cuts::Copy(TObject &c) const
227{
228 //
229 // Copy function
230 //
231
232 AliESDv0Cuts& target = (AliESDv0Cuts &) c;
233
234 target.Init();
235
236 target.fCutMinDcaPosToVertex = fCutMinDcaPosToVertex;
237 target.fCutMinDcaNegToVertex = fCutMinDcaNegToVertex;
238 target.fCutMaxChi2 = fCutMaxChi2;
239 target.fCutMaxDcaV0Daughters = fCutMaxDcaV0Daughters;
240 target.fCutMinRadius = fCutMinRadius;
241 target.fCutMaxRadius = fCutMaxRadius;
242 target.fCutMinCosinePointingAngle = fCutMinCosinePointingAngle;
243 target.fCutRequireOnFlyStatus = fCutRequireOnFlyStatus;
244 target.fCutMaxDcaV0ToVertex = fCutMaxDcaV0ToVertex;
245
246 target.fPMin = fPMin;
247 target.fPMax = fPMax;
248 target.fPtMin = fPtMin;
249 target.fPtMax = fPtMax;
250 target.fPxMin = fPxMin;
251 target.fPxMax = fPxMax;
252 target.fPyMin = fPyMin;
253 target.fPyMax = fPyMax;
254 target.fPzMin = fPzMin;
255 target.fPzMax = fPzMax;
256
257 target.fHistogramsOn = fHistogramsOn;
258
259 for (Int_t i=0; i<2; ++i)
260 {
261 if (fhDcaPosToVertex[i]) target.fhDcaPosToVertex[i] = (TH1F*) fhDcaPosToVertex[i]->Clone();
262 if (fhDcaNegToVertex[i]) target.fhDcaNegToVertex[i] = (TH1F*) fhDcaNegToVertex[i]->Clone();
263 if (fhChi2[i]) target.fhChi2[i] = (TH1F*) fhChi2[i]->Clone();
264 if (fhDcaV0Daughters[i]) target.fhDcaV0Daughters[i] = (TH1F*) fhDcaV0Daughters[i]->Clone();
265 if (fhRadius[i]) target.fhRadius[i] = (TH1F*) fhRadius[i]->Clone();
266 if (fhCosinePointingAngle[i]) target.fhCosinePointingAngle[i] = (TH1F*) fhCosinePointingAngle[i]->Clone();
267 if (fhOnFlyStatus[i]) target.fhOnFlyStatus[i] = (TH1F*) fhOnFlyStatus[i]->Clone();
268 if (fhDcaV0ToVertex[i]) target.fhDcaV0ToVertex[i] = (TH1F*) fhDcaV0ToVertex[i]->Clone();
269
270 if (fhPt[i]) target.fhPt[i] = (TH1F*) fhPt[i]->Clone();
271 }
272 if (fhCutStatistics) target.fhCutStatistics = (TH1F*) fhCutStatistics->Clone();
273 if (fhCutCorrelation) target.fhCutCorrelation = (TH2F*) fhCutCorrelation->Clone();
274
275 TNamed::Copy(c);
276}
277
278//_____________________________________________________________________________
279Long64_t AliESDv0Cuts::Merge(TCollection* list) {
280 // Merge a list of AliESDv0Cuts objects with this (needed for PROOF)
281 // Returns the number of merged objects (including this)
282
283 if (!list)
284 return 0;
285
286 if (list->IsEmpty())
287 return 1;
288
289 if (!fHistogramsOn)
290 return 0;
291
292 TIterator* iter = list->MakeIterator();
293 TObject* obj;
294
295
296 // collection of measured and generated histograms
297 Int_t count = 0;
298 while ((obj = iter->Next())) {
299
300 AliESDv0Cuts* entry = dynamic_cast<AliESDv0Cuts*>(obj);
301 if (entry == 0)
302 continue;
303
304 if (!entry->fHistogramsOn)
305 continue;
306
307 for (Int_t i=0; i<2; i++) {
308
309 fhDcaPosToVertex[i] ->Add(entry->fhDcaPosToVertex[i] );
310 fhDcaNegToVertex[i] ->Add(entry->fhDcaNegToVertex[i] );
311 fhChi2[i] ->Add(entry->fhChi2[i] );
312 fhDcaV0Daughters[i] ->Add(entry->fhDcaV0Daughters[i] );
313 fhRadius[i] ->Add(entry->fhRadius[i] );
314 fhCosinePointingAngle[i]->Add(entry->fhCosinePointingAngle[i]);
315 fhOnFlyStatus[i] ->Add(entry->fhOnFlyStatus[i] );
316 fhDcaV0ToVertex[i] ->Add(entry->fhDcaV0ToVertex[i] );
317
318 fhPt[i] ->Add(entry->fhPt[i]);
319 }
320 fhCutStatistics ->Add(entry->fhCutStatistics);
321 fhCutCorrelation ->Add(entry->fhCutCorrelation);
322
323 count++;
324 }
325
326 return count+1;
327}
328
329void AliESDv0Cuts::EnableNeededBranches(TTree* tree)
330{
331 // enables the branches needed by AcceptV0, for a list see comment of AcceptV0
332
333 tree->SetBranchStatus("fV0s.fDcaV0Daughters", 1);
334 tree->SetBranchStatus("fV0s.fChi2V0", 1);
335 tree->SetBranchStatus("fV0s.fPos*", 1);
336 tree->SetBranchStatus("fV0s.fNmom*", 1);
337 tree->SetBranchStatus("fV0s.fPmom*", 1);
338 tree->SetBranchStatus("fV0s.fRr", 1);
339 tree->SetBranchStatus("fV0s.fPointAngle*", 1);
340 tree->SetBranchStatus("fV0s.fOnFlyStatus", 1);
341}
342
343//____________________________________________________________________
344Bool_t
345AliESDv0Cuts::IsSelected(TList* listObj) {
346 if(listObj->GetSize()!=4) return kFALSE;
347 AliESDv0 *esdV0 = (AliESDv0*)listObj->At(0);
348 AliESDtrack *trackPos = (AliESDtrack*)listObj->At(1);
349 AliESDtrack *trackNeg = (AliESDtrack*)listObj->At(2);
350 const AliESDVertex *esdVertex = (AliESDVertex*)listObj->At(3);
351 return AcceptV0(esdV0,trackPos,trackNeg,esdVertex);
352}
353
354//____________________________________________________________________
355Bool_t
356AliESDv0Cuts::AcceptV0(AliESDv0* esdV0, AliESDtrack* trackPos, AliESDtrack* trackNeg, const AliESDVertex* esdVertex) {
357 //
358 // figure out if the v0s survives all the v0 cuts defined
359 //
360 // the different quality parameter and kinematic values are first
361 // retrieved from the v0. then it is found out what cuts the
362 // v0 did not survive and finally the cuts are imposed.
363
364 // this function needs the following branches... but this is not enough
365 // fV0s.fDcaV0Daughters
366 // fV0s.fChi2V0
367 // fV0s.fPos*
368 // fV0s.fNmom*
369 // fV0s.fPmom*
370 // fV0s.fRr
371 // fV0s.fPointAngle
372 // fV0s.fOnFlyStatus
373
374 Float_t dcaPosToVertex = 0, dcaNegToVertex = 0;
375 Float_t tdcaPosToVertex[2], tdcaNegToVertex[2];
376
377 if (trackPos) trackPos->GetImpactParameters(tdcaPosToVertex[0],tdcaPosToVertex[1]);
378 if (trackNeg) trackNeg->GetImpactParameters(tdcaNegToVertex[0],tdcaNegToVertex[1]);
379
380 dcaPosToVertex = TMath::Sqrt(tdcaPosToVertex[0]*tdcaPosToVertex[0]+tdcaPosToVertex[1]*tdcaPosToVertex[1]);
381 dcaNegToVertex = TMath::Sqrt(tdcaNegToVertex[0]*tdcaNegToVertex[0]+tdcaNegToVertex[1]*tdcaNegToVertex[1]);
382
383 UInt_t status = esdV0->GetOnFlyStatus();
384 Float_t chi2 = esdV0->GetChi2V0();
385
386 Double_t dcaV0Daughters = esdV0->GetDcaV0Daughters();
387
388 Double_t vtxPosition[3]; esdVertex->GetXYZ(vtxPosition);
389 Double_t dcaV0ToVertex = esdV0->GetD(vtxPosition[0],vtxPosition[1],vtxPosition[2]);
390
391 Double_t v0Position[3];
392 esdV0->GetXYZ(v0Position[0],v0Position[1],v0Position[2]);
393 Double_t radius = TMath::Sqrt(TMath::Power(v0Position[0],2) + TMath::Power(v0Position[1],2));
394 Double_t v0CosinePointingAngle = esdV0->GetV0CosineOfPointingAngle();
395
396 // getting the kinematic variables of the v0
397 Double_t p[3];
398 esdV0->GetPxPyPz(p[0],p[1],p[2]);
399 Float_t momentum = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2) + TMath::Power(p[2],2));
400 Float_t pt = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2));
401
402 //########################################################################
403 // cut the v0?
404
405 Bool_t cuts[kNCuts];
406 for (Int_t i=0; i<kNCuts; i++) cuts[i]=kFALSE;
407
408 // v0 quality cuts
409 if (dcaPosToVertex < fCutMinDcaPosToVertex)
410 cuts[0]=kTRUE;
411 if (dcaNegToVertex < fCutMinDcaNegToVertex)
412 cuts[1]=kTRUE;
413 if (chi2 > fCutMaxChi2)
414 cuts[2]=kTRUE;
415 if (dcaV0Daughters > fCutMaxDcaV0Daughters)
416 cuts[3]=kTRUE;
417 if (radius < fCutMinRadius)
418 cuts[4]=kTRUE;
419 if (radius > fCutMaxRadius)
420 cuts[5]=kTRUE;
421 if (v0CosinePointingAngle < fCutMinCosinePointingAngle)
422 cuts[6]=kTRUE;
423 if (fCutRequireOnFlyStatus && !status)
424 cuts[7]=kTRUE;
425 if (dcaV0ToVertex > fCutMaxDcaV0ToVertex)
426 cuts[8] = kTRUE;
427
428 // v0 kinematics cut
429 if((momentum < fPMin) || (momentum > fPMax))
430 cuts[9]=kTRUE;
431 if((pt < fPtMin) || (pt > fPtMax))
432 cuts[10] = kTRUE;
433 if((p[0] < fPxMin) || (p[0] > fPxMax))
434 cuts[11] = kTRUE;
435 if((p[1] < fPyMin) || (p[1] > fPyMax))
436 cuts[12] = kTRUE;
437 if((p[2] < fPzMin) || (p[2] > fPzMax))
438 cuts[13] = kTRUE;
439
440 Bool_t cut=kFALSE;
441 for (Int_t i=0; i<kNCuts; i++)
442 if (cuts[i]) cut = kTRUE;
443
444 //########################################################################
445 // filling histograms
446 if (fHistogramsOn) {
447 fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n v0s")));
448
449 if (cut)
450 fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n cut v0s")));
451
452 for (Int_t i=0; i<kNCuts; i++) {
453 if (cuts[i])
454 fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin(fgkCutNames[i])));
455
456 for (Int_t j=i; j<kNCuts; j++) {
457 if (cuts[i] && cuts[j]) {
458 Float_t x = fhCutCorrelation->GetXaxis()->GetBinCenter(fhCutCorrelation->GetXaxis()->FindBin(fgkCutNames[i]));
459 Float_t y = fhCutCorrelation->GetYaxis()->GetBinCenter(fhCutCorrelation->GetYaxis()->FindBin(fgkCutNames[j]));
460 fhCutCorrelation->Fill(x,y);
461 }
462 }
463 }
464
465 fhDcaPosToVertex[0]->Fill(dcaPosToVertex);
466 fhDcaNegToVertex[0]->Fill(dcaNegToVertex);
467 fhChi2[0]->Fill(chi2);
468 fhDcaV0Daughters[0]->Fill(dcaV0Daughters);
469 fhRadius[0]->Fill(radius);
470 fhCosinePointingAngle[0]->Fill(v0CosinePointingAngle);
471 fhOnFlyStatus[0]->Fill(status);
472 fhDcaV0ToVertex[0]->Fill(dcaV0ToVertex);
473
474 fhPt[0]->Fill(pt);
475 }
476
477 //########################################################################
478 // cut the v0!
479 if (cut) return kFALSE;
480
481 //########################################################################
482 // filling histograms after cut
483 if (fHistogramsOn) {
484 fhDcaPosToVertex[1]->Fill(dcaPosToVertex);
485 fhDcaNegToVertex[1]->Fill(dcaNegToVertex);
486 fhChi2[1]->Fill(chi2);
487 fhDcaV0Daughters[1]->Fill(dcaV0Daughters);
488 fhRadius[1]->Fill(radius);
489 fhCosinePointingAngle[1]->Fill(v0CosinePointingAngle);
490 fhOnFlyStatus[1]->Fill(status);
491 fhDcaV0ToVertex[1]->Fill(dcaV0ToVertex);
492
493 fhPt[1]->Fill(pt);
494 }
495
496 return kTRUE;
497}
498
499//____________________________________________________________________
500TObjArray* AliESDv0Cuts::GetAcceptedV0s(AliESD* esd)
501{
502 //
503 // returns an array of all v0s that pass the cuts
504 //
505
506 TObjArray* acceptedV0s = new TObjArray();
507 // const AliESDVertex *spdVertex = esd->GetVertex();
508 const AliESDVertex *primaryVertex = esd->GetPrimaryVertex();
509 Int_t lIndexTrackPos = 0, lIndexTrackNeg = 0;
510
511 // loop over esd v0s
512 for (Int_t iV0 = 0; iV0 < esd->GetNumberOfV0s(); iV0++) {
513 AliESDv0* v0 = esd->GetV0(iV0);
514
515 lIndexTrackPos = TMath::Abs(v0->GetPindex());
516 lIndexTrackNeg = TMath::Abs(v0->GetNindex());
517 AliESDtrack *trackPos = esd->GetTrack(lIndexTrackPos);
518 AliESDtrack *trackNeg = esd->GetTrack(lIndexTrackNeg);
519
520 if (AcceptV0(v0,trackPos,trackNeg,primaryVertex))
521 acceptedV0s->Add(v0);
522 }
523
524 return acceptedV0s;
525}
526
527//____________________________________________________________________
528Int_t AliESDv0Cuts::CountAcceptedV0s(AliESD* esd)
529{
530 //
531 // returns an the number of v0s that pass the cuts
532 //
533
534 Int_t count = 0;
535 // const AliESDVertex *spdVertex = esd->GetVertex();
536 const AliESDVertex *primaryVertex = esd->GetPrimaryVertex();
537 Int_t lIndexTrackPos = 0, lIndexTrackNeg = 0;
538
539 // loop over esd v0s
540 for (Int_t iV0 = 0; iV0 < esd->GetNumberOfV0s(); iV0++) {
541 AliESDv0* v0 = esd->GetV0(iV0);
542
543 lIndexTrackPos = TMath::Abs(v0->GetPindex());
544 lIndexTrackNeg = TMath::Abs(v0->GetNindex());
545 AliESDtrack *trackPos = esd->GetTrack(lIndexTrackPos);
546 AliESDtrack *trackNeg = esd->GetTrack(lIndexTrackNeg);
547
548 if (AcceptV0(v0,trackPos,trackNeg,primaryVertex))
549 count++;
550 }
551
552 return count;
553}
554
555//____________________________________________________________________
556TObjArray* AliESDv0Cuts::GetAcceptedV0s(AliESDEvent* esd)
557{
558 //
559 // returns an array of all v0s that pass the cuts
560 //
561
562 TObjArray* acceptedV0s = new TObjArray();
563 // const AliESDVertex *spdVertex = esd->GetVertex();
564 const AliESDVertex *primaryVertex = esd->GetPrimaryVertex();
565 Int_t lIndexTrackPos = 0, lIndexTrackNeg = 0;
566
567 // loop over esd v0s
568 for (Int_t iV0 = 0; iV0 < esd->GetNumberOfV0s(); iV0++) {
569 AliESDv0* v0 = esd->GetV0(iV0);
570
571 lIndexTrackPos = TMath::Abs(v0->GetPindex());
572 lIndexTrackNeg = TMath::Abs(v0->GetNindex());
573 AliESDtrack *trackPos = esd->GetTrack(lIndexTrackPos);
574 AliESDtrack *trackNeg = esd->GetTrack(lIndexTrackNeg);
575
576 if (AcceptV0(v0,trackPos,trackNeg,primaryVertex))
577 acceptedV0s->Add(v0);
578 }
579
580 return acceptedV0s;
581}
582
583//____________________________________________________________________
584Int_t AliESDv0Cuts::CountAcceptedV0s(AliESDEvent* esd)
585{
586 //
587 // returns an the number of v0s that pass the cuts
588 //
589
590 Int_t count = 0;
591 // const AliESDVertex *spdVertex = esd->GetVertex();
592 const AliESDVertex *primaryVertex = esd->GetPrimaryVertex();
593 Int_t lIndexTrackPos = 0, lIndexTrackNeg = 0;
594
595 // loop over esd v0s
596 for (Int_t iV0 = 0; iV0 < esd->GetNumberOfV0s(); iV0++) {
597 AliESDv0* v0 = esd->GetV0(iV0);
598
599 lIndexTrackPos = TMath::Abs(v0->GetPindex());
600 lIndexTrackNeg = TMath::Abs(v0->GetNindex());
601 AliESDtrack *trackPos = esd->GetTrack(lIndexTrackPos);
602 AliESDtrack *trackNeg = esd->GetTrack(lIndexTrackNeg);
603
604 if (AcceptV0(v0,trackPos,trackNeg,primaryVertex))
605 count++;
606 }
607
608 return count;
609}
610
611//____________________________________________________________________
612 void AliESDv0Cuts::DefineHistograms(Int_t color) {
613 //
614 // diagnostics histograms are defined
615 //
616
617 fHistogramsOn=kTRUE;
618
619 Bool_t oldStatus = TH1::AddDirectoryStatus();
620 TH1::AddDirectory(kFALSE);
621
622 //###################################################################################
623 // defining histograms
624
625 fhCutStatistics = new TH1F("cut_statistics","cut statistics",kNCuts+4,-0.5,kNCuts+3.5);
626
627 fhCutStatistics->GetXaxis()->SetBinLabel(1,"n v0s");
628 fhCutStatistics->GetXaxis()->SetBinLabel(2,"n cut v0s");
629
630 fhCutCorrelation = new TH2F("cut_correlation","cut correlation",kNCuts,-0.5,kNCuts-0.5,kNCuts,-0.5,kNCuts-0.5);;
631
632 for (Int_t i=0; i<kNCuts; i++) {
633 fhCutStatistics->GetXaxis()->SetBinLabel(i+4,fgkCutNames[i]);
634 fhCutCorrelation->GetXaxis()->SetBinLabel(i+1,fgkCutNames[i]);
635 fhCutCorrelation->GetYaxis()->SetBinLabel(i+1,fgkCutNames[i]);
636 }
637
638 fhCutStatistics ->SetLineColor(color);
639 fhCutCorrelation ->SetLineColor(color);
640 fhCutStatistics ->SetLineWidth(2);
641 fhCutCorrelation ->SetLineWidth(2);
642
643 Char_t str[256];
644 for (Int_t i=0; i<2; i++) {
645 if (i==0) sprintf(str," ");
646 else sprintf(str,"_cut");
647
648 fhDcaPosToVertex[i] = new TH1F(Form("dcaPosToVertex%s",str),"",120,0,3);
649 fhDcaNegToVertex[i] = new TH1F(Form("dcaNegToVertex%s",str),"",120,0,3);
650 fhChi2[i] = new TH1F(Form("chi2%s",str),"",50,0,50);
651 fhDcaV0Daughters[i] = new TH1F(Form("dcaV0Daughters%s",str),"",200,0,5);
652 fhRadius[i] = new TH1F(Form("decayRadius%s",str),"",300,0,150);
653 fhCosinePointingAngle[i] = new TH1F(Form("cosinePointingAngle%s",str),"",100,-1,1);
654 fhOnFlyStatus[i] = new TH1F(Form("onflystatus%s",str),"",5,0,5);
655 fhDcaV0ToVertex[i] = new TH1F(Form("dcaV0ToVertex%s",str),"",100,0,5);
656
657 fhPt[i] = new TH1F(Form("pt%s",str) ,"p_{T} distribution;p_{T} (GeV/c)",500,0.0,100.0);
658
659 fhDcaPosToVertex[i]->SetTitle("Dca of positive daughter to parent vertex");
660 fhDcaNegToVertex[i]->SetTitle("Dca of negative daughter to parent vertex");
661 fhChi2[i]->SetTitle("#Chi^{2} for v0");
662 fhDcaV0Daughters[i]->SetTitle("Dca between daughter tracks");
663 fhRadius[i]->SetTitle("Decay radius of the v0");
664 fhCosinePointingAngle[i]->SetTitle("Cosine of the Pointing Angle");
665 fhOnFlyStatus[i]->SetTitle("On-the-Fly Status");
666 fhDcaV0ToVertex[i]->SetTitle("Dca of v0 to parent vertex");
667
668 fhDcaPosToVertex[i]->SetLineColor(color); fhDcaPosToVertex[i]->SetLineWidth(2);
669 fhDcaNegToVertex[i]->SetLineColor(color); fhDcaNegToVertex[i]->SetLineWidth(2);
670 fhChi2[i]->SetLineColor(color); fhChi2[i]->SetLineWidth(2);
671 fhDcaV0Daughters[i]->SetLineColor(color); fhDcaV0Daughters[i]->SetLineWidth(2);
672 fhRadius[i]->SetLineColor(color); fhRadius[i]->SetLineWidth(2);
673 fhCosinePointingAngle[i]->SetLineColor(color); fhCosinePointingAngle[i]->SetLineWidth(2);
674 fhOnFlyStatus[i]->SetLineColor(color); fhOnFlyStatus[i]->SetLineWidth(2);
675 fhDcaV0ToVertex[i]->SetLineColor(color); fhDcaV0ToVertex[i]->SetLineWidth(2);
676 }
677
678 TH1::AddDirectory(oldStatus);
679}
680
681//____________________________________________________________________
682Bool_t AliESDv0Cuts::LoadHistograms(const Char_t* dir)
683{
684 //
685 // loads the histograms from a file
686 // if dir is empty a directory with the name of this object is taken (like in SaveHistogram)
687 //
688
689 if (!dir)
690 dir = GetName();
691
692 if (!gDirectory->cd(dir))
693 return kFALSE;
694
695 fhCutStatistics = dynamic_cast<TH1F*> (gDirectory->Get("cut_statistics"));
696 fhCutCorrelation = dynamic_cast<TH2F*> (gDirectory->Get("cut_correlation"));
697
698 Char_t str[5];
699 for (Int_t i=0; i<2; i++) {
700 if (i==0)
701 {
702 gDirectory->cd("before_cuts");
703 str[0] = 0;
704 }
705 else
706 {
707 gDirectory->cd("after_cuts");
708 sprintf(str,"_cut");
709 }
710
711 fhDcaPosToVertex[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("dcaPosToVertex%s",str) ));
712 fhDcaNegToVertex[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("dcaNegToVertex%s",str) ));
713 fhChi2[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("chi2%s",str)));
714 fhDcaV0Daughters[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("dcaV0Daughters%s",str)));
715 fhRadius[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("decayRadius%s",str)));
716 fhCosinePointingAngle[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("cosinepointingangle%s",str)));
717 fhOnFlyStatus[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("onflystatus%s",str)));
718 fhDcaV0ToVertex[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("dcaV0ToVertex%s",str)));
719
720 fhPt[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("pt%s",str)));
721
722 gDirectory->cd("../");
723 }
724
725 gDirectory->cd("..");
726
727 return kTRUE;
728}
729
730//____________________________________________________________________
731void AliESDv0Cuts::SaveHistograms(const Char_t* dir) {
732 //
733 // saves the histograms in a directory (dir)
734 //
735
736 if (!fHistogramsOn) {
737 AliDebug(0, "Histograms not on - cannot save histograms!!!");
738 return;
739 }
740
741 if (!dir)
742 dir = GetName();
743
744 gDirectory->mkdir(dir);
745 gDirectory->cd(dir);
746
747 gDirectory->mkdir("before_cuts");
748 gDirectory->mkdir("after_cuts");
749
750 fhCutStatistics->Write();
751 fhCutCorrelation->Write();
752
753 for (Int_t i=0; i<2; i++) {
754 if (i==0)
755 gDirectory->cd("before_cuts");
756 else
757 gDirectory->cd("after_cuts");
758
759 fhDcaPosToVertex[i] ->Write();
760 fhDcaNegToVertex[i] ->Write();
761 fhChi2[i] ->Write();
762 fhDcaV0Daughters[i] ->Write();
763 fhRadius[i] ->Write();
764 fhCosinePointingAngle[i] ->Write();
765 fhOnFlyStatus[i] ->Write();
766 fhDcaV0ToVertex[i] ->Write();
767
768 fhPt[i] ->Write();
769
770 gDirectory->cd("../");
771 }
772
773 gDirectory->cd("../");
774}
775
776//____________________________________________________________________
777void AliESDv0Cuts::DrawHistograms()
778{
779 // draws some histograms
780
781 TCanvas* canvas1 = new TCanvas(Form("%s_1", GetName()), "V0 Quality Results1", 800, 800);
782 canvas1->Divide(2, 2);
783
784 canvas1->cd(1);
785 fhDcaPosToVertex[0]->SetStats(kFALSE);
786 fhDcaPosToVertex[0]->Draw();
787
788 canvas1->cd(2);
789 fhChi2[0]->SetStats(kFALSE);
790 fhChi2[0]->Draw();
791
792 canvas1->cd(3);
793 fhDcaV0ToVertex[0]->SetStats(kFALSE);
794 fhDcaV0ToVertex[0]->Draw();
795
796 canvas1->SaveAs(Form("%s_%s.gif", GetName(), canvas1->GetName()));
797
798 TCanvas* canvas2 = new TCanvas(Form("%s_2", GetName()), "V0 Quality Results2", 1200, 800);
799 canvas2->Divide(3, 2);
800
801 canvas2->cd(1);
802 fhDcaV0Daughters[0]->SetStats(kFALSE);
803 gPad->SetLogy();
804 fhDcaV0Daughters[0]->Draw();
805
806 canvas2->cd(2);
807 fhRadius[0]->SetStats(kFALSE);
808 gPad->SetLogy();
809 fhRadius[0]->Draw();
810
811
812 canvas2->cd(4);
813 fhCosinePointingAngle[0]->SetStats(kFALSE);
814 gPad->SetLogy();
815 fhCosinePointingAngle[0]->Draw();
816
817 canvas2->cd(5);
818 fhOnFlyStatus[0]->SetStats(kFALSE);
819 gPad->SetLogy();
820 fhOnFlyStatus[0]->Draw();
821
822 canvas2->SaveAs(Form("%s_%s.gif", GetName(), canvas2->GetName()));
823
824 TCanvas* canvas3 = new TCanvas(Form("%s_4", GetName()), "V0 Quality Results3", 800, 500);
825 canvas3->Divide(2, 1);
826
827 canvas3->cd(1);
828 fhCutStatistics->SetStats(kFALSE);
829 fhCutStatistics->LabelsOption("v");
830 gPad->SetBottomMargin(0.3);
831 fhCutStatistics->Draw();
832
833 canvas3->cd(2);
834 fhCutCorrelation->SetStats(kFALSE);
835 fhCutCorrelation->LabelsOption("v");
836 gPad->SetBottomMargin(0.3);
837 gPad->SetLeftMargin(0.3);
838 fhCutCorrelation->Draw("COLZ");
839
840 canvas3->SaveAs(Form("%s_%s.gif", GetName(), canvas3->GetName()));
841}
842