Warning, file /include/Rivet/Math/Vector3.hh was not indexed
or was modified since last indexation (in which case cross-reference links may be missing, inaccurate or erroneous).
0001 #ifndef RIVET_MATH_VECTOR3
0002 #define RIVET_MATH_VECTOR3
0003
0004 #include "Rivet/Tools/TypeTraits.hh"
0005 #include "Rivet/Math/MathConstants.hh"
0006 #include "Rivet/Math/MathUtils.hh"
0007 #include "Rivet/Math/VectorN.hh"
0008
0009 namespace Rivet {
0010
0011
0012 class Vector3;
0013 typedef Vector3 ThreeVector;
0014 typedef Vector3 V3;
0015 Vector3 multiply(const double, const Vector3&);
0016 Vector3 multiply(const Vector3&, const double);
0017 Vector3 add(const Vector3&, const Vector3&);
0018 Vector3 operator*(const double, const Vector3&);
0019 Vector3 operator*(const Vector3&, const double);
0020 Vector3 operator/(const Vector3&, const double);
0021 Vector3 operator+(const Vector3&, const Vector3&);
0022 Vector3 operator-(const Vector3&, const Vector3&);
0023
0024 class ThreeMomentum;
0025 typedef ThreeMomentum P3;
0026 ThreeMomentum multiply(const double, const ThreeMomentum&);
0027 ThreeMomentum multiply(const ThreeMomentum&, const double);
0028 ThreeMomentum add(const ThreeMomentum&, const ThreeMomentum&);
0029 ThreeMomentum operator*(const double, const ThreeMomentum&);
0030 ThreeMomentum operator*(const ThreeMomentum&, const double);
0031 ThreeMomentum operator/(const ThreeMomentum&, const double);
0032 ThreeMomentum operator+(const ThreeMomentum&, const ThreeMomentum&);
0033 ThreeMomentum operator-(const ThreeMomentum&, const ThreeMomentum&);
0034
0035 class Matrix3;
0036
0037
0038
0039
0040 class Vector3 : public Vector<3> {
0041
0042 friend class Matrix3;
0043 friend Vector3 multiply(const double, const Vector3&);
0044 friend Vector3 multiply(const Vector3&, const double);
0045 friend Vector3 add(const Vector3&, const Vector3&);
0046 friend Vector3 subtract(const Vector3&, const Vector3&);
0047
0048 public:
0049 Vector3() : Vector<3>() { }
0050
0051 template<typename V3TYPE>
0052 Vector3(const V3TYPE& other) {
0053 this->setX(other.x());
0054 this->setY(other.y());
0055 this->setZ(other.z());
0056 }
0057
0058 Vector3(const Vector<3>& other) {
0059 this->setX(other.get(0));
0060 this->setY(other.get(1));
0061 this->setZ(other.get(2));
0062 }
0063
0064 Vector3(double x, double y, double z) {
0065 this->setX(x);
0066 this->setY(y);
0067 this->setZ(z);
0068 }
0069
0070 ~Vector3() { }
0071
0072
0073 public:
0074
0075 static Vector3 mkX() { return Vector3(1,0,0); }
0076 static Vector3 mkY() { return Vector3(0,1,0); }
0077 static Vector3 mkZ() { return Vector3(0,0,1); }
0078
0079
0080 public:
0081
0082 double x() const { return get(0); }
0083 double x2() const { return sqr(x()); }
0084 Vector3& setX(double x) { set(0, x); return *this; }
0085
0086 double y() const { return get(1); }
0087 double y2() const { return sqr(y()); }
0088 Vector3& setY(double y) { set(1, y); return *this; }
0089
0090 double z() const { return get(2); }
0091 double z2() const { return sqr(z()); }
0092 Vector3& setZ(double z) { set(2, z); return *this; }
0093
0094
0095
0096 double dot(const Vector3& v) const {
0097 return _vec.dot(v._vec);
0098 }
0099
0100
0101 Vector3 cross(const Vector3& v) const {
0102 Vector3 result;
0103 result._vec = _vec.cross(v._vec);
0104 return result;
0105 }
0106
0107
0108 double angle(const Vector3& v) const {
0109 const double localDotOther = unit().dot(v.unit());
0110 if (localDotOther > 1.0) return 0.0;
0111 if (localDotOther < -1.0) return M_PI;
0112 return acos(localDotOther);
0113 }
0114
0115
0116
0117 Vector3 unitVec() const {
0118 double md = mod();
0119 if ( fuzzyLessEquals(md, 0.0) ) return Vector3();
0120 else return *this * 1.0/md;
0121 }
0122
0123
0124 Vector3 unit() const {
0125 return unitVec();
0126 }
0127
0128
0129 Vector3 polarVec() const {
0130 Vector3 rtn = *this;
0131 rtn.setZ(0.);
0132 return rtn;
0133 }
0134
0135 Vector3 perpVec() const {
0136 return polarVec();
0137 }
0138
0139 Vector3 rhoVec() const {
0140 return polarVec();
0141 }
0142
0143
0144 double polarRadius2() const {
0145 return x()*x() + y()*y();
0146 }
0147
0148 double perp2() const {
0149 return polarRadius2();
0150 }
0151
0152 double rho2() const {
0153 return polarRadius2();
0154 }
0155
0156
0157 double polarRadius() const {
0158 return sqrt(polarRadius2());
0159 }
0160
0161 double perp() const {
0162 return polarRadius();
0163 }
0164
0165 double rho() const {
0166 return polarRadius();
0167 }
0168
0169
0170
0171
0172
0173 double azimuthalAngle(const PhiMapping mapping = ZERO_2PI) const {
0174
0175
0176 if (x() == 0 && y() == 0) return 0.0;
0177
0178 const double value = atan2( y(), x() );
0179 return mapAngle(value, mapping);
0180 }
0181
0182 double phi(const PhiMapping mapping = ZERO_2PI) const {
0183 return azimuthalAngle(mapping);
0184 }
0185
0186
0187 double tanTheta() const {
0188 return polarRadius()/z();
0189 }
0190
0191
0192 double sinTheta() const {
0193 return sqrt(polarRadius2()/mod2());
0194 }
0195
0196
0197 double cosTheta() const {
0198 return z()/mod2();
0199 }
0200
0201
0202 double polarAngle() const {
0203
0204 const double polarangle = atan2(polarRadius(), z());
0205 return mapAngle0ToPi(polarangle);
0206 }
0207
0208
0209 double theta() const {
0210 return polarAngle();
0211 }
0212
0213
0214
0215
0216
0217
0218
0219
0220
0221 double pseudorapidity() const {
0222 if (mod() == 0.0) return 0.0;
0223 if (mod() == fabs(z()) ) return std::copysign(INF, z());
0224 const double eta = std::log((mod() + fabs(z())) / perp());
0225 return std::copysign(eta, z());
0226 }
0227
0228
0229 double eta() const {
0230 return pseudorapidity();
0231 }
0232
0233
0234 double abseta() const {
0235 return fabs(eta());
0236 }
0237
0238
0239 public:
0240
0241
0242 Vector3& operator *= (const double a) {
0243 _vec = multiply(a, *this)._vec;
0244 return *this;
0245 }
0246
0247
0248 Vector3& operator /= (const double a) {
0249 _vec = multiply(1.0/a, *this)._vec;
0250 return *this;
0251 }
0252
0253
0254 Vector3& operator += (const Vector3& v) {
0255 _vec = add(*this, v)._vec;
0256 return *this;
0257 }
0258
0259
0260 Vector3& operator -= (const Vector3& v) {
0261 _vec = subtract(*this, v)._vec;
0262 return *this;
0263 }
0264
0265
0266 Vector3 operator - () const {
0267 Vector3 rtn;
0268 rtn._vec = -_vec;
0269 return rtn;
0270 }
0271
0272 };
0273
0274
0275
0276
0277 inline double dot(const Vector3& a, const Vector3& b) {
0278 return a.dot(b);
0279 }
0280
0281
0282 inline Vector3 cross(const Vector3& a, const Vector3& b) {
0283 return a.cross(b);
0284 }
0285
0286
0287 inline Vector3 multiply(const double a, const Vector3& v) {
0288 Vector3 result;
0289 result._vec = a * v._vec;
0290 return result;
0291 }
0292
0293
0294 inline Vector3 multiply(const Vector3& v, const double a) {
0295 return multiply(a, v);
0296 }
0297
0298
0299 inline Vector3 operator * (const double a, const Vector3& v) {
0300 return multiply(a, v);
0301 }
0302
0303
0304 inline Vector3 operator * (const Vector3& v, const double a) {
0305 return multiply(a, v);
0306 }
0307
0308
0309 inline Vector3 operator / (const Vector3& v, const double a) {
0310 return multiply(1.0/a, v);
0311 }
0312
0313
0314 inline Vector3 add(const Vector3& a, const Vector3& b) {
0315 Vector3 result;
0316 result._vec = a._vec + b._vec;
0317 return result;
0318 }
0319
0320
0321 inline Vector3 subtract(const Vector3& a, const Vector3& b) {
0322 Vector3 result;
0323 result._vec = a._vec - b._vec;
0324 return result;
0325 }
0326
0327
0328 inline Vector3 operator + (const Vector3& a, const Vector3& b) {
0329 return add(a, b);
0330 }
0331
0332
0333 inline Vector3 operator - (const Vector3& a, const Vector3& b) {
0334 return subtract(a, b);
0335 }
0336
0337
0338
0339
0340 inline double angle(const Vector3& a, const Vector3& b) {
0341 return a.angle(b);
0342 }
0343
0344
0345
0346
0347
0348
0349 class ThreeMomentum : public ThreeVector {
0350 public:
0351 ThreeMomentum() { }
0352
0353 template<typename V3TYPE, typename std::enable_if<HasXYZ<V3TYPE>::value, int>::type DUMMY=0>
0354 ThreeMomentum(const V3TYPE& other) {
0355 this->setPx(other.x());
0356 this->setPy(other.y());
0357 this->setPz(other.z());
0358 }
0359
0360 ThreeMomentum(const Vector<3>& other)
0361 : ThreeVector(other) { }
0362
0363 ThreeMomentum(const double px, const double py, const double pz) {
0364 this->setPx(px);
0365 this->setPy(py);
0366 this->setPz(pz);
0367 }
0368
0369 ~ThreeMomentum() {}
0370
0371 public:
0372
0373
0374
0375
0376
0377
0378 ThreeMomentum& setPx(double px) {
0379 setX(px);
0380 return *this;
0381 }
0382
0383
0384 ThreeMomentum& setPy(double py) {
0385 setY(py);
0386 return *this;
0387 }
0388
0389
0390 ThreeMomentum& setPz(double pz) {
0391 setZ(pz);
0392 return *this;
0393 }
0394
0395
0396
0397
0398
0399
0400
0401
0402 double px() const { return x(); }
0403
0404 double px2() const { return x2(); }
0405
0406
0407 double py() const { return y(); }
0408
0409 double py2() const { return y2(); }
0410
0411
0412 double pz() const { return z(); }
0413
0414 double pz2() const { return z2(); }
0415
0416
0417
0418 double p() const { return mod(); }
0419
0420 double p2() const { return mod2(); }
0421
0422
0423
0424 ThreeMomentum pTvec() const {
0425 return polarVec();
0426 }
0427
0428 ThreeMomentum ptvec() const {
0429 return pTvec();
0430 }
0431
0432
0433 double pT2() const {
0434 return polarRadius2();
0435 }
0436
0437 double pt2() const {
0438 return polarRadius2();
0439 }
0440
0441
0442 double pT() const {
0443 return sqrt(pT2());
0444 }
0445
0446 double pt() const {
0447 return sqrt(pT2());
0448 }
0449
0450
0451
0452
0453
0454
0455
0456
0457
0458
0459
0460 ThreeMomentum& operator *= (double a) {
0461 _vec = multiply(a, *this)._vec;
0462 return *this;
0463 }
0464
0465
0466 ThreeMomentum& operator /= (double a) {
0467 _vec = multiply(1.0/a, *this)._vec;
0468 return *this;
0469 }
0470
0471
0472 ThreeMomentum& operator += (const ThreeMomentum& v) {
0473 _vec = add(*this, v)._vec;
0474 return *this;
0475 }
0476
0477
0478 ThreeMomentum& operator -= (const ThreeMomentum& v) {
0479 _vec = add(*this, -v)._vec;
0480 return *this;
0481 }
0482
0483
0484 ThreeMomentum operator - () const {
0485 ThreeMomentum result;
0486 result._vec = -_vec;
0487 return result;
0488 }
0489
0490
0491
0492
0493
0494
0495
0496
0497 };
0498
0499
0500 inline ThreeMomentum multiply(const double a, const ThreeMomentum& v) {
0501 ThreeMomentum result;
0502 result._vec = a * v._vec;
0503 return result;
0504 }
0505
0506 inline ThreeMomentum multiply(const ThreeMomentum& v, const double a) {
0507 return multiply(a, v);
0508 }
0509
0510 inline ThreeMomentum operator*(const double a, const ThreeMomentum& v) {
0511 return multiply(a, v);
0512 }
0513
0514 inline ThreeMomentum operator*(const ThreeMomentum& v, const double a) {
0515 return multiply(a, v);
0516 }
0517
0518 inline ThreeMomentum operator/(const ThreeMomentum& v, const double a) {
0519 return multiply(1.0/a, v);
0520 }
0521
0522 inline ThreeMomentum add(const ThreeMomentum& a, const ThreeMomentum& b) {
0523 ThreeMomentum result;
0524 result._vec = a._vec + b._vec;
0525 return result;
0526 }
0527
0528 inline ThreeMomentum operator+(const ThreeMomentum& a, const ThreeMomentum& b) {
0529 return add(a, b);
0530 }
0531
0532 inline ThreeMomentum operator-(const ThreeMomentum& a, const ThreeMomentum& b) {
0533 return add(a, -b);
0534 }
0535
0536
0537
0538
0539 inline Vector3 operator+(const ThreeMomentum& a, const Vector3& b) {
0540 return add(static_cast<const Vector3&>(a), b);
0541 }
0542 inline Vector3 operator+(const Vector3& a, const ThreeMomentum& b) {
0543 return add(a, static_cast<const Vector3&>(b));
0544 }
0545
0546 inline Vector3 operator-(const ThreeMomentum& a, const Vector3& b) {
0547 return add(static_cast<const Vector3&>(a), -b);
0548 }
0549 inline Vector3 operator-(const Vector3& a, const ThreeMomentum& b) {
0550 return add(a, -static_cast<const Vector3&>(b));
0551 }
0552
0553
0554
0555
0556
0557
0558
0559
0560
0561
0562 inline double deltaEta(const Vector3& a, const Vector3& b, bool sign=false) {
0563 return deltaEta(a.pseudorapidity(), b.pseudorapidity(), sign);
0564 }
0565
0566
0567 inline double deltaEta(const Vector3& v, double eta2, bool sign=false) {
0568 return deltaEta(v.pseudorapidity(), eta2, sign);
0569 }
0570
0571
0572 inline double deltaEta(double eta1, const Vector3& v, bool sign=false) {
0573 return deltaEta(eta1, v.pseudorapidity(), sign);
0574 }
0575
0576
0577
0578
0579
0580
0581
0582
0583 inline double deltaPhi(const Vector3& a, const Vector3& b, bool sign=false) {
0584 return deltaPhi(a.azimuthalAngle(), b.azimuthalAngle(), sign);
0585 }
0586
0587
0588 inline double deltaPhi(const Vector3& v, double phi2, bool sign=false) {
0589 return deltaPhi(v.azimuthalAngle(), phi2, sign);
0590 }
0591
0592
0593 inline double deltaPhi(double phi1, const Vector3& v, bool sign=false) {
0594 return deltaPhi(phi1, v.azimuthalAngle(), sign);
0595 }
0596
0597
0598
0599
0600
0601
0602
0603
0604 inline double deltaR2(const Vector3& a, const Vector3& b) {
0605 return deltaR2(a.pseudorapidity(), a.azimuthalAngle(),
0606 b.pseudorapidity(), b.azimuthalAngle());
0607 }
0608
0609
0610 inline double deltaR(const Vector3& a, const Vector3& b) {
0611 return sqrt(deltaR2(a,b));
0612 }
0613
0614
0615 inline double deltaR2(const Vector3& v, double eta2, double phi2) {
0616 return deltaR2(v.pseudorapidity(), v.azimuthalAngle(), eta2, phi2);
0617 }
0618
0619
0620 inline double deltaR(const Vector3& v, double eta2, double phi2) {
0621 return sqrt(deltaR2(v, eta2, phi2));
0622 }
0623
0624
0625 inline double deltaR2(double eta1, double phi1, const Vector3& v) {
0626 return deltaR2(eta1, phi1, v.pseudorapidity(), v.azimuthalAngle());
0627 }
0628
0629
0630 inline double deltaR(double eta1, double phi1, const Vector3& v) {
0631 return sqrt(deltaR2(eta1, phi1, v));
0632 }
0633
0634
0635
0636
0637
0638
0639
0640
0641
0642
0643 inline double mT(const Vector3& vis, const Vector3& invis) {
0644
0645 return mT(vis.perp(), invis.perp(), deltaPhi(vis, invis));
0646 }
0647
0648
0649 inline double pT(const Vector3& a, const Vector3& b) {
0650 return (a+b).perp();
0651 }
0652
0653
0654
0655
0656 }
0657
0658 #endif