StRoot  1
HelicityMatrixElements.h
1 // HelicityMatrixElements.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 physics classes used in tau decays.
7
8 #ifndef Pythia8_HelicityMatrixElements_H
9 #define Pythia8_HelicityMatrixElements_H
10
11 #include "Basics.h"
12 #include "Event.h"
13 #include "HelicityBasics.h"
14 #include "PythiaComplex.h"
15 #include "PythiaStdlib.h"
16 #include "StandardModel.h"
17
18 namespace Pythia8 {
19
20 //==========================================================================
21
22 // The helicity matrix element class.
23
25
26 public:
27
28  // Constructor and destructor.
30  virtual ~HelicityMatrixElement() {};
31
32  // Initialize the physics matrices and pointers.
33  virtual void initPointers(ParticleData*, Couplings*);
34
35  // Initialize the channel.
36  virtual HelicityMatrixElement* initChannel(vector<HelicityParticle>&);
37
38  // Calculate the matrix element weight for a decay.
39  virtual double decayWeight(vector<HelicityParticle>&);
40
41  // Calculate the maximum matrix element decay weight.
42  virtual double decayWeightMax(vector<HelicityParticle>&)
43  {return DECAYWEIGHTMAX;}
44
45  // Calculate the helicity matrix element.
46  virtual complex calculateME(vector<int>){return complex(0,0);}
47
48  // Calculate the decay matrix for a particle.
49  virtual void calculateD(vector<HelicityParticle>&);
50
51  // Calculate the density matrix for a particle.
52  virtual void calculateRho(unsigned int, vector<HelicityParticle>&);
53
54  // Set a fermion line.
55  void setFermionLine(int, HelicityParticle&, HelicityParticle&);
56
57  // Calculate Breit-Wigner's with running widths.
58  virtual complex sBreitWigner(double m0, double m1, double s,
59  double M, double G);
60  virtual complex pBreitWigner(double m0, double m1, double s,
61  double M, double G);
62  virtual complex dBreitWigner(double m0, double m1, double s,
63  double M, double G);
64
65 protected:
66
67  // Maximum decay weight. WARNING: hardcoded constant.
68  double DECAYWEIGHTMAX;
69
70  // Physics matrices.
71  vector< GammaMatrix > gamma;
72
73  // Particle map vector.
74  vector< int > pMap;
75
76  // Particle ID vector.
77  vector< int > pID;
78
79  // Particle mass vector.
80  vector< double > pM;
81
82  // Wave functions.
83  vector< vector< Wave4 > > u;
84
85  // Initialize the constants for the matrix element (called by initChannel).
86  virtual void initConstants() {};
87
88  // Initialize the wave functions (called by decayWeight and calculateRho/D).
89  virtual void initWaves(vector<HelicityParticle>&) {};
90
91  // Pointer to particle data.
92  ParticleData* particleDataPtr;
93
94  // Pointer to Standard Model constants.
95  Couplings* couplingsPtr;
96
97 private:
98
99  // Recursive sub-method to calculate the density matrix for a particle.
100  void calculateRho(unsigned int, vector<HelicityParticle>&,
101  vector<int>&, vector<int>&, unsigned int);
102
103  // Recursive sub-method to calculate the decay matrix for a particle.
104  void calculateD(vector<HelicityParticle>&, vector<int>&, vector<int>&,
105  unsigned int);
106
107  // Recursive sub-method to calculate the matrix element weight for a decay.
108  void decayWeight(vector<HelicityParticle>&, vector<int>&, vector<int>&,
109  complex&, unsigned int);
110
111  // Calculate the product of the decay matrices for a hard process.
112  complex calculateProductD(unsigned int, unsigned int,
113  vector<HelicityParticle>&, vector<int>&, vector<int>&);
114
115  // Calculate the product of the decay matrices for a decay process.
116  complex calculateProductD(vector<HelicityParticle>&,
117  vector<int>&, vector<int>&);
118
119 };
120
121 //==========================================================================
122
123 // Helicity matrix element for the hard process of two fermions -> W ->
124 // two fermions.
125
127
128 public:
129
130  void initWaves(vector<HelicityParticle>&);
131
132  complex calculateME(vector<int>);
133
134 };
135
136 //==========================================================================
137
138 // Helicity matrix element for the hard process of two fermions ->
139 // photon -> two fermions.
140
142
143 public:
144
145  void initWaves(vector<HelicityParticle>&);
146
147  complex calculateME(vector<int>);
148
149 private:
150
151  // Fermion line charge and interaction energy.
152  double p0Q, p2Q, s;
153
154 };
155
156 //==========================================================================
157
158 // Helicity matrix element for the hard process of two fermions ->
159 // Z -> two fermions.
160
162
163 public:
164
165  void initConstants();
166
167  void initWaves(vector<HelicityParticle>&);
168
169  complex calculateME(vector<int>);
170
171 private:
172
173  // Vector and axial couplings.
174  double p0CA, p2CA, p0CV, p2CV;
175
176  // Weinberg angle, Z width (on-shell), Z mass (on-shell), and s.
177  double cos2W, sin2W, zG, zM, s;
178
179  // Bool whether the incoming fermions are oriented with the z-axis.
180  bool zaxis;
181
182 };
183
184 //==========================================================================
185
186 // Helicity matrix element for the hard process of two fermions ->
187 // photon/Z -> two fermions with full interference.
188
189 // In general the initPointers and initChannel methods should not need to be
190 // redeclared, but in this case each matrix element needs to be initialized.
191
193
194 public:
195
196  HelicityMatrixElement* initChannel(vector<HelicityParticle>&);
197
198  void initPointers(ParticleData*, Couplings*);
199
200  void initWaves(vector<HelicityParticle>&);
201
202  complex calculateME(vector<int>);
203
204 private:
205
208
209 };
210
211 //==========================================================================
212
213 // Helicity matrix element for the decay of a CP even Higgs -> two fermions.
214
215 // Because the Higgs is spin zero the Higgs production mechanism is not
216 // needed for calculating helicity density matrices.
217
219
220 public:
221
222  void initWaves(vector<HelicityParticle>&);
223
224  complex calculateME(vector<int>);
225
226 private:
227
228  // Coupling constants of the fermions with the Higgs.
229  double p2CA, p2CV;
230
231 };
232
233 //==========================================================================
234
235 // Helicity matrix element for the decay of a CP odd Higgs -> two fermions.
236
238
239 public:
240
241  void initWaves(vector<HelicityParticle>&);
242
243  complex calculateME(vector<int>);
244
245 private:
246
247  // Coupling constants of the fermions with the Higgs.
248  double p2CA, p2CV;
249
250 };
251
252 //==========================================================================
253
254 // Helicity matrix element for the decay of a charged Higgs -> two fermions.
255
257
258 public:
259
260  void initWaves(vector<HelicityParticle>&);
261
262  complex calculateME(vector<int>);
263
264 private:
265
266  // Coupling constants of the fermions with the Higgs.
267  double p2CA, p2CV;
268
269 };
270
271 //==========================================================================
272
273 // Helicity matrix element which provides an unpolarized on-diagonal helicity
274 // density matrix. Used for unknown hard processes.
275
277
278 public:
279
280  void calculateRho(unsigned int, vector<HelicityParticle>&);
281
282 };
283
284 //==========================================================================
285
286 // Base class for all tau decay helicity matrix elements.
287
289
290 public:
291
292  virtual void initWaves(vector<HelicityParticle>&);
293
294  virtual complex calculateME(vector<int>);
295
296  virtual double decayWeightMax(vector<HelicityParticle>&);
297
298 protected:
299
301
302  virtual void calculateResonanceWeights(vector<double>&, vector<double>&,
303  vector<complex>&);
304
305 };
306
307 //==========================================================================
308
309 // Helicity matrix element for a tau decaying into a single scalar meson.
310
311 class HMETau2Meson : public HMETauDecay {
312
313 public:
314
315  void initConstants();
316
318
319 };
320
321 //==========================================================================
322
323 // Helicity matrix element for a tau decaying into two leptons.
324
326
327 public:
328
329  void initConstants();
330
331  void initWaves(vector<HelicityParticle>&);
332
333  complex calculateME(vector<int>);
334
335 };
336
337 //==========================================================================
338
339 // Helicity matrix element for a tau decaying into two mesons through a
340 // vector meson resonance.
341
343
344 public:
345
346  void initConstants();
347
349
350 private:
351
352  // Resonance masses, widths, and weights.
353  vector<double> vecM, vecG, vecP, vecA;
354  vector<complex> vecW;
355
356 };
357
358 //==========================================================================
359
360 // Helicity matrix element for a tau decay into two mesons through a vector
361 // or scalar meson resonance.
362
364
365 public:
366
367  void initConstants();
368
370
371 private:
372
373  // Coupling to vector and scalar resonances.
374  double scaC, vecC;
375
376  // Resonance masses, widths, and weights.
377  vector<double> scaM, scaG, scaP, scaA, vecM, vecG, vecP, vecA;
378  vector<complex> scaW, vecW;
379
380 };
381
382 //==========================================================================
383
384 // Helicity matrix element for a tau decay into three pions.
385
387
388 public:
389
390  void initConstants();
391
393
394 private:
395
396  // Resonance masses, widths, and weights.
397  vector<double> rhoM, rhoG, rhoPp, rhoAp, rhoPd, rhoAd;
398  double f0M, f0G, f0P, f0A, f2M, f2G, f2P, f2A;
399  double sigM, sigG, sigP, sigA;
400  vector<complex> rhoWp, rhoWd;
401  complex f0W, f2W, sigW;
402
403  // Center of mass energies.
404  double s1, s2, s3, s4;
405
406  // Form factors.
407  complex F1();
408  complex F2();
409  complex F3();
410
411  // Running width and Breit-Wigner for the a1.
412  double a1Width(double);
413  complex a1BreitWigner(double);
414
415 };
416
417 //==========================================================================
418
419 // Helicity matrix element for a tau decay into four pions.
420
422
423 public:
424
425  void initConstants();
426
428
429 private:
430
431  // G-function form factors (fits).
432  double G(int i, double s);
433
434  // T-vector functions.
435  Wave4 t1(int, int, int, int);
436  Wave4 t2(int, int, int, int);
437  Wave4 t3(int, int, int, int);
438
439  // Breit-Wigner denominators for the intermediate mesons.
440  complex a1D(double s);
441  complex rhoD(double s);
442  complex sigD(double s);
443  complex omeD(double s);
444
445  // Form factors needed for the a1, rho, and omega.
446  double a1FormFactor(double s);
447  double rhoFormFactor1(double s);
448  double rhoFormFactor2(double s);
449  double omeFormFactor(double s);
450
451  // Masses and widths of the intermediate mesons.
452  double a1M, a1G, rhoM, rhoG, sigM, sigG, omeM, omeG;
453
454  // Masses for the pions (charged and neutral).
455  double picM, pinM;
456
457  // Amplitude, phases, and weights for mixing.
458  double sigA, sigP, omeA, omeP;
459  complex sigW, omeW;
460
461  // Cut-off for a1 form factor.
462  double lambda2;
463
464 };
465
466 //==========================================================================
467
468 // Helicity matrix element for a tau decay into flat phase space.
469
471
472 public:
473
474  void initWaves(vector<HelicityParticle>&) {};
475
476  complex calculateME(vector<int>) {return 1;}
477
478  void calculateD(vector<HelicityParticle>&) {};
479
480  void calculateRho(unsigned int, vector<HelicityParticle>&) {};
481
482  double decayWeight(vector<HelicityParticle>&) {return 1.0;}
483
484  double decayWeightMax(vector<HelicityParticle>&) {return 1.0;}
485
486 };
487
488 //==========================================================================
489
490 } // end namespace Pythia8
491
492 #endif // end Pythia8_HelicityMatrixElements_H