StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
ResonanceWidths.h
1 // ResonanceWidths.h is a part of the PYTHIA event generator.
2 // Copyright (C) 2014 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 resonance properties: dynamical widths etc.
7 // ResonanceWidths: base class for all resonances.
8 // ResonanceGmZ, ...: derived classes for individual resonances.
9 
10 #ifndef Pythia8_ResonanceWidths_H
11 #define Pythia8_ResonanceWidths_H
12 
13 #include "Pythia8/Basics.h"
14 #include "Pythia8/Info.h"
15 #include "Pythia8/PythiaStdlib.h"
16 #include "Pythia8/Settings.h"
17 #include "Pythia8/StandardModel.h"
18 
19 namespace Pythia8 {
20 
21 //==========================================================================
22 
23 // Forward references to ParticleData and StandardModel classes.
24 class DecayChannel;
25 class ParticleData;
26 class ParticleDataEntry;
27 class Couplings;
28 
29 //==========================================================================
30 
31 // The ResonanceWidths is the base class. Also used for generic resonaces.
32 
33 class ResonanceWidths {
34 
35 public:
36 
37  // Destructor.
38  virtual ~ResonanceWidths() {}
39 
40  // Set up standard properties.
41  void initBasic(int idResIn, bool isGenericIn = false) {
42  idRes = idResIn; isGeneric = isGenericIn;}
43 
44  // Calculate and store partial and total widths at the nominal mass.
45  virtual bool init(Info* infoPtrIn, Settings* settingsPtrIn,
46  ParticleData* particleDataPtrIn, Couplings* couplingsPtrIn);
47 
48  // Return identity of particle species.
49  int id() const {return idRes;}
50 
51  // Calculate the total/open width for given mass, charge and instate.
52  double width(int idSgn, double mHatIn, int idInFlavIn = 0,
53  bool openOnly = false, bool setBR = false, int idOutFlav1 = 0,
54  int idOutFlav2 = 0);
55 
56  // Special case to calculate open final-state width.
57  double widthOpen(int idSgn, double mHatIn, int idIn = 0) {
58  return width( idSgn, mHatIn, idIn, true, false);}
59 
60  // Special case to store open final-state widths for channel selection.
61  double widthStore(int idSgn, double mHatIn, int idIn = 0) {
62  return width( idSgn, mHatIn, idIn, true, true);}
63 
64  // Return fraction of width open for particle and antiparticle.
65  double openFrac(int idSgn) {return (idSgn > 0) ? openPos : openNeg;}
66 
67  // Return forced rescaling factor of resonance width.
68  double widthRescaleFactor() {return forceFactor;}
69 
70  // Special case to calculate one final-state width.
71  // Currently only used for Higgs -> qqbar, g g or gamma gamma.
72  double widthChan(double mHatIn, int idOutFlav1, int idOutFlav2) {
73  return width( 1, mHatIn, 0, false, false, idOutFlav1, idOutFlav2);}
74 
75 protected:
76 
77  // Constructor.
78  ResonanceWidths() {}
79 
80  // Constants: could only be changed in the code itself.
81  static const int NPOINT;
82  static const double MASSMIN, MASSMARGIN;
83 
84  // Particle properties always present.
85  int idRes, hasAntiRes;
86  bool doForceWidth, isGeneric, allowCalcWidth;
87  double minWidth, minThreshold, mRes, GammaRes, m2Res, GamMRat,
88  openPos, openNeg, forceFactor;
89 
90  // Properties for currently studied decay channel(s).
91  int iChannel, onMode, meMode, mult, id1, id2, id3, id1Abs,
92  id2Abs, id3Abs, idInFlav;
93  double widNow, mHat, mf1, mf2, mf3, mr1, mr2, mr3, ps, kinFac,
94  alpEM, alpS, colQ, preFac;
95 
96  // Pointer to properties of the particle species.
97  ParticleDataEntry* particlePtr;
98 
99  // Pointer to various information on the generation.
100  Info* infoPtr;
101 
102  // Pointer to the settings database.
103  Settings* settingsPtr;
104 
105  // Pointer to the particle data table.
106  ParticleData* particleDataPtr;
107 
108  // Pointers to Standard Model and SUSY couplings.
109  Couplings* couplingsPtr;
110 
111  // Initialize constants.
112  virtual void initConstants() {}
113 
114  // Virtual methods to handle model-specific (non-SM) part of initialization
115  // for use by derived classes that implement additional models (eg SUSY).
116  virtual bool initBSM() {return true;}
117  virtual bool allowCalc() {return true;}
118 
119  // Calculate various common prefactors for the current mass.
120  // Optional argument calledFromInit only used for Z0.
121  virtual void calcPreFac(bool = false) {}
122 
123  // Calculate width for currently considered channel.
124  // Optional argument calledFromInit only used for Z0.
125  virtual void calcWidth(bool = false) {}
126 
127  // Simple routines for matrix-element integration over Breit-Wigners.
128  double numInt1BW(double mHatIn, double m1, double Gamma1, double mMin1,
129  double m2, int psMode = 1);
130  double numInt2BW(double mHatIn, double m1, double Gamma1, double mMin1,
131  double m2, double Gamma2, double mMin2, int psMode = 1);
132 
133 };
134 
135 //==========================================================================
136 
137 // The ResonanceGeneric class handles a generic resonance.
138 // Only needs a constructor and allowCalc = false; for the rest uses
139 // defaults in base class.
140 
141 class ResonanceGeneric : public ResonanceWidths {
142 
143 public:
144 
145  // Constructor.
146  ResonanceGeneric(int idResIn) {initBasic(idResIn, true);}
147 
148  // By default, assume no dedicated code exists to compute width.
149  virtual bool allowCalc() {return false;}
150 
151 };
152 
153 //==========================================================================
154 
155 // The ResonanceGmZ class handles the gamma*/Z0 resonance.
156 
157 class ResonanceGmZ : public ResonanceWidths {
158 
159 public:
160 
161  // Constructor.
162  ResonanceGmZ(int idResIn) {initBasic(idResIn);}
163 
164 private:
165 
166  // Locally stored properties and couplings.
167  int gmZmode;
168  double thetaWRat, ei2, eivi, vi2ai2, gamNorm, intNorm, resNorm;
169 
170  // Initialize constants.
171  virtual void initConstants();
172 
173  // Calculate various common prefactors for the current mass.
174  virtual void calcPreFac(bool = false);
175 
176  // Caclulate width for currently considered channel.
177  virtual void calcWidth(bool calledFromInit = false);
178 
179 };
180 
181 //==========================================================================
182 
183 // The ResonanceW class handles the W+- resonance.
184 
185 class ResonanceW : public ResonanceWidths {
186 
187 public:
188 
189  // Constructor.
190  ResonanceW(int idResIn) {initBasic(idResIn);}
191 
192 private:
193 
194  // Locally stored properties and couplings.
195  double thetaWRat, alpEM;
196 
197  // Initialize constants.
198  virtual void initConstants();
199 
200  // Calculate various common prefactors for the current mass.
201  virtual void calcPreFac(bool = false);
202 
203  // Caclulate width for currently considered channel.
204  virtual void calcWidth(bool = false);
205 
206 };
207 
208 //==========================================================================
209 
210 // The ResonanceTop class handles the top/antitop resonance.
211 
212 class ResonanceTop : public ResonanceWidths {
213 
214 public:
215 
216  // Constructor.
217  ResonanceTop(int idResIn) {initBasic(idResIn);}
218 
219 private:
220 
221  // Locally stored properties and couplings.
222  double thetaWRat, m2W, tanBeta, tan2Beta, mbRun;
223 
224  // Initialize constants.
225  virtual void initConstants();
226 
227  // Calculate various common prefactors for the current mass.
228  virtual void calcPreFac(bool = false);
229 
230  // Caclulate width for currently considered channel.
231  virtual void calcWidth(bool = false);
232 
233 };
234 
235 //==========================================================================
236 
237 // The ResonanceFour class handles fourth-generation resonances.
238 
239 class ResonanceFour : public ResonanceWidths {
240 
241 public:
242 
243  // Constructor.
244  ResonanceFour(int idResIn) {initBasic(idResIn);}
245 
246 private:
247 
248  // Locally stored properties and couplings.
249  double thetaWRat, m2W;
250 
251  // Initialize constants.
252  virtual void initConstants();
253 
254  // Calculate various common prefactors for the current mass.
255  virtual void calcPreFac(bool = false);
256 
257  // Caclulate width for currently considered channel.
258  virtual void calcWidth(bool = false);
259 
260 };
261 
262 //==========================================================================
263 
264 // The ResonanceH class handles the SM and BSM Higgs resonance.
265 // higgsType = 0 : SM H; = 1: h^0/H_1; = 2 : H^0/H_2; = 3 : A^0/A_3.
266 
267 class ResonanceH : public ResonanceWidths {
268 
269 public:
270 
271  // Constructor.
272  ResonanceH(int higgsTypeIn, int idResIn) : higgsType(higgsTypeIn)
273  {initBasic(idResIn);}
274 
275 private:
276 
277  // Constants: could only be changed in the code itself.
278  static const double MASSMINWZ, MASSMINT, GAMMAMARGIN;
279 
280  // Higgs type in current instance.
281  int higgsType;
282 
283  // Locally stored properties and couplings.
284  bool useCubicWidth, useRunLoopMass, useNLOWidths;
285  double sin2tW, cos2tW, mT, mZ, mW, mHchg, GammaT, GammaZ, GammaW,
286  rescAlpS, rescColQ, coup2d, coup2u, coup2l, coup2Z, coup2W,
287  coup2Hchg, coup2H1H1, coup2A3A3, coup2H1Z, coup2A3Z, coup2A3H1,
288  coup2HchgW, mLowT, mStepT, mLowZ, mStepZ, mLowW, mStepW,
289  kinFacT[101], kinFacZ[101], kinFacW[101];
290 
291  // Initialize constants.
292  virtual void initConstants();
293 
294  // Calculate various common prefactors for the current mass.
295  virtual void calcPreFac(bool = false);
296 
297  // Caclulate width for currently considered channel.
298  virtual void calcWidth(bool = false);
299 
300  // Sum up loop contributions in Higgs -> g + g.
301  double eta2gg();
302 
303  // Sum up loop contributions in Higgs -> gamma + gamma.
304  double eta2gaga();
305 
306  // Sum up loop contributions in Higgs -> gamma + Z0.
307  double eta2gaZ();
308 
309 };
310 
311 //==========================================================================
312 
313 // The ResonanceHchg class handles the H+- resonance.
314 
315 class ResonanceHchg : public ResonanceWidths {
316 
317 public:
318 
319  // Constructor.
320  ResonanceHchg(int idResIn) {initBasic(idResIn);}
321 
322 private:
323 
324  // Locally stored properties and couplings.
325  bool useCubicWidth;
326  double thetaWRat, mW, tanBeta, tan2Beta, coup2H1W;
327 
328  // Initialize constants.
329  virtual void initConstants();
330 
331  // Calculate various common prefactors for the current mass.
332  virtual void calcPreFac(bool = false);
333 
334  // Caclulate width for currently considered channel.
335  virtual void calcWidth(bool = false);
336 
337 };
338 
339 //==========================================================================
340 
341 // The ResonanceZprime class handles the gamma*/Z0 /Z'^0 resonance.
342 
343 class ResonanceZprime : public ResonanceWidths {
344 
345 public:
346 
347  // Constructor.
348  ResonanceZprime(int idResIn) {initBasic(idResIn);}
349 
350 private:
351 
352  // Locally stored properties and couplings.
353  int gmZmode;
354  double sin2tW, cos2tW, thetaWRat, mZ, GammaZ, m2Z, GamMRatZ, afZp[20],
355  vfZp[20], coupZpWW, ei2, eivi, vai2, eivpi, vaivapi, vapi2,
356  gamNorm, gamZNorm, ZNorm, gamZpNorm, ZZpNorm, ZpNorm;
357 
358  // Initialize constants.
359  virtual void initConstants();
360 
361  // Calculate various common prefactors for the current mass.
362  virtual void calcPreFac(bool = false);
363 
364  // Caclulate width for currently considered channel.
365  virtual void calcWidth(bool calledFromInit = false);
366 
367 };
368 
369 //==========================================================================
370 
371 // The ResonanceWprime class handles the W'+- resonance.
372 
373 class ResonanceWprime : public ResonanceWidths {
374 
375 public:
376 
377  // Constructor.
378  ResonanceWprime(int idResIn) {initBasic(idResIn);}
379 
380 private:
381 
382  // Locally stored properties and couplings.
383  double thetaWRat, cos2tW, alpEM, aqWp, vqWp, alWp, vlWp, coupWpWZ;
384 
385  // Initialize constants.
386  virtual void initConstants();
387 
388  // Calculate various common prefactors for the current mass.
389  virtual void calcPreFac(bool = false);
390 
391  // Caclulate width for currently considered channel.
392  virtual void calcWidth(bool = false);
393 
394 };
395 
396 //==========================================================================
397 
398 // The ResonanceRhorizontal class handles the R^0 resonance.
399 
400 class ResonanceRhorizontal : public ResonanceWidths {
401 
402 public:
403 
404  // Constructor.
405  ResonanceRhorizontal(int idResIn) {initBasic(idResIn);}
406 
407 private:
408 
409  // Locally stored properties and couplings.
410  double thetaWRat;
411 
412  // Initialize constants.
413  virtual void initConstants();
414 
415  // Calculate various common prefactors for the current mass.
416  virtual void calcPreFac(bool = false);
417 
418  // Caclulate width for currently considered channel.
419  virtual void calcWidth(bool = false);
420 
421 };
422 
423 //==========================================================================
424 
425 // The ResonanceExcited class handles excited-fermion resonances.
426 
427 class ResonanceExcited : public ResonanceWidths {
428 
429 public:
430 
431  // Constructor.
432  ResonanceExcited(int idResIn) {initBasic(idResIn);}
433 
434 private:
435 
436  // Locally stored properties and couplings.
437  double Lambda, coupF, coupFprime, coupFcol, sin2tW, cos2tW;
438 
439  // Initialize constants.
440  virtual void initConstants();
441 
442  // Calculate various common prefactors for the current mass.
443  virtual void calcPreFac(bool = false);
444 
445  // Caclulate width for currently considered channel.
446  virtual void calcWidth(bool = false);
447 
448 };
449 
450 //==========================================================================
451 
452 // The ResonanceGraviton class handles the excited Graviton resonance.
453 
454 class ResonanceGraviton : public ResonanceWidths {
455 
456 public:
457 
458  // Constructor.
459  ResonanceGraviton(int idResIn) {initBasic(idResIn);}
460 
461 private:
462 
463  // Locally stored properties and couplings.
464  bool eDsmbulk, eDvlvl;
465  double kappaMG;
466 
467  // Couplings between graviton and SM (map from particle id to coupling).
468  double eDcoupling[27];
469 
470  // Initialize constants.
471  virtual void initConstants();
472 
473  // Calculate various common prefactors for the current mass.
474  virtual void calcPreFac(bool = false);
475 
476  // Caclulate width for currently considered channel.
477  virtual void calcWidth(bool = false);
478 
479 };
480 
481 //==========================================================================
482 
483 // The ResonanceKKgluon class handles the g^*/KK-gluon^* resonance.
484 
485 class ResonanceKKgluon : public ResonanceWidths {
486 
487 public:
488 
489  // Constructor.
490  ResonanceKKgluon(int idResIn) {initBasic(idResIn);}
491 
492 private:
493 
494  // Locally stored properties.
495  double normSM, normInt, normKK;
496 
497  // Couplings between kk gluon and SM (indexed by particle id).
498  // Helicity dependent couplings. Use vector/axial-vector
499  // couplings internally, gv/ga = 0.5 * (gL +/- gR).
500  double eDgv[10], eDga[10];
501 
502  // Interference parameter.
503  int interfMode;
504 
505  // Initialize constants.
506  virtual void initConstants();
507 
508  // Calculate various common prefactors for the current mass.
509  virtual void calcPreFac(bool calledFromInit = false);
510 
511  // Caclulate width for currently considered channel.
512  virtual void calcWidth(bool calledFromInit = false);
513 
514 };
515 
516 //==========================================================================
517 
518 // The ResonanceLeptoquark class handles the LQ/LQbar resonance.
519 
520 class ResonanceLeptoquark : public ResonanceWidths {
521 
522 public:
523 
524  // Constructor.
525  ResonanceLeptoquark(int idResIn) {initBasic(idResIn);}
526 
527 private:
528 
529  // Locally stored properties and couplings.
530  double kCoup;
531 
532  // Initialize constants.
533  virtual void initConstants();
534 
535  // Calculate various common prefactors for the current mass.
536  virtual void calcPreFac(bool = false);
537 
538  // Caclulate width for currently considered channel.
539  virtual void calcWidth(bool = false);
540 
541 };
542 
543 //==========================================================================
544 
545 // The ResonanceNuRight class handles righthanded Majorana neutrinos.
546 
547 class ResonanceNuRight : public ResonanceWidths {
548 
549 public:
550 
551  // Constructor.
552  ResonanceNuRight(int idResIn) {initBasic(idResIn);}
553 
554 private:
555 
556  // Locally stored properties and couplings.
557  double thetaWRat, mWR;
558 
559  // Initialize constants.
560  virtual void initConstants();
561 
562  // Calculate various common prefactors for the current mass.
563  virtual void calcPreFac(bool = false);
564 
565  // Caclulate width for currently considered channel.
566  virtual void calcWidth(bool = false);
567 
568 };
569 
570 //==========================================================================
571 
572 // The ResonanceZRight class handles the Z_R^0 resonance.
573 
574 class ResonanceZRight : public ResonanceWidths {
575 
576 public:
577 
578  // Constructor.
579  ResonanceZRight(int idResIn) {initBasic(idResIn);}
580 
581 private:
582 
583  // Locally stored properties and couplings.
584  double sin2tW, thetaWRat;
585 
586  // Initialize constants.
587  virtual void initConstants();
588 
589  // Calculate various common prefactors for the current mass.
590  virtual void calcPreFac(bool = false);
591 
592  // Caclulate width for currently considered channel.
593  virtual void calcWidth(bool = false);
594 
595 };
596 
597 //==========================================================================
598 
599 // The ResonanceWRight class handles the W_R+- resonance.
600 
601 class ResonanceWRight : public ResonanceWidths {
602 
603 public:
604 
605  // Constructor.
606  ResonanceWRight(int idResIn) {initBasic(idResIn);}
607 
608 private:
609 
610  // Locally stored properties and couplings.
611  double thetaWRat;
612 
613  // Initialize constants.
614  virtual void initConstants();
615 
616  // Calculate various common prefactors for the current mass.
617  virtual void calcPreFac(bool = false);
618 
619  // Caclulate width for currently considered channel.
620  virtual void calcWidth(bool = false);
621 
622 };
623 
624 //==========================================================================
625 
626 // The ResonanceHchgchgLeft class handles the H++/H-- (left) resonance.
627 
628 class ResonanceHchgchgLeft : public ResonanceWidths {
629 
630 public:
631 
632  // Constructor.
633  ResonanceHchgchgLeft(int idResIn) {initBasic(idResIn);}
634 
635 private:
636 
637  // Locally stored properties and couplings.
638  double yukawa[4][4], gL, vL, mW;
639 
640  // Initialize constants.
641  virtual void initConstants();
642 
643  // Calculate various common prefactors for the current mass.
644  virtual void calcPreFac(bool = false);
645 
646  // Caclulate width for currently considered channel.
647  virtual void calcWidth(bool = false);
648 
649 };
650 
651 //==========================================================================
652 
653 // The ResonanceHchgchgRight class handles the H++/H-- (right) resonance.
654 
655 class ResonanceHchgchgRight : public ResonanceWidths {
656 
657 public:
658 
659  // Constructor.
660  ResonanceHchgchgRight(int idResIn) {initBasic(idResIn);}
661 
662 private:
663 
664  // Locally stored properties and couplings.
665  int idWR;
666  double yukawa[4][4], gR;
667 
668  // Initialize constants.
669  virtual void initConstants();
670 
671  // Calculate various common prefactors for the current mass.
672  virtual void calcPreFac(bool = false);
673 
674  // Caclulate width for currently considered channel.
675  virtual void calcWidth(bool = false);
676 
677 };
678 
679 //==========================================================================
680 
681 } // end namespace Pythia8
682 
683 #endif // Pythia8_ResonanceWidths_H