DPMJET/user3.0-5F.f

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