[cvs] / MG_ME / Models / sm / couplings.f Repository:
ViewVC logotype

Annotation of /MG_ME/Models/sm/couplings.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.13 - (view) (download)

1 : madgraph 1.1 c----------------------------------------------------------------------
2 :     C couplings.f
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.
7 :     c
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.
12 :     c
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.
17 :     c
18 :     c This file contains the following routines:
19 :     c
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.
23 :     c subroutine setpara
24 :     c
25 :     c-routine to read the LH parameters in
26 :     c subroutine lh_readin
27 :     c
28 :     c-to set
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----------------------------------------------------------------------
33 :    
34 :    
35 : madgraph 1.8 subroutine setpara(param_name,readlha)
36 : madgraph 1.1 c***********************************************************************
37 :     c This subroutine sets up the HELAS couplings of the STANDARD MODEL.
38 :     c***********************************************************************
39 :     implicit none
40 :     c
41 :     c local
42 :     c
43 : madgraph 1.7 character*(*) param_name
44 : madgraph 1.8 logical readlha
45 : madgraph 1.1 integer i
46 :     real*8 dum
47 :     c
48 :     c common file with the couplings
49 :     c
50 :     include 'coupl.inc'
51 :     c
52 :     c local
53 :     c
54 :     double precision v
55 : madgraph 1.10 double precision ee, ee2, ez, ey, sw, cw, sc2, sin2w, wm
56 : madgraph 1.1 double precision gwne, gwud, lambda, lam4, xt, rew, rqcd
57 :     double precision alphas, alfa, alfaw, mfrun
58 :     external alphas, alfa, alfaw, mfrun
59 :     c
60 :     c Common to lh_readin and printout
61 :     c
62 : madgraph 1.10 double precision alpha, gfermi, alfas
63 : madgraph 1.1 double precision mtMS,mbMS,mcMS,mtaMS!MSbar masses
64 :     double precision Vud,Vus !CKM matrix elements
65 : madgraph 1.10 common/values/ alpha,gfermi,alfas,
66 : madgraph 1.1 & mtMS,mbMS,mcMS,mtaMS,
67 :     & Vud
68 :     c
69 :     c constants
70 :     c
71 :     double complex ci
72 :     parameter( ci = ( 0.0d0, 1.0d0 ) )
73 :     double precision Zero, One, Two, Three, Four, Half, Rt2
74 :     parameter( Zero = 0.0d0, One = 1.0d0, Two = 2.0d0 )
75 :     parameter( Three = 3.0d0, Four = 4.0d0, Half = 0.5d0 )
76 :     parameter( Rt2 = 1.414213562d0 )
77 :     double precision Pi, Fourpi
78 :     parameter( Pi = 3.14159265358979323846d0 )
79 :     parameter( Fourpi = Four * Pi )
80 :     c
81 :     c------------------------------------------
82 :     c Start calculating the couplings for HELAS
83 :     c------------------------------------------
84 :     c
85 : madgraph 1.8 if(readlha) then
86 :     call lh_readin(param_name)
87 :     G = DSQRT(4d0*PI*ALFAS) ! use setting of the param_card.dat @ NLO
88 : madgraph 1.1 c
89 :     c auxiliary local values
90 :     c
91 : madgraph 1.10 wm = sqrt(zmass**2/Two+
92 :     $ sqrt(zmass**4/Four-Pi/Rt2*alpha/gfermi*zmass**2))
93 :     sin2w = One-(wm/zmass)**2
94 : madgraph 1.1 cw = sqrt( One - sin2w )
95 :     ee2 = alpha * Fourpi
96 :     sw = sqrt( sin2w )
97 :     ee = sqrt( ee2 )
98 :     ez = ee/(sw*cw)
99 :     ey = ee*(sw/cw)
100 :     sc2 = sin2w*( One - sin2w )
101 : madgraph 1.10 v = Two*wm*sw/ee ! the wmass is used to calculate v
102 : madgraph 1.1 lambda = hmass**2 / (Two * v**2)
103 :     c
104 :     c vector boson couplings
105 :     c
106 :     gw = ee/sw
107 :     gwwa = ee
108 :     gwwz = ee*cw/sw
109 :     c
110 :     c gauge & higgs boson coupling constants
111 :     c
112 :     gwwh = dcmplx( ee2/sin2w*Half*v, Zero )
113 :     gzzh = dcmplx( ee2/sc2*Half*v, Zero )
114 :     ghhh = dcmplx( -hmass**2/v*Three, Zero )
115 :     gwwhh = dcmplx( ee2/sin2w*Half, Zero )
116 :     gzzhh = dcmplx( ee2/sc2*Half, Zero)
117 :     ghhhh = ghhh/v
118 :     c
119 :     c fermion-fermion-vector couplings
120 :     c
121 :     gal(1) = dcmplx( ee , Zero )
122 :     gal(2) = dcmplx( ee , Zero )
123 :     gau(1) = dcmplx( -ee*Two/Three, Zero )
124 :     gau(2) = dcmplx( -ee*Two/Three, Zero )
125 :     gad(1) = dcmplx( ee/Three , Zero )
126 :     gad(2) = dcmplx( ee/Three , Zero )
127 :    
128 :     gwf(1) = dcmplx( -ee/sqrt(Two*sin2w), Zero )
129 :     gwf(2) = dcmplx( Zero , Zero )
130 :    
131 :     gzn(1) = dcmplx( -ez*Half , Zero )
132 :     gzn(2) = dcmplx( Zero , Zero )
133 :     gzl(1) = dcmplx( -ez*(-Half + sin2w) , Zero )
134 :     gzl(2) = dcmplx( -ey , Zero )
135 :     gzu(1) = dcmplx( -ez*( Half - sin2w*Two/Three), Zero )
136 :     gzu(2) = dcmplx( ey*Two/Three , Zero )
137 :     gzd(1) = dcmplx( -ez*(-Half + sin2w/Three) , Zero )
138 :     gzd(2) = dcmplx( -ey/Three , Zero )
139 :     c
140 :     c fermion-fermion-Higgs couplings (complex) hff(2)
141 :     c
142 :     c NOTE: the running mass is evaluated @ the same order
143 :     c nloop of alpha_s set by the PDF choice
144 :     c
145 :    
146 :     if(mtMS.gt.1d0) then
147 :     ghtop(1) = dcmplx( -mtMS/v, Zero )
148 :     else
149 :     ghtop(1) = dcmplx( Zero,Zero)
150 :     endif
151 :     ghtop(2) = ghtop(1)
152 :    
153 :     if(mbMS.gt.1d0) then
154 :     ghbot(1) = dcmplx( -mbMS/v, Zero )
155 :     else
156 :     ghbot(1) = dcmplx( Zero, Zero )
157 :     endif
158 :     ghbot(2) = ghbot(1)
159 :    
160 :     if(mcMS.gt.1d0) then
161 : madgraph 1.6 ghcha(1) = dcmplx( -mcMS/v, Zero )
162 : madgraph 1.1 else
163 :     ghcha(1) = dcmplx( Zero, Zero )
164 :     endif
165 :     ghcha(2) = ghcha(1)
166 :    
167 :     ghtau(1) = dcmplx( -mtaMS/v, Zero )
168 :     ghtau(2) = ghtau(1)
169 :     c
170 :     c CKM matrix:
171 :     c symmetric 3x3 matrix, Vud=Vcs, Vus=Vcd Vcb=Vub=0
172 :     c
173 :     c >>>>>>>>>>>>>>>***** NOTE****<<<<<<<<<<<<<<<<<<<<<<<<<
174 :     c these couplings matter only when interaction_CKM.dat
175 :     c is used to generate all the diagrams with off-diagonal
176 :     c couplings. The default of MadEvent is a diagonal
177 :     c CKM matrix.
178 :    
179 :     Vus=DSQRT(1d0-Vud**2)
180 :     do i=1,2
181 :     gwfc(i) = gwf(i)*Vud
182 :     gwfs(i) = gwf(i)*Vus
183 :     gwfm(i) =-gwf(i)*Vus
184 :     enddo
185 :    
186 :     c---------------------------------------------------------
187 :     c Set Photon Width to Zero, used by symmetry optimization
188 :     c---------------------------------------------------------
189 :    
190 :     awidth = 0d0
191 :    
192 : madgraph 1.12 endif
193 :     c
194 :     c Strong coupling
195 :     c
196 :     c As a rule we first check if a pdf has been chosen in the
197 :     c run_card.dat (which has been already read at this stage).
198 :     c If there pdfs in the initial state, then the alpha_s(MZ) used
199 :     c is set to the corresponding value.
200 :    
201 :     GG(1) = -G
202 :     GG(2) = -G
203 :    
204 : madgraph 1.1 c----------------------------
205 :     c end subroutine coupsm
206 :     c----------------------------
207 :    
208 :    
209 :     return
210 :     end
211 :    
212 : madgraph 1.7 subroutine lh_readin(param_name)
213 : madgraph 1.1 c----------------------------------------------------------------------
214 :     c Read the parameters from the lh file
215 :     c
216 :     c 1. Input values for the EW sector
217 :     c 2. Higgs mass and width
218 :     c 3. Fermion masses (pole and MSbar) and widths
219 :     c----------------------------------------------------------------------
220 :     implicit none
221 : madgraph 1.7
222 : madgraph 1.1 c
223 :     c parameters
224 :     c
225 :     integer maxpara
226 :     parameter (maxpara=100)
227 : madgraph 1.10 double precision Two, Four, Rt2, Pi
228 :     parameter( Two = 2.0d0, Four = 4.0d0 )
229 :     parameter( Rt2 = 1.414213562d0 )
230 :     parameter( Pi = 3.14159265358979323846d0 )
231 : madgraph 1.1 c
232 :     c local
233 :     c
234 : madgraph 1.7 character*(*) param_name
235 : madgraph 1.1 integer npara,l1,l2,id
236 :     character*132 buff
237 :     real*8 real_value
238 :     real*8 value(maxpara)
239 :     integer ivalue(maxpara),n
240 :     character*20 name(maxpara),bn
241 :     logical block_found,done,fopened
242 :     integer iunit,i,name_length,idum
243 :     c
244 :     c block info
245 :     c
246 :     character*20 block_name
247 :     c
248 :     c Common
249 :     c
250 :     include 'coupl.inc'
251 :     c
252 :     c Common to lh_readin and printout
253 :     c
254 : madgraph 1.10 double precision alpha, gfermi, alfas
255 : madgraph 1.1 double precision mtMS,mbMS,mcMS,mtaMS!MSbar masses
256 :     double precision Vud,Vus !CKM matrix elements
257 : madgraph 1.10 common/values/ alpha,gfermi,alfas,
258 : madgraph 1.1 & mtMS,mbMS,mcMS,mtaMS,
259 :     & Vud
260 :     c
261 :     c----------
262 :     c start
263 :     c----------
264 :     c
265 :     c open file
266 :     c
267 :     iunit=14
268 : madgraph 1.9 call open_file_mdl(iunit,param_name,fopened)
269 : madgraph 1.1 done=.false.
270 :    
271 :     n=0
272 :     do while(.not.done)
273 :     block_found=.false.
274 :     c
275 :     c looks for the blocks or for decay
276 :     c
277 :     do while(.not.block_found)
278 :     read(iunit,'(a132)',end=99,err=99) buff
279 :     c-- change to lower case
280 :     call case_trap(buff,20)
281 :     if(buff(1:5).eq.'block') then
282 :     if(index(buff,"#").ne.0) l1=index(buff,"#")-1
283 :     block_name=buff(6:min(l1,26))
284 :     call no_spaces(block_name,name_length)
285 :     c write(*,*) block_name(1:name_length)
286 :     block_found=.true.
287 :     elseif(buff(1:5).eq.'decay') then
288 :     n=n+1
289 :     l1=30
290 :     if(index(buff,"#").ne.0) l1=index(buff,"#")-1 ! ignore comments
291 :     read(buff(6:l1),*) ivalue(n),value(n)
292 :     name(n)="decay"
293 :     endif
294 :     end do
295 :     c
296 :     c
297 :     if(block_found) then
298 :     do while(.true.)
299 :     read(iunit,'(a132)',end=99,err=99) buff
300 :     call case_trap(buff,20)
301 :     if(buff(1:1).eq.'b'.or.buff(1:1).eq.'d') then
302 :     backspace iunit
303 :     exit
304 :     endif
305 :     if(buff(1:1).ne.'#') then !if it not a comment
306 :     n=n+1
307 :     l1=30
308 :     if(index(buff,"#").ne.0) l1=index(buff,"#")-1 ! ignore comments
309 :     c
310 :     c WARNING:... not all blocks have the same sintax!! You need to change it
311 :     c depending on the block you are reading
312 :     c
313 :     if(block_name(1:5).eq."mgckm") then
314 :     read(buff(1:l1),*) ivalue(n),idum,value(n)
315 :     else
316 :     read(buff(1:l1),*) ivalue(n),value(n)
317 :     endif
318 :     name(n)=block_name(1:name_length)
319 :     c write(*,"(1x,i2,2x,e16.8,1x,a)")
320 :     c & ivalue(n),value(n),name(n)
321 :     endif
322 :     end do ! do while in the block
323 :     else
324 :     done=.true.
325 :     endif
326 :     end do ! do while the entire file
327 :    
328 :    
329 :     99 continue
330 :    
331 :     bn="sminputs"
332 :     call set_it(n,ivalue,value,name,1,bn,alpha,128.9d0)
333 :     alpha=1d0/alpha
334 :     call set_it(n,ivalue,value,name,2,bn,gfermi,0.1166d-4)
335 :     call set_it(n,ivalue,value,name,3,bn,alfas,0.119d0)
336 :     call set_it(n,ivalue,value,name,4,bn,zmass,91.188d0)
337 :     call set_it(n,ivalue,value,name,6,bn,tmass,174.3d0)
338 :     call set_it(n,ivalue,value,name,7,bn,lmass,1.777d0)
339 :     bn="mgyukawa"
340 :     call set_it(n,ivalue,value,name,4,bn,mcMS,1.25d0)
341 :     call set_it(n,ivalue,value,name,5,bn,mbMS,4.2d0)
342 :     call set_it(n,ivalue,value,name,6,bn,mtMS,174d0)
343 :     call set_it(n,ivalue,value,name,15,bn,mtaMS,1.777d0)
344 :     bn="mgckm"
345 :     call set_it(n,ivalue,value,name,1,bn,vud,1d0)
346 :     bn="mass"
347 :     call set_it(n,ivalue,value,name,4,bn,cmass,1.4d0)
348 :     call set_it(n,ivalue,value,name,5,bn,bmass,4.7d0)
349 :     call set_it(n,ivalue,value,name,6,bn,tmass,tmass*1d0)
350 :     call set_it(n,ivalue,value,name,15,bn,lmass,lmass*1d0)
351 :     call set_it(n,ivalue,value,name,25,bn,hmass,120d0)
352 :     call set_it(n,ivalue,value,name,23,bn,zmass,zmass*1d0)
353 : madgraph 1.10 call set_it(n,ivalue,value,name,24,bn,wmass,sqrt(zmass**2/Two+
354 :     $ sqrt(zmass**4/Four-Pi/Rt2*alpha/gfermi*zmass**2)))
355 : madgraph 1.1 bn="decay"
356 :     call set_it(n,ivalue,value,name,6,bn,twidth,1.5083d0)
357 :     call set_it(n,ivalue,value,name,25,bn,hwidth,0.0037d0)
358 :     call set_it(n,ivalue,value,name,23,bn,zwidth,2.441d0)
359 :     call set_it(n,ivalue,value,name,24,bn,wwidth,2.0476d0)
360 :    
361 :     return
362 :     end
363 :    
364 :    
365 :     subroutine set_it(npara,ivalue,value,name,id,
366 :     & block_name,var,def_value)
367 :     c----------------------------------------------------------------------------------
368 :     c finds the parameter value in block_name and associate var to it.
369 :     c If it is not found a default is given.
370 :     c----------------------------------------------------------------------------------
371 :     implicit none
372 :    
373 :     c
374 :     c parameters
375 :     c
376 :     integer maxpara
377 :     parameter (maxpara=100)
378 :     c
379 :     c arguments
380 :     c
381 :     integer npara,ivalue(maxpara),id
382 :     character*20 block_name,name(maxpara)
383 :     real*8 var,def_value,value(maxpara)
384 :     c
385 :     c local
386 :     c
387 :     logical found
388 :     integer i
389 :     c
390 :     c start
391 :     c
392 :     found=.false.
393 :     do i=1,npara
394 :     found = (id.eq.ivalue(i)).and.(name(i).eq.block_name)
395 :     if(found) then
396 :     var=value(i)
397 :     exit
398 :     endif
399 :     enddo
400 :    
401 :     if (.not.found) then
402 :     c write (*,*) "Warning: parameter not found"
403 :     c write (*,*) " setting it to default value ",def_value
404 :     var=def_value
405 :     endif
406 :     return
407 :    
408 :     end
409 :    
410 :    
411 :     subroutine case_trap(string,length)
412 :     c**********************************************************
413 :     c change string to lowercase if the input is not
414 :     c**********************************************************
415 :     implicit none
416 :     c
417 :     c ARGUMENT
418 :     c
419 :     character*(*) string
420 :     integer length
421 :     c
422 :     c LOCAL
423 :     c
424 :     integer i,k
425 :    
426 :     do i=1,length
427 :     k=ichar(string(i:i))
428 :     if(k.ge.65.and.k.le.90) then !upper case A-Z
429 :     k=ichar(string(i:i))+32
430 :     string(i:i)=char(k)
431 :     endif
432 :     enddo
433 :    
434 :     return
435 :     end
436 :    
437 :    
438 :     subroutine no_spaces(buff,nchars)
439 :     c**********************************************************************
440 :     c Given buff a buffer of words separated by spaces
441 :     c returns it where all space are moved to the right
442 :     c returns also the length of the single word.
443 :     c maxlength is the length of the buffer
444 :     c**********************************************************************
445 :     implicit none
446 :     c
447 :     c Constants
448 :     c
449 :     integer maxline
450 :     parameter (maxline=20)
451 :     character*1 null
452 :     parameter (null=' ')
453 :     c
454 :     c Arguments
455 :     c
456 :     character*(maxline) buff
457 :     integer nchars,maxlength
458 :     c
459 :     c Local
460 :     c
461 :     integer i,j
462 :     character*(maxline) temp
463 :     c-----
464 :     c Begin Code
465 :     c-----
466 :     nchars=0
467 :     c write (*,*) "buff=",buff(1:maxlength)
468 :     do i=1,maxline
469 :     if(buff(i:i).ne.null) then
470 :     nchars=nchars+1
471 :     temp(nchars:nchars)=buff(i:i)
472 :     endif
473 :     c write(*,*) i,":",buff(1:maxlength),":",temp(1:nchars),":"
474 :     enddo
475 :     buff=temp
476 :     end
477 :    
478 : madgraph 1.13 subroutine open_file_mdl(lun,filename,fopened)
479 : madgraph 1.4 c***********************************************************************
480 :     c opens file input-card.dat in current directory or above
481 :     c***********************************************************************
482 :     implicit none
483 :     c
484 :     c Arguments
485 :     c
486 :     integer lun
487 :     logical fopened
488 :     character*(*) filename
489 : madgraph 1.13 character*90 tempname
490 : madgraph 1.4 integer fine
491 : madgraph 1.13 integer dirup,i
492 : madgraph 1.4
493 :     c-----
494 :     c Begin Code
495 :     c-----
496 :     c
497 :     c first check that we will end in the main directory
498 :     c
499 :     tempname=filename
500 :     fine=index(tempname,' ')
501 :     if(fine.eq.0) fine=len(tempname)
502 : madgraph 1.13 tempname=tempname(1:fine)
503 : madgraph 1.4 c
504 :     c if I have to read a card
505 :     c
506 : madgraph 1.13 if(index(filename,"_card").gt.0) then
507 :     tempname='./Cards/'//tempname
508 : madgraph 1.4 endif
509 :    
510 :    
511 : madgraph 1.13 fopened=.false.
512 :     do i=0,5
513 :     open(unit=lun,file=tempname,status='old',ERR=30)
514 :     fopened=.true.
515 :     write(*,*) 'read model file',tempname
516 :     exit
517 :     30 tempname='../'//tempname
518 :     if (i.eq.5)then
519 :     print 'Warning: file ',tempname,' is not correct'
520 :     stop
521 :     endif
522 :     enddo
523 :    
524 : madgraph 1.4
525 :     return
526 :     end
527 :    
528 : madgraph 1.13
529 :    

Mail to administrator
ViewVC Help
Powered by ViewVC 1.0.0