2 #ifndef RIVET_FastJets_HH
3 #define RIVET_FastJets_HH
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"
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"
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"
29 return Vector3(pj.px(), pj.py(), pj.pz());
55 CDFJETCLU, CDFMIDPOINT, D0ILCONE,
56 JADE, DURHAM, TRACKJET };
67 double rparameter,
double seed_threshold=1.0)
68 :
JetAlg(fsp), _adef(0) { _init1(alg, rparameter, seed_threshold); }
72 fastjet::RecombinationScheme recom,
double rparameter)
73 :
JetAlg(fsp), _adef(0) { _init2(type, recom, rparameter); }
77 :
JetAlg(fsp), _adef(0) { _init3(plugin); }
80 :
JetAlg(fsp), _adef(0) { _init3(&plugin); }
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); }
118 size_t numJets(
double ptmin = 0.0)
const;
126 Jets _jets(
double ptmin = 0.0)
const;
154 if (_adef == 0)
return (fastjet::ClusterSequenceArea*) 0;
155 return dynamic_cast<fastjet::ClusterSequenceArea*
>(_cseq.get());
159 const fastjet::JetDefinition&
jetDef()
const {
164 const fastjet::AreaDefinition*
areaDef()
const {
169 vector<double>
ySubJet(
const fastjet::PseudoJet& jet)
const;
173 fastjet::PseudoJet
splitJet(fastjet::PseudoJet jet,
double& last_R)
const;
177 fastjet::PseudoJet
filterJet(fastjet::PseudoJet jet,
double& stingy_R,
const double def_R)
const;
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);
206 fastjet::JetDefinition _jdef;
209 fastjet::AreaDefinition* _adef;
212 shared_ptr<fastjet::ClusterSequence> _cseq;
215 shared_ptr<fastjet::JetDefinition::Plugin> _plugin;
218 mutable map<int, vector<double> > _yscales;
222 map<int, Particle> _particles;
const fastjet::ClusterSequence * clusterSeq() const
Return the cluster sequence (FastJet-specific).
Definition: FastJets.hh:147
FastJets(const FinalState &fsp, fastjet::JetDefinition::Plugin &plugin)
Explicitly pass in an externally-constructed plugin (must be heap-allocated, Rivet will delete) ...
Definition: FastJets.hh:79
FastJets(JetAlgName alg, double rparameter, double seed_threshold=1.0)
Same thing as above, but without an FS (for when we want to pass the particles directly to the calc m...
Definition: FastJets.hh:84
FastJets(fastjet::JetDefinition::Plugin &plugin)
Same thing as above, but without an FS (for when we want to pass the particles directly to the calc m...
Definition: FastJets.hh:93
const fastjet::AreaDefinition * areaDef() const
Return the area definition (FastJet-specific). May be null.
Definition: FastJets.hh:164
size_t numJets(double ptmin=0.0) const
Number of jets above the cut.
Definition: FastJets.cc:145
const fastjet::ClusterSequenceArea * clusterSeqArea() const
Return the cluster sequence (FastJet-specific).
Definition: FastJets.hh:152
PseudoJets pseudoJetsByE(double ptmin=0.0) const
Get the pseudo jets, ordered by .
Definition: FastJets.hh:137
void project(const Event &e)
Perform the projection on the Event.
Definition: FastJets.cc:82
void useJetArea(fastjet::AreaDefinition *adef)
Use provided jet area definition.
Definition: FastJets.hh:113
PseudoJets pseudoJets(double ptmin=0.0) const
Get the pseudo jets (unordered).
Definition: FastJets.cc:175
vector< fastjet::PseudoJet > PseudoJets
Typedef for a collection of PseudoJets.
Definition: FastJets.hh:39
std::vector< Jet > Jets
Typedef for a collection of Jet objects.
Definition: Jet.hh:177
FastJets(const FinalState &fsp, JetAlgName alg, double rparameter, double seed_threshold=1.0)
Definition: FastJets.hh:66
vector< double > ySubJet(const fastjet::PseudoJet &jet) const
Get the subjet splitting variables for the given jet.
Definition: FastJets.cc:184
Vector3 momentum3(const fastjet::PseudoJet &pj)
Make a 3-momentum vector from a FastJet pseudo-jet.
Definition: FastJets.hh:28
Project out jets found using the FastJet package jet algorithms.
Definition: FastJets.hh:48
FastJets(const FinalState &fsp, fastjet::JetAlgorithm type, fastjet::RecombinationScheme recom, double rparameter)
Native argument constructor, using FastJet alg/scheme enums.
Definition: FastJets.hh:71
const fastjet::JetDefinition & jetDef() const
Return the jet definition (FastJet-specific).
Definition: FastJets.hh:159
fastjet::PseudoJet splitJet(fastjet::PseudoJet jet, double &last_R) const
Split a jet a la PRL100,242001(2008). Based on code from G.Salam, A.Davison.
Definition: FastJets.cc:200
PseudoJets pseudoJetsByPt(double ptmin=0.0) const
Get the pseudo jets, ordered by .
Definition: FastJets.hh:132
FastJets(fastjet::JetAlgorithm type, fastjet::RecombinationScheme recom, double rparameter)
Same thing as above, but without an FS (for when we want to pass the particles directly to the calc m...
Definition: FastJets.hh:87
PseudoJets pseudoJetsByRapidity(double ptmin=0.0) const
Get the pseudo jets, ordered by rapidity.
Definition: FastJets.hh:142
virtual const Projection * clone() const
Clone on the heap.
Definition: FastJets.hh:98
FastJets(const FinalState &fsp, fastjet::JetDefinition::Plugin *plugin)
Explicitly pass in an externally-constructed plugin (must be heap-allocated, Rivet will delete) ...
Definition: FastJets.hh:76
FastJets(fastjet::JetDefinition::Plugin *plugin)
Same thing as above, but without an FS (for when we want to pass the particles directly to the calc m...
Definition: FastJets.hh:90
int compare(const Projection &p) const
Compare projections.
Definition: FastJets.cc:69
fastjet::PseudoJet filterJet(fastjet::PseudoJet jet, double &stingy_R, const double def_R) const
Filter a jet a la PRL100,242001(2008). Based on code from G.Salam, A.Davison.
Definition: FastJets.cc:239
Project out all final-state particles in an event. Probably the most important projection in Rivet! ...
Definition: FinalState.hh:14
std::vector< Particle > ParticleVector
Typedef for a vector of Particle objects.
Definition: Particle.fhh:14
size_t size() const
Number of jets.
Definition: FastJets.hh:121
void reset()
Reset the projection. Jet def, etc. are unchanged.
Definition: FastJets.cc:138
Base class for all Rivet projections.
Definition: Projection.hh:28
Three-dimensional specialisation of Vector.
Definition: Vector3.hh:26
FourMomentum momentum(const fastjet::PseudoJet &pj)
Make a 4-momentum vector from a FastJet pseudo-jet.
Definition: FastJets.hh:33
JetAlgName
Wrapper enum for selected Fastjet jet algorithms.
Definition: FastJets.hh:52
Specialized version of the FourVector with momentum/energy functionality.
Definition: Vector4.hh:324
Abstract base class for projections which can return a set of Jets.
Definition: JetAlg.hh:65
void calc(const ParticleVector &ps)
Do the calculation locally (no caching).
Definition: FastJets.cc:93