AliPhysics  8bb951a (8bb951a)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliTrackContainer.cxx
Go to the documentation of this file.
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 #include <TClonesArray.h>
16 
17 #include "AliVEvent.h"
18 #include "AliLog.h"
19 #include "AliVCuts.h"
20 #include "AliESDtrack.h"
21 
22 #include "AliTLorentzVector.h"
25 #include "AliTrackContainer.h"
26 
30 
31 TString AliTrackContainer::fgDefTrackCutsPeriod = "";
32 
38  fTrackFilterType(AliEmcalTrackSelection::kHybridTracks),
39  fListOfCuts(0),
40  fSelectionModeAny(kFALSE),
41  fAODFilterBits(0),
42  fTrackCutsPeriod(),
43  fEmcalTrackSelection(0),
44  fFilteredTracks(0),
45  fTrackTypes(5000)
46 {
47  fBaseClassName = "AliVTrack";
48  SetClassName("AliVTrack");
49  fMassHypothesis = 0.139;
50 }
51 
57 AliTrackContainer::AliTrackContainer(const char *name, const char *period):
59  fTrackFilterType(AliEmcalTrackSelection::kHybridTracks),
60  fListOfCuts(0),
61  fSelectionModeAny(kFALSE),
62  fAODFilterBits(0),
63  fTrackCutsPeriod(period),
64  fEmcalTrackSelection(0),
65  fFilteredTracks(0),
66  fTrackTypes(5000)
67 {
68  fBaseClassName = "AliVTrack";
69  SetClassName("AliVTrack");
70 
72  AliInfo(Form("Default track cuts period is %s", AliTrackContainer::fgDefTrackCutsPeriod.Data()));
74  }
75  fMassHypothesis = 0.139;
76 }
77 
84 void AliTrackContainer::SetArray(AliVEvent *event)
85 {
87 
91  }
92  else {
94 
95  AliInfo("Using custom track cuts");
96 
97  if (fLoadedClass) {
98  if (fLoadedClass->InheritsFrom("AliAODTrack")) {
99  AliInfo(Form("Objects are of type %s: AOD track selection will be done.", fLoadedClass->GetName()));
101  }
102  else if (fLoadedClass->InheritsFrom("AliESDtrack")) {
103  AliInfo(Form("Objects are of type %s: ESD track selection will be done.", fLoadedClass->GetName()));
105  }
106  else {
107  AliWarning(Form("Objects are of type %s: no track filtering will be done!!", fLoadedClass->GetName()));
108  }
109  }
110 
111  if (fEmcalTrackSelection) {
112  if (fSelectionModeAny) {
114  }
115  else {
117  }
118 
120  }
121  }
122  else {
123  if (!fTrackCutsPeriod.IsNull()) {
124  AliInfo(Form("Using track cuts %d for period %s", fTrackFilterType, fTrackCutsPeriod.Data()));
125  }
126  else {
127  AliInfo(Form("Using track cuts %d (no data period was provided!)", fTrackFilterType));
128  }
129 
130  if (fLoadedClass->InheritsFrom("AliAODTrack")) {
131  AliInfo(Form("Objects are of type %s: AOD track selection will be done.", fLoadedClass->GetName()));
133  }
134  else if (fLoadedClass->InheritsFrom("AliESDtrack")) {
135  AliInfo(Form("Objects are of type %s: ESD track selection will be done.", fLoadedClass->GetName()));
137  }
138  else {
139  AliWarning(Form("Objects are of type %s: no track filtering will be done!!", fLoadedClass->GetName()));
140  }
141  }
142  }
143 }
144 
151 {
152  fTrackTypes.Reset(kUndefined);
153  if (fEmcalTrackSelection) {
155 
156  const TClonesArray* trackBitmaps = fEmcalTrackSelection->GetAcceptedTrackBitmaps();
157  TIter nextBitmap(trackBitmaps);
158  TBits* bits = 0;
159  Int_t i = 0;
160  while ((bits = static_cast<TBits*>(nextBitmap()))) {
161  if (i >= fTrackTypes.GetSize()) fTrackTypes.Set((i+1)*2);
162  AliVTrack* vTrack = static_cast<AliVTrack*>(fFilteredTracks->At(i));
163  if (!vTrack) {
164  fTrackTypes[i] = kRejected;
165  }
167  if (bits->FirstSetBit() == 0) {
169  }
170  else if (bits->FirstSetBit() == 1) {
171  if ((vTrack->GetStatus()&AliVTrack::kITSrefit) != 0) {
173  }
174  else {
176  }
177  }
178  }
179  i++;
180  }
181  }
182  else {
184  }
185 }
186 
192 AliVTrack* AliTrackContainer::GetTrack(Int_t i) const
193 {
194  //Get i^th jet in array
195 
196  if (i == -1) i = fCurrentID;
197  if (i < 0 || i >= fFilteredTracks->GetEntriesFast()) return 0;
198  AliVTrack *vp = static_cast<AliVTrack*>(fFilteredTracks->At(i));
199  return vp;
200 }
201 
207 AliVTrack* AliTrackContainer::GetAcceptTrack(Int_t i) const
208 {
209  UInt_t rejectionReason;
210  if (i == -1) i = fCurrentID;
211  if (AcceptTrack(i, rejectionReason)) {
212  return GetTrack(i);
213  }
214  else {
215  AliDebug(2,"Track not accepted.");
216  return 0;
217  }
218 }
219 
226 {
227  const Int_t n = GetNEntries();
228  AliVTrack *p = 0;
229  do {
230  fCurrentID++;
231  if (fCurrentID >= n) break;
233  } while (!p);
234 
235  return p;
236 }
237 
244 {
245  const Int_t n = GetNEntries();
246  AliVTrack *p = 0;
247  do {
248  fCurrentID++;
249  if (fCurrentID >= n) break;
250  p = GetTrack(fCurrentID);
251  } while (!p);
252 
253  return p;
254 }
255 
262 Char_t AliTrackContainer::GetTrackType(const AliVTrack* track) const
263 {
264  Int_t id = fFilteredTracks->IndexOf(track);
265  if (id >= 0) {
266  return fTrackTypes[id];
267  }
268  else {
269  return kUndefined;
270  }
271 }
272 
282 Bool_t AliTrackContainer::GetMomentumFromTrack(TLorentzVector &mom, const AliVTrack* track, Double_t mass) const
283 {
284  if (track) {
285  if (mass < 0) mass = track->M();
286 
287  Bool_t useConstrainedParams = kFALSE;
288  if (fLoadedClass->InheritsFrom("AliESDtrack") && fTrackFilterType == AliEmcalTrackSelection::kHybridTracks) {
289  Char_t trackType = GetTrackType(track);
290  if (trackType == kHybridConstrained || trackType == kHybridConstrainedNoITSrefit) {
291  useConstrainedParams = kTRUE;
292  }
293  }
294 
295  if (useConstrainedParams) {
296  const AliESDtrack *esdtrack = static_cast<const AliESDtrack*>(track);
297  mom.SetPtEtaPhiM(esdtrack->GetConstrainedParam()->Pt(), esdtrack->GetConstrainedParam()->Eta(), esdtrack->GetConstrainedParam()->Phi(), mass);
298  }
299  else {
300  mom.SetPtEtaPhiM(track->Pt(), track->Eta(), track->Phi(), mass);
301  }
302  return kTRUE;
303  }
304  else {
305  mom.SetPtEtaPhiM(0, 0, 0, 0);
306  return kFALSE;
307  }
308 }
309 
317 Bool_t AliTrackContainer::GetMomentumFromTrack(TLorentzVector &mom, const AliVTrack* part) const
318 {
319  return GetMomentumFromTrack(mom,part,fMassHypothesis);
320 }
321 
331 Bool_t AliTrackContainer::GetMomentum(TLorentzVector &mom, Int_t i) const
332 {
333  Double_t mass = fMassHypothesis;
334 
335  if (i == -1) i = fCurrentID;
336  AliVTrack *vp = GetTrack(i);
337  if (vp) {
338  if (mass < 0) mass = vp->M();
339 
340  if (fLoadedClass->InheritsFrom("AliESDtrack") && fTrackFilterType == AliEmcalTrackSelection::kHybridTracks &&
341  (fTrackTypes[i] == kHybridConstrained || fTrackTypes[i] == kHybridConstrainedNoITSrefit)) {
342  AliESDtrack *track = static_cast<AliESDtrack*>(vp);
343  mom.SetPtEtaPhiM(track->GetConstrainedParam()->Pt(), track->GetConstrainedParam()->Eta(), track->GetConstrainedParam()->Phi(), mass);
344  }
345  else {
346  mom.SetPtEtaPhiM(vp->Pt(), vp->Eta(), vp->Phi(), mass);
347  }
348  return kTRUE;
349  }
350  else {
351  mom.SetPtEtaPhiM(0, 0, 0, 0);
352  return kFALSE;
353  }
354 }
355 
365 Bool_t AliTrackContainer::GetNextMomentum(TLorentzVector &mom)
366 {
367  Double_t mass = fMassHypothesis;
368 
369  AliVTrack *vp = GetNextTrack();
370  if (vp) {
371  if (mass < 0) mass = vp->M();
372 
373  if (fLoadedClass->InheritsFrom("AliESDtrack") && fTrackFilterType == AliEmcalTrackSelection::kHybridTracks &&
375  AliESDtrack *track = static_cast<AliESDtrack*>(vp);
376  mom.SetPtEtaPhiM(track->GetConstrainedParam()->Pt(), track->GetConstrainedParam()->Eta(), track->GetConstrainedParam()->Phi(), mass);
377  }
378  else {
379  mom.SetPtEtaPhiM(vp->Pt(), vp->Eta(), vp->Phi(), mass);
380  }
381  return kTRUE;
382  }
383  else {
384  mom.SetPtEtaPhiM(0, 0, 0, 0);
385  return kFALSE;
386  }
387 }
388 
399 Bool_t AliTrackContainer::GetAcceptMomentum(TLorentzVector &mom, Int_t i) const
400 {
401 
402  Double_t mass = fMassHypothesis;
403 
404  if (i == -1) i = fCurrentID;
405  AliVTrack *vp = GetAcceptTrack(i);
406  if (vp) {
407  if (mass < 0) mass = vp->M();
408 
409  if (fLoadedClass->InheritsFrom("AliESDtrack") && fTrackFilterType == AliEmcalTrackSelection::kHybridTracks &&
410  (fTrackTypes[i] == kHybridConstrained || fTrackTypes[i] == kHybridConstrainedNoITSrefit)) {
411  AliESDtrack *track = static_cast<AliESDtrack*>(vp);
412  mom.SetPtEtaPhiM(track->GetConstrainedParam()->Pt(), track->GetConstrainedParam()->Eta(), track->GetConstrainedParam()->Phi(), mass);
413  }
414  else {
415  mom.SetPtEtaPhiM(vp->Pt(), vp->Eta(), vp->Phi(), mass);
416  }
417 
418  return kTRUE;
419  }
420  else {
421  mom.SetPtEtaPhiM(0, 0, 0, 0);
422  return kFALSE;
423  }
424 }
425 
435 Bool_t AliTrackContainer::GetNextAcceptMomentum(TLorentzVector &mom)
436 {
437  Double_t mass = fMassHypothesis;
438 
439  AliVTrack *vp = GetNextAcceptTrack();
440  if (vp) {
441  if (mass < 0) mass = vp->M();
442 
443  if (fLoadedClass->InheritsFrom("AliESDtrack") && fTrackFilterType == AliEmcalTrackSelection::kHybridTracks &&
445  AliESDtrack *track = static_cast<AliESDtrack*>(vp);
446  mom.SetPtEtaPhiM(track->GetConstrainedParam()->Pt(), track->GetConstrainedParam()->Eta(), track->GetConstrainedParam()->Phi(), mass);
447  }
448  else {
449  mom.SetPtEtaPhiM(vp->Pt(), vp->Eta(), vp->Phi(), mass);
450  }
451 
452  return kTRUE;
453  }
454  else {
455  mom.SetPtEtaPhiM(0, 0, 0, 0);
456  return kFALSE;
457  }
458 }
459 
470 Bool_t AliTrackContainer::AcceptTrack(const AliVTrack *vp, UInt_t &rejectionReason) const
471 {
472  Bool_t r = ApplyTrackCuts(vp, rejectionReason);
473  if (!r) return kFALSE;
474 
475  AliTLorentzVector mom;
476  GetMomentumFromTrack(mom, vp);
477 
478  return ApplyKinematicCuts(mom, rejectionReason);
479 }
480 
492 Bool_t AliTrackContainer::AcceptTrack(Int_t i, UInt_t &rejectionReason) const
493 {
494  Bool_t r = ApplyTrackCuts(GetTrack(i), rejectionReason);
495  if (!r) return kFALSE;
496 
497  AliTLorentzVector mom;
498  GetMomentum(mom, i);
499 
500  return ApplyKinematicCuts(mom, rejectionReason);
501 }
502 
512 Bool_t AliTrackContainer::ApplyTrackCuts(const AliVTrack* vp, UInt_t &rejectionReason) const
513 {
514  return ApplyParticleCuts(vp, rejectionReason);
515 }
516 
522 {
523  if (!fListOfCuts) {
524  fListOfCuts = new TObjArray;
525  fListOfCuts->SetOwner(true);
526  }
527  fListOfCuts->Add(cuts);
528 }
529 
535 {
536  if (!fListOfCuts) return 0;
537  return fListOfCuts->GetEntries();
538 }
539 
545 AliVCuts* AliTrackContainer::GetTrackCuts(Int_t icut)
546 {
547  if (!fListOfCuts) return NULL;
548  if (icut < fListOfCuts->GetEntries()) {
549  return static_cast<AliVCuts *>(fListOfCuts->At(icut));
550  }
551  return NULL;
552 }
553 
560 const char* AliTrackContainer::GetTitle() const
561 {
562  static TString trackString;
563 
564  if (GetMinPt() == 0) {
565  trackString = TString::Format("%s_pT0000", GetArrayName().Data());
566  }
567  else if (GetMinPt() < 1.0) {
568  trackString = TString::Format("%s_pT0%3.0f", GetArrayName().Data(), GetMinPt()*1000.0);
569  }
570  else {
571  trackString = TString::Format("%s_pT%4.0f", GetArrayName().Data(), GetMinPt()*1000.0);
572  }
573 
574  return trackString.Data();
575 }
Interface for virtual track selection.
ClassImp(AliAnalysisTaskTriggerRates) AliAnalysisTaskTriggerRates
TObjArray * fFilteredTracks
! tracks filtered using fEmcalTrackSelection
static TString fgDefTrackCutsPeriod
! default period string used to generate track cuts
Double_t mass
Container with name, TClonesArray and cuts for particles.
Track selected under the constrained hybrid track cuts.
virtual Bool_t ApplyParticleCuts(const AliVParticle *vp, UInt_t &rejectionReason) const
UInt_t fAODFilterBits
track filter bits
virtual void SetArray(AliVEvent *event)
const TClonesArray * GetAcceptedTrackBitmaps() const
AliVCuts * GetTrackCuts(Int_t icut)
virtual Bool_t GetNextMomentum(TLorentzVector &mom)
AliEmcalTrackSelection * fEmcalTrackSelection
! track selection object
TClass * fLoadedClass
! Class of the objects contained in the TClonesArray
virtual Bool_t GetMomentum(TLorentzVector &mom, Int_t i) const
ETrackFilterType_t fTrackFilterType
track filter type
TString fTrackCutsPeriod
period string used to generate track cuts
void SetArray(AliVEvent *event)
Container for particles within the EMCAL framework.
Track selected under the constrained hybrid track cuts without ITS refit.
virtual Bool_t GetAcceptMomentum(TLorentzVector &mom, Int_t i) const
virtual Bool_t GetMomentumFromTrack(TLorentzVector &mom, const AliVTrack *track, Double_t mass) const
virtual Bool_t ApplyTrackCuts(const AliVTrack *vp, UInt_t &rejectionReason) const
virtual Bool_t AcceptTrack(const AliVTrack *vp, UInt_t &rejectionReason) const
Double_t GetMinPt() const
Int_t GetNumberOfCutObjects() const
TArrayC fTrackTypes
! track types
TObjArray * fListOfCuts
list of track cut objects
void AddTrackCuts(AliVCuts *cuts)
Implement virtual track selection for AOD analysis.
void SetClassName(const char *clname)
virtual AliVTrack * GetNextTrack()
Bool_t Data(TH1F *h, Double_t *rangefit, Bool_t writefit, Double_t &sgn, Double_t &errsgn, Double_t &bkg, Double_t &errbkg, Double_t &sgnf, Double_t &errsgnf, Double_t &sigmafit, Int_t &status)
Bool_t fSelectionModeAny
accept track if any of the cuts is fulfilled
const TString & GetArrayName() const
Double_t fMassHypothesis
if < 0 it will use a PID mass when available
virtual AliVTrack * GetTrack(Int_t i=-1) const
Track selected under the global hybrid track cuts.
TObjArray * GetAcceptedTracks(const TClonesArray *const tracks)
virtual AliVTrack * GetNextAcceptTrack()
Char_t GetTrackType(const AliVTrack *track) const
void AddTrackCuts(AliVCuts *cuts)
TString fBaseClassName
name of the base class that this container can handle
virtual AliVTrack * GetAcceptTrack(Int_t i=-1) const
virtual Bool_t GetNextAcceptMomentum(TLorentzVector &mom)
const char * GetTitle() const
virtual Bool_t ApplyKinematicCuts(const AliTLorentzVector &mom, UInt_t &rejectionReason) const
TClonesArray * fClArray
! Pointer to array in input event
Int_t GetNEntries() const
Int_t fCurrentID
! current ID for automatic loops
Implementation of virtual track selection for ESDs.
Track status undefined.