HepMC3 event record library
|
Go to the documentation of this file.
6 #ifndef HEPMC3_FOURVECTOR_H
7 #define HEPMC3_FOURVECTOR_H
16 #define M_PI 3.14159265358979323846264338327950288
54 void set(
double x1,
double x2,
double x3,
double x4) {
64 if (i==0) {
m_v1=
x;
return; }
65 if (i==1) {
m_v2=
x;
return; }
66 if (i==2) {
m_v3=
x;
return; }
67 if (i==3) {
m_v4=
x;
return; }
72 if (i==0)
return m_v1;
73 if (i==1)
return m_v2;
74 if (i==2)
return m_v3;
75 if (i==3)
return m_v4;
81 double x()
const {
return m_v1; }
88 double y()
const {
return m_v2; }
95 double z()
const {
return m_v3; }
110 double px()
const {
return x(); }
117 double py()
const {
return y(); }
124 double pz()
const {
return z(); }
131 double e()
const {
return t(); }
165 double m()
const {
return (
m2() > 0.0) ? std::sqrt(
m2()) : -std::sqrt(-
m2()); }
168 double phi()
const {
return std::atan2(
y(),
x() ); }
174 double rap()
const {
return (
e() == 0.0 ) ? 0.0: (0.5*std::log( (
e() +
pz()) / (
e() -
pz()) )); }
191 bool is_zero()
const {
return x() == 0 &&
y() == 0 &&
z() == 0 &&
t() == 0; }
195 double dphi =
phi() - v.
phi();
196 if (dphi != dphi)
return dphi;
197 while (dphi >=
M_PI) dphi -= 2.*
M_PI;
198 while (dphi < -
M_PI) dphi += 2.*
M_PI;
236 return x() == rhs.
x() &&
y() == rhs.
y() &&
z() == rhs.
z() &&
t() == rhs.
t();
double delta_phi(const FourVector &a, const FourVector &b)
Signed azimuthal angle separation in [-pi, pi] between vecs a and b.
double delta_eta(const FourVector &a, const FourVector &b)
Pseudorapidity separation between vecs a and b.
double m2() const
Squared invariant mass m^2 = E^2 - px^2 - py^2 - pz^2.
void operator-=(const FourVector &rhs)
Arithmetic operator -=.
FourVector operator*(const double rhs) const
Arithmetic operator * by scalar.
double delta_r2_eta(const FourVector &v) const
R_eta^2-distance separation dR^2 = dphi^2 + deta^2.
double p3mod2() const
Squared magnitude of p3 = (px, py, pz) vector.
double delta_rap(const FourVector &v) const
Rapidity separation.
double x() const
x-component of position/displacement
void set_x(double xx)
Set x-component of position/displacement.
double pz() const
z-component of momentum
double theta() const
Polar angle w.r.t. z direction.
void set(double x1, double x2, double x3, double x4)
Set all FourVector fields, in order x,y,z,t.
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.
void operator/=(const double rhs)
Arithmetic operator /= by scalar.
double length2() const
Squared magnitude of (x, y, z) 3-vector.
void set_z(double zz)
Set z-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_v3
pz or z. Interpretation depends on accessors used
double m_v1
px or x. Interpretation depends on accessors used
double delta_eta(const FourVector &v) const
Pseudorapidity separation.
double m_v4
e or t. Interpretation depends on accessors used
#define M_PI
Definition of PI. Needed on some platforms.
void set_pz(double pzz)
Set z-component of momentum.
FourVector(const FourVector &v)
Copy constructor.
void set_py(double pyy)
Set y-component of momentum.
double pt2() const
Squared transverse momentum px^2 + py^2.
double y() const
y-component of position/displacement
double t() const
Time component of position/displacement.
double length() const
Magnitude of spatial (x, y, z) 3-vector.
double p3mod() const
Magnitude of p3 = (px, py, pz) vector.
bool operator==(const FourVector &rhs) const
Equality.
double interval() const
Spacetime invariant interval s^2 = t^2 - x^2 - y^2 - z^2.
static const FourVector & ZERO_VECTOR()
Static null FourVector = (0,0,0,0)
double e() const
Energy component of momentum.
void operator*=(const double rhs)
Arithmetic operator *= by scalar.
Feature< Feature_type > abs(const Feature< Feature_type > &input)
Obtain the absolute value of a Feature. This works as you'd expect. If foo is a valid Feature,...
double delta_r_eta(const FourVector &v) const
R_eta-distance separation dR = sqrt(dphi^2 + deta^2)
double delta_r2_rap(const FourVector &v) const
R_rap^2-distance separation dR^2 = dphi^2 + drap^2.
bool operator!=(const FourVector &rhs) const
Inequality.
double phi() const
Azimuthal angle.
double rap() const
Rapidity.
double m() const
Invariant mass. Returns -sqrt(-m) if e^2 - P^2 is negative.
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.
void set_px(double pxx)
Set x-component of momentum.
double perp() const
Magnitude of (x, y) vector.
double delta_rap(const FourVector &a, const FourVector &b)
Rapidity separation between vecs a and b.
void set_component(const int i, const double x)
set component of position/displacement
double get_component(const int i) const
get component of position/displacement
void set_y(double yy)
Set y-component of position/displacement.
double delta_r_rap(const FourVector &v) const
R-rap-distance separation dR = sqrt(dphi^2 + drap^2)
FourVector operator/(const double rhs) const
Arithmetic operator / by scalar.
double pseudoRapidity() const
double pt() const
Transverse momentum.
FourVector operator-(const FourVector &rhs) const
Arithmetic operator -.
bool is_zero() const
Check if the length of this vertex is zero.
double py() const
y-component of momentum
double abs_rap() const
Absolute rapidity.
FourVector()
Default constructor.
double eta() const
Pseudorapidity.
double delta_phi(const FourVector &v) const
Signed azimuthal angle separation in [-pi, pi].
double px() const
x-component of momentum
FourVector(double xx, double yy, double zz, double ee)
Sets all FourVector fields.
double abs_eta() const
Absolute pseudorapidity.
double perp2() const
Squared magnitude of (x, y) vector.
FourVector operator+(const FourVector &rhs) const
Arithmetic operator +.
double m_v2
py or y. Interpretation depends on accessors used
void operator+=(const FourVector &rhs)
Arithmetic operator +=.
void set_e(double ee)
Set energy component of momentum.
void set_t(double tt)
Set time component of position/displacement.
double z() const
z-component of position/displacement