HepMC event record
include/HepMC/FourVector.h
1 // -*- C++ -*-
2 //
3 // This file is part of HepMC
4 // Copyright (C) 2014-2015 The HepMC collaboration (see AUTHORS for details)
5 //
6 #ifndef HEPMC_FOURVECTOR_H
7 #define HEPMC_FOURVECTOR_H
8 /**
9  * @file FourVector.h
10  * @brief Definition of \b class FourVector
11  */
12 #include "HepMC/Common.h"
13 #include <cmath>
14 
15 namespace HepMC {
16 
17 
18 /**
19  * @brief Generic 4-vector
20  *
21  * Interpretation of its content depends on accessors used: it's much simpler to do this
22  * than to distinguish between space and momentum vectors via the type system (especially
23  * given the need for backward compatibility with HepMC2). Be sensible and don't call
24  * energy functions on spatial vectors! To avoid duplication, most definitions are only
25  * implemented on the spatial function names, with the energy-momentum functions as aliases.
26  *
27  * This is @a not intended to be a fully featured 4-vector, but does contain the majority
28  * of common non-boosting functionality, as well as a few support operations on
29  * 4-vectors.
30  *
31  * The implementations in this class are fully inlined.
32  */
33 class FourVector {
34 public:
35 
36  /** @brief Default constructor */
38  : m_v1(0.0), m_v2(0.0), m_v3(0.0), m_v4(0.0) {}
39  /** @brief Sets all FourVector fields */
40  FourVector(double xx, double yy, double zz, double ee)
41  : m_v1(xx), m_v2(yy), m_v3(zz), m_v4(ee) {}
42  /** @brief Copy constructor */
44  : m_v1(v.m_v1), m_v2(v.m_v2), m_v3(v.m_v3), m_v4(v.m_v4) {}
45 
46 
47  /// @name Component accessors
48  //@{
49 
50  /** @brief Set all FourVector fields, in order x,y,z,t */
51  void set(double x1, double x2, double x3, double x4) {
52  m_v1 = x1;
53  m_v2 = x2;
54  m_v3 = x3;
55  m_v4 = x4;
56  }
57 
58 
59  /// x-component of position/displacement
60  double x() const { return m_v1; }
61  /// Set x-component of position/displacement
62  void setX(double xx) { m_v1 = xx; }
63 
64  /// y-component of position/displacement
65  double y() const { return m_v2; }
66  /// Set y-component of position/displacement
67  void setY(double yy) { m_v2 = yy; }
68 
69  /// z-component of position/displacement
70  double z() const { return m_v3; }
71  /// Set z-component of position/displacement
72  void setZ(double zz) { m_v3 = zz; }
73 
74  /// Time component of position/displacement
75  double t() const { return m_v4; }
76  /// Set time component of position/displacement
77  void setT(double tt) { m_v4 = tt; }
78 
79 
80  /// x-component of momentum
81  double px() const { return x(); }
82  /// Set x-component of momentum
83  void setPx(double pxx) { setX(pxx); }
84 
85  /// y-component of momentum
86  double py() const { return y(); }
87  /// Set y-component of momentum
88  void setPy(double pyy) { setY(pyy); }
89 
90  /// z-component of momentum
91  double pz() const { return z(); }
92  /// Set z-component of momentum
93  void setPz(double pzz) { setZ(pzz); }
94 
95  /// Energy component of momentum
96  double e() const { return t(); }
97  /// Set energy component of momentum
98  void setE (double ee ) { setT(ee); }
99 
100  //@}
101 
102 
103  /// @name Computed properties
104  //@{
105 
106  /// Squared magnitude of (x, y, z) 3-vector
107  double length2() const { return sqr(x()) + sqr(y()) + sqr(z()); }
108  /// Magnitude of spatial (x, y, z) 3-vector
109  double length() const { return sqrt(length2()); }
110  /// Squared magnitude of (x, y) vector
111  double perp2() const { return sqr(x()) + sqr(y()); }
112  /// Magnitude of (x, y) vector
113  double perp() const { return sqrt(perp2()); }
114  /// Spacetime invariant interval s^2 = t^2 - x^2 - y^2 - z^2
115  double interval() const { return sqr(t()) - length2(); }
116 
117  /// Squared magnitude of p3 = (px, py, pz) vector
118  double p3mod2() const { return length2(); }
119  /// Magnitude of p3 = (px, py, pz) vector
120  double p3mod() const { return length(); }
121  /// Squared transverse momentum px^2 + py^2
122  double pt2() const { return perp2(); }
123  /// Transverse momentum
124  double pt() const { return perp(); }
125  /// Squared invariant mass m^2 = E^2 - px^2 - py^2 - pz^2
126  double m2() const { return interval(); }
127  /// Invariant mass. Returns -sqrt(-m) if e^2 - P^2 is negative
128  double m() const { return (m2() > 0.0) ? sqrt(m2()) : -sqrt(-m2()); }
129 
130  /// Azimuthal angle
131  double phi() const { return atan2( y(), x() ); }
132  /// Polar angle w.r.t. z direction
133  double theta() const { return atan2( perp(), z() ); }
134  /// Pseudorapidity
135  /// @todo Improve numerical stability
136  double eta() const { return 0.5*log( (p3mod() + pz()) / (p3mod() - pz()) ); }
137  /// Rapidity
138  /// @todo Improve numerical stability
139  double rap() const { return 0.5*log( (e() + pz()) / (e() - pz()) ); }
140  /// Absolute pseudorapidity
141  double abs_eta() const { return std::abs( eta() ); }
142  /// Absolute rapidity
143  double abs_rap() const { return std::abs( rap() ); }
144 
145  #ifndef HEPMC_NO_DEPRECATED
146  /// Same as eta
147  double pseudoRapidity() const { return eta(); }
148  #endif
149 
150  //@}
151 
152 
153  /// @name Comparisons to another FourVector
154  //@{
155 
156  /// Check if the length of this vertex is zero
157  bool is_zero() const { return x() == 0 && y() == 0 && z() == 0 && t() == 0; }
158 
159  /// Signed azimuthal angle separation in [-pi, pi]
160  double delta_phi(const FourVector &v) const {
161  double dphi = phi() - v.phi();
162  if (dphi != dphi) return dphi;
163  while (dphi >= M_PI) dphi -= 2.*M_PI;
164  while (dphi < -M_PI) dphi += 2.*M_PI;
165  return dphi;
166  }
167 
168  /// Pseudorapidity separation
169  double delta_eta(const FourVector &v) const { return eta() - v.eta(); }
170 
171  /// Rapidity separation
172  double delta_rap(const FourVector &v) const { return rap() - v.rap(); }
173 
174  /// R_eta^2-distance separation dR^2 = dphi^2 + deta^2
175  double delta_r2_eta(const FourVector &v) const {
176  return sqr(delta_phi(v)) + sqr(delta_eta(v));
177  }
178 
179  /// R_eta-distance separation dR = sqrt(dphi^2 + deta^2)
180  double delta_r_eta(const FourVector &v) const {
181  return sqrt( delta_r2_eta(v) );
182  }
183 
184  /// R_rap^2-distance separation dR^2 = dphi^2 + drap^2
185  double delta_r2_rap(const FourVector &v) const {
186  return sqr(delta_phi(v)) + sqr(delta_rap(v));
187  }
188 
189  /// R-rap-distance separation dR = sqrt(dphi^2 + drap^2)
190  double delta_r_rap(const FourVector &v) const {
191  return sqrt( delta_r2_rap(v) );
192  }
193 
194  //@}
195 
196 
197  /// @name Operators
198  //@{
199 
200  /// Equality
201  bool operator==(const FourVector& rhs) const {
202  return x() == rhs.x() && y() == rhs.y() && z() == rhs.z() && t() == rhs.t();
203  }
204  /// Inequality
205  bool operator!=(const FourVector& rhs) const { return !(*this == rhs); }
206 
207  /// Arithmetic operator +
208  FourVector operator+ (const FourVector& rhs) const {
209  return FourVector( x() + rhs.x(), y() + rhs.y(), z() + rhs.z(), t() + rhs.t() );
210  }
211  /// Arithmetic operator -
212  FourVector operator- (const FourVector& rhs) const {
213  return FourVector( x() - rhs.x(), y() - rhs.y(), z() - rhs.z(), t() - rhs.t() );
214  }
215  /// Arithmetic operator * by scalar
216  FourVector operator* (const double rhs) const {
217  return FourVector( x()*rhs, y()*rhs, z()*rhs, t()*rhs );
218  }
219  /// Arithmetic operator / by scalar
220  FourVector operator/ (const double rhs) const {
221  return FourVector( x()/rhs, y()/rhs, z()/rhs, t()/rhs );
222  }
223 
224  /// Arithmetic operator +=
225  void operator += (const FourVector& rhs) {
226  setX(x() + rhs.x());
227  setY(y() + rhs.y());
228  setZ(z() + rhs.z());
229  setT(t() + rhs.t());
230  }
231  /// Arithmetic operator -=
232  void operator -= (const FourVector& rhs) {
233  setX(x() - rhs.x());
234  setY(y() - rhs.y());
235  setZ(z() - rhs.z());
236  setT(t() - rhs.t());
237  }
238  /// Arithmetic operator *= by scalar
239  void operator *= (const double rhs) {
240  setX(x()*rhs);
241  setY(y()*rhs);
242  setZ(z()*rhs);
243  setT(t()*rhs);
244  }
245  /// Arithmetic operator /= by scalar
246  void operator /= (const double rhs) {
247  setX(x()/rhs);
248  setY(y()/rhs);
249  setZ(z()/rhs);
250  setT(t()/rhs);
251  }
252 
253  //@}
254 
255 
256  /// Static null FourVector = (0,0,0,0)
257  static const FourVector& ZERO_VECTOR() {
258  static const FourVector v;
259  return v;
260  }
261 
262 
263 private:
264 
265  double m_v1; ///< px or x. Interpretation depends on accessors used
266  double m_v2; ///< py or y. Interpretation depends on accessors used
267  double m_v3; ///< pz or z. Interpretation depends on accessors used
268  double m_v4; ///< e or t. Interpretation depends on accessors used
269 
270 };
271 
272 
273 /// @name Unbound vector comparison functions
274 //@{
275 
276 /// Signed azimuthal angle separation in [-pi, pi] between vecs @c a and @c b
277 inline double delta_phi(const FourVector &a, const FourVector &b) { return b.delta_phi(a); }
278 
279 /// Pseudorapidity separation between vecs @c a and @c b
280 inline double delta_eta(const FourVector &a, const FourVector &b) { return b.delta_eta(a); }
281 
282 /// Rapidity separation between vecs @c a and @c b
283 inline double delta_rap(const FourVector &a, const FourVector &b) { return b.delta_rap(a); }
284 
285 /// R_eta^2-distance separation dR^2 = dphi^2 + deta^2 between vecs @c a and @c b
286 inline double delta_r2_eta(const FourVector &a, const FourVector &b) { return b.delta_r2_eta(a); }
287 
288 /// R_eta-distance separation dR = sqrt(dphi^2 + deta^2) between vecs @c a and @c b
289 inline double delta_r_eta(const FourVector &a, const FourVector &b) { return b.delta_r_eta(a); }
290 
291 /// R_rap^2-distance separation dR^2 = dphi^2 + drap^2 between vecs @c a and @c b
292 inline double delta_r2_rap(const FourVector &a, const FourVector &b) { return b.delta_r2_rap(a); }
293 
294 /// R_rap-distance separation dR = sqrt(dphi^2 + drap^2) between vecs @c a and @c b
295 inline double delta_r_rap(const FourVector &a, const FourVector &b) { return b.delta_r_rap(a); }
296 
297 //@}
298 
299 
300 } // namespace HepMC
301 
302 
303 #endif
double phi() const
Azimuthal angle.
void operator+=(const FourVector &rhs)
Arithmetic operator +=.
bool operator!=(const FourVector &rhs) const
Inequality.
void operator/=(const double rhs)
Arithmetic operator /= by scalar.
double px() const
x-component of momentum
double y() const
y-component of position/displacement
FourVector(double xx, double yy, double zz, double ee)
Sets all FourVector fields.
void operator-=(const FourVector &rhs)
Arithmetic operator -=.
double pt2() const
Squared transverse momentum px^2 + py^2.
void operator*=(const double rhs)
Arithmetic operator *= by scalar.
FourVector operator+(const FourVector &rhs) const
Arithmetic operator +.
double p3mod2() const
Squared magnitude of p3 = (px, py, pz) vector.
void setE(double ee)
Set energy component of momentum.
double pseudoRapidity() const
Same as eta.
double perp() const
Magnitude of (x, y) vector.
double abs_rap() const
Absolute rapidity.
static const FourVector & ZERO_VECTOR()
Static null FourVector = (0,0,0,0)
bool operator==(const FourVector &rhs) const
Equality.
double delta_r2_eta(const FourVector &v) const
R_eta^2-distance separation dR^2 = dphi^2 + deta^2.
double m_v4
e or t. Interpretation depends on accessors used
FourVector(const FourVector &v)
Copy constructor.
bool is_zero() const
Check if the length of this vertex is zero.
double pt() const
Transverse momentum.
double interval() const
Spacetime invariant interval s^2 = t^2 - x^2 - y^2 - z^2.
double delta_rap(const FourVector &v) const
Rapidity separation.
double delta_eta(const FourVector &v) const
Pseudorapidity separation.
FourVector operator*(const double rhs) const
Arithmetic operator * by scalar.
double x() const
x-component of position/displacement
double perp2() const
Squared magnitude of (x, y) vector.
void setT(double tt)
Set time component of position/displacement.
void setX(double xx)
Set x-component of position/displacement.
FourVector()
Default constructor.
double delta_r2_rap(const FourVector &a, const FourVector &b)
R_rap^2-distance separation dR^2 = dphi^2 + drap^2 between vecs a and b.
double delta_phi(const FourVector &a, const FourVector &b)
Signed azimuthal angle separation in [-pi, pi] between vecs a and b.
void setZ(double zz)
Set z-component of position/displacement.
void setPy(double pyy)
Set y-component of momentum.
double theta() const
Polar angle w.r.t. z direction.
NUM sqr(NUM x)
Handy number squaring function.
double delta_r_eta(const FourVector &v) const
R_eta-distance separation dR = sqrt(dphi^2 + deta^2)
double delta_eta(const FourVector &a, const FourVector &b)
Pseudorapidity separation between vecs a and b.
double t() const
Time component of position/displacement.
double delta_r_eta(const FourVector &a, const FourVector &b)
R_eta-distance separation dR = sqrt(dphi^2 + deta^2) between vecs a and b.
double delta_r_rap(const FourVector &a, const FourVector &b)
R_rap-distance separation dR = sqrt(dphi^2 + drap^2) between vecs a and b.
double m() const
Invariant mass. Returns -sqrt(-m) if e^2 - P^2 is negative.
double length2() const
Squared magnitude of (x, y, z) 3-vector.
double m_v2
py or y. Interpretation depends on accessors used
FourVector operator-(const FourVector &rhs) const
Arithmetic operator -.
double py() const
y-component of momentum
double m2() const
Squared invariant mass m^2 = E^2 - px^2 - py^2 - pz^2.
double delta_phi(const FourVector &v) const
Signed azimuthal angle separation in [-pi, pi].
double e() const
Energy component of momentum.
void setPx(double pxx)
Set x-component of momentum.
double delta_r2_rap(const FourVector &v) const
R_rap^2-distance separation dR^2 = dphi^2 + drap^2.
double delta_r_rap(const FourVector &v) const
R-rap-distance separation dR = sqrt(dphi^2 + drap^2)
double length() const
Magnitude of spatial (x, y, z) 3-vector.
FourVector operator/(const double rhs) const
Arithmetic operator / by scalar.
double delta_r2_eta(const FourVector &a, const FourVector &b)
R_eta^2-distance separation dR^2 = dphi^2 + deta^2 between vecs a and b.
double p3mod() const
Magnitude of p3 = (px, py, pz) vector.
double m_v3
pz or z. Interpretation depends on accessors used
double pz() const
z-component of momentum
void setPz(double pzz)
Set z-component of momentum.
Definition of template class SmartPointer.
double abs_eta() const
Absolute pseudorapidity.
void setY(double yy)
Set y-component of position/displacement.
double m_v1
px or x. Interpretation depends on accessors used
double z() const
z-component of position/displacement
double delta_rap(const FourVector &a, const FourVector &b)
Rapidity separation between vecs a and b.