Rivet  1.8.0
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  // /// The charge of this Particle.
83  // double charge() const {
84  // return PID::charge(*this);
85  // }
86 
87  // /// Three times the charge of this Particle (i.e. integer multiple of smallest quark charge).
88  // int threeCharge() const {
89  // return PID::threeCharge(*this);
90  // }
91 
92 
94  bool hasAncestor(PdgId pdg_id) const;
95 
96 
97  private:
98 
100  const GenParticle* _original;
101 
103  PdgId _id;
104 
106  FourMomentum _momentum;
107  };
108 
109 
111 
112 
114  inline std::string toString(const ParticlePair& pair) {
115  stringstream out;
116  out << "["
117  << toParticleName(pair.first.pdgId()) << " @ "
118  << pair.first.momentum().E()/GeV << " GeV, "
119  << toParticleName(pair.second.pdgId()) << " @ "
120  << pair.second.momentum().E()/GeV << " GeV]";
121  return out.str();
122  }
123 
125  inline std::ostream& operator<<(std::ostream& os, const ParticlePair& pp) {
126  os << toString(pp);
127  return os;
128  }
129 
131 
132 
134 
135 
136  inline bool cmpParticleByPt(const Particle& a, const Particle& b) {
137  return a.momentum().pT() > b.momentum().pT();
138  }
140  inline bool cmpParticleByAscPt(const Particle& a, const Particle& b) {
141  return a.momentum().pT() < b.momentum().pT();
142  }
144  inline bool cmpParticleByP(const Particle& a, const Particle& b) {
145  return a.momentum().vector3().mod() > b.momentum().vector3().mod();
146  }
148  inline bool cmpParticleByAscP(const Particle& a, const Particle& b) {
149  return a.momentum().vector3().mod() < b.momentum().vector3().mod();
150  }
152  inline bool cmpParticleByEt(const Particle& a, const Particle& b) {
153  return a.momentum().Et() > b.momentum().Et();
154  }
156  inline bool cmpParticleByAscEt(const Particle& a, const Particle& b) {
157  return a.momentum().Et() < b.momentum().Et();
158  }
160  inline bool cmpParticleByE(const Particle& a, const Particle& b) {
161  return a.momentum().E() > b.momentum().E();
162  }
164  inline bool cmpParticleByAscE(const Particle& a, const Particle& b) {
165  return a.momentum().E() < b.momentum().E();
166  }
168  inline bool cmpParticleByDescPseudorapidity(const Particle& a, const Particle& b) {
169  return a.momentum().pseudorapidity() > b.momentum().pseudorapidity();
170  }
172  inline bool cmpParticleByAscPseudorapidity(const Particle& a, const Particle& b) {
173  return a.momentum().pseudorapidity() < b.momentum().pseudorapidity();
174  }
176  inline bool cmpParticleByDescAbsPseudorapidity(const Particle& a, const Particle& b) {
177  return fabs(a.momentum().pseudorapidity()) > fabs(b.momentum().pseudorapidity());
178  }
180  inline bool cmpParticleByAscAbsPseudorapidity(const Particle& a, const Particle& b) {
181  return fabs(a.momentum().pseudorapidity()) < fabs(b.momentum().pseudorapidity());
182  }
184  inline bool cmpParticleByDescRapidity(const Particle& a, const Particle& b) {
185  return a.momentum().rapidity() > b.momentum().rapidity();
186  }
188  inline bool cmpParticleByAscRapidity(const Particle& a, const Particle& b) {
189  return a.momentum().rapidity() < b.momentum().rapidity();
190  }
192  inline bool cmpParticleByDescAbsRapidity(const Particle& a, const Particle& b) {
193  return fabs(a.momentum().rapidity()) > fabs(b.momentum().rapidity());
194  }
196  inline bool cmpParticleByAscAbsRapidity(const Particle& a, const Particle& b) {
197  return fabs(a.momentum().rapidity()) < fabs(b.momentum().rapidity());
198  }
200 
201  inline double deltaR(const Particle& p1, const Particle& p2,
202  RapScheme scheme = PSEUDORAPIDITY) {
203  return deltaR(p1.momentum(), p2.momentum(), scheme);
204  }
205 
206  inline double deltaR(const Particle& p, const FourMomentum& v,
207  RapScheme scheme = PSEUDORAPIDITY) {
208  return deltaR(p.momentum(), v, scheme);
209  }
210 
211  inline double deltaR(const Particle& p, const FourVector& v,
212  RapScheme scheme = PSEUDORAPIDITY) {
213  return deltaR(p.momentum(), v, scheme);
214  }
215 
216  inline double deltaR(const Particle& p, const Vector3& v) {
217  return deltaR(p.momentum(), v);
218  }
219 
220  inline double deltaR(const Particle& p, double eta, double phi) {
221  return deltaR(p.momentum(), eta, phi);
222  }
223 
224  inline double deltaR(const FourMomentum& v, const Particle& p,
225  RapScheme scheme = PSEUDORAPIDITY) {
226  return deltaR(v, p.momentum(), scheme);
227  }
228 
229  inline double deltaR(const FourVector& v, const Particle& p,
230  RapScheme scheme = PSEUDORAPIDITY) {
231  return deltaR(v, p.momentum(), scheme);
232  }
233 
234  inline double deltaR(const Vector3& v, const Particle& p) {
235  return deltaR(v, p.momentum());
236  }
237 
238  inline double deltaR(double eta, double phi, const Particle& p) {
239  return deltaR(eta, phi, p.momentum());
240  }
241 
242 
243  inline double deltaPhi(const Particle& p1, const Particle& p2) {
244  return deltaPhi(p1.momentum(), p2.momentum());
245  }
246 
247  inline double deltaPhi(const Particle& p, const FourMomentum& v) {
248  return deltaPhi(p.momentum(), v);
249  }
250 
251  inline double deltaPhi(const Particle& p, const FourVector& v) {
252  return deltaPhi(p.momentum(), v);
253  }
254 
255  inline double deltaPhi(const Particle& p, const Vector3& v) {
256  return deltaPhi(p.momentum(), v);
257  }
258 
259  inline double deltaPhi(const Particle& p, double phi) {
260  return deltaPhi(p.momentum(), phi);
261  }
262 
263  inline double deltaPhi(const FourMomentum& v, const Particle& p) {
264  return deltaPhi(v, p.momentum());
265  }
266 
267  inline double deltaPhi(const FourVector& v, const Particle& p) {
268  return deltaPhi(v, p.momentum());
269  }
270 
271  inline double deltaPhi(const Vector3& v, const Particle& p) {
272  return deltaPhi(v, p.momentum());
273  }
274 
275  inline double deltaPhi(double phi, const Particle& p) {
276  return deltaPhi(phi, p.momentum());
277  }
278 
279 
280  inline double deltaEta(const Particle& p1, const Particle& p2) {
281  return deltaEta(p1.momentum(), p2.momentum());
282  }
283 
284  inline double deltaEta(const Particle& p, const FourMomentum& v) {
285  return deltaEta(p.momentum(), v);
286  }
287 
288  inline double deltaEta(const Particle& p, const FourVector& v) {
289  return deltaEta(p.momentum(), v);
290  }
291 
292  inline double deltaEta(const Particle& p, const Vector3& v) {
293  return deltaEta(p.momentum(), v);
294  }
295 
296  inline double deltaEta(const Particle& p, double eta) {
297  return deltaEta(p.momentum(), eta);
298  }
299 
300  inline double deltaEta(const FourMomentum& v, const Particle& p) {
301  return deltaEta(v, p.momentum());
302  }
303 
304  inline double deltaEta(const FourVector& v, const Particle& p) {
305  return deltaEta(v, p.momentum());
306  }
307 
308  inline double deltaEta(const Vector3& v, const Particle& p) {
309  return deltaEta(v, p.momentum());
310  }
311 
312  inline double deltaEta(double eta, const Particle& p) {
313  return deltaEta(eta, p.momentum());
314  }
315 
316 }
317 
318 #endif