rivet is hosted by Hepforge, IPPP Durham
Rivet 4.1.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#include "fastjet/contrib/VariableRPlugin.hh"
21
22#include "Rivet/Projections/PxConePlugin.hh"
23#include "Rivet/Tools/TypeTraits.hh"
24
25namespace Rivet {
26
27 typedef std::shared_ptr<fastjet::JetDefinition::Plugin> FJPluginPtr;
28
30 class FastJets : public JetFinder {
31 public:
32
33 using JetFinder::operator =;
34
35
38
42 FastJets(const FinalState& fsp,
43 const fastjet::JetDefinition& jdef,
44 JetMuons usemuons=JetMuons::ALL,
45 JetInvisibles useinvis=JetInvisibles::NONE,
46 fastjet::AreaDefinition* adef=nullptr)
47 : JetFinder(fsp, usemuons, useinvis), _jdef(jdef),
48 _adef(adef), _cuts(Cuts::OPEN)
49 {
50 _initBase();
51 }
52
56 FastJets(const FinalState& fsp,
57 const fastjet::JetDefinition& jdef,
58 const Cut& c,
59 JetMuons usemuons=JetMuons::ALL,
60 JetInvisibles useinvis=JetInvisibles::NONE,
61 fastjet::AreaDefinition* adef=nullptr)
62 : JetFinder(fsp, usemuons, useinvis), _jdef(jdef),
63 _adef(adef), _cuts(c)
64 {
65 _initBase();
66 }
67
71 FastJets(const FinalState& fsp,
72 const fastjet::JetDefinition& jdef,
73 fastjet::AreaDefinition* adef,
74 JetMuons usemuons=JetMuons::ALL,
75 JetInvisibles useinvis=JetInvisibles::NONE)
76 : FastJets(fsp, jdef, usemuons, useinvis, adef)
77 { }
78
83 FastJets(const FinalState& fsp,
84 const fastjet::JetDefinition& jdef,
85 fastjet::AreaDefinition* adef,
86 const Cut& c,
87 JetMuons usemuons=JetMuons::ALL,
88 JetInvisibles useinvis=JetInvisibles::NONE)
89 : FastJets(fsp, jdef, c, usemuons, useinvis, adef)
90 { }
91
92
96 FastJets(const FinalState& fsp,
97 fastjet::JetAlgorithm type,
98 fastjet::RecombinationScheme recom, double rparameter,
99 JetMuons usemuons=JetMuons::ALL,
100 JetInvisibles useinvis=JetInvisibles::NONE,
101 fastjet::AreaDefinition* adef=nullptr)
102 : FastJets(fsp, fastjet::JetDefinition(type, rparameter, recom), usemuons, useinvis, adef)
103 { }
104
109 fastjet::JetAlgorithm type,
110 fastjet::RecombinationScheme recom, double rparameter,
111 const Cut& c,
112 JetMuons usemuons=JetMuons::ALL,
113 JetInvisibles useinvis=JetInvisibles::NONE,
114 fastjet::AreaDefinition* adef=nullptr)
115 : FastJets(fsp, fastjet::JetDefinition(type, rparameter, recom), c, usemuons, useinvis, adef)
116 { }
117
122 fastjet::JetAlgorithm type,
123 fastjet::RecombinationScheme recom, double rparameter,
124 fastjet::AreaDefinition* adef,
125 JetMuons usemuons=JetMuons::ALL,
126 JetInvisibles useinvis=JetInvisibles::NONE)
127 : FastJets(fsp, type, recom, rparameter, usemuons, useinvis, adef)
128 { }
129
134 fastjet::JetAlgorithm type,
135 fastjet::RecombinationScheme recom, double rparameter,
136 fastjet::AreaDefinition* adef,
137 const Cut& c,
138 JetMuons usemuons=JetMuons::ALL,
139 JetInvisibles useinvis=JetInvisibles::NONE)
140 : FastJets(fsp, type, recom, rparameter, c, usemuons, useinvis, adef)
141 { }
142
147 FJPluginPtr plugin,
148 JetMuons usemuons=JetMuons::ALL,
149 JetInvisibles useinvis=JetInvisibles::NONE,
150 fastjet::AreaDefinition* adef=nullptr)
151 : FastJets(fsp, fastjet::JetDefinition(plugin.get()), usemuons, useinvis, adef)
152 {
153 _plugin = plugin;
154 }
155
160 FJPluginPtr plugin,
161 const Cut& c,
162 JetMuons usemuons=JetMuons::ALL,
163 JetInvisibles useinvis=JetInvisibles::NONE,
164 fastjet::AreaDefinition* adef=nullptr)
165 : FastJets(fsp, fastjet::JetDefinition(plugin.get()), c, usemuons, useinvis, adef)
166 {
167 _plugin = plugin;
168 }
169
174 FJPluginPtr plugin,
175 fastjet::AreaDefinition* adef,
176 JetMuons usemuons=JetMuons::ALL,
177 JetInvisibles useinvis=JetInvisibles::NONE)
178 : FastJets(fsp, plugin, usemuons, useinvis, adef)
179 { }
180
185 FJPluginPtr plugin,
186 fastjet::AreaDefinition* adef,
187 const Cut& c,
188 JetMuons usemuons=JetMuons::ALL,
189 JetInvisibles useinvis=JetInvisibles::NONE)
190 : FastJets(fsp, plugin, c, usemuons, useinvis, adef)
191 { }
192
201 JetAlg alg, double rparameter=-1,
202 JetMuons usemuons=JetMuons::ALL,
203 JetInvisibles useinvis=JetInvisibles::NONE,
204 fastjet::AreaDefinition* adef=nullptr)
205 : JetFinder(fsp, usemuons, useinvis),
206 _adef(adef), _cuts(Cuts::OPEN)
207 {
208 _initBase();
209 _initJdef(alg, rparameter);
210 }
211
220 JetAlg alg, double rparameter,
221 const Cut& c,
222 JetMuons usemuons=JetMuons::ALL,
223 JetInvisibles useinvis=JetInvisibles::NONE,
224 fastjet::AreaDefinition* adef=nullptr)
225 : JetFinder(fsp, usemuons, useinvis), _adef(adef), _cuts(c)
226 {
227 _initBase();
228 _initJdef(alg, rparameter);
229 }
230
231
232
235
237
238
240 using Projection::operator =;
241
242
245
247 static fastjet::JetDefinition mkJetDef(JetAlg alg, double rparameter);
248
251 template<JetAlg JETALG, typename... Args>
252 static FJPluginPtr mkPlugin(Args&&... args){
253 return std::make_shared<mapJetAlg2Plugin_t<JETALG>>(std::forward<Args>(args)...);
254 }
255
257 static FJPluginPtr mkPlugin(JetAlg alg, double rparameter=-1);
258
260 static PseudoJets mkClusterInputs(const Particles& fsparticles, const Particles& tagparticles=Particles());
262 static Jet mkJet(const PseudoJet& pj, const Particles& fsparticles, const Particles& tagparticles=Particles());
264 static Jets mkJets(const PseudoJets& pjs, const Particles& fsparticles=Particles(), const Particles& tagparticles=Particles());
265
267
268
271
272 typedef std::pair<PseudoJets, Particles> PJetsParts;
273
284 template <
285 typename CONTAINER,
286 typename = std::enable_if_t<
287 is_citerable_v<CONTAINER>,
288 Jet
289 >
290 >
291 static CONTAINER reclusterJets(const CONTAINER &jetsIn, const fastjet::JetDefinition &jDef){
292 PJetsParts reclusteredConsts = reclusterJetsParts(jetsIn, jDef);
293 return mkTaggedJets(jetsIn, reclusteredConsts);
294 }
295
304 template <
305 typename CONTAINER,
306 typename = std::enable_if_t<
307 is_citerable_v<CONTAINER>,
308 Jet
309 >
310 >
311 static CONTAINER reclusterJets(const CONTAINER &jetsIn, const fastjet::JetDefinition &jDef, const fastjet::Filter &filter){
312 PJetsParts reclusteredConsts = reclusterJetsParts(jetsIn, jDef);
313 ifilterPseudoJets(reclusteredConsts.first, filter);
314 return mkTaggedJets(jetsIn, reclusteredConsts);
315 }
316
320 template <
321 typename CONTAINER,
322 typename = std::enable_if_t<
323 is_citerable_v<CONTAINER>,
324 Jet
325 >
326 >
327 static CONTAINER reclusterJets(const CONTAINER &jetsIn, const FJPluginPtr &jetAlg){
328 const fastjet::JetDefinition jDef(jetAlg.get());
329 return reclusterJets(jetsIn, jDef);
330 }
331
338 template <
339 typename CONTAINER,
340 typename = std::enable_if_t<
341 is_citerable_v<CONTAINER>,
342 Jet
343 >
344 >
345 static CONTAINER reclusterJets(const CONTAINER &jetsIn, const FJPluginPtr &jetAlg, const fastjet::Filter &filter){
346 const fastjet::JetDefinition jDef(jetAlg.get());
347 return reclusterJets(jetsIn, jDef, filter);
348 }
349
353 template <
354 JetAlg JETALG, typename... Args, typename CONTAINER,
355 typename = std::enable_if_t<
356 is_citerable_v<CONTAINER>,
357 Jet
358 >
359 >
360 static CONTAINER reclusterJets(const CONTAINER &jetsIn, Args&&... args){
361 if constexpr (JETALG<JetAlg::SISCONE){ // JetDefinition
362 const fastjet::JetDefinition jDef = mkJetDef(JETALG, std::forward<Args>(args)...);
363 return reclusterJets(jetsIn, jDef);
364 } else { // Plugin
365 const FJPluginPtr plugin = mkPlugin<JETALG>(std::forward<Args>(args)...);
366 return reclusterJets(jetsIn, plugin);
367 }
368 throw std::invalid_argument( "Unknown jet algorithm: "+to_string(int(JETALG)) );
369 }
370
377 template <
378 JetAlg JETALG, typename... Args, typename CONTAINER,
379 typename = std::enable_if_t<
380 is_citerable_v<CONTAINER>,
381 Jet
382 >
383 >
384 static CONTAINER reclusterJets(const CONTAINER &jetsIn, const fastjet::Filter &filter, Args&&... args){
385 if constexpr (JETALG<JetAlg::SISCONE){ // JetDefinition
386 const fastjet::JetDefinition jDef = mkJetDef(JETALG, std::forward<Args>(args)...);
387 return reclusterJets(jetsIn, jDef, filter);
388 } else { // Plugin
389 const FJPluginPtr plugin = mkPlugin<JETALG>(std::forward<Args>(args)...);
390 return reclusterJets(jetsIn, plugin, filter);
391 }
392 throw std::invalid_argument( "Unknown jet algorithm: "+to_string(int(JETALG)) );
393 }
394
398 template <typename T, typename U, typename... Args>
399 static std::map<T, U> reclusterJets(const std::map<T, U> &jetsMap, Args&&... args){
400 std::map<T, U> rtn;
401 for ( auto const &[key, jetsIn] : jetsMap ) rtn[key] = reclusterJets(jetsIn, std::forward<Args>(args)...);
402 return rtn;
403 }
404
408 template <JetAlg JETALG, typename T, typename U, typename... Args>
409 static std::map<T, U> reclusterJets(const std::map<T, U> &jetsMap, Args&&... args){
410 std::map<T, U> rtn;
411 for ( auto const &[key, jetsIn] : jetsMap ) rtn[key] = reclusterJets<JETALG>(jetsIn, std::forward<Args>(args)...);
412 return rtn;
413 }
414
416
417
419 void reset();
420
421
424
429 void useJetArea(fastjet::AreaDefinition* adef) {
430 _adef.reset(adef);
431 }
432
435 _adef.reset();
436 }
437
439
440
443
448 void addTrf(fastjet::Transformer* trf) {
449 _trfs.push_back(shared_ptr<fastjet::Transformer>(trf));
450 }
451
456 template<typename TRFS, typename TRF=typename TRFS::value_type>
457 typename std::enable_if<Derefable<TRF>::value, void>::type
458 addTrfs(const TRFS& trfs) {
459 for (auto& trf : trfs) addTrf(trf);
460 }
461
465 void addFilter(fastjet::Filter* filter) {
466 addTrf(filter);
467 }
468
472 template<typename FILTERS, typename FILTER=typename FILTERS::value_type>
473 typename std::enable_if<Derefable<FILTER>::value, void>::type
474 addFilters(const FILTERS& filters) {
475 addTrfs(filters);
476 }
477
479 void clearTrfs() {
480 _trfs.clear();
481 }
482
484
485
488
490 Jets _jets() const;
491
494 PseudoJets pseudojets(double ptmin=0.0) const;
495
497 PseudoJets pseudojetsByPt(double ptmin=0.0) const {
498 return sorted_by_pt(pseudojets(ptmin));
499 }
500
502 PseudoJets pseudojetsByE(double ptmin=0.0) const {
503 return sorted_by_E(pseudojets(ptmin));
504 }
505
507 PseudoJets pseudojetsByRapidity(double ptmin=0.0) const {
508 return sorted_by_rapidity(pseudojets(ptmin));
509 }
510
512
513
516
519 const shared_ptr<fastjet::ClusterSequence> clusterSeq() const {
520 return _cseq;
521 }
522
525 const shared_ptr<fastjet::ClusterSequenceArea> clusterSeqArea() const {
526 return areaDef() ? dynamic_pointer_cast<fastjet::ClusterSequenceArea>(_cseq) : nullptr;
527 }
528
530 const fastjet::JetDefinition& jetDef() const {
531 return _jdef;
532 }
533
538 const shared_ptr<fastjet::AreaDefinition> areaDef() const {
539 return _adef;
540 }
541
543
544
545 protected:
546
548 void _initBase();
549
550
551 void _initJdef(JetAlg alg, double rparameter){
552 MSG_DEBUG("JetAlg = " << static_cast<int>(alg));
553 MSG_DEBUG("R parameter = " << rparameter);
554 if ( alg < JetAlg::SISCONE ){ // fastjet::JetDefinitions
555 _jdef = mkJetDef(alg, rparameter);
556 } else { // fastjet::JetDefinition::Plugins
557 _plugin = mkPlugin(alg, rparameter);
558 _jdef = fastjet::JetDefinition(_plugin.get());
559 }
560 }
561
563 static Jets mkTaggedJets(const Jets &jetsIn, const PJetsParts &pJetsParts);
564
566 static PJetsParts reclusterJetsParts(const Jets &jetsIn, const fastjet::JetDefinition &jDef);
567
568 protected:
569
571 void project(const Event& e);
572
574 CmpState compare(const Projection& p) const;
575
576 public:
577
579 void calc(const Particles& fsparticles, const Particles& tagparticles=Particles());
580
581
582 protected:
583
585 fastjet::JetDefinition _jdef;
586
588 std::shared_ptr<fastjet::AreaDefinition> _adef;
589
591 std::shared_ptr<fastjet::ClusterSequence> _cseq;
592
594 Cut _cuts;
595
597 FJPluginPtr _plugin;
598
600 std::vector< std::shared_ptr<fastjet::Transformer> > _trfs;
601
603 static const std::map<JetAlg, std::pair<fastjet::JetAlgorithm, fastjet::RecombinationScheme>> jetAlgMap;
604
607 template<JetAlg, typename X> struct mapJetAlg2Plugin;
608
610 template<typename X> struct mapJetAlg2Plugin<JetAlg::SISCONE, X>{ using type = fastjet::SISConePlugin; };
611 template<typename X> struct mapJetAlg2Plugin<JetAlg::PXCONE, X>{ using type = Rivet::PxConePlugin; };
612 template<typename X> struct mapJetAlg2Plugin<JetAlg::CDFJETCLU, X>{ using type = fastjet::CDFJetCluPlugin; };
613 template<typename X> struct mapJetAlg2Plugin<JetAlg::CDFMIDPOINT, X>{ using type = fastjet::CDFMidPointPlugin; };
614 template<typename X> struct mapJetAlg2Plugin<JetAlg::D0ILCONE, X>{ using type = fastjet::D0RunIIConePlugin; };
615 template<typename X> struct mapJetAlg2Plugin<JetAlg::JADE, X>{ using type = fastjet::JadePlugin; };
616 template<typename X> struct mapJetAlg2Plugin<JetAlg::TRACKJET, X>{ using type = fastjet::TrackJetPlugin; };
617 template<typename X> struct mapJetAlg2Plugin<JetAlg::VARIABLER, X>{ using type = fastjet::contrib::VariableRPlugin; };
618
620 template<JetAlg JETALG>
622
624 mutable std::map<int, vector<double> > _yscales;
625
627 Particles _fsparticles, _tagparticles;
628
629 };
630
631}
632
633#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:30
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.
FastJets(const FinalState &fsp, FJPluginPtr 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:173
static PJetsParts reclusterJetsParts(const Jets &jetsIn, const fastjet::JetDefinition &jDef)
Aux function for reclustering.
PseudoJets pseudojetsByE(double ptmin=0.0) const
Get the pseudo jets, ordered by .
Definition FastJets.hh:502
FastJets(const FinalState &fsp, JetAlg alg, double rparameter=-1, JetMuons usemuons=JetMuons::ALL, JetInvisibles useinvis=JetInvisibles::NONE, fastjet::AreaDefinition *adef=nullptr)
Convenience constructor using Rivet enums for most common jet algs (including some plugins).
Definition FastJets.hh:200
void clearTrfs()
Don't apply any jet transformers.
Definition FastJets.hh:479
FastJets(const FinalState &fsp, const fastjet::JetDefinition &jdef, fastjet::AreaDefinition *adef, JetMuons usemuons=JetMuons::ALL, JetInvisibles useinvis=JetInvisibles::NONE)
Definition FastJets.hh:71
std::enable_if< Derefable< FILTER >::value, void >::type addFilters(const FILTERS &filters)
Add a list of grooming filters.
Definition FastJets.hh:474
PseudoJets pseudojetsByRapidity(double ptmin=0.0) const
Get the pseudo jets, ordered by rapidity.
Definition FastJets.hh:507
FastJets(const FinalState &fsp, const fastjet::JetDefinition &jdef, JetMuons usemuons=JetMuons::ALL, JetInvisibles useinvis=JetInvisibles::NONE, fastjet::AreaDefinition *adef=nullptr)
Definition FastJets.hh:42
FastJets(const FinalState &fsp, const fastjet::JetDefinition &jdef, const Cut &c, JetMuons usemuons=JetMuons::ALL, JetInvisibles useinvis=JetInvisibles::NONE, fastjet::AreaDefinition *adef=nullptr)
Definition FastJets.hh:56
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:96
FastJets(const FinalState &fsp, fastjet::JetAlgorithm type, fastjet::RecombinationScheme recom, double rparameter, const Cut &c, JetMuons usemuons=JetMuons::ALL, JetInvisibles useinvis=JetInvisibles::NONE, fastjet::AreaDefinition *adef=nullptr)
Definition FastJets.hh:108
std::enable_if< Derefable< TRF >::value, void >::type addTrfs(const TRFS &trfs)
Add a list of grooming transformers.
Definition FastJets.hh:458
FastJets(const FinalState &fsp, const fastjet::JetDefinition &jdef, fastjet::AreaDefinition *adef, const Cut &c, JetMuons usemuons=JetMuons::ALL, JetInvisibles useinvis=JetInvisibles::NONE)
Definition FastJets.hh:83
void clearJetArea()
Don't calculate a jet area.
Definition FastJets.hh:434
const shared_ptr< fastjet::ClusterSequence > clusterSeq() const
Definition FastJets.hh:519
void addFilter(fastjet::Filter *filter)
Add a grooming filter.
Definition FastJets.hh:465
static FJPluginPtr mkPlugin(JetAlg alg, double rparameter=-1)
Non-templated version only allowing to set rparameter.
FastJets(const FinalState &fsp, FJPluginPtr plugin, const Cut &c, JetMuons usemuons=JetMuons::ALL, JetInvisibles useinvis=JetInvisibles::NONE, fastjet::AreaDefinition *adef=nullptr)
Explicitly pass in an externally-constructed plugin.
Definition FastJets.hh:159
FastJets(const FinalState &fsp, FJPluginPtr plugin, fastjet::AreaDefinition *adef, const Cut &c, 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:184
PseudoJets pseudojets(double ptmin=0.0) const
FastJets(const FinalState &fsp, JetAlg alg, double rparameter, const Cut &c, JetMuons usemuons=JetMuons::ALL, JetInvisibles useinvis=JetInvisibles::NONE, fastjet::AreaDefinition *adef=nullptr)
Convenience constructor using Cut argument and Rivet enums for most common jet algs (including some p...
Definition FastJets.hh:219
FastJets(const FinalState &fsp, FJPluginPtr plugin, JetMuons usemuons=JetMuons::ALL, JetInvisibles useinvis=JetInvisibles::NONE, fastjet::AreaDefinition *adef=nullptr)
Explicitly pass in an externally-constructed plugin.
Definition FastJets.hh:146
static fastjet::JetDefinition mkJetDef(JetAlg alg, double rparameter)
Make a FastJet JetDefinition according to enum JetAlg.
static Jets mkTaggedJets(const Jets &jetsIn, const PJetsParts &pJetsParts)
Aux function for reclustering.
void calc(const Particles &fsparticles, const Particles &tagparticles=Particles())
Do the calculation locally (no caching).
typename mapJetAlg2Plugin< JETALG, void >::type mapJetAlg2Plugin_t
Make usage of "map" a bit more convenient.
Definition FastJets.hh:621
static const std::map< JetAlg, std::pair< fastjet::JetAlgorithm, fastjet::RecombinationScheme > > jetAlgMap
Map of Rivet::JetAlg to targeted fastjet::JetAlgorithm and fastjet::RecombinationScheme.
Definition FastJets.hh:603
const shared_ptr< fastjet::ClusterSequenceArea > clusterSeqArea() const
Definition FastJets.hh:525
void project(const Event &e)
Perform the projection on the Event.
RIVET_DEFAULT_PROJ_CLONE(FastJets)
Clone on the heap.
static FJPluginPtr mkPlugin(Args &&... args)
Definition FastJets.hh:252
void useJetArea(fastjet::AreaDefinition *adef)
Use provided jet area definition.
Definition FastJets.hh:429
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:530
void addTrf(fastjet::Transformer *trf)
Add a grooming transformer (base class of fastjet::Filter, etc.)
Definition FastJets.hh:448
static Jets mkJets(const PseudoJets &pjs, const Particles &fsparticles=Particles(), const Particles &tagparticles=Particles())
Convert a whole list of PseudoJets to a list of Jets, with mkJet-style unpacking.
PseudoJets pseudojetsByPt(double ptmin=0.0) const
Get the pseudo jets, ordered by .
Definition FastJets.hh:497
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:121
void reset()
Reset the projection. Jet def, etc. are unchanged.
const shared_ptr< fastjet::AreaDefinition > areaDef() const
Return the area definition.
Definition FastJets.hh:538
FastJets(const FinalState &fsp, fastjet::JetAlgorithm type, fastjet::RecombinationScheme recom, double rparameter, fastjet::AreaDefinition *adef, const Cut &c, JetMuons usemuons=JetMuons::ALL, JetInvisibles useinvis=JetInvisibles::NONE)
Definition FastJets.hh:133
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:41
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
const PROJ & get(const std::string &name) const
Definition ProjectionApplier.hh:83
Base class for all Rivet projections.
Definition Projection.hh:29
PseudoJet & ifilterPseudoJets(PseudoJet &pj, const fastjet::Filter &filter)
Apply given FastJet::Filter.
static CONTAINER reclusterJets(const CONTAINER &jetsIn, const fastjet::JetDefinition &jDef, const fastjet::Filter &filter)
Recluster Rivet::Jets.
Definition FastJets.hh:311
static CONTAINER reclusterJets(const CONTAINER &jetsIn, Args &&... args)
Apply the.
Definition FastJets.hh:360
static CONTAINER reclusterJets(const CONTAINER &jetsIn, const fastjet::Filter &filter, Args &&... args)
Apply the.
Definition FastJets.hh:384
static CONTAINER reclusterJets(const CONTAINER &jetsIn, const fastjet::JetDefinition &jDef)
Recluster Rivet::Jets.
Definition FastJets.hh:291
static CONTAINER reclusterJets(const CONTAINER &jetsIn, const FJPluginPtr &jetAlg)
Apply the.
Definition FastJets.hh:327
static CONTAINER reclusterJets(const CONTAINER &jetsIn, const FJPluginPtr &jetAlg, const fastjet::Filter &filter)
Apply the.
Definition FastJets.hh:345
static std::map< T, U > reclusterJets(const std::map< T, U > &jetsMap, Args &&... args)
Apply the.
Definition FastJets.hh:409
static std::map< T, U > reclusterJets(const std::map< T, U > &jetsMap, Args &&... args)
Apply the.
Definition FastJets.hh:399
#define MSG_DEBUG(x)
Debug messaging, not enabled by default, using MSG_LVL.
Definition Logging.hh:182
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
Definition FastJets.hh:607