MINT2
polVector.h
Go to the documentation of this file.
1 #ifndef POLARIZATION_VECTOR_HH
2 #define POLARIZATION_VECTOR_HH
3 // author: Philippe d'Argent (p.dargent@cern.ch)
4 // status: Tu 6 May 2014
5 
6 // Spin 1 polarization vectors taken from:
7 // S.U. Chung, Spin Formalisms (CERN Yellow Report No. CERN 71-8, Geneva, Switzerland, 1971)
8 
9 #include "TLorentzVectorC.h"
10 #include "TLorentzVector.h"
11 
12 using namespace std;
13 
14 class polVector : public TLorentzVectorC{
15 
16  public:
17 
18  polVector(const TLorentzVector p,const double m, const int lambda) :TLorentzVectorC(){
19 
20 
21  if(m>0.){
22  if(lambda==-1){
23  _v[0].SetXYZT(m+p.Px()*p.Px()/(p.E()+m),p.Px()*p.Py()/(p.E()+m),p.Px()*p.Pz()/(p.E()+m),p.Px());
24  _v[0]*=1./(m*sqrt(2));
25  _v[1].SetXYZT(p.Px()*p.Py()/(p.E()+m),m+p.Py()*p.Py()/(p.E()+m),p.Py()*p.Pz()/(p.E()+m),p.Py());
26  _v[1]*=-1./(m*sqrt(2));
27  }
28 
29  else if(lambda==0){
30  _v[0].SetXYZT(p.Pz()*p.Px()/(p.E()+m),p.Pz()*p.Py()/(p.E()+m),m+p.Pz()*p.Pz()/(p.E()+m),p.Pz());
31  _v[0]*=1./m;
32  _v[1].SetXYZT(0,0,0,0);
33  }
34  else if(lambda==1){
35  _v[0].SetXYZT(m+p.Px()*p.Px()/(p.E()+m),p.Px()*p.Py()/(p.E()+m),p.Px()*p.Pz()/(p.E()+m),p.Px());
36  _v[0]*=-1./(m*sqrt(2));
37  _v[1].SetXYZT(p.Px()*p.Py()/(p.E()+m),m+p.Py()*p.Py()/(p.E()+m),p.Py()*p.Pz()/(p.E()+m),p.Py());
38  _v[1]*=-1./(m*sqrt(2));
39  }
40  else if(lambda==-999){
41  TLorentzVectorC tmp = polVector(p,m,-1) + polVector(p,m,0) + polVector(p,m,1);
42  _v[0]= tmp.Re();
43  _v[1]= tmp.Im();
44  }
45 
46  else {
47  std::cout << "I can't handle spin > 1 particles. I'll set everything to 0. " << std::endl;
48  _v[0].SetXYZT(0,0,0,0);
49  _v[1].SetXYZT(0,0,0,0);
50  }
51 
52  }
53  else{
54 
55  double px = p.Px()/p.E();
56  double py = p.Py()/p.E();
57  double pz = p.Pz()/p.E();
58 
59  if((abs(px) > 1.E-4)||(abs(py) > 1.E-4)){
60  double x = -py/sqrt(px*px + py*py);
61  double y = px/sqrt(px*px + py*py);
62  double c = pz;
63  double s = sin(acos(c));
64 
65  if(lambda==-1){
66  _v[0].SetXYZT((x*x*(1-c)+c),x*y*(1-c),-y*s,0 );
67  _v[0]*=1./sqrt(2);
68  _v[1].SetXYZT(-x*y*(1-c),-(y*y*(1-c)+c),-x*s,0);
69  _v[1]*=1./sqrt(2);
70  }
71 
72  else if(lambda==1){
73  _v[0].SetXYZT(-(x*x*(1-c)+c),-x*y*(1-c),y*s,0 );
74  _v[0]*=1./sqrt(2);
75  _v[1].SetXYZT(-x*y*(1-c),-(y*y*(1-c)+c),-x*s,0);
76  _v[1]*=1./sqrt(2);
77  }
78 
79  else if(lambda==-999){
80  TLorentzVectorC tmp = polVector(p,m,-1) + polVector(p,m,1);
81  _v[0]= tmp.Re();
82  _v[1]= tmp.Im();
83  }
84 
85  else {
86  _v[0].SetXYZT(0,0,0,0);
87  _v[1].SetXYZT(0,0,0,0);
88  }
89 
90 
91  }
92 
93  else{
94 
95  if(lambda==-1){
96  _v[0].SetXYZT(1/sqrt(2),0,0,0);
97  _v[1].SetXYZT(0,-1/sqrt(2),0,0);
98  }
99 
100  else if(lambda==1){
101  _v[0].SetXYZT(-1/sqrt(2),0,0,0);
102  _v[1].SetXYZT(0,-1/sqrt(2),0,0);
103  }
104 
105  else if(lambda==-999){
106  TLorentzVectorC tmp = polVector(p,m,-1) + polVector(p,m,1);
107  _v[0]= tmp.Re();
108  _v[1]= tmp.Im();
109  }
110 
111  else {
112  _v[0].SetXYZT(0,0,0,0);
113  _v[1].SetXYZT(0,0,0,0);
114  }
115 
116 
117  }
118 
119 
120  }
121 
122 
123  }
124 
125  polVector(const polVector& other, bool conj = false) :TLorentzVectorC(other,conj){}
126 
127 
128 };
129 #endif
130 //
s
static const double s
Definition: CLHEPSystemOfUnits.h:126
TLorentzVectorC::Re
const TLorentzVector & Re() const
Definition: TLorentzVectorC.h:43
polVector::polVector
polVector(const polVector &other, bool conj=false)
Definition: polVector.h:125
m
static const double m
Definition: CLHEPSystemOfUnits.h:86
TLorentzVectorC::Im
const TLorentzVector & Im() const
Definition: TLorentzVectorC.h:44
lambda
double lambda(double x, double y, double z)
Definition: lambda.h:8
polVector::polVector
polVector(const TLorentzVector p, const double m, const int lambda)
Definition: polVector.h:18
Generated by   doxygen 1.8.15