Changes to make code work for pp collisions.
[u/mrichter/AliRoot.git] / DPMJET / user3.0-4.f
1 *
2 *===program dpmjet=====================================================*
3 *
4       PROGRAM DPMJET
5
6       IMPLICIT DOUBLE PRECISION (A-H,O-Z)
7       SAVE
8
9       EXTERNAL BDANDI, BDEMFL, BDEVAP, BDMULS, BDNOPT, BDPHFL, BDPWXS,
10      &         BDRNDM, BDSQZI, BLKDT1, BLKDT2, BLKDT3, BLKDT4, BLKDT5,
11      &         BLKDT6, BLKDT7, COUNTR, ELPROI, HADINI, BDPREE
12
13 C     EXTERNAL PYDATA
14
15 * event flag
16       COMMON /DTEVNO/ NEVENT,ICASCA
17
18 *-----------------------------------------------------------------------
19 * initialization
20
21 *   the following statement provides a call to DT_USRHIS(MODE=1) for
22 *   histogram initialization etc.
23       CALL DT_DTUINI(NEVTS,EPN,NPMASS,NPCHAR,NTMASS,NTCHAR,IDP,IEMU)
24
25 *-----------------------------------------------------------------------
26 * generation of events
27
28       DO 1 IEVT=1,NEVTS
29
30 *   some defaults, do not change!
31          NEVENT = IEVT
32          KKMAT  = -1
33          ELAB   = EPN
34 *   uncomment if dpmjet3 is linked to particle transport code
35 C        ICASCA = 1
36
37 ************************************************************************
38 * The following lines show how to select the target nucleus for runs
39 * with composite targets (and fixed projectile and energy!).
40 *
41 *   Sampling of the target nucleus (mass number NTMASS, charge NTCHAR)
42 *   according to the fractions defined with EMULSION input-cards.
43 *   The different nuclei are numbered as KKMAT = 1,2,3,...  according to
44 *   their appearance in the input-file.
45          IF (IEMU.GT.0) THEN
46 *   Replace this selection by your own one if needed.
47             CALL DT_GETEMU(NTMASS,NTCHAR,KKMAT,0)
48 *   Kkmat has to be negative for composite targets!
49             KKMAT = -KKMAT
50          ENDIF
51 ************************************************************************
52
53 ************************************************************************
54 * The following lines show how to define projectile, target and energy
55 * for this event in runs with Glauber-data file pre-initialized for a
56 * certain range of projectiles, targets and energies. The definitions
57 * have to be within the pre-initialized parameter range.
58 *
59 *   projectile-id (for hadron projectiles)
60 C        IDP    = 1
61 *   projectile mass and charge numbers
62 C        NPMASS = 12
63 C        NPCHAR = 6
64 *   target mass and charge numbers
65 C        NTMASS = 16
66 C        NTCHAR = 8
67 *   lab energy
68 C        ELAB = 200.0D0
69 ************************************************************************
70
71 ************************************************************************
72 * If an energy-range has been defined with the ENERGY input-card the
73 * laboratory energy ELAB can be set to any value within that range. For
74 * example:
75 C        ELO  = 10.0D0
76 C        EHI  = 1000.0D0
77 C        ELAB = DT_RNDM(ELAB)*(EHI-ELO)+ELO
78 ************************************************************************
79
80 *   sampling of one event
81          CALL DT_KKINC(NPMASS,NPCHAR,NTMASS,NTCHAR,IDP,ELAB,KKMAT,IREJ)
82          IF (IREJ.NE.0) GOTO 1
83
84 *   the following statement provides a call to DT_USRHIS(MODE=2) from
85 *   where the final state particles can be obtained
86
87          CALL PHO_PHIST(2000,DUM)
88
89     1 CONTINUE
90
91 *-----------------------------------------------------------------------
92 * output, statistics etc.
93
94 *   the following statement provides a call to DT_USRHIS(MODE=3) in
95 *   order to calculate histograms etc.
96       CALL DT_DTUOUT
97
98       END
99 *
100 *===usrhis=============================================================*
101 *
102 CDECK  ID>, DT_USRHIS
103       SUBROUTINE DT_USRHIS(MODE)
104
105       IMPLICIT DOUBLE PRECISION (A-H,O-Z)
106       SAVE
107 *
108 * COMMON /DTEVT1/ :
109 *                   NHKK         number of entries in common block
110 *                   NEVHKK       number of the event
111 *                   ISTHKK(i)    status code for entry i
112 *                   IDHKK(i)     identifier for the entry
113 *                                (for particles: identifier according
114 *                                 to the PDG numbering scheme)
115 *                   JMOHKK(1,i)  pointer to the entry of the first mother
116 *                                of entry i
117 *                   JMOHKK(2,i)  pointer to the entry of the second mother
118 *                                of entry i
119 *                   JDAHKK(1,i)  pointer to the entry of the first daughter
120 *                                of entry i
121 *                   JDAHKK(2,i)  pointer to the entry of the second daughter
122 *                                of entry i
123 *                   PHKK(1..3,i) 3-momentum
124 *                   PHKK(4,i)    energy
125 *                   PHKK(5,i)    mass
126 *
127 * event history
128
129       PARAMETER (NMXHKK=200000)
130
131       COMMON /DTEVT1/ NHKK,NEVHKK,ISTHKK(NMXHKK),IDHKK(NMXHKK),
132      &                JMOHKK(2,NMXHKK),JDAHKK(2,NMXHKK),
133      &                PHKK(5,NMXHKK),VHKK(4,NMXHKK),WHKK(4,NMXHKK)
134 * extended event history
135       COMMON /DTEVT2/ IDRES(NMXHKK),IDXRES(NMXHKK),NOBAM(NMXHKK),
136      &                IDBAM(NMXHKK),IDCH(NMXHKK),NPOINT(10),
137      &                IHIST(2,NMXHKK)
138
139       GOTO (1,2,3) MODE
140
141 *------------------------------------------------------------------
142 *
143     1 CONTINUE
144 *
145 * initializations
146 *
147 *  Called with MODE=1 once at the beginning of the run.
148 *
149       RETURN
150 *
151 *------------------------------------------------------------------
152 *
153     2 CONTINUE
154 *
155 * scoring of the present event
156 *
157 *  Called with MODE=2 every time one event has been finished.
158 *
159 *  The final state particles from the actual event (number NEVHKK)
160 *  can be found in DTEVT1 and identified by their status:
161 *
162 *     ISTHKK(i) = 1    final state particle produced in
163 *                      photon-/hadron-/nucleon-nucleon collisions or
164 *                      in intranuclear cascade processes
165 *                -1    nucleons, deuterons, H-3, He-3, He-4 evaporated
166 *                      from excited nucleus and
167 *                      photons produced in nuclear deexcitation processes
168 *                1001  residual nucleus (ground state)
169 *
170 *  The types of these particles/nuclei are given in IDHKK as follows
171 *
172 *     all final state part. except nuclei :
173 *       IDHKK(i)=particle identifier according to PDG numbering scheme
174 *     nuclei (evaporation products, and residual nucleus) :
175 *       IDHKK(i)=80000, IDRES(i)=mass number, IDXRES(i)=charge number
176 *
177 *  The 4-momenta and masses can be found in PHKK (target nucleus rest frame):
178 *                   PHKK(1..3,i) 3-momentum (p_x,p_y,p_z)
179 *                   PHKK(4,i)    energy
180 *                   PHKK(5,i)    mass
181 *
182 *
183 *
184 *  Pick out the final state particles from DTEVT1 in each event for
185 *  instance by the following loop (NHKK=number of entries in the present
186 *  event) and fill your histograms
187       NPART = 0
188       DO 19 I=1,NHKK
189 C         IF (ABS(ISTHKK(I)).NE. -1 .AND. ABS(ISTHKK(I)) .NE. 1001 ) THEN
190          IF (ABS(ISTHKK(I)).EQ. 1) THEN
191             NPART=NPART+1
192          ENDIF
193  19   CONTINUE  
194       write(16,*) npart
195       DO 20 I=1,NHKK
196 C         IF (ABS(ISTHKK(I)).NE. -1 .AND. ABS(ISTHKK(I)) .NE. 1001 ) THEN
197          IF (ABS(ISTHKK(I)).EQ. 1) THEN
198             WRITE(16,100) 
199      &           i, isthkk(i), idhkk(i),
200      &           phkk(1,i), phkk(2,i), phkk(3,i), phkk(4,i), phkk(5,i),
201      &           jmohkk(1,i), jmohkk(2,i), jdahkk(1,i), jdahkk(2,i)
202          ELSEIF (ABS(ISTHKK(I)).EQ.1001) THEN
203          ENDIF
204  20   CONTINUE
205  100  FORMAT(3i6, 5e13.5, 4i7)
206 *  At any time during the run a list of the actual entries in DTEVT1 and
207 *  DTEVT2 can be obtained (output unit 6) by the following statement:
208 C     CALL DT_EVTOUT(4)
209
210       RETURN
211 *
212 *------------------------------------------------------------------
213 *
214     3 CONTINUE
215 *
216 * output/statistics/histograms etc.
217 *
218 *  Called with MODE=3 once after all events have been sampled.
219 *
220       RETURN
221
222       END