]>
Commit | Line | Data |
---|---|---|
63ba5337 | 1 | // HelicityBasics.h is a part of the PYTHIA event generator. |
2 | // Copyright (C) 2012 Philip Ilten, Torbjorn Sjostrand. | |
3 | // PYTHIA is licenced under the GNU GPL version 2, see COPYING for details. | |
4 | // Please respect the MCnet Guidelines, see GUIDELINES for details. | |
5 | ||
6 | // Header file for a number of helper classes used in tau decays. | |
7 | ||
8 | #ifndef Pythia8_HelicityBasics_H | |
9 | #define Pythia8_HelicityBasics_H | |
10 | ||
11 | #include "Basics.h" | |
12 | #include "Event.h" | |
13 | #include "PythiaComplex.h" | |
14 | #include "PythiaStdlib.h" | |
15 | ||
16 | namespace Pythia8 { | |
17 | ||
18 | //========================================================================== | |
19 | ||
20 | // The Wave4 class provides a class for complex four-vector wave functions. | |
21 | // The Wave4 class can be multiplied with the GammaMatrix class to allow | |
22 | // for the writing of helicity matrix elements. | |
23 | ||
24 | class Wave4 { | |
25 | ||
26 | public: | |
27 | ||
28 | // Constructors and destructor. | |
29 | Wave4() {}; | |
30 | Wave4(complex v0, complex v1, complex v2, complex v3) {val[0] = v0; | |
31 | val[1] = v1; val[2] = v2; val[3] = v3;} | |
32 | Wave4(Vec4 v) {val[0] = v.e(); val[1] = v.px(); val[2] = v.py(); | |
33 | val[3] = v.pz();} | |
34 | ~Wave4() {}; | |
35 | ||
36 | // Access an element of the wave vector. | |
37 | complex& operator() (int i) {return val[i];} | |
38 | ||
39 | // Wave4 + Wave4. | |
40 | Wave4 operator+(Wave4 w) {return Wave4( val[0] + w.val[0], | |
41 | val[1] + w.val[1], val[2] + w.val[2], val[3] + w.val[3]);} | |
42 | ||
43 | // Wave4 - Wave4. | |
44 | Wave4 operator-(Wave4 w) {return Wave4( val[0] - w.val[0], | |
45 | val[1] - w.val[1], val[2] - w.val[2], val[3] - w.val[3]);} | |
46 | ||
47 | // - Wave4. | |
48 | Wave4 operator-() {return Wave4(-val[0], -val[1], -val[2], -val[3]);} | |
49 | ||
50 | // Wave4 * Wave4. | |
51 | complex operator*(Wave4 w) {return val[0] * w.val[0] | |
52 | + val[1] * w.val[1] + val[2] * w.val[2] + val[3] * w.val[3];} | |
53 | ||
54 | // Wave4 * complex. | |
55 | Wave4 operator*(complex s) {return Wave4(val[0] * s, val[1] * s, | |
56 | val[2] * s, val[3] * s);} | |
57 | ||
58 | // complex * Wave4. | |
59 | friend Wave4 operator*(complex s, const Wave4& w); | |
60 | ||
61 | // Wave4 * double. | |
62 | Wave4 operator*(double s) {return Wave4(val[0] * s, val[1] * s, | |
63 | val[2] * s, val[3] * s);} | |
64 | ||
65 | // double * Wave4. | |
66 | friend Wave4 operator*(double s, const Wave4& w); | |
67 | ||
68 | // Wave4 / complex. | |
69 | Wave4 operator/(complex s) {return Wave4(val[0] / s, val[1] / s, | |
70 | val[2] / s, val[3] / s);} | |
71 | ||
72 | // Wave4 / double. | |
73 | Wave4 operator/(double s) {return Wave4(val[0] / s, val[1] / s, | |
74 | val[2]/s, val[3]/s);} | |
75 | ||
76 | // Complex conjugate. | |
77 | friend Wave4 conj(Wave4 w); | |
78 | ||
79 | // Permutation operator. | |
80 | friend Wave4 epsilon(Wave4 w1, Wave4 w2, Wave4 w3); | |
81 | ||
82 | // Invariant squared mass for REAL Wave4 (to save time). | |
83 | friend double m2(Wave4 w); | |
84 | friend double m2(Wave4 w1, Wave4 w2); | |
85 | ||
86 | // Wave4 * GammaMatrix multiplication is defined in the GammaMatrix class. | |
87 | ||
88 | // Print a Wave4 vector. | |
89 | friend ostream& operator<<(ostream& output, Wave4 w); | |
90 | ||
91 | protected: | |
92 | ||
93 | complex val[4]; | |
94 | ||
95 | }; | |
96 | ||
97 | //-------------------------------------------------------------------------- | |
98 | ||
99 | // Namespace function declarations; friends of Wave4 class. | |
100 | Wave4 operator*(complex s, const Wave4& w); | |
101 | Wave4 conj(Wave4 w); | |
102 | Wave4 epsilon(Wave4 w1, Wave4 w2, Wave4 w3); | |
103 | double m2(Wave4 w); | |
104 | double m2(Wave4 w1, Wave4 w2); | |
105 | ostream& operator<< (ostream& os, Wave4 w); | |
106 | ||
107 | //========================================================================== | |
108 | ||
109 | // The GammaMatrix class is a special sparse matrix class used to write | |
110 | // helicity matrix elements in conjuction with the Wave4 class. Note that | |
111 | // only left to right multplication of Wave4 vectors with the GammaMatrix | |
112 | // class is allowed. Additionally, subtracting a scalar from a GammaMatrix | |
113 | // (or subtracting a GammaMatrix from a scalar) subtracts the scalar from | |
114 | //each non-zero element of the GammaMatrix. This is designed specifically | |
115 | // with the (1 - gamma^5) structure of matrix elements in mind. | |
116 | ||
117 | class GammaMatrix { | |
118 | ||
119 | public: | |
120 | ||
121 | // Constructors and destructor. | |
122 | GammaMatrix() {}; | |
123 | GammaMatrix(int mu); | |
124 | ~GammaMatrix() {}; | |
125 | ||
126 | // Access an element of the matrix. | |
127 | complex& operator() (int I, int J) {if (index[J] == I) return val[J]; | |
128 | else return COMPLEXZERO; } | |
129 | ||
130 | // Wave4 * GammaMatrix. | |
131 | friend Wave4 operator*(Wave4 w, GammaMatrix g); | |
132 | ||
133 | // GammaMatrix * Scalar. | |
134 | GammaMatrix operator*(complex s) {val[0] = s*val[0]; val[1] = s*val[1]; | |
135 | val[2] = s*val[2]; val[3] = s*val[3]; return *this;} | |
136 | ||
137 | // Scalar * GammaMatrix. | |
138 | friend GammaMatrix operator*(complex s, GammaMatrix g); | |
139 | ||
140 | // Gamma5 - I * Scalar. | |
141 | GammaMatrix operator-(complex s) {val[0] = val[0] - s; val[1] = val[1] - s; | |
142 | val[2] = val[2] - s; val[3] = val[3] - s; return *this;} | |
143 | ||
144 | // I * Scalar - Gamma5. | |
145 | friend GammaMatrix operator-(complex s, GammaMatrix g); | |
146 | ||
147 | // Gamma5 + I * Scalar | |
148 | GammaMatrix operator+(complex s) {val[0] = val[0] + s; val[1] = val[1] + s; | |
149 | val[2] = val[2] + s; val[3] = val[3] + s; return *this;} | |
150 | ||
151 | // I * Scalar + Gamma5 | |
152 | friend GammaMatrix operator+(complex s, GammaMatrix g); | |
153 | ||
154 | // << GammaMatrix. | |
155 | friend ostream& operator<< (ostream& os, GammaMatrix g); | |
156 | ||
157 | protected: | |
158 | ||
159 | complex val[4]; | |
160 | int index[4]; | |
161 | ||
162 | // Need to define complex 0 as a variable for operator() to work. | |
163 | complex COMPLEXZERO; | |
164 | ||
165 | }; | |
166 | ||
167 | //-------------------------------------------------------------------------- | |
168 | ||
169 | // Namespace function declarations; friends of GammaMatrix class. | |
170 | Wave4 operator*(Wave4 w, GammaMatrix g); | |
171 | GammaMatrix operator*(complex s, GammaMatrix g); | |
172 | GammaMatrix operator-(complex s, GammaMatrix g); | |
173 | GammaMatrix operator+(complex s, GammaMatrix g); | |
174 | ostream& operator<< (ostream& os, GammaMatrix g); | |
175 | ||
176 | //========================================================================== | |
177 | ||
178 | // Helicity particle class containing helicity information, derived from | |
179 | // particle base class. | |
180 | ||
181 | class HelicityParticle : public Particle { | |
182 | ||
183 | public: | |
184 | ||
185 | // Constructors. | |
186 | HelicityParticle() : Particle() { direction = 1;} | |
187 | HelicityParticle(int idIn, int statusIn = 0, int mother1In = 0, | |
188 | int mother2In = 0, int daughter1In = 0, int daughter2In = 0, | |
189 | int colIn = 0, int acolIn = 0, double pxIn = 0., | |
190 | double pyIn = 0., double pzIn = 0., double eIn = 0., | |
191 | double mIn = 0., double scaleIn = 0., ParticleData* ptr = 0) | |
192 | : Particle(idIn, statusIn, mother1In, mother2In, daughter1In, daughter2In, | |
193 | colIn, acolIn, pxIn, pyIn, pzIn, eIn, mIn, scaleIn) { | |
194 | if (ptr) { setPDTPtr(ptr); setPDEPtr(); } | |
195 | rho = vector< vector<complex> >(spinStates(), | |
196 | vector<complex>(spinStates(), 0)); | |
197 | D = vector< vector<complex> >(spinStates(), | |
198 | vector<complex>(spinStates(), 0)); | |
199 | for (int i = 0; i < spinStates(); i++) { rho[i][i] = 0.5; D[i][i] = 1.;} | |
200 | direction = 1; } | |
201 | HelicityParticle(int idIn, int statusIn, int mother1In, int mother2In, | |
202 | int daughter1In, int daughter2In, int colIn, int acolIn, Vec4 pIn, | |
203 | double mIn = 0., double scaleIn = 0., ParticleData* ptr = 0) | |
204 | : Particle(idIn, statusIn, mother1In, mother2In, daughter1In, daughter2In, | |
205 | colIn, acolIn, pIn, mIn, scaleIn) { | |
206 | if (ptr) { setPDTPtr(ptr); setPDEPtr();} | |
207 | rho = vector< vector<complex> >(spinStates(), | |
208 | vector<complex>(spinStates(), 0)); | |
209 | D = vector< vector<complex> >(spinStates(), | |
210 | vector<complex>(spinStates(), 0)); | |
211 | for (int i = 0; i < spinStates(); i++) { rho[i][i] = 0.5; D[i][i] = 1;} | |
212 | direction = 1; } | |
213 | HelicityParticle(const Particle& ptIn, ParticleData* ptr = 0) | |
214 | : Particle(ptIn) { | |
215 | if (ptr) { setPDTPtr(ptr); setPDEPtr();} | |
216 | rho = vector< vector<complex> >(spinStates(), | |
217 | vector<complex>(spinStates(), 0)); | |
218 | D = vector< vector<complex> >(spinStates(), | |
219 | vector<complex>(spinStates(), 0)); | |
220 | for (int i = 0; i < spinStates(); i++) { rho[i][i] = 0.5; D[i][i] = 1;} | |
221 | direction = 1; } | |
222 | ||
223 | // Methods. | |
224 | Wave4 wave(int h); | |
225 | Wave4 waveBar(int h); | |
226 | void normalize(vector< vector<complex> >& m); | |
227 | int spinStates(); | |
228 | ||
229 | // Event record position. | |
230 | int idx; | |
231 | ||
232 | // Flag for whether particle is incoming (-1) or outgoing (1). | |
233 | int direction; | |
234 | ||
235 | // Helicity density matrix. | |
236 | vector< vector<complex> > rho; | |
237 | ||
238 | // Decay matrix. | |
239 | vector< vector<complex> > D; | |
240 | ||
241 | private: | |
242 | ||
243 | // Constants: could only be changed in the code itself. | |
244 | static const double TOLERANCE; | |
245 | ||
246 | }; | |
247 | ||
248 | //========================================================================== | |
249 | ||
250 | } // end namespace Pythia8 | |
251 | ||
252 | #endif // end Pythia8_HelicityBasics_H |