Changes to make code work for pp collisions.
[u/mrichter/AliRoot.git] / DPMJET / dpmjet3.0-4.f
1 *
2 *    +-------------------------------------------------------------+
3 *    |                                                             |
4 *    |                                                             |
5 *    |                        DPMJET 3.0                           |
6 *    |                                                             |
7 *    |                                                             |
8 *    |         S. Roesler+), R. Engel#), J. Ranft*)                |
9 *    |                                                             |
10 *    |         +) CERN, TIS-RP                                     |
11 *    |            CH-1211 Geneva 23, Switzerland                   |
12 *    |            Email: Stefan.Roesler@cern.ch                    |
13 *    |                                                             |
14 *    |         #) University of Delaware, BRI                      |
15 *    |            Newark, DE 19716, USA                            |
16 *    |                                                             |
17 *    |         *) University of Siegen, Dept. of Physics           |
18 *    |            D-57068 Siegen, Germany                          |
19 *    |                                                             |
20 *    |                                                             |
21 *    |       http://home.cern.ch/sroesler/dpmjet3.html             |
22 *    |                                                             |
23 *    |                                                             |
24 *    |       Monte Carlo models used for event generation:         |
25 *    |          PHOJET 1.12, JETSET 7.4 and LEPTO 6.5.1            |
26 *    |                                                             |
27 *    +-------------------------------------------------------------+
28 *
29 *
30 *===init===============================================================*
31 *
32 CDECK  ID>, DT_INIT
33       SUBROUTINE DT_INIT(NCASES,EPN,NPMASS,NPCHAR,NTMASS,NTCHAR,
34      &                                             IDP,IGLAU)
35
36 ************************************************************************
37 * Initialization of event generation                                   *
38 * This version dated  7.4.98  is written by S. Roesler.                *
39 ************************************************************************
40
41       IMPLICIT DOUBLE PRECISION (A-H,O-Z)
42       SAVE
43
44       PARAMETER ( LINP = 5 ,
45      &            LOUT = 6 ,
46      &            LDAT = 9 )
47
48       PARAMETER (ZERO=0.0D0,ONE=1.0D0)
49
50 * particle properties (BAMJET index convention)
51       CHARACTER*8  ANAME
52       COMMON /DTPART/ ANAME(210),AAM(210),GA(210),TAU(210),
53      &                IICH(210),IIBAR(210),K1(210),K2(210)
54 * names of hadrons used in input-cards
55       CHARACTER*8 BTYPE
56       COMMON /DTPAIN/ BTYPE(30)
57
58       INCLUDE '(DIMPAR)'
59       INCLUDE '(PAREVT)'
60       INCLUDE '(EVAPAR)'
61       INCLUDE '(FRBKCM)'
62
63       PARAMETER (NCOMPX=20,NEB=8,NQB= 5,KSITEB=50)
64
65 * emulsion treatment
66       COMMON /DTCOMP/ EMUFRA(NCOMPX),IEMUMA(NCOMPX),IEMUCH(NCOMPX),
67      &                NCOMPO,IEMUL
68 * Glauber formalism: parameters
69       COMMON /DTGLAM/ RASH(NCOMPX),RBSH(NCOMPX),
70      &                BMAX(NCOMPX),BSTEP(NCOMPX),
71      &                SIGSH,ROSH,GSH,BSITE(0:NEB,NQB,NCOMPX,KSITEB),
72      &                NSITEB,NSTATB
73 * Glauber formalism: cross sections
74       COMMON /DTGLXS/ ECMNN(NEB),Q2G(NQB),ECMNOW,Q2,
75      &                XSTOT(NEB,NQB,NCOMPX),XSELA(NEB,NQB,NCOMPX),
76      &                XSQEP(NEB,NQB,NCOMPX),XSQET(NEB,NQB,NCOMPX),
77      &                XSQE2(NEB,NQB,NCOMPX),XSPRO(NEB,NQB,NCOMPX),
78      &                XSDEL(NEB,NQB,NCOMPX),XSDQE(NEB,NQB,NCOMPX),
79      &                XETOT(NEB,NQB,NCOMPX),XEELA(NEB,NQB,NCOMPX),
80      &                XEQEP(NEB,NQB,NCOMPX),XEQET(NEB,NQB,NCOMPX),
81      &                XEQE2(NEB,NQB,NCOMPX),XEPRO(NEB,NQB,NCOMPX),
82      &                XEDEL(NEB,NQB,NCOMPX),XEDQE(NEB,NQB,NCOMPX),
83      &                BSLOPE,NEBINI,NQBINI
84 * interface HADRIN-DPM
85       COMMON /HNTHRE/ EHADTH,EHADLO,EHADHI,INTHAD,IDXTA
86 * central particle production, impact parameter biasing
87       COMMON /DTIMPA/ BIMIN,BIMAX,XSFRAC,ICENTR
88 * parameter for intranuclear cascade
89       LOGICAL LPAULI
90       COMMON /DTFOTI/ TAUFOR,KTAUGE,ITAUVE,INCMOD,LPAULI
91 * various options for treatment of partons (DTUNUC 1.x)
92 * (chain recombination, Cronin,..)
93       LOGICAL LCO2CR,LINTPT
94       COMMON /DTCHAI/ SEASQ,CRONCO,CUTOF,MKCRON,ISICHA,IRECOM,
95      &                LCO2CR,LINTPT
96 * threshold values for x-sampling (DTUNUC 1.x)
97       COMMON /DTXCUT/ XSEACU,UNON,UNOM,UNOSEA,CVQ,CDQ,CSEA,SSMIMA,
98      &                SSMIMQ,VVMTHR
99 * flags for input different options
100       LOGICAL LEMCCK,LHADRO,LSEADI,LEVAPO
101       COMMON /DTFLG1/ IFRAG(2),IRESCO,IMSHL,IRESRJ,IOULEV(6),
102      &                LEMCCK,LHADRO(0:9),LSEADI,LEVAPO,IFRAME,ITRSPT
103 * nuclear potential
104       LOGICAL LFERMI
105       COMMON /DTNPOT/ PFERMP(2),PFERMN(2),FERMOD,
106      &                EBINDP(2),EBINDN(2),EPOT(2,210),
107      &                ETACOU(2),ICOUL,LFERMI
108 * n-n cross section fluctuations
109       PARAMETER (NBINS = 1000)
110       COMMON /DTXSFL/ FLUIXX(NBINS),IFLUCT
111 * flags for particle decays
112       COMMON /DTFRPA/ MSTUX(20),PARUX(20),MSTJX(20),PARJX(20),
113      &                IMSTU(20),IPARU(20),IMSTJ(20),IPARJ(20),
114      &                NMSTU,NPARU,NMSTJ,NPARJ,PDB,PDBSEA(3),ISIG0,IPI0
115 * diquark-breaking mechanism
116       COMMON /DTDIQB/ DBRKR(3,8),DBRKA(3,8),CHAM1,CHAM3,CHAB1,CHAB3
117 * nucleon-nucleon event-generator
118       CHARACTER*8 CMODEL
119       LOGICAL LPHOIN
120       COMMON /DTMODL/ CMODEL(4),ELOJET,MCGENE,LPHOIN
121 * properties of interacting particles
122       COMMON /DTPRTA/ IT,ITZ,IP,IPZ,IJPROJ,IBPROJ,IJTARG,IBTARG
123 * properties of photon/lepton projectiles
124       COMMON /DTGPRO/ VIRT,PGAMM(4),PLEPT0(4),PLEPT1(4),PNUCL(4),IDIREC
125 * flags for diffractive interactions (DTUNUC 1.x)
126       COMMON /DTFLG3/ ISINGD,IDOUBD,IFLAGD,IDIFF
127 * parameters for hA-diffraction
128       COMMON /DTDIHA/ DIBETA,DIALPH
129 * Lorentz-parameters of the current interaction
130       COMMON /DTLTRA/ GACMS(2),BGCMS(2),GALAB,BGLAB,BLAB,
131      &                UMO,PPCM,EPROJ,PPROJ
132 * kinematical cuts for lepton-nucleus interactions
133       COMMON /DTLCUT/ ECMIN,ECMAX,XBJMIN,ELMIN,EGMIN,EGMAX,YMIN,YMAX,
134      &                Q2MIN,Q2MAX,THMIN,THMAX,Q2LI,Q2HI,ECMLI,ECMHI
135 * VDM parameter for photon-nucleus interactions
136       COMMON /DTVDMP/ RL2,EPSPOL,INTRGE(2),IDPDF,MODEGA,ISHAD(3)
137 * Glauber formalism: flags and parameters for statistics
138       LOGICAL LPROD
139       CHARACTER*8 CGLB
140       COMMON /DTGLGP/ JSTATB,JBINSB,CGLB,IOGLB,LPROD
141 * cuts for variable energy runs
142       COMMON /DTVARE/ VARELO,VAREHI,VARCLO,VARCHI
143 * flags for activated histograms
144       COMMON /DTHIS3/ IHISPP(50),IHISXS(50),IXSTBL
145
146       COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
147
148       COMMON/PYDAT3/MDCY(500,3),MDME(4000,2),BRAT(4000),KFDP(4000,5)
149
150 * LEPTO
151 **LUND single / double precision
152       REAL CUT,PARL,TMPX,TMPY,TMPW2,TMPQ2,TMPU
153       COMMON /LEPTOU/ CUT(14),LST(40),PARL(30),
154      &                TMPX,TMPY,TMPW2,TMPQ2,TMPU
155 * LEPTO
156       REAL RPPN
157       COMMON /LEPTOI/ RPPN,LEPIN,INTER
158 * steering flags for qel neutrino scattering modules
159       COMMON /QNEUTO/ DSIGSU,DSIGMC,NDSIG,NEUTYP,NEUDEC
160 * event flag
161       COMMON /DTEVNO/ NEVENT,ICASCA
162
163       INTEGER PYCOMP
164
165 C     DIMENSION XPARA(5)
166       DIMENSION XDUMB(40),IPRANG(5)
167
168       PARAMETER (MXCARD=58)
169       CHARACTER*78 CLINE,CTITLE
170       CHARACTER*60 CWHAT
171       CHARACTER*8  BLANK,SDUM
172       CHARACTER*10 CODE,CODEWD
173       CHARACTER*72 HEADER
174       LOGICAL LSTART,LEINP,LXSTAB
175       DIMENSION WHAT(6),CODE(MXCARD)
176       DATA CODE/
177      &   'TITLE     ','PROJPAR   ','TARPAR    ','ENERGY    ',
178      &   'MOMENTUM  ','CMENERGY  ','EMULSION  ','FERMI     ',
179      &   'TAUFOR    ','PAULI     ','COULOMB   ','HADRIN    ',
180      &   'EVAP      ','EMCCHECK  ','MODEL     ','PHOINPUT  ',
181      &   'GLAUBERI  ','FLUCTUAT  ','CENTRAL   ','RECOMBIN  ',
182      &   'COMBIJET  ','XCUTS     ','INTPT     ','CRONINPT  ',
183      &   'SEADISTR  ','SEASU3    ','DIQUARKS  ','RESONANC  ',
184      &   'DIFFRACT  ','SINGLECH  ','NOFRAGME  ','HADRONIZE ',
185      &   'POPCORN   ','PARDECAY  ','BEAM      ','LUND-MSTU ',
186      &   'LUND-MSTJ ','LUND-MDCY ','LUND-PARJ ','LUND-PARU ',
187      &   'OUTLEVEL  ','FRAME     ','L-TAG     ','L-ETAG    ',
188      &   'ECMS-CUT  ','VDM-PAR1  ','HISTOGRAM ','XS-TABLE  ',
189      &   'GLAUB-PAR ','GLAUB-INI ','VDM-PAR2  ','XS-QELPRO ',
190      &   'RNDMINIT  ','LEPTO-CUT ','LEPTO-LST ','LEPTO-PARL',
191      &   'START     ','STOP      '/
192       DATA BLANK /'        '/
193
194       DATA LSTART,LXSTAB,IFIRST /.TRUE.,.FALSE.,1/
195       DATA CMEOLD /0.0D0/
196
197 * --- Added by Chiara
198       
199       CHARACTER*100  ALIROOT  
200       CHARACTER*100 FILNAM
201       INTEGER*4 LNROOT
202       LOGICAL EXISTS
203       ALIROOT=' '
204
205 *---------------------------------------------------------------------
206 * at the first call of INIT: initialize event generation
207       EPNSAV = EPN
208       IF (LSTART) THEN
209          CALL DT_TITLE
210 *   initialization and test of the random number generator
211          IF (ITRSPT.NE.1) THEN
212
213             CALL FL48UT (ISRM48,ISEED1,ISEED2)
214             CALL FL48IN (54217137,ISEED1,ISEED2)
215
216          ENDIF
217 *   initialization of BAMJET, DECAY and HADRIN
218          CALL DT_DDATAR
219          CALL DT_DHADDE
220          CALL DT_DCHANT
221          CALL DT_DCHANH
222 *   set default values for input variables
223          CALL DT_DEFAUL(EPN,PPN)
224          IGLAU  = 0
225          IXSQEL = 0
226 *   flag for collision energy input
227          LEINP  = .FALSE.
228          LSTART = .FALSE.
229       ENDIF
230
231 *---------------------------------------------------------------------
232    10 CONTINUE
233
234 * bypass reading input cards (e.g. for use with Fluka)
235 *  in this case Epn is expected to carry the beam momentum
236       IF (NCASES.EQ.-1) THEN
237          IP      = NPMASS
238          IPZ     = NPCHAR
239          PPN     = EPNSAV
240          EPN     = ZERO
241          CMENER  = ZERO
242          LEINP   = .TRUE.
243          MKCRON  = 0
244          WHAT(1) = 1
245          WHAT(2) = 0
246          CODEWD  = 'START     '
247          GOTO 900
248       ENDIF
249
250 * read control card from input-unit LINP
251 C      READ(LINP,'(A78)',END=9999) CLINE
252 * ###   Read control card from specified file 
253 * ### Changed by Chiara (original version LINP=5)
254 *      OPEN(UNIT=7,
255 *     + FILE='/home/oppedisa/AliRoot/new/DPMJET/inp/PbPbLHC.inp',
256 *     + STATUS='OLD')
257
258       CALL GETENVF('ALICE_ROOT',ALIROOT)
259       LNROOT = LNBLNK(ALIROOT)
260
261       FILNAM=ALIROOT(1:LNROOT)//'/DPMJET/inp/ppLHC.inp'
262       OPEN(UNIT=7,FILE=FILNAM,STATUS='OLD')
263       OPEN(UNIT=14,FILE="nuclear.bin",STATUS='OLD')
264 *     OPEN(UNIT=6,FILE="dpm.out",STATUS='UNKNOWN')
265
266       READ(7,'(A78)',END=9999) CLINE
267
268       IF (CLINE(1:1).EQ.'*') THEN
269 * comment-line
270 C         WRITE(LOUT,'(A78)') CLINE
271          GOTO 10
272       ENDIF
273 C     READ(CLINE,1000,END=9999) CODEWD,(WHAT(I),I=1,6),SDUM
274 C1000 FORMAT(A10,6E10.0,A8)
275       DO 1008 I=1,6
276          WHAT(I) = ZERO
277  1008 CONTINUE
278       READ(CLINE,1006,END=9999) CODEWD,CWHAT,SDUM
279  1006 FORMAT(A10,A60,A8)
280       READ(CWHAT,*,END=1007) (WHAT(I),I=1,6)
281  1007 CONTINUE
282       WRITE(LOUT,1001) CODEWD,(WHAT(I),I=1,6),SDUM
283  1001 FORMAT(A10,6G10.3,A8)
284
285   900 CONTINUE
286
287 * check for valid control card and get card index
288       ICW = 0
289       DO 11 I=1,MXCARD
290          IF (CODEWD.EQ.CODE(I)) ICW = I
291    11 CONTINUE
292       IF (ICW.EQ.0) THEN
293          WRITE(LOUT,1002) CODEWD
294  1002    FORMAT(/,1X,'---> ',A10,': invalid control-card !',/)
295          GOTO 10
296       ENDIF
297       GOTO(
298 *------------------------------------------------------------
299 *       TITLE   ,  PROJPAR ,  TARPAR  ,  ENERGY  ,  MOMENTUM,
300      &  100     ,  110     ,  120     ,  130     ,  140     ,
301 *
302 *------------------------------------------------------------
303 *       CMENERGY,  EMULSION,  FERMI   ,  TAUFOR  ,  PAULI   ,
304      &  150     ,  160     ,  170     ,  180     ,  190     ,
305 *
306 *------------------------------------------------------------
307 *       COULOMB ,  HADRIN  ,  EVAP    ,  EMCCHECK,  MODEL   ,
308      &  200     ,  210     ,  220     ,  230     ,  240     ,
309 *
310 *------------------------------------------------------------
311 *       PHOINPUT,  GLAUBERI,  FLUCTUAT,  CENTRAL ,  RECOMBIN,
312      &  250     ,  260     ,  270     ,  280     ,  290     ,
313 *
314 *------------------------------------------------------------
315 *       COMBIJET,  XCUTS   ,  INTPT   ,  CRONINPT,  SEADISTR,
316      &  300     ,  310     ,  320     ,  330     ,  340     ,
317 *
318 *------------------------------------------------------------
319 *       SEASU3  ,  DIQUARKS,  RESONANC,  DIFFRACT,  SINGLECH,
320      &  350     ,  360     ,  370     ,  380     ,  390     ,
321 *
322 *------------------------------------------------------------
323 *       NOFRAGME, HADRONIZE,  POPCORN ,  PARDECAY,  BEAM    ,
324      &  400     ,  410     ,  420     ,  430     ,  440     ,
325 *
326 *------------------------------------------------------------
327 *      LUND-MSTU, LUND-MSTJ, LUND-MDCY, LUND-PARJ, LUND-PARU,
328      &  450     ,  451     ,  452     ,  460     ,  470     ,
329 *
330 *------------------------------------------------------------
331 *       OUTLEVEL,  FRAME   , L-TAG    ,  L-ETAG  ,  ECMS-CUT,
332      &  480     ,  490     ,  500     ,  510     ,  520     ,
333 *
334 *------------------------------------------------------------
335 *       VDM-PAR1, HISTOGRAM, XS-TABLE , GLAUB-PAR, GLAUB-INI,
336      &  530     ,  540     ,  550     ,  560     ,  565     ,
337 *
338 *------------------------------------------------------------
339 *               ,          ,  VDM-PAR2, XS-QELPRO, RNDMINIT ,
340      &                        570     ,  580     ,  590     ,
341 *
342 *------------------------------------------------------------
343 *      LEPTO-CUT, LEPTO-LST,LEPTO-PARL,  START   ,  STOP    )
344      &  600     ,  610     ,  620     ,  630     ,  640     ) , ICW
345 *
346 *------------------------------------------------------------
347
348       GOTO 10
349
350 *********************************************************************
351 *                                                                   *
352 *               control card:  codewd = TITLE                       *
353 *                                                                   *
354 *       what (1..6), sdum   no meaning                              *
355 *                                                                   *
356 *       Note:  The control-card following this must consist of      *
357 *              a string of characters usually giving the title of   *
358 *              the run.                                             *
359 *                                                                   *
360 *********************************************************************
361
362   100 CONTINUE
363 C      READ(LINP,'(A78)') CTITLE
364 * ###   Read control card from specified file 
365 * ### Changed by Chiara (original version LINP=5)
366       READ(7,'(A78)') CTITLE
367
368       WRITE(LOUT,'(//,5X,A78,//)') CTITLE
369       GOTO 10
370
371 *********************************************************************
372 *                                                                   *
373 *               control card:  codewd = PROJPAR                     *
374 *                                                                   *
375 *       what (1) =  mass number of projectile nucleus  default: 1   *
376 *       what (2) =  charge of projectile nucleus       default: 1   *
377 *       what (3..6)   no meaning                                    *
378 *       sdum        projectile particle code word                   *
379 *                                                                   *
380 *       Note: If sdum is defined what (1..2) have no meaning.       *
381 *                                                                   *
382 *********************************************************************
383
384   110 CONTINUE
385       IF (SDUM.EQ.BLANK) THEN
386          IP     = INT(WHAT(1))
387          IPZ    = INT(WHAT(2))
388          IJPROJ = 1
389          IBPROJ = 1
390       ELSE
391          IJPROJ = 0
392          DO 111 II=1,30
393             IF (SDUM.EQ.BTYPE(II)) THEN
394                IP     = 1
395                IPZ    = 1
396                IF (II.EQ.26) THEN
397                   IJPROJ = 135
398                ELSEIF (II.EQ.27) THEN
399                   IJPROJ = 136
400                ELSEIF (II.EQ.28) THEN
401                   IJPROJ = 133
402                ELSEIF (II.EQ.29) THEN
403                   IJPROJ = 134
404                ELSE
405                   IJPROJ = II
406                ENDIF
407                IBPROJ = IIBAR(IJPROJ)
408 * photon
409                IF ((IJPROJ.EQ.7).AND.(WHAT(1).GT.ZERO)) VIRT = WHAT(1)
410 * lepton
411                IF (((IJPROJ.EQ. 3).OR.(IJPROJ.EQ. 4).OR.
412      &              (IJPROJ.EQ.10).OR.(IJPROJ.EQ.11)).AND.
413      &                              (WHAT(1).GT.ZERO)) Q2HI = WHAT(1)
414             ENDIF
415   111    CONTINUE
416          IF (IJPROJ.EQ.0) THEN
417             WRITE(LOUT,1110)
418  1110       FORMAT(/,1X,'invalid PROJPAR card !',/)
419             GOTO 9999
420          ENDIF
421       ENDIF
422       GOTO 10
423
424 *********************************************************************
425 *                                                                   *
426 *               control card:  codewd = TARPAR                      *
427 *                                                                   *
428 *       what (1) =  mass number of target nucleus      default: 1   *
429 *       what (2) =  charge of target nucleus           default: 1   *
430 *       what (3..6)   no meaning                                    *
431 *       sdum        target particle code word                       *
432 *                                                                   *
433 *       Note: If sdum is defined what (1..2) have no meaning.       *
434 *                                                                   *
435 *********************************************************************
436
437   120 CONTINUE
438       IF (SDUM.EQ.BLANK) THEN
439          IT     = INT(WHAT(1))
440          ITZ    = INT(WHAT(2))
441          IJTARG = 1
442          IBTARG = 1
443       ELSE
444          IJTARG = 0
445          DO 121 II=1,30
446             IF (SDUM.EQ.BTYPE(II)) THEN
447                IT     = 1
448                ITZ    = 1
449                IJTARG = II
450                IBTARG = IIBAR(IJTARG)
451             ENDIF
452   121    CONTINUE
453          IF (IJTARG.EQ.0) THEN
454             WRITE(LOUT,1120)
455  1120       FORMAT(/,1X,'invalid TARPAR card !',/)
456             GOTO 9999
457          ENDIF
458       ENDIF
459       GOTO 10
460
461 *********************************************************************
462 *                                                                   *
463 *               control card:  codewd = ENERGY                      *
464 *                                                                   *
465 *       what (1) =  energy (GeV) of projectile in Lab.              *
466 *                   if what(1) < 0:  |what(1)| = kinetic energy     *
467 *                                                default: 200 GeV   *
468 *                   if |what(2)| > 0: min. energy for variable      *
469 *                                     energy runs                   *
470 *       what (2) =  max. energy for variable energy runs            *
471 *                   if what(2) < 0:  |what(2)| = kinetic energy     *
472 *                                                                   *
473 *********************************************************************
474
475   130 CONTINUE
476       EPN    = WHAT(1)
477       PPN    = ZERO
478       CMENER = ZERO
479       IF ((ABS(WHAT(2)).GT.ZERO).AND.
480      &    (ABS(WHAT(2)).GT.ABS(WHAT(1)))) THEN
481          VARELO = WHAT(1)
482          VAREHI = WHAT(2)
483          EPN    = VAREHI
484       ENDIF
485       LEINP  = .TRUE.
486       GOTO 10
487
488 *********************************************************************
489 *                                                                   *
490 *               control card:  codewd = MOMENTUM                    *
491 *                                                                   *
492 *       what (1) =  momentum (GeV/c) of projectile in Lab.          *
493 *                                                default: 200 GeV/c *
494 *       what (2..6), sdum   no meaning                              *
495 *                                                                   *
496 *********************************************************************
497
498   140 CONTINUE
499       EPN    = ZERO
500       PPN    = WHAT(1)
501       CMENER = ZERO
502       LEINP  = .TRUE.
503       GOTO 10
504
505 *********************************************************************
506 *                                                                   *
507 *               control card:  codewd = CMENERGY                    *
508 *                                                                   *
509 *       what (1) =  energy in nucleon-nucleon cms.                  *
510 *                                                default: none      *
511 *       what (2..6), sdum   no meaning                              *
512 *                                                                   *
513 *********************************************************************
514
515   150 CONTINUE
516       EPN    = ZERO
517       PPN    = ZERO
518       CMENER = WHAT(1)
519       LEINP  = .TRUE.
520       GOTO 10
521
522 *********************************************************************
523 *                                                                   *
524 *               control card:  codewd = EMULSION                    *
525 *                                                                   *
526 *               definition of nuclear emulsions                     *
527 *                                                                   *
528 *     what(1)      mass number of emulsion component                *
529 *     what(2)      charge of emulsion component                     *
530 *     what(3)      fraction of events in which a scattering on a    *
531 *                  nucleus of this properties is performed          *
532 *     what(4,5,6)  as what(1,2,3) but for another component         *
533 *                                             default: no emulsion  *
534 *     sdum         no meaning                                       *
535 *                                                                   *
536 *     Note: If this input-card is once used with valid parameters   *
537 *           TARPAR is obsolete.                                     *
538 *           Not the absolute values of the fractions are important  *
539 *           but only the ratios of fractions of different comp.     *
540 *           This control card can be repeatedly used to define      *
541 *           emulsions consisting of up to 10 elements.              *
542 *                                                                   *
543 *********************************************************************
544
545   160 CONTINUE
546       IF ((WHAT(1).GT.ZERO).AND.(WHAT(2).GT.ZERO)
547      &                     .AND.(ABS(WHAT(3)).GT.ZERO)) THEN
548          NCOMPO = NCOMPO+1
549          IF (NCOMPO.GT.NCOMPX) THEN
550             WRITE(LOUT,1600)
551             STOP
552          ENDIF
553          IEMUMA(NCOMPO) = INT(WHAT(1))
554          IEMUCH(NCOMPO) = INT(WHAT(2))
555          EMUFRA(NCOMPO) = WHAT(3)
556          IEMUL = 1
557 C        CALL SHMAKF(IDUM,IDUM,IEMUMA(NCOMPO),IEMUCH(NCOMPO))
558       ENDIF
559       IF ((WHAT(4).GT.ZERO).AND.(WHAT(5).GT.ZERO)
560      &                     .AND.(ABS(WHAT(6)).GT.ZERO)) THEN
561          NCOMPO = NCOMPO+1
562          IF (NCOMPO.GT.NCOMPX) THEN
563             WRITE(LOUT,1001)
564             STOP
565          ENDIF
566          IEMUMA(NCOMPO) = INT(WHAT(4))
567          IEMUCH(NCOMPO) = INT(WHAT(5))
568          EMUFRA(NCOMPO) = WHAT(6)
569 C        CALL SHMAKF(IDUM,IDUM,IEMUMA(NCOMPO),IEMUCH(NCOMPO))
570       ENDIF
571  1600 FORMAT(1X,'too many emulsion components - program stopped')
572       GOTO 10
573
574 *********************************************************************
575 *                                                                   *
576 *               control card:  codewd = FERMI                       *
577 *                                                                   *
578 *       what (1) = -1 Fermi-motion of nucleons not treated          *
579 *                                                 default: 1        *
580 *       what (2) =    scale factor for Fermi-momentum               *
581 *                                                 default: 0.75     *
582 *       what (3..6), sdum   no meaning                              *
583 *                                                                   *
584 *********************************************************************
585
586   170 CONTINUE
587       IF (WHAT(1).EQ.-1.0D0) THEN
588          LFERMI = .FALSE.
589       ELSE
590          LFERMI = .TRUE.
591       ENDIF
592       XMOD = WHAT(2)
593       IF (XMOD.GE.ZERO) FERMOD = XMOD
594       GOTO 10
595
596 *********************************************************************
597 *                                                                   *
598 *               control card:  codewd = TAUFOR                      *
599 *                                                                   *
600 *          formation time supressed intranuclear cascade            *
601 *                                                                   *
602 *    what (1)      formation time (in fm/c)                         *
603 *                  note: what(1)=10. corresponds roughly to an      *
604 *                        average formation time of 1 fm/c           *
605 *                                                 default: 5. fm/c  *
606 *    what (2)      number of generations followed                   *
607 *                                                 default: 25       *
608 *    what (3) = 1. p_t-dependent formation zone                     *
609 *             = 2. constant formation zone                          *
610 *                                                 default: 1        *
611 *    what (4)      modus of selection of nucleus where the          *
612 *                  cascade if followed first                        *
613 *             = 1.  proj./target-nucleus with probab. 1/2           *
614 *             = 2.  nucleus with highest mass                       *
615 *             = 3.  proj. nucleus if particle is moving in pos. z   *
616 *                   targ. nucleus if particle is moving in neg. z   *
617 *                                                 default: 1        *
618 *    what (5..6), sdum   no meaning                                 *
619 *                                                                   *
620 *********************************************************************
621
622   180 CONTINUE
623       TAUFOR = WHAT(1)
624       KTAUGE = INT(WHAT(2))
625       INCMOD = 1
626       IF ((WHAT(3).GE.1.0D0).AND.(WHAT(3).LE.2.0D0))
627      &                                    ITAUVE = INT(WHAT(3))
628       IF ((WHAT(4).GE.1.0D0).AND.(WHAT(4).LE.3.0D0))
629      &                                    INCMOD = INT(WHAT(4))
630       GOTO 10
631
632 *********************************************************************
633 *                                                                   *
634 *               control card:  codewd = PAULI                       *
635 *                                                                   *
636 *       what (1) =  -1  Pauli's principle for secondary             *
637 *                       interactions not treated                    *
638 *                                                    default: 1     *
639 *       what (2..6), sdum   no meaning                              *
640 *                                                                   *
641 *********************************************************************
642
643   190 CONTINUE
644       IF (WHAT(1).EQ.-1.0D0) THEN
645          LPAULI = .FALSE.
646       ELSE
647          LPAULI = .TRUE.
648       ENDIF
649       GOTO 10
650
651 *********************************************************************
652 *                                                                   *
653 *               control card:  codewd = COULOMB                     *
654 *                                                                   *
655 *       what (1) = -1. Coulomb-energy treatment switched off        *
656 *                                                    default: 1     *
657 *       what (2..6), sdum   no meaning                              *
658 *                                                                   *
659 *********************************************************************
660
661   200 CONTINUE
662       ICOUL = 1
663       IF (WHAT(1).EQ.-1.0D0) THEN
664          ICOUL = 0
665       ELSE
666          ICOUL = 1
667       ENDIF
668       GOTO 10
669
670 *********************************************************************
671 *                                                                   *
672 *               control card:  codewd = HADRIN                      *
673 *                                                                   *
674 *                       HADRIN module                               *
675 *                                                                   *
676 *    what (1) = 0. elastic/inelastic interactions with probab.      *
677 *                  as defined by cross-sections                     *
678 *             = 1. inelastic interactions forced                    *
679 *             = 2. elastic interactions forced                      *
680 *                                                 default: 1        *
681 *    what (2)      upper threshold in total energy (GeV) below      *
682 *                  which interactions are sampled by HADRIN         *
683 *                                                 default: 5. GeV   *
684 *    what (3..6), sdum   no meaning                                 *
685 *                                                                   *
686 *********************************************************************
687
688   210 CONTINUE
689       IWHAT = INT(WHAT(1))
690       IF ((IWHAT.GE.0).AND.(IWHAT.LE.2)) INTHAD = IWHAT
691       IF ((WHAT(2).GT.ZERO).AND.(WHAT(2).LT.15.0D0)) EHADTH = WHAT(2)
692       GOTO 10
693
694 *********************************************************************
695 *                                                                   *
696 *               control card:  codewd = EVAP                        *
697 *                                                                   *
698 *                    evaporation module                             *
699 *                                                                   *
700 *  what (1) =< -1 ==> evaporation is switched off                   *
701 *           >=  1 ==> evaporation is performed                      *
702 *                                                                   *
703 *         what (1) = i1 + i2*10 + i3*100 + i4*10000                 *
704 *                    (i1, i2, i3, i4 >= 0 )                         *
705 *                                                                   *
706 *   i1 is the flag for selecting the T=0 level density option used  *
707 *      =  1: standard EVAP level densities with Cook pairing        *
708 *            energies                                               *
709 *      =  2: Z,N-dependent Gilbert & Cameron level densities        *
710 *                                                        (default)  *
711 *      =  3: Julich A-dependent level densities                     *
712 *      =  4: Z,N-dependent Brancazio & Cameron level densities      *
713 *                                                                   *
714 *   i2 >= 1: high energy fission activated                          *
715 *            (default high energy fission activated)                *
716 *                                                                   *
717 *   i3 =  0: No energy dependence for level densities               *
718 *      =  1: Standard Ignyatuk (1975, 1st) energy dependence        *
719 *            for level densities (default)                          *
720 *      =  2: Standard Ignyatuk (1975, 1st) energy dependence        *
721 *            for level densities with NOT used set of parameters    *
722 *      =  3: Standard Ignyatuk (1975, 1st) energy dependence        *
723 *            for level densities with NOT used set of parameters    *
724 *      =  4: Second   Ignyatuk (1975, 2nd) energy dependence        *
725 *            for level densities                                    *
726 *      =  5: Second   Ignyatuk (1975, 2nd) energy dependence        *
727 *            for level densities with fit 1 Iljinov & Mebel set of  *
728 *            parameters                                             *
729 *      =  6: Second   Ignyatuk (1975, 2nd) energy dependence        *
730 *            for level densities with fit 2 Iljinov & Mebel set of  *
731 *            parameters                                             *
732 *      =  7: Second   Ignyatuk (1975, 2nd) energy dependence        *
733 *            for level densities with fit 3 Iljinov & Mebel set of  *
734 *            parameters                                             *
735 *      =  8: Second   Ignyatuk (1975, 2nd) energy dependence        *
736 *            for level densities with fit 4 Iljinov & Mebel set of  *
737 *            parameters                                             *
738 *                                                                   *
739 *   i4 >= 1: Original Gilbert and Cameron pairing energies used     *
740 *            (default Cook's modified pairing energies)             *
741 *                                                                   *
742 *  what (2) = ig + 10 * if   (ig and if must have the same sign)    *
743 *                                                                   *
744 *   ig =< -1 ==> deexcitation gammas are not produced               *
745 *                (if the evaporation step is not performed          *
746 *                 they are never produced)                          *
747 *   if =< -1 ==> Fermi Break Up is not invoked                      *
748 *                (if the evaporation step is not performed          *
749 *                 it is never invoked)                              *
750 *   The default is: deexcitation gamma produced and Fermi break up  *
751 *                   activated for the new  preequilibrium, not      *
752 *                   activated otherwise.                            *
753 *  what (3..6), sdum   no meaning                                   *
754 *                                                                   *
755 *********************************************************************
756
757  220  CONTINUE
758
759       IF (WHAT(1).LE.-1.0D0) THEN
760          LEVPRT = .FALSE.
761          LDEEXG = .FALSE.
762          LHEAVY = .FALSE.
763          GOTO 10
764       ENDIF
765       WHTSAV = WHAT (1)
766       IF ( NINT (WHAT (1)) .GE. 10000 ) THEN
767          LLVMOD   = .FALSE.
768          JLVHLP   = NINT (WHAT (1)) / 10000
769          WHAT (1) = WHAT (1) - 10000.D+00 * JLVHLP
770       END IF
771       IF ( NINT (WHAT (1)) .GE. 100 ) THEN
772          JLVMOD   = NINT (WHAT (1)) / 100
773          WHAT (1) = WHAT (1) - 100.D+00 * JLVMOD
774       END IF
775       IF ( NINT (WHAT (1)) .GE. 10  ) THEN
776          IFISS    = 1
777          JLVHLP   = NINT (WHAT (1)) / 10
778          WHAT (1) = WHAT (1) - 10.D+00 * JLVHLP
779       ELSE IF ( NINT (WHTSAV) .NE. 0 ) THEN
780          IFISS    = 0
781       END IF
782       IF ( NINT (WHAT (1)) .GE. 0 ) THEN
783          LEVPRT = .TRUE.
784          ILVMOD = NINT (WHAT(1))
785          IF ( ABS (NINT (WHAT (2))) .GE. 10  ) THEN
786             LFRMBK   = .TRUE.
787             JLVHLP   = NINT (WHAT (2)) / 10
788             WHAT (2) = WHAT (2) - 10.D+00 * JLVHLP
789          ELSE IF ( NINT (WHAT (2)) .NE. 0 ) THEN
790             LFRMBK   = .FALSE.
791          END IF
792          IF ( NINT (WHAT (2)) .GE. 0 ) THEN
793             LDEEXG = .TRUE.
794          ELSE
795             LDEEXG = .FALSE.
796          END IF
797 **sr heavies are always put to /FKFHVY/
798 C        IF ( NINT (WHAT(3)) .GE. 1 ) THEN
799 C           LHEAVY = .TRUE.
800 C        ELSE
801 C           LHEAVY = .FALSE.
802 C        END IF
803          LHEAVY = .TRUE.
804       ELSE
805          LEVPRT = .FALSE.
806          LDEEXG = .FALSE.
807          LHEAVY = .FALSE.
808       END IF
809
810       LOLDEV = .FALSE.
811
812       GOTO 10
813
814 *********************************************************************
815 *                                                                   *
816 *               control card:  codewd = EMCCHECK                    *
817 *                                                                   *
818 *    extended energy-momentum / quantum-number conservation check   *
819 *                                                                   *
820 *       what (1) = -1   extended check not performed                *
821 *                                                    default: 1.    *
822 *       what (2..6), sdum   no meaning                              *
823 *                                                                   *
824 *********************************************************************
825
826   230 CONTINUE
827       IF (WHAT(1).EQ.-1) THEN
828          LEMCCK = .FALSE.
829       ELSE
830          LEMCCK = .TRUE.
831       ENDIF
832       GOTO 10
833
834 *********************************************************************
835 *                                                                   *
836 *               control card:  codewd = MODEL                       *
837 *                                                                   *
838 *     Model to be used to treat nucleon-nucleon interactions        *
839 *                                                                   *
840 *       sdum = DTUNUC    two-chain model                            *
841 *            = PHOJET    multiple chains including minijets         *
842 *            = LEPTO     DIS                                        *
843 *            = QNEUTRIN  quasi-elastic neutrino scattering          *
844 *                                                  default: PHOJET  *
845 *                                                                   *
846 *       if sdum = LEPTO:                                            *
847 *       what (1)         (variable INTER)                           *
848 *                        = 1  gamma exchange                        *
849 *                        = 2  W+-   exchange                        *
850 *                        = 3  Z0    exchange                        *
851 *                        = 4  gamma/Z0 exchange                     *
852 *                                                                   *
853 *       if sdum = QNEUTRIN:                                         *
854 *       what (1)         = 0  elastic scattering on nucleon and     *
855 *                             tau does not decay (default)          *
856 *                        = 1  decay of tau into mu..                *
857 *                        = 2  decay of tau into e..                 *
858 *                        = 10 CC events on p and n                  *
859 *                        = 11 NC events on p and n                  *
860 *                                                                   *
861 *       what (2..6)      no meaning                                 *
862 *                                                                   *
863 *********************************************************************
864
865   240 CONTINUE
866       IF (SDUM.EQ.CMODEL(1)) THEN
867          MCGENE = 1
868       ELSEIF (SDUM.EQ.CMODEL(2)) THEN
869          MCGENE = 2
870       ELSEIF (SDUM.EQ.CMODEL(3)) THEN
871          MCGENE = 3
872          IF ((WHAT(1).GE.1.0D0).AND.(WHAT(1).LE.4.0D0))
873      &      INTER = INT(WHAT(1))
874       ELSEIF (SDUM.EQ.CMODEL(4)) THEN
875          MCGENE = 4
876          IWHAT  = INT(WHAT(1))
877          IF ((IWHAT.EQ.1 ).OR.(IWHAT.EQ.2 ).OR.
878      &       (IWHAT.EQ.10).OR.(IWHAT.EQ.11))
879      &      NEUDEC = IWHAT
880       ELSE
881          STOP ' Unknown model !'
882       ENDIF
883       GOTO 10
884
885 *********************************************************************
886 *                                                                   *
887 *               control card:  codewd = PHOINPUT                    *
888 *                                                                   *
889 *       Start of input-section for PHOJET-specific input-cards      *
890 *       Note:  This section will not be finished before giving      *
891 *              ENDINPUT-card                                        *
892 *       what (1..6), sdum   no meaning                              *
893 *                                                                   *
894 *********************************************************************
895
896   250 CONTINUE
897       IF (LPHOIN) THEN
898
899 C         CALL PHO_INIT(LINP,IREJ1)
900 * ###   Read control card from specified file 
901 * ### Changed by Chiara (original version LINP=5)
902          CALL PHO_INIT(7,IREJ1)
903
904          IF (IREJ1.NE.0) THEN
905             WRITE(LOUT,'(1X,A)')'INIT:   reading PHOJET-input failed'
906             STOP
907          ENDIF
908          LPHOIN = .FALSE.
909       ENDIF
910       GOTO 10
911
912 *********************************************************************
913 *                                                                   *
914 *               control card:  codewd = GLAUBERI                    *
915 *                                                                   *
916 *        Pre-initialization of impact parameter selection           *
917 *                                                                   *
918 *        what (1..6), sdum   no meaning                             *
919 *                                                                   *
920 *********************************************************************
921
922   260 CONTINUE
923       IF (IFIRST.NE.99) THEN
924          CALL DT_RNDMST(12,34,56,78)
925          CALL DT_RNDMTE(1)
926          OPEN(40,FILE='shm.out',STATUS='UNKNOWN')
927 C        OPEN(11,FILE='shm.dbg',STATUS='UNKNOWN')
928          IFIRST = 99
929       ENDIF
930
931       IPPN = 8
932       PLOW = 10.0D0
933 C     IPPN = 1
934 C     PLOW = 100.0D0
935       PHI  = 1.0D5
936       APLOW = LOG10(PLOW)
937       APHI  = LOG10(PHI)
938       ADP   = (APHI-APLOW)/DBLE(IPPN)
939
940       IPLOW = 1
941       IDIP  = 1
942       IIP   = 5
943 C     IPLOW = 1
944 C     IDIP  = 1
945 C     IIP   = 1
946       IPRANG(1) = 1
947       IPRANG(2) = 2
948       IPRANG(3) = 5
949       IPRANG(4) = 10
950       IPRANG(5) = 20
951
952       ITLOW = 30
953       IDIT  = 3
954       IIT   = 60
955 C     IDIT  = 10
956 C     IIT   = 21
957
958       DO 473 NCIT=1,IIT
959          IT   = ITLOW+(NCIT-1)*IDIT
960 C        IPHI = IT
961 C        IDIP = 10
962 C        IIP  = (IPHI-IPLOW)/IDIP
963 C        IF (IIP.EQ.0) IIP = 1
964 C        IF (IT.EQ.IPLOW) IIP = 0
965
966          DO 472 NCIP=1,IIP
967             IP = IPRANG(NCIP)
968 CC           IF (NCIP.LE.IIP) THEN
969 C               IP = IPLOW+(NCIP-1)*IDIP
970 CC           ELSE
971 CC              IP = IT
972 CC           ENDIF
973             IF (IP.GT.IT) GOTO 472
974
975             DO 471 NCP=1,IPPN+1
976                APPN = APLOW+DBLE(NCP-1)*ADP
977                PPN  = 10**APPN
978
979                OPEN(12,FILE='shm.sta',STATUS='UNKNOWN')
980                WRITE(12,'(1X,2I5,E15.3)') IP,IT,PPN
981                CLOSE(12)
982
983                XLIM1 = 0.0D0
984                XLIM2 = 50.0D0
985                XLIM3 = ZERO
986                IBIN  = 50
987                CALL DT_NEWHGR(XDUM,XDUM,XDUM,XDUMB,-1,IHDUM)
988                CALL DT_NEWHGR(XLIM1,XLIM2,XLIM3,XDUMB,IBIN,IHSHMA)
989
990                NEVFIT = 5
991 C              IF ((IP.GT.10).OR.(IT.GT.10)) THEN
992 C                 NEVFIT = 5
993 C              ELSE
994 C                 NEVFIT = 10
995 C              ENDIF
996                SIGAV  = 0.0D0
997
998                DO 478 I=1,NEVFIT
999                   CALL DT_SHMAKI(IP,IDUM1,IT,IDUM1,IJPROJ,PPN,99)
1000                   SIGAV = SIGAV+XSPRO(1,1,1)
1001                   DO 479 J=1,50
1002                      XC = DBLE(J)
1003                      CALL DT_FILHGR(XC,BSITE(1,1,1,J),IHSHMA,I)
1004   479             CONTINUE
1005   478          CONTINUE
1006
1007                CALL DT_EVTHIS(IDUM)
1008                HEADER = ' BSITE'
1009 C              CALL OUTGEN(IHSHMA,0,0,0,0,0,HEADER,0,NEVFIT,ONE,0,1,-1)
1010
1011 C              CALL GENFIT(XPARA)
1012 C              WRITE(40,'(2I4,E11.3,F6.0,5E11.3)')
1013 C    &              IP,IT,PPN,SIGAV/DBLE(NEVFIT),XPARA
1014
1015   471       CONTINUE
1016
1017   472    CONTINUE
1018
1019   473 CONTINUE
1020
1021       STOP
1022
1023 *********************************************************************
1024 *                                                                   *
1025 *               control card:  codewd = FLUCTUAT                    *
1026 *                                                                   *
1027 *           Treatment of cross section fluctuations                 *
1028 *                                                                   *
1029 *       what (1) = 1  treat cross section fluctuations              *
1030 *                                                    default: 0.    *
1031 *       what (1..6), sdum   no meaning                              *
1032 *                                                                   *
1033 *********************************************************************
1034
1035  270  CONTINUE
1036       IFLUCT = 0
1037       IF (WHAT(1).EQ.ONE) THEN
1038          IFLUCT = 1
1039          CALL DT_FLUINI
1040       ENDIF
1041       GOTO 10
1042
1043 *********************************************************************
1044 *                                                                   *
1045 *               control card:  codewd = CENTRAL                     *
1046 *                                                                   *
1047 *       what (1) = 1.  central production forced     default: 0     *
1048 *  if what (1) < 0 and > -100                                       *
1049 *       what (2) = min. impact parameter             default: 0     *
1050 *       what (3) = max. impact parameter             default: b_max *
1051 *  if what (1) < -99                                                *
1052 *       what (2) = fraction of cross section         default: 1     *
1053 *  if what (1) = -1 : evaporation/fzc suppressed                    *
1054 *  if what (1) < -1 : evaporation/fzc allowed                       *
1055 *                                                                   *
1056 *       what (4..6), sdum   no meaning                              *
1057 *                                                                   *
1058 *********************************************************************
1059
1060   280 CONTINUE
1061       ICENTR = INT(WHAT(1))
1062       IF (ICENTR.LT.0) THEN
1063          IF (ICENTR.GT.-100) THEN
1064             BIMIN = WHAT(2)
1065             BIMAX = WHAT(3)
1066          ELSE
1067             XSFRAC = WHAT(2)
1068          ENDIF
1069       ENDIF
1070       GOTO 10
1071
1072 *********************************************************************
1073 *                                                                   *
1074 *               control card:  codewd = RECOMBIN                    *
1075 *                                                                   *
1076 *                     Chain recombination                           *
1077 *        (recombine S-S and V-V chains to V-S chains)               *
1078 *                                                                   *
1079 *       what (1) = -1. recombination switched off    default: 1     *
1080 *       what (2..6), sdum   no meaning                              *
1081 *                                                                   *
1082 *********************************************************************
1083
1084   290 CONTINUE
1085       IRECOM = 1
1086       IF (WHAT(1).EQ.-1.0D0) IRECOM = 0
1087       GOTO 10
1088
1089 *********************************************************************
1090 *                                                                   *
1091 *               control card:  codewd = COMBIJET                    *
1092 *                                                                   *
1093 *               chain fusion (2 q-aq --> qq-aqaq)                   *
1094 *                                                                   *
1095 *       what (1) = 1   fusion treated                               *
1096 *                                                    default: 0.    *
1097 *       what (2)       minimum number of uncombined chains from     *
1098 *                      single projectile or target nucleons         *
1099 *                                                    default: 0.    *
1100 *       what (3..6), sdum   no meaning                              *
1101 *                                                                   *
1102 *********************************************************************
1103
1104   300 CONTINUE
1105       LCO2CR = .FALSE.
1106       IF (INT(WHAT(1)).EQ.1) LCO2CR = .TRUE.
1107       IF (WHAT(2).GE.ZERO) CUTOF = WHAT(2)
1108       GOTO 10
1109
1110 *********************************************************************
1111 *                                                                   *
1112 *               control card:  codewd = XCUTS                       *
1113 *                                                                   *
1114 *                 thresholds for x-sampling                         *
1115 *                                                                   *
1116 *    what (1)    defines lower threshold for val.-q x-value (CVQ)   *
1117 *                                                 default: 1.       *
1118 *    what (2)    defines lower threshold for val.-qq x-value (CDQ)  *
1119 *                                                 default: 2.       *
1120 *    what (3)    defines lower threshold for sea-q x-value (CSEA)   *
1121 *                                                 default: 0.2      *
1122 *    what (4)    sea-q x-values in S-S chains (SSMIMA)              *
1123 *                                                 default: 0.14     *
1124 *    what (5)    not used                                           *
1125 *                                                 default: 2.       *
1126 *    what (6), sdum   no meaning                                    *
1127 *                                                                   *
1128 *    Note: Lower thresholds (what(1..3)) are def. as x_thr=CXXX/ECM *
1129 *                                                                   *
1130 *********************************************************************
1131
1132   310 CONTINUE
1133       IF (WHAT(1).GE.0.5D0) CVQ    = WHAT(1)
1134       IF (WHAT(2).GE.ONE)   CDQ    = WHAT(2)
1135       IF (WHAT(3).GE.0.1D0) CSEA   = WHAT(3)
1136       IF (WHAT(4).GE.ZERO) THEN
1137          SSMIMA = WHAT(4)
1138          SSMIMQ = SSMIMA**2
1139       ENDIF
1140       IF (WHAT(5).GT.2.0D0) VVMTHR = WHAT(5)
1141       GOTO 10
1142
1143 *********************************************************************
1144 *                                                                   *
1145 *               control card:  codewd = INTPT                       *
1146 *                                                                   *
1147 *     what (1) = -1   intrinsic transverse momenta of partons       *
1148 *                     not treated                default: 1         *
1149 *     what (2..6), sdum   no meaning                                *
1150 *                                                                   *
1151 *********************************************************************
1152
1153   320 CONTINUE
1154       IF (WHAT(1).EQ.-1.0D0) THEN
1155          LINTPT = .FALSE.
1156       ELSE
1157          LINTPT = .TRUE.
1158       ENDIF
1159       GOTO 10
1160
1161 *********************************************************************
1162 *                                                                   *
1163 *               control card:  codewd = CRONINPT                    *
1164 *                                                                   *
1165 *    Cronin effect (multiple scattering of partons at chain ends)   *
1166 *                                                                   *
1167 *       what (1) = -1  Cronin effect not treated     default: 1     *
1168 *       what (2) = 0   scattering parameter          default: 0.64  *
1169 *       what (3..6), sdum   no meaning                              *
1170 *                                                                   *
1171 *********************************************************************
1172
1173   330 CONTINUE
1174       IF (WHAT(1).EQ.-1.0D0) THEN
1175          MKCRON = 0
1176       ELSE
1177          MKCRON = 1
1178       ENDIF
1179       CRONCO = WHAT(2)
1180       GOTO 10
1181
1182 *********************************************************************
1183 *                                                                   *
1184 *               control card:  codewd = SEADISTR                    *
1185 *                                                                   *
1186 *     what (1)  (XSEACO)  sea(x) prop. 1/x**what (1)   default: 1.  *
1187 *     what (2)  (UNON)                                 default: 2.  *
1188 *     what (3)  (UNOM)                                 default: 1.5 *
1189 *     what (4)  (UNOSEA)                               default: 5.  *
1190 *                        qdis(x) prop. (1-x)**what (1)  etc.        *
1191 *     what (5..6), sdum   no meaning                                *
1192 *                                                                   *
1193 *********************************************************************
1194
1195   340 CONTINUE
1196       XSEACO = WHAT(1)
1197       XSEACU = 1.05D0-XSEACO
1198       UNON   = WHAT(2)
1199       IF (UNON.LT.0.1D0) UNON = 2.0D0
1200       UNOM   = WHAT(3)
1201       IF (UNOM.LT.0.1D0) UNOM = 1.5D0
1202       UNOSEA = WHAT(4)
1203       IF (UNOSEA.LT.0.1D0) UNOSEA = 5.0D0
1204       GOTO 10
1205
1206 *********************************************************************
1207 *                                                                   *
1208 *               control card:  codewd = SEASU3                      *
1209 *                                                                   *
1210 *          Treatment of strange-quarks at chain ends                *
1211 *                                                                   *
1212 *       what (1)   (SEASQ)  strange-quark supression factor         *
1213 *                  iflav = 1.+rndm*(2.+SEASQ)                       *
1214 *                                                    default: 1.    *
1215 *       what (2..6), sdum   no meaning                              *
1216 *                                                                   *
1217 *********************************************************************
1218
1219   350 CONTINUE
1220       SEASQ = WHAT(1)
1221       GOTO 10
1222
1223 *********************************************************************
1224 *                                                                   *
1225 *               control card:  codewd = DIQUARKS                    *
1226 *                                                                   *
1227 *     what (1) = -1.  sea-diquark/antidiquark-pairs not treated     *
1228 *                                                    default: 1.    *
1229 *     what (2..6), sdum   no meaning                                *
1230 *                                                                   *
1231 *********************************************************************
1232
1233  360  CONTINUE
1234       IF (WHAT(1).EQ.-1.0D0) THEN
1235          LSEADI = .FALSE.
1236       ELSE
1237          LSEADI = .TRUE.
1238       ENDIF
1239       GOTO 10
1240
1241 *********************************************************************
1242 *                                                                   *
1243 *               control card:  codewd = RESONANC                    *
1244 *                                                                   *
1245 *                 treatment of low mass chains                      *
1246 *                                                                   *
1247 *    what (1) = -1 low chain masses are not corrected for resonance *
1248 *                  masses (obsolete for BAMJET-fragmentation)       *
1249 *                                       default: 1.                 *
1250 *    what (2) = -1 massless partons     default: 1. (massive)       *
1251 *                                       default: 1. (massive)       *
1252 *    what (3) = -1 chain-system containing chain of too small       *
1253 *                  mass is rejected (note: this does not fully      *
1254 *                  apply to S-S chains) default: 0.                 *
1255 *    what (4..6), sdum   no meaning                                 *
1256 *                                                                   *
1257 *********************************************************************
1258
1259   370 CONTINUE
1260       IRESCO = 1
1261       IMSHL  = 1
1262       IRESRJ = 0
1263       IF (WHAT(1).EQ.-ONE) IRESCO = 0
1264       IF (WHAT(2).EQ.-ONE) IMSHL  = 0
1265       IF (WHAT(3).EQ.-ONE) IRESRJ = 1
1266       GOTO 10
1267
1268 *********************************************************************
1269 *                                                                   *
1270 *               control card:  codewd = DIFFRACT                    *
1271 *                                                                   *
1272 *                Treatment of diffractive events                    *
1273 *                                                                   *
1274 *     what (1) = (ISINGD) 0  no single diffraction                  *
1275 *                         1  single diffraction included            *
1276 *                       +-2  single diffractive events only         *
1277 *                       +-3  projectile single diffraction only     *
1278 *                       +-4  target single diffraction only         *
1279 *                        -5  double pomeron exchange only           *
1280 *                      (neg. sign applies to PHOJET events)         *
1281 *                                                     default: 0.   *
1282 *                                                                   *
1283 *     what (2) = (IDOUBD) 0  no double diffraction                  *
1284 *                         1  double diffraction included            *
1285 *                         2  double diffractive events only         *
1286 *                                                     default: 0.   *
1287 *     what (3) = 1 projectile diffraction treated (2-channel form.) *
1288 *                                                     default: 0.   *
1289 *     what (4) = alpha-parameter in projectile diffraction          *
1290 *                                                     default: 0.   *
1291 *     what (5..6), sdum   no meaning                                *
1292 *                                                                   *
1293 *********************************************************************
1294
1295   380 CONTINUE
1296       IF (ABS(WHAT(1)).GT.ZERO) ISINGD = INT(WHAT(1))
1297       IF (ABS(WHAT(2)).GT.ZERO) IDOUBD = INT(WHAT(2))
1298       IF ((ISINGD.GT.1).AND.(IDOUBD.GT.1)) THEN
1299          WRITE(LOUT,1380)
1300  1380    FORMAT(1X,'INIT:   inconsistent DIFFRACT - input !',/,
1301      &          11X,'IDOUBD is reset to zero')
1302          IDOUBD = 0
1303       ENDIF
1304       IF (WHAT(3).GT.ZERO) DIBETA = WHAT(3)
1305       IF (WHAT(4).GT.ZERO) DIALPH = WHAT(4)
1306       GOTO 10
1307
1308 *********************************************************************
1309 *                                                                   *
1310 *               control card:  codewd = SINGLECH                    *
1311 *                                                                   *
1312 *       what (1) = 1.  Regge contribution (one chain) included      *
1313 *                                                   default: 0.     *
1314 *       what (2..6), sdum   no meaning                              *
1315 *                                                                   *
1316 *********************************************************************
1317
1318  390  CONTINUE
1319       ISICHA = 0
1320       IF (WHAT(1).EQ.ONE) ISICHA = 1
1321       GOTO 10
1322
1323 *********************************************************************
1324 *                                                                   *
1325 *               control card:  codewd = NOFRAGME                    *
1326 *                                                                   *
1327 *                 biased chain hadronization                        *
1328 *                                                                   *
1329 *       what (1..6) = -1  no of hadronizsation of S-S chains        *
1330 *                   = -2  no of hadronizsation of D-S chains        *
1331 *                   = -3  no of hadronizsation of S-D chains        *
1332 *                   = -4  no of hadronizsation of S-V chains        *
1333 *                   = -5  no of hadronizsation of D-V chains        *
1334 *                   = -6  no of hadronizsation of V-S chains        *
1335 *                   = -7  no of hadronizsation of V-D chains        *
1336 *                   = -8  no of hadronizsation of V-V chains        *
1337 *                   = -9  no of hadronizsation of comb. chains      *
1338 *                                  default:  complete hadronization *
1339 *       sdum   no meaning                                           *
1340 *                                                                   *
1341 *********************************************************************
1342
1343   400 CONTINUE
1344       DO 401 I=1,6
1345          ICHAIN = INT(WHAT(I))
1346          IF ((ICHAIN.LE.-1).AND.(ICHAIN.GE.-9))
1347      &      LHADRO(ABS(ICHAIN)) = .FALSE.
1348   401 CONTINUE
1349       GOTO 10
1350
1351 *********************************************************************
1352 *                                                                   *
1353 *               control card:  codewd = HADRONIZE                   *
1354 *                                                                   *
1355 *           hadronization model and parameter switch                *
1356 *                                                                   *
1357 *       what (1) = 1    hadronization via BAMJET                    *
1358 *                = 2    hadronization via JETSET                    *
1359 *                                                    default: 2     *
1360 *       what (2) = 1..3 parameter set to be used                    *
1361 *                       JETSET: 3 sets available                    *
1362 *                               ( = 3 default JETSET-parameters)    *
1363 *                       BAMJET: 1 set available                     *
1364 *                                                    default: 1     *
1365 *       what (3..6), sdum   no meaning                              *
1366 *                                                                   *
1367 *********************************************************************
1368
1369   410 CONTINUE
1370       IWHAT1 = INT(WHAT(1))
1371       IWHAT2 = INT(WHAT(2))
1372       IF ((IWHAT1.EQ.1).OR.(IWHAT1.EQ.2)) IFRAG(1) = IWHAT1
1373       IF ((IWHAT1.EQ.2).AND.(IWHAT2.GE.1).AND.(IWHAT2.LE.3))
1374      &                                    IFRAG(2) = IWHAT2
1375       GOTO 10
1376
1377 *********************************************************************
1378 *                                                                   *
1379 *               control card:  codewd = POPCORN                     *
1380 *                                                                   *
1381 *  "Popcorn-effect" in fragmentation and diquark breaking diagrams  *
1382 *                                                                   *
1383 *   what (1) = (PDB) frac. of diquark fragmenting directly into     *
1384 *                    baryons (PYTHIA/JETSET fragmentation)          *
1385 *                    (JETSET: = 0. Popcorn mechanism switched off)  *
1386 *                                                    default: 0.5   *
1387 *   what (2) = probability for accepting a diquark breaking         *
1388 *              diagram involving the generation of a u/d quark-     *
1389 *              antiquark pair                        default: 0.0   *
1390 *   what (3) = same a what (2), here for s quark-antiquark pair     *
1391 *                                                    default: 0.0   *
1392 *   what (4..6), sdum   no meaning                                  *
1393 *                                                                   *
1394 *********************************************************************
1395
1396   420 CONTINUE
1397       IF (WHAT(1).GE.0.0D0) PDB = WHAT(1)
1398       IF (WHAT(2).GE.0.0D0) THEN
1399          PDBSEA(1) = WHAT(2)
1400          PDBSEA(2) = WHAT(2)
1401       ENDIF
1402       IF (WHAT(3).GE.0.0D0) PDBSEA(3) = WHAT(3)
1403       DO 421 I=1,8
1404          DBRKA(1,I) = DBRKR(1,I)*PDBSEA(1)/(1.D0-PDBSEA(1))
1405          DBRKA(2,I) = DBRKR(2,I)*PDBSEA(2)/(1.D0-PDBSEA(2))
1406          DBRKA(3,I) = DBRKR(3,I)*PDBSEA(3)/(1.D0-PDBSEA(3))
1407   421 CONTINUE
1408       GOTO 10
1409
1410 *********************************************************************
1411 *                                                                   *
1412 *               control card:  codewd = PARDECAY                    *
1413 *                                                                   *
1414 *      what (1) = 1.  Sigma0/Asigma0 are decaying within JETSET     *
1415 *               = 2.  pion^0 decay after intranucl. cascade         *
1416 *                                                default: no decay  *
1417 *      what (2..6), sdum   no meaning                               *
1418 *                                                                   *
1419 *********************************************************************
1420
1421  430  CONTINUE
1422       IF (WHAT(1).EQ.ONE)  ISIG0 = 1
1423       IF (WHAT(1).EQ.2.0D0) IPI0 = 1
1424       GOTO 10
1425
1426 *********************************************************************
1427 *                                                                   *
1428 *               control card:  codewd = BEAM                        *
1429 *                                                                   *
1430 *              definition of beam parameters                        *
1431 *                                                                   *
1432 *      what (1/2)  > 0 : energy of beam 1/2 (GeV)                   *
1433 *                  < 0 : abs(what(1/2)) energy per charge of        *
1434 *                        beam 1/2 (GeV)                             *
1435 *                  (beam 1 is directed into positive z-direction)   *
1436 *      what (3)    beam crossing angle, defined as 2x angle between *
1437 *                  one beam and the z-axis (micro rad)              *
1438 *      what (4)    angle with x-axis defining the collision plane   *
1439 *      what (5..6), sdum   no meaning                               *
1440 *                                                                   *
1441 *      Note: this card requires previously defined projectile and   *
1442 *            target identities (PROJPAR, TARPAR)                    *
1443 *                                                                   *
1444 *********************************************************************
1445
1446   440 CONTINUE
1447       CALL DT_BEAMPR(WHAT,PPN,1)
1448       EPN    = ZERO
1449       CMENER = ZERO
1450       LEINP  = .TRUE.
1451       GOTO 10
1452
1453 *********************************************************************
1454 *                                                                   *
1455 *               control card:  codewd = LUND-MSTU                   *
1456 *                                                                   *
1457 *          set parameter MSTU in JETSET-common /LUDAT1/             *
1458 *                                                                   *
1459 *       what (1) =  index according to LUND-common block            *
1460 *       what (2) =  new value of MSTU( int(what(1)) )               *
1461 *       what (3), what(4) and what (5), what(6) further             *
1462 *                   parameter in the same way as what (1) and       *
1463 *                   what (2)                                        *
1464 *                        default: default-Lund or corresponding to  *
1465 *                                 the set given in HADRONIZE        *
1466 *                                                                   *
1467 *********************************************************************
1468
1469   450 CONTINUE
1470       IF (WHAT(1).GT.ZERO) THEN
1471          NMSTU = NMSTU+1
1472          IMSTU(NMSTU) = INT(WHAT(1))
1473          MSTUX(NMSTU) = INT(WHAT(2))
1474       ENDIF
1475       IF (WHAT(3).GT.ZERO) THEN
1476          NMSTU = NMSTU+1
1477          IMSTU(NMSTU) = INT(WHAT(3))
1478          MSTUX(NMSTU) = INT(WHAT(4))
1479       ENDIF
1480       IF (WHAT(5).GT.ZERO) THEN
1481          NMSTU = NMSTU+1
1482          IMSTU(NMSTU) = INT(WHAT(5))
1483          MSTUX(NMSTU) = INT(WHAT(6))
1484       ENDIF
1485       GOTO 10
1486
1487 *********************************************************************
1488 *                                                                   *
1489 *               control card:  codewd = LUND-MSTJ                   *
1490 *                                                                   *
1491 *          set parameter MSTJ in JETSET-common /LUDAT1/             *
1492 *                                                                   *
1493 *       what (1) =  index according to LUND-common block            *
1494 *       what (2) =  new value of MSTJ( int(what(1)) )               *
1495 *       what (3), what(4) and what (5), what(6) further             *
1496 *                   parameter in the same way as what (1) and       *
1497 *                   what (2)                                        *
1498 *                        default: default-Lund or corresponding to  *
1499 *                                 the set given in HADRONIZE        *
1500 *                                                                   *
1501 *********************************************************************
1502
1503   451 CONTINUE
1504       IF (WHAT(1).GT.ZERO) THEN
1505          NMSTJ = NMSTJ+1
1506          IMSTJ(NMSTJ) = INT(WHAT(1))
1507          MSTJX(NMSTJ) = INT(WHAT(2))
1508       ENDIF
1509       IF (WHAT(3).GT.ZERO) THEN
1510          NMSTJ = NMSTJ+1
1511          IMSTJ(NMSTJ) = INT(WHAT(3))
1512          MSTJX(NMSTJ) = INT(WHAT(4))
1513       ENDIF
1514       IF (WHAT(5).GT.ZERO) THEN
1515          NMSTJ = NMSTJ+1
1516          IMSTJ(NMSTJ) = INT(WHAT(5))
1517          MSTJX(NMSTJ) = INT(WHAT(6))
1518       ENDIF
1519       GOTO 10
1520
1521 *********************************************************************
1522 *                                                                   *
1523 *               control card:  codewd = LUND-MDCY                   *
1524 *                                                                   *
1525 *  set parameter MDCY(I,1) for particle decays in JETSET-common     *
1526 *                                                      /LUDAT3/     *
1527 *                                                                   *
1528 *       what (1-6) = PDG particle index of particle which should    *
1529 *                    not decay                                      *
1530 *                        default: default-Lund or forced in         *
1531 *                                 DT_INITJS                         *
1532 *                                                                   *
1533 *********************************************************************
1534
1535   452 CONTINUE
1536       DO 4521 I=1,6
1537          IF (WHAT(I).NE.ZERO) THEN
1538
1539             KC = PYCOMP(INT(WHAT(I)))
1540
1541             MDCY(KC,1) = 0
1542          ENDIF
1543  4521 CONTINUE
1544       GOTO 10
1545
1546 *********************************************************************
1547 *                                                                   *
1548 *               control card:  codewd = LUND-PARJ                   *
1549 *                                                                   *
1550 *          set parameter PARJ in JETSET-common /LUDAT1/             *
1551 *                                                                   *
1552 *       what (1) =  index according to LUND-common block            *
1553 *       what (2) =  new value of PARJ( int(what(1)) )               *
1554 *       what (3), what(4) and what (5), what(6) further             *
1555 *                   parameter in the same way as what (1) and       *
1556 *                   what (2)                                        *
1557 *                        default: default-Lund or corresponding to  *
1558 *                                 the set given in HADRONIZE        *
1559 *                                                                   *
1560 *********************************************************************
1561
1562   460 CONTINUE
1563       IF (WHAT(1).NE.ZERO) THEN
1564          NPARJ = NPARJ+1
1565          IPARJ(NPARJ) = INT(WHAT(1))
1566          PARJX(NPARJ) = WHAT(2)
1567       ENDIF
1568       IF (WHAT(3).NE.ZERO) THEN
1569          NPARJ = NPARJ+1
1570          IPARJ(NPARJ) = INT(WHAT(3))
1571          PARJX(NPARJ) = WHAT(4)
1572       ENDIF
1573       IF (WHAT(5).NE.ZERO) THEN
1574          NPARJ = NPARJ+1
1575          IPARJ(NPARJ) = INT(WHAT(5))
1576          PARJX(NPARJ) = WHAT(6)
1577       ENDIF
1578       GOTO 10
1579
1580 *********************************************************************
1581 *                                                                   *
1582 *               control card:  codewd = LUND-PARU                   *
1583 *                                                                   *
1584 *          set parameter PARJ in JETSET-common /LUDAT1/             *
1585 *                                                                   *
1586 *       what (1) =  index according to LUND-common block            *
1587 *       what (2) =  new value of PARU( int(what(1)) )               *
1588 *       what (3), what(4) and what (5), what(6) further             *
1589 *                   parameter in the same way as what (1) and       *
1590 *                   what (2)                                        *
1591 *                        default: default-Lund or corresponding to  *
1592 *                                 the set given in HADRONIZE        *
1593 *                                                                   *
1594 *********************************************************************
1595
1596   470 CONTINUE
1597       IF (WHAT(1).GT.ZERO) THEN
1598          NPARU = NPARU+1
1599          IPARU(NPARU) = INT(WHAT(1))
1600          PARUX(NPARU) = WHAT(2)
1601       ENDIF
1602       IF (WHAT(3).GT.ZERO) THEN
1603          NPARU = NPARU+1
1604          IPARU(NPARU) = INT(WHAT(3))
1605          PARUX(NPARU) = WHAT(4)
1606       ENDIF
1607       IF (WHAT(5).GT.ZERO) THEN
1608          NPARU = NPARU+1
1609          IPARU(NPARU) = INT(WHAT(5))
1610          PARUX(NPARU) = WHAT(6)
1611       ENDIF
1612       GOTO 10
1613
1614 *********************************************************************
1615 *                                                                   *
1616 *               control card:  codewd = OUTLEVEL                    *
1617 *                                                                   *
1618 *                    output control switches                        *
1619 *                                                                   *
1620 *       what (1) =  internal rejection informations  default: 0     *
1621 *       what (2) =  energy-momentum conservation check output       *
1622 *                                                    default: 0     *
1623 *       what (3) =  internal warning messages        default: 0     *
1624 *       what (4..6), sdum    not yet used                           *
1625 *                                                                   *
1626 *********************************************************************
1627
1628   480 CONTINUE
1629       DO 481 K=1,6
1630          IOULEV(K) = INT(WHAT(K))
1631   481 CONTINUE
1632       GOTO 10
1633
1634 *********************************************************************
1635 *                                                                   *
1636 *               control card:  codewd = FRAME                       *
1637 *                                                                   *
1638 *          frame in which final state is given in DTEVT1            *
1639 *                                                                   *
1640 *       what (1) = 1  target rest frame (laboratory)                *
1641 *                = 2  nucleon-nucleon cms                           *
1642 *                                                    default: 1     *
1643 *                                                                   *
1644 *********************************************************************
1645
1646   490 CONTINUE
1647       KFRAME = INT(WHAT(1))
1648       IF ((KFRAME.GE.1).AND.(KFRAME.LE.2)) IFRAME = KFRAME
1649       GOTO 10
1650
1651 *********************************************************************
1652 *                                                                   *
1653 *               control card:  codewd = L-TAG                       *
1654 *                                                                   *
1655 *                        lepton tagger:                             *
1656 *   definition of kinematical cuts for radiated photon and          *
1657 *   outgoing lepton detection in lepton-nucleus interactions        *
1658 *                                                                   *
1659 *       what (1) = y_min                                            *
1660 *       what (2) = y_max                                            *
1661 *       what (3) = Q^2_min                                          *
1662 *       what (4) = Q^2_max                                          *
1663 *       what (5) = theta_min  (Lab)                                 *
1664 *       what (6) = theta_max  (Lab)                                 *
1665 *                                       default: no cuts            *
1666 *       sdum    no meaning                                          *
1667 *                                                                   *
1668 *********************************************************************
1669
1670   500 CONTINUE
1671       YMIN  = WHAT(1)
1672       YMAX  = WHAT(2)
1673       Q2MIN = WHAT(3)
1674       Q2MAX = WHAT(4)
1675       THMIN = WHAT(5)
1676       THMAX = WHAT(6)
1677       GOTO 10
1678
1679 *********************************************************************
1680 *                                                                   *
1681 *               control card:  codewd = L-ETAG                      *
1682 *                                                                   *
1683 *                        lepton tagger:                             *
1684 *       what (1) = min. outgoing lepton energy  (in Lab)            *
1685 *       what (2) = min. photon energy           (in Lab)            *
1686 *       what (3) = max. photon energy           (in Lab)            *
1687 *                                       default: no cuts            *
1688 *       what (2..6), sdum    no meaning                             *
1689 *                                                                   *
1690 *********************************************************************
1691
1692   510 CONTINUE
1693       ELMIN = MAX(WHAT(1),ZERO)
1694       EGMIN = MAX(WHAT(2),ZERO)
1695       EGMAX = MAX(WHAT(3),ZERO)
1696       GOTO 10
1697
1698 *********************************************************************
1699 *                                                                   *
1700 *               control card:  codewd = ECMS-CUT                    *
1701 *                                                                   *
1702 *     what (1) = min. c.m. energy to be sampled                     *
1703 *     what (2) = max. c.m. energy to be sampled                     *
1704 *     what (3) = min x_Bj         to be sampled                     *
1705 *                                       default: no cuts            *
1706 *     what (3..6), sdum    no meaning                               *
1707 *                                                                   *
1708 *********************************************************************
1709
1710   520 CONTINUE
1711       ECMIN  = WHAT(1)
1712       ECMAX  = WHAT(2)
1713       IF (ECMIN.GT.ECMAX) ECMIN = ECMAX
1714       XBJMIN = MAX(WHAT(3),ZERO)
1715       GOTO 10
1716
1717 *********************************************************************
1718 *                                                                   *
1719 *               control card:  codewd = VDM-PAR1                    *
1720 *                                                                   *
1721 *      parameters in gamma-nucleus cross section calculation        *
1722 *                                                                   *
1723 *       what (1) =  Lambda^2                       default: 2.      *
1724 *       what (2)    lower limit in M^2 integration                  *
1725 *                =  1  (3m_pi)^2                                    *
1726 *                =  2  (m_rho0)^2                                   *
1727 *                =  3  (m_phi)^2                   default: 1       *
1728 *       what (3)    upper limit in M^2 integration                  *
1729 *                =  1   s/2                                         *
1730 *                =  2   s/4                                         *
1731 *                =  3   s                          default: 3       *
1732 *       what (4)    CKMT F_2 structure function                     *
1733 *                =  2212  proton                                    *
1734 *                =  100   deuteron                 default: 2212    *
1735 *       what (5)    calculation of gamma-nucleon xsections          *
1736 *                =  1  according to CKMT-parametrization of F_2     *
1737 *                =  2  integrating SIGVP over M^2                   *
1738 *                =  3  using SIGGA                                  *
1739 *                =  4  PHOJET cross sections       default:  4      *
1740 *                                                                   *
1741 *       what (6), sdum    no meaning                                *
1742 *                                                                   *
1743 *********************************************************************
1744
1745   530 CONTINUE
1746       IF (WHAT(1).GE.ZERO) RL2 = WHAT(1)
1747       IF ((WHAT(2).GE.1).AND.(WHAT(2).LE.3)) INTRGE(1) = INT(WHAT(2))
1748       IF ((WHAT(3).GE.1).AND.(WHAT(3).LE.3)) INTRGE(2) = INT(WHAT(3))
1749       IF ((WHAT(4).EQ.2212).OR.(WHAT(4).EQ.100)) IDPDF = INT(WHAT(4))
1750       IF ((WHAT(5).GE.1).AND.(WHAT(5).LE.4)) MODEGA = INT(WHAT(5))
1751       GOTO 10
1752
1753 *********************************************************************
1754 *                                                                   *
1755 *               control card:  codewd = HISTOGRAM                   *
1756 *                                                                   *
1757 *           activate different classes of histograms                *
1758 *                                                                   *
1759 *                                default: no histograms             *
1760 *                                                                   *
1761 *********************************************************************
1762
1763   540 CONTINUE
1764       DO 541 J=1,6
1765          IF ((WHAT(J).GE.100).AND.(WHAT(J).LE.150)) THEN
1766             IHISPP(INT(WHAT(J))-100) = 1
1767          ELSEIF ((ABS(WHAT(J)).GE.200).AND.(ABS(WHAT(J)).LE.250)) THEN
1768             IHISXS(INT(ABS(WHAT(J)))-200) = 1
1769             IF (WHAT(J).LT.ZERO) IXSTBL = 1
1770          ENDIF
1771   541 CONTINUE
1772       GOTO 10
1773
1774 *********************************************************************
1775 *                                                                   *
1776 *               control card:  codewd = XS-TABLE                    *
1777 *                                                                   *
1778 *    output of cross section table for requested interaction        *
1779 *              - particle production deactivated ! -                *
1780 *                                                                   *
1781 *       what (1)      lower energy limit for tabulation             *
1782 *                > 0  Lab. frame                                    *
1783 *                < 0  nucleon-nucleon cms                           *
1784 *       what (2)      upper energy limit for tabulation             *
1785 *                > 0  Lab. frame                                    *
1786 *                < 0  nucleon-nucleon cms                           *
1787 *       what (3) > 0  # of equidistant lin. bins in E               *
1788 *                < 0  # of equidistant log. bins in E               *
1789 *       what (4)      lower limit of particle virtuality (photons)  *
1790 *       what (5)      upper limit of particle virtuality (photons)  *
1791 *       what (6) > 0  # of equidistant lin. bins in Q^2             *
1792 *                < 0  # of equidistant log. bins in Q^2             *
1793 *                                                                   *
1794 *********************************************************************
1795
1796   550 CONTINUE
1797       IF (WHAT(1).EQ.99999.0D0) THEN
1798          IRATIO = INT(WHAT(2))
1799          GOTO 10
1800       ENDIF
1801       CMENER = ABS(WHAT(2))
1802       IF (.NOT.LXSTAB) THEN
1803
1804          CALL BERTTP
1805          CALL INCINI
1806
1807       ENDIF
1808       IF ((.NOT.LXSTAB).OR.(CMENER.NE.CMEOLD)) THEN
1809          CMEOLD = CMENER
1810          IF (WHAT(2).GT.ZERO)
1811      &      CMENER = SQRT(2.0D0*AAM(1)**2+2.0D0*WHAT(2)*AAM(1))
1812          EPN = ZERO
1813          PPN = ZERO
1814 C        WRITE(LOUT,*) 'CMENER = ',CMENER
1815          CALL DT_LTINI(IJPROJ,IJTARG,EPN,PPN,CMENER,1)
1816          CALL DT_PHOINI
1817       ENDIF
1818       CALL DT_XSTABL(WHAT,IXSQEL,IRATIO)
1819       IXSQEL = 0
1820       LXSTAB = .TRUE.
1821       GOTO 10
1822
1823 *********************************************************************
1824 *                                                                   *
1825 *               control card:  codewd = GLAUB-PAR                   *
1826 *                                                                   *
1827 *                parameters in Glauber-formalism                    *
1828 *                                                                   *
1829 *    what (1)  # of nucleon configurations sampled in integration   *
1830 *              over nuclear desity                default: 1000     *
1831 *    what (2)  # of bins for integration over impact-parameter and  *
1832 *              for profile-function calculation   default: 49       *
1833 *    what (3)  = 1 calculation of tot., el. and qel. cross sections *
1834 *                                                 default: 0        *
1835 *    what (4)  = 1   read pre-calculated impact-parameter distrib.  *
1836 *                    from "sdum".glb                                *
1837 *              =-1   dump pre-calculated impact-parameter distrib.  *
1838 *                    into "sdum".glb                                *
1839 *              = 100 read pre-calculated impact-parameter distrib.  *
1840 *                    for variable projectile/target/energy runs     *
1841 *                    from "sdum".glb                                *
1842 *                                                 default: 0        *
1843 *    what (5..6)   no meaning                                       *
1844 *    sdum      if |what (4)| = 1 name of in/output-file (sdum.glb)  *
1845 *                                                                   *
1846 *********************************************************************
1847
1848   560 CONTINUE
1849       IF (WHAT(1).GT.ZERO) JSTATB = INT(WHAT(1))
1850       IF (WHAT(2).GT.ZERO) JBINSB = INT(WHAT(2))
1851       IF (WHAT(3).EQ.ONE) LPROD = .FALSE.
1852       IF ((ABS(WHAT(4)).EQ.ONE).OR.(WHAT(4).EQ.100)) THEN
1853          IOGLB = INT(WHAT(4))
1854          CGLB  = SDUM
1855       ENDIF
1856       GOTO 10
1857
1858 *********************************************************************
1859 *                                                                   *
1860 *               control card:  codewd = GLAUB-INI                   *
1861 *                                                                   *
1862 *             pre-initialization of profile function                *
1863 *                                                                   *
1864 *       what (1)      lower energy limit for initialization         *
1865 *                > 0  Lab. frame                                    *
1866 *                < 0  nucleon-nucleon cms                           *
1867 *       what (2)      upper energy limit for initialization         *
1868 *                > 0  Lab. frame                                    *
1869 *                < 0  nucleon-nucleon cms                           *
1870 *       what (3) > 0  # of equidistant lin. bins in E               *
1871 *                < 0  # of equidistant log. bins in E               *
1872 *       what (4)      maximum projectile mass number for which the  *
1873 *                     Glauber data are initialized for each         *
1874 *                     projectile mass number                        *
1875 *                     (if <= mass given with the PROJPAR-card)      *
1876 *                                              default: 18          *
1877 *       what (5)      steps in mass number starting from what (4)   *
1878 *                     up to mass number defined with PROJPAR-card   *
1879 *                     for which Glauber data are initialized        *
1880 *                                              default: 5           *
1881 *       what (6)      no meaning                                    *
1882 *       sdum          no meaning                                    *
1883 *                                                                   *
1884 *********************************************************************
1885
1886   565 CONTINUE
1887       IOGLB = -100
1888       CALL DT_GLBINI(WHAT)
1889       GOTO 10
1890
1891 *********************************************************************
1892 *                                                                   *
1893 *               control card:  codewd = VDM-PAR2                    *
1894 *                                                                   *
1895 *      parameters in gamma-nucleus cross section calculation        *
1896 *                                                                   *
1897 *      what (1) = 0 no suppression of shadowing by direct photon    *
1898 *                   processes                                       *
1899 *               = 1 suppression ..                   default: 1     *
1900 *      what (2) = 0 no suppression of shadowing by anomalous        *
1901 *                   component if photon-F_2                         *
1902 *               = 1 suppression ..                   default: 1     *
1903 *      what (3) = 0 no suppression of shadowing by coherence        *
1904 *                   length of the photon                            *
1905 *               = 1 suppression ..                   default: 1     *
1906 *      what (4) = 1 longitudinal polarized photons are taken into   *
1907 *                   account                                         *
1908 *                   eps*R*Q^2/M^2 = what(4)*Q^2/M^2  default: 0     *
1909 *      what (5..6), sdum    no meaning                              *
1910 *                                                                   *
1911 *********************************************************************
1912
1913   570 CONTINUE
1914       IF ((WHAT(1).EQ.ZERO).OR.(WHAT(1).EQ.ONE)) ISHAD(1) = INT(WHAT(1))
1915       IF ((WHAT(2).EQ.ZERO).OR.(WHAT(2).EQ.ONE)) ISHAD(2) = INT(WHAT(2))
1916       IF ((WHAT(3).EQ.ZERO).OR.(WHAT(3).EQ.ONE)) ISHAD(3) = INT(WHAT(3))
1917       EPSPOL  = WHAT(4)
1918       GOTO 10
1919
1920 *********************************************************************
1921 *                                                                   *
1922 *               control card:  XS-QELPRO                            *
1923 *                                                                   *
1924 *     what (1..6), sdum    no meaning                               *
1925 *                                                                   *
1926 *********************************************************************
1927
1928   580 CONTINUE
1929       IXSQEL = ABS(WHAT(1))
1930       GOTO 10
1931
1932 *********************************************************************
1933 *                                                                   *
1934 *               control card:  RNDMINIT                             *
1935 *                                                                   *
1936 *           initialization of random number generator               *
1937 *                                                                   *
1938 *     what (1..4)    values for initialization (= 1..168)           *
1939 *     what (5..6), sdum    no meaning                               *
1940 *                                                                   *
1941 *********************************************************************
1942
1943   590 CONTINUE
1944       IF ((WHAT(1).LT.1.0D0).OR.(WHAT(1).GT.168.0D0)) THEN
1945          NA1 = 22
1946       ELSE
1947          NA1 = WHAT(1)
1948       ENDIF
1949       IF ((WHAT(2).LT.1.0D0).OR.(WHAT(2).GT.168.0D0)) THEN
1950          NA2 = 54
1951       ELSE
1952          NA2 = WHAT(2)
1953       ENDIF
1954       IF ((WHAT(3).LT.1.0D0).OR.(WHAT(3).GT.168.0D0)) THEN
1955          NA3 = 76
1956       ELSE
1957          NA3 = WHAT(3)
1958       ENDIF
1959       IF ((WHAT(4).LT.1.0D0).OR.(WHAT(4).GT.168.0D0)) THEN
1960          NA4 = 92
1961       ELSE
1962          NA4 = WHAT(4)
1963       ENDIF
1964       CALL DT_RNDMST(NA1,NA2,NA3,NA4)
1965       GOTO 10
1966
1967 *********************************************************************
1968 *                                                                   *
1969 *               control card:  codewd = LEPTO-CUT                   *
1970 *                                                                   *
1971 *          set parameter CUT in LEPTO-common /LEPTOU/               *
1972 *                                                                   *
1973 *       what (1) =  index in CUT-array                              *
1974 *       what (2) =  new value of CUT( int(what(1)) )                *
1975 *       what (3), what(4) and what (5), what(6) further             *
1976 *                   parameter in the same way as what (1) and       *
1977 *                   what (2)                                        *
1978 *                        default: default-LEPTO parameters          *
1979 *                                                                   *
1980 *********************************************************************
1981
1982   600 CONTINUE
1983       IF (WHAT(1).GT.ZERO) CUT(INT(WHAT(1))) = WHAT(2)
1984       IF (WHAT(3).GT.ZERO) CUT(INT(WHAT(3))) = WHAT(4)
1985       IF (WHAT(5).GT.ZERO) CUT(INT(WHAT(5))) = WHAT(6)
1986       GOTO 10
1987
1988 *********************************************************************
1989 *                                                                   *
1990 *               control card:  codewd = LEPTO-LST                   *
1991 *                                                                   *
1992 *          set parameter LST in LEPTO-common /LEPTOU/               *
1993 *                                                                   *
1994 *       what (1) =  index in LST-array                              *
1995 *       what (2) =  new value of LST( int(what(1)) )                *
1996 *       what (3), what(4) and what (5), what(6) further             *
1997 *                   parameter in the same way as what (1) and       *
1998 *                   what (2)                                        *
1999 *                        default: default-LEPTO parameters          *
2000 *                                                                   *
2001 *********************************************************************
2002
2003   610 CONTINUE
2004       IF (WHAT(1).GT.ZERO) LST(INT(WHAT(1))) = INT(WHAT(2))
2005       IF (WHAT(3).GT.ZERO) LST(INT(WHAT(3))) = INT(WHAT(4))
2006       IF (WHAT(5).GT.ZERO) LST(INT(WHAT(5))) = INT(WHAT(6))
2007       GOTO 10
2008
2009 *********************************************************************
2010 *                                                                   *
2011 *               control card:  codewd = LEPTO-PARL                  *
2012 *                                                                   *
2013 *          set parameter PARL in LEPTO-common /LEPTOU/              *
2014 *                                                                   *
2015 *       what (1) =  index in PARL-array                             *
2016 *       what (2) =  new value of PARL( int(what(1)) )               *
2017 *       what (3), what(4) and what (5), what(6) further             *
2018 *                   parameter in the same way as what (1) and       *
2019 *                   what (2)                                        *
2020 *                        default: default-LEPTO parameters          *
2021 *                                                                   *
2022 *********************************************************************
2023
2024   620 CONTINUE
2025       IF (WHAT(1).GT.ZERO) PARL(INT(WHAT(1))) = WHAT(2)
2026       IF (WHAT(3).GT.ZERO) PARL(INT(WHAT(3))) = WHAT(4)
2027       IF (WHAT(5).GT.ZERO) PARL(INT(WHAT(5))) = WHAT(6)
2028       GOTO 10
2029
2030 *********************************************************************
2031 *                                                                   *
2032 *               control card:  codewd = START                       *
2033 *                                                                   *
2034 *       what (1) =   number of events                default: 100.  *
2035 *       what (2) = 0 Glauber initialization follows                 *
2036 *                = 1 Glauber initialization supressed, fitted       *
2037 *                    results are used instead                       *
2038 *                    (this does not apply if emulsion-treatment     *
2039 *                     is requested)                                 *
2040 *                = 2 Glauber initialization is written to           *
2041 *                    output-file shmakov.out                        *
2042 *                = 3 Glauber initialization is read from input-file *
2043 *                    shmakov.out                     default: 0     *
2044 *       what (3..6)  no meaning                                     *
2045 *       what (3..6)  no meaning                                     *
2046 *                                                                   *
2047 *********************************************************************
2048
2049   630 CONTINUE
2050
2051 * check for cross-section table output only
2052       IF (LXSTAB) STOP
2053
2054       NCASES = INT(WHAT(1))
2055       IF (NCASES.LE.0) NCASES = 100
2056       IGLAU = INT(WHAT(2))
2057       IF ((IGLAU.NE.1).AND.(IGLAU.NE.2).AND.(IGLAU.NE.3))
2058      &                                            IGLAU = 0
2059
2060       NPMASS = IP
2061       NPCHAR = IPZ
2062       NTMASS = IT
2063       NTCHAR = ITZ
2064       IDP    = IJPROJ
2065       IDT    = IJTARG
2066       IF (IDP.LE.0) IDP = 1
2067 * muon neutrinos: temporary (missing index)
2068 * (new patch in projpar: therefore the following this is probably not
2069 *  necessary anymore..)
2070 C     IF (IDP.EQ.26) IDP = 5
2071 C     IF (IDP.EQ.27) IDP = 6
2072
2073 * redefine collision energy
2074       IF (LEINP) THEN
2075          IF (ABS(VAREHI).GT.ZERO) THEN
2076             PDUM = ZERO
2077             IF (VARELO.LT.EHADLO) VARELO = EHADLO
2078             CALL DT_LTINI(IDP,IDT,VARELO,PDUM,VARCLO,1)
2079             PDUM = ZERO
2080             CALL DT_LTINI(IDP,IDT,VAREHI,PDUM,VARCHI,1)
2081          ENDIF
2082          CALL DT_LTINI(IDP,IDT,EPN,PPN,CMENER,1)
2083       ELSE
2084          WRITE(LOUT,1003)
2085  1003    FORMAT(1X,'INIT:   collision energy not defined!',/,
2086      &          1X,'              -program stopped-      ')
2087          STOP
2088       ENDIF
2089
2090 * switch off evaporation (even if requested) if central coll. requ.
2091       IF ((ICENTR.EQ.-1).OR.(ICENTR.GT.0).OR.(XSFRAC.LT.0.5D0)) THEN
2092          IF (LEVPRT) THEN
2093             WRITE(LOUT,1004)
2094  1004       FORMAT(1X,/,'Warning!  Evaporation request rejected since',
2095      &             ' central collisions forced.')
2096             LEVPRT = .FALSE.
2097             LDEEXG = .FALSE.
2098             LHEAVY = .FALSE.
2099          ENDIF
2100       ENDIF
2101
2102 * initialization of evaporation-module
2103
2104 *  initialize evaporation if the code is not used as Fluka event generator
2105       IF (ITRSPT.NE.1) THEN
2106 *         CALL BERTTP
2107 *         CALL INCINI
2108       ENDIF
2109       IF (LEVPRT) LHEAVY = .TRUE.
2110
2111
2112 * save the default JETSET-parameter
2113       CALL DT_JSPARA(0)
2114
2115 * force use of phojet for g-A
2116       IF ((IDP.EQ.7).AND.(MCGENE.NE.3)) MCGENE = 2
2117 * initialization of nucleon-nucleon event generator
2118       IF (MCGENE.EQ.2) CALL DT_PHOINI
2119 * initialization of LEPTO event generator
2120       IF (MCGENE.EQ.3) THEN
2121
2122          STOP ' This version does not contain LEPTO !'
2123
2124       ENDIF
2125
2126 * initialization of quasi-elastic neutrino scattering
2127       IF (MCGENE.EQ.4) THEN
2128          IF (IJPROJ.EQ.5) THEN
2129             NEUTYP = 1
2130          ELSEIF (IJPROJ.EQ.6) THEN
2131             NEUTYP = 2
2132          ELSEIF (IJPROJ.EQ.135) THEN
2133             NEUTYP = 3
2134          ELSEIF (IJPROJ.EQ.136) THEN
2135             NEUTYP = 4
2136          ELSEIF (IJPROJ.EQ.133) THEN
2137             NEUTYP = 5
2138          ELSEIF (IJPROJ.EQ.134) THEN
2139             NEUTYP = 6
2140          ENDIF
2141       ENDIF
2142
2143 * normalize fractions of emulsion components
2144       IF (NCOMPO.GT.0) THEN
2145          SUMFRA = ZERO
2146          DO 491 I=1,NCOMPO
2147             SUMFRA = SUMFRA+EMUFRA(I)
2148   491    CONTINUE
2149          IF (SUMFRA.GT.ZERO) THEN
2150             DO 492 I=1,NCOMPO
2151                EMUFRA(I) = EMUFRA(I)/SUMFRA
2152   492       CONTINUE
2153          ENDIF
2154       ENDIF
2155
2156 * disallow Cronin's multiple scattering for nucleus-nucleus interactions
2157       IF ((IP.GT.1).AND.(MKCRON.GT.0)) THEN
2158          WRITE(LOUT,1005)
2159  1005    FORMAT(/,1X,'INIT:  multiple scattering disallowed',/)
2160          MKCRON = 0
2161       ENDIF
2162
2163 * initialization of Glauber-formalism (moved to xAEVT, sr 26.3.96)
2164 C     IF (NCOMPO.LE.0) THEN
2165 C        CALL DT_SHMAKI(IP,IPZ,IT,ITZ,IDP,PPN,IGLAU)
2166 C     ELSE
2167 C        DO 493 I=1,NCOMPO
2168 C           CALL DT_SHMAKI(IP,IPZ,IEMUMA(I),IEMUCH(I),IDP,PPN,0)
2169 C 493    CONTINUE
2170 C     ENDIF
2171
2172 * pre-tabulation of elastic cross-sections
2173       CALL DT_SIGTBL(JDUM,JDUM,DUM,DUM,-1)
2174
2175       CALL DT_XTIME
2176
2177       RETURN
2178
2179 *********************************************************************
2180 *                                                                   *
2181 *               control card:  codewd = STOP                        *
2182 *                                                                   *
2183 *               stop of the event generation                        *
2184 *                                                                   *
2185 *       what (1..6)  no meaning                                     *
2186 *                                                                   *
2187 *********************************************************************
2188
2189  9999 CONTINUE
2190       WRITE(LOUT,9000)
2191  9000 FORMAT(1X,'---> unexpected end of input !')
2192
2193   640 CONTINUE
2194       STOP
2195
2196       END
2197 *
2198 *===kkinc==============================================================*
2199 *
2200 CDECK  ID>, DT_KKINC
2201       SUBROUTINE DT_KKINC(NPMASS,NPCHAR,NTMASS,NTCHAR,IDP,EPN,KKMAT,
2202      &                                                         IREJ)
2203
2204 ************************************************************************
2205 * Treatment of complete nucleus-nucleus or hadron-nucleus scattering   *
2206 * This subroutine is an update of the previous version written         *
2207 * by J. Ranft/ H.-J. Moehring.                                         *
2208 * This version dated 19.11.95 is written by S. Roesler                 *
2209 ************************************************************************
2210
2211       IMPLICIT DOUBLE PRECISION (A-H,O-Z)
2212       SAVE
2213
2214       PARAMETER ( LINP = 5 ,
2215      &            LOUT = 6 ,
2216      &            LDAT = 9 )
2217
2218       PARAMETER (ZERO=0.0D0,ONE=1.0D0,TINY5=1.0D-5,
2219      &           TINY2=1.0D-2,TINY3=1.0D-3)
2220
2221       LOGICAL LFZC
2222
2223 * event history
2224
2225       PARAMETER (NMXHKK=200000)
2226
2227       COMMON /DTEVT1/ NHKK,NEVHKK,ISTHKK(NMXHKK),IDHKK(NMXHKK),
2228      &                JMOHKK(2,NMXHKK),JDAHKK(2,NMXHKK),
2229      &                PHKK(5,NMXHKK),VHKK(4,NMXHKK),WHKK(4,NMXHKK)
2230 * extended event history
2231       COMMON /DTEVT2/ IDRES(NMXHKK),IDXRES(NMXHKK),NOBAM(NMXHKK),
2232      &                IDBAM(NMXHKK),IDCH(NMXHKK),NPOINT(10),
2233      &                IHIST(2,NMXHKK)
2234 * particle properties (BAMJET index convention)
2235       CHARACTER*8  ANAME
2236       COMMON /DTPART/ ANAME(210),AAM(210),GA(210),TAU(210),
2237      &                IICH(210),IIBAR(210),K1(210),K2(210)
2238 * properties of interacting particles
2239       COMMON /DTPRTA/ IT,ITZ,IP,IPZ,IJPROJ,IBPROJ,IJTARG,IBTARG
2240 * Lorentz-parameters of the current interaction
2241       COMMON /DTLTRA/ GACMS(2),BGCMS(2),GALAB,BGLAB,BLAB,
2242      &                UMO,PPCM,EPROJ,PPROJ
2243 * flags for input different options
2244       LOGICAL LEMCCK,LHADRO,LSEADI,LEVAPO
2245       COMMON /DTFLG1/ IFRAG(2),IRESCO,IMSHL,IRESRJ,IOULEV(6),
2246      &                LEMCCK,LHADRO(0:9),LSEADI,LEVAPO,IFRAME,ITRSPT
2247 * flags for particle decays
2248       COMMON /DTFRPA/ MSTUX(20),PARUX(20),MSTJX(20),PARJX(20),
2249      &                IMSTU(20),IPARU(20),IMSTJ(20),IPARJ(20),
2250      &                NMSTU,NPARU,NMSTJ,NPARJ,PDB,PDBSEA(3),ISIG0,IPI0
2251 * cuts for variable energy runs
2252       COMMON /DTVARE/ VARELO,VAREHI,VARCLO,VARCHI
2253 * Glauber formalism: flags and parameters for statistics
2254       LOGICAL LPROD
2255       CHARACTER*8 CGLB
2256       COMMON /DTGLGP/ JSTATB,JBINSB,CGLB,IOGLB,LPROD
2257
2258       DIMENSION WHAT(6)
2259
2260       IREJ  = 0
2261       ILOOP = 0
2262   100 CONTINUE
2263       IF (ILOOP.EQ.4) THEN
2264          WRITE(LOUT,1000) NEVHKK
2265  1000    FORMAT(1X,'KKINC: event ',I8,' rejected!')
2266          GOTO 9999
2267       ENDIF
2268       ILOOP = ILOOP+1
2269
2270 * variable energy-runs, recalculate parameters for LT's
2271       IF ((ABS(VAREHI).GT.ZERO).OR.(IOGLB.EQ.100)) THEN
2272          PDUM = ZERO
2273          CDUM = ZERO
2274          CALL DT_LTINI(IDP,1,EPN,PDUM,CDUM,1)
2275       ENDIF
2276       IF (EPN.GT.EPROJ) THEN
2277          WRITE(LOUT,'(A,E9.3,2A,E9.3,A)')
2278      &      ' Requested energy (',EPN,'GeV) exceeds',
2279      &      ' initialization energy (',EPROJ,'GeV) !'
2280          STOP
2281       ENDIF
2282
2283 * re-initialize /DTPRTA/
2284       IP  = NPMASS
2285       IPZ = NPCHAR
2286       IT  = NTMASS
2287       ITZ = NTCHAR
2288       IJPROJ = IDP
2289       IBPROJ = IIBAR(IJPROJ)
2290
2291 * calculate nuclear potentials (common /DTNPOT/)
2292       CALL DT_NCLPOT(IPZ,IP,ITZ,IT,ZERO,ZERO,0)
2293
2294 * initialize treatment for residual nuclei
2295       CALL DT_RESNCL(EPN,NLOOP,1)
2296
2297 * sample hadron/nucleus-nucleus interaction
2298       CALL DT_KKEVNT(KKMAT,IREJ1)
2299       IF (IREJ1.GT.0) THEN
2300          IF (IOULEV(1).GT.0) WRITE(LOUT,*) 'rejected 1 in KKINC'
2301          GOTO 9999
2302       ENDIF
2303
2304       IF ((NPMASS.GT.1).OR.(NTMASS.GT.1)) THEN
2305
2306 * intranuclear cascade of final state particles for KTAUGE generations
2307 * of secondaries
2308          CALL DT_FOZOCA(LFZC,IREJ1)
2309          IF (IREJ1.GT.0) THEN
2310             IF (IOULEV(1).GT.0) WRITE(LOUT,*) 'rejected 2 in KKINC'
2311             GOTO 9999
2312          ENDIF
2313
2314 * baryons unable to escape the nuclear potential are treated as
2315 * excited nucleons (ISTHKK=15,16)
2316          CALL DT_SCN4BA
2317
2318 * decay of resonances produced in intranuclear cascade processes
2319 **sr 15-11-95 should be obsolete
2320 C        IF (LFZC) CALL DT_DECAY1
2321
2322   101    CONTINUE
2323 * treatment of residual nuclei
2324          CALL DT_RESNCL(EPN,NLOOP,2)
2325
2326 * evaporation / fission / fragmentation
2327 * (if intranuclear cascade was sampled only)
2328          IF (LFZC) THEN
2329             CALL DT_FICONF(IJPROJ,IP,IPZ,IT,ITZ,NLOOP,IREJ1)
2330             IF (IREJ1.GT.1) GOTO 101
2331             IF (IREJ1.EQ.1) GOTO 100
2332          ENDIF
2333
2334       ENDIF
2335
2336 * transform finale state into Lab.
2337       IFLAG = 2
2338       CALL DT_BEAMPR(WHAT,DUM,IFLAG)
2339       IF ((IFRAME.EQ.1).AND.(IFLAG.EQ.-1)) CALL DT_LT2LAB
2340
2341       IF (IPI0.EQ.1) CALL DT_DECPI0
2342
2343 C     IF (NEVHKK.EQ.5) CALL DT_EVTOUT(4)
2344
2345       RETURN
2346  9999 CONTINUE
2347       IREJ = 1
2348       RETURN
2349       END
2350 *
2351 *===defaul=============================================================*
2352 *
2353 CDECK  ID>, DT_DEFAUL
2354       SUBROUTINE DT_DEFAUL(EPN,PPN)
2355
2356 ************************************************************************
2357 * Variables are set to default values.                                 *
2358 * This version dated 8.5.95 is written by S. Roesler.                  *
2359 ************************************************************************
2360
2361       IMPLICIT DOUBLE PRECISION (A-H,O-Z)
2362       SAVE
2363       PARAMETER (ZERO=0.0D0,ONE=1.0D0,TINY10=1.0D-10)
2364       PARAMETER (TWOPI  = 6.283185307179586454D+00)
2365
2366 * particle properties (BAMJET index convention)
2367       CHARACTER*8  ANAME
2368       COMMON /DTPART/ ANAME(210),AAM(210),GA(210),TAU(210),
2369      &                IICH(210),IIBAR(210),K1(210),K2(210)
2370 * nuclear potential
2371       LOGICAL LFERMI
2372       COMMON /DTNPOT/ PFERMP(2),PFERMN(2),FERMOD,
2373      &                EBINDP(2),EBINDN(2),EPOT(2,210),
2374      &                ETACOU(2),ICOUL,LFERMI
2375 * interface HADRIN-DPM
2376       COMMON /HNTHRE/ EHADTH,EHADLO,EHADHI,INTHAD,IDXTA
2377 * central particle production, impact parameter biasing
2378       COMMON /DTIMPA/ BIMIN,BIMAX,XSFRAC,ICENTR
2379 * properties of interacting particles
2380       COMMON /DTPRTA/ IT,ITZ,IP,IPZ,IJPROJ,IBPROJ,IJTARG,IBTARG
2381 * properties of photon/lepton projectiles
2382       COMMON /DTGPRO/ VIRT,PGAMM(4),PLEPT0(4),PLEPT1(4),PNUCL(4),IDIREC
2383
2384       PARAMETER (NCOMPX=20,NEB=8,NQB= 5,KSITEB=50)
2385
2386 * emulsion treatment
2387       COMMON /DTCOMP/ EMUFRA(NCOMPX),IEMUMA(NCOMPX),IEMUCH(NCOMPX),
2388      &                NCOMPO,IEMUL
2389 * parameter for intranuclear cascade
2390       LOGICAL LPAULI
2391       COMMON /DTFOTI/ TAUFOR,KTAUGE,ITAUVE,INCMOD,LPAULI
2392 * various options for treatment of partons (DTUNUC 1.x)
2393 * (chain recombination, Cronin,..)
2394       LOGICAL LCO2CR,LINTPT
2395       COMMON /DTCHAI/ SEASQ,CRONCO,CUTOF,MKCRON,ISICHA,IRECOM,
2396      &                LCO2CR,LINTPT
2397 * threshold values for x-sampling (DTUNUC 1.x)
2398       COMMON /DTXCUT/ XSEACU,UNON,UNOM,UNOSEA,CVQ,CDQ,CSEA,SSMIMA,
2399      &                SSMIMQ,VVMTHR
2400 * flags for input different options
2401       LOGICAL LEMCCK,LHADRO,LSEADI,LEVAPO
2402       COMMON /DTFLG1/ IFRAG(2),IRESCO,IMSHL,IRESRJ,IOULEV(6),
2403      &                LEMCCK,LHADRO(0:9),LSEADI,LEVAPO,IFRAME,ITRSPT
2404 * n-n cross section fluctuations
2405       PARAMETER (NBINS = 1000)
2406       COMMON /DTXSFL/ FLUIXX(NBINS),IFLUCT
2407 * flags for particle decays
2408       COMMON /DTFRPA/ MSTUX(20),PARUX(20),MSTJX(20),PARJX(20),
2409      &                IMSTU(20),IPARU(20),IMSTJ(20),IPARJ(20),
2410      &                NMSTU,NPARU,NMSTJ,NPARJ,PDB,PDBSEA(3),ISIG0,IPI0
2411 * diquark-breaking mechanism
2412       COMMON /DTDIQB/ DBRKR(3,8),DBRKA(3,8),CHAM1,CHAM3,CHAB1,CHAB3
2413 * nucleon-nucleon event-generator
2414       CHARACTER*8 CMODEL
2415       LOGICAL LPHOIN
2416       COMMON /DTMODL/ CMODEL(4),ELOJET,MCGENE,LPHOIN
2417 * flags for diffractive interactions (DTUNUC 1.x)
2418       COMMON /DTFLG3/ ISINGD,IDOUBD,IFLAGD,IDIFF
2419 * VDM parameter for photon-nucleus interactions
2420       COMMON /DTVDMP/ RL2,EPSPOL,INTRGE(2),IDPDF,MODEGA,ISHAD(3)
2421 * Glauber formalism: flags and parameters for statistics
2422       LOGICAL LPROD
2423       CHARACTER*8 CGLB
2424       COMMON /DTGLGP/ JSTATB,JBINSB,CGLB,IOGLB,LPROD
2425 * kinematical cuts for lepton-nucleus interactions
2426       COMMON /DTLCUT/ ECMIN,ECMAX,XBJMIN,ELMIN,EGMIN,EGMAX,YMIN,YMAX,
2427      &                Q2MIN,Q2MAX,THMIN,THMAX,Q2LI,Q2HI,ECMLI,ECMHI
2428 * flags for activated histograms
2429       COMMON /DTHIS3/ IHISPP(50),IHISXS(50),IXSTBL
2430 * cuts for variable energy runs
2431       COMMON /DTVARE/ VARELO,VAREHI,VARCLO,VARCHI
2432 * parameters for hA-diffraction
2433       COMMON /DTDIHA/ DIBETA,DIALPH
2434 * LEPTO
2435       REAL RPPN
2436       COMMON /LEPTOI/ RPPN,LEPIN,INTER
2437 * steering flags for qel neutrino scattering modules
2438       COMMON /QNEUTO/ DSIGSU,DSIGMC,NDSIG,NEUTYP,NEUDEC
2439 * event flag
2440       COMMON /DTEVNO/ NEVENT,ICASCA
2441
2442       DATA POTMES /0.002D0/
2443
2444 * common /DTNPOT/
2445       DO 10 I=1,2
2446          PFERMP(I) = ZERO
2447          PFERMN(I) = ZERO
2448          EBINDP(I) = ZERO
2449          EBINDN(I) = ZERO
2450          DO 11 J=1,210
2451             EPOT(I,J) = ZERO
2452    11    CONTINUE
2453 * nucleus independent meson potential
2454          EPOT(I,13) = POTMES
2455          EPOT(I,14) = POTMES
2456          EPOT(I,15) = POTMES
2457          EPOT(I,16) = POTMES
2458          EPOT(I,23) = POTMES
2459          EPOT(I,24) = POTMES
2460          EPOT(I,25) = POTMES
2461    10 CONTINUE
2462 **sr 7.4.98: changed after corrected B-sampling
2463 C     FERMOD    = 0.55D0
2464       FERMOD    = 0.68D0
2465       ETACOU(1) = ZERO
2466       ETACOU(2) = ZERO
2467       ICOUL     = 1
2468       LFERMI    = .TRUE.
2469
2470 * common /HNTHRE/
2471       EHADTH = -99.0D0
2472       EHADLO = 4.06D0
2473       EHADHI = 6.0D0
2474       INTHAD = 1
2475       IDXTA  = 2
2476
2477 * common /DTIMPA/
2478       ICENTR = 0
2479       BIMIN  = ZERO
2480       BIMAX  = 1.0D10
2481       XSFRAC = 1.0D0
2482
2483 * common /DTPRTA/
2484       IP  = 1
2485       IPZ = 1
2486       IT  = 1
2487       ITZ = 1
2488       IJPROJ = 1
2489       IBPROJ = 1
2490       IJTARG = 1
2491       IBTARG = 1
2492 * common /DTGPRO/
2493       VIRT = ZERO
2494       DO 14 I=1,4
2495          PGAMM(I)  = ZERO
2496          PLEPT0(I) = ZERO
2497          PLEPT1(I) = ZERO
2498          PNUCL(I)  = ZERO
2499    14 CONTINUE
2500       IDIREC   = 0
2501
2502 * common /DTFOTI/
2503 **sr 7.4.98: changed after corrected B-sampling
2504 C     TAUFOR = 4.4D0
2505       TAUFOR = 3.1D0
2506       KTAUGE = 25
2507       ITAUVE = 1
2508       INCMOD = 1
2509       LPAULI = .TRUE.
2510
2511 * common /DTCHAI/
2512       SEASQ  = ONE
2513       MKCRON = 1
2514       CRONCO = 0.64D0
2515       ISICHA = 0
2516       CUTOF  = 100.0D0
2517       LCO2CR = .FALSE.
2518       IRECOM = 1
2519       LINTPT = .TRUE.
2520
2521 * common /DTXCUT/
2522 *  definition of soft quark distributions
2523       XSEACU = 0.05D0
2524       UNON   = 2.0D0
2525       UNOM   = 1.5D0
2526       UNOSEA = 5.0D0
2527 *  cutoff parameters for x-sampling
2528       CVQ    = 1.0D0
2529       CDQ    = 2.0D0
2530 C     CSEA   = 0.3D0
2531       CSEA   = 0.1D0
2532       SSMIMA = 1.2D0
2533       SSMIMQ = SSMIMA**2
2534       VVMTHR = 2.0D0
2535
2536 * common /DTXSFL/
2537       IFLUCT = 0
2538
2539 * common /DTFRPA/
2540       PDB = 0.15D0
2541       PDBSEA(1) = 0.0D0
2542       PDBSEA(2) = 0.0D0
2543       PDBSEA(3) = 0.0D0
2544       ISIG0 = 0
2545       IPI0  = 0
2546       NMSTU = 0
2547       NPARU = 0
2548       NMSTJ = 0
2549       NPARJ = 0
2550
2551 * common /DTDIQB/
2552       DO 15 I=1,8
2553          DBRKR(1,I) = 5.0D0
2554          DBRKR(2,I) = 5.0D0
2555          DBRKR(3,I) = 10.0D0
2556          DBRKA(1,I) = ZERO
2557          DBRKA(2,I) = ZERO
2558          DBRKA(3,I) = ZERO
2559    15 CONTINUE
2560       CHAM1 = 0.2D0
2561       CHAM3 = 0.5D0
2562       CHAB1 = 0.7D0
2563       CHAB3 = 1.0D0
2564
2565 * common /DTFLG3/
2566       ISINGD = 0
2567       IDOUBD = 0
2568       IFLAGD = 0
2569       IDIFF  = 0
2570
2571 * common /DTMODL/
2572       MCGENE    = 2
2573       CMODEL(1) = 'DTUNUC  '
2574       CMODEL(2) = 'PHOJET  '
2575       CMODEL(3) = 'LEPTO   '
2576       CMODEL(4) = 'QNEUTRIN'
2577       LPHOIN    = .TRUE.
2578       ELOJET    = 5.0D0
2579
2580 * common /DTLCUT/
2581       ECMIN  = 3.5D0
2582       ECMAX  = 1.0D10
2583       XBJMIN = ZERO
2584       ELMIN = ZERO
2585       EGMIN = ZERO
2586       EGMAX = 1.0D10
2587       YMIN  = TINY10
2588       YMAX  = 0.999D0
2589       Q2MIN = TINY10
2590       Q2MAX = 10.0D0
2591       THMIN = ZERO
2592       THMAX = TWOPI
2593       Q2LI  = ZERO
2594       Q2HI  = 1.0D10
2595       ECMLI = ZERO
2596       ECMHI = 1.0D10
2597
2598 * common /DTVDMP/
2599       RL2       = 2.0D0
2600       INTRGE(1) = 1
2601       INTRGE(2) = 3
2602       IDPDF     = 2212
2603       MODEGA    = 4
2604       ISHAD(1)  = 1
2605       ISHAD(2)  = 1
2606       ISHAD(3)  = 1
2607       EPSPOL    = ZERO
2608
2609 * common /DTGLGP/
2610       JSTATB = 1000
2611       JBINSB = 49
2612       CGLB   = '        '
2613       IF (ITRSPT.EQ.1) THEN
2614          IOGLB  = 100
2615       ELSE
2616          IOGLB  = 0
2617       ENDIF
2618       LPROD  = .TRUE.
2619
2620 * common /DTHIS3/
2621       DO 16 I=1,50
2622          IHISPP(I) = 0
2623          IHISXS(I) = 0
2624    16 CONTINUE
2625       IXSTBL = 0
2626
2627 * common /DTVARE/
2628       VARELO = ZERO
2629       VAREHI = ZERO
2630       VARCLO = ZERO
2631       VARCHI = ZERO
2632
2633 * common /DTDIHA/
2634       DIBETA = -1.0D0
2635       DIALPH = ZERO
2636
2637 * common /LEPTOI/
2638       RPPN  = 0.0
2639       LEPIN = 0
2640       INTER = 0
2641
2642 * common /QNEUTO/
2643       NEUTYP = 1
2644       NEUDEC = 0
2645
2646 * common /DTEVNO/
2647       NEVENT = 1
2648       IF (ITRSPT.EQ.1) THEN
2649          ICASCA = 1
2650       ELSE
2651          ICASCA = 0
2652       ENDIF
2653
2654 * default Lab.-energy
2655       EPN = 200.0D0
2656       PPN = SQRT((EPN-AAM(IJPROJ))*(EPN+AAM(IJPROJ)))
2657
2658       RETURN
2659       END
2660 *
2661 *===aaevt==============================================================*
2662 *
2663 CDECK  ID>, DT_AAEVT
2664       SUBROUTINE DT_AAEVT(NEVTS,EPN,NPMASS,NPCHAR,NTMASS,NTCHAR,
2665      &                                             IDP,IGLAU)
2666
2667 ************************************************************************
2668 * This version dated 22.03.96 is written by S. Roesler.                *
2669 ************************************************************************
2670
2671       IMPLICIT DOUBLE PRECISION (A-H,O-Z)
2672       SAVE
2673
2674       PARAMETER ( LINP = 5 ,
2675      &            LOUT = 6 ,
2676      &            LDAT = 9 )
2677
2678       PARAMETER (NCOMPX=20,NEB=8,NQB= 5,KSITEB=50)
2679
2680 * emulsion treatment
2681       COMMON /DTCOMP/ EMUFRA(NCOMPX),IEMUMA(NCOMPX),IEMUCH(NCOMPX),
2682      &                NCOMPO,IEMUL
2683 * event flag
2684       COMMON /DTEVNO/ NEVENT,ICASCA
2685
2686       CHARACTER*8 DATE,HHMMSS
2687       DIMENSION IDMNYR(3)
2688
2689       KKMAT  = 1
2690       NMSG   = MAX(NEVTS/100,1)
2691
2692 * initialization of run-statistics and histograms
2693       CALL DT_STATIS(1)
2694
2695       CALL PHO_PHIST(1000,DUM)
2696
2697 * initialization of Glauber-formalism
2698       IF (NCOMPO.LE.0) THEN
2699          CALL DT_SHMAKI(NPMASS,NPCHAR,NTMASS,NTCHAR,IDP,EPN,IGLAU)
2700       ELSE
2701          DO 1 I=1,NCOMPO
2702             CALL DT_SHMAKI(NPMASS,NPCHAR,IEMUMA(I),IEMUCH(I),IDP,EPN,0)
2703     1    CONTINUE
2704       ENDIF
2705       CALL DT_SIGEMU
2706
2707       CALL IDATE(IDMNYR)
2708       WRITE(DATE,'(I2,''/'',I2,''/'',I2)')
2709      &   IDMNYR(1),IDMNYR(2),MOD(IDMNYR(3),100)
2710       CALL ITIME(IDMNYR)
2711       WRITE(HHMMSS,'(I2,'':'',I2,'':'',I2)')
2712      &   IDMNYR(1),IDMNYR(2),IDMNYR(3)
2713       WRITE(LOUT,1001) DATE,HHMMSS
2714  1001 FORMAT(/,' DT_AAEVT: Initialisation finished. ( Date: ',A8,
2715      &       '   Time: ',A8,' )')
2716
2717 * generate NEVTS events
2718       DO 2 IEVT=1,NEVTS
2719
2720 *  print run-status message
2721          IF (MOD(IEVT,NMSG).EQ.0) THEN
2722             CALL IDATE(IDMNYR)
2723             WRITE(DATE,'(I2,''/'',I2,''/'',I2)')
2724      &         IDMNYR(1),IDMNYR(2),MOD(IDMNYR(3),100)
2725             CALL ITIME(IDMNYR)
2726             WRITE(HHMMSS,'(I2,'':'',I2,'':'',I2)')
2727      &         IDMNYR(1),IDMNYR(2),IDMNYR(3)
2728             WRITE(LOUT,1000) IEVT-1,NEVTS,DATE,HHMMSS
2729  1000       FORMAT(/,1X,I8,' out of ',I8,' events sampled ( Date: ',A,
2730      &             '   Time: ',A,' )',/)
2731 C           WRITE(LOUT,1000) IEVT-1
2732 C1000       FORMAT(1X,I8,' events sampled')
2733          ENDIF
2734          NEVENT = IEVT
2735 *  treat nuclear emulsions
2736          IF (IEMUL.GT.0) CALL DT_GETEMU(NTMASS,NTCHAR,KKMAT,0)
2737 *  composite targets only
2738          KKMAT = -KKMAT
2739 *  sample this event
2740          CALL DT_KKINC(NPMASS,NPCHAR,NTMASS,NTCHAR,IDP,EPN,KKMAT,IREJ)
2741
2742          CALL PHO_PHIST(2000,DUM)
2743
2744     2 CONTINUE
2745
2746 * print run-statistics and histograms to output-unit 6
2747
2748       CALL PHO_PHIST(3000,DUM)
2749
2750       CALL DT_STATIS(2)
2751
2752       RETURN
2753       END
2754 *
2755 *===laevt==============================================================*
2756 *
2757 CDECK  ID>, DT_LAEVT
2758       SUBROUTINE DT_LAEVT(NEVTS,EPN,NPMASS,NPCHAR,NTMASS,NTCHAR,
2759      &                                             IDP,IGLAU)
2760
2761 ************************************************************************
2762 * Interface to run DPMJET for lepton-nucleus interactions.             *
2763 * Kinematics is sampled using the equivalent photon approximation      *
2764 * Based on GPHERA-routine by R. Engel.                                 *
2765 * This version dated 23.03.96 is written by S. Roesler.                *
2766 ************************************************************************
2767
2768       IMPLICIT DOUBLE PRECISION (A-H,O-Z)
2769       SAVE
2770
2771       PARAMETER ( LINP = 5 ,
2772      &            LOUT = 6 ,
2773      &            LDAT = 9 )
2774
2775       PARAMETER (TINY10=1.0D-10,TINY4=1.0D-4,
2776      &           ZERO=0.0D0,ONE=1.0D0,TWO=2.0D0,THREE=3.0D0)
2777       PARAMETER (TWOPI  = 6.283185307179586454D+00,
2778      &           PI     = TWOPI/TWO,
2779      &           ALPHEM = ONE/137.0D0)
2780
2781 C     CHARACTER*72 HEADER
2782
2783 * particle properties (BAMJET index convention)
2784       CHARACTER*8  ANAME
2785       COMMON /DTPART/ ANAME(210),AAM(210),GA(210),TAU(210),
2786      &                IICH(210),IIBAR(210),K1(210),K2(210)
2787 * event history
2788
2789       PARAMETER (NMXHKK=200000)
2790
2791       COMMON /DTEVT1/ NHKK,NEVHKK,ISTHKK(NMXHKK),IDHKK(NMXHKK),
2792      &                JMOHKK(2,NMXHKK),JDAHKK(2,NMXHKK),
2793      &                PHKK(5,NMXHKK),VHKK(4,NMXHKK),WHKK(4,NMXHKK)
2794 * extended event history
2795       COMMON /DTEVT2/ IDRES(NMXHKK),IDXRES(NMXHKK),NOBAM(NMXHKK),
2796      &                IDBAM(NMXHKK),IDCH(NMXHKK),NPOINT(10),
2797      &                IHIST(2,NMXHKK)
2798 * kinematical cuts for lepton-nucleus interactions
2799       COMMON /DTLCUT/ ECMIN,ECMAX,XBJMIN,ELMIN,EGMIN,EGMAX,YMIN,YMAX,
2800      &                Q2MIN,Q2MAX,THMIN,THMAX,Q2LI,Q2HI,ECMLI,ECMHI
2801 * properties of interacting particles
2802       COMMON /DTPRTA/ IT,ITZ,IP,IPZ,IJPROJ,IBPROJ,IJTARG,IBTARG
2803 * properties of photon/lepton projectiles
2804       COMMON /DTGPRO/ VIRT,PGAMM(4),PLEPT0(4),PLEPT1(4),PNUCL(4),IDIREC
2805 * kinematics at lepton-gamma vertex
2806       COMMON /DTLGVX/ PPL0(4),PPL1(4),PPG(4),PPA(4)
2807 * flags for activated histograms
2808       COMMON /DTHIS3/ IHISPP(50),IHISXS(50),IXSTBL
2809
2810       PARAMETER (NCOMPX=20,NEB=8,NQB= 5,KSITEB=50)
2811
2812 * emulsion treatment
2813       COMMON /DTCOMP/ EMUFRA(NCOMPX),IEMUMA(NCOMPX),IEMUCH(NCOMPX),
2814      &                NCOMPO,IEMUL
2815 * Glauber formalism: cross sections
2816       COMMON /DTGLXS/ ECMNN(NEB),Q2G(NQB),ECMNOW,Q2,
2817      &                XSTOT(NEB,NQB,NCOMPX),XSELA(NEB,NQB,NCOMPX),
2818      &                XSQEP(NEB,NQB,NCOMPX),XSQET(NEB,NQB,NCOMPX),
2819      &                XSQE2(NEB,NQB,NCOMPX),XSPRO(NEB,NQB,NCOMPX),
2820      &                XSDEL(NEB,NQB,NCOMPX),XSDQE(NEB,NQB,NCOMPX),
2821      &                XETOT(NEB,NQB,NCOMPX),XEELA(NEB,NQB,NCOMPX),
2822      &                XEQEP(NEB,NQB,NCOMPX),XEQET(NEB,NQB,NCOMPX),
2823      &                XEQE2(NEB,NQB,NCOMPX),XEPRO(NEB,NQB,NCOMPX),
2824      &                XEDEL(NEB,NQB,NCOMPX),XEDQE(NEB,NQB,NCOMPX),
2825      &                BSLOPE,NEBINI,NQBINI
2826 * nucleon-nucleon event-generator
2827       CHARACTER*8 CMODEL
2828       LOGICAL LPHOIN
2829       COMMON /DTMODL/ CMODEL(4),ELOJET,MCGENE,LPHOIN
2830 * flags for input different options
2831       LOGICAL LEMCCK,LHADRO,LSEADI,LEVAPO
2832       COMMON /DTFLG1/ IFRAG(2),IRESCO,IMSHL,IRESRJ,IOULEV(6),
2833      &                LEMCCK,LHADRO(0:9),LSEADI,LEVAPO,IFRAME,ITRSPT
2834 * event flag
2835       COMMON /DTEVNO/ NEVENT,ICASCA
2836
2837       DIMENSION XDUMB(40),BGTA(4)
2838
2839 * LEPTO
2840       IF (MCGENE.EQ.3) THEN
2841
2842          STOP ' This version does not contain LEPTO !'
2843
2844       ENDIF
2845
2846       KKMAT  = 1
2847       NMSG   = MAX(NEVTS/10,1)
2848
2849 * mass of incident lepton
2850       AMLPT  = AAM(IDP)
2851       AMLPT2 = AMLPT**2
2852       IDPPDG = IDT_IPDGHA(IDP)
2853
2854 * consistency of kinematical limits
2855       Q2MIN  = MAX(Q2MIN,TINY10)
2856       Q2MAX  = MAX(Q2MAX,TINY10)
2857       YMIN   = MIN(MAX(YMIN,TINY10),0.999D0)
2858       YMAX   = MIN(MAX(YMAX,TINY10),0.999D0)
2859
2860 * total energy of the lepton-nucleon system
2861       PTOTLN = SQRT( (PLEPT0(1)+PNUCL(1))**2+(PLEPT0(2)+PNUCL(2))**2
2862      &                                      +(PLEPT0(3)+PNUCL(3))**2 )
2863       ETOTLN = PLEPT0(4)+PNUCL(4)
2864       ECMLN  = SQRT((ETOTLN-PTOTLN)*(ETOTLN+PTOTLN))
2865       ECMAX  = MIN(ECMAX,ECMLN)
2866       WRITE(LOUT,1003) ECMIN,ECMAX,YMIN,YMAX,Q2MIN,Q2MAX,EGMIN,
2867      &                 THMIN,THMAX,ELMIN
2868  1003 FORMAT(1X,'LAEVT:',16X,'kinematical cuts',/,22X,
2869      &       '------------------',/,9X,'W (min)   =',
2870      &       F7.1,' GeV    (max) =',F7.1,' GeV',/,9X,'y (min)   =',
2871      &       F7.3,8X,'(max) =',F7.3,/,9X,'Q^2 (min) =',F7.1,
2872      &       ' GeV^2  (max) =',F7.1,' GeV^2',/,' (Lab)   E_g (min) ='
2873      &       ,F7.1,' GeV',/,' (Lab) theta (min) =',F7.4,8X,'(max) =',
2874      &       F7.4,'   for E_lpt >',F7.1,' GeV',/)
2875
2876 * Lorentz-parameter for transf. into Lab
2877       BGTA(1) = PNUCL(1)/AAM(1)
2878       BGTA(2) = PNUCL(2)/AAM(1)
2879       BGTA(3) = PNUCL(3)/AAM(1)
2880       BGTA(4) = PNUCL(4)/AAM(1)
2881 * LT of incident lepton into Lab and dump it in DTEVT1
2882       CALL DT_DALTRA(BGTA(4),-BGTA(1),-BGTA(2),-BGTA(3),
2883      &            PLEPT0(1),PLEPT0(2),PLEPT0(3),PLEPT0(4),
2884      &            PLTOT,PPL0(1),PPL0(2),PPL0(3),PPL0(4))
2885       CALL DT_DALTRA(BGTA(4),-BGTA(1),-BGTA(2),-BGTA(3),
2886      &            PNUCL(1),PNUCL(2),PNUCL(3),PNUCL(4),
2887      &            PLTOT,PPA(1),PPA(2),PPA(3),PPA(4))
2888 * maximum energy of photon nucleon system
2889       PTOTGN = SQRT((YMAX*PPL0(1)+PPA(1))**2+(YMAX*PPL0(2)+PPA(2))**2
2890      &                                      +(YMAX*PPL0(3)+PPA(3))**2)
2891       ETOTGN = YMAX*PPL0(4)+PPA(4)
2892       EGNMAX = SQRT((ETOTGN-PTOTGN)*(ETOTGN+PTOTGN))
2893       EGNMAX = MIN(EGNMAX,ECMAX)
2894 * minimum energy of photon nucleon system
2895       PTOTGN = SQRT((YMIN*PPL0(1)+PPA(1))**2+(YMIN*PPL0(2)+PPA(2))**2
2896      &                                      +(YMIN*PPL0(3)+PPA(3))**2)
2897       ETOTGN = YMIN*PPL0(4)+PPA(4)
2898       EGNMIN = SQRT((ETOTGN-PTOTGN)*(ETOTGN+PTOTGN))
2899       EGNMIN = MAX(EGNMIN,ECMIN)
2900
2901 * limits for Glauber-initialization
2902       Q2LI  = Q2MIN
2903       Q2HI  = MAX(Q2LI,MIN(Q2HI,Q2MAX))
2904       ECMLI = MAX(EGNMIN,THREE)
2905       ECMHI = EGNMAX
2906       WRITE(LOUT,1004) EGNMIN,EGNMAX,ECMLI,ECMHI,Q2LI,Q2HI
2907  1004 FORMAT(1X,'resulting limits:',/,9X,'W (min)   =',F7.1,
2908      &       ' GeV    (max) =',F7.1,' GeV',/,/,' limits for ',
2909      &       'Glauber-initialization:',/,9X,'W (min)   =',F7.1,
2910      &       ' GeV    (max) =',F7.1,' GeV',/,9X,'Q^2 (min) =',F7.1,
2911      &       ' GeV^2  (max) =',F7.1,' GeV^2',/)
2912 * initialization of Glauber-formalism
2913       IF (NCOMPO.LE.0) THEN
2914          CALL DT_SHMAKI(NPMASS,NPCHAR,NTMASS,NTCHAR,IDP,EPN,IGLAU)
2915       ELSE
2916          DO 9 I=1,NCOMPO
2917             CALL DT_SHMAKI(NPMASS,NPCHAR,IEMUMA(I),IEMUCH(I),IDP,EPN,0)
2918     9    CONTINUE
2919       ENDIF
2920       CALL DT_SIGEMU
2921
2922 * initialization of run-statistics and histograms
2923       CALL DT_STATIS(1)
2924
2925       CALL PHO_PHIST(1000,DUM)
2926
2927 * maximum photon-nucleus cross section
2928       I1  = 1
2929       I2  = 1
2930       RAT = ONE
2931       IF (EGNMAX.GE.ECMNN(NEBINI)) THEN
2932          I1  = NEBINI
2933          I2  = NEBINI
2934      &nbs