00001
00002
00003
00004
00005
00006
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
00048
00049
00050
00051
00052
00053
00054
00055
00056
00057
00058
00059
00060
00061
00062
00063
00064
00065
00066
00067
00068
00069
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 SUBROUTINE PHOINI
00095
00096
00097
00098
00099
00100
00101
00102
00103
00104
00105
00106
00107
00108
00109
00110
00111 IMPLICIT NONE
00112 INTEGER INIT,IDUM,IPHQRK,IPHEKL
00113 SAVE INIT
00114 DATA INIT/ 0/
00115
00116
00117 IF (INIT.NE.0) RETURN
00118 INIT=1
00119
00120
00121
00122
00123
00124
00125
00126
00127 CALL PHCORK(1)
00128
00129
00130
00131 IDUM= IPHQRK(1)
00132
00133
00134 IDUM= IPHEKL(2)
00135
00136
00137 CALL PHOCIN
00138
00139
00140 CALL PHOINF
00141
00142
00143
00144 CALL PHORIN
00145 RETURN
00146 END
00147 SUBROUTINE PHOCIN
00148
00149
00150
00151
00152
00153
00154
00155
00156
00157
00158
00159
00160
00161
00162
00163 IMPLICIT NONE
00164 INTEGER d_h_NMXHEP
00165 PARAMETER (d_h_NMXHEP=10000)
00166 LOGICAL QEDRAD
00167 COMMON/PHOQED/QEDRAD(d_h_NMXHEP)
00168 INTEGER PHLUN
00169 COMMON/PHOLUN/PHLUN
00170 REAL*8 ALPHA,XPHCUT
00171 COMMON/PHOCOP/ALPHA,XPHCUT
00172 REAL*8 PI,TWOPI
00173 COMMON/PHPICO/PI,TWOPI
00174 INTEGER ISEED,I97,J97
00175 REAL*8 URAN,CRAN,CDRAN,CMRAN
00176 COMMON/PHSEED/ISEED(2),I97,J97,URAN(97),CRAN,CDRAN,CMRAN
00177 INTEGER PHOMES
00178 PARAMETER (PHOMES=10)
00179 INTEGER STATUS
00180 COMMON/PHOSTA/STATUS(PHOMES)
00181 LOGICAL INTERF,ISEC,ITRE,IEXP,IFTOP,IFW
00182 REAL*8 FINT,FSEC,EXPEPS
00183 COMMON /PHOKEY/ FSEC,FINT,EXPEPS,INTERF,ISEC,ITRE,IEXP,IFTOP,IFW
00184 INTEGER INIT,I
00185 SAVE INIT
00186 DATA INIT/ 0/
00187
00188
00189 IF (INIT.NE.0) RETURN
00190 INIT=1
00191
00192
00193
00194 DO 10 I=1,d_h_NMXHEP
00195 10 QEDRAD(I)=.TRUE.
00196
00197
00198 PHLUN=6
00199
00200
00201 XPHCUT=0.01 D0
00202
00203
00204
00205 ALPHA=0.00729735039D0
00206 PI=3.14159265358979324D0
00207 TWOPI=6.28318530717958648D0
00208
00209
00210 ISEED(1)=1802
00211 ISEED(2)=9373
00212
00213
00214
00215
00216 INTERF=.TRUE.
00217
00218
00219 ISEC=.TRUE.
00220
00221
00222 ITRE=.FALSE.
00223
00224 IEXP=.FALSE.
00225 IF (IEXP) THEN
00226 ISEC=.FALSE.
00227 ITRE=.FALSE.
00228 CALL PHCORK(5)
00229
00230 XPHCUT=0.000 000 1
00231 EXPEPS=1D-4
00232 ENDIF
00233
00234
00235
00236 IFTOP=.TRUE.
00237
00238
00239
00240
00241
00242
00243
00244 IF (INTERF) THEN
00245
00246
00247
00248 FINT=2.0D0
00249 ELSE
00250 FINT=1.0D0
00251 ENDIF
00252
00253
00254
00255
00256
00257
00258
00259
00260
00261
00262
00263
00264
00265
00266 IFW=.TRUE.
00267
00268 DO 20 I=1,PHOMES
00269 20 STATUS(I)=0
00270 RETURN
00271 END
00272 SUBROUTINE PHOINF
00273
00274
00275
00276
00277
00278
00279
00280
00281
00282
00283
00284
00285
00286
00287 IMPLICIT NONE
00288 INTEGER IV1,IV2,IV3
00289 INTEGER PHOVN1,PHOVN2
00290 COMMON/PHOVER/PHOVN1,PHOVN2
00291 INTEGER PHLUN
00292 COMMON/PHOLUN/PHLUN
00293 LOGICAL INTERF,ISEC,ITRE,IEXP,IFTOP,IFW
00294 REAL*8 FINT,FSEC,EXPEPS
00295 COMMON /PHOKEY/ FSEC,FINT,EXPEPS,INTERF,ISEC,ITRE,IEXP,IFTOP,IFW
00296 REAL*8 ALPHA,XPHCUT
00297 COMMON/PHOCOP/ALPHA,XPHCUT
00298
00299
00300 PHOVN1=215
00301 PHOVN2=111005
00302
00303
00304 WRITE(PHLUN,9000)
00305 WRITE(PHLUN,9020)
00306 WRITE(PHLUN,9010)
00307 WRITE(PHLUN,9030)
00308 IV1=PHOVN1/100
00309 IV2=PHOVN1-IV1*100
00310 WRITE(PHLUN,9040) IV1,IV2
00311 IV1=PHOVN2/10000
00312 IV2=(PHOVN2-IV1*10000)/100
00313 IV3=PHOVN2-IV1*10000-IV2*100
00314 WRITE(PHLUN,9050) IV1,IV2,IV3
00315 WRITE(PHLUN,9030)
00316 WRITE(PHLUN,9010)
00317 WRITE(PHLUN,9060)
00318 WRITE(PHLUN,9010)
00319 WRITE(PHLUN,9070)
00320 WRITE(PHLUN,9010)
00321 WRITE(PHLUN,9020)
00322 WRITE(PHLUN,9010)
00323 WRITE(PHLUN,9064) INTERF,ISEC,ITRE,IEXP,IFTOP,IFW,ALPHA,XPHCUT
00324 WRITE(PHLUN,9010)
00325 IF (INTERF) WRITE(PHLUN,9061)
00326 IF (ISEC) WRITE(PHLUN,9062)
00327 IF (ITRE) WRITE(PHLUN,9066)
00328 IF (IEXP) WRITE(PHLUN,9067) EXPEPS
00329 IF (IFTOP) WRITE(PHLUN,9063)
00330 IF (IFW) WRITE(PHLUN,9065)
00331 WRITE(PHLUN,9080)
00332 WRITE(PHLUN,9010)
00333 WRITE(PHLUN,9020)
00334 RETURN
00335 9000 FORMAT(1H1)
00336 9010 FORMAT(1H ,'*',T81,'*')
00337 9020 FORMAT(1H ,80('*'))
00338 9030 FORMAT(1H ,'*',26X,26('='),T81,'*')
00339 9040 FORMAT(1H ,'*',28X,'PHOTOS, Version: ',I2,'.',I2,T81,'*')
00340 9050 FORMAT(1H ,'*',28X,'Released at: ',I2,'/',I2,'/',I2,T81,'*')
00341 9060 FORMAT(1H ,'*',18X,'PHOTOS QED Corrections in Particle Decays',
00342 &T81,'*')
00343 9061 FORMAT(1H ,'*',18X,'option with interference is active ',
00344 &T81,'*')
00345 9062 FORMAT(1H ,'*',18X,'option with double photons is active ',
00346 &T81,'*')
00347 9066 FORMAT(1H ,'*',18X,'option with triple/quatric photons is active',
00348 &T81,'*')
00349 9067 FORMAT(1H ,'*',18X,'option with exponentiation is active EPSEXP=',
00350 &E10.4,T81,'*')
00351 9063 FORMAT(1H ,'*',18X,'emision in t tbar production is active ',
00352 &T81,'*')
00353 9064 FORMAT(1H ,'*',18X,'Internal input parameters:',T81,'*'
00354 &,/, 1H ,'*',T81,'*'
00355 &,/, 1H ,'*',18X,'INTERF=',L2,' ISEC=',L2,' ITRE=',L2,
00356 & ' IEXP=',L2,' IFTOP=',L2,
00357 & ' IFW=',L2,T81,'*'
00358 &,/, 1H ,'*',18X,'ALPHA_QED=',F8.5,' XPHCUT=',E8.3,T81,'*')
00359 9065 FORMAT(1H ,'*',18X,'correction wt in decay of W is active ',
00360 &T81,'*')
00361 9070 FORMAT(1H ,'*',9X,
00362 &'Monte Carlo Program - by E. Barberio, B. van Eijk and Z. Was',
00363 & T81,'*',/,
00364 & 1H ,'*',9X,'Version 2.09 - by P. Golonka and Z.W.',T81,'*')
00365 9080 FORMAT( 1H ,'*',9X,' ',T81,'*',/,
00366 & 1H ,'*',9X,
00367 & ' WARNING (1): /HEPEVT/ is not anymore the standard common block'
00368 & ,T81,'*',/,
00369 & 1H ,'*',9X,' ',T81,'*',/,
00370 & 1H ,'*',9X,
00371 & ' PHOTOS expects /HEPEVT/ to have REAL*8 variables. To change to'
00372 & ,T81,'*',/, 1H ,'*',9X,
00373 & ' REAL*4 modify its declaration in subr. PHOTOS_GET PHOTOS_SET:'
00374 & ,T81,'*',/, 1H ,'*',9X,
00375 & ' REAL*8 d_h_phep, d_h_vhep'
00376 & ,T81,'*',/, 1H ,'*',9X,
00377 & ' WARNING (2): check dims. of /hepevt/ /phoqed/ /ph_hepevt/.'
00378 & ,T81,'*',/, 1H ,'*',9X,
00379 & ' HERE: d_h_nmxhep=10000 and NMXHEP=10000'
00380 & ,T81,'*')
00381 END
00382 SUBROUTINE PHOTOS(ID)
00383
00384
00385
00386
00387
00388
00389
00390
00391
00392
00393
00394
00395
00396
00397
00398
00399
00400 COMMON /PHLUPY/ IPOIN,IPOINM
00401 INTEGER IPOIN,IPOINM
00402 COMMON /PHNUM/ IEV
00403 INTEGER IEV
00404 INTEGER PHLUN
00405 COMMON/PHOLUN/PHLUN
00406
00407 IF (1.GT.IPOINM.AND.1.LT.IPOIN ) THEN
00408 WRITE(PHLUN,*) 'EVENT NR=',IEV,
00409 $ 'WE ARE TESTING /HEPEVT/ at IPOINT=1 (input)'
00410 CALL PHODMP
00411 ENDIF
00412 CALL PHOTOS_GET
00413 CALL PHOTOS_MAKE(ID)
00414 CALL PHOTOS_SET
00415 IF (1.GT.IPOINM.AND.1.LT.IPOIN ) THEN
00416 WRITE(PHLUN,*) 'EVENT NR=',IEV,
00417 $ 'WE ARE TESTING /HEPEVT/ at IPOINT=1 (output)'
00418 CALL PHODMP
00419 ENDIF
00420
00421 END
00422
00423 SUBROUTINE PHOTOS_GET
00424
00425
00426
00427
00428
00429
00430
00431
00432
00433
00434
00435
00436
00437
00438
00439
00440 IMPLICIT NONE
00441 INTEGER d_h_nmxhep
00442 PARAMETER (d_h_NMXHEP=10000)
00443 REAL*8 d_h_phep, d_h_vhep
00444 INTEGER d_h_nevhep,d_h_nhep,d_h_isthep,d_h_idhep,d_h_jmohep,
00445 $ d_h_jdahep
00446 COMMON /hepevt/
00447 $ d_h_nevhep,
00448 $ d_h_nhep,
00449 $ d_h_isthep(d_h_nmxhep),
00450 $ d_h_idhep(d_h_nmxhep),
00451 $ d_h_jmohep(2,d_h_nmxhep),
00452 $ d_h_jdahep(2,d_h_nmxhep),
00453 $ d_h_phep(5,d_h_nmxhep),
00454 $ d_h_vhep(4,d_h_nmxhep)
00455
00456 LOGICAL d_h_qedrad
00457 COMMON /phoqed/
00458 $ d_h_qedrad(d_h_nmxhep)
00459 INTEGER NMXHEP
00460 PARAMETER (NMXHEP=10000)
00461 INTEGER IDHEP,ISTHEP,JDAHEP,JMOHEP,NEVHEP,NHEP
00462 REAL*8 PHEP,VHEP
00463 COMMON/PH_HEPEVT/NEVHEP,NHEP,ISTHEP(NMXHEP),IDHEP(NMXHEP),
00464 &JMOHEP(2,NMXHEP),JDAHEP(2,NMXHEP),PHEP(5,NMXHEP),VHEP(4,NMXHEP)
00465 LOGICAL QEDRAD
00466 COMMON/PH_PHOQED/QEDRAD(NMXHEP)
00467 integer k,l
00468 nevhep= d_h_nevhep
00469 nhep = d_h_nhep
00470 DO K=1,nhep
00471 isthep(k) =d_h_isthep(k)
00472 idhep(k) =d_h_idhep(k)
00473 jmohep(1,k) =d_h_jmohep(1,k)
00474 jdahep(1,k) =d_h_jdahep(1,k)
00475 jmohep(2,k) =d_h_jmohep(2,k)
00476 jdahep(2,k) =d_h_jdahep(2,k)
00477 DO l=1,4
00478 phep(l,k) =d_h_phep(l,k)
00479 vhep(l,k) =d_h_vhep(l,k)
00480 ENDDO
00481 phep(5,k) =d_h_phep(5,k)
00482 qedrad(k) =d_h_qedrad(k)
00483 ENDDO
00484 END
00485
00486
00487 SUBROUTINE PHOTOS_SET
00488
00489
00490
00491
00492
00493
00494
00495
00496
00497
00498
00499
00500
00501
00502
00503 IMPLICIT NONE
00504 INTEGER d_h_nmxhep
00505 PARAMETER (d_h_NMXHEP=10000)
00506 REAL*8 d_h_phep, d_h_vhep
00507 INTEGER d_h_nevhep,d_h_nhep,d_h_isthep,d_h_idhep,d_h_jmohep,
00508 $ d_h_jdahep
00509 COMMON /hepevt/
00510 $ d_h_nevhep,
00511 $ d_h_nhep,
00512 $ d_h_isthep(d_h_nmxhep),
00513 $ d_h_idhep(d_h_nmxhep),
00514 $ d_h_jmohep(2,d_h_nmxhep),
00515 $ d_h_jdahep(2,d_h_nmxhep),
00516 $ d_h_phep(5,d_h_nmxhep),
00517 $ d_h_vhep(4,d_h_nmxhep)
00518
00519 LOGICAL d_h_qedrad
00520 COMMON /phoqed/
00521 $ d_h_qedrad(d_h_nmxhep)
00522 INTEGER NMXHEP
00523 PARAMETER (NMXHEP=10000)
00524 INTEGER IDHEP,ISTHEP,JDAHEP,JMOHEP,NEVHEP,NHEP
00525 REAL*8 PHEP,VHEP
00526 COMMON/PH_HEPEVT/NEVHEP,NHEP,ISTHEP(NMXHEP),IDHEP(NMXHEP),
00527 &JMOHEP(2,NMXHEP),JDAHEP(2,NMXHEP),PHEP(5,NMXHEP),VHEP(4,NMXHEP)
00528 LOGICAL QEDRAD
00529 COMMON/PH_PHOQED/QEDRAD(NMXHEP)
00530 INTEGER K,L
00531
00532 d_h_nevhep= nevhep
00533 d_h_nhep = nhep
00534 DO K=1,nhep
00535 d_h_isthep(k) =isthep(k)
00536 d_h_idhep(k) =idhep(k)
00537 d_h_jmohep(1,k) =jmohep(1,k)
00538 d_h_jdahep(1,k) =jdahep(1,k)
00539 d_h_jmohep(2,k) =jmohep(2,k)
00540 d_h_jdahep(2,k) =jdahep(2,k)
00541 DO l=1,4
00542 d_h_phep(l,k) =phep(l,k)
00543 d_h_vhep(l,k) =vhep(l,k)
00544 ENDDO
00545 d_h_phep(5,k) =phep(5,k)
00546 d_h_qedrad(k) =qedrad(k)
00547 ENDDO
00548 END
00549 SUBROUTINE PHOTOS_MAKE(IPARR)
00550
00551
00552
00553
00554
00555
00556
00557
00558
00559
00560
00561
00562
00563
00564
00565
00566
00567
00568
00569
00570 IMPLICIT NONE
00571 REAL*8 PHOTON(5)
00572 INTEGER IP,IPARR,IPPAR,I,J,K,L,NLAST
00573 DOUBLE PRECISION DATA
00574 INTEGER MOTHER,POSPHO
00575 LOGICAL CASCAD
00576 INTEGER NMXHEP
00577 PARAMETER (NMXHEP=10000)
00578 INTEGER IDHEP,ISTHEP,JDAHEP,JMOHEP,NEVHEP,NHEP
00579 REAL*8 PHEP,VHEP
00580 COMMON/PH_HEPEVT/NEVHEP,NHEP,ISTHEP(NMXHEP),IDHEP(NMXHEP),
00581 &JMOHEP(2,NMXHEP),JDAHEP(2,NMXHEP),PHEP(5,NMXHEP),VHEP(4,NMXHEP)
00582 LOGICAL QEDRAD
00583 COMMON/PH_PHOQED/QEDRAD(NMXHEP)
00584 INTEGER NMXPHO
00585 PARAMETER (NMXPHO=10000)
00586 INTEGER ISTACK(0:NMXPHO),NUMIT,NTRY,KK,LL,II,NA,FIRST,LAST
00587 INTEGER FIRSTA,LASTA,IPP,IDA1,IDA2,MOTHER2,IDPHO,ISPHO
00588 REAL*8 PORIG(5,NMXPHO)
00589
00590 IPPAR=ABS(IPARR)
00591
00592 IP=IPPAR
00593 NLAST=NHEP
00594 CASCAD=.FALSE.
00595
00596
00597 IF ((JDAHEP(1,IP).EQ.0).OR.(JMOHEP(1,JDAHEP(1,IP)).NE.IP)) RETURN
00598
00599
00600
00601
00602 ISTACK(0)=IPPAR
00603
00604 NUMIT=0
00605
00606
00607 NTRY=0
00608
00609 IF (IPARR.GT.0) THEN
00610 30 CONTINUE
00611 DO I=JDAHEP(1,IP),JDAHEP(2,IP)
00612 IF (JDAHEP(1,I).NE.0.AND.JMOHEP(1,JDAHEP(1,I)).EQ.I) THEN
00613 NUMIT=NUMIT+1
00614 IF (NUMIT.GT.NMXPHO) THEN
00615 DATA=NUMIT
00616 CALL PHOERR(7,'PHOTOS',DATA)
00617 ENDIF
00618 ISTACK(NUMIT)=I
00619 ENDIF
00620 ENDDO
00621 IF(NUMIT.GT.NTRY) THEN
00622 NTRY=NTRY+1
00623 IP=ISTACK(NTRY)
00624 GOTO 30
00625 ENDIF
00626 ENDIF
00627
00628 DO 25 KK=0,NUMIT
00629 NA=NHEP
00630 FIRST=JDAHEP(1,ISTACK(KK))
00631 LAST=JDAHEP(2,ISTACK(KK))
00632 DO II=1,LAST-FIRST+1
00633 DO LL=1,5
00634 PORIG(LL,II)=PHEP(LL,FIRST+II-1)
00635 ENDDO
00636 ENDDO
00637
00638 CALL PHTYPE(ISTACK(KK))
00639
00640
00641 IF(NHEP.GT.NA) THEN
00642 DO II=1,LAST-FIRST+1
00643 IPP=FIRST+II-1
00644 FIRSTA=JDAHEP(1,IPP)
00645 LASTA=JDAHEP(2,IPP)
00646 IF(JMOHEP(1,IPP).EQ.ISTACK(KK))
00647 $ CALL PHOBOS(IPP,PORIG(1,II),PHEP(1,IPP),FIRSTA,LASTA)
00648 ENDDO
00649 ENDIF
00650 25 CONTINUE
00651
00652
00653 IF (NHEP.GT.NLAST) THEN
00654 DO 160 I=NLAST+1,NHEP
00655
00656
00657 MOTHER=JMOHEP(1,I)
00658 POSPHO=JDAHEP(2,MOTHER)+1
00659
00660 DO 90 J=1,5
00661 90 PHOTON(J)=PHEP(J,I)
00662 ISPHO =ISTHEP(I)
00663 IDPHO =IDHEP(I)
00664 MOTHER2 =JMOHEP(2,I)
00665 IDA1 =JDAHEP(1,I)
00666 IDA2 =JDAHEP(2,I)
00667
00668
00669 IF (POSPHO.NE.NHEP) THEN
00670
00671
00672
00673 DO 120 K=I,POSPHO+1,-1
00674 ISTHEP(K)=ISTHEP(K-1)
00675 QEDRAD(K)=QEDRAD(K-1)
00676 IDHEP(K)=IDHEP(K-1)
00677 DO 100 L=1,2
00678 JMOHEP(L,K)=JMOHEP(L,K-1)
00679 100 JDAHEP(L,K)=JDAHEP(L,K-1)
00680 DO 110 L=1,5
00681 110 PHEP(L,K)=PHEP(L,K-1)
00682 DO 120 L=1,4
00683 120 VHEP(L,K)=VHEP(L,K-1)
00684
00685
00686 DO 130 K=1,NHEP
00687 DO 130 L=1,2
00688 IF ((JMOHEP(L,K).NE.0).AND.(JMOHEP(L,K).GE.
00689 & POSPHO)) JMOHEP(L,K)=JMOHEP(L,K)+1
00690 IF ((JDAHEP(L,K).NE.0).AND.(JDAHEP(L,K).GE.
00691 & POSPHO)) JDAHEP(L,K)=JDAHEP(L,K)+1
00692 130 CONTINUE
00693
00694
00695 DO 140 J=1,5
00696 140 PHEP(J,POSPHO)=PHOTON(J)
00697 ENDIF
00698
00699
00700 JDAHEP(2,MOTHER)=POSPHO
00701 ISTHEP(POSPHO)=ISPHO
00702 IDHEP(POSPHO)=IDPHO
00703 JMOHEP(1,POSPHO)=MOTHER
00704 JMOHEP(2,POSPHO)=MOTHER2
00705 JDAHEP(1,POSPHO)=IDA1
00706 JDAHEP(2,POSPHO)=IDA2
00707
00708
00709 DO 150 J=1,4
00710 150 VHEP(J,POSPHO)=VHEP(J,POSPHO-1)
00711 160 CONTINUE
00712 ENDIF
00713 RETURN
00714 END
00715 SUBROUTINE PHOBOS(IP,PBOOS1,PBOOS2,FIRST,LAST)
00716
00717
00718
00719
00720
00721
00722
00723
00724
00725
00726
00727
00728
00729
00730
00731
00732
00733
00734
00735
00736
00737
00738 IMPLICIT NONE
00739 DOUBLE PRECISION BET1(3),BET2(3),GAM1,GAM2,PB,DATA
00740 INTEGER I,J,FIRST,LAST,MAXSTA,NSTACK,IP
00741 PARAMETER (MAXSTA=10000)
00742 INTEGER STACK(MAXSTA)
00743 REAL*8 PBOOS1(5),PBOOS2(5)
00744 INTEGER NMXHEP
00745 PARAMETER (NMXHEP=10000)
00746 INTEGER IDHEP,ISTHEP,JDAHEP,JMOHEP,NEVHEP,NHEP
00747 REAL*8 PHEP,VHEP
00748 COMMON/PH_HEPEVT/NEVHEP,NHEP,ISTHEP(NMXHEP),IDHEP(NMXHEP),
00749 &JMOHEP(2,NMXHEP),JDAHEP(2,NMXHEP),PHEP(5,NMXHEP),VHEP(4,NMXHEP)
00750 IF ((LAST.EQ.0).OR.(LAST.LT.FIRST)) RETURN
00751 NSTACK=0
00752 DO 10 J=1,3
00753 BET1(J)=-PBOOS1(J)/PBOOS1(5)
00754 10 BET2(J)=PBOOS2(J)/PBOOS2(5)
00755 GAM1=PBOOS1(4)/PBOOS1(5)
00756 GAM2=PBOOS2(4)/PBOOS2(5)
00757
00758
00759 20 DO 50 I=FIRST,LAST
00760 PB=BET1(1)*PHEP(1,I)+BET1(2)*PHEP(2,I)+BET1(3)*PHEP(3,I)
00761 IF (JMOHEP(1,I).EQ.IP) THEN
00762 DO 30 J=1,3
00763 30 PHEP(J,I)=PHEP(J,I)+BET1(J)*(PHEP(4,I)+PB/(GAM1+1.D0))
00764 PHEP(4,I)=GAM1*PHEP(4,I)+PB
00765
00766
00767 PB=BET2(1)*PHEP(1,I)+BET2(2)*PHEP(2,I)+BET2(3)*PHEP(3,I)
00768 DO 40 J=1,3
00769 40 PHEP(J,I)=PHEP(J,I)+BET2(J)*(PHEP(4,I)+PB/(GAM2+1.D0))
00770 PHEP(4,I)=GAM2*PHEP(4,I)+PB
00771 IF (JDAHEP(1,I).NE.0) THEN
00772 NSTACK=NSTACK+1
00773
00774
00775 IF (NSTACK.GT.MAXSTA) THEN
00776 DATA=NSTACK
00777 CALL PHOERR(7,'PHOBOS',DATA)
00778 ENDIF
00779 STACK(NSTACK)=I
00780 ENDIF
00781 ENDIF
00782 50 CONTINUE
00783 IF (NSTACK.NE.0) THEN
00784
00785
00786 FIRST=JDAHEP(1,STACK(NSTACK))
00787 LAST=JDAHEP(2,STACK(NSTACK))
00788 IP=STACK(NSTACK)
00789 NSTACK=NSTACK-1
00790 GOTO 20
00791 ENDIF
00792 RETURN
00793 END
00794 SUBROUTINE PHOIN(IP,BOOST,NHEP0)
00795
00796
00797
00798
00799
00800
00801
00802
00803
00804
00805
00806
00807
00808
00809
00810
00811
00812
00813 IMPLICIT NONE
00814 INTEGER NMXHEP
00815 PARAMETER (NMXHEP=10000)
00816 INTEGER IDHEP,ISTHEP,JDAHEP,JMOHEP,NEVHEP,NHEP
00817 REAL*8 PHEP,VHEP
00818 COMMON/PH_HEPEVT/NEVHEP,NHEP,ISTHEP(NMXHEP),IDHEP(NMXHEP),
00819 &JMOHEP(2,NMXHEP),JDAHEP(2,NMXHEP),PHEP(5,NMXHEP),VHEP(4,NMXHEP)
00820 INTEGER NMXPHO
00821 PARAMETER (NMXPHO=10000)
00822 INTEGER IDPHO,ISTPHO,JDAPHO,JMOPHO,NEVPHO,NPHO
00823 REAL*8 PPHO,VPHO
00824 COMMON/PHOEVT/NEVPHO,NPHO,ISTPHO(NMXPHO),IDPHO(NMXPHO),
00825 &JMOPHO(2,NMXPHO),JDAPHO(2,NMXPHO),PPHO(5,NMXPHO),VPHO(4,NMXPHO)
00826 INTEGER IP,IP2,I,FIRST,LAST,LL,NA
00827 LOGICAL BOOST
00828 INTEGER J,NHEP0
00829 DOUBLE PRECISION BET(3),GAM,PB
00830 COMMON /PHOCMS/ BET,GAM
00831 LOGICAL INTERF,ISEC,ITRE,IEXP,IFTOP,IFW
00832 REAL*8 FINT,FSEC,EXPEPS
00833 COMMON /PHOKEY/ FSEC,FINT,EXPEPS,INTERF,ISEC,ITRE,IEXP,IFTOP,IFW
00834
00835
00836 FIRST=JDAHEP(1,IP)
00837 LAST =JDAHEP(2,IP)
00838 NPHO=3+LAST-FIRST+NHEP-NHEP0
00839 NEVPHO=NPHO
00840
00841 IDPHO(1)=IDHEP(IP)
00842 JDAPHO(1,1)=3
00843 JDAPHO(2,1)=3+LAST-FIRST
00844 DO I=1,5
00845 PPHO(I,1)=PHEP(I,IP)
00846 ENDDO
00847
00848 IP2=JMOHEP(2,JDAHEP(1,IP))
00849 IF((IP2.NE.0).AND.(IP2.NE.IP)) THEN
00850 IDPHO(2)=IDHEP(IP2)
00851 JDAPHO(1,2)=3
00852 JDAPHO(2,2)=3+LAST-FIRST
00853 DO I=1,5
00854 PPHO(I,2)=PHEP(I,IP2)
00855 ENDDO
00856 ELSE
00857 IDPHO(2)=0
00858 DO I=1,5
00859 PPHO(I,2)=0.0D0
00860 ENDDO
00861 ENDIF
00862
00863 DO LL=0,LAST-FIRST
00864 IDPHO(3+LL)=IDHEP(FIRST+LL)
00865 JMOPHO(1,3+LL)=JMOHEP(1,FIRST+LL)
00866 IF (JMOHEP(1,FIRST+LL).EQ.IP) JMOPHO(1,3+LL)=1
00867 DO I=1,5
00868 PPHO(I,3+LL)=PHEP(I,FIRST+LL)
00869 ENDDO
00870 ENDDO
00871 IF (NHEP.GT.NHEP0) THEN
00872
00873 NA=3+LAST-FIRST
00874 DO LL=1,NHEP-NHEP0
00875 IDPHO(NA+LL)=IDHEP(NHEP0+LL)
00876 JMOPHO(1,NA+LL)=JMOHEP(1,NHEP0+LL)
00877 IF (JMOHEP(1,NHEP0+LL).EQ.IP) JMOPHO(1,NA+LL)=1
00878 DO I=1,5
00879 PPHO(I,NA+LL)=PHEP(I,NHEP0+LL)
00880 ENDDO
00881 ENDDO
00882
00883 JDAPHO(2,1)=3+LAST-FIRST+NHEP-NHEP0
00884 ENDIF
00885 IF(IDPHO(NPHO).EQ.22)CALL PHLUPA(100001)
00886
00887 CALL PHCORK(0)
00888 IF(IDPHO(NPHO).EQ.22)CALL PHLUPA(100002)
00889
00890 IF(IFTOP) CALL PHOTWO(0)
00891 BOOST=.FALSE.
00892
00893
00894
00895
00896
00897
00898
00899
00900
00901
00902
00903
00904 IF ((ABS(PPHO(1,1)+ABS(PPHO(2,1))+ABS(PPHO(3,1))).GT.
00905 $ PPHO(5,1)*1.D-8).AND.(PPHO(5,1).NE.0)) THEN
00906
00907 BOOST=.TRUE.
00908
00909
00910
00911 DO 10 J=1,3
00912 10 BET(J)=-PPHO(J,1)/PPHO(5,1)
00913 GAM=PPHO(4,1)/PPHO(5,1)
00914 DO 30 I=JDAPHO(1,1),JDAPHO(2,1)
00915 PB=BET(1)*PPHO(1,I)+BET(2)*PPHO(2,I)+BET(3)*PPHO(3,I)
00916 DO 20 J=1,3
00917 20 PPHO(J,I)=PPHO(J,I)+BET(J)*(PPHO(4,I)+PB/(GAM+1.D0))
00918 30 PPHO(4,I)=GAM*PPHO(4,I)+PB
00919
00920 I=1
00921 PB=BET(1)*PPHO(1,I)+BET(2)*PPHO(2,I)+BET(3)*PPHO(3,I)
00922 DO J=1,3
00923 PPHO(J,I)=PPHO(J,I)+BET(J)*(PPHO(4,I)+PB/(GAM+1.D0))
00924 ENDDO
00925 PPHO(4,I)=GAM*PPHO(4,I)+PB
00926 ENDIF
00927
00928 IF(IFTOP) CALL PHOTWO(1)
00929 CALL PHLUPA(2)
00930 IF(IDPHO(NPHO).EQ.22) CALL PHLUPA(10000)
00931
00932 END
00933 SUBROUTINE PHOTWO(MODE)
00934
00935
00936
00937
00938
00939
00940
00941
00942
00943
00944
00945
00946
00947
00948
00949 IMPLICIT NONE
00950 INTEGER NMXPHO
00951 PARAMETER (NMXPHO=10000)
00952 INTEGER IDPHO,ISTPHO,JDAPHO,JMOPHO,NEVPHO,NPHO
00953 REAL*8 PPHO,VPHO
00954 COMMON/PHOEVT/NEVPHO,NPHO,ISTPHO(NMXPHO),IDPHO(NMXPHO),
00955 &JMOPHO(2,NMXPHO),JDAPHO(2,NMXPHO),PPHO(5,NMXPHO),VPHO(4,NMXPHO)
00956 DOUBLE PRECISION BET(3),GAM
00957 COMMON /PHOCMS/ BET,GAM
00958 INTEGER I,MODE
00959 REAL*8 MPASQR
00960 LOGICAL IFRAD
00961
00962
00963
00964
00965
00966
00967 IF(MODE.EQ.0) THEN
00968 IFRAD=(IDPHO(1).EQ.21).AND.(IDPHO(2).EQ.21)
00969 IFRAD=IFRAD.OR.(IDPHO(1).EQ.-IDPHO(2).AND.ABS(IDPHO(1)).LE.6)
00970 IFRAD=IFRAD
00971 & .AND.(ABS(IDPHO(3)).EQ.6).AND.(ABS(IDPHO(4)).EQ.6)
00972 MPASQR= (PPHO(4,1)+PPHO(4,2))**2-(PPHO(3,1)+PPHO(3,2))**2
00973 & -(PPHO(2,1)+PPHO(2,2))**2-(PPHO(1,1)+PPHO(1,2))**2
00974 IFRAD=IFRAD.AND.(MPASQR.GT.0.0D0)
00975 IF(IFRAD) THEN
00976
00977 DO I=1,4
00978 PPHO(I,1)=PPHO(I,1)+PPHO(I,2)
00979 ENDDO
00980 PPHO(5,1)=SQRT(MPASQR)
00981
00982 DO I=1,5
00983 PPHO(I,2)=0.0D0
00984 ENDDO
00985 ENDIF
00986 ELSE
00987
00988
00989
00990 ENDIF
00991 END
00992 SUBROUTINE PHOOUT(IP,BOOST,NHEP0)
00993
00994
00995
00996
00997
00998
00999
01000
01001
01002
01003
01004
01005
01006
01007
01008
01009
01010
01011 IMPLICIT NONE
01012 INTEGER NMXHEP
01013 PARAMETER (NMXHEP=10000)
01014 INTEGER IDHEP,ISTHEP,JDAHEP,JMOHEP,NEVHEP,NHEP
01015 REAL*8 PHEP,VHEP
01016 COMMON/PH_HEPEVT/NEVHEP,NHEP,ISTHEP(NMXHEP),IDHEP(NMXHEP),
01017 &JMOHEP(2,NMXHEP),JDAHEP(2,NMXHEP),PHEP(5,NMXHEP),VHEP(4,NMXHEP)
01018 INTEGER NMXPHO
01019 PARAMETER (NMXPHO=10000)
01020 INTEGER IDPHO,ISTPHO,JDAPHO,JMOPHO,NEVPHO,NPHO
01021 REAL*8 PPHO,VPHO
01022 COMMON/PHOEVT/NEVPHO,NPHO,ISTPHO(NMXPHO),IDPHO(NMXPHO),
01023 &JMOPHO(2,NMXPHO),JDAPHO(2,NMXPHO),PPHO(5,NMXPHO),VPHO(4,NMXPHO)
01024 INTEGER IP,LL,FIRST,LAST,I
01025 LOGICAL BOOST
01026 INTEGER NN,J,K,NHEP0,NA
01027 DOUBLE PRECISION BET(3),GAM,PB
01028 COMMON /PHOCMS/ BET,GAM
01029 IF(NPHO.EQ.NEVPHO) RETURN
01030
01031 CALL PHLUPA(10)
01032 IF (BOOST) THEN
01033 DO 110 J=JDAPHO(1,1),JDAPHO(2,1)
01034 PB=-BET(1)*PPHO(1,J)-BET(2)*PPHO(2,J)-BET(3)*PPHO(3,J)
01035 DO 100 K=1,3
01036 100 PPHO(K,J)=PPHO(K,J)-BET(K)*(PPHO(4,J)+PB/(GAM+1.D0))
01037 110 PPHO(4,J)=GAM*PPHO(4,J)+PB
01038
01039 DO NN=NEVPHO+1,NPHO
01040 PB=-BET(1)*PPHO(1,NN)-BET(2)*PPHO(2,NN)-BET(3)*PPHO(3,NN)
01041 DO 120 K=1,3
01042 120 PPHO(K,NN)=PPHO(K,NN)-BET(K)*(PPHO(4,NN)+PB/(GAM+1.D0))
01043 PPHO(4,NN)=GAM*PPHO(4,NN)+PB
01044 ENDDO
01045 ENDIF
01046 FIRST=JDAHEP(1,IP)
01047 LAST =JDAHEP(2,IP)
01048
01049 DO LL=0,LAST-FIRST
01050 IDHEP(FIRST+LL) = IDPHO(3+LL)
01051 DO I=1,5
01052 PHEP(I,FIRST+LL) = PPHO(I,3+LL)
01053 ENDDO
01054 ENDDO
01055
01056 NA=3+LAST-FIRST
01057 DO LL=1,NPHO-NA
01058 IDHEP(NHEP0+LL) = IDPHO(NA+LL)
01059 ISTHEP(NHEP0+LL)=ISTPHO(NA+LL)
01060 JMOHEP(1,NHEP0+LL)=IP
01061 JMOHEP(2,NHEP0+LL)=JMOHEP(2,JDAHEP(1,IP))
01062 JDAHEP(1,NHEP0+LL)=0
01063 JDAHEP(2,NHEP0+LL)=0
01064 DO I=1,5
01065 PHEP(I,NHEP0+LL) = PPHO(I,NA+LL)
01066 ENDDO
01067 ENDDO
01068 NHEP=NHEP+NPHO-NEVPHO
01069 CALL PHLUPA(20)
01070 END
01071 SUBROUTINE PHOCHK(JFIRST)
01072
01073
01074
01075
01076
01077
01078
01079
01080
01081
01082
01083
01084
01085
01086
01087 IMPLICIT NONE
01088 INTEGER NMXPHO
01089 PARAMETER (NMXPHO=10000)
01090 INTEGER IDPHO,ISTPHO,JDAPHO,JMOPHO,NEVPHO,NPHO
01091 REAL*8 PPHO,VPHO
01092 COMMON/PHOEVT/NEVPHO,NPHO,ISTPHO(NMXPHO),IDPHO(NMXPHO),
01093 &JMOPHO(2,NMXPHO),JDAPHO(2,NMXPHO),PPHO(5,NMXPHO),VPHO(4,NMXPHO)
01094 LOGICAL CHKIF
01095 COMMON/PHOIF/CHKIF(NMXPHO)
01096 INTEGER NMXHEP
01097 PARAMETER (NMXHEP=10000)
01098 LOGICAL QEDRAD
01099 COMMON/PH_PHOQED/QEDRAD(NMXHEP)
01100 INTEGER JFIRST
01101 LOGICAL F
01102 INTEGER IDABS,NLAST,I,IPPAR
01103 LOGICAL INTERF,ISEC,ITRE,IEXP,IFTOP,IFW,IFNPI0,IFKL
01104 REAL*8 FINT,FSEC,EXPEPS
01105 COMMON /PHOKEY/ FSEC,FINT,EXPEPS,INTERF,ISEC,ITRE,IEXP,IFTOP,IFW
01106 LOGICAL IFRAD
01107 INTEGER IDENT,K,IQRK,IPHQRK,IEKL,IPHEKL
01108
01109 F(IDABS)=
01110 & ( ((IDABS.GT.9.OR.IQRK.NE.1).AND.(IDABS.LE.40))
01111 & .OR.(IDABS.GT.100) )
01112 & .AND.(IDABS.NE.21)
01113 $ .AND.(IDABS.NE.2101).AND.(IDABS.NE.3101).AND.(IDABS.NE.3201)
01114 & .AND.(IDABS.NE.1103).AND.(IDABS.NE.2103).AND.(IDABS.NE.2203)
01115 & .AND.(IDABS.NE.3103).AND.(IDABS.NE.3203).AND.(IDABS.NE.3303)
01116
01117 IQRK=IPHQRK(0)
01118 IEKL=IPHEKL(0)
01119 NLAST = NPHO
01120
01121 IPPAR=1
01122
01123 IFNPI0=.TRUE.
01124 IF (IEKL.GT.1) THEN
01125
01126 IFNPI0= (IDPHO(1).NE.111)
01127 IFKL = ((IDPHO(1).EQ.130).AND.
01128 $ ((IDPHO(3).EQ.22).OR.(IDPHO(4).EQ.22).OR.
01129 $ (IDPHO(5).EQ.22)).AND.
01130 $ ((IDPHO(3).EQ.11).OR.(IDPHO(4).EQ.11).OR.
01131 $ (IDPHO(5).EQ.11)) )
01132
01133 IFNPI0=(IFNPI0.AND.(.NOT.IFKL))
01134 ENDIF
01135 DO 10 I=IPPAR,NLAST
01136 IDABS = ABS(IDPHO(I))
01137
01138 CHKIF(I)= F(IDABS) .AND.F(ABS(IDPHO(1)))
01139 & .AND. (IDPHO(2).EQ.0)
01140 IF(I.GT.2) CHKIF(I)=CHKIF(I).AND.QEDRAD(JFIRST+I-IPPAR-2)
01141 & .AND.IFNPI0
01142 10 CONTINUE
01143
01144
01145
01146 IF(IFTOP) THEN
01147
01148 DO K=JDAPHO(2,1),JDAPHO(1,1),-1
01149 IF(IDPHO(K).NE.22) THEN
01150 IDENT=K
01151 GOTO 15
01152 ENDIF
01153 ENDDO
01154 15 CONTINUE
01155 IFRAD=((IDPHO(1).EQ.21).AND.(IDPHO(2).EQ.21))
01156 & .OR. ((ABS(IDPHO(1)).LE.6).AND.((IDPHO(2)).EQ.(-IDPHO(1))))
01157 IFRAD=IFRAD
01158 & .AND.(ABS(IDPHO(3)).EQ.6).AND.((IDPHO(4)).EQ.(-IDPHO(3)))
01159 & .AND.(IDENT.EQ.4)
01160 IF(IFRAD) THEN
01161 DO 20 I=IPPAR,NLAST
01162 CHKIF(I)= .TRUE.
01163 IF(I.GT.2) CHKIF(I)=CHKIF(I).AND.QEDRAD(JFIRST+I-IPPAR-2)
01164 20 CONTINUE
01165 ENDIF
01166 ENDIF
01167
01168
01169 IF(IFTOP) THEN
01170
01171 DO K=JDAPHO(2,1),JDAPHO(1,1),-1
01172 IF(IDPHO(K).NE.22) THEN
01173 IDENT=K
01174 GOTO 25
01175 ENDIF
01176 ENDDO
01177 25 CONTINUE
01178 IFRAD=((ABS(IDPHO(1)).EQ.6).AND.(IDPHO(2).EQ.0))
01179 IFRAD=IFRAD
01180 & .AND.((ABS(IDPHO(3)).EQ.24).AND.(ABS(IDPHO(4)).EQ.5)
01181 & .OR.(ABS(IDPHO(3)).EQ.5).AND.(ABS(IDPHO(4)).EQ.24))
01182 & .AND.(IDENT.EQ.4)
01183 IF(IFRAD) THEN
01184 DO 30 I=IPPAR,NLAST
01185 CHKIF(I)= .TRUE.
01186 IF(I.GT.2) CHKIF(I)=CHKIF(I).AND.QEDRAD(JFIRST+I-IPPAR-2)
01187 30 CONTINUE
01188 ENDIF
01189 ENDIF
01190
01191
01192 END
01193 SUBROUTINE PHTYPE(ID)
01194
01195
01196
01197
01198
01199
01200
01201
01202
01203
01204
01205
01206
01207
01208
01209
01210 IMPLICIT NONE
01211 INTEGER NMXHEP
01212 PARAMETER (NMXHEP=10000)
01213 INTEGER IDHEP,ISTHEP,JDAHEP,JMOHEP,NEVHEP,NHEP
01214 REAL*8 PHEP,VHEP
01215 COMMON/PH_HEPEVT/NEVHEP,NHEP,ISTHEP(NMXHEP),IDHEP(NMXHEP),
01216 &JMOHEP(2,NMXHEP),JDAHEP(2,NMXHEP),PHEP(5,NMXHEP),VHEP(4,NMXHEP)
01217 LOGICAL INTERF,ISEC,ITRE,IEXP,IFTOP,IFW
01218 REAL*8 FINT,FSEC,EXPEPS
01219 COMMON /PHOKEY/ FSEC,FINT,EXPEPS,INTERF,ISEC,ITRE,IEXP,IFTOP,IFW
01220 LOGICAL EXPINI
01221 INTEGER NX,K,NCHAN
01222 PARAMETER (NX=10)
01223 REAL*8 PRO,PRSUM,ESU
01224 COMMON /PHOEXP/ PRO(NX),NCHAN,EXPINI
01225
01226 INTEGER ID,NHEP0
01227 LOGICAL IPAIR
01228 REAL*8 RN,PHORAN,SUM
01229 INTEGER WTDUM
01230 LOGICAL IFOUR
01231
01232 IFOUR=(.TRUE.).AND.(ITRE)
01233
01234 IPAIR=.TRUE.
01235
01236 IF (JDAHEP(1,ID).EQ.0) RETURN
01237
01238
01239 NHEP0=NHEP
01240
01241 IF (IEXP) THEN
01242 EXPINI=.TRUE.
01243 DO NCHAN=1,NX
01244 PRO(NCHAN)=0.D0
01245 ENDDO
01246 NCHAN=0
01247
01248 FSEC=1.0D0
01249 CALL PHOMAK(ID,NHEP0)
01250
01251 EXPINI=.FALSE.
01252 RN=PHORAN(WTDUM)
01253 PRSUM=0
01254 DO K=1,NX
01255 PRSUM=PRSUM+PRO(K)
01256 ENDDO
01257 ESU=EXP(-PRSUM)
01258
01259
01260 SUM=ESU
01261
01262 DO K=1,100
01263 IF(RN.LT.SUM) GOTO 100
01264 ESU=ESU*PRSUM/K
01265 SUM=SUM+ESU
01266 NCHAN=0
01267 CALL PHOMAK(ID,NHEP0)
01268 IF(SUM.GT.1D0-EXPEPS) GOTO 100
01269 ENDDO
01270 100 CONTINUE
01271 ELSEIF(IFOUR) THEN
01272
01273 FSEC=1.0D0
01274 RN=PHORAN(WTDUM)
01275 IF (RN.GE.23.D0/24D0) THEN
01276 CALL PHOMAK(ID,NHEP0)
01277 CALL PHOMAK(ID,NHEP0)
01278 CALL PHOMAK(ID,NHEP0)
01279 CALL PHOMAK(ID,NHEP0)
01280 ELSEIF (RN.GE.17.D0/24D0) THEN
01281 CALL PHOMAK(ID,NHEP0)
01282 CALL PHOMAK(ID,NHEP0)
01283 ELSEIF (RN.GE.9.D0/24D0) THEN
01284 CALL PHOMAK(ID,NHEP0)
01285 ENDIF
01286 ELSEIF(ITRE) THEN
01287
01288 FSEC=1.0D0
01289 RN=PHORAN(WTDUM)
01290 IF (RN.GE.5.D0/6D0) THEN
01291 CALL PHOMAK(ID,NHEP0)
01292 CALL PHOMAK(ID,NHEP0)
01293 CALL PHOMAK(ID,NHEP0)
01294 ELSEIF (RN.GE.2.D0/6D0) THEN
01295 CALL PHOMAK(ID,NHEP0)
01296 ENDIF
01297 ELSEIF(ISEC) THEN
01298
01299 FSEC=1.0D0
01300 RN=PHORAN(WTDUM)
01301 IF (RN.GE.0.5D0) THEN
01302 CALL PHOMAK(ID,NHEP0)
01303 CALL PHOMAK(ID,NHEP0)
01304 ENDIF
01305 ELSE
01306
01307 FSEC=1.0D0
01308 CALL PHOMAK(ID,NHEP0)
01309 ENDIF
01310
01311
01312
01313 END
01314 SUBROUTINE PHOMAK(IPPAR,NHEP0)
01315
01316
01317
01318
01319
01320
01321
01322
01323
01324
01325
01326
01327
01328
01329
01330
01331
01332
01333
01334 IMPLICIT NONE
01335 DOUBLE PRECISION DATA
01336 REAL*8 PHORAN
01337 INTEGER IP,IPPAR,NCHARG
01338 INTEGER WTDUM,IDUM,NHEP0
01339 INTEGER NCHARB,NEUDAU
01340 REAL*8 RN,WT,PHINT
01341 LOGICAL BOOST
01342 INTEGER NMXHEP
01343 PARAMETER (NMXHEP=10000)
01344 INTEGER IDHEP,ISTHEP,JDAHEP,JMOHEP,NEVHEP,NHEP
01345 REAL*8 PHEP,VHEP
01346 COMMON/PH_HEPEVT/NEVHEP,NHEP,ISTHEP(NMXHEP),IDHEP(NMXHEP),
01347 &JMOHEP(2,NMXHEP),JDAHEP(2,NMXHEP),PHEP(5,NMXHEP),VHEP(4,NMXHEP)
01348 LOGICAL INTERF,ISEC,ITRE,IEXP,IFTOP,IFW
01349 REAL*8 FINT,FSEC,EXPEPS
01350 COMMON /PHOKEY/ FSEC,FINT,EXPEPS,INTERF,ISEC,ITRE,IEXP,IFTOP,IFW
01351
01352 IP=IPPAR
01353 IDUM=1
01354 NCHARG=0
01355
01356 CALL PHOIN(IP,BOOST,NHEP0)
01357 CALL PHOCHK(JDAHEP(1,IP))
01358 WT=0.0D0
01359 CALL PHOPRE(1,WT,NEUDAU,NCHARB)
01360
01361 IF (WT.EQ.0.0D0) RETURN
01362 RN=PHORAN(WTDUM)
01363
01364 CALL PHODO(1,NCHARB,NEUDAU)
01365
01366 IF (INTERF) WT=WT*PHINT(IDUM) /FINT
01367 IF (IFW) CALL PHOBW(WT)
01368 DATA=WT
01369 IF (WT.GT.1.0D0) CALL PHOERR(3,'WT_INT',DATA)
01370
01371 IF (RN.LE.WT) THEN
01372 CALL PHOOUT(IP,BOOST,NHEP0)
01373 ENDIF
01374 RETURN
01375 END
01376 FUNCTION PHINT1(IDUM)
01377
01378
01379
01380
01381
01382
01383
01384
01385
01386
01387
01388
01389
01390
01391
01392
01393
01394
01395
01396 IMPLICIT NONE
01397 REAL*8 PHINT,phint1
01398 REAL*8 PHOCHA
01399 INTEGER IDUM
01400 INTEGER NMXPHO
01401 PARAMETER (NMXPHO=10000)
01402 INTEGER IDPHO,ISTPHO,JDAPHO,JMOPHO,NEVPHO,NPHO
01403 REAL*8 PPHO,VPHO
01404 COMMON/PHOEVT/NEVPHO,NPHO,ISTPHO(NMXPHO),IDPHO(NMXPHO),
01405 &JMOPHO(2,NMXPHO),JDAPHO(2,NMXPHO),PPHO(5,NMXPHO),VPHO(4,NMXPHO)
01406 DOUBLE PRECISION MCHSQR,MNESQR
01407 REAL*8 PNEUTR
01408 COMMON/PHOMOM/MCHSQR,MNESQR,PNEUTR(5)
01409 DOUBLE PRECISION COSTHG,SINTHG
01410 REAL*8 XPHMAX,XPHOTO
01411 COMMON/PHOPHS/XPHMAX,XPHOTO,COSTHG,SINTHG
01412 REAL*8 MPASQR,XX,BETA
01413 LOGICAL IFINT
01414 INTEGER K,IDENT
01415
01416 DO K=JDAPHO(2,1),JDAPHO(1,1),-1
01417 IF(IDPHO(K).NE.22) THEN
01418 IDENT=K
01419 GOTO 20
01420 ENDIF
01421 ENDDO
01422 20 CONTINUE
01423
01424 IFINT= NPHO.GT.IDENT
01425
01426 IFINT= IFINT.AND.(IDENT-JDAPHO(1,1)).EQ.1
01427
01428 IFINT= IFINT.AND.IDPHO(JDAPHO(1,1)).EQ.-IDPHO(IDENT)
01429
01430 IFINT= IFINT.AND.PHOCHA(IDPHO(IDENT)).NE.0
01431
01432 IF(IFINT) THEN
01433 MPASQR = PPHO(5,1)**2
01434 XX=4.D0*MCHSQR/MPASQR*(1.D0-XPHOTO)/(1.D0-XPHOTO+(MCHSQR-MNESQR)
01435 & /MPASQR)**2
01436 BETA=SQRT(1.D0-XX)
01437 PHINT = 2D0/(1D0+COSTHG**2*BETA**2)
01438 ELSE
01439 PHINT = 1D0
01440 ENDIF
01441 phint1=1
01442 END
01443
01444 FUNCTION PHINT2(IDUM)
01445
01446
01447
01448
01449
01450
01451
01452
01453
01454
01455
01456
01457
01458
01459
01460
01461
01462
01463
01464 IMPLICIT NONE
01465 REAL*8 PHINT,PHINT1,PHINT2
01466 REAL*8 PHOCHA
01467 INTEGER IDUM
01468 INTEGER NMXPHO
01469 PARAMETER (NMXPHO=10000)
01470 INTEGER IDPHO,ISTPHO,JDAPHO,JMOPHO,NEVPHO,NPHO
01471 REAL*8 PPHO,VPHO
01472 COMMON/PHOEVT/NEVPHO,NPHO,ISTPHO(NMXPHO),IDPHO(NMXPHO),
01473 &JMOPHO(2,NMXPHO),JDAPHO(2,NMXPHO),PPHO(5,NMXPHO),VPHO(4,NMXPHO)
01474 DOUBLE PRECISION MCHSQR,MNESQR
01475 REAL*8 PNEUTR
01476 COMMON/PHOMOM/MCHSQR,MNESQR,PNEUTR(5)
01477 DOUBLE PRECISION COSTHG,SINTHG
01478 REAL*8 XPHMAX,XPHOTO
01479 COMMON/PHOPHS/XPHMAX,XPHOTO,COSTHG,SINTHG
01480 REAL*8 MPASQR,XX,BETA,PQ1(4),PQ2(4),PPHOT(4)
01481 REAL*8 SS,PP2,PP,E1,E2,Q1,Q2,COSTHE
01482 LOGICAL IFINT
01483 INTEGER K,IDENT
01484
01485 DO K=JDAPHO(2,1),JDAPHO(1,1),-1
01486 IF(IDPHO(K).NE.22) THEN
01487 IDENT=K
01488 GOTO 20
01489 ENDIF
01490 ENDDO
01491 20 CONTINUE
01492
01493 IFINT= NPHO.GT.IDENT
01494
01495 IFINT= IFINT.AND.(IDENT-JDAPHO(1,1)).EQ.1
01496
01497
01498
01499 IFINT= IFINT.AND.abs(PHOCHA(IDPHO(IDENT))).GT.0.01D0
01500
01501 IFINT= IFINT.AND.
01502 $ abs(PHOCHA(IDPHO(JDAPHO(1,1)))).gt.0.01D0
01503
01504 IF(IFINT) THEN
01505 MPASQR = PPHO(5,1)**2
01506 XX=4.D0*MCHSQR/MPASQR*(1.-XPHOTO)/(1.-XPHOTO+(MCHSQR-MNESQR)/
01507 & MPASQR)**2
01508 BETA=SQRT(1.D0-XX)
01509 PHINT = 2D0/(1D0+COSTHG**2*BETA**2)
01510 SS =MPASQR*(1.D0-XPHOTO)
01511 PP2=((SS-MCHSQR-MNESQR)**2-4*MCHSQR*MNESQR)/SS/4
01512 PP =SQRT(PP2)
01513 E1 =SQRT(PP2+MCHSQR)
01514 E2 =SQRT(PP2+MNESQR)
01515 PHINT= (E1+E2)**2/((E2+COSTHG*PP)**2+(E1-COSTHG*PP)**2)
01516
01517 q1=PHOCHA(IDPHO(JDAPHO(1,1)))
01518 q2=PHOCHA(IDPHO(IDENT))
01519 do k=1,4
01520 pq1(k)=ppho(k,JDAPHO(1,1))
01521 pq2(k)=ppho(k,JDAPHO(1,1)+1)
01522 pphot(k)=ppho(k,npho)
01523 enddo
01524 costhe=(pphot(1)*pq1(1)+pphot(2)*pq1(2)+pphot(3)*pq1(3))
01525 costhe=costhe/sqrt(pq1(1)**2+pq1(2)**2+pq1(3)**2)
01526 costhe=costhe/sqrt(pphot(1)**2+pphot(2)**2+pphot(3)**2)
01527
01528
01529
01530
01531
01532 IF (costhg*costhe.GT.0) then
01533
01534 PHINT= (q1*(E2+COSTHG*PP)-q2*(E1-COSTHG*PP))**2
01535 & /(q1**2*(E2+COSTHG*PP)**2+q2**2*(E1-COSTHG*PP)**2)
01536 ELSE
01537
01538 PHINT= (q1*(E1-COSTHG*PP)-q2*(E2+COSTHG*PP))**2
01539 & /(q1**2*(E1-COSTHG*PP)**2+q2**2*(E2+COSTHG*PP)**2)
01540 ENDIF
01541 ELSE
01542 PHINT = 1D0
01543 ENDIF
01544 phint1=1
01545 phint2=1
01546 END
01547
01548
01549 SUBROUTINE PHOPRE(IPARR,WT,NEUDAU,NCHARB)
01550
01551
01552
01553
01554
01555
01556
01557
01558
01559
01560
01561
01562
01563
01564
01565
01566
01567
01568
01569
01570
01571
01572 IMPLICIT NONE
01573 DOUBLE PRECISION MINMAS,MPASQR,MCHREN
01574 DOUBLE PRECISION BETA,EPS,DEL1,DEL2,DATA,BIGLOG
01575 REAL*8 PHOCHA,PHOSPI,PHORAN,PHOCOR,MASSUM
01576 INTEGER IP,IPARR,IPPAR,I,J,ME,NCHARG,NEUPOI,NLAST,THEDUM
01577 INTEGER IDABS,IDUM
01578 INTEGER NCHARB,NEUDAU
01579 REAL*8 WT,WGT
01580 INTEGER NMXPHO
01581 PARAMETER (NMXPHO=10000)
01582 INTEGER IDPHO,ISTPHO,JDAPHO,JMOPHO,NEVPHO,NPHO
01583 REAL*8 PPHO,VPHO
01584 COMMON/PHOEVT/NEVPHO,NPHO,ISTPHO(NMXPHO),IDPHO(NMXPHO),
01585 &JMOPHO(2,NMXPHO),JDAPHO(2,NMXPHO),PPHO(5,NMXPHO),VPHO(4,NMXPHO)
01586 LOGICAL CHKIF
01587 COMMON/PHOIF/CHKIF(NMXPHO)
01588 INTEGER CHAPOI(NMXPHO)
01589 DOUBLE PRECISION MCHSQR,MNESQR
01590 REAL*8 PNEUTR
01591 COMMON/PHOMOM/MCHSQR,MNESQR,PNEUTR(5)
01592 DOUBLE PRECISION COSTHG,SINTHG
01593 REAL*8 XPHMAX,XPHOTO
01594 COMMON/PHOPHS/XPHMAX,XPHOTO,COSTHG,SINTHG
01595 REAL*8 ALPHA,XPHCUT
01596 COMMON/PHOCOP/ALPHA,XPHCUT
01597 INTEGER IREP
01598 REAL*8 PROBH,CORWT,XF
01599 COMMON/PHOPRO/PROBH,CORWT,XF,IREP
01600
01601 LOGICAL INTERF,ISEC,ITRE,IEXP,IFTOP,IFW
01602 REAL*8 FINT,FSEC,EXPEPS
01603 COMMON /PHOKEY/ FSEC,FINT,EXPEPS,INTERF,ISEC,ITRE,IEXP,IFTOP,IFW
01604
01605
01606 IPPAR=IPARR
01607
01608 IP=IPPAR
01609 NLAST=NPHO
01610 IDUM=1
01611
01612
01613 IF (JDAPHO(1,IP).EQ.0) RETURN
01614
01615
01616 10 NCHARG=0
01617 IREP=0
01618 MINMAS=0.D0
01619 MASSUM=0.D0
01620 DO 20 I=JDAPHO(1,IP),JDAPHO(2,IP)
01621
01622
01623
01624 IDABS=ABS(IDPHO(I))
01625 IF (CHKIF(I-JDAPHO(1,IP)+3)) THEN
01626 IF (PHOCHA(IDPHO(I)).NE.0) THEN
01627 NCHARG=NCHARG+1
01628 IF (NCHARG.GT.NMXPHO) THEN
01629 DATA=NCHARG
01630 CALL PHOERR(1,'PHOTOS',DATA)
01631 ENDIF
01632 CHAPOI(NCHARG)=I
01633 ENDIF
01634 MINMAS=MINMAS+PPHO(5,I)**2
01635 ENDIF
01636 MASSUM=MASSUM+PPHO(5,I)
01637 20 CONTINUE
01638 IF (NCHARG.NE.0) THEN
01639
01640
01641 IF ((PPHO(5,IP)-MASSUM)/PPHO(5,IP).GT.2.D0*XPHCUT) THEN
01642
01643
01644
01645 IF (NCHARG.GT.1) CALL PHOOMA(1,NCHARG,CHAPOI)
01646
01647 30 CONTINUE
01648 DO 70 J=1,3
01649 70 PNEUTR(J)=-PPHO(J,CHAPOI(NCHARG))
01650 PNEUTR(4)=PPHO(5,IP)-PPHO(4,CHAPOI(NCHARG))
01651
01652
01653 MPASQR=PPHO(5,IP)**2
01654 MCHSQR=PPHO(5,CHAPOI(NCHARG))**2
01655 IF ((JDAPHO(2,IP)-JDAPHO(1,IP)).EQ.1) THEN
01656 NEUPOI=JDAPHO(1,IP)
01657 IF (NEUPOI.EQ.CHAPOI(NCHARG)) NEUPOI=JDAPHO(2,IP)
01658 MNESQR=PPHO(5,NEUPOI)**2
01659 PNEUTR(5)=PPHO(5,NEUPOI)
01660 ELSE
01661 MNESQR=PNEUTR(4)**2-PNEUTR(1)**2-PNEUTR(2)**2-PNEUTR(3)**2
01662 MNESQR=MAX(MNESQR,MINMAS-MCHSQR)
01663 PNEUTR(5)=SQRT(MNESQR)
01664 ENDIF
01665
01666
01667 XPHMAX=(MPASQR-(PNEUTR(5)+PPHO(5,CHAPOI(NCHARG)))**2)/MPASQR
01668
01669
01670 CALL PHOENE(MPASQR,MCHREN,BETA,BIGLOG,IDPHO(CHAPOI(NCHARG)))
01671
01672 IF (XPHOTO.LT.-4D0) THEN
01673 NCHARG=0
01674 XPHOTO=0d0
01675
01676 ELSEIF ((XPHOTO.LT.XPHCUT).OR.(XPHOTO.GT.XPHMAX)) THEN
01677
01678
01679
01680 NCHARG=NCHARG-1
01681 IF (NCHARG.GT.0) THEN
01682 IREP=IREP+1
01683 GOTO 30
01684 ENDIF
01685 ELSE
01686
01687
01688
01689 EPS=MCHREN/(1.D0+BETA)
01690
01691
01692 DEL1=(2.D0-EPS)*(EPS/(2.D0-EPS))**PHORAN(THEDUM)
01693 DEL2=2.D0-DEL1
01694
01695
01696
01697
01698
01699
01700
01701
01702
01703
01704
01705
01706
01707
01708
01709
01710
01711
01712
01713
01714
01715
01716
01717
01718
01719 COSTHG=(1.D0-DEL1)/BETA
01720 SINTHG=SQRT(DEL1*DEL2-MCHREN)/BETA
01721 WGT=1D0
01722
01723
01724
01725
01726 ME=2.D0*PHOSPI(IDPHO(CHAPOI(NCHARG)))+1.D0
01727
01728
01729
01730
01731 DO I=JDAPHO(1,IP),JDAPHO(2,IP)
01732 IF (I.NE.CHAPOI(NCHARG)) THEN
01733 NEUDAU=I
01734 GOTO 51
01735 ENDIF
01736 ENDDO
01737
01738
01739 DATA=NCHARG
01740 CALL PHOERR(5,'PHOKIN',DATA)
01741 51 CONTINUE
01742 NCHARB=CHAPOI(NCHARG)
01743 NCHARB=NCHARB-JDAPHO(1,IP)+3
01744 NEUDAU=NEUDAU-JDAPHO(1,IP)+3
01745 WT=PHOCOR(MPASQR,MCHREN,ME)*WGT
01746 ENDIF
01747 ELSE
01748 DATA=PPHO(5,IP)-MASSUM
01749 CALL PHOERR(10,'PHOTOS',DATA)
01750 ENDIF
01751 ENDIF
01752
01753 RETURN
01754 END
01755 SUBROUTINE PHOOMA(IFIRST,ILAST,POINTR)
01756
01757
01758
01759
01760
01761
01762
01763
01764
01765
01766
01767
01768
01769
01770
01771
01772
01773
01774
01775 IMPLICIT NONE
01776 INTEGER NMXPHO
01777 PARAMETER (NMXPHO=10000)
01778 INTEGER IDPHO,ISTPHO,JDAPHO,JMOPHO,NEVPHO,NPHO
01779 REAL*8 PPHO,VPHO
01780 COMMON/PHOEVT/NEVPHO,NPHO,ISTPHO(NMXPHO),IDPHO(NMXPHO),
01781 &JMOPHO(2,NMXPHO),JDAPHO(2,NMXPHO),PPHO(5,NMXPHO),VPHO(4,NMXPHO)
01782 INTEGER IFIRST,ILAST,I,J,BUFPOI,POINTR(NMXPHO)
01783 REAL*8 BUFMAS,MASS(NMXPHO)
01784 IF (IFIRST.EQ.ILAST) RETURN
01785
01786
01787 DO 10 I=IFIRST,ILAST
01788 10 MASS(I)=PPHO(5,POINTR(I))
01789
01790
01791 DO 30 I=IFIRST,ILAST-1
01792 DO 20 J=I+1,ILAST
01793 IF (MASS(J).LE.MASS(I)) GOTO 20
01794 BUFPOI=POINTR(J)
01795 POINTR(J)=POINTR(I)
01796 POINTR(I)=BUFPOI
01797 BUFMAS=MASS(J)
01798 MASS(J)=MASS(I)
01799 MASS(I)=BUFMAS
01800 20 CONTINUE
01801 30 CONTINUE
01802 RETURN
01803 END
01804 SUBROUTINE PHOENE(MPASQR,MCHREN,BETA,BIGLOG,IDENT)
01805
01806
01807
01808
01809
01810
01811
01812
01813
01814
01815
01816
01817
01818
01819
01820
01821
01822
01823
01824
01825
01826 IMPLICIT NONE
01827 DOUBLE PRECISION MPASQR,MCHREN,BIGLOG,BETA,DATA
01828 INTEGER IWT1,IRN,IWT2
01829 REAL*8 PRSOFT,PRHARD,PHORAN,PHOFAC
01830 DOUBLE PRECISION MCHSQR,MNESQR
01831 REAL*8 PNEUTR
01832 INTEGER IDENT
01833 REAL*8 PHOCHA,PRKILL,RRR
01834 COMMON/PHOMOM/MCHSQR,MNESQR,PNEUTR(5)
01835 DOUBLE PRECISION COSTHG,SINTHG
01836 REAL*8 XPHMAX,XPHOTO
01837 COMMON/PHOPHS/XPHMAX,XPHOTO,COSTHG,SINTHG
01838 REAL*8 ALPHA,XPHCUT
01839 COMMON/PHOCOP/ALPHA,XPHCUT
01840 REAL*8 PI,TWOPI
01841 COMMON/PHPICO/PI,TWOPI
01842 INTEGER IREP
01843 REAL*8 PROBH,CORWT,XF
01844 COMMON/PHOPRO/PROBH,CORWT,XF,IREP
01845 LOGICAL INTERF,ISEC,ITRE,IEXP,IFTOP,IFW
01846 REAL*8 FINT,FSEC,EXPEPS
01847 COMMON /PHOKEY/ FSEC,FINT,EXPEPS,INTERF,ISEC,ITRE,IEXP,IFTOP,IFW
01848 INTEGER NX,NCHAN,K
01849 PARAMETER (NX=10)
01850 LOGICAL EXPINI
01851 REAL*8 PRO,PRSUM
01852 COMMON /PHOEXP/ PRO(NX),NCHAN,EXPINI
01853
01854 IF (XPHMAX.LE.XPHCUT) THEN
01855 BETA=PHOFAC(-1)
01856 XPHOTO=0.0D0
01857 RETURN
01858 ENDIF
01859
01860 MCHREN=4.D0*MCHSQR/MPASQR/(1.D0+MCHSQR/MPASQR)**2
01861 BETA=SQRT(1.D0-MCHREN)
01862
01863
01864
01865
01866
01867
01868
01869
01870
01871
01872
01873
01874 BIGLOG=LOG(MPASQR/MCHSQR*(1.D0+BETA)**2/4.D0*
01875 & (1.D0+MCHSQR/MPASQR)**2)
01876 PRHARD=ALPHA/PI*(1D0/BETA*BIGLOG)*
01877 &(LOG(XPHMAX/XPHCUT)-.75D0+XPHCUT/XPHMAX-.25D0*XPHCUT**2/XPHMAX**2)
01878 PRHARD=PRHARD*PHOCHA(IDENT)**2*FSEC*FINT
01879
01880 IF (IREP.EQ.0) PROBH=0.D0
01881 PRKILL=0d0
01882 IF (IEXP) THEN
01883 NCHAN=NCHAN+1
01884 IF (EXPINI) THEN
01885 PRO(NCHAN)=PRHARD+0.05*(1.0+FINT)
01886
01887 PRHARD=0D0
01888 PROBH=PRHARD
01889 ELSE
01890 PRSUM=0
01891 DO K=NCHAN,NX
01892 PRSUM=PRSUM+PRO(K)
01893 ENDDO
01894 PRHARD=PRHARD/PRSUM
01895
01896
01897
01898 PRKILL=PRO(NCHAN)/PRSUM-PRHARD
01899
01900 ENDIF
01901 PRSOFT=1.D0-PRHARD
01902 ELSE
01903 PRHARD=PRHARD*PHOFAC(0)
01904
01905
01906 PROBH=PRHARD
01907 ENDIF
01908 PRSOFT=1.D0-PRHARD
01909
01910
01911 IF (IEXP) THEN
01912 IF (PRSOFT.LT.-5.0D-8) THEN
01913 DATA=PRSOFT
01914 CALL PHOERR(2,'PHOENE',DATA)
01915 ENDIF
01916 ELSE
01917 IF (PRSOFT.LT.0.1D0) THEN
01918 DATA=PRSOFT
01919 CALL PHOERR(2,'PHOENE',DATA)
01920 ENDIF
01921 ENDIF
01922
01923 RRR=PHORAN(IWT1)
01924 IF (RRR.LT.PRSOFT) THEN
01925
01926
01927 XPHOTO=0.D0
01928 IF (RRR.LT.PRKILL) XPHOTO=-5d0
01929 ELSE
01930
01931
01932
01933 10 XPHOTO=EXP(PHORAN(IRN)*LOG(XPHCUT/XPHMAX))
01934 XPHOTO=XPHOTO*XPHMAX
01935 IF (PHORAN(IWT2).GT.((1.D0+(1.D0-XPHOTO/XPHMAX)**2)/2.D0))
01936 & GOTO 10
01937 ENDIF
01938
01939
01940 XF=4.D0*MCHSQR*MPASQR/(MPASQR+MCHSQR-MNESQR)**2
01941 RETURN
01942 END
01943 FUNCTION PHOCOR(MPASQR,MCHREN,ME)
01944
01945
01946
01947
01948
01949
01950
01951
01952
01953
01954
01955
01956
01957
01958
01959
01960
01961
01962
01963 IMPLICIT NONE
01964 DOUBLE PRECISION MPASQR,MCHREN,BETA,XX,YY,DATA
01965 INTEGER ME
01966 REAL*8 PHOCOR,PHOFAC,WT1,WT2,WT3
01967 DOUBLE PRECISION MCHSQR,MNESQR
01968 REAL*8 PNEUTR
01969 COMMON/PHOMOM/MCHSQR,MNESQR,PNEUTR(5)
01970 DOUBLE PRECISION COSTHG,SINTHG
01971 REAL*8 XPHMAX,XPHOTO
01972 COMMON/PHOPHS/XPHMAX,XPHOTO,COSTHG,SINTHG
01973 INTEGER IREP
01974 REAL*8 PROBH,CORWT,XF
01975 COMMON/PHOPRO/PROBH,CORWT,XF,IREP
01976
01977
01978 XX=4.D0*MCHSQR/MPASQR*(1.D0-XPHOTO)/(1.D0-XPHOTO+(MCHSQR-MNESQR)/
01979 &MPASQR)**2
01980 IF (ME.EQ.1) THEN
01981 YY=1.D0
01982 WT3=(1.D0-XPHOTO/XPHMAX)/((1.D0+(1.D0-XPHOTO/XPHMAX)**2)/2.D0)
01983 ELSEIF (ME.EQ.2) THEN
01984 YY=0.5D0*(1.D0-XPHOTO/XPHMAX+1.D0/(1.D0-XPHOTO/XPHMAX))
01985 WT3=1.D0
01986 ELSEIF ((ME.EQ.3).OR.(ME.EQ.4).OR.(ME.EQ.5)) THEN
01987 YY=1.D0
01988 WT3=(1.D0+(1.D0-XPHOTO/XPHMAX)**2-(XPHOTO/XPHMAX)**3)/
01989 & (1.D0+(1.D0-XPHOTO/XPHMAX)** 2)
01990 ELSE
01991 DATA=(ME-1.D0)/2.D0
01992 CALL PHOERR(6,'PHOCOR',DATA)
01993 YY=1.D0
01994 WT3=1.D0
01995 ENDIF
01996 BETA=SQRT(1.D0-XX)
01997 WT1=(1.D0-COSTHG*SQRT(1.D0-MCHREN))/(1.D0-COSTHG*BETA)
01998 WT2=(1.D0-XX/YY/(1.D0-BETA**2*COSTHG**2))*(1.D0+COSTHG*BETA)/2.D0
01999 WT2=WT2*PHOFAC(1)
02000 PHOCOR=WT1*WT2*WT3
02001 CORWT=PHOCOR
02002 IF (PHOCOR.GT.1.D0) THEN
02003 DATA=PHOCOR
02004 CALL PHOERR(3,'PHOCOR',DATA)
02005 ENDIF
02006 RETURN
02007 END
02008 FUNCTION PHOFAC(MODE)
02009
02010
02011
02012
02013
02014
02015
02016
02017
02018
02019
02020
02021
02022
02023
02024
02025
02026
02027
02028
02029
02030
02031
02032
02033
02034
02035
02036
02037
02038 IMPLICIT NONE
02039 REAL*8 PHOFAC,FF,PRX
02040 INTEGER MODE
02041 INTEGER IREP
02042 REAL*8 PROBH,CORWT,XF
02043 COMMON/PHOPRO/PROBH,CORWT,XF,IREP
02044 LOGICAL INTERF,ISEC,ITRE,IEXP,IFTOP,IFW
02045 REAL*8 FINT,FSEC,EXPEPS
02046 COMMON /PHOKEY/ FSEC,FINT,EXPEPS,INTERF,ISEC,ITRE,IEXP,IFTOP,IFW
02047 SAVE PRX,FF
02048 DATA PRX,FF/ 0.D0, 0.D0/
02049 IF (IEXP) THEN
02050 PHOFAC=1
02051 RETURN
02052 ENDIF
02053 IF (MODE.EQ.-1) THEN
02054 PRX=1.D0
02055 FF=1.D0
02056 PROBH=0.0
02057 ELSEIF (MODE.EQ.0) THEN
02058 IF (IREP.EQ.0) PRX=1.D0
02059 PRX=PRX/(1.D0-PROBH)
02060 FF=1.D0
02061
02062
02063
02064
02065
02066
02067 PHOFAC=FF*PRX
02068 ELSE
02069 PHOFAC=1.D0/FF
02070 ENDIF
02071 END
02072 SUBROUTINE PHOBW(WT)
02073
02074
02075
02076
02077
02078
02079
02080
02081
02082
02083
02084
02085
02086
02087
02088
02089
02090
02091
02092
02093
02094
02095 IMPLICIT NONE
02096 DOUBLE PRECISION WT
02097 INTEGER NMXPHO
02098 PARAMETER (NMXPHO=10000)
02099 INTEGER IDPHO,ISTPHO,JDAPHO,JMOPHO,NEVPHO,NPHO
02100 REAL*8 PPHO,VPHO
02101 COMMON/PHOEVT/NEVPHO,NPHO,ISTPHO(NMXPHO),IDPHO(NMXPHO),
02102 &JMOPHO(2,NMXPHO),JDAPHO(2,NMXPHO),PPHO(5,NMXPHO),VPHO(4,NMXPHO)
02103 INTEGER I
02104 DOUBLE PRECISION EMU,MCHREN,BETA,COSTHG,MPASQR,XPH
02105
02106 IF (ABS(IDPHO(1)).EQ.24.AND.
02107 $ ABS(IDPHO(JDAPHO(1,1) )).GE.11.AND.
02108 $ ABS(IDPHO(JDAPHO(1,1) )).LE.16.AND.
02109 $ ABS(IDPHO(JDAPHO(1,1)+1)).GE.11.AND.
02110 $ ABS(IDPHO(JDAPHO(1,1)+1)).LE.16 ) THEN
02111
02112 IF(
02113 $ ABS(IDPHO(JDAPHO(1,1) )).EQ.11.OR.
02114 $ ABS(IDPHO(JDAPHO(1,1) )).EQ.13.OR.
02115 $ ABS(IDPHO(JDAPHO(1,1) )).EQ.15 ) THEN
02116 I=JDAPHO(1,1)
02117 ELSE
02118 I=JDAPHO(1,1)+1
02119 ENDIF
02120 EMU=PPHO(4,I)
02121 MCHREN=ABS(PPHO(4,I)**2-PPHO(3,I)**2
02122 $ -PPHO(2,I)**2-PPHO(1,I)**2)
02123 BETA=SQRT(1- MCHREN/ PPHO(4,I)**2)
02124 COSTHG=(PPHO(3,I)*PPHO(3,NPHO)+PPHO(2,I)*PPHO(2,NPHO)
02125 $ +PPHO(1,I)*PPHO(1,NPHO))/
02126 $ SQRT(PPHO(3,I)**2+PPHO(2,I)**2+PPHO(1,I)**2) /
02127 $ SQRT(PPHO(3,NPHO)**2+PPHO(2,NPHO)**2+PPHO(1,NPHO)**2)
02128 MPASQR=PPHO(4,1)**2
02129 XPH=PPHO(4,NPHO)
02130 WT=WT*(1-8*EMU*XPH*(1-COSTHG*BETA)*
02131 $ (MCHREN+2*XPH*SQRT(MPASQR))/
02132 $ MPASQR**2/(1-MCHREN/MPASQR)/(4-MCHREN/MPASQR))
02133 ENDIF
02134
02135
02136
02137 END
02138 SUBROUTINE PHODO(IP,NCHARB,NEUDAU)
02139
02140
02141
02142
02143
02144
02145
02146
02147
02148
02149
02150
02151
02152
02153
02154
02155
02156
02157
02158
02159
02160
02161 IMPLICIT NONE
02162 DOUBLE PRECISION PHOAN1,PHOAN2,ANGLE,FI1,FI3,FI4,FI5,TH1,TH3,TH4
02163 DOUBLE PRECISION PARNE,QNEW,QOLD,DATA
02164 INTEGER IP,FI3DUM,I,J,NEUDAU,FIRST,LAST
02165 INTEGER NCHARB
02166 REAL*8 EPHOTO,PMAVIR,PHOTRI
02167 REAL*8 GNEUT,PHORAN,CCOSTH,SSINTH,PVEC(4)
02168 INTEGER NMXPHO
02169 PARAMETER (NMXPHO=10000)
02170 INTEGER IDPHO,ISTPHO,JDAPHO,JMOPHO,NEVPHO,NPHO
02171 REAL*8 PPHO,VPHO
02172 COMMON/PHOEVT/NEVPHO,NPHO,ISTPHO(NMXPHO),IDPHO(NMXPHO),
02173 &JMOPHO(2,NMXPHO),JDAPHO(2,NMXPHO),PPHO(5,NMXPHO),VPHO(4,NMXPHO)
02174 DOUBLE PRECISION MCHSQR,MNESQR
02175 REAL*8 PNEUTR
02176 COMMON/PHOMOM/MCHSQR,MNESQR,PNEUTR(5)
02177 DOUBLE PRECISION COSTHG,SINTHG
02178 REAL*8 XPHMAX,XPHOTO
02179 COMMON/PHOPHS/XPHMAX,XPHOTO,COSTHG,SINTHG
02180 REAL*8 PI,TWOPI
02181 COMMON/PHPICO/PI,TWOPI
02182
02183 EPHOTO=XPHOTO*PPHO(5,IP)/2.D0
02184 PMAVIR=SQRT(PPHO(5,IP)*(PPHO(5,IP)-2.D0*EPHOTO))
02185
02186
02187 FI1=PHOAN1(PNEUTR(1),PNEUTR(2))
02188
02189
02190
02191 TH1=PHOAN2(PNEUTR(3),SQRT(PNEUTR(1)**2+PNEUTR(2)**2))
02192 CALL PHORO3(-FI1,PNEUTR(1))
02193 CALL PHORO2(-TH1,PNEUTR(1))
02194
02195
02196
02197
02198
02199
02200
02201 QNEW=PHOTRI(PMAVIR,PNEUTR(5),PPHO(5,NCHARB))
02202 QOLD=PNEUTR(3)
02203 GNEUT=(QNEW**2+QOLD**2+MNESQR)/(QNEW*QOLD+SQRT((QNEW**2+MNESQR)*
02204 &(QOLD**2+MNESQR)))
02205 IF (GNEUT.LT.1.D0) THEN
02206 DATA=0.D0
02207 CALL PHOERR(4,'PHOKIN',DATA)
02208 ENDIF
02209 PARNE=GNEUT-SQRT(MAX(GNEUT**2-1.0D0,0.D0))
02210
02211
02212 CALL PHOBO3(PARNE,PNEUTR)
02213
02214
02215 NPHO=NPHO+1
02216 ISTPHO(NPHO)=1
02217 IDPHO(NPHO) =22
02218
02219 JMOPHO(1,NPHO)=IP
02220 JMOPHO(2,NPHO)=0
02221 JDAPHO(1,NPHO)=0
02222 JDAPHO(2,NPHO)=0
02223 PPHO(4,NPHO)=EPHOTO*PPHO(5,IP)/PMAVIR
02224
02225
02226 CCOSTH=-COSTHG
02227 SSINTH=SINTHG
02228 TH3=PHOAN2(CCOSTH,SSINTH)
02229 FI3=TWOPI*PHORAN(FI3DUM)
02230 PPHO(1,NPHO)=PPHO(4,NPHO)*SINTHG*COS(FI3)
02231 PPHO(2,NPHO)=PPHO(4,NPHO)*SINTHG*SIN(FI3)
02232
02233
02234 PPHO(3,NPHO)=-PPHO(4,NPHO)*COSTHG
02235 PPHO(5,NPHO)=0.D0
02236
02237
02238 CALL PHORO3(-FI3,PNEUTR(1))
02239 CALL PHORO3(-FI3,PPHO(1,NPHO))
02240 CALL PHORO2(-TH3,PNEUTR(1))
02241 CALL PHORO2(-TH3,PPHO(1,NPHO))
02242 ANGLE=EPHOTO/PPHO(4,NPHO)
02243
02244
02245 CALL PHOBO3(ANGLE,PNEUTR(1))
02246 CALL PHOBO3(ANGLE,PPHO(1,NPHO))
02247
02248
02249 FI4=PHOAN1(PNEUTR(1),PNEUTR(2))
02250 TH4=PHOAN2(PNEUTR(3),SQRT(PNEUTR(1)**2+PNEUTR(2)**2))
02251 CALL PHORO3(FI4,PNEUTR(1))
02252 CALL PHORO3(FI4,PPHO(1,NPHO))
02253
02254 DO 60 I=2,4
02255 60 PVEC(I)=0.D0
02256 PVEC(1)=1.D0
02257 CALL PHORO3(-FI3,PVEC)
02258 CALL PHORO2(-TH3,PVEC)
02259 CALL PHOBO3(ANGLE,PVEC)
02260 CALL PHORO3(FI4,PVEC)
02261 CALL PHORO2(-TH4,PNEUTR)
02262 CALL PHORO2(-TH4,PPHO(1,NPHO))
02263 CALL PHORO2(-TH4,PVEC)
02264 FI5=PHOAN1(PVEC(1),PVEC(2))
02265
02266
02267 CALL PHORO3(-FI5,PNEUTR)
02268 CALL PHORO3(-FI5,PPHO(1,NPHO))
02269 CALL PHORO2(TH1,PNEUTR(1))
02270 CALL PHORO2(TH1,PPHO(1,NPHO))
02271 CALL PHORO3(FI1,PNEUTR)
02272 CALL PHORO3(FI1,PPHO(1,NPHO))
02273
02274 IF ((JDAPHO(2,IP)-JDAPHO(1,IP)).GT.1) THEN
02275
02276
02277 FIRST=NEUDAU
02278 LAST=JDAPHO(2,IP)
02279 DO 70 I=FIRST,LAST
02280 IF (I.NE.NCHARB.AND.(JMOPHO(1,I).EQ.IP)) THEN
02281
02282
02283 CALL PHORO3(-FI1,PPHO(1,I))
02284 CALL PHORO2(-TH1,PPHO(1,I))
02285
02286
02287 CALL PHOBO3(PARNE,PPHO(1,I))
02288
02289
02290 CALL PHORO3(-FI3,PPHO(1,I))
02291 CALL PHORO2(-TH3,PPHO(1,I))
02292
02293
02294 CALL PHOBO3(ANGLE,PPHO(1,I))
02295
02296
02297 CALL PHORO3(FI4,PPHO(1,I))
02298 CALL PHORO2(-TH4,PPHO(1,I))
02299
02300
02301 CALL PHORO3(-FI5,PPHO(1,I))
02302 CALL PHORO2(TH1,PPHO(1,I))
02303 CALL PHORO3(FI1,PPHO(1,I))
02304 ENDIF
02305 70 CONTINUE
02306 ELSE
02307
02308
02309 DO 80 J=1,4
02310 80 PPHO(J,NEUDAU)=PNEUTR(J)
02311 ENDIF
02312
02313
02314 DO 90 J=1,3
02315 90 PPHO(J,NCHARB)=-(PPHO(J,NPHO)+PNEUTR(J))
02316 PPHO(4,NCHARB)=PPHO(5,IP)-(PPHO(4,NPHO)+PNEUTR(4))
02317
02318 END
02319 FUNCTION PHOTRI(A,B,C)
02320
02321
02322
02323
02324
02325
02326
02327
02328
02329
02330
02331
02332
02333
02334
02335 IMPLICIT NONE
02336 DOUBLE PRECISION DA,DB,DC,DAPB,DAMB,DTRIAN
02337 REAL*8 A,B,C,PHOTRI
02338 DA=A
02339 DB=B
02340 DC=C
02341 DAPB=DA+DB
02342 DAMB=DA-DB
02343 DTRIAN=SQRT((DAMB-DC)*(DAPB+DC)*(DAMB+DC)*(DAPB-DC))
02344 PHOTRI=DTRIAN/(DA+DA)
02345 RETURN
02346 END
02347 FUNCTION PHOAN1(X,Y)
02348
02349
02350
02351
02352
02353
02354
02355
02356
02357
02358
02359
02360
02361
02362 IMPLICIT NONE
02363 DOUBLE PRECISION PHOAN1
02364 REAL*8 X,Y
02365 REAL*8 PI,TWOPI
02366 COMMON/PHPICO/PI,TWOPI
02367 IF (ABS(Y).LT.ABS(X)) THEN
02368 PHOAN1=ATAN(ABS(Y/X))
02369 IF (X.LE.0.D0) PHOAN1=PI-PHOAN1
02370 ELSE
02371 PHOAN1=ACOS(X/SQRT(X**2+Y**2))
02372 ENDIF
02373 IF (Y.LT.0.D0) PHOAN1=TWOPI-PHOAN1
02374 RETURN
02375 END
02376 FUNCTION PHOAN2(X,Y)
02377
02378
02379
02380
02381
02382
02383
02384
02385
02386
02387
02388
02389
02390
02391 IMPLICIT NONE
02392 DOUBLE PRECISION PHOAN2
02393 REAL*8 X,Y
02394 REAL*8 PI,TWOPI
02395 COMMON/PHPICO/PI,TWOPI
02396 IF (ABS(Y).LT.ABS(X)) THEN
02397 PHOAN2=ATAN(ABS(Y/X))
02398 IF (X.LE.0.D0) PHOAN2=PI-PHOAN2
02399 ELSE
02400 PHOAN2=ACOS(X/SQRT(X**2+Y**2))
02401 ENDIF
02402 RETURN
02403 END
02404 SUBROUTINE PHOBO3(ANGLE,PVEC)
02405
02406
02407
02408
02409
02410
02411
02412
02413
02414
02415
02416
02417
02418
02419
02420 IMPLICIT NONE
02421 DOUBLE PRECISION QPL,QMI,ANGLE
02422 REAL*8 PVEC(4)
02423 QPL=(PVEC(4)+PVEC(3))*ANGLE
02424 QMI=(PVEC(4)-PVEC(3))/ANGLE
02425 PVEC(3)=(QPL-QMI)/2.D0
02426 PVEC(4)=(QPL+QMI)/2.D0
02427 RETURN
02428 END
02429 SUBROUTINE PHORO2(ANGLE,PVEC)
02430
02431
02432
02433
02434
02435
02436
02437
02438
02439
02440
02441
02442
02443
02444
02445 IMPLICIT NONE
02446 DOUBLE PRECISION CS,SN,ANGLE
02447 REAL*8 PVEC(4)
02448 CS=COS(ANGLE)*PVEC(1)+SIN(ANGLE)*PVEC(3)
02449 SN=-SIN(ANGLE)*PVEC(1)+COS(ANGLE)*PVEC(3)
02450 PVEC(1)=CS
02451 PVEC(3)=SN
02452 RETURN
02453 END
02454 SUBROUTINE PHORO3(ANGLE,PVEC)
02455
02456
02457
02458
02459
02460
02461
02462
02463
02464
02465
02466
02467
02468
02469
02470 IMPLICIT NONE
02471 DOUBLE PRECISION CS,SN,ANGLE
02472 REAL*8 PVEC(4)
02473 CS=COS(ANGLE)*PVEC(1)-SIN(ANGLE)*PVEC(2)
02474 SN=SIN(ANGLE)*PVEC(1)+COS(ANGLE)*PVEC(2)
02475 PVEC(1)=CS
02476 PVEC(2)=SN
02477 RETURN
02478 END
02479 SUBROUTINE PHORIN
02480
02481
02482
02483
02484
02485
02486
02487
02488
02489
02490
02491
02492
02493
02494
02495
02496 IMPLICIT NONE
02497 DOUBLE PRECISION DATA
02498 REAL*8 S,T
02499 INTEGER I,IS1,IS2,IS3,IS4,IS5,J
02500 INTEGER ISEED,I97,J97
02501 REAL*8 URAN,CRAN,CDRAN,CMRAN
02502 COMMON/PHSEED/ISEED(2),I97,J97,URAN(97),CRAN,CDRAN,CMRAN
02503
02504
02505 IF ((ISEED(1).LT.0).OR.(ISEED(1).GE.31328)) THEN
02506 DATA=ISEED(1)
02507 CALL PHOERR(8,'PHORIN',DATA)
02508 ENDIF
02509 IF ((ISEED(2).LT.0).OR.(ISEED(2).GE.30081)) THEN
02510 DATA=ISEED(2)
02511 CALL PHOERR(9,'PHORIN',DATA)
02512 ENDIF
02513
02514
02515 IS1=MOD(ISEED(1)/177,177)+2
02516 IS2=MOD(ISEED(1),177)+2
02517 IS3=MOD(ISEED(2)/169,178)+1
02518 IS4=MOD(ISEED(2),169)
02519 DO 20 I=1,97
02520 S=0.D0
02521 T=0.5D0
02522 DO 10 J=1,24
02523 IS5=MOD (MOD(IS1*IS2,179)*IS3,179)
02524 IS1=IS2
02525 IS2=IS3
02526 IS3=IS5
02527 IS4=MOD(53*IS4+1,169)
02528 IF (MOD(IS4*IS5,64).GE.32) S=S+T
02529 10 T=0.5D0*T
02530 20 URAN(I)=S
02531 CRAN=362436.D0/16777216.D0
02532 CDRAN=7654321.D0/16777216.D0
02533 CMRAN=16777213.D0/16777216.D0
02534 I97=97
02535 J97=33
02536 RETURN
02537 END
02538 FUNCTION PHORAN(IDUM)
02539
02540
02541
02542
02543
02544
02545
02546
02547
02548
02549
02550
02551
02552
02553
02554
02555
02556
02557
02558
02559 IMPLICIT NONE
02560 REAL*8 PHORAN
02561 INTEGER IDUM
02562 INTEGER ISEED,I97,J97
02563 REAL*8 URAN,CRAN,CDRAN,CMRAN
02564 COMMON/PHSEED/ISEED(2),I97,J97,URAN(97),CRAN,CDRAN,CMRAN
02565 10 PHORAN=URAN(I97)-URAN(J97)
02566 IF (PHORAN.LT.0.D0) PHORAN=PHORAN+1.D0
02567 URAN(I97)=PHORAN
02568 I97=I97-1
02569 IF (I97.EQ.0) I97=97
02570 J97=J97-1
02571 IF (J97.EQ.0) J97=97
02572 CRAN=CRAN-CDRAN
02573 IF (CRAN.LT.0.D0) CRAN=CRAN+CMRAN
02574 PHORAN=PHORAN-CRAN
02575 IF (PHORAN.LT.0.D0) PHORAN=PHORAN+1.D0
02576 IF (PHORAN.LE.0.D0) GOTO 10
02577 RETURN
02578 END
02579 FUNCTION PHOCHA(IDHEP)
02580
02581
02582
02583
02584
02585
02586
02587
02588
02589
02590
02591
02592
02593
02594
02595
02596
02597 IMPLICIT NONE
02598 REAL*8 PHOCHA
02599 INTEGER IDHEP,IDABS,Q1,Q2,Q3
02600
02601
02602
02603 REAL*8 CHARGE(0:100)
02604 DATA CHARGE/ 0.D0,
02605 &-0.3333333333D0, 0.6666666667D0, -0.3333333333D0, 0.6666666667D0,
02606 &-0.3333333333D0, 0.6666666667D0, -0.3333333333D0, 0.6666666667D0,
02607 & 2*0.D0, -1.D0, 0.D0, -1.D0, 0.D0, -1.D0, 0.D0, -1.D0, 6*0.D0,
02608 & 1.D0, 12*0.D0, 1.D0, 63*0.D0/
02609 IDABS=ABS(IDHEP)
02610 IF (IDABS.LE.100) THEN
02611
02612
02613 PHOCHA = CHARGE(IDABS)
02614 ELSE
02615
02616
02617 Q3=MOD(IDABS/1000,10)
02618 Q2=MOD(IDABS/100,10)
02619 Q1=MOD(IDABS/10,10)
02620 IF (Q3.EQ.0) THEN
02621
02622
02623 IF(MOD(Q2,2).EQ.0) THEN
02624 PHOCHA=CHARGE(Q2)-CHARGE(Q1)
02625 ELSE
02626 PHOCHA=CHARGE(Q1)-CHARGE(Q2)
02627 ENDIF
02628 ELSE
02629
02630
02631 PHOCHA=CHARGE(Q1)+CHARGE(Q2)+CHARGE(Q3)
02632 ENDIF
02633 ENDIF
02634
02635
02636 IF (IDHEP.LT.0.D0) PHOCHA=-PHOCHA
02637 IF (PHOCHA**2.lt.1d-6) PHOCHA=0.D0
02638 RETURN
02639 END
02640 FUNCTION PHOSPI(IDHEP)
02641
02642
02643
02644
02645
02646
02647
02648
02649
02650
02651
02652
02653
02654
02655
02656
02657
02658
02659 IMPLICIT NONE
02660 REAL*8 PHOSPI
02661 INTEGER IDHEP,IDABS
02662
02663
02664
02665 REAL*8 SPIN(100)
02666 DATA SPIN/ 8*.5D0, 1.D0, 0.D0, 8*.5D0, 2*0.D0, 4*1.D0, 76*0.D0/
02667 IDABS=ABS(IDHEP)
02668
02669
02670 IF (IDABS.LE.100) THEN
02671 PHOSPI=SPIN(IDABS)
02672 ELSE
02673
02674
02675 PHOSPI=(MOD(IDABS,10)-1.D0)/2.D0
02676
02677
02678 PHOSPI=MAX(PHOSPI,0.D0)
02679 ENDIF
02680 RETURN
02681 END
02682 SUBROUTINE PHOERR(IMES,TEXT,DATA)
02683
02684
02685
02686
02687
02688
02689
02690
02691
02692
02693
02694
02695
02696
02697
02698 IMPLICIT NONE
02699 DOUBLE PRECISION DATA
02700 INTEGER IMES,IERROR
02701 REAL*8 SDATA
02702 INTEGER PHLUN
02703 COMMON/PHOLUN/PHLUN
02704 INTEGER PHOMES
02705 PARAMETER (PHOMES=10)
02706 INTEGER STATUS
02707 COMMON/PHOSTA/STATUS(PHOMES)
02708 CHARACTER TEXT*(*)
02709 SAVE IERROR
02710
02711 LOGICAL ISEC
02712 SAVE ISEC
02713 DATA ISEC /.TRUE./
02714 DATA IERROR/ 0/
02715 IF (IMES.LE.PHOMES) STATUS(IMES)=STATUS(IMES)+1
02716
02717
02718 IF ((IMES.EQ. 6).AND.(STATUS(IMES).GE.2)) RETURN
02719 IF ((IMES.EQ.10).AND.(STATUS(IMES).GE.2)) RETURN
02720 SDATA=DATA
02721 WRITE(PHLUN,9000)
02722 WRITE(PHLUN,9120)
02723 GOTO (10,20,30,40,50,60,70,80,90,100),IMES
02724 WRITE(PHLUN,9130) IMES
02725 GOTO 120
02726 10 WRITE(PHLUN,9010) TEXT,INT(SDATA)
02727 GOTO 110
02728 20 WRITE(PHLUN,9020) TEXT,SDATA
02729 GOTO 110
02730 30 WRITE(PHLUN,9030) TEXT,SDATA
02731 GOTO 110
02732 40 WRITE(PHLUN,9040) TEXT
02733 GOTO 110
02734 50 WRITE(PHLUN,9050) TEXT,INT(SDATA)
02735 GOTO 110
02736 60 WRITE(PHLUN,9060) TEXT,SDATA
02737 GOTO 130
02738 70 WRITE(PHLUN,9070) TEXT,INT(SDATA)
02739 GOTO 110
02740 80 WRITE(PHLUN,9080) TEXT,INT(SDATA)
02741 GOTO 110
02742 90 WRITE(PHLUN,9090) TEXT,INT(SDATA)
02743 GOTO 110
02744 100 WRITE(PHLUN,9100) TEXT,SDATA
02745 GOTO 130
02746 110 CONTINUE
02747 WRITE(PHLUN,9140)
02748 WRITE(PHLUN,9120)
02749 WRITE(PHLUN,9000)
02750 IF (ISEC) THEN
02751 STOP
02752 ELSE
02753 GOTO 130
02754 ENDIF
02755 120 IERROR=IERROR+1
02756 IF (IERROR.GE.10) THEN
02757 WRITE(PHLUN,9150)
02758 WRITE(PHLUN,9120)
02759 WRITE(PHLUN,9000)
02760 IF (ISEC) THEN
02761 STOP
02762 ELSE
02763 GOTO 130
02764 ENDIF
02765 ENDIF
02766 130 WRITE(PHLUN,9120)
02767 WRITE(PHLUN,9000)
02768 RETURN
02769 9000 FORMAT(1H ,80('*'))
02770 9010 FORMAT(1H ,'* ',A,': Too many charged Particles, NCHARG =',I6,T81,
02771 &'*')
02772 9020 FORMAT(1H ,'* ',A,': Too much Bremsstrahlung required, PRSOFT = ',
02773 &F15.6,T81,'*')
02774 9030 FORMAT(1H ,'* ',A,': Combined Weight is exceeding 1., Weight = ',
02775 &F15.6,T81,'*')
02776 9040 FORMAT(1H ,'* ',A,
02777 &': Error in Rescaling charged and neutral Vectors',T81,'*')
02778 9050 FORMAT(1H ,'* ',A,
02779 &': Non matching charged Particle Pointer, NCHARG = ',I5,T81,'*')
02780 9060 FORMAT(1H ,'* ',A,
02781 &': Do you really work with a Particle of Spin: ',F4.1,' ?',T81,
02782 &'*')
02783 9070 FORMAT(1H ,'* ',A, ': Stack Length exceeded, NSTACK = ',I5 ,T81,
02784 &'*')
02785 9080 FORMAT(1H ,'* ',A,
02786 &': Random Number Generator Seed(1) out of Range: ',I8,T81,'*')
02787 9090 FORMAT(1H ,'* ',A,
02788 &': Random Number Generator Seed(2) out of Range: ',I8,T81,'*')
02789 9100 FORMAT(1H ,'* ',A,
02790 &': Available Phase Space below Cut-off: ',F15.6,' GeV/c^2',T81,
02791 &'*')
02792 9120 FORMAT(1H ,'*',T81,'*')
02793 9130 FORMAT(1H ,'* Funny Error Message: ',I4,' ! What to do ?',T81,'*')
02794 9140 FORMAT(1H ,'* Fatal Error Message, I stop this Run !',T81,'*')
02795 9150 FORMAT(1H ,'* 10 Error Messages generated, I stop this Run !',T81,
02796 &'*')
02797 END
02798 SUBROUTINE PHOREP
02799
02800
02801
02802
02803
02804
02805
02806
02807
02808
02809
02810
02811
02812
02813
02814 IMPLICIT NONE
02815 INTEGER PHLUN
02816 COMMON/PHOLUN/PHLUN
02817 INTEGER PHOMES
02818 PARAMETER (PHOMES=10)
02819 INTEGER STATUS
02820 COMMON/PHOSTA/STATUS(PHOMES)
02821 INTEGER I
02822 LOGICAL ERROR
02823 ERROR=.FALSE.
02824 WRITE(PHLUN,9000)
02825 WRITE(PHLUN,9010)
02826 WRITE(PHLUN,9020)
02827 WRITE(PHLUN,9030)
02828 WRITE(PHLUN,9040)
02829 WRITE(PHLUN,9030)
02830 WRITE(PHLUN,9020)
02831 DO 10 I=1,PHOMES
02832 IF (STATUS(I).EQ.0) GOTO 10
02833 IF ((I.EQ.6).OR.(I.EQ.10)) THEN
02834 WRITE(PHLUN,9050) I,STATUS(I)
02835 ELSE
02836 ERROR=.TRUE.
02837 WRITE(PHLUN,9060) I,STATUS(I)
02838 ENDIF
02839 10 CONTINUE
02840 IF (.NOT.ERROR) WRITE(PHLUN,9070)
02841 WRITE(PHLUN,9020)
02842 WRITE(PHLUN,9010)
02843 RETURN
02844 9000 FORMAT(1H1)
02845 9010 FORMAT(1H ,80('*'))
02846 9020 FORMAT(1H ,'*',T81,'*')
02847 9030 FORMAT(1H ,'*',26X,25('='),T81,'*')
02848 9040 FORMAT(1H ,'*',30X,'PHOTOS Run Summary',T81,'*')
02849 9050 FORMAT(1H ,'*',22X,'Warning #',I2,' occured',I6,' times',T81,'*')
02850 9060 FORMAT(1H ,'*',23X,'Error #',I2,' occured',I6,' times',T81,'*')
02851 9070 FORMAT(1H ,'*',16X,'PHOTOS Execution has successfully terminated',
02852 &T81,'*')
02853 END
02854 SUBROUTINE PHLUPA(IPOINT)
02855 IMPLICIT NONE
02856
02857
02858
02859
02860
02861
02862
02863
02864
02865
02866
02867
02868
02869
02870
02871
02872 INTEGER NMXPHO
02873 PARAMETER (NMXPHO=10000)
02874 INTEGER IDPHO,ISTPHO,JDAPHO,JMOPHO,NEVPHO,NPHO,I,J,IPOINT
02875 INTEGER IPOIN,IPOIN0,IPOINM,IEV
02876 INTEGER IOUT
02877 REAL*8 PPHO,VPHO,SUM
02878 COMMON/PHOEVT/NEVPHO,NPHO,ISTPHO(NMXPHO),IDPHO(NMXPHO),
02879 &JMOPHO(2,NMXPHO),JDAPHO(2,NMXPHO),PPHO(5,NMXPHO),VPHO(4,NMXPHO)
02880 COMMON /PHNUM/ IEV
02881 INTEGER PHLUN
02882 COMMON/PHOLUN/PHLUN
02883 DIMENSION SUM(5)
02884 DATA IPOIN0/ -5/
02885 COMMON /PHLUPY/ IPOIN,IPOINM
02886 SAVE IPOIN0
02887 IF (IPOIN0.LT.0) THEN
02888 IPOIN0=300 000
02889 IPOIN =IPOIN0
02890 IPOINM=300 001
02891 ENDIF
02892 IF (IPOINT.LE.IPOINM.OR.IPOINT.GE.IPOIN ) RETURN
02893 IOUT=56
02894 IF (IEV.LT.1000) THEN
02895 DO I=1,5
02896 SUM(I)=0.0D0
02897 ENDDO
02898 WRITE(PHLUN,*) 'EVENT NR=',IEV,
02899 $ 'WE ARE TESTING /PHOEVT/ at IPOINT=',IPOINT
02900 WRITE(PHLUN,10)
02901 I=1
02902 WRITE(PHLUN,20) IDPHO(I),PPHO(1,I),PPHO(2,I),PPHO(3,I),
02903 $ PPHO(4,I),PPHO(5,I),JDAPHO(1,I),JDAPHO(2,I)
02904 I=2
02905 WRITE(PHLUN,20) IDPHO(I),PPHO(1,I),PPHO(2,I),PPHO(3,I),
02906 $ PPHO(4,I),PPHO(5,I),JDAPHO(1,I),JDAPHO(2,I)
02907 WRITE(PHLUN,*) ' '
02908 DO I=3,NPHO
02909 WRITE(PHLUN,20) IDPHO(I),PPHO(1,I),PPHO(2,I),PPHO(3,I),
02910 $ PPHO(4,I),PPHO(5,I),JMOPHO(1,I),JMOPHO(2,I)
02911 DO J=1,4
02912 SUM(J)=SUM(J)+PPHO(J,I)
02913 ENDDO
02914 ENDDO
02915 SUM(5)=SQRT(ABS(SUM(4)**2-SUM(1)**2-SUM(2)**2-SUM(3)**2))
02916 WRITE(PHLUN,30) SUM
02917 10 FORMAT(1X,' ID ','p_x ','p_y ','p_z ',
02918 $ 'E ','m ',
02919 $ 'ID-MO_DA1','ID-MO DA2' )
02920 20 FORMAT(1X,I4,5(F14.9),2I9)
02921 30 FORMAT(1X,' SUM',5(F14.9))
02922 ENDIF
02923 END
02924
02925
02926
02927 FUNCTION IPHQRK(MODCOR)
02928 implicit none
02929
02930
02931
02932
02933
02934
02935
02936
02937
02938
02939
02940
02941
02942
02943
02944 INTEGER IPHQRK,MODCOR,MODOP
02945 INTEGER PHLUN
02946 COMMON/PHOLUN/PHLUN
02947
02948 SAVE MODOP
02949 DATA MODOP /0/
02950 IF (MODCOR.NE.0) THEN
02951
02952 MODOP=MODCOR
02953
02954 WRITE(PHLUN,*)
02955 $ 'Message from PHOTOS: IPHQRK(MODCOR):: (re)initialization'
02956 IF (MODOP.EQ.1) THEN
02957 WRITE(PHLUN,*)
02958 $ 'MODOP=1 -- blocks emission from light quarks: DEFAULT'
02959 ELSEIF (MODOP.EQ.2) THEN
02960 WRITE(PHLUN,*)
02961 $ 'MODOP=2 -- enables emission from light quarks: TEST '
02962 ELSE
02963 WRITE(PHLUN,*) 'IPHQRK wrong MODCOR=',MODCOR
02964 STOP
02965 ENDIF
02966 RETURN
02967 ENDIF
02968
02969 IF (MODOP.EQ.0.AND.MODCOR.EQ.0) THEN
02970 WRITE(PHLUN,*) 'IPHQRK lack of initialization'
02971 STOP
02972 ENDIF
02973 IPHQRK=MODOP
02974 END
02975
02976
02977 FUNCTION IPHEKL(MODCOR)
02978 implicit none
02979
02980
02981
02982
02983
02984
02985
02986
02987
02988
02989
02990
02991
02992
02993
02994 INTEGER IPHEKL,MODCOR,MODOP
02995 INTEGER PHLUN
02996 COMMON/PHOLUN/PHLUN
02997
02998 SAVE MODOP
02999 DATA MODOP /0/
03000
03001 IF (MODCOR.NE.0) THEN
03002
03003 MODOP=MODCOR
03004
03005 WRITE(PHLUN,*)
03006 $ 'Message from PHOTOS: IPHEKL(MODCOR):: (re)initialization'
03007 IF (MODOP.EQ.2) THEN
03008 WRITE(PHLUN,*)
03009 $ 'MODOP=2 -- blocks emission in pi0 to gamma e+e-: DEFAULT'
03010 WRITE(PHLUN,*)
03011 $ 'MODOP=2 -- blocks emission in Kl to gamma e+e-: DEFAULT'
03012 ELSEIF (MODOP.EQ.1) THEN
03013 WRITE(PHLUN,*)
03014 $ 'MODOP=1 -- enables emission in pi0 to gamma e+e- : TEST '
03015 WRITE(PHLUN,*)
03016 $ 'MODOP=1 -- enables emission in Kl to gamma e+e- : TEST '
03017 ELSE
03018 WRITE(PHLUN,*) 'IPHEKL wrong MODCOR=',MODCOR
03019 STOP
03020 ENDIF
03021 RETURN
03022 ENDIF
03023
03024 IF (MODOP.EQ.0.AND.MODCOR.EQ.0) THEN
03025 WRITE(PHLUN,*) 'IPHELK lack of initialization'
03026 STOP
03027 ENDIF
03028 IPHEKL=MODOP
03029 END
03030
03031 SUBROUTINE PHCORK(MODCOR)
03032 implicit none
03033
03034
03035
03036
03037
03038
03039
03040
03041
03042
03043
03044
03045
03046
03047
03048
03049
03050
03051
03052
03053
03054
03055 INTEGER NMXPHO
03056 PARAMETER (NMXPHO=10000)
03057
03058 REAL*8 M,P2,PX,PY,PZ,E,EN,MCUT,XMS
03059 INTEGER MODCOR,MODOP,I,IEV,IPRINT,K
03060 INTEGER IDPHO,ISTPHO,JDAPHO,JMOPHO,NEVPHO,NPHO
03061 REAL*8 PPHO,VPHO
03062 COMMON/PHOEVT/NEVPHO,NPHO,ISTPHO(NMXPHO),IDPHO(NMXPHO),
03063 &JMOPHO(2,NMXPHO),JDAPHO(2,NMXPHO),PPHO(5,NMXPHO),VPHO(4,NMXPHO)
03064
03065 INTEGER PHLUN
03066 COMMON/PHOLUN/PHLUN
03067
03068 COMMON /PHNUM/ IEV
03069 SAVE MODOP
03070 DATA MODOP /0/
03071 SAVE IPRINT
03072 DATA IPRINT /0/
03073 SAVE MCUT
03074 IF (MODCOR.NE.0) THEN
03075
03076 MODOP=MODCOR
03077
03078 WRITE(PHLUN,*) 'Message from PHCORK(MODCOR):: initialization'
03079 IF (MODOP.EQ.1) THEN
03080 WRITE(PHLUN,*) 'MODOP=1 -- no corrections on event: DEFAULT'
03081 ELSEIF (MODOP.EQ.2) THEN
03082 WRITE(PHLUN,*) 'MODOP=2 -- corrects Energy from mass'
03083 ELSEIF (MODOP.EQ.3) THEN
03084 WRITE(PHLUN,*) 'MODOP=3 -- corrects mass from Energy'
03085 ELSEIF (MODOP.EQ.4) THEN
03086 WRITE(PHLUN,*) 'MODOP=4 -- corrects Energy from mass to Mcut'
03087 WRITE(PHLUN,*) ' and mass from energy above Mcut '
03088 MCUT=0.4
03089 WRITE(PHLUN,*) 'Mcut=',MCUT,'GeV'
03090 ELSEIF (MODOP.EQ.5) THEN
03091 WRITE(PHLUN,*) 'MODOP=5 -- corrects Energy from mass+flow'
03092
03093 ELSE
03094 WRITE(PHLUN,*) 'PHCORK wrong MODCOR=',MODCOR
03095 STOP
03096 ENDIF
03097 RETURN
03098 ENDIF
03099
03100 IF (MODOP.EQ.0.AND.MODCOR.EQ.0) THEN
03101 WRITE(PHLUN,*) 'PHCORK lack of initialization'
03102 STOP
03103 ENDIF
03104
03105
03106
03107
03108
03109
03110 PX=0
03111 PY=0
03112 PZ=0
03113 E =0
03114
03115 IF (MODOP.EQ.1) THEN
03116
03117
03118 RETURN
03119 ELSEIF(MODOP.EQ.2) THEN
03120
03121
03122
03123
03124 DO I=3,NPHO
03125
03126 PX=PX+PPHO(1,I)
03127 PY=PY+PPHO(2,I)
03128 PZ=PZ+PPHO(3,I)
03129
03130 P2=PPHO(1,I)**2+PPHO(2,I)**2+PPHO(3,I)**2
03131
03132 EN=SQRT( PPHO(5,I)**2 + P2)
03133
03134 IF (IPRINT.EQ.1)
03135 & WRITE(PHLUN,*) "CORRECTING ENERGY OF ",I,":",
03136 & PPHO(4,I),"=>",EN
03137
03138 PPHO(4,I)=EN
03139 E = E+PPHO(4,I)
03140
03141 ENDDO
03142
03143 ELSEIF(MODOP.EQ.5) THEN
03144
03145
03146
03147
03148 DO I=3,NPHO
03149
03150 PX=PX+PPHO(1,I)
03151 PY=PY+PPHO(2,I)
03152 PZ=PZ+PPHO(3,I)
03153
03154 P2=PPHO(1,I)**2+PPHO(2,I)**2+PPHO(3,I)**2
03155
03156 EN=SQRT( PPHO(5,I)**2 + P2)
03157
03158 IF (IPRINT.EQ.1)
03159 & WRITE(PHLUN,*) "CORRECTING ENERGY OF ",I,":",
03160 & PPHO(4,I),"=>",EN
03161
03162 PPHO(4,I)=EN
03163 E = E+PPHO(4,I)
03164
03165 ENDDO
03166 DO K=1,4
03167 PPHO(K,1)=0d0
03168 DO I=3,NPHO
03169 PPHO(K,1)=PPHO(K,1)+PPHO(K,I)
03170 ENDDO
03171 ENDDO
03172 XMS=SQRT(PPHO(4,1)**2-PPHO(3,1)**2-PPHO(2,1)**2-PPHO(1,1)**2)
03173 PPHO(5,1)=XMS
03174 ELSEIF(MODOP.EQ.3) THEN
03175
03176
03177
03178
03179
03180 DO I=3,NPHO
03181
03182 PX=PX+PPHO(1,I)
03183 PY=PY+PPHO(2,I)
03184 PZ=PZ+PPHO(3,I)
03185 E = E+PPHO(4,I)
03186
03187 P2=PPHO(1,I)**2+PPHO(2,I)**2+PPHO(3,I)**2
03188
03189 M=SQRT(ABS( PPHO(4,I)**2 - P2))
03190
03191 IF (IPRINT.EQ.1)
03192 & WRITE(PHLUN,*) "CORRECTING MASS OF ",I,":",
03193 & PPHO(5,I),"=>",M
03194
03195 PPHO(5,I)=M
03196
03197 ENDDO
03198
03199
03200 ELSEIF(MODOP.EQ.4) THEN
03201
03202
03203
03204
03205
03206 DO I=3,NPHO
03207
03208 PX=PX+PPHO(1,I)
03209 PY=PY+PPHO(2,I)
03210 PZ=PZ+PPHO(3,I)
03211
03212 P2=PPHO(1,I)**2+PPHO(2,I)**2+PPHO(3,I)**2
03213
03214 M=SQRT(ABS( PPHO(4,I)**2 - P2))
03215
03216 IF (M.GT.MCUT) THEN
03217 IF (IPRINT.EQ.1)
03218 & WRITE(PHLUN,*) "CORRECTING MASS OF ",I,":",
03219 & PPHO(5,I),"=>",M
03220 PPHO(5,I)=M
03221 E = E+PPHO(4,I)
03222 ELSE
03223
03224 EN=SQRT( PPHO(5,I)**2 + P2)
03225
03226 IF (IPRINT.EQ.1)
03227 & WRITE(PHLUN,*) "CORRECTING ENERGY OF ",I,":",
03228 & PPHO(4,I),"=>",EN
03229
03230 PPHO(4,I)=EN
03231 E = E+PPHO(4,I)
03232 ENDIF
03233
03234 ENDDO
03235 ENDIF
03236
03237
03238 IF (IPRINT.EQ.1) THEN
03239 WRITE(PHLUN,*) "CORRECTING MOTHER"
03240 WRITE(PHLUN,*) "PX:",PPHO(1,1),"=>",PX-PPHO(1,2)
03241 WRITE(PHLUN,*) "PY:",PPHO(2,1),"=>",PY-PPHO(2,2)
03242 WRITE(PHLUN,*) "PZ:",PPHO(3,1),"=>",PZ-PPHO(3,2)
03243 WRITE(PHLUN,*) " E:",PPHO(4,1),"=>",E-PPHO(4,2)
03244 ENDIF
03245
03246 PPHO(1,1)=PX-PPHO(1,2)
03247 PPHO(2,1)=PY-PPHO(2,2)
03248 PPHO(3,1)=PZ-PPHO(3,2)
03249 PPHO(4,1)=E -PPHO(4,2)
03250
03251 P2=PPHO(1,1)**2+PPHO(2,1)**2+PPHO(3,1)**2
03252
03253 IF (PPHO(4,1)**2.GT.P2) THEN
03254 M=SQRT( PPHO(4,1)**2 - P2 )
03255 IF (IPRINT.EQ.1)
03256 & WRITE(PHLUN,*) " M:",PPHO(5,1),"=>",M
03257 PPHO(5,1)=M
03258 ENDIF
03259
03260 CALL PHLUPA(25)
03261
03262 END
03263
03264
03265
03266 FUNCTION PHINT(IDUM)
03267
03268
03269
03270
03271
03272
03273
03274
03275
03276
03277
03278
03279
03280
03281
03282
03283 IMPLICIT NONE
03284 REAL*8 PHINT,PHINT2
03285 INTEGER IDUM
03286 INTEGER NMXPHO
03287 PARAMETER (NMXPHO=10000)
03288 INTEGER IDPHO,ISTPHO,JDAPHO,JMOPHO,NEVPHO,NPHO
03289 REAL*8 PPHO,VPHO
03290 COMMON/PHOEVT/NEVPHO,NPHO,ISTPHO(NMXPHO),IDPHO(NMXPHO),
03291 &JMOPHO(2,NMXPHO),JDAPHO(2,NMXPHO),PPHO(5,NMXPHO),VPHO(4,NMXPHO)
03292 INTEGER I,K,L
03293 DOUBLE PRECISION EMU,MCHREN,BETA,COSTHG,MPASQR,XPH, XC1, XC2,XDENO
03294 DOUBLE PRECISION XNUM1,XNUM2
03295 DOUBLE PRECISION EPS1(4),EPS2(4),PH(4),PL(4)
03296 REAL*8 PHOCHA
03297
03298
03299
03300
03301 DO K=1,4
03302 PH(K)=PPHO(K,NPHO)
03303 EPS2(K)=1D0
03304 ENDDO
03305
03306 CALL PHOEPS(PH,EPS2,EPS1)
03307 CALL PHOEPS(PH,EPS1,EPS2)
03308
03309
03310 XNUM1=0D0
03311 XNUM2=0D0
03312 XDENO=0D0
03313
03314 DO K=JDAPHO(1,1),NPHO-1
03315
03316
03317 DO L=1,4
03318 PL(L)=PPHO(L,K)
03319 ENDDO
03320
03321
03322 XC1 = - PHOCHA(IDPHO(K)) *
03323 & ( PL(1)*EPS1(1) + PL(2)*EPS1(2) + PL(3)*EPS1(3) ) /
03324 & ( PH(4)*PL(4) - PH(1)*PL(1) - PH(2)*PL(2) - PH(3)*PL(3) )
03325
03326 XC2 = - PHOCHA(IDPHO(K)) *
03327 & ( PL(1)*EPS2(1) + PL(2)*EPS2(2) + PL(3)*EPS2(3) ) /
03328 & ( PH(4)*PL(4) - PH(1)*PL(1) - PH(2)*PL(2) - PH(3)*PL(3) )
03329
03330
03331
03332 XNUM1 = XNUM1+XC1
03333 XNUM2 = XNUM2+XC2
03334
03335 XDENO = XDENO + XC1**2 + XC2**2
03336
03337 ENDDO
03338
03339 PHINT=(XNUM1**2 + XNUM2**2) / XDENO
03340 PHINT2=PHINT
03341
03342 END
03343
03344
03345 SUBROUTINE PHOEPS (VEC1, VEC2, EPS)
03346
03347
03348
03349
03350
03351
03352
03353
03354
03355
03356
03357
03358
03359
03360
03361
03362
03363
03364 DOUBLE PRECISION VEC1(4), VEC2(4), EPS(4),XN
03365
03366 EPS(1)=VEC1(2)*VEC2(3) - VEC1(3)*VEC2(2)
03367 EPS(2)=VEC1(3)*VEC2(1) - VEC1(1)*VEC2(3)
03368 EPS(3)=VEC1(1)*VEC2(2) - VEC1(2)*VEC2(1)
03369 EPS(4)=0D0
03370
03371 XN=SQRT( EPS(1)**2 +EPS(2)**2 +EPS(3)**2)
03372
03373 EPS(1)=EPS(1)/XN
03374 EPS(2)=EPS(2)/XN
03375 EPS(3)=EPS(3)/XN
03376
03377
03378 END
03379 SUBROUTINE PHODMP
03380
03381
03382
03383
03384
03385
03386
03387
03388
03389
03390
03391
03392
03393
03394
03395 DOUBLE PRECISION SUMVEC(5)
03396 INTEGER I,J
03397
03398 INTEGER NMXHEP
03399 PARAMETER (NMXHEP=10000)
03400 REAL*8 phep, vhep
03401 INTEGER nevhep,nhep,isthep,idhep,jmohep,
03402 $ jdahep
03403 COMMON /hepevt/
03404 $ nevhep,
03405 $ nhep,
03406 $ isthep(nmxhep),
03407 $ idhep(nmxhep),
03408 $ jmohep(2,nmxhep),
03409 $ jdahep(2,nmxhep),
03410 $ phep(5,nmxhep),
03411 $ vhep(4,nmxhep)
03412
03413 LOGICAL qedrad
03414 COMMON /phoqed/
03415 $ qedrad(nmxhep)
03416
03417 SAVE hepevt,phoqed
03418 INTEGER PHLUN
03419 COMMON/PHOLUN/PHLUN
03420 DO 10 I=1,5
03421 10 SUMVEC(I)=0.
03422
03423
03424 WRITE(PHLUN,9000)
03425 WRITE(PHLUN,9010) NEVHEP
03426 WRITE(PHLUN,9080)
03427 WRITE(PHLUN,9020)
03428 DO 30 I=1,NHEP
03429
03430
03431 IF (JDAHEP(1,I).EQ.0) THEN
03432 DO 20 J=1,4
03433 20 SUMVEC(J)=SUMVEC(J)+PHEP(J,I)
03434 IF (JMOHEP(2,I).EQ.0) THEN
03435 WRITE(PHLUN,9030) I,IDHEP(I),JMOHEP(1,I),(PHEP(J,I),J=1,5)
03436 ELSE
03437 WRITE(PHLUN,9040) I,IDHEP(I),JMOHEP(1,I),JMOHEP(2,I),(PHEP
03438 & (J,I),J=1,5)
03439 ENDIF
03440 ELSE
03441 IF (JMOHEP(2,I).EQ.0) THEN
03442 WRITE(PHLUN,9050) I,IDHEP(I),JMOHEP(1,I),JDAHEP(1,I),
03443 & JDAHEP(2,I),(PHEP(J,I),J=1,5)
03444 ELSE
03445 WRITE(PHLUN,9060) I,IDHEP(I),JMOHEP(1,I),JMOHEP(2,I),
03446 & JDAHEP(1,I),JDAHEP(2,I),(PHEP(J,I),J=1,5)
03447 ENDIF
03448 ENDIF
03449 30 CONTINUE
03450 SUMVEC(5)=SQRT(SUMVEC(4)**2-SUMVEC(1)**2-SUMVEC(2)**2-
03451 &SUMVEC(3)**2)
03452 WRITE(PHLUN,9070) (SUMVEC(J),J=1,5)
03453 RETURN
03454 9000 FORMAT(1H0,80('='))
03455 9010 FORMAT(1H ,29X,'Event No.:',I10)
03456 9020 FORMAT(1H0,1X,'Nr',3X,'Type',3X,'Parent(s)',2X,'Daughter(s)',6X,
03457 &'Px',7X,'Py',7X,'Pz',7X,'E',4X,'Inv. M.')
03458 9030 FORMAT(1H ,I4,I7,3X,I4,9X,'Stable',2X,5F9.2)
03459 9040 FORMAT(1H ,I4,I7,I4,' - ',I4,5X,'Stable',2X,5F9.2)
03460 9050 FORMAT(1H ,I4,I7,3X,I4,6X,I4,' - ',I4,5F9.2)
03461 9060 FORMAT(1H ,I4,I7,I4,' - ',I4,2X,I4,' - ',I4,5F9.2)
03462 9070 FORMAT(1H0,23X,'Vector Sum: ', 5F9.2)
03463 9080 FORMAT(1H0,6X,'Particle Parameters')
03464 END