]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ANALYSIS/AliESDv0Cuts.cxx
Using detector quality flag (taken from ALICE logbook) to decide whether to rpodcue...
[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) {
f690bf48 346// Selection cuts
94b50f5f 347 if(listObj->GetSize()!=4) return kFALSE;
348 AliESDv0 *esdV0 = (AliESDv0*)listObj->At(0);
349 AliESDtrack *trackPos = (AliESDtrack*)listObj->At(1);
350 AliESDtrack *trackNeg = (AliESDtrack*)listObj->At(2);
351 const AliESDVertex *esdVertex = (AliESDVertex*)listObj->At(3);
352 return AcceptV0(esdV0,trackPos,trackNeg,esdVertex);
353}
354
355//____________________________________________________________________
356Bool_t
f690bf48 357AliESDv0Cuts::AcceptV0(AliESDv0* const esdV0, AliESDtrack* const trackPos, AliESDtrack* const trackNeg, const AliESDVertex* esdVertex) {
94b50f5f 358 //
359 // figure out if the v0s survives all the v0 cuts defined
360 //
361 // the different quality parameter and kinematic values are first
362 // retrieved from the v0. then it is found out what cuts the
363 // v0 did not survive and finally the cuts are imposed.
364
365 // this function needs the following branches... but this is not enough
366 // fV0s.fDcaV0Daughters
367 // fV0s.fChi2V0
368 // fV0s.fPos*
369 // fV0s.fNmom*
370 // fV0s.fPmom*
371 // fV0s.fRr
372 // fV0s.fPointAngle
373 // fV0s.fOnFlyStatus
374
375 Float_t dcaPosToVertex = 0, dcaNegToVertex = 0;
dd7faf60 376 Float_t tdcaPosToVertex[2]={999,999};
377 Float_t tdcaNegToVertex[2]={999,999};
94b50f5f 378
379 if (trackPos) trackPos->GetImpactParameters(tdcaPosToVertex[0],tdcaPosToVertex[1]);
380 if (trackNeg) trackNeg->GetImpactParameters(tdcaNegToVertex[0],tdcaNegToVertex[1]);
381
382 dcaPosToVertex = TMath::Sqrt(tdcaPosToVertex[0]*tdcaPosToVertex[0]+tdcaPosToVertex[1]*tdcaPosToVertex[1]);
383 dcaNegToVertex = TMath::Sqrt(tdcaNegToVertex[0]*tdcaNegToVertex[0]+tdcaNegToVertex[1]*tdcaNegToVertex[1]);
384
385 UInt_t status = esdV0->GetOnFlyStatus();
386 Float_t chi2 = esdV0->GetChi2V0();
387
388 Double_t dcaV0Daughters = esdV0->GetDcaV0Daughters();
389
390 Double_t vtxPosition[3]; esdVertex->GetXYZ(vtxPosition);
391 Double_t dcaV0ToVertex = esdV0->GetD(vtxPosition[0],vtxPosition[1],vtxPosition[2]);
392
393 Double_t v0Position[3];
394 esdV0->GetXYZ(v0Position[0],v0Position[1],v0Position[2]);
395 Double_t radius = TMath::Sqrt(TMath::Power(v0Position[0],2) + TMath::Power(v0Position[1],2));
396 Double_t v0CosinePointingAngle = esdV0->GetV0CosineOfPointingAngle();
397
398 // getting the kinematic variables of the v0
399 Double_t p[3];
400 esdV0->GetPxPyPz(p[0],p[1],p[2]);
401 Float_t momentum = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2) + TMath::Power(p[2],2));
402 Float_t pt = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2));
403
404 //########################################################################
405 // cut the v0?
406
407 Bool_t cuts[kNCuts];
408 for (Int_t i=0; i<kNCuts; i++) cuts[i]=kFALSE;
409
410 // v0 quality cuts
411 if (dcaPosToVertex < fCutMinDcaPosToVertex)
412 cuts[0]=kTRUE;
413 if (dcaNegToVertex < fCutMinDcaNegToVertex)
414 cuts[1]=kTRUE;
415 if (chi2 > fCutMaxChi2)
416 cuts[2]=kTRUE;
417 if (dcaV0Daughters > fCutMaxDcaV0Daughters)
418 cuts[3]=kTRUE;
419 if (radius < fCutMinRadius)
420 cuts[4]=kTRUE;
421 if (radius > fCutMaxRadius)
422 cuts[5]=kTRUE;
423 if (v0CosinePointingAngle < fCutMinCosinePointingAngle)
424 cuts[6]=kTRUE;
425 if (fCutRequireOnFlyStatus && !status)
426 cuts[7]=kTRUE;
427 if (dcaV0ToVertex > fCutMaxDcaV0ToVertex)
428 cuts[8] = kTRUE;
429
430 // v0 kinematics cut
431 if((momentum < fPMin) || (momentum > fPMax))
432 cuts[9]=kTRUE;
433 if((pt < fPtMin) || (pt > fPtMax))
434 cuts[10] = kTRUE;
435 if((p[0] < fPxMin) || (p[0] > fPxMax))
436 cuts[11] = kTRUE;
437 if((p[1] < fPyMin) || (p[1] > fPyMax))
438 cuts[12] = kTRUE;
439 if((p[2] < fPzMin) || (p[2] > fPzMax))
440 cuts[13] = kTRUE;
441
442 Bool_t cut=kFALSE;
443 for (Int_t i=0; i<kNCuts; i++)
444 if (cuts[i]) cut = kTRUE;
445
446 //########################################################################
447 // filling histograms
448 if (fHistogramsOn) {
449 fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n v0s")));
450
451 if (cut)
452 fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n cut v0s")));
453
454 for (Int_t i=0; i<kNCuts; i++) {
455 if (cuts[i])
456 fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin(fgkCutNames[i])));
457
458 for (Int_t j=i; j<kNCuts; j++) {
459 if (cuts[i] && cuts[j]) {
460 Float_t x = fhCutCorrelation->GetXaxis()->GetBinCenter(fhCutCorrelation->GetXaxis()->FindBin(fgkCutNames[i]));
461 Float_t y = fhCutCorrelation->GetYaxis()->GetBinCenter(fhCutCorrelation->GetYaxis()->FindBin(fgkCutNames[j]));
462 fhCutCorrelation->Fill(x,y);
463 }
464 }
465 }
466
467 fhDcaPosToVertex[0]->Fill(dcaPosToVertex);
468 fhDcaNegToVertex[0]->Fill(dcaNegToVertex);
469 fhChi2[0]->Fill(chi2);
470 fhDcaV0Daughters[0]->Fill(dcaV0Daughters);
471 fhRadius[0]->Fill(radius);
472 fhCosinePointingAngle[0]->Fill(v0CosinePointingAngle);
473 fhOnFlyStatus[0]->Fill(status);
474 fhDcaV0ToVertex[0]->Fill(dcaV0ToVertex);
475
476 fhPt[0]->Fill(pt);
477 }
478
479 //########################################################################
480 // cut the v0!
481 if (cut) return kFALSE;
482
483 //########################################################################
484 // filling histograms after cut
485 if (fHistogramsOn) {
486 fhDcaPosToVertex[1]->Fill(dcaPosToVertex);
487 fhDcaNegToVertex[1]->Fill(dcaNegToVertex);
488 fhChi2[1]->Fill(chi2);
489 fhDcaV0Daughters[1]->Fill(dcaV0Daughters);
490 fhRadius[1]->Fill(radius);
491 fhCosinePointingAngle[1]->Fill(v0CosinePointingAngle);
492 fhOnFlyStatus[1]->Fill(status);
493 fhDcaV0ToVertex[1]->Fill(dcaV0ToVertex);
494
495 fhPt[1]->Fill(pt);
496 }
497
498 return kTRUE;
499}
500
501//____________________________________________________________________
f690bf48 502TObjArray* AliESDv0Cuts::GetAcceptedV0s(const AliESD* esd)
94b50f5f 503{
504 //
505 // returns an array of all v0s that pass the cuts
506 //
507
508 TObjArray* acceptedV0s = new TObjArray();
509 // const AliESDVertex *spdVertex = esd->GetVertex();
510 const AliESDVertex *primaryVertex = esd->GetPrimaryVertex();
511 Int_t lIndexTrackPos = 0, lIndexTrackNeg = 0;
512
513 // loop over esd v0s
514 for (Int_t iV0 = 0; iV0 < esd->GetNumberOfV0s(); iV0++) {
515 AliESDv0* v0 = esd->GetV0(iV0);
516
517 lIndexTrackPos = TMath::Abs(v0->GetPindex());
518 lIndexTrackNeg = TMath::Abs(v0->GetNindex());
519 AliESDtrack *trackPos = esd->GetTrack(lIndexTrackPos);
520 AliESDtrack *trackNeg = esd->GetTrack(lIndexTrackNeg);
521
522 if (AcceptV0(v0,trackPos,trackNeg,primaryVertex))
523 acceptedV0s->Add(v0);
524 }
525
526 return acceptedV0s;
527}
528
529//____________________________________________________________________
f690bf48 530Int_t AliESDv0Cuts::CountAcceptedV0s(const AliESD* esd)
94b50f5f 531{
532 //
533 // returns an the number of v0s that pass the cuts
534 //
535
536 Int_t count = 0;
537 // const AliESDVertex *spdVertex = esd->GetVertex();
538 const AliESDVertex *primaryVertex = esd->GetPrimaryVertex();
539 Int_t lIndexTrackPos = 0, lIndexTrackNeg = 0;
540
541 // loop over esd v0s
542 for (Int_t iV0 = 0; iV0 < esd->GetNumberOfV0s(); iV0++) {
543 AliESDv0* v0 = esd->GetV0(iV0);
544
545 lIndexTrackPos = TMath::Abs(v0->GetPindex());
546 lIndexTrackNeg = TMath::Abs(v0->GetNindex());
547 AliESDtrack *trackPos = esd->GetTrack(lIndexTrackPos);
548 AliESDtrack *trackNeg = esd->GetTrack(lIndexTrackNeg);
549
550 if (AcceptV0(v0,trackPos,trackNeg,primaryVertex))
551 count++;
552 }
553
554 return count;
555}
556
557//____________________________________________________________________
f690bf48 558TObjArray* AliESDv0Cuts::GetAcceptedV0s(const AliESDEvent* esd)
94b50f5f 559{
560 //
561 // returns an array of all v0s that pass the cuts
562 //
563
564 TObjArray* acceptedV0s = new TObjArray();
565 // const AliESDVertex *spdVertex = esd->GetVertex();
566 const AliESDVertex *primaryVertex = esd->GetPrimaryVertex();
567 Int_t lIndexTrackPos = 0, lIndexTrackNeg = 0;
568
569 // loop over esd v0s
570 for (Int_t iV0 = 0; iV0 < esd->GetNumberOfV0s(); iV0++) {
571 AliESDv0* v0 = esd->GetV0(iV0);
572
573 lIndexTrackPos = TMath::Abs(v0->GetPindex());
574 lIndexTrackNeg = TMath::Abs(v0->GetNindex());
575 AliESDtrack *trackPos = esd->GetTrack(lIndexTrackPos);
576 AliESDtrack *trackNeg = esd->GetTrack(lIndexTrackNeg);
577
578 if (AcceptV0(v0,trackPos,trackNeg,primaryVertex))
579 acceptedV0s->Add(v0);
580 }
581
582 return acceptedV0s;
583}
584
585//____________________________________________________________________
f690bf48 586Int_t AliESDv0Cuts::CountAcceptedV0s(const AliESDEvent* esd)
94b50f5f 587{
588 //
589 // returns an the number of v0s that pass the cuts
590 //
591
592 Int_t count = 0;
593 // const AliESDVertex *spdVertex = esd->GetVertex();
594 const AliESDVertex *primaryVertex = esd->GetPrimaryVertex();
595 Int_t lIndexTrackPos = 0, lIndexTrackNeg = 0;
596
597 // loop over esd v0s
598 for (Int_t iV0 = 0; iV0 < esd->GetNumberOfV0s(); iV0++) {
599 AliESDv0* v0 = esd->GetV0(iV0);
600
601 lIndexTrackPos = TMath::Abs(v0->GetPindex());
602 lIndexTrackNeg = TMath::Abs(v0->GetNindex());
603 AliESDtrack *trackPos = esd->GetTrack(lIndexTrackPos);
604 AliESDtrack *trackNeg = esd->GetTrack(lIndexTrackNeg);
605
606 if (AcceptV0(v0,trackPos,trackNeg,primaryVertex))
607 count++;
608 }
609
610 return count;
611}
612
613//____________________________________________________________________
614 void AliESDv0Cuts::DefineHistograms(Int_t color) {
615 //
616 // diagnostics histograms are defined
617 //
618
619 fHistogramsOn=kTRUE;
620
621 Bool_t oldStatus = TH1::AddDirectoryStatus();
622 TH1::AddDirectory(kFALSE);
623
624 //###################################################################################
625 // defining histograms
626
627 fhCutStatistics = new TH1F("cut_statistics","cut statistics",kNCuts+4,-0.5,kNCuts+3.5);
628
629 fhCutStatistics->GetXaxis()->SetBinLabel(1,"n v0s");
630 fhCutStatistics->GetXaxis()->SetBinLabel(2,"n cut v0s");
631
632 fhCutCorrelation = new TH2F("cut_correlation","cut correlation",kNCuts,-0.5,kNCuts-0.5,kNCuts,-0.5,kNCuts-0.5);;
633
634 for (Int_t i=0; i<kNCuts; i++) {
635 fhCutStatistics->GetXaxis()->SetBinLabel(i+4,fgkCutNames[i]);
636 fhCutCorrelation->GetXaxis()->SetBinLabel(i+1,fgkCutNames[i]);
637 fhCutCorrelation->GetYaxis()->SetBinLabel(i+1,fgkCutNames[i]);
638 }
639
640 fhCutStatistics ->SetLineColor(color);
641 fhCutCorrelation ->SetLineColor(color);
642 fhCutStatistics ->SetLineWidth(2);
643 fhCutCorrelation ->SetLineWidth(2);
644
645 Char_t str[256];
646 for (Int_t i=0; i<2; i++) {
accb0197 647 if (i==0) snprintf(str,256, " ");
648 else snprintf(str,256, "_cut");
94b50f5f 649
650 fhDcaPosToVertex[i] = new TH1F(Form("dcaPosToVertex%s",str),"",120,0,3);
651 fhDcaNegToVertex[i] = new TH1F(Form("dcaNegToVertex%s",str),"",120,0,3);
652 fhChi2[i] = new TH1F(Form("chi2%s",str),"",50,0,50);
653 fhDcaV0Daughters[i] = new TH1F(Form("dcaV0Daughters%s",str),"",200,0,5);
654 fhRadius[i] = new TH1F(Form("decayRadius%s",str),"",300,0,150);
655 fhCosinePointingAngle[i] = new TH1F(Form("cosinePointingAngle%s",str),"",100,-1,1);
656 fhOnFlyStatus[i] = new TH1F(Form("onflystatus%s",str),"",5,0,5);
657 fhDcaV0ToVertex[i] = new TH1F(Form("dcaV0ToVertex%s",str),"",100,0,5);
658
659 fhPt[i] = new TH1F(Form("pt%s",str) ,"p_{T} distribution;p_{T} (GeV/c)",500,0.0,100.0);
660
661 fhDcaPosToVertex[i]->SetTitle("Dca of positive daughter to parent vertex");
662 fhDcaNegToVertex[i]->SetTitle("Dca of negative daughter to parent vertex");
663 fhChi2[i]->SetTitle("#Chi^{2} for v0");
664 fhDcaV0Daughters[i]->SetTitle("Dca between daughter tracks");
665 fhRadius[i]->SetTitle("Decay radius of the v0");
666 fhCosinePointingAngle[i]->SetTitle("Cosine of the Pointing Angle");
667 fhOnFlyStatus[i]->SetTitle("On-the-Fly Status");
668 fhDcaV0ToVertex[i]->SetTitle("Dca of v0 to parent vertex");
669
670 fhDcaPosToVertex[i]->SetLineColor(color); fhDcaPosToVertex[i]->SetLineWidth(2);
671 fhDcaNegToVertex[i]->SetLineColor(color); fhDcaNegToVertex[i]->SetLineWidth(2);
672 fhChi2[i]->SetLineColor(color); fhChi2[i]->SetLineWidth(2);
673 fhDcaV0Daughters[i]->SetLineColor(color); fhDcaV0Daughters[i]->SetLineWidth(2);
674 fhRadius[i]->SetLineColor(color); fhRadius[i]->SetLineWidth(2);
675 fhCosinePointingAngle[i]->SetLineColor(color); fhCosinePointingAngle[i]->SetLineWidth(2);
676 fhOnFlyStatus[i]->SetLineColor(color); fhOnFlyStatus[i]->SetLineWidth(2);
677 fhDcaV0ToVertex[i]->SetLineColor(color); fhDcaV0ToVertex[i]->SetLineWidth(2);
678 }
679
680 TH1::AddDirectory(oldStatus);
681}
682
683//____________________________________________________________________
684Bool_t AliESDv0Cuts::LoadHistograms(const Char_t* dir)
685{
686 //
687 // loads the histograms from a file
688 // if dir is empty a directory with the name of this object is taken (like in SaveHistogram)
689 //
690
691 if (!dir)
692 dir = GetName();
693
694 if (!gDirectory->cd(dir))
695 return kFALSE;
696
697 fhCutStatistics = dynamic_cast<TH1F*> (gDirectory->Get("cut_statistics"));
698 fhCutCorrelation = dynamic_cast<TH2F*> (gDirectory->Get("cut_correlation"));
699
700 Char_t str[5];
701 for (Int_t i=0; i<2; i++) {
702 if (i==0)
703 {
704 gDirectory->cd("before_cuts");
705 str[0] = 0;
706 }
707 else
708 {
709 gDirectory->cd("after_cuts");
accb0197 710 snprintf(str,5, "_cut");
94b50f5f 711 }
712
713 fhDcaPosToVertex[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("dcaPosToVertex%s",str) ));
714 fhDcaNegToVertex[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("dcaNegToVertex%s",str) ));
715 fhChi2[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("chi2%s",str)));
716 fhDcaV0Daughters[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("dcaV0Daughters%s",str)));
717 fhRadius[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("decayRadius%s",str)));
718 fhCosinePointingAngle[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("cosinepointingangle%s",str)));
719 fhOnFlyStatus[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("onflystatus%s",str)));
720 fhDcaV0ToVertex[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("dcaV0ToVertex%s",str)));
721
722 fhPt[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("pt%s",str)));
723
724 gDirectory->cd("../");
725 }
726
727 gDirectory->cd("..");
728
729 return kTRUE;
730}
731
732//____________________________________________________________________
733void AliESDv0Cuts::SaveHistograms(const Char_t* dir) {
734 //
735 // saves the histograms in a directory (dir)
736 //
737
738 if (!fHistogramsOn) {
739 AliDebug(0, "Histograms not on - cannot save histograms!!!");
740 return;
741 }
742
743 if (!dir)
744 dir = GetName();
745
746 gDirectory->mkdir(dir);
747 gDirectory->cd(dir);
748
749 gDirectory->mkdir("before_cuts");
750 gDirectory->mkdir("after_cuts");
751
752 fhCutStatistics->Write();
753 fhCutCorrelation->Write();
754
755 for (Int_t i=0; i<2; i++) {
756 if (i==0)
757 gDirectory->cd("before_cuts");
758 else
759 gDirectory->cd("after_cuts");
760
761 fhDcaPosToVertex[i] ->Write();
762 fhDcaNegToVertex[i] ->Write();
763 fhChi2[i] ->Write();
764 fhDcaV0Daughters[i] ->Write();
765 fhRadius[i] ->Write();
766 fhCosinePointingAngle[i] ->Write();
767 fhOnFlyStatus[i] ->Write();
768 fhDcaV0ToVertex[i] ->Write();
769
770 fhPt[i] ->Write();
771
772 gDirectory->cd("../");
773 }
774
775 gDirectory->cd("../");
776}
777
778//____________________________________________________________________
779void AliESDv0Cuts::DrawHistograms()
780{
781 // draws some histograms
782
783 TCanvas* canvas1 = new TCanvas(Form("%s_1", GetName()), "V0 Quality Results1", 800, 800);
784 canvas1->Divide(2, 2);
785
786 canvas1->cd(1);
787 fhDcaPosToVertex[0]->SetStats(kFALSE);
788 fhDcaPosToVertex[0]->Draw();
789
790 canvas1->cd(2);
791 fhChi2[0]->SetStats(kFALSE);
792 fhChi2[0]->Draw();
793
794 canvas1->cd(3);
795 fhDcaV0ToVertex[0]->SetStats(kFALSE);
796 fhDcaV0ToVertex[0]->Draw();
797
798 canvas1->SaveAs(Form("%s_%s.gif", GetName(), canvas1->GetName()));
799
800 TCanvas* canvas2 = new TCanvas(Form("%s_2", GetName()), "V0 Quality Results2", 1200, 800);
801 canvas2->Divide(3, 2);
802
803 canvas2->cd(1);
804 fhDcaV0Daughters[0]->SetStats(kFALSE);
805 gPad->SetLogy();
806 fhDcaV0Daughters[0]->Draw();
807
808 canvas2->cd(2);
809 fhRadius[0]->SetStats(kFALSE);
810 gPad->SetLogy();
811 fhRadius[0]->Draw();
812
813
814 canvas2->cd(4);
815 fhCosinePointingAngle[0]->SetStats(kFALSE);
816 gPad->SetLogy();
817 fhCosinePointingAngle[0]->Draw();
818
819 canvas2->cd(5);
820 fhOnFlyStatus[0]->SetStats(kFALSE);
821 gPad->SetLogy();
822 fhOnFlyStatus[0]->Draw();
823
824 canvas2->SaveAs(Form("%s_%s.gif", GetName(), canvas2->GetName()));
825
826 TCanvas* canvas3 = new TCanvas(Form("%s_4", GetName()), "V0 Quality Results3", 800, 500);
827 canvas3->Divide(2, 1);
828
829 canvas3->cd(1);
830 fhCutStatistics->SetStats(kFALSE);
831 fhCutStatistics->LabelsOption("v");
832 gPad->SetBottomMargin(0.3);
833 fhCutStatistics->Draw();
834
835 canvas3->cd(2);
836 fhCutCorrelation->SetStats(kFALSE);
837 fhCutCorrelation->LabelsOption("v");
838 gPad->SetBottomMargin(0.3);
839 gPad->SetLeftMargin(0.3);
840 fhCutCorrelation->Draw("COLZ");
841
842 canvas3->SaveAs(Form("%s_%s.gif", GetName(), canvas3->GetName()));
843}
844