Merging the VirtualMC branch to the main development branch (HEAD)
[u/mrichter/AliRoot.git] / HIJING / hipyset1_35 / pyinki_hijing.F
1 * $Id$
2     
3 C*********************************************************************  
4     
5       SUBROUTINE PYINKI_HIJING(CHFRAM,CHBEAM,CHTARG,WIN)   
6     
7 C...Identifies the two incoming particles and sets up kinematics,   
8 C...including rotations and boosts to/from CM frame.    
9 #include "lujets_hijing.inc"
10 #include "ludat1_hijing.inc"
11 #include "pysubs_hijing.inc"
12 #include "pypars_hijing.inc"
13 #include "pyint1_hijing.inc"
14       CHARACTER CHFRAM*8,CHBEAM*8,CHTARG*8,CHCOM(3)*8,CHALP(2)*26,  
15      &CHIDNT(3)*8,CHTEMP*8,CHCDE(18)*8,CHINIT*76    
16       DIMENSION LEN(3),KCDE(18) 
17       DATA CHALP/'abcdefghijklmnopqrstuvwxyz',  
18      &'ABCDEFGHIJKLMNOPQRSTUVWXYZ'/ 
19       DATA CHCDE/'e-      ','e+      ','nue     ','nue~    ',   
20      &'mu-     ','mu+     ','numu    ','numu~   ','tau-    ',   
21      &'tau+    ','nutau   ','nutau~  ','pi+     ','pi-     ',   
22      &'n       ','n~      ','p       ','p~      '/  
23       DATA KCDE/11,-11,12,-12,13,-13,14,-14,15,-15,16,-16,  
24      &211,-211,2112,-2112,2212,-2212/   
25     
26 C...Convert character variables to lowercase and find their length. 
27       CHCOM(1)=CHFRAM   
28       CHCOM(2)=CHBEAM   
29       CHCOM(3)=CHTARG   
30       DO 120 I=1,3  
31       LEN(I)=8  
32       DO 100 LL=8,1,-1  
33       IF(LEN(I).EQ.LL.AND.CHCOM(I)(LL:LL).EQ.' ') LEN(I)=LL-1   
34       DO 100 LA=1,26    
35   100 IF(CHCOM(I)(LL:LL).EQ.CHALP(2)(LA:LA)) CHCOM(I)(LL:LL)=   
36      &CHALP(1)(LA:LA)   
37       CHIDNT(I)=CHCOM(I)    
38       DO 110 LL=1,6 
39       IF(CHIDNT(I)(LL:LL+2).EQ.'bar') THEN  
40         CHTEMP=CHIDNT(I)    
41         CHIDNT(I)=CHTEMP(1:LL-1)//'~'//CHTEMP(LL+3:8)//'  ' 
42       ENDIF 
43   110 CONTINUE  
44       DO 120 LL=1,8 
45       IF(CHIDNT(I)(LL:LL).EQ.'_') THEN  
46         CHTEMP=CHIDNT(I)    
47         CHIDNT(I)=CHTEMP(1:LL-1)//CHTEMP(LL+1:8)//' '   
48       ENDIF 
49   120 CONTINUE  
50     
51 C...Set initial state. Error for unknown codes. Reset variables.    
52       N=2   
53       DO 140 I=1,2  
54       K(I,2)=0  
55       DO 130 J=1,18 
56   130 IF(CHIDNT(I+1).EQ.CHCDE(J)) K(I,2)=KCDE(J)    
57       P(I,5)=ULMASS_HIJING(K(I,2)) 
58       MINT(40+I)=1  
59       IF(IABS(K(I,2)).GT.100) MINT(40+I)=2  
60       DO 140 J=1,5  
61   140 V(I,J)=0. 
62       IF(K(1,2).EQ.0) WRITE(MSTU(11),1000) CHBEAM(1:LEN(2)) 
63       IF(K(2,2).EQ.0) WRITE(MSTU(11),1100) CHTARG(1:LEN(3)) 
64       IF(K(1,2).EQ.0.OR.K(2,2).EQ.0) STOP   
65       DO 150 J=6,10 
66   150 VINT(J)=0.    
67       CHINIT=' '    
68     
69 C...Set up kinematics for events defined in CM frame.   
70       IF(CHCOM(1)(1:2).EQ.'cm') THEN    
71         IF(CHCOM(2)(1:1).NE.'e') THEN   
72           LOFFS=(34-(LEN(2)+LEN(3)))/2  
73           CHINIT(LOFFS+1:76)='PYTHIA will be initialized for a '//  
74      &    CHCOM(2)(1:LEN(2))//'-'//CHCOM(3)(1:LEN(3))//' collider'//' ' 
75         ELSE    
76           LOFFS=(33-(LEN(2)+LEN(3)))/2  
77           CHINIT(LOFFS+1:76)='PYTHIA will be initialized for an '// 
78      &    CHCOM(2)(1:LEN(2))//'-'//CHCOM(3)(1:LEN(3))//' collider'//' ' 
79         ENDIF   
80 C        WRITE(MSTU(11),1200) CHINIT 
81 C        WRITE(MSTU(11),1300) WIN    
82         S=WIN**2    
83         P(1,1)=0.   
84         P(1,2)=0.   
85         P(2,1)=0.   
86         P(2,2)=0.   
87         P(1,3)=SQRT(((S-P(1,5)**2-P(2,5)**2)**2-(2.*P(1,5)*P(2,5))**2)/ 
88      &  (4.*S)) 
89         P(2,3)=-P(1,3)  
90         P(1,4)=SQRT(P(1,3)**2+P(1,5)**2)    
91         P(2,4)=SQRT(P(2,3)**2+P(2,5)**2)    
92     
93 C...Set up kinematics for fixed target events.  
94       ELSEIF(CHCOM(1)(1:3).EQ.'fix') THEN   
95         LOFFS=(29-(LEN(2)+LEN(3)))/2    
96         CHINIT(LOFFS+1:76)='PYTHIA will be initialized for '//  
97      &  CHCOM(2)(1:LEN(2))//' on '//CHCOM(3)(1:LEN(3))//    
98      &  ' fixed target'//' '    
99 C        WRITE(MSTU(11),1200) CHINIT 
100 C        WRITE(MSTU(11),1400) WIN    
101         P(1,1)=0.   
102         P(1,2)=0.   
103         P(2,1)=0.   
104         P(2,2)=0.   
105         P(1,3)=WIN  
106         P(1,4)=SQRT(P(1,3)**2+P(1,5)**2)    
107         P(2,3)=0.   
108         P(2,4)=P(2,5)   
109         S=P(1,5)**2+P(2,5)**2+2.*P(2,4)*P(1,4)  
110         VINT(10)=P(1,3)/(P(1,4)+P(2,4)) 
111         CALL LUROBO_HIJING(0.,0.,0.,0.,-VINT(10))  
112 C        WRITE(MSTU(11),1500) SQRT(S)    
113     
114 C...Set up kinematics for events in user-defined frame. 
115       ELSEIF(CHCOM(1)(1:3).EQ.'use') THEN   
116         LOFFS=(13-(LEN(1)+LEN(2)))/2    
117         CHINIT(LOFFS+1:76)='PYTHIA will be initialized for '//  
118      &  CHCOM(2)(1:LEN(2))//' on '//CHCOM(3)(1:LEN(3))//    
119      &  'user-specified configuration'//' ' 
120 C        WRITE(MSTU(11),1200) CHINIT 
121 C        WRITE(MSTU(11),1600)    
122 C        WRITE(MSTU(11),1700) CHCOM(2),P(1,1),P(1,2),P(1,3)  
123 C        WRITE(MSTU(11),1700) CHCOM(3),P(2,1),P(2,2),P(2,3)  
124         P(1,4)=SQRT(P(1,1)**2+P(1,2)**2+P(1,3)**2+P(1,5)**2)    
125         P(2,4)=SQRT(P(2,1)**2+P(2,2)**2+P(2,3)**2+P(2,5)**2)    
126         DO 160 J=1,3    
127   160   VINT(7+J)=(DBLE(P(1,J))+DBLE(P(2,J)))/DBLE(P(1,4)+P(2,4))   
128         CALL LUROBO_HIJING(0.,0.,-VINT(8),-VINT(9),-VINT(10))  
129         VINT(7)=ULANGL_HIJING(P(1,1),P(1,2))   
130         CALL LUROBO_HIJING(0.,-VINT(7),0.,0.,0.)   
131         VINT(6)=ULANGL_HIJING(P(1,3),P(1,1))   
132         CALL LUROBO_HIJING(-VINT(6),0.,0.,0.,0.)   
133         S=P(1,5)**2+P(2,5)**2+2.*(P(1,4)*P(2,4)-P(1,3)*P(2,3))  
134 C        WRITE(MSTU(11),1500) SQRT(S)    
135     
136 C...Unknown frame. Error for too low CM energy. 
137       ELSE  
138         WRITE(MSTU(11),1800) CHFRAM(1:LEN(1))   
139         STOP    
140       ENDIF 
141       IF(S.LT.PARP(2)**2) THEN  
142         WRITE(MSTU(11),1900) SQRT(S)    
143         STOP    
144       ENDIF 
145     
146 C...Save information on incoming particles. 
147       MINT(11)=K(1,2)   
148       MINT(12)=K(2,2)   
149       MINT(43)=2*MINT(41)+MINT(42)-2    
150       VINT(1)=SQRT(S)   
151       VINT(2)=S 
152       VINT(3)=P(1,5)    
153       VINT(4)=P(2,5)    
154       VINT(5)=P(1,3)    
155     
156 C...Store constants to be used in generation.   
157       IF(MSTP(82).LE.1) VINT(149)=4.*PARP(81)**2/S  
158       IF(MSTP(82).GE.2) VINT(149)=4.*PARP(82)**2/S  
159     
160 C...Formats for initialization and error information.   
161  1000 FORMAT(1X,'Error: unrecognized beam particle ''',A,'''.'/ 
162      &1X,'Execution stopped!')  
163  1100 FORMAT(1X,'Error: unrecognized target particle ''',A,'''.'/   
164      &1X,'Execution stopped!')  
165  1200 FORMAT(/1X,78('=')/1X,'I',76X,'I'/1X,'I',A76,'I') 
166  1300 FORMAT(1X,'I',18X,'at',1X,F10.3,1X,'GeV center-of-mass energy',   
167      &19X,'I'/1X,'I',76X,'I'/1X,78('='))    
168  1400 FORMAT(1X,'I',22X,'at',1X,F10.3,1X,'GeV/c lab-momentum',22X,'I')  
169  1500 FORMAT(1X,'I',76X,'I'/1X,'I',11X,'corresponding to',1X,F10.3,1X,  
170      &'GeV center-of-mass energy',12X,'I'/1X,'I',76X,'I'/1X,78('='))    
171  1600 FORMAT(1X,'I',76X,'I'/1X,'I',24X,'px (GeV/c)',3X,'py (GeV/c)',3X, 
172      &'pz (GeV/c)',16X,'I') 
173  1700 FORMAT(1X,'I',15X,A8,3(2X,F10.3,1X),15X,'I')  
174  1800 FORMAT(1X,'Error: unrecognized coordinate frame ''',A,'''.'/  
175      &1X,'Execution stopped!')  
176  1900 FORMAT(1X,'Error: too low CM energy,',F8.3,' GeV for event ', 
177      &'generation.'/1X,'Execution stopped!')    
178     
179       RETURN    
180       END