STARLIGHT/starlight/dpmjet/dpmjetint.f

   1 *
   2 *===program crint======================================================*
   3 *
   4 C      OPTIONS/ EXTEND_SOURCE
   5 C      SUBROUTINE CRINT
   6       SUBROUTINE DT_PRODUCEEVENT(ENERGY_SL, NPARTICLES)
   7
   8
   9       IMPLICIT DOUBLE PRECISION (A-H,O-Z)
  10       REAL ENERGY_SL
  11       INTEGER INIT
  12       REAL ne,etest,prob,slump
  13       SAVE
  14
  15 * Call the init sub routine in the first event
  16       DATA INIT /0/
  17
  18       PARAMETER (NMXHKK=200000)
  19
  20       COMMON /DTIONT/ LINP,LOUT,LDAT
  21
  22       COMMON /DTEVT1/ NHKK,NEVHKK,ISTHKK(NMXHKK),IDHKK(NMXHKK),
  23      &                JMOHKK(2,NMXHKK),JDAHKK(2,NMXHKK),
  24      &                PHKK(5,NMXHKK),VHKK(4,NMXHKK),WHKK(4,NMXHKK)
  25
  26 *     event flag
  27       COMMON /DTEVNO/ NEVENT, ICASCA
  28
  29       IF(INIT.EQ.0) THEN
  30          OPEN (UNIT = 50, file = "my.input")
  31          LINP = 50
  32          CALL DT_DTUINI(NEVTS,EPN,NPMASS,NPCHAR,NTMASS,NTCHAR,IDP,IEMU)
  33 *        Init called, make sure it's not called again
  34          INIT = 1
  35       ENDIF
  36 *-----------------------------------------------------------------------
  37 *     generation of one event
  38       NEVENT = 1
  39       KKMAT = -1
  40
  41 *   If an energy-range has been defined with the ENERGY input-card the
  42 *   laboratory energy ELAB can be set to any value within that range,..
  43 C        ELAB = DT_RNDM(EPN)*(EPN-0.5D7)+0.5D7
  44
  45 *   ..otherwise it has to coincide with EPN.
  46 C        ELAB = EPN
  47
  48       ELAB = ENERGY_SL
  49
  50 *   sampling of one event
  51
  52 *     TEST
  53
  54       CALL DT_KKINC(NPMASS,NPCHAR,NTMASS,NTCHAR,IDP,ELAB,KKMAT,IREJ)
  55
  56       IF (IREJ.NE.0) RETURN
  57
  58 c     Return the number of particles produced
  59
  60 c     Fill the particle info
  61       CALL DT_GETPARTICLES(NPARTICLES)
  62
  63       END
  64
  65
  66       SUBROUTINE DT_GETPARTICLES(NPARTICLES)
  67
  68       IMPLICIT DOUBLE PRECISION (A-H,O-Z)
  69       INTEGER pid,qch,q_sum,Ntpc,Nfinal,NACCEPT,IPART,RES
  70       DOUBLE PRECISION yrap,pt,mass,mt,etot
  71       DOUBLE PRECISION pt_cut_tpc
  72       PARAMETER(pt_cut_tpc=0.050)
  73
  74       SAVE
  75 *
  76 * COMMON /DTEVT1/ :
  77 *                   NHKK         number of entries in common block
  78 *                   NEVHKK       number of the event
  79 *                   ISTHKK(i)    status code for entry i
  80 *                   IDHKK(i)     identifier for the entry
  81 *                                (for particles: identifier according
  82 *                                 to the PDG numbering scheme)
  83 *                   JMOHKK(1,i)  pointer to the entry of the first mother
  84 *                                of entry i
  85 *                   JMOHKK(2,i)  pointer to the entry of the second mother
  86 *                                of entry i
  87 *                   JDAHKK(1,i)  pointer to the entry of the first daughter
  88 *                                of entry i
  89 *                   JDAHKK(2,i)  pointer to the entry of the second daughter
  90 *                                of entry i
  91 *                   PHKK(1..3,i) 3-momentum
  92 *                   PHKK(4,i)    energy
  93 *                   PHKK(5,i)    mass
  94 *
  95 * event history
  96
  97       PARAMETER (NMXHKK=200000)
  98
  99       COMMON /DTEVT1/ NHKK,NEVHKK,ISTHKK(NMXHKK),IDHKK(NMXHKK),
 100      &                JMOHKK(2,NMXHKK),JDAHKK(2,NMXHKK),
 101      &                PHKK(5,NMXHKK),VHKK(4,NMXHKK),WHKK(4,NMXHKK)
 102
 103 * extended event history
 104       COMMON /DTEVT2/ IDRES(NMXHKK),IDXRES(NMXHKK),NOBAM(NMXHKK),
 105      &                IDBAM(NMXHKK),IDCH(NMXHKK),NPOINT(10),
 106      &                IHIST(2,NMXHKK)
 107
 108       DOUBLE PRECISION SLPX, SLPY, SLPZ, SLE, SLM
 109       INTEGER SLPID, SLCHARGE
 110       COMMON /DPMJETPARTICLE/ SLPX(NMXHKK), SLPY(NMXHKK), SLPZ(NMXHKK),
 111      &       SLE(NMXHKK), SLM(NMXHKK), SLPID(NMXHKK), SLCHARGE(NMXHKK)
 112
 113
 114 C     >> Set Counter to Zero
 115
 116       Nfinal=0
 117
 118       DO 42 I=1, NHKK
 119 c      I = IPART
 120
 121 CC       >> Remove all non-final-state particles
 122         IF(.not.(ISTHKK(I).eq.1.or.ISTHKK(I).eq.-1.or.
 123      $ISTHKK(I).eq.1001)) GOTO 42
 124
 125 C       >> Find Particle Charge, qch
 126         IF((ABS(ISTHKK(I)).eq.1).and.(IDHKK(I).ne.80000))THEN
 127 C         >> final state ptcles except nuclei
 128
 129           qch=IPHO_CHR3(IDHKK(I),1)/3
 130         ELSEIF(IDHKK(I).eq.80000)THEN
 131 C         >> final state nuclei
 132           qch=IDXRES(I)
 133         ELSE
 134 C         >> not a final state particle, qch not interesting
 135           qch=-999
 136         ENDIF
 137
 138         Nfinal = Nfinal + 1
 139         SLPX(Nfinal) = PHKK(1,I)
 140         SLPY(Nfinal) = PHKK(2,I)
 141         SLPZ(Nfinal) = PHKK(3,I)
 142         SLE(Nfinal) = PHKK(4,I)
 143         SLM(Nfinal) = PHKK(5,I)
 144         SLPID(Nfinal) = IDHKK(I)
 145         SLCHARGE(Nfinal) = qch
 146
 147  42     CONTINUE
 148         NPARTICLES = Nfinal
 149
 150       END
 151
 152       SUBROUTINE DT_USRHIS(MODE)
 153 c Dummy to make the linker happy
 154       END
 155