AliPhysics  8bb951a (8bb951a)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliEmcalJetTask.cxx
Go to the documentation of this file.
1 //
2 // Emcal jet finder task.
3 //
4 // Authors: C.Loizides, S.Aiola, M. Verweij
5 
6 #include <vector>
7 #include "AliEmcalJetTask.h"
8 
9 #include <TChain.h>
10 #include <TClonesArray.h>
11 #include <TList.h>
12 #include <TMath.h>
13 #include <TRandom3.h>
14 #include <TClass.h>
15 
16 #include "AliTLorentzVector.h"
17 #include "AliAnalysisManager.h"
18 #include "AliCentrality.h"
19 #include "AliEMCALGeometry.h"
20 #include "AliESDEvent.h"
21 #include "AliEmcalJet.h"
22 #include "AliEmcalParticle.h"
23 #include "AliFJWrapper.h"
24 #include "AliVCluster.h"
25 #include "AliVEvent.h"
26 #include "AliVParticle.h"
27 #include "AliEmcalJetUtility.h"
28 #include "AliParticleContainer.h"
29 #include "AliClusterContainer.h"
30 
31 using std::cout;
32 using std::endl;
33 using std::cerr;
34 
36 
37 const Int_t AliEmcalJetTask::fgkConstIndexShift = 100000;
38 
39 //________________________________________________________________________
42  fJetsTag(),
43  fJetAlgo(AliJetContainer::antikt_algorithm),
44  fJetType(AliJetContainer::kFullJet),
45  fRecombScheme(AliJetContainer::pt_scheme),
46  fRadius(0.4),
47  fMinJetArea(0.001),
48  fMinJetPt(1.0),
49  fJetPhiMin(-10),
50  fJetPhiMax(+10),
51  fJetEtaMin(-1),
52  fJetEtaMax(+1),
53  fGhostArea(0.005),
54  fTrackEfficiency(1.),
55  fUtilities(0),
56  fLocked(0),
57  fJetsName(),
58  fIsInit(0),
59  fIsPSelSet(0),
60  fIsEmcPart(0),
61  fLegacyMode(kFALSE),
62  fFillGhost(kFALSE),
63  fJets(0),
64  fFastJetWrapper("AliEmcalJetTask","AliEmcalJetTask")
65 {
66  // Default constructor.
67 }
68 
69 //________________________________________________________________________
72  fJetsTag("Jets"),
73  fJetAlgo(AliJetContainer::antikt_algorithm),
74  fJetType(AliJetContainer::kFullJet),
75  fRecombScheme(AliJetContainer::pt_scheme),
76  fRadius(0.4),
77  fMinJetArea(0.001),
78  fMinJetPt(1.0),
79  fJetPhiMin(-10),
80  fJetPhiMax(+10),
81  fJetEtaMin(-1),
82  fJetEtaMax(+1),
83  fGhostArea(0.005),
84  fTrackEfficiency(1.),
85  fUtilities(0),
86  fLocked(0),
87  fJetsName(),
88  fIsInit(0),
89  fIsPSelSet(0),
90  fIsEmcPart(0),
91  fLegacyMode(kFALSE),
92  fFillGhost(kFALSE),
93  fJets(0),
94  fFastJetWrapper(name,name)
95 {
96  // Standard constructor.
97 
98  fBranchNames="ESD:AliESDRun.,AliESDHeader.,PrimaryVertex.";
99 }
100 
101 //________________________________________________________________________
103 {
104  // Destructor
105 }
106 
107 //________________________________________________________________________
109 {
110  // Addition of utilities.
111 
112  if (!fUtilities) fUtilities = new TObjArray();
113  if (fUtilities->FindObject(utility)) {
114  Error("AddSupply", "Jet utility %s already connected.", utility->GetName());
115  return utility;
116  }
117  fUtilities->Add(utility);
118  utility->SetJetTask(this);
119 
120  return utility;
121 }
122 
123 //________________________________________________________________________
125 {
126  TIter next(fUtilities);
127  AliEmcalJetUtility *utility = 0;
128  while ((utility=static_cast<AliEmcalJetUtility*>(next()))) utility->Init();
129 }
130 
131 //________________________________________________________________________
133 {
134  TIter next(fUtilities);
135  AliEmcalJetUtility *utility = 0;
136  while ((utility=static_cast<AliEmcalJetUtility*>(next()))) utility->Prepare(fFastJetWrapper);
137 }
138 
139 //________________________________________________________________________
141 {
142  TIter next(fUtilities);
143  AliEmcalJetUtility *utility = 0;
144  while ((utility=static_cast<AliEmcalJetUtility*>(next()))) utility->ProcessJet(jet, ij, fFastJetWrapper);
145 }
146 
147 //________________________________________________________________________
149 {
150  TIter next(fUtilities);
151  AliEmcalJetUtility *utility = 0;
152  while ((utility=static_cast<AliEmcalJetUtility*>(next()))) utility->Terminate(fFastJetWrapper);
153 }
154 
155 //________________________________________________________________________
157 {
158  // Main loop, called for each event.
159 
160  // clear the jet array (normally a null operation)
161  fJets->Delete();
162 
163  Int_t n = FindJets();
164 
165  if (n == 0) return kFALSE;
166 
167  FillJetBranch();
168 
169  return kTRUE;
170 }
171 
172 //________________________________________________________________________
174 {
175  // Find jets.
176 
177 
178  if (fParticleCollArray.GetEntriesFast() == 0 && fClusterCollArray.GetEntriesFast() == 0){
179  AliError("No tracks or clusters, returning.");
180  return 0;
181  }
182 
184 
185  AliDebug(2,Form("Jet type = %d", fJetType));
186 
187  AliTLorentzVector mom;
188 
189  Int_t iColl = 1;
190  TIter nextPartColl(&fParticleCollArray);
191  AliParticleContainer* tracks = 0;
192  while ((tracks = static_cast<AliParticleContainer*>(nextPartColl()))) {
193  AliDebug(2,Form("Tracks from collection %d: '%s'.", iColl-1, tracks->GetName()));
194  tracks->ResetCurrentID();
195  AliVParticle* t = 0;
196  while ((t = tracks->GetNextAcceptParticle())) {
197  tracks->GetMomentum(mom, tracks->GetCurrentID());
198  if (((fJetType & AliJetContainer::kChargedJet) != 0) && (t->Charge() == 0)) {
199  AliDebug(2,Form("Skipping track %d because it is neutral.", tracks->GetCurrentID()));
200  continue;
201  }
202 
203  if (((fJetType & AliJetContainer::kNeutralJet) != 0) && (t->Charge() != 0)) {
204  AliDebug(2,Form("Skipping track %d because it is charged.", tracks->GetCurrentID()));
205  continue;
206  }
207 
208  // artificial inefficiency
209  if (fTrackEfficiency < 1.) {
210  Double_t rnd = gRandom->Rndm();
211  if (fTrackEfficiency < rnd) {
212  AliDebug(2,Form("Track %d rejected due to artificial tracking inefficiency", tracks->GetCurrentID()));
213  continue;
214  }
215  }
216 
217  AliDebug(2,Form("Track %d accepted (label = %d, pt = %f)", tracks->GetCurrentID(), t->GetLabel(), t->Pt()));
218  Int_t uid = tracks->GetCurrentID() + fgkConstIndexShift * iColl;
219  fFastJetWrapper.AddInputVector(mom.Px(), mom.Py(), mom.Pz(), mom.E(), uid);
220  }
221  iColl++;
222  }
223 
224  iColl = 1;
225  TIter nextClusColl(&fClusterCollArray);
226  AliClusterContainer* clusters = 0;
227  while ((clusters = static_cast<AliClusterContainer*>(nextClusColl()))) {
228  AliDebug(2,Form("Clusters from collection %d: '%s'.", iColl-1, clusters->GetName()));
229  clusters->ResetCurrentID();
230  AliVCluster* c = 0;
231  while ((c = clusters->GetNextAcceptCluster())) {
232  clusters->GetMomentum(mom, clusters->GetCurrentID());
233  Double_t cEta = mom.Eta();
234  Double_t cPhi = mom.Phi_0_2pi();
235  Double_t cPt = mom.Pt();
236  Double_t cPx = mom.Px();
237  Double_t cPy = mom.Py();
238  Double_t cPz = mom.Pz();
239 
240  Double_t e = TMath::Sqrt(cPx*cPx+cPy*cPy+cPz*cPz);
241 
242  AliDebug(2,Form("Cluster %d accepted (label = %d, energy = %.3f)", clusters->GetCurrentID(), c->GetLabel(), e));
243  Int_t uid = -clusters->GetCurrentID() - fgkConstIndexShift * iColl;
244  fFastJetWrapper.AddInputVector(cPx, cPy, cPz, e, uid);
245  }
246  iColl++;
247  }
248 
249  if (fFastJetWrapper.GetInputVectors().size() == 0) return 0;
250 
251  // run jet finder
253 
254  return fFastJetWrapper.GetInclusiveJets().size();
255 }
256 
257 //________________________________________________________________________
259 {
260  // Fill jet branch with fastjet jets.
261 
263 
264  // loop over fastjet jets
265  std::vector<fastjet::PseudoJet> jets_incl = fFastJetWrapper.GetInclusiveJets();
266  // sort jets according to jet pt
267  static Int_t indexes[9999] = {-1};
268  GetSortedArray(indexes, jets_incl);
269 
270  AliDebug(1,Form("%d jets found", (Int_t)jets_incl.size()));
271  for (UInt_t ijet = 0, jetCount = 0; ijet < jets_incl.size(); ++ijet) {
272  Int_t ij = indexes[ijet];
273  AliDebug(3,Form("Jet pt = %f, area = %f", jets_incl[ij].perp(), fFastJetWrapper.GetJetArea(ij)));
274 
275  if (jets_incl[ij].perp() < fMinJetPt) continue;
276  if (fFastJetWrapper.GetJetArea(ij) < fMinJetArea) continue;
277  if ((jets_incl[ij].eta() < fJetEtaMin) || (jets_incl[ij].eta() > fJetEtaMax) ||
278  (jets_incl[ij].phi() < fJetPhiMin) || (jets_incl[ij].phi() > fJetPhiMax))
279  continue;
280 
281  AliEmcalJet *jet = new ((*fJets)[jetCount])
282  AliEmcalJet(jets_incl[ij].perp(), jets_incl[ij].eta(), jets_incl[ij].phi(), jets_incl[ij].m());
283  jet->SetLabel(ij);
284 
285  fastjet::PseudoJet area(fFastJetWrapper.GetJetAreaVector(ij));
286  jet->SetArea(area.perp());
287  jet->SetAreaEta(area.eta());
288  jet->SetAreaPhi(area.phi());
289  jet->SetAreaE(area.E());
290 
291  // Fill constituent info
292  std::vector<fastjet::PseudoJet> constituents(fFastJetWrapper.GetJetConstituents(ij));
293  FillJetConstituents(jet, constituents, constituents);
294 
295  if (fGeom) {
296  if ((jet->Phi() > fGeom->GetArm1PhiMin() * TMath::DegToRad()) &&
297  (jet->Phi() < fGeom->GetArm1PhiMax() * TMath::DegToRad()) &&
298  (jet->Eta() > fGeom->GetArm1EtaMin()) &&
299  (jet->Eta() < fGeom->GetArm1EtaMax()))
300  jet->SetAxisInEmcal(kTRUE);
301  }
302 
303  ExecuteUtilities(jet, ij);
304 
305  AliDebug(2,Form("Added jet n. %d, pt = %f, area = %f, constituents = %d", jetCount, jet->Pt(), jet->Area(), jet->GetNumberOfConstituents()));
306  jetCount++;
307  }
308 
310 }
311 
312 //________________________________________________________________________
313 Bool_t AliEmcalJetTask::GetSortedArray(Int_t indexes[], std::vector<fastjet::PseudoJet> array) const
314 {
315  // Get the leading jets.
316 
317  static Float_t pt[9999] = {0};
318 
319  const Int_t n = (Int_t)array.size();
320 
321  if (n < 1)
322  return kFALSE;
323 
324  for (Int_t i = 0; i < n; i++)
325  pt[i] = array[i].perp();
326 
327  TMath::Sort(n, pt, indexes);
328 
329  return kTRUE;
330 }
331 
332 //________________________________________________________________________
334 {
335  // Init the task.
336 
337  SetNeedEmcalGeom(kFALSE);
338 
339  if (fTrackEfficiency < 1.) {
340  if (gRandom) delete gRandom;
341  gRandom = new TRandom3(0);
342  }
343 
345 
346  // add jets to event if not yet there
347  if (!(InputEvent()->FindListObject(fJetsName))) {
348  fJets = new TClonesArray("AliEmcalJet");
349  fJets->SetName(fJetsName);
350  ::Info("AliEmcalJetTask::ExecOnce", "Jet collection with name '%s' has been added to the event.", fJetsName.Data());
351  InputEvent()->AddObject(fJets);
352  }
353  else {
354  AliError(Form("%s: Object with name %s already in event! Returning", GetName(), fJetsName.Data()));
355  return;
356  }
357 
358  // setup fj wrapper
359  fFastJetWrapper.SetAreaType(fastjet::active_area_explicit_ghosts);
365 
366  // setting legacy mode
367  if (fLegacyMode) {
369  }
370 
371  InitUtilities();
372 
374 }
375 
376 //___________________________________________________________________________________________________________________
377 void AliEmcalJetTask::FillJetConstituents(AliEmcalJet *jet, std::vector<fastjet::PseudoJet>& constituents,
378  std::vector<fastjet::PseudoJet>& constituents_unsub, Int_t flag, TClonesArray *particles_sub)
379 {
380  Int_t nt = 0;
381  Int_t nc = 0;
382  Double_t neutralE = 0.;
383  Double_t maxCh = 0.;
384  Double_t maxNe = 0.;
385  Int_t gall = 0;
386  Int_t gemc = 0;
387  Int_t cemc = 0;
388  Int_t ncharged = 0;
389  Int_t nneutral = 0;
390  Double_t mcpt = 0.;
391  Double_t emcpt = 0.;
392 
393  Int_t uid = -1;
394 
395  jet->SetNumberOfTracks(constituents.size());
396  jet->SetNumberOfClusters(constituents.size());
397 
398  for (UInt_t ic = 0; ic < constituents.size(); ++ic) {
399 
400  if (flag == 0) {
401  uid = constituents[ic].user_index();
402  }
403  else {
404  if (constituents[ic].perp()<1.e-10) {
405  uid=-1;
406  }
407  else {
408  uid = GetIndexSub(constituents[ic].phi(), constituents_unsub);
409  }
410  if (uid==0) {
411  AliError("correspondence between un/subtracted constituent not found");
412  continue;
413  }
414  }
415 
416  AliDebug(3,Form("Processing constituent %d", uid));
417  if (uid == -1) { //ghost particle
418  ++gall;
419  if (fGeom) {
420  Double_t gphi = constituents[ic].phi();
421  if (gphi < 0) gphi += TMath::TwoPi();
422  gphi *= TMath::RadToDeg();
423  Double_t geta = constituents[ic].eta();
424  if ((gphi > fGeom->GetArm1PhiMin()) && (gphi < fGeom->GetArm1PhiMax()) &&
425  (geta > fGeom->GetArm1EtaMin()) && (geta < fGeom->GetArm1EtaMax()))
426  ++gemc;
427  }
428 
429  if (fFillGhost) jet->AddGhost(constituents[ic].px(),
430  constituents[ic].py(),
431  constituents[ic].pz(),
432  constituents[ic].e());
433  }
434  else if (uid >= fgkConstIndexShift) { // track constituent
435  Int_t iColl = uid / fgkConstIndexShift;
436  Int_t tid = uid - iColl * fgkConstIndexShift;
437  iColl--;
438  AliDebug(3,Form("Constituent %d is a track from collection %d and with ID %d", uid, iColl, tid));
439  AliParticleContainer* partCont = GetParticleContainer(iColl);
440  if (!partCont) {
441  AliError(Form("Could not find particle container %d",iColl));
442  continue;
443  }
444  AliVParticle *t = partCont->GetParticle(tid);
445  if (!t) {
446  AliError(Form("Could not find track %d",tid));
447  continue;
448  }
449 
450  Double_t cEta = t->Eta();
451  Double_t cPhi = t->Phi();
452  Double_t cPt = t->Pt();
453  Double_t cP = t->P();
454  if (t->Charge() == 0) {
455  neutralE += cP;
456  ++nneutral;
457  if (cPt > maxNe) maxNe = cPt;
458  } else {
459  ++ncharged;
460  if (cPt > maxCh) maxCh = cPt;
461  }
462 
463  // check if MC particle
464  if (TMath::Abs(t->GetLabel()) > fMinMCLabel) mcpt += cPt;
465 
466  if (fGeom) {
467  if (cPhi < 0) cPhi += TMath::TwoPi();
468  cPhi *= TMath::RadToDeg();
469  if ((cPhi > fGeom->GetArm1PhiMin()) && (cPhi < fGeom->GetArm1PhiMax()) &&
470  (cEta > fGeom->GetArm1EtaMin()) && (cEta < fGeom->GetArm1EtaMax())) {
471  emcpt += cPt;
472  ++cemc;
473  }
474  }
475 
476  if (flag == 0 || particles_sub == 0) {
477  jet->AddTrackAt(tid, nt);
478  }
479  else {
480  Int_t part_sub_id = particles_sub->GetEntriesFast();
481  AliEmcalParticle* part_sub = new ((*particles_sub)[part_sub_id]) AliEmcalParticle(dynamic_cast<AliVTrack*>(t)); // SA: probably need to be fixed!!
482  part_sub->SetPtEtaPhiM(constituents[ic].perp(),constituents[ic].eta(),constituents[ic].phi(),constituents[ic].m());
483  jet->AddTrackAt(part_sub_id, nt);
484  }
485 
486  ++nt;
487  }
488  else if (uid <= -fgkConstIndexShift) { // cluster constituent
489  Int_t iColl = -uid / fgkConstIndexShift;
490  Int_t cid = -uid - iColl * fgkConstIndexShift;
491  iColl--;
492  AliClusterContainer* clusCont = GetClusterContainer(iColl);
493  AliVCluster *c = clusCont->GetCluster(cid);
494 
495  if (!c) continue;
496 
497  AliTLorentzVector nP;
498  clusCont->GetMomentum(nP, cid);
499 
500  Double_t cEta = nP.Eta();
501  Double_t cPhi = nP.Phi_0_2pi();
502  Double_t cPt = nP.Pt();
503  Double_t cP = nP.P();
504 
505  neutralE += cP;
506  if (cPt > maxNe) maxNe = cPt;
507 
508  // MC particle
509  if (TMath::Abs(c->GetLabel()) > fMinMCLabel) mcpt += c->GetMCEnergyFraction() > 1e-6 ? cPt * c->GetMCEnergyFraction() : cPt;
510 
511  if (fGeom) {
512  if (cPhi<0) cPhi += TMath::TwoPi();
513  cPhi *= TMath::RadToDeg();
514  if ((cPhi > fGeom->GetArm1PhiMin()) && (cPhi < fGeom->GetArm1PhiMax()) &&
515  (cEta > fGeom->GetArm1EtaMin()) && (cEta < fGeom->GetArm1EtaMax())) {
516  emcpt += cPt;
517  ++cemc;
518  }
519  }
520 
521  if (flag == 0 || particles_sub == 0) {
522  jet->AddClusterAt(cid, nc);
523  }
524  else {
525  Int_t part_sub_id = particles_sub->GetEntriesFast();
526  AliEmcalParticle* part_sub = new ((*particles_sub)[part_sub_id]) AliEmcalParticle(c);
527  part_sub->SetPtEtaPhiM(constituents[ic].perp(),constituents[ic].eta(),constituents[ic].phi(),constituents[ic].m());
528  jet->AddTrackAt(part_sub_id, nt);
529  }
530 
531  ++nc;
532  ++nneutral;
533  }
534  else {
535  AliError(Form("%s: No logical way to end up here (uid = %d).", GetName(), uid));
536  continue;
537  }
538  }
539 
540  jet->SetNumberOfTracks(nt);
541  jet->SetNumberOfClusters(nc);
542  jet->SetNEF(neutralE / jet->E());
543  jet->SetMaxChargedPt(maxCh);
544  jet->SetMaxNeutralPt(maxNe);
545  if (gall > 0) jet->SetAreaEmc(jet->Area() * gemc / gall);
546  else jet->SetAreaEmc(-1);
547  jet->SetNEmc(cemc);
548  jet->SetNumberOfCharged(ncharged);
549  jet->SetNumberOfNeutrals(nneutral);
550  jet->SetMCPt(mcpt);
551  jet->SetPtEmc(emcpt);
552  jet->SortConstituents();
553 }
554 
555 //______________________________________________________________________________________
556 Int_t AliEmcalJetTask::GetIndexSub(Double_t phi_sub, std::vector<fastjet::PseudoJet>& constituents_unsub)
557 {
558  Double_t dphi=0;
559  for (UInt_t ii = 0; ii < constituents_unsub.size(); ii++) {
560  Double_t phi_unsub = constituents_unsub[ii].phi();
561  dphi=phi_unsub-phi_sub;
562  if (dphi < -1*TMath::Pi()) dphi += (2*TMath::Pi());
563  else if (dphi > TMath::Pi()) dphi -= (2*TMath::Pi());
564  if (TMath::Abs(dphi)<0.1 && constituents_unsub[ii].user_index()!=-1) return constituents_unsub[ii].user_index();
565  }
566 
567  return 0;
568 }
569 
570 //______________________________________________________________________________________
572 {
573  if (fLocked) {
574  AliFatal("Jet finder task is locked! Changing properties is not allowed.");
575  return kTRUE;
576  }
577  else {
578  return kFALSE;
579  }
580 }
581 
582 //______________________________________________________________________________________
583 void AliEmcalJetTask::SelectCollisionCandidates(UInt_t offlineTriggerMask)
584 {
585  if(!fIsPSelSet) {
586  fIsPSelSet = kTRUE;
587  fOfflineTriggerMask = offlineTriggerMask;
588  }
589  else {
590  AliWarning("Manually setting the event selection for jet finders is not allowed! Using trigger=old_trigger | your_trigger");
591  fOfflineTriggerMask = fOfflineTriggerMask | offlineTriggerMask;
592  }
593 }
594 
595 //______________________________________________________________________________________
596 void AliEmcalJetTask::SetEtaRange(Double_t emi, Double_t ema)
597 {
598  if (IsLocked()) return;
599 
600  TIter nextPartColl(&fParticleCollArray);
601  AliParticleContainer* tracks = 0;
602  while ((tracks = static_cast<AliParticleContainer*>(nextPartColl()))) {
603  tracks->SetParticleEtaLimits(emi, ema);
604  }
605 }
606 
607 //______________________________________________________________________________________
609 {
610  if (IsLocked()) return;
611 
612  TIter nextClusColl(&fClusterCollArray);
613  AliClusterContainer* clusters = 0;
614  while ((clusters = static_cast<AliClusterContainer*>(nextClusColl()))) {
615  clusters->SetClusPtCut(min);
616  }
617 }
618 
619 //______________________________________________________________________________________
621 {
622  if (IsLocked()) return;
623 
624  TIter nextClusColl(&fClusterCollArray);
625  AliClusterContainer* clusters = 0;
626  while ((clusters = static_cast<AliClusterContainer*>(nextClusColl()))) {
627  clusters->SetClusECut(min);
628  }
629 }
630 
631 //______________________________________________________________________________________
633 {
634  if (IsLocked()) return;
635 
636  TIter nextPartColl(&fParticleCollArray);
637  AliParticleContainer* tracks = 0;
638  while ((tracks = static_cast<AliParticleContainer*>(nextPartColl()))) {
639  tracks->SetParticlePtCut(min);
640  }
641 }
642 
643 //______________________________________________________________________________________
644 void AliEmcalJetTask::SetPhiRange(Double_t pmi, Double_t pma)
645 {
646  if (IsLocked()) return;
647 
648  TIter nextPartColl(&fParticleCollArray);
649  AliParticleContainer* tracks = 0;
650  while ((tracks = static_cast<AliParticleContainer*>(nextPartColl()))) {
651  tracks->SetParticlePhiLimits(pmi, pma);
652  }
653 }
654 
655 //______________________________________________________________________________________
656 fastjet::JetAlgorithm AliEmcalJetTask::ConvertToFJAlgo(EJetAlgo_t algo)
657 {
658  switch(algo) {
660  return fastjet::kt_algorithm;
662  return fastjet::antikt_algorithm;
664  return fastjet::cambridge_algorithm;
666  return fastjet::genkt_algorithm;
668  return fastjet::cambridge_for_passive_algorithm;
670  return fastjet::genkt_for_passive_algorithm;
672  return fastjet::plugin_algorithm;
674  return fastjet::undefined_jet_algorithm;
675 
676  default:
677  ::Error("AliEmcalJetTask::ConvertToFJAlgo", "Jet algorithm %d not recognized!!!", algo);
678  return fastjet::undefined_jet_algorithm;
679  }
680 }
681 
682 //______________________________________________________________________________________
683 fastjet::RecombinationScheme AliEmcalJetTask::ConvertToFJRecoScheme(ERecoScheme_t reco)
684 {
685  switch(reco) {
687  return fastjet::E_scheme;
688 
690  return fastjet::pt_scheme;
691 
693  return fastjet::pt2_scheme;
694 
696  return fastjet::Et_scheme;
697 
699  return fastjet::Et2_scheme;
700 
702  return fastjet::BIpt_scheme;
703 
705  return fastjet::BIpt2_scheme;
706 
708  return fastjet::external_scheme;
709 
710  default:
711  ::Error("AliEmcalJetTask::ConvertToFJRecoScheme", "Recombination scheme %d not recognized!!!", reco);
712  return fastjet::external_scheme;
713  }
714 }
void SetMaxNeutralPt(Double32_t t)
Definition: AliEmcalJet.h:132
fastjet::PseudoJet GetJetAreaVector(UInt_t idx) const
ClassImp(AliAnalysisTaskTriggerRates) AliAnalysisTaskTriggerRates
Double_t Area() const
Definition: AliEmcalJet.h:69
TObjArray fClusterCollArray
cluster collection array
void AddTrackAt(Int_t track, Int_t idx)
Definition: AliEmcalJet.h:115
void SetParticlePtCut(Double_t cut)
virtual AliVParticle * GetNextAcceptParticle()
TClonesArray * fJets
=true ghost particles will be filled in AliEmcalJet obj
virtual Int_t Run()
void SetEtaRange(Double_t emi, Double_t ema)
Bool_t IsLocked() const
virtual Bool_t GetMomentum(TLorentzVector &mom, Int_t i) const
Double_t Eta() const
Definition: AliEmcalJet.h:60
Base task in the EMCAL framework.
EJetType_t fJetType
Double_t Phi() const
Definition: AliEmcalJet.h:55
void FillJetConstituents(AliEmcalJet *jet, std::vector< fastjet::PseudoJet > &constituents, std::vector< fastjet::PseudoJet > &constituents_sub, Int_t flag=0, TClonesArray *particles_sub=0)
void AddClusterAt(Int_t clus, Int_t idx)
Definition: AliEmcalJet.h:114
Bool_t fFillGhost
=true to enable FJ 2.x behavior
AliEmcalJetUtility * AddUtility(AliEmcalJetUtility *utility)
virtual void ProcessJet(AliEmcalJet *jet, Int_t ij, AliFJWrapper &fjw)=0
void SetArea(Double_t a)
Definition: AliEmcalJet.h:126
void ExecuteUtilities(AliEmcalJet *jet, Int_t ij)
void SetMaxRap(Double_t maxrap)
Definition: AliFJWrapper.h:86
virtual void Terminate(AliFJWrapper &fjw)=0
void SetAreaE(Double_t a)
Definition: AliEmcalJet.h:129
ERecoScheme_t fRecombScheme
void SetRecombScheme(const fastjet::RecombinationScheme &scheme)
Definition: AliFJWrapper.h:82
Double_t fTrackEfficiency
void SetMCPt(Double_t p)
Definition: AliEmcalJet.h:139
Double_t E() const
Definition: AliEmcalJet.h:58
virtual ~AliEmcalJetTask()
TRandom * gRandom
static FJRecoScheme ConvertToFJRecoScheme(ERecoScheme_t reco)
void SetLabel(Int_t l)
Definition: AliEmcalJet.h:125
static TString GenerateJetName(EJetType_t jetType, EJetAlgo_t jetAlgo, ERecoScheme_t recoScheme, Double_t radius, AliParticleContainer *partCont, AliClusterContainer *clusCont, TString tag)
UShort_t GetNumberOfConstituents() const
Definition: AliEmcalJet.h:84
Container for particles within the EMCAL framework.
Bool_t fIsPSelSet
=true if already initialized
void SetR(Double_t r)
Definition: AliFJWrapper.h:87
TObjArray fParticleCollArray
particle/track collection array
AliParticleContainer * GetParticleContainer(Int_t i=0) const
const std::vector< fastjet::PseudoJet > & GetInclusiveJets() const
Definition: AliFJWrapper.h:32
Double_t GetJetArea(UInt_t idx) const
Int_t GetCurrentID() const
AliEMCALGeometry * fGeom
!emcal geometry
void SetLegacyMode(Bool_t mode)
Definition: AliFJWrapper.h:97
virtual AliVParticle * GetParticle(Int_t i=-1) const
TObjArray * fUtilities
void SetPtEtaPhiM(Double_t pt, Double_t eta, Double_t phi, Double_t m)
void SetAlgorithm(const fastjet::JetAlgorithm &algor)
Definition: AliFJWrapper.h:81
AliClusterContainer * GetClusterContainer(Int_t i=0) const
void SetMaxChargedPt(Double32_t t)
Definition: AliEmcalJet.h:133
virtual void Clear(const Option_t *="")
void SetNEF(Double_t nef)
Definition: AliEmcalJet.h:134
Bool_t GetSortedArray(Int_t indexes[], std::vector< fastjet::PseudoJet > array) const
void AddGhost(const Double_t dPx, const Double_t dPy, const Double_t dPz, const Double_t dE)
Definition: AliEmcalJet.h:269
void SetNEmc(Int_t n)
Definition: AliEmcalJet.h:142
AliVCluster * GetCluster(Int_t i) const
virtual void Prepare(AliFJWrapper &fjw)=0
Int_t fMinMCLabel
minimum MC label value for the tracks/clusters being considered MC particles
void SetJetTask(AliEmcalJetTask *jetTask)
void SetMinJetClusE(Double_t min)
Double_t Pt() const
Definition: AliEmcalJet.h:47
void SelectCollisionCandidates(UInt_t offlineTriggerMask=AliVEvent::kMB)
void SetGhostArea(Double_t gharea)
Definition: AliFJWrapper.h:85
EJetAlgo_t fJetAlgo
void SetNumberOfCharged(Int_t n)
Definition: AliEmcalJet.h:137
std::vector< fastjet::PseudoJet > GetJetConstituents(UInt_t idx) const
const char * GetName() const
void SetAreaEta(Double_t a)
Definition: AliEmcalJet.h:127
static FJJetAlgo ConvertToFJAlgo(EJetAlgo_t algo)
void SetParticleEtaLimits(Double_t min, Double_t max)
const std::vector< fastjet::PseudoJet > & GetInputVectors() const
Definition: AliFJWrapper.h:30
void SetMinJetTrackPt(Double_t min)
void SetNeedEmcalGeom(Bool_t n)
Bool_t GetMomentum(TLorentzVector &mom, const AliVCluster *vc, Double_t mass) const
void SetNumberOfTracks(Int_t n)
Definition: AliEmcalJet.h:136
Bool_t fLegacyMode
=true if emcal particles are given as input (for clusters)
void SetClusPtCut(Double_t cut)
virtual void AddInputVector(Double_t px, Double_t py, Double_t pz, Double_t E, Int_t index=-99999)
void SortConstituents()
void SetMinJetClusPt(Double_t min)
Int_t GetIndexSub(Double_t phi_sub, std::vector< fastjet::PseudoJet > &constituents_unsub)
void SetPhiRange(Double_t pmi, Double_t pma)
void SetClusECut(Double_t cut)
void SetAreaPhi(Double_t a)
Definition: AliEmcalJet.h:128
void SetAreaType(const fastjet::AreaType &atype)
Definition: AliFJWrapper.h:83
void SetParticlePhiLimits(Double_t min, Double_t max)
void SetPtEmc(Double_t pt)
Definition: AliEmcalJet.h:143
Container structure for EMCAL clusters.
AliVCluster * GetNextAcceptCluster()
void ResetCurrentID(Int_t i=-1)
AliFJWrapper fFastJetWrapper
jet collection
void SetNumberOfClusters(Int_t n)
Definition: AliEmcalJet.h:135
virtual void Init()=0
void SetAxisInEmcal(Bool_t b)
Definition: AliEmcalJet.h:131
static const Int_t fgkConstIndexShift
fastjet wrapper
void SetAreaEmc(Double_t a)
Definition: AliEmcalJet.h:130
void SetNumberOfNeutrals(Int_t n)
Definition: AliEmcalJet.h:138