PYTHIA  8.312
HelicityMatrixElements.h
1 // HelicityMatrixElements.h is a part of the PYTHIA event generator.
2 // Copyright (C) 2024 Philip Ilten, Torbjorn Sjostrand.
3 // PYTHIA is licenced under the GNU GPL v2 or later, 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 "Pythia8/Basics.h"
12 #include "Pythia8/Event.h"
13 #include "Pythia8/HelicityBasics.h"
14 #include "Pythia8/PythiaComplex.h"
15 #include "Pythia8/PythiaStdlib.h"
16 #include "Pythia8/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  coupSMPtr(), settingsPtr() {};
31  virtual ~HelicityMatrixElement() {};
32 
33  // Initialize the physics matrices and pointers.
34  virtual void initPointers(ParticleData*, CoupSM* , Settings* = 0);
35 
36  // Initialize the channel.
37  virtual HelicityMatrixElement* initChannel(vector<HelicityParticle>&);
38 
39  // Calculate the matrix element weight for a decay.
40  virtual double decayWeight(vector<HelicityParticle>&);
41 
42  // Calculate the maximum matrix element decay weight.
43  virtual double decayWeightMax(vector<HelicityParticle>&)
44  {return DECAYWEIGHTMAX;}
45 
46  // Calculate the helicity matrix element.
47  virtual complex calculateME(vector<int>){return complex(0,0);}
48 
49  // Calculate the decay matrix for a particle.
50  virtual void calculateD(vector<HelicityParticle>&);
51 
52  // Calculate the density matrix for a particle.
53  virtual void calculateRho(unsigned int, vector<HelicityParticle>&);
54 
55  // Set a fermion line.
57 
58  // Calculate Breit-Wigner's with running widths and fixed.
59  virtual complex breitWigner(double s, double M, double G);
60  virtual complex sBreitWigner(double m0, double m1, double s,
61  double M, double G);
62  virtual complex pBreitWigner(double m0, double m1, double s,
63  double M, double G);
64  virtual complex dBreitWigner(double m0, double m1, double s,
65  double M, double G);
66 
67 protected:
68 
69  // Maximum decay weight. WARNING: hardcoded constant.
71 
72  // Physics matrices.
73  vector< GammaMatrix > gamma;
74 
75  // Particle map vector.
76  vector< int > pMap;
77 
78  // Particle ID vector.
79  vector< int > pID;
80 
81  // Particle mass vector.
82  vector< double > pM;
83 
84  // Wave functions.
85  vector< vector< Wave4 > > u;
86 
87  // Initialize the constants for the matrix element (called by initChannel).
88  virtual void initConstants() {};
89 
90  // Initialize the wave functions (called by decayWeight and calculateRho/D).
91  virtual void initWaves(vector<HelicityParticle>&) {};
92 
93  // Pointer to particle data.
95 
96  // Pointer to Standard Model constants.
98 
99  // Pointer to Settings.
101 
102 private:
103 
104  // Recursive sub-method to calculate the density matrix for a particle.
105  void calculateRho(unsigned int, vector<HelicityParticle>&,
106  vector<int>&, vector<int>&, unsigned int);
107 
108  // Recursive sub-method to calculate the decay matrix for a particle.
109  void calculateD(vector<HelicityParticle>&, vector<int>&, vector<int>&,
110  unsigned int);
111 
112  // Recursive sub-method to calculate the matrix element weight for a decay.
113  void decayWeight(vector<HelicityParticle>&, vector<int>&, vector<int>&,
114  complex&, unsigned int);
115 
116  // Calculate the product of the decay matrices for a hard process.
117  complex calculateProductD(unsigned int, unsigned int,
118  vector<HelicityParticle>&, vector<int>&, vector<int>&);
119 
120  // Calculate the product of the decay matrices for a decay process.
121  complex calculateProductD(vector<HelicityParticle>&,
122  vector<int>&, vector<int>&);
123 
124 };
125 
126 //==========================================================================
127 
128 // Helicity matrix element for the hard process of two fermions -> W/W' ->
129 // two fermions.
130 
132 
133 public:
134 
135  HMETwoFermions2W2TwoFermions() : p0CA(), p2CA(), p0CV(), p2CV() {}
136 
137  void initConstants();
138 
139  void initWaves(vector<HelicityParticle>&);
140 
141  complex calculateME(vector<int>);
142 
143 private:
144 
145  // Vector and axial couplings.
146  double p0CA, p2CA, p0CV, p2CV;
147 
148 };
149 
150 //==========================================================================
151 
152 // Helicity matrix element for the hard process of two fermions ->
153 // photon/Z/Z' -> two fermions.
154 
156 
157 public:
158 
159  HMETwoFermions2GammaZ2TwoFermions() : p0CAZ(), p2CAZ(), p0CVZ(), p2CVZ(),
160  p0CAZp(), p2CAZp(), p0CVZp(), p2CVZp(), cos2W(), sin2W(), zG(), zM(),
161  zpG(), zpM(), s(), p0Q(), p2Q(), sMin(), zaxis(), includeGamma(),
162  includeZ(), includeZp() {}
163 
164  void initConstants();
165 
166  void initWaves(vector<HelicityParticle>&);
167 
168  complex calculateME(vector<int>);
169 
170 private:
171 
172  // Return gamma element for the helicity matrix element.
173  complex calculateGammaME(vector<int>);
174 
175  // Return Z/Z' element for helicity matrix element.
176  complex calculateZME(vector<int>, double, double, double, double,
177  double, double);
178 
179  // Return Z/Z' element, assuming massless fermions.
180  complex calculateZMEMasslessFermions(vector<int>, double, double, double,
181  double, double, double);
182 
183  // Return the Z' vector or axial coupling for a fermion.
184  double zpCoupling(int id, string type);
185 
186  // Vector and axial couplings.
187  double p0CAZ, p2CAZ, p0CVZ, p2CVZ, p0CAZp, p2CAZp, p0CVZp, p2CVZp;
188 
189  // Weinberg angle, Z/Z' width (on-shell), Z/Z' mass (on-shell), and s.
190  double cos2W, sin2W, zG, zM, zpG, zpM, s;
191 
192  // Fermion line charge.
193  double p0Q, p2Q;
194 
195  // Minimum s to use massless approximation.
196  double sMin;
197 
198  // Bool whether the incoming fermions are oriented with the z-axis.
199  bool zaxis, includeGamma, includeZ, includeZp;
200 
201 };
202 
203 //==========================================================================
204 
205 // Helicity matrix element for the hard process of two photons ->
206 // two fermions.
207 
209 
210 public:
211 
212  void initWaves(vector<HelicityParticle>&);
213 
214  complex calculateME(vector<int>);
215 
216 private:
217 
218  double tm, um, m;
219  Vec4 q0, q1;
220 
221 };
222 
223 //==========================================================================
224 
225 // Helicity matrix element for the hard process of X -> two fermions.
226 
228 
229 public:
230 
231  void initWaves(vector<HelicityParticle>&);
232 
233 };
234 
235 //==========================================================================
236 
237 // Helicity matrix element for the hard process of W/W' -> two fermions.
238 
240 
241 public:
242 
243  HMEW2TwoFermions() : p2CA(), p2CV() {}
244 
245  void initConstants();
246 
247  complex calculateME(vector<int>);
248 
249 private:
250 
251  // Vector and axial couplings.
252  double p2CA, p2CV;
253 
254 };
255 
256 //==========================================================================
257 
258 // Helicity matrix element for the hard process of photon -> two fermions.
259 
261 
262 public:
263 
264  complex calculateME(vector<int>);
265 
266 };
267 
268 //==========================================================================
269 
270 // Helicity matrix element for the hard process of Z/Z' -> two fermions.
271 
273 
274 public:
275 
276  HMEZ2TwoFermions() : p2CA(), p2CV() {}
277 
278  void initConstants();
279 
280  complex calculateME(vector<int>);
281 
282 private:
283 
284  // Return the Z' vector or axial coupling for a fermion.
285  double zpCoupling(int id, string type);
286 
287  // Vector and axial couplings.
288  double p2CA, p2CV;
289 
290 };
291 
292 //==========================================================================
293 
294 // Helicity matrix element for the decay of a Higgs -> two fermions.
295 
296 // Because the Higgs is spin zero the Higgs production mechanism is not
297 // needed for calculating helicity density matrices. However, the CP mixing
298 // is needed.
299 
301 
302 public:
303 
304  void initConstants();
305 
306  void initWaves(vector<HelicityParticle>&);
307 
308  complex calculateME(vector<int>);
309 
310 private:
311 
312  // Coupling constants of the fermions with the Higgs.
313  complex p2CA, p2CV;
314 
315 };
316 
317 //==========================================================================
318 
319 // Base class for all tau decay helicity matrix elements.
320 
322 
323 public:
324 
325  virtual void initWaves(vector<HelicityParticle>&);
326 
327  virtual complex calculateME(vector<int>);
328 
329  virtual double decayWeightMax(vector<HelicityParticle>&);
330 
331 protected:
332 
333  virtual void initHadronicCurrent(vector<HelicityParticle>&) {};
334 
335  virtual void calculateResonanceWeights(vector<double>&, vector<double>&,
336  vector<complex>&);
337 
338 };
339 
340 //==========================================================================
341 
342 // Helicity matrix element for a tau decaying into a single scalar meson.
343 
344 class HMETau2Meson : public HMETauDecay {
345 
346 public:
347 
348  void initConstants();
349 
350  void initHadronicCurrent(vector<HelicityParticle>&);
351 
352 };
353 
354 //==========================================================================
355 
356 // Helicity matrix element for a tau decaying into two leptons.
357 
359 
360 public:
361 
362  void initConstants();
363 
364  void initWaves(vector<HelicityParticle>&);
365 
366  complex calculateME(vector<int>);
367 
368 };
369 
370 //==========================================================================
371 
372 // Helicity matrix element for a tau decaying into two mesons through a
373 // vector meson resonance.
374 
376 
377 public:
378 
379  void initConstants();
380 
381  void initHadronicCurrent(vector<HelicityParticle>&);
382 
383 private:
384 
385  // Resonance masses, widths, and weights.
386  vector<double> vecM, vecG, vecP, vecA;
387  vector<complex> vecW;
388 
389 };
390 
391 //==========================================================================
392 
393 // Helicity matrix element for a tau decay into two mesons through a vector
394 // or scalar meson resonance.
395 
397 
398 public:
399 
400  HMETau2TwoMesonsViaVectorScalar() : scaC(), vecC() {}
401 
402  void initConstants();
403 
404  void initHadronicCurrent(vector<HelicityParticle>&);
405 
406 private:
407 
408  // Coupling to vector and scalar resonances.
409  double scaC, vecC;
410 
411  // Resonance masses, widths, and weights.
412  vector<double> scaM, scaG, scaP, scaA, vecM, vecG, vecP, vecA;
413  vector<complex> scaW, vecW;
414 
415 };
416 
417 //==========================================================================
418 
419 // Helicity matrix element for a tau decay into three mesons (base class).
420 
422 
423 public:
424 
425  HMETau2ThreeMesons() = default;
426 
427  void initConstants();
428 
429  void initHadronicCurrent(vector<HelicityParticle>&);
430 
431 protected:
432 
433  // Decay mode of the tau.
434  enum Mode{Pi0Pi0Pim, PimPimPip, Pi0PimK0b, PimPipKm, Pi0PimEta, PimKmKp,
435  Pi0K0Km, KlPimKs, Pi0Pi0Km, KlKlPim, PimKsKs, PimK0bK0, Uknown};
436  Mode mode{};
437 
438  // Initialize decay mode and resonance constants (called by initConstants).
439  virtual void initMode();
440  virtual void initResonances() {;}
441 
442  // Initialize the momenta.
443  virtual void initMomenta(vector<HelicityParticle>&);
444 
445  // Center of mass energies and momenta.
446  double s1{}, s2{}, s3{}, s4{};
447  Wave4 q{}, q2{}, q3{}, q4{};
448 
449  // Stored a1 Breit-Wigner (for speed).
450  complex a1BW{};
451 
452  // Form factors.
453  virtual complex F1() {return complex(0, 0);}
454  virtual complex F2() {return complex(0, 0);}
455  virtual complex F3() {return complex(0, 0);}
456  virtual complex F4() {return complex(0, 0);}
457 
458  // Phase space and Breit-Wigner for the a1.
459  virtual double a1PhaseSpace(double);
460  virtual complex a1BreitWigner(double);
461 
462  // Sum running p and fixed width Breit-Wigner resonances.
463  complex T(double m0, double m1, double s,
464  vector<double>& M, vector<double>& G, vector<double>& W);
465  complex T(double s, vector<double>& M, vector<double>& G, vector<double>& W);
466 
467 };
468 
469 //==========================================================================
470 
471 // Helicity matrix element for a tau decay into three pions.
472 
474 
475 public:
476 
477  HMETau2ThreePions() : f0M(), f0G(), f0P(), f0A(), f2M(), f2G(), f2P(),
478  f2A(), sigM(), sigG(), sigP(), sigA() {}
479 
480 private:
481 
482  void initResonances();
483 
484  // Resonance masses, widths, and weights.
485  vector<double> rhoM, rhoG, rhoPp, rhoAp, rhoPd, rhoAd;
486  double f0M, f0G, f0P, f0A, f2M, f2G, f2P, f2A;
487  double sigM, sigG, sigP, sigA;
488  vector<complex> rhoWp, rhoWd;
489  complex f0W, f2W, sigW;
490 
491  // Form factors.
492  complex F1();
493  complex F2();
494  complex F3();
495 
496  // Running width and Breit-Wigner for the a1.
497  double a1PhaseSpace(double);
498  complex a1BreitWigner(double);
499 
500 };
501 
502 //==========================================================================
503 
504 // Helicity matrix element for a tau decay into three mesons with kaons.
505 
507 
508 public:
509 
510  HMETau2ThreeMesonsWithKaons() : kM(), piM(), piW() {}
511 
512 private:
513 
514  void initResonances();
515 
516  // Resonance masses, widths, and weights.
517  vector<double> rhoMa, rhoGa, rhoWa, rhoMv, rhoGv, rhoWv;
518  vector<double> kstarMa, kstarGa, kstarWa, kstarMv, kstarGv, kstarWv;
519  vector<double> k1Ma, k1Ga, k1Wa, k1Mb, k1Gb, k1Wb;
520  vector<double> omegaM, omegaG, omegaW;
521  double kM, piM, piW;
522 
523  // Form factors.
524  complex F1();
525  complex F2();
526  complex F4();
527 
528 };
529 
530 //==========================================================================
531 
532 // Helicity matrix element for a tau decay into generic three mesons.
533 
535 
536 public:
537 
538  HMETau2ThreeMesonsGeneric() : kM(), piM(), piW() {}
539 
540 private:
541 
542  void initResonances();
543 
544  // Resonance masses, widths, and weights.
545  vector<double> rhoMa, rhoGa, rhoWa, rhoMv, rhoGv, rhoWv;
546  vector<double> kstarM, kstarG, kstarW, k1M, k1G, k1W;
547  double kM, piM, piW;
548 
549  // Form factors.
550  complex F1();
551  complex F2();
552  complex F4();
553 
554 };
555 
556 //==========================================================================
557 
558 // Helicity matrix element for a tau decay into two pions and a photon.
559 
561 
562 public:
563 
564  HMETau2TwoPionsGamma() : piM() {}
565 
566  void initConstants();
567 
568  void initWaves(vector<HelicityParticle>&);
569 
570  complex calculateME(vector<int>);
571 
572 protected:
573 
574  // Resonance masses, widths, and weights.
575  vector<double> rhoM, rhoG, rhoW, omegaM, omegaG, omegaW;
576  double piM;
577 
578  // Form factor.
579  complex F(double s, vector<double> M, vector<double> G, vector<double> W);
580 
581 };
582 
583 //==========================================================================
584 
585 // Helicity matrix element for a tau decay into four pions.
586 
588 
589 public:
590 
591  HMETau2FourPions() : a1M(), a1G(), rhoM(), rhoG(), sigM(), sigG(), omeM(),
592  omeG(), picM(), pinM(), sigA(), sigP(), omeA(), omeP(), lambda2() {}
593 
594  void initConstants();
595 
596  void initHadronicCurrent(vector<HelicityParticle>& p);
597 
598 private:
599 
600  // G-function form factors (fits).
601  double G(int i, double s);
602 
603  // T-vector functions.
604  Wave4 t1(Wave4&, Wave4&, Wave4&, Wave4&, Wave4&);
605  Wave4 t2(Wave4&, Wave4&, Wave4&, Wave4&, Wave4&);
606  Wave4 t3(Wave4&, Wave4&, Wave4&, Wave4&, Wave4&);
607 
608  // Breit-Wigner denominators for the intermediate mesons.
609  complex a1D(double s);
610  complex rhoD(double s);
611  complex sigD(double s);
612  complex omeD(double s);
613 
614  // Form factors needed for the a1, rho, and omega.
615  double a1FormFactor(double s);
616  double rhoFormFactor1(double s);
617  double rhoFormFactor2(double s);
618  double omeFormFactor(double s);
619 
620  // Masses and widths of the intermediate mesons.
621  double a1M, a1G, rhoM, rhoG, sigM, sigG, omeM, omeG;
622 
623  // Masses for the pions (charged and neutral).
624  double picM, pinM;
625 
626  // Amplitude, phases, and weights for mixing.
627  double sigA, sigP, omeA, omeP;
628  complex sigW, omeW;
629 
630  // Cut-off for a1 form factor.
631  double lambda2;
632 
633 };
634 
635 //==========================================================================
636 
637 // Helicity matrix element for a tau decaying into five pions.
638 
640 
641 public:
642 
643  HMETau2FivePions() : a1M(), a1G(), rhoM(), rhoG(), omegaM(), omegaG(),
644  omegaW(), sigmaM(), sigmaG(), sigmaW() {}
645 
646  void initConstants();
647 
648  void initHadronicCurrent(vector<HelicityParticle>&);
649 
650 private:
651 
652  // Hadronic currents.
653  Wave4 Ja(Wave4 &q, Wave4 &q1, Wave4 &q2, Wave4 &q3, Wave4 &q4, Wave4 &q5);
654  Wave4 Jb(Wave4 &q, Wave4 &q1, Wave4 &q2, Wave4 &q3, Wave4 &q4, Wave4 &q5);
655 
656  // Simplified s-wave Breit-Wigner assuming massless products.
657  complex breitWigner(double s, double M, double G);
658 
659  // Masses and widths of intermediates.
660  double a1M, a1G, rhoM, rhoG, omegaM, omegaG, omegaW, sigmaM, sigmaG, sigmaW;
661 
662 };
663 
664 //==========================================================================
665 
666 // Helicity matrix element for a tau decay into flat phase space.
667 
669 
670 public:
671 
672  void initWaves(vector<HelicityParticle>&) {};
673 
674  complex calculateME(vector<int>) {return 1;}
675 
676  void calculateD(vector<HelicityParticle>&) {};
677 
678  void calculateRho(unsigned int, vector<HelicityParticle>&) {};
679 
680  double decayWeight(vector<HelicityParticle>&) {return 1.0;}
681 
682  double decayWeightMax(vector<HelicityParticle>&) {return 1.0;}
683 
684 };
685 
686 //==========================================================================
687 
688 } // end namespace Pythia8
689 
690 #endif // end Pythia8_HelicityMatrixElements_H
vector< vector< Wave4 > > u
Wave functions.
Definition: HelicityMatrixElements.h:85
std::complex< double > complex
Convenient typedef for double precision complex numbers.
Definition: PythiaComplex.h:17
Helicity matrix element for a tau decay into three mesons (base class).
Definition: HelicityMatrixElements.h:421
Helicity matrix element for a tau decay into three mesons with kaons.
Definition: HelicityMatrixElements.h:506
Definition: HelicityMatrixElements.h:208
vector< GammaMatrix > gamma
Physics matrices.
Definition: HelicityMatrixElements.h:73
Helicity matrix element for a tau decay into three pions.
Definition: HelicityMatrixElements.h:473
void initWaves(vector< HelicityParticle > &)
Initialize wave functions for the helicity matrix element.
Definition: HelicityMatrixElements.h:672
Helicity matrix element for a tau decay into generic three mesons.
Definition: HelicityMatrixElements.h:534
Definition: HelicityMatrixElements.h:155
ParticleData * particleDataPtr
Pointer to particle data.
Definition: HelicityMatrixElements.h:91
virtual complex sBreitWigner(double m0, double m1, double s, double M, double G)
Return an s-wave BreitWigner.
Definition: HelicityMatrixElements.cc:274
void calculateRho(unsigned int, vector< HelicityParticle > &)
Calculate the density matrix for a particle.
Definition: HelicityMatrixElements.h:678
virtual void initPointers(ParticleData *, CoupSM *, Settings *=0)
Initialize the physics matrices and pointers.
Definition: HelicityMatrixElements.cc:21
Definition: HelicityBasics.h:24
vector< int > pID
Particle ID vector.
Definition: HelicityMatrixElements.h:79
Helicity matrix element for the hard process of Z/Z&#39; -> two fermions.
Definition: HelicityMatrixElements.h:272
Helicity matrix element for a tau decaying into a single scalar meson.
Definition: HelicityMatrixElements.h:344
vector< int > pMap
Particle map vector.
Definition: HelicityMatrixElements.h:76
vector< double > pM
Particle mass vector.
Definition: HelicityMatrixElements.h:82
Helicity matrix element for a tau decay into four pions.
Definition: HelicityMatrixElements.h:587
Helicity matrix element for the hard process of W/W&#39; -> two fermions.
Definition: HelicityMatrixElements.h:239
Definition: HelicityMatrixElements.h:375
double decayWeightMax(vector< HelicityParticle > &)
Return the maximum decay weight for the helicity matrix element.
Definition: HelicityMatrixElements.h:682
virtual complex dBreitWigner(double m0, double m1, double s, double M, double G)
Return a d-wave BreitWigner.
Definition: HelicityMatrixElements.cc:300
CoupSM * coupSMPtr
Pointer to Standard Model constants.
Definition: HelicityMatrixElements.h:97
double decayWeight(vector< HelicityParticle > &)
Calculate the matrix element weight for a decay.
Definition: HelicityMatrixElements.h:680
Definition: HelicityBasics.h:182
virtual complex F1()
Form factors.
Definition: HelicityMatrixElements.h:453
virtual void initConstants()
Initialize the constants for the matrix element (called by initChannel).
Definition: HelicityMatrixElements.h:88
Definition: StandardModel.h:135
Helicity matrix element for a tau decay into two pions and a photon.
Definition: HelicityMatrixElements.h:560
Helicity matrix element for the hard process of X -> two fermions.
Definition: HelicityMatrixElements.h:227
vector< double > rhoM
Resonance masses, widths, and weights.
Definition: HelicityMatrixElements.h:575
virtual complex pBreitWigner(double m0, double m1, double s, double M, double G)
Return a p-wave BreitWigner.
Definition: HelicityMatrixElements.cc:287
virtual double decayWeightMax(vector< HelicityParticle > &)
Calculate the maximum matrix element decay weight.
Definition: HelicityMatrixElements.h:43
Helicity matrix element for a tau decaying into five pions.
Definition: HelicityMatrixElements.h:639
virtual double decayWeight(vector< HelicityParticle > &)
Calculate the matrix element weight for a decay.
Definition: HelicityMatrixElements.cc:163
virtual void calculateD(vector< HelicityParticle > &)
Calculate the decay matrix for a particle.
Definition: HelicityMatrixElements.cc:54
Helicity matrix element for a tau decaying into two leptons.
Definition: HelicityMatrixElements.h:358
double DECAYWEIGHTMAX
Maximum decay weight. WARNING: hardcoded constant.
Definition: HelicityMatrixElements.h:70
virtual complex calculateME(vector< int >)
Calculate the helicity matrix element.
Definition: HelicityMatrixElements.h:47
virtual complex breitWigner(double s, double M, double G)
Calculate Breit-Wigner&#39;s with running widths and fixed.
Definition: HelicityMatrixElements.cc:264
HelicityMatrixElement()
Constructor and destructor.
Definition: HelicityMatrixElements.h:29
virtual HelicityMatrixElement * initChannel(vector< HelicityParticle > &)
Initialize the channel.
Definition: HelicityMatrixElements.cc:36
Header for classes to set beam momentum and interaction vertex spread.
Definition: Analysis.h:20
double m(const Vec4 &v1)
Invariant mass and its square.
Definition: Basics.cc:595
void setFermionLine(int, HelicityParticle &, HelicityParticle &)
Set a fermion line.
Definition: HelicityMatrixElements.cc:239
Settings * settingsPtr
Pointer to Settings.
Definition: HelicityMatrixElements.h:100
Helicity matrix element for the hard process of photon -> two fermions.
Definition: HelicityMatrixElements.h:260
virtual void calculateRho(unsigned int, vector< HelicityParticle > &)
Calculate the density matrix for a particle.
Definition: HelicityMatrixElements.cc:103
complex calculateME(vector< int >)
Return element for the helicity matrix element.
Definition: HelicityMatrixElements.h:674
Helicity matrix element for a tau decay into flat phase space.
Definition: HelicityMatrixElements.h:668
virtual void initWaves(vector< HelicityParticle > &)
Initialize the wave functions (called by decayWeight and calculateRho/D).
Definition: HelicityMatrixElements.h:91
Mode
Decay mode of the tau.
Definition: HelicityMatrixElements.h:434
This class holds a map of all ParticleDataEntries.
Definition: ParticleData.h:422
void calculateD(vector< HelicityParticle > &)
Calculate the decay matrix for a particle.
Definition: HelicityMatrixElements.h:676
The helicity matrix element class.
Definition: HelicityMatrixElements.h:24
Definition: HelicityMatrixElements.h:396
Definition: HelicityMatrixElements.h:131
Definition: Basics.h:32
Helicity matrix element for the decay of a Higgs -> two fermions.
Definition: HelicityMatrixElements.h:300
Base class for all tau decay helicity matrix elements.
Definition: HelicityMatrixElements.h:321
Definition: Settings.h:195