[u/mrichter/AliRoot.git] / HIJING / hipyset1_35 / pyinit_hijing.F

* $Id$
      SUBROUTINE PYINIT_HIJING(FRAME,BEAM,TARGET,WIN)  
    
C...Initializes the generation procedure; finds maxima of the   
C...differential cross-sections to be used for weighting.   
#include "ludat1_hijing.inc"
#include "ludat2_hijing.inc"
#include "ludat3_hijing.inc"
#include "ludat4_hijing.inc"
#include "pysubs_hijing.inc"
#include "pypars_hijing.inc"
#include "pyint1_hijing.inc"
#include "pyint2_hijing.inc"
#include "pyint5_hijing.inc"
      CHARACTER*(*) FRAME,BEAM,TARGET   
      CHARACTER CHFRAM*8,CHBEAM*8,CHTARG*8,CHMO(12)*3,CHLH(2)*6 
      DATA CHMO/'Jan','Feb','Mar','Apr','May','Jun','Jul','Aug','Sep',  
     &'Oct','Nov','Dec'/, CHLH/'lepton','hadron'/   
    
C...Write headers.  
C      IF(MSTP(122).GE.1) WRITE(MSTU(11),1000) MSTP(181),MSTP(182),  
C     &MSTP(185),CHMO(MSTP(184)),MSTP(183)   
      CALL LULIST_HIJING(0)
C      IF(MSTP(122).GE.1) WRITE(MSTU(11),1100)  
    
C...Identify beam and target particles and initialize kinematics.   
      CHFRAM=FRAME//' ' 
      CHBEAM=BEAM//' '  
      CHTARG=TARGET//' '    
      CALL PYINKI_HIJING(CHFRAM,CHBEAM,CHTARG,WIN) 
    
C...Select partonic subprocesses to be included in the simulation.  
      IF(MSEL.NE.0) THEN    
        DO 100 I=1,200  
  100   MSUB(I)=0   
      ENDIF 
      IF(MINT(43).EQ.1.AND.(MSEL.EQ.1.OR.MSEL.EQ.2)) THEN   
C...Lepton+lepton -> gamma/Z0 or W. 
        IF(MINT(11)+MINT(12).EQ.0) MSUB(1)=1    
        IF(MINT(11)+MINT(12).NE.0) MSUB(2)=1    
      ELSEIF(MSEL.EQ.1) THEN    
C...High-pT QCD processes:  
        MSUB(11)=1  
        MSUB(12)=1  
        MSUB(13)=1  
        MSUB(28)=1  
        MSUB(53)=1  
        MSUB(68)=1  
        IF(MSTP(82).LE.1.AND.CKIN(3).LT.PARP(81)) MSUB(95)=1    
        IF(MSTP(82).GE.2.AND.CKIN(3).LT.PARP(82)) MSUB(95)=1    
      ELSEIF(MSEL.EQ.2) THEN    
C...All QCD processes:  
        MSUB(11)=1  
        MSUB(12)=1  
        MSUB(13)=1  
        MSUB(28)=1  
        MSUB(53)=1  
        MSUB(68)=1  
        MSUB(91)=1  
        MSUB(92)=1  
        MSUB(93)=1  
        MSUB(95)=1  
      ELSEIF(MSEL.GE.4.AND.MSEL.LE.8) THEN  
C...Heavy quark production. 
        MSUB(81)=1  
        MSUB(82)=1  
        DO 110 J=1,MIN(8,MDCY(21,3))    
  110   MDME(MDCY(21,2)+J-1,1)=0    
        MDME(MDCY(21,2)+MSEL-1,1)=1 
      ELSEIF(MSEL.EQ.10) THEN   
C...Prompt photon production:   
        MSUB(14)=1  
        MSUB(18)=1  
        MSUB(29)=1  
      ELSEIF(MSEL.EQ.11) THEN   
C...Z0/gamma* production:   
        MSUB(1)=1   
      ELSEIF(MSEL.EQ.12) THEN   
C...W+/- production:    
        MSUB(2)=1   
      ELSEIF(MSEL.EQ.13) THEN   
C...Z0 + jet:   
        MSUB(15)=1  
        MSUB(30)=1  
      ELSEIF(MSEL.EQ.14) THEN   
C...W+/- + jet: 
        MSUB(16)=1  
        MSUB(31)=1  
      ELSEIF(MSEL.EQ.15) THEN   
C...Z0 & W+/- pair production:  
        MSUB(19)=1  
        MSUB(20)=1  
        MSUB(22)=1  
        MSUB(23)=1  
        MSUB(25)=1  
      ELSEIF(MSEL.EQ.16) THEN   
C...H0 production:  
        MSUB(3)=1   
        MSUB(5)=1   
        MSUB(8)=1   
        MSUB(102)=1 
      ELSEIF(MSEL.EQ.17) THEN   
C...H0 & Z0 or W+/- pair production:    
        MSUB(24)=1  
        MSUB(26)=1  
      ELSEIF(MSEL.EQ.21) THEN   
C...Z'0 production: 
        MSUB(141)=1 
      ELSEIF(MSEL.EQ.22) THEN   
C...H+/- production:    
        MSUB(142)=1 
      ELSEIF(MSEL.EQ.23) THEN   
C...R production:   
        MSUB(143)=1 
      ENDIF 
    
C...Count number of subprocesses on.    
      MINT(44)=0    
      DO 120 ISUB=1,200 
      IF(MINT(43).LT.4.AND.ISUB.GE.91.AND.ISUB.LE.96.AND.   
     &MSUB(ISUB).EQ.1) THEN 
        WRITE(MSTU(11),1200) ISUB,CHLH(MINT(41)),CHLH(MINT(42)) 
        STOP    
      ELSEIF(MSUB(ISUB).EQ.1.AND.ISET(ISUB).EQ.-1) THEN 
        WRITE(MSTU(11),1300) ISUB   
        STOP    
      ELSEIF(MSUB(ISUB).EQ.1.AND.ISET(ISUB).LE.-2) THEN 
        WRITE(MSTU(11),1400) ISUB   
        STOP    
      ELSEIF(MSUB(ISUB).EQ.1) THEN  
        MINT(44)=MINT(44)+1 
      ENDIF 
  120 CONTINUE  
      IF(MINT(44).EQ.0) THEN    
        WRITE(MSTU(11),1500)    
        STOP    
      ENDIF 
      MINT(45)=MINT(44)-MSUB(91)-MSUB(92)-MSUB(93)-MSUB(94) 
    
C...Maximum 4 generations; set maximum number of allowed flavours.  
      MSTP(1)=MIN(4,MSTP(1))    
      MSTU(114)=MIN(MSTU(114),2*MSTP(1))    
      MSTP(54)=MIN(MSTP(54),2*MSTP(1))  
    
C...Sum up Cabibbo-Kobayashi-Maskawa factors for each quark/lepton. 
      DO 140 I=-20,20   
      VINT(180+I)=0.    
      IA=IABS(I)    
      IF(IA.GE.1.AND.IA.LE.2*MSTP(1)) THEN  
        DO 130 J=1,MSTP(1)  
        IB=2*J-1+MOD(IA,2)  
        IPM=(5-ISIGN(1,I))/2    
        IDC=J+MDCY(IA,2)+2  
  130   IF(MDME(IDC,1).EQ.1.OR.MDME(IDC,1).EQ.IPM) VINT(180+I)= 
     &  VINT(180+I)+VCKM((IA+1)/2,(IB+1)/2) 
      ELSEIF(IA.GE.11.AND.IA.LE.10+2*MSTP(1)) THEN  
        VINT(180+I)=1.  
      ENDIF 
  140 CONTINUE  
    
C...Choose Lambda value to use in alpha-strong. 
      MSTU(111)=MSTP(2) 
      IF(MSTP(3).GE.1) THEN 
        ALAM=PARP(1)    
        IF(MSTP(51).EQ.1) ALAM=0.2  
        IF(MSTP(51).EQ.2) ALAM=0.29 
        IF(MSTP(51).EQ.3) ALAM=0.2  
        IF(MSTP(51).EQ.4) ALAM=0.4  
        IF(MSTP(51).EQ.11) ALAM=0.16    
        IF(MSTP(51).EQ.12) ALAM=0.26    
        IF(MSTP(51).EQ.13) ALAM=0.36    
        PARP(1)=ALAM    
        PARP(61)=ALAM   
        PARU(112)=ALAM  
        PARJ(81)=ALAM   
      ENDIF 
    
C...Initialize widths and partial widths for resonances.    
      CALL PYINRE_HIJING   
    
C...Reset variables for cross-section calculation.  
      DO 150 I=0,200    
      DO 150 J=1,3  
      NGEN(I,J)=0   
  150 XSEC(I,J)=0.  
      VINT(108)=0.  
    
C...Find parametrized total cross-sections. 
      IF(MINT(43).EQ.4) CALL PYXTOT_HIJING 
    
C...Maxima of differential cross-sections.  
      IF(MSTP(121).LE.0) CALL PYMAXI_HIJING    
    
C...Initialize possibility of overlayed events. 
      IF(MSTP(131).NE.0) CALL PYOVLY_HIJING(1) 
    
C...Initialize multiple interactions with variable impact parameter.    
      IF(MINT(43).EQ.4.AND.(MINT(45).NE.0.OR.MSTP(131).NE.0).AND.   
     &MSTP(82).GE.2) CALL PYMULT_HIJING(1) 
C      IF(MSTP(122).GE.1) WRITE(MSTU(11),1600)  
    
C...Formats for initialization information. 
 1000 FORMAT(///20X,'The Lund Monte Carlo - PYTHIA version ',I1,'.',I1/ 
     &20X,'**  Last date of change:  ',I2,1X,A3,1X,I4,'  **'/)  
 1100 FORMAT('1',18('*'),1X,'PYINIT_HIJING: initialization of PYTHIA ',    
     &'routines',1X,17('*'))    
 1200 FORMAT(1X,'Error: process number ',I3,' not meaningful for ',A6,  
     &'-',A6,' interactions.'/1X,'Execution stopped!')  
 1300 FORMAT(1X,'Error: requested subprocess',I4,' not implemented.'/   
     &1X,'Execution stopped!')  
 1400 FORMAT(1X,'Error: requested subprocess',I4,' not existing.'/  
     &1X,'Execution stopped!')  
 1500 FORMAT(1X,'Error: no subprocess switched on.'/    
     &1X,'Execution stopped.')  
 1600 FORMAT(/1X,22('*'),1X,'PYINIT_HIJING: initialization completed',1X
     $     ,22('*'))  
    
      RETURN    
      END
Commit	Line	Data
e74335a4	1	* $Id$
	2	SUBROUTINE PYINIT_HIJING(FRAME,BEAM,TARGET,WIN)
	3
	4	C...Initializes the generation procedure; finds maxima of the
	5	C...differential cross-sections to be used for weighting.
	6	#include "ludat1_hijing.inc"
	7	#include "ludat2_hijing.inc"
	8	#include "ludat3_hijing.inc"
	9	#include "ludat4_hijing.inc"
	10	#include "pysubs_hijing.inc"
	11	#include "pypars_hijing.inc"
	12	#include "pyint1_hijing.inc"
	13	#include "pyint2_hijing.inc"
	14	#include "pyint5_hijing.inc"
	15	CHARACTER() FRAME,BEAM,TARGET
	16	CHARACTER CHFRAM8,CHBEAM8,CHTARG8,CHMO(12)3,CHLH(2)*6
	17	DATA CHMO/'Jan','Feb','Mar','Apr','May','Jun','Jul','Aug','Sep',
	18	&'Oct','Nov','Dec'/, CHLH/'lepton','hadron'/
	19
	20	C...Write headers.
	21	C IF(MSTP(122).GE.1) WRITE(MSTU(11),1000) MSTP(181),MSTP(182),
	22	C &MSTP(185),CHMO(MSTP(184)),MSTP(183)
	23	CALL LULIST_HIJING(0)
	24	C IF(MSTP(122).GE.1) WRITE(MSTU(11),1100)
	25
	26	C...Identify beam and target particles and initialize kinematics.
	27	CHFRAM=FRAME//' '
	28	CHBEAM=BEAM//' '
	29	CHTARG=TARGET//' '
	30	CALL PYINKI_HIJING(CHFRAM,CHBEAM,CHTARG,WIN)
	31
	32	C...Select partonic subprocesses to be included in the simulation.
	33	IF(MSEL.NE.0) THEN
	34	DO 100 I=1,200
	35	100 MSUB(I)=0
	36	ENDIF
	37	IF(MINT(43).EQ.1.AND.(MSEL.EQ.1.OR.MSEL.EQ.2)) THEN
	38	C...Lepton+lepton -> gamma/Z0 or W.
	39	IF(MINT(11)+MINT(12).EQ.0) MSUB(1)=1
	40	IF(MINT(11)+MINT(12).NE.0) MSUB(2)=1
	41	ELSEIF(MSEL.EQ.1) THEN
	42	C...High-pT QCD processes:
	43	MSUB(11)=1
	44	MSUB(12)=1
	45	MSUB(13)=1
	46	MSUB(28)=1
	47	MSUB(53)=1
	48	MSUB(68)=1
	49	IF(MSTP(82).LE.1.AND.CKIN(3).LT.PARP(81)) MSUB(95)=1
	50	IF(MSTP(82).GE.2.AND.CKIN(3).LT.PARP(82)) MSUB(95)=1
	51	ELSEIF(MSEL.EQ.2) THEN
	52	C...All QCD processes:
	53	MSUB(11)=1
	54	MSUB(12)=1
	55	MSUB(13)=1
	56	MSUB(28)=1
	57	MSUB(53)=1
	58	MSUB(68)=1
	59	MSUB(91)=1
	60	MSUB(92)=1
	61	MSUB(93)=1
	62	MSUB(95)=1
	63	ELSEIF(MSEL.GE.4.AND.MSEL.LE.8) THEN
	64	C...Heavy quark production.
65	MSUB(81)=1
66	MSUB(82)=1
67	DO 110 J=1,MIN(8,MDCY(21,3))
68	110 MDME(MDCY(21,2)+J-1,1)=0
69	MDME(MDCY(21,2)+MSEL-1,1)=1
70	ELSEIF(MSEL.EQ.10) THEN
71	C...Prompt photon production:
72	MSUB(14)=1
73	MSUB(18)=1
74	MSUB(29)=1
75	ELSEIF(MSEL.EQ.11) THEN
76	C...Z0/gamma* production:
77	MSUB(1)=1
78	ELSEIF(MSEL.EQ.12) THEN
79	C...W+/- production:
80	MSUB(2)=1
81	ELSEIF(MSEL.EQ.13) THEN
82	C...Z0 + jet:
83	MSUB(15)=1
84	MSUB(30)=1
85	ELSEIF(MSEL.EQ.14) THEN
86	C...W+/- + jet:
87	MSUB(16)=1
88	MSUB(31)=1
89	ELSEIF(MSEL.EQ.15) THEN
90	C...Z0 & W+/- pair production:
91	MSUB(19)=1
92	MSUB(20)=1
93	MSUB(22)=1
94	MSUB(23)=1
95	MSUB(25)=1
96	ELSEIF(MSEL.EQ.16) THEN
97	C...H0 production:
98	MSUB(3)=1
99	MSUB(5)=1
100	MSUB(8)=1
101	MSUB(102)=1
102	ELSEIF(MSEL.EQ.17) THEN
103	C...H0 & Z0 or W+/- pair production:
104	MSUB(24)=1
105	MSUB(26)=1
106	ELSEIF(MSEL.EQ.21) THEN
107	C...Z'0 production:
108	MSUB(141)=1
109	ELSEIF(MSEL.EQ.22) THEN
110	C...H+/- production:
111	MSUB(142)=1
112	ELSEIF(MSEL.EQ.23) THEN
113	C...R production:
114	MSUB(143)=1
115	ENDIF
116
117	C...Count number of subprocesses on.
118	MINT(44)=0
119	DO 120 ISUB=1,200
120	IF(MINT(43).LT.4.AND.ISUB.GE.91.AND.ISUB.LE.96.AND.
121	&MSUB(ISUB).EQ.1) THEN
122	WRITE(MSTU(11),1200) ISUB,CHLH(MINT(41)),CHLH(MINT(42))
123	STOP
124	ELSEIF(MSUB(ISUB).EQ.1.AND.ISET(ISUB).EQ.-1) THEN
125	WRITE(MSTU(11),1300) ISUB
126	STOP
127	ELSEIF(MSUB(ISUB).EQ.1.AND.ISET(ISUB).LE.-2) THEN
128	WRITE(MSTU(11),1400) ISUB
129	STOP
130	ELSEIF(MSUB(ISUB).EQ.1) THEN
131	MINT(44)=MINT(44)+1
132	ENDIF
133	120 CONTINUE
134	IF(MINT(44).EQ.0) THEN
135	WRITE(MSTU(11),1500)
136	STOP
137	ENDIF
138	MINT(45)=MINT(44)-MSUB(91)-MSUB(92)-MSUB(93)-MSUB(94)
139
140	C...Maximum 4 generations; set maximum number of allowed flavours.
141	MSTP(1)=MIN(4,MSTP(1))
142	MSTU(114)=MIN(MSTU(114),2*MSTP(1))
143	MSTP(54)=MIN(MSTP(54),2*MSTP(1))
144
145	C...Sum up Cabibbo-Kobayashi-Maskawa factors for each quark/lepton.
146	DO 140 I=-20,20
147	VINT(180+I)=0.
148	IA=IABS(I)
149	IF(IA.GE.1.AND.IA.LE.2*MSTP(1)) THEN
150	DO 130 J=1,MSTP(1)
151	IB=2*J-1+MOD(IA,2)
152	IPM=(5-ISIGN(1,I))/2
153	IDC=J+MDCY(IA,2)+2
154	130 IF(MDME(IDC,1).EQ.1.OR.MDME(IDC,1).EQ.IPM) VINT(180+I)=
155	& VINT(180+I)+VCKM((IA+1)/2,(IB+1)/2)
156	ELSEIF(IA.GE.11.AND.IA.LE.10+2*MSTP(1)) THEN
157	VINT(180+I)=1.
158	ENDIF
159	140 CONTINUE
160
161	C...Choose Lambda value to use in alpha-strong.
162	MSTU(111)=MSTP(2)
163	IF(MSTP(3).GE.1) THEN
164	ALAM=PARP(1)
165	IF(MSTP(51).EQ.1) ALAM=0.2
166	IF(MSTP(51).EQ.2) ALAM=0.29
167	IF(MSTP(51).EQ.3) ALAM=0.2
168	IF(MSTP(51).EQ.4) ALAM=0.4
169	IF(MSTP(51).EQ.11) ALAM=0.16
170	IF(MSTP(51).EQ.12) ALAM=0.26
171	IF(MSTP(51).EQ.13) ALAM=0.36
172	PARP(1)=ALAM
173	PARP(61)=ALAM
174	PARU(112)=ALAM
175	PARJ(81)=ALAM
176	ENDIF
177
178	C...Initialize widths and partial widths for resonances.
179	CALL PYINRE_HIJING
180
181	C...Reset variables for cross-section calculation.
182	DO 150 I=0,200
183	DO 150 J=1,3
184	NGEN(I,J)=0
185	150 XSEC(I,J)=0.
186	VINT(108)=0.
187
188	C...Find parametrized total cross-sections.
189	IF(MINT(43).EQ.4) CALL PYXTOT_HIJING
190
191	C...Maxima of differential cross-sections.
192	IF(MSTP(121).LE.0) CALL PYMAXI_HIJING
193
194	C...Initialize possibility of overlayed events.
195	IF(MSTP(131).NE.0) CALL PYOVLY_HIJING(1)
196
197	C...Initialize multiple interactions with variable impact parameter.
198	IF(MINT(43).EQ.4.AND.(MINT(45).NE.0.OR.MSTP(131).NE.0).AND.
199	&MSTP(82).GE.2) CALL PYMULT_HIJING(1)
200	C IF(MSTP(122).GE.1) WRITE(MSTU(11),1600)
201
202	C...Formats for initialization information.
203	1000 FORMAT(///20X,'The Lund Monte Carlo - PYTHIA version ',I1,'.',I1/
204	&20X,' Last date of change: ',I2,1X,A3,1X,I4,' '/)
205	1100 FORMAT('1',18('*'),1X,'PYINIT_HIJING: initialization of PYTHIA ',
206	&'routines',1X,17('*'))
207	1200 FORMAT(1X,'Error: process number ',I3,' not meaningful for ',A6,
208	&'-',A6,' interactions.'/1X,'Execution stopped!')
209	1300 FORMAT(1X,'Error: requested subprocess',I4,' not implemented.'/
210	&1X,'Execution stopped!')
211	1400 FORMAT(1X,'Error: requested subprocess',I4,' not existing.'/
212	&1X,'Execution stopped!')
213	1500 FORMAT(1X,'Error: no subprocess switched on.'/
214	&1X,'Execution stopped.')
215	1600 FORMAT(/1X,22('*'),1X,'PYINIT_HIJING: initialization completed',1X
216	$ ,22('*'))
217
218	RETURN
219	END