00001
00002
00003
00004
00005
00006 SUBROUTINE JAKER(JAK)
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044
00045
00046
00047 COMMON / TAUBRA / GAMPRT(30),JLIST(30),NCHAN
00048
00049
00050
00051
00052 REAL CUMUL(30),RRR(1)
00053
00054 IF(NCHAN.LE.0.OR.NCHAN.GT.30) GOTO 902
00055 CALL RANMAR(RRR,1)
00056 SUM=0
00057 DO 20 I=1,NCHAN
00058 SUM=SUM+GAMPRT(I)
00059 20 CUMUL(I)=SUM
00060 DO 25 I=NCHAN,1,-1
00061 IF(RRR(1).LT.CUMUL(I)/CUMUL(NCHAN)) JI=I
00062 25 CONTINUE
00063 JAK=JLIST(JI)
00064 RETURN
00065 902 PRINT 9020
00066 9020 FORMAT(' ----- JAKER: WRONG NCHAN')
00067 STOP
00068 END
00069 SUBROUTINE DEKAY(KTO,HX)
00070
00071
00072
00073
00074
00075
00076
00077
00078
00079
00080
00081
00082
00083
00084
00085
00086
00087
00088
00089
00090
00091
00092
00093
00094
00095
00096
00097
00098
00099 REAL H(4)
00100 REAL*8 HX(4)
00101 COMMON / JAKI / JAK1,JAK2,JAKP,JAKM,KTOM
00102
00103
00104
00105 COMMON / IDFC / IDF
00106
00107 COMMON /TAUPOS/ NP1,NP2
00108 COMMON / TAUBMC / GAMPMC(30),GAMPER(30),NEVDEC(30)
00109 REAL*4 GAMPMC ,GAMPER
00110 PARAMETER (NMODE=15,NM1=0,NM2=1,NM3=8,NM4=2,NM5=1,NM6=3)
00111
00112
00113
00114 COMMON / TAUDCD /IDFFIN(9,NMODE),MULPIK(NMODE)
00115
00116 & ,NAMES
00117 CHARACTER NAMES(NMODE)*31
00118 COMMON / INOUT / INUT,IOUT
00119 REAL PDUM1(4),PDUM2(4),PDUM3(4),PDUM4(4),PDUM5(4),HDUM(4)
00120 REAL PDUMX(4,9)
00121 DATA IWARM/0/
00122 KTOM=KTO
00123
00124
00125
00126 IF(KTO.EQ.-1) THEN
00127
00128
00129
00130 NP1=3
00131 NP2=4
00132 KTOM=1
00133 IF (IWARM.EQ.1) X=5/(IWARM-1)
00134 IWARM=1
00135 WRITE(IOUT,7001) JAK1,JAK2
00136 NEVTOT=0
00137 NEV1=0
00138 NEV2=0
00139 IF(JAK1.NE.-1.OR.JAK2.NE.-1) THEN
00140 CALL DADMEL(-1,IDUM,HDUM,PDUM1,PDUM2,PDUM3,PDUM4,PDUM5)
00141 CALL DADMMU(-1,IDUM,HDUM,PDUM1,PDUM2,PDUM3,PDUM4,PDUM5)
00142 CALL DADMPI(-1,IDUM,PDUM,PDUM1,PDUM2)
00143 CALL DADMRO(-1,IDUM,HDUM,PDUM1,PDUM2,PDUM3,PDUM4)
00144 CALL DADMAA(-1,IDUM,HDUM,PDUM1,PDUM2,PDUM3,PDUM4,PDUM5,JDUM)
00145 CALL DADMKK(-1,IDUM,PDUM,PDUM1,PDUM2)
00146 CALL DADMKS(-1,IDUM,HDUM,PDUM1,PDUM2,PDUM3,PDUM4,JDUM)
00147 CALL DADNEW(-1,IDUM,HDUM,PDUM1,PDUM2,PDUMX,JDUM)
00148 ENDIF
00149 DO 21 I=1,30
00150 NEVDEC(I)=0
00151 GAMPMC(I)=0
00152 21 GAMPER(I)=0
00153 ELSEIF(KTO.EQ.1) THEN
00154
00155
00156 NEVTOT=NEVTOT+1
00157 IF(IWARM.EQ.0) GOTO 902
00158 ISGN= IDF/IABS(IDF)
00159
00160
00161
00162
00163 CALL TAURDF(KTO)
00164
00165 CALL DEKAY1(0,H,ISGN)
00166 ELSEIF(KTO.EQ.2) THEN
00167
00168
00169 NEVTOT=NEVTOT+1
00170 IF(IWARM.EQ.0) GOTO 902
00171 ISGN=-IDF/IABS(IDF)
00172
00173
00174
00175
00176 CALL TAURDF(KTO)
00177
00178 CALL DEKAY2(0,H,ISGN)
00179 ELSEIF(KTO.EQ.11) THEN
00180
00181
00182 NEV1=NEV1+1
00183 ISGN= IDF/IABS(IDF)
00184 CALL DEKAY1(1,H,ISGN)
00185 ELSEIF(KTO.EQ.12) THEN
00186
00187
00188 NEV2=NEV2+1
00189 ISGN=-IDF/IABS(IDF)
00190 CALL DEKAY2(1,H,ISGN)
00191 ELSEIF(KTO.EQ.100) THEN
00192
00193 IF(JAK1.NE.-1.OR.JAK2.NE.-1) THEN
00194 CALL DADMEL( 1,IDUM,HDUM,PDUM1,PDUM2,PDUM3,PDUM4,PDUM5)
00195 CALL DADMMU( 1,IDUM,HDUM,PDUM1,PDUM2,PDUM3,PDUM4,PDUM5)
00196 CALL DADMPI( 1,IDUM,PDUM,PDUM1,PDUM2)
00197 CALL DADMRO( 1,IDUM,HDUM,PDUM1,PDUM2,PDUM3,PDUM4)
00198 CALL DADMAA( 1,IDUM,HDUM,PDUM1,PDUM2,PDUM3,PDUM4,PDUM5,JDUM)
00199 CALL DADMKK( 1,IDUM,PDUM,PDUM1,PDUM2)
00200 CALL DADMKS( 1,IDUM,HDUM,PDUM1,PDUM2,PDUM3,PDUM4,JDUM)
00201 CALL DADNEW( 1,IDUM,HDUM,PDUM1,PDUM2,PDUMX,JDUM)
00202 WRITE(IOUT,7010) NEV1,NEV2,NEVTOT
00203 WRITE(IOUT,7011) (NEVDEC(I),GAMPMC(I),GAMPER(I),I= 1,7)
00204 WRITE(IOUT,7012)
00205 $ (NEVDEC(I),GAMPMC(I),GAMPER(I),NAMES(I-7),I=8,7+NMODE)
00206 WRITE(IOUT,7013)
00207 ENDIF
00208 ELSE
00209
00210 GOTO 910
00211 ENDIF
00212
00213 DO 78 K=1,4
00214 78 HX(K)=H(K)
00215 RETURN
00216 7001 FORMAT(///1X,15(5H*****)
00217 $ /,' *', 25X,'*****TAUOLA LIBRARY: VERSION 2.7 ******',9X,1H*,
00218 $ /,' *', 25X,'***********August 1995***************',9X,1H*,
00219 $ /,' *', 25X,'**AUTHORS: S.JADACH, Z.WAS*************',9X,1H*,
00220 $ /,' *', 25X,'**R. DECKER, M. JEZABEK, J.H.KUEHN*****',9X,1H*,
00221 $ /,' *', 25X,'**AVAILABLE FROM: WASM AT CERNVM ******',9X,1H*,
00222 $ /,' *', 25X,'***** PUBLISHED IN COMP. PHYS. COMM.***',9X,1H*,
00223 $ /,' *', 25X,' Physics initialization by CLEO collab ',9X,1H*,
00224 $ /,' *', 25X,' see Alain Weinstein www home page: ',9X,1H*,
00225 $ /,' *', 25X,'http://www.cithep.caltech.edu/~ajw/ ',9X,1H*,
00226 $ /,' *', 25X,'/korb_doc.html#files ',9X,1H*,
00227 $ /,' *', 25X,'*******CERN TH-6793 NOVEMBER 1992*****',9X,1H*,
00228 $ /,' *', 25X,'**5 or more pi dec.: precision limited ',9X,1H*,
00229 $ /,' *', 25X,'****DEKAY ROUTINE: INITIALIZATION******',9X,1H*,
00230 $ /,' *',I20 ,5X,'JAK1 = DECAY MODE TAU+ ',9X,1H*,
00231 $ /,' *',I20 ,5X,'JAK2 = DECAY MODE TAU- ',9X,1H*,
00232 $ /,1X,15(5H*****)/)
00233 7010 FORMAT(///1X,15(5H*****)
00234 $ /,' *', 25X,'*****TAUOLA LIBRARY: VERSION 2.7 ******',9X,1H*,
00235 $ /,' *', 25X,'***********August 1995***************',9X,1H*,
00236 $ /,' *', 25X,'**AUTHORS: S.JADACH, Z.WAS*************',9X,1H*,
00237 $ /,' *', 25X,'**R. DECKER, M. JEZABEK, J.H.KUEHN*****',9X,1H*,
00238 $ /,' *', 25X,'**AVAILABLE FROM: WASM AT CERNVM ******',9X,1H*,
00239 $ /,' *', 25X,'***** PUBLISHED IN COMP. PHYS. COMM.***',9X,1H*,
00240 $ /,' *', 25X,'*******CERN-TH-5856 SEPTEMBER 1990*****',9X,1H*,
00241 $ /,' *', 25X,'*******CERN-TH-6195 SEPTEMBER 1991*****',9X,1H*,
00242 $ /,' *', 25X,'*******CERN TH-6793 NOVEMBER 1992*****',9X,1H*,
00243 $ /,' *', 25X,'*****DEKAY ROUTINE: FINAL REPORT*******',9X,1H*,
00244 $ /,' *',I20 ,5X,'NEV1 = NO. OF TAU+ DECS. ACCEPTED ',9X,1H*,
00245 $ /,' *',I20 ,5X,'NEV2 = NO. OF TAU- DECS. ACCEPTED ',9X,1H*,
00246 $ /,' *',I20 ,5X,'NEVTOT = SUM ',9X,1H*,
00247 $ /,' *',' NOEVTS ',
00248 $ ' PART.WIDTH ERROR ROUTINE DECAY MODE ',9X,1H*)
00249 7011 FORMAT(1X,'*'
00250 $ ,I10,2F12.7 ,' DADMEL ELECTRON ',9X,1H*
00251 $ /,' *',I10,2F12.7 ,' DADMMU MUON ',9X,1H*
00252 $ /,' *',I10,2F12.7 ,' DADMPI PION ',9X,1H*
00253 $ /,' *',I10,2F12.7, ' DADMRO RHO (->2PI) ',9X,1H*
00254 $ /,' *',I10,2F12.7, ' DADMAA A1 (->3PI) ',9X,1H*
00255 $ /,' *',I10,2F12.7, ' DADMKK KAON ',9X,1H*
00256 $ /,' *',I10,2F12.7, ' DADMKS K* ',9X,1H*)
00257 7012 FORMAT(1X,'*'
00258 $ ,I10,2F12.7,A31 ,8X,1H*)
00259 7013 FORMAT(1X,'*'
00260 $ ,20X,'THE ERROR IS RELATIVE AND PART.WIDTH ',10X,1H*
00261 $ /,' *',20X,'IN UNITS GFERMI**2*MASS**5/192/PI**3 ',10X,1H*
00262 $ /,1X,15(5H*****)/)
00263 902 PRINT 9020
00264 9020 FORMAT(' ----- DEKAY: LACK OF INITIALISATION')
00265 STOP
00266 910 PRINT 9100
00267 9100 FORMAT(' ----- DEKAY: WRONG VALUE OF KTO ')
00268 STOP
00269 END
00270 SUBROUTINE DEKAY1(IMOD,HH,ISGN)
00271
00272
00273 COMMON / DECP4 / PP1(4),PP2(4),KF1,KF2
00274 COMMON / JAKI / JAK1,JAK2,JAKP,JAKM,KTOM
00275 COMMON / TAUBMC / GAMPMC(30),GAMPER(30),NEVDEC(30)
00276 REAL*4 GAMPMC ,GAMPER
00277 REAL HH(4)
00278 REAL HV(4),PNU(4),PPI(4)
00279 REAL PWB(4),PMU(4),PNM(4)
00280 REAL PRHO(4),PIC(4),PIZ(4)
00281 REAL PAA(4),PIM1(4),PIM2(4),PIPL(4)
00282 REAL PKK(4),PKS(4)
00283 REAL PNPI(4,9)
00284 REAL PHOT(4)
00285 REAL PDUM(4)
00286 DATA NEV,NPRIN/0,10/
00287 KTO=1
00288 IF(JAK1.EQ.-1) RETURN
00289 IMD=IMOD
00290 IF(IMD.EQ.0) THEN
00291
00292 JAK=JAK1
00293 IF(JAK1.EQ.0) CALL JAKER(JAK)
00294 IF(JAK.EQ.1) THEN
00295 CALL DADMEL(0, ISGN,HV,PNU,PWB,PMU,PNM,PHOT)
00296 ELSEIF(JAK.EQ.2) THEN
00297 CALL DADMMU(0, ISGN,HV,PNU,PWB,PMU,PNM,PHOT)
00298 ELSEIF(JAK.EQ.3) THEN
00299 CALL DADMPI(0, ISGN,HV,PPI,PNU)
00300 ELSEIF(JAK.EQ.4) THEN
00301 CALL DADMRO(0, ISGN,HV,PNU,PRHO,PIC,PIZ)
00302 ELSEIF(JAK.EQ.5) THEN
00303 CALL DADMAA(0, ISGN,HV,PNU,PAA,PIM1,PIM2,PIPL,JAA)
00304 ELSEIF(JAK.EQ.6) THEN
00305 CALL DADMKK(0, ISGN,HV,PKK,PNU)
00306 ELSEIF(JAK.EQ.7) THEN
00307 CALL DADMKS(0, ISGN,HV,PNU,PKS ,PKK,PPI,JKST)
00308 ELSE
00309 CALL DADNEW(0, ISGN,HV,PNU,PWB,PNPI,JAK-7)
00310 ENDIF
00311 DO 33 I=1,3
00312 33 HH(I)=HV(I)
00313 HH(4)=1.0
00314
00315 ELSEIF(IMD.EQ.1) THEN
00316
00317 NEV=NEV+1
00318 IF (JAK.LT.31) THEN
00319 NEVDEC(JAK)=NEVDEC(JAK)+1
00320 ENDIF
00321 DO 34 I=1,4
00322 34 PDUM(I)=.0
00323 IF(JAK.EQ.1) THEN
00324 CALL DWLUEL(1,ISGN,PNU,PWB,PMU,PNM)
00325 CALL DWRPH(KTOM,PHOT)
00326 DO 10 I=1,4
00327 10 PP1(I)=PMU(I)
00328
00329 ELSEIF(JAK.EQ.2) THEN
00330 CALL DWLUMU(1,ISGN,PNU,PWB,PMU,PNM)
00331 CALL DWRPH(KTOM,PHOT)
00332 DO 20 I=1,4
00333 20 PP1(I)=PMU(I)
00334
00335 ELSEIF(JAK.EQ.3) THEN
00336 CALL DWLUPI(1,ISGN,PPI,PNU)
00337 DO 30 I=1,4
00338 30 PP1(I)=PPI(I)
00339
00340 ELSEIF(JAK.EQ.4) THEN
00341 CALL DWLURO(1,ISGN,PNU,PRHO,PIC,PIZ)
00342 DO 40 I=1,4
00343 40 PP1(I)=PRHO(I)
00344
00345 ELSEIF(JAK.EQ.5) THEN
00346 CALL DWLUAA(1,ISGN,PNU,PAA,PIM1,PIM2,PIPL,JAA)
00347 DO 50 I=1,4
00348 50 PP1(I)=PAA(I)
00349 ELSEIF(JAK.EQ.6) THEN
00350 CALL DWLUKK(1,ISGN,PKK,PNU)
00351 DO 60 I=1,4
00352 60 PP1(I)=PKK(I)
00353 ELSEIF(JAK.EQ.7) THEN
00354 CALL DWLUKS(1,ISGN,PNU,PKS,PKK,PPI,JKST)
00355 DO 70 I=1,4
00356 70 PP1(I)=PKS(I)
00357 ELSE
00358
00359 CALL DWLNEW(1,ISGN,PNU,PWB,PNPI,JAK)
00360 DO 80 I=1,4
00361 80 PP1(I)=PWB(I)
00362 ENDIF
00363
00364 ENDIF
00365
00366 END
00367 SUBROUTINE DEKAY2(IMOD,HH,ISGN)
00368
00369
00370 COMMON / DECP4 / PP1(4),PP2(4),KF1,KF2
00371 COMMON / JAKI / JAK1,JAK2,JAKP,JAKM,KTOM
00372 COMMON / TAUBMC / GAMPMC(30),GAMPER(30),NEVDEC(30)
00373 REAL*4 GAMPMC ,GAMPER
00374 REAL HH(4)
00375 REAL HV(4),PNU(4),PPI(4)
00376 REAL PWB(4),PMU(4),PNM(4)
00377 REAL PRHO(4),PIC(4),PIZ(4)
00378 REAL PAA(4),PIM1(4),PIM2(4),PIPL(4)
00379 REAL PKK(4),PKS(4)
00380 REAL PNPI(4,9)
00381 REAL PHOT(4)
00382 REAL PDUM(4)
00383 DATA NEV,NPRIN/0,10/
00384 KTO=2
00385 IF(JAK2.EQ.-1) RETURN
00386 IMD=IMOD
00387 IF(IMD.EQ.0) THEN
00388
00389 JAK=JAK2
00390 IF(JAK2.EQ.0) CALL JAKER(JAK)
00391 IF(JAK.EQ.1) THEN
00392 CALL DADMEL(0, ISGN,HV,PNU,PWB,PMU,PNM,PHOT)
00393 ELSEIF(JAK.EQ.2) THEN
00394 CALL DADMMU(0, ISGN,HV,PNU,PWB,PMU,PNM,PHOT)
00395 ELSEIF(JAK.EQ.3) THEN
00396 CALL DADMPI(0, ISGN,HV,PPI,PNU)
00397 ELSEIF(JAK.EQ.4) THEN
00398 CALL DADMRO(0, ISGN,HV,PNU,PRHO,PIC,PIZ)
00399 ELSEIF(JAK.EQ.5) THEN
00400 CALL DADMAA(0, ISGN,HV,PNU,PAA,PIM1,PIM2,PIPL,JAA)
00401 ELSEIF(JAK.EQ.6) THEN
00402 CALL DADMKK(0, ISGN,HV,PKK,PNU)
00403 ELSEIF(JAK.EQ.7) THEN
00404 CALL DADMKS(0, ISGN,HV,PNU,PKS ,PKK,PPI,JKST)
00405 ELSE
00406 CALL DADNEW(0, ISGN,HV,PNU,PWB,PNPI,JAK-7)
00407 ENDIF
00408 DO 33 I=1,3
00409 33 HH(I)=HV(I)
00410 HH(4)=1.0
00411 ELSEIF(IMD.EQ.1) THEN
00412
00413 NEV=NEV+1
00414 IF (JAK.LT.31) THEN
00415 NEVDEC(JAK)=NEVDEC(JAK)+1
00416 ENDIF
00417 DO 34 I=1,4
00418 34 PDUM(I)=.0
00419 IF(JAK.EQ.1) THEN
00420 CALL DWLUEL(2,ISGN,PNU,PWB,PMU,PNM)
00421 CALL DWRPH(KTOM,PHOT)
00422 DO 10 I=1,4
00423 10 PP2(I)=PMU(I)
00424
00425 ELSEIF(JAK.EQ.2) THEN
00426 CALL DWLUMU(2,ISGN,PNU,PWB,PMU,PNM)
00427 CALL DWRPH(KTOM,PHOT)
00428 DO 20 I=1,4
00429 20 PP2(I)=PMU(I)
00430
00431 ELSEIF(JAK.EQ.3) THEN
00432 CALL DWLUPI(2,ISGN,PPI,PNU)
00433 DO 30 I=1,4
00434 30 PP2(I)=PPI(I)
00435
00436 ELSEIF(JAK.EQ.4) THEN
00437 CALL DWLURO(2,ISGN,PNU,PRHO,PIC,PIZ)
00438 DO 40 I=1,4
00439 40 PP2(I)=PRHO(I)
00440
00441 ELSEIF(JAK.EQ.5) THEN
00442 CALL DWLUAA(2,ISGN,PNU,PAA,PIM1,PIM2,PIPL,JAA)
00443 DO 50 I=1,4
00444 50 PP2(I)=PAA(I)
00445 ELSEIF(JAK.EQ.6) THEN
00446 CALL DWLUKK(2,ISGN,PKK,PNU)
00447 DO 60 I=1,4
00448 60 PP1(I)=PKK(I)
00449 ELSEIF(JAK.EQ.7) THEN
00450 CALL DWLUKS(2,ISGN,PNU,PKS,PKK,PPI,JKST)
00451 DO 70 I=1,4
00452 70 PP1(I)=PKS(I)
00453 ELSE
00454
00455 CALL DWLNEW(2,ISGN,PNU,PWB,PNPI,JAK)
00456 DO 80 I=1,4
00457 80 PP1(I)=PWB(I)
00458 ENDIF
00459
00460 ENDIF
00461
00462 END
00463 SUBROUTINE DEXAY(KTO,POL)
00464
00465
00466
00467
00468
00469
00470
00471
00472
00473
00474
00475
00476
00477 COMMON / TAUBMC / GAMPMC(30),GAMPER(30),NEVDEC(30)
00478 REAL*4 GAMPMC ,GAMPER
00479 COMMON / JAKI / JAK1,JAK2,JAKP,JAKM,KTOM
00480 COMMON / IDFC / IDFF
00481 COMMON /TAUPOS/ NP1,NP2
00482 PARAMETER (NMODE=15,NM1=0,NM2=1,NM3=8,NM4=2,NM5=1,NM6=3)
00483 COMMON / TAUDCD /IDFFIN(9,NMODE),MULPIK(NMODE)
00484 & ,NAMES
00485 CHARACTER NAMES(NMODE)*31
00486 COMMON / INOUT / INUT,IOUT
00487 REAL POL(4)
00488 REAL PDUM1(4),PDUM2(4),PDUM3(4),PDUM4(4),PDUM5(4)
00489 REAL PDUM(4)
00490 REAL PDUMI(4,9)
00491 DATA IWARM/0/
00492 KTOM=KTO
00493
00494 IF(KTO.EQ.-1) THEN
00495
00496
00497
00498
00499 NP1=3
00500 NP2=4
00501 IWARM=1
00502 WRITE(IOUT, 7001) JAK1,JAK2
00503 NEVTOT=0
00504 NEV1=0
00505 NEV2=0
00506 IF(JAK1.NE.-1.OR.JAK2.NE.-1) THEN
00507 CALL DEXEL(-1,IDUM,PDUM,PDUM1,PDUM2,PDUM3,PDUM4,PDUM5)
00508 CALL DEXMU(-1,IDUM,PDUM,PDUM1,PDUM2,PDUM3,PDUM4,PDUM5)
00509 CALL DEXPI(-1,IDUM,PDUM,PDUM1,PDUM2)
00510 CALL DEXRO(-1,IDUM,PDUM,PDUM1,PDUM2,PDUM3,PDUM4)
00511 CALL DEXAA(-1,IDUM,PDUM,PDUM1,PDUM2,PDUM3,PDUM4,PDUM5,IDUM)
00512 CALL DEXKK(-1,IDUM,PDUM,PDUM1,PDUM2)
00513 CALL DEXKS(-1,IDUM,PDUM,PDUM1,PDUM2,PDUM3,PDUM4,IDUM)
00514 CALL DEXNEW(-1,IDUM,PDUM,PDUM1,PDUM2,PDUMI,IDUM)
00515 ENDIF
00516 DO 21 I=1,30
00517 NEVDEC(I)=0
00518 GAMPMC(I)=0
00519 21 GAMPER(I)=0
00520 ELSEIF(KTO.EQ.1) THEN
00521
00522
00523 NEVTOT=NEVTOT+1
00524 NEV1=NEV1+1
00525 IF(IWARM.EQ.0) GOTO 902
00526 ISGN=IDFF/IABS(IDFF)
00527
00528 CALL DEXAY1(KTO,JAK1,JAKP,POL,ISGN)
00529 ELSEIF(KTO.EQ.2) THEN
00530
00531
00532 NEVTOT=NEVTOT+1
00533 NEV2=NEV2+1
00534 IF(IWARM.EQ.0) GOTO 902
00535 ISGN=-IDFF/IABS(IDFF)
00536
00537 CALL DEXAY1(KTO,JAK2,JAKM,POL,ISGN)
00538 ELSEIF(KTO.EQ.100) THEN
00539
00540 IF(JAK1.NE.-1.OR.JAK2.NE.-1) THEN
00541 CALL DEXEL( 1,IDUM,PDUM,PDUM1,PDUM2,PDUM3,PDUM4,PDUM5)
00542 CALL DEXMU( 1,IDUM,PDUM,PDUM1,PDUM2,PDUM3,PDUM4,PDUM5)
00543 CALL DEXPI( 1,IDUM,PDUM,PDUM1,PDUM2)
00544 CALL DEXRO( 1,IDUM,PDUM,PDUM1,PDUM2,PDUM3,PDUM4)
00545 CALL DEXAA( 1,IDUM,PDUM,PDUM1,PDUM2,PDUM3,PDUM4,PDUM5,IDUM)
00546 CALL DEXKK( 1,IDUM,PDUM,PDUM1,PDUM2)
00547 CALL DEXKS( 1,IDUM,PDUM,PDUM1,PDUM2,PDUM3,PDUM4,IDUM)
00548 CALL DEXNEW( 1,IDUM,PDUM,PDUM1,PDUM2,PDUMI,IDUM)
00549 WRITE(IOUT,7010) NEV1,NEV2,NEVTOT
00550 WRITE(IOUT,7011) (NEVDEC(I),GAMPMC(I),GAMPER(I),I= 1,7)
00551 WRITE(IOUT,7012)
00552 $ (NEVDEC(I),GAMPMC(I),GAMPER(I),NAMES(I-7),I=8,7+NMODE)
00553 WRITE(IOUT,7013)
00554 ENDIF
00555 ELSE
00556 GOTO 910
00557 ENDIF
00558 RETURN
00559 7001 FORMAT(///1X,15(5H*****)
00560 $ /,' *', 25X,'*****TAUOLA LIBRARY: VERSION 2.7 ******',9X,1H*,
00561 $ /,' *', 25X,'***********August 1995***************',9X,1H*,
00562 $ /,' *', 25X,'**AUTHORS: S.JADACH, Z.WAS*************',9X,1H*,
00563 $ /,' *', 25X,'**R. DECKER, M. JEZABEK, J.H.KUEHN*****',9X,1H*,
00564 $ /,' *', 25X,'**AVAILABLE FROM: WASM AT CERNVM ******',9X,1H*,
00565 $ /,' *', 25X,'***** PUBLISHED IN COMP. PHYS. COMM.***',9X,1H*,
00566 $ /,' *', 25X,' Physics initialization by CLEO collab ',9X,1H*,
00567 $ /,' *', 25X,' see Alain Weinstein www home page: ',9X,1H*,
00568 $ /,' *', 25X,'http://www.cithep.caltech.edu/~ajw/ ',9X,1H*,
00569 $ /,' *', 25X,'/korb_doc.html#files ',9X,1H*,
00570 $ /,' *', 25X,'*******CERN TH-6793 NOVEMBER 1992*****',9X,1H*,
00571 $ /,' *', 25X,'**5 or more pi dec.: precision limited ',9X,1H*,
00572 $ /,' *', 25X,'*******CERN-TH-6793 NOVEMBER 1992*****',9X,1H*,
00573 $ /,' *', 25X,'**5 or more pi dec.: precision limited ',9X,1H*,
00574 $ /,' *', 25X,'******DEXAY ROUTINE: INITIALIZATION****',9X,1H*
00575 $ /,' *',I20 ,5X,'JAK1 = DECAY MODE FERMION1 (TAU+) ',9X,1H*
00576 $ /,' *',I20 ,5X,'JAK2 = DECAY MODE FERMION2 (TAU-) ',9X,1H*
00577 $ /,1X,15(5H*****)/)
00578
00579
00580 7010 FORMAT(///1X,15(5H*****)
00581 $ /,' *', 25X,'*****TAUOLA LIBRARY: VERSION 2.7 ******',9X,1H*,
00582 $ /,' *', 25X,'***********August 1995***************',9X,1H*,
00583 $ /,' *', 25X,'**AUTHORS: S.JADACH, Z.WAS*************',9X,1H*,
00584 $ /,' *', 25X,'**R. DECKER, M. JEZABEK, J.H.KUEHN*****',9X,1H*,
00585 $ /,' *', 25X,'**AVAILABLE FROM: WASM AT CERNVM ******',9X,1H*,
00586 $ /,' *', 25X,'***** PUBLISHED IN COMP. PHYS. COMM.***',9X,1H*,
00587 $ /,' *', 25X,'*******CERN-TH-5856 SEPTEMBER 1990*****',9X,1H*,
00588 $ /,' *', 25X,'*******CERN-TH-6195 SEPTEMBER 1991*****',9X,1H*,
00589 $ /,' *', 25X,'*******CERN-TH-6793 NOVEMBER 1992*****',9X,1H*,
00590 $ /,' *', 25X,'******DEXAY ROUTINE: FINAL REPORT******',9X,1H*
00591 $ /,' *',I20 ,5X,'NEV1 = NO. OF TAU+ DECS. ACCEPTED ',9X,1H*
00592 $ /,' *',I20 ,5X,'NEV2 = NO. OF TAU- DECS. ACCEPTED ',9X,1H*
00593 $ /,' *',I20 ,5X,'NEVTOT = SUM ',9X,1H*
00594 $ /,' *',' NOEVTS ',
00595 $ ' PART.WIDTH ERROR ROUTINE DECAY MODE ',9X,1H*)
00596 7011 FORMAT(1X,'*'
00597 $ ,I10,2F12.7 ,' DADMEL ELECTRON ',9X,1H*
00598 $ /,' *',I10,2F12.7 ,' DADMMU MUON ',9X,1H*
00599 $ /,' *',I10,2F12.7 ,' DADMPI PION ',9X,1H*
00600 $ /,' *',I10,2F12.7, ' DADMRO RHO (->2PI) ',9X,1H*
00601 $ /,' *',I10,2F12.7, ' DADMAA A1 (->3PI) ',9X,1H*
00602 $ /,' *',I10,2F12.7, ' DADMKK KAON ',9X,1H*
00603 $ /,' *',I10,2F12.7, ' DADMKS K* ',9X,1H*)
00604 7012 FORMAT(1X,'*'
00605 $ ,I10,2F12.7,A31 ,8X,1H*)
00606 7013 FORMAT(1X,'*'
00607 $ ,20X,'THE ERROR IS RELATIVE AND PART.WIDTH ',10X,1H*
00608 $ /,' *',20X,'IN UNITS GFERMI**2*MASS**5/192/PI**3 ',10X,1H*
00609 $ /,1X,15(5H*****)/)
00610 902 WRITE(IOUT, 9020)
00611 9020 FORMAT(' ----- DEXAY: LACK OF INITIALISATION')
00612 STOP
00613 910 WRITE(IOUT, 9100)
00614 9100 FORMAT(' ----- DEXAY: WRONG VALUE OF KTO ')
00615 STOP
00616 END
00617 SUBROUTINE DEXAY1(KTO,JAKIN,JAK,POL,ISGN)
00618
00619
00620
00621
00622
00623 COMMON / TAUBMC / GAMPMC(30),GAMPER(30),NEVDEC(30)
00624 REAL*4 GAMPMC ,GAMPER
00625 COMMON / INOUT / INUT,IOUT
00626 REAL POL(4),POLAR(4)
00627 REAL PNU(4),PPI(4)
00628 REAL PRHO(4),PIC(4),PIZ(4)
00629 REAL PWB(4),PMU(4),PNM(4)
00630 REAL PAA(4),PIM1(4),PIM2(4),PIPL(4)
00631 REAL PKK(4),PKS(4)
00632 REAL PNPI(4,9)
00633 REAL PHOT(4)
00634 REAL PDUM(4)
00635
00636 IF(JAKIN.EQ.-1) RETURN
00637 DO 33 I=1,3
00638 33 POLAR(I)=POL(I)
00639 POLAR(4)=0.
00640 DO 34 I=1,4
00641 34 PDUM(I)=.0
00642 JAK=JAKIN
00643 IF(JAK.EQ.0) CALL JAKER(JAK)
00644
00645 IF(JAK.EQ.1) THEN
00646 CALL DEXEL(0, ISGN,POLAR,PNU,PWB,PMU,PNM,PHOT)
00647 CALL DWLUEL(KTO,ISGN,PNU,PWB,PMU,PNM)
00648 CALL DWRPH(KTO,PHOT )
00649 ELSEIF(JAK.EQ.2) THEN
00650 CALL DEXMU(0, ISGN,POLAR,PNU,PWB,PMU,PNM,PHOT)
00651 CALL DWLUMU(KTO,ISGN,PNU,PWB,PMU,PNM)
00652 CALL DWRPH(KTO,PHOT )
00653 ELSEIF(JAK.EQ.3) THEN
00654 CALL DEXPI(0, ISGN,POLAR,PPI,PNU)
00655 CALL DWLUPI(KTO,ISGN,PPI,PNU)
00656 ELSEIF(JAK.EQ.4) THEN
00657 CALL DEXRO(0, ISGN,POLAR,PNU,PRHO,PIC,PIZ)
00658 CALL DWLURO(KTO,ISGN,PNU,PRHO,PIC,PIZ)
00659 ELSEIF(JAK.EQ.5) THEN
00660 CALL DEXAA(0, ISGN,POLAR,PNU,PAA,PIM1,PIM2,PIPL,JAA)
00661 CALL DWLUAA(KTO,ISGN,PNU,PAA,PIM1,PIM2,PIPL,JAA)
00662 ELSEIF(JAK.EQ.6) THEN
00663 CALL DEXKK(0, ISGN,POLAR,PKK,PNU)
00664 CALL DWLUKK(KTO,ISGN,PKK,PNU)
00665 ELSEIF(JAK.EQ.7) THEN
00666 CALL DEXKS(0, ISGN,POLAR,PNU,PKS,PKK,PPI,JKST)
00667 CALL DWLUKS(KTO,ISGN,PNU,PKS,PKK,PPI,JKST)
00668 ELSE
00669 JNPI=JAK-7
00670 CALL DEXNEW(0, ISGN,POLAR,PNU,PWB,PNPI,JNPI)
00671 CALL DWLNEW(KTO,ISGN,PNU,PWB,PNPI,JAK)
00672 ENDIF
00673 NEVDEC(JAK)=NEVDEC(JAK)+1
00674 END
00675 SUBROUTINE DEXEL(MODE,ISGN,POL,PNU,PWB,Q1,Q2,PH)
00676
00677
00678
00679
00680
00681
00682 REAL POL(4),HV(4),PWB(4),PNU(4),Q1(4),Q2(4),PH(4),RN(1)
00683 DATA IWARM/0/
00684
00685 IF(MODE.EQ.-1) THEN
00686
00687 IWARM=1
00688 CALL DADMEL( -1,ISGN,HV,PNU,PWB,Q1,Q2,PH)
00689
00690
00691 ELSEIF(MODE.EQ. 0) THEN
00692
00693 300 CONTINUE
00694 IF(IWARM.EQ.0) GOTO 902
00695 CALL DADMEL( 0,ISGN,HV,PNU,PWB,Q1,Q2,PH)
00696 WT=(1+POL(1)*HV(1)+POL(2)*HV(2)+POL(3)*HV(3))/2.
00697
00698 CALL RANMAR(RN,1)
00699 IF(RN(1).GT.WT) GOTO 300
00700
00701 ELSEIF(MODE.EQ. 1) THEN
00702
00703 CALL DADMEL( 1,ISGN,HV,PNU,PWB,Q1,Q2,PH)
00704
00705 ENDIF
00706
00707 RETURN
00708 902 PRINT 9020
00709 9020 FORMAT(' ----- DEXEL: LACK OF INITIALISATION')
00710 STOP
00711 END
00712 SUBROUTINE DEXMU(MODE,ISGN,POL,PNU,PWB,Q1,Q2,PH)
00713
00714
00715
00716
00717
00718
00719
00720
00721 COMMON / INOUT / INUT,IOUT
00722 REAL POL(4),HV(4),PWB(4),PNU(4),Q1(4),Q2(4),PH(4),RN(1)
00723 DATA IWARM/0/
00724
00725 IF(MODE.EQ.-1) THEN
00726
00727 IWARM=1
00728 CALL DADMMU( -1,ISGN,HV,PNU,PWB,Q1,Q2,PH)
00729
00730
00731 ELSEIF(MODE.EQ. 0) THEN
00732
00733 300 CONTINUE
00734 IF(IWARM.EQ.0) GOTO 902
00735 CALL DADMMU( 0,ISGN,HV,PNU,PWB,Q1,Q2,PH)
00736 WT=(1+POL(1)*HV(1)+POL(2)*HV(2)+POL(3)*HV(3))/2.
00737
00738 CALL RANMAR(RN,1)
00739 IF(RN(1).GT.WT) GOTO 300
00740
00741 ELSEIF(MODE.EQ. 1) THEN
00742
00743 CALL DADMMU( 1,ISGN,HV,PNU,PWB,Q1,Q2,PH)
00744
00745 ENDIF
00746
00747 RETURN
00748 902 WRITE(IOUT, 9020)
00749 9020 FORMAT(' ----- DEXMU: LACK OF INITIALISATION')
00750 STOP
00751 END
00752 SUBROUTINE DADMEL(MODE,ISGN,HHV,PNU,PWB,Q1,Q2,PHX)
00753
00754
00755
00756
00757 COMMON / DECPAR / GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
00758 REAL*4 GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
00759 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
00760 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
00761 * ,AMK,AMKZ,AMKST,GAMKST
00762
00763 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
00764 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
00765 * ,AMK,AMKZ,AMKST,GAMKST
00766 COMMON / TAUBMC / GAMPMC(30),GAMPER(30),NEVDEC(30)
00767 REAL*4 GAMPMC ,GAMPER
00768 REAL*4 PHX(4)
00769 COMMON / INOUT / INUT,IOUT
00770 REAL HHV(4),HV(4),PWB(4),PNU(4),Q1(4),Q2(4)
00771 REAL PDUM1(4),PDUM2(4),PDUM3(4),PDUM4(4),PDUM5(4)
00772 REAL*4 RRR(3)
00773 REAL*8 SWT, SSWT
00774 DATA PI /3.141592653589793238462643/
00775 DATA IWARM/0/
00776
00777 IF(MODE.EQ.-1) THEN
00778
00779 IWARM=1
00780 NEVRAW=0
00781 NEVACC=0
00782 NEVOVR=0
00783 SWT=0
00784 SSWT=0
00785 WTMAX=1E-20
00786 DO 15 I=1,500
00787 CALL DPHSEL(WT,HV,PDUM1,PDUM2,PDUM3,PDUM4,PDUM5)
00788 IF(WT.GT.WTMAX/1.2) WTMAX=WT*1.2
00789 15 CONTINUE
00790
00791
00792 ELSEIF(MODE.EQ. 0) THEN
00793
00794 300 CONTINUE
00795 IF(IWARM.EQ.0) GOTO 902
00796 NEVRAW=NEVRAW+1
00797 CALL DPHSEL(WT,HV,PNU,PWB,Q1,Q2,PHX)
00798
00799 SWT=SWT+WT
00800 SSWT=SSWT+WT**2
00801 CALL RANMAR(RRR,3)
00802 RN=RRR(1)
00803 IF(WT.GT.WTMAX) NEVOVR=NEVOVR+1
00804 IF(RN*WTMAX.GT.WT) GOTO 300
00805
00806 RR2=RRR(2)
00807 COSTHE=-1.+2.*RR2
00808 THET=ACOS(COSTHE)
00809 RR3=RRR(3)
00810 PHI =2*PI*RR3
00811 CALL ROTOR2(THET,PNU,PNU)
00812 CALL ROTOR3( PHI,PNU,PNU)
00813 CALL ROTOR2(THET,PWB,PWB)
00814 CALL ROTOR3( PHI,PWB,PWB)
00815 CALL ROTOR2(THET,Q1,Q1)
00816 CALL ROTOR3( PHI,Q1,Q1)
00817 CALL ROTOR2(THET,Q2,Q2)
00818 CALL ROTOR3( PHI,Q2,Q2)
00819 CALL ROTOR2(THET,HV,HV)
00820 CALL ROTOR3( PHI,HV,HV)
00821 CALL ROTOR2(THET,PHX,PHX)
00822 CALL ROTOR3( PHI,PHX,PHX)
00823 DO 44,I=1,3
00824 44 HHV(I)=-ISGN*HV(I)
00825 NEVACC=NEVACC+1
00826
00827 ELSEIF(MODE.EQ. 1) THEN
00828
00829 IF(NEVRAW.EQ.0) RETURN
00830 PARGAM=SWT/FLOAT(NEVRAW+1)
00831 ERROR=0
00832 IF(NEVRAW.NE.0) ERROR=SQRT(SSWT/SWT**2-1./FLOAT(NEVRAW))
00833 RAT=PARGAM/GAMEL
00834 WRITE(IOUT, 7010) NEVRAW,NEVACC,NEVOVR,PARGAM,RAT,ERROR
00835
00836 GAMPMC(1)=RAT
00837 GAMPER(1)=ERROR
00838
00839 ENDIF
00840
00841 RETURN
00842 7010 FORMAT(///1X,15(5H*****)
00843 $ /,' *', 25X,'******** DADMEL FINAL REPORT ******** ',9X,1H*
00844 $ /,' *',I20 ,5X,'NEVRAW = NO. OF EL DECAYS TOTAL ',9X,1H*
00845 $ /,' *',I20 ,5X,'NEVACC = NO. OF EL DECS. ACCEPTED ',9X,1H*
00846 $ /,' *',I20 ,5X,'NEVOVR = NO. OF OVERWEIGHTED EVENTS ',9X,1H*
00847 $ /,' *',E20.5,5X,'PARTIAL WTDTH ( ELECTRON) IN GEV UNITS ',9X,1H*
00848 $ /,' *',F20.9,5X,'IN UNITS GFERMI**2*MASS**5/192/PI**3 ',9X,1H*
00849 $ /,' *',F20.9,5X,'RELATIVE ERROR OF PARTIAL WIDTH ',9X,1H*
00850 $ /,' *',25X, 'COMPLETE QED CORRECTIONS INCLUDED ',9X,1H*
00851 $ /,' *',25X, 'BUT ONLY V-A CUPLINGS ',9X,1H*
00852 $ /,1X,15(5H*****)/)
00853 902 WRITE(IOUT, 9020)
00854 9020 FORMAT(' ----- DADMEL: LACK OF INITIALISATION')
00855 STOP
00856 END
00857 SUBROUTINE DADMMU(MODE,ISGN,HHV,PNU,PWB,Q1,Q2,PHX)
00858
00859 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
00860 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
00861 * ,AMK,AMKZ,AMKST,GAMKST
00862
00863 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
00864 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
00865 * ,AMK,AMKZ,AMKST,GAMKST
00866 COMMON / DECPAR / GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
00867 REAL*4 GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
00868 COMMON / TAUBMC / GAMPMC(30),GAMPER(30),NEVDEC(30)
00869 REAL*4 GAMPMC ,GAMPER
00870 COMMON / INOUT / INUT,IOUT
00871 REAL*4 PHX(4)
00872 REAL HHV(4),HV(4),PNU(4),PWB(4),Q1(4),Q2(4)
00873 REAL PDUM1(4),PDUM2(4),PDUM3(4),PDUM4(4),PDUM5(4)
00874 REAL*4 RRR(3)
00875 REAL*8 SWT, SSWT
00876 DATA PI /3.141592653589793238462643/
00877 DATA IWARM /0/
00878
00879 IF(MODE.EQ.-1) THEN
00880
00881 IWARM=1
00882 NEVRAW=0
00883 NEVACC=0
00884 NEVOVR=0
00885 SWT=0
00886 SSWT=0
00887 WTMAX=1E-20
00888 DO 15 I=1,500
00889 CALL DPHSMU(WT,HV,PDUM1,PDUM2,PDUM3,PDUM4,PDUM5)
00890 IF(WT.GT.WTMAX/1.2) WTMAX=WT*1.2
00891 15 CONTINUE
00892
00893
00894 ELSEIF(MODE.EQ. 0) THEN
00895
00896 300 CONTINUE
00897 IF(IWARM.EQ.0) GOTO 902
00898 NEVRAW=NEVRAW+1
00899 CALL DPHSMU(WT,HV,PNU,PWB,Q1,Q2,PHX)
00900
00901 SWT=SWT+WT
00902 SSWT=SSWT+WT**2
00903 CALL RANMAR(RRR,3)
00904 RN=RRR(1)
00905 IF(WT.GT.WTMAX) NEVOVR=NEVOVR+1
00906 IF(RN*WTMAX.GT.WT) GOTO 300
00907
00908 COSTHE=-1.+2.*RRR(2)
00909 THET=ACOS(COSTHE)
00910 PHI =2*PI*RRR(3)
00911 CALL ROTOR2(THET,PNU,PNU)
00912 CALL ROTOR3( PHI,PNU,PNU)
00913 CALL ROTOR2(THET,PWB,PWB)
00914 CALL ROTOR3( PHI,PWB,PWB)
00915 CALL ROTOR2(THET,Q1,Q1)
00916 CALL ROTOR3( PHI,Q1,Q1)
00917 CALL ROTOR2(THET,Q2,Q2)
00918 CALL ROTOR3( PHI,Q2,Q2)
00919 CALL ROTOR2(THET,HV,HV)
00920 CALL ROTOR3( PHI,HV,HV)
00921 CALL ROTOR2(THET,PHX,PHX)
00922 CALL ROTOR3( PHI,PHX,PHX)
00923 DO 44,I=1,3
00924 44 HHV(I)=-ISGN*HV(I)
00925 NEVACC=NEVACC+1
00926
00927 ELSEIF(MODE.EQ. 1) THEN
00928
00929 IF(NEVRAW.EQ.0) RETURN
00930 PARGAM=SWT/FLOAT(NEVRAW+1)
00931 ERROR=0
00932 IF(NEVRAW.NE.0) ERROR=SQRT(SSWT/SWT**2-1./FLOAT(NEVRAW))
00933 RAT=PARGAM/GAMEL
00934 WRITE(IOUT, 7010) NEVRAW,NEVACC,NEVOVR,PARGAM,RAT,ERROR
00935
00936 GAMPMC(2)=RAT
00937 GAMPER(2)=ERROR
00938
00939 ENDIF
00940
00941 RETURN
00942 7010 FORMAT(///1X,15(5H*****)
00943 $ /,' *', 25X,'******** DADMMU FINAL REPORT ******** ',9X,1H*
00944 $ /,' *',I20 ,5X,'NEVRAW = NO. OF MU DECAYS TOTAL ',9X,1H*
00945 $ /,' *',I20 ,5X,'NEVACC = NO. OF MU DECS. ACCEPTED ',9X,1H*
00946 $ /,' *',I20 ,5X,'NEVOVR = NO. OF OVERWEIGHTED EVENTS ',9X,1H*
00947 $ /,' *',E20.5,5X,'PARTIAL WTDTH (MU DECAY) IN GEV UNITS ',9X,1H*
00948 $ /,' *',F20.9,5X,'IN UNITS GFERMI**2*MASS**5/192/PI**3 ',9X,1H*
00949 $ /,' *',F20.9,5X,'RELATIVE ERROR OF PARTIAL WIDTH ',9X,1H*
00950 $ /,' *',25X, 'COMPLETE QED CORRECTIONS INCLUDED ',9X,1H*
00951 $ /,' *',25X, 'BUT ONLY V-A CUPLINGS ',9X,1H*
00952 $ /,1X,15(5H*****)/)
00953 902 WRITE(IOUT, 9020)
00954 9020 FORMAT(' ----- DADMMU: LACK OF INITIALISATION')
00955 STOP
00956 END
00957 SUBROUTINE DPHSEL(DGAMX,HVX,XNX,PAAX,QPX,XAX,PHX)
00958
00959
00960
00961
00962 REAL*4 PHX(4)
00963 REAL*4 HVX(4),PAAX(4),XAX(4),QPX(4),XNX(4)
00964 REAL*8 HV(4),PH(4),PAA(4),XA(4),QP(4),XN(4)
00965 REAL*8 DGAMT
00966 IELMU=1
00967 CALL DRCMU(DGAMT,HV,PH,PAA,XA,QP,XN,IELMU)
00968 DO 7 K=1,4
00969 HVX(K)=HV(K)
00970 PHX(K)=PH(K)
00971 PAAX(K)=PAA(K)
00972 XAX(K)=XA(K)
00973 QPX(K)=QP(K)
00974 XNX(K)=XN(K)
00975 7 CONTINUE
00976 DGAMX=DGAMT
00977 END
00978 SUBROUTINE DPHSMU(DGAMX,HVX,XNX,PAAX,QPX,XAX,PHX)
00979
00980
00981
00982
00983 REAL*4 PHX(4)
00984 REAL*4 HVX(4),PAAX(4),XAX(4),QPX(4),XNX(4)
00985 REAL*8 HV(4),PH(4),PAA(4),XA(4),QP(4),XN(4)
00986 REAL*8 DGAMT
00987 IELMU=2
00988 CALL DRCMU(DGAMT,HV,PH,PAA,XA,QP,XN,IELMU)
00989 DO 7 K=1,4
00990 HVX(K)=HV(K)
00991 PHX(K)=PH(K)
00992 PAAX(K)=PAA(K)
00993 XAX(K)=XA(K)
00994 QPX(K)=QP(K)
00995 XNX(K)=XN(K)
00996 7 CONTINUE
00997 DGAMX=DGAMT
00998 END
00999 SUBROUTINE DRCMU(DGAMT,HV,PH,PAA,XA,QP,XN,IELMU)
01000 IMPLICIT REAL*8 (A-H,O-Z)
01001
01002
01003
01004
01005 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
01006 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
01007 * ,AMK,AMKZ,AMKST,GAMKST
01008
01009 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
01010 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
01011 * ,AMK,AMKZ,AMKST,GAMKST
01012 COMMON / DECPAR / GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
01013 REAL*4 GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
01014 COMMON / INOUT / INUT,IOUT
01015 COMMON / TAURAD / XK0DEC,ITDKRC
01016 REAL*8 XK0DEC
01017 REAL*8 HV(4),PT(4),PH(4),PAA(4),XA(4),QP(4),XN(4)
01018 REAL*8 PR(4)
01019 REAL*4 RRR(6)
01020 LOGICAL IHARD
01021 DATA PI /3.141592653589793238462643D0/
01022
01023
01024
01025
01026
01027 PHSPAC=1./2**17/PI**8
01028 AMTAX=AMTAU
01029
01030 PT(1)=0.D0
01031 PT(2)=0.D0
01032 PT(3)=0.D0
01033 PT(4)=AMTAX
01034
01035 CALL RANMAR(RRR,6)
01036
01037 IF (IELMU.EQ.1) THEN
01038 AMU=AMEL
01039 ELSE
01040 AMU=AMMU
01041 ENDIF
01042
01043 PRHARD=0.30D0
01044 IF ( ITDKRC.EQ.0) PRHARD=0D0
01045 PRSOFT=1.-PRHARD
01046 IF(PRSOFT.LT.0.1) THEN
01047 PRINT *, 'ERROR IN DRCMU; PRSOFT=',PRSOFT
01048 STOP
01049 ENDIF
01050
01051 RR5=RRR(5)
01052 IHARD=(RR5.GT.PRSOFT)
01053 IF (IHARD) THEN
01054
01055 RR1=RRR(1)
01056 AMS1=(AMU+AMNUTA)**2
01057 AMS2=(AMTAX)**2
01058 XK1=1-AMS1/AMS2
01059 XL1=LOG(XK1/2/XK0DEC)
01060 XL0=LOG(2*XK0DEC)
01061 XK=EXP(XL1*RR1+XL0)
01062 AM3SQ=(1-XK)*AMS2
01063 AM3 =SQRT(AM3SQ)
01064 PHSPAC=PHSPAC*AMS2*XL1*XK
01065 PHSPAC=PHSPAC/PRHARD
01066 ELSE
01067 AM3=AMTAX
01068 PHSPAC=PHSPAC*2**6*PI**3
01069 PHSPAC=PHSPAC/PRSOFT
01070 ENDIF
01071
01072 RR2=RRR(2)
01073 AMS1=(AMNUTA)**2
01074 AMS2=(AM3-AMU)**2
01075
01076
01077
01078 AM2SQ=AMS1+ RR2*(AMS2-AMS1)
01079 AM2 =SQRT(AM2SQ)
01080 PHSPAC=PHSPAC*(AMS2-AMS1)
01081
01082 ENQ1=(AM2SQ+AMNUTA**2)/(2*AM2)
01083 ENQ2=(AM2SQ-AMNUTA**2)/(2*AM2)
01084 PPI= ENQ1**2-AMNUTA**2
01085 PPPI=SQRT(ABS(ENQ1**2-AMNUTA**2))
01086 PHSPAC=PHSPAC*(4*PI)*(2*PPPI/AM2)
01087
01088 CALL SPHERD(PPPI,XN)
01089 XN(4)=ENQ1
01090
01091 DO 30 I=1,3
01092 30 XA(I)=-XN(I)
01093 XA(4)=ENQ2
01094
01095
01096 PR(1)=0
01097 PR(2)=0
01098 PR(4)=1.D0/(2*AM3)*(AM3**2+AM2**2-AMU**2)
01099 PR(3)= SQRT(ABS(PR(4)**2-AM2**2))
01100 PPI = PR(4)**2-AM2**2
01101
01102 QP(1)=0
01103 QP(2)=0
01104 QP(4)=1.D0/(2*AM3)*(AM3**2-AM2**2+AMU**2)
01105 QP(3)=-PR(3)
01106 PHSPAC=PHSPAC*(4*PI)*(2*PR(3)/AM3)
01107
01108 EXE=(PR(4)+PR(3))/AM2
01109 CALL BOSTD3(EXE,XN,XN)
01110 CALL BOSTD3(EXE,XA,XA)
01111 RR3=RRR(3)
01112 RR4=RRR(4)
01113 IF (IHARD) THEN
01114 EPS=4*(AMU/AMTAX)**2
01115 XL1=LOG((2+EPS)/EPS)
01116 XL0=LOG(EPS)
01117 ETA =EXP(XL1*RR3+XL0)
01118 CTHET=1+EPS-ETA
01119 THET =ACOS(CTHET)
01120 PHSPAC=PHSPAC*XL1/2*ETA
01121 PHI = 2*PI*RR4
01122 CALL ROTPOX(THET,PHI,XN)
01123 CALL ROTPOX(THET,PHI,XA)
01124 CALL ROTPOX(THET,PHI,QP)
01125 CALL ROTPOX(THET,PHI,PR)
01126
01127
01128
01129 PAA(1)=0
01130 PAA(2)=0
01131 PAA(4)=1/(2*AMTAX)*(AMTAX**2+AM3**2)
01132 PAA(3)= SQRT(ABS(PAA(4)**2-AM3**2))
01133 PPI = PAA(4)**2-AM3**2
01134 PHSPAC=PHSPAC*(4*PI)*(2*PAA(3)/AMTAX)
01135
01136 PH(1)=0
01137 PH(2)=0
01138 PH(4)=PAA(3)
01139 PH(3)=-PAA(3)
01140
01141
01142 EXE=(PAA(4)+PAA(3))/AM3
01143 CALL BOSTD3(EXE,XN,XN)
01144 CALL BOSTD3(EXE,XA,XA)
01145 CALL BOSTD3(EXE,QP,QP)
01146 CALL BOSTD3(EXE,PR,PR)
01147 ELSE
01148 THET =ACOS(-1.+2*RR3)
01149 PHI = 2*PI*RR4
01150 CALL ROTPOX(THET,PHI,XN)
01151 CALL ROTPOX(THET,PHI,XA)
01152 CALL ROTPOX(THET,PHI,QP)
01153 CALL ROTPOX(THET,PHI,PR)
01154
01155
01156
01157 PAA(1)=0
01158 PAA(2)=0
01159 PAA(4)=AMTAX
01160 PAA(3)=0
01161
01162 PH(1)=0
01163 PH(2)=0
01164 PH(4)=0
01165 PH(3)=0
01166 ENDIF
01167
01168 CALL DAMPRY(ITDKRC,XK0DEC,PH,XA,QP,XN,AMPLIT,HV)
01169 DGAMT=1/(2.*AMTAX)*AMPLIT*PHSPAC
01170 END
01171 SUBROUTINE DAMPRY(ITDKRC,XK0DEC,XK,XA,QP,XN,AMPLIT,HV)
01172 IMPLICIT REAL*8 (A-H,O-Z)
01173
01174
01175
01176
01177
01178 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
01179 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
01180 * ,AMK,AMKZ,AMKST,GAMKST
01181
01182 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
01183 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
01184 * ,AMK,AMKZ,AMKST,GAMKST
01185 REAL*8 HV(4),QP(4),XN(4),XA(4),XK(4)
01186
01187 HV(4)=1.D0
01188 AK0=XK0DEC*AMTAU
01189 IF(XK(4).LT.0.1D0*AK0) THEN
01190 AMPLIT=THB(ITDKRC,QP,XN,XA,AK0,HV)
01191 ELSE
01192 AMPLIT=SQM2(ITDKRC,QP,XN,XA,XK,AK0,HV)
01193 ENDIF
01194 RETURN
01195 END
01196 FUNCTION SQM2(ITDKRC,QP,XN,XA,XK,AK0,HV)
01197
01198
01199
01200
01201
01202
01203
01204
01205
01206
01207
01208
01209 IMPLICIT REAL*8(A-H,O-Z)
01210 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
01211 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
01212 * ,AMK,AMKZ,AMKST,GAMKST
01213
01214 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
01215 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
01216 * ,AMK,AMKZ,AMKST,GAMKST
01217 COMMON / DECPAR / GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
01218 REAL*4 GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
01219 COMMON / QEDPRM /ALFINV,ALFPI,XK0
01220 REAL*8 ALFINV,ALFPI,XK0
01221 REAL*8 QP(4),XN(4),XA(4),XK(4)
01222 REAL*8 R(4)
01223 REAL*8 HV(4)
01224 REAL*8 S0(3),RXA(3),RXK(3),RQP(3)
01225 DATA PI /3.141592653589793238462643D0/
01226
01227 TMASS=AMTAU
01228 GF=GFERMI
01229 ALPHAI=ALFINV
01230 TMASS2=TMASS**2
01231 EMASS2=QP(4)**2-QP(1)**2-QP(2)**2-QP(3)**2
01232 R(4)=TMASS
01233
01234 DO 7 I=1,3
01235 R(1)=0.D0
01236 R(2)=0.D0
01237 R(3)=0.D0
01238 R(I)=TMASS
01239 RXA(I)=R(4)*XA(4)-R(1)*XA(1)-R(2)*XA(2)-R(3)*XA(3)
01240
01241 RXK(I)=R(4)*XK(4)-R(1)*XK(1)-R(2)*XK(2)-R(3)*XK(3)
01242 RQP(I)=R(4)*QP(4)-R(1)*QP(1)-R(2)*QP(2)-R(3)*QP(3)
01243 7 CONTINUE
01244 QPXN=QP(4)*XN(4)-QP(1)*XN(1)-QP(2)*XN(2)-QP(3)*XN(3)
01245 QPXA=QP(4)*XA(4)-QP(1)*XA(1)-QP(2)*XA(2)-QP(3)*XA(3)
01246 QPXK=QP(4)*XK(4)-QP(1)*XK(1)-QP(2)*XK(2)-QP(3)*XK(3)
01247
01248 XNXK=XN(4)*XK(4)-XN(1)*XK(1)-XN(2)*XK(2)-XN(3)*XK(3)
01249 XAXK=XA(4)*XK(4)-XA(1)*XK(1)-XA(2)*XK(2)-XA(3)*XK(3)
01250 TXN=TMASS*XN(4)
01251 TXA=TMASS*XA(4)
01252 TQP=TMASS*QP(4)
01253 TXK=TMASS*XK(4)
01254
01255 X= XNXK/QPXN
01256 Z= TXK/TQP
01257 A= 1+X
01258 B= 1+ X*(1+Z)/2+Z/2
01259 S1= QPXN*TXA*( -EMASS2/QPXK**2*A + 2*TQP/(QPXK*TXK)*B-
01260 $TMASS2/TXK**2) +
01261 $QPXN/TXK**2* ( TMASS2*XAXK - TXA*TXK+ XAXK*TXK) -
01262 $TXA*TXN/TXK - QPXN/(QPXK*TXK)* (TQP*XAXK-TXK*QPXA)
01263 CONST4=256*PI/ALPHAI*GF**2
01264 IF (ITDKRC.EQ.0) CONST4=0D0
01265 SQM2=S1*CONST4
01266 DO 5 I=1,3
01267 S0(I) = QPXN*RXA(I)*(-EMASS2/QPXK**2*A + 2*TQP/(QPXK*TXK)*B-
01268 $ TMASS2/TXK**2) +
01269 $ QPXN/TXK**2* (TMASS2*XAXK - TXA*RXK(I)+ XAXK*RXK(I))-
01270 $ RXA(I)*TXN/TXK - QPXN/(QPXK*TXK)*(RQP(I)*XAXK- RXK(I)*QPXA)
01271 5 HV(I)=S0(I)/S1-1.D0
01272 RETURN
01273 END
01274 FUNCTION THB(ITDKRC,QP,XN,XA,AK0,HV)
01275
01276
01277
01278
01279
01280
01281
01282
01283
01284
01285
01286
01287
01288 IMPLICIT REAL*8(A-H,O-Z)
01289 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
01290 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
01291 * ,AMK,AMKZ,AMKST,GAMKST
01292
01293 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
01294 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
01295 * ,AMK,AMKZ,AMKST,GAMKST
01296 COMMON / DECPAR / GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
01297 REAL*4 GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
01298 COMMON / QEDPRM /ALFINV,ALFPI,XK0
01299 REAL*8 ALFINV,ALFPI,XK0
01300 DIMENSION QP(4),XN(4),XA(4)
01301 REAL*8 HV(4)
01302 DIMENSION R(4)
01303 REAL*8 RXA(3),RXN(3),RQP(3)
01304 REAL*8 BORNPL(3),AM3POL(3),XM3POL(3)
01305 DATA PI /3.141592653589793238462643D0/
01306
01307 TMASS=AMTAU
01308 GF=GFERMI
01309 ALPHAI=ALFINV
01310
01311 TMASS2=TMASS**2
01312 R(4)=TMASS
01313 DO 7 I=1,3
01314 R(1)=0.D0
01315 R(2)=0.D0
01316 R(3)=0.D0
01317 R(I)=TMASS
01318 RXA(I)=R(4)*XA(4)-R(1)*XA(1)-R(2)*XA(2)-R(3)*XA(3)
01319 RXN(I)=R(4)*XN(4)-R(1)*XN(1)-R(2)*XN(2)-R(3)*XN(3)
01320
01321 RQP(I)=R(4)*QP(4)-R(1)*QP(1)-R(2)*QP(2)-R(3)*QP(3)
01322 7 CONTINUE
01323
01324 U0=QP(4)/TMASS
01325 U3=SQRT(QP(1)**2+QP(2)**2+QP(3)**2)/TMASS
01326 W3=U3
01327 W0=(XN(4)+XA(4))/TMASS
01328 UP=U0+U3
01329 UM=U0-U3
01330 WP=W0+W3
01331 WM=W0-W3
01332 YU=LOG(UP/UM)/2
01333 YW=LOG(WP/WM)/2
01334 EPS2=U0**2-U3**2
01335 EPS=SQRT(EPS2)
01336 Y=W0**2-W3**2
01337 AL=AK0/TMASS
01338
01339 F0=2*U0/U3*( DILOGT(1-(UM*WM/(UP*WP)))- DILOGT(1-WM/WP) +
01340 $DILOGT(1-UM/UP) -2*YU+ 2*LOG(UP)*(YW+YU) ) +
01341 $1/Y* ( 2*U3*YU + (1-EPS2- 2*Y)*LOG(EPS) ) +
01342 $ 2 - 4*(U0/U3*YU -1)* LOG(2*AL)
01343 FP= YU/(2*U3)*(1 + (1-EPS2)/Y ) + LOG(EPS)/Y
01344 FM= YU/(2*U3)*(1 - (1-EPS2)/Y ) - LOG(EPS)/Y
01345 F3= EPS2*(FP+FM)/2
01346
01347 QPXN=QP(4)*XN(4)-QP(1)*XN(1)-QP(2)*XN(2)-QP(3)*XN(3)
01348 QPXA=QP(4)*XA(4)-QP(1)*XA(1)-QP(2)*XA(2)-QP(3)*XA(3)
01349 XNXA=XN(4)*XA(4)-XN(1)*XA(1)-XN(2)*XA(2)-XN(3)*XA(3)
01350 TXN=TMASS*XN(4)
01351 TXA=TMASS*XA(4)
01352 TQP=TMASS*QP(4)
01353
01354 CONST3=1/(2*ALPHAI*PI)*64*GF**2
01355 IF (ITDKRC.EQ.0) CONST3=0D0
01356 XM3= -( F0* QPXN*TXA + FP*EPS2* TXN*TXA +
01357 $FM* QPXN*QPXA + F3* TMASS2*XNXA )
01358 AM3=XM3*CONST3
01359
01360 BRAK= (GV+GA)**2*TQP*XNXA+(GV-GA)**2*TXA*QPXN
01361 & -(GV**2-GA**2)*TMASS*AMNUTA*QPXA
01362 BORN= 32*(GFERMI**2/2.)*BRAK
01363 DO 5 I=1,3
01364 XM3POL(I)= -( F0* QPXN*RXA(I) + FP*EPS2* TXN*RXA(I) +
01365 $ FM* QPXN* (QPXA + (RXA(I)*TQP-TXA*RQP(I))/TMASS2 ) +
01366 $ F3* (TMASS2*XNXA +TXN*RXA(I) -RXN(I)*TXA) )
01367 AM3POL(I)=XM3POL(I)*CONST3
01368
01369 BORNPL(I)=BORN+(
01370 & (GV+GA)**2*TMASS*XNXA*QP(I)
01371 & -(GV-GA)**2*TMASS*QPXN*XA(I)
01372 & +(GV**2-GA**2)*AMNUTA*TXA*QP(I)
01373 & -(GV**2-GA**2)*AMNUTA*TQP*XA(I) )*
01374 & 32*(GFERMI**2/2.)
01375 5 HV(I)=(BORNPL(I)+AM3POL(I))/(BORN+AM3)-1.D0
01376 THB=BORN+AM3
01377 IF (THB/BORN.LT.0.1D0) THEN
01378 PRINT *, 'ERROR IN THB, THB/BORN=',THB/BORN
01379 THB=0.D0
01380 ENDIF
01381 RETURN
01382 END
01383 SUBROUTINE DEXPI(MODE,ISGN,POL,PPI,PNU)
01384
01385
01386
01387
01388
01389
01390 REAL POL(4),HV(4),PNU(4),PPI(4),RN(1)
01391
01392 IF(MODE.EQ.-1) THEN
01393
01394 CALL DADMPI(-1,ISGN,HV,PPI,PNU)
01395
01396
01397 ELSEIF(MODE.EQ. 0) THEN
01398
01399 300 CONTINUE
01400 CALL DADMPI( 0,ISGN,HV,PPI,PNU)
01401 WT=(1+POL(1)*HV(1)+POL(2)*HV(2)+POL(3)*HV(3))/2.
01402
01403 CALL RANMAR(RN,1)
01404 IF(RN(1).GT.WT) GOTO 300
01405
01406 ELSEIF(MODE.EQ. 1) THEN
01407
01408 CALL DADMPI( 1,ISGN,HV,PPI,PNU)
01409
01410 ENDIF
01411
01412 RETURN
01413 END
01414 SUBROUTINE DADMPI(MODE,ISGN,HV,PPI,PNU)
01415
01416 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
01417 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
01418 * ,AMK,AMKZ,AMKST,GAMKST
01419
01420 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
01421 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
01422 * ,AMK,AMKZ,AMKST,GAMKST
01423 COMMON / DECPAR / GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
01424 REAL*4 GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
01425 COMMON / TAUBMC / GAMPMC(30),GAMPER(30),NEVDEC(30)
01426 REAL*4 GAMPMC ,GAMPER
01427 COMMON / INOUT / INUT,IOUT
01428 REAL PPI(4),PNU(4),HV(4)
01429 DATA PI /3.141592653589793238462643/
01430
01431 IF(MODE.EQ.-1) THEN
01432
01433 NEVTOT=0
01434 ELSEIF(MODE.EQ. 0) THEN
01435
01436 NEVTOT=NEVTOT+1
01437 EPI= (AMTAU**2+AMPI**2-AMNUTA**2)/(2*AMTAU)
01438 ENU= (AMTAU**2-AMPI**2+AMNUTA**2)/(2*AMTAU)
01439 XPI= SQRT(EPI**2-AMPI**2)
01440
01441 CALL SPHERA(XPI,PPI)
01442 PPI(4)=EPI
01443
01444 DO 30 I=1,3
01445 30 PNU(I)=-PPI(I)
01446 PNU(4)=ENU
01447 PXQ=AMTAU*EPI
01448 PXN=AMTAU*ENU
01449 QXN=PPI(4)*PNU(4)-PPI(1)*PNU(1)-PPI(2)*PNU(2)-PPI(3)*PNU(3)
01450 BRAK=(GV**2+GA**2)*(2*PXQ*QXN-AMPI**2*PXN)
01451 & +(GV**2-GA**2)*AMTAU*AMNUTA*AMPI**2
01452 DO 40 I=1,3
01453 40 HV(I)=-ISGN*2*GA*GV*AMTAU*(2*PPI(I)*QXN-PNU(I)*AMPI**2)/BRAK
01454 HV(4)=1
01455
01456 ELSEIF(MODE.EQ. 1) THEN
01457
01458 IF(NEVTOT.EQ.0) RETURN
01459 FPI=0.1284
01460
01461
01462
01463
01464 GAMM=(GFERMI*FPI)**2/(16.*PI)*AMTAU**3*
01465 $ (BRAK/AMTAU**4)*
01466 $ SQRT((AMTAU**2-AMPI**2-AMNUTA**2)**2
01467 $ -4*AMPI**2*AMNUTA**2 )/AMTAU**2
01468 ERROR=0
01469 RAT=GAMM/GAMEL
01470 WRITE(IOUT, 7010) NEVTOT,GAMM,RAT,ERROR
01471 GAMPMC(3)=RAT
01472 GAMPER(3)=ERROR
01473
01474 ENDIF
01475
01476 RETURN
01477 7010 FORMAT(///1X,15(5H*****)
01478 $ /,' *', 25X,'******** DADMPI FINAL REPORT ******** ',9X,1H*
01479 $ /,' *',I20 ,5X,'NEVTOT = NO. OF PI DECAYS TOTAL ',9X,1H*
01480 $ /,' *',E20.5,5X,'PARTIAL WTDTH ( PI DECAY) IN GEV UNITS ',9X,1H*
01481 $ /,' *',F20.9,5X,'IN UNITS GFERMI**2*MASS**5/192/PI**3 ',9X,1H*
01482 $ /,' *',F20.8,5X,'RELATIVE ERROR OF PARTIAL WIDTH (STAT.)',9X,1H*
01483 $ /,1X,15(5H*****)/)
01484 END
01485 SUBROUTINE DEXRO(MODE,ISGN,POL,PNU,PRO,PIC,PIZ)
01486
01487
01488
01489
01490
01491
01492
01493
01494 COMMON / INOUT / INUT,IOUT
01495 REAL POL(4),HV(4),PRO(4),PNU(4),PIC(4),PIZ(4),RN(1)
01496 DATA IWARM/0/
01497
01498 IF(MODE.EQ.-1) THEN
01499
01500 IWARM=1
01501 CALL DADMRO( -1,ISGN,HV,PNU,PRO,PIC,PIZ)
01502
01503
01504
01505 ELSEIF(MODE.EQ. 0) THEN
01506
01507 300 CONTINUE
01508 IF(IWARM.EQ.0) GOTO 902
01509 CALL DADMRO( 0,ISGN,HV,PNU,PRO,PIC,PIZ)
01510 WT=(1+POL(1)*HV(1)+POL(2)*HV(2)+POL(3)*HV(3))/2.
01511
01512
01513
01514 CALL RANMAR(RN,1)
01515 IF(RN(1).GT.WT) GOTO 300
01516
01517 ELSEIF(MODE.EQ. 1) THEN
01518
01519 CALL DADMRO( 1,ISGN,HV,PNU,PRO,PIC,PIZ)
01520
01521
01522 ENDIF
01523
01524 RETURN
01525 902 WRITE(IOUT, 9020)
01526 9020 FORMAT(' ----- DEXRO: LACK OF INITIALISATION')
01527 STOP
01528 END
01529 SUBROUTINE DADMRO(MODE,ISGN,HHV,PNU,PRO,PIC,PIZ)
01530
01531 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
01532 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
01533 * ,AMK,AMKZ,AMKST,GAMKST
01534
01535 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
01536 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
01537 * ,AMK,AMKZ,AMKST,GAMKST
01538 COMMON / DECPAR / GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
01539 REAL*4 GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
01540 COMMON / TAUBMC / GAMPMC(30),GAMPER(30),NEVDEC(30)
01541 REAL*4 GAMPMC ,GAMPER
01542 COMMON / INOUT / INUT,IOUT
01543 REAL HHV(4)
01544 REAL HV(4),PRO(4),PNU(4),PIC(4),PIZ(4)
01545 REAL PDUM1(4),PDUM2(4),PDUM3(4),PDUM4(4)
01546 REAL*4 RRR(3)
01547 REAL*8 SWT, SSWT
01548 DATA PI /3.141592653589793238462643/
01549 DATA IWARM/0/
01550
01551 IF(MODE.EQ.-1) THEN
01552
01553 IWARM=1
01554 NEVRAW=0
01555 NEVACC=0
01556 NEVOVR=0
01557 SWT=0
01558 SSWT=0
01559 WTMAX=1E-20
01560 DO 15 I=1,500
01561 CALL DPHSRO(WT,HV,PDUM1,PDUM2,PDUM3,PDUM4)
01562 IF(WT.GT.WTMAX/1.2) WTMAX=WT*1.2
01563 15 CONTINUE
01564
01565
01566
01567 ELSEIF(MODE.EQ. 0) THEN
01568
01569 300 CONTINUE
01570 IF(IWARM.EQ.0) GOTO 902
01571 CALL DPHSRO(WT,HV,PNU,PRO,PIC,PIZ)
01572
01573 NEVRAW=NEVRAW+1
01574 SWT=SWT+WT
01575 SSWT=SSWT+WT**2
01576 CALL RANMAR(RRR,3)
01577 RN=RRR(1)
01578 IF(WT.GT.WTMAX) NEVOVR=NEVOVR+1
01579 IF(RN*WTMAX.GT.WT) GOTO 300
01580
01581 COSTHE=-1.+2.*RRR(2)
01582 THET=ACOS(COSTHE)
01583 PHI =2*PI*RRR(3)
01584 CALL ROTOR2(THET,PNU,PNU)
01585 CALL ROTOR3( PHI,PNU,PNU)
01586 CALL ROTOR2(THET,PRO,PRO)
01587 CALL ROTOR3( PHI,PRO,PRO)
01588 CALL ROTOR2(THET,PIC,PIC)
01589 CALL ROTOR3( PHI,PIC,PIC)
01590 CALL ROTOR2(THET,PIZ,PIZ)
01591 CALL ROTOR3( PHI,PIZ,PIZ)
01592 CALL ROTOR2(THET,HV,HV)
01593 CALL ROTOR3( PHI,HV,HV)
01594 DO 44 I=1,3
01595 44 HHV(I)=-ISGN*HV(I)
01596 NEVACC=NEVACC+1
01597
01598 ELSEIF(MODE.EQ. 1) THEN
01599
01600 IF(NEVRAW.EQ.0) RETURN
01601 PARGAM=SWT/FLOAT(NEVRAW+1)
01602 ERROR=0
01603 IF(NEVRAW.NE.0) ERROR=SQRT(SSWT/SWT**2-1./FLOAT(NEVRAW))
01604 RAT=PARGAM/GAMEL
01605 WRITE(IOUT, 7010) NEVRAW,NEVACC,NEVOVR,PARGAM,RAT,ERROR
01606
01607 GAMPMC(4)=RAT
01608 GAMPER(4)=ERROR
01609
01610 ENDIF
01611
01612 RETURN
01613 7003 FORMAT(///1X,15(5H*****)
01614 $ /,' *', 25X,'******** DADMRO INITIALISATION ********',9X,1H*
01615 $ /,' *',E20.5,5X,'WTMAX = MAXIMUM WEIGHT ',9X,1H*
01616 $ /,1X,15(5H*****)/)
01617 7010 FORMAT(///1X,15(5H*****)
01618 $ /,' *', 25X,'******** DADMRO FINAL REPORT ******** ',9X,1H*
01619 $ /,' *',I20 ,5X,'NEVRAW = NO. OF RHO DECAYS TOTAL ',9X,1H*
01620 $ /,' *',I20 ,5X,'NEVACC = NO. OF RHO DECS. ACCEPTED ',9X,1H*
01621 $ /,' *',I20 ,5X,'NEVOVR = NO. OF OVERWEIGHTED EVENTS ',9X,1H*
01622 $ /,' *',E20.5,5X,'PARTIAL WTDTH (RHO DECAY) IN GEV UNITS ',9X,1H*
01623 $ /,' *',F20.9,5X,'IN UNITS GFERMI**2*MASS**5/192/PI**3 ',9X,1H*
01624 $ /,' *',F20.8,5X,'RELATIVE ERROR OF PARTIAL WIDTH ',9X,1H*
01625 $ /,1X,15(5H*****)/)
01626 902 WRITE(IOUT, 9020)
01627 9020 FORMAT(' ----- DADMRO: LACK OF INITIALISATION')
01628 STOP
01629 END
01630 SUBROUTINE DPHSRO(DGAMT,HV,PN,PR,PIC,PIZ)
01631
01632
01633
01634
01635 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
01636 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
01637 * ,AMK,AMKZ,AMKST,GAMKST
01638
01639 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
01640 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
01641 * ,AMK,AMKZ,AMKST,GAMKST
01642 COMMON / DECPAR / GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
01643 REAL*4 GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
01644 REAL HV(4),PT(4),PN(4),PR(4),PIC(4),PIZ(4),QQ(4),RR1(1)
01645 DATA PI /3.141592653589793238462643/
01646 DATA ICONT /0/
01647
01648
01649 PHSPAC=1./2**11/PI**5
01650
01651 PT(1)=0.
01652 PT(2)=0.
01653 PT(3)=0.
01654 PT(4)=AMTAU
01655
01656 AMS1=(AMPI+AMPIZ)**2
01657 AMS2=(AMTAU-AMNUTA)**2
01658
01659
01660
01661
01662
01663 ALP1=ATAN((AMS1-AMRO**2)/AMRO/GAMRO)
01664 ALP2=ATAN((AMS2-AMRO**2)/AMRO/GAMRO)
01665
01666 100 CONTINUE
01667 CALL RANMAR(RR1,1)
01668 ALP=ALP1+RR1(1)*(ALP2-ALP1)
01669 AMX2=AMRO**2+AMRO*GAMRO*TAN(ALP)
01670 AMX=SQRT(AMX2)
01671 IF(AMX.LT.2.*AMPI) GO TO 100
01672
01673 PHSPAC=PHSPAC*((AMX2-AMRO**2)**2+(AMRO*GAMRO)**2)/(AMRO*GAMRO)
01674 PHSPAC=PHSPAC*(ALP2-ALP1)
01675
01676
01677 PN(1)=0
01678 PN(2)=0
01679 PN(4)=1./(2*AMTAU)*(AMTAU**2+AMNUTA**2-AMX**2)
01680 PN(3)=-SQRT((PN(4)-AMNUTA)*(PN(4)+AMNUTA))
01681
01682 PR(1)=0
01683 PR(2)=0
01684 PR(4)=1./(2*AMTAU)*(AMTAU**2-AMNUTA**2+AMX**2)
01685 PR(3)=-PN(3)
01686 PHSPAC=PHSPAC*(4*PI)*(2*PR(3)/AMTAU)
01687
01688
01689 ENQ1=(AMX2+AMPI**2-AMPIZ**2)/(2.*AMX)
01690 ENQ2=(AMX2-AMPI**2+AMPIZ**2)/(2.*AMX)
01691 PPPI=SQRT((ENQ1-AMPI)*(ENQ1+AMPI))
01692 PHSPAC=PHSPAC*(4*PI)*(2*PPPI/AMX)
01693
01694 CALL SPHERA(PPPI,PIC)
01695 PIC(4)=ENQ1
01696
01697 DO 20 I=1,3
01698 20 PIZ(I)=-PIC(I)
01699 PIZ(4)=ENQ2
01700 EXE=(PR(4)+PR(3))/AMX
01701
01702 CALL BOSTR3(EXE,PIC,PIC)
01703 CALL BOSTR3(EXE,PIZ,PIZ)
01704 DO 30 I=1,4
01705 30 QQ(I)=PIC(I)-PIZ(I)
01706
01707 PRODPQ=PT(4)*QQ(4)
01708 PRODNQ=PN(4)*QQ(4)-PN(1)*QQ(1)-PN(2)*QQ(2)-PN(3)*QQ(3)
01709 PRODPN=PT(4)*PN(4)
01710 QQ2= QQ(4)**2-QQ(1)**2-QQ(2)**2-QQ(3)**2
01711 BRAK=(GV**2+GA**2)*(2*PRODPQ*PRODNQ-PRODPN*QQ2)
01712 & +(GV**2-GA**2)*AMTAU*AMNUTA*QQ2
01713 AMPLIT=(GFERMI*CCABIB)**2*BRAK*2*FPIRHO(AMX)
01714 DGAMT=1/(2.*AMTAU)*AMPLIT*PHSPAC
01715 DO 40 I=1,3
01716 40 HV(I)=2*GV*GA*AMTAU*(2*PRODNQ*QQ(I)-QQ2*PN(I))/BRAK
01717 RETURN
01718 END
01719 SUBROUTINE DEXAA(MODE,ISGN,POL,PNU,PAA,PIM1,PIM2,PIPL,JAA)
01720
01721
01722
01723
01724
01725
01726
01727
01728
01729
01730 COMMON / INOUT / INUT,IOUT
01731 REAL POL(4),HV(4),PAA(4),PNU(4),PIM1(4),PIM2(4),PIPL(4),RN(1)
01732 DATA IWARM/0/
01733
01734 IF(MODE.EQ.-1) THEN
01735
01736 IWARM=1
01737 CALL DADMAA( -1,ISGN,HV,PNU,PAA,PIM1,PIM2,PIPL,JAA)
01738
01739
01740 ELSEIF(MODE.EQ. 0) THEN
01741
01742 300 CONTINUE
01743 IF(IWARM.EQ.0) GOTO 902
01744 CALL DADMAA( 0,ISGN,HV,PNU,PAA,PIM1,PIM2,PIPL,JAA)
01745 WT=(1+POL(1)*HV(1)+POL(2)*HV(2)+POL(3)*HV(3))/2.
01746
01747 CALL RANMAR(RN,1)
01748 IF(RN(1).GT.WT) GOTO 300
01749
01750 ELSEIF(MODE.EQ. 1) THEN
01751
01752 CALL DADMAA( 1,ISGN,HV,PNU,PAA,PIM1,PIM2,PIPL,JAA)
01753
01754 ENDIF
01755
01756 RETURN
01757 902 WRITE(IOUT, 9020)
01758 9020 FORMAT(' ----- DEXAA: LACK OF INITIALISATION')
01759 STOP
01760 END
01761 SUBROUTINE DADMAA(MODE,ISGN,HHV,PNU,PAA,PIM1,PIM2,PIPL,JAA)
01762
01763
01764
01765 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
01766 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
01767 * ,AMK,AMKZ,AMKST,GAMKST
01768
01769 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
01770 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
01771 * ,AMK,AMKZ,AMKST,GAMKST
01772 COMMON / DECPAR / GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
01773 REAL*4 GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
01774 COMMON / TAUBMC / GAMPMC(30),GAMPER(30),NEVDEC(30)
01775 REAL*4 GAMPMC ,GAMPER
01776 COMMON / INOUT / INUT,IOUT
01777 REAL HHV(4)
01778 REAL HV(4),PAA(4),PNU(4),PIM1(4),PIM2(4),PIPL(4)
01779 REAL PDUM1(4),PDUM2(4),PDUM3(4),PDUM4(4),PDUM5(4)
01780 REAL*4 RRR(3)
01781 REAL*8 SWT, SSWT
01782 DATA PI /3.141592653589793238462643/
01783 DATA IWARM/0/
01784
01785 IF(MODE.EQ.-1) THEN
01786
01787 IWARM=1
01788 NEVRAW=0
01789 NEVACC=0
01790 NEVOVR=0
01791 SWT=0
01792 SSWT=0
01793 WTMAX=1E-20
01794 DO 15 I=1,500
01795 CALL DPHSAA(WT,HV,PDUM1,PDUM2,PDUM3,PDUM4,PDUM5,JAA)
01796 IF(WT.GT.WTMAX/1.2) WTMAX=WT*1.2
01797 15 CONTINUE
01798
01799
01800 ELSEIF(MODE.EQ. 0) THEN
01801
01802 300 CONTINUE
01803 IF(IWARM.EQ.0) GOTO 902
01804 CALL DPHSAA(WT,HV,PNU,PAA,PIM1,PIM2,PIPL,JAA)
01805
01806 NEVRAW=NEVRAW+1
01807 SWT=SWT+WT
01808
01809
01810 SSWT=SSWT+dble(WT)**2
01811
01812 CALL RANMAR(RRR,3)
01813 RN=RRR(1)
01814 IF(WT.GT.WTMAX) NEVOVR=NEVOVR+1
01815 IF(RN*WTMAX.GT.WT) GOTO 300
01816
01817 COSTHE=-1.+2.*RRR(2)
01818 THET=ACOS(COSTHE)
01819 PHI =2*PI*RRR(3)
01820 CALL ROTPOL(THET,PHI,PNU)
01821 CALL ROTPOL(THET,PHI,PAA)
01822 CALL ROTPOL(THET,PHI,PIM1)
01823 CALL ROTPOL(THET,PHI,PIM2)
01824 CALL ROTPOL(THET,PHI,PIPL)
01825 CALL ROTPOL(THET,PHI,HV)
01826 DO 44 I=1,3
01827 44 HHV(I)=-ISGN*HV(I)
01828 NEVACC=NEVACC+1
01829
01830 ELSEIF(MODE.EQ. 1) THEN
01831
01832 IF(NEVRAW.EQ.0) RETURN
01833 PARGAM=SWT/FLOAT(NEVRAW+1)
01834 ERROR=0
01835 IF(NEVRAW.NE.0) ERROR=SQRT(SSWT/SWT**2-1./FLOAT(NEVRAW))
01836 RAT=PARGAM/GAMEL
01837 WRITE(IOUT, 7010) NEVRAW,NEVACC,NEVOVR,PARGAM,RAT,ERROR
01838
01839 GAMPMC(5)=RAT
01840 GAMPER(5)=ERROR
01841
01842 ENDIF
01843
01844 RETURN
01845 7003 FORMAT(///1X,15(5H*****)
01846 $ /,' *', 25X,'******** DADMAA INITIALISATION ********',9X,1H*
01847 $ /,' *',E20.5,5X,'WTMAX = MAXIMUM WEIGHT ',9X,1H*
01848 $ /,1X,15(5H*****)/)
01849 7010 FORMAT(///1X,15(5H*****)
01850 $ /,' *', 25X,'******** DADMAA FINAL REPORT ******** ',9X,1H*
01851 $ /,' *',I20 ,5X,'NEVRAW = NO. OF A1 DECAYS TOTAL ',9X,1H*
01852 $ /,' *',I20 ,5X,'NEVACC = NO. OF A1 DECS. ACCEPTED ',9X,1H*
01853 $ /,' *',I20 ,5X,'NEVOVR = NO. OF OVERWEIGHTED EVENTS ',9X,1H*
01854 $ /,' *',E20.5,5X,'PARTIAL WTDTH (A1 DECAY) IN GEV UNITS ',9X,1H*
01855 $ /,' *',F20.9,5X,'IN UNITS GFERMI**2*MASS**5/192/PI**3 ',9X,1H*
01856 $ /,' *',F20.8,5X,'RELATIVE ERROR OF PARTIAL WIDTH ',9X,1H*
01857 $ /,1X,15(5H*****)/)
01858 902 WRITE(IOUT, 9020)
01859 9020 FORMAT(' ----- DADMAA: LACK OF INITIALISATION')
01860 STOP
01861 END
01862 SUBROUTINE DPHSAA(DGAMT,HV,PN,PAA,PIM1,PIM2,PIPL,JAA)
01863
01864
01865
01866
01867 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
01868 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
01869 * ,AMK,AMKZ,AMKST,GAMKST
01870
01871 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
01872 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
01873 * ,AMK,AMKZ,AMKST,GAMKST
01874 COMMON / TAUKLE / BRA1,BRK0,BRK0B,BRKS
01875 REAL*4 BRA1,BRK0,BRK0B,BRKS
01876 REAL HV(4),PN(4),PAA(4),PIM1(4),PIM2(4),PIPL(4)
01877
01878
01879 REAL*4 RRR(1)
01880
01881 MNUM=0
01882
01883 KEYT=1
01884 CALL RANMAR(RRR,1)
01885 RMOD=RRR(1)
01886 IF (RMOD.LT.BRA1) THEN
01887 JAA=1
01888 AMP1=AMPI
01889 AMP2=AMPI
01890 AMP3=AMPI
01891 ELSE
01892 JAA=2
01893 AMP1=AMPIZ
01894 AMP2=AMPIZ
01895 AMP3=AMPI
01896 ENDIF
01897 CALL
01898 $ DPHTRE(DGAMT,HV,PN,PAA,PIM1,AMP1,PIM2,AMP2,PIPL,AMP3,KEYT,MNUM)
01899 END
01900 SUBROUTINE DEXKK(MODE,ISGN,POL,PKK,PNU)
01901
01902
01903
01904
01905
01906
01907 REAL POL(4),HV(4),PNU(4),PKK(4),RN(1)
01908
01909 IF(MODE.EQ.-1) THEN
01910
01911 CALL DADMKK(-1,ISGN,HV,PKK,PNU)
01912
01913
01914 ELSEIF(MODE.EQ. 0) THEN
01915
01916 300 CONTINUE
01917 CALL DADMKK( 0,ISGN,HV,PKK,PNU)
01918 WT=(1+POL(1)*HV(1)+POL(2)*HV(2)+POL(3)*HV(3))/2.
01919
01920 CALL RANMAR(RN,1)
01921 IF(RN(1).GT.WT) GOTO 300
01922
01923 ELSEIF(MODE.EQ. 1) THEN
01924
01925 CALL DADMKK( 1,ISGN,HV,PKK,PNU)
01926
01927 ENDIF
01928
01929 RETURN
01930 END
01931 SUBROUTINE DADMKK(MODE,ISGN,HV,PKK,PNU)
01932
01933
01934 COMMON / DECPAR / GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
01935 REAL*4 GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
01936 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
01937 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
01938 * ,AMK,AMKZ,AMKST,GAMKST
01939
01940 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
01941 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
01942 * ,AMK,AMKZ,AMKST,GAMKST
01943 COMMON / TAUBMC / GAMPMC(30),GAMPER(30),NEVDEC(30)
01944 REAL*4 GAMPMC ,GAMPER
01945 COMMON / INOUT / INUT,IOUT
01946 REAL PKK(4),PNU(4),HV(4)
01947 DATA PI /3.141592653589793238462643/
01948
01949 IF(MODE.EQ.-1) THEN
01950
01951 NEVTOT=0
01952 ELSEIF(MODE.EQ. 0) THEN
01953
01954 NEVTOT=NEVTOT+1
01955 EKK= (AMTAU**2+AMK**2-AMNUTA**2)/(2*AMTAU)
01956 ENU= (AMTAU**2-AMK**2+AMNUTA**2)/(2*AMTAU)
01957 XKK= SQRT(EKK**2-AMK**2)
01958
01959 CALL SPHERA(XKK,PKK)
01960 PKK(4)=EKK
01961
01962 DO 30 I=1,3
01963 30 PNU(I)=-PKK(I)
01964 PNU(4)=ENU
01965 PXQ=AMTAU*EKK
01966 PXN=AMTAU*ENU
01967 QXN=PKK(4)*PNU(4)-PKK(1)*PNU(1)-PKK(2)*PNU(2)-PKK(3)*PNU(3)
01968 BRAK=(GV**2+GA**2)*(2*PXQ*QXN-AMK**2*PXN)
01969 & +(GV**2-GA**2)*AMTAU*AMNUTA*AMK**2
01970 DO 40 I=1,3
01971 40 HV(I)=-ISGN*2*GA*GV*AMTAU*(2*PKK(I)*QXN-PNU(I)*AMK**2)/BRAK
01972 HV(4)=1
01973
01974 ELSEIF(MODE.EQ. 1) THEN
01975
01976 IF(NEVTOT.EQ.0) RETURN
01977 FKK=0.0354
01978
01979
01980
01981
01982
01983 GAMM=(GFERMI*FKK)**2/(16.*PI)*AMTAU**3*
01984 $ (BRAK/AMTAU**4)*
01985 $ SQRT((AMTAU**2-AMK**2-AMNUTA**2)**2
01986 $ -4*AMK**2*AMNUTA**2 )/AMTAU**2
01987 ERROR=0
01988
01989 ERROR=0
01990 RAT=GAMM/GAMEL
01991 WRITE(IOUT, 7010) NEVTOT,GAMM,RAT,ERROR
01992 GAMPMC(6)=RAT
01993 GAMPER(6)=ERROR
01994
01995 ENDIF
01996
01997 RETURN
01998 7010 FORMAT(///1X,15(5H*****)
01999 $ /,' *', 25X,'******** DADMKK FINAL REPORT ********',9X,1H*
02000 $ /,' *',I20 ,5X,'NEVTOT = NO. OF K DECAYS TOTAL ',9X,1H*,
02001 $ /,' *',E20.5,5X,'PARTIAL WTDTH ( K DECAY) IN GEV UNITS ',9X,1H*,
02002 $ /,' *',F20.9,5X,'IN UNITS GFERMI**2*MASS**5/192/PI**3 ',9X,1H*
02003 $ /,' *',F20.8,5X,'RELATIVE ERROR OF PARTIAL WIDTH (STAT.)',9X,1H*
02004 $ /,1X,15(5H*****)/)
02005 END
02006 SUBROUTINE DEXKS(MODE,ISGN,POL,PNU,PKS,PKK,PPI,JKST)
02007
02008
02009
02010
02011
02012
02013
02014
02015
02016
02017
02018 COMMON / INOUT / INUT,IOUT
02019 REAL POL(4),HV(4),PKS(4),PNU(4),PKK(4),PPI(4),RN(1)
02020 DATA IWARM/0/
02021
02022 IF(MODE.EQ.-1) THEN
02023
02024 IWARM=1
02025
02026 CALL DADMKS( -1,ISGN,HV,PNU,PKS,PKK,PPI,JKST)
02027
02028
02029
02030 ELSEIF(MODE.EQ. 0) THEN
02031
02032 300 CONTINUE
02033 IF(IWARM.EQ.0) GOTO 902
02034 CALL DADMKS( 0,ISGN,HV,PNU,PKS,PKK,PPI,JKST)
02035 WT=(1+POL(1)*HV(1)+POL(2)*HV(2)+POL(3)*HV(3))/2.
02036
02037
02038
02039 CALL RANMAR(RN,1)
02040 IF(RN(1).GT.WT) GOTO 300
02041
02042 ELSEIF(MODE.EQ. 1) THEN
02043
02044 CALL DADMKS( 1,ISGN,HV,PNU,PKS,PKK,PPI,JKST)
02045
02046
02047 ENDIF
02048
02049 RETURN
02050 902 WRITE(IOUT, 9020)
02051 9020 FORMAT(' ----- DEXKS: LACK OF INITIALISATION')
02052 STOP
02053 END
02054 SUBROUTINE DADMKS(MODE,ISGN,HHV,PNU,PKS,PKK,PPI,JKST)
02055
02056 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
02057 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
02058 * ,AMK,AMKZ,AMKST,GAMKST
02059
02060 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
02061 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
02062 * ,AMK,AMKZ,AMKST,GAMKST
02063 COMMON / DECPAR / GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
02064 REAL*4 GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
02065 COMMON / TAUBMC / GAMPMC(30),GAMPER(30),NEVDEC(30)
02066 REAL*4 GAMPMC ,GAMPER
02067 COMMON / TAUKLE / BRA1,BRK0,BRK0B,BRKS
02068 REAL*4 BRA1,BRK0,BRK0B,BRKS
02069 COMMON / INOUT / INUT,IOUT
02070 REAL HHV(4)
02071 REAL HV(4),PKS(4),PNU(4),PKK(4),PPI(4)
02072 REAL PDUM1(4),PDUM2(4),PDUM3(4),PDUM4(4)
02073 REAL*4 RRR(3),RMOD(1)
02074 REAL*8 SWT, SSWT
02075 DATA PI /3.141592653589793238462643/
02076 DATA IWARM/0/
02077
02078 IF(MODE.EQ.-1) THEN
02079
02080 IWARM=1
02081 NEVRAW=0
02082 NEVACC=0
02083 NEVOVR=0
02084 SWT=0
02085 SSWT=0
02086 WTMAX=1E-20
02087 DO 15 I=1,500
02088
02089 JKST=10
02090 CALL DPHSKS(WT,HV,PDUM1,PDUM2,PDUM3,PDUM4,JKST)
02091 IF(WT.GT.WTMAX/1.2) WTMAX=WT*1.2
02092 15 CONTINUE
02093
02094
02095
02096 ELSEIF(MODE.EQ. 0) THEN
02097
02098 IF(IWARM.EQ.0) GOTO 902
02099
02100
02101 DEC1=BRKS
02102 400 CONTINUE
02103 CALL RANMAR(RMOD,1)
02104 IF(RMOD(1).LT.DEC1) THEN
02105 JKST=10
02106 ELSE
02107 JKST=20
02108 ENDIF
02109 CALL DPHSKS(WT,HV,PNU,PKS,PKK,PPI,JKST)
02110 CALL RANMAR(RRR,3)
02111 RN=RRR(1)
02112 IF(WT.GT.WTMAX) NEVOVR=NEVOVR+1
02113 NEVRAW=NEVRAW+1
02114 SWT=SWT+WT
02115 SSWT=SSWT+WT**2
02116 IF(RN*WTMAX.GT.WT) GOTO 400
02117
02118 COSTHE=-1.+2.*RRR(2)
02119 THET=ACOS(COSTHE)
02120 PHI =2*PI*RRR(3)
02121 CALL ROTOR2(THET,PNU,PNU)
02122 CALL ROTOR3( PHI,PNU,PNU)
02123 CALL ROTOR2(THET,PKS,PKS)
02124 CALL ROTOR3( PHI,PKS,PKS)
02125 CALL ROTOR2(THET,PKK,PKK)
02126 CALL ROTOR3(PHI,PKK,PKK)
02127 CALL ROTOR2(THET,PPI,PPI)
02128 CALL ROTOR3( PHI,PPI,PPI)
02129 CALL ROTOR2(THET,HV,HV)
02130 CALL ROTOR3( PHI,HV,HV)
02131 DO 44 I=1,3
02132 44 HHV(I)=-ISGN*HV(I)
02133 NEVACC=NEVACC+1
02134
02135 ELSEIF(MODE.EQ. 1) THEN
02136
02137 IF(NEVRAW.EQ.0) RETURN
02138 PARGAM=SWT/FLOAT(NEVRAW+1)
02139 ERROR=0
02140 IF(NEVRAW.NE.0) ERROR=SQRT(SSWT/SWT**2-1./FLOAT(NEVRAW))
02141 RAT=PARGAM/GAMEL
02142 WRITE(IOUT, 7010) NEVRAW,NEVACC,NEVOVR,PARGAM,RAT,ERROR
02143
02144 GAMPMC(7)=RAT
02145 GAMPER(7)=ERROR
02146
02147 ENDIF
02148
02149 RETURN
02150 7003 FORMAT(///1X,15(5H*****)
02151 $ /,' *', 25X,'******** DADMKS INITIALISATION ********',9X,1H*
02152 $ /,' *',E20.5,5X,'WTMAX = MAXIMUM WEIGHT ',9X,1H*
02153 $ /,1X,15(5H*****)/)
02154 7010 FORMAT(///1X,15(5H*****)
02155 $ /,' *', 25X,'******** DADMKS FINAL REPORT ********',9X,1H*
02156 $ /,' *',I20 ,5X,'NEVRAW = NO. OF K* DECAYS TOTAL ',9X,1H*,
02157 $ /,' *',I20 ,5X,'NEVACC = NO. OF K* DECS. ACCEPTED ',9X,1H*,
02158 $ /,' *',I20 ,5X,'NEVOVR = NO. OF OVERWEIGHTED EVENTS ',9X,1H*
02159 $ /,' *',E20.5,5X,'PARTIAL WTDTH (K* DECAY) IN GEV UNITS ',9X,1H*,
02160 $ /,' *',F20.9,5X,'IN UNITS GFERMI**2*MASS**5/192/PI**3 ',9X,1H*
02161 $ /,' *',F20.8,5X,'RELATIVE ERROR OF PARTIAL WIDTH ',9X,1H*
02162 $ /,1X,15(5H*****)/)
02163 902 WRITE(IOUT, 9020)
02164 9020 FORMAT(' ----- DADMKS: LACK OF INITIALISATION')
02165 STOP
02166 END
02167 SUBROUTINE DPHSKS(DGAMT,HV,PN,PKS,PKK,PPI,JKST)
02168
02169
02170
02171
02172
02173
02174 COMMON / DECPAR / GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
02175 REAL*4 GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
02176 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
02177 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
02178 * ,AMK,AMKZ,AMKST,GAMKST
02179
02180 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
02181 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
02182 * ,AMK,AMKZ,AMKST,GAMKST
02183 REAL HV(4),PT(4),PN(4),PKS(4),PKK(4),PPI(4),QQ(4),RR1(1)
02184 COMPLEX BWIGS
02185 DATA PI /3.141592653589793238462643/
02186
02187 DATA ICONT /0/
02188
02189 PHSPAC=1./2**11/PI**5
02190
02191 PT(1)=0.
02192 PT(2)=0.
02193 PT(3)=0.
02194 PT(4)=AMTAU
02195 CALL RANMAR(RR1,1)
02196
02197 IF(JKST.EQ.10)THEN
02198
02199
02200 AMS1=(AMPI+AMKZ)**2
02201 AMS2=(AMTAU-AMNUTA)**2
02202
02203
02204
02205
02206
02207 ALP1=ATAN((AMS1-AMKST**2)/AMKST/GAMKST)
02208 ALP2=ATAN((AMS2-AMKST**2)/AMKST/GAMKST)
02209 ALP=ALP1+RR1(1)*(ALP2-ALP1)
02210 AMX2=AMKST**2+AMKST*GAMKST*TAN(ALP)
02211 AMX=SQRT(AMX2)
02212 PHSPAC=PHSPAC*((AMX2-AMKST**2)**2+(AMKST*GAMKST)**2)
02213 & /(AMKST*GAMKST)
02214 PHSPAC=PHSPAC*(ALP2-ALP1)
02215
02216
02217 PN(1)=0
02218 PN(2)=0
02219 PN(4)=1./(2*AMTAU)*(AMTAU**2+AMNUTA**2-AMX**2)
02220 PN(3)=-SQRT((PN(4)-AMNUTA)*(PN(4)+AMNUTA))
02221
02222
02223 PKS(1)=0
02224 PKS(2)=0
02225 PKS(4)=1./(2*AMTAU)*(AMTAU**2-AMNUTA**2+AMX**2)
02226 PKS(3)=-PN(3)
02227 PHSPAC=PHSPAC*(4*PI)*(2*PKS(3)/AMTAU)
02228
02229
02230 ENPI=( AMX**2+AMPI**2-AMKZ**2 ) / ( 2*AMX )
02231 PPPI=SQRT((ENPI-AMPI)*(ENPI+AMPI))
02232 PHSPAC=PHSPAC*(4*PI)*(2*PPPI/AMX)
02233
02234 CALL SPHERA(PPPI,PPI)
02235 PPI(4)=ENPI
02236
02237 DO 20 I=1,3
02238 20 PKK(I)=-PPI(I)
02239 PKK(4)=( AMX**2+AMKZ**2-AMPI**2 ) / ( 2*AMX )
02240 EXE=(PKS(4)+PKS(3))/AMX
02241
02242 CALL BOSTR3(EXE,PPI,PPI)
02243 CALL BOSTR3(EXE,PKK,PKK)
02244 DO 30 I=1,4
02245 30 QQ(I)=PPI(I)-PKK(I)
02246
02247 PKSD =PKS(4)*PKS(4)-PKS(3)*PKS(3)-PKS(2)*PKS(2)-PKS(1)*PKS(1)
02248 QQPKS=PKS(4)* QQ(4)-PKS(3)* QQ(3)-PKS(2)* QQ(2)-PKS(1)* QQ(1)
02249 DO 31 I=1,4
02250 31 QQ(I)=QQ(I)-PKS(I)*QQPKS/PKSD
02251
02252 PRODPQ=PT(4)*QQ(4)
02253 PRODNQ=PN(4)*QQ(4)-PN(1)*QQ(1)-PN(2)*QQ(2)-PN(3)*QQ(3)
02254 PRODPN=PT(4)*PN(4)
02255 QQ2= QQ(4)**2-QQ(1)**2-QQ(2)**2-QQ(3)**2
02256 BRAK=(GV**2+GA**2)*(2*PRODPQ*PRODNQ-PRODPN*QQ2)
02257 & +(GV**2-GA**2)*AMTAU*AMNUTA*QQ2
02258
02259 FKS=CABS(BWIGS(AMX2,AMKST,GAMKST))**2
02260 AMPLIT=(GFERMI*SCABIB)**2*BRAK*2*FKS
02261 DGAMT=1/(2.*AMTAU)*AMPLIT*PHSPAC
02262 DO 40 I=1,3
02263 40 HV(I)=2*GV*GA*AMTAU*(2*PRODNQ*QQ(I)-QQ2*PN(I))/BRAK
02264
02265
02266 ELSEIF(JKST.EQ.20)THEN
02267
02268
02269 AMS1=(AMPIZ+AMK)**2
02270 AMS2=(AMTAU-AMNUTA)**2
02271
02272
02273
02274
02275
02276 ALP1=ATAN((AMS1-AMKST**2)/AMKST/GAMKST)
02277 ALP2=ATAN((AMS2-AMKST**2)/AMKST/GAMKST)
02278 ALP=ALP1+RR1(1)*(ALP2-ALP1)
02279 AMX2=AMKST**2+AMKST*GAMKST*TAN(ALP)
02280 AMX=SQRT(AMX2)
02281 PHSPAC=PHSPAC*((AMX2-AMKST**2)**2+(AMKST*GAMKST)**2)
02282 & /(AMKST*GAMKST)
02283 PHSPAC=PHSPAC*(ALP2-ALP1)
02284
02285
02286 PN(1)=0
02287 PN(2)=0
02288 PN(4)=1./(2*AMTAU)*(AMTAU**2+AMNUTA**2-AMX**2)
02289 PN(3)=-SQRT((PN(4)-AMNUTA)*(PN(4)+AMNUTA))
02290
02291 PKS(1)=0
02292 PKS(2)=0
02293 PKS(4)=1./(2*AMTAU)*(AMTAU**2-AMNUTA**2+AMX**2)
02294 PKS(3)=-PN(3)
02295 PHSPAC=PHSPAC*(4*PI)*(2*PKS(3)/AMTAU)
02296
02297
02298 ENPI=( AMX**2+AMPIZ**2-AMK**2 ) / ( 2*AMX )
02299 PPPI=SQRT((ENPI-AMPIZ)*(ENPI+AMPIZ))
02300 PHSPAC=PHSPAC*(4*PI)*(2*PPPI/AMX)
02301
02302 CALL SPHERA(PPPI,PPI)
02303 PPI(4)=ENPI
02304
02305 DO 50 I=1,3
02306 50 PKK(I)=-PPI(I)
02307 PKK(4)=( AMX**2+AMK**2-AMPIZ**2 ) / ( 2*AMX )
02308 EXE=(PKS(4)+PKS(3))/AMX
02309
02310 CALL BOSTR3(EXE,PPI,PPI)
02311 CALL BOSTR3(EXE,PKK,PKK)
02312 DO 60 I=1,4
02313 60 QQ(I)=PKK(I)-PPI(I)
02314
02315 PKSD =PKS(4)*PKS(4)-PKS(3)*PKS(3)-PKS(2)*PKS(2)-PKS(1)*PKS(1)
02316 QQPKS=PKS(4)* QQ(4)-PKS(3)* QQ(3)-PKS(2)* QQ(2)-PKS(1)* QQ(1)
02317 DO 61 I=1,4
02318 61 QQ(I)=QQ(I)-PKS(I)*QQPKS/PKSD
02319
02320 PRODPQ=PT(4)*QQ(4)
02321 PRODNQ=PN(4)*QQ(4)-PN(1)*QQ(1)-PN(2)*QQ(2)-PN(3)*QQ(3)
02322 PRODPN=PT(4)*PN(4)
02323 QQ2= QQ(4)**2-QQ(1)**2-QQ(2)**2-QQ(3)**2
02324 BRAK=(GV**2+GA**2)*(2*PRODPQ*PRODNQ-PRODPN*QQ2)
02325 & +(GV**2-GA**2)*AMTAU*AMNUTA*QQ2
02326
02327 FKS=CABS(BWIGS(AMX2,AMKST,GAMKST))**2
02328 AMPLIT=(GFERMI*SCABIB)**2*BRAK*2*FKS
02329 DGAMT=1/(2.*AMTAU)*AMPLIT*PHSPAC
02330 DO 70 I=1,3
02331 70 HV(I)=2*GV*GA*AMTAU*(2*PRODNQ*QQ(I)-QQ2*PN(I))/BRAK
02332 ENDIF
02333 RETURN
02334 END
02335
02336
02337
02338 SUBROUTINE DPHNPI(DGAMT,HVX,PNX,PRX,PPIX,JNPI)
02339
02340
02341
02342
02343 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
02344 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
02345 * ,AMK,AMKZ,AMKST,GAMKST
02346
02347 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
02348 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
02349 * ,AMK,AMKZ,AMKST,GAMKST
02350 COMMON / DECPAR / GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
02351 REAL*4 GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
02352 PARAMETER (NMODE=15,NM1=0,NM2=1,NM3=8,NM4=2,NM5=1,NM6=3)
02353 COMMON / TAUDCD /IDFFIN(9,NMODE),MULPIK(NMODE)
02354 & ,NAMES
02355 CHARACTER NAMES(NMODE)*31
02356 REAL*8 WETMAX(20)
02357
02358 REAL*8 PN(4),PR(4),PPI(4,9),HV(4)
02359 REAL*4 PNX(4),PRX(4),PPIX(4,9),HVX(4)
02360 REAL*8 PV(5,9),PT(4),UE(3),BE(3)
02361 REAL*8 PAWT,AMX,AMS1,AMS2,PA,PHS,PHSMAX,PMIN,PMAX
02362
02363 REAL*8 PHSPAC
02364
02365 REAL*8 GAM,BEP,PHI,A,B,C
02366 REAL*8 AMPIK
02367 REAL*4 RRR(9),RRX(2),RN(1),RR2(1)
02368
02369 DATA PI /3.141592653589793238462643/
02370 DATA WETMAX /20*1D-15/
02371
02372
02373
02374 PAWT(A,B,C)=
02375 $ SQRT(MAX(0.D0,(A**2-(B+C)**2)*(A**2-(B-C)**2)))/(2.D0*A)
02376
02377 AMPIK(I,J)=DCDMAS(IDFFIN(I,J))
02378
02379
02380 IF ((JNPI.LE.0).OR.JNPI.GT.20) THEN
02381 WRITE(6,*) 'JNPI OUTSIDE RANGE DEFINED BY WETMAX; JNPI=',JNPI
02382 STOP
02383 ENDIF
02384
02385
02386 PT(1)=0.
02387 PT(2)=0.
02388 PT(3)=0.
02389 PT(4)=AMTAU
02390
02391 500 CONTINUE
02392
02393 ND=MULPIK(JNPI)
02394 PS=0.
02395 PHSPAC = 1./2.**5 /PI**2
02396 DO 4 I=1,ND
02397 4 PS =PS+AMPIK(I,JNPI)
02398 CALL RANMAR(RR2,1)
02399 AMS1=PS**2
02400 AMS2=(AMTAU-AMNUTA)**2
02401
02402
02403 AMX2=AMS1+ RR2(1)*(AMS2-AMS1)
02404 AMX =SQRT(AMX2)
02405 AMW =AMX
02406 PHSPAC=PHSPAC * (AMS2-AMS1)
02407
02408
02409 PN(1)=0
02410 PN(2)=0
02411 PN(4)=1./(2*AMTAU)*(AMTAU**2+AMNUTA**2-AMX2)
02412 PN(3)=-SQRT((PN(4)-AMNUTA)*(PN(4)+AMNUTA))
02413
02414 PR(1)=0
02415 PR(2)=0
02416 PR(4)=1./(2*AMTAU)*(AMTAU**2-AMNUTA**2+AMX2)
02417 PR(3)=-PN(3)
02418 PHSPAC=PHSPAC * (4.*PI) * (2.*PR(3)/AMTAU)
02419
02420
02421
02422
02423 PXQ=AMTAU*PR(4)
02424 PXN=AMTAU*PN(4)
02425 QXN=PR(4)*PN(4)-PR(1)*PN(1)-PR(2)*PN(2)-PR(3)*PN(3)
02426
02427
02428 BRAK=2*(GV**2+GA**2)*(2*PXQ*QXN+AMX2*PXN)
02429 & -6*(GV**2-GA**2)*AMTAU*AMNUTA*AMX2
02430
02431
02432 AMPLIT=CCABIB**2*GFERMI**2/2. * BRAK * AMX2*SIGEE(AMX2,JNPI)
02433 DGAMT=1./(2.*AMTAU)*AMPLIT*PHSPAC
02434
02435
02436 PHSMAX = 1.
02437 DO 200 I=1,4
02438 200 PV(I,1)=PR(I)
02439 PV(5,1)=AMW
02440 PV(5,ND)=AMPIK(ND,JNPI)
02441
02442 PMAX=AMW-PS+AMPIK(ND,JNPI)
02443 PMIN=.0
02444 DO 220 IL=ND-1,1,-1
02445 PMAX=PMAX+AMPIK(IL,JNPI)
02446 PMIN=PMIN+AMPIK(IL+1,JNPI)
02447 220 PHSMAX=PHSMAX*PAWT(PMAX,PMIN,AMPIK(IL,JNPI))/PMAX
02448
02449
02450 AMX=AMW
02451 DO 222 IL=1,ND-2
02452 AMS1=.0
02453 DO 223 JL=IL+1,ND
02454 223 AMS1=AMS1+AMPIK(JL,JNPI)
02455 AMS1=AMS1**2
02456 AMX =(AMX-AMPIK(IL,JNPI))
02457 AMS2=(AMX)**2
02458 PHSMAX=PHSMAX * (AMS2-AMS1)
02459 222 CONTINUE
02460 NCONT=0
02461 100 CONTINUE
02462 NCONT=NCONT+1
02463
02464 PHS=1.D0
02465 PHSPAC = 1./2.**(6*ND-7) /PI**(3*ND-4)
02466 AMX=AMW
02467 CALL RANMAR(RRR,ND-2)
02468 DO 230 IL=1,ND-2
02469 AMS1=.0D0
02470 DO 231 JL=IL+1,ND
02471 231 AMS1=AMS1+AMPIK(JL,JNPI)
02472 AMS1=AMS1**2
02473 AMS2=(AMX-AMPIK(IL,JNPI))**2
02474 RR1=RRR(IL)
02475 AMX2=AMS1+ RR1*(AMS2-AMS1)
02476 AMX=SQRT(AMX2)
02477 PV(5,IL+1)=AMX
02478 PHSPAC=PHSPAC * (AMS2-AMS1)
02479
02480 PHS=PHS* (AMS2-AMS1)
02481 PA=PAWT(PV(5,IL),PV(5,IL+1),AMPIK(IL,JNPI))
02482 PHS =PHS *PA/PV(5,IL)
02483 230 CONTINUE
02484 PA=PAWT(PV(5,ND-1),AMPIK(ND-1,JNPI),AMPIK(ND,JNPI))
02485 PHS =PHS *PA/PV(5,ND-1)
02486 CALL RANMAR(RN,1)
02487 WETMAX(JNPI)=1.2D0*MAX(WETMAX(JNPI)/1.2D0,PHS/PHSMAX)
02488 IF (NCONT.EQ.500 000) THEN
02489 XNPI=0.0
02490 DO KK=1,ND
02491 XNPI=XNPI+AMPIK(KK,JNPI)
02492 ENDDO
02493 WRITE(6,*) 'ROUNDING INSTABILITY IN DPHNPI ?'
02494 WRITE(6,*) 'AMW=',AMW,'XNPI=',XNPI
02495 WRITE(6,*) 'IF =AMW= IS NEARLY EQUAL =XNPI= THAT IS IT'
02496 WRITE(6,*) 'PHS=',PHS,'PHSMAX=',PHSMAX
02497 GOTO 500
02498 ENDIF
02499 IF(RN(1)*PHSMAX*WETMAX(JNPI).GT.PHS) GO TO 100
02500
02501 280 DO 300 IL=1,ND-1
02502 PA=PAWT(PV(5,IL),PV(5,IL+1),AMPIK(IL,JNPI))
02503 CALL RANMAR(RRX,2)
02504 UE(3)=2.*RRX(1)-1.
02505 PHI=2.*PI*RRX(2)
02506 UE(1)=SQRT(1.D0-UE(3)**2)*COS(PHI)
02507 UE(2)=SQRT(1.D0-UE(3)**2)*SIN(PHI)
02508 DO 290 J=1,3
02509 PPI(J,IL)=PA*UE(J)
02510 290 PV(J,IL+1)=-PA*UE(J)
02511 PPI(4,IL)=SQRT(PA**2+AMPIK(IL,JNPI)**2)
02512 PV(4,IL+1)=SQRT(PA**2+PV(5,IL+1)**2)
02513 PHSPAC=PHSPAC *(4.*PI)*(2.*PA/PV(5,IL))
02514 300 CONTINUE
02515
02516 DO 310 J=1,4
02517 310 PPI(J,ND)=PV(J,ND)
02518 DO 340 IL=ND-1,1,-1
02519 DO 320 J=1,3
02520 320 BE(J)=PV(J,IL)/PV(4,IL)
02521 GAM=PV(4,IL)/PV(5,IL)
02522 DO 340 I=IL,ND
02523 BEP=BE(1)*PPI(1,I)+BE(2)*PPI(2,I)+BE(3)*PPI(3,I)
02524 DO 330 J=1,3
02525 330 PPI(J,I)=PPI(J,I)+GAM*(GAM*BEP/(1.D0+GAM)+PPI(4,I))*BE(J)
02526 PPI(4,I)=GAM*(PPI(4,I)+BEP)
02527 340 CONTINUE
02528
02529 HV(4)=1.
02530 HV(3)=0.
02531 HV(2)=0.
02532 HV(1)=0.
02533 DO K=1,4
02534 PNX(K)=PN(K)
02535 PRX(K)=PR(K)
02536 HVX(K)=HV(K)
02537 DO L=1,ND
02538 PPIX(K,L)=PPI(K,L)
02539 ENDDO
02540 ENDDO
02541 RETURN
02542 END
02543 FUNCTION SIGEE(Q2,JNP)
02544
02545
02546
02547
02548
02549
02550
02551
02552
02553
02554
02555
02556
02557
02558
02559
02560
02561
02562 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
02563 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
02564 * ,AMK,AMKZ,AMKST,GAMKST
02565
02566 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
02567 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
02568 * ,AMK,AMKZ,AMKST,GAMKST
02569 REAL*4 DATSIG(17,6)
02570
02571 DATA DATSIG/
02572 1 7.40,12.00,16.15,21.25,24.90,29.55,34.15,37.40,37.85,37.40,
02573 2 36.00,33.25,30.50,27.70,24.50,21.25,18.90,
02574 3 1.24, 2.50, 3.70, 5.40, 7.45,10.75,14.50,18.20,22.30,28.90,
02575 4 29.35,25.60,22.30,18.60,14.05,11.60, 9.10,
02576 5 17*.0,
02577 6 17*.0,
02578 7 9*.0,.65,1.25,2.20,3.15,5.00,5.75,7.80,8.25,
02579 8 17*.0/
02580 DATA SIG0 / 86.8 /
02581 DATA PI /3.141592653589793238462643/
02582 DATA INIT / 0 /
02583
02584 JNPI=JNP
02585 IF(JNP.EQ.4) JNPI=3
02586 IF(JNP.EQ.3) JNPI=4
02587 IF(INIT.EQ.0) THEN
02588 INIT=1
02589
02590 IF (AMPI.LT.0.139) AMPI = 0.1395675
02591 AMPI2=AMPI**2
02592 FPI = .943*AMPI
02593 DO 100 I=1,17
02594 DATSIG(I,2) = DATSIG(I,2)/2.
02595 DATSIG(I,1) = DATSIG(I,1) + DATSIG(I,2)
02596 S = 1.025+(I-1)*.05
02597 FACT=0.
02598 S2=S**2
02599 DO 200 J=1,17
02600 T= 1.025+(J-1)*.05
02601 IF(T . GT. S-AMPI ) GO TO 201
02602 T2=T**2
02603 FACT=(T2/S2)**2*SQRT((S2-T2-AMPI2)**2-4.*T2*AMPI2)/S2 *2.*T*.05
02604 FACT = FACT * (DATSIG(J,1)+DATSIG(J+1,1))
02605 200 DATSIG(I,3) = DATSIG(I,3) + FACT
02606 201 DATSIG(I,3) = DATSIG(I,3) /(2*PI*FPI)**2
02607 DATSIG(I,4) = DATSIG(I,3)
02608 DATSIG(I,6) = DATSIG(I,5)
02609 100 CONTINUE
02610
02611 1000 FORMAT(///1X,' EE SIGMA USED IN MULTIPI DECAYS'/
02612 % (17F7.2/))
02613 ENDIF
02614 Q=SQRT(Q2)
02615 QMIN=1.
02616 IF(Q.LT.QMIN) THEN
02617 SIGEE=DATSIG(1,JNPI)+
02618 & (DATSIG(2,JNPI)-DATSIG(1,JNPI))*(Q-1.)/.05
02619 ELSEIF(Q.LT.1.8) THEN
02620 DO 1 I=1,16
02621 QMAX = QMIN + .05
02622 IF(Q.LT.QMAX) GO TO 2
02623 QMIN = QMIN + .05
02624 1 CONTINUE
02625 2 SIGEE=DATSIG(I,JNPI)+
02626 & (DATSIG(I+1,JNPI)-DATSIG(I,JNPI)) * (Q-QMIN)/.05
02627 ELSEIF(Q.GT.1.8) THEN
02628 SIGEE=DATSIG(17,JNPI)+
02629 & (DATSIG(17,JNPI)-DATSIG(16,JNPI)) * (Q-1.8)/.05
02630 ENDIF
02631 IF(SIGEE.LT..0) SIGEE=0.
02632
02633 SIGEE = SIGEE/(6.*PI**2*SIG0)
02634
02635 RETURN
02636 END
02637
02638 FUNCTION SIGOLD(Q2,JNPI)
02639
02640
02641
02642
02643
02644
02645
02646
02647
02648
02649
02650
02651
02652
02653
02654
02655
02656
02657 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
02658 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
02659 * ,AMK,AMKZ,AMKST,GAMKST
02660
02661 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
02662 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
02663 * ,AMK,AMKZ,AMKST,GAMKST
02664 REAL*4 DATSIG(17,4)
02665
02666 DATA DATSIG/
02667 1 7.40,12.00,16.15,21.25,24.90,29.55,34.15,37.40,37.85,37.40,
02668 2 36.00,33.25,30.50,27.70,24.50,21.25,18.90,
02669 3 1.24, 2.50, 3.70, 5.40, 7.45,10.75,14.50,18.20,22.30,28.90,
02670 4 29.35,25.60,22.30,18.60,14.05,11.60, 9.10,
02671 5 17*.0,
02672 6 9*.0,.65,1.25,2.20,3.15,5.00,5.75,7.80,8.25/
02673 DATA SIG0 / 86.8 /
02674 DATA PI /3.141592653589793238462643/
02675 DATA INIT / 0 /
02676
02677 IF(INIT.EQ.0) THEN
02678 INIT=1
02679 AMPI2=AMPI**2
02680 FPI = .943*AMPI
02681 DO 100 I=1,17
02682 DATSIG(I,2) = DATSIG(I,2)/2.
02683 DATSIG(I,1) = DATSIG(I,1) + DATSIG(I,2)
02684 S = 1.025+(I-1)*.05
02685 FACT=0.
02686 S2=S**2
02687 DO 200 J=1,17
02688 T= 1.025+(J-1)*.05
02689 IF(T . GT. S-AMPI ) GO TO 201
02690 T2=T**2
02691 FACT=(T2/S2)**2*SQRT((S2-T2-AMPI2)**2-4.*T2*AMPI2)/S2 *2.*T*.05
02692 FACT = FACT * (DATSIG(J,1)+DATSIG(J+1,1))
02693 200 DATSIG(I,3) = DATSIG(I,3) + FACT
02694 201 DATSIG(I,3) = DATSIG(I,3) /(2*PI*FPI)**2
02695 100 CONTINUE
02696
02697 1000 FORMAT(///1X,' EE SIGMA USED IN MULTIPI DECAYS'/
02698 % (17F7.2/))
02699 ENDIF
02700 Q=SQRT(Q2)
02701 QMIN=1.
02702 IF(Q.LT.QMIN) THEN
02703 SIGEE=DATSIG(1,JNPI)+
02704 & (DATSIG(2,JNPI)-DATSIG(1,JNPI))*(Q-1.)/.05
02705 ELSEIF(Q.LT.1.8) THEN
02706 DO 1 I=1,16
02707 QMAX = QMIN + .05
02708 IF(Q.LT.QMAX) GO TO 2
02709 QMIN = QMIN + .05
02710 1 CONTINUE
02711 2 SIGEE=DATSIG(I,JNPI)+
02712 & (DATSIG(I+1,JNPI)-DATSIG(I,JNPI)) * (Q-QMIN)/.05
02713 ELSEIF(Q.GT.1.8) THEN
02714 SIGEE=DATSIG(17,JNPI)+
02715 & (DATSIG(17,JNPI)-DATSIG(16,JNPI)) * (Q-1.8)/.05
02716 ENDIF
02717 IF(SIGEE.LT..0) SIGEE=0.
02718
02719 SIGEE = SIGEE/(6.*PI**2*SIG0)
02720 SIGOLD=SIGEE
02721
02722 RETURN
02723 END
02724 SUBROUTINE DPHSPK(DGAMT,HV,PN,PAA,PNPI,JAA)
02725
02726
02727
02728
02729 PARAMETER (NMODE=15,NM1=0,NM2=1,NM3=8,NM4=2,NM5=1,NM6=3)
02730 COMMON / TAUDCD /IDFFIN(9,NMODE),MULPIK(NMODE)
02731 & ,NAMES
02732 CHARACTER NAMES(NMODE)*31
02733
02734 REAL HV(4),PN(4),PAA(4),PIM1(4),PIM2(4),PIPL(4),PNPI(4,9)
02735
02736 MNUM=JAA
02737
02738 KEYT=4
02739 IF(JAA.EQ.7) KEYT=3
02740
02741 AMP1=DCDMAS(IDFFIN(1,JAA+NM4+NM5+NM6))
02742 AMP2=DCDMAS(IDFFIN(2,JAA+NM4+NM5+NM6))
02743 AMP3=DCDMAS(IDFFIN(3,JAA+NM4+NM5+NM6))
02744 CALL
02745 $ DPHTRE(DGAMT,HV,PN,PAA,PIM1,AMP1,PIM2,AMP2,PIPL,AMP3,KEYT,MNUM)
02746 DO I=1,4
02747 PNPI(I,1)=PIM1(I)
02748 PNPI(I,2)=PIM2(I)
02749 PNPI(I,3)=PIPL(I)
02750 ENDDO
02751 END
02752
02753
02754
02755
02756 SUBROUTINE
02757 $ DPHTRE(DGAMT,HV,PN,PAA,PIM1,AMPA,PIM2,AMPB,PIPL,AMP3,KEYT,MNUM)
02758
02759
02760
02761
02762
02763
02764
02765
02766
02767
02768
02769
02770
02771
02772
02773
02774 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
02775 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
02776 * ,AMK,AMKZ,AMKST,GAMKST
02777
02778 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
02779 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
02780 * ,AMK,AMKZ,AMKST,GAMKST
02781 COMMON / DECPAR / GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
02782 REAL*4 GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
02783 REAL HV(4),PT(4),PN(4),PAA(4),PIM1(4),PIM2(4),PIPL(4)
02784 REAL PR(4)
02785 REAL*4 RRR(5)
02786 DATA PI /3.141592653589793238462643/
02787 DATA ICONT /0/
02788 XLAM(X,Y,Z)=SQRT(ABS((X-Y-Z)**2-4.0*Y*Z))
02789
02790
02791
02792
02793 PHSPAC=1./2**17/PI**8
02794
02795 PT(1)=0.
02796 PT(2)=0.
02797 PT(3)=0.
02798 PT(4)=AMTAU
02799
02800 CALL RANMAR(RRR,5)
02801 RR=RRR(5)
02802
02803 CALL CHOICE(MNUM,RR,ICHAN,PROB1,PROB2,PROB3,
02804 $ AMRX,GAMRX,AMRA,GAMRA,AMRB,GAMRB)
02805 IF (ICHAN.EQ.1) THEN
02806 AMP1=AMPB
02807 AMP2=AMPA
02808 ELSEIF (ICHAN.EQ.2) THEN
02809 AMP1=AMPA
02810 AMP2=AMPB
02811 ELSE
02812 AMP1=AMPB
02813 AMP2=AMPA
02814 ENDIF
02815
02816 RR1=RRR(1)
02817 AMS1=(AMP1+AMP2+AMP3)**2
02818 AMS2=(AMTAU-AMNUTA)**2
02819
02820 ALP1=ATAN((AMS1-AMRX**2)/AMRX/GAMRX)
02821 ALP2=ATAN((AMS2-AMRX**2)/AMRX/GAMRX)
02822 ALP=ALP1+RR1*(ALP2-ALP1)
02823 AM3SQ =AMRX**2+AMRX*GAMRX*TAN(ALP)
02824 AM3 =SQRT(AM3SQ)
02825 PHSPAC=PHSPAC*((AM3SQ-AMRX**2)**2+(AMRX*GAMRX)**2)/(AMRX*GAMRX)
02826 PHSPAC=PHSPAC*(ALP2-ALP1)
02827
02828 RR2=RRR(2)
02829 AMS1=(AMP2+AMP3)**2
02830 AMS2=(AM3-AMP1)**2
02831 IF (ICHAN.LE.2) THEN
02832
02833 ALP1=ATAN((AMS1-AMRA**2)/AMRA/GAMRA)
02834 ALP2=ATAN((AMS2-AMRA**2)/AMRA/GAMRA)
02835 ALP=ALP1+RR2*(ALP2-ALP1)
02836 AM2SQ =AMRA**2+AMRA*GAMRA*TAN(ALP)
02837 AM2 =SQRT(AM2SQ)
02838
02839
02840
02841
02842 ELSE
02843
02844 AM2SQ=AMS1+ RR2*(AMS2-AMS1)
02845 AM2 =SQRT(AM2SQ)
02846 PHF0=(AMS2-AMS1)
02847 ENDIF
02848
02849 ENQ1=(AM2SQ-AMP2**2+AMP3**2)/(2*AM2)
02850 ENQ2=(AM2SQ+AMP2**2-AMP3**2)/(2*AM2)
02851 PPI= ENQ1**2-AMP3**2
02852 PPPI=SQRT(ABS(ENQ1**2-AMP3**2))
02853
02854 PHF1=(4*PI)*(2*PPPI/AM2)
02855
02856 CALL SPHERA(PPPI,PIPL)
02857 PIPL(4)=ENQ1
02858
02859 DO 30 I=1,3
02860 30 PIM1(I)=-PIPL(I)
02861 PIM1(4)=ENQ2
02862
02863
02864 PR(1)=0
02865 PR(2)=0
02866 PR(4)=1./(2*AM3)*(AM3**2+AM2**2-AMP1**2)
02867 PR(3)= SQRT(ABS(PR(4)**2-AM2**2))
02868 PPI = PR(4)**2-AM2**2
02869
02870 PIM2(1)=0
02871 PIM2(2)=0
02872 PIM2(4)=1./(2*AM3)*(AM3**2-AM2**2+AMP1**2)
02873 PIM2(3)=-PR(3)
02874 PHF2=(4*PI)*(2*PR(3)/AM3)
02875
02876 EXE=(PR(4)+PR(3))/AM2
02877 CALL BOSTR3(EXE,PIPL,PIPL)
02878 CALL BOSTR3(EXE,PIM1,PIM1)
02879 RR3=RRR(3)
02880 RR4=RRR(4)
02881
02882 THET =ACOS(-1.+2*RR3)
02883 PHI = 2*PI*RR4
02884 CALL ROTPOL(THET,PHI,PIPL)
02885 CALL ROTPOL(THET,PHI,PIM1)
02886 CALL ROTPOL(THET,PHI,PIM2)
02887 CALL ROTPOL(THET,PHI,PR)
02888
02889
02890
02891 PAA(1)=0
02892 PAA(2)=0
02893 PAA(4)=1./(2*AMTAU)*(AMTAU**2-AMNUTA**2+AM3**2)
02894 PAA(3)= SQRT(ABS(PAA(4)**2-AM3**2))
02895 PPI = PAA(4)**2-AM3**2
02896 PHSPAC=PHSPAC*(4*PI)*(2*PAA(3)/AMTAU)
02897
02898 PN(1)=0
02899 PN(2)=0
02900 PN(4)=1./(2*AMTAU)*(AMTAU**2+AMNUTA**2-AM3**2)
02901 PN(3)=-PAA(3)
02902
02903
02904 AMS1=(AMP2+AMP3)**2
02905 AMS2=(AM3-AMP1)**2
02906 ALP1=ATAN((AMS1-AMRA**2)/AMRA/GAMRA)
02907 ALP2=ATAN((AMS2-AMRA**2)/AMRA/GAMRA)
02908 XPRO = (PIM1(3)+PIPL(3))**2
02909 $ +(PIM1(2)+PIPL(2))**2+(PIM1(1)+PIPL(1))**2
02910 AM2SQ=-XPRO+(PIM1(4)+PIPL(4))**2
02911
02912 FF1 = ((AM2SQ-AMRA**2)**2+(AMRA*GAMRA)**2)/(AMRA*GAMRA)
02913 FF1 =FF1 *(ALP2-ALP1)
02914
02915 GG1 = (4*PI)*(XLAM(AM2SQ,AMP2**2,AMP3**2)/AM2SQ)
02916
02917 GG1 =GG1 *(4*PI)*SQRT(4*XPRO/AM3SQ)
02918 XJAJE=GG1*(AMS2-AMS1)
02919
02920 AMS1=(AMP1+AMP3)**2
02921 AMS2=(AM3-AMP2)**2
02922 ALP1=ATAN((AMS1-AMRB**2)/AMRB/GAMRB)
02923 ALP2=ATAN((AMS2-AMRB**2)/AMRB/GAMRB)
02924 XPRO = (PIM2(3)+PIPL(3))**2
02925 $ +(PIM2(2)+PIPL(2))**2+(PIM2(1)+PIPL(1))**2
02926 AM2SQ=-XPRO+(PIM2(4)+PIPL(4))**2
02927 FF2 = ((AM2SQ-AMRB**2)**2+(AMRB*GAMRB)**2)/(AMRB*GAMRB)
02928 FF2 =FF2 *(ALP2-ALP1)
02929 GG2 = (4*PI)*(XLAM(AM2SQ,AMP1**2,AMP3**2)/AM2SQ)
02930 GG2 =GG2 *(4*PI)*SQRT(4*XPRO/AM3SQ)
02931 XJADW=GG2*(AMS2-AMS1)
02932
02933 A1=0.0
02934 A2=0.0
02935 A3=0.0
02936 XJAC1=FF1*GG1
02937 XJAC2=FF2*GG2
02938 IF (ICHAN.EQ.2) THEN
02939 XJAC3=XJADW
02940 ELSE
02941 XJAC3=XJAJE
02942 ENDIF
02943 IF (XJAC1.NE.0.0) A1=PROB1/XJAC1
02944 IF (XJAC2.NE.0.0) A2=PROB2/XJAC2
02945 IF (XJAC3.NE.0.0) A3=PROB3/XJAC3
02946
02947 IF (A1+A2+A3.NE.0.0) THEN
02948 PHSPAC=PHSPAC/(A1+A2+A3)
02949 ELSE
02950 PHSPAC=0.0
02951 ENDIF
02952 IF(ICHAN.EQ.2) THEN
02953 DO 70 I=1,4
02954 X=PIM1(I)
02955 PIM1(I)=PIM2(I)
02956 70 PIM2(I)=X
02957 ENDIF
02958
02959
02960 EXE=(PAA(4)+PAA(3))/AM3
02961 CALL BOSTR3(EXE,PIPL,PIPL)
02962 CALL BOSTR3(EXE,PIM1,PIM1)
02963 CALL BOSTR3(EXE,PIM2,PIM2)
02964 CALL BOSTR3(EXE,PR,PR)
02965
02966 IF (MNUM.EQ.8) THEN
02967 CALL DAMPOG(PT,PN,PIM1,PIM2,PIPL,AMPLIT,HV)
02968
02969
02970 ELSE
02971 CALL DAMPPK(MNUM,PT,PN,PIM1,PIM2,PIPL,AMPLIT,HV)
02972 ENDIF
02973 IF (KEYT.EQ.1.OR.KEYT.EQ.2) THEN
02974
02975
02976 PHSPAC=PHSPAC*2.0
02977 PHSPAC=PHSPAC/2.
02978 ENDIF
02979 DGAMT=1/(2.*AMTAU)*AMPLIT*PHSPAC
02980 END
02981 SUBROUTINE DAMPAA(PT,PN,PIM1,PIM2,PIPL,AMPLIT,HV)
02982
02983
02984
02985
02986
02987
02988
02989
02990
02991 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
02992 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
02993 * ,AMK,AMKZ,AMKST,GAMKST
02994
02995 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
02996 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
02997 * ,AMK,AMKZ,AMKST,GAMKST
02998 COMMON / DECPAR / GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
02999 REAL*4 GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
03000 COMMON /TESTA1/ KEYA1
03001 REAL HV(4),PT(4),PN(4),PIM1(4),PIM2(4),PIPL(4)
03002 REAL PAA(4),VEC1(4),VEC2(4)
03003 REAL PIVEC(4),PIAKS(4),HVM(4)
03004 COMPLEX BWIGN,HADCUR(4),FPIK
03005 DATA ICONT /1/
03006
03007
03008
03009 DATA FPI /93.3E-3/
03010
03011 BWIGN(XM,AM,GAMMA)=1./CMPLX(XM**2-AM**2,GAMMA*AM)
03012
03013
03014 DO 10 I=1,4
03015 10 PAA(I)=PIM1(I)+PIM2(I)+PIPL(I)
03016
03017 XMAA =SQRT(ABS(PAA(4)**2-PAA(3)**2-PAA(2)**2-PAA(1)**2))
03018 XMRO1 =SQRT(ABS((PIPL(4)+PIM1(4))**2-(PIPL(1)+PIM1(1))**2
03019 $ -(PIPL(2)+PIM1(2))**2-(PIPL(3)+PIM1(3))**2))
03020 XMRO2 =SQRT(ABS((PIPL(4)+PIM2(4))**2-(PIPL(1)+PIM2(1))**2
03021 $ -(PIPL(2)+PIM2(2))**2-(PIPL(3)+PIM2(3))**2))
03022
03023 PROD1 =PAA(4)*(PIM1(4)-PIPL(4))-PAA(1)*(PIM1(1)-PIPL(1))
03024 $ -PAA(2)*(PIM1(2)-PIPL(2))-PAA(3)*(PIM1(3)-PIPL(3))
03025 PROD2 =PAA(4)*(PIM2(4)-PIPL(4))-PAA(1)*(PIM2(1)-PIPL(1))
03026 $ -PAA(2)*(PIM2(2)-PIPL(2))-PAA(3)*(PIM2(3)-PIPL(3))
03027 DO 40 I=1,4
03028 VEC1(I)= PIM1(I)-PIPL(I) -PAA(I)*PROD1/XMAA**2
03029 40 VEC2(I)= PIM2(I)-PIPL(I) -PAA(I)*PROD2/XMAA**2
03030
03031 IF (KEYA1.EQ.1) THEN
03032 FA1=9.87
03033 FAROPI=1.0
03034 FRO2PI=1.0
03035 FNORM=FA1/SQRT(2.)*FAROPI*FRO2PI
03036 DO 45 I=1,4
03037 HADCUR(I)= CMPLX(FNORM) *AMA1**2*BWIGN(XMAA,AMA1,GAMA1)
03038 $ *(CMPLX(VEC1(I))*AMRO**2*BWIGN(XMRO1,AMRO,GAMRO)
03039 $ +CMPLX(VEC2(I))*AMRO**2*BWIGN(XMRO2,AMRO,GAMRO))
03040 45 CONTINUE
03041 ELSE
03042 FNORM=2.0*SQRT(2.)/3.0/FPI
03043 GAMAX=GAMA1*GFUN(XMAA**2)/GFUN(AMA1**2)
03044 DO 46 I=1,4
03045 HADCUR(I)= CMPLX(FNORM) *AMA1**2*BWIGN(XMAA,AMA1,GAMAX)
03046 $ *(CMPLX(VEC1(I))*FPIK(XMRO1)
03047 $ +CMPLX(VEC2(I))*FPIK(XMRO2))
03048 46 CONTINUE
03049 ENDIF
03050
03051
03052 CALL CLVEC(HADCUR,PN,PIVEC)
03053 CALL CLAXI(HADCUR,PN,PIAKS)
03054 CALL CLNUT(HADCUR,BRAKM,HVM)
03055
03056 BRAK= (GV**2+GA**2)*PT(4)*PIVEC(4) +2.*GV*GA*PT(4)*PIAKS(4)
03057 & +2.*(GV**2-GA**2)*AMNUTA*AMTAU*BRAKM
03058 AMPLIT=(GFERMI*CCABIB)**2*BRAK/2.
03059
03060
03061
03062 DO 90 I=1,3
03063 HV(I)=-(AMTAU*((GV**2+GA**2)*PIAKS(I)+2.*GV*GA*PIVEC(I)))
03064 & +(GV**2-GA**2)*AMNUTA*AMTAU*HVM(I)
03065
03066 HV(I)=-HV(I)/BRAK
03067 90 CONTINUE
03068 END
03069
03070 FUNCTION GFUN(QKWA)
03071
03072
03073
03074 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
03075 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
03076 * ,AMK,AMKZ,AMKST,GAMKST
03077
03078 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
03079 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
03080 * ,AMK,AMKZ,AMKST,GAMKST
03081
03082 IF (QKWA.LT.(AMRO+AMPI)**2) THEN
03083 GFUN=4.1*(QKWA-9*AMPIZ**2)**3
03084 $ *(1.-3.3*(QKWA-9*AMPIZ**2)+5.8*(QKWA-9*AMPIZ**2)**2)
03085 ELSE
03086 GFUN=QKWA*(1.623+10.38/QKWA-9.32/QKWA**2+0.65/QKWA**3)
03087 ENDIF
03088 END
03089 COMPLEX FUNCTION BWIGS(S,M,G)
03090
03091
03092
03093 REAL S,M,G
03094 REAL PI,PIM,QS,QM,W,GS,MK
03095
03096 REAL PM, PG, PBETA
03097 COMPLEX BW,BWP
03098 DATA INIT /0/
03099 P(A,B,C)=SQRT(ABS(ABS(((A+B-C)**2-4.*A*B)/4./A)
03100 $ +(((A+B-C)**2-4.*A*B)/4./A))/2.0)
03101
03102 IF (INIT.EQ.0) THEN
03103 INIT=1
03104 PI=3.141592654
03105 PIM=.139
03106 MK=.493667
03107
03108 PM = PKORB(1,16)
03109 PG = PKORB(2,16)
03110 PBETA = PKORB(3,16)
03111
03112 ENDIF
03113 QS=P(S,PIM**2,MK**2)
03114 QM=P(M**2,PIM**2,MK**2)
03115 W=SQRT(S)
03116 GS=G*(M/W)*(QS/QM)**3
03117 BW=M**2/CMPLX(M**2-S,-M*GS)
03118 QPM=P(PM**2,PIM**2,MK**2)
03119 G1=PG*(PM/W)*(QS/QPM)**3
03120 BWP=PM**2/CMPLX(PM**2-S,-PM*G1)
03121 BWIGS= (BW+PBETA*BWP)/(1+PBETA)
03122 RETURN
03123 END
03124 COMPLEX FUNCTION BWIG(S,M,G)
03125
03126
03127
03128 REAL S,M,G
03129 REAL PI,PIM,QS,QM,W,GS
03130 DATA INIT /0/
03131
03132 IF (INIT.EQ.0) THEN
03133 INIT=1
03134 PI=3.141592654
03135 PIM=.139
03136
03137 ENDIF
03138 IF (S.GT.4.*PIM**2) THEN
03139 QS=SQRT(ABS(ABS(S/4.-PIM**2)+(S/4.-PIM**2))/2.0)
03140 QM=SQRT(M**2/4.-PIM**2)
03141 W=SQRT(S)
03142 GS=G*(M/W)*(QS/QM)**3
03143 ELSE
03144 GS=0.0
03145 ENDIF
03146 BWIG=M**2/CMPLX(M**2-S,-M*GS)
03147 RETURN
03148 END
03149 COMPLEX FUNCTION FPIK(W)
03150
03151
03152
03153 COMPLEX BWIG
03154 REAL ROM,ROG,ROM1,ROG1,BETA1,PI,PIM,S,W
03155 EXTERNAL BWIG
03156 DATA INIT /0/
03157
03158
03159 IF (INIT.EQ.0 ) THEN
03160 INIT=1
03161 PI=3.141592654
03162 PIM=.140
03163 ROM=PKORB(1,9)
03164 ROG=PKORB(2,9)
03165 ROM1=PKORB(1,15)
03166 ROG1=PKORB(2,15)
03167 BETA1=PKORB(3,15)
03168 ENDIF
03169
03170 S=W**2
03171 FPIK= (BWIG(S,ROM,ROG)+BETA1*BWIG(S,ROM1,ROG1))
03172 & /(1+BETA1)
03173 RETURN
03174 END
03175 FUNCTION FPIRHO(W)
03176
03177
03178
03179 COMPLEX FPIK
03180 FPIRHO=CABS(FPIK(W))**2
03181 END
03182 SUBROUTINE CLVEC(HJ,PN,PIV)
03183
03184
03185
03186
03187
03188
03189 REAL PIV(4),PN(4)
03190 COMPLEX HJ(4),HN
03191
03192 HN= HJ(4)*CMPLX(PN(4))-HJ(3)*CMPLX(PN(3))
03193 HH= REAL(HJ(4)*CONJG(HJ(4))-HJ(3)*CONJG(HJ(3))
03194 $ -HJ(2)*CONJG(HJ(2))-HJ(1)*CONJG(HJ(1)))
03195 DO 10 I=1,4
03196 10 PIV(I)=4.*REAL(HN*CONJG(HJ(I)))-2.*HH*PN(I)
03197 RETURN
03198 END
03199 SUBROUTINE CLAXI(HJ,PN,PIA)
03200
03201
03202
03203
03204
03205
03206 COMMON / JAKI / JAK1,JAK2,JAKP,JAKM,KTOM
03207 COMMON / IDFC / IDFF
03208 REAL PIA(4),PN(4)
03209 COMPLEX HJ(4),HJC(4)
03210
03211
03212 DET2(I,J)=AIMAG(HJC(I)*HJ(J)-HJC(J)*HJ(I))
03213
03214
03215
03216
03217 IF (KTOM.EQ.1.OR.KTOM.EQ.-1) THEN
03218 SIGN= IDFF/ABS(IDFF)
03219 ELSEIF (KTOM.EQ.2) THEN
03220 SIGN=-IDFF/ABS(IDFF)
03221 ELSE
03222 PRINT *, 'STOP IN CLAXI: KTOM=',KTOM
03223 STOP
03224 ENDIF
03225
03226 DO 10 I=1,4
03227 10 HJC(I)=CONJG(HJ(I))
03228 PIA(1)= -2.*PN(3)*DET2(2,4)+2.*PN(4)*DET2(2,3)
03229 PIA(2)= -2.*PN(4)*DET2(1,3)+2.*PN(3)*DET2(1,4)
03230 PIA(3)= 2.*PN(4)*DET2(1,2)
03231 PIA(4)= 2.*PN(3)*DET2(1,2)
03232
03233 DO 20 I=1,4
03234 20 PIA(I)=PIA(I)*SIGN
03235 END
03236 SUBROUTINE CLNUT(HJ,B,HV)
03237
03238
03239
03240
03241
03242
03243 COMPLEX HJ(4)
03244 REAL HV(4),P(4)
03245 DATA P /3*0.,1.0/
03246
03247 CALL CLAXI(HJ,P,HV)
03248 B=REAL( HJ(4)*AIMAG(HJ(4)) - HJ(3)*AIMAG(HJ(3))
03249 & - HJ(2)*AIMAG(HJ(2)) - HJ(1)*AIMAG(HJ(1)) )
03250 RETURN
03251 END
03252 SUBROUTINE DAMPOG(PT,PN,PIM1,PIM2,PIPL,AMPLIT,HV)
03253
03254
03255
03256
03257
03258
03259
03260
03261
03262 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
03263 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
03264 * ,AMK,AMKZ,AMKST,GAMKST
03265
03266 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
03267 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
03268 * ,AMK,AMKZ,AMKST,GAMKST
03269 COMMON / DECPAR / GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
03270 REAL*4 GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
03271 COMMON /TESTA1/ KEYA1
03272 REAL HV(4),PT(4),PN(4),PIM1(4),PIM2(4),PIPL(4)
03273 REAL PAA(4),VEC1(4),VEC2(4)
03274 REAL PIVEC(4),PIAKS(4),HVM(4)
03275 COMPLEX BWIGN,HADCUR(4),FNORM,FORMOM
03276 DATA ICONT /1/
03277
03278
03279
03280
03281 DO 10 I=1,4
03282 VEC1(I)=0.0
03283 VEC2(I)=0.0
03284 HV(I) =0.0
03285 10 PAA(I)=PIM1(I)+PIM2(I)+PIPL(I)
03286 VEC1(1)=1.0
03287
03288 XMAA =SQRT(ABS(PAA(4)**2-PAA(3)**2-PAA(2)**2-PAA(1)**2))
03289 XMOM =SQRT(ABS( (PIM2(4)+PIPL(4))**2-(PIM2(3)+PIPL(3))**2
03290 $ -(PIM2(2)+PIPL(2))**2-(PIM2(1)+PIPL(1))**2 ))
03291 XMRO2 =(PIPL(1))**2 +(PIPL(2))**2 +(PIPL(3))**2
03292
03293 PROD1 =VEC1(1)*PIPL(1)
03294 PROD2 =VEC2(2)*PIPL(2)
03295 P12 =PIM1(4)*PIM2(4)-PIM1(1)*PIM2(1)
03296 $ -PIM1(2)*PIM2(2)-PIM1(3)*PIM2(3)
03297 P1PL =PIM1(4)*PIPL(4)-PIM1(1)*PIPL(1)
03298 $ -PIM1(2)*PIPL(2)-PIM1(3)*PIPL(3)
03299 P2PL =PIPL(4)*PIM2(4)-PIPL(1)*PIM2(1)
03300 $ -PIPL(2)*PIM2(2)-PIPL(3)*PIM2(3)
03301 DO 40 I=1,3
03302 VEC1(I)= (VEC1(I)-PROD1/XMRO2*PIPL(I))
03303 40 CONTINUE
03304 GNORM=SQRT(VEC1(1)**2+VEC1(2)**2+VEC1(3)**2)
03305 DO 41 I=1,3
03306 VEC1(I)= VEC1(I)/GNORM
03307 41 CONTINUE
03308 VEC2(1)=(VEC1(2)*PIPL(3)-VEC1(3)*PIPL(2))/SQRT(XMRO2)
03309 VEC2(2)=(VEC1(3)*PIPL(1)-VEC1(1)*PIPL(3))/SQRT(XMRO2)
03310 VEC2(3)=(VEC1(1)*PIPL(2)-VEC1(2)*PIPL(1))/SQRT(XMRO2)
03311 P1VEC1 =PIM1(4)*VEC1(4)-PIM1(1)*VEC1(1)
03312 $ -PIM1(2)*VEC1(2)-PIM1(3)*VEC1(3)
03313 P2VEC1 =VEC1(4)*PIM2(4)-VEC1(1)*PIM2(1)
03314 $ -VEC1(2)*PIM2(2)-VEC1(3)*PIM2(3)
03315 P1VEC2 =PIM1(4)*VEC2(4)-PIM1(1)*VEC2(1)
03316 $ -PIM1(2)*VEC2(2)-PIM1(3)*VEC2(3)
03317 P2VEC2 =VEC2(4)*PIM2(4)-VEC2(1)*PIM2(1)
03318 $ -VEC2(2)*PIM2(2)-VEC2(3)*PIM2(3)
03319
03320 FNORM=FORMOM(XMAA,XMOM)
03321 BRAK=0.0
03322 DO 120 JJ=1,2
03323 DO 45 I=1,4
03324 IF (JJ.EQ.1) THEN
03325 HADCUR(I) = FNORM *(
03326 $ VEC1(I)*(AMPI**2*P1PL-P2PL*(P12-P1PL))
03327 $ -PIM2(I)*(P2VEC1*P1PL-P1VEC1*P2PL)
03328 $ +PIPL(I)*(P2VEC1*P12 -P1VEC1*(AMPI**2+P2PL)) )
03329 ELSE
03330 HADCUR(I) = FNORM *(
03331 $ VEC2(I)*(AMPI**2*P1PL-P2PL*(P12-P1PL))
03332 $ -PIM2(I)*(P2VEC2*P1PL-P1VEC2*P2PL)
03333 $ +PIPL(I)*(P2VEC2*P12 -P1VEC2*(AMPI**2+P2PL)) )
03334 ENDIF
03335 45 CONTINUE
03336
03337
03338 CALL CLVEC(HADCUR,PN,PIVEC)
03339 CALL CLAXI(HADCUR,PN,PIAKS)
03340 CALL CLNUT(HADCUR,BRAKM,HVM)
03341
03342 BRAK=BRAK+(GV**2+GA**2)*PT(4)*PIVEC(4) +2.*GV*GA*PT(4)*PIAKS(4)
03343 & +2.*(GV**2-GA**2)*AMNUTA*AMTAU*BRAKM
03344 DO 90 I=1,3
03345 HV(I)=HV(I)-(AMTAU*((GV**2+GA**2)*PIAKS(I)+2.*GV*GA*PIVEC(I)))
03346 & +(GV**2-GA**2)*AMNUTA*AMTAU*HVM(I)
03347 90 CONTINUE
03348
03349 120 CONTINUE
03350 AMPLIT=(GFERMI*CCABIB)**2*BRAK/2.
03351
03352
03353
03354 DO 91 I=1,3
03355 HV(I)=-HV(I)/BRAK
03356 91 CONTINUE
03357
03358 END
03359 SUBROUTINE DAMPPK(MNUM,PT,PN,PIM1,PIM2,PIM3,AMPLIT,HV)
03360
03361
03362
03363
03364
03365
03366
03367
03368
03369 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
03370 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
03371 * ,AMK,AMKZ,AMKST,GAMKST
03372
03373 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
03374 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
03375 * ,AMK,AMKZ,AMKST,GAMKST
03376 COMMON / DECPAR / GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
03377 REAL*4 GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
03378 REAL HV(4),PT(4),PN(4),PIM1(4),PIM2(4),PIM3(4)
03379 REAL PAA(4),VEC1(4),VEC2(4),VEC3(4),VEC4(4),VEC5(4)
03380 REAL PIVEC(4),PIAKS(4),HVM(4)
03381 REAL FNORM(0:7),COEF(1:5,0:7)
03382 COMPLEX HADCUR(4),FORM1,FORM2,FORM3,FORM4,FORM5,UROJ
03383 COMPLEX F1,F2,F3,F4,F5
03384 EXTERNAL FORM1,FORM2,FORM3,FORM4,FORM5
03385 DATA PI /3.141592653589793238462643/
03386 DATA ICONT /0/
03387
03388 DATA FPI /93.3E-3/
03389 IF (ICONT.EQ.0) THEN
03390 ICONT=1
03391 UROJ=CMPLX(0.0,1.0)
03392 DWAPI0=SQRT(2.0)
03393 FNORM(0)=CCABIB/FPI
03394 FNORM(1)=CCABIB/FPI
03395 FNORM(2)=CCABIB/FPI
03396 FNORM(3)=CCABIB/FPI
03397 FNORM(4)=SCABIB/FPI/DWAPI0
03398 FNORM(5)=SCABIB/FPI
03399 FNORM(6)=SCABIB/FPI
03400 FNORM(7)=CCABIB/FPI
03401
03402 COEF(1,0)= 2.0*SQRT(2.)/3.0
03403 COEF(2,0)=-2.0*SQRT(2.)/3.0
03404
03405 COEF(3,0)= 2.0*SQRT(2.)/3.0
03406 COEF(4,0)= FPI
03407 COEF(5,0)= 0.0
03408
03409 COEF(1,1)=-SQRT(2.)/3.0
03410 COEF(2,1)= SQRT(2.)/3.0
03411 COEF(3,1)= 0.0
03412 COEF(4,1)= FPI
03413 COEF(5,1)= SQRT(2.)
03414
03415 COEF(1,2)=-SQRT(2.)/3.0
03416 COEF(2,2)= SQRT(2.)/3.0
03417 COEF(3,2)= 0.0
03418 COEF(4,2)= 0.0
03419 COEF(5,2)=-SQRT(2.)
03420
03421
03422 COEF(1,3)= 1./3.
03423 COEF(2,3)=-2./3.
03424 COEF(3,3)= 2./3.
03425 COEF(4,3)= 0.0
03426 COEF(5,3)= 0.0
03427
03428 COEF(1,4)= 1.0/SQRT(2.)/3.0
03429 COEF(2,4)=-1.0/SQRT(2.)/3.0
03430 COEF(3,4)= 0.0
03431 COEF(4,4)= 0.0
03432 COEF(5,4)= 0.0
03433
03434 COEF(1,5)=-SQRT(2.)/3.0
03435 COEF(2,5)= SQRT(2.)/3.0
03436 COEF(3,5)= 0.0
03437 COEF(4,5)= 0.0
03438 COEF(5,5)=-SQRT(2.)
03439
03440
03441 COEF(1,6)= 1./3.
03442 COEF(2,6)=-2./3.
03443 COEF(3,6)= 2./3.
03444 COEF(4,6)= 0.0
03445 COEF(5,6)=-2.0
03446
03447 COEF(1,7)= 0.0
03448 COEF(2,7)= 0.0
03449 COEF(3,7)= 0.0
03450 COEF(4,7)= 0.0
03451 COEF(5,7)=-SQRT(2.0/3.0)
03452
03453 ENDIF
03454
03455 DO 10 I=1,4
03456 10 PAA(I)=PIM1(I)+PIM2(I)+PIM3(I)
03457 XMAA =SQRT(ABS(PAA(4)**2-PAA(3)**2-PAA(2)**2-PAA(1)**2))
03458 XMRO1 =SQRT(ABS((PIM3(4)+PIM2(4))**2-(PIM3(1)+PIM2(1))**2
03459 $ -(PIM3(2)+PIM2(2))**2-(PIM3(3)+PIM2(3))**2))
03460 XMRO2 =SQRT(ABS((PIM3(4)+PIM1(4))**2-(PIM3(1)+PIM1(1))**2
03461 $ -(PIM3(2)+PIM1(2))**2-(PIM3(3)+PIM1(3))**2))
03462 XMRO3 =SQRT(ABS((PIM1(4)+PIM2(4))**2-(PIM1(1)+PIM2(1))**2
03463 $ -(PIM1(2)+PIM2(2))**2-(PIM1(3)+PIM2(3))**2))
03464
03465 PROD1 =PAA(4)*(PIM2(4)-PIM3(4))-PAA(1)*(PIM2(1)-PIM3(1))
03466 $ -PAA(2)*(PIM2(2)-PIM3(2))-PAA(3)*(PIM2(3)-PIM3(3))
03467 PROD2 =PAA(4)*(PIM3(4)-PIM1(4))-PAA(1)*(PIM3(1)-PIM1(1))
03468 $ -PAA(2)*(PIM3(2)-PIM1(2))-PAA(3)*(PIM3(3)-PIM1(3))
03469 PROD3 =PAA(4)*(PIM1(4)-PIM2(4))-PAA(1)*(PIM1(1)-PIM2(1))
03470 $ -PAA(2)*(PIM1(2)-PIM2(2))-PAA(3)*(PIM1(3)-PIM2(3))
03471 DO 40 I=1,4
03472 VEC1(I)= PIM2(I)-PIM3(I) -PAA(I)*PROD1/XMAA**2
03473 VEC2(I)= PIM3(I)-PIM1(I) -PAA(I)*PROD2/XMAA**2
03474 VEC3(I)= PIM1(I)-PIM2(I) -PAA(I)*PROD3/XMAA**2
03475 40 VEC4(I)= PIM1(I)+PIM2(I)+PIM3(I)
03476 CALL PROD5(PIM1,PIM2,PIM3,VEC5)
03477
03478
03479
03480 F1 = CMPLX(COEF(1,MNUM))*FORM1(MNUM,XMAA**2,XMRO1**2,XMRO2**2)
03481 F2 = CMPLX(COEF(2,MNUM))*FORM2(MNUM,XMAA**2,XMRO2**2,XMRO1**2)
03482 F3 = CMPLX(COEF(3,MNUM))*FORM3(MNUM,XMAA**2,XMRO3**2,XMRO1**2)
03483 F4 = (-1.0*UROJ)*
03484 $CMPLX(COEF(4,MNUM))*FORM4(MNUM,XMAA**2,XMRO1**2,XMRO2**2,XMRO3**2)
03485 F5 = (-1.0)*UROJ/4.0/PI**2/FPI**2*
03486 $ CMPLX(COEF(5,MNUM))*FORM5(MNUM,XMAA**2,XMRO1**2,XMRO2**2)
03487
03488 DO 45 I=1,4
03489 HADCUR(I)= CMPLX(FNORM(MNUM)) * (
03490 $ CMPLX(VEC1(I))*F1+CMPLX(VEC2(I))*F2+CMPLX(VEC3(I))*F3+
03491 $ CMPLX(VEC4(I))*F4+CMPLX(VEC5(I))*F5)
03492 45 CONTINUE
03493
03494
03495 CALL CLVEC(HADCUR,PN,PIVEC)
03496 CALL CLAXI(HADCUR,PN,PIAKS)
03497 CALL CLNUT(HADCUR,BRAKM,HVM)
03498
03499 BRAK= (GV**2+GA**2)*PT(4)*PIVEC(4) +2.*GV*GA*PT(4)*PIAKS(4)
03500 & +2.*(GV**2-GA**2)*AMNUTA*AMTAU*BRAKM
03501 AMPLIT=(GFERMI)**2*BRAK/2.
03502 IF (MNUM.GE.9) THEN
03503 PRINT *, 'MNUM=',MNUM
03504 ZNAK=-1.0
03505 XM1=0.0
03506 XM2=0.0
03507 XM3=0.0
03508 DO 77 K=1,4
03509 IF (K.EQ.4) ZNAK=1.0
03510 XM1=ZNAK*PIM1(K)**2+XM1
03511 XM2=ZNAK*PIM2(K)**2+XM2
03512 XM3=ZNAK*PIM3(K)**2+XM3
03513 77 PRINT *, 'PIM1=',PIM1(K),'PIM2=',PIM2(K),'PIM3=',PIM3(K)
03514 PRINT *, 'XM1=',SQRT(XM1),'XM2=',SQRT(XM2),'XM3=',SQRT(XM3)
03515 PRINT *, '************************************************'
03516 ENDIF
03517
03518 DO 90 I=1,3
03519 HV(I)=-(AMTAU*((GV**2+GA**2)*PIAKS(I)+2.*GV*GA*PIVEC(I)))
03520 & +(GV**2-GA**2)*AMNUTA*AMTAU*HVM(I)
03521
03522 HV(I)=-HV(I)/BRAK
03523 90 CONTINUE
03524 END
03525 SUBROUTINE PROD5(P1,P2,P3,PIA)
03526
03527
03528
03529
03530
03531 COMMON / JAKI / JAK1,JAK2,JAKP,JAKM,KTOM
03532 COMMON / IDFC / IDFF
03533 REAL PIA(4),P1(4),P2(4),P3(4)
03534 DET2(I,J)=P1(I)*P2(J)-P2(I)*P1(J)
03535
03536 IF (KTOM.EQ.1.OR.KTOM.EQ.-1) THEN
03537 SIGN= IDFF/ABS(IDFF)
03538 ELSEIF (KTOM.EQ.2) THEN
03539 SIGN=-IDFF/ABS(IDFF)
03540 ELSE
03541 PRINT *, 'STOP IN PROD5: KTOM=',KTOM
03542 STOP
03543 ENDIF
03544
03545
03546
03547 PIA(1)= -P3(3)*DET2(2,4)+P3(4)*DET2(2,3)+P3(2)*DET2(3,4)
03548 PIA(2)= -P3(4)*DET2(1,3)+P3(3)*DET2(1,4)-P3(1)*DET2(3,4)
03549 PIA(3)= P3(4)*DET2(1,2)-P3(2)*DET2(1,4)+P3(1)*DET2(2,4)
03550 PIA(4)= P3(3)*DET2(1,2)-P3(2)*DET2(1,3)+P3(1)*DET2(2,3)
03551
03552 DO 20 I=1,4
03553 20 PIA(I)=PIA(I)*SIGN
03554 END
03555
03556 SUBROUTINE DEXNEW(MODE,ISGN,POL,PNU,PAA,PNPI,JNPI)
03557
03558
03559
03560
03561
03562
03563
03564
03565
03566
03567 COMMON / INOUT / INUT,IOUT
03568 REAL POL(4),HV(4),PAA(4),PNU(4),PNPI(4,9),RN(1)
03569 DATA IWARM/0/
03570
03571 IF(MODE.EQ.-1) THEN
03572
03573 IWARM=1
03574 CALL DADNEW( -1,ISGN,HV,PNU,PAA,PNPI,JDUMM)
03575
03576
03577 ELSEIF(MODE.EQ. 0) THEN
03578
03579 300 CONTINUE
03580 IF(IWARM.EQ.0) GOTO 902
03581 CALL DADNEW( 0,ISGN,HV,PNU,PAA,PNPI,JNPI)
03582 WT=(1+POL(1)*HV(1)+POL(2)*HV(2)+POL(3)*HV(3))/2.
03583
03584 CALL RANMAR(RN,1)
03585 IF(RN(1).GT.WT) GOTO 300
03586
03587 ELSEIF(MODE.EQ. 1) THEN
03588
03589 CALL DADNEW( 1,ISGN,HV,PNU,PAA,PNPI,JDUMM)
03590
03591 ENDIF
03592
03593 RETURN
03594 902 WRITE(IOUT, 9020)
03595 9020 FORMAT(' ----- DEXNEW: LACK OF INITIALISATION')
03596 STOP
03597 END
03598 SUBROUTINE DADNEW(MODE,ISGN,HV,PNU,PWB,PNPI,JNPI)
03599
03600 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
03601 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
03602 * ,AMK,AMKZ,AMKST,GAMKST
03603
03604 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
03605 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
03606 * ,AMK,AMKZ,AMKST,GAMKST
03607 COMMON / DECPAR / GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
03608 REAL*4 GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
03609 COMMON / TAUBMC / GAMPMC(30),GAMPER(30),NEVDEC(30)
03610 REAL*4 GAMPMC ,GAMPER
03611 COMMON / INOUT / INUT,IOUT
03612 PARAMETER (NMODE=15,NM1=0,NM2=1,NM3=8,NM4=2,NM5=1,NM6=3)
03613 COMMON / TAUDCD /IDFFIN(9,NMODE),MULPIK(NMODE)
03614 & ,NAMES
03615 CHARACTER NAMES(NMODE)*31
03616
03617 REAL*4 PNU(4),PWB(4),PNPI(4,9),HV(4),HHV(4)
03618 REAL*4 PDUM1(4),PDUM2(4),PDUMI(4,9)
03619 REAL*4 RRR(3)
03620 REAL*4 WTMAX(NMODE)
03621 REAL*8 SWT(NMODE),SSWT(NMODE)
03622 DIMENSION NEVRAW(NMODE),NEVOVR(NMODE),NEVACC(NMODE)
03623
03624 DATA PI /3.141592653589793238462643/
03625 DATA IWARM/0/
03626
03627 IF(MODE.EQ.-1) THEN
03628
03629
03630 NMOD=NMODE
03631 IWARM=1
03632
03633 DO 1 JNPI=1,NMOD
03634 NEVRAW(JNPI)=0
03635 NEVACC(JNPI)=0
03636 NEVOVR(JNPI)=0
03637 SWT(JNPI)=0
03638 SSWT(JNPI)=0
03639 WTMAX(JNPI)=-1.
03640
03641
03642 NTRIALS = 500
03643 IF (JNPI.LE.NM4) THEN
03644 A = PKORB(3,37+JNPI)
03645 IF (A.LT.0.) THEN
03646 NTRIALS = 10000
03647 ELSE
03648 NTRIALS = 0
03649 WTMAX(JNPI)=A
03650 END IF
03651 END IF
03652 DO I=1,NTRIALS
03653 IF (JNPI.LE.0) THEN
03654 GOTO 903
03655 ELSEIF(JNPI.LE.NM4) THEN
03656 CALL DPH4PI(WT,HV,PDUM1,PDUM2,PDUMI,JNPI)
03657 ELSEIF(JNPI.LE.NM4+NM5) THEN
03658 CALL DPH5PI(WT,HV,PDUM1,PDUM2,PDUMI,JNPI)
03659 ELSEIF(JNPI.LE.NM4+NM5+NM6) THEN
03660 CALL DPHNPI(WT,HV,PDUM1,PDUM2,PDUMI,JNPI)
03661 ELSEIF(JNPI.LE.NM4+NM5+NM6+NM3) THEN
03662 INUM=JNPI-NM4-NM5-NM6
03663 CALL DPHSPK(WT,HV,PDUM1,PDUM2,PDUMI,INUM)
03664 ELSEIF(JNPI.LE.NM4+NM5+NM6+NM3+NM2) THEN
03665 INUM=JNPI-NM4-NM5-NM6-NM3
03666 CALL DPHSRK(WT,HV,PDUM1,PDUM2,PDUMI,INUM)
03667 ELSE
03668 GOTO 903
03669 ENDIF
03670 IF(WT.GT.WTMAX(JNPI)/1.2) WTMAX(JNPI)=WT*1.2
03671 ENDDO
03672
03673
03674
03675 1 CONTINUE
03676 WRITE(IOUT,7005)
03677
03678 ELSEIF(MODE.EQ. 0) THEN
03679
03680 IF(IWARM.EQ.0) GOTO 902
03681
03682 300 CONTINUE
03683 IF (JNPI.LE.0) THEN
03684 GOTO 903
03685 ELSEIF(JNPI.LE.NM4) THEN
03686 CALL DPH4PI(WT,HHV,PNU,PWB,PNPI,JNPI)
03687 ELSEIF(JNPI.LE.NM4+NM5) THEN
03688 CALL DPH5PI(WT,HHV,PNU,PWB,PNPI,JNPI)
03689 ELSEIF(JNPI.LE.NM4+NM5+NM6) THEN
03690 CALL DPHNPI(WT,HHV,PNU,PWB,PNPI,JNPI)
03691 ELSEIF(JNPI.LE.NM4+NM5+NM6+NM3) THEN
03692 INUM=JNPI-NM4-NM5-NM6
03693 CALL DPHSPK(WT,HHV,PNU,PWB,PNPI,INUM)
03694 ELSEIF(JNPI.LE.NM4+NM5+NM6+NM3+NM2) THEN
03695 INUM=JNPI-NM4-NM5-NM6-NM3
03696 CALL DPHSRK(WT,HHV,PNU,PWB,PNPI,INUM)
03697 ELSE
03698 GOTO 903
03699 ENDIF
03700 DO I=1,4
03701 HV(I)=-ISGN*HHV(I)
03702 ENDDO
03703
03704 NEVRAW(JNPI)=NEVRAW(JNPI)+1
03705 SWT(JNPI)=SWT(JNPI)+WT
03706
03707
03708 SSWT(JNPI)=SSWT(JNPI)+dble(WT)**2
03709
03710 CALL RANMAR(RRR,3)
03711 RN=RRR(1)
03712 IF(WT.GT.WTMAX(JNPI)) NEVOVR(JNPI)=NEVOVR(JNPI)+1
03713 IF(RN*WTMAX(JNPI).GT.WT) GOTO 300
03714
03715 COSTHE=-1.+2.*RRR(2)
03716 THET=ACOS(COSTHE)
03717 PHI =2*PI*RRR(3)
03718 CALL ROTOR2(THET,PNU,PNU)
03719 CALL ROTOR3( PHI,PNU,PNU)
03720 CALL ROTOR2(THET,PWB,PWB)
03721 CALL ROTOR3( PHI,PWB,PWB)
03722 CALL ROTOR2(THET,HV,HV)
03723 CALL ROTOR3( PHI,HV,HV)
03724 ND=MULPIK(JNPI)
03725 DO 301 I=1,ND
03726 CALL ROTOR2(THET,PNPI(1,I),PNPI(1,I))
03727 CALL ROTOR3( PHI,PNPI(1,I),PNPI(1,I))
03728 301 CONTINUE
03729 NEVACC(JNPI)=NEVACC(JNPI)+1
03730
03731 ELSEIF(MODE.EQ. 1) THEN
03732
03733 DO 500 JNPI=1,NMOD
03734 IF(NEVRAW(JNPI).EQ.0) GOTO 500
03735 PARGAM=SWT(JNPI)/FLOAT(NEVRAW(JNPI)+1)
03736 ERROR=0
03737 IF(NEVRAW(JNPI).NE.0)
03738 & ERROR=SQRT(SSWT(JNPI)/SWT(JNPI)**2-1./FLOAT(NEVRAW(JNPI)))
03739 RAT=PARGAM/GAMEL
03740 WRITE(IOUT, 7010) NAMES(JNPI),
03741 & NEVRAW(JNPI),NEVACC(JNPI),NEVOVR(JNPI),PARGAM,RAT,ERROR
03742
03743 GAMPMC(8+JNPI-1)=RAT
03744 GAMPER(8+JNPI-1)=ERROR
03745
03746 500 CONTINUE
03747 ENDIF
03748
03749 RETURN
03750 7003 FORMAT(///1X,15(5H*****)
03751 $ /,' *', 25X,'******** DADNEW INITIALISATION ********',9X,1H*
03752 $ )
03753 7004 FORMAT(' *',E20.5,5X,'WTMAX = MAXIMUM WEIGHT ',9X,1H*/)
03754 7005 FORMAT(
03755 $ /,1X,15(5H*****)/)
03756 7010 FORMAT(///1X,15(5H*****)
03757 $ /,' *', 25X,'******** DADNEW FINAL REPORT ******** ',9X,1H*
03758 $ /,' *', 25X,'CHANNEL:',A31 ,9X,1H*
03759 $ /,' *',I20 ,5X,'NEVRAW = NO. OF DECAYS TOTAL ',9X,1H*
03760 $ /,' *',I20 ,5X,'NEVACC = NO. OF DECAYS ACCEPTED ',9X,1H*
03761 $ /,' *',I20 ,5X,'NEVOVR = NO. OF OVERWEIGHTED EVENTS ',9X,1H*
03762 $ /,' *',E20.5,5X,'PARTIAL WTDTH IN GEV UNITS ',9X,1H*
03763 $ /,' *',F20.9,5X,'IN UNITS GFERMI**2*MASS**5/192/PI**3 ',9X,1H*
03764 $ /,' *',F20.8,5X,'RELATIVE ERROR OF PARTIAL WIDTH ',9X,1H*
03765 $ /,1X,15(5H*****)/)
03766 902 WRITE(IOUT, 9020)
03767 9020 FORMAT(' ----- DADNEW: LACK OF INITIALISATION')
03768 STOP
03769 903 WRITE(IOUT, 9030) JNPI,MODE
03770 9030 FORMAT(' ----- DADNEW: WRONG JNPI',2I5)
03771 STOP
03772 END
03773
03774
03775 SUBROUTINE DPH4PI(DGAMT,HV,PN,PAA,PMULT,JNPI)
03776
03777
03778
03779
03780 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
03781 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
03782 * ,AMK,AMKZ,AMKST,GAMKST
03783
03784 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
03785 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
03786 * ,AMK,AMKZ,AMKST,GAMKST
03787 COMMON / DECPAR / GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
03788 REAL*4 GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
03789 REAL HV(4),PT(4),PN(4),PAA(4),PIM1(4),PIM2(4),PIPL(4),PMULT(4,9)
03790 REAL PR(4),PIZ(4)
03791 REAL*4 RRR(9)
03792 REAL*8 UU,FF,FF1,FF2,FF3,FF4,GG1,GG2,GG3,GG4,RR
03793 DATA PI /3.141592653589793238462643/
03794 DATA ICONT /0/
03795 XLAM(X,Y,Z)=SQRT(ABS((X-Y-Z)**2-4.0*Y*Z))
03796
03797
03798
03799
03800 PHSPAC=1./2**23/PI**11
03801 PHSP=1./2**5/PI**2
03802 IF (JNPI.EQ.1) THEN
03803 PREZ=0.7
03804 AMP1=AMPI
03805 AMP2=AMPI
03806 AMP3=AMPI
03807 AMP4=AMPIZ
03808 AMRX=PKORB(1,14)
03809 GAMRX=PKORB(2,14)
03810
03811
03812
03813
03814 AMROP =1.2
03815 GAMROP=.46
03816 ELSE
03817 PREZ=0.0
03818 AMP1=AMPIZ
03819 AMP2=AMPIZ
03820 AMP3=AMPIZ
03821 AMP4=AMPI
03822 AMRX=1.4
03823 GAMRX=.6
03824 AMROP =AMRX
03825 GAMROP=GAMRX
03826
03827 ENDIF
03828 RRB=0.3
03829 CALL CHOICE(100+JNPI,RRB,ICHAN,PROB1,PROB2,PROB3,
03830 $ AMROP,GAMROP,AMRX,GAMRX,AMRB,GAMRB)
03831 PREZ=PROB1+PROB2
03832
03833 PT(1)=0.
03834 PT(2)=0.
03835 PT(3)=0.
03836 PT(4)=AMTAU
03837
03838 CALL RANMAR(RRR,9)
03839
03840
03841
03842
03843 RR1=RRR(6)
03844 AMS1=(AMP1+AMP2+AMP3+AMP4)**2
03845 AMS2=(AMTAU-AMNUTA)**2
03846 ALP1=ATAN((AMS1-AMROP**2)/AMROP/GAMROP)
03847 ALP2=ATAN((AMS2-AMROP**2)/AMROP/GAMROP)
03848 ALP=ALP1+RR1*(ALP2-ALP1)
03849 AM4SQ =AMROP**2+AMROP*GAMROP*TAN(ALP)
03850 AM4 =SQRT(AM4SQ)
03851 PHSPAC=PHSPAC*
03852 $ ((AM4SQ-AMROP**2)**2+(AMROP*GAMROP)**2)/(AMROP*GAMROP)
03853 PHSPAC=PHSPAC*(ALP2-ALP1)
03854
03855
03856 RR1=RRR(1)
03857 AMS1=(AMP2+AMP3+AMP4)**2
03858 AMS2=(AM4-AMP1)**2
03859 IF (RRR(9).GT.PREZ) THEN
03860 AM3SQ=AMS1+ RR1*(AMS2-AMS1)
03861 AM3 =SQRT(AM3SQ)
03862
03863 FF1=AMS2-AMS1
03864 ELSE
03865
03866 ALP1=ATAN((AMS1-AMRX**2)/AMRX/GAMRX)
03867 ALP2=ATAN((AMS2-AMRX**2)/AMRX/GAMRX)
03868 ALP=ALP1+RR1*(ALP2-ALP1)
03869 AM3SQ =AMRX**2+AMRX*GAMRX*TAN(ALP)
03870 AM3 =SQRT(AM3SQ)
03871
03872 FF1=((AM3SQ-AMRX**2)**2+(AMRX*GAMRX)**2)/(AMRX*GAMRX)
03873 FF1=FF1*(ALP2-ALP1)
03874 ENDIF
03875
03876 RR2=RRR(2)
03877 AMS1=(AMP3+AMP4)**2
03878 AMS2=(AM3-AMP2)**2
03879
03880 AM2SQ=AMS1+ RR2*(AMS2-AMS1)
03881 AM2 =SQRT(AM2SQ)
03882
03883 FF2=(AMS2-AMS1)
03884
03885 ENQ1=(AM2SQ-AMP3**2+AMP4**2)/(2*AM2)
03886 ENQ2=(AM2SQ+AMP3**2-AMP4**2)/(2*AM2)
03887 PPI= ENQ1**2-AMP4**2
03888 PPPI=SQRT(ABS(ENQ1**2-AMP4**2))
03889 PHSPAC=PHSPAC*(4*PI)*(2*PPPI/AM2)
03890
03891 CALL SPHERA(PPPI,PIZ)
03892 PIZ(4)=ENQ1
03893
03894 DO 30 I=1,3
03895 30 PIPL(I)=-PIZ(I)
03896 PIPL(4)=ENQ2
03897
03898
03899 PR(1)=0
03900 PR(2)=0
03901 PR(4)=1./(2*AM3)*(AM3**2+AM2**2-AMP2**2)
03902 PR(3)= SQRT(ABS(PR(4)**2-AM2**2))
03903 PPI = PR(4)**2-AM2**2
03904
03905 PIM1(1)=0
03906 PIM1(2)=0
03907 PIM1(4)=1./(2*AM3)*(AM3**2-AM2**2+AMP2**2)
03908 PIM1(3)=-PR(3)
03909
03910 FF3=(4*PI)*(2*PR(3)/AM3)
03911
03912 EXE=(PR(4)+PR(3))/AM2
03913 CALL BOSTR3(EXE,PIZ,PIZ)
03914 CALL BOSTR3(EXE,PIPL,PIPL)
03915 RR3=RRR(3)
03916 RR4=RRR(4)
03917 THET =ACOS(-1.+2*RR3)
03918 PHI = 2*PI*RR4
03919 CALL ROTPOL(THET,PHI,PIPL)
03920 CALL ROTPOL(THET,PHI,PIM1)
03921 CALL ROTPOL(THET,PHI,PIZ)
03922 CALL ROTPOL(THET,PHI,PR)
03923
03924
03925 PR(1)=0
03926 PR(2)=0
03927 PR(4)=1./(2*AM4)*(AM4**2+AM3**2-AMP1**2)
03928 PR(3)= SQRT(ABS(PR(4)**2-AM3**2))
03929 PPI = PR(4)**2-AM3**2
03930
03931 PIM2(1)=0
03932 PIM2(2)=0
03933 PIM2(4)=1./(2*AM4)*(AM4**2-AM3**2+AMP1**2)
03934 PIM2(3)=-PR(3)
03935
03936 FF4=(4*PI)*(2*PR(3)/AM4)
03937
03938 EXE=(PR(4)+PR(3))/AM3
03939 CALL BOSTR3(EXE,PIZ,PIZ)
03940 CALL BOSTR3(EXE,PIPL,PIPL)
03941 CALL BOSTR3(EXE,PIM1,PIM1)
03942 RR3=RRR(7)
03943 RR4=RRR(8)
03944 THET =ACOS(-1.+2*RR3)
03945 PHI = 2*PI*RR4
03946 CALL ROTPOL(THET,PHI,PIPL)
03947 CALL ROTPOL(THET,PHI,PIM1)
03948 CALL ROTPOL(THET,PHI,PIM2)
03949 CALL ROTPOL(THET,PHI,PIZ)
03950 CALL ROTPOL(THET,PHI,PR)
03951
03952
03953
03954 PAA(1)=0
03955 PAA(2)=0
03956 PAA(4)=1./(2*AMTAU)*(AMTAU**2-AMNUTA**2+AM4**2)
03957 PAA(3)= SQRT(ABS(PAA(4)**2-AM4**2))
03958 PPI = PAA(4)**2-AM4**2
03959 PHSPAC=PHSPAC*(4*PI)*(2*PAA(3)/AMTAU)
03960 PHSP=PHSP*(4*PI)*(2*PAA(3)/AMTAU)
03961
03962 PN(1)=0
03963 PN(2)=0
03964 PN(4)=1./(2*AMTAU)*(AMTAU**2+AMNUTA**2-AM4**2)
03965 PN(3)=-PAA(3)
03966
03967 IF(RRR(9).LE.0.5*PREZ) THEN
03968 DO 72 I=1,4
03969 X=PIM1(I)
03970 PIM1(I)=PIM2(I)
03971 72 PIM2(I)=X
03972 ENDIF
03973
03974
03975
03976 AM3SQ=(PIM1(4)+PIZ(4)+PIPL(4))**2-(PIM1(3)+PIZ(3)+PIPL(3))**2
03977 $ -(PIM1(2)+PIZ(2)+PIPL(2))**2-(PIM1(1)+PIZ(1)+PIPL(1))**2
03978 AMS2=(AM4-AMP2)**2
03979 AMS1=(AMP1+AMP3+AMP4)**2
03980 FF1=(AMS2-AMS1)
03981 AMS1=(AMP3+AMP4)**2
03982 AMS2=(SQRT(AM3SQ)-AMP1)**2
03983 FF2=AMS2-AMS1
03984 FF3=(4*PI)*(XLAM(AM2**2,AMP1**2,AM3SQ)/AM3SQ)
03985 FF4=(4*PI)*(XLAM(AM3SQ,AMP2**2,AM4**2)/AM4**2)
03986 UU=FF1*FF2*FF3*FF4
03987
03988 AM3SQ=(PIM1(4)+PIZ(4)+PIPL(4))**2-(PIM1(3)+PIZ(3)+PIPL(3))**2
03989 $ -(PIM1(2)+PIZ(2)+PIPL(2))**2-(PIM1(1)+PIZ(1)+PIPL(1))**2
03990 AMS2=(AM4-AMP2)**2
03991 AMS1=(AMP1+AMP3+AMP4)**2
03992 ALP1=ATAN((AMS1-AMRX**2)/AMRX/GAMRX)
03993 ALP2=ATAN((AMS2-AMRX**2)/AMRX/GAMRX)
03994 FF1=((AM3SQ-AMRX**2)**2+(AMRX*GAMRX)**2)/(AMRX*GAMRX)
03995 FF1=FF1*(ALP2-ALP1)
03996 AMS1=(AMP3+AMP4)**2
03997 AMS2=(SQRT(AM3SQ)-AMP1)**2
03998 FF2=AMS2-AMS1
03999 FF3=(4*PI)*(XLAM(AM2**2,AMP1**2,AM3SQ)/AM3SQ)
04000 FF4=(4*PI)*(XLAM(AM3SQ,AMP2**2,AM4**2)/AM4**2)
04001 FF=FF1*FF2*FF3*FF4
04002
04003 AM3SQ=(PIM2(4)+PIZ(4)+PIPL(4))**2-(PIM2(3)+PIZ(3)+PIPL(3))**2
04004 $ -(PIM2(2)+PIZ(2)+PIPL(2))**2-(PIM2(1)+PIZ(1)+PIPL(1))**2
04005 AMS2=(AM4-AMP1)**2
04006 AMS1=(AMP2+AMP3+AMP4)**2
04007 ALP1=ATAN((AMS1-AMRX**2)/AMRX/GAMRX)
04008 ALP2=ATAN((AMS2-AMRX**2)/AMRX/GAMRX)
04009 GG1=((AM3SQ-AMRX**2)**2+(AMRX*GAMRX)**2)/(AMRX*GAMRX)
04010 GG1=GG1*(ALP2-ALP1)
04011 AMS1=(AMP3+AMP4)**2
04012 AMS2=(SQRT(AM3SQ)-AMP2)**2
04013 GG2=AMS2-AMS1
04014 GG3=(4*PI)*(XLAM(AM2**2,AMP2**2,AM3SQ)/AM3SQ)
04015 GG4=(4*PI)*(XLAM(AM3SQ,AMP1**2,AM4**2)/AM4**2)
04016 GG=GG1*GG2*GG3*GG4
04017
04018 IF ( ( (FF+GG)*UU+FF*GG ).GT.0.0D0) THEN
04019 RR=FF*GG*UU/(0.5*PREZ*(FF+GG)*UU+(1.0-PREZ)*FF*GG)
04020 PHSPAC=PHSPAC*RR
04021 ELSE
04022 PHSPAC=0.0
04023 ENDIF
04024
04025 IF (JNPI.EQ.1) THEN
04026 RR5= RRR(5)
04027 IF(RR5.LE.0.5) THEN
04028 DO 70 I=1,4
04029 X=PIM1(I)
04030 PIM1(I)=PIM2(I)
04031 70 PIM2(I)=X
04032 ENDIF
04033 PHSPAC=PHSPAC/2.
04034 ELSE
04035
04036 RR5= RRR(5)
04037 IF(RR5.LE.0.5) THEN
04038 DO 71 I=1,4
04039 X=PIM1(I)
04040 PIM1(I)=PIM2(I)
04041 71 PIM2(I)=X
04042 ENDIF
04043 PHSPAC=PHSPAC/6.
04044 ENDIF
04045
04046
04047 EXE=(PAA(4)+PAA(3))/AM4
04048 CALL BOSTR3(EXE,PIZ,PIZ)
04049 CALL BOSTR3(EXE,PIPL,PIPL)
04050 CALL BOSTR3(EXE,PIM1,PIM1)
04051 CALL BOSTR3(EXE,PIM2,PIM2)
04052 CALL BOSTR3(EXE,PR,PR)
04053
04054
04055
04056
04057
04058
04059
04060 IF (JNPI.EQ.1) THEN
04061 CALL DAM4PI(JNPI,PT,PN,PIM1,PIM2,PIZ,PIPL,AMPLIT,HV)
04062 ELSEIF (JNPI.EQ.2) THEN
04063 CALL DAM4PI(JNPI,PT,PN,PIM1,PIM2,PIPL,PIZ,AMPLIT,HV)
04064 ENDIF
04065 DGAMT=1/(2.*AMTAU)*AMPLIT*PHSPAC
04066
04067
04068 DO 77 K=1,4
04069 PMULT(K,1)=PIM1(K)
04070 PMULT(K,2)=PIM2(K)
04071 PMULT(K,3)=PIPL(K)
04072 PMULT(K,4)=PIZ (K)
04073 77 CONTINUE
04074 END
04075 SUBROUTINE DAM4PI(MNUM,PT,PN,PIM1,PIM2,PIM3,PIM4,AMPLIT,HV)
04076
04077
04078
04079
04080
04081
04082
04083
04084
04085 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
04086 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
04087 * ,AMK,AMKZ,AMKST,GAMKST
04088
04089 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
04090 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
04091 * ,AMK,AMKZ,AMKST,GAMKST
04092 COMMON / DECPAR / GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
04093 REAL*4 GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
04094 REAL HV(4),PT(4),PN(4),PIM1(4),PIM2(4),PIM3(4),PIM4(4)
04095 REAL PIVEC(4),PIAKS(4),HVM(4)
04096 COMPLEX HADCUR(4),FORM1,FORM2,FORM3,FORM4,FORM5
04097 EXTERNAL FORM1,FORM2,FORM3,FORM4,FORM5
04098 DATA PI /3.141592653589793238462643/
04099 DATA ICONT /0/
04100
04101 CALL CURR_CLEO(MNUM,PIM1,PIM2,PIM3,PIM4,HADCUR)
04102
04103
04104 CALL CLVEC(HADCUR,PN,PIVEC)
04105 CALL CLAXI(HADCUR,PN,PIAKS)
04106 CALL CLNUT(HADCUR,BRAKM,HVM)
04107
04108 BRAK= (GV**2+GA**2)*PT(4)*PIVEC(4) +2.*GV*GA*PT(4)*PIAKS(4)
04109 & +2.*(GV**2-GA**2)*AMNUTA*AMTAU*BRAKM
04110 AMPLIT=(CCABIB*GFERMI)**2*BRAK/2.
04111
04112 DO 90 I=1,3
04113 HV(I)=-(AMTAU*((GV**2+GA**2)*PIAKS(I)+2.*GV*GA*PIVEC(I)))
04114 & +(GV**2-GA**2)*AMNUTA*AMTAU*HVM(I)
04115
04116 IF (BRAK.NE.0.0)
04117 &HV(I)=-HV(I)/BRAK
04118 90 CONTINUE
04119 END
04120 SUBROUTINE DPH5PI(DGAMT,HV,PN,PAA,PMULT,JNPI)
04121
04122
04123
04124
04125 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
04126 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
04127 * ,AMK,AMKZ,AMKST,GAMKST
04128
04129 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
04130 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
04131
04132
04133 * ,AMK,AMKZ,AMKST,GAMKST
04134 COMMON / DECPAR / GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
04135 REAL*4 GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
04136 PARAMETER (NMODE=15,NM1=0,NM2=1,NM3=8,NM4=2,NM5=1,NM6=3)
04137 COMMON / TAUDCD /IDFFIN(9,NMODE),MULPIK(NMODE)
04138 & ,NAMES
04139 CHARACTER NAMES(NMODE)*31
04140 REAL HV(4),PT(4),PN(4),PAA(4),PMULT(4,9)
04141 REAL*4 PR(4),PI1(4),PI2(4),PI3(4),PI4(4),PI5(4)
04142 REAL*8 AMP1,AMP2,AMP3,AMP4,AMP5,ams1,ams2,amom,gamom
04143 REAL*8 AM5SQ,AM4SQ,AM3SQ,AM2SQ,AM5,AM4,AM3
04144 REAL*4 RRR(10)
04145 REAL*8 gg1,gg2,gg3,ff1,ff2,ff3,ff4,alp,alp1,alp2
04146 REAL*8 XM,AM,GAMMA
04147
04148 real*8 phspac
04149
04150 DATA PI /3.141592653589793238462643/
04151 DATA ICONT /0/
04152 data fpi /93.3e-3/
04153
04154 COMPLEX BWIGN
04155
04156 BWIGN(XM,AM,GAMMA)=XM**2/CMPLX(XM**2-AM**2,GAMMA*AM)
04157
04158
04159 AMOM=.782
04160 GAMOM=0.0085
04161
04162
04163
04164 PHSPAC=1./2**29/PI**14
04165
04166
04167 AMP1=DCDMAS(IDFFIN(1,JNPI))
04168 AMP2=DCDMAS(IDFFIN(2,JNPI))
04169 AMP3=DCDMAS(IDFFIN(3,JNPI))
04170 AMP4=DCDMAS(IDFFIN(4,JNPI))
04171 AMP5=DCDMAS(IDFFIN(5,JNPI))
04172
04173
04174 PT(1)=0.
04175 PT(2)=0.
04176 PT(3)=0.
04177 PT(4)=AMTAU
04178
04179 CALL RANMAR(RRR,10)
04180
04181
04182
04183
04184
04185 rr1=rrr(10)
04186 ams1=(amp1+amp2+amp3+amp4+amp5)**2
04187 ams2=(amtau-amnuta)**2
04188 am5sq=ams1+ rr1*(ams2-ams1)
04189 am5 =sqrt(am5sq)
04190 phspac=phspac*(ams2-ams1)
04191
04192
04193
04194 rr1=rrr(9)
04195 ams1=(amp2+amp3+amp4+amp5)**2
04196 ams2=(am5-amp1)**2
04197 am4sq=ams1+ rr1*(ams2-ams1)
04198 am4 =sqrt(am4sq)
04199 gg1=ams2-ams1
04200
04201
04202
04203 rr1=rrr(1)
04204 ams1=(amp2+amp3+amp4)**2
04205 ams2=(am4-amp5)**2
04206 alp1=atan((ams1-amom**2)/amom/gamom)
04207 alp2=atan((ams2-amom**2)/amom/gamom)
04208 alp=alp1+rr1*(alp2-alp1)
04209 am3sq =amom**2+amom*gamom*tan(alp)
04210 am3 =sqrt(am3sq)
04211
04212 gg2=((am3sq-amom**2)**2+(amom*gamom)**2)/(amom*gamom)
04213 gg2=gg2*(alp2-alp1)
04214
04215
04216
04217
04218
04219
04220
04221 rr2=rrr(2)
04222 ams1=(amp3+amp4)**2
04223 ams2=(am3-amp2)**2
04224
04225 am2sq=ams1+ rr2*(ams2-ams1)
04226 am2 =sqrt(am2sq)
04227
04228 gg3=ams2-ams1
04229
04230
04231 enq1=(am2sq+amp3**2-amp4**2)/(2*am2)
04232 enq2=(am2sq-amp3**2+amp4**2)/(2*am2)
04233 ppi= enq1**2-amp3**2
04234 pppi=sqrt(abs(enq1**2-amp3**2))
04235 ff1=(4*pi)*(2*pppi/am2)
04236
04237 call sphera(pppi,pi3)
04238 pi3(4)=enq1
04239
04240 do 30 i=1,3
04241 30 pi4(i)=-pi3(i)
04242 pi4(4)=enq2
04243
04244
04245
04246 pr(1)=0
04247 pr(2)=0
04248 pr(4)=1./(2*am3)*(am3**2+am2**2-amp2**2)
04249 pr(3)= sqrt(abs(pr(4)**2-am2**2))
04250 ppi = pr(4)**2-am2**2
04251
04252 pi2(1)=0
04253 pi2(2)=0
04254 pi2(4)=1./(2*am3)*(am3**2-am2**2+amp2**2)
04255 pi2(3)=-pr(3)
04256
04257 ff2=(4*pi)*(2*pr(3)/am3)
04258
04259 exe=(pr(4)+pr(3))/am2
04260 call bostr3(exe,pi3,pi3)
04261 call bostr3(exe,pi4,pi4)
04262 rr3=rrr(3)
04263 rr4=rrr(4)
04264 thet =acos(-1.+2*rr3)
04265 phi = 2*pi*rr4
04266 call rotpol(thet,phi,pi2)
04267 call rotpol(thet,phi,pi3)
04268 call rotpol(thet,phi,pi4)
04269
04270
04271
04272 pr(1)=0
04273 pr(2)=0
04274 pr(4)=1./(2*am4)*(am4**2+am3**2-amp5**2)
04275 pr(3)= sqrt(abs(pr(4)**2-am3**2))
04276 ppi = pr(4)**2-am3**2
04277
04278 pi5(1)=0
04279 pi5(2)=0
04280 pi5(4)=1./(2*am4)*(am4**2-am3**2+amp5**2)
04281 pi5(3)=-pr(3)
04282
04283 ff3=(4*pi)*(2*pr(3)/am4)
04284
04285 exe=(pr(4)+pr(3))/am3
04286 call bostr3(exe,pi2,pi2)
04287 call bostr3(exe,pi3,pi3)
04288 call bostr3(exe,pi4,pi4)
04289 rr3=rrr(5)
04290 rr4=rrr(6)
04291 thet =acos(-1.+2*rr3)
04292 phi = 2*pi*rr4
04293 call rotpol(thet,phi,pi2)
04294 call rotpol(thet,phi,pi3)
04295 call rotpol(thet,phi,pi4)
04296 call rotpol(thet,phi,pi5)
04297
04298
04299
04300 pr(1)=0
04301 pr(2)=0
04302 pr(4)=1./(2*am5)*(am5**2+am4**2-amp1**2)
04303 pr(3)= sqrt(abs(pr(4)**2-am4**2))
04304 ppi = pr(4)**2-am4**2
04305
04306 pi1(1)=0
04307 pi1(2)=0
04308 pi1(4)=1./(2*am5)*(am5**2-am4**2+amp1**2)
04309 pi1(3)=-pr(3)
04310
04311 ff4=(4*pi)*(2*pr(3)/am5)
04312
04313 exe=(pr(4)+pr(3))/am4
04314 call bostr3(exe,pi2,pi2)
04315 call bostr3(exe,pi3,pi3)
04316 call bostr3(exe,pi4,pi4)
04317 call bostr3(exe,pi5,pi5)
04318 rr3=rrr(7)
04319 rr4=rrr(8)
04320 thet =acos(-1.+2*rr3)
04321 phi = 2*pi*rr4
04322 call rotpol(thet,phi,pi1)
04323 call rotpol(thet,phi,pi2)
04324 call rotpol(thet,phi,pi3)
04325 call rotpol(thet,phi,pi4)
04326 call rotpol(thet,phi,pi5)
04327
04328
04329
04330 paa(1)=0
04331 paa(2)=0
04332
04333
04334
04335 paa(4)=1./(2*amtau)*(amtau**2-amnuta**2+am5sq)
04336 paa(3)= sqrt(abs(paa(4)**2-am5sq))
04337 ppi = paa(4)**2-am5sq
04338 phspac=phspac*(4*pi)*(2*paa(3)/amtau)
04339
04340 pn(1)=0
04341 pn(2)=0
04342 pn(4)=1./(2*amtau)*(amtau**2+amnuta**2-am5**2)
04343 pn(3)=-paa(3)
04344
04345 phspac=phspac * gg1*gg2*gg3*ff1*ff2*ff3*ff4
04346
04347
04348
04349 exe=(paa(4)+paa(3))/am5
04350 call bostr3(exe,pi1,pi1)
04351 call bostr3(exe,pi2,pi2)
04352 call bostr3(exe,pi3,pi3)
04353 call bostr3(exe,pi4,pi4)
04354 call bostr3(exe,pi5,pi5)
04355
04356
04357
04358
04359
04360 PXQ=AMTAU*PAA(4)
04361 PXN=AMTAU*PN(4)
04362 QXN=PAA(4)*PN(4)-PAA(1)*PN(1)-PAA(2)*PN(2)-PAA(3)*PN(3)
04363 BRAK=2*(GV**2+GA**2)*(2*PXQ*QXN+AM5SQ*PXN)
04364 & -6*(GV**2-GA**2)*AMTAU*AMNUTA*AM5SQ
04365 fompp = cabs(bwign(am3,amom,gamom))**2
04366
04367
04368 fnorm = 1/fpi**6
04369
04370 AMPLIT=CCABIB**2*GFERMI**2/2. * BRAK
04371 amplit = amplit * fompp * fnorm
04372
04373
04374 DGAMT=1/(2.*AMTAU)*AMPLIT*PHSPAC
04375
04376 DO 40 I=1,3
04377 40 HV(I)=0.
04378
04379 do 77 k=1,4
04380 pmult(k,1)=pi1(k)
04381 pmult(k,2)=pi2(k)
04382 pmult(k,3)=pi3(k)
04383 pmult(k,4)=pi4(k)
04384 pmult(k,5)=pi5(k)
04385 77 continue
04386 return
04387
04388
04389
04390
04391 end
04392 SUBROUTINE DPHSRK(DGAMT,HV,PN,PR,PMULT,INUM)
04393
04394
04395
04396
04397
04398 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
04399 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
04400 * ,AMK,AMKZ,AMKST,GAMKST
04401
04402 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
04403 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
04404 * ,AMK,AMKZ,AMKST,GAMKST
04405 COMMON / DECPAR / GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
04406 REAL*4 GFERMI,GV,GA,CCABIB,SCABIB,GAMEL
04407 REAL HV(4),PT(4),PN(4),PR(4),PKC(4),PKZ(4),QQ(4),PMULT(4,9)
04408 REAL RR1(1)
04409 DATA PI /3.141592653589793238462643/
04410 DATA ICONT /0/
04411
04412
04413 PHSPAC=1./2**11/PI**5
04414
04415 PT(1)=0.
04416 PT(2)=0.
04417 PT(3)=0.
04418 PT(4)=AMTAU
04419
04420 AMS1=(AMK+AMKZ)**2
04421 AMS2=(AMTAU-AMNUTA)**2
04422
04423 CALL RANMAR(RR1,1)
04424 AMX2=AMS1+ RR1(1)*(AMS2-AMS1)
04425 AMX=SQRT(AMX2)
04426 PHSPAC=PHSPAC*(AMS2-AMS1)
04427
04428
04429
04430
04431 100 CONTINUE
04432
04433
04434
04435
04436
04437
04438
04439
04440
04441
04442 PN(1)=0
04443 PN(2)=0
04444 PN(4)=1./(2*AMTAU)*(AMTAU**2+AMNUTA**2-AMX**2)
04445 PN(3)=-SQRT((PN(4)-AMNUTA)*(PN(4)+AMNUTA))
04446
04447 PR(1)=0
04448 PR(2)=0
04449 PR(4)=1./(2*AMTAU)*(AMTAU**2-AMNUTA**2+AMX**2)
04450 PR(3)=-PN(3)
04451 PHSPAC=PHSPAC*(4*PI)*(2*PR(3)/AMTAU)
04452
04453
04454 ENQ1=(AMX2+AMK**2-AMKZ**2)/(2.*AMX)
04455 ENQ2=(AMX2-AMK**2+AMKZ**2)/(2.*AMX)
04456 PPPI=SQRT((ENQ1-AMK)*(ENQ1+AMK))
04457 PHSPAC=PHSPAC*(4*PI)*(2*PPPI/AMX)
04458
04459 CALL SPHERA(PPPI,PKC)
04460 PKC(4)=ENQ1
04461
04462 DO 20 I=1,3
04463 20 PKZ(I)=-PKC(I)
04464 PKZ(4)=ENQ2
04465 EXE=(PR(4)+PR(3))/AMX
04466
04467 CALL BOSTR3(EXE,PKC,PKC)
04468 CALL BOSTR3(EXE,PKZ,PKZ)
04469 DO 30 I=1,4
04470 30 QQ(I)=PKC(I)-PKZ(I)
04471
04472 PKSD =PR(4)*PR(4)-PR(3)*PR(3)-PR(2)*PR(2)-PR(1)*PR(1)
04473 QQPKS=PR(4)* QQ(4)-PR(3)* QQ(3)-PR(2)* QQ(2)-PR(1)* QQ(1)
04474 DO 31 I=1,4
04475 31 QQ(I)=QQ(I)-PR(I)*QQPKS/PKSD
04476
04477 PRODPQ=PT(4)*QQ(4)
04478 PRODNQ=PN(4)*QQ(4)-PN(1)*QQ(1)-PN(2)*QQ(2)-PN(3)*QQ(3)
04479 PRODPN=PT(4)*PN(4)
04480 QQ2= QQ(4)**2-QQ(1)**2-QQ(2)**2-QQ(3)**2
04481 BRAK=(GV**2+GA**2)*(2*PRODPQ*PRODNQ-PRODPN*QQ2)
04482 & +(GV**2-GA**2)*AMTAU*AMNUTA*QQ2
04483 AMPLIT=(GFERMI*CCABIB)**2*BRAK*2*FPIRK(AMX)
04484 DGAMT=1/(2.*AMTAU)*AMPLIT*PHSPAC
04485 DO 40 I=1,3
04486 40 HV(I)=2*GV*GA*AMTAU*(2*PRODNQ*QQ(I)-QQ2*PN(I))/BRAK
04487 do 77 k=1,4
04488 pmult(k,1)=pkc(k)
04489 pmult(k,2)=pkz(k)
04490 77 continue
04491 RETURN
04492 END
04493 FUNCTION FPIRK(W)
04494
04495
04496
04497 COMMON / PARMAS / AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
04498 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
04499 * ,AMK,AMKZ,AMKST,GAMKST
04500
04501 REAL*4 AMTAU,AMNUTA,AMEL,AMNUE,AMMU,AMNUMU
04502 * ,AMPIZ,AMPI,AMRO,GAMRO,AMA1,GAMA1
04503 * ,AMK,AMKZ,AMKST,GAMKST
04504
04505 COMPLEX FPIKM
04506 FPIRK=CABS(FPIKM(W,AMK,AMKZ))**2
04507
04508 END
04509 COMPLEX FUNCTION FPIKMK(W,XM1,XM2)
04510
04511
04512
04513 COMPLEX BWIGM
04514 REAL ROM,ROG,ROM1,ROG1,BETA1,PI,PIM,S,W
04515 EXTERNAL BWIG
04516 DATA INIT /0/
04517
04518
04519 IF (INIT.EQ.0 ) THEN
04520 INIT=1
04521 PI=3.141592654
04522 PIM=.140
04523 ROM=0.773
04524 ROG=0.145
04525 ROM1=1.570
04526 ROG1=0.510
04527
04528 BETA1=-0.221
04529 ENDIF
04530
04531 S=W**2
04532 FPIKMK=(BWIGM(S,ROM,ROG,XM1,XM2)+BETA1*BWIGM(S,ROM1,ROG1,XM1,XM2))
04533 & /(1+BETA1)
04534 RETURN
04535 END
04536 SUBROUTINE RESLUX
04537
04538
04539 NHEP=0
04540 END
04541 SUBROUTINE DWRPH(KTO,PHX)
04542
04543
04544
04545 IMPLICIT REAL*8 (A-H,O-Z)
04546 REAL*4 PHX(4)
04547 REAL*4 QHOT(4)
04548
04549 DO 9 K=1,4
04550 QHOT(K) =0.0
04551 9 CONTINUE
04552
04553
04554 DO 1002 I=1,4
04555 1002 QHOT(I)=PHX(I)
04556 IF (QHOT(4).GT.1.E-5) CALL DWLUPH(KTO,QHOT)
04557 RETURN
04558 END
04559 SUBROUTINE DWLUPH(KTO,PHOT)
04560
04561
04562
04563
04564
04565
04566
04567
04568
04569 REAL PHOT(4)
04570 COMMON /TAUPOS/ NP1,NP2
04571
04572
04573 IF (PHOT(4).LE.0.0) RETURN
04574
04575
04576 IF((KTO.EQ. 1).OR.(KTO.EQ.11)) THEN
04577 NPS=NP1
04578 ELSE
04579 NPS=NP2
04580 ENDIF
04581
04582 KTOS=KTO
04583 IF(KTOS.GT.10) KTOS=KTOS-10
04584
04585 CALL TRALO4(KTOS,PHOT,PHOT,AM)
04586 CALL FILHEP(0,1,22,NPS,NPS,0,0,PHOT,0.0,.TRUE.)
04587
04588 RETURN
04589 END
04590
04591 SUBROUTINE DWLUEL(KTO,ISGN,PNU,PWB,PEL,PNE)
04592
04593
04594
04595
04596
04597
04598
04599
04600
04601 REAL PNU(4),PWB(4),PEL(4),PNE(4)
04602 COMMON /TAUPOS/ NP1,NP2
04603
04604
04605 IF(KTO.EQ. 1) THEN
04606 NPS=NP1
04607 ELSE
04608 NPS=NP2
04609 ENDIF
04610
04611
04612 CALL TRALO4(KTO,PNU,PNU,AM)
04613 CALL FILHEP(0,1,16*ISGN,NPS,NPS,0,0,PNU,AM,.TRUE.)
04614
04615
04616 CALL TRALO4(KTO,PWB,PWB,AM)
04617
04618
04619
04620 CALL TRALO4(KTO,PEL,PEL,AM)
04621 CALL FILHEP(0,1,11*ISGN,NPS,NPS,0,0,PEL,AM,.FALSE.)
04622
04623
04624 CALL TRALO4(KTO,PNE,PNE,AM)
04625 CALL FILHEP(0,1,-12*ISGN,NPS,NPS,0,0,PNE,AM,.TRUE.)
04626
04627 RETURN
04628 END
04629 SUBROUTINE DWLUMU(KTO,ISGN,PNU,PWB,PMU,PNM)
04630
04631
04632
04633
04634
04635
04636
04637
04638
04639 REAL PNU(4),PWB(4),PMU(4),PNM(4)
04640 COMMON /TAUPOS/ NP1,NP2
04641
04642
04643 IF(KTO.EQ. 1) THEN
04644 NPS=NP1
04645 ELSE
04646 NPS=NP2
04647 ENDIF
04648
04649
04650 CALL TRALO4(KTO,PNU,PNU,AM)
04651 CALL FILHEP(0,1,16*ISGN,NPS,NPS,0,0,PNU,AM,.TRUE.)
04652
04653
04654 CALL TRALO4(KTO,PWB,PWB,AM)
04655
04656
04657
04658 CALL TRALO4(KTO,PMU,PMU,AM)
04659 CALL FILHEP(0,1,13*ISGN,NPS,NPS,0,0,PMU,AM,.FALSE.)
04660
04661
04662 CALL TRALO4(KTO,PNM,PNM,AM)
04663 CALL FILHEP(0,1,-14*ISGN,NPS,NPS,0,0,PNM,AM,.TRUE.)
04664
04665 RETURN
04666 END
04667 SUBROUTINE DWLUPI(KTO,ISGN,PPI,PNU)
04668
04669
04670
04671
04672
04673
04674
04675
04676
04677 REAL PNU(4),PPI(4)
04678 COMMON /TAUPOS/ NP1,NP2
04679
04680
04681 IF(KTO.EQ. 1) THEN
04682 NPS=NP1
04683 ELSE
04684 NPS=NP2
04685 ENDIF
04686
04687
04688 CALL TRALO4(KTO,PNU,PNU,AM)
04689 CALL FILHEP(0,1,16*ISGN,NPS,NPS,0,0,PNU,AM,.TRUE.)
04690
04691
04692 CALL TRALO4(KTO,PPI,PPI,AM)
04693 CALL FILHEP(0,1,-211*ISGN,NPS,NPS,0,0,PPI,AM,.TRUE.)
04694
04695 RETURN
04696 END
04697 SUBROUTINE DWLURO(KTO,ISGN,PNU,PRHO,PIC,PIZ)
04698
04699
04700
04701
04702
04703
04704
04705
04706
04707 REAL PNU(4),PRHO(4),PIC(4),PIZ(4)
04708 COMMON /TAUPOS/ NP1,NP2
04709
04710
04711 IF(KTO.EQ. 1) THEN
04712 NPS=NP1
04713 ELSE
04714 NPS=NP2
04715 ENDIF
04716
04717
04718 CALL TRALO4(KTO,PNU,PNU,AM)
04719 CALL FILHEP(0,1,16*ISGN,NPS,NPS,0,0,PNU,AM,.TRUE.)
04720
04721
04722 CALL TRALO4(KTO,PRHO,PRHO,AM)
04723 CALL FILHEP(0,2,-213*ISGN,NPS,NPS,0,0,PRHO,AM,.TRUE.)
04724
04725
04726 CALL TRALO4(KTO,PIC,PIC,AM)
04727 CALL FILHEP(0,1,-211*ISGN,-1,-1,0,0,PIC,AM,.TRUE.)
04728
04729
04730 CALL TRALO4(KTO,PIZ,PIZ,AM)
04731 CALL FILHEP(0,1,111,-2,-2,0,0,PIZ,AM,.TRUE.)
04732
04733 RETURN
04734 END
04735 SUBROUTINE DWLUAA(KTO,ISGN,PNU,PAA,PIM1,PIM2,PIPL,JAA)
04736
04737
04738
04739
04740
04741
04742
04743
04744
04745
04746 REAL PNU(4),PAA(4),PIM1(4),PIM2(4),PIPL(4)
04747 COMMON /TAUPOS/ NP1,NP2
04748
04749
04750 IF(KTO.EQ. 1) THEN
04751 NPS=NP1
04752 ELSE
04753 NPS=NP2
04754 ENDIF
04755
04756
04757 CALL TRALO4(KTO,PNU,PNU,AM)
04758 CALL FILHEP(0,1,16*ISGN,NPS,NPS,0,0,PNU,AM,.TRUE.)
04759
04760
04761 CALL TRALO4(KTO,PAA,PAA,AM)
04762 CALL FILHEP(0,1,-20213*ISGN,NPS,NPS,0,0,PAA,AM,.TRUE.)
04763
04764
04765 IF(JAA.EQ.1) THEN
04766
04767
04768
04769
04770 CALL TRALO4(KTO,PIM2,PIM2,AM)
04771 CALL FILHEP(0,1,-211*ISGN,-1,-1,0,0,PIM2,AM,.TRUE.)
04772
04773
04774 CALL TRALO4(KTO,PIM1,PIM1,AM)
04775 CALL FILHEP(0,1,-211*ISGN,-2,-2,0,0,PIM1,AM,.TRUE.)
04776
04777
04778 CALL TRALO4(KTO,PIPL,PIPL,AM)
04779 CALL FILHEP(0,1, 211*ISGN,-3,-3,0,0,PIPL,AM,.TRUE.)
04780
04781 ELSE IF (JAA.EQ.2) THEN
04782
04783
04784
04785
04786 CALL TRALO4(KTO,PIM2,PIM2,AM)
04787 CALL FILHEP(0,1,111,-1,-1,0,0,PIM2,AM,.TRUE.)
04788
04789
04790 CALL TRALO4(KTO,PIM1,PIM1,AM)
04791 CALL FILHEP(0,1,111,-2,-2,0,0,PIM1,AM,.TRUE.)
04792
04793
04794 CALL TRALO4(KTO,PIPL,PIPL,AM)
04795 CALL FILHEP(0,1,-211*ISGN,-3,-3,0,0,PIPL,AM,.TRUE.)
04796
04797 ENDIF
04798
04799 RETURN
04800 END
04801 SUBROUTINE DWLUKK (KTO,ISGN,PKK,PNU)
04802
04803
04804
04805
04806
04807
04808
04809
04810 REAL PKK(4),PNU(4)
04811 COMMON /TAUPOS/ NP1,NP2
04812
04813
04814 IF (KTO.EQ.1) THEN
04815 NPS=NP1
04816 ELSE
04817 NPS=NP2
04818 ENDIF
04819
04820
04821 CALL TRALO4 (KTO,PNU,PNU,AM)
04822 CALL FILHEP(0,1,16*ISGN,NPS,NPS,0,0,PNU,AM,.TRUE.)
04823
04824
04825 CALL TRALO4 (KTO,PKK,PKK,AM)
04826 CALL FILHEP(0,1,-321*ISGN,NPS,NPS,0,0,PKK,AM,.TRUE.)
04827
04828 RETURN
04829 END
04830 SUBROUTINE DWLUKS(KTO,ISGN,PNU,PKS,PKK,PPI,JKST)
04831 COMMON / TAUKLE / BRA1,BRK0,BRK0B,BRKS
04832 REAL*4 BRA1,BRK0,BRK0B,BRKS
04833
04834
04835
04836
04837
04838
04839
04840
04841
04842 REAL PNU(4),PKS(4),PKK(4),PPI(4),XIO(1)
04843 COMMON /TAUPOS/ NP1,NP2
04844
04845
04846 IF(KTO.EQ. 1) THEN
04847 NPS=NP1
04848 ELSE
04849 NPS=NP2
04850 ENDIF
04851
04852
04853 CALL TRALO4(KTO,PNU,PNU,AM)
04854 CALL FILHEP(0,1,16*ISGN,NPS,NPS,0,0,PNU,AM,.TRUE.)
04855
04856
04857 CALL TRALO4(KTO,PKS,PKS,AM)
04858 CALL FILHEP(0,1,-323*ISGN,NPS,NPS,0,0,PKS,AM,.TRUE.)
04859
04860
04861 IF(JKST.EQ.10) THEN
04862
04863
04864
04865
04866 CALL TRALO4(KTO,PPI,PPI,AM)
04867 CALL FILHEP(0,1,-211*ISGN,-1,-1,0,0,PPI,AM,.TRUE.)
04868
04869 BRAN=BRK0B
04870 IF (ISGN.EQ.-1) BRAN=BRK0
04871
04872 CALL RANMAR(XIO,1)
04873 IF(XIO(1).GT.BRAN) THEN
04874 K0TYPE = 130
04875 ELSE
04876 K0TYPE = 310
04877 ENDIF
04878
04879 CALL TRALO4(KTO,PKK,PKK,AM)
04880 CALL FILHEP(0,1,K0TYPE,-2,-2,0,0,PKK,AM,.TRUE.)
04881
04882 ELSE IF(JKST.EQ.20) THEN
04883
04884
04885
04886
04887 CALL TRALO4(KTO,PPI,PPI,AM)
04888 CALL FILHEP(0,1,111,-1,-1,0,0,PPI,AM,.TRUE.)
04889
04890
04891 CALL TRALO4(KTO,PKK,PKK,AM)
04892 CALL FILHEP(0,1,-321*ISGN,-2,-2,0,0,PKK,AM,.TRUE.)
04893
04894 ENDIF
04895
04896 RETURN
04897 END
04898 SUBROUTINE DWLNEW(KTO,ISGN,PNU,PWB,PNPI,MODE)
04899
04900
04901
04902
04903
04904
04905
04906
04907
04908 PARAMETER (NMODE=15,NM1=0,NM2=1,NM3=8,NM4=2,NM5=1,NM6=3)
04909 COMMON / TAUDCD /IDFFIN(9,NMODE),MULPIK(NMODE)
04910 & ,NAMES
04911 COMMON /TAUPOS/ NP1,NP2
04912 CHARACTER NAMES(NMODE)*31
04913 REAL PNU(4),PWB(4),PNPI(4,9)
04914 REAL PPI(4)
04915
04916 JNPI=MODE-7
04917
04918 IF(KTO.EQ. 1) THEN
04919 NPS=NP1
04920 ELSE
04921 NPS=NP2
04922 ENDIF
04923
04924
04925 CALL TRALO4(KTO,PNU,PNU,AM)
04926 CALL FILHEP(0,1,16*ISGN,NPS,NPS,0,0,PNU,AM,.TRUE.)
04927
04928
04929 CALL TRALO4(KTO,PWB,PWB,AM)
04930 CALL FILHEP(0,1,-24*ISGN,NPS,NPS,0,0,PWB,AM,.TRUE.)
04931
04932
04933
04934
04935 ND=MULPIK(JNPI)
04936 DO I=1,ND
04937 KFPI=LUNPIK(IDFFIN(I,JNPI),-ISGN)
04938
04939
04940 DO J=1,4
04941 PPI(J)=PNPI(J,I)
04942 END DO
04943 CALL TRALO4(KTO,PPI,PPI,AM)
04944 CALL FILHEP(0,1,KFPI,-I,-I,0,0,PPI,AM,.TRUE.)
04945 END DO
04946
04947 RETURN
04948 END
04949 FUNCTION AMAST(PP)
04950
04951
04952
04953
04954
04955 IMPLICIT REAL*8 (A-H,O-Z)
04956 REAL*8 PP(4)
04957 AAA=PP(4)**2-PP(3)**2-PP(2)**2-PP(1)**2
04958
04959 IF(AAA.NE.0.0) AAA=AAA/SQRT(ABS(AAA))
04960 AMAST=AAA
04961 RETURN
04962 END
04963 FUNCTION AMAS4(PP)
04964
04965
04966
04967
04968
04969
04970 REAL PP(4)
04971 AAA=PP(4)**2-PP(3)**2-PP(2)**2-PP(1)**2
04972 IF(AAA.NE.0.0) AAA=AAA/SQRT(ABS(AAA))
04973 AMAS4=AAA
04974 RETURN
04975 END
04976 FUNCTION ANGXY(X,Y)
04977
04978
04979
04980
04981 IMPLICIT DOUBLE PRECISION (A-H,O-Z)
04982 DATA PI /3.141592653589793238462643D0/
04983
04984 IF(ABS(Y).LT.ABS(X)) THEN
04985 THE=ATAN(ABS(Y/X))
04986 IF(X.LE.0D0) THE=PI-THE
04987 ELSE
04988 THE=ACOS(X/SQRT(X**2+Y**2))
04989 ENDIF
04990 ANGXY=THE
04991 RETURN
04992 END
04993 FUNCTION ANGFI(X,Y)
04994
04995
04996
04997
04998
04999 IMPLICIT DOUBLE PRECISION (A-H,O-Z)
05000 DATA PI /3.141592653589793238462643D0/
05001
05002 IF(ABS(Y).LT.ABS(X)) THEN
05003 THE=ATAN(ABS(Y/X))
05004 IF(X.LE.0D0) THE=PI-THE
05005 ELSE
05006 THE=ACOS(X/SQRT(X**2+Y**2))
05007 ENDIF
05008 IF(Y.LT.0D0) THE=2D0*PI-THE
05009 ANGFI=THE
05010 END
05011 SUBROUTINE ROTOD1(PH1,PVEC,QVEC)
05012
05013
05014
05015
05016 IMPLICIT DOUBLE PRECISION (A-H,O-Z)
05017 DIMENSION PVEC(4),QVEC(4),RVEC(4)
05018
05019 PHI=PH1
05020 CS=COS(PHI)
05021 SN=SIN(PHI)
05022 DO 10 I=1,4
05023 10 RVEC(I)=PVEC(I)
05024 QVEC(1)=RVEC(1)
05025 QVEC(2)= CS*RVEC(2)-SN*RVEC(3)
05026 QVEC(3)= SN*RVEC(2)+CS*RVEC(3)
05027 QVEC(4)=RVEC(4)
05028 RETURN
05029 END
05030 SUBROUTINE ROTOD2(PH1,PVEC,QVEC)
05031
05032
05033
05034
05035 IMPLICIT DOUBLE PRECISION (A-H,O-Z)
05036 DIMENSION PVEC(4),QVEC(4),RVEC(4)
05037
05038 PHI=PH1
05039 CS=COS(PHI)
05040 SN=SIN(PHI)
05041 DO 10 I=1,4
05042 10 RVEC(I)=PVEC(I)
05043 QVEC(1)= CS*RVEC(1)+SN*RVEC(3)
05044 QVEC(2)=RVEC(2)
05045 QVEC(3)=-SN*RVEC(1)+CS*RVEC(3)
05046 QVEC(4)=RVEC(4)
05047 RETURN
05048 END
05049 SUBROUTINE ROTOD3(PH1,PVEC,QVEC)
05050
05051
05052
05053
05054 IMPLICIT DOUBLE PRECISION (A-H,O-Z)
05055
05056 DIMENSION PVEC(4),QVEC(4),RVEC(4)
05057 PHI=PH1
05058 CS=COS(PHI)
05059 SN=SIN(PHI)
05060 DO 10 I=1,4
05061 10 RVEC(I)=PVEC(I)
05062 QVEC(1)= CS*RVEC(1)-SN*RVEC(2)
05063 QVEC(2)= SN*RVEC(1)+CS*RVEC(2)
05064 QVEC(3)=RVEC(3)
05065 QVEC(4)=RVEC(4)
05066 END
05067 SUBROUTINE BOSTR3(EXE,PVEC,QVEC)
05068
05069
05070
05071
05072
05073 REAL*4 PVEC(4),QVEC(4),RVEC(4)
05074
05075 DO 10 I=1,4
05076 10 RVEC(I)=PVEC(I)
05077 RPL=RVEC(4)+RVEC(3)
05078 RMI=RVEC(4)-RVEC(3)
05079 QPL=RPL*EXE
05080 QMI=RMI/EXE
05081 QVEC(1)=RVEC(1)
05082 QVEC(2)=RVEC(2)
05083 QVEC(3)=(QPL-QMI)/2
05084 QVEC(4)=(QPL+QMI)/2
05085 END
05086 SUBROUTINE BOSTD3(EXE,PVEC,QVEC)
05087
05088
05089
05090
05091
05092 IMPLICIT DOUBLE PRECISION (A-H,O-Z)
05093 DIMENSION PVEC(4),QVEC(4),RVEC(4)
05094
05095 DO 10 I=1,4
05096 10 RVEC(I)=PVEC(I)
05097 RPL=RVEC(4)+RVEC(3)
05098 RMI=RVEC(4)-RVEC(3)
05099 QPL=RPL*EXE
05100 QMI=RMI/EXE
05101 QVEC(1)=RVEC(1)
05102 QVEC(2)=RVEC(2)
05103 QVEC(3)=(QPL-QMI)/2
05104 QVEC(4)=(QPL+QMI)/2
05105 RETURN
05106 END
05107 SUBROUTINE ROTOR1(PH1,PVEC,QVEC)
05108
05109
05110
05111
05112 REAL*4 PVEC(4),QVEC(4),RVEC(4)
05113
05114 PHI=PH1
05115 CS=COS(PHI)
05116 SN=SIN(PHI)
05117 DO 10 I=1,4
05118 10 RVEC(I)=PVEC(I)
05119 QVEC(1)=RVEC(1)
05120 QVEC(2)= CS*RVEC(2)-SN*RVEC(3)
05121 QVEC(3)= SN*RVEC(2)+CS*RVEC(3)
05122 QVEC(4)=RVEC(4)
05123 END
05124 SUBROUTINE ROTOR2(PH1,PVEC,QVEC)
05125
05126
05127
05128
05129 IMPLICIT REAL*4(A-H,O-Z)
05130 REAL*4 PVEC(4),QVEC(4),RVEC(4)
05131
05132 PHI=PH1
05133 CS=COS(PHI)
05134 SN=SIN(PHI)
05135 DO 10 I=1,4
05136 10 RVEC(I)=PVEC(I)
05137 QVEC(1)= CS*RVEC(1)+SN*RVEC(3)
05138 QVEC(2)=RVEC(2)
05139 QVEC(3)=-SN*RVEC(1)+CS*RVEC(3)
05140 QVEC(4)=RVEC(4)
05141 END
05142 SUBROUTINE ROTOR3(PHI,PVEC,QVEC)
05143
05144
05145
05146
05147 REAL*4 PVEC(4),QVEC(4),RVEC(4)
05148
05149 CS=COS(PHI)
05150 SN=SIN(PHI)
05151 DO 10 I=1,4
05152 10 RVEC(I)=PVEC(I)
05153 QVEC(1)= CS*RVEC(1)-SN*RVEC(2)
05154 QVEC(2)= SN*RVEC(1)+CS*RVEC(2)
05155 QVEC(3)=RVEC(3)
05156 QVEC(4)=RVEC(4)
05157 END
05158 SUBROUTINE SPHERD(R,X)
05159
05160
05161
05162
05163 REAL*8 R,X(4),PI,COSTH,SINTH
05164 REAL*4 RRR(2)
05165 DATA PI /3.141592653589793238462643D0/
05166
05167 CALL RANMAR(RRR,2)
05168 COSTH=-1+2*RRR(1)
05169 SINTH=SQRT(1 -COSTH**2)
05170 X(1)=R*SINTH*COS(2*PI*RRR(2))
05171 X(2)=R*SINTH*SIN(2*PI*RRR(2))
05172 X(3)=R*COSTH
05173 RETURN
05174 END
05175 SUBROUTINE ROTPOX(THET,PHI,PP)
05176 IMPLICIT REAL*8 (A-H,O-Z)
05177
05178
05179
05180 DIMENSION PP(4)
05181
05182 CALL ROTOD2(THET,PP,PP)
05183 CALL ROTOD3( PHI,PP,PP)
05184 RETURN
05185 END
05186 SUBROUTINE SPHERA(R,X)
05187
05188
05189
05190
05191
05192 REAL X(4)
05193 REAL*4 RRR(2)
05194 DATA PI /3.141592653589793238462643/
05195
05196 CALL RANMAR(RRR,2)
05197 COSTH=-1.+2.*RRR(1)
05198 SINTH=SQRT(1.-COSTH**2)
05199 X(1)=R*SINTH*COS(2*PI*RRR(2))
05200 X(2)=R*SINTH*SIN(2*PI*RRR(2))
05201 X(3)=R*COSTH
05202 RETURN
05203 END
05204 SUBROUTINE ROTPOL(THET,PHI,PP)
05205
05206
05207
05208
05209 REAL PP(4)
05210
05211 CALL ROTOR2(THET,PP,PP)
05212 CALL ROTOR3( PHI,PP,PP)
05213 RETURN
05214 END
05215 SUBROUTINE RANMAR(RVEC,LENV)
05216
05217
05218
05219
05220
05221
05222
05223
05224
05225
05226
05227
05228
05229
05230
05231 COMMON / INOUT / INUT,IOUT
05232 DIMENSION RVEC(*)
05233 COMMON/RASET1/U(97),C,I97,J97
05234 PARAMETER (MODCNS=1000000000)
05235 DATA NTOT,NTOT2,IJKL/-1,0,0/
05236
05237 IF (NTOT .GE. 0) GO TO 50
05238
05239
05240 IJKL = 54217137
05241 NTOT = 0
05242 NTOT2 = 0
05243 KALLED = 0
05244 GO TO 1
05245
05246 ENTRY RMARIN(IJKLIN, NTOTIN,NTOT2N)
05247
05248
05249
05250
05251
05252
05253
05254 IJKL = IJKLIN
05255 NTOT = MAX(NTOTIN,0)
05256 NTOT2= MAX(NTOT2N,0)
05257 KALLED = 1
05258
05259 1 CONTINUE
05260 IJ = IJKL/30082
05261 KL = IJKL - 30082*IJ
05262 I = MOD(IJ/177, 177) + 2
05263 J = MOD(IJ, 177) + 2
05264 K = MOD(KL/169, 178) + 1
05265 L = MOD(KL, 169)
05266 WRITE(IOUT,201) IJKL,NTOT,NTOT2
05267 201 FORMAT(1X,' RANMAR INITIALIZED: ',I10,2X,2I10)
05268 DO 2 II= 1, 97
05269 S = 0.
05270 T = .5
05271 DO 3 JJ= 1, 24
05272 M = MOD(MOD(I*J,179)*K, 179)
05273 I = J
05274 J = K
05275 K = M
05276 L = MOD(53*L+1, 169)
05277 IF (MOD(L*M,64) .GE. 32) S = S+T
05278 3 T = 0.5*T
05279 2 U(II) = S
05280 TWOM24 = 1.0
05281 DO 4 I24= 1, 24
05282 4 TWOM24 = 0.5*TWOM24
05283 C = 362436.*TWOM24
05284 CD = 7654321.*TWOM24
05285 CM = 16777213.*TWOM24
05286 I97 = 97
05287 J97 = 33
05288
05289
05290 DO 45 LOOP2= 1, NTOT2+1
05291 NOW = MODCNS
05292 IF (LOOP2 .EQ. NTOT2+1) NOW=NTOT
05293 IF (NOW .GT. 0) THEN
05294 WRITE (IOUT,'(A,I15)') ' RMARIN SKIPPING OVER ',NOW
05295 DO 40 IDUM = 1, NTOT
05296 UNI = U(I97)-U(J97)
05297 IF (UNI .LT. 0.) UNI=UNI+1.
05298 U(I97) = UNI
05299 I97 = I97-1
05300 IF (I97 .EQ. 0) I97=97
05301 J97 = J97-1
05302 IF (J97 .EQ. 0) J97=97
05303 C = C - CD
05304 IF (C .LT. 0.) C=C+CM
05305 40 CONTINUE
05306 ENDIF
05307 45 CONTINUE
05308 IF (KALLED .EQ. 1) RETURN
05309
05310
05311 50 CONTINUE
05312 DO 100 IVEC= 1, LENV
05313 UNI = U(I97)-U(J97)
05314 IF (UNI .LT. 0.) UNI=UNI+1.
05315 U(I97) = UNI
05316 I97 = I97-1
05317 IF (I97 .EQ. 0) I97=97
05318 J97 = J97-1
05319 IF (J97 .EQ. 0) J97=97
05320 C = C - CD
05321 IF (C .LT. 0.) C=C+CM
05322 UNI = UNI-C
05323 IF (UNI .LT. 0.) UNI=UNI+1.
05324
05325 IF (UNI .EQ. 0.) THEN
05326 UNI = TWOM24*U(2)
05327
05328 IF (UNI .EQ. 0.) UNI= TWOM24*TWOM24
05329 ENDIF
05330 RVEC(IVEC) = UNI
05331 100 CONTINUE
05332 NTOT = NTOT + LENV
05333 IF (NTOT .GE. MODCNS) THEN
05334 NTOT2 = NTOT2 + 1
05335 NTOT = NTOT - MODCNS
05336 ENDIF
05337 RETURN
05338
05339 ENTRY RMARUT(IJKLUT,NTOTUT,NTOT2T)
05340 IJKLUT = IJKL
05341 NTOTUT = NTOT
05342 NTOT2T = NTOT2
05343 RETURN
05344 END
05345 FUNCTION DILOGT(X)
05346
05347 IMPLICIT REAL*8(A-H,O-Z)
05348
05349 Z=-1.64493406684822
05350 IF(X .LT.-1.0) GO TO 1
05351 IF(X .LE. 0.5) GO TO 2
05352 IF(X .EQ. 1.0) GO TO 3
05353 IF(X .LE. 2.0) GO TO 4
05354 Z=3.2898681336964
05355 1 T=1.0/X
05356 S=-0.5
05357 Z=Z-0.5* LOG(ABS(X))**2
05358 GO TO 5
05359 2 T=X
05360 S=0.5
05361 Z=0.
05362 GO TO 5
05363 3 DILOGT=1.64493406684822
05364 RETURN
05365 4 T=1.0-X
05366 S=-0.5
05367 Z=1.64493406684822 - LOG(X)* LOG(ABS(T))
05368 5 Y=2.66666666666666 *T+0.66666666666666
05369 B= 0.00000 00000 00001
05370 A=Y*B +0.00000 00000 00004
05371 B=Y*A-B+0.00000 00000 00011
05372 A=Y*B-A+0.00000 00000 00037
05373 B=Y*A-B+0.00000 00000 00121
05374 A=Y*B-A+0.00000 00000 00398
05375 B=Y*A-B+0.00000 00000 01312
05376 A=Y*B-A+0.00000 00000 04342
05377 B=Y*A-B+0.00000 00000 14437
05378 A=Y*B-A+0.00000 00000 48274
05379 B=Y*A-B+0.00000 00001 62421
05380 A=Y*B-A+0.00000 00005 50291
05381 B=Y*A-B+0.00000 00018 79117
05382 A=Y*B-A+0.00000 00064 74338
05383 B=Y*A-B+0.00000 00225 36705
05384 A=Y*B-A+0.00000 00793 87055
05385 B=Y*A-B+0.00000 02835 75385
05386 A=Y*B-A+0.00000 10299 04264
05387 B=Y*A-B+0.00000 38163 29463
05388 A=Y*B-A+0.00001 44963 00557
05389 B=Y*A-B+0.00005 68178 22718
05390 A=Y*B-A+0.00023 20021 96094
05391 B=Y*A-B+0.00100 16274 96164
05392 A=Y*B-A+0.00468 63619 59447
05393 B=Y*A-B+0.02487 93229 24228
05394 A=Y*B-A+0.16607 30329 27855
05395 A=Y*A-B+1.93506 43008 6996
05396 DILOGT=S*T*(A-B)+Z
05397 RETURN
05398
05399
05400
05401 END
05402
05403
05404
05405
05406
05407
05408
05409 SUBROUTINE INIRCHL(FLAG)
05410 INTEGER FLAG
05411
05412 IF(FLAG.NE.0) THEN
05413 WRITE(*,25) FLAG
05414 25 FORMAT(1X, "TAUOLA IniRChL: Fatal error, FLAG=",I2," but RChL currents missing")
05415 WRITE(*,*) " in loaded version of TAUOLA-FORTRAN library."
05416 STOP
05417 ENDIF
05418
05419 RETURN
05420 END
