3 // Standalone jet finder
7 #include "AliFJWrapper.h"
8 #include "AliEmcalJet.h"
9 #include "AliEmcalJetFinder.h"
14 //________________________________________________________________________
15 AliEmcalJetFinder::AliEmcalJetFinder() :
16 TNamed("EmcalJetFinder","EmcalJetFinder"), fFastjetWrapper(0), fInputVectorIndex(0), fJetCount(0), fJetArray(), fGhostArea(0.005), fRadius(0.4), fJetAlgorithm(0), fRecombScheme(-1), fTrackMaxEta(0.9), fJetMaxEta(0.5), fJetMinPt(0), fJetMinArea(0)
19 fFastjetWrapper = new AliFJWrapper("FJWrapper", "FJWrapper");
22 //________________________________________________________________________
23 AliEmcalJetFinder::AliEmcalJetFinder(const char* name) :
24 TNamed(name, name), fFastjetWrapper(0), fInputVectorIndex(0), fJetCount(0), fJetArray(), fGhostArea(0.005), fRadius(0.4), fJetAlgorithm(0), fRecombScheme(-1), fTrackMaxEta(0.9), fJetMaxEta(0.5), fJetMinPt(0), fJetMinArea(0)
27 fFastjetWrapper = new AliFJWrapper("FJWrapper", "FJWrapper");
30 //________________________________________________________________________
31 AliEmcalJetFinder::~AliEmcalJetFinder()
34 delete fFastjetWrapper;
38 //________________________________________________________________________
39 Bool_t AliEmcalJetFinder::FindJets()
41 // Tidy up and check input
44 if(!fInputVectorIndex)
46 AliError("No input vectors added to jet finder!");
50 // Pass settings to fastjet
51 fFastjetWrapper->SetAreaType(fastjet::active_area_explicit_ghosts);
52 fFastjetWrapper->SetGhostArea(fGhostArea);
53 fFastjetWrapper->SetR(fRadius);
54 if(fJetAlgorithm == 0)
55 fFastjetWrapper->SetAlgorithm(fastjet::antikt_algorithm);
56 if(fJetAlgorithm == 1)
57 fFastjetWrapper->SetAlgorithm(fastjet::kt_algorithm);
59 fFastjetWrapper->SetRecombScheme(static_cast<fastjet::RecombinationScheme>(fRecombScheme));
61 fFastjetWrapper->SetMaxRap(fTrackMaxEta);
64 fFastjetWrapper->Run();
66 // Save the found jets as light-weight objects
67 std::vector<fastjet::PseudoJet> fastjets = fFastjetWrapper->GetInclusiveJets();
68 fJetArray.resize(fastjets.size());
70 for (UInt_t i=0; i<fastjets.size(); i++)
73 if (fastjets[i].perp()<fJetMinPt)
75 if (fFastjetWrapper->GetJetArea(i)<fJetMinArea)
77 if (TMath::Abs(fastjets[i].eta())>fJetMaxEta)
80 AliEmcalJet* jet = new AliEmcalJet(fastjets[i].perp(), fastjets[i].eta(), fastjets[i].phi(), fastjets[i].m());
82 // Set the most important properties of the jet
83 Int_t nConstituents(fFastjetWrapper->GetJetConstituents(i).size());
84 jet->SetNumberOfTracks(nConstituents);
85 jet->SetNumberOfClusters(nConstituents);
86 fJetArray[fJetCount] = jet;
90 fJetArray.resize(fJetCount);
93 fFastjetWrapper->Clear();
94 fInputVectorIndex = 0;
98 //________________________________________________________________________
99 void AliEmcalJetFinder::AddInputVector(Double_t px, Double_t py, Double_t pz)
101 fFastjetWrapper->AddInputVector(px, py, pz, TMath::Sqrt(px*px + py*py + pz*pz), fInputVectorIndex + 100);
105 //________________________________________________________________________
106 void AliEmcalJetFinder::AddInputVector(Double_t px, Double_t py, Double_t pz, Double_t E)
108 fFastjetWrapper->AddInputVector(px, py, pz, E, fInputVectorIndex + 100);
112 //________________________________________________________________________
113 void AliEmcalJetFinder::FillPtHistogram(TH1* histogram)
117 for (Int_t i=0; i<fJetArray.size(); i++)
119 histogram->Fill(fJetArray[i]->Pt());
123 //________________________________________________________________________
124 void AliEmcalJetFinder::FillPhiHistogram(TH1* histogram)
128 for (Int_t i=0; i<fJetArray.size(); i++)
130 histogram->Fill(fJetArray[i]->Phi());
134 //________________________________________________________________________
135 void AliEmcalJetFinder::FillEtaHistogram(TH1* histogram)
139 for (Int_t i=0; i<fJetArray.size(); i++)
141 histogram->Fill(fJetArray[i]->Eta());