Belle II Software  light-2205-abys
B2Vector3.h
1 /**************************************************************************
2  * basf2 (Belle II Analysis Software Framework) *
3  * Author: The Belle II Collaboration *
4  * *
5  * See git log for contributors and copyright holders. *
6  * This file is licensed under LGPL-3.0, see LICENSE.md. *
7  **************************************************************************/
8 
9 #pragma once
10 
11 
12 #include <framework/logging/Logger.h>
13 
14 #include <TVector3.h>
15 #include <Math/Vector3D.h>
16 #include <string>
17 #include <iostream> // std::cout, std::fixed
18 #include <iomanip> // std::setprecision
19 #include <typeinfo>
20 #include <cmath>
21 
22 
23 namespace Belle2 {
41  template<typename DataType>
42  class B2Vector3 {
43  protected:
45  static_assert(std::is_floating_point<DataType>::value, "B2Vector3 only works with floating point types");
47  DataType m_coordinates[3];
48  public:
50  typedef DataType value_type;
51 
53  B2Vector3(void) : m_coordinates {static_cast<DataType>(0), static_cast<DataType>(0), static_cast<DataType>(0)} {};
55  B2Vector3(const DataType xVal, const DataType yVal, const DataType zVal): m_coordinates {xVal, yVal, zVal} {};
57  explicit B2Vector3(const DataType(& coords)[3]): m_coordinates {coords[0], coords[1], coords[2]} {};
59  explicit B2Vector3(const DataType(* coords)[3]): m_coordinates {(*coords)[0], (*coords)[1], (*coords)[2]} {};
61  // cppcheck-suppress noExplicitConstructor
62  B2Vector3(const TVector3& tVec3): m_coordinates {static_cast<DataType>(tVec3.X()), static_cast<DataType>(tVec3.Y()), static_cast<DataType>(tVec3.Z())} {};
64  // cppcheck-suppress noExplicitConstructor
65  B2Vector3(const TVector3* tVec3): m_coordinates {static_cast<DataType>(tVec3->X()), static_cast<DataType>(tVec3->Y()), static_cast<DataType>(tVec3->Z())} {};
67  explicit B2Vector3(const B2Vector3<DataType>& b2Vec3): m_coordinates {b2Vec3.X(), b2Vec3.Y(), b2Vec3.Z()} {};
69  explicit B2Vector3(const B2Vector3<DataType>* b2Vec3): m_coordinates {b2Vec3->X(), b2Vec3->Y(), b2Vec3->Z()} {};
71  // cppcheck-suppress noExplicitConstructor
72  template <typename OtherType> B2Vector3(const B2Vector3<OtherType>& b2Vec3):
73  m_coordinates {static_cast<DataType>(b2Vec3.X()), static_cast<DataType>(b2Vec3.Y()), static_cast<DataType>(b2Vec3.Z())} {};
75  template <typename OtherType> explicit B2Vector3(const B2Vector3<OtherType>* b2Vec3):
76  m_coordinates {static_cast<DataType>(b2Vec3->X()), static_cast<DataType>(b2Vec3->Y()), static_cast<DataType>(b2Vec3->Z())} {};
78  // cppcheck-suppress noExplicitConstructor
79  B2Vector3(const ROOT::Math::XYZVector& xyzVec): m_coordinates {static_cast<DataType>(xyzVec.X()), static_cast<DataType>(xyzVec.Y()), static_cast<DataType>(xyzVec.Z())} {};
81  // cppcheck-suppress noExplicitConstructor
82  B2Vector3(const ROOT::Math::XYZVector* xyzVec): m_coordinates {static_cast<DataType>(xyzVec->X()), static_cast<DataType>(xyzVec->Y()), static_cast<DataType>(xyzVec->Z())} {};
83 
85  DataType operator()(unsigned i) const { return m_coordinates[i]; }
87  DataType operator[](unsigned i) const { return m_coordinates[i]; }
89  DataType& operator()(unsigned i) { return m_coordinates[i]; }
91  DataType& operator[](unsigned i) { return m_coordinates[i]; }
92 
96  B2Vector3<DataType>& operator = (const TVector3& b);
98  B2Vector3<DataType>& operator = (const ROOT::Math::XYZVector& b);
99 
101  operator TVector3() const { return GetTVector3(); }
102 
104  bool operator == (const B2Vector3<DataType>& b) const { return X() == b.X() && Y() == b.Y() && Z() == b.Z(); }
106  bool operator == (const TVector3& b) const { return X() == b.X() && Y() == b.Y() && Z() == b.Z(); }
108  bool operator != (const B2Vector3<DataType>& b) const { return !(*this == b); }
110  bool operator != (const TVector3& b) const { return !(*this == b); }
111 
119  B2Vector3<DataType> operator - () const { return B2Vector3<DataType>(-X(), -Y(), -Z()); }
122  {
123  return B2Vector3<DataType>(X() + b.X(), Y() + b.Y(), Z() + b.Z());
124  }
127  {
128  return B2Vector3<DataType>(X() - b.X(), Y() - b.Y(), Z() - b.Z());
129  }
132  {
133  return B2Vector3<DataType>(a * X(), a * Y(), a * Z());
134  }
137  {
138  return B2Vector3<DataType>(X() / a, Y() / a, Z() / a);
139  }
141  DataType operator * (const B2Vector3<DataType>& b) const { return Dot(b); }
142 
143 
145  DataType Phi() const { return X() == 0 && Y() == 0 ? 0 : atan2(Y(), X()); }
147  DataType Theta() const { return X() == 0 && Y() == 0 && Z() == 0 ? 0 : atan2(Perp(), Z()); }
149  DataType CosTheta() const { const double pTot = Mag(); return pTot == 0 ? 1 : Z() / pTot; }
151  DataType Mag2() const { return X() * X() + Y() * Y() + Z() * Z(); }
153  DataType Mag() const { return std::hypot((double)Perp(), (double)Z()); }
154 
156  void SetPhi(DataType phi)
157  {
158  const double perp = Perp();
159  SetX(perp * cos((double)phi));
160  SetY(perp * sin((double)phi));
161  }
162 
164  void SetTheta(DataType theta)
165  {
166  const double ma = Mag();
167  const double ph = Phi();
168  const double ctheta = std::cos((double) theta);
169  const double stheta = std::sin((double) theta);
170  SetX(ma * stheta * std::cos(ph));
171  SetY(ma * stheta * std::cos(ph));
172  SetZ(ma * ctheta);
173  }
174 
176  void SetMag(DataType mag)
177  {
178  double factor = Mag();
179  if (factor == 0) {
180  B2WARNING(name() << "::SetMag: zero vector can't be stretched");
181  } else {
182  factor = mag / factor;
183  SetX(X()*factor);
184  SetY(Y()*factor);
185  SetZ(Z()*factor);
186  }
187  }
188 
190  DataType Perp2() const { return X() * X() + Y() * Y(); }
192  DataType Pt() const { return Perp(); }
194  DataType Perp() const { return std::hypot((double)X(), (double)Y()); }
195 
197  void SetPerp(DataType r)
198  {
199  const double p = Perp();
200  if (p != 0.0) {
201  m_coordinates[0] *= r / p;
202  m_coordinates[1] *= r / p;
203  }
204  }
205 
207  DataType Perp2(const B2Vector3<DataType>& axis) const
208  {
209  const double tot = axis.Mag2();
210  const double ss = Dot(axis);
211  double per = Mag2();
212  if (tot > 0.0) per -= ss * ss / tot;
213  if (per < 0) per = 0;
214  return per;
215  }
216 
218  DataType Pt(const B2Vector3<DataType>& axis) const { return Perp(axis); }
220  DataType Perp(const B2Vector3<DataType>& axis) const { return std::sqrt(Perp2(axis)); }
222  DataType DeltaPhi(const B2Vector3<DataType>& v) const { return Mpi_pi(Phi() - v.Phi()); }
223 
224 
226  static DataType Mpi_pi(DataType angle)
227  {
228  if (std::isnan(angle)) {
229  B2ERROR(name() << "::Mpi_pi: function called with NaN");
230  return angle;
231  }
232  angle = std::remainder(angle, 2 * M_PI);
233  //for compatibility with ROOT we flip the sign for exactly pi
234  if (angle == M_PI) angle = -M_PI;
235  return angle;
236  }
237 
239  DataType DeltaR(const B2Vector3<DataType>& v) const
240  {
241  const double deta = Eta() - v.Eta();
242  const double dphi = DeltaPhi(v);
243  return std::hypot(deta, dphi);
244  }
245 
247  DataType DrEtaPhi(const B2Vector3<DataType>& v) const
248  {
249  return DeltaR(v);
250  }
251 
253  void SetMagThetaPhi(DataType mag, DataType theta, DataType phi)
254  {
255  const double amag = std::abs(mag);
256  const double sinTheta = std::sin((double)theta);
257  m_coordinates[0] = amag * sinTheta * std::cos((double)phi);
258  m_coordinates[1] = amag * sinTheta * std::sin((double)phi);
259  m_coordinates[2] = amag * std::cos((double)theta);
260  }
261 
264  {
265  const double tot = Mag2();
266  B2Vector3<DataType> p(X(), Y(), Z());
267  return tot > 0.0 ? p *= (1.0 / std::sqrt(tot)) : p;
268  }
269 
272  {
273  const double xVal = std::abs((double)X());
274  const double yVal = std::abs((double)Y());
275  const double zVal = std::abs((double)Z());
276  if (xVal < yVal) {
277  return xVal < zVal ? B2Vector3<DataType>(0, Z(), -Y()) : B2Vector3<DataType>(Y(), -X(), 0);
278  } else {
279  return yVal < zVal ? B2Vector3<DataType>(-Z(), 0, X()) : B2Vector3<DataType>(Y(), -X(), 0);
280  }
281  }
282 
284  DataType Dot(const B2Vector3<DataType>& p) const
285  {
286  return X() * p.X() + Y() * p.Y() + Z() * p.Z();
287  }
288 
291  {
292  return B2Vector3<DataType>(Y() * p.Z() - p.Y() * Z(), Z() * p.X() - p.Z() * X(), X() * p.Y() - p.X() * Y());
293  }
294 
296  DataType Angle(const B2Vector3<DataType>& q) const
297  {
298  const double ptot2 = Mag2() * q.Mag2();
299  if (ptot2 <= 0) {
300  return 0.0;
301  } else {
302  double arg = Dot(q) / std::sqrt(ptot2);
303  if (arg > 1.0) arg = 1.0;
304  if (arg < -1.0) arg = -1.0;
305  return std::acos(arg);
306  }
307  }
308 
313  DataType PseudoRapidity() const
314  {
315  const double cosTheta = CosTheta();
316  if (std::abs(cosTheta) < 1) return -0.5 * std::log((1.0 - cosTheta) / (1.0 + cosTheta));
317  if (Z() == 0) return 0;
318  //B2WARNING(name() << "::PseudoRapidity: transverse momentum = 0! return +/- 10e10");
319  if (Z() > 0) return 10e10;
320  else return -10e10;
321  }
322 
323 
325  DataType Eta() const { return PseudoRapidity(); }
326 
327 
329  void RotateX(DataType angle)
330  {
331  //rotate vector around X
332  const double s = std::sin((double)angle);
333  const double c = std::cos((double)angle);
334  const double yOld = Y();
335  m_coordinates[1] = c * yOld - s * Z();
336  m_coordinates[2] = s * yOld + c * Z();
337  }
338 
339 
341  void RotateY(DataType angle)
342  {
343  //rotate vector around Y
344  const double s = std::sin((double)angle);
345  const double c = std::cos((double)angle);
346  const double zOld = Z();
347  m_coordinates[0] = s * zOld + c * X();
348  m_coordinates[2] = c * zOld - s * X();
349  }
350 
351 
353  void RotateZ(DataType angle)
354  {
355  //rotate vector around Z
356  const double s = std::sin((double)angle);
357  const double c = std::cos((double)angle);
358  const double xOld = X();
359  m_coordinates[0] = c * xOld - s * Y();
360  m_coordinates[1] = s * xOld + c * Y();
361  }
362 
364  void RotateUz(const B2Vector3<DataType>& NewUzVector)
365  {
366  // NewUzVector must be normalized !
367 
368  const double u1 = NewUzVector.X();
369  const double u2 = NewUzVector.Y();
370  const double u3 = NewUzVector.Z();
371  double up = u1 * u1 + u2 * u2;
372 
373  if (up) {
374  up = std::sqrt(up);
375  DataType px = X(), py = Y(), pz = Z();
376  m_coordinates[0] = (u1 * u3 * px - u2 * py + u1 * up * pz) / up;
377  m_coordinates[1] = (u2 * u3 * px + u1 * py + u2 * up * pz) / up;
378  m_coordinates[2] = (u3 * u3 * px - px + u3 * up * pz) / up;
379  } else if (u3 < 0.) {
380  m_coordinates[0] = -m_coordinates[0];
381  m_coordinates[2] = -m_coordinates[2];
382  }
383  }
384 
393  void Rotate(DataType alpha, const B2Vector3<DataType>& v)
394  {
395  B2Vector3<DataType> n = v.Unit();
396  *this = (n * (n.Dot(*this)) + cos(alpha) * ((n.Cross(*this)).Cross(n)) + sin(alpha) * (n.Cross(*this)));
397  }
398 
400  void Abs()
401  {
402  m_coordinates[0] = std::abs(m_coordinates[0]);
403  m_coordinates[1] = std::abs(m_coordinates[1]);
404  m_coordinates[2] = std::abs(m_coordinates[2]);
405  }
406 
408  void Sqrt()
409  {
410  Abs();
411  m_coordinates[0] = std::sqrt(m_coordinates[0]);
412  m_coordinates[1] = std::sqrt(m_coordinates[1]);
413  m_coordinates[2] = std::sqrt(m_coordinates[2]);
414  }
415 
417  DataType at(unsigned i) const;
419  DataType x() const { return m_coordinates[0]; }
421  DataType y() const { return m_coordinates[1]; }
423  DataType z() const { return m_coordinates[2]; }
425  DataType X() const { return x(); }
427  DataType Y() const { return y(); }
429  DataType Z() const { return z(); }
431  DataType Px() const { return x(); }
433  DataType Py() const { return y(); }
435  DataType Pz() const { return z(); }
436 
438  void GetXYZ(Double_t* carray) const;
440  void GetXYZ(Float_t* carray) const;
442  void GetXYZ(TVector3* tVec) const;
444  TVector3 GetTVector3() const;
445 
447  void SetX(DataType x) { m_coordinates[0] = x; }
449  void SetY(DataType y) { m_coordinates[1] = y; }
451  void SetZ(DataType z) { m_coordinates[2] = z; }
452 
454  void SetXYZ(DataType x, DataType y, DataType z)
455  {
456  SetX(x); SetY(y); SetZ(z);
457  }
459  void SetXYZ(const TVector3& tVec);
461  void SetXYZ(const TVector3* tVec);
462 
464  static std::string name();
465 
467  std::string PrintString(unsigned precision = 4) const
468  {
469  return name() + " " + PrintStringXYZ(precision) + " " + PrintStringCyl(precision);
470  }
471 
473  std::string PrintStringXYZ(unsigned precision = 4) const
474  {
475  std::ostringstream output;
476  output << "(x,y,z)=("
477  << std::fixed << std::setprecision(precision)
478  << X() << "," << Y() << "," << Z() << ")";
479  return output.str();
480  }
481 
483  std::string PrintStringCyl(unsigned precision = 4) const
484  {
485  std::ostringstream output;
486  output << "(rho, theta, phi)=("
487  << std::fixed << std::setprecision(precision)
488  << Mag() << "," << Theta() * 180. / M_PI << "," << Phi() * 180. / M_PI << ")";
489  return output.str();
490  }
491 
493  void Print()
494  {
495  //print vector parameters
496  Print(PrintString().c_str());
497  }
498 
499  };
500 
503 
506 
508  template <typename DataType>
509  Bool_t operator == (const TVector3& a, const B2Vector3<DataType>& b)
510  {
511  return (a.X() == b.X() && a.Y() == b.Y() && a.Z() == b.Z());
512  }
513 
515  template < typename DataType>
516  Bool_t operator != (const TVector3& a, const B2Vector3<DataType>& b)
517  {
518  return !(a == b);
519  }
520 
522  template < typename DataType>
524  {
525  return B2Vector3<DataType>(a * p.X(), a * p.Y(), a * p.Z());
526  }
527 
529  template < typename DataType>
531  {
532  return B2Vector3<DataType>(a.X() + b.X(), a.Y() + b.Y(), a.Z() + b.Z());
533  }
534 
536  template < typename DataType>
538  {
539  return B2Vector3<DataType>(a.X() - b.X(), a.Y() - b.Y(), a.Z() - b.Z());
540  }
541 
543  template < typename DataType>
545  {
546  return B2Vector3<DataType>(a.X() + b.X(), a.Y() + b.Y(), a.Z() + b.Z());
547  }
548 
550  template < typename DataType>
552  {
553  return B2Vector3<DataType>(a.X() - b.X(), a.Y() - b.Y(), a.Z() - b.Z());
554  }
555 
556 
558  template< typename DataType >
560  {
561  m_coordinates[0] = b.X();
562  m_coordinates[1] = b.Y();
563  m_coordinates[2] = b.Z();
564  return *this;
565  }
566 
568  template< typename DataType >
570  {
571  m_coordinates[0] = b.X();
572  m_coordinates[1] = b.Y();
573  m_coordinates[2] = b.Z();
574  return *this;
575  }
576 
578  template< typename DataType >
580  {
581  m_coordinates[0] += b.X();
582  m_coordinates[1] += b.Y();
583  m_coordinates[2] += b.Z();
584  return *this;
585  }
586 
587 
589  template< typename DataType >
591  {
592  m_coordinates[0] -= b.X();
593  m_coordinates[1] -= b.Y();
594  m_coordinates[2] -= b.Z();
595  return *this;
596  }
597 
599  template< typename DataType >
601  {
602  m_coordinates[0] *= a;
603  m_coordinates[1] *= a;
604  m_coordinates[2] *= a;
605  return *this;
606  }
607 
609  template< typename DataType >
610  void B2Vector3<DataType>::SetXYZ(const TVector3& tVec)
611  {
612  m_coordinates[0] = static_cast<Double_t>(tVec.X());
613  m_coordinates[1] = static_cast<Double_t>(tVec.Y());
614  m_coordinates[2] = static_cast<Double_t>(tVec.Z());
615  }
616 
618  template< typename DataType >
619  void B2Vector3<DataType>::SetXYZ(const TVector3* tVec)
620  {
621  m_coordinates[0] = static_cast<Double_t>(tVec->X());
622  m_coordinates[1] = static_cast<Double_t>(tVec->Y());
623  m_coordinates[2] = static_cast<Double_t>(tVec->Z());
624  }
625 
626  template< typename DataType >
627  void B2Vector3<DataType>::GetXYZ(double* carray) const
628  {
629  carray[0] = X();
630  carray[1] = Y();
631  carray[2] = Z();
632  }
633 
635  template< typename DataType >
636  void B2Vector3<DataType>::GetXYZ(TVector3* tVec) const
637  {
638  tVec->SetXYZ(static_cast<Double_t>(X()),
639  static_cast<Double_t>(Y()),
640  static_cast<Double_t>(Z()));
641  }
642 
643 
645  template< typename DataType >
647  {
648  return
649  TVector3(
650  static_cast<Double_t>(X()),
651  static_cast<Double_t>(Y()),
652  static_cast<Double_t>(Z())
653  );
654  }
655 
656 
658  template < typename DataType>
659  DataType B2Vector3<DataType>::at(unsigned i) const
660  {
661  switch (i) {
662  case 0:
664  case 1:
666  case 2:
668  }
669  B2FATAL(this->name() << "::access operator: given index (i=" << i << ") is out of bounds!");
670  return 0.;
671  }
672 
674  template < typename DataType>
676  {
677  return std::string("B2Vector3<") + typeid(DataType).name() + std::string(">");
678  }
679 
681 } // end namespace Belle2
A fast and root compatible alternative to TVector3.
Definition: B2Vector3.h:42
DataType Phi() const
The azimuth angle.
Definition: B2Vector3.h:145
B2Vector3< DataType > operator+(const B2Vector3< DataType > &b) const
Addition of 3-vectors.
Definition: B2Vector3.h:121
B2Vector3< DataType > Cross(const B2Vector3< DataType > &p) const
Cross product.
Definition: B2Vector3.h:290
void SetMag(DataType mag)
Set magnitude keeping theta and phi constant.
Definition: B2Vector3.h:176
DataType Pz() const
access variable Z (= .at(2) without boundary check)
Definition: B2Vector3.h:435
DataType m_coordinates[3]
Make sure that we only have floating point vectors.
Definition: B2Vector3.h:45
B2Vector3(const B2Vector3< DataType > &b2Vec3)
Constructor expecting a B2Vector3 of same type.
Definition: B2Vector3.h:67
B2Vector3< DataType > operator/(DataType a) const
Scaling of 3-vectors with a real number.
Definition: B2Vector3.h:136
DataType & operator()(unsigned i)
member access without boundary check
Definition: B2Vector3.h:89
DataType operator()(unsigned i) const
member access without boundary check
Definition: B2Vector3.h:85
B2Vector3< DataType > operator*(DataType a) const
Scaling of 3-vectors with a real number.
Definition: B2Vector3.h:131
DataType Z() const
access variable Z (= .at(2) without boundary check)
Definition: B2Vector3.h:429
DataType Perp2() const
The transverse component squared (R^2 in cylindrical coordinate system).
Definition: B2Vector3.h:190
void SetPerp(DataType r)
Set the transverse component keeping phi and z constant.
Definition: B2Vector3.h:197
DataType y() const
access variable Y (= .at(1) without boundary check)
Definition: B2Vector3.h:421
void SetX(DataType x)
set X/1st-coordinate
Definition: B2Vector3.h:447
DataType Theta() const
The polar angle.
Definition: B2Vector3.h:147
DataType CosTheta() const
Cosine of the polar angle.
Definition: B2Vector3.h:149
DataType & operator[](unsigned i)
member access without boundary check
Definition: B2Vector3.h:91
void GetXYZ(Float_t *carray) const
directly copies coordinates to an array of float
DataType z() const
access variable Z (= .at(2) without boundary check)
Definition: B2Vector3.h:423
B2Vector3(const B2Vector3< DataType > *b2Vec3)
Constructor expecting a pointer to a B2Vector3.
Definition: B2Vector3.h:69
void SetMagThetaPhi(DataType mag, DataType theta, DataType phi)
setter with mag, theta, phi
Definition: B2Vector3.h:253
DataType X() const
access variable X (= .at(0) without boundary check)
Definition: B2Vector3.h:425
std::string PrintStringXYZ(unsigned precision=4) const
create a string containing vector in cartesian coordinates
Definition: B2Vector3.h:473
DataType Eta() const
Returns the pseudo-rapidity.
Definition: B2Vector3.h:325
DataType DeltaPhi(const B2Vector3< DataType > &v) const
returns phi in the interval [-PI,PI)
Definition: B2Vector3.h:222
void RotateY(DataType angle)
Rotates the B2Vector3 around the y-axis.
Definition: B2Vector3.h:341
DataType Y() const
access variable Y (= .at(1) without boundary check)
Definition: B2Vector3.h:427
DataType operator[](unsigned i) const
member access without boundary check
Definition: B2Vector3.h:87
std::string PrintString(unsigned precision=4) const
create a string containing vector in cartesian and spherical coordinates
Definition: B2Vector3.h:467
DataType Mag() const
The magnitude (rho in spherical coordinate system).
Definition: B2Vector3.h:153
DataType value_type
storage type of the vector
Definition: B2Vector3.h:50
B2Vector3(const ROOT::Math::XYZVector &xyzVec)
Constructor expecting a XYZVector.
Definition: B2Vector3.h:79
bool operator==(const B2Vector3< DataType > &b) const
Comparison for equality with a B2Vector3.
Definition: B2Vector3.h:104
DataType x() const
access variable X (= .at(0) without boundary check)
Definition: B2Vector3.h:419
B2Vector3(const DataType(*coords)[3])
Constructor using a pointer.
Definition: B2Vector3.h:59
B2Vector3< DataType > operator-() const
unary minus
Definition: B2Vector3.h:119
DataType DeltaR(const B2Vector3< DataType > &v) const
return deltaR with respect to input-vector
Definition: B2Vector3.h:239
DataType Perp(const B2Vector3< DataType > &axis) const
The transverse component w.r.t.
Definition: B2Vector3.h:220
void RotateX(DataType angle)
Rotates the B2Vector3 around the x-axis.
Definition: B2Vector3.h:329
void Abs()
calculates the absolute value of the coordinates element-wise
Definition: B2Vector3.h:400
B2Vector3(const DataType(&coords)[3])
Constructor using a reference.
Definition: B2Vector3.h:57
B2Vector3(const ROOT::Math::XYZVector *xyzVec)
Constructor expecting a pointer to a XYZVector.
Definition: B2Vector3.h:82
DataType Pt(const B2Vector3< DataType > &axis) const
The transverse component w.r.t.
Definition: B2Vector3.h:218
DataType Mag2() const
The magnitude squared (rho^2 in spherical coordinate system).
Definition: B2Vector3.h:151
void SetTheta(DataType theta)
Set theta keeping mag and phi constant.
Definition: B2Vector3.h:164
B2Vector3(const TVector3 &tVec3)
Constructor expecting a TVector3.
Definition: B2Vector3.h:62
void RotateZ(DataType angle)
Rotates the B2Vector3 around the z-axis.
Definition: B2Vector3.h:353
void Print()
just for backward compatibility, should not be used with new code
Definition: B2Vector3.h:493
void Sqrt()
calculates the square root of the absolute values of the coordinates element-wise
Definition: B2Vector3.h:408
B2Vector3< DataType > Unit() const
Unit vector parallel to this.
Definition: B2Vector3.h:263
void SetZ(DataType z)
set Z/3rd-coordinate
Definition: B2Vector3.h:451
DataType Dot(const B2Vector3< DataType > &p) const
Scalar product.
Definition: B2Vector3.h:284
B2Vector3< DataType > Orthogonal() const
Vector orthogonal to this one.
Definition: B2Vector3.h:271
void SetY(DataType y)
set Y/2nd-coordinate
Definition: B2Vector3.h:449
B2Vector3(const B2Vector3< OtherType > *b2Vec3)
Constructor expecting a pointer to a B2Vector3 of different type.
Definition: B2Vector3.h:75
DataType PseudoRapidity() const
Returns the pseudo-rapidity, i.e.
Definition: B2Vector3.h:313
DataType DrEtaPhi(const B2Vector3< DataType > &v) const
return DrEtaPhi with respect to input-vector
Definition: B2Vector3.h:247
bool operator!=(const B2Vector3< DataType > &b) const
Comparison != with a B2Vector3.
Definition: B2Vector3.h:108
DataType Perp() const
The transverse component (R in cylindrical coordinate system).
Definition: B2Vector3.h:194
B2Vector3(void)
empty Constructor sets everything to 0
Definition: B2Vector3.h:53
DataType Perp2(const B2Vector3< DataType > &axis) const
The transverse component w.r.t.
Definition: B2Vector3.h:207
B2Vector3(const B2Vector3< OtherType > &b2Vec3)
Constructor expecting a B2Vector3 of different type.
Definition: B2Vector3.h:72
void GetXYZ(Double_t *carray) const
directly copies coordinates to an array of double
void RotateUz(const B2Vector3< DataType > &NewUzVector)
Rotates reference frame from Uz to newUz (unit vector).
Definition: B2Vector3.h:364
DataType Py() const
access variable Y (= .at(1) without boundary check)
Definition: B2Vector3.h:433
static DataType Mpi_pi(DataType angle)
returns given angle in the interval [-PI,PI)
Definition: B2Vector3.h:226
DataType Pt() const
The transverse component (R in cylindrical coordinate system).
Definition: B2Vector3.h:192
B2Vector3(const TVector3 *tVec3)
Constructor expecting a pointer to a TVector3.
Definition: B2Vector3.h:65
DataType Px() const
access variable X (= .at(0) without boundary check)
Definition: B2Vector3.h:431
void SetXYZ(DataType x, DataType y, DataType z)
set all coordinates using data type
Definition: B2Vector3.h:454
DataType Angle(const B2Vector3< DataType > &q) const
The angle w.r.t.
Definition: B2Vector3.h:296
std::string PrintStringCyl(unsigned precision=4) const
create a string containing vector in spherical coordinates
Definition: B2Vector3.h:483
void Rotate(DataType alpha, const B2Vector3< DataType > &v)
Rotation around an arbitrary axis v with angle alpha.
Definition: B2Vector3.h:393
void SetPhi(DataType phi)
Set phi keeping mag and theta constant.
Definition: B2Vector3.h:156
B2Vector3(const DataType xVal, const DataType yVal, const DataType zVal)
Constructor expecting 3 coordinates.
Definition: B2Vector3.h:55
bool operator==(const DecayNode &node1, const DecayNode &node2)
Compare two Decay Nodes: They are equal if All daughter decay nodes are equal or one of the daughter ...
Definition: DecayNode.cc:48
bool operator!=(const DecayNode &node1, const DecayNode &node2)
Not equal: See operator==.
Definition: DecayNode.cc:65
void SetXYZ(const TVector3 &tVec)
set all coordinates using a reference to TVector3
Definition: B2Vector3.h:610
B2Vector3< DataType > & operator-=(const B2Vector3< DataType > &b)
subtraction
Definition: B2Vector3.h:590
TVector3 GetTVector3() const
returns a TVector3 containing the same coordinates
Definition: B2Vector3.h:646
B2Vector3< float > B2Vector3F
typedef for common usage with float
Definition: B2Vector3.h:505
void GetXYZ(TVector3 *tVec) const
directly copies coordinates to a TVector3
Definition: B2Vector3.h:636
void SetXYZ(const TVector3 *tVec)
set all coordinates using a pointer to TVector3
Definition: B2Vector3.h:619
B2Vector3< DataType > & operator+=(const B2Vector3< DataType > &b)
addition
Definition: B2Vector3.h:579
B2Vector3< DataType > operator*(DataType a, const B2Vector3< DataType > &p)
non-memberfunction Scaling of 3-vectors with a real number
Definition: B2Vector3.h:523
B2Vector3< DataType > & operator*=(DataType a)
scaling with real numbers
Definition: B2Vector3.h:600
B2Vector3< DataType > & operator=(const B2Vector3< DataType > &b)
Assignment via B2Vector3.
Definition: B2Vector3.h:559
B2Vector3< double > B2Vector3D
typedef for common usage with double
Definition: B2Vector3.h:502
DataType at(unsigned i) const
safe member access (with boundary check!)
Definition: B2Vector3.h:659
static std::string name()
Returns the name of the B2Vector.
Definition: B2Vector3.h:675
B2Vector3< DataType > operator-(const TVector3 &a, const B2Vector3< DataType > &b)
non-memberfunction for substracting a TVector3 from a B2Vector3
Definition: B2Vector3.h:537
B2Vector3< DataType > operator+(const TVector3 &a, const B2Vector3< DataType > &b)
non-memberfunction for adding a TVector3 to a B2Vector3
Definition: B2Vector3.h:530
Abstract base class for different kinds of events.
Definition: ClusterUtils.h:23