AliPhysics  75b74d3 (75b74d3)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
AliTrackContainer.cxx
Go to the documentation of this file.
1 //
2 // Container with name, TClonesArray and cuts for particles
3 //
4 // Author: M. Verweij, S. Aiola
5 
6 #include <TClonesArray.h>
7 
8 #include "AliVEvent.h"
9 #include "AliLog.h"
10 #include "AliVCuts.h"
11 #include "AliESDtrack.h"
12 
13 #include "AliTLorentzVector.h"
16 #include "AliTrackContainer.h"
17 
19 
20 TString AliTrackContainer::fgDefTrackCutsPeriod = "";
21 
22 //________________________________________________________________________
25  fTrackFilterType(AliEmcalTrackSelection::kHybridTracks),
26  fListOfCuts(0),
27  fSelectionModeAny(kFALSE),
28  fAODFilterBits(0),
29  fTrackCutsPeriod(),
30  fEmcalTrackSelection(0),
31  fFilteredTracks(0),
32  fTrackTypes(5000)
33 {
34  // Default constructor.
35 
36  fClassName = "AliVTrack";
37  fMassHypothesis = 0.139;
38 }
39 
40 //________________________________________________________________________
41 AliTrackContainer::AliTrackContainer(const char *name, const char *period):
43  fTrackFilterType(AliEmcalTrackSelection::kHybridTracks),
44  fListOfCuts(0),
45  fSelectionModeAny(kFALSE),
46  fAODFilterBits(0),
47  fTrackCutsPeriod(period),
48  fEmcalTrackSelection(0),
49  fFilteredTracks(0),
50  fTrackTypes(5000)
51 {
52  // Standard constructor.
53 
54  fClassName = "AliVTrack";
55 
57  AliInfo(Form("Default track cuts period is %s", AliTrackContainer::fgDefTrackCutsPeriod.Data()));
59  }
60  fMassHypothesis = 0.139;
61 }
62 
63 //________________________________________________________________________
64 void AliTrackContainer::SetArray(AliVEvent *event)
65 {
66  // Get array from event.
67 
69 
73  }
74  else {
76 
77  AliInfo("Using custom track cuts");
78 
79  if (fLoadedClass) {
80  if (fLoadedClass->InheritsFrom("AliAODTrack")) {
81  AliInfo(Form("Objects are of type %s: AOD track selection will be done.", fLoadedClass->GetName()));
83  }
84  else if (fLoadedClass->InheritsFrom("AliESDtrack")) {
85  AliInfo(Form("Objects are of type %s: ESD track selection will be done.", fLoadedClass->GetName()));
87  }
88  else {
89  AliWarning(Form("Objects are of type %s: no track filtering will be done!!", fLoadedClass->GetName()));
90  }
91  }
92 
94  if (fSelectionModeAny) {
96  }
97  else {
99  }
100 
102  }
103  }
104  else {
105  if (!fTrackCutsPeriod.IsNull()) {
106  AliInfo(Form("Using track cuts %d for period %s", fTrackFilterType, fTrackCutsPeriod.Data()));
107  }
108  else {
109  AliInfo(Form("Using track cuts %d (no data period was provided!)", fTrackFilterType));
110  }
111 
112  if (fLoadedClass->InheritsFrom("AliAODTrack")) {
113  AliInfo(Form("Objects are of type %s: AOD track selection will be done.", fLoadedClass->GetName()));
115  }
116  else if (fLoadedClass->InheritsFrom("AliESDtrack")) {
117  AliInfo(Form("Objects are of type %s: ESD track selection will be done.", fLoadedClass->GetName()));
119  }
120  else {
121  AliWarning(Form("Objects are of type %s: no track filtering will be done!!", fLoadedClass->GetName()));
122  }
123  }
124  }
125 }
126 
127 //________________________________________________________________________
129 {
130  fTrackTypes.Reset(kUndefined);
131  if (fEmcalTrackSelection) {
133 
134  const TClonesArray* trackBitmaps = fEmcalTrackSelection->GetAcceptedTrackBitmaps();
135  TIter nextBitmap(trackBitmaps);
136  TBits* bits = 0;
137  Int_t i = 0;
138  while ((bits = static_cast<TBits*>(nextBitmap()))) {
139  if (i >= fTrackTypes.GetSize()) fTrackTypes.Set((i+1)*2);
140  AliVTrack* vTrack = static_cast<AliVTrack*>(fFilteredTracks->At(i));
141  if (!vTrack) {
142  fTrackTypes[i] = kRejected;
143  }
145  if (bits->FirstSetBit() == 0) {
147  }
148  else if (bits->FirstSetBit() == 1) {
149  if ((vTrack->GetStatus()&AliVTrack::kITSrefit) != 0) {
151  }
152  else {
154  }
155  }
156  }
157  i++;
158  }
159  }
160  else {
162  }
163 }
164 
165 //________________________________________________________________________
166 AliVTrack* AliTrackContainer::GetTrack(Int_t i) const
167 {
168  //Get i^th jet in array
169 
170  if (i == -1) i = fCurrentID;
171  if (i < 0 || i >= fFilteredTracks->GetEntriesFast()) return 0;
172  AliVTrack *vp = static_cast<AliVTrack*>(fFilteredTracks->At(i));
173  return vp;
174 }
175 
176 //________________________________________________________________________
178 {
179  //return pointer to particle if particle is accepted
180 
181  if (i == -1) i = fCurrentID;
182  if (AcceptTrack(i)) {
183  return GetTrack(i);
184  }
185  else {
186  AliDebug(2,"Track not accepted.");
187  return 0;
188  }
189 }
190 
191 //________________________________________________________________________
193 {
194  //Get next accepted particle
195 
196  const Int_t n = GetNEntries();
197  AliVTrack *p = 0;
198  do {
199  fCurrentID++;
200  if (fCurrentID >= n) break;
202  } while (!p);
203 
204  return p;
205 }
206 
207 //________________________________________________________________________
209 {
210  //Get next particle
211 
212  const Int_t n = GetNEntries();
213  AliVTrack *p = 0;
214  do {
215  fCurrentID++;
216  if (fCurrentID >= n) break;
217  p = GetTrack(fCurrentID);
218  } while (!p);
219 
220  return p;
221 }
222 
223 //________________________________________________________________________
224 Bool_t AliTrackContainer::GetMomentum(TLorentzVector &mom, const AliVTrack* part, Double_t mass)
225 {
226  if (part) {
227  if (mass < 0) mass = part->M();
228  mom.SetPtEtaPhiM(part->Pt(), part->Eta(), part->Phi(), mass);
229  return kTRUE;
230  }
231  else {
232  mom.SetPtEtaPhiM(0, 0, 0, 0);
233  return kFALSE;
234  }
235 }
236 
237 //________________________________________________________________________
238 Bool_t AliTrackContainer::GetMomentum(TLorentzVector &mom, const AliVTrack* part)
239 {
240  return GetMomentum(mom,part,fMassHypothesis);
241 }
242 
243 //________________________________________________________________________
244 Bool_t AliTrackContainer::GetMomentum(TLorentzVector &mom, Int_t i)
245 {
246  //Get momentum of the i^th particle in array
247 
248  Double_t mass = fMassHypothesis;
249 
250  if (i == -1) i = fCurrentID;
251  AliVTrack *vp = GetTrack(i);
252  if (vp) {
253  if (mass < 0) mass = vp->M();
254 
255  if (fLoadedClass->InheritsFrom("AliESDtrack") && fTrackFilterType == AliEmcalTrackSelection::kHybridTracks &&
256  (fTrackTypes[i] == kHybridConstrained || fTrackTypes[i] == kHybridConstrainedNoITSrefit)) {
257  AliESDtrack *track = static_cast<AliESDtrack*>(vp);
258  mom.SetPtEtaPhiM(track->GetConstrainedParam()->Pt(), track->GetConstrainedParam()->Eta(), track->GetConstrainedParam()->Phi(), mass);
259  }
260  else {
261  mom.SetPtEtaPhiM(vp->Pt(), vp->Eta(), vp->Phi(), mass);
262  }
263  return kTRUE;
264  }
265  else {
266  mom.SetPtEtaPhiM(0, 0, 0, 0);
267  return kFALSE;
268  }
269 }
270 
271 //________________________________________________________________________
272 Bool_t AliTrackContainer::GetNextMomentum(TLorentzVector &mom)
273 {
274  //Get momentum of the next particle in array
275 
276  Double_t mass = fMassHypothesis;
277 
278  AliVTrack *vp = GetNextTrack();
279  if (vp) {
280  if (mass < 0) mass = vp->M();
281 
282  if (fLoadedClass->InheritsFrom("AliESDtrack") && fTrackFilterType == AliEmcalTrackSelection::kHybridTracks &&
284  AliESDtrack *track = static_cast<AliESDtrack*>(vp);
285  mom.SetPtEtaPhiM(track->GetConstrainedParam()->Pt(), track->GetConstrainedParam()->Eta(), track->GetConstrainedParam()->Phi(), mass);
286  }
287  else {
288  mom.SetPtEtaPhiM(vp->Pt(), vp->Eta(), vp->Phi(), mass);
289  }
290  return kTRUE;
291  }
292  else {
293  mom.SetPtEtaPhiM(0, 0, 0, 0);
294  return kFALSE;
295  }
296 }
297 
298 //________________________________________________________________________
299 Bool_t AliTrackContainer::GetAcceptMomentum(TLorentzVector &mom, Int_t i)
300 {
301  //Get momentum of the i^th particle in array
302 
303  Double_t mass = fMassHypothesis;
304 
305  if (i == -1) i = fCurrentID;
306  AliVTrack *vp = GetAcceptTrack(i);
307  if (vp) {
308  if (mass < 0) mass = vp->M();
309 
310  if (fLoadedClass->InheritsFrom("AliESDtrack") && fTrackFilterType == AliEmcalTrackSelection::kHybridTracks &&
311  (fTrackTypes[i] == kHybridConstrained || fTrackTypes[i] == kHybridConstrainedNoITSrefit)) {
312  AliESDtrack *track = static_cast<AliESDtrack*>(vp);
313  mom.SetPtEtaPhiM(track->GetConstrainedParam()->Pt(), track->GetConstrainedParam()->Eta(), track->GetConstrainedParam()->Phi(), mass);
314  }
315  else {
316  mom.SetPtEtaPhiM(vp->Pt(), vp->Eta(), vp->Phi(), mass);
317  }
318 
319  return kTRUE;
320  }
321  else {
322  mom.SetPtEtaPhiM(0, 0, 0, 0);
323  return kFALSE;
324  }
325 }
326 
327 //________________________________________________________________________
328 Bool_t AliTrackContainer::GetNextAcceptMomentum(TLorentzVector &mom)
329 {
330  //Get momentum of the next accepted particle in array
331 
332  Double_t mass = fMassHypothesis;
333 
334  AliVTrack *vp = GetNextAcceptTrack();
335  if (vp) {
336  if (mass < 0) mass = vp->M();
337 
338  if (fLoadedClass->InheritsFrom("AliESDtrack") && fTrackFilterType == AliEmcalTrackSelection::kHybridTracks &&
340  AliESDtrack *track = static_cast<AliESDtrack*>(vp);
341  mom.SetPtEtaPhiM(track->GetConstrainedParam()->Pt(), track->GetConstrainedParam()->Eta(), track->GetConstrainedParam()->Phi(), mass);
342  }
343  else {
344  mom.SetPtEtaPhiM(vp->Pt(), vp->Eta(), vp->Phi(), mass);
345  }
346 
347  return kTRUE;
348  }
349  else {
350  mom.SetPtEtaPhiM(0, 0, 0, 0);
351  return kFALSE;
352  }
353 }
354 
355 //________________________________________________________________________
356 Bool_t AliTrackContainer::AcceptTrack(const AliVTrack *vp)
357 {
358  // Return true if vp is accepted.
359  Bool_t r = ApplyTrackCuts(vp);
360  if (!r) return kFALSE;
361 
362  AliTLorentzVector mom;
363 
364  Int_t id = fFilteredTracks->IndexOf(vp);
365  if (id >= 0) {
366  GetMomentum(mom, id);
367  }
368  else {
369  GetMomentum(mom, vp);
370  }
371 
372  return ApplyKinematicCuts(mom);
373 }
374 
375 //________________________________________________________________________
377 {
378  // Return true if vp is accepted.
379  Bool_t r = ApplyTrackCuts(GetTrack(i));
380  if (!r) return kFALSE;
381 
382  AliTLorentzVector mom;
383  GetMomentum(mom, i);
384 
385  return ApplyKinematicCuts(mom);
386 }
387 
388 //________________________________________________________________________
389 Bool_t AliTrackContainer::ApplyTrackCuts(const AliVTrack* vp)
390 {
391  // Return true if i^th particle is accepted.
392 
393  return ApplyParticleCuts(vp);
394 }
395 
396 //________________________________________________________________________
397 void AliTrackContainer::SetClassName(const char *clname)
398 {
399  // Set the class name
400 
401  TClass cls(clname);
402  if (cls.InheritsFrom("AliVTrack")) fClassName = clname;
403  else AliError(Form("Unable to set class name %s for a AliTrackContainer, it must inherits from AliVTrack!",clname));
404 }
405 
406 //________________________________________________________________________
408 {
409  if (!fListOfCuts) {
410  fListOfCuts = new TObjArray;
411  fListOfCuts->SetOwner(true);
412  }
413  fListOfCuts->Add(cuts);
414 }
415 
416 //________________________________________________________________________
418 {
419  if (!fListOfCuts) return 0;
420  return fListOfCuts->GetEntries();
421 }
422 
423 //________________________________________________________________________
424 AliVCuts* AliTrackContainer::GetTrackCuts(Int_t icut)
425 {
426  if (!fListOfCuts) return NULL;
427  if (icut < fListOfCuts->GetEntries()) {
428  return static_cast<AliVCuts *>(fListOfCuts->At(icut));
429  }
430  return NULL;
431 }
432 
433 //________________________________________________________________________
434 const char* AliTrackContainer::GetTitle() const
435 {
436  static TString trackString;
437 
438  if (GetMinPt() == 0) {
439  trackString = TString::Format("%s_pT0000", GetArrayName().Data());
440  }
441  else if (GetMinPt() < 1.0) {
442  trackString = TString::Format("%s_pT0%3.0f", GetArrayName().Data(), GetMinPt()*1000.0);
443  }
444  else {
445  trackString = TString::Format("%s_pT%4.0f", GetArrayName().Data(), GetMinPt()*1000.0);
446  }
447 
448  return trackString.Data();
449 }
Interface for virtual track selection.
ClassImp(AliAnalysisTaskTriggerRates) AliAnalysisTaskTriggerRates
TObjArray * fFilteredTracks
track selection object
void SetClassName(const char *clname)
virtual Bool_t ApplyTrackCuts(const AliVTrack *vp)
static TString fgDefTrackCutsPeriod
Double_t mass
virtual Bool_t AcceptTrack(const AliVTrack *vp)
TString fClassName
name of branch
virtual void SetArray(AliVEvent *event)
const TClonesArray * GetAcceptedTrackBitmaps() const
virtual Bool_t ApplyParticleCuts(const AliVParticle *vp)
AliVCuts * GetTrackCuts(Int_t icut)
virtual Bool_t GetNextMomentum(TLorentzVector &mom)
AliEmcalTrackSelection * fEmcalTrackSelection
TClass * fLoadedClass
!Class of teh objects contained in the TClonesArray
ETrackFilterType_t fTrackFilterType
default period string used to generate track cuts
void SetArray(AliVEvent *event)
Container for particles within the EMCAL framework.
Double_t GetMinPt() const
Int_t GetNumberOfCutObjects() const
TArrayC fTrackTypes
tracks filtered using fEmcalTrackSelection
TObjArray * fListOfCuts
void AddTrackCuts(AliVCuts *cuts)
Implement virtual track selection for AOD analysis.
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)
const TString & GetArrayName() const
Double_t fMassHypothesis
maximum MC label
virtual AliVTrack * GetTrack(Int_t i=-1) const
TObjArray * GetAcceptedTracks(const TClonesArray *const tracks)
virtual AliVTrack * GetNextAcceptTrack()
virtual Bool_t GetMomentum(TLorentzVector &mom, const AliVTrack *part, Double_t mass)
virtual Bool_t GetAcceptMomentum(TLorentzVector &mom, Int_t i)
virtual Bool_t ApplyKinematicCuts(const AliTLorentzVector &mom)
void AddTrackCuts(AliVCuts *cuts)
virtual Bool_t GetNextAcceptMomentum(TLorentzVector &mom)
const char * GetTitle() const
virtual AliVTrack * GetAcceptTrack(Int_t i=-1)
TClonesArray * fClArray
if < 0 it will use a PID mass when available
Int_t GetNEntries() const
Int_t fCurrentID
!current ID for automatic loops
Implementation of virtual track selection for ESDs.