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

Mail to administrator
ViewVC Help
Powered by ViewVC 1.0.0