Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-01-18 10:06:30

0001 // SigmaLeftRightSym.h is a part of the PYTHIA event generator.
0002 // Copyright (C) 2024 Torbjorn Sjostrand.
0003 // PYTHIA is licenced under the GNU GPL v2 or later, see COPYING for details.
0004 // Please respect the MCnet Guidelines, see GUIDELINES for details.
0005 
0006 // Header file for left-rights-symmetry differential cross sections.
0007 // Contains classes derived from SigmaProcess via Sigma(1/2/3)Process.
0008 
0009 #ifndef Pythia8_SigmaLeftRightSym_H
0010 #define Pythia8_SigmaLeftRightSym_H
0011 
0012 #include "Pythia8/SigmaProcess.h"
0013 
0014 namespace Pythia8 {
0015 
0016 //==========================================================================
0017 
0018 // A derived class for f fbar -> Z_R^0 (righthanded gauge boson).
0019 
0020 class Sigma1ffbar2ZRight : public Sigma1Process {
0021 
0022 public:
0023 
0024   // Constructor.
0025   Sigma1ffbar2ZRight() : idZR(), mRes(), GammaRes(), m2Res(), GamMRat(),
0026     sin2tW(), sigma0(), ZRPtr() {}
0027 
0028   // Initialize process.
0029   virtual void initProc();
0030 
0031   // Calculate flavour-independent parts of cross section.
0032   virtual void sigmaKin();
0033 
0034   // Evaluate sigmaHat(sHat).
0035   virtual double sigmaHat();
0036 
0037   // Select flavour, colour and anticolour.
0038   virtual void setIdColAcol();
0039 
0040   // Evaluate weight for G* decay angle.
0041   virtual double weightDecay( Event& process, int iResBeg, int iResEnd);
0042 
0043   // Info on the subprocess.
0044   virtual string name()       const {return "f fbar -> Z_R^0";}
0045   virtual int    code()       const {return 3101;}
0046   virtual string inFlux()     const {return "ffbarSame";}
0047   virtual int    resonanceA() const {return idZR;}
0048 
0049 private:
0050 
0051   // Parameters set at initialization or for current kinematics.
0052   int    idZR;
0053   double mRes, GammaRes, m2Res, GamMRat, sin2tW, sigma0;
0054 
0055   // Pointer to properties of the particle species, to access decay channels.
0056   ParticleDataEntryPtr ZRPtr;
0057 
0058 };
0059 
0060 //==========================================================================
0061 
0062 // A derived class for f fbar' -> W_R^+- (righthanded gauge boson).
0063 
0064 class Sigma1ffbar2WRight : public Sigma1Process {
0065 
0066 public:
0067 
0068   // Constructor.
0069   Sigma1ffbar2WRight() : idWR(), mRes(), GammaRes(), m2Res(), GamMRat(),
0070     thetaWRat(), sigma0Pos(), sigma0Neg(), particlePtr() {}
0071 
0072   // Initialize process.
0073   virtual void initProc();
0074 
0075   // Calculate flavour-independent parts of cross section.
0076   virtual void sigmaKin();
0077 
0078   // Evaluate sigmaHat(sHat).
0079   virtual double sigmaHat();
0080 
0081   // Select flavour, colour and anticolour.
0082   virtual void setIdColAcol();
0083 
0084   // Evaluate weight for W decay angle.
0085   virtual double weightDecay( Event& process, int iResBeg, int iResEnd);
0086 
0087   // Info on the subprocess.
0088   virtual string name()       const {return "f fbar' -> W_R^+-";}
0089   virtual int    code()       const {return 3102;}
0090   virtual string inFlux()     const {return "ffbarChg";}
0091   virtual int    resonanceA() const {return idWR;}
0092 
0093 private:
0094 
0095   // Parameters set at initialization.
0096   int    idWR;
0097   double mRes, GammaRes, m2Res, GamMRat, thetaWRat, sigma0Pos, sigma0Neg;
0098 
0099   // Pointer to properties of the particle species, to access decay channels.
0100   ParticleDataEntryPtr particlePtr;
0101 
0102 };
0103 
0104 //==========================================================================
0105 
0106 // A derived class for l l -> H_L^++-- or H_R^++-- (doubly charged Higgs).
0107 
0108 class Sigma1ll2Hchgchg : public Sigma1Process {
0109 
0110 public:
0111 
0112   // Constructor.
0113   Sigma1ll2Hchgchg(int leftRightIn ) : leftRight(leftRightIn), idHLR(),
0114     codeSave(), mRes(), GammaRes(), m2Res(), GamMRat(), yukawa(),
0115     particlePtr() {}
0116 
0117   // Initialize process.
0118   virtual void initProc();
0119 
0120   // Evaluate sigmaHat(sHat).
0121   virtual double sigmaHat();
0122 
0123   // Select flavour, colour and anticolour.
0124   virtual void setIdColAcol();
0125 
0126   // Evaluate weight for W decay angle.
0127   virtual double weightDecay( Event& process, int iResBeg, int iResEnd);
0128 
0129   // Info on the subprocess.
0130   virtual string name()       const {return nameSave;}
0131   virtual int    code()       const {return codeSave;}
0132   virtual string inFlux()     const {return "ff";}
0133   virtual int    resonanceA() const {return idHLR;}
0134 
0135 private:
0136 
0137   // Parameters set at initialization.
0138   int    leftRight, idHLR, codeSave;
0139   string nameSave;
0140   double mRes, GammaRes, m2Res, GamMRat, yukawa[4][4];
0141 
0142   // Pointer to properties of the particle species, to access decay channels.
0143   ParticleDataEntryPtr particlePtr;
0144 
0145 };
0146 
0147 //==========================================================================
0148 
0149 // A derived class for l- gamma -> H_(L/R)^-- l+  (doubly charged Higgs).
0150 
0151 class Sigma2lgm2Hchgchgl : public Sigma2Process {
0152 
0153 public:
0154 
0155   // Constructor.
0156   Sigma2lgm2Hchgchgl(int leftRightIn, int idLepIn ) : leftRight(leftRightIn),
0157     idHLR(), idLep(idLepIn), codeSave(), yukawa(), openFracPos(),
0158     openFracNeg() {}
0159 
0160   // Initialize process.
0161   virtual void initProc();
0162 
0163   // Evaluate sigmaHat(sHat).
0164   virtual double sigmaHat();
0165 
0166   // Select flavour, colour and anticolour.
0167   virtual void setIdColAcol();
0168 
0169   // Evaluate weight for W decay angle.
0170   virtual double weightDecay( Event& process, int iResBeg, int iResEnd);
0171 
0172   // Info on the subprocess.
0173   virtual string name()       const {return nameSave;}
0174   virtual int    code()       const {return codeSave;}
0175   virtual string inFlux()     const {return "fgm";}
0176   virtual int    resonanceA() const {return idHLR;}
0177 
0178 private:
0179 
0180   // Parameters set at initialization.
0181   int    leftRight, idHLR, idLep, codeSave;
0182   string nameSave;
0183   double yukawa[4], openFracPos, openFracNeg;
0184 
0185 };
0186 
0187 //==========================================================================
0188 
0189 // A derived class for f_1 f_2 -> H_(L/R)^++-- f_3 f_4 (W+- W+- fusion).
0190 
0191 class Sigma3ff2HchgchgfftWW : public Sigma3Process {
0192 
0193 public:
0194 
0195   // Constructor.
0196   Sigma3ff2HchgchgfftWW(int leftRightIn) : leftRight(leftRightIn), idHLR(),
0197     codeSave(), mWS(), prefac(), sigma0TU(), sigma0T(), openFracPos(),
0198     openFracNeg() {}
0199 
0200   // Initialize process.
0201   virtual void initProc();
0202 
0203   // Calculate flavour-independent parts of cross section.
0204   virtual void sigmaKin();
0205 
0206   // Evaluate sigmaHat(sHat).
0207   virtual double sigmaHat();
0208 
0209   // Select flavour, colour and anticolour.
0210   virtual void setIdColAcol();
0211 
0212   // Evaluate weight for decay angles.
0213   virtual double weightDecay( Event& process, int iResBeg, int iResEnd);
0214 
0215   // Info on the subprocess.
0216   virtual string name()    const {return nameSave;}
0217   virtual int    code()    const {return codeSave;}
0218   virtual string inFlux()  const {return "ff";}
0219   virtual int    id3Mass() const {return idHLR;}
0220 
0221   // Instructions for 3-body phase space with t-channel propagators.
0222   virtual int    idTchan1()        const {return 9900024;}
0223   virtual int    idTchan2()        const {return 9900024;}
0224   virtual double tChanFracPow1()   const {return 0.05;}
0225   virtual double tChanFracPow2()   const {return 0.9;}
0226   virtual bool   useMirrorWeight() const {return true;}
0227 
0228 private:
0229 
0230   // Store standard prefactor.
0231   int    leftRight, idHLR, codeSave;
0232   string nameSave;
0233   double mWS, prefac, sigma0TU, sigma0T, openFracPos, openFracNeg;
0234 
0235 };
0236 
0237 //==========================================================================
0238 
0239 // A derived class for f fbar -> H_(L/R)^++ H_(L/R)^--  (doubly charged Higgs).
0240 
0241 class Sigma2ffbar2HchgchgHchgchg : public Sigma2Process {
0242 
0243 public:
0244 
0245   // Constructor.
0246   Sigma2ffbar2HchgchgHchgchg(int leftRightIn) : leftRight(leftRightIn),
0247     idHLR(), codeSave(), mRes(), GammaRes(), m2Res(), GamMRat(), sin2tW(),
0248     preFac(), yukawa(), openFrac() {}
0249 
0250   // Initialize process.
0251   virtual void initProc();
0252 
0253   // Evaluate sigmaHat(sHat).
0254   virtual double sigmaHat();
0255 
0256   // Select flavour, colour and anticolour.
0257   virtual void setIdColAcol();
0258 
0259   // Evaluate weight for W decay angle.
0260   virtual double weightDecay( Event& process, int iResBeg, int iResEnd);
0261 
0262   // Info on the subprocess.
0263   virtual string name()       const {return nameSave;}
0264   virtual int    code()       const {return codeSave;}
0265   virtual string inFlux()     const {return "ffbarSame";}
0266   virtual int    id3Mass()    const {return idHLR;}
0267   virtual int    id4Mass()    const {return idHLR;}
0268   virtual int    resonanceA() const {return 23;}
0269 
0270 private:
0271 
0272   // Parameters set at initialization.
0273   int    leftRight, idHLR, codeSave;
0274   string nameSave;
0275   double mRes, GammaRes, m2Res, GamMRat, sin2tW, preFac, yukawa[4][4],
0276          openFrac;
0277 
0278 };
0279 
0280 //==========================================================================
0281 
0282 } // end namespace Pythia8
0283 
0284 #endif // Pythia8_SigmaLeftRightSym_H