1
c----------------------------------------------------------------------
3
c----------------------------------------------------------------------
4
c This files takes the inputs of the standard model from a Les Houches
5
c file (param_card.dat) and calculates all the couplings that end up
6
c in the Feynman rules, i.e., in the HELAS routines.
8
c With the convention of the New Web Generation MadEvent in this
9
c part no non-trivial computation is made. The SM model secondary
10
c parameters have been all calculated by the SM calculator, SMCalc
11
c which produced param_card.dat.
13
c The only exception to the above rule is for alpha_S. In the case
14
c where a pdf is used in the initial state, then the values as(MZ)
15
c set in the les houches card is superseeded.
16
c Scales and pdf's are set in the run_card.dat.
18
c This file contains the following routines:
20
c- madgraph original call to set the parameters
21
c- lh_readin is called from here.
22
c It talks to the lh_readin through the common block values.
25
c-routine to read the LH parameters in
26
c subroutine lh_readin
29
c subroutine set_it(npara,ivalue,value,name,id,
30
c subroutine case_trap(string,length)
31
c subroutine no_spaces(buff,nchars)
32
c----------------------------------------------------------------------
36
c***********************************************************************
37
c This subroutine sets up the HELAS couplings of the STANDARD MODEL.
38
c***********************************************************************
46
c common file with the couplings
50
include 'calc_values.inc'
55
double precision ee, ee2, ez, ey, sw, cw, sc2
56
double precision gwne, gwud, lambda, lam4, xt, rew, rqcd
57
double precision alphas, alfa, alfaw, mfrun
58
external alphas, alfa, alfaw, mfrun
60
c Common to lh_readin and printout
62
double precision alpha, sin2w, gfermi, alfas
63
double precision mtMS,mbMS,mcMS,mtaMS!MSbar masses
64
double precision Vud,Vus !CKM matrix elements
65
common/values/ alpha,sin2w,gfermi,alfas,
66
& mtMS,mbMS,mcMS,mtaMS,
72
parameter( ci = ( 0.0d0, 1.0d0 ) )
84
BETA(XX) = DSQRT(One-Four*XX)
85
HFF(XX,Y) = One/(Two*Four*PI)*XX*(BETA(Y))**3
86
HEAVY(XX,Y,Z)= ( One - Half*(y**2+z**2) - Half*(y**2-z**2)**2
87
& + Three*y*z*(XX**2 - One)/(XX**2 + One) )
88
& * sqrt( (One-y**2-z**2)**2 - Four * y**2 * z**2 )
90
c------------------------------------------
91
c Start calculating the couplings for HELAS
92
c------------------------------------------
98
c As a rule we first check if a pdf has been chosen in the
99
c run_card.dat (which has been already read at this stage).
100
c If there pdfs in the initial state, then the alpha_s(MZ) used
101
c is set to the corresponding value.
103
if(scale.le.1d0) scale=zmass
105
if(lpp(1).ne.0.or.lpp(2).ne.0) then
106
if(asmz .le.0.01d0) asmz =0.118d0
108
asmz=alfas !value read from the param_card.dat
111
if(nloop.eq.0) nloop=1
113
G = DSQRT(4d0*PI*ALPHAS(SCALE)) ! use setting of the param_card.dat @ NLO
117
c auxiliary local values
119
cw = sqrt( One - sin2w )
125
sc2 = sin2w*( One - sin2w )
126
v = Two*wmass*sw/ee ! the wmass is used to calculate v
127
lambda = hmass**2 / (Two * v**2)
129
c vector boson couplings
135
c gauge & higgs boson coupling constants
137
gwwh = dcmplx( ee2/sin2w*Half*v, Zero )
138
gzzh = dcmplx( ee2/sc2*Half*v, Zero )
139
ghhh = dcmplx( -hmass**2/v*Three, Zero )
140
gwwhh = dcmplx( ee2/sin2w*Half, Zero )
141
gzzhh = dcmplx( ee2/sc2*Half, Zero)
144
c fermion-fermion-vector couplings
146
gal(1) = dcmplx( ee , Zero )
147
gal(2) = dcmplx( ee , Zero )
148
gau(1) = dcmplx( -ee*Two/Three, Zero )
149
gau(2) = dcmplx( -ee*Two/Three, Zero )
150
gad(1) = dcmplx( ee/Three , Zero )
151
gad(2) = dcmplx( ee/Three , Zero )
153
gwf(1) = dcmplx( -ee/sqrt(Two*sin2w), Zero )
154
gwf(2) = dcmplx( Zero , Zero )
156
gzn(1) = dcmplx( -ez*Half , Zero )
157
gzn(2) = dcmplx( Zero , Zero )
158
gzl(1) = dcmplx( -ez*(-Half + sin2w) , Zero )
159
gzl(2) = dcmplx( -ey , Zero )
160
gzu(1) = dcmplx( -ez*( Half - sin2w*Two/Three), Zero )
161
gzu(2) = dcmplx( ey*Two/Three , Zero )
162
gzd(1) = dcmplx( -ez*(-Half + sin2w/Three) , Zero )
163
gzd(2) = dcmplx( -ey/Three , Zero )
165
c fermion-fermion-Higgs couplings (complex) hff(2)
167
c NOTE: the running mass is evaluated @ the same order
168
c nloop of alpha_s set by the PDF choice
172
ghtop(1) = dcmplx( -mtMS/v, Zero )
174
ghtop(1) = dcmplx( Zero,Zero)
179
ghbot(1) = dcmplx( -mbMS/v, Zero )
181
ghbot(1) = dcmplx( Zero, Zero )
186
ghcha(1) = dcmplx( -mcMS/v, Zero )
188
ghcha(1) = dcmplx( Zero, Zero )
192
ghtau(1) = dcmplx( -mtaMS/v, Zero )
196
c symmetric 3x3 matrix, Vud=Vcs, Vus=Vcd Vcb=Vub=0
198
c >>>>>>>>>>>>>>>***** NOTE****<<<<<<<<<<<<<<<<<<<<<<<<<
199
c these couplings matter only when interaction_CKM.dat
200
c is used to generate all the diagrams with off-diagonal
201
c couplings. The default of MadEvent is a diagonal
204
Vus=DSQRT(1d0-Vud**2)
211
c---------------------------------------------------------
212
c Set Photon Width to Zero, used by symmetry optimization
213
c---------------------------------------------------------
217
c Z boson partial widths
219
decz = zmass / ( 24.0d0 * Pi )
220
w_z_nn = decz * ( gzn(1)**2 + gzn(2)**2 )
221
w_z_ll = decz * ( gzl(1)**2 + gzl(2)**2 )
223
w_z_uu = decz * ( gzu(1)**2 + gzu(2)**2 )
224
w_z_dd = decz * ( gzd(1)**2 + gzd(2)**2 )
225
dum = dble( (gzl(2)+gzl(1))/(gzl(2)-gzl(1)) )
226
w_z_tau = w_z_ll * heavy( dum, lmass/zmass, lmass/zmass )
227
dum = dble( (gzu(2)+gzu(1))/(gzu(2)-gzu(1)) )
228
w_z_cc = w_z_uu * heavy( dum, cmass/zmass, cmass/zmass )
229
dum = dble( (gzd(2)+gzd(1))/(gzd(2)-gzd(1)) )
230
w_z_bb = w_z_dd * heavy( dum, bmass/zmass, bmass/zmass )
232
zwidth = Three*w_z_nn + Two*w_z_ll + w_z_tau
233
& + Two*w_z_dd + w_z_uu + w_z_cc + w_z_bb
236
c W boson partial widths
238
decw = wmass / ( 24.0d0 * Pi )
239
w_w_nl = decw * ( gwf(1)**2 + gwf(2)**2 )
240
dum = dble( (gwf(2)+gwf(1))/(gwf(2)-gwf(1)) )
241
w_w_tau = w_w_nl * heavy( dum, lmass/wmass, Zero )
242
w_w_ud = w_w_nl * Three
243
w_w_cs = w_w_ud * heavy( dum, cmass/wmass, Zero )
245
wwidth = Two*w_w_nl + w_w_tau + w_w_ud + w_w_cs
250
call topwid(tmass,wmass,bmass,wwidth,gw,twidth)
255
lwidth = 2.36d-12 !tau width, PDG value
257
c LO withds of the Higgs into tt~,bb~,tau tau~.
259
if(hmass.gt.2d0*tmass) then
260
w_h_tt =3d0*cdabs(ghtop(1)**2)*hff(hmass,(tmass/hmass)**2)
265
w_h_bb =3d0*cdabs(ghbot(1)**2)*hff(hmass,(bmass/hmass)**2)
266
w_h_tau=cdabs(ghtau(1)**2)*hff(hmass,(lmass/hmass)**2)
267
w_h_bb =3d0*cdabs(ghbot(1)**2)*hff(hmass,(bmass/hmass)**2)
268
w_h_cc =3d0*cdabs(ghcha(1)**2)*hff(hmass,(cmass/hmass)**2)
269
w_h_tau=cdabs(ghtau(1)**2)*hff(hmass,(lmass/hmass)**2)
272
if(hmass.gt.2d0*wmass) then
273
w_h_ww=gfermi/8d0/pi/rt2*hmass**3*
274
& dsqrt(1-4d0*(wmass/hmass)**2)*
275
& (12d0*(wmass/hmass)**4 -4d0*(wmass/hmass)**2+1d0)
276
w_h_WLWL=gw**2/64d0/pi*hmass**3/wmass**2* !longitudinal W's
277
& dsqrt(1-4d0*(wmass/hmass)**2)*
278
& (1-2d0*(wmass/hmass)**2)**2
284
if(hmass.gt.2d0*zmass) then
285
w_h_zz=gfermi/16d0/pi/rt2*hmass**3*
286
& dsqrt(1-4d0*(zmass/hmass)**2)*
287
& (12d0*(zmass/hmass)**4 -4d0*(zmass/hmass)**2+1d0)
288
w_h_ZLZL=gw**2/128d0/pi*hmass**3/wmass**2* !longitudinal Z's
289
& dsqrt(1-4d0*(zmass/hmass)**2)*
290
& (1-2d0*(zmass/hmass)**2)**2
296
c--------------------------
297
c start interface to HDECAY
298
c--------------------------
299
c do not change things here unless you exactly know what you are doing
310
c-- do not set masses to zero here
311
ams = 0.190d0 !strange pole mass
313
if(cmass.gt.0d0) then
316
amc = 1.42d0 !charm pole mass
319
if(bmass.gt.0d0) then
322
amb = 4.7d0 !bottom pole mass
325
if(tmass.gt.0d0) then
328
amt = 174.3d0 !top pole mass
331
if(lmass.gt.0d0) then
334
almass = 1.777d0 !tau mass
337
ammuon = 0.105658389d0 !muon mass
339
alph = 137.0359895D0 !alpha_EM
342
if(wwidth.gt.0d0) then
348
if(zwidth.gt.0d0) then
361
c this calculates the branching ratios of the Higgs
362
c at the best of our knowledge
367
c------------------------
368
c end interface to HDECAY
369
c------------------------
372
c----------------------------
373
c end subroutine coupsm
374
c----------------------------
382
c----------------------------------------------------------------------
383
c Read the parameters from the lh file
385
c 1. Input values for the EW sector
386
c 2. Higgs mass and width
387
c 3. Fermion masses (pole and MSbar) and widths
388
c----------------------------------------------------------------------
394
parameter (maxpara=1000)
401
real*8 value(maxpara)
402
integer ivalue(maxpara),n
403
character*20 name(maxpara),bn,dumstring
404
logical block_found,done,fopened
405
integer i,name_length,idum
410
character*20 block_name
417
c Common to lh_readin and printout
419
double precision alpha, sin2w, gfermi, alfas
420
double precision mtMS,mbMS,mcMS,mtaMS!MSbar masses
421
double precision Vud,Vus !CKM matrix elements
422
common/values/ alpha,sin2w,gfermi,alfas,
423
& mtMS,mbMS,mcMS,mtaMS,
437
c looks for the blocks or for decay
439
do while(.not.block_found)
440
read(lunp,'(a132)',end=99,err=99) buff
441
c-- change to lower case
442
call case_trap(buff,20)
443
if(buff(1:5).eq.'block') then
444
if(index(buff,"#").ne.0) l1=index(buff,"#")-1
445
block_name=buff(6:min(l1,26))
446
call no_spaces(block_name,name_length)
447
c write(*,*) block_name(1:name_length)
449
elseif(buff(1:5).eq.'decay') then
452
if(index(buff,"#").ne.0) l1=index(buff,"#")-1 ! ignore comments
453
read(buff(6:l1),*) ivalue(n),value(n)
462
read(lunp,'(a132)',end=99,err=99) buff
463
call case_trap(buff,20)
464
if(buff(1:1).eq.'b'.or.buff(1:1).eq.'d') then
468
if(buff(1:1).ne.'#') then !if it not a comment
471
if(index(buff,"#").ne.0) l1=index(buff,"#")-1 ! ignore comments
473
c WARNING:... not all blocks have the same sintax!! You need to change it
474
c depending on the block you are reading
477
if(block_name(1:5).eq."mgckm") then
478
read(buff(1:l1),*) ivalue(n),idum,value(n)
479
elseif (block_name(1:6).eq."spinfo".or.
480
&block_name(1:6).eq."dcinfo") then
481
read(buff(1:l1),*) ivalue(n),dumstring
483
read(buff(1:l1),*) ivalue(n),value(n)
485
name(n)=block_name(1:name_length)
486
c write(*,"(1x,i2,2x,e16.8,1x,a)")
487
c & ivalue(n),value(n),name(n)
489
end do ! do while in the block
493
end do ! do while the entire file
499
call set_it(n,ivalue,value,name,1,bn,alpha,128.9d0)
501
call set_it(n,ivalue,value,name,2,bn,gfermi,0.1166d-4)
502
call set_it(n,ivalue,value,name,3,bn,alfas,0.119d0)
503
call set_it(n,ivalue,value,name,4,bn,zmass,91.188d0)
504
call set_it(n,ivalue,value,name,6,bn,tmass,174.3d0)
505
call set_it(n,ivalue,value,name,7,bn,lmass,1.777d0)
507
call set_it(n,ivalue,value,name,1,bn,sin2w,0.2312d0)
508
call set_it(n,ivalue,value,name,2,bn,wmass,80.419d0)
510
call set_it(n,ivalue,value,name,4,bn,mcMS,1.25d0)
511
call set_it(n,ivalue,value,name,5,bn,mbMS,4.2d0)
512
call set_it(n,ivalue,value,name,6,bn,mtMS,174d0)
513
call set_it(n,ivalue,value,name,15,bn,mtaMS,1.777d0)
515
call set_it(n,ivalue,value,name,1,bn,vud,1d0)
517
call set_it(n,ivalue,value,name,4,bn,cmass,1.4d0)
518
call set_it(n,ivalue,value,name,5,bn,bmass,4.7d0)
519
call set_it(n,ivalue,value,name,6,bn,tmass,tmass*1d0)
520
call set_it(n,ivalue,value,name,15,bn,lmass,lmass*1d0)
521
call set_it(n,ivalue,value,name,25,bn,hmass,120d0)
522
call set_it(n,ivalue,value,name,23,bn,zmass,zmass*1d0)
523
call set_it(n,ivalue,value,name,24,bn,wmass,wmass*1d0)
526
c call set_it(n,ivalue,value,name,6,bn,twidth,1.5083d0)
527
c call set_it(n,ivalue,value,name,25,bn,hwidth,0.0037d0)
528
c call set_it(n,ivalue,value,name,23,bn,zwidth,2.441d0)
529
c call set_it(n,ivalue,value,name,24,bn,wwidth,2.0476d0)
532
write(*,*) ' >>> Widths in param_card.dat are ignored <<<'
539
subroutine set_it(npara,ivalue,value,name,id,
540
& block_name,var,def_value)
541
c----------------------------------------------------------------------------------
542
c finds the parameter value in block_name and associate var to it.
543
c If it is not found a default is given.
544
c----------------------------------------------------------------------------------
551
parameter (maxpara=100)
555
integer npara,ivalue(maxpara),id
556
character*20 block_name,name(maxpara)
557
real*8 var,def_value,value(maxpara)
568
found = (id.eq.ivalue(i)).and.(name(i).eq.block_name)
576
c write (*,*) "Warning: parameter not found"
577
c write (*,*) " setting it to default value ",def_value
585
subroutine case_trap(string,length)
586
c**********************************************************
587
c change string to lowercase if the input is not
588
c**********************************************************
602
if(k.ge.65.and.k.le.90) then !upper case A-Z
603
k=ichar(string(i:i))+32
612
subroutine no_spaces(buff,nchars)
613
c**********************************************************************
614
c Given buff a buffer of words separated by spaces
615
c returns it where all space are moved to the right
616
c returns also the length of the single word.
617
c maxlength is the length of the buffer
618
c**********************************************************************
624
parameter (maxline=20)
630
character*(maxline) buff
631
integer nchars,maxlength
636
character*(maxline) temp
641
c write (*,*) "buff=",buff(1:maxlength)
643
if(buff(i:i).ne.null) then
645
temp(nchars:nchars)=buff(i:i)
647
c write(*,*) i,":",buff(1:maxlength),":",temp(1:nchars),":"
653
SUBROUTINE TOPWID(RMT,RMW,RMB,RGW,GW,RGT)
654
c*************************************************************************
655
c THE TOTAL WEAK DECAY WIDTH OF THE TOP QUARK, INCLUDING
656
c THE EFFECTS OF BOTTOM MASS AND, IF IGW=1, A FINITE W WIDTH.
657
c From James Stirling 6-10-94
667
c*************************************************************************
668
IMPLICIT COMPLEX*16(A-H,O-Z)
669
REAL*8 RMT,RMB,RMW,XW,XB,RGW,RGT,GW
672
XGW=dcmplx(GW/2d0/dsqrt(2d0))
678
OM=1.+RM-DCMPLX(RMW**2,RMW*RGW)/RMT**2
679
Y1=OM+CDSQRT(OM*OM-4.*RM)
680
Y0=OM-CDSQRT(OM*OM-4.*RM)
684
D0=(-Y0**8+3.*Y0**7*RM+3.*Y0**7-8.*Y0**6*RM-12.*Y0**5*RM**
685
. 2-12.*Y0**5*RM+96.*Y0**4*RM**2-48.*Y0**3*RM**3-48.*Y0**3*
686
. RM**2-128.*Y0**2*RM**3+192.*Y0*RM**4+192.*Y0*RM**3-256.*
687
. RM**4)/(24.*Y0**4*(Y1-Y0))
688
D1=(-Y1**8+3.*Y1**7*RM+3.*Y1**7-8.*Y1**6*RM-12.*Y1**5*RM**
689
. 2-12.*Y1**5*RM+96.*Y1**4*RM**2-48.*Y1**3*RM**3-48.*Y1**3*
690
. RM**2-128.*Y1**2*RM**3+192.*Y1*RM**4+192.*Y1*RM**3-256.*
691
. RM**4)/(24.*Y1**4*(Y1-Y0))
692
A4=(32.*RM**4*(Y1-Y0))/(3.*Y1*Y0*(Y1-Y0))
693
A3=(8.*RM**3*(-3.*Y1**2*Y0*RM-3.*Y1**2*Y0+4.*Y1**2*RM+3.*
694
. Y1*Y0**2*RM+3.*Y1*Y0**2-4.*Y0**2*RM))/(3.*Y1**2*Y0**2*(Y1
696
A2=(8.*RM**3*(2.*Y1**3*Y0**2-3.*Y1**3*Y0*RM-3.*Y1**3*Y0+4.
697
. *Y1**3*RM-2.*Y1**2*Y0**3+3.*Y1*Y0**3*RM+3.*Y1*Y0**3-4.*Y0
698
. **3*RM))/(3.*Y1**3*Y0**3*(Y1-Y0))
699
A1=(2.*RM**2*(3.*Y1**4*Y0**3*RM+3.*Y1**4*Y0**3+8.*Y1**4*Y0
700
. **2*RM-12.*Y1**4*Y0*RM**2-12.*Y1**4*Y0*RM+16.*Y1**4*RM**2
701
. -3.*Y1**3*Y0**4*RM-3.*Y1**3*Y0**4-8.*Y1**2*Y0**4*RM+12.*
702
. Y1*Y0**4*RM**2+12.*Y1*Y0**4*RM-16.*Y0**4*RM**2))/(3.*Y1**
704
B0=(Y1**3-3.*Y1**2*RM-3.*Y1**2+8.*Y1*RM-Y0**3+3.*Y0**2*RM+
705
. 3.*Y0**2-8.*Y0*RM)/(24.*(Y1-Y0))
706
B1=(Y1+Y0-3.*RM-3.)/24.
709
RINT=D0*CDLOG((Z1-Y0)/(Z0-Y0))
710
. -D1*CDLOG((Y1-Z1)/(Y1-Z0))
711
. -A4/3.*(1./Z1**3-1./Z0**3)
712
. -A3/2.*(1./Z1**2-1./Z0**2)
713
. -A2 *(1./Z1 -1./Z0 )
716
. +B1/2.*(Z1**2-Z0**2)
717
. +B2/3.*(Z1**3-Z0**3)
721
* TOTAL WIDTH INCLUDES FLAVOUR & COLOUR FACTORS
722
RGT=RMT**3/(RMW*RGW)*XGW4/(8.*PI**3)*DIMAG(RINT)