File indexing completed on 2023-03-17 11:18:21
0001 #ifndef RecoJets_JetAlgorithms_QJets_h
0002 #define RecoJets_JetAlgorithms_QJets_h
0003 #include <queue>
0004 #include <vector>
0005 #include <list>
0006 #include <algorithm>
0007 #include "fastjet/JetDefinition.hh"
0008 #include "fastjet/PseudoJet.hh"
0009 #include "fastjet/ClusterSequence.hh"
0010 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0011 #include "FWCore/Utilities/interface/isFinite.h"
0012 #include "FWCore/Utilities/interface/RandomNumberGenerator.h"
0013 #include "FWCore/ServiceRegistry/interface/Service.h"
0014 #include "CLHEP/Random/RandomEngine.h"
0015
0016 struct JetDistance {
0017 double dij;
0018 int j1;
0019 int j2;
0020 };
0021
0022 class JetDistanceCompare {
0023 public:
0024 JetDistanceCompare(){};
0025 bool operator()(const JetDistance& lhs, const JetDistance& rhs) const { return lhs.dij > rhs.dij; };
0026 };
0027
0028 class Qjets {
0029 private:
0030 double omega;
0031 bool _rand_seed_set;
0032 unsigned int _seed;
0033 double _zcut, _dcut, _dcut_fctr, _exp_min, _exp_max, _rigidity, _truncation_fctr;
0034 std::map<int, bool> _merged_jets;
0035 std::priority_queue<JetDistance, std::vector<JetDistance>, JetDistanceCompare> _distances;
0036 CLHEP::HepRandomEngine* _rnEngine;
0037
0038 double d_ij(const fastjet::PseudoJet& v1, const fastjet::PseudoJet& v2) const;
0039 void computeDCut(fastjet::ClusterSequence& cs);
0040
0041 double Rand();
0042 bool Prune(JetDistance& jd, fastjet::ClusterSequence& cs);
0043 bool JetsUnmerged(const JetDistance& jd) const;
0044 bool JetUnmerged(int num) const;
0045 void ComputeNewDistanceMeasures(fastjet::ClusterSequence& cs, unsigned int new_jet);
0046 void ComputeAllDistances(const std::vector<fastjet::PseudoJet>& inp);
0047 double ComputeMinimumDistance();
0048 double ComputeNormalization(double dmin);
0049 JetDistance GetNextDistance();
0050 bool Same(const JetDistance& lhs, const JetDistance& rhs);
0051
0052 public:
0053 Qjets(double zcut,
0054 double dcut_fctr,
0055 double exp_min,
0056 double exp_max,
0057 double rigidity,
0058 double truncation_fctr,
0059 CLHEP::HepRandomEngine* rnEngine)
0060 : _rand_seed_set(false),
0061 _zcut(zcut),
0062 _dcut(-1.),
0063 _dcut_fctr(dcut_fctr),
0064 _exp_min(exp_min),
0065 _exp_max(exp_max),
0066 _rigidity(rigidity),
0067 _truncation_fctr(truncation_fctr),
0068 _rnEngine(rnEngine){};
0069
0070 void Cluster(fastjet::ClusterSequence& cs);
0071 void SetRandSeed(unsigned int seed);
0072 };
0073
0074 class QjetsBaseExtras : public fastjet::ClusterSequence::Extras {
0075 public:
0076 QjetsBaseExtras() : _wij(-1.) {}
0077 ~QjetsBaseExtras() override {}
0078 virtual double weight() const { return _wij; }
0079 friend class Qjets;
0080
0081 protected:
0082 double _wij;
0083 };
0084
0085 #endif