herwig is hosted by Hepforge, IPPP Durham
Herwig  7.2.1
MEPP2VVPowheg.h
1 // -*- C++ -*-
2 #ifndef HERWIG_MEPP2VVPowheg_H
3 #define HERWIG_MEPP2VVPowheg_H
4 //
5 // This is the declaration of the MEPP2VVPowheg class.
6 //
7 
8 #include "Herwig/MatrixElement/Hadron/MEPP2VV.h"
9 #include "Herwig/MatrixElement/Powheg/VVKinematics.h"
10 #include "Herwig/Utilities/Maths.h"
11 #include "Herwig/Models/StandardModel/StandardCKM.h"
12 #include "Herwig/Shower/ShowerAlpha.h"
13 
14 namespace Herwig {
15 using namespace ThePEG;
16 using Math::ReLi2;
17 using Constants::pi;
18 
25 class MEPP2VVPowheg: public MEPP2VV {
26 
27 public:
28 
32  MEPP2VVPowheg();
33 
39  virtual POWHEGType hasPOWHEGCorrection() {return ISR;}
40 
44  virtual RealEmissionProcessPtr generateHardest(RealEmissionProcessPtr,
45  ShowerInteraction inter);
47 
48 public:
49 
58  virtual bool generateKinematics(const double * r);
59 
64  virtual int nDim() const;
65 
73  virtual double me2() const;
74 
78  virtual Energy2 scale() const;
79 
86  bool sanityCheck() const;
87 
91  Complex CKM(int ix,int iy) const { return ckm_[ix][iy]; }
92 
93 public:
94 
101  void persistentOutput(PersistentOStream & os) const;
102 
108  void persistentInput(PersistentIStream & is, int version);
110 
117  static void Init();
118 
119 protected:
120 
128  virtual void doinit();
130 
131 public:
132 
136  void getKinematics(double xt, double y, double theta2);
137 
142  double NLOweight() const;
143 
148  double Lhat_ab(tcPDPtr a, tcPDPtr b, realVVKinematics Kinematics) const;
149 
153  double Vtilde_universal(realVVKinematics S) const;
154 
159  double Ctilde_Ltilde_qq_on_x(tcPDPtr a,tcPDPtr b,realVVKinematics C) const;
160 
165  double Ctilde_Ltilde_gq_on_x(tcPDPtr a,tcPDPtr b,realVVKinematics C) const;
166 
171  double Rtilde_Ltilde_qqb_on_x(tcPDPtr a,tcPDPtr b) const;
172 
177  double Rtilde_Ltilde_qg_on_x(tcPDPtr a,tcPDPtr b) const;
178 
183  double Rtilde_Ltilde_gqb_on_x(tcPDPtr a,tcPDPtr b) const;
184 
191  double M_V_regular(realVVKinematics S) const;
192 
200  mutable double M_V_regular_;
201 
205  Energy2 t_u_M_R_qqb(realVVKinematics R) const;
206 
210  mutable Energy2 t_u_M_R_qqb_;
211 
215  Energy2 t_u_M_R_qg(realVVKinematics R) const;
216 
220  mutable Energy2 t_u_M_R_qg_;
221 
225  Energy2 t_u_M_R_gqb(realVVKinematics R) const;
226 
230  mutable Energy2 t_u_M_R_gqb_;
231 
235  Energy2 t_u_M_R_qqb_hel_amp(realVVKinematics R) const;
236 
241  mutable Energy2 t_u_M_R_qqb_hel_amp_;
242 
246  Energy2 t_u_M_R_qg_hel_amp(realVVKinematics R) const;
247 
252  mutable Energy2 t_u_M_R_qg_hel_amp_;
253 
257  Energy2 t_u_M_R_gqb_hel_amp(realVVKinematics R) const;
258 
263  mutable Energy2 t_u_M_R_gqb_hel_amp_;
264 
268  double lo_me() const;
269 
275  double M_Born_WZ(bornVVKinematics B) const;
276 
283  mutable double M_Born_;
284 
290  double M_Born_ZZ(bornVVKinematics B) const;
291 
298  double M_V_regular_ZZ(realVVKinematics S) const;
299 
306  Energy2 t_u_M_R_qqb_ZZ(realVVKinematics R) const;
307 
313  double M_Born_WW(bornVVKinematics B) const;
314 
321  double M_V_regular_WW(realVVKinematics S) const;
322 
329  Energy2 t_u_M_R_qqb_WW(realVVKinematics R) const;
330 
334  Energy2 mu_F2() const;
335 
339  Energy2 mu_UV2() const;
340 
341 protected:
342 
349  inline virtual IBPtr clone() const { return new_ptr(*this); }
350 
355  inline virtual IBPtr fullclone() const { return new_ptr(*this); }
357 
358 private:
359 
364  MEPP2VVPowheg & operator=(const MEPP2VVPowheg &) = delete;
365 
366 private:
367 
372 
377  double tiny;
378 
383  tcBeamPtr hadron_B_;
384 
389 
394 
399 
404 
409 
414 
419 
423  tcPDPtr ab_, bb_;
424 
429  tcPDPtr quark_, antiquark_;
430 
434  double lo_lumi_;
435 
439  mutable double lo_me2_;
440 
444  double CF_;
445 
449  double TR_;
450 
454  double NC_;
455 
459  mutable double gW_, sin2ThetaW_;
460 
464  mutable double guL_, gdL_;
465 
469  mutable double guR_, gdR_;
470 
474  mutable double eZ_;
475 
482  mutable double eZ2_;
483 
487  mutable double Fij2_;
488 
492  unsigned int contrib_;
493 
498  unsigned int channels_;
499 
503  unsigned int nlo_alphaS_opt_;
504 
509 
513  unsigned int removebr_;
514 
518  unsigned int scaleopt_;
519 
523  Energy mu_F_;
524 
528  Energy mu_UV_;
529 
533  Energy2 k1r_perp2_lab_;
534 
538  Energy2 k2r_perp2_lab_;
539 
543  vector< vector<Complex> > ckm_;
544 
549 
554 
559 
565 
573  AbstractFFVVertexPtr FFPvertex_;
574 
578  AbstractFFVVertexPtr FFWvertex_;
579 
583  AbstractFFVVertexPtr FFZvertex_;
584 
588  AbstractVVVVertexPtr WWWvertex_;
589 
593  AbstractFFVVertexPtr FFGvertex_;
595 
599  mutable double alphaS_;
600 
601 protected:
602 
611  double getResult(int emis_type, realVVKinematics R, Energy pT);
612 
619  bool getEvent(vector<Lorentz5Momentum> & pnew,unsigned int & emissiontype);
620 
625  void setTheScales(Energy pT);
626 
630  Energy2 t_u_M_R_qqb_hel_amp(realVVKinematics R, bool getMatrix) const;
631 
632 
636  Energy2 t_u_M_R_qg_hel_amp(realVVKinematics R, bool getMatrix) const;
637 
641  Energy2 t_u_M_R_gqb_hel_amp(realVVKinematics R, bool getMatrix) const;
642 
650  double lo_me(bool getMatrix) const;
651 
655  void recalculateVertex();
656 
661  bool isotropicDecayer();
662 
666  Energy2 triangleFn(Energy2,Energy2,Energy2);
667 
668 private:
669 
676 
680  double lo_me_;
681 
686 
693  unsigned int channel_;
694 
700 
704  ShowerAlphaPtr showerAlphaS_;
705 
714  double power_;
715 
719  double preqqbar_;
720 
724  double preqg_;
725 
729  double pregqbar_;
730 
734  double b0_;
735 
742  Energy LambdaQCD_;
743 
747  vector<double> prefactor_;
749 
758  PPtr qbProgenitor_;
759 
764  PPtr showerAntiquark_;
765 
770  tcBeamPtr qbHadron_;
772 
781  PPtr V1_;
782  PPtr V2_;
783  vector<PPtr> children_;
784 
789  bool flipped_;
790 
794  double Yk_;
795 
799  Energy pT_;
801 
805  Energy min_pT_;
806 
807  // Work out the scales we want to use in the matrix elements and the pdfs:
811  Energy2 QCDScale_;
812 
816  Energy2 PDFScale_;
817 
821  Energy2 EWScale_;
822 
826  mutable Complex productionMatrix_[3][3][3][3];
827 
828 };
829 
830 }
831 
832 #endif /* HERWIG_MEPP2VVPowheg_H */
Ptr< BeamParticleData >::transient_const_pointer tcBeamPtr
A typedef for the BeamParticleData.
Definition: HwMEBase.h:20
ThePEG::Ptr< Particle >::pointer PPtr
realVVKinematics H_
The resolved 2->3 real emission kinematics:
Energy2 QCDScale_
Scale for alpha_S: pT^2 of the diboson system.
double lo_lumi_
Values of the PDF&#39;s before radiation.
ProductionMatrixElement gqb_hel_amps_
The g + qb -> v1 + v2 + qb helicity amplitudes.
std::complex< double > Complex
Energy2 k1r_perp2_lab_
The pT of V1 in a radiative event in the lab frame (for scale setting only)
PPtr showerQuark_
Pointers to the Shower particle objects for the partons.
unsigned int nlo_alphaS_opt_
Whether to use a fixed or a running QCD coupling for the NLO weight.
Energy2 PDFScale_
Scale for real emission PDF:
ShowerInteraction
Handy header file to be included in all Shower classes.
virtual IBPtr fullclone() const
Make a clone of this object, possibly modifying the cloned object to make it sane.
realVVKinematics R_
The resolved 2->3 real emission kinematics.
bool realMESpinCorrelations_
If this boolean is true the n+1 body helicity amplitudes will be used to calculate a hard vertex base...
realVVKinematics SCm_
The collinear limit of the 2->3 kinematics (emission in -z direction).
double tiny
Parameters for the NLO weight.
bornVVKinematics B_
Born / virtual 2->2 kinematics.
int fermionNumberOfMother_
Identifies the space-like mother of the branching as quark (+1) or antiquark (-1): ...
ThePEG::Ptr< InterfacedBase >::pointer IBPtr
virtual IBPtr clone() const
Make a simple clone of this object.
Energy mu_UV_
The renormalization scale.
double TR_
The TR_ colour factor.
Energy2 EWScale_
Scale of electroweak vertices: mVV^2 the invariant mass of the diboson system.
double M_V_regular_
Member variable to store the regular part of the virtual correction matrix element(s).
ProductionMatrixElement qqb_hel_amps_
The q + qb -> v1 + v2 + g helicity amplitudes.
AbstractFFVVertexPtr FFPvertex_
The vertices.
unsigned int channels_
Whether to generate all channels contributions or just qqb or just qg+gqb contributions.
unsigned int contrib_
Whether to generate the positive, negative or leading order contribution.
double gW_
The weak coupling and the sin (squared) of the Weinberg angle.
unsigned int removebr_
Flag to remove or multiply in MCFM branching fractions for testing.
unsigned int scaleopt_
Selects a dynamic (sHat) or fixed factorization scale.
realVVKinematics SCp_
Soft-collinear limit of the 2->3 kinematics (emission in +z direction).
unsigned int channel_
This specifies the emitting configuration: 1: q + qbar -> V1 + V2 + g 2: q + g -> V1 + V2 + q 3: g + ...
tcPDPtr quark_
The ParticleData object for the quark and antiquark (which can be in a different order to ab_ and bb_...
Energy2 t_u_M_R_qqb_
Member variable to store the matrix element q + qb -> n + g times tk*uk.
double NC_
The number of colours.
ShowerAlphaPtr showerAlphaS_
Pointer to the object calculating the strong coupling.
Energy pT_
The transverse momentum of the jet.
double b0_
The QCD beta function divided by 4pi, (11-2/3*nf)/4/pi, with nf = 5.
bool helicityConservation_
Option to impose helicity conservation on the real NLO ME&#39;s (greatly improves evaluation time)...
double Yk_
the rapidity of the jet
virtual POWHEGType hasPOWHEGCorrection()
Has a POWHEG style correction.
Definition: MEPP2VVPowheg.h:39
tcPDPtr ab_
The ParticleData object for the plus and minus lo partons.
The storage of the helicity amplitude expression for the matrix element of a hard process...
double guR_
The up and down, right handed, quark-boson couplings (for WW & ZZ)
PPtr qProgenitor_
Properties of the incoming particles.
Here is the documentation of the MEPP2VVPowheg class.
Definition: MEPP2VVPowheg.h:25
double pregqbar_
The prefactor, for the channel.
Complex CKM(int ix, int iy) const
Return the CKM matrix elements.
Definition: MEPP2VVPowheg.h:91
Energy2 t_u_M_R_qqb_hel_amp_
Member variable to store the matrix element q + qb -> n + g times (tk*uk)^2 - using helicity amplitud...
double preqg_
The prefactor, for the channel.
realVVKinematics S_
Soft limit of the 2->3 real emission kinematics.
Energy mu_F_
The factorization scale.
double CF_
The CF_ colour factor.
vector< double > prefactor_
The prefactors as a vector for easy use.
ThePEG::Ptr< ParticleData >::transient_const_pointer tcPDPtr
realVVKinematics Cm_
The collinear limit of the 2->3 kinematics (emission in -z direction).
double fixed_alphaS_
The value of alphaS to use for the nlo weight if nlo_alphaS_opt_=1.
The realVVKinematics class is used to store information on the the kinematics of the real emission pr...
Definition: VVKinematics.h:261
Energy2 t_u_M_R_qg_hel_amp_
Member variable to store the matrix element q + g -> n + q times (tk*uk)^2 - using helicity amplitude...
PPtr gluon_
Properties of the boson and jets.
Energy2 t_u_M_R_gqb_hel_amp_
Member variable to store the matrix element g + qb -> n + qb times (tk*uk)^2 - using helicity amplitu...
double power_
Constants for the sampling.
POWHEGType
Virtual members to be overridden by inheriting classes which implement hard corrections.
Definition: HwMEBase.h:87
Energy LambdaQCD_
The fundamental QCD scale in the one-loop alpha_{S} used for the crude (not the very crude) overestim...
Energy2 t_u_M_R_gqb_
Member variable to store the matrix element g + qb -> n + q times tk*uk.
double lo_me2_
The value of the leading order qqbar->VV matrix element.
tcBeamPtr qHadron_
Pointers to the BeamParticleData objects.
Energy2 k2r_perp2_lab_
The pT of V2 in a radiative event in the lab frame (for scale setting only)
-*- C++ -*-
Energy min_pT_
The transverse momentum of the jet.
Energy2 t_u_M_R_qg_
Member variable to store the matrix element q + g -> n + q times tk*uk.
vector< vector< Complex > > ckm_
The ckm matrix elements (unsquared, to allow interference)
double M_Born_
Member variable to store the Born matrix element as given in Equation 3.1 - 3.3 in NPB 383 (1992) ***...
AbstractFFVVertexPtr FFGvertex_
The quark-antiquark gluon vertex.
double eZ2_
The TGC coupling squared.
double alphaS_
The value of used for the calculation.
The bornVVKinematics class is used to store information on the the kinematics of the real emission pr...
Definition: VVKinematics.h:18
long double ReLi2(long double)
The real part of the dilog function taken from FORTRAN Herwig.
bool flipped_
Flag indicating if the q & qbar are flipped or not i.e.
AbstractFFVVertexPtr FFZvertex_
The Z fermion-antifermionvertex.
tcBeamPtr hadron_A_
The BeamParticleData object for the plus and minus direction hadrons.
double Fij2_
The CKM factor (Fij^2)
AbstractVVVVertexPtr WWWvertex_
The triple electroweak gauge boson vertex.
double preqqbar_
The prefactor, for the channel.
The MEPP2VV class implements the production of , and in hadron-hadron collisions.
Definition: MEPP2VV.h:24
double lo_me_
The colour & spin averaged n-body (leading order) matrix element squared.
realVVKinematics Cp_
The collinear limit of the 2->3 kinematics (emission in +z direction).
double eZ_
The TGC coupling.
AbstractFFVVertexPtr FFWvertex_
The W fermion-antifermion vertex.
ProductionMatrixElement qg_hel_amps_
The q + g -> v1 + v2 + q helicity amplitudes.
double guL_
The up and down, left handed, quark-boson couplings.