Logo ROOT  
Reference Guide
LorentzVector.h
Go to the documentation of this file.
1// @(#)root/mathcore:$Id$
2// Authors: W. Brown, M. Fischler, L. Moneta 2005
3
4/**********************************************************************
5 * *
6 * Copyright (c) 2005 , LCG ROOT MathLib Team *
7 * *
8 * *
9 **********************************************************************/
10
11// Header file for class LorentzVector
12//
13// Created by: moneta at Tue May 31 17:06:09 2005
14// Major mods by: fischler at Wed Jul 20 2005
15//
16// Last update: $Id$
17//
18#ifndef ROOT_Math_GenVector_LorentzVector
19#define ROOT_Math_GenVector_LorentzVector 1
20
22
24
26
27#include <cmath>
28
29namespace ROOT {
30
31 namespace Math {
32
33//__________________________________________________________________________________________
34 /**
35 Class describing a generic LorentzVector in the 4D space-time,
36 using the specified coordinate system for the spatial vector part.
37 The metric used for the LorentzVector is (-,-,-,+).
38 In the case of LorentzVector we don't distinguish the concepts
39 of points and displacement vectors as in the 3D case,
40 since the main use case for 4D Vectors is to describe the kinematics of
41 relativistic particles. A LorentzVector behaves like a
42 DisplacementVector in 4D. The Minkowski components could be viewed as
43 v and t, or for kinematic 4-vectors, as p and E.
44
45 ROOT provides specialisations and aliases to them of the ROOT::Math::LorentzVector template:
46 - ROOT::Math::PtEtaPhiMVector based on pt (rho),eta,phi and M (t) coordinates in double precision
47 - ROOT::Math::PtEtaPhiEVector based on pt (rho),eta,phi and E (t) coordinates in double precision
48 - ROOT::Math::PxPyPzMVector based on px,py,pz and M (mass) coordinates in double precision
49 - ROOT::Math::PxPyPzEVector based on px,py,pz and E (energy) coordinates in double precision
50 - ROOT::Math::XYZTVector based on x,y,z,t coordinates (cartesian) in double precision (same as PxPyPzEVector)
51 - ROOT::Math::XYZTVectorF based on x,y,z,t coordinates (cartesian) in float precision (same as PxPyPzEVector but float)
52
53More details about the GenVector package can be found [here](Vector.html).
54
55
56 @ingroup GenVector
57 */
58 template< class CoordSystem >
60
61 public:
62
63 // ------ ctors ------
64
65 typedef typename CoordSystem::Scalar Scalar;
66 typedef CoordSystem CoordinateType;
67
68 /**
69 default constructor of an empty vector (Px = Py = Pz = E = 0 )
70 */
72
73 /**
74 generic constructors from four scalar values.
75 The association between values and coordinate depends on the
76 coordinate system. For PxPyPzE4D,
77 \param a scalar value (Px)
78 \param b scalar value (Py)
79 \param c scalar value (Pz)
80 \param d scalar value (E)
81 */
83 const Scalar & b,
84 const Scalar & c,
85 const Scalar & d) :
86 fCoordinates(a , b, c, d) { }
87
88 /**
89 constructor from a LorentzVector expressed in different
90 coordinates, or using a different Scalar type
91 */
92 template< class Coords >
95
96 /**
97 Construct from a foreign 4D vector type, for example, HepLorentzVector
98 Precondition: v must implement methods x(), y(), z(), and t()
99 */
100 template<class ForeignLorentzVector>
101 explicit LorentzVector( const ForeignLorentzVector & v) :
102 fCoordinates(PxPyPzE4D<Scalar>( v.x(), v.y(), v.z(), v.t() ) ) { }
103
104#ifdef LATER
105 /**
106 construct from a generic linear algebra vector implementing operator []
107 and with a size of at least 4. This could be also a C array
108 In this case v[0] is the first data member
109 ( Px for a PxPyPzE4D base)
110 \param v LA vector
111 \param index0 index of first vector element (Px)
112 */
113 template< class LAVector >
114 explicit LorentzVector(const LAVector & v, size_t index0 ) {
115 fCoordinates = CoordSystem ( v[index0], v[index0+1], v[index0+2], v[index0+3] );
116 }
117#endif
118
119
120 // ------ assignment ------
121
122 /**
123 Assignment operator from a lorentz vector of arbitrary type
124 */
125 template< class OtherCoords >
127 fCoordinates = v.Coordinates();
128 return *this;
129 }
130
131 /**
132 assignment from any other Lorentz vector implementing
133 x(), y(), z() and t()
134 */
135 template<class ForeignLorentzVector>
136 LorentzVector & operator = ( const ForeignLorentzVector & v) {
137 SetXYZT( v.x(), v.y(), v.z(), v.t() );
138 return *this;
139 }
140
141#ifdef LATER
142 /**
143 assign from a generic linear algebra vector implementing operator []
144 and with a size of at least 4
145 In this case v[0] is the first data member
146 ( Px for a PxPyPzE4D base)
147 \param v LA vector
148 \param index0 index of first vector element (Px)
149 */
150 template< class LAVector >
151 LorentzVector & AssignFrom(const LAVector & v, size_t index0=0 ) {
152 fCoordinates.SetCoordinates( v[index0], v[index0+1], v[index0+2], v[index0+3] );
153 return *this;
154 }
155#endif
156
157 // ------ Set, Get, and access coordinate data ------
158
159 /**
160 Retrieve a const reference to the coordinates object
161 */
162 const CoordSystem & Coordinates() const {
163 return fCoordinates;
164 }
165
166 /**
167 Set internal data based on an array of 4 Scalar numbers
168 */
170 fCoordinates.SetCoordinates(src);
171 return *this;
172 }
173
174 /**
175 Set internal data based on 4 Scalar numbers
176 */
178 fCoordinates.SetCoordinates(a, b, c, d);
179 return *this;
180 }
181
182 /**
183 Set internal data based on 4 Scalars at *begin to *end
184 */
185 template< class IT >
187 IT a = begin; IT b = ++begin; IT c = ++begin; IT d = ++begin;
188 (void)end;
189 assert (++begin==end);
190 SetCoordinates (*a,*b,*c,*d);
191 return *this;
192 }
193
194 /**
195 get internal data into 4 Scalar numbers
196 */
197 void GetCoordinates( Scalar& a, Scalar& b, Scalar& c, Scalar & d ) const
198 { fCoordinates.GetCoordinates(a, b, c, d); }
199
200 /**
201 get internal data into an array of 4 Scalar numbers
202 */
203 void GetCoordinates( Scalar dest[] ) const
204 { fCoordinates.GetCoordinates(dest); }
205
206 /**
207 get internal data into 4 Scalars at *begin to *end
208 */
209 template <class IT>
210 void GetCoordinates( IT begin, IT end ) const
211 { IT a = begin; IT b = ++begin; IT c = ++begin; IT d = ++begin;
212 (void)end;
213 assert (++begin==end);
214 GetCoordinates (*a,*b,*c,*d);
215 }
216
217 /**
218 get internal data into 4 Scalars at *begin
219 */
220 template <class IT>
221 void GetCoordinates( IT begin ) const {
222 Scalar a,b,c,d = 0;
223 GetCoordinates (a,b,c,d);
224 *begin++ = a;
225 *begin++ = b;
226 *begin++ = c;
227 *begin = d;
228 }
229
230 /**
231 set the values of the vector from the cartesian components (x,y,z,t)
232 (if the vector is held in another coordinates, like (Pt,eta,phi,m)
233 then (x, y, z, t) are converted to that form)
234 */
236 fCoordinates.SetPxPyPzE(xx,yy,zz,tt);
237 return *this;
238 }
240 fCoordinates.SetPxPyPzE(xx,yy,zz,ee);
241 return *this;
242 }
243
244 // ------------------- Equality -----------------
245
246 /**
247 Exact equality
248 */
249 bool operator==(const LorentzVector & rhs) const {
250 return fCoordinates==rhs.fCoordinates;
251 }
252 bool operator!= (const LorentzVector & rhs) const {
253 return !(operator==(rhs));
254 }
255
256 // ------ Individual element access, in various coordinate systems ------
257
258 // individual coordinate accessors in various coordinate systems
259
260 /**
261 spatial X component
262 */
263 Scalar Px() const { return fCoordinates.Px(); }
264 Scalar X() const { return fCoordinates.Px(); }
265 /**
266 spatial Y component
267 */
268 Scalar Py() const { return fCoordinates.Py(); }
269 Scalar Y() const { return fCoordinates.Py(); }
270 /**
271 spatial Z component
272 */
273 Scalar Pz() const { return fCoordinates.Pz(); }
274 Scalar Z() const { return fCoordinates.Pz(); }
275 /**
276 return 4-th component (time, or energy for a 4-momentum vector)
277 */
278 Scalar E() const { return fCoordinates.E(); }
279 Scalar T() const { return fCoordinates.E(); }
280 /**
281 return magnitude (mass) squared M2 = T**2 - X**2 - Y**2 - Z**2
282 (we use -,-,-,+ metric)
283 */
284 Scalar M2() const { return fCoordinates.M2(); }
285 /**
286 return magnitude (mass) using the (-,-,-,+) metric.
287 If M2 is negative (space-like vector) a GenVector_exception
288 is suggested and if continuing, - sqrt( -M2) is returned
289 */
290 Scalar M() const { return fCoordinates.M();}
291 /**
292 return the spatial (3D) magnitude ( sqrt(X**2 + Y**2 + Z**2) )
293 */
294 Scalar R() const { return fCoordinates.R(); }
295 Scalar P() const { return fCoordinates.R(); }
296 /**
297 return the square of the spatial (3D) magnitude ( X**2 + Y**2 + Z**2 )
298 */
299 Scalar P2() const { return P() * P(); }
300 /**
301 return the square of the transverse spatial component ( X**2 + Y**2 )
302 */
303 Scalar Perp2( ) const { return fCoordinates.Perp2();}
304
305 /**
306 return the transverse spatial component sqrt ( X**2 + Y**2 )
307 */
308 Scalar Pt() const { return fCoordinates.Pt(); }
309 Scalar Rho() const { return fCoordinates.Pt(); }
310
311 /**
312 return the transverse mass squared
313 \f[ m_t^2 = E^2 - p{_z}^2 \f]
314 */
315 Scalar Mt2() const { return fCoordinates.Mt2(); }
316
317 /**
318 return the transverse mass
319 \f[ \sqrt{ m_t^2 = E^2 - p{_z}^2} X sign(E^ - p{_z}^2) \f]
320 */
321 Scalar Mt() const { return fCoordinates.Mt(); }
322
323 /**
324 return the transverse energy squared
325 \f[ e_t = \frac{E^2 p_{\perp}^2 }{ |p|^2 } \f]
326 */
327 Scalar Et2() const { return fCoordinates.Et2(); }
328
329 /**
330 return the transverse energy
331 \f[ e_t = \sqrt{ \frac{E^2 p_{\perp}^2 }{ |p|^2 } } X sign(E) \f]
332 */
333 Scalar Et() const { return fCoordinates.Et(); }
334
335 /**
336 azimuthal Angle
337 */
338 Scalar Phi() const { return fCoordinates.Phi();}
339
340 /**
341 polar Angle
342 */
343 Scalar Theta() const { return fCoordinates.Theta(); }
344
345 /**
346 pseudorapidity
347 \f[ \eta = - \ln { \tan { \frac { \theta} {2} } } \f]
348 */
349 Scalar Eta() const { return fCoordinates.Eta(); }
350
351 /**
352 get the spatial components of the Vector in a
353 DisplacementVector based on Cartesian Coordinates
354 */
356 return ::ROOT::Math::DisplacementVector3D<Cartesian3D<Scalar> >( X(), Y(), Z() );
357 }
358
359 // ------ Operations combining two Lorentz vectors ------
360
361 /**
362 scalar (Dot) product of two LorentzVector vectors (metric is -,-,-,+)
363 Enable the product using any other LorentzVector implementing
364 the x(), y() , y() and t() member functions
365 \param q any LorentzVector implementing the x(), y() , z() and t()
366 member functions
367 \return the result of v.q of type according to the base scalar type of v
368 */
369
370 template< class OtherLorentzVector >
371 Scalar Dot(const OtherLorentzVector & q) const {
372 return t()*q.t() - x()*q.x() - y()*q.y() - z()*q.z();
373 }
374
375 /**
376 Self addition with another Vector ( v+= q )
377 Enable the addition with any other LorentzVector
378 \param q any LorentzVector implementing the x(), y() , z() and t()
379 member functions
380 */
381 template< class OtherLorentzVector >
382 inline LorentzVector & operator += ( const OtherLorentzVector & q)
383 {
384 SetXYZT( x() + q.x(), y() + q.y(), z() + q.z(), t() + q.t() );
385 return *this;
386 }
387
388 /**
389 Self subtraction of another Vector from this ( v-= q )
390 Enable the addition with any other LorentzVector
391 \param q any LorentzVector implementing the x(), y() , z() and t()
392 member functions
393 */
394 template< class OtherLorentzVector >
395 LorentzVector & operator -= ( const OtherLorentzVector & q) {
396 SetXYZT( x() - q.x(), y() - q.y(), z() - q.z(), t() - q.t() );
397 return *this;
398 }
399
400 /**
401 addition of two LorentzVectors (v3 = v1 + v2)
402 Enable the addition with any other LorentzVector
403 \param v2 any LorentzVector implementing the x(), y() , z() and t()
404 member functions
405 \return a new LorentzVector of the same type as v1
406 */
407 template<class OtherLorentzVector>
408 LorentzVector operator + ( const OtherLorentzVector & v2) const
409 {
411 v3 += v2;
412 return v3;
413 }
414
415 /**
416 subtraction of two LorentzVectors (v3 = v1 - v2)
417 Enable the subtraction of any other LorentzVector
418 \param v2 any LorentzVector implementing the x(), y() , z() and t()
419 member functions
420 \return a new LorentzVector of the same type as v1
421 */
422 template<class OtherLorentzVector>
423 LorentzVector operator - ( const OtherLorentzVector & v2) const {
425 v3 -= v2;
426 return v3;
427 }
428
429 //--- scaling operations ------
430
431 /**
432 multiplication by a scalar quantity v *= a
433 */
435 fCoordinates.Scale(a);
436 return *this;
437 }
438
439 /**
440 division by a scalar quantity v /= a
441 */
443 fCoordinates.Scale(1/a);
444 return *this;
445 }
446
447 /**
448 product of a LorentzVector by a scalar quantity
449 \param a scalar quantity of type a
450 \return a new mathcoreLorentzVector q = v * a same type as v
451 */
453 LorentzVector tmp(*this);
454 tmp *= a;
455 return tmp;
456 }
457
458 /**
459 Divide a LorentzVector by a scalar quantity
460 \param a scalar quantity of type a
461 \return a new mathcoreLorentzVector q = v / a same type as v
462 */
465 tmp /= a;
466 return tmp;
467 }
468
469 /**
470 Negative of a LorentzVector (q = - v )
471 \return a new LorentzVector with opposite direction and time
472 */
474 //LorentzVector<CoordinateType> v(*this);
475 //v.Negate();
476 return operator*( Scalar(-1) );
477 }
479 return *this;
480 }
481
482 // ---- Relativistic Properties ----
483
484 /**
485 Rapidity relative to the Z axis: .5 log [(E+Pz)/(E-Pz)]
486 */
487 Scalar Rapidity() const {
488 // TODO - It would be good to check that E > Pz and use the Throw()
489 // mechanism or at least load a NAN if not.
490 // We should then move the code to a .cpp file.
491 const Scalar ee = E();
492 const Scalar ppz = Pz();
493 return Scalar(0.5) * log((ee + ppz) / (ee - ppz));
494 }
495
496 /**
497 Rapidity in the direction of travel: atanh (|P|/E)=.5 log[(E+P)/(E-P)]
498 */
500 // TODO - It would be good to check that E > P and use the Throw()
501 // mechanism or at least load a NAN if not.
502 const Scalar ee = E();
503 const Scalar pp = P();
504 return Scalar(0.5) * log((ee + pp) / (ee - pp));
505 }
506
507 /**
508 Determine if momentum-energy can represent a physical massive particle
509 */
510 bool isTimelike( ) const {
511 Scalar ee = E(); Scalar pp = P(); return ee*ee > pp*pp;
512 }
513
514 /**
515 Determine if momentum-energy can represent a massless particle
516 */
517 bool isLightlike( Scalar tolerance
519 Scalar ee = E(); Scalar pp = P(); Scalar delta = ee-pp;
520 if ( ee==0 ) return pp==0;
521 return delta*delta < tolerance * ee*ee;
522 }
523
524 /**
525 Determine if momentum-energy is spacelike, and represents a tachyon
526 */
527 bool isSpacelike( ) const {
528 Scalar ee = E(); Scalar pp = P(); return ee*ee < pp*pp;
529 }
530
532
533 /**
534 The beta vector for the boost that would bring this vector into
535 its center of mass frame (zero momentum)
536 */
538 if (E() == 0) {
539 if (P() == 0) {
540 return BetaVector();
541 } else {
542 // TODO - should attempt to Throw with msg about
543 // boostVector computed for LorentzVector with t=0
544 return -Vect()/E();
545 }
546 }
547 if (M2() <= 0) {
548 // TODO - should attempt to Throw with msg about
549 // boostVector computed for a non-timelike LorentzVector
550 }
551 return -Vect()/E();
552 }
553
554 /**
555 The beta vector for the boost that would bring this vector into
556 its center of mass frame (zero momentum)
557 */
558 template <class Other4Vector>
559 BetaVector BoostToCM(const Other4Vector& v ) const {
560 Scalar eSum = E() + v.E();
562 if (eSum == 0) {
563 if (vecSum.Mag2() == 0) {
564 return BetaVector();
565 } else {
566 // TODO - should attempt to Throw with msg about
567 // boostToCM computed for two 4-vectors with combined t=0
568 return BetaVector(vecSum/eSum);
569 }
570 // TODO - should attempt to Throw with msg about
571 // boostToCM computed for two 4-vectors with combined e=0
572 }
573 return BetaVector (vecSum * (-1./eSum));
574 }
575
576 //beta and gamma
577
578 /**
579 Return beta scalar value
580 */
581 Scalar Beta() const {
582 if ( E() == 0 ) {
583 if ( P2() == 0)
584 // to avoid Nan
585 return 0;
586 else {
587 GenVector::Throw ("LorentzVector::Beta() - beta computed for LorentzVector with t = 0. Return an Infinite result");
588 return 1./E();
589 }
590 }
591 if ( M2() <= 0 ) {
592 GenVector::Throw ("LorentzVector::Beta() - beta computed for non-timelike LorentzVector . Result is physically meaningless" );
593 }
594 return P() / E();
595 }
596 /**
597 Return Gamma scalar value
598 */
599 Scalar Gamma() const {
600 const Scalar v2 = P2();
601 const Scalar t2 = E() * E();
602 if (E() == 0) {
603 if ( P2() == 0) {
604 return 1;
605 } else {
606 GenVector::Throw ("LorentzVector::Gamma() - gamma computed for LorentzVector with t = 0. Return a zero result");
607
608 }
609 }
610 if ( t2 < v2 ) {
611 GenVector::Throw ("LorentzVector::Gamma() - gamma computed for a spacelike LorentzVector. Imaginary result");
612 return 0;
613 }
614 else if ( t2 == v2 ) {
615 GenVector::Throw ("LorentzVector::Gamma() - gamma computed for a lightlike LorentzVector. Infinite result");
616 }
617 return Scalar(1) / sqrt(Scalar(1) - v2 / t2);
618 } /* gamma */
619
620
621 // Method providing limited backward name compatibility with CLHEP ----
622
623 Scalar x() const { return fCoordinates.Px(); }
624 Scalar y() const { return fCoordinates.Py(); }
625 Scalar z() const { return fCoordinates.Pz(); }
626 Scalar t() const { return fCoordinates.E(); }
627 Scalar px() const { return fCoordinates.Px(); }
628 Scalar py() const { return fCoordinates.Py(); }
629 Scalar pz() const { return fCoordinates.Pz(); }
630 Scalar e() const { return fCoordinates.E(); }
631 Scalar r() const { return fCoordinates.R(); }
632 Scalar theta() const { return fCoordinates.Theta(); }
633 Scalar phi() const { return fCoordinates.Phi(); }
634 Scalar rho() const { return fCoordinates.Rho(); }
635 Scalar eta() const { return fCoordinates.Eta(); }
636 Scalar pt() const { return fCoordinates.Pt(); }
637 Scalar perp2() const { return fCoordinates.Perp2(); }
638 Scalar mag2() const { return fCoordinates.M2(); }
639 Scalar mag() const { return fCoordinates.M(); }
640 Scalar mt() const { return fCoordinates.Mt(); }
641 Scalar mt2() const { return fCoordinates.Mt2(); }
642
643
644 // Methods requested by CMS ---
645 Scalar energy() const { return fCoordinates.E(); }
646 Scalar mass() const { return fCoordinates.M(); }
647 Scalar mass2() const { return fCoordinates.M2(); }
648
649
650 /**
651 Methods setting a Single-component
652 Work only if the component is one of which the vector is represented.
653 For example SetE will work for a PxPyPzE Vector but not for a PxPyPzM Vector.
654 */
655 LorentzVector<CoordSystem>& SetE ( Scalar a ) { fCoordinates.SetE (a); return *this; }
657 LorentzVector<CoordSystem>& SetM ( Scalar a ) { fCoordinates.SetM (a); return *this; }
659 LorentzVector<CoordSystem>& SetPt ( Scalar a ) { fCoordinates.SetPt (a); return *this; }
660 LorentzVector<CoordSystem>& SetPx ( Scalar a ) { fCoordinates.SetPx (a); return *this; }
661 LorentzVector<CoordSystem>& SetPy ( Scalar a ) { fCoordinates.SetPy (a); return *this; }
662 LorentzVector<CoordSystem>& SetPz ( Scalar a ) { fCoordinates.SetPz (a); return *this; }
663
664 private:
665
666 CoordSystem fCoordinates; // internal coordinate system
667
668
669 }; // LorentzVector<>
670
671
672
673 // global nethods
674
675 /**
676 Scale of a LorentzVector with a scalar quantity a
677 \param a scalar quantity of typpe a
678 \param v mathcore::LorentzVector based on any coordinate system
679 \return a new mathcoreLorentzVector q = v * a same type as v
680 */
681 template< class CoordSystem >
683 ( const typename LorentzVector<CoordSystem>::Scalar & a,
686 tmp *= a;
687 return tmp;
688 }
689
690 // ------------- I/O to/from streams -------------
691
692 template< class char_t, class traits_t, class Coords >
693 inline
694 std::basic_ostream<char_t,traits_t> &
695 operator << ( std::basic_ostream<char_t,traits_t> & os
696 , LorentzVector<Coords> const & v
697 )
698 {
699 if( !os ) return os;
700
701 typename Coords::Scalar a, b, c, d;
702 v.GetCoordinates(a, b, c, d);
703
706 // TODO: call MF's bitwise-accurate functions on each of a, b, c, d
707 }
708 else {
709 os << detail::get_manip( os, detail::open ) << a
710 << detail::get_manip( os, detail::sep ) << b
711 << detail::get_manip( os, detail::sep ) << c
712 << detail::get_manip( os, detail::sep ) << d
714 }
715
716 return os;
717
718 } // op<< <>()
719
720
721 template< class char_t, class traits_t, class Coords >
722 inline
723 std::basic_istream<char_t,traits_t> &
724 operator >> ( std::basic_istream<char_t,traits_t> & is
726 )
727 {
728 if( !is ) return is;
729
730 typename Coords::Scalar a, b, c, d;
731
734 // TODO: call MF's bitwise-accurate functions on each of a, b, c
735 }
736 else {
738 detail::require_delim( is, detail::sep ); is >> b;
739 detail::require_delim( is, detail::sep ); is >> c;
740 detail::require_delim( is, detail::sep ); is >> d;
742 }
743
744 if( is )
745 v.SetCoordinates(a, b, c, d);
746 return is;
747
748 } // op>> <>()
749
750
751
752 } // end namespace Math
753
754} // end namespace ROOT
755
756#include <sstream>
757namespace cling
758{
759template<typename CoordSystem>
760std::string printValue(const ROOT::Math::LorentzVector<CoordSystem> *v)
761{
762 std::stringstream s;
763 s << *v;
764 return s.str();
765}
766
767} // end namespace cling
768
769#endif
770
771//#include "Math/GenVector/LorentzVectorOperations.h"
772
773
774
#define d(i)
Definition: RSha256.hxx:102
#define b(i)
Definition: RSha256.hxx:100
#define c(i)
Definition: RSha256.hxx:101
float * q
Definition: THbookFile.cxx:87
double log(double)
typedef void((*Func_t)())
Class describing a generic displacement vector in 3 dimensions.
Scalar Mag2() const
Magnitute squared ( r^2 in spherical coordinate)
Scalar E() const
return 4-th component (time, or energy for a 4-momentum vector)
Scalar Et() const
return the transverse energy
LorentzVector< CoordSystem > & SetPz(Scalar a)
void GetCoordinates(Scalar dest[]) const
get internal data into an array of 4 Scalar numbers
BetaVector BoostToCM() const
The beta vector for the boost that would bring this vector into its center of mass frame (zero moment...
bool operator==(const LorentzVector &rhs) const
Exact equality.
Scalar M() const
return magnitude (mass) using the (-,-,-,+) metric.
LorentzVector< CoordSystem > & SetPy(Scalar a)
Scalar Pt() const
return the transverse spatial component sqrt ( X**2 + Y**2 )
Scalar Mt2() const
return the transverse mass squared
Scalar Rapidity() const
Rapidity relative to the Z axis: .5 log [(E+Pz)/(E-Pz)].
Scalar Beta() const
Return beta scalar value.
LorentzVector< CoordSystem > & SetPx(Scalar a)
Scalar Perp2() const
return the square of the transverse spatial component ( X**2 + Y**2 )
void GetCoordinates(IT begin) const
get internal data into 4 Scalars at *begin
LorentzVector< CoordSystem > operator/(const Scalar &a) const
Divide a LorentzVector by a scalar quantity.
LorentzVector< CoordSystem > & SetPxPyPzE(Scalar xx, Scalar yy, Scalar zz, Scalar ee)
LorentzVector & operator+=(const OtherLorentzVector &q)
Self addition with another Vector ( v+= q ) Enable the addition with any other LorentzVector.
Scalar Eta() const
pseudorapidity
Scalar Py() const
spatial Y component
LorentzVector(const Scalar &a, const Scalar &b, const Scalar &c, const Scalar &d)
generic constructors from four scalar values.
Definition: LorentzVector.h:82
bool isLightlike(Scalar tolerance=100 *std::numeric_limits< Scalar >::epsilon()) const
Determine if momentum-energy can represent a massless particle.
::ROOT::Math::DisplacementVector3D< Cartesian3D< Scalar > > Vect() const
get the spatial components of the Vector in a DisplacementVector based on Cartesian Coordinates
Scalar ColinearRapidity() const
Rapidity in the direction of travel: atanh (|P|/E)=.5 log[(E+P)/(E-P)].
Scalar Pz() const
spatial Z component
LorentzVector(const ForeignLorentzVector &v)
Construct from a foreign 4D vector type, for example, HepLorentzVector Precondition: v must implement...
LorentzVector operator-() const
Negative of a LorentzVector (q = - v )
LorentzVector< CoordSystem > & SetPt(Scalar a)
LorentzVector operator+() const
const CoordSystem & Coordinates() const
Retrieve a const reference to the coordinates object.
LorentzVector< CoordSystem > & SetEta(Scalar a)
LorentzVector operator*(const Scalar &a) const
product of a LorentzVector by a scalar quantity
Scalar Phi() const
azimuthal Angle
LorentzVector< CoordSystem > & SetCoordinates(const Scalar src[])
Set internal data based on an array of 4 Scalar numbers.
LorentzVector & operator/=(Scalar a)
division by a scalar quantity v /= a
LorentzVector & operator*=(Scalar a)
multiplication by a scalar quantity v *= a
Scalar M2() const
return magnitude (mass) squared M2 = T**2 - X**2 - Y**2 - Z**2 (we use -,-,-,+ metric)
Scalar Dot(const OtherLorentzVector &q) const
scalar (Dot) product of two LorentzVector vectors (metric is -,-,-,+) Enable the product using any ot...
void GetCoordinates(Scalar &a, Scalar &b, Scalar &c, Scalar &d) const
get internal data into 4 Scalar numbers
Scalar Gamma() const
Return Gamma scalar value.
bool isSpacelike() const
Determine if momentum-energy is spacelike, and represents a tachyon.
Scalar R() const
return the spatial (3D) magnitude ( sqrt(X**2 + Y**2 + Z**2) )
CoordSystem::Scalar Scalar
Definition: LorentzVector.h:65
LorentzVector(const LorentzVector< Coords > &v)
constructor from a LorentzVector expressed in different coordinates, or using a different Scalar type
Definition: LorentzVector.h:93
Scalar P2() const
return the square of the spatial (3D) magnitude ( X**2 + Y**2 + Z**2 )
Scalar Et2() const
return the transverse energy squared
LorentzVector()
default constructor of an empty vector (Px = Py = Pz = E = 0 )
Definition: LorentzVector.h:71
bool isTimelike() const
Determine if momentum-energy can represent a physical massive particle.
LorentzVector & operator-=(const OtherLorentzVector &q)
Self subtraction of another Vector from this ( v-= q ) Enable the addition with any other LorentzVect...
Scalar Px() const
spatial X component
Scalar Mt() const
return the transverse mass
LorentzVector< CoordSystem > & SetM(Scalar a)
LorentzVector< CoordSystem > & SetE(Scalar a)
Methods setting a Single-component Work only if the component is one of which the vector is represent...
LorentzVector< CoordSystem > & SetCoordinates(IT begin, IT end)
Set internal data based on 4 Scalars at *begin to *end.
LorentzVector< CoordSystem > & SetXYZT(Scalar xx, Scalar yy, Scalar zz, Scalar tt)
set the values of the vector from the cartesian components (x,y,z,t) (if the vector is held in anothe...
DisplacementVector3D< Cartesian3D< Scalar > > BetaVector
bool operator!=(const LorentzVector &rhs) const
void GetCoordinates(IT begin, IT end) const
get internal data into 4 Scalars at *begin to *end
Scalar Theta() const
polar Angle
BetaVector BoostToCM(const Other4Vector &v) const
The beta vector for the boost that would bring this vector into its center of mass frame (zero moment...
LorentzVector & operator=(const LorentzVector< OtherCoords > &v)
Assignment operator from a lorentz vector of arbitrary type.
LorentzVector< CoordSystem > & SetCoordinates(Scalar a, Scalar b, Scalar c, Scalar d)
Set internal data based on 4 Scalar numbers.
LorentzVector< CoordSystem > & SetPhi(Scalar a)
Class describing a 4D cartesian coordinate system (x, y, z, t coordinates) or momentum-energy vectors...
Definition: PxPyPzE4D.h:42
Namespace for new Math classes and functions.
void Throw(const char *)
function throwing exception, by creating internally a GenVector_exception only when needed
char_t get_manip(std::basic_ios< char_t, traits_t > &ios, manip_t m)
Definition: GenVectorIO.h:54
std::basic_istream< char_t, traits_t > & require_delim(std::basic_istream< char_t, traits_t > &is, manip_t m)
Definition: GenVectorIO.h:113
void set_manip(std::basic_ios< char_t, traits_t > &ios, manip_t m, char_t ch)
Definition: GenVectorIO.h:74
std::ostream & operator<<(std::ostream &os, const AxisAngle &a)
Stream Output and Input.
Definition: AxisAngle.cxx:91
std::basic_istream< char_t, traits_t > & operator>>(std::basic_istream< char_t, traits_t > &is, DisplacementVector2D< T, U > &v)
VecExpr< UnaryOp< Sqrt< T >, VecExpr< A, T, D >, T >, T, D > sqrt(const VecExpr< A, T, D > &rhs)
Rotation3D::Scalar Scalar
VSD Structures.
Definition: StringConv.hxx:21
static constexpr double s
auto * tt
Definition: textangle.C:16
auto * a
Definition: textangle.C:12
REAL epsilon
Definition: triangle.c:617
#define dest(otri, vertexptr)
Definition: triangle.c:1040