STARLIGHT code and interface

[u/mrichter/AliRoot.git] / STARLIGHT / starlight / dpmjet / .svn / text-base / dpmjetint.f.svn-base

*\r
*===program crint======================================================*\r
*\r
C      OPTIONS/ EXTEND_SOURCE\r
C      SUBROUTINE CRINT\r
      SUBROUTINE DT_PRODUCEEVENT(ENERGY_SL, NPARTICLES)\r
\r
\r
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)\r
      REAL ENERGY_SL\r
      INTEGER INIT\r
      REAL ne,etest,prob,slump\r
      SAVE\r
\r
* Call the init sub routine in the first event\r
      DATA INIT /0/\r
\r
      PARAMETER (NMXHKK=200000)\r
\r
      COMMON /DTIONT/ LINP,LOUT,LDAT\r
\r
      COMMON /DTEVT1/ NHKK,NEVHKK,ISTHKK(NMXHKK),IDHKK(NMXHKK),\r
     &                JMOHKK(2,NMXHKK),JDAHKK(2,NMXHKK),\r
     &                PHKK(5,NMXHKK),VHKK(4,NMXHKK),WHKK(4,NMXHKK)\r
\r
*     event flag\r
      COMMON /DTEVNO/ NEVENT, ICASCA\r
\r
      IF(INIT.EQ.0) THEN\r
         OPEN (UNIT = 50, file = "my.input")    \r
	 LINP = 50\r
         CALL DT_DTUINI(NEVTS,EPN,NPMASS,NPCHAR,NTMASS,NTCHAR,IDP,IEMU)\r
*        Init called, make sure it's not called again\r
         INIT = 1\r
      ENDIF\r
*-----------------------------------------------------------------------\r
*     generation of one event\r
      NEVENT = 1\r
      KKMAT = -1\r
\r
*   If an energy-range has been defined with the ENERGY input-card the\r
*   laboratory energy ELAB can be set to any value within that range,..\r
C        ELAB = DT_RNDM(EPN)*(EPN-0.5D7)+0.5D7\r
\r
*   ..otherwise it has to coincide with EPN.\r
C        ELAB = EPN\r
\r
      ELAB = ENERGY_SL\r
\r
*   sampling of one event\r
\r
*     TEST\r
\r
      CALL DT_KKINC(NPMASS,NPCHAR,NTMASS,NTCHAR,IDP,ELAB,KKMAT,IREJ)\r
\r
      IF (IREJ.NE.0) RETURN\r
\r
c     Return the number of particles produced\r
      \r
c     Fill the particle info \r
      CALL DT_GETPARTICLES(NPARTICLES)\r
\r
      END\r
\r
\r
      SUBROUTINE DT_GETPARTICLES(NPARTICLES)\r
\r
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)\r
      INTEGER pid,qch,q_sum,Ntpc,Nfinal,NACCEPT,IPART,RES\r
      DOUBLE PRECISION yrap,pt,mass,mt,etot\r
      DOUBLE PRECISION pt_cut_tpc\r
      PARAMETER(pt_cut_tpc=0.050)\r
\r
      SAVE\r
*\r
* COMMON /DTEVT1/ :\r
*                   NHKK         number of entries in common block\r
*                   NEVHKK       number of the event\r
*                   ISTHKK(i)    status code for entry i\r
*                   IDHKK(i)     identifier for the entry\r
*                                (for particles: identifier according\r
*                                 to the PDG numbering scheme)\r
*                   JMOHKK(1,i)  pointer to the entry of the first mother\r
*                                of entry i\r
*                   JMOHKK(2,i)  pointer to the entry of the second mother\r
*                                of entry i\r
*                   JDAHKK(1,i)  pointer to the entry of the first daughter\r
*                                of entry i\r
*                   JDAHKK(2,i)  pointer to the entry of the second daughter\r
*                                of entry i\r
*                   PHKK(1..3,i) 3-momentum\r
*                   PHKK(4,i)    energy\r
*                   PHKK(5,i)    mass\r
*\r
* event history\r
\r
      PARAMETER (NMXHKK=200000)\r
\r
      COMMON /DTEVT1/ NHKK,NEVHKK,ISTHKK(NMXHKK),IDHKK(NMXHKK),\r
     &                JMOHKK(2,NMXHKK),JDAHKK(2,NMXHKK),\r
     &                PHKK(5,NMXHKK),VHKK(4,NMXHKK),WHKK(4,NMXHKK)\r
\r
* extended event history\r
      COMMON /DTEVT2/ IDRES(NMXHKK),IDXRES(NMXHKK),NOBAM(NMXHKK),\r
     &                IDBAM(NMXHKK),IDCH(NMXHKK),NPOINT(10),\r
     &                IHIST(2,NMXHKK)\r
\r
      DOUBLE PRECISION SLPX, SLPY, SLPZ, SLE, SLM\r
      INTEGER SLPID, SLCHARGE\r
      COMMON /DPMJETPARTICLE/ SLPX(NMXHKK), SLPY(NMXHKK), SLPZ(NMXHKK),\r
     &       SLE(NMXHKK), SLM(NMXHKK), SLPID(NMXHKK), SLCHARGE(NMXHKK)\r
\r
\r
C     >> Set Counter to Zero\r
\r
      Nfinal=0\r
      \r
      DO 42 I=1, NHKK\r
c      I = IPART\r
\r
CC       >> Remove all non-final-state particles\r
        IF(.not.(ISTHKK(I).eq.1.or.ISTHKK(I).eq.-1.or.\r
     $ISTHKK(I).eq.1001)) GOTO 42\r
\r
C	>> Find Particle Charge, qch\r
        IF((ABS(ISTHKK(I)).eq.1).and.(IDHKK(I).ne.80000))THEN\r
C         >> final state ptcles except nuclei\r
\r
          qch=IPHO_CHR3(IDHKK(I),1)/3\r
        ELSEIF(IDHKK(I).eq.80000)THEN\r
C         >> final state nuclei\r
          qch=IDXRES(I)\r
        ELSE\r
C         >> not a final state particle, qch not interesting\r
          qch=-999\r
        ENDIF\r
\r
	Nfinal = Nfinal + 1\r
	SLPX(Nfinal) = PHKK(1,I)\r
        SLPY(Nfinal) = PHKK(2,I)\r
        SLPZ(Nfinal) = PHKK(3,I)\r
        SLE(Nfinal) = PHKK(4,I)\r
        SLM(Nfinal) = PHKK(5,I)\r
        SLPID(Nfinal) = IDHKK(I)\r
        SLCHARGE(Nfinal) = qch\r
\r
 42     CONTINUE\r
        NPARTICLES = Nfinal\r
  \r
      END\r
\r
      SUBROUTINE DT_USRHIS(MODE)\r
c Dummy to make the linker happy\r
      END\r
\r