Rivet  3.1.5
FastJets.hh
1 // -*- C++ -*-
2 #ifndef RIVET_FastJets_HH
3 #define RIVET_FastJets_HH
4 
5 #include "Rivet/Jet.hh"
6 #include "Rivet/Particle.hh"
7 #include "Rivet/Projection.hh"
8 #include "Rivet/Projections/JetAlg.hh"
9 #include "Rivet/Projections/FinalState.hh"
10 #include "Rivet/Tools/RivetFastJet.hh"
11 
12 #include "fastjet/SISConePlugin.hh"
13 #include "fastjet/ATLASConePlugin.hh"
14 #include "fastjet/CMSIterativeConePlugin.hh"
15 #include "fastjet/CDFJetCluPlugin.hh"
16 #include "fastjet/CDFMidPointPlugin.hh"
17 #include "fastjet/D0RunIIConePlugin.hh"
18 #include "fastjet/TrackJetPlugin.hh"
19 #include "fastjet/JadePlugin.hh"
20 
21 #include "Rivet/Projections/PxConePlugin.hh"
22 #include "Rivet/Tools/TypeTraits.hh"
23 
24 namespace Rivet {
25 
26 
28  class FastJets : public JetAlg {
29  public:
30 
33  enum Algo { KT=0,
34  AKT=1, ANTIKT=1,
35  CA=2, CAM=2, CAMBRIDGE=2,
36  SISCONE, PXCONE,
37  ATLASCONE, CMSCONE,
38  CDFJETCLU, CDFMIDPOINT, D0ILCONE,
39  JADE, DURHAM, TRACKJET, GENKTEE ,
40  KTET, ANTIKTET };
41 
42 
45 
49  FastJets(const FinalState& fsp,
50  const fastjet::JetDefinition& jdef,
51  JetAlg::Muons usemuons=JetAlg::Muons::ALL,
52  JetAlg::Invisibles useinvis=JetAlg::Invisibles::NONE,
53  fastjet::AreaDefinition* adef=nullptr)
54  : JetAlg(fsp, usemuons, useinvis), _jdef(jdef), _adef(adef)
55  {
56  _initBase();
57  }
58 
62  FastJets(const FinalState& fsp,
63  const fastjet::JetDefinition& jdef,
64  fastjet::AreaDefinition* adef,
65  JetAlg::Muons usemuons=JetAlg::Muons::ALL,
66  JetAlg::Invisibles useinvis=JetAlg::Invisibles::NONE)
67  : FastJets(fsp, jdef, usemuons, useinvis, adef)
68  { }
69 
73  FastJets(const FinalState& fsp,
74  fastjet::JetAlgorithm type,
75  fastjet::RecombinationScheme recom, double rparameter,
76  JetAlg::Muons usemuons=JetAlg::Muons::ALL,
77  JetAlg::Invisibles useinvis=JetAlg::Invisibles::NONE,
78  fastjet::AreaDefinition* adef=nullptr)
79  : FastJets(fsp, fastjet::JetDefinition(type, rparameter, recom), usemuons, useinvis, adef)
80  { }
81 
85  FastJets(const FinalState& fsp,
86  fastjet::JetAlgorithm type,
87  fastjet::RecombinationScheme recom, double rparameter,
88  fastjet::AreaDefinition* adef,
89  JetAlg::Muons usemuons=JetAlg::Muons::ALL,
90  JetAlg::Invisibles useinvis=JetAlg::Invisibles::NONE)
91  : FastJets(fsp, type, recom, rparameter, usemuons, useinvis, adef)
92  { }
93 
97  FastJets(const FinalState& fsp,
98  fastjet::JetDefinition::Plugin* plugin,
99  JetAlg::Muons usemuons=JetAlg::Muons::ALL,
100  JetAlg::Invisibles useinvis=JetAlg::Invisibles::NONE,
101  fastjet::AreaDefinition* adef=nullptr)
102  : FastJets(fsp, fastjet::JetDefinition(plugin), usemuons, useinvis, adef)
103  {
104  _plugin.reset(plugin);
105  }
106 
110  FastJets(const FinalState& fsp,
111  fastjet::JetDefinition::Plugin* plugin,
112  fastjet::AreaDefinition* adef,
113  JetAlg::Muons usemuons=JetAlg::Muons::ALL,
114  JetAlg::Invisibles useinvis=JetAlg::Invisibles::NONE)
115  : FastJets(fsp, plugin, usemuons, useinvis, adef)
116  { }
117 
125  FastJets(const FinalState& fsp,
126  Algo alg, double rparameter,
127  JetAlg::Muons usemuons=JetAlg::Muons::ALL,
128  JetAlg::Invisibles useinvis=JetAlg::Invisibles::NONE,
129  fastjet::AreaDefinition* adef=nullptr,
130  double seed_threshold=1.0)
131  : JetAlg(fsp, usemuons, useinvis)
132  {
133  _initBase();
134  _initJdef(alg, rparameter, seed_threshold);
135  }
136 
137 
138  // /// Same thing as above, but without an FS (for when we want to pass the particles directly to the calc method)
139  // /// @todo Does this work properly, without internal HeavyQuarks etc.?
140  // FastJets(Algo alg, double rparameter, double seed_threshold=1.0) { _initJdef(alg, rparameter, seed_threshold); }
141  // /// Same thing as above, but without an FS (for when we want to pass the particles directly to the calc method)
142  // /// @todo Does this work properly, without internal HeavyQuarks etc.?
143  // FastJets(fastjet::JetAlgorithm type, fastjet::RecombinationScheme recom, double rparameter) { _initJdef(type, recom, rparameter); }
144  // /// Same thing as above, but without an FS (for when we want to pass the particles directly to the calc method)
145  // /// @todo Does this work properly, without internal HeavyQuarks etc.?
146  // FastJets(fastjet::JetDefinition::Plugin* plugin) : _jdef(plugin), _plugin(plugin) {
147  // // _plugin.reset(plugin);
148  // // _jdef = fastjet::JetDefinition(plugin);
149  // }
150 
151 
154 
156 
157 
160 
162  static PseudoJets mkClusterInputs(const Particles& fsparticles, const Particles& tagparticles=Particles());
164  static Jet mkJet(const PseudoJet& pj, const Particles& fsparticles, const Particles& tagparticles=Particles());
166  static Jets mkJets(const PseudoJets& pjs, const Particles& fsparticles, const Particles& tagparticles=Particles());
167 
169 
170 
172  void reset();
173 
174 
177 
182  void useJetArea(fastjet::AreaDefinition* adef) {
183  _adef.reset(adef);
184  }
185 
187  void clearJetArea() {
188  _adef.reset();
189  }
190 
192 
193 
196 
201  void addTrf(fastjet::Transformer* trf) {
202  _trfs.push_back(shared_ptr<fastjet::Transformer>(trf));
203  }
204 
209  template<typename TRFS, typename TRF=typename TRFS::value_type>
210  typename std::enable_if<Derefable<TRF>::value, void>::type
211  addTrfs(const TRFS& trfs) {
212  for (auto& trf : trfs) addTrf(trf);
213  }
214 
216  void clearTrfs() {
217  _trfs.clear();
218  }
219 
223  Jet trimJet(const Jet& input, const fastjet::Filter& trimmer) const;
224 
226 
227 
229 
230 
232  Jets _jets() const;
233 
235  PseudoJets pseudoJets(double ptmin=0.0) const;
237  PseudoJets pseudojets(double ptmin=0.0) const { return pseudoJets(ptmin); }
238 
240  PseudoJets pseudoJetsByPt(double ptmin=0.0) const {
241  return sorted_by_pt(pseudoJets(ptmin));
242  }
244  PseudoJets pseudojetsByPt(double ptmin=0.0) const { return pseudoJetsByPt(ptmin); }
245 
247  PseudoJets pseudoJetsByE(double ptmin=0.0) const {
248  return sorted_by_E(pseudoJets(ptmin));
249  }
251  PseudoJets pseudojetsByE(double ptmin=0.0) const { return pseudoJetsByE(ptmin); }
252 
254  PseudoJets pseudoJetsByRapidity(double ptmin=0.0) const {
255  return sorted_by_rapidity(pseudoJets(ptmin));
256  }
258  PseudoJets pseudojetsByRapidity(double ptmin=0.0) const { return pseudoJetsByRapidity(ptmin); }
259 
261 
262 
264 
265 
268  const shared_ptr<fastjet::ClusterSequence> clusterSeq() const {
269  return _cseq;
270  }
271 
274  const shared_ptr<fastjet::ClusterSequenceArea> clusterSeqArea() const {
275  return areaDef() ? dynamic_pointer_cast<fastjet::ClusterSequenceArea>(_cseq) : nullptr;
276  }
277 
279  const fastjet::JetDefinition& jetDef() const {
280  return _jdef;
281  }
282 
287  const shared_ptr<fastjet::AreaDefinition> areaDef() const {
288  return _adef;
289  }
290 
292 
293 
294  private:
295 
297  void _initBase();
298  void _initJdef(Algo alg, double rparameter, double seed_threshold);
299 
300  protected:
301 
303  void project(const Event& e);
304 
306  CmpState compare(const Projection& p) const;
307 
308  public:
309 
311  void calc(const Particles& fsparticles, const Particles& tagparticles=Particles());
312 
313 
314  private:
315 
317  fastjet::JetDefinition _jdef;
318 
320  std::shared_ptr<fastjet::AreaDefinition> _adef;
321 
323  std::shared_ptr<fastjet::ClusterSequence> _cseq;
324 
326  std::shared_ptr<fastjet::JetDefinition::Plugin> _plugin;
327 
329  std::vector< std::shared_ptr<fastjet::Transformer> > _trfs;
330 
332  mutable std::map<int, vector<double> > _yscales;
333 
335  Particles _fsparticles, _tagparticles;
336 
337  };
338 
339 }
340 
341 #endif
Definition: MC_Cent_pPb.hh:10
void clearJetArea()
Don&#39;t calculate a jet area.
Definition: FastJets.hh:187
const shared_ptr< fastjet::ClusterSequenceArea > clusterSeqArea() const
Definition: FastJets.hh:274
PseudoJets pseudojetsByE(double ptmin=0.0) const
Alias.
Definition: FastJets.hh:251
DEFAULT_RIVET_PROJ_CLONE(FastJets)
Clone on the heap.
PseudoJets pseudojetsByRapidity(double ptmin=0.0) const
Alias.
Definition: FastJets.hh:258
const shared_ptr< fastjet::ClusterSequence > clusterSeq() const
Definition: FastJets.hh:268
FastJets(const FinalState &fsp, fastjet::JetDefinition::Plugin *plugin, fastjet::AreaDefinition *adef, JetAlg::Muons usemuons=JetAlg::Muons::ALL, JetAlg::Invisibles useinvis=JetAlg::Invisibles::NONE)
Explicitly pass in an externally-constructed plugin, with reordered args for easier specification of ...
Definition: FastJets.hh:110
void project(const Event &e)
Perform the projection on the Event.
FastJets(const FinalState &fsp, const fastjet::JetDefinition &jdef, JetAlg::Muons usemuons=JetAlg::Muons::ALL, JetAlg::Invisibles useinvis=JetAlg::Invisibles::NONE, fastjet::AreaDefinition *adef=nullptr)
Definition: FastJets.hh:49
void useJetArea(fastjet::AreaDefinition *adef)
Use provided jet area definition.
Definition: FastJets.hh:182
PseudoJets pseudoJets(double ptmin=0.0) const
Get the pseudo jets (unordered).
const shared_ptr< fastjet::AreaDefinition > areaDef() const
Return the area definition.
Definition: FastJets.hh:287
Project out jets found using the FastJet package jet algorithms.
Definition: FastJets.hh:28
static PseudoJets mkClusterInputs(const Particles &fsparticles, const Particles &tagparticles=Particles())
Make PseudoJets for input to a ClusterSequence, with user_index codes for constituent- and tag-partic...
FastJets(const FinalState &fsp, fastjet::JetAlgorithm type, fastjet::RecombinationScheme recom, double rparameter, JetAlg::Muons usemuons=JetAlg::Muons::ALL, JetAlg::Invisibles useinvis=JetAlg::Invisibles::NONE, fastjet::AreaDefinition *adef=nullptr)
Definition: FastJets.hh:73
Representation of a HepMC event, and enabler of Projection caching.
Definition: Event.hh:22
PseudoJets pseudoJetsByE(double ptmin=0.0) const
Get the pseudo jets, ordered by .
Definition: FastJets.hh:247
Algo
Definition: FastJets.hh:33
Abstract base class for projections which can return a set of Jets.
Definition: JetFinder.hh:15
Definition: RivetFastJet.hh:14
PseudoJets pseudojets(double ptmin=0.0) const
Alias.
Definition: FastJets.hh:237
CmpState compare(const Projection &p) const
Compare projections.
Invisibles
Enum for the treatment of invisible particles: whether to include all, some, or none in jet-finding...
Definition: JetFinder.hh:22
static Jets mkJets(const PseudoJets &pjs, const Particles &fsparticles, const Particles &tagparticles=Particles())
Convert a whole list of PseudoJets to a list of Jets, with mkJet-style unpacking. ...
FastJets(const FinalState &fsp, const fastjet::JetDefinition &jdef, fastjet::AreaDefinition *adef, JetAlg::Muons usemuons=JetAlg::Muons::ALL, JetAlg::Invisibles useinvis=JetAlg::Invisibles::NONE)
Definition: FastJets.hh:62
Project out all final-state particles in an event. Probably the most important projection in Rivet! ...
Definition: FinalState.hh:12
static Jet mkJet(const PseudoJet &pj, const Particles &fsparticles, const Particles &tagparticles=Particles())
Make a Rivet Jet from a PseudoJet holding a user_index code for lookup of Rivet fsparticle or tagpart...
PseudoJets pseudojetsByPt(double ptmin=0.0) const
Alias.
Definition: FastJets.hh:244
std::vector< PseudoJet > PseudoJets
Definition: RivetFastJet.hh:29
const fastjet::JetDefinition & jetDef() const
Return the jet definition.
Definition: FastJets.hh:279
void reset()
Reset the projection. Jet def, etc. are unchanged.
void calc(const Particles &fsparticles, const Particles &tagparticles=Particles())
Do the calculation locally (no caching).
Representation of a clustered jet of particles.
Definition: Jet.hh:18
void clearTrfs()
Don&#39;t apply any jet transformers.
Definition: FastJets.hh:216
PseudoJets pseudoJetsByRapidity(double ptmin=0.0) const
Get the pseudo jets, ordered by rapidity.
Definition: FastJets.hh:254
std::enable_if< Derefable< TRF >::value, void >::type addTrfs(const TRFS &trfs)
Add a list of grooming transformers.
Definition: FastJets.hh:211
Base class for all Rivet projections.
Definition: Projection.hh:29
PseudoJets pseudoJetsByPt(double ptmin=0.0) const
Get the pseudo jets, ordered by .
Definition: FastJets.hh:240
Jet trimJet(const Jet &input, const fastjet::Filter &trimmer) const
Trim (filter) a jet, keeping tag and constituent info in the resulting jet.
FastJets(const FinalState &fsp, Algo alg, double rparameter, JetAlg::Muons usemuons=JetAlg::Muons::ALL, JetAlg::Invisibles useinvis=JetAlg::Invisibles::NONE, fastjet::AreaDefinition *adef=nullptr, double seed_threshold=1.0)
Convenience constructor using Rivet enums for most common jet algs (including some plugins)...
Definition: FastJets.hh:125
FastJets(const FinalState &fsp, fastjet::JetAlgorithm type, fastjet::RecombinationScheme recom, double rparameter, fastjet::AreaDefinition *adef, JetAlg::Muons usemuons=JetAlg::Muons::ALL, JetAlg::Invisibles useinvis=JetAlg::Invisibles::NONE)
Definition: FastJets.hh:85
void addTrf(fastjet::Transformer *trf)
Add a grooming transformer (base class of fastjet::Filter, etc.)
Definition: FastJets.hh:201
Muons
Enum for the treatment of muons: whether to include all, some, or none in jet-finding.
Definition: JetFinder.hh:19
FastJets(const FinalState &fsp, fastjet::JetDefinition::Plugin *plugin, JetAlg::Muons usemuons=JetAlg::Muons::ALL, JetAlg::Invisibles useinvis=JetAlg::Invisibles::NONE, fastjet::AreaDefinition *adef=nullptr)
Explicitly pass in an externally-constructed plugin.
Definition: FastJets.hh:97