C++InterfacetoTauola
initwksw.f
1 
2  SUBROUTINE initwkswdelt(mode,IDEX,IDFX,SVAR,SWSQEFF, DELTSQ, DeltV, GMU, ALPHAINV, AMZi, GAMMZi, KEYGSW,
3  &regsw1,cimgsw1,regsw2,cimgsw2,regsw3,cimgsw3,regsw4,cimgsw4,regsw6,cimgsw6 )
4 
5 
6 ! initialization routine coupling masses etc., explicitly varying SWSQ
7  IMPLICIT REAL*8 (a-h,o-z)
8  COMMON / t_beampm / ene ,amin,amfin,ide,idf
9  REAL*8 ene ,amin,amfin
10  COMMON / t_gauspm /ss,poln,t3e,qe,t3f,qf
11  & ,xupgi0 ,xupzi0 ,xupgf0 ,xupzf0
12  & ,ndiag0,ndiaga,keya,keyz
13  & ,itce,jtce,itcf,jtcf,kolor
14  REAL*8 ss,poln,t3e,qe,t3f,qf
15  & ,xupgi0(2),xupzi0(2),xupgf0(2),xupzf0(2)
16  COMMON / t_gauspm1/vvcor, zetvpi, gamvpi
17  & ,xupgi ,xupzi ,xupgf ,xupzf
18  COMPLEX*16 vvcor, zetvpi, gamvpi
19  COMPLEX*16 xupgi(2),xupzi(2),xupgf(2),xupzf(2)
20 
21  COMMON / t_gswprmn /swsq,amw,amz,amh,amtop,gammz
22  REAL*8 swsq,amw,amz,amh,amtop,gammz
23  COMMON / t_ewn / gmun, alphainvn
24  REAL*8 gmun, alphainvn
25  COMPLEX *16 gsw(10)
26  REAL*8 pi
27  DATA pi /3.141592653589793238462643d0/
28  gsw(1) = dcmplx(regsw1,cimgsw1)
29  gsw(2) = dcmplx(regsw2,cimgsw2)
30  gsw(3) = dcmplx(regsw3,cimgsw3)
31  gsw(4) = dcmplx(regsw4,cimgsw4)
32  ! GSW(5) out
33  gsw(6) = dcmplx(regsw6,cimgsw6)
34 
35 C PRINT *, ' initwksw GSW = ', SWSQEFF, ReGSW1, CImGSW1, ReGSW2, CImGSW2, ReGSW6, CImGSW6
36 
37 C SWSQ = sin2 (theta Weinberg)
38 C AMW,AMZ = W & Z boson masses respectively
39 C AMH = the Higgs mass
40 C AMTOP = the top mass
41 C GAMMZ = Z0 width
42 C
43  ene=sqrt(svar)/2
44  amin=0.511d-3
45  swsq=swsqeff
46  amz=amzi !91.1887
47  gammz=gammzi !2.4952
48  gmun=gmu
49  alphainvn=alphainv
50 
51 
52 C Gfermi=1.16639d-5
53  gfermi=gmu
54 
55  zetvpi = gfermi *amz**2 *alphainv /(dsqrt(2.d0)*8.d0*pi)
56  $ *(swsq*(1d0-swsq)) *16d0
57  $ * gsw(1)
58 C updated following KK2f_defaults
59 C IF( KEYGSW.NE.0) THEN
60 C GAMMZ=2.50072032
61 C ENDIF
62 
63 
64 
65 
66  gamvpi = 1d0 /(2d0-gsw(6))
67 
68 C PRINT *, ' initwksw ZetVPi, GamVPi = ', GSW(1), ZetVPi, GamVPi
69 
70 
71  IF (idfx.EQ. 11) then
72  idf=2 ! denotes tau +2 tau-
73  amfin=0.511d-3 !this mass is irrelevant if small, used in ME only
74  ELSEIF (idfx.EQ.-11) then
75  idf=-2 ! denotes tau -2 tau-
76  amfin=0.511d-3 !this mass is irrelevant if small, used in ME only
77  ELSEIF (idfx.EQ. 15) then
78  idf=2 ! denotes tau +2 tau-
79  amfin=1.77703 !this mass is irrelevant if small, used in ME only
80  ELSEIF (idfx.EQ.-15) then
81  idf=-2 ! denotes tau -2 tau-
82  amfin=1.77703 !this mass is irrelevant if small, used in ME only
83  ELSE
84  WRITE(*,*) 'INITWKSW: WRONG IDFX'
85  stop
86  ENDIF
87 
88  IF (idex.EQ. 11) then !electron
89  ide= 2
90  amin=0.511d-3
91  ELSEIF (idex.EQ.-11) then !positron
92  ide=-2
93  amin=0.511d-3
94  ELSEIF (idex.EQ. 13) then !mu+
95  ide= 2
96  amin=0.105659
97  ELSEIF (idex.EQ.-13) then !mu-
98  ide=-2
99  amin=0.105659
100  ELSEIF (idex.EQ. 1) then !d
101  ide= 4
102  amin=0.05
103  ELSEIF (idex.EQ.- 1) then !d~
104  ide=-4
105  amin=0.05
106  ELSEIF (idex.EQ. 2) then !u
107  ide= 3
108  amin=0.02
109  ELSEIF (idex.EQ.- 2) then !u~
110  ide=-3
111  amin=0.02
112  ELSEIF (idex.EQ. 3) then !s
113  ide= 4
114  amin=0.3
115  ELSEIF (idex.EQ.- 3) then !s~
116  ide=-4
117  amin=0.3
118  ELSEIF (idex.EQ. 4) then !c
119  ide= 3
120  amin=1.3
121  ELSEIF (idex.EQ.- 4) then !c~
122  ide=-3
123  amin=1.3
124  ELSEIF (idex.EQ. 5) then !b
125  ide= 4
126  amin=4.5
127  ELSEIF (idex.EQ.- 5) then !b~
128  ide=-4
129  amin=4.5
130  ELSEIF (idex.EQ. 12) then !nu_e
131  ide= 1
132  amin=0.1d-3
133  ELSEIF (idex.EQ.- 12) then !nu_e~
134  ide=-1
135  amin=0.1d-3
136  ELSEIF (idex.EQ. 14) then !nu_mu
137  ide= 1
138  amin=0.1d-3
139  ELSEIF (idex.EQ.- 14) then !nu_mu~
140  ide=-1
141  amin=0.1d-3
142  ELSEIF (idex.EQ. 16) then !nu_tau
143  ide= 1
144  amin=0.1d-3
145  ELSEIF (idex.EQ.- 16) then !nu_tau~
146  ide=-1
147  amin=0.1d-3
148 
149  ELSE
150  WRITE(*,*) 'INITWKSW: WRONG IDEX'
151  stop
152  ENDIF
153 
154 C ----------------------------------------------------------------------
155 C
156 C INITIALISATION OF COUPLING CONSTANTS AND FERMION-GAMMA / Z0 VERTEX
157 C
158 C called by : KORALZ
159 C ----------------------------------------------------------------------
160  itce=ide/iabs(ide)
161  jtce=(1-itce)/2
162  itcf=idf/iabs(idf)
163  jtcf=(1-itcf)/2
164  CALL t_givizo( ide, 1,aizor,qe,kdumm)
165  CALL t_givizo( ide,-1,aizol,qe,kdumm)
166  xupgi(1)=qe
167  xupgi(2)=qe
168  t3e = (aizol+aizor)/2.
169  xupzi(1)=(aizor-qe*(swsq+deltsq)*gsw(3)-qe*deltv)/sqrt(swsq*(1-swsq))
170  xupzi(2)=(aizol-qe*(swsq+deltsq)*gsw(3)-qe*deltv)/sqrt(swsq*(1-swsq))
171  ve =(xupzi(1)+xupzi(2))/2.
172  CALL t_givizo( idf, 1,aizor,qf,kolor)
173  CALL t_givizo( idf,-1,aizol,qf,kolor)
174  xupgf(1)=qf
175  xupgf(2)=qf
176  t3f = (aizol+aizor)/2.
177  xupzf(1)=(aizor-qf*(swsq+deltsq)*gsw(2)-qf*deltv)/sqrt(swsq*(1-swsq))
178  xupzf(2)=(aizol-qf*(swsq+deltsq)*gsw(2)-qf*deltv)/sqrt(swsq*(1-swsq))
179  vf =(xupzf(1)+xupzf(2))/2.
180 
181 * Coupling costants times EW form-factors
182  deno = dsqrt(swsq*(1d0-swsq))
183  ! Ve = (2*T3e -4*Qe*m_Sw2*CorEle)/Deno
184  ! Vf = (2*T3f -4*Qf*m_Sw2*CorFin)/Deno
185  ! Ae = 2*T3e /Deno
186  ! Af = 2*T3f /Deno
187 * Angle dependent double-vector extra-correction
188  vvcef = ( (t3e) *(t3f)
189  $ -(qe*swsq+deltsq) *(t3f) *gsw(3) -qe*(t3f)*deltv
190  $ -(qf*swsq+deltsq) *(t3e) *gsw(2) -qf*(t3e)*deltv
191  $ + (qe*swsq) *(qf*swsq) *gsw(4)
192  $ + 2*qe*qf*deltsq*swsq + 2*qe*qf*deltv*swsq )/deno**2
193 
194  vvcor = 1d0
195  IF(keygsw.NE.0.AND.keygsw.NE.4) THEN
196  vvcor = vvcef/(ve*vf)
197  ENDIF
198 C
199 C PRINT *,' initwksw VVCor = ', VVCor
200  ndiag0=2
201  ndiaga=11
202  keya = 1
203  keyz = 1
204 C
205 C
206  RETURN
207  END
208  FUNCTION t_bornew(MODE,KEYGSW,SVAR,COSTHE,TA,TB)
209 C ----------------------------------------------------------------------
210 C THIS ROUTINE PROVIDES BORN CROSS SECTION. IT HAS THE SAME
211 C STRUCTURE AS FUNTIS AND FUNTIH, THUS CAN BE USED AS SIMPLER
212 C EXAMPLE OF THE METHOD APPLIED THERE
213 C INPUT PARAMETERS ARE: SVAR -- transfer
214 C COSTHE -- cosine of angle between tau+ and 1st beam
215 C TA,TB -- helicity states of tau+ tau-
216 C mode -- parameter for mass terms; 1 means mass terms are on.
217 C keyGSW -- keyGSW=0 gamma propagator is off
218 C keyGSW=10 running Z width
219 C
220 C called by : BORNY, BORAS, BORNV, WAGA, WEIGHT
221 C ----------------------------------------------------------------------
222  IMPLICIT REAL*8(a-h,o-z)
223  COMMON / t_beampm / ene ,amin,amfin,ide,idf
224  REAL*8 ene ,amin,amfin
225  COMMON / t_gauspm /ss,poln,t3e,qe,t3f,qf
226  & ,xupgi0 ,xupzi0 ,xupgf0 ,xupzf0
227  & ,ndiag0,ndiaga,keya,keyz
228  & ,itce,jtce,itcf,jtcf,kolor
229  REAL*8 ss,poln,t3e,qe,t3f,qf
230  & ,xupgi0(2),xupzi0(2),xupgf0(2),xupzf0(2)
231  COMMON / t_gauspm1/vvcor, zetvpi, gamvpi
232  & ,xupgi ,xupzi ,xupgf ,xupzf
233  COMPLEX*16 vvcor, zetvpi, gamvpi
234  COMPLEX*16 xupgi(2),xupzi(2),xupgf(2),xupzf(2)
235  COMMON / t_ewn / gmun, alphainvn
236  REAL*8 gmun, alphainvn
237 
238 
239  REAL*8 seps1,seps2
240 C=====================================================================
241  COMMON / t_gswprmn /swsq,amw,amz,amh,amtop,gammz
242  REAL*8 swsq,amw,amz,amh,amtop,gammz
243 C SWSQ = sin2 (theta Weinberg)
244 C AMW,AMZ = W & Z boson masses respectively
245 C AMH = the Higgs mass
246 C AMTOP = the top mass
247 C GAMMZ = Z0 width
248  COMPLEX*16 aborn(2,2),aphot(2,2),azett(2,2)
249  COMPLEX*16 xupzfp(2),xupzip(2),xupzif(2,2)
250  COMPLEX*16 abornm(2,2),aphotm(2,2),azettm(2,2)
251  COMPLEX*16 propa,propz
252  COMPLEX*16 xr,xi
253  COMPLEX*16 xupf,xupi
254  COMPLEX*16 xthing
255  DATA xi/(0.d0,1.d0)/,xr/(1.d0,0.d0)/
256  DATA mode0 /-5/
257  DATA ide0 /-55/
258  DATA svar0,cost0 /-5.d0,-6.d0/
259  DATA pi /3.141592653589793238462643d0/
260  DATA seps1,seps2 /0d0,0d0/
261 
262 C
263 C MEMORIZATION =========================================================
264  IF ( mode.NE.mode0.OR.svar.NE.svar0.OR.costhe.NE.cost0
265  $ .OR.ide0.NE.ide)THEN
266 C
267 
268  ! PRINT *,' T_BORN EW loop ( ',sqrt(svar),XUPGI(1),')= ', VVcor, ZetVPi!, GamVPi
269  ! PRINT *,' T_BORN new( ',mode,')= ',SWSQ,AMW,AMZ,AMH,AMTOP,GAMMZ
270 C ** SWITCH OF MEMORISATION
271 C IDE0=IDE
272 C MODE0=MODE
273 C SVAR0=SVAR
274 C COST0=COSTHE
275 C ** PROPAGATORS
276  sinthe=sqrt(1.d0-costhe**2)
277  beta=sqrt(max(0d0,1d0-4d0*amfin**2/svar))
278 ! BETA=1.D0! Dec 10, 2019 mass term may need to be killed for EW tests
279 C I MULTIPLY AXIAL COUPLING BY BETA FACTOR.
280  xupzfp(1)=0.5d0*(xupzf(1)+xupzf(2))+0.5d0*beta*(xupzf(1)-xupzf(2))
281  xupzfp(2)=0.5d0*(xupzf(1)+xupzf(2))-0.5d0*beta*(xupzf(1)-xupzf(2))
282  xupzip(1)=0.5d0*(xupzi(1)+xupzi(2))+0.5d0*(xupzi(1)-xupzi(2))
283  xupzip(2)=0.5d0*(xupzi(1)+xupzi(2))-0.5d0*(xupzi(1)-xupzi(2))
284  xupzif(1,1)=(0.5d0*(xupzi(1)+xupzi(2))+0.5d0*(xupzi(1)-xupzi(2)))*(0.5d0*(xupzf(1)+xupzf(2))+0.5d0*beta*(xupzf(1)-xupzf(2)))
285  $ +(0.5d0*(xupzi(1)+xupzi(2)))*(0.5d0*(xupzf(1)+xupzf(2)))*(vvcor-1)
286  xupzif(1,2)=(0.5d0*(xupzi(1)+xupzi(2))+0.5d0*(xupzi(1)-xupzi(2)))*(0.5d0*(xupzf(1)+xupzf(2))-0.5d0*beta*(xupzf(1)-xupzf(2)))
287  $ +(0.5d0*(xupzi(1)+xupzi(2)))*(0.5d0*(xupzf(1)+xupzf(2)))*(vvcor-1)
288  xupzif(2,1)=(0.5d0*(xupzi(1)+xupzi(2))-0.5d0*(xupzi(1)-xupzi(2)))*(0.5d0*(xupzf(1)+xupzf(2))+0.5d0*beta*(xupzf(1)-xupzf(2)))
289  $ +(0.5d0*(xupzi(1)+xupzi(2)))*(0.5d0*(xupzf(1)+xupzf(2)))*(vvcor-1)
290  xupzif(2,2)=(0.5d0*(xupzi(1)+xupzi(2))-0.5d0*(xupzi(1)-xupzi(2)))*(0.5d0*(xupzf(1)+xupzf(2))-0.5d0*beta*(xupzf(1)-xupzf(2)))
291  $ +(0.5d0*(xupzi(1)+xupzi(2)))*(0.5d0*(xupzf(1)+xupzf(2)))*(vvcor-1)
292 
293 C FINAL STATE VECTOR COUPLING
294  xupf =0.5d0*(xupzf(1)+xupzf(2))
295  xupi =0.5d0*(xupzi(1)+xupzi(2))
296  xthing =0d0
297 
298 
299  propa =1d0/svar*gamvpi
300 C use running width
301  propz =1d0/dcmplx(svar-amz**2,svar/amz*gammz)*zetvpi
302 
303 
304  IF( keygsw. eq. 2) THEN
305  gfermi=gmun
306  alphainv=alphainvn
307  zetv = gfermi *amz**2 *alphainv /(dsqrt(2.d0)*8.d0*pi)
308  $ *(swsq*(1d0-swsq)) *16d0
309 
310 ! variants of the Z propagators for the non-ew case
311 ! ==1==
312  ! PROPZ =1D0/DCMPLX(SVAR-AMZ**2,AMZ*GAMMZ)*ZetV !default
313 ! ==2==
314 ! PROPZ =1D0/DCMPLX(SVAR-AMZ**2/(1+GAMMZ**2/AMZ**2), ! alternative as
315 ! $ AMZ*GAMMZ /(1+GAMMZ**2/AMZ**2) ) ! running width
316 ! $ *ZetV
317 ! PROPZ =PROPZ*DCMPLX(1,-GAMMZ/AMZ/(1+GAMMZ**2/AMZ**2))
318 ! ==3==
319  propz =1d0/dcmplx(svar-amz**2 , ! running
320  $ gammz*svar/amz )*zetv
321  ENDIF
322 
323 C use fixed width
324  IF( keygsw. eq. 10) THEN
325 ! PROPZ =1D0/DCMPLX(SVAR-AMZ**2,AMZ*GAMMZ)*ZetVPi ! this form need redefined M_Z and G_Z
326 ! below variant with this rescaling implemented
327  propz =1d0/dcmplx(svar-amz**2/(1+gammz**2/amz**2), ! alternative as
328  $ amz*gammz /(1+gammz**2/amz**2) ) ! running width
329  $ *zetv
330  propz =propz*dcmplx(1,-gammz/amz/(1+gammz**2/amz**2))
331 
332  ENDIF
333  IF (keygsw.EQ.0) propa=0.d0
334  DO 50 i=1,2
335  DO 50 j=1,2
336  regula= (3-2*i)*(3-2*j) + costhe
337  regulm=-(3-2*i)*(3-2*j) * sinthe *2.d0*amfin/sqrt(svar)
338  aphot(i,j)=propa*(xupgi(i)*xupgf(j)*regula)
339  azett(i,j)=propz*(xupzip(i)*xupzfp(j)+xthing)*regula
340  azett(i,j)=propz*(xupzif(i,j)+xthing)*regula ! with electroweak effects in.
341  aborn(i,j)=aphot(i,j)+azett(i,j)
342  aphotm(i,j)=propa*dcmplx(0d0,1d0)*xupgi(i)*xupgf(j)*regulm
343  azettm(i,j)=propz*dcmplx(0d0,1d0)*(xupzip(i)*xupf+xthing)*regulm
344  abornm(i,j)=aphotm(i,j)+azettm(i,j)
345  50 CONTINUE
346  ENDIF
347 C
348 C******************
349 C* IN CALCULATING CROSS SECTION ONLY DIAGONAL ELEMENTS
350 C* OF THE SPIN DENSITY MATRICES ENTER (LONGITUD. POL. ONLY.)
351 C* HELICITY CONSERVATION EXPLICITLY OBEYED
352  polar1= (seps1)
353  polar2= (-seps2)
354  born=0d0
355  DO 150 i=1,2
356  helic= 3-2*i
357  DO 150 j=1,2
358  helit=3-2*j
359  factor=kolor*(1d0+helic*polar1)*(1d0-helic*polar2)/4d0
360  factom=factor*(1+helit*ta)*(1-helit*tb)
361  factor=factor*(1+helit*ta)*(1+helit*tb)
362 
363  born=born+cdabs(aborn(i,j))**2*factor
364 C MASS TERM IN BORN
365  IF (mode.GE.1) THEN
366  born=born+cdabs(abornm(i,j))**2*factom
367  ENDIF
368 
369  150 CONTINUE
370 C************
371  funt=born
372  IF(funt.LT.0.d0) funt=born
373 
374 C
375  IF (svar.GT.4d0*amfin**2) THEN
376 C PHASE SPACE THRESHOLD FACTOR
377  thresh=sqrt(1-4d0*amfin**2/svar)
378  t_bornew= funt*svar**2*thresh
379  ELSE
380  thresh=0.d0
381  t_bornew=0.d0
382  ENDIF
383  END