2 #ifndef RIVET_FastJets_HH 3 #define RIVET_FastJets_HH 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" 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" 21 #include "Rivet/Projections/PxConePlugin.hh" 32 enum Algo { KT, CAM, SISCONE, ANTIKT, PXCONE,
34 CDFJETCLU, CDFMIDPOINT, D0ILCONE,
35 JADE, DURHAM, TRACKJET, GENKTEE };
45 const fastjet::JetDefinition& jdef,
48 fastjet::AreaDefinition* adef=
nullptr)
49 :
JetAlg(fsp, usemuons, useinvis), _jdef(jdef), _adef(adef)
58 const fastjet::JetDefinition& jdef,
59 fastjet::AreaDefinition* adef,
62 :
FastJets(fsp, jdef, usemuons, useinvis, adef)
69 fastjet::JetAlgorithm type,
70 fastjet::RecombinationScheme recom,
double rparameter,
73 fastjet::AreaDefinition* adef=
nullptr)
74 :
FastJets(fsp,
fastjet::JetDefinition(type, rparameter, recom), usemuons, useinvis, adef)
81 fastjet::JetAlgorithm type,
82 fastjet::RecombinationScheme recom,
double rparameter,
83 fastjet::AreaDefinition* adef,
86 :
FastJets(fsp, type, recom, rparameter, usemuons, useinvis, adef)
93 fastjet::JetDefinition::Plugin* plugin,
96 fastjet::AreaDefinition* adef=
nullptr)
97 :
FastJets(fsp,
fastjet::JetDefinition(plugin), usemuons, useinvis, adef)
99 _plugin.
reset(plugin);
106 fastjet::JetDefinition::Plugin* plugin,
107 fastjet::AreaDefinition* adef,
110 :
FastJets(fsp, plugin, usemuons, useinvis, adef)
121 Algo alg,
double rparameter,
124 fastjet::AreaDefinition* adef=
nullptr,
125 double seed_threshold=1.0)
126 :
JetAlg(fsp, usemuons, useinvis)
129 _initJdef(alg, rparameter, seed_threshold);
159 static Jet mkJet(
const PseudoJet& pj,
const Particles& fsparticles,
const Particles& tagparticles=Particles());
161 static Jets
mkJets(
const PseudoJets& pjs,
const Particles& fsparticles,
const Particles& tagparticles=Particles());
212 Jet trimJet(
const Jet& input,
const fastjet::Filter& trimmer)
const;
222 const shared_ptr<fastjet::ClusterSequence>
clusterSeq()
const {
229 return areaDef() ? dynamic_pointer_cast<fastjet::ClusterSequenceArea>(_cseq) :
nullptr;
233 const fastjet::JetDefinition&
jetDef()
const {
241 const shared_ptr<fastjet::AreaDefinition>
areaDef()
const {
252 void _initJdef(
Algo alg,
double rparameter,
double seed_threshold);
265 void calc(
const Particles& fsparticles,
const Particles& tagparticles=Particles());
271 fastjet::JetDefinition _jdef;
274 std::shared_ptr<fastjet::AreaDefinition> _adef;
277 std::shared_ptr<fastjet::ClusterSequence> _cseq;
280 std::shared_ptr<fastjet::JetDefinition::Plugin> _plugin;
283 mutable std::map<int, vector<double> > _yscales;
286 Particles _fsparticles, _tagparticles;
Definition: MC_Cent_pPb.hh:10
const shared_ptr< fastjet::ClusterSequenceArea > clusterSeqArea() const
Definition: FastJets.hh:228
PseudoJets pseudojetsByE(double ptmin=0.0) const
Alias.
Definition: FastJets.hh:202
DEFAULT_RIVET_PROJ_CLONE(FastJets)
Clone on the heap.
PseudoJets pseudojetsByRapidity(double ptmin=0.0) const
Alias.
Definition: FastJets.hh:209
const shared_ptr< fastjet::ClusterSequence > clusterSeq() const
Definition: FastJets.hh:222
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:105
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:44
void useJetArea(fastjet::AreaDefinition *adef)
Use provided jet area definition.
Definition: FastJets.hh:174
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:241
Project out jets found using the FastJet package jet algorithms.
Definition: FastJets.hh:27
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:68
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:198
Algo
Definition: FastJets.hh:32
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:188
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:57
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:195
std::vector< PseudoJet > PseudoJets
Definition: RivetFastJet.hh:29
const fastjet::JetDefinition & jetDef() const
Return the jet definition.
Definition: FastJets.hh:233
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
PseudoJets pseudoJetsByRapidity(double ptmin=0.0) const
Get the pseudo jets, ordered by rapidity.
Definition: FastJets.hh:205
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:191
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:120
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:80
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:92