Rivet  1.8.3
Jet.hh
1 // -*- C++ -*-
2 #ifndef RIVET_Jet_HH
3 #define RIVET_Jet_HH
4 
5 #include "Rivet/Rivet.hh"
6 #include <numeric>
7 
8 namespace Rivet {
9 
10 
12  class Jet : public ParticleBase {
13  public:
14 
16 
17 
18  Jet() : ParticleBase() { clear(); }
19 
21  Jet(const vector<Particle>& particles, const FourMomentum& pjet)
22  : ParticleBase() {
23  setState(particles, pjet);
24  }
25 
26  // /// Set all the jet data, without particle ID information.
27  // Jet(const vector<FourMomentum>& momenta, const FourMomentum& pjet)
28  // : ParticleBase() {
29  // setState(momenta, pjet);
30  // }
31 
33 
34 
36 
37 
39  size_t size() const { return _particles.size(); }
40 
41  // /// Define a Jet::iterator via a typedef.
42  // typedef vector<FourMomentum>::iterator iterator;
43 
44  // /// Define a Jet::const_iterator via a typedef.
45  // typedef vector<FourMomentum>::const_iterator const_iterator;
46 
47  // /// Get a begin iterator over the particle/track four-momenta in this jet.
48  // iterator begin() {
49  // return _momenta.begin();
50  // }
51 
52  // /// Get an end iterator over the particle/track four-momenta in this jet.
53  // iterator end() {
54  // return _momenta.end();
55  // }
56 
57  // /// Get a const begin iterator over the particle/track four-momenta in this jet.
58  // const_iterator begin() const {
59  // return _momenta.begin();
60  // }
61 
62  // /// Get a const end iterator over the particle/track four-momenta in this jet.
63  // const_iterator end() const {
64  // return _momenta.end();
65  // }
66 
67  // /// Get the track momenta in this jet.
68  // vector<FourMomentum>& momenta() {
69  // return _momenta;
70  // }
71 
72  // /// Get the track momenta in this jet (const version).
73  // const vector<FourMomentum>& momenta() const {
74  // return _momenta;
75  // }
76 
78  vector<Particle>& particles() { return _particles; }
79 
81  const vector<Particle>& particles() const { return _particles; }
82 
84  bool containsParticle(const Particle& particle) const;
85 
87  bool containsParticleId(PdgId pid) const;
88 
90  bool containsParticleId(const vector<PdgId>& pids) const;
91 
93  bool containsCharm() const;
94 
96  bool containsBottom() const;
97 
99 
100 
102 
103 
105  const FourMomentum& momentum() const { return _momentum; }
106 
108  double eta() const { return momentum().eta(); }
109 
111  double phi() const { return momentum().phi(); }
112 
114  double totalEnergy() const { return momentum().E(); }
115 
117  double neutralEnergy() const;
118 
120  double hadronicEnergy() const;
121 
123  double ptSum() const { return momentum().pT(); }
124 
126  double EtSum() const { return momentum().Et(); }
127 
129 
130 
132 
133 
135  Jet& setState(const vector<Particle>& particles, const FourMomentum& pjet);
136 
137  // /// Set all the jet data, without particle ID information.
138  // Jet& setState(const vector<FourMomentum>& momenta, const FourMomentum& pjet);
139 
142 
144  Jet& setParticles(const vector<Particle>& particles);
145 
146  // /// Set the particles collection with momentum information only.
147  // Jet& setParticles(const vector<FourMomentum>& momenta);
148 
149  // /// Add a particle/track to this jet.
150  // Jet& addParticle(const FourMomentum& particle);
151 
152  // /// Add a particle/track to this jet.
153  // Jet& addParticle(const Particle& particle);
154 
156  Jet& clear();
157 
159 
160 
161  private:
162 
164  ParticleVector _particles;
165 
166  // /// The particle momenta.
167  // /// @todo Eliminate this to ensure consistency.
168  // std::vector<FourMomentum> _momenta;
169 
171  FourMomentum _momentum;
172 
173  };
174 
175 
177  typedef std::vector<Jet> Jets;
178 
179 
181 
182 
185  inline bool cmpJetsByPt(const Jet& a, const Jet& b) {
186  return a.ptSum() > b.ptSum();
187  }
190  inline bool cmpJetsByAscPt(const Jet& a, const Jet& b) {
191  return a.ptSum() < b.ptSum();
192  }
193 
195  inline bool cmpJetsByP(const Jet& a, const Jet& b) {
196  return a.momentum().vector3().mod() > b.momentum().vector3().mod();
197  }
199  inline bool cmpJetsByAscP(const Jet& a, const Jet& b) {
200  return a.momentum().vector3().mod() < b.momentum().vector3().mod();
201  }
202 
203  // @brief Compare jets by \f$ E_\perp \f$ (descending - usual sorting for HEP)
205  inline bool cmpJetsByEt(const Jet& a, const Jet& b) {
206  return a.EtSum() > b.EtSum();
207  }
208  // @brief Compare jets by \f$ E_\perp \f$ (ascending)
210  inline bool cmpJetsByEtDesc(const Jet& a, const Jet& b) {
211  return a.EtSum() < b.EtSum();
212  }
213 
216  inline bool cmpJetsByE(const Jet& a, const Jet& b) {
217  return a.momentum().E() > b.momentum().E();
218  }
221  inline bool cmpJetsByAscE(const Jet& a, const Jet& b) {
222  return a.momentum().E() < b.momentum().E();
223  }
224 
227  inline bool cmpJetsByDescPseudorapidity(const Jet& a, const Jet& b) {
228  return a.momentum().pseudorapidity() > b.momentum().pseudorapidity();
229  }
232  inline bool cmpJetsByAscPseudorapidity(const Jet& a, const Jet& b) {
233  return a.momentum().pseudorapidity() < b.momentum().pseudorapidity();
234  }
235 
238  inline bool cmpJetsByDescAbsPseudorapidity(const Jet& a, const Jet& b) {
239  return fabs(a.momentum().pseudorapidity()) > fabs(b.momentum().pseudorapidity());
240  }
243  inline bool cmpJetsByAscAbsPseudorapidity(const Jet& a, const Jet& b) {
244  return fabs(a.momentum().pseudorapidity()) < fabs(b.momentum().pseudorapidity());
245  }
246 
249  inline bool cmpJetsByDescRapidity(const Jet& a, const Jet& b) {
250  return a.momentum().rapidity() > b.momentum().rapidity();
251  }
254  inline bool cmpJetsByAscRapidity(const Jet& a, const Jet& b) {
255  return a.momentum().rapidity() < b.momentum().rapidity();
256  }
257 
260  inline bool cmpJetsByDescAbsRapidity(const Jet& a, const Jet& b) {
261  return fabs(a.momentum().rapidity()) > fabs(b.momentum().rapidity());
262  }
265  inline bool cmpJetsByAscAbsRapidity(const Jet& a, const Jet& b) {
266  return fabs(a.momentum().rapidity()) < fabs(b.momentum().rapidity());
267  }
268 
270 
271  inline double deltaR(const Jet& j1, const Jet& j2,
272  RapScheme scheme = PSEUDORAPIDITY) {
273  return deltaR(j1.momentum(), j2.momentum(), scheme);
274  }
275 
276  inline double deltaR(const Jet& j, const Particle& p,
277  RapScheme scheme = PSEUDORAPIDITY) {
278  return deltaR(j.momentum(), p.momentum(), scheme);
279  }
280 
281  inline double deltaR(const Particle& p, const Jet& j,
282  RapScheme scheme = PSEUDORAPIDITY) {
283  return deltaR(p.momentum(), j.momentum(), scheme);
284  }
285 
286  inline double deltaR(const Jet& j, const FourMomentum& v,
287  RapScheme scheme = PSEUDORAPIDITY) {
288  return deltaR(j.momentum(), v, scheme);
289  }
290 
291  inline double deltaR(const Jet& j, const FourVector& v,
292  RapScheme scheme = PSEUDORAPIDITY) {
293  return deltaR(j.momentum(), v, scheme);
294  }
295 
296  inline double deltaR(const Jet& j, const Vector3& v) {
297  return deltaR(j.momentum(), v);
298  }
299 
300  inline double deltaR(const Jet& j, double eta, double phi) {
301  return deltaR(j.momentum(), eta, phi);
302  }
303 
304  inline double deltaR(const FourMomentum& v, const Jet& j,
305  RapScheme scheme = PSEUDORAPIDITY) {
306  return deltaR(v, j.momentum(), scheme);
307  }
308 
309  inline double deltaR(const FourVector& v, const Jet& j,
310  RapScheme scheme = PSEUDORAPIDITY) {
311  return deltaR(v, j.momentum(), scheme);
312  }
313 
314  inline double deltaR(const Vector3& v, const Jet& j) {
315  return deltaR(v, j.momentum());
316  }
317 
318  inline double deltaR(double eta, double phi, const Jet& j) {
319  return deltaR(eta, phi, j.momentum());
320  }
321 
322 
323  inline double deltaPhi(const Jet& j1, const Jet& j2) {
324  return deltaPhi(j1.momentum(), j2.momentum());
325  }
326 
327  inline double deltaPhi(const Jet& j, const Particle& p) {
328  return deltaPhi(j.momentum(), p.momentum());
329  }
330 
331  inline double deltaPhi(const Particle& p, const Jet& j) {
332  return deltaPhi(p.momentum(), j.momentum());
333  }
334 
335  inline double deltaPhi(const Jet& j, const FourMomentum& v) {
336  return deltaPhi(j.momentum(), v);
337  }
338 
339  inline double deltaPhi(const Jet& j, const FourVector& v) {
340  return deltaPhi(j.momentum(), v);
341  }
342 
343  inline double deltaPhi(const Jet& j, const Vector3& v) {
344  return deltaPhi(j.momentum(), v);
345  }
346 
347  inline double deltaPhi(const Jet& j, double phi) {
348  return deltaPhi(j.momentum(), phi);
349  }
350 
351  inline double deltaPhi(const FourMomentum& v, const Jet& j) {
352  return deltaPhi(v, j.momentum());
353  }
354 
355  inline double deltaPhi(const FourVector& v, const Jet& j) {
356  return deltaPhi(v, j.momentum());
357  }
358 
359  inline double deltaPhi(const Vector3& v, const Jet& j) {
360  return deltaPhi(v, j.momentum());
361  }
362 
363  inline double deltaPhi(double phi, const Jet& j) {
364  return deltaPhi(phi, j.momentum());
365  }
366 
367 
368  inline double deltaEta(const Jet& j1, const Jet& j2) {
369  return deltaEta(j1.momentum(), j2.momentum());
370  }
371 
372  inline double deltaEta(const Jet& j, const Particle& p) {
373  return deltaEta(j.momentum(), p.momentum());
374  }
375 
376  inline double deltaEta(const Particle& p, const Jet& j) {
377  return deltaEta(p.momentum(), j.momentum());
378  }
379 
380  inline double deltaEta(const Jet& j, const FourMomentum& v) {
381  return deltaEta(j.momentum(), v);
382  }
383 
384  inline double deltaEta(const Jet& j, const FourVector& v) {
385  return deltaEta(j.momentum(), v);
386  }
387 
388  inline double deltaEta(const Jet& j, const Vector3& v) {
389  return deltaEta(j.momentum(), v);
390  }
391 
392  inline double deltaEta(const Jet& j, double eta) {
393  return deltaEta(j.momentum(), eta);
394  }
395 
396  inline double deltaEta(const FourMomentum& v, const Jet& j) {
397  return deltaEta(v, j.momentum());
398  }
399 
400  inline double deltaEta(const FourVector& v, const Jet& j) {
401  return deltaEta(v, j.momentum());
402  }
403 
404  inline double deltaEta(const Vector3& v, const Jet& j) {
405  return deltaEta(v, j.momentum());
406  }
407 
408  inline double deltaEta(double eta, const Jet& j) {
409  return deltaEta(eta, j.momentum());
410  }
411 
412 }
413 
414 #endif
bool cmpJetsByEt(const Jet &a, const Jet &b)
Use this so that highest is at the front of the list.
Definition: Jet.hh:205
Definition: MC_JetAnalysis.hh:9
size_t size() const
Number of particles in this jet.
Definition: Jet.hh:39
double pT() const
Calculate the transverse momentum .
Definition: Vector4.hh:410
bool cmpJetsByAscRapidity(const Jet &a, const Jet &b)
Compare jets by (ascending) Use this so that lowest is at the front of the list.
Definition: Jet.hh:254
double eta() const
Synonym for pseudorapidity.
Definition: Vector4.hh:134
bool cmpJetsByAscE(const Jet &a, const Jet &b)
Compare jets by (ascending) Use this so that lowest is at the front of the list.
Definition: Jet.hh:221
bool cmpJetsByE(const Jet &a, const Jet &b)
Compare jets by (descending - usual sorting for HEP) Use this so that highest is at the front of th...
Definition: Jet.hh:216
bool cmpJetsByP(const Jet &a, const Jet &b)
Compare jets by descending momentum, .
Definition: Jet.hh:195
double Et() const
Calculate the transverse energy .
Definition: Vector4.hh:420
bool containsParticleId(PdgId pid) const
Check whether this jet contains a certain particle type.
Definition: Jet.cc:67
double neutralEnergy() const
Get the energy carried in this jet by neutral particles.
Definition: Jet.cc:88
bool cmpJetsByDescAbsRapidity(const Jet &a, const Jet &b)
Compare jets by (descending) Use this so that highest is at the front of the list.
Definition: Jet.hh:260
bool containsBottom() const
Check whether this jet contains a bottom-flavoured hadron (or decay products from one)...
Definition: Jet.cc:129
std::vector< Jet > Jets
Typedef for a collection of Jet objects.
Definition: Jet.hh:177
double phi() const
Get the unweighted average for this jet. (caches)
Definition: Jet.hh:111
vector< Particle > & particles()
Get the particles in this jet.
Definition: Jet.hh:78
Jet(const vector< Particle > &particles, const FourMomentum &pjet)
Set all the jet data, with full particle information.
Definition: Jet.hh:21
double rapidity() const
Calculate the rapidity.
Definition: Vector4.hh:400
Representation of particles from a HepMC::GenEvent.
Definition: Particle.hh:16
double pseudorapidity() const
Pseudorapidity (defined purely by the 3-vector components)
Definition: Vector4.hh:129
Jet & setState(const vector< Particle > &particles, const FourMomentum &pjet)
Set all the jet data, with full particle information.
Definition: Jet.cc:10
bool cmpJetsByEtDesc(const Jet &a, const Jet &b)
Use this so that lowest is at the front of the list.
Definition: Jet.hh:210
bool cmpJetsByAscAbsRapidity(const Jet &a, const Jet &b)
Compare jets by (ascending) Use this so that lowest is at the front of the list.
Definition: Jet.hh:265
bool cmpJetsByAscPt(const Jet &a, const Jet &b)
Compare jets by (ascending) Use this so that lowest is at the front of the list.
Definition: Jet.hh:190
Jet & clear()
Reset this jet as empty.
Definition: Jet.cc:146
bool cmpJetsByAscP(const Jet &a, const Jet &b)
Compare jets by ascending momentum, .
Definition: Jet.hh:199
bool cmpJetsByPt(const Jet &a, const Jet &b)
Compare jets by (descending - usual sorting for HEP) Use this so that highest is at the front of th...
Definition: Jet.hh:185
Jet & setMomentum(const FourMomentum &momentum)
Set the effective 4-momentum of the jet.
Definition: Jet.cc:24
const FourMomentum & momentum() const
Get equivalent single momentum four-vector.
Definition: Jet.hh:105
double eta() const
Get the unweighted average for this jet. (caches)
Definition: Jet.hh:108
Specialisation of VectorN to a general (non-momentum) Lorentz 4-vector.
Definition: Vector4.hh:20
bool cmpJetsByDescRapidity(const Jet &a, const Jet &b)
Compare jets by (descending) Use this so that highest is at the front of the list.
Definition: Jet.hh:249
const FourMomentum & momentum() const
The momentum of this Particle.
Definition: Particle.hh:68
bool containsParticle(const Particle &particle) const
Check whether this jet contains a particular particle.
Definition: Jet.cc:58
double EtSum() const
Get the sum of the values of the constituent tracks. (caches)
Definition: Jet.hh:126
bool cmpJetsByDescPseudorapidity(const Jet &a, const Jet &b)
Compare jets by (descending) Use this so that highest is at the front of the list.
Definition: Jet.hh:227
Base class for particle-like things like Particle and Jet.
Definition: ParticleBase.hh:11
const vector< Particle > & particles() const
Get the particles in this jet (const version)
Definition: Jet.hh:81
double E() const
Get energy (time component of momentum).
Definition: Vector4.hh:355
RapScheme
Enum for rapidity variable to be used in calculating , applying rapidity cuts, etc.
Definition: MathHeader.hh:60
Vector3 vector3() const
Get the spatial part of the 4-vector as a 3-vector.
Definition: Vector4.hh:139
double ptSum() const
Get the sum of the values of the constituent tracks. (caches)
Definition: Jet.hh:123
Representation of a clustered jet of particles.
Definition: Jet.hh:12
double phi(const PhiMapping mapping=ZERO_2PI) const
Synonym for azimuthalAngle.
Definition: Vector4.hh:114
Three-dimensional specialisation of Vector.
Definition: Vector3.hh:26
Jet & setParticles(const vector< Particle > &particles)
Set the particles collection with full particle information.
Definition: Jet.cc:30
bool containsCharm() const
Check whether this jet contains a charm-flavoured hadron (or decay products from one).
Definition: Jet.cc:112
bool cmpJetsByAscAbsPseudorapidity(const Jet &a, const Jet &b)
Compare jets by (ascending) Use this so that lowest is at the front of the list.
Definition: Jet.hh:243
bool cmpJetsByAscPseudorapidity(const Jet &a, const Jet &b)
Compare jets by (ascending) Use this so that lowest is at the front of the list.
Definition: Jet.hh:232
Specialized version of the FourVector with momentum/energy functionality.
Definition: Vector4.hh:324
double hadronicEnergy() const
Get the energy carried in this jet by hadrons.
Definition: Jet.cc:100
double totalEnergy() const
Get the total energy of this jet.
Definition: Jet.hh:114
double mod() const
Calculate the modulus of a vector. .
Definition: VectorN.hh:87
bool cmpJetsByDescAbsPseudorapidity(const Jet &a, const Jet &b)
Compare jets by (descending) Use this so that highest is at the front of the list.
Definition: Jet.hh:238