| 0795afa3 |
1 | #include "isajet/pilot.h" |
| 2 | C ********************************************************************* |
| 3 | C *** *** |
| 4 | C *** coded by H. Murayama & I. Watanabe *** |
| 5 | C *** For the formalism and notations, see the following reference: *** |
| 6 | C *** H. Murayama, I. Watanabe and K. Hagiwara *** |
| 7 | C *** "HELAS: HELicity Amplitude Subroutines *** |
| 8 | C *** for Feynman diagram evaluation" *** |
| 9 | C *** KEK Report 91-11, December 1991 *** |
| 10 | C *** *** |
| 11 | C ********************************************************************* |
| 12 | C |
| 13 | C Converted to double precision by W. Long and T. Seltzer for MadGraph. |
| 14 | C |
| 15 | C Minor changes for portability by FEP, July 1999. The code is not ANSI |
| 16 | C standard, but that cannot be helped if MadGraph compatibility is to |
| 17 | C be maintained. |
| 18 | C |
| 19 | C ====================================================================== |
| 20 | C |
| 21 | SUBROUTINE BOOSTX(P,Q , PBOOST) |
| 22 | C |
| 23 | C this subroutine performs the lorentz boost of a four-momentum. the |
| 24 | C momentum p is assumed to be given in the rest frame of q. pboost is |
| 25 | C the momentum p boosted to the frame in which q is given. q must be a |
| 26 | C timelike momentum. |
| 27 | C |
| 28 | C input: |
| 29 | C real p(0:3) : four-momentum p in the q rest frame |
| 30 | C real q(0:3) : four-momentum q in the boosted frame |
| 31 | C |
| 32 | C output: |
| 33 | C real pboost(0:3) : four-momentum p in the boosted frame |
| 34 | C |
| 35 | #if defined(CERNLIB_IMPNONE) |
| 36 | IMPLICIT NONE |
| 37 | #endif |
| 38 | REAL*8 P(0:3),Q(0:3),PBOOST(0:3),PQ,QQ,M,LF |
| 39 | REAL*8 RXZERO |
| 40 | PARAMETER( RXZERO=0.0D0 ) |
| 41 | C |
| 42 | QQ=Q(1)**2+Q(2)**2+Q(3)**2 |
| 43 | C |
| 44 | IF ( QQ .NE. RXZERO ) THEN |
| 45 | PQ=P(1)*Q(1)+P(2)*Q(2)+P(3)*Q(3) |
| 46 | M=SQRT(Q(0)**2-QQ) |
| 47 | LF=((Q(0)-M)*PQ/QQ+P(0))/M |
| 48 | PBOOST(0) = (P(0)*Q(0)+PQ)/M |
| 49 | PBOOST(1) = P(1)+Q(1)*LF |
| 50 | PBOOST(2) = P(2)+Q(2)*LF |
| 51 | PBOOST(3) = P(3)+Q(3)*LF |
| 52 | ELSE |
| 53 | PBOOST(0)=P(0) |
| 54 | PBOOST(1)=P(1) |
| 55 | PBOOST(2)=P(2) |
| 56 | PBOOST(3)=P(3) |
| 57 | ENDIF |
| 58 | C |
| 59 | RETURN |
| 60 | END |
| 61 | C |
| 62 | C ********************************************************************** |
| 63 | C |
| 64 | SUBROUTINE COUP1X(SW2 , GW,GWWA,GWWZ) |
| 65 | C |
| 66 | C this subroutine sets up the coupling constants of the gauge bosons in |
| 67 | C the standard model. |
| 68 | C |
| 69 | C input: |
| 70 | C real sw2 : square of sine of the weak angle |
| 71 | C |
| 72 | C output: |
| 73 | C real gw : weak coupling constant |
| 74 | C real gwwa : dimensionless coupling of w-,w+,a |
| 75 | C real gwwz : dimensionless coupling of w-,w+,z |
| 76 | C |
| 77 | #if defined(CERNLIB_IMPNONE) |
| 78 | IMPLICIT NONE |
| 79 | #endif |
| 80 | REAL*8 SW2,GW,GWWA,GWWZ,ALPHA,FOURPI,EE,SW,CW |
| 81 | REAL*8 RXONE, RXFOUR, RXOTE, RXPI, RIALPH |
| 82 | PARAMETER( RXONE=1.0D0, RXFOUR=4.0D0, RXOTE=128.0D0 ) |
| 83 | PARAMETER( RXPI=3.14159265358979323846D0, RIALPH=137.0359895D0 ) |
| 84 | C |
| 85 | ALPHA = RXONE / RXOTE |
| 86 | C alpha = r_one / r_ialph |
| 87 | FOURPI = RXFOUR * RXPI |
| 88 | EE=SQRT( ALPHA * FOURPI ) |
| 89 | SW=SQRT( SW2 ) |
| 90 | CW=SQRT( RXONE - SW2 ) |
| 91 | C |
| 92 | GW = EE/SW |
| 93 | GWWA = EE |
| 94 | GWWZ = EE*CW/SW |
| 95 | C |
| 96 | RETURN |
| 97 | END |
| 98 | C |
| 99 | C ---------------------------------------------------------------------- |
| 100 | C |
| 101 | SUBROUTINE COUP2X(SW2 , GAL,GAU,GAD,GWF,GZN,GZL,GZU,GZD,G1) |
| 102 | C |
| 103 | C this subroutine sets up the coupling constants for the fermion- |
| 104 | C fermion-vector vertices in the standard model. the array of the |
| 105 | C couplings specifies the chirality of the flowing-in fermion. g??(1) |
| 106 | C denotes a left-handed coupling, and g??(2) a right-handed coupling. |
| 107 | C |
| 108 | C input: |
| 109 | C real sw2 : square of sine of the weak angle |
| 110 | C |
| 111 | C output: |
| 112 | C real gal(2) : coupling with a of charged leptons |
| 113 | C real gau(2) : coupling with a of up-type quarks |
| 114 | C real gad(2) : coupling with a of down-type quarks |
| 115 | C real gwf(2) : coupling with w-,w+ of fermions |
| 116 | C real gzn(2) : coupling with z of neutrinos |
| 117 | C real gzl(2) : coupling with z of charged leptons |
| 118 | C real gzu(2) : coupling with z of up-type quarks |
| 119 | C real gzd(2) : coupling with z of down-type quarks |
| 120 | C real g1(2) : unit coupling of fermions |
| 121 | C |
| 122 | #if defined(CERNLIB_IMPNONE) |
| 123 | IMPLICIT NONE |
| 124 | #endif |
| 125 | REAL*8 GAL(2),GAU(2),GAD(2),GWF(2),GZN(2),GZL(2),GZU(2),GZD(2), |
| 126 | & G1(2),SW2,ALPHA,FOURPI,EE,SW,CW,EZ,EY |
| 127 | C |
| 128 | REAL*8 RXZERO, RXHALF, RXONE, RXTWO, RTHREE, RXFOUR, RXOTE |
| 129 | REAL*8 RXPI, RIALPH |
| 130 | PARAMETER( RXZERO=0.0D0, RXHALF=0.5D0, RXONE=1.0D0, RXTWO=2.0D0, |
| 131 | $ RTHREE=3.0D0 ) |
| 132 | PARAMETER( RXFOUR=4.0D0, RXOTE=128.0D0 ) |
| 133 | PARAMETER( RXPI=3.14159265358979323846D0, RIALPH=137.0359895D0 ) |
| 134 | C |
| 135 | ALPHA = RXONE / RXOTE |
| 136 | C alpha = r_one / r_ialph |
| 137 | FOURPI = RXFOUR * RXPI |
| 138 | EE=SQRT( ALPHA * FOURPI ) |
| 139 | SW=SQRT( SW2 ) |
| 140 | CW=SQRT( RXONE - SW2 ) |
| 141 | EZ=EE/(SW*CW) |
| 142 | EY=EE*(SW/CW) |
| 143 | C |
| 144 | GAL(1) = EE |
| 145 | GAL(2) = EE |
| 146 | GAU(1) = -EE*RXTWO/RTHREE |
| 147 | GAU(2) = -EE*RXTWO/RTHREE |
| 148 | GAD(1) = EE /RTHREE |
| 149 | GAD(2) = EE /RTHREE |
| 150 | GWF(1) = -EE/SQRT(RXTWO*SW2) |
| 151 | GWF(2) = RXZERO |
| 152 | GZN(1) = -EZ* RXHALF |
| 153 | GZN(2) = RXZERO |
| 154 | GZL(1) = -EZ*(-RXHALF+SW2) |
| 155 | GZL(2) = -EY |
| 156 | GZU(1) = -EZ*( RXHALF-SW2*RXTWO/RTHREE) |
| 157 | GZU(2) = EY* RXTWO/RTHREE |
| 158 | GZD(1) = -EZ*(-RXHALF+SW2 /RTHREE) |
| 159 | GZD(2) = -EY /RTHREE |
| 160 | G1(1) = RXONE |
| 161 | G1(2) = RXONE |
| 162 | C |
| 163 | RETURN |
| 164 | END |
| 165 | C |
| 166 | C ---------------------------------------------------------------------- |
| 167 | C |
| 168 | SUBROUTINE COUP3X(SW2,ZMASS,HMASS , |
| 169 | & GWWH,GZZH,GHHH,GWWHH,GZZHH,GHHHH) |
| 170 | C |
| 171 | C this subroutine sets up the coupling constants of the gauge bosons and |
| 172 | C higgs boson in the standard model. |
| 173 | C |
| 174 | C input: |
| 175 | C real sw2 : square of sine of the weak angle |
| 176 | C real zmass : mass of z |
| 177 | C real hmass : mass of higgs |
| 178 | C |
| 179 | C output: |
| 180 | C real gwwh : dimensionful coupling of w-,w+,h |
| 181 | C real gzzh : dimensionful coupling of z, z, h |
| 182 | C real ghhh : dimensionful coupling of h, h, h |
| 183 | C real gwwhh : dimensionful coupling of w-,w+,h, h |
| 184 | C real gzzhh : dimensionful coupling of z, z, h, h |
| 185 | C real ghhhh : dimensionless coupling of h, h, h, h |
| 186 | C |
| 187 | #if defined(CERNLIB_IMPNONE) |
| 188 | IMPLICIT NONE |
| 189 | #endif |
| 190 | REAL*8 SW2,ZMASS,HMASS,GWWH,GZZH,GHHH,GWWHH,GZZHH,GHHHH, |
| 191 | & ALPHA,FOURPI,EE2,SC2,V |
| 192 | C |
| 193 | REAL*8 RXHALF, RXONE, RXTWO, RTHREE, RXFOUR, RXOTE |
| 194 | REAL*8 RXPI, RIALPH |
| 195 | PARAMETER( RXHALF=0.5D0, RXONE=1.0D0, RXTWO=2.0D0, RTHREE=3.0D0 ) |
| 196 | PARAMETER( RXFOUR=4.0D0, RXOTE=128.0D0 ) |
| 197 | PARAMETER( RXPI=3.14159265358979323846D0, RIALPH=137.0359895D0 ) |
| 198 | C |
| 199 | ALPHA = RXONE / RXOTE |
| 200 | C alpha = r_one / r_ialph |
| 201 | FOURPI = RXFOUR * RXPI |
| 202 | EE2=ALPHA*FOURPI |
| 203 | SC2=SW2*( RXONE - SW2 ) |
| 204 | V = RXTWO * ZMASS*SQRT(SC2)/SQRT(EE2) |
| 205 | C |
| 206 | GWWH = EE2/SW2*RXHALF*V |
| 207 | GZZH = EE2/SC2*RXHALF*V |
| 208 | GHHH = -HMASS**2/V*RTHREE |
| 209 | GWWHH = EE2/SW2*RXHALF |
| 210 | GZZHH = EE2/SC2*RXHALF |
| 211 | GHHHH = -(HMASS/V)**2*RTHREE |
| 212 | C |
| 213 | RETURN |
| 214 | END |
| 215 | C |
| 216 | C ---------------------------------------------------------------------- |
| 217 | C |
| 218 | SUBROUTINE COUP4X(SW2,ZMASS,FMASS , GCHF) |
| 219 | C |
| 220 | C This subroutine sets up the coupling constant for the fermion-fermion- |
| 221 | C Higgs vertex in the STANDARD MODEL. The coupling is COMPLEX and the |
| 222 | C array of the coupling specifies the chirality of the flowing-IN |
| 223 | C fermion. GCHF(1) denotes a left-handed coupling, and GCHF(2) a right- |
| 224 | C handed coupling. |
| 225 | C |
| 226 | C INPUT: |
| 227 | C real SW2 : square of sine of the weak angle |
| 228 | C real ZMASS : Z mass |
| 229 | C real FMASS : fermion mass |
| 230 | C |
| 231 | C OUTPUT: |
| 232 | C complex GCHF(2) : coupling of fermion and Higgs |
| 233 | C |
| 234 | #if defined(CERNLIB_IMPNONE) |
| 235 | IMPLICIT NONE |
| 236 | #endif |
| 237 | COMPLEX*16 GCHF(2) |
| 238 | REAL*8 SW2,ZMASS,FMASS,ALPHA,FOURPI,EZ,G |
| 239 | C |
| 240 | ALPHA=1.D0/128.D0 |
| 241 | C ALPHA=1./REAL(137.0359895) |
| 242 | FOURPI=4.D0*3.14159265358979323846D0 |
| 243 | EZ=SQRT(ALPHA*FOURPI)/SQRT(SW2*(1.D0-SW2)) |
| 244 | G=EZ*FMASS*0.5D0/ZMASS |
| 245 | C |
| 246 | GCHF(1) = DCMPLX( -G ) |
| 247 | GCHF(2) = DCMPLX( -G ) |
| 248 | C |
| 249 | RETURN |
| 250 | END |
| 251 | C |
| 252 | C ====================================================================== |
| 253 | C |
| 254 | SUBROUTINE EAIXXX(EB,EA,SHLF,CHLF,PHI,NHE,NHA , EAI) |
| 255 | C |
| 256 | C This subroutine computes an off-shell electron wavefunction after |
| 257 | C emitting a photon from the electron beam, with a special care for the |
| 258 | C small angle region. The momenta are measured in the laboratory frame, |
| 259 | C where the e- beam is along the positive z axis. |
| 260 | C |
| 261 | C INPUT: |
| 262 | C real EB : energy (GeV) of beam e- |
| 263 | C real EA : energy (GeV) of final photon |
| 264 | C real SHLF : sin(theta/2) of final photon |
| 265 | C real CHLF : cos(theta/2) of final photon |
| 266 | C real PHI : azimuthal angle of final photon |
| 267 | C integer NHE = -1 or 1 : helicity of beam e- |
| 268 | C integer NHA = -1 or 1 : helicity of final photon |
| 269 | C |
| 270 | C OUTPUT: |
| 271 | C complex EAI(6) : off-shell electron |e',A,e> |
| 272 | C |
| 273 | #if defined(CERNLIB_IMPNONE) |
| 274 | IMPLICIT NONE |
| 275 | #endif |
| 276 | COMPLEX*16 EAI(6),PHS |
| 277 | REAL*8 EB,EA,SHLF,CHLF,PHI,ME,ALPHA,GAL,RNHE,X,C,S,D,COEFF, |
| 278 | & XNNP,XNNM,SNP,CSP |
| 279 | INTEGER NHE,NHA,NN |
| 280 | C |
| 281 | ME = 0.51099906D-3 |
| 282 | ALPHA=1./128. |
| 283 | GAL =SQRT(ALPHA*4.*3.14159265D0) |
| 284 | C |
| 285 | NN=NHA*NHE |
| 286 | RNHE=NHE |
| 287 | X=EA/EB |
| 288 | C=(CHLF+SHLF)*(CHLF-SHLF) |
| 289 | S=2.*CHLF*SHLF |
| 290 | D=-1./(EA*EB*(4.*SHLF**2+(ME/EB)**2*C)) |
| 291 | COEFF=-NN*GAL*SQRT(EB)*D |
| 292 | XNNP=X*(1+NN) |
| 293 | XNNM=X*(1-NN) |
| 294 | SNP=SIN(PHI) |
| 295 | CSP=COS(PHI) |
| 296 | PHS=DCMPLX( CSP , RNHE*SNP ) |
| 297 | C |
| 298 | EAI((5-3*NHE)/2) = -RNHE*COEFF*ME*S*(1.+XNNP*.5) |
| 299 | EAI((5-NHE)/2) = XNNP*COEFF*ME*CHLF**2*PHS |
| 300 | EAI((5+NHE)/2) = RNHE*COEFF*EB*S*(-2.+XNNM) |
| 301 | EAI((5+3*NHE)/2) = XNNM*COEFF*EB*SHLF**2*PHS*2. |
| 302 | C |
| 303 | EAI(5) = EB*DCMPLX( 1.-X , 1.-X*C ) |
| 304 | EAI(6) = -EB*X*S*DCMPLX( CSP , SNP ) |
| 305 | C |
| 306 | RETURN |
| 307 | END |
| 308 | C |
| 309 | C ---------------------------------------------------------------------- |
| 310 | C |
| 311 | SUBROUTINE EAOXXX(EB,EA,SHLF,CHLF,PHI,NHE,NHA , EAO) |
| 312 | C |
| 313 | C This subroutine computes an off-shell positron wavefunction after |
| 314 | C emitting a photon from the positron beam, with a special care for the |
| 315 | C small angle region. The momenta are measured in the laboratory frame, |
| 316 | C where the e+ beam is along the negative z axis. |
| 317 | C |
| 318 | C INPUT: |
| 319 | C real EB : energy (GeV) of beam e+ |
| 320 | C real EA : energy (GeV) of final photon |
| 321 | C real SHLF : sin(theta/2) of final photon |
| 322 | C real CHLF : cos(theta/2) of final photon |
| 323 | C real PHI : azimuthal angle of final photon |
| 324 | C integer NHE = -1 or 1 : helicity of beam e+ |
| 325 | C integer NHA = -1 or 1 : helicity of final photon |
| 326 | C |
| 327 | C OUTPUT: |
| 328 | C complex EAO(6) : off-shell positron <e,A,e'| |
| 329 | C |
| 330 | #if defined(CERNLIB_IMPNONE) |
| 331 | IMPLICIT NONE |
| 332 | #endif |
| 333 | COMPLEX*16 EAO(6),PHS |
| 334 | REAL*8 EB,EA,SHLF,CHLF,PHI,ME,ALPHA,GAL,RNHE,X,C,S,D,COEFF, |
| 335 | & XNNP,XNNM,SNP,CSP |
| 336 | INTEGER NHE,NHA,NN |
| 337 | C |
| 338 | ME = 0.51099906D-3 |
| 339 | ALPHA=1./128. |
| 340 | GAL =SQRT(ALPHA*4.*3.14159265D0) |
| 341 | C |
| 342 | NN=NHA*NHE |
| 343 | RNHE=NHE |
| 344 | X=EA/EB |
| 345 | C=(CHLF+SHLF)*(CHLF-SHLF) |
| 346 | S=2.*CHLF*SHLF |
| 347 | D=-1./(EA*EB*(4.*CHLF**2-(ME/EB)**2*C)) |
| 348 | COEFF=NN*GAL*SQRT(EB)*D |
| 349 | XNNP=X*(1+NN) |
| 350 | XNNM=X*(1-NN) |
| 351 | SNP=SIN(PHI) |
| 352 | CSP=COS(PHI) |
| 353 | PHS=DCMPLX( CSP ,-RNHE*SNP ) |
| 354 | C |
| 355 | EAO((5-3*NHE)/2) = COEFF*ME*S*(1.+XNNP*.5) |
| 356 | EAO((5-NHE)/2) = RNHE*XNNP *COEFF*ME*SHLF**2*PHS |
| 357 | EAO((5+NHE)/2) = COEFF*EB*S*(-2.+XNNM) |
| 358 | EAO((5+3*NHE)/2) = REAL(NHA-NHE)*COEFF*EB*X*CHLF**2*PHS*2. |
| 359 | C |
| 360 | EAO(5) = EB*DCMPLX( X-1. , X*C+1. ) |
| 361 | EAO(6) = EB*X*S*DCMPLX( CSP , SNP ) |
| 362 | C |
| 363 | RETURN |
| 364 | END |
| 365 | C |
| 366 | C ---------------------------------------------------------------------- |
| 367 | C |
| 368 | SUBROUTINE FSIXXX(FI,SC,GC,FMASS,FWIDTH , FSI) |
| 369 | C |
| 370 | C this subroutine computes an off-shell fermion wavefunction from a |
| 371 | C flowing-in external fermion and a vector boson. |
| 372 | C |
| 373 | C input: |
| 374 | C complex*16 fi(6) : flow-in fermion |fi> |
| 375 | C complex*16 sc(3) : input scalar s |
| 376 | C complex*16 gc(2) : coupling constants gchf |
| 377 | C real*8 fmass : mass of output fermion f' |
| 378 | C real*8 fwidth : width of output fermion f' |
| 379 | C |
| 380 | C output: |
| 381 | C complex fsi(6) : off-shell fermion |f',s,fi> |
| 382 | C |
| 383 | #if defined(CERNLIB_IMPNONE) |
| 384 | IMPLICIT NONE |
| 385 | #endif |
| 386 | COMPLEX*16 FI(6),SC(3),FSI(6),GC(2),SL1,SL2,SR1,SR2,DS |
| 387 | REAL*8 PF(0:3),FMASS,FWIDTH,PF2,P0P3,P0M3 |
| 388 | C |
| 389 | FSI(5) = FI(5)-SC(2) |
| 390 | FSI(6) = FI(6)-SC(3) |
| 391 | C |
| 392 | PF(0)=DBLE( FSI(5)) |
| 393 | PF(1)=DBLE( FSI(6)) |
| 394 | PF(2)=DIMAG(FSI(6)) |
| 395 | PF(3)=DIMAG(FSI(5)) |
| 396 | PF2=PF(0)**2-(PF(1)**2+PF(2)**2+PF(3)**2) |
| 397 | C |
| 398 | DS=-SC(1)/DCMPLX(PF2-FMASS**2,MAX(DSIGN(FMASS*FWIDTH ,PF2),0D0)) |
| 399 | P0P3=PF(0)+PF(3) |
| 400 | P0M3=PF(0)-PF(3) |
| 401 | SL1=GC(1)*(P0P3*FI(1)+DCONJG(FSI(6))*FI(2)) |
| 402 | SL2=GC(1)*(P0M3*FI(2) +FSI(6) *FI(1)) |
| 403 | SR1=GC(2)*(P0M3*FI(3)-DCONJG(FSI(6))*FI(4)) |
| 404 | SR2=GC(2)*(P0P3*FI(4) -FSI(6) *FI(3)) |
| 405 | C |
| 406 | FSI(1) = ( GC(1)*FMASS*FI(1) + SR1 )*DS |
| 407 | FSI(2) = ( GC(1)*FMASS*FI(2) + SR2 )*DS |
| 408 | FSI(3) = ( GC(2)*FMASS*FI(3) + SL1 )*DS |
| 409 | FSI(4) = ( GC(2)*FMASS*FI(4) + SL2 )*DS |
| 410 | C |
| 411 | RETURN |
| 412 | END |
| 413 | C |
| 414 | C ---------------------------------------------------------------------- |
| 415 | C |
| 416 | SUBROUTINE FSOXXX(FO,SC,GC,FMASS,FWIDTH , FSO) |
| 417 | C |
| 418 | C this subroutine computes an off-shell fermion wavefunction from a |
| 419 | C flowing-out external fermion and a vector boson. |
| 420 | C |
| 421 | C input: |
| 422 | C complex*16 fo(6) : flow-out fermion <fo| |
| 423 | C complex*16 sc(6) : input scalar s |
| 424 | C complex*16 gc(2) : coupling constants gchf |
| 425 | C real*8 fmass : mass of output fermion f' |
| 426 | C real*8 fwidth : width of output fermion f' |
| 427 | C |
| 428 | C output: |
| 429 | C complex fso(6) : off-shell fermion <fo,s,f'| |
| 430 | C |
| 431 | #if defined(CERNLIB_IMPNONE) |
| 432 | IMPLICIT NONE |
| 433 | #endif |
| 434 | COMPLEX*16 FO(6),SC(6),FSO(6),GC(2),SL1,SL2,SR1,SR2,DS |
| 435 | REAL*8 PF(0:3),FMASS,FWIDTH,PF2,P0P3,P0M3 |
| 436 | C |
| 437 | FSO(5) = FO(5)+SC(2) |
| 438 | FSO(6) = FO(6)+SC(3) |
| 439 | C |
| 440 | PF(0)=DBLE( FSO(5)) |
| 441 | PF(1)=DBLE( FSO(6)) |
| 442 | PF(2)=DIMAG(FSO(6)) |
| 443 | PF(3)=DIMAG(FSO(5)) |
| 444 | PF2=PF(0)**2-(PF(1)**2+PF(2)**2+PF(3)**2) |
| 445 | C |
| 446 | DS=-SC(1)/DCMPLX(PF2-FMASS**2,MAX(DSIGN(FMASS*FWIDTH ,PF2),0D0)) |
| 447 | P0P3=PF(0)+PF(3) |
| 448 | P0M3=PF(0)-PF(3) |
| 449 | SL1=GC(2)*(P0P3*FO(3) +FSO(6) *FO(4)) |
| 450 | SL2=GC(2)*(P0M3*FO(4)+DCONJG(FSO(6))*FO(3)) |
| 451 | SR1=GC(1)*(P0M3*FO(1) -FSO(6) *FO(2)) |
| 452 | SR2=GC(1)*(P0P3*FO(2)-DCONJG(FSO(6))*FO(1)) |
| 453 | C |
| 454 | FSO(1) = ( GC(1)*FMASS*FO(1) + SL1 )*DS |
| 455 | FSO(2) = ( GC(1)*FMASS*FO(2) + SL2 )*DS |
| 456 | FSO(3) = ( GC(2)*FMASS*FO(3) + SR1 )*DS |
| 457 | FSO(4) = ( GC(2)*FMASS*FO(4) + SR2 )*DS |
| 458 | C |
| 459 | RETURN |
| 460 | END |
| 461 | C |
| 462 | C ---------------------------------------------------------------------- |
| 463 | C |
| 464 | SUBROUTINE FVIXXX(FI,VC,G,FMASS,FWIDTH , FVI) |
| 465 | C |
| 466 | C this subroutine computes an off-shell fermion wavefunction from a |
| 467 | C flowing-in external fermion and a vector boson. |
| 468 | C |
| 469 | C input: |
| 470 | C complex fi(6) : flow-in fermion |fi> |
| 471 | C complex vc(6) : input vector v |
| 472 | C real g(2) : coupling constants gvf |
| 473 | C real fmass : mass of output fermion f' |
| 474 | C real fwidth : width of output fermion f' |
| 475 | C |
| 476 | C output: |
| 477 | C complex fvi(6) : off-shell fermion |f',v,fi> |
| 478 | C |
| 479 | #if defined(CERNLIB_IMPNONE) |
| 480 | IMPLICIT NONE |
| 481 | #endif |
| 482 | COMPLEX*16 FI(6),VC(6),FVI(6),SL1,SL2,SR1,SR2,D |
| 483 | REAL*8 G(2),PF(0:3),FMASS,FWIDTH,PF2 |
| 484 | REAL*8 RXZERO, RXONE |
| 485 | PARAMETER( RXZERO=0.0D0, RXONE=1.0D0 ) |
| 486 | COMPLEX*16 CXIMAG |
| 487 | C |
| 488 | LOGICAL FIRST |
| 489 | SAVE CXIMAG,FIRST |
| 490 | DATA FIRST/.TRUE./ |
| 491 | C |
| 492 | C Fix compilation with g77 |
| 493 | IF(FIRST) THEN |
| 494 | FIRST=.FALSE. |
| 495 | CXIMAG=DCMPLX( RXZERO, RXONE ) |
| 496 | ENDIF |
| 497 | C |
| 498 | FVI(5) = FI(5)-VC(5) |
| 499 | FVI(6) = FI(6)-VC(6) |
| 500 | C |
| 501 | PF(0)=DBLE( FVI(5)) |
| 502 | PF(1)=DBLE( FVI(6)) |
| 503 | PF(2)=DIMAG(FVI(6)) |
| 504 | PF(3)=DIMAG(FVI(5)) |
| 505 | PF2=PF(0)**2-(PF(1)**2+PF(2)**2+PF(3)**2) |
| 506 | C |
| 507 | D=-RXONE/DCMPLX( PF2-FMASS**2,MAX(SIGN(FMASS*FWIDTH,PF2),RXZERO)) |
| 508 | SL1= (VC(1)+ VC(4))*FI(1) |
| 509 | & +(VC(2)-CXIMAG*VC(3))*FI(2) |
| 510 | SL2= (VC(2)+CXIMAG*VC(3))*FI(1) |
| 511 | & +(VC(1)- VC(4))*FI(2) |
| 512 | C |
| 513 | IF ( G(2) .NE. RXZERO ) THEN |
| 514 | SR1= (VC(1)- VC(4))*FI(3) |
| 515 | & -(VC(2)-CXIMAG*VC(3))*FI(4) |
| 516 | SR2=-(VC(2)+CXIMAG*VC(3))*FI(3) |
| 517 | & +(VC(1)+ VC(4))*FI(4) |
| 518 | C |
| 519 | FVI(1) = ( G(1)*((PF(0)-PF(3))*SL1 -DCONJG(FVI(6))*SL2) |
| 520 | & +G(2)*FMASS*SR1)*D |
| 521 | FVI(2) = ( G(1)*( -FVI(6)*SL1 +(PF(0)+PF(3))*SL2) |
| 522 | & +G(2)*FMASS*SR2)*D |
| 523 | FVI(3) = ( G(2)*((PF(0)+PF(3))*SR1 +DCONJG(FVI(6))*SR2) |
| 524 | & +G(1)*FMASS*SL1)*D |
| 525 | FVI(4) = ( G(2)*( FVI(6)*SR1 +(PF(0)-PF(3))*SR2) |
| 526 | & +G(1)*FMASS*SL2)*D |
| 527 | C |
| 528 | ELSE |
| 529 | FVI(1) = G(1)*((PF(0)-PF(3))*SL1 -DCONJG(FVI(6))*SL2)*D |
| 530 | FVI(2) = G(1)*( -FVI(6)*SL1 +(PF(0)+PF(3))*SL2)*D |
| 531 | FVI(3) = G(1)*FMASS*SL1*D |
| 532 | FVI(4) = G(1)*FMASS*SL2*D |
| 533 | END IF |
| 534 | C |
| 535 | RETURN |
| 536 | END |
| 537 | C |
| 538 | C ---------------------------------------------------------------------- |
| 539 | C |
| 540 | SUBROUTINE FVOXXX(FO,VC,G,FMASS,FWIDTH , FVO) |
| 541 | C |
| 542 | C this subroutine computes an off-shell fermion wavefunction from a |
| 543 | C flowing-out external fermion and a vector boson. |
| 544 | C |
| 545 | C input: |
| 546 | C complex fo(6) : flow-out fermion <fo| |
| 547 | C complex vc(6) : input vector v |
| 548 | C real g(2) : coupling constants gvf |
| 549 | C real fmass : mass of output fermion f' |
| 550 | C real fwidth : width of output fermion f' |
| 551 | C |
| 552 | C output: |
| 553 | C complex fvo(6) : off-shell fermion <fo,v,f'| |
| 554 | C |
| 555 | #if defined(CERNLIB_IMPNONE) |
| 556 | IMPLICIT NONE |
| 557 | #endif |
| 558 | COMPLEX*16 FO(6),VC(6),FVO(6),SL1,SL2,SR1,SR2,D |
| 559 | REAL*8 G(2),PF(0:3),FMASS,FWIDTH,PF2 |
| 560 | REAL*8 RXZERO, RXONE |
| 561 | PARAMETER( RXZERO=0.0D0, RXONE=1.0D0 ) |
| 562 | COMPLEX*16 CXIMAG |
| 563 | LOGICAL FIRST |
| 564 | SAVE CXIMAG,FIRST |
| 565 | DATA FIRST/.TRUE./ |
| 566 | C |
| 567 | C Fix compilation with g77 |
| 568 | IF(FIRST) THEN |
| 569 | FIRST=.FALSE. |
| 570 | CXIMAG=DCMPLX( RXZERO, RXONE ) |
| 571 | ENDIF |
| 572 | C |
| 573 | FVO(5) = FO(5)+VC(5) |
| 574 | FVO(6) = FO(6)+VC(6) |
| 575 | C |
| 576 | PF(0)=DBLE( FVO(5)) |
| 577 | PF(1)=DBLE( FVO(6)) |
| 578 | PF(2)=DIMAG(FVO(6)) |
| 579 | PF(3)=DIMAG(FVO(5)) |
| 580 | PF2=PF(0)**2-(PF(1)**2+PF(2)**2+PF(3)**2) |
| 581 | C |
| 582 | D=-RXONE/DCMPLX( PF2-FMASS**2,MAX(SIGN(FMASS*FWIDTH,PF2),RXZERO)) |
| 583 | SL1= (VC(1)+ VC(4))*FO(3) |
| 584 | & +(VC(2)+CXIMAG*VC(3))*FO(4) |
| 585 | SL2= (VC(2)-CXIMAG*VC(3))*FO(3) |
| 586 | & +(VC(1)- VC(4))*FO(4) |
| 587 | C |
| 588 | IF ( G(2) .NE. RXZERO ) THEN |
| 589 | SR1= (VC(1)- VC(4))*FO(1) |
| 590 | & -(VC(2)+CXIMAG*VC(3))*FO(2) |
| 591 | SR2=-(VC(2)-CXIMAG*VC(3))*FO(1) |
| 592 | & +(VC(1)+ VC(4))*FO(2) |
| 593 | C |
| 594 | FVO(1) = ( G(2)*( (PF(0)+PF(3))*SR1 +FVO(6)*SR2) |
| 595 | & +G(1)*FMASS*SL1)*D |
| 596 | FVO(2) = ( G(2)*( DCONJG(FVO(6))*SR1 +(PF(0)-PF(3))*SR2) |
| 597 | & +G(1)*FMASS*SL2)*D |
| 598 | FVO(3) = ( G(1)*( (PF(0)-PF(3))*SL1 -FVO(6)*SL2) |
| 599 | & +G(2)*FMASS*SR1)*D |
| 600 | FVO(4) = ( G(1)*(-DCONJG(FVO(6))*SL1 +(PF(0)+PF(3))*SL2) |
| 601 | & +G(2)*FMASS*SR2)*D |
| 602 | C |
| 603 | ELSE |
| 604 | FVO(1) = G(1)*FMASS*SL1*D |
| 605 | FVO(2) = G(1)*FMASS*SL2*D |
| 606 | FVO(3) = G(1)*( (PF(0)-PF(3))*SL1 -FVO(6)*SL2)*D |
| 607 | FVO(4) = G(1)*(-DCONJG(FVO(6))*SL1 +(PF(0)+PF(3))*SL2)*D |
| 608 | END IF |
| 609 | C |
| 610 | RETURN |
| 611 | END |
| 612 | C |
| 613 | C ---------------------------------------------------------------------- |
| 614 | C |
| 615 | SUBROUTINE GGGGXX(WM,W31,WP,W32,G, VERTEX) |
| 616 | C |
| 617 | C this subroutine computes an amplitude of the four-point coupling of |
| 618 | C the w-, w+ and two w3/z/a. the amplitude includes the contributions |
| 619 | C of w exchange diagrams. the internal w propagator is given in unitary |
| 620 | C gauge. if one sets wmass=0.0, then the gggg vertex is given (see sect |
| 621 | C 2.9.1 of the manual). |
| 622 | C |
| 623 | C input: |
| 624 | C complex wm(0:3) : flow-out w- wm |
| 625 | C complex w31(0:3) : first w3/z/a w31 |
| 626 | C complex wp(0:3) : flow-out w+ wp |
| 627 | C complex w32(0:3) : second w3/z/a w32 |
| 628 | C real g : coupling of w31 with w-/w+ |
| 629 | C (see the table below) |
| 630 | C |
| 631 | C the possible sets of the inputs are as follows: |
| 632 | C ------------------------------------------- |
| 633 | C | wm | w31 | wp | w32 | g31 | g32 | |
| 634 | C ------------------------------------------- |
| 635 | C | w- | w3 | w+ | w3 | gw | gw | |
| 636 | C | w- | w3 | w+ | z | gw | gwwz | |
| 637 | C | w- | w3 | w+ | a | gw | gwwa | |
| 638 | C | w- | z | w+ | z | gwwz | gwwz | |
| 639 | C | w- | z | w+ | a | gwwz | gwwa | |
| 640 | C | w- | a | w+ | a | gwwa | gwwa | |
| 641 | C ------------------------------------------- |
| 642 | C where all the bosons are defined by the flowing-out quantum number. |
| 643 | C |
| 644 | C output: |
| 645 | C complex vertex : amplitude gamma(wm,w31,wp,w32) |
| 646 | C |
| 647 | #if defined(CERNLIB_IMPNONE) |
| 648 | IMPLICIT NONE |
| 649 | #endif |
| 650 | COMPLEX*16 WM(6),W31(6),WP(6),W32(6),VERTEX |
| 651 | COMPLEX*16 DV1(0:3),DV2(0:3),DV3(0:3),DV4(0:3), |
| 652 | & DVERTX,V12,V13,V14,V23,V24,V34 |
| 653 | REAL*8 PWM(0:3),PW31(0:3),PWP(0:3),PW32(0:3),G |
| 654 | REAL*8 DP1(0:3),DP2(0:3),DP3(0:3),DP4(0:3) |
| 655 | REAL*8 RXZERO, RXONE |
| 656 | PARAMETER( RXZERO=0.0D0, RXONE=1.0D0 ) |
| 657 | C |
| 658 | PWM(0)=DBLE( WM(5)) |
| 659 | PWM(1)=DBLE( WM(6)) |
| 660 | PWM(2)=DIMAG(WM(6)) |
| 661 | PWM(3)=DIMAG(WM(5)) |
| 662 | PWP(0)=DBLE( WP(5)) |
| 663 | PWP(1)=DBLE( WP(6)) |
| 664 | PWP(2)=DIMAG(WP(6)) |
| 665 | PWP(3)=DIMAG(WP(5)) |
| 666 | PW31(0)=DBLE( W31(5)) |
| 667 | PW31(1)=DBLE( W31(6)) |
| 668 | PW31(2)=DIMAG(W31(6)) |
| 669 | PW31(3)=DIMAG(W31(5)) |
| 670 | PW32(0)=DBLE( W32(5)) |
| 671 | PW32(1)=DBLE( W32(6)) |
| 672 | PW32(2)=DIMAG(W32(6)) |
| 673 | PW32(3)=DIMAG(W32(5)) |
| 674 | C |
| 675 | DV1(0)=DCMPLX(WM(1)) |
| 676 | DV1(1)=DCMPLX(WM(2)) |
| 677 | DV1(2)=DCMPLX(WM(3)) |
| 678 | DV1(3)=DCMPLX(WM(4)) |
| 679 | DP1(0)=DBLE(PWM(0)) |
| 680 | DP1(1)=DBLE(PWM(1)) |
| 681 | DP1(2)=DBLE(PWM(2)) |
| 682 | DP1(3)=DBLE(PWM(3)) |
| 683 | DV2(0)=DCMPLX(W31(1)) |
| 684 | DV2(1)=DCMPLX(W31(2)) |
| 685 | DV2(2)=DCMPLX(W31(3)) |
| 686 | DV2(3)=DCMPLX(W31(4)) |
| 687 | DP2(0)=DBLE(PW31(0)) |
| 688 | DP2(1)=DBLE(PW31(1)) |
| 689 | DP2(2)=DBLE(PW31(2)) |
| 690 | DP2(3)=DBLE(PW31(3)) |
| 691 | DV3(0)=DCMPLX(WP(1)) |
| 692 | DV3(1)=DCMPLX(WP(2)) |
| 693 | DV3(2)=DCMPLX(WP(3)) |
| 694 | DV3(3)=DCMPLX(WP(4)) |
| 695 | DP3(0)=DBLE(PWP(0)) |
| 696 | DP3(1)=DBLE(PWP(1)) |
| 697 | DP3(2)=DBLE(PWP(2)) |
| 698 | DP3(3)=DBLE(PWP(3)) |
| 699 | DV4(0)=DCMPLX(W32(1)) |
| 700 | DV4(1)=DCMPLX(W32(2)) |
| 701 | DV4(2)=DCMPLX(W32(3)) |
| 702 | DV4(3)=DCMPLX(W32(4)) |
| 703 | DP4(0)=DBLE(PW32(0)) |
| 704 | DP4(1)=DBLE(PW32(1)) |
| 705 | DP4(2)=DBLE(PW32(2)) |
| 706 | DP4(3)=DBLE(PW32(3)) |
| 707 | C |
| 708 | V12= DV1(0)*DV2(0)-DV1(1)*DV2(1)-DV1(2)*DV2(2)-DV1(3)*DV2(3) |
| 709 | V13= DV1(0)*DV3(0)-DV1(1)*DV3(1)-DV1(2)*DV3(2)-DV1(3)*DV3(3) |
| 710 | V14= DV1(0)*DV4(0)-DV1(1)*DV4(1)-DV1(2)*DV4(2)-DV1(3)*DV4(3) |
| 711 | V23= DV2(0)*DV3(0)-DV2(1)*DV3(1)-DV2(2)*DV3(2)-DV2(3)*DV3(3) |
| 712 | V24= DV2(0)*DV4(0)-DV2(1)*DV4(1)-DV2(2)*DV4(2)-DV2(3)*DV4(3) |
| 713 | V34= DV3(0)*DV4(0)-DV3(1)*DV4(1)-DV3(2)*DV4(2)-DV3(3)*DV4(3) |
| 714 | |
| 715 | DVERTX = V14*V23 -V13*V24 |
| 716 | C |
| 717 | VERTEX = DCMPLX( DVERTX ) * (G*G) |
| 718 | C |
| 719 | RETURN |
| 720 | END |
| 721 | C |
| 722 | C ====================================================================== |
| 723 | C |
| 724 | SUBROUTINE GGGXXX(WM,WP,W3,G , VERTEX) |
| 725 | C |
| 726 | C this subroutine computes an amplitude of the three-point coupling of |
| 727 | C the gauge bosons. |
| 728 | C |
| 729 | C input: |
| 730 | C complex wm(6) : vector flow-out w- |
| 731 | C complex wp(6) : vector flow-out w+ |
| 732 | C complex w3(6) : vector j3 or a or z |
| 733 | C real g : coupling constant gw or gwwa or gwwz |
| 734 | C |
| 735 | C output: |
| 736 | C complex vertex : amplitude gamma(wm,wp,w3) |
| 737 | C |
| 738 | #if defined(CERNLIB_IMPNONE) |
| 739 | IMPLICIT NONE |
| 740 | #endif |
| 741 | COMPLEX*16 WM(6),WP(6),W3(6),VERTEX, |
| 742 | & XV1,XV2,XV3,V12,V23,V31,P12,P13,P21,P23,P31,P32 |
| 743 | REAL*8 PWM(0:3),PWP(0:3),PW3(0:3),G |
| 744 | REAL*8 RXZERO, RTENTH |
| 745 | PARAMETER( RXZERO=0.0D0, RTENTH=0.1D0 ) |
| 746 | C |
| 747 | PWM(0)=DBLE( WM(5)) |
| 748 | PWM(1)=DBLE( WM(6)) |
| 749 | PWM(2)=DIMAG(WM(6)) |
| 750 | PWM(3)=DIMAG(WM(5)) |
| 751 | PWP(0)=DBLE( WP(5)) |
| 752 | PWP(1)=DBLE( WP(6)) |
| 753 | PWP(2)=DIMAG(WP(6)) |
| 754 | PWP(3)=DIMAG(WP(5)) |
| 755 | PW3(0)=DBLE( W3(5)) |
| 756 | PW3(1)=DBLE( W3(6)) |
| 757 | PW3(2)=DIMAG(W3(6)) |
| 758 | PW3(3)=DIMAG(W3(5)) |
| 759 | C |
| 760 | V12=WM(1)*WP(1)-WM(2)*WP(2)-WM(3)*WP(3)-WM(4)*WP(4) |
| 761 | V23=WP(1)*W3(1)-WP(2)*W3(2)-WP(3)*W3(3)-WP(4)*W3(4) |
| 762 | V31=W3(1)*WM(1)-W3(2)*WM(2)-W3(3)*WM(3)-W3(4)*WM(4) |
| 763 | XV1=RXZERO |
| 764 | XV2=RXZERO |
| 765 | XV3=RXZERO |
| 766 | IF ( ABS(WM(1)) .NE. RXZERO ) THEN |
| 767 | IF (ABS(WM(1)).GE.MAX(ABS(WM(2)),ABS(WM(3)),ABS(WM(4))) |
| 768 | $ *RTENTH) |
| 769 | & XV1=PWM(0)/WM(1) |
| 770 | ENDIF |
| 771 | IF ( ABS(WP(1)) .NE. RXZERO) THEN |
| 772 | IF (ABS(WP(1)).GE.MAX(ABS(WP(2)),ABS(WP(3)),ABS(WP(4))) |
| 773 | $ *RTENTH) |
| 774 | & XV2=PWP(0)/WP(1) |
| 775 | ENDIF |
| 776 | IF ( ABS(W3(1)) .NE. RXZERO) THEN |
| 777 | IF ( ABS(W3(1)).GE.MAX(ABS(W3(2)),ABS(W3(3)),ABS(W3(4))) |
| 778 | $ *RTENTH) |
| 779 | & XV3=PW3(0)/W3(1) |
| 780 | ENDIF |
| 781 | P12= (PWM(0)-XV1*WM(1))*WP(1)-(PWM(1)-XV1*WM(2))*WP(2) |
| 782 | & -(PWM(2)-XV1*WM(3))*WP(3)-(PWM(3)-XV1*WM(4))*WP(4) |
| 783 | P13= (PWM(0)-XV1*WM(1))*W3(1)-(PWM(1)-XV1*WM(2))*W3(2) |
| 784 | & -(PWM(2)-XV1*WM(3))*W3(3)-(PWM(3)-XV1*WM(4))*W3(4) |
| 785 | P21= (PWP(0)-XV2*WP(1))*WM(1)-(PWP(1)-XV2*WP(2))*WM(2) |
| 786 | & -(PWP(2)-XV2*WP(3))*WM(3)-(PWP(3)-XV2*WP(4))*WM(4) |
| 787 | P23= (PWP(0)-XV2*WP(1))*W3(1)-(PWP(1)-XV2*WP(2))*W3(2) |
| 788 | & -(PWP(2)-XV2*WP(3))*W3(3)-(PWP(3)-XV2*WP(4))*W3(4) |
| 789 | P31= (PW3(0)-XV3*W3(1))*WM(1)-(PW3(1)-XV3*W3(2))*WM(2) |
| 790 | & -(PW3(2)-XV3*W3(3))*WM(3)-(PW3(3)-XV3*W3(4))*WM(4) |
| 791 | P32= (PW3(0)-XV3*W3(1))*WP(1)-(PW3(1)-XV3*W3(2))*WP(2) |
| 792 | & -(PW3(2)-XV3*W3(3))*WP(3)-(PW3(3)-XV3*W3(4))*WP(4) |
| 793 | C |
| 794 | VERTEX = -(V12*(P13-P23)+V23*(P21-P31)+V31*(P32-P12))*G |
| 795 | C |
| 796 | RETURN |
| 797 | END |
| 798 | SUBROUTINE HIOXXX(FI,FO,GC,SMASS,SWIDTH , HIO) |
| 799 | C |
| 800 | C this subroutine computes an off-shell scalar current from an external |
| 801 | C fermion pair. |
| 802 | C |
| 803 | C input: |
| 804 | C complex fi(6) : flow-in fermion |fi> |
| 805 | C complex fo(6) : flow-out fermion <fo| |
| 806 | C complex gc(2) : coupling constants gchf |
| 807 | C real smass : mass of output scalar s |
| 808 | C real swidth : width of output scalar s |
| 809 | C |
| 810 | C output: |
| 811 | C complex hio(3) : scalar current j(<fi|s|fo>) |
| 812 | C |
| 813 | #if defined(CERNLIB_IMPNONE) |
| 814 | IMPLICIT NONE |
| 815 | #endif |
| 816 | COMPLEX*16 FI(6),FO(6),HIO(3),GC(2),DN |
| 817 | REAL*8 Q(0:3),SMASS,SWIDTH,Q2 |
| 818 | C |
| 819 | HIO(2) = FO(5)-FI(5) |
| 820 | HIO(3) = FO(6)-FI(6) |
| 821 | C |
| 822 | Q(0)=DBLE( HIO(2)) |
| 823 | Q(1)=DBLE( HIO(3)) |
| 824 | Q(2)=DIMAG(HIO(3)) |
| 825 | Q(3)=DIMAG(HIO(2)) |
| 826 | Q2=Q(0)**2-(Q(1)**2+Q(2)**2+Q(3)**2) |
| 827 | C |
| 828 | DN=-DCMPLX(Q2-SMASS**2,DMAX1(DSIGN(SMASS*SWIDTH,Q2),0.D0)) |
| 829 | C |
| 830 | HIO(1) = ( GC(1)*(FO(1)*FI(1)+FO(2)*FI(2)) |
| 831 | & +GC(2)*(FO(3)*FI(3)+FO(4)*FI(4)) )/DN |
| 832 | C |
| 833 | RETURN |
| 834 | END |
| 835 | C ---------------------------------------------------------------------- |
| 836 | C |
| 837 | SUBROUTINE HSSSXX(S1,S2,S3,G,SMASS,SWIDTH , HSSS) |
| 838 | C |
| 839 | C This subroutine computes an off-shell scalar current from the four- |
| 840 | C scalar coupling. |
| 841 | C |
| 842 | C INPUT: |
| 843 | C complex S1(3) : first scalar S1 |
| 844 | C complex S2(3) : second scalar S2 |
| 845 | C complex S3(3) : third scalar S3 |
| 846 | C real G : coupling constant GHHHH |
| 847 | C real SMASS : mass of OUTPUT scalar S' |
| 848 | C real SWIDTH : width of OUTPUT scalar S' |
| 849 | C |
| 850 | C OUTPUT: |
| 851 | C complex HSSS(3) : scalar current J(S':S1,S2,S3) |
| 852 | C |
| 853 | #if defined(CERNLIB_IMPNONE) |
| 854 | IMPLICIT NONE |
| 855 | #endif |
| 856 | COMPLEX*16 S1(3),S2(3),S3(3),HSSS(3),DG |
| 857 | REAL*8 Q(0:3),G,SMASS,SWIDTH,Q2 |
| 858 | C |
| 859 | HSSS(2) = S1(2)+S2(2)+S3(2) |
| 860 | HSSS(3) = S1(3)+S2(3)+S3(3) |
| 861 | C |
| 862 | Q(0)=DBLE( HSSS(2)) |
| 863 | Q(1)=DBLE( HSSS(3)) |
| 864 | Q(2)=DIMAG(HSSS(3)) |
| 865 | Q(3)=DIMAG(HSSS(2)) |
| 866 | Q2=Q(0)**2-(Q(1)**2+Q(2)**2+Q(3)**2) |
| 867 | C |
| 868 | DG=-G/DCMPLX( Q2-SMASS**2,MAX(SIGN(SMASS*SWIDTH ,Q2),0.D0)) |
| 869 | C |
| 870 | HSSS(1) = DG * S1(1)*S2(1)*S3(1) |
| 871 | C |
| 872 | RETURN |
| 873 | END |
| 874 | C ---------------------------------------------------------------------- |
| 875 | C |
| 876 | SUBROUTINE HSSXXX(S1,S2,G,SMASS,SWIDTH , HSS) |
| 877 | C |
| 878 | C This subroutine computes an off-shell scalar current from the three- |
| 879 | C scalar coupling. |
| 880 | C |
| 881 | C INPUT: |
| 882 | C complex S1(3) : first scalar S1 |
| 883 | C complex S2(3) : second scalar S2 |
| 884 | C real G : coupling constant GHHH |
| 885 | C real SMASS : mass of OUTPUT scalar S' |
| 886 | C real SWIDTH : width of OUTPUT scalar S' |
| 887 | C |
| 888 | C OUTPUT: |
| 889 | C complex HSS(3) : scalar current J(S':S1,S2) |
| 890 | C |
| 891 | #if defined(CERNLIB_IMPNONE) |
| 892 | IMPLICIT NONE |
| 893 | #endif |
| 894 | COMPLEX*16 S1(3),S2(3),HSS(3),DG |
| 895 | REAL*8 Q(0:3),G,SMASS,SWIDTH,Q2 |
| 896 | C |
| 897 | HSS(2) = S1(2)+S2(2) |
| 898 | HSS(3) = S1(3)+S2(3) |
| 899 | C |
| 900 | Q(0)=DBLE( HSS(2)) |
| 901 | Q(1)=DBLE( HSS(3)) |
| 902 | Q(2)=DIMAG(HSS(3)) |
| 903 | Q(3)=DIMAG(HSS(2)) |
| 904 | Q2=Q(0)**2-(Q(1)**2+Q(2)**2+Q(3)**2) |
| 905 | C |
| 906 | DG=-G/DCMPLX( Q2-SMASS**2, MAX(SIGN(SMASS*SWIDTH ,Q2),0.D0)) |
| 907 | C |
| 908 | HSS(1) = DG*S1(1)*S2(1) |
| 909 | C |
| 910 | RETURN |
| 911 | END |
| 912 | C |
| 913 | C ====================================================================== |
| 914 | C ---------------------------------------------------------------------- |
| 915 | C |
| 916 | SUBROUTINE HVSXXX(VC,SC,G,SMASS,SWIDTH , HVS) |
| 917 | C |
| 918 | C this subroutine computes an off-shell scalar current from the vector- |
| 919 | C scalar-scalar coupling. the coupling is absent in the minimal sm in |
| 920 | C unitary gauge. |
| 921 | C |
| 922 | C input: |
| 923 | C complex vc(6) : input vector v |
| 924 | C complex sc(3) : input scalar s |
| 925 | C complex g : coupling constant (s charge) |
| 926 | C real smass : mass of output scalar s' |
| 927 | C real swidth : width of output scalar s' |
| 928 | C |
| 929 | C examples of the coupling constant g for susy particles are as follows: |
| 930 | C ----------------------------------------------------------- |
| 931 | C | s1 | (q,i3) of s1 || v=a | v=z | v=w | |
| 932 | C ----------------------------------------------------------- |
| 933 | C | nu~_l | ( 0 , +1/2) || --- | gzn(1) | gwf(1) | |
| 934 | C | e~_l | ( -1 , -1/2) || gal(1) | gzl(1) | gwf(1) | |
| 935 | C | u~_l | (+2/3 , +1/2) || gau(1) | gzu(1) | gwf(1) | |
| 936 | C | d~_l | (-1/3 , -1/2) || gad(1) | gzd(1) | gwf(1) | |
| 937 | C ----------------------------------------------------------- |
| 938 | C | e~_r-bar | ( +1 , 0 ) || -gal(2) | -gzl(2) | -gwf(2) | |
| 939 | C | u~_r-bar | (-2/3 , 0 ) || -gau(2) | -gzu(2) | -gwf(2) | |
| 940 | C | d~_r-bar | (+1/3 , 0 ) || -gad(2) | -gzd(2) | -gwf(2) | |
| 941 | C ----------------------------------------------------------- |
| 942 | C where the sc charge is defined by the flowing-out quantum number. |
| 943 | C |
| 944 | C output: |
| 945 | C complex hvs(3) : scalar current j(s':v,s) |
| 946 | C |
| 947 | #if defined(CERNLIB_IMPNONE) |
| 948 | IMPLICIT NONE |
| 949 | #endif |
| 950 | COMPLEX*16 VC(6),SC(3),HVS(3),DG,QVV,QPV,G |
| 951 | REAL*8 QV(0:3),QP(0:3),QA(0:3),SMASS,SWIDTH,Q2 |
| 952 | C |
| 953 | HVS(2) = VC(5)+SC(2) |
| 954 | HVS(3) = VC(6)+SC(3) |
| 955 | C |
| 956 | QV(0)=DBLE( VC(5)) |
| 957 | QV(1)=DBLE( VC(6)) |
| 958 | QV(2)=DIMAG( VC(6)) |
| 959 | QV(3)=DIMAG( VC(5)) |
| 960 | QP(0)=DBLE( SC(2)) |
| 961 | QP(1)=DBLE( SC(3)) |
| 962 | QP(2)=DIMAG( SC(3)) |
| 963 | QP(3)=DIMAG( SC(2)) |
| 964 | QA(0)=DBLE( HVS(2)) |
| 965 | QA(1)=DBLE( HVS(3)) |
| 966 | QA(2)=DIMAG(HVS(3)) |
| 967 | QA(3)=DIMAG(HVS(2)) |
| 968 | Q2=QA(0)**2-(QA(1)**2+QA(2)**2+QA(3)**2) |
| 969 | C |
| 970 | DG=-G/DCMPLX( Q2-SMASS**2 , MAX(DSIGN( SMASS*SWIDTH ,Q2),0D0) ) |
| 971 | QVV=QV(0)*VC(1)-QV(1)*VC(2)-QV(2)*VC(3)-QV(3)*VC(4) |
| 972 | QPV=QP(0)*VC(1)-QP(1)*VC(2)-QP(2)*VC(3)-QP(3)*VC(4) |
| 973 | C |
| 974 | HVS(1) = DG*(2D0*QPV+QVV)*SC(1) |
| 975 | C |
| 976 | RETURN |
| 977 | END |
| 978 | C |
| 979 | C ---------------------------------------------------------------------- |
| 980 | C |
| 981 | SUBROUTINE HVVXXX(V1,V2,G,SMASS,SWIDTH , HVV) |
| 982 | C |
| 983 | C this subroutine computes an off-shell scalar current from the vector- |
| 984 | C vector-scalar coupling. |
| 985 | C |
| 986 | C input: |
| 987 | C complex v1(6) : first vector v1 |
| 988 | C complex v2(6) : second vector v2 |
| 989 | C real g : coupling constant gvvh |
| 990 | C real smass : mass of output scalar s |
| 991 | C real swidth : width of output scalar s |
| 992 | C |
| 993 | C output: |
| 994 | C complex hvv(3) : off-shell scalar current j(s:v1,v2) |
| 995 | C |
| 996 | #if defined(CERNLIB_IMPNONE) |
| 997 | IMPLICIT NONE |
| 998 | #endif |
| 999 | COMPLEX*16 V1(6),V2(6),HVV(3),DG |
| 1000 | REAL*8 Q(0:3),G,SMASS,SWIDTH,Q2 |
| 1001 | REAL*8 RXZERO |
| 1002 | PARAMETER( RXZERO=0.0D0 ) |
| 1003 | C |
| 1004 | HVV(2) = V1(5)+V2(5) |
| 1005 | HVV(3) = V1(6)+V2(6) |
| 1006 | C |
| 1007 | Q(0)=DBLE( HVV(2)) |
| 1008 | Q(1)=DBLE( HVV(3)) |
| 1009 | Q(2)=DIMAG(HVV(3)) |
| 1010 | Q(3)=DIMAG(HVV(2)) |
| 1011 | Q2=Q(0)**2-(Q(1)**2+Q(2)**2+Q(3)**2) |
| 1012 | C |
| 1013 | DG=-G/DCMPLX( Q2-SMASS**2 , MAX(SIGN( SMASS*SWIDTH ,Q2),RXZERO) ) |
| 1014 | C |
| 1015 | HVV(1) = DG*(V1(1)*V2(1)-V1(2)*V2(2)-V1(3)*V2(3)-V1(4)*V2(4)) |
| 1016 | C |
| 1017 | RETURN |
| 1018 | END |
| 1019 | C |
| 1020 | C ====================================================================== |
| 1021 | C |
| 1022 | SUBROUTINE IOSXXX(FI,FO,SC,GC , VERTEX) |
| 1023 | C |
| 1024 | C This subroutine computes an amplitude of the fermion-fermion-scalar |
| 1025 | C coupling. |
| 1026 | C |
| 1027 | C INPUT: |
| 1028 | C complex FI(6) : flow-in fermion |FI> |
| 1029 | C complex FO(6) : flow-out fermion <FO| |
| 1030 | C complex SC(3) : input scalar S |
| 1031 | C complex GC(2) : coupling constants GCHF |
| 1032 | C |
| 1033 | C OUTPUT: |
| 1034 | C complex VERTEX : amplitude <FO|S|FI> |
| 1035 | C |
| 1036 | #if defined(CERNLIB_IMPNONE) |
| 1037 | IMPLICIT NONE |
| 1038 | #endif |
| 1039 | COMPLEX*16 FI(6),FO(6),SC(3),GC(2),VERTEX |
| 1040 | C |
| 1041 | VERTEX = SC(1)*( GC(1)*(FI(1)*FO(1)+FI(2)*FO(2)) |
| 1042 | & +GC(2)*(FI(3)*FO(3)+FI(4)*FO(4)) ) |
| 1043 | C |
| 1044 | RETURN |
| 1045 | END |
| 1046 | C |
| 1047 | C ====================================================================== |
| 1048 | C |
| 1049 | SUBROUTINE IOVXXX(FI,FO,VC,G , VERTEX) |
| 1050 | C |
| 1051 | C this subroutine computes an amplitude of the fermion-fermion-vector |
| 1052 | C coupling. |
| 1053 | C |
| 1054 | C input: |
| 1055 | C complex fi(6) : flow-in fermion |fi> |
| 1056 | C complex fo(6) : flow-out fermion <fo| |
| 1057 | C complex vc(6) : input vector v |
| 1058 | C real g(2) : coupling constants gvf |
| 1059 | C |
| 1060 | C output: |
| 1061 | C complex vertex : amplitude <fo|v|fi> |
| 1062 | C |
| 1063 | #if defined(CERNLIB_IMPNONE) |
| 1064 | IMPLICIT NONE |
| 1065 | #endif |
| 1066 | COMPLEX*16 FI(6),FO(6),VC(6),VERTEX |
| 1067 | REAL*8 G(2) |
| 1068 | REAL*8 RXZERO, RXONE |
| 1069 | PARAMETER( RXZERO=0.0D0, RXONE=1.0D0 ) |
| 1070 | COMPLEX*16 CXIMAG |
| 1071 | LOGICAL FIRST |
| 1072 | SAVE CXIMAG,FIRST |
| 1073 | DATA FIRST/.TRUE./ |
| 1074 | C |
| 1075 | C Fix compilation with g77 |
| 1076 | IF(FIRST) THEN |
| 1077 | FIRST=.FALSE. |
| 1078 | CXIMAG=DCMPLX( RXZERO, RXONE ) |
| 1079 | ENDIF |
| 1080 | C |
| 1081 | |
| 1082 | VERTEX = G(1)*( (FO(3)*FI(1)+FO(4)*FI(2))*VC(1) |
| 1083 | & +(FO(3)*FI(2)+FO(4)*FI(1))*VC(2) |
| 1084 | & -(FO(3)*FI(2)-FO(4)*FI(1))*VC(3)*CXIMAG |
| 1085 | & +(FO(3)*FI(1)-FO(4)*FI(2))*VC(4) ) |
| 1086 | C |
| 1087 | IF ( G(2) .NE. RXZERO ) THEN |
| 1088 | VERTEX = VERTEX |
| 1089 | & + G(2)*( (FO(1)*FI(3)+FO(2)*FI(4))*VC(1) |
| 1090 | & -(FO(1)*FI(4)+FO(2)*FI(3))*VC(2) |
| 1091 | & +(FO(1)*FI(4)-FO(2)*FI(3))*VC(3)*CXIMAG |
| 1092 | & -(FO(1)*FI(3)-FO(2)*FI(4))*VC(4) ) |
| 1093 | END IF |
| 1094 | C |
| 1095 | RETURN |
| 1096 | END |
| 1097 | C |
| 1098 | C Subroutine returns the desired fermion or |
| 1099 | C anti-fermion spinor. ie., |f> |
| 1100 | C A replacement for the HELAS routine IXXXXX |
| 1101 | C |
| 1102 | C Adam Duff, 1992 August 31 |
| 1103 | C <duff@phenom.physics.wisc.edu> |
| 1104 | C |
| 1105 | SUBROUTINE IXXXXX(P,FMASS,NHEL,NSF,FI) |
| 1106 | C P IN: FOUR VECTOR MOMENTUM |
| 1107 | C FMASS IN: FERMION MASS |
| 1108 | C NHEL IN: SPINOR HELICITY, -1 OR 1 |
| 1109 | C NSF IN: -1=ANTIFERMION, 1=FERMION |
| 1110 | C FI OUT: FERMION WAVEFUNCTION |
| 1111 | C |
| 1112 | C declare input/output variables |
| 1113 | C |
| 1114 | #if defined(CERNLIB_IMPNONE) |
| 1115 | IMPLICIT NONE |
| 1116 | #endif |
| 1117 | COMPLEX*16 FI(6) |
| 1118 | INTEGER*4 NHEL, NSF |
| 1119 | REAL*8 P(0:3), FMASS |
| 1120 | REAL*8 RXZERO, RXONE, RXTWO |
| 1121 | PARAMETER( RXZERO=0.0D0, RXONE=1.0D0, RXTWO=2.0D0 ) |
| 1122 | REAL*8 PLAT, PABS, OMEGAP, OMEGAM, RS2PA, SPAZ |
| 1123 | COMPLEX*16 CXZERO |
| 1124 | C |
| 1125 | C declare local variables |
| 1126 | C |
| 1127 | LOGICAL FIRST |
| 1128 | SAVE CXZERO,FIRST |
| 1129 | DATA FIRST/.TRUE./ |
| 1130 | C |
| 1131 | C Fix compilation with g77 |
| 1132 | IF(FIRST) THEN |
| 1133 | FIRST=.FALSE. |
| 1134 | CXZERO=DCMPLX( RXZERO, RXZERO ) |
| 1135 | ENDIF |
| 1136 | C |
| 1137 | C define kinematic parameters |
| 1138 | C |
| 1139 | FI(5) = DCMPLX( P(0), P(3) ) * NSF |
| 1140 | FI(6) = DCMPLX( P(1), P(2) ) * NSF |
| 1141 | PLAT = SQRT( P(1)**2 + P(2)**2 ) |
| 1142 | PABS = SQRT( P(1)**2 + P(2)**2 + P(3)**2 ) |
| 1143 | OMEGAP = SQRT( P(0) + PABS ) |
| 1144 | C |
| 1145 | C do massive fermion case |
| 1146 | C |
| 1147 | IF ( FMASS .NE. RXZERO ) THEN |
| 1148 | OMEGAM = FMASS / OMEGAP |
| 1149 | IF ( NSF .EQ. 1 ) THEN |
| 1150 | IF ( NHEL .EQ. 1 ) THEN |
| 1151 | IF ( P(3) .GE. RXZERO ) THEN |
| 1152 | IF ( PLAT .EQ. RXZERO ) THEN |
| 1153 | FI(1) = DCMPLX( OMEGAM, RXZERO ) |
| 1154 | FI(2) = CXZERO |
| 1155 | FI(3) = DCMPLX( OMEGAP, RXZERO ) |
| 1156 | FI(4) = CXZERO |
| 1157 | ELSE |
| 1158 | RS2PA = RXONE / SQRT( RXTWO * PABS ) |
| 1159 | SPAZ = SQRT( PABS + P(3) ) |
| 1160 | FI(1) = OMEGAM * RS2PA |
| 1161 | & * DCMPLX( SPAZ, RXZERO ) |
| 1162 | FI(2) = OMEGAM * RS2PA / SPAZ |
| 1163 | & * DCMPLX( P(1), P(2) ) |
| 1164 | FI(3) = OMEGAP * RS2PA |
| 1165 | & * DCMPLX( SPAZ, RXZERO ) |
| 1166 | FI(4) = OMEGAP * RS2PA / SPAZ |
| 1167 | & * DCMPLX( P(1), P(2) ) |
| 1168 | END IF |
| 1169 | ELSE |
| 1170 | IF ( PLAT .EQ. RXZERO ) THEN |
| 1171 | FI(1) = CXZERO |
| 1172 | FI(2) = DCMPLX( OMEGAM, RXZERO ) |
| 1173 | FI(3) = CXZERO |
| 1174 | FI(4) = DCMPLX( OMEGAP, RXZERO ) |
| 1175 | ELSE |
| 1176 | RS2PA = RXONE / SQRT( RXTWO * PABS ) |
| 1177 | SPAZ = SQRT( PABS - P(3) ) |
| 1178 | FI(1) = OMEGAM * RS2PA / SPAZ |
| 1179 | & * DCMPLX( PLAT, RXZERO ) |
| 1180 | FI(2) = OMEGAM * RS2PA * SPAZ / PLAT |
| 1181 | & * DCMPLX( P(1), P(2) ) |
| 1182 | FI(3) = OMEGAP * RS2PA / SPAZ |
| 1183 | & * DCMPLX( PLAT, RXZERO ) |
| 1184 | FI(4) = OMEGAP * RS2PA * SPAZ / PLAT |
| 1185 | & * DCMPLX( P(1), P(2) ) |
| 1186 | END IF |
| 1187 | END IF |
| 1188 | ELSE IF ( NHEL .EQ. -1 ) THEN |
| 1189 | IF ( P(3) .GE. RXZERO ) THEN |
| 1190 | IF ( PLAT .EQ. RXZERO ) THEN |
| 1191 | FI(1) = CXZERO |
| 1192 | FI(2) = DCMPLX( OMEGAP, RXZERO ) |
| 1193 | FI(3) = CXZERO |
| 1194 | FI(4) = DCMPLX( OMEGAM, RXZERO ) |
| 1195 | ELSE |
| 1196 | RS2PA = RXONE / SQRT( RXTWO * PABS ) |
| 1197 | SPAZ = SQRT( PABS + P(3) ) |
| 1198 | FI(1) = OMEGAP * RS2PA / SPAZ |
| 1199 | & * DCMPLX( -P(1), P(2) ) |
| 1200 | FI(2) = OMEGAP * RS2PA |
| 1201 | & * DCMPLX( SPAZ, RXZERO ) |
| 1202 | FI(3) = OMEGAM * RS2PA / SPAZ |
| 1203 | & * DCMPLX( -P(1), P(2) ) |
| 1204 | FI(4) = OMEGAM * RS2PA |
| 1205 | & * DCMPLX( SPAZ, RXZERO ) |
| 1206 | END IF |
| 1207 | ELSE |
| 1208 | IF ( PLAT .EQ. RXZERO ) THEN |
| 1209 | FI(1) = DCMPLX( -OMEGAP, RXZERO ) |
| 1210 | FI(2) = CXZERO |
| 1211 | FI(3) = DCMPLX( -OMEGAM, RXZERO ) |
| 1212 | FI(4) = CXZERO |
| 1213 | ELSE |
| 1214 | RS2PA = RXONE / SQRT( RXTWO * PABS ) |
| 1215 | SPAZ = SQRT( PABS - P(3) ) |
| 1216 | FI(1) = OMEGAP * RS2PA * SPAZ / PLAT |
| 1217 | & * DCMPLX( -P(1), P(2) ) |
| 1218 | FI(2) = OMEGAP * RS2PA / SPAZ |
| 1219 | & * DCMPLX( PLAT, RXZERO ) |
| 1220 | FI(3) = OMEGAM * RS2PA * SPAZ / PLAT |
| 1221 | & * DCMPLX( -P(1), P(2) ) |
| 1222 | FI(4) = OMEGAM * RS2PA / SPAZ |
| 1223 | & * DCMPLX( PLAT, RXZERO ) |
| 1224 | END IF |
| 1225 | END IF |
| 1226 | ELSE |
| 1227 | STOP 'IXXXXX: FERMION HELICITY MUST BE +1,-1' |
| 1228 | END IF |
| 1229 | ELSE IF ( NSF .EQ. -1 ) THEN |
| 1230 | IF ( NHEL .EQ. 1 ) THEN |
| 1231 | IF ( P(3) .GE. RXZERO ) THEN |
| 1232 | IF ( PLAT .EQ. RXZERO ) THEN |
| 1233 | FI(1) = CXZERO |
| 1234 | FI(2) = DCMPLX( -OMEGAP, RXZERO ) |
| 1235 | FI(3) = CXZERO |
| 1236 | FI(4) = DCMPLX( OMEGAM, RXZERO ) |
| 1237 | ELSE |
| 1238 | RS2PA = RXONE / SQRT( RXTWO * PABS ) |
| 1239 | SPAZ = SQRT( PABS + P(3) ) |
| 1240 | FI(1) = -OMEGAP * RS2PA / SPAZ |
| 1241 | & * DCMPLX( -P(1), P(2) ) |
| 1242 | FI(2) = -OMEGAP * RS2PA |
| 1243 | & * DCMPLX( SPAZ, RXZERO ) |
| 1244 | FI(3) = OMEGAM * RS2PA / SPAZ |
| 1245 | & * DCMPLX( -P(1), P(2) ) |
| 1246 | FI(4) = OMEGAM * RS2PA |
| 1247 | & * DCMPLX( SPAZ, RXZERO ) |
| 1248 | END IF |
| 1249 | ELSE |
| 1250 | IF ( PLAT .EQ. RXZERO ) THEN |
| 1251 | FI(1) = DCMPLX( OMEGAP, RXZERO ) |
| 1252 | FI(2) = CXZERO |
| 1253 | FI(3) = DCMPLX( -OMEGAM, RXZERO ) |
| 1254 | FI(4) = CXZERO |
| 1255 | ELSE |
| 1256 | RS2PA = RXONE / SQRT( RXTWO * PABS ) |
| 1257 | SPAZ = SQRT( PABS - P(3) ) |
| 1258 | FI(1) = -OMEGAP * RS2PA * SPAZ / PLAT |
| 1259 | & * DCMPLX( -P(1), P(2) ) |
| 1260 | FI(2) = -OMEGAP * RS2PA / SPAZ |
| 1261 | & * DCMPLX( PLAT, RXZERO ) |
| 1262 | FI(3) = OMEGAM * RS2PA * SPAZ / PLAT |
| 1263 | & * DCMPLX( -P(1), P(2) ) |
| 1264 | FI(4) = OMEGAM * RS2PA / SPAZ |
| 1265 | & * DCMPLX( PLAT, RXZERO ) |
| 1266 | END IF |
| 1267 | END IF |
| 1268 | ELSE IF ( NHEL .EQ. -1 ) THEN |
| 1269 | IF ( P(3) .GE. RXZERO ) THEN |
| 1270 | IF ( PLAT .EQ. RXZERO ) THEN |
| 1271 | FI(1) = DCMPLX( OMEGAM, RXZERO ) |
| 1272 | FI(2) = CXZERO |
| 1273 | FI(3) = DCMPLX( -OMEGAP, RXZERO ) |
| 1274 | FI(4) = CXZERO |
| 1275 | ELSE |
| 1276 | RS2PA = RXONE / SQRT( RXTWO * PABS ) |
| 1277 | SPAZ = SQRT( PABS + P(3) ) |
| 1278 | FI(1) = OMEGAM * RS2PA |
| 1279 | & * DCMPLX( SPAZ, RXZERO ) |
| 1280 | FI(2) = OMEGAM * RS2PA / SPAZ |
| 1281 | & * DCMPLX( P(1), P(2) ) |
| 1282 | FI(3) = -OMEGAP * RS2PA |
| 1283 | & * DCMPLX( SPAZ, RXZERO ) |
| 1284 | FI(4) = -OMEGAP * RS2PA / SPAZ |
| 1285 | & * DCMPLX( P(1), P(2) ) |
| 1286 | END IF |
| 1287 | ELSE |
| 1288 | IF ( PLAT .EQ. RXZERO ) THEN |
| 1289 | FI(1) = CXZERO |
| 1290 | FI(2) = DCMPLX( OMEGAM, RXZERO ) |
| 1291 | FI(3) = CXZERO |
| 1292 | FI(4) = DCMPLX( -OMEGAP, RXZERO ) |
| 1293 | ELSE |
| 1294 | RS2PA = RXONE / SQRT( RXTWO * PABS ) |
| 1295 | SPAZ = SQRT( PABS - P(3) ) |
| 1296 | FI(1) = OMEGAM * RS2PA / SPAZ |
| 1297 | & * DCMPLX( PLAT, RXZERO ) |
| 1298 | FI(2) = OMEGAM * RS2PA * SPAZ / PLAT |
| 1299 | & * DCMPLX( P(1), P(2) ) |
| 1300 | FI(3) = -OMEGAP * RS2PA / SPAZ |
| 1301 | & * DCMPLX( PLAT, RXZERO ) |
| 1302 | FI(4) = -OMEGAP * RS2PA * SPAZ / PLAT |
| 1303 | & * DCMPLX( P(1), P(2) ) |
| 1304 | END IF |
| 1305 | END IF |
| 1306 | ELSE |
| 1307 | STOP 'IXXXXX: FERMION HELICITY MUST BE +1,-1' |
| 1308 | END IF |
| 1309 | ELSE |
| 1310 | STOP 'IXXXXX: FERMION TYPE MUST BE +1,-1' |
| 1311 | END IF |
| 1312 | C |
| 1313 | C do massless fermion case |
| 1314 | C |
| 1315 | ELSE |
| 1316 | IF ( NSF .EQ. 1 ) THEN |
| 1317 | IF ( NHEL .EQ. 1 ) THEN |
| 1318 | IF ( P(3) .GE. RXZERO ) THEN |
| 1319 | IF ( PLAT .EQ. RXZERO ) THEN |
| 1320 | FI(1) = CXZERO |
| 1321 | FI(2) = CXZERO |
| 1322 | FI(3) = DCMPLX( OMEGAP, RXZERO ) |
| 1323 | FI(4) = CXZERO |
| 1324 | ELSE |
| 1325 | SPAZ = SQRT( PABS + P(3) ) |
| 1326 | FI(1) = CXZERO |
| 1327 | FI(2) = CXZERO |
| 1328 | FI(3) = DCMPLX( SPAZ, RXZERO ) |
| 1329 | FI(4) = RXONE / SPAZ |
| 1330 | & * DCMPLX( P(1), P(2) ) |
| 1331 | END IF |
| 1332 | ELSE |
| 1333 | IF ( PLAT .EQ. RXZERO ) THEN |
| 1334 | FI(1) = CXZERO |
| 1335 | FI(2) = CXZERO |
| 1336 | FI(3) = CXZERO |
| 1337 | FI(4) = DCMPLX( OMEGAP, RXZERO ) |
| 1338 | ELSE |
| 1339 | SPAZ = SQRT( PABS - P(3) ) |
| 1340 | FI(1) = CXZERO |
| 1341 | FI(2) = CXZERO |
| 1342 | FI(3) = RXONE / SPAZ |
| 1343 | & * DCMPLX( PLAT, RXZERO ) |
| 1344 | FI(4) = SPAZ / PLAT |
| 1345 | & * DCMPLX( P(1), P(2) ) |
| 1346 | END IF |
| 1347 | END IF |
| 1348 | ELSE IF ( NHEL .EQ. -1 ) THEN |
| 1349 | IF ( P(3) .GE. RXZERO ) THEN |
| 1350 | IF ( PLAT .EQ. RXZERO ) THEN |
| 1351 | FI(1) = CXZERO |
| 1352 | FI(2) = DCMPLX( OMEGAP, RXZERO ) |
| 1353 | FI(3) = CXZERO |
| 1354 | FI(4) = CXZERO |
| 1355 | ELSE |
| 1356 | SPAZ = SQRT( PABS + P(3) ) |
| 1357 | FI(1) = RXONE / SPAZ |
| 1358 | & * DCMPLX( -P(1), P(2) ) |
| 1359 | FI(2) = DCMPLX( SPAZ, RXZERO ) |
| 1360 | FI(3) = CXZERO |
| 1361 | FI(4) = CXZERO |
| 1362 | END IF |
| 1363 | ELSE |
| 1364 | IF ( PLAT .EQ. RXZERO ) THEN |
| 1365 | FI(1) = DCMPLX( -OMEGAP, RXZERO ) |
| 1366 | FI(2) = CXZERO |
| 1367 | FI(3) = CXZERO |
| 1368 | FI(4) = CXZERO |
| 1369 | ELSE |
| 1370 | SPAZ = SQRT( PABS - P(3) ) |
| 1371 | FI(1) = SPAZ / PLAT |
| 1372 | & * DCMPLX( -P(1), P(2) ) |
| 1373 | FI(2) = RXONE / SPAZ |
| 1374 | & * DCMPLX( PLAT, RXZERO ) |
| 1375 | FI(3) = CXZERO |
| 1376 | FI(4) = CXZERO |
| 1377 | END IF |
| 1378 | END IF |
| 1379 | ELSE |
| 1380 | STOP 'IXXXXX: FERMION HELICITY MUST BE +1,-1' |
| 1381 | END IF |
| 1382 | ELSE IF ( NSF .EQ. -1 ) THEN |
| 1383 | IF ( NHEL .EQ. 1 ) THEN |
| 1384 | IF ( P(3) .GE. RXZERO ) THEN |
| 1385 | IF ( PLAT .EQ. RXZERO ) THEN |
| 1386 | FI(1) = CXZERO |
| 1387 | FI(2) = DCMPLX( -OMEGAP, RXZERO ) |
| 1388 | FI(3) = CXZERO |
| 1389 | FI(4) = CXZERO |
| 1390 | ELSE |
| 1391 | SPAZ = SQRT( PABS + P(3) ) |
| 1392 | FI(1) = -RXONE / SPAZ |
| 1393 | & * DCMPLX( -P(1), P(2) ) |
| 1394 | FI(2) = DCMPLX( -SPAZ, RXZERO ) |
| 1395 | FI(3) = CXZERO |
| 1396 | FI(4) = CXZERO |
| 1397 | END IF |
| 1398 | ELSE |
| 1399 | IF ( PLAT .EQ. RXZERO ) THEN |
| 1400 | FI(1) = DCMPLX( OMEGAP, RXZERO ) |
| 1401 | FI(2) = CXZERO |
| 1402 | FI(3) = CXZERO |
| 1403 | FI(4) = CXZERO |
| 1404 | ELSE |
| 1405 | SPAZ = SQRT( PABS - P(3) ) |
| 1406 | FI(1) = -SPAZ / PLAT |
| 1407 | & * DCMPLX( -P(1), P(2) ) |
| 1408 | FI(2) = -RXONE / SPAZ |
| 1409 | & * DCMPLX( PLAT, RXZERO ) |
| 1410 | FI(3) = CXZERO |
| 1411 | FI(4) = CXZERO |
| 1412 | END IF |
| 1413 | END IF |
| 1414 | ELSE IF ( NHEL .EQ. -1 ) THEN |
| 1415 | IF ( P(3) .GE. RXZERO ) THEN |
| 1416 | IF ( PLAT .EQ. RXZERO ) THEN |
| 1417 | FI(1) = CXZERO |
| 1418 | FI(2) = CXZERO |
| 1419 | FI(3) = DCMPLX( -OMEGAP, RXZERO ) |
| 1420 | FI(4) = CXZERO |
| 1421 | ELSE |
| 1422 | SPAZ = SQRT( PABS + P(3) ) |
| 1423 | FI(1) = CXZERO |
| 1424 | FI(2) = CXZERO |
| 1425 | FI(3) = DCMPLX( -SPAZ, RXZERO ) |
| 1426 | FI(4) = -RXONE / SPAZ |
| 1427 | & * DCMPLX( P(1), P(2) ) |
| 1428 | END IF |
| 1429 | ELSE |
| 1430 | IF ( PLAT .EQ. RXZERO ) THEN |
| 1431 | FI(1) = CXZERO |
| 1432 | FI(2) = CXZERO |
| 1433 | FI(3) = CXZERO |
| 1434 | FI(4) = DCMPLX( -OMEGAP, RXZERO ) |
| 1435 | ELSE |
| 1436 | SPAZ = SQRT( PABS - P(3) ) |
| 1437 | FI(1) = CXZERO |
| 1438 | FI(2) = CXZERO |
| 1439 | FI(3) = -RXONE / SPAZ |
| 1440 | & * DCMPLX( PLAT, RXZERO ) |
| 1441 | FI(4) = -SPAZ / PLAT |
| 1442 | & * DCMPLX( P(1), P(2) ) |
| 1443 | END IF |
| 1444 | END IF |
| 1445 | ELSE |
| 1446 | STOP 'IXXXXX: FERMION HELICITY MUST BE +1,-1' |
| 1447 | END IF |
| 1448 | ELSE |
| 1449 | STOP 'IXXXXX: FERMION TYPE MUST BE +1,-1' |
| 1450 | END IF |
| 1451 | END IF |
| 1452 | C |
| 1453 | C done |
| 1454 | C |
| 1455 | RETURN |
| 1456 | END |
| 1457 | C |
| 1458 | C ---------------------------------------------------------------------- |
| 1459 | C |
| 1460 | SUBROUTINE J3XXXX(FI,FO,GAF,GZF,ZMASS,ZWIDTH , J3) |
| 1461 | C |
| 1462 | C this subroutine computes the sum of photon and z currents with the |
| 1463 | C suitable weights ( j(w3) = cos(theta_w) j(z) + sin(theta_w) j(a) ). |
| 1464 | C the output j3 is useful as an input of vvvxxx, jvvxxx or w3w3xx. |
| 1465 | C the photon propagator is given in feynman gauge, and the z propagator |
| 1466 | C is given in unitary gauge. |
| 1467 | C |
| 1468 | C input: |
| 1469 | C complex fi(6) : flow-in fermion |fi> |
| 1470 | C complex fo(6) : flow-out fermion <fo| |
| 1471 | C real gaf(2) : fi couplings with a gaf |
| 1472 | C real gzf(2) : fi couplings with z gzf |
| 1473 | C real zmass : mass of z |
| 1474 | C real zwidth : width of z |
| 1475 | C |
| 1476 | C output: |
| 1477 | C complex j3(6) : w3 current j^mu(<fo|w3|fi>) |
| 1478 | C |
| 1479 | #if defined(CERNLIB_IMPNONE) |
| 1480 | IMPLICIT NONE |
| 1481 | #endif |
| 1482 | COMPLEX*16 FI(6),FO(6),J3(6), |
| 1483 | & C0L,C1L,C2L,C3L,CSL,C0R,C1R,C2R,C3R,CSR,DZ,DDIF |
| 1484 | REAL*8 GAF(2),GZF(2),Q(0:3),ZMASS,ZWIDTH,ZM2,ZMW,Q2,DA,WW, |
| 1485 | & CW,SW,GN,GZ3L,GA3L |
| 1486 | C |
| 1487 | REAL*8 RXZERO, RXONE |
| 1488 | PARAMETER( RXZERO=0.0D0, RXONE=1.0D0 ) |
| 1489 | COMPLEX*16 CXIMAG |
| 1490 | LOGICAL FIRST |
| 1491 | SAVE CXIMAG,FIRST |
| 1492 | DATA FIRST/.TRUE./ |
| 1493 | C |
| 1494 | C Fix compilation with g77 |
| 1495 | IF(FIRST) THEN |
| 1496 | FIRST=.FALSE. |
| 1497 | CXIMAG=DCMPLX( RXZERO, RXONE ) |
| 1498 | ENDIF |
| 1499 | C |
| 1500 | J3(5) = FO(5)-FI(5) |
| 1501 | J3(6) = FO(6)-FI(6) |
| 1502 | C |
| 1503 | Q(0)=-DBLE( J3(5)) |
| 1504 | Q(1)=-DBLE( J3(6)) |
| 1505 | Q(2)=-DIMAG(J3(6)) |
| 1506 | Q(3)=-DIMAG(J3(5)) |
| 1507 | Q2=Q(0)**2-(Q(1)**2+Q(2)**2+Q(3)**2) |
| 1508 | ZM2=ZMASS**2 |
| 1509 | ZMW=ZMASS*ZWIDTH |
| 1510 | C |
| 1511 | DA=RXONE/Q2 |
| 1512 | WW=MAX(DSIGN( ZMW ,Q2),RXZERO) |
| 1513 | DZ=RXONE/DCMPLX( Q2-ZM2 , WW ) |
| 1514 | DDIF=DCMPLX( -ZM2 , WW )*DA*DZ |
| 1515 | C |
| 1516 | C ddif is the difference : ddif=da-dz |
| 1517 | C for the running width, use below instead of the above ww,dz and ddif. |
| 1518 | C ww=max( zwidth*q2/zmass ,r_zero) |
| 1519 | C dz=r_one/dcmplx( q2-zm2 , ww ) |
| 1520 | C ddif=dcmplx( -zm2 , ww )*da*dz |
| 1521 | C |
| 1522 | CW=RXONE/SQRT(RXONE+(GZF(2)/GAF(2))**2) |
| 1523 | SW=SQRT((RXONE-CW)*(RXONE+CW)) |
| 1524 | GN=GAF(2)*SW |
| 1525 | GZ3L=GZF(1)*CW |
| 1526 | GA3L=GAF(1)*SW |
| 1527 | C0L= FO(3)*FI(1)+FO(4)*FI(2) |
| 1528 | C0R= FO(1)*FI(3)+FO(2)*FI(4) |
| 1529 | C1L=-(FO(3)*FI(2)+FO(4)*FI(1)) |
| 1530 | C1R= FO(1)*FI(4)+FO(2)*FI(3) |
| 1531 | C2L= (FO(3)*FI(2)-FO(4)*FI(1))*CXIMAG |
| 1532 | C2R=(-FO(1)*FI(4)+FO(2)*FI(3))*CXIMAG |
| 1533 | C3L= -FO(3)*FI(1)+FO(4)*FI(2) |
| 1534 | C3R= FO(1)*FI(3)-FO(2)*FI(4) |
| 1535 | CSL=(Q(0)*C0L-Q(1)*C1L-Q(2)*C2L-Q(3)*C3L)/ZM2 |
| 1536 | CSR=(Q(0)*C0R-Q(1)*C1R-Q(2)*C2R-Q(3)*C3R)/ZM2 |
| 1537 | C |
| 1538 | J3(1) = GZ3L*DZ*(C0L-CSL*Q(0))+GA3L*C0L*DA |
| 1539 | & + GN*(C0R*DDIF-CSR*Q(0)*DZ) |
| 1540 | J3(2) = GZ3L*DZ*(C1L-CSL*Q(1))+GA3L*C1L*DA |
| 1541 | & + GN*(C1R*DDIF-CSR*Q(1)*DZ) |
| 1542 | J3(3) = GZ3L*DZ*(C2L-CSL*Q(2))+GA3L*C2L*DA |
| 1543 | & + GN*(C2R*DDIF-CSR*Q(2)*DZ) |
| 1544 | J3(4) = GZ3L*DZ*(C3L-CSL*Q(3))+GA3L*C3L*DA |
| 1545 | & + GN*(C3R*DDIF-CSR*Q(3)*DZ) |
| 1546 | C |
| 1547 | RETURN |
| 1548 | END |
| 1549 | C |
| 1550 | C ---------------------------------------------------------------------- |
| 1551 | C |
| 1552 | SUBROUTINE JEEXXX(EB,EF,SHLF,CHLF,PHI,NHB,NHF,NSF , JEE) |
| 1553 | C |
| 1554 | C This subroutine computes an off-shell photon wavefunction emitted from |
| 1555 | C the electron or positron beam, with a special care for the small angle |
| 1556 | C region. The momenta are measured in the laboratory frame, where the |
| 1557 | C e- (e+) beam is along the positive (negative) z axis. |
| 1558 | C |
| 1559 | C INPUT: |
| 1560 | C real EB : energy (GeV) of beam e-/e+ |
| 1561 | C real EF : energy (GeV) of final e-/e+ |
| 1562 | C real SHLF : sin(theta/2) of final e-/e+ |
| 1563 | C real CHLF : cos(theta/2) of final e-/e+ |
| 1564 | C real PHI : azimuthal angle of final e-/e+ |
| 1565 | C integer NHB = -1 or 1 : helicity of beam e-/e+ |
| 1566 | C integer NHF = -1 or 1 : helicity of final e-/e+ |
| 1567 | C integer NSF = -1 or 1 : +1 for electron, -1 for positron |
| 1568 | C |
| 1569 | C OUTPUT: |
| 1570 | C complex JEE(6) : off-shell photon J^mu(<e|A|e>) |
| 1571 | C |
| 1572 | #if defined(CERNLIB_IMPNONE) |
| 1573 | IMPLICIT NONE |
| 1574 | #endif |
| 1575 | COMPLEX*16 JEE(6),COEFF |
| 1576 | REAL*8 CS(2),EB,EF,SHLF,CHLF,PHI,ME,ALPHA,GAL,HI,SF,SFH,X,ME2,Q2, |
| 1577 | & RFP,RFM,SNP,CSP,RXC,C,S |
| 1578 | INTEGER NHB,NHF,NSF |
| 1579 | C |
| 1580 | ME =0.51099906D-3 |
| 1581 | ALPHA=1./128. |
| 1582 | GAL =SQRT(ALPHA*4.*3.14159265D0) |
| 1583 | C |
| 1584 | HI =NHB |
| 1585 | SF =NSF |
| 1586 | SFH=NHB*NSF |
| 1587 | CS((3+NSF)/2)=SHLF |
| 1588 | CS((3-NSF)/2)=CHLF |
| 1589 | C CS(1)=CHLF and CS(2)=SHLF for electron |
| 1590 | C CS(1)=SHLF and CS(2)=CHLF for positron |
| 1591 | X=EF/EB |
| 1592 | ME2=ME**2 |
| 1593 | Q2=-4.*CS(2)**2*(EF*EB-ME2) |
| 1594 | & +SF*(1.-X)**2/X*(SHLF+CHLF)*(SHLF-CHLF)*ME2 |
| 1595 | RFP=(1+NSF) |
| 1596 | RFM=(1-NSF) |
| 1597 | SNP=SIN(PHI) |
| 1598 | CSP=COS(PHI) |
| 1599 | C |
| 1600 | IF (NHB.EQ.NHF) THEN |
| 1601 | RXC=2.*X/(1.-X)*CS(1)**2 |
| 1602 | COEFF= GAL*2.*EB*SQRT(X)*CS(2)/Q2 |
| 1603 | & *(DCMPLX( RFP )-RFM*DCMPLX( CSP ,-SNP*HI ))*.5 |
| 1604 | JEE(1) = DCMPLX( 0.D0 ) |
| 1605 | JEE(2) = COEFF*DCMPLX( (1.+RXC)*CSP ,-SFH*SNP ) |
| 1606 | JEE(3) = COEFF*DCMPLX( (1.+RXC)*SNP , SFH*CSP ) |
| 1607 | JEE(4) = COEFF*(-SF*RXC/CS(1)*CS(2)) |
| 1608 | ELSE |
| 1609 | COEFF= GAL*ME/Q2/SQRT(X) |
| 1610 | & *(DCMPLX( RFP )+RFM*DCMPLX( CSP , SNP*HI ))*.5*HI |
| 1611 | JEE(1) = -COEFF*(1.+X)*CS(2)*DCMPLX( CSP , SFH*SNP ) |
| 1612 | JEE(2) = COEFF*(1.-X)*CS(1) |
| 1613 | JEE(3) = JEE(2)*DCMPLX( 0.D0 , SFH ) |
| 1614 | JEE(4) = JEE(1)*SF*(1.-X)/(1.+X) |
| 1615 | ENDIF |
| 1616 | C |
| 1617 | C=(CHLF+SHLF)*(CHLF-SHLF) |
| 1618 | S=2.*CHLF*SHLF |
| 1619 | C |
| 1620 | JEE(5) = -EB*DCMPLX( 1.-X , SF-X*C ) |
| 1621 | JEE(6) = EB*X*S*DCMPLX( CSP , SNP ) |
| 1622 | C |
| 1623 | RETURN |
| 1624 | END |
| 1625 | C |
| 1626 | C |
| 1627 | C ---------------------------------------------------------------------- |
| 1628 | C |
| 1629 | SUBROUTINE JGGGXX(W1,W2,W3,G, JW3W) |
| 1630 | C |
| 1631 | C this subroutine computes an off-shell w+, w-, w3, z or photon current |
| 1632 | C from the four-point gauge boson coupling, including the contributions |
| 1633 | C of w exchange diagrams. the vector propagator is given in feynman |
| 1634 | C gauge for a photon and in unitary gauge for w and z bosons. if one |
| 1635 | C sets wmass=0.0, then the ggg-->g current is given (see sect 2.9.1 of |
| 1636 | C the manual). |
| 1637 | C |
| 1638 | C input: |
| 1639 | C complex w1(6) : first vector w1 |
| 1640 | C complex w2(6) : second vector w2 |
| 1641 | C complex w3(6) : third vector w3 |
| 1642 | C real g : first coupling constant |
| 1643 | C (see the table below) |
| 1644 | C |
| 1645 | C output: |
| 1646 | C complex jw3w(6) : w current j^mu(w':w1,w2,w3) |
| 1647 | C |
| 1648 | #if defined(CERNLIB_IMPNONE) |
| 1649 | IMPLICIT NONE |
| 1650 | #endif |
| 1651 | COMPLEX*16 W1(6),W2(6),W3(6),JW3W(6) |
| 1652 | COMPLEX*16 DW1(0:3),DW2(0:3),DW3(0:3), |
| 1653 | & JJ(0:3),DV,W32,W13 |
| 1654 | REAL*8 P1(0:3),P2(0:3),P3(0:3),Q(0:3),G,DG2,Q2 |
| 1655 | C |
| 1656 | REAL*8 RXZERO |
| 1657 | PARAMETER( RXZERO=0.0D0 ) |
| 1658 | C |
| 1659 | JW3W(5) = W1(5)+W2(5)+W3(5) |
| 1660 | JW3W(6) = W1(6)+W2(6)+W3(6) |
| 1661 | C |
| 1662 | DW1(0)=DCMPLX(W1(1)) |
| 1663 | DW1(1)=DCMPLX(W1(2)) |
| 1664 | DW1(2)=DCMPLX(W1(3)) |
| 1665 | DW1(3)=DCMPLX(W1(4)) |
| 1666 | DW2(0)=DCMPLX(W2(1)) |
| 1667 | DW2(1)=DCMPLX(W2(2)) |
| 1668 | DW2(2)=DCMPLX(W2(3)) |
| 1669 | DW2(3)=DCMPLX(W2(4)) |
| 1670 | DW3(0)=DCMPLX(W3(1)) |
| 1671 | DW3(1)=DCMPLX(W3(2)) |
| 1672 | DW3(2)=DCMPLX(W3(3)) |
| 1673 | DW3(3)=DCMPLX(W3(4)) |
| 1674 | P1(0)=DBLE( W1(5)) |
| 1675 | P1(1)=DBLE( W1(6)) |
| 1676 | P1(2)=DBLE(DIMAG(W1(6))) |
| 1677 | P1(3)=DBLE(DIMAG(W1(5))) |
| 1678 | P2(0)=DBLE( W2(5)) |
| 1679 | P2(1)=DBLE( W2(6)) |
| 1680 | P2(2)=DBLE(DIMAG(W2(6))) |
| 1681 | P2(3)=DBLE(DIMAG(W2(5))) |
| 1682 | P3(0)=DBLE( W3(5)) |
| 1683 | P3(1)=DBLE( W3(6)) |
| 1684 | P3(2)=DBLE(DIMAG(W3(6))) |
| 1685 | P3(3)=DBLE(DIMAG(W3(5))) |
| 1686 | Q(0)=-(P1(0)+P2(0)+P3(0)) |
| 1687 | Q(1)=-(P1(1)+P2(1)+P3(1)) |
| 1688 | Q(2)=-(P1(2)+P2(2)+P3(2)) |
| 1689 | Q(3)=-(P1(3)+P2(3)+P3(3)) |
| 1690 | |
| 1691 | Q2 =Q(0)**2 -(Q(1)**2 +Q(2)**2 +Q(3)**2) |
| 1692 | |
| 1693 | DG2=DBLE(G)*DBLE(G) |
| 1694 | C |
| 1695 | DV = 1.0D0/DCMPLX( Q2 ) |
| 1696 | |
| 1697 | C for the running width, use below instead of the above dv. |
| 1698 | C dv = 1.0d0/dcmplx( q2 -mv2 , dmax1(dwv*q2/dmv,0.d0) ) |
| 1699 | C |
| 1700 | W32=DW3(0)*DW2(0)-DW3(1)*DW2(1)-DW3(2)*DW2(2)-DW3(3)*DW2(3) |
| 1701 | C |
| 1702 | C |
| 1703 | W13=DW1(0)*DW3(0)-DW1(1)*DW3(1)-DW1(2)*DW3(2)-DW1(3)*DW3(3) |
| 1704 | C |
| 1705 | JJ(0)=DG2*( DW1(0)*W32 - DW2(0)*W13 ) |
| 1706 | JJ(1)=DG2*( DW1(1)*W32 - DW2(1)*W13 ) |
| 1707 | JJ(2)=DG2*( DW1(2)*W32 - DW2(2)*W13 ) |
| 1708 | JJ(3)=DG2*( DW1(3)*W32 - DW2(3)*W13 ) |
| 1709 | C |
| 1710 | JW3W(1) = DCMPLX( JJ(0)*DV ) |
| 1711 | JW3W(2) = DCMPLX( JJ(1)*DV ) |
| 1712 | JW3W(3) = DCMPLX( JJ(2)*DV ) |
| 1713 | JW3W(4) = DCMPLX( JJ(3)*DV ) |
| 1714 | C |
| 1715 | RETURN |
| 1716 | END |
| 1717 | C |
| 1718 | C ---------------------------------------------------------------------- |
| 1719 | C |
| 1720 | SUBROUTINE JGGXXX(V1,V2,G, JVV) |
| 1721 | C |
| 1722 | C this subroutine computes an off-shell vector current from the three- |
| 1723 | C point gauge boson coupling. the vector propagator is given in feynman |
| 1724 | C gauge for a massless vector and in unitary gauge for a massive vector. |
| 1725 | C |
| 1726 | C input: |
| 1727 | C complex v1(6) : first vector v1 |
| 1728 | C complex v2(6) : second vector v2 |
| 1729 | C real g : coupling constant (see the table below) |
| 1730 | C |
| 1731 | C output: |
| 1732 | C complex jvv(6) : vector current j^mu(v:v1,v2) |
| 1733 | C |
| 1734 | #if defined(CERNLIB_IMPNONE) |
| 1735 | IMPLICIT NONE |
| 1736 | #endif |
| 1737 | COMPLEX*16 V1(6),V2(6),JVV(6),J12(0:3), |
| 1738 | & SV1,SV2,V12 |
| 1739 | REAL*8 P1(0:3),P2(0:3),Q(0:3),G,GS,S |
| 1740 | C |
| 1741 | REAL*8 RXZERO |
| 1742 | PARAMETER( RXZERO=0.0D0 ) |
| 1743 | C |
| 1744 | JVV(5) = V1(5)+V2(5) |
| 1745 | JVV(6) = V1(6)+V2(6) |
| 1746 | C |
| 1747 | P1(0)=DBLE( V1(5)) |
| 1748 | P1(1)=DBLE( V1(6)) |
| 1749 | P1(2)=DIMAG(V1(6)) |
| 1750 | P1(3)=DIMAG(V1(5)) |
| 1751 | P2(0)=DBLE( V2(5)) |
| 1752 | P2(1)=DBLE( V2(6)) |
| 1753 | P2(2)=DIMAG(V2(6)) |
| 1754 | P2(3)=DIMAG(V2(5)) |
| 1755 | Q(0)=-DBLE( JVV(5)) |
| 1756 | Q(1)=-DBLE( JVV(6)) |
| 1757 | Q(2)=-DIMAG(JVV(6)) |
| 1758 | Q(3)=-DIMAG(JVV(5)) |
| 1759 | S=Q(0)**2-(Q(1)**2+Q(2)**2+Q(3)**2) |
| 1760 | C |
| 1761 | V12=V1(1)*V2(1)-V1(2)*V2(2)-V1(3)*V2(3)-V1(4)*V2(4) |
| 1762 | SV1= (P2(0)-Q(0))*V1(1) -(P2(1)-Q(1))*V1(2) |
| 1763 | & -(P2(2)-Q(2))*V1(3) -(P2(3)-Q(3))*V1(4) |
| 1764 | SV2=-(P1(0)-Q(0))*V2(1) +(P1(1)-Q(1))*V2(2) |
| 1765 | & +(P1(2)-Q(2))*V2(3) +(P1(3)-Q(3))*V2(4) |
| 1766 | J12(0)=(P1(0)-P2(0))*V12 +SV1*V2(1) +SV2*V1(1) |
| 1767 | J12(1)=(P1(1)-P2(1))*V12 +SV1*V2(2) +SV2*V1(2) |
| 1768 | J12(2)=(P1(2)-P2(2))*V12 +SV1*V2(3) +SV2*V1(3) |
| 1769 | J12(3)=(P1(3)-P2(3))*V12 +SV1*V2(4) +SV2*V1(4) |
| 1770 | C |
| 1771 | GS=-G/S |
| 1772 | C |
| 1773 | JVV(1) = GS*J12(0) |
| 1774 | JVV(2) = GS*J12(1) |
| 1775 | JVV(3) = GS*J12(2) |
| 1776 | JVV(4) = GS*J12(3) |
| 1777 | C |
| 1778 | RETURN |
| 1779 | END |
| 1780 | C |
| 1781 | C ---------------------------------------------------------------------- |
| 1782 | C |
| 1783 | SUBROUTINE JIOXXX(FI,FO,G,VMASS,VWIDTH , JIO) |
| 1784 | C |
| 1785 | C this subroutine computes an off-shell vector current from an external |
| 1786 | C fermion pair. the vector boson propagator is given in feynman gauge |
| 1787 | C for a massless vector and in unitary gauge for a massive vector. |
| 1788 | C |
| 1789 | C input: |
| 1790 | C complex fi(6) : flow-in fermion |fi> |
| 1791 | C complex fo(6) : flow-out fermion <fo| |
| 1792 | C real g(2) : coupling constants gvf |
| 1793 | C real vmass : mass of output vector v |
| 1794 | C real vwidth : width of output vector v |
| 1795 | C |
| 1796 | C output: |
| 1797 | C complex jio(6) : vector current j^mu(<fo|v|fi>) |
| 1798 | C |
| 1799 | #if defined(CERNLIB_IMPNONE) |
| 1800 | IMPLICIT NONE |
| 1801 | #endif |
| 1802 | COMPLEX*16 FI(6),FO(6),JIO(6),C0,C1,C2,C3,CS,D |
| 1803 | REAL*8 G(2),Q(0:3),VMASS,VWIDTH,Q2,VM2,DD |
| 1804 | C |
| 1805 | REAL*8 RXZERO, RXONE |
| 1806 | PARAMETER( RXZERO=0.0D0, RXONE=1.0D0 ) |
| 1807 | COMPLEX*16 CXIMAG |
| 1808 | LOGICAL FIRST |
| 1809 | SAVE CXIMAG,FIRST |
| 1810 | DATA FIRST/.TRUE./ |
| 1811 | C |
| 1812 | C Fix compilation with g77 |
| 1813 | IF(FIRST) THEN |
| 1814 | FIRST=.FALSE. |
| 1815 | CXIMAG=DCMPLX( RXZERO, RXONE ) |
| 1816 | ENDIF |
| 1817 | C |
| 1818 | JIO(5) = FO(5)-FI(5) |
| 1819 | JIO(6) = FO(6)-FI(6) |
| 1820 | C |
| 1821 | Q(0)=DBLE( JIO(5)) |
| 1822 | Q(1)=DBLE( JIO(6)) |
| 1823 | Q(2)=DIMAG(JIO(6)) |
| 1824 | Q(3)=DIMAG(JIO(5)) |
| 1825 | Q2=Q(0)**2-(Q(1)**2+Q(2)**2+Q(3)**2) |
| 1826 | VM2=VMASS**2 |
| 1827 | C |
| 1828 | IF (VMASS.NE.RXZERO) THEN |
| 1829 | C |
| 1830 | D=RXONE/DCMPLX( Q2-VM2 , MAX(SIGN( VMASS*VWIDTH ,Q2),RXZERO) ) |
| 1831 | C for the running width, use below instead of the above d. |
| 1832 | C d=r_one/dcmplx( q2-vm2 , max( vwidth*q2/vmass ,r_zero) ) |
| 1833 | C |
| 1834 | IF (G(2).NE.RXZERO) THEN |
| 1835 | C |
| 1836 | C0= G(1)*( FO(3)*FI(1)+FO(4)*FI(2)) |
| 1837 | & +G(2)*( FO(1)*FI(3)+FO(2)*FI(4)) |
| 1838 | C1= -G(1)*( FO(3)*FI(2)+FO(4)*FI(1)) |
| 1839 | & +G(2)*( FO(1)*FI(4)+FO(2)*FI(3)) |
| 1840 | C2=( G(1)*( FO(3)*FI(2)-FO(4)*FI(1)) |
| 1841 | & +G(2)*(-FO(1)*FI(4)+FO(2)*FI(3)))*CXIMAG |
| 1842 | C3= G(1)*(-FO(3)*FI(1)+FO(4)*FI(2)) |
| 1843 | & +G(2)*( FO(1)*FI(3)-FO(2)*FI(4)) |
| 1844 | ELSE |
| 1845 | C |
| 1846 | D=D*G(1) |
| 1847 | C0= FO(3)*FI(1)+FO(4)*FI(2) |
| 1848 | C1= -FO(3)*FI(2)-FO(4)*FI(1) |
| 1849 | C2=( FO(3)*FI(2)-FO(4)*FI(1))*CXIMAG |
| 1850 | C3= -FO(3)*FI(1)+FO(4)*FI(2) |
| 1851 | END IF |
| 1852 | C |
| 1853 | CS=(Q(0)*C0-Q(1)*C1-Q(2)*C2-Q(3)*C3)/VM2 |
| 1854 | C |
| 1855 | JIO(1) = (C0-CS*Q(0))*D |
| 1856 | JIO(2) = (C1-CS*Q(1))*D |
| 1857 | JIO(3) = (C2-CS*Q(2))*D |
| 1858 | JIO(4) = (C3-CS*Q(3))*D |
| 1859 | C |
| 1860 | ELSE |
| 1861 | DD=RXONE/Q2 |
| 1862 | C |
| 1863 | IF (G(2).NE.RXZERO) THEN |
| 1864 | JIO(1) = ( G(1)*( FO(3)*FI(1)+FO(4)*FI(2)) |
| 1865 | & +G(2)*( FO(1)*FI(3)+FO(2)*FI(4)) )*DD |
| 1866 | JIO(2) = (-G(1)*( FO(3)*FI(2)+FO(4)*FI(1)) |
| 1867 | & +G(2)*( FO(1)*FI(4)+FO(2)*FI(3)) )*DD |
| 1868 | JIO(3) = ( G(1)*( FO(3)*FI(2)-FO(4)*FI(1)) |
| 1869 | & +G(2)*(-FO(1)*FI(4)+FO(2)*FI(3))) |
| 1870 | $ *DCMPLX(RXZERO,DD) |
| 1871 | JIO(4) = ( G(1)*(-FO(3)*FI(1)+FO(4)*FI(2)) |
| 1872 | & +G(2)*( FO(1)*FI(3)-FO(2)*FI(4)) )*DD |
| 1873 | C |
| 1874 | ELSE |
| 1875 | DD=DD*G(1) |
| 1876 | C |
| 1877 | JIO(1) = ( FO(3)*FI(1)+FO(4)*FI(2))*DD |
| 1878 | JIO(2) = -( FO(3)*FI(2)+FO(4)*FI(1))*DD |
| 1879 | JIO(3) = ( FO(3)*FI(2)-FO(4)*FI(1))*DCMPLX(RXZERO,DD) |
| 1880 | JIO(4) = (-FO(3)*FI(1)+FO(4)*FI(2))*DD |
| 1881 | END IF |
| 1882 | END IF |
| 1883 | C |
| 1884 | RETURN |
| 1885 | END |
| 1886 | C ---------------------------------------------------------------------- |
| 1887 | C |
| 1888 | SUBROUTINE JSSXXX(S1,S2,G,VMASS,VWIDTH , JSS) |
| 1889 | C |
| 1890 | C This subroutine computes an off-shell vector current from the vector- |
| 1891 | C scalar-scalar coupling. The coupling is absent in the minimal SM in |
| 1892 | C unitary gauge. The propagator is given in Feynman gauge for a |
| 1893 | C massless vector and in unitary gauge for a massive vector. |
| 1894 | C |
| 1895 | C INPUT: |
| 1896 | C complex S1(3) : first scalar S1 |
| 1897 | C complex S2(3) : second scalar S2 |
| 1898 | C real G : coupling constant (S1 charge) |
| 1899 | C real VMASS : mass of OUTPUT vector V |
| 1900 | C real VWIDTH : width of OUTPUT vector V |
| 1901 | C |
| 1902 | C Examples of the coupling constant G for SUSY particles are as follows: |
| 1903 | C ----------------------------------------------------------- |
| 1904 | C | S1 | (Q,I3) of S1 || V=A | V=Z | V=W | |
| 1905 | C ----------------------------------------------------------- |
| 1906 | C | nu~_L | ( 0 , +1/2) || --- | GZN(1) | GWF(1) | |
| 1907 | C | e~_L | ( -1 , -1/2) || GAL(1) | GZL(1) | GWF(1) | |
| 1908 | C | u~_L | (+2/3 , +1/2) || GAU(1) | GZU(1) | GWF(1) | |
| 1909 | C | d~_L | (-1/3 , -1/2) || GAD(1) | GZD(1) | GWF(1) | |
| 1910 | C ----------------------------------------------------------- |
| 1911 | C | e~_R-bar | ( +1 , 0 ) || -GAL(2) | -GZL(2) | -GWF(2) | |
| 1912 | C | u~_R-bar | (-2/3 , 0 ) || -GAU(2) | -GZU(2) | -GWF(2) | |
| 1913 | C | d~_R-bar | (+1/3 , 0 ) || -GAD(2) | -GZD(2) | -GWF(2) | |
| 1914 | C ----------------------------------------------------------- |
| 1915 | C where the S1 charge is defined by the flowing-OUT quantum number. |
| 1916 | C |
| 1917 | C OUTPUT: |
| 1918 | C complex JSS(6) : vector current J^mu(V:S1,S2) |
| 1919 | C |
| 1920 | #if defined(CERNLIB_IMPNONE) |
| 1921 | IMPLICIT NONE |
| 1922 | #endif |
| 1923 | COMPLEX*16 S1(3),S2(3),JSS(6),DG,ADG |
| 1924 | REAL*8 PP(0:3),PA(0:3),Q(0:3),G,VMASS,VWIDTH,Q2,VM2,MP2,MA2,M2D |
| 1925 | C |
| 1926 | JSS(5) = S1(2)+S2(2) |
| 1927 | JSS(6) = S1(3)+S2(3) |
| 1928 | C |
| 1929 | Q(0)=DBLE( JSS(5)) |
| 1930 | Q(1)=DBLE( JSS(6)) |
| 1931 | Q(2)=DIMAG(JSS(6)) |
| 1932 | Q(3)=DIMAG(JSS(5)) |
| 1933 | Q2=Q(0)**2-(Q(1)**2+Q(2)**2+Q(3)**2) |
| 1934 | VM2=VMASS**2 |
| 1935 | C |
| 1936 | IF (VMASS.EQ.0.) GOTO 10 |
| 1937 | C |
| 1938 | DG=G/DCMPLX( Q2-VM2, MAX(SIGN( VMASS*VWIDTH ,Q2),0.D0)) |
| 1939 | C For the running width, use below instead of the above DG. |
| 1940 | C DG=G/dCMPLX( Q2-VM2 , MAX( VWIDTH*Q2/VMASS ,0.) ) |
| 1941 | C |
| 1942 | ADG=DG*S1(1)*S2(1) |
| 1943 | C |
| 1944 | PP(0)=DBLE( S1(2)) |
| 1945 | PP(1)=DBLE( S1(3)) |
| 1946 | PP(2)=DIMAG(S1(3)) |
| 1947 | PP(3)=DIMAG(S1(2)) |
| 1948 | PA(0)=DBLE( S2(2)) |
| 1949 | PA(1)=DBLE( S2(3)) |
| 1950 | PA(2)=DIMAG(S2(3)) |
| 1951 | PA(3)=DIMAG(S2(2)) |
| 1952 | MP2=PP(0)**2-(PP(1)**2+PP(2)**2+PP(3)**2) |
| 1953 | MA2=PA(0)**2-(PA(1)**2+PA(2)**2+PA(3)**2) |
| 1954 | M2D=MP2-MA2 |
| 1955 | C |
| 1956 | JSS(1) = ADG*( (PP(0)-PA(0)) - Q(0)*M2D/VM2) |
| 1957 | JSS(2) = ADG*( (PP(1)-PA(1)) - Q(1)*M2D/VM2) |
| 1958 | JSS(3) = ADG*( (PP(2)-PA(2)) - Q(2)*M2D/VM2) |
| 1959 | JSS(4) = ADG*( (PP(3)-PA(3)) - Q(3)*M2D/VM2) |
| 1960 | C |
| 1961 | RETURN |
| 1962 | C |
| 1963 | 10 ADG=G*S1(1)*S2(1)/Q2 |
| 1964 | C |
| 1965 | JSS(1) = ADG*DBLE( S1(2)-S2(2)) |
| 1966 | JSS(2) = ADG*DBLE( S1(3)-S2(3)) |
| 1967 | JSS(3) = ADG*DIMAG(S1(3)-S2(3)) |
| 1968 | JSS(4) = ADG*DIMAG(S1(2)-S2(2)) |
| 1969 | C |
| 1970 | RETURN |
| 1971 | END |
| 1972 | C |
| 1973 | C |
| 1974 | C ---------------------------------------------------------------------- |
| 1975 | C |
| 1976 | SUBROUTINE JTIOXX(FI,FO,G , JIO) |
| 1977 | C |
| 1978 | C this subroutine computes an off-shell vector current from an external |
| 1979 | C fermion pair. the vector boson propagator is not included in this |
| 1980 | C routine. |
| 1981 | C |
| 1982 | C input: |
| 1983 | C complex fi(6) : flow-in fermion |fi> |
| 1984 | C complex fo(6) : flow-out fermion <fo| |
| 1985 | C real g(2) : coupling constants gvf |
| 1986 | C |
| 1987 | C output: |
| 1988 | C complex jio(6) : vector current j^mu(<fo|v|fi>) |
| 1989 | C |
| 1990 | #if defined(CERNLIB_IMPNONE) |
| 1991 | IMPLICIT NONE |
| 1992 | #endif |
| 1993 | COMPLEX*16 FI(6),FO(6),JIO(6) |
| 1994 | REAL*8 G(2) |
| 1995 | C |
| 1996 | REAL*8 RXZERO, RXONE |
| 1997 | PARAMETER( RXZERO=0.0D0, RXONE=1.0D0 ) |
| 1998 | COMPLEX*16 CXIMAG |
| 1999 | LOGICAL FIRST |
| 2000 | SAVE CXIMAG,FIRST |
| 2001 | DATA FIRST/.TRUE./ |
| 2002 | C |
| 2003 | C Fix compilation with g77 |
| 2004 | IF(FIRST) THEN |
| 2005 | FIRST=.FALSE. |
| 2006 | CXIMAG=DCMPLX( RXZERO, RXONE ) |
| 2007 | ENDIF |
| 2008 | C |
| 2009 | JIO(5) = FO(5)-FI(5) |
| 2010 | JIO(6) = FO(6)-FI(6) |
| 2011 | C |
| 2012 | IF ( G(2) .NE. RXZERO ) THEN |
| 2013 | JIO(1) = ( G(1)*( FO(3)*FI(1)+FO(4)*FI(2)) |
| 2014 | & +G(2)*( FO(1)*FI(3)+FO(2)*FI(4)) ) |
| 2015 | JIO(2) = (-G(1)*( FO(3)*FI(2)+FO(4)*FI(1)) |
| 2016 | & +G(2)*( FO(1)*FI(4)+FO(2)*FI(3)) ) |
| 2017 | JIO(3) = ( G(1)*( FO(3)*FI(2)-FO(4)*FI(1)) |
| 2018 | & +G(2)*(-FO(1)*FI(4)+FO(2)*FI(3)) )*CXIMAG |
| 2019 | JIO(4) = ( G(1)*(-FO(3)*FI(1)+FO(4)*FI(2)) |
| 2020 | & +G(2)*( FO(1)*FI(3)-FO(2)*FI(4)) ) |
| 2021 | C |
| 2022 | ELSE |
| 2023 | JIO(1) = ( FO(3)*FI(1)+FO(4)*FI(2))*G(1) |
| 2024 | JIO(2) = -( FO(3)*FI(2)+FO(4)*FI(1))*G(1) |
| 2025 | JIO(3) = ( FO(3)*FI(2)-FO(4)*FI(1))*DCMPLX(RXZERO,G(1)) |
| 2026 | JIO(4) = (-FO(3)*FI(1)+FO(4)*FI(2))*G(1) |
| 2027 | END IF |
| 2028 | C |
| 2029 | RETURN |
| 2030 | END |
| 2031 | C ---------------------------------------------------------------------- |
| 2032 | C |
| 2033 | SUBROUTINE JVSSXX(VC,S1,S2,G,VMASS,VWIDTH , JVSS) |
| 2034 | C |
| 2035 | C This subroutine computes an off-shell vector current from the vector- |
| 2036 | C vector-scalar-scalar coupling. The vector propagator is given in |
| 2037 | C Feynman gauge for a massless vector and in unitary gauge for a massive |
| 2038 | C vector. |
| 2039 | C |
| 2040 | C INPUT: |
| 2041 | C complex VC(6) : input vector V |
| 2042 | C complex S1(3) : first scalar S1 |
| 2043 | C complex S2(3) : second scalar S2 |
| 2044 | C real G : coupling constant GVVHH |
| 2045 | C real VMASS : mass of OUTPUT vector V' |
| 2046 | C real VWIDTH : width of OUTPUT vector V' |
| 2047 | C |
| 2048 | C OUTPUT: |
| 2049 | C complex JVSS(6) : vector current J^mu(V':V,S1,S2) |
| 2050 | C |
| 2051 | #if defined(CERNLIB_IMPNONE) |
| 2052 | IMPLICIT NONE |
| 2053 | #endif |
| 2054 | COMPLEX*16 VC(6),S1(3),S2(3),JVSS(6),DG |
| 2055 | REAL*8 Q(0:3),G,VMASS,VWIDTH,Q2,VK,VM2 |
| 2056 | C |
| 2057 | JVSS(5) = VC(5)+S1(2)+S2(2) |
| 2058 | JVSS(6) = VC(6)+S1(3)+S2(3) |
| 2059 | C |
| 2060 | Q(0)=DBLE( JVSS(5)) |
| 2061 | Q(1)=DBLE( JVSS(6)) |
| 2062 | Q(2)=DIMAG(JVSS(6)) |
| 2063 | Q(3)=DIMAG(JVSS(5)) |
| 2064 | Q2=Q(0)**2-(Q(1)**2+Q(2)**2+Q(3)**2) |
| 2065 | VM2=VMASS**2 |
| 2066 | C |
| 2067 | IF (VMASS.EQ.0.) GOTO 10 |
| 2068 | C |
| 2069 | DG=G*S1(1)*S2(1)/DCMPLX( Q2-VM2,MAX(SIGN( VMASS*VWIDTH,Q2),0.D0)) |
| 2070 | C For the running width, use below instead of the above DG. |
| 2071 | C DG=G*S1(1)*S2(1)/CMPLX( Q2-VM2 , MAX( VWIDTH*Q2/VMASS ,0.)) |
| 2072 | C |
| 2073 | VK=(Q(0)*VC(1)-Q(1)*VC(2)-Q(2)*VC(3)-Q(3)*VC(4))/VM2 |
| 2074 | C |
| 2075 | JVSS(1) = DG*(VC(1)-VK*Q(0)) |
| 2076 | JVSS(2) = DG*(VC(2)-VK*Q(1)) |
| 2077 | JVSS(3) = DG*(VC(3)-VK*Q(2)) |
| 2078 | JVSS(4) = DG*(VC(4)-VK*Q(3)) |
| 2079 | C |
| 2080 | RETURN |
| 2081 | C |
| 2082 | 10 DG= G*S1(1)*S2(1)/Q2 |
| 2083 | C |
| 2084 | JVSS(1) = DG*VC(1) |
| 2085 | JVSS(2) = DG*VC(2) |
| 2086 | JVSS(3) = DG*VC(3) |
| 2087 | JVSS(4) = DG*VC(4) |
| 2088 | C |
| 2089 | RETURN |
| 2090 | END |
| 2091 | C |
| 2092 | C |
| 2093 | C ---------------------------------------------------------------------- |
| 2094 | C |
| 2095 | SUBROUTINE JVSXXX(VC,SC,G,VMASS,VWIDTH , JVS) |
| 2096 | C |
| 2097 | C this subroutine computes an off-shell vector current from the vector- |
| 2098 | C vector-scalar coupling. the vector propagator is given in feynman |
| 2099 | C gauge for a massless vector and in unitary gauge for a massive vector. |
| 2100 | C |
| 2101 | C input: |
| 2102 | C complex vc(6) : input vector v |
| 2103 | C complex sc(3) : input scalar s |
| 2104 | C real g : coupling constant gvvh |
| 2105 | C real vmass : mass of output vector v' |
| 2106 | C real vwidth : width of output vector v' |
| 2107 | C |
| 2108 | C output: |
| 2109 | C complex jvs(6) : vector current j^mu(v':v,s) |
| 2110 | C |
| 2111 | #if defined(CERNLIB_IMPNONE) |
| 2112 | IMPLICIT NONE |
| 2113 | #endif |
| 2114 | COMPLEX*16 VC(6),SC(3),JVS(6),DG,VK |
| 2115 | REAL*8 Q(0:3),VMASS,VWIDTH,Q2,VM2,G |
| 2116 | C |
| 2117 | JVS(5) = VC(5)+SC(2) |
| 2118 | JVS(6) = VC(6)+SC(3) |
| 2119 | C |
| 2120 | Q(0)=DBLE( JVS(5)) |
| 2121 | Q(1)=DBLE( JVS(6)) |
| 2122 | Q(2)=DIMAG(JVS(6)) |
| 2123 | Q(3)=DIMAG(JVS(5)) |
| 2124 | Q2=Q(0)**2-(Q(1)**2+Q(2)**2+Q(3)**2) |
| 2125 | VM2=VMASS**2 |
| 2126 | C |
| 2127 | IF (VMASS.EQ.0.) GOTO 10 |
| 2128 | C |
| 2129 | DG=G*SC(1)/DCMPLX( Q2-VM2 , MAX(DSIGN( VMASS*VWIDTH ,Q2),0.D0) ) |
| 2130 | C for the running width, use below instead of the above dg. |
| 2131 | C dg=g*sc(1)/dcmplx( q2-vm2 , max( vwidth*q2/vmass ,0.) ) |
| 2132 | C |
| 2133 | VK=(-Q(0)*VC(1)+Q(1)*VC(2)+Q(2)*VC(3)+Q(3)*VC(4))/VM2 |
| 2134 | C |
| 2135 | JVS(1) = DG*(Q(0)*VK+VC(1)) |
| 2136 | JVS(2) = DG*(Q(1)*VK+VC(2)) |
| 2137 | JVS(3) = DG*(Q(2)*VK+VC(3)) |
| 2138 | JVS(4) = DG*(Q(3)*VK+VC(4)) |
| 2139 | C |
| 2140 | RETURN |
| 2141 | C |
| 2142 | 10 DG=G*SC(1)/Q2 |
| 2143 | C |
| 2144 | JVS(1) = DG*VC(1) |
| 2145 | JVS(2) = DG*VC(2) |
| 2146 | JVS(3) = DG*VC(3) |
| 2147 | JVS(4) = DG*VC(4) |
| 2148 | C |
| 2149 | RETURN |
| 2150 | END |
| 2151 | |
| 2152 | |
| 2153 | C |
| 2154 | C ---------------------------------------------------------------------- |
| 2155 | C |
| 2156 | SUBROUTINE JVVXXX(V1,V2,G,VMASS,VWIDTH , JVV) |
| 2157 | C |
| 2158 | C this subroutine computes an off-shell vector current from the three- |
| 2159 | C point gauge boson coupling. the vector propagator is given in feynman |
| 2160 | C gauge for a massless vector and in unitary gauge for a massive vector. |
| 2161 | C |
| 2162 | C input: |
| 2163 | C complex v1(6) : first vector v1 |
| 2164 | C complex v2(6) : second vector v2 |
| 2165 | C real g : coupling constant (see the table below) |
| 2166 | C real vmass : mass of output vector v |
| 2167 | C real vwidth : width of output vector v |
| 2168 | C |
| 2169 | C the possible sets of the inputs are as follows: |
| 2170 | C ------------------------------------------------------------------ |
| 2171 | C | v1 | v2 | jvv | g | vmass | vwidth | |
| 2172 | C ------------------------------------------------------------------ |
| 2173 | C | w- | w+ | a/z | gwwa/gwwz | 0./zmass | 0./zwidth | |
| 2174 | C | w3/a/z | w- | w+ | gw/gwwa/gwwz | wmass | wwidth | |
| 2175 | C | w+ | w3/a/z | w- | gw/gwwa/gwwz | wmass | wwidth | |
| 2176 | C ------------------------------------------------------------------ |
| 2177 | C where all the bosons are defined by the flowing-out quantum number. |
| 2178 | C |
| 2179 | C output: |
| 2180 | C complex jvv(6) : vector current j^mu(v:v1,v2) |
| 2181 | C |
| 2182 | #if defined(CERNLIB_IMPNONE) |
| 2183 | IMPLICIT NONE |
| 2184 | #endif |
| 2185 | COMPLEX*16 V1(6),V2(6),JVV(6),J12(0:3),JS,DG, |
| 2186 | & SV1,SV2,S11,S12,S21,S22,V12 |
| 2187 | REAL*8 P1(0:3),P2(0:3),Q(0:3),G,VMASS,VWIDTH,GS,S,VM2,M1,M2 |
| 2188 | C |
| 2189 | REAL*8 RXZERO |
| 2190 | PARAMETER( RXZERO=0.0D0 ) |
| 2191 | C |
| 2192 | JVV(5) = V1(5)+V2(5) |
| 2193 | JVV(6) = V1(6)+V2(6) |
| 2194 | C |
| 2195 | P1(0)=DBLE( V1(5)) |
| 2196 | P1(1)=DBLE( V1(6)) |
| 2197 | P1(2)=DIMAG(V1(6)) |
| 2198 | P1(3)=DIMAG(V1(5)) |
| 2199 | P2(0)=DBLE( V2(5)) |
| 2200 | P2(1)=DBLE( V2(6)) |
| 2201 | P2(2)=DIMAG(V2(6)) |
| 2202 | P2(3)=DIMAG(V2(5)) |
| 2203 | Q(0)=-DBLE( JVV(5)) |
| 2204 | Q(1)=-DBLE( JVV(6)) |
| 2205 | Q(2)=-DIMAG(JVV(6)) |
| 2206 | Q(3)=-DIMAG(JVV(5)) |
| 2207 | S=Q(0)**2-(Q(1)**2+Q(2)**2+Q(3)**2) |
| 2208 | C |
| 2209 | V12=V1(1)*V2(1)-V1(2)*V2(2)-V1(3)*V2(3)-V1(4)*V2(4) |
| 2210 | SV1= (P2(0)-Q(0))*V1(1) -(P2(1)-Q(1))*V1(2) |
| 2211 | & -(P2(2)-Q(2))*V1(3) -(P2(3)-Q(3))*V1(4) |
| 2212 | SV2=-(P1(0)-Q(0))*V2(1) +(P1(1)-Q(1))*V2(2) |
| 2213 | & +(P1(2)-Q(2))*V2(3) +(P1(3)-Q(3))*V2(4) |
| 2214 | J12(0)=(P1(0)-P2(0))*V12 +SV1*V2(1) +SV2*V1(1) |
| 2215 | J12(1)=(P1(1)-P2(1))*V12 +SV1*V2(2) +SV2*V1(2) |
| 2216 | J12(2)=(P1(2)-P2(2))*V12 +SV1*V2(3) +SV2*V1(3) |
| 2217 | J12(3)=(P1(3)-P2(3))*V12 +SV1*V2(4) +SV2*V1(4) |
| 2218 | C |
| 2219 | IF ( VMASS .NE. RXZERO ) THEN |
| 2220 | VM2=VMASS**2 |
| 2221 | M1=P1(0)**2-(P1(1)**2+P1(2)**2+P1(3)**2) |
| 2222 | M2=P2(0)**2-(P2(1)**2+P2(2)**2+P2(3)**2) |
| 2223 | S11=P1(0)*V1(1)-P1(1)*V1(2)-P1(2)*V1(3)-P1(3)*V1(4) |
| 2224 | S12=P1(0)*V2(1)-P1(1)*V2(2)-P1(2)*V2(3)-P1(3)*V2(4) |
| 2225 | S21=P2(0)*V1(1)-P2(1)*V1(2)-P2(2)*V1(3)-P2(3)*V1(4) |
| 2226 | S22=P2(0)*V2(1)-P2(1)*V2(2)-P2(2)*V2(3)-P2(3)*V2(4) |
| 2227 | JS=(V12*(-M1+M2) +S11*S12 -S21*S22)/VM2 |
| 2228 | C |
| 2229 | DG=-G/DCMPLX( S-VM2 , MAX(SIGN( VMASS*VWIDTH ,S),RXZERO) ) |
| 2230 | C |
| 2231 | C for the running width, use below instead of the above dg. |
| 2232 | C dg=-g/dcmplx( s-vm2 , max( vwidth*s/vmass ,r_zero) ) |
| 2233 | C |
| 2234 | JVV(1) = DG*(J12(0)-Q(0)*JS) |
| 2235 | JVV(2) = DG*(J12(1)-Q(1)*JS) |
| 2236 | JVV(3) = DG*(J12(2)-Q(2)*JS) |
| 2237 | JVV(4) = DG*(J12(3)-Q(3)*JS) |
| 2238 | C |
| 2239 | ELSE |
| 2240 | GS=-G/S |
| 2241 | C |
| 2242 | JVV(1) = GS*J12(0) |
| 2243 | JVV(2) = GS*J12(1) |
| 2244 | JVV(3) = GS*J12(2) |
| 2245 | JVV(4) = GS*J12(3) |
| 2246 | END IF |
| 2247 | C |
| 2248 | RETURN |
| 2249 | END |
| 2250 | C |
| 2251 | C ---------------------------------------------------------------------- |
| 2252 | C |
| 2253 | SUBROUTINE JW3WXX(W1,W2,W3,G1,G2,WMASS,WWIDTH,VMASS,VWIDTH , JW3W) |
| 2254 | C |
| 2255 | C this subroutine computes an off-shell w+, w-, w3, z or photon current |
| 2256 | C from the four-point gauge boson coupling, including the contributions |
| 2257 | C of w exchange diagrams. the vector propagator is given in feynman |
| 2258 | C gauge for a photon and in unitary gauge for w and z bosons. if one |
| 2259 | C sets wmass=0.0, then the ggg-->g current is given (see sect 2.9.1 of |
| 2260 | C the manual). |
| 2261 | C |
| 2262 | C input: |
| 2263 | C complex w1(6) : first vector w1 |
| 2264 | C complex w2(6) : second vector w2 |
| 2265 | C complex w3(6) : third vector w3 |
| 2266 | C real g1 : first coupling constant |
| 2267 | C real g2 : second coupling constant |
| 2268 | C (see the table below) |
| 2269 | C real wmass : mass of internal w |
| 2270 | C real wwidth : width of internal w |
| 2271 | C real vmass : mass of output w' |
| 2272 | C real vwidth : width of output w' |
| 2273 | C |
| 2274 | C the possible sets of the inputs are as follows: |
| 2275 | C ------------------------------------------------------------------- |
| 2276 | C | w1 | w2 | w3 | g1 | g2 |wmass|wwidth|vmass|vwidth || jw3w | |
| 2277 | C ------------------------------------------------------------------- |
| 2278 | C | w- | w3 | w+ | gw |gwwz|wmass|wwidth|zmass|zwidth || z | |
| 2279 | C | w- | w3 | w+ | gw |gwwa|wmass|wwidth| 0. | 0. || a | |
| 2280 | C | w- | z | w+ |gwwz|gwwz|wmass|wwidth|zmass|zwidth || z | |
| 2281 | C | w- | z | w+ |gwwz|gwwa|wmass|wwidth| 0. | 0. || a | |
| 2282 | C | w- | a | w+ |gwwa|gwwz|wmass|wwidth|zmass|zwidth || z | |
| 2283 | C | w- | a | w+ |gwwa|gwwa|wmass|wwidth| 0. | 0. || a | |
| 2284 | C ------------------------------------------------------------------- |
| 2285 | C | w3 | w- | w3 | gw | gw |wmass|wwidth|wmass|wwidth || w+ | |
| 2286 | C | w3 | w+ | w3 | gw | gw |wmass|wwidth|wmass|wwidth || w- | |
| 2287 | C | w3 | w- | z | gw |gwwz|wmass|wwidth|wmass|wwidth || w+ | |
| 2288 | C | w3 | w+ | z | gw |gwwz|wmass|wwidth|wmass|wwidth || w- | |
| 2289 | C | w3 | w- | a | gw |gwwa|wmass|wwidth|wmass|wwidth || w+ | |
| 2290 | C | w3 | w+ | a | gw |gwwa|wmass|wwidth|wmass|wwidth || w- | |
| 2291 | C | z | w- | z |gwwz|gwwz|wmass|wwidth|wmass|wwidth || w+ | |
| 2292 | C | z | w+ | z |gwwz|gwwz|wmass|wwidth|wmass|wwidth || w- | |
| 2293 | C | z | w- | a |gwwz|gwwa|wmass|wwidth|wmass|wwidth || w+ | |
| 2294 | C | z | w+ | a |gwwz|gwwa|wmass|wwidth|wmass|wwidth || w- | |
| 2295 | C | a | w- | a |gwwa|gwwa|wmass|wwidth|wmass|wwidth || w+ | |
| 2296 | C | a | w+ | a |gwwa|gwwa|wmass|wwidth|wmass|wwidth || w- | |
| 2297 | C ------------------------------------------------------------------- |
| 2298 | C where all the bosons are defined by the flowing-out quantum number. |
| 2299 | C |
| 2300 | C output: |
| 2301 | C complex jw3w(6) : w current j^mu(w':w1,w2,w3) |
| 2302 | C |
| 2303 | #if defined(CERNLIB_IMPNONE) |
| 2304 | IMPLICIT NONE |
| 2305 | #endif |
| 2306 | COMPLEX*16 W1(6),W2(6),W3(6),JW3W(6) |
| 2307 | COMPLEX*16 DW1(0:3),DW2(0:3),DW3(0:3), |
| 2308 | & JJ(0:3),J4(0:3), |
| 2309 | & DV,W12,W32,W13, |
| 2310 | & JQ |
| 2311 | REAL*8 G1,G2,WMASS,WWIDTH,VMASS,VWIDTH |
| 2312 | REAL*8 P1(0:3),P2(0:3),P3(0:3),Q(0:3), |
| 2313 | & DG2,DMV,DWV,MV2,Q2 |
| 2314 | C |
| 2315 | REAL*8 RXZERO |
| 2316 | PARAMETER( RXZERO=0.0D0 ) |
| 2317 | C |
| 2318 | JW3W(5) = W1(5)+W2(5)+W3(5) |
| 2319 | JW3W(6) = W1(6)+W2(6)+W3(6) |
| 2320 | C |
| 2321 | DW1(0)=DCMPLX(W1(1)) |
| 2322 | DW1(1)=DCMPLX(W1(2)) |
| 2323 | DW1(2)=DCMPLX(W1(3)) |
| 2324 | DW1(3)=DCMPLX(W1(4)) |
| 2325 | DW2(0)=DCMPLX(W2(1)) |
| 2326 | DW2(1)=DCMPLX(W2(2)) |
| 2327 | DW2(2)=DCMPLX(W2(3)) |
| 2328 | DW2(3)=DCMPLX(W2(4)) |
| 2329 | DW3(0)=DCMPLX(W3(1)) |
| 2330 | DW3(1)=DCMPLX(W3(2)) |
| 2331 | DW3(2)=DCMPLX(W3(3)) |
| 2332 | DW3(3)=DCMPLX(W3(4)) |
| 2333 | P1(0)=DBLE( W1(5)) |
| 2334 | P1(1)=DBLE( W1(6)) |
| 2335 | P1(2)=DBLE(DIMAG(W1(6))) |
| 2336 | P1(3)=DBLE(DIMAG(W1(5))) |
| 2337 | P2(0)=DBLE( W2(5)) |
| 2338 | P2(1)=DBLE( W2(6)) |
| 2339 | P2(2)=DBLE(DIMAG(W2(6))) |
| 2340 | P2(3)=DBLE(DIMAG(W2(5))) |
| 2341 | P3(0)=DBLE( W3(5)) |
| 2342 | P3(1)=DBLE( W3(6)) |
| 2343 | P3(2)=DBLE(DIMAG(W3(6))) |
| 2344 | P3(3)=DBLE(DIMAG(W3(5))) |
| 2345 | Q(0)=-(P1(0)+P2(0)+P3(0)) |
| 2346 | Q(1)=-(P1(1)+P2(1)+P3(1)) |
| 2347 | Q(2)=-(P1(2)+P2(2)+P3(2)) |
| 2348 | Q(3)=-(P1(3)+P2(3)+P3(3)) |
| 2349 | |
| 2350 | |
| 2351 | Q2 =Q(0)**2 -(Q(1)**2 +Q(2)**2 +Q(3)**2) |
| 2352 | DG2=DBLE(G1)*DBLE(G2) |
| 2353 | DMV=DBLE(VMASS) |
| 2354 | DWV=DBLE(VWIDTH) |
| 2355 | MV2=DMV**2 |
| 2356 | IF (VMASS.EQ. RXZERO) THEN |
| 2357 | DV = 1.0D0/DCMPLX( Q2 ) |
| 2358 | ELSE |
| 2359 | DV = 1.0D0/DCMPLX( Q2 -MV2 , DMAX1(DSIGN(DMV*DWV,Q2 ),0.D0) ) |
| 2360 | ENDIF |
| 2361 | C for the running width, use below instead of the above dv. |
| 2362 | C dv = 1.0d0/dcmplx( q2 -mv2 , dmax1(dwv*q2/dmv,0.d0) ) |
| 2363 | C |
| 2364 | W12=DW1(0)*DW2(0)-DW1(1)*DW2(1)-DW1(2)*DW2(2)-DW1(3)*DW2(3) |
| 2365 | W32=DW3(0)*DW2(0)-DW3(1)*DW2(1)-DW3(2)*DW2(2)-DW3(3)*DW2(3) |
| 2366 | C |
| 2367 | IF ( WMASS .NE. RXZERO ) THEN |
| 2368 | W13=DW1(0)*DW3(0)-DW1(1)*DW3(1)-DW1(2)*DW3(2)-DW1(3)*DW3(3) |
| 2369 | C |
| 2370 | J4(0)=DG2*( DW1(0)*W32 + DW3(0)*W12 - 2.D0*DW2(0)*W13 ) |
| 2371 | J4(1)=DG2*( DW1(1)*W32 + DW3(1)*W12 - 2.D0*DW2(1)*W13 ) |
| 2372 | J4(2)=DG2*( DW1(2)*W32 + DW3(2)*W12 - 2.D0*DW2(2)*W13 ) |
| 2373 | J4(3)=DG2*( DW1(3)*W32 + DW3(3)*W12 - 2.D0*DW2(3)*W13 ) |
| 2374 | C |
| 2375 | JJ(0)=J4(0) |
| 2376 | JJ(1)=J4(1) |
| 2377 | JJ(2)=J4(2) |
| 2378 | JJ(3)=J4(3) |
| 2379 | |
| 2380 | ELSE |
| 2381 | C |
| 2382 | W12=DW1(0)*DW2(0)-DW1(1)*DW2(1)-DW1(2)*DW2(2)-DW1(3)*DW2(3) |
| 2383 | W32=DW3(0)*DW2(0)-DW3(1)*DW2(1)-DW3(2)*DW2(2)-DW3(3)*DW2(3) |
| 2384 | W13=DW1(0)*DW3(0)-DW1(1)*DW3(1)-DW1(2)*DW3(2)-DW1(3)*DW3(3) |
| 2385 | C |
| 2386 | J4(0)=DG2*( DW1(0)*W32 - DW2(0)*W13 ) |
| 2387 | J4(1)=DG2*( DW1(1)*W32 - DW2(1)*W13 ) |
| 2388 | J4(2)=DG2*( DW1(2)*W32 - DW2(2)*W13 ) |
| 2389 | J4(3)=DG2*( DW1(3)*W32 - DW2(3)*W13 ) |
| 2390 | C |
| 2391 | JJ(0)=J4(0) |
| 2392 | JJ(1)=J4(1) |
| 2393 | JJ(2)=J4(2) |
| 2394 | JJ(3)=J4(3) |
| 2395 | |
| 2396 | END IF |
| 2397 | C |
| 2398 | IF ( VMASS .NE. RXZERO ) THEN |
| 2399 | C |
| 2400 | JQ=(JJ(0)*Q(0)-JJ(1)*Q(1)-JJ(2)*Q(2)-JJ(3)*Q(3))/MV2 |
| 2401 | C |
| 2402 | JW3W(1) = DCMPLX( (JJ(0)-JQ*Q(0))*DV ) |
| 2403 | JW3W(2) = DCMPLX( (JJ(1)-JQ*Q(1))*DV ) |
| 2404 | JW3W(3) = DCMPLX( (JJ(2)-JQ*Q(2))*DV ) |
| 2405 | JW3W(4) = DCMPLX( (JJ(3)-JQ*Q(3))*DV ) |
| 2406 | C |
| 2407 | ELSE |
| 2408 | C |
| 2409 | JW3W(1) = DCMPLX( JJ(0)*DV ) |
| 2410 | JW3W(2) = DCMPLX( JJ(1)*DV ) |
| 2411 | JW3W(3) = DCMPLX( JJ(2)*DV ) |
| 2412 | JW3W(4) = DCMPLX( JJ(3)*DV ) |
| 2413 | END IF |
| 2414 | C |
| 2415 | RETURN |
| 2416 | END |
| 2417 | C |
| 2418 | C ---------------------------------------------------------------------- |
| 2419 | C |
| 2420 | SUBROUTINE JWWWXX(W1,W2,W3,GWWA,GWWZ,ZMASS,ZWIDTH,WMASS,WWIDTH , |
| 2421 | & JWWW) |
| 2422 | C |
| 2423 | C this subroutine computes an off-shell w+/w- current from the four- |
| 2424 | C point gauge boson coupling, including the contributions of photon and |
| 2425 | C z exchanges. the vector propagators for the output w and the internal |
| 2426 | C z bosons are given in unitary gauge, and that of the internal photon |
| 2427 | C is given in feynman gauge. |
| 2428 | C |
| 2429 | C input: |
| 2430 | C complex w1(6) : first vector w1 |
| 2431 | C complex w2(6) : second vector w2 |
| 2432 | C complex w3(6) : third vector w3 |
| 2433 | C real gwwa : coupling constant of w and a gwwa |
| 2434 | C real gwwz : coupling constant of w and z gwwz |
| 2435 | C real zmass : mass of internal z |
| 2436 | C real zwidth : width of internal z |
| 2437 | C real wmass : mass of output w |
| 2438 | C real wwidth : width of output w |
| 2439 | C |
| 2440 | C the possible sets of the inputs are as follows: |
| 2441 | C ------------------------------------------------------------------- |
| 2442 | C | w1 | w2 | w3 |gwwa|gwwz|zmass|zwidth|wmass|wwidth || jwww | |
| 2443 | C ------------------------------------------------------------------- |
| 2444 | C | w- | w+ | w- |gwwa|gwwz|zmass|zwidth|wmass|wwidth || w+ | |
| 2445 | C | w+ | w- | w+ |gwwa|gwwz|zmass|zwidth|wmass|wwidth || w- | |
| 2446 | C ------------------------------------------------------------------- |
| 2447 | C where all the bosons are defined by the flowing-out quantum number. |
| 2448 | C |
| 2449 | C output: |
| 2450 | C complex jwww(6) : w current j^mu(w':w1,w2,w3) |
| 2451 | C |
| 2452 | #if defined(CERNLIB_IMPNONE) |
| 2453 | IMPLICIT NONE |
| 2454 | #endif |
| 2455 | COMPLEX*16 W1(6),W2(6),W3(6),JWWW(6) |
| 2456 | COMPLEX*16 DW1(0:3),DW2(0:3),DW3(0:3), |
| 2457 | & JJ(0:3),JS(0:3),JT(0:3),J4(0:3), |
| 2458 | & JT12(0:3),JT32(0:3),J12(0:3),J32(0:3), |
| 2459 | & DZS,DZT,DW,W12,W32,W13,P1W2,P2W1,P3W2,P2W3, |
| 2460 | & JK12,JK32,JSW3,JTW1,P3JS,KSW3,P1JT,KTW1,JQ |
| 2461 | REAL*8 GWWA,GWWZ,ZMASS,ZWIDTH,WMASS,WWIDTH |
| 2462 | REAL*8 P1(0:3),P2(0:3),P3(0:3),Q(0:3),KS(0:3),KT(0:3), |
| 2463 | & DGWWA2,DGWWZ2,DGW2,DMZ,DWZ,DMW,DWW,MZ2,MW2,Q2,KS2,KT2, |
| 2464 | & DAS,DAT |
| 2465 | C |
| 2466 | JWWW(5) = W1(5)+W2(5)+W3(5) |
| 2467 | JWWW(6) = W1(6)+W2(6)+W3(6) |
| 2468 | C |
| 2469 | DW1(0)=DCMPLX(W1(1)) |
| 2470 | DW1(1)=DCMPLX(W1(2)) |
| 2471 | DW1(2)=DCMPLX(W1(3)) |
| 2472 | DW1(3)=DCMPLX(W1(4)) |
| 2473 | DW2(0)=DCMPLX(W2(1)) |
| 2474 | DW2(1)=DCMPLX(W2(2)) |
| 2475 | DW2(2)=DCMPLX(W2(3)) |
| 2476 | DW2(3)=DCMPLX(W2(4)) |
| 2477 | DW3(0)=DCMPLX(W3(1)) |
| 2478 | DW3(1)=DCMPLX(W3(2)) |
| 2479 | DW3(2)=DCMPLX(W3(3)) |
| 2480 | DW3(3)=DCMPLX(W3(4)) |
| 2481 | P1(0)=DBLE( W1(5)) |
| 2482 | P1(1)=DBLE( W1(6)) |
| 2483 | P1(2)=DBLE(DIMAG(W1(6))) |
| 2484 | P1(3)=DBLE(DIMAG(W1(5))) |
| 2485 | P2(0)=DBLE( W2(5)) |
| 2486 | P2(1)=DBLE( W2(6)) |
| 2487 | P2(2)=DBLE(DIMAG(W2(6))) |
| 2488 | P2(3)=DBLE(DIMAG(W2(5))) |
| 2489 | P3(0)=DBLE( W3(5)) |
| 2490 | P3(1)=DBLE( W3(6)) |
| 2491 | P3(2)=DBLE(DIMAG(W3(6))) |
| 2492 | P3(3)=DBLE(DIMAG(W3(5))) |
| 2493 | Q(0)=-(P1(0)+P2(0)+P3(0)) |
| 2494 | Q(1)=-(P1(1)+P2(1)+P3(1)) |
| 2495 | Q(2)=-(P1(2)+P2(2)+P3(2)) |
| 2496 | Q(3)=-(P1(3)+P2(3)+P3(3)) |
| 2497 | KS(0)=P1(0)+P2(0) |
| 2498 | KS(1)=P1(1)+P2(1) |
| 2499 | KS(2)=P1(2)+P2(2) |
| 2500 | KS(3)=P1(3)+P2(3) |
| 2501 | KT(0)=P2(0)+P3(0) |
| 2502 | KT(1)=P2(1)+P3(1) |
| 2503 | KT(2)=P2(2)+P3(2) |
| 2504 | KT(3)=P2(3)+P3(3) |
| 2505 | Q2 =Q(0)**2 -(Q(1)**2 +Q(2)**2 +Q(3)**2) |
| 2506 | KS2=KS(0)**2-(KS(1)**2+KS(2)**2+KS(3)**2) |
| 2507 | KT2=KT(0)**2-(KT(1)**2+KT(2)**2+KT(3)**2) |
| 2508 | DGWWA2=DBLE(GWWA)**2 |
| 2509 | DGWWZ2=DBLE(GWWZ)**2 |
| 2510 | DGW2 =DGWWA2+DGWWZ2 |
| 2511 | DMZ=DBLE(ZMASS) |
| 2512 | DWZ=DBLE(ZWIDTH) |
| 2513 | DMW=DBLE(WMASS) |
| 2514 | DWW=DBLE(WWIDTH) |
| 2515 | MZ2=DMZ**2 |
| 2516 | MW2=DMW**2 |
| 2517 | C |
| 2518 | DAS=-DGWWA2/KS2 |
| 2519 | DAT=-DGWWA2/KT2 |
| 2520 | DZS=-DGWWZ2/DCMPLX( KS2-MZ2 , DMAX1(DSIGN(DMZ*DWZ,KS2),0.D0) ) |
| 2521 | DZT=-DGWWZ2/DCMPLX( KT2-MZ2 , DMAX1(DSIGN(DMZ*DWZ,KT2),0.D0) ) |
| 2522 | DW =-1.0D0/DCMPLX( Q2 -MW2 , DMAX1(DSIGN(DMW*DWW,Q2 ),0.D0) ) |
| 2523 | C for the running width, use below instead of the above dw. |
| 2524 | C dw =-1.0d0/dcmplx( q2 -mw2 , dmax1(dww*q2/dmw,0.d0) ) |
| 2525 | C |
| 2526 | W12=DW1(0)*DW2(0)-DW1(1)*DW2(1)-DW1(2)*DW2(2)-DW1(3)*DW2(3) |
| 2527 | W32=DW3(0)*DW2(0)-DW3(1)*DW2(1)-DW3(2)*DW2(2)-DW3(3)*DW2(3) |
| 2528 | C |
| 2529 | P1W2= (P1(0)+KS(0))*DW2(0)-(P1(1)+KS(1))*DW2(1) |
| 2530 | & -(P1(2)+KS(2))*DW2(2)-(P1(3)+KS(3))*DW2(3) |
| 2531 | P2W1= (P2(0)+KS(0))*DW1(0)-(P2(1)+KS(1))*DW1(1) |
| 2532 | & -(P2(2)+KS(2))*DW1(2)-(P2(3)+KS(3))*DW1(3) |
| 2533 | P3W2= (P3(0)+KT(0))*DW2(0)-(P3(1)+KT(1))*DW2(1) |
| 2534 | & -(P3(2)+KT(2))*DW2(2)-(P3(3)+KT(3))*DW2(3) |
| 2535 | P2W3= (P2(0)+KT(0))*DW3(0)-(P2(1)+KT(1))*DW3(1) |
| 2536 | & -(P2(2)+KT(2))*DW3(2)-(P2(3)+KT(3))*DW3(3) |
| 2537 | C |
| 2538 | JT12(0)= (P1(0)-P2(0))*W12 + P2W1*DW2(0) - P1W2*DW1(0) |
| 2539 | JT12(1)= (P1(1)-P2(1))*W12 + P2W1*DW2(1) - P1W2*DW1(1) |
| 2540 | JT12(2)= (P1(2)-P2(2))*W12 + P2W1*DW2(2) - P1W2*DW1(2) |
| 2541 | JT12(3)= (P1(3)-P2(3))*W12 + P2W1*DW2(3) - P1W2*DW1(3) |
| 2542 | JT32(0)= (P3(0)-P2(0))*W32 + P2W3*DW2(0) - P3W2*DW3(0) |
| 2543 | JT32(1)= (P3(1)-P2(1))*W32 + P2W3*DW2(1) - P3W2*DW3(1) |
| 2544 | JT32(2)= (P3(2)-P2(2))*W32 + P2W3*DW2(2) - P3W2*DW3(2) |
| 2545 | JT32(3)= (P3(3)-P2(3))*W32 + P2W3*DW2(3) - P3W2*DW3(3) |
| 2546 | C |
| 2547 | JK12=(JT12(0)*KS(0)-JT12(1)*KS(1)-JT12(2)*KS(2)-JT12(3)*KS(3))/MZ2 |
| 2548 | JK32=(JT32(0)*KT(0)-JT32(1)*KT(1)-JT32(2)*KT(2)-JT32(3)*KT(3))/MZ2 |
| 2549 | C |
| 2550 | J12(0)=JT12(0)*(DAS+DZS)-KS(0)*JK12*DZS |
| 2551 | J12(1)=JT12(1)*(DAS+DZS)-KS(1)*JK12*DZS |
| 2552 | J12(2)=JT12(2)*(DAS+DZS)-KS(2)*JK12*DZS |
| 2553 | J12(3)=JT12(3)*(DAS+DZS)-KS(3)*JK12*DZS |
| 2554 | J32(0)=JT32(0)*(DAT+DZT)-KT(0)*JK32*DZT |
| 2555 | J32(1)=JT32(1)*(DAT+DZT)-KT(1)*JK32*DZT |
| 2556 | J32(2)=JT32(2)*(DAT+DZT)-KT(2)*JK32*DZT |
| 2557 | J32(3)=JT32(3)*(DAT+DZT)-KT(3)*JK32*DZT |
| 2558 | C |
| 2559 | JSW3=J12(0)*DW3(0)-J12(1)*DW3(1)-J12(2)*DW3(2)-J12(3)*DW3(3) |
| 2560 | JTW1=J32(0)*DW1(0)-J32(1)*DW1(1)-J32(2)*DW1(2)-J32(3)*DW1(3) |
| 2561 | C |
| 2562 | P3JS= (P3(0)-Q(0))*J12(0)-(P3(1)-Q(1))*J12(1) |
| 2563 | & -(P3(2)-Q(2))*J12(2)-(P3(3)-Q(3))*J12(3) |
| 2564 | KSW3= (KS(0)-Q(0))*DW3(0)-(KS(1)-Q(1))*DW3(1) |
| 2565 | & -(KS(2)-Q(2))*DW3(2)-(KS(3)-Q(3))*DW3(3) |
| 2566 | P1JT= (P1(0)-Q(0))*J32(0)-(P1(1)-Q(1))*J32(1) |
| 2567 | & -(P1(2)-Q(2))*J32(2)-(P1(3)-Q(3))*J32(3) |
| 2568 | KTW1= (KT(0)-Q(0))*DW1(0)-(KT(1)-Q(1))*DW1(1) |
| 2569 | & -(KT(2)-Q(2))*DW1(2)-(KT(3)-Q(3))*DW1(3) |
| 2570 | C |
| 2571 | JS(0)= (KS(0)-P3(0))*JSW3 + P3JS*DW3(0) - KSW3*J12(0) |
| 2572 | JS(1)= (KS(1)-P3(1))*JSW3 + P3JS*DW3(1) - KSW3*J12(1) |
| 2573 | JS(2)= (KS(2)-P3(2))*JSW3 + P3JS*DW3(2) - KSW3*J12(2) |
| 2574 | JS(3)= (KS(3)-P3(3))*JSW3 + P3JS*DW3(3) - KSW3*J12(3) |
| 2575 | JT(0)= (KT(0)-P1(0))*JTW1 + P1JT*DW1(0) - KTW1*J32(0) |
| 2576 | JT(1)= (KT(1)-P1(1))*JTW1 + P1JT*DW1(1) - KTW1*J32(1) |
| 2577 | JT(2)= (KT(2)-P1(2))*JTW1 + P1JT*DW1(2) - KTW1*J32(2) |
| 2578 | JT(3)= (KT(3)-P1(3))*JTW1 + P1JT*DW1(3) - KTW1*J32(3) |
| 2579 | C |
| 2580 | W13=DW1(0)*DW3(0)-DW1(1)*DW3(1)-DW1(2)*DW3(2)-DW1(3)*DW3(3) |
| 2581 | C |
| 2582 | J4(0)=DGW2*( DW1(0)*W32 + DW3(0)*W12 - 2.D0*DW2(0)*W13 ) |
| 2583 | J4(1)=DGW2*( DW1(1)*W32 + DW3(1)*W12 - 2.D0*DW2(1)*W13 ) |
| 2584 | J4(2)=DGW2*( DW1(2)*W32 + DW3(2)*W12 - 2.D0*DW2(2)*W13 ) |
| 2585 | J4(3)=DGW2*( DW1(3)*W32 + DW3(3)*W12 - 2.D0*DW2(3)*W13 ) |
| 2586 | C |
| 2587 | C jj(0)=js(0)+jt(0)+j4(0) |
| 2588 | C jj(1)=js(1)+jt(1)+j4(1) |
| 2589 | C jj(2)=js(2)+jt(2)+j4(2) |
| 2590 | C jj(3)=js(3)+jt(3)+j4(3) |
| 2591 | |
| 2592 | JJ(0)=J4(0) |
| 2593 | JJ(1)=J4(1) |
| 2594 | JJ(2)=J4(2) |
| 2595 | JJ(3)=J4(3) |
| 2596 | C |
| 2597 | JQ=(JJ(0)*Q(0)-JJ(1)*Q(1)-JJ(2)*Q(2)-JJ(3)*Q(3))/MW2 |
| 2598 | C |
| 2599 | |
| 2600 | JWWW(1) = DCMPLX( (JJ(0)-JQ*Q(0))*DW ) |
| 2601 | JWWW(2) = DCMPLX( (JJ(1)-JQ*Q(1))*DW ) |
| 2602 | JWWW(3) = DCMPLX( (JJ(2)-JQ*Q(2))*DW ) |
| 2603 | JWWW(4) = DCMPLX( (JJ(3)-JQ*Q(3))*DW ) |
| 2604 | C |
| 2605 | RETURN |
| 2606 | END |
| 2607 | |
| 2608 | C |
| 2609 | C ---------------------------------------------------------------------- |
| 2610 | C |
| 2611 | SUBROUTINE MOM2CX(ESUM,MASS1,MASS2,COSTH1,PHI1 , P1,P2) |
| 2612 | C |
| 2613 | C This subroutine sets up two four-momenta in the two particle rest |
| 2614 | C frame. |
| 2615 | C |
| 2616 | C INPUT: |
| 2617 | C real ESUM : energy sum of particle 1 and 2 |
| 2618 | C real MASS1 : mass of particle 1 |
| 2619 | C real MASS2 : mass of particle 2 |
| 2620 | C real COSTH1 : cos(theta) of particle 1 |
| 2621 | C real PHI1 : azimuthal angle of particle 1 |
| 2622 | C |
| 2623 | C OUTPUT: |
| 2624 | C real P1(0:3) : four-momentum of particle 1 |
| 2625 | C real P2(0:3) : four-momentum of particle 2 |
| 2626 | C |
| 2627 | #if defined(CERNLIB_IMPNONE) |
| 2628 | IMPLICIT NONE |
| 2629 | #endif |
| 2630 | REAL*8 P1(0:3),P2(0:3), |
| 2631 | & ESUM,MASS1,MASS2,COSTH1,PHI1,MD2,ED,PP,SINTH1 |
| 2632 | C |
| 2633 | MD2=(MASS1-MASS2)*(MASS1+MASS2) |
| 2634 | ED=MD2/ESUM |
| 2635 | IF (MASS1*MASS2.EQ.0.) THEN |
| 2636 | PP=(ESUM-ABS(ED))*0.5D0 |
| 2637 | C |
| 2638 | ELSE |
| 2639 | PP=SQRT((MD2/ESUM)**2-2.0D0*(MASS1**2+MASS2**2)+ESUM**2)*0.5D0 |
| 2640 | ENDIF |
| 2641 | SINTH1=SQRT((1.0D0-COSTH1)*(1.0D0+COSTH1)) |
| 2642 | C |
| 2643 | P1(0) = MAX((ESUM+ED)*0.5D0,0.D0) |
| 2644 | P1(1) = PP*SINTH1*COS(PHI1) |
| 2645 | P1(2) = PP*SINTH1*SIN(PHI1) |
| 2646 | P1(3) = PP*COSTH1 |
| 2647 | C |
| 2648 | P2(0) = MAX((ESUM-ED)*0.5D0,0.D0) |
| 2649 | P2(1) = -P1(1) |
| 2650 | P2(2) = -P1(2) |
| 2651 | P2(3) = -P1(3) |
| 2652 | C |
| 2653 | RETURN |
| 2654 | END |
| 2655 | C ********************************************************************** |
| 2656 | C |
| 2657 | SUBROUTINE MOMNTX(ENERGY,MASS,COSTH,PHI , P) |
| 2658 | C |
| 2659 | C This subroutine sets up a four-momentum from the four inputs. |
| 2660 | C |
| 2661 | C INPUT: |
| 2662 | C real ENERGY : energy |
| 2663 | C real MASS : mass |
| 2664 | C real COSTH : cos(theta) |
| 2665 | C real PHI : azimuthal angle |
| 2666 | C |
| 2667 | C OUTPUT: |
| 2668 | C real P(0:3) : four-momentum |
| 2669 | C |
| 2670 | #if defined(CERNLIB_IMPNONE) |
| 2671 | IMPLICIT NONE |
| 2672 | #endif |
| 2673 | REAL*8 P(0:3),ENERGY,MASS,COSTH,PHI,PP,SINTH |
| 2674 | C |
| 2675 | P(0) = ENERGY |
| 2676 | IF (ENERGY.EQ.MASS) THEN |
| 2677 | P(1) = 0. |
| 2678 | P(2) = 0. |
| 2679 | P(3) = 0. |
| 2680 | ELSE |
| 2681 | PP=SQRT((ENERGY-MASS)*(ENERGY+MASS)) |
| 2682 | SINTH=SQRT((1.-COSTH)*(1.+COSTH)) |
| 2683 | P(3) = PP*COSTH |
| 2684 | IF (PHI.EQ.0.) THEN |
| 2685 | P(1) = PP*SINTH |
| 2686 | P(2) = 0. |
| 2687 | ELSE |
| 2688 | P(1) = PP*SINTH*COS(PHI) |
| 2689 | P(2) = PP*SINTH*SIN(PHI) |
| 2690 | ENDIF |
| 2691 | ENDIF |
| 2692 | RETURN |
| 2693 | END |
| 2694 | C |
| 2695 | C |
| 2696 | C |
| 2697 | C Subroutine returns the desired fermion or |
| 2698 | C anti-fermion anti-spinor. ie., <f| |
| 2699 | C A replacement for the HELAS routine OXXXXX |
| 2700 | C |
| 2701 | C Adam Duff, 1992 August 31 |
| 2702 | C <duff@phenom.physics.wisc.edu> |
| 2703 | C |
| 2704 | SUBROUTINE OXXXXX(P,FMASS,NHEL,NSF,FO) |
| 2705 | C |
| 2706 | C P IN: FOUR VECTOR MOMENTUM |
| 2707 | C FMASS IN: FERMION MASS |
| 2708 | C NHEL IN: ANTI-SPINOR HELICITY, -1 OR 1 |
| 2709 | C NSF IN: -1=ANTIFERMION, 1=FERMION |
| 2710 | C FO OUT: FERMION WAVEFUNCTION |
| 2711 | C |
| 2712 | C declare input/output variables |
| 2713 | C |
| 2714 | #if defined(CERNLIB_IMPNONE) |
| 2715 | IMPLICIT NONE |
| 2716 | #endif |
| 2717 | COMPLEX*16 FO(6) |
| 2718 | INTEGER*4 NHEL, NSF |
| 2719 | REAL*8 P(0:3), FMASS |
| 2720 | C |
| 2721 | C declare local variables |
| 2722 | C |
| 2723 | REAL*8 RXZERO, RXONE, RXTWO |
| 2724 | PARAMETER( RXZERO=0.0D0, RXONE=1.0D0, RXTWO=2.0D0 ) |
| 2725 | REAL*8 PLAT, PABS, OMEGAP, OMEGAM, RS2PA, SPAZ |
| 2726 | COMPLEX*16 CXZERO |
| 2727 | LOGICAL FIRST |
| 2728 | SAVE CXZERO,FIRST |
| 2729 | DATA FIRST/.TRUE./ |
| 2730 | C |
| 2731 | C Fix compilation with g77 |
| 2732 | IF(FIRST) THEN |
| 2733 | FIRST=.FALSE. |
| 2734 | CXZERO=DCMPLX( RXZERO, RXZERO ) |
| 2735 | ENDIF |
| 2736 | C |
| 2737 | C define kinematic parameters |
| 2738 | C |
| 2739 | FO(5) = DCMPLX( P(0), P(3) ) * NSF |
| 2740 | FO(6) = DCMPLX( P(1), P(2) ) * NSF |
| 2741 | PLAT = SQRT( P(1)**2 + P(2)**2 ) |
| 2742 | PABS = SQRT( P(1)**2 + P(2)**2 + P(3)**2 ) |
| 2743 | OMEGAP = SQRT( P(0) + PABS ) |
| 2744 | C |
| 2745 | C do massive fermion case |
| 2746 | C |
| 2747 | IF ( FMASS .NE. RXZERO ) THEN |
| 2748 | OMEGAM = FMASS / OMEGAP |
| 2749 | IF ( NSF .EQ. 1 ) THEN |
| 2750 | IF ( NHEL .EQ. 1 ) THEN |
| 2751 | IF ( P(3) .GE. RXZERO ) THEN |
| 2752 | IF ( PLAT .EQ. RXZERO ) THEN |
| 2753 | FO(1) = DCMPLX( OMEGAP, RXZERO ) |
| 2754 | FO(2) = CXZERO |
| 2755 | FO(3) = DCMPLX( OMEGAM, RXZERO ) |
| 2756 | FO(4) = CXZERO |
| 2757 | ELSE |
| 2758 | RS2PA = RXONE / SQRT( RXTWO * PABS ) |
| 2759 | SPAZ = SQRT( PABS + P(3) ) |
| 2760 | FO(1) = OMEGAP * RS2PA |
| 2761 | & * DCMPLX( SPAZ, RXZERO ) |
| 2762 | FO(2) = OMEGAP * RS2PA / SPAZ |
| 2763 | & * DCMPLX( P(1), -P(2) ) |
| 2764 | FO(3) = OMEGAM * RS2PA |
| 2765 | & * DCMPLX( SPAZ, RXZERO ) |
| 2766 | FO(4) = OMEGAM * RS2PA / SPAZ |
| 2767 | & * DCMPLX( P(1), -P(2) ) |
| 2768 | END IF |
| 2769 | ELSE |
| 2770 | IF ( PLAT .EQ. RXZERO ) THEN |
| 2771 | FO(1) = CXZERO |
| 2772 | FO(2) = DCMPLX( OMEGAP, RXZERO ) |
| 2773 | FO(3) = CXZERO |
| 2774 | FO(4) = DCMPLX( OMEGAM, RXZERO ) |
| 2775 | ELSE |
| 2776 | RS2PA = RXONE / SQRT( RXTWO * PABS ) |
| 2777 | SPAZ = SQRT( PABS - P(3) ) |
| 2778 | FO(1) = OMEGAP * RS2PA / SPAZ |
| 2779 | & * DCMPLX( PLAT, RXZERO ) |
| 2780 | FO(2) = OMEGAP * RS2PA * SPAZ / PLAT |
| 2781 | & * DCMPLX( P(1), -P(2) ) |
| 2782 | FO(3) = OMEGAM * RS2PA / SPAZ |
| 2783 | & * DCMPLX( PLAT, RXZERO ) |
| 2784 | FO(4) = OMEGAM * RS2PA * SPAZ / PLAT |
| 2785 | & * DCMPLX( P(1), -P(2) ) |
| 2786 | END IF |
| 2787 | END IF |
| 2788 | ELSE IF ( NHEL .EQ. -1 ) THEN |
| 2789 | IF ( P(3) .GE. RXZERO ) THEN |
| 2790 | IF ( PLAT .EQ. RXZERO ) THEN |
| 2791 | FO(1) = CXZERO |
| 2792 | FO(2) = DCMPLX( OMEGAM, RXZERO ) |
| 2793 | FO(3) = CXZERO |
| 2794 | FO(4) = DCMPLX( OMEGAP, RXZERO ) |
| 2795 | ELSE |
| 2796 | RS2PA = RXONE / SQRT( RXTWO * PABS ) |
| 2797 | SPAZ = SQRT( PABS + P(3) ) |
| 2798 | FO(1) = OMEGAM * RS2PA / SPAZ |
| 2799 | & * DCMPLX( -P(1), -P(2) ) |
| 2800 | FO(2) = OMEGAM * RS2PA |
| 2801 | & * DCMPLX( SPAZ, RXZERO ) |
| 2802 | FO(3) = OMEGAP * RS2PA / SPAZ |
| 2803 | & * DCMPLX( -P(1), -P(2) ) |
| 2804 | FO(4) = OMEGAP * RS2PA |
| 2805 | & * DCMPLX( SPAZ, RXZERO ) |
| 2806 | END IF |
| 2807 | ELSE |
| 2808 | IF ( PLAT .EQ. RXZERO ) THEN |
| 2809 | FO(1) = DCMPLX( -OMEGAM, RXZERO ) |
| 2810 | FO(2) = CXZERO |
| 2811 | FO(3) = DCMPLX( -OMEGAP, RXZERO ) |
| 2812 | FO(4) = CXZERO |
| 2813 | ELSE |
| 2814 | RS2PA = RXONE / SQRT( RXTWO * PABS ) |
| 2815 | SPAZ = SQRT( PABS - P(3) ) |
| 2816 | FO(1) = OMEGAM * RS2PA * SPAZ / PLAT |
| 2817 | & * DCMPLX( -P(1), -P(2) ) |
| 2818 | FO(2) = OMEGAM * RS2PA / SPAZ |
| 2819 | & * DCMPLX( PLAT, RXZERO ) |
| 2820 | FO(3) = OMEGAP * RS2PA * SPAZ / PLAT |
| 2821 | & * DCMPLX( -P(1), -P(2) ) |
| 2822 | FO(4) = OMEGAP * RS2PA / SPAZ |
| 2823 | & * DCMPLX( PLAT, RXZERO ) |
| 2824 | END IF |
| 2825 | END IF |
| 2826 | ELSE |
| 2827 | STOP 'OXXXXX: FERMION HELICITY MUST BE +1,-1' |
| 2828 | END IF |
| 2829 | ELSE IF ( NSF .EQ. -1 ) THEN |
| 2830 | IF ( NHEL .EQ. 1 ) THEN |
| 2831 | IF ( P(3) .GE. RXZERO ) THEN |
| 2832 | IF ( PLAT .EQ. RXZERO ) THEN |
| 2833 | FO(1) = CXZERO |
| 2834 | FO(2) = DCMPLX( OMEGAM, RXZERO ) |
| 2835 | FO(3) = CXZERO |
| 2836 | FO(4) = DCMPLX( -OMEGAP, RXZERO ) |
| 2837 | ELSE |
| 2838 | RS2PA = RXONE / SQRT( RXTWO * PABS ) |
| 2839 | SPAZ = SQRT( PABS + P(3) ) |
| 2840 | FO(1) = OMEGAM * RS2PA / SPAZ |
| 2841 | & * DCMPLX( -P(1), -P(2) ) |
| 2842 | FO(2) = OMEGAM * RS2PA |
| 2843 | & * DCMPLX( SPAZ, RXZERO ) |
| 2844 | FO(3) = -OMEGAP * RS2PA / SPAZ |
| 2845 | & * DCMPLX( -P(1), -P(2) ) |
| 2846 | FO(4) = -OMEGAP * RS2PA |
| 2847 | & * DCMPLX( SPAZ, RXZERO ) |
| 2848 | END IF |
| 2849 | ELSE |
| 2850 | IF ( PLAT .EQ. RXZERO ) THEN |
| 2851 | FO(1) = DCMPLX( -OMEGAM, RXZERO ) |
| 2852 | FO(2) = CXZERO |
| 2853 | FO(3) = DCMPLX( OMEGAP, RXZERO ) |
| 2854 | FO(4) = CXZERO |
| 2855 | ELSE |
| 2856 | RS2PA = RXONE / SQRT( RXTWO * PABS ) |
| 2857 | SPAZ = SQRT( PABS - P(3) ) |
| 2858 | FO(1) = OMEGAM * RS2PA * SPAZ / PLAT |
| 2859 | & * DCMPLX( -P(1), -P(2) ) |
| 2860 | FO(2) = OMEGAM * RS2PA / SPAZ |
| 2861 | & * DCMPLX( PLAT, RXZERO ) |
| 2862 | FO(3) = -OMEGAP * RS2PA * SPAZ / PLAT |
| 2863 | & * DCMPLX( -P(1), -P(2) ) |
| 2864 | FO(4) = -OMEGAP * RS2PA / SPAZ |
| 2865 | & * DCMPLX( PLAT, RXZERO ) |
| 2866 | END IF |
| 2867 | END IF |
| 2868 | ELSE IF ( NHEL .EQ. -1 ) THEN |
| 2869 | IF ( P(3) .GE. RXZERO ) THEN |
| 2870 | IF ( PLAT .EQ. RXZERO ) THEN |
| 2871 | FO(1) = DCMPLX( -OMEGAP, RXZERO ) |
| 2872 | FO(2) = CXZERO |
| 2873 | FO(3) = DCMPLX( OMEGAM, RXZERO ) |
| 2874 | FO(4) = CXZERO |
| 2875 | ELSE |
| 2876 | RS2PA = RXONE / SQRT( RXTWO * PABS ) |
| 2877 | SPAZ = SQRT( PABS + P(3) ) |
| 2878 | FO(1) = -OMEGAP * RS2PA |
| 2879 | & * DCMPLX( SPAZ, RXZERO ) |
| 2880 | FO(2) = -OMEGAP * RS2PA / SPAZ |
| 2881 | & * DCMPLX( P(1), -P(2) ) |
| 2882 | FO(3) = OMEGAM * RS2PA |
| 2883 | & * DCMPLX( SPAZ, RXZERO ) |
| 2884 | FO(4) = OMEGAM * RS2PA / SPAZ |
| 2885 | & * DCMPLX( P(1), -P(2) ) |
| 2886 | END IF |
| 2887 | ELSE |
| 2888 | IF ( PLAT .EQ. RXZERO ) THEN |
| 2889 | FO(1) = CXZERO |
| 2890 | FO(2) = DCMPLX( -OMEGAP, RXZERO ) |
| 2891 | FO(3) = CXZERO |
| 2892 | FO(4) = DCMPLX( OMEGAM, RXZERO ) |
| 2893 | ELSE |
| 2894 | RS2PA = RXONE / SQRT( RXTWO * PABS ) |
| 2895 | SPAZ = SQRT( PABS - P(3) ) |
| 2896 | FO(1) = -OMEGAP * RS2PA / SPAZ |
| 2897 | & * DCMPLX( PLAT, RXZERO ) |
| 2898 | FO(2) = -OMEGAP * RS2PA * SPAZ / PLAT |
| 2899 | & * DCMPLX( P(1), -P(2) ) |
| 2900 | FO(3) = OMEGAM * RS2PA / SPAZ |
| 2901 | & * DCMPLX( PLAT, RXZERO ) |
| 2902 | FO(4) = OMEGAM * RS2PA * SPAZ / PLAT |
| 2903 | & * DCMPLX( P(1), -P(2) ) |
| 2904 | END IF |
| 2905 | END IF |
| 2906 | ELSE |
| 2907 | STOP 'OXXXXX: FERMION HELICITY MUST BE +1,-1' |
| 2908 | END IF |
| 2909 | ELSE |
| 2910 | STOP 'OXXXXX: FERMION TYPE MUST BE +1,-1' |
| 2911 | END IF |
| 2912 | C |
| 2913 | C do massless case |
| 2914 | C |
| 2915 | ELSE |
| 2916 | IF ( NSF .EQ. 1 ) THEN |
| 2917 | IF ( NHEL .EQ. 1 ) THEN |
| 2918 | IF ( P(3) .GE. RXZERO ) THEN |
| 2919 | IF ( PLAT .EQ. RXZERO ) THEN |
| 2920 | FO(1) = DCMPLX( OMEGAP, RXZERO ) |
| 2921 | FO(2) = CXZERO |
| 2922 | FO(3) = CXZERO |
| 2923 | FO(4) = CXZERO |
| 2924 | ELSE |
| 2925 | SPAZ = SQRT( PABS + P(3) ) |
| 2926 | FO(1) = DCMPLX( SPAZ, RXZERO ) |
| 2927 | FO(2) = RXONE / SPAZ |
| 2928 | & * DCMPLX( P(1), -P(2) ) |
| 2929 | FO(3) = CXZERO |
| 2930 | FO(4) = CXZERO |
| 2931 | END IF |
| 2932 | ELSE |
| 2933 | IF ( PLAT .EQ. RXZERO ) THEN |
| 2934 | FO(1) = CXZERO |
| 2935 | FO(2) = DCMPLX( OMEGAP, RXZERO ) |
| 2936 | FO(3) = CXZERO |
| 2937 | FO(4) = CXZERO |
| 2938 | ELSE |
| 2939 | SPAZ = SQRT( PABS - P(3) ) |
| 2940 | FO(1) = RXONE / SPAZ |
| 2941 | & * DCMPLX( PLAT, RXZERO ) |
| 2942 | FO(2) = SPAZ / PLAT |
| 2943 | & * DCMPLX( P(1), -P(2) ) |
| 2944 | FO(3) = CXZERO |
| 2945 | FO(4) = CXZERO |
| 2946 | END IF |
| 2947 | END IF |
| 2948 | ELSE IF ( NHEL .EQ. -1 ) THEN |
| 2949 | IF ( P(3) .GE. RXZERO ) THEN |
| 2950 | IF ( PLAT .EQ. RXZERO ) THEN |
| 2951 | FO(1) = CXZERO |
| 2952 | FO(2) = CXZERO |
| 2953 | FO(3) = CXZERO |
| 2954 | FO(4) = DCMPLX( OMEGAP, RXZERO ) |
| 2955 | ELSE |
| 2956 | SPAZ = SQRT( PABS + P(3) ) |
| 2957 | FO(1) = CXZERO |
| 2958 | FO(2) = CXZERO |
| 2959 | FO(3) = RXONE / SPAZ |
| 2960 | & * DCMPLX( -P(1), -P(2) ) |
| 2961 | FO(4) = DCMPLX( SPAZ, RXZERO ) |
| 2962 | END IF |
| 2963 | ELSE |
| 2964 | IF ( PLAT .EQ. RXZERO ) THEN |
| 2965 | FO(1) = CXZERO |
| 2966 | FO(2) = CXZERO |
| 2967 | FO(3) = DCMPLX( -OMEGAP, RXZERO ) |
| 2968 | FO(4) = CXZERO |
| 2969 | ELSE |
| 2970 | SPAZ = SQRT( PABS - P(3) ) |
| 2971 | FO(1) = CXZERO |
| 2972 | FO(2) = CXZERO |
| 2973 | FO(3) = SPAZ / PLAT |
| 2974 | & * DCMPLX( -P(1), -P(2) ) |
| 2975 | FO(4) = RXONE / SPAZ |
| 2976 | & * DCMPLX( PLAT, RXZERO ) |
| 2977 | END IF |
| 2978 | END IF |
| 2979 | ELSE |
| 2980 | STOP 'OXXXXX: FERMION HELICITY MUST BE +1,-1' |
| 2981 | END IF |
| 2982 | ELSE IF ( NSF .EQ. -1 ) THEN |
| 2983 | IF ( NHEL .EQ. 1 ) THEN |
| 2984 | IF ( P(3) .GE. RXZERO ) THEN |
| 2985 | IF ( PLAT .EQ. RXZERO ) THEN |
| 2986 | FO(1) = CXZERO |
| 2987 | FO(2) = CXZERO |
| 2988 | FO(3) = CXZERO |
| 2989 | FO(4) = DCMPLX( -OMEGAP, RXZERO ) |
| 2990 | ELSE |
| 2991 | SPAZ = SQRT( PABS + P(3) ) |
| 2992 | FO(1) = CXZERO |
| 2993 | FO(2) = CXZERO |
| 2994 | FO(3) = -RXONE / SPAZ |
| 2995 | & * DCMPLX( -P(1), -P(2) ) |
| 2996 | FO(4) = DCMPLX( -SPAZ, RXZERO ) |
| 2997 | END IF |
| 2998 | ELSE |
| 2999 | IF ( PLAT .EQ. RXZERO ) THEN |
| 3000 | FO(1) = CXZERO |
| 3001 | FO(2) = CXZERO |
| 3002 | FO(3) = DCMPLX( OMEGAP, RXZERO ) |
| 3003 | FO(4) = CXZERO |
| 3004 | ELSE |
| 3005 | SPAZ = SQRT( PABS - P(3) ) |
| 3006 | FO(1) = CXZERO |
| 3007 | FO(2) = CXZERO |
| 3008 | FO(3) = -SPAZ / PLAT |
| 3009 | & * DCMPLX( -P(1), -P(2) ) |
| 3010 | FO(4) = -RXONE / SPAZ |
| 3011 | & * DCMPLX( PLAT, RXZERO ) |
| 3012 | END IF |
| 3013 | END IF |
| 3014 | ELSE IF ( NHEL .EQ. -1 ) THEN |
| 3015 | IF ( P(3) .GE. RXZERO ) THEN |
| 3016 | IF ( PLAT .EQ. RXZERO ) THEN |
| 3017 | FO(1) = DCMPLX( -OMEGAP, RXZERO ) |
| 3018 | FO(2) = CXZERO |
| 3019 | FO(3) = CXZERO |
| 3020 | FO(4) = CXZERO |
| 3021 | ELSE |
| 3022 | SPAZ = SQRT( PABS + P(3) ) |
| 3023 | FO(1) = DCMPLX( -SPAZ, RXZERO ) |
| 3024 | FO(2) = -RXONE / SPAZ |
| 3025 | & * DCMPLX( P(1), -P(2) ) |
| 3026 | FO(3) = CXZERO |
| 3027 | FO(4) = CXZERO |
| 3028 | END IF |
| 3029 | ELSE |
| 3030 | IF ( PLAT .EQ. RXZERO ) THEN |
| 3031 | FO(1) = CXZERO |
| 3032 | FO(2) = DCMPLX( -OMEGAP, RXZERO ) |
| 3033 | FO(3) = CXZERO |
| 3034 | FO(4) = CXZERO |
| 3035 | ELSE |
| 3036 | SPAZ = SQRT( PABS - P(3) ) |
| 3037 | FO(1) = -RXONE / SPAZ |
| 3038 | & * DCMPLX( PLAT, RXZERO ) |
| 3039 | FO(2) = -SPAZ / PLAT |
| 3040 | & * DCMPLX( P(1), -P(2) ) |
| 3041 | FO(3) = CXZERO |
| 3042 | FO(4) = CXZERO |
| 3043 | END IF |
| 3044 | END IF |
| 3045 | ELSE |
| 3046 | STOP 'OXXXXX: FERMION HELICITY MUST BE +1,-1' |
| 3047 | END IF |
| 3048 | ELSE |
| 3049 | STOP 'OXXXXX: FERMION TYPE MUST BE +1,-1' |
| 3050 | END IF |
| 3051 | END IF |
| 3052 | C |
| 3053 | C done |
| 3054 | C |
| 3055 | RETURN |
| 3056 | END |
| 3057 | C |
| 3058 | C ---------------------------------------------------------------------- |
| 3059 | C |
| 3060 | SUBROUTINE ROTXXX(P,Q , PROT) |
| 3061 | C |
| 3062 | C this subroutine performs the spacial rotation of a four-momentum. |
| 3063 | C the momentum p is assumed to be given in the frame where the spacial |
| 3064 | C component of q points the positive z-axis. prot is the momentum p |
| 3065 | C rotated to the frame where q is given. |
| 3066 | C |
| 3067 | C input: |
| 3068 | C real p(0:3) : four-momentum p in q(1)=q(2)=0 frame |
| 3069 | C real q(0:3) : four-momentum q in the rotated frame |
| 3070 | C |
| 3071 | C output: |
| 3072 | C real prot(0:3) : four-momentum p in the rotated frame |
| 3073 | C |
| 3074 | #if defined(CERNLIB_IMPNONE) |
| 3075 | IMPLICIT NONE |
| 3076 | #endif |
| 3077 | REAL*8 P(0:3),Q(0:3),PROT(0:3),QT2,QT,PSGN,QQ,P1 |
| 3078 | C |
| 3079 | REAL*8 RXZERO, RXONE |
| 3080 | PARAMETER( RXZERO=0.0D0, RXONE=1.0D0 ) |
| 3081 | C |
| 3082 | PROT(0) = P(0) |
| 3083 | C |
| 3084 | QT2=Q(1)**2+Q(2)**2 |
| 3085 | C |
| 3086 | IF ( QT2 .EQ. RXZERO ) THEN |
| 3087 | IF ( Q(3) .EQ. RXZERO ) THEN |
| 3088 | PROT(1) = P(1) |
| 3089 | PROT(2) = P(2) |
| 3090 | PROT(3) = P(3) |
| 3091 | ELSE |
| 3092 | PSGN=DSIGN(RXONE,Q(3)) |
| 3093 | PROT(1) = P(1)*PSGN |
| 3094 | PROT(2) = P(2)*PSGN |
| 3095 | PROT(3) = P(3)*PSGN |
| 3096 | ENDIF |
| 3097 | ELSE |
| 3098 | QQ=SQRT(QT2+Q(3)**2) |
| 3099 | QT=SQRT(QT2) |
| 3100 | P1=P(1) |
| 3101 | PROT(1) = Q(1)*Q(3)/QQ/QT*P1 -Q(2)/QT*P(2) +Q(1)/QQ*P(3) |
| 3102 | PROT(2) = Q(2)*Q(3)/QQ/QT*P1 +Q(1)/QT*P(2) +Q(2)/QQ*P(3) |
| 3103 | PROT(3) = -QT/QQ*P1 +Q(3)/QQ*P(3) |
| 3104 | ENDIF |
| 3105 | C |
| 3106 | RETURN |
| 3107 | END |
| 3108 | C ====================================================================== |
| 3109 | C |
| 3110 | SUBROUTINE SSSSXX(S1,S2,S3,S4,G , VERTEX) |
| 3111 | C |
| 3112 | C This subroutine computes an amplitude of the four-scalar coupling. |
| 3113 | C |
| 3114 | C INPUT: |
| 3115 | C complex S1(3) : first scalar S1 |
| 3116 | C complex S2(3) : second scalar S2 |
| 3117 | C complex S3(3) : third scalar S3 |
| 3118 | C complex S4(3) : fourth scalar S4 |
| 3119 | C real G : coupling constant GHHHH |
| 3120 | C |
| 3121 | C OUTPUT: |
| 3122 | C complex VERTEX : amplitude Gamma(S1,S2,S3,S4) |
| 3123 | C |
| 3124 | #if defined(CERNLIB_IMPNONE) |
| 3125 | IMPLICIT NONE |
| 3126 | #endif |
| 3127 | COMPLEX*16 S1(3),S2(3),S3(3),S4(3),VERTEX |
| 3128 | REAL*8 G |
| 3129 | C |
| 3130 | VERTEX = G*S1(1)*S2(1)*S3(1)*S4(1) |
| 3131 | C |
| 3132 | RETURN |
| 3133 | END |
| 3134 | C |
| 3135 | C ====================================================================== |
| 3136 | C |
| 3137 | SUBROUTINE SSSXXX(S1,S2,S3,G , VERTEX) |
| 3138 | C |
| 3139 | C This subroutine computes an amplitude of the three-scalar coupling. |
| 3140 | C |
| 3141 | C INPUT: |
| 3142 | C complex S1(3) : first scalar S1 |
| 3143 | C complex S2(3) : second scalar S2 |
| 3144 | C complex S3(3) : third scalar S3 |
| 3145 | C real G : coupling constant GHHH |
| 3146 | C |
| 3147 | C OUTPUT: |
| 3148 | C complex VERTEX : amplitude Gamma(S1,S2,S3) |
| 3149 | C |
| 3150 | #if defined(CERNLIB_IMPNONE) |
| 3151 | IMPLICIT NONE |
| 3152 | #endif |
| 3153 | COMPLEX*16 S1(3),S2(3),S3(3),VERTEX |
| 3154 | REAL*8 G |
| 3155 | C |
| 3156 | VERTEX = G*S1(1)*S2(1)*S3(1) |
| 3157 | C |
| 3158 | RETURN |
| 3159 | END |
| 3160 | C |
| 3161 | C |
| 3162 | C ---------------------------------------------------------------------- |
| 3163 | C |
| 3164 | SUBROUTINE SXXXXX(P,NSS , SC) |
| 3165 | C |
| 3166 | C This subroutine computes a complex SCALAR wavefunction. |
| 3167 | C |
| 3168 | C INPUT: |
| 3169 | C real P(0:3) : four-momentum of scalar boson |
| 3170 | C integer NSS = -1 or 1 : +1 for final, -1 for initial |
| 3171 | C |
| 3172 | C OUTPUT: |
| 3173 | C complex SC(3) : scalar wavefunction S |
| 3174 | C |
| 3175 | #if defined(CERNLIB_IMPNONE) |
| 3176 | IMPLICIT NONE |
| 3177 | #endif |
| 3178 | COMPLEX*16 SC(3) |
| 3179 | REAL*8 P(0:3) |
| 3180 | INTEGER NSS |
| 3181 | C |
| 3182 | SC(1) = DCMPLX( 1.0 ) |
| 3183 | SC(2) = DCMPLX(P(0),P(3))*NSS |
| 3184 | SC(3) = DCMPLX(P(1),P(2))*NSS |
| 3185 | C |
| 3186 | RETURN |
| 3187 | END |
| 3188 | C |
| 3189 | C ====================================================================== |
| 3190 | C |
| 3191 | SUBROUTINE VSSXXX(VC,S1,S2,G , VERTEX) |
| 3192 | C |
| 3193 | C this subroutine computes an amplitude from the vector-scalar-scalar |
| 3194 | C coupling. the coupling is absent in the minimal sm in unitary gauge. |
| 3195 | C |
| 3196 | C complex vc(6) : input vector v |
| 3197 | C complex s1(3) : first scalar s1 |
| 3198 | C complex s2(3) : second scalar s2 |
| 3199 | C complex g : coupling constant (s1 charge) |
| 3200 | C |
| 3201 | C examples of the coupling constant g for susy particles are as follows: |
| 3202 | C ----------------------------------------------------------- |
| 3203 | C | s1 | (q,i3) of s1 || v=a | v=z | v=w | |
| 3204 | C ----------------------------------------------------------- |
| 3205 | C | nu~_l | ( 0 , +1/2) || --- | gzn(1) | gwf(1) | |
| 3206 | C | e~_l | ( -1 , -1/2) || gal(1) | gzl(1) | gwf(1) | |
| 3207 | C | u~_l | (+2/3 , +1/2) || gau(1) | gzu(1) | gwf(1) | |
| 3208 | C | d~_l | (-1/3 , -1/2) || gad(1) | gzd(1) | gwf(1) | |
| 3209 | C ----------------------------------------------------------- |
| 3210 | C | e~_r-bar | ( +1 , 0 ) || -gal(2) | -gzl(2) | -gwf(2) | |
| 3211 | C | u~_r-bar | (-2/3 , 0 ) || -gau(2) | -gzu(2) | -gwf(2) | |
| 3212 | C | d~_r-bar | (+1/3 , 0 ) || -gad(2) | -gzd(2) | -gwf(2) | |
| 3213 | C ----------------------------------------------------------- |
| 3214 | C where the s1 charge is defined by the flowing-out quantum number. |
| 3215 | C |
| 3216 | C output: |
| 3217 | C complex vertex : amplitude gamma(v,s1,s2) |
| 3218 | C |
| 3219 | #if defined(CERNLIB_IMPNONE) |
| 3220 | IMPLICIT NONE |
| 3221 | #endif |
| 3222 | COMPLEX*16 VC(6),S1(3),S2(3),VERTEX,G |
| 3223 | REAL*8 P(0:3) |
| 3224 | C |
| 3225 | P(0)=DBLE( S1(2)-S2(2)) |
| 3226 | P(1)=DBLE( S1(3)-S2(3)) |
| 3227 | P(2)=DIMAG(S1(3)-S2(3)) |
| 3228 | P(3)=DIMAG(S1(2)-S2(2)) |
| 3229 | C |
| 3230 | VERTEX = G*S1(1)*S2(1) |
| 3231 | & *(VC(1)*P(0)-VC(2)*P(1)-VC(3)*P(2)-VC(4)*P(3)) |
| 3232 | C |
| 3233 | RETURN |
| 3234 | END |
| 3235 | C |
| 3236 | SUBROUTINE VVSSXX(V1,V2,S1,S2,G , VERTEX) |
| 3237 | C |
| 3238 | C This subroutine computes an amplitude of the vector-vector-scalar- |
| 3239 | C scalar coupling. |
| 3240 | C |
| 3241 | C INPUT: |
| 3242 | C complex V1(6) : first vector V1 |
| 3243 | C complex V2(6) : second vector V2 |
| 3244 | C complex S1(3) : first scalar S1 |
| 3245 | C complex S2(3) : second scalar S2 |
| 3246 | C real G : coupling constant GVVHH |
| 3247 | C |
| 3248 | C OUTPUT: |
| 3249 | C complex VERTEX : amplitude Gamma(V1,V2,S1,S2) |
| 3250 | C |
| 3251 | #if defined(CERNLIB_IMPNONE) |
| 3252 | IMPLICIT NONE |
| 3253 | #endif |
| 3254 | COMPLEX*16 V1(6),V2(6),S1(3),S2(3),VERTEX |
| 3255 | REAL*8 G |
| 3256 | C |
| 3257 | VERTEX = G*S1(1)*S2(1) |
| 3258 | & *(V1(1)*V2(1)-V1(2)*V2(2)-V1(3)*V2(3)-V1(4)*V2(4)) |
| 3259 | C |
| 3260 | RETURN |
| 3261 | END |
| 3262 | C |
| 3263 | C |
| 3264 | C ====================================================================== |
| 3265 | C |
| 3266 | SUBROUTINE VVSXXX(V1,V2,SC,G , VERTEX) |
| 3267 | C |
| 3268 | C this subroutine computes an amplitude of the vector-vector-scalar |
| 3269 | C coupling. |
| 3270 | C |
| 3271 | C input: |
| 3272 | C complex v1(6) : first vector v1 |
| 3273 | C complex v2(6) : second vector v2 |
| 3274 | C complex sc(3) : input scalar s |
| 3275 | C real g : coupling constant gvvh |
| 3276 | C |
| 3277 | C output: |
| 3278 | C complex vertex : amplitude gamma(v1,v2,s) |
| 3279 | C |
| 3280 | #if defined(CERNLIB_IMPNONE) |
| 3281 | IMPLICIT NONE |
| 3282 | #endif |
| 3283 | COMPLEX*16 V1(6),V2(6),SC(3),VERTEX |
| 3284 | REAL*8 G |
| 3285 | C |
| 3286 | VERTEX = G*SC(1)*(V1(1)*V2(1)-V1(2)*V2(2)-V1(3)*V2(3)-V1(4)*V2(4)) |
| 3287 | C |
| 3288 | RETURN |
| 3289 | END |
| 3290 | C |
| 3291 | C ====================================================================== |
| 3292 | C |
| 3293 | SUBROUTINE VVVXXX(WM,WP,W3,G , VERTEX) |
| 3294 | C |
| 3295 | C this subroutine computes an amplitude of the three-point coupling of |
| 3296 | C the gauge bosons. |
| 3297 | C |
| 3298 | C input: |
| 3299 | C complex wm(6) : vector flow-out w- |
| 3300 | C complex wp(6) : vector flow-out w+ |
| 3301 | C complex w3(6) : vector j3 or a or z |
| 3302 | C real g : coupling constant gw or gwwa or gwwz |
| 3303 | C |
| 3304 | C output: |
| 3305 | C complex vertex : amplitude gamma(wm,wp,w3) |
| 3306 | C |
| 3307 | #if defined(CERNLIB_IMPNONE) |
| 3308 | IMPLICIT NONE |
| 3309 | #endif |
| 3310 | COMPLEX*16 WM(6),WP(6),W3(6),VERTEX, |
| 3311 | & XV1,XV2,XV3,V12,V23,V31,P12,P13,P21,P23,P31,P32 |
| 3312 | REAL*8 PWM(0:3),PWP(0:3),PW3(0:3),G |
| 3313 | C |
| 3314 | REAL*8 RXZERO, RTENTH |
| 3315 | PARAMETER( RXZERO=0.0D0, RTENTH=0.1D0 ) |
| 3316 | C |
| 3317 | PWM(0)=DBLE( WM(5)) |
| 3318 | PWM(1)=DBLE( WM(6)) |
| 3319 | PWM(2)=DIMAG(WM(6)) |
| 3320 | PWM(3)=DIMAG(WM(5)) |
| 3321 | PWP(0)=DBLE( WP(5)) |
| 3322 | PWP(1)=DBLE( WP(6)) |
| 3323 | PWP(2)=DIMAG(WP(6)) |
| 3324 | PWP(3)=DIMAG(WP(5)) |
| 3325 | PW3(0)=DBLE( W3(5)) |
| 3326 | PW3(1)=DBLE( W3(6)) |
| 3327 | PW3(2)=DIMAG(W3(6)) |
| 3328 | PW3(3)=DIMAG(W3(5)) |
| 3329 | C |
| 3330 | V12=WM(1)*WP(1)-WM(2)*WP(2)-WM(3)*WP(3)-WM(4)*WP(4) |
| 3331 | V23=WP(1)*W3(1)-WP(2)*W3(2)-WP(3)*W3(3)-WP(4)*W3(4) |
| 3332 | V31=W3(1)*WM(1)-W3(2)*WM(2)-W3(3)*WM(3)-W3(4)*WM(4) |
| 3333 | XV1=RXZERO |
| 3334 | XV2=RXZERO |
| 3335 | XV3=RXZERO |
| 3336 | IF ( ABS(WM(1)) .NE. RXZERO ) THEN |
| 3337 | IF (ABS(WM(1)).GE.MAX(ABS(WM(2)),ABS(WM(3)),ABS(WM(4))) |
| 3338 | $ *RTENTH) |
| 3339 | & XV1=PWM(0)/WM(1) |
| 3340 | ENDIF |
| 3341 | IF ( ABS(WP(1)) .NE. RXZERO) THEN |
| 3342 | IF (ABS(WP(1)).GE.MAX(ABS(WP(2)),ABS(WP(3)),ABS(WP(4))) |
| 3343 | $ *RTENTH) |
| 3344 | & XV2=PWP(0)/WP(1) |
| 3345 | ENDIF |
| 3346 | IF ( ABS(W3(1)) .NE. RXZERO) THEN |
| 3347 | IF ( ABS(W3(1)).GE.MAX(ABS(W3(2)),ABS(W3(3)),ABS(W3(4))) |
| 3348 | $ *RTENTH) |
| 3349 | & XV3=PW3(0)/W3(1) |
| 3350 | ENDIF |
| 3351 | P12= (PWM(0)-XV1*WM(1))*WP(1)-(PWM(1)-XV1*WM(2))*WP(2) |
| 3352 | & -(PWM(2)-XV1*WM(3))*WP(3)-(PWM(3)-XV1*WM(4))*WP(4) |
| 3353 | P13= (PWM(0)-XV1*WM(1))*W3(1)-(PWM(1)-XV1*WM(2))*W3(2) |
| 3354 | & -(PWM(2)-XV1*WM(3))*W3(3)-(PWM(3)-XV1*WM(4))*W3(4) |
| 3355 | P21= (PWP(0)-XV2*WP(1))*WM(1)-(PWP(1)-XV2*WP(2))*WM(2) |
| 3356 | & -(PWP(2)-XV2*WP(3))*WM(3)-(PWP(3)-XV2*WP(4))*WM(4) |
| 3357 | P23= (PWP(0)-XV2*WP(1))*W3(1)-(PWP(1)-XV2*WP(2))*W3(2) |
| 3358 | & -(PWP(2)-XV2*WP(3))*W3(3)-(PWP(3)-XV2*WP(4))*W3(4) |
| 3359 | P31= (PW3(0)-XV3*W3(1))*WM(1)-(PW3(1)-XV3*W3(2))*WM(2) |
| 3360 | & -(PW3(2)-XV3*W3(3))*WM(3)-(PW3(3)-XV3*W3(4))*WM(4) |
| 3361 | P32= (PW3(0)-XV3*W3(1))*WP(1)-(PW3(1)-XV3*W3(2))*WP(2) |
| 3362 | & -(PW3(2)-XV3*W3(3))*WP(3)-(PW3(3)-XV3*W3(4))*WP(4) |
| 3363 | C |
| 3364 | VERTEX = -(V12*(P13-P23)+V23*(P21-P31)+V31*(P32-P12))*G |
| 3365 | C |
| 3366 | RETURN |
| 3367 | END |
| 3368 | C |
| 3369 | C |
| 3370 | C Subroutine returns the value of evaluated |
| 3371 | C helicity basis boson polarisation wavefunction. |
| 3372 | C Replaces the HELAS routine VXXXXX |
| 3373 | C |
| 3374 | C Adam Duff, 1992 September 3 |
| 3375 | C <duff@phenom.physics.wisc.edu> |
| 3376 | C |
| 3377 | SUBROUTINE VXXXXX(P,VMASS,NHEL,NSV,VC) |
| 3378 | C |
| 3379 | C P IN: BOSON FOUR MOMENTUM |
| 3380 | C VMASS IN: BOSON MASS |
| 3381 | C NHEL IN: BOSON HELICITY |
| 3382 | C NSV IN: INCOMING (-1) OR OUTGOING (+1) |
| 3383 | C VC OUT: BOSON WAVEFUNCTION |
| 3384 | C |
| 3385 | C declare input/output variables |
| 3386 | C |
| 3387 | #if defined(CERNLIB_IMPNONE) |
| 3388 | IMPLICIT NONE |
| 3389 | #endif |
| 3390 | COMPLEX*16 VC(6) |
| 3391 | INTEGER*4 NHEL, NSV |
| 3392 | REAL*8 P(0:3), VMASS |
| 3393 | C |
| 3394 | C declare local variables |
| 3395 | C |
| 3396 | REAL*8 RXZERO, RXONE, RXTWO |
| 3397 | PARAMETER( RXZERO=0.0D0, RXONE=1.0D0, RXTWO=2.0D0 ) |
| 3398 | REAL*8 PLAT, PABS, RS2, RPLAT, RPABS, RDEN |
| 3399 | COMPLEX*16 CXZERO |
| 3400 | LOGICAL FIRST |
| 3401 | SAVE CXZERO,FIRST |
| 3402 | DATA FIRST/.TRUE./ |
| 3403 | C |
| 3404 | C Fix compilation with g77 |
| 3405 | IF(FIRST) THEN |
| 3406 | FIRST=.FALSE. |
| 3407 | CXZERO=DCMPLX( RXZERO, RXZERO ) |
| 3408 | ENDIF |
| 3409 | C |
| 3410 | C define internal/external momenta |
| 3411 | C |
| 3412 | IF ( NSV**2 .NE. 1 ) THEN |
| 3413 | STOP 'VXXXXX: NSV IS NOT ONE OF -1, +1' |
| 3414 | END IF |
| 3415 | C |
| 3416 | RS2 = SQRT( RXONE / RXTWO ) |
| 3417 | VC(5) = DCMPLX( P(0), P(3) ) * NSV |
| 3418 | VC(6) = DCMPLX( P(1), P(2) ) * NSV |
| 3419 | PLAT = SQRT( P(1)**2 + P(2)**2 ) |
| 3420 | PABS = SQRT( P(1)**2 + P(2)**2 + P(3)**2 ) |
| 3421 | C |
| 3422 | C calculate polarisation four vectors |
| 3423 | C |
| 3424 | IF ( NHEL**2 .EQ. 1 ) THEN |
| 3425 | IF ( (PABS .EQ. RXZERO) .OR. (PLAT .EQ. RXZERO) ) THEN |
| 3426 | VC(1) = CXZERO |
| 3427 | VC(2) = DCMPLX( -NHEL * RS2 * DSIGN( RXONE, P(3) ), RXZERO ) |
| 3428 | VC(3) = DCMPLX( RXZERO, NSV * RS2 ) |
| 3429 | VC(4) = CXZERO |
| 3430 | ELSE |
| 3431 | RPLAT = RXONE / PLAT |
| 3432 | RPABS = RXONE / PABS |
| 3433 | VC(1) = CXZERO |
| 3434 | VC(2) = DCMPLX( -NHEL * RS2 * RPABS * RPLAT * P(1) * P(3), |
| 3435 | & -NSV * RS2 * RPLAT * P(2) ) |
| 3436 | VC(3) = DCMPLX( -NHEL * RS2 * RPABS * RPLAT * P(2) * P(3), |
| 3437 | & NSV * RS2 * RPLAT * P(1) ) |
| 3438 | VC(4) = DCMPLX( NHEL * RS2 * RPABS * PLAT, |
| 3439 | & RXZERO ) |
| 3440 | END IF |
| 3441 | ELSE IF ( NHEL .EQ. 0 ) THEN |
| 3442 | IF ( VMASS .GT. RXZERO ) THEN |
| 3443 | IF ( PABS .EQ. RXZERO ) THEN |
| 3444 | VC(1) = CXZERO |
| 3445 | VC(2) = CXZERO |
| 3446 | VC(3) = CXZERO |
| 3447 | VC(4) = DCMPLX( RXONE, RXZERO ) |
| 3448 | ELSE |
| 3449 | RDEN = P(0) / ( VMASS * PABS ) |
| 3450 | VC(1) = DCMPLX( PABS / VMASS, RXZERO ) |
| 3451 | VC(2) = DCMPLX( RDEN * P(1), RXZERO ) |
| 3452 | VC(3) = DCMPLX( RDEN * P(2), RXZERO ) |
| 3453 | VC(4) = DCMPLX( RDEN * P(3), RXZERO ) |
| 3454 | END IF |
| 3455 | ELSE |
| 3456 | STOP 'VXXXXX: NHEL = 0 IS ONLY FOR MASSIVE BOSONS' |
| 3457 | END IF |
| 3458 | ELSE IF ( NHEL .EQ. 4 ) THEN |
| 3459 | IF ( VMASS .GT. RXZERO ) THEN |
| 3460 | RDEN = RXONE / VMASS |
| 3461 | VC(1) = DCMPLX( RDEN * P(0), RXZERO ) |
| 3462 | VC(2) = DCMPLX( RDEN * P(1), RXZERO ) |
| 3463 | VC(3) = DCMPLX( RDEN * P(2), RXZERO ) |
| 3464 | VC(4) = DCMPLX( RDEN * P(3), RXZERO ) |
| 3465 | ELSEIF (VMASS .EQ. RXZERO) THEN |
| 3466 | RDEN = RXONE / P(0) |
| 3467 | VC(1) = DCMPLX( RDEN * P(0), RXZERO ) |
| 3468 | VC(2) = DCMPLX( RDEN * P(1), RXZERO ) |
| 3469 | VC(3) = DCMPLX( RDEN * P(2), RXZERO ) |
| 3470 | VC(4) = DCMPLX( RDEN * P(3), RXZERO ) |
| 3471 | ELSE |
| 3472 | STOP 'VXXXXX: NHEL = 4 IS ONLY FOR M>=0' |
| 3473 | END IF |
| 3474 | ELSE |
| 3475 | STOP 'VXXXXX: NHEL IS NOT ONE OF -1, 0, 1 OR 4' |
| 3476 | END IF |
| 3477 | C |
| 3478 | C done |
| 3479 | C |
| 3480 | RETURN |
| 3481 | END |
| 3482 | C |
| 3483 | C ---------------------------------------------------------------------- |
| 3484 | C |
| 3485 | SUBROUTINE W3W3XX(WM,W31,WP,W32,G31,G32,WMASS,WWIDTH , VERTEX) |
| 3486 | C |
| 3487 | C this subroutine computes an amplitude of the four-point coupling of |
| 3488 | C the w-, w+ and two w3/z/a. the amplitude includes the contributions |
| 3489 | C of w exchange diagrams. the internal w propagator is given in unitary |
| 3490 | C gauge. if one sets wmass=0.0, then the gggg vertex is given (see sect |
| 3491 | C 2.9.1 of the manual). |
| 3492 | C |
| 3493 | C input: |
| 3494 | C complex wm(0:3) : flow-out w- wm |
| 3495 | C complex w31(0:3) : first w3/z/a w31 |
| 3496 | C complex wp(0:3) : flow-out w+ wp |
| 3497 | C complex w32(0:3) : second w3/z/a w32 |
| 3498 | C real g31 : coupling of w31 with w-/w+ |
| 3499 | C real g32 : coupling of w32 with w-/w+ |
| 3500 | C (see the table below) |
| 3501 | C real wmass : mass of w |
| 3502 | C real wwidth : width of w |
| 3503 | C |
| 3504 | C the possible sets of the inputs are as follows: |
| 3505 | C ------------------------------------------- |
| 3506 | C | wm | w31 | wp | w32 | g31 | g32 | |
| 3507 | C ------------------------------------------- |
| 3508 | C | w- | w3 | w+ | w3 | gw | gw | |
| 3509 | C | w- | w3 | w+ | z | gw | gwwz | |
| 3510 | C | w- | w3 | w+ | a | gw | gwwa | |
| 3511 | C | w- | z | w+ | z | gwwz | gwwz | |
| 3512 | C | w- | z | w+ | a | gwwz | gwwa | |
| 3513 | C | w- | a | w+ | a | gwwa | gwwa | |
| 3514 | C ------------------------------------------- |
| 3515 | C where all the bosons are defined by the flowing-out quantum number. |
| 3516 | C |
| 3517 | C output: |
| 3518 | C complex vertex : amplitude gamma(wm,w31,wp,w32) |
| 3519 | C |
| 3520 | #if defined(CERNLIB_IMPNONE) |
| 3521 | IMPLICIT NONE |
| 3522 | #endif |
| 3523 | COMPLEX*16 WM(6),W31(6),WP(6),W32(6),VERTEX |
| 3524 | COMPLEX*16 DV1(0:3),DV2(0:3),DV3(0:3),DV4(0:3),DVERTX, |
| 3525 | & V12,V13,V14,V23,V24,V34 |
| 3526 | REAL*8 G31,G32,WMASS,WWIDTH |
| 3527 | C |
| 3528 | REAL*8 RXZERO, RXONE |
| 3529 | PARAMETER( RXZERO=0.0D0, RXONE=1.0D0 ) |
| 3530 | |
| 3531 | DV1(0)=DCMPLX(WM(1)) |
| 3532 | DV1(1)=DCMPLX(WM(2)) |
| 3533 | DV1(2)=DCMPLX(WM(3)) |
| 3534 | DV1(3)=DCMPLX(WM(4)) |
| 3535 | DV2(0)=DCMPLX(W31(1)) |
| 3536 | DV2(1)=DCMPLX(W31(2)) |
| 3537 | DV2(2)=DCMPLX(W31(3)) |
| 3538 | DV2(3)=DCMPLX(W31(4)) |
| 3539 | DV3(0)=DCMPLX(WP(1)) |
| 3540 | DV3(1)=DCMPLX(WP(2)) |
| 3541 | DV3(2)=DCMPLX(WP(3)) |
| 3542 | DV3(3)=DCMPLX(WP(4)) |
| 3543 | DV4(0)=DCMPLX(W32(1)) |
| 3544 | DV4(1)=DCMPLX(W32(2)) |
| 3545 | DV4(2)=DCMPLX(W32(3)) |
| 3546 | DV4(3)=DCMPLX(W32(4)) |
| 3547 | C |
| 3548 | IF ( DBLE(WMASS) .NE. RXZERO ) THEN |
| 3549 | C dm2inv = r_one / dmw2 |
| 3550 | C |
| 3551 | V12= DV1(0)*DV2(0)-DV1(1)*DV2(1)-DV1(2)*DV2(2)-DV1(3)*DV2(3) |
| 3552 | V13= DV1(0)*DV3(0)-DV1(1)*DV3(1)-DV1(2)*DV3(2)-DV1(3)*DV3(3) |
| 3553 | V14= DV1(0)*DV4(0)-DV1(1)*DV4(1)-DV1(2)*DV4(2)-DV1(3)*DV4(3) |
| 3554 | V23= DV2(0)*DV3(0)-DV2(1)*DV3(1)-DV2(2)*DV3(2)-DV2(3)*DV3(3) |
| 3555 | V24= DV2(0)*DV4(0)-DV2(1)*DV4(1)-DV2(2)*DV4(2)-DV2(3)*DV4(3) |
| 3556 | V34= DV3(0)*DV4(0)-DV3(1)*DV4(1)-DV3(2)*DV4(2)-DV3(3)*DV4(3) |
| 3557 | C |
| 3558 | DVERTX = V12*V34 +V14*V23 -2.D0*V13*V24 |
| 3559 | C |
| 3560 | VERTEX = DCMPLX( DVERTX ) * (G31*G32) |
| 3561 | C |
| 3562 | ELSE |
| 3563 | V12= DV1(0)*DV2(0)-DV1(1)*DV2(1)-DV1(2)*DV2(2)-DV1(3)*DV2(3) |
| 3564 | V13= DV1(0)*DV3(0)-DV1(1)*DV3(1)-DV1(2)*DV3(2)-DV1(3)*DV3(3) |
| 3565 | V14= DV1(0)*DV4(0)-DV1(1)*DV4(1)-DV1(2)*DV4(2)-DV1(3)*DV4(3) |
| 3566 | V23= DV2(0)*DV3(0)-DV2(1)*DV3(1)-DV2(2)*DV3(2)-DV2(3)*DV3(3) |
| 3567 | V24= DV2(0)*DV4(0)-DV2(1)*DV4(1)-DV2(2)*DV4(2)-DV2(3)*DV4(3) |
| 3568 | V34= DV3(0)*DV4(0)-DV3(1)*DV4(1)-DV3(2)*DV4(2)-DV3(3)*DV4(3) |
| 3569 | C |
| 3570 | |
| 3571 | DVERTX = V14*V23 -V13*V24 |
| 3572 | C |
| 3573 | VERTEX = DCMPLX( DVERTX ) * (G31*G32) |
| 3574 | END IF |
| 3575 | C |
| 3576 | RETURN |
| 3577 | END |
| 3578 | C |
| 3579 | C ====================================================================== |
| 3580 | C |
| 3581 | SUBROUTINE WWWWXX(WM1,WP1,WM2,WP2,GWWA,GWWZ,ZMASS,ZWIDTH , VERTEX) |
| 3582 | C |
| 3583 | C this subroutine computes an amplitude of the four-point w-/w+ |
| 3584 | C coupling, including the contributions of photon and z exchanges. the |
| 3585 | C photon propagator is given in feynman gauge and the z propagator is |
| 3586 | C given in unitary gauge. |
| 3587 | C |
| 3588 | C input: |
| 3589 | C complex wm1(0:3) : first flow-out w- wm1 |
| 3590 | C complex wp1(0:3) : first flow-out w+ wp1 |
| 3591 | C complex wm2(0:3) : second flow-out w- wm2 |
| 3592 | C complex wp2(0:3) : second flow-out w+ wp2 |
| 3593 | C real gwwa : coupling constant of w and a gwwa |
| 3594 | C real gwwz : coupling constant of w and z gwwz |
| 3595 | C real zmass : mass of z |
| 3596 | C real zwidth : width of z |
| 3597 | C |
| 3598 | C output: |
| 3599 | C complex vertex : amplitude gamma(wm1,wp1,wm2,wp2) |
| 3600 | C |
| 3601 | #if defined(CERNLIB_IMPNONE) |
| 3602 | IMPLICIT NONE |
| 3603 | #endif |
| 3604 | COMPLEX*16 WM1(6),WP1(6),WM2(6),WP2(6),VERTEX |
| 3605 | COMPLEX*16 DV1(0:3),DV2(0:3),DV3(0:3),DV4(0:3), |
| 3606 | & J12(0:3),J34(0:3),J14(0:3),J32(0:3),DVERTX, |
| 3607 | & SV1,SV2,SV3,SV4,TV1,TV2,TV3,TV4,DZS,DZT, |
| 3608 | & V12,V13,V14,V23,V24,V34,JS12,JS34,JS14,JS32,JS,JT |
| 3609 | REAL*8 PWM1(0:3),PWP1(0:3),PWM2(0:3),PWP2(0:3), |
| 3610 | & GWWA,GWWZ,ZMASS,ZWIDTH |
| 3611 | REAL*8 Q(0:3),K(0:3),DP1(0:3),DP2(0:3),DP3(0:3),DP4(0:3), |
| 3612 | & DGWWA2,DGWWZ2,DGW2,DMZ,DWIDTH,S,T,DAS,DAT |
| 3613 | C |
| 3614 | REAL*8 RXZERO, RXONE, RXTWO |
| 3615 | PARAMETER( RXZERO=0.0D0, RXONE=1.0D0, RXTWO=2.0D0 ) |
| 3616 | C |
| 3617 | PWM1(0)=DBLE( WM1(5)) |
| 3618 | PWM1(1)=DBLE( WM1(6)) |
| 3619 | PWM1(2)=DIMAG(WM1(6)) |
| 3620 | PWM1(3)=DIMAG(WM1(5)) |
| 3621 | PWP1(0)=DBLE( WP1(5)) |
| 3622 | PWP1(1)=DBLE( WP1(6)) |
| 3623 | PWP1(2)=DIMAG(WP1(6)) |
| 3624 | PWP1(3)=DIMAG(WP1(5)) |
| 3625 | PWM2(0)=DBLE( WM2(5)) |
| 3626 | PWM2(1)=DBLE( WM2(6)) |
| 3627 | PWM2(2)=DIMAG(WM2(6)) |
| 3628 | PWM2(3)=DIMAG(WM2(5)) |
| 3629 | PWP2(0)=DBLE( WP2(5)) |
| 3630 | PWP2(1)=DBLE( WP2(6)) |
| 3631 | PWP2(2)=DIMAG(WP2(6)) |
| 3632 | PWP2(3)=DIMAG(WP2(5)) |
| 3633 | C |
| 3634 | DV1(0)=DCMPLX(WM1(1)) |
| 3635 | DV1(1)=DCMPLX(WM1(2)) |
| 3636 | DV1(2)=DCMPLX(WM1(3)) |
| 3637 | DV1(3)=DCMPLX(WM1(4)) |
| 3638 | DP1(0)=DBLE(PWM1(0)) |
| 3639 | DP1(1)=DBLE(PWM1(1)) |
| 3640 | DP1(2)=DBLE(PWM1(2)) |
| 3641 | DP1(3)=DBLE(PWM1(3)) |
| 3642 | DV2(0)=DCMPLX(WP1(1)) |
| 3643 | DV2(1)=DCMPLX(WP1(2)) |
| 3644 | DV2(2)=DCMPLX(WP1(3)) |
| 3645 | DV2(3)=DCMPLX(WP1(4)) |
| 3646 | DP2(0)=DBLE(PWP1(0)) |
| 3647 | DP2(1)=DBLE(PWP1(1)) |
| 3648 | DP2(2)=DBLE(PWP1(2)) |
| 3649 | DP2(3)=DBLE(PWP1(3)) |
| 3650 | DV3(0)=DCMPLX(WM2(1)) |
| 3651 | DV3(1)=DCMPLX(WM2(2)) |
| 3652 | DV3(2)=DCMPLX(WM2(3)) |
| 3653 | DV3(3)=DCMPLX(WM2(4)) |
| 3654 | DP3(0)=DBLE(PWM2(0)) |
| 3655 | DP3(1)=DBLE(PWM2(1)) |
| 3656 | DP3(2)=DBLE(PWM2(2)) |
| 3657 | DP3(3)=DBLE(PWM2(3)) |
| 3658 | DV4(0)=DCMPLX(WP2(1)) |
| 3659 | DV4(1)=DCMPLX(WP2(2)) |
| 3660 | DV4(2)=DCMPLX(WP2(3)) |
| 3661 | DV4(3)=DCMPLX(WP2(4)) |
| 3662 | DP4(0)=DBLE(PWP2(0)) |
| 3663 | DP4(1)=DBLE(PWP2(1)) |
| 3664 | DP4(2)=DBLE(PWP2(2)) |
| 3665 | DP4(3)=DBLE(PWP2(3)) |
| 3666 | DGWWA2=DBLE(GWWA)**2 |
| 3667 | DGWWZ2=DBLE(GWWZ)**2 |
| 3668 | DGW2 =DGWWA2+DGWWZ2 |
| 3669 | DMZ =DBLE(ZMASS) |
| 3670 | DWIDTH=DBLE(ZWIDTH) |
| 3671 | C |
| 3672 | V12= DV1(0)*DV2(0)-DV1(1)*DV2(1)-DV1(2)*DV2(2)-DV1(3)*DV2(3) |
| 3673 | V13= DV1(0)*DV3(0)-DV1(1)*DV3(1)-DV1(2)*DV3(2)-DV1(3)*DV3(3) |
| 3674 | V14= DV1(0)*DV4(0)-DV1(1)*DV4(1)-DV1(2)*DV4(2)-DV1(3)*DV4(3) |
| 3675 | V23= DV2(0)*DV3(0)-DV2(1)*DV3(1)-DV2(2)*DV3(2)-DV2(3)*DV3(3) |
| 3676 | V24= DV2(0)*DV4(0)-DV2(1)*DV4(1)-DV2(2)*DV4(2)-DV2(3)*DV4(3) |
| 3677 | V34= DV3(0)*DV4(0)-DV3(1)*DV4(1)-DV3(2)*DV4(2)-DV3(3)*DV4(3) |
| 3678 | C |
| 3679 | Q(0)=DP1(0)+DP2(0) |
| 3680 | Q(1)=DP1(1)+DP2(1) |
| 3681 | Q(2)=DP1(2)+DP2(2) |
| 3682 | Q(3)=DP1(3)+DP2(3) |
| 3683 | K(0)=DP1(0)+DP4(0) |
| 3684 | K(1)=DP1(1)+DP4(1) |
| 3685 | K(2)=DP1(2)+DP4(2) |
| 3686 | K(3)=DP1(3)+DP4(3) |
| 3687 | C |
| 3688 | S=Q(0)**2-Q(1)**2-Q(2)**2-Q(3)**2 |
| 3689 | T=K(0)**2-K(1)**2-K(2)**2-K(3)**2 |
| 3690 | C |
| 3691 | DAS=-RXONE/S |
| 3692 | DAT=-RXONE/T |
| 3693 | DZS=-RXONE/DCMPLX( S-DMZ**2 , DMAX1(DSIGN(DMZ*DWIDTH,S),RXZERO) ) |
| 3694 | DZT=-RXONE/DCMPLX( T-DMZ**2 , DMAX1(DSIGN(DMZ*DWIDTH,T),RXZERO) ) |
| 3695 | C |
| 3696 | SV1= (DP2(0)+Q(0))*DV1(0) -(DP2(1)+Q(1))*DV1(1) |
| 3697 | & -(DP2(2)+Q(2))*DV1(2) -(DP2(3)+Q(3))*DV1(3) |
| 3698 | SV2=-(DP1(0)+Q(0))*DV2(0) +(DP1(1)+Q(1))*DV2(1) |
| 3699 | & +(DP1(2)+Q(2))*DV2(2) +(DP1(3)+Q(3))*DV2(3) |
| 3700 | SV3= (DP4(0)-Q(0))*DV3(0) -(DP4(1)-Q(1))*DV3(1) |
| 3701 | & -(DP4(2)-Q(2))*DV3(2) -(DP4(3)-Q(3))*DV3(3) |
| 3702 | SV4=-(DP3(0)-Q(0))*DV4(0) +(DP3(1)-Q(1))*DV4(1) |
| 3703 | & +(DP3(2)-Q(2))*DV4(2) +(DP3(3)-Q(3))*DV4(3) |
| 3704 | C |
| 3705 | TV1= (DP4(0)+K(0))*DV1(0) -(DP4(1)+K(1))*DV1(1) |
| 3706 | & -(DP4(2)+K(2))*DV1(2) -(DP4(3)+K(3))*DV1(3) |
| 3707 | TV2=-(DP3(0)-K(0))*DV2(0) +(DP3(1)-K(1))*DV2(1) |
| 3708 | & +(DP3(2)-K(2))*DV2(2) +(DP3(3)-K(3))*DV2(3) |
| 3709 | TV3= (DP2(0)-K(0))*DV3(0) -(DP2(1)-K(1))*DV3(1) |
| 3710 | & -(DP2(2)-K(2))*DV3(2) -(DP2(3)-K(3))*DV3(3) |
| 3711 | TV4=-(DP1(0)+K(0))*DV4(0) +(DP1(1)+K(1))*DV4(1) |
| 3712 | & +(DP1(2)+K(2))*DV4(2) +(DP1(3)+K(3))*DV4(3) |
| 3713 | C |
| 3714 | J12(0)=(DP1(0)-DP2(0))*V12 +SV1*DV2(0) +SV2*DV1(0) |
| 3715 | J12(1)=(DP1(1)-DP2(1))*V12 +SV1*DV2(1) +SV2*DV1(1) |
| 3716 | J12(2)=(DP1(2)-DP2(2))*V12 +SV1*DV2(2) +SV2*DV1(2) |
| 3717 | J12(3)=(DP1(3)-DP2(3))*V12 +SV1*DV2(3) +SV2*DV1(3) |
| 3718 | J34(0)=(DP3(0)-DP4(0))*V34 +SV3*DV4(0) +SV4*DV3(0) |
| 3719 | J34(1)=(DP3(1)-DP4(1))*V34 +SV3*DV4(1) +SV4*DV3(1) |
| 3720 | J34(2)=(DP3(2)-DP4(2))*V34 +SV3*DV4(2) +SV4*DV3(2) |
| 3721 | J34(3)=(DP3(3)-DP4(3))*V34 +SV3*DV4(3) +SV4*DV3(3) |
| 3722 | C |
| 3723 | J14(0)=(DP1(0)-DP4(0))*V14 +TV1*DV4(0) +TV4*DV1(0) |
| 3724 | J14(1)=(DP1(1)-DP4(1))*V14 +TV1*DV4(1) +TV4*DV1(1) |
| 3725 | J14(2)=(DP1(2)-DP4(2))*V14 +TV1*DV4(2) +TV4*DV1(2) |
| 3726 | J14(3)=(DP1(3)-DP4(3))*V14 +TV1*DV4(3) +TV4*DV1(3) |
| 3727 | J32(0)=(DP3(0)-DP2(0))*V23 +TV3*DV2(0) +TV2*DV3(0) |
| 3728 | J32(1)=(DP3(1)-DP2(1))*V23 +TV3*DV2(1) +TV2*DV3(1) |
| 3729 | J32(2)=(DP3(2)-DP2(2))*V23 +TV3*DV2(2) +TV2*DV3(2) |
| 3730 | J32(3)=(DP3(3)-DP2(3))*V23 +TV3*DV2(3) +TV2*DV3(3) |
| 3731 | C |
| 3732 | JS12=Q(0)*J12(0)-Q(1)*J12(1)-Q(2)*J12(2)-Q(3)*J12(3) |
| 3733 | JS34=Q(0)*J34(0)-Q(1)*J34(1)-Q(2)*J34(2)-Q(3)*J34(3) |
| 3734 | JS14=K(0)*J14(0)-K(1)*J14(1)-K(2)*J14(2)-K(3)*J14(3) |
| 3735 | JS32=K(0)*J32(0)-K(1)*J32(1)-K(2)*J32(2)-K(3)*J32(3) |
| 3736 | C |
| 3737 | JS=J12(0)*J34(0)-J12(1)*J34(1)-J12(2)*J34(2)-J12(3)*J34(3) |
| 3738 | JT=J14(0)*J32(0)-J14(1)*J32(1)-J14(2)*J32(2)-J14(3)*J32(3) |
| 3739 | C |
| 3740 | DVERTX = (V12*V34 +V14*V23 -RXTWO*V13*V24)*DGW2 |
| 3741 | |
| 3742 | C & +(dzs*dgwwz2+das*dgwwa2)*js -dzs*dgwwz2*js12*js34/dmz**2 |
| 3743 | C & +(dzt*dgwwz2+dat*dgwwa2)*jt -dzt*dgwwz2*js14*js32/dmz**2 |
| 3744 | C |
| 3745 | VERTEX = -DCMPLX( DVERTX ) |
| 3746 | C |
| 3747 | RETURN |
| 3748 | END |
git://git.uio.no / u / mrichter / AliRoot.git / blame