ThePEG/Vectors/Lorentz5Vector.h

0001 // -*- C++ -*-
0002 //
0003 // Lorentz5Vector.h is a part of ThePEG - Toolkit for HEP Event Generation
0004 // Copyright (C) 1999-2019 Leif Lonnblad
0005 //
0006 // ThePEG is licenced under version 3 of the GPL, see COPYING for details.
0007 // Please respect the MCnet academic guidelines, see GUIDELINES for details.
0008 //
0009 #ifndef ThePEG_Lorentz5Vector_H
0010 #define ThePEG_Lorentz5Vector_H
0011
0012 // This is the declaration of the Lorentz5vector class.
0013
0014 #include "LorentzVector.h"
0015 #include "Lorentz5Vector.fh"
0016 #include "ThePEG/Utilities/Maths.h"
0017 #include "ThePEG/Utilities/Direction.h"
0018 #include "ThePEG/Utilities/UnitIO.h"
0019 #include "LorentzRotation.h"
0020
0021 namespace ThePEG {
0022
0023 template <typename Value>
0024 /**
0025  * The Lorentz5Vector inherits from the
0026  * <code>LorentzVector</code> class. It is templated on the
0027  * type of the member variables. The <code>Lorentz5Vector</code> class is a
0028  * <code>LorentzVector</code> with an extra member for the invariant
0029  * length/mass of the vector. Note that an object of the
0030  * <code>Lorentz5Vector</code> class may be internally inconsistent in
0031  * that the invariant length/mass of the <code>LorentzVector</code>
0032  * class need not be the same as the member variable representing the
0033  * invariant length/mass. The degree of inconsistency can be accessed
0034  * with the <code>massError()</code>, <code>energyError()</code> and
0035  * <code>rhoError()</code> methods and an object can be made consistent
0036  * using the <code>rescaleMass()</code>, <code>rescaleEnergy()</code> or
0037  * <code>rescaleRho()</code> methods.
0038  *
0039  * @see Math
0040  *
0041  */
0042 class Lorentz5Vector: public LorentzVector<Value> {
0043
0044 public:
0045
0046   /** Template argument typedef. */
0047   using Value2 = decltype(sqr(std::declval<Value>()));
0048
0049 public:
0050   /// Component access.
0051   //@{
0052   Value x() const { return LorentzVector<Value>::x(); }
0053   Value y() const { return LorentzVector<Value>::y(); }
0054   Value z() const { return LorentzVector<Value>::z(); }
0055   Value t() const { return LorentzVector<Value>::t(); }
0056   //@}
0057
0058 public:
0059
0060   /** @name Constructors and destructor. */
0061   //@{
0062   /**
0063    * Constructor giving the null vector.
0064    */
0065   Lorentz5Vector() : mm() {}
0066
0067   /**
0068    * Constructor giving the invariant length.
0069    */
0070   Lorentz5Vector(Value m)
0071     : LorentzVector<Value>(Value(), Value(), Value(), m), mm(m) {}
0072
0073   /**
0074    * Constructor giving the components x, y, z, t. The invariant
0075    * length is set to LorentzVector::mag().
0076    */
0077   Lorentz5Vector(Value x, Value y, Value z, Value t = Value())
0078     : LorentzVector<Value>(x, y, z, t) { rescaleMass(); }
0079
0080   /**
0081    * Constructor giving the components x, y, z, t and invariant length.
0082    * May result in an inconsistent Lorentz5Vector.
0083    */
0084   Lorentz5Vector(Value x, Value y, Value z, Value t, Value tau)
0085     : LorentzVector<Value>(x, y, z, t), mm(tau) {}
0086
0087   /**
0088    * Constructor giving a 3-Vector and a time component. The invariant
0089    * length is set to LorentzVector::mag().
0090    */
0091   Lorentz5Vector(const ThreeVector<Value> & p, Value e)
0092     : LorentzVector<Value>(p, e) { rescaleMass(); }
0093
0094   /**
0095    * Constructor giving an invariant length and a 3-Vector
0096    * component. The time component is set to the corresponding value.
0097    */
0098   Lorentz5Vector(Value m, const ThreeVector<Value> & p)
0099     : LorentzVector<Value>(p, sqrt(p.mag2() + m*m)), mm(m) {}
0100
0101   /**
0102    * Constructor giving a 3-Vector, a time component and an invariant
0103    * length. May result in an inconsistent Lorentz5Vector.
0104    */
0105   Lorentz5Vector(const ThreeVector<Value> & p, Value t, Value tau)
0106     : LorentzVector<Value>(p, t), mm(tau) {}
0107
0108   /**
0109    * Constructor giving a LorentzVector and an invariant length.
0110    * May result in an inconsistent Lorentz5Vector.
0111    */
0112   Lorentz5Vector(const LorentzVector<Value> & p, Value m)
0113     : LorentzVector<Value>(p), mm(m) {}
0114
0115   /**
0116    * Copy from HepLorentzVector constructor. The invariant
0117    * length is set to LorentzVector::mag().
0118    */
0119   Lorentz5Vector(const LorentzVector<Value> & p)
0120     : LorentzVector<Value>(p) { rescaleMass(); }
0121
0122   /**
0123    * Construct from value type U convertible to Value.
0124    */
0125   template<class U>
0126   Lorentz5Vector(const Lorentz5Vector<U> & p)
0127     : LorentzVector<Value>(p), mm(p.m) {}
0128   //@}
0129
0130   /** @name Assignment and set functions. */
0131   //@{
0132   /**
0133    * Set invariant length/mass.
0134    */
0135   void setTau(Value a) { mm = a; }
0136
0137   /**
0138    * Set invariant length/mass.
0139    */
0140   void setMass(Value a) { mm = a; }
0141
0142   /**
0143    * Assignment. The invariant length is kept fixed. May result in an
0144    * inconsistent Lorentz5Vector.
0145    */
0146   Lorentz5Vector & operator=(const LorentzVector<Value> & q) {
0147     LorentzVector<Value>::operator=(q);
0148     return *this;
0149   }
0150   //@}
0151
0152   /** @name Rescale functions to make consistent. */
0153   //@{
0154   /**
0155    * Rescale energy, so that the invariant length/mass of the
0156    * LorentzVector agrees with the current one.
0157    */
0158   void rescaleEnergy() {
0159     LorentzVector<Value>::setT(sqrt(LorentzVector<Value>::vect().mag2() + mass2()));
0160   }
0161
0162   /**
0163    * Rescale spatial component, so that the invariant length/mass of
0164    * the LorentzVector agrees with the current one.
0165    */
0166   void rescaleRho() {
0167     LorentzVector<Value>::setRho(sqrt(t()*t() - mass2()));
0168   }
0169
0170   /**
0171    * Set the invariant length/mass member, so that it agrees with the
0172    * invariant length/mass of the LorentzVector.
0173    */
0174   void rescaleMass() {
0175     mm = LorentzVector<Value>::m();
0176   }
0177   //@}
0178
0179   /** @name Check consistency. */
0180   //@{
0181   /**
0182    * Return the relative inconsistency in the mass component.
0183    */
0184   double massError() const {
0185     return sqrt(abs(Math::relativeError(mass2(),
0186                     LorentzVector<Value>::m2())));
0187   }
0188
0189   /**
0190    * Return the relative inconsistency in the energy component.
0191    */
0192   double energyError() const {
0193     return sqrt(abs(Math::relativeError(t()*t(), mass2()
0194                     + LorentzVector<Value>::vect().mag2())));
0195   }
0196
0197   /**
0198    * Return the relative inconsistency in the spatial components.
0199    */
0200   double rhoError() const {
0201     return sqrt(abs(Math::relativeError(LorentzVector<Value>::vect().mag2(),
0202                     t()*t() - mass2())));
0203   }
0204   //@}
0205
0206   /** @name Access components. */
0207   //@{
0208   /**
0209    * Mass/invariant length component squared. m2() gives
0210    * the same calculated from the LorentzVector
0211    */
0212   Value2 mass2() const { return mm > Value() ? mm*mm: -mm*mm; }
0213
0214   /**
0215    * Mass/invariant length component squared. m2() gives
0216    * the same calculated from the LorentzVector
0217    */
0218   Value2 tau2() const { return mass2(); }
0219
0220   /**
0221    * Mass/invariant length component. m() gives the same
0222    * calculated from the LorentzVector
0223    */
0224   Value mass() const { return mm; }
0225
0226
0227   /**
0228    * Mass/invariant length component. m() gives the same
0229    * calculated from the LorentzVector
0230    */
0231   Value tau() const { return mass(); }
0232
0233   /**
0234    * Return the positive negative light-cone components (depending on
0235    * the value of Direction<0>.
0236    */
0237   Value dirPlus() const {
0238     return Direction<0>::pos() ?
0239       LorentzVector<Value>::plus()
0240       :
0241       LorentzVector<Value>::minus();
0242   }
0243
0244   /**
0245    * Return the positive negative light-cone components (depending on
0246    * the value of Direction<0>.
0247    */
0248   Value dirMinus() const {
0249     return Direction<0>::neg() ?
0250       LorentzVector<Value>::plus()
0251       :
0252       LorentzVector<Value>::minus();
0253   }
0254   //@}
0255
0256   /**
0257    *  Perform a Lorentz transformation
0258    */
0259   Lorentz5Vector & transform(const LorentzRotation & r)
0260   {
0261     LorentzVector<Value>::transform(r.one());
0262     return *this;
0263   }
0264
0265 private:
0266
0267   /** The invariant mass/length member. */
0268   Value mm;
0269
0270 };
0271
0272 /** Output a Lorentz5Vector to a stream. */
0273 template <typename OStream, typename T, typename UT>
0274 void ounitstream(OStream & os, const Lorentz5Vector<T> & p, UT & u) {
0275   os << ounit(p.x(), u) << ounit(p.y(), u) << ounit(p.z(), u)
0276      << ounit(p.e(), u) << ounit(p.mass(), u);
0277 }
0278
0279 /** Input a Lorentz5Vector from a stream. */
0280 template <typename IStream, typename T, typename UT>
0281 void iunitstream(IStream & is, Lorentz5Vector<T> & p, UT & u) {
0282   T x, y, z, e, mass;
0283   is >> iunit(x, u) >> iunit(y, u) >> iunit(z, u) >> iunit(e, u)
0284      >> iunit(mass, u);
0285   p = Lorentz5Vector<T>(x, y, z, e, mass);
0286 }
0287
0288
0289 /// @name Dot product overloads.
0290 //@{
0291 template <typename ValueA, typename ValueB>
0292 inline auto
0293 operator*(const Lorentz5Vector<ValueA> & a, const Lorentz5Vector<ValueB> & b)
0294 -> decltype(a.dot(b))
0295 {
0296   return a.dot(b);
0297 }
0298
0299 template <typename ValueA, typename ValueB>
0300 inline auto
0301 operator*(const LorentzVector<ValueA> & a, const Lorentz5Vector<ValueB> & b)
0302 -> decltype(a.dot(b))
0303 {
0304   return a.dot(b);
0305 }
0306
0307 template <typename ValueA, typename ValueB>
0308 inline auto
0309 operator*(const Lorentz5Vector<ValueA> & a, const LorentzVector<ValueB> & b)
0310 -> decltype(a.dot(b))
0311 {
0312   return a.dot(b);
0313 }
0314
0315 template <typename Value>
0316 inline auto
0317 operator*(const Lorentz5Vector<Value> & a, const Lorentz5Vector<Value> & b)
0318 -> decltype(a.dot(b))
0319 {
0320   return a.dot(b);
0321 }
0322 //@}
0323 }
0324
0325 #endif /* ThePEG_Particle_H */