Rivet
1.8.0
|
00001 // -*- C++ -*- 00002 #ifndef RIVET_FastJets_HH 00003 #define RIVET_FastJets_HH 00004 00005 #include "Rivet/Projection.hh" 00006 #include "Rivet/Projections/JetAlg.hh" 00007 #include "Rivet/Projections/FinalState.hh" 00008 #include "Rivet/Particle.hh" 00009 #include "Rivet/Jet.hh" 00010 00011 #include "fastjet/JetDefinition.hh" 00012 #include "fastjet/AreaDefinition.hh" 00013 #include "fastjet/ClusterSequence.hh" 00014 #include "fastjet/ClusterSequenceArea.hh" 00015 #include "fastjet/PseudoJet.hh" 00016 00017 #include "fastjet/SISConePlugin.hh" 00018 #include "fastjet/CMSIterativeConePlugin.hh" 00019 #include "fastjet/D0RunIIConePlugin.hh" 00020 #include "fastjet/TrackJetPlugin.hh" 00021 #include "fastjet/JadePlugin.hh" 00022 //#include "fastjet/PxConePlugin.hh" 00023 00024 namespace Rivet { 00025 00026 00028 inline Vector3 momentum3(const fastjet::PseudoJet& pj) { 00029 return Vector3(pj.px(), pj.py(), pj.pz()); 00030 } 00031 00033 inline FourMomentum momentum(const fastjet::PseudoJet& pj) { 00034 return FourMomentum(pj.E(), pj.px(), pj.py(), pj.pz()); 00035 } 00036 00037 00039 typedef vector<fastjet::PseudoJet> PseudoJets; 00040 00041 00042 00044 00045 00046 00048 class FastJets : public JetAlg { 00049 00050 public: 00052 enum JetAlgName { KT, CAM, SISCONE, ANTIKT, 00053 PXCONE, 00054 ATLASCONE, CMSCONE, 00055 CDFJETCLU, CDFMIDPOINT, D0ILCONE, 00056 JADE, DURHAM, TRACKJET }; 00057 00058 00060 00061 00066 FastJets(const FinalState& fsp, JetAlgName alg, 00067 double rparameter, double seed_threshold=1.0) 00068 : JetAlg(fsp), _adef(0) { _init1(alg, rparameter, seed_threshold); } 00069 00071 FastJets(const FinalState& fsp, fastjet::JetAlgorithm type, 00072 fastjet::RecombinationScheme recom, double rparameter) 00073 : JetAlg(fsp), _adef(0) { _init2(type, recom, rparameter); } 00074 00076 FastJets(const FinalState& fsp, fastjet::JetDefinition::Plugin* plugin) 00077 : JetAlg(fsp), _adef(0) { _init3(plugin); } 00079 FastJets(const FinalState& fsp, fastjet::JetDefinition::Plugin& plugin) 00080 : JetAlg(fsp), _adef(0) { _init3(&plugin); } 00081 00082 00084 FastJets(JetAlgName alg, double rparameter, double seed_threshold=1.0) 00085 : _adef(0) { _init1(alg, rparameter, seed_threshold); } 00087 FastJets(fastjet::JetAlgorithm type, fastjet::RecombinationScheme recom, double rparameter) 00088 : _adef(0) { _init2(type, recom, rparameter); } 00090 FastJets(fastjet::JetDefinition::Plugin* plugin) 00091 : _adef(0) { _init3(plugin); } 00093 FastJets(fastjet::JetDefinition::Plugin& plugin) 00094 : _adef(0) { _init3(&plugin); } 00095 00096 00098 virtual const Projection* clone() const { 00099 return new FastJets(*this); 00100 } 00101 00103 00104 00105 public: 00106 00108 void reset(); 00109 00113 void useJetArea(fastjet::AreaDefinition* adef) { 00114 _adef = adef; 00115 } 00116 00118 size_t numJets(double ptmin = 0.0) const; 00119 00121 size_t size() const { 00122 return numJets(); 00123 } 00124 00126 Jets _jets(double ptmin = 0.0) const; 00127 00129 PseudoJets pseudoJets(double ptmin = 0.0) const; 00130 00132 PseudoJets pseudoJetsByPt(double ptmin = 0.0) const { 00133 return sorted_by_pt(pseudoJets(ptmin)); 00134 } 00135 00137 PseudoJets pseudoJetsByE(double ptmin = 0.0) const { 00138 return sorted_by_E(pseudoJets(ptmin)); 00139 } 00140 00142 PseudoJets pseudoJetsByRapidity(double ptmin = 0.0) const { 00143 return sorted_by_rapidity(pseudoJets(ptmin)); 00144 } 00145 00147 const fastjet::ClusterSequence* clusterSeq() const { 00148 return _cseq.get(); 00149 } 00150 00152 const fastjet::ClusterSequenceArea* clusterSeqArea() const { 00154 if (_adef == 0) return (fastjet::ClusterSequenceArea*) 0; 00155 return dynamic_cast<fastjet::ClusterSequenceArea*>(_cseq.get()); 00156 } 00157 00159 const fastjet::JetDefinition& jetDef() const { 00160 return _jdef; 00161 } 00162 00164 const fastjet::AreaDefinition* areaDef() const { 00165 return _adef; 00166 } 00167 00169 vector<double> ySubJet(const fastjet::PseudoJet& jet) const; 00170 00173 fastjet::PseudoJet splitJet(fastjet::PseudoJet jet, double& last_R) const; 00174 00177 fastjet::PseudoJet filterJet(fastjet::PseudoJet jet, double& stingy_R, const double def_R) const; 00178 00179 private: 00180 00181 Jets _pseudojetsToJets(const PseudoJets& pjets) const; 00182 00184 void _init1(JetAlgName alg, double rparameter, double seed_threshold); 00185 void _init2(fastjet::JetAlgorithm type, 00186 fastjet::RecombinationScheme recom, double rparameter); 00187 void _init3(fastjet::JetDefinition::Plugin* plugin); 00188 00189 protected: 00190 00192 void project(const Event& e); 00193 00195 int compare(const Projection& p) const; 00196 00197 public: 00198 00200 void calc(const ParticleVector& ps); 00201 00202 00203 private: 00204 00206 fastjet::JetDefinition _jdef; 00207 00209 fastjet::AreaDefinition* _adef; 00210 00212 shared_ptr<fastjet::ClusterSequence> _cseq; 00213 00215 shared_ptr<fastjet::JetDefinition::Plugin> _plugin; 00216 00218 mutable map<int, vector<double> > _yscales; 00219 00221 //set<Particle, ParticleBase::byPTAscending> _particles; 00222 map<int, Particle> _particles; 00223 00224 }; 00225 00226 } 00227 00228 #endif