standalone jet finder class (Ruediger)

diff --git a/PWGJE/CMakelibPWGJEEMCALJetTasks.pkg b/PWGJE/CMakelibPWGJEEMCALJetTasks.pkg
index 793fa673d321adf9de38770b1a51de917db7333a..317002a902d81093e312bad70acd34210c56567d 100644 (file)
--- a/PWGJE/CMakelibPWGJEEMCALJetTasks.pkg
+++ b/PWGJE/CMakelibPWGJEEMCALJetTasks.pkg
@@ -50,6 +50,7 @@ set ( SRCS
  EMCALJetTasks/AliJetRandomizerTask.cxx
  EMCALJetTasks/AliJetResponseMaker.cxx
  EMCALJetTasks/AliRhoParameter.cxx
+ EMCALJetTasks/AliEmcalJetFinder.cxx
  EMCALJetTasks/UserTasks/AliAnalysisTaskCLQA.cxx
  EMCALJetTasks/UserTasks/AliAnalysisTaskChargedJetsPA.cxx
  EMCALJetTasks/UserTasks/AliAnalysisTaskEmcalJetHMEC.cxx

diff --git a/PWGJE/EMCALJetTasks/AliEmcalJetFinder.cxx b/PWGJE/EMCALJetTasks/AliEmcalJetFinder.cxx
new file mode 100644 (file)
index 0000000..eea371d
--- /dev/null
+++ b/PWGJE/EMCALJetTasks/AliEmcalJetFinder.cxx
@@ -0,0 +1,137 @@
+// $Id$
+//
+// Standalone jet finder
+//
+// Authors: R.Haake
+
+#ifndef ALIEMCALJETFINDER_CXX
+#define ALIEMCALJETFINDER_CXX
+
+#include "AliFJWrapper.h"
+#include "AliEmcalJet.h"
+#include "AliEmcalJetFinder.h"
+#include "AliLog.h"
+#include <vector>
+#include "TH1.h"
+
+//________________________________________________________________________
+AliEmcalJetFinder::AliEmcalJetFinder() :
+  TNamed("EmcalJetFinder","EmcalJetFinder"), fFastjetWrapper(0), fInputVectorIndex(0), fJetCount(0), fJetArray(), fGhostArea(0.005), fRadius(0.4), fJetAlgorithm(0), fTrackMaxEta(0.9), fJetMaxEta(0.5), fJetMinPt(0), fJetMinArea(0)
+{
+  // Constructor
+  fFastjetWrapper = new AliFJWrapper("FJWrapper", "FJWrapper");
+}
+
+//________________________________________________________________________
+AliEmcalJetFinder::AliEmcalJetFinder(const char* name) :
+  TNamed(name, name), fFastjetWrapper(0), fInputVectorIndex(0), fJetCount(0), fJetArray(), fGhostArea(0.005), fRadius(0.4), fJetAlgorithm(0), fTrackMaxEta(0.9), fJetMaxEta(0.5), fJetMinPt(0), fJetMinArea(0)
+{
+  // Constructor
+  fFastjetWrapper = new AliFJWrapper("FJWrapper", "FJWrapper");
+}
+
+//________________________________________________________________________
+AliEmcalJetFinder::~AliEmcalJetFinder()
+{
+  if(fFastjetWrapper)
+    delete fFastjetWrapper;
+}
+
+
+//________________________________________________________________________
+Bool_t AliEmcalJetFinder::FindJets()
+{
+  // Tidy up and check input
+  fJetArray.clear();
+  fJetCount = 0;
+  if(!fInputVectorIndex)
+  {
+    AliError("No input vectors added to jet finder!");
+    return kFALSE;
+  }
+
+  // Pass settings to fastjet
+  fFastjetWrapper->SetAreaType(fastjet::active_area_explicit_ghosts);
+  fFastjetWrapper->SetGhostArea(fGhostArea);
+  fFastjetWrapper->SetR(fRadius);
+  if(fJetAlgorithm == 0)
+    fFastjetWrapper->SetAlgorithm(fastjet::antikt_algorithm);  
+  if(fJetAlgorithm == 1)
+    fFastjetWrapper->SetAlgorithm(fastjet::kt_algorithm);  
+
+  fFastjetWrapper->SetMaxRap(fTrackMaxEta);
+
+  // Run jet finding
+  fFastjetWrapper->Run();
+
+  // Save the found jets as light-weight objects
+  std::vector<fastjet::PseudoJet> fastjets = fFastjetWrapper->GetInclusiveJets();
+  fJetArray.resize(fastjets.size());
+
+  for (UInt_t i=0; i<fastjets.size(); i++)
+  {
+    // Apply jet cuts
+    if (fastjets[i].perp()<fJetMinPt) 
+      continue;
+    if (fFastjetWrapper->GetJetArea(i)<fJetMinArea)
+      continue;
+    if (TMath::Abs(fastjets[i].eta())>fJetMaxEta)
+      continue;
+
+    AliEmcalJet* jet = new AliEmcalJet(fastjets[i].perp(), fastjets[i].eta(), fastjets[i].phi(), fastjets[i].m());
+
+    // Set the most important properties of the jet
+    Int_t nConstituents(fFastjetWrapper->GetJetConstituents(i).size());
+    jet->SetNumberOfTracks(nConstituents);
+    jet->SetNumberOfClusters(nConstituents);
+    fJetArray[fJetCount] = jet;
+    fJetCount++;
+  }
+
+  fJetArray.resize(fJetCount);
+
+  fFastjetWrapper->Clear();
+  fInputVectorIndex = 0;
+  return kTRUE;
+}
+    
+//________________________________________________________________________
+void AliEmcalJetFinder::AddInputVector(Double_t px, Double_t py, Double_t pz)
+{
+  fFastjetWrapper->AddInputVector(px, py, pz, TMath::Sqrt(px*px + py*py + pz*pz), fInputVectorIndex + 100);
+  fInputVectorIndex++;
+}
+
+//________________________________________________________________________
+void AliEmcalJetFinder::FillPtHistogram(TH1* histogram)
+{
+  if(!histogram)
+    return;
+  for (Int_t i=0; i<fJetArray.size(); i++)
+  {
+    histogram->Fill(fJetArray[i]->Pt());
+  }
+}
+
+//________________________________________________________________________
+void AliEmcalJetFinder::FillPhiHistogram(TH1* histogram)
+{
+  if(!histogram)
+    return;
+  for (Int_t i=0; i<fJetArray.size(); i++)
+  {
+    histogram->Fill(fJetArray[i]->Phi());
+  }
+}
+
+//________________________________________________________________________
+void AliEmcalJetFinder::FillEtaHistogram(TH1* histogram)
+{
+  if(!histogram)
+    return;
+  for (Int_t i=0; i<fJetArray.size(); i++)
+  {
+    histogram->Fill(fJetArray[i]->Eta());
+  }
+}
+#endif

diff --git a/PWGJE/EMCALJetTasks/AliEmcalJetFinder.h b/PWGJE/EMCALJetTasks/AliEmcalJetFinder.h
new file mode 100644 (file)
index 0000000..95dda54
--- /dev/null
+++ b/PWGJE/EMCALJetTasks/AliEmcalJetFinder.h
@@ -0,0 +1,58 @@
+#ifndef ALIEMCALJETFINDER_H
+#define ALIEMCALJETFINDER_H
+
+// $Id$
+
+namespace fastjet {
+  class PseudoJet;
+}
+
+class AliEmcalJet;
+class AliFJWrapper;
+class TNamed;
+class TH1;
+
+class AliEmcalJetFinder : public TNamed
+{
+  public:
+    AliEmcalJetFinder();
+    AliEmcalJetFinder(const char* name);
+    ~AliEmcalJetFinder();
+    
+    Bool_t                        FindJets();
+    void                          AddInputVector(Double_t px, Double_t py, Double_t pz);
+    void                          FillPtHistogram(TH1* histogram);
+    void                          FillEtaHistogram(TH1* histogram);
+    void                          FillPhiHistogram(TH1* histogram);
+
+    Int_t                         GetJets(AliEmcalJet*& jets)     {jets = fJetArray[0]; return fJetCount;}
+    AliEmcalJet*                  GetJet(Int_t index)             {return fJetArray[index];}
+    void                          SetGhostArea(Double_t val)      {fGhostArea = val;}
+    void                          SetRadius(Double_t val)         {fRadius = val;}
+    void                          SetJetAlgorithm(Int_t val)      {fJetAlgorithm = val;}
+    void                          SetTrackMaxEta(Double_t val)    {fTrackMaxEta = val;}
+    void                          SetJetMaxEta(Double_t val)      {fJetMaxEta = val;}
+    void                          SetJetMinPt(Double_t val)       {fJetMinPt = val;}
+    void                          SetJetMinArea(Double_t val)     {fJetMinArea = val;}
+
+    void                          SetManualIndex(Int_t val)       {fInputVectorIndex = val;}
+  private:
+    // General properties
+    AliFJWrapper*                 fFastjetWrapper;
+    Int_t                         fInputVectorIndex;
+    Int_t                         fJetCount;
+    std::vector<AliEmcalJet*>     fJetArray;
+    // Settings for fastjet
+    Double_t                      fGhostArea;
+    Double_t                      fRadius;
+    Int_t                         fJetAlgorithm;
+    Double_t                      fTrackMaxEta;
+    // Jet cuts
+    Double_t                      fJetMaxEta;
+    Double_t                      fJetMinPt;
+    Double_t                      fJetMinArea;
+
+    ClassDef(AliEmcalJetFinder, 1); // Lightweight fastjet implementation outside analysis tasks
+};
+
+#endif

diff --git a/PWGJE/PWGJEEMCALJetTasksLinkDef.h b/PWGJE/PWGJEEMCALJetTasksLinkDef.h
index d6a2f7a33e575d2e038955c3a5db55394478902d..0a92f044954e6c349d25560acd5def7072df07b3 100644 (file)
--- a/PWGJE/PWGJEEMCALJetTasksLinkDef.h
+++ b/PWGJE/PWGJEEMCALJetTasksLinkDef.h
@@ -28,6 +28,7 @@
 #pragma link C++ class AliJetConstituentTagCopier+;
 #pragma link C++ class AliJetResponseMaker+;
 #pragma link C++ class AliRhoParameter+;
+#pragma link C++ class AliEmcalJetFinder+;
 
 // user tasks
 #pragma link C++ class AliAnalysisTaskCLQA+;