1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 // The class AliCFTrackIsPrimaryCut is designed to select reconstructed tracks
17 // with a small impact parameter and tracks which are (not) daughters of kink
18 // decays and to provide corresponding QA histograms.
19 // This class inherits from the Analysis' Framework abstract base class
20 // AliAnalysisCuts and is a part of the Correction Framework.
21 // This class acts on single, reconstructed tracks, it is applicable on
23 // It mainly consists of a IsSelected function that returns a boolean.
24 // This function checks whether the considered track passes a set of cuts:
25 // - min. and max. distance to main vertex in transverse plane (xy)
26 // - min. and max. longitudinal distance to main vertex (z)
27 // - min. and max. distance to main vertex as ellpise in xy - z plane
28 // - all above cuts on absolute values or in units of sigma (resolution)
29 // - min. and max. distance to main vertex in units of sigma (resolution)
30 // - max. transverse (xy) and longitudinal (z) impact parameter resolution
31 // - require that the dca calculation doesn't fail
32 // - accept or not accept daughter tracks of kink decays
34 // The cut values for these cuts are set with the corresponding set functions.
35 // All cut classes provided by the correction framework are supposed to be
36 // added in the Analysis Framwork's class AliAnalysisFilter and applied by
37 // the filter via a loop.
39 // author: I. Kraus (Ingrid.Kraus@cern.ch)
41 // AliESDtrackCuts writte by Jan Fiete Grosse-Oetringhaus and
42 // AliRsnDaughterCut class written by A. Pulvirenti.
45 #include <TDirectory.h>
49 #include <AliESDtrack.h>
51 #include "AliCFTrackIsPrimaryCuts.h"
53 ClassImp(AliCFTrackIsPrimaryCuts)
55 //__________________________________________________________________________________
56 AliCFTrackIsPrimaryCuts::AliCFTrackIsPrimaryCuts() :
64 fNSigmaToVertexMin(0),
65 fNSigmaToVertexMax(0),
68 fRequireSigmaToVertex(0),
69 fAODType(AliAODTrack::kUndef),
70 fAcceptKinkDaughters(0),
75 fhNBinsRequireSigma(0),
84 fhBinLimRequireSigma(0x0),
85 fhBinLimAcceptKink(0x0),
88 fhBinLimDcaXYnorm(0x0),
89 fhBinLimDcaZnorm(0x0),
90 fhBinLimSigmaDcaXY(0x0),
91 fhBinLimSigmaDcaZ(0x0)
94 // Default constructor
98 //__________________________________________________________________________________
99 AliCFTrackIsPrimaryCuts::AliCFTrackIsPrimaryCuts(Char_t* name, Char_t* title) :
100 AliCFCutBase(name,title),
101 fMinDCAToVertexXY(0),
103 fMaxDCAToVertexXY(0),
107 fNSigmaToVertexMin(0),
108 fNSigmaToVertexMax(0),
111 fRequireSigmaToVertex(0),
112 fAODType(AliAODTrack::kUndef),
113 fAcceptKinkDaughters(0),
118 fhNBinsRequireSigma(0),
119 fhNBinsAcceptKink(0),
124 fhNBinsSigmaDcaXY(0),
127 fhBinLimRequireSigma(0x0),
128 fhBinLimAcceptKink(0x0),
131 fhBinLimDcaXYnorm(0x0),
132 fhBinLimDcaZnorm(0x0),
133 fhBinLimSigmaDcaXY(0x0),
134 fhBinLimSigmaDcaZ(0x0)
141 //__________________________________________________________________________________
142 AliCFTrackIsPrimaryCuts::AliCFTrackIsPrimaryCuts(const AliCFTrackIsPrimaryCuts& c) :
144 fMinDCAToVertexXY(c.fMinDCAToVertexXY),
145 fMinDCAToVertexZ(c.fMinDCAToVertexZ),
146 fMaxDCAToVertexXY(c.fMaxDCAToVertexXY),
147 fMaxDCAToVertexZ(c.fMaxDCAToVertexZ),
148 fDCAToVertex2D(c.fDCAToVertex2D),
149 fAbsDCAToVertex(c.fAbsDCAToVertex),
150 fNSigmaToVertexMin(c.fNSigmaToVertexMin),
151 fNSigmaToVertexMax(c.fNSigmaToVertexMax),
152 fSigmaDCAxy(c.fSigmaDCAxy),
153 fSigmaDCAz(c.fSigmaDCAz),
154 fRequireSigmaToVertex(c.fRequireSigmaToVertex),
155 fAODType(c.fAODType),
156 fAcceptKinkDaughters(c.fAcceptKinkDaughters),
157 fhCutStatistics(c.fhCutStatistics),
158 fhCutCorrelation(c.fhCutCorrelation),
160 fhNBinsNSigma(c.fhNBinsNSigma),
161 fhNBinsRequireSigma(c.fhNBinsRequireSigma),
162 fhNBinsAcceptKink(c.fhNBinsAcceptKink),
163 fhNBinsDcaXY(c.fhNBinsDcaXY),
164 fhNBinsDcaZ(c.fhNBinsDcaZ),
165 fhNBinsDcaXYnorm(c.fhNBinsDcaXYnorm),
166 fhNBinsDcaZnorm(c.fhNBinsDcaZnorm),
167 fhNBinsSigmaDcaXY(c.fhNBinsSigmaDcaXY),
168 fhNBinsSigmaDcaZ(c.fhNBinsSigmaDcaZ),
169 fhBinLimNSigma(c.fhBinLimNSigma),
170 fhBinLimRequireSigma(c.fhBinLimRequireSigma),
171 fhBinLimAcceptKink(c.fhBinLimAcceptKink),
172 fhBinLimDcaXY(c.fhBinLimDcaXY),
173 fhBinLimDcaZ(c.fhBinLimDcaZ),
174 fhBinLimDcaXYnorm(c.fhBinLimDcaXYnorm),
175 fhBinLimDcaZnorm(c.fhBinLimDcaZnorm),
176 fhBinLimSigmaDcaXY(c.fhBinLimSigmaDcaXY),
177 fhBinLimSigmaDcaZ(c.fhBinLimSigmaDcaZ)
182 ((AliCFTrackIsPrimaryCuts &) c).Copy(*this);
184 //__________________________________________________________________________________
185 AliCFTrackIsPrimaryCuts& AliCFTrackIsPrimaryCuts::operator=(const AliCFTrackIsPrimaryCuts& c)
188 // Assignment operator
191 AliCFCutBase::operator=(c) ;
192 fMinDCAToVertexXY = c.fMinDCAToVertexXY;
193 fMinDCAToVertexZ = c.fMinDCAToVertexZ;
194 fMaxDCAToVertexXY = c.fMaxDCAToVertexXY;
195 fMaxDCAToVertexZ = c.fMaxDCAToVertexZ;
196 fDCAToVertex2D = c.fDCAToVertex2D;
197 fAbsDCAToVertex = c.fAbsDCAToVertex;
198 fNSigmaToVertexMin = c.fNSigmaToVertexMin ;
199 fNSigmaToVertexMax = c.fNSigmaToVertexMax ;
200 fSigmaDCAxy = c.fSigmaDCAxy ;
201 fSigmaDCAz = c.fSigmaDCAz ;
202 fRequireSigmaToVertex = c.fRequireSigmaToVertex ;
203 fAODType = c.fAODType ;
204 fAcceptKinkDaughters = c.fAcceptKinkDaughters ;
205 fhCutStatistics = c.fhCutStatistics ;
206 fhCutCorrelation = c.fhCutCorrelation ;
208 fhNBinsNSigma = c.fhNBinsNSigma;
209 fhNBinsRequireSigma = c.fhNBinsRequireSigma;
210 fhNBinsAcceptKink = c.fhNBinsAcceptKink;
211 fhNBinsDcaXY = c.fhNBinsDcaXY;
212 fhNBinsDcaZ = c.fhNBinsDcaZ;
213 fhNBinsDcaXYnorm = c.fhNBinsDcaXYnorm;
214 fhNBinsDcaZnorm = c.fhNBinsDcaZnorm;
215 fhNBinsSigmaDcaXY = c.fhNBinsSigmaDcaXY;
216 fhNBinsSigmaDcaZ = c.fhNBinsSigmaDcaZ;
217 fhBinLimNSigma = c.fhBinLimNSigma;
218 fhBinLimRequireSigma = c.fhBinLimRequireSigma;
219 fhBinLimAcceptKink = c.fhBinLimAcceptKink;
220 fhBinLimDcaXY = c.fhBinLimDcaXY;
221 fhBinLimDcaZ = c.fhBinLimDcaZ;
222 fhBinLimDcaXYnorm = c.fhBinLimDcaXYnorm;
223 fhBinLimDcaZnorm = c.fhBinLimDcaZnorm;
224 fhBinLimSigmaDcaXY = c.fhBinLimSigmaDcaXY;
225 fhBinLimSigmaDcaZ = c.fhBinLimSigmaDcaZ;
227 for (Int_t j=0; j<c.kNStepQA; j++){
228 if(c.fhDcaXYvsDcaZ[j]) fhDcaXYvsDcaZ[j] = (TH2F*)c.fhDcaXYvsDcaZ[j]->Clone();
229 if(c.fhDcaXYvsDcaZnorm[j]) fhDcaXYvsDcaZnorm[j] = (TH2F*)c.fhDcaXYvsDcaZnorm[j]->Clone();
230 for (Int_t i=0; i<c.kNHist; i++){
231 if(c.fhQA[i][j]) fhQA[i][j] = (TH1F*)c.fhQA[i][j]->Clone();
234 ((AliCFTrackIsPrimaryCuts &) c).Copy(*this);
238 //__________________________________________________________________________________
239 AliCFTrackIsPrimaryCuts::~AliCFTrackIsPrimaryCuts()
244 if (fhCutStatistics) delete fhCutStatistics;
245 if (fhCutCorrelation) delete fhCutCorrelation;
247 for (Int_t j=0; j<kNStepQA; j++){
248 if(fhDcaXYvsDcaZ[j]) delete fhDcaXYvsDcaZ[j];
249 if(fhDcaXYvsDcaZnorm[j]) delete fhDcaXYvsDcaZnorm[j];
250 for (Int_t i=0; i<kNHist; i++)
251 if(fhQA[i][j]) delete fhQA[i][j];
253 if(fBitmap) delete fBitmap;
254 if(fhBinLimNSigma) delete fhBinLimNSigma;
255 if(fhBinLimRequireSigma) delete fhBinLimRequireSigma;
256 if(fhBinLimAcceptKink) delete fhBinLimAcceptKink;
257 if(fhBinLimDcaXY) delete fhBinLimDcaXY;
258 if(fhBinLimDcaZ) delete fhBinLimDcaZ;
259 if(fhBinLimDcaXYnorm) delete fhBinLimDcaXYnorm;
260 if(fhBinLimDcaZnorm) delete fhBinLimDcaZnorm;
261 if(fhBinLimSigmaDcaXY) delete fhBinLimSigmaDcaXY;
262 if(fhBinLimSigmaDcaZ) delete fhBinLimSigmaDcaZ;
264 //__________________________________________________________________________________
265 void AliCFTrackIsPrimaryCuts::Initialise()
268 // sets everything to zero
270 fMinDCAToVertexXY = 0;
271 fMinDCAToVertexZ = 0;
272 fMaxDCAToVertexXY = 0;
273 fMaxDCAToVertexZ = 0;
276 fNSigmaToVertexMin = 0;
277 fNSigmaToVertexMax = 0;
280 fRequireSigmaToVertex = 0;
281 fAcceptKinkDaughters = 0;
282 fAODType = AliAODTrack::kUndef;
284 SetMinDCAToVertexXY();
285 SetMinDCAToVertexZ();
286 SetMaxDCAToVertexXY();
287 SetMaxDCAToVertexZ();
290 SetMinNSigmaToVertex();
291 SetMaxNSigmaToVertex();
294 SetRequireSigmaToVertex();
295 SetAcceptKinkDaughters();
298 for (Int_t j=0; j<kNStepQA; j++) {
299 fhDcaXYvsDcaZ[j] = 0x0;
300 fhDcaXYvsDcaZnorm[j] = 0x0;
301 for (Int_t i=0; i<kNHist; i++)
305 fhCutCorrelation = 0;
306 fBitmap=new TBits(0);
308 //set default bining for QA histograms
309 SetHistogramBins(kCutNSigmaToVertex,500,0.,50.);
310 SetHistogramBins(kCutRequireSigmaToVertex,5,-0.75,1.75);
311 SetHistogramBins(kCutAcceptKinkDaughters,5,-0.75,1.75);
312 SetHistogramBins(kDcaXY,500,-10.,10.);
313 SetHistogramBins(kDcaZ,500,-10.,10.);
314 SetHistogramBins(kDcaXYnorm,500,-10.,10.);
315 SetHistogramBins(kDcaZnorm,500,-10.,10.);
316 SetHistogramBins(kSigmaDcaXY,500,-0.1,0.9);
317 SetHistogramBins(kSigmaDcaZ,500,-0.1,0.9);
319 //__________________________________________________________________________________
320 void AliCFTrackIsPrimaryCuts::Copy(TObject &c) const
325 AliCFTrackIsPrimaryCuts& target = (AliCFTrackIsPrimaryCuts &) c;
329 if (fhCutStatistics) target.fhCutStatistics = (TH1F*) fhCutStatistics->Clone();
330 if (fhCutCorrelation) target.fhCutCorrelation = (TH2F*) fhCutCorrelation->Clone();
332 for (Int_t j=0; j<kNStepQA; j++){
333 if(fhDcaXYvsDcaZ[j]) target.fhDcaXYvsDcaZ[j] = (TH2F*)fhDcaXYvsDcaZ[j]->Clone();
334 if(fhDcaXYvsDcaZnorm[j]) target.fhDcaXYvsDcaZnorm[j] = (TH2F*)fhDcaXYvsDcaZnorm[j]->Clone();
335 for (Int_t i=0; i<kNHist; i++)
336 if(fhQA[i][j]) target.fhQA[i][j] = (TH1F*)fhQA[i][j]->Clone();
340 //__________________________________________________________________________________
341 void AliCFTrackIsPrimaryCuts::SelectionBitMap(TObject* obj)
344 // test if the track passes the single cuts
345 // and store the information in a bitmap
348 // bitmap stores the decision of each single cut
349 for(Int_t i=0; i<kNCuts; i++)fBitmap->SetBitNumber(i,kFALSE);
351 // check TObject and cast into ESDtrack
353 if (!obj->InheritsFrom("AliVParticle")) {
354 AliError("object must derived from AliVParticle !");
358 Bool_t isESDTrack = strcmp(obj->ClassName(),"AliESDtrack") == 0 ? kTRUE : kFALSE ;
359 Bool_t isAODTrack = strcmp(obj->ClassName(),"AliAODTrack") == 0 ? kTRUE : kFALSE ;
361 AliESDtrack * esdTrack = 0x0 ;
362 AliAODTrack * aodTrack = 0x0 ;
363 if (isESDTrack) esdTrack = dynamic_cast<AliESDtrack*>(obj);
364 if (isAODTrack) aodTrack = dynamic_cast<AliAODTrack*>(obj);
366 // get the track to vertex parameter for ESD track
367 Float_t nSigmaToVertex = 0.;
368 if (isESDTrack) nSigmaToVertex = AliESDtrackCuts::GetSigmaToVertex(esdTrack);
373 Float_t b2Dmin = 0, b2Dmax = 0;
375 esdTrack->GetImpactParameters(b,bCov);
376 if (bCov[0]<=0 || bCov[2]<=0) {
377 // // // AliDebugClass(1, "Estimated b resolution lower or equal zero!");
378 bCov[0]=0; bCov[2]=0;
380 b[0] = TMath::Abs(b[0]);
381 b[1] = TMath::Abs(b[1]);
382 bRes[0] = TMath::Sqrt(bCov[0]);
383 bRes[1] = TMath::Sqrt(bCov[2]);
384 if (!fAbsDCAToVertex) {
385 if (bRes[0] > 0) b[0] = b[0]/bRes[0];
386 // // // else AliDebugClass(1, "Estimated b resolution equal zero!");
387 if (bRes[1] > 0) b[1] = b[1]/bRes[1];
388 // // // else AliDebugClass(1, "Estimated b resolution equal zero!");
390 if (fMinDCAToVertexXY>0 && fMinDCAToVertexZ>0)
391 b2Dmin = b[0]*b[0]/fMinDCAToVertexXY/fMinDCAToVertexXY + b[1]*b[1]/fMinDCAToVertexZ/fMinDCAToVertexZ;
392 if (fMaxDCAToVertexXY>0 && fMaxDCAToVertexZ>0)
393 b2Dmax = b[0]*b[0]/fMaxDCAToVertexXY/fMaxDCAToVertexXY + b[1]*b[1]/fMaxDCAToVertexZ/fMaxDCAToVertexZ;
400 if (fDCAToVertex2D || (!fDCAToVertex2D && b[0] >= fMinDCAToVertexXY && b[0] <= fMaxDCAToVertexXY)) {
401 fBitmap->SetBitNumber(iCutBit,kTRUE);
404 else fBitmap->SetBitNumber(iCutBit,kTRUE);
409 if (fDCAToVertex2D || (!fDCAToVertex2D && b[1] >= fMinDCAToVertexZ && b[1] <= fMaxDCAToVertexZ)) {
410 fBitmap->SetBitNumber(iCutBit,kTRUE);
413 else fBitmap->SetBitNumber(iCutBit,kTRUE);
418 if (!fDCAToVertex2D || (fDCAToVertex2D && b2Dmin > 1 && b2Dmax < 1)) {
419 fBitmap->SetBitNumber(iCutBit,kTRUE);
422 else fBitmap->SetBitNumber(iCutBit,kTRUE);
427 if (nSigmaToVertex >= fNSigmaToVertexMin && nSigmaToVertex <= fNSigmaToVertexMax) {
428 fBitmap->SetBitNumber(iCutBit,kTRUE);
431 else fBitmap->SetBitNumber(iCutBit,kTRUE);
436 if (bRes[0] < fSigmaDCAxy) {
437 fBitmap->SetBitNumber(iCutBit,kTRUE);
440 else fBitmap->SetBitNumber(iCutBit,kTRUE);
445 if (bRes[1] < fSigmaDCAz) {
446 fBitmap->SetBitNumber(iCutBit,kTRUE);
449 else fBitmap->SetBitNumber(iCutBit,kTRUE);
454 if (!fRequireSigmaToVertex || (nSigmaToVertex>=0 && fRequireSigmaToVertex)) {
455 fBitmap->SetBitNumber(iCutBit,kTRUE);
458 else fBitmap->SetBitNumber(iCutBit,kTRUE);
463 if (fAcceptKinkDaughters || (!fAcceptKinkDaughters && esdTrack->GetKinkIndex(0)<=0)) {
464 fBitmap->SetBitNumber(iCutBit,kTRUE);
467 else fBitmap->SetBitNumber(iCutBit,kTRUE);
472 if (fAODType==AliAODTrack::kUndef || fAODType == aodTrack->GetType()) {
473 fBitmap->SetBitNumber(iCutBit,kTRUE);
476 else fBitmap->SetBitNumber(iCutBit,kTRUE);
480 //__________________________________________________________________________________
481 Bool_t AliCFTrackIsPrimaryCuts::IsSelected(TObject* obj) {
483 // loops over decisions of single cuts and returns if the track is accepted
485 SelectionBitMap(obj);
487 if (fIsQAOn) FillHistograms(obj,0);
488 Bool_t isSelected = kTRUE;
490 for (UInt_t icut=0; icut<fBitmap->GetNbits();icut++) {
491 if(!fBitmap->TestBitNumber(icut)) {
496 if (!isSelected) return kFALSE ;
497 if (fIsQAOn) FillHistograms(obj,1);
500 //__________________________________________________________________________________
501 void AliCFTrackIsPrimaryCuts::SetHistogramBins(Int_t index, Int_t nbins, Double_t *bins)
509 case kCutNSigmaToVertex:
510 fhNBinsNSigma=nbins+1;
511 fhBinLimNSigma=new Double_t[nbins+1];
512 for(Int_t i=0;i<nbins+1;i++)fhBinLimNSigma[i]=bins[i];
515 case kCutRequireSigmaToVertex:
516 fhNBinsRequireSigma=nbins+1;
517 fhBinLimRequireSigma=new Double_t[nbins+1];
518 for(Int_t i=0;i<nbins+1;i++)fhBinLimRequireSigma[i]=bins[i];
521 case kCutAcceptKinkDaughters:
522 fhNBinsAcceptKink=nbins+1;
523 fhBinLimAcceptKink=new Double_t[nbins+1];
524 for(Int_t i=0;i<nbins+1;i++)fhBinLimAcceptKink[i]=bins[i];
528 fhNBinsDcaXY=nbins+1;
529 fhBinLimDcaXY=new Double_t[nbins+1];
530 for(Int_t i=0;i<nbins+1;i++)fhBinLimDcaXY[i]=bins[i];
535 fhBinLimDcaZ=new Double_t[nbins+1];
536 for(Int_t i=0;i<nbins+1;i++)fhBinLimDcaZ[i]=bins[i];
540 fhNBinsDcaXYnorm=nbins+1;
541 fhBinLimDcaXYnorm=new Double_t[nbins+1];
542 for(Int_t i=0;i<nbins+1;i++)fhBinLimDcaXYnorm[i]=bins[i];
546 fhNBinsDcaZnorm=nbins+1;
547 fhBinLimDcaZnorm=new Double_t[nbins+1];
548 for(Int_t i=0;i<nbins+1;i++)fhBinLimDcaZnorm[i]=bins[i];
552 fhNBinsSigmaDcaXY=nbins+1;
553 fhBinLimSigmaDcaXY=new Double_t[nbins+1];
554 for(Int_t i=0;i<nbins+1;i++)fhBinLimSigmaDcaXY[i]=bins[i];
558 fhNBinsSigmaDcaZ=nbins+1;
559 fhBinLimSigmaDcaZ=new Double_t[nbins+1];
560 for(Int_t i=0;i<nbins+1;i++)fhBinLimSigmaDcaZ[i]=bins[i];
564 //__________________________________________________________________________________
565 void AliCFTrackIsPrimaryCuts::SetHistogramBins(Int_t index, Int_t nbins, Double_t xmin, Double_t xmax)
571 case kCutNSigmaToVertex:
572 fhNBinsNSigma=nbins+1;
573 fhBinLimNSigma=new Double_t[nbins+1];
574 for(Int_t i=0;i<nbins+1;i++)fhBinLimNSigma[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
577 case kCutRequireSigmaToVertex:
578 fhNBinsRequireSigma=nbins+1;
579 fhBinLimRequireSigma=new Double_t[nbins+1];
580 for(Int_t i=0;i<nbins+1;i++)fhBinLimRequireSigma[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
583 case kCutAcceptKinkDaughters:
584 fhNBinsAcceptKink=nbins+1;
585 fhBinLimAcceptKink=new Double_t[nbins+1];
586 for(Int_t i=0;i<nbins+1;i++)fhBinLimAcceptKink[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
590 fhNBinsDcaXY=nbins+1;
591 fhBinLimDcaXY=new Double_t[nbins+1];
592 for(Int_t i=0;i<nbins+1;i++)fhBinLimDcaXY[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
597 fhBinLimDcaZ=new Double_t[nbins+1];
598 for(Int_t i=0;i<nbins+1;i++)fhBinLimDcaZ[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
602 fhNBinsDcaXYnorm=nbins+1;
603 fhBinLimDcaXYnorm=new Double_t[nbins+1];
604 for(Int_t i=0;i<nbins+1;i++)fhBinLimDcaXYnorm[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
608 fhNBinsDcaZnorm=nbins+1;
609 fhBinLimDcaZnorm=new Double_t[nbins+1];
610 for(Int_t i=0;i<nbins+1;i++)fhBinLimDcaZnorm[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
614 fhNBinsSigmaDcaXY=nbins+1;
615 fhBinLimSigmaDcaXY=new Double_t[nbins+1];
616 for(Int_t i=0;i<nbins+1;i++)fhBinLimSigmaDcaXY[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
620 fhNBinsSigmaDcaZ=nbins+1;
621 fhBinLimSigmaDcaZ=new Double_t[nbins+1];
622 for(Int_t i=0;i<nbins+1;i++)fhBinLimSigmaDcaZ[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
626 //__________________________________________________________________________________
627 void AliCFTrackIsPrimaryCuts::DefineHistograms() {
629 // histograms for cut variables, cut statistics and cut correlations
633 // book cut statistics and cut correlation histograms
634 fhCutStatistics = new TH1F(Form("%s_cut_statistics",GetName()), Form("%s cut statistics",GetName()), kNCuts,0.5,kNCuts+0.5);
635 fhCutStatistics->SetLineWidth(2);
636 fhCutStatistics->GetXaxis()->SetBinLabel(1,"dca xy");
637 fhCutStatistics->GetXaxis()->SetBinLabel(2,"dca z");
638 fhCutStatistics->GetXaxis()->SetBinLabel(3,"dca ellipse");
639 fhCutStatistics->GetXaxis()->SetBinLabel(4,"n dca");
640 fhCutStatistics->GetXaxis()->SetBinLabel(5,"sigma dca xy");
641 fhCutStatistics->GetXaxis()->SetBinLabel(6,"sigma dca z");
642 fhCutStatistics->GetXaxis()->SetBinLabel(7,"require dca");
643 fhCutStatistics->GetXaxis()->SetBinLabel(8,"kink daughter");
644 fhCutStatistics->GetXaxis()->SetBinLabel(9,"AOD type");
646 fhCutCorrelation = new TH2F(Form("%s_cut_correlation",GetName()), Form("%s cut correlation",GetName()), kNCuts,0.5,kNCuts+0.5,kNCuts,0.5,kNCuts+0.5);
647 fhCutCorrelation->SetLineWidth(2);
648 fhCutCorrelation->GetXaxis()->SetBinLabel(1,"dca xy");
649 fhCutCorrelation->GetXaxis()->SetBinLabel(2,"dca z");
650 fhCutCorrelation->GetXaxis()->SetBinLabel(3,"dca ellipse");
651 fhCutCorrelation->GetXaxis()->SetBinLabel(4,"n dca");
652 fhCutCorrelation->GetXaxis()->SetBinLabel(5,"sigma dca xy");
653 fhCutCorrelation->GetXaxis()->SetBinLabel(6,"sigma dca z");
654 fhCutCorrelation->GetXaxis()->SetBinLabel(7,"require dca");
655 fhCutCorrelation->GetXaxis()->SetBinLabel(8,"kink daughter");
656 fhCutCorrelation->GetXaxis()->SetBinLabel(9,"AOD type");
658 fhCutCorrelation->GetYaxis()->SetBinLabel(1,"dca xy");
659 fhCutCorrelation->GetYaxis()->SetBinLabel(2,"dca z");
660 fhCutCorrelation->GetYaxis()->SetBinLabel(3,"dca ellipse");
661 fhCutCorrelation->GetYaxis()->SetBinLabel(4,"n dca");
662 fhCutCorrelation->GetYaxis()->SetBinLabel(5,"sigma dca xy");
663 fhCutCorrelation->GetYaxis()->SetBinLabel(6,"sigma dca z");
664 fhCutCorrelation->GetYaxis()->SetBinLabel(7,"require dca");
665 fhCutCorrelation->GetYaxis()->SetBinLabel(8,"kink daughter");
666 fhCutCorrelation->GetYaxis()->SetBinLabel(9,"AOD type");
668 // book QA histograms
670 for (Int_t i=0; i<kNStepQA; i++) {
671 if (i==0) sprintf(str," ");
672 else sprintf(str,"_cut");
674 fhDcaXYvsDcaZ[i] = new TH2F(Form("%s_dcaXYvsDcaZ%s",GetName(),str),"",200,-10,10,200,-10,10);
675 fhDcaXYvsDcaZnorm[i] = new TH2F(Form("%s_dcaXYvsDcaZnorm%s",GetName(),str),"",200,-10,10,200,-10,10);
676 fhQA[kCutNSigmaToVertex][i] = new TH1F(Form("%s_nSigmaToVertex%s",GetName(),str),"",fhNBinsNSigma-1,fhBinLimNSigma);
677 fhQA[kCutRequireSigmaToVertex][i] = new TH1F(Form("%s_requireSigmaToVertex%s",GetName(),str),"",fhNBinsRequireSigma-1,fhBinLimRequireSigma);
678 fhQA[kCutAcceptKinkDaughters][i] = new TH1F(Form("%s_acceptKinkDaughters%s",GetName(),str),"",fhNBinsAcceptKink-1,fhBinLimAcceptKink);
679 fhQA[kDcaXY][i] = new TH1F(Form("%s_dcaXY%s",GetName(),str),"",fhNBinsDcaXY-1,fhBinLimDcaXY);
680 fhQA[kDcaZ][i] = new TH1F(Form("%s_dcaZ%s",GetName(),str),"",fhNBinsDcaZ-1,fhBinLimDcaZ);
681 fhQA[kDcaXYnorm][i] = new TH1F(Form("%s_dcaXYnorm%s",GetName(),str),"",fhNBinsDcaXYnorm-1,fhBinLimDcaXYnorm);
682 fhQA[kDcaZnorm][i] = new TH1F(Form("%s_dcaZnorm%s",GetName(),str),"",fhNBinsDcaZnorm-1,fhBinLimDcaZnorm);
683 fhQA[kSigmaDcaXY][i] = new TH1F(Form("%s_sigmaDcaXY%s",GetName(),str),"",fhNBinsSigmaDcaXY-1,fhBinLimSigmaDcaXY);
684 fhQA[kSigmaDcaZ][i] = new TH1F(Form("%s_sigmaDcaZ%s",GetName(),str),"",fhNBinsSigmaDcaZ-1,fhBinLimSigmaDcaZ);
686 fhDcaXYvsDcaZ[i]->SetXTitle("impact par. d_{z}");
687 fhDcaXYvsDcaZ[i]->SetYTitle("impact par. d_{xy}");
688 fhDcaXYvsDcaZnorm[i]->SetXTitle("norm. impact par. d_{z} / #sigma_{z}");
689 fhDcaXYvsDcaZnorm[i]->SetYTitle("norm. impact par. d_{xy} / #sigma_{xy}");
691 fhQA[kCutNSigmaToVertex][i]->SetXTitle("n #sigma to vertex");
692 fhQA[kCutRequireSigmaToVertex][i]->SetXTitle("require #sigma to vertex");
693 fhQA[kCutAcceptKinkDaughters][i]->SetXTitle("accept kink daughters");
694 fhQA[kDcaXY][i]->SetXTitle("impact par. d_{xy}");
695 fhQA[kDcaZ][i]->SetXTitle("impact par. d_{z}");
696 fhQA[kDcaXYnorm][i]->SetXTitle("norm. impact par. d_{xy} / #sigma_{xy}");
697 fhQA[kDcaZnorm][i]->SetXTitle("norm. impact par. d_{z} / #sigma_{z}");
698 fhQA[kSigmaDcaXY][i]->SetXTitle("impact par. resolution #sigma_{xy}");
699 fhQA[kSigmaDcaZ][i]->SetXTitle("impact par. resolution #sigma_{z}");
701 for(Int_t i=0; i<kNHist; i++) fhQA[i][1]->SetLineColor(color);
703 //__________________________________________________________________________________
704 void AliCFTrackIsPrimaryCuts::FillHistograms(TObject* obj, Bool_t f)
707 // fill the QA histograms
712 Bool_t isESDTrack = strcmp(obj->ClassName(),"AliESDtrack") == 0 ? kTRUE : kFALSE ;
713 Bool_t isAODTrack = strcmp(obj->ClassName(),"AliAODTrack") == 0 ? kTRUE : kFALSE ;
715 AliESDtrack * esdTrack = 0x0 ;
716 AliAODTrack * aodTrack = 0x0 ;
717 if (isESDTrack) esdTrack = dynamic_cast<AliESDtrack*>(obj);
718 if (isAODTrack) aodTrack = dynamic_cast<AliAODTrack*>(obj);
720 // f = 0: fill histograms before cuts
721 // f = 1: fill histograms after cuts
723 // get the track to vertex parameter for ESD track
724 Float_t nSigmaToVertex = 0.;
725 if (isESDTrack) nSigmaToVertex = AliESDtrackCuts::GetSigmaToVertex(esdTrack);
730 esdTrack->GetImpactParameters(b,bCov);
731 if (bCov[0]<=0 || bCov[2]<=0) {
732 AliDebug(1, "Estimated b resolution lower or equal zero!");
733 bCov[0]=0; bCov[2]=0;
735 bRes[0] = TMath::Sqrt(bCov[0]);
736 bRes[1] = TMath::Sqrt(bCov[2]);
738 fhQA[kDcaZ][f]->Fill(b[1]);
739 fhQA[kDcaXY][f]->Fill(b[0]);
740 fhDcaXYvsDcaZ[f]->Fill(b[1],b[0]);
741 fhQA[kSigmaDcaXY][f]->Fill(bRes[0]);
742 fhQA[kSigmaDcaZ][f]->Fill(bRes[1]);
744 if (bRes[0]!=0 && bRes[1]!=0) {
745 fhQA[kDcaZnorm][f]->Fill(b[1]/bRes[1]);
746 fhQA[kDcaXYnorm][f]->Fill(b[0]/bRes[0]);
747 fhDcaXYvsDcaZnorm[f]->Fill(b[1]/bRes[1], b[0]/bRes[0]);
750 fhQA[kCutNSigmaToVertex][f]->Fill(nSigmaToVertex);
751 if (nSigmaToVertex<0 && fRequireSigmaToVertex) fhQA[kCutRequireSigmaToVertex][f]->Fill(0.);
752 if (!(nSigmaToVertex<0 && fRequireSigmaToVertex)) fhQA[kCutRequireSigmaToVertex][f]->Fill(1.);
754 if (!fAcceptKinkDaughters && esdTrack->GetKinkIndex(0)>0) fhQA[kCutAcceptKinkDaughters][f]->Fill(0.);
755 if (!(!fAcceptKinkDaughters && esdTrack->GetKinkIndex(0)>0)) fhQA[kCutAcceptKinkDaughters][f]->Fill(0.);
758 // fill cut statistics and cut correlation histograms with information from the bitmap
761 // Get the bitmap of the single cuts
763 SelectionBitMap(obj);
765 // Number of single cuts in this class
766 UInt_t ncuts = fBitmap->GetNbits();
767 for(UInt_t bit=0; bit<ncuts;bit++) {
768 if (!fBitmap->TestBitNumber(bit)) {
769 fhCutStatistics->Fill(bit+1);
770 for (UInt_t bit2=bit; bit2<ncuts;bit2++) {
771 if (!fBitmap->TestBitNumber(bit2))
772 fhCutCorrelation->Fill(bit+1,bit2+1);
777 //__________________________________________________________________________________
778 void AliCFTrackIsPrimaryCuts::SaveHistograms(const Char_t* dir) {
780 // saves the histograms in a directory (dir)
787 gDirectory->mkdir(dir);
790 gDirectory->mkdir("before_cuts");
791 gDirectory->mkdir("after_cuts");
793 fhCutStatistics->Write();
794 fhCutCorrelation->Write();
796 for (Int_t j=0; j<kNStepQA; j++) {
798 gDirectory->cd("before_cuts");
800 gDirectory->cd("after_cuts");
802 fhDcaXYvsDcaZ[j] ->Write();
803 fhDcaXYvsDcaZnorm[j]->Write();
805 for(Int_t i=0; i<kNHist; i++) fhQA[i][j]->Write();
807 gDirectory->cd("../");
809 gDirectory->cd("../");
811 //__________________________________________________________________________________
812 void AliCFTrackIsPrimaryCuts::DrawHistograms()
815 // draws some histograms
820 Float_t right = 0.03;
821 Float_t left = 0.175;
823 Float_t bottom = 0.175;
825 TCanvas* canvas1 = new TCanvas("Track_QA_Primary_1", "Track QA Primary 1", 800, 500);
826 canvas1->Divide(2, 1);
829 fhCutStatistics->SetStats(kFALSE);
830 fhCutStatistics->LabelsOption("v");
831 gPad->SetLeftMargin(left);
832 gPad->SetBottomMargin(0.25);
833 gPad->SetRightMargin(right);
834 gPad->SetTopMargin(0.1);
835 fhCutStatistics->Draw();
838 fhCutCorrelation->SetStats(kFALSE);
839 fhCutCorrelation->LabelsOption("v");
840 gPad->SetLeftMargin(0.30);
841 gPad->SetRightMargin(bottom);
842 gPad->SetTopMargin(0.1);
843 gPad->SetBottomMargin(0.25);
844 fhCutCorrelation->Draw("COLZ");
846 canvas1->SaveAs(Form("%s.eps", canvas1->GetName()));
847 canvas1->SaveAs(Form("%s.ps", canvas1->GetName()));
851 TCanvas* canvas2 = new TCanvas("Track_QA_Primary_2", "Track QA Primary 2", 800, 800);
852 canvas2->Divide(2, 2);
855 gPad->SetRightMargin(right);
856 gPad->SetLeftMargin(left);
857 gPad->SetTopMargin(top);
858 gPad->SetBottomMargin(bottom);
859 fhQA[kDcaXY][0]->SetStats(kFALSE);
860 fhQA[kDcaXY][0]->Draw();
861 fhQA[kDcaXY][1]->Draw("same");
864 gPad->SetRightMargin(right);
865 gPad->SetLeftMargin(left);
866 gPad->SetTopMargin(top);
867 gPad->SetBottomMargin(bottom);
868 fhQA[kDcaZ][0]->SetStats(kFALSE);
869 fhQA[kDcaZ][0]->Draw();
870 fhQA[kDcaZ][1]->Draw("same");
873 // fhDXYvsDZ[0]->SetStats(kFALSE);
875 gPad->SetLeftMargin(bottom);
876 gPad->SetTopMargin(0.1);
877 gPad->SetBottomMargin(bottom);
878 gPad->SetRightMargin(0.2);
879 fhDcaXYvsDcaZ[0]->Draw("COLZ");
882 // fhDXYvsDZ[1]->SetStats(kFALSE);
884 gPad->SetLeftMargin(bottom);
885 gPad->SetTopMargin(0.1);
886 gPad->SetBottomMargin(bottom);
887 gPad->SetRightMargin(0.2);
888 fhDcaXYvsDcaZ[1]->Draw("COLZ");
890 canvas2->SaveAs(Form("%s.eps", canvas2->GetName()));
891 canvas2->SaveAs(Form("%s.ps", canvas2->GetName()));
895 TCanvas* canvas3 = new TCanvas("Track_QA_Primary_3", "Track QA Primary 3", 800, 800);
896 canvas3->Divide(2, 2);
899 gPad->SetRightMargin(right);
900 gPad->SetLeftMargin(left);
901 gPad->SetTopMargin(top);
902 gPad->SetBottomMargin(bottom);
903 fhQA[kDcaXYnorm][0]->SetStats(kFALSE);
904 fhQA[kDcaXYnorm][0]->Draw();
905 fhQA[kDcaXYnorm][1]->Draw("same");
908 gPad->SetRightMargin(right);
909 gPad->SetLeftMargin(left);
910 gPad->SetTopMargin(top);
911 gPad->SetBottomMargin(bottom);
912 fhQA[kDcaZnorm][0]->SetStats(kFALSE);
913 fhQA[kDcaZnorm][0]->Draw();
914 fhQA[kDcaZnorm][1]->Draw("same");
917 // fhDXYvsDZ[0]->SetStats(kFALSE);
919 gPad->SetLeftMargin(bottom);
920 gPad->SetTopMargin(0.1);
921 gPad->SetBottomMargin(bottom);
922 gPad->SetRightMargin(0.2);
923 fhDcaXYvsDcaZnorm[0]->Draw("COLZ");
926 // fhDXYvsDZ[1]->SetStats(kFALSE);
928 gPad->SetLeftMargin(bottom);
929 gPad->SetTopMargin(0.1);
930 gPad->SetBottomMargin(bottom);
931 gPad->SetRightMargin(0.2);
932 fhDcaXYvsDcaZnorm[1]->Draw("COLZ");
934 canvas3->SaveAs(Form("%s.eps", canvas3->GetName()));
935 canvas3->SaveAs(Form("%s.ps", canvas3->GetName()));
939 TCanvas* canvas4 = new TCanvas("Track_QA_Primary_4", "Track QA Primary 4", 1200, 800);
940 canvas4->Divide(3, 2);
943 gPad->SetRightMargin(right);
944 gPad->SetLeftMargin(left);
945 gPad->SetTopMargin(top);
946 gPad->SetBottomMargin(bottom);
947 fhQA[kSigmaDcaXY][0]->SetStats(kFALSE);
948 fhQA[kSigmaDcaXY][0]->Draw();
949 fhQA[kSigmaDcaXY][1]->Draw("same");
952 gPad->SetRightMargin(right);
953 gPad->SetLeftMargin(left);
954 gPad->SetTopMargin(top);
955 gPad->SetBottomMargin(bottom);
956 fhQA[kSigmaDcaZ][0]->SetStats(kFALSE);
957 fhQA[kSigmaDcaZ][0]->Draw();
958 fhQA[kSigmaDcaZ][1]->Draw("same");
961 gPad->SetRightMargin(right);
962 gPad->SetLeftMargin(left);
963 gPad->SetTopMargin(top);
964 gPad->SetBottomMargin(bottom);
965 fhQA[kCutNSigmaToVertex][0]->SetStats(kFALSE);
966 fhQA[kCutNSigmaToVertex][0]->Draw();
967 fhQA[kCutNSigmaToVertex][1]->Draw("same");
970 gPad->SetRightMargin(right);
971 gPad->SetLeftMargin(left);
972 gPad->SetTopMargin(top);
973 gPad->SetBottomMargin(bottom);
974 fhQA[kCutRequireSigmaToVertex][0]->SetStats(kFALSE);
975 fhQA[kCutRequireSigmaToVertex][0]->Draw();
976 fhQA[kCutRequireSigmaToVertex][1]->Draw("same");
979 gPad->SetRightMargin(right);
980 gPad->SetLeftMargin(left);
981 gPad->SetTopMargin(top);
982 gPad->SetBottomMargin(bottom);
983 fhQA[kCutAcceptKinkDaughters][0]->SetStats(kFALSE);
984 fhQA[kCutAcceptKinkDaughters][0]->Draw();
985 fhQA[kCutAcceptKinkDaughters][1]->Draw("same");
987 canvas4->SaveAs(Form("%s.eps", canvas4->GetName()));
988 canvas4->SaveAs(Form("%s.ps", canvas4->GetName()));
990 //__________________________________________________________________________________
991 void AliCFTrackIsPrimaryCuts::AddQAHistograms(TList *qaList) {
993 // saves the histograms in a TList
997 qaList->Add(fhCutStatistics);
998 qaList->Add(fhCutCorrelation);
1000 for (Int_t j=0; j<kNStepQA; j++) {
1001 qaList->Add(fhDcaXYvsDcaZ[j]);
1002 qaList->Add(fhDcaXYvsDcaZnorm[j]);
1003 for(Int_t i=0; i<kNHist; i++)
1004 qaList->Add(fhQA[i][j]);