Rivet  1.8.3
Particle.hh
1 // -*- C++ -*-
2 #ifndef RIVET_Particle_HH
3 #define RIVET_Particle_HH
4 
5 #include "Rivet/Rivet.hh"
6 #include "Rivet/Particle.fhh"
7 #include "Rivet/ParticleBase.hh"
8 #include "Rivet/ParticleName.hh"
9 #include "Rivet/Math/Vectors.hh"
10 #include "Rivet/Tools/Logging.hh"
11 
12 namespace Rivet {
13 
14 
16  class Particle : public ParticleBase {
17  public:
18 
22  : ParticleBase(),
23  _original(0), _id(0), _momentum()
24  { }
25 
27  Particle(PdgId pid, const FourMomentum& mom)
28  : ParticleBase(),
29  _original(0), _id(pid), _momentum(mom)
30  { }
31 
33  Particle(const GenParticle& gp)
34  : ParticleBase(),
35  _original(&gp), _id(gp.pdg_id()),
36  _momentum(gp.momentum())
37  { }
38 
39 
40  public:
41 
43  const GenParticle& genParticle() const {
44  assert(_original);
45  return *_original;
46  }
47 
48 
50  bool hasGenParticle() const {
51  return bool(_original);
52  }
53 
54 
56  PdgId pdgId() const {
57  return _id;
58  }
59 
60 
63  _momentum = momentum;
64  return *this;
65  }
66 
68  const FourMomentum& momentum() const {
69  return _momentum;
70  }
71 
73  double energy() const {
74  return momentum().E();
75  }
76 
78  double mass() const {
79  return momentum().mass();
80  }
81 
82 
84  // /// The charge of this Particle.
85  // double charge() const {
86  // return PID::charge(*this);
87  // }
88 
90  // /// Three times the charge of this Particle (i.e. integer multiple of smallest quark charge).
91  // int threeCharge() const {
92  // return PID::threeCharge(*this);
93  // }
94 
100  bool hasAncestor(PdgId pdg_id) const;
101 
110  bool fromDecay() const;
111 
113 
114 
115  private:
116 
118  const GenParticle* _original;
119 
121  PdgId _id;
122 
124  FourMomentum _momentum;
125  };
126 
127 
129 
130 
132  inline std::string toString(const ParticlePair& pair) {
133  stringstream out;
134  out << "["
135  << toParticleName(pair.first.pdgId()) << " @ "
136  << pair.first.momentum().E()/GeV << " GeV, "
137  << toParticleName(pair.second.pdgId()) << " @ "
138  << pair.second.momentum().E()/GeV << " GeV]";
139  return out.str();
140  }
141 
143  inline std::ostream& operator<<(std::ostream& os, const ParticlePair& pp) {
144  os << toString(pp);
145  return os;
146  }
147 
149 
150 
152 
153 
154  inline bool cmpParticleByPt(const Particle& a, const Particle& b) {
155  return a.momentum().pT() > b.momentum().pT();
156  }
158  inline bool cmpParticleByAscPt(const Particle& a, const Particle& b) {
159  return a.momentum().pT() < b.momentum().pT();
160  }
162  inline bool cmpParticleByP(const Particle& a, const Particle& b) {
163  return a.momentum().vector3().mod() > b.momentum().vector3().mod();
164  }
166  inline bool cmpParticleByAscP(const Particle& a, const Particle& b) {
167  return a.momentum().vector3().mod() < b.momentum().vector3().mod();
168  }
170  inline bool cmpParticleByEt(const Particle& a, const Particle& b) {
171  return a.momentum().Et() > b.momentum().Et();
172  }
174  inline bool cmpParticleByAscEt(const Particle& a, const Particle& b) {
175  return a.momentum().Et() < b.momentum().Et();
176  }
178  inline bool cmpParticleByE(const Particle& a, const Particle& b) {
179  return a.momentum().E() > b.momentum().E();
180  }
182  inline bool cmpParticleByAscE(const Particle& a, const Particle& b) {
183  return a.momentum().E() < b.momentum().E();
184  }
186  inline bool cmpParticleByDescPseudorapidity(const Particle& a, const Particle& b) {
187  return a.momentum().pseudorapidity() > b.momentum().pseudorapidity();
188  }
190  inline bool cmpParticleByAscPseudorapidity(const Particle& a, const Particle& b) {
191  return a.momentum().pseudorapidity() < b.momentum().pseudorapidity();
192  }
194  inline bool cmpParticleByDescAbsPseudorapidity(const Particle& a, const Particle& b) {
195  return fabs(a.momentum().pseudorapidity()) > fabs(b.momentum().pseudorapidity());
196  }
198  inline bool cmpParticleByAscAbsPseudorapidity(const Particle& a, const Particle& b) {
199  return fabs(a.momentum().pseudorapidity()) < fabs(b.momentum().pseudorapidity());
200  }
202  inline bool cmpParticleByDescRapidity(const Particle& a, const Particle& b) {
203  return a.momentum().rapidity() > b.momentum().rapidity();
204  }
206  inline bool cmpParticleByAscRapidity(const Particle& a, const Particle& b) {
207  return a.momentum().rapidity() < b.momentum().rapidity();
208  }
210  inline bool cmpParticleByDescAbsRapidity(const Particle& a, const Particle& b) {
211  return fabs(a.momentum().rapidity()) > fabs(b.momentum().rapidity());
212  }
214  inline bool cmpParticleByAscAbsRapidity(const Particle& a, const Particle& b) {
215  return fabs(a.momentum().rapidity()) < fabs(b.momentum().rapidity());
216  }
218 
219  inline double deltaR(const Particle& p1, const Particle& p2,
220  RapScheme scheme = PSEUDORAPIDITY) {
221  return deltaR(p1.momentum(), p2.momentum(), scheme);
222  }
223 
224  inline double deltaR(const Particle& p, const FourMomentum& v,
225  RapScheme scheme = PSEUDORAPIDITY) {
226  return deltaR(p.momentum(), v, scheme);
227  }
228 
229  inline double deltaR(const Particle& p, const FourVector& v,
230  RapScheme scheme = PSEUDORAPIDITY) {
231  return deltaR(p.momentum(), v, scheme);
232  }
233 
234  inline double deltaR(const Particle& p, const Vector3& v) {
235  return deltaR(p.momentum(), v);
236  }
237 
238  inline double deltaR(const Particle& p, double eta, double phi) {
239  return deltaR(p.momentum(), eta, phi);
240  }
241 
242  inline double deltaR(const FourMomentum& v, const Particle& p,
243  RapScheme scheme = PSEUDORAPIDITY) {
244  return deltaR(v, p.momentum(), scheme);
245  }
246 
247  inline double deltaR(const FourVector& v, const Particle& p,
248  RapScheme scheme = PSEUDORAPIDITY) {
249  return deltaR(v, p.momentum(), scheme);
250  }
251 
252  inline double deltaR(const Vector3& v, const Particle& p) {
253  return deltaR(v, p.momentum());
254  }
255 
256  inline double deltaR(double eta, double phi, const Particle& p) {
257  return deltaR(eta, phi, p.momentum());
258  }
259 
260 
261  inline double deltaPhi(const Particle& p1, const Particle& p2) {
262  return deltaPhi(p1.momentum(), p2.momentum());
263  }
264 
265  inline double deltaPhi(const Particle& p, const FourMomentum& v) {
266  return deltaPhi(p.momentum(), v);
267  }
268 
269  inline double deltaPhi(const Particle& p, const FourVector& v) {
270  return deltaPhi(p.momentum(), v);
271  }
272 
273  inline double deltaPhi(const Particle& p, const Vector3& v) {
274  return deltaPhi(p.momentum(), v);
275  }
276 
277  inline double deltaPhi(const Particle& p, double phi) {
278  return deltaPhi(p.momentum(), phi);
279  }
280 
281  inline double deltaPhi(const FourMomentum& v, const Particle& p) {
282  return deltaPhi(v, p.momentum());
283  }
284 
285  inline double deltaPhi(const FourVector& v, const Particle& p) {
286  return deltaPhi(v, p.momentum());
287  }
288 
289  inline double deltaPhi(const Vector3& v, const Particle& p) {
290  return deltaPhi(v, p.momentum());
291  }
292 
293  inline double deltaPhi(double phi, const Particle& p) {
294  return deltaPhi(phi, p.momentum());
295  }
296 
297 
298  inline double deltaEta(const Particle& p1, const Particle& p2) {
299  return deltaEta(p1.momentum(), p2.momentum());
300  }
301 
302  inline double deltaEta(const Particle& p, const FourMomentum& v) {
303  return deltaEta(p.momentum(), v);
304  }
305 
306  inline double deltaEta(const Particle& p, const FourVector& v) {
307  return deltaEta(p.momentum(), v);
308  }
309 
310  inline double deltaEta(const Particle& p, const Vector3& v) {
311  return deltaEta(p.momentum(), v);
312  }
313 
314  inline double deltaEta(const Particle& p, double eta) {
315  return deltaEta(p.momentum(), eta);
316  }
317 
318  inline double deltaEta(const FourMomentum& v, const Particle& p) {
319  return deltaEta(v, p.momentum());
320  }
321 
322  inline double deltaEta(const FourVector& v, const Particle& p) {
323  return deltaEta(v, p.momentum());
324  }
325 
326  inline double deltaEta(const Vector3& v, const Particle& p) {
327  return deltaEta(v, p.momentum());
328  }
329 
330  inline double deltaEta(double eta, const Particle& p) {
331  return deltaEta(eta, p.momentum());
332  }
333 
334 }
335 
336 #endif