Rivet  1.8.3
FinalState.hh
1 // -*- C++ -*-
2 #ifndef RIVET_FinalState_HH
3 #define RIVET_FinalState_HH
4 
5 #include "Rivet/Projection.hh"
6 #include "Rivet/Particle.hh"
7 #include "Rivet/Event.hh"
8 
9 namespace Rivet {
10 
11 
14  class FinalState : public Projection {
15  public:
16 
18 
19  FinalState(double mineta = -MAXRAPIDITY,
22  double maxeta = MAXRAPIDITY,
23  double minpt = 0.0*GeV);
24 
27  FinalState(const vector<pair<double, double> >& etaRanges,
28  double minpt = 0.0*GeV);
29 
31  virtual const Projection* clone() const {
32  return new FinalState(*this);
33  }
34 
36 
37 
39  virtual const ParticleVector& particles() const { return _theParticles; }
40 
42  template <typename F>
43  const ParticleVector& particles(F sorter) const {
44  std::sort(_theParticles.begin(), _theParticles.end(), sorter);
45  return _theParticles;
46  }
47 
49  const ParticleVector& particlesByPt() const {
50  return particles(cmpParticleByPt);
51  }
52 
54  const ParticleVector& particlesByP() const {
55  return particles(cmpParticleByP);
56  }
57 
59  const ParticleVector& particlesByE() const {
60  return particles(cmpParticleByE);
61  }
62 
64  const ParticleVector& particlesByEt() const {
65  return particles(cmpParticleByEt);
66  }
67 
69  const ParticleVector& particlesByEta() const {
71  }
72 
74  const ParticleVector& particlesByModEta() const {
76  }
77 
79  const ParticleVector& particlesByRapidity() const {
81  }
82 
84  const ParticleVector& particlesByModRapidity() const {
86  }
87 
89  virtual size_t size() const { return _theParticles.size(); }
90 
92  virtual bool empty() const { return _theParticles.empty(); }
94  virtual bool isEmpty() const { return _theParticles.empty(); }
95 
97  virtual double ptMin() const { return _ptmin; }
98 
99 
100  public:
101 
102  typedef Particle entity_type;
103  typedef ParticleVector collection_type;
104 
106  const collection_type& entities() const {
107  return particles();
108  }
109 
110 
111  protected:
112 
114  virtual void project(const Event& e);
115 
117  virtual int compare(const Projection& p) const;
118 
120  bool accept(const Particle& p) const;
121 
122 
123  protected:
124 
126  vector<pair<double,double> > _etaRanges;
127 
129  double _ptmin;
130 
132  mutable ParticleVector _theParticles;
133 
134  };
135 
136 
137 }
138 
139 #endif
Definition: MC_JetAnalysis.hh:9
const ParticleVector & particlesByE() const
Get the final-state particles, ordered by decreasing .
Definition: FinalState.hh:59
const ParticleVector & particlesByEta() const
Get the final-state particles, ordered by increasing .
Definition: FinalState.hh:69
bool cmpParticleByEt(const Particle &a, const Particle &b)
Sort by descending transverse energy, .
Definition: Particle.hh:170
const ParticleVector & particlesByModRapidity() const
Get the final-state particles, ordered by increasing .
Definition: FinalState.hh:84
bool cmpParticleByAscAbsRapidity(const Particle &a, const Particle &b)
Sort by ascending absolute rapidity, .
Definition: Particle.hh:214
bool cmpParticleByAscRapidity(const Particle &a, const Particle &b)
Sort by ascending rapidity, .
Definition: Particle.hh:206
const ParticleVector & particlesByEt() const
Get the final-state particles, ordered by decreasing .
Definition: FinalState.hh:64
virtual const ParticleVector & particles() const
Get the final-state particles.
Definition: FinalState.hh:39
Representation of particles from a HepMC::GenEvent.
Definition: Particle.hh:16
virtual bool isEmpty() const
Definition: FinalState.hh:94
const ParticleVector & particlesByPt() const
Get the final-state particles, ordered by decreasing .
Definition: FinalState.hh:49
const ParticleVector & particlesByP() const
Get the final-state particles, ordered by decreasing .
Definition: FinalState.hh:54
bool accept(const Particle &p) const
Decide if a particle is to be accepted or not.
Definition: FinalState.cc:95
Definition: Event.hh:22
bool cmpParticleByAscPseudorapidity(const Particle &a, const Particle &b)
Sort by ascending pseudorapidity, .
Definition: Particle.hh:190
virtual size_t size() const
Access the projected final-state particles.
Definition: FinalState.hh:89
Project out all final-state particles in an event. Probably the most important projection in Rivet! ...
Definition: FinalState.hh:14
const collection_type & entities() const
Template-usable interface common to JetAlg.
Definition: FinalState.hh:106
virtual double ptMin() const
Minimum- requirement.
Definition: FinalState.hh:97
bool cmpParticleByE(const Particle &a, const Particle &b)
Sort by descending energy, .
Definition: Particle.hh:178
virtual const Projection * clone() const
Clone on the heap.
Definition: FinalState.hh:31
FinalState(double mineta=-MAXRAPIDITY, double maxeta=MAXRAPIDITY, double minpt=0.0 *GeV)
Definition: FinalState.cc:9
static const double MAXRAPIDITY
A sensible default maximum value of rapidity for Rivet analyses to use.
Definition: Rivet.hh:24
virtual int compare(const Projection &p) const
Compare projections.
Definition: FinalState.cc:44
const ParticleVector & particlesByRapidity() const
Get the final-state particles, ordered by increasing .
Definition: FinalState.hh:79
virtual void project(const Event &e)
Apply the projection to the event.
Definition: FinalState.cc:63
bool cmpParticleByP(const Particle &a, const Particle &b)
Sort by descending momentum, .
Definition: Particle.hh:162
Base class for all Rivet projections.
Definition: Projection.hh:28
bool cmpParticleByPt(const Particle &a, const Particle &b)
Sort by descending transverse momentum, .
Definition: Particle.hh:154
virtual bool empty() const
Is this final state empty?
Definition: FinalState.hh:92
const ParticleVector & particles(F sorter) const
Get the final-state particles, ordered by supplied sorting function object.
Definition: FinalState.hh:43
const ParticleVector & particlesByModEta() const
Get the final-state particles, ordered by increasing .
Definition: FinalState.hh:74
bool cmpParticleByAscAbsPseudorapidity(const Particle &a, const Particle &b)
Sort by ascending absolute pseudorapidity, .
Definition: Particle.hh:198