FastJet wrapper and CMakefile
authorfronchet <fronchet@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 8 Aug 2011 12:30:51 +0000 (12:30 +0000)
committerfronchet <fronchet@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 8 Aug 2011 12:30:51 +0000 (12:30 +0000)
HLT/FJWrapper/AliFJWrapper.cxx [new file with mode: 0644]
HLT/FJWrapper/AliFJWrapper.h [new file with mode: 0644]
HLT/FJWrapper/CMakejet.conf [new file with mode: 0644]

diff --git a/HLT/FJWrapper/AliFJWrapper.cxx b/HLT/FJWrapper/AliFJWrapper.cxx
new file mode 100644 (file)
index 0000000..c29f124
--- /dev/null
@@ -0,0 +1,514 @@
+#include "AliFJWrapper.h"
+
+#include "TString.h"
+
+#include "AliLog.h"
+#include "AliFastJetHeader.h"
+
+#include <cmath>
+#include <iostream>
+/*
+@
+@
+@
+@
+@
+@
+@
+ */
+using namespace std;
+namespace fj = fastjet;
+
+AliFJWrapper::AliFJWrapper(const char *name, const char *title)
+  : TNamed(name, title)
+    
+  , fInputVectors      ( )
+  , fInclusiveJets     ( )
+
+  , fSubtractedJetsPt  ( )
+
+  , fAreaDef           (0)
+  , fVorAreaSpec       (0)
+  , fGhostedAreaSpec   (0)
+  , fJetDef            (0)
+  , fPlugin            (0)
+  , fRange             (0)
+  , fClustSeq          (0)
+
+  , fStrategy          (fj::Best)
+  , fAlgor             (fj::kt_algorithm)
+  , fScheme            (fj::BIpt_scheme)
+  , fAreaType          (fj::active_area)
+
+  , fNGhostRepeats     (1)
+  , fGhostArea         (0.01)
+  , fMaxRap            (1.)
+  , fR                 (0.4)
+
+  , fGridScatter       (1.0)
+  , fKtScatter         (1.0)
+  , fMeanGhostKt       (1e-100)
+  , fPluginAlgor       (0)
+
+  , fMedianUsedForBgSubtraction (0)
+  , fUseArea4Vector    (false)
+{
+  ;
+}
+
+void AliFJWrapper::CopySettingsFrom (const AliFJWrapper& wrapper)
+{
+  //
+  // copy some settings - avoid double typing...
+  // you very often want to keep most of the settings 
+  // but change only the algorithm or R - do it after call to this function
+  //
+
+  fStrategy       = wrapper.fStrategy;
+  fAlgor          = wrapper.fAlgor;
+  fScheme         = wrapper.fScheme;
+  fAreaType       = wrapper.fAreaType;
+  
+  fNGhostRepeats  = wrapper.fNGhostRepeats;
+  fGhostArea      = wrapper.fGhostArea;
+  fMaxRap         = wrapper.fMaxRap;
+  fR              = wrapper.fR;
+  
+  fGridScatter    = wrapper.fGridScatter;
+  fKtScatter      = wrapper.fKtScatter;
+  fMeanGhostKt    = wrapper.fMeanGhostKt;
+  fPluginAlgor    = wrapper.fPluginAlgor;
+  
+  fUseArea4Vector = wrapper.fUseArea4Vector;
+  
+}
+
+AliFJWrapper::~AliFJWrapper()
+{  
+  //
+  // destructor
+  //
+  delete fAreaDef;
+  delete fVorAreaSpec;
+  delete fGhostedAreaSpec;
+  delete fJetDef;
+  delete fPlugin;
+  delete fRange;
+  delete fClustSeq;  
+}
+
+void AliFJWrapper::Clear(const Option_t *opt)
+{
+  //
+  // simply clear the input vectors
+  // make sure done on every event if the instance is reused
+  // reset the median to zero
+  //
+
+  fInputVectors.clear();
+  fMedianUsedForBgSubtraction = 0;
+
+  // for the moment brute force delete everything
+  delete fAreaDef;
+  delete fVorAreaSpec;
+  delete fGhostedAreaSpec;
+  delete fJetDef;
+  delete fPlugin;
+  delete fRange;
+  delete fClustSeq;  
+
+  TNamed::Clear(opt);
+}
+
+void AliFJWrapper::SetupFromHeader(const AliFastJetHeader *header)
+{
+  //
+  // feed the wrapper from the parameters of the header
+  //
+  
+  fR             = header->GetRparam();
+  fStrategy      = header->GetStrategy();
+  fScheme        = header->GetRecombScheme();
+  fAlgor         = header->GetAlgorithm(); 
+  fAreaType      = header->GetAreaType();
+  fMaxRap        = header->GetGhostEtaMax(); 
+  fGhostArea     = header->GetGhostArea(); 
+  fNGhostRepeats = header->GetActiveAreaRepeats(); 
+
+}
+
+void AliFJWrapper::AddInputVector(double px, double py, double pz, double E, int index)
+{
+  // make the input pseudojet 
+  fastjet::PseudoJet inVec(px, py, pz, E);
+  
+  if (index > -99999)
+    {
+      inVec.set_user_index(index);
+    }
+  else
+    {
+      inVec.set_user_index(fInputVectors.size());
+    }
+
+  // add to the fj container of input vectors
+  fInputVectors.push_back(inVec);
+}
+
+void AliFJWrapper::AddInputVector(const fj::PseudoJet& vec, int index)
+{
+  // add an input pseudojet 
+  fj::PseudoJet inVec = vec;
+  
+  if (index > -99999)
+    {
+      inVec.set_user_index(index);
+    }
+  else
+    {
+      inVec.set_user_index(fInputVectors.size());
+    }
+
+  // add to the fj container of input vectors
+  fInputVectors.push_back(inVec);
+}
+
+void AliFJWrapper::AddInputVectors(const vector<fastjet::PseudoJet>& vecs, int offsetIndex)
+{
+  //
+  // add the input from vector of pseudojets
+  //
+
+  for (unsigned int i = 0; i < vecs.size(); i++)
+    {
+      fj::PseudoJet inVec = vecs[i];
+      
+      if (offsetIndex > -99999)
+       inVec.set_user_index(fInputVectors.size() + offsetIndex);
+      // add to the fj container of input vectors
+      
+      fInputVectors.push_back(inVec);
+    }
+}
+
+void AliFJWrapper::GetMedianAndSigma(double &median, double &sigma, int remove)
+{
+  //
+  // get the median and sigma from fastjet
+  // user can also do it on his own because the cluster sequence is exposed (via a getter)
+  //
+
+  if (!fClustSeq)
+    {
+      AliError("[e] Run the jfinder first.");
+    }
+
+  double mean_area = 0;
+  try 
+    {
+      if(0 == remove)
+       {
+         fClustSeq->get_median_rho_and_sigma(*fRange, fUseArea4Vector, median, sigma, mean_area);
+       }
+      else
+       {
+         std::vector<fastjet::PseudoJet> input_jets = sorted_by_pt(fClustSeq->inclusive_jets());
+         input_jets.erase(input_jets.begin(), input_jets.begin() + remove);
+         fClustSeq->get_median_rho_and_sigma(input_jets, *fRange, false, median, sigma, mean_area);
+         input_jets.clear();
+       }
+    }
+  
+  catch (fj::Error)
+    {
+      AliError(" [w] FJ Exception caught.");
+      median = -1.;
+      sigma = -1;
+    }
+}
+
+int AliFJWrapper::Run()
+{
+  //
+  // run the actual jet finder
+  //
+
+  if (fAreaType == fj::voronoi_area)
+    {
+      // Rfact - check dependence - default is 1.
+      // NOTE: hardcoded variable!
+      fVorAreaSpec = new fj::VoronoiAreaSpec(1.); 
+      fAreaDef     = new fj::AreaDefinition(*fVorAreaSpec);      
+    }
+  else
+    {
+      fGhostedAreaSpec = new fj::GhostedAreaSpec(fMaxRap,
+                                                 fNGhostRepeats, 
+                                                 fGhostArea,
+                                                 fGridScatter,
+                                                 fKtScatter,
+                                                 fMeanGhostKt);
+
+      fAreaDef = new fj::AreaDefinition(*fGhostedAreaSpec, fAreaType);
+    }
+
+  // this is acceptable by fastjet:
+  fRange = new fj::RangeDefinition(fMaxRap - 0.95 * fR);
+
+  if (fAlgor == fj::plugin_algorithm) 
+    {
+      if (fPluginAlgor == 0)
+       {
+         // SIS CONE ALGOR
+         // NOTE: hardcoded split parameter
+         double overlap_threshold = 0.75; // NOTE: this actually splits a lot: thr/min(pt1,pt2)
+         fPlugin = new fj::SISConePlugin(fR, 
+                                          overlap_threshold,
+                                          0, //search of stable cones - zero = until no more
+                                          1.0); // this should be seed effectively for proto jets
+         //__INFO("Define siscone");
+         fJetDef = new fastjet::JetDefinition(fPlugin);
+         //__INFO("SIS cone initialized");
+       }
+      else
+       {
+         AliError("[e] Unrecognized plugin number!");
+       }
+    }
+  else
+    {
+      //fJetDef = new fj::JetDefinition(fj::kt_algorithm, 0.4, fj::pt_scheme, fj::Best);
+      fJetDef = new fj::JetDefinition(fAlgor, fR, fScheme, fStrategy);
+    }
+  
+  try 
+    {
+      fClustSeq = new fj::ClusterSequenceArea(fInputVectors, *fJetDef, *fAreaDef);
+    }
+  catch (fj::Error)
+    {
+      AliError(" [w] FJ Exception caught.");
+      return -1;
+    }
+
+  // inclusive jets:
+  fInclusiveJets.clear();
+  fInclusiveJets = fClustSeq->inclusive_jets(0.0); //sorted_by_pt(fClustSeq->inclusive_jets(0.0));
+
+  return 0;
+}
+
+void AliFJWrapper::SubtractBackground(const double median_pt)
+{
+  //
+  // subtract the background (specify the value - see below the meaning)
+  // negative argument means the bg will be determined with the current algorithm
+  // this is the default behavior
+  // zero is allowed
+  // Note: user may set the switch for area4vector based subtraction
+  //
+
+  double median = 0;
+  double sigma = 0;
+  double mean_area = 0;
+
+  // clear the subtracted jet pt's vector<double>
+  fSubtractedJetsPt.clear();
+      
+  // check what was specified (default is -1)
+  if (median_pt < 0)
+    {
+      try 
+       {
+         fClustSeq->get_median_rho_and_sigma(*fRange, false, median, sigma, mean_area);
+       }
+      
+      catch (fj::Error)
+       {
+         AliError(" [w] FJ Exception caught.");
+         median = -9999.;
+         sigma = -1;
+         fMedianUsedForBgSubtraction = median;
+         return;
+       }
+      fMedianUsedForBgSubtraction = median;
+    }
+  else
+    {
+      // we do not know the sigma in this case
+      sigma = -1;
+      if (0.0 == median_pt)
+       {
+         // AliWarning(Form(" [w] Median specified for bg subtraction is ZERO: %f", median_pt )) << endl;
+         fMedianUsedForBgSubtraction = 0.;
+       }
+      else
+       {
+         fMedianUsedForBgSubtraction = median_pt;
+       }
+    }
+
+  // subtract:
+  for (unsigned i = 0; i < fInclusiveJets.size(); i++) 
+    {
+      //fj::PseudoJet jet_sub = fClustSeq->subtracted_jet(fInclusiveJets[i], fMedianUsedForBgSubtraction_);
+      //cout << fInclusiveJets[i].perp() << " (" << fMedianUsedForBgSubtraction_ << ") -> " << jet_sub.perp() << endl;
+      
+      if (true == fUseArea4Vector)
+       {
+         // sutract the background using the area4vector
+         fj::PseudoJet area4v = fClustSeq->area_4vector(fInclusiveJets[i]);
+         fj::PseudoJet jet_sub = fInclusiveJets[i] - area4v * fMedianUsedForBgSubtraction;
+         fSubtractedJetsPt.push_back(jet_sub.perp()); // here we put only the pt of the jet - note: this can be negative
+       }
+      else
+       {
+         // subtract the background using scalars
+         // fj::PseudoJet jet_sub = fInclusiveJets[i] - area * fMedianUsedForBgSubtraction_;
+
+         double area = fClustSeq->area(fInclusiveJets[i]);
+
+         // standard subtraction
+         double pt_sub = fInclusiveJets[i].perp() - fMedianUsedForBgSubtraction * area;
+         fSubtractedJetsPt.push_back(pt_sub); // here we put only the pt of the jet - note: this can be negative
+         //cout << "phi incl : " << fInclusiveJets[i].phi() << " subtracted : " << jet_sub.phi() << endl;
+       }
+    }
+}
+
+double AliFJWrapper::GetJetArea(const unsigned int idx)
+{
+  //
+  // get the jet area
+  //
+  double retval = -1; // really wrong area..
+  if ( idx < fInclusiveJets.size() )
+    {
+      retval = fClustSeq->area(fInclusiveJets[idx]);
+    }
+  else
+    {
+      cerr << "[e] ::GetJetArea wrong index: " << idx << endl;
+    }
+  return retval;
+}
+
+std::vector<double> AliFJWrapper::GetSubtractedJetsPts(const double median_pt, const bool sorted)
+{ 
+  //
+  // get subtracted jets' pTs returns vector
+  //
+  SubtractBackground(median_pt);
+  
+  if (true == sorted)
+    {
+      std::sort(fSubtractedJetsPt.begin(), fSubtractedJetsPt.begin());
+    }
+  return fSubtractedJetsPt;
+}
+
+double AliFJWrapper::GetJetSubtractedPt(const unsigned int idx)
+{
+  //
+  // get subtracted jets' pTs returns double (pT)
+  //
+  double retval = -99999.; // really wrong pt..
+  if ( idx < fSubtractedJetsPt.size() )
+    {
+      retval = fSubtractedJetsPt[idx];
+    }
+  return retval;
+}
+
+std::vector<fastjet::PseudoJet> 
+AliFJWrapper::GetJetConstituents(const unsigned int idx)
+{
+  //
+  // get jets' constituents
+  //
+  std::vector<fastjet::PseudoJet> retval;
+  
+  if ( idx < fInclusiveJets.size() )
+    {
+      retval = fClustSeq->constituents(fInclusiveJets[idx]);
+    }  
+  else
+    {
+      cerr << "[e] ::GetJetConstituents wrong index: " << idx << endl;
+    }
+  
+  return retval;
+}
+
+void AliFJWrapper::SetupAlgorithmfromOpt(const char *option)
+{
+  //
+  // setup algorithm from char
+  //
+  string opt(option);
+  
+  if (!opt.compare("kt"))                fAlgor    = fj::kt_algorithm;
+  if (!opt.compare("antikt"))            fAlgor    = fj::antikt_algorithm;
+  if (!opt.compare("cambridge"))         fAlgor    = fj::cambridge_algorithm;
+  if (!opt.compare("genkt"))             fAlgor    = fj::genkt_algorithm;
+  if (!opt.compare("cambridge_passive")) fAlgor    = fj::cambridge_for_passive_algorithm;
+  if (!opt.compare("genkt_passive"))     fAlgor    = fj::genkt_for_passive_algorithm;
+  if (!opt.compare("ee_kt"))             fAlgor    = fj::ee_kt_algorithm;
+  if (!opt.compare("ee_genkt"))          fAlgor    = fj::ee_genkt_algorithm;
+  if (!opt.compare("plugin"))            fAlgor    = fj::plugin_algorithm;
+
+
+}
+void AliFJWrapper::SetupStrategyfromOpt(const char *option)
+{
+  //
+  // setup strategy from char
+  //
+  string opt(option);
+
+  if (!opt.compare("Best"))            fStrategy = fj::Best;
+  if (!opt.compare("N2MinHeapTiled"))  fStrategy = fj::N2MinHeapTiled;
+  if (!opt.compare("N2Tiled"))         fStrategy = fj::N2Tiled;
+  if (!opt.compare("N2PoorTiled"))     fStrategy = fj::N2PoorTiled;
+  if (!opt.compare("N2Plain"))         fStrategy = fj::N2Plain;
+  if (!opt.compare("N3Dumb"))          fStrategy = fj::N3Dumb;
+  if (!opt.compare("NlnN"))            fStrategy = fj::NlnN;
+  if (!opt.compare("NlnN3pi"))         fStrategy = fj::NlnN3pi;
+  if (!opt.compare("NlnN4pi"))         fStrategy = fj::NlnN4pi;
+  if (!opt.compare("NlnNCam4pi"))      fStrategy = fj::NlnNCam4pi;
+  if (!opt.compare("NlnNCam2pi2R"))    fStrategy = fj::NlnNCam2pi2R;
+  if (!opt.compare("NlnNCam"))         fStrategy = fj::NlnNCam;
+  if (!opt.compare("plugin"))          fStrategy = fj::plugin_strategy;
+
+}
+void AliFJWrapper::SetupSchemefromOpt(const char *option)
+{
+  //
+  // setup scheme from char
+  //
+  string opt(option);
+
+  if (!opt.compare("BIpt"))   fScheme   = fj::BIpt_scheme;
+  if (!opt.compare("BIpt2"))  fScheme   = fj::BIpt2_scheme;
+  if (!opt.compare("E"))      fScheme   = fj::E_scheme;
+  if (!opt.compare("pt"))     fScheme   = fj::pt_scheme;
+  if (!opt.compare("pt2"))    fScheme   = fj::pt2_scheme;
+  if (!opt.compare("Et"))     fScheme   = fj::Et_scheme;
+  if (!opt.compare("Et2"))    fScheme   = fj::Et2_scheme;
+
+}
+void AliFJWrapper::SetupAreaTypefromOpt(const char *option)
+{
+  //
+  // setup area type from char
+  //
+  string opt(option);
+
+  if (!opt.compare("active"))                      fAreaType = fj::active_area;
+  if (!opt.compare("invalid"))                     fAreaType = fj::invalid_area;
+  if (!opt.compare("active_area_explicit_ghosts")) fAreaType = fj::active_area_explicit_ghosts;
+  if (!opt.compare("one_ghost_passive"))           fAreaType = fj::one_ghost_passive_area;
+  if (!opt.compare("passive"))                     fAreaType = fj::passive_area;
+  if (!opt.compare("voronoi"))                     fAreaType = fj::voronoi_area;
+
+}
diff --git a/HLT/FJWrapper/AliFJWrapper.h b/HLT/FJWrapper/AliFJWrapper.h
new file mode 100644 (file)
index 0000000..61d8870
--- /dev/null
@@ -0,0 +1,125 @@
+#ifndef AliFJWrapper_HH
+#define AliFJWrapper_HH
+/*
+@Comments:
+@
+@
+@
+@
+@
+@
+@
+*/
+#include <vector>
+#include <fastjet/PseudoJet.hh>
+#include <fastjet/JetDefinition.hh>
+#include <fastjet/ClusterSequence.hh>
+#include <fastjet/ClusterSequenceArea.hh>
+#include <fastjet/AreaDefinition.hh>
+#include <fastjet/SISConePlugin.hh>
+
+#include "TNamed.h"
+
+class AliFastJetHeader;
+
+class AliFJWrapper : public TNamed
+{
+ public:
+
+  AliFJWrapper(const char *name, const char *title);
+
+  virtual ~AliFJWrapper();
+  
+  virtual void CopySettingsFrom (const AliFJWrapper&     wrapper);
+  virtual void SetupFromHeader  (const AliFastJetHeader* header);
+
+  virtual void AddInputVector (double px, double py, double pz, double E, int index = -99999);
+
+  virtual void AddInputVector (const fastjet::PseudoJet& vec,                int index = -99999);
+  virtual void AddInputVectors(const std::vector<fastjet::PseudoJet>& vecs,  int offsetIndex = -99999);
+
+  virtual Int_t  Run();
+  virtual void   Clear(const Option_t* opt = "");
+
+  void SetStrategy     (const fastjet::Strategy strat)             { fStrategy = strat;  }
+  void SetAlgorithm    (const fastjet::JetAlgorithm algor)         { fAlgor    = algor;  }
+  void SetRecombScheme (const fastjet::RecombinationScheme scheme) { fScheme   = scheme; }
+  void SetAreaType     (const fastjet::AreaType atype)             { fAreaType = atype;  }
+
+  void SetupAlgorithmfromOpt (const char *option);
+  void SetupStrategyfromOpt  (const char *option);
+  void SetupSchemefromOpt    (const char *option);
+  void SetupAreaTypefromOpt  (const char *option);
+
+  void SetNRepeats      (const int    nrepeat) { fNGhostRepeats  = nrepeat; }
+  void SetGhostArea     (const double gharea)  { fGhostArea      = gharea;  }
+  void SetMaxRap        (const double maxrap)  { fMaxRap         = maxrap;  }
+  void SetR             (const double r)       { fR              = r;       }
+  void SetGridScatter   (const double gridSc)  { fGridScatter    = gridSc;  }
+  void SetKtScatter     (const double ktSc)    { fKtScatter      = ktSc;    }
+  void SetMeanGhostKt   (const double meankt)  { fMeanGhostKt    = meankt;  }
+  void SetPluginAlgor   (const int    plugin)  { fPluginAlgor    = plugin;  }
+  void SetUseArea4Vector(const bool   useA4v)  { fUseArea4Vector = useA4v; }
+
+  virtual void   GetMedianAndSigma(double& median, double& sigma, int remove = 0);
+  virtual double GetMedianUsedForBgSubtraction() {return fMedianUsedForBgSubtraction;}
+
+  const std::vector<fastjet::PseudoJet> GetInputVectors()    { return fInputVectors;}
+  const std::vector<fastjet::PseudoJet> GetInclusiveJets()   { return fInclusiveJets; }
+
+  std::vector<fastjet::PseudoJet> GetJetConstituents(const unsigned int idx);
+
+  fastjet::ClusterSequenceArea*   GetClusterSequence() { return fClustSeq;      }
+
+  std::vector<double>             GetSubtractedJetsPts(const double median_pt = -1, const bool sorted = false);
+
+  double                          GetJetArea         (const unsigned int idx);
+  double                          GetJetSubtractedPt (const unsigned int idx);
+
+ protected:
+
+  AliFJWrapper();
+  AliFJWrapper(const AliFJWrapper& wrapper);
+  AliFJWrapper& operator = (const AliFJWrapper& wrapper);// {return *this;}
+
+  std::vector<fastjet::PseudoJet> fInputVectors; //!
+  std::vector<fastjet::PseudoJet> fInclusiveJets; //!
+
+  std::vector<double> fSubtractedJetsPt; //!
+  
+  fastjet::AreaDefinition        *fAreaDef;//!
+  fastjet::VoronoiAreaSpec       *fVorAreaSpec;//!
+  fastjet::GhostedAreaSpec       *fGhostedAreaSpec;//!
+  fastjet::JetDefinition         *fJetDef;//!
+  fastjet::JetDefinition::Plugin *fPlugin;//!
+  fastjet::RangeDefinition       *fRange;//!
+  fastjet::ClusterSequenceArea   *fClustSeq;//!
+
+  // same as in the alifjheader
+  fastjet::Strategy               fStrategy; //!
+  fastjet::JetAlgorithm           fAlgor;    //!
+  fastjet::RecombinationScheme    fScheme;   //!
+  fastjet::AreaType               fAreaType; //!
+
+  int                             fNGhostRepeats; //!
+  double                          fGhostArea;    //!
+  double                          fMaxRap;       //!
+  double                          fR;             //!
+
+  // not included in the header (within AliFastJetHeader)
+  // no setters for the moment - used default values in the constructor
+  double                          fGridScatter; //!
+  double                          fKtScatter;  //!
+  double                          fMeanGhostKt;        //!
+  int                             fPluginAlgor; //!
+
+  // extra parameters
+  double                          fMedianUsedForBgSubtraction; //!
+  bool                            fUseArea4Vector; //!
+
+  virtual void   SubtractBackground(const double median_pt = -1);
+
+  ClassDef(AliFJWrapper, 0)
+};
+
+#endif
diff --git a/HLT/FJWrapper/CMakejet.conf b/HLT/FJWrapper/CMakejet.conf
new file mode 100644 (file)
index 0000000..b5d88b1
--- /dev/null
@@ -0,0 +1,35 @@
+#--------------------------------------------------------------------------------#
+# Package File for jet                                                           #
+# Author : Johny Jose (johny.jose@cern.ch)                                       #
+# Variables Defined :                                                            #
+#                                                                                #
+# SRCS - C++ source files                                                        #
+# HDRS - C++ header files                                                        #
+# DHDR - ROOT Dictionary Linkdef header file                                     #
+# CSRCS - C source files                                                         #
+# CHDRS - C header files                                                         #
+# EINCLUDE - Include directories                                                 #
+# EDEFINE - Compiler definitions                                                 #
+# ELIBS - Extra libraries to link                                                #
+# ELIBSDIR - Extra library directories                                           #
+# PACKFFLAGS - Fortran compiler flags for package                                #
+# PACKCXXFLAGS - C++ compiler flags for package                                  #
+# PACKCFLAGS - C compiler flags for package                                      #
+# PACKSOFLAGS - Shared library linking flags                                     #
+# PACKLDFLAGS - Module linker flags                                              #
+# PACKBLIBS - Libraries to link (Executables only)                               #
+# EXPORT - Header files to be exported                                           #
+# CINTHDRS - Dictionary header files                                             #
+# CINTAUTOLINK - Set automatic dictionary generation                             #
+# ARLIBS - Archive Libraries and objects for linking (Executables only)          #
+# SHLIBS - Shared Libraries and objects for linking (Executables only)           #
+#--------------------------------------------------------------------------------#
+
+if(FASTJET)
+  set ( CLASS_HDRS ${CLASS_HDRS}       ${CLASS_HDRS_FJ})
+  set ( HLTCXXFLAGS_FASTJET    "-I${FASTJET}/include")
+  set ( HLTDEF_FASTJET         "-DHAVE_FASTJET -DWITHFASTJET")
+  set ( HLTSOFLAGS_FASTJET     "-L${FASTJET}/lib -lCGAL -lfastjet")
+
+  set ( CINTFLAGS ${CINTFLAGS}   ${HLTCXXFLAGS_FASTJET} ${HLTDEF_FASTJET})
+endif(FASTJET)