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