herwig is hosted by Hepforge, IPPP Durham
Herwig 7.3.0
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
14namespace Herwig {
15using namespace ThePEG;
16using Math::ReLi2;
17using Constants::pi;
18
25class MEPP2VVPowheg: public MEPP2VV {
26
27public:
28
33
39 virtual POWHEGType hasPOWHEGCorrection() {return ISR;}
40
44 virtual RealEmissionProcessPtr generateHardest(RealEmissionProcessPtr,
45 ShowerInteraction inter);
47
48public:
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
93public:
94
102
108 void persistentInput(PersistentIStream & is, int version);
110
117 static void Init();
118
119protected:
120
128 virtual void doinit();
130
131public:
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
154
160
166
172
178
184
192
200 mutable double M_V_regular_;
201
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
226
230 mutable Energy2 t_u_M_R_gqb_;
231
236
241 mutable Energy2 t_u_M_R_qqb_hel_amp_;
242
247
252 mutable Energy2 t_u_M_R_qg_hel_amp_;
253
258
263 mutable Energy2 t_u_M_R_gqb_hel_amp_;
264
268 double lo_me() const;
269
276
283 mutable double M_Born_;
284
291
299
307
314
322
330
334 Energy2 mu_F2() const;
335
339 Energy2 mu_UV2() const;
340
341protected:
342
349 inline virtual IBPtr clone() const { return new_ptr(*this); }
350
355 inline virtual IBPtr fullclone() const { return new_ptr(*this); }
357
358private:
359
365
366private:
367
372
377 double tiny;
378
383 tcBeamPtr hadron_B_;
384
389
394
399
404
409
414
419
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
534
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
601protected:
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
656
662
666 Energy2 triangleFn(Energy2,Energy2,Energy2);
667
668private:
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
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
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 */
POWHEGType
Virtual members to be overridden by inheriting classes which implement hard corrections.
Definition: HwMEBase.h:87
Here is the documentation of the MEPP2VVPowheg class.
Definition: MEPP2VVPowheg.h:25
double M_V_regular_
Member variable to store the regular part of the virtual correction matrix element(s).
AbstractFFVVertexPtr FFPvertex_
The vertices.
double Yk_
the rapidity of the jet
virtual POWHEGType hasPOWHEGCorrection()
Has a POWHEG style correction.
Definition: MEPP2VVPowheg.h:39
double alphaS_
The value of used for the calculation.
PPtr gluon_
Properties of the boson and jets.
Energy2 t_u_M_R_gqb(realVVKinematics R) const
The matrix element g + qb -> n + q times tk*uk.
Energy2 t_u_M_R_gqb_hel_amp(realVVKinematics R) const
The matrix element g + qb -> n + qb times (tk*uk)^2 - using helicity amplitudes.
Energy2 t_u_M_R_qqb_WW(realVVKinematics R) const
t_u_M_R_qqb_WW is the q + qb -> n + g times tk*uk real emission matrix element as defined in Eq.
double Lhat_ab(tcPDPtr a, tcPDPtr b, realVVKinematics Kinematics) const
Calculate the ratio of the NLO luminosity to the LO luminosity function for the initiated channel.
double tiny
Parameters for the NLO weight.
void recalculateVertex()
Recalculate hard vertex to include spin correlations for radiative events.
double M_Born_WZ(bornVVKinematics B) const
The Born matrix element as given in Equation 3.1 - 3.3 in NPB 383 (1992) *** modulo a factor 1/(2s) *...
realVVKinematics SCm_
The collinear limit of the 2->3 kinematics (emission in -z direction).
Energy2 t_u_M_R_qqb_ZZ(realVVKinematics R) const
t_u_M_R_qqb_ZZ is the q + qb -> n + g times tk*uk real emission matrix element as defined in Eq.
Energy2 t_u_M_R_qqb_hel_amp(realVVKinematics R) const
The matrix element q + qb -> n + g times (tk*uk)^2 - using helicity amplitudes.
double Ctilde_Ltilde_gq_on_x(tcPDPtr a, tcPDPtr b, realVVKinematics C) const
Function for calculation of the initiated real contribution.
double preqg_
The prefactor, for the channel.
double pregqbar_
The prefactor, for the channel.
Complex productionMatrix_[3][3][3][3]
A matrix to hold the home-grown production matrix element.
static void Init()
The standard Init function used to initialize the interfaces.
Energy2 triangleFn(Energy2, Energy2, Energy2)
The triangle function lambda(x,y,z)=sqrt(x^2+y^2+z^2-2*x*y-2*y*z-2*x*z)
Energy2 k1r_perp2_lab_
The pT of V1 in a radiative event in the lab frame (for scale setting only)
double Rtilde_Ltilde_qqb_on_x(tcPDPtr a, tcPDPtr b) const
Function for calculation of the initiated real contribution.
Energy2 mu_UV2() const
Return the renormalisation scale squared.
ShowerAlphaPtr showerAlphaS_
Pointer to the object calculating the strong coupling.
Energy2 mu_F2() const
Return the factorion scale squared.
double TR_
The TR_ colour factor.
double Fij2_
The CKM factor (Fij^2)
unsigned int removebr_
Flag to remove or multiply in MCFM branching fractions for testing.
Energy mu_UV_
The renormalization scale.
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...
unsigned int contrib_
Whether to generate the positive, negative or leading order contribution.
PPtr qProgenitor_
Properties of the incoming particles.
double Rtilde_Ltilde_qg_on_x(tcPDPtr a, tcPDPtr b) const
Function for calculation of the initiated real contribution.
ProductionMatrixElement qg_hel_amps_
The q + g -> v1 + v2 + q helicity amplitudes
AbstractVVVVertexPtr WWWvertex_
The triple electroweak gauge boson vertex.
virtual Energy2 scale() const
Return the scale associated with the last set phase space point.
Energy2 QCDScale_
Scale for alpha_S: pT^2 of the diboson system.
bool flipped_
Flag indicating if the q & qbar are flipped or not i.e.
vector< vector< Complex > > ckm_
The ckm matrix elements (unsquared, to allow interference)
virtual bool generateKinematics(const double *r)
Generate internal degrees of freedom given nDim() uniform random numbers in the interval .
double M_V_regular(realVVKinematics S) const
The regular part of the virtual correction matrix element(s).
tcBeamPtr hadron_A_
The BeamParticleData object for the plus and minus direction hadrons.
virtual IBPtr clone() const
Make a simple clone of this object.
AbstractFFVVertexPtr FFGvertex_
The quark-antiquark gluon vertex.
AbstractFFVVertexPtr FFWvertex_
The W fermion-antifermion vertex.
Energy2 k2r_perp2_lab_
The pT of V2 in a radiative event in the lab frame (for scale setting only)
unsigned int nlo_alphaS_opt_
Whether to use a fixed or a running QCD coupling for the NLO weight.
Energy mu_F_
The factorization scale.
virtual RealEmissionProcessPtr generateHardest(RealEmissionProcessPtr, ShowerInteraction inter)
Apply the POWHEG style correction.
unsigned int channels_
Whether to generate all channels contributions or just qqb or just qg+gqb contributions.
Energy2 t_u_M_R_qqb_
Member variable to store the matrix element q + qb -> n + g times tk*uk.
virtual int nDim() const
The number of internal degrees of freedom used in the matrix element.
void getKinematics(double xt, double y, double theta2)
Function to set the born variables.
Complex CKM(int ix, int iy) const
Return the CKM matrix elements.
Definition: MEPP2VVPowheg.h:91
double guR_
The up and down, right handed, quark-boson couplings (for WW & ZZ)
ProductionMatrixElement qqb_hel_amps_
The q + qb -> v1 + v2 + g helicity amplitudes
unsigned int scaleopt_
Selects a dynamic (sHat) or fixed factorization scale.
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...
PPtr showerQuark_
Pointers to the Shower particle objects for the partons.
realVVKinematics H_
The resolved 2->3 real emission kinematics:
ProductionMatrixElement gqb_hel_amps_
The g + qb -> v1 + v2 + qb helicity amplitudes
tcPDPtr ab_
The ParticleData object for the plus and minus lo partons.
Energy2 t_u_M_R_gqb_hel_amp(realVVKinematics R, bool getMatrix) const
The matrix element g + qb -> n + q times tk*uk.
double M_Born_WW(bornVVKinematics B) const
The Born matrix element as given in Equation 3.2 - 3.8 in NPB 410 (1993) *** modulo a factor 1/(2s) *...
bool helicityConservation_
Option to impose helicity conservation on the real NLO ME's (greatly improves evaluation time).
Energy2 t_u_M_R_qg(realVVKinematics R) const
The matrix element q + g -> n + q times tk*uk.
bornVVKinematics B_
Born / virtual 2->2 kinematics.
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.
void setTheScales(Energy pT)
sets the QCD, EW and PDF scales
double M_Born_
Member variable to store the Born matrix element as given in Equation 3.1 - 3.3 in NPB 383 (1992) ***...
double Rtilde_Ltilde_gqb_on_x(tcPDPtr a, tcPDPtr b) const
Function for calculation of the initiated real contribution.
double Vtilde_universal(realVVKinematics S) const
Calculate the universal soft-virtual contribution to the NLO weight.
void persistentInput(PersistentIStream &is, int version)
Function used to read in object persistently.
double lo_me(bool getMatrix) const
The matrix element for the kinematical configuration previously provided by the last call to setKinem...
double fixed_alphaS_
The value of alphaS to use for the nlo weight if nlo_alphaS_opt_=1.
double eZ_
The TGC coupling.
double lo_me2_
The value of the leading order qqbar->VV matrix element.
double preqqbar_
The prefactor, for the channel.
unsigned int channel_
This specifies the emitting configuration: 1: q + qbar -> V1 + V2 + g 2: q + g -> V1 + V2 + q 3: g + ...
Energy2 t_u_M_R_qg_hel_amp(realVVKinematics R, bool getMatrix) const
The matrix element q + g -> n + q times tk*uk.
MEPP2VVPowheg()
The default constructor.
bool sanityCheck() const
This member check the collinear limits of the real emission matrix elements are equal to the appropri...
double M_Born_ZZ(bornVVKinematics B) const
The Born matrix element as given in Equation 2.18 - 2.19 in NPB 357 (1991) *** modulo a factor 1/(2s)...
bool isotropicDecayer()
Member which selects a two body decay mode for each vector boson and distributes decay products isotr...
double M_V_regular_WW(realVVKinematics S) const
M_V_regular_WW is the regular part of the one-loop WW matrix element exactly as defined in Eqs.
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...
double Ctilde_Ltilde_qq_on_x(tcPDPtr a, tcPDPtr b, realVVKinematics C) const
Function for calculation of the initiated real contribution.
double b0_
The QCD beta function divided by 4pi, (11-2/3*nf)/4/pi, with nf = 5.
tcPDPtr quark_
The ParticleData object for the quark and antiquark (which can be in a different order to ab_ and bb_...
virtual double me2() const
The matrix element for the kinematical configuration previously provided by the last call to setKinem...
bool getEvent(vector< Lorentz5Momentum > &pnew, unsigned int &emissiontype)
generates the hardest emission (yj,p)
Energy pT_
The transverse momentum of the jet.
double lo_me() const
The leading order matrix element - using helicity amplitudes.
double gW_
The weak coupling and the sin (squared) of the Weinberg angle.
AbstractFFVVertexPtr FFZvertex_
The Z fermion-antifermionvertex.
virtual void doinit()
Initialize this object after the setup phase before saving an EventGenerator to disk.
double NC_
The number of colours.
realVVKinematics S_
Soft limit of the 2->3 real emission kinematics.
realVVKinematics R_
The resolved 2->3 real emission kinematics.
double NLOweight() const
Calculate the correction weight with which leading-order configurations are re-weighted.
bool realMESpinCorrelations_
If this boolean is true the n+1 body helicity amplitudes will be used to calculate a hard vertex base...
double lo_me_
The colour & spin averaged n-body (leading order) matrix element squared.
Energy2 t_u_M_R_gqb_
Member variable to store the matrix element g + qb -> n + q times tk*uk.
double getResult(int emis_type, realVVKinematics R, Energy pT)
Returns the matrix element for a given type of process, rapidity of the jet and transverse momentum ...
double CF_
The CF_ colour factor.
void persistentOutput(PersistentOStream &os) const
Function used to write out object persistently.
tcBeamPtr qHadron_
Pointers to the BeamParticleData objects.
double power_
Constants for the sampling.
realVVKinematics SCp_
Soft-collinear limit of the 2->3 kinematics (emission in +z direction).
double guL_
The up and down, left handed, quark-boson couplings.
Energy2 t_u_M_R_qg_hel_amp(realVVKinematics R) const
The matrix element q + g -> n + q times (tk*uk)^2 - using helicity amplitudes.
Energy2 t_u_M_R_qqb(realVVKinematics R) const
The matrix element q + qb -> n + g times tk*uk.
Energy2 PDFScale_
Scale for real emission PDF:
MEPP2VVPowheg & operator=(const MEPP2VVPowheg &)=delete
The assignment operator is private and must never be called.
vector< double > prefactor_
The prefactors as a vector for easy use.
virtual IBPtr fullclone() const
Make a clone of this object, possibly modifying the cloned object to make it sane.
Energy2 EWScale_
Scale of electroweak vertices: mVV^2 the invariant mass of the diboson system.
realVVKinematics Cp_
The collinear limit of the 2->3 kinematics (emission in +z direction).
double M_V_regular_ZZ(realVVKinematics S) const
M_V_regular_ZZ is the regular part of the one-loop ZZ matrix element exactly as defined in Eqs.
double eZ2_
The TGC coupling squared.
int fermionNumberOfMother_
Identifies the space-like mother of the branching as quark (+1) or antiquark (-1):
Energy LambdaQCD_
The fundamental QCD scale in the one-loop alpha_{S} used for the crude (not the very crude) overestim...
double lo_lumi_
Values of the PDF's before radiation.
realVVKinematics Cm_
The collinear limit of the 2->3 kinematics (emission in -z direction).
Energy2 t_u_M_R_qqb_hel_amp(realVVKinematics R, bool getMatrix) const
The matrix element q + qb -> n + g times tk*uk.
The MEPP2VV class implements the production of , and in hadron-hadron collisions.
Definition: MEPP2VV.h:24
The storage of the helicity amplitude expression for the matrix element of a hard process.
The bornVVKinematics class is used to store information on the the kinematics of the real emission pr...
Definition: VVKinematics.h:18
The realVVKinematics class is used to store information on the the kinematics of the real emission pr...
Definition: VVKinematics.h:261
ShowerInteraction
Handy header file to be included in all Shower classes.
long double ReLi2(long double)
The real part of the dilog function taken from FORTRAN Herwig.
-*- C++ -*-
Ptr< BeamParticleData >::transient_const_pointer tcBeamPtr
A typedef for the BeamParticleData.
Definition: HwMEBase.h:20
constexpr double pi
std::complex< double > Complex
ThePEG::Ptr< InterfacedBase >::pointer IBPtr
ThePEG::Ptr< Particle >::pointer PPtr
ThePEG::Ptr< ParticleData >::transient_const_pointer tcPDPtr