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 // By default, the distance to 'vertex calculated from tracks' is used.
35 // Optionally the SPD (tracklet based) or TPC (TPC only tracks based) vertex
37 // Note: the distance to the TPC-vertex is already stored in the ESD,
38 // the distance to the SPD-vertex has to be re-calculated by propagating each
39 // track while executing this cut.
41 // The cut values for these cuts are set with the corresponding set functions.
42 // All cut classes provided by the correction framework are supposed to be
43 // added in the Analysis Framwork's class AliAnalysisFilter and applied by
44 // the filter via a loop.
46 // author: I. Kraus (Ingrid.Kraus@cern.ch)
48 // AliESDtrackCuts writte by Jan Fiete Grosse-Oetringhaus and
49 // AliRsnDaughterCut class written by A. Pulvirenti.
52 #include <TDirectory.h>
56 #include <AliESDtrack.h>
57 #include <AliAODTrack.h>
58 #include <AliESDEvent.h>
59 #include <AliAODEvent.h>
61 #include "AliCFTrackIsPrimaryCuts.h"
63 ClassImp(AliCFTrackIsPrimaryCuts)
65 //__________________________________________________________________________________
66 AliCFTrackIsPrimaryCuts::AliCFTrackIsPrimaryCuts() :
77 fNSigmaToVertexMin(0),
78 fNSigmaToVertexMax(0),
81 fRequireSigmaToVertex(0),
82 fAODType(AliAODTrack::kUndef),
83 fAcceptKinkDaughters(0),
88 fhNBinsRequireSigma(0),
97 fhBinLimRequireSigma(0x0),
98 fhBinLimAcceptKink(0x0),
101 fhBinLimDcaXYnorm(0x0),
102 fhBinLimDcaZnorm(0x0),
103 fhBinLimSigmaDcaXY(0x0),
104 fhBinLimSigmaDcaZ(0x0)
107 // Default constructor
111 //__________________________________________________________________________________
112 AliCFTrackIsPrimaryCuts::AliCFTrackIsPrimaryCuts(Char_t* name, Char_t* title) :
113 AliCFCutBase(name,title),
117 fMinDCAToVertexXY(0),
119 fMaxDCAToVertexXY(0),
123 fNSigmaToVertexMin(0),
124 fNSigmaToVertexMax(0),
127 fRequireSigmaToVertex(0),
128 fAODType(AliAODTrack::kUndef),
129 fAcceptKinkDaughters(0),
134 fhNBinsRequireSigma(0),
135 fhNBinsAcceptKink(0),
140 fhNBinsSigmaDcaXY(0),
143 fhBinLimRequireSigma(0x0),
144 fhBinLimAcceptKink(0x0),
147 fhBinLimDcaXYnorm(0x0),
148 fhBinLimDcaZnorm(0x0),
149 fhBinLimSigmaDcaXY(0x0),
150 fhBinLimSigmaDcaZ(0x0)
157 //__________________________________________________________________________________
158 AliCFTrackIsPrimaryCuts::AliCFTrackIsPrimaryCuts(const AliCFTrackIsPrimaryCuts& c) :
161 fUseSPDvertex(c.fUseSPDvertex),
162 fUseTPCvertex(c.fUseTPCvertex),
163 fMinDCAToVertexXY(c.fMinDCAToVertexXY),
164 fMinDCAToVertexZ(c.fMinDCAToVertexZ),
165 fMaxDCAToVertexXY(c.fMaxDCAToVertexXY),
166 fMaxDCAToVertexZ(c.fMaxDCAToVertexZ),
167 fDCAToVertex2D(c.fDCAToVertex2D),
168 fAbsDCAToVertex(c.fAbsDCAToVertex),
169 fNSigmaToVertexMin(c.fNSigmaToVertexMin),
170 fNSigmaToVertexMax(c.fNSigmaToVertexMax),
171 fSigmaDCAxy(c.fSigmaDCAxy),
172 fSigmaDCAz(c.fSigmaDCAz),
173 fRequireSigmaToVertex(c.fRequireSigmaToVertex),
174 fAODType(c.fAODType),
175 fAcceptKinkDaughters(c.fAcceptKinkDaughters),
176 fhCutStatistics(c.fhCutStatistics),
177 fhCutCorrelation(c.fhCutCorrelation),
179 fhNBinsNSigma(c.fhNBinsNSigma),
180 fhNBinsRequireSigma(c.fhNBinsRequireSigma),
181 fhNBinsAcceptKink(c.fhNBinsAcceptKink),
182 fhNBinsDcaXY(c.fhNBinsDcaXY),
183 fhNBinsDcaZ(c.fhNBinsDcaZ),
184 fhNBinsDcaXYnorm(c.fhNBinsDcaXYnorm),
185 fhNBinsDcaZnorm(c.fhNBinsDcaZnorm),
186 fhNBinsSigmaDcaXY(c.fhNBinsSigmaDcaXY),
187 fhNBinsSigmaDcaZ(c.fhNBinsSigmaDcaZ),
188 fhBinLimNSigma(c.fhBinLimNSigma),
189 fhBinLimRequireSigma(c.fhBinLimRequireSigma),
190 fhBinLimAcceptKink(c.fhBinLimAcceptKink),
191 fhBinLimDcaXY(c.fhBinLimDcaXY),
192 fhBinLimDcaZ(c.fhBinLimDcaZ),
193 fhBinLimDcaXYnorm(c.fhBinLimDcaXYnorm),
194 fhBinLimDcaZnorm(c.fhBinLimDcaZnorm),
195 fhBinLimSigmaDcaXY(c.fhBinLimSigmaDcaXY),
196 fhBinLimSigmaDcaZ(c.fhBinLimSigmaDcaZ)
201 ((AliCFTrackIsPrimaryCuts &) c).Copy(*this);
203 //__________________________________________________________________________________
204 AliCFTrackIsPrimaryCuts& AliCFTrackIsPrimaryCuts::operator=(const AliCFTrackIsPrimaryCuts& c)
207 // Assignment operator
210 AliCFCutBase::operator=(c) ;
212 fUseSPDvertex = c.fUseSPDvertex;
213 fUseTPCvertex = c.fUseTPCvertex;
214 fMinDCAToVertexXY = c.fMinDCAToVertexXY;
215 fMinDCAToVertexZ = c.fMinDCAToVertexZ;
216 fMaxDCAToVertexXY = c.fMaxDCAToVertexXY;
217 fMaxDCAToVertexZ = c.fMaxDCAToVertexZ;
218 fDCAToVertex2D = c.fDCAToVertex2D;
219 fAbsDCAToVertex = c.fAbsDCAToVertex;
220 fNSigmaToVertexMin = c.fNSigmaToVertexMin ;
221 fNSigmaToVertexMax = c.fNSigmaToVertexMax ;
222 fSigmaDCAxy = c.fSigmaDCAxy ;
223 fSigmaDCAz = c.fSigmaDCAz ;
224 fRequireSigmaToVertex = c.fRequireSigmaToVertex ;
225 fAODType = c.fAODType ;
226 fAcceptKinkDaughters = c.fAcceptKinkDaughters ;
227 fhCutStatistics = c.fhCutStatistics ;
228 fhCutCorrelation = c.fhCutCorrelation ;
230 fhNBinsNSigma = c.fhNBinsNSigma;
231 fhNBinsRequireSigma = c.fhNBinsRequireSigma;
232 fhNBinsAcceptKink = c.fhNBinsAcceptKink;
233 fhNBinsDcaXY = c.fhNBinsDcaXY;
234 fhNBinsDcaZ = c.fhNBinsDcaZ;
235 fhNBinsDcaXYnorm = c.fhNBinsDcaXYnorm;
236 fhNBinsDcaZnorm = c.fhNBinsDcaZnorm;
237 fhNBinsSigmaDcaXY = c.fhNBinsSigmaDcaXY;
238 fhNBinsSigmaDcaZ = c.fhNBinsSigmaDcaZ;
239 fhBinLimNSigma = c.fhBinLimNSigma;
240 fhBinLimRequireSigma = c.fhBinLimRequireSigma;
241 fhBinLimAcceptKink = c.fhBinLimAcceptKink;
242 fhBinLimDcaXY = c.fhBinLimDcaXY;
243 fhBinLimDcaZ = c.fhBinLimDcaZ;
244 fhBinLimDcaXYnorm = c.fhBinLimDcaXYnorm;
245 fhBinLimDcaZnorm = c.fhBinLimDcaZnorm;
246 fhBinLimSigmaDcaXY = c.fhBinLimSigmaDcaXY;
247 fhBinLimSigmaDcaZ = c.fhBinLimSigmaDcaZ;
249 for (Int_t j=0; j<6; j++) fDCA[j] = c.fDCA[j];
250 for (Int_t j=0; j<c.kNStepQA; j++){
251 if(c.fhDcaXYvsDcaZ[j]) fhDcaXYvsDcaZ[j] = (TH2F*)c.fhDcaXYvsDcaZ[j]->Clone();
252 for (Int_t i=0; i<c.kNHist; i++){
253 if(c.fhQA[i][j]) fhQA[i][j] = (TH1F*)c.fhQA[i][j]->Clone();
256 ((AliCFTrackIsPrimaryCuts &) c).Copy(*this);
260 //__________________________________________________________________________________
261 AliCFTrackIsPrimaryCuts::~AliCFTrackIsPrimaryCuts()
266 if (fhCutStatistics) delete fhCutStatistics;
267 if (fhCutCorrelation) delete fhCutCorrelation;
269 for (Int_t j=0; j<kNStepQA; j++){
270 if(fhDcaXYvsDcaZ[j]) delete fhDcaXYvsDcaZ[j];
271 for (Int_t i=0; i<kNHist; i++)
272 if(fhQA[i][j]) delete fhQA[i][j];
274 if(fEvt) delete fEvt;
275 if(fBitmap) delete fBitmap;
276 if(fhBinLimNSigma) delete fhBinLimNSigma;
277 if(fhBinLimRequireSigma) delete fhBinLimRequireSigma;
278 if(fhBinLimAcceptKink) delete fhBinLimAcceptKink;
279 if(fhBinLimDcaXY) delete fhBinLimDcaXY;
280 if(fhBinLimDcaZ) delete fhBinLimDcaZ;
281 if(fhBinLimDcaXYnorm) delete fhBinLimDcaXYnorm;
282 if(fhBinLimDcaZnorm) delete fhBinLimDcaZnorm;
283 if(fhBinLimSigmaDcaXY) delete fhBinLimSigmaDcaXY;
284 if(fhBinLimSigmaDcaZ) delete fhBinLimSigmaDcaZ;
286 //__________________________________________________________________________________
287 void AliCFTrackIsPrimaryCuts::Initialise()
290 // sets everything to zero
294 fMinDCAToVertexXY = 0;
295 fMinDCAToVertexZ = 0;
296 fMaxDCAToVertexXY = 0;
297 fMaxDCAToVertexZ = 0;
300 fNSigmaToVertexMin = 0;
301 fNSigmaToVertexMax = 0;
304 fRequireSigmaToVertex = 0;
305 fAcceptKinkDaughters = 0;
306 fAODType = AliAODTrack::kUndef;
308 SetMinDCAToVertexXY();
309 SetMinDCAToVertexZ();
310 SetMaxDCAToVertexXY();
311 SetMaxDCAToVertexZ();
314 SetMinNSigmaToVertex();
315 SetMaxNSigmaToVertex();
318 SetRequireSigmaToVertex();
319 SetAcceptKinkDaughters();
322 for (Int_t j=0; j<6; j++) fDCA[j] = 0.;
323 for (Int_t j=0; j<kNStepQA; j++) {
324 fhDcaXYvsDcaZ[j] = 0x0;
325 for (Int_t i=0; i<kNHist; i++)
329 fhCutCorrelation = 0;
330 fBitmap=new TBits(0);
332 //set default bining for QA histograms
333 SetHistogramBins(kCutNSigmaToVertex,100,0.,10.);
334 SetHistogramBins(kCutRequireSigmaToVertex,5,-0.75,1.75);
335 SetHistogramBins(kCutAcceptKinkDaughters,5,-0.75,1.75);
336 SetHistogramBins(kDcaXY,500,-10.,10.);
337 SetHistogramBins(kDcaZ,500,-10.,10.);
338 SetHistogramBins(kDcaXYnorm,500,-10.,10.);
339 SetHistogramBins(kDcaZnorm,500,-10.,10.);
340 SetHistogramBins(kSigmaDcaXY,500,-0.1,0.9);
341 SetHistogramBins(kSigmaDcaZ,500,-0.1,0.9);
343 //__________________________________________________________________________________
344 void AliCFTrackIsPrimaryCuts::Copy(TObject &c) const
349 AliCFTrackIsPrimaryCuts& target = (AliCFTrackIsPrimaryCuts &) c;
353 if (fhCutStatistics) target.fhCutStatistics = (TH1F*) fhCutStatistics->Clone();
354 if (fhCutCorrelation) target.fhCutCorrelation = (TH2F*) fhCutCorrelation->Clone();
356 for (Int_t j=0; j<6; j++) target.fDCA[j] = fDCA[j];
357 for (Int_t j=0; j<kNStepQA; j++){
358 if(fhDcaXYvsDcaZ[j]) target.fhDcaXYvsDcaZ[j] = (TH2F*)fhDcaXYvsDcaZ[j]->Clone();
359 for (Int_t i=0; i<kNHist; i++)
360 if(fhQA[i][j]) target.fhQA[i][j] = (TH1F*)fhQA[i][j]->Clone();
364 //__________________________________________________________________________________
365 void AliCFTrackIsPrimaryCuts::SetRecEventInfo(const TObject* evt) {
367 // Sets pointer to event information (AliESDEvent or AliAODEvent)
370 AliError("Pointer to AliVEvent !");
373 TString className(evt->ClassName());
374 if (! (className.CompareTo("AliESDEvent")==0 || className.CompareTo("AliAODEvent")==0)) {
375 AliError("argument must point to an AliESDEvent or AliAODEvent !");
378 fEvt = (AliVEvent*) evt;
380 //__________________________________________________________________________________
381 void AliCFTrackIsPrimaryCuts::UseSPDvertex(Bool_t b) {
383 if(fUseTPCvertex && fUseSPDvertex) {
384 fUseSPDvertex = kFALSE;
385 AliError("SPD and TPC vertex chosen. TPC vertex is preferred.");
388 //__________________________________________________________________________________
389 void AliCFTrackIsPrimaryCuts::UseTPCvertex(Bool_t b) {
391 if(fUseTPCvertex) fUseSPDvertex = kFALSE;
393 //__________________________________________________________________________________
394 void AliCFTrackIsPrimaryCuts::GetDCA(AliESDtrack* esdTrack)
396 if (!esdTrack) return;
398 Float_t b[2] = {0.,0.};
399 Float_t bCov[3] = {0.,0.,0.};
400 if(!fUseSPDvertex && !fUseTPCvertex) esdTrack->GetImpactParameters(b,bCov);
401 if( fUseTPCvertex) esdTrack->GetImpactParametersTPC(b,bCov);
405 AliESDEvent * evt = 0x0 ;
406 evt = dynamic_cast<AliESDEvent*>(fEvt);
408 AliError("event not found");
411 const AliESDVertex *vtx = evt->GetVertex();
412 const Double_t Bz = evt->GetMagneticField();
413 AliExternalTrackParam *cParam = 0;
414 Bool_t success = esdTrack->RelateToVertex(vtx, Bz, kVeryBig, cParam);
415 if (success) esdTrack->GetImpactParameters(b,bCov);
418 if (bCov[0]<=0 || bCov[2]<=0) {
419 bCov[0]=0; bCov[2]=0;
421 fDCA[0] = b[0]; // impact parameter xy
422 fDCA[1] = b[1]; // impact parameter z
423 fDCA[2] = TMath::Sqrt(bCov[0]); // resolution xy
424 fDCA[3] = TMath::Sqrt(bCov[2]); // resolution z
426 if (!fAbsDCAToVertex) {
427 if (fDCA[2] > 0) fDCA[0] = fDCA[0]/fDCA[2]; // normalised impact parameter xy
428 if (fDCA[3] > 0) fDCA[1] = fDCA[1]/fDCA[3]; // normalised impact parameter z
432 if(!fUseSPDvertex && !fUseTPCvertex)
433 fDCA[5] = AliESDtrackCuts::GetSigmaToVertex(esdTrack);
437 if (fDCA[2]==0 || fDCA[3]==0)
440 Float_t d = TMath::Sqrt(TMath::Power(b[0]/fDCA[2],2) + TMath::Power(b[1]/fDCA[3],2));
441 if (TMath::Exp(-d * d / 2) < 1e-15)
443 fDCA[5] = TMath::ErfInverse(1 - TMath::Exp(-d * d / 2)) * TMath::Sqrt(2);
447 //__________________________________________________________________________________
448 void AliCFTrackIsPrimaryCuts::GetDCA(AliAODTrack* aodTrack)
450 if (!aodTrack) return;
452 Double_t p[3] = {0.,0.,0.};
453 Double_t cov[21] = {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.};
455 aodTrack->XYZAtDCA(p); // position at DCA
456 aodTrack->GetCovarianceXYZPxPyPz(cov);
459 fDCA[5] = -1; // n_sigma
460 if (p[0]==-999. || p[1]==-999. || p[2]==-999.) {
461 AliError("dca info not available !");
462 fDCA[0] = -999.; // impact parameter xy
463 fDCA[1] = -999.; // impact parameter z
467 AliAODEvent * evt = 0x0;
468 evt = dynamic_cast<AliAODEvent*>(fEvt);
471 // primary vertex is the "best": tracks, SPD or TPC vertex
472 AliAODVertex * primaryVertex = evt->GetVertex(0);
473 // dca = track postion - primary vertex position
474 p[0] = p[0] - primaryVertex->GetX();
475 p[1] = p[1] - primaryVertex->GetY();
476 p[2] = p[2] - primaryVertex->GetZ();
477 // calculate dca in transverse plane
478 Float_t b[2] = {0.,0.};
479 b[0] = TMath::Sqrt(p[0]*p[0]+p[1]*p[1]);
482 Double_t bCov[3] = {0.,0.,0.};
483 // how to calculate the resoultion in the transverse plane ?
484 bCov[0] = 0.; // to do: calculate cov in transverse plane
485 bCov[2] = 0.; // from cov in x and y, need to know correlation
487 if (bCov[0]<=0 || bCov[2]<=0) {
488 bCov[0]=0; bCov[2]=0;
490 fDCA[0] = b[0]; // impact parameter xy
491 fDCA[1] = b[1]; // impact parameter z
492 fDCA[2] = TMath::Sqrt(bCov[0]); // resolution xy
493 fDCA[3] = TMath::Sqrt(bCov[2]); // resolution z
495 if (!fAbsDCAToVertex) {
496 AliError("resolution of impact parameter not available, use absolute dca cut instead !");
497 if (fDCA[2] > 0) fDCA[0] = fDCA[0]/fDCA[2]; // normalised impact parameter xy
498 if (fDCA[3] > 0) fDCA[1] = fDCA[1]/fDCA[3]; // normalised impact parameter z
502 if (fDCA[2]==0 || fDCA[3]==0)
505 Float_t d = TMath::Sqrt(TMath::Power(b[0]/fDCA[2],2) + TMath::Power(b[1]/fDCA[3],2));
506 if (TMath::Exp(-d * d / 2) < 1e-15)
508 fDCA[5] = TMath::ErfInverse(1 - TMath::Exp(-d * d / 2)) * TMath::Sqrt(2);
512 //__________________________________________________________________________________
513 void AliCFTrackIsPrimaryCuts::SelectionBitMap(TObject* obj)
516 // test if the track passes the single cuts
517 // and store the information in a bitmap
520 // bitmap stores the decision of each single cut
521 for(Int_t i=0; i<kNCuts; i++)fBitmap->SetBitNumber(i,kFALSE);
523 // check TObject and cast into ESDtrack
525 if (!obj->InheritsFrom("AliVParticle")) {
526 AliError("object must derived from AliVParticle !");
530 AliESDtrack * esdTrack = dynamic_cast<AliESDtrack*>(obj);
531 AliAODTrack * aodTrack = dynamic_cast<AliAODTrack*>(obj);
533 if (!(esdTrack || aodTrack)) {
534 AliError("object must be an ESDtrack or an AODtrack !");
538 Bool_t isESDTrack = kFALSE;
539 Bool_t isAODTrack = kFALSE;
541 if (esdTrack) isESDTrack = strcmp(obj->ClassName(),"AliESDtrack") == 0 ? kTRUE : kFALSE ;
542 if (aodTrack) isAODTrack = strcmp(obj->ClassName(),"AliAODTrack") == 0 ? kTRUE : kFALSE ;
544 // get the track to vertex parameter for ESD track
545 if (isESDTrack) GetDCA(esdTrack);
546 if (isAODTrack) GetDCA(aodTrack);
548 // check whether dca info is filled
550 if (fDCA[0]>-990. && fDCA[1]>-990.) dcaInfo = 1;
552 Float_t bxy = 0, bz = 0;
553 bxy = TMath::Abs(fDCA[0]);
554 bz = TMath::Abs(fDCA[1]);
556 Float_t b2Dmin = 0, b2Dmax = 0;
557 if (fMinDCAToVertexXY>0 && fMinDCAToVertexZ>0)
558 b2Dmin = fDCA[0]*fDCA[0]/fMinDCAToVertexXY/fMinDCAToVertexXY + fDCA[1]*fDCA[1]/fMinDCAToVertexZ/fMinDCAToVertexZ;
559 if (fMaxDCAToVertexXY>0 && fMaxDCAToVertexZ>0)
560 b2Dmax = fDCA[0]*fDCA[0]/fMaxDCAToVertexXY/fMaxDCAToVertexXY + fDCA[1]*fDCA[1]/fMaxDCAToVertexZ/fMaxDCAToVertexZ;
566 if (!dcaInfo || fDCAToVertex2D || (!fDCAToVertex2D && bxy >= fMinDCAToVertexXY && bxy <= fMaxDCAToVertexXY))
567 fBitmap->SetBitNumber(iCutBit,kTRUE);
570 if (!dcaInfo || fDCAToVertex2D || (!fDCAToVertex2D && bz >= fMinDCAToVertexZ && bz <= fMaxDCAToVertexZ))
571 fBitmap->SetBitNumber(iCutBit,kTRUE);
574 if (!dcaInfo || !fDCAToVertex2D || (fDCAToVertex2D && TMath::Sqrt(b2Dmin) > 1 && TMath::Sqrt(b2Dmax) < 1))
575 fBitmap->SetBitNumber(iCutBit,kTRUE);
578 if (!dcaInfo || (fDCA[5] >= fNSigmaToVertexMin && fDCA[5] <= fNSigmaToVertexMax))
579 fBitmap->SetBitNumber(iCutBit,kTRUE);
582 if (!dcaInfo || fDCA[2] < fSigmaDCAxy)
583 fBitmap->SetBitNumber(iCutBit,kTRUE);
586 if (!dcaInfo || fDCA[3] < fSigmaDCAz)
587 fBitmap->SetBitNumber(iCutBit,kTRUE);
590 if (!dcaInfo || !fRequireSigmaToVertex || (fDCA[5]>=0 && fRequireSigmaToVertex))
591 fBitmap->SetBitNumber(iCutBit,kTRUE);
594 if (!dcaInfo || fAcceptKinkDaughters || (!fAcceptKinkDaughters && esdTrack->GetKinkIndex(0)<=0))
595 fBitmap->SetBitNumber(iCutBit,kTRUE);
599 if (fAODType==AliAODTrack::kUndef || fAODType == aodTrack->GetType()) {
600 fBitmap->SetBitNumber(iCutBit,kTRUE);
603 else fBitmap->SetBitNumber(iCutBit,kTRUE);
607 //__________________________________________________________________________________
608 Bool_t AliCFTrackIsPrimaryCuts::IsSelected(TObject* obj) {
610 // loops over decisions of single cuts and returns if the track is accepted
612 SelectionBitMap(obj);
614 if (fIsQAOn) FillHistograms(obj,0);
615 Bool_t isSelected = kTRUE;
617 for (UInt_t icut=0; icut<fBitmap->GetNbits();icut++) {
618 if(!fBitmap->TestBitNumber(icut)) {
624 if (!isSelected) return kFALSE ;
625 if (fIsQAOn) FillHistograms(obj,1);
628 //__________________________________________________________________________________
629 void AliCFTrackIsPrimaryCuts::SetHistogramBins(Int_t index, Int_t nbins, Double_t *bins)
637 case kCutNSigmaToVertex:
638 fhNBinsNSigma=nbins+1;
639 fhBinLimNSigma=new Double_t[nbins+1];
640 for(Int_t i=0;i<nbins+1;i++)fhBinLimNSigma[i]=bins[i];
643 case kCutRequireSigmaToVertex:
644 fhNBinsRequireSigma=nbins+1;
645 fhBinLimRequireSigma=new Double_t[nbins+1];
646 for(Int_t i=0;i<nbins+1;i++)fhBinLimRequireSigma[i]=bins[i];
649 case kCutAcceptKinkDaughters:
650 fhNBinsAcceptKink=nbins+1;
651 fhBinLimAcceptKink=new Double_t[nbins+1];
652 for(Int_t i=0;i<nbins+1;i++)fhBinLimAcceptKink[i]=bins[i];
656 fhNBinsDcaXY=nbins+1;
657 fhBinLimDcaXY=new Double_t[nbins+1];
658 for(Int_t i=0;i<nbins+1;i++)fhBinLimDcaXY[i]=bins[i];
663 fhBinLimDcaZ=new Double_t[nbins+1];
664 for(Int_t i=0;i<nbins+1;i++)fhBinLimDcaZ[i]=bins[i];
668 fhNBinsDcaXYnorm=nbins+1;
669 fhBinLimDcaXYnorm=new Double_t[nbins+1];
670 for(Int_t i=0;i<nbins+1;i++)fhBinLimDcaXYnorm[i]=bins[i];
674 fhNBinsDcaZnorm=nbins+1;
675 fhBinLimDcaZnorm=new Double_t[nbins+1];
676 for(Int_t i=0;i<nbins+1;i++)fhBinLimDcaZnorm[i]=bins[i];
680 fhNBinsSigmaDcaXY=nbins+1;
681 fhBinLimSigmaDcaXY=new Double_t[nbins+1];
682 for(Int_t i=0;i<nbins+1;i++)fhBinLimSigmaDcaXY[i]=bins[i];
686 fhNBinsSigmaDcaZ=nbins+1;
687 fhBinLimSigmaDcaZ=new Double_t[nbins+1];
688 for(Int_t i=0;i<nbins+1;i++)fhBinLimSigmaDcaZ[i]=bins[i];
692 //__________________________________________________________________________________
693 void AliCFTrackIsPrimaryCuts::SetHistogramBins(Int_t index, Int_t nbins, Double_t xmin, Double_t xmax)
699 case kCutNSigmaToVertex:
700 fhNBinsNSigma=nbins+1;
701 fhBinLimNSigma=new Double_t[nbins+1];
702 for(Int_t i=0;i<nbins+1;i++)fhBinLimNSigma[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
705 case kCutRequireSigmaToVertex:
706 fhNBinsRequireSigma=nbins+1;
707 fhBinLimRequireSigma=new Double_t[nbins+1];
708 for(Int_t i=0;i<nbins+1;i++)fhBinLimRequireSigma[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
711 case kCutAcceptKinkDaughters:
712 fhNBinsAcceptKink=nbins+1;
713 fhBinLimAcceptKink=new Double_t[nbins+1];
714 for(Int_t i=0;i<nbins+1;i++)fhBinLimAcceptKink[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
718 fhNBinsDcaXY=nbins+1;
719 fhBinLimDcaXY=new Double_t[nbins+1];
720 for(Int_t i=0;i<nbins+1;i++)fhBinLimDcaXY[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
725 fhBinLimDcaZ=new Double_t[nbins+1];
726 for(Int_t i=0;i<nbins+1;i++)fhBinLimDcaZ[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
730 fhNBinsDcaXYnorm=nbins+1;
731 fhBinLimDcaXYnorm=new Double_t[nbins+1];
732 for(Int_t i=0;i<nbins+1;i++)fhBinLimDcaXYnorm[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
736 fhNBinsDcaZnorm=nbins+1;
737 fhBinLimDcaZnorm=new Double_t[nbins+1];
738 for(Int_t i=0;i<nbins+1;i++)fhBinLimDcaZnorm[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
742 fhNBinsSigmaDcaXY=nbins+1;
743 fhBinLimSigmaDcaXY=new Double_t[nbins+1];
744 for(Int_t i=0;i<nbins+1;i++)fhBinLimSigmaDcaXY[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
748 fhNBinsSigmaDcaZ=nbins+1;
749 fhBinLimSigmaDcaZ=new Double_t[nbins+1];
750 for(Int_t i=0;i<nbins+1;i++)fhBinLimSigmaDcaZ[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
754 //__________________________________________________________________________________
755 void AliCFTrackIsPrimaryCuts::DefineHistograms() {
757 // histograms for cut variables, cut statistics and cut correlations
761 // book cut statistics and cut correlation histograms
762 fhCutStatistics = new TH1F(Form("%s_cut_statistics",GetName()), Form("%s cut statistics",GetName()), kNCuts,0.5,kNCuts+0.5);
763 fhCutStatistics->SetLineWidth(2);
764 fhCutStatistics->GetXaxis()->SetBinLabel(1,"dca xy");
765 fhCutStatistics->GetXaxis()->SetBinLabel(2,"dca z");
766 fhCutStatistics->GetXaxis()->SetBinLabel(3,"dca ellipse");
767 fhCutStatistics->GetXaxis()->SetBinLabel(4,"n dca");
768 fhCutStatistics->GetXaxis()->SetBinLabel(5,"sigma dca xy");
769 fhCutStatistics->GetXaxis()->SetBinLabel(6,"sigma dca z");
770 fhCutStatistics->GetXaxis()->SetBinLabel(7,"require dca");
771 fhCutStatistics->GetXaxis()->SetBinLabel(8,"kink daughter");
772 fhCutStatistics->GetXaxis()->SetBinLabel(9,"AOD type");
774 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);
775 fhCutCorrelation->SetLineWidth(2);
776 fhCutCorrelation->GetXaxis()->SetBinLabel(1,"dca xy");
777 fhCutCorrelation->GetXaxis()->SetBinLabel(2,"dca z");
778 fhCutCorrelation->GetXaxis()->SetBinLabel(3,"dca ellipse");
779 fhCutCorrelation->GetXaxis()->SetBinLabel(4,"n dca");
780 fhCutCorrelation->GetXaxis()->SetBinLabel(5,"sigma dca xy");
781 fhCutCorrelation->GetXaxis()->SetBinLabel(6,"sigma dca z");
782 fhCutCorrelation->GetXaxis()->SetBinLabel(7,"require dca");
783 fhCutCorrelation->GetXaxis()->SetBinLabel(8,"kink daughter");
784 fhCutCorrelation->GetXaxis()->SetBinLabel(9,"AOD type");
786 fhCutCorrelation->GetYaxis()->SetBinLabel(1,"dca xy");
787 fhCutCorrelation->GetYaxis()->SetBinLabel(2,"dca z");
788 fhCutCorrelation->GetYaxis()->SetBinLabel(3,"dca ellipse");
789 fhCutCorrelation->GetYaxis()->SetBinLabel(4,"n dca");
790 fhCutCorrelation->GetYaxis()->SetBinLabel(5,"sigma dca xy");
791 fhCutCorrelation->GetYaxis()->SetBinLabel(6,"sigma dca z");
792 fhCutCorrelation->GetYaxis()->SetBinLabel(7,"require dca");
793 fhCutCorrelation->GetYaxis()->SetBinLabel(8,"kink daughter");
794 fhCutCorrelation->GetYaxis()->SetBinLabel(9,"AOD type");
796 // book QA histograms
798 for (Int_t i=0; i<kNStepQA; i++) {
799 if (i==0) snprintf(str,5," ");
800 else snprintf(str,5,"_cut");
802 fhDcaXYvsDcaZ[i] = new TH2F(Form("%s_dcaXYvsDcaZ%s",GetName(),str),"",200,-10,10,200,-10,10);
803 fhQA[kCutNSigmaToVertex][i] = new TH1F(Form("%s_nSigmaToVertex%s",GetName(),str),"",fhNBinsNSigma-1,fhBinLimNSigma);
804 fhQA[kCutRequireSigmaToVertex][i] = new TH1F(Form("%s_requireSigmaToVertex%s",GetName(),str),"",fhNBinsRequireSigma-1,fhBinLimRequireSigma);
805 fhQA[kCutAcceptKinkDaughters][i] = new TH1F(Form("%s_acceptKinkDaughters%s",GetName(),str),"",fhNBinsAcceptKink-1,fhBinLimAcceptKink);
806 fhQA[kDcaXY][i] = new TH1F(Form("%s_dcaXY%s",GetName(),str),"",fhNBinsDcaXY-1,fhBinLimDcaXY);
807 fhQA[kDcaZ][i] = new TH1F(Form("%s_dcaZ%s",GetName(),str),"",fhNBinsDcaZ-1,fhBinLimDcaZ);
808 fhQA[kDcaXYnorm][i] = new TH1F(Form("%s_dcaXYnorm%s",GetName(),str),"",fhNBinsDcaXYnorm-1,fhBinLimDcaXYnorm);
809 fhQA[kDcaZnorm][i] = new TH1F(Form("%s_dcaZnorm%s",GetName(),str),"",fhNBinsDcaZnorm-1,fhBinLimDcaZnorm);
810 fhQA[kSigmaDcaXY][i] = new TH1F(Form("%s_sigmaDcaXY%s",GetName(),str),"",fhNBinsSigmaDcaXY-1,fhBinLimSigmaDcaXY);
811 fhQA[kSigmaDcaZ][i] = new TH1F(Form("%s_sigmaDcaZ%s",GetName(),str),"",fhNBinsSigmaDcaZ-1,fhBinLimSigmaDcaZ);
813 fhDcaXYvsDcaZ[i]->SetXTitle("impact par. d_{z}");
814 fhDcaXYvsDcaZ[i]->SetYTitle("impact par. d_{xy}");
816 fhQA[kCutNSigmaToVertex][i]->SetXTitle("n #sigma to vertex");
817 fhQA[kCutRequireSigmaToVertex][i]->SetXTitle("require #sigma to vertex");
818 fhQA[kCutAcceptKinkDaughters][i]->SetXTitle("accept kink daughters");
819 fhQA[kDcaXY][i]->SetXTitle("impact par. d_{xy}");
820 fhQA[kDcaZ][i]->SetXTitle("impact par. d_{z}");
821 fhQA[kDcaXYnorm][i]->SetXTitle("norm. impact par. d_{xy} / #sigma_{xy}");
822 fhQA[kDcaZnorm][i]->SetXTitle("norm. impact par. d_{z} / #sigma_{z}");
823 fhQA[kSigmaDcaXY][i]->SetXTitle("impact par. resolution #sigma_{xy}");
824 fhQA[kSigmaDcaZ][i]->SetXTitle("impact par. resolution #sigma_{z}");
826 for(Int_t i=0; i<kNHist; i++) fhQA[i][1]->SetLineColor(color);
828 //__________________________________________________________________________________
829 void AliCFTrackIsPrimaryCuts::FillHistograms(TObject* obj, Bool_t f)
832 // fill the QA histograms
836 Bool_t isESDTrack = strcmp(obj->ClassName(),"AliESDtrack") == 0 ? kTRUE : kFALSE ;
837 Bool_t isAODTrack = strcmp(obj->ClassName(),"AliAODTrack") == 0 ? kTRUE : kFALSE ;
839 AliESDtrack * esdTrack = 0x0 ;
840 AliAODTrack * aodTrack = 0x0 ;
841 if (isESDTrack) esdTrack = dynamic_cast<AliESDtrack*>(obj);
842 if (isAODTrack) aodTrack = dynamic_cast<AliAODTrack*>(obj);
844 // f = 0: fill histograms before cuts
845 // f = 1: fill histograms after cuts
847 // get the track to vertex parameter for ESD track
850 fhQA[kDcaZ][f]->Fill(fDCA[1]);
851 fhQA[kDcaXY][f]->Fill(fDCA[0]);
852 fhDcaXYvsDcaZ[f]->Fill(fDCA[1],fDCA[0]);
853 fhQA[kSigmaDcaXY][f]->Fill(fDCA[2]);
854 fhQA[kSigmaDcaZ][f]->Fill(fDCA[3]);
855 // // // // // // // delete histograms
856 fhQA[kDcaZnorm][f]->Fill(fDCA[1]);
857 fhQA[kDcaXYnorm][f]->Fill(fDCA[0]);
859 fhQA[kCutNSigmaToVertex][f]->Fill(fDCA[5]);
860 if (fDCA[5]<0 && fRequireSigmaToVertex)
861 fhQA[kCutRequireSigmaToVertex][f]->Fill(0.);
863 fhQA[kCutRequireSigmaToVertex][f]->Fill(1.);
865 if (!fAcceptKinkDaughters && esdTrack->GetKinkIndex(0)>0)
866 fhQA[kCutAcceptKinkDaughters][f]->Fill(0.);
868 fhQA[kCutAcceptKinkDaughters][f]->Fill(1.);
871 // fill cut statistics and cut correlation histograms with information from the bitmap
874 // Get the bitmap of the single cuts
875 SelectionBitMap(obj);
877 // Number of single cuts in this class
878 UInt_t ncuts = fBitmap->GetNbits();
879 for(UInt_t bit=0; bit<ncuts;bit++) {
880 if (!fBitmap->TestBitNumber(bit)) {
881 fhCutStatistics->Fill(bit+1);
882 for (UInt_t bit2=bit; bit2<ncuts;bit2++) {
883 if (!fBitmap->TestBitNumber(bit2))
884 fhCutCorrelation->Fill(bit+1,bit2+1);
889 //__________________________________________________________________________________
890 void AliCFTrackIsPrimaryCuts::SaveHistograms(const Char_t* dir) {
892 // saves the histograms in a directory (dir)
899 gDirectory->mkdir(dir);
902 gDirectory->mkdir("before_cuts");
903 gDirectory->mkdir("after_cuts");
905 fhCutStatistics->Write();
906 fhCutCorrelation->Write();
908 for (Int_t j=0; j<kNStepQA; j++) {
910 gDirectory->cd("before_cuts");
912 gDirectory->cd("after_cuts");
914 fhDcaXYvsDcaZ[j] ->Write();
916 for(Int_t i=0; i<kNHist; i++) fhQA[i][j]->Write();
918 gDirectory->cd("../");
920 gDirectory->cd("../");
922 //__________________________________________________________________________________
923 void AliCFTrackIsPrimaryCuts::DrawHistograms()
926 // draws some histograms
931 Float_t right = 0.03;
932 Float_t left = 0.175;
934 Float_t bottom = 0.175;
936 TCanvas* canvas1 = new TCanvas("Track_QA_Primary_1", "Track QA Primary 1", 800, 500);
937 canvas1->Divide(2, 1);
940 fhCutStatistics->SetStats(kFALSE);
941 fhCutStatistics->LabelsOption("v");
942 gPad->SetLeftMargin(left);
943 gPad->SetBottomMargin(0.25);
944 gPad->SetRightMargin(right);
945 gPad->SetTopMargin(0.1);
946 fhCutStatistics->Draw();
949 fhCutCorrelation->SetStats(kFALSE);
950 fhCutCorrelation->LabelsOption("v");
951 gPad->SetLeftMargin(0.30);
952 gPad->SetRightMargin(bottom);
953 gPad->SetTopMargin(0.1);
954 gPad->SetBottomMargin(0.25);
955 fhCutCorrelation->Draw("COLZ");
957 canvas1->SaveAs(Form("%s.eps", canvas1->GetName()));
958 canvas1->SaveAs(Form("%s.ps", canvas1->GetName()));
962 TCanvas* canvas2 = new TCanvas("Track_QA_Primary_2", "Track QA Primary 2", 800, 800);
963 canvas2->Divide(2, 2);
966 gPad->SetRightMargin(right);
967 gPad->SetLeftMargin(left);
968 gPad->SetTopMargin(top);
969 gPad->SetBottomMargin(bottom);
970 fhQA[kDcaXY][0]->SetStats(kFALSE);
971 fhQA[kDcaXY][0]->Draw();
972 fhQA[kDcaXY][1]->Draw("same");
975 gPad->SetRightMargin(right);
976 gPad->SetLeftMargin(left);
977 gPad->SetTopMargin(top);
978 gPad->SetBottomMargin(bottom);
979 fhQA[kDcaZ][0]->SetStats(kFALSE);
980 fhQA[kDcaZ][0]->Draw();
981 fhQA[kDcaZ][1]->Draw("same");
984 // fhDXYvsDZ[0]->SetStats(kFALSE);
986 gPad->SetLeftMargin(bottom);
987 gPad->SetTopMargin(0.1);
988 gPad->SetBottomMargin(bottom);
989 gPad->SetRightMargin(0.2);
990 fhDcaXYvsDcaZ[0]->Draw("COLZ");
993 // fhDXYvsDZ[1]->SetStats(kFALSE);
995 gPad->SetLeftMargin(bottom);
996 gPad->SetTopMargin(0.1);
997 gPad->SetBottomMargin(bottom);
998 gPad->SetRightMargin(0.2);
999 fhDcaXYvsDcaZ[1]->Draw("COLZ");
1001 canvas2->SaveAs(Form("%s.eps", canvas2->GetName()));
1002 canvas2->SaveAs(Form("%s.ps", canvas2->GetName()));
1006 TCanvas* canvas3 = new TCanvas("Track_QA_Primary_3", "Track QA Primary 3", 800, 400);
1007 canvas3->Divide(2, 1);
1010 gPad->SetRightMargin(right);
1011 gPad->SetLeftMargin(left);
1012 gPad->SetTopMargin(top);
1013 gPad->SetBottomMargin(bottom);
1014 fhQA[kDcaXYnorm][0]->SetStats(kFALSE);
1015 fhQA[kDcaXYnorm][0]->Draw();
1016 fhQA[kDcaXYnorm][1]->Draw("same");
1019 gPad->SetRightMargin(right);
1020 gPad->SetLeftMargin(left);
1021 gPad->SetTopMargin(top);
1022 gPad->SetBottomMargin(bottom);
1023 fhQA[kDcaZnorm][0]->SetStats(kFALSE);
1024 fhQA[kDcaZnorm][0]->Draw();
1025 fhQA[kDcaZnorm][1]->Draw("same");
1027 canvas3->SaveAs(Form("%s.eps", canvas3->GetName()));
1028 canvas3->SaveAs(Form("%s.ps", canvas3->GetName()));
1032 TCanvas* canvas4 = new TCanvas("Track_QA_Primary_4", "Track QA Primary 4", 1200, 800);
1033 canvas4->Divide(3, 2);
1036 gPad->SetRightMargin(right);
1037 gPad->SetLeftMargin(left);
1038 gPad->SetTopMargin(top);
1039 gPad->SetBottomMargin(bottom);
1040 fhQA[kSigmaDcaXY][0]->SetStats(kFALSE);
1041 fhQA[kSigmaDcaXY][0]->Draw();
1042 fhQA[kSigmaDcaXY][1]->Draw("same");
1045 gPad->SetRightMargin(right);
1046 gPad->SetLeftMargin(left);
1047 gPad->SetTopMargin(top);
1048 gPad->SetBottomMargin(bottom);
1049 fhQA[kSigmaDcaZ][0]->SetStats(kFALSE);
1050 fhQA[kSigmaDcaZ][0]->Draw();
1051 fhQA[kSigmaDcaZ][1]->Draw("same");
1054 gPad->SetRightMargin(right);
1055 gPad->SetLeftMargin(left);
1056 gPad->SetTopMargin(top);
1057 gPad->SetBottomMargin(bottom);
1058 fhQA[kCutNSigmaToVertex][0]->SetStats(kFALSE);
1059 fhQA[kCutNSigmaToVertex][0]->Draw();
1060 fhQA[kCutNSigmaToVertex][1]->Draw("same");
1063 gPad->SetRightMargin(right);
1064 gPad->SetLeftMargin(left);
1065 gPad->SetTopMargin(top);
1066 gPad->SetBottomMargin(bottom);
1067 fhQA[kCutRequireSigmaToVertex][0]->SetStats(kFALSE);
1068 fhQA[kCutRequireSigmaToVertex][0]->Draw();
1069 fhQA[kCutRequireSigmaToVertex][1]->Draw("same");
1072 gPad->SetRightMargin(right);
1073 gPad->SetLeftMargin(left);
1074 gPad->SetTopMargin(top);
1075 gPad->SetBottomMargin(bottom);
1076 fhQA[kCutAcceptKinkDaughters][0]->SetStats(kFALSE);
1077 fhQA[kCutAcceptKinkDaughters][0]->Draw();
1078 fhQA[kCutAcceptKinkDaughters][1]->Draw("same");
1080 canvas4->SaveAs(Form("%s.eps", canvas4->GetName()));
1081 canvas4->SaveAs(Form("%s.ps", canvas4->GetName()));
1083 //__________________________________________________________________________________
1084 void AliCFTrackIsPrimaryCuts::AddQAHistograms(TList *qaList) {
1086 // saves the histograms in a TList
1090 qaList->Add(fhCutStatistics);
1091 qaList->Add(fhCutCorrelation);
1093 for (Int_t j=0; j<kNStepQA; j++) {
1094 qaList->Add(fhDcaXYvsDcaZ[j]);
1095 for(Int_t i=0; i<kNHist; i++)
1096 qaList->Add(fhQA[i][j]);