1 #ifndef RIVET_MATH_VECTOR4 2 #define RIVET_MATH_VECTOR4 4 #include "Rivet/Math/MathConstants.hh" 5 #include "Rivet/Math/MathUtils.hh" 6 #include "Rivet/Math/VectorN.hh" 7 #include "Rivet/Math/Vector3.hh" 14 class LorentzTransform;
15 typedef FourVector Vector4;
16 FourVector
transform(
const LorentzTransform& lt,
const FourVector& v4);
34 this->setT(other.t());
35 this->setX(other.x());
36 this->setY(other.y());
37 this->setZ(other.z());
43 FourVector(
const double t,
const double x,
const double y,
const double z) {
50 virtual ~FourVector() { }
54 double t()
const {
return get(0); }
55 double t2()
const {
return sqr(t()); }
56 FourVector& setT(
const double t) {
set(0, t);
return *
this; }
58 double x()
const {
return get(1); }
59 double x2()
const {
return sqr(x()); }
60 FourVector& setX(
const double x) {
set(1, x);
return *
this; }
62 double y()
const {
return get(2); }
63 double y2()
const {
return sqr(y()); }
64 FourVector& setY(
const double y) {
set(2, y);
return *
this; }
66 double z()
const {
return get(3); }
67 double z2()
const {
return sqr(z()); }
68 FourVector& setZ(
const double z) {
set(3, z);
return *
this; }
70 double invariant()
const {
72 return (t() + z())*(t() - z()) - x()*x() - y()*y();
80 double angle(
const FourVector& v)
const {
163 return Vector3(
get(1),
get(2),
get(3));
174 const double result = t()*v.t() - x()*v.x() - y()*v.y() - z()*v.z();
179 double dot(
const FourVector& v)
const {
190 _vec = multiply(a, *
this)._vec;
196 _vec = multiply(1.0/a, *
this)._vec;
202 _vec = add(*
this, v)._vec;
208 _vec = add(*
this, -v)._vec;
221 FourVector result = -*
this;
222 result.setT(-result.t());
241 result._vec = a * v._vec;
246 return multiply(a, v);
250 return multiply(a, v);
254 return multiply(a, v);
258 return multiply(1.0/a, v);
263 result._vec = a._vec + b._vec;
278 return lv.invariant();
310 template<
typename V4>
312 this->setE(other.t());
313 this->setPx(other.x());
314 this->setPy(other.y());
315 this->setPz(other.z());
321 FourMomentum(
const double E,
const double px,
const double py,
const double pz) {
362 FourMomentum&
setPE(
double px,
double py,
double pz,
double E) {
364 throw std::invalid_argument(
"Negative energy given as argument: " +
to_str(E));
365 setPx(px); setPy(py); setPz(pz); setE(E);
369 FourMomentum&
setXYZE(
double px,
double py,
double pz,
double E) {
370 return setPE(px, py, pz, E);
383 FourMomentum&
setPM(
double px,
double py,
double pz,
double mass) {
385 throw std::invalid_argument(
"Negative mass given as argument: " +
to_str(mass));
386 const double E = sqrt(
sqr(mass) +
sqr(px) +
sqr(py) +
sqr(pz) );
388 return setPE(px, py, pz, E);
391 FourMomentum&
setXYZM(
double px,
double py,
double pz,
double mass) {
392 return setPM(px, py, pz, mass);
402 throw std::invalid_argument(
"Negative mass given as argument");
404 throw std::invalid_argument(
"Negative energy given as argument");
405 const double theta = 2 * atan(exp(-eta));
406 if (theta < 0 || theta > M_PI)
407 throw std::domain_error(
"Polar angle outside 0..pi in calculation");
408 setThetaPhiME(theta, phi, mass, E);
418 throw std::invalid_argument(
"Negative mass given as argument");
420 throw std::invalid_argument(
"Negative transverse momentum given as argument");
421 const double theta = 2 * atan(exp(-eta));
422 if (theta < 0 || theta > M_PI)
423 throw std::domain_error(
"Polar angle outside 0..pi in calculation");
424 const double p = pt / sin(theta);
425 const double E = sqrt(
sqr(p) +
sqr(mass) );
426 setThetaPhiME(theta, phi, mass, E);
440 throw std::invalid_argument(
"Negative mass given as argument");
442 throw std::invalid_argument(
"Negative energy given as argument");
443 const double sqrt_pt2_m2 = E / cosh(y);
444 const double pt = sqrt(
sqr(sqrt_pt2_m2) -
sqr(mass) );
446 throw std::domain_error(
"Negative transverse momentum in calculation");
447 const double pz = sqrt_pt2_m2 * sinh(y);
448 const double px = pt * cos(phi);
449 const double py = pt * sin(phi);
450 setPE(px, py, pz, E);
460 throw std::invalid_argument(
"Negative mass given as argument");
462 throw std::invalid_argument(
"Negative transverse mass given as argument");
463 const double E = sqrt(
sqr(pt) +
sqr(mass) ) * cosh(y);
465 throw std::domain_error(
"Negative energy in calculation");
466 setRapPhiME(y, phi, mass, E);
476 if (theta < 0 || theta > M_PI)
477 throw std::invalid_argument(
"Polar angle outside 0..pi given as argument");
479 throw std::invalid_argument(
"Negative mass given as argument");
481 throw std::invalid_argument(
"Negative energy given as argument");
482 const double p = sqrt(
sqr(E) -
sqr(mass) );
483 const double pz = p * cos(theta);
484 const double pt = p * sin(theta);
486 throw std::invalid_argument(
"Negative transverse momentum in calculation");
487 const double px = pt * cos(phi);
488 const double py = pt * sin(phi);
489 setPE(px, py, pz, E);
499 if (theta < 0 || theta > M_PI)
500 throw std::invalid_argument(
"Polar angle outside 0..pi given as argument");
502 throw std::invalid_argument(
"Negative mass given as argument");
504 throw std::invalid_argument(
"Negative transverse momentum given as argument");
505 const double p = pt / sin(theta);
506 const double px = pt * cos(phi);
507 const double py = pt * sin(phi);
508 const double pz = p * cos(theta);
509 const double E = sqrt(
sqr(p) +
sqr(mass) );
510 setPE(px, py, pz, E);
519 throw std::invalid_argument(
"Negative transverse momentum given as argument");
521 throw std::invalid_argument(
"Negative mass given as argument");
523 throw std::invalid_argument(
"Negative energy given as argument");
524 const double px = pt * cos(phi);
525 const double py = pt * sin(phi);
526 const double pz = sqrt(
sqr(E) -
sqr(mass) -
sqr(pt));
527 setPE(px, py, pz, E);
538 double E()
const {
return t(); }
540 double E2()
const {
return t2(); }
543 double px()
const {
return x(); }
545 double px2()
const {
return x2(); }
548 double py()
const {
return y(); }
550 double py2()
const {
return y2(); }
553 double pz()
const {
return z(); }
555 double pz2()
const {
return z2(); }
568 return sign(mass2()) * sqrt(fabs(mass2()));
593 return 0.5 * std::log( (E() + pz()) / (E() - pz()) );
654 return sqrt(E2()/mass2());
660 return gamma() * p3().
unit();
686 struct byEAscending {
687 bool operator()(
const FourMomentum& left,
const FourMomentum& right)
const{
688 const double pt2left = left.
E();
689 const double pt2right = right.
E();
690 return pt2left < pt2right;
693 bool operator()(
const FourMomentum* left,
const FourMomentum* right)
const{
694 return (*
this)(*left, *right);
701 struct byEDescending {
702 bool operator()(
const FourMomentum& left,
const FourMomentum& right)
const{
703 return byEAscending()(right, left);
706 bool operator()(
const FourMomentum* left,
const FourVector* right)
const{
707 return (*
this)(*left, *right);
722 _vec = multiply(a, *
this)._vec;
728 _vec = multiply(1.0/a, *
this)._vec;
734 _vec = add(*
this, v)._vec;
740 _vec = add(*
this, -v)._vec;
753 FourMomentum result = -*
this;
754 result.
setE(-result.
E());
768 static FourMomentum
mkXYZE(
double px,
double py,
double pz,
double E) {
769 return FourMomentum().
setPE(px, py, pz, E);
773 static FourMomentum
mkXYZM(
double px,
double py,
double pz,
double mass) {
774 return FourMomentum().
setPM(px, py, pz, mass);
779 return FourMomentum().
setEtaPhiME(eta, phi, mass, E);
788 static FourMomentum
mkRapPhiME(
double y,
double phi,
double mass,
double E) {
789 return FourMomentum().
setRapPhiME(y, phi, mass, E);
808 static FourMomentum
mkPtPhiME(
double pt,
double phi,
double mass,
double E) {
809 return FourMomentum().
setPtPhiME(pt, phi, mass, E);
821 result._vec = a * v._vec;
826 return multiply(a, v);
830 return multiply(a, v);
834 return multiply(a, v);
838 return multiply(1.0/a, v);
843 result._vec = a._vec + b._vec;
873 case PSEUDORAPIDITY :
880 string err =
"deltaR with scheme RAPIDITY can only be called with FourMomentum objects, not FourVectors";
881 throw std::runtime_error(err);
883 return deltaR2(*ma, *mb, scheme);
886 throw std::runtime_error(
"The specified deltaR scheme is not yet implemented");
900 return sqrt(
deltaR2(a, b, scheme));
912 double eta2,
double phi2,
915 case PSEUDORAPIDITY :
921 string err =
"deltaR with scheme RAPIDITY can only be called with FourMomentum objects, not FourVectors";
922 throw std::runtime_error(err);
924 return deltaR2(*mv, eta2, phi2, scheme);
927 throw std::runtime_error(
"The specified deltaR scheme is not yet implemented");
938 double eta2,
double phi2,
940 return sqrt(
deltaR2(v, eta2, phi2, scheme));
950 inline double deltaR2(
double eta1,
double phi1,
954 case PSEUDORAPIDITY :
960 string err =
"deltaR with scheme RAPIDITY can only be called with FourMomentum objects, not FourVectors";
961 throw std::runtime_error(err);
963 return deltaR2(eta1, phi1, *mv, scheme);
966 throw std::runtime_error(
"The specified deltaR scheme is not yet implemented");
976 inline double deltaR(
double eta1,
double phi1,
979 return sqrt(
deltaR2(eta1, phi1, v, scheme));
997 throw std::runtime_error(
"The specified deltaR scheme is not yet implemented");
1009 return sqrt(
deltaR2(a, b, scheme));
1019 double eta2,
double phi2,
1022 case PSEUDORAPIDITY:
1027 throw std::runtime_error(
"The specified deltaR scheme is not yet implemented");
1037 double eta2,
double phi2,
1039 return sqrt(
deltaR2(v, eta2, phi2, scheme));
1052 case PSEUDORAPIDITY:
1057 throw std::runtime_error(
"The specified deltaR scheme is not yet implemented");
1066 inline double deltaR(
double eta1,
double phi1,
1069 return sqrt(
deltaR2(eta1, phi1, v, scheme));
1081 case PSEUDORAPIDITY:
1086 throw std::runtime_error(
"The specified deltaR scheme is not yet implemented");
1097 return sqrt(
deltaR2(a, b, scheme));
1118 return deltaR(b, a, scheme);
1343 return a.
pt() > b.
pt();
1347 return a.
pt() < b.
pt();
1361 return a.
Et() > b.
Et();
1365 return a.
Et() < b.
Et();
1370 return a.
E() > b.
E();
1374 return a.
E() < b.
E();
1388 return a.
eta() < b.
eta();
1398 return fabs(a.
eta()) < fabs(b.
eta());
1403 return fabs(a.
eta()) > fabs(b.
eta());
1430 template<
typename MOMS,
typename CMP>
1432 std::sort(pbs.begin(), pbs.end(),
cmp);
1436 template<
typename MOMS,
typename CMP>
1439 std::sort(rtn.begin(), rtn.end(),
cmp);
1444 template<
typename MOMS>
1449 template<
typename MOMS>
1455 template<
typename MOMS>
1460 template<
typename MOMS>
1466 template<
typename MOMS>
1471 template<
typename MOMS>
1490 return mT(vis.
p3(), invis.
p3());
1495 return mT(vis.
p3(), invis);
1500 return mT(vis, invis.
p3());
1514 std::ostringstream out;
1515 out <<
"(" << (fabs(lv.t()) < 1E-30 ? 0.0 : lv.t())
1516 <<
"; " << (fabs(lv.x()) < 1E-30 ? 0.0 : lv.x())
1517 <<
", " << (fabs(lv.y()) < 1E-30 ? 0.0 : lv.y())
1518 <<
", " << (fabs(lv.z()) < 1E-30 ? 0.0 : lv.z())
1534 typedef std::vector<FourMomentum> FourMomenta;
FourMomentum & setPy(double py)
Set y-component of momentum .
Definition: Vector4.hh:349
FourVector operator-() const
Multiply all components (space and time) by -1.
Definition: Vector4.hh:213
bool cmpMomByAscE(const FourMomentum &a, const FourMomentum &b)
Comparison to give a sorting by increasing energy.
Definition: Vector4.hh:1373
Definition: MC_Cent_pPb.hh:10
double abseta() const
Get the directly (alias).
Definition: Vector4.hh:159
double py2() const
Get y-squared .
Definition: Vector4.hh:550
FourMomentum & setXYZM(double px, double py, double pz, double mass)
Alias for setPM.
Definition: Vector4.hh:391
Vector3 rhoVec() const
Synonym for polarVec.
Definition: Vector3.hh:121
double pT() const
Calculate the transverse momentum .
Definition: Vector4.hh:628
FourMomentum & operator+=(const FourMomentum &v)
Add to this 4-vector. NB time as well as space components are added.
Definition: Vector4.hh:733
Vector3 betaVec() const
Definition: Vector4.hh:671
double perp() const
Synonym for polarRadius.
Definition: Vector4.hh:108
double rapidity(double E, double pz)
Calculate a rapidity value from the supplied energy E and longitudinal momentum pz.
Definition: MathUtils.hh:617
double px() const
Get x-component of momentum .
Definition: Vector4.hh:543
double mod2() const
Calculate the modulus-squared of a vector. .
Definition: VectorN.hh:84
double eta() const
Synonym for pseudorapidity.
Definition: Vector4.hh:152
static FourMomentum mkRapPhiMPt(double y, double phi, double mass, double pt)
Make a vector from (y,phi,pT) coordinates and the mass.
Definition: Vector4.hh:793
FourMomentum & setPz(double pz)
Set z-component of momentum .
Definition: Vector4.hh:355
FourVector & operator*=(double a)
Multiply by a scalar.
Definition: Vector4.hh:189
double polarRadius() const
Polar radius.
Definition: Vector3.hh:139
double absrap() const
Absolute rapidity.
Definition: Vector4.hh:605
bool cmpMomByMass(const FourMomentum &a, const FourMomentum &b)
Comparison to give a sorting by decreasing mass.
Definition: Vector4.hh:1378
double pt2() const
Calculate the squared transverse momentum .
Definition: Vector4.hh:623
bool cmpMomByDescAbsRap(const FourMomentum &a, const FourMomentum &b)
Comparison to give a sorting by decreasing absolute rapidity.
Definition: Vector4.hh:1422
MOMS & sortByPt(MOMS &pbs)
Sort a container of momenta by pT (decreasing) and return by reference for non-const inputs...
Definition: Vector4.hh:1445
Vector3 polarVec() const
Projection of 3-vector on to the plane.
Definition: Vector4.hh:117
double azimuthalAngle(const PhiMapping mapping=ZERO_2PI) const
Angle subtended by the 3-vector's projection in x-y and the x-axis.
Definition: Vector4.hh:130
double Et() const
Calculate the transverse energy .
Definition: Vector4.hh:641
bool cmpMomByEt(const FourMomentum &a, const FourMomentum &b)
Comparison to give a sorting by decreasing transverse energy.
Definition: Vector4.hh:1360
Vector3 ptvec() const
Synonym for pTvec.
Definition: Vector4.hh:614
double rho() const
Synonym for polarRadius.
Definition: Vector3.hh:147
Vector3 gammaVec() const
Definition: Vector4.hh:659
double absrapidity() const
Absolute rapidity.
Definition: Vector4.hh:601
double p2() const
Get the modulus-squared of the 3-momentum.
Definition: Vector4.hh:586
FourMomentum & operator*=(double a)
Multiply by a scalar.
Definition: Vector4.hh:721
double perp() const
Synonym for polarRadius.
Definition: Vector3.hh:143
bool cmpMomByEta(const FourMomentum &a, const FourMomentum &b)
Comparison to give a sorting by increasing eta (pseudorapidity)
Definition: Vector4.hh:1387
static FourMomentum mkThetaPhiME(double theta, double phi, double mass, double E)
Make a vector from (theta,phi,energy) coordinates and the mass.
Definition: Vector4.hh:798
FourMomentum & setE(double E)
Set energy (time component of momentum).
Definition: Vector4.hh:337
double contract(const FourVector &v) const
Contract two 4-vectors, with metric signature (+ - - -).
Definition: Vector4.hh:173
bool cmpMomByPt(const FourMomentum &a, const FourMomentum &b)
Comparison to give a sorting by decreasing pT.
Definition: Vector4.hh:1342
static FourMomentum mkEtaPhiMPt(double eta, double phi, double mass, double pt)
Make a vector from (eta,phi,pT) coordinates and the mass.
Definition: Vector4.hh:783
static FourMomentum mkEtaPhiME(double eta, double phi, double mass, double E)
Make a vector from (eta,phi,energy) coordinates and the mass.
Definition: Vector4.hh:778
double mass2() const
Get the squared mass (the Lorentz self-invariant).
Definition: Vector4.hh:572
double abspseudorapidity() const
Get the directly.
Definition: Vector4.hh:157
bool cmpMomByDescEta(const FourMomentum &a, const FourMomentum &b)
Comparison to give a sorting by decreasing eta (pseudorapidity)
Definition: Vector4.hh:1392
double rapidity() const
Calculate the rapidity.
Definition: Vector4.hh:592
double pz2() const
Get z-squared .
Definition: Vector4.hh:555
FourMomentum & operator-=(const FourMomentum &v)
Subtract from this 4-vector. NB time as well as space components are subtracted.
Definition: Vector4.hh:739
Vector3 p3() const
Get 3-momentum part, .
Definition: Vector4.hh:578
double pseudorapidity() const
Pseudorapidity (defined purely by the 3-vector components)
Definition: Vector4.hh:148
bool cmpMomByE(const FourMomentum &a, const FourMomentum &b)
Comparison to give a sorting by decreasing energy.
Definition: Vector4.hh:1369
double polarRadius() const
Magnitude of projection of 3-vector on to the plane.
Definition: Vector4.hh:104
double E2() const
Get energy-squared .
Definition: Vector4.hh:540
double deltaEta(double eta1, double eta2, bool sign=false)
Definition: MathUtils.hh:590
FourMomentum & setThetaPhiME(double theta, double phi, double mass, double E)
Definition: Vector4.hh:475
bool cmpMomByAscPt(const FourMomentum &a, const FourMomentum &b)
Comparison to give a sorting by increasing pT.
Definition: Vector4.hh:1346
double angle(const Vector3 &v) const
Angle in radians to another vector.
Definition: Vector3.hh:89
FourMomentum & setEtaPhiMPt(double eta, double phi, double mass, double pt)
Definition: Vector4.hh:416
MOMS & sortBy(MOMS &pbs, const CMP &cmp)
Sort a container of momenta by cmp and return by reference for non-const inputs.
Definition: Vector4.hh:1431
bool cmpMomByP(const FourMomentum &a, const FourMomentum &b)
Comparison to give a sorting by decreasing 3-momentum magnitude |p|.
Definition: Vector4.hh:1351
Vector3 polarVec() const
Polar projection of this vector into the x-y plane.
Definition: Vector3.hh:111
FourMomentum & setPtPhiME(double pt, double phi, double mass, double E)
Definition: Vector4.hh:517
FourMomentum & setThetaPhiMPt(double theta, double phi, double mass, double pt)
Definition: Vector4.hh:498
static FourMomentum mkPtPhiME(double pt, double phi, double mass, double E)
Make a vector from (pT,phi,energy) coordinates and the mass.
Definition: Vector4.hh:808
string to_str(const T &x)
Convert any object to a string.
Definition: Utils.hh:75
FourMomentum & setPx(double px)
Set x-component of momentum .
Definition: Vector4.hh:343
double polarAngle() const
Angle subtended by the vector and the z-axis.
Definition: Vector3.hh:174
PhiMapping
Enum for range of to be mapped into.
Definition: MathConstants.hh:49
double gamma() const
Definition: Vector4.hh:653
double dot(const FourVector &v) const
Contract two 4-vectors, with metric signature (+ - - -).
Definition: Vector4.hh:179
double pseudorapidity() const
Definition: Vector3.hh:190
double deltaPhi(double phi1, double phi2, bool sign=false)
Calculate the difference between two angles in radians.
Definition: MathUtils.hh:582
double angle(const Vector3 &v3) const
Angle between this vector and another (3-vector)
Definition: Vector4.hh:84
double rho2() const
Synonym for polarRadius2.
Definition: Vector3.hh:134
double mass() const
Get the mass (the Lorentz self-invariant).
Definition: Vector4.hh:561
double beta() const
Definition: Vector4.hh:665
double rho2() const
Synonym for polarRadius2.
Definition: Vector4.hh:99
FourMomentum operator-() const
Multiply all components (time and space) by -1.
Definition: Vector4.hh:745
bool cmpMomByDescAbsEta(const FourMomentum &a, const FourMomentum &b)
Comparison to give a sorting by increasing absolute eta (pseudorapidity)
Definition: Vector4.hh:1402
bool cmpMomByAscMass(const FourMomentum &a, const FourMomentum &b)
Comparison to give a sorting by increasing mass.
Definition: Vector4.hh:1382
Specialisation of VectorN to a general (non-momentum) Lorentz 4-vector.
Definition: Vector4.hh:22
double mT(double pT1, double pT2, double dphi)
Definition: MathUtils.hh:634
double p() const
Get the modulus of the 3-momentum.
Definition: Vector4.hh:581
MOMS & sortByE(MOMS &pbs)
Sort a container of momenta by E (decreasing) and return by reference for non-const inputs...
Definition: Vector4.hh:1456
double invariant(const FourVector &lv)
Definition: Vector4.hh:277
static FourMomentum mkXYZE(double px, double py, double pz, double E)
Make a vector from (px,py,pz,E) coordinates.
Definition: Vector4.hh:768
std::vector< FourVector > FourVectors
Definition: Vector4.hh:1533
double deltaRap(double y1, double y2, bool sign=false)
Definition: MathUtils.hh:598
A minimal base class for -dimensional vectors.
Definition: VectorN.hh:13
double deltaR2(double rap1, double phi1, double rap2, double phi2)
Definition: MathUtils.hh:605
double polarRadius2() const
Mod-square of the projection of the 3-vector on to the plane This is a more efficient function than ...
Definition: Vector4.hh:91
double deltaR(double rap1, double phi1, double rap2, double phi2)
Definition: MathUtils.hh:612
double E() const
Get energy (time component of momentum).
Definition: Vector4.hh:538
const CONTAINER2 & transform(const CONTAINER1 &in, CONTAINER2 &out, const FN &f)
A single-container-arg version of std::transform, aka map.
Definition: Utils.hh:385
Vector3 perpVec() const
Synonym for polarVec.
Definition: Vector4.hh:121
bool cmpMomByAscP(const FourMomentum &a, const FourMomentum &b)
Comparison to give a sorting by increasing 3-momentum magnitude |p|.
Definition: Vector4.hh:1355
RapScheme
Enum for rapidity variable to be used in calculating , applying rapidity cuts, etc.
Definition: MathConstants.hh:46
FourMomentum & operator/=(double a)
Divide by a scalar.
Definition: Vector4.hh:727
Vector3 perpVec() const
Synonym for polarVec.
Definition: Vector3.hh:117
double operator*(const FourVector &v) const
Contract two 4-vectors, with metric signature (+ - - -).
Definition: Vector4.hh:184
FourMomentum & setRapPhiMPt(double y, double phi, double mass, double pt)
Definition: Vector4.hh:458
double perp2() const
Synonym for polarRadius2.
Definition: Vector4.hh:95
double pt() const
Calculate the transverse momentum .
Definition: Vector4.hh:632
Vector3 vector3() const
Get the spatial part of the 4-vector as a 3-vector.
Definition: Vector4.hh:162
Vector3 unit() const
Synonym for unitVec.
Definition: Vector3.hh:105
FourMomentum & setXYZE(double px, double py, double pz, double E)
Alias for setPE.
Definition: Vector4.hh:369
std::string toString(const AnalysisInfo &ai)
String representation.
static FourMomentum mkThetaPhiMPt(double theta, double phi, double mass, double pt)
Make a vector from (theta,phi,pT) coordinates and the mass.
Definition: Vector4.hh:803
FourMomentum & setRapPhiME(double y, double phi, double mass, double E)
Definition: Vector4.hh:438
double angle(const FourVector &v) const
Angle between this vector and another.
Definition: Vector4.hh:80
double theta() const
Synonym for polarAngle.
Definition: Vector4.hh:143
std::enable_if< std::is_floating_point< NUM >::value, bool >::type isZero(NUM val, double tolerance=1e-8)
Compare a number to zero.
Definition: MathUtils.hh:21
double polarAngle() const
Angle subtended by the 3-vector and the z-axis.
Definition: Vector4.hh:139
MOMS & sortByEt(MOMS &pbs)
Sort a container of momenta by Et (decreasing) and return by reference for non-const inputs...
Definition: Vector4.hh:1467
Vector3 pTvec() const
Calculate the transverse momentum vector .
Definition: Vector4.hh:610
FourMomentum & setEtaPhiME(double eta, double phi, double mass, double E)
Definition: Vector4.hh:400
double pz() const
Get z-component of momentum .
Definition: Vector4.hh:553
double phi(const PhiMapping mapping=ZERO_2PI) const
Synonym for azimuthalAngle.
Definition: Vector4.hh:134
double py() const
Get y-component of momentum .
Definition: Vector4.hh:548
double azimuthalAngle(const PhiMapping mapping=ZERO_2PI) const
Definition: Vector3.hh:155
FourVector & operator-=(const FourVector &v)
Subtract from this 4-vector. NB time as well as space components are subtracted.
Definition: Vector4.hh:207
double px2() const
Get x-squared .
Definition: Vector4.hh:545
FourMomentum & setPE(double px, double py, double pz, double E)
Set the p coordinates and energy simultaneously.
Definition: Vector4.hh:362
double theta() const
Synonym for polarAngle.
Definition: Vector3.hh:181
double phi(const PhiMapping mapping=ZERO_2PI) const
Synonym for azimuthalAngle.
Definition: Vector3.hh:164
bool cmpMomByRap(const FourMomentum &a, const FourMomentum &b)
Comparison to give a sorting by increasing rapidity.
Definition: Vector4.hh:1407
Three-dimensional specialisation of Vector.
Definition: Vector3.hh:26
bool cmpMomByAscEt(const FourMomentum &a, const FourMomentum &b)
Comparison to give a sorting by increasing transverse energy.
Definition: Vector4.hh:1364
FourVector & operator+=(const FourVector &v)
Add to this 4-vector.
Definition: Vector4.hh:201
Vector3 rhoVec() const
Synonym for polarVec.
Definition: Vector4.hh:125
std::ostream & operator<<(std::ostream &os, const AnalysisInfo &ai)
Stream an AnalysisInfo as a text description.
Definition: AnalysisInfo.hh:365
FourVector & operator/=(double a)
Divide by a scalar.
Definition: Vector4.hh:195
Specialized version of the FourVector with momentum/energy functionality.
Definition: Vector4.hh:301
double eta() const
Synonym for pseudorapidity.
Definition: Vector3.hh:200
double Et2() const
Calculate the transverse energy .
Definition: Vector4.hh:637
static FourMomentum mkRapPhiME(double y, double phi, double mass, double E)
Make a vector from (y,phi,energy) coordinates and the mass.
Definition: Vector4.hh:788
static FourMomentum mkXYZM(double px, double py, double pz, double mass)
Make a vector from (px,py,pz) coordinates and the mass.
Definition: Vector4.hh:773
FourMomentum reverse() const
Multiply space components only by -1.
Definition: Vector4.hh:752
FourMomentum & setPM(double px, double py, double pz, double mass)
Set the p coordinates and mass simultaneously.
Definition: Vector4.hh:383
std::enable_if< std::is_arithmetic< NUM >::value, NUM >::type sqr(NUM a)
Named number-type squaring operation.
Definition: MathUtils.hh:198
double rap() const
Alias for rapidity.
Definition: Vector4.hh:596
double polarRadius2() const
Square of the polar radius (.
Definition: Vector3.hh:126
double mod() const
Calculate the modulus of a vector. .
Definition: VectorN.hh:95
double rho() const
Synonym for polarRadius.
Definition: Vector4.hh:112
double perp2() const
Synonym for polarRadius2.
Definition: Vector3.hh:130
bool cmpMomByAbsEta(const FourMomentum &a, const FourMomentum &b)
Comparison to give a sorting by increasing absolute eta (pseudorapidity)
Definition: Vector4.hh:1397
FourVector reverse() const
Multiply space components only by -1.
Definition: Vector4.hh:220
std::enable_if< std::is_arithmetic< NUM >::value, int >::type sign(NUM val)
Find the sign of a number.
Definition: MathUtils.hh:245
bool cmpMomByDescRap(const FourMomentum &a, const FourMomentum &b)
Comparison to give a sorting by decreasing rapidity.
Definition: Vector4.hh:1412
double pT2() const
Calculate the squared transverse momentum .
Definition: Vector4.hh:619
Cmp< T > cmp(const T &t1, const T &t2)
Global helper function for easy creation of Cmp objects.
Definition: Cmp.hh:255
bool cmpMomByAbsRap(const FourMomentum &a, const FourMomentum &b)
Comparison to give a sorting by increasing absolute rapidity.
Definition: Vector4.hh:1417