Rivet  1.8.0
FastJets.hh
1 // -*- C++ -*-
2 #ifndef RIVET_FastJets_HH
3 #define RIVET_FastJets_HH
4 
5 #include "Rivet/Projection.hh"
6 #include "Rivet/Projections/JetAlg.hh"
7 #include "Rivet/Projections/FinalState.hh"
8 #include "Rivet/Particle.hh"
9 #include "Rivet/Jet.hh"
10 
11 #include "fastjet/JetDefinition.hh"
12 #include "fastjet/AreaDefinition.hh"
13 #include "fastjet/ClusterSequence.hh"
14 #include "fastjet/ClusterSequenceArea.hh"
15 #include "fastjet/PseudoJet.hh"
16 
17 #include "fastjet/SISConePlugin.hh"
18 #include "fastjet/CMSIterativeConePlugin.hh"
19 #include "fastjet/D0RunIIConePlugin.hh"
20 #include "fastjet/TrackJetPlugin.hh"
21 #include "fastjet/JadePlugin.hh"
22 //#include "fastjet/PxConePlugin.hh"
23 
24 namespace Rivet {
25 
26 
28  inline Vector3 momentum3(const fastjet::PseudoJet& pj) {
29  return Vector3(pj.px(), pj.py(), pj.pz());
30  }
31 
33  inline FourMomentum momentum(const fastjet::PseudoJet& pj) {
34  return FourMomentum(pj.E(), pj.px(), pj.py(), pj.pz());
35  }
36 
37 
39  typedef vector<fastjet::PseudoJet> PseudoJets;
40 
41 
42 
44 
45 
46 
48  class FastJets : public JetAlg {
49 
50  public:
52  enum JetAlgName { KT, CAM, SISCONE, ANTIKT,
53  PXCONE,
54  ATLASCONE, CMSCONE,
55  CDFJETCLU, CDFMIDPOINT, D0ILCONE,
56  JADE, DURHAM, TRACKJET };
57 
58 
60 
61 
66  FastJets(const FinalState& fsp, JetAlgName alg,
67  double rparameter, double seed_threshold=1.0)
68  : JetAlg(fsp), _adef(0) { _init1(alg, rparameter, seed_threshold); }
69 
71  FastJets(const FinalState& fsp, fastjet::JetAlgorithm type,
72  fastjet::RecombinationScheme recom, double rparameter)
73  : JetAlg(fsp), _adef(0) { _init2(type, recom, rparameter); }
74 
76  FastJets(const FinalState& fsp, fastjet::JetDefinition::Plugin* plugin)
77  : JetAlg(fsp), _adef(0) { _init3(plugin); }
79  FastJets(const FinalState& fsp, fastjet::JetDefinition::Plugin& plugin)
80  : JetAlg(fsp), _adef(0) { _init3(&plugin); }
81 
82 
84  FastJets(JetAlgName alg, double rparameter, double seed_threshold=1.0)
85  : _adef(0) { _init1(alg, rparameter, seed_threshold); }
87  FastJets(fastjet::JetAlgorithm type, fastjet::RecombinationScheme recom, double rparameter)
88  : _adef(0) { _init2(type, recom, rparameter); }
90  FastJets(fastjet::JetDefinition::Plugin* plugin)
91  : _adef(0) { _init3(plugin); }
93  FastJets(fastjet::JetDefinition::Plugin& plugin)
94  : _adef(0) { _init3(&plugin); }
95 
96 
98  virtual const Projection* clone() const {
99  return new FastJets(*this);
100  }
101 
103 
104 
105  public:
106 
108  void reset();
109 
113  void useJetArea(fastjet::AreaDefinition* adef) {
114  _adef = adef;
115  }
116 
118  size_t numJets(double ptmin = 0.0) const;
119 
121  size_t size() const {
122  return numJets();
123  }
124 
126  Jets _jets(double ptmin = 0.0) const;
127 
129  PseudoJets pseudoJets(double ptmin = 0.0) const;
130 
132  PseudoJets pseudoJetsByPt(double ptmin = 0.0) const {
133  return sorted_by_pt(pseudoJets(ptmin));
134  }
135 
137  PseudoJets pseudoJetsByE(double ptmin = 0.0) const {
138  return sorted_by_E(pseudoJets(ptmin));
139  }
140 
142  PseudoJets pseudoJetsByRapidity(double ptmin = 0.0) const {
143  return sorted_by_rapidity(pseudoJets(ptmin));
144  }
145 
147  const fastjet::ClusterSequence* clusterSeq() const {
148  return _cseq.get();
149  }
150 
152  const fastjet::ClusterSequenceArea* clusterSeqArea() const {
154  if (_adef == 0) return (fastjet::ClusterSequenceArea*) 0;
155  return dynamic_cast<fastjet::ClusterSequenceArea*>(_cseq.get());
156  }
157 
159  const fastjet::JetDefinition& jetDef() const {
160  return _jdef;
161  }
162 
164  const fastjet::AreaDefinition* areaDef() const {
165  return _adef;
166  }
167 
169  vector<double> ySubJet(const fastjet::PseudoJet& jet) const;
170 
173  fastjet::PseudoJet splitJet(fastjet::PseudoJet jet, double& last_R) const;
174 
177  fastjet::PseudoJet filterJet(fastjet::PseudoJet jet, double& stingy_R, const double def_R) const;
178 
179  private:
180 
181  Jets _pseudojetsToJets(const PseudoJets& pjets) const;
182 
184  void _init1(JetAlgName alg, double rparameter, double seed_threshold);
185  void _init2(fastjet::JetAlgorithm type,
186  fastjet::RecombinationScheme recom, double rparameter);
187  void _init3(fastjet::JetDefinition::Plugin* plugin);
188 
189  protected:
190 
192  void project(const Event& e);
193 
195  int compare(const Projection& p) const;
196 
197  public:
198 
200  void calc(const ParticleVector& ps);
201 
202 
203  private:
204 
206  fastjet::JetDefinition _jdef;
207 
209  fastjet::AreaDefinition* _adef;
210 
212  shared_ptr<fastjet::ClusterSequence> _cseq;
213 
215  shared_ptr<fastjet::JetDefinition::Plugin> _plugin;
216 
218  mutable map<int, vector<double> > _yscales;
219 
221  //set<Particle, ParticleBase::byPTAscending> _particles;
222  map<int, Particle> _particles;
223 
224  };
225 
226 }
227 
228 #endif