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

Mail to administrator
ViewVC Help
Powered by ViewVC 1.0.0