hijing/src/pyhiinki.f

0001
0002 C*********************************************************************
0003
0004       SUBROUTINE PYHIINKI(CHFRAM,CHBEAM,CHTARG,WIN)
0005
0006 C...Identifies the two incoming particles and sets up kinematics,
0007 C...including rotations and boosts to/from CM frame.
0008       COMMON/LUJETS/N,K(9000,5),P(9000,5),V(9000,5)
0009       SAVE /LUJETS/
0010       COMMON/LUDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
0011       SAVE /LUDAT1/
0012       COMMON/PYHISUBS/MSEL,MSUB(200),KFIN(2,-40:40),CKIN(200)
0013       SAVE /PYHISUBS/
0014       COMMON/PYHIPARS/MSTP(200),PARP(200),MSTI(200),PARI(200)
0015       SAVE /PYHIPARS/
0016       COMMON/PYHIINT1/MINT(400),VINT(400)
0017       SAVE /PYHIINT1/
0018       CHARACTER CHFRAM*8,CHBEAM*8,CHTARG*8,CHCOM(3)*8,CHALP(2)*26,
0019      &CHIDNT(3)*8,CHTEMP*8,CHCDE(18)*8,CHINIT*76
0020       DIMENSION LEN(3),KCDE(18)
0021       DATA CHALP/'abcdefghijklmnopqrstuvwxyz',
0022      &'ABCDEFGHIJKLMNOPQRSTUVWXYZ'/
0023       DATA CHCDE/'e-      ','e+      ','nue     ','nue~    ',
0024      &'mu-     ','mu+     ','numu    ','numu~   ','tau-    ',
0025      &'tau+    ','nutau   ','nutau~  ','pi+     ','pi-     ',
0026      &'n       ','n~      ','p       ','p~      '/
0027       DATA KCDE/11,-11,12,-12,13,-13,14,-14,15,-15,16,-16,
0028      &211,-211,2112,-2112,2212,-2212/
0029
0030 C...Convert character variables to lowercase and find their length.
0031       CHCOM(1)=CHFRAM
0032       CHCOM(2)=CHBEAM
0033       CHCOM(3)=CHTARG
0034       DO 120 I=1,3
0035       LEN(I)=8
0036       DO 100 LL=8,1,-1
0037       IF(LEN(I).EQ.LL.AND.CHCOM(I)(LL:LL).EQ.' ') LEN(I)=LL-1
0038       DO 100 LA=1,26
0039   100 IF(CHCOM(I)(LL:LL).EQ.CHALP(2)(LA:LA)) CHCOM(I)(LL:LL)=
0040      &CHALP(1)(LA:LA)
0041       CHIDNT(I)=CHCOM(I)
0042       DO 110 LL=1,6
0043       IF(CHIDNT(I)(LL:LL+2).EQ.'bar') THEN
0044         CHTEMP=CHIDNT(I)
0045         CHIDNT(I)=CHTEMP(1:LL-1)//'~'//CHTEMP(LL+3:8)//'  '
0046       ENDIF
0047   110 CONTINUE
0048       DO 120 LL=1,8
0049       IF(CHIDNT(I)(LL:LL).EQ.'_') THEN
0050         CHTEMP=CHIDNT(I)
0051         CHIDNT(I)=CHTEMP(1:LL-1)//CHTEMP(LL+1:8)//' '
0052       ENDIF
0053   120 CONTINUE
0054
0055 C...Set initial state. Error for unknown codes. Reset variables.
0056       N=2
0057       DO 140 I=1,2
0058       K(I,2)=0
0059       DO 130 J=1,18
0060   130 IF(CHIDNT(I+1).EQ.CHCDE(J)) K(I,2)=KCDE(J)
0061       P(I,5)=ULMASS(K(I,2))
0062       MINT(40+I)=1
0063       IF(IABS(K(I,2)).GT.100) MINT(40+I)=2
0064       DO 140 J=1,5
0065   140 V(I,J)=0.
0066       IF(K(1,2).EQ.0) WRITE(MSTU(11),1000) CHBEAM(1:LEN(2))
0067       IF(K(2,2).EQ.0) WRITE(MSTU(11),1100) CHTARG(1:LEN(3))
0068       IF(K(1,2).EQ.0.OR.K(2,2).EQ.0) STOP
0069       DO 150 J=6,10
0070   150 VINT(J)=0.
0071       CHINIT=' '
0072
0073 C...Set up kinematics for events defined in CM frame.
0074       IF(CHCOM(1)(1:2).EQ.'cm') THEN
0075         IF(CHCOM(2)(1:1).NE.'e') THEN
0076           LOFFS=(34-(LEN(2)+LEN(3)))/2
0077           CHINIT(LOFFS+1:76)='PYTHIA will be initialized for a '//
0078      &    CHCOM(2)(1:LEN(2))//'-'//CHCOM(3)(1:LEN(3))//' collider'//' '
0079         ELSE
0080           LOFFS=(33-(LEN(2)+LEN(3)))/2
0081           CHINIT(LOFFS+1:76)='PYTHIA will be initialized for an '//
0082      &    CHCOM(2)(1:LEN(2))//'-'//CHCOM(3)(1:LEN(3))//' collider'//' '
0083         ENDIF
0084 C        WRITE(MSTU(11),1200) CHINIT
0085 C        WRITE(MSTU(11),1300) WIN
0086         S=WIN**2
0087         P(1,1)=0.
0088         P(1,2)=0.
0089         P(2,1)=0.
0090         P(2,2)=0.
0091         P(1,3)=SQRT(((S-P(1,5)**2-P(2,5)**2)**2-(2.*P(1,5)*P(2,5))**2)/
0092      &  (4.*S))
0093         P(2,3)=-P(1,3)
0094         P(1,4)=SQRT(P(1,3)**2+P(1,5)**2)
0095         P(2,4)=SQRT(P(2,3)**2+P(2,5)**2)
0096
0097 C...Set up kinematics for fixed target events.
0098       ELSEIF(CHCOM(1)(1:3).EQ.'fix') THEN
0099         LOFFS=(29-(LEN(2)+LEN(3)))/2
0100         CHINIT(LOFFS+1:76)='PYTHIA will be initialized for '//
0101      &  CHCOM(2)(1:LEN(2))//' on '//CHCOM(3)(1:LEN(3))//
0102      &  ' fixed target'//' '
0103 C        WRITE(MSTU(11),1200) CHINIT
0104 C        WRITE(MSTU(11),1400) WIN
0105         P(1,1)=0.
0106         P(1,2)=0.
0107         P(2,1)=0.
0108         P(2,2)=0.
0109         P(1,3)=WIN
0110         P(1,4)=SQRT(P(1,3)**2+P(1,5)**2)
0111         P(2,3)=0.
0112         P(2,4)=P(2,5)
0113         S=P(1,5)**2+P(2,5)**2+2.*P(2,4)*P(1,4)
0114         VINT(10)=P(1,3)/(P(1,4)+P(2,4))
0115         CALL LUROBO(0.,0.,0.,0.,-VINT(10))
0116 C        WRITE(MSTU(11),1500) SQRT(S)
0117
0118 C...Set up kinematics for events in user-defined frame.
0119       ELSEIF(CHCOM(1)(1:3).EQ.'use') THEN
0120         LOFFS=(13-(LEN(1)+LEN(2)))/2
0121         CHINIT(LOFFS+1:76)='PYTHIA will be initialized for '//
0122      &  CHCOM(2)(1:LEN(2))//' on '//CHCOM(3)(1:LEN(3))//
0123      &  'user-specified configuration'//' '
0124 C        WRITE(MSTU(11),1200) CHINIT
0125 C        WRITE(MSTU(11),1600)
0126 C        WRITE(MSTU(11),1700) CHCOM(2),P(1,1),P(1,2),P(1,3)
0127 C        WRITE(MSTU(11),1700) CHCOM(3),P(2,1),P(2,2),P(2,3)
0128         P(1,4)=SQRT(P(1,1)**2+P(1,2)**2+P(1,3)**2+P(1,5)**2)
0129         P(2,4)=SQRT(P(2,1)**2+P(2,2)**2+P(2,3)**2+P(2,5)**2)
0130         DO 160 J=1,3
0131   160   VINT(7+J)=(DBLE(P(1,J))+DBLE(P(2,J)))/DBLE(P(1,4)+P(2,4))
0132         CALL LUROBO(0.,0.,-VINT(8),-VINT(9),-VINT(10))
0133         VINT(7)=ULANGL(P(1,1),P(1,2))
0134         CALL LUROBO(0.,-VINT(7),0.,0.,0.)
0135         VINT(6)=ULANGL(P(1,3),P(1,1))
0136         CALL LUROBO(-VINT(6),0.,0.,0.,0.)
0137         S=P(1,5)**2+P(2,5)**2+2.*(P(1,4)*P(2,4)-P(1,3)*P(2,3))
0138 C        WRITE(MSTU(11),1500) SQRT(S)
0139
0140 C...Unknown frame. Error for too low CM energy.
0141       ELSE
0142         WRITE(MSTU(11),1800) CHFRAM(1:LEN(1))
0143         STOP
0144       ENDIF
0145       IF(S.LT.PARP(2)**2) THEN
0146         WRITE(MSTU(11),1900) SQRT(S)
0147         STOP
0148       ENDIF
0149
0150 C...Save information on incoming particles.
0151       MINT(11)=K(1,2)
0152       MINT(12)=K(2,2)
0153       MINT(43)=2*MINT(41)+MINT(42)-2
0154       VINT(1)=SQRT(S)
0155       VINT(2)=S
0156       VINT(3)=P(1,5)
0157       VINT(4)=P(2,5)
0158       VINT(5)=P(1,3)
0159
0160 C...Store constants to be used in generation.
0161       IF(MSTP(82).LE.1) VINT(149)=4.*PARP(81)**2/S
0162       IF(MSTP(82).GE.2) VINT(149)=4.*PARP(82)**2/S
0163
0164 C...Formats for initialization and error information.
0165  1000 FORMAT(1X,'Error: unrecognized beam particle ''',A,'''.'/
0166      &1X,'Execution stopped!')
0167  1100 FORMAT(1X,'Error: unrecognized target particle ''',A,'''.'/
0168      &1X,'Execution stopped!')
0169  1200 FORMAT(/1X,78('=')/1X,'I',76X,'I'/1X,'I',A76,'I')
0170  1300 FORMAT(1X,'I',18X,'at',1X,F10.3,1X,'GeV center-of-mass energy',
0171      &19X,'I'/1X,'I',76X,'I'/1X,78('='))
0172  1400 FORMAT(1X,'I',22X,'at',1X,F10.3,1X,'GeV/c lab-momentum',22X,'I')
0173  1500 FORMAT(1X,'I',76X,'I'/1X,'I',11X,'corresponding to',1X,F10.3,1X,
0174      &'GeV center-of-mass energy',12X,'I'/1X,'I',76X,'I'/1X,78('='))
0175  1600 FORMAT(1X,'I',76X,'I'/1X,'I',24X,'px (GeV/c)',3X,'py (GeV/c)',3X,
0176      &'pz (GeV/c)',16X,'I')
0177  1700 FORMAT(1X,'I',15X,A8,3(2X,F10.3,1X),15X,'I')
0178  1800 FORMAT(1X,'Error: unrecognized coordinate frame ''',A,'''.'/
0179      &1X,'Execution stopped!')
0180  1900 FORMAT(1X,'Error: too low CM energy,',F8.3,' GeV for event ',
0181      &'generation.'/1X,'Execution stopped!')
0182
0183       RETURN
0184       END