Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
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
24#include "Math/Vector2D.h"
25
26#include "TMath.h"
27
28#include <cmath>
29#include <string>
30
31namespace ROOT {
32
33 namespace Math {
34
35//__________________________________________________________________________________________
36/** @ingroup GenVector
37
38Class describing a generic LorentzVector in the 4D space-time,
39using the specified coordinate system for the spatial vector part.
40The metric used for the LorentzVector is (-,-,-,+).
41In the case of LorentzVector we don't distinguish the concepts
42of points and displacement vectors as in the 3D case,
43since the main use case for 4D Vectors is to describe the kinematics of
44relativistic particles. A LorentzVector behaves like a
45DisplacementVector in 4D. The Minkowski components could be viewed as
46v and t, or for kinematic 4-vectors, as p and E.
47
48ROOT provides specialisations (and aliases to them) of the ROOT::Math::LorentzVector template:
49- ROOT::Math::PtEtaPhiMVector based on pt (rho),eta,phi and M (t) coordinates in double precision
50- ROOT::Math::PtEtaPhiEVector based on pt (rho),eta,phi and E (t) coordinates in double precision
51- ROOT::Math::PxPyPzMVector based on px,py,pz and M coordinates in double precision
52- ROOT::Math::PxPyPzEVector based on px,py,pz and E coordinates in double precision
53- ROOT::Math::XYZTVector based on x,y,z,t coordinates (cartesian) in double precision (same as PxPyPzEVector)
54- ROOT::Math::XYZTVectorF based on x,y,z,t coordinates (cartesian) in float precision (same as PxPyPzEVector but float)
55Pt (or rho) refers to transverse momentum, whereas eta refers to pseudorapidity. M is mass, E is energy, p is momentum.
56
57@see GenVector
58*/
59
60 template< class CoordSystem >
62
63 public:
64
65 // ------ ctors ------
66
67 typedef typename CoordSystem::Scalar Scalar;
68 typedef CoordSystem CoordinateType;
69
70 /**
71 default constructor of an empty vector (Px = Py = Pz = E = 0 )
72 */
74
75 /**
76 generic constructors from four scalar values.
77 The association between values and coordinate depends on the
78 coordinate system. For PxPyPzE4D,
79 \param a scalar value (Px)
80 \param b scalar value (Py)
81 \param c scalar value (Pz)
82 \param d scalar value (E)
83 */
85 const Scalar & b,
86 const Scalar & c,
87 const Scalar & d) :
88 fCoordinates(a , b, c, d) { }
89
90 /**
91 constructor from a LorentzVector expressed in different
92 coordinates, or using a different Scalar type
93 */
94 template< class Coords >
95 explicit constexpr LorentzVector(const LorentzVector<Coords> & v ) :
97
98 /**
99 Construct from a foreign 4D vector type, for example, HepLorentzVector
100 Precondition: v must implement methods x(), y(), z(), and t()
101 */
102 template<class ForeignLorentzVector,
103 typename = decltype(std::declval<ForeignLorentzVector>().x()
104 + std::declval<ForeignLorentzVector>().y()
105 + std::declval<ForeignLorentzVector>().z()
106 + std::declval<ForeignLorentzVector>().t())>
107 explicit constexpr LorentzVector( const ForeignLorentzVector & v) :
108 fCoordinates(PxPyPzE4D<Scalar>( v.x(), v.y(), v.z(), v.t() ) ) { }
109
110#ifdef LATER
111 /**
112 construct from a generic linear algebra vector implementing operator []
113 and with a size of at least 4. This could be also a C array
114 In this case v[0] is the first data member
115 ( Px for a PxPyPzE4D base)
116 \param v LA vector
117 \param index0 index of first vector element (Px)
118 */
119 template< class LAVector >
120 explicit constexpr LorentzVector(const LAVector & v, size_t index0 ) {
121 fCoordinates = CoordSystem ( v[index0], v[index0+1], v[index0+2], v[index0+3] );
122 }
123#endif
124
125
126 // ------ assignment ------
127
128 /**
129 Assignment operator from a lorentz vector of arbitrary type
130 */
131 template< class OtherCoords >
133 fCoordinates = v.Coordinates();
134 return *this;
135 }
136
137 /**
138 assignment from any other Lorentz vector implementing
139 x(), y(), z() and t()
140 */
141 template<class ForeignLorentzVector,
142 typename = decltype(std::declval<ForeignLorentzVector>().x()
143 + std::declval<ForeignLorentzVector>().y()
144 + std::declval<ForeignLorentzVector>().z()
145 + std::declval<ForeignLorentzVector>().t())>
146 LorentzVector & operator = ( const ForeignLorentzVector & v) {
147 SetXYZT( v.x(), v.y(), v.z(), v.t() );
148 return *this;
149 }
150
151#ifdef LATER
152 /**
153 assign from a generic linear algebra vector implementing operator []
154 and with a size of at least 4
155 In this case v[0] is the first data member
156 ( Px for a PxPyPzE4D base)
157 \param v LA vector
158 \param index0 index of first vector element (Px)
159 */
160 template< class LAVector >
161 LorentzVector & AssignFrom(const LAVector & v, size_t index0=0 ) {
162 fCoordinates.SetCoordinates( v[index0], v[index0+1], v[index0+2], v[index0+3] );
163 return *this;
164 }
165#endif
166
167 // ------ Set, Get, and access coordinate data ------
168
169 /**
170 Retrieve a const reference to the coordinates object
171 */
172 const CoordSystem & Coordinates() const {
173 return fCoordinates;
174 }
175
176 /**
177 Set internal data based on an array of 4 Scalar numbers
178 */
180 fCoordinates.SetCoordinates(src);
181 return *this;
182 }
183
184 /**
185 Set internal data based on 4 Scalar numbers
186 */
188 fCoordinates.SetCoordinates(a, b, c, d);
189 return *this;
190 }
191
192 /**
193 Set internal data based on 4 Scalars at *begin to *end
194 */
195 template< class IT >
197 IT a = begin; IT b = ++begin; IT c = ++begin; IT d = ++begin;
198 (void)end;
199 assert (++begin==end);
200 SetCoordinates (*a,*b,*c,*d);
201 return *this;
202 }
203
204 /**
205 get internal data into 4 Scalar numbers
206 \note Alternatively, you may use structured bindings: `auto const [a, b, c, d] = v`.
207 */
208 void GetCoordinates( Scalar& a, Scalar& b, Scalar& c, Scalar & d ) const
209 { fCoordinates.GetCoordinates(a, b, c, d); }
210
211 /**
212 get internal data into an array of 4 Scalar numbers
213 */
214 void GetCoordinates( Scalar dest[] ) const
215 { fCoordinates.GetCoordinates(dest); }
216
217 /**
218 get internal data into 4 Scalars at *begin to *end
219 */
220 template <class IT>
221 void GetCoordinates( IT begin, IT end ) const
222 { IT a = begin; IT b = ++begin; IT c = ++begin; IT d = ++begin;
223 (void)end;
224 assert (++begin==end);
225 GetCoordinates (*a,*b,*c,*d);
226 }
227
228 /**
229 get internal data into 4 Scalars at *begin
230 */
231 template <class IT>
232 void GetCoordinates( IT begin ) const {
233 Scalar a,b,c,d = 0;
234 GetCoordinates (a,b,c,d);
235 *begin++ = a;
236 *begin++ = b;
237 *begin++ = c;
238 *begin = d;
239 }
240
241 /**
242 set the values of the vector from the cartesian components (x,y,z,t)
243 (if the vector is held in another coordinates, like (Pt,eta,phi,m)
244 then (x, y, z, t) are converted to that form)
245 */
247 fCoordinates.SetPxPyPzE(xx,yy,zz,tt);
248 return *this;
249 }
251 fCoordinates.SetPxPyPzE(xx,yy,zz,ee);
252 return *this;
253 }
254
255 // ------------------- Equality -----------------
256
257 /**
258 Exact equality
259 */
260 bool operator==(const LorentzVector & rhs) const {
261 return fCoordinates==rhs.fCoordinates;
262 }
263 bool operator!= (const LorentzVector & rhs) const {
264 return !(operator==(rhs));
265 }
266
267 // ------ Individual element access, in various coordinate systems ------
268
269 /**
270 dimension
271 */
272 unsigned int Dimension() const
273 {
274 return fDimension;
275 };
276
277 // individual coordinate accessors in various coordinate systems
278
279 /**
280 spatial X component
281 */
282 Scalar Px() const { return fCoordinates.Px(); }
283 Scalar X() const { return fCoordinates.Px(); }
284 /**
285 spatial Y component
286 */
287 Scalar Py() const { return fCoordinates.Py(); }
288 Scalar Y() const { return fCoordinates.Py(); }
289 /**
290 spatial Z component
291 */
292 Scalar Pz() const { return fCoordinates.Pz(); }
293 Scalar Z() const { return fCoordinates.Pz(); }
294 /**
295 return 4-th component (time, or energy for a 4-momentum vector)
296 */
297 Scalar E() const { return fCoordinates.E(); }
298 Scalar T() const { return fCoordinates.E(); }
299 /**
300 return magnitude (mass) squared M2 = T**2 - X**2 - Y**2 - Z**2
301 (we use -,-,-,+ metric)
302 */
303 Scalar M2() const { return fCoordinates.M2(); }
304 /**
305 return magnitude (mass) using the (-,-,-,+) metric.
306 If M2 is negative (space-like vector) a GenVector_exception
307 is suggested and if continuing, - sqrt( -M2) is returned
308 */
309 Scalar M() const { return fCoordinates.M();}
310 /**
311 return the spatial (3D) magnitude ( sqrt(X**2 + Y**2 + Z**2) )
312 */
313 Scalar R() const { return fCoordinates.R(); }
314 Scalar P() const { return fCoordinates.R(); }
315 /**
316 return the square of the spatial (3D) magnitude ( X**2 + Y**2 + Z**2 )
317 */
318 Scalar P2() const { return P() * P(); }
319 /**
320 return the square of the transverse spatial component ( X**2 + Y**2 )
321 */
322 Scalar Perp2( ) const { return fCoordinates.Perp2();}
323
324 /**
325 return the transverse spatial component sqrt ( X**2 + Y**2 )
326 */
327 Scalar Pt() const { return fCoordinates.Pt(); }
328 Scalar Rho() const { return fCoordinates.Pt(); }
329
330 /**
331 return the transverse mass squared
332 \f[ m_t^2 = E^2 - p{_z}^2 \f]
333 */
334 Scalar Mt2() const { return fCoordinates.Mt2(); }
335
336 /**
337 return the transverse mass
338 \f[ \sqrt{ m_t^2 = E^2 - p{_z}^2} X sign(E^ - p{_z}^2) \f]
339 */
340 Scalar Mt() const { return fCoordinates.Mt(); }
341
342 /**
343 return the transverse energy squared
344 \f[ e_t = \frac{E^2 p_{\perp}^2 }{ |p|^2 } \f]
345 */
346 Scalar Et2() const { return fCoordinates.Et2(); }
347
348 /**
349 return the transverse energy
350 \f[ e_t = \sqrt{ \frac{E^2 p_{\perp}^2 }{ |p|^2 } } X sign(E) \f]
351 */
352 Scalar Et() const { return fCoordinates.Et(); }
353
354 /**
355 azimuthal Angle
356 */
357 Scalar Phi() const { return fCoordinates.Phi();}
358
359 /**
360 polar Angle
361 */
362 Scalar Theta() const { return fCoordinates.Theta(); }
363
364 /**
365 pseudorapidity
366 \f[ \eta = - \ln { \tan { \frac { \theta} {2} } } \f]
367 */
368 Scalar Eta() const { return fCoordinates.Eta(); }
369
370 /**
371 deltaRapidity between this and vector v
372 \f[ \Delta R = \sqrt { \Delta \eta ^2 + \Delta \phi ^2 } \f]
373 \param useRapidity true to use Rapidity(), false to use Eta()
374 */
375 template <class OtherLorentzVector>
376 Scalar DeltaR(const OtherLorentzVector &v, const bool useRapidity = false) const
377 {
378 const double delta = useRapidity ? Rapidity() - v.Rapidity() : Eta() - v.Eta();
379 double dphi = Phi() - v.Phi();
380 // convert dphi angle to the interval (-PI,PI]
381 if (dphi > TMath::Pi() || dphi <= -TMath::Pi()) {
382 if (dphi > 0) {
383 int n = static_cast<int>(dphi / TMath::TwoPi() + 0.5);
384 dphi -= TMath::TwoPi() * n;
385 } else {
386 int n = static_cast<int>(0.5 - dphi / TMath::TwoPi());
387 dphi += TMath::TwoPi() * n;
388 }
389 }
390 return std::sqrt(delta * delta + dphi * dphi);
391 }
392
393 /**
394 get the spatial components of the Vector in a
395 DisplacementVector based on Cartesian Coordinates
396 */
398 return ::ROOT::Math::DisplacementVector3D<Cartesian3D<Scalar> >( X(), Y(), Z() );
399 }
400
401 // ------ Operations combining two Lorentz vectors ------
402
403 /**
404 scalar (Dot) product of two LorentzVector vectors (metric is -,-,-,+)
405 Enable the product using any other LorentzVector implementing
406 the x(), y() , y() and t() member functions
407 \param q any LorentzVector implementing the x(), y() , z() and t()
408 member functions
409 \return the result of v.q of type according to the base scalar type of v
410 */
411
412 template< class OtherLorentzVector >
414 return t()*q.t() - x()*q.x() - y()*q.y() - z()*q.z();
415 }
416
417 /**
418 Self addition with another Vector ( v+= q )
419 Enable the addition with any other LorentzVector
420 \param q any LorentzVector implementing the x(), y() , z() and t()
421 member functions
422 */
423 template< class OtherLorentzVector >
425 {
426 SetXYZT( x() + q.x(), y() + q.y(), z() + q.z(), t() + q.t() );
427 return *this;
428 }
429
430 /**
431 Self subtraction of another Vector from this ( v-= q )
432 Enable the addition with any other LorentzVector
433 \param q any LorentzVector implementing the x(), y() , z() and t()
434 member functions
435 */
436 template< class OtherLorentzVector >
438 SetXYZT( x() - q.x(), y() - q.y(), z() - q.z(), t() - q.t() );
439 return *this;
440 }
441
442 /**
443 addition of two LorentzVectors (v3 = v1 + v2)
444 Enable the addition with any other LorentzVector
445 \param v2 any LorentzVector implementing the x(), y() , z() and t()
446 member functions
447 \return a new LorentzVector of the same type as v1
448 */
449 template<class OtherLorentzVector>
451 {
453 v3 += v2;
454 return v3;
455 }
456
457 /**
458 subtraction of two LorentzVectors (v3 = v1 - v2)
459 Enable the subtraction of any other LorentzVector
460 \param v2 any LorentzVector implementing the x(), y() , z() and t()
461 member functions
462 \return a new LorentzVector of the same type as v1
463 */
464 template<class OtherLorentzVector>
467 v3 -= v2;
468 return v3;
469 }
470
471 //--- scaling operations ------
472
473 /**
474 multiplication by a scalar quantity v *= a
475 */
477 fCoordinates.Scale(a);
478 return *this;
479 }
480
481 /**
482 division by a scalar quantity v /= a
483 */
485 fCoordinates.Scale(1/a);
486 return *this;
487 }
488
489 /**
490 product of a LorentzVector by a scalar quantity
491 \param a scalar quantity of type a
492 \return a new mathcoreLorentzVector q = v * a same type as v
493 */
495 LorentzVector tmp(*this);
496 tmp *= a;
497 return tmp;
498 }
499
500 /**
501 Divide a LorentzVector by a scalar quantity
502 \param a scalar quantity of type a
503 \return a new mathcoreLorentzVector q = v / a same type as v
504 */
507 tmp /= a;
508 return tmp;
509 }
510
511 /**
512 Negative of a LorentzVector (q = - v )
513 \return a new LorentzVector with opposite direction and time
514 */
516 //LorentzVector<CoordinateType> v(*this);
517 //v.Negate();
518 return operator*( Scalar(-1) );
519 }
521 return *this;
522 }
523
524 // ---- Relativistic Properties ----
525
526 /**
527 Rapidity relative to the Z axis: .5 log [(E+Pz)/(E-Pz)]
528 */
529 Scalar Rapidity() const {
530 // TODO - It would be good to check that E > Pz and use the Throw()
531 // mechanism or at least load a NAN if not.
532 // We should then move the code to a .cpp file.
533 const Scalar ee = E();
534 const Scalar ppz = Pz();
535 using std::log;
536 return Scalar(0.5) * log((ee + ppz) / (ee - ppz));
537 }
538
539 /**
540 Rapidity in the direction of travel: atanh (|P|/E)=.5 log[(E+P)/(E-P)]
541 */
543 // TODO - It would be good to check that E > P and use the Throw()
544 // mechanism or at least load a NAN if not.
545 const Scalar ee = E();
546 const Scalar pp = P();
547 using std::log;
548 return Scalar(0.5) * log((ee + pp) / (ee - pp));
549 }
550
551 /**
552 Determine if momentum-energy can represent a physical massive particle
553 */
554 bool isTimelike( ) const {
555 Scalar ee = E(); Scalar pp = P(); return ee*ee > pp*pp;
556 }
557
558 /**
559 Determine if momentum-energy can represent a massless particle
560 */
562 = 100*std::numeric_limits<Scalar>::epsilon() ) const {
563 Scalar ee = E(); Scalar pp = P(); Scalar delta = ee-pp;
564 if ( ee==0 ) return pp==0;
565 return delta*delta < tolerance * ee*ee;
566 }
567
568 /**
569 Determine if momentum-energy is spacelike, and represents a tachyon
570 */
571 bool isSpacelike( ) const {
572 Scalar ee = E(); Scalar pp = P(); return ee*ee < pp*pp;
573 }
574
576
577 /**
578 The beta vector for the boost that would bring this vector into
579 its center of mass frame (zero momentum)
580 */
582 if (E() == 0) {
583 if (P() == 0) {
584 return BetaVector();
585 } else {
586 // TODO - should attempt to Throw with msg about
587 // boostVector computed for LorentzVector with t=0
588 return -Vect()/E();
589 }
590 }
591 if (M2() <= 0) {
592 // TODO - should attempt to Throw with msg about
593 // boostVector computed for a non-timelike LorentzVector
594 }
595 return -Vect()/E();
596 }
597
598 /**
599 The beta vector for the boost that would bring this vector into
600 its center of mass frame (zero momentum)
601 */
602 template <class Other4Vector>
604 Scalar eSum = E() + v.E();
606 if (eSum == 0) {
607 if (vecSum.Mag2() == 0) {
608 return BetaVector();
609 } else {
610 // TODO - should attempt to Throw with msg about
611 // boostToCM computed for two 4-vectors with combined t=0
612 return BetaVector(vecSum/eSum);
613 }
614 // TODO - should attempt to Throw with msg about
615 // boostToCM computed for two 4-vectors with combined e=0
616 }
617 return BetaVector (vecSum * (-1./eSum));
618 }
619
620 //beta and gamma
621
622 /**
623 Return beta scalar value
624 */
625 Scalar Beta() const {
626 if ( E() == 0 ) {
627 if ( P2() == 0)
628 // to avoid Nan
629 return 0;
630 else {
632 "LorentzVector::Beta() - beta computed for LorentzVector with t = 0. Return an Infinite result");
633 return 1./E();
634 }
635 }
636 if ( M2() <= 0 ) {
637 GenVector_Throw("LorentzVector::Beta() - beta computed for non-timelike LorentzVector . Result is "
638 "physically meaningless");
639 }
640 return P() / E();
641 }
642 /**
643 Return Gamma scalar value
644 */
645 Scalar Gamma() const {
646 const Scalar v2 = P2();
647 const Scalar t2 = E() * E();
648 if (E() == 0) {
649 if ( P2() == 0) {
650 return 1;
651 } else {
653 "LorentzVector::Gamma() - gamma computed for LorentzVector with t = 0. Return a zero result");
654 }
655 }
656 if ( t2 < v2 ) {
657 GenVector_Throw("LorentzVector::Gamma() - gamma computed for a spacelike LorentzVector. Imaginary result");
658 return 0;
659 }
660 else if ( t2 == v2 ) {
661 GenVector_Throw("LorentzVector::Gamma() - gamma computed for a lightlike LorentzVector. Infinite result");
662 }
663 using std::sqrt;
664 return Scalar(1) / sqrt(Scalar(1) - v2 / t2);
665 } /* gamma */
666
667
668 // Method providing limited backward name compatibility with CLHEP ----
669
670 Scalar x() const { return fCoordinates.Px(); }
671 Scalar y() const { return fCoordinates.Py(); }
672 Scalar z() const { return fCoordinates.Pz(); }
673 Scalar t() const { return fCoordinates.E(); }
674 Scalar px() const { return fCoordinates.Px(); }
675 Scalar py() const { return fCoordinates.Py(); }
676 Scalar pz() const { return fCoordinates.Pz(); }
677 Scalar e() const { return fCoordinates.E(); }
678 Scalar r() const { return fCoordinates.R(); }
679 Scalar theta() const { return fCoordinates.Theta(); }
680 Scalar phi() const { return fCoordinates.Phi(); }
681 Scalar rho() const { return fCoordinates.Rho(); }
682 Scalar eta() const { return fCoordinates.Eta(); }
683 Scalar pt() const { return fCoordinates.Pt(); }
684 Scalar perp2() const { return fCoordinates.Perp2(); }
685 Scalar mag2() const { return fCoordinates.M2(); }
686 Scalar mag() const { return fCoordinates.M(); }
687 Scalar mt() const { return fCoordinates.Mt(); }
688 Scalar mt2() const { return fCoordinates.Mt2(); }
689
690
691 // Methods requested by CMS ---
692 Scalar energy() const { return fCoordinates.E(); }
693 Scalar mass() const { return fCoordinates.M(); }
694 Scalar mass2() const { return fCoordinates.M2(); }
695
696
697 /**
698 Methods setting a Single-component
699 Work only if the component is one of which the vector is represented.
700 For example SetE will work for a PxPyPzE Vector but not for a PxPyPzM Vector.
701 */
702 LorentzVector<CoordSystem>& SetE ( Scalar a ) { fCoordinates.SetE (a); return *this; }
704 LorentzVector<CoordSystem>& SetM ( Scalar a ) { fCoordinates.SetM (a); return *this; }
706 LorentzVector<CoordSystem>& SetPt ( Scalar a ) { fCoordinates.SetPt (a); return *this; }
707 LorentzVector<CoordSystem>& SetPx ( Scalar a ) { fCoordinates.SetPx (a); return *this; }
708 LorentzVector<CoordSystem>& SetPy ( Scalar a ) { fCoordinates.SetPy (a); return *this; }
709 LorentzVector<CoordSystem>& SetPz ( Scalar a ) { fCoordinates.SetPz (a); return *this; }
710
711 private:
712
713 CoordSystem fCoordinates; // internal coordinate system
714 static constexpr unsigned int fDimension = CoordinateType::Dimension;
715
716 }; // LorentzVector<>
717
718
719
720 // global methods
721
722 /**
723 Scale of a LorentzVector with a scalar quantity a
724 \param a scalar quantity of type a
725 \param v LorentzVector based on any coordinate system
726 \return a new mathcoreLorentzVector q = v * a same type as v
727 */
728 template <class CoordSystem>
731 {
733 tmp *= a;
734 return tmp;
735 }
736
737 /**
738 pair (p+ p-) acoplanarity `alpha = 1 - |phi+ - phi-|/pi`.
739 \param pp p+, LorentzVector based on any coordinate system
740 \param pm p-, LorentzVector based on any coordinate system
741 \return a scalar
742 \ingroup GenVector
743 \see http://doi.org/10.1103/PhysRevLett.121.212301
744 */
745 template <class CoordSystem>
748 {
749 auto dphi = pp.Phi() - pm.Phi();
750 // convert dphi angle to the interval (-PI,PI]
751 if (dphi > TMath::Pi() || dphi <= -TMath::Pi()) {
752 if (dphi > 0) {
753 int n = static_cast<int>(dphi / TMath::TwoPi() + 0.5);
754 dphi -= TMath::TwoPi() * n;
755 } else {
756 int n = static_cast<int>(0.5 - dphi / TMath::TwoPi());
757 dphi += TMath::TwoPi() * n;
758 }
759 }
760 return 1 - std::abs(dphi) / TMath::Pi();
761 }
762
763 /**
764 pair (p+ p-) vectorial asymmetry `Av = |Pt+ - Pt-|/|Pt+ + Pt-|`.
765 In an experimental setting, it reflects a convolution of the experimental resolutions
766 of particle energy and azimuthal angle measurement.
767 \param pp p+, LorentzVector based on any coordinate system
768 \param pm p-, LorentzVector based on any coordinate system
769 \return a scalar. Returns -1 if both momenta are exactly mirrored.
770 \ingroup GenVector
771 \see http://doi.org/10.1103/PhysRevLett.121.212301, https://doi.org/10.1103/PhysRevD.99.093013
772 */
773 template <class CoordSystem>
776 {
777 ROOT::Math::XYVector vp(pp.Px(), pp.Py());
778 ROOT::Math::XYVector vm(pm.Px(), pm.Py());
779 auto denom = (vp + vm).R();
780 if (denom == 0.)
781 return -1;
782 return (vp - vm).R() / denom;
783 }
784
785 /**
786 pair (p+ p-) scalar asymmetry `As = ||Pt+| - |Pt-|/||Pt+| + |Pt-||`.
787 Measures the relative difference in transverse momentum of the pair, e.g. two photons,
788 and would be ideally zero for two back-to-back photons.
789 \param pp p+, LorentzVector based on any coordinate system
790 \param pm p-, LorentzVector based on any coordinate system
791 \return a scalar. Returns 0 if both transverse momenta are zero
792 \ingroup GenVector
793 \see http://doi.org/10.1103/PhysRevLett.121.212301, https://doi.org/10.1103/PhysRevD.99.093013
794 */
795 template< class CoordSystem >
797 {
798 auto ptp = pp.Pt();
799 auto ptm = pm.Pt();
800 if (ptp == 0 && ptm == 0)
801 return 0;
802 return std::abs(ptp - ptm) / (ptp + ptm);
803 }
804
805 // ------------- I/O to/from streams -------------
806
807 template< class char_t, class traits_t, class Coords >
808 inline
809 std::basic_ostream<char_t,traits_t> &
810 operator << ( std::basic_ostream<char_t,traits_t> & os
811 , LorentzVector<Coords> const & v
812 )
813 {
814 if( !os ) return os;
815
816 typename Coords::Scalar a, b, c, d;
817 v.GetCoordinates(a, b, c, d);
818
821 // TODO: call MF's bitwise-accurate functions on each of a, b, c, d
822 }
823 else {
824 os << detail::get_manip( os, detail::open ) << a
825 << detail::get_manip( os, detail::sep ) << b
826 << detail::get_manip( os, detail::sep ) << c
827 << detail::get_manip( os, detail::sep ) << d
829 }
830
831 return os;
832
833 } // op<< <>()
834
835
836 template< class char_t, class traits_t, class Coords >
837 inline
838 std::basic_istream<char_t,traits_t> &
839 operator >> ( std::basic_istream<char_t,traits_t> & is
841 )
842 {
843 if( !is ) return is;
844
845 typename Coords::Scalar a, b, c, d;
846
849 // TODO: call MF's bitwise-accurate functions on each of a, b, c
850 }
851 else {
857 }
858
859 if( is )
860 v.SetCoordinates(a, b, c, d);
861 return is;
862
863 } // op>> <>()
864
865 // Structured bindings
866 template <std::size_t I, class CoordSystem>
867 typename CoordSystem::Scalar get(LorentzVector<CoordSystem> const& p)
868 {
869 static_assert(I < 4);
870 if constexpr (I == 0) {
871 return p.X();
872 } else if constexpr (I == 1) {
873 return p.Y();
874 } else if constexpr (I == 2) {
875 return p.Z();
876 } else {
877 return p.E();
878 }
879 }
880
881 } // end namespace Math
882
883} // end namespace ROOT
884
885// Structured bindings
886#include <tuple>
887namespace std {
888 template <class CoordSystem>
889 struct tuple_size<ROOT::Math::LorentzVector<CoordSystem>> : integral_constant<size_t, 4> {};
890 template <size_t I, class CoordSystem>
891 struct tuple_element<I, ROOT::Math::LorentzVector<CoordSystem>> {
892 static_assert(I < 4);
893 using type = typename CoordSystem::Scalar;
894 };
895}
896
897#include <sstream>
898namespace cling
899{
900template<typename CoordSystem>
901std::string printValue(const ROOT::Math::LorentzVector<CoordSystem> *v)
902{
903 std::stringstream s;
904 s << *v;
905 return s.str();
906}
907
908} // end namespace cling
909
910#endif
911
912//#include "Math/GenVector/LorentzVectorOperations.h"
#define d(i)
Definition RSha256.hxx:102
#define b(i)
Definition RSha256.hxx:100
#define c(i)
Definition RSha256.hxx:101
#define a(i)
Definition RSha256.hxx:99
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
winID h TVirtualViewer3D TVirtualGLPainter p
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t dest
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t src
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h Atom_t Int_t ULong_t ULong_t unsigned char prop_list Atom_t Atom_t Atom_t Time_t type
float * q
Class describing a generic displacement vector in 3 dimensions.
Class describing a generic LorentzVector in the 4D space-time, using the specified coordinate system ...
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
constexpr LorentzVector(const LorentzVector< Coords > &v)
constructor from a LorentzVector expressed in different coordinates, or using a different Scalar type
Scalar Py() const
spatial Y component
constexpr LorentzVector(const ForeignLorentzVector &v)
Construct from a foreign 4D vector type, for example, HepLorentzVector Precondition: v must implement...
LorentzVector(const Scalar &a, const Scalar &b, const Scalar &c, const Scalar &d)
generic constructors from four scalar values.
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 operator-() const
Negative of a LorentzVector (q = - v )
static constexpr unsigned int fDimension
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
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 )
Scalar DeltaR(const OtherLorentzVector &v, const bool useRapidity=false) const
deltaRapidity between this and vector v
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
unsigned int Dimension() const
dimension
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)
LorentzVector< CoordSystem >::Scalar AsymmetryScalar(LorentzVector< CoordSystem > const &pp, LorentzVector< CoordSystem > const &pm)
pair (p+ p-) scalar asymmetry As = ||Pt+| - |Pt-|/||Pt+| + |Pt-||.
LorentzVector< CoordSystem >::Scalar Acoplanarity(LorentzVector< CoordSystem > const &pp, LorentzVector< CoordSystem > const &pm)
pair (p+ p-) acoplanarity alpha = 1 - |phi+ - phi-|/pi.
LorentzVector< CoordSystem >::Scalar AsymmetryVectorial(LorentzVector< CoordSystem > const &pp, LorentzVector< CoordSystem > const &pm)
pair (p+ p-) vectorial asymmetry Av = |Pt+ - Pt-|/|Pt+ + Pt-|.
Double_t y[n]
Definition legend1.C:17
Double_t x[n]
Definition legend1.C:17
const Int_t n
Definition legend1.C:16
#define I(x, y, z)
Namespace for new Math classes and functions.
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)
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
void GenVector_Throw(const char *)
function throwing exception, by creating internally a GenVector_exception only when needed
std::basic_istream< char_t, traits_t > & operator>>(std::basic_istream< char_t, traits_t > &is, DisplacementVector2D< T, U > &v)
CoordSystem::Scalar get(DisplacementVector2D< CoordSystem, Tag > const &p)
AxisAngle operator*(RotationX const &r1, AxisAngle const &r2)
Multiplication of an axial rotation by an AxisAngle.
VecExpr< UnaryOp< Sqrt< T >, VecExpr< A, T, D >, T >, T, D > sqrt(const VecExpr< A, T, D > &rhs)
constexpr Double_t Pi()
Definition TMath.h:40
constexpr Double_t TwoPi()
Definition TMath.h:47
auto * tt
Definition textangle.C:16