rivet is hosted by Hepforge, IPPP Durham
Rivet 4.0.0
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/JetFinder.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
24namespace Rivet {
25
26
28 class FastJets : public JetFinder {
29 public:
30
31 using JetFinder::operator =;
32
33
36
40 FastJets(const FinalState& fsp,
41 const fastjet::JetDefinition& jdef,
42 JetMuons usemuons=JetMuons::ALL,
43 JetInvisibles useinvis=JetInvisibles::NONE,
44 fastjet::AreaDefinition* adef=nullptr)
45 : JetFinder(fsp, usemuons, useinvis), _jdef(jdef), _adef(adef)
46 {
47 _initBase();
48 }
49
53 FastJets(const FinalState& fsp,
54 const fastjet::JetDefinition& jdef,
55 fastjet::AreaDefinition* adef,
56 JetMuons usemuons=JetMuons::ALL,
57 JetInvisibles useinvis=JetInvisibles::NONE)
58 : FastJets(fsp, jdef, usemuons, useinvis, adef)
59 { }
60
64 FastJets(const FinalState& fsp,
65 fastjet::JetAlgorithm type,
66 fastjet::RecombinationScheme recom, double rparameter,
67 JetMuons usemuons=JetMuons::ALL,
68 JetInvisibles useinvis=JetInvisibles::NONE,
69 fastjet::AreaDefinition* adef=nullptr)
70 : FastJets(fsp, fastjet::JetDefinition(type, rparameter, recom), usemuons, useinvis, adef)
71 { }
72
76 FastJets(const FinalState& fsp,
77 fastjet::JetAlgorithm type,
78 fastjet::RecombinationScheme recom, double rparameter,
79 fastjet::AreaDefinition* adef,
80 JetMuons usemuons=JetMuons::ALL,
81 JetInvisibles useinvis=JetInvisibles::NONE)
82 : FastJets(fsp, type, recom, rparameter, usemuons, useinvis, adef)
83 { }
84
88 FastJets(const FinalState& fsp,
89 fastjet::JetDefinition::Plugin* plugin,
90 JetMuons usemuons=JetMuons::ALL,
91 JetInvisibles useinvis=JetInvisibles::NONE,
92 fastjet::AreaDefinition* adef=nullptr)
93 : FastJets(fsp, fastjet::JetDefinition(plugin), usemuons, useinvis, adef)
94 {
95 _plugin.reset(plugin);
96 }
97
102 fastjet::JetDefinition::Plugin* plugin,
103 fastjet::AreaDefinition* adef,
104 JetMuons usemuons=JetMuons::ALL,
105 JetInvisibles useinvis=JetInvisibles::NONE)
106 : FastJets(fsp, plugin, usemuons, useinvis, adef)
107 { }
108
117 JetAlg alg, double rparameter,
118 JetMuons usemuons=JetMuons::ALL,
119 JetInvisibles useinvis=JetInvisibles::NONE,
120 fastjet::AreaDefinition* adef=nullptr,
121 double seed_threshold=1.0)
122 : JetFinder(fsp, usemuons, useinvis), _adef(adef)
123 {
124 _initBase();
125 _initJdef(alg, rparameter, seed_threshold);
126 }
127
128
131
133
134
136 using Projection::operator =;
137
138
141
143 static PseudoJets mkClusterInputs(const Particles& fsparticles, const Particles& tagparticles=Particles());
145 static Jet mkJet(const PseudoJet& pj, const Particles& fsparticles, const Particles& tagparticles=Particles());
147 static Jets mkJets(const PseudoJets& pjs, const Particles& fsparticles, const Particles& tagparticles=Particles());
148
150
151
153 void reset();
154
155
158
163 void useJetArea(fastjet::AreaDefinition* adef) {
164 _adef.reset(adef);
165 }
166
169 _adef.reset();
170 }
171
173
174
177
182 void addTrf(fastjet::Transformer* trf) {
183 _trfs.push_back(shared_ptr<fastjet::Transformer>(trf));
184 }
185
190 template<typename TRFS, typename TRF=typename TRFS::value_type>
191 typename std::enable_if<Derefable<TRF>::value, void>::type
192 addTrfs(const TRFS& trfs) {
193 for (auto& trf : trfs) addTrf(trf);
194 }
195
197 void clearTrfs() {
198 _trfs.clear();
199 }
200
202
203
206
208 Jets _jets() const;
209
212 PseudoJets pseudojets(double ptmin=0.0) const;
213
215 PseudoJets pseudojetsByPt(double ptmin=0.0) const {
216 return sorted_by_pt(pseudojets(ptmin));
217 }
218
220 PseudoJets pseudojetsByE(double ptmin=0.0) const {
221 return sorted_by_E(pseudojets(ptmin));
222 }
223
225 PseudoJets pseudojetsByRapidity(double ptmin=0.0) const {
226 return sorted_by_rapidity(pseudojets(ptmin));
227 }
228
230
231
234
237 const shared_ptr<fastjet::ClusterSequence> clusterSeq() const {
238 return _cseq;
239 }
240
243 const shared_ptr<fastjet::ClusterSequenceArea> clusterSeqArea() const {
244 return areaDef() ? dynamic_pointer_cast<fastjet::ClusterSequenceArea>(_cseq) : nullptr;
245 }
246
248 const fastjet::JetDefinition& jetDef() const {
249 return _jdef;
250 }
251
256 const shared_ptr<fastjet::AreaDefinition> areaDef() const {
257 return _adef;
258 }
259
261
262
263 protected:
264
266 void _initBase();
267 void _initJdef(JetAlg alg, double rparameter, double seed_threshold);
268
269 protected:
270
272 void project(const Event& e);
273
275 CmpState compare(const Projection& p) const;
276
277 public:
278
280 void calc(const Particles& fsparticles, const Particles& tagparticles=Particles());
281
282
283 protected:
284
286 fastjet::JetDefinition _jdef;
287
289 std::shared_ptr<fastjet::AreaDefinition> _adef;
290
292 std::shared_ptr<fastjet::ClusterSequence> _cseq;
293
295 std::shared_ptr<fastjet::JetDefinition::Plugin> _plugin;
296
298 std::vector< std::shared_ptr<fastjet::Transformer> > _trfs;
299
301 mutable std::map<int, vector<double> > _yscales;
302
304 Particles _fsparticles, _tagparticles;
305
306 };
307
308}
309
310#endif
Representation of a HepMC event, and enabler of Projection caching.
Definition Event.hh:22
Find jets using jet algorithms via the FastJet package.
Definition FastJets.hh:28
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...
CmpState compare(const Projection &p) const
Compare projections.
PseudoJets pseudojetsByE(double ptmin=0.0) const
Get the pseudo jets, ordered by .
Definition FastJets.hh:220
FastJets(const FinalState &fsp, fastjet::JetDefinition::Plugin *plugin, fastjet::AreaDefinition *adef, JetMuons usemuons=JetMuons::ALL, JetInvisibles useinvis=JetInvisibles::NONE)
Explicitly pass in an externally-constructed plugin, with reordered args for easier specification of ...
Definition FastJets.hh:101
void clearTrfs()
Don't apply any jet transformers.
Definition FastJets.hh:197
FastJets(const FinalState &fsp, const fastjet::JetDefinition &jdef, fastjet::AreaDefinition *adef, JetMuons usemuons=JetMuons::ALL, JetInvisibles useinvis=JetInvisibles::NONE)
Definition FastJets.hh:53
PseudoJets pseudojetsByRapidity(double ptmin=0.0) const
Get the pseudo jets, ordered by rapidity.
Definition FastJets.hh:225
FastJets(const FinalState &fsp, const fastjet::JetDefinition &jdef, JetMuons usemuons=JetMuons::ALL, JetInvisibles useinvis=JetInvisibles::NONE, fastjet::AreaDefinition *adef=nullptr)
Definition FastJets.hh:40
FastJets(const FinalState &fsp, fastjet::JetDefinition::Plugin *plugin, JetMuons usemuons=JetMuons::ALL, JetInvisibles useinvis=JetInvisibles::NONE, fastjet::AreaDefinition *adef=nullptr)
Explicitly pass in an externally-constructed plugin.
Definition FastJets.hh:88
FastJets(const FinalState &fsp, fastjet::JetAlgorithm type, fastjet::RecombinationScheme recom, double rparameter, JetMuons usemuons=JetMuons::ALL, JetInvisibles useinvis=JetInvisibles::NONE, fastjet::AreaDefinition *adef=nullptr)
Definition FastJets.hh:64
std::enable_if< Derefable< TRF >::value, void >::type addTrfs(const TRFS &trfs)
Add a list of grooming transformers.
Definition FastJets.hh:192
void clearJetArea()
Don't calculate a jet area.
Definition FastJets.hh:168
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.
const shared_ptr< fastjet::ClusterSequence > clusterSeq() const
Definition FastJets.hh:237
PseudoJets pseudojets(double ptmin=0.0) const
FastJets(const FinalState &fsp, JetAlg alg, double rparameter, JetMuons usemuons=JetMuons::ALL, JetInvisibles useinvis=JetInvisibles::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:116
void calc(const Particles &fsparticles, const Particles &tagparticles=Particles())
Do the calculation locally (no caching).
const shared_ptr< fastjet::ClusterSequenceArea > clusterSeqArea() const
Definition FastJets.hh:243
void project(const Event &e)
Perform the projection on the Event.
RIVET_DEFAULT_PROJ_CLONE(FastJets)
Clone on the heap.
void useJetArea(fastjet::AreaDefinition *adef)
Use provided jet area definition.
Definition FastJets.hh:163
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...
const fastjet::JetDefinition & jetDef() const
Return the jet definition.
Definition FastJets.hh:248
void addTrf(fastjet::Transformer *trf)
Add a grooming transformer (base class of fastjet::Filter, etc.)
Definition FastJets.hh:182
PseudoJets pseudojetsByPt(double ptmin=0.0) const
Get the pseudo jets, ordered by .
Definition FastJets.hh:215
FastJets(const FinalState &fsp, fastjet::JetAlgorithm type, fastjet::RecombinationScheme recom, double rparameter, fastjet::AreaDefinition *adef, JetMuons usemuons=JetMuons::ALL, JetInvisibles useinvis=JetInvisibles::NONE)
Definition FastJets.hh:76
void reset()
Reset the projection. Jet def, etc. are unchanged.
const shared_ptr< fastjet::AreaDefinition > areaDef() const
Return the area definition.
Definition FastJets.hh:256
Project out all final-state particles in an event. Probably the most important projection in Rivet!
Definition FinalState.hh:12
Abstract base class for projections which can return a set of Jets.
Definition JetFinder.hh:23
Representation of a clustered jet of particles.
Definition Jet.hh:42
Specialised vector of Jet objects.
Definition Jet.hh:21
Specialised vector of Particle objects.
Definition Particle.hh:21
Base class for all Rivet projections.
Definition Projection.hh:29
Definition MC_CENT_PPB_Projections.hh:10
JetInvisibles
Enum for the treatment of invisible particles: whether to include all, some, or none in jet-finding.
Definition JetFinder.hh:18
std::vector< PseudoJet > PseudoJets
Definition RivetFastJet.hh:30
JetMuons
Enum for the treatment of muons: whether to include all, some, or none in jet-finding.
Definition JetFinder.hh:15